基于NORTA-MCS方法的香港秀茂坪巖坡失穩(wěn)概率分析

孫元元，周生通，杜曉鵬

(華東交通大學機電與車輛工程學院，江西 南昌 330013)

0 引言

香港秀茂坪巖坡是巖土工程中一個典型的邊坡穩(wěn)定性案例，由于受到隨機因素的影響，其穩(wěn)定可靠性問題一直是國內(nèi)外學者關(guān)注的重點之一[1-4]。目前二階矩法和Monte Carlo仿真(MCS)方法是分析邊坡可靠性的常用方法。其中，前者屬于近似解法，雖然計算效率高，但往往會過高或過低地估計可靠度結(jié)果。相反，MCS方法雖然計算成本高、效率低，但當選取足夠樣本數(shù)時結(jié)果將趨近精確解，是當前工程中獲得近似真實失穩(wěn)概率的通用方法，而且該方法適用于任何規(guī)模、復雜度和非線性程度的可靠性問題。

不過，產(chǎn)生符合巖坡隨機因素概率信息的樣本隨機數(shù)是基于MCS方法的巖坡可靠性評估正確實施的關(guān)鍵環(huán)節(jié)之一。在影響巖質(zhì)邊坡平面穩(wěn)定性的因素中，結(jié)構(gòu)面的抗剪強度、結(jié)構(gòu)面傾角、張裂縫深度、張裂縫充水程度以及地震作用等都可能是不確定性的[3-4]。這些不確定性因素在實際中的概率分布類型往往是任意的，而且因素之間還常存在相關(guān)關(guān)系，但由于受實際條件限制這些隨機因素的完全概率信息多是無法獲得的，尤其是在隨機因素間完整依賴結(jié)構(gòu)的函數(shù)關(guān)系表達方面。目前，在眾多存在各種限制要求的相關(guān)隨機數(shù)生成方法[5-6]中，NORTA (NORmal To Anything)方法正是滿足這一情況的可生成具有任意邊緣分布和可行相關(guān)矩陣的多維隨機數(shù)生成算法。它最初由CARIO和NELSON[5]提出，而后GHOSH和HENDERSON[7]討論了其在高維情況下的性能。從本質(zhì)上講，NORTA方法是Nataf變換[8]在相關(guān)隨機數(shù)生成方面的應(yīng)用，其利用邊緣分布和高斯Copula假設(shè)[9]來構(gòu)造物理空間隨機因素的偽完備概率信息。

為此，本文將深入探討聯(lián)合NORTA方法和MCS方法解決相關(guān)隨機因素作用下的香港秀茂坪巖坡可靠性問題，給出基于NORTA-MCS方法的可靠度計算基本原理和實現(xiàn)方法，并討論典型工況下的香港秀茂坪巖坡穩(wěn)定可靠度及影響規(guī)律等問題。

1 NORTA-MCS方法原理

基于MCS方法的可靠度計算關(guān)鍵在于生成給定概率信息的足夠樣本，之后利用下式估算失穩(wěn)概率：

式中：E[·]——數(shù)學期望運算符；

Df={X|g(X)<0}——失穩(wěn)域；

g(X)——極限狀態(tài)函數(shù)；

IDf(x)——失穩(wěn)域的指示函數(shù)，當x∈Df時，IDf(x)=1，否則，IDf(x)=0；

N——總樣本點個數(shù)；

(2)

式中：Φ(·)——標準正態(tài)分布的累積分布函數(shù)；

m——隨機變量個數(shù)。

那么即通過式可將一個維數(shù)為m、相關(guān)系數(shù)矩陣為R0的多維標準正態(tài)向量Y=[Y1,Y2,…,Ym]T變換為實際物理隨機向量X，而基于此一組服從Y的高斯相關(guān)隨機數(shù)即可轉(zhuǎn)換為服從X的相關(guān)隨機數(shù)。不過需要注意的是，依據(jù)Copula理論可知上述變換式等同于在X中引入了高斯Copula假設(shè)[9]。由式(2)得X的概率密度函數(shù)為：

式中：φm(Y,R0)——m維標準正態(tài)聯(lián)合密度函數(shù)；

R0——向量Y的相關(guān)系數(shù)矩陣，可由物理向量X的相關(guān)系數(shù)矩陣RX得到；

φ(·)——標準正態(tài)分布的概率密度函數(shù)；

fi(·)——變量Xi的概率密度函數(shù)。

