CC BY
23
ФИЗИКА ЗЕМЛИ, АТМОСФЕРЫ И ГИДРОСФЕРЫ
Численное моделирование гравитационных волн, возбуждаемых в океане низкочастотными поверхностными сейсмическими волнами,
на основе записей GPS-станций
К. А. Семенцов1а, М. А. Носов ^, С. В. Колесов1, Ю. Ву3
1 Московский государственный университет имени М. В. Ломоносова, физический факультет, кафедра физики моря и вод суши. Россия, 119991, Москва, Ленинские горы, д. 1, стр. 2.
2 Институт морской геологии и геофизики Дальневосточного отделения РАН. Россия, 693022, Южно-Сахалинск, ул. Науки, д. 1Б.
3 Университет Токио, Институт исследования землетрясений.
Япония, 113-0032, Токио, 1-1-1 Яйои, Банкио.
E-mail: a sebbest@yandex.ru, b nosov@phys.msu.ru Статья поступила 07.11.2016, подписана в печать 24.11.2016.
Описан численный эксперимент по воспроизведению процесса генерации свободных гравитационных волн в океане пробегающими по дну низкочастотными поверхностными сейсмическими волнами. Динамика движения дна восстановлена на основе реальных записей GPS-станций, выполненных во время катастрофического землетрясения Тохоку 11 марта 2011 г. Результаты численного моделирования показывают, что ключевую роль в генерации свободных гравитационных волн играют горизонтальные колебания подводных склонов.
Ключевые слова: цунами, предшественники цунами, волны Рэлея, волны Лява, поверхностные сейсмические волны, гравитационные волны в океане, потенциальная теория волн.
УДК: 551.466. PACS: 92.10.H, 92.10.hl, 91.30.Nw.
Введение
Колебания уровня моря, предшествующие приходу волн цунами и значительно уступающие им по амплитуде и периоду, неоднократно наблюдались в последние полвека [1-4]. Так, после чилийского землетрясения 24.05.1960 на записях японских мареографов перед приходом волн цунами (высоты заплеска 4-6 м) проявляются относительно слабые колебания водной поверхности с амплитудой около 0.5 м и периодом от 17 до 48 мин в зависимости от региона [1]. Записи донного датчика давления, выполненные во время землетрясения и цунами на Аляске 06.03.1988 (М¥ 7.6, амплитуда волн цунами в точке постановки датчика — 2 см), показывают, что вскоре после прохождения волн Рэ-лея возникают колебания водной поверхности с амплитудой порядка 0.5 см и периодами 100-200 с, опережающие приход лидирующей волны цунами на 40 мин [2].
Во время катастрофического землетрясения Тохоку 11.03.2011 (М¥ 9.0) предшественники цунами были зафиксированы как в Тихом океане, так и в Японском море. Анализ записей глубоководных станций ООЫБТ, расположенных в Тихом океане вблизи п-ова Кии, позволил выявить существование слабых колебаний уровня моря с амплитудами порядка 2-3 см и периодами порядка 50-150 с. Эти колебания наблюдались вскоре после прохождения по дну океана поверхностных сейсмических волн и опережали приход лидирующей волны цунами более
чем на час [3]. Российские и японские мареографы, расположенные на побережье Японского моря, также зафиксировали волны амплитудой до нескольких десятков сантиметров с периодом до 10 мин, наблюдавшиеся сразу же после землетрясения, т. е. за 1-2 ч до расчетного времени прихода лидирующей волны цунами [4].
В работах [3] и [4] высказываются различные гипотезы о механизме формирования предшественников цунами. Авторы работы [3] предполагают динамический механизм, связанный с колебаниями дна при прохождении длиннопериодных сейсмических волн. Авторы работы [4] высказываются в пользу статического механизма, предполагая, что предшественники цунами возбуждаются аналогично основной волне цунами — косейсмическими (остаточными) деформациями дна, но не в очаге, а в непосредственной близости от точек наблюдения.
Известно, что величина косейсмической деформации дна быстро убывает при удалении от источника (обратнопропорционально квадрату расстояния) [5, 6]. В этой связи статический механизм должен иметь ограниченную область действия. В противоположность статическим деформациям амплитуда поверхностных сейсмических волн убывает с расстоянием по существенно более слабому закону (обратно пропорционально корню квадратному из расстояния) [7]. Поэтому динамический механизм генерации предшественников цунами должен быть работоспособен и на значительном удалении от источника.
В случае когда поверхностные сейсмические волны распространяются по ровному горизонтальному дну океана, за формирование предшественников цунами могут быть ответственны только волны Рэлея, которые вызывают вертикальные колебания дна. Но сейсмические волны распространяются со скоростью, которая существенно превосходит скорость длинных гравитационных волн в океане. Из теории длинных волн известно, что быстро бегущее по ровному дну возмущение не способно возбуждать свободные гравитационные волны [3, 6]. В такой ситуации движение водного слоя представляет собой вынужденное возмущение, повторяющее движения дна, т. е. смещение поверхности океана оказывается эквивалентно вертикальному смещению дна в волне Рэлея. Для океана переменной глубины ситуация принципиально меняется. Теперь вклад в генерацию предшественников цунами могут обеспечивать горизонтальные движения дна, ассоциирующиеся не только с волнами Рэлея, но и с волнами Лява. Вклад горизонтальных движений отличается тем, что он может проявляться только в узких областях крутых подводных склонов. Поэтому источник гравитационных волн в водном слое уже не имеет характер бегущего возмущения, а скорее напоминает отдельные колеблющиеся участки дна.
В работе [3] мы представили приближенные теоретические оценки амплитуды предшественников цунами, основанные на динамическом механизме их генерации. Эти оценки позволили нам предположить, что наиболее вероятной причиной возникновения предшественников цунами являются горизонтальные колебания подводных склонов, вызванные прохождением по дну океана пакета длиннопери-одных сейсмических волн. Проверить это предположение можно путем численного моделирования процесса генерации гравитационных волн в океане сейсмическими движениями дна.
Цели настоящей работы: 1) воспроизвести процесс генерации предшественников цунами поверхностными сейсмическими волнами методом численного моделирования на основе фактических (измеренных) данных о движении земной поверхности; 2) оценить параметры предшественников цунами при катастрофическом землетрясении Тохоку 11.03.2011 по данным численного моделирования и сопоставить их с данными наблюдений.
1. Математическая модель
Предшественники цунами, описанные в работах [3, 4], — относительно короткие гравитационные волны, и применение для их описания теории длинных волн было бы не совсем корректно. В этой связи мы будем использовать негидростатическую трехмерную численную модель, построенную в рамках линейной потенциальной теории [8, 9].
Опишем физическую постановку задачи и уравнения, лежащие в основе численной модели. Рас-
смотрим слой идеальной однородной сжимаемой жидкости переменной глубины Н(х, у), ограниченный сверху свободной поверхностью, а снизу — подвижным, но непроницаемым и абсолютно жестким дном. Разделим расчетную область на глубоководную и мелководную части изобатой Н0. В глубоководной части задача о генерации и распространении волн решается в рамках линейной потенциальной теории сжимаемой жидкости. В мелководной части рассматривается уравнение линейной теории длинных волн.
В линейном приближении уравнения гидродинамики сводятся к волновому уравнению относительно потенциала скорости течения Р(х, у, г, 1) [10]:
d2F
- c2AF = 0.
(1)
Дополним это уравнение граничными условиями на поверхности (при г = 0):
д2Р дР
-312 = дг (2)
и на дне (при г = —Н(х, у)): дР
дп =(", и). (3)
Здесь п = (пх, пу, пг) — единичный вектор нормали к поверхности дна, g — ускорение свободного падения, и = (их, иу, иг) — вектор скорости сейсмических движений дна [11]. Решив поставленную задачу, найдем потенциал скорости течения Р(х, у, г, 1), после чего по известным формулам рассчитаем поле вариаций давления в любой момент времени р(х, у, г, 1) и смещение свободной поверхности £(х, у, 1):
дР
p(x, y, г, 1) = —р—,
*(х,y,t) = -g %
z=0
Сжимаемость жидкости не влияет на конечный результат моделирования, так как для любой точки Мирового океана гравитационные и акустические волны лежат в различных и непересекающихся частотных диапазонах [6].
Уравнения и граничные условия глубоководного блока (1)—(3) записываются в криволинейных (сигма-сферических) координатах. Переход к сигма-координатам отображает расчетную область сложной формы, которая ограничена снизу неровным дном, на прямоугольный параллелепипед, что существенно упрощает постановку граничного условия на дне. Кроме того, по сравнению с традиционной «z-le-veled» моделью возникает возможность проводить трехмерное моделирование как угодно близко к берегу, не увеличивая при этом число узлов сетки по вертикали в глубоководной области. Связь сигма-сферической системы координат с декартовой осуществляется согласно следующим формулам:
х = (Re - аИ(ф, в)) cos в cos ф,
у = (RE — оИ(ф, в)) cos в sin ф, z = (Re — оИ(ф, в)) sin в,
где Re — радиус Земли, ф — долгота, в — широта, о — «сигма-координата», величина которой варьируется от 0 на поверхности океана до 1 на дне.
В мелководной области решается уравнение линейной теории длинных волн, записанное также в терминах потенциала скорости течения F(x, y, z, t) [6]:
= V(gH VF) — gW,
где V — дифференциальный оператор, действующий в горизонтальной плоскости, W — скорость движения поверхности дна в вертикальном направлении. Скалярная величина W связана с вектором скорости сейсмических движений дна U = (Ux, Uy, Uz) и распределением глубин следующей формулой [6, 12]:
дИ . г, дИ
ду
W = Ux ^ + Uy— + Uz.
2. Восстановление динамики движения дна
Для моделирования процесса генерации предшественников цунами необходимо обладать данными о динамике движения дна. Во время землетрясения Тохоку 11.03.2011 сейсмические волны были зарегистрированы большим количеством ОРБ-станций на всей территории Японских островов (http://www.gps-solutions.com/data_2011_ ^оки_ея). Для региона (138-143° в.д., 34.5-
42.5° с.ш.) частота дискретизации данных составляла 1 Гц. Эти данные открывают уникальную возможность для восстановления динамики движения дна. Однако стоит отметить, что прямая «математическая» экстраполяция записей наземных ОРБ-стан-ций на прибрежную область океанического дна не вполне корректна. В связи с этим в качестве области моделирования желательно выбирать заливы, так как в этом случае ОРБ-станции частично окружают участок океана, что повышает корректность экстраполяции. Также необходимо дополнить сугубо математическую экстраполяцию физическими представлениями о распространении сейсмических волн.
Мы осуществляли моделирование в области 138-140° в.д., 34-35.5° с.ш. (см. рис. 1, а). Размеры области составляли примерно 250 х 150 км, что сопоставимо с длинами рассматриваемых сейсмических волн. В качестве входных данных использовалась реальная батиметрия (ОЕВСО) и пространственно-временной закон движения дна, восстановленный с использованием записей 79 ОРБ-станций, расположенных в указанной области (см. рис. 1, б). Примеры записей сейсмических волн ОРБ-станция-ми показаны на рис. 2.
Известно, что поверхностные гравитационные волны в океане возбуждаются только при колебаниях дна с достаточно низкими частотами, удовлетворяющими соотношению
f < fg = «J H,
(4)
Рис. 1. а — Расположение эпицентров основного землетрясения 11.03.2011 и афтершока (звездочки). Красным прямоугольником показана расчетная область численного моделирования гравитационных волн. Стрелка обозначает направление распространения поверхностных сейсмических волн в области моделирования; б — расчетная область численного моделирования генерации гравитационных волн. Треугольниками отмечены ОРБ-станции, расположенные в этой области. Подписями снабжены две станции, записи которых приводятся для примера на рис. 2. Звездочками показаны точки, для которых представлены результаты
моделирования (см. рис. 3)
Рис. 2. Смещения земной поверхности, записанные ОРБ-станциями (слева направо: восточная, северная и вертикальная компоненты соответственно). Верхний ряд: записи ОРБ-станции 1102 (138.19° в.д., 35.32° с.ш.), нижний ряд: записи ОРБ-станции 3053 (138.99° в.д., 34.75° с.ш.). Расположение указанных
станций показано на рис. 1, б
где а — числовой коэффициент, точное значение которого можно вычислить из аналитического решения задачи о генерации гравитационных волн колебаниями участка дна [6, 13]. Примем для данного коэффициента значение а « 0.366, соответствующее стократному ослаблению амплитуды волны по сравнению с амплитудой колебаний дна.
Из формулы (4) следует, что чем глубже океан, тем ниже предельная частота гравитационных волн, которые могут возникнуть при колебаниях дна. Ограничиваясь волнами, возбуждаемыми на глубинах H ^ 500 м, согласно формуле (4) получаем следующее значение пороговой частоты: fg = 0.05 Гц. Именно эта частота бралась в качестве частоты отсечки при удалении высокочастотных компонент сигналов из исходных записей ОРБ-станций.
Указанной критической частоте fg соответствует критическое значение периода T = 1/fg = 20 с. Анализ движений ОРБ-станций, вызванных прохождением длиннопериодных компонент (T >20 с) сейсмических волн, показал, что волны в рассматриваемой области можно считать плоскими и распространяющимися вдоль некоторого направления ф с определенной скоростью v (на рис. 1, а направление распространения сейсмических волн показано стрелкой). Выбрав в качестве «эталонного» профиля волны запись одной из ОРБ-станций, мы можем восстановить динамику движения дна, предварительно вычислив направление и скорость распространения плоских волн.
Для того чтобы вычислить направление и скорость распространения плоских сейсмических волн, воспользуемся следующим алгоритмом: 1) найдем
промежутки времени Д^-, на которые сдвинут сигнал, записанный I-й станцией, относительно сигнала, записанного --й станцией (записи станций имеют сходный вид, следовательно, ДЦ можно найти, используя кросс-корреляционный анализ записей I-й и ] -й станций); 2) рассчитаем проекции расстояния между I-й и - -й станциями Дг- (ф) на произвольное направление ф; 3) для каждой станции рассчитаем среднюю «мнимую» (зависящую от выбранного направления) скорость распространения волны по формуле vi(ф) = ДД—)) N; 4) найдем
такое значение ф, при котором ^ VI (ф)2 принимает
i
минимальное значение. Это значение и будет соответствовать направлению распространения плоской волны; 5) скорость распространения плоской волны рассчитаем, используя полученное значение ф и формулу из пункта 3. Для I возьмем значение, соответствующее ранее выбранной ОРБ-стан-ции с «эталонным» профилем волны.
Применяя описанный выше алгоритм для каждой из трех компонент записей ОРБ-станций, мы рассчитали значения скоростей распространения вертикальной и горизонтальной компонент движения земной поверхности: vv = 3.14 км/с и vH = 3.61 км/с соответственно. Известно, что в отличие от волн Лява, при прохождении волны Рэлея частицы земной поверхности колеблются как в горизонтальной, так и в вертикальной плоскостях. При этом амплитуда горизонтальных колебаний, вызванных прохождением волны Рэлея, меньше, чем амплитуда вертикальных колебаний (примерно в 2 раза
в зависимости от типа грунта) [14]. Кроме того, в нашем случае амплитуда вертикальных колебаний, записанных ОРБ-станциями, сама по себе уступает амплитуде записанных горизонтальных колебаний (см. рис. 2). Это позволяет нам ассоциировать найденное значение vV со скоростью распространения волн Рэлея, а vH — со скоростью распространения волн Лява.
Зная направление и скорость распространения плоской волны, а также ее профиль, мы можем найти пространственно-временной закон движения дна, т. е. смещение дна в любой момент времени в любой точке. В нашей гидродинамической модели в качестве граничного условия на дне используется условие непротекания: скорость жидкости по нормали к дну должна совпадать с проекцией скорости самого дна на то же направление (см. формулу (3)). Необходимое для моделирования векторное поле скорости движения дна легко найти из пространственно-временного закона движения дна путем дифференцирования по времени.
3. Результаты численного моделирования и их обсуждение
Рассмотрим результаты численного моделирования предшественников цунами. На рис. 3 показаны зависимости возвышения свободной поверхности от времени в точках с координатами 139.75 ° в.д., 35.15° с.ш. (слева) и с координатами 139.85 ° в.д., 34.7° с.ш. (справа). Расположение этих точек на карте показано на рис. 1, б. В обеих точках расчеты проводились для двух случаев: а) учитывались и горизонтальная, и вертикальная компоненты движений дна; б) учитывалась только горизонтальная компонента движений дна. Видно, что вынужденные колебания водной поверхности, возникающие непосредственно над пробегающей по дну сейсмической волной (на рис. 3 — ярко выраженные всплески амплитудой около 15 см, наблюдающиеся около
5:50 иТС), воспроизводятся только при учете вертикальной компоненты движений дна. В то же время свободные гравитационные волны (на графиках — колебания с амплитудой 3-5 см, наблюдающиеся вскоре после вынужденных колебаний) возбуждаются почти исключительно горизонтальными компонентами движений дна, а вертикальная компонента не вносит существенного вклада в их генерацию. Это подтверждает высказанное в работе [3] предположение о том, что основным механизмом генерации указанных гравитационных волн (предшественников цунами) являются горизонтальные колебания подводных склонов.
Отметим, что полученные в результате численного моделирования параметры предшественников цунами (амплитуда — 3-5 см, период — 50-100 с) в целом совпадают с величинами, измеренными станциями ЭОЫЕТ [3]. Тот факт, что амплитуда синтезированных гравитационных волн на 1-2 см превосходит амплитуду наблюдаемых, по-видимому, связан с тем, что область моделирования расположена существенно ближе к эпицентру землетрясения, чем область постановки станций ЭОЫЕТ. А так как амплитуда поверхностных сейсмических волн убывает по мере удаления от эпицентра, то убывает и амплитуда гравитационных волн в океане, порожденных этими сейсмическими волнами.
Заключение
На основе записей ОРБ-станций, выполненных во время катастрофического японского землетрясения 11.03.2011, восстановлена динамика движения океанического дна в прибрежной зоне в области с координатами 138-140° в.д., 34-35.5° с.ш. Восстановление динамики движения дна производилось в рамках приближения плоской волны. На основе восстановленной динамики движения дна было выполнено численное моделирование колебаний водно-
Рис. 3. Зависимости возвышения свободной поверхности от времени в точках с координатами 138.25° в.д., 34.55° с.ш. (слева) и с координатами 139.85° в.д., 34.70° с.ш. (справа). Расположение данных точек показано звездочками на рис. 1, б. Красным цветом изображены результаты расчетов с учетом и вертикальной и горизонтальной компонент движения дна, синим — результаты расчетов с учетом только горизонтальной
компоненты
го слоя. Результаты численного моделирования воспроизводят процесс формирования свободных гравитационных волн в океане (предшественников цунами), возникающих после прохождения по дну поверхностных сейсмических волн и опережающих приход волны цунами более чем на час. Параметры воспроизведенных гравитационных волн по порядку величины совпадают с наблюдаемыми значениями. Также результаты численного моделирования подтверждают предположение о том, что ключевую роль в генерации указанных гравитационных волн играют горизонтальные колебания подводных склонов.
Авторы признательны компаниям GSI (Geospatial Information Authority, Japan), NGDS (Nippon GPS Data Service, Japan), Hitz (Hitachi Zosen Co., Japan), GPSS (GPS Solutions, Boulder, Co, USA) и VERIPOS за предоставленные данные.
Работа выполнена при финансовой поддержке РФФИ (гранты 16-35-00232, 16-55-50018).
Список литературы
1. Nakamura K., Watanabe И. Report on the Chilean Tsunami / Ed. by Field Investigation Committee for the Chilean Tsunami. 1961. P. 82.
2. Milburn H.B., Nakamura A.I., Gonzalez F.I. // Proc. of the Oceans 96 MTS/IEEE Conference, 23-26 September 1996, Fort Lauderdale, FL. P. 390.
3. Носов М.А., Семенцов К.А., Колесов С.В. и др. // ДАН. 2015. 461, № 5. C. 593. (Nosov M.A., Semen-tsov K.A., Kolesov S.V. et al. // Doklady. Earth Sciences. 2015, 461, N 2. P. 408.)
4. Murotani S., Iwai M., Satake K. et. al. // Pure and Applied Geophysics. 2015. 172. P. 683.
5. Okada, Y. // J. Phys. Earth. 1995. 43. P. 697.
6. Levin B.W., Nosov M.A. Physics of Tsunamis. Second Edition. Springer, 2016.
7. Аки K., Ричардс П. Количественная сейсмология. Теория и методы. T. 1. М.: Мир, 1983. (Aki K., Richards P. Quantitative Seismology. Theory and Methods. Vol. 1. San Francisco: W.H. Freeman and Company, 1980.)
8. Nosov M.A., Kolesov S.V. // NHESS. 2007. 7. P. 243.
9. Колесов С.В., Носов М.А. // Ученые записки физического ф-та Московского ун-та. 2016. № 3. С. 163904.
10. Ландау Л.Д., Лившиц Е.М. Теоретическая физика. T. 4. Гидродинамика. М.: Наука, 1988.
11. Nosov M.A., Moshenceva A.V., Kolesov S.V. // Pure Appl. Geophys. 2013. 170, N 9-10. P. 1647.
12. Nosov M.A., Bolshakova A.V., Kolesov S.V. // Pure Appl. Geophys. 2014. 171. P. 3515.
13. Носов М.А. // Вестн. Моск. ун-та. Физ. Астрон. 1992. № 1. C. 109.
14. Love A.E.H. Some Problems of Geodynamics. Cambridge University Press, 1911.
Numerical simulation of gravity waves excited in the ocean by low-frequency surface seismic waves based on GPS recordings
K.A. Sementsov1a, M.A. Nosov12b, S.V. Kolesov1, Y. Wu3
1 Department of Marine and Inland Water Physics, Faculty of Physics, Lomonosov Moscow State University. Moscow 119991, Russia.
2 Institute of Marine Geology and Geophysics, Far Eastern Branch of Russian Academy of Sciences. Yuzhno-Sakhalinsk 693022, Russia
4 Earthquake Research Institute, The University of Tokyo. 1-1-1 Yayoi, Bunkyo-ku, Tokyo, 113-0032 Japan. E-mail: a sebbest@yandex.ru, b nosov@phys.msu.ru.
A numerical experiment for reproducing the generation of free gravity waves in the ocean by lowfrequency surface seismic waves passing across the bottom is described. The dynamics of the bottom movement is reconstructed based on the real GPS data recorded during the disastrous Tohoku earthquake of March 11, 2011. Results of the numerical simulation show that horizontal movements of underwater slopes play a key role in the generation of free gravity waves.
Keywords: tsunami, tsunami forerunners, Rayleigh waves, Love waves, surface seismic waves, gravity waves in the ocean, potential wave theory. PACS: 92.10.H, 92.10.hl, 91.30.Nw.
Received 7 November 2016.
English version: Moscow University Physics Bulletin. 2017. 72, No. 6. Pp. 614-619.
Сведения об авторах
1. Семенцов Кирилл Александрович — канд. физ.-мат. наук, ст. лаборант; тел.: (495) 939-36-98, e-mail: sebbest@yandex.ru.
2. Носов Михаил Александрович — доктор физ.-мат. наук, профессор; тел.: (495) 939-36-98, e-mail: nosov@phys.msu.ru.
3. Колесов Сергей Владимирович — канд. физ.-мат. наук, мл. науч. сотрудник; тел.: (495) 939-36-98, e-mail: kolesov@ocean.phys.msu.ru.
4. Ву Юфей — аспирант; e-mail: wuyifei@eri.u-tokyo.ac.jp.
