我国黄河流域大部分区域均处在干旱半干旱地区，该地区雨量稀少而集中，水资源严重不足^[1-2]，且水土流失十分严重，致使黄河成为世界上含沙量最高的河流^[3-4]。为解决水资源短缺问题，该地区一些灌区结合自身流域特点及生产实际，开展了引浑淤灌^[5-6]。多年来，引浑淤灌的实践证明，浑水灌溉不仅缓减了该地区的旱情，而且显著提高了土壤肥力，改善土壤结构，促进作物增产提质^[7-8]。浑水灌溉与清水灌溉的本质区别在于灌溉水中挟带泥沙，浑水入渗过程中泥沙颗粒易在土壤表层沉积形成致密层使得土壤结构发生变化，从而改变了水分入渗的上边界条件，影响了土壤的入渗机制^[9-11]，因此现有的清水入渗条件下的有关灌溉研究成果已不再适用于浑水条件下应用，浑水入渗机理亟待研究。

浑水入渗过程中致密层的形成，将土壤变成了上实下松的双层结构，其中土壤上层致密层是整个水分入渗的控制层，对水分入渗起到关键作用，即浑水灌溉条件下，土壤的入渗能力将主要受浑水泥沙致密层的影响^[12]。致密层的形成包含两种性质不同的现象：一种是泥沙颗粒沉积形成落淤层，落淤层对水分入渗的影响增加了入渗介质层，延长了水分进入土体的路径；另一种是泥沙颗粒滞留形成滞留层，入渗滞留改变的是介质层内颗粒级配，使得土壤孔隙发生变化，减小了水分入渗通道^[13]。关于浑水入渗的研究，目前国外相关报道较少，目前有文献报道的仅有Bouwer等^[14]曾针对大池塘渗漏问题对浑水入渗规律进行了探讨，主要研究了浑水泥沙沉积中的泥沙颗粒运动状况对渗流的影响特性。国内学者对浑水入渗问题的研究重点在于落淤层对入渗的影响，没有考虑滞留层的影响，曹惠提等^[15]通过一维垂直入渗试验，探讨了浑水土壤一维垂直入渗特性及其影响因素，指出浑水降低了土壤入渗能力，这主要与浑水含沙率和泥沙中物理性黏粒含量有关；费良军和王文焰^[16]基于不同的已知条件和影响因素，提出了3个浑水间歇入渗模型；汪志荣等^[17]提出了基于不同含沙率的浑水波涌灌入渗的Green-Ampt模型；费良军等^[18]、钟韵等^[19]、刘利华等^[20]、张全菊等^[21]对浑水膜孔灌入渗特性及影响因素进行了一系列研究，指出浑水中泥沙对土壤入渗能力影响显著。以上针对浑水入渗的研究均只简单阐述浑水入渗特性及其影响因素，而未涉及到田面地表土壤孔隙、泥沙沉积致密层形成与发展规律等方面，而浑水灌溉的入渗机制较清水更为复杂。因此，本文拟基于现有研究成果，开展浑水含沙率对一维垂直入渗特性及致密层形成特性的影响研究，以期为进一步深入探讨浑水灌溉提供理论依据。

1 材料与方法 1.1 试验材料

试验在西安理工大学西北旱区生态水利工程国家重点实验室进行，试验土样为粉砂壤土(美国制标准USDA)，采自西安市灞桥区，风干碾碎后过2 mm筛，土壤初始质量含水量为0.03 kg/kg，土壤饱和质量含水量为0.37 kg/kg，物理性黏粒含量为25.06%，饱和导水率为0.017 2 cm/min，田间持水量为0.26 kg/kg，土壤颗粒组成通过Mastersizer-2000型激光粒度分析仪测定，粒径d在d≤0.002mm，0.002mm < d≤0.05mm，0.05mm < d≤2 mm分别占12.58%、74.98%、12.44%，各粒径含量为体积分数。将从黄河流域泾惠渠灌区渠底取得的泥沙风干后过1 mm筛后，通过称重法人工调配成质量含沙率分别为3%、6%、9% 和12% 的浑水，其泥沙颗粒级配组成见表 1。

1.2 试验方法

浑水一维垂直入渗试验装置由两部分组成：土柱和马氏瓶。试验土柱采用内径为12 cm、高60 cm的有机玻璃制成，壁厚10 mm。试验前将土柱内壁清理干净，土样按照干容重1.40 g/cm³分层(5 cm)称重均匀装入土柱，采用夯实工具将土壤压实到预定高度后，在土壤表面覆膜，静置24 h。试验过程中通过加装有自动搅拌系统的马氏瓶调整水室水头为2.5 cm，使浑水的含沙率维持在稳定的范围。入渗结束后每隔1 cm取土30 g左右，烘干碾碎过2 mm筛后，测定各层的土壤颗粒组成。每组试验均进行3次重复，取其平均值进行结果分析。

