基于蒙特卡洛方法的降雨条件下非饱和土边坡稳定性研究

王东

（中国铁路上海局集团有限公司工务部，上海 200071）

我国南方铁路路堤多为非饱和土质边坡。该类型路基边坡在降雨过程中，随着土体含水率和孔隙水压力增大，土层基质吸力丧失，抗剪强度不断降低，抗滑力减弱，易发生边坡溜坍、浅层滑坡等失稳性破坏，给铁路运输安全带来严重威胁。

降雨入渗过程受降雨模式和边坡土体水力特性的影响。非饱和土体中土壤的渗透率与含水率成正相关，降雨强度超过土壤饱和渗透率时，土壤的渗透率限制了雨水的入渗量，部分雨水以坡面径流的形式排出，另一部分入渗到土体中的水分增强了土体含水率和孔隙水压力，滑坡体重量增加，下滑力增大，同时伴随这一过程土体的抗剪强度不断降低，抗滑力减弱，进而导致边坡失稳。国内外学者对此开展了大量研究。何忠明等［1］通过室内模拟实验发现不同降雨模式下，路堤边坡中的水分分布特征有明显差异。李丞等［2］研究发现降雨强度对渣土场边坡的渗透特征和稳定性有显著影响。一些学者将非饱和土强度理论与极限平衡法结合求解边坡非饱和土边坡稳定性。卢玉林等［3］通过渗流、地震与极限平衡法的耦合研究了砂土边坡稳定性；刘子振等［4］讨论了非饱和土基质吸力和渗流力对边坡稳定性的影响。然而自然边坡岩土体离散性大［5］，依靠单一土体参数的现场实验和数值模拟研究难以获得边坡的失效概率，为解决这一问题，一些学者利用基于随机抽样的蒙特卡洛方法求解边坡的失效概率［6］。李侃等［7］通过Matlab 和Geoslope实现了边坡可靠度评价；张丽波等［8］利用蒙特卡洛法对边坡内摩擦角、黏聚力进行抽样，并结合极限平衡法计算了边坡的可靠性。

目前关于降雨条件下边坡稳定性问题已经开展了大量研究，但针对不同降雨强度与土壤水力特性下非饱和土边坡失稳问题的研究较少，本文选取三种不同水力参数的非饱和边坡进行数值模拟，研究其在不同降雨条件下的水分传递规律及边坡稳定性，并通过蒙特卡洛方法研究边坡体在不同强度参数下的失效概率，最后对土体黏聚力、内摩擦角和土体重度进行敏感性分析。

1 试验方法

1.1 非饱和土理论

Richards 方程是描述饱和-非饱和渗流时最常用的方程之一，基于单元体的连续性方程和广义Darcy定律建立。二维渗流的微分方程，计算公式为

式中：H为总水头；kx为x方向渗透系数；ky为y方向渗透系数；Q为边界流量；θ为单位体积水量；t为时间。

非饱和土坡中地下水位之上的孔隙水压力相对于大气压为负值，其强度不同于饱和土，本文稳定性计算采用目前被广泛使用的Fredlund理论。

Fredlund 抗剪强度公式采用双应力状态变量（σ-ua）和（ua-uw）来描述非饱和土抗剪强度τf，计算式为

式中：c′为有效黏聚力，为Mohr-Coulomb破坏包线的延伸线与剪应力轴的截距；σf为破坏时在破坏面上的法向总应力；ua，uw分别为破坏时破坏面上的孔隙气压力和孔隙水压力，在实际分析中，通常可认为ua= 0；φ′为与净法向应力状态变量(σf-ua)f有关的内摩擦角；φb为抗剪强度随基质吸力而增加的速率。

1.2 蒙特卡洛可靠性理论

蒙特卡洛法是基于随机数的概率方法，其在Slope/w 中实现过程为：①已知输入变量的概率分布，程序根据输入变量的概率分布进行随机抽样得到样本S=（S1，S2，…，SN），N为样本数量；②将样本S带入边坡安全系数计算模型F（S）中计算边坡安全系数；③统计所有安全系数小于1的样本n，计算边坡的失效概率Pf。

1.3 数值计算模型

1.3.1 数值模型参数

GEO-Studio 是一种适用岩土工程渗流与稳定性分析的有限元软件，其内置的GEO-Seep 模块能够较好地模拟降雨条件下边坡体内部渗流场的动态变化，配合GEO-Slope模块可以定量分析边坡稳定性。

