吳敏強(qiáng)，蔡宜洲，周 豪

(三峽大學(xué) 水利與環(huán)境學(xué)院，湖北 宜昌 443002)

1 概 述

極限平衡法于1927年由Fellenius首次提出，是邊坡穩(wěn)定性分析中最常用的方法[1]。當(dāng)今應(yīng)用較為廣泛的極限平衡法有Fellenius法[2]、簡化Bishop法[3]、Janbu法[4]、Sarma法[5]、Morgenstern-Price法[6]等。極限平衡法的研究領(lǐng)域分布廣泛，條間力的假設(shè)與滑動面的搜索是其主要應(yīng)用方向，而且在三維擴(kuò)展方面也成果顯著。

關(guān)于剛體極限平衡法的研究，國內(nèi)學(xué)界大致有如下幾種觀點(diǎn)。朱大勇等[7]通過假設(shè)一組合式為零的條間剪力分布，使整個滑體所有平衡條件得到滿足，驗證了簡化Bishop法計算公式中沒有條塊間剪力，是因某種組合式為零，而不是條間剪力實際為零的設(shè)想。鄧東平等[8]找到了一種用隨機(jī)角來搜索隨機(jī)滑動面、簡化Janbu計算安全系數(shù)的新方法，通過隨機(jī)模擬曲線，發(fā)現(xiàn)隨機(jī)角物理意義明確，優(yōu)勢獨(dú)特。鄭穎人等[9]針對Sarma法條間剪力方程存在的兩個問題：不滿足邊坡穩(wěn)定性分析的合理性要求與不能正確表示任意條塊的條間剪力，給出了新的條間剪力方程。

在整個發(fā)展過程中，學(xué)者們一直將多余未知量的假設(shè)選在條間作用力或作用位置上。對于嚴(yán)格的二維極限平衡法而言，這種處理方法尚能夠得到令人比較滿意的安全系數(shù)結(jié)果，但有時也存在收斂性問題。對于三維極限平衡法而言，由于假設(shè)較多，例如假設(shè)合理性、力學(xué)嚴(yán)格性，這種研究思路已經(jīng)漸顯不足。

Bell、朱大勇和鄭宏等學(xué)者將研究視線轉(zhuǎn)到滑面正應(yīng)力的分布模式上來。這種極限平衡分析方法具有力學(xué)嚴(yán)格性，與其他嚴(yán)格的條分法相比，求解過程更便捷、易于編程實現(xiàn)且不存在收斂性問題。

嚴(yán)格極限平衡法發(fā)展?jié)摿薮?，如采用此種方法對地震力作用下的邊坡穩(wěn)定性進(jìn)行研究，將會產(chǎn)生可靠有效的理論成果，也有望取得不錯的實際效益。如若研究成果的合理性與適用性被證明，并被廣泛接受，大量應(yīng)用，便可在工程邊坡的穩(wěn)定性分析中，獲得省時省力、節(jié)約成本、提升速率的良好效果。

2 基本原理

在滑面正應(yīng)力修正法計算分析中，若能準(zhǔn)確構(gòu)建滑面上的正應(yīng)力分布模式，就可以求解出準(zhǔn)確的安全系數(shù)。對于實際計算過程，滑面正應(yīng)力分布基本無法直接獲得。但若是假設(shè)一個合理的分布模式模擬真實狀態(tài)，計算精度達(dá)到要求也是允許的。對于滑面正應(yīng)力的分布模式，研究工作主要存在于兩個部分：①滑面正應(yīng)力假設(shè)模式；②求解安全系數(shù)。

2.1 滑面正應(yīng)力假設(shè)模式

Bell提出將滑面正應(yīng)力函數(shù)假設(shè)成含有兩個待定參數(shù)的方程，這樣就與安全系數(shù)這個未知數(shù)一起構(gòu)成3個未知數(shù)，3個未知數(shù)3個方程理論上一定能求出數(shù)值解，這樣安全系數(shù)在該假設(shè)模式下一定可以求出。假設(shè)形式見式(1)。

