能量法結合徑向條分法求邊坡穩定性

莫智豪，劉保東，朱劍宏，饒平平

（1.上海理工大學 環境與建筑學院，上海 200093；2.南寧城建管廊建設投資有限公司，南寧 530219；3.廣西路橋工程集團有限公司，南寧 530011）

邊坡是一種自然地質體，在受到外界因素影響時，邊坡土體會沿著一些不穩定的結構面滑動，導致邊坡的失穩。邊坡穩定性分析是巖土工程的重要研究內容，已經形成一個應用研究課題，其穩定性問題涉及采礦工程、水利和水電工程等許多工程領域，近年來受到越來越多的關注。邊坡穩定性分析方法有很多種，其中基于定性判斷的工程地質分析方法得到廣泛應用，還有基于剛性定量計算的極限平衡法和基于數值計算的強度折減法[1-2]。工程地質分析方法主要分析了邊坡的成因歷史，研究了邊坡所處的地質環境、形成的地質歷史，以及影響邊坡穩定性的變形破壞軌跡之間的關系，從而對其穩定狀況作出宏觀評價[3]。極限平衡法基于Mohr-Coulomb 屈服準則和靜力平衡，根據邊坡潛在破壞面上抗滑力與滑動力的比值來確定邊坡穩定系數，然后合理評估斜坡的穩定性[4]。強度折減法是有限元中引入的計算安全系數的方法，通過對邊坡土體參數進行折減，使邊坡達到最終的極限平衡狀態，根據該方法可以得到安全系數和潛在破壞面的位置[5]。

以上的邊坡穩定性分析方法均在工程上有廣泛的應用，且獲得了工程界的認可，但各自均存在明顯的不足。工程地質法類比條件因地而異，經驗性強，沒有數量界限，需要具有豐富的工程實踐經驗的地質工作者才能作出較為準確的判斷。通過極限平衡法難以確定破壞面的形狀和位置，這給邊坡的穩定性分析帶來了很大的困難。此外，極限平衡法沒有考慮土的應力-應變關系，計算的土體條間力或滑塊底部的反作用力不能代表滑坡滑動時的真實力，這給邊坡穩定性分析帶來了很大的不確定性，無法判斷邊坡的變形破壞模式。在強度折減法中，降低巖土強度參數c和φ的方法缺乏理論依據，并且以相同的比例減小，還是以不同的比例折減等問題仍缺乏合理依據，但目前還是得到了廣泛的應用與認可。

本文使用的能量法從極限分析的上限定理演化而來。極限分析方法最初由Drucker 等[6]提出，為斜坡的穩定性分析開辟了一條新的道路。Chen[7]首次討論了巖土結構的穩定性分析，并討論了極限分析在巖土邊坡穩定性中的應用。Michalowski[8-9]使用楔形體和豎直條塊平動破壞機構以及轉動機構對邊坡進行了上限分析。王根龍等[10-12]采用塑性極限分析法結合條分法和等分圓弧面法對土坡進行了分析。趙煉恒[13]、陳靜瑜等[14]基于極限分析上限定理，考慮孔隙水壓力的影響，提出折線形滑面邊坡穩定性計算模型。易慶林等[15]根據虛擬滑動位移，結合地下水的影響，提出滑坡穩定性計算模型。

本文從能量守恒的角度出發，建立了邊坡穩定性計算模型。針對滑動面為對數螺旋面的邊坡，采用徑向條分法，虛設虛擬轉角位移，考慮塊體間相對位移產生的內部耗能，計算滑坡重力勢能和滑面摩擦內能的變化，使用強度折減法計算邊坡穩定性。

1 基本原理和物理模型

1.1 基本原理

利用滑坡體和滑動面作為能量守恒法判斷邊坡穩定狀態的研究系統。假設滑動面上的土體沿對數螺旋滑動面的極點O產生非常小的角位移dθ，將其稱為虛擬轉角位移。在此過程中，重力勢能的變化量為ΔU，抗滑力（滑動面上黏聚阻力與滑動摩擦力）做功ΔW1，塊體分界面上由于土體自適應變形而相對運動產生的內能耗散（黏聚阻力所做的功）為ΔW2。

采用強度折減法[16]對土體參數c和φ進行折減，安全系數為

式中：c，tanφ為實際土體強度參數；cm，tanφm為維持邊坡穩定的極限強度參數。土體參數折減后邊坡達到極限平衡狀態，即

1.2 物理概化模型

本文所研究的邊坡穩定性計算物理模型如圖1所示。圖中：滑動面處于對數螺旋曲線AB上；點O是對數螺旋線的極點；點A是破壞面與邊坡上部平面的交點；點B為坡腳位置；θA和θB分別為A點和B點對應的角度；r0為A點所對應的初始極徑；H為邊坡高度；α為坡角。假設邊坡土體服從相關聯流動法則和Mohr-Coulomb 屈服準則的理想塑性變形假設，根據Chen[7]所述，采用對數螺旋面作為邊坡的滑動面。

