雙參數強度折減法研究中存在的問題分析

陳子玉，宋彥輝, 2，嚴 豪，陳康達

(1. 長安大學地質工程與測繪學院，陜西 西安 710054； 2. 西部礦產資源與地質工程教育部重點實驗室，陜西 西安 710054)

強度折減法自Zienkiewicz于1975年提出以來，在邊坡穩定性分析中廣泛應用[1]。相對于極限平衡法，強度折減法不需要事先假定邊坡失穩滑動面，通過降低巖土體強度參數的方法，結合有限元方法，得到邊坡失穩破壞時的應變和應力圖。強度折減法研究大都基于庫倫-摩爾強度準則進行強度折減計算[2-5]，在參數折減時對c，φ同步折減 (Strength Reduction Method,以下簡稱為SRM)。近年來有學者根據巖土體在變形破壞時強度參數c，φ的變化規律，采用雙參數折減方式進行強度折減法研究[6-11](Double Reduction Method,以下簡稱為DRM)，在具體計算時，巖土體的c，φ值不同步折減。但是這些結果常常沒有足夠的理論支撐，同時眾多研究的計算理念互不相通，強度參數c，φ數值的衰減規律和DRM中對c，φ的不同步折減不同步，使得DRM研究在“百花齊放”的同時，其結果互相矛盾。DRM是否合理，計算結果是否適合工程應用，尚沒有充分的理論研究和實驗驗證。

相對SRM，DRM由于存在2個折減系數Rc，Rφ，其研究過程中需要確定雙折減系數Rc，Rφ之間的關系，然后通過一定的數學運算求得綜合穩定性系數F(Rc，Rφ)。DRM研究過程中面臨的主要問題包括三個方面：(1)折減路徑選擇；(2)綜合安全系數確定；(3)結果驗證。這些問題給DRM的研究和應用帶來了很多質疑和限制。強度折減法已經在國家規范[12]中得以應用，我國工程界開始正式將其納入設計和工程治理中，關系到萬千工程的安全穩定性，DRM作為一種新興的邊坡穩定性分析方法，受到越來越多的關注，因此有必要進行研究和探討DRM應用的合理性。本文從分析強度參數c0，φ0的物理意義著手，就強度折減法本質進行討論，同時逐一分析DRM存在的問題，以討論和研究DRM應用的合理性和工程適用性。

1 DRM原理及強度折減本質討論

1.1 強度折減本質討論

有限元強度折減法對邊坡穩定性系數的定義與簡化畢肖普法相同[2]，兩者都采用強度儲備的概念計算邊坡穩定性系數，通過對邊坡巖土體強度參數折減迭代計算，邊坡處于極限平衡狀態時折減系數的大小等同于穩定性系數。通常在描述巖土體強度參數時有抗拉、抗壓、抗剪等不同的強度參數，在邊坡穩定性分析中，抗剪強度是主要控制因素。破裂面上，巖土體的抗剪強度值τf是一個關于σn，c0和φ0的函數：

τf=τ(c0,φ0,σn) (1)

傳統有限元強度折減法在采用莫爾-庫倫強度準則來表示抗剪強度大小時，對巖土體強度參數c0，φ0折減，莫爾-庫倫強度理論公式如下：

τf=c0+σntanφ0(2)

式中：τf——莫爾庫侖式確定的抗剪強度；

c0，φ0——巖土體黏聚力和內摩擦角；

σn——邊坡滑動面上正應力。

從式(2)中可知，巖土體抗剪強度由黏聚力和摩擦強度組成。巖土體黏聚強度c0由靜電引力、范德華力、顆粒間膠結、顆粒間接觸點的化合價健、表觀黏聚力等構成，摩擦強度σn·tanφ0由滑動摩擦和咬合摩擦構成。c0，φ0值是巖土體根據同一實驗手段得到的，兩者之間互相關聯，但又不是簡單的疊加，在不同實驗手段中，如常規的三軸實驗和直剪實驗中，不同的應力水平和排水條件下測量出的結果不同，而且兩者互相隱含，既有區別的一面，同樣有相同的一面[13]。因此c0，φ0是計算巖土體抗剪強度的2個參數，且在物理實質上難以區分[14]。

此外，有限元強度折減法在計算過程邊坡變形破壞的過程中，應力場不斷變化[15]，即式(2)中的應力值σn在不斷調整。而巖土體邊坡中c0，φ0值受到邊坡應力狀態的影響，即不同的圍壓、軸向應力、不同的中間主應力[16]作用下，測得的c0，φ0數值不同。因此強度參數值隨著應力場在不斷調整。

因此對于一個確定的邊坡穩定性問題，邊坡內部應力場確定，其強度參數c0，φ0大小為定值，在計算時根據抗剪強度τf的折減進行計算，更符合強度折減法的本質。從式(2)中可以看出，SRM對式(2)左側進行折減時，式(2)右側正好對應了c0，tanφ0進行同步折減，這就表明SRM的本質是對抗剪強度進行折減，折減系數與強度儲備系數相同，此時折減系數即為邊坡的穩定性系數。而DRM計算中并沒有對抗剪強度進行等比例折減，此時的折減系數與強度儲備系數不同。

1.2 DRM原理

DRM通過巖土體初始強度參數c0，φ0采取不同的折減系數進行強度折減計算，具體計算如下：

(3)

(4)

式中：τs——剪切應力；

L——邊坡滑動面長度；

Rc，Rφ——c0，φ0的強度折減系數。

當Rc=Rφ時，DRM變化為SRM。式(5)為1時，表示邊坡處于臨界破壞狀態，此時邊坡的穩定性系數為1。

2 DRM存在的問題

2.1 折減路徑選擇問題

在DRM中，巖土體強度參數c0，φ0各自的折減系數Rc，Rφ之間的比例關系是需要首先確定的關鍵問題。通過折減比K=Rφ/Rc的選擇，式(3)、(4)中雙折減系數Rc，Rφ簡化成1個，此時DRM折減路徑的選擇問題變成折減比K的確定問題。

朱彥鵬等[17]、康校森[18]、唐芬等[19]采用固定比例關系，根據強度參數在衰減過程中變化規律，使得K=Rφ/Rc保持常數不變，在具體計算時根據不同土體取不同的K值。文獻[10]中根據c-tanφ臨界狀態曲線，采取統計學的方法得到固定的折減比K。薛海濱等[11]基于巖土體強度的弱化規律，根據簡單的線性軟化模型，采用非等比例關系，即Rc，Rφ之間呈一定的函數關系：

(6)

式中：cr，φr——巖土體殘余階段強度參數。

以上方法中對雙折減系數的比例確定多基于浸水軟化[6,18-19]、巖土體應變軟化[11]或巖土體自然劣化規律[8]等，或采用統計學的方法[10]。巖土體抗剪強度參數的變化規律受到多種因素影響，這些研究成果受到樣本值及不同巖土體之間差異性等影響，僅根據以上這些因素并不能準確得出c0，φ0的變化規律。式(6)采用巖土體應變軟化模型得出的折減比例，是在邊坡巖土體在穩定到極限平衡的過程中，巖土體強度參數由c0，φ0線性降低到cr，φr的假定上得出，但是有限元強度折減法的計算過程一般基于有限元計算不收斂得出，其計算終止時并不一定對應于強度參數值cr，φr，且巖土體的軟化模型不局限于線性模型。因此以上對折減路徑的選擇不具有普遍適用性。

2.2 DRM計算原理

(1) 在有限元計算時，首先建立分析時步與場變量之間關系：

m=a+bt(7)

式中：m——場變量；

t——分析時步，且0≤t≤1；

a——場變量初始值，計算時一般設為1；

b——系數，根據邊坡算例實際情況試算或經驗取值。

(2) 將m與tanφ0設為反比例關系，即將DRM中內摩擦角φ0的折減系數Rφ作為場變量的值：

(8)

此時黏聚力c0與m之間關系如下：

(9)

(3) DRM計算時折減系數與場變量關系如下：

Rφ=m=a+bt(10)

由上述可知，雙折減系數Rc，Rφ均為場變量的函數。DRM計算時，通過場變量m的變化，實現了強度參數c0，φ0的折減，設置最小時步為△tmin，假設經過n步計算有限元計算收斂，此時場變量為mt，內摩擦角折減系數Rφ=mt，黏聚力折減系數Rc=mt/K。當K=1時，上述DRM計算過程和SRM計算過程相同。

對比SRM和DRM計算過程，結合式(10)和式(11)可知，在每一個時步計算時，DRM可以看成是以Kc0，tanφ0為初始強度參數進行SRM計算。為了進一步說明DRM與SRM之間的內在聯系，研究DRM的本質，增加來自文獻[19]中經典邊坡算例如圖1所示，算例中邊坡坡高20 m，坡度45°，巖土體材料參數取值見表1。邊坡計算時以有限元計算不收斂作為失穩判據，采用非相關流動法則，采用莫爾-庫倫屈服準則，三角形網格單元，模型兩側及底面固定，模型頂部自由。

圖1 算例有限元模型Fig.1 FEM model

γ/(kN·m-3)c/kPaφ/(°)μE/kPa2042170.3105

