东海海啸反问题预报模式研究——以日本“3·11”海啸为例

应超，于普兵，穆锦斌，周华民

（浙江省水利河口研究院，浙江杭州310020）

东海海啸反问题预报模式研究——以日本“3·11”海啸为例

应超，于普兵，穆锦斌，周华民

（浙江省水利河口研究院，浙江杭州310020）

摘要：采用COMCOT数值模式，将日本东海岸划分为18个单位震源，建立了日本东海岸海啸数据库并用非负约束的最小二乘法建立反问题预报模式。将本模式应用于日本“3·11”海啸，计算所得的海啸初始水位有10 m的抬高与3 m的降低，与前人研究结果基本一致，预报的浮标水位与实测资料符合良好。对比浙江省近岸潮位站实测海啸波高，预报值与实测值偏差较大。若采用反问题反演的震源，通过COMCOT非线性模式求解近岸水位，可以大幅提高预报精度。

关键词：COMCOT；海啸；反问题；最小二乘法

1　引言

海啸是一种具有强大破坏力的海浪，具有传播速度快、破坏力大等特点。准确及时的海啸预报是减小海啸灾害的重要举措。海啸预报中，海啸波的到达时间、波高与预报速度是需要重点关心的问题。

反问题预报方法是海啸预报中的一种重要方法，是指用预先确定好的地震源参数把地震源划分为不同的断层，用不同断层单位活动产生的水位过程建立数据库，然后回归分析记录的海啸波形确定断层滑动分布和地震震级。国内外很多学者对反问题预报模式进行过研究，并取得了丰硕的成果。Titov[1]等利用断层滑动分布通过反问题算法确定越洋海啸实时模型的初始条件，Wei[2]等采用反问题方法通过预先计算好的潮位过程的数据库确定远场海啸波形预测海啸。在Johnson[3]地震源参数基础上，Wei[2]等用反问题方法准确确定滑层分布，并准确计算夏威夷各潮位站的波高。Yamazaki[4]等用此方法很好地模拟了1944年Tonakai和1994年Kuril地震海啸通过太平洋在日本产生的海啸过程。李林燕[5]等采用反问题方法对马尼拉海沟建立海啸数据库并建立预报模式，但并未采用实测海啸资料进行验证，仅采用正反问题的计算结果进行相互验证。以往的这些研究，均未对中国东海附近海域进行研究，COMCOT模型及反问题方法在中国东海海域的适用性均需进行验证与探讨，另外反问题方法在近岸预报的适用性问题也应进行进一步分析。

北京时间2011年3月11日13时46分，日本东北部海域发生了9.0级特大地震，从而诱发了灾害性海啸，造成了巨大的人员伤亡与财产损失。此次海啸已经对我国浙东沿海造成影响，地震发生后6—9 h，浙江沈家门、大陈、坎门、石坪、石浦、健跳、福建东山及广东汕头、汕尾等潮位站先后监测到振幅为10—55 cm的海啸波[6]。因此研究日本“3·11”海啸反问题预报模式对我国东海海啸预警机制的建立与完善具有很大的借鉴意义。

2　反问题预报方法

反问题预报方法如下：

（1）设置海啸浮标。海啸浮标一般设置在地形平坦水深较大的海域，距震源有一定距离但不宜太远，从而保证浮标既能较快采集数据又不被海啸破坏；

（2）建立海啸数据库。将关心区域内，所有可能发生地震的俯冲带均划分为若干个小单元，称为单位震源。计算每个单位震源发生单位错动时所产生的海啸数据，并存为海啸数据库；

（3）海啸数据记录与震源反演。当地震发生后，浮标会监测海啸数据。当首波数据监测完成后，即可以利用海啸数据库与反问题算法，得到计算域内任意位置的海啸波振幅与到达时间，并反演震源信息。

3　海啸数据库建立

3.1COMCOT数值模式简介

COMCOT（Cornell Multi-grid Coupled Tsunami Model）模式由美国Cornell大学开发，可以计算海啸的整个生命周期，包括海啸生成、传播、爬高与淹没的整个过程。根据研究尺度不同，模型可以灵活选择球面与直角坐标系下的线性与非线性浅水方程，所有方程均采用显式蛙跳有限差分法进行离散。具体方程与数植格式详见文献[7]。COMCOT数值模型是目前较为成熟的海啸模型，已经成功模拟了2004印尼苏门答腊海啸[8]等著名海啸。

图1　计算区域与水深条件图

3.2网格布置与计算参数

模型共布设3层嵌套网格，计算范围包括东海、台湾海峡、南中国海、日本及部分西太平洋海域，网格大小为2 min，首层网格地形数据采用ETOPO1数据，第二、三层网格地形数据采用海图与实测地形插值得到。计算范围及地形见图1，网格布置见表1。计算时间步长应满足CFL稳定条件，本文取固定时间步长为2 s。建立海啸数据库是基于线性假定，所以仅计算首层网格。

3.3初始条件与边界条件

