基于耦合水动力模型的广州市东濠涌流域洪涝模拟

2019-08-06罗海婉陈文杰李志威黄国如

水资源与水工程学报 2019年3期

罗海婉，陈文杰，李志威，潘健，高强，黄国如,3,4

(1.华南理工大学土木与交通学院，广东广州 510640； 2. 广州市水务科学研究所广东广州 510220；3.华南理工大学亚热带建筑科学国家重点实验室，广东广州 510640；4. 广东省水利工程安全与绿色水利工程技术研究中心，广东广州 510640)

1 研究背景

模拟和预报城市洪涝是防洪减灾的重要基础，也是当今城市水文学研究的热点课题[1-3]。国内外研究学者对城市区域的产汇流特性和洪水计算方法做了大量研究并与计算机技术结合，基于现有的水文水动力知识基础，利用计算机技术建立数值模型对城市洪涝过程的地表产汇流、管道汇流、河道行洪和地表淹没进行模拟和计算，分析淹没区域和水流水位等情况，可为城市防洪排涝、抢险救灾等提供决策依据[4-5]。

城市洪涝模型可分为一维模型、二维模型和一二维耦合模型3种。一维、二维模型均可对城市洪涝过程进行水文水动力过程模拟，但一维、二维模型均有各自的优势与局限性，前者计算效率较高但不具备模拟洪水在地表演进过程的能力[6]，后者可模拟二维地表水流却不具备模拟地下管网汇流过程能力[7]。因此，将一维、二维模型进行耦合可克服两者的缺点、发挥两者的优势[8-9]。目前，常用于一二维城市流域洪涝模拟的模型有华霖富水力研究有限公司(Wallingford)的Infoworks、丹麦水利研究所(DHI)的MIKE URBAN等[10-11]，以上两个模型均为使用成本较高的商业软件。随着国内外学者的研究，近年来一维、二维耦合模型有了较大的发展。Schmitt等[12]通过耦合一、二维模型研究城市洪水特征。耿艳芬[13]基于水动力学方程组，建立可准确模拟城市流域地下管网模拟与地面淹没模拟的一、二维模型及其耦合模型。Seyoum等[14]将SWMM模型与自主开发的二维地表模型耦合，克服了地表和地下管网水流交换的困难，实现了真正意义上的一、二维水流交换。李传奇等[15]基于一维非恒定流方程与二维自由水面流动模型构建一二维水动力耦合模型，并应用于济南市地区洪水模拟中。陈鹏宇[16]基于SWMM与自主开发的二维漫流模型，构建城市流域洪涝一维、二维耦合模型系统，并通过济南市洪涝算例验证模型。栾慕等[17]构建SWMM-MIKE11耦合模型，对桐庐县的管网排水能力和城市内涝风险进行了评估。王庆平等[18]对一、二维水动力模型进行耦合并改进，应用于山区半山区流域的洪涝灾害预警。

由上可知，将一维、二维水动力模型耦合可充分发挥两者的优势，对一维水动力模型进行二次开发也是研究城市洪涝模拟问题的首选技术路线[19]，本文将一维模型SWMM(Storm Water Management Model)与二维水动力模型进行耦合构建耦合水动力模型，克服SWMM无法模拟地面水流淹没与流动的局限性，将验证后的模型应用于广州市东濠涌流域的洪涝过程模拟并分析该模型在城市流域的适用性，取得的模拟结果表明模型可靠性和精度较好。

2 研究区域概况及积水原因

东濠涌是珠江广州段主要河涌之一，总长约4.4 km，发源于白云山麓湖，流经越秀区汇入珠江。东濠涌流域水系主要由麓湖、东濠涌和新河浦涌组成，全流域面积12.4 km2。麓湖的洪水和流域的暴雨实行错峰泄洪，东濠涌流域防洪排涝主要依靠东濠涌及其支涌新河浦涌两条河道和地下排水管网，因此根据研究区域实际汇流情况对东濠涌流域边界进行调整，扣除麓湖及其汇水面积后，东濠涌流域总面积为10.3 km2，如图1所示。

近几十年的社会发展中，广州市的土地利用复杂性不断上升[20]，城市化的快速进程使城市的水循环、水生态和水安全受到一定影响。近年来广州市区年降雨量逐年增加，暴雨中心逐渐向市中心转移，降雨时段更为集中、空间分布不均匀，广州市的暴雨洪涝问题越发严重[21]。城市暴雨洪涝问题是由多种因素共同作用而产生的灾害现象[22]，东濠涌流域位于广州市典型老城区越秀区内，具有截污、防洪、排涝等功能。降雨期间，流域内地面径流通过地下管道、地表汇流进入东濠涌。因近年来流域内建筑物密度变大，管网规划混乱且设计标准偏低等因素，加以东濠涌断面宽度变小、河道淤积情况严重等原因，造成河涌过流能力下降，无法满足日益上升的防洪排涝要求，流域内溢流污染和积水问题严峻。当降雨量较大时，东濠涌流域局部地区容易产生积水淹没现象，轻则影响居民的日常生活，重则造成经济财产损失。通过相关部门提供的信息整理与实地调研，调查积水点的分布并分析其成因，流域内积水点分布与成因分别如图1、表1所示。

3 耦合水动力模型

本文所构建的耦合水动力模型的一维管道模型基础是美国环保署(EPA)资助下研发的SWMM，其广泛应用于城市区域和非城市区域的排水系统的规划、分析和设计[6,23]。SWMM的局限性在于无法模拟二维地表积水的流动，但其模型代码开源易于实现模型耦合。SWMM提供恒定流法、运动波法和动力波法3种方法用于管道流量传输计算，本模型计算考虑回水、蓄变、有压流、逆流等情况，采用动力波法进行求解，即完整求解圣维南方程以求得管道流量。圣维南方程包括连续方程和动量方程，如公式(1)、(2)所示：

(1)

(2)

式中：Q为流量，m3/s；A为过水断面面积，m2；x为距离，m；t为时间，s；H为静压水头，m；v为平均流速，m/s；g为重力加速度，可取9.8 m/s2；S0为重力项；Sf为摩擦力。

图1 东濠涌流域边界和积水点分布图

表1 研究区积水点位置及其成因分析

本文所构建的耦合水动力模型的二维模型基于二维浅水方程，如公式(3)所示。结合Godunov型有限体积方法、隐式双时间步方法以及HLLC格式黎曼近似求解器等，建立一个高效、稳定、时空均具有二阶精度的二维非恒定流水动力学模型，改进后控制方程如式公(4)所示。

(3)

(4)

式中：h为水深，m；u和v分别为x、y方向平均流速，m/s；b为水底高程，m；Sox和Soy分别为底坡源项Sb在x、y方向的分量；Sfx和Sfy分别为Sf摩阻源项在x、y方向分量。

双时间步法具有计算效率高、稳定性好等优点，自提出后便得到广泛的应用[24-28]，本文将此算法应用于二维水动力模型并建立隐式高效的有限体积法模型，可优化模型计算效率与稳定性。二维水动力模型构建完成后，利用恒定流算例、Stocker和Malpasset溃坝算例、斜水跃算例，对该模型处理恒定流、间断流、实际地形等模拟能力进行检验与验证，并将隐式双时间步法结果与显式方法结果进行对比，对比结果表明隐式双时间步法的计算耗时和稳定性明显优于显式方法[29]。模型具有处理实际洪水的能力以及良好的计算精度，能够应用于城市洪涝数值模拟计算[30]。

本文的耦合水动力模型以SWMM为一维模型基础，通过动态链接库文件(DLL)方式耦合二维模型。为建立一维模型和二维模型在水平和垂向方向上的连接，分别采用经验公式、相互提供边界条件、堰流公式和孔口流量公式的方式建立地表侧向连接、正向连接和垂向连接，并在每个连接方向上采用理想算例对连接方式进行分析与校验，计算结果表明模型水量守恒、计算结果精度良好、连接方式合理可行，能很好地预测城市洪涝的发生过程[31]。耦合模型结构示意图如图2所示。

