层状土毛细水上升过程中Lucas-Washburn模型评价及修正

2018-12-13曹鼎峰孙梦雅魏广庆

水文地质工程地质 2018年6期

郝瑞，施斌，曹鼎峰，孙梦雅，魏广庆

(1. 南京大学地球科学与工程学院，江苏南京 210023；2. 苏州南智传感科技有限公司，江苏苏州 215123)

受淋溶-淀积、剥蚀-沉积、生物扰动等作用，层状结构土壤在自然界广泛分布[1]。层状土中水力特性的不连续，导致其水分运移和赋存特征与均质土壤中的情况完全不同。在毛细水上升过程中，土壤中的分层结构可能会阻碍(或加速)水分垂直方向的运动[2～4]。有关层状土中界面效应对水分运移规律影响的研究成果已较多。Hill等[5]发现土壤层状结构对水分入渗起抑制减渗作用。Baker等[6]指出水分在穿越细砂-粗砂的界面时会出现滞留。付志文等[7]将毛细水运移过程分为四个阶段，即惯性力作用阶段，黏性力-惯性力作用阶段，黏性力作用阶段和黏性力-重力作用阶段。Zettl等[8]比较了加拿大阿尔伯达省油砂地区田间深1 m的土壤剖面持水量，结果表明具有分层结构的土壤田间持水量高于均匀土壤，并有利于作物高产。因此，建立层状土中描述毛细水上升规律的模型非常重要。

目前描述毛细水运移最常用的渗吸模型有：Lucas-Washburn模型、Terzaghi模型、Handy模型、Bosanquet模型和纯惯性力模型等[9～12]。19世纪初，Lucas等[13]分析了等截面毛细管和多孔材料界面中水分自吸的作用因素，并提出了经典的Lucas-Washburn(LW)模型。后来很多关于自吸的理论分析和实验研究大多数在LW模型的基础上进行。蔡建超等[1]综述了常用的自发渗吸理论模型，总结了近十年来自发渗吸理论的研究进展和现状，分析了渗吸机理判别参数，简述了数值模拟研究及渗吸率影响机理的实验研究现状；Cai等[14]基于弯曲流线的分形特征，通过引入弯曲毛细管的弯曲度和分形维数，获得毛细管上升的高度/重量随时间演变的解析表达式。然而，目前大部分模型都是用于描述均质土中毛细水上升，关于层状土的研究还很少，主要原因是受测试技术的限制，层状土界面影响很难实际测定，从而难以进一步验证相关模型。

目前，对于毛细水运移与赋存特征的测试方法都是点式的，如取样烘干法、干湿计法和时域反射法，这些方法对均质土含水率的测试效果较好，但是空间分辨率低，难以更加精确测定层状土中水分信息；遥测法、地面热辐射测量法、卫星遥感法等测量精度较低，主要用于测量地表土层含水率[15～16]。而近几年发展起来的主动加热光纤法(Actively heated fiber optic method, AHFO)具有测试精度高、速度快、操作简便、可以连续分布式实时监测土中水分场等优势，曹鼎峰等[17]、严珺凡等[18]已对AHFO法测试水分场的可行性进行了验证。

本文将基于AHFO测试结果对LW模型进行评价，进一步对其进行改进，并通过室内试验测定改进前后LW模型的精度。

1 毛细水运移与赋存原理

1.1 Lucas-Washburn(LW)渗吸模型

毛细水主要受到基质吸力Fcap、毛细管侧壁的黏性阻力Fvisco和毛细水自身重力Fgrav作用。毛细水上升过程中液体充分发展的稳态阶段在整个上升过程中占有重要地位，此时的流体可视为层流，尤其是通过一个等截面圆管的流体流动被认为是Hagen-Poiseuille流动[19]，其流动的控制方程可直接由纳维-斯托克斯方程(Navier-Stokes)[20]求得。根据计算有：

(1)

式中：γ——毛细水液面张力；

φ——接触角；

R——平均空隙半径；

η——流体黏度；

ρ——毛细水密度；

g——重力加速度；

t——时间；

h——毛细水上升高度。

在黏性力作用阶段，Lucas和Wash burn忽略式(1)中惯性力和重力的作用，提出了该阶段的简化模型(LW)：

(2)

对于黏性力-重力作用阶段，Washburn在考虑了毛细水重力、黏性阻力和基质吸力作用，求解出式(1)的隐形解：

(3)

当基质吸力和毛细水重力相平衡时，毛细水上升到最大高度，可以求得最大高度heq为：