式(3)即為基于X的邊緣密度函數(shù)fi(xi),i=1,2,L,m及其相關(guān)系數(shù)矩陣RX，并結(jié)合高斯Copula假設(shè)得到的偽完備聯(lián)合概率密度函數(shù)，亦稱為Nataf分布[8]。

但NOTRA方法的難點也是其應(yīng)用的關(guān)鍵在于能否順利求解R0并對其進行Cholesky分解。若用ρx,ij表示變量Xi和Xi在相關(guān)矩陣RX中對應(yīng)的元素，那么其與R0中元素ρ0,ij的映射關(guān)系為：

式中：μi,σi——分別是變量Xi的均值和標準差。

φ(y)dy-μiμj)(5)

φ(y)dy-μiμj)(6)

綜上可知，基于式(4)、(5)和(6)即可計算得到R0，然后在檢驗并確認R0非負定的情況下，即可利用式(2)中的變換原理生成具有給定邊緣分布和相關(guān)系數(shù)的相關(guān)隨機數(shù)，最后利用這些隨機數(shù)和式(1)表征的MCS求失穩(wěn)概率的方法，即可實現(xiàn)結(jié)構(gòu)可靠度計算。本文將這一可靠性分析過程稱為NORTA-MCS方法。

2 秀茂坪巖坡的穩(wěn)定可靠性模型

秀茂坪巖坡的巖體是帶有席狀節(jié)理的未風化花崗巖，屬于平面滑動失穩(wěn)機理，可用圖1所示帶有張力裂縫的巖坡模型描述其穩(wěn)定性?；谠撃Ｐ涂赏茖С鲂忝簬r坡的安全系數(shù)S函數(shù)表達式為[4]：

圖1 帶張力裂縫的秀茂坪巖坡模型Fig.1 Sau Mau Ping rock slope model with tension crack

(7)

式中：A——滑動面單位寬度面積/m2；

W——滑動面上方巖體的重力/t；

U——滑動面上水壓產(chǎn)生的上舉力/t；

V——張力裂縫中水壓產(chǎn)生的水平力/t。

上述各量的計算式分別為：

A=(H-z)/sinψp(8)

W=0.5γH2[(1-(z/H)2)cotψp-cotψf](9)

U=0.5γwrzA(10)

V=0.5γwr2z2(11)

在計算時，式中參數(shù)的取值通常為：

H=60 m——總坡高/m；

ψf=50°——坡度/(°)；

ψp=35——滑面傾角/(°)；

γ=2.6t/m3——巖石比重/(t/m3)；

γw=1.0t/m3——水比重/(t/m3)；

c=10t/m3——黏聚力/(t/m3)；

φ=35°——內(nèi)摩擦角/(°)；

α=0.08——水平地震加速度系數(shù)；

T——錨索拉力/t，依據(jù)實際情況取值；

θ——拉力傾角/(°)，依據(jù)實際情況取值，在自然工況下(即未采取錨固措施)，T=0、θ=0°；

r=zw/z——為張力裂縫的充水程度；

zw——裂縫水深/m，與當?shù)赜晁鞖庥嘘P(guān)，計算時一般取為0.5。

從公式中可以看出，該巖坡穩(wěn)定模型除考慮了巖體重力、滑面摩擦阻力等內(nèi)部因素的影響外，還考慮了地震、裂縫水壓以及人工錨固等外部因素的作用。

在模型參數(shù)取為上述定值情況下，可得自然工況下秀茂坪巖坡的安全系數(shù)約為1.22。按照傳統(tǒng)觀點，該巖坡是穩(wěn)定的且有一定裕度。但現(xiàn)實是模型中的部分參數(shù)通常是無法準確測定或確定的，而更應(yīng)該以隨機量來表征這些參數(shù)，例如：材料不確定參數(shù)c和φ，幾何不確定參數(shù)z以及外界環(huán)境或載荷不確定參數(shù)zw、r和α等，研究還發(fā)現(xiàn)諸如黏聚力和摩擦角φ之間、充水程度r與張裂縫深度z之間還存在有明顯的負相關(guān)性[1,4]。按照不確定性傳播原理，受不確定性因素影響的秀茂坪巖坡穩(wěn)定安全系數(shù)也將具有不確定性，這就會使得確定條件下認定為安全的巖坡在實際不確定性因素影響下的安全性受到質(zhì)疑。

為此基于結(jié)構(gòu)可靠性理論，本文以c,φ,z,r,α五個參數(shù)為隨機參量，構(gòu)建秀茂坪巖坡穩(wěn)定可靠性模型，即在式基礎(chǔ)上定義極限狀態(tài)函數(shù)：

