許曉亮 張家富 曾林風 徐健文 史為政

(三峽庫區(qū)地質災害教育部重點實驗室(三峽大學)， 湖北 宜昌 443002)

邊坡作為常見的巖體之一，具有顯著的不確定性，其中一個主要的表現(xiàn)便是巖土參數(shù)的空間變異性[1]．為了更好地考慮巖土體參數(shù)空間變異特性，部分學者利用隨機場[2]來進行刻畫，也有學者通過機器學習[3]與貝葉斯壓縮感知[4]等方法利用有限現(xiàn)場實測數(shù)據對參數(shù)進下概率反演，從而更好地表征巖土參數(shù)變異性．巖土參數(shù)空間變異性與可靠性分析結合受到越來越學者的關注和探究[5-6]．

邊坡確定性分析方法主要包括極限平衡(LE)、極限分析法(LA)和有限元法(FEM)等．LE法需要提前假設一個破壞機制，需人為地假定一系列平衡條件，同時沒有考慮巖土體的應力應變關系，故理論嚴密性不足，而LA法將真實解限制在下限解[7]和上限解[8]之間，可信度要高于LE法，但是對于復雜的巖土工程問題，傳統(tǒng)LA法很難人為地構建靜止許可應力場和機動許可速度場，這導致了其應用受到較大限制．直到有限元和極限分析法相結合形成了有限元極限分析法[9](finite element limit analysis，F(xiàn)ELA)，該方法既克服了LA法在復雜工程下無法人為構建應力場與速度場問題，也避免了LE法理論嚴密性不足等問題．

如前所述，將隨機場理論與可靠性分析相結合，是考慮巖土體參數(shù)空間變異性影響的一個重要途徑．在邊坡工程的研究中，Cho[10]將整個邊坡區(qū)域生成隨機場，并進行極限平衡分析，但當參數(shù)空間變異性較大時，得到的滑移面不盡合理；Griffiths等[11]較早地將FEM與隨機場有效結合，研究了空間變異性對邊坡穩(wěn)定性的影響．同時，基于FELA的優(yōu)勢，將FELA與隨機場相結合的研究也自然受到了重視．Krabbenhoft等[12]將強度折減法應用于FELA，并結合隨機場理論得到了基于FELA的隨機場方法，實現(xiàn)了考慮參數(shù)變異性下的邊坡FELA強度折減穩(wěn)定性計算；Huang等[13]將有限元極限分析應用于滑坡問題，通過Monte-Carlo模擬建立了一種新的滑坡概率分析框架；陳朝暉等[14]基于FELA分析了參數(shù)空間變異性對邊坡穩(wěn)定性的影響，并與經典極限平衡法進行比較，發(fā)現(xiàn)通過傳統(tǒng)極限平衡法計算安全系數(shù)再進行可靠性分析，會高估邊坡失效概率．上述研究表明，除在邊坡穩(wěn)定性的確定性分析外，F(xiàn)ELA在邊坡可靠性分析中也逐漸得到了應用．然而，在FELA計算中，為了得到更為接近真解的上下解，需要增大單元數(shù)量，壓縮上下限區(qū)間至極小狀態(tài)，進而要開展大量的網格剖分，導致計算規(guī)模大幅增加，因此部分學者[15]通過網格自適應加密，通過對塑性區(qū)等應變率較大區(qū)域繼續(xù)加密處理，從而提高計算效率得到更小的上下限區(qū)間．目前，考慮網格自適應的FELA法在邊坡確定性分析中已經取得了一定的進展[16-17]，但在考慮空間變異性的隨機場條件下，采用網格自適應的FELA對邊坡穩(wěn)定性及可靠性分析結果的影響并不清楚．

為此，本文結合隨機場理論和FELA方法，利用Matlab編制了基于Karhunen-Loeve隨機場離散的二維隨機場程序以及Matlab與極限平衡、有限元軟件交互程序，提出了基于網格自適應有限元極限分析邊坡可靠性分析方法，針對一邊坡算例，將網格自適應前后的FELA法以及傳統(tǒng)極限平衡法進行對比分析，并基于網格自適FELA邊坡可靠性分析結果，探討了巖土參數(shù)變異系數(shù)和波動范圍對安全系數(shù)均值、可靠度指標以及邊坡潛在滑移面等因素的影響．

