基于蒙特卡洛方法的降雨條件下非飽和土邊坡穩(wěn)定性研究

王東

（中國鐵路上海局集團有限公司工務(wù)部，上海 200071）

我國南方鐵路路堤多為非飽和土質(zhì)邊坡。該類型路基邊坡在降雨過程中，隨著土體含水率和孔隙水壓力增大，土層基質(zhì)吸力喪失，抗剪強度不斷降低，抗滑力減弱，易發(fā)生邊坡溜坍、淺層滑坡等失穩(wěn)性破壞，給鐵路運輸安全帶來嚴(yán)重威脅。

降雨入滲過程受降雨模式和邊坡土體水力特性的影響。非飽和土體中土壤的滲透率與含水率成正相關(guān)，降雨強度超過土壤飽和滲透率時，土壤的滲透率限制了雨水的入滲量，部分雨水以坡面徑流的形式排出，另一部分入滲到土體中的水分增強了土體含水率和孔隙水壓力，滑坡體重量增加，下滑力增大，同時伴隨這一過程土體的抗剪強度不斷降低，抗滑力減弱，進而導(dǎo)致邊坡失穩(wěn)。國內(nèi)外學(xué)者對此開展了大量研究。何忠明等［1］通過室內(nèi)模擬實驗發(fā)現(xiàn)不同降雨模式下，路堤邊坡中的水分分布特征有明顯差異。李丞等［2］研究發(fā)現(xiàn)降雨強度對渣土場邊坡的滲透特征和穩(wěn)定性有顯著影響。一些學(xué)者將非飽和土強度理論與極限平衡法結(jié)合求解邊坡非飽和土邊坡穩(wěn)定性。盧玉林等［3］通過滲流、地震與極限平衡法的耦合研究了砂土邊坡穩(wěn)定性；劉子振等［4］討論了非飽和土基質(zhì)吸力和滲流力對邊坡穩(wěn)定性的影響。然而自然邊坡巖土體離散性大［5］，依靠單一土體參數(shù)的現(xiàn)場實驗和數(shù)值模擬研究難以獲得邊坡的失效概率，為解決這一問題，一些學(xué)者利用基于隨機抽樣的蒙特卡洛方法求解邊坡的失效概率［6］。李侃等［7］通過Matlab 和Geoslope實現(xiàn)了邊坡可靠度評價；張麗波等［8］利用蒙特卡洛法對邊坡內(nèi)摩擦角、黏聚力進行抽樣，并結(jié)合極限平衡法計算了邊坡的可靠性。

目前關(guān)于降雨條件下邊坡穩(wěn)定性問題已經(jīng)開展了大量研究，但針對不同降雨強度與土壤水力特性下非飽和土邊坡失穩(wěn)問題的研究較少，本文選取三種不同水力參數(shù)的非飽和邊坡進行數(shù)值模擬，研究其在不同降雨條件下的水分傳遞規(guī)律及邊坡穩(wěn)定性，并通過蒙特卡洛方法研究邊坡體在不同強度參數(shù)下的失效概率，最后對土體黏聚力、內(nèi)摩擦角和土體重度進行敏感性分析。

1 試驗方法

1.1 非飽和土理論

Richards 方程是描述飽和-非飽和滲流時最常用的方程之一，基于單元體的連續(xù)性方程和廣義Darcy定律建立。二維滲流的微分方程，計算公式為

式中：H為總水頭；kx為x方向滲透系數(shù)；ky為y方向滲透系數(shù)；Q為邊界流量；θ為單位體積水量；t為時間。

非飽和土坡中地下水位之上的孔隙水壓力相對于大氣壓為負值，其強度不同于飽和土，本文穩(wěn)定性計算采用目前被廣泛使用的Fredlund理論。

Fredlund 抗剪強度公式采用雙應(yīng)力狀態(tài)變量（σ-ua）和（ua-uw）來描述非飽和土抗剪強度τf，計算式為