(1)

式中：λ1，λ2為假設(shè)的未知參數(shù)；σ(x)為滑面正應(yīng)力分布函數(shù)；σ0(x)為瑞典法的分布函數(shù)；a，b為滑體兩端點(diǎn)的橫坐標(biāo)。

鄭宏教授提議將滑體看作積分域，通過數(shù)學(xué)原理將其轉(zhuǎn)變?yōu)檫吔绶e分的研究方法不再局限于對滑體進(jìn)行條分，可直接獲得其安全系數(shù)。他建立了條塊的平衡條件后發(fā)現(xiàn)，滑面正應(yīng)力由兩部分構(gòu)成，其中滑體體積力和條間作用力各占一部分。通過進(jìn)一步的研究使滑面正應(yīng)力的具體力學(xué)組成結(jié)構(gòu)關(guān)系更加清晰明朗[10]。同時，在他的研究中提出了一種新的滑面正應(yīng)力分布的假設(shè)模式：

(2)

2.2 求解安全系數(shù)

圖1為具有一般形狀滑動面的邊坡，y=s(x)為任意形狀滑動面；y=g(x)為坡表地形線；單位寬度條塊上重力w(x)；單位寬度條塊上水平地震力為kcw(x)；表面水平方向上的分布荷載為qx(x)、豎直方向上的分布荷載為qy(x)；滑面上孔隙水壓力為u(x)；總正應(yīng)力為σ(x)；剪應(yīng)力為τ(x)；y=yt(x)為邊坡體推力線；滑面形心點(diǎn)坐標(biāo)為(xc,yc)。向坡體設(shè)置預(yù)應(yīng)力錨索并視為一個外部集中力對邊坡整體進(jìn)行極限平衡分析，其中P為每根錨桿(預(yù)應(yīng)力錨索)所提供的錨固力，n為錨桿(預(yù)應(yīng)力錨索)的數(shù)量，m為錨桿(預(yù)應(yīng)力錨索)編號，i為錨桿(預(yù)應(yīng)力錨索)在錨固過程中與水平線的夾角，(xpm,ypm)為錨桿形心點(diǎn)坐標(biāo)。

圖1 邊坡滑體及受力分布

設(shè)邊坡安全系數(shù)為Fs，滑面摩擦系數(shù)為tanψ(x)，有效黏聚力為c(x)，根據(jù)摩爾-庫侖準(zhǔn)則，可得：

(3)

安全系數(shù)的平衡方程如下：

(4)

(5)

(6)

式中：s′=ds/dx=tanα;α為計算點(diǎn)滑面的傾角。若能確定滑面正應(yīng)力的分布函數(shù)，就可得到安全系數(shù)。

對邊坡滑體的任一條塊進(jìn)行受力分析，見圖2。

圖2 邊坡滑體受力模型

圖2中，E為條間水平向作用力；T為條間豎向作用力；h為條塊高度；ht為條間力的作用點(diǎn)高度；Δx為條塊寬度。

(7)

(8)

根據(jù)摩爾-庫侖強(qiáng)度準(zhǔn)則：

(9)

由式(9)可見，滑面正應(yīng)力由滑體體積力和坡面外力對滑面正應(yīng)力的貢獻(xiàn)(記為σ0(x))和土條間作用力對滑面正應(yīng)力的貢獻(xiàn)(記為σ1(x))兩個部分構(gòu)成，則可簡記為：

σ(x)=σ0(x)+σ1(x)

(10)

(11)

(12)

在滑面已知前提下，σ0(x)屬于已知函數(shù)。為了方程組能夠求解，可以采用以下形式來對滑面正應(yīng)力進(jìn)行假設(shè)：

σ(x)=σ0(x)+f(x)

(13)

式中：f(x)為逼近函數(shù)，f(x)構(gòu)造為含有兩個待定參數(shù)的函數(shù)才能保證平衡方程能夠順利求解。逼近函數(shù)對滑面正應(yīng)力的貢獻(xiàn)是非常有限的，不需過分追求與真實分布之間的吻合度，可以選用線性函數(shù)逼近：

