非饱和湿陷性黄土二维渗流分析

第一章绪论
第一节黄土与黄土的湿陷性
一、黄土和湿陷性黄土的概念
黄土【loess】指的是在干燥气候条件下形成的多孔性具有柱状节理的黄色粉性土，湿陷性黄土受水浸湿后会产生较大的沉陷。

黄土主要分布于世界大陆比较干燥的中纬度地带。

我国的黄土的分布广泛，西起甘肃祁连山脉的东端，东至山西、河南、河北交接处的太行山脉，南抵陕西秦岭，北到长城，包括陕西、陕西、宁夏、甘肃、青海等五个省区的220多个县市，面积达54万平方公里，占全国土地面积的百分之六。

湿陷性黄土是一种特殊性质的土，在一定的压力下，下沉稳定后，受水浸湿，土结构迅速破坏，并产生显著附加下沉，故在润陷性黄土场地上进行建设，应根据建筑物的重要性、地基受水浸湿可能性的大小和在使用期间对不均匀沉降限制的严格程度，采取以地基处理为主的综合措施，防止地基湿陷对建筑产生危害。

二、湿陷性黄土的分布
在我国，黄土和黄土状土广泛分布在华北、西北等地，且地层多、厚度大。

在这些地区，一般气候干燥、降雨量少，蒸发量大，属于干旱、半干旱气候类型。

黄土分布地区年平均降水量在250～500 mm之间。

黄土在自重或一定荷重作用下受水浸湿后其结构迅速破坏而发生显著的附加下沉，以至在其上的建筑物遭到破坏。

这是黄土的一种特殊性质。

我国湿陷性黄土分布约占黄土分布面积的60％，大部分在黄河中游地区。

由于各地黄土堆积环境、地理和气候条件不同，致使其在堆积黄土的物理、力学性质等方面都具有明显的差别，湿陷性有自西向东、自北向南逐渐减弱的规律。

第二节国内外湿陷性黄土渗流研究现状及发展趋势
半个世纪以来，黄土对水的特殊敏感性及其在黄土变形、强度和本构关系等力学特性上所表现出来的影响规律一直是黄土力学特性研究的中心，研究中出现了一系列重要的发展变化：由单轴压缩试验发展到三轴压缩试验，由常规三轴应力路径发展到多种复杂的应力路径，由对浸水湿陷量的研究发展到对湿陷敏感性的研究，由狭义的浸水饱和湿陷发展到广义的浸水增湿湿陷，由单调的增湿湿陷变形发展到增湿、减湿
1
和间歇性湿陷变形，由将黄土作为直线变形体、非线性弹性体发展到作为弹塑性体，由宏观力学特性分析发展到宏、微观结合的力学特性分析，由研究原状黄土发展到研究击实黄土、挤密黄土和饱和黄土，由假定各向同性发展到研究各向异性，由对一些点上黄土的研究发展到黄土地区面上分布规律和经验关系的研究，由静力学特性的研究发展到动力学特性的研究，由对黄土力学特性的描述发展到对黄土力学特性的改造，由曾、应力研究发展到吸力及有效应力研究等。

一、国外湿陷性黄土渗流的研究进展
土体是自然生态环境的一个重要组成部分，它与人类的生活和生产活动有着极为密切的关系。

在土体的固、液、气三相组成中，土体中的水分是最活跃的因素。

州门对饱和土体中水分运动问题的研究有较长的历史，其科学理论的形成始于130年前法国工程师Darcy所进行的试验，经不断完善和发展，形成了现在的地下水动力学。

非饱和土体中水分运动问题一直是土体学研究的重要课题之一。

上世纪下半叶以来，人们对非饱和土体中水分运动问题的研究主要采用形态学观点，定性地分析和描述土体中水分的保持和运动。

随着科学技术的进步，本世纪初Buckinghma(1907)首次将“毛管势”应用于土体水，出现了土体水分的能态学观点。

经30年代和50年代的发展和逐步完善，加之电子计算技术的迅猛发展和各学科的相互渗透，非饱和土体水分运动研究出现了由经验到理论、从定性到定量的深刻变化，从而对土体中水分运动规律的描述更加准确，土体水分特征参数的测定开始变得简单易行。

