江汉平原北部黏土层土壤水分特征曲线的测定与模拟

石浩楠，陈植华，胡 成，黄 琨，刘添文，王 清

(1.中国地质大学(武汉)环境学院，湖北 武汉 430074；2.中国地质调查局武汉地质调查中心，湖北 武汉 430205)

非饱和带是连接地表水、土壤水和地下水的纽带，存在着气体、液体等流体的流动和多种物质与成分之间的迁移转化过程，是人类活动的最根本的载体[1]。包气带是大气水和地表水同饱和带地下水进行水分和能量交换的枢纽，也是地表污染物进入地下水的通道，该带内的水分运动对整个水循环活动及土壤污染物迁移起着至关重要的作用[2-3]。土壤水分特征曲线是描述土体吸力h与土壤含水量θ之间关系的曲线(h-θ曲线)，反映了土壤水能量和数量的关系，在研究非饱和带土壤水分流动、溶质运移以及土壤污染质迁移转化过程中有着非常重要的作用[4]，也是研究“地表水—土壤水—地下水”三水转化机理以及非饱和带水汽热耦合的一个基本参数。非饱和带土壤水分特征曲线反映了土壤的持水性能及保水性能，在农业农田灌溉中有着非常关键的指导作用[5]。研究大别山区-江汉平原三水转换关键带黏土层典型剖面的土壤水分特征曲线，可以为研究区域内包气带水分和溶质的迁移、转化以及数值模拟提供理论和数据支持，可进一步探明非饱和带水分运移规律，对当地水资源的开发利用与保护也有着重要意义。

近些年来，学者们通过大量的研究建立了许多经验公式来拟合土壤水分特征曲线，比较常用的有：van Genuchten模型(1980)[6]、Gardner模型(1970)[7]、Gardner-Russo模型(1988)[8]、Brooks-Corey(1964)模型[9]、Frdlund and Xing模型(1994)[10]，其中van Genuchten模型、Gardner模型应用较为广泛[11]。一些学者通过评价van Genuchten模型和Gardner模型对不同地区不同条件下土壤水分特征曲线的拟合效果，结果发现van Genuchten模型的拟合效果稳定，适应范围更广，吻合度和拟合精度均较高[12-16]。

1 研究区概况与研究方法

1.1 研究区概况

研究区位于湖北省孝感市肖港镇，北纬31°00′～31°10′，东经113°45′～114°00′，地处江汉平原东北部，北接大别山，南连江汉平原，位于大别山-江汉平原过渡带。该地区属亚热带大陆季风性气候，雨量充沛，年平均气温为16.2℃，年平均降水量为1 146 mm，海拔高度为20～50 m。

图1 研究区域水文地质图

1.2 土壤样品的采集与测定

在研究区黏土层典型剖面处分别采集0～0.2 m、0.2～0.5 m、0.5～0.9 m、0.9～1.4 m、1.4～2.0 m、2.0～2.5 m、2.5～3.0 m、3.0～3.5 m、3.5～4.0 m、4.0～5.0 m、5.0～6.0 m、6.0～6.5 m深度的土壤样品，采用取土环刀取原状土样，同时在环刀取样位置的周边取1.5 kg扰动土样，并进行密封、遮光处理，防止土壤水分损失。本次共采集研究区原状土壤样品24个，用于测定土壤的水分特征曲线和物理性质，并选取不同层位的原状土样进行物理性质测试，得到研究区不同深度土壤样品的物理性质参数见表1。

表1 研究区不同深度土壤样品的物理性质参数一览表

1.3 土壤水分特征曲线的测定方法及拟合模型

本试验所用的仪器为美国SOILMOISTURE公司生产的5bar压力膜仪，采取压力膜法测定土壤的水分特征曲线，其测定方法为：将采取的原状土样处理后，水中浸泡24 h直至饱和，设定10个压力水头，分别为100 cm、300 cm、500 cm、700 cm、1 000 cm、1 500 cm、2 000 cm、3 000 cm、4 000 cm、5 000 cm，对不同层位土壤样品土壤水分特征曲线进行测定，同一层位同时进行3次重复试验。

由于van Genuchten模型和Gardner模型有着广泛的适用性和较好的拟合性，因此本次采用van Genuchten模型和Gardner模型分别对研究区不同层位土壤的水分特征曲线进行拟合。

Gardner模型的表达式为[7]

h=aθ-b

(1)

van Genuchten模型的表达式为[6]

(2)

式中：θ为土壤体积含水量(cm3/cm3)；θs为土壤饱和体积含水率(cm3/cm3)；θr为土壤残余体积含水率(cm3/cm3)；h为压力水头(cm)；α为与土壤进气值有关的参数；a、b为大于零的正常数；m、n为曲线形状参数，其中m=1-1/n。

