大家好，欢迎来到IT知识分享网。

摘 要：

【目的】城市下垫面透水面与不透水面之间错综复杂的空间关系使得城市地表降水径流关系变得更加复杂，与此同时，现阶段城市产汇流数值模型的产流计算模块仅停留在不透水地面产流量和透水地面产流量简单叠加的阶段。因此，有必要深入开展城市复杂下垫面径流机理研究，以便能够为城市内涝防治与治理提供理论基础和科技支撑。【方法】通过耦合一维垂向土壤水分运动方程和二维地表浅水方程，构建了适用于城市复杂下垫面的产汇流数值模型，结合室外水文观测试验对数值模型参数进行率定与验证。在此基础上，以汇水单元为计算尺度，就不同汇流模式下城市不透水面积比对城市地表产流系数的影响进行了定量描述。【结果】采用多场场地降雨径流监测结果对所建模型参数进行率定，且模型模拟精度验证结果表明，不同降雨径流过程的纳什效率系数分布在0.71～0.88之间，模型模拟效果较好。【结论】研究表明，在不同汇流模式下，不透水面积占比与径流系数呈线性函数关系，“有效不透水型”汇水单元的场次径流系数要略高于“无效不透水型”汇水单元的场次径流系数。

关键词：

城市洪涝;数值模型;Richard方程;浅水方程;

作者简介：

任梅芳(1987—),女，工程师，博士，主要从事气候变化、城市水文模拟、城市洪涝成因分析、海绵城市建设碳减排等研究。

*宋利祥(1985—),男，高级工程师，博士，主要从事水动力学模型研发等研究。

基金：

国家重点研发计划“城市内涝演化规律与多尺度内涝精准模拟研究”(2021YFC3001401);

引用：

任梅芳, 宋利祥, 庞博, 等. 基于物理机制的城市汇水单元降雨径流特性研究[J]. 水利水电技术(中英文), 2023, 54(9): 37- 47.

REN Meifang, SONG Lixiang, PANG Bo, et al. Physical mechanism-based study on rainfall and runoff characteristics of urban catchment unit[J]. Water Resources and Hydropower Engineering, 2023, 54(9): 37- 47.

0 引 言

随着全球气候变暖和城市化飞速发展，我国城市洪涝事件频繁发生，例如2021年7月20日，郑州市普降大雨，遭遇罕见特大暴雨，发生了严重的城市洪涝事件，造成了严重的人员伤亡和经济损失。

面对逐渐加剧的城市暴雨洪涝灾害，国内外开展了大量的城市洪涝模拟模型和预警系统的研发、应急响应关键技术、城市洪涝承载体脆弱性评价及城市排水管理机制等相关领域的研究。然而，揭示城市产汇流物理机制及其调控机理是城市雨洪灾害成因分析的前提条件；对于城市流域降水径流形成机制方面的研究，尽管国内外开展了诸多的室内外小区及流域试验、数值模型模拟试验等来探讨不同类型城市下垫面的产汇流分布规律，但由于城市区域的下垫面包括自然地表和人工地表，地表覆盖物种类复杂，不同地表的产汇流机制差异性较大，城市地表不透水面与透水面空间分布错综复杂，尤其是现阶段我国对城市复杂下垫面上产流规律的认识还存在一定的缺陷。因此，前人学者研究指出，深入开展城市复杂下垫面径流机理研究，在科学和实践方面都具有十分重要的意义。大量的研究表明降雨强度、土壤含水量、土地利用形式以及不透水地表的面积比例是影响城市降雨径流的主要因素,例如大量研究采用SWMM等模型模拟，说明城市化的发展大幅地增加了城市区域径流系数，与此同时，城市汇流路径的连通性也是影响城市降雨径流的重要影响要素，例如宋晓猛等研究表明不透水面的空间位置对城市流域的流量有不容忽视的作用，石树兰等采用SWMM模型研究表明，考虑不透水面有效性将有效提高雨洪模型的精度，改变汇流路径、减少有效不透水面是减轻雨洪灾害的重要措施。

