基于Hydrus-1D模型的LID措施雨水径流控制效应研究

麦叶鹏，黄国如，解河海，曾碧球，马兴华

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

1 研究背景

试验和监测是探究低影响开发(low impact development，LID)措施对雨水径流控制规律和机理的必不可少的环节，但是进行大量工况下的试验和监测会耗费大量财力物力，因此基于物理模型试验结果，构建数值模型，并通过数值模型来模拟分析更加经济合理。模型模拟是优化LID措施设计和预测运行结果的有效手段，目前有不少研究[1-4]通过模型模拟或现场监测探究LID措施在地块或片区尺度的雨水径流控制效果，也有部分研究基于试验结果运用模型模拟分析LID单项措施的雨水径流控制效应。

对于LID单项措施模拟模型而言，常用的模型有SWC、Hydrus-1D、DRAINMOD和RECARGA。SWC模型可用于场次和长序列降雨条件下LID措施的雨水径流水量分析，一般适用于选用LID措施类型及其面积[5]；Hydrus-1D模型可用于场次和长序列降雨下可变饱和多孔介质中的水流运动和溶质运移分析[6]，因此可用于分析LID措施径流水量变化和径流污染物在填料中的迁移机理；DRAINMOD模型可通过长序列降雨数据来预测地下水位、排水速率和排水总量，以及模拟氮素的转化和盐分运移，但不能进行场次降雨分析[7]；RECARGA模型适用于模拟生物滞留设施水量变化过程，但其降雨间隔最小为1 h[8]。由于Hydrus-1D模型不仅可以模拟LID措施的径流水量和径流污染物变化过程，而且还可用于场次降雨分析，因此本研究选用Hydrus-1D模型展开LID单项措施雨水径流控制效应模拟研究。

李家科等[9]利用Hydrus-1D模型对不同情景下生物滞留池的水量削减效果进行模拟分析，结果表明人工填料层为粉煤灰+沙的生物滞留池对水量削减效果最佳，进入水量的增大会降低生物滞留池的水量削减效果。Jiang等[10]运用试验观测数据来构建生物滞留设施的Hydrus-1D模型，发现实测和模拟值之间的误差是由入流和入渗特性引起的。Hilten等[11]较早地将Hydrus-1D模型运用于绿色屋顶的模拟，以确定绿色屋顶在24 h设计降雨下的峰值流量和产流时间。Hakimdavar等[12]评估了Hydrus-1D模型预测相同绿色屋顶安装在不同空间尺度上的差异，结果表明其预测的绝对误差随降雨量增大而增加。Xie等[13]构建了绿色屋顶的Hydrus-1D模型和SWMM模型，并比较两者的产流量，结果表明Hydrus-1D的模拟结果总是大于SWMM模型。Qin等[14]基于Hydrus-1D模型开发了含有存水层的绿色屋顶数值模型，并用于绿色屋顶的土壤水分变化模拟，结果表明，在不进行灌溉的情况下，绿色屋顶一年有35%以上的时间遭受水分胁迫，且随着土层深度的增加，水分胁迫天数减少，但灌溉总量增加。Zhang等[15]通过Hydrus-1D模型模拟分析表明绿色屋顶的水量削减率随降雨量周期呈指数下降，但随基质厚度增加而增加。相对而言，运用Hydrus-1D模型进行下凹式绿地和透水铺装模型模拟的研究相对较少，比如张建锋等[16]仅用Hydrus-1D模型模拟分析了下凹式绿地的雨水径流污染物在土壤中的迁移消减变化规律，Brunetti等[17]运用Hydrus-1D模型模拟研究了透水路面的水力性能。此外，影响LID措施雨水径流控制效应的参数众多，仍需进一步结合区域土壤特性和降雨特征探究LID措施在不同情形下的雨水径流控制效应。因此，本文基于试验研究成果[18-19]，运用Hydrus-1D模型进一步探究生物滞留池、下凹式绿地和透水铺装在多种工况下的雨水径流控制效果和规律，以期为LID措施的规划和设计提供参考。

2 试验研究