式中：c′為有效黏聚力，為Mohr-Coulomb破壞包線的延伸線與剪應(yīng)力軸的截距；σf為破壞時在破壞面上的法向總應(yīng)力；ua，uw分別為破壞時破壞面上的孔隙氣壓力和孔隙水壓力，在實際分析中，通常可認為ua= 0；φ′為與凈法向應(yīng)力狀態(tài)變量(σf-ua)f有關(guān)的內(nèi)摩擦角；φb為抗剪強度隨基質(zhì)吸力而增加的速率。

1.2 蒙特卡洛可靠性理論

蒙特卡洛法是基于隨機數(shù)的概率方法，其在Slope/w 中實現(xiàn)過程為：①已知輸入變量的概率分布，程序根據(jù)輸入變量的概率分布進行隨機抽樣得到樣本S=（S1，S2，…，SN），N為樣本數(shù)量；②將樣本S帶入邊坡安全系數(shù)計算模型F（S）中計算邊坡安全系數(shù)；③統(tǒng)計所有安全系數(shù)小于1的樣本n，計算邊坡的失效概率Pf。

1.3 數(shù)值計算模型

1.3.1 數(shù)值模型參數(shù)

GEO-Studio 是一種適用巖土工程滲流與穩(wěn)定性分析的有限元軟件，其內(nèi)置的GEO-Seep 模塊能夠較好地模擬降雨條件下邊坡體內(nèi)部滲流場的動態(tài)變化，配合GEO-Slope模塊可以定量分析邊坡穩(wěn)定性。

數(shù)值計算模型如圖1 所示：長度50 m，高度16 m，模型邊坡坡高8 m，坡角45°。AH，BC分別為10 m 和3 m 的水頭邊界，底部AB、左右側(cè)的GH，DC為不透水邊界，頂部DE，EF及FG設(shè)置為隨時間變化的流量邊界模擬降雨。同時在長度20.5 m 和32.5 m 處布置了兩個監(jiān)測斷面分析研究水分的變化。網(wǎng)格類型采用軟件內(nèi)置的四邊形和三角形單元，網(wǎng)格尺寸為0.5 m，共計2561節(jié)點，2437個單元。

圖1 模型邊坡網(wǎng)格劃分

根據(jù)文獻選取了三種不同滲透特性的土體：黏性土、粉質(zhì)砂土［9］和壤質(zhì)砂土，其土壤水土特征曲線SWCC 和相對滲透率曲線見圖2，黏性土滲透系數(shù)最低，壤質(zhì)砂土最大。

圖2 模型土壤水土特征曲線和相對滲透率曲線

為討論滲透特性對邊坡穩(wěn)定性的影響，模型邊坡土體被設(shè)置為相同的初始力學(xué)參數(shù)（表1）。模型計算中非飽和土體強度使用體積含水率進行計算得到。

表1 土壤物理參數(shù)

1.3.2 計算工況

綜合我國氣象降雨分級和邊坡失穩(wěn)案例對模型邊坡施加10，20，40 mm/h 三種不同強度的降雨，降雨持時選取為6 h，降雨量分別為60，120，240 mm，對應(yīng)氣象分級分別為暴雨、大暴雨和特大暴雨。

表2 模擬工況

1.3.3 基于正態(tài)分布的參數(shù)分布

通過將粉質(zhì)砂土的內(nèi)摩擦角、黏聚力和飽和重度作為服從正態(tài)分布的隨機變量進行邊坡穩(wěn)定性分析，利用蒙特卡洛法研究土體參數(shù)對邊坡可靠性的影響，其參數(shù)取值見表3，模擬計算次數(shù)為2000。

表3 粉質(zhì)砂土取值參數(shù)

2 結(jié)果與討論

2.1 非飽和土邊坡穩(wěn)定性

降雨條件下模型邊坡穩(wěn)定性見圖3?？芍吼ば酝吝吰略谌N類型降雨條件下安全系數(shù)的變化一致，初始安全系數(shù)最大為1.162，最低安全系數(shù)1.147出現(xiàn)在720 h 處。粉質(zhì)砂土邊坡在暴雨和大暴雨工況下安全系數(shù)的變化與黏性土類似，安全系數(shù)隨時間不斷降低，最小值1.09出現(xiàn)在720 h處，特大暴雨工況下安全系數(shù)最小值1.03出現(xiàn)在降雨結(jié)束時刻，隨后安全系數(shù)有所回升。壤質(zhì)砂土邊坡的安全系數(shù)在三種降雨下都表現(xiàn)出先快速降低后升高的變化，其中特大暴雨工況中邊坡最低安全系數(shù)為0.877，低于穩(wěn)定性限值1，表明邊坡發(fā)生失穩(wěn)。

為驗證本文所采用的Geostudio 程序的有效性，使用Plaxis 軟件采用的強度折減法對案例邊坡的安全系數(shù)和滑動面位置進行計算。本文使用的極限平衡法所計算的邊坡臨界滑動面與Plaxis 軟件確定的滑動面位置基本吻合，見圖4。

圖4 工況6位移及極限平衡法確定的臨界滑動面

降雨6 h 時各工況安全系數(shù)見表4?？芍?，本文采用極限平衡法所計算的安全系數(shù)和強度折減法所得到的安全系數(shù)基本一致，誤差小于2%。

表4 降雨6 h時各工況安全系數(shù)

為更細致地分析邊坡的失穩(wěn)模式，本文提取了各種工況下的臨界滑動面。如圖5 和圖6 所示，臨界滑動面與邊坡土質(zhì)相關(guān)，滑動面最大埋深出現(xiàn)在坡頂與坡面相交處，壤質(zhì)砂土的滑動面最深，粉質(zhì)砂土滑動面最淺。壤質(zhì)砂土邊坡滑動面的范圍延伸至坡腳平臺處，伴隨降雨強度的增加滑動體的截面積從52.63 m2增加到61.03 m2，最大深度從6.510 m增至6.907 m。

圖5 降雨6 h 后邊坡臨界滑動面深度

2.2 降雨條件下非飽和土邊坡水分分布

圖7 模型邊坡1-1斷面土壤飽和度

本文在坡頂和坡腳處設(shè)置兩個監(jiān)測斷面以獲得降雨作用下邊坡土壤水分的變化。邊坡頂部監(jiān)測斷面結(jié)果見圖7，不同土質(zhì)邊坡體水分差異較大，黏性土邊坡中土壤滲透性最差，水分的改變集中于土壤表層，6，24 h 水分的濕潤峰分別位于0.84，0.93 m 深度處，由于土壤滲透率小于降雨強度，水分入滲受到土壤滲透率的控制，所以三種不同強度的降雨沒有造成明顯入滲差異；與黏性土相比，粉質(zhì)砂土具有更大的土壤滲透率，暴雨和大暴雨工況下6 h 時濕潤峰的推進深度達1.78 m，特大暴雨濕潤峰達到了2.26 m 深度，降雨結(jié)束后伴隨著水分的運移，坡面處土壤水分不斷降低，地下水水位升高；降雨過程中壤質(zhì)砂土邊坡的濕潤峰已經(jīng)達到了地下水位處，不同等級的降雨顯著影響6 h 時的水分剖面，降雨強度提高了坡頂飽和度，不同于黏性土和粉質(zhì)砂土的中部小兩端大的水分分布，壤質(zhì)砂土呈現(xiàn)出上小下大的分布模式。

圖8為位于坡腳處的2-2斷面的水分剖面，可見三種降雨條件下黏性土邊坡坡腳水分分布沒有明顯的差異，6 h時水分的改變集中在土體表層，濕潤峰的深度為0.548 m；暴雨和大暴雨作用下粉質(zhì)砂土中濕潤峰在6 h 時達到1.069 m 深度，其下水分受到地下水位抬升的影響有稍微升高，而特大暴雨作用下坡腳土壤基本完全飽和；壤質(zhì)砂土邊坡坡腳土體基本保持飽和狀態(tài)。6 h至24 h這段時間，黏性土和粉質(zhì)黏土邊坡中都出現(xiàn)了地下水位的抬升，粉質(zhì)黏土中這種抬升表現(xiàn)得更加明顯。