COMCOT模式利用海床位移量来估算地震引起的初始水面高度，这种做法的前提是假定海水不可压缩，海床为刚性介质，水面变动与地震引起的地层错动同时发生[9]。海啸的初始条件的确定是通过输入地震断层参数，由弹性断层模型计算得到。

本文采用Okada断层模型，共需输入9个参数，分别为：震中经度、纬度、震源深度、断裂长度、断裂宽度、滑动量、走向角、倾角、滑移角。断层模型示意图见图2。

边界条件的设定中，水边界设为开边界；选择线性浅水波方程时，海陆边界设为垂直反射边界，选择非线性浅水波方程时，采用移动边界方案。

3.4日本东海岸单位滑块布置

图2　断层模型示意图

NOAA的太平洋海洋环境实验室（PMEL）将太平洋潜在地震海啸由15个板块804个海啸源组成，每个单位海啸源的大小均为100 km×50 km，并预先计算每个海啸源在深海传播的结果组成太平洋海啸数据库[10]。本文针对日本“3·11”海啸进行研究，所以仅选取震源附近的18个单位震源建立海啸数据库，单位震源布置见图3，震源参数见表2。

图3　单位震源分布图

4　反问题模式验证

4.1反问题模式的数学方法

假设震源区域被分割为NS个单位震源。第i个震源在地震发生后t时刻，在计算域内任意一点(x,y )所产生的水面波动可记为Gi(x,y,t )，那么t时刻处的总水面波动可表示为：

表1　网格布置与模式设置

表2　单位震源参数表

式（1）中，ci表示第i个单位震源的权重系数，需要通过反问题算法来确定。

地震发生后，一旦海啸波到达浮标处，浮标就会记录水面波动。假设浮标所处位置为(x0,y0) ，那么浮标所记录的数据集可以表示为Zk(x0,y0,tk) , { k=1,Nt}，式中，Nt代表浮标记录的数据总数。从而权重系数ci可由非负约束的最小二乘法求得，即式（2）：

式中，A为预先建立的海啸数据库，c为待求的权重系数，b为浮标记录数据集。

值得注意的是，浮标记录个数Nt需大于等于参于反演的单位震源个数Ns，这样ci才有唯一解。

4.2日本“3·11”海啸反问题模式验证

“3·11”海啸发生后，浮标21401在地震发生后78 min时记录了完整的海啸首波波形。本文采用表2所示的单位震源建立海啸数据库，并根据21401浮标前78 min的海啸实测数据进行反演，求解各单位震源的权重系数，计算结果表明主要发生错动的单元为3b—7b，其余单元的权重系数均为0，单位震源权重系数分布见表3，并根据OKADA模型计算地震发生后的初始水位，见图4。“3·11”地震海啸发生后，Tang[11]等利用浮标站压力数据，反演出震源可由六块沿着日本海沟的100 km×50 km的单位震源组合而成，断层的详细参数如表4所示，利用这些断层参数与OKADA模型作为初始条件，在海啸到岸前5 h计算出了美国海岸线30多个城市的海啸波爬高和淹没情况，模拟值与后来的观测值较为符合。Wei[12]等也利用同样的断层设置，模拟计算了海啸在日本近岸的爬高和淹没，结果表明，模拟计算的茨城县和青森县的淹没范围准确度达到85.5%。利用Tang等反演的震源计算地震发生后的初始水位，见图5。对比初始水位图可以看出，两者初始水位分布形态较为相似，且量级上也非常接近，均反演出初始海面位移有10 m的抬升和3 m的下沉。

表3　单位震源权重系数表

为进一步验证模型，采用21401，21419与21413 这3个浮标的实测数据对本文的预报结果进行验证（见图6）。由浮标水位验证图可以看出，反问题模式的预报结果在首波到达时间与首波波高上的计算上基本与实测值相符。虽然反问题模式在首波波谷处的计算上有所偏差，但是从预报的角度来说，本模式已经可以在极短的时间内回答预报所关心的问题。

5　反问题模式近岸预报准确度分析

表4　Tang等反演的震源信息

图4　本文反演震源所产生初始水位（单位：m）

图5　Tang反演震源所产生初始水位（单位：m）

图6　浮标实测数据验证图

图7　浙江沿海潮位站验证站布置图

表5　近岸潮位站验证表

图8　浙江沿海潮位站验证图

由于反问题模式的基本假定是海啸波在大洋中为线性传播，计算海啸数据库时也仅采用了线性浅水方程，所以对于近岸点直接采用式（1）进行预报难免会存在误差，可以考虑首先利用反问题法反演震源信息，再用COMCOT非线性方程求解近岸海啸波传播过程。本文采用前述两种方法分别计算“3·11”海啸期坎门、石浦、沈家门潮位站的海啸参数，COMCOT非线模拟时网格布置同3.2节，三层网格均参于计算。各潮位站布置见图7，各方法计算值与实测值的对比见表5与图8。由图与表可见，各方法所计算的首波波峰到达时间与实测资料相差均在20 min以内，具有较高的精度。反问题法直接预报的首波波峰值，较实测值有明显的偏小，这是因为，反问题模式并无考虑海啸波的近海非线性效应与海啸波的爬高过程。若采用反问题法反演的震源，采用COMCOT非线性模式进行计算，首波波峰有了明显的增大，与实测值的误差基本在20%以内，仅沈家门站偏差较大，这可能是地形精度的影响所致。本节计算说明，在近岸浅水区域，反问题预报法仅能对海啸波的到达时间作初步估计，若要取得较准确的预报结果，需利用反问题法反演的震源，采用非线性模型进行计算得到。

