蔣水華，李典慶

（1. 南昌大學 建筑工程學院，江西 南昌 330031；2. 武漢大學 水資源與水電工程科學國家重點實驗室，湖北 武漢 430072）

1 引 言

天然土體受沉積、后沉積、化學風化、物理降解、水熱變化和搬運等作用以及不同荷載歷史的影響，在巖土工程實際中多層土坡十分常見，Ireland[1]通過對芝加哥國會街切坡失穩(wěn)事故調查，發(fā)現(xiàn)土體框架呈現(xiàn)明顯的分層現(xiàn)象。美國密蘇里州東北部的Clarence Cannon壩[2]和加拿大魁北克省詹姆斯壩地基土[3]都存在明顯的層狀分布特征。2010年2月15日意大利南部城鎮(zhèn)馬耶拉托城附近一場暴雨引發(fā)的山體滑坡，揭露出來的新鮮連續(xù)地層中灰黃褐色和淺色地層縱橫交錯，呈現(xiàn)明顯的非均質土體框架[4]。Suchomel等[5]通過試驗發(fā)現(xiàn)，捷克波希米亞南部Trebon盆地南邊砂性土也呈現(xiàn)明顯的層狀分布特征。上述案例均表明，不僅土體參數(shù)存在一定的空間變異性，而且天然土體框架呈現(xiàn)明顯的層狀分布特征。目前已對忽略土體參數(shù)空間變異性的多層土坡可靠度問題進行了大量有益的研究[6－9]，同時考慮土體參數(shù)空間變異性的單層土坡可靠度分析也取得了可喜的進展[10－15]。然而，對考慮土體參數(shù)空間變異性的多層（2層以上）土坡可靠度問題研究的卻遠遠不夠。

另一方面，邊坡穩(wěn)定性分析通常采用有限元方法或者極限平衡方法，所獲得邊坡安全系數(shù)一般是土體參數(shù)（如黏聚力和內摩擦角）的非線性隱式函數(shù)，比如采用直接蒙特卡洛模擬（MCS）方法計算邊坡可靠度，需要進行成千上萬次邊坡穩(wěn)定性分析獲取臨界安全系數(shù)，計算量非常龐大，尤其對于低失效概率水平復雜邊坡可靠度問題。為了提高計算效率，通常先采用多項式展開、人工神經網(wǎng)絡、支持向量機或克里金等代理模型建立邊坡安全系數(shù)與輸入?yún)?shù)間的近似顯式函數(shù)關系，再進行邊坡可靠度分析，這對于忽略土體參數(shù)空間變異性隨機變量數(shù)目較少且失效模式較為單一的邊坡可靠度問題非常有效。然而對于考慮土體參數(shù)空間變異性，參數(shù)隨機場需離散為大量的隨機變量，并存在多個空間變異潛在失效模式的多層土坡可靠度問題，這種方法的適用性較差[10]。因此，亟需發(fā)展高效的代理模型分析考慮參數(shù)空間變異性的多層土坡可靠度問題。

本文提出了基于多重響應面的考慮參數(shù)空間變異性邊坡系統(tǒng)可靠度分析的蒙特卡洛模擬方法，系統(tǒng)地研究了考慮土體參數(shù)空間變異性的多層土坡系統(tǒng)可靠度問題。

2 邊坡系統(tǒng)可靠度蒙特卡洛模擬

天然土體參數(shù)存在一定的空間變異性，同時呈現(xiàn)明顯的層狀分布特征，相應地邊坡存在多個潛在滑動面[7－8]，并且這些潛在滑動面分布也存在一定的空間變異性[12]。對于含多條潛在滑動面的邊坡穩(wěn)定性問題，通常只要有任何一條潛在滑動面的安全系數(shù)小于1.0，邊坡就會失穩(wěn)，因此邊坡穩(wěn)定性可定義為一個串聯(lián)系統(tǒng)可靠度問題[8，13]。為了有效地分析含多失效模式的多層土坡系統(tǒng)可靠度問題，本文提出基于多重響應面的蒙特卡洛模擬方法。根據(jù)串聯(lián)系統(tǒng)失效概率定義，可得邊坡系統(tǒng)失效概率pf,s：

