基于HYDRUS-2D的雨水集聚深层入渗系统土壤水分运移模拟

张 伟，赵西宁，高晓东，吴普特，潘岱立，宋小林，，杨世伟，姚 杰

(1.西北农林科技大学水利与建筑工程学院，陕西 杨凌 712100；2.西北农林科技大学中国旱区节水农业研究院，陕西 杨凌 712100；3.国家节水灌溉杨凌工程技术研究中心，陕西 杨凌 712100；4.延安市宝塔区果业局，陕西 延安 716000)

黄土高原由于独特的气候条件和疏松的土壤结构类型，被公认为优质苹果的生产区[1-3]，目前苹果产业已成为该区支柱性产业和农民的重要经济来源[4-5]。干旱缺水和水土流失是限制黄土高原地区生态和经济发展的主要因素[6-7]。降雨是该区旱作苹果园唯一水分输入[8-9]，但该区降雨稀少且分布不均，加之苹果树需水量大，果树水分供需矛盾十分突出[4-5,10]。已有研究表明，苹果树在生育期内缺水会显著降低产量和果实品质[11]，从而影响当地农业经济的发展。因此，增加果园的降雨利用、减少水土流失和土壤蒸发对旱作果园产业的可持续发展具有重要意义。

多年来，国内外学者针对旱作及灌溉果园的水分高效利用技术进行了大量研究。研究内容主要包括：不同覆盖方式下经济作物与果树间套作措施[12-14]；滴灌、微灌、喷灌等新型节水技术的应用；梯田[15]和鱼鳞坑[16]等集雨工程措施。这些研究均取得了一定的成果，有效增加了果园土壤表面降雨拦蓄量，显著减少果园地表蒸发。然而，鉴于黄土塬区土壤水分点尺度的入渗机制[17-19]，降雨入渗补给存在一定的滞后性[20-22]，到达土壤深层需要一定的时间且水量较少，难以对深层土壤水分形成有效补给[23]。由此，有学者在陕北旱作果园推广中应用了一种具有蓄水、保水和水肥一体化等优点的雨水集聚深层入渗系统(rainwater collection and infiltration systems, RWCI)，发现RWCI系统能够显著增加土壤含水率低值区土壤平均含水率[24-25]，能够显著提高0.2～1.0 m土壤的土壤含水量，且能够提高0～2 m土层的果树根系分布[26-27]。该系统对雨养果园土壤水分的影响机制与果树生长过程、根系分布特征以及不同生长阶段耗水特征等因素密切相关，但RWCI系统的水分空间运动尚不明晰，亟需进一步研究。

本文通过在陕西省延安市宝塔区万庄村试验点开展室内不同灌水量和不同RWCI设计深度条件下土壤水分入渗试验，旨在揭示RWCI系统土壤水分空间运动规律；在此基础上，采用HYDRUS-2D模型建立RWCI系统的土壤水分二维入渗模型[28-30]，分析土壤湿润锋运移过程和含水率的动态变化，旨在为RWCI技术的进一步改进和发展提供理论依据。

1 材料与方法

1.1 雨水集聚入渗系统介绍

RWCI系统是一种中深层尺度立体集雨灌溉技术[24-27]，该技术是在树冠下距树干一定距离和方位，挖掘面积为80 cm×80 cm，深度为40～80 cm范围的立方体集雨坑，在其中心安置一根多孔集雨管，集雨坑中用粉碎秸秆等有机类物质分层填充压实，至坑口处修成凹面状，覆盖黑色集雨膜(中间钻孔，用集雨管盖固定在集雨管口处)。当降雨时，雨水通过黑色集雨膜进入多孔集雨管，渗入坑内基质，最终沿渗水坑壁渗入根区土壤。其田间布设图及具体的操作技术详见文献[26]。该技术最大特点是通过收集雨水将其直接输送到作物根区供作物吸收利用，从而减少水分的蒸发损失，提高降水利用效率。

1.2 室内模拟试验

1.2.1 试验材料 试验在陕西省延安市宝塔区万庄村实验点进行。试验装置由有机玻璃土箱和供水系统两部分组成。有机玻璃土箱规格为100 cm×30 cm×120 cm，底部设若干排气孔，以防止气阻，试验采用马氏瓶进行定流量灌水。以RWCI系统为研究对象，截取以集水管为中心的半土体剖面进行试验，试验装置如图1所示。供试土壤采用当地苹果园的10～50 cm土层黄绵土，土壤质地颗粒组成为粘粒16.3%、粉粒25.65%、砂粒58.04%。集水管为直径10 mm的PVC管，为了保证水分均匀入渗，PVC管壁四周间隔2 cm均匀开Ф5 mm的孔隙，集水管周围填充有机材料(粉碎的秸秆)。

