北部湾台风风暴潮数值模拟及重点区域风险分析

蒋昌波，赵兵兵，邓 斌，伍志元

（1.长沙理工大学水利工程学院，湖南长沙 410004；2.水科学与水灾害防治湖南省重点实验室，湖南长沙 410004）

北部湾台风风暴潮数值模拟及重点区域风险分析

蒋昌波1,2，赵兵兵1，邓 斌1,2，伍志元1

（1.长沙理工大学水利工程学院，湖南长沙 410004；2.水科学与水灾害防治湖南省重点实验室，湖南长沙 410004）

基于Delft3D模型建立了适用于北部湾海域的台风风暴潮数学模型，同时根据进入北部湾海域67 a（1949—2015年）历史台风资料，采用蒙特卡罗方法随机构造55场台风进行风暴增水计算，并选取受风暴潮影响最严重的铁山港与最频繁的涠洲岛进行了风险分析。结果表明：（1）从铁山湾下侧进入北部湾海域，中心气压为952 hPa的台风会造成石头埠站风暴增水为3.68 m，达到1 000 a一遇的级别；（2）广西沿岸的风暴增水对经过北部湾海域的台风存在滞后效应；（3）涠洲岛风暴增水普遍较低，但进入北部湾的台风均会在涠洲岛产生风暴潮，因此要做好涠洲岛的防护工程以防止海岸线继续后退。

Delft3D；风暴潮；蒙特卡罗；风险分析

1 引言

台风是沿海地区主要的自然灾害，当风暴潮与天文大潮高潮叠加时可能产生极端高水位，导致海岸防御建筑物的破坏、局部地区出现沿岸洪水、岸滩侵蚀与海岸线后退，造成大量人员伤亡与财产损失。据1986—2010年数据统计，北部湾广西沿海风暴潮灾害造成的直接经济损失达94.70亿元，受灾人数1053.73万人，冲毁海岸工程476.57 km[1]。为使北部湾海域的经济得到可持续发展，开展数值模拟工作及重点地区的风险分析显得尤为重要。

国外风暴潮数值模拟研究起始于20世纪50年代，主要通过对二维流体动力学方程积分计算出风暴增水的极值。近年来，国外对风暴潮的数值模拟研究主要集中在台风特性对风暴增水的影响[2-3]、风暴潮与天文潮的耦合以及波浪与风暴潮的耦合[4]3个方面，均取得了一定的进展。国内对风暴潮的研究开始较晚，20世纪70年代，冯士筰[5]从封闭、半封闭海域以及开敞海域的风暴潮问题出发，建立了我国风暴潮研究的理论体系，为我国风暴潮研究奠定了基础。李树华等[6-7]对广西沿岸主要港口风暴潮进行了研究，初步建立了各港口台风风暴潮的预报模型。陈波等[8-9]对北部湾海域风暴潮进行了整体性研究，包括风暴潮的形成与地形、台风路径的关系以及增减水的分布规律等。上述研究多关注于不同条件下风暴增水的计算，对北部湾区域进行极端情况的风险分析相对较少。

本文首先基于Delft3D水动力模型以及Holland风场模型建立了适用于北部湾海域的风暴潮数学模型，随后依据1949—2015年间北部湾海域的台风发展规律建立了台风数学模型，最后，对模拟台风产生的风暴潮进行计算，并对铁山港与涠洲岛风暴增水情况进行风险分析，为广西沿岸的风暴潮预报及防灾减灾工作提供帮助。

2 风暴增水数值模拟计算

2.1 水动力模型

水动力模型采用Delft3D，数值方式采用ADI法。Delft3D模型计算台风风暴潮在国际上得到广泛认可[10-11]。

在球坐标系中，连续性方程为：

式中：ζ是参考平面以上自由水面，d是参考平面以下水深，λ是经度，ϕ是纬度，R是地球半径（6 378.137 km，WGS84），U、V是水深平均流速，Q代表源和汇的作用，如降水、蒸发等：