自从1856年法国工程师达西(H(.Dalcy)通过均质砂土中液体渗流试验的研究总结出渗流的基本公式一一达西定律以来，渗流理论及其应用研究得到了迅速发展。

虽然达西的试验研究是在均质砂土中液体作均匀流的情况下进行的，但是这个研究成果己被后来的学者推广到整个渗流计算中去。

1863年裘布依以达西定律为基础，研究了单向和平面径向地下水流稳定运动。

1889年俄国数学力学家H.E.儒可夫斯基首先推导出了渗流微分方程。

1922年，H.E巴甫洛夫斯基正式提出了求解渗流场的电模拟法，为解决比较复杂的渗流问题提出了一个有效的工具。

1931年，Richards将Darcy的线性渗流理论推广应用到非饱和渗流中以后，人们才开始了非饱和渗流的研究。

水相流所满足的控制方程很快便建立起来，通常称之为Richards方程。

基于Richards控制方程的饱和一非饱和渗流后来得到了深入的研究，并成功地应用到许多实际工程中[1.2]。

Buekley&Leverett(1942)在忽略毛管力时给出T一维两相渗流方程的解。

Coleman和Bodllman(1944，1945)较早的研究了最简单最典型的垂直入渗问题一一干土在积水条
2
件下的入渗。

他们将土体含水量剖面分为四个区:饱和区、过渡区、传导区和湿润区。

这使人布门对入渗过程有了初步的定性认识。

Philip(1957年)研究T一维垂直入渗的级数解。

Parlange(1971，1972)提出了一种半解析迭代方法，求解地表处为第一类边界(含水量已知)和第二类边界(供水强度已知)的一维垂直入渗问题。

Rubin(1968)研究了二维饱和一非饱和土中的非稳定流，他用有限差分法给出二维Richards方程的数值解。

Freeze(1971)研究了三维地下水含水层饱和一非饱和非稳态流，并给出数值解法。

Warren&Root(1963)提出了双重介质渗流问题。

Closmann(1975)提出了三重介质渗流模型。

Lam和Frdelund(1987)对饱和一非饱和土渗流问题作了比较完整的论述，把非饱和土体水运动理论与非饱和土固结理论相结合，得到了符合岩土工程师使用习惯的饱和一非饱和渗流控制方程，并运用二维有限元方法对复杂地下水流动系统的几个暂态渗流实例问题进行了数值模拟。

Remson等(1971)出版了《地下水文学的数值法》专著，他们对70年代以前的饱和一非饱和渗流有限差分数值模拟作了比较全面的回顾和介绍[3,4]。

Neuman< 1973,1974)最早将有限元方法应用到解饱和一非饱和渗流问题。

他用Galerkin法对Richards方程进行空间域的离散，用Crank-Nicolson有限差分格式对时间域进行离散。

不过为了解决边界条件突变时Crank-Nicolson格式不能求解的困难，他建议用隐式向后差分格式对时间域进行离散。

为了解决求解时可能出现的解振荡问题，他建议在时间域离散中采用集中质量矩阵。

他还对非线性迭代时初值的选取作了具体的指导。

Neuman这方面的研究成果后来被人们广泛采用，他的文献因此也成为饱和一非饱和土研究方面的经典之作[5]。

赤井浩一[6]等(1977)在Neuman的基础上研究了考虑土水特征曲线吸湿与脱湿不同情形影响的饱和一非饱和渗流。

他们做了砂槽模型试验，并用有限元方法模拟了砂槽试验，其数值模拟结果与试验结果基本吻合。

Samnori和Tsnboyanma[7](1991)采用Galerkin有限元法模拟了恒定降雨强度下边坡暂态渗流方程，并对边坡的稳定性进行了参数研究，其中考虑了斜坡长度、土层深度、横截面形状和土的性质等因素的影响。

结果表明对斜坡稳定性不利的参数取值情况是较低的导水率，较长的斜坡，较浅的土层和凹形斜坡表面。

Shimada等[8](1995)用Galerkin有限元法模拟了不同降雨强度和不同土的类型条件下的二维非饱和渗流。

并采用刚体弹簧模型进行了斜坡稳定性数值分析。