城市雨洪数值模拟模型是研究城市水安全问题的基本工具，自上世纪70年代起，城市雨洪模型得到了迅速的发展。早期的研究成果中，多为水文学城市雨洪模型，其中具有代表性的有：美国环保署(Environmental Protection Agency, EPA)推出的雨水管理模型(Storm Water Management Model, SWMM),美国伊利诺大学研制的伊利诺伊城市排水区域模拟模型(ILLinois Urban Drainage Area Simulation, ILLUDAS),美国陆军工程兵团研发的暴雨径流模型(Storage Treament Overflow Runoff Model, STORM)等;20世纪90年代以后，国内学者陆续开始进行城市雨洪模型研究。目前，大多数雨洪模型的汇流过程采用圣维南方程组来计算，采用圣维南方程组求解的城市汇流计算模型能够反映水流运动的物理过程，可以描述较为详尽的地表汇流过程。然而，大多数雨洪模型在产流过程的计算中，对不透水地表均采用相同的处理方法，即产流量为降雨量扣除填洼、截留、蒸发等损失量；而对透水部分采用简易的计算方法，例如径流系数法、Horton公式法、SCS法以及Green-Ampt法等；总产流量常采用将透水面和不透水地表的产流量简单叠加形式。上述产流经验公式不仅不能系统准确地反映透水地表的产流物理过程，经验下渗公式的选择和参数确定也多采用试算或者率定的方法得到，这对流域雨洪观测资料的依赖性较强。据此，与城市雨洪模型中的汇流计算过程相比较，现阶段产流数值模型的计算存在缺乏一定的物理机制、计算精度相对较低的问题。因此，研究城市降水径流过程的物理机制，构建物理性较强的产流计算数值模型，对准确、系统地认识城市洪涝致灾机理，提高产汇流数值模拟模型精度都具有十分重要的现实意义。

本文耦合一维垂向非饱和土壤水分运动方程和二维地表浅水方程，构建物理性较强的数值计算模型，形成适用于城市复杂下垫面并基于物理机制的产汇流数值模型，在此基础上，本文以汇水单元为研究尺度，就不同汇流路径下城市化率对城市地表产流系数的影响进行了较为详细的探讨，并进行了定量的描述，以期为城市复杂下垫面径流机理研究奠定理论基础，对准确、系统地认识城市洪涝致灾机理、乃至为我国城市水安全保障研究提供理论基础和科技支撑。

1 数值模型构建

1.1 一维垂向土壤水分运动方程

1.1.1 非饱和土壤水分运动控制方程

非饱和土壤水分运动达西定律

土壤水分运动连续方程为

一维垂向土壤水流运动方程可表示为

式中，t为时间(s);θa为均匀分布的初始含水量(cm3/cm3),θb为地表因湿润条件而维持不变的含水量(cm3/cm3);K(θ)为非饱和土壤导水率(cm/h);Z为土层厚度(m),向下为正。

由于一维垂向土壤水分运动方程为二阶非线性偏微分方程，求解该方程需要不同土壤类型的土壤水分特征曲线，通常情况下，采用直接测定法和经验公式法来得到土壤水分特征曲线，其中，经验公式法是通过测定土壤体积含水率与基质吸力的关系来得到土壤水分特征曲线，并采用某些数学函数拟合这一关系,经验公式有Gardner指数模型;Brooks-Corey模型;Van-Genuchten模型;Kosugi模型,大量前人研究表明Van-Genuchten模型的拟合效果较好，因此，本研究模型中采用Van-Genuchten模型对一维垂向土壤水分运动方程进行求解计算，Van-Genuchten模型可以表示为

式中，θ为土壤体积含水量(cm3/cm3);θs为饱和含水量(cm3/cm3);θr为残余含水量(cm3/cm3);Ks为饱和土壤导水率(cm/h);h为土壤吸力(cm);a、m、n为拟合参数，m=1-1/n。

1.1.2 土壤下渗强度计算模型

土壤累计入渗量的表达式为

式中，I(t)为累计入渗量(cm);θ(z,0)为初始含水率分布(cm3/cm3);L为土层厚度(cm)。

t时刻的表层土壤入渗强度为(不考虑地下水交换)

式中，i(t)为表层土壤入渗强度(cm/s)。

