台州温岭市风暴潮淹没危险性分析
2017-01-09郑国诞谢亚力胡金春陈韬霄
郑国诞,谢亚力,胡金春,陈韬霄
台州温岭市风暴潮淹没危险性分析
郑国诞1,2,3,谢亚力1,2,3,胡金春1,2,3,陈韬霄1,2,3
(1.浙江省水利河口研究院,浙江杭州310020;2.浙江省海洋规划设计研究院,浙江杭州310020;3.浙江省河口海岸重点实验室,浙江杭州310020)
利用MIKE21FM模块建立了适用于台州温岭市的高分辨率风暴潮漫滩数值模式,以9711号台风路径为基准,分别以5 km为间隔向两侧平移,构造各种设计台风路径,找到最不利登陆位置进行风暴潮计算,再利用SWAN计算不利路径的堤前台风浪要素,通过经验公式计算越浪量,进行溃堤判断。最后将各路径的淹没计算结果取最大包络,画出6档台风的风暴潮淹没范围图。结果表明:在最高等级的915 hPa时,由于沿海堤防几乎全溃,因此处于沿海平原区的乡镇基本受淹,随台风强度降低淹没水深及范围逐渐减小,最低等级的965 hPa时无淹没,其他等级风暴潮淹没介于两者之间。
风暴潮;台风浪;淹没;越浪量;溃堤
1 引言
台州温岭市地处浙东南沿海台州湾以南,东临东海,南接玉环县,西邻乐清市,北连路桥区,岛屿较多,岸线较长,海岸线总长约316.3 km,拥有丰富的海洋自然资源。温岭市属中亚热带季风气候,受海洋气象影响明显,是热带气旋频繁活动区域,易受台风风暴潮、近岸波浪的影响和侵袭[1],导致水位暴涨、岸堤决口,淹没农田、围区、厂房等,给温岭市沿海经济发展、人民生命财产安全造成严重损失,同时也制约了当地海洋经济发展。如9417号、9608号、9711号和0414号等台风都造成了水利设施的损坏、人员的伤亡和大量的经济损失,因此开展风暴潮淹没分析工作对台州温岭市海洋经济发展和防灾减灾工作有着十分重要的意义。
目前,国内外很多学者对风暴潮漫滩、漫堤等进行了大量的研究工作,且取得了一定的成果。吴玮等[2]利用ADCIRC(A(Parallel)Advanced Circulation Model for Oceanic,Coastal and Estuarine Waters)模式计算了不同台风强度等级的风暴潮对台州市的影响,在模式中充分考虑了一维海堤情况,并对不同强度的淹没情况进行了评估分析;殷杰等[3]采用了高精度洪水数值模型(FloodMap)构建了两处溃堤点6种重现期台风风暴潮溃堤情景,开展上海台风风暴潮淹没情景模拟,结果显示上海受台风风暴潮灾害影响有限;朱军政等[4]利用二维浅水方程建立钱塘江河口地区和杭州湾的风暴潮溢流模型,探讨风暴潮出现溢流的计算方法,模拟了风暴潮漫溢堤防的过程。国内外还有众多学者对风暴潮灾害情况进行了计算分析研究[5-9],大部分以漫堤为主,但是随着极端气候的频繁发生,超强台风发生的概率增大,海堤溃堤风险也随之增大,不容忽视。本文在前人研究的基础上,建立了适用于台州温岭市的风暴潮和台风浪计算模型,利用越浪量进行了海堤溃堤判断,最后开展了不利路径不同强度台风影响下的台州温岭市风暴潮溃堤淹没情况计算。
2 风暴潮模型简介
2.1 控制方程
本次计算模型采用丹麦水力学研究所的无结构网格平面二维模型MIKE21软件FM模块。采用三角形网格和四边形网格的耦合计算模型,该模型网格节点布设灵活,与岸线吻合良好,利于拟合堤防等复杂边界线,便于局部加密,具有算法可靠、计算稳定的优点,已在很多工程中得到成功应用。
基本方程包括一个连续性方程和两个动量方程:
式中:ς为潮位;g为重力加速度;p、q分别为x、y方向上的垂线平均单宽流量;h为水深;Ω为科氏力参数;ρw为水密度;C为谢才系数;Pa为大气压力;f为风摩擦系数;V,Vx,Vy分别为风速及其在x、y方向的分量;Ex、Ey为涡动粘性系数。
初始条件:
边界条件:
(1)水边界:ς(x,y,t)=ς*(x,y,t)带“*”表示已知值;
(2)陆边界:法向通量为0。
法线方向流速为零。
有了以上条件,就可用一定的离散格式求出方程的解。目前求解上述方程的数值计算方法很多,较为流行的有控制体积法、交替方向法(Alternating Direction implicit method,ADI)、破开算子法、直接差分法、特征线法和有限单元法等,MIKE21软件FM模块采用的是控制体积法显式迎风格式。
图1 风暴潮计算模型示意图
图2 研究区域网格图
2.2 计算范围及网格布置
模型计算范围及网格如图1和图2所示,采用三角形网格和四边形网格耦合的方法剖分计算域,其中对研究区域的堤坝及道路采用四边形网格,其他区域采用三角形网格,由480 418个单元和246 395个节点组成,台州近海及陆上网格边长约50 m,外海边界基本囊括台风第一警戒线,涵盖渤海、黄海、东海、日本海、朝鲜海峡、台湾海峡、长江口、杭州湾及钱塘江,并对浙江沿海、台州附近海域局部加密网格,椒江上游边界至临海上游7 km的三江村。坐标投影统一采用北京54坐标6度带(中央经线123°)。
2.3 边界条件
外海边界取静压水位叠加天文潮位,该潮位由全球潮波模型TPXO.6(网址为:http://csdms. colorado.edu/wiki/Data:TPXO6.2)推算求得,其中包含8个主要分潮M2、S2、K1、O1、N2、P1、K2、Q1,以及两个长周期分潮Mf和Mn。该模型的全球网格数为1 440×721,分辨率为0.25°,基本能够构造出外海深水处真实的天文潮过程,其式如下:
式中:ζ0为边界处的潮位,ζp为边界处静压水位,i为1—10,分别对应上述10个分潮,Ai、αi分别为分潮在3条边界处的振幅和迟角,ωi为分潮的角频率。
2.4 风场和气压场
在风暴潮计算中,台风场和气压场的计算是一个重要的环节。本次工作选用Jelesnianski模型风场和气压场[10],其式如下:
式中:R为最大风速半径;r为计算点到台风中心的距离;V0为台风移动速度;WR为台风最大风速;A=-[(x-xc)sinθ+(y-yc)cosθ]、B=(x-xc)cosθ-(y-yc)sinθ;(x,y)、(xc,yc)分别为计算点坐标和台风中心坐标;θ为流入角(计算中当r≤R时θ取10°,当r>1.2R时θ取25°,其余的θ在10°和25°之间线性内插);P0为台风中心气压,P∞为无穷远处的大气压(计算中取1 010 hPa);β为台风风速距离衰减系数。
最大风速计算使用Atkinson Hollidy提出的风-压关系式:
式中:P0为中心气压;R为最大风速半径;Rk为经验常数,推荐值为40,也可由气压或风速的拟合精度予以调节。
3 台风浪模型
3.1 基本方程
第三代浅海海浪数值模式(Simulating WAves Nearshore,SWAN)能够描述在浅水区特定的风、流和水下地形条件下的波浪场的演化。控制方程如下:
式中:N为波作用量;σ为波的相对频率;θ为波向;S项为源汇项。式中左端第一项为波作用量的局地变化;第二、三项为波作用量在地理空间上的传播,其中Cx和Cy分别为波作用量在x和y方向上的传播速度;第四项为由于水深和水流变化造成的相对频率变化,其中Cσ为波作用量在频率空间中的传播速度;第五项为水深和流产生的波浪折射,其中Cθ为波作用量在波向空间中的传播速度。波作用量方程的离散基于迪卡尔坐标下的隐式差分,对于波的传播和源项采用固定的时间步长。
3.2 模型范围
台风浪模型计算域西至广东汕头,南至台湾岛南端,东至日本,北至山东半岛,大致范围为21.5°—37°N、116.5°—129.5°E,计算域面积为2.21×106km2。计算域采用矩形网格并用4层嵌套逐步加密,其中大范围的网格尺寸为8 000 m×8 000 m;第二层网格尺寸为2 000 m×2 000 m,覆盖了杭州湾以南的浙江沿海;第三层网格尺寸为400 m×400 m,覆盖了健跳至坎门的台州外海;第四层网格尺寸为100 m× 100 m,覆盖临海至温岭沿海(见图3),各层粗网格计算结果同时为下一层细网格计算提供边界条件。
台风浪计算采用的风场模型与风暴潮计算所用风场一致,风场计算采用的台风路径、中心气压、最大风速、最大风速半径等台风参数的取值均与风暴潮计算相同。
图3 台风浪计算网格
图4 9711号台风风暴潮位验证图
4 模型验证
作为对一种大范围海洋水体运动现象的模拟,模型的验证必须建立在对较大范围内观测点的观测值进行模拟的基础上,同时,模型验证所选择的台风也必须有代表性和典型性。本次主要验证在浙江或者福建登陆的对台州影响较大的台风,如8707号、9417号、0414号“云娜”、0608号“桑美”、0908号“莫拉克”以及1323号台风“菲特”等,由于篇幅关系,本文仅罗列9711号台风的验证情况为代表,潮位验证见图4,台风浪验证见图5。
5 计算条件的确定
5.1 路径选取
9711号台风是1949年以来影响温岭市最严重的、风暴增水最显著的典型台风,其过程最低气压920 hPa,登陆气压960 hPa,登陆地点在温岭南部,该台风登陆时恰逢天文大潮高潮位,形成了沿海历史最高潮位。0608号“桑美”台风是1949年以来登陆浙江气压最低,风速最强的超强台风,其过程最低气压915 hPa,登陆气压920 hPa,登陆地点为距温岭约170 km的苍南,虽然由于其登陆时遭遇天文潮低潮位,且台风半径较小,登陆点离温岭较远,对温岭造成的影响不是很严重,但如果将该台风登陆点北移至玉环、温州一带,且在大潮高潮位时登陆必然会形成温岭市极端高潮位.为了在同等条件下比较9711号台风和0608号台风路径对温岭增水的影响程度,设计以下台风参数:设定两条路径台风气压均为915 hPa,且沿程不变;台风半径均取35 km;将0608号台风平移至9711号台风相同的登陆位置。结果显示9711型路径和0608型路径海门的最大增水分别为6.00 m、5.64 m,温岭石塘的最大增水分别为2.14 m、1.93 m可见9711型路径增水相对较高,因此选定9711型路径作为设计台风路径。
5.2 台风强度
将台风强度划分为6档,由于台风依据中心附近地面最大风速来划分等级,而风暴潮模式的输入参数为台风中心最低气压,因此依据统计结果,建立台风中心最低气压和中心附近最大风速的对应关系,如表1所示。
图5 9711号台风期波浪验证
表1 台风强度等级表
5.3 天文潮的叠加
一次风暴潮过程的最大值可以发生在天文潮的任何时段,但是最大淹没风险是发生在风暴增水叠加到当地天文高潮时。选取代表站(石塘验潮站)6—10月份19 a的逐月高潮平均值2.10 m作为天文潮控制条件。
5.4 溃堤判断
溃堤处理在本模型中是一个难点,溃堤过程,溃堤何时发生很难准确估计,它取决于潮位和浪高的综合作用,同时溃堤宽度较难确定。一般把海堤结构、风暴潮位及堤前波浪情况三者作为是否溃堤的主要因子,较为理想的手段是针对研究区域内各代表海堤开展波浪模型试验,用各级台风对应的计算潮浪条件检验海堤受损情况,研究其溃决过程,但由于模型试验工作量十分巨大,实际工作中很难做到。
海堤越浪量可体现海堤结构、风暴潮位及堤前波浪这三者的综合效应,越浪量过大是超强台风中海堤溃决的主因,越浪量与海堤溃决具有较明确的对应关系,因此本次工作采用计算海堤越浪量作为是否溃堤的判据。
对于带胸墙的斜坡堤采用《海港水文规范》(JTS145-2-2013)中的经验式计算越浪量:
对于直立堤采用《浙江省海塘工程技术规定》中的经验式计算越浪量:
无风条件下:
有风的越浪量为无风条件下的越浪量乘风校正因子K′:
式中:Wf取决于风速的系数,具体参数详见《海港水文规范》(JTS145-2-2013)和《浙江省海塘工程技术规定》。根据堤前各代表点的台风浪过程及同步风暴潮位,按上述公式计算研究区域内海堤的越浪量,并对照阀值确定是否溃堤,考虑到保护区重要程度,参考《浙江省海塘技术规定》,确定本次溃堤阀值越浪量取0.05 m3/(m·s)。
一旦判定溃堤,则设定溃堤在一瞬间完成,不考虑溃堤过程。溃堤后的高程取堤内地面高程,溃堤长度由分布于该段堤防的代表点的越浪量决定,每个代表点代表一段堤防长度。
6 计算结果
6.1 台风路径比选
针对6档台风强度,以选取的9711号台风实际路径作为基准位置,并分别以5 km为间隔向两侧平移,构造各种设计台风路径,进行海平面条件下台风纯增水计算,以寻求最不利(增水最大)的登陆位置。在温岭沿海布设金清、东海塘、松门、石塘、东浦、担屿、江厦共7个代表点(见图6)分析其增水大小。
计算结果见表2,对于面向东海的前6个代表点基本是在台风登陆点右侧最大风速半径处出现增水最大值,等级最高的915 hPa台风时其最大增水3.5—5.8 m,等级最低的965 hPa台风时其最大增水1.7—3 m。对于面向乐清湾的江厦代表点,其增水是由外海风暴潮经乐清湾口传入后在乐清湾沿岸形成,故增水最大值受控于乐清湾口的增水最大值,因此台风路径需进一步南移该点才出现增水最大值,由于位于乐清湾顶水体堆积,该处增水最大可达.1 m,即使在965 hPa台风时增水仍超过3.3 m。
6.2 两潮耦合计算
通过路径比选,每档台风均筛选出了最不利的台风路径。将这些不利路径台风叠加当地台风期平均高潮位(石塘2.10 m)的天文潮过程,并使最大增水遭遇天文高潮位,进行两潮耦合风暴潮计算。各海堤前沿的计算结果见表3,海堤位置见图7。
封闭海湾的湾顶附近,水体易于堆积,而后者位于湾口附近。在面向乐清湾的海堤中,由于靠近湾顶,风暴高潮位较高,在915 hPa,台风时达到8.15 m,在965 hPa,台风时为5.81 m,但该地风浪较小,且浪向基本离岸或顺岸方向,因此乐清湾中不考虑溃堤,只计算当潮位高于堤顶后的漫堤淹没。
6.3 台风浪计算
根据风暴潮中确定的针对不同区域的不利台风路径分别计算温岭沿海台风浪,所采用台风风场与风暴潮的相同,台风浪计算时叠加风暴潮位背景场。
图6 增水代表点位置图
图7 海堤位置图
根据台风浪计算结果,选取受向岸浪作用的海堤进行溃堤判断,其他海堤按漫堤情形进行淹没计算。各堤位置见图7,堤前最大有效波高列于表4。由表可见:黄礁门、担屿堤前波高相对较大,因该海堤向东面朝开敞海区,无岛屿或岸线掩护,其有效波高在915 hPa台风时达到3.6 m;上箬堤因其朝向南,且堤线内凹波浪不易传入,其最大有效波高小于1.7 m;其他海堤的波高介于两者之间。
6.4 风暴潮淹没计算结果
针对6档台风强度及各区域的不利路径进行风暴潮淹没计算,其中受向岸浪作用的海堤按溃堤考虑,非向岸浪作用的海堤或无堤防的岸段以及经计算不溃的海堤均按漫堤考虑。溃堤计算时,将网格地形设为镇压层或堤内地面高程,通过设定“门”在某一时刻开启来模拟溃堤现象。漫堤计算时,将网格地形设为海堤堤顶高程作为防潮高程,不考虑挡浪墙的防潮功能,即假设挡浪墙在风暴潮来时已经
损毁,当外海水位高于堤顶时自然漫溢。其中溃堤判断则根据之前所述按越浪量0.05 m3/(m·s)的阀值判断,如根据前面越浪量计算公式计算,如果堤前波浪较小,根据计算越浪量没达到0.05 m3/(m·s)或者为离岸浪,则不溃堤。在915 hPa超强台风作用下,除了受岛屿或岸线掩护较好的神址、上箬海堤,以及堤防标准较高的黄礁门海堤外,其余海堤均全溃或大部溃决;随台风强度降低溃堤岸段逐渐减少,955 hPa台风时已不存在全溃海堤,不溃岸段的比例远高于溃堤岸段;965 hPa台风时所有海堤均不溃,溃堤判断结果见表5。
表2 各档台风最不利路径(单位:km)及其最大增水(单位:m)
表3 各级台风各不利路径时的最大风暴潮位(单位:m)
表4 各级台风不利路径时海堤前沿的最大有效波高(单位:m)
表5 各海堤溃堤判断结果
图8 不同气压下台风淹没范围图
将各路径的淹没计算结果取最大包络,6档台风的风暴潮淹没范围分别见图8。由图可见,在最高等级的915 hPa时,由于沿海堤防几乎全溃,因此处于沿海平原区的乡镇基本受淹,其中东海塘围区、东浦镇、观岙、担屿围区、乐清湾沿岸受淹严重,横河、滨海、箬横、新街、松门由于东海塘后的二线塘保护受淹相对较轻。随台风强度降低淹没水深及范围逐渐减小,最低等级的965 hPa时无淹没,其他等级风暴潮淹没介于两者之间。
7 结论与不足
7.1 结论
(1)本文以遭受风暴潮灾严重的台州温岭市为例,以9711号台风路径为基础,利用MIKE21FM模块进行风暴潮计算,通过平移路径计算不同堤段的最大增水,构造出对台州温岭市各堤段最不利的不同强度等级台风路径,在方案计算中,利用SWAN模式构建台风浪模型计算不同强度不同路径下堤前波要素,以越浪量因子作为溃堤判断依据,最后计算分析了不同强度台风下台州温岭市的最大淹没范围,对当地的海洋防灾以及产业布局有着重要的指导意义。其中利用台风浪模型和风暴潮模型,综合计算越浪量,来进行溃堤判断,并以设定“门”在某一时刻开启来模拟溃堤现象为本文的创新之处;
(2)通过众多方案的风暴潮和台风浪计算,可知在最高等级的915 hPa时,由于沿海堤防几乎全溃,因此处于沿海平原区的乡镇基本受淹,其中东海塘围区、东浦镇、观岙、担屿围区、乐清湾沿岸受淹严重,横河、滨海、箬横、新街、松门由于东海塘后的二线塘保护受淹相对较轻。随台风强度降低淹没水深及范围逐渐减小,最低等级的965 hPa时无淹没,其他等级风暴潮淹没介于两者之间。
7.2 不足
(1)本文构造的台风路径偏不利且增水是叠加至当地天文高潮,为“不利”加“不利”双重叠加,结果可能会偏保守;
(2)越浪量的计算采用规范所列的经验公式计算,而实际中影响越浪量的因素较多,公式计算存在误差,且目前溃堤过程概化较为简略,这些不足将在今后研究中逐步改善。
[1]朱业,丁骏,卢美,等.1949-2009年登陆和影响浙江的热带气旋分析[J].海洋预报,2012,29(2):8-13.
[2]吴玮,刘秋兴,于福江,等.台州沿海地区台风风暴潮淹没风险分析[J].海洋预报,2012,29(2):25-31.
[3]殷杰,尹占娥,于大鹏,等.基于情景的上海台风风暴潮淹没模拟研究[J].地理科学,2013,33(1):110-115.
[4]朱军政,于普兵.钱塘江河口杭州湾风暴潮溢流计算方法研究[J].水科学进展,2009,20(2):269-274.
[5]何佩东,左军成,顾云碧,等.普陀沿海风暴潮淹没危险性评估[J].海洋湖沼通报,2015(1):1-8.
[6]莎日娜,尹宝树,杨德周,等.天津近岸台风暴潮漫滩数值模式研究[J].海洋科学,2007,31(7):63-67.
[7]傅赐福,于福江,王培涛,等.滨海新区温带风暴潮灾害风险评估研究[J].海洋学报,2013,35(1):55-62.
[8]梁海燕,邹欣庆.海口湾沿岸风暴潮漫滩风险计算[J].海洋通报, 2004,23(3):20-26.
[9]聂会,孙志林,黄森军,等.新岸线下近期台风暴潮模拟[J].浙江大学学报(理学版),2015,42(2):205-211.
Inundation risk assessment of typhoon storm surge along Taizhou Wenling City
ZHENG Guo-dan1,2,3,XIE Ya-li1,2,3,HU Jin-chun1,2,3,CHEN Tao-xiao1,2,3
(1.Zhejiang Institute of Hydraulics and Estuary,Hangzhou 310020 China; 2.Zhejiang Institute of Marine Planning and Design,Hangzhou 310020 China; 3.Zhejiang Provincial Key Laboratory of Estuary and Coast,Hangzhou 310020 China)
The typhoon storm surge is the main type of natural disasters in Wenling city.So it has a very important significance to analyze the storm surge inundation risk.In this paper,the MIKE21FM module was used to establish a high-resolution storm surge floodplain numerical model which was applicable to Wenling City. Based on the Typhoon“9711”,the value 5 km as the intervals on both sides of pan was used to construct various design typhoon path,in order to find the most unfavorable landing position to calculate the storm surge,and then the SWAN model was used to calculate the wave elements of the most unfavorable landing position.The empirical formula was used to calculate overtopping discharge to judge dike break.Finally,the largest envelope value from the calculated result was used and drew the storm surge submerged area map under six-speed typhoon.At the highest level 915 hPa,dike of coastal was almost broken,and the rural communities along the coastal were basically flooded,while there was no flooding at the lowest level 965 hPa,and the flooding extent at other grades of the storm surge were in between.The typhoon intensity was decreased with depth and scope of the drowned gradually.
storm surge;storm waves;flooding;overtopping;dike break
P731.23
:A
:1003-0239(2016)06-0040-11
10.11737/j.issn.1003-0239.2016.06.005
2016-05-26
水利部公益性行业科研专项(201401010);浙江省科技计划项目(2014C33057)。
郑国诞(1986-),男,工程师,硕士,主要从事河口海岸水动力学及海洋灾害防洪等研究。E-mail:12971414@qq.com