(14)

將式(11)、式(13)、式(14)和摩爾-庫侖準(zhǔn)則帶入平衡方程，經(jīng)整理可得：

(15)

(16)

(17)

解平衡方程，可以得到：

(18)

其中，

(19)

(20)

(21)

其中，

T0=A3B2-A2B3

(22)

S0=A1B3-A3B1

(23)

G0=A1B2-A2B1

(24)

則式(18)可以寫成：

(25)

(26)

解式(25)，取實根得：

(27)

方程具有解析意義，是邊坡嚴(yán)格極限平衡法顯式解答，適用一般形狀滑裂面。

3 地震力作用下邊坡穩(wěn)定性計算

3.1 均質(zhì)邊坡穩(wěn)定性計算

取某均質(zhì)土質(zhì)邊坡，邊坡長100 m，邊坡高20 m，邊坡坡角為45°，邊坡的斷面尺寸見圖3。該坡體重度為19 kN/m3，黏聚力為30 kPa，內(nèi)摩擦角為35°。

圖3 土質(zhì)邊坡斷面

均質(zhì)邊坡的穩(wěn)定性分別采用Slide軟件和Geostudio軟件來計算。由于邊坡的坡比為1∶1，可知該邊坡最大水平動力系數(shù)為2，最大豎向動力系數(shù)為1，地震烈度為Ⅶ度，參照《水利水電工程邊坡設(shè)計規(guī)范》(SL 386-2007)取地震作用綜合系數(shù)為 0.25，取地震作用重要性系數(shù)為 1.0。Slide軟件Bishop法計算得到安全系數(shù)為1.571，Janbu法計算得到安全系數(shù)為1.488，見圖4。

3.1.1 數(shù)值應(yīng)力場獲取

通過對均質(zhì)邊坡的二維建模以及材料參數(shù)的賦值，進(jìn)行分析計算，以便從計算軟件中提取出滑動面上各節(jié)點(diǎn)的有用應(yīng)力信息。土坡的穩(wěn)定性計算分析采用理正巖土計算軟件，采用軟件中內(nèi)置的地震烈度考慮地震效應(yīng)。

圖4 均質(zhì)邊坡穩(wěn)定性計算

3.1.2 有限元法與滑面正應(yīng)力修正法結(jié)合

初始正應(yīng)力分布函數(shù)σ0(x)通過軟件得取可按照以下步驟操作來進(jìn)行：

1) 條分整個滑體，讀取滑面穿過的條塊相對應(yīng)位置的坐標(biāo)。

2) 讀取相應(yīng)坐標(biāo)點(diǎn)的應(yīng)力信息(σx,σy,τxy)，參照式(28)可以換算成對應(yīng)點(diǎn)的初始正應(yīng)力信息：

(28)

式中：條塊的底面中點(diǎn)與數(shù)值方向的夾角為β。單元的應(yīng)力信息可直接由計算軟件導(dǎo)出,當(dāng)正應(yīng)力初始分布形式通過σ0(x)數(shù)值方法得到后，再結(jié)合線性修正方法，對任意邊坡安全系數(shù)的求解會更加方便、準(zhǔn)確。

在Slide軟件中，根據(jù)斷面圖所提供的信息建立計算模型，由提供的已知條件賦予相應(yīng)的荷載，地震烈度采用7級，可由計算軟件根據(jù)極限平衡理論自動搜索安全系數(shù)最小的滑面。在ANDIA軟件中提取滑面上的應(yīng)力信息，經(jīng)過公式轉(zhuǎn)換擬合出對應(yīng)的初始滑面正應(yīng)力分布函數(shù)，見圖5。

圖5 初始滑面正應(yīng)力函數(shù)分布曲線

將得到的初始滑面正應(yīng)力分布函數(shù)進(jìn)行修正，計算得出邊坡在地震力作用下嚴(yán)格極限平衡法的安全系數(shù)。

