陆面-河网耦合模型在极端风暴潮淹没危险性分析中的应用

2021-06-03贝竹园邱桔斐章晓洁

华东师范大学学报（自然科学版） 2021年3期

滕飞, 李路, 贝竹园, 王军, 邱桔斐, 章晓洁

（1. 上海市水利工程设计研究院有限公司, 上海 200061; 2. 上海滩涂海岸工程技术研究中心,上海 200061; 3. 上海市海洋管理事务中心, 上海 200050; 4. 国家海洋局东海海洋环境调查勘察中心, 上海 200137）

0 引言

上海市面临着较为严重的台风风暴潮威胁, 是我国风暴潮灾害的易发和频发区域之一, 历史上多次出现过造成较大经济损失的强风暴潮灾害[1−2]. 随着沿海风暴潮易发地区土地资源的不断开发, 这些区域的风暴潮灾害的“放大效应”不断凸显[3]. 在国家海洋局国家尺度风暴潮灾害危险性评价结果中,上海市的5个沿江沿海区域被评为Ⅰ级(高), 当地人民的生命财产、城市运行、社会稳定和经济发展受到严重的威胁, 对风暴潮灾害危险性研究的需求较为迫切.

由于风暴潮灾害可能与暴雨、上游洪水等灾害叠加发生, 形成风(台风)、暴(暴雨)、潮(高潮)、洪(流域洪水)双碰头甚至多碰头的情景[4], 因此双因子、多因子影响下的复合风暴潮灾害(海−气耦合[5]、浪−潮耦合[6−8]、潮−洪耦合[9−10]等)成为沿海地区风暴潮灾害脆弱性评估的主要研究对象. 风暴潮水通过漫堤和溃口进入陆域后, 河网对潮水进行调蓄和运输, 在减少靠近海塘区域地面积水量的同时, 也使部分潮水通过河网向上游输运, 造成远离海塘的区域地面积水. 国内外已开展了大量对风暴潮灾害风险的评估研究, 并逐渐将研究重点放在多种灾害复合作用情景下的风暴潮过程模拟分析. 对于上海市面临的风暴潮灾害造成的海塘溃决风险, 现有研究常采用MIKE21、FloodMAP、Delft3D等数值模型进行二维水动力模拟[11−13], 分析海塘溃决后高潮叠加暴雨情形下的区域洪水演进过程和内涝风险, 评估各研究区域在面临上述情境灾害下的脆弱性. 在较大时空尺度的风暴潮灾害模拟中, 现有研究往往对河网的刻画较为简化, 保留区域内主干河道[4,14]的同时忽略了支级河道对潮水的调蓄和运输功能,从而易出现对漫溢进入陆面潮水的时空分布估计不准确问题, 影响受灾区域危险性评估的准确性. 作为城市最大的“海绵体”[15], 河网对区域内涝的调蓄作用不可忽视[16−17], 在上海市易受风暴潮影响的区域建立陆面−河网耦合模型, 精细化模拟风暴潮水进入陆域后的演进过程, 对于准确评估当地的风暴潮灾害具有重要的意义.

在上海沿江沿海5个行政区中, 位于西南部的金山区的社会经济面对风暴潮灾害时脆弱性相对较高[18]. 金山区南高北低，地坪高程为2.5 ～ 4.5 m, 较低的地坪高程使金山区在承受超标“风、潮、暴、洪”等灾害时面临着更大的受灾风险. 上海市有记录的历史数据中, 发生在金山区沿岸的风暴潮灾害占据较高的比例, 自唐代至民国有记录的风暴潮灾害达百起以上[18]. 新中国成立后如5612号台风、7413号台风、8114号台风等引起的较大的风暴潮灾害期间, 当地最高潮位均达到5.5 m以上, 对金山区岸线造成了不同程度的破坏[19]. 根据上海市相关工程建设规范[20]中的上海地区外围主要水文站点设计高潮位, 金山区金山嘴水文站200年一遇设计高潮位为6.83 m, 在5个沿江沿海行政区中最高,故而在相同超标条件下超高水位绝对值大, 风险相对较高. 与此同时, 金山区沿海具有多处重点保护目标, 其中重要工业目标包括上海石化和上海化工区两个重要工业区, 多达100家危险化学品企业,这些承灾体一旦受灾, 极有可能引发危险化学品泄漏等次生灾害. 综上所述, 本研究选取金山区作为典型区域, 基于MIKE11和MIKE21分别建立金山区一维河网模型和二维陆面模型, 并在MIKE FLOOD模块中进行一、二维模型耦合水动力计算, 模拟分析在风暴潮、台风、区间暴雨和流域洪水的综合影响下金山区的积水变化过程, 为当地的风暴潮灾害危险性分析提供依据.

