邊坡凸集模型非概率可靠度求解方法

高樂星，梁 斌，吳 政

(湖南工業大學 土木工程學院，湖南 株洲 412007)

0 引言

在實際的邊坡工程中，巖土體的基本參數具有很強的復雜性和變異性，試驗雖能獲取其物理力學參數，但可獲得的試驗樣本數量有限，不足以表明其分布類型和特征，而獲取大批量的樣本數據卻因為環境或成本等因素的限制變得困難。此時，難以再用例如矩法等常規的可靠度方法來進行概率可靠性分析。盡管試驗獲得的原始樣本數據較少，但它還是能夠確定參數的一個大致區間，即表現出“未知而有界”特征。為解決此類問題，非概率可靠性的概念[1]便應運而生，最早是在應用力學領域由Ben-Haim和Elishakoff提出。1994年，穩健可靠性理論[2]被提出，作為穩健可靠性的數學基礎，凸集合的概念被引入到非概率可靠性模型中，認為當已知的樣本信息匱乏時，采用集合模型來描述不確定性的波動范圍較合適。郭書祥等[3]則將用來表述不確定性參數的變量區間化，從而定義了一種新的非概率可靠性指標。王曉軍等[4]用集合表示結構的不確定性信息，并將結構安全域與區間變量域體積之比用來度量結構的非概率可靠性。喬心州等[5]在提出了一種基于橢球凸集的非概率可靠性模型的同時，還證明了概率模型與非概率模型的相容性。周凌等[6]則針對已有的超橢球凸集非概率可靠性指標的片面性，提出一種綜合指標并給出求解算法。

除了可靠度指標的求解途徑外，邊坡功能函數的建立也是邊坡結構穩定可靠度研究的重點。眾多不確定性因素導致邊坡的極限狀態方程一般都頗為復雜，呈隱含形式且非線性程度較高。因此，將邊坡功能函數顯式化后再進行相應的可靠度計算分析是本研究旨在解決的問題之一。運用代理模型不失為一種方法，響應面模型[7]、支持向量機模型[8]、神經網絡模型[9]等都是目前應用比較廣泛的代理模型。穆雪峰等[10]曾指出地質統計學中的Kriging模型對未知信息的預測能力要優于一些常用的代理模型。

本研究將按照如下思路展開研究：首先，由有限的樣本數據確定一個大致區間，在參數的區間描述下構建邊坡超橢球凸集模型；然后，將拉丁超立方試驗獲得的區間內樣本信息代入邊坡模型，得其響應值后利用Kriging代理模型實現邊坡功能函數的顯式化。最后，通過換算得到在標準正態空間內的非概率指標η。當η>1和0<η<1時，邊坡穩定性評價分別采用非概率指標和概率指標。本研究將Kriging代理模型與已有的超橢球凸集模型非概率可靠性理論結合起來，為求解樣本模糊特征下邊坡可靠度提供一種方法。

1 邊坡超橢球凸集模型的構建

目前在工程領域中，描述不確定性參數的凸集模型有不少，大體分為能量有界模型、包容有界模型、斜率有界模型和有界Fourier模型[11]這4類。有界Fourier模型因其在幾何上表示為多維超橢球也被稱為橢球凸集模型。描述不確定性參數的凸集模型中，區間凸集模型和橢球凸集模型這兩種凸集模型被應用得最為廣泛，且區間凸集模型適用于各個參數互不相關的情況，而當參數具有相關性時則適用橢球凸集模型。因此本研究將采用超橢球模型來界定邊坡參數的不確定量。

將結構的隨機參數和功能函數分別設為X=(X1,X2,…,Xn)T,Z=gX(X)。根據有限的樣本數據取Xi∈[Ximin,Ximax]，則隨機參數Xi的變差δi為：

(1)

(2)

故有:

|δi|≤ΔXi。

(3)

上式區間凸模型在幾何空間中表示一超長方體，該長方體的外接等效單位橢球用式子表示如下:

(4)

式中θei為橢球半軸。

