竖直微润灌土壤水分运动数值模拟与验证

范严伟,赵 彤,赵廷红

(兰州理工大学能源与动力工程学院,甘肃 兰州 730050)

微润灌是利用高分子半透膜制成微润管,将半透膜技术引入灌溉领域,通过微润管渗出水分进行地下灌溉的一种低定额局部灌水技术[1-3]。微润灌技术结构简单，动力消耗少，运行费用低，对温室、大田作物具有良好的节水增产效果。近年来，在山丘、盐碱荒地和沙漠地带等复杂环境均得到应用。

微润管铺设方式主要分为水平铺设和竖直铺设两种,水平铺设为微润管沿行水平铺埋安装,适宜于密植、浅根系植株灌溉,如蔬菜、玉米等作物;竖直铺设为微润管在根系附近竖直插入,适宜于宽距、深根系植株灌溉,如果树、林木等植物。目前,对微润灌土壤水分运动规律的研究主要针对水平铺设情况,且以室内试验为主[4-8]。采用数值模拟方法对微润灌土壤水分运动规律的研究相对较少[9]。随着计算机技术和数值分析方法的发展,数值模拟方法越来越多地应用于农田灌水技术的研究中[10-14]。数值模拟方法具有简便、快捷、灵活等特点,可在给定的初始、边界条件下,模拟不同影响因素对土壤水分运动规律的影响,为确定适宜的灌水参数和设计灌溉系统提供了实用而方便的手段[15]。本文以非饱和土壤水动力学理论为基础,针对微润管竖直铺设特点,建立微润灌土壤水分运动数学模型,用SWMS-2D软件求解,通过室内试验对模拟结果进行分析验证,以期借助数值模拟方法,进一步认识竖直微润灌土壤水分运动机理,为后期模拟分析土壤特性和灌水技术要素对微润灌土壤水分运动规律的影响提供技术支撑。

1 试验材料与方法

1.1 供试土壤

供试土壤取自甘肃省白银市景泰县柏林山庄的风沙土和兰州市黄峪镇祁家营的粉壤土,取土深度为0～40 cm。将取得的土壤风干、碾压,均匀混合,过2 mm筛后制成试验土样。室内测定土壤基本物理特性,结果如表1所示。

表1 试验土壤基本物理特性

1.2 试验装置与方法

试验装置由土箱、微润管、马氏瓶和橡胶软管四部分组成,如图1所示。土箱采用10 mm厚有机玻璃制成,长60 cm、宽60 cm、高100 cm,为防止水分入渗过程中产生气阻,土箱底部设置通气孔,为便于灌水结束后取土,土箱壁每隔5 cm开直径为2 cm的取土孔。将供试土壤按设定密度(风沙土1.56 g/cm3,粉壤土1.33 g/cm3)分层均匀装入土箱。为便于观察湿润锋形状和运移情况,根据竖直线源入渗的对称性,将微润管紧靠土箱壁放置,微润管埋深设置为20 cm。采用深圳市微润灌溉技术有限公司的白色微润管,其管径为2 cm,孔径为10～900 nm。风沙土中管长为40 cm,粉壤土中管长为30 cm。试验过程中,马氏瓶为微润管提供恒压水头,马氏瓶直径为80 cm,高度为100 cm。

图1 试验装置示意图

试验开始后,定时记录马氏瓶读数和绘制湿润锋轮廓线,入渗达70 h后,停止供水,从取土孔取土测定土壤体积含水率。为尽量消除试验误差,每种土壤设置2个重复,试验结果取两者平均值。

2 数学模型

2.1 基本方程

竖直微润灌土壤水分运动属于三维问题,若假定土壤为各向同性的均质体,根据对称性,则可视为轴对称条件下的二维问题。其土壤水分运动基本方程为

(1)

式中:x为横向坐标;z为垂向坐标,规定z向上为正;θ为土壤含水率,cm3/cm3;φ为基质势,cm;t为入渗时间,min;K(θ)为非饱和导水率,cm/min。

式(1)涉及θ、φ和K(θ)三者之间的关系,模拟中采用van Genuchten-Mualem模型拟合[16-17]，即

(2)

(3)

由于微润管内外管径较小,为简化研究工作,微润管按线源处理[13]。考虑到竖直微润灌土壤水分运动的轴对称特性,选取图2所示模拟区域(EFGH为求解区域,AB为微润管,C点为微润管中点)。

图2 求解区域示意图

2.2 定解条件

a. 初始条件：

(4)

式中:θ0(x,z)为土壤初始含水率,cm3/cm3；xF和zH为模拟区域边界(装置物理边界)在x和z方向的坐标。

b. 边界条件。不考虑土壤蒸发和降雨,上边界EF按零通量面处理:

(5)

考虑地下水埋深较大,灌溉水量未到达,下边界GH按零通量面处理:

(6)

左边界EA和BH为微润管中心入渗面,由于微润管为轴对称,可按零通量面处理:

(7)

式中:zA和zB为微润管边界在z方向的坐标。

微润管通过管内水压调节出流量,室内试验结果表明,灌溉过程中,微润管渗水速率基本恒定[4-9]。故左边界AB可采用第二类边界条件处理:

(8)

式中:q0为微润管渗水速率,cm/min;Q(t)为单位长度微润管入渗速率,mL/(cm·min);D为微润管管径,cm。

由于对称性,右边界FG为交汇界面,微润管之间关于交汇界面对称,可按零通量面处理:

(9)

2.3 模型求解

利用二维有限元土壤水分运动模拟软件SWMS-2D[18]进行数值求解。将模拟计算区域剖分为长方形单元,考虑到田间实际和计算精度要求,有限元计算区域的深度为100 cm,宽度为30 cm,深度和宽度间隔为1 cm。数值模拟中供试土壤的van Genuchten-Mualem模型参数见表2。

表2 不同土质水力特性参数

2.4 微润管出流特性

根据试验数据,计算不同时刻微润管单位长度入渗量,得到单位长度入渗量与时间的关系曲线如图3所示。

图3 单位长度入渗量与时间的关系曲线

由图3可看出,微润管单位长度入渗量与时间具有良好的线性关系(R2=0.999 5、0.998 9)。对曲线进行线性回归分析,可拟合得出风沙土和粉壤土的单位长度入渗速率分别为0.025 6 mL/(cm·min)和0.029 7 mL/(cm·min)。

3 数值模拟与验证

用SWMS-2D软件求解竖直微润灌土壤水分运动数学方程,获得微润管单位长度入渗量、土壤湿润锋运移值、土壤剖面含水率等指标,将模拟结果与试验结果进行对比验证。

3.1 微润管单位长度入渗量

求解过程中,单位长度微润管渗水速率采用实测资料拟合值,模拟获得不同时刻单位长度入渗量,与实测值进行比较,结果如表3所示。

由表3可知,入渗初期(9 h内),单位长度入渗量的模拟值与实测值差别较大,相对误差绝对值最大为16.66%,平均为11.66%,主要是由于入渗初期,微润管壁土壤含水率急剧增加,土壤基质势对微润管影响明显;入渗后期(9 h后),单位长度入渗量的模拟值与实测值差别较小,相对误差绝对值最大为2.54%,平均为1.19%,主要是由于入渗后期,微润管壁土壤含水率基本恒定,土壤基质势对微润管影响微弱。由于微润灌溉为线源续灌入渗方式,其入渗后期的稳定性更有利于模拟结果的准确性。对模拟值与实测值进行统计性分析可得,风沙土模拟值与实测值决定系数R2=0.999,均方根误差(RMSE)为1.130 mL/cm,t检验参数(Sig)为0.09>0.05;粉壤土模拟值与实测值决定系数R2=0.999,均方根误差为0.935 mL/cm,t检验参数为0.19>0.05。说明模拟值与实测值无显著性差异,模拟过程中,采用实测单位长度渗水速率进行数值求解是可行的。

表3 微润管单位长度入渗量模拟值与实测值对比

3.2 土壤湿润锋运移值

模拟获得土壤湿润锋随时间的动态变化,如图4所示,图中实线代表实测值,虚线代表模拟值。

图4 湿润锋动态变化对比

由图4可看出,不同观测历时的土壤湿润体形状均为围绕微润管的椭球形。入渗初期,各方向运移距离从大到小的顺序为水平向湿润锋运移距离、垂直向下运移距离、垂直向上运移距离,主要是由于竖直放置微润管,土壤水分运动受基质势和重力势影响,入渗初期,土壤水分运动主要受基质势影响,重力势影响较小,而水平向基质势水力梯度大于垂直向基质势水力梯度,导致水平向湿润锋运移较快。随着灌水时间的延长,湿润体各个方向上的湿润锋不断扩展,但扩展速度逐渐减慢,湿润锋垂直向下运移距离有超过水平向湿润锋运移距离的趋势,如风沙土入渗达到70 h时,主要是入渗后期,风沙土重力势影响逐渐增大,基质势影响逐渐减弱。

微润管竖直埋设时,合理的微润管埋深可降低表层土壤含水量,减少土壤水分的无效蒸发,提高土壤水分利用率,此外,竖直微润灌土壤湿润锋运移规律特点,利于根系的伸展和深扎,显著提高深根性植物吸水能力,符合固沙植物高效灌溉要求。

定量分析A点垂直向上及水平方向、B点垂直向下及水平方向和C点水平方向湿润锋运移值,并与实测值比较,如表4所示。

