基于洪水演进模拟的溃堤洪灾损失评估

李文欢,支欢乐,蒋水华*,雷 声,黄 河

(1.南昌大学工程建设学院,江西 南昌 330031;2.中铁水利水电规划设计集团有限公司,江西 南昌 330029;3.江西省水利科学院,江西 南昌 330029)

洪水作为世界上最常见且最具威胁的自然灾害之一,严重影响了人们生产生活,制约了中国经济社会的发展。中国水旱灾害公报资料表明[1],2019年,中国因洪涝灾害共造成4 766.6万人受灾,658人死亡失踪,直接经济损失1 922.7亿元。因此,开展洪水灾害风险评估,进而建立洪水灾害预警系统,制定有效的防灾减灾措施具有重要的工程应用价值。

目前,众多学者在洪水演进分析方面开展了大量有益的研究[2-4],并发现相较于Delft3D[2]、HEC-RAS[3]和Delft3D-Flow[4]等软件,MIKE 21软件不仅可以考虑风速和入渗等多种因素的影响,而且模拟结果更准确,计算效率更高。例如,刘卫林等[5]针对水文资料缺乏的中小河流通过小流域水文计算和MIKE 11水动力模型进行研究,提出了临界水位预警指标确定方法。在水动力学模型方面,詹明强等[6]研究了溃坝对下游淹没范围和程度的影响,建立了一维非恒定流水动力学模型。

洪涝灾害损失评估主要分为生命损失、经济损失和生态环境损失这3类[7]。对于生命损失评估,目前主要参考了国外经验并结合国内情况估算洪灾人口死亡率和生命损失[8]。郭磊等[9]通过夜间灯光的遥感数据和土地利用类型推求人口空间展布数据,再叠加历史洪水模拟得到受灾人口数量。经济损失评估常分为直接和间接损失评估,如Pariartha等[10]采用MIKE FLOOD计算洪水淹没水深和淹没面积,叠加分析住宅地理信息,再根据住房损失曲线和平均房价来评估洪水造成的经济损失。目前关于洪水造成的生态环境损失研究较为欠缺,这方面的成果不多,国外则主要以定性方法评估洪水对于生态环境的影响。近年来发展迅猛的地理信息技术受到了学者们的青睐,与洪水灾害损失评估结合进行洪水淹没空间分布数据分析和承灾体空间信息灾害评估[8-10]。

采用MIKE 21软件模拟溃堤洪水演进过程,根据土地利用类型不同合理确定下垫面糙率值,并基于Python对ArcGIS进行二次开发,在此基础上建立具有明确物理意义的损失定量评估模型,进而量化经济财产价值和生态系统服务价值,结合溃堤洪水时空分布数据,建立洪水灾害损失评估模型。最后,以鄱阳湖区某重点实际工程为例验证方法的有效性,并评估溃堤洪水造成的生命、经济和生态环境损失。

1 溃堤洪水演进模拟方法

水动力学模型常被用于模拟洪水过程演进过程。MIKE 21软件是目前国际上流行的水动力学模拟软件[11],常用于模拟溃堤/溃坝洪水演进过程。Symonds等[2]分析对比Delft3D和MIKE 21的模拟效果和计算效率,表明两者的模拟结果都与实测数据吻合良好,但是MIKE 21在多核运行时使用非结构化模型,计算效率更高。Shrestha等[3]分别采用HEC-RAS和MIKE 21对比分析了布伦特伍德地区Deer Creek河流的洪水演进过程,表明两者都可以精确模拟洪水演进过程,但是相对于HEC-RAS,MIKE 21还可以考虑风速、入渗等因素,更真实地模拟实际洪水过程。Parsapour等[4]分别采用Delft3D-Flow和MIKE 21进行水动力建模,在同等网格分辨率的情况下,MIKE 21会获得更准确的模拟结果。MIKE 21采用的二维水流模型基本方程组包括二维流水流连续方程与水流动量方程,其中二维流水流连续方程见式(1):

(1)

式中ξ——水面到基准面的高度,m;t——计算时间,s;u——流速在x方向上的分量,m/s;v——流速在y方向上的分量,m/s。

x和y方向上的二维水流动量方程见式(2)、(3):

(2)

(3)

式中h——静止水深,m;g——重力加速度,m2/s;V——风速;Vx——x方向的风速矢量,m/s;Vy——y方向的风速矢量,m/s;C——阻力系数;f——风摩擦因素,f=γa2ρa;γa——风应力系数;ρa——空气密度;Ω——科氏力系数,与纬度相关。

2 溃堤洪水灾害损失评估方法

基于研究区域洪灾承灾体类别、价值和空间分布情况,构建空间人群分布、经济价值量化、生态系统服务价值分布数据空间化拓扑数据库,发展可融合洪水时空分布数据的灾害损失评估方法,并结合洪水淹没的时空分布特征,动态评估溃堤洪水所造成的生命、经济和生态环境损失。对于生命损失,一方面需要考虑人群的空间分布数据,另一方面还需要考虑洪灾的生命损失率,采用的生命损失LOL计算见式(4)[12]:

(4)

式中 PAR——风险人口数量,即研究区域内居民受洪水影响的总数量;a——人群分布数据的像元网格总数量;IRi——第i个像元网格上的风险个体生命损失率,计算见式(5)、(6)[12]:

IR=f0αβ

(5)

α=qm1+(1-q)m2

(6)

式中f0——中国风险人口死亡率;α——洪灾严重程度系数;β——修正系数,一般取值为1.4;m1——灾难程度主要影响因子;m2——灾难程度次要影响因子;q——权重系数,取值为0.8。

其中,m1可视为各个主要影响因素对风险人口死亡率的影响程度与其权重系数的乘积,见式(7)[13]:

(7)

式中si——第i个主要影响因素对风险人口死亡率的影响程度,主要包含风险人口数量、溃堤洪水严重性SD、预警时间WT以及对洪水严重性的理解程度UD;θi——权重系数,分别取值为θ1=θ2=0.2,θ3=θ4=0.3。

m2为各个次要影响因素对风险人口死亡率的影响程度与其权重系数的乘积,见式(8)[13]:

(8)

式中ti——第i个次要影响因素对风险人口死亡率的影响程度;ni——ti对应的权重系数值。

经济损失评估主要包括直接经济损失和间接经济损失两方面。其中,直接经济损失采取损失率方法计算,通过建立各类财产洪水损失率,结合洪水淹没特征、经济损失统计指标以及淹没水深等级,进行叠加计算见式(9)[14]:

(9)

式中D——总直接经济损失;Sp——第p类财产的直接经济损失;e——研究区域财产分布在二维空间上的像元网格数量;h——第j个像元网格上的财产种类;l——不同淹没水深等级;βijk——第i个像元网格上第j类财产在第k种淹没水深对应的损失率;Wijk——第i个像元网格上第j类财产对应的财产价值。

目前间接经济损失R主要通过洪灾间接系数进行确定,根据直接经济损失与间接经济损失的关系计算见式(10)[14]:

(10)

式中Ki——第p类财产对应的间接系数,农业取值为15%～28%,住宅区取值为15%,公路取值为25%。

对于生态环境损失计算,改进了文献[15]方法,见式(11):

E=c×ESV

(11)

式中c——经验损失率,取c为0.1;ESV——被洪水淹没区域的生态系统服务价值;E——按经验系数法计算的生态环境损失。

3 洪灾损失评估系统实现

为定量评估溃堤洪水灾害损失,借助面向对象的计算机程序语言——Python,调用Arcpy站点包实现内嵌于ArcGIS的内置插件条,基于堤防保护区域人群空间分布数据、经济财产价值量化空间分布数据和生态系统服务价值空间分布数据,有机结合洪水时空分布特征,建立洪水灾害损失动态评估模型,实现对生命、经济以及生态环境损失的动态定量评估。

在损失评估模型的实现过程中,洪水分布数据和承载体空间分布数据的类型及精度是不一致的。在后续的损失计算过程中,要求数据可以相互叠加计算,需要保证洪水数据和承载体数据的投影坐标系、栅格像元大小以及数据类型等一致,所以需要提前进行数据栅格标准化处理。具体步骤如下:①根据溃堤洪水演进模型,提取未被标准化的洪水数据;②将洪水数据进行面转栅格处理,然后将栅格化的洪水数据与承载体数据空间分布数据进一步进行标准化处理。

对于生命损失评估模块,在数据标准化处理模块的基础上,实现动态生命损失评估的具体步骤如下:①输入洪水的时间尺度,本文输入时间是溃堤洪水演进过程中的模型起始时间和末尾时间,设定的间隔时间为2 h,生命损失模块可根据间隔时间循环读取淹没水深数据和流速数据;②设定模型参数,计算个体风险生命损失率IR分布情况,累加计算所有像元网格的个体风险生命损失率总值,叠加计算人群数据,将获得的生命损失值统计情况存入表格中。

