SWMM产汇流模型参数的设计条件等价优化

刘子龙，周玉文，王 强，叶婉璐

(1.北京工业大学 建筑工程学院，100124北京;2.北京城市规划设计研究院，100037北京)

中国雨水系统设计长期以来均是采用推理公式法结合恒定均匀流理论进行管网设计计算，该方法由于其自身适用条件的限制，仅适用于较小区域的设计，应用于较大汇水流域时计算精度偏低［1］.因此，新的室外排水设计规范规定当汇水面积超过2 km2时，宜采用数学模型法校核并调整其雨水设计方案［2］.!影响模型设计校核是否有效的关键因素在于整个建模过程中，模型相关参数是否等价地反映了初始设计条件，即“设计条件等价”［3］.模型中，子汇水区产汇流计算和管网汇流计算是相互分开的，因此，“设计条件等价”包括了子汇水区参数等价和管网参数等价.管网参数的等价均较易实现，如管径、管长、粗糙系数和其他水力构筑物的结构参数等，是设计阶段和模型校核阶段一一对应的确定性内容，而子汇水区参数包括汇水时间和洪峰径流系数的等价较难实现.

时间面积法是最简单的实现子汇水区参数等价的方法，因为其参数是一致的，但模拟超标降雨情景时，固定径流系数将造成较大误差.暴雨洪水管理模型(storm water management model，SWMM)作为一种更先进的模拟技术，已在国内外广泛应用，其地表计算模块分别采用前损后损法和非线性水库模型进行产流和汇流过程模拟［4］.与时间面积法不同，SWMM模型应用于设计方案校核时需要大量基础数据支持，而国内对其产汇流计算参数的研究较少，部分成果也仅仅是提供了模型参数一个合理化取值的范围［5-6］，在该范围内参数的不同取值对最终模拟结果仍具有较大影响，即参数并不具有良好的稳健性.过去往往由于采用了不合适的概化参数，使设计方案得到了错误的评价结果.如何在无法获取全部参数的情况下，通过限定参数的合理化取值范围，得到符合设计条件的参数组合是本研究的内容.

1 设计条件等价优化原理

根据推理公式法原理，考虑流域地表前期湿润条件下，认为其地表产流强度r在时空分布上为恒定常数，在汇水面积a、单位时间dt内形成的出流流量 dQ也为常数，即［7］

式中:r为产流强度;q为暴雨强度(mm/h);μ为损失强度(mm/h);dQ/dt为单位时间的出流流量(m3/s);0.278为单位换算系数;a为汇水面积(km2).

开始阶段，地表径流随汇水面积的不断增加而增加，直至全部面积参与汇流(全部汇流)，此时产生流域洪峰流量Qm.此后，如果降雨强度和损失强度保持恒定，则出流量Q不再继续增加，达到出流稳定阶段.至降雨停止，出流流量随历时逐渐减少，推理公式法的流量过程线假设如图1所示.

图1 推理公式法的流量过程线假设

推理公式法假设设计降雨历时等于全面汇流时间，即设计汇流时间tc时，形成地表洪峰流量Qm，令ψ=1-μ/q，则洪峰流量

式中ψ为设计洪峰径流系数，中国当前采用地表综合径流系数近似代替，暴雨强度q采用某设计重现期P下的暴雨强度公式计算，即

式中:a为与重现期相关系数;b和c为暴雨强度公式常数;tc为汇水区设计汇水时间.

以推理公式法的流量过程线假设为基础，SWMM产汇流模型等价的实现，以恒定设计暴雨为模型输入条件，在合理参数范围内，调整产汇流模型参数使子汇水区出流流量过程线与推理公式法的流量过程线特征相吻合.推理公式法的流量过程线特征包括设计汇水时间，即出流流量达到稳定的时间，以及设计洪峰流量.因此，判断产汇流模型与设计条件是否等价的要素同样包括了径流流量达到稳定的时间及计算洪峰径流量.

由于需模拟流量稳定阶段，SWMM产汇流模型的输入降雨过程线时长需至少大于设计汇水时间，采用2tc为输入降雨时长，同时考虑地表前期湿润条件，不计前损的洼蓄以及后损的初始土壤下渗速率及其衰减过程，仅采用霍顿(Horton)模型的稳定下渗速率［8］.SWMM产汇流模型需调整参数包括汇水区宽度Width(m)、汇水区坡度PcntSlop(%)，不透水面积所占比例 PcntImperv(%)、不透水区地表曼宁系数N-Imperv、透水区的地表曼宁系数N-Perv以及霍顿模型的土壤稳定渗透速率MinRate(mm/h)共6项参数.

通常情况下，模型参数的调整过程可以通过人工和机算两种途径实现.需调整参数较少时，可根据经验进行人工调整，但这通常需要多次试算以及巨大的工作量，本文采用粒子群全局优化算法，研究其应用于模型参数的自动优化调整过程.

2 非线性权值递减粒子群算法及实现

