基于DEM-SRM和盲數(shù)理論的邊坡地震穩(wěn)定可靠性分析

劉煥玉，魏汝明，胡文奎，程愛華，陳訓(xùn)龍

(1.濟(jì)南市勘察測繪研究院，山東 濟(jì)南 250000； 2.華中科技大學(xué)巖土與地下工程研究所，湖北 武漢 430074)

1 引 言

我國是一個山脈眾多和地震頻發(fā)的國家，因地震而產(chǎn)生的滑坡災(zāi)害給國家和人民的財(cái)產(chǎn)和人身安全造成了巨大的損失[1]。因此，對邊坡工程的抗震分析和研究具有重要的意義。

現(xiàn)有的邊坡抗震分析主要是基于擬靜力法和動力時(shí)程分析法[2～4]。擬靜力法的基本思想[5]是在靜力計(jì)算的基礎(chǔ)上，將地震作用簡化為一個慣性力系附加在研究對象上，其核心是設(shè)計(jì)地震加速度的確定問題。由于擬靜力法的計(jì)算簡單，工程物理意義明確，且計(jì)算結(jié)果偏于安全，已成為邊坡工程中最常用的計(jì)算分析方法[6]。動力時(shí)程分析法是考慮地震的動態(tài)過程性，將地震力以地震加速度時(shí)程的方式施加于研究對象上，以此來分析研究對象的動力響應(yīng)及穩(wěn)定性，其核心問題是邊坡模型確定及地震輸入的準(zhǔn)確性[7，8]。隨著計(jì)算機(jī)技術(shù)及數(shù)值分析方法的發(fā)展，動力時(shí)程分析法可以實(shí)現(xiàn)復(fù)雜邊坡模型的建立和地震動力響應(yīng)的快速計(jì)算和結(jié)果輸出。因此，對于復(fù)雜邊坡模型的地震穩(wěn)定性分析，動力時(shí)程分析法是一個優(yōu)選方法[9]。然而，現(xiàn)有基于以上兩種分析方法的邊坡地震穩(wěn)定可靠性研究往往僅考慮了邊坡巖土體參數(shù)的一種或兩種不確定性，即隨機(jī)性和模糊性[10，11]。事實(shí)上，邊坡工程存有諸多不確定性因素，使得邊坡巖土體參數(shù)也具有多種不確定性，如隨機(jī)性、模糊性、灰色性和未確知性等。單獨(dú)僅考慮其中的一種或兩種不確定性易造成對邊坡穩(wěn)定可靠程度的高估，從而影響邊坡的抗震分析[12]。盲數(shù)理論可以同時(shí)考慮參數(shù)的多種不確定性，但目前采用盲數(shù)理論分析邊坡穩(wěn)定性問題的研究還很少[13，14]，而將盲數(shù)理論應(yīng)用于邊坡抗震穩(wěn)定可靠性分析的研究更是未見報(bào)道。

因此，本文基于盲數(shù)理論進(jìn)行邊坡的地震穩(wěn)定可靠性分析，并提出了一種分析方法。首先，基于離散元強(qiáng)度折減法和擬靜力法計(jì)算分析地震作用下邊坡的安全穩(wěn)定狀態(tài)，并結(jié)合響應(yīng)面法和蒙特卡洛法建立地震邊坡的極限狀態(tài)方程，其次，運(yùn)用盲數(shù)理論對邊坡巖土體參數(shù)進(jìn)行盲數(shù)化處理，并進(jìn)一步建立地震邊坡的安全系數(shù)盲數(shù)模型，最后，通過MATLAB編寫相應(yīng)的計(jì)算程序?qū)崿F(xiàn)安全系數(shù)盲數(shù)解的快速輸出和分析。

2 方法與基本原理

2.1 離散元強(qiáng)度折減法基本原理

離散元強(qiáng)度折減法的基本思想是在采用離散元法分析邊坡穩(wěn)定性時(shí)，采用強(qiáng)度折減法計(jì)算出邊坡的安全系數(shù)。