对于经济损失评估模块,洪水淹没是1个动态变化的过程,某一瞬时的洪水数据仅能分析这个时间点的经济损失空间分布情况。经济损失评估具体步骤如下:①为实现经济损失动态计算,输入溃堤洪水演进过程中起始、末尾时间以及间隔时间,表示每过一次间隔时间计算一次经济损失;②输入经济分布数据和淹没水深数据,而财产种类数据和损失率关系可根据蓄滞洪区或者防洪保护区的实际社会情况进行修改;③输出各类财产损失的统计表以及经济损失的分布数据。

对于生态环境损失评估模块,根据式(11)可知,需要输入生态系统服务价值空间分布数据。生态环境损失评估具体步骤如下:①输入始末以及间隔时间;②通过叠加洪水淹没的面积筛选受影响的生态系统服务价值空间分布数据,设定经验损失率,计算生态损失;③输出生态损失统计表。

本文建立的洪灾损失动态评估模型有效利用了溃堤洪水演进模拟获取的洪水空间分布数据,其主要优势如下:①基于ArcGIS平台构建空间化拓扑数据库,并且与国家数据统计年鉴内容紧密结合,保证了数据有效性,同时解决了数据的及时更新问题;②考虑了洪水和承载体数据的空间分布差异特征,有效提升了洪水损失评估结果准确性;③通过有限的输入和输出以及模型参数修改反映了洪灾损失动态变化情况;④可以根据研究对象的洪灾损失典型资料,调查分析和调整洪灾损失率关系,兼顾了模型的普适性;⑤实现了溃堤洪水演进模拟、空间化拓扑数据库动态连接以及洪灾损失的实时评估。

4 工程应用

鄱阳湖区某重点堤防地处江西省九江市永修县境内,圩堤起于永修县城郊的观音岩,止于杨公脑,全长33.57 km,保护面积56.28 km2,保护耕地5.03×104亩(1亩约等于666.67 m2,下同),常住人口2.34万人。据2005年现场勘察,三角联圩各类险工险段多达6处,累计堤线长30多km,大多堤段的断面较为单薄且内外坡均较陡,汛期水流逐年冲刷,极易导致堤身渗透破坏,存在部分堤段填筑压实不均匀且高程不足,脱坡险情突出。2020年7月12日19时40分,溃口宽度达20余m,溃堤后,圩内平均水位在0.5 h内上涨0.5 m多,由于洪水冲刷,溃口迅速扩大,截至7月13日10点,三角联圩溃口宽度已经扩大到200余m,并于7月14日14时启动封堵作业,7月16日21时43分,三角联圩完成合拢。淹没耕地面积5.54×104亩,水产业面积1.05×104亩,淹没5 100余栋房屋,其中倒塌了40余栋,溃口经过55 h的抢修后合拢。

4.1 模型计算范围确定

基于中国大地坐标系,建立堤防保护区域数字高程模型见图1。采用三角形网格进行模型剖分见图2,将东南部堤防溃口(红色边界)设置为入流边界,北部三角联圩和西部陆地(绿色边界)设置为闭边界。

图1 堤防保护区域数字高程模型

图2 堤防水动力学模型网格剖分

本数值模型共存在4个片区,分别为堤防南部片区、北部片区、东部片区和西部片区,在每个片区内分别设置4个特征点(即自然村)。假设决口发生在三角联圩南部片区下方见图3,其具体地理位置为九江市永修县三角乡、永丰垦殖场和南昌市新建区大塘坪乡交界,桩号27k+880～28k+010 m处。

图3 三角堤防区域特征点和溃口位置

4.2 模型参数设置及有效性

提出的洪水演进模拟方法涉及的物理参数(包括边界条件、网格形式、初始条件、涡黏系数)见表1。

表1 洪水演进数值模型主要物理参数取值[16]

本模型中土地利用类型主要包括旱地、水田、草地、滩地、居民区和水面6种,其相应糙率参数按照下垫面分布情况,并参照相关经验进行率定选取[19],见表2。为了验证模型的有效性,图4比较了采用所建模型模拟的溃口淹没水位过程线和实测的溃口区域内淹没水位过程线。由图4可知,在水位基本不变后,水位模拟值为22.65 m,水位实测值为22.58 m,两条水位过程线较为吻合,表明基于MIKE 21软件能够合理模拟溃堤洪水在堤防保护区域内的演进过程。此外,所建模型模拟的淹没面积为56.15 km2,与堤防溃口合拢后的实测淹没面积(56.28 km2)相差仅0.13 km2,进一步说明了本模型的有效性,可以满足计算精度要求。