上式中a、b、α、m、n为拟合模型参数，本文通过非线性最小二乘法并借助RETC 软件求解拟合模型的各参数。

2 研究结果与分析

2.1 土壤水分特征曲线的测定

土壤水分特征曲线是土壤质地、结构和孔隙等物理性质综合作用的结果，是土壤物理特性的表征[18]。研究区不同深度土壤水分特征曲线的实测结果，见图2。

由图2可以看出：

(1) 研究区不同深度土壤的体积含水量(即排水量)伴随土体吸力的变化呈规律性变化：在低吸力(<1 000 cm)阶段，各个层位土壤体积含水量整体快速降低，土壤水分特征曲线较缓，比水容量较大，排水主要在相对大孔隙中进行，随着土体吸力的增大，该曲线逐渐变陡，由大孔隙排水转变为中小孔隙排水，比水容量减小，在此阶段5.0～6.0 m层位土壤的持水能力较弱；在中吸力(1 000～5 000 cm)阶段，土壤排水主要来自毛细孔隙内水分的变化，此阶段由于土体对水分的强吸持能力，水分变化幅度较小，土壤水分特征曲线较陡，在此阶段2.0～2.5 m层位土壤的持水能力较弱。

图2 研究区不同深度土壤水分特征曲线的实测结果

(2) 研究区不同深度土壤的持水能力整体表现为1.4～2.0 m层位最强和2.0～2.5 m、5.0～6.0 m层位较弱，这可能是由于2.0～2.5 m、5.0～6.0 m层位的土壤颗粒之间排列相对较松散，孔隙结构较疏松，且大孔隙较多，而大孔隙具有较小的基质势，脱湿状态下使得土壤水从大孔隙优先排出，土样在极小的吸力条件下就能快速失水；相反，细质地土壤的孔隙结构密实，大孔隙减少，中小孔隙增多且连通性变差，使得土壤的水分运移速率较粗质地土壤慢，土壤水分特征曲线较陡，从而表现出较好的土壤持水能力，这与赵雅琼等[19]的研究结论一致。

(3) 研究区0～0.2 m、0.2～0.5 m、1.4～2.0 m层位的土壤，在吸力水柱为200～1 200 cm阶段，随着土样吸力的增长，土壤没有排水现象，表明此阶段土壤的持水能力较强；在吸力水柱大于1 200 cm阶段，随着土样吸力的增长，土壤排水相对增多，土壤水分特征曲线较缓，这可能与土样吸力到达某一临界值，土壤内孔隙发生突变，导致之前较低吸力集聚的水分突然释放出来有关。

上述试验结果表明，研究区不同深度土壤的水分特征曲线表现形态各异，说明土壤的持水性能各不相同：在低吸力阶段，除了0.2～0.5 m层位的土壤几乎没有排水外，其余各个层位土壤快速排水，比水容量大，土壤水分特征曲线较缓，表明饱和状态的土壤变为非饱和状态土壤的最初阶段，土体基质吸力对土壤含水量的变化非常敏感；在中吸力阶段，土壤整体排水变缓，比水容量较小，土壤水分特征曲线较陡，表明随着土壤含水量的减少，其饱和度降低，土体基质吸力对土壤含水量变化的敏感性逐渐降低，当土壤饱和度降至所能维持的最低饱和度时，土体基质吸力达到最大。

2.2 土壤当量孔径的分布

土壤水分特征曲线可以间接地反映土壤孔隙大小的分布，而孔隙的形态和分布情况直接影响着土壤中水分运动和溶质迁移[1]。为了更好地考察研究区不同深度土壤孔隙的分布情况，对研究区不同层位土壤的孔隙分布情况进行了分析。将土壤中的孔隙设想为各种孔径的圆形毛管，则土壤吸力S与毛管直径d的关系可表示如下：

S=4δ/d

(3)

式中:δ为水的表面张力系数(N·cm-1)，在室温条件为75×10-5N·cm-1。

如果土壤孔隙直径d的单位为mm，土壤吸力S的单位为Pa，则土壤孔隙直径d与土壤吸力S的关系式为d= 4δ/S，通过此式计算的土壤孔径为土壤当量孔径。雷志栋等[1]研究发现：如果土壤吸力为S1时对应的土壤含水量和土壤当量孔径为θ1和d1，表明土壤中小于或等于d1的孔隙中充水，其所占体积与孔隙总体积之比为θ1；如果土壤吸力为S2时对应的土壤含水量和土壤当量孔径为θ2和d2，表明土壤中小于或等于d2的孔隙中充水，其所占体积与孔隙总体积之比为θ2，则土壤中孔径在d1和d2之间的孔隙所占的体积与孔隙总体积之比为θ1-θ2。

