Содержание к диссертации
Введение
ГЛАВА 1. Изученность вопроса 7
ГЛАВА 2. Построение "почти" обратного оператора линеаризованной обратной задачи сейсмки в скалярном случае 12
2.1. Постановка обратной задачи и её линеаризация 12
2.2. Асимптотическое решение обратной задачи 15
2.2.1. Ключевая интегральная формула 16
2.2.2. Формула обращения, и её свойства 18
2.3. Численные примеры построения изображений для простой модели 25
ГЛАВА 3. Численное построение изображений на примере синтетических данных SIGSBEE2A 30
3.1. Описание модели Sigsbee2a 30
3.2. Построение изображений в простой части модели Sigsbee2a... 33
3.3. Построение изображений в сложной части модели Sigsbee2a...41
3.3.1. Трассирование гауссовых пучков в сглаженной модели...42
3.3.2. Построение изображений в сглаженной модели 47
3.4. Численное исследование устойчивости метода 53
ГЛАВА 4. Построение изображений в анизотропных средах 63
2.1. Гауссовы пучки в анизотропной среде 63
2.2. Построение изображений рассеивающих объектов 68
2.3. Численные примеры построения изображений рассеивающих объектов в трансверсально-изотропной среде 72
Заключение 74
Литература
- Асимптотическое решение обратной задачи
- Формула обращения, и её свойства
- Построение изображений в простой части модели Sigsbee2a...
- Построение изображений рассеивающих объектов
Введение к работе
Объектом исследований в данной работе являются волновые сейсмические поля многократного перекрытия.
Актуальность. Одним из основных этапов обработки сейсмических данных является построение волновых изображений. Входной информацией при их реализации, кроме самих данных, является наличие макроскоростной модели. Ниже под этим выражением, следуя устоявшейся к настоящему времени терминологии, будет пониматься такое распределение скорости, которое предписывает правильные времена распространения волн, но не приводит к возникновению сколько-нибудь значимых отражений. Естественно, это верно для достаточно больших расстояний - порядка нескольких длин волн. Построение изображений, таким образом, сводится к отображению в заданной скоростной модели объектов, вызывающих изменение направления распространения сейсмической энергии и ее возвращение на свободную поверхность, к расположенной там системе наблюдений.
Используемые в настоящее время методы построения изображений в заданной макроскоростной модели весьма разнообразны и имеют довольно давнюю историю. Наиболее исследованным среди них является сейсмическая миграция, восходящая еще к предложенному Ю.В. Ризниченко [28] методу полей - времен. К настоящему времени это направление является одним из наиболее развиваемых в современной геофизике. Здесь необходимо подчеркнуть, что подавляющее большинство используемых в настоящее время миграционных процедур в качестве "скелета" используют прямое и обратное "кинематическое отображение" элементарных рассеивающих объектов. Такое отображение переводит элементарный рассеивающий объект (рассеивающая "точка") из пространственной области в поверхность Гюйгенса во временной области (годограф дифрагированной волны ), а точку из временной области - в поверхность равных времен пробега дифрагированных волн или изохрону [58],[77]. Подавляющее количество используемых в настоящее время алгоритмов миграции сейсмических данных базируется на построении изображений путем суммирования исходных данных вдоль годографа дифрагированной волны. Естественно, что при использовании корректной макроскоростной модели будет получено "кинематически" корректное изображение изучаемой среды. Однако его динамическая достоверность при этом никоим образом не гарантируется, так как не принимается в учет изменчивость амплитуды волны в процессе ее распространения. Поэтому, при прочих равных условиях, интенсивность, или "яркость", более глубоко залегающих границ на получаемом волновом изображении, будет, как правило, меньше, чем у вышележащих. В то же время, получение изображений, яркость которых отражает распределение физических свойств горных пород, слагающих разрез, является чрезвычайно важной проблемой современной сейсморазведки и необходимым условием для проведения правильной интерпретации результатов обработки сейсмических данных. В настоящее время для получения таких изображений осуществляют вышеупомянутое отображение с некоторым весом, предназначенным компенсировать изменчивость амплитуды волн при их распространении в среде.
Приемы, используемые для построения таких весов берут свое начало от работ [39], [63], [41] и основаны на применении различных видов асимптотического разложения функции Грина в неоднородной среде. Наиболее употребительным при этом является использование нулевого члена лучевого разложения. Однако, его область применения ограничена только теми средами, для которых центральное поле лучей оказывается регулярным. Предлагаемый в настоящей работе подход к построению изображений в истинных амплитудах свободен от такого недостатка, так как основан на асимптотическом представлении функции Грина в виде суперпозиции Гауссовых пучков, являющихся регулярными во всем пространстве. Замечательным фактом при этом оказалась и возможность получения с помощью этих же Гауссовых пучков так называемых "селективных" изображений. Селективные изображения - это такие изображения среды, на которых представлены только лишь определенным образом расположенные (наклоненные) в пространстве отражающие объекты. Они несут в себе чрезвычайно полезную информацию при изучении распределения в среде различного рода сингулярных объектов, таких как разломы, зоны повышенной трещшговатости и другие. Используемые в настоящее время приемы для построения изображений таких зон никоим образом не позволяют получать "истинные" амплитуды и, следовательно, получать количественные оценки концентрации и степени неоднородности этих объектов. Предлагаемый в данной работе подход позволяет это сделать и, тем самым, существенно повысить информативность результатов обработки.
Цель исследований - обеспечить повышение информативности результатов обработки сейсмических данных за счет получения изображений (включая и селективные изображения) отражающих/рассеивающих объектов в "истинных" амплитудах для макроско-ростных моделей сколь угодно сложного строения.
Научная задача - разработать численные методы решения линеаризованной обратной задачи сейсмики и на этой основе создать алгоритмы обработки сейсмических данных многократного перекрытия, обеспечивающие достоверное восстановление расположения и физических свойств отражающих/рассеивающих объектов, помещенных в макроскорост-ную модель сколь угодно сложного строения.
Основные этапы исследований: получить основные интегральные тождества, связывающие точные решения уравнения Гельмгольца и уравнения эластодинамики с их асимптотическими решениями (Гауссовы пучки); с помощью конечно-разностного моделирования изучить основные особенности процессов распространения Гауссовых пучков в сложноустроенных средах на примере модели SIGSBEE2A и установить на этой основе границы применимости полученных интегральных тождеств; на основе полученных интегральных тождеств и теории псевдодифференциальных операторов построить "почти" обратный оператор для отыскания приближенного решения линеаризованной обратной задачи сейсмики; провести серию численных экспериментов для верификации полученных соотношений.
Методы исследований и фактический материал.
Теоретической основой решения поставленной научной задачи являются: результаты современной теории распространения волновых полей в сложнопостро-енных средах, особенно относящиеся к построению приближенных, асимптотических решений, обладающих свойством локальной сосредоточенности; теория псевдодифференциальных операторов, позволившая получить представление действия сконструированного "почти" обратного оператора по степеням гладкости; конечно-разностные методы, позволившие верифицировать полученные результаты.
Численные алгоритмы создавались и тестировались в процессе обработки представительной серии синтетических данных многократного перекрытия, среди которых ключевую роль занимал международно признанный набор данных для модели SIGSBEE2A.
Верификация результатов осуществлялась на наборе данных SIGSBEE2A, представляющем из себя данные многократного перекрытия для тонкослоистой среды с серией малоамплитудных разломов, содержащей в себе массивное соляное тело с границей чрезвычайно сложной формы.
Защищаемые научные положения:
Построен "почти" обратный оператор линеаризованной обратной задачи сейсми-ки, восстанавливающий амплитуды "разрывов" контраста скорости среды.
Построены селективные изображения рассеивающих объектов и отражающих площадок, на основе "почти" обратного оператора.
Новизна работы. Личный вклад.
Выполненные автором лично преобразования тождества, полученного с использова нием теоремы Грина, и анализ полученных соотношений, основывающиеся на исполь зовании теории псевдодифференциальных операторов, позволили построить "почти" обратный оператор для линеаризованной задачи сейсмики. - Проведенный автором лично анализ данного оператора, позволил построить селективные изображений среды, на которых преимущественно представлены только лишь определенным образом ориентированные отражающие и рассеивающие элементы. — Соискателем лично разработано программное обеспечение, реализующее действие данного "почти" обратного оператора.
Соискателем лично получены тождества, посредством использования теоремы Гри на, на основе которых построены изображения рассеивающих объектов в анизотроп ных средах.
Теоретическая и практическая значимость результатов.
Впервые предложен и обоснован метод построения изображений отражающих и рассеивающих объектов, сохраняющий их амплитудные характеристики (коэффициенты отражения и интенсивность рассеяния), без каких-либо ограничений на строение макроско-ростной модели. Проведенная серия численных экспериментов на международно признанном тестовом наборе синтетических данных SIGSBEE2A позволила верифицировать предложенный метод и уточнить его границы применимости. Разработанное программное обеспечение может служить прототипом для создания промышленных обрабатывающих систем, обеспечивающих построение изображений отражающих/рассеивающих объектов в истинных амплитудах.
Апробация работы и публикации.
Основные положения и результаты докладывались на Международной научной студенческой конференции "МНСК - 2004" (Новосибирск, 2004) [24], Международной конференции по вычислительной математике "МКВМ 2004" (Новосибирск, 2004), Международной конференции "Сейсмические волны в латералыю-неоднородных средах" (Прага(Чехия), 2005) [71], международной конференции "Дни дифракции 2005" (Санкт-Петербург, 2005) [70], а также на семинарах Института геофизики СО РАН, Института вычислительной математики и математической геофизики СО РАН и Института математики СО РАН.
Результаты исследований по теме диссертации изложены в 8 опубликованных работах. Из них 4 работы - это тезисы докладов [24], [25], [26], [70] на российских международных конференциях. Одна работа - это тезисы доклада [71] на Международной конференции, проводимой раз в четыре года Институтом геофизики Чешской АН. Одна работа - это расширенные тезисы на Международной конференции, проводимой EAGE в Мадриде в июне 2005 r[G9]. Ещё одна работа - это статья в рецензируемых материалах международной конференции "Дни дифракции 2005" [72]. И последняя работа - это статья в журнале Докл. РАН [27].
Диссертация выполнена в лаборатории динамических проблем сейсмики Института геофизики СО РАН.
Автор выражает искреннюю признательность научному руководителю, к.ф.-м.н. В.А. Чеверде за всестороннюю поддержку и постоянное внимание. Хочется также поблагодарить к.ф.-м.н. A.M. Айзенберга и к.ф.-м.н. А.А. Дучкова за то, что они взялись за труд ознакомиться с работой и высказать о ней свое мнение.
Асимптотическое решение обратной задачи
Излагаемый ниже метод построения изображений в истинных амплитудах основывается на работе [33]. В ней было предложено использовать разложение волнового поля по гауссовым пучкам для картирования слабоконтрастиых рассеивающих объектов, расположенных вблизи от резких границ раздела. Для этого с помощью взвешенного суммирования данных многократного перекрытия формируется гауссов пучок с заранее заданными параметрами. После осуществления такого суммирования исходные данные трансформируются в следы на системе наблюдения от соответствующего Гауссова пучка, претерпевшего отражения и рассеяния на локальных объектах внутри среды. Последующая миграция набора таких просуммированных данных вдоль подходящих направлений и позволяет получать изображения рассеивающих объектов за счёт существенного подавления интенсивности регулярных отражений. Развитие этого метода осуществляется в два этапа: построение селективных изображений в истинных амплитудах; суммирование селективных изображений и получение полного изображения в истинных амплитудах.
Прежде чем начать изложение подхода, получим интегральную формулу, которая нам понадобится в дальнейшем. Давайте рассмотрим краевую задачу и cl(x,z) dw(x,z;u ) ,2 Aw + л w = 0; (2.14) dz = f(x;u); (2.15) г=0 -Ьусловие на бесконечности (см. [50],[4]). (2.16) Применение теоремы Грина (см. [5]) к функции w(x,z;u ) и функции Грина G(x,z; T];u ): [ . , » dG(x,z;,ri;u ) . . 8W(X,Z;UJ) J [w{x,z;u) Q Z G(x,z;,r);u) — }dx = Я ыг „ ,. J JAG + 4M0) -w-(Aw + 4MW) G)]dxiz приводит к следующему соотношению: і f(x; u)G(x, z\ , rj; u)dx = -w{, rj; w). (2.17) z=0
Давайте из точки ( ,) исследуемой области в макроскоростной модели выпустим 2 луча в направлении свободной поверхности. Назовём один луч "s - лучом"1, а другой луч "g - лучом"2. Вдоль этих специфических лучей давайте сконструируем Гауссовы пучки: uf \X,Z;OJ), щ9 \x,z;u). Теперь рассмотрим краевую задачу: &ws + " .w3 = 0; (2.18) dws(x,z;u) cg(x, z) Oz = т?ь {х;и); (2.19) z=0 4-условие на бесконечности (см. [50],[4]), (2.20) 1м8мозначает "source -англ. источник", в некотором смысле этот луч должен быть протрассирован в направлении расположения источников. 2как и "s - луч", этот луч может быть назван как "geophone - луч" еорЬопе- англ. приёмник), и он должен быть протрассирован в направлении расположения приемников. здесь Ts (х;и)— это производная от гауссова пучка щ9 (х,z;u) по нормали к свободной поверхности, вычисленная на ней. Потому как гауссов пучок является асимптотическим решением уравнения Гельмгольца, решение краевой задачи (2.18)-(2.20) в окрестности "s-луча" совпадает с построенным выше гауссовым пучком, а вне окрестности "s-луча" решение этой краевой задачи отличается пренебрежимо мало от построенного гауссова пучка(см. [2]) (здесь важно, что макроскоростная модель является достаточно гладкой). Поэтому можно утверждать, что решение краевой задачи (2.18)-(2.20) практически совпадает с функцией щ3 (X,Z;LJ): uf (X,Z;UJ) = WS(X,Z;U ) + 0(UJ 2). Аналогично, и/ (x,z;u ): щ \x,z ,u) = wg(x,z;u ) +0(w 2), удовлетворяет краевой задаче: Au)9 + WJf9 = 0; (2 21) dwg(x,z;ui) +условие на бесконечности (см. [50],[4]), (2.23) Тд(х;иі)— это производная от гауссова пучка щ9 \x,z\u) по нормали к свободной поверхности, вычисленная на ней. Теперь, двойное применение теоремы Грина к гауссовым пучкам и функции Грина, приводит к соотношениям, аналогичным (2.17): - / т (х-,ш) G(x,z ,rr,u)dx = u b\ rj-,u) + 0(и У, (2.24) - / T b\x-,u)-G(x,z]t,V-,u,)dx = ugb\t rJ-,u) + 0(u). (2.25) Jz=0 Применение соотношений (2.24), (2.25) к данным многократного перекрытия (2.13), описываемых интегральным приближением (2.12): G(xg,0;f,г)\u )ddr)dxgdxs = / т 9Ь)(х3;и) т зЬ\хд;и) ф{хд;xs;uj)dxgdxs Xu J Хід приводит нас к ключевому соотношению rX2s rX-z J Х\. J Х\а 2ш2 [ [ К&тш).Щ$г= [ " [ 2\w\xs;u) (xa]u)-cP(xg-xs-,u)dxgdxs + 0(u,-i), J Jx CQ( V) Jxla Jxlg (2.26) где т,Я-и) = - .uy\t,r,]W)uW(t,ri;u ). (2.27) В дальнейшем при использовании соотношения (2.26) добавка 0(и 2") будет опущена. Как можно видеть, задача сводится к разрешению линейного интегрального уравнения первого рода относительно функции . Теперь, чтобы упростить анализ интегрального уравнения (2.26), введем лучевые координаты для обоих лучей и используем представление гауссовых пучков в лучевых координатах: W4t( .,q.)M .,QM - Jф ош) ехр т(з »схр уг ( 2}- (2-28) Ugb\aSs,Qs)A s,qM = )} %}. (2.29) Здесь индексы "s" и "g" обозначают параметры, связанные с "s-лучом" и "g-лучом" соответственно. В частности, текущая точка интегрирования у = (,т])(см. Рис.2.3) имеет координаты (sg,qg) и (ss,qs) в лучевых координат, связанных с g-лучом и s-лучом соответственно. Длины лучей ss, sg равны нулю на свободной поверхности в точках (XQ3, 0) и (х0д, 0) соответственно. ss0, sgo это длины лучей между точками (хоз, 0) и (ХІ, Zi), а также {х0д,0) и (xi,Zi) соответственно. Напомним, что система кинематической трассировки луча выглядит следующим образом: .-& = - )Vco(s); (2-30) Ш(0) = х0; р{0) = р0. Функция Га(за) : r„(ss) = Pa(s3)Qj1(ss), где параметры Qs(ss) и P3(s3) являются решениями задачи Коши для следующей системы дифференциальных уравнений: —-CQ{SS,0)PS; — = ——T -Z Q3; (2.31) dss dss c (s3,0) oqj qa=o Qs(0) = QQ; P,(0) = P0. Начальные данные Q0 и P0 выбираются таким образом, что минимальная ширина гауссова пучка находится в точке ss = SSQ, которая в то же время является исследуемой точкой (ХІ, zi) (см. [68],[43]). Уравнения для "g-луча" и для динамических параметров Qg, Рд соответствующего гауссова пучка остаются такими же как (2.30), (2.31), с точностью до замены индекса а на индекс д. Следует подчеркнуть, что для того, чтобы избежать разрушения гауссова пучка, нужно следить за тем, чтобы его ширина была больше доминирующей длины волны. Так как гауссовы пучки зависят от лучей, в окрестности которых они сосредоточе-ны(т.е. от параметров задающих лучи: координат (xi,Zi),Xog,xos), то будем писать:
Формула обращения, и её свойства
Теперь несложно показать, что оператор Т содержится в классе псевдодифференциальных операторов bj"o(X х X х R2)(cu. [34]). Действительно, очевидно, что ao = 1 Є 5j0(X х X х Л2). Это означает, что TQ Є LQlfi{X х X х R2). Мы показали, что функции А2 Є 5j0, #2 Є Зі 0. Кроме того А2 - 0(\х - у\1), Н2 = 0(\х - у\2). Это означает, что А2 Щ 1 Є «S o1, кроме того А2 Щ 1 = 0(\х - у\2п 1)- Отсюда следует, что ап 6 Si0(X х X х R2), однако пользуясь свойством (см. [34]), приведенном ниже: Лемма 1. Пусть А - псевдодифференциалъпый оператор, заданный формулой Аи{х) = f f е«х-у а{х, у, e)u{y)dyd9. (2.34)
Функция а(х,у,в) Є SS(X х X х R2). Тогда, если а(х,у,в) имеет при х = у нуль порядка к, то оператор А можно записать в виде , где вместо а(х, у, в) стоит Ь(х, у, в) Є S-kip-s)(X xXxR2). заключаем, что оператор Т% Є Ь\ (Х х X х R2).
Таким образом, мы приходим к тому, что оператор М13 представляется в виде суммы всё более гладких псевдодифференциальных операторов:
Стоит отметить, что оператор из Lj g увеличивает на п количество производных функции, к которой он применен, что означает, что функция из Ск(Х) переходит в функцию из Ск+п(Х). В частности, результатом применения оператора Т к, искомой нами, кусочно-непрерывной функции а, является п — 1 раз непрерывно-дифференцируемая функция. Из этого следует, что результатом применения оператора М — 1рат к функции а является непрерывная функция. Т.е. мы можем сказать, что оператор М0 эквивалентен оператору частичного восстановления Ij]ar: Мр « I" .
Эквивалентность « следует понимать в том смысле, что действие операторов Мр и 1&аг на разрывы функции практически одинаково. Подводя итог, мы можем сказать, что по формуле обращения: /аг гиг рХгз Дг9 / / / / T \xs;u;a;(3) b\xg;u;a;p) (x9;xs;u) ОСІ J W\ J Xls JXlg .e-Mr(xo9;x)+T(x0s;x)] . HdXgdXadujda = M0 l (Ш) J0 l (ї) (2.35) W Co Co мы можем восстанавливать амплитуды разрывов функции -. Теперь, давайте обратим внимание на действие оператора 1 аг. Он является "почти" единичным оператором. Действительно он может быть интерпретирован как суперпозиция прямого и "почти" обратного преобразования Фурье: /r SL (г) = /" / [ [ Щ- eW- »dgdp = [ [ e - S/S Jdp. 0 У JXparix) J JX СОКУ) J JXpaAx)
To, что может быть восстановлено, определяется множеством ХраГ. В частности, этот оператор не меняет функцию, если носитель фурье-образа этой функции принадлежит множеству Храг. И наоборот, если носитель фурье-образа этой функции не принадлежит множеству Храг, 0п лежит внутри ядра этого оператора. Например, если в среде имеется плоская граница: -І--Т- = аь 5(т] - ) -» ci/co(p) = аь 5{рх + рг) -» supp CI/CQ = {рх,рг : рх = -#»}. со (У) И например, если в среде имеется точечный рассеиватель: col = ая s(У - Узе) - ci/co(p) = asc -+ supp а/со = {р Є Я2}.
Если теперь вспомнить структуру множества Храг (2.32), то становится ясно, что для восстановления амплитуды границы достаточно взять 45 — є а 45 + є. Это означает, что мы должны ориентировать начальные направления "s" и "g" лучей на эту границу(т.е. в соответствии с законом Снеллиуса). Амплитуда же рассеивателя восстанавливается при любой ориентации начальных направлений лучей. Однако амплитуда рассеивателя восстанавливается почти всегда частично, да и бывают ситуации, когда амплитуда рассеивателя меньше амплитуды границы, и поэтому их не "видно". В этом случае можно построить селективные изображения с разной ориентацией начальных направлений, на которых например при а = 0 проявятся только лишь рассеиватель.
Нужно отметить, что параметр (3, который является углом раствора между начальными направлениями остался "свободным", т.е. можно сказать, что построен набор операторов по параметру /3. На самом деле, выбор этого параметра (соответственно оператора) будет определяться системой наблюдения и геометрией лучей в макроскоростной модели. Основными условиями, которые должны быть выполнены при выборе 0, являются "s(g) y4HM, и гауссовы пучки, построенные вдоль этих лучей, идущие из целевой области должны быть, такими, чтобы из набора источников(приемников) можно было сформировать гауссов пучок.
В дальнейшем выбор /3 будет производиться исходя из этого условия. Ясно, что этот выбор будет определяться системой наблюдения и макроскоростной моделью (условия освещения).
Прежде чем применять метод к синтетическим, а уж тем более к реальным данным, нужно его опробовать на простых примерах. Убедиться в том, что теория подтверждается практикой на элементарных объектах, и какие могут возникать практические трудности. Именно этому и посвящен настоящий параграф.
В качестве простой тестовой модели была выбрана среда изображенная на Рис.2.5. Она заполнена водой, плюс имеется несколько слабоконтрастных слоев, которые в природе соответствуют слоям, образованным осадочными породами.
Целевая область, в которой будем строить изображение приведена на Рис.2.7. Эта область содержит несколько слабоконтрастных горизонтальных слоев. Полная и макроско-ростная модель в целевой области, а также, искомый нами, "рефлектор" и контраст отображены на Рис.2.8.
Прежде чем строить изображение, нужно определиться с тем, какие углы а и /3 (см. Рис.2.3) выбрать. С одной стороны, диапазон этих углов определяется системой наблюдения. С другой стороны, если вспомнить структуру множества Храг (2.32), то становится ясно, что для восстановления амплитуды разрывов, соответствующих горизонтально залегающим границам, достаточно взять 0 — є а 0 + е. Это означает, что мы должны ориентировать начальные направления "s" и "g" лучей на эту границу. А угол /3 выбирается из тех соображений, чтобы в данных можно было сформировать гауссов пучок. В нашем случае (3 выбран равным 18. Тем самым мы строим селективное изображение среды, которое в данном конкретном случае одновременно является и полным. Нужно заметить, что при обработке реальных данных, когда неизвестно, как ориентированы границы, ясно, что нужно выбирать, как можно более больший интервал углов а.
Построение изображений в простой части модели Sigsbee2a...
Напомним, что алгоритм построения изображения завязан на трассировании гауссовых пучков из исследуемой точки до свободной поверхности. Поэтому прежде чем переходить к построению изображений в простой части, давайте сначала разберемся с трассированием в этой области.
Трассирование гауссовых пучков подразумевает решение двух систем обыкновенных дифференциальных уравнений (2.30) и (2.31). Стандартными наиболее употребительными численными методами решения задачи копій для обыкновенных дифференциальных уравнений являются методы Руиге-Кутта [3], которые и были использованы в численных экспериментах. Напомним, что алгоритм построения изображения предполагает трассирование гауссова пучка изнутри среды и вычисление производной на поверхности. Также напомним, что при изложении теории, утверждалось, что построенный таким способом Гауссов пучок является решением соответствующей краевой задачи (2.18)-(2.19). Поэтому при проведении численных экспериментов разумно было практически проверить это утверждение. Для этого данная краевая задача (2.18)-(2.19) была решена во временной области конечно-разностными методами. Решение, полученное таким способом, можно считать в определенной степени точным. Т.е. если верно был построен Гауссов пучок, то конечно-разностное решение краевой задачи должно совпасть с ним. На Рис.(3.5)приведены примеры для простой части среды, где изображены луч и ширина гауссова пучка, протрассированного изнутри среды. И в тоже время отображено как распространяется Гауссов пучок, построенный как конечно-разностное решение краевой задачи. Видно, что максимум амплитуды конечно-разностоного решения действительно находится на луче, и распространение происходит действительно внутри ширины гауссова пучка. Кроме того фокусируется конечно-разностное решение в том же месте, что и протрассированный Гауссов пучок. А также кривизна волнового фронта конечно-разностного решения меняет знак при переходе через точку фокусировки, как и должно быть.
Давайте теперь перейдём непосредственно к построению изображения. Хотя для восстановления полного изображения контрастности требуется суммирование по биссектрисе угла между начальными направлениями лучей, для начала мы ограничились восстановлением селективных изображений, соответствующих фиксированным ориентациям лучей. Напомним, что для восстановления контраста плоской границы достаточно сориентировать гауссовы пучки на эту границу. Стратиграфическая модель в целевой области содержит в основном горизонтальные границы, поэтому в первом эксперименте гауссовы пучки были сориентированы на горизонтально залегающие границы, т.е. а = 0. Угол раствора между начальными направлениями был задан равным 2(5 = 36. Построенное таким образом селективное изображение приведено на Рис. З.б. Видно, что верно восстановлено расположение границ. Также видно, что амплитуды границ затухают по мере того, как нормаль к границе начинает уклоняться от биссектрисы угла между начальными направлениями гауссовых пучков, чего и следовало ожидать. На изображении также частично проявились разломы и рассеивающие объекты. Разлом но сути состоит из набора "углов", на которых происходит рассеяние, поэтому они проявились, как набор рассеивателей. Амплитуда же рассеивателей восстанавливаются при любой ориентации начальных направлений лучей. Однако, это частичное восстановление. Для того, чтобы получить более полное изображение среды, были построены ещё два селективных изображения. В одном из них гауссовы пучки были ориентированы на границы, которые находятся в левой нижней части целевой области и одновременно на разлом, который находится в серединной части целевой области, т.е. а = 18. Построенное таким образом селективное изображение и ориентация гауссовых пучков отображены на Рис.3.7. Видно, что проявились амплитуды границ, па которые были сориентированы пучки, а также снова проявились разломы и рассеиватели. В другом случае гауссовы пучки были ориентированы на разлом, который находится в правой срединной части, т.е. а = —18. Угол раствора между начальными направлениями и в том, и в другом случае по-прежнему равен 2/3 = 36. Построенное таким образом селективное изображение и ориентация гауссовых пучков отображены на Рис.3.8. Видно, как и ранее, что проявились амплитуды границ, на которые были сориентированы пучки, а также снова проявились разломы и рассеиватели.
После того, как все эти полученные селективные изображения были просуммированы в соответствии с формулой обращения (2.35), мы получили более полное изображение среды, которое можно видеть на Рис.3.9. Видно, что все особенности среды восстановлены. О качестве восстанавливаемых амплитуд разрывов контраста можно судить по Рис.3.10. Видно, что амплитуды восстановлены верно.
Полученные численные результаты позволяют нам говорить, что для несложных реалистичных моделей сред метод не содержит ограничений, кроме тех, которые накладывает на него геометрия системы наблюдения. Более того, продемонстрирована возможность получать селективные изображения, которые тоже могут оказаться практически значимыми.
Построение изображений рассеивающих объектов
При построении изображений как для скалярного случая, так и для упругости важным инструментом являются гауссовы пучки. Гауссовы пучки для изотропной упругой среды уже достаточно давно получены [13] и долгое время используются. Формальные асимптотические решения, сосредоточенные в окрестности луча для анизотропной теории упругости впервые построены в работе [19]. После публикации этой работы прошло более 20 лет, и за это время появились работы, в которых обнаружены трудности, связанные с так называемыми сингулярностями или сингулярными направлениями. Имея ввиду эти трудности, были использованы только квазипродольные гауссовы пучки. В этом случае не возникает сингулярностей. Далее, следуя работе [19] получим формулу для квазипродольных гауссовых пучков. Будем искать решение уравнений эластодинамики: (C«« ) + / AJ«I = 0, j = 1,2,3; (4.1) i,l,k=l г в виде : иі(х,ш) = еіш х)Ц(х), / = 1,2,3. (4.2)
Подставим (4.2) в (4.1), приведем подобные при одинаковых степенях и и будем удовлетворять уравнению (4.1) до первой степени и включительно.
Чтобы уравнение (4.3) имело нетривиальное решение необходимо, чтобы detN = 0, а это возможно, когда хотя бы одно из собственных значений тензора Грина-Кристоффеля: OS dS - " дХідХк равно 1. Нужно отметить, что три различных собственных значения соответствуют трём разным типам волн: квазипродольной и двум квазипоперечным. Введём обозначения для этих собственных значений: Gqp(x,VS), Gq3l(x,VS), GqS2(x,VS). Напомню, что мы ограничились отысканием лишь квазипродольного решения, поэтому решать будем уравнение Gqp(x,VS) = l. (4.5) Теоретически возможен случай, когда Gqp(x, VS) = 1 = Gq31(x, V5), однако в природе это практически не встречается, поэтому мы рассматриваем случай Gqp(x,VS) ф Gqsi(x,VS),Gqp(x,VS) ф GqS2(x,VS). Это означает, что квазипродольная волна отделена от квазипоперечных, т.е. квазипродольная всегда быстрее квазипоперечных. Введём обозначение: H{x,VS)=l-Gqp{x,VS). Н(х, VS1) называется функцией Гамильтона или Гамильтонианом. Уравнение (4.5) относительно функции S, известное как уравнение Эйконала, принимает вид: Н(х, V5) = 1/2. (4.6) В дальнейшем второй аргумент функции Гамильтона будем обозначать через р: Н(х,р). Давайте теперь V будем искать в виде: v = 1/V: If {р і ЇЇ) = 1) 9 Є KerN (dim(KerN) = 1). Т.е. мы получили следующее уравнение на плотностнонормализованиый вектор поляризации: N t = 0, t- Ш, t) = 1- (4.7) Тогда уравнение (4.4) примет следующий вид: М ifip = Tf = (М 7f у? , if) = 0. /.ж _ ч v dcy« 95 95 дді OS dgi , d2S . {M g p ,g) = —--- рт + с..ы .( + + - - )+ i,j,l,k=l dxi dxk Dxidxk OxkOxi ОХІОХ dS dip dS dp sr дсцы dS dS dg, dS dg{ d2S , dS dip d{cijkl-gl9j) dS dip dH. 0Я В изложенной выше цепочке преобразований были использованы следующие свойства: E dS dgi v-v dS dgj - . dH dS cHki -4—9i= У. ciiki—---gi\ Cijkigigjpi:=—-,pi= —-.
В результате мы получили уравнение относительно функции р, известное как уравнение переноса: 2(V ,—) + (V, - -)=0. (4.8)
Теперь для того, чтобы построить решение уравнения (4.1), мы должны решить уравнение Эйконала (4.6) и уравнение переноса (4.8). Далее мы будем искать комплексные решения этих уравнений в окрестности фиксированного луча. Давайте построим луч х(т) или, что то же самое, характеристику уравнения Эйконала, которая является решением системы уравнений Гамильтона-Якоби [17]: dx _ dH(x,p) dp _ ЭН(х,р) dr dp dr dx Теперь перейдём к лучевым координатам (хі,Х2,х$) —» (г, q\,q2), (pi,P2,Pi) — (рТіРчх,Р42)-Будем иногда обозначать г = g3- Базисные векторы лучевых координат удовлетворяют небезызвестным формулам Френе [12].
Далее, следуя работе [13] будем предполагать, что qt = 0(u;-1//2), і = 1,2, и, как и ранее, будем удовлетворять уравнению (4.1) по всем степеням о; до 1 степени включительно. Уравнения на (7-,01,02) и gi(T,qi,g2) мы получили при второй степени со, а уравнение на ір(т, 1,) при первой степени и: 84 и2 : N -$ = 0; det(N) = 0 & H(q, —) = 1/2; ш : М iplf, 3 = 0.
Т.к. решение мы ищем в виде ряда по степеням q\ и q2, то для того, чтобы удовлетворить уравнению (4.1) по всем степеням w до 1 степени включительно мы должны найти квадратичное приближение эйконала S(r,qi,q2) и вектора поляризации gi(T,qi,q2) и нулевое приближение функции р(т, qi,q2) H(q, Щ) = Ъ + 0(ш- 2) = S = r + 1(Го,о) + 0(с"3/2); М р?, = 0 + 0{ 2) = ф, q) = + 0(Ш- У, y/c{T)det{Q{T)) N = N0 + N.u-"2 + N2u l + 0(L -3/2) =» ф,gi,q2) = g? + g}u 2 + g2u l + 0(uT3/2). Нужно отметить, что зная найденную функцию S = г + (Tq,q), мы можем вычислить N NQ + AW"1/2 + N2u l = [ aijkl - - - Sjt]. , , OXx OXk Окончательно, получаем формулу для гауссова пучка в анизотропной неоднородной среде: ИдЬ(т,Яг,Я2,и) = е 1 -л ==Щг,диЯ2). у/с(т) det Q(T)
Прежде чем перейти к построению изображений рассеивающих объектов, давайте сначала рассмотрим обратную задачу на определение расположения источников в упругой анизотропной среде.
Пусть в упругом неоднородном анизотропном полупространстве X = {(х\,х2, хЗ) : %з 0} в точке х3 расположен источник типа сосредоточенной силы, действующей в направлении / (см. Рис.4.1). Процесс формирования и распространения волнового поля, порожденного этим источником описывается краевой задачей:
