山区河流倒灌引发管网溢流洪水数值模拟研究

靳 伟，赵军伟，孙 健，黄鹏程

（1.中国电建集团贵阳勘测设计研究院有限公司，贵州贵阳 550081；2.贵州省建筑信息模型（BIM）工程技术研究中心，贵州贵阳 550081）

0 引言

在全球气候变化的大背景下，洪水灾害在世界各地发生越来越普遍［1］。我国目前正处于城镇化高速发展的阶段，城市区域下垫面条件发生了巨大变化，硬化路面增加、河湖水系等自然水体萎缩，城市排水高度依赖地下管网系统，城市“韧性”或对大中型洪水的抵御及灾后迅速恢复能力尚且较弱［2］。

对于内陆城市，常见的洪水类型包括汛期河道洪水、降雨导致的内涝或山洪。非常规的洪水如因上游水库泄洪或溃坝导致的下游洪水事件在近年来也有报道，如2018年山东省寿光县水灾、2020 年8 月重庆主城区大洪水等。因河道水位短时间内暴涨，水流倒灌入管网造成内陆检查井溢流而淹没街道的现象，在山区河道雨季或上游水库泄洪、溃坝期间极有可能发生，若溢流洪水发生在河道洪水冲毁或漫过河堤之前，会给公众的逃生和财产转移带来巨大安全隐患和不利影响，但此类洪水受到的关注较少。

洪水的时空要素信息对于防洪减灾和指导应急疏散具有重要意义。目前，学术界和工业界应用水文水力模型在城市河道防洪规划和雨洪内涝管理方面已有较多成功经验，模型模拟成为一种重要的洪水研究方法［3］。在内涝防治方面，暴雨径流管理模型SWMM［4］是一个被全世界广泛应用的开源免费城市暴雨径流管理模型，其扩展运输模块（EXTRAN）是计算排水管网水力过程的有力工具，径流模块基于半分布式水文学方法，但由于子流域的边界通常不易明确，且该模型难以直接对流域内特定位置的流量和水深进行描述。当基于水文学方法的雨洪模型较难满足城市防洪精细化要求时，部分学者利用水动力模型或者耦合水动力及水文过程来模拟雨洪过程［5-7］，比如使用以MIKE、InfoWorks 系列产品为代表的商业软件来计算洪水演进过程［8-10］。相较于水文学方法，水动力学方法通过直接求解圣维南方程组或其简化形式来描述洪水的地表演进过程，从而获得更加详尽精确的结果，但求解过程通常过于复杂，并且随着计算区域及地形复杂程度增大，计算时间成本会增加。为解决这类问题，一些学者开始探索使用图形处理器（GPU，Graphics Processing Unit）来加速水动力模型［11，12］。在河道洪水模拟和预测方面，大部分学者的关注点主要集中在对河道洪水自身演进规律的描述和对该类洪水常规影响因素的探讨，如上游溃坝、暴雨、下游顶托、闸坝的调度管理等。例如，马利平［13］等采用集成HLLC 近似黎曼求解器的二维水动力模型，模拟研究了支沟清水沟水库溃坝洪水对主河道大理河行洪过程的影响；魏乾坤［14］等利用MIKE 软件对平原感潮河网地区河道洪水对村镇内涝的影响进行了系统研究；刘业森［15］等以深圳市布吉河流域为例，研究了降雨中心和降雨移动两个时空特征对河道洪水过程的影响。对于河道倒灌洪水，既往的研究大多围绕沟渠倒灌或污水截流系统的截流井和泵站溢流管倒灌进行分析，例如，岳志春［16］等基于一二维水动力学数值模拟方法及沟渠倒灌影响评估方法，建立了考虑排水沟渠倒灌影响的溃堤风险分析模型，模拟了黄河宁夏卫宁段遭遇百年一遇洪水时溃堤造成洪水淹没风险，评估了引排水沟渠倒灌影响；张楠［17］等采用MIKE 软件对截流井倒灌与泵站溢流管溢流进行了模拟诊断，分析了晴天河水倒灌量与雨天溢流量，提出景观河道的合理控制水位。

城市排水管网的终端往往是河湖等自然水体，如何有效地对因河道洪水导致的内陆管网溢流进行模拟研究是一个亟待解决的问题。本文针对此类洪水，建立了二维地表漫流模型和地表与管道水量交换模型，为模拟预测和防范该类洪水风险提供有力支撑和借鉴。

1 管网溢流与地表漫流耦合模型构建

本研究涉及的管道溢水、地表洪水漫流分别由包含检查井溢流模型的一维管道流模型和二维地表洪水演进模型构成。

1.1 二维地表漫流模型

二维地表洪水演进模型以非线性长波方程［18-21］为基础公式：

式中：h为地表洪水淹没深度，m；M，N分别为x和y轴方向的线流量，m∕s；Q为检查井内的溢流量，负值时表示地表洪水流入井内流量，m3∕s；∆x，∆y为计算网格的长和宽，m；z为地表高程，m；g为重力加速度；n为曼宁粗糙系数。该控制方程的离散化采用Leap-frog方法。