(4)

毛细水达到最大高度后，其运移过程趋于平衡，此时其表现出的特征称为毛细水赋存特征，常用土壤水分特征曲线表征。

1.2 土壤水分特征曲线Van Genuchten模型

土壤水分特征曲线又名土壤持水曲线(soil water retention curve，SWRC)，表征土壤中含水率与势能之间的关系，是研究非饱和土壤水分的保持能力和运移规律所用到的基本特性曲线，本文选用Van Genuchten模型以描述土水特征曲线[21]：

(5)

(6)

(7)

(8)

式中：θ(h)——土壤体积含水率；

θr，θs——残余含水率、饱和含水率/(cm3·cm-3)；

α,m,n——土壤水分曲线参数；

h——土壤水吸力/cm；

K(h)，KS——饱和渗透系数；

Se——有效饱和度[22～23]。

在毛细水上升过程中，当砂土层覆盖黏土层时，水分在界面处运移速率取决于两层土壤的土壤持水曲线(soil water retention curve，SWRC)，该曲线反映了土壤水的基质势(或土壤水吸力)与土壤含水率的变化关系，主要受到土壤颗粒粒径、矿物成分、孔隙率、孔径大小、容重等因素影响。砂土和黏土初始含水率相近，但下层黏土的基质吸力明显高于上层砂土的基质吸力[24]。因此，当下部毛细水湿润锋上升到界面处时，不能立即进入上部砂层，水流会在界面处滞留。直至黏土基质吸力随含水率升高而降低，并低于上层砂土的基质吸力时，毛细水流才会越过界面继续向上运动[24～27]，该现象即为“毛细屏障作用”。

图1显示了在毛细水湿润锋穿过砂土层达到平衡状态后层状土壤剖面中含水率的分布。当毛细水进入砂土层后，由于砂土的基质势低于黏土，因此，在毛细水运移过程稳定后，整个土壤剖面的含水率不连续，砂土层的含水率会出现低于黏土层的突变。

图1 平衡状态下层状土中土壤含水率示意图Fig.1 Moisture content of layered soils at equilibrium

1.3 AHFO技术测试原理

本文采用AHFO法测量毛细水运移过程与赋存特征。主动加热光纤法(Actively heated fiber optic method, 简称AHFO)，是利用热耗散原理测量土壤中含水率[28]。通过AHFO法测得的温度值满足：

(9)

式中：T——AHFO法测得土壤温度/ ℃；

P——单位长度加热功率/(J·m-1·s1)；

λ——土壤有效导热系数；

D——线性热源长度，本文中即为测管长度/m;

rb——线性热源半径，本文中即为测管半径/m;

Rb——土壤热阻系数/(K·m·W-1)；

T0——土壤初始温度/℃;

t——加热时间/s；

α——所测土壤热扩散系数/(m2·s1)。

式(9)可以简化为[28]:

(10)

土壤的导热系数与其含水率有关，含水率越高，土壤的导热能力越强，导热系数越大，在相同功率下加热相同时间时所达到的温度越低。Sayde等[29]研究建议，将温度特征值和土壤的含水率建立函数关系，根据所测温度特征值直接推算土壤含水率，可以大大降低求解误差。碳纤维加热光缆周边形成的温度场升温梯度稳定时，选取某个时间段[t1,t2]计算单位时间平均温度升高值，该平均温度升高值即为温度特征值Tt[30]。如图2所示，测量时通过碳纤维或金属材料对土体施加温度场，利用分布式光纤监测周围土体温度特征值Tt，根据温度特征值Tt与体积含水率θ的函数关系就可以计算出土壤的体积含水率θ。根据土壤含水率θ分布曲线，箭头所指含水率突增点为湿润锋所达高度，即毛细水上升高度。

图2 AHFO技术测试原理Fig.2 Principle of AHFO technology

曹鼎峰等[17]通过测试证明，单用指数函数、幂函数或对数函数模型不能满足整个时刻的需求，而分段函数模型具有测试精度高、标定简单、工程领域实用性强的优点。所以本文采用分段函数模型，根据文献[17]研究结果，含水率较低时采用对数函数，含水率较高时采用线性函数，其函数表达式为：

(11)

式中：θc——界限含水率，即Tt-θ-h曲线上曲率半径最小点所对应的含水率；

θs——饱和含水率；

a,b,c,k,d——常参数。

