两种方法对土壤水分特征曲线的拟合及比较

梁晨璟+李春光+赵文娟

摘要：土壤水分特征曲线是土壤水吸力与含水率之间的关系曲线。此次试验利用压力膜仪测定银川北部盐渍土的土壤水吸力和含水率。对土壤水分特征曲线的Van-Genuchten模型采用Matlab软件和Microsoft Excel软件进行拟合。经过比较和分析后发现，用Microsoft Excel软件与Matlab软件拟合土壤水分特征曲线精度近似相同，但是Microsoft Excel软件对于不懂任何编程语言的人来说操作更为简便。

关键词：土壤水分特征曲线；Van-Genuchten模型；拟合方法；比较

中图分类号：S152.7+1 文献标识码：A 文章编号：0439-8114（2014）01-0056-03

Comparative Analysis of Two Methods Fitting Soil Moisture Characteristic Curve

LIANG Chen-jing1，LI Chun-guang2，ZHAO Wen-juan1

（1. Civil Engineering and Water Conservancy， Ningxia University， Yinchuan 750021，China； 2. Institute of Numerical Computation and Engineering Applications， Beifang University of Nationalities，Yinchuan 750021，China）

Abstract： Soil moisture characteristic curve describes the relationship between water suction and moisture content in soil. This paper is aimed to measure the suction and moisture content by using pressure membrane analyzer in the northern area of Yin Chuan saline soil. Van - Genuchten model is used to fit soil moisture characteristic curve use by Matlab and Microsoft Excel. By comparative analysis， it was found that it had approximately the same accuracy of the soil moisture characteristic curve fitted by using Microsoft Excel and Matlab software. But it is more simple and convenient for people not knowing any programming language by using Microsoft Excel software.

Key words： soil moisture characteristic curve； Van-Genuchten model；fitting method； comparison

收稿日期：2013-03-29

基金项目：国家自然科学基金重大研究计划培养项目（91230111）

作者简介：梁晨璟（1988-），女，河北邢台人，在读硕士研究生，研究方向为旱区河流泥沙动力学理论与数值模拟，（电话）15909510629（电子信箱）

liangmeng88-12@163.com；通讯作者，李春光（1964-），男，教授，博士，主要从事计算机数学和流体力学研究，（电子信箱）

cglizd@hotmail.com。

土壤水分特征曲线是土壤水的基质势或土壤水吸力随土壤含水率变化的关系曲线，是研究土壤水分的保持和运动所用到的反映土壤水分基本特征的曲线。土壤水分特征曲线表示了土壤的一个基本特性，有重要的实用价值：①可以利用它进行土壤水吸力和含水率之间的换算；②它可以间接地反映出土壤中空隙大小的分布；③它可以用来分析不同质地土壤的持水性和土壤水分的有效性；④应用数学物理方法对土壤中的水分运动进行定量分析[1]。因而，对土壤水分特征曲线的研究是土壤物理学工作者关注的重点问题。

1 材料与方法

此次试验是针对宁夏银川北部地区的盐渍土进行研究，因此选取银川市西夏区双渠口村、银川市永宁县望远乡望远三队、石嘴山市大武口区和石嘴山市平罗县4个地方进行采样并依次编号为1、2、3、4号。经测定它们的土壤物理参数如表1。

测定土壤水分特征曲线的方法有很多种[2]。此次试验采用的是直接测定法中的压力膜仪法。压力膜仪法与其他方法相比，不仅可完整地测定出整个土壤水分特征曲线，而且可测定原装土或是扰动土。

1.1 研究选取的模型

由于土壤水分特征曲线受土壤质地、土壤结构、土壤容重、温度及土壤的膨胀和收缩等影响，所以很难定性拟合土壤水分特征曲线。为了计算和分析方便，常把水分特征曲线近似处理为经验公式和理论模型，如目前比较常用的经验公式有Brooks-Corey模型[3]、Gardner模型[4，5]、Van-Genuchten模型[6]和Gardner-Russo模型[7]。

此次的曲线拟合采用Van-Genuchten模型。Van-Genuchten模型既连续又有连续斜率，得到的曲线光滑，对绝大多数土壤在相当宽的水势或含水量范围内具有普遍适用性，并可得到相对导水率的解析解，因而在土壤物理领域得到了最为广泛的应用[8]，其公式表示为：

■=■■ （1）

式中，θ为体积含水率（cm3/cm3），θr为残留含水率（cm3/cm3），θs为饱和含水率（cm3/cm3），h为土壤水吸力（cm），m，n，α为拟合参数，其中m=1-（1/n）。

1.2 拟合参数初值的确定

在拟合时，Van-Genuchten模型中含水率θ和土壤水吸力h为试验实测数值，而残留含水率θr，饱和含水率θs，α，m，n在拟合过程中都视为非线性曲线拟合参数。无论用什么方法都需要对这些拟合参数设置初值。

参数θr为残留含水率，它的物理意义为土壤水分特征曲线导数为零时的土壤含水率。在实际的应用中，一般认为凋萎点处的含水率为残留含水率，即在土壤水吸力为1 500 cm时的含水率可视为凋萎含水率。参数θs为饱和含水率，也就是说在土壤水吸力等于零时的含水率。θr和θs可以通过试验测得，一般土壤的取值范围为0<θr<θs<0.5 cm3/cm3。参数n没有明显的物理意义但是决定了土壤水分特征曲线的坡度。当n大时曲线缓，当n小时曲线陡。参数α一般认为是进气值。吸力小时土壤无水排出，土壤含水率维持饱和；当吸力超过某一临界值，土壤孔隙中水分开始外排，这一临界负压值称为进气值。重质黏性土壤进气值较大，轻质土或结构良好土壤进气值较小[9]。

1.3 利用Matlab软件拟合土壤水分曲线

Matlab软件是用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言。在Matlab中内置许多命令可以直接调用来进行数据处理，其中命令lsqcurvefit可用于非线性曲线拟合。命令lsqcurvefit的原理是利用最小二乘法，即∑■■[F（x，xdatai）-ydatai]2最小，式中，xdata和ydata为向量，F（x，xdata）为向量函数值，x为拟合参数向量。根据实测得到的数据xdata和ydata，找到与函数F（x，xdata）最佳拟合参数x[10]。

想要使用命令lsqcurvefit，事先应定义M-函数文件，而后在脚本函数中调用命令lsqcurvefit。命令lsqcurvefit用法为a=lsqcurvefit[@fun，a0，xdata，ydata]，式中，a为非线性曲线拟合参数值，fun为定义的文件名称，a0为非线性曲线拟合参数初值，xdata，ydata为试验实测数值。

1.4 利用Microsoft Excel 2007拟合土壤水分特征曲线

在Microsoft Excel 2007中内置有“规划求解”这一功能，可利用这一方法进行土壤水分特征曲线的拟合。Microsoft Excel 2007的“规划求解”为非线性最优化代码。Microsoft Excel 2007 Solver程序代码是以宏的方式提供调用的。使用时不需要关心其具体的实现方法，只需要和它的对话框进行交互即可[11]。

使用Microsoft Excel 2007具体方法如下：

1）在Microsoft Excel 2007中的一列输入测得的土壤水吸力值，在另一列输入含水率，同时把残留含水率θr，饱和含水率θs，α，m，n在程序中都视为参数，单写一行并假定初值。

2）利用Microsoft Excel 2007的计算功能编写Van-Genuchten模型关系式，根据现有的数据计算出含水率，并求得残差平方和。

3）利用Microsoft Excel 2007的“规划求解”功能计算模型中的参数。“设置目标单元格”为所算出的残差平方和并使它等于最小值。“可变单元格”选择残留含水率θr，饱和含水率θs，α，m，n所设的初值，“约束”为模型中的拟合参数并使其设定为非负。

4）在右侧的“选项（Q）”可设置最长运算时间、迭代次数、精度等有关具体运算的相关信息，待选择好之后确定。选择“求解（s）”，便可求出拟合曲线的拟合参数。

5）根据上述计算结果，利用Microsoft Excel 2007的画图功能便可画出Van-Genuchten模型的土壤水分特征曲线。

2 Matlab 与Microsoft Excel 2007拟合结果比较

图1至图4、表2为1～4号土壤的水分特征曲线拟合的结果。

