土的饱和－非饱和渗流研究现状

渗流力学最先是作为许多学科的边缘学科而发展的，其既是水力学和岩土力学所不可缺少的组成部分，同时也是水工结构、水文地质、地下水文、农田水利、灌溉排水、地下水资源开发利用以及石油开采等学科中的部分内容。渗流力学作为渗流学科的核心课程，在20世纪更是受到国际学术界和工程界的高度重视，随着各有关学科的不断发展和生产实践中提出的渗流问题日益广泛复杂，渗流己逐渐发展成为具有自己的理论、方法和应用范围的独立学科。目前，已有渗流水力学、岩土水力学、地下水力学、渗流力学、地下水动力学、多孔介质中流体动力学等多部渗流专著出版（柴军瑞等，1970；庄礼贤等，2009；陈崇希等，1999；孔祥言，2010）。

渗流发展始于18世纪中期，1852—1855年，法国工程师达西通过垂直圆管中的砂土透水性渗流试验提出了线性渗流定律（达西定律），为渗流理论的发展奠定了基础。1889年茹可夫斯基首先推导出了渗流微分方程，但仅求出了适用于均质渗流介质和简单边界条件的解析解。1922年巴甫洛夫斯基提出了求解渗流场的电拟法，为解决比较复杂的渗流问题提供了一个有效的工具，但它无法模拟非均值各向异性渗流介质，对复杂的地质和边界条件也不尽适应，为了解决更复杂的渗流问题，逐步发展和研究了电网模拟法（毛昶熙，2003；杜延龄等，1992）。

1931年Richards将Darcy推广应用到非饱和渗流中以后，人们才开始了非饱和渗流的研究。因此，直到数值模拟应用到Richards控制方程中以后，使得饱和－非饱和的渗流场获得合理的数值解成为可能。早期主要采用有限差分的方法求解Richards方程，随着有限元技术的发展和逐渐成熟，有限元方法成为求解饱和－非饱和渗流问题的主要方法。

20世纪60年代，随着计算机的发展，数值方法，即有限差分法、有限单元法和边界元法在渗流分析中得到了愈来愈广泛的应用（唐辉明等，2002；唐辉明等，2000；胡新丽等，2006；朱伯芳，2000；王元汉等，2001；Richard等，1971；Neuman，1973，1974）。以往用数值法分析渗流问题，一般是以自由水面为边界，在饱和区内进行计算研究，但这种分析方法在计算的每个时段都要试求自由水面边界，过程比较繁琐，而且也没有考虑非饱和区的孔隙水压力状况，因此不能全面真实地反映地下水的渗流动态。为此，国外从20世纪70年代开始考虑非饱和区域水的流动，即将饱和区与非饱和区一起进行整体分析。压力水头在饱和区为正值，而在非饱和区为负值，零压力面就是自由水面，即饱和区与非饱和区的分界面。这样计算域内不再有自由水面边界，使得计算简化，程序处理也比较容易。

在有限差分数值模拟饱和－非饱和渗流问题方面，Rubin对二维饱和－非饱和土进行了研究，用有限差分的方法给出了Richards控制方程的数值解；弗雷兹（Freeze，1979）研究了地下含水层饱和－非饱和非稳态流，给出了数值解法，并对以前的饱和－非饱和渗流有限差分数值模拟进行了比较系统的阐述和总结。

有限单元法的基本思想最早在1943年由柯朗提出，1965年津克维茨等提出了有限元法适用于所有可按变分形式进行计算的场问题，为该方法在渗流分析中的应用提供了理论基础。有限单元法最早应用于分析稳定渗流场问题，20世纪70年代初才被应用于非稳定渗流场问题，随后用来解决各种具体的工程问题（唐辉明等，2002；唐辉明等，2000；胡新丽等，2006；丁秀丽等，2004；汪斌等，2006）。津克威茨（Zienkiewicz O C）第一次具体地将有限单元法应用于稳定渗流的分析，德赛（Desai，1985）将有限单元法应引入非稳定渗流，Remson等（1971）出版了《地下水文学的数值法》专著，他们对20世纪70年代以前的饱和－非饱和渗流有限差分数值模拟做了比较全面的回顾和介绍。

最早应用有限元方法求解饱和非饱和渗流问题的是纽曼（Neuman，1974，1977），他采用Galerkin法对Richards方程进行空间域的离散，用Crank－Nicolson有限差分格式对时间域进行离散，在1973年提出了一维非饱和渗流的有限元法之后，又提出了解决二维饱和－非饱和渗流问题的有限元法。同时对非线性迭代时初值的选取做了具体的指导。纽曼以后，许多专家和技术人员对饱和非饱和渗流问题进行了更加广泛深入的研究和探讨，对Richards方程从不同的方面进行了大量的数值计算工作，也取得了一些可喜的成果和经验。