2.1 试验装置生物滞留池试验装置长6.35 m、宽1.3 m、深1 m，填料主要由蓄水层(200 mm)、种植土层(250 mm)、透水土工布、石屑填料层(300 mm)、砾石排水层(250 mm)和基土层等构成。下凹式绿地试验装置长和宽均为2 m，填料主要由土壤基层(赤红壤，500 mm)和种植土层(200 mm)组成，雨水口高度为50 mm。下凹式绿地装置共设置4个等大的隔间，其中3个隔间分别种植玉龙草、马尼拉草、大叶油草，另一个隔间为裸土，本文主要根据大叶油草下凹式绿地的试验结果展开模型模拟研究。透水铺装试验装置长2 m、宽1 m，主要由砾石排水层(35 mm)、土壤基层(430 mm)、粗砂垫层(36 mm)、级配碎石层(364 mm)、找平层(20 mm)、透水砖(50 mm)等组成。各LID措施的试验装置填料结构图如图1所示。在试验过程中，生物滞留池和下凹式绿地采用水塔提供模拟雨水径流，而透水铺装采用人工模拟降雨系统提供模拟降雨。

图1 各LID措施试验装置填料结构

2.2 试验方案对于生物滞留池，主要探究在汇水面积比为8%时4个不同降雨重现期(即0.2 a、0.5 a、1 a和2 a)下的雨水径流控制效果；对于下凹式绿地，主要探究在汇水面积比为20%时4个不同降雨重现期(即0.5 a、1 a、2 a和5 a)下的雨水径流控制效果；对于透水铺装，主要探究小雨、中雨、大雨及暴雨等4种不同降雨强度下的雨水径流控制效果。不同LID措施试验对应的入流流量或降雨强度计算方法有所不同。对于生物滞留池，主要基于广州市降雨强度公式(式(1))计算模拟降雨，并利用芝加哥雨型进行雨型分配，然后结合雨水流量公式(式(2))计算出不同降雨重现期下的入流流量过程，其中生物滞留池的综合径流系数根据《室外排水设计标准》(GB50014—2021)确定为0.9。对于下凹式绿地，也基于广州市降雨强度公式计算降雨强度，并按照均匀雨型进行雨型分配，然后结合雨水流量公式计算出不同降雨重现期下的入流流量过程，并根据《室外排水设计标准》确定下凹式绿地的综合径流系数为0.85。对于透水铺装，则按照人工模拟降雨系统的不同等级降雨进行模拟降雨。不同LID措施试验过程中对应的入流流量或降雨强度如表1所示。

(1)

式中：q为降雨强度，L/(s·hm2)；P为重现期，a；t为降雨历时，min。

Q=ψqF

(2)

式中：Q为入流流量，L/s；ψ为综合径流系数；F为汇水面积，hm2。

表1 不同LID措施试验对应的入流流量或降雨强度

2.3 试验结果根据不同LID措施的试验方案展开试验，降雨历时或入流时长都为1 h。生物滞留池的入流过程按照芝加哥雨型(雨峰系数0.4)进行入流，下凹式绿地按照均匀雨型进行入流，而透水铺装按照均匀雨型进行降雨，并对生物滞留池和下凹式绿地试验过程中的入流和溢流以及透水铺装试验过程中的降雨强度和溢流进行测量。各LID措施在不同降雨情形下的溢流过程如图2所示。

图2 各LID措施在不同降雨情形下的溢流过程

3 Hydrus-1D模型原理与构建

Hydrus-1D是基于土壤水动力学原理的一维入渗模型，可模拟多孔介质中变饱和情况下的水分、溶质运移和热量传递过程[20]。该模型可以考虑多种汇源项的影响及边界条件，模拟最小时间间隔为1 s。构建Hydrus-1D模型需要输入包括土壤几何和时间信息、土壤水力参数、植被特性、初始和边界条件以及土壤剖面信息等。由于本文研究的降雨时长较短，因此忽略植物生长及根系吸水的影响，而只运用Hydrus-1D模型的水分运移模块来进行LID单项措施雨水径流水量控制效应模拟评估。模型的长度单位设定为cm，时间单位为min，质量单位为g，迭代条件采用默认值。

3.1 模型参数设定在Hydrus-1D中，多孔介质中的一维垂直水分运动采用Richards方程进行计算，如下：

(3)