从上述实例可以发现，利用不同的软件拟合的结果是不同的。这主要归结于它们各自算法的不同。用Microsoft Excel 2007软件进行曲线拟合解决了大多数人对程序语言不了解的问题。操作人员仅需要了解Microsoft Excel 2007表格公式的计算方法和对“规划求解”对话框内容的填充。待所有计算完毕，计算结果又会表现在Microsoft Excel 2007表格内。这样使得操作人员可以清晰地看到拟合点的具体数值。但是在计算过程中会发现，用Microsoft Excel 2007拟合土壤水分特征曲线时，有时会将残留含水率θr计算为零，这显然是不对的。出现这种情况的原因是在设置约束条件时，发现约束条件只能设置5种情况，即大于等于、小于等于、等于、取整和二进制。但是在Van-Genuchten模型中要求拟合参数都大于零，这样Microsoft Excel 2007的约束条件便无法完全满足计算需要。因而用Microsoft Excel 2007只能适用大多数情况。除去上述缺点，Microsoft Excel 2007在土壤水的低吸力和高吸力下都能很好地拟合出曲线。

利用Matlab软件拟合土壤特征曲线时，需要操作人对编程语言了解且能熟练应用。但是经上述实例发现，Matlab软件只是在土壤水的低吸力状态下拟合很好。通过与Microsoft Excel 2007软件比较各自所算得的残差平方和和相关系数，可以发现用Microsoft Excel 2007软件拟合曲线与Matlab软件拟合曲线精度相近。

3 结论

无论用哪种软件拟合曲线，都要注意初值的设定。Matlab与Microsoft Excel 2007相比，两种软件各有利弊。总的来说虽然用Microsoft Excel 2007软件拟合曲线的精度近似于用Matlab软件拟合曲线的精度，且Microsoft Excel 2007的操作不需要有任何语言基础，但是Microsoft Excel 2007软件在拟合上操作较繁琐，在设置约束条件上有局限性，有可能使计算结果为零。而Matlab对懂编程语言的人来说更为简单，仅需调用非线性曲线拟合命令lsqcurvefit便可精确拟合出曲线。

参考文献：

[1] 雷志栋，杨师秀，谢森传. 土壤水动力学[M].北京：清华大学出版社，1988.

[2] 徐绍辉，刘建立.土壤水力性质确定方法研究进展[J].水科学进展，2003，14（4）：494-501.

[3] MILLY P C D. Estimation of the Brooks-Corey parameters from water retention data[J].Water Resour，1987，23：1085-1089.

[4] GARDNER W R， HILLEL D， BENYAMIN Y. Post irrigation of soil water I. Redistribution[J]. Water Resour，1970，6：851-861.

[5] GARDNER W R，HILLEL D， BENYAMIN Y. Post irrigation of soil water II. Simultaneous redistribution and evaporation[J]. Water Resour，1970，6：1148-1153.

[6] VAN GENUCHTEN M TH. A closed-from equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Sci Soc Am，1980，44：892-898.

[7] RUSSO D. Determining soil hydraulic properties by parameter estimation： on the selection of a model for the hydraulic properties[J].Water Resour，1988，24：453-459.

[8] 朱蔚利，肖自幸，牛健植，等. 两种模型对土壤水分特征曲线拟合的比较分析[J].湖南农业科学，2011（17）：47-51.

[9] 李法虎，傅建平，孙雪峰. 土壤水分运动参数的确定及其灵敏性能分析[J].灌溉排水，1993，12（2）：6-14.

[10] 彭建平，邵爱军.基于Matlab方法确定VG模型参数[J].水文地质工程地质，2006，33（6）：25-28.

[11] 余 亮. 利用Microsoft Excel软件进行非线性拟合的非编程方法[J].微型机与应用，2000（5）：16-17.

（责任编辑 郑 威）

3 结论

无论用哪种软件拟合曲线，都要注意初值的设定。Matlab与Microsoft Excel 2007相比，两种软件各有利弊。总的来说虽然用Microsoft Excel 2007软件拟合曲线的精度近似于用Matlab软件拟合曲线的精度，且Microsoft Excel 2007的操作不需要有任何语言基础，但是Microsoft Excel 2007软件在拟合上操作较繁琐，在设置约束条件上有局限性，有可能使计算结果为零。而Matlab对懂编程语言的人来说更为简单，仅需调用非线性曲线拟合命令lsqcurvefit便可精确拟合出曲线。