试验前通过取样烘干测得的含水率与温度特征值拟合出a,b,c,k,d等5个参数。

2 室内试验

为了测量均质土和层状土中毛细水运移规律，评价分析LW模型在毛细水运移表述中的误差，设计了均质土和层状土的毛细水上升模型试验。试验装置如图3所示。试验装置主要包括铁桶和水槽。铁桶高度100 cm，直径40 cm，上部有桶盖防止水分蒸发，下部在距桶底3 cm处铁桶侧面均匀设置8个半径2 cm小孔，保证铁桶内部和水槽内保持相同水位高度。铁桶放置于1个深度为10 cm、直径60 cm水槽中。为防止土壤堵塞下部小孔，在桶底铺设厚10 cm、粒径4～8 mm砾石层作为反滤层。在砾石层上覆1层纱布防止上部黏土进入。

图3 试验装置Fig.3 Experimental apparatus

为了提高光纤测量的空间分辨率，将碳纤维加热光缆紧密缠绕直径50 mm的PVC管，制成碳纤维感测光缆测管，简称测管，设备与测管如图4所示。并在测管外侧紧套热缩管防止水分沿光缆上移产生误差。光纤的空间分辨率为1 m，因此沿测管方向空间分辨率达到20 mm，满足室内试验要求。

图4 碳纤维测管与仪器实物图Fig.4 Carbon fiber tube and instrument physical map

试验所采用的工业砂和高岭土初始含水率分别为7.7%、4.7%，工业砂干密度为1.53 g/cm3，高岭土干密度为2.53 g/cm3，工业砂与高岭土按照9∶1、7∶3比例均匀混合后土壤颗粒级配曲线如图5所示。本次所用DTS解调仪主要技术指标如表1所示。

图5 土壤颗粒分析结果Fig.5 Particle analysis results of sand material

光纤类型测量距离/km空间分辨率/m测温范围/ ℃测温精度/ ℃功率/W62.5/125101-40～1200.1300

试验时将工业砂和高岭土分别按照比例9∶1和7∶3均匀混合，得到砂土和黏土材料。在桶底砾石层上覆盖一层纱布，将光纤测管埋置于铁桶正中心固定，再将试验材料分层装满铁桶，每次装填10 cm后击实，并做刨毛处理，减小分层装填对模型整体性的影响。对于均质土模型，分别直接用按9∶1和7∶3比例均匀混合的材料分层填装满铁桶；对于层状土模型，先填装厚度20 cm 按7∶3比例均匀混合的黏土，再填装厚度10 cm按9∶1比例均匀混合的砂土，再用黏土将铁桶填装满。在初始状态下测得模型中土壤的初始含水率，然后向水槽中持续加水保持满水状态。前24 h通过AHFO法每1 h测试1次含水率，24 h后每隔12 h左右测1次含水率并取样。

3 试验结果

3.1 含水率与温度特征值标定

已知温度特征值Tt与土壤含水率θ之间的函数表达形式，参照文献[17]取样烘干分段拟合的标定方式，分别拟合出工业砂与高岭土按照9∶1、7∶3均匀混合土函数表达式中的常参数a,b,c,k,d，其表达式为：

9∶1土样

Tt=

(12)

7∶3土样

Tt=

(13)

3.2 土壤含水率运移与分布

基于AHFO技术测得不同时刻均质土与层状土剖面中土壤含水率分布如图6所示，图6(a)为工业砂和高岭土按照比例7∶3均匀混合所得均质黏性土中不同时刻体积含水率分布曲线，图6(b)为层状土中不同时刻体积含水率分布曲线，图中水平虚线为分层界面，0～20 cm为7∶3均质黏土，20～30 cm为9∶1均质砂土，30～100 cm为7∶3均质黏土。图中含水率突增点所对应高度即为毛细水湿润锋所达到高度。可以看出，随着时间延长，两种土毛细水上升高度逐渐增高，湿润锋以下土壤含水率也逐渐增大。而湿润锋上升速率和含水率增长速率逐渐降低。主要是因为随着毛细水量的增加，重力势增大，基质势不变，毛细水所受向上驱动力减小。当含水率增加到重力势与基质势相同时，向上驱动力为零，毛细水达到最大上升高度。在毛细水湿润锋达到20 cm前，两者运移规律相同，但达到20 cm后，层状土砂土层中含水率急剧下降，远低于均质土中同层位土壤含水率。这是因为砂土基质吸力突然变小，其持水能力降低。而且湿润锋超过30 cm后，从下到上均质土中含水率缓慢下降，但层状土中含水率却先升高后降低。这是因为上层黏土层基质吸力变大，持水能力增强。

