comsol在非饱和土渗流的应用

1、基于 comso 的非饱和土渗流研究/comso 祢岩土工程渗流的应用摘要：岩土工程的核心难点即解决地下水问题，一般岩土工程事故都是由于对地下水的影响重视不够而造成的，然而解决这一难点关键在于解决地下水渗流问题。目前对于非饱和土渗流研究的理论仍相对落后，本文结合非饱和土渗流场基本方程以及由水土特征曲线得到的相关渗流参数(渗透系数，体积含水量)，阐明了如何解决渗透模型要求渗流场方程的连续性与现场实测数据的非连续性之间的矛盾，并利用 comsolMultiphysics 软件对某工程中非饱和土渗流问题进行了模拟，并验证了 Fredlund 和 xing(1994)土水特征曲线方程的正确性。 这种解

2、决非饱和土渗流问题的思想可供学者参考。关键词：非饱和土；渗流场；渗流参数；连续性矛盾；ComsolMultiphysicsStudyonseepageofunsaturatedsoilseepagebasedoncomsolAbstract:Thecoredifficultyofgeotechnicalengineeringistosolvegroundwaterproblems,thegeneralgeotechnicalengineeringaccidentsareduetotheimpactofgroundwatercausedbyinsufficientattention,howeve

3、r,thekeytosolvethisdifficultproblemistodealwiththegroundwaterflow.Atpresent,forthestudyofunsaturatedsoilseepagetheoryisstillrelativelybackward,thispapercombinesbasicequationofunsaturatedsoilseepagewithsoil-watercharacteristicscurveandobtainstherelevantflowparameters(hydraulicconductivity,volumetricw

4、atercontent)fromthem,andillustrateshowtosolvetheconflictbetweentheseepagefieldpenetrationmodelrequiringEquationofcontinuityandthemeasureddataofnon-continuity,andusingthesoftwarecomsolMultiphysicstosimulateunsaturatedsoilseepageproblemsinoneprojectandverifiedtherightofFredlundandxing(1994)soil-waterc

5、haracteristiccurveequation.Theideaofsolvingunsaturatedsoilseepageproblemsmaybereferredbysimilarprojects.Keywords:unsaturatedsoil;seepagefield;seepageparameters;continuousconflict;ComsolMultiphysics1 引言岩土工程设计与施工的难点在于解决地下水问题，一般岩土工程事故都是由于对地下水的影响重视不够而造成的， 像2003年7月14日上海轨道交通4号线工程事故； 2007年8月17日山东新汶煤矿透水事故；

6、2008年11月15日杭州地铁工地塌陷事故以及2011年 1 月 1 日杭州余杭区-工地土方坍塌事故等等都是由于忽视地下水的影响而造成的。 然而解决这一难题的关键在于解决地下水渗流问题。虽然众多国内外学者对土的渗流问题做了大量的研究，但是目前对于非饱和土渗流研究的理论以及实践应用仍相对落后。一般来说，解决非饱和土渗流设计的问题以及与其相关的工程实践问题，可以归结于就具体的非饱和土渗流工程概况而建立渗流场基本方程，然后解这一渗流场基本方程，从而得出相关的渗流流线（水位）分布、水流渗流力矢量分布、流速矢量分布和相关的趋势，最后以此来指导实践施工。然而在求解非饱和土渗流场基本方程时，首要要解决两个重

7、要未知参数，即体积含水量 8 和渗透系数 k,这两个参数在实际工程中是通过实验得到的，试验得到的是一系列孤立的点，然而这与渗流场基本方程建立于连续性模型相悖，这就给求解渗流场基本方程带来了很大的困难，于是国内外很多学者对此进行了大量的研究。为了解决非饱和土的体积含水量 9 与渗透系数 k 测量的耗时以及模型的连续性问题，Gardner;Brooks&Corey,VanGenuchten,McKee&Bumb 和 Fredlund&Xing 等人1先后通过试验建立了一系列的土水特征方程即体积含水量 9 与吸力甲的拟合连续性表达式；但是渗透系数 k 的连续性仍未解决，后来