圖8 模型邊坡2-2斷面土壤飽和度

2.3 降雨條件下非飽和土邊坡孔隙水壓力變化

孔隙水壓力直接影響非飽和土強度，通過對邊坡穩(wěn)定性的計算得到邊坡各個時段臨界滑動面，選取特大降雨工況提取滑動面上孔隙水壓力繪制得到圖9。黏性土邊坡臨界滑動面孔隙水壓力的變化主要發(fā)生在坡頂和坡腳處，坡頂和坡腳處負孔隙水壓力分別增加了32.3，19.4 kPa，其安全系數(shù)降低是因為滑動體重量增加和坡頂及坡腳處非飽和土剪切強度降低；粉質(zhì)黏土中滑動面孔壓的變化主要分為兩個階段：第一階段坡頂和坡腳處的孔隙水壓力快速升高，降雨停止后進入第二階段，滑動面頂部和底部孔隙水壓力逐步降低，中部孔壓不斷升高，最終滑動面孔隙水壓力穩(wěn)定在-25～-15 kPa。伴隨著雨水的入滲，壤質(zhì)砂土邊坡臨界滑動面上孔隙水壓力快速升高，并在滑動面中部出現(xiàn)正孔隙水壓力，造成土體有效應(yīng)力降低，誘發(fā)邊坡失穩(wěn)，孔隙水壓力極值25.45 kPa 出現(xiàn)在降雨停止時，伴隨著水分的入滲滑動面上的孔隙水壓力不斷降低。

圖9 臨界滑動面孔隙水壓力

2.4 非飽和土邊坡可靠性及敏感性分析

本文選取特大暴雨工況下粉質(zhì)黏土邊坡進行可靠性分析和參數(shù)敏感性實驗。邊坡在不同土體參數(shù)條件下的安全系數(shù)概率分布見圖10，黏聚力、內(nèi)摩擦角和飽和重度參數(shù)改變誘發(fā)的邊坡失效概率分別為5.75%，0.60%，0.05%。

圖10 概率分布

為量化不同土體參數(shù)對邊坡穩(wěn)定性的影響，本文計算了不同土體參數(shù)下邊坡穩(wěn)定性安全系數(shù)，結(jié)果見圖11。黏聚力和內(nèi)摩擦角與邊坡穩(wěn)定性安全系數(shù)成正相關(guān)，其中黏聚力從0.5 kPa 增加到5.5 kPa 的過程中邊坡安全系數(shù)從0.878 增加到了1.217。土體飽和重度的增加對邊坡的安全性造成消極影響，飽和重度為20.5 kN/m3時，邊坡安全系數(shù)為0.98。

圖11 安全系數(shù)敏感度

3 結(jié)論

本文通過數(shù)值模型試驗對不同降雨和土質(zhì)條件下非飽和土邊坡穩(wěn)定性進行分析，得出如下結(jié)論：

1）在降雨條件下非飽和土質(zhì)邊坡穩(wěn)定性變化有顯著差異，高滲透率的壤質(zhì)砂土邊坡安全系數(shù)的極小值出現(xiàn)在降雨結(jié)束時刻，而低滲透率的黏性土和粉質(zhì)黏土邊坡的穩(wěn)定性安全系數(shù)在降雨之后仍持續(xù)降低。

2）黏性土邊坡滲透性較低，降雨過程中滑動面上孔隙水壓力的增長主要發(fā)生在滑動面頂部和底部，該部分基質(zhì)吸力的喪失降低了土體邊坡穩(wěn)定性。

3）高滲透性土體的邊坡失穩(wěn)模式是降雨導(dǎo)致的地下水位抬升誘發(fā)邊坡失穩(wěn)，其規(guī)模和影響深度大于低滲透性邊坡。

4）邊坡土體物理力學(xué)參數(shù)影響邊坡穩(wěn)定性，通過可靠性分析發(fā)現(xiàn)，土體的黏聚力和內(nèi)摩擦角與邊坡穩(wěn)定性安全系數(shù)成正相關(guān)，與土體飽和重度成負相關(guān)。