Расчёт синтетических сейсмограмм методом гауссовых пучков с заданной шириной затухания Текст научной статьи по специальности «Физика»

М. А. Гейер

РАСЧЁТ СИНТЕТИЧЕСКИХ СЕЙСМОГРАММ МЕТОДОМ ГАУССОВЫХ ПУЧКОВ С ЗАДАННОЙ ШИРИНОЙ ЗАТУХАНИЯ

Введение. Расчёт сейсмического волнового поля в сложнопостроенных средах на относительно небольших расстояниях от источника возмущения необходим в задачах сейсмологии и сейсморазведки. Естественным и физически оправданным для этих целей представляется использовать лучевой метод [1]. Однако в чистом виде использование лучевого метода сталкивается с рядом трудностей из-за наличия каустик и скачков геометрического расхождения, обусловленных способом аппроксимации скоростного разреза, и для описания волнового поля приходится привлекать специальные функции. Применение конечно-разностного метода [8] позволяет избавиться от этих недостатков, однако требует больших объёмов вычислений, что приводит к значительным затратам оперативной памяти и занимает много времени.

М. М. Поповым в работе [4] был предложен метод суммирования гауссовых пучков для расчёта волновых полей в высокочастотном приближении. Данный метод независимо от того, попадает ли луч на каустику или нет, позволяет единым образом вычислять волновое поле в области с регулярным и нерегулярным поведением поля лучей и не связан с расчётами специальных функций. Поэтому предлагается использовать метод суммирования гауссовых пучков для вычисления волнового поля на небольших расстояниях от источника. Однако, несмотря на то, что он свободен от трудностей, присущих лучевому и конечно-разностному методам, до настоящего времени он не нашёл широкого распространения в задачах вычисления волновых полей в моделях реальных сред.

Суть метода, между тем, состоит в том, что в некоторой окрестности точки наблюдения М строится веер лучей. Для каждого данного луча вычисляется простейшее сосредоточенное около луча решение исходных уравнений в частных производных - так называемый гауссов пучок. Тогда для вычисления поля в точке М (в главном члене высокочастотной асимптотики) нужно просуммировать по всем лучам веера вклад от каждого гауссова пучка в этой точке. При этом каждый гауссов пучок берётся с некоторой начальной амплитудой, которая определяется источником, порождающим лучевое поле. Алгоритм вычисления волнового поля не зависит от положения точки наблюдения и не связан с необходимостью вычисления специальных функций. Поскольку гауссовы пучки при больших частотах быстро убывают при удалении от луча в ортогональной к нему плоскости, вклад в поле в точке М дают лишь те лучи веера, которые проходят не слишком далеко от М. При этом нет необходимости иметь луч, проходящий ровно через М, но точность вычисления поля, разумеется, зависит от того, сколь плотный веер лучей построен в её окрестности. Это выгодное обстоятельство избавляет нас от затратной процедуры трассирования лучей между двумя точками.

Однако построение гауссовых пучков требует задания некоторых свободных параметров, которые определяют их вид (ширину затухания). При их неудачном выборе гауссовы пучки могут иметь довольно сложную форму - например, они могут медленно затухать и/или иметь большие осцилляции. Использование таких не оптимально построенных пучков приводит не только к потере точности вычислений, но и значительно увеличивает время, необходимое для расчёта волнового поля. До сих пор выбор

© М. А. Гейер, 2010

свободных параметров основывался в большей степени на их итерационном предрасчёт-ном подборе, а выбранные параметры являлись оптимальными только для рассматриваемой среды, в которой необходимо произвести расчёты, и для данного источника колебаний [6, 11]. Форма таким способом построенных гауссовых пучков крайне неустойчива к различного рода изменениям свойств среды, в случае сложных сред гауссов пучок как бы пускается в свободное плавание, его форма уже не контролируется. Эта проблема значительно затрудняет применение метода гауссовых пучков для расчёта синтетических сейсмограмм в моделях реальных сред, что и объясняет тот факт, что метод суммирования гауссовых пучков не нашёл до сих пор широкого применения в сейсмологии и сейсморазведке.

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