1.2.2 试验设计 试验方案设计见表1，采用正交试验理论设置不同坑深和灌水量，共8个处理，其中坑内集水管外无填充材料为对照处理(处理4和处理5)。试验前，将土样风干过2 mm筛，人工配置初始含水率为0.07 cm3·cm-3的土壤基质，箱内装土容重按设计容重1.35 g·cm-3每10 cm分层填筑，层间进行打毛处理，防止光滑面对土壤水分入渗产生影响；土壤表面用塑料薄膜覆盖，故不考虑室内土壤的蒸发，装土过程中埋设土壤水分传感器，装土完成后静置1天使土壤含水率分布均匀。

试验开始时，开挖设计深度的土槽装填粉碎秸秆，调节马氏瓶到一定高度进行定流量灌水，当灌水量增大时，灌水时间随之增大。试验开始灌水后3 min描绘第一次湿润锋，其后每间隔5 min描绘一次，试验进行一段时间后根据湿润锋运移距离和入渗速率适当延长湿润锋描绘时间。试验结束后测量湿润锋径向和垂向侧渗距离，试验过程中采用EM50每间隔1 min进行土壤水分的动态监测。

图1 试验装置示意图(单位: cm)Fig.1 Diagram of experiment device (unit: cm)

表1 试验设计

1.3 数值模拟

1.3.1 基本方程 设计试验以截取集雨管为中心的半土体剖面为研究对象，其土壤水分入渗过程方程可以简化为二维入渗过程，采用Richards方程进行描述[19]。

(1)

式中,θ为体积含水率(cm3·cm-3);t为入渗时间(min);h为基质势(cm);r、z分别为径向和垂直坐标值(cm);K(h)为土壤非饱和导水率(cm·min-1)。

土壤基质势h、非饱和导水率K(h)与含水率的关系采用Van Genuchten-Mualem模型[31]进行描述。

(2)

(3)

(4)

式中,θr为残余含水率(cm3·cm-3);θs为饱和含水率(cm3·cm-3);Ks为饱和导水率(cm·min-1);m,n,α为经验系数,m=1-1/n。

1.3.2 初始和边界条件 模型中初始含水率设置如图2：根据试验实测值，秸秆初始含水率设置为0.03 cm3·cm-3，初始土壤含水率设置为0.07 cm3·cm-3；上边界AH、HG和GF用塑料薄膜覆盖，CD和EF为有机玻璃，均无水量交换，因此设置为零通量边界；DE为若干排气孔的下边界，设置为自由排水边界。定流量进行灌水时，多孔集雨管边界AB的水头由零增大到峰值a后灌水停止，水头随着时间逐渐又变为零的变边界面和变水头这样一个复杂过程，当灌水量增大，灌水时间延长，峰值a也相应增大。模型为简化这一过程，将模拟中多孔集雨管AB边界的长度假设为最大峰值0.5a的变水头边界；BC边设置为变水头边界。土壤水分传感器的布设如图2所示，图中每个格子间隔10 cm，探头主要集中在水分分布区域。

1.3.3 模型参数获取 根据确定的土壤水分运动模型，在HYDRUS-2D中对所建模型进行数值模拟。土壤水力学参数获取是利用土壤颗粒级配和容重，采用ROSETTA预测初始的土壤水力参数，结合离心机测得的实际土壤水分特征曲线，用RETC软件拟合得到的土壤的水力学参数，得出优化后的土壤水力学参数，以文献[32]中秸秆的水力学参数为初始的秸秆水力学参数(表2)，利用inverse模块反演求得所建模型中秸秆的水力参数。

1.3.4 模型评价标准和统计分析 本研究用相对均方差RE[33]、平均绝对误差MAE和纳什系数NE[34-37]3个指标评价模型模拟效果，各指标计算方程如下：

(5)

(6)

(7)

式中,Pi和Oi分别为实测值和模拟值;n为实测值总数;Pm为实测值均值。RE和MAE(评价湿润锋和土壤含水率单位分别为cm和cm3·cm-3)的数值越接近于0,NE约接近1, 表示模拟值与实测值差异越小, 模拟效果越好[33,38]。

1.4 试验数据处理

