hsm2 (Расчет сейсмограмм для плоской однородной земли)

РАСЧЕТ СЕЙСМОГРАММ ДЛЯПЛОСКОЙ ОДНОРОДНОЙ ЗЕМЛИД. В. ПерегудовАннотация. Рассматривается проблема вычисления сейсмограмм в однородном полупространстве для источника с произвольным тензором сейсмического момента. Результат выражается через эллиптические интегралы. Приведены примеры расчетасейсмограмм для разных источников по разным направлениям.ВведениеПредлагаемая статья является продолжением работы [1], в которой детальнорассмотрено решение двумерной задачи Лэмба.

Там же приведен обзор литературы [2]—[10] (не только по двумерным, но и по трехмерным задачам для однородного полупространства), который мы здесь не будем повторять. Настоящая работа“ближе” к реальности тем, что исследуется трехмерный случай, к тому же источник выбирается вполне физическим и соответствует современным представлениямо механизме тектонического землетрясения.

Некоторые сведения о сейсмическомисточнике приведены в приложении к статье.Отметим сразу же особенности предлагаемой работы по сравнению с работами [2]—[10]. Во-первых, в предлагаемой работе рассматривается источник, расположенный внутри полупространства, тогда как в предыдущих работах рассматривался источник в виде внешней силы, приложенной к поверхности. Во-вторых,в предлагаемой работе рассматривается источник-тензор момента, чего ранее неделалось.

В-третьих, в предлагаемой работе подробно исследованы особенностиэталонных интегралов (фактически — функции Грина), соответствующие различным волнам. Это исследование основано на специальной замене переменных приинтегрировании по горизонтальной медленности, при этом оказывается, что всекомпоненты функции Грина для заглубленного источника (а не только компоненты функции Грина для источника на поверхности, как показано в [9]) выражаютсячерез эллиптические интегралы.Постановка задачиПростейшей моделью Земли является однородное изотропное упругое полупространство. Направим ось z декартовой системы координат внутрь полупространства, а оси x и y — параллельно границе (см.

рис. 1). Уравнения теории упругостиимеют видρüi = cijkl ∂j ∂l uk + fi ,Typeset by AMS-TEX12где ρ — плотность, ui — вектор смещения, cijkl = λδij δkl +µ(δik δjl +δil δjk ) — тензорГука. Плотность внешних сил fi выберем в виде сейсмического источника-тензорамомента∂fi (t, r) = −Mik (t)δ(r − r ),∂rkгде r — место расположения источника. (Более детальное обсуждение сейсмического источника смотри в приложении “Замечание о сейсмическом источнике”.)До землетрясения (при t < 0) тензор сейсмического момента равен нулю, а после землетрясения — некоторой постоянной величине, так что изменение Mik (t)происходит в течение короткого промежутка времени вблизи t = 0 (см. рис.

2).Для упрощения формул выберем начало координат в источнике, глубину которого обозначим через h. Границу полупространства будем считать свободной отнапряженийσiz (t, x, y, −h) = 0,где σij = cijkl ∂k ul . Предполагается также, что в начальный момент среда находилась в равновесии и покоилась ui (0, r) = 0, u̇i (0, r) = 0. Требуется вычислитьсейсмограмму, то есть зависимость от времени смещения в заданной (произвольной) точке поверхности.BMhxyAtzРис. 1. К постановке задачиРис. 2.

Характерная зависимость тензорасейсмического момента от времениСоберем все предположения вместе, но сперва введем скорости продольных и2= µ/ρ, а также их отношение γ = v /v⊥ ипоперечных волн v2 = (λ + 2µ)/ρ, v⊥новое время τ = v t. Система уравнений, начальных и граничных условий приметвидγ 2 üi = (γ 2 − 1)∂i ∂j uj + ∂j2 ui −(1)ui (0, r) = 0,Mik (τ )∂k δ(r),µu̇i (0, r) = 0,σiz (τ, x, y, −h) = 0.Сделаем несколько замечаний по постановке задачи. В силу линейности уравнений решение имеет следующий общий вид τds uijk (τ − s, x)Mjk (s).(2)ui (τ, x) =0Однако методом Каньяра, который мы применим для построения решения, естественно определяются не коэффициенты uijk , а некоторые функции Γijk , такие что3Γ̈ijk = uijk . Удобно выполнить в формуле (2) двукратное интегрирование по частям, учитывая, что внеинтегральные члены выпадают (на нижнем пределе этообеспечивается тензором сейсмического момента, а на верхнем — коэффициентамиuijk , которые соответствуют запаздывающему решению) τds Γijk (τ − s, x)M̈jk (s).(3)ui (τ, x) =0Оказывается, что даже коэффициенты Γijk (не говоря уж об uijk ) имеет довольносильные особенности (самая сильная, порядка ln τ — на фронте S-волны, следующей за головной волной, более подробное обсуждение смотри ниже в разделе“Особенности интегралов”).

Поскольку решение будет получено в виде интеграла, который предполагается вычислять численно, нужно по возможности избегатьбесконечностей. Этого можно добиться, еще раз проинтегрировав в формуле (3) почастям τ...ds Gijk (τ − s, x)M jk (s).ui (τ, x) =0Фукнции Gijk уже непрерывны, однако нужно накладывать дополнительное условие Gijk (τ, x) = 0 при τ < 0 для исчезновения внеинтегральных членов.В силу симметрии среды азимутальная зависимость смещения имеет очень простой вид и определяется трансформационными свойствами тензора сейсмическогомомента.

