棚内微咸水覆膜滴灌下带间土壤水盐分布特征模拟

陆培榕 邢玮麟 杨玉杰 刘文龙 罗 纨

(扬州大学水利科学与工程学院, 扬州 225009)

0 引言

覆膜滴灌是将具有保水保墒作用的地表覆膜与局部控制性供水的滴灌方式相结合的高效节水技术[1]。由于覆膜限制了土壤与外部大气的接触,降低了灌后土壤的水分蒸发及溶质上溯,使得即便在湿润淋洗程度有限的滴灌模式下,微咸水也能作为有效的灌水资源[2-3]。然而,滴灌属于典型的点源入渗,灌水受土壤基质势的影响将以湿润体的形式分布于土壤内,在蒸发过程中,微咸水输入的盐分将滞留于湿润体外侧,并随灌水次数增加而逐渐积聚[4],导致土壤中水分和盐分的动态分布过程存在明显的不一致性,当涉及地表覆膜时,土表蒸发与非蒸发面交替存在,水盐动态运移的复杂程度将进一步加剧[5],尤其在缺乏雨水和漫灌淋洗的设施大棚内,盐分的积聚将存在局部叠加,对作物生长造成胁迫[6]。

现阶段,随着农业种植模式的多样化以及设施大棚内滴灌布设模式的集约化,针对覆膜滴灌条件下的水盐运移研究由单点源入渗扩展至多点源交汇。研究的对象以多点源入渗形成的交汇型湿润体几何形态、收缩及扩张速率、脱盐及积盐区的划分为主。虎胆·吐马尔白等[7]以室内砂壤土土槽试验为基础,研究了滴头间距为30 cm时不同滴头流量和灌水量下湿润体交汇区的盐分分布特性,发现滴头下方的脱盐深度大于交汇区的脱盐深度,且脱盐深度的差异随滴头流量的减小而增大。王卫华等[8]以采用膜下滴灌的粉壤土棉田为研究区域,发现膜下湿润体交汇区的含水率随滴头流量的增大而增大,且两滴头间湿润体交汇区径向宽度与灌水历时呈幂函数关系。齐智娟等[9]在轻度盐渍化的玉米田中对比了滴头间全膜覆盖和半膜覆盖的处理方式,发现使用微咸水滴灌时全膜覆盖处理中两滴头间0～30 cm的土层内形成了明显的脱盐带,而在半膜覆盖的条件下土壤脱盐区仅存在于膜下局部区域。为进一步可视化、定量化地描述土壤水盐的分布特征,数值模拟方法伴随计算机技术的发展得以在滴灌湿润体交汇下土壤水盐动态变化的研究中发挥作用,并以HYDRUS数值模型为主取得了一系列成果。例如,王维娟等[10]运用HYDRUS模型,对不同间距下的双点源滴灌湿润体交汇轮廓进行了三维状态下的动态模拟,结果显示,交汇区湿润体形状从初期的1个半球体逐步演变为2个分离的近似半球体。CHEN等[11]模拟研究了间作模式下采用膜下滴灌时土壤水盐的运移特征,分别对不同作物的膜下根系区以及膜间裸地进行子区划分,分析了各子区内的水分及溶质的动态形式以及各子区间的通量交互,发现瞬时的水平向通量交互多发生于灌水阶段,而缓慢的垂向通量交互主要存在于灌后蒸发阶段。GUO等[12]模拟了全膜覆盖形式下湿润体存在交汇时棉田内土壤水盐的分布变化,结果表明,土壤盐分的垂向分布差异随灌水次数的增加而逐步显现,而改变滴头流量造成的影响在不同生育期内表现不一。上述模拟研究表明,覆膜滴灌形成的湿润体交汇下水盐的分布变化具有瞬时性和局部性,在灌后蒸发过程中的变化则较为缓慢且一致。然而,个别的灌水蒸发过程对土壤水盐分布的改变不足以对长期的作物生长造成显著影响。因此,进行针对湿润体交汇区域内水分总体变化趋势以及盐分累积分布特征的研究将更具实用意义。