8、E.C.Leong 和 H.Rahardjo 等人通过大量的试验发现非饱和土的体积含水量 8 和渗透系数 k 都与吸力甲存在密切的联系，并且依赖于土水特征曲线，于是建立了体积含水量 8 与渗透系数k 拟合关系。于此解决了求解非饱和土渗流场基本方程的首要问题。本文以以上内容为核心思想，利用 ComsolMultiphysics 软件，对具体工程中非饱和土渗流问题进行了数值f 真拟，通过反馈验证了 Fredlund 和 Xing4等人士水特征曲线方程的正确性。以上解决非饱和土渗流问题的思想可供学者参考。2 模型的建立 2.1非饱和土渗流场基本方程一般来说，解决非饱和土渗流设计的问题以及与其相关的工

9、程实践问题，特别是数值分析，都归结于就具体的非饱和土渗流工程概况而建立渗流场基本方程，然后解这一渗流场基本方程，从而得出相关的渗流流线(水位)分布、水流渗流力矢量分布、流速矢量分布和相关的趋势，最后以此来指导实践施工。二维非饱和土渗流场基本方程为：式中：H 一总水头；kxx 方向的渗透系数(变量)；kyy 方向的渗透系数(变量)；Q应用边界流；0体积含水量；t 一时间。由于孔隙气压力通常为常量 ua,对体积含水量 9 的改变没有影响，改变的体积含水量由下式决定：.F-mwfuw式中：mwuw-G)关系曲线的斜率。H 的表达式为：H=亍+y(3)w式中：uw孔隙水压力；M 水的重度；y海拔。(3

10、)式可得：uw=(H-y早w(4)由(1)(2)(4)得：已f小)，已己H)FH-、一.kx十一ky一十Q=mw，w(5)exvexJcycyJct在求解上述非饱和土渗流场基本方程(1)需要考虑边界条件和初始条件，边界条件可以分为三类3:第一类边界条件为水头边界条件，即 h|n=H1(x,y,-=1%=H2(x,y,t),其中 qn为单位二 n(D(2)t);第二类边界条件为流量边界，即：2面积边界上穿过的已知流量；第三类边界条件为混合边界条件，即=p，其中&B为参数。初始条件表示为：h(x,y,t)|t=0=h0(x,y,t)。2.2 土水特征曲线方程在求解非饱和土渗流场基本方程(1

11、)时，要知道体积含水量 8 和渗透系数 k,而这两个参数在实际工程中是通过实验得到的，试验得到的是一系列孤立的点，然而这与渗流场基本方程建立于连续性模型相悖。 为了解决非饱和土的体积含水量 9 与渗透系数 k 测量的耗时以及模型的连续性问题，Gardner,Brooks&Corey,VanGenuchten,McKee&Bumb 和 Fredlund&Xing 等人1先后通过试验建立了一系列的土水特征方程即体积含水量 9 与吸力甲的拟合连续性表达式。但是没有形成统一的体系，后来 E.C.Leong 和 H.Rahardjo 通过大量实验对比在1中建议我们采用 Fredl

12、und 和 Xing(1994)的土水特征方程，其为:式中：8w 为体积含水量；8s 饱和体积含水量；里，a,b,c 为未知量。通过试验测得非饱和土的体积含水量8与吸力甲关系的数据点组， 可以拟合成如(6)式的连续性表达式，从而解决 9 的连续性。需要注意的是 E.C.Leong 和 H.Rahardjo 在1中提到 Fredlund 和 Xing(1994)拟合效果最差，11 个点的数据集，而不是 6 个点的数据集。2.3 渗透参数试验方程上面解决了体积含水量 8 非连续性的矛盾， 渗透系数 k 的非连续仍未解决， 后来很多学者人通过大量的试验发现非饱和土的体积含水量 8 和渗透系数 k 都

13、与吸力里存在密切的联系， 并且依赖于土水特征曲线， 于是建立了体积含水量 8 与渗透系数 k 拟合关系。E.C.Leong 和 H.Rahardjo 在2中提到的拟合关系如kr-Op式中：kr 一相对渗透系数，kr=kw/ks,kw 为渗透系数，ks 为饱和渗透系数；0 一标准体积含水量，=(9w-9r)/(9s-9r),8r 为残余体积含水量；p为常量(一般在 2.524.5 之间)。由土水特征曲线得到饱和渗透系数 ks 和残余体积含水量 Or,以及由(6)式从而解决了非饱和土渗透系数 k 的连续性劭h+ctaF或4ab(6)甲为吸力，与孔隙气压力和孔隙水压力有关;133 工程实例模拟3.1

14、 工程概况某水库如图 1 所示，具体物理参数在图中已给出，其中有 8m 的透水层；试验测得的体积含水量 8 以及渗透系数 k 与吸力甲关系的数据如表 1,水的重度旭=9.8kN/m3,初始水位高度为，H0=10m。44m8m图 1 某水库的物理参数Fig. 1 Physicalparametersofareservoir表 1 土水特征参数Table1Soil-watercharacteristicparameters实验组号23456渗透参数、1吸力中(kPa)020406080100体积含水量 e0.38250.35000.20030.09960.07150.0627渗透系数 k(m/s)

15、0.00790.00320.000520X10-644-61.30644-61.0063.2 模拟参数的处理(1)体积含水量根据表 1 中测得的数据绘制出吸力甲与体积含水量 9(48 关系图)如图 2 所示。图 2 实测甲和 9 关系图Fig. 2 MeasuredWand8diagramsFredlund 和 xing（1994）提到（6）式中的系数 C（由系数接近于 1,简化计算，取C（甲）=1,饱和体积含水量8s=0.383,根据图 2,拟合的含水量连续性方程（6）式，可得（6） 式中的参数为 a=30.72;b=3.851;c=1.257 （一般可以用 matlab 或 origin

16、软件来拟合，本文通过 origin 拟合）。（2）渗透参数根据表 1 和图 2,（7）式中饱和渗透系数 ks=0.008m/s；残余体积含水量8r=0.062;通过拟合可得到 p=4.52 再由（6）式就可以得到了渗透系数 k 的连续性方程。本工程虽然已测出了渗透系数 k 和吸力甲的关系，仅供验证的目的。（3） mw处理将已求出未知参数的含水量连续性方程（6）式中的吸力甲换成-uw,并对 Uw求一阶导数即得到 mw。3.3comsol 的模拟在 comsol 中利用非饱和土渗流场基本方程（5）式，输入已求出的连续性表达式，以及边界条件（图 1 中，BC 为第一类边界条件；CD 为第二类边界条件

17、，DE 为第三类边 界 条 件 ； 在 comsol 中 需 要 转 化 相 应 的 Dirichletboundarycondition 和Neumannboundarycondition 或这两类边界条件的组合），本例的模拟过程和模拟结果如图 3图 110.0501200.1020406080100吸力用（kPa）O O5 54,3,3SO.SO.05050505,32,2V,32,2VSSO.SSSO.S%图 3 水库划分网格Fig.3Reservoirmesh图 4 水库初始应力分布图Fig.4Initialstressdistributionofthereservoir图 5 一年后

18、水库应力分布图Fig. 5 Reservoirstressdistributionafteroneyear图 6 初始浸润线位置Fig. 6 Theinitiallocationofsaturationlineii广&ingy.ki;!_i.J图 7 一年后浸润线位置Fig.7Locationofsaturationlineafteroneyear图 8 初始流速矢量Fig.8Initialvelocityvector图 9 一年后流速矢量Fig. 9 Velocityvectorafteroneyear图 10 初始流线图Fig. 10 Initialflowchart1g.fML；.LlJ.Tj.iiJiJLjlMavtridiiu7.UHdll图 11 一年后流线图Fig. 11 Flowchartafteroneyear4 结论本文提供了一种解决非饱和土渗流问题的方法，即结合非饱和土渗流场基本方程以及由水土特征曲线得到的相关渗流参数（渗透系数，体积含水量），以此解决渗透模型要求渗流场方程的连续性与现场实测数据的非连续性之间的矛盾，并利用数值软件对具体工程中非饱和土渗流问题进行模拟，