降雨入渗下膨胀性黄土隧道围岩破坏演化

张万志 徐帮树 曾仲毅 田斌华 石伟航

(1山东大学岩土与结构工程研究中心, 济南 250061)(2贵州省交通规划勘察设计研究院股份有限公司, 贵阳 550081)

膨胀性黄土是黄土地区自然地质变化过程中形成的一种对建筑工程具有特殊危害性的地质体,兼具黄土和膨胀土的特性[1-2],在我国中西部地区分布广泛.膨胀性黄土具有的胀缩性和裂隙性对建筑工程具有严重危害,较为突出的是持续降雨作用下穿越膨胀性黄土地层的隧道施工遇到的大变形、塌方、洞口滑坡等问题[3].

降雨入渗是诱发膨胀性黄土地层隧道工程失稳破坏的主要因素,其影响可概括为3个方面：① 雨水渗入土体后,导致土体原有非饱和区内部的负孔隙水压力急剧降低,根据非饱和抗剪强度理论,负孔隙水压力的降低甚至消失将引起土体抗剪强度的降低[4]；② 入渗的雨水将导致土体不同部位出现暂态饱和区,增大区域土体自重,使隧道洞顶的下滑力增加,且土体体积含水率增大也将使土体强度显著降低；③ 膨胀土会随降雨增湿过程而吸水膨胀,其膨胀变形在受到初期支护的约束作用后,对其产生一定的膨胀力作用,从而对支护效果产生较大影响[5].以上3个方面与土体的渗流特征密切相关,因此研究膨胀性黄土隧道在降雨条件下的渗流特征是研究其稳定性的必要前提.

目前,许多学者已对降雨入渗作用下隧道和边坡工程的稳定性开展了相关研究.王军等[6]基于有限差分方法,通过自编程序,对降雨入渗下隧道边坡力学-流变-渗流特性进行数值模拟分析,研究了隧道围岩潜在滑动带和破坏影响区,确定了隧道围岩重点加固范围.陈伟等[7]利用渗流-应力耦合方法对降雨入渗条件下隧道有无排水系统情况的围岩体力学特性进行数值分析,研究表明隧道内部不设排水条件时,降雨入渗后围岩体的塑性区扩展明显,对隧道稳定性有极大影响.冯丙阳[8]基于热膨胀与增湿膨胀外在形式的相似性,以热膨胀模拟增湿膨胀,开展了降雨入渗下某膨胀性黄土隧道在不同膨胀应力作用下的围岩体变形研究,并提出了“钢拱架+钢格栅”联合支护技术.综上可知,现有成果已从降雨入渗条件下非饱和渗流分析、力学计算、热膨胀变形等角度进行了隧道围岩稳定性的研究,但在综合考虑降雨入渗条件下“非饱和渗流-围岩软化-膨胀力”共同作用的膨胀性黄土隧道围岩稳定性方面尚缺乏系统性分析.

鉴于此,本文提出了降雨入渗下膨胀性黄土隧道围岩大变形破坏的演化规律研究方法.通过FISH语言改进FLAC3D渗流计算模块,使其更适用于隧道工程的非饱和渗流计算；等效模拟渗流过程中土体抗剪强度参数(内摩擦角和黏聚力)随含水率增大而下降的关系,并在每个计算时间步中考虑膨胀性黄土的增湿膨胀作用,以此来分析降雨入渗下膨胀性黄土隧道围岩稳定性问题.研究成果在膨胀性黄土隧道围岩稳定性分析和变形控制方面具有一定的理论意义和应用价值.

1 FLAC3D膨胀土非饱和降雨入渗分析功能二次开发

1.1 饱和非饱和渗流理论[9]

土体饱和非饱和渗流偏微分方程形式如下：

(1)

式中,Vw为土单元体中水的体积;vx,vy,vz分别为x,y,z方向流速;t为时间.等式左边是单元体流进流出水的体积差,等式右边为单位时间单元体水体积变化量.土体饱和-非饱和渗流服从达西定律：

(2)

(3)