1 隨機場理論

1.1 土體空間變異性

為了更真實地表征土體參數(shù)的空間特征，需在對邊坡進行穩(wěn)定性分析時考慮參數(shù)空間變異性．根據隨機場理論，隨機場可通過參數(shù)均值、變異系數(shù)、波動范圍來表征空間變異性，為生成滿足要求的平穩(wěn)隨機場，采用被廣泛使用的指數(shù)型自相關函數(shù)來表示土體的空間變異性[18]，表達式如下：

(1)

式中：(x，y)分別為兩點空間坐標(x1，y1)與(x2，y2)；δh和δv分別為水平波動范圍與豎向波動范圍，波動范圍越大，表明土性參數(shù)分布越平穩(wěn)．

1.2 隨機場離散

隨機場模型一般不能在有限元分析中直接使用，需將隨機場離散從而與有限元相結合．隨機場離散的方法有協(xié)方差法、局部平均細分、Karhunen-Loeve(K-L)級數(shù)展開法等，本文采用K-L展開法[19]生成滿足要求的平穩(wěn)正態(tài)隨機場如下：

(2)

(3)

由于巖土體物理參數(shù)往往不能取負值[20]，需要將標準正態(tài)隨機場轉化成對數(shù)正態(tài)、Beta和截尾指數(shù)等分布的非高斯隨機場來表征參數(shù)的空間特征．考慮到巖土體參數(shù)大多服從正態(tài)和對數(shù)正態(tài)分布[20]，本文采用對數(shù)正態(tài)隨機場，將任一點的參數(shù)X作為隨機變量并服從對數(shù)正態(tài)分布，其K-L離散結果表達式如下：

(4)

式中：ulnXi和σlnXi分別為相應高斯參數(shù)隨機場中l(wèi)nXi均值和方差，由對數(shù)正態(tài)變換式(5)～(6)計算得：

(5)

(6)

上述二維隨機場進行離散時，對于可分離型自相關函數(shù)，二維的特征值、特征函數(shù)可簡化為一維方向上的特征值乘積、特征函數(shù)乘積[21]；

λi=λxi·λyi

fi(x，y)=fi(x)·fi(y)

(7)

一維方向上特征值與特征函數(shù)求解步驟如下，以x方向為例(y方向同理)：

1)設x區(qū)域隨機場范圍為Ω=[xmin，xmax]，用[-1，1]區(qū)間的標準Legendre多項式函數(shù)Lj-1(x)表示[xmin，xmax]區(qū)間的正交基函數(shù)：

(8)

其中：a和b分別為縮放參數(shù)和平移參數(shù)，a=(xmax-xmin)/2，b=(xmax+xmin)/2．

2)分別由式(7)生成x方向矩陣Ax(Ay同理)．

(9)

式中：ρ(x，y)為前述指數(shù)型自相關函數(shù)；通過積分換元以及參數(shù)的縮放平移可進一步改進為：

Lj-1(u)Lk-1(v)dudv

(10)

3)利用矩陣方程求解兩個方向上特征值λx(λy對應Dy)

AxDx=DxΛ

(11)

式中：Dx為矩陣Ax特征矩陣；Λ為對角矩陣，將其主對角元素從大到小排列即為λxi．

4)利用正交基函數(shù)hi(x)[22]求解特征函數(shù)fi(x)(fi(y)對應hk(y))

(12)

上述特征值以及特征函數(shù)的大小均取決于展開項數(shù)M，M越大離散越精確，為了控制離散精度和計算量，引入期望能比率因子ε[22]≥0.95，來控制離散誤差從而確定合適的展開項數(shù)M．

(13)

式中：S為二維離散區(qū)域面積；ε為期望能比率因子．

2 網格自適應有限元極限分析

2.1 有限元極限分析原理