数值计算模型如图1 所示：长度50 m，高度16 m，模型边坡坡高8 m，坡角45°。AH，BC分别为10 m 和3 m 的水头边界，底部AB、左右侧的GH，DC为不透水边界，顶部DE，EF及FG设置为随时间变化的流量边界模拟降雨。同时在长度20.5 m 和32.5 m 处布置了两个监测断面分析研究水分的变化。网格类型采用软件内置的四边形和三角形单元，网格尺寸为0.5 m，共计2561节点，2437个单元。

图1 模型边坡网格划分

根据文献选取了三种不同渗透特性的土体：黏性土、粉质砂土［9］和壤质砂土，其土壤水土特征曲线SWCC 和相对渗透率曲线见图2，黏性土渗透系数最低，壤质砂土最大。

图2 模型土壤水土特征曲线和相对渗透率曲线

为讨论渗透特性对边坡稳定性的影响，模型边坡土体被设置为相同的初始力学参数（表1）。模型计算中非饱和土体强度使用体积含水率进行计算得到。

表1 土壤物理参数

1.3.2 计算工况

综合我国气象降雨分级和边坡失稳案例对模型边坡施加10，20，40 mm/h 三种不同强度的降雨，降雨持时选取为6 h，降雨量分别为60，120，240 mm，对应气象分级分别为暴雨、大暴雨和特大暴雨。

表2 模拟工况

1.3.3 基于正态分布的参数分布

通过将粉质砂土的内摩擦角、黏聚力和饱和重度作为服从正态分布的随机变量进行边坡稳定性分析，利用蒙特卡洛法研究土体参数对边坡可靠性的影响，其参数取值见表3，模拟计算次数为2000。

表3 粉质砂土取值参数

2 结果与讨论

2.1 非饱和土边坡稳定性

降雨条件下模型边坡稳定性见图3。可知：黏性土边坡在三种类型降雨条件下安全系数的变化一致，初始安全系数最大为1.162，最低安全系数1.147出现在720 h 处。粉质砂土边坡在暴雨和大暴雨工况下安全系数的变化与黏性土类似，安全系数随时间不断降低，最小值1.09出现在720 h处，特大暴雨工况下安全系数最小值1.03出现在降雨结束时刻，随后安全系数有所回升。壤质砂土边坡的安全系数在三种降雨下都表现出先快速降低后升高的变化，其中特大暴雨工况中边坡最低安全系数为0.877，低于稳定性限值1，表明边坡发生失稳。

为验证本文所采用的Geostudio 程序的有效性，使用Plaxis 软件采用的强度折减法对案例边坡的安全系数和滑动面位置进行计算。本文使用的极限平衡法所计算的边坡临界滑动面与Plaxis 软件确定的滑动面位置基本吻合，见图4。

图4 工况6位移及极限平衡法确定的临界滑动面

降雨6 h 时各工况安全系数见表4。可知，本文采用极限平衡法所计算的安全系数和强度折减法所得到的安全系数基本一致，误差小于2%。

表4 降雨6 h时各工况安全系数

为更细致地分析边坡的失稳模式，本文提取了各种工况下的临界滑动面。如图5 和图6 所示，临界滑动面与边坡土质相关，滑动面最大埋深出现在坡顶与坡面相交处，壤质砂土的滑动面最深，粉质砂土滑动面最浅。壤质砂土边坡滑动面的范围延伸至坡脚平台处，伴随降雨强度的增加滑动体的截面积从52.63 m2增加到61.03 m2，最大深度从6.510 m增至6.907 m。

图5 降雨6 h 后边坡临界滑动面深度

2.2 降雨条件下非饱和土边坡水分分布

图7 模型边坡1-1断面土壤饱和度

本文在坡顶和坡脚处设置两个监测断面以获得降雨作用下边坡土壤水分的变化。边坡顶部监测断面结果见图7，不同土质边坡体水分差异较大，黏性土边坡中土壤渗透性最差，水分的改变集中于土壤表层，6，24 h 水分的湿润峰分别位于0.84，0.93 m 深度处，由于土壤渗透率小于降雨强度，水分入渗受到土壤渗透率的控制，所以三种不同强度的降雨没有造成明显入渗差异；与黏性土相比，粉质砂土具有更大的土壤渗透率，暴雨和大暴雨工况下6 h 时湿润峰的推进深度达1.78 m，特大暴雨湿润峰达到了2.26 m 深度，降雨结束后伴随着水分的运移，坡面处土壤水分不断降低，地下水水位升高；降雨过程中壤质砂土边坡的湿润峰已经达到了地下水位处，不同等级的降雨显著影响6 h 时的水分剖面，降雨强度提高了坡顶饱和度，不同于黏性土和粉质砂土的中部小两端大的水分分布，壤质砂土呈现出上小下大的分布模式。

