福建沿海浪潮耦合漫堤风险评估：以台风天兔为例*

王 凯 侯一筠,5① 冯兴如 李水清 傅赐福

(1.中国科学院海洋研究所 青岛 266071；2.中国科学院海洋环流与波动重点实验室 青岛 266071；3.中国科学院大学北京 100049；4.中国科学院海洋大科学研究中心 青岛 266071；5.青岛海洋科学与技术试点国家实验室 青岛 266237；6.国家海洋环境预报中心 北京 100081)

福建沿海位于西北太平洋台风的主要移动路径上,台风及其引起的风暴潮灾害是福建沿海遭受的主要的海洋灾害之一。风暴潮灾害是由于台风风暴潮、天文大潮和台风浪等共同作用引起的沿岸水位异常上升造成的。根据国家海洋局2018年统计公报,2018年,风暴潮灾害在各类海洋灾害中导致的直接经济损失最严重,占直接经济损失总额的93%；海浪灾害导致的死亡人数(含失踪)最多,占总死亡人数(含失踪)的96%。其中,风暴潮灾害给福建造成的直接经济损失为11.41亿元,占风暴潮灾害直接经济损失总额的25.6%。海浪灾害对福建造成的直接经济损失为0.13亿元,占海浪灾害直接经济损失总额的37.1%；并且海浪灾害导致福建29人死亡(含失踪),占海浪灾害总死亡(含失踪)人数的41.4%。

海堤是在河口和沿海地区修建的一种专门用来挡水的建筑物,其目的是为了防止天文大潮的高潮位和台风风暴潮以及台风浪的侵袭,保障沿海城市的工业设施和农田的安全。近年来,为了提高沿海城市的防洪抗灾能力,不少沿海城市都修建了不同等级的海堤。海堤在预防和减少海洋灾害方面发挥了重要的作用。但是由于海水的侵蚀,且有些海堤年久失修,当沿海地区遭遇的天文大潮和超强台风极易引发漫堤甚至溃堤灾害,给人民的生命和财产带来重大损失。因此,有必要研究风暴潮期间漫堤风险等级评估,以分析其风险性和危害程度。及时并且准确的漫堤风险等级评估会把人员伤亡和经济损失减小到最低程度。

漫堤风险等级评估一般以堤前水位是否超过堤顶高程来判断。20世纪50年代,国内外学者在物理模型试验和现场试验的基础上,提出了一系列越浪和波浪爬高的计算公式。国内许多学者对风暴潮漫堤开展了大量卓有成效的工作,尹宝树等(2006)针对黄河三角洲示范区,根据越浪的程度提出了漫堤灾害风险等级标准,并且根据其建立的浪潮耦合数值模式以及长期预报的结果,提出了风险评估的方法和步骤；傅赐福等(2014)利用ADCIRC+SWAN耦合模型构建福建沿海精细化风暴潮漫堤数值预报系统,利用预报潮位与近岸海浪波高数据采用通用爬高公式计算各海堤堤前爬高与海堤高程对比,并按爬高的大小分为四个漫堤等级；张莉等(2016)基于SWAN海浪模式和MATLAB软件,建立了福建沿海天文潮-风暴潮-台风浪耦合漫堤风险评估系统。

福建沿海东临广阔的西北太平洋,由于其特殊的地理位置,每年夏季和秋季经常遭受到台风的侵袭,是我国易受海洋灾害影响的区域之一,海洋灾害形势复杂严峻。近年来,随着全球变暖和海平面上升,极端气候事件不断增多,热带气旋中产生台风的比例上升,登陆或者严重影响福建沿海的强台风频率不断上升,强台风甚至超强台风频繁袭击福建沿海。本文拟针对台风风暴潮灾害较为严重的福建海域建立浪潮耦合漫堤系统,设计四种不同强度的台风等级对福建沿海206条海堤进行漫堤风险等级评估研究,对比了在不同台风强度下风暴潮增水与波浪爬高对于漫堤风险的影响程度,为我国的海洋防灾减灾提供参考。

1 高分辨率风暴潮漫堤风险评估系统的建立

为了进行台风过程下的漫堤风险评估,本研究首先使用ADCIRC+SWAN耦合模式针对福建沿海建立了高分辨率的天文潮-风暴潮-台风浪耦合数值系统。

1.1 台风风场和气压场

要模拟台风引起的风暴潮和海浪,需要将台风过程中的气压场和风场作为模式的强迫场,风场的准确性直接影响模式模拟的准确性。经过比较和验证,本研究采用了在福建沿海地区应用效果较好的Jelesnianski台风模型(Jelesnianski,1965),该模型考虑了台风移动速度对于台风风场的影响,其对风场和气压场的计算公式如下：