首先對算例邊坡進行DRM計算，分別給定K值為0.2，0.6，1.0，1.2，1.6，2.0，2.2，2.6，3.0，共9組，計算結果如表2所示。然后針對不同K值，改變算例邊坡中的強度參數，將黏聚力c0乘以K，內摩擦角φ0不變，以Kc0,φ0為初始狀態參數進行SRM計算，K取值和DRM相同，分別為0.2，0.6，1.0，1.2，1.6，2.0，2.2，2.6，3.0共9組，SRM計算采用的初始狀態參數和邊坡破壞時臨界狀態參數及折減系數如表3所示。可以看出DRM計算出折減系數Rφ與SRM折減系數R完全相同，且其臨界狀態參數相同。因此DRM在折減計算啟動時，可以看成是以(Kc0、tanφ0)為參數的SRM，即DRM的本質還是SRM。DRM相對SRM并沒有計算原理上的改變，僅相當于改變了原始狀態參數的SRM，依據文獻[9]提出的參照邊坡概念，將臨界狀態邊坡作為參照邊坡，待評價的邊坡作為原始邊坡，依據上述分析可知，DRM中參照的原始邊坡出現了改變，使得雙參數強度折減法對邊坡的數值模擬出現了重大失真。

表2 不同K值下DRM計算結果Table 2 DRM results of various K

表3 不同初始狀態參數下SRM計算結果Table 3 SRM results of various initial state parameters

2.3 DRM對邊坡潛在滑動面的影響

取出邊坡潛在滑動面上一個單元點(圖2)，單元點應力按照主應力形式表示如下：

圖2 單元受力情況Fig.2 Stress of unit

式中：σ1，σ3——破壞面上一點的最大主應力、最小主應力；

α——最危險滑動面AB與σ3之間的夾角。

由式(12)通過變換可以得到最危險滑移面角度α計算公式：

(13)

由式(12)與式(13)可以看出，對于一個確定的邊坡，邊坡巖土體強度參數c0，φ0已知，邊坡內部應力場確定時，邊坡潛在滑動面是確定不變的。對式(13)右側的強度參數按照SRM進行等比例折減，得到式(14)，雖然式(14)右側進行了參數折減，但是式(14)與式(13)計算結果相同。式(14)表明SRM在計算時，等比例折減強度參數c0，φ0并未改變邊坡滑動面。

(14)

對式(13)右側c0，φ0分別以Rc，Rφ進行折減，得到式(15)，可以看出式(13)與式(15)的結果明顯不同，即DRM計算時邊坡滑動面的形狀發生了變化。

(15)

對前文中邊坡算例采用DRM計算，在不同K值下，有限元計算臨界狀態最大剪應變帶云圖如圖3所示。

圖3 不同K值下臨界狀態最大剪應變云圖Fig.3 Maximum shear strains of the critical state under different values of K

由圖3可以看出，當K=0.2時，最大剪應變帶靠近坡面，剪出口位于坡腳上方；K=1.0時對應于SRM計算結果；隨著K值增大，K=2.0、K=3.0時最大剪應變帶加深,最大剪應變帶向坡內移動。因此DRM在計算時，隨著K取值變化，邊坡潛在滑動面出現變化，K值增大，滑動面加深，K值減小，滑動面靠近坡面，剪出口上移。

袁維等[20]指出DRM在計算時由于對c0，φ0的非等比例折減，其滑動面破壞形式和SRM不同，并以此作為DRM優于SRM的證據。對于一個穩定邊坡，其穩定性系數大于1，其破壞雖然存在無數種可能，而根據潘家錚最大、最小值原理[21]，其破壞面會沿著抗滑力最小的滑動面發生，即邊坡在物理力學參數和外部邊界確定時，在不受外界影響時，其潛在破壞面唯一確定。而根據林杭[22]、王少華等[23]的研究，邊坡巖土體強度參數c0，φ0的相對大小，控制著邊坡剪切破壞面的形式：在一個外部邊界固定的邊坡中，內摩擦角φ0值減小，滑動面加深；黏聚力c0減小，臨界滑動面靠近坡面。綜合以上分析，結合式(15)可知，DRM由于在折減計算時，由于折減比K的作用，改變了邊坡的原始參數，從而導致了計算出的滑動面并非邊坡潛在滑動面。

3 綜合穩定性系數及結果驗證

在工程實踐中，一般采用固定的穩定性系數來評價邊坡穩定性，SRM將有限元計算不收斂時的折減系數作為邊坡的穩定性系數，而DRM在計算時強度折減系數有兩個，對一個確定的邊坡，邊坡穩定性系數應當如何確定？不同的研究者根據不同的研究方法提出了邊坡穩定性系數的確定方法。

唐芬等[6]、陳依琳等[24]采用平均值確定綜合穩定性系數F:

(16)

朱彥鵬等[17]、袁維[20]采用Rc，Rφ乘積開方確定綜合穩定性系數F:

(17)

YUAN等[25]采用統計學的方法，通過利用大量不同形態特征的邊坡數值算例，采用數據擬合的辦法確定邊坡的綜合安全系數F：

(18)