并且此時超長方體的頂點在球殼上，則有:

(5)

由此看來，如何確定式(3)的最小外接橢球問題就轉換成了一個最值問題，即在已知式(5)的情況下求外接橢球的體積最小值。可利用拉格朗日乘數法求解，則半徑為1的超橢球凸集模型[12]表示為:

(6)

式中將W視為一各元素均大于0的對角矩陣，即意味著橢球各主軸分別平行于坐標軸，則有:

(7)

引進向量ui:

(8)

聯合式(8)和式(6)，原定超橢球凸集模型就轉化成:

Ec(u)={u|uTu≤1}。

(9)

而Ec為空間U中的一個單位超球集合，由式(8)可得:

(10)

具體求解時，要實現原始空間向標準向量空間的轉換[12]，只需將式(10)代入結構功能函數即可。

2 基于凸集模型的非概率綜合指標可靠度求解

2.1 邊坡功能函數

在邊坡穩定性分析中,邊坡的功能函數通常被定義為g(X)=g(X1,X2,…,Xn)。一般來說，當g(X)>0時,邊坡處于穩定狀態；當g(X)<0時,邊坡處于破壞狀態；當g(X)=0時,邊坡處于極限平衡的臨界狀態。具體計算時,通常利用滑動力S與抗滑力T之間的關系來表示邊坡的極限狀態方程[13],即

g(X)=g(X1,X2,…,Xn)=T-S。

(11)

此時傳統的邊坡穩定性分析方法可能需要進行一些變換來適應這種表示形式。為此，以確定性分析中的安全系數計算方法為基礎，一些專家學者重新定義了邊坡穩定極限狀態方程,即

g(X)=Fs-1，

(12)

式中安全系數Fs是關于隨機變量X的函數,X=(X1,X2,…,Xn)。

而確定性分析中安全系數的計算普遍采用極限平衡法，極限平衡法因其簡便適用至今已形成多種計算方法，通常有瑞典圓弧法、Bishop法、Janbu法、薩爾瑪法、摩根斯坦-普瑞斯法[14]等。本研究基于極限平衡法中的簡化Bishop法[15]來建立邊坡穩定功能函數。

簡化Bishop法安全系數Fs計算公式如下：

(13)

邊坡功能函數即為：

(14)

mαi=cosαi+sinαitanφi/Fs，

(15)

式中，Wi為計算參數；bi為第i個條塊重力、寬度、網弧底面傾角；ui，Ni，Ti為第i個條塊孔隙水壓力、法向作用力、圓弧底面剪力；ci為黏聚力；φi為滑面內摩擦角，見圖1。

圖1 滑體及條塊間作用力示意圖Fig.1 Schematic diagram of force between sliding body and bar

顯而易見，式(14)和(15)中都含有Fs和mαi，同時Fs又是關于ci，φi，Wi的函數，因此式(14)邊坡的功能函數Z是一個高度非線性的隱式函數。需尋求與要求偏導的傳統可靠度方法(如矩法)不同的方法以解決邊坡的可靠度分析問題。

2.2 Kriging代理模型

Kriging方法通過部分已知信息去模擬某點的未知信息，是一種更具統計性的半參數化插值方法，可以在某點一定范圍內求其線性最優無偏估計值。Kriging模型一般由回歸部分加隨機過程組成，其關系表示成如下形式[16]:

Y(X)=fT(X)ω+Z(X)，

(16)

式中，ω為回歸系數；f(X)為回歸模型，一般是關于變量X的多項式函數；Z(X)為一隨機過程。

(17)

式中，R為一對稱矩陣，由R(λ:S)構成，大小為m×m，對角線上的元素為1；λ為相關參數；r(Xnew)為待測點和訓練樣本間的相關向量；F是m維向量，由m個樣本點處的回歸模型組成；f(Xnew)為回歸多項式，視實際情況確定，一般可采用多項式的形式(2階或以下)。

在正態隨機過程的假設下，構造最優Kriging模型需要通過求解未知量λ，λ則通過求解式(19)所示優化問題獲得。于是通過最大似然估計得:

(18)

(19)

根據上述公式推導，將式(17)作為邊坡功能函數的近似顯示表達式，將邊坡穩定性的影響因素設為隨機變量形式X=(X1,X2,…,Xn)，則邊坡功能函數Z的近似表達式為：

(20)

首先將邊坡參數訓練樣本S=[X(1),X(2),…,X(m)]代入簡化Bishop模型中獲得樣本的真實響應值Y=[Y(1),Y(2),…,Y(m)]，然后由式(17)～式 (19)得到功能函數近似表達式(20)，求得各待定系數后取回歸多項式f(x)為0階，得邊坡功能函數的Kriging代理模型:

(21)

2.3 非概率綜合指標可靠度求解原理

文獻[6]所定義的一種多個超橢球凸集合描述時的可靠性綜合指標k，其表達式為：

(22)

式中，Rset為非概率可靠度；η為非概率可靠性指標。

圖2為二維空間中橢球凸集模型與極限狀態曲面的干涉關系，如圖所示，不確定變量所對應的凸域由單位圓表示。將標準向量空間的坐標原點到極限狀態曲面的最短距離定義為η。功能函數經變換由原始空間轉化至標準向量空間U中，極限狀態曲線將U劃分為兩部分，一個g(u)<0的區域，一個g(u)>0的區域，分別對應失效域和可靠域。由此可得當η>1時，結構總是穩定的，此時可用非概率可靠性指標η，即坐標原點到極限狀態曲面的最短距離來評估結構可靠度；而當0<η<1時，η指標便不再適用，因此對式(22)中的非概率可靠度Rset進行求解，以此作為可靠度評估指標，將功能函數變換到標準正態空間后再用蒙特卡洛法利用式(23)求解失效概率和可靠度指標[17]。

圖2 不同情況下非概率可靠性指標示意圖Fig.2 Schematic diagram of non-probabilistic reliability indicators under different conditions

(23)

式中，Z為邊坡功能函數；NMC為抽樣總次數，NZ≤0為功能函數Z≤0的樣本數量。

綜上，邊坡凸集模型非概率綜合指標可靠度求解方法的操作流程為:

(1)根據已有的樣本參數信息確定樣本區間范圍，設影響邊坡穩定性的隨機變量為X=(X1,X2,…,Xn)。

(2) 構建邊坡超橢球凸集模型,得到實現邊坡功能函數由原始空間轉換至標準正態空間的橋梁式(10)。

(3) 對樣本區間采用拉丁超立方試驗來構造樣本點，代入簡化Bishop模型式(14)獲得其響應值。

(4) 通過Kriging代理模型使邊坡功能函數顯式化。

(5) 引入基于超橢球凸集的非概率可靠性綜合指標η，當η>1時，通過迭代計算得到坐標原點到極限狀態曲面的最短距離，以此作為評估標準；當0<η<1時，將式(10)代入Kriging代理模型，運用蒙特卡洛法根據式(23)求解失效概率Pf。

3 實例分析

一雙層邊坡[18]，坡比為1∶2，高15.24 m。其土層參數信息見表1。

表1 參數樣本信息Tab.1 Parameter sample information

其中黏聚力c1，c2和內摩擦角φ2未知其具體數值與服從的分布類型，根據試驗經驗，把c1，c2，φ2作為隨機參數，表示為X=[X1,X2,X3]。

根據上述的可靠度求解步驟，所建立的半徑為1的超橢球凸集模型為：

(24)

表和ΔXi計算結果Tab.2 Calculated result of and ΔXi

求得加權矩陣W:

W=diag(6.31×10-4,1.62×10-3,0.025 7)。

(25)

對區間采用拉丁超立方試驗構造樣本點，代入邊坡功能函數式(14)得響應值后，由Kriging代理模型擬合獲得邊坡功能函數。

然后由式(10)，進行標準化變換的矩陣：

(26)