1 研究区域概况

金山区位于长江三角洲南翼, 太湖流域碟形洼地东南端, 东邻奉贤区, 西与浙江省平湖市、嘉善县交界, 南濒杭州湾, 北与松江区、青浦区接壤. 区域东西长约44 km, 南北宽约26 km, 行政区划面积为586.05 km2. 全境地势低平, 地面高程自西北至东南略有升高, 东南地面高程为4 ～ 4.5 m, 中部高程为3.6 ～ 4 m, 西北部大多在3.5 m以下, 最低处在2.5 m以下. 金山区水系西临杭嘉湖浙水下游, 受浙江客水和黄浦江潮水共同影响, 河渠交织成网. 金山区现状河道长度为2 473.64 km, 河湖面积为41.74 km2, 河湖水面率约为7.12%[21]. 根据金山区河道蓝线规划, 当地规划河道长约2 125 km, 其中骨干河道长约336 km, 支级河道长约1 789 km, 河面积53.98 km2, 河湖水面率约为9.21%. 规划河湖布局如图1所示.

图 1 金山区规划河湖布局示意图Fig. 1 Distributions of plan rivers and lakes in Jinshan District

2 模型构建及验证

本研究采用丹麦水力研究所(DHI)研发的MIKE 21模块进行平面二维水动力计算, 用于模拟地面积水变化过程. MIKE 21模块基于Boussinesq静压假定后的简化Navier−Stokes方程, 通过有限体积法求解垂向平均后的二维质量和动量守恒方程.

采用MIKE 11构建一维河网, 用于区域河道水动力计算. MIKE11为一维一层(垂向均质)的水力学模型, 其差分格式采用了Abbott−Ionescu六点中心隐式格式, 其数值计算采用传统的“追赶法”, 即“双扫”算法.

采用MIKE FLOOD模块进行MIKE11与MIKE21一、二维模型动态耦合模拟计算, 可以模拟河道水位超过堤岸漫溢至周边地坪以及地坪积水通过河道堤岸进入河道的动态过程.

各模型统一采用上海城建平面坐标和上海吴淞高程坐标进行计算.

2.1 一维模型基本方程

MIKE11水动力模块的基本方程采用一维非恒定流的圣维南方程组:

式(1)和式(2)中:A为过水断面面积;Q为过水断面流量;x为沿程坐标;t为时间坐标;h为断面水位;C为谢才系数;R为水力半径;g为重力加速度.

2.2 二维模型基本方程

MIKE 21水动力模型基于Boussinesq静压假定后的简化Navier−Stokes方程, 其基本方程具体如下:

2.3 一、二维模型耦合

一、二维模型之间通过MIKE FLOOD模块建立连接关系, 计算中通过双扫描耦合计算方法进行自动耦合. 具体为: MIKE 11模型从MIKE 21模型中提取n时步的总水深h(n), 计算出n时步的流量Q(n); MIKE 11模型内部预测器再根据n时步的流量Q(n)及总水深h(n), 应用公式(6)预测出n+1/2时步的流量Q(n+1/2), 并作为源项提供给MIKE 21模型; MIKE 21模型再根据n+ 1/2时步的流量Q(n+1/2)计算出n+ 1时步的总水深h(n+1), 依此类推.

2.4 一、二维耦合模型设置

一维河网计算范围包括金山区河道及邻近的青浦区、松江区和奉贤区部分河道, 未包括浙江省平湖市范围内的河道. 模型中包括所有规划骨干河道及主要支级河道共计339条, 并对部分支级河道进行概化, 确保金山区范围内规划河湖水面率为9.21%(见图2(a)). 区域内各级河道采用规划河道断面,断面为梯形, 河口顶高取金山区平均地坪高程4.0 m, 骨干河道河底高程设置为 –1.0 m, 支级河道底高程设为0.0 m.