有限元極限分析方法[12-13]基于塑性極限理論，結合廣義原理和混合有限元理論．有限元極限分析通過有限單元將巖土體內應力場和機動場離散化，從而在離散后的場內構建滿足條件的約束方程，并將巖土體內總的能量耗散(上限分析)或外力荷載(下限分析)作為上、下限目標函數(shù)，最終通過數(shù)學規(guī)劃的手段求解目標函數(shù)，從而得到上、下限解，基本思想與方法如下：

對于彈塑性材料，其廣義變分原理如下式：

s.t.F(α，κ)≤0

(14)

式中：α為荷載因子；σ為任一點應力狀態(tài)；b為重力矢量；u為位移矢量；t代表邊界應力．F(α，κ)為屈服函數(shù)，當κ=0即為理想彈塑性材料．

FELA下限解(LB)假定應力場內滿足平衡條件，應力邊界條件和屈服條件且理想彈塑性材料不破壞，此時的應力場即為靜力許可應力場，對應的外部荷載即為破壞荷載，下限解即滿足靜力許可條件的最大荷載，對應的目標函數(shù)優(yōu)化問題可如下表述：

maxα

PTσ=αtonsσ

f(σ)≤0

(15)

FELA上限解(UB)假定機動場滿足應變-位移關系、關聯(lián)流動法則和位移邊界條件，在彈塑性材料破壞時，內力功大于等于外力功．則對應的目標函數(shù)優(yōu)化問題可表述為：

(16)

2.2 網格自適應加密

在上述優(yōu)化問題求解過程中，為了縮小上下解區(qū)間，需增加網格單元數(shù)，若對網格均勻劃分的話，要達到滿足的精度需要極大的計算規(guī)模，但對塑性區(qū)外的網格進行劃分并不會提高精度，因此采用網格自適應加密，只在塑性區(qū)等應變率較大的區(qū)域進行加密，計算規(guī)模較小并能提高計算精度．

文中上限自適應采用一種基于能量耗散的自適應加密方法，即計算出每個單元的能量耗散，然后通過反復加密能量耗散較大的三角單元來降低破壞區(qū)域內單元間能量耗散的差異，從而在破壞區(qū)形成能量耗散密集的單元，而在其他區(qū)域單元呈疏松狀態(tài)，具體加密過程見文獻[16]．

下限自適采用一種基于節(jié)點應力的自適應加密策略[15]，將Mohr-Coulomb準則改寫成一個等式約束和一個二階錐約束，將前面優(yōu)化問題轉化為錐優(yōu)化的形式求解，最終通過衡量屈服準則殘余以及等效變形來指導網格進行剖分，從而使破壞區(qū)域網格更密集．

3 考慮網格自適應的邊坡可靠性隨機有限元極限分析方法

3.1 參數(shù)隨機場模擬

通過前述K-L離散方法生成平穩(wěn)隨機場，為驗證所生成隨機場的合理性，以黏聚力為例，生成不同均值度下不同波動范圍隨機場，如圖1所示．

圖1 不同波動范圍下的隨機場

由圖1中(a)、(b)可知，當δh/δv=1時，所模擬的黏聚力呈塊狀的分布，且隨著波動范圍的增大，塊狀分布更大，低黏聚力以及高黏聚力分布越集中；由圖中(c)、(d)可知，當δh/δv=10時，由于水平波動范圍較大，豎向波動范圍較小，不同波動范圍下的隨機場中，黏聚力分布呈水平層狀，且隨著波動范圍的增大，水平層狀帶更長，這一規(guī)律與已有研究一致[23]．

3.2 考慮網格自適應的邊坡可靠性FELA流程

本文通過Matlab生成N組對數(shù)正態(tài)參數(shù)隨機場，結合Monte-Carlo法，將每組參數(shù)隨機場作為隨機變量樣本，結合OPTUM-G2對邊坡進行N次確定性分析并結合強度折減[12]得到N組邊坡上、下限安全系數(shù)，代入式(17)所表征的功能函數(shù)g(X)，并由式(18)可計算得出邊坡失效概率Pf與可靠度指標β．

g(X)=Fs(X)-1

(17)

β=Φ-1(1-Pf)

(18)

實現(xiàn)流程如圖2所示，該過程由4個步驟組成：