(1)離散元法基本原理

離散元法[15](Distinct Element Method，簡稱DEM)是20世紀(jì)70年代由Cundall提出的一種處理非連續(xù)介質(zhì)問題的數(shù)值模擬方法，其離散單元能更真實(shí)地表達(dá)巖土體的幾何特征，便于處理所有非線性變形等巖土體破壞問題，因此，被廣泛應(yīng)用于模擬邊坡滑坡等大變形的力學(xué)過程。該方法的理論基礎(chǔ)是結(jié)合不同本構(gòu)關(guān)系的牛頓第二定律，采用動態(tài)松弛法求解方程，其基本方程由塊體運(yùn)動方程和物理方程組成。

在數(shù)值計(jì)算分析中，將塊體離散化為有限差分三角形，網(wǎng)點(diǎn)由三角形的頂點(diǎn)組成，各網(wǎng)點(diǎn)的運(yùn)動方程如式(1)表示：

(1)

物理方程表示塊體間接觸點(diǎn)的力和位移的關(guān)系。各塊體之間的接觸可歸納為角-邊接觸和邊-邊接觸。力與位移的關(guān)系采用增量表示，具體如式(2)和式(3)所示：

△Fn=Kn△un

(2)

△Fs=Ks△us

(3)

其中：△Fn和△Fs分別為法向力增量和切向力增量；△un和△us分別為其對應(yīng)的法向位移增量和切向位移增量；Kn和Ks分別為接觸面上的法向剛度和切向剛度。

根據(jù)牛頓第二定律及中心差分法，可進(jìn)一步得到：

(4)

由此，可計(jì)算出每一時(shí)步的變形、振動及其對應(yīng)的變形速度、轉(zhuǎn)動速度等參數(shù)。

(2)強(qiáng)度折減法基本原理

強(qiáng)度折減法[16]又稱為抗剪強(qiáng)度折減系數(shù)法，其基本原理為：對邊坡巖土體的抗剪強(qiáng)度指標(biāo)，即黏聚力C和內(nèi)摩擦角φ，同時(shí)除以一個折減系數(shù)K，如式(5)和(6)所示的形式進(jìn)行折減，然后將這組折減后的值C′和φ′作為新的參數(shù)代入進(jìn)行計(jì)算。當(dāng)計(jì)算剛好不收斂時(shí)，此時(shí)坡體被認(rèn)為達(dá)到極限穩(wěn)定狀態(tài)，對應(yīng)的值稱為邊坡的安全系數(shù)。

C′=C/K

(5)

φ′=arctan(tanφ/K)

(6)

2.2 盲數(shù)理論

(1)盲數(shù)[12]的概念

(2)盲數(shù)的運(yùn)算

設(shè)*為運(yùn)算符號+，-，×，÷中的一種，A與B為2個盲數(shù)，其表達(dá)式分別為：

(7)

(8)

則A和B間的運(yùn)算法則A*B可通過可能值和可信度值兩個矩陣來實(shí)現(xiàn)。需要注意的是，在進(jìn)行盲數(shù)運(yùn)算時(shí)，A*B中的相同元素須進(jìn)行合并，當(dāng)成一個元素重新進(jìn)行排列，同時(shí)，在可信度矩陣中須將相同元素對應(yīng)的可信度值相加，所得之和即為這個元素的最終可信度值。

(3)盲數(shù)可信度

設(shè)A與B均為盲數(shù)，A*B在r處的可信度可表示為P(A*B≥r)，其表達(dá)式為：

(9)

其中：r為具體問題中要求值或標(biāo)準(zhǔn)值。

3 地震邊坡穩(wěn)定可靠性分析的盲數(shù)模型

3.1 地震邊坡極限狀態(tài)方程的建立

以一單層均質(zhì)土坡為例，考慮地震力作用，選用擬靜力法施加地震作用力，水平地震作用力Fe的表達(dá)式如式(10)所示：

Fe=Wαhξαi/g

(10)

