劉廸森， 莫慶炎

（香港天文臺(tái)， 中國(guó) 香港）

0 引言

海嘯可以由海底地震、 海底巨大山泥傾瀉、海底火山爆發(fā)和隕石撞擊海洋而引發(fā)， 其中海底地震是引發(fā)海嘯最常見的原因， 占所有海嘯個(gè)案90%以上。 過去10年， 全球發(fā)生的大海嘯引起了世界關(guān)注， 如2004年12月印度尼西亞蘇門答臘北部以西的印度洋海底9 級(jí)地震引發(fā)跨海洋海嘯，及2011年3月在日本東北的太平洋發(fā)生9 級(jí)地震所產(chǎn)生的海嘯均造成嚴(yán)重的人命傷亡及破壞。 香港位處珠江口， 面對(duì)南海， 來自太平洋的海嘯有時(shí)會(huì)進(jìn)入南海， 而在馬尼拉海溝所產(chǎn)生的海嘯也會(huì)對(duì)華南沿岸地區(qū)以至香港構(gòu)成威脅。

香港天文臺(tái)于2005年底通過政府間海洋學(xué)委員會(huì)（IOC）的海嘯淹浸模型交換計(jì)劃（TIME）引進(jìn)了近場(chǎng)海嘯數(shù)值模型， 利用南海及香港的水深和地形數(shù)據(jù)， 模擬海嘯抵達(dá)香港時(shí)香港沿岸的海嘯高度①黃永德， 楊繼興， 李健威. 香港海嘯警報(bào)服務(wù)的新發(fā)展和數(shù)值模式的應(yīng)用. 第二屆粵港澳地震科技研討會(huì)論文集.中國(guó)澳門， 2006.。 受當(dāng)時(shí)技術(shù)限制， 該數(shù)值模型模擬海嘯所需的運(yùn)算長(zhǎng)達(dá)數(shù)十小時(shí)， 只能用作建立預(yù)設(shè)地震情景之下的模擬海嘯結(jié)果數(shù)據(jù)庫(kù)， 在需要時(shí)提取最接近地震情景的海嘯模擬數(shù)據(jù)以作參考， 而不能在發(fā)生地震后利用實(shí)時(shí)分析所得的地震參數(shù)進(jìn)行運(yùn)算和提供模擬海嘯的結(jié)果。

隨著海嘯數(shù)值模型技術(shù)不斷改進(jìn)及計(jì)算機(jī)運(yùn)算速度的不斷提升， 香港天文臺(tái)在2012年引進(jìn)了一套可應(yīng)用平行運(yùn)算技術(shù)操作的海嘯數(shù)值模型-Cornell Multigrid COupled Tsunami（COMCOT）model[1]。 該模型由美國(guó)康乃爾大學(xué)土木與環(huán)境工程學(xué)系劉立方教授所領(lǐng)導(dǎo)的團(tuán)隊(duì)研發(fā)， 大大縮短模擬海嘯的運(yùn)算時(shí)間， 并可與香港天文臺(tái)在2010年開始運(yùn)作的地震數(shù)據(jù)處理及分析系統(tǒng)結(jié)合， 實(shí)時(shí)把系統(tǒng)探測(cè)和分析所得的太平洋或南海地震參數(shù)， 輸入模型以模擬海嘯傳播過程， 并計(jì)算海嘯在海面上和抵岸時(shí)的情況， 為應(yīng)急決策提供及時(shí)和貼近實(shí)際情景的模擬結(jié)果。 本文利用日本311地震所引發(fā)的海嘯和過去香港曾經(jīng)錄得的海嘯記錄以驗(yàn)證COMCOT 模擬海嘯的能力， 并探討COMCOT 在香港天文臺(tái)海嘯預(yù)警工作上的應(yīng)用。

1 COMCOT 海嘯數(shù)值模型概述

COMCOT 是基于淺水波理論方程而發(fā)展的數(shù)值模型， 現(xiàn)時(shí)在全球多個(gè)海嘯預(yù)報(bào)中心及學(xué)術(shù)機(jī)構(gòu)的業(yè)務(wù)運(yùn)作或?qū)W術(shù)研究上被廣泛地應(yīng)用。

