含軟弱夾層邊坡穩定性上限分析法研究與評價

付一鳴 王 珍 萬 峰 張洪清

（1.遼寧工程技術大學礦業學院，遼寧 阜新 123000；2.本溪鋼鐵（集團）礦業有限責任公司歪頭山鐵礦，遼寧 本溪 117000；3.黃河科技學院建筑工程學院，河南 鄭州 450063；4.扎賚諾爾煤業有限責任公司，內蒙古 滿洲里 021410）

0 引言

軟弱夾層是指在邊坡巖土體內廣泛存在的不連續結構面，由于其力學性質薄弱，是導致邊坡發生滑坡災害的重要因素。

軟弱夾層對邊坡穩定性影響已成為當下巖土工程中重要課題。現階段，邊坡工程領域常用的分析方法包括極限平衡法、極限分析法等方法［1］。上世紀20年代，學者Fellenius提出瑞典圓弧法［2］，該方法被認為是極限平衡法的開端。之后，簡化Bishop法［3］、Morgenstern-Price法［4］、Sarma法［5］等極限平衡法相繼提出并廣泛應用。極限平衡法從力與力矩平衡角度進行穩定性求解，具有完善、成熟的理論體系。但極限平衡法弊端也是十分明顯的，該方法由于不能有效考慮巖土體應力應變關系，須引入大量假設條件才可解決邊坡穩定性問題，而引入的假定條件會使該方法嚴密性不足［6］。自上世紀中期以來，隨著力學理論的快速發展，基于塑性力學中關聯及非關聯流動法則衍生的極限分析法成為邊坡工程領域研究熱點。其中上限分析法從能量與速度的角度進行邊坡穩定性分析，計算過程嚴謹，可解決邊坡巖土體破壞中的運動學問題。學者Chen［7］率先將上限分析法應用至邊坡工程領域，并詳細說明了上限分析法在邊坡穩定性分析中的計算過程及程序實現。學者Donald［8］將對邊坡破壞區域內的巖土體進行斜條分劃，提出斜條分法與平動運動場相結合的上限分析法，有效擴寬了該方法的適用范圍。顯然，上限分析法在穩定性分析中具有一定優點，但由于計算的限制，該方法在含軟弱夾層邊坡穩定性分析中應用較少。

為解決該問題，本項目應用力學理論分析的方法建立含軟弱夾層邊坡破壞機構，基于最小二乘法擬合推導邊坡巖土體外力功率及內能耗散率計算公式，并結合上限分析定理以及強度折減原理，編制含軟弱夾層邊坡穩定性分析循環程序，最終達到含軟弱夾層邊坡穩定性上限分析法求解的目標。

1 邊坡穩定性分析

1.1 破壞機構的構建

由于軟弱夾層力學性質薄弱，對于含軟弱夾層邊坡，其潛在破壞模式為以軟弱夾層為底界面的剪切—順層破壞。由于含軟弱夾層邊坡破壞時部分巖土體沿弱層面滑動，因此邊坡破壞機構由轉動破壞機構和平動破壞機構兩部分共同構成。對于含軟弱夾層邊坡，巖土體黏聚力為c1，內摩擦角為?1，容重為γ1，軟弱夾層黏聚力為c2，內摩擦角為?2，容重為γ2，則該邊坡破壞機構如圖1所示。

圖1中，邊坡高度為h，坡面角為α，軟弱夾層傾角為β，破壞機構由轉動破壞機構（對數螺旋線AM）與平動破壞機構（軟弱夾層線CM）共同構成。

1.2 計算模型的建立

計算模型的建立是進行邊坡穩定性分析的重要基礎，上限分析法是從能量和速度的角度進行邊坡穩定性分析，因此分別建立邊坡幾何模型與運動學模型，并將2種模型結合進行穩定性分析。

1.2.1 幾何模型的建立

圖1中，r為轉動破壞機構與旋轉中心連線的線段距離，m；θ為轉動破壞機構與旋轉中心連線的線段水平向夾角，rad；則該邊坡轉動破壞機構方程可表示為

式中，a為方程待定參數。

以旋轉中心O點為原點，分別以水平方向以及豎直方向為坐標軸的X軸、Y軸，建立直角坐標系。則根據極坐標與直角坐標的轉換關系，式（1）可進一步表示為

顯然，式（2）是只關于x、y2個未知數的隱函數，同時該隱函數無法進行有效的顯化處理。因此，本研究采用最小二乘法的方式對該隱函數進行擬合。由于邊坡破壞機構為嚴格增函數，因此只對增區間內進行多項式擬合，令擬合結果為x=f（a，y）。

同時，在直角坐標系下，令該邊坡的坡底線方程為y=y0，則該邊坡的坡頂線方程為y=y0+h。邊坡坡面角為α，軟弱夾層傾角為β，則坡面方程f（y）以及軟弱夾層方程g（y）可分別表示為

式中，b、c分別為坡面方程及軟弱夾層方程中的待定參數。