1.3 测定项目

浑水灌溉后土壤的结构见图 1所示，本文将致密层分成落淤层和滞留层两个部分来分析。

1) 累积入渗量：按照先密后疏的时间间隔通过马氏瓶读取。

2) 湿润锋运移距离：由于所用风干土壤颜色较淡，湿润体颜色显著加深，因此，在透明有机玻璃外壁上定时描绘出不同入渗时刻所对应的湿润锋的位置及形状，试验结束后进行测量。

3) 落淤层厚度：土样均匀装入试验土柱后，在透明有机玻璃土柱外壁标记上表土位置，将土柱表土处的标记作为表土“本体”和落淤层“客体”的分界线；入渗开始后，每隔5 min使用精度为0.01 mm的游标卡尺测量从表土标记处到落淤层表面的距离即为落淤层厚度，待厚度变化稳定后逐渐加大时间间隔，直至入渗结束。

4) 落淤层、滞留层颗粒组成：采用Mastersizer- 2000型激光泥沙粒度仪测定。

1.4 数据处理与分析

利用统计学中均方根误差(RMSE)、偏差百分比(PBIAS)和平均绝对误差(MAE)3个指标，对计算值与试验观测值之间的符合度进行评价分析。其中，RMSE和MAE越接近0、PBIAS < ±10% 时，表示计算值与试验观测值差异越小，两者吻合越好。各指标定义如下：

$ {\text{RMSE}} = \sqrt {\frac{1}{{n - 1}}\sum\limits_{i = 1}^n {{{\left( {{x_i} - {y_i}} \right)}^2}} } $

(1)

$ {\text{PBIAS}} = \frac{{\sum\limits_{i = 1}^n {\left( {{x_i} - {y_i}} \right)} }}{{\sum\limits_{i = 1}^n {{x_i}} }} \times 100\% $

(2)

$ {\text{MAE}} = \frac{{\sum\limits_{i = 1}^n {\left| {{x_i} - {y_i}} \right|} }}{n} $

(3)

式中：x_i为第i个观测值；y_i为第i个计算值；n为数据总个数。

2 结果与讨论 2.1 含沙率对入渗特性的影响 2.1.1 含沙率对入渗能力的影响

反映土壤入渗能力的指标为累积入渗量和入渗率。由图 2、图 3可以看出，不同浑水含沙率的累积入渗量随着入渗历时的延长不断增大，而入渗率则不断变小，最终趋于稳定。入渗历时相同时，累积入渗量和土壤入渗率随浑水含沙率的增加而减小。当浑水含沙率由3% 增大至12% 时，入渗结束时其累积入渗量由7.90 cm减少至4.90 cm，降幅为37.97%；相对应的稳定入渗率由0.021 6 cm/min降低至0.010 3 cm/min，降幅为52.31%，表明含沙率对浑水一维垂直入渗能力影响显著。这与文献[22-23]的研究结论类似，其中王全九等^[23]将浑水泥沙对一维垂直入渗能力的减渗原因归咎于浑水泥沙对湿润锋处吸力的减小，并且研究表明湿润锋处平均吸力与浑水含沙率呈负相关关系。土壤浑水入渗特性与土壤自身特性和浑水含沙率密切相关^{[20, 24]}，浑水入渗时，泥沙颗粒逐渐在土壤表面挂淤落於，沉积在土壤表层(落淤层)，并伴随着入渗水将细颗粒带入土壤孔隙中，浑水泥沙滞留减少了表层土壤的孔隙率(滞留层)，形成一个密度较大的致密层^[23]。这个致密层阻渗作用强，而浑水中含沙率越大，该致密层形成越快、越密实，对土壤的入渗能力影响越大。

研究表明，不同浑水含沙率的累积入渗量与入渗历时的关系曲线均符合Kostiakov入渗模型。

$ I = K{t^a} $

(4)

式中：I为累积入渗量(cm)；K为入渗系数(cm/min)，主要反映入渗开始后第一个单位时段内的入渗量；a为入渗指数，主要反映了土壤入渗能力的衰减速度，a值越大，土壤入渗能力衰减的速度越快；$ t $为入渗历时(min)。