根据Neumann边界条件，可得

式中，q为边界通量(cm/s)。

对于向下为正的垂向坐标，可得表层土壤入渗强度为

式中，i0为表层土壤入渗强度(cm/s);下标1、2分别代表第1层、第2层土壤的值。

同理，底层土壤下渗强度为

式中，iN为底层土壤下渗强度(cm/s);下标N-1、N分别代表第N-1层、第N层土壤的值；N为土壤的分层总数。

1.2 二维地表浅水方程

采用守恒形式的二维浅水方程

式中，U为守恒向量；Eadv、Gadv分别为x和y方向的对流通量向量；S为源项向量。公式为

式中，h为水深(m);u、ν分别为x、y方向流速(m/s);b为底高程(m);r为降雨强度(m/s);i为入渗强度(m/s);g为重力加速度(m/s2);

为摩阻斜率；S0x=-∂b(x,y)/∂x、S0y=-∂b(x,y)/∂y为底坡斜率；n为曼宁糙率系数。

采用二维有限体积法进行地表浅水方程求解，具体算法见文献[18,19]。

1.3 模型耦合

本研究采用松散耦合的方式，实现地表与土壤水动力过程的耦合计算，模型中地表与土壤水动力模型通过下渗进行耦合，地表模型的降雨强度和积水深度为土壤水动力模型提供边界条件，土壤水动力模型求解得到下渗强度，下渗强度作为地表水动力模型圣维南方程的“源汇项”,参与方程求解，进而实现两个模型之间的耦合，具体耦合计算流程如图1所示。

图1 地表-土壤水分运动过程耦合计算流程

2 模型率定与验证

2.1 场地观测试验布置与试验方案

2.1.1 实验场地设计与布置

为了率定本模型参数并验证模型在城市地表的适用性及稳定性，本研究开展场地降水径流模拟观测试验。室外水文观测试验在北京泰宁科创雨水利用技术股份有限公司园区内进行，园区位于北京市昌平区。整个实验场地呈矩形，场地下垫面主要组成部分为不透水屋面和透水地面(自然绿地),场地内透水面地表土壤类型属于沙壤土。

在城市雨洪模拟模型中，汇水单元的地表径流进入排水口的汇流方式被分为三种类型，通常将与排水系统直接相连接的不透水面(Directly Connected Impervious Area, DCIA)称为有效不透水面，与排水系统非直接连接的不透水面(Unconnected impervious area, UIA)称为无效不透水面,相应地，在试验设置中将透水面到不透水(不透水面直接连接到排水系统中)的汇流方式定义为“有效不透水型”,将不透水面到透水面(不透水面间接连接到排水系统中)的汇流方式定义为“无效不透水型”,而将透水面和不透水面并列分布的模式称为“并联型”。三类典型汇水单元的拓扑示意如图2所示;对于“并联型”和“有效不透水型”排水模式，可认为不透水面完全产流，其汇水单元产流量可对透水面产流量和不透水面产流量进行错峰叠加；而“无效不透水型”排水模式中，不透水面的产流量首先汇集到透水面区域，通过透水面之后再汇集到雨水口，该类排水模式中，不透水面所产生的径流在经过透水面时伴随着填洼、下渗等过程，其产流过程不仅和不透水面本身的产流量有关，还和与它相连接透水面的特性有关。因此，考虑到 “无效不透水型”排水模式产流过程的复杂性，同时为了深入研究 “无效不透水型”排水模式中不透水面的产流特性，本研究将降雨径流观测试验场地设计为典型的“无效不透水型”汇水单元模式。

图2 典型汇水单元汇水路径示意

场地试验装置主要由供水装置、降水模拟装置、土壤含水量及土壤水分下渗监测装置、地表径流收集、测量装置，以及排水装置组成。试验布置如图3所示。

图3 实验场地布置

2.1.2 实验场地设计与布置

为了率定一维垂向土壤水分运动方程的相关参数及验证产汇流模型的模拟效果，本研究将降雨强度、土壤前期含水量、雨峰系数及不透水面积比作为变量，开展不同方案的室外场地观测试验。

