東海海嘯反問題預(yù)報(bào)模式研究——以日本“3·11”海嘯為例

應(yīng)超，于普兵，穆錦斌，周華民

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

東海海嘯反問題預(yù)報(bào)模式研究——以日本“3·11”海嘯為例

應(yīng)超，于普兵，穆錦斌，周華民

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

摘要：采用COMCOT數(shù)值模式，將日本東海岸劃分為18個(gè)單位震源，建立了日本東海岸海嘯數(shù)據(jù)庫并用非負(fù)約束的最小二乘法建立反問題預(yù)報(bào)模式。將本模式應(yīng)用于日本“3·11”海嘯，計(jì)算所得的海嘯初始水位有10 m的抬高與3 m的降低，與前人研究結(jié)果基本一致，預(yù)報(bào)的浮標(biāo)水位與實(shí)測資料符合良好。對比浙江省近岸潮位站實(shí)測海嘯波高，預(yù)報(bào)值與實(shí)測值偏差較大。若采用反問題反演的震源，通過COMCOT非線性模式求解近岸水位，可以大幅提高預(yù)報(bào)精度。

關(guān)鍵詞：COMCOT；海嘯；反問題；最小二乘法

1　引言

海嘯是一種具有強(qiáng)大破壞力的海浪，具有傳播速度快、破壞力大等特點(diǎn)。準(zhǔn)確及時(shí)的海嘯預(yù)報(bào)是減小海嘯災(zāi)害的重要舉措。海嘯預(yù)報(bào)中，海嘯波的到達(dá)時(shí)間、波高與預(yù)報(bào)速度是需要重點(diǎn)關(guān)心的問題。

反問題預(yù)報(bào)方法是海嘯預(yù)報(bào)中的一種重要方法，是指用預(yù)先確定好的地震源參數(shù)把地震源劃分為不同的斷層，用不同斷層單位活動(dòng)產(chǎn)生的水位過程建立數(shù)據(jù)庫，然后回歸分析記錄的海嘯波形確定斷層滑動(dòng)分布和地震震級(jí)。國內(nèi)外很多學(xué)者對反問題預(yù)報(bào)模式進(jìn)行過研究，并取得了豐碩的成果。Titov[1]等利用斷層滑動(dòng)分布通過反問題算法確定越洋海嘯實(shí)時(shí)模型的初始條件，Wei[2]等采用反問題方法通過預(yù)先計(jì)算好的潮位過程的數(shù)據(jù)庫確定遠(yuǎn)場海嘯波形預(yù)測海嘯。在Johnson[3]地震源參數(shù)基礎(chǔ)上，Wei[2]等用反問題方法準(zhǔn)確確定滑層分布，并準(zhǔn)確計(jì)算夏威夷各潮位站的波高。Yamazaki[4]等用此方法很好地模擬了1944年Tonakai和1994年Kuril地震海嘯通過太平洋在日本產(chǎn)生的海嘯過程。李林燕[5]等采用反問題方法對馬尼拉海溝建立海嘯數(shù)據(jù)庫并建立預(yù)報(bào)模式，但并未采用實(shí)測海嘯資料進(jìn)行驗(yàn)證，僅采用正反問題的計(jì)算結(jié)果進(jìn)行相互驗(yàn)證。以往的這些研究，均未對中國東海附近海域進(jìn)行研究，COMCOT模型及反問題方法在中國東海海域的適用性均需進(jìn)行驗(yàn)證與探討，另外反問題方法在近岸預(yù)報(bào)的適用性問題也應(yīng)進(jìn)行進(jìn)一步分析。

北京時(shí)間2011年3月11日13時(shí)46分，日本東北部海域發(fā)生了9.0級(jí)特大地震，從而誘發(fā)了災(zāi)害性海嘯，造成了巨大的人員傷亡與財(cái)產(chǎn)損失。此次海嘯已經(jīng)對我國浙東沿海造成影響，地震發(fā)生后6—9 h，浙江沈家門、大陳、坎門、石坪、石浦、健跳、福建東山及廣東汕頭、汕尾等潮位站先后監(jiān)測到振幅為10—55 cm的海嘯波[6]。因此研究日本“3·11”海嘯反問題預(yù)報(bào)模式對我國東海海嘯預(yù)警機(jī)制的建立與完善具有很大的借鑒意義。

2　反問題預(yù)報(bào)方法

反問題預(yù)報(bào)方法如下：