式中：θ为土壤含水率；h为压力水头，cm；K为土壤导水率，cm/min；z为土壤深度向下为正，cm；s(h)为源汇项。

Hydrus-1D提供了多种模型来计算土壤水分特性，主要分为单孔模型和双孔模型两大类。本文选择目前应用较为广泛的单孔模型即van Genuchten模型来计算土壤水分特性，同时不考虑土壤中水流运动的滞后效应。van Genuchten模型基于土壤水分特征参数来预测土壤含水率和非饱和渗透系数，其计算公式如式(4)—(6)所示。

(4)

(5)

(6)

式中：θs为土壤饱和含水率；θr为土壤残余含水率；α为与土壤进气吸力有关的参数；n为土壤水分特征曲线的形状参数；Ks为土壤饱和渗透系数，cm/min；Se为相对饱和度；l为孔隙连通性系数，一般取值为0.5，无量纲。

采用ku-pF非饱和导水率测量系统测得石屑、赤红壤及基土等主要填料的土壤水分特征曲线，并拟合得到各填料及基土的土壤水分特征曲线van Genuchten模型拟合参数(如表2)，进而为各LID措施的Hydrus-1D模型提供土壤水力特征参数的初始数据。各LID措施水流模块的初始土壤含水率为试验前测得的各层填料含水率，上边界均为具有表层的大气边界，下边界均为自由排水边界。

表2 填料及基土的土壤水分特征曲线

图3 各LID措施的填料剖面图

3.2 填料剖面划分在Hydrus-1D的Soil Profile-Graphical Editor模块中依据LID措施各层填料的厚度进行一维剖面划分。各LID措施水流模块的下边界设定为自由排水边界，生物滞留池的剖面划分为赤红壤层、生物炭层、石屑层、砾石层和基土层；下凹式绿地的剖面则划分为种植土层和土壤基层；透水铺装的剖面则划分为透水砖、找平层、级配碎石层、粗砂垫层和土壤基层；下凹式绿地和透水铺装的土壤基层采用的是赤红壤进行铺填。各LID措施的填料剖面在Hydrus-1D模型中的几何信息如图3所示。

3.3 模型参数率定与验证根据试验测量数据进行模型参数的率定与验证，并选用纳什效率系数(NSE)、相对误差系数(Re)和决定系数(R2)来评价模型模拟值和实测值之间的拟合情况，各评价指标的计算公式如式(7)—(9)所示。NSE常用于水文模型参数的率定与验证，其反映了模型模拟值和实测值之间统计差异。NSE的取值范围为负无穷至1，当NSE≤0时为模型不可信，当0

(7)

(8)

(9)

根据试验过程中采用的降雨类型，生物滞留池的率定期工况为0.5 a和2 a，验证期工况为0.2 a和1 a；下凹式绿地的率定期工况为1 a和5 a，验证期工况为0.5 a和2 a；透水铺装的率定期工况为中雨和暴雨，验证期工况为小雨和大雨。采用试错法调整土壤水力参数，使得率定期时各项LID措施的NSE都高于0.65，相对误差绝对值控制在20%以内。经参数调整后，得到了符合要求的模拟结果，各LID措施在率定期不同降雨情形下的溢流的模拟与实测结果如图4所示。根据式(7)—(9)计算分析得到不同LID措施在率定期与验证期时实测值与模拟值之间的评价结果(表3)，进而确定不同LID措施的模型参数如表4所示。

图4 LID措施在不同降雨情形下的溢流的模拟与实测结果

表3 不同LID措施的参数率定与验证评价结果

表4 不同LID措施的模型参数

由表3可知，生物滞留池的模型模拟值与实测值之间的NSE在率定期与验证期均大于0.73，Re的绝对值小于17.9%，R2大于0.76；下凹式绿地的模型模拟值与实测值之间的NSE在率定期与验证期均大于0.69，Re的绝对值小于13.2%，R2大于0.79；透水铺装的模型模拟值与实测值之间的NSE在率定期大于0.75，在验证期大于0.57，Re的绝对值小于16.49%，R2大于0.74。以上结果表明所构建LID措施Hydrus-1D模型的模拟精度较高，可用于进一步的模拟评估。

4 雨水径流控制效应分析

