蓄滞洪区洪水演进过程数值模拟及洪灾损失分析
周洁 董增川 朱振业 陈序
摘要：基于实测断面资料建立了研究区的一维水动力模型，基于高精度DEM以及1：10000地形图建立了研究区的二维水动力模型，并用MIKE FLOOD将一维模型和二维模型进行耦合，构建了洪泽湖周边滞洪区一、二维耦合的洪水演进数学模型。利用2003年历史洪水资料对模型参数进行了率定，并以2007年历史洪水资料进行了验证。以洪泽湖百年一遇设计洪水为模型上边界，二河闸、三河闸以及高良涧闸的现行调度方案的水位-流量关系为模型下边界，对洪泽湖百年一遇设计洪水方案进行模拟计算，当蒋坝水位达到14 33 m时，洪泽湖周边滞洪区开始滞洪，得到开始滞洪后不同时段研究区内各类洪水风险要素的动态分布情况以及最大淹没水深、淹没历时，验证了模型的合理性，可用于蓄滞洪区洪水演算分析。
关键词：蓄滞洪区；洪水演算；洪泽湖；风险要素；MIKE FLOOD
中圖分类号：TV 122 文献标识码：A 文章编号：1672 1683（2017）05-0056-07
蓄滞洪区是我国防洪减灾的重点对象，适时启用蓄滞洪区削减洪峰、分蓄洪水是保障下游和防洪重点保护区域安全的有效措施。为做好防洪减灾工作，我国高度重视防洪规划、防洪减灾工程建设与非工程措施部署。其中，洪水演进模拟是非工程措施的重要组成部分，是洪水风险分析、洪水风险图绘制的重要依据。洪水演进模拟的研究方法主要有差分法、有限元法、有限体积法等。对于洪水数值模拟，国内外已有不少专家利用有限差分法和矩形网格建立了二维洪水演进模型，但是独立的二维模型无法计算获得河道中某个断面、某个时刻的水位和流量等水文要素，尤其在研究区内水工建筑物众多、控制调度复杂的情况下，更是难以高精度地模拟水库调度、溃口等对河道洪水演进的影响。
为了更好更精确的反映决堤后洪水演进的情况，本文拟采用以一、二维耦合的水力学方法，依据水量守恒及水动量守恒原理建立洪水数值模型，模拟洪水的恒定及非恒定演进过程，分为一维水力学法和二维水力学法。其中一维水力学法是指采用一维非恒定流微分方程即圣维南方程组对河道非恒定流的水动力学进行模拟，该方法能够获得每个模拟时刻河道内的水面线、流量等情况，便于检查和验证，也便于过程控制。二维水力学方法采用浅水方程来模拟洪水在淹没区域内耕地、林地、山丘区、道路等演进的情况，与一维数学模型相比，二维数学模型能够提供更加详细的水情信息，如：淹没范围、淹没水深、淹没历时等情况。目前世界上最为广泛应用的丹麦水利研究所（Danish Hydraulic Institute）研发的MIKE软件，包含MIKE11、MIKE21、MIKEFLOOD等多个功能强大的模块，其中MIKEFLOOD模块可以模拟一维河网水动力学系统以及二维的洪泛区，实现一维和二维区域之间自由的水体交换，适用于宏观上的流域控制性工程规模论证分析和流域洪水调度研究等。本文以洪泽湖地区为例，利用MIKE FLOOD建立一、二维耦合的洪水演进数学模型进行模拟研究。
1研究区概况
洪泽湖是淮河流域最大的湖泊型水库，地处苏北平原中部偏西，位于淮河中下游结合部，承泄淮河上、中游15.8万km3的来水，总库容135亿m3。地貌上，湖西、湖南为低山岗阜，湖西岗陇和洼地宽窄不等、“三洼四岗”高低相间；湖南区蒋坝至盱眙县城是连绵的低山，湖东、湖北由河湖冲刷堆积而成的平原，地势低下。洪泽湖整体呈西北高东南低的形态。洪泽湖周边滞洪区位于洪泽湖大堤以西，废黄河以南，泗洪县西南高地以东，包含盱眙县的沿湖、沿淮地区，地面高程在16.83 m（1985国家高程，下同）以下，涉及淮安、宿迁两市的淮阴区、洪泽县、盱眙县和宿城区、泗洪县、泗阳县六县区部分地区，还有省属洪泽湖、三河农场，滞洪区面积2132 k m3。
洪泽湖入湖河流主要在湖西部，有淮河、怀洪新河、新老濉河、新老汴河和徐洪河、安东河等，在湖北侧入成子湖的河流有古山河、五河、肖河等，南侧主要入湖河道为维桥河、高桥河，淮河入量占流入总量的70%以上。区域内地势外高内低，唯一排水出路即洪泽湖湖区，洪泽湖洪水由三河闸、二河闸、高良涧闸分别经入江水道、入海水道、分淮入沂、灌溉总渠入江入海。区域防洪工程包括外围防洪、内部防洪工程，主要为迎湖挡洪堤和河道堤防。入湖控制建筑物主要有团结闸、高松河闸、古山河闸、西民便河闸、安东河闸、濉河闸、老汴河闸等。出湖控制建筑物主要是三河闸、二河闸、高良涧闸；区域内防洪工程主要有滞洪区间内部河道堤防和排涝涵闸、泵站。详见图1。
2模型构建
2.1模型控制方程
描述一维水流运动的Saint-Venant方程组是建立在质量和能量守恒的基础上的，以水位和流量为研究对象，其表达式为：
2.2河道一维模型
一维河网的洪水运动用Saint-Venant方程组描述，其上、下游边界的控制条件一般采用水位过程控制、流量过程控制、流量～水位关系控制等形式。由基本方程Saint-Venant方程、边界条件和初始条件共同组成一维洪水运动的定解问题。另外，在天然河道中普遍存在着诸如支流交汇、集中分（入）流、洼地蓄水、断面突扩（缩）、堰、闸等，在这些局部地区，由于水流受固体边壁的影响，水流流态急变，Saint-Venantt方程组不再适用，必须重新根据守恒定律，补充必要的计算条件，这类计算条件是位于域内的物理条件，称之为内边界条件。
洪泽湖周边滞洪区内河网纵横交错、水系复杂，着重提取对洪水分析比较重要的河流沟渠来构造河网模型，包括主要通湖河道、出湖通道以及圩内河道等。收集了128条主要河道的实测断面资料，对圩内河道进行一定的概化，概化后的一维河网共有500多个河段，见图2。模型以淮河、怀洪新河、徐洪河、新汴河、新濉河、老濉河等主要通湖河道的入湖流量作为河道的上边界，以三河闸、二河闸和高良涧闸的水位或者水位一流量关系作为下边界条件。endprint
2.3洪泛区二维模型
2.3.1线状地物处理
区内线状地物处理主要为圩区堤防与重要道路的处理。由于原始道路节点间距不规则，划分网格时容易产生小网格，对模型计算不利，为此，对道路节点进行抽稀处理。道路抽稀长度控制为200 m一个节点，然后根据道路测量高程点对道路进行打断，并将高程点作为分段道路的高程，之后对交汇点进行修正，使所有交汇点都完全拟合，作区内挡水建筑物处理。本次研究将县道及县道以上道路均纳入模型中，合计49条，对县道以下有实测资料的道路经过修正后纳入模型，区域共概化道路1634段。洪泽农场是研究区内面积最大的圩区，其内部道路概化情况见图3。
2.3.2网格剖分
研究区二维地形采用1：10000地形图以及5 m×5 m的高精度DEM资料，采用非结构化三角形网格，网格的大小随地形地势和阻水建筑物的分布灵活确定，充分反映计算区域的地形特征。对区域内的典型的线性阻水建筑物，如堤防、公路等，经合理概化，并对网格适当加密，在二维地形中充分反映其特征。最大网格面积不超过0.1km2，其中湖面网格适当放大，最大网格面积不超过1km2。共剖分网格85 894个，见图4。
2.4一二维耦合模型
MIKE FLOOD是把一维模型和二维模型连接在一起，进行动态耦合的模型系统，耦合模型既利用了一维模型和二维模型的优点，又避免采用单一模型时遇到的网格精度和准确性方面的问题。将一、二维模型通过MIKE FLOOD进行耦合计算，耦合点使用标准连接或侧向连接进行计算。其中，标准连接是将连接线映射到一个或多个二维网格的耦合线上，一维模型为二维模型提供流量值Q作为二维模型的边界条件，将Q值分布到二维计算单元的各节点上；在连接处二维计算网格的水位值并不相等因此取各个计算网格的平均水位值Z返回给一维模型，以进行下一时段的计算。侧向连接允许二维模型的网格单元从侧面连接到一维模型的部分河段甚至是整个河段，利用建筑物的流量公式来计算通过侧向连接的水流。
本次模型的耦合方式主要采用MIKE FLOOD提供的标准连接和侧向连接，其中，通湖河道与洪泽湖湖区、进洪口门与洪泛区主要采用了标准连接；实测河道与洪泛区主要采用了侧向连接。
3模型参数率定验证
3.1參数选取
糙率是表征河道底部、岸坡和洪泛区地表影响水流阻力的综合系数，是水力计算的重要灵敏参数，也是水动力数学模型中最重要的参数，包括一维河道糙率和二维洪泛区糙率。根据水普下垫面信息，确定不同区域的糙率值，将下垫面数据导入模型，创建糙率分区。根据项目区地形、地貌和河道实际情况确定初始糙率值。
3.2模型率定
2003年淮河流域发生了1954年以来最大洪水，洪泽湖蒋坝水位最高达14.20 m，采用2003年6月26日-7月31日的实测降雨和洪泽湖的出湖水位及蒋坝等四个水位站的水位资料对模型进行率定，模型计算范围内水位站分布见图5。将实测入流过程和洪泽湖出湖水位过程输入模型进行模拟计算，将计算的洪泽湖水位过程与实测水位进行对比分析，来率定参数。糙率率定结果见表1。
根据实测资料和模型计算值对比可知，计算的蒋坝、临淮头、香城庄、尚嘴水位及变化过程与2003年实测水位基本一致，见图6。
3.3模型验证
根据之前率定确定的模型参数，采用2007年的实测入流和洪泽湖的出湖水位及蒋坝、临淮头、香城庄、尚嘴四个水位站的水位资料对模型进行验证。选取2007年7月1日-8月3日作为计算时段，将计算的洪泽湖水位过程与实测水位进行对比分析，来检验模型的准确性与合理性。
根据实测资料和模型计算值对比可知，计算的蒋坝、临淮头站、香城庄站、尚嘴站水位及变化过程与2007年实测水位基本一致，见表2。综合以上率定验证成果可知，所建模型较好地重现了2003年及2007年研究区内洪水演进过程。
4滞洪模拟计算
以现状工况遇洪泽湖百年一遇入流情况进行洪水模拟计算，模型的上边界为1996年淮委规划设计研究院的入湖设计洪水流量过程，下边界为洪泽湖出口的水位一流量关系曲线，见表3。根据相关调度规则，当洪泽湖蒋坝水位达到14.33 m且有上涨趋势时，洪泽湖周边滞洪区一次性破圩滞洪。利用之前构建的一、二维耦合模型对洪泽湖周边滞洪区进行滞洪模拟，不仅可以得到最终时刻的淹没要素，还能得到洪水演进过程中，任意网格任意时刻的洪水淹没要素，如淹没水深、流速、流向等。
洪泽湖百年一遇入流情况及模型模拟计算得到的总入流、总出流及蒋坝水位过程见图7，模型模拟时长为42 d，最大入湖流量21703 m3/s，最大出湖流量为18 160 m3/s，总计入湖水量为546 3亿m3，总计出湖为500 9亿m3/s，蒋坝初始水位为13.33m，最高水位达15.16 m；当模型运行至第7 d时，蒋坝水位超过14 33 m时，滞洪区开始滞洪；当模型运行至第21 d时，滞洪区达到最大滞洪量，最大滞洪量为17.63亿m3。
不同时刻研究区内淹没水深图见图8。由于研究区内圩区众多，每个圩区都有独立的进洪口门，概化的314个圩区在达到滞洪条件，即洪泽湖蒋坝水位达到14.33 m时，同时开始滞洪，在整个研究区最初会出现多个淹没区域，随着时间的推移，区内被洪水淹没的面积逐步变大，淹没水深也逐步上升。当t=170h时，蒋坝水位14.33 m，达到滞洪条件，开始滞洪。图8（b）是开始滞洪12 h后的淹没情况，由于314个圩区同时破圩滞洪，滞洪量迅速增加，蒋坝水位不仅减缓了上升趋势，还有些许下降，此时蒋坝水位14.27 m，滞洪量3.27亿m3，淹没面积3.11 km2；图8（c）是开始滞洪后24 h的淹没情况，此时蒋坝水位14.25 m，滞洪量达到4.80亿m3，淹没面积400.50 km2；随着入湖流量的持续增加，滞洪区水位不断上升，滞洪量渐渐达到饱和，削峰效果已不如刚破圩滞洪时显著，洪泽湖水位持续升高。当t=504 h时，蒋坝水位15.08 m，滞洪区达到最大滞洪量17.63亿m3，淹没面积940.02 km2。
通过模拟计算，该区域最大淹没面积940.02km3，最大淹没历时35 d。不同淹没水深对应的淹没面积、不同淹没历时对应的淹没面积见表4、表5。由表可知，随着淹没水深、淹没历时的增大，淹没面积均相应增大，淹没水深、淹没历时在研究区内分布合理。
5结论
为了更准确的反映决堤后洪水演进的情况，采用一、二维耦合的水力学方法，建立洪水演进数学模型。其中，河道一维模型能够获得每个模拟时刻河道内的水面线、流量等情况，洪泛区二维模型能够提供更加详细的水情信息。
利用2003年、2007年历史洪水资料对模型中主要参数糙率进行率定验证，将实测资料和模型计算值对比分析合理性，结果表明，水位变化趋势、最高洪水位到达时间基本一致，具有较高的模拟精度，能够反映模型模拟的准确性。
对现状工况遇洪泽湖百年一遇入流的设计方案进行模拟计算，得到各风险要素在研究区内的分布情况，验证了模型的合理性，可用于蓄滞洪区洪水模拟及风险分析。endprint