COMCOT 采用平行運(yùn)算技術(shù)以加快運(yùn)算速度，在發(fā)生地震后可以利用實(shí)時(shí)地震參數(shù)進(jìn)行運(yùn)算，并在短時(shí)間內(nèi)提供運(yùn)算結(jié)果。 此外， COMCOT 可以根據(jù)波高與水深的比值以及模擬范圍的大小，適當(dāng)?shù)剡x擇使用球坐標(biāo)系統(tǒng)或卡式坐標(biāo)系統(tǒng)、 網(wǎng)格大小、 線性或非線性波方程進(jìn)行模擬運(yùn)算， 在基本上不影響計(jì)算準(zhǔn)確度的情況下盡量縮短運(yùn)算的時(shí)間。 在大范圍運(yùn)算的深水區(qū)中， 由于海嘯在深水傳播時(shí)， 網(wǎng)格分辨率對(duì)模式的準(zhǔn)確度影響不大， 模式采用較大的網(wǎng)格及使用線性方程， 因應(yīng)模擬范圍大和需要考慮地球面弧度而采用球坐標(biāo)系統(tǒng)。 到近岸及淺水地區(qū)， 因需要考慮海洋底部磨擦項(xiàng)的影響及模擬海嘯造成的淹浸情況， 便采用非線性淺水波方程及高分辨率的網(wǎng)格以獲取更詳盡的海嘯高度。 由于模擬范圍縮小， 地球曲率效應(yīng)影響相對(duì)地輕微， 模式改為采用卡式坐標(biāo)系統(tǒng)以減少運(yùn)算時(shí)間[2]。

網(wǎng)格方面， COMCOT 采用嵌套方法， 網(wǎng)格層的數(shù)目視乎最外層與最內(nèi)層分辨率的差別而定。要確保模式運(yùn)行的穩(wěn)定性， 時(shí)間步長(zhǎng)的選取、 嵌套網(wǎng)格區(qū)的設(shè)定和水深的變化必須配合， 因此依據(jù)目標(biāo)地的地理情況融合水深和地形資料是成功運(yùn)行模式的重要工作。 根據(jù)香港的地理情況， 天文臺(tái)在南海制作了一套五層的嵌套網(wǎng)格（圖1）， 最外層的網(wǎng)格覆蓋南海（5° ～25°N 及105° ～125°E），而最內(nèi)層的網(wǎng)格包括整個(gè)香港（22°10′ ～22°35′ N及113°50′ ～114°30′ E）， 各層網(wǎng)格的分辨率從內(nèi)至外分別為1、 3、 9、 27 和81 弧秒（1 弧秒約等于30 m)。 各網(wǎng)格層之間的分辨率比例為1： 3。至于南海以外的太平洋則采用美國(guó) NOAA 的National Geophysical Data Centre（NGDC）所提供的ETOP1 地形及水深數(shù)據(jù)[3]， 分辨率為2 弧分。

南海水深和地形資料以美國(guó)國(guó)家航空航天局（NASA）和美國(guó)影像與制圖局（NIMA）合作通過航天飛機(jī)雷達(dá)地形測(cè)量任務(wù)制作的SRTM30 PLUS 數(shù)位高程資料（DEM）為基礎(chǔ)[4]， 并在近岸處根據(jù)較精確的區(qū)域資料作出修定。 SRTM30 PLUS 的最高分辨率為30 弧秒， 即約1 km。 而最內(nèi)層網(wǎng)格則采用香港地政總署的香港地形資料及香港海事處的香港境內(nèi)水深資料， 最高分辨率達(dá)30 m。

圖1 迭加五層嵌套網(wǎng)格Fig 1 Overlap of 5 layer nested grids