Поэтому фактически вычисление сейсмограмм для заданного горизонтального расстояния до источника и произвольного азимута сводится к вычислению коэффициентов разложения смещения по компонентам тензора момента длякакого-то одного (удобного) азимута.Нетрудно определить количество независимых коэффициентов. Функции Gijk ,очевидно, представляют собой тензор третьего ранга. Имеется всего два вектора,из которых его можно построить: это вектор нормали к границе n и вектор r,направленный из источника в точку наблюдения.

В нашем распоряжении имеетсятакже единичный тензор δij . В силу симметрии тензора сейсмического моментатензор Gijk можно считать симметричным по двум последним индексам. Оказывается, что можно составить 10 комбинаций с указанными свойствами, это ri δjk ,ni δjk , rj δik +rk δij , nj δik +nk δij , ri rj rk , ni nj nk , ni rj rk , ri (nj rk +nk rj ), ni (nj rk +nk rj ),ri nj nk . Кроме того, имеется остаточная симметрия по отношению к вращениювокруг оси z.

Эта симметрия дает возможность сократить число независимых коэффициентов до девяти. Действительно, свертка Gijj есть линейная комбинациявекторов n и r. Выбрав систему координат так, чтобы точка наблюдения лежалана оси y, получим Gxjj = 0, то есть линейное соотношение между коэффициентами.В дальнейшем будем вычислять коэффициенты именно для случая, когда точканаблюдения расположена на оси y.

Если мы хотим вычислить сейсмограмму не наоси y, а по какому-то другому азимуту, то вычисления можно проделать по тойже формуле (2), нужно только преобразовать коэффициенты uijk . Именно, введемсистему координат x y z, такую что ось y лежит по заданному азимуту. Тогдапреобразования вектора смещения и тензора сейсмического момента имеют вид(см. рис. 3)sin ϕcos ϕ 0, O =  − cos ϕ sin ϕ 0  .(4)ui = Oij uj , Mkl = Okn Olp Mnp0014yyeϕerBϕxxРис. 3. Преобразование координатИспользуя ортогональность матрицы O, найдем, что смещения по выбранномуазимуту вычисляются по формуле (2) с заменой uijk наuijk = Oik Ojl Okn ukln .(5)Более наглядными являются компоненты вектора смещения в цилиндрическихкоординатах, поэтому практически удобнее пользоваться не формулой (5), а преобразовывать тензор момента по формуле (4) и вычислять ui по формулам (2) споследующим отождествлением ux → −uϕ , uy → ur .Решение задачи в образахДля решения задачи (1) применим метод интегральных преобразований.

Сделаем преобразование Лапласа по времени и преобразования Фурье по x и yUi (ω, kx , ky , z) =+∞−i(kx x+ky y)dx dy e−∞0∞dτ e−ωτ ui (τ, r).Для образов получим систему обыкновенных дифференциальных уравнений. Источник входит только в точке z = 0, поэтому его наличие мы учтем, вводя условияскачка. При z = 0 получаем систему однородных уравнений, которая в явнойпокомпонентной записи имеет видγ 2 ω 2 Ux = (γ 2 − 1)ikx (ikx Ux + iky Uy + ∂z Uz ) + (−k 2 + ∂z2 )Ux ,γ 2 ω 2 Uy = (γ 2 − 1)iky (ikx Ux + iky Uy + ∂z Uz ) + (−k 2 + ∂z2 )Uy ,γ 2 ω 2 Uz = (γ 2 − 1)∂z (ikx Ux + iky Uy + ∂z Uz ) + (−k 2 + ∂z2 )Uz .Подстановка U ∼ e−αz приводит к системе алгебраических уравнений(6) 2 2γ (ω + kx2 ) + ky2 − α2(γ 2 − 1)kx kyi(γ 2 − 1)kx αUx(γ 2 − 1)kx kyγ 2 (ω 2 + ky2 ) + kx2 − α2i(γ 2 − 1)ky α   Uy  = 0.Uzi(γ 2 − 1)kx αi(γ 2 − 1)ky αγ 2 ω 2 + k 2 − γ 2 α25Как известно, линейная однородная система алгебраических уравнений имеет нетривиальное решение, если определитель ее матрицы равен нулю∆ = γ 2 (α2 − ω 2 γ 2 − k 2 )2 (α2 − ω 2 − k 2 ) = 0.Таким образом, нетривиальное решение существует приα1 =ω 2 γ 2 + k2 ,α2 =ω 2 + k2 ,α3 = −α1 ,α4 = −α2 .причем корни α1 и α3 двукратные.

Подставляя найденные корни в уравнение (6),можно найти его решения, причем корням α1 и α3 соответствуют по два линейнонезависимых решения. Считая, что α1 и α2 имеют положительные вещественныечасти (более подробное обсуждение смотри в [1]), видим, что решения с α1 и α2убывают при z → +∞, а решения с α3 и α4 — при z → −∞. Общее решение,убывающее при z → +∞ имеет видky−iα1 kxkxUx Uy  = a+  −kx  e−zα1 + b+  −iα1 ky  e−zα1 + c+  ky  e−zα2 ,Uz0k2iα2(7)причем линейно независимые решения для α1 выбраны так, что первое и второеслагаемые представляют собой вклады SH- и SV -волн, третье же представляетвклад P -волн.