Step1：建立邊坡計算模型，讀取單元中心坐標．首先在OPTUM-G2建立邊坡模型并進行網格劃分，導出對應單元信息．

Step2：生成參數(shù)隨機場．結合巖土參數(shù)統(tǒng)計信息，選取自相關函數(shù)與波動范圍，利用Matlab編制對數(shù)正態(tài)隨機場程序，基于Step1中單元節(jié)點信息生成N組巖土參數(shù)隨機場．

Step3：批量確定性分析．將上述隨機場參數(shù)批量導入原始邊坡模型中，生成N組邊坡計算文件，利用Matlab調用OPTUM-G2進行批量網格自適應FELA確定性分析并結合強度折減計算上、下限安全系數(shù)．

Step4：計算失效概率及可靠度指標．統(tǒng)計Step3中N組確定性分析計算結果，結合式(18)計算出失效概率與可靠度指標．

圖2 基于網格自適應FELA邊坡可靠性分析流程圖

4 算例分析

以一簡單土坡為算例，如圖3所示，坡高為8 m，邊坡頂部寬度為4 m，臨空坡面的坡比為1∶2．邊坡模型均勻劃分為1 000個三角單元，單元尺寸為0.46 m，單元尺寸與豎直波動范圍之比[24]為0.153<0.25，滿足要求，邊坡土體采用使用摩爾-庫侖準則和非關聯(lián)流動法則的理想彈塑性本構模型，計算模型底部采取水平與豎直約束，左側采取水平約束，右側無約束；模型加自重荷載．

圖3 邊坡模型圖

確定性分析時，不考慮參數(shù)空間變異特征，黏聚力c、內摩擦角φ、彈性模量E、重度γ、泊松比ν等參數(shù)采用均值見表1．

表1 土性參數(shù)取值表

在可靠性分析時，僅考慮c和φ的空間變異特征，其余參數(shù)視為定值，彈性模量、泊松比、剪脹角等參數(shù)對邊坡穩(wěn)定性狀態(tài)分析結果的影響可忽略不計[25]．算例中，c和φ均服從對數(shù)正態(tài)分布，變異系數(shù)、波動范圍等參數(shù)參考文獻[26-27]選取，初步設立隨機場參數(shù)取值如下，δh=30 m，δv=3 m．

4.1 確定性分析結果

當不考慮空間變異特征進行確定性分析時，計算參數(shù)采用均值，見表1．借助GEO-SLOPE采用Morgenstern-Price(M-P)法得出的安全系數(shù)為1.412．由FELA計算的塑性乘子云圖如圖4(a)、(c)所示，上限解(UB)安全系數(shù)為1.453，下限解(LB)為1.393；而采用網格自適應3次迭代后FELA計算云圖如圖4(b)、(d)所示，上限安全系數(shù)為1.419，下限為1.395．需要指出的是，隨著迭代次數(shù)的增大，計算所需時間逐漸增加，但經試算發(fā)現(xiàn)，當?shù)螖?shù)大于3以后，迭代次數(shù)的增加帶來的上下限區(qū)間縮小幅度逐漸減小，為提高計算效率，文中采取3次迭代的計算結果．

圖4 網格自適應前后FELA法確定性分析圖

通過塑性乘子云圖對比分析，自適應后的等效塑性區(qū)潛在滑移面網格分布較未自適應更加密集，網格加密區(qū)所反映的破壞形態(tài)更為清晰，滑動帶明顯更窄；且FELA法在網格自適應后計算的上限安全系數(shù)較小，下限也有所增大，上下限的區(qū)間范圍明顯更??；若要增加單元數(shù)量要達到相同效果，由表2可知，以上限解分析為例，需要增加單元至8 000個，計算時間也將由網格自適應所需的8.91 s增加至19.65 s．

表2 時間對比

可見，若不考慮網格自適應，在一次計算中需通過增加單元數(shù)來縮小上下限區(qū)間，其用時更長，而對于后續(xù)可靠性分析，成千上萬次的計算將大幅增加時長，而考慮網格自適應則能在保障上、下限解逼近真解的同時顯著提高計算效率．

4.2 隨機性分析結果

