蔣水華，李典慶

（1. 南昌大學(xué) 建筑工程學(xué)院，江西 南昌 330031；2. 武漢大學(xué) 水資源與水電工程科學(xué)國家重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430072）

1 引 言

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

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

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

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

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

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

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

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

3 多重響應(yīng)面構(gòu)建

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

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

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

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

4 算例分析

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

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

圖3 2 816條隨機(jī)產(chǎn)生的潛在滑動面Fig.3 Slope with 2 816 randomly generated potential slip surfaces

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

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

圖5 驗(yàn)證多重響應(yīng)面的有效性Fig.5 Validation of multiple response surfaces

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

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

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

5 結(jié) 語

本文提出的基于多重響應(yīng)面的考慮參數(shù)空間變異性邊坡系統(tǒng)可靠度分析的蒙特卡洛模擬方法，能夠有效地分析考慮參數(shù)空間變異性低失效概率水平的多層土坡系統(tǒng)可靠度問題，并且具有較高的參數(shù)敏感性分析計(jì)算效率，為研究土體參數(shù)統(tǒng)計(jì)特征（均值、變異系數(shù)、邊緣概率分布、自相關(guān)距離和自相關(guān)函數(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)兵, 等. 考慮自相關(guān)函數(shù)影響的邊坡可靠度分析[J]. 巖土工程學(xué)報(bào), 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]. 應(yīng)用基礎(chǔ)與工程科學(xué)學(xué)報(bào),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.