Ю. Б. Ашмарина, Т. Б. Яновская
ФУНКЦИЯ НАПРАВЛЕННОСТИ ПРИПОВЕРХНОСТНОГО ИСТОЧНИКА ТИПА ДВОЙНОЙ ПАРЫ СИЛ
Введение. Исследования, посвящённые определению механизма очага землетрясения, основываются на представлении точечной сдвиговой дислокации системой двойной пары сил. Ориентация этих сил и моментов определяет основные геометрические параметры разрыва в очаге: падение и простирание плоскости разрыва и направление подвижки. Вывод, что именно такая система сил оказывается эквивалентна элементарной дислокации сдвига, основан на теореме представления, в которой сейсмическое волновое поле в рассматриваемом объёме среды, ограниченном некоторой поверхностью, выражается через смещения и напряжения на этой поверхности и тензорную функцию Грина для упругой среды [1]. При этом предполагается, что «объём» ограничен поверхностями разрыва и поверхностью на бесконечности. Интеграл по поверхности, отнесённой на бесконечность, равен нулю, а тензорная функция Грина принимается такой же, как для неограниченного пространства.
Общее выражение для n-й компоненты смещения
^п(^, t) JJ j Cijpq * ~7ÿ~Gnpd^i (1)
s q
где Gnp - компоненты тензора Грина, S - поверхность разрыва, v - нормаль к поверхности разрыва, ^ - направление подвижки, Oijpq - тензор упругих постоянных, [м] - скачок смещения на разрыве.
Подстановка в эту формулу части функции Грина, соответствующей продольной волне в однородном изотропном пространстве на больших расстояниях от источника r = || — x|,
Gnp(4,x,t)= УпУр2 b(t-r/a), (2)
4npa2r
где Yi = Xi/r, даёт поле продольной волны в дальней зоне в виде
wn(x,t) = qrVjCijpq j j [ùi(t - r/a)]dE, (3)
s
где [iii(t — r/a)] - производная скачка смещений на S. При малых размерах площадки
разрыва по сравнению с расстоянием от источника до точки наблюдения эта формула
приводится к виду
м™СМ) = lllPlq vicijpq[Mt ~ r/a)\S = ]l'{plq Мря(г ~ r/a)> (4)
4npa3r 4npa3r
где Mpq (t) = Vjcijpq [Mi(t)]5 - тензор сейсмического момента. Это выражение представляет собой смещение от двойной пары сил.
В результате применения такого подхода к построению поля продольных волн оказывается, что смещения в первой фазе волнового цуга в зависимости от направления
© Ю. Б. Ашмарина, Т. Б. Яновская, 2010
выхода волны из источника меняют знак на нодальных линиях, представляющих собой два взаимно перпендикулярных экватора на окружающей очаг сфере [1]. Положение таких линий, найденное по данным сейсмических наблюдений, даёт возможность определить механизм очага (ориентацию плоскости разрыва и направление подвижки). Знание механизма очага применяется при оценке ориентации главных осей тектонических напряжений, а также в задачах, использующих динамические характеристики волн, для введения поправки за функцию направленности излучения.
В случае же, когда источник (точечная дислокация) расположен в полупространстве вблизи свободной поверхности (близость к поверхности означает, что глубина очага меньше длины волны), необходимо либо в выражении для поля (1) учитывать вклад, вносимый смещениями на свободной поверхности, либо использовать функцию Грина, отвечающую отсутствию напряжений на свободной поверхности.
Очевидно, первый способ нереализуем, поскольку смещение на свободной поверхности неизвестно. Поэтому для определения системы сил, эквивалентной сдвиговой дислокации вблизи свободной поверхности, следует в выражениях для поля использовать функцию Грина для полупространства со свободной поверхностью. В этом случае оказывается, что, если разрыв не выходит на поверхность, система сил будет такой же, как и в безграничном пространстве, но волновое поле, вызванное этой системой сил, должно учитывать наличие свободной поверхности.
В настоящей работе выводится формула для поля продольной волны в дальней зоне от источника типа двойной пары сил, приложенных к поверхности полупространства, и исследуется различие функций направленности в случае расположения источника на поверхности полупространства и в безграничном пространстве путём численного моделирования. Для однотипных источников сравниваются значения функций направленности, характеризующих интенсивность излучения в заданном направлении, и положение нодальных линий на сфере, окружающей очаг.
Функция Грина для полупространства. Поставленная задача построения поля от источника типа двойной пары сил, приложенных к поверхности полупространства, сводится к нахождению тензорной функции Грина для случая, когда источник расположен на поверхности полупространства.
Начиная с пионерских исследований Лэмба (1904), функции Грина в различных средах изучались множеством авторов (Peceris 1955; Helmberger 1968; Johnson 1974; Luco, Apsel 1983; Yao, Harkrider 1983; Chen 1999; Zhang, Chen 2006). В нашем случае задача упрощается благодаря тому, что поле будет рассматриваться на больших (по сравнению с длиной волны) расстояниях от источника.
Функцию Грина для источника на поверхности полупространства можно было бы построить, если определить поле от источника на некоторой глубине h под поверхностью как сумму прямой и отражённых от поверхности волн [2, 7] в приближении дальней зоны и в пределе устремить к нулю. Однако такой подход не очевиден: при приближении источника к поверхности падающие на поверхность и отражённые волны нельзя рассматривать в лучевом приближении. Поэтому следует решить задачу построения поля в полупространстве, возбуждаемого приложенными к поверхности сосредоточенными простыми силами, направленными по осям X, Y,, Z.
Выведем вначале выражение для поля, вызванного вертикальной сосредоточенной силой. Поскольку источник осесимметричный, то поле представляется в виде
(Jo(kr)Wo(k, z)ez + J0(kr)Uo(k, z)er)dk,
(5)
где .1о(кт)- функции Бесселя, к - волновой вектор, а функции Шо(к, х),Ио(к, х) для однородного полупространства имеют вид
Ш0(к, х) = Аехр(-гатх/е) + В ехр(—г^ых/с),
А
11о(к,г) = — — ехр(—шсох/с) — *|ЗВ ехр(—*|Зсо,г/с), га
где
^-1
а2 1
cos dp с2
sin Фр ’ V Ь2
cos ds а . b m
-----—, - = sm§p, - = sm§s, к = —,
sin с с с
Фр, - углы падения на свободную границу Р- и Б-волн.
Коэффициенты А, В определяются из граничных условий на поверхности полупространства х = 0:
6(г)
ПГ
-2jt J
Откуда
Pa2~J~ ~ кР(а2 ~ 2b2)t/o =
dz
fdUp
pb2 [^+kWo ) =0.
Поскольку мы далее будем рассматривать только волну Р, то достаточно определить лишь коэффициент А:
А= —
g2 cos 2'fr.s sin
i2df
4пра2 g2 cos2 2ds + sin 2ds sin 2dp ’
где д = а/Ь - отношение скоростей Р- и Б-волн.
Поле продольной волны определится интегралом
(6)
Up
/(Aexp(—ikaz)Jo(kr)ez----ехр(—ikaz) J'0(kr)e.
\ га
dk,
(7)
который оценим далее методом стационарной фазы при больших kr:
U
(Z)
= ieRA
e-imR/a
R sin dp
= eR-
2g2 cos 2ds
cos dp
4jtpa2R g2 eos2 2'fr.s + sin 2'fr.s sin 2Qp
-imR/a
(8)
где eR = ez cos dp + er sin dp, R = r/ sin dp. Это выражение определяет поле волны, распространяющейся со скоростью а, на больших расстояниях от источника. Фронт этой волны сферический, амплитуда волны изменяется вдоль фронта, и в отличие от поля волны в безграничном пространстве выражение для амплитуды содержит дополнительный множитель
2g2 cos 2ds
Cz =
g2 cos2 2ds + sin 2ds sin 2dp
(9)
а
a
zz
т
xz
e
0, град.
Рис. 1. Зависимость коэффициентов Сг и Сн от угла падения продольной волны
Нетрудно показать, что этот множитель равен
g2 cos dp sin ds
cos dp (1 — Крр) +
cos ds sin dp
dp
sin ds Ksp,
(10)
где Крр, к^р - коэффициенты отражения от свободной поверхности. Это доказывает справедливость предположения, что поле продольной волны от источника на поверхности может быть на больших расстояниях представлено как суперпозиция прямой волны Р и отражённых от поверхности рР иэР в рамках лучевого представления. Аналогичным образом выводится выражение для поля волны Р, вызванного горизонтальной силой (для конкретности направленной вдоль оси X):
,(*)
ед-
cos ф sin d р
4g cos ds cos dp
4яра2Д \g2 cos2 2ds + sin2ds sin2dp
exp(—iwR/a).
(11)
Как и в случае вертикальной силы, это выражение отличается от выражения для поля в безграничном пространстве коэффициентом
Ch
4g cos ds cos dp
g2 cos
2 2ds
s + sin 2ds sin 2dp
(12)
Точно такое же выражение получается, если интерпретировать поле на большом расстоянии от источника как суперпозицию прямой и отражённых рР- и вР волн.
Поле продольной волны от двойной пары сил. Зависимость коэффициентов С я,Сн от угла падения продольной волны изображена на рис. 1. При углах меньше ~ 60° оба коэффициента равны приблизительно 2. Удвоение амплитуды волны по сравнению с безграничным пространством объясняется тем, что вся энергия волны распределяется только в сторону нижнего полупространства.
Возрастание коэффициента Ся при больших углах падения определяется увеличением вклада волны вР, а уменьшение Сн тем, что волны Р и рР гасят друг друга при приближении угла падения к 90°, а поперечная волна, излучаемая горизонтальной силой, исчезает в направлении силы.
Теперь определим смещение в продольной волне, вызванное двойной парой сил.
Пусть t - единичный вектор в направлении силы в какой-то из двух пар, а единич-
ный вектор n определяет направление плеча (или, что то же самое, направление силы в другой паре), у - единичный вектор, направленный из источника в точку наблюдения. Введём матрицу
(Ch 0 0 \
C = I 0 Ch 0 I . (13)
V 0 0 Cz)
Тогда смещение от единичной силы будет
Соответственно для смещения от пары сил («с» - от «couple», англ.) (с точностью до коэффициента) получим
Uc = (Ct, y)(y, n)Y, (15)
а смещение от двойной пары («DC» - аббревиатура от «double couple», англ.)
UDC = [(tTcY)^ n) + (ïiTcY)^ t)]Y. (16)
Обозначим
г = Cy, (17)
где Гт = (Ch sin d cos ф, Ch sin d sin ф, Cz cos d).
Амплитуда смещения соответственно
ti Y'ijj nj + nj Tj Yiti = ti(TiYj + y¿r¿ )nj = tT Sn, (18)
где элементы матрицы S имеют вид
Sij = TiYj + Yirj. (19)
Учитывая выражения для компонент вектора Y, зависящих от угла между лучом и вертикалью d и азимутального угла ф, получаем
/ 2Cz sin2 d cos2 ф Cz sin2 d sin 2ф cz+cH sjn 2-9. cos ф\
S= Сz sin2 d sin 2ф 2Сг8т2д8т2ф Сг+Сд sin 2d sin ф . (20)
\Cz+cH sin 2d cos Ф Çz±Çh_ gin 2d sin ф “2-CH cos2 d )
В случае безграничного пространства амплитуда смещения также может быть выражена формулой (18), но при этом следует принять S = YYT.
Результаты численного моделирования. Для сопоставления диаграмм излучения поверхностного и заглублённого источников мы далее нормируем амплитуду волны от поверхностного источника на коэффициент 2, поскольку, как уже выше упоминалось, излучение от поверхностного источника распределяется на полусферу, тогда как
от заглублённого той же интенсивности - на полную сферу.
Результаты численного моделирования - стереографические проекции амплитуд смещения на нижней полусфере, окружающей источник, для разных механизмов очага.
Под механизмом очага понимают ориентацию плоскости разрыва и направление подвижки, которые можно охарактеризовать в пространстве тремя углами (рис. 2).
Ориентация плоскости разрыва в пространстве определяется направлением нормали п. Вектор подвижки Ь лежит в плоскости разрыва и, как правило, направлен от движущегося края борта к неподвижному. В такой системе координат ось Х\ направлена вдоль направления простирания линии пересечения плоскости разрыва с поверхностью Земли. Ось Х3 направлена вверх, а Х2 перпендикулярна выбранным осям соответственно. Угол падения 6 определяет ориентацию плоскости разрыва относительно горизонтальной плоскости. Направление подвижки определяется углом скольжения X - углом между направлением простирания и вектором подвижки Ь. Чтобы перейти к географической системе, вводят угол простирания ф^, который в плоскости земной поверхности простирание составляет с направлением на север.
Чтобы оптимально охватить все возможные механизмы очага, угол простирания можно фиксировать: ф^ = 0о, для угла скольжения взять три значения: X = 0о, 45°, 90о, а угол падения проварьировать в промежутке 0° ^ 6 ^ 90°, для X = 0° и X = 45°, а для X = 90°, учитывая симметричность диаграмм, достаточно было варьировать в промежутке 0° ^ 6 ^ 45°. В работе 6 варьировался через каждые 5°, на рис. 3 приведён через 15°.
Для всех механизмов очага оказывается, что чем ближе нодальные линии к центру проекции (углы выхода близки к вертикали), тем меньше они различаются для поверхностного и заглублённого источников, а у чистого сброса (взброса) 6 = 45° (внизу) нодальные линии почти идентичны.
Для механизмов очага, близких к сбросу (взбросу), наибольшее различие в положении нодальных линий достигается у первой нодальной линии при наклоне плоскости разрыва 6 = 10° и составляет Д6 = 10°, а у второй при наклоне 6 = 35° и составляет Д6 = 2°.
Рис. 3. Поля смещений для приповерхностного источника (слева) и заглублённого (справа) при ф^ = 0°, X = 0°:
нодальные линии выделены белым цветом
¡-1
1
0
Рис. 4- Поля смещений для припо- Рис. 5. Поля смещений для приповерхностного источника (сле- верхностного источника (слева) и заглублённого (справа) ва) и заглублённого (справа)
при фз = 0°, X = 90°: при фз = 0°, X = 45°:
нодальные линии выделены белым цветом
Для сдвига по простиранию 6 = 90о и близких к нему механизмов нодальные линии различаются слабо (рис. 4), хотя мы видим различия в амплитудах смещения для углов 6 > 60о.
Нодальные линии могут не пересекаться, например, у таких сложных механизмов, как на рис. 3 и рис. 5 для углов падения 6 = 15о, 30о.
Что касается амплитуд смещения, то различия удобно наблюдать, построив модули отношения амплитуд смещений поверхностного и заглублённого источников 1^=0/С4=о| (рис. 6).
Для всех типов механизмов очага видны различия в амплитудах на краях стереографической проекции (т. е. для углов выхода из источника, близких к 90о), где отношение амплитуд находится в промежутке от 0 до 0,9, а также вблизи нодальных линий, где оно достигает больших значений. На практике для определения механизма
I
> 1,2
1,1
1,0
0,9
< 0,8
5 = 30
5 = 45
5 = 30
5 = 75
5 = 45
5 = 90
\ / \ / :::
I I I I
5 = 30
5 = 75
5 = 45
5 = 90
X = 15
Рис. 6. Модуль отношения амплитуд смещений поверхностного и заглублённого источника для разных механизмов очага
очага используют знаки вступления Р-волн в промежутке эпицентральных расстояний 20о ^ А ^ 90о, что соответствует примерно углам выхода из источника 16о ^ г ^ 35о. Для таких углов выхода (почти в центре стереографической проекции) отношение амплитуд смещений составляет 1 ± 0,1.
Выводы.
1. Смещение в волновом поле поверхностного источника в дальней зоне с точностью до коэффициента выражается формулой
и = е [(Су)т у + Ут (Су)]п.
2. Из-за близости источника к свободной поверхности его функция направленности Р-волн отличается от функции направленности аналогичного заглублённого источника, так как излучение от поверхностного источника распределяется на полусферу и амплитуда удваивается.
3. При обработке данных по механизмам очага для приповерхностного источника следует учитывать изменение положения нодальных линий, изменение же амплитуд смещения по сравнению с «удвоенным» заглублённым источником незначительно (10 %) в центре стереографической проекции (для углов выхода из источника 16о ^ г ^ 35о), однако в окрестности нодальных линий различие в амплитудах увеличивается. Это означает, что в задачах, основанных на интерпретации амплитуд (или спектров) продольных волн, не следует использовать данные, отвечающие выходу лучей вблизи нодальных линий.
1. Аки К., Ричардс П. Количественная сейсмология: Теория и методы. М., 1983. Т. 1.
2. Ben-Menahem A., Smith S. W., Teng T. L. (1965) A Procedure for Source Studies from Spectrums of Long-period Seismic Waves // Bull. Seismol. Soc. Am. 55, 203-235.
3. Chen X. E. Seismogram synthesis in multi-layered half-space. Part I. Theoretical formulations // Earth. Res. China, 13, 150-174.
4. Helmberger D. V. (1968) The crust-mantle transition in the Bering Sea // Bull. Seism. Soc. Am. 58, 179-214.
5. Johnson L. R. Green’s function for Lamb’s problem // Geophys. J. R. Astron. Soc. 37, 99-131.
6. Lamb H. On the propagation of tremors over the surface of an elastic solid // Phil. Trans. Roy. Soc. Lond. A 203 (1904), 1-42.
7. Langston C. A. Body wave synthesis for shallow earthquake sources: inversion for source and earth structure parameters, Dissertation // California Institute of Technology.
8. Luco J. E.,Apsel R. J. (1983) On the Green’s functions for a layered half-space. Part I // Bull. Seism. Soc. Am. 73, 909-929.
9. Pekeris C. The seismic surface pulse // Proc. Nat. acad. Sci., 41, 469-480.
10. Yao Z. X., Harkrider D.C. (1983) A generalized reflection-transmission coefficient matrix and discrete wavenumber method for synthetic seismograms // Bull. Seismol. Soc. Am. 73, 1685-1699.
11. Zhang H., Chen X. Dynamic rupture on a planar fault in three-dimensional half-space - I. Theory. // Geophys. J. Int. 164, 633-652.
Статья поступила в редакцию 16 февраля 2010 г.