式中,V0=dxdydz.研究表明,土体单元中水体积变化可通过单位时间单元体净法向应力σv-ua和基质吸力ua-uw的变化确定,关系如下[10]：

(4)

非饱和渗流过程中,通常假设总应力不随时间变化,即∂σv/∂t=0.假设土孔隙与外界大气相通,则∂ua/∂t=0.总水头h取为重力水头和基质吸力水头之和,即

h=y-uw/(ρwg)

(5)

式中,ρw为水的密度；y为重力水头；g为重力加速度.

将式(4)、(5)代入式(3),可得到土体非恒定渗流偏微分方程：

(6)

在已知初始条件和边界条件的情况下,可通过式(6)得到孔隙水压力分布随时间的变化情况.

1.2 非饱和土土水特征函数与渗透系数函数

在饱和非饱和渗流数值分析中,土体的孔隙水压力、渗透系数与单元饱和度有关.为实现FLAC3D土体非饱和渗流计算功能,需定义非饱和土体孔隙水压力和渗透系数与土体含水率的关系，即土-水特征曲线和渗透系数关系曲线.

含水量可表示为含水率或饱和度.Fredlund等[11]通过大量实验,利用数理统计方法推导了土体含水率与基质吸力的函数关系方程：

(7)

式中,ψ为基质吸力；ψr为残余含水率对应的基质吸力；ω,ωsat分别为土含水率和饱和含水率；a,m,n为试验拟合参数.

Nguyen等[12]研究表明土体非饱和渗透系数kw与土饱和度s的函数关系为

kw=ks(s)ksat

(8)

式中,ksat为饱和渗透系数；kw为非饱和渗透系数；ks(s)为比例因子,且ks(s)=s3,因此式(8)简化为

kw=s3ksat

(9)

1.3 非饱和土抗剪强度函数

文献[13]研究表明,土的黏聚力、摩擦角分别与含水率呈近似指数和线性关系：

c=Ae-Bω

(10)

式中,c为黏聚力；A,B为与土样初始条件有关的系数.

φ=F-Eω

(11)

式中,φ为内摩擦角；E,F为与土样初始条件有关的系数.

1.4 非饱和膨胀土膨胀力函数

文献[14]研究了增湿过程中膨胀土膨胀力随含水率变化的关系,指出在一定范围内,含水率与膨胀力基本呈线性变化.若膨胀土初始含水率为15%,增湿程度至23%时,膨胀力发挥90%以上,该范围内膨胀力近似为线性变化.因此,在数值计算中可假定膨胀力与含水率的关系为线性关系,表达如下：

(12)

式中,ΔP为含水率ω对应的膨胀力；ω0为初始含水率；ω1为最终含水率；Ps为最终含水率ω1对应的膨胀力.

综合上述分析,利用FISH语言在FLAC3D中实现降雨入渗作用下膨胀性黄土隧道围岩大变形破坏演化机理分析的方法如下：

① 模型建立及初始化.建立三维模型,输入各初始材料参数(体积模量、剪切模量、内摩擦角、黏聚力、渗透系数等),赋予初始边界条件(初始饱和度分布、孔隙水压力分布、重力场分布等).

② 计算一个渗流时间步.计算后节点饱和度(或含水率)将会发生变化.

③ 更新基质吸力(负孔隙水压力).检测节点饱和度,若饱和度不等于1(不饱和),则依据土-水特征曲线函数,通过节点含水率求得基质吸力,并赋值给节点孔隙水压力；若饱和度等于1(饱和),则直接进行下一步计算.

④ 将节点饱和度转化为单元饱和度.由于FLAC3D中饱和度是节点变量,而渗透系数等是单元变量,因此需要将节点饱和度转化为单元饱和度,才能通过单元饱和度求取单元渗透系数.其中单元饱和度的转化可以通过节点饱和度的加权插值方法实现.

⑤ 通过单元饱和度求得对应的渗透系数,并更新单元渗透系数值.

⑥ 通过单元含水率更新求得对应的内摩擦角和黏聚力,且赋给单元.

⑦ 通过单元含水率增量计算得到膨胀力增量,假定三向膨胀力增量大小相等,施加给相应单元.

