1.4 土壤水分运动参数

土壤水分运动参数包括土壤水分特征曲线、非饱和导水率和土壤水分扩散率,是利用土壤水分运动基本方程分析土壤水分运动特征不可缺少的参数。

1.4.1 土壤水分特征曲线

1.4.1.1 土壤水分特征曲线的定义与影响因素

土壤水的基质势或土壤水吸力是随土壤含水率而变化的,其关系曲线称为土壤水分特征曲线,也称为土壤持水曲线,反映了土壤水分静态能量特征。由于土壤水基质势一般为负值,为了应用方便起见,将负的基质势称为土壤水吸力。因此土壤水分特征曲线也可解释为土壤水吸力与含水量之间函数关系。随着土壤含水量增加,土壤水吸力降低。该曲线表示土壤水的能量和数量之间的关系,反映了土壤水分的基本特性,如图1.1所示。

由图1.1可以看出,当土壤中的水分处于饱和状态时,含水率为饱和含水率θs,而吸力h为零。

图1.1 土壤水分特征曲线示意图

土壤水分特征曲线斜率的倒数,即单位基质势的变化所引起的含水量变化,称为比水容量,用c表示。c值随土壤含水率θ或土壤水吸力h而变化,所以也表示为cθ)或ch)。比水容量c是分析土壤水分运动的重要参数之一。

土壤水分特征曲线受到众多因素的影响,主要包括土壤质地、容重、结构、温度、有机质含量、湿润方式等。一般而言土壤黏粒含量愈高,同一吸力对应的土壤含水量愈高。土壤质地不仅影响土壤含水量所对应的土壤吸力大小,而且影响土壤水分特征曲线形状。对于砂质土壤而言,一般在高吸力区曲线比较陡,而在低吸力部分土壤水分特征曲线比较平缓。随着容重增加,同一吸力相应土壤含水量升高,土壤饱和含水量减少,土壤水分特征曲线变陡。随着土壤有机质含量增加,土壤团粒结构增加,高含水量段变缓。土壤温度通过改变土壤水分黏滞性和表面张力而影响土壤水分特征曲线。一般随着温度升高,基质势增加。土壤水分特征曲线与土壤水分变化过程也存在密切关系。

图1.2 土壤水分特征曲线滞后现象示意图

在其他条件一致的情况下,在脱湿和吸湿过程中测定的土壤水分特征曲线也存在明显差别。在同样含水量条件下,脱湿过程的吸力较吸湿过程要大,将这一现象定义为滞后效应,如图1.2所示。滞后现象在轻质土壤中表现较为明显,其形成机制较为复杂。目前对滞后效应的解释存在三种理论,即瓶颈理论、接触角理论和弯月面延迟形成理论。

1.4.1.2 土壤水分特征曲线的函数形式

由于土壤水分特征曲线的影响因素复杂,至今尚没有从理论上建立土壤含水量和土壤基质势之间的关系,通常用经验公式来描述。常用的公式有如下几种形式:

(1)Gardner模型

式中:h为土壤吸力,L;θ为土壤含水率;ab为经验参数。

(2)Brooks-Corey模型

式中:θ为土壤含水率;θs为土壤饱和含水率;θr为土壤滞留含水率;h为土壤吸力,L;hd为进气吸力,L;N为形状系数。

(3)van Genuchten模型

式中:θ为土壤含水率;θs为土壤饱和含水率;θr为土壤滞留含水率;h为土壤吸力,L;α为与进气吸力相关的参数;mn为形状系数,一般m=1-1/n

Gardner模型是较早采用的土壤水分特征曲线的经验模型,由于参数少、形式简单,在早期的研究中应用较多。但由于其没有考虑滞留含水量,对高吸力段的描述可能存在偏差;曲线在近饱和段没有弯点,饱和点也与实际不符,所以该模型不能很好地描述实测土壤水分特征曲线,在目前对土壤水力特性精度要求较高的土壤水分运动定量模拟中采用不多。

从Brooks-Corey模型可以看出,当土壤处于饱和状态时,土壤吸力等于进气吸力,因此该模型描述了脱湿过程的土壤水分特征曲线。

在van Genuchten模型中,当土壤含水量处于饱和状态时,土壤水吸力为零。因此,该模型描述了土壤吸湿过程的土壤水分特征曲线,同时由于该模型能够配合大部分土壤水分特征曲线的形状,因此得到广泛应用。为了比较式(1.24)和式(1.25)所包含参数间关系,将式(1.25)近似简化,忽略式(1.25)右边括号中分母的1,则变为

若令α=1/hdmn=N,则式(1.26)与式(1.24)相同,当然这种转化存在一定误差。因此,可以认为van Genuchten模型与Brooks-Corey模型之间存在着一定的联系。

1.4.1.3 土壤水分特征曲线确定方法

土壤水分特征曲线确定方法主要有两大类,即直接测定方法和间接推求方法。直接测定方法可以直接准确测定出土壤水分特征曲线,但费时费力,需要特殊设备;而间接确定方法相对比较简单,但准确性和参数不唯一性影响该方法的应用。

1.直接测定方法

实验室内测定土壤水分特征曲线的方法主要有张力计法、砂性漏斗法、压力膜法、离心机法等,田间原位测定土壤水分特征曲线大都用张力计法。在测定过程中,土壤水分含量的测定方法主要有取土烘干法、频域反射仪法(Frequency-Domain Reflectometry,FDR)、时域反射仪法(Time-Domain Reflectometry,TDR)等。

通过试验测定各含水率及相对应的吸力,便可依据最小二乘法等方法对于上述Gardner模型、Brooks-Corey模型、van Genuchten模型进行拟合,从而得到模型参数。以下是几种常用的实验测定方法。

(1)张力计法:该方法是在一定含水量下测定土壤的水吸力,其原理是:当陶土头插入被测土壤后,管中的纯自由水便通过多孔陶土壁与土壤水建立了水力联系。由于仪器中自由水的势值总是高于非饱和土壤水的势值,所以张力计管中的水很快流向土壤,同时在管中形成负压。当仪器内外的势值达到平衡时,由与管相连通的真空表或负压传感器测得的负压就是土壤水(陶土头处)的吸力值。该方法既可以用于室内扰动土样和原状土样水吸力测定,也可以用于田间土壤水吸力测定,试验设备和操作都比较简单。但是,此方法只能测定0~80kPa吸力范围(低吸力)内的水分特征曲线,而且实验时间较长。

(2)压力膜法:与砂性漏斗法的实验过程相反,砂性漏斗法是对土样施加吸力,压力仪法则是对土样施以压力。压力室内气压增加,陶土板上土样与板下水室的自由水相联系,当压力室气压增加,陶土板上土样中的总土水势高于板下自由水水势时,土样便开始排水,直到陶土板上的土壤水和陶土板自由水水势相等为止,此即为平衡状态。土样的基质势等于压力室内压力的负值。在已知土样初始含水率的情况下,根据排出的水量就可以算出相应的土壤含水率。该方法仅用于室内实验,其优点是能得到完整的水分特征曲线,缺点是不易得到吸湿曲线(土壤水分逐渐增加的曲线)。

(3)离心机法:用离心率法测定土壤水分特征曲线,实际是将重力场的装置移至离心力场。将粉碎筛分后的土样装入底部有滤纸的离心管中压实,放在浸水中的纱布上约2~4h,即达到毛管饱和,将此离心管中土样放在离心机中旋转过一定时间后称重一次,记录离心机转速及土样失水量即可计算出土壤水吸力与含水率之间的关系曲线。这种方法测定范围较大,可测定高达20个大气压的负压值,但只能得到脱湿曲线(水分逐渐减少)。

2.间接推求方法

直接方法存在着耗时长、工作量大等缺点,因此人们也在积极探索估计土壤水分特征曲线的间接方法。

间接估计土壤水力性质的方法主要包括:土壤转换函数方法(Pedotransfer functions,PTFS)、分形方法(Fractal method)、土壤形态学方法(Morphological method)、数值反演方法(Numerical inverse method)等。此处仅对该4种方法进行简要介绍。

土壤转换函数法:该方法是利用容易获得的土壤物理性质,如土壤颗粒大小分布、容重和有机质含量等,来估计土壤的水力性质。目前,已有研究得到了van Genuchten模型中的参数与土壤颗粒大小分布、容重和有机质含量等的分布多元线性回归关系。