式中：P(?)為某一事件或系統(tǒng)的失效概率；E[Sj]為邊坡沿第j條滑動面Sj失穩(wěn)的事件；Ns為潛在滑動面（失效模式）數(shù)目。式（1）可進一步表示為

式中：f(X)為輸入隨機向量X的聯(lián)合概率密度函數(shù)；g(X)為邊坡穩(wěn)定性分析的功能函數(shù)，g(X)=FS(X)- 1 .0。由于巖土工程實際中現(xiàn)場實測數(shù)據(jù)非常有限，一般難以獲得較完整的土體參數(shù)概率分布等信息，同時式（2）積分區(qū)域是一個復雜的函數(shù)，故對式（2）直接進行積分的計算難度較大。為提高計算效率，常采用MCS方法將式（2）簡化[8，13]為

式中：Nt為MCS抽樣次數(shù)，k=1， 2，··，Nt；FSmin表示對于給定的任意一組隨機樣本通過邊坡穩(wěn)定性分析所獲得的臨界安全系數(shù)；I[?]為邊坡失效區(qū)域的指示性函數(shù)[13，15]。

3 多重響應面構建

采用蒙特卡洛模擬進行邊坡系統(tǒng)可靠度分析最關鍵的一步是計算式（3）中的邊坡臨界安全系數(shù)FSmin，如果直接獲取FSmin，一般需要進行成千上萬次邊坡穩(wěn)定性分析計算量非常龐大。為提高FSmin的計算效率，本文也先采用代理模型建立邊坡安全系數(shù)與輸入?yún)?shù)間的近似顯式函數(shù)關系。由于考慮土體參數(shù)空間變異性的多層土坡存在多條潛在滑動面，故采用二階不含交叉項多項式展開建立每條潛在滑動面安全系數(shù)與原始空間隨機變量X間的響應面函數(shù)關系[8]，如對于第j條潛在滑動面，其安全系數(shù)響應面函數(shù)表達式為

式中：FSj(X)為邊坡第j條潛在滑動面安全系數(shù)，j= 1， 2，··，Ns；X為原始空間隨機向量；ai,j為待定系數(shù)，i=1， 2，··，Nc，Nc為待定系數(shù)數(shù)目，Nc=2N+1，N為隨機變量數(shù)目；Ψi,j(?)為多項式展開。需要說明的是，考慮參數(shù)空間變異性邊坡可靠度分析中原始空間隨機變量指的是空間任意點處的土體參數(shù)（如黏聚力和內摩擦角），為建立邊坡安全系數(shù)的多重響應面函數(shù)，重要一步是計算多項式展開系數(shù)。本文采用 2N+1組合樣本設計方法[8]，分別在每個原始空間隨機變量均值附近以其均值加減2倍標準差μi+2σi和μi-2σi產生樣本點，并將這 2N組樣本點和隨機變量均值一起作為邊坡穩(wěn)定性分析模型輸入?yún)?shù)，計算得到相應的邊坡安全系數(shù)，然后基于式（4）建立線性代數(shù)方程組，采用回歸分析方法計算二階多項式展開系數(shù)，進而得到第j條潛在滑動面安全系數(shù)響應面。類似地，對于Ns條潛在滑動面，只需對式（4）重復計算Ns次，便可獲得Ns條潛在滑動面安全系數(shù)的多（Ns）重響應面。需要指出，多重響應面的構建需要保證能夠獲得邊坡每條潛在滑動面安全系數(shù)，并且其計算效率隨著邊坡潛在滑動面數(shù)目的增加而降低。