由表4可知,各时段的模拟值与实测值基本一致,相对误差绝对值最大为15.38%,平均为3.80%,误差较大点可能是装土不均匀所致。对模拟值与实测值进行统计分析可得,模拟值与实测值决定系数R2=0.997,均方根误差为0.291 cm,t检验参数为0.77>0.05,说明模拟值与实测值无显著性差异,所建模型是可靠的,SWMS-2D软件能较准确地模拟获得竖直微润灌土壤湿润锋运移规律。

3.3 土壤剖面含水率

图5显示了灌水结束时(70 h)土壤剖面含水率的模拟值和实测值(图中曲线为模拟值,标记点为实测值,土壤含水率均为体积含水率)。

由图5可以看出,土壤质地越黏重，湿润体体积越小。土壤含水率等值线围绕微润管呈近似椭球形分布,土壤含水率从微润管附近向四周逐渐降低,微润管附近土壤含水率未达到饱和状态。主要是由于微润灌灌水流量极小,灌水时段内水分主要依靠土壤吸力迅速向四周扩散,微润管周围未产生积水,不存在入渗水头。

分析图5中灌水结束后土壤剖面含水率的模拟值分布规律,可知,除微润管附近误差较大外,其他处误差均较小。模拟值与实测值之间误差小于5%的占66%,误差在5%～10%之间的占18%,误差在10%～20%之间的占13%,误差在20%以上的占3%。误差较大点可能是土体不够均匀、取土量较少等原因造成的。对模拟值与实测值进行统计分析可得,模拟值与实测值决定系数R2=0.993,均方根误差为0.010 cm3/cm3,t检验参数为0.10>0.05,说明模拟值与实测值无显著性差异,表明所建模型是正确的,SWMS-2D软件能较好地模拟竖直微润灌土壤剖面含水率的分布规律。

4 结 语

为分析竖直微润灌土壤水分的运动状况,依据非饱和土壤水分运动理论,借助计算机数值模拟方法,应用SWMS-2D软件对竖直微润灌土壤水分运动进行数值模拟。采用微润管单位长度入渗量、土壤湿润锋运移值及土壤剖面含水率等指标的实测值与模拟值对模型进行了分析验证,结果表明,数值计算结果与实测结果具有良好的一致性,所建立的数学模型适用于竖直微润灌土壤水分运动状况。

表4 湿润锋运移值模拟值与实测值对比

图5 土壤剖面含水率数值模拟值与实测值对比

试验和模拟均表明,不同观测历时的土壤湿润体形状及含水率等值线均为围绕微润管的椭球形;土壤湿润锋运移规律主要表现为水平向运移距离>垂直向下运移距离>垂直向上运移距离,该规律利于根系的伸展和深扎,符合固沙植物高效灌溉要求。

文中所建立数学模型及求解方法为进一步研究不同土壤特性(土壤质地、密度、含水率)和灌水要素(埋深、压力水头、流量)等条件下的竖直微润灌土壤水分运动规律提供理论依据和技术手段。

[1] PETTY J D,HUCKINS J N,MARTIN D B,et al.Use of semipermeable membrane devices (SPMDS) to determine bioavailable organochlorine pesticide residues in streams receiving irrigation drainwater[J].Chemosphere,1995,30(10):1891-1903.

[2] 杨文君,田磊,杜太生,等.半透膜节水灌溉技术的研究进展[J].水资源与水工程学报,2008,19(6):60-63.(YANG Wenjun,TIAN Lei,DU Taisheng,et al.Research prospect of the water-saving irrigation by semi-permeable film[J].Journal of Water Resources & Water Engineering,2008,19(6): 60-63.(in Chinese))

[3] YANG Q.Irrigation method with semi permeable membrane,irrigation container and irrigation system made of semi permeable membrane and their application: EP,EP 2153714 A1[P].2010.

[4] 张俊,牛文全,张琳琳,等.微润灌溉线源入渗湿润体特性试验研究[J].中国水土保持科学,2012,10(6): 32-38.(ZHANG Jun,NIU Wenquan,ZHANG Linlin,et al.Experimental study on characters of wetted soil in moistube irrigation[J].Science of Soil and Water Conservation,2012,10(6): 32-38.(in Chinese))

[5] 牛文全,张俊,张琳琳,等.埋深与压力对微润灌湿润体水分运移的影响[J].农业机械学报,2013,44(12):128-134.(NIU Wenquan,ZHANG Jun,ZHANG Linlin,et al.Effects of buried depth and pressure head on water movement of wetted soil during moistube-irrigation[J].Transactions of the Chinese Society for Agricultural Machinery,2013,44(12):128-134.(in Chinese))

