劉 亞 坤

(西南交通大學土木工程學院，四川 成都 610000)

0 引言

巖土工程中分析邊坡穩(wěn)定性方法主要可分為傳統(tǒng)的極限平衡法和有限元法。傳統(tǒng)的極限平衡法以剛體極限平衡理論為基礎，分析結果物理意義明確，但需事先假設滑動面的形狀、位置，然后不斷求解對應的安全系數，直到最小安全系數，此時對應的滑動面為臨界滑動面。該方法不僅計算量大，而且難以保證臨界滑動面的準確性，尤其是對復雜地質邊坡的實用價值大打折扣；而有限元法不僅能夠模擬邊坡的變形過程及其滑動面的形狀，搜索臨界滑動面時無需事先假定滑動面形狀和位置[1-3]。

Davidon[4]和CHEN Z Y & SHAO C M[5]基于最優(yōu)控制理論，提出了能最快速搜索目標的變度量法，這為邊坡臨界滑動面快速搜索提供了新思路，但其計算精度必須建立在準確確定初始滑動面的基礎之上。已有研究表明，有限元強度折減法能較為精確的確定初始滑動面。

為此，本文借鑒有限元強度折減法獲取初始滑動面、局部優(yōu)化法的單變量精確計算臨界滑動面和最快速獲得搜索目標的邊度量法的優(yōu)點，提出了一種能編程簡單、搜索效率較高且能搜索任意形狀邊坡的臨界滑動面交替變度量法搜索方法。

1 交替變度量法基本原理

1959年，D’esopo[4]把多變量求最優(yōu)值問題轉化成單變量輪流求最優(yōu)值問題，從而提出了坐標輪換法——以節(jié)點坐標為變量，每次搜索只改變一個變量，其余變量保持不變，即沿坐標方向輪流進行搜索的尋優(yōu)方法；1959年Davidon[4]提出變度量法。CHEN Zuyu, SHAO Changming[5]于1988年將最優(yōu)控制理論中變度量法改進，成功的用于邊坡滑動面的搜索中。本文借鑒坐標輪換法及變度量法的原理，將多變量問題轉化成單變量問題，利用變度量法逐一搜索每個節(jié)點的最優(yōu)值，從而提出了交替變度量法。

為確定臨界滑動面，利用有限元強度折減法確定等效塑性應變區(qū)域，利用樣條函數沿最大塑性應變區(qū)域擬合一條光滑曲線(最大塑性應變值連線為波浪形[6])，作為初始滑動面。在初始滑動面上確定若干節(jié)點，先令其余節(jié)點固定不變，使第一個節(jié)點(坡腳節(jié)點)沿邊界方向進行搜索，獲取最小安全系數，如圖1a)所示；然后沿著最優(yōu)化方向移動另一節(jié)點，獲取最小安全系數，其余節(jié)點仍然保持固定，見圖1b)；滑動面上所有節(jié)點均按照這種方法搜索最小函數值，直至得到一次循環(huán)的最小安全系數。

利用交替變度量法確定邊坡最小安全系數和臨界滑動面的具體步驟如下：

1)通過樣條函數沿最大塑性應變值擬合的滑動面作為初始滑動面。

2)i為循環(huán)變量，第一次循環(huán)記為i=1。

3)在初始滑動面上確定若干節(jié)點，節(jié)點位置需分布均勻且滿足運動學規(guī)律，即節(jié)點之間距離盡可能相等。節(jié)點間的距離過小則無法滿足搜索靈活性，距離過大，則可能無法滿足運動學規(guī)律。節(jié)點之間用直線段連接，初始滑動面上節(jié)點的坐標為Zi=(x1,y1,x2,y2,…,xm,ym)。

繼續(xù)對第二個節(jié)點進行沿最優(yōu)化方向搜索，其余節(jié)點不變，直至獲得最小安全系數，最優(yōu)化方向記為Sk，按式(1)確定。

(1)

其中，Hk為對稱矩陣。本文節(jié)點移動的步長λ采用“黃金分割”法沿Sk方向進行一維搜索。第二次搜索后臨界滑動面上的節(jié)點坐標變?yōu)椋?/p>

5)按步驟4)對除最后一個節(jié)點外的剩余節(jié)點依次進行搜索，直至第(m-1)個節(jié)點搜索完畢，此時臨界滑動面的坐標變?yōu)椋?/p>

6)對最后一個節(jié)點(第m節(jié)點)沿邊界移動搜索，直到獲得對應的最小安全系數值。

7)重復4)～7)步，當最后兩次循環(huán)的結果差值滿足式(2)時，循環(huán)搜索結束，最后一次循環(huán)所得的坐標連線即為臨界滑動面，所得最小函數值即為邊坡安全系數。

2 算例分析

對澳大利亞計算機應用協會(ACADS)提供的均質邊坡算例[7]進行計算驗證，并與幾種常用的極限平衡法計算結果進行了比較，驗證了上述方法的正確性及優(yōu)越性。

