CC BY
12
УДК 550.34.013.4
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ КИНЕМАТИКИ СЕЙСМИЧЕСКИХ ВОЛН В СЛОЖНОПОСТРОЕННЫХ СРЕДАХ
Павел Юрьевич Степанов1 Юлия Александровна Гоманюк2
1 Московский государственный университет имени М.В. Ломоносова, Москва, Россия; [email protected], https://orcid.org/0000-0002-8131-8998
2 Московский государственный университет имени М.В. Ломоносова, Москва, Россия; [email protected], https://orcid.org/0000-0002-1253-4689
Аннотация. Представлен вариационный метод математического моделирования кинематики сейсмических волн, разработанный на кафедре сейсмометрии и геоакустики геологического факультета МГУ имени М.В. Ломоносова для исследования особенностей кинематики сейсмических волн разных типов в двумерных изотропных средах (градиентных и слоистых). Задача определения траекторий сейсмических лучей решалась путем интегрирования методом Рунге-Кутта системы дифференциальных уравнений с заданными начальными условиями. Выполнено исследование алгоритма с целью проверки точности и корректности получаемых решений, а также его опробование на ряде теоретических моделей неоднородных сред. Разработанный метод лучевого трассирования применен для изучения влияния градиента скорости и геометрии сейсмических границ на кинематику отраженных волн в многослойных средах. По результатам математического моделирования кинематики отраженных волн сделаны выводы о пределах применимости упрощенных моделей горизонтально-слоистых сред, которыми часто аппроксимируют сложнопостроенные неоднородные среды.
Ключевые слова: лучевое трассирование, сейсмический луч, неоднородные среды, вариационное исчисление, кинематика сейсмических волн
Для цитирования: Степанов П.Ю., Гоманюк Ю.А. Математическое моделирование кинематики сейсмических волн в сложнопостроенных средах // Вестн. Моск. ун-та. Сер. 4. Геология. 2022. № 6. С. 167-178.
MATHEMATICAL MODELING OF SEISMIC WAVE KINEMATICS IN COMPLEX MEDIA
Pavel Yu. Stepanov1^, Julia A. Gomanyuk2
1 Lomonosov Moscow State University, Moscow, Russia; [email protected], https://orcid.org/0000-0002-8131-8998
2 Lomonosov Moscow State University, Moscow, Russia; [email protected], https://orcid.org/0000-0002-1253-4689
Abstract. A variational method of mathematical modeling of seismic wave kinematics has been developed at the Department of Seismometry and Geoacoustics of the Faculty of Geology of Lomonosov Moscow State University for studying the kinematics of seismic waves of different types in two-dimensional isotropic media (gradient and layered) is presented. The problem of determining the trajectories of seismic rays was solved by integrating using the Runge-Kutta method a system of differential equations with given initial conditions. The algorithm was studied in order to verify the accuracy and correctness of the solutions obtained, as well as its testing on a number of theoretical models of heterogeneous media. The developed ray tracing method was used to study the effect of the velocity gradient and the geometry of seismic boundaries on the kinematics of reflected waves in multilayer media. Based on the results of mathematical modeling of the kinematics of reflected waves, conclusions were drawn about the limits of applicability of simplified models of horizontally layered media, which often approximate complex inhomogeneous media.
Key words: ray tracing, seismic ray, inhomogeneous media, variational theory, kinematics of seismic waves For citation: Stepanov P.Yu., Gomanyuk J.A. Mathematical modeling of seismic wave kinematics in complex media. Moscow University Geol. Bull. 2022; 6: 167-178. (In Russ.).
Введение. Лучевое трассирование (ray tracing) — один из основных подходов к решению прямых кинематических задач сейсморазведки. Прослеживание сейсмических лучей волн разных типов позволяет детально изучить характер распространения упругих волн в сложнопостроенных средах, построить систему волновых фронтов и, самое главное, с высокой точностью рассчитать время пробега волны от источника до любой точки геологической среды, вовлеченной в волновой процесс.
Методы решения прямых кинематических задач сейсморазведки в области лучевого трассирования, развитие которых началось в последние десятилетия ХХ в., не утратило своей актуальности и в настоящее время [Virieux, Farra, 1991; Cerveny, 2001; Gpystdal et. al., 2002; Облогина et. al., 2004; Rawlinson et. al., 2007]. В задачи сейсморазведки, которые решаются с применением методов лучевого трассирования, входят расчет синтетических сейсмограмм, решение томографических задач и задач скоростной инвер-
сии, проведение кинематическом миграции данных МОВ по опорным отражающим горизонтам, глубинно-скоростное моделирование и др.
Все методы лучевого трассирования делятся на два основных типа: лучевое трассирование с заданием начальных условий (initial-value ray tracing) и лучевое трассирование с заданием граничных условий (boundary-value ray tracing) [Cerveny, 2001]. Первый тип алгоритмов подразумевает, что при прослеживании луча известны его направление и текущая координата, т. е. известны начальные условия для решения системы дифференциальных уравнений, описывающих траекторию луча.
В качестве начальных условий для методов лучевого трассирования с заданием граничных условий (boundary-value ray tracing) задаются точки или поверхности, через которые должен пройти сейсмический луч, при этом направление луча в каждой его точке неизвестно. Наиболее распространенные методы решения задачи лучевого трассирования с заданием граничных условий — методы пристрелки [Yang, 2003; Галактионова, Белоносов, 2020] и методы изгиба луча [Moser et. al., 1992; Xu et. al., 2010].
Кроме того, существует ряд методов, в которых рассчитываются время прихода волн и волновые фронты, однако не рассчитываются (или рассчитываются в качестве вспомогательных вычислений) сейсмические лучи. К этой группе методов относятся сеточные методы кратчайшего пути (network shortest-path ray tracing) [Klimes, Kvasnicka, 1994], конечно-разностные методы решения уравнения эйконала [Cao, Greenhalgh, 1994; Тульчинский, 2009], метод быстрой маршировки [Sethian, Popovici, 1999] и др.
Нами был использован вариационный метод математического моделирования кинематики сейсмических волн, разработанный на кафедре сейсмометрии и геоакустики геологического факультета МГУ имени М.В. Ломоносова, который относится к методам лучевого трассирования с заданием начальных условий [Облогина, 1998; Облогина и др., 2004].
Цель работы — применение разработанного вариационного алгоритма лучевого трассирования для исследования особенностей кинематики сейсмических волн разных типов в двумерных изотропных средах (градиентных и слоистых).
Математическая постановка задачи. В основе геометрической сейсмики лежит известный принцип Ферма, который гласит: действительным путем распространения волны из точки А в точку В является тот путь, вдоль которого по сравнению со всевозможными другими путями, соединяющими точки А и В, время пробега волны минимально. Время пробега волны из точки А в точку В можно представить в виде
<•« ds )
-1' v(x, y)' (1)
t =
где ds — дифференциал дуги кривой, соединяющей точки А и В; v(x, у) — скорость распространения
волны в среде. Время пробега можно рассматривать как минимальное значение функционала Ферма
t[У(x)] = dx,
v( x> У)
(2)
где XQ и XQ — абсциссы точек А^0, у0) и B(x1, у1).
Определение траекторий сейсмических лучей в среде с заданной скоростью v(x, у) сводится в соответствии с принципом Ферма к нахождению тех функций у(ж), которые реализуют минимум функционала [Облогина, 1998].
Представим функционал Ферма в виде
Xj _
t[У(х)] = J a(x, y)yjj + y'2, a(x, y):
J
v(x> У)
(3)
с граничными условиями
УС^ = Уo, У^О = У^ (4)
означающими, что кривые у(%) проходят через заданные граничные точки А^0, у0) и B(x1, у1) (например, через источник и точку наблюдения).
Из теории вариационного исчисления [Эль-сгольц, 1969; Облогина и др., 2004] известно, что по необходимому условию экстремума функционала кривые у^) должны удовлетворять уравнению Эйлера:
Fy — Fxy — Fyyy' — Fyyy" = 0. (5)
Интегральные кривые этого уравнения называются экстремалями. Таким образом, сейсмические лучи можно рассматривать как экстремали, реализующие минимум функционала Ферма.
Преобразуя левую часть уравнения Эйлера и полагая при этом, что
F (X, y (x), y' (x)) = a (x, y)^J + y
(6)
получим
Fy" diF"= "Л77
(ax - "уУ)у'
ay
VT
1 + y (J + y")2 После простых преобразований будем иметь
(7)
y' (x) = (J + y' 2( x))
д ln a д ln a
dy
dx
У'(x)
(8)
Полученное уравнение — дифференциальное уравнение лучей. Его общее решение у = у^, С1, С2) зависит от двух произвольных постоянных, которые определяются из граничных условий y(xQ)=yQ, У(xl)=Уl.
Дифференциальное уравнение лучей интегрируется в конечном виде только для частных простых
видов скоростной функции у(х, у), когда скорость зависит только от х или только от у. Это обстоятельство вынуждает прибегать к численным методам интегрирования обыкновенных дифференциальных уравнений.
Явная форма представления кривых у=у(х) неудобна при применении численных методов интегрирования. Чтобы обойти указанные затруднения, перейдем к параметрической форме задания кривых сравнения. В качестве параметра возьмем t — время пробега волны вдоль луча. Такой параметр удобен еще и тем, что при интегрировании вместе с координатами луча х и у одновременно определяется и время пробега t от источника до рассматриваемой точки наблюдения [Ермаков, Степанов, 2018].
Для применения широко известного метода Рунге-Кутта для интегрирования дифференциального уравнения лучей необходимо уравнение (8) свести к эквивалентной системе дифференциальных уравнений 1-го порядка. Эта система имеет вид
dx(t )
dt
= v cos У,
dy • о — = v sin 9, dt
d9(t) . 0 0
-= vx sin 9 - vy cos 9.
dt x y
(12)
Начальные условия для этих уравнений
x(t)|^ = x» y(t)|^ = Ус e(t)l^ =°о (13)
задают координаты x0, y0 источника и начальный угол 90 выхода луча из источника.
Таким образом, задача определения лучей в двумерно-неоднородных средах может быть поставлена и решена как задача Коши для приведенной выше системы дифференциальных уравнений 1-го порядка (12). Если вычисления вести с постоянным шагом по t, то решение получается в узлах сетки через равные промежутки At, что можно использовать для построения изохрон. Годограф строится по временам выхода лучей на земную поверхность у = 0 [Ермаков, Степанов, 2018].
Программная реализация алгоритма. Программа Ray Tracing, реализующая вариационный алгоритм лучевого трассирования на основе интегрирования системы дифференциальных уравнений лучей, разработана в среде ОС Windows 10 с помощью Microsoft Visual Studio Community на программной платформе Microsoft.NET Framework 4.7. Программа написана на языке С#.
Исследование вариационного алгоритма. Прямая кинематическая задача сейсморазведки может быть решена аналитическим и численным способами. Аналитический метод расчета траекторий сейсмических лучей представляет собой точное решение уравнений лучей. Численный метод за-
ключается в интегрировании уравнения Эйлера для функционала Ферма, определяющего время пробега луча в среде. Аналитические расчеты выполнены для однородной и вертикально-неоднородной среды с линейным изменением скорости с глубиной. Результаты расчетов аналитическим и численным методами совпадают с высокой степенью точности. Расчеты подтверждают, что рассматриваемый в работе численный алгоритм решения прямой задачи позволяет получить достаточно точные и корректные решения.
Существует несколько способов задания скоростной модели среды: скорость в произвольной точке среды может быть определена по аналитической зависимости, описывающей распределение скорости в среде, или скорость задана таблично в узлах сетки. Последний способ подразумевает использование интерполяционных алгоритмов для расчета значений скорости в любой точке модели среды. Решение прямой кинематической задачи при разных способах задания скорости в модели и сравнение результатов расчетов позволяет сделать вывод о точности и корректности решения, получаемого при табличном способе задания скоростей, по сравнению с точным аналитическим заданием значений скорости и их частных производных. Нами сделаны расчеты для горизонтально- и вертикально-неоднородной, а также двумерно-неоднородной моделей сред. Полученные лучевые траектории оказались полностью идентичны, что подтверждает корректность табличного способа задания значений скорости и используемого интерполяционного алгоритма.
Опробование алгоритма решения прямой кинематической задачи на теоретических моделях. Рассмотрим результаты применения вариационного алгоритма для расчетов лучей и годографов сейсмических волн в различных неоднородных средах, включая среды с границами 1-го и 2-го рода.
Вертикально-неоднородная среда. Вертикально-неоднородная модель среды (ВНС) — одна из базовых моделей сред, используемых в сейсморазведке. В некоторых ситуациях горизонтально-слоистая среда, которая чаще всего рассматривается в сейсморазведке, может быть описана вертикально-неоднородной моделью.
На рис. 1, а представлен результат решения прямой кинематической задачи для модели ВНС с линейным законом возрастания скорости с глубиной. Полученные результаты полностью согласуются с базовыми свойствами лучей и годографов в вертикально-неоднородных средах, известными из кинематической теории сейсморазведки:
1) луч в вертикально-неоднородной среде симметричен относительно вертикальной прямой, проходящей через точку максимального проникновения луча;
2) лучевые поля рефрагированных волн и соответствующие им годографы симметричны относительно вертикальной прямой, проходящей через центр расстановки;
Рис. 1. Результат решения прямой кинематической задачи для вертикально-неоднородной (а) и двумерно-неоднородной (б) сред (годографы и лучи рефрагированных волн)
3) кажущиеся значения скорости во взаимных точках встречных годографов совпадают.
Двумерно-неоднородная среда. Двумерно-неоднородная модель среды — естественное продолжение вертикально-неоднородной среды и отличается от нее наличием в среде горизонтального градиента скорости (который, как правило, существенно меньше, чем вертикальный).
Рассмотрим модель, скорость в которой задается следующим уравнением:
У = VQ (1 + к1х + к2х),
где к1 — горизонтальный градиент скорости, к2 — вертикальный градиент скорости. На рис. 1, б представлен результат решения прямой кинематической задачи для этой модели. Анализ результатов позволяет выявить основные отличия двумерно-неоднородной среды от вертикально-неоднородной (для фланговой системы наблюдения):
1) лучевые поля и годографы для прямого и встречного пунктов возбуждения не симметричны относительно центра расстановки;
2) точки пересечения лучей и годографов для прямого и встречного пунктов возбуждения смещены относительно центра расстановки;
3) времена во взаимных точках встречных годографов равны в соответствии с принципом взаимности, а значения кажущейся скорости различаются.
Вертикально-неоднородная среда с границей 2-го рода и повышенным градиентом скорости во втором слое. На рис. 2 представлены результаты решения прямой кинематической задачи для модели с одной границей 2-го рода, в которой градиент скорости во втором слое значительно превышает
градиент скорости в первом. Лучи рефрагированных волн, распространяющиеся только в первом слое, и соответствующие им участки годографов показаны на рис. 2 черным, лучи, проходящие через оба слоя, и их участки годографов — бордовым. Скорость в модели задается по линейному закону.
Лучи рефрагированных волн в случае постоянного градиента скорости (в первом слое) с увеличением угла выхода луча с поверхности наблюдений проникают на большую глубину и возвращаются на поверхность наблюдений на больших удалениях от пункта возбуждения, чем лучи с меньшим углом выхода. Однако из-за значительной разницы в градиентах значений скорости в слоях и высокого градиента во втором слое лучи рефрагированных волн, вошедшие в среду под большими углами, возвращаются из второго слоя на поверхность наблюдений на меньшем удалении от пункта возбуждения, чем лучи из первого слоя с меньшими углами выхода. В результате на годографе наблюдается петля возврата. На дальних удалениях лучи рефрагированной волны из второго слоя выходят в первые вступления, так как градиент скорости в этом слое значительно выше, чем в первом (рис. 2).
Вертикально-неоднородная среда с границей 1-го и 2-го рода и постоянной скоростью во втором слое. Рассмотрим примеры вертикально-неоднородных сред с меняющимися на границах значениями скорости и градиентами скорости, т. е. среды с границами 1-го и 2-го рода. На рис. 3 представлен результат решения прямой кинематической задачи для такой модели среды. Скорость в первом слое задается по линейному закону, во втором слое скорость постоянна (рис. 3, в).
а
3600 4200 6000
Координата приемника, м
Рис. 2. Результат решения прямой кинематической задачи для вертикально-неоднородной среды с границей второго рода: а — годограф рефрагированной волны; б — лучи сейсмических волн; в — скоростная модель среды
6
егас! У(2)
§гасНУг>
В первом слое градиент скорости постоянный, скорость на границе увеличивается скачком, во втором слое скорость постоянна, как следствие лучи представляют собой прямые линии, показанные на лучевых построениях зеленым. В первом слое при малых углах выхода луча с поверхности распространяются лучи рефрагированных волн (на рисунке черные). Существует «критический» луч, который касается границы и совпадает с «первым» лучом отраженной волны (при возрастании угла выхода лучей из источника). При больших углах выхода луча в среде распространяются отраженные волны, показанные на лучевых картинах бордовым. Годограф отраженной волны на больших удалениях от источника стремится к годографу рефрагиро-ванной волны и совпадает с ним в последней точке каждого годографа — точке прихода к поверхности
наблюдений луча, общего для отраженной и рефра-гированной волн.
При достижении критического угла падения луча на границу вдоль нее распространяется головная волна, первый луч которой совпадает с лучом отраженной волны. Поэтому на годографах отраженной и головной волн есть общая точка. Начиная с некоторых удалений головная волна выходит в первые вступления, на годографе этому моменту соответствует точка пересечения годографов реф-рагированной и головной волны (рис. 3).
Вертикально-неоднородная среда с волноводом. Эта модель содержит две границы 1-го и 2-го рода, скорость во втором слое постоянная и пониженная по сравнению с соседними слоями, т. е. второй слой представляет собой однородный волновод.
а
1000
се
5 ю
PÎ
[—I
1500
2000-
2500
2100 2800 3500
Координата приемника, м
G
erad V(z)
grad V(z)
Рис. 3. Результат решения прямой кинематической задачи для вертикально-неоднородной среды с границами 1-го и 2-го рода и постоянной скоростью во втором слое: а — годографы рефрагированной, отраженной и головной волн; б — лучи сейсмических волн; в — скоростная модель среды
На рис. 4 представлен пример решения прямой кинематической задачи для модели с однородным волноводом. Скорость в модели задается по линейному закону. В первом и третьем слоях распространяются лучи рефрагированных волн. Скорость в волноводе постоянна, поэтому лучи во втором слое представляют собой прямые линии. Понижение скорости в волноводе приводит к возникновению зоны тени — области, которая не освещена лучами и в которой наблюдается разрыв годографа.
Лучи с меньшими углами выхода из источника проходят в волноводе больший путь, чем лучи с большими углами выхода, поэтому приходят к поверхности наблюдений позже лучей, прошедших малый путь в волноводе. Последние в свою очередь быстрее возвращаются к поверхности наблюдений из-за того, что большую часть пути проходят в третьем слое с высоким градиентом скорости. В результате на годографе образуется ветвь возврата, верхняя часть которой соответствует лучам с малыми углами
выхода, прошедшим больший путь в волноводе, а нижняя — лучам с большими углами выхода, быстро вернувшимся на поверхность наблюдений из третьего слоя (рис. 4). Конфигурация годографа с волноводом может быть разной в зависимости от заданного скоростного поля — значений скорости, величины градиентов скорости и типа волновода.
Влияние градиента скорости и геометрии сейсмических границ на кинематику отраженных волн в многослойных средах. Математическое моделирование кинематики отраженных волн — еще одна задача, которую можно решить с помощью рассматриваемого вариационного алгоритма решения прямой кинематической задачи. Качественно и количественно проанализируем влияние разных факторов среды, таких, как наклон и форма отражающих и преломляющих границ, а также величины градиента скорости в слоях, на результаты решения прямой кинематической задачи. Эти факторы будут определять отличие рассчитанного годографа отра-
а
О
§ 14-<и
^ ю I
СО 12Н 1086420
^16
-14 -12 -10 -8 -6 4 \-2
4000
8000
12000 18000 20000 Координата приемника, м
Ъ
Рис. 4. Результат решения прямой кинематической задачи для вертикально-неоднородной среды с однородным волноводом: а — годограф; б — лучи сейсмических волн; в — скоростная модель среды
V
§гас! У(2)
женной волны от «эталонного» годографа для горизонтально-слоистой среды. Анализируя полученные результаты, можно сделать выводы о том, когда влиянием кривизны границ и неоднородностью среды можно пренебречь и работать в рамках упрощенных моделей горизонтально-слоистых сред, которые наиболее часто используются в сейсморазведке.
Для расчетов использована четырехслойная модель, приближенно описывающая геологическое строение верхней 100-метровой толщи разреза в районе МГУ имени М.В. Ломоносова (рис. 5, а). Траектории лучей и годографы волн, отраженных от трех границ, представлены на рис. 5, б, в. Рассмотрим последовательно влияние на форму годографов отраженных волн каждого из перечисленных выше факторов.
Влияние градиента скорости. Добавим в первые три слоя исходной модели вертикальный градиент скорости. Рассмотрим четыре варианта модели с разными градиентами скорости в слоях (рис. 6, а). Как видно на рис. 6, б, невязки годографов для исходной горизонтально-слоистой и вертикально-неоднородных моделей составляют от 7-10 мс (для минимального градиента скорости) до 40-50 мс (для максимального градиента).
Рассмотрим теперь влияние горизонтального градиента скорости на кинематику отраженных волн. На рис. 7, а представлены четыре двумерно-неоднородные модели с разными градиентами скорости (соотношение вертикального градиента
к горизонтальному составляет приблизительно 2 : 1).
Значения времени прихода отраженных волн для рассматриваемых моделей отличаются от исходной модели гораздо существеннее, чем в предыдущем эксперименте (от 20-30 мс для модели с минимальным градиентом до 70-140 мс для модели с максимальным градиентом) (рис. 7, б). Таким образом, видно, что даже при минимальном градиенте скорости величина невязки в значениях времени прихода отраженных волн слишком велика, чтобы наличием градиента скорости в модели можно было пренебречь.
Наклон отражающей границы. Оценим влияние наклона отражающих границ на форму годографов. На рис. 8 представлены лучевые поля для случая, когда третья граница модели наклонена на 5° (против часовой стрелки), годографы для вариантов наклона границы на 2 и 5° против часовой стрелки и на 10 и 20° по часовой стрелке. Отклонение рассчитанных годографов от «эталонного» (для горизонтально-слоистой модели) изменяется от 10-20 мс для наклона 2° до 50-150 мс для наклона границы на 20°. Если наклонять обе границы, то эффект от наклона каждой границы суммируется. Анализируя степень отличия рассчитанных годографов от эталонного горизонтально-слоистого варианта среды, можно сделать вывод, что среду можно рассматривать как горизонтально-слоистую для наклонов границ не более 2-3°.
О 20 40
I I
60
80
а
100
120
140
160
180
200
20
40-
60-
80
1004
120-
350 м/с 400 м/с Суглинок
Песок
55 м
500 м/с Глина
и™
Рис. 5. Исходная модель (а); годографы, соответствующие трем отражающим границам (б), траектории лучей волн, отраженных от третьей границы (в), для исходной модели; 1-3 отражающие границы
_1_I.
3(0 Вертикально-неоднородная модель 2 И» (1+0,0023 р
■ ■
Рис. 6. Сравнение годографов (б) отраженных волн для четырех вертикально-неоднородных моделей (а) с разным вертикальным градиентом скорости (номер годографа соответствует номеру модели). Черный годограф соответствует исходному годографу для горизонтально-слоистой модели
о
2 540 Н
6 520-1
си рр
500480460440420400
50
100 150 200
Коорднната приемника, м
\
-540 -520 -500 -480 -460 -440 -420 -400
Форма отражающей границы. Рассмотрим влияние формы отражающих границ на времена прихода отраженных волн. На рис. 9 изображены годографы для вариантов моделей с разной формой отражающей границы. Отклонение рассчитанных годографов от эталонного (для горизонтально-слоистой модели) составляет от 10 до 25-30 мс. Таким образом, пренебречь кривизной границы возможно лишь для перепадов высоты не более нескольких метров вдоль границы.
Влияние геометрии преломляющих границ. Изучим влияние наклона и кривизны преломляющих границ на времена прихода отраженных волн. Рассмотрим сначала влияние наклона одной из преломляющих границ (рис. 10, а). Согласно расчетам, чем больше глубина границы, тем больше годограф отраженной волны отличается от эталонного годографа для горизонтально-слоистой модели. Так, для наклонной первой границы невязки
годографов составляют не больше 10 мс, для второй — до 20 мс, для третьей — до 40-45 мс. В случае, когда наклонены обе преломляющие границы, эффект от наклона каждой границы суммируется, т. е. начинает работать эффект накопления ошибок.
Также были выполнены расчеты лучевых полей и значений времени прихода отраженных волн для вариантов модели с измененной формой преломляющих границ (рис. 10, б). Отклонение рассчитанных годографов от эталонного (для горизонтально-слоистой модели) не превышает 10 мс. При изменении формы обеих границ эффекты влияния формы каждой из границ на форму годографа отраженной волны суммируются.
Численные эксперименты по расчету лучевых полей и времен прихода отраженных волн в многослойных средах позволили сделать следующие выводы:
1) наличие скоростного градиента в слоях — наиболее существенный фактор, определяющий от-
а
?оа м
?оа м
Рис. 7. Сравнение годографов (б) отраженных волн для четырех двумерно-неоднородных моделей (а) (номер годографа соответствует номеру модели). Черный годограф соответствует исходному годографу для горизонтально-слоистой модели
550 Ф
530510-
55
«5 й> а,
РР 4904
470450430410390370350
50
клонение годографа отраженной волны от годографа в горизонтально-слоистой модели с постоянными значениями скорости в слоях;
2) при наличии в среде горизонтального градиента расхождение годографов еще существеннее и может достигать нескольких сотен миллисекунд;
3) наклоном границ можно пренебречь для углов наклона до нескольких градусов;
4)аппроксимация неровной отражающей границы плоскостью приводит к возникновению ошибок во времени прихода отраженных волн в несколько десятков миллисекунд;
5) при наклоне или изменении формы нескольких границ одновременно эффект влияния каждой из них суммируется;
6) влияние изменения формы преломляющих границ на времена пробега лучей существенно меньше, чем отражающих. С увеличением расстояния, пройденного лучами, влияние формы границ увеличивается.
^ -550 -530 -510 -490 -470 -450 -430 -410 -390 -370 350
100 150 200
Коорднната приемника, м
Заключение. На основе методов вариационного исчисления разработан алгоритм решения прямой кинематической задачи, который относится к методам лучевого трассирования с заданием начальных условий. Исследован разработанный алгоритм и выполнено его тестирование на моделях сред разной сложности: вертикально- и двумерно-неоднородных, с границами 1-го и 2-го рода. Показано, что рассматриваемый алгоритм решения прямой кинематической задачи работает корректно и дает результаты, практически не отличающиеся от точных аналитических решений, позволяя получать корректные решения прямой кинематической задачи для сложнопостроенных сред. Также изучено влияние градиента скорости и геометрии сейсмических границ на кинематику отраженных волн в многослойных средах. В результате исследований сделаны выводы о пределах применимости упрощенных моделей горизонтально-слоистых сред, которыми часто аппроксимируют сложнопостроенные неоднородные среды.
Рис. 8. Сравнение годографов
(а) отраженных волн для моделей с наклонной третьей границей. Лучевые построения
(б) для волны, отраженной от третьей границы, наклоненной на 5° против часовой стрелки. Черный годограф соответствует исходному годографу для горизонтально-слоистой модели
100 150 200
Координата приемника, м
100 150 200
Координата приемника, м
Рис. 9. Сравнение годографов отраженных волн для моделей с неровной формой отражающих границ, заданных формулами г = Н + №т0,05х (а) и г = 100 - к(х - 100)2 (б). Черный годограф соответствует исходному годографу
для горизонтально-слоистой модели
Рис.10. Сравнение годографов отраженных волн для моделей с измененной формой преломляющих границ: наклонной (а) и «волнистой» (заданной формулой г = Н + №т0,05х) (б). Черный годограф соответствует исходному годографу для горизонтально-слоистой модели
СПИСОК ЛИТЕРАТУРЫ
Галактионова А.А., Белоносов А.С. Алгоритм решения прямой кинематической задачи сейсмики в трехмерных неоднородных изотропных средах // Математические заметки СВФУ 2020. Т. 27, № 1. С. 53-68.
Ермаков А.П., Степанов П.Ю. Сейсморазведка неоднородных сред. М.: КДУ, Университетская книга, 2018. 122 с.
Облогина Т.И. Кинематическая теория сейсмических волн в неоднородных анизотропных средах // Вестн. Моск. ун-та. 1998. № 6. С. 52-59.
Облогина Т.И., Степанов П.Ю. Ефимова Е.А., Пийп В.Б. Сейсмика неоднородных сред в Московском университете // Тр. конф. «Ломоносовская школа МГУ по геофизическим методам исследования земных недр: прошлое, настоящее, будущее». М.: Изд-во Моск. ун-та, 2004. С. 51-53.
Тульчинский П.Г. Трассирование луча по конечно-разностной двумерной сейсмической модели // Компьютерная математика. 2009. № 1. С. 29-36.
Эльсгольц Л.Э. Дифференциальные уравнения и вариационное исчисление. М.: Наука, 1969. 424 с.
Cao S., Greenhalgh S. Finite-difference solution of the eikonal equation using an efficient, first-arrival, wavefront tracking scheme // Geophysics. 1994. Vol. 59, N 4. P. 632-643.
Cerveny V. Seismic Ray Theory. 1st Ed. Cambridge: Cambridge Univ. Press, 2001. 713 p.
Gjoystdal H., Iversen E., Laurain R. et al. K. Review of ray theory applications in modelling and imaging of seismic data // Stud. Geophys. Geod. 2002. Vol. 46. P. 113-164.
Klimes L., Kvasnicka M. 3D network ray tracing // Geophys. J. Int. 1994. Vol. 116. P. 726-738.
Moser T.J., Nolet G., Snieder R. Ray bending revisted // Bull. Seismol. Soc. Amer. 1992. Vol. 82, N 1. P. 259-288.
Rawlinson N., Hauser J., Sambridge M. Seismic ray tracing and wavefront tracking in laterally heterogeneous media // Advances in Geophysics. 2007. Vol. 49. P. 203-267.
Sethian J.A., Popovici A.M. 3-D traveltime computation using the fast marching method // Geophysics. 1999. Vol. 64, N 2. P. 516-523.
Virieux J., Farra V. Ray tracing in 3D complex isotropic media: An analysis of the problem // Geophysics. 1991. Vol. 56, N 12. P. 2057-2069.
Xu T., Li F., Wu Z., Wu C., Gao E., Zhou B., Zhang Z. A successive three-point perturbation method for fast ray tracing in complex 2D and 3D geological models // Tectono-physics. 2014. Vol. 627. 10 p.
Yang W. A Basical study on two-point seismic ray tracing // 2003. 8 p.
Статья поступила в редакцию 04.10.2022, одобрена после рецензирования 02.11.2022, принята к публикации 26.01.2023