分形方法:用分形维数表征土壤孔隙弯曲程度,建立描述结构性土壤的孔隙和团聚体大小分布的分形模型,于是可以根据颗粒大小分布资料和分形理论得到土壤水分特征曲线的模型。

土壤形态学方法:其基本思想是通过对高分辨率的系列土壤剖面图像分析处理,确定土壤的孔隙大小分布和连通情况,在此基础上建立网格模型,预测土壤的水力性质。

数值反演方法:该方法根据可控制的初、边值条件,通过求解Richards方程,并与水分分布的实测结果相比较,从而确定Gardner模型、Brooks-Corey模型、van Genuchten模型的参数。

1.4.2 非饱和导水率

1.4.2.1 非饱和导水率的概念及其表达式

非饱和导水率是指非饱和土壤在单位水势梯度作用下单位断面面积的流量(即水流通量)。与饱和土壤导水率(渗透系数)不同,非饱和导水率是土壤含水率的函数,随着土壤含水率的增大而增大,饱和导水率为其极限值;当然非饱和导水率也是土壤吸力的函数,随着土壤吸力的增大而减小。一般来说,粗颗粒土壤的非饱和导水率随含水率的变化,要比细颗粒土壤来得更加剧烈。

常用的非饱和导水率表达式有

式中:Kθ)、Kh)分别为以土壤含水率、吸力为自变量的非饱和导水率,LT-1h为土壤吸力,L;θ为土壤含水率;θs为土壤饱和含水率;θr为土壤滞留含水率;αcβNmM为经验参数。

1.4.2.2 瞬时剖面法测定非饱和导水率

瞬时剖面法是在室内进行均质土壤的一维上渗或下渗试验时,测定不同时刻土壤剖面的含水率和吸力分布,通过计算求得非饱和导水率Kθ)。瞬时剖面法对扰动土和原状土均适用,可测定吸湿和脱湿过程,试验和计算都比较方便,因此应用较为普遍。

1.基本原理

对非饱和垂直一维水分运动,当取z坐标向上为正时,由达西定律可导出非饱和导水率的计算式:

式中:Kθ)为非饱和导水率,LT-1h为土壤基质吸力,L;q为水流通量,LT-1

由式(1.33)可知,只要知道了某一点的土壤含水率θ、水流通量q和吸力梯度∂h/∂z,就可计算出其相应于土壤含水率θ的非饱和导水率Kθ)。为此,试验需要测得两个时刻t1t2的土壤含水量和吸力剖面(图1.3)。吸力梯度可直接由吸力剖面计算近似得到,并取两个时刻的均值;各个位置的水流通量可以通过水量平衡原理、水流连续性原理获得。

2.试验过程

瞬时剖面法测定非饱和导水率常用垂直土柱上渗试验装置,如图1.4所示。土柱顶端不密封,但要防止蒸发。利用马氏瓶维持水位不变,并测出补给水量,可以得到入渗表面(z=0)的水流通量q0

试验过程中,设定两个观测时刻t1t2,通常采用TDR或FDR观测土壤含水量得到相应时刻的剖面含水量分布,同时读取张力计读数得到相应时刻的吸力分布。根据水量平衡原理,土柱内含水量的增加量应该等于马氏瓶的补给水量,由此可以计算任意位置z处的水流通量qz),由此可进行水量平衡验证。由水流连续性原理得

图1.3 土壤含水率分布曲线和水吸力分布曲线

(a)含水率分布曲线;(b)水吸力分布曲线

对式(1.34)在[0,z]范围内积分得

于是有

图1.4 垂直土柱上渗试验装置示意图

式中的后两个积分项分别是t1t2时刻[0,z]范围土柱的土壤储水量,可根据观测的土壤水分分布采用近似法计算得到。

1.4.2.3 确定Brooks-Corey模型参数的水平吸渗法

水平一维土壤水分运动的达西定律为

式中:q为土壤水分通量,cm/min;Kh)为非饱和导水率,cm/min;h为土壤水吸力,cm;x为横坐标,cm。

土壤水分特征曲线和非饱和导水率仍采用Brooks-Corey(1964)模型