结果表明较高降雨强度引起安全因数显著降低，渗透性函数对斜坡安全因数降低影响较小，接近饱和时，基质吸力较小改变可引起斜坡安全因数大幅度降低。

Forsyth[21] (1988)对地下水非饱和单相流与两相流模型进行了数值模拟的对比
3
研究，他认为两者存在比较大的差别，两相流的分析方法更接近实际一些。

Fredlund和Hasan在假定气相是连续的，Fick定律和Darcy定律分别适用于气相和水相的流动，水相和气相的渗透系数都是土的基质吸力或某一体积一质量的函数的条件下，提出了求解非饱和土固结过程中孔隙气压力、孔隙水压力的偏微分方程，用Galerkin加权余量原理推导了二维稳态渗流的有限元形式。

Alonso (1995)进行了土坡二维非饱和渗流和极限平衡法的联合分析，渗流分析中采用了考虑空气压力变化的耦合型控制方程，考虑的影响因素包括土的类型、降雨持时、降雨强度、水分保持曲线的形状和土的渗透性。

二、国内湿陷性黄土渗流的研究进展：
国内近十几年来在这方面也开展了许多工作。

李信(1992)等应用伽辽金有限元法对三维饱和一非饱和土渗流问题进行计算研究，介绍了数值方法和主要计算公式，给出了典型算例的计算结果，并与前人的试验资料进行比较。

该研究成果表明，在进行渗流分析时，该方法比只在饱和区内进行饱和渗流分析更接近实际。

另外，将饱和区与非饱和区耦合在一起分析可以避开难以处理的自由水面边界问题。

因而使水位升降、降雨和蒸发等引起的饱和一非饱和渗流问题得到较好的解决。

黄俊，苏向明等[9],(1990)采用伽辽金有限元法隐式向后差分法，编制出简单合理、切实可并以唐山陡河土坝为计算对象，对该坝进行了饱和一非饱和渗流数值计算，计算结果规律正常，令人满意。

吴宏伟((1999)等针对香港地区一种典型非饱和土斜坡，用有限元法模拟雨水入渗引起的暂态渗流场，然后将计算得到的暂态孔隙水压力分布用于斜坡的极限平衡分析。

研究了降雨特征、水文地质条件及坡面防渗处理等因素对暂态渗流场和边坡安全因素的影响。

但是该文在研究中没有考虑雨水入渗随着土体入渗能力的变化而变化的特性，而只是把入渗量按降雨量的一定比例降低来大致确定。

张家发[10] (1997)研究了土坝三维饱和一非饱和隐态渗流。

他编制了三维饱和一非饱和稳定和非稳态渗流程序，该程序曾应用于三峡船闸高边坡渗流场分析中。

陈虹等[11]（1997）采用有限单元法对飞来峡水利枢纽工程纵向导流围堰典型剖面进行渗流数值模拟，分析降雨及非稳态渗流对堰体及高喷防渗治理的影响，反映了堤坝中流体的实际物理现象和运动规律。

马博恒(1997)对露天开采边坡渗流实例进行有限元分析;李爱兵(1994)用边界元法对边坡中的地下水渗流进行了分析。

4
刘洁、毛艇熙[12] (1997)年对饱和与非饱和渗流计算的有限单元法进行了推导，并建立了数学模型。

通过计算，比较了一般的饱和渗流与饱和一非饱和渗流两种模型的计算结果，并引用模型试验资料进行验证。

其结果表明，饱和非饱和渗流计算方法具有不进行自由面调整等优点，在物理概念上更为完善。

陈善雄，陈守义[13.14] (1998. 2001)用积分有限差分方法模拟了降雨条件下土体中水分的运动情况，并对降雨条件下非饱和土坡的稳定性的分析方法做了研究。

朱伟等[15] (1999)通过大型降雨渗透试验实测了土堤内浸润线的变化和水分移动，并利用有限元方法对饱和一非饱和渗流进行数值模拟，进而分析了降雨时土堤内饱和一非饱和的渗透机制。

吴梦喜、高莲士(1999)对饱和一非饱和土体非稳定渗流作了数值分析，对一般的非饱和渗流有限元计算方法加以改进，以消除非饱和渗流数值计算存在的数值弥散现象。

同时还提出了一种逸出面处理新方法，并给出了非饱和月既急定渗流计算的实例。

邵龙潭、王助贫[24] C2000)等采用孔隙介质力学分析方法，把土体骨架、孔隙水和孔隙气分别作为独立的研知寸象，结合孔隙水和孔隙气在气液交界面上满足的力学引牛建立藕合方程，求解非饱和土中孔隙水的入渗和孔隙气体的排出过程。