图2 模型耦合结构示意图

4 东濠涌流域洪涝模型构建

4.1 一维模型构建

本文的一二维耦合模型的一维模型基础为SWMM，因此一维排水模型的构建主要包括以下几个步骤：

(1)地下管网和河道数据的处理。CAD图中的地下管网系统较复杂，检查井、管道和河道种类较多，因此需根据构建城市雨洪模型的需要对排水系统进行概化，将雨水篦子、个别短管道等对模型计算精度和可靠性有影响的管网排除。此外，河道有较大的行洪作用，是一维排水模型的重要组成部分，因此需对河道断面数据进行整理和概化。完成上述概化工作后，提取管网和河道信息并构建排水系统数据库，整理、分析、检查和修正管道和河道的数据，最终简化为3 928条管道、河道以及3 471个检查井，如图3所示。

(2)子汇水区划分及其数据提取。东濠涌流域面积约10.3 km2，检查井和管线数量较多，实际产汇流过程复杂，手动划分子汇水区的工作量较大，故需利用GIS通过地形数据初步提取分水岭，再结合检查井的位置采用泰森多边形法划分子汇水区，利用Arc GIS空间分析功能，计算子汇水区的面积、坡度、不透水率等模型构建所需特征值并构建子汇水区数据库，最终共划分7 700个子汇水区，如图3所示。

图3 研究区域一维模型结构图

(3)模型参数设置。模型运行还需设置模型运行参数，如透水区曼宁系数、不透水区曼宁系数及霍顿产流参数等，则需根据SWMM使用手册和邻近地区相关研究成果确定[32-33]，一维模型参数取值如表2所示。

表2 一维模型参数取值

4.2 一维、二维模型耦合

一维、二维模型耦合主要分成以下几个步骤：

(1)确定二维模拟研究范围。根据东濠涌流域的实际情况对研究区域进行分析，在数据允许的范围内尽量与一维排水模型保持一致，因此二维模型构建的范围与一维模型范围相同，均为东濠涌流域，见图1所示。

(2)地表建筑物及河道概化。东濠涌流域为高度城市化的区域，流域内建筑密集，建筑物在实际中对水流有阻挡作用，在模型中的处理思路有两种，一是将建筑物高程抬高到一定的高程，使其在模拟时不出现积水；二是在模拟范围中排除建筑物而进行模拟。本研究采用第二种方法，根据土地利用分布图划出建筑物轮廓并在模拟范围中将建筑物排除在外，只考虑建筑物对水流的阻挡作用，不考虑建筑物内部淹没作用。此外，河道作为重要的泄洪通道，实际中不可积水，因此也将河道排除在二维模型模拟范围之外，建筑物及河道的分布如图4所示。

图4 研究区域建筑物、河道分布与监测点布置图

(3)边界条件确定。由于东濠涌流域为闭合的城市流域，本研究设置二维模型的边界为固边界，雨水只能通过地下管网和东濠涌排到流域外。此外，流域内存在需要进行侧向耦合的河道，须对耦合边界进行标记，便于模型识别耦合边界。

(4)网格划分。二维模型网格划分的最大尺寸为100 m2，在建筑物密集区域、道路区域等重要区域的网格根据建筑物边界和道路边界自动加密，网格最小尺寸可达0.1 m2。以上的网格划分方法既可提高了道路、建筑物等重要区域的计算精度，又保证了模型整体的计算效率。

(5)构建地表高程模型。建立二维模型须引入地面高程模型并将地面高程赋值到划分好的网格之中。本模型的高程点数据由地形图、管网图等包含地形信息的资料中获取，因此高程点数据较为密集，容易出现高程点缺失或错误等高程误差现象，造成局部地区坡度变化过大影响模拟结果的可靠性和精度。研究区域内部高程范围较大，须利用Arc GIS的3D分析功能建立TIN模型，对高程误差点进行修正。

