强降雨入渗作用下土质边坡失稳灾变研究

,　　,　

(1.广西大学 机械工程学院, 广西 南宁　530007；2.广西交通科学研究院，广西 南宁　530007；3.广西翔路建设有限责任公司, 广西 南宁　530029；4.广州大学 工程抗震研究中心，广东 广州　510405)

0　引言

降雨入渗是诱发土质边坡失稳最主要的原因之一[1]。对于风化的深度较浅的边坡，降雨时大多失效仅发生在表层的风化残积土[2]。由于造价的因素，大量土质路堑边坡在公路施工期及运营期均裸露在外或者以简单的植草进行防护。在较长时间阳光、风、降雨等自然因素的影响下，其表层土体风化程度较大，细粒结构因渗流因素大量流失，土体性状相较深层土体有较大差异，形成典型的薄壳双层土体结构。在遭遇连续降雨后，容易发生以浅层滑坡为主的坡体失稳，滑动类型一般为小型的平面或非圆弧滑动和崩塌，风化表层土与深层土体结合部经常构成滑动面，“湿润峰”深度一般不超过2 m[3]。这种层状土入渗特性不同于均质土的入渗特性，层状土的降雨入渗分析对研究多层土的入渗分析具有重要意义。

国内外许多学者对边坡的降雨入渗过程进行大量研究，建立相应的求解方法。朱伟等[4]研究设计室内降雨入渗土柱试验，在初步揭示降雨入渗过程和规律的基础上，探讨准确反映降雨入渗量的有限元计算方法；张社荣等[5]研究了强降雨特性对饱和-非饱和边坡失稳的影响；李宁等[6]用Python语言对ABAQUS软件的降雨入渗边界进行二次开发，开发出基于ABAQUS软件的降雨入渗模块；吴李泉等[7]考虑降雨入渗引起土体质量和渗透力增大及抗剪强度对边坡稳定性的影响，开发USLOPEFEM计算程序，研究强降雨入渗条件下浙江武义平头村山体高边坡的瞬态渗流场及稳定性；W.T.Oh和S.K.Vanapalli[8]研究降雨入渗对土体抗剪强度及边坡稳定性的影响；Guan， J. et al[9]综合分析了地表径流、渗流作用与边坡稳定的关系；Oh S和Lu N[10]研究发现尽管边坡几何、流体力学的特性、抗剪强度和降水时间会有差异，但边坡失效时其安全系数都小于1.0。上述研究表明降雨诱发滑坡的机制主要是由于降雨入渗造成岩土体基质吸力降低、孔隙水压力变化及地下水位的上升，并增加岩土体的容重、降低岩土体的抗剪强度。但是，降雨对于深层滑坡，主要地下水位的抬升，以及水的弱化作用对滑坡的影响较大；而对于浅层滑坡，主要受瞬时孔隙水压力的变化控制，同时伴随着冲刷和侵蚀作用。

上述研究大多假设边坡由均质土构成，对层状边坡的降雨入渗研究甚少。因此以广西百色一典型公路边坡为例，利用公路运营前设计阶段的边坡勘察资料和公路运营期边坡补充地质勘查和土工试验资料，分析了原状土体和公路运营后边坡土体的变化，并利用有限元程序ABAQUS分析了双层土质高边坡降雨入渗全过程，以更好地理解降雨对滑坡的作用过程和诱发机制，为边坡病害防治和预测预报提供科学依据。

1　研究边坡的选择与概况

图1　目标高边坡工程地质图(单位：cm)

选择广西百色至隆林高速公路田林县境某路堑高边坡作为研究对象。图1为高边坡的工程地质图，2006年设计阶段进行了一系列的现场勘查工作,钻孔分别为ZK-1与ZK-2，得到了相应岩土物理力学参数。施工开挖后只剩下强风化泥质粉砂岩的单一均质土层。2014年为研究边坡失稳，再次对边坡进行了补充勘查，钻孔分别为ZK-3与ZK-4，重点研究边坡的风化程度和风化后土层的力学参数。钻探结果显示，边坡表层平均深2 m左右的土，由于风化作用，表现出与2006年钻探结果中全风化泥质粉砂岩III级黏土类似性质，而深层硬土的力学参数变化很小。具体岩土物理力学参数如表1所示，因土的物理参数随深度而发生细微变化，为简化计算，表中数值经过均值处理。

表1　各地层岩土物理力学参数