降雨强度：根据北京市降雨特性，将观测试验中的降雨强度设置为10 mm/h、25 mm/h、50 mm/h和100 mm/h, 降雨历时控制在1 h, 表层土壤前期含水量控制在0.35 m3/m3左右，开展不同降雨强度条件下的降水径流观测试验。

土壤前期含水量：根据对不同土层厚度的土壤体积含水量长时间监测结果(见图4),场次降雨对距离地表10 cm土壤的体积含水量影响较大，其次为距离地表30 cm的土层；对距离地表50 cm和100 cm的土层的体积含水量影响不是很明显；据此，本研究只对距离地表10 cm和30 cm土层的前期含水量进行调节控制。观测试验将降雨强度固定为50 mm/h, 降雨历时固定在1 h, 通过前期预处理方式，将距离地表10 cm处土壤的前期含水量调整为0.36 m3/m3左右，距离地表30 cm处土壤的前期含水量调整为0.39 m3/m3左右；再将距离地表10 cm处土壤的前期含水量调整为0.45 m3/m3左右，将距离地表30 cm处土壤的前期含水量调整为0.42 m3/m3左右进行观测试验，据此，开展两组不同土壤前期含水量条件下的降雨径流水文观测试验。

图4 不同土层土壤含水量变化观测值示意

不同雨峰系数：本研究采用芝加哥雨型，试验采用5 a一遇设计降雨过程，将雨峰系数设置为0.3、0.4和0.5,降雨历时控制在1 h, 表层土壤前期含水量控制在0.35 m3/m3左右，开展不同雨峰系数条件下降水径流观测试验。

不透水面积比：将降雨强度固定为50 mm/h, 降雨历时控制在1 h, 表层土壤前期含水量控制在0.35 m3/m3左右，根据试验场地特性，将不透水面积比固定为15.9%、27.5%、36.6%、42.9%、85%和全不透水地表(100%),开展不同不透水面积比条件下降水径流观测试验。

2.2 模型验证与结果

为了率定模型中土壤水分运行方程的参数取值并验证本次数值模型的计算精度，构建实验场地数值模型(见图5,图示中不透水面积比为42.9%),并采用纳西效率系数(Nash),计算观测与模拟径流值之间的Nash效率系数，Nash效率系数越接近1,说明模拟效果越接近观测值，Nash效率系数计算公式为

式中，Nash为纳西效率系数；Qobs为试验观测径流值；Qsim为数值试验模拟径流值；Qobs为观测径流平均值。

图5 实验场地数值高程及糙率插场示意

表1给出用于模型率定的3场实测降雨径流过程的各项参数，Ⅰ号降雨径流过程降雨强度为50 mm/h, 降雨历时为1 h, 场地内不透水地表占比为42.09%;Ⅱ号降雨径流过程降雨强度为50 mm, 降雨历时为1 h, 场地为全透水面地表；Ⅲ号降雨径流过程降雨强度为100 mm/h, 降雨历时为1 h, 场地为全透水面地表。观测试验场地透水面土质为砂壤土，经率定，模型中透水面入渗率取值范围为1.04～7.20 cm/h, 根据不同土质入渗率经验值，本次率定结果透水面土质入渗率介于壤土和砂壤土之间，根据场地土壤入渗率实测结果，可认为入渗率取值与实际基本吻合。

表2给出用于验证的3场实测降雨径流过程的各项参数，Ⅳ号降雨径流过程为5 a一遇设计降雨，降雨历时为1 h, 雨峰系数为0.3,场地为全透水面地表；Ⅴ号降雨径流过程为5 a一遇设计降雨，降雨历时为1 h, 雨峰系数为0.4,场地为全透水面地表；Ⅵ号降雨径流过程为5 a一遇设计降雨，降雨历时为1 h, 雨峰系数为0.5,场地为全透水面地表。表2给出了场次降雨径流观测过程和模拟过程的洪峰值及降雨径流过程的Nash效率系数(见表2),经计算，不同降雨径流过程的Nash效率系数分布在0.71～0.88之间，模型模拟效果较好。