模型验证：试验各处理中土壤含水率的测点均有15个，本研究选取4个具有代表性的测点(测点2、6、8和14)进行模拟与实测值评价分析；湿润锋对比分析：选取两种不同设计坑深(40 cm、60 cm)、不同灌水量(21 L、36 L)和有无填充材料条件下土壤水分随时间的入渗过程，分析不同时间下土壤水分湿润锋动态运移；土壤含水率变化分析：由于试验各处理土壤水分观测点相对较多，取处理2中的4个具有代表性的土壤水分测点进行实测值和模拟值对比分析。

数据分析采用Excel 2010，利用Sigmaplot 12.5软件进行图像绘制。

2 结果与分析

2.1 模型效果验证

RWCI系统下，室内模拟试验与HYDRUS-2D模型模拟得到的湿润锋和土壤含水率对比如图3所示。从图3(a)、3(b)可以看出，湿润锋在垂向和径向的模拟与实测行进速度和趋势的吻合度均相对较好，垂向和径向湿润锋的RE、MAE和NE分别为0.019、0.011 cm、0.994和0.018、0.851 cm、0.977；径向湿润锋相对于垂向湿润锋的MAE明显偏大，模型低估了径向湿润锋的距离，原因可能是模型秸秆填充材料的水力学参数设置与实际情况存在一定差异；RE均小于0.02，NE均大于0.95，故认为模型对湿润锋的模拟效果较好。

图2 模型边界条件及传感器布设示意图(单位：cm)Fig.2 Schematic diagram of model boundary conditions and sensor layout (unit: cm)

表2 土壤水分特征参数

图3 湿润锋及土壤含水率实测与模拟对比关系图Fig.3 Comparison of observed and simulated wetting front and soil water content

将试验土壤含水率模拟值与实测值评价分析，得出土壤含水率的RE、MAE和NE分别为0.188、0.016 cm3·cm-3和0.916。从图3(c)可以看出模型对土壤含水率的模拟效果相对较差，在某些时刻低估或者高估了土壤含水率值，原因可能是模型中未考虑土壤水分特征曲线的滞后现象和温度、空气对土壤水分运动的影响，相关研究同样出现了土壤水分的模拟值和实测值存在较大差异的现象[39-40]。总体而言，所建模型能够较好地对湿润锋和土壤含水率进行模拟，能够用来描述RWCI系统的土壤水分运移变化规律。

2.2 不同处理条件下湿润锋的变化

图4为不同设计坑深(40 cm、60 cm)、不同灌水量(21 L、36 L)和有无填充材料条件下土壤水分随时间的入渗过程，湿润锋曲线运移变化的时间为3, 5, 10, 15, 20, 25, 30, 40, 60, 80, 100, 160, 280 min。

从图4(a)、(d)、(c)、(f)可以看出，当灌水量相同时，秸秆填充处理能够明显地减小径向和垂向湿润锋的运移距离，说明秸秆中储存了一定的水量，由于秸秆的吸水效应，减少了土壤净入渗水量，缩小了水分的扩展空间。在果园的实际应用中，秸秆中的水分有利于果园培肥土壤，因此一定的秸秆施用也是有必要的。在相同灌水量下，60 cm设计深度比40 cm土壤水分入渗更深，但湿润体的体积范围与40 cm之间无明显的差异，因此，较深的集雨坑更有利于土壤水分入渗至深层土壤，增加深层土壤水分含量，最大程度增加了果树在深层的根系有效生长，进而可以提高深层土壤的根系分布。在相同的设计深度(40 cm)下，增加灌水量能够增加土壤湿润体的体积，当灌水量为36 L时，湿润体在垂向上的运移已经扩散到地表，而设计深度60 cm的湿润体垂直向上的运移距离距地表仍有20 cm左右，说明设计坑深较浅的RWCI系统在降雨量较大时，与设计深度较深(60 cm)的RWCI系统相比，明显加大了土壤水分的无效蒸发，相应减小水分有效利用率，而果树根系主要分布区域在0～100 cm[41]，因此，设计深度为60 cm RWCI系统较设计深度为40 cm RWCI系统在不同灌水量下均能促进果树根系对土壤水分的有效利用。

图4 各处理不同时刻湿润锋分布Fig.4 Wetting front distribution at different treatment time

2.3 垂直和水平湿润锋变化规律