（1）設(shè)置海嘯浮標(biāo)。海嘯浮標(biāo)一般設(shè)置在地形平坦水深較大的海域，距震源有一定距離但不宜太遠(yuǎn)，從而保證浮標(biāo)既能較快采集數(shù)據(jù)又不被海嘯破壞；

（2）建立海嘯數(shù)據(jù)庫。將關(guān)心區(qū)域內(nèi)，所有可能發(fā)生地震的俯沖帶均劃分為若干個(gè)小單元，稱為單位震源。計(jì)算每個(gè)單位震源發(fā)生單位錯(cuò)動(dòng)時(shí)所產(chǎn)生的海嘯數(shù)據(jù)，并存為海嘯數(shù)據(jù)庫；

（3）海嘯數(shù)據(jù)記錄與震源反演。當(dāng)?shù)卣鸢l(fā)生后，浮標(biāo)會(huì)監(jiān)測海嘯數(shù)據(jù)。當(dāng)首波數(shù)據(jù)監(jiān)測完成后，即可以利用海嘯數(shù)據(jù)庫與反問題算法，得到計(jì)算域內(nèi)任意位置的海嘯波振幅與到達(dá)時(shí)間，并反演震源信息。

3　海嘯數(shù)據(jù)庫建立

3.1COMCOT數(shù)值模式簡介

COMCOT（Cornell Multi-grid Coupled Tsunami Model）模式由美國Cornell大學(xué)開發(fā)，可以計(jì)算海嘯的整個(gè)生命周期，包括海嘯生成、傳播、爬高與淹沒的整個(gè)過程。根據(jù)研究尺度不同，模型可以靈活選擇球面與直角坐標(biāo)系下的線性與非線性淺水方程，所有方程均采用顯式蛙跳有限差分法進(jìn)行離散。具體方程與數(shù)植格式詳見文獻(xiàn)[7]。COMCOT數(shù)值模型是目前較為成熟的海嘯模型，已經(jīng)成功模擬了2004印尼蘇門答臘海嘯[8]等著名海嘯。

圖1　計(jì)算區(qū)域與水深條件圖

3.2網(wǎng)格布置與計(jì)算參數(shù)

模型共布設(shè)3層嵌套網(wǎng)格，計(jì)算范圍包括東海、臺(tái)灣海峽、南中國海、日本及部分西太平洋海域，網(wǎng)格大小為2 min，首層網(wǎng)格地形數(shù)據(jù)采用ETOPO1數(shù)據(jù)，第二、三層網(wǎng)格地形數(shù)據(jù)采用海圖與實(shí)測地形插值得到。計(jì)算范圍及地形見圖1，網(wǎng)格布置見表1。計(jì)算時(shí)間步長應(yīng)滿足CFL穩(wěn)定條件，本文取固定時(shí)間步長為2 s。建立海嘯數(shù)據(jù)庫是基于線性假定，所以僅計(jì)算首層網(wǎng)格。

3.3初始條件與邊界條件

COMCOT模式利用海床位移量來估算地震引起的初始水面高度，這種做法的前提是假定海水不可壓縮，海床為剛性介質(zhì)，水面變動(dòng)與地震引起的地層錯(cuò)動(dòng)同時(shí)發(fā)生[9]。海嘯的初始條件的確定是通過輸入地震斷層參數(shù)，由彈性斷層模型計(jì)算得到。

本文采用Okada斷層模型，共需輸入9個(gè)參數(shù)，分別為：震中經(jīng)度、緯度、震源深度、斷裂長度、斷裂寬度、滑動(dòng)量、走向角、傾角、滑移角。斷層模型示意圖見圖2。

邊界條件的設(shè)定中，水邊界設(shè)為開邊界；選擇線性淺水波方程時(shí)，海陸邊界設(shè)為垂直反射邊界，選擇非線性淺水波方程時(shí)，采用移動(dòng)邊界方案。

3.4日本東海岸單位滑塊布置

圖2　斷層模型示意圖

NOAA的太平洋海洋環(huán)境實(shí)驗(yàn)室（PMEL）將太平洋潛在地震海嘯由15個(gè)板塊804個(gè)海嘯源組成，每個(gè)單位海嘯源的大小均為100 km×50 km，并預(yù)先計(jì)算每個(gè)海嘯源在深海傳播的結(jié)果組成太平洋海嘯數(shù)據(jù)庫[10]。本文針對日本“3·11”海嘯進(jìn)行研究，所以僅選取震源附近的18個(gè)單位震源建立海嘯數(shù)據(jù)庫，單位震源布置見圖3，震源參數(shù)見表2。