赤井浩一等（1977）日本学者采用了纽曼的数值模型和有限元法进行了试验与数值计算，研究了考虑土－水特征曲线吸湿与脱湿不同情形影响的饱和－非饱和渗流，做了砂槽模型试验，用有限元方法对砂槽试验进行了数值模拟，而且数值模拟的结果与砂槽模型试验结果比较接近。Fredlund（1998）用Galerkin加权余量原理推导了二维稳态渗流的有限元形式；Mille（1962）提出了非饱和介质的渗透系数是含水量或水压力水头的函数即kr（θ）和kr（h），这种函数关系随多孔介质性质的变化而变化。Genuchten（1984）提出了一个土体非饱和参数的表达式，在这个表达式中存在未知的参数，可以根据土体的非饱和特征曲线krθ和θ－h来确定。这种方法在实际的计算中应用比较广泛。

近年来，由于非饱和土理论的发展，国内外学者越来越注意非饱和土的性质研究和对边坡或滑坡稳定性的影响研究。如加拿大弗雷德隆德（Fredlund，1998），应用其所发展的非饱和土渗流理论，对香港地区花岗岩残积土边坡进行了研究，强调基质吸力对边坡稳定性的作用和影响，认为许多非饱和土组成的稳定边坡，在长时间连续暴雨期间，有可能由于基质吸力消失而发生灾难性的破坏。Lam和Fredlund在1984年应用饱和－非饱和渗流分析程序Trasee对一些渗流问题进行了求解；对饱和－非饱和土体的渗流问题做了较完整的论述，把非饱和土壤水运动理论与非饱和土固结理论相结合，得到了岩土工程师使用习惯的饱和－非饱和渗流控制方程，并用二维有限元方法对复杂的地下水流动系统的几个暂态渗流实例进行了数值模拟。

我国也有不少专家对饱和－非饱和渗流做过研究。毛昶熙（1997）对饱和与非饱和渗流计算的有限单元法进行了推导，并建立了数学模型。通过计算，比较了一般的饱和渗流与饱和－非饱和渗流两种模型的计算结果，并引用模型试验资料进行验证。其结果表明，饱和－非饱和渗流计算方法具有不进行自由面调整等优点，在物理概念上更为完善。汪自力、李莉等（1991）提出了堤坝饱和－非饱和有限元分析的高斯点法，该方法计算量比较小，迭代的格式也比较简单。如黄俊（1990）采用伽辽金有限元法编制了土坝饱和非饱和渗流程序；张家发等（1998）研究了土坝三维饱和－非饱和非稳定渗流，并编制了三维饱和－非饱和稳定和非稳定渗流程序，曾应用于三峡船闸高边坡渗流场分析中；陈虹（1997）采用有限元法分析了降雨条件下的非稳定渗流；陈守义（1997）用积分有限差分法模拟了降雨条件下土体中水分的运动情况；吴梦喜（1999）也对饱和－非饱和渗流分析进行了改进。朱岳明等（2003）提出了三维饱和－非饱和渗流场求解及其逸出面边界条件处理方法。

库水位波动过程中斜坡内的地下水问题属非稳定渗流问题，由于斜坡地下水的排泄条件随库水位的下降处于不断变化当中，求解这样一个非稳定渗流场是非常复杂的问题，稳定性分析和工程设计中，常采用简化模型。受传统的地下水模式的影响，一般将滑坡地下水概化为具有连续的浸润面的潜水，以渗流计算得到的水位来分析滑坡的稳定性。因此，国内很多学者把求解斜坡内的非稳定渗流场问题归结为求解自由面问题，引入使问题简化的假设条件，提出了各种估计自由面位置的方法。

无论采用哪种数学模型或哪种数值计算方法来求解渗流场，均存在一个共同的问题，即渗流自由面位置事先未知，这使得边值问题具有未定边界，因此只能采用逐次逼近的试算方法求解。早期学者所做的研究工作就是基于这种思想，用试算法不断改变自由面位置，直到计算结果满足自由面上的边界条件为止，从而获得近似解答。

彭良泉认为，自由面下降的位置在库岸滑坡稳定性分析中有着决定性的意义。他根据Neuman（1974）提出的采用Galerkin加权残值法求解二维饱和－非饱和渗流的有限元法，并根据自由面下降的快慢提出了滑坡稳定性计算的危险水力条件。

目前有限元方法已经成为渗流计算的主要的数值计算方法，主要是因为有限差分法处理不规则的边界比较困难，有限单元法在模拟曲线边界和各向异性方面比有限差分法具有更大的灵活性（唐辉明等，2002；汪斌等，2006；刘洁等，1997）。

在国内，对斜坡地下水渗流的研究尚不够深入，也不够系统，大多停留在定性分析和评价的层面上，有待于研究理论和分析方法等方面不断的补充和完善，特别是对受地下水渗流影响显著的库岸边坡进行数值模拟与分析，进一步阐明库岸边坡的失稳机理和成灾规律。