2　边坡计算条件

2.1　计算模型及参数

目标边坡数值计算模型如图2所示，边坡坡度为42°，表层风化土为厚2 m的全风化泥质粉砂岩，下部为强风化泥质粉砂岩，土体的物理力学参数如表1所示。A、B为监测单元，C、D、F为监控点，利用ABAQUS有限元程序强大的渗流与变形耦合分析功能，进行边坡稳态渗流场与应力场的耦合仿真分析，以探求不同强降雨作用下边坡失稳的内在规律。分析中边坡的非饱和水力参数参考相关研究[11]类比确定，计算所采用的孔压-相对渗透系数的关系如图3所示。位移边界条件为模型左右两侧施加法向约束，边坡上部为自由边界，下部为全约束边界；流量边界条件为模型两侧、地下水位以下考虑为定水头边界，地下水位以上设为零流量边界，模型底面为不透水边界。模型采用孔压/位移耦合的CPE4P四节点平面应变单元，单元数为1 186个。本构模型和屈服准则分别采用理想弹塑性本构和Mohr-Coulomb强度准则，计算模型内所有节点的应力、位移及单元的塑性应变。

2.2　降雨入渗计算条件

参考勘查报告，模型两侧初始水头边界为8 m。参照对百色市1961—1990降雨统计资料。确定模拟降雨的过程为：模拟总时间为3 d，降雨总量为288 mm。因实际降雨情况非常复杂，难以逐一讨论。本文降雨形式采用全局平均方法，将实际降雨理想化地划分为3次直线上升下降型、直线上升型、单次直线上升下降型，以表示降雨强度在降雨过程中的变化情况。则3种型式的降雨分布如图4所示。

图2　边坡有限元模型

图3　边坡孔压与相对渗透系数的关系

图4　降雨雨型

3　计算结果与分析

3.1　降雨前后孔压与应力对比分析

降雨前后边坡的孔隙水压力(POR)分布如图5所示。图5(a)中，边坡的水压力呈线性分布，这与本文给出的初始假设孔压随深度线性分布是相符合的。对比图5(a)与图5(b)、(c)可知，降雨入渗后的孔压分布图与初始状态有明显的区别，边坡顶部以下的吸力区范围减小。由图5(b)、(c)可知，随着降雨时间的延长，边坡的饱和度增大，孔隙水压力增大，土体浅层的基质吸力则减小或消失。验证了模拟强降雨方法的可靠性，能在此基础上进行降雨特性对边坡破坏的深入研究。

图5　孔隙水压力分布图

降雨前后边坡的竖向有效应力(S22)分布如图6所示。图6(a)中，降雨前边坡的顶部竖向有效应力并不为零，是由于Abaqus/Standard有效应力中考虑了吸力的影响。另外竖向有效应力分布呈现从坡面向里逐渐增加的特点，这与水平地基的应力分布情况是截然不同的。对比图6(a)与图6(b)、(c)可知，降雨入渗后的边坡的竖向应力发生了较大的改变，边坡内部最大竖向应力逐渐增加，边坡表层与土层接触位置的竖向应力都随着降雨时间增加逐渐减小。这也正好说明了，随着降雨时间的增长，边坡可能在浅层发生破坏。

图6　竖向有效应力分布图

3.2　降雨雨型对边坡响应的影响

图7为单元A在雨型1～3作用下塑性应变(ε)随时间的变化曲线，A1表示单元A在雨型1作用下，A2、A3类似。单元A在雨型1～3作用下时，分别当t=7.6 h、t=43.7 h与t=17.2 h时，坡脚出现塑性区，之后都随着时间的增长而快速增加，大约到t=12 h、36 h和72 h时，此时3种降雨强度首次达到最大值，坡脚的塑性应变也分别达到最大值0.067、0.069和0.051，对于雨型1和2塑性应变达到峰值之后不再随时间增长而增加。雨型1作用下，边坡很快就出现塑性应变。雨型2作用时，虽塑性应变出现的时间较晚，但峰值最大。雨型3作用下，边坡的塑性应变峰值最小。由此可知，对于不同雨型，边坡塑性应变峰值都是在降雨达到最大时出现，当降雨强度达到最大后，塑性应变不再增加；降雨总量相同的情况下，雨型会影响塑性应变出现的时间与峰值。