1.2 地表与管道水量交换

为综合考虑管道内的流动情况和从检查井溢出的水量，对检查井和地表连接处的流量交换以及管道内水头变化进行了单独建模，如图1所示，其特征是考虑了检查井和管道连接处的局部压力差。A为检查井的横截面积，m2；D为管道横截面积，m2；η为检查井内水位，m；L为管道长度，m；v为管道内流速，m∕s；w为检查井内水流垂直流速，m∕s；i为管道和检查井的编号；h为地表洪水淹没深度，m；z为地表高程，m；ρ为水的密度，kg∕m3。当检查井内水位η高于z时，发生溢流；反之，则未无溢流情况。

图1 管道-检查井-地表连接水量交换示意图Fig.1 Water volume exchange of Pipe-Manhole-Land system

1.2.1 检查井内水位溢流

如图1（a）所示，此时检查井上游管道断面处［图1（a）中的P1］的压力水头等于检查井内的全水头压力，断面P1和P2的压力计算公式分别用式（4）和（5）表示：

式中：f0为流速损失系数（取1.0）；fd为分流损失系数。

右侧第二项表示从检查井内溢出水量带来的压力损失。式（5）右侧第三项表示检查井上下游管道断面因压力差造成的流速水头差。管道内流速由式（6）进行计算：

式中：as为分流损失外的形状损失和摩擦损失。

溢流量可由下式计算：

将式（7）所得溢流量Q代入式（1）右侧，即可将管道水量和地表水量联系起来，并求得溢流量导致的地表洪水淹没情况。

1.2.2 检查井内水位未溢流

如图1（b）所示，此时断面P1和P2的压力计算公式分别用式（8）和（9）表示：

管道内流速用式（6）计算，检查井内水位可用连续方程求得：

式（10）中Qi为地表流入检查井内流量，m3∕s，其值可用堰流公式计算：

式中：μ为流量系数；S为检查井周长，m。

2 模型验证

2.1 理想算例实验

为了验证所构建模型能够模拟水流从河流水体通过管道检查井溢流或反向流动，应用图2 所示的装置构建数值模型。水池T1 和水池T2 的横截面面积分别为6 000 m2和3 000 m2，二者通过管道连接，初始条件为水池T1的水位从压力平衡水面线下降20 m，处于静止状态。为了与解析解进行比较，该实验忽略了检查井中垂向流速对水池中压力场的影响，水池底部的压力近似于根据水池内水位计算的静水压力。模型中，水平方向的网格宽度设为15 m，管道计算网格为2 m。模拟结果如图3所示，可以看出水池T1的水位变动在振幅和周期上都可以准确地得到再现。但考虑到验证和对比方便，此算例中模型未考虑摩擦力，因此水位变动在振幅和周期上并无变化。

图2 模拟装置图Fig.2 Modeling configuration diagram

图3 水池T1中水位数值解与解析解比较Fig.3 Comparison between modeling results and analytical solution of water level in Tank T1

2.2 小尺寸管网溢流实验

为验证模型对实际洪水的模拟能力，进行了小尺寸模型实验，实验用的排水管网和陆地区域模型为实际的1∕30 尺寸，如图4 所示。排水管网模型由聚乙烯管制作，管道及检查井布置和设计参考了最新版《室外排水设计规范》GB50014-2021，排水管段M1～M2～M3～M4～M5 为主干管，M6～M7（M7’）～M2～M8～M9和M10～M11～M3为支管。三条管道分别包含4条、5条和2条管段，各管段管长和管径均在图4（a）中标出，检查井均设置为边长5 cm 的矩形截面立方体。垂向布置上，管道距底部距离为0.18 m。为方便通过颜色辨识地表溢流洪水分布，陆地地表用浅色聚氯乙烯板制作，高程为0.7 m，共设置有9栋建筑物。

本实验设置有4 种河道洪水入侵情景，河堤前的水位变化（包括由河堤引起的反射洪水波）如图5 所示。4 种洪水情景的区别主要是洪水位初期上升方式，即前10 s 的最大洪水位和上涨速度，Case1、2、3 洪水位以较陡峭的方式上升，最高水位依次递减，Case4的水位则相对上升缓慢，而4种情景在15 s后，洪水位的变化区别不大。陆地区域发生管网溢流后，在图4（a）所示的10个淹水深度测量点进行水位测量，由于本实验淹水深度很小，因此在测量点附近均安装有刻度尺，通过对其视频拍摄，读取水深深度变化进行测量。图4（b）显示了装置沿主干管的纵断面示意图和高程设置。

图4 小尺寸实验装置图Fig.4 Schematic of the small-scale experiment

图5 河堤前水位随时间变化图Fig.5 Time-series of water level in front of river bank

3 结果及讨论

3.1 溢流量