当0≤r≤R时

当r＞R时,

其中,Vx,Vy分别为台风中心的移动速度；(x0,y0)为台风中心坐标；(x,y)为计算点的位置；θ为入流角,取20°；β为衰减系数,取0.4；r为计算点到台风中心的距离；R为最大风速半径；(Wx,Wy)为计算点的风速；Pa为计算点的气压；中心气压P0使用Atkinson等(1977)提出的经验公式,如公式(3)所示,其中WR表示最大风速,P∞表示无穷远处的大气压(取1013.25hPa)。

最大风速半径是确定台风影响范围的关键参数,本研究中最大风速半径R使用MEF经验公式(许灵静等,2016)：

其中,P∞为无穷远处的大气压(取 1013.25hPa)；P0为台风的中心气压(单位hPa)；φ为台风中心的纬度(单位为°)；v为台风中心的移动速度；M为起算半径,取45km。台风的位置、中心气压等基本信息来自中国气象局热带气旋资料中心最佳路径数据集(Yinget al,2014)(http：//tcdata.typhoon.org.cn/)。

1.2 风暴潮模式

目前,可用于风暴潮模拟的数值模式有很多,风暴潮模式不仅要使网格的计算区域足够大,能够描述海水的运动,而且要考虑沿海地区复杂地形对模拟结果的影响。ADCIRC型具有计算效率高和稳定性好的优点,应用较为广泛。因此,本研究选取由美国北卡罗来纳州大学Luettich博士和美国诺特丹大学Westerink博士等人共同主持开发的ADCIRC模式(Westerinket al,1992)。ADCIRC采用三角形网格可以在水深变化剧烈的区域使用更高的分辨率,能很好的拟合复杂的岸线,而在地形变化缓慢的地方采用低分辨率,既满足了计算要求,又节省了计算资源和计算时间(Blainet al,1998)。

为了避免数值解法引起的错误,ADCIRC求解通用波连续方程(GWCE)计算水位的变化,GWCE由连续方程带入动量方程得到,球坐标系下垂向平均的二维连续方程形式如下：

其中,t代表时间；Re代表地球的半径；λ和φ分别代表经度和纬度；ξ代表自由海面距离大地水准面的高度；H代表海面到海底的总深度；U和V分别代表经向和纬向的垂直平均流速。

ADCIRC求解基于深度平均的动量方程获得流场数据,动量方程在球坐标中的形式如下所示：

其中,g代表重力加速度；Ps代表海表大气压；Re代表地球半径；ρ0代表海水的密度；f代表柯氏参数；λ、φ、z分别代表经度、纬度和垂向深度；u、v、w分别代表流速的在球坐标3个坐标上的分量；为球坐标下对时间t的全微分；代表垂向应力；mλ,mφ代表径向和纬向的水平应力梯度；bλ,bφ代表径向和纬向的斜压压强梯度；η代表牛顿平衡潮势；α是常数。更多有关ADCIRC模型的详细描述可参看其用户手册(http：//adcirc.org/home/documentation/users-manual-v52/)。

1.3 海浪模式

本研究采用的目前国际上通用的第三代海浪模式SWAN(Booijet al,1999),SWAN模式是由荷兰Delft科技大学在WAM模式基础上研发的第三代浅水海浪模式,经过多年不断的改进,已经趋于成熟。SWAN支持三角形网格,避免了精细模拟时的嵌套计算,可以高效地对包含远海和近岸整个区域的海浪进行模拟。

当波浪在背景流场中传播时,波浪的动谱能量密度守恒,但是能量密度不守恒。因此,在SWAN模型中,用二维动谱能量密度N(σ,θ)来描述随机海浪场E(σ,θ),两者之间的关系有：

其中,σ为相对频率,θ为波向。

SWAN海浪模式使用波作用平衡方程,笛卡尔坐标系和球坐标系都可以使用,笛卡尔坐标系下的方程形式：

其中,x,y为地理坐标；t为时间；θ为传播方向；σ为相对频率；左边第一项代表动谱能量密度N随着时间的变化,第二项和第三项代表动谱能量密度N在几何空间的传播(变化率为Cx和Cy),第四项和第五项代表地形和流场导致的频移和折射对动谱能量密度N的作用。Stot是波能源函数项,代表系统中能量的输入和输出,包括风浪的生成、发展,底面摩擦、白浪、波浪破碎导致的能量损失和波-波相互作用。

SWAN中的源项由如下6个过程组成：

其中Sin代表风能输入项,Snl3代表三波非线性相互作用,Snl4代表四波非线性相互作用引起的能量传递项,Sds,w代表波陡过高引起的白帽耗散,Sds,b代表底摩擦引起的能量耗散,Sds,br代表水深变浅导致海浪破碎造成的能量耗散。

1.4 耦合模式

海洋中的波动和流动现象是相互影响的,水位和流场会影响波浪的传播以及破碎位置；波浪传播所产生的辐射应力会影响水位和流场。因此,波流过程的耦合十分必要。本研究使用美国北卡罗来纳州立大学(North Carolina State University)的Dietrich博士开发的ADCIRC和SWAN的耦合方法(Dietrichet al,2011)。

ADCIRC和SWAN是由同一个风场驱动的。由于ADCIRC的数值格式是发散的,所以根据Courant条件,其运行的时间步长为2s；但是SWAN采取全隐式绝对稳定的差分格式是无条件收敛的,所以其计算的时间步长设置为1800s,本耦合模式数据交换的时间间隔SWAN计算时间步长相同。

漫堤风险等级评估系统根据耦合模式计算结果(总水位、波浪爬高)对福建海域206条海堤进行可视化风险等级评估,海堤数据来源于海洋公益性行业科研专项(2009513)。Feng等(2016)对这个耦合模式进行了详细的检验,与实测结果吻合的很好,为可靠的漫堤风险评估打下了良好的基础

2 数值实验方案

本研究使用的计算网格引用自(Fenget al,2016),该网格从远海到近岸逐渐加密,沿岸网格分辨率最高达50m,对福建沿海复杂地形和海堤信息实现了高精度的刻画,同时又大大提高了计算效率。模式在开边界处加入了8个主要天文分潮(M2、S2、N2、K2、K1、O1、P1和Q1)的调和常数进行驱动,调和常数来自于NAO99潮汐模型(Matsumotoet al,2000).

图1 计算区域水深(a)和模式计算网格(b)Fig.1 Water depth of computational domain(a)and computational grid of the model(b)

2.1 漫堤风险等级评估方案

海堤的漫堤和越浪程度主要受到海堤的高度、水位和波高的控制,在不同的水位和波浪条件下,同一海堤形成越浪和漫堤的可能性差别很大。因此,海堤的越浪漫堤风险,应当依据实际情况进行具体的评估。本文采用尹宝树等(2006)提出的漫堤越浪风险程度标准,将海堤风险等级根据堤前水位与波浪爬高之和与海堤高程的对比分成五个等级：

一级：最大波浪爬高远未达到堤顶(距离堤顶大于0.5m)；

二级：最大波浪爬高接近堤顶(距离堤顶小于0.5m)；

三级：部分海浪越过堤顶,越浪率小于或等于13%；

四级：越浪率大于13%,至总水位与堤顶齐高；

五级：总水位高于堤顶高度。

总水位为风暴潮与天文潮的耦合水位。

2.2 波浪爬高计算

波浪爬高是波浪向岸传播遇到海堤时,海堤坡面上的水体的爬高程与静水高程之差。波浪爬高的准确计算在漫堤灾害预报中起着重要的作用。《海港水文规范》推荐的直立堤波浪在爬高公式为Hc/d=B(HW/d)A,其中Hc为爬高,d为堤前水深,HW为入射波高,A、B为系数,分别为：

式中,Ts为有效波周期；为平均波周期；T*为无因次波周期。

3 台风天兔漫堤风险等级评估

北京时间2013年9月17日02时,台风天兔(1319)在菲律宾以东的西北太平洋上形成,最大风速18m/s(8级),台风中心最低气压1000hPa；天兔缓慢向西移动并且逐渐增强,18日20时增强为台风,19日11时增强为强台风,9月19日17时增强为超强台风并且维持了30h,期间最大风速为52m/s(16级),中心最低气压为930hPa；之后台风继续往西北方向移动,于21日20时减弱为强台风,9月22日19时40分左右在广东省汕尾市南部沿海地区登陆,登陆时最大风速为45m/s(14级),中心最低气压为940hPa。台风天兔期间受风暴潮和台风浪的共同影响,福建省受灾人口9.29万人,紧急转移安置14.46万人,水产养殖受灾面积达3.03千公顷,40艘渔船遭受毁坏,264艘渔船损坏,码头损毁1.26km,防波堤损毁7.36km,海堤和护岸损毁4.86km,直接经济损失达6.36亿元。天兔在登陆时恰逢天文大潮,风暴潮增水与天文大潮叠加,加重了灾情。本文对超强台风天兔期间,福建沿海的风暴潮漫堤风险等级进行评估,同时以超强台风天兔为基础,设计四种不同台风强度,研究不同台风强度下福建沿岸海堤的漫堤风险等级。

图2 台风天兔路径Fig.2 Trajectory of the typhoon Usagi

将风险等级评估结果与本次台风期间福建沿海有实况调查的两条海堤的漫堤情况进行比较验证,实况调查数据来自(张莉等,2016)。系统评估结果显示,东山县沃角海堤和东山县康美海堤漫堤风险等级为五级,即总水位高于堤顶高度,说明这两处海堤会发生漫堤、漫滩甚至溃堤灾害,在灾后调查中发现,这两处海堤被不同程度的毁坏,说明本系统模拟得到的漫堤风险等级评估基本准确。

表1 台风强度等级Tab.1 Typhoon intensity scale

表2 台风天兔过程漫堤风险等级预报及实况对比Tab.2 Comparison in risk level of overtopping seawall between the hindcast and the actual situation during Typhoon Usagi

在不同台风强度的情况下,对海堤逐条进行了风险等级划分,并在此基础上做了统计分析,并将结果展示在(图3)中。得出在强度一的情况下,有35处海堤的漫堤风险等级达到了四级,10处海堤的漫堤风险等级达到了五级；在强度二的情况下,有35处海堤的漫堤风险等级达到了四级,7处海堤的漫堤风险等级达到了五级；在强度三的情况下,有48处海堤的漫堤风险等级达到了四级,10处海堤的漫堤风险等级达到了五级；在强度四的情况下,有55处海堤的漫堤风险等级达到了四级,10处海堤的漫堤风险等级达到了五级。

以下是在4种强度下漫堤风险等级均达到了五级的7处海堤,有必要做好防范措施。

选取了其中四条典型的海堤做出堤前风暴潮增水和波浪爬高过程曲线。

由图4可知,在台风强度为最低级别时,四条海堤的风暴潮增水普遍低于1m。随着台风强度的增强,风暴潮增水也在逐渐增加,当台风强度达到最大时,风暴潮增水最大可达1.3m。但是风暴潮增水的值相对于福建沿岸海堤的高度较小,台风强度增强导致的风暴潮水位的增长对海漫堤风险等级影响程度很小。在福建省现有的海堤条件下,漫堤风险等级达到五级的海堤一般都是高度较低的海堤,除了强度二由于台风强度较低,只有7条海堤漫堤风险等级达到五级以外(图3),其余台风强度条件下,都是漳浦县杜浔海堤等10条海堤的漫堤风险等级达到五级(图3),风暴潮水位的升高没有使更多的海堤的堤前水位超过海堤高度。

表3 海堤信息Tab.3 Seawall information

图3 不同台风强度下的漫堤风险等级Fig.3 Risk level of overtopping seawall under different typhoon intensities

图4 四种不同强度情况下风暴潮增水时间序列Fig.4 Time series of storm surge under four typhoon intensity scales

由图5可知,即使在最低台风强度情况下,波浪爬高普遍都超过了2m,在强度四的情况下波浪爬高最大可达4.7m,波浪爬高的值比风暴潮增水的值大很多,考虑波浪爬高后,海堤由于越浪导致的漫堤风险会大大增加,从强度二到强度四,漫堤风险等级为四级的海堤从35条迅速上升到55条(图3)。台风浪引起的4m高的越浪对海堤的冲击力会很大,会导致海堤有溃堤甚至毁堤的风险,需要增强沿岸海堤对于台风浪的防范措施。

图5 四种不同强度情况波浪爬高时间序列Fig.5 Time series of wave run-up under four typhoon intensity scales

4 结论

本文建立了福建沿海地区天文潮-风暴潮-台风浪耦合漫堤风险等级评估系统。通过波浪爬高所及高程与海堤高程进行比较的方式来判断漫堤风险等级；选取2013年超强台风天兔,并且制定了4种不同的台风强度等级,通过数值模拟对福建沿海206条海堤进行了漫堤风险等级评估,并给出了漫堤风险等级分布图。结果表明：波浪爬高对漫堤风险的影响高于单纯的风暴潮增水；风暴潮增水随台风强度的增强增量较小,对于漫堤的风险影响较小；福建沿海波浪爬高普遍较高,随着台风强度的增强,波浪爬高会显著增加漫堤的风险等级,并且应该重视台风浪对海堤造成的冲击所导致的溃堤甚至毁堤灾害,本研究可为沿海防灾减灾提供科学依据。

致谢本文的数值模拟是在中国科学院海洋研究所的高性能计算平台上完成的,特此致谢。