考慮土體參數(shù)的空間變異特征進行可靠性分析時，參數(shù)取值見表1，采用網格自適應的FELA方法，利用Monte Carlo法模擬2 000次．結果表明，土體參數(shù)變異性會導致邊坡滑移面發(fā)生顯著變化，同一安全系數(shù)往往會對應多個滑移面形式，對應于4.1節(jié)中自適應下限解為1.395的情況，在可靠性分析中共出現(xiàn)11個滑面(如圖5所示)，滑移面與圖4(b)中確定性分析潛在滑移面有較大區(qū)別，部分滑移面不再經過坡腳而經過坡面，滑移面積差異性也較大，其中最小滑移面積28.43 m2，最大滑移面面積64.87 m2，這表明參數(shù)空間變異性下，坡體潛在滑移面位置和規(guī)模也呈現(xiàn)出較為顯著的不確定性．

圖5 同一安全系數(shù)滑移面分布

為了進一步研究參數(shù)隨機場分布與潛在滑面的關系，以黏聚力c的隨機場為例，圖6(a)、圖6(b)分別給出了安全系數(shù)最小與最大時對應的潛在滑移面以及隨機場，圖6(c)、圖6(d)為滑移面積最小與最大時對應的潛在滑移面和隨機場．可看出潛在滑移面更靠近且穿過土體參數(shù)較小的區(qū)域，且基于M-P法得出的邊坡安全系數(shù)及潛在滑移面積均大于網格自適應FELA的結果．

圖6 典型樣本滑移面及隨機場

對比圖7中M-P法、FELA法和自適應后的FELA所得的安全系數(shù)概率密度曲線可知，采用FELA法，不論其上限或下限解，在考慮網格自適應和未考慮時其安全系數(shù)變異性均比較接近，且都明顯小于M-P法的計算結果．同時，F(xiàn)ELA法所得安全系數(shù)概率密度曲線均比M-P法“高瘦”，無論是否采用網格自適應，F(xiàn)ELA法所得失效概率均小于M-P法的失效概率，故由傳統(tǒng)極限平衡法進行可靠性分析往往會高估邊坡的失效概率．

圖7 概率密度曲線規(guī)律對比

分別對FELA法和自適應后的FELA生成的2 000次安全系數(shù)進行統(tǒng)計分析，圖8給出了相應的累計分布曲線和頻率分布直方圖．由圖8可知，自適應前后得到的邊坡安全系數(shù)均服從對數(shù)正態(tài)分布，自適應后，安全系數(shù)上限均值從1.408減小至1.378，下限均值從1.331增大至1.338，上下限均值區(qū)間明顯減??；同時，上限可靠度指標從1.903減小至1.805，下限可靠度指標從1.594增加至1.630，上下限更為趨近，由此計算出來的可靠度指標明顯更接近與真實狀態(tài)．且由圖8中的概率密度曲線可知，無論是否采用網格自適應，2 000次FELA法計算出來的安全系數(shù)中有約65%以上都低于確定性分析安全系數(shù)，這表明在該邊坡算例中，忽略巖土體參數(shù)空間變異性將會得出較高的安全系數(shù)，造成對邊坡穩(wěn)定性的不保守估計．

圖8 網格自適應前后FELA法安全系數(shù)統(tǒng)計規(guī)律

4.3 隨機場參數(shù)敏感性分析

變異系數(shù)與波動范圍是表征隨機場的重要參數(shù)，為研究其對計算結果的敏感性，仍采用網格自適應FELA方法，基于已有研究[26-27]設計計算參數(shù)，見表3．計算中，改變黏聚力變異系數(shù)ccov，內摩擦角變異系數(shù)φcov，水平波動范圍δh和豎向波動范圍δv，研究隨機場參數(shù)對邊坡的安全系數(shù)均值、可靠度指標、安全系數(shù)概率密度曲線以及滑移面分布范圍的影響規(guī)律．

表3 敏感性分析參數(shù)取值表

4.3.1 變異系數(shù)影響