图6 不同时刻均质土的含水率剖面Fig.6 Volumetric moisture content of the homogeneous clay

为了便于观察分析毛细水湿润锋随时间升高规律，由图6中含水率分布曲线绘制出均质砂土、均质黏土和层状土湿润锋高度随时间变化曲线，即毛细水上升高度随时间变化关系曲线(图7)。

图7 毛细水上升高度随时间变化曲线Fig.7 Curve of the rise of capillary water with time

由图7可以发现，当毛细水上升到20 cm时，毛细水湿润锋到达黏土层(下部)与砂土层(上部)界面高度，均质黏土和均质砂土中毛细水持续稳定上升，20 h时毛细水分别上升到25 cm和30 cm，但层状土中毛细水发生滞留与减缓，20 h时毛细水仍滞留在20 cm高度。原因是下层黏土的基质吸力大于上层砂土，毛细水不能立即进入砂土层。但是界面以下黏土中水分持续增加，基质吸力不断减小，直至下层黏土中基质吸力随着含水率的升高而下降到与上层砂土的基质吸力相同时，毛细水才能越过界面继续向上运移进入砂土层，40 h时毛细水达到25 cm。该现象即“毛细屏障作用”。

80 h后，均质黏土和均质砂土中毛细水高度均已超过40 cm，上升速率缓慢，而层状土中毛细水达到30 cm进入黏土层，其上升速率加快，与均值土毛细水高度差距逐渐减小。根据田町正誉[31]的解释，当小孔径毛细管连通在大孔径毛细管之上时，毛细水的上升速率明显加快，是因为砂土层中大孔径毛细管促进了水分在黏土小孔径毛细管中的运动。究其本质，是由于在接近地下水位处，土壤剖面的含水量较大，砂土的导水率明显大于黏土所致[32]。

对比均质砂土与均质黏土上升曲线，85 h前黏土中毛细水上升速率低于砂土，根据张志权[33]的研究，毛细水前期的上升速度与土壤毛细管通畅性相关，黏粒含量越高，毛细管截面内被黏粒强烈吸附的结合水越多，阻碍了自由水向上运动，所以毛细水上升速度越慢，与本试验实际测得的结果相同。

根据Van Genuchten模型理论，在黏土层与砂土层界面处含水率会发生突变。而图6(b)中实际测得的结果却是连续缓变的。选取120 h实测界面处含水率与理论值绘制在图8中，可以看出在层内含水率实测值与理论值很吻合，但是在界面处吻合度较低。这是由AHFO技术测量原理导致的。已知测管空间分辨率为2 cm，即测管每一点监测值为上下各1 cm内监测值的平均值，因此在含水率的突变处由光纤测管测得的突变程度是低于实际突变程度的。

图8 界面处含水率实测值与理论值Fig.8 Measured and theoretical values of moisture content at the interface

3.3 Lucas-Washburn渗吸模型预测与修正

通过Lucas-Washburn模型可以计算预测毛细水湿润锋上升高度，根据文献[7]中所给出ρ,γ,φ,η的经验值，并通过文献[34]中给出的平均空隙直径D0与d20和孔隙率之间的关系D0=0.63nd20，代入毛细水在黏性力作用阶段和黏性力-重力作用阶段的表达式(2)和式(3)，绘制出Lucas-Washburn渗吸模型预测曲线和层状土湿润锋高度实测值曲线如图9所示。

根据图9层状土中毛细水湿润锋实测高度和Lucas-Washburn模型预测值，在高度达到20 cm前模型预测值较准确，但毛细水进入砂土层后该模型无法预测湿润锋高度。一是因为砂土层减缓了毛细水上升速度，二是因为毛细水穿过砂土层重新进入上层黏土层后，下部砂土层导水率加大，有利于毛细水运移。据此，为了描述毛细水在上层黏土中运移规律，本文对Lucas-Washburn模型进行修正。由于上述两点原因的主要影响因素是毛细管径，因此，定义穿过砂土层后上层黏土的等效毛细管半径R为：

(14)

式中：R1，R2——砂土和黏土的平均空隙半径/m；

D1，D2——砂土层厚度和毛细水进入上部黏土层高度/cm。

则Lucas-Washburn模型黏性力-重力作用阶段修正后Improved Lucas-Washburn模型(ILW模型)表达式为：

(15)