h<hd时,S=1。其中,S为有效饱和度,nM为参数。

水平一维土壤水分运动基本方程为

式中:θi为初始含水量,cm3/cm3θs为饱和含水量,cm3/cm3;其余符号意义同前。

假定任意时刻土壤水吸力分布可以表示为

式中:hx)为任意坐标点x处的吸力,cm;xf为湿润锋距离,cm;a为系数。

这样土壤水分通量表示为

x=0,土壤水分通量q(0)即是土壤入渗率i,即

对式(1.40)进行积分,得

累积入渗量与湿润锋关系可表示为

根据土壤水分特征曲线和土壤水吸力分布表达式,土壤水分分布表示为

土壤累积入渗量表示为

为了简单起见,将式(1.54)、式(1.43)和式(1.48)分别表示为

其中:

因此,参数nMhd可表示为

Brooks-Corey模型中包括6个参数(θrθsKshdMn),其中θr可根据风干土壤含水量进行确定,θs可根据土壤容重进行计算或在通过实验方法直接测定,Ks可利用实验方法直接测定。其他3个参数可以利用累积入渗量和入渗率与湿润锋距离关系以及湿润锋随时间关系进行计算。具体为利用累积入渗量与湿润锋距离关系计算参数n,利用入渗率与湿润锋深度关系计算参数hd,利用湿润锋深度与时间关系可以计算参数M

1.4.2.4 根据土壤水分特征曲线估算土壤非饱和导水率

虽然人们提出了多种测定土壤水分运动参数的室内试验方法,但这些方法在田间都存在一些问题,所以难以获得非饱和土壤水分运动模型参数是限制水分运动数值模拟技术应用于田间的重要因素。鉴于土壤水分特征曲线相对容易获得,土壤水分特征曲线、非饱和导水率均与土壤孔隙分布之间存在明确的关系,于是一些学者通过对土壤孔隙分布特征的概化,建立了土壤水分特征曲线与非饱和导水率间的函数关系。最具代表性的是Burdine和Mualem建立的由土壤水分特征曲线预报非饱和导水率的模式。基于这些模式,Brooks-Corey和van Genuchten先后提出根据不同的土壤水分特征曲线形式获得非饱和导水率、扩散率的计算方法。

1.Brooks-Corey模型

Burdine建立的非饱和导水率预测模型为

其中,Θ为有效饱和度,即

式中:hx)为土壤水分特征曲线;x为积分符后的代换变量;l为孔隙弯曲度。

Brooks-Corey 采用的土壤水分特征曲线模型为

式中:h为土壤吸力,L;hd为进气吸力,L;N为形状系数。

将式(1.59)代入式(1.57),积分可得非饱和导水率的表达式为

Brooks-Corey采用孔隙弯曲度l=2,再令

于是,根据土壤水分特征曲线(1.57)得到非饱和导水率

根据土壤水分扩散率的定义,由上述土壤水分特征曲线、导水率公式可得到:

2.van Genuchten模型

Mualem建立的非饱和导水率预测模型为

式中:各符号意义同前。

van Genuchten采用的土壤水分特征曲线表达式为

式中:θ为土壤含水率;θs为土壤饱和含水率;θr为土壤滞留含水率;h为土壤吸力,L;α为与进气吸力相关的参数;mn为形状系数。

van Genuchten取m=1-1/n、孔隙弯曲度l=1/2,将该土壤水分特征曲线表达式代入Mualem建立的非饱和导水率预测模型中,积分得到非饱和导水率表达式:

根据土壤水分扩散率的定义,由上述土壤水分特征曲线、导水率公式可得

1.4.3 非饱和土壤水分扩散率

为了在以土壤含水率为因变量和以基质势(或吸力)为因变量的土壤水分运动方程的转换中便于表达,定义了非饱和土壤扩散率,其定义式为

其中, cθ)为土壤容水度,即

式中:Dθ)为非饱扩散率,L2T-1Kθ)为非饱和导水率,LT-1θ为土壤体积含水率;ψ为土壤基质势,L;h为以水头表示的土壤吸力,L。

1.4.3.1 水平入渗法测定非饱和土壤水分扩散率

