曹永港 劉長建 劉同木馮硯青 劉愉強 廖世智

1）國家海洋局南海調查技術中心，廣州市海珠區(qū)新港西路155號 510300

2）天津大學水利工程仿真與安全國家重點實驗室，天津市南開區(qū)衛(wèi)津路92號 300072

3）國家海洋局南海維權技術與應用重點實驗室，廣州市新港中路353號 510310

4）水沙科學與水災害防治湖南省重點實驗室，長沙市天心區(qū)萬家麗南路2段960號 410114

0 引言

海嘯是由水下地震、火山爆發(fā)、水下塌陷及滑坡等大地活動或氣象變化產生的破壞性海浪，其中，由水下地震產生的地震海嘯較為普遍，表現為海底地形急劇升降變動引起的海水強烈擾動。根據海嘯傳播至岸邊時間的長短，它可分為近海地震海嘯和遠洋地震海嘯。海嘯的發(fā)生頻率雖然較低，可一旦發(fā)生，破壞力是巨大的。2004年12月26日印度尼西亞蘇門答臘島西北附近海域發(fā)生了MW9.1強烈地震，地震在印度洋沿岸引發(fā)罕見的地震海嘯并傳播到全球各大洋（Yeh et al，2007）。2011年3月11日日本東北部海域發(fā)生了MW9.0強烈地震，這是日本近海發(fā)生的最強烈的一次大地震，并引發(fā)了太平洋范圍的海嘯（王培濤等，2012）。據統(tǒng)計，每2年全球發(fā)生1次局地破壞性海嘯，每10年發(fā)生1次越洋大海嘯（于福江等，2011）。Jing等（2013）研究發(fā)現，地震引發(fā)海嘯的概率很高，而在過去約100年里，震級≥8.8的地震都觸發(fā)了海嘯。

我國處于太平洋西部，海區(qū)遼闊，海岸線長達18000km，近海海域大多位于環(huán)太平洋地震帶上，歷史上地震海嘯時有發(fā)生，是世界上最早記錄到地震海嘯的國家。據不完全統(tǒng)計，公元47～2004年，中國沿海共發(fā)生29次地震海嘯，其中，8～9次為破壞性海嘯。雖然自1949年新中國成立以來我國并未發(fā)生破壞性海嘯，但近期南海地震頻發(fā)，一旦誘發(fā)海嘯將給我國南部、東南部發(fā)達城市造成巨大損失（包澄瀾等，2005；楊馬陵等，2005；Liu et al，2007）。我國歷史上曾經發(fā)生過多次海嘯，其中，臺灣地區(qū)周圍是海嘯的高發(fā)區(qū)域，其次是大陸架區(qū)域，低發(fā)區(qū)是渤海區(qū)域。在各省份中，按照有歷史記錄的海嘯分析，浙江的發(fā)生次數最多，高達45次，其次是江蘇、山東、上海、福建、臺灣、廣東。按照已確定的海嘯記錄分析，臺灣發(fā)生的次數最多，有8次，其次是山東、廣東和浙江。不論哪個數據都表明，我國是一個海嘯發(fā)生的危險區(qū)，而南海區(qū)域是中國受海嘯威脅最大的地方，這是因為呂宋島西側和馬尼拉海溝一帶是地震活躍帶（潘文亮等，2009）。Li等（2016）認為，常用的均勻滑動模型明顯低估了海嘯災害。因此，研究建立廣東沿海海嘯預警機制，在海嘯發(fā)生后迅速預報海嘯到達時間、海嘯最大波高，可為海嘯防災減災提供有力支持，以期最大化地減少海嘯造成的損失。

現有的海嘯計算模式中，一般不含有潮汐模塊。這是由于通常不考慮潮汐與海嘯的非線性作用，將計算出的海嘯波高與潮位值進行線性疊加來得到關心站點的水位。但是，對于近岸淺水域，潮汐與海嘯的非線性作用會明顯增大，故影響海嘯波的波幅及到達時間。本文利用COMCOT模型和全球潮汐模式TPXO7.2，建立廣東沿海天文潮與海嘯的耦合數學模型，針對馬尼拉地震帶海嘯，通過設計海嘯震源計算，分析廣東省沿海海嘯風險。

1 數值模型

