高 峰, 朱川林, 李 昊, 唐友剛, 黃杰斯, 張慈珩, 劉 針, 姜云鵬
(1. 交通運輸部天津水運工程科學研究所,港口水工建筑技術國家工程實驗室,工程泥沙交通行業重點實驗室,天津 300456;2. 天津大學建筑工程學院,天津 300072; 3. 海南金海灣投資開發有限公司,海南 海口 570000)
中國是世界上海洋災害最嚴重、最頻發的少數國家之一,而四面環海的海南省又是最易受到海洋災害影響的地區,對于地處海南省經濟發展前沿的最大城市——海口市而言,其瀕臨海口灣水域所面臨的海洋災害形勢十分嚴峻,超強臺風、風暴潮、極端大浪、海嘯、海平面持續上升、岸灘侵蝕加劇以及海岸生態環境變化加速等,均對沿海經濟造成威脅相當嚴重,已成為沿海開發規劃實施和經濟建設的制約因素之一。項目用海災害風險分析主要是對項目災害性風險的種類、項目用海對海洋功能的風險、項目用海對周邊利益相關者的風險進行論證和分析[1],進而提出避免災害性風險發生的方案,或是災害性風險發生后的應急措施。因此,開展沿海工程海洋災害風險評價分析,是有效提高災害防御能力、控制風險和減少災害損失的關鍵技術方法之一,對該地區近海涉水工程的建設具有重要意義。
目前,一般風險分析的程序是:風險識別、風險估計、風險評價、風險處理和風險決策[2]。來自海洋水動力的風險主要包括兩個方面:(1)一方面是用海項目自身引起的突發或緩發事件導致對海域動力環境造成的危害,發生于施工期居多,如波浪增加、流速加快、懸浮物擴散、船舶溢油等;(2)另一方面是由于自然災害對海域使用項目造成的危害,發生于營運期居多,如遭受熱帶氣旋、風暴潮、海嘯等海洋災害對工程造成的損壞。海南省四面環海、瀕臨南海,地理位置特殊,海口灣位于海南島北部的海岸線上,分布于瓊北中部岬間海灣內,隔瓊州海峽與大陸的雷州半島相望,背依海口市及澄邁縣的馬村(老城)工業開發區。其所在海域的海洋水動力災害主要包括:風暴潮、海浪、海嘯等。按照一般風險分析的程序,就近岸海洋動力風險分析研究來說,重點在于上述災害的風險識別與估計,本文即通過實測資料分析結合數值模擬計算,確定上述各主要風險因子發生的強度大小和破壞情況。
登陸我國的臺風和強溫帶天氣過程往往造成風暴潮災害,其成災頻率高、致災強度大,造成的人員和經濟損失慘重。風暴潮災害主要是由異常的風暴增水使潮位大幅升高而形成,致災因子不僅包括風暴潮,還包括天文大潮、近岸浪及其三者之間的耦合作用。影響海南島北部的熱帶氣旋較多[3],是受風暴潮影響較為嚴重的區域之一,海水堆積引發的潮災是海口灣水域主要海洋災害之一。根據相關資料統計,海口年均發生風暴潮約3次,而嚴重或特大風暴潮約每2年1次[4]。海口灣每年6—12月均有可能發生風暴潮,且集中發生在7—10月,尤以10月份最多。在對歷史風暴潮資料進行搜集整理的基礎上,本文研究基于現有的實測臺風資料通過將Monte Carlo法[5](又稱為隨機統計模型)與水動力模型相結合的方法實現對可能出現的臺風風暴潮進行數值模擬研究,從而推求出重現期為50、100以及500 a的極端高水位預測值。
1.1.1 計算方法 采用假定水體上下層均勻一致的垂向平均二維風暴潮數值模型作為近海臺風風暴潮的水動力模型,且二維模型突出了風暴潮預報與計算的最重要變量—風暴潮位,這也是預測臺風風暴潮引起的極端高水位的核心變量。模型采用大小模型嵌套方式,其中大模型的北側和西側以大陸岸線為陸地邊界,南面以北緯15.0°為外海邊界,東面以東經117.9°為外海邊界,模型范圍覆蓋了北部灣、瓊州海峽、瓊州海峽以東及南海部分海域,這種規模的大尺度計算網格可以很好地對一次臺風影響過程進行數值模擬,其網格尺寸約為1 km,共有420×808個網格點,339 360個單元。同時,小模型南北兩側以大陸岸線為陸地邊界,南北向跨度最大處約180.0 km,東西向從108.7°E~110.9°E,跨度約為229.2 km。計算范圍及網格劃分如圖1(a)~(b)所示。