Квазидвумерный гауссов пучок. В данной работе рассматривался следующий частный случай (который можно назвать квазидвумерным): параметры среды зависят от двух пространственных переменных, т. е. не зависят, например, от координаты у, а источник и точка наблюдения М находятся в одной плоскости у = 0. Впервые такой гауссов пучок был подробно рассмотрен Качаловым и Поповым в [3]. Приведём некоторые основные идеи и формулы из этой работы, необходимые нам в дальнейшем.

Уравнению упругости будет соответствовать центральное поле лучей с центром в начале координат, лучи которого г(в, 0, ф) параметризируются углами сферической системы 0, ф координат (в - длина дуги луча, отсчитываемая от источника). Лучи г(в, 0, 0), для которых ф = 0, будут оставаться в плоскости у = 0 при всех в ^ 0. В окрестности каждого луча г(в, 0, 0), который будем называть центральным, вводятся локальные координаты в,п [4, 5], где п - расстояние по нормали от точки наблюдения М до луча в плоскости у = 0.

Рассмотрим высокочастотное стационарное волновое поле. Гауссов пучок О0(в, п; ю) такого поля, распространяющийся вдоль центрального луча г(в, 0,0), или, иными словами, асимптотическое решение стационарного уравнения упругости в некоторой окрестности данного луча строится явно через соответствующие решения <7(в),р(в) и д(в),р(в) линейных систем уравнений Якоби [4-7], описывающих в первом приближении лучи, близкие к центральному. В гамильтоновой форме эти системы имеют вид:

£Ф) = У(в)р(в) а^00 = °

(1)

и

где соответственно V(в) и д2У(в)/дп2 |п=о - скорость и вторая производная от скорости по направлению п, вычисляемые на центральном луче.

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

В работе [3] показано, что в качестве решения системы (1) необходимо выбрать такое, которое описывает расходимость лучевой трубки в плоскости, ортогональной ф = 0, а именно:

/ 5

$00 = <300 = Л =

1р(*) = Р(3) = ц$.

Для решения системы (2) вводится её фундаментальная матрица

Ж (в) =

с начальными условиями в виде

Ж (0)

ді(в) ®(в)

Рі(в) Р2(в)

1 0 0 1

(3)

(4)

Введём далее функции Q(в) или Р(в), задаваемые следующими выражениями

|д(в) = еді(в) + д2(в)

(в) = ері(в) +Р2(«)

№) = ді (в) +ед2(в) \Р (в) = рі(в) + ер2(в),

или аналогично

(5)

(5*)

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

1т