由上可知，多重響應面構建過程中沒有涉及隨機場離散等信息，一旦多重響應面建立之后，便可根據(jù)土體參數(shù)隨機場統(tǒng)計特征（均值、變異系數(shù)、邊緣概率分布、互相關系數(shù)、自相關距離和自相關函數(shù)），采用喬列斯基分解技術[14]或 Karhunen-Loève（K-L）展開方法[15]等隨機場離散方法生成原始空間土體參數(shù)隨機場實現(xiàn)值。對于給定的第k組MCS樣本點，本文采用喬列斯基分解技術生成原始空間土體參數(shù)隨機場的第k次實現(xiàn)Xk(x,y)，(x,y)為二維計算區(qū)域內的任意點坐標，其詳細計算步驟見文獻[14]。最后，基于多重響應面通過式（5）計算第k組MCS樣本點處的邊坡FSmin(k)：

圖1 基于多重響應面蒙特卡洛模擬方法的計算流程圖Fig.1 Flow chart of multiple-response surface based Monte-Carlo simulation approach

4 算例分析

以文獻[9，16]中的3層不排水飽和黏土邊坡為例，驗證本文提出方法分析考慮土體參數(shù)空間變異性多層土坡系統(tǒng)可靠度問題的有效性。3層土坡計算模型如圖 2所示，坡高 6 m，坡度為 1︰3。3個黏土層的不排水抗剪強度均值分別為 18、20、25 kPa，3層土體的單位重度γsat均為18 kN/m3。采用簡化畢肖普法進行邊坡穩(wěn)定性分析得到臨界安全系數(shù)FSmin=1.285，自動搜索的臨界確定性滑動面如圖2虛線所示，同時經過黏土層1和2。與文獻[16]采用簡化畢肖普法計算的臨界安全系數(shù)1.282和臨界滑動面位置均非常吻合。

圖2 邊坡計算模型及穩(wěn)定性分析結果Fig.2 Analysis results of slope model and stability

圖3 2 816條隨機產生的潛在滑動面Fig.3 Slope with 2 816 randomly generated potential slip surfaces

為了驗證所構建的安全系數(shù)多重響應面能否有效地代替確定性邊坡穩(wěn)定性分析計算邊坡臨界安全系數(shù)FSmin，本文隨機產生了100組MCS樣本點，在此基礎上采用喬列斯基分解技術得到隨機場的100次實現(xiàn)。途徑1：將3參數(shù)隨機場的100次實現(xiàn)分別直接代入所構建的多重響應面中計算FSmin。途徑2：將3參數(shù)隨機場的100次實現(xiàn)分別賦給確定性邊坡穩(wěn)定性分析模型，采用簡化畢肖普法直接計算FSmin。其中隨機場的2次典型實現(xiàn)如圖4所示。圖中，顏色較深部分表示不排水抗剪強度參數(shù)值較大區(qū)域，顏色較淺部分表示不排水抗剪強度參數(shù)值較小區(qū)域。途徑2中對3參數(shù)隨機場的2次實現(xiàn)采用簡化畢肖普法進行邊坡穩(wěn)定性分析，得到的FSmin及對應的臨界滑動面分別如圖4（a）和4（b）所示。圖5比較了由途徑1多重響應面計算的邊坡FSmin與由途徑2簡化畢肖普法直接計算的FSmin。由圖5可知，以上兩種途徑計算的FSmin非常吻合，表明所構建的多重響應面可以有效地代替確定性邊坡穩(wěn)定性分析準確地計算FSmin。因此，本文提出方法一旦獲得土體參數(shù)隨機場實現(xiàn)值之后，只需通過簡單的數(shù)學表達式計算邊坡臨界安全系數(shù)，無需再進行確定性邊坡穩(wěn)定性分析，極大地提高了邊坡系統(tǒng)可靠度計算效率。

圖4 3參數(shù)隨機場的2次典型實現(xiàn)Fig.4 Two typical realizations of random fields