表2 糙率参数参考取值

图4 溃口区域内模拟水位与实测水位的比较

4.3 溃堤洪水演进过程分析

图5给出了三角联圩溃口流量随溃决时间的变化情况。2020年7月12日19时为图5中的起始时刻,在7月12日19时40分堤防溃决,随后堤防溃口宽度逐渐扩大至最大值200 m,溃堤后的21 h内流量从0急剧增加到最大值1 724 m3/s,随后流量开始下降,于7月14日14时启动封堵作业,7月16日21时43分,三角联圩溃口完成合拢,流量接近于零,整个堤防保护区域进洪过程在98 h内完成。

图5 三角堤防溃口流量过程线

为了进一步展示溃堤洪水演进的淹没过程,图6给出了6个典型进洪时段的淹没水深分布。洪水在决堤后10 h内到达堤防南部片区,虽然红旗村处于溃口附近,但是因红旗村的地势较高,在溃决40 h内受灾情况较轻。由于堤防北部片区地势较低,洪水在20 h内向周家湖、三角圩等村庄挺进,溃堤40 h后,堤防北部片区基本被淹没,洪水逐渐向东西部蔓延,堤防保护区域内水深高达5.2 m。溃堤60 h后,堤防北部片区淹没水深最高达6.0 m;溃堤80 h后,因堤防西部片区比东部片区地势较高,故东部区域淹没深度最深达到6.5 m,而西部区域仍旧有部分地区未被淹没。进洪100 h后,堤防保护区域绝大部分被淹没,此时需要尽快对圩内居民实施避险转移抢险救援等举措。

a)10 h

4.4 溃堤洪水损失评估

根据收集的洪水影响对象统计数据,对人群数据进行空间展布,并量化经济财产和生态系统服务价值。然后,搭建空间化拓扑数据库来反映生命、经济、生态环境在空间上的分布差异,接着在溃堤洪水演进模型中提取洪水淹没时空分布数据,并将数据栅格进行标准化处理。最后,进行包括生命、经济和生态环境损失在内的综合损失动态评估。图7给出了溃堤造成的生命损失随溃口历时的变化曲线。由图7可知,随着预警时间的延长,生命损失不断减少。同时,居民对洪水严重性的理解程度也对生命损失具有重要的影响因素。例如,在预警时间不超过0.25 h且居民对洪水严重性理解程度为模糊的情况下,生命损失高达454人,见图7a;相较之下,如果预警时间足够长,且居民对洪水严重性的理解程度为明确的情况下,生命损失骤降至2人,见图7b。

a)居民模糊理解洪水严重性情况下生命损失曲线

图8给出了溃堤造成的经济损失随溃口历时的变化关系曲线。由图8可知,直接经济损失随着溃口历时的增加而增大,其中居民住房直接经济损失最大,高达2.86亿元,其次为家庭财产直接经济损失。同时,图9给出了堤防保护区域内总经济损失分布情况,不难看出居民用地损失情况最为严重。

a)耕地-家庭财产-住房直接经济损失

图9 溃堤造成的总经济损失分布

5 结论

以江西省鄱阳湖区某重点堤防为例,建立了基于MIKE 21软件的洪水演进数值模型,获得了溃堤洪水下淹没水深和流速等洪水空间分布数据,丰富了溃堤洪水演进模拟方法。同时,考虑生命、经济、生态环境的空间分布特性,基于Python平台建立了一种具有明确物理意义的,可考虑生命、经济和环境损失的定量动态评估模型。主要结论如下。

a)根据研究区下垫面土地利用情况设置不同的糙率值,以及与实测的溃口处淹没水位过程线进行对比,验证了洪水演进数值模型的有效性,表明所建模型具有较高的计算精度,可有效用于获取淹没水深、流速等洪水空间分布数据。

b)发展了溃堤洪水灾害损失评估方法,并采用Python语言编写了基于ArcGIS的嵌入式耦合插件工具,建立了洪灾损失动态评估模型,可以有效评估不同工况下溃堤洪水造成的生命、经济以及生态环境损失,为洪水灾害防控提供了有效途径。