⑧ 计算一个力学时间步.若不满足平衡条件或小于力学计算子步数,将重复过程③～⑧,否则进入下一步计算.

⑨ 检验是否达到渗流计算终止条件,若未达到则重复②～⑨继续计算,若是则结束计算.

2 工程实例

2.1 工程概况

太兴铁路TXXS-2标段小河沟隧道位于山西省太原市娄烦县向阳村,隧道全长1 803 m,隧道最大埋深约82 m,最小埋深约5 m.该隧道地貌为黄土梁、峁地貌,其间冲沟发育,多呈“V”字形.地表分布砂质黄土,呈黄褐色,稍密至中密；地表以下20～60 m为上第三系(N)膨胀性黄土,呈红褐色,硬塑,结构致密.2011年7月施工至DK75+190附近发生大塌方,塌方段里程DK75+190～DK75+260,长约70 m.塌方段洞身位于膨胀性黄土地层中,隧道埋深约为28 m,位于偏压地带.隧道地质纵剖面图及塌方段地质断面图如图1所示.

2.2 基本物理力学参数

根据地质勘探和土工试验测定塌方段膨胀性黄土膨胀潜势分级为中等,试验结果见表1.塌方段土体物理力学参数见表2,其中砂质黄土初始含水率为17.9%,饱和含水率为28.5%.膨胀性黄土的初始含水率为16.7%,塌方后含水率为22.8%.

(a) 隧道地质纵剖面图(单位：m)

(b) 隧道塌方段地质断面图

3 数值计算方法

3.1 基质吸力和渗透系数确定

基质吸力采用滤纸法测定,通过对滤纸吸力的率定结果推算出土的基质吸力大小[15].将试验结果通过式(7)拟合,膨胀性黄土的土-水特征曲线如图2所示.

数值计算中膨胀性黄土渗透系数取值通过式(9)确定.饱和渗透系数ksat取值见表2.

3.2 黏聚力和内摩擦角的确定

通过室内固结慢剪试验获得膨胀性黄土不同含水率下的黏聚力和内摩擦角,试验结果如图3所示.应用式(10)和(11),得黏聚力c、内摩擦角φ与含水率ω的关系式：

lgc=-0.118 6ω+4.223 1

(13)

φ=-2.264 2ω+65.092 4

(14)

3.3 膨胀力的确定

根据不同初始含水率配制膨胀性黄土试样,采用恒体积试验法研究膨胀性黄土的膨胀力特性.依据该膨胀性黄土含水率分布范围16.7%～22.8%,制作不同初始含水率试样5组,每组2个试样.按照《土工试验规程》(SL237—1999)中膨胀力试验操作方法测定试样膨胀力,试验结果如表3所示.数值计算中,假定该膨胀性黄土含水率达22.8%时,其膨胀力发挥了100%.膨胀力与含水率的函数关系由式(12)确定,式中Ps=180 kPa.

表1 膨胀土膨胀性试验结果统计

表2 各土层物理力学性质参数

图2 滤纸法测定的膨胀性黄土的土-水特征曲线

图3 膨胀性黄土内摩擦角、黏聚力与含水率关系

含水率/%膨胀力/kPa试样1试样2均值14.14100.65102.80101.7316.51112.47111.28111.8818.95126.83127.11126.9720.33147.24146.37151.8023.21176.88179.63178.26

3.4 数值模型与边界条件

以隧道塌方段为对象建立三维数值模型,模型的左、右边界以实际地貌的山脊和沟谷为界,距离洞身约为60 m,下边界取5倍隧道的开挖直径,距离洞身约为50 m,上边界为地表,数值模型如图4所示.坐标原点位于仰拱横断面中点,x轴沿水平方向指向坡底,y轴沿竖直方向向上,z轴沿隧道轴向向内.围岩力学模型采用摩尔-库仑弹塑性本构模型.

图4 数值模型