經(jīng)過曲線擬合，邊坡圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.0165x3-1.0566x2+17.38x+5.1878

計算在地震工況下邊坡的穩(wěn)定性，將得到的滑面正應(yīng)力的初始正應(yīng)力分布函數(shù)σ0(x)，根據(jù)式(29)對滑面正應(yīng)力進(jìn)行修正，其中λ和ξ為假定未知參數(shù)。

σ(x)=σ0(x)+λζ1(x)+ξζ2(x)

(29)

(30)

(31)

修正后的滑面正應(yīng)力分布函數(shù)見圖6。

圖6 修正后滑面正應(yīng)力函數(shù)分布曲線

經(jīng)過曲線擬合，修正后邊坡圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.0177x3-1.1266x2+15.366x+2.405

根據(jù)式(29)對滑面正應(yīng)力進(jìn)行修正，根據(jù)式(3)可以得出安全系數(shù)的顯示解。選取Geostudio軟件與嚴(yán)格極限平衡法結(jié)果進(jìn)行比較，具體結(jié)果見表1。

表1 均質(zhì)邊坡安全系數(shù)匯總

根據(jù)計算結(jié)果可以看出，基于極限平衡理論的滑面正應(yīng)力修正法與數(shù)值方法相結(jié)合的穩(wěn)定性計算，可以直接從數(shù)值模擬軟件中了解滑面底部的正應(yīng)力分布情況，通過簡單的數(shù)值計算即可獲得比較準(zhǔn)確的初始正應(yīng)力，再對初始正應(yīng)力進(jìn)行修正，所得滑面正應(yīng)力函數(shù)比較接近實際情況，所得安全系數(shù)的值與傳統(tǒng)極限平衡法所得結(jié)果的相對誤差非常小，從而驗證了嚴(yán)格極限平衡法的合理性，便于實際應(yīng)用。

3.2 非均質(zhì)邊坡穩(wěn)定性計算

選取具有一般形狀典型的非均質(zhì)邊坡，其剖面形狀見圖7。邊坡坡高15 m，邊坡比為1∶2，共有4個土層，各土層物理力學(xué)參數(shù)見表2。取地震峰值加速度為1.63 m/s2，對照中國地震烈度表取對應(yīng)的地震烈度為Ⅶ度；參照《水利水電工程邊坡設(shè)計規(guī)范》(SL 386-2007)取地震作用綜合系數(shù)為 0.25，取地震作用重要性系數(shù)為1.0。

圖7 非均質(zhì)邊坡剖面圖

表2 土層物理力學(xué)參數(shù)

基于數(shù)值應(yīng)力場的非均質(zhì)邊坡穩(wěn)定性計算，可參照以下步驟進(jìn)行：

1) 在Geostudio軟件中，根據(jù)剖面圖所提供的信息建立計算模型。

2) 可由計算軟件根據(jù)極限平衡理論自動搜索出安全系數(shù)最小的圓弧滑面。

3) 在大型有限元計算軟件ANSYS中，建立該剖面的二維有限元計算模型，其二維模型見圖8，采用D-P彈塑性本構(gòu)模型。

4) 在有限元軟件中提取兩種滑動面上的應(yīng)力信息，經(jīng)過公式的轉(zhuǎn)換擬合出對應(yīng)的初始滑面正應(yīng)力分布函數(shù)。

5) 將得到的初始滑面正應(yīng)力分布函數(shù)進(jìn)行修正，計算得出安全系數(shù)。

6) 對比嚴(yán)格極限平衡法與傳統(tǒng)極限平衡法所得結(jié)果。

圖8 非均質(zhì)邊坡二維有限元模型

在進(jìn)行材料參數(shù)賦值后對該邊坡進(jìn)行穩(wěn)定性分析，對應(yīng)極限平衡法搜索出的滑面位置在相應(yīng)位置導(dǎo)出其應(yīng)力信息，通過式(28)的轉(zhuǎn)換，可得到選取單元上的初始正應(yīng)力值。經(jīng)過正應(yīng)力修正后便可以擬合出修正后滑面正應(yīng)力分布函數(shù)，見圖9。