[6] 张俊,牛文全,张琳琳,等.初始含水率对微润灌溉线源入渗特征的影响[J].排灌机械工程学报,2014,32(1): 72-79.(ZHANG Jun,NIU Wenquan,ZHANG Linlin,et al.Effects of soil initial water content on line-source infiltration characteristic in moistube irrigation[J].Journal of Drainage & Irrigation Machinery Engineering,2014,32(1): 72-79.(in Chinese))

[7] 吴卫熊,何令祖,张廷强,等.滴灌和微润灌条件下桂西北山区典型土壤水分运移规律分析[J].节水灌溉,2016(9): 18-25.(WU Weixiong,HE Lingzu,ZHANG Tingqiang,et al.Research on soil-water movement characteristics of Guangxi main soil under the condition of moisture irrigation and drip irrigation [J].Water Saving Irrigation,2016 (9): 18-25.(in Chinese))

[8] 祁世磊,谢香文,邱秀云,等.低压微润带出流与入渗试验研究[J].灌溉排水学报,2013,32(2): 90-93.(QI Shilei,XIE Xiangwen,QIU Xiuyun,et al.Research on flow and infiltration of low-pressure moistube[J].Journal of Irrigation & Drainage,2013,32(2): 90-93.(in Chinese))

[9] 陈高听,郭凤台,王利书,等.基于HYDRUS的微润灌溉线源入渗数值模拟[J].人民黄河,2016,38(4): 144-148.(CHEN Gaoting,GUO FengTai,WANG Lishu,et al.Numerical simulation of line source infiltration under moisture irrigation by HYDRUS[J].Yellow River,2016,38(4): 144-148.(in Chinese))

[10] 池宝亮,黄学芳,张冬梅,等.点源地下滴灌土壤水分运动数值模拟及验证[J].农业工程学报,2005,21(3): 56-59.(CHI Baoliang,HUANG Xuefang,ZHANG Dongmei,et al.Numerical simulation and validation of soil water movement undersubsurface drip irrigation with point-source emitter[J].Transactions of the Chinese Society of Agricultural Engineering,2005,21(3): 56-59.(in Chinese))

[11] 范严伟,马孝义,王波雷,等.膜孔灌土壤湿润体水分分布与入渗特性数值模拟[J].农业机械学报,2008,39(11): 35-41.(FAN Yanyanwei,MA Xiaoyi WANG Bolei,et al.Numerical simulation on soil wetting pattern moisture distribution and infiltration characteristics for film hole irrigation[J].Transactions of the Chinese Society for Agricultural Machinery,2008,39(11): 35-41.(in Chinese))

[12] 郭向红,孙西欢,马娟娟.降雨灌溉蒸发条件下苹果园土壤水分运动数值模拟[J].农业机械学报,2009,40(11):68-73.(GUO Xianghong,SUN Xihuan,MA Juanjuan.Numerical simulation for root zone soil moisture movement of apple orchard under rainfall-irrigation-evaporation[J].Transactions of the Chinese Society for Agricultural Machinery,2009,40(11):68-73.(in Chinese))

[13] 李淑芹,王全九.垂直线源入渗土壤水分分布特性模拟[J].农业机械学报,2011,42(3): 51-57.(LI Shuqing,WANG Quanjiu.Simulation of soil water distribution under vertical line source infiltration[J].Transactions of the Chinese Society for Agricultural Machinery,2011,42(3): 51-57.(in Chinese))

[14] 程勤波,陈喜,赵玲玲,等.饱和与非饱和带土壤水动力耦合模拟及入渗试验[J].河海大学学报(自然科学版),2009,37(3):284-289.(CHENG Qinbo,CHEN Xi,ZHAO Lingling,et al.Dynamic coupled modeling and infiltration experiment on soil water in saturated and unsaturated zones[J].Journal of Hohai University(Natural Science),2009,37(3):284-289.(in Chinese))

[15] 雷志栋,胡和平.土壤水研究进展与评述[J].水科学进展,1999,10(3): 311-318.(LEI Zhidong,HU Heping.A review of soil water research[J].Advances in Water Science,1999,10(3): 311-318.(in Chinese))

[16] VAN GENUCHTEN M T.A closed-form equation for predicting the hydraulic conductivity of unsaturated soils[J].Soil Science Society of America Journal,1980,44(5): 892-898.

[17] 杨红,邵孝侯,张林,等.天津滨海地区3种典型土壤水分物理特性[J].河海大学学报(自然科学版),2011,39(4): 431-435.(YANG Hong,SHAO Xiaohou,ZHANG Lin,et al.Moisture-physical characteristics of three types of soils in Tianjin Binhai area[J].Journal of Hohai University(Natural Science),2011,39(4): 431-435.(in Chinese))

[18] SIMUNEK J,VOGEL T,VAN GENUCHTEN M T.The SWMS-2D code for simulating water flow and solute transport in two-dimensional variably saturated media[R].Riverside:US Salinity Laboratory,1992.