本文以设施大棚内田间小区试验为基础,结合HYDRUS模型对微咸水膜下滴灌时土壤二维剖面内的水盐动态进行数值模拟。同时,为充分考虑滴头间湿润体交汇程度在灌水和蒸发过程中的空间多变性,采用剖面子区分析的方式,将存在滴头的滴灌带带间区域作为主要的分析对象,探究在不同覆膜宽度及滴头流量下带间剖面水盐分布的动态变化及累积特征,以期为设施环境微咸水覆膜滴灌模式下创造适宜作物根系生长的低盐环境提供理论支撑与方案参考。

1 材料与方法

1.1 试验区概况与试验设计

试验于江苏省扬州市沙头镇现代农业产业园(32°17′N,119°29′E)内的设施大棚中进行。棚内0～100 cm土层的质地为砂壤土,平均砂粒(粒径0.05～2 mm)占比为55.36%,粉粒(粒径0.002～0.05 mm)占比31.78%,粘粒(粒径0～0.002 mm)占比为12.86%;平均土壤容重为1.41 g/cm3,饱和导水率为47.28 cm/d,饱和含水率和田间持水率分别为0.42 cm3/cm3和0.25 cm3/cm3。试验期间棚内地下水埋深的波动范围为1.63～1.97 m。该设施大棚长期抽取地下水(平均电导率1.37 dS/m)进行滴灌,且缺乏雨水或漫灌的淋洗,导致棚内土壤呈轻度盐渍化,试验前0～100 cm土层的平均含盐量为1.52 g/kg。

以田间小区试验为研究基础,在棚内划分 2 m×2 m的地块,采用“两膜两行”滴灌布置模式(图1),膜宽20 cm,滴头间距为35 cm,两条滴灌带间距 80 cm,由外接马氏瓶对小区内的滴灌系统进行恒压供水,并采用压力补偿式灌水滴头将流量统一控制为1.6 L/h。同时,由于取水井内地下水含盐量存在波动,故使用产业园内提供的自来水(平均电导率0.22～0.30 dS/m)与氯化钠晶体(分析纯)配比质量浓度固定为3 g/L的盐溶液进行灌水。各小区中均嵌入张力计,用于监测地表以下深度15 cm处的土壤水势,试验设计两组灌水处理T20和T40,分别对应上一轮灌水后土壤基质势降至-20 kPa和-40 kPa时进行下一轮滴灌的灌水制度,两组处理每轮的灌水持续时间均为2 h。此外,每组处理重复3次,在棚内共计划分6块滴灌小区。试验于2021年7月中旬开始至10月下旬结束,总时长为90 d,为避免作物根系随机性的扩张对土壤水盐动态分布的影响,试验期内未涉及作物的种植。试验过程中,设施大棚仅开启两侧部分的通风口,无降雨或其他形式的灌水。

图1 试验小区覆膜滴灌布设形式示意图

1.2 观测指标与换算方法

在每轮灌水前选用内径为2 cm的土钻对与滴头水平距离10、20、30、40 cm处,深度10、20、30、40 cm处的点位进行取样,用于测定含水率及电导率。由于小区中各滴头对应的流量及灌水量相同,故每轮取样时选取不重复的滴头对相应的位置进行取样。每次取样完毕后,使用试验小区内的土壤对取土孔进行回填,防止大孔隙优先流对灌水及蒸发过程中土壤水分及溶质运移的影响。采用干燥法测定不同取样点位的土壤质量含水率,并根据土壤容重换算成相应的体积含水率。部分取回的土样经风干、碾磨、过筛(孔径1 mm)后,与蒸馏水混合(土水质量比1∶5),制备土壤浸提液。采用电导率仪(DDBJ-350型,上海雷磁创益仪有限公司)室温(20℃)下测定浸提液的电导率,并与蒸发结晶法进行比对,得出土壤含盐量与电导率的关系式为

S=4.12EC1∶5+0.23

(1)

式中S——土壤含盐量,g/kg

EC1∶5——土水质量比为1∶5的土壤浸提液电导率,dS/m