圖5 驗證多重響應面的有效性Fig.5 Validation of multiple response surfaces

最后，根據(jù)文獻[17]建議的方法模擬 3參數(shù)空間變異性，即將3層土體參數(shù)隨機場考慮為全局非平穩(wěn)隨機場，當參數(shù)自相關距離小于每層土體厚度，同一土層內參數(shù)隨機場遵循平穩(wěn)假設或者準平穩(wěn)假設，不同土層內任意兩點處土體參數(shù)之間的空間自相關性為0。隨機產生10萬組MCS樣本點，采用喬列斯基分解技術得到三參數(shù)隨機場的 10×104次實現(xiàn)，分別代入式（4）得到每條潛在滑動面的 10×104個安全系數(shù)，基于式（5）得到10×104個邊坡臨界安全系數(shù)，再采用式（3）計算得到邊坡系統(tǒng)失效概率pf,s。表 1為本文提出方法（MRSM+ MCS）的計算結果，包括安全系數(shù)統(tǒng)計矩（均值μFS、標準差σFS、偏度δFS和峰度κFS）和邊坡系統(tǒng)失效概率pf,s。為了進一步驗證本文提出方法進行邊坡系統(tǒng)可靠度分析的有效性，同時采用直接拉丁超立方抽樣（LHS）結合K-L展開方法進行邊坡可靠度分析，1 000和1×104次LHS方法的計算結果也列入表 1。由表中可知，本文提出方法與1×104次LHS方法的計算結果非常吻合，并且1 000次LHS方法也滿足計算精度要求，似乎沒有必要事先構建多重響應面，對于表1中工況1的邊坡失效概率水平(pf,s=0.136)較高的情況確實如此。然而，對于低失效概率水平邊坡系統(tǒng)可靠度問題(如pf,s<0.001)，這在巖土工程實際中更為常見，如當?shù)淖儺愊禂?shù)分別取為 0.3、0.1、0.1，其計算結果見表1中工況2。邊坡失效概率水平大大降低，此時1 000次LHS方法的計算精度顯然不夠，當計算次數(shù)增加至 3×104時，計算的pf,s為 7.33×10－4，與本文提出方法計算結果(pf,s=9.9×10－4)也較為吻合，說明本文提出方法此時同樣可以滿足計算要求，但是需要指出的是本文提出方法計算效率卻遠遠高于直接LHS方法，前者只需進行3 137次邊坡穩(wěn)定性分析，而后者則需要進行 3×104次邊坡穩(wěn)定性分析。證明本文方法能夠有效地分析考慮參數(shù)空間變異性的低失效概率水平多層土坡系統(tǒng)可靠度問題。

表1 3層土坡可靠度結果Table 1 Reliability analysis results for three-layered clay slope

對于表1工況2，采用本文提出方法重新計算邊坡系統(tǒng)失效概率時，無需重新進行邊坡穩(wěn)定性分析構建安全系數(shù)多重響應面，因為多重響應面構建過程中沒有涉及隨機場離散等信息，一旦獲得安全系數(shù)多重響應面之后，采用直接 MCS方法計算邊坡系統(tǒng)失效概率不再依賴輸入?yún)?shù)的統(tǒng)計特征（均值、變異系數(shù)、邊緣概率分布、自相關距離和自相關函數(shù)），可見本文提出方法還具有較高的參數(shù)敏感性分析計算效率，為研究土體參數(shù)統(tǒng)計特征對邊坡可靠度的影響提供了技術支持。

5 結 語

本文提出的基于多重響應面的考慮參數(shù)空間變異性邊坡系統(tǒng)可靠度分析的蒙特卡洛模擬方法，能夠有效地分析考慮參數(shù)空間變異性低失效概率水平的多層土坡系統(tǒng)可靠度問題，并且具有較高的參數(shù)敏感性分析計算效率，為研究土體參數(shù)統(tǒng)計特征（均值、變異系數(shù)、邊緣概率分布、自相關距離和自相關函數(shù)）對邊坡可靠度的影響提供了技術支持。