图 2 一维河网模型与二维陆面模型计算范围及河网分布Fig. 2 Range and framework of the 1−D river network model and the 2−D surface model

二维模型基于《上海市海洋灾害风险评估与区划(金山区)》成果[22]提取建立. 模型计算范围为: 东边界至东经121°47′附近, 西边界至东经121°20′附近, 北边界至北纬30°57′附近, 南边界至北纬30°23′附近. 计算区域包括整个金山区、部分奉贤区和浙江省平湖市以及上述区域外侧杭州湾近岸水域, 共70 611个网格单元, 35 607个网格节点, 较好地拟合了岸线. 模型的陆面地形数据中包含了具有挡水作用的海塘及道路工程(如沪杭线), 其他尺度较小的地形信息则进行了概化处理. 杭州湾水域开边界处的模型网格尺寸较大, 约1 000 m, 海塘附近的网格分辨率较高, 网格尺度最小约为90 m(见图2(b)). 研究区域溃口以及海表面风场、气压的设置与《上海市海洋灾害风险评估与区划(金山区)》成果[22]一致. 金山区海堤模拟溃口共设置两处, 分别位于戚家墩至龙泉港出海闸间单塘(位置1)和城市沙滩(位置2), 如图3所示. 溃口宽度为90 m, 溃决方式设置为瞬间全溃, 溃口的底高程取溃口所在海堤堤内的地面高程. 计算风场取典型5612号台风路径平移至金山区的不利条件, 设计台风中心气压设置为900 hpa, 约为该区域200年一遇的低气压. 二维模型水位边界由上述评估成果选取的计算风场叠加金山嘴潮位站连续19年的月最大天文潮的10%超越高潮位数确定, 此时金山区沿岸过程最高潮位可达9.4 m以上.

由于在本研究设置的风暴潮溃堤漫溢工况下, 陆面、河网水位已较高, 城市排水管网已基本无法正常发挥其功能, 且管网本身容量有限, 故在本文中不考虑城市排水管网在模型中的作用. 模型综合金山区土地利用、植被类型等数据在建模过程中对陆面底摩擦系数进行了率定, 并使用前述评估成果[22]进行了验证, 总体上结果良好.

在一、二维模型的重叠区域, 通过链接一维模型河道左右两侧堤岸与二维模型相应的陆域网格,形成侧向流的交互, 完成一、二维模型的耦合建模.

2.5 计算工况设置

为比较有无河网调蓄对金山区风暴潮漫溢积水分布的影响, 设置风暴潮影响下无河网调蓄(工况0)和有河网调蓄(工况1)两种计算工况.

图 3 一、二维耦合模型计算示意图Fig. 3 Construction of the 1−D and 2−D coupled model

为分析降雨的影响, 在工况1的基础上, 设置叠加区间暴雨影响的工况(工况2), 其中模型的降雨根据上海市治涝标准[23]中对设计暴雨重现期的相关要求, 选取金山区20年一遇最大24 h降雨(24 h总雨量191 mm), 雨型为“菲特”雨型, 如图4所示. 考虑到风暴潮通过漫堤、溃口进入陆面时, 地表陆续出现大范围积水, 植被等基本已失去调蓄雨水的能力, 因此设定所有降雨全部转化为地表径流. 按照不利因素组合考虑, 模型中设定最大单小时降雨与近岸最高风暴潮位遭遇.

图 4 “菲特”雨型20年一遇24 h降雨逐时雨量分布Fig. 4 Hourly distribution of precipitation with “Fate” pattern and 20−year return period over 24 hours