本文利用COMCOT 在模擬海嘯對(duì)香港近岸的影響時(shí)， 網(wǎng)格分辨率達(dá)30 m 及90 m 的最內(nèi)兩層網(wǎng)格層將采用非線性方程運(yùn)算， 而其余的網(wǎng)格層則采用線性淺方程。 一般在水深少于100 m 的情況下， 淺水非線性效應(yīng)比較明顯， 利用非線性淺方程運(yùn)算可以得到較佳的計(jì)算精度。 而天文臺(tái)在香港境內(nèi)制作的最內(nèi)兩層網(wǎng)格亦剛好在水深少于100 m 的范圍內(nèi)而分辨率也是最高。 但這兩層在運(yùn)算的過程中所需要的時(shí)間會(huì)相對(duì)較長(zhǎng)。 在水深大于100 m 的海洋上， 利用線性方程可節(jié)省運(yùn)算時(shí)間。

以一個(gè)南海馬尼拉海溝發(fā)生的地震所引發(fā)的海嘯為例， 使用全部五層嵌套網(wǎng)格作運(yùn)算的時(shí)間為45 min。 但如果減少最內(nèi)一層而只使用四層嵌套網(wǎng)格， 運(yùn)算時(shí)間則可縮短至7 min。

2 模擬海嘯及驗(yàn)證

COMCOT 在一臺(tái)內(nèi)有兩顆八核中央處理器的服務(wù)器運(yùn)行， 運(yùn)算速度為每秒 294 兆浮點(diǎn)（GFLOPS）。 模擬海嘯之前， 需要把地震數(shù)據(jù)處理及分析系統(tǒng)所計(jì)算的各項(xiàng)參數(shù)和模擬時(shí)間的長(zhǎng)度、嵌套網(wǎng)格等數(shù)據(jù)輸入COMCOT 內(nèi)作好準(zhǔn)備。

圖2 斷層模型Fig 2 Fault model

初始條件方面， 模型假設(shè)水面初始垂直位移幅度與海底因斷層而造成的變形幅度相等， 并采用了Manshinha and Smylie 斷層模型（1971）[5]及配合日本氣象廳提供的經(jīng)驗(yàn)方程計(jì)算斷層在海底的變形幅度， 例如長(zhǎng)度（L）、 寬度（W）及滑動(dòng)距離（D）（圖2）。 而震中位置、 矩震級(jí)（Mw）、 震源深度、 走向角、 傾角及滑動(dòng)角可透過香港天文臺(tái)的地震數(shù)據(jù)處理及分析系統(tǒng)在地震發(fā)生后短時(shí)間內(nèi)自動(dòng)計(jì)算出來。

我們利用日本311 地震在日本及香港錄得的海嘯資料和其它香港曾經(jīng)錄得的海嘯個(gè)案以驗(yàn)證COMCOT 模擬海嘯的能力。

2.1 2011年3月11日在日本東北發(fā)生九級(jí)地震所引發(fā)的海嘯

根據(jù)美國(guó)地質(zhì)勘察局（USGS）的地震參數(shù)（表1）及利用日本氣象廳提供的經(jīng)驗(yàn)方程求得的斷層模型變形幅度參數(shù)（表2）， COMCOT 模擬海嘯影響日本近岸的情況， 模擬時(shí)間長(zhǎng)度為20 h（圖3）。

表1 2011年3月11日在日本東北發(fā)生9 級(jí)地震的地震參數(shù)Table 1 Earthquake parameters for M9 Tohoku earthquake on 11 Mar, 2011

表2 2011年3月11日在日本東北發(fā)生九級(jí)地震的斷層模型變形幅度參數(shù)Table 2 Fault model parameters for M9 Tohoku earthquake on 11 Mar，2011

模擬結(jié)果顯示影響相馬及石卷的第一個(gè)海嘯波高度分別為9.5 m 和7.0 m（圖4）， 與當(dāng)?shù)貎蓚€(gè)驗(yàn)潮站分別錄得的9.3 m 和7.6 m 接近（圖5）。 在模擬最大海嘯波高度圖顯示（圖6）， 其他日本東北沿岸的海嘯波高度跟驗(yàn)潮站錄得的海嘯波高度亦大致相約。

圖3 利用COMCOT 模擬2011年3月11日在日本東北發(fā)生的海嘯Fig.3 Simulation of Tohoku tsunami on 11 March，2011 by using COMCOT