圖3　單位震源分布圖

4　反問題模式驗(yàn)證

4.1反問題模式的數(shù)學(xué)方法

假設(shè)震源區(qū)域被分割為NS個(gè)單位震源。第i個(gè)震源在地震發(fā)生后t時(shí)刻，在計(jì)算域內(nèi)任意一點(diǎn)(x,y )所產(chǎn)生的水面波動(dòng)可記為Gi(x,y,t )，那么t時(shí)刻處的總水面波動(dòng)可表示為：

表1　網(wǎng)格布置與模式設(shè)置

表2　單位震源參數(shù)表

式（1）中，ci表示第i個(gè)單位震源的權(quán)重系數(shù)，需要通過反問題算法來確定。

地震發(fā)生后，一旦海嘯波到達(dá)浮標(biāo)處，浮標(biāo)就會(huì)記錄水面波動(dòng)。假設(shè)浮標(biāo)所處位置為(x0,y0) ，那么浮標(biāo)所記錄的數(shù)據(jù)集可以表示為Zk(x0,y0,tk) , { k=1,Nt}，式中，Nt代表浮標(biāo)記錄的數(shù)據(jù)總數(shù)。從而權(quán)重系數(shù)ci可由非負(fù)約束的最小二乘法求得，即式（2）：

式中，A為預(yù)先建立的海嘯數(shù)據(jù)庫，c為待求的權(quán)重系數(shù)，b為浮標(biāo)記錄數(shù)據(jù)集。

值得注意的是，浮標(biāo)記錄個(gè)數(shù)Nt需大于等于參于反演的單位震源個(gè)數(shù)Ns，這樣ci才有唯一解。

4.2日本“3·11”海嘯反問題模式驗(yàn)證

“3·11”海嘯發(fā)生后，浮標(biāo)21401在地震發(fā)生后78 min時(shí)記錄了完整的海嘯首波波形。本文采用表2所示的單位震源建立海嘯數(shù)據(jù)庫，并根據(jù)21401浮標(biāo)前78 min的海嘯實(shí)測數(shù)據(jù)進(jìn)行反演，求解各單位震源的權(quán)重系數(shù)，計(jì)算結(jié)果表明主要發(fā)生錯(cuò)動(dòng)的單元為3b—7b，其余單元的權(quán)重系數(shù)均為0，單位震源權(quán)重系數(shù)分布見表3，并根據(jù)OKADA模型計(jì)算地震發(fā)生后的初始水位，見圖4?！?·11”地震海嘯發(fā)生后，Tang[11]等利用浮標(biāo)站壓力數(shù)據(jù)，反演出震源可由六塊沿著日本海溝的100 km×50 km的單位震源組合而成，斷層的詳細(xì)參數(shù)如表4所示，利用這些斷層參數(shù)與OKADA模型作為初始條件，在海嘯到岸前5 h計(jì)算出了美國海岸線30多個(gè)城市的海嘯波爬高和淹沒情況，模擬值與后來的觀測值較為符合。Wei[12]等也利用同樣的斷層設(shè)置，模擬計(jì)算了海嘯在日本近岸的爬高和淹沒，結(jié)果表明，模擬計(jì)算的茨城縣和青森縣的淹沒范圍準(zhǔn)確度達(dá)到85.5%。利用Tang等反演的震源計(jì)算地震發(fā)生后的初始水位，見圖5。對比初始水位圖可以看出，兩者初始水位分布形態(tài)較為相似，且量級(jí)上也非常接近，均反演出初始海面位移有10 m的抬升和3 m的下沉。

表3　單位震源權(quán)重系數(shù)表

為進(jìn)一步驗(yàn)證模型，采用21401，21419與21413 這3個(gè)浮標(biāo)的實(shí)測數(shù)據(jù)對本文的預(yù)報(bào)結(jié)果進(jìn)行驗(yàn)證（見圖6）。由浮標(biāo)水位驗(yàn)證圖可以看出，反問題模式的預(yù)報(bào)結(jié)果在首波到達(dá)時(shí)間與首波波高上的計(jì)算上基本與實(shí)測值相符。雖然反問題模式在首波波谷處的計(jì)算上有所偏差，但是從預(yù)報(bào)的角度來說，本模式已經(jīng)可以在極短的時(shí)間內(nèi)回答預(yù)報(bào)所關(guān)心的問題。

5　反問題模式近岸預(yù)報(bào)準(zhǔn)確度分析

表4　Tang等反演的震源信息

圖4　本文反演震源所產(chǎn)生初始水位（單位：m）

