第二章 土壤水分运动基本方程2

1、第二章 土壤水分运动基本方程如前所述，达西定律是由达西（Darcy，Henry 1856）通过饱和砂柱渗透试验得出，后由Richards（1931）将其扩伸至非饱和水流中，并规定导水率为土壤负压h的函数，即 （2-2-1）式中：为水势梯度； k（h）为导水率，是土壤负压h的函数； q 为水流通量或流速。Richards方程垂向一维方程为注意：H=hz，垂直坐标向上为“+”；向下时为“”。由于k（h）受滞后影响较大，上式仅适用于单纯的吸湿或脱湿过程。若将导水率作为容积含水率函数，即以k（）代替人k（h），则可避免滞后作用的影响。一般说来达西定律对饱和与非饱和水流均可适用，即水流通量与势能梯度成正

2、比。但在饱和土壤中，压力为正值，其总水头包括了由该点在地下水面以下深度来确定的静水压力（正值）和相对于基准面高度来确定的位置水头，总水头为压力水头和位置水头之和，水由总水头高处向低处流动。在非饱和土壤中，基质势为负值，土水势在不考虑溶质势、温度势及气压势时，只包括重力势和基质势。因此，总水头常以负压水头和位置水头之和来表示。一维Richards方程的几种形式：根据（K=CD）得： 第一节 直角坐标系中土壤水分运动基本方程一、基本方程的推导土壤水分运动一般遵循达西定律，且符合质量守恒的连续性原理。土壤水分运动基本方程可通过达西定律和连续方程进行推导。如图2-2-1所示，从土壤中取出微分单元体ab

3、cdefgh，其体积为，由于该立方体很小，在各个面上的每一点流速可以看成是相等的，设其流速为，在tt+t时段内，流入立方体的质量为（3个面流入）： （2-2-2）流出立方体的质量为（3个面流出）： （2-2-3）式中：水的密度；分别表示微分体x、y、z方向长度；，分别表示水流经微分体后，其流速在x、y、z方向的变化值。由式(2一22）、式（223）之差可求得流入和流出立方体的质量差： （224）设为立方体内土壤含水率，则在t时间内立方体内质量变化又可写为 （225）根据质量平衡原理（流入量流出量储存量变化量），式（324）、式（325）应相等，即 （2-2-6）根据达西定律得：， （2-2-7

4、）式中k（）土壤水力传导度，为含水率的函数；H总土水势，为基质势与重力势之和（Hhz）。因此，式（226）可以写作以下形式： （2-2-8）上式可以简写为 （2-2-9）式（228）或式（229）为土壤水分运动基本方程。在饱和土壤中，含水量和基质势均为常量。水力传导度也为常量，常称渗透系数，则方程（228）可写为 （2210）或写作 （2210） （2211）式中：2拉普拉斯算子。式（2210）或式（2210）为饱和土壤水流的拉普拉斯方程。二、基本方程的不同形式为运用基本方程分析各种实际问题的方便，可将基本方程改写为多种表达形式。为简便起见，以下均以一维垂向土壤水分运动为例，给出基本方程的不同

5、表达形式。（一）以含水率为变量的基本方程由式（228）可得一维垂向土壤水分运动的基本方程为 （2212）式中：H总土水势； z为水流方向坐标，取z向上为正。因为H=h十z，所以上式可写作 （2213）式(2213）为以为变量的基本方程，将代入式（2213）得：令，则式（2213）可以写成（一维垂向土壤水分运动方程）： （2214）在水平运动的情况下，重力项等于0，所以，其形式与Fick扩散定律相同。式（2214）具有扩散方程的形式，故将D（）称为扩散度。 （2214）Fick定律：自由水中溶质的分子扩散通量符合Fick定律：式中：J为溶质的扩散通量； D为溶质的扩散系数； 为溶质的浓度梯度。（

6、二）以基质势h为变量的基本方程由于 ，则式（2214）可以写成： （2215）式中：c（h）比水容量（也称容水度），c（h）=，表示单位基质势变化时含水率变化。（三）以参数v为因变量的基本方程采用Kirchhoff变换，令则 由式（2-2-15）得： （2216）式中hc土壤的进气值，即土壤含水率开始小于饱和含水率时的负压值。另外，；在非饱和区： 在饱和区： 且因为 ，所以 ；则方程式（2216）为：（四）以位置坐标z为变量的土壤水运动方程以z为变量，则z为、t的函数，z（，t）为未知函数。已知=（z，t），当处，可以解出z= z（，t），即14对z，t分别求导数：，于是 及将以上式子代入方程

7、（2214）得： （2217）（五）以参数u为因变量的土壤水运动方程定义式中：初始含水率； ； 饱和含水率。由式（2-2-14）得： 将代入上式得： 所以 （2-2-18）以上各式中式（2214）、式（2215）是二种经常采用的形式，形式的选定取决于要解决问题的边界条件和初始条件。以含水率为因变量的基本方程常用于求解均质土层或全剖面为非饱和流动问题，这种方程形式对于层状土壤或求解饱和非饱和流问题不适用；以负压水头h为因变量的基本方程是应用较多的一种形式，可适用于饱和非饱和水流求解及层状土壤的水分运动分析计算，但由于非饱和土壤水的导水率k（h）及容水度c(h)，受滞后影响较大，计算中参数选取不当

8、会造成较大误差；以v，u为因变量基本方程实际上分别相当于以负压水头h和含水率为因变量的基本方程，在某些情况下由于经代换后方程较为简单，易于求解；以坐标为因变量的基本方程根据定解条件需要求解较简单的土壤水分运动问题。以上为直角坐标系中土壤水分运动的基本方程，求解某些土壤水分运动问题时，采用柱坐标系可能更方便。第二节 柱坐标系中土壤水分运动基本方程在推导柱坐标系中的基本方程时，方法同直角坐标系，同样可用达西定律与连续方程相结合的方法导出。若以z轴为轴的柱坐标系，根据达西定律，在此坐标系中可表示为：式中：r、z分别表示柱半径，角坐标和垂直坐标（图222） q r、q、qz相应于r、z三个方向的通量；

9、 H总水势。下面利用质量守恒来推导连续方程。t时段内，在r方向的流入量为zt，流出量为，则流入与流出量之差（忽略高阶无穷小量）为 （2-2-19）同理，在方向流入流出量之差为 （2-2-20）在z方向土壤水分流入流出量之差为 （2-2-21）上述三个方向流入和流出单元体的水量差总计为 （2-2-22）单元体体积应为 ，略去高阶无穷小量后为rzt，在t时间内单元体内水分增量为 （2-2-23）根据质量守恒原理武（3222）应与式（3223）相等，即 （2-2-24）式（2224）为柱坐标系中土壤水分运动的连续方程。将式（2218）代入上式，即得柱坐标系中土壤水分运动基本方程： （2-2-25）以

10、总水势H=h+z，水容度c=，以及导水率k（），扩散度D（）等代入，基本方程可表示为 （2-2-26）对于平面轴时称问题，上式可改写为 （2-2-26）同理可推得以x（或y）轴为轴的柱坐标系的基本方程： （2-2-27）关于球坐标系中的基本方程应用较少，推导方法同上，这里不再论述。第三节 土壤水分运动基本方程的定解条件土壤水运动基本方程的定解条件包括初始条件和边界条件，为了简单起见，将以直角坐标系中基本方程常用形式为例进行论述。（一）初始条件相应于前式（1）、式（15）的初始条件分别以下式表达： （2-2-28） （2-2-29）脚标“i”表示初始已知量。初始条件：t=0时剖面上、h的分布已知

11、。（二）边界条件边界条件一般有一类边界、二类边界、三类边界三种。1一类边界条件（变量已知的边界1）对干式（2214）、式（2215）的一类边界的表达式为 （2-2-30） （2-2-31）脚标“0”均表示一类边界上的值；z0为一类边界的坐标。在一维垂向土壤水分运动中，一类边界的情况发生在压力入渗（地表形成水层）时，地表含水率达到饱和含水率，或当强烈蒸发时，表土达到风干土含水率的情况。2二类边界条件（边界2上水流通量已知的情况）相应于式（2214）、式（2215）表达式为 （2-2-32） （2-2-33）式（2232）及式（2233）中，均假设垂直坐标z向下为正。在一维垂向土壤水分运动中，这种

12、情况常发生在降雨、灌水入渗或蒸发强度已知的边界。在降雨或灌水入渗时，（t）为正值，在蒸发时（t）为负值。在不透水边界和无蒸发入渗的边界，（t）=0，则式（2232）、式（2233）分别为 （2-2-32，） （2-2-33，）3三类运界条件相当于水流通量随边界3上的变量（或h）值而变化的情况三类边界的一般形式为 （2-2-34）式中，f为变量。在土壤蒸发强度为表土含水率或表土负压的函数的情况下，式（2214）、式（2215）的三类边界条件表达式为 （2-2-35） （2-2-36）方程式（2235）右端的a+b表示三类边界上水流通量为表土含水率的线性函数。方程式（2236）右端的a(h)+b表

13、示三类边界上水流通量为表土负压的函数。例如、土壤的下部有弱透水层阻隔，边界受顶托补给，其补给强度决定于下部弱透水层的导水率k2，弱透水层上、下的压力h1、h2，其厚度为，方程的三类边界条件可写成： （2-2-36，）上述的二种边界条件是经常遇到的情况。在野外实际情况下，有时还存在地下水位为已知的运动边界，此时可将地下水面处h=0作为边界条件。如在任一时间，地下水埋深为d（t），则 （2-2-37）d（t）表示 t 时刻的地下水面所在位置，如地下水位是等速下降的，则 （2-2-38）式中：v地下水位下降的速度。如地下水位是由于水井抽水引起下降的，则 （2-2-39）式中：d0下水初始埋深； Q井

14、的抽水流量； T地下水含水层的导水系数； a地下水含水层的导压系数； 井函数（见泰斯公式）； r计算点离抽水井若地下水位变化规律未知，不能作动边界处理。第四节 土壤水分运动参数土壤水分运动中的主要参数有土壤水力传导度k（又称导水率），比水容量c（也称容水度）及扩散度D（也称扩散率）等。这些参数的变化决定了土壤水分运动状态，所以了解和掌握这些参数的特性及其变化规律是十分重要的。一、土壤水力传导度k土壤水力传导度是反映土壤水分在压力水头差作用下流动的性能。一般在饱和土壤中导水率称为渗透系数。土壤水力传导度为在单位水头差作用下，单位断面面积上流过的水流通量。它是土壤含水率或土壤负压的函数。饱和土壤孔

15、隙中都充满水，导水率达到最大值，且为常量。在非饱和土壤中，因土壤孔隙中部分充气，导水孔隙相应减少，导水率低于饱和土壤水情况，而且导水率是负压或含水率的函数，随着含水率降低而减小。由于在吸力作用下，土壤水首先从大孔隙中排除，随着吸力增加，水流仅能在小孔隙中流动。所以，土壤从饱和到非饱和将引起导水率的急剧降低。当吸力由零增至1105 Pa时，导水率可能降低好几个数量级，有时降低到饱和导水率的1100000。对于不同结构土壤，饱和与非饱和土壤水导水性能的相对关系是不同的。饱和土壤导水性能最好的是粗粒砂性土壤，导水最差的土壤是细质粘土，但非饱和土壤在较大负压情况下则情况可能相反。具有大孔隙粗质土壤，在

16、吸力作用下孔隙中水分很快排除，导水率迅速下降；而粘质细颗粒土壤，在较高吸力下，许多小孔隙仍充满水，仍具有一定的导水性能，因此，导水率下降较缓慢。所以，细颗粒粘质土壤在同一吸力条件下可能较大孔隙粗质土壤具有较高导水率15。所以，导水率与含水率（或负压）关系是较复杂的，目前还不能用理论分析方法推导它们之间关系式，需通过试验测定。图223为非饱和土壤水在负压水头作用下流动的模型。在水平土柱两端有多孔板，分别由平水箱保持一定水位，使其负压为h1和h2，在负压梯度h/x的作用下，立柱中土壤水从1端向2端运移。土壤水通量q可由1端补给量或2端溢出量测得，两者相等时，水流处于稳定状态。非饱和土壤水力传导度可

17、由达西定律求得。由于水平土柱沿程负压（或含水率）是变化的，求得的导水率k也应是变化的，若距离较小，可以平均负压（或含水率）确定平均土壤水力传导度。在不同的平均负压（吸力）值下，通量与负压梯度成正比，如图224所示15。两者呈直线关系，但其斜率（即水力传导度）随平均负压而变。此外，土壤水力传导度还与土质有关，如图225所示，砂性土壤饱和导水率高于粘性土壤，随着土壤吸力增加，砂性土壤导水率降低速率较粘性土壤快，所以吸力增大时，粘质土壤导水率反大于砂质土壤。非饱和土壤水力传导度k与土壤负压h或含水率的关系通常由试验资料拟合成经验公式，一般有以下几种形式。（1）土壤水力传导度与负压（吸力）h的关系式：

18、 （2-2-40） (2-2-41) (2-2-42)式中：ks饱和土壤水导水率，或称渗透系数； a，b经验常数； c，n经验指数。（2）土壤水力传导度与土壤含水率的关系式： (2-2-43) (2-2-44) (2-2-45) (2-2-46)式中： r某一特征含水率，通常采用最大分子持水率； s饱和含水率； c，m，n均为经验指数，在式（2246）中， 其他符号意义同前。k值的量纲单位为LT-1。由于土壤负压与含水率的关系曲线土壤水分特征曲线有滞后现象，所以，土壤水力传导度随负压的变化同样也存在滞后现象，即在同一负压下，干燥过程中的土壤水力传导度高于湿润过程中的土壤水力传导度。但土壤水力传

19、导度与含水率关系受滞后作用的影响较小。二、土壤水分扩散度D土壤水扩散度为单位含水率梯度下，通过单位面积的土壤水流量，其值为土壤含水率的函数，即。扩散度与土壤含水率的关系如图226所示，这种关系有时可用以下经验公式表示： (2-2-47)式中：a，b经验常数，D的量纲为L2T-1。以上公式仅能用于含水率较高的阶段。在土壤含水率很低时，由于土壤水汽扩散速度增大，使扩散度随土壤含水率降低而增大。在土壤含水率很高的情况下，土壤接近于饱和，扩散率趋向于无限大。三、容水度c容水度表示压力水头减小一个单位时，自单位体积土壤中所能释放出来的水体体积，量纲为L。容水度可以用下式表示：它是负压的函数，为水分特征曲

20、线上任一特定含水率值时的斜率（导数），并随土壤水分特征曲线而变化，所以它取决于土壤含水率和土壤质地等。第五节 考虑水汽热耦合关系的土壤水分运动基本方程长期以来，人们都采用等温模型研究土壤水分运动。在自然条件下，日夜温差很大，地表以下不同深度处温度的差异和变化影响土壤水分的转化和运移，用等温模型来研究土壤水分运动常带来一定误差。一些学者根据能量平衡和热传导理论，提出了考虑水、汽、热耦合关系的土壤水运动基本方程。根据Philip、De Vries和Milly的理论，多孔介质中的水分质量通量可以表示为 (2-2-48)式中：水分总质量通量（kg/s）； 水流质量通量（kg/s）； 水汽质量通量(kg

21、/s)。在入渗速率不很大的情况下，土壤中水和汽之间存在局部热动力学平衡时，及可表示为 (2-2-49) (2-2-50)式中：土壤水密度（kgcm3）； h土壤水压力水头（cmH2O，1cmH2O=98103Pa）； T绝对温度（C）； DTV热蒸汽扩散系数； DHV蒸汽传导率； 水蒸汽气体常数，R4615erg（g）（1erg=10J）。在水汽和多孔介质中水体达到局部平衡时，两者之间的自由能相等，则 （2-2-51）式中：v绝对湿度（gcm3）； 0（T）水汽饱和状态下的绝对湿度（gcm3），根据Camilo等19研究表明： (2-2-52)Constantz等对水力传导度与温度的关系进行了

22、研究，提出了水力传导度可表示为 (2-2-53)若忽略温度变化对流体密度影响，则k（h，T）可表示为 (2-2-54)以上各式中：k多孔介质内渗透率（cm2）；kr(h）相对非饱和渗透率；g重力加速度（cms2）；k(h, T0）在参考温度T。时的水力传导率（cm/d）流体动力粘滞系数。式（2254）中（T0）（T）可表示为 (2-2-55)根据Philip和De Vries等研究，热蒸汽扩散系数DTV和蒸汽传导率DHV分别可表示为， (2-2-56) (2-2-57) (2-2-58) (2-2-59) (2-2-60)式中：f在Philip和 De Vries公式中对水汽扩散引入的修正因子

23、，f=f，经验系数=1332；k液体水流流动盯以忽略时的含水率；Da空气中分于扩散系数（cm2/s）；由气体所充填孔隙的弯曲率；L孔隙中蒸汽的体积与土壤体积比；H为相对湿度，H=exp（hgRT）。根据质量守恒定律得； (2-2-61)将式（2249）、式（2250）、式（2251）代人上式，得到水分运动方程： (2-2-62)其中 在单位体积多孔介质中的热量为 (2-2-63)式中：C=为水汽介质的热容量（calcm3，lcal41868J）；Ci，i土壤中水、汽、石英和其他矿物质及有机质的体积热容量及占土体百分数；L0在参考温度T下的蒸发潜热，calg，在T=20时，L0=585calg。

24、蒸发潜热L可以下式表示：式中：Cv常压下水蒸汽比热；Ci水的比热。W为微分吸湿热（calg），由热动力学原理，可表示为17： (2-2-64)式中：j热功当量（ergcal）。多孔介质内任一点处的热通量为 (2-2-65) (2-2-66)式中：土壤的热传导系数cal（cms），它是通过土壤中各种组分的热传导系数加权平均而求得的。 ki固体颗粒温度梯度与水体温度梯度之比。Chung20提出了一个热传导系数的经验式 (2-2-66,)将式（2263）及式（2265）代入热量连续方程得： （2-2-67）其中 C为湿土壤的热容量cal/（cm3）。式（2262）和式（2267）及其相应定解条件组成

25、了水、汽、热耦合求解模型。对一维问题，方程可简化为： （2-2-68） （2-2-69）由于方程中系数也为待求函数的函数，这些非线性方程需通过迭代求解。在求解过程中，常用的压力和热量（能量）单位有以下几种。（1）压力单位常采用帕（Pa）、巴（bar）、标准大气压（atm）、毫米汞柱（mmHg）、毫个水柱（mmH2O）等。（2）热量（能量）常用单位有焦尔（J）、卡路里（cal），尔格（erg）。以上单位的相互转换关系见下表所示。第六节 土壤水分通量法一、土壤水分通量法基本原理土壤水分通量法是直接利用达西定律和质量守恒原理分析计算土壤水通量及潜水的入渗量或蒸发量的一种方法，这种方法在有定位点负压h

26、（基膜势）和含水率观测资料地区应用简便。根据质量守恒原理，一维垂向土壤水流连续方程（见2-2-6）可写作： （2-2-70）上式由z1至z2积分得： （2-2-71）式中：q（z1）、q（z2）分别表示高度为z1和z2处土壤水运动通量。在t1至t2时段内（t），上式可写作： （2-2-72）式中：Q（z1）、Q（z2）分别表示高度为z1和z2处通过单位断面面积的水量（t1至t2时段内）。式（2272）表明，当已知时段前后两个瞬时土壤剖面上含水率分布时，仅需已知一断面上土壤水通量即可求得任一断面的通量或水量。因此，称该方法为土壤水通量法。由于这种方法是根据时段前后两个瞬时含水率剖面确定水流通量和

27、水量的，在某些情况下，称之为瞬时剖面法。通量法可分为零通量面法和已知通量法两种。二、零通量面法由前述可知，当不考虑溶质势、气压势和温度势时，土壤水的总水势为基膜势m及重力势g之和，常用负压水头h和位置水头z之和表示。在测定土壤剖面上基质势和重力势后，可计算出土壤剖面上总水势分布曲线，如图227所示。由达西定律，土壤水通量为：q，当0时，q0，即为零通量面。图2中A、B两断面均为零通量面，但A、B两断面的水流状况是不同的。对A断面而言，水流向上、下断面运移，也可称为发散型零通量面。而B断面的上、下断面水流分别向断面B汇集，故称聚合型零通量面。图227中有两个零通量面，这表明降雨（或灌水）入渗与蒸发是交替发生的。一般土壤较长时间处于单一的蒸发或入渗状态，剖面上可能不存在零通量面。当剖面上存在零通量面时，可根据式（2272）求任一断面z处土壤水通量，即： （2