粒子群优化算法是全局集群寻优算法的一种，以“粒子”表征一组待求解的未知解，根据未知解的数量，单个粒子可能具有多维的搜索空间.各个粒子的适应度通过一定的计算方法进行评估，从而影响下一代粒子种群的生成.影响因素包括两个方面，一是某个粒子自身的历史最优组合，二是整个群体的历史最优组合.通过该继承策略，算法反复迭代，直至寻找到最适宜的参数组合.

依据算法原理，粒子群单个粒子的多维搜索空间即为需优化的未知参数数量(本实例共6项参数).设粒子群落中某一粒子为i，其当前位置为xi，同时，各粒子还包括一个多维速度项决定下一代粒子变化的方向及幅度，记为飞翔速度v.为了保证群体的多样性，根据限定条件，即参数的合理化取值范围，随机生成粒子的初代群体.

式中:xmin、xmax分别为参数取值范围的上下限值;r1和r2为0～1间的均布随机数;vmax为飞翔速度的最大限值.从初代开始，群体通过不断搜索适应度最高的解，从而使最适宜的粒子迭代得到更新.在每一次迭代中，各粒子均根据当前其自身的历史最优解pbest以及群体的历史最优解gbest计算各自的飞翔速度，并得到下一代的位置.设第n次迭代中，粒子i的位置及速度分别为xi(n)和vi(n)，则n+1次迭代中，其速度以及位置计算如下:

式中:w为粒子的惯性系数，代表粒子保持其当前速度的能力;c1、c2为置信系数，通常均取为 2.r1、r2均为(0，1)范围内的均布随机数，采用随机数可使算法更易实现全局的概率搜索.根据文献［9］的研究结果，采用非线性的惯性系数权值递减策略，可确保初始阶段不易陷入局部最优，而后期随着权值的递减，可辅助算法得到更高精度的优化解，本文采用其指数递减函数计算w值，即［9］

式中:ws、we分别为初始及终止时刻的惯性系数;n为当前迭代次数;N为算法终止迭代次数;cN为指数因子，取为 10.依据 Shi等［10］试验结果，ws和we分别取为0.95和0.4时，算法能够获得较快的收敛速度且求解精度高.

粒子群优化计算步骤如下:

步骤一，根据求解未知解的个数，确定粒子的维度，初始化生成初代粒子群及各粒子的初始速度，计算初代惯性系数，设置群体规模和最大迭代次数.

步骤二，根据设计汇水时间和重现期，通过暴雨强度公式计算平均降雨强度.根据粒子各个维度取值，构建SWMM模型文件(INP文件)，并调用水力分析计算引擎完成延时模拟计算.模拟完成后，通过调用模型结果文件(OUT文件)获取汇水区径流(runoff)过程线，按目标评价函数计算粒子适应度.

步骤三，比较各粒子当前适应度和其个体历史最优解并对其进行更新.对比每个粒子的个体历史最优解，选择其中最优适应度粒子为群体历史最优解.

步骤四，根据式(6)更新各粒子的飞翔速度，并检验是否超过速度限值.依据式(7)更新粒子位置，同样检验是否超过设计变量的取值范围.

步骤五，满足迭代终止条件则算法终止，返回最后一次迭代的群体历史最优值，否则返回步骤二，根据式(8)计算递减惯性系数，继续算法运行.

单个汇水区的优化程序流程如图2所示.

图2 单个汇水区优化程序流程

3 实例及分析

3.1 实例一

对单个子汇水区进行人工试算调整，验证提出的设计条件等价原理实施的可行性.人工试算调整需首先定义未知参数的合理化取值范围，避免调整结果偏离实际.其次进行参数的灵敏度分析，根据目标明确参数的调整方向及需调整的大小.

以某设计子汇水区为例，其面积为1.45 hm2，设计重现期和设计汇水时间分别为P=2 a以及tc=10 min.按当地暴雨强度公式

计算其设计降雨强度为1.51 mm/min，设计洪峰径流系数为 0.65，按式(2)计算设计洪峰流量为0.23 m3/s.

汇水区模型宽度和坡度可结合实际情况采用合理的取值范围.本文宽度范围取100～1 000 m，地表坡度值范围取0～1%.根据文献［5-6］的研究成果，总结了模型6项参数的合理化取值范围，如表1所示.

表1 模型参数合理化取值范围

模型参数的灵敏度分析采用扰动分析法，即在参数基准值附近给定一个人工干扰值(本文采用0.4、0.7、1.3 和 1.6 的变化倍数)，而其他参数保持不变，计算该参数在范围内波动产生的模型结果变化率.需统计的模型结果包括出流流量达到稳定的时间(即汇流时间)以及洪峰径流量，出流稳定阶段以计算径流量达到洪峰径流量95%为判定标准.各参数的灵敏度分析结果如图3所示.