土壤潜在蒸发量通过小型蒸发皿(DF-AM3型,内径20 cm,北京东方鑫鸿科技有限公司)的蒸发量进行估测。根据已有研究理论,当土壤含水率较高(大于70%田间持水率)时,土壤处于非限制性蒸发状态,相应的蒸发量与大气蒸发力成正比,可用水面蒸发量乘以相应的系数进行估测[13-14]。然而,在实际田块中,维持大体积土壤的高含水率,并持续地监测蒸发造成的水分耗散量,实现起来较为困难。故采用100 mL圆柱形环刀进行土壤取样(20组),在环刀内饱和后底部加盖,埋置于田块中,并保证顶部与表土齐平。每隔3 h对环刀进行称量,并换算成相应的蒸发量,直至环刀内土壤含水率小于田间持水率时停止观测。随后,将观测时段内土样的蒸发量与棚内的小型蒸发皿对应的蒸发量建立拟合关系式

Es=0.63Epan

(2)

式中Es——棚内土壤潜在蒸发量,mm

Epan——小型蒸发皿蒸发量,mm

试验期内,棚内蒸发皿的日蒸发量及T20和T40两组处理的灌水量如图2所示。

图2 试验期棚内蒸发皿蒸发量及T20和T40处理灌水量

1.3 HYDRUS-2D模型构建

1.3.1模型原理及基本方程

HYDRUS-2D是以有限元计算为基础的数值模型,可用于模拟二维条件下,变饱和介质中水分、热量及溶质的动态分布状况及边界处的瞬时及累积通量[15]。由于滴灌对应的点源入渗将在土壤中形成趋于轴对称的半椭球状湿润体[16],故滴灌对应的含水率分布模拟可简化至二维垂向剖面内。基于此,以土壤均匀且各项同性为基本假设,并忽略土壤水分变化的滞后效应,模型采用修正后的Richards方程[17]描述二维模式下的土壤水分运动,即

(3)

式中θ——土壤体积含水率,cm3/cm3

h——压力水头,cm

x——横向坐标,cm

z——纵向坐标,向上为正,cm

t——时间,h

K(h)——当土壤内的压力水头为h时的土壤导水率,cm/h

S(h)——根系吸水项,h-1,本研究没有涉及作物栽培,故该项忽略

同时,模型采用van Genuchten模型[18]来描述土壤含水率、基质势以及导水率的相互关系,即

(4)

(5)

(6)

(7)

式中θs——土壤饱和含水率,cm3/cm3

θr——土壤残余含水率,cm3/cm3

α——进气吸力的倒数,cm-1

Ks——土壤饱和导水率,cm/h

Se——土壤相对饱和度

l——孔隙关联度系数,取0.5[19]

m、n——土壤水分特征曲线形状系数

本研究中,灌溉用水均为氯化钠溶液,未涉及到具有反应性溶质的化肥或农药的输入。因此,在模拟过程中,将土壤中的溶质统一简化为保守性溶质,模型选用的保守性溶质条件下描述溶质运移的对流弥散方程的二维形式为[5]

(8)

式中C——土壤溶液质量浓度,g/L

qx、qz——横向、纵向通量,cm/h

DT、DL——横向、纵向弥散系数,cm

1.3.2边界及初始条件设置

以灌水滴头正下方0～100 cm的土壤剖面为对象建立二维模拟域。如图3所示,模拟域水平长度与试验小区一致设为200 cm,垂向深度为 100 cm。模拟域由约30 000个三角形有限单元组成,对存在蒸发和灌水的上边界加密节点的分布,并按照试验取样点位在模拟域上设置观测点位。模拟域的上边界可分为覆膜区、无膜裸地以及滴头进水区。模拟时假设地膜完全隔绝土壤表面与外部空气的接触,故将覆膜区设置为零通量边界。无膜裸地设置为大气边界,由于棚内试验无降雨及作物蒸腾,因此,模拟域的大气边界上仅输入估测所得的土壤潜在蒸发量Es。滴头进水区的范围取决于滴头的流量及土壤的水分特征,已有研究表明,滴灌时较快的出流会在土壤表面形成趋于圆形的稳定进水区域,采用该进水区域作为入渗边界,能较为精确地描述二维形式下滴灌形成的土壤水分动态过程[20-21],相应的进水区域半径可近似表示为