代入Kriging模型將功能函數轉化至標準正態空間。最后通過判斷，該邊坡有失穩的可能性，為0<η<1的情況，再運用蒙特卡洛法求得失效概率Pf和可靠度指標β，結果如表3所示。

表3 不同組數據可靠度計算結果Tab.3 Calculated reliability of different sets of data

由表可知，當取30組或40組樣本數據時求得的失效概率和可靠度指標相對誤差較小，且40組數據抽樣結果：Pf=1.47%，β=2.177 5，與文獻[18]中已知隨機變量的具體數值和已知分布類型的情況下求得的結果：Pf=1.46%，β=2.180 0非常接近。由此證明本研究提出的模糊樣本特征條件下的基于超橢球凸集模型和Kriging代理模型的非概率邊坡穩定可靠度求解方法是可行的，并且精確度很高，較于蒙特卡洛法其效率也明顯提高。

4 工程案例分析

湖南省某高速公路是國家規劃“七縱九橫”的第4縱南益高速公路的重要一段，全長86.724 km，從G56杭瑞高速程家山樞紐開始到S7101益陽繞城高速迎豐橋樞紐止。其中一段路塹邊坡高度約17.5 m，坡比為1∶3，層狀構造，其巖土層主要由(1)雜填土(2)粉質黏土(3)強風化泥質砂巖組成。

實地勘探得到其有限的物理力學參數樣本數據如表4所示。

表4 巖土層參數樣本數據Tab.4 Sample data of rock and soil layer parameters

表中黏聚力c1，c2，c3，內摩擦角φ1，φ2，φ3和重度γ1，γ2，γ3都只知其大致范圍，根據試驗經驗，把c1，c2，c3，φ1，φ2，φ3，γ1，γ2，γ3都作為隨機變量，表示為X=[X1,X2,…,X9]。

根據2.3節可靠度求解步驟,首先構建半徑為1的超橢球凸集模型：

(27)

表和ΔXi計算結果Tab.5 Calculated result of and ΔXi

(28)

然后通過區間拉丁超立方試驗得到40組樣本，代入邊坡功能函數式(14)得到響應值，如表6所示，將表6中數據代入Kriging代理模型計算得功能函數。

表6 工程實例樣本點及功能函數計算結果Tab.6 Project example sample points and function function calculation result

續表6

再由式(10)，進行標準化變換的矩陣：

(29)

代入Kriging模型將功能函數轉化至標準正態空間。

最后判斷該邊坡屬于0<η<1的情況，運用蒙特卡洛法求得失效概率和可靠度指標Pf=4.0×10-5，β=3.944 4，參考《公路工程結構可靠性設計統一標準》(JTG 2120—2020)中高速公路路基目標可靠度指標：3～2.5，結果表明該邊坡失穩概率非常低，邊坡比較穩定，同時對照表6中40組數據，由簡化Bishop法計算得到的40組安全系數平均值FS=2.033，遠大于邊坡安全系數臨界值1，進一步驗證了該邊坡的穩定性。

5 結論

研究表明，本研究主要取得以下幾個方面的成果：

(1)由于邊坡巖土體參數信息的不完備性，利用其有限的原始樣本數據構建了邊坡超橢球凸集模型，為邊坡功能函數的空間轉換提供了便利。

(2)利用LHS抽樣和具備超強學習能力的Kriging代理模型實現邊坡功能函數顯式化，然后將其變換至標準正態空間，求解失效概率和可靠度指標。

(3)求解過程引入基于超橢球凸集的非概率可靠性綜合指標η，當η>1時，通過坐標原點到極限狀態曲面的最短距離來衡量邊坡穩定性；當0<η<1時，運用蒙特卡洛法即可求得邊坡失效概率Pf與可靠度指標β。

(4)算例分析表明本研究方法可行且選取40組樣本數據時精確度較高；對照《公路工程結構可靠性設計統一標準》(JTG 2120—2020)，采用本研究方法分析工程案例得出的結論與運用簡化Bishop法的計算結果一致。算例分析與工程案例分析充分證明了本研究方法的可操作性、計算精度及工程實用性。