不同观测与模拟降雨径流过程及其对比如图6所示，研究表明，模型对峰值的计算能力较好，但对达到峰值之前径流过程的捕捉能力稍差，数值模型对峰值之前的径流过程拟合偏差的主要原因可能有以下几个方面：(1)在场地观测试验过程中，降雨模拟系统从开启到稳定需要一段时间，这一段时间内降雨设备并没有达到预设的稳定降雨值，因此，导致数值模拟结果与观测结果之间存在一定的偏差；(2)场地降雨径流过程地表水层较薄，与明渠均匀流水力学计算存在一定的差异，这也是导致观测过程与模拟过程存在一定差异的原因。

图6 观测试验与数值模拟试验径流过程对比

3 不透水面面积对径流系数的影响探讨

为了能够充分体现出城市下垫面透水面与不透水面空间分布对城市地表径流系数的影响，本研究分别对“无效不透水型”和“有效不透水型”汇流路径条件下，不透水地表面积占比与径流系数之间的响应关系进行模拟和初步探讨。

3.1 无效不透水型汇流路径条件

为了能够探讨“无效不透水型”汇流模式下，不透水面面积占比对径流系数的影响，数值模型试验以纵坡为1‰,边长为100 m×100 m的矩形区域作为假想汇水单元，分别以1 a一遇、5 a一遇、10 a一遇和20 a一遇设计降雨(芝加哥雨型，降雨历时为1 h, 雨峰系数为0.4)作为降雨条件，将透水面导水率固定为3.6 cm/h, 土壤前期含水量设置为0.35 cm3/cm3,对汇水单元内不透水面积的比例进行调整，从而探讨不透水面积占比与径流系数的响应关系。

图7(a)给出不同不透水面积占比条件下场次降雨量与径流量的关系曲线，图7(b)给出不同不透水面积占比条件下径流系数的分布情况及不透水面积占比与径流系数的响应关系曲线；取径流系数分布的中位数，拟合不透水面积占比与径流系数的关系曲线，经拟合计算，不透水面积占比与径流系数呈线性函数关系。

图7 不同不透水面积占比条件下降雨-径流关系及不透水面积比与径流系数关系(无效不透水型)

3.2 有效不透水型汇流路径条件下

为了能够探讨“有效不透水型”汇流模式下，不透水面面积占比对径流系数的影响，数值模拟试验同样以纵坡为1‰、边长为100 m×100 m的矩形区域作为假想汇水单元，分别以1 a一遇、5 a一遇、10 a一遇和20 a一遇设计降雨(芝加哥雨型，降雨历时为1 h, 雨峰系数为0.4)作为降雨条件，将透水面导水率固定为3.6 cm/h, 只对汇水单元内不透水面积的比例进行改变，可得到不透水面积占比与径流系数的响应关系。

图8(a)给出在有效不透水型汇水单元上，在不同的不透水面积比条件下，场次降水径流关系曲线及不透水面积占比与径流系数的关系曲线，该曲线与“无效不透水型”汇水单元上得到的关系曲线基本相同；与“无效不透水型”汇水单元计算结果进行对比分析，在同等边界条件下，“有效不透水型”汇水单元的场次径流系数要略高于“无效不透水型”汇水单元的场次径流系数，对不同不透水面积占比条件下场次径流系数组的中位数进行拟合计算结果表明[见图8(b)],在“有效不透水型”汇水单元上，不透水面积占比与径流系数也呈线性函数关系。

图8 不同不透水面积占比条件下降雨-径流关系及不透水面积比与径流系数关系(有效不透水型)

4 结果讨论

目前，城市地表不透水面的空间分布及有效性问题已经受到了城市水文研究的广泛关注，其中班玉龙等借助于SWMM模型，在不同汇流路径下的管网径流进行模拟分析，研究结果表明透水面对不透水面进行阻隔能会产生较为明显的水文效应，最高可使研究区总下渗量增加1倍；石树兰等采用SWMM模型，以北京市凉水河流域的大红门排水片区为研究区，通过设置不同重现期暴雨情景，分析不透水面空间组合方式对洪涝的影响，结果表明考虑不透水面有效性后，模型模拟效果产生了较为显著的提升；SHUSTER等研究表明，由于城市化地区不透水地表空间分布的复杂性，使用总不透水面积进行城市雨洪模拟结果会高于实测径流量。本次研究结果表明城市汇水单元降水径流过程会受到城市地表不透水面的空间分布的影响，这与前人研究结果一致，本研究在前人研究的基础上，采用一维垂向土壤水分运动方程对场地下渗过程进行描述，对实现基于物理机制的产汇流数值模型研究城市降水径流关系进行了初探。