模型左、右、前、后边界和下边界均为法向约束,上边界为自由边界.非饱和渗流数值计算中合理的水力边界条件决定模型的可靠性和准确性.小河沟隧道地处黄土山地地区,属温带半干旱季风性气候，年均降雨量约为428 mm，地下水位较低.因此数值模拟中不考虑地下水的影响,模型左右和底部设为不透水边界；自2011年7月份进入雨季以来,本地连续降雨,其中28～29日日降雨量较大，达到37.0～49.3 mm,至30日DK75+220～DK75+235里程段洞内初期支护纵向开裂变形,并有水滴.31日凌晨,DK75+190～DK75+260段发生大体积塌方.数值模拟塌方前60 h内隧道围岩破坏演化规律,此时塌方段已经有较大降雨量,因此数值模拟过程中假定地表一直处于饱和状态,故将地表节点饱和度设置为1,并将地表处的孔隙水压力值设置为0.流体模型采用均质各向同性流体本构模型.

隧道该塌方段支护形式采用V级围岩加强(黄土)复合式衬砌.数值计算中,将由钢拱架、钢筋网和喷射混凝土构成的初期支护和临时仰拱按照弯矩等效原理,模拟为相应厚度壳单元材料[16],等效的壳单元力学参数如表4所示.

表4 等效的壳单元物理力学参数

3.5 数值计算

数值计算中首先进行初始地应力平衡,然后进行隧道开挖模拟.开挖平衡后的竖向地应力场如图5所示,图中S2表示土体竖向应力.

隧道开挖施作初期支护至稳定后,开展降雨入渗增湿条件下膨胀性黄土隧道围岩变形演化过程模拟,采用流固耦合模式,在模型上边界设置相应的降雨入渗边界条件.数值模拟的降雨入渗过程是一个瞬态变化过程,不同的时间塌方段土层具有不同的含水率分布,也对应着不同应力状态.由于土层中初始的含水率分布未知,假设每层土中的含水均匀分布.图6～图8分别为降雨持续作用2, 12, 36, 60 h时数值计算结果对应的含水率、合位移、等效应力分布云图.图7中D表示合位移,图8中S表示等效应力.分析图6～图8,降雨持续至60 h过程中,隧道洞身及以上膨胀性黄土和砂质黄土含水率逐步增大,膨胀性黄土和砂质黄土自重增加且膨胀性黄土增湿膨胀,此外洞身处于偏压地带,左拱腰以上位置岩体自重较大,右拱肩以上位置埋深较浅,致使对应位置的围岩变形和洞周等效应力不断增大,最大位移约为22 cm,最大等效应力约为1.2 MPa.

图5 开挖后竖向地应力云图

(a) 降雨持续2 h

(b) 降雨持续12 h

(c) 降雨持续36 h

(d) 降雨持续60 h

图6 降雨持续作用下含水率分布变化

(a) 降雨持续2 h

(b) 降雨持续12 h

(c) 降雨持续36 h

(d) 降雨持续60 h

图7 降雨持续作用下合位移分布变化

(a) 降雨持续2 h

(b) 降雨持续12 h

(c) 降雨持续36 h

(d) 降雨持续60 h

图8 降雨持续作用下等效应力分布变化

4 结果分析与讨论

4.1 隧道围岩变形情况对比

隧道围岩变形随含水率变化情况的数值计算结果如图9(a)所示,隧道围岩变形现场实测结果如图9(b)所示.

(a) 数值计算结果

(b) 现场实测结果

分析图9可知,自隧道开挖完成至围岩变形稳定阶段,隧道拱顶沉降及拱腰水平收敛值实测与数值计算结果相差不大,二者拱顶沉降均约为5 cm,拱腰水平收敛位移均约为2.5 cm.该结果表明数值模型可有效模拟隧道开挖,为考虑降雨入渗作用下膨胀性黄土隧道围岩破坏演化研究奠定了基础.