图8为位于坡脚处的2-2断面的水分剖面，可见三种降雨条件下黏性土边坡坡脚水分分布没有明显的差异，6 h时水分的改变集中在土体表层，湿润峰的深度为0.548 m；暴雨和大暴雨作用下粉质砂土中湿润峰在6 h 时达到1.069 m 深度，其下水分受到地下水位抬升的影响有稍微升高，而特大暴雨作用下坡脚土壤基本完全饱和；壤质砂土边坡坡脚土体基本保持饱和状态。6 h至24 h这段时间，黏性土和粉质黏土边坡中都出现了地下水位的抬升，粉质黏土中这种抬升表现得更加明显。

图8 模型边坡2-2断面土壤饱和度

2.3 降雨条件下非饱和土边坡孔隙水压力变化

孔隙水压力直接影响非饱和土强度，通过对边坡稳定性的计算得到边坡各个时段临界滑动面，选取特大降雨工况提取滑动面上孔隙水压力绘制得到图9。黏性土边坡临界滑动面孔隙水压力的变化主要发生在坡顶和坡脚处，坡顶和坡脚处负孔隙水压力分别增加了32.3，19.4 kPa，其安全系数降低是因为滑动体重量增加和坡顶及坡脚处非饱和土剪切强度降低；粉质黏土中滑动面孔压的变化主要分为两个阶段：第一阶段坡顶和坡脚处的孔隙水压力快速升高，降雨停止后进入第二阶段，滑动面顶部和底部孔隙水压力逐步降低，中部孔压不断升高，最终滑动面孔隙水压力稳定在-25～-15 kPa。伴随着雨水的入渗，壤质砂土边坡临界滑动面上孔隙水压力快速升高，并在滑动面中部出现正孔隙水压力，造成土体有效应力降低，诱发边坡失稳，孔隙水压力极值25.45 kPa 出现在降雨停止时，伴随着水分的入渗滑动面上的孔隙水压力不断降低。

图9 临界滑动面孔隙水压力

2.4 非饱和土边坡可靠性及敏感性分析

本文选取特大暴雨工况下粉质黏土边坡进行可靠性分析和参数敏感性实验。边坡在不同土体参数条件下的安全系数概率分布见图10，黏聚力、内摩擦角和饱和重度参数改变诱发的边坡失效概率分别为5.75%，0.60%，0.05%。

图10 概率分布

为量化不同土体参数对边坡稳定性的影响，本文计算了不同土体参数下边坡稳定性安全系数，结果见图11。黏聚力和内摩擦角与边坡稳定性安全系数成正相关，其中黏聚力从0.5 kPa 增加到5.5 kPa 的过程中边坡安全系数从0.878 增加到了1.217。土体饱和重度的增加对边坡的安全性造成消极影响，饱和重度为20.5 kN/m3时，边坡安全系数为0.98。

图11 安全系数敏感度

3 结论

本文通过数值模型试验对不同降雨和土质条件下非饱和土边坡稳定性进行分析，得出如下结论：

1）在降雨条件下非饱和土质边坡稳定性变化有显著差异，高渗透率的壤质砂土边坡安全系数的极小值出现在降雨结束时刻，而低渗透率的黏性土和粉质黏土边坡的稳定性安全系数在降雨之后仍持续降低。

2）黏性土边坡渗透性较低，降雨过程中滑动面上孔隙水压力的增长主要发生在滑动面顶部和底部，该部分基质吸力的丧失降低了土体边坡稳定性。

3）高渗透性土体的边坡失稳模式是降雨导致的地下水位抬升诱发边坡失稳，其规模和影响深度大于低渗透性边坡。

4）边坡土体物理力学参数影响边坡稳定性，通过可靠性分析发现，土体的黏聚力和内摩擦角与边坡稳定性安全系数成正相关，与土体饱和重度成负相关。