COMCOT模型（Cornell Multi-grid Coupled Tsunamimodel）是由Cornell大學土木與環(huán)境工程系Philip Liu研究組開發(fā)基于淺水長波方程的海嘯數值計算模式，該模式已被多次用來模擬歷史海嘯事件，如1992年印度尼西亞Flores Islands海嘯（Liu et al，1995）、2003年Algeria海嘯（Wand et al，2005）、2004年印度洋大海嘯（Wand et al，2006）以及2011年日本海嘯（應超等，2015；Jing et al，2012），均有效模擬了海嘯生成、傳播、爬高和淹沒的整個過程。模式可根據研究范圍，靈活選擇球面與直角坐標系下的線性與非線性淺水方程，所有方程均采用顯式蛙跳有限差分法進行離散。

1.1 基本方程

球坐標下，考慮科氏力忽略海底摩擦，線性淺水方程為

式中，η為相對于平均海平面的自由表面位移；R為地球半徑；H＝η＋h為總的水深，h為靜水深；P為沿緯度單位寬度的通量；Q為沿經度單位寬度的通量；φ、ψ分別為經度和緯度；f為科氏力系數；g為重力加速度。海嘯傳播至近岸，采用笛卡爾坐標下非線性淺水波方程，并考慮底摩擦效應，其連續(xù)方程和動量方程為

式中，P、Q分別為x和y方向上的體積通量；n為曼寧粗糙系數。

1.2 數值格式

對式（1）～（5）采用了蛙跳有限差分格式進行離散，非線性對流項采用迎風向上格式離散。圖1為交錯格式示意圖。由圖1可見，波高和體積通量在空間和時間上都是交錯的，波高及水深位于網格中心，體積通量位于網格線上。數值格式的截斷誤差為o（Δx2，Δy2，Δt2），在時空上均有二階精度。COMCOT是比較成熟的海嘯數值模型，具體算法詳見文獻（Wang，2009）。

圖1 交錯格式示意圖

1.3 邊界條件

海嘯初始條件的確定是通過輸入地震斷層參數，由彈性斷層模型計算得到。模式中可選擇Mansinha和Smylie的彈性半空間錯移模型和Okada的理論模型，這也是現階段常用的2種彈性斷層模型。彈性斷層模型基于彈性錯移理論建立，主要利用斷層錯動資料計算特定點處的位移量。例如，斷層在j方向上錯動Δui，對整個矩形斷層范圍積分，即得到

式中，δjk為克羅內克函數變量；vk為∑面向外的垂直向量；μ、λ為拉梅常數；uji為受到破裂面在j方向的單位應力作用而在地表i方向產生的位移。

COMCOT模型中的斷層模型共輸入9個參數，分別為震中經度、緯度、震源深度、斷裂長度、斷裂寬度、滑動量、走向角、傾角、滑移角。斷層模型示意圖見圖2。

邊界條件的設定中，水邊界設為開邊界。選擇線性淺水波方程時，海陸邊界設為垂直反射邊界；選擇非線性淺水波方程時，采用移動邊界方案，即考慮岸灘的干濕變化。本文通過源碼修改在模型中加入潮汐強迫邊界條件（李林燕等，2012；趙鑫等，2016），外海潮位邊界由全球潮汐模型（TPXO7.2）求得，該模型通過10個分潮推算天文潮位，包含8個主要分潮M2、S2、K1、O1、N2、P1、K2、Q1以及2個長周期分潮Mf和Mn，基本能夠構造出外海深水處真實的天文潮過程

圖2 斷層示意圖

式中，ζ0為邊界處的潮位；ζp為邊界處靜壓水位；i為1～10，分別對應上述分潮；Ai、αi分別為分潮在3條邊界處的振幅和遲角；ωi為分潮的角頻率。

2 海嘯模型驗證

2.1 計算區(qū)域與網格布置

利用日本“311”海嘯對COMCOT模型進行驗證。模型采用2層嵌套網格，第1層網格的范圍是5°S～62°N、105.0°～179.5°E，計算范圍包括東海、臺灣海峽、南中國海、日本及部分西太平洋海域，網格大小為4′，控制方程選擇球坐標下的線性淺水長波方程，忽略摩擦作用；第2層網格的范圍是14.5°～28.0°N、110°～126°E，網格大小為0.8′（圖3）。

圖3 網格范圍

2.2 震源參數設置