并对标准砂进行了一维有压水流入渗的试验和计算。

介玉新等(2001)研究了墙后填土中存在二维渗流时对作用在墙背上的土压力和总压力的影响，比较了库仑理论和朗肯理论计算结果的差别.
李兆平(2001)以土体体系含水量作为控制变量，应用非饱和土水分运动基本理论建立了降雨入渗过程中土体瞬态含水量的计算模型，采用非饱和土强度理论和极限平衡方法，得出了可以考虑基质吸力影响的边坡安全系数计算公式，最后通过具体的算例，讨论了降雨入渗对土质边坡稳定性的影响。

但在对渗流场进行模拟时，采用的是一维垂直入渗模型，没有考虑水平方向的影响。

朱文彬、刘宝探(2002)等利用饱和一非饱和土的渗流理论，并运用有限元法和有限差分法对公路边坡在降雨过程中的渗流规律进行了一个实例分析，模拟了雨水在土坡中的渗流过程。

陈力等[16] (2001)采用运动波理论和两次改进后的Green-Ampt入渗模型建立了坡面降雨入渗产流的动力学模型，并得到了试验资料的良好验证。

运用该模型分析研究了简单坡面上降雨入渗产流的动力学规律。

朱军[17] (2001)对饱和一非饱和的渗流的计算做了系统的研究，提出了求解渗透
5
流量的一种简单的方法。

并对随机渗流场的数值分析作了深入的研究。

姚海林、郑少河(2002)等专门针对非饱和膨胀土边坡在考虑暂态饱和一非饱和渗流的情况下进行了参数研究，在研究中着重考虑了裂隙对边坡中孔隙水压力和体积含水量分布的影响，认为对于膨胀土而言，渗透系数越低，越应注意裂隙的作用。