图3 HYDRUS-2D构建的试验小区土壤剖面模拟域示意图

(9)

式中Rs——进水区域半径,cm

q——滴头流量,cm3/h

进水区设置为变通量边界,长度为2Rs,非灌水时通量为零,灌水时通量由滴头流量和相应的通量边界尺寸确定[11]。试验过程中,棚内各小区相对独立,除滴灌外无其他形式的水分输入,故模拟域的左右两侧设置为零通量边界。此外,在试验过程中,地下水的观测时间间隔较长,而模拟中以“小时”为时间步长,使得边界的赋值难以动态匹配。结合现有的研究结论,当地下水埋深超过1 m时,裸地相较于作物耕地的潜水上溯蒸发量将显著下降[22-23]。本研究未考虑作物生长,且观测深度集中在0～40 cm,故模拟过程中忽略地下水波动对近地表土壤水盐分布的影响,将模拟域底边设置为自由出流边界。

模拟域剖面的初始水分及盐分根据试验前的实测值分层对应设置。在输入含盐量时,土壤含盐量转换为HYDRUS-2D模型采用的溶质液相浓度[24],关系式为

(10)

式中Sc——土壤体积含水率为θ时土壤的盐分质量浓度,g/L

ρ——土壤干容重,g/cm3

初始的土壤水分参数除实测所得的饱和含水率θs及饱和导水率Ks,其余参数均根据设施大棚内土壤的颗粒级配及容重通过HYDRUS内置的Rosetta软件进行推演。溶质运移相关的纵向弥散系数DL及横向弥散系数DT与模拟域对应的空间尺寸相关。参照有关HYDRUS-2D对土壤剖面水盐运移的模拟研究[25-26],将DL的初始值设为模拟域深度的1/10,且DT=DL/10。本研究以小时(h)为时间步长模拟土壤剖面内水分及溶质在灌水及蒸发过程中的运移状况,模拟时长为90 d,共计2 160 h。

1.3.3模型检验

通过对比不同取样点处水盐含量的实测值与模拟值的差异来分析所建模型的模拟精度,进而对相关土壤水分及溶质运移参数进行率定和验证,并采用均方根误差(Root mean square error,RMSE)、纳什效率系数(Nash-sutcliffe efficiency coefficient,NSE)、平均相对误差(Mean relative error,MRE)作为量化模拟值与实测值差异程度的评价指标。RMSE和MRE越接近0,说明模拟值与实测值的差异越低;NSE的波动范围为-∞～1,该值越接近1说明模拟值与实测值的一致性越高,通常当NSE大于0.5时,可认为模拟结果满足精度要求[27]。

2 结果与分析

2.1 模型率定与验证