其中：W為巖土體重力；αh為水平向設(shè)計(jì)地震加速度系數(shù)；ξ為地震作用效應(yīng)的折減系數(shù)；αi為巖土體質(zhì)心處水平地震加速度分布系數(shù)，可根據(jù)圖1內(nèi)插算得[17]，其中αm為分布系數(shù)最大值，當(dāng)?shù)卣鹪O(shè)計(jì)烈度為Ⅶ、Ⅷ、Ⅸ時(shí)，對應(yīng)αm分別取3.0、2.5和2.0；g為重力加速度，一般取 9.8 m/s2。

圖1 水平地震加速度分布系數(shù)參考取值

選取重度γ、黏聚力C和內(nèi)摩擦角φ作為基本變量，則根據(jù)響應(yīng)面法原理，建立地震作用下邊坡極限狀態(tài)的響應(yīng)面函數(shù)為：

g(γ，C，φ)=Fs-1=α+?1γ+?2C+?3φ+χ1γ2+χ2C2+χ3φ2

(11)

其中：Fs為地震邊坡安全系數(shù)；α，?1，?2，?3，χ1，χ2，χ3為響應(yīng)面的待定系數(shù)，需要7個方程確定求解，其求解步驟如下[12]：

(1)假定初始迭代點(diǎn)Y(1)={γ(1)，C(1)，φ(1)}，一般選取均值點(diǎn)。

(2)采用擬靜力法和離散元強(qiáng)度折減法，分別計(jì)算出地震邊坡的極限狀態(tài)函數(shù)值g{γ(1)，C(1)，φ(1)}、g{γ(1)±fσγ，C(1)，φ(1)}、g{γ(1)，C(1)±fσC，φ(1)}及g{γ(1)，C(1)，φ(1)±fσφ}，其中f為響應(yīng)面參數(shù)，一般選取固定值為2，σγ，σC，σφ分別為對應(yīng)參數(shù)γ，C，φ的標(biāo)準(zhǔn)差。

(3)將步驟(2)得到的迭代點(diǎn)和函數(shù)值代入式(11)中，計(jì)算得到7個待定系數(shù)，從而確定響應(yīng)面函數(shù)。

(4)采用蒙特卡洛法，計(jì)算步驟(3)中得到的極限狀態(tài)響應(yīng)面函數(shù)的可靠度指標(biāo)βk、可靠度Pk和驗(yàn)算點(diǎn)Y*(k)。

(5)計(jì)算驗(yàn)算點(diǎn)Y*(k)的極限狀態(tài)函數(shù)值g{γ*(k)，C*(k)，φ*(k)}，并利用下式(12)計(jì)算下一步迭代點(diǎn)Y(k+1)。

(12)

(6)根據(jù)步驟(5)得到的迭代點(diǎn)Y(k+1)，重復(fù)步驟(2)～(5)進(jìn)行下一步迭代計(jì)算，如滿足下式(13)的條件，可認(rèn)為計(jì)算收斂，則結(jié)束迭代計(jì)算，輸出此時(shí)的可靠度指標(biāo)βk以及可靠度Pk，此時(shí)得到的極限狀態(tài)響應(yīng)面函數(shù)即為所求。

|βk-βk-1|<0.001

(13)

3.2 地震邊坡安全系數(shù)盲數(shù)模型的建立

(14)

其中：γmin為γ的下限，γmax為γ的上限，f1(γ)為γ的盲數(shù)；Cmin為C的下限，Cmax為C的上限，f2(C)為C的盲數(shù)；φmin為φ的下限，φmax為φ的上限，f3(φ)為φ的盲數(shù)。

由此，根據(jù)式(11)，易建立并得到地震邊坡安全系數(shù)的盲數(shù)模型為：

Fs=α+1+?1{[γmin，γmax]，f1(γ)}+?2{[Cmin，Cmax]，f2(C)}+?3{[φmin，φmax]，f3(φ)}

(15)

3.3 地震邊坡安全系數(shù)盲數(shù)模型的快速計(jì)算