圖4 COMCOT 模擬日本沿岸相馬及石卷的海嘯波的高度時(shí)間序列Fig.4 Time series of tsunami wave height of Sōma， Fukushima and Ishinomaki， Miyagi of the coastal area of Japan simulated by COMCOT

海嘯波在地震發(fā)生后約8 h 到達(dá)香港（圖7）。香港三個(gè)驗(yàn)潮站（圖8）所錄得的海嘯到達(dá)時(shí)間及第一個(gè)海嘯波的高度和COMCOT 模擬結(jié)果的比較列于表3。

結(jié)果顯示各潮站錄得的第一個(gè)海嘯波的高度與COMCOT 的模擬頗為接近， 相差不多于0.1 m，而錄得的海嘯到達(dá)時(shí)間與模擬所得也大致相約，差距少于15 min。 此外， 在驗(yàn)潮站的時(shí)間序列圖中顯示（圖9、 10 及11)， 海嘯到達(dá)后5 h 內(nèi)的海嘯波幅變化與模擬的波幅變化接近同步， 顯示模式模擬海嘯拍岸時(shí)的變化表現(xiàn)相當(dāng)理想。

2.2 香港錄得的其它海嘯記錄

除日本311 地震所引發(fā)的越洋海嘯， 香港自1950年開始運(yùn)作自動(dòng)驗(yàn)潮站以來另外曾錄得6 次海嘯， 全部屬于輕微水位變化（少于0.5 m）。 我們利用USGS 的地震參數(shù)模擬這6 次海嘯事件， 模擬結(jié)果與香港鲗魚涌驗(yàn)潮站的觀測(cè)數(shù)據(jù)作比較（表4）。 結(jié)果顯示其中5 次的模擬結(jié)果與實(shí)測(cè)的數(shù)值頗為接近， 只有1988年6月24日呂宋海峽所發(fā)生的5.7 級(jí)地震， COMCOT 沒有顯示任何海嘯到達(dá)香港， 似乎模式的海嘯波在抵達(dá)香港前已完全減弱。

圖5日本沿岸驗(yàn)潮站在2011年3月11日錄得的最大海嘯波高度圖（數(shù)據(jù)來源：日本氣象廳）Fig.5 Maximum tsunami wave height observed by tide gauging system along the coast of Japan on 11 March，2011

圖6 COMCOT 模擬日本沿岸在2011年3月11日的最大海嘯波高度圖Fig.6 Maximum tsunami wave height simulated by COMCOT along the coast of Japan on 11 March，2011

圖7 COMCOT 模擬日本311 海嘯到達(dá)香港的情況(時(shí)間序列：左至右，上至下)Fig.7 Arrival of tsunami at Hong Kong for 311 Japan tsunami simulated by COMCOT（Time sequence from left to right,up to down）

圖8 香港其中三個(gè)驗(yàn)潮站橫瀾島、鲗魚涌和石壁的位置Fig.8 Positions of three tide gauge stations of Hong Kong (Wagland Island，Quarry Bay, Shek Pik)

表3 香港驗(yàn)潮站在日本311 地震錄得的海嘯到達(dá)時(shí)間及海嘯高度和COMCOT 模擬結(jié)果的比較Table 3 Comparison of arrival time and wave height recorded by Hong Kong tide gauge stations and simulated by COMCOT for 311 Japan earthquake

圖9 石壁驗(yàn)潮站的實(shí)測(cè)及模擬海嘯波的時(shí)間序列Fig.9 Time series of measured and simulated tsunami wave at Shek Pik tide gauge station

圖11 橫瀾島驗(yàn)潮站的實(shí)測(cè)及模擬海嘯波的時(shí)間序列Fig.11 Time series of measured and simulated tsunami wave at Waglan Island tide gauge station

表4 鲗魚涌驗(yàn)潮站曾錄得的六次海嘯與COMCOT 的模擬比較Table 4 Comparison of the six recorded tsunamis at Quarry Bay and COMCOT simulation

3 COMCOT 的應(yīng)用