[1] IRELAND H O. Stability analysis of the Congress Street open cut in Chicago[J]. Géotechnique, 1954, 4(4): 163－168.

[2] HASSAN A M, WOLFF T F. Search algorithm for minimum reliability index of earth slopes[J]. Journal of Geotechnical and Geoenvironmental Engineering,1999, 125 (4): 301－308.

[3] EL-RAMLY H, MORGENSTERN N R, CRUDEN D M.Probabilistic slope stability analysis for practice[J].Canadian Geotechnical Journal, 2002, 39(3): 665－683.

[4] GATTINONI P, SCESI L, ARIENI L, et al. The February 2010 large landslide at Maierato, Vibo Valentia, Southern Italy[J]. Landslides, 2012, 9(2): 255－261.

[5] SUCHOMEL R, MA?IN D. Probabilistic analyses of a strip footing on horizontally stratified sandy deposit using advanced constitutive model[J]. Computers and Geotechnics, 2011, 38(3): 363－374.

[6] OKA Y, WU T H. System reliability of slope stability[J].Journal of Geotechnical Engineering, 1990, 116(8):1185－1189.

[7] CHOWDHURY R N, XU D W. Geotechnical system reliability of slopes[J]. Reliability Engineering and System Safety, 1995, 47(3): 141－151.

[8] ZHANG J, ZHANG L M, TANG W H. New methods for system reliability analysis of soil slopes[J]. Canadian Geotechnical Journal, 2011, 48(7): 1138－1148.

[9] ZHANG J, HUANG H W, JUANG C H, et al. Extension of Hassan and Wolff method for system reliability analysis of soil slopes[J]. Engineering Geology, 2013,160: 81－88.

[10] HUANG J S, GRIFFITHS D V, FENTON G A. System reliability of slopes by RFEM[J]. Soils and Foundations,2010, 50(3): 345－355.

[11] CHO S E. Probabilistic assessment of slope stability that considers the spatial variability of soil properties[J].Journal of Geotechnical and Geoenvironmental Engineering, 2010, 136(7): 975－984.

[12] WANG Y, CAO Z, AU S K. Practical reliability analysis of slope stability by advanced Monte Carlo simulations in a spreadsheet[J]. Canadian Geotechnical Journal, 2011,48(1): 162－172.

[13] LI L, WANG Y, CAO Z J, et al. Risk de-aggregation and system reliability analysis of slope stability using representative slip surfaces[J]. Computers and Geotechnics, 2013, 53: 95－105.

[14] 蔣水華, 李典慶, 周創(chuàng)兵, 等. 考慮自相關函數(shù)影響的邊坡可靠度分析[J]. 巖土工程學報, 2014, 36(3): 508－518.JIANG Shui-hua, LI Dian-qing, ZHOU Chuang-bing,et al. Slope reliability analysis considering effect of autocorrelation functions[J]. Chinese Journal of Geotechnical Engineering, 2014, 36(3): 508－518.

[15] 蔣水華, 李典慶, 曹子君, 等. 考慮參數(shù)空間變異性的邊坡系統(tǒng)可靠度分析[J]. 應用基礎與工程科學學報,2014, 22(5): 841－855.JIANG Shui-hua, LI Dian-qing, CAO Zi-jun, et al.System reliability analysis of slopes considering spatial variability of soil properties[J]. Chinese Journal of Basic Science and Engineering, 2014, 22(5): 841－855.

[16] KANG F, HAN S, SALGADO R, et al. System probabilistic stability analysis of soil slopes using Gaussian process regression with Latin hypercube sampling[J]. Computers and Geotechnics, 2015, 63: 13－25.

[17] LI D Q, JIANG S H, CAO Z J, et al. A multiple response-surface method for slope reliability analysis considering spatial variability of soil properties[J].Engineering Geology, 2015, 187: 60－72.