日本“311”地震海嘯發(fā)生后，Tang等（2012）利用D21418浮標站在震后0.42～0.58h和D21401在0.80～1.23h之間的監(jiān)測數據，反演出初始海面位移有10m的抬升和3m的下沉，震源可由6塊沿著日本海溝的100km×50km的單位震源組合而成，在海嘯到岸前5h計算出了美國海岸線30多個城市的海嘯波爬高和淹沒情況，模擬值與后來的觀測值較為符合。Wei等（2013）也利用同樣的斷層設置，模擬計算了海嘯在日本近岸的爬高和淹沒，結果表明，模擬計算的茨城縣和青森縣的淹沒范圍準確率達到85.5%。以上述National Oceanic and Atmospheric Administration的斷層設置作為海嘯模型的輸入條件，模擬此次海嘯波的產生和傳播，并輸出在日本近岸的DART浮標點以及我國近海驗潮站的波高序列。地震斷層參數設計見表1。

表1 地震斷層參數

圖4 日本東部海嘯浮標分布

2.3 浮標站驗證

本研究收集海嘯發(fā)生后日本東部編號為21401＃（42.617°N，152.583°E）、21413＃（30.515°N，152.117°E）、21418＃（38.688°N，148.769°E）、21419＃（44.455°N，155.736°E）的海嘯浮標站水位資料進行驗證，浮標分布如圖4所示，模型驗證如圖5所示。在4個浮標點處，首波到岸時間的模擬值與觀測值之差均小于2min，21401＃和21419＃浮標站，波高的模擬值與觀測值之差小于2cm，浮標點處的模擬值與實測值較為符合；在21418＃h和21413＃浮標站，波高的計算值與觀測值的誤差分別為36、21cm，誤差相對較大，這與海底地形誤差、驗潮站局部海浪初值累積誤差等有關。

圖5 浮標點處模擬值與觀測值對比

3 計算結果分析

在我國南海海域附近，3大俯沖帶被認為具有地震激發(fā)海嘯的潛在危險，分別為馬尼拉俯沖帶、琉球群島俯沖帶和北蘇拉威西島俯沖帶。馬尼拉海溝位于亞歐大陸板塊和菲律賓海板塊交界處，南起菲律賓巴拉望島北端，沿菲律賓呂宋西部邊緣向北發(fā)展，北至臺灣島，總長度約1000km。亞歐大陸板塊以70mm／a的速度向菲律賓海板塊之下俯沖，俯沖板塊與上覆板塊之間的匯聚擠壓作用持續(xù)了相當長的一段時間，當集聚的應力釋放，便會發(fā)生海底地震，地震激發(fā)的海嘯波會使臺灣、福建、廣東、海南等沿海地區(qū)毫無遮攔地暴露在海嘯波的威脅之下。

本文通過使用COMCOT模型模擬了馬尼拉海溝地震后產生的海嘯對廣東省附近海域的影響，模型采用2層嵌套網格，網格范圍如圖5所示。美國地質調查局根據海溝方位角和斷層幾何學將馬尼拉海溝分為6個斷裂帶，震源位置如圖6所示，地震震級為MW9.3（Phuong et al，2013），震源斷層參數見表2。

圖7為馬尼拉海溝發(fā)生地震后海嘯波傳播和先導波到達時間圖。馬尼拉海溝因地震引發(fā)的海嘯波3h之后傳播到廣東省沿岸海域附近。COMCOT可以得到外海島嶼附近海域更精細的模擬結果，由COMCOT第2層模型結果可以得到廣東省附近海嘯引發(fā)的最大增水值。

圖8為馬尼拉海溝地震引發(fā)的海嘯對廣東省附近海域影響的最大增水圖。由圖8可知，馬尼拉海溝處地震引發(fā)的海嘯由深海區(qū)傳入淺水區(qū)后，在淺水效應的作用下，在東沙群島附近海嘯波有明顯的增高，最大增水達到11.95m，海嘯波傳到廣東省沿岸后最大增水約為4.87m。海嘯波的速度近似為（g H）1／2，當海嘯波與高潮位重合時，耦合計算總水深H增大，海嘯傳播速度增大，傳播時間變短。若不考慮潮汐、海嘯的計算結果，在時間上會有所延遲，約15min左右，沿海各地受地形的影響，會略有不同。表3為在馬尼拉海溝地震影響下廣東各海洋站增、減水值及對應時間。

圖6 馬尼拉海溝斷裂帶位置

表2 馬尼拉海溝斷層參數