采用T20处理条件下与滴头水平距离为10、20、30、40 cm处0～40 cm土层内每轮灌水前土壤的体积含水率和电导率EC1∶5对模型参数进行率定,随后,采用T40处理下相应的观测值对率定后的参数进行验证。如图4所示,在率定和验证过程中,T20和T40处理对应的土壤体积含水率及电导率的模拟值和实测值基本集中分布在1∶1线附近。T20和T40处理的含水率模拟值和实测值的拟合线斜率均大于1,而土壤电导率模拟值和实测值的拟合线斜率略小于1,说明两组处理中含水率实测值(T20:0.18～0.27 cm3/cm3;T40:0.19～0.23 cm3/cm3)的波动范围均大于模拟值(T20:0.20～0.25 cm3/cm3;T40:0.15～0.26 cm3/cm3),而土壤电导率实测值(T20:0.17～0.62 dS/m;T40:0.19～0.67 dS/m)的波动范围略高于模拟值(T20:0.16～0.56 dS/m;T40:0.21～0.62 dS/m)。表1中,率定和验证阶段土壤含水率RMSE波动范围为0.01～0.02 cm3/cm3,NSE为0.52～0.63,MRE为4.39%～8.53%;土壤电导率对应的RMSE、NSE和MRE在相同取样位置处较含水率均有所提升,波动范围分别为0.03～0.07 dS/m,0.53～0.82以及8.36%～13.42%。同时,RMSE和MRE随与滴头水平距离的增加而有所提升,且NSE随之减小,说明土壤水分和盐分的模拟精度随与滴头水平距离的增加而逐渐下降。总体而言,RMSE接近于0,MRE较低,NSE均大于0.5,表明在构建的二维模拟域中,所输入的土壤水分及溶质运移参数能够较为准确地模拟微咸水膜下滴灌时土壤水盐的动态变化。经率定及验证,选定的土壤残余含水率θr及饱和含水率θs分别为0.047 cm3/cm3和0.415 cm3/cm3,土壤进气吸力倒数α为0.02 cm-1,土壤水分特征曲线形状系数n为1.5,饱和导水率Ks为1.97 cm/h,孔隙关联度系数l为0.5,纵向和横向弥散系数DL和DT分别为 30 cm和3 cm。

表1 滴头不同水平位置处深度0～40 cm内土壤含水率及电导率模拟精度评价指标

图4 T20处理和T40处理下距滴头不同水平位置处深度0～40 cm内土壤含水率及电导率模拟值与实测值对比

2.2 模型应用

2.2.1土壤剖面内水分及盐分的变化特征

在试验及模拟过程中,T20和T40处理单次灌水水量一致,但总灌水次数不同(图2)。因此,本文以两组处理在第5次灌水后的1、3、6 d的模拟剖面含水率分布为对象,对比分析在带间区域和滴灌带外侧区域(滴头和试验小区边界侧之间,以下简称“带外区”)土壤水分的变化特征。如图5所示,由于T40处理的灌水时间间隔长,剖面上含水率的差异较T20处理更大。在灌水后1 d,两组处理均在滴头下方形成了含水率较高的圆形湿润体,且湿润体的边缘存在交汇,使得膜间裸地下方的含水率相应提升。在灌后第3天和第6天,受土表蒸发的影响,覆膜区和裸地的含水率差异逐渐体现,同时,土壤剖面内的局部含水率较高的湿润体逐渐消失,剖面内水平向的含水率趋于一致。

图5 T20和T40处理灌水后1、3、6 d土壤剖面含水率分布

在土壤含盐量变化方面,本文通过HYDRUS-2D模拟的土壤盐分质量浓度来反映土壤溶质的分布及运移状况,并针对T20和T40处理在模拟过程中第30、60、90天的剖面盐分分布进行分析。由 图6 可知,由于T20处理在相同时段内的淋洗次数多于T40处理,使得T20处理在整个剖面的盐分质量浓度均低于T40处理;两组处理受灌水淋洗作用,盐分质量浓度在10 g/L以下的区域自滴头向下逐渐扩张;在土表裸地处,土壤内溶质随水分蒸发而逐渐表聚,且带间区盐分的表聚程度和范围均低于模拟域两侧的带外区。可见,土壤剖面内盐分的分布状况受滴灌后湿润体交汇及扩散的影响,虽然灌水次数较多的T20处理中剖面输入的盐分更多,但灌水次数较低的T40处理条件下盐分的表聚程度更为明显。

图6 T20和T40处理在第30、60、90天土壤剖面盐分分布

此外,为进一步比对带间区与带外区内沿水平方向的盐分动态变化特征,对T20和T40处理土壤剖面分别距滴头两侧10、20、30、40 cm水平位置处0～40 cm的常规耕作土层的水盐含量变化进行模拟分析。如图7所示,在土壤含水率变化方面,与滴头的水平距离越小,灌水前后土壤含水率的波动范围越大,但带间区和带外区相同水平位置处含水率的差异越小;同时,在不同水平位置处,各处理中每轮灌水前后土壤含水率的变化幅度趋于一致,没有明显的总体上升或下降趋势。此外,距滴头10、20、30 cm水平位置处带外区与带间区在灌水后含水率的峰值基本相同,而在40 cm位置处T20和T40处理带间区的含水率峰值分别比带外区高11.55%～17.78%和13.34%～16.52%,说明在40 cm处两滴头形成的湿润体发生交汇重合,造成含水率上升。