為研究土體參數(shù)變異系數(shù)對邊坡的影響，取表3中組合1、2，分別研究c、φ變異系數(shù)的影響規(guī)律．以組合1為例，φcov=0.2，ccov取值范圍為[0.2,0.8]，波動范圍取值為δh=30 m和δv=3 m．

圖9和圖10分別為不同變異系數(shù)ccov、φcov下的可靠性指標β和安全系數(shù)均值uF值，結果表明β和uF有著相同的變化趨勢，均隨著變異系數(shù)的增大而減小．當ccov從0.2增加到0.8時，相應的uF上、下限解降幅分別達到8.1%和10%；β上、下限解降幅分別為68.6%和72.4%；而當φcov從0.05增加到0.25時，相應的uF上、下限解降幅為3.1%和3.9%；β上、下限解降幅則為45.9%和45.4%．對比數(shù)據，β比uF變化幅度明顯更大，說明相對于安全系數(shù)，可靠度指標更為敏感，這與類似研究的結論一致[27]．

圖9 β和uF隨ccov變化曲線

圖10 β和uF隨φcov變化曲線

以下限法為例，圖11和圖12分別給出了安全系數(shù)概率密度曲線隨ccov和φcov的變化趨勢．

圖11 不同ccov下安全系數(shù)概率密度曲線

圖12 不同φcov下安全系數(shù)概率密度曲線

不難看出，隨著參數(shù)變異系數(shù)增大，坡體物理力學性質的離散性變大，薄弱區(qū)更集中，更容易出現(xiàn)局部強弱化區(qū)，在概率密度曲線上整體呈現(xiàn)出向左偏移的同時逐漸扁平趨于“矮胖”狀的特點；相應的安全系數(shù)均值下降但變異性增大，坡體失效概率逐漸增大，可靠度指標減?。畬Ρ萩cov和φcov對概率密度曲線的影響規(guī)律可見，兩者雖然有相同趨勢，但后者相較于前者的影響更為顯著，即參數(shù)φ的變異性更為敏感．

仍以下限法為例，圖13和14為不同ccov和φcov下的坡體潛在滑移面分布范圍，圖中同一顏色的兩條滑面分別為某一變異系數(shù)下計算得到的滑面上下界限．當ccov和φcov增加時，淺層滑面位置上移，最小滑移面面積減小，深層滑面下移，最大滑移面面積增大，即坡體潛在滑面的分布范圍更大，滑移面分布更為離散．

圖13 不同ccov下滑移面分布圖

圖14 不同φcov下滑移面分布圖

4.3.2 波動范圍影響

為研究隨機場參數(shù)波動范圍對邊坡可靠性的影響，取表3中參數(shù)組合3～5進行分析．其中組合3，控制豎向波動范圍為3 m，水平波動范圍取值范圍[10 m，50 m]，主要研究水平向波動范圍的影響；組合4中，控制水平波動范圍為30 m，豎向波動范圍取值范圍[1 m，5 m]，主要研究豎向波動范圍影響；由于一般水平波動范圍比豎向波動范圍高一個數(shù)量級，為研究兩者共同影響，控制δh/δv=10，即表3中組合5．圖15和圖16反映了隨機場波動范圍δh、δv對β和uF的影響，與前述變異系數(shù)的影響有所不同，隨著波動范圍的增大，β呈現(xiàn)減小而uF呈現(xiàn)增大趨勢．當δh從10 m增加50 m時，相應的uF上、下限解增幅分別為0.9%和0.6%，β上下限降幅分別為11.9%和13.4%；當δv從0.05增加到0.25時，uF上、下限解增幅分別為3.0%和5.7%；而β上、下限解降幅分別為27.9%和28.5%．

圖15 β和uF隨δh變化曲線

圖16 β和uF隨δv變化曲線

從圖15～16可以看出，相對安全系數(shù)，可靠度指標對波動范圍更敏感，且豎向波動距離δv比水平向波動距離δh對安全系數(shù)和可靠性指標影響更大；此外，值得注意的是，隨著波動范圍增大，坡體安全系數(shù)均值變大的同時，其可靠性指標減小，這與一般地“越安全，越可靠”的認知相斥．