不同深度土壤的孔隙分布状况可以反映出土层之间持水能力的变化，土壤当量孔径大小与土壤持水能力呈负相关，即土壤当量孔径越大，土壤的持水能力越弱[20]。根据雷志栋等[1]提出的将土壤水的数量和形态联系起来的特征含水量，将土壤有效含水率的基质势吸力区间定为330～15 000 cm压力水头，根据土壤当量孔径的计算公式(3)可计算得到相应的土壤当量孔径区间为0.000 2～0.009 mm。将土壤当量孔径大于0.009 mm、0.000 2～0.009 mm、小于0.000 2 mm分别定义为大孔隙、有效孔隙、微小孔隙，其孔隙内所含的水分分为重力水、有效水、无效水三种类型。其中，土壤重力水能够在重力作用下排走，土壤有效水能够被植物所利用，土壤无效水难以被植物利用。为了便于对研究区不同层位的土壤水分进行分析，将土壤当量孔径按土壤水吸力水头设定值分为田间持水吸力段(<330 cm)、低吸力段(330～1 000 cm)、中吸力段(1 000～15 000 cm)和高吸力段(>15 000 cm)，相应地将土壤当量孔径分为4个区间，即大于0.009 mm、0.003～0.009 mm、0.000 2～0.003 mm、小于0.000 2 mm。受试验条件限制，本文将吸力水头最高设定为5 000 cm，高于5 000 cm则采用拟合数据。研究区不同深度土壤当量孔径分布比例的对比，见图3。

图3 研究区不同深度土壤当量孔径分布比例的对比图

由图3可见，当土壤的当量孔径大于0.009 mm时，研究区不同层位由浅及深土壤的当量孔径所占的比例变化明显：在浅部0～0.2 m、0.2～0.5 m层位土壤的当量孔径所占的比例最低，可能由于浅部易受到人类活动影响，改变了土壤的部分性质，导致土壤孔隙发生改变；在2.0～2.5 m、3.0～3.5 m、3.5～4.0 m层位土壤的大孔隙所占的比例明显高于其他层位，表明此层位土壤相对疏松，且此层位与邻近层位存在孔隙差异，可能存在优先流；当土壤的当量孔径在0.003～0.009 mm范围内时，研究区整体呈现出0.5 m以上浅部层位土壤的当量孔径所占的比例小于0.5 m以下层位；当土壤的当量孔径在0.000 2～0.003 mm范围内时，研究区整体表现出2.0 m以上层位土壤的当量孔径所占的比例大于之下的层位，且随着层位加深土壤的当量孔径所占的比例降低，其中0.2～0.5 m、1.4～2.0 m层位土壤的当量孔径所占的比例分别高达49.75%、41.71%，此层位主要以有效孔隙为主，土壤孔隙内的水分可被植物吸收；当土壤的当量孔径小于0.000 2 mm时，研究区总体表现出层位由浅及深土壤的当量孔径所占的比例由小变大的趋势，其中2.5～3.0 m层位微小孔隙土壤的当量孔径所占的比例为38.88%。

综上分析可见，研究区2.0 m以上层位土壤主要以有效孔隙为主，2.0 m以下层位土壤的大孔隙较其他层位增多，同时微小孔隙所占的比例较大。浅部土壤大孔隙所占的比例较低，土壤含有的重力水相对较少，主要以毛细水为主，土壤水分向下传导缓慢，土壤的持水能力相对较好；土壤水分继续向下传导至2.0 m以下层位，由于孔隙变化，此层位土壤的大孔隙增多，土壤所含的重力水相对增多，水分向下传导相对增快，可能在此层位存在优先流。

依据土壤水分特征曲线的形态，通过土壤当量孔径的分析，可以简单地判断研究区不同深度土壤的分选情况，在一定程度上可反映出研究区第四系的形成条件。本文结合图2和图3分析发现，研究区0.2～0.5 m层位土壤的水分特征曲线较平缓，土壤当量孔径类别单一且占比高，其分选性最好。综合分析发现，研究区随着层位由浅及深，土壤的水分特征曲线整体逐渐变陡，土壤当量孔径类别复杂且占比不一，其分选性变差。

2.3 土壤水分特征曲线的拟合

本文选取van Genuchten模型和Gardner模型分别对研究区不同深度土壤样品的水分特征曲线进行了拟合，并对拟合效果进行评价。研究区不同深度土壤水分特征曲线的模型拟合结果见图4。土壤水分特征曲线是高度非线性的函数，本文采用非线性最小二乘法，并借助RETC 软件来拟合模型的各参数值，得到拟合模型的相关参数见表2。