从图6、图9(a)可看出,在降雨开始后约12 h内,隧道拱顶含水率几乎无变化.12 h后其含水率随时间呈现持续增大趋势,降雨至25 h后,拱顶含水率增大至21%,且增速开始减缓；此时,隧道拱顶沉降和拱腰水平收敛位移不断上升到达稳定值,并在较长时间内保持稳定.降雨持续50 h后,隧道拱腰至拱顶位置含水率增至约23%,且几乎不再增大；此时,隧道拱顶沉降和拱腰水平收敛位移急剧增大.降雨持续至60 h,拱顶沉降最大值约为21.7 cm,拱腰水平收敛最大值约为17.4 cm；对比分析图9(b)可看出,2011-07-31隧道塌方前,实测围岩拱顶沉降和水平收敛值也急剧增大,拱顶沉降最大值约为18.1 mm,拱腰水平收敛最大值约为14.8 mm.计算值与实测值偏差分别为16.6%和15.4%.综合以上结果,隧道围岩变形情况数值计算结果较实测值稍偏大,但总体规律吻合较好.原因是：① 数值计算模拟持续较大降雨情况,未考虑降雨量的变化；② 隧道塌方前一段时间内,测试条件恶劣,无法测得该时段围岩变形值,因而最终隧道围岩变形实测值偏小.

4.2 隧道围岩受力情况对比

隧道开挖后,为防止由围岩应力释放引起的大变形甚至坍塌,应及时施作初期支护,使初期支护与围岩共同受力.因此,为准确获得隧道围岩受力情况,可在隧道初期支护内设置应力监测点.图10为数值模型与隧道初期支护实测测点布置情况,图中SW表示钢架外侧测点,SN表示钢架内侧测点.

(a) 数值模型监测点

(b) 钢拱架应力实测测点

图10 隧道围岩应力测点布设形式

图11为隧道开挖稳定和降雨持续至60 h后初期支护的弯矩计算结果.分析可知,隧道开挖平衡后,周边最大正弯矩(弯矩以洞内侧受拉为正)主要分布在直边墙的中部位置,拱部弯矩在左右拱肩处较大,而在拱顶位置弯矩较小,左拱肩位置为正弯矩,拱顶及右拱肩位置为负弯矩.降雨持续作用60 h后,与开挖平衡阶段相比,左右边墙及右拱腰弯矩增大明显,其中右边墙正弯矩增大约2倍.图12为隧道开挖稳定和降雨持续至60 h后初期支护的轴力计算结果,图中轴力以压为正.分析可知,2种工况下初期支护各点轴力增量相差不大,右边墙和左拱腰位置的轴力较大.

图11 隧道初期支护弯矩计算值分布曲线

图12 隧道初期支护轴力计算值分布曲线

现场测试位置为2个断面,断面Ⅰ为塌方前正常里程段DK75+382,断面Ⅱ为塌方区域里程段DK75+220,两断面钢拱架应力监测结果如图13所示(钢架应力以压为负,拉为正).分析可知,在监测初期,断面Ⅱ实测钢架压应力与断面Ⅰ相差不大,随着降雨入渗影响,断面Ⅱ钢架压应力明显增大,隧道拱顶、拱肩及左右边墙压应力均增大3～4倍.从图13(b)可看出,2011-07-31隧道塌方前,拱顶及右拱肩位置内侧压应力逐渐超过外侧压应力,出现负弯矩区,钢架有向外挤出的趋势；左拱肩及左右边墙位置外侧压应力逐渐超过内侧压应力,出现正弯矩区,有向内挤出的趋势.钢架弯矩分布情况计算结果与实测结果相吻合.

(a) 隧道正常里程段断面Ⅰ

(b) 隧道塌方里程段断面Ⅱ

4.3 隧道围岩破坏机理分析

图14为持续降雨条件下隧道塌方前洞内初期支护破坏情况,隧道两侧边墙出现沿隧道轴线方向的裂缝、初支混凝土剥离等现象.

图14 隧道工程初期支护破坏情况

在降雨入渗条件下,膨胀性黄土围岩体增湿膨胀,膨胀力增大,围岩体普遍发生剪切破坏,因此剪应变增量变化在一定程度上可以表征围岩体变形发展情况,剪切应变增量大的位置更容易形成裂隙,甚至破坏.图15为降雨48 h与60 h后围岩体剪应变增量分布情况,图中Si表示剪切应变增量.可看出,降雨持续48 h隧道围岩剪切应变增量主要分布在左右边墙处和拱腰处;随着围岩体增湿程度进一步增大,剪应变增量迅速增大并向外延伸扩展;到降雨60 h后,围岩体出现较大的剪切滑动面,在左边墙中部和右边墙上侧剪切应变均十分显著,是造成围岩体塌方失稳的主要因素.这也较好地解释了图14中隧道初支破坏形式发生的原因.