基于所构建的各项LID措施Hydrus-1D模型，进行各项LID措施在不同情景方案下的雨水径流控制效应模拟评估，探索其雨水径流控制效应变化规律。对于生物滞留池，研究的影响因素主要为降雨重现期、汇水面积比(生物滞留池面积与汇水面面积之比)、蓄水层厚度和种植土厚度等关键参数。对于下凹式绿地，研究的影响因素主要为降雨重现期、雨水口高度、汇水面积比(下凹式绿地面积与汇水面面积之比)和种植土厚度等关键参数。对于透水铺装，研究的影响因素主要为降雨类型、降雨历时、雨峰系数和透水砖渗透系数等。模型模拟分析所用降雨数据由广州市暴雨强度公式(式(1))和芝加哥雨型计算得到。对于生物滞留池和下凹式绿地，雨峰系数取值为0.4，降雨历时取值为60 min，同时还需结合雨水流量公式(式(2))计算其入流流量。

4.1 生物滞留池基于已验证的生物滞留池Hydrus-1D模型，进行不同降雨重现期下，不同汇水面积比、不同蓄水层厚度和不同种植土厚度等情景下生物滞留池雨水径流控制效果的模型模拟。生物滞留池的种植土包含赤红壤层和生物炭层，在本研究中通过改变生物炭层下的赤红壤层厚度来调整种植土厚度，而表层的赤红壤层和生物炭层的厚度保持不变。《海绵城市建设技术指南——低影响开发雨水系统构建(试行)》(简称《指南》)[21]中规定生物滞留池的汇水面积比一般为5%～10%，蓄水层厚度一般为20～30 cm。此外，为满足植物生长，需要种植土厚度不宜小于20 cm。因此，生物滞留池的汇水面积比、蓄水层厚度及种植土厚度应设置在规定范围内。在模拟分析生物滞留池汇水面积比对径流量削减能力影响时，蓄水层厚度为20 cm，种植土厚度为25 cm；在模拟分析蓄水层厚度对径流量削减能力影响时，汇水面积比为10%，种植土厚度为25 cm；在模拟分析种植土厚度对径流量削减能力影响时，汇水面积比为10%，蓄水层厚度为20 cm。经模型模拟与计算分析，生物滞留池在不同降雨重现期、汇水面积比、蓄水层厚度和种植土厚度时的径流量削减率如图5所示。

图5 生物滞留池在不同情景下的径流量削减率

由图5可知，生物滞留池的径流量削减率随降雨重现期增大而减小。由图5(a)可知，生物滞留池的径流量削减率随汇水面积比增大而增大，因为汇水面积比越大，单位面积的生物滞留池所需处理的雨水径流就越少，径流量削减率也就越高。由图5(b)可知，生物滞留池的径流量削减率随蓄水层厚度增加而增加，其原因是蓄水层厚度越厚则生物滞留池所蓄滞和削减的水量越多。由图5(c)可知，在种植土厚度为20～40 cm时生物滞留池的径流量削减率随种植土厚度增加而减小，究其原因，虽然种植土越厚会使得土壤蓄水性能有所提升，但生物滞留池的种植土厚度增加使得其渗径长度也在增加，渗径长度增加则会使得渗透速度下降，从而导致下渗水量减少，最终使得径流量削减率随着种植土厚度增大而减小。为探究各降雨重现期下生物滞留池的汇水面积比、蓄水层厚度和种植土厚度与径流量削减率之间的函数关系，将模型模拟结果进行拟合分析，得到其拟合曲线方程如表5所示。

表5 生物滞留池的汇水面积比、蓄水层厚度和种植土厚度与径流量削减率之间的拟合曲线方程