图7 与滴头不同水平位置处0～40 cm土层平均土壤含水率与盐分质量浓度变化曲线

在土壤含盐量变化方面,HYDRUS-2D模型通过土壤盐分质量浓度波动来反映观测点位含盐量的变化。T20的带间区与带外区的含盐量均低于T40处理,且差异随水平距离的增加而逐渐扩大;两组处理内各水平位置处外侧区域的含盐量均高于带间区,且差异随灌水次数的增加而扩大。距滴头10、20、30、40 cm水平位置处T20处理带间区的土壤含盐量总体呈下降趋势,模拟时段末的土壤含盐量均低于初始值。相比之下,在T40处理下带间区内距滴头水平位置30、40 cm处土壤盐分呈总体上升趋势,尤其在距离40 cm处灌水淋洗后土壤含盐量的最低值与初始值几乎一致。上述现象表明,相比于带外区,湿润体的交汇能够抑制土壤盐分的集聚程度,但当灌水频率降低时带间区内膜间裸地的积盐程度随之上升。

2.2.2情景模拟

为进一步探究带间区内土壤水盐变化受滴灌强度和覆膜宽度的影响,本研究以率定及验证后的模型参数和T40处理的灌水制度为基础,利用HYDRUS-2D模型模拟了不同设计情景下带间区深度0～40 cm范围内土壤水盐的变化过程。共设计了6组滴头流量(0.5、1.0、1.5、2.0、2.5、3.0 L/h)递增的情景,且各流量情景对应的单次灌水量(3.2 L)、模拟期内总灌水量(28.8 L)以及灌水含盐量(3 g/L)均与T40处理保持一致。覆膜宽度的情景同样为6组,在带间区一侧的膜宽分别为15、20、25、30、35、40 cm,相应的带间区内膜间裸地间距分别为50、40、30、20、10、0 cm(图8)。

图8 不同覆膜宽度的设计情景

(1)带间区土壤水盐总体变化趋势

图9为设计的36组模拟情景下带间区内深度0～40 cm范围内土壤含水率及盐分质量浓度在 90 d 内的变化过程。图9中,滴头流量对带间区内的含水率影响程度较低,不同覆膜宽度下带间区内的含水率在每轮灌水前后的变化量基本一致,仅有0.02～0.04 cm3/cm3的差异。带间区内土壤含水率随膜间裸地间距的缩小而略有上升,平均含水率最高值(28.76 cm3/cm3)在裸地间距为0 cm的处理中,最低值(25.12 cm3/cm3)在膜间裸地间距50 cm的处理中。相比之下,带间区的含盐量随灌水次数的增加呈逐步下降趋势,且下降程度随膜间裸地间距的减小而增大,当膜间裸地间距由50 cm减至0 cm时,平均土壤盐分质量浓度由9.53 g/L下降至6.25 g/L。此外,不同的滴头流量仅在灌水间歇期内土壤含水率和盐分质量浓度呈轻微差异,所有情景中平均土壤含水率和盐分质量浓度的最大差异仅为0.14 cm3/cm3和0.22 g/L,出现在膜间裸地间距50 cm时流量为0.5 L/h和3.0 L/h的处理之间。上述模拟结果表明,带间区在多轮滴灌后均能得到盐分淋洗,且淋洗程度随地膜间距的增大而增加,但滴头流量对水盐分布的影响相对较弱。

图9 不同设计情景下带间区深度0～40 cm内土壤含水率及盐分质量浓度变化曲线

(2)膜间裸地及滴头流量对带间区土壤盐分分布的影响