經(jīng)過曲線擬合，修正后邊坡圓弧滑面正應(yīng)力分布函數(shù)為：

y=-0.0071x3+0.2277x2+8.5385x+41.095

按照均質(zhì)邊坡的計算過程，可以得出非均質(zhì)邊坡在地震作用下嚴(yán)格極限平衡法的安全系數(shù)，其結(jié)果匯總見表3。

圖9 修正后滑面正應(yīng)力函數(shù)分布曲線

表3 非均質(zhì)邊坡安全系數(shù)匯總

根據(jù)計算結(jié)果可以看出，在地震作用下嚴(yán)格極限平衡法計算得出的結(jié)果小于Bishop法、M-P法和瑞典圓弧法的安全系數(shù)，但是大于Janbu法和Bishop法折線滑面的計算結(jié)果。與Bishop法相比，結(jié)果相對誤差為1.9%；與Bishop法折線滑面相比，結(jié)果相對誤差為2%；但是與Janbu法相比，相對誤差為11.1%?？傮w來說，誤差范圍在1.9%～11.1%之間。造成誤差的原因主要在于正應(yīng)力提取的不完整和正應(yīng)力修正造成一定的誤差，但在容許的范圍內(nèi)，說明嚴(yán)格極限平衡法在非均質(zhì)邊坡應(yīng)用的合理性與適用性。

4 工程應(yīng)用

4.1 工程概況

牙根水電站在雅礱江干流上，位于四川省甘孜州雅江縣境內(nèi)，是我國重要水電能源基地雅礱江流域的大型電站工程?？碧浇衣秹沃穮^(qū)右岸Ⅰ-Ⅴ勘探線之間發(fā)育拉裂松動巖體，其地表后緣邊界高程自Ⅰ-Ⅴ線從2 640 m→2 667 m→2 690 m→2 697 m→2 597 m呈弧形變化，上游從PD17硐下部基巖突出處延伸至江邊，下游側(cè)從PD15硐延伸至江邊。此外，在Ⅶ線高高程PD31平硐(硐口高程2 612.67 m)硐深0～57 m，Ⅸ線低高程PD33平硐(硐口高程2 513.74 m)硐深0～129.5 m，中高程PD35平硐(硐口高程2 564.66m)等處亦發(fā)現(xiàn)有拉裂松動巖體分布。通過勘探，壩址區(qū)右岸存在的拉裂松動巖體對大壩的整體布局、溢洪道等水工建筑物全面布置，以及對施工和運(yùn)行期的安全都存在隱患，評價其在不同工況下的穩(wěn)定性，提出合適的加固處理措施。選取2號拉裂松動巖體為對象進(jìn)行穩(wěn)定性研究。

2號拉裂松動巖體分布于Ⅲ區(qū)，即Ⅶ-Ⅸ線下游側(cè)，拉裂松動巖體自Ⅶ線高高程斜向下延伸至Ⅸ線下游低高程江邊，呈現(xiàn)出從上游至下游自高向低的斜河向形態(tài)分布，橫向上受NE向陡傾角斷層的控制，其地表出露面積約6.3×104m2，體積約為300×104m3。該拉裂松動巖體內(nèi)地表發(fā)育的斷層有：N20-40°E/SE∠80-85°的近橫斜河陡傾的f11-f14共4條Ⅲ級結(jié)構(gòu)斷面(圖10)。

圖10 2號拉裂松動巖體平面地質(zhì)

選取橫Ⅶ剖面對2號拉裂松動巖體的穩(wěn)定性作計算。其計算剖面簡圖見圖11。

各地層物理力學(xué)參數(shù)見表4。

圖11 2號拉裂松動巖體橫Ⅶ剖面邊坡計算簡圖

表4 巖體物理力學(xué)參數(shù)表

4.2 計算成果分析