然而，城市降水径流关系不仅受到不透水地表面积占比、透水地表与不透水地表空间分布的影响，同时还会受到汇水单元坡度、汇水单元透水地表土壤前期含水量、透水地表入渗率及汇水单元面积等不同因素的影响，本研究仅从不透水地表的空间分布及面积占比的角度初步探讨了城市化对径流系数的影响，存在一定的缺陷，但对城市降雨径流响应关系的研究奠定了理论基础，未来在此模型的基础上拟进一步开展城市地表不同驱动要素与径流系数响应关系的研究。

5 结 论

城市下垫面高度的空间变异性使得城市降雨径流之间的物理机制成为城市水文科学研究的难点，而相比于城市雨洪模型中的汇流模型，现阶段产流数值模型的计算存在缺乏一定的物理机制、计算精度相对较低等问题。

鉴于一维垂向土壤水运动方程可以较准确地描述不同边界条件下土壤水分的运动过程，本研究采用一维垂向土壤水运动基本方程作为产流模型，耦合二维地表浅水方程，形成地表-土壤水运动过程耦合的城市地表产汇流计算模型，并基于该模型，在汇水单元尺度上，对城市降水径流响应关系进行了初探研究。主要研究结论如下。

(1)采用室外场地水文试验观测结果对所建模型参数进行率定并对模型模拟精度进行验证，经计算，不同降雨径流过程的Nash效率系数分布在0.71～0.88之间，模型模拟效果较好。研究表明，模型对峰值的计算能力较好，但对达到峰值之前径流过程的捕捉能力稍差，数值模型对峰值之前的径流过程拟合偏差的主要原因可能有以下几个方面：第一，在场地观测试验过程中，降雨模拟系统从开启到稳定需要一段时间，此段时间内降雨设备并没有达到预设的稳定降雨值，导致数值模拟结果与观测结果之间存在一定的偏差；第二，场地降雨径流过程地表水层较薄，与明渠均匀流水力学计算存在一定的差异，这也是导致观测过程与模拟过程存在一定差异的原因。

(2)本研究在汇水单元尺度上，就不同汇流模式下城市化率对城市地表产流系数的影响进行了定量描述，研究表明，在“有效不透水型”和“无效不透水型”两种汇流模式下，不透水面积占比与径流系数呈线性函数关系，同时，“有效不透水型”汇水单元的场次径流系数要略高于“无效不透水型”汇水单元的场次径流系数。

[18] SONG L X,ZHOU J Z,LI Q Q,et al.An unstructured finite volume model for dam-break floods with wet/dry fronts over complex topography[J].International Journal for Numerical Methods in Fluids,2011,67(8):960-980.

[19]宋利祥，周建中，王光谦，等.溃坝水流数值计算的非结构有限体积模型[J].水科学进展，2011,22(3),373-381.SONG Lixiang,ZHOU Jianzhong,WANG Guangqian,et al.Unstructured finite volume model for numerical simulation of dam-break flow[J].Advances in Water Science,2011,22(3):373-381.

水利水电技术（中英文）

水利部《水利水电技术（中英文）》杂志是中国水利水电行业的综合性技术期刊（月刊），为全国中文核心期刊，面向国内外公开发行。本刊以介绍我国水资源的开发、利用、治理、配置、节约和保护，以及水利水电工程的勘测、设计、施工、运行管理和科学研究等方面的技术经验为主，同时也报道国外的先进技术。期刊主要栏目有:水文水资源、水工建筑、工程施工、工程基础、水力学、机电技术、泥沙研究、水环境与水生态、运行管理、试验研究、工程地质、金属结构、水利经济、水利规划、防汛抗旱、建设管理、新能源、城市水利、农村水利、水土保持、水库移民、水利现代化、国际水利等。