单元A与单元B的平均应力路径如图8所示，横坐标与纵坐标分别为单元的平均应力(P)与等效应力(σ)。图8中，单元A与单元B的平均应力路径出现明显不同的特征，雨型1～3作用下单元A的平均应力的路径基本一致，整体为开口向右下的斜U字型。这是因为坡脚的A单元，在降雨入渗的作用下，孔压增加，平均应力是减小的。当减小到一定程度时，平均应力路径达到屈服面，此时应力路径沿着屈服面向左下方移动。但这之后对于雨型1与雨型3，由于降雨强度的减少，边坡吸力增加，孔压减小，平均应力增大后，会逐渐偏移屈服面。而雨型2，降雨强度不会减少，因而不会偏移屈服面。雨型1、3作用下单元B的平均应力的路径也基本一致，平均应力与等效应力先增大，到一定值后，又分别减小。这是由于单元B处于土坡内部，上方的土吸水后容重增加，导致单元B的平均应力与等效应力都会增加，但对雨型1、3接近降雨结束时，单元B的平均应力与等效应力会有所下降。但雨型2作用时，单元B的平均应力与等效应力一值是上升的。由于该双层结构土质路堑高边坡表层A单元更容易进入屈服，这正好说明，降雨入渗作用下边坡的表层可能先出现失稳现象。由此可见，应对实际工程中的高边坡，采取一定的措施，增加边坡的稳定。

图7　单元A塑性应变随时间的变化

图8　单元A与单元B平均应力路径

图9为雨型1～3作用下边坡的位移时程曲线，U11表示雨型1作用下边坡的水平位移U1，U21表示雨型1作用下边坡的竖向位移U2。由图9可知，对于该双层结构土质路堑高边坡，相同降水总量的情况下，降雨形式对边坡表层监测点的位移响应影响较大。边坡在各雨型作用初期，降雨强度较低时，边坡竖向位移都会略微向上，水平位移会略微水平向左，这是由于初期边坡土体吸水而发生膨胀。对于边坡表层监测点C而言(图9(a))，各雨型降水强度首次接近峰值时，边坡的水平与竖向位移都急速增大。这是由于降雨强度的增大，导致双层结构土质路堑高边坡的坡体内部的土体孔压增大，吸力降低，土体内部发生软化，导致土体的抗剪强度降低，加之坡体本身自重作用下，边坡的变形会迅速加大。图9(a)中，雨型1、3作用下时，边坡水平与竖向位移达到峰值后，随着时间的增长，位移的变化很小。监测点C在雨型1作用下的竖向位移最大，为-0.34 m，雨型2次之，雨型3作用下最小。监测点C在雨型2作用下的水平位移最大，为-0.22 m，雨型1次之，雨型3作用下最小。对于边坡表层监测点D而言(图9(b))，各雨型降水强度首次接近峰值时，边坡的水平与竖向位移都急速增大。原因同监测点C。图9(b)中，雨型1～3作用下时，边坡水平与竖向位移变化趋势同监测点C。监测点D在雨型1作用下的竖向位移最大，为-0.31 m，雨型2次之，雨型3作用下最小。监测点D在雨型2作用下的水平位移最大，为-0.40 m，雨型1次之，雨型3作用下最小。这里监测点C与D的位移都较大，此时边坡位移较大，边坡已处于濒临失稳的边缘。图9(c)中，监测点F位移响应与雨量的变化较为一致；边坡的竖向位移都为向上，水平位移都是向左，这是由于边坡在降雨的作用下坡体软化和边坡自重的增加，进而边坡的滑动力将增大，该滑动力传递至边坡下部并挤压坡脚土体，从而导致该点会抬升，同时向左滑移。总体而言，不同雨型作用下，双层土边坡的不同部位的位移响应差别较大，降雨对边坡中部的水平位移影响最大，上部次之，底部最小；降雨对边坡上部的竖向位移影响最大，中部次之，底部最小。可见对双层结构土质路堑高边坡在连续降雨作用下的失稳，应予以重视，在实际工程中应采取相应的措施，增强边坡的稳定性。

图9　雨型1～3作用下边坡的位移时程曲线

3.3　降雨对边坡失稳的影响