g(X)=SX-1(12)

式中：S(·)——安全系數(shù)函數(shù)，即式(7)；

X=[c,φ,z,r,α]T——物理空間不確定參量組成的隨機向量。

進而，巖坡失穩(wěn)概率表達為：

(13)

式中：f(X)——隨機因素的聯(lián)合概率密度函數(shù)。

式(7)、(12)和(13)共同構(gòu)成了香港秀茂坪巖坡穩(wěn)定可靠性模型。而利用上節(jié)給出的NORTA-MCS方法即可實現(xiàn)對該可靠性模型進行求解。

3 結(jié)果分析與討論

表1是秀茂坪巖坡穩(wěn)定可靠性模型中的隨機參數(shù)信息，有關(guān)取值參照文獻[4]，同時考慮實際情況將各個隨機參量取為截尾分布類型。在此基礎(chǔ)上，本節(jié)將基于NORTA-MCS方法分別對自然工況、相關(guān)性、錨固措施以及地震動影響下的秀茂坪巖坡穩(wěn)定性進行分析和對比，最后還討論了巖坡滑面傾角不確定時的巖坡穩(wěn)定性問題。

表1 秀茂坪邊坡穩(wěn)定可靠性模型的隨機參數(shù)

3.1 自然工況時的失穩(wěn)概率

對比可以發(fā)現(xiàn)，由于c和φ服從的截尾正態(tài)分布與正態(tài)分布差別不大，故其的可行域非常接近[-1, 1]區(qū)間，且相應(yīng)Y中的相關(guān)系數(shù)也與X中的相關(guān)系數(shù)僅存微小差別。而z和r由于兩者分布類型不同，且r服從的截尾指數(shù)分布與正態(tài)分布差別較大，故可行域縮小明顯且變換前后的相關(guān)系數(shù)差別也比較明顯。

圖2 相關(guān)系數(shù)ρ0,zr和ρx,zr的映射關(guān)系Fig.2 The mapping relationship between ρ0,zr and ρx,zr

3.2 相關(guān)性對失穩(wěn)概率的影響

秀茂坪巖坡中黏聚力c與內(nèi)摩擦角φ、張力裂縫深度z和充水程度r之間具有明顯相關(guān)性，如果忽略這些相關(guān)性或不恰當?shù)刂付ㄏ嚓P(guān)系數(shù)就很可能對巖坡失穩(wěn)概率估算結(jié)果造成較大影響。為此，本節(jié)分三種情形討論這種影響：

情形1：ρzr=-0.5，ρcφ=-0.5:0.01:0.5；

情形2：ρcφ=-0.5，ρzr=-0.5:0.01:0.5；

情形3：ρcφ=ρzr=-0.5:0.01:0.5。

其中，情形1考慮了c與φ相關(guān)性的變化，由于實際中它們多為負相關(guān)性，因而0～0.5的正相關(guān)關(guān)系是虛設(shè)的；情形2考慮了z和r的相關(guān)性變化，其中0～0.5的正相關(guān)系數(shù)可以反映Low指出的由流行地下水位控制的正相關(guān)工況；情形3則是一種虛設(shè)情形，當兩者同時為0時，代表不考慮相關(guān)關(guān)系的獨立工況。由圖3可以看到，在[-0.5, 0.5]的相關(guān)系數(shù)范圍內(nèi)，隨著相關(guān)系數(shù)代數(shù)值的變大，三種情況下的巖坡失穩(wěn)概率均逐步增大且近似呈線性變化。在情形3的獨立工況時失穩(wěn)概率(6.510 8%)約是自然工況(2.21%)下的3倍，可見忽略相關(guān)性對秀茂坪巖坡失穩(wěn)的評估影響非常大。從情形2中的z和r間正相關(guān)關(guān)系(0～0.5)可以看到，流行地下水位控制下的巖坡張力裂縫充水工況將使得巖坡失穩(wěn)概率增大，因而實際中非常有必要弄清巖坡張力裂縫充水的原因，正確判斷z和r間的相關(guān)性質(zhì)。另外，從變化程度來看，c和φ間相關(guān)性的變化對失穩(wěn)概率的影響程度要比z和r的大。例如：忽略c和φ間相關(guān)性(情形1：ρzr=-0.5，ρcφ=0)時，失穩(wěn)概率由原來自然狀態(tài)下的2.21%變化為5.075 6%；而忽略z和r間相關(guān)性(情形2：ρzr=0,ρcφ=-0.5)時，失穩(wěn)概率變化為3.633 6%。因而實際中要更注重對巖坡材料相關(guān)性的精確測定。需要注意到，圖4中失穩(wěn)概率變化曲線并非光滑曲線，這是MCS仿真結(jié)果的隨機性導致的。

