王浅宁，邹德昊，王清正，彭勇，吴剑
(1.大连理工大学建设工程学部，辽宁 大连 116024；2.珠江水利委员会珠江水利科学研究院，广州 510610)
近年来，随着全球气候的极端变化，我国极端暴雨多发频发[1]，城市内涝形势严峻。城市暴雨内涝模拟计算作为城市防洪排涝核心技术之一，其模型参数的选取对模拟结果的影响不可忽视。地表径流速度作为坡面水文过程重要参数之一，对流域水文学计算、城市暴雨内涝模拟计算等至关重要。地表径流速度即地表坡面水流汇集运动过程时水流的流速[2]。在流域产汇流方面，地表径流速度控制降雨向径流和下渗的分配[3]，直接影响产流量，因此为计算符合实际的径流系数，需要准确估算地表径流速度[4-5]。在水文水动力学原理方面，地表径流速度常用于推导雷诺数、曼宁系数等水文参数[6-7]。在城市暴雨内涝模拟计算方面，地表径流速度的准确估计对于计算集水区径流过程至关重要[8]。Li等[9]在海绵城市洪涝及水质模拟中明确指出地表径流速度对城市洪水总径流量、洪峰流量和水质均有较大影响。
目前，国内外学者对地表径流速度的研究大多围绕山区流域或高原土壤展开[10-12]，侧重于对产沙过程及土壤侵蚀影响进行分析[13-15]，鲜有涉及地表径流速度对城市内涝灾害模拟的影响。地表径流速度作为MIKE-URBAN模型中集水区产汇流模块的重要参数之一，可通过修改URBAN模型中对应参数来设置地表径流速度取值。因此，本文设置系列地表径流速度方案，模拟研究区在不同频率暴雨条件下的内涝积水过程，并量化分析地表径流速度对研究区内涝灾害的影响规律，为城市暴雨内涝模拟分析以及洪涝灾害预防、治理提供支撑。
研究区位于大连市沙河口区北部，见图1，面积为0.91 km2。该区域城市化程度较高，排水体系相对独立，下垫面属性复杂，产汇流特征不明显，地形起伏较大，是沙河口区内涝积水多发区域。
图1 研究区
1.2.1设计降雨
根据住建部《试点城市内涝积水点分布图绘制说明》，对1、2、5 a一遇的场景计算时推荐采用短历时(2 h或3 h)设计降雨雨型，对50 a一遇及以上的场景计算时推荐采用24 h长历时设计降雨雨型。本文设置1、2、3、5、10、20 a一遇短历时2 h以及50、100、200、500 a一遇长历时24 h共10种设计降雨作为模型输入。Qin等[16]的研究指出，芝加哥雨型在中国城市排水设计中有较好的应用效果。根据历史暴雨统计资料，大连市中心城区暴雨历时短且以单峰型为主。因此，本文选取芝加哥雨型作为短历时2 h设计降雨雨型，雨峰系数(峰现时间与暴雨历时的比值)r=0.5。根据《大连市中心城区暴雨强度公式及查算图表》，采用式(1)计算研究区暴雨强度。
(1)
式中：q为降雨强度，mm/min；A=1 230.157；C=0.724；P为重现期，a；t为降雨历时，min；b=5.783；n=0.661。图2为大连市中心城区1、2、3、5、10、20 a一遇短历时2 h降雨过程线，时段长为5 min。其中，不同重现期2 h总降雨量为36.5～70.8 mm。
图2 短历时2 h设计降雨过程线
图3 长历时设计降雨过程线
表1 长历时24 h设计降雨特征值
1.2.2地形、土地利用及管网现状
研究区内地形高程空间分布见图4。基于Landsat TM遥感影像反演得到研究区域内建筑、道路及人造表面等3种土地利用类型。研究区地下排水管网排水设计标准一般为1至5 a一遇，总长9.76 km，管网空间布置见图5。
图4 研究区DEM高程
图5 研究区土地利用及管网分布
采用丹麦水利研究所(Danish hydraulic institute,DHI)开发的水文水动力数值模型，模型具有较高的模拟精度和计算速度，已被广泛应用于城市内涝灾害模拟并得到了较好的验证[17-20]。研究区主要包含地下排水管网和地表地形，构建一维排水管网模型和二维地表模型的耦合数值模型，对研究区暴雨积水内涝过程进行模拟。
2.1.1一维管网模型
研究区一维管网模型构建主要包含产流模块、汇流模块以及排水管网水动力模块。首先，产流模块采用泰森多边形法将研究区划分成与人孔(Manhole)同等数量的子集水区(Catchment)，计算子集水区综合不渗透系数，降雨经产流模块转化成地表径流量；接着，汇流模块将子集水区产生的径流汇入与之相连的人孔，水流汇集时间受地表径流速度参数的影响，参数越大汇集时间越短，从而影响了人孔处水流的汇集过程。研究区各子集水区水流汇集时间计算公式为
Ti=Li/V
(2)
最后，排水管网水动力模块运用隐式有限差分法求解圣维南方程组，并用自适应步长法推进模拟时间，高效准确的计算多连通的分支环管网络的水动力过程。经概化，研究区一维地下排水管网模型共包含181个人孔、4个排放口(Outlet)、193条排水管道(Link)以及181个子集水区，其中，子集水区最大面积为2.29万m2，最小面积为116 m2。
2.1.2二维地表模型
研究区二维地表模型构建主要包含模型边界、剖分网格、特殊地物处理等内容。二维地表模型边界是人孔处的雨水溢流。模型采用有限体积法对原始方程进行空间离散，并将连续介质细分为不重叠的单元，使空间域离散化。二维地形模型采用矩形结构化网格进行剖分划分为4 032个矩形网格，网格大小15 m×15 m。同时，为真实反映研究区域实际地形，对房屋、道路等特殊地物进行处理，将建筑物处地形高程抬高0.7 m，道路处地形高程降低0.15 m。
构建的模型参数分为确定性参数和不确定性参数。模型中管道长度、管径、管底埋深以及管道曼宁系数等确定性参数根据管网实测资料确定。其中，曼宁系数是影响排水管道中水流水动力特征的重要因素，曼宁系数的改变势必会引起管道内水流状态的变化。郭永鑫等[21]指出混凝土排水管曼宁系数为0.011～0.015，考虑到研究区管道磨损情况，模型排水管道曼宁系数取0.012。模型中不确定性参数则需根据历史实测水文数据率定得到，包括地表糙率、不透水系数、折减系数等。根据文献[22-26]确定模型的不确定性参数，具体见表2。
表2 模型参数取值
表3 积水点验证信息汇总
模型中地面径流速度取值直接影响区域内涝模拟结果，其取值目前仍没有统一的经验公式和方法。其中：吴淑芳等[27]指出在不同径流调控措施下坡面径流流速大致范围为0.03～0.23 m/s；杨大明等[28]模拟在自然降雨条件下不同坡度对应的坡面流运动情况，得出流速分布为0.08～0.40 m/s。因此，本文设置地表径流速度为0.02、0.05、0.10和0.30 m/s共4种场景。按穷举法设计不同重现期降雨、不同流速的多种组合工况，最终共确定40个模拟方案。
以管道充满度评估管道承载状态，该指标表示管道内水流的充满程度，以管道水深高程与管径之比表示[29]。在不同重现期设计降雨条件下，管道充满度大于1时表示管道内水位超过管道顶端高程，管道已经变成有压管道，不再是重力流排水，严重影响了周围相邻管道雨水的正常排放，很容易导致管道漫溢，不满足当前降雨的排水需求。管道充满度小于1时表示管道可满足当前降雨的排水需求，故选取F<1作为管道是否达标的判定标准。具体公式为
(3)
式中：F为管道充满度；Wlevel为水位高程，m；Pinertlevel为管道底高程，m；Phight为管道高度，m。
采用人孔溢流量表征研究区域的内涝积水情况。为进一步评估研究区地表径流速度对研究区内涝风险的影响，需对研究区内涝淹没进行风险等级划分。在借鉴赵燕霞[30]、戴晶晶等[31]以及侯精明等[32]对城区的内涝风险等级划分方法的基础上，综合考虑积水深度和积水历时两方面因素，并按如下标准划分研究区内涝风险等级：低风险区为路面积水深度>0.03～≤0.25 m且积水历时大于30 min；中风险区为路面积水深度>0.25～≤0.50 m且积水历时大于60 min；高风险区为路面积水深度超过0.50 m且积水历时大于60 min。
基于6种短历时2 h设计降雨条件模型模拟结果，统计在不同重现期设计降雨和地表径流速度条件下达标管道长度及其与管道总长的比值，结果见表4。
表4 各重现期降雨条件下管道达标统计
由表4可知，地表径流速度越低，管道排水负荷越小，管道达标率越高，同时随着降雨重现期变短，地表径流流速对管道排水负荷的削减效果变得明显。比如：在20 a一遇设计降雨条件下，地表径流速度0.02、0.05、0.10 m/s比地表径流速度0.30 m/s的管道达标百分率提高7.6%、6.4%、2.1%；在5 a一遇设计降雨条件下对应提高13.0%、5.7%、2.1%；在1 a一遇设计降雨条件下对应提高36.2%、6.3%、2.6%。降雨沿集水区坡面汇流至雨水井口后流入排水管道，当地表径流速度变小时，集水区汇流时间相应变长，进而导致管道排水过程延长，排水负荷变小。当下渗能力相同时，设计降雨重现期越短，集水区产流量越小，进入排水管道的水量越少，故降低地表径流速度可明显缓解管道排水负荷，同时设计降雨重现期越短，地表径流流速对管道排水负荷的削减效果越明显。
模型中人孔的溢流量决定了内涝积水量，对模型中各人孔溢流量叠加得到研究区总溢流量，其随时间变化过程见图7。在不同重现期设计降雨条件下溢流量峰值及峰现迟滞时间见表5。由表5和图7可知，地表径流速度越低，检查井溢流量峰值越小，峰现时间滞后越长，随着降雨重现期增长，溢流量的削减效果减弱。比如：在20 a一遇设计降雨条件下，当地表径流速度为0.02、0.05、0.10 m/s时，比地表径流速度0.30 m/s的溢流量峰值分别降低了2 750、893、312 m3，峰现时间分别滞后56、14、5 min；在5 a一遇设计降雨条件下，地表径流速度为0.02、0.05、0.10 m/s时，比地表径流速度0.30 m/s的溢流总量分别减少了1 776、538、207 m3，峰现时间分别提前53、11、4 min。地表径流速度降低导致汇流时间增长，当下渗能力相同时，下渗水量增多，进而对检查井溢流量起削减和延缓的作用。因此，在应对城市内涝积水时，可通过改变下垫面条件等方式降低地表径流速度，达到对洪量的削减、滞留，进而减少内涝积水造成的损失。
图7 不同重现期设计降雨溢流量变化过程
采用上述内涝风险等级划分方法，得到研究区各等级风险区淹没面积随地表径流速度及设计降雨重现期变化，见图8，可见地表径流速度对研究区内涝风险影响较大。表6为在不同地表径流速度下各等级风险区淹没面积。
图8 研究区各等级风险区淹没面积汇总
表6 研究区各等级风险区评估指标
由表6可知，研究区高风险区和低风险区淹没面积变化趋势规律一致，在相同重现期设计降雨条件下，研究区低、高风险区淹没面积随地表径流速度降低明显减小，在地表径流速度为0.02 m/s时淹没面积最小，且设计降雨重现期越长淹没面积差距越明显。在50 a一遇设计降雨条件下，地表径流速度0.02 m/s比地表径流速度0.30 m/s的高风险区淹没面积降低了0.9万m2，低风险区淹没面积降低了7.04万m2；在500 a一遇设计降雨条件下，地表径流速度0.02 m/s比地表径流速度0.30 m/s的高风险区淹没面积降低了1.64万m2，低风险区淹没面积降低了8.37万m2；研究区中风险区淹没面积随地表径流速度及设计降雨重现期变化反复。在50 a、100 a和500 a一遇设计降雨条件下，地表径流速度0.02 m/s的中风险区淹没面积最小，分别为4.03和4.14万m2；在200 a一遇设计降雨条件下，地表径流速度0.05 m/s的中风险区淹没面积最小，为4.37万m2。为进一步分析中风险区淹没面积变化反复原因，给出研究区不同地表径流速度场景下内涝风险区分布，见图9。
图9 研究区不同场景内涝风险区分布
由图9可知，中、高风险区淹没面积主要位于研究区地势低洼区域。在200 a一遇设计降雨条件下，随着地表径流速度由0.30 m/s降低至0.02 m/s，高风险区淹没面积有1.35万m2降为中风险区，中风险区淹没面积有0.97万m2降为低风险区，低风险区淹没面积有9.32万m2降为无风险区。结合内涝积水结果可知，地表径流速度变小或设计降雨重现期变短均会减少研究区的内涝积水，从而使部分高等级风险区转为低等级风险区，且低风险区转为无风险区易于中风险区转为低风险区。因此，当地表径流流速降低、设计降雨重现期变短时，低、高风险区淹没面积大幅降低，中风险区淹没面积变化反复。
通过设置系列地表径流速度模拟方案，分析研究区不同重现期设计降雨条件下的内涝积水过程，探讨了地表径流速度对研究区管网排水能力、内涝积水以及内涝风险区淹没面积等的影响规律，主要结论如下：
在短历时2 h设计降雨条件下，地表径流速度越低、降雨重现期越短，降雨后管道的排水负荷越小，管道达标率越高，地表径流流速对管道排水负荷的削减效果越明显。随着地表径流速度降低，达标管道长度和管道达标率最大可分别提高4.7 km和48.05%。
地表径流流速对内涝积水影响显著，地表径流速度越低，降雨产生的内涝积水峰值越小，峰现时间滞后越长。随着降雨重现期增长，内涝积水的削减效果减弱;随着地表径流速度降低，溢流总量峰值最大可降低2 750 m3，峰现时间最长可滞后56 min。
在长历时24 h设计降雨条件下，研究区低、高风险区面积随地表径流速度降低明显减小，且设计降雨重现期越长风险区面积削减量越大。随着地表径流速度降低，低、高风险区面积最高可分别降低1.64万、8.37万m2，但中风险区面积变化反复。
本文设置的地表径流速度场景虽然是针对本文研究区域的特殊地形和土地利用状况的，但对其他影响内涝模拟结果的因素的研究有较好的参考价值，可为城市规划、建设提供基础支撑依据。此外，本文采用穷举法的方式设计了不同重现期降雨、不同流速的多种组合工况，避免了降雨重现期与流速相关性的考量。同时，本文也未考虑坡度对地表水流汇集速度大小的影响，仅探究了在均匀流速下改变地表径流速度大小对研究区内涝模拟的影响，以上不足将是作者今后研究当中的重要工作之一。