((a)大模型; (b)小模型; (c)1409號臺風模擬; (d)驗證結果。(a)Large scale model; (b) Local nested model; (c)Typhoon 1409 simulation; (d)Verification results. )
圖1 模型網格嵌套分布與1409號臺風模擬驗證情況
Fig.1 Distribution of model nested grid and verification result of No.1409 typhoon’s simulation
氣壓分布模型采用Holland模型,環流風速模型采用梯度風速公式,移行風速模型采用宮崎正衛模型。由氣壓分布模型和風速模型組合成臺風風場的模擬方案,簡稱Holl模型[6]。典型的實測臺風依據美國FEMA使用的臺風歷史記錄數據處理方法—模擬圓法,提取在一定時期內途徑以要模擬的城市為模擬點,以設定的距離為半徑的圓的所有臺風,以這些臺風的記錄作為研究對象[7-8]。在1949—2015年期間曾有151個臺風經過以海口灣為中心,200 km為半徑的圓形研究區域,篩選登陸點中心最低氣壓在0~980 hPa范圍內共計60個。在所有的60個實測臺風中,2014年7月18日經過的“威馬遜”臺風強度最大,故利用鐵山港海區英羅灣附近氣象觀測站實測資料對該場臺風進行了風場驗證,“威馬遜”登陸后整個研究區域于2014年7月18日18:00的氣壓場與風場模擬情況可分別見下圖1(c)所示,模擬與實測的最大風速擬合良好圖1(d)。進而采用Monte Carlo模型生成100年內可能出現的模擬臺風數據。由于臺風歷史資料采用的是1949—2015年間的60個典型臺風數據,則500年內應當產生的臺風數為448個。對這些個可能出現臺風的參數進行模擬,可得所需的模擬臺風數據。采用隨機模型(Monte Carlo)隨機模型和水動力模型(Delft3D)相結合的方法得到水位過程線,并得出極值水位。首先采用Monte Carlo隨機模型,將篩選出的67年里60條臺風的氣壓和路徑分別作為獨立的元素隨機組合得到448個重組臺風。將重組臺風要素帶入已驗證的風暴潮水動力模型進行計算,其中臺風過程與天文潮位隨機組合,最后得到448個水位過程線,統計每條水位過程線中的最大水位值,用于后面的極端水位分析。
1.1.2 計算結果 通過以上過程,統計每條水位過程線中的最大水位值,部分較大水位過程線如圖2(a)所示。將統計出的448個最大水位值進行降序排列,計算得秀英站處最大水位序列的最高值達3.75 m。極端水位的回歸周期TR可按水位的排序參照式:TR=500/M估計,其中,M=8是448個最高水位的從高到低的排序,TR為重現期。例如,當M=1(對應的是448個極端水位的最大值)時,這個水位的最大值對的重現期是TR1=500/1=500 a。相應的極端水位所對應的重現期可點繪成圖,且可由該圖確定出任一重現期所對應的極端水位值。第五高水位3.02 m所對應的是100年一遇的極端高水位值,重現期為500、100和50年所對應的秀英站極端高水位預測值分別為3.75、3.02和2.68 m。
另外,通過采用重現期100年一遇氣壓P0=883 hPa,最大風速半徑R為33 km,臺風移速為35 km/h,采用假想路徑,對可能最大風暴潮增水(PMSS-ProbableMaximum Storm Surge)[9]進行了計算,計算結果表明最大增水達到4.42 m,其對應增水過程見圖2(b)。