圖1 邊坡模型Fig.1 Slope model

對數螺旋滑動面曲線方程為

根據極限分析上限法的原理和許多大型滑坡的運動特征研究，土體在運動過程中，下部土體的運動速度通常大于上部土體的運動速度，在此基礎上，如果仍然采用垂直條分法或斜條分法顯然不合適。本文采用徑向條分法研究滑動面上的土體，即從極點O中提取出一系列射線，將滑動面上的土體根據角度等分成n個小滑塊，每個小塊所占角度，如圖2 所示。

圖2 徑向條分后的邊坡模型Fig.2 Slope model after radial slicing

本文對能量法邊坡穩定性計算模型進行了以下簡化假設：a.假設各塊體為剛體；b.所研究的邊坡滑動面為對數螺旋滑動面；c.邊坡在產生虛擬轉角位移dθ的過程中，由于自適應變形各塊體間產生相對運動；d.不考慮相鄰塊體間土體自重所產生的摩擦耗能。

2 穩定性計算模型

2.1 重力勢能改變量

如圖2 所示，采用徑向條分法將邊坡對數螺旋滑動面上部的土體根據極角等分為n個塊體。根據塊體重心的變化來計算塊體重力做功，對于扇形結構其重心在其角平分線上，所以分析中取塊體中線進行計算。假設每個塊體的重心位于塊體中心線（塊體對應角度的角平分線）的中點，第i個塊體的質量為Gi。當某一時刻產生虛擬轉角位移dθ時，第i個滑塊的重心的位移為dSi，其對應的豎向位移為dyi，則重力勢能改變量為

將圖形簡化為如圖3 所示，計算第i個塊體的重力勢能改變量。∠QOA=θA，第i個塊體重心線與水平線OQ之間的夾角為位移后新重心線與水平線OQ之間的角度為

圖3 滑塊的重心位移示意圖Fig.3 Sketch of slider center under gravity changes

設初始重心Pi到極點O的長度為′，在轉動角度dθ后，Pi與之間的距離為

轉動的角度為dθ，由圖中幾何關系可知：

其中

所以

2.2 滑動面抗滑力做功

由于滑動面上部土體產生位移，且存在黏聚力和滑動摩擦力，則滑動面上會產生能量耗散。

a.滑動摩擦力做功

式中，F為滑塊自重在滑動面上的法向分力與摩擦系數μ的乘積，dsi=ridθ，r為滑塊中線的長度。

b.黏聚力做功

式中，c為土體的黏聚力，Δs=riΔθ（當n足夠大時，滑塊在滑動面上的部分可近似為圓弧形曲線）。

對第i個塊體進行計算，如圖4 所示。

對其運動幾何圖進行分析，如圖5 所示。

由幾何關系可得

圖4 滑動面上內能耗散示意圖Fig.4 Energy dissipation on the sliding surface

圖5 滑塊摩擦耗能計算幾何分解圖Fig.5 Geometric diagram of the calculated energy dissipation due to slider friction

2.3 相鄰塊體接觸面上的能量耗散

由于采用對數螺旋滑動面，其上部土體會產生變形，塊體間會產生能量耗散，如圖6 所示。

圖6 塊體徑向位移Fig.6 Radial displacement of the slider

所以第i個塊體與第i+1 個塊體接觸面相對位移

由于n個塊體間存在n-1 個接觸面，所以

3 算例分析

3.1 算例概況

為了說明本文的合理性，采用工程算例進行驗證。如圖7 所示邊坡，設邊坡土體均質且服從Mohr-Coulomb 屈服準則，土坡高度H=20 m，坡角α=60°，土體重度γ=19 kN/m3，摩擦角φ=24°，黏聚力c=35 kPa。

圖7 穩定性計算模型Fig.7 Stability calculation model

3.2 計算方法

采用徑向條分法，將對數螺旋滑動面上部土體根據角度等分為20 份，每份對應極角為（θA-θB）/20。假設虛擬轉角位移dθ=0.01°。將推導出來的表達式使用Matlab 軟件進行編程，使用fmincon優化程序對表達式中的變量θA和θB進行優化，得到最小安全系數FS，根據其對應的參數θA和θB即可得到最危險滑動面。

3.3 計算結果

通過本文介紹的能量法，計算求得穩定性系數FS=1.228。使用極限分析法對該算例進行建模和計算，得出穩定性系數FS=1.240，能量法和極限分析法的計算結果之間的差異為0.97%。通過本文方法得到的計算結果與常用的邊坡穩定性極限分析方法獲得的結果差別不大且更加準確，說明了本文方法的正確性，對邊坡的穩定性評估及設計施工有一定的參考價值。

4 參數影響分析

對于給定的邊坡模型，不同的參數對計算結果有很大的影響。基于算例1 中的邊坡穩定性計算模型，分析了條塊數n、虛擬轉角位移dθ、黏聚力c、摩擦角φ和坡角α對穩定系數的影響，并將計算結果與極限分析上限法得出的結果進行了對比分析。

4.1 條塊數n 對穩定系數的影響