以试验结束的时间(第90天)为分析对象,图10(图中不同颜色的实线为不同滴头流量下相同膜间裸地间距对应的平均值,实线周边的阴影为相应的标准偏差)为各模拟情景下土壤盐分质量浓度在深度0～40 cm处带间区盐分的垂向分布形式。由图10可知,在相同的膜间裸地间距内滴头流量的改变对盐分垂向分布造成的差异较小,所有模拟情景中最大偏差出现在深度5 cm处裸地间距为50 cm的处理中,仅为0.046 g/L。此外,除膜间裸地间距为0 cm的处理外,其余情景中土壤的含盐量均随深度增加而降低,随膜间裸地间距的扩大而增大。各类覆膜宽度条件下土壤盐分质量浓度的差异随深度的加深而逐渐缩小,在深度 5 cm 处每增加10 cm的膜间裸地,土壤盐分质量浓度将提升约1.73 g/L,而在深度40 cm时,提升幅度降为1.34 g/L。

图10 不同膜间裸地间距下带间区内土壤深度0～40 cm土壤盐分质量浓度变化曲线

以带间区水平向的中点为原点,分析从中点至两侧滴头处深度0～40 cm土壤盐分质量浓度在各类情景下的分布状况(图11)。由图11可知,不同水平位置处的土壤盐分质量浓度随膜间裸地间距的增加而降低;除膜间裸地间距为0 cm处理外,土壤盐分质量浓度在不同的覆膜宽度处理下的最大值均出现在中点位置,并向两侧滴头处逐渐降低。此外,不同滴头流量下土壤盐分的分布趋势基本一致,相应的差异主要存在于中点处,且盐分浓度随滴头流量的增加而增加。例如,在膜间裸地间距为50 cm处理中,中点处的盐分质量浓度从流量3.0 h/L对应的12.82 g/L逐渐降至流量 0.5 h/L 对应的12.15 g/L。随着膜间裸地间距的降低,含盐量将在靠近两侧滴头处呈现抬升,并在距两侧滴头 10～20 cm 处回落。该现象主要是由带外区的覆膜宽度较窄(图8),集聚的盐分逐步扩散,并推进至带间区所致。

图11 滴头流量和膜间裸地间距对带间区内不同水平位置处深度0～40 cm土壤盐分质量浓度的影响

3 讨论

采用数值模型的方式对土壤水盐动态进行模拟是当前分析及评估滴灌实施效益的重要方式[28],而所建模型的模拟精度保障是决定模拟研究充分可行的前提条件。在本研究的模型率定与验证过程中,含水率模拟值的变化量小于实测值的变化量,原因是模型假设土壤的剖面为均质且各向同性,而实地条件下,由于取样或土壤干缩裂缝等因素容易在土体局部形成优先流[29],导致土壤中模拟与实测的含水率差异扩大。同时,浅层地下水波动会引起上层土壤水分分布的变化,而在模拟域构建时将下边界简化为自由出流边界,这将导致模拟时水分的下渗速率过高,影响了对水分分布的模拟效果。模拟所得的土壤含盐量略高于实际观测值,原因可能是在模拟域构建过程中,含盐量假设为自上而下均匀线性分布的状态,与实际取样所得的观测值之间必然存在差异,进而影响后续的模拟精度。另外,在膜下滴灌时,土壤的盐分迁移与水分运动并不同步[8],灌水和蒸发造成的湿润体交汇与收缩,加剧了土壤水盐运移及分布的复杂性,使得越接近带间区中间部位的模拟结果误差越大(表1)。考虑到相应的误差分析指标(RMSE、NSE、MRE)均处于精度可接受的范围内,表明验证后的模型能够有效描述微咸水膜下滴灌过程中的水盐分布状况。