利用式(4)拟合累积入渗量与入渗历时的关系曲线结果见表 2，拟合的决定系数均接近于1，表明利用Kostiakov入渗模型建立的相关关系显著。从表 2可以看出，随着浑水含沙率的增大，入渗系数K减小，表明随着浑水含沙率的增大，入渗开始后第一个单位时段内浑水单位面积累积入渗量减少；而入渗指数a随之逐渐增大，表明随着入渗时间的延长，浑水含沙率越大，其入渗率衰减的速率较快，最终达到了稳渗阶段。研究表明，入渗系数K与入渗指数a随着含沙率ρ的变化十分明显，其关系为^{[22, 25]}：

$ K/{K_0} = 1 - e{\rho ^f} $

(5)

$ \alpha /{\alpha _0} = 1 + m{\rho ^n} $

(6)

式中：K₀、a₀分别为清水的入渗系数和入渗指数；ρ为浑水含沙率，%；e、f、m、n为拟合参数。

对表 2数据进行拟合得：

$ K/{K_0} = 1 - 0.102\;5{\rho ^{0.6851}}\;\;\;\;\;\;\; R^{2}=0.997 5 $

(7)

$ \alpha /{\alpha _0} = 1 + 0.000\;3{\rho ^{1.8625}}\;\;\;\;\;\;\; R^{2}=0.995 5 $

(8)

经分析，入渗系数K和入渗指数a与浑水含沙率关系拟合的决定系数均接近1。将式(7)和式(8)代入式(4)得到以浑水含沙率和入渗历时为自变量的浑水一维垂直入渗的累积入渗量模型，即：

$ I\left( t \right) = {K_0}(1 - 0.102\;5{\rho ^{0.685\;1}}){t^{{\alpha _0}(1 + 0.000\;3{\rho ^{1.862\;5}})}} $

(9)

为评价该经验模型的精确度，利用相同试验条件下浑水含沙率ρ为5% 的试验资料验证，通过式(9)计算该浑水含沙率的累积入渗量，试验实测值与模型计算值对比结果见表 3。可以看出，试验实测值与模型计算值的相对误差均在±10% 以内，两者的RMSE、PBIAS和MAE分别为0.14 cm、-2.40% 和0.11 cm，RMSE和MAE都接近0，|PBIAS| < 10%，表明该经验模型计算精度较高。

2.1.2 含沙率对湿润锋运移的影响

由图 4可以看出，入渗历时相同时，湿润锋运移距离随浑水含沙率的增加而减小，如浑水含沙率由3%增大至12% 时，入渗结束时其湿润锋运移距离由21.15 cm减少至17.75 cm，降幅为16.08%，表明浑水含沙率对湿润锋运移距离影响显著，这与白瑞等^[11]的研究结果相似。这主要是由于浑水入渗时，浑水中泥沙逐渐在土壤表面挂淤落於，形成一个密度较大的致密层，而浑水含沙率越大，该致密层形成越快越密实，对水分的阻渗作用越强，因此导致其湿润锋运移距离越小。

经分析，浑水入渗湿润锋运移距离与入渗历时之间符合幂函数规律，即：

$ Z = b{t^c} $

(10)

式中：Z为湿润锋运移距离(cm)；b为湿润锋拟合系数；c为湿润锋拟合指数。

利用式(10)对图 4中实测数据拟合的结果见表 4。经研究表明，拟合系数b、拟合指数c与浑水含沙率ρ关系为^{[22, 25]}：

$ b/{b_0} = g{{\text{e}}^{h\rho }} $

(11)

$ c/{c_0} = j{{\text{e}}^{k\rho }} $

(12)

式中：b₀、c₀分别为清水入渗湿润锋的拟合系数和指数；g、h、j、k为拟合参数。

对表 4数据拟合得：

$ b/{b_0} = 0.906\;3{{\text{e}}^{ - 0.008\;9\rho }} \;\;\;\;\;\;\;R^{2}=0.981 7 $

(13)

$ c/{c_0} = 1.002\;2{{\text{e}}^{ - 0.002\;0\rho }} \;\;\;\;\;\;\;R^{2}=0.991 0 $

(14)

经分析，拟合系数b和拟合指数$ c $与浑水含沙率关系拟合的决定系数均接近于1。将式(13)和式(14)代入式(10)中得到以浑水含沙率和入渗历时为自变量的浑水一维垂直入渗的湿润锋运移距离模型，即：

$ Z\left( t \right) = {b_0}(0.906\;3{{\text{e}}^{ - 0.008\;9\rho }}){t^{{c_0}(1.002\;2{{\text{e}}^{ - 0.002\;0\rho }})}} \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;0≤t≤210 $

(15)