为了探讨强降雨对边坡失稳的影响，在雨型2的基础上，降雨总量考虑为720 mm的极端情况。图10(a)为边坡接近破坏时位移等值云图，此时边坡的最大位移为0.463 m，在边坡坡脚表层土位置，然而，由于位移中包含了由重力荷载所引起的位移，根据总位移等值线图无法判断滑动面的位置。但考虑到边坡失稳的趋势应该在增量位移中得到反映，因而，可以通过计算终止的最后一个增量步的增量位移来判断出边坡滑动面的位置。边坡最后一个增量步的增量位移如图10(b)所示，最大位移为0.025 m，该位置边坡的变形速度为0.054 m/h，边坡已发生紧急变动，说明对于该双层结构土质路堑高边坡，在降雨的持续作用下，最先发生失稳是该位置，实际中应采取相应的措施进行加固处理。其它雨型作用下，边坡的失稳的过程与雨型2类似。

图10　边坡位移等值云图

4　结论与建议

基于饱和非饱和降雨入渗模拟方法，采用ABAQUS软件，结合地质勘查资料和土工试验成果对不同雨型作用下的典型双层土质高边坡降雨入渗全过程进行了数值模拟，得到了如下结论：

(1) 在强降雨条件下，地表水入渗将导致土体浅层的基质吸力减小或消失，土体表层竖向有效应力减小。

(2) 不同雨型作用下，边坡表层土塑性应变峰值是在降雨达到最大时出现。在降雨入渗的作用下，边坡表层土的平均应力是减小的，当减小到一定程度时，平均应力路径达到屈服面，此时应力路径沿着屈服面向左下方移动，后期降雨强度减少，平均应力增大后，会逐渐偏移屈服面，但边坡内部平均应力与等效应力都会增加。不同雨型作用下，双层土边坡的不同部位的位移响应差别较大，降雨对边坡中部的水平位移影响最大，上部次之，底部最小；降雨对边坡上部的竖向位移影响最大，中部次之，底部最小。

(3) 双层土边坡失稳过程的塑性区向上延伸的方向有两个，一个是沿着两层土的过渡区向上，另一个是沿着弧线向上。但是，此类边坡主要体现在表层失稳，深层沿滑动面滑移的可能性较小。

(4) 对于这类土质边坡，在运营过程中如果能够加强山背截水沟、坡面排水等构造设计，加快排水的速度，减少降雨入渗的强度，对于防止土体表层失稳有比较明显的效果。

[1]Anagnostopoulos G G, Fatichi S, Burlando P. An advanced process-based distributed model for the investigation of rainfall-induced landslides: The effect of process representation and boundary conditions[J]. Water Resources Research, 2015, 51(9): 7501-7523.

[2]Cho S E. Probabilistic stability analysis of rainfall-induced landslides considering spatial variability of permeability[J]. Engineering Geology, 2014, 171: 11-20.

[3]Tu X B, Kwong A K，Tu X B, et al. Field monitoring of rainfall infiltration in a loess slope and analysis of failure mechanism of rainfall-induced landslides[J]. Engineering Geology, 2009, 105(1): 134-150.

[4]朱伟, 陈学东, 钟小春. 降雨入渗规律的实测与分析[J]. 岩土力学, 2006, 27(11)：1873-1879.

[5]张社荣, 谭尧升, 王超, 等. 强降雨特性对饱和-非饱和边坡失稳破坏的影响[J]. 岩石力学与工程学报, 2014, 33(S2):4102-4112.

[6]李宁,许建聪. 基于ABAQUS 的三维边坡降雨入渗模块的开发及其应用[J]. 岩土工程学报, 2015, 37(4): 667-674.

[7]吴李泉, 张峰, 凌贤长, 等. 强降雨条件下浙江武义平头村山体高边坡稳定性分析[J]. 岩石力学与工程学报, 2009, 28(6):1193-1199.

[8]Oh W T, Vanapalli S K. Influence of rain infiltration on the stability of compacted soil slopes[J]. Computers and Geotechnics, 2010, 37(5):649-657.

[9]Guan J, Mok C M, Yeung A T. Integrated analysis framework for predicting surface runoff, infiltration, and slope stability[C]//Geo-Congress 2014 Technical Papers@ Geo-characterization and Modeling for Sustainability. ASCE, 2014: 2588-2599.

[10]Oh S, Lu N. Slope stability analysis under unsaturated conditions: Case studies of rainfall-induced failure of cut slopes[J]. Engineering Geology, 2015, 184: 96-103.

[11]费康, 张建伟. Abaqus 在岩土工程中的应用[M]. 北京:中国水利水电出版社, 2010:222-233.