2_非饱和水流运动基本方程

1、非饱和土壤水流的基本方程非饱和土壤水流的基本方程1 非饱和土壤水流的达西定律非饱和土壤水流的达西定律2 土壤水分运动参数及其测定方法土壤水分运动参数及其测定方法3 非饱和土壤水运动的基本方程非饱和土壤水运动的基本方程4土壤水分运动基本方程的定解条件土壤水分运动基本方程的定解条件1.非饱和土壤水流的达西定律非饱和土壤水流的达西定律非饱和土壤水分运动和饱和土壤水分运动一样，水分从水势高非饱和土壤水分运动和饱和土壤水分运动一样，水分从水势高处向水势低处运动。一般认为，适用于饱和水流动的达西定律处向水势低处运动。一般认为，适用于饱和水流动的达西定律在很多情况下也同样适用于非饱和土壤水分流动。在很多情况

2、下也同样适用于非饱和土壤水分流动。19311931年，年，RichardsRichards最早将达西定律引入非饱和土壤水流动。非最早将达西定律引入非饱和土壤水流动。非饱和土壤水分流动的达西定律饱和土壤水分流动的达西定律: :()( )mqKqK 或sqKH 饱和土壤水分流动的达西定律: 达西定律的推导达西定律的推导流体中由动量和连续方程可以推导出流体中由动量和连续方程可以推导出渗流中与渗透率成反比渗流中与渗透率成反比速度速度v由流量由流量/孔隙度替换孔隙度替换忽略可压缩性忽略可压缩性得得 的处理的处理渗流边界上速度不为渗流边界上速度不为 达西定律有速度上限，速度太高时需要达西定律有速度上限，速

3、度太高时需要 考虑考虑惯性和湍流效应。而对于多孔介质，骨架会阻惯性和湍流效应。而对于多孔介质，骨架会阻止流体运动，止流体运动， 这一项的存在就不合理了。这一项的存在就不合理了。去除去除 项得项得对于稳态得达西定律对于稳态得达西定律ks，渗透系数=const. 非饱和流饱和流水势组成： =g +p =g +m流动准则： 高 低土水势总水头高 低k( ) 为土壤含水率的函数ks k() k( ) g :p:m = 0相对参考平面之高度至地下水面的高度p= 0m :m 取决于土壤的干湿程度K:2.土壤水分运动参数及其测定方法土壤水分运动参数及其测定方法非饱和非饱和/饱和水力传导度饱和水力传导度容水度

4、容水度土壤水分扩散度土壤水分扩散度l是反映土壤水分在水势梯度作用下流动的性能。一般在饱是反映土壤水分在水势梯度作用下流动的性能。一般在饱和土壤中导水率称为渗透系数，为常量。和土壤中导水率称为渗透系数，为常量。l定义定义：在水势梯度作用下，单位断面面积上流过的水流通：在水势梯度作用下，单位断面面积上流过的水流通量。量。非饱和水力传导度的概念及特征非饱和水力传导度的概念及特征l在非饱和土壤中，在非饱和土壤中，导水率是负压或含水率的函数，随着含导水率是负压或含水率的函数，随着含水率降低而减小水率降低而减小。Hqhk)(Hhkq)(Hkq)(Hqk)(非饱和水力传导度及其测定K(h)nhahk)(hc

5、sekhk)(K()nrsrskk)(bhahkn)(mskk)(scsekk)(K() K(h)多种近似多种近似van Genuchten-Mualem :Brooks and Corey (BC) :HH1H212x平水装置多孔板非饱和水力传导度的测定非饱和水力传导度的测定在水平土柱两端有多孔板，分别在水平土柱两端有多孔板，分别由平水箱保持一定水位，使其负由平水箱保持一定水位，使其负压为压为h h1 1和和h h2 2，在梯度作用下，土在梯度作用下，土柱中土壤水从柱中土壤水从l l端向端向2 2端运移。土端运移。土壤水通量壤水通量q q可由可由l l端补给量或端补给量或2 2端端溢出量测得

6、，两者相等时，水流溢出量测得，两者相等时，水流处于稳定状态。处于稳定状态。非饱和土壤水力传导度可由达西非饱和土壤水力传导度可由达西定律求得。定律求得。HXqXHqHqhk)( 计算的计算的k k是平均的是平均的k k 试样中各点的试样中各点的 是不相同的是不相同的l非饱和达西实验非饱和达西实验)(hkhl在不同的平均负压在不同的平均负压( (吸力吸力) )值下，通量与值下，通量与负压梯度负压梯度成正比，两者成正比，两者呈直线关系，但其斜率呈直线关系，但其斜率( (即水力传导度即水力传导度) )随平均负压而变。随平均负压而变。HXqXHqHqhk)(通量q负压梯度h/xh=-50cmh=-30c

7、mh=-10cmh=0,k=ks容水度容水度(或比水容量或比水容量)单位基膜势单位基膜势( (负压值负压值) )变化所引起土壤含水率的变化，一变化所引起土壤含水率的变化，一般称为容水度或比水容量（般称为容水度或比水容量（C)C)，可以下式表示：可以下式表示： dhdhC表示在单位压力水头降低时自单位体积土壤中所释放表示在单位压力水头降低时自单位体积土壤中所释放出来的水的体积，它与饱和土壤的给水度相似。出来的水的体积，它与饱和土壤的给水度相似。用测水分特征曲线的方法来测定用测水分特征曲线的方法来测定土壤水分扩散度土壤水分扩散度D土壤水分扩散度为单位含水率梯度下，通过单位面积的土壤水分扩散度为单位