為了探究上述安全系數(shù)與可靠性指標變化趨勢互異的原因，圖17和圖18給出了不同δh和δv下的安全系數(shù)概率密度曲線．

圖17 不同δv下安全系數(shù)概率密度曲線

圖18 不同δh下安全系數(shù)概率密度曲線

可以看出，隨著波動范圍的增大，安全系數(shù)概率密度曲線整體呈現(xiàn)出向右偏移且趨于“矮胖”狀的特點，且隨豎向波動范圍δv的變化中更為明顯．概率密度曲線整體右移將使得安全系數(shù)均值有所增大，趨于“矮胖”狀意味著安全系數(shù)變異性增大，進而坡體失效概率增大，相應的可靠度指標減小．由于波動范圍表示隨機場中任意兩點力學參數(shù)不在相關的臨界距離，不難理解，波動范圍越大，表明土性參數(shù)分布越平穩(wěn)，其變異程度減弱，可進一步解釋安全系數(shù)均值有所增大的原因．

此外，結合圖15～16可以發(fā)現(xiàn)，無論δv還是δh，當增加到一定值時，安全系數(shù)均值與可靠度指標變化趨于平緩，受波動范圍影響漸漸減弱．這是由于波動范圍較小時，邊坡土體參數(shù)可能更容易出現(xiàn)局部弱強度區(qū)，當波動范圍增大，弱強度區(qū)范圍也逐漸增大，但是增加到一定值時，土體參數(shù)整體分布更趨于均值，對坡體可靠性指標等影響不再顯著．

圖19給出了δh和δv同步增大時2 000次蒙特卡羅模擬中下限解對應的潛在滑移面分布范圍．與坡體抗剪強度參數(shù)變異系數(shù)增大引起滑移面分布更為離散相反，當波動范圍整體增大時，土性參數(shù)分布越平穩(wěn)，其變異程度呈減弱趨勢，最小滑移面與最大滑移面之間區(qū)間范圍減小，滑移面分布范圍越來越集中．

圖19 不同波動范圍下滑移面分布圖

5 結 論

提出了考慮網格自適應的邊坡可靠性隨機有限元極限分析方法，通過均值邊坡算例，探討了網格自適應前后的FELA法的確定性和可靠性分析結果，同時在網格自適應FELA可靠性分析的基礎上，從安全系數(shù)均值、可靠度指標、安全系數(shù)分布規(guī)律、潛在滑移面分布范圍等角度探討了參數(shù)變異系數(shù)和波動范圍的影響，主要得出以下結論：

1)同等計算精度下，若不考慮網格自適應，需通過大幅增加單元數(shù)來縮小FELA的上下限區(qū)間，其用時更長，而考慮網格自適應的FELA上下限更為接近，潛在滑面位置更加明確，在更易于逼近真解的同時可顯著提高計算效率，尤其對于需大量抽樣計算的可靠性分析，其優(yōu)勢更突出．

2)對比基于M-P和網格自適應的FELA的隨機性分析結果，前者更易得出更為深層的潛在滑面和偏大的失效概率；如忽略巖土體參數(shù)空間變異性，采用確定性分析則會得出更大的安全系數(shù)，進而引起對邊坡穩(wěn)定性的過高估計．

3)參數(shù)變異系數(shù)對邊坡可靠性有顯著影響，隨著ccov和φcov的增大，安全系數(shù)概率密度曲線在整體向左偏移的同時逐漸略趨于扁平，安全系數(shù)均值下降但變異性增大，可靠度指標減小，同時，坡體潛在滑面的分布范圍更大，更為離散；相對于黏聚力c，可靠性分析結果及潛在滑面位置對內摩擦角φ的變異性更為敏感．

4)波動范圍越大，安全系數(shù)概率密度曲線在整體向右偏移的同時逐漸略趨于扁平，坡體安全系數(shù)均值與變異性均增大，進而可靠性指標減小，且豎向波動范圍的影響更顯著；但當波動范圍增加到一定程度時，波動范圍對可靠度指標的影響也漸漸減弱，主要是由于波動范圍越大土性參數(shù)分布越平穩(wěn)，其變異程度減弱，相應得出的潛在滑移面分布則更加集中．