在对耦合模型进行验证分析前，先对管网溢流模型溢流量计算的可靠性进行确认。图6 所示为Case1 溢流量计算值与实测值比较，二者的结果吻合度较好。Case1～Case4 的溢流量计算值与实测值比较统计如图7 所示，Case1 至Case4 的溢流量依次降低，4 种情景的溢流量与洪水波初始水位及上升速度成正相关。Case1 的洪水位在实验开始10 s 内最先达到并超过20 cm（检查井顶端距河道静止水面的高度），而后Case2 至Case4依次达到和超过该高度值，洪水位持续超过该值的时间越长，造成管网溢流的水量也越大。在比较Case1 至Case3 溢流过程时，发现即使不同情景下总的溢流量变化不大，但由于洪水位初始上升速度（洪水波陡峭程度）不同，也会导致溢流过程中短时间内溢流流量有较大差别，洪水波上升陡峭的情景下，检查井单位时间的溢流量较大。4 种情景下管网溢流量分布特点是，主干管各检查井M2、M3、M4 溢流水量大，临近河堤并靠近主干管的M7、M8检查井溢流量也相对其他检查井较大，这些溢流量较大的检查井也是管网中管径前后变化的关键节点。

图6 Case1的溢流量计算值与实测值比较Fig.6 Comparison between simulated and measured overflow volume for Case 1

图7 Case1～Case4溢流量统计Fig.7 Overflow volume for Case1～Case4

3.2 地表淹没

通过耦合模型对溢流洪水在地表漫流进行模拟分析，模型中地表网格精度为5 cm，曼宁粗糙系数根据聚氯乙烯板材质，取为0.010。图8 所示为Case1 情景下，各淹没水深测量点处的水位计算值与实测值的比较，二者基本一致。但在大约20 s 之后，模拟水深值均略低于测量值，这是因为模型计算中高估了地表返流入检查井的水量。在实验中，由于检查井本身尺寸所限以及地表洪水淹没水深较小（大部分情况下低于0.5 cm），导致当其返流入检查井时，受到地表和检查井井壁表面张力影响，抑制了水流流动，阻碍其流入检查井及在井内水位下降。在模型中，并未考虑该表面张力的影响，因此在水位下降期，模拟水位值比测量值普遍偏小。

图8 Case1地表淹没水深计算值与测量值随时间变化比较Fig.8 Time-series comparison of the simulated and measured inundation depth

此外，在溢流发生的初期，P1 和P2 测量点处，可以观察到溢流和地表淹没水深急剧增大，二者的淹没范围大致一样，但是从检查井溢水的方式略有不同，Case1 相对Case4，其溢流状态更接近于喷涌。因此，河道洪水位上升的速度（洪水波的陡峭程度）对总的溢流水量体积影响较小，但对短期的溢流流量有较大影响。Case1 的初期洪水位上升快速，短期从检查井溢出的水量大，若发生在实际洪水中，大量洪水短时间内从街道检查井喷出，会造成严重后果，妨碍人员疏散避难，甚至对行人和车辆产生直接安全威胁。所以，为正确把握对象区域的管网溢流洪水风险，除考虑洪水位的高低之外，应对洪水波的形态、陡峭程度等因素进行系统分析，在模型中均应予以考虑。

图9对比了不同时刻模型计算和实测的地表淹水范围，选择颜料着色效果最好的Case2 进行比较。M2 和M4 处的溢流淹水较显著，并朝四周扩散，在模型计算结果中得到较好地再现。

图9 地表淹水范围比较（Case2）（单位：cm）Fig.9 Comparison of flooding area of the simulated and small-scale experiment

4 结论与展望

针对山区河道汛期雨季或因上游水库泄洪造成的水位快速上涨，导致临江管网发生溢流洪水的现象，建立了定量评价检查井溢流与地表洪水淹没分布的数值模型。通过理想算例和小尺寸管网溢流实验，验证了模型计算检查井溢流和地表洪水淹没分布的可靠性。主管干和临江管线的检查井在河道水位暴涨初期均有较大的溢流风险。河道洪水的最高水位、上涨速度和持续时间均对管网溢流产生影响，波形陡峭的洪水可使临江检查井产生喷涌，对街道行人和车辆避险有极大安全隐患。对有此类洪水风险的地区，建议可通过在管道或排水道末端合理安装防潮门进行规避。使用本文提出的模拟方法，可有效地对管网溢流进行预测和安全风险评价，对区域洪水安全管理与应急疏散路线的合理制定可提供较强的技术支持，具有一定的工程借鉴和科学意义。

目前的模型当中，尚未考虑地面和检查井井壁表面张力对水位下降期间水流流动的抑制作用，对小尺寸模型试验验证结果产生了不利影响，因此在未来会在模型中加以改进。此外，管道的长期使用与缺少维护会导致出现淤积、堵塞等情况，影响其水力性能，尤其是在管道下游或管道坡度平缓处，由于流速降低，泥沙或污染物固体大颗粒下沉，旱季时附着固结在管壁，对汛期或雨季及时排水不利，在本研究中未对此加以考虑。在后期研究中，还会针对工程中可能出现的防潮门设置情况进行进一步探讨，提出更符合实际工程背景和需求的有效建议。