现有的试验及数值模拟研究均证实,滴灌过程中湿润体交汇程度越高,对应的带间区内土壤含水率的提升程度越大[10,30-31]。本研究中,不同覆膜宽度下,每轮灌水过程中带间区内含水率的抬升程度趋于一致,说明覆膜宽度的变化对湿润体交汇程度的影响较小。原因可能是灌水历时较短,期间的水分蒸发量不足以体现地膜抑制蒸发的效果。在盐分分布方面,SELIM等[4]以两滴头间距为40 cm的壤砂土剖面为对象,同样利用HYDRUS模型模拟了微咸水(0～2 dS/m)滴灌下水盐的分布状况,发现以滴头流量1.0 L/h灌水后,将在湿润体交汇区以下10～40 cm范围内形成了超过初始含盐量2～4倍的锥形积盐区。然而,在本研究中,滴灌形成的湿润体交汇区下方并未出现明显的积盐区(图6),主要原因可能是设置的土壤初始含水率较低,滴灌使用的微咸水浓度低于土壤的盐分质量浓度;另一方面,在灌水后的蒸发阶段,地表覆膜抑制了土壤水分的散失[32],限制了溶质的上溯,减缓了湿润体外围积盐区的扩张。在模拟末期,由于带外区覆膜宽度较窄,外侧裸地集聚的盐分向膜内扩散,使得部分模拟情景中,含盐量的最低值出现在距滴头10～20 cm的范围内。该现象说明覆膜虽能够抑制蒸发导致的垂向盐分上溯,但并不能限制盐分在土壤剖面内的侧向扩散,尤其在淋洗强度和范围有限的滴灌模式下,盐分的空间分布由早期的“膜外表聚型”逐渐向后期的“膜内底聚型”转变[33]。

灌水量一致的情况下,不同滴头流量下带间区内的总体含水率变化幅度的差异较低,且水分的差异将随膜间裸地面积的降低而进一步缩小;土壤盐分的差异则主要集中在膜间裸地内,且积聚程度随滴头流量的增大而增加(图11),原因是滴头流量越大形成的湿润体深度范围越小,但水平向范围越大[34],使得滴头流量较大的处理中土壤盐分的淋洗深度有限。同时,在点源入渗过程中,土壤中湿润体的扩张速率和含水率随与滴头距离的增加而降低[35],即便当两滴头形成的湿润体存在交汇时土壤水平向的淋洗程度也相对有限,导致灌水间歇期内滴头流量较高处理的蒸发返盐程度高于滴头流量较低的处理。另外,本研究的灌水制度是根据设计的土壤含水率阈值而制定的,由于研究中没有考虑根系吸水的情况,土壤达到含水率阈值主要依靠灌后的表层蒸发,这将导致灌水间歇期相对延长,加剧了蒸发造成的土壤盐分集聚,且模拟过程中土壤的潜在蒸发是根据小型蒸发皿估算所得,造成与实际耕种时土壤水盐环境的偏差。同时,带间区内垂向及水平向的水分及溶质通量交换有待进一步细化,并加强对浅层地下水波动形式的高频动态监测,从而提升模型的实用性以及棚内微咸水覆膜滴灌方案设计的合理性。

4 结论

(1)以HYDRUS-2D软件为基础构建了“两膜两行”的微咸水覆膜滴灌二维土壤剖面有限元模拟域,并引用设施大棚内的蒸发状况作为大气边界条件模拟棚内环境。率定及验证的结果表明,构建的土壤水盐模型精度可靠,在滴灌带带间区内土壤含水率及含盐量的均方根误差、纳什效率系数以及平均相对误差均在合理的范围内波动,模拟精度随与滴头水平距离的减小而增加,且含水率的模拟精度较含盐量高。

(2)针对不同膜间裸地间距及滴头流量展开的情景模拟结果表明,随灌水次数的增加,带间区域的土壤含盐量总体呈下降趋势,裸地间距越小下降程度越大;当灌水量不变时,不同滴头流量下带间区内土壤水盐含量的差异较低,且膜间裸地间距越小差异将进一步减小。

(3)多轮微咸水滴灌后,膜间裸地间距越大土壤的近地表积盐现象越突出,不同滴头流量造成的土壤盐分含量差异在膜间裸地内较为明显,且提升滴头流量将加剧积盐程度;在膜间裸地距离为0～30 cm处理中,带外区积聚程度较高的土壤盐分将向含盐量较低的带间区扩散,使得水平向的含盐量最低值并非出现在淋洗程度最大的滴头处,而存在于与滴头水平距离10～20 cm范围内。