水平入渗法是测定土壤水分扩散率Dθ)的非稳定流法,最早由Bruce和Klute提出。该方法利用半无限长水平土柱吸渗试验资料,结合解析公式计算出土壤水分扩散率Dθ)。

图1.5 水平土柱渗吸试验装置示意图

试验装置如图1.5所示。试验前按一定容重装土,保证土柱初始含水量和各点密度均匀。为保证进水端含水率保持不变,采用马氏瓶供水。在试验结束时从湿润锋开始向进水端等间隔迅速测出从进水端到湿润锋之间的n段各断面土壤含水量分布,即θii=0,1,2,…,n),其中i为断面编号。

水平入渗法测定土壤水分扩散率的基础是一维水平土壤水分运动的Boltzmann变换解。问题的控制方程和定解条件为

式中:θ0为初始土壤含水量;θb为进水端边界土壤含水量(接近饱和含水率)。

式(1.68)中,后两式分别为左边界条件和右边界条件。左边界为进水端,含水率接近饱和含水率;右边界为渗吸试验过程影响不到(湿润锋未达)之处,含水率仍然保持初始含水率。该方程为非线性偏微分方程,采用Boltzmann变换得到解析解:

其中,λ为Boltzmann变换参数,即

式中:Dθ)为土壤水分扩散率,L2T-1θ0为初始土壤含水量;θb为进水端边界土壤含水量;x为自进水端沿水平土柱渗水方向的坐标,L;t为渗吸结束取土时的入渗历时,T。

由式(1.69)可以得到渗吸试验结束时(t时刻)参数λ沿x坐标的分布,同时由试验观测得到同时刻的土壤含水率θ沿x坐标的分布,于是可得到λ-θ关系曲线。将式(1.69)离散化,便可得到土壤水分扩散率的计算式

式中:Dθ)为土壤水分扩散率,L2T-1n为从进水端到湿润锋之间的分段数;i为从进水端到湿润锋之间的观测断面编号;θi为第i断面的含水率;λi为第i断面的Boltzmann变换参数,LT-1/2

1.4.3.2 确定土壤水分扩散率的水平一维分析法

目前常采用Bruce和Klute(1956)方法进行土壤水分扩散率的测定。该方法主要依据土壤水分运动基本方程的Boltzmann变换结果,利用水平一维土壤吸湿过程中的含水量剖面计算土壤的扩散率。因此该方法需要准确测定土壤含水量剖面,并对含水量剖面进行人为光滑,增加人为误差。在土壤入渗过程中,湿润锋距离比较容易测定,下面介绍相关理论。

水平一维土壤水分运动的达西定律表示为

式中:q为水分通量,cm/min,kh)为非饱和导水率,cm/min,h为土壤吸力,cm,x为坐标,cm,Dθ)为扩散率;θ为含水量,cm3/cm3

非饱和导水率可表示为(Brooks-Corey,1964):

式中:Ks为饱和导水率,cm/min;hd为进气吸力,cm;M为参数。

土壤水分特征曲线表示为(Brooks-Corey,1964):

式中:θs为饱和含水量,cm3/cm3θr为滞留含水量,cm3/cm3n为参数。

水平一维吸渗条件下土壤水分运动的基本方程为

式中:θi为初始含水量,cm3/cm3;其余符号意义同前。

王全九(2004)对式(1.75)进行分析,得到了入渗率i表达式为

式(1.76)表示了入渗率与湿润锋距离间关系,即入渗率与湿润锋距离的倒数呈线形关系。

土壤吸力分布可表示为

土壤含水量分布为

土壤累积入渗量表示为

如初始含水量比较低,且认为滞留含水量与初始含水量相等时,累积入渗量简化为

由于入渗率是累积入渗量的导数,即

这样可得到湿润锋距离与时间关系:

式(1.82)显示湿润锋距离与时间平方根呈线形关系。

扩散率是土壤水分特征曲线和非饱和导水率的函数,即

并且有

在式(1.85)和式(1.86)中的(1-M)/n可以利用式(1.80)进行计算。kshd/(M-1)可利用式(1.76)进行计算。因此,当测定累积入渗量、入渗率湿润锋距离和随时间变化过程就可计算土壤水分扩散率。