在研究区域遭遇风暴潮与区间暴雨的基础上, 设置黄浦江上游同时发生流域洪水产生溃堤, 即遭遇“风、暴、潮、洪”四碰头的共同影响(工况3). 设置黄浦江沿岸各排水闸门因上游洪水而全部关闭,并在行洪期间出现溃堤情况. 黄浦江上游红旗塘–大蒸塘–圆泄泾河道共设置3处溃口(见图5), 溃口宽度均为90 m, 溃决方式设置为瞬间全溃, 溃口底高程取溃口所在河道堤防外地面高程. 考虑到“菲特”台风期间, 黄浦江上游米市渡水文站达到了历史最高水位, 故黄浦江上游采用“菲特”台风期间的实况洪水. 按照不利因素组合考虑, 设置黄浦江最高水位出现在区域内地面过程积水量最大时刻(第140小时), 并在该时刻发生溃堤.

各计算工况考虑的因素具体见表1.

3 风暴潮作用下金山区积水分布情况

不考虑河网(工况0)和考虑陆域河网(工况1)条件下计算的金山区界范围内总进水量、积水量、河网调蓄水量的变化过程如图6所示, 模型设置在模拟的第134小时开始出现溃口进水. 在工况0中,地面最大积水出现在第150小时; 工况1中, 地面最大积水出现在第140小时. 各工况计算结束时刻均为第240 小时. 未考虑河道时, 金山区范围内风暴潮进水量先增后减, 过程最大积水量为2.81亿m3(第150小时), 计算结束时的总积水量为2.39亿m3. 考虑河道后, 风暴潮进水量变化过程与无河网时略有不同, 总体上呈递增趋势. 地面最大积水量为1.54亿m3(第140小时), 显著小于无河道条件; 河道调蓄量为0.97亿m3,总水量为2.51亿m3, 该工况下风暴潮进水量等于地面积水量与河道调蓄量之和. 计算结束时地面积水量为1.22亿m3, 河道调蓄量为1.45亿m3, 总水量为2.67亿m3. 总体上, 考虑陆域河网后, 金山区风暴潮过程最大进水量从2.81 亿m3减少至2.67 亿m3, 减幅约5%, 由于河道的调蓄作用, 地面最大积水量从2.81 亿m3减少至1.54 亿m3(见表2), 减幅约45%.

图 5 黄浦江上游红旗塘–大蒸塘–圆泄泾溃口位置分布Fig. 5 Location of the Hongqitang to Yuanxiejing breach at the upper Huangpu River

表 1 计算工况表Tab. 1 Calculation conditions

图 6 无河网和有河网条件下金山区风暴潮进水量、地面积水量、河道调蓄量变化过程Fig. 6 The quantity of total inflow, surface waterlogging, and stream channel storage during a storm surge event with (and without) considering the river network in Jinshan district

有无河网条件下陆域最大积水深度分布如图7所示. 考虑河网调蓄后, 由于部分积水进入河网,海塘附近地面积水深度显著减小, 金山区内3 m以上积水区域仅出现在沿海黄姑塘、卫城南河、老龙泉港以南地区, 1.5 m以上积水影响范围由工况0下的到达金山区中部缩小至工况1下的红旗港以南区域; 金山区红旗港以北区域积水深度逐渐降低, 至中运河附近已无较为明显的积水. 有无河网工况下陆域最大积水分布差异如图8所示, 河网调蓄显著减小了金山区中部及南部区域特别是邻近海塘内侧区域的积水深度, 积水深度的减小幅度最大超过1 m. 由于河道的连通性及金山区地形南高北低的特点, 从近海塘一侧地面进入河道的风暴潮水被迅速向河道上游输运, 导致金山区西北角局部区域积水程度略有增加, 增加幅度最大约0.1 m.

表 2 不同计算工况下特征时刻金山区内积水量统计Tab. 2 Statistics on water accumulation in Jinshan District at a characteristic time under different calculation conditions

图 7 风暴潮条件下无河网(左)和有河网(右)陆域最大积水深度分布Fig. 7 Distribution of maximum surface waterlogging during a storm surge event under calculation condition 0 (left) and calculation condition 1 (right)

工况0与工况1下金山区地面最大积水面积和平均最大积水深度统计见表3. 考虑河网调蓄功能后, 由于河网的调蓄能力, 最大积水面积从561.74 km2减至426.39 km2, 减幅约24%, 平均最大积水深度从115.65 cm减至71.29 cm, 减幅约38%, 最大积水深度从7.71 m减至7.70 m, 地面积水程度大幅度减小.

4 风暴潮叠加区间降雨和流域洪水作用下金山区积水分布情况