圖5　Tang反演震源所產(chǎn)生初始水位（單位：m）

圖6　浮標(biāo)實(shí)測數(shù)據(jù)驗(yàn)證圖

圖7　浙江沿海潮位站驗(yàn)證站布置圖

表5　近岸潮位站驗(yàn)證表

圖8　浙江沿海潮位站驗(yàn)證圖

由于反問題模式的基本假定是海嘯波在大洋中為線性傳播，計(jì)算海嘯數(shù)據(jù)庫時(shí)也僅采用了線性淺水方程，所以對于近岸點(diǎn)直接采用式（1）進(jìn)行預(yù)報(bào)難免會(huì)存在誤差，可以考慮首先利用反問題法反演震源信息，再用COMCOT非線性方程求解近岸海嘯波傳播過程。本文采用前述兩種方法分別計(jì)算“3·11”海嘯期坎門、石浦、沈家門潮位站的海嘯參數(shù)，COMCOT非線模擬時(shí)網(wǎng)格布置同3.2節(jié)，三層網(wǎng)格均參于計(jì)算。各潮位站布置見圖7，各方法計(jì)算值與實(shí)測值的對比見表5與圖8。由圖與表可見，各方法所計(jì)算的首波波峰到達(dá)時(shí)間與實(shí)測資料相差均在20 min以內(nèi)，具有較高的精度。反問題法直接預(yù)報(bào)的首波波峰值，較實(shí)測值有明顯的偏小，這是因?yàn)?，反問題模式并無考慮海嘯波的近海非線性效應(yīng)與海嘯波的爬高過程。若采用反問題法反演的震源，采用COMCOT非線性模式進(jìn)行計(jì)算，首波波峰有了明顯的增大，與實(shí)測值的誤差基本在20%以內(nèi)，僅沈家門站偏差較大，這可能是地形精度的影響所致。本節(jié)計(jì)算說明，在近岸淺水區(qū)域，反問題預(yù)報(bào)法僅能對海嘯波的到達(dá)時(shí)間作初步估計(jì)，若要取得較準(zhǔn)確的預(yù)報(bào)結(jié)果，需利用反問題法反演的震源，采用非線性模型進(jìn)行計(jì)算得到。

6　結(jié)論

本文采用COMCOT數(shù)值模式建立了日本東海岸的海嘯數(shù)據(jù)庫，并采用非負(fù)約束的最小二乘法構(gòu)建海嘯預(yù)報(bào)模式。應(yīng)用本模式對日本“3·11”海嘯進(jìn)行反演，計(jì)算結(jié)果表明，本模式可在2 min內(nèi)完成反演與預(yù)報(bào)，計(jì)算的初始水位與前人研究成果基本一致，計(jì)算的浮標(biāo)水位波動(dòng)與實(shí)測資料較為符合。在近岸處，反問題模式計(jì)算的首波波峰到達(dá)時(shí)間仍有一定的精度，但在首波波峰值的計(jì)算上有較大偏差。若需獲得較準(zhǔn)確的近岸預(yù)報(bào)結(jié)果，應(yīng)采用反問題法反演的震源信息，應(yīng)用非線性方程進(jìn)行計(jì)算。

參考文獻(xiàn)

[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]李林燕,毛獻(xiàn)忠.南海海嘯反問題預(yù)報(bào)模式[J].水動(dòng)力學(xué)研究與進(jìn)展:A輯, 2012, 27(1): 62-67.

[6]王培濤,于福江,趙聯(lián)大,等. 2011年3月11日日本地震海嘯越洋傳播及對中國影響的數(shù)值分析[J].地球物理學(xué)報(bào), 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數(shù)值模式的介紹和應(yīng)用[J].海洋預(yù)報(bào), 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

作者簡介：應(yīng)超（1988-），男，工程師，碩士，主要從事海岸防災(zāi)減災(zāi)及海洋水環(huán)境研究。E-mail:xinqing928@126.com

基金項(xiàng)目：浙江省自然基金（LY13E090001）；浙江省創(chuàng)新人才培養(yǎng)項(xiàng)目（2012F20049）；浙江省公益基金（2014C33057）；浙江省公益技術(shù)研究社會(huì)發(fā)展項(xiàng)目（2014C33057）；浙江省科技計(jì)劃（2015F50064）

收稿日期：2014-07-31

DOI:10.11737/j.issn.1003-0239.2015.03.005

中圖分類號(hào)：P731.36

文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1003-0239（2015）03-0036-07