根據(jù)已知相關(guān)參數(shù)值及建立的盲數(shù)模型，通過MATLAB編寫相應(yīng)的計(jì)算程序，可以實(shí)現(xiàn)地震邊坡安全系數(shù)盲數(shù)解的快速輸出和分析，程序的具體計(jì)算步驟如下：

(1)根據(jù)式(15)，計(jì)算地震邊坡安全系數(shù)盲數(shù)解，輸出地震邊坡的總體可信度概率分布圖，同時(shí)，將計(jì)算結(jié)果以盲數(shù)標(biāo)準(zhǔn)形式輸出，如式(16)所示的三段式盲數(shù)解，實(shí)際可根據(jù)具體需要劃分區(qū)間段數(shù)：

(16)

其中：Fmin，F(xiàn)1，F(xiàn)2，F(xiàn)max為安全系數(shù)；P1，P2，P3為可信度。

(2)根據(jù)步驟(1)中盲數(shù)解以及前述響應(yīng)面法算得的地震邊坡可靠度Pk，進(jìn)一步計(jì)算出邊坡盲數(shù)安全系數(shù)Fps[12]。

(3)根據(jù)《建筑邊坡工程技術(shù)規(guī)范》(GB50330-2013)[18]中規(guī)定的邊坡安全最低安全系數(shù)值及式(16)，快速計(jì)算出邊坡最小安全要求可靠度Ps[12]。

4 算例分析

4.1 算例概況及模型建立

選用文獻(xiàn)[12]中的算例進(jìn)行驗(yàn)證和分析，已知該邊坡為一級邊坡，邊坡巖土體的物理力學(xué)參數(shù)如表1所示，其中γ，C，φ服從正態(tài)分布。該邊坡的場區(qū)抗震設(shè)防烈度為Ⅶ度，地震的相關(guān)參數(shù)如下：αh=0.10g；ξ=0.25；αm=3.0。

邊坡巖土體物理力學(xué)參數(shù) 表1

根據(jù)該邊坡的剖面尺寸，采用離散元法建立邊坡的三維模型如圖2所示，其中Z軸方向長度取為 1 m。邊界條件采用左右兩側(cè)水平約束，底部固定，上部自由。為簡化計(jì)算，固定Z軸約束和變形，將邊坡模型分析轉(zhuǎn)化為平面應(yīng)變問題，采用莫爾-庫倫屈服準(zhǔn)則的彈塑性條件進(jìn)行計(jì)算和分析。

圖2 邊坡的離散元三維數(shù)值模型

4.2 盲數(shù)計(jì)算與分析

根據(jù)已知條件，采用本文方法首先計(jì)算出地震作用下邊坡的極限狀態(tài)方程，同時(shí)，得到響應(yīng)面法迭代計(jì)算最終的可靠度Pk為58.32%。然后，進(jìn)一步求解出地震邊坡安全系數(shù)盲數(shù)解，具體結(jié)果如表2所示。

地震邊坡安全系數(shù)盲數(shù)解 表2

表2給出了地震邊坡安全系數(shù)取值區(qū)間及對應(yīng)的可信度大小，其中可信度可直接理解為可靠度。由此，可以計(jì)算出任一安全系數(shù)值所對應(yīng)的可靠度值，從而可以更全面地分析地震邊坡的安全穩(wěn)定狀態(tài)，并根據(jù)實(shí)際工程及抗震要求來進(jìn)行邊坡抗震設(shè)計(jì)與加固處理。根據(jù)表2的結(jié)果，可進(jìn)一步采用MATLAB軟件快速計(jì)算并輸出地震邊坡的總體可信度概率分布，結(jié)果如圖3所示。

圖3 地震邊坡總體可信度概率分布圖

從圖3中可以很直觀地看出各安全系數(shù)值所對應(yīng)的地震邊坡可信度的大小。當(dāng)安全系數(shù)取值為0.432時(shí)，對應(yīng)地震邊坡的可信度最大，此時(shí)可信度值為1，可認(rèn)為邊坡在地震作用下取此安全系數(shù)值最保險(xiǎn)；而當(dāng)安全系數(shù)取值為1.714時(shí)，對應(yīng)地震邊坡的可信度最小，可信度值為0，可認(rèn)為邊坡在地震作用下取此安全系數(shù)風(fēng)險(xiǎn)最大，風(fēng)險(xiǎn)概率為100%。