該例來自ACADS考核題EX1(a)，土體重度γ=20.0 kN/m3，黏聚力c=3.0 kPa，內摩擦角φ=19.6°，彈性模量E=1.0×104kPa，泊松比v=0.25，具體材料參數、幾何參數見I.B.Donald和P.Giam[7]的研究。分別采用屬于圓弧法的瑞典法、簡化Bishop法以及屬于非圓弧法的Morgenstern-Price法與本文方法進行最危險滑動面搜索和安全系數計算。

圖2為初始滑動面位置及計算節(jié)點的布置。由文獻[6]知，臨界滑面是由沿深部方向上的等效塑性應變的極大點所組成，但最大塑性應變值的連線成波浪形，故本文采用樣條函數沿最大塑性應變區(qū)域擬合一條光滑曲線作為初始滑動面。然后等間距選取8個計算節(jié)點，通過交替變度量法搜索每個節(jié)點最優(yōu)值(最小安全系數值對應的節(jié)點坐標)。

圖3為交替變度量法及幾種常用傳統(tǒng)極限平衡法搜索的臨界滑面位置，可以看出，本文方法確定的臨界滑動面總體上要比圓弧法確定的臨界滑動面更深，但比Morgenstern-Price確定的臨界滑動面更淺。本文所得結果為折線型滑動面，計算節(jié)點僅為8個，但已然很接近于圓弧狀滑面，若增加滑動面計算節(jié)點數，則最終搜索結果基本等同于圓弧狀滑面。

表1給出了幾種傳統(tǒng)極限平衡法及本文方法計算得到的安全系數，雖然本文方法計算結果小于算例給出的答案，但只相差1.2%，比瑞典法、簡化Bishop法都接近算例答案，幾乎與Morgenstern-Price法結果相同。

表1 安全系數計算成果匯總表

表1還總結了不同方法嘗試搜索滑動面的次數?？梢钥闯?，傳統(tǒng)極限平衡法中瑞典法搜索次數較少，為140次，Bishop法搜索次數為400次，Morgenstern-Price法搜索次數明顯大于其他傳統(tǒng)極限平衡法的搜索次數，達到了667次。從表1搜索時間還可以看出，本文方法僅用時3 min，相較瑞典法、Bishop法和Morgenstern-Price法分別降低了9.33，14.00和25.33倍。結合表1的安全系數可知，瑞典法的精度較差，安全系數為0.942，Bishop法精度適中，安全系數約為0.986，Morgenstern-Price法精度最高，達到了0.989，而采用本文方法，計算程序只需循環(huán)搜索15次，最后兩次計算結果的差值就可滿足誤差要求|Fs(Zi+1*-Fs(Zi*)|<1×10-5，此時安全系數達到0.988。可見，與傳統(tǒng)極限平衡法相比，本文方法在保證精度的同時能大大提高搜索效率。

傳統(tǒng)極限平衡法假設在整個臨界滑動面上安全系數相等，這并不符合實際情況，邊坡的破壞并不是整體同時破壞而引起的，而是局部破壞的逐步擴展從而引起整體失穩(wěn)，本文方法恰好能反映這一點。圖4給出了臨界滑動面上安全系數的變化，每點位置與圖3坐標表示一致。由圖4可以看出，邊坡坡腳處較為穩(wěn)定，沿著坡腳到坡頂安全系數逐漸降低，接近坡頂時穩(wěn)定性又逐漸增加，與邊坡實際破壞機理基本一致。

3 結語

1)通過有限元強度折減法確定初始滑動面，再利用交替變度量法可以快速、有效的搜索臨界滑動面。通過ACADS的邊坡算例驗證，其在邊坡滑動面的搜索中的應用是可行的，且搜索效率明顯高于傳統(tǒng)極限平衡方法。

2)利用交替變度量法搜索結果雖然為折線型滑動面，但只要增加計算節(jié)點，就能很好的擬合為圓弧狀滑面。因此利用交替變度量法可以搜索任意形狀的滑動面。為了提高計算效率，在保證準確性的情況下，計算節(jié)點數量以8個～15個為宜。

3)本文方法以有限元強度折減法作為計算基礎，可以真實反映滑動面不同位置安全系數的變化，符合邊坡破壞的實際機理。

參考文獻：

[1] 鄭 宏,李春光,葛修潤.求解安全系數的有限元法[J].巖土工程學報,2002,24(5):626-628.

[2] 李春光,朱宇飛,劉 豐,等.基于下限原理有限元的強度折減法[J].巖土力學,2012,33(6):1816-1821.

[3] 楊光華,鐘志輝,張玉成,等.用局部強度折減法進行邊坡穩(wěn)定性分析[J].巖土力學,2010,31(S2):53-58.

[4] 馬昌鳳.最優(yōu)化方法及其Matlab程序設計[M].北京:科學出版社,2010.

[5] CHEN Zuyu,SHAO Changming.Evaluation of factor of safety in slope stability analysis[J].Canadian Geotechnical Journal,1988,25(4):735-748.

[6] 孫冠華,鄭 宏,李春光.基于等效塑性應變的邊坡滑動面搜索[J].巖土力學,2008,29(5):1159-1163.

[7] I.B.Donald, P.Giam. The ACADS slope stability programs review. Proc[A].The 6th International Symposium on Landslides[C].1992:1665-1670.