参考文献：

[1] 雷志栋，杨师秀，谢森传. 土壤水动力学[M].北京：清华大学出版社，1988.

[2] 徐绍辉，刘建立.土壤水力性质确定方法研究进展[J].水科学进展，2003，14（4）：494-501.

[3] MILLY P C D. Estimation of the Brooks-Corey parameters from water retention data[J].Water Resour，1987，23：1085-1089.

[4] GARDNER W R， HILLEL D， BENYAMIN Y. Post irrigation of soil water I. Redistribution[J]. Water Resour，1970，6：851-861.

[5] GARDNER W R，HILLEL D， BENYAMIN Y. Post irrigation of soil water II. Simultaneous redistribution and evaporation[J]. Water Resour，1970，6：1148-1153.

[6] VAN GENUCHTEN M TH. A closed-from equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Sci Soc Am，1980，44：892-898.

[7] RUSSO D. Determining soil hydraulic properties by parameter estimation： on the selection of a model for the hydraulic properties[J].Water Resour，1988，24：453-459.

[8] 朱蔚利，肖自幸，牛健植，等. 两种模型对土壤水分特征曲线拟合的比较分析[J].湖南农业科学，2011（17）：47-51.

[9] 李法虎，傅建平，孙雪峰. 土壤水分运动参数的确定及其灵敏性能分析[J].灌溉排水，1993，12（2）：6-14.

[10] 彭建平，邵爱军.基于Matlab方法确定VG模型参数[J].水文地质工程地质，2006，33（6）：25-28.

[11] 余 亮. 利用Microsoft Excel软件进行非线性拟合的非编程方法[J].微型机与应用，2000（5）：16-17.

（责任编辑 郑 威）

3 结论

无论用哪种软件拟合曲线，都要注意初值的设定。Matlab与Microsoft Excel 2007相比，两种软件各有利弊。总的来说虽然用Microsoft Excel 2007软件拟合曲线的精度近似于用Matlab软件拟合曲线的精度，且Microsoft Excel 2007的操作不需要有任何语言基础，但是Microsoft Excel 2007软件在拟合上操作较繁琐，在设置约束条件上有局限性，有可能使计算结果为零。而Matlab对懂编程语言的人来说更为简单，仅需调用非线性曲线拟合命令lsqcurvefit便可精确拟合出曲线。

参考文献：

[1] 雷志栋，杨师秀，谢森传. 土壤水动力学[M].北京：清华大学出版社，1988.

[2] 徐绍辉，刘建立.土壤水力性质确定方法研究进展[J].水科学进展，2003，14（4）：494-501.

[3] MILLY P C D. Estimation of the Brooks-Corey parameters from water retention data[J].Water Resour，1987，23：1085-1089.

[4] GARDNER W R， HILLEL D， BENYAMIN Y. Post irrigation of soil water I. Redistribution[J]. Water Resour，1970，6：851-861.

[5] GARDNER W R，HILLEL D， BENYAMIN Y. Post irrigation of soil water II. Simultaneous redistribution and evaporation[J]. Water Resour，1970，6：1148-1153.

[6] VAN GENUCHTEN M TH. A closed-from equation for predicting the hydraulic conductivity of unsaturated soils[J]. Soil Sci Soc Am，1980，44：892-898.

[7] RUSSO D. Determining soil hydraulic properties by parameter estimation： on the selection of a model for the hydraulic properties[J].Water Resour，1988，24：453-459.

[8] 朱蔚利，肖自幸，牛健植，等. 两种模型对土壤水分特征曲线拟合的比较分析[J].湖南农业科学，2011（17）：47-51.

[9] 李法虎，傅建平，孙雪峰. 土壤水分运动参数的确定及其灵敏性能分析[J].灌溉排水，1993，12（2）：6-14.

[10] 彭建平，邵爱军.基于Matlab方法确定VG模型参数[J].水文地质工程地质，2006，33（6）：25-28.

[11] 余 亮. 利用Microsoft Excel软件进行非线性拟合的非编程方法[J].微型机与应用，2000（5）：16-17.

（责任编辑 郑 威）