為了更加詳細地了解海嘯波對廣東省沿岸的影響，選取廣東省沿岸的海門、甲子、汕尾、惠州港、香港、桂山島、廣州、大萬山、珠海港、海陵山島（閘坡）、湛江、下泊等處，得到其水位隨時間的變化序列，從而分析海嘯波對廣東省沿岸海域的影響。敏感點位置見圖9，受影響的時間序列見圖10。在各海洋站的統(tǒng)計中，東沙島和大萬山的最大增水相對較為顯著，分別為11.95、4.87m；其次是在珠海港和甲子附近，最大增水達到4m左右，惠州港附近3.76m，海門、汕尾、井岸、香港最大增水為2.0～2.5m，海陵山島（閘坡）、汕頭、茂石化（水東港）、東澳島、桂山島、西葛、潮州港（三百門）最大增水為1.5～2.0m，大鵬灣（鹽田港）、大亞灣、湛江、南澳島（云澳灣）、馬鞭洲（廣石化）、橫山、三灶島、燈籠山、流沙最大增水為1.0～1.5m；隨后較弱的為烏石、下泊、下港、博賀、電白、硵洲島（北港）、北津、海安，最大增水為0.5～1.0m；澳門、珠海（九洲港）、內伶仃島、黃埔、橫門、珠海（香洲）、北街、南沙（水牛頭）、廣州、海沁、深圳機場（油碼頭）、上川島（三洲灣）、舢舨洲、蛇口（赤灣）最大增水為0.19～0.50m，相對稍弱。由圖10可知，馬尼拉海溝地震引發(fā)的海嘯約2.5h后到達廣東沿岸，海嘯波到達敏感點區(qū)附近后，影響在1h內達到最大，隨后影響逐漸變小。內伶仃洋因外海星羅棋布的島嶼阻攔和淺水地形影響，受海嘯波的影響微弱，海嘯波高相對較小。不到5h，最大增水覆蓋了廣東大部分沿海地區(qū)，而對于廣州、湛江以西的影響則較為微弱。

圖7 馬尼拉海溝地震海嘯波傳播及先導波到達時間

圖8 馬尼拉海溝地震引發(fā)的海嘯對廣東省沿岸海域最大增水圖

表3 在馬尼拉海溝地震影響下廣東各海洋站增、減水值及對應時間

圖9 站位示意圖

圖10 馬尼拉海溝地震后典型站位水位時間序列

4 結語

本文基于COMCOT模型建立了廣東沿海天文潮與海嘯耦合數學模型，對日本“311”海嘯進行了良好地驗證，同時模擬了馬尼拉海溝地震后產生的海嘯對廣東附近海域的影響，對廣東沿岸由東向西的海洋站（海門至下泊等處）數據進行了分析，主要研究結果如下。

（1）本文建立的海嘯耦合數學模型能夠較好地模擬日本“311”海嘯，首波到岸時間的模擬值與觀測值之差均小于2min，波高的計算結果與實測值對比良好，兩者的差異與海底地形誤差、驗潮站局部海浪初值累積誤差等有關。本模型為準確評估馬尼拉海溝地震后產生的海嘯對廣東沿海的影響提供了研究基礎。

（2）馬尼拉海溝地震引發(fā)的海嘯約2.5h后到達廣東沿岸，海嘯波到達敏感點區(qū)附近后，影響在1h內達到最大，隨后影響逐漸變小。內伶仃洋因外海星羅棋布的島嶼阻攔和淺水地形影響，受海嘯波的影響微弱，海嘯波高相對較小。

（3）基本上5h以內，最大增水覆蓋了廣東大部分沿海地區(qū)，對近海工程將造成影響。東沙島和大萬山的最大增水相對較為顯著，分別為11.95、4.87m；其次是在珠海港和甲子附近最大增水達到4m左右，惠州港附近3.76m，海門、汕尾、井岸、香港最大增水為2.0～2.5m。而對于廣州，湛江以西的影響則較為微弱。

（4）由于南海海嘯發(fā)生的歷史記錄尚少，地震數據亦有所欠缺，震源機制解不能有效精確，故計算時間段無法準確判定。線性疊加潮汐則比耦合潮汐的海嘯波計算結果稍安全。同時，潮汐對海嘯的影響敘述篇幅太長，這也是后續(xù)的研究工作，海嘯波抵達時間并無多大差異，主要差異體現在海嘯最大波高，將來可考慮多年最大、最小潮差時間段，一年當中不同月份的最大、最小潮差時間段，以及受到東北季風和西南季風影響，甚至受到極端天氣影響情況下的模擬，以構建相應的海嘯源數據庫，建立更加完善的南海海嘯耦合數學模型。