由表5可知，不同降雨重现期下，汇水面积比和蓄水层厚度与径流量削减率之间的拟合曲线方程均为二元一次方程，其斜率均大于0且R2均大于0.999，说明生物滞留池的汇水面积比和蓄水层厚度与径流量削减率之间为正线性相关关系。在0.5～2 a降雨重现期，汇水面积比与径流量削减率拟合曲线的斜率为5.3066～6.9206，而蓄水层厚度与径流量削减率拟合曲线的斜率为2.058～2.6154，说明生物滞留池的径流量削减率随汇水面积比增大而增大的幅度比蓄水层厚度时的大。而种植土厚度与径流量削减率之间的拟合曲线方程均为二元一次方程且R2大于0.941，两者为负线性相关关系。上述结果说明生物滞留池的径流量削减率随着汇水面积比和蓄水层厚度增大而增大，但随着种植土厚度增大而减小。此外，由于考虑到植物耐淹性能和土壤渗透性能，《指南》[21]中建议生物滞留池的蓄水层厚度一般为20～30 cm，广州市《城市绿化工程施工和验收规范》[22]则规定草本植被的种植土厚度不小于30 cm。结合本文研究结果，建议生物滞留池的蓄水层厚度取值为30 cm，种植土厚度取值为30 cm。

4.2 下凹式绿地基于已验证的下凹式绿地Hydrus-1D模型，进一步探究不同降雨重现期、汇水面积比、雨水口高度和种植土厚度对下凹式绿地雨水径流控制效应的影响。在模拟分析下凹式绿地汇水面积比对径流量削减能力影响时，雨水口高度为7.5 cm，种植土厚度为20 cm；在模拟分析雨水口高度对径流量削减能力影响时，汇水面积比为15%，种植土厚度为20 cm；在模拟分析种植土厚度对径流量削减能力影响时，汇水面积比为15%，雨水口高度为7.5 cm。经模型模拟与计算分析，下凹式绿地在不同降雨重现期、汇水面积比、雨水口高度和种植土厚度时的径流量削减率如图6所示。

图6 下凹式绿地在不同情景下的径流量削减率

由图6可知，下凹式绿地的径流量削减率随降雨重现期增大而减小，随着汇水面积比、雨水口高度及种植土厚度增大而增大。由图6(a)可知，在低汇水面积比时下凹式绿地的径流量削减率较小，而在高汇水面积比时下凹式绿地的径流量削减率就很高。由图6(b)可知，下凹式绿地的径流量削减率随雨水口高度增加而增加，其原因是下凹式绿地的雨水口高度越高，则其所蓄滞的径流水量就越多，所以径流量削减率就越高。由图6(c)可知，下凹式绿地的径流量削减率随种植土厚度增加而增加，但增加的幅度较小。为探究各降雨重现期下下凹式绿地的汇水面积比、雨水口高度和种植土厚度与径流量削减率之间的函数关系，将模型模拟分析结果进行拟合分析，得到其拟合曲线方程如表6所示。

表6 下凹式绿地的汇水面积比、雨水口高度和种植土厚度与径流量削减率之间的拟合曲线方程

由表6可知，不同降雨重现期下，下凹式绿地的汇水面积比、雨水口高度和种植土厚度与径流量削减率之间的拟合曲线方程均为二元一次方程，其斜率均大于0且R2均大于0.998，说明下凹式绿地的汇水面积比、雨水口高度和种植土厚度与径流量削减率之间为正线性相关关系。在0.5～2 a降雨重现期下汇水面积比与径流量削减率拟合曲线方程的斜率为3.557～4.634，雨水口高度与径流量削减率拟合曲线方程的斜率为2.9364～3.8364，而种植土厚度与径流量削减率拟合曲线方程的斜率为0.324～0.4219，说明提高汇水面积比和雨水口高度可以较好的提高径流量削减率，而种植土厚度的提高效果不明显。此外，《指南》[21]中建议下凹式绿地的雨水口高度一般为5～10 cm，广州市《城市绿化工程施工和验收规范》[22]规定草本植被的种植土厚度不小于30 cm。因此结合本文研究结果，建议下凹式绿地的汇水面积比不宜小于15%，雨水口高度取值为10 cm，种植土厚度取值为30 cm。

4.3 透水铺装基于已验证的透水铺装Hydrus-1D模型，模拟评估不同降雨重现期、降雨历时、雨峰系数和透水砖渗透系数对透水铺装雨水径流控制效果的影响。在模拟分析降雨历时对透水铺装径流量削减能力影响时，雨峰系数为0.4，透水砖渗透系数为0.018 cm/min；在模拟分析雨峰系数对径流量削减能力影响时，降雨历时为60 min，透水砖渗透系数为0.018 cm/min；在模拟分析透水砖渗透系数对径流量削减能力影响时，降雨历时为60 min，雨峰系数为0.4。经模型模拟与计算分析，透水铺装在不同降雨重现期、降雨历时、雨峰系数和透水砖渗透系数时的径流量削减率如图7所示。