4.3 地震邊坡穩(wěn)定可靠性分析

根據(jù)表2和圖3的結(jié)果，可進(jìn)一步計(jì)算地震作用下邊坡的盲數(shù)安全系數(shù)Fps和最小安全要求可靠度值Ps。根據(jù)響應(yīng)面法迭代計(jì)算最終得到的可靠度Pk=58.32%，可計(jì)算出地震邊坡盲數(shù)安全系數(shù)Fps=1.073。已知該邊坡為一級邊坡，其對應(yīng)的邊坡穩(wěn)定安全系數(shù)要求值為1.30，由此，可計(jì)算出地震邊坡最小安全要求可靠度值Ps=33.08%。由Fps=1.073<1.30，可知該邊坡未達(dá)到一級邊坡穩(wěn)定安全系數(shù)的最低要求值，是偏不安全的，同時(shí)，由Ps=33.08%，可知該邊坡在地震作用下具有較高的失效風(fēng)險(xiǎn)。

未驗(yàn)證本文方法的可行性和有效性，與現(xiàn)有的一般方法進(jìn)行了對比，各方法的計(jì)算結(jié)果如表3所示，其中方法(1)為采用擬靜力法考慮地震作用，M-P法計(jì)算均值安全系數(shù)以及蒙特卡洛法計(jì)算可靠度；方法(2)為采用擬靜力法考慮地震作用，有限元強(qiáng)度折減法計(jì)算均值安全系數(shù)以及響應(yīng)面法計(jì)算可靠度；方法(3)為不考慮地震作用，基于盲數(shù)理論計(jì)算的盲數(shù)安全系數(shù)和最小安全要求可靠度[12]。

邊坡穩(wěn)定可靠性計(jì)算結(jié)果比較 表3

從表3中可以看出：①在考慮地震作用的情況下，本文方法是有效可行的，且計(jì)算結(jié)果較現(xiàn)有方法(1)和(2)偏低，由于本文方法考慮了邊坡巖土體參數(shù)的多種不確定性，分析因素更為全面，從而使得計(jì)算結(jié)果更為可靠且更具參考價(jià)值，對工程而言也更安全；②本文方法與方法(3)對比可知，地震作用對邊坡的整體穩(wěn)定可靠性具有明顯影響，這與現(xiàn)有已知結(jié)論一致，因此，在邊坡工程分析、加固處理及抗震設(shè)計(jì)時(shí)，應(yīng)考慮地震作用。

5 結(jié) 論

本文采用擬靜力法考慮地震作用，并采用離散元強(qiáng)度折減法、響應(yīng)面法和盲數(shù)理論提出了一種邊坡地震穩(wěn)定可靠性分析方法。結(jié)合具體算例分析，得出如下結(jié)論：

(1)本文方法考慮了地震作用和邊坡巖土體參數(shù)的多種不確定性，計(jì)算結(jié)果更為可靠且更具參考價(jià)值，對工程而言也更安全。

(2)根據(jù)本文方法計(jì)算出的地震邊坡安全系數(shù)盲數(shù)解和邊坡總體可信度概率分布，可以計(jì)算出任一安全系數(shù)值所對應(yīng)的可靠度值，可以更全面地分析地震邊坡的安全穩(wěn)定狀態(tài)，并根據(jù)實(shí)際工程等級情況及抗震要求來進(jìn)行邊坡抗震設(shè)計(jì)與加固方案設(shè)計(jì)。

(3)地震作用對邊坡的整體穩(wěn)定性具有明顯影響，在邊坡工程分析中應(yīng)考慮地震的因素，同時(shí)，本文方法也為邊坡地震可靠度分析提供了一個新的研究思路和參考。