從理論上來說，將滑動土體進行分塊，考慮其塊體間摩擦耗能，是對計算結果的進一步優化，且條塊數n取值越大，結果越精確。對算例1 中的邊坡穩定性計算模型，分別取條塊數n為5,10,15,20,25,30，模型中其他參數不變，計算邊坡穩定性系數，其計算結果如圖8 所示。從圖中可知，隨著條塊數n的增大，邊坡穩定性系數FS逐漸減小，且當n大于一定值后，邊坡穩定性系數FS趨向于收斂。由于計算過程比較復雜，因此在后面的例子中，模型的徑向條塊數均為20。

圖8 條塊數對穩定系數的影響Fig.8 Influence of number of slider blocks on the stability factor

4.2 虛擬轉角位移dθ 對穩定系數的影響

從理論上來說，虛擬轉角位移dθ數值越小，計算結果越精確。采用算例1 中的計算模型，令dθ取值0.001°,0.01°,0.1°,1°,10°，且其他參數保持不變，計算邊坡穩定性，結果如圖9 所示。從計算結果可以看出，在一定范圍內，dθ取值越小，邊坡穩定性系數FS越小，但當dθ小于某一數量級后，FS基本上趨于穩定。根據上述計算所得的結果，當dθ＜0.1°后FS值基本趨于穩定，所以后面的算例中均按照dθ=0.01°進行計算。

圖9 轉角位移對穩定系數的影響Fig.9 Influence of virtual displacement on the stability factor

4.3 物理參數c,φ,α 對穩定系數的影響

采用算例中的數據，分別對c,φ,α取不同值計算穩定性系數，并與極限分析法的計算結果進行了對比分析。如圖10～12 所示。

圖10 黏聚力對穩定系數的影響Fig.10 Influence of cohesive force on the stability factor

圖11 內摩擦角對穩定系數的影響Fig.11 Influence of friction angel on the stability factor

根據圖10 中的曲線，能量法與極限分析法計算出來的穩定系數相差不大，但能量法計算出來的結果更小。其變化規律基本一致，穩定系數隨著c值的增大而增大。根據圖11 中的曲線，邊坡穩定系數隨著土體內摩擦角φ值的增大而增大，但能量法的計算結果要小于極限分析法的計算結果。根據圖12 中的曲線，能量法與極限分析法穩定系數均隨坡角的增大而減小，且其減小規律基本一致，但能量法的計算結果更優。

圖12 坡角對穩定系數的影響Fig.12 Influence of slope angel on the stability factor

4.4 能量法與極限分析法的對比

本文的能量法與極限分析法有很多異同點。相同點是兩者滑動面均為對數螺旋滑動面，且將滑動面上部劃分的土體視為剛性塊體。不同點主要在于建模方法的不同，在建立邊坡物理計算模型過程中，能量法采用徑向條分法。從極點引出一系列射線，將滑動面上部土體劃分為若干徑向條塊，考慮塊體與塊體間的能量損耗，使安全系數更為精確。

從本質上來講，極限分析法是通過構建虛擬速度場，從功率的角度來求解穩定系數，而能量法是假設虛擬角度，從能量的變化角度來求解。上述兩種方法都是考慮滑塊體的重力和滑動面上的抗滑力來進行計算，故而兩種方法的計算結果應該是較為接近的。不同的是能量法通過虛擬轉角位移dθ，計算力與位移的乘積，位移前后邊坡狀態不同。極限分析法通過虛設角速度ω，計算內外做功功率，邊坡狀態未發生變化。

5 結論

本文根據邊坡失穩破壞過程中土體的運動模式，采用徑向條分法和對數螺旋滑動面建立物理模型，虛設滑動面上部土體轉角位移，采用能量法對邊坡的穩定性進行計算，總結出以下幾點：

a.通過建立能量法邊坡穩定模型，考慮相鄰塊體間的能量損耗，推導出滑動土體重力做功公式ΔU、滑動面上內能耗散公式ΔW1以及各塊體間的內能耗散公式ΔW2。結合強度折減法，計算邊坡穩定系數FS。通過算例計算，并與極限分析法的計算結果進行對比驗證。結果表明考慮內部能量耗散的能量法計算結果比極限分析法計算結果更優。

b.通過采用不同的土體參數與邊坡幾何參數，運用能量法與極限分析法分別計算穩定性系數，對比分析了各參數對邊坡穩定性系數的影響規律，可為工程實踐以及設計施工提供參考。

c.采用能量法計算邊坡穩定性時，條塊數n越大，計算結果越精確，虛擬轉角位移dθ越小，計算結果越精確，邊坡安全系數隨著內黏聚力與內摩擦角增大而增大，隨著坡角的增大而減小。

d.本文分析結果對工程實踐有一定參考價值，但也存在不足。本文所研究的滑動面為對數螺旋滑動面，但實際中滑動面更為復雜。在計算塊體間摩擦耗能時，只考慮了土體黏聚力的影響，未考慮各塊體內部土壓力造成的內能耗散。