8、含水率梯度下，通过单位面积的土壤水流量，其值为土壤含水率的函数，即土壤水流量，其值为土壤含水率的函数，即 baeD CkhkD土壤水分扩散度与土壤的关土壤水分扩散度与土壤的关系可用以下经验公式表示系可用以下经验公式表示扩散率扩散率D D土壤体积含水量土壤体积含水量直角坐标非饱和水分基本方程直角坐标非饱和水分基本方程Hhkq)(3.非饱和土壤水运动的基本方程l理论基础：理论基础：达西定律达西定律质量守恒定律（水流连续原理）质量守恒定律（水流连续原理）直角坐标非饱和水分基本方程直角坐标非饱和水分基本方程)(zzvvzz基本方程推导基本方程推导)(yyvvyy；)(xxvvxx同理同理(y+ y)

9、、(z+ z) 面流速为：面流速为：取微分单元体，体积：取微分单元体，体积： x y z设设沿沿x、y和和z方向方向流速流速分别分别为：为：vx，vy，vz，则则 (x+ x)面流速为面流速为xxyzzyvxvzvy设六面体土壤含水量为设六面体土壤含水量为 ，则，则 t t内内六面体内土壤水质量变六面体内土壤水质量变化量为：化量为：根据质量守恒原理有根据质量守恒原理有pmx y z tt tzyxzvyvxvzyxoutinm-mtzyxtzvyvxvtzyx即即非饱和土壤水运动基本方程，可简写为：非饱和土壤水运动基本方程，可简写为：Hkt)(zHkzyHkyxHkxtzyx)()()(xHk

10、vxx)(根据达西定律有根据达西定律有: :yHkvyy)(zHkvzz)(将上式代入将上式代入zvyvxvtzyx假定土壤各向同性，则有：假定土壤各向同性，则有： )()()(kkkkzyxzHkzyHkyxHkxt)()()( zhkzhhkzyhhkyxhhkxthhc)()()()(1zhzH hhcthhcthht)(基本方程的不同形式基本方程的不同形式zhHxhxHyhyHl用基质势用基质势h h为变量的基本方程为变量的基本方程对于非饱和土壤水，总水头对于非饱和土壤水，总水头H H由负压水头由负压水头h h和重力水头和重力水头z z组成：组成：c（h）表示比水容量（也称容水度）表示

11、比水容量（也称容水度）zHkzyHkyxHkxt)()()( zhKzhhkzyhhkyxhhkxt)()()(令令则有则有对上式求偏导对上式求偏导, ,则有则有故故 zhkzhhkzxhhkxthhc)()()(l用基质势用基质势h h为变量的基本方程为变量的基本方程 zhkzhhkzthhc)()(剖面二维：剖面二维：垂向一维：垂向一维：l以基质势以基质势h为变量的基本方程，最突出的优点是适用于饱和为变量的基本方程，最突出的优点是适用于饱和-非饱和问题的求解，也可用于分层土壤的水分运动的计算，但非饱和问题的求解，也可用于分层土壤的水分运动的计算，但非饱和土壤的导水率和容水度受滞后影响较大，

12、计算中参数选非饱和土壤的导水率和容水度受滞后影响较大，计算中参数选取不当会造成较大误差。取不当会造成较大误差。zkzhkzyhkyxhkxt)()()()( hkDxhxh基本方程的不同形式基本方程的不同形式l用含水量用含水量为变量的基本方程为变量的基本方程D D（）为土壤水的扩散率，）为土壤水的扩散率，yhyhzhzhzkzDzyDyxDxt)()()()( ckhkD zKzhkzyhkyxhkxt)()()(令令则有则有上式中上式中代入上式有代入上式有l用用含水量含水量为变量的基本方程为变量的基本方程zkzDzxDxt)()()(剖面二维：剖面二维：垂向一维：垂向一维：l以含水量以含水量

13、为变量的基本方程常用于求解均质土层或非饱和流为变量的基本方程常用于求解均质土层或非饱和流问题，但不适宜层状土壤或求解饱和问题，但不适宜层状土壤或求解饱和-非饱和问题。非饱和问题。zkzDzt)()(这里达西定律其实就是流体的动量定理的转变，当不考虑温度变化时与连续方程耦合。若要考虑能量方程，就得给出温度条件，再和温度耦合。4.土壤水分运动方程的定解条件土壤水分运动方程的定解条件初始条件（初始条件（t）)()0 ,();()0 ,(zhzhzzii边界条件边界条件)(),();(),(0000thtzhttz以垂向一维流动为例：以垂向一维流动为例：l一类边界（变量已知边界）：一类边界（变量已知边

14、界）：在一维垂向土壤水分运动中，一类边界的情况发生在：在一维垂向土壤水分运动中，一类边界的情况发生在：地表形成积水时；地表形成积水时；地表含水率达到饱和含水率；地表含水率达到饱和含水率；当强烈蒸发时，表土达到风干土含水率。当强烈蒸发时，表土达到风干土含水率。 )(),(00thtzh)(),(0ttzs)(),(0ttzrl二类边界条件二类边界条件( (边界上水流通量已知边界上水流通量已知) )在一维垂向土壤水分运动中，这种情况常发生在降雨、灌水入在一维垂向土壤水分运动中，这种情况常发生在降雨、灌水入渗或蒸发强度已知的边界上。渗或蒸发强度已知的边界上。在降雨或灌水入渗时，在降雨或灌水入渗时，

15、( (t)t)为负值，在蒸发时为负值，在蒸发时 ( (t)t) 为正值。为正值。在不透水边界和无蒸发入渗的边界，在不透水边界和无蒸发入渗的边界， ( (t)t) =0 =0，则，则 hkD)()()(2tkzDB)0)(;()()(2tBzhkzhhk)0)(;()()(2tBzkzD)()()(2tzzhhkB tkzhkB2 thkzhhkB2l三类边界条件三类边界条件：相当于水流通量随边界上的变量：相当于水流通量随边界上的变量( (含水率或压力含水率或压力) )值而变化的情况值而变化的情况在土壤蒸发强度为表土含水率或表土负压的函数的情在土壤蒸发强度为表土含水率或表土负压的函数的情况下，三

16、类边界条件表达式为：况下，三类边界条件表达式为：为变量ffzf321bhahkzhhkbakzD)()()()()(三类边界的一般形式为三类边界的一般形式为Comsol上的算例上的算例Richards模型方程及模型方程及van Genuchten型土壤型土壤Richards模型方程及模型方程及van Genuchten型土壤型土壤储储水系数水系数容水度容水度相对水利传导率相对水利传导率几何几何参数设置参数设置10008.25e-5 or 5.83e-54.4e-101e-80.43 or 0.600.04514.5m2.68根据已有的算例设置土壤参数边界条件边界条件21对称边界流量为0确定流量

17、边界条件边界条件对称边界：对称边界：下边界流量为零：下边界流量为零：上边界流量确定：上边界流量确定： 分段函数分段函数编号编号日降雨量日降雨量（mm/d）降雨时间降雨时间（h）再渗透时间再渗透时间（h）总时间（总时间（h） 降雨类别降雨类别124121224中雨22481624中雨3100121224暴雨410081624暴雨5200121224大暴雨620081624大暴雨初始条件初始条件在不考虑重力的情况下，初始土在不考虑重力的情况下，初始土体的水分应该均匀分布体的水分应该均匀分布设置设置 ，可以使有效饱，可以使有效饱和度和度计算计算一：单种土壤降雨入渗一：单种土壤降雨入渗24mm、100mm、200mm 8小时降雨16小时再渗透过程有效饱和度的变化计算计算一：单种土壤降雨入渗一：单种土壤降雨入渗24mm、100mm、200mm 8小时降雨16小时再渗透过程有效饱和度-深度关系取取 的等值线的深度为降雨入渗深的等值线的深度为降雨入渗深度，统计不同降雨量、不同降雨时间（编度，统计不同降雨量、不同降雨时间（编号号1-6）条件下入渗深度，画出入渗深度）条件下入渗深度，画出入渗深度-时间曲线时间曲线计算计算一：单种土壤降雨入渗一：单种土壤降雨入渗计算计算一：单种土壤降雨入渗一：单种土壤降雨