圖3 不同情形下巖坡失穩(wěn)概率隨相關(guān)系數(shù)變化的規(guī)律Fig.3 Instability probability of rock slope with the changed correlation coefficients for different cases

3.3 錨固措施對失穩(wěn)概率影響

作為最常見的巖坡失穩(wěn)預防措施，本節(jié)探討錨固措施(T≠0時)對巖坡失穩(wěn)概率的影響規(guī)律。取錨索拉力變化范圍為T=[0,300](t)，拉力角變化范圍為θ=[0,55](°)(即拉力方向按順時針方向由垂直于坡面變化為水平方向，見圖1)，分析秀茂坪巖坡失穩(wěn)概率隨著錨索拉力T及拉力角θ的變化規(guī)律，結(jié)果如圖4所示。

圖4 巖坡失穩(wěn)概率隨錨索拉力T及拉力角θ的變化曲線Fig.4 Instable probability of rock slope with the changed anchor tension force T and the inclination θ

3.4 地震動對失穩(wěn)概率的影響

3.5 滑面傾角不確定時的最大失穩(wěn)概率

通常復雜巖坡的結(jié)構(gòu)面產(chǎn)狀較為離散且有限條件下常無法確定其確切分布，此時巖坡穩(wěn)定模型中的滑面傾角ψP便無法準確確定。這時一種常見的分析方法是：將ψP看作變量尋找?guī)r坡失穩(wěn)概率最大的滑面，即轉(zhuǎn)化為一個概率最優(yōu)化問題。本節(jié)采用文獻[3]中的方法，即：通過在一系列確定性的ψP值之下，得到一系列對應(yīng)的巖坡失穩(wěn)概率值；然后擬合得到最大失穩(wěn)概率所對應(yīng)的滑面傾角，以該滑面傾角計算巖坡最大失穩(wěn)概率值。一系列確定性的ψP值確定方法可在HOEK和BRAY給出的滑面臨界傾角附近取值，該臨界傾角的估算公式為[3]：

αcr=(ψF+φ)/2=(50°+35°)/2=42.5°

式中：ψF——坡度；

一是構(gòu)建項目成本核算體系。首先，依照作業(yè)成本法理論，建立以病人為中心的項目作業(yè)庫。融入醫(yī)院自身特色，進一步優(yōu)化資源動因和作業(yè)動因，規(guī)范作業(yè)流程。建立時間技術(shù)難度風險系數(shù)算法模型（時間×工作量×技術(shù)難度×風險系數(shù)）計算作業(yè)成本更為合理，也更能體現(xiàn)勞動技術(shù)價值，對鼓勵醫(yī)院績效考核具有指導意義，綜合現(xiàn)實情況，更具說服力。

φ——內(nèi)摩擦角。

圖5 不同滑動傾角時的巖坡失穩(wěn)概率Fig.5 Instable probability of rock slope for different failure surface angles

4 結(jié)論

(1)在高斯Copula假設(shè)下，NORTA方法可僅基于隨機參量的邊緣分布和線性相關(guān)系數(shù)產(chǎn)生所需的相關(guān)隨機數(shù)，在結(jié)合MCS方法的基礎(chǔ)上實現(xiàn)巖坡有限概率信息下的穩(wěn)定可靠性評估。

(2)忽略相關(guān)關(guān)系的秀茂坪巖坡失穩(wěn)概率要比自然工況下的高出約3倍；且整個巖坡的失穩(wěn)概率會隨著黏聚力c與內(nèi)摩擦角φ以及張力裂縫深度z與充水程度r間相關(guān)依賴性的增大而增大。因而，正確考慮和評估各變量間的相關(guān)關(guān)系和相關(guān)程度是準確評估秀茂坪巖坡失穩(wěn)概率的前提。

(3)相比自然工況下的秀茂坪巖坡失穩(wěn)概率，合理的人工錨固措施可有效地降低巖坡失穩(wěn)概率，而忽略地震動作用的巖坡穩(wěn)定性評估結(jié)果則會遠遠高估現(xiàn)實秀茂坪巖坡的安全性，從而埋下安全隱患。在另一方面，當滑動傾角難以確定時，應(yīng)該將其作為變量來考查秀茂坪巖坡的穩(wěn)定可靠性，求出最大的失穩(wěn)概率作為巖坡穩(wěn)定性預防的決策指標。