((a)水位過程線; (b)最大可能增水過程。(a)Stage hydrograph; (b)Maximum possible water increase process.)
圖2 各重現期秀英站臺風暴潮水位過程線和最大可能風暴增水過程圖
Fig.2 Storm surge water level process of different return period for Xiuying and the greatest possible storm surge elevation
海口灣屬臺風作用活躍區域,受臺風和大陸季風作用明顯,但該水域及周邊缺乏長期觀測站,現僅有玉苞角、白沙門兩個臨時站。影響該海域的臺風主要來自南海自東或東南經廣東、海南登陸或從海南島南部海域向偏西方向移動的臺風。海口灣水域全年涌浪較少、以風浪為主,風浪頻率為76%~85%,涌浪頻率占14%~23%,常浪向為ENE,頻率為30.1%,次常浪向為NE,頻率為22.9%。波浪出現最少的方向為S~WSW。受季風風向影響,工程海域常浪向隨季節風而變化,東北季風期(11 月—次年3 月),常浪向為ENE,次常浪向為NE,季風轉換期(4—5 月,9—10 月)常浪向、次常浪向與東北季風期一樣,但常浪向出現頻率略低于冬季,次常浪向出現頻率略高于冬季;西南季風期(6—8 月)常浪向為NE。目前已知最大波浪由1409 號威馬遜臺風超強臺風產生,該超強臺風在海口灣的東北側通過,極有利于西北向波浪的成長和傳入海口灣,在海口灣出現最大波浪,形成有效波高6.8 m,周期約7.8 s,該波高為偏西波向的異常特大值,在N-E-S方位的最大有效波高3.9 m,周期約6.3 s。
1.2.1 計算方法 該海域主要成風原因有兩種,一為季風、二為氣旋。因此,研究主要通過風場模擬計算波浪場,風場資料需綜合考慮季風和臺風兩種影響。具體研究思路如下:
(1)采用中國氣象數據網(www.data.cma.cn)公布的臺風路徑,結合歐洲中長期在南海海域的后報風場資料,選取每年影響工程海區的臺風及季風過程,對風場進行分析。
(2)根據多年風場數據,采用SWAN模式建立大區域波浪計算模型,得到每年工程區附近外海深水區的年波浪極值,進而采用Pearson-Ⅲ型曲線適線法推算不同重現期波要素。
(3)根據外海不同重現期波浪條件,在自然地形條件下即工程建設前,利用SWAN模式建立小區域波浪計算模型,推算工程區設計波浪條件。
1.2.2 主要結果與分析 根據工程海區最近30年(1986—2015年)臺風和季風資料,通過臺風場模型和季風場數據分析計算工程海區大風過程,進而采用大范圍風浪模型計算工程海區風浪場,從中選取主要計算浪向ENE、NE、NNE、N、NNW與NW向的波浪年極值,分析工程外海深水區不同重現期波浪。提取海口灣北側外海計算點(水深約-70 m)處的主要計算方向波浪重現期結果見下表1,其中前文所提到的60個典型臺風中,分別選取1522號“啟德”臺風與1117號“飛燕”臺風過境瓊州海峽時波高分布圖(見圖3),其中重現期100年有效波高為7.14 m,重現期50年有效波高為6.80 m。

表1 海口灣不同重現期有效波高結果Table 1 Results of different return period of significant wave height in Haikou-Bay /m
另外,依據上述波高結果,通過HUDSON公式對海口灣近岸工程可能采用斜坡式護岸或者防波堤的護面塊體重量進行了核算,如目前海口灣內各人工島如普遍采用30 t護面塊體時,其對應破壞波高可達8.78 m,當于重現期為1 400年波浪(年概率為0.071%),若以設計年限100年考慮,這一等級的塊體穩定性保障率可達93%。同時,根據爬高和越浪的計算分析,極端高水位條件下有效波高超過7.9 m就會出現危害堤頂鋪面結構安全的越浪量,該波浪條件對應的重現期約為400年,既年出現頻率為0.25%。

圖3 典型臺風工程海域作用過程Fig.3 Functional process of typical typhoon in engineering sea area
隨著波浪由外海向近岸傳播,充分考慮了灣內自然水深地形條件下波浪分布。采用SWAN模型對工程區的波浪分布情況進行計算分析。計算方向為ENE、NE、NNE、N、NNW與NW向。包含100年一遇高水位、設計高水位與設計低水位。其中,ENE最大有效波高為6.7 m(見圖4(a)),NE向最大有效波高5.8 m(見圖4(b)),NNE向為5.5 m,N向為3.8 m,NNW向為4.4 m,NW向為4.9 m。對海口灣近岸影響較大的2個方向為ENE與NE向,重現期100年一遇高水位時H13%波高分布見圖4所示。

圖4 100年一遇高水位重現期100年有效波高分布結果Fig.4 Distribution of significant wave height for 100 years return period by a high water of 100 years return period
依據美國地質調查局(USGS)提供的一份潛在海嘯地震斷裂帶風險圖[10],我國南海海域附近有三大俯沖帶被認為具有地震激發海嘯的潛在危險,分別為馬尼拉俯沖帶、琉球群島俯沖帶和北蘇拉威西島俯沖帶。其中,馬尼拉俯沖帶(又稱為馬尼拉海溝)對南海周圍的國家和地區會帶來很大的海嘯威脅,其位于亞歐大陸板塊和菲律賓海板塊交界處,南起菲律賓巴拉望島北端,沿菲律賓呂宋西部邊緣向北發展,北至臺灣島,總長度約1 000 km[11]。USGS根據馬尼拉海溝的方位特征和斷層幾何學將其劃分為6個斷裂帶[12],沿馬尼拉海溝分布如圖5(b)所示。
上述區域地震激發的海嘯波會使海南等我國南部海域暴露在海嘯波威脅之下,雖然歷史統計資料顯示幾次大的跨洋海嘯對我國東南沿海地區影響有限,但南海馬尼拉海溝地震帶形成的海嘯仍會對臺灣、海南、廣東和福建等沿海地區曾造成不同程度的影響,南海曾發生過約15次局地海嘯,造成了一定影響[13]。由于近幾年全球地質運動較為活躍,特大地震海嘯發生頻率較高,馬尼拉海溝被認為是最有可能引發海嘯的潛在地震源。基于上述分析,仍不可忽視海嘯波的潛在威脅。
2.3.1 計算方法 本文采用淺水波數值模型GeoClaw,基于有限體積法對非線性淺水方程進行數值離散,考慮了海嘯波在近岸傳播的非線性作用,同時嵌入了自適應網格技術可追蹤海嘯波高判斷網格加密與否。通過對馬尼拉海溝斷層特征的分析,Lau等[14]根據歷史上曾經在附近區域發生過的地震和經驗關系公式,確定了7.0、8.0級地震時E1~E3板塊的數值模擬參數(見表2),并基于Nguyen等[15]在USGS震源參數基礎上給出的震級為9.3級的震源參數,此時對應各板塊震源參數見表3,從而依此計算相應的傳播對海南島周邊海域的影響。

( (a)南海附近俯沖帶;(b)馬尼拉海溝不同板塊分布。(a)Subduction zone near the south China sea。(b)Distribution of different plates in the Manila trench.)
圖5 中國南海附近的俯沖帶與馬尼拉海溝不同板塊產生的海底變形及位置
Fig.5 Sea bottom deformation and location of subduction zone and different tectonic plates from the Manila trench nearby China South Sea
利用Okada彈性半空間錯位理論模型,根據前文所述六塊斷裂帶參數生成相應的海底變形可分析出海南島受E1、E2和E3板塊的影響較大。據此將對上述各板塊(震源深度15 km)分別發生7.0級和8.0級地震形成的海嘯進行數值模擬,進而估算對海口灣水域的影響情況,還考慮了整個斷裂帶(E1~E6)發生破裂形成的特大地震海嘯時的極端情況。計算范圍從99°E~130°E,1°S~33°N,包含了我國東南沿海至南海海域以及菲律賓、越南等國家。

表2 震級分別在7.8和8.0時馬尼拉海溝E1~E3斷裂帶震源模擬參數Table 2 Source modeling parameters of E1 ~ E3 fault zone in Mw=7.8 and 8.0 respectively in the Manila trench

表3 Mw=9.3時,馬尼拉海溝E1~E6斷裂帶基本參數Table 3 Source modeling parameters of E1 ~ E6 fault zone in Mw=9.3 in the Manila trench

表4 各計算工況下工程區海嘯波特征值 Table 4 Characteristic value of tsunami under the different calculation conditions
2.3.2 結果與分析 數值模擬結果如圖6所示,分別給出了馬尼拉海溝北部3個板塊(E1、E2、E3)分別發生里氏7.0、8.0級海底地震形成的海嘯波。各工況下工程區的海嘯波特征值如表4所示。從走時圖中可以看到,海嘯波大概經歷4~5 h即可到達海口灣附近海域。對于單一板塊的一般震級地震,海嘯波到達海口附近波幅都小于4 cm,波長在50 km左右,周期1.6~2 h,海嘯持續時間(第一個大波的作用時間)約為0.7~1.1 h,海嘯波行進速度為8.1 m/s。對于整個斷裂帶發生9.3級大地震引發的極端海嘯,海嘯波到達工程區附近海域時最大波幅為2.25 m,波長115 km,周期3.4 h,海嘯持續時間1.3 h,海嘯波行進速度為9.4 m/s,該情況將對海口灣會產生影響,特別是岬灣形態的海口灣近岸,更易形成能量集中,應予以關注。

圖6 各板塊分別在7.0和8.0級地震生成的海嘯波走勢圖Fig.6 The tsunami wave distribution of each plate in Mw=7.0 and 8.0 generated by earth quake
海洋水動力風險是用海風險的重要組成部分,關鍵在于風險識別和風險估測。本文根據對海口灣水域可能出現的風暴潮、極端波浪、海嘯等問題的研究,對該海域水動力風險因素進行了研究分析,得到風暴潮不同重現期增水水位,并在確定外海重現期100 a有效波高為7.14 m,重現期50 a有效波高為6.8 m的基礎上,計算了海口灣近岸不同波向的最大波高,對灣內近岸斜坡式防護工程采用護面塊體的重量校核。針對潛在海嘯威脅最大的馬尼拉海溝,進行了多種工況數值模擬,除9.3級大地震引發極端海嘯導致海嘯波達2.25 m之外,7、8級海嘯對于海口灣所受影響較小。這些計算結果反映海口灣沿岸在不同水位下的漫灘情況,對海口灣沿岸重要的工程項目以及防災規劃的制定,具有參考意義。