Численное решение уравнения двойного корня в криволинейных координатах Текст научной статьи по специальности «Науки о Земле и смежные экологические науки»

ЧИСЛЕННОЕ РЕШЕНИЕ УРАВНЕНИЯ ДВОЙНОГО КОРНЯ В КРИВОЛИНЕЙНЫХ КООРДИНАТАХ

Александр Сергеевич Сердюков

Новосибирский государственный университет, 630090, Россия, г. Новосибирск, ул. Пирогова, 2, старший преподаватель кафедры Высшей математики Физического факультета, к.ф.-м.н., тел. +79137679513, e-mail: AleksanderSerdyukov@yandex.ru

Антон Альбертович Дучков

Институт нефтегазовой геологии и геофизики им. А.А. Трофимука СО РАН, 630090, Россия, г. Новосибирск, пр. Ак. Коптюга, 3, к.ф.-м.н., зав. лаб., e-mail: DuchkovAA@ipgg.sbras.ru

Уравнение двойного корня описывает продолжение времен пробега отраженных волн по глубине (пересчет кинематических данных по глубине). Соответствующее эволюционное уравнение по глубине и его аппроксимации используются для реализации сейсмической миграции до суммирования. Само уравнение двойного корня является уравнением типа Г а-мильтона-Якоби и может быть решено численно для расчета поля времен отраженных волн сразу для всей системы наблюдений многократных перекрытий. В настоящей работе получена и программно реализована численная схема (WENO-RK) решения уравнения двойного корня в криволинейных координатах. Она позволяет рассчитывать поле времен отраженных волн на основе обобщенного принципа взрывающихся границ. Полученные результаты будут полезны для решения различных кинематических задач в случае наличия топографии и резких отражающих границ.

Ключевые слова: уравнение двойного корня, поле времен, принцип взрывающихся границ, WENO-схема.

NUMERICAL SOLVER FOR DOUBLE SQUARE ROOT EQUATION IN CURVILINEAR COORDINATES

Alexandr S. Serdyukov

Novosibirsk State University, 630090, Russia, Novosibirsk, Pirogova str. 2, assistant lecturer, department of Physics, tel. +79137679513, e-mail: AleksanderSerdyukov@yandex.ru

Anton A. Duchkov

A.A. Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, 630090, Russia, Novosibirsk, Prospect Akad. Koptyuga 3, head of the laboratory, e-mail: DuchkovAA@ipgg.sbras.ru

Double-square-root (DSR) equation describes kinematics of downward data continuation in depth. Corresponding evolution equation in depth and its approximations are used for implementing prestack depth migration. The DSR equation is of Hamilton-Jacobi type and can be solved numerically for the whole acquisition system simultaneously. We have developed WENO-RK type numerical solver for DSR equation in curvilinear coordinates. These scheme allows us to compute reflection traveltimes, using the exploding reflector concept. The obtained results could be used for solving different kinematic problems in the presence of topography and high-contrast reflectors.

Key words: double-square-root equation, reflected wave traveltimes, exploding reflector concept, WENO-scheme.

Введение

Уравнение двойного корня (double-square-root, DSR) описывает продолжение поля времен отраженных волн по глубине [1]. Соответствующее эволюционное уравнение по глубине и его аппроксимации используются для сейсмической миграции до суммирования [2].

Уравнение двойного корня может рассматриваться, как уравнение Г амиль-тона-Якоби (характеристическое уравнение), описывающее распространение слабых разрывов (волн). От классического уравнения эйконала его отличает то, что по глубине пересчитываются времена отраженных волн, т.е. учитывается набег времен как по лучу падающей, так и отраженной волны. Таким образом, применение уравнения двойного корня позволяет пересчитывать поле времен с одной глубины на другую при условии, что лучи не имеют горизонтальных сегментов (условие DSR).

По аналогии со стандартным уравнением эйконала для численного решения уравнения Гамильтона-Якоби можно использовать WENO-RK схему [5]. В случае гладких криволинейных границ расчетной области возможно переписать DSR уравнение и соответствующую численную схему в криволинейных координатах. Это позволяет пересчитывать поле времен с одной криволинейной границы на другую.

В работе рассмотрен обобщенный принцип взрывающихся границ, который позволяет рассчитывать поле времен отраженных и головных волн сразу для всей системы наблюдений многократных перекрытий. Приведены примеры расчета отраженных и головных волн.

DSR кинематика для криволинейных границ

Уравнение двойного корня в декартовых координатах имеет вид [1, 2]:

где т(хя, хг, г) - время двойного пробега отраженной волны между источником с горизонтальной координатой х3 и приемником с горизонтальной координатой хг, находящимися на одной глубине г, т г = дт / дг, тхг = дт / дхг, тх5 = дт/ дх3. Таким образом, уравнение (1) полезно для пересчета времен с одной горизонтальной плоскости на другую, причем вычисления проводятся в расширенном (не физическом) пространстве (х5, хг, г, £) .

Во многих задачах возникает необходимость расчета времен отраженных и головных волн от криволинейной границы на криволинейную поверхность наблюдения (учет топографии). Рассмотрим криволинейную отражающую границу, г = Н0 + f (х) , и поверхность наблюдения г г = д (х). Применим следующее преобразование координат [4] (см. рис. 1):

4V (xr , z)~2 ~Tl) , (i)

z - g (x)

(2)

Уравнение (1) может быть переписано для криволинейных координат (2), пересчет времени идет в пространстве (%3, %г , у, Ь) между поверхностями у = с о ш Ь. Новое уравнение имеет вид:

Ту = н°БН ( % 3, %г , у, т ^, т^) = (? ( % 3, у, т(з) + (? ( %г , у, т ^г), (3)

где Ту = д т / д у, т = д т / д %г , т ^3 = д т / д %3

-т^С(^у) + V 2({, у)(С2(\', У) + Г2(\')) - т}р2({)

<? (%,у,т() =------------------------------------------------------------------------І (4)

с2^,г) + Р2(0 здесь V - скорость в новых координатах:

V (%,у) = V (%,уН°+Г У'9 (х) + д (х) ).

Но

(5)

Функции ¥(£), 0(%,у) - это компоненты матрицы Якоби преобразования (2)

«о , „ а'«)(2-/ет-н0)-2гет

р(0 =

-,С(?,У) = Н0:

Но + ґ(х) - 0(хУ ~" "и (Но + Г(О-0(())2

где г для заданных % и у может быть найдено из (2).

(6)

Рис. 1. Исходная расчетная область и криволинейная сетка в декартовой системе координат (слева) и соответствующая регулярная сетка в новых координатах (справа)

WENO схема для криволинейного DSR уравнения

’^ЕКО схема для уравнения (3) с переменной у в качестве эволюционного параметра конструируется следующим образом [5].

• Аппроксимируются «горизонтальные» производные т^,т^ в каждой точке сетки для рассматриваемого у. Эти производные могут иметь разрыв, так

где НI есть частная производная Н по /-му аргументу и мы опускаем зависимость от % и у в обозначениях. В данной работе мы использовали функцию потока Лакса-Фридрихса [5].

• Модифицированное уравнение с функцией потока в правой части (

Для расчета кинематики отраженных (головных) волн для заданной отражающей границы можно использовать концепцию взрывающейся границы для ББЯ уравнения [3]. При этом происходит инициализация численной схемы для ББЯ уравнения в форме линейного распределенного источника в пространстве (^ %г> У > С), форма которого совпадает с формой границы.

Отраженные волны. Пусть отражающая граница задана уравнением . Тогда в криволинейных координатах граница соответствует . для расчета времен пробега отраженных волн мы можем использовать функцию источника в терминах взрывающейся границы (при нулевом выносе х8 = хг и времени 1 = 0 ):

Физический смысл здесь состоит в том, что в нулевой момент времени двойного пробега (С = 0 ) падающие и отраженные лучи совпадают (х8 = хг) на отражающей границе. Из выражения (7) мы так же можем получить начальные условия для DSR уравнения.

Головные волны. Для задания головных волн нужно отдельно вычислите времена пробега в нижнем слое между различными точками границы. Затем может быть использована следующая формула взрывающейся границы:

что вычисляются левая и правая производные с использованием техни-

ки ’^ЕКО-аппроксимации.

• Вычисленные производные затем используются для вычисления

ту, при этом гамильтониан Нвш (т^,т%г) , заменяется функцией потока

Н(т+ т~ т+ т~ ) =

ту = Н(ті, , т^г ,т^г)) интегрируется по у методом Рунге-Кутта.

Принцип взрывающихся границ

/ (%з,%г,у,0 = 8 (у-уо ) 8 (%3~%Л 8 ( Ь) .

(7)

Примеры расчета времен пробега отраженной и головной волн

Пример расчета поля времен отраженных волн показан на рис. 2. Слева приведена скоростная модель, представляющая собой градиентную среду с криволинейной отражающей границей. Справа показан результат расчета поля времен (двойного времени пробега) отраженных волн для системы наблюдений многократных перекрытий с помощью WENO - схемы. Поле времен (рис. 2, справа) показано для 2=0 в координатах средней точки и половины удаления источник-приемник ( х = (хг+хх)/2, к = (хг - х^) / 2 ).

На рис 3. приведен пример расчета времен пробега головных волн при помощи прдлагаемого алгоритма для неровной поверхности наблюдений. Слева показана скоростная модель и лучи головных волн. Справа показано поле времен головной волны в координатах ( ).

Рис. 2. Скоростная модель (слева) и времена пробега отраженной волны (справа)

Времена двойного пробега (в секундах)

Скоростная модель и лучи

Рис. 3. Скоростная модель и лучи (слева) и рассчитанные времена пробега головной волны (источник-граница-приемник) (справа)

Заключение

В работе получено уравнение двойного корня в случае криволинейных координат и предложена численная схема его решения. Данные результаты могут быть использованы в ряде кинематических задач, в частности для скоростного анализа верхней части разреза или для задач сейсмологии.

Авторы выражают благодарность T. Alkhalifah, J. Virieux, F. Andersson и А.М. Айзенбергу за конструктивное обсуждение результатов. Работа была поддержана в рамках Соглашения с KAUST № 1029-3 и Шведским фондом по международному сотрудничеству в науке и высшем образовании (the Swedish Foundation for International Cooperation in Research and Higher Education).

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Белоносова А.В., Алексеев А.С. Об одной постановке обратной кинематической задачи сейсмики для двумерной неоднородной среды // Некоторые методы и алгоритмы интерпретации геофизических данных. М. 1967. С. 137-154.

2. Claerbout, J., 1985, Imaging the earth's interior: Blackwell Scientific Publishing.Clayton, R., 1978, Common midpoint migration: Technical Report, Stanford University, SEP-14.

3. Duchkov, A. and M. De Hoop, 2010, Extended isochron rays in prestack depth \map migration: Geophysics, 75 (4), S139-S150.

4. Fornberg, B, 1988, The pseudospectral method; accurate representation of interfaces in elastic wave calculations; Geophysics 53, 625-637.

5. Osher, S. and C. Shu, 1991, High-order essentially nonoscillatory schemes for Hamilton-Jacobi equations; SIAM Journal on Numerical Analysis, 907-922.

© А.С. Сердюков, А.Д. Дучков, 2013