(a) 48 h

4.4 地表变形破坏机理分析

图16为持续降雨条件下隧道塌方诱发地表山体滑塌情况.可见,由于黄土为垂直节理发育,隧道洞内塌方引起地表山体出现局部开裂、滑塌等灾害.

降雨60 h后,降雨入渗作用下膨胀性黄土隧道围岩大变形引起的地表山体塑性区分布情况如图17所示.分析可知降雨入渗作用下,下覆膨胀性黄土的山体边坡地表发育出大面积塑性区,且塑性区分布沿山坡走势连成一体,形成塑性贯通区.这有效揭示了降雨入渗作用下隧道塌方引起的地表开裂及塌陷等现象.

(a) 塌方段地表错动开裂

(b) 塌方段地表滑塌

图16 塌方里程段地表山体滑坡

图17 围岩体塑性破坏区数值计算结果

5 结论

1) 通过改进FLAC3D软件非饱和渗流计算模块,及渗流过程中土体抗剪强度参数和膨胀力随含水率变化关系,可较准确实现降雨入渗下膨胀性黄土隧道围岩大变形破坏演化分析.

2) 降雨持续作用下隧道围岩位移随含水率呈台阶形增长.降雨至25 h后,拱顶含水率增长至21%,此时,增速开始减缓,相应地隧道拱顶沉降和水平收敛位移出现极大值且趋于相对稳定.待降雨持续至50 h,隧道拱腰至拱顶含水率约为23%,此时,隧道位移急剧增大,最终引发围岩土体破坏,故在膨胀土隧道工程实践中应特别注意和重视含水率增长转折点,可以此作为安全施工的一项监控指标.

3) 降雨持续作用60 h后,隧道数值计算结果左右边墙弯矩增大明显、初期支护变形向内挤出,围岩体剪切应变发育显著,山体边坡发育出大面积塑性贯通区,有效揭示了降雨持续作用下洞内发生大塌方,且诱发洞外山体滑坡等现象.鉴于此,对于膨胀土隧道工程,支护结构建议采用刚度大、柔性好的支护方式；隧道洞顶地表及洞外边坡做好排水,必要时需做加固措施.

DOI:10.3969/j.issn.1000-7598.2009.01.001. Chen Zhenghan, Fang Xiangwei, Zhu Yuanqing, et al. Research on meso-structures and their evolution laws of expansive soil and loess [J].RockandSoilMechanics, 2009, 30(1): 1-11.

DOI:10.3969/j.issn.1000-7598.2009.01.001. (in Chinese)

[2]Rao S M, Reddy B V V, Muttharam M. The impact of cyclic wetting and drying on the swelling behaviour of stabilized expansive soils[J].EngineeringGeology, 2001, 60(1): 223-233. DOI:10.1016/s0013-7952(00)00103-4.

[3]曾仲毅. 降雨入渗下膨胀性黄土隧道围岩力学特性及稳定性分析[D].济南：山东大学土建与水利学院, 2014.

[4]曾铃, 史振宁, 付宏渊, 等. 降雨入渗对边坡暂态饱和区分布特征的影响[J].中国公路学报, 2017,30(1): 25-34. Zeng Ling, Shi Zhenning, Fu Hongyuan, et al. Influence of rainfall infiltration on distribution characteristics of slope transient saturated zone [J].ChinaJournalofHighwayandTransport, 2017, 30(1): 25-34. (in Chinese)

[5]李树忱, 晏勤, 谢璨, 等. 膨胀性黄土隧道钢拱架-格栅联合支护力学特性研究[J]. 岩石力学与工程学报, 2017, 36(2): 446-456.DOI:10.13722/j.cnki.jrme.2016.0599. Li Shuchen, Yan Qin, Xie Can, et al. The mechanical behavior of composite supports of steel-grid in loess tunnel [J].ChineseJournalofRockMechanicsandEngineering, 2017, 36(2): 446-456. DOI：10.13722/j.cnki.jrme.2016.0599.(in Chinese)