以参数变化倍数为横坐标，扰动结果的变化率为纵坐标，进行线性回归曲线拟合，各参数的拟合曲线斜率如表2所示.曲线斜率大小即反映各参数的灵敏程度，斜率正负反映了模拟结果随参数的变化方向.以参数基准值作为模型初始值，设计条件等价的人工调试过程如表3所示.

图3 SWMM产汇流参数灵敏度分析

表2 扰动结果变化率的拟合曲线斜率

表3 SWMM模型参数的人工调试过程

通过人工调试的模型计算径流量过程线变化如图4所示.可以看出，通过模型参数灵敏度分析，人为控制参数的变化方向，并在合理范围内进行调整，可使模型计算流量过程线与推理公式法的假设流量过程线特征不断趋于吻合，并最终实现与设计条件的等价.

图4 人工调试的计算径流过程线变化

3.2 实例二

采用本文改进粒子群优化算法计算实例一汇水区的产汇流参数，各项参数的合理化取值范围与实例一相同.考虑计算汇水时间和峰值流量的设计条件等价，采用其相对误差之和为粒子适应度评价指标，构建算法寻优目标函数

式中:TcS为模型计算汇水时间(模拟径流量达到模拟洪峰径流量95%的时间);Tcd为设计汇水时间;QmS为模拟洪峰流量;Qmd为设计洪峰流量.

设置粒子群群体规模20，算法终止条件为迭代50次，各维度粒子的最大飞翔速度为±0.1×(最大限值-最小限值)，非线性权值递减粒子群算法优化过程如图5所示.可以看出，随着迭代次数的增加，粒子群适应度不断趋近于零，模拟汇水时间和洪峰流量存在一定的波动，迭代20次后达到稳定，得到的最终优化参数组合如表4所示.

图5 粒子群算法优化迭代过程

表4 粒子群最终优化参数组合

最终优化计算汇水时间为10 min，相对误差为0%，优化计算洪峰流量为0.237 m3/s，相对误差为3%.

基于流域的实际情况，当部分参数存在更明确的取值范围时，可通过固定值或更小的寻优空间限制其优化值.假定采用宽度限值100～200 m，不透水率30%～40%，再次进行寻优计算.增加限值条件后，最终优化参数组合如表5所示.

表5 增加限值条件的粒子群最终优化参数组合

增加限值条件后的最终优化计算汇水时间为11 min，计算洪峰流量为 0.236 m3/s.可见，通过增加参数的限值可使优化参数组合更符合实际.

对比实例一和实例二的参数调试过程，采用粒子群优化算法，无需设计人员有较多的建模工程经验，即可快速获得适宜的等价模型参数组合，减少设计人员工作量且无需对参数的影响进行大量分析.

4 结语

采用模型技术对雨水管网设计方案校核的前提是必须实现模型参数的“设计条件等价”.本文依据推理公式法的基本假设提出SWMM产汇流模型参数的等价化原理，并通过实例验证了等价原理实施的可行性.采用非线性权值递减粒子群优化算法可快速实现参数的等价优化，在无法获取全部产汇流参数的情况下，辅助模型技术在雨水管网设计方案校核中的应用，具有一定的实际价值.

值得注意的是，虽然当前推理公式法设计参数可通过经验或理论公式计算［1］，对于较小的子汇水区，地表产汇流计算精度较高，但由于设计条件的不确定性，模型计算结果与设计条件之间允许存在一定的误差，采用模型参数的等价优化同样可验证设计参数的合理性，即当优化模型计算结果与设计条件偏差较大时，往往是由于方案选取的设计参数不符合工程实际，而通过模型计算可辅助设计参数的选取.

［1］周玉文.排水管网理论与计算［M］.1版.北京:中国建筑工业出版社，2000:126-173.

［2］GB 50014—2006室外排水设计规范［S］.北京:中国建筑工业出版社，2014.

［3］王磊.基于模型的城市排水管网积水灾害评价与防治研究［D］.北京:北京工业大学，2010.

［4］ROSSMAN A.Storm water management model user's manual［M］.Washington D C:USEPA，2010:33-55.

［5］ENGMAN E T.Roughness coefficients for routing surface runoff［J］.Journal of Irrigation and Drainage Engineering，1986，112(1):35-46.

［6］周玉文，余永琦，李阳，等.城市雨水管网系统地面径流损失规律研究［J］.沈阳建筑工程学院学报，1995(2):133-137.

［7］DURRANS S R.Durrans stormwater conveyance modeling and design［M］.Waterbury:Heastad Press，2003:477-527.

［8］HORTON R E.An approach toward a physical interpretation of infiltration-capacity［J］.Soil Science Society of America Proceedings，1940，5:399-417.

［9］陈贵敏，贾建援，韩琪.粒子群优化算法的惯性权值递减策略研究［J］.西安交通大学学报，2006(1):61-64.

［10］SHI Y，EBERHART R.Empirical study of particle swarm optimization［M］.Washington:Evolutionary Computation，1999:1945-1950.