在考虑陆域河道调蓄功能时, 在无降雨(工况1)及降雨(工况2)条件下, 金山区范围内风暴潮进水量、降雨产流量、地面积水量、河道调蓄量等统计见表2. 有降雨时, 由于风暴潮通过漫堤和溃堤进入陆域之前, 地面已存在一部分积水, 使得进入地面的风暴潮水较无降雨时略有减少, 在地面最大积水时刻(第140 h), “菲特”雨型下风暴潮进水量约为2.47亿m3, 较无降雨时减小0.04亿m3.

图 8 风暴潮条件下有无河网工况下陆域最大积水深度差异(有河网工况减无河网工况)Fig. 8 Difference in maximum waterlogging depth with the river network included and excluded during a storm surge event (i.e., calculation condition 1 minus calculation condition 0)

表 3 不同工况下金山区地面最大积水面积和积水深度Tab. 3 Maximum area and depth of waterlogging in Jinshan District under different calculation conditions

考虑降雨时, 在第140小时, 降雨的产流量为1.14亿m3, 至第240小时为1.17亿m3(见表2). 考虑降雨时, 金山区范围内河道调蓄量在第140小时和第240小时均出现0.06亿m3的增加, 地面积水量在第140小时和第240小时的增幅分别为1.04亿m3和1.14亿m3. 降雨导致区域内总水量大幅提高,在第140小时, 区域内总水量(3.61亿m3)较无降雨时(2.51亿m3)增加约44%; 在第240小时, 区域总水量(3.87亿m3)较无降雨时(2.67亿m3)增加约45%.

“菲特”雨型降雨条件下, 陆域最大积水深度分布如图9(a)所示, 有降雨情况下积水的空间分布总体趋势上与无降雨条件下相同, 积水深度总体增加. “菲特”雨型降雨条件与无降雨条件的最大积水差异如图10(a)所示, 降雨主要造成金山区红旗港以北的陆域积水量增加, 最大积水深度增加超过0.5 m;红旗港以南特别是黄姑塘–卫城河–运石河以南的沿海区域, 在降雨条件下积水深度基本无变化.

不同工况下金山区最大积水面积和平均最大积水深度统计如表3所示. 考虑降雨时, 雨水降至地面导致整个区域全部出现积水, 最大积水面积为603.87 km2, 较无降雨时(426.39 km2)增加46%, “菲特”雨型降雨条件下平均最大积水深度为72.95 cm, 较无降雨时(72.29 cm)增加0.9%, 最大积水深度为7.70 m, 与无降雨时一致. 考虑降雨后区域内总积水程度大于无降雨条件.

在考虑“菲特”雨型降雨前提下, 黄浦江上游红旗塘–大蒸塘–圆泄泾无溃堤和发生溃堤两种情况下的金山区最大积水深度分布如图9(b)所示. 由积水分布可见, 黄浦江上游发生溃堤对金山区内总体积水淹没分布无明显的影响, 仅在溃口邻近局部区域增加了积水深度. 有无溃堤两种条件下模型计算得出的积水差值如图10(b)所示, 溃口2和溃口3造成了向荡港以北、掘石港以东及茹港以西所包含的东西长约10 km、南北宽约3 km的陆域内积水增加0.2 ～ 0.5 m; 溃口1在降雨同期发生的情况下, 使溃口附近陆域积水增加约0.1 m. 由表3统计的积水数据可知, 考虑流域洪水溃决后, 金山区内平均最大积水深度由72.95 cm提高到73.46 cm, 增幅为0.70%, 金山区最大积水面积及最大积水深度较无溃口工况均无明显变化.

图 9 风暴潮仅叠加降雨与同时叠加降雨和溃堤下陆域最大积水范围和最大积水深度分布比较Fig. 9 (a) Maximum area and depth of surface waterlogging during a storm surge event with regional precipitation(calculation condition 2); (b) Maximum area and depth of surface waterlogging during a storm surge event with regional precipitation and flooding (calculation condition 3)