[6]王军, 曹平, 李江腾, 等. 降雨入渗对流变介质隧道边坡稳定性的分析[J]. 岩土力学, 2009, 30(7): 2158-2162, 2167.DOI:10.3969/j.issn.1000-7598.2009.07.050. Wang Jun, Cao Ping, Li Jiangteng, et al. Analysis of stability of tunnel-slope with rheological medium under rainfall infiltration [J].RockandSoilMechanics, 2009, 30(7): 2158-2162, 2167. DOI:10.3969/j.issn.1000-7598.2009.07.050. (in Chinese)

[7]陈伟, 阮怀宁, 张辉. 降雨入渗对浅埋偏压隧道及其支护系统的影响[J].现代隧道技术, 2008, 45(1): 30-34, 38.DOI:10.3969/j.issn.1009-6582.2008.01.006. Chen Wei, Ruan Huaining, Zhang Hui. Influence of rainfall seepage on the behavior of an evenly-pressured shallow tunnel and its support system [J].ModernTunnellingTechnology, 2008, 45(1): 30-34, 38. DOI:10.3969/j.issn.1009-6582.2008.01.006. (in Chinese)

[8]冯丙阳. 膨胀性黄土隧道大变形演化特征及支护对策研究[D]. 济南：山东大学土建与水利学院, 2014.

[9]毛昶熙. 渗流计算分析与控制[M]. 北京: 中国水利水电出版社, 2003:7-18.

[10]廖红建, 王铁行. 岩土工程数值分析[M].北京:机械工业出版社, 2006:225-230.

[11]Fredlund D G, Xing A Q. Equations for the soil-water characteristic curve[J].CanadianGeotechnicalJournal, 1994, 31(4): 521-532. DOI:10.1139/t94-061.

[12]Nguyen H V, Durso D F. Absorption of water by fiber webs: An illustration of diffusion transport [J].ProceedTechnicalAssocofthePulp&PaperIndustry, 1983, 66(12): 76-79.

[13]贾文聪, 李永红, 党进谦, 等. 膨胀土强度特性的试验研究[J]. 西北农林科技大学学报(自然科学版),2015(5): 217-221,228. DOI:10.13207/j.cnki.jnwafu.2015.05.016. Jia Wencong, Li Yonghong, Dang Jinqian, et al. Strength characteristics of expansive soil [J].JournalofNorthwestA&FUniversity(NaturalScienceEdition), 2015(5): 217-221, 228.DOI：10.13207/j.cnki.jnwafu.2015.05.016. (in Chinese)

[14]丁振洲, 郑颖人, 李利晟. 膨胀力变化规律试验研究[J]. 岩土力学, 2007, 28(7): 1328-1332.DOI:10.3969/j.issn.1000-7598.2007.07.008. Ding Zhenzhou, Zheng Yingren, Li Lisheng. Trial study on variation regularity of swelling force [J].RockandSoilMechanics, 2007, 28(7): 1328-1332. DOI:10.3969/j.issn.1000-7598.2007.07.008. (in Chinese)

[15]白福青, 刘斯宏, 袁骄. 滤纸法测定南阳中膨胀土土水特征曲线试验研究[J]. 岩土工程学报, 2011, 33(6):928-933. Bai Fuqing, Liu Sihong, Yuan Jiao. Measurement of SWCC of Nanyang expansive soil using the filter paper method [J].ChineseJournalofGeotechnicalEngineering, 2011, 33(6): 928-933.(in Chinese)

[16]曾仲毅, 徐帮树, 胡世权, 等. 增湿条件下膨胀土隧道衬砌破坏数值分析[J]. 岩土力学, 2014, 35(3): 871-880. Zeng Zhongyi, Xu Bangshu, Hu Shiquan, et al. Numerical analysis of tunnel liner failure mechanism in expansive soil considering water-increased state [J].RockandSoilMechanics, 2014,35(3): 871-880.(in Chinese)