Я(в) = °, (Р (в)

VQ(s)

> °.

(6)

(7)

Тогда функция, описывающая квазидвумерный гауссов пучок О0(в, п; ю), определяется формулой

^ £(в, 0,0) Ф(0; ю) (. ію Р(в) 2

С?е(в, п; со) = —. =—. = ехр *сот(в) Н-г^гп

ву,,) V 2д(.)

(8)

где і(в, 0,0) - единичный вектор, касательный к центральному лучу, р(в) и т(в) - соответственно плотность среды и эйконал, вычисленные на центральном луче, т. е. когда

п = 0. Функция Ф(6, ю) определяется конкретным типом источника. В частности, для центра расширения [3] Р(и) = V(Ь(х)Ь(у)Ь(г)) имеем

,, . 1 / со \2 1 1о((ЛэтВ / л\

«*> = 2У

Вследствие того, что волновое уравнение является линейным, полное поле О (г; ю) в точке М вычисляется суммированием вкладов от каждого гауссова пучка Об(в,п; ю) по всем лучам веера, построенного в окрестности М, а именно

П

О (г; ю) = ^ Об(в,п; ю)d6. (9)

0

Для перехода от стационарного случая к нестационарному используется преобразование Фурье

СЮ

и(г',Ь) = — 11е J /р(<х>)С(г; со) ехр(—i^x>t)d^x>, (10)

0

где /р(ю) - спектральная плотность временного воздействия /(Ь) в источнике. Интегрирование происходит только по положительной полуоси ю, так как подразумевается, что / (t) и и (г; t) вещественные, обращающиеся в нуль при Ь < 0.

Построение гауссовых пучков любой заранее заданной ширины затухания. Для использования метода гауссовых пучков необходимо прежде всего определить форму этих суммируемых пучков. Из формулы (8) функции Об(в,п; ю) следует, что вид гауссовых пучков задаётся функциями Q(в) и Р(в), представляющими линейную комбинацию (5) или (5*) фундаментальных решений системы (2). Однако пара функций Q(в) и Р(в) является многозначной, так как зависит от выбора свободного параметра е, ограниченного лишь условиями (6) и (7). Выбирая параметр е, можно значительно менять форму гауссовых пучков, а следовательно, точность и скорость вычисления волнового поля. Так, при суммировании очень узких пучков метод гауссовых пучков переходит в лучевой [1], при очень широких - в метод Чепмана-Маслова [13], а использование гауссовых пучков конечной ширины делает рассчитываемое волновое поле неразрывным [6].

В такой ситуации выбор неоптимальных параметров е может привести, например, к неоправданно узкому или, наоборот, к слишком широкому, с точки зрения физики, и/или имеющему сильные осцилляции гауссову пучку, что, безусловно, приведёт к потере точности счёта волнового поля. К сожалению, до настоящего времени не имелось надёжного и простого метода выбора этого параметра. Обычно использовалось угадывание некоторых качественных критериев выбора формы пучка, но не было формализованной удобной процедуры выразить их количественно, чтобы использовать в программах расчёта. Конечно, хорошо известно, как выбрать форму пучка в некоторых простых структурах, но для сложных сред такие угадывания невозможны. Это создавало серьёзную проблему при использовании метода гауссовых пучков для расчётов синтетических сейсмограмм в моделях реальных сред. Несколькими авторами, однако, были предложены некоторые алгоритмы выбора формы гауссовых пучков [4, 6, 9-11], но ни один из них не нашёл широкого применения. Поэтому ниже предлагается новый, простой и универсальный алгоритм построения гауссовых пучков любой заранее

заданной ширины затухания Ьо для любых сред и одновременно с минимальным для данной ширины количеством осцилляций функции Об.

Перепишем формулу (8) в виде

Из физических соображений следует, что К (в) отвечает за кривизну волнового фронта, а Ь(в) есть эффективная ширина затухания гауссова пучка [6, 11]. При вычислении интеграла в формуле (9), очевидно, удобно использовать узкие гауссовы пучки, так как это позволяет ограничиться лишь достаточно близкими к точке наблюдения М лучами. Кроме того, удобны такие пучки, у которых кривизны волнового фронта в окрестности точки М наименьшие - это уменьшает осцилляции функции Об(в,п; ю). В некоторой степени оба эти требования можно удовлетворить специфическим подбором свободного параметра е.

Пусть функции Q и Р (для простоты записи опустим в соответствующих функциях явную зависимость от в) задаются формулами в виде (5). Представим теперь е в виде

где А и В - вещественные функции. Наша задача определить их таким образом, чтобы гауссов пучок получился шириной Ьо.

Определяя отношение Р/( с помощью формул (5) и (13) и используя тот факт, что величина дір2 — рід2 есть определитель фундаментальной матрицы Ш (3), получаем мнимую часть отношения Р/(( в виде

Потребуем теперь, чтобы ширина гауссова пучка Ь(в) была равна Ь0, и тогда из (11) и (14) получаем

Выберем теперь А таким образом, чтобы последнее слагаемое в (15) обращалось в нуль:

(11)

где

(12)

е = А + ІВ,

(13)

(14)

Отсюда

В\1 Н-----------^ В (іеі \¥ + (АЧ1 + д2)2 = 0.

(15)

(Аді + <й) = 0 => А = —— (ді ф 0).

(16)

ді

Тогда уравнение (15) может быть записано в форме

Данному уравнению соответствуют два решения: В і =0 и

(17)

Из них следует выбрать решение В2, так как первое решение В\ соответствует бесконечно широким гауссовым пучкам.

Используя теперь эти выражения, преобразуем формулы для Q(в) и Р(в)/Q(в):

Отсюда функция К (в), отвечающая за кривизну волнового фронта, представляется в виде

Таким образом, определяя комплексную постоянную е в виде (13), а параметры А и В - исходя из уравнений (16) и (17), мы можем задавать гауссов пучок с любой необходимой нам величиной затухания Ь(в) = Ьо. Если на пути распространения гауссова пучка присутствуют границы раздела сред, то на каждой из данных границ он претерпевает отражения/преломления, оставаясь при этом гауссовым [4, 5]. Для дальнейшего построения отражённого/преломлённого гауссова пучка необходимо отразить/преломить соответствующее комплексное решение Ш(в) на границе по закону геометрической оптики лишь в первом приближении [5, 7], т. е. с точностью до линейных по д(в) и р(в) членов. Последующее описание отражённого/преломлённого гауссова пучка осуществляется аналогично вышеописанному начальному падающему пучку с той лишь разницей, что величины д(в), р(в) во всех формулах будут отвечать уже соответствующим отражённым/преломлённым фундаментальным решениям.

Продифференцируем теперь определитель матрицы:

Заменяя производные соответственными величинами из уравнения (2), получаем

Это выражение справедливо и для отражённых/преломлённых лучей, так как их соответствующая матрица трансформации [5, 7] является линейной. Отсюда следует, что определитель det Ш(в) остаётся постоянным вдоль всего луча и определяется только начальными данными матрицы Ш(в). В частности, для фундаментальной системы уравнений (2) с начальными условиями (4) получим

выкладки делаются аналогично, так что в конечных формулах (18) и (19) функции (в) и р1 (в) заменяются на соответствующие ®(в) и Р2(в), а знак в формуле (18) меняется на противоположный.

(18)

Р(в) Р1(в) , 2

= —ГГ + *—7

(19)

^ (в) = д'і (в)Р2(в) + ді(в)р2(в) — д2 (в)рі(в) — д2 (в)р1(в).

(20)

det 'Ш (в) = 0.

(21)

Тогда

det Ш(в) = det Ш(0) = со^ .

(22)

det Ш (в) = det Ш (0) = 1.

В том случае, если функции Q(в) и Р(в) задаются формулой (5*), то все математические

Резюмируя вышеприведённые соображения, получаем, что для построения гауссовых пучков (8) с заранее заданной шириной затухания Ьо величины Q(в) и Р(в)/Q(в) должны задаваться на основании следующих формул:

(23)

где верхний знак соответствует формуле (5) (к = 1), а нижний - формуле (5*) (к = 2).

Из соотношений (23) вытекает, что для описания гауссова пучка, построенного таким способом, достаточно знать только одно фундаментальное решение системы (2). Важно также отметить, что для системы (23) тождественно выполняются условия (6) и (7) выбора свободного параметра е, т. е. Q(в) ни при каком в не обращается в нуль, что обеспечивает отсутствие сингулярностей гауссовых пучков, а мнимая часть отношения Р(в)/Q(в) всегда больше нуля, что обеспечивает сосредоточенность гауссова пучка в окрестности центрального луча.

Согласно свойству фундаментальной системы уравнений (3), (в) и ®(в) не могут

одновременно обращаться в нуль, следовательно, система (23) всегда может быть определена. В том случае, если оба эти решения далеки от нуля, то для заданной ширины затухания Ьо можно выбрать такую пару ^(в), Р1(в) или д2(в), Р2(в), чтобы функция гауссова пучка давала наименьшее количество осцилляций в интеграле (9). Как следует из (11), это достигается путём выбора минимального из множителей В1 =0:

Таким образом, предложенная методика выбора условий для свободного параметра е позволяет строить решение уравнения упругости, состоящее из суперпозиции гауссовых пучков Об(в, п; ю) любой ширины затухания Ьо и минимальным для данной ширины количеством осцилляций. Свойства таких гауссовых пучков позволяют, во-первых, заранее ограничиться лишь достаточно близкими к точке наблюдения М лучами, лежащими в небольшой окрестности от неё, например \п\ ^ 3Ьо, так как очевидно, что основной вклад в точке М дают именно гауссовы пучки от этих лучей, в то время как вкладами лучей с \п\ > 3Ьо можно пренебречь, что позволяет уменьшить время расчёта волнового поля. Во-вторых, минимизация количества осцилляций функции Об(в,п; ю) для данной ширины Ьо повышает точность вычисления интеграла (9) для заданной плотности лучей.

Из физических соображений разумно задавать ширину гауссова пучка Ьо пропорционально длине волны - чем выше частота, тем уже пучок. Такой выбор ширины гауссова пучка хорошо согласуется и с его математической теорией [4, 6], согласно которой п ~ 0(ю-1/2). В самом простом случае ширина Ьо может быть задана равной длине волны X:

тогда необходимо учитывать только те лучи, для которых выполняется условие

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

г . 1 Р 1<» г . 1 Р2(з)

2т(в) </1(в) ’ 2т(в)д2(«)’

Ьо =

(24)

(25)

Поэтому необходимо следить, чтобы количество лучей, проходящих в данной окрестности, всегда было достаточно для корректного вычисления интеграла (9). Этого можно добиться, например, уменьшением шага выхода лучей из источника. Для расчёта численных моделей волнового поля, представленных ниже, ширина гауссова пучка Ьо задавалась постоянной и равной длине волны в точке наблюдения (25).

Алгоритм расчёта. Алгоритм вычисления волнового поля по предлагаемому методу для квазидвумерного случая естественным образом распадается на четыре этапа.

Первый этап заключается в построении веера лучей, проходящих в некоторой окрестности точки наблюдения М. Это реализуется с помощью стандартных для лучевого метода алгоритмов. Подразумевается, что все лучи лежат в плоскости у = 0. Если лучи веера ги = г (в, 6и), которые попадают в окрестность точки М, занумерованы целым индексом к = 1, 2, 3,..., где 6и - значение углов выхода лучей из источника, то далее необходимо определить координаты в и ,пи точки наблюдения М в локальной системе координат, связанной с центральным лучом г и. Координата в и определяется как точка пересечения луча ги с перпендикуляром, лежащим в плоскости у = 0 и проходящим через точку М. Координата пи находится как длина перпендикуляра из точки М в точку в и.

На втором этапе вычисляется стационарный гауссов пучок Об(в,п; ю), связанный с лучом г и. Для этого необходимо для точки в = в и знать два линейно независимых решения системы уравнений (2) или, другими словами, фундаментальную матрицу решений Ш(ви) = Ши. Это достигается численным интегрированием вдоль луча линейной системы уравнений (2). Если луч ги попадает на границы раздела сред, то через точки отражения/преломления матрица Ши продолжается в соответствии с работами [5, 7], а амплитуда отражённого/преломлённого гауссова пучка должна быть домноже-на на соответствующий коэффициент отражения/преломления. На основании равенств системы (23) находятся значения величин Q(в) и Р(в)/Q(в), необходимые для построения гауссова пучка (8) с заранее заданной шириной Ьо, равной длине волны в точке наблюдения.

Далее на основании формулы (9) вычисляется стационарное волновое поле О(М; ю) в точке М (точнее, главный член высокочастотной асимптотики О(М)) путём суммирования всех вкладов в данную точку от гауссовых пучков лучей ги веера в окрестности \пи\ < 3Ьо.

Для перехода к нестационарному случаю применяется быстрое преобразование Фурье по частоте на основе формулы (10).

Модельные расчёты. Вышеописанный метод был численно реализован в вычислительном пакете, позволяющем рассчитывать волновые поля в неоднородных упругих средах, имеющих слоистую структуру. Для подтверждения его работоспособности в примерах 1 и 2 использовались те же модели, что и в работах [2, 6], поэтому рассчитанные стационарные волновые поля сравнивались с соответствующими результатами, полученными в указанных статьях. В качестве дополнения к стационарным задачам этих примеров здесь приводятся также синтетические сейсмограммы. В примере 3 демонстрируются теоретические сейсмограммы отраженных и обменных волн, образующих многочисленные пересекающиеся друг с другом каустики.

В рассматриваемых ниже примерах точечный источник помещён в начало координат и генерирует только продольные волны. Функция ](Ь) принимается в виде

f (Ь) = |0, 1 < 0

1 at exp-|3^ 8т(ю£), Ь > 0,

2

и

Эпицентральное расстояние X, км

Рис. 1. Лучи в упругом волноводе и образуемая ими каустика (пример 1): треугольниками отмечены положения приёмников

Рис. 2. Модули амплитуды волнового поля в зависимости от эпицентрального расстояния для стационарной задачи (пример 1):

нижняя сплощная кривая соответствует частоте 10п Гц, пунктирная — 20п Гц, верхняя сплошная — 30п Гц

где а = 200 с-1, (3 = 40 с-1. Шаг выхода лучей из источника равен 0,03°. Для всех сред плотности принимаются за единицу, а скорости поперечных волны в \/3 раза меньше, чем соответствующие продольные.

Пример 1. Для расчётов был взят волновод продольных волн из работы [2]. Скорость среды V(г) в зависимости от нормального удаления г от его оси аппроксимируется функцией V(г) = (3,0 + 0,2г2) км/с. Лучи выходят из источника, располагающегося на оси волновода (г = 0), под углами в диапазоне от —80 до +80° и образуют симметричные, состоящие из двух ветвей каустики, которые пересекаются на оси волновода. Соответствующая картина лучей показана на рис. 1. Волновое поле вычисляется в окрестности первой точки возврата на оси волновода в интервале от 6 до 12 км. Синтетические сейсмограммы рассчитываются для приёмников, распределённых внутри вышеуказанного интервала с шагом 500 м.

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

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

Рис. 3. Синтетические сейсмограммы (пример 1)

характерной для него фазой. При последовательном смещении точки М вдоль оси волновода от источника к острию каустики происходит, во-первых, увеличение количества лучей, дающих вклады в точке М с большими амплитудами, а, во-вторых, сдвиг фаз вкладов пучков в точке наблюдения М по мере её приближении к фокусу лучей будет стремиться к нулю. Очевидно, что глобальный амплитудный максимум соответствует их конструктивной суперпозиции. Дальнейшее осцилляционное поведение амплитуды соответствует периодическому синфазному и противофазному суммированию пучков. Графики хорошо согласуются с результатами, полученными в работе [2].

Теоретические сейсмограммы, рассчитанные для приёмников, располагающихся в местах, отмеченных треугольниками (см. рис. 1), показаны на рис. 3. Здесь видно, что первые семь приёмников регистрируют только одну волну. Если последнюю рассматривать как интерференцию множества нестационарных гауссовых пучков, то регистрация только одной волны объясняется тем, что все гауссовы пучки, дающие заметные вклады в точку наблюдения, достигают приёмников в небольшой интервал времени, сравнимый с периодом колебания волны. В результате их интерференции формируется сосредоточенный волновой пакет, так что его можно интерпретировать как одну волну. При приближении приёмников к острию каустики амплитуда регистрируемых волн последовательно увеличивается. Приёмники же, которые расположены за этим остриём, будут регистрировать уже две волны, так как гауссовы пучки лучей, проходящих в окрестности этих приёмников, могут условно быть разделены на две группы, которые достигают приёмников в два разных временных интервала. Первая регистрируемая волна соответствует интерференции группы гауссовых пучков от симметричных параксиальных лучей, а вторая - от приосевых. Вследствие скоростных свойств среды эти параксиальные гауссовы пучки достигают приёмников, находящихся за точкой возврата, быстрей, чем приосевые гауссовы пучки, проходя при этом большее расстояние, чем последние, что вызывает более сильное затухание амплитуд параксиальных пучков.

Скорость Т(2), км/с Эпицентральное расстояние X, км

Рис. 4- Скоростная модель среды (а); лучи в данной среде и образуемая ими каустика, пересекающая дневную поверхность (б) (в примере 2): треугольниками отмечены положения приёмников

Рис. 5. Модули амплитуды волнового поля в зависимости от эпицентрального расстояния для стационарной задачи (пример 2):

нижняя сплошная кривая соответствует частоте 10п Гц, пунктирная — 20п Гц, верхняя сплошная — 30п Гц

Эпицентральное расстояние X, км

Пример 2. Рассматривается модель, аналогичная той, которая приводится в работе [6]. Она состоит из однородного слоя со свободной верхней границей, располагающейся над слоем с постоянным градиентом скорости, равным 0,2 с-1 (рис. 4а). Лучи выходят из источника в диапазоне углов от —50 до 0°. Соответствующие им волны можно разделить на три вспомогательных группы: прямые, неглубоко рефрагирован-ные и глубоко рефрагированные волны. Глубоко рефрагированные волны образуют простую каустику, пересекающую дневную поверхность (рис. 4б). Волновое поле вычисляется в окрестности этой каустики на эпицентральном расстоянии в интервале от 30 до 45 км. Синтетические сейсмограммы рассчитываются для приёмников, распределённых внутри вышеуказанного интервала с шагом 500 м. Отражённые от свободной поверхности волны в этом примере не рассматриваются.

На рис. 5 показаны кривые зависимости амплитуды стационарных (ю = 10п Гц, 20п Гц, 30п Гц) волновых полей от удаления для рефрагированных волн —50 ^—10°.

10,0

10.5 11,0

11.5 12,0

с

мя, 12,5 е р

Щ

13,0

13.5 14,0

14.5

15,0

15,5

Рис. 6. Синтетические сейсмограммы (пример 2):

по горизонтальной оси отложено эпицентральное расстояние приёмников (км), по вертикальной оси — время (с)

Как и в случае с волноводом, максимум амплитуды поля достигается в окрестности каустики, после чего она начинает осциллировать. Данные амплитудные кривые хорошо согласуются с теми, которые были получены в работе [6].

На рис. 6 показаны рассчитанные теоретические сейсмограммы. Видно, что первые несколько приёмников, расположенных до каустики, регистрируют только одну прямую волну, а приёмники, размещённые за каустикой, начинают регистрировать ещё и рефрагированные волны. В дальнейшем при увеличении эпицентрального расстояния рефрагированные волны расщепляются на две волны: приходящая первой из них будет соответствовать глубоко рефрагированной, а приходящая второй - неглубоко ре-фрагированной волнам.

Пример 3. Внесём изменения в модель из примера 2 так, что во втором слое градиент скорости будет равен 0,4 с-1. Точечный источник генерирует продольные волны в диапазоне углов от —90 до 0° градусов. Эти волны из-за наличия градиента скорости во втором слое рефрагируют в нём, образуя простую каустику, и возвращаются в первый слой, где отражаются от свободной поверхности. При этом формируются однотипные и обменные отражённые волны, которые опять рефрагируют во второй среде.

0 20 40 60

Эпицентральное расстояние X, км

Рис. 7. Лучи в диапазоне их углов выхода из источника — 75 ^—35° (пример 3): треугольниками отмечены положения приёмников

В результате картина волнового поля будет содержать множество периодически образующихся и взаимно пересекающихся простых каустик, соответствующих продольным и поперечным волнам. Для удобства визуализации лучевая картина волнового поля показана на рис. 7, 8 последовательно для лучей, выходящих из источника под углами —75 ^—35° и -35 ^—0° соответственно.

На этих рисунках отчётливо видно множество различных каустик, которые пересекаются как друг с другом, так и с линией свободной поверхности. Виртуальные приёмники, регистрирующие волны, располагаются на дневной поверхности с равным шагом 100 м в интервале от 20 до 60 км. Соответствующие им сейсмограммы (рис. 9), рассчитанные методом гауссовых пучков, представляют собой комбинацию записей волн различных типов.

Сравнивая ход лучей на рисунках и отвечающие им сейсмограммы, проведём их простейшую качественную интерпретацию. Годограф 1 соответствует прямой продольной волне Р, годографы 2 и 3 - глубоко и неглубоко рефрагированным продольным волнам Р, годографы 4 и 5 - глубоко и неглубоко рефрагированным продольным волнам РР, отражённым один раз от свободной поверхности, и, наконец, годографы 6 и 7 - глубоко и неглубоко рефрагированным обменным поперечным волнам РБ.

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

Глубина 2, км

“I 1 Г

20 40

Эпицентральное расстояние X, км

Рис. 8. Лучи в диапазоне их углов выхода из источника —35 ^ 0° (пример 3): треугольниками отмечены положения приёмников

25,0 30,0 35,0 40,0 45,0 50,0 55,0

Рис. 9. Синтетические сейсмограммы (пример 3):

по горизонтальной оси отложено эпицентральное расстояние приёмников (км), по вертикальной оси — время (с)

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

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

Автор выражает благодарность научному руководителю профессору Татьяне Борисовне Яновской за многочисленные советы при выполнении данной работы и за полезные критические замечания при подготовке данной статьи, а также Михаилу Михайловичу Попову за консультации по методу гауссовых пучков.

Литература

1. Бабич В. М. Лучевой метод вычисления интенсивности волновых фронтов // Докл. АН

СССР. 1956. Т. 110. № 3. С. 551-568.

2. Качалов А. П., Попов М. М. Применение метода суммирования гауссовых пучков для расчёта высокочастотных волновых полей // Докл. АН СССР. 1981. Т. 258. № 5. С. 1097-1100.

3. Качалов А. П., Попов М. М. Применение метода суммирования гауссовых пучков для расчёта теоретических сейсмограмм // Зап. науч. семин. ЛОМИ. 1986. Т. 156. С. 73-97.

4. Попов М. М. Новый метод расчёта волновых полей в высокочастотном приближении // Зап. науч. семин. ЛОМИ. 1981. Т. 104. С. 195-216.

5. Попов М. М. Метод суммирования гауссовых пучков в изотропной теории упругости // Изв. АН СССР. Физика Земли. 1983. № 9. С. 39-50.

6. Cerveny V., Popov M. M., Psencik I. Computation of wave fields in inhomogeneous media - Gaussian beam approach // Geophys. J. R. Astr. Soc. 1982. Vol. 70. P. 109-128.

7. Cerveny V., Popov M. M., Psencik I. Gaussian beams in two-dimensional laterally varying

layered structures // Geophys. J. R. Astr. Soc. 1984. Vol. 78. P. 65-91.

8. Kelly K. R., Ward R. W., Treitel S., Alford R. M. Synthetic seismograms: a finite-difference approach // Geophysics. 1976. Vol. 41. N 1. P. 2-27.

9. Kim W., Garmany J. Optimum Beamwidth for Gaussian beams in three dimensional problems // EOS. 1985. Vol. 66. P. 980.

10. Klimes L. Optimization of the shape of Gaussian beams of a fixed length // Studia geo-physica et geodaetica 1989. Vol. 33. P. 146-163.

11. Konopaskova J., Cerveny V. Numerical modeling of the time-harmonic seismic wave fields in simple structures by the Gaussian beam method. Part I // Studia geophysica et geodaetica. 1984. Vol. 28. P. 19-35.

12. Popov M. M, Psencik I. Ray amplitudes in inhomogeneous media with Curved Interfaces // Geofysikalni Sbornik XXIV. 1976. P. 111-129.

13. Thompson C. J., Chapman C. H. An introduction to Maslov’s asymptotic method // Geophys. J. R. Astr. Soc. 1985. Vol. 83. P. 143.

Статья поступила в редакцию 11 мая 2010 г.