图 10 风暴潮叠加降雨或区间洪水后陆域积水差异Fig. 10 (a) Difference in maximum waterlogging depth with precipitation included and excluded during a storm surge event (i.e., calculation condition 2 minus calculation condition 1); (b) Difference in maximum waterlogging depth with flood included and excluded during a storm surge event(calculation condition 3 minus calculation condition 2)

综上所述, 黄浦江上游红旗塘–大蒸塘–圆泄泾上发生溃堤时, 在“菲特”台风期间实况水位过程条件下, 仅在溃口周边局部区域增加了积水程度, 对金山区范围内的积水程度基本无影响.

5 风暴潮叠加区间降雨和流域洪水对危险性等级的影响

根据《风暴潮灾害风险评估和区划技术导则》[24]的规定, 对区域淹没水深危险性的评价按照不同积水深度分为4级, 如表4所示. 按此划分标准, 根据不同工况下金山区最大积水分布计算了相应的风险等级, 用红(Ⅰ级)、橙(Ⅱ级)、黄(Ⅲ级)、蓝(Ⅳ级)四色标识表征危险性等级大小.

表 4 市(县)尺度淹没水深危险性等级划分标准Tab. 4 Risk level of waterlogging depth at the city or county level

仅考虑风暴潮(工况0)与综合风暴潮、区间降雨和流域洪水, 并考虑河道调蓄(工况3)条件下金山区危险性等级分布如图11所示. 模型在叠加20年一遇区间降雨及“菲特”台风期间黄浦江上游的洪水溃堤后, 由于考虑了河道的调蓄作用, 地面积水较仅考虑风暴潮淹没的工况显著减小, 除靠近海塘的区域危险性等级未出现变化外, 金山区中部和北部大部分区域危险性等级总体减小. 两种工况下不同等级区域的面积统计结果如表5所示. 工况0和工况3相比, 高危险性等级(Ⅰ–Ⅲ级)区域面积减小, 低危险性等级(Ⅳ级)区域面积增加, 金山区危险性等级大幅下降, 其中危险性等级Ⅰ的面积从35.04 km2降至17.99 km2, 降幅约49%; 危险性等级Ⅱ的面积从144.84 km2降至66.92 km2, 降幅约54%; 危险性等级Ⅲ的面积从275.74 km2降至226.41 km2, 降幅约18%; 危险性等级Ⅳ的面积大幅增加, 从95.76 km2增至281.17 km2, 增幅约194%.

图 11 工况0(左)和工况3(右)条件下金山区危险性等级分布Fig. 11 Risk level distribution in Jinshan District under calculation condition 0(left) and calculation condition 3 (right)

表 5 不同计算工况下金山区不同危险性等级面积Tab. 5 Waterlogging area for different risk levels in Jinshan District under calculation condition 0 and calculation condition 3

工况0与工况3条件下金山区危险性等级变化分布如图12所示. 总体上金山区南部、中部和北部大部分区域积水危险性等级无变化或降低1级, 局部降低2级; 西北角部分区域积水危险性等级提高1级, 局部增加2级; 金山化工区由于降雨作用积水危险性等级提高1级. 区域的积水危险性不同等级变化如表6所示, 危险性等级无变化、降低1级和提高1级的区域面积分别为309.23 km2、243.33 km2和23.75 km2, 三者占总面积的95%以上; 变化等级超过2级以上的区域较少, 面积占比不足5%.

图 12 工况0和工况3条件下金山区危险性等级差异(工况3减工况0)Fig. 12 Change in waterlogging risk level from calculation condition 0 to calculation condition 3 in Jinshan District (+ denotes rise, – denotes fall)

表 6 金山区工况3较工况0危险性等级变化面积及占比Tab. 6 Statistics on changes in risk level for waterlogging area in Jinshan District from calculation condition 0 to calculation condition 3

综合考虑河网调蓄及风暴潮叠加区间降雨和流域洪水的情形相比仅考虑风暴潮的情形, 金山区危险性等级不变或降低1个等级的区域约占85%, 危险性等级提升1个等级的区域约占全区总面积的10%. 由结果可见, 考虑河网调蓄的影响后, 区域内整体淹没分布出现一定程度改变, 金山区西北角危险性等级提高, 中部和北部大部分区域危险性等级降低, 结果较仅考虑风暴潮情形所得到的危险性等级更符合实际.