朱彥鵬等[26]還提出了采用Rc，Rφ中最小值來確定綜合穩定性系數F：

F=min(Rc，Rφ) (19)

Isakov等[27]、趙煉恒等[28]認為，綜合穩定性系數F為Rc，Rφ兩者的隱含表達式，通過c0，φ0的不同折減系數Rc，Rφ之間的函數關系式，利用雙參數強度折減最短路徑來定義綜合穩定性系數F：

(20)

薛海濱等[11]、王強志等[29]將綜合穩定性系數F定義為邊坡土體強度參數折減前滑動面提供的抗滑力與極限狀態下滑動面提供的抗滑力比值，從而根據黏聚力與摩擦強度的權重確定出綜合穩定性系數。

但是上述研究方法在理論上差別較大，在對最終穩定性系數確定上理論依據不充分。如式(16)、式(17)在確定邊坡綜合穩定性系數時，其計算式并沒有明確的物理意義，并且缺乏理論證明。式(19)由于采用最小值來作為邊坡穩定性系數，忽略了另一個參數的強度貢獻。式(18)、式(20)雖然有公式的推導，但是其穩定性系數F值的物理意義并不明確，而且基于數據擬合的公式推導，導致其適用性受到樣本值的影響而難以推廣。DRM認為黏聚力c0和摩擦系數tanφ0發揮作用不同，并因此采用不同的折減系數，但是在確定綜合穩定性系數時，并沒有體現出黏聚力c0和摩擦系數tanφ0發揮作用不同的影響，各種綜合穩定性系數確定方法在原理和結果上不能統一。

另外，DRM在綜合穩定性系數的確定方法上難以應用在多層巖土體邊坡，已有文獻資料對這一問題均沒有提及。在實際工程實踐中，純粹的均質土邊坡基本不存在，即使對于單一巨厚沉積層的邊坡，也會由于地下水位、卸荷裂隙、風化等原因，呈現出分層性。在雙參數強度折減法計算過程中，對于多層土或含裂隙、節理的巖質邊坡，多層土之間的強度參數ci，φi各不相同，含節理、裂隙巖質邊坡巖塊和節理面強度參數不同。按照DRM理論，假設某邊坡內有n層巖土類型，n種類的巖土體都會對應2n個折減系數，此時邊坡的穩定性系數無法用式(16)～(21)所給定的方法確定出一個單一的綜合穩定系數F。

SRM一般在計算機有限元軟件中進行應力、應變分析，并根據三類判據的不同在一定精度范圍內和極限平衡法做對比(一般取5%為界限)。這是由于極限平衡法已經在實際工程建設領域得到了檢驗，而5%的計算誤差在工程允許的范圍內。但在某些研究中，將DRM的綜合穩定性系數與SRM計算結果進行對比，或者與其他研究者的不同綜合穩定性系數確定方法之間作對比，其精度誤差在2%以內[8]，雖然在一定程度上計算結果與極限平衡法和傳統強度折減法一致，但是依前述分析可知，兩者之間并無本質上的不同，偏離了邊坡穩定性計算的本質。

4 討論

DRM相比SRM固然可以反映巖土體強度參數c0，φ0的變化規律，但是在具體的研究和應用上偏離了強度折減的本質，存在著物理意義不明確、概念混淆的問題。在不同的數值分析軟件中，SRM只需要確定出數值計算收斂、位移突變、塑性區貫通三種破壞判據中的某一種即可求得邊坡的穩定性系數，目前已經在國內外的工程實踐和理論研究中得到了大量應用，而DRM相較于SRM多了折減比K的確定、綜合穩定性系數確定等兩個步驟，此外由前文分析可知，DRM不適用于多層邊坡的穩定性分析。SRM和DRM作為一種數值計算手段，相較于傳統極限平衡法，可以對邊坡內部應力及應變場進行分析無疑是一種進步，但是同樣地由于建模及數值分析手段的限制，在具體計算時參數的輸入及結果輸出上引起了很多爭論。在邊坡穩定性分析的工程實踐中，首要的是確定出邊坡的穩定性系數和潛在滑動面，而DRM由于本身計算原理和特性，其計算結果受到折減比K的影響而呈現多樣性，使得該方法更加復雜和難以推廣，強度折減法(SRM)的研究應該更多地針對強度準則選用、輸入參數的選取、破壞判據的選擇及輸出結果分析等內容。

5 結論

(1)雙參數強度折減法與傳統強度折減法無本質上不同，有限元計算時兩者計算原理相同，但雙參數強度折減法與基于強度儲備的穩定性系數計算方法不符。

(2)折減比K的存在改變了邊坡潛在滑動面的形態，不同的K值會得到不同的邊坡臨界滑動面，與邊坡實際不符。

(3)已有的雙參數強度折減法研究對多層巖土體邊坡不適用，無法確定出多層巖土體邊坡的穩定性系數。