为评价该经验模型的精确度，利用相同试验条件下浑水含沙率ρ为5% 的试验资料验证发现，试验实测值与模型计算值的相对误差均在±11% 以内，两者的RMSE、PBIAS和MAE分别为0.48 cm、-3.14%和0.39 cm，RMSE和MAE都接近0，|PBIAS| < 10%，表明该经验模型计算精度较高。

2.2 含沙率对致密层形成特性的影响 2.2.1 含沙率对落淤层厚度的影响

由图 5可以看出，不同含沙率对应的落淤层厚度均随入渗历时的增大而不断增大；相同入渗历时，落淤层厚度随浑水含沙率的增加而增大，入渗结束时，含沙率为3%、6%、9%、12% 的落淤层厚度分别为3.02、3.56、4.51、5.29 mm。入渗初期(0 ~ 20 min)，落淤层厚度较小，这主要是由于浑水中泥沙颗粒伴随着入渗水将进入到土壤孔隙中，主要以滞留层形式进入到土壤表面，所以落淤层厚度小得无法测量。入渗中期(20 ~ 130 min)的落淤层厚度增加较快，随着入渗历时的增大，厚度增加速率逐渐变小，这是由于随着浑水累积入渗量的增加，土壤孔隙逐渐被泥沙颗粒填充，浑水中泥沙小颗粒继续洗入到土壤孔隙中，以滞留形式进入到土壤，而浑水中较大泥沙颗粒则沉积在土体表面形成落淤层，入渗过程中落淤层逐渐加厚，阻塞入渗通道，并延长了水分进入土体的路径，使得土壤导水能力变弱，进一步阻碍浑水入渗，由于浑水入渗量减小，供水速度减慢，导致落淤层厚度增加速率逐渐变小。入渗后期(130 ~ 210 min)的落淤层厚度随着入渗时间稳定增加，此时土壤表面致密层已经形成，泥沙颗粒主要以落淤层形式沉积于土壤表面，浑水入渗能力随入渗时间减小的幅度较小，所以落淤层厚度增幅也逐渐趋于基本稳定。

研究表明，浑水一维垂直入渗落淤层厚度随入渗时间的变化曲线符合幂函数关系，根据两者的关系构造如下幂函数模型：

$ H = g\;{(t - 20)^j} $

(16)

式中：H为落淤层厚度(mm)；g、j为拟合系数和拟合指数。

利用式(16)对图 5中实测数据拟合的结果见表 5。可以看出，拟合的决定系数均接近于1，拟合系数g和拟合指数j均随着浑水含沙率增大而不断增大，经分析发现拟合系数g和拟合指数j与浑水含沙率之间均呈线性正相关关系。

$ g = 0.0274\;\rho + 0.2940\;\;\;\;\;\; R^{2}=0.990 5 $

(17)

$ j = 0.0004\rho + 0.4011 \;\;\;\;\;\;R^{2}=0.991 7 $

(18)

将式(17)和式(18)代入式(16)可得到以浑水含沙率和入渗历时为自变量的浑水一维垂直入渗的落淤层厚度计算模型，即：

$ H = (0.027\;4\;\rho + 0.294\;0){(t - 20)^{0.000\;4\;\;\rho + 0.401\;1}} \\ 0≤t≤210 $

(19)

利用相同试验条件下浑水含沙率ρ为5% 的试验资料对该经验模型进行验证发现，试验实测值与模型计算值两者的相对误差均在±15% 以内，两者的RMSE、PBIAS和MAE分别为0.17 mm、–6.14% 和0.15 mm，RMSE和MAE都接近0，|PBIAS| < 10%，表明该经验模型计算精度较高。

2.2.2 含沙率对致密层颗粒组成的影响

由图 6A的落淤层土壤颗粒组成含量曲线可知，在灌水结束时，不同浑水含沙率的落淤层中土壤细颗粒含量均低于浑水泥沙，大颗粒含量均高于浑水泥沙，表明浑水入渗过程中存在滞留现象，细颗粒进入了原土壤中，即损失了部分细颗粒的泥沙沉积在表层土壤形成落淤层。以浑水泥沙的颗粒组成为对照，含沙率为3%、6%、9%、12% 的落淤层对应的黏粒相对含量减幅依次为3.82%、25.95%、29.01% 和48.09%，表明随着浑水含沙率的增大，滞留现象越明显，落淤层细颗粒损失越多使得其相对含量越少；沙粒相对含量增幅依次为17.04%、37.53%、42.89% 和52.56%，表明随着浑水含沙率的增大，滞留现象越明显，落淤层粗颗粒相对含量越多。落淤层土壤颗粒分析结果见表 6，可以看出，不同浑水含沙率的落淤层土壤物理性黏粒含量显著低于浑水泥沙，以浑水泥沙的颗粒组成为对照，浑水含沙率3%、6%、9%、12% 的落淤层对应的物理性黏粒含量d_0.01减幅依次为16.72%、42.55%、46.30% 和61.50%，表明浑水中泥沙含量越大，细颗粒伴随着入渗水洗入到土体的越多，其沉积在土壤表面形成落淤层的相对含量越少。