5 结果与分析

为提高模型的精度和可靠性，本文选取东濠涌流域内所布置的雨量站监测所得的实测降雨数据、检查井变化过程、地表淹没水深等实际情况对东濠涌流域洪涝模型参数进行率定和验证。通过整理雨量站监测资料、检查井监测和水深调研数据，本文选定两场实测降雨场次，即采用20170715场次实测降雨进行模拟，根据该场暴雨造成检查井水位变化数据对东濠涌模型进行参数率定并采用20170905场次实测降雨验证该模型。

5.1 实测降雨数据

东濠涌流域范围较大，内含3个雨量站。由于降雨在时空上具有较大的不均匀性，因此在具备多个雨量站的流域内，常采用泰森多边形对流域进行划分[34-35]。本文根据雨量站所在位置，采用泰森多边形法将东濠涌流域划分3个雨量站的监测区域，雨量站的详细信息如表3所示。检查井监测点布置、雨量站布置及其监测区域如图4所示。

表3 雨量站信息表

通过整理雨量站的降雨数据，选出20170715(18:24-18:44)、20170905(13:01-14:00)两场实测降雨用于模型的参数率定及验证，两场实测降雨的降雨过程见图5(a)、5(b)所示。

5.2 实测降雨模拟结果分析

本文所构建的耦合水动力模型可模拟地下管网一维水力特征值、子汇水区蒸发和下渗等水文特征值以及二维地面积水深度和流速等水力特征值等一系列城市洪涝过程特征值的变化过程。由于20170715场次降雨历时较短、降雨强度较小，流域内检查井几乎不发生溢流，即无法利用地表淹没范围及最大淹没水深数据与模型结果进行对比，因此20170715场次降雨选取2个检查井内液位计实测水位与模型模拟水位进行对比研究，检查井分别是位于区庄立交附近的J1和位于麓景路与恒安路路口附近的J2。以20170715场次降雨数据输入模型进行模拟，模拟结果显示流域内仅小部分区域发生内涝，与实际情况相符。检查井的实测水位与模拟水位变化过程如图6(a)、6(b)，由图6可知本文所构建的东濠涌流域洪涝模型的一维模拟结果与实测数据较相符。

图5 实测降雨过程图

图6 20170715场次降雨中检查井J1和J2水位变化过程

与20170715场次降雨不同，20170905场次降雨历时较长且降雨强度较大，流域内较多检查井产生溢流进而发生淹没现象，但检查井内的液位计只能测得检查井地面以下的水位，无法测得检查井溢流后地表淹没深度，即无法利用检查井实测水位与模拟水位结果进行比较。因此20170905场次降雨选取该场次暴雨下调研所得的主要受淹街道的最大淹没深度与模拟所得最大淹没深度进行比较。以20170905场次降雨数据输入模型进行模拟，模拟结果显示流域较大范围发生淹没现象，与实际情况相符。统计本场次降雨受淹街道的调研结果与模拟结果并进行对比分析，如表4所示。由表4可看出，模拟结果总体比调研结果偏大，这可能是由于在二维地表建筑物概化时，将建筑物排除在模拟范围之内，忽略了地表水流对建筑物的实际淹没，致使模拟水位比实际调研水位偏大。但两者差值较小，因此可说明本文所构建的东濠涌洪涝模型的二维模拟结果与实测数据较为相符。

选用两场实测降雨作为降雨情景对所构建的东濠涌洪涝模型进行模拟，通过两场降雨的实测水位、调研淹没深度与模拟结果进行对比分析，可发现无论在检查井的水位变化过程还是地表淹没的范围和深度上，模型模拟的结果与实际情况较为相符，能够反映出在降雨过程中东濠涌流域的一维排水系统的运行状态和二维地表淹没情况，可认为该模型在洪涝过程模拟上具有良好的精度和可靠性。