图4 研究区不同深度土壤水分特征曲线的模型拟合结果

表2 研究区不同深度土壤水分特征曲线拟合模型的相关参数表

对比研究区不同深度土壤水分特征曲线拟合模型的相关参数，并结合相关参数的物理意义，通过分析表2可知，研究区各个层位土壤的饱和体积含水率(θs)分布区间为0.370 4～0.545 9 cm3/cm3，变化幅度较大，其中2.5～3.0 m层位土壤的饱和体积含水率最高，为0.545 9 cm3/cm3，0.5～0.9 m层位土壤的饱和体积含水率最低，为0.370 4 cm3/cm3；整体上看，1.4～2.0 m以上层位土壤的饱和体积含水率小于之下层位，与各层位的土壤质地、土壤孔隙分布情况不同有关。另外，由表1和图3可知，1.4～2.0 m以下层位土壤的孔隙比大，大孔隙所占的比例较大，土壤饱和状态下能够储存较多水分，其饱和体积含水率较高。参数α值表示的是土壤的进气值的倒数[21]。由表2可知，整体上看，研究区1.4～2.0 m以上层位土壤的α值小于之下层位，表明1.4～2.0 m以上层位土壤的进气值大于以下层位，同时结合研究区不同层位土壤当量孔径的分布比例分析发现，土壤进气值与土壤大孔隙所占的比例呈负相关，即土壤大孔隙所占的比例越高，则土壤进气值越小，这与实际土壤的水分特征曲线相符。

本文采用均方根误差(RMSE)和确定性系数(R2)评价两个模型的拟合效果，即R2值越大、RMSE值越小，表明模型的拟合精度越高，拟合效果越好。van Genuchten模型和Gardner模型拟合误差的结果见表2。

由表2可知，van Genuchten模型对研究区不同深度土壤水分特征曲线拟合结果的确定性系数(R2)值，除了0～0.2 m层位的土壤外，其余层位的土壤R2值均大于0.96，均值大于0.968，且均方根误差(RMSE)均控制在0.01左右，浮动小，表明该模型的拟合精度高；Gardner模型的R2均值为0.872，远小于van Genuchten模型的R2均值(0.968)，且其标准差为0.127，远大于van Genuchten模型的标准差0.025。单独对比Gardner模型对研究区不同深度土壤水分特征曲线的拟合效果，发现深部层位土壤的拟合效果整体好于浅部层位。此外，由图4可见，van Genuchten模型对研究区不同深度土壤水分特征曲线的拟合结果与实测曲线基本重合，很好地拟合出研究区各层位土体吸力与土壤体积含水量之间的关系，拟合得到的模型各参数值也符合其表达的实际物理意义，表明van Genuchten模型是研究区土壤水分特征曲线较为理想的拟合模型。

3 结 论

本文通过测定江汉平原北部黏土层不同深度土壤的水分特征曲线，对比分析了不同层位土壤水分特征曲线的变化规律以及土壤当量孔径的分布规律，并利用常用的van Genuchten模型和Gardner模型对研究区不同深度土壤的水分特征曲线进行了拟合，综合评价了模型的拟合效果，主要得出以下结论：

(1)研究区不同深度土壤的水分特征曲线表现出不同的形态，反映出土壤不同的持水能力，且1.4～2.0 m层位土壤的持水性能好于其他土层。在低吸力阶段，整体上研究区各个层位土壤快速排水，土壤水分特征曲线较缓；在高吸力阶段，整体上研究区土壤排水变缓，土壤水分特征曲线较陡，表明饱和态的土壤转为非饱和状态土壤的最初阶段，土壤基质吸力对土壤含水量的变化非常敏感，随着土壤饱和度的降低，其敏感性逐渐降低。

(2) 研究区不同深度土壤的水分特征曲线特征可反映土体孔隙分布的状况。浅层2.0 m以上层位的土壤主要以有效孔隙为主，连通性相对较差，水分向下传导相对缓慢，土壤的进气值较高，分选性较好；2.0 m以下层位的土壤大孔隙较其他层位增多，土壤在较小的吸力条件下就开始快速排水，土壤的进气值相对较低，失水速率较高，分选性较差，此层位土壤可能存在优先流。

(3) 选取van Genuchten模型和Gardner模型分别对研究区不同深度土壤水分特征曲线进行了拟合，获得了模型的相关参数，通过对比发现：2.0 m以下层位土壤的饱和体积含水率θs较高；参数α值表明2.0 m以上层位土壤的进气值大于以下层位。通过综合评价，结果表明van Genuchten模型对研究区不同深度土壤的水分特征曲线的拟合效果好、精度高，是研究区土壤水分特征曲线较为理想的拟合模型。