图 6B ~ 6D分别为滞留层土壤颗粒组成含量曲线。可以看出，在灌水结束时，不同浑水含沙率的滞留层中土壤细颗粒相对含量均高于原土壤，大颗粒相对含量均低于原土壤，这再一次证明浑水入渗过程中，泥沙不仅会沉积在土壤表面，还存在滞留现象，部分细颗粒会伴随着入渗水洗入土体滞留在表层土壤孔隙中，改变了表层土壤的颗粒级配组成，因此土壤表层的滞留层细颗粒含量高于原土壤，这与姚雷^[10]的研究结果类似。当前学者们对浑水入渗机理的研究重点在于落淤层对入渗的影响，没有考虑滞留层的影响，均认为浑水灌溉时一部分泥沙将在流动过程中不断地挂淤、落淤在粗糙的田面上^[23]，最终浑水中的泥沙也将全部“滞留”于田面；这种说法实际上指的都是落淤层，其中所提到的“滞留”并不是说泥沙伴随水分入渗，而是浑水中水分入渗后，在土壤上方的水分含沙量增加，由此引起泥沙直接沉积在土壤表面^[10]。浑水含沙率对滞留层土壤颗粒组成影响较大，随着浑水含沙率的增加，泥沙细颗粒伴随着入渗水洗入到土壤越多，使滞留层细颗粒相对含量增多，粗颗粒含量相对减少，特别是在入渗土壤深度为0 ~ 1 cm处。以原土壤颗粒组成为对照，浑水含沙率为3%、6%、9% 和12% 的滞留层0 ~ 1 cm入渗深度处对应的土壤黏粒含量增幅依次为39.98%、49.36%、82.75% 和101.75%，沙粒含量减幅依次为16.73%、23.96%、32.48% 和39.07%；在入渗深度为2 ~ 3 cm处土壤颗粒组成与原土壤差异不明显，说明泥沙细颗粒通过入渗进入土体后，其主要滞留分布在0 ~ 2 cm深度的土体中，2 cm以下入渗深度处土壤受浑水泥沙细颗粒影响较小。滞留层土壤颗粒分析结果详见表 6，可以看出，不同浑水含沙率的滞留层土壤物理性黏粒含量$ {d_{0.01}} $显著高于原土壤；以原土壤为对照，浑水含沙率3%、6%、9%、12% 的滞留层0 ~ 1 cm入渗深度处对应的物理性黏粒含量$ {d_{0.01}} $增幅依次为16.49%、27.03%、38.09% 和47.44%，1 ~ 2 cm入渗土壤深度处增幅依次为5.34%、10.23%、18.39% 和27.17%，表明浑水中泥沙含量越大，物理性黏粒伴随着入渗水被洗入到土体的越多，因此其滞留在土壤中形成滞留层的相对含量越多。

3 结论

1) 浑水累积入渗量和入渗率随含沙率的增加而减小，不同含沙率下的累积入渗量与入渗历时的关系均符合Kostiakov模型，其中入渗系数与含沙率呈幂函数负相关关系，入渗指数与含沙率呈幂函数正相关关系。

2) 湿润锋运移距离随浑水含沙率的增加而减小，不同含沙率下的湿润锋运移距离与入渗历时之间符合幂函数关系，其中拟合系数和拟合指数与含沙率均呈指数函数负相关关系。

3) 不同浑水含沙率对应的落淤层厚度均随入渗历时的增大而不断增大，两者符合幂函数关系，拟合系数和拟合参数均随着含沙率增大而不断增大，两者与含沙率之间均呈线性正相关关系；相同入渗历时，落淤层厚度随浑水含沙率的增加而增大。

4) 随着浑水含沙率的增大，滞留现象越明显，落淤层细颗粒相对含量越少，物理性黏粒含量均显著低于浑水泥沙，粗颗粒相对含量越多；滞留层细颗粒相对含量随着含沙率的增加而增多，其物理性黏粒含量显著高于原土壤，特别是在入渗深度为0 ~ 1 cm处。