图5为有填充材料、坑深分别为40 cm和60 cm下三种不同灌水量时湿润锋在水平和垂直方向运移距离实测值和模拟值对比。从图5可以看出在相同灌水量下，设计深度对径向湿润锋运移影响不显著(P>0.05)，但是随着设计深度的增加，湿润锋在垂直方向上的分布范围相应下降，说明设计坑深对土壤水分在垂直方向上的运移存在显著影响(P<0.05)。随着灌水量的增大，湿润锋在垂向与径向的运移距离差异逐渐增大，这主要是由于土壤水分运动在垂直方向受基质吸力和重力势的共同作用，而在径向方向只受基质吸力的作用，所以湿润锋的垂向下移速度大于径向运移速度，且在垂向与径向的运移距离差异也逐渐增大。

2.4 土壤含水率随时间的变化

图6表示处理2中观测点2、6、8和14的实测和模拟土壤水分变化情况。从图6中可以看出在不同观测点中土壤含水率的模拟值与实测值曲线的拟合度较好，测点2、6和8的土壤含水率均是增加到一定的值后产生水分的消退过程；由于测点14距离灌水面相对较远，土壤水分增加后没有产生消退过程。观测点2在20～90 min对土壤水分的模拟效果相对较差，高估了这段时间的水分入渗，可能是因为模拟开始时模型对变水头边界条件的假设水头值大于试验开始时实际的水头值，但是总体的模拟趋势相近，因此可认为模型能够较好地反映湿润体内土壤含水率的分布规律。

2.5 不同坑深条件下土壤剖面含水率分布

图7表示径向距离50 cm处垂直剖面模拟土壤含水率变化。从图中可以看出RWCI系统剖面的土壤含水率分布规律主要受设计坑深和灌水量的影响。剖面土壤含水率随着灌水量的增加而递增，坑深为40 cm，水平距离50 cm处的垂直剖面土壤含水率最大峰值出现在垂向40 cm左右处；而坑深为60 cm，水平距离50 cm处的垂直剖面土壤含水率最大峰值出现在垂向60 cm左右处，随着坑深的增加，土壤含水率峰值位置也相应下延。

图5 实测与模拟土壤湿润锋变化对比Fig.5 Comparisonof observed and simulated soil wetting distances

图8表示垂直距离70 cm处径向剖面模拟土壤含水率变化。从图中看出垂向70 cm处的各处理径向方向的土壤含水率最大值均在径向距离起始点，随着径向距离的增加呈递减的趋势。在垂直距离70 cm处径向剖面的土壤水分含量均受设计深度和灌水量的共同影响，在相同的设计深度下，灌水量越大，剖面的土壤含水率越高。

图6 不同测点土壤含水率实测值与模拟值比较Fig.6 Comparisonof observed and simulated soil water contents at different observing points

图7 径向距离50 cm处垂直剖面模拟土壤含水率变化Fig.7 Changes of simulated soil water content in a vertical profile at radial distance of 50 cm

图8 垂直距离70 cm处径向剖面模拟土壤含水率变化Fig.8 Changes of simulated soil water content in a radial profile at vertical distance of 70 cm

4 结 论

(1)采用HYDRUS-2D模拟RWCI系统土壤水分入渗规律，得出模拟湿润锋与土壤含水率和实测值结果均较理想：垂向湿润锋相对均方差(RE)、平均绝对误差(MAE)和纳什系数(NE)分别为0.019、0.011 cm和0.994，径向湿润锋RE、MAE和NE分别为0.018、0.851 cm和0.977，土壤含水率RE、MAE和NE分别为0.188、0.016 cm3·cm-3和0.916。构建的模型较为合理，能够用于模拟RWCI系统在不同灌水量的变水头边界条件下的土壤水分运动情况。

(2)设计坑深较深的RWCI系统在灌水量较大时较设计深度较浅的RWCI系统明显使得土壤水分入渗深度增加，这为减少土壤的无效蒸发，提高水分利用效率创造了有利条件；当灌水量较低时，两种设计深度RWCI系统均能较好地将土壤水分保持在果树主要根系分布层，设计深度为40 cm RWCI系统更加经济合理，而设计深度为60 cm RWCI系统则适用于较大灌水量时土壤水分的有效利用。相同灌水量下，设计深度对径向湿润锋运移影响较小，但随着设计深度增加，湿润锋在垂直方向上的运移存在明显差异；且当灌水量增大时，湿润锋的垂向下移速度大于径向运移速度，且在垂向与径向的运移距离差异也逐渐增大。

以上研究结果可为RWCI系统在旱作果园的田间布设提供一定的理论依据，但试验设计中没有考虑果树根系吸水和蒸发等条件的影响，后期研究将结合不同果树品种的根系分布深度和RWCI系统在实际应用中的工程造价，确定不同条件下适宜的RWCI技术参数。