汪自力((2002)等在饱和一非饱和渗流不动网格有限元计算的基础上，寻求用土体单元所受的渗透力代替其周边的孔隙水压力，以达到利用渗流计算时的剖分网格和计算结果，直接连续进行渗流作用下的边坡稳定分析的目的。

但在对非饱和土土骨架进行受力分析时，仍然沿用了饱和土的概念。

吴良骥、[18]雷光耀[19]、张培文[20]、刘德富等对饱和非饱和岩土月稳定渗流问题进行了有限元的分析研究，结合工程实例对堤坝饱和非饱和渗流进行了数值模拟。

上述得到广泛研究的饱和一非饱和渗流都是基于Richards控制方程，它忽略了空气相的流动，仅仅考虑水相在土中的流动。

然而严格意义上的非饱和渗流理论应当考虑空气和水的两相流。

但是在土坡稳定分析工作中，一般忽略空气流动和体变的影响，采用简化的Richards方程进行计算。

张华、陈善雄、陈守义[22]推导了考虑土骨架体变、气相流动的非饱和渗流一般性控制方程，并给出了几种简化形式。

清华大学朱京义、张丙印、张彦等[23]讨论了非饱和土中水和气体两相不相容、不可压缩孔隙流体的流动定律、不同流体间的相互作用以及反映土体孔隙持水特性的土水特征曲线等，据此推导了非饱和土两相流体耦合渗流的控制方程。

邵龙潭，王助贫等[24]采用孔隙介质力学分析方法，把土体骨架、孔隙水和孔隙气分别作为独立的研究对象，结合孔隙水和孔隙气在气液交界面上满足的力学条件建立藕合方程，求解非饱和土中孔隙水的入渗和孔隙气体的排出过程。

并对标准砂进行了一维有压水流入渗的试验和计算，并用数值模拟方法进行对比。

得出在研究水流入渗问题时，对一些导气率较低的土类，考虑气相的压缩和运动的影响是必要的结论。

邓英尔[25]、刘慈群等基于三参数非线性渗流运动定律、质量守恒定律及椭圆渗流的概念，建立了低渗透介质中两相流体椭圆非线性渗流数学模型，运用有限差分法与外推法求解，导出了两相流体椭圆非线性渗流条件下油井见水前后开发指标的计算公式，进行了实例分析。

骆祖江、王增辉((1999)以多相渗流理论为基础，从水气二相渗流的连续性方程和达西定律出发，推导了非饱和带水气二相渗流的耦合动力学模型，讨论了模型的求解
6
方法的原理和步骤。

目前虽然许多学者对非饱和土体的水分迁移训于了大量的研究，但是对.水分迁移机理问题还存在比较大的争议。

综合考虑温度、密度及含水量对水分迁移规律的影响研究比较少，这些都制约了非饱和土水分迁移研究的进展。

第三节本文的研究内容、意义和技术路线
一、研究内容
二、研究意义
三、本文研究的技术路线
第二章湿陷性黄土的工程性质
第一节湿陷性黄土的成因及主要特征
一、湿陷性黄土的定义
黄土是一种第四纪沉积物，具有一系列内部物质成分和外部特征，不同于同时期的其它沉积物。

黄土具有以下全部特征：
1.颜色以黄色、褐黄色为主，有时呈灰黄色；
2.颗粒组成以粉粒(0.05-0.005mm)为主，含量一般在60%以上，几乎没有粒径大于0.25mm的颗粒；
3.孔隙比较大，一般在1.0左右；
4.富含碳酸钙盐类；
5.垂直节理发育；
6.一般有肉眼可见的大孔隙。

当缺少其中一项或几项特征的称黄土状土。

二、黄土的分类
黄土安成因分为原生（或典型）黄土和次生黄土。

一般认为不具层理的风成黄土为原生黄土。

原生黄土经过流水冲刷、搬运重新沉积而形成的具有层理含较多砂粒以至细粒的黄土称次生黄土。

7
次生黄土的结构强度比原生黄土低，湿陷性较高，地质界常将原生黄土称为黄土，次生黄土称为黄土状土。

从工程建设角度看，主要是土的物理力学性质，因此不区别黄土或黄土状土，而通称为黄土。

三、黄土的成因
1、风成说
在当时干旱的大陆性气候作用下，高度风化的黄土物质受到强大的反旋风作用，从中部呈离心状吹向荒漠边缘地区，当遇到异向风或降雨沉落于地面，经风化作用形成黄土。

一般认为我国黄土材料是从中亚西亚搬运来的。

形成黄土的自然环境是干旱或搬干旱地荒漠草原。

我国黄土的颗粒组成、矿物和化学成分都有自南向北逐渐变化的趋势，黄土颗粒自东南向西北逐渐变粗，在同一区域内，高山和低地都有黄土分布，而且黄土堆积常有坡向性，在迎风面堆积量大，背风面少。

黄土地形于下伏基岩有密切关系，常随古地形起伏而起伏，且与下伏基岩性质无关，成分复杂，呈不整合或假整合接触，在构造上，黄土无层理，柱状节理发育，这与具有层理的冲、洪积形成的黄土状土有明显差别。

2、水成说
认为黄土是由冲积、洪积、湖积等的假说。

认为黄土的整个堆积过程于整个地形地貌的发展过程密切相关，早期随盆地四周山坡降水下流而汇集于山间或三角洲处的黄土冲积物堆积成黄土高原，在一定盆地内有一定的分布高度，称为黄土线，黄土线就代表着过去河流淤积的最高地面，超过这一高度就没有黄土分布。

陕西高原黄土来自上游大小盆地，晚期在新构造运动作用下地层上升切割形成的河谷中，黄土沉积物堆积成阶地性状，在大陆性干旱气候下，这些沉积物在风化和成土作用下形成黄土。

该学说将原生黄土与次生黄土的沉积过程与地貌形成过程统一起来，将黄土堆积物概括为两大类：较古老的高原黄土，年青的阶地黄土。

3、多种成因说
认为各地区黄土形成的地质地理环境以及这些环境的演化历史是不同的。

利用地层、地貌、岩性、古地理、地球化学、土的物质成分等分析方法的假说。

该假说将物质来源、堆积方式、黄土的形成以及演化看做是一个统一的过程。

划分成因类型的主要标志是地质因素。

根据之一原则，将我国黄土划分为十种类型：冲
8
9
积、洪积、坡积、风积、冰水沉积、湖积、洪积－坡积、冲积－洪积、残积－坡积、冲积－坡积。

河谷平原类型内的黄土状土的成因主要是冲积类型。

另外，还有土壤残积假说。

黄土成因虽有不用的假说，但分歧主要是对分布于高原和高分水岭上的黄土，而分布于河谷地带的黄土则意见较为一致，都认为是冲积类型，也间有洪积－坡积类型。

目前，我国大部分建筑工程多兴建在河谷地带的一、二级阶地上或是靠山、近山的洪－坡积地带，在高原上较少。

四、黄土的分布
黄土分布很广，面积达1300万平方公里，约占陆地总面积的9.3%。

世界各大洲黄土覆盖面积占总面积的比例为：欧洲7%，北美5%，南美10%，亚洲3%，此外，在澳大利亚、北非也有零星分布。

我国黄土分布面积635280平方公里，占世界黄土分布总面积的4.9%左右，主要分布在北纬33-47度，以34-45度之间最为发育，属于干旱、半干旱气候类型。

我国湿陷性黄土分布面积约占我国黄土分布总面积的60%左右，为27万平方公里，大部分在黄土中游地区，北起长城附近，南达秦岭，西自乌鞘岭，东至太行山，即北纬34-41度，东经102-114度之间。

湿陷性黄土一般都覆盖在下卧的非湿陷性黄土层上，厚度以六盘山以西地区较大，达30米，六盘山以东地区稍薄，如汾渭河谷多为几米至十几米，再向东至河南西部则更少，并且有非湿陷性黄土位于湿陷性黄土层之间。

我国黄土地层的化分 地质时代
地层名称 说明 全新世4Q
近期Q 24
新 黄 土 新近堆积黄土 一般有湿陷性、高压缩性 早期Q 1
4
一般湿陷性黄土 有湿陷性 晚更新世Ｑ3
马兰黄土 中更新世Ｑ2
离石黄土 老 黄 土 一般无湿陷性
早更新世Ｑ1 午城黄土 由于我国各地黄土堆积环境、地理、地质和气候条件不同，致使其在堆积厚度、
土的物理、力学性质等方面都有明显的差别，如湿陷性具有自西向东，自北向南逐渐减弱的规律。

第二节湿陷性黄土的物理化学性质
一、湿陷性黄土的物理性质
黄土是由固、液、气三相组成的，其三相组成间重量和体积的比例关系，可以反映出一系列物理性质，这些性质常用以下一系列指标表示：颗粒组成，土粒比重，含水量，重度，孔隙比，饱和度，液限，塑性指数。

1、颗粒组成:
土的颗粒是指那些岩石、矿物和非晶体化合物的零散碎片或碎屑。

颗粒本身可以是矿物的结晶构造，也可以是非结晶构造，如二氧化硅和氢氧化铁。

不同地质时代的黄土，其颗粒组成不同，第四纪早期的黄土比晚期黄土中的粘粒含量高，而砂含量则低。

湿陷性黄土以粉粒为主，含量达60%以上，其中细粉粒(0.05-0.01mm)占7-9%，粗粉粒(0.01-0.005mm)占45-65%。

我国由西北向东南方向，砂粒减少而粘粒增多，这与我国湿陷性黄土由西北向东南递减的趋势大体相关。

2、土粒比重和天然重度:
黄土的土粒比重一般为2.51-2.84，平原区为2.62-2.76，比重大小与土的颗粒组成有关，当粗粉粒和砂粒含量较多时，比重常在2.69以下，粘粒含量多，则在2.72以上。

黄土的颗粒组成与其液限、塑限有一定关系。

Ip G
<7 2.67
7-10 2.69
10-13 2.71
13-17 2.72
>17 2.73-2.74 湿陷性黄土天然重度一般在13.3-18.1kN/m3之间，它不仅取决于颗粒的大小和含量的多寡，还与土的含水量有关，一般工程常用干重度和孔隙比反映土的密实度。

10。