利用嚴(yán)格極限平衡法，計算邊坡在天然、暴雨、地震及地震+暴雨4種工況下的安全系數(shù)。暴雨工況表層參數(shù)c，f折減至90%，容重取飽和容重，地震工況的加速度為0.173 g，采用擬靜力法來計算，利用Slide軟件自帶的Bishop法和Janbu法計算不同工況下滑面的安全系數(shù)，結(jié)果見表5。

表5 計算結(jié)果匯總

注：表5中數(shù)值表示差值相對于Janbu法結(jié)果的百分比。

按照前述的方法計算不同工況下滑面的應(yīng)力分布函數(shù)，按照式(29)進(jìn)行修正并進(jìn)行擬合，結(jié)果見圖12-圖15。

天然工況下圓弧滑面正應(yīng)力分布見圖12。

則天然狀態(tài)圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.00008x3-0.051x2+8.1193x+187.36

帶入平衡方程中，解出該工況下的安全系數(shù)為1.193。

暴雨工況下圓弧滑面正應(yīng)力分布見圖13。

則暴雨工況下圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.00007x3-0.0528x2+9.3091x+64.77

帶入平衡方程中，解出該工況下的安全系數(shù)為1.126。

地震工況下圓弧滑面正應(yīng)力分布見圖14。

則地震工況下圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.0001x3-0.0678x2+10.537x+45.782

帶入平衡方程中，解出該工況下的安全系數(shù)為0.972。

地震+暴雨工況下圓弧滑面正應(yīng)力分布見圖15。

則地震+暴雨工況下圓弧滑面正應(yīng)力分布函數(shù)為：

y=0.00006x3-0.0481x2+8.9016x+53.2541

帶入平衡方程中，解出該工況下的安全系數(shù)為0.897。

圖12 天然工況下圓弧滑面正應(yīng)力分布

圖13 暴雨工況下圓弧滑面正應(yīng)力分布

圖14 地震工況下圓弧滑面正應(yīng)力分布

圖15 地震+暴雨工況下圓弧滑面正應(yīng)力分布

不同方法、不同工況計算結(jié)果匯總見表5。

通過對牙根水電站2號拉裂松動巖體采用基于數(shù)值應(yīng)力場的滑面正應(yīng)力修正法，分析邊坡在不同工況下的穩(wěn)定性以及得出相對應(yīng)的安全系數(shù)，并與Bishop法、Janbu法的安全系數(shù)作對比，可得出如下結(jié)論：3種極限平衡方法的計算結(jié)果相差在3%以內(nèi)，三者得出的結(jié)論基本吻合。將滑面正應(yīng)力修正法與有限元法相結(jié)合，理論明確，操作簡單，計算結(jié)果精度較高，適合大規(guī)模推廣。

5 結(jié) 語

本文采用滑面正應(yīng)力修正法與數(shù)值應(yīng)力場的分析方法對邊坡進(jìn)行穩(wěn)定性分析，其結(jié)果與傳統(tǒng)剛體極限平衡法的分析結(jié)果相比較相對誤差較小，進(jìn)一步證明該方法可行性。通過對不同工況下的邊坡進(jìn)行穩(wěn)定性計算發(fā)現(xiàn)，邊坡在天然狀態(tài)下安全系數(shù)較小，當(dāng)遇到暴雨或地震時安全系數(shù)下降幅度大，需要及時治理。通過對牙根水電站2號拉裂松動巖體采用基于數(shù)值應(yīng)力場的滑面正應(yīng)力修正法，分析邊坡在不同工況下的穩(wěn)定性以及得出相對應(yīng)的安全系數(shù)，并與Bishop法、Janbu法的安全系數(shù)作對比，可得出結(jié)論：3種極限平衡方法的計算結(jié)果相差在3%以內(nèi)，三者得出的結(jié)論基本吻合。將滑面正應(yīng)力修正法與有限元法相結(jié)合，理論明確，操作簡單，計算結(jié)果精度較高，適合大規(guī)模推廣。