图7 透水铺装在不同情景下的径流量削减率

由图7可知，透水铺装的径流量削减率随降雨重现期和雨峰系数增大而减小，随着降雨历时和透水砖渗透系数增加而增加。由图7(a)可知，当降雨历时较长时，径流量削减率则较大。这是因为降雨历时越长，平均雨强则越小，入渗水量更多，所以径流量削减率更大。由图7(b)可知，透水铺装径流量削减率随着雨峰系数增大而降低，但降低幅度相对较小。其可能原因是前期的降雨使得透水铺装的填料逐渐饱和，而雨峰系数越大则降雨峰值越靠后，从而导致透水铺装的径流水量削减能力越差。由图7(c)可知，透水砖渗透系数在0.018～0.109 cm/min时，透水铺装径流量削减率随透水砖渗透系数增大而迅速增加。当透水砖渗透系数大于0.109 cm/min后，径流量削减率的增加趋势逐渐趋缓。当透水砖渗透系数为0.291 cm/min时，0.5～2 a重现期下的透水铺装削减率为98.5%～100%，说明当透水砖渗透系数大于0.291 cm/min时，透水铺装产生雨水径流较少。为探究各降雨重现期下透水铺装的降雨历时、雨峰系数和透水砖渗透系数与径流量削减率之间的函数关系，将模型模拟分析结果进行拟合分析，得到其拟合曲线方程如表7所示。

表7 降雨历时、雨峰系数和透水砖渗透系数与透水铺装径流量削减率之间的拟合曲线方程

由表7可知，降雨历时和透水铺装径流量削减率之间的拟合曲线方程均为二元一次方程，其斜率大于0且R2均大于0.99，说明降雨历时与透水铺装径流量削减率之间为正线性相关关系。雨峰系数和透水铺装径流量削减率之间的拟合曲线为二元一次方程，其斜率小于0且R2均大于0.965，说明两者之间为负线性相关关系。渗透系数和透水铺装径流量削减率之间的拟合曲线方程均为二元三次方程且R2均大于0.996，两者之间为正相关关系。此外，《透水砖路面技术规程》[23]规定透水砖的渗透系数不应小于0.6 cm/min。结合本文研究结果，建议透水砖的渗透系数设计值亦不应小于0.291 cm/min。

5 结论

(1)生物滞留池的径流量削减率随降雨重现期增大而减小，汇水面积比和蓄水层厚度与径流量削减率之间为正线性相关关系，在种植土厚度为20～40 cm时种植土厚度与径流量削减率之间为负相关关系。提高汇水面积比和蓄水层厚度可以明显地提高径流量削减率，而在20～40 cm时种植土厚度越小则径流量削减率越大。此外，根据研究结果和相关规定，建议生物滞留池的蓄水层厚度取值为30 cm，种植土厚度取值为30 cm。

(2)下凹式绿地的径流量削减率随降雨重现期增大而减小，随汇水面积比、雨水口高度和种植土厚度增加而增加。下凹式绿地的汇水面积比、雨水口高度和种植土厚度与径流量削减率之间为正线性相关关系。提高下凹式绿地汇水面积比和雨水口高度可以较好的提高其径流量削减率，而增加种植土厚度的提高效果不明显。此外，为更好地控制雨水径流，建议下凹式绿地的汇水面积比不宜小于15%，雨水口高度取值为10 cm，种植土厚度取值为30 cm。

(3)透水铺装的径流量削减率随降雨重现期和雨峰系数增大而减小，随着降雨历时和透水砖渗透系数增加而增加。降雨历时与透水铺装径流量削减率之间为正线性相关关系，雨峰系数和透水铺装径流量削减率之间为对数函数关系且为负相关关系，透水砖渗透系数和透水铺装径流量削减率之间为正相关关系。当透水砖渗透系数大于0.29 cm/min时，透水铺装产生的雨水径流较少。此外，根据本研究结果，建议透水砖的渗透系数设计值不应小于0.291 cm/min。