應(yīng)用平行運(yùn)算技術(shù)， COMCOT 在一臺(tái)內(nèi)有兩顆八核中央處理器， 運(yùn)算速度為每秒294 兆浮點(diǎn)（GFLOPS）的服務(wù)器運(yùn)行， 模擬日本311 地震后從太平洋到香港長(zhǎng)約八小時(shí)的海嘯過程， 運(yùn)算時(shí)間只需約75 min。 而一個(gè)在南海馬尼拉海溝發(fā)生的地震所引發(fā)的海嘯， 到達(dá)香港的時(shí)間約需3 h，COMCOT 的運(yùn)算時(shí)間也只不過是45 min 左右， 可見利用COMCOT 配合香港天文臺(tái)的實(shí)時(shí)地震數(shù)據(jù)處理及分析系統(tǒng)， 足以支持香港天文臺(tái)海嘯預(yù)警系統(tǒng)的業(yè)務(wù)應(yīng)用。

香港天文臺(tái)的地震數(shù)據(jù)處理及分析系統(tǒng)一般可以在6 級(jí)或以上的地震發(fā)生后大約6～8 min 計(jì)算出地震發(fā)生的時(shí)間、 位置和震級(jí)。 而地震的震源機(jī)制亦可在發(fā)震后大約20 min 透過香港天文臺(tái)的地震數(shù)據(jù)處理及分析系統(tǒng)計(jì)算出來， 或經(jīng)由USGS 提供。 如果馬尼拉海溝發(fā)生地震， COMCOT 在地震發(fā)生后大約一個(gè)小時(shí)便可以提供數(shù)值模式的計(jì)算結(jié)果， 為海嘯預(yù)警的應(yīng)急工作帶來莫大裨益。

4 總結(jié)

本文利用日本311 海嘯個(gè)案及其他六次香港曾錄得的海嘯記錄驗(yàn)證了COMCOT 模擬海嘯的能力。 在311 海嘯個(gè)案的模擬結(jié)果顯示， 影響日本東北沿岸的海嘯波高度跟驗(yàn)潮站錄得的海嘯波高度大致相約， 而香港各潮站錄得的第一個(gè)海嘯波的高度與模擬頗為接近， 相差不多于0.1 m， 而錄得的海嘯到達(dá)時(shí)間與模擬所得也大致相約， 差距少于15 min。 此外， 由于COMCOT 應(yīng)用平行運(yùn)算技術(shù)， 縮短模擬所需的運(yùn)算時(shí)間， 使其可以有效地應(yīng)用于海嘯應(yīng)急的業(yè)務(wù)運(yùn)作， 大大提升香港天文臺(tái)海嘯預(yù)警方面的能力。 再者， 一個(gè)在南海馬尼拉海溝發(fā)生的地震所引發(fā)的海嘯， 由于到達(dá)香港的時(shí)間約需3 h， 如果只使用四層嵌套網(wǎng)格模擬， 運(yùn)算的時(shí)間則大大縮短為只需7 min， 對(duì)初步評(píng)估海嘯的影響起了很大的參考價(jià)值。

[1] Xiaoming Wang. Technical document for COMCOT v.1.7[M]. New Zealand： Institute of Geological & Nuclear Science， 2009.

[2] Philip L.F. Liu. Seung -Buhm Woo， Yong -Sik Cho.Computer Programs for Tsunami Propagation and Inundation [M]. Ithaca， NY： School of Civil and Environmental Engineering Cornell University， 1998.

[3] Amante, C. B. W. Eakins, ETOPO1 1 Arc-Minute Global Relief Model： Procedures, Data Sources and Analysis [C].NOAA Technical Memorandum NESDIS NGDC-24, 19 pp,March 2009.

[4] Farr， T.G.， M. Kobrick. Shuttle Radar Topography Mission produces a wealth of data [J]. Amer. Geophys. Union Eos， 2000， 81： 583-585.

[5] Manshinha， L.， and D. Smylie. The displacement field of inclined faults [J]. Bulletin of the Seismological Society of America， 1971， 61（5）： 1433-1440.