добавить поле дислокации, расположенной на расстоянии а над свободной поверхностью (при у = 2а). Тогда для интересующей нас сдвиговой дислокации на свободной поверхности (у = а) получим:
е у=а = (Дю/ 2п)[ у/ (X2 + у 2) +
+ (2a - y)/((2a- y)2 + x2)]y=a = = Д rna/[rc (a2 + y2)].
(1)
Первый член описывает деформации, вызванные самой дислокацией, а второй учитывает влияние свободной поверхности (изображения дислокации). Смещения на поверхности получаются интегрированием по х (нужно помнить о скачке смещений при х = 0):
J е dy = (Дю/ 2)[1 - (2/ n)arctg( x/a)].
(2)
Рассмотрим модель смещений во время землетрясений при образовании разломов со смещением по протиранию (модель 2, рис. 2, ). Попробуем определить деформации, возникающие в результате смещения на разломе. Использование модели со смещением по плоскости вытекает из простого соображения. Измерения методом GPS проводились через несколько месяцев после основного события и главной серии афтершоков, т.е. это интегральный эффект смещения для большой серии афтершоков по линии разлома (нодальной плоскости, 130°N). В данном случае не следует рассматривать отдельно смещения при каждом событии.
В рамках теории упругости рассмотрим полупространство y > 0, в котором первоначально имеется однородное сдвиговое напряжение axz,0 (рис. 2, ). Начало осей координат лежит на поверхности полупространства, будущая трещина простирается по оси z. Далее в плоскости x = 0 возникает двумерная трещина, идущая от поверхности до глубины y = a. При появлении трещины возникают смещения. Возникающее поле деформаций может служить моделью деформации, вызванной землетрясением на разломе со смещением по простиранию. Предположим, что имеется только одна ненулевая компонента rnz смещения в направлении оси z. После образования трещины появляются две ненулевые компоненты напряжения: axz (azx) и ayz (azy). Пренебрегаем напряжением трения на разломе после землетрясения и положим nxz= 0 при x = 0, 0 < y < a. Поля смещений и напряжений в данной двумерной задаче не зависят от z. Изменение поля смещений при землетрясении происходят быстро, т.е. можно использовать реологию Гука.
Для определения напряжений, смещений и деформаций в полупространстве y > 0 нужно вывести дифференциальное уравнение, выражающее равновесие сил, действующих на элемент среды. Уравнение равновесия для сил в направлении z можно записать в следующей форме:
dojdx + dyyz/dy = 0. (3)
Величины охг и оу2 можно выразить через производные смещения, где G — модуль сдвига земной коры:
0х2 = Одю2/дх, 0уг = Од<о2/ду. (4)
Подставляя выражение (4) в (3), получаем:
д 2ь^1дх2 +д 2ю2/ ду2 = 0. (5)
Таким образом, смещение ю2 удовлетворяет уравнению Лапласа.
Решение, описывающее смещение перед образованием трещины, имеет простой вид:
Ю =(до »,«>/G) х- (6)
После образования трещины смещение описывается решением уравнения Лапласа со следующими граничными условиями:
дю2/ду = 0 при у = 0, (7)
дю2/ дх = 0 при х = 0, 0 < у < а, (8)
ю2 = 0 при х = 0, а < у, (9)
ю2 ^ (ОО)х при (10)
Решение уравнения Лапласа с граничными условиями (7)—(10) находится с помощью методов теории функций комплексного переменного. Решение для напряжений в плоскости трещины получается в следующем виде:
30 при 0 < у < а,
°хг =4 Л 2 2ч1/2 ^ (11)
[-Ох2,0 у/(у - а У пРи у > а
= 2-0х2,0 у/(а2 - у1)1'2 при 0 < у < а,
^ yz =
(12)
[0 при у > а.
В вершине трещины оба напряжения обращаются в бесконечность. На поверхности у = 0 напряжения равны:
°х2 =±°х2,0 X (X2 + а2)12, (13)
О у2 = 0. (14)
При х ^ ^ получаем: ох2 = о 0.
Смещение поверхности, обусловленное образованием трещины, равно:
Дю2 = ±0х2,0[(х2 + а2)^2 - 1x1] О. (15)
Отсюда относительное смещение бортов разлома имеет вид:
Дю20 = 2а0х2,0 / О. (16)
Сопровождающее землетрясение смещение поверхности в функции расстояния от разлома можно записать как:
Дю2 = ±Дюг0[(1 + х21а1 )^2 -1хЦа]2. (17)
Используем полученные соотношения 2D-модели для определения параметров Чуйского землетрясения. Рассмотрим изменение составляющей горизонтальной компоненты смещения, параллельной нодальной плоскости 130°^ в зависимости от расстояния пункта до нодальной плоскости. Из решения (17) для центральных станций КЦКА и CHAG получаем величины для смеще-
ния на разрыве 1.8 м, для глубины 16.5 км. Глубина разрыва из экспериментальных данных с использованием 2D-мoдели при смещении на 2 м на разломе для станций вдоль разрыва составит: ULAG — 6.4 км, КиИА — 14.4 км, CHAG — 12.4 км, иКОК — 8.2 км. Таким образом, при постоянной величине скачка смещений по разрыву глубина разрыва уменьшается по краям. Для сброса напряжений при значениях модуля G от 30 до 55 ГПа, глубине от 9 до 16 км и смещении на 2 м получаем оценки от 2 до 6 МПа и среднее значение— 4 МПа. Для описания распределения глубины разрыва и определения длины разрыва использована линейная зависимость распределения hi = ALi + B, при этом северо-западная часть простирается на 71.0 км от точки главного события, юго-восточная часть — на 101.3 км, длина разрыва составит 172.3 км. В рамках этой модели разрыва глубина в центре равна 17.5 ± ±0.9км при смещении на 2 м. При рассмотрении квадратичной зависимости распределения hi = aL^ + Ь получаем для северо-западной части — 61.4 км, для юго-восточной части — 80.8 км, при этом длина разрыва составит 142.2 км, глубина в центре — 14.6 ± 0.8 км при смещении на 2 м. Последнее распределение использовалось при 3D-мoделирoвании. Моделирование смещений проводилось с помощью вычислений по аналитическим выражениям [14] при следующих параметрах: смещение — сдвиг на 2 м; площадки смещений — для глубины от 0 до 5 км длина плоскости 130 км (от -57 до +73 км), от 5 до 10 км — 97 км (от -46 до +51 км); от 10 до 15 км — 42 км (от -20 до +22 км); коэффициент Пуассона—0.25. В этом случае для краев разрыва получаем вертикальные смещения, сходные с данными эксперимента. Например, получено для пункта иКОК по осям у, x, z соответственно 246.9, 96.2, 25.6 мм и модуль 265 мм, из наблюдений — 204.0, 192.5, 28.7 мм и модуль 280 мм. Расхождения в горизонтальных движениях говорят о нелинейности линии разрыва по краям разлома. Можно предположить, что плоскость разрыва имеет форму «пропеллера», а в смещениях по краям разрыва появляется надвиговая составляющая.
Итак, по экспериментальным GPS-данным получена модель разрыва и определена величина смещения (2 м). Используя модуль упругости коры 3.3 • 1010 Па, получаем величину сейсмического момента 0.9 • 1020 Н • м. Таким образом, магнитуда землетрясения равна MW = 7.2. Оценки по сейсмологическим данным составляют (0.91.0) • 1020 Н-м для момента, 7.2-7.3 для магнитуды MW, 7.2-7.5 для MS .
5. Постсейсмические смещения в зоне Чуйского землетрясения и тектонические движения Горного Алтая
Проведение ежегодных измерений в 2004-2007 гг. позволяет разделить тектоническую составляющую
Горного Алтая и постсейсмические движения в эпи-центральной зоне. Исключив пункты эпицентральной зоны Чуйского землетрясения, было получено поле смещений пунктов Алтайской сети с 2000 по 2006 гг. Вычисления проведены относительно пункта Ануй, расположенного на гранитном батолите на равнине южнее г. Бийск. Получено достаточно однородное смещение западной части Горного Алтая на северо-запад со скоростью около 2 мм в год. Далее к северо-западу на равнине расположена зона тектонического торошения, возможно связанная с отмеченной выше тектонической нагрузкой (рис. 3). Результаты определения тектонического эффекта подтверждаются, например, данными станции Усть-Кан, где существуют два пункта наблюдений, заложенные в 2000 и 2004 гг. Данные по обоим пунктам очень хорошо совпадают. Также на этом пункте с 2000 по 2006 г. практически отсутствует изменение высоты пункта (± 1 мм) и зарегистрировано стабильное абсолютное значение силы тяжести (± 1 мкГал).
Постсейсмические процессы в эпицентральной зоне изучались в эпоху 2004-2007 гг. Из смещений, полученных в геоцентрической системе координат, вычитались модельные движения Евроазиатской плиты, вычисленные по модели AR-IR-2006 [7]. Поле постсейсмических смещений представлено на рис. 4. Как видим, правостороннее движение повторяет косейсмический скачок, при этом скорости смещений уменьшились на два порядка. Осредненная скорость горизонтальных смещений для пунктов, расположенных к северо-востоку от линии разрыва, составила 5-7 мм в год. Рассмотрим существующие модели постсейсмических горизонтальных движений для землетрясений сдвигового типа с вертикальным разрывом [15-19].
6. Постсейсмические движения и вязкоупругая модель явления
Косейсмические и постсейсмические распределения смещений точек на поверхности отличаются введением в решение времени, т.е. от решения с упругим полупространством переходим к двухслойной модели (упругий и вязкоупругий слои [16]). Рассматривается решение при t = 0 и далее вводится время вплоть до t =
Рассмотрим двухмерный случай. Пусть имеется упругий слой с модулем сдвига ц1, объемным модулем K1 и толщиной H, лежащий на упругом полупространстве с модулем сдвига и объемным модулем K1. Поверхность у3 = 0 является свободной. Рассматриваются горизонтальный сдвиг по вертикальной границе и надвиг по границе, расположенной под углом к горизонтальной плоскости. Для них находится упругое решение.
Рассмотрим два типа разломов и смещений: разломы с горизонтальным и вертикальным смещением по плоскости разрыва. В первом случае сдвиг на разломе моделируется введением на глубине у3 = ^ сдвиговой
55° N
53° N
51°Ы
49° N
Рис. 3. Горизонтальные скорости тектонических смещений Горного Алтая ( ), исключена эпицентральная зона Чуйского землетрясения. Карта рельефа земной поверхности и направление смещений Горного Алтая ( )
дислокации со смещением по разлому А и. Плоскость, по которой происходит смещение, — это поверхность Ух = 0. Деформация, связанная с этой дислокацией, пол-
ностью описывается отдельной компонентой смещения и1 (и двумя сдвиговыми напряжениями ст12 и ст,2) и на свободной поверхности у3 = 0 смещение
50°48'
50°24'
50°00'
49°36'
86°12'
Ч
86°36'
87°00'
87-24'
87°48'
88°12'
88°36'
89°00'
89°24'
Рис. 4. Косейсмические смещения пунктов в эпицентральной зоне, в правом верхнем углу показано направление тектонических смещений в северо-западной части Горного Алтая
и2 = &и/п {ап^ (у1 /£>) + Х[ап^ (ух!(1пН + £>)) -
-агс1ё(Л/(2иЯ-/)))]}. (18)
Во втором случае моделируется надвиг по разлому с краевой дислокацией, расположенной на глубине Б с постоянным смещением А и. Плоскость разрыва простирается вверх от линии дислокации под углом ср в положительном направлении по оси 3. Здесь деформация описывается смещениями (/, и (/, и напряжениями ст33, стп и ст31. Вертикальное смещение на свободной поверхности ( з — 0) приблизительно описывается соотношением:
и3 = Аи {(1 - ц2)/(1 + ц2 /И! )] + + [/2 (82 - ■Кх \12/^ )/(82 + ц2 /^)] + /з}, (19) где ^ — геометрический фактор и 8г- = (ЪК1 +7цг)х х(3 К{ + М-,)1. Скачок в смещении А и для обоих разломов происходит в направлении, ортогональном плоскости разлома. Если }Х1 = (х2, и преобразуется в выражении для края или винтовой дислокации в упругом полупространстве. Когда и2 = 0, £/описывает деформацию поверхности упругого слоя с внутренней винтовой или краевой дислокацией.
Предположим, что нижнее полупространство вязко-упругое. Используем соответствующие подходы, чтобы получить медленную, зависимую от времени деформацию как отклик на мгновенную деформацию, вызванную смещением (дислокацией) в момент г = 0. Рассмотрим вязкое поведение материала, вызванное сдвиговым процессом, таким как крип, или релаксацией вязкого флюида в зонах плавления. В этом случае модуль сдвига и2 необходимо заменить на соответствующий оператор Лапласа ,и для стандартного вязкоупругого твердого тела:
= 1/2 [(Ь0 + Ь^)/(а0 + о,*)], (20)
где а0, ах, Ь0 и Ъх — постоянные коэффициенты в генерализированном соотношении между сдвиговым напряжением ст у и сдвиговой деформацией е^ (г Ф у):
/ /
20. Savage J.C., Prescott W.H. Asthenosphere readjustment and the earthquake cycle // J. Geophys. Res. - 1978. - V. 83. - No. B7. -P. 3369-3376.
21. Turcott D.L., Schubert G. Geodynamics. Applications of Continuum Physics to Geological Problems. - New York: John Wiley & Sons, 1982. - 730 p.
22. Calais E., Vergnolle M., Déverchère J., San'kov V., Lukhnev A., Amarjargal S. Are post-seismic effects of the M = 8.4 Bolnay earthquake (1905 July 23) still influencing GPS velocities in the Mongolia-
Baikal area? // Geophys. J. Int. - 2002. - V. 149. - No. 1. - P. 157168.
23. Segall P. Integrating geologic and geodetic estimates of slip rate on the San Andreas fault system // Int. Geol. Rev. - 2002. - V. 44. -P. 62-82.
24. Johnson K.M., Segall P. Viscoelastic earthquake cycle models with deep stress-driven creep along the San Andreas fault system // J. Geophys. Res. - 2005. - V. 109. - P. B10403.
Поступила в редакцию 02.10.2008 г.
С
Тимофеев Владимир Юрьевич, д.ф.-м.н., с.н.с., заведующий лабораторией ИНГГ СО РАН, TimofeevVY@ipgg.nsc.ru. Ардюков Дмитрий Геннадьевич, научный сотрудник ИНГГ СО РАН, ArdyukovDG@ipgg.nsc.ru Бойко Елена Валерьевна, инженер ИНГГ СО РАН, BoykoEV@ipgg.nsc.ru