圖1中，對數螺旋線與軟弱夾層相交點為M，令M點坐標為（xm，ym），則平面直角坐標系下直線OM方程可表示為

直線OM將邊坡破壞區域劃分為兩部分，直線OM上部（曲變形ABPM）為轉動破壞區域，該區域應用幾何公式可表示為

直線OM下部（三角形CPM）為平動破壞區域，該區域應用幾何公式可表示為

式（5）和式（6）中，yp為P點的縱坐標，可通過邊坡坡面方程f（y）與直線OM方程m（y）結合求得。

1.2.2 運動學模型的建立

曲變形ABPM為轉動破壞區域，該區域內滑體以固定角速度ω繞旋轉中心O轉動運動，則該區域內滑體速度可表示為

三角形CPM為平動破壞區域，該區域內滑體以固定線速度平動運動。邊坡破壞區域內形成多塊體的破壞模式，相鄰的塊體之間應滿足速度相容關系［9］。M點速度相容關系如圖2所示。

圖2中vm為轉動破壞區域滑體速度，v2為平動破壞區域滑體速度，v'為OM上塊體相對速度。顯然，線段OM、軟弱夾層CM均為速度間斷面。根據圖2中的速度相容關系圖，可確定平動破壞區域滑體速度v2可表示為

1.3 上限分析計算

上限分析法是從能量與速度的角度進行邊坡穩定性分析，其理論基礎為當邊坡坡體內的巖土體處于極限平衡狀態時，外力功率應該等于其內部能量耗散率。因此，下文對含軟弱夾層邊坡破壞區域巖土體外力功率以及內能耗散率分別進行計算。

1.3.1 外力功率計算方法

自然狀態下，含軟弱夾層邊坡坡體內外力功率由重力提供，因此：

式中，E為含軟弱夾層邊坡坡體內的外力功率，kW；dS為任意取出的面積微元，m2；v為該微元對應的重力方向分速度，m/s。

由于轉動破壞區域與平動破壞區域速度求解方式不同，因此分別對轉動破壞區域與平動破壞區域外力功率進行求解。

對于轉動破壞區域，面積微元在重力方向上分速度為v1cosθ，根據式（7）重力方向分速度可簡化為ωx，則轉動破壞區域外力功率E1可表示為

對于平動破壞區域，面積微元在重力方向上分速度為v2sin(β??2)，則平動破壞區域外力功率E2可表示為

該邊坡破壞區域外力功率可通過轉動破壞區域外力功率E1與平動破壞區域外力功率E2二者相加求得，即E=E1+E2。

1.3.2 內能耗散率計算方法

傳統的力學理論認為，摩擦能耗與剪脹能耗恰好可以抵消，因此含軟弱夾層邊坡破壞區域內能耗散全部由黏聚力提供，且全部集中于速度間斷面上。圖1中邊坡速度間斷面由對數螺旋線AM、軟弱夾層CM以及線段OM三部分構成，因此，分別對上述三部分內能耗散率求解。

式中，θ0為線段OA的水平方向夾角，rad；r0為線段OA的長度，m；θ1為線段OM的水平方向夾角，rad；r1為線段OM的長度，m。

根據2點之間的距離計算公式，對圖1中AM上的內能耗散率可表示為

對于軟弱夾層CM，平動破壞區域滑體在軟弱夾層方向上分速度為v2cos?2，因此圖1中軟弱夾層CM上的內能耗散率D2可表示為

對于速度間斷面OM，轉動破壞區域滑體在OM方向上分速度為0，平動破壞區域滑體在OM方向上分速度為v2sin2?2。因此OM上內能耗散率D3可表示為

則該邊坡內能耗散率可通過三部分內能耗散率相加求得，即D=D1+D2+D3。

1.3.3 強度折減運算

上文中給出了含軟弱夾層邊坡功率的計算方法，為使含軟弱夾層邊坡逐漸過渡至極限平衡狀態，同時考慮到利于計算機編程的需求，本研究結合強度折減法進行穩定性分析，其公式為

式中，F為強度折減系數。顯然，通過人工計算的手段將含軟弱夾層邊坡折減到極限平衡狀態十分復雜，本項目借助計算機編程技術實現對于含軟弱夾層邊坡穩定性分析的強度折減。

1.3.4 計算機循環編程

強度折減規律表明，當巖土體強度折減系數等于邊坡穩定系數時，含軟弱夾層邊坡破壞區域外力功率與內能耗散率二者恰好相等。同時，強度折減系數與二者比值呈現嚴格相關性。基于上述特性，本項目構建函數ε=（D/E-1）min，顯然函數ε在［0，+∞）上只有一個零點，可采用二分法進行穩定系數求解。同時為使程序可有效終止，本項目將滿足∣ε∣≤0.001近似為計算零點，并作為邊坡穩定性分析程序的終止條件。

經過邊坡穩定性分析二分法計算流程可輸出滿足∣ε∣≤0.001的折減系數。同時當D/E-1取得最小值時，根據各參數值可繪制邊坡滑裂面。該折減系數及滑裂面形態即為邊坡穩定系數及最危險滑面，因此可最終實現含軟弱夾層邊坡穩定性分析的目標。