6　结论

本文采用COMCOT数值模式建立了日本东海岸的海啸数据库，并采用非负约束的最小二乘法构建海啸预报模式。应用本模式对日本“3·11”海啸进行反演，计算结果表明，本模式可在2 min内完成反演与预报，计算的初始水位与前人研究成果基本一致，计算的浮标水位波动与实测资料较为符合。在近岸处，反问题模式计算的首波波峰到达时间仍有一定的精度，但在首波波峰值的计算上有较大偏差。若需获得较准确的近岸预报结果，应采用反问题法反演的震源信息，应用非线性方程进行计算。

参考文献

[1] Titov V V, Mofjeld H O, Gonzalez F I, et al. Offshore Forecasting of Alaskan Tsunamis in Hawaii[M]//Tsunami Research at the End of a Critical Decade. Netherlands: Springer, 2001: 75-90.

[2] Wei Y, Cheung K F, Curtis G D, et al. Inverse Algorithm for Tsunami Forecasts[J]. Journal of Waterway, Port, Coastal, and Ocean Engineering, 2003, 129(2): 60-69.

[3] Johnson J M. Heterogeneous Coupling Along Alaska-Aleutians as Inferred From Tsunami, Seismic, and Geodetic Inversions[J]. Advances in Geophysics, 1998, 39: 1-116.

[4] Yamazaki Y, Wei Y, Cheung K F, et al. Forecast of Tsunamis from the Japan-Kuril-Kamchatka Source Region[J]. Natural Hazards, 2006, 38(3): 411-435.

[5]李林燕,毛献忠.南海海啸反问题预报模式[J].水动力学研究与进展:A辑, 2012, 27(1): 62-67.

[6]王培涛,于福江,赵联大,等. 2011年3月11日日本地震海啸越洋传播及对中国影响的数值分析[J].地球物理学报, 2012, 55(9): 3088-3096.

[7] Liu P L F, Woo S B, Cho Y S. Computer Programs for Tsunami Propagation and Inundation[J]. 1998.

[8] Wang X M, Liu P L F. An analysis of 2004 Sumatra earthquake fault plane mechanisms and Indian Ocean tsunami[J]. Journal of Hydraulic Research, 2006, 44(2): 147-154.

[9]潘文亮,王盛安. COMCOT数值模式的介绍和应用[J].海洋预报, 2009, 26(3): 45-52.

[10] Gica E, Spillane M C, Titov V V, et al. Development of the Forecast Propagation Database for NOAA's Short-Term Inundation Forecast for Tsunamis(SIFT)[R]. 2008.

[11] Tang L J, Titov V V, Bernard E N, et al. Direct Energy Estimation of the 2011 Japan Tsunami Using Deep-Ocean Pressure Measurements[J]. Journal of Geophysical Research: Oceans (1978-2012), 2012, 117(C8): C08008.

[12] Wei Y, Chamberlin C, Titov V V, et al. Modeling of the 2011 Japan Tsunami: Lessons for Near-Field Forecast[J]. Pure and Applied Geophysics, 2013, 170(6-8): 1309-1331.

Study on inversion forecasting model for East China Sea——A case study of Japan“3·11”tsunami

YING Chao，YU Pu-bing，MU Jin-bin，ZHOU Hua-min
（Zhejiang Inst．of Hydraulics＆Estuary，Hangzhou 310020 China）

Abstract：The Japan’s east coast was divided into 18 unit source to set up a tsunami database by using COMCOT numerical model. An inversion forecasting model was established by using least non-negative square method based on the database. The model was applied to 2011 Tohoku tsuami, the initial tsunami water level with 10m increase and 3m decrease calculated by the model was basically the same as previous research, the buoy level of prediction is in good agreement with measured data. Comparing with tsunami heights measured by tidal stations at coastal area of Zhejiang Province, the deviation of forecasted and measured value is large. But the prediction accuracy.can be greatly improved by solving COMCOT nonlinear equations with source parameters inversed by the forecasting model.

Key words：COMCOT；tsunami；inverse problem epicenter；least square method

作者简介：应超（1988-），男，工程师，硕士，主要从事海岸防灾减灾及海洋水环境研究。E-mail:xinqing928@126.com

基金项目：浙江省自然基金（LY13E090001）；浙江省创新人才培养项目（2012F20049）；浙江省公益基金（2014C33057）；浙江省公益技术研究社会发展项目（2014C33057）；浙江省科技计划（2015F50064）

收稿日期：2014-07-31

DOI:10.11737/j.issn.1003-0239.2015.03.005

中图分类号：P731.36

文献标识码：A

文章编号：1003-0239（2015）03-0036-07