Содержание к диссертации
Введение
ГЛАВА 1. Изученность вопроса 8
ГЛАВА 2. Динамически обоснованный кинематический оператор (ко) 14
2.1. Построение динамически обоснованного кинематического оператора 14
2.2. Сопряженный КО. Миграция временных невязок 28
2.3. КО в случае вертикально-неоднородной вмещающей среды .41
ГЛАВА 3. Совместное обращение кинематического и динамического операторов в случае вертикально -неоднородной вмещающей среды 52
3.1. Динамический оператор для вертикально-неоднородной вмещающей среды 52
3.2. Составной оператор для вертикально - неоднородной вмещающей среды ' 56
3.3. Численная реализация совместного обращения 58
Заключение , 65
Литература
- Сопряженный КО. Миграция временных невязок
- КО в случае вертикально-неоднородной вмещающей среды
- Составной оператор для вертикально - неоднородной вмещающей среды
- Численная реализация совместного обращения
Введение к работе
В настоящее время в обработке сейсмических данных можно выделить две основные процедуры - определение скоростного строения среды и миграция данных сейсмических наблюдений в среду с известным скоростным строением. Относительно миграционных процедур к настоящему времени сложилось ясное понимание их природы и унифицированное описание как результата действия сопряженного оператора для соответствующей прямой динамической задачи. Миграционные процедуры позволяют получать достаточно подробные изображения исследуемой среды, которые характеризуют положение и форму отражающих объектов, но не позволяют судить ни о распределении скоростей над отражающими границами, ни об абсолютных значениях локальных скоростных неоднородностей. Все попытки итерационного применения этих процедур для полного восстановления среды не увенчались успехом. Как было показано в ряде работ за последние десять лет, это связано со спецификой строения сингулярного спектра упомянутого оператора - оказалось, что сингулярные векторы, отвечающие за плавное изменение распределения скоростей в среде, лежат в недоступной для устойчивых вычислений области. В то же время, из практики обработки сейсмических данных известно, что использование только кинематической составляющей волнового поля позволяет устойчиво определять именно плавную (трендовую) компоненту скоростного разреза. В некотором роде сложилась парадоксальная ситуация -кинематическая информация, безусловно, содержится в полных динамических данных, но, фактически, не работает при их обращении. В связи с этим представляется актуальным предложить такой подход при обработке динамической информации, который нацелен именно на выделение и дальнейшее использование кинематической информации.
Объектом исследований в данной работе является один из возможных подходов к решению линеаризованной обратной кинематической задачи сейсмики, который позволяет учитывать динамику распространения сейсмических волн.
Цель исследований - разработка способа обращения невязок времен пробега сейсмических волн, учитывающего динамические характеристики волновых полей, а также спектральную ограниченность сигналов, используемых в сейсморазведке.
Задачи исследований: достоверное восстановление скоростного строения среды с помощью обращения невязок времен пробега, которые определяются с учетом динамики распространения сейсмических волн.
Основные этапы исследований:
получить выражение для линеаризованного кинематического оператора (т.е. оператора, переводящего возмущения скоростного строения вмещающей среды в невязки времен пробега), возникающего при реализации динамического подхода и изучить его основные свойства;
на примере синтетических и реальных данных показать, что применение динамического подхода при обращении временных невязок позволяет более полно учитывать эффекты, связанные с распространением реальных сейсмических волн и позволяет восстанавливать трендовую составляющую скоростного разреза;
показать, что в случае вертикально - неоднородной вмещающей среды существует возможность совместного обращения кинематических (невязки времен пробега) и динамических (невязки волновых полей) данных; при этом восстанавливается более полная информация о скоростном строении, чем при независимом обращении.
Методы исследований и фактический материал
Основным методом исследования является математическое моделирование волнового поля от точечного источника в двумерной среде с помощью начально-краевой задачи для скалярного волнового уравнения. Численное моделирование включает задание двумерной модели, источника, системы регистрации и выполняется с помощью конечно-разностных методов. Теоретической основой решения поставленных задач является построение сингулярного спектра компактного оператора и изучение структуры его устойчивых подпространств,
то есть подпространств, натянутых на старшие сингулярные векторы. Обработка данных предусматривает получение матричного представления оператора и сводится к реализации процедуры устойчивого решения систем линейных уравнений, основанной на выполнении сингулярного разложения.
Основным фактическим материалом являются результаты обработки данных численного моделирования процессов распространения сейсмических волн. Для опробования алгоритмов и ключевых процедур использовались как синтетические данные, так и данные реальных сейсмических наблюдений по методу вертикального сейсмического профилирования, полученные в Западной Сибири (предоставлены СБ. Горшкалевым в 2001 г.).
Сформулированы и защищаются следующие научные результаты
1. Линеаризованный кинематический оператор
строится на основе применения динамического подхода к определению временных невязок, который заключается в том, что задержка целевой волны в реальной среде относительно волны того же типа, рассчитанной для референтной модели, определяется как аргумент максимума функции взаимной корреляции этих двух сигналов. Рассеянная компонента волнового поля описывается с помощью борцовского приближения в спектральной области
сводится к линейному интегральному уравнению Фредгольма первого рода относительно возмущения скорости.
2. Решение обратной кинематической задачи
основано на аппроксимации интегрального оператора конечномерной системой линейных алгебраических уравнений
строится с привлечением регуляризирующей процедуры, которая заключается в построении приближенного решения на основе усеченного SVD - разложения матричного представления компактного оператора
позволяет устойчиво восстанавливать треидовую компоненту скоростного разреза.
3. Совместное обращение кинематических и динамических данных в случае
вертикально - неоднородной вмещающей среды
заключается в построении динамического и кинематического операторов, действующих на один и тот же элемент пространства моделей
решении совместной системы линейных алгебраических уравнений, составленной из матриц, которые возникают при конечномерной аппроксимации динамического и кинематического интегральных операторов.
Новизна работы. Личный вклад.
Предложен оригинальный подход к решению линеаризованной обратной кинематической задачи сеисмики, который существенно опирается на использование динамической информации, содержащейся в сейсмических записях.
С использованием метода возмущений получено линейное приближение для кинематического оператора; изучены свойства ядра чувствительности этого оператора; в случае вертикально-неоднородной вмещающей среды получено его представление в виде однопарамстрического семейства одномерных интегральных операторов.
Проанализировано строение сингулярного спектра данного оператора в сравнении с сингулярным спектром линеаризованного оператора динамического обращения и показано, что устойчивые области этих спектров дополняют друг друга. На основании этого предложен новый подход к совместному обращению кинематических и динамических данных (для вертикально-неоднородной вмещающей среды).
Проведена представительная серия численных экспериментов, позволившая выяснить разрешающую способность метода динамического обращения временных невязок и определить границы его применимости. Выполнена обработка реальных данных, позволившая уточнить скоростное строение околоскважинного пространства и существенно улучшить результаты миграции сейсмических данных. На примере синтетических данных показано, что предложенный подход к совместному обращению кинематических и динамических данных позволяет получать более полную информацию о строении среды.
Теоретическая и практическая значимость результатов
Предложенный динамический подход к обращению временных невязок может применяться при обработке реальных данных и позволит в определенных ситуациях (особенно в тех, когда нельзя пренебречь эффектами, связанными с ограниченностью спектров реальных сигналов) получать более достоверную информацию о распределении скоростей в среде, чем стандартные методы лучевой томографии. В случаях поверхностных систем наблюдения, когда в качестве референтной модели используется вертикально-неоднородная вмещающая среда, существует возможность совместного обращения кинематических и динамических данных, что существенно улучшает информативность получаемых результатов.
Апробация работы и публикации
Основные положения и результаты докладывались па международной конференции молодых ученых, специалистов и студентов "Геофизика - 2001" (Новосибирск, 2001), международной конференции "Математические методы в геофизике" (Новосибирск, 2003).
Результаты исследований по теме диссертации изложены в 4 опубликованных работах.
Диссертация выполнена в Лаборатории динамических проблем сейсмики Института геофизики СО РАН.
Автор выражает искреннюю признательность научному руководителю, к.ф.-м.и. В.А. Чсверде за всестороннюю поддержку и постоянное внимание. Хочется поблагодарить сотрудников Лаборатории динамических проблем сейсмики Института геофизики СО РАН к.ф.-м.н. В.И. Костина, к.ф.-м.н. В.Г. Хайдукова, Д.М. Вишневского за помощь в работе, а также д.ф.-м.н. СИ. Кабанихина, к.т.н. Ю.Л. Орлова, д.ф.-м.н. Б.П. Сибирякова зато, что они взялись за труд ознакомиться с работой и высказать о ней свое мнение.
Сопряженный КО. Миграция временных невязок
Численное решение линейного интегрального уравнения (2.1.18) в случае произвольной вмещающей среды приводит к необходимости обращения больших, плохо обусловленных матриц. Для решения подобных задач успешно применяются итерационные алгоритмы. Одним из наиболее эффективных является алгоритм LSQR [48]. Первый шаг этого алгоритма - действие сопряженным оператором на правую часть уравнения, в данном случае на временные невязки. Получим выражение для сопряженного оператора.
Прежде всего отметим, что выражение, стоящее в знаменателе (2.1.18) является известной функцией и может быть перенесено в левую часть: Щг8,?,) = Л(Я{г))= frf -dl J(/ii)3- (rr,rf;iv). (rlf iv).G( r,;w)Av (2.2.1) D V0 \h) -НІ
Определим скалярное произведение в пространстве моделей M = L2(D). Пусть w,(г) є М, т2(?) є М , (т{ (г),т2 (г))м = J/и, m2dr (2.2.2) В пространстве данных Л/, = L2 (D,) . Пусть рх (rs, rg) є М,, р2 (г,, rg) є Л/, {P\irs rg)tp2{rl9rg)) = \\px-p2drsdrg (2.2.3) По определению сопряженного оператора имеем: (А( ),5Г{г,,г,))иі =(Л(г) -(Щг„г,))) = = \\sf(rt,rg)drgdrs Hg2«/ \( y u0(fg,r/tW).u0(rs, v)-G(C,rg;w)dw = D, D V0 \Q) -w2 = f d? \\ST(r,9rg)drgdr, 1( .5-0 , ).110( , ).0( ) D V0 \h) D, -w} A\W(rs,rg)) = \\tf{rs,rg)drgdrs (iw)3.i/0(rf,rf;W).i/0(r, ;W).C(C, ;iv)rfw (2.2.4) Перепишем это выражение в несколько другом виде: Обозначим выражение в квадратных скобках через W. W удовлетворяет уравнению Гельмгольца AW(r,w) + k(r)2W(r,w) = iw-u0(rs,rg;w)-3T(rs,rg)-S(r-?g) (2.2.5) с условием излучения, соответствующим приходящей волне. Используя равенство Парсеваля и переставляя интеграл по приемникам получим: A\ rgJs)) = -— \drs\-W{rsJj) u0{rs,r,t)dt (2.2.6) где W является решением задачи 1 d2W A„r г д v02(r) dt2 " J" V1" dt dW W\ = = AW+ jsf(fs;rg)-—u0(rsjg,t)drg (2.2.7) = 0 t=T dt
Здесь T - время регистрации данных, достаточно большое, чтобы считать наблюдаемое волновое поле существенно ослабленным при t Т.
Таким образом, функция W(rs,r,t) представляет собой "продолжение в обратном времени" данных ( ,» »0 = ( "Т о О изо всех приемников внутрь среды с известным скоростным распределением v0(r). Напомним, что величина 5T(r5,rg)связана с невязками времен 5T(rs,rg) = ST(rs,rg)- jw2\u0(rs9rg;w)\ dw (2.2.8)
Значение функции 5v{r) є М: 5v(r) = А 1ST) для каждой фиксированной точки из D представляет из себя сумму по всем положениям источника временных корреляций, "продолженных в обратном времени" данных Ф(гх гв 0 и поля, вызванного точечным источником, расположенным в точке rs. Действие сопряженного оператора на временные невязки напоминает процедуру миграции до суммирования, поэтому далее будем называть его миграцией временных невязок. Миграция временных невязок использовалась для решения ряда задач. Миграция временных невязок
Локализация геологического объекта типа кимберлитовой трубки В качестве исходных данных были взяты синтетические волновые поля, рассчитанные для сложно построенной среды, которая моделировала геологические условия Восточной Сибири. (На Рис.14 представлен фрагмент модели - та часть, которая использовалась в эксперименте). Применялась следующая система наблюдения. В скважине глубиной 200 м. находятся 20 источников с шагом 10 м. Сигнал от каждого из источников (импульс Риксра, и0 =70Hz) регистрируется в четырех скважинах, расположенных симметрично относительно первой на расстоянии 200 и 400 м., в каждой из которых находятся по 20 приемников с шагом 10 м. Эта система наблюдения перемещается вдоль профиля с шагом 200 м. (Рис.13.). В качестве фоновой была взята двухслойная модель. Граница сглаживалась, чтобы не возникала мешающая отраженная волна. Для миграции используются невязки времен пробега прямой волны.
КО в случае вертикально-неоднородной вмещающей среды
Рассмотрим теперь случай, когда вмещающая среда является вертикально-неоднородной, причем предполагается, что функция v0(z) непрерывна вплоть до целевой отражающей границы z = Н. Источники и приемники заполняют всю прямую z = zs. Возмущение вмещающей среды Sv(x, z) есть финитная функция с носителем, сосредоточенным в квадрате [-L,L]x[h,H], h zs, которая должна быть найдена по данным системы многократного перекрытия u(xg,zg;xs,zs;t) = g(xg,x/,t) (2.3.1) Функция u(x,z;xs,zs;t) в (2.3.1) есть решение волнового уравнения 1 я2 г , л, 2TT = Au + f(0S(x-xs,z-zs) (2.3.2) [v0(z) + 3v(x,z)] дґ с нулевыми начальными данными. Будем считать, что функция f(t) достаточно гладкая, удовлетворяет условию /(0) = f\Q) = 0 и имеет ограниченный спектр, т.е. F(u ) = 0, wg[w,,Yv2 -oo
Введем следующее преобразование u(x,z;xs,zs ,t) (в 2-D случае это преобразование Радона [16]): u(x,z;a;r) = ju(x,zxs\T-a-xs)dxt (2.3.3) -00 где u(x,z;a\r) будет решением задачи і я2 - —г- = Au+f(t-a-x)- S(z - z.) [v0(z) + Sv(x,z)]2 dt2 5« "I „ = — (2.3.4) = 0 /=0
Эта функция описывает процесс распространения волны, которая изначально была плоской. Действительно, если положить v(x, z) = v0 = const (такое предположение справедливо, например, вблизи источника), то решение (2.3.4) запишется как u(x,z;a;T) = -Q(T-a-x-f]-\z-zs\), /3 = ]\-а2 (2.3.5) Р Yvo где б(г) = )f{r)dz о
Определение временной невязки, данное в (2.1.3) будет справедливо и для сигналов, преобразованных с помощью (2.3.3). Пусть uQ{x,z\a\t) - решение (2.3.4) для фоновой скорости v0(z), a Su(x,z;a;t)- его возмущение, вызванное неоднородностью Sv(x,z). Тогда по аналогии с (2.1.11) получим выражение для временной невязки преобразованных сигналов: hrv ,ч 5м0(дг ,z ;a;/) = д2її (x z «-,) 2 3-6) —00
Воспользовавшись равенством Парсеваля, с учетом ограниченности спектра функции /(/), запишем (2.3.6) в частотной области: jiw SU(xg, zg; a; w) U0 (xg ,zg\a; w) dw ST{Xg,a) = --- (2.3.7) jw2 -p0(xg,zg;a\\v)\2dw Функция U0(x,z;a;w) удовлетворяет неоднородному уравнению Гельмгольца AUQ + - -UQ = - eiwaxF{w)-8{z-zs) (2.3.8) и может быть записана как U0(x,z;a;w) = ei axF(w)-G(z1zs;a-iw) (2.3.9) где G(z,zs;a ,w)- функция Грина для обыкновенного дифференциального уравнения dlG г dz -a .2 2 W ) G = -S(z-zs) (2.3.10)
Подставляя (2.3.9) в (2.3.7) получим jVw SU(xg, ze; a; w) F(w) -G(zg,zs; a; w) e"waXg dw ЯХ ,. ) = -— : (2.3.11) ]w2.\F(w)\2.\G(zg,zs;a;w)\2dw Функция SU(x,z;a;w) является решением ДЯ7 + --т и7 = 2-dv{x,z) - e,waxF{w)-G{z,zs\a\w) (2.3.12) v000 v03(z) Выполним преобразование Фурье для задержки ST(xg;a) по переменной xg . Из (2.3.11) следует \iw 8U{k8x -a-w,zg;a;w)-F(\v)-G(zg,zs;a;w)dw 57W,a) = -— (2.3.13) jw2 -\F(wf .\G(zg,zs;a;w)\2dw SU(kx,z;a;w) - преобразование Фурье $U(x,z;a;w) по переменной xg, которое удовлетворяет обыкновенному дифференциальному уравнению: + \ -k2\8U = -Sv(kx+a-w-,zyF(w)-G(z,z/,a-w) (2.3.14) dz \y0{z) ) v0(z) Sv(kx,z) - преобразование Фурье Sv(x,z) по x. Отсюда, с помощью определенной в (2.3.10) функции Грина, получим следующее представление для 5U(kx,z;a;w): SO(kx,z;a;w) = 2w2F(w) Г(к" \" К .G( z,;a;w)-G(z; ; -;w)dC (2.3.15) vo(0 w
Теперь выражение (2.3.13) может быть записано как одноиараметрическое семейство (кх -параметр) линейных интегральных уравнений, связывающих преобразование Фурье от известных невязок времен ST(kx;a) с преобразованием Фурье искомого возмущения вмещающей среды 5v(&f \ z)
Составной оператор для вертикально - неоднородной вмещающей среды
Расчеты проводились для вмещающей среды, представляющей собой слой мощностью Я = 2 со скоростью v0 = 2.5, лежащий на однородном полупространстве с v,=5. Область реконструкции D выбрана как D = {(x,z):xe[-\;l],z є [0.4; 1.6]}. Линия наблюдения расположена при z = 0. Функция источника - импульс Рикера с доминирующей частотой и0 = 30 Hz. В качестве базиса в пространстве моделей используем систему из 301 ступенчатой функции с шагом Иг = 0.004. В область D вносилось возмущение в виде эллипса с центром в точке х0 = 0 , z0 = 1 с полуосями а = 0.8 , Ъ = 0.25 и перепадом скорости 5v = 0.2. Правые части рассчитывались как результат действия соответствующих матриц на вектор Svj(z), являющийся коэффициентом Фурье в разложении (3.1.7). При обращении матриц использовался алгоритм усеченного SVD - разложения. На Рис.43-46 представлены результаты независимого обращения операторов Aj и Aj , для j = 3,6 а также результат их совместного обращения в сравнении с точной гармоникой. (При построении этих решений использовались усечения, которые обеспечивали обусловленность систем порядка 1000). На Рис.47 приведены результаты суммирования ряда (3.1.7) с N = 6 точных и восстановленных гармоник (при совместном обращении) (без нулевых). Отмстим, как точно определяются границы восстановленной неоднородности. В то же время и ее абсолютная величина весьма близка к той, которая дается рядом (3.1.7) с у = 1,...,6 (Рис.48). К сожалению, в случае однородной вмещающей среды нельзя воспользоваться совместным обращением для восстановления нулевой гармоники (см. Гл.2 п.З).
Выясним, как располагаются подпространства, натянутые на старшие правые сингулярные векторы рассматриваемых операторов, относительно друг друга. Для этого вычислим раствор между ними. (Раствор между подпространствами - совокупность углов, которые составляют друг с другом ортогональные базисные векторы из каждого подпространства [4]).
В силу того, что обусловленность систем бралась равной 1000, решение Л3(Г) лежит в span iT)j,i = \...,\0, а решение A\D) лежит в span I, j, і = 1,...,80. Раствор для этих подпространств представлен на Рис.49 Раствор между подпространствами, в которых отыскиваются решения для А6 и Л , представлен на Рис.50. Как видно, в этих подпространствах имеется некоторое количество взаимно ортогональных элементов. Таким образом, при совместном обращении происходит пополнение дополнительными базисными векторами того подпространства, в котором отыскивается решение. inaex Рис.50 Раствор между spanfyT)\ / = 1,...,14 и spanffiD)}, / = 1,...,80. Рассмотрим пример с более сложной вмещающей средой. Пусть имеется двухслойная среда Рис.51. Внутри первого слоя имеется возмущение в виде эллипса с центром в точке х0 = 0 , z0 = 0.35 с полуосями а = 0.8 , b = 0.08. Линия наблюдений находится при z = 0.005. Для расчета исходных данных всего использовался 41 источник (функция в источниках - импульс Рикера с доминирующей частотой и0 =30 Hz) и 401 приемник, которые располагались с равномерным шагом на всем интервале [-2;2]. Фоновая модель представлена на Рис.52. При расчете кинематического оператора граница при z = 0.7 сглаживалась; в качестве целевой волны использовалась волна, отраженная от второй границы при z = 1.
Область реконструкции D = {(x,z) :х є [-1.5; 1.5], z є [0.1;0.9]}. На Рис.53-56 представлены результаты независимого обращения операторов Aj и Aj , для у = 2,7 а также результат их совместного обращения в сравнении с точной гармоникой. При построении этих решений использовались усечения, которые обеспечивали обусловленность систем порядка 10000. При использовании больших значений числа обусловленности решение начинает разрушаться вследствие неизбежных ошибок в исходных данных. На Рис.57 приведены результаты суммирования ряда (3.1.7) с N = 1 точных и восстановленных гармоник (при совместном обращении) (с нулевой гармоникой).
Основным результатом работы является разработка оригинального подхода к обращению невязок времен пробега сейсмических волн, учитывающего динамику их распространения. Его характерная особенность заключается в том, что для поверхностной системы наблюдений в случае использования вертикально-неоднородной вмещающей среды удается свести исходное интегральное уравнение (в общем случае по трехмерной области) к распадающейся системе одномерных линейных интегральных уравнений Фредгольма первого рода. Представленные результаты численных экспериментов по анализу сингулярного спектра линейных конечномерных операторов, возникающих при дискретизации этой системы доказывают возможность устойчивого восстановления трендовой компоненты скоростного разреза при наличии значительных помех во входных данных. Сравнение сингулярного спектра вышеупомянутого кинематического оператора с сингулярным спектром линеаризованного оператора динамического обращения показывает, что устойчивые области этих спектров дополняют друг друга. Этот факт служит теоретическим обоснованием для предложенного подхода к совместному обращению кинематических и динамических данных. Он заключается в построении динамического и кинематического операторов, действующих на один и тот же элемент пространства моделей и решении совместной системы линейных алгебраических уравнений, составленной из матриц, которые возникают при конечномерной аппроксимации этих операторов. Численные эксперименты показывают, что при реализации данного подхода происходит значительное улучшение получаемых результатов по сравнению с результами независимых обращений.
Численная реализация совместного обращения
Построение кинематического оператора
Для последующих выкладок нам будет необходимо переписать (2.1.11) в спектральной области. Но прежде сделаем несколько предположений. Если используется волна, отраженная от некоторой известной границы, будем считать, что v0(r)- гладкая функция, не имеющая разрывов вплоть до этой границы. В случае, когда используется прямая волна (наблюдения в скважине) будем считать v0(г) гладкой. Далее, полагаем, что носитель dv(r) содержится в некоторой области D, целиком лежащей в верхней полуплоскости, а источники и приемники располагаются вне этой области. Теперь, воспользовавшись равенством Парсеваля, запишем J/w SU(rg ;rs;w)-U0 (rg ;rs;w) dw W(rg;rs) = - (2112) jw1 ]UQ(rg;rs;w)\2dw Формула (2.1.12) устанавливает связь между невязкой времени пробега волны из точки rs в точку rg с возмущением волнового поля в точке fg. Чтобы получить окончательное представление кинематического оператора, необходимо выразить невязку ST через возмущение скорости 5v(r). Для этого найдем выражение для возмущения волнового поля du(r,rs;t). В качестве математической модели для описания процесса распространения волн выберем начально - краевую задачу для скалярного волнового уравнения. Наблюдаемое волновое поле uobs(. s iO удовлетворяет следующей начально- краевой задаче: Z\\2 д,2 (v0(r) + 5V(r))2[u0+du] = A[uQ+du] + f(t)d(?-rs) (2.1.13)
Пренебрегая членами второго порядка малости, с учетом (2.1.13а) получим, что возмущение du(rtrs\t), вызванное наличием скоростной неоднородности, удовлетворяет следующей начально-краевой задаче:
Отмстим непрерывность ядра оператора в области D, так как особенности у функций Грина, входящих в выражение (2.1.18), могут возникнуть только в точках г = rs и г = rg. Кинематический оператор для однородной вмещающей среды
Выражение (2.1.18) представляет собой интегральное уравнение первого рода с непрерывным ядром. Как известно, задача, связанная с его решением, относится к классу некорректных. Поэтому для численного решения (2.1.18) необходимо привлечение регуляризирующих процедур. В связи с этим полезно рассмотреть свойства ядра оператора, а также провести SVD - анализ его матричного представления. Это позволит выяснить структуру решения, и то, какая часть этого решения может быть восстановлена при наличии ошибок во входных данных и погрешности аппроксимации оператора конечномерной системой линейных уравнений. Однако для получения матричного представления оператора (2.1.18) в случае произвольной вмещающей среды необходимо вычислять функции Грина для каждого фиксированного положения источника, что весьма трудоемко. Поэтому рассмотрим здесь простейший случай однородной вмещающей среды.
Функция Грина для двумерного уравнения Гельмгольца [3]: G(f,w) = HJ,2)(k-\r\) (2.1.19) Функция Грина выбирается в зависимости от знака при экспоненте прямого преобразования Фурье ди(г\rs\t). Воспользуемся асимптотическим выражением для функции Ханкеля: Н«\х) = Л\— ехр 7DC Л -IX- X -» (2.1.20) Тогда для функций Грина, входящих в (2.1.18), получаем следующее представление: G(r;f;w)= (br-f) = i (тг---.ехр (2.1.21) 4Н?Чк-\г-ё\) = - I . . -ехв -/"" Ч v0 н J J Вводя обозначения Л = - , ,=rf-r, /?2 = rg-л получим следующее представление кинематического оператора для однородной вмещающей среды: ST(rg,rs) = \KTtf,?s ,rt )Sv(OdC (2.1.22) D с ядром K(\rs;rs): Кт(Л,г,) = I W V2 »\vo X; R Л /?2 f F(w)f sin R-Rx-R2)-« »2 Q\F(WfdW л dw (2.1.22a) Зависимость временных невязок от положения возмущения относительно луча
Особый интерес представляет поведение ядра уравнения (2.1.22а) в окрестности луча. Для его изучения был проведен следующий численный эксперимент (Рис.2). В окрестность луча при фиксированном расстоянии 1 от источника вносилась скоростная неоднородность — = 10% в виде круга с радиусом R —.
После чего конечно-разностным методом считалась прямая динамическая задача и находилось поле иоЫ в точке rg. По формуле (2.1.3) определялась временная задержка, которая приписывалась соответствующему положению возмущения. В качестве функции источника использовался импульс Рикера /"(/) = (l-2;r2Uo/2)-exp(-;r2L o/2) с доминирующей частотой и0 =30 Hz. Остальные параметры: v0=2km/scc, 5v = +0.2 /ти/sec L = lkm, R = 0.03 km. Полученные зависимости для двух значений 1 представлены на Рис.3.
Кривые зависимости ST(s). /", - радиус первой френелевской зоны для частоты v = SO Hz; гг - радиус первой френелевской зоны для частоты v = ЗОЯг
Максимальная задержка возникает не в том случае, когда возмущение находится непосредственно на луче, а в некоторой его окрестности, в тот момент, когда возмущение начинает пересекать первые френелевские зоны, рассчитанные для частот из верхней части значимого диапазона (Рис.4). Весьма неожиданно то, что при определенных положениях возмущения невязка меняет знак. Подобные зависимости получаются и при использовании других функций источников /(/), например, если f(t) = s m(2nv0t)/2KU0t или сигнала, зарегистрированного вблизи взрывного источника при проведении ВСП. Этот результат отличается от того, что получалось в [42,43], где возмущение непосредственно на луче вообще не вызывало задержки, что по-видимому объясняется использованием в упомянутой работе низкочастотных сигналов с очень ограниченным спектром (частотный диапазон 0.01-0.06 Hz). Зависимость временной задержки от положения возмущения в пространстве описывается функцией чувствительности ядра. Она определяется при вычислении ядра уравнения (2.1.22) в каждой точке и достигает максимума в тех точках, в которых внесение скоростного возмущения вызывает наибольшую задержку (таким образом, графики ST(s) на Рис.3, представляют собой сечения функции чувствительности на заданном расстоянии от источника). Функция чувствительности для рассматриваемого примера представлена на Рис.5. Отметим повышение чувствительности в окрестности приемника и источника (этот факт описывался в [59,60]), в то время как в лучевой томографии чувствительность постоянна вдоль луча.