2 算例分析

某露天礦山形成的邊坡坡高為120 m，坡面角α=25°，巖性為炭質泥巖，粘聚力c1=150 kPa，內摩擦角?1=30°，容重γ1=20 kN/m3。邊坡巖土體內發育有通過邊坡坡趾的軟弱夾層，傾角β=5°，黏聚力c2=30 kPa，內摩擦角?2=10°，容重γ2=12 kN/m3，邊坡形態及破壞機構如圖3所示。

2.1 邊坡穩定性計算

根據計算流程，選取初始值a1=0，b1=3，則在第一次的循環程序中，中值F1=1.5。經過強度折減后c1'=100 kPa，?1'=21.05°，c2=20 kPa，?2'=6.7°，則對數螺旋線AM的表達式為

函數增區間為y∈（-0.417a，-1.831a），對區間內函數表達式進行多項式擬合，擬合結果為

同時，在平面直角坐標系下，邊坡坡底線方程為y=y0，坡頂線方程為y=y0+120，邊坡坡面方程為f（y）=2.145y+b，軟弱夾層方程為g（y）=11.43y+c，直線OM方程為m（y）=xm×y/ym。則含軟弱夾層邊坡功率計算公式可分別表示為

為確定內能耗散率與外力功率的比值，應針對邊坡工程特性給出相應的約束條件。根據邊坡工程實踐，各參數應滿足以下方面約束條件：①邊坡坡面與軟弱夾層相交于邊坡坡趾C點，則約束條件為f（y0）=g（y0）；②對數螺旋線與軟弱夾層相交于點M，則約束條件為xm=g（ym）=f（a，ym）；③邊坡坡面與直線OM相交于點P，則約束條件為xp=f（yp）=f（a，yp）；④對數螺旋線位于第四象限增區間，則約束條件為 -0.417a≤ym≤-1.831a，-0.417a≤y0+120≤-1.831a；⑤邊坡破壞機構位于邊坡坡面內部，則約束條件為f（a，y）-f（y）≥0，y∈（ym，y0）。

根據內能耗散率與外力功率計算公式以及上述約束條件，可求得D/E-1的最小值ε。第一次循環中，強度折減系數F=1.5，求得ε<-0.001，不滿足∣ε∣≤0.001的程序終止要求，因此進入下一次循環。反復執行運算可最終輸出滿足∣ε∣≤0.001的折減系數，循環計算結果如表1所示。

表1表明：經過11次循環，輸出的計算結果滿足∣ε∣≤0.001要求，此時折減系數為1.195，則該邊坡穩定系數Fs=1.195。同時，根據各參數值可繪制含軟弱夾層邊坡最危險滑面，上限分析法穩定性分析結果見后文圖4。

2.2 計算結果評價

根據前文中分析可知，上限分析法和極限平衡法分析角度完全不同，但對于同一邊坡，其穩定性分析結果應是一致的。考慮到極限平衡法中簡化Janbu法對于含軟弱夾層邊坡具有良好適用性，因此選取簡化Janbu法對上限分析法計算結果進行準確性評價。簡化Janbu法穩定性分析結果如圖4所示。

從穩定系數的角度分析，兩方法計算偏差小于5%，表明兩方法都具備較高的準確性。簡化Janbu法雖然對于含軟弱夾層邊坡穩定性具有良好適用性，但由于其求解結果并非嚴格的上限解或下限解，計算誤差無法估算。而上限分析法對于含軟弱夾層邊坡穩定系數的計算結果為嚴格上限解，計算誤差可通過工程經驗進行消除［10］。

從最危險滑面形態的角度分析，兩方法最危險滑面在軟弱夾層處交點位置近似相同，表明兩方法求得的最危險滑面平動破壞機構形態近似相同。相比于簡化Janbu法，上限分析法最危險滑面起始點更靠近于邊坡后部，表明兩方法最危險滑面轉動破壞機構形態存在一定差異。上限分析法求得最危險滑面轉動破壞機構形態為對數螺旋線，而剩余推力法求得最危險滑面轉動破壞機構形態為圓弧線，根據Mohr-Coulomb屈服準則，上限分析法求得的最危險滑面形態可有效滿足巖土體材料的速度分離要求。因此從最危險滑面形態的角度分析，上限分析法計算結果具有更強的工程實踐意義。

3 結 論

（1）應用力學理論分析構建了滿足關聯流動法則及速度分離要求的含軟弱夾層邊坡破壞機構。

（2）建立了含軟弱夾層邊坡破壞機構幾何模型和運動學模型，并基于最小二乘擬合的方式推導了含軟弱夾層邊坡破壞機構外力功率及內能耗散率計算公式。

（3）綜合應用上限分析定理和強度折減原理，編制了含軟弱夾層邊坡穩定性上限分析法計算機循環程序，并通過工程算例對上限分析法計算結果進行了準確性評價，充分驗證了該方法計算結果的可靠性。