孫元元，周生通，杜曉鵬

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

0 引言

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

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

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

1 NORTA-MCS方法原理

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

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

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

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

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

N——總樣本點(diǎn)個(gè)數(shù)；

(2)

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

m——隨機(jī)變量個(gè)數(shù)。

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

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

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

φ(·)——標(biāo)準(zhǔn)正態(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方法的難點(diǎn)也是其應(yīng)用的關(guān)鍵在于能否順利求解R0并對其進(jìn)行Cholesky分解。若用ρx,ij表示變量Xi和Xi在相關(guān)矩陣RX中對應(yīng)的元素，那么其與R0中元素ρ0,ij的映射關(guān)系為：

式中：μi,σi——分別是變量Xi的均值和標(biāo)準(zhǔn)差。

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

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

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

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

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

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

(7)

式中：A——滑動(dòng)面單位寬度面積/m2；

W——滑動(dòng)面上方巖體的重力/t；

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

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

上述各量的計(jì)算式分別為：

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)

在計(jì)算時(shí)，式中參數(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ù)實(shí)際情況取值；

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

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

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

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

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

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

g(X)=SX-1(12)

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

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

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

(13)

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

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

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

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

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

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

圖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)預(yù)防措施，本節(jié)探討錨固措施(T≠0時(shí))對巖坡失穩(wěn)概率的影響規(guī)律。取錨索拉力變化范圍為T=[0,300](t)，拉力角變化范圍為θ=[0,55](°)(即拉力方向按順時(shí)針方向由垂直于坡面變化為水平方向，見圖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 地震動(dòng)對失穩(wěn)概率的影響

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

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

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

式中：ψF——坡度；

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

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

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

4 結(jié)論

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

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

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