式中：H=d+ζ，P是降雨，E是蒸发，qin是水的来源，qout是水的渗流。

λ、ϕ方向动量守恒方程为：

式中：Pλ、Pϕ是压强梯度，Fλ、Fϕ是雷诺应力项，f是科氏力系数，Mλ、Mϕ分别是外部源和汇的分量。

2.2 风场和气压场

风场和气压场采用Holland台风模型[12]：

距离台风中心r处的风速V(r)：

式中：Pn是外围表面气压，Pc是台风中心气压，Rmax是最大风速半径，ρa是空气密度，B是台风形状参数。

2.3 天文潮与风暴潮的耦合

图1 计算区域及网格

对天文潮和风暴潮的耦合，是在计算区域内和深水开边界处，即要考虑天文潮的作用，又要考虑台风的作用。天文潮的作用是通过在开边界上提供M2、N2、S2、K2、K1、O1、P1、Q1等8个主要分潮作为驱动力。台风的作用是通过上述的风场和气压场提供计算区域内和开边界上的风和气压，按照静压假设和自由表面边界条件作用在水体上。

2.4 模型设置与验证

为提高模型精度，本研究采用嵌套网格进行计算。网格一的计算范围为：105.5°～117.8°E，16.0°～23.9°N；网格二的计算范围为北部湾广西沿海。两套网格均为正交曲线网格（见图1）。铁山港海域水深由数字化航海保证部海图得到，其他位置采用GEBCO_14的30"水深数据。外海潮汐调和常数由全球潮流模型TPXO7.2提供。验证模型采用1409号“威马逊”台风，计算时间为2014年7月10日00时（世界时，下同）—2014年7月20日18时。

石头埠站、临时站风暴增水验证见图2。石头埠站最大增水比实测小19.1 cm，最大减水比实测大43.1 cm；临时测站最大增水比实测小21.1 cm，最大减水比实测大28.5 cm，两测站的风暴增减水趋势与实测较为一致，证明该风暴潮数值模型准确合理。

3 台风模型的建立

3.1 北部湾海域台风统计

进入北部湾海域的台风分3种类型[9]：一是斜穿雷州半岛和海南岛东北部进入北部湾，在广西沿海或越南北部登陆，该类台风引起广西沿海港湾强烈的增水；二是横穿海南岛或雷州半岛进入北部湾，在越南北部沿海登陆，该类台风引起的港湾增水程度和范围要小于一类路径；三是绕过海南岛向北发展，在广西沿岸登陆，该类台风也会引起广西沿海港湾的水位变化。进入北部湾海域的台风路径类型如图3所示。本文统计了1949—2015年进入北部湾海域的共100场台风（见表1）。

由于Ⅰ类台风进入北部湾的频率最高，并且Ⅰ类台风在北部湾海域内会产生大于Ⅱ、Ⅲ类台风增水，因此下面对Ⅰ类台风进行模拟。

3.2 Ⅰ类台风生成方法

文中对62场Ⅰ类台风的台风路径（经纬度）、最大风速以及台风中心气压资料进行统计整理，台风资料来源于中国气象局热带气旋资料中心，随后依据蒙特卡罗方法随机构造进入北部湾海域的台风。蒙特卡罗方法是一种纯概率统计方法，依据此方法构造台风的真实性与历史台风样本的大小相关，同时，此方法生成的台风不能反映出地形对台风强度的影响，如台风经过雷州半岛或海南岛时强度会有所减弱。但蒙塔卡罗方法具有计算快，可同时计算多个台风方案等优点，在增加路径转折点、风速/气压转折点以及风速/气压最大值等限定条件后，可较好的反应进入北部湾海域台风的特征，因此本文采用蒙特卡罗方法构造台风。

表1 进入北部湾台风类型统计

图2 风暴潮增、减水验证

3.2.1 台风起点的生成

台风经度起点在127.5°～145.5°E间随机生成，纬度起点在8°～10°N间随机生成。

3.2.2 台风路径的生成

统计历史数据，Ⅰ类台风到达路径转折点（120°～125°E）会发生转折，偏向西方前进（见图3）。转折前、后台风单位统计时间（6 h）经、纬度前进距离分布见图4和5。

图3 进入北部湾海域的台风路径分类示意图

转折前经度前进距离分布：

转折后经度前进距离分布：

转折前纬度前进距离分布：

转折后纬度前进距离分布：

图4 单位统计时间内台风前进经度分布拟合

图5 单位统计时间内台风前进纬度分布拟合（——为高斯分布拟合曲线）

在路径转折点前，由式（7）、（9）分别生成单位统计时间的经、纬度前进距离；一旦到达路径转折点，由式（8）、（10）分别生成单位统计时间的经、纬度前进距离。

3.2.3 台风风速与中心气压的统计与生成

台风起点风速在10～20 m/s之间，起点中心气压在990～1 005 hPa之间。去除掉极端情况，在台风到达风速/气压转折点（108.5°～111.5°E）前，最大风速在0～5 m/s之间增大，当风速达到70 m/s时，则会维持在70 m/s不变；台风中心气压在0～5 hPa之间减小，当中心气压达到920 hPa时，则维持在920 hPa不变。在台风行进至风速/气压转折点后，台风强度降低，最大风速在0～10 m/s之间减小，中心气压在0～10 hPa之间增大。最大风速变化统计见表2，中心气压变化统计见表3。

3.2.4 台风最大风速半径

由于缺乏实测资料，台风最大风速半径是最难确定的参数。而台风最大风速半径与台风中心最大风速有密切关系，因此本文在计算过程中，采用经验公式进行计算[13]：

式中：Rmax是最大风速半径；Vmax是近中心最大风速；ϕ为纬度。

3.2.5 台风终点的生成

为避免台风不切实际的深入内陆，规定当台风风速小于15 m/s或路径经度小于经度终止点（102°～106°E）时，台风终止。

表2 最大风速变化统计

3.3 模拟台风结果

I类台风多发生于7月、8月，规定本文中的台风起点在7月15日00时生成，运行程序，生成的55场台风见图6，通过与实测台风资料对比，排除掉5场风速与气压拟合效果不佳的台风，模拟结果能够反应进入北部湾海域的台风特征。

4 结果与分析

历史上北部湾海域各个测站的增水最大值中，铁山港内石头埠站的值最大，为2.71 m（1409号台风“威马逊”），涠洲站的值最小，为1.78 m（0312号台风“科罗旺”）。因此，本文对石头埠站与涠洲站的风暴增水情况进行分析，模拟台风风暴增水频次见图7。

与历史资料相比，石头埠站产生2 m以上风暴增水的频率较大，这是由于模拟台风均为I类台风。其中增水在0～1 m之间发生28次，1～2 m之间发生15次，2～3 m之间发生6次，增水在3 m以上的发生1次。石头埠站产生增水在2 m以上的7场台风均斜穿雷州半岛，最大风速在40 m/s以上，中心气压在945～960 hPa之间。产生增水在1.5～2 m之间的台风有两种类型，第一种类型台风，路径与产生2 m以上增水台风路径类似，最大风速在30～40 m/s之间，中心气压在965～980 hPa之间；第二种类型台风，路径在涠洲岛至琼州海峡之间，最大风速在50 m/s以上，中心气压在935 hPa以下，此类强台风在历史上并未发生过，但铁山港海域对此类台风仍应加强防范。产生增水在1～1.5 m之间的台风多数经过海南岛中部，最大风速达40 m/s，中心气压低于965 hPa。历史上此路径台风，风速多在18～30 m/s之间，达到40 m/s的较少，仅有3场台风（9106号“Zeke”、6403号“Winnie”、5513号“Kate”）。因此通过对台风特征参数的分析，可以初步判断进入北部湾的台风产生的风暴潮对铁山港的危险性。

表3 中心气压变化统计

图6 模拟台风路径

图7 模拟台风最大增水频次

纪燕新[14]采用耿贝尔方法，推算出石头埠站重现期是1 000 a的增水为3.49 m。本次模拟中石头埠站的天文潮历时曲线和36号台风引起的风暴增减水曲线见图8，最大增水为3.68 m，达到了千年一遇的级别。此次台风路径斜穿雷州半岛，于2014年7月21日00时登陆雷州半岛，台风风速达到45 m/s，中心气压为952 hPa，于7月21日06时进入铁山湾，最大风速及中心气压不变，按照《热带气旋等级》[15]中的划分，达到了强台风级别。

从图9可以看出，路径在北部湾海域的台风（以36号台风为例）产生的风暴增水在广西沿岸存在滞后效应。以铁山港为例，当台风中心未到达铁山港时，铁山港海域最大风速的风向从港内指向外海，即离岸风，所以首先出现一个较强的减水过程；当台风中心行进至铁山港时，最大风速的风向与港湾方向垂直，产生的增水与减水较小；当台风中心穿越铁山港、行进至廉州湾海域时，铁山港海域处于台风的右半圆区域，最大风速的风向径直指向港内，即向岸风，海水在强风作用下进入港区内，所以产生较强的风暴增水。因此，应注意这类台风过境后的港区内增水情况，做好灾害防护工作。

图8 石头埠站天文潮位及风暴增减水示意图

铁山港由于其地理形状近似口袋型，水体容易堆积且不易流出，从而产生较大的增水，而涠洲岛是一个孤岛，四周相对开阔，因此产生的风暴增水高度较低。据统计，1956—2014年，涠洲站增水＞0.5 m的风暴潮共发生83次，＞1 m的风暴潮共发生3次[16]。

图9 36号台风经过北部湾时的风场及水位分布

本文模拟的台风造成涠洲站最大增水为1.40 m，其路径与0312号台风“科罗旺”类似，均为从涠洲岛下侧穿越北部湾，最后在越南登陆并消失。增水在0.4～0.9 m的台风共发生41次，占本次模拟的82%。虽然风暴增水较低，但经过北部湾海域的台风均会在涠洲岛产生风暴增水。涠洲岛的潮间带为沙质岸滩，台风过境产生的风暴潮及近岸浪会导致较为严重的岸滩侵蚀与海岸线后退，2006—2013年间，涠洲岛西南部岸滩年均下蚀达0.18 m，在2013年6月风暴潮过后，北部岸滩最大下蚀达0.4 m[17]。因此需在岛上易受侵蚀区域种植固沙植被，建造防护性工程以阻止沙滩的后退。

5 结论

本文的台风模型是在统计北部湾海域历史台风数据基础上建立的，目的是为了生成产生潜在性灾难的台风，而并不是历史上真实的台风。通过对该台风模型与风暴潮模型的结合应用，对铁山港与涠洲岛的风暴增水情况进行了风险分析，得出了以下结论：

（1）基于Delft3D水动力模型和Holland台风模型建立的二维数学模型可以较好的模拟北部湾海域的天文潮与风暴潮耦合的情况；

（2）通过对台风特征参数的分析，可以初步判断铁山港海域风暴潮的危险性，为防灾减灾工作提供帮助；

（3）广西沿岸风暴增水对北部湾海域的台风存在滞后效应，即台风过境前会在沿岸产生减水，台风过境后产生增水；

（4）涠洲岛风暴增水较低，但为防止海岸线继续后退，需要采用适当的防护工程，如人工补沙、栽种固沙植被等。

[1]陈宪云,刘晖,董德信,等.广西主要海洋灾害风险分析[J].广西科学,2013,20(3):248-253.

[2]Irish J L,Resio D T,Ratcliff J J.The influence of storm size on hurricane surge[J].Journal of Physical Oceanography,2008,38(9): 2003-2013.

[3]Lima R J,Li C Y.On the importance of the forward speed of hurricanes in storm surge forecasting:A numerical study[J]. Geophysical Research Letters,2009,36(7):L07609.

[4]Xie L,Wu K J,Pietrafesa L,et al.A numerical study of wavecurrent interaction through surface and bottom stresses:Winddriven circulation in the South Atlantic Bight under uniform winds [J].Journal of Geophysical Research:Oceans,2001,106(C8): 16841-16855.

[5]冯士筰.风暴潮导论[M].北京:科学出版社,1982.

[6]李树华,陈文广,陈波,等.广西沿海台风暴潮数值模拟试验[J].海洋学报,1992,14(5):15-25.

[7]李树华.珍珠港台风暴潮特征及其预报的初步研究[J].海洋预报,1986,3(4):17-24.

[8]陈波,魏更生.广西沿海风暴潮的数值计算研究[J].海洋湖沼通报,2002(1):1-8.

[9]陈波,邱绍芳.广西沿海港湾风暴潮增减水与台风路径和地形效应的关系[J].广西科学,2000,7(4):282-285.

[10]Gomes M P,Pinho J L,Do Carmo J S A,et al.Hazard assessment of storm events for The Battery,New York[J].Ocean&Coastal Management,2015,118:22-31.

[11]黄潘阳,叶银灿,韦雁机,等.“威马逊”台风暴潮增水及水动力响应数值模拟[J].海洋预报,2012,29(2):32-38.

[12]Deltares.Wind Enhanced Scheme(WES)for Cyclone Modelling [M].2014.

[13]刘永玲,冯建龙,江文胜,等.热带气旋资料长度对风暴潮危险性评估结果的影响[J].海洋学报,2016,38(3):60-70.

[14]纪燕新.北部湾广西沿海风暴潮灾害及防灾减灾研究[D].南宁:广西大学,2007.

[15]中国国际标准化管理委员会.GB/T 19201-2006热带气旋等级[S].北京:中国标准出版社,2006.

[16]李明杰,吴少华,刘秋兴,等.风暴潮、大潮对广西涠洲岛西南沙滩侵蚀的影响分析[J].海洋学报,2015,37(9):126-137.

[17]姚子恒,高伟,高珊,等.广西北海涠洲岛海岸侵蚀特征[J].海岸工程,2013,32(4):31-40.

Numerical simulation of typhoon storm surge in the Beibu Gulf and hazardous analysis at key areas

JIANG Chang-bo1,2，ZHAO Bing-bing1，DENG Bin1,2，WU Zhi-yuan1
（1.School of Hydraulic Engineering,Changsha University of Science&Technology,Changsha 410004 China;2.Hu nan Province Key Laboratory of Water,Sediment Sciences&Flood Hazard Prevention,Changsha 410004 China）

A mathematic model of typhoon storm surge for the Beibu Gulf was established based on Delft3D hydrodynamic models,and the Monte Carlo method was used to generate 55 typhoons randomly according to the historical typhoon data of 67 a(1949—2015)in the Beibu Gulf.A study on hazard research of the most severe area(Tieshan Bay)and the most frequent area(Weizhou Island)affected by storm surge has been carried out by storm surge calculation with the generated typhoons.The result showed that(1)the surge at Shitoubu station affected by a typhoon with central pressure of 952 hPa is 3.68 m which return period is 1 000 a.(2)The storm surge at the coast of Guangxi has lagged effect on the typhoon passing through the Beibu Gulf.(3)The storm surge at Weizhou Island is generally low,but the typhoon into the Beibu Gulf will generate a storm surge here. Hence,the protection project at Weizhou Island should be prepared to prevent the coastline continues to retreat.

Delft3D;storm surge;Monte Carlo;hazard research

P731.23

A

1003-0239（2017）03-0032-09

10.11737/j.issn.1003-0239.2017.03.005

2016-12-11；

2017-02-10。

国家自然科学基金资助项目（51239001）；交通运输部应用基础研究项目（2015319825080）；湖南省研究生科研创新项目（CX2015B348）。

蒋昌波（1970-），男，教授，博士，主要从事河流、海岸动力过程及其模拟技术研究。E-mail:jcb36@vip.163.com