復(fù)雜渠坡穩(wěn)定性分析方法比較

岑威鈞，鄧同春，石從浩，王 蒙

(河海大學(xué) 水利水電學(xué)院，江蘇 南京 210098)

目前邊坡抗滑穩(wěn)定性分析方法主要有3類:經(jīng)典的剛體極限平衡法、有限元極限平衡法和Zienkiewicz等提出的有限元強(qiáng)度折減法.經(jīng)典的剛體極限平衡法視邊坡土體為剛體，不考慮土體本身變形對邊坡穩(wěn)定的影響，同時在土條受力平衡分析時引入一些假設(shè)和簡化，不同假設(shè)下得到的安全系數(shù)有一定的差異.但由于原理簡單，易于編程實(shí)現(xiàn)，在工程應(yīng)用中最為廣泛.有限元極限平衡法［1］在有限元應(yīng)力變形計(jì)算成果的基礎(chǔ)上對不同滑動體進(jìn)行土體的極限平衡分析，滑動體的確定可以采用與剛體極限平衡法中相同的優(yōu)化算法，也可采用一些專門的滑動軌跡搜索算法［2－3］.在計(jì)算土體應(yīng)力變形時可以采用不同的土體本構(gòu)模型，以反映土體的非線性、非彈性、剪脹(縮)性和各向異性等土體基本力學(xué)特性，故在理論上較剛體極限平衡法更為先進(jìn)合理，得到的安全系數(shù)及滑動軌跡具有明確的物理意義.目前這類方法在工程邊坡分析中也已得到了較為廣泛的應(yīng)用.O.C.Zienkiewicz等提出的有限元強(qiáng)度折減法［4－5］完全拋棄了前2種方法要不斷搜索可能破壞面的基本假定，對土體抗剪強(qiáng)度指標(biāo)進(jìn)行折減后直接進(jìn)行有限元彈塑性計(jì)算，最終能夠得到一個土體恰好“破壞”時的折減系數(shù).該折減系數(shù)具有明確的安全系數(shù)物理意義.判斷邊坡失穩(wěn)“破壞”的依據(jù)主要有塑性區(qū)是否貫通、迭代計(jì)算是否收斂及特征部位變形是否發(fā)生突變等判據(jù)［6］.邊坡彈塑性計(jì)算時的屈服準(zhǔn)則以Mohr-Coulomb(M-C)準(zhǔn)則和Drucker-Prager(D-P)系列準(zhǔn)則應(yīng)用最為廣泛.有限元強(qiáng)度折減法既可以用于土質(zhì)邊坡穩(wěn)定分析，也可以用于巖質(zhì)邊坡穩(wěn)定分析.一般含M-C準(zhǔn)則或D-P準(zhǔn)則的基本理想彈性有限元程序就可直接使用，不需要額外增加搜索破壞面的程序代碼，因此較前2種方法更為方便.

本文以南水北調(diào)工程河南境內(nèi)某典型復(fù)雜渠道邊坡為例，采用剛體極限平衡法中的Spencer法［7］、有限元極限平衡法以及有限元強(qiáng)度折減法進(jìn)行比較計(jì)算，比較這3種方法以及不同屈服準(zhǔn)則對渠坡穩(wěn)定性分析結(jié)果的影響.

1 3種穩(wěn)定計(jì)算方法的基本理論

1.1 剛體極限平衡法

剛體極限平衡法將可能滑動體視為剛體，分別計(jì)算滑動體所受的抗滑力(矩)和滑動力(矩)，以抗滑力(矩)與滑動力(矩)之比作為衡量邊坡穩(wěn)定的依據(jù)，其基本公式如下:

式中:τf為土體的抗剪強(qiáng)度;τ為剪應(yīng)力.

根據(jù)可能滑動體各土條是否考慮條間力、條間力的作用位置以及滑弧形狀等假定的不同，剛體極限平衡法又分為瑞典圓弧法、Bishop法、Morgenstern-Price法、Spencer法等.其中，Spencer法引入了條間法向力與切向力夾角恒定以及切向力合力作用點(diǎn)位于土條底面的假設(shè)，使得各土條能夠滿足全部靜力學(xué)平衡條件，最為嚴(yán)格，故本文選擇Spencer法作為剛體極限平衡法的典型計(jì)算方法.

1.2 有限元極限平衡法

1.2.1 滑動面上的應(yīng)力 有限元極限平衡法中，假定滑動面上某點(diǎn)的應(yīng)力狀態(tài)為σx，σy和τxy，則該點(diǎn)的法向應(yīng)力 σn和切向應(yīng)力 τ［8－9］為:

式中:α為滑動面與水平面的夾角.

又，該點(diǎn)的抗剪強(qiáng)度為:

式中:c為黏聚力;φ為內(nèi)摩擦角.則整個滑動面的安全系數(shù)為:

為了計(jì)算滑弧面上各點(diǎn)的應(yīng)力，首先要確定滑弧附近單元結(jié)點(diǎn)的應(yīng)力狀態(tài).結(jié)點(diǎn)處應(yīng)力需通過單元高斯點(diǎn)處應(yīng)力映射(外推)而來，其計(jì)算公式［10］如下:

式中:f為單元結(jié)點(diǎn)應(yīng)力;N為映射矩陣;F為高斯點(diǎn)的應(yīng)力值.

對于被滑弧穿過的每個單元都執(zhí)行以上映射過程，通過繞結(jié)點(diǎn)平均的方法求得滑弧附近每一結(jié)點(diǎn)處的σx，σy和τxy，則按式(2)即可求得滑弧上相應(yīng)點(diǎn)的法向和切向應(yīng)力.

1.2.2 滑動面搜索方法 如果滑動面的形式已知，如常見的圓弧或復(fù)式滑動面等，則最危險滑動面的確定可以采用與剛體極限平衡法完全相同的搜索方法.但是對于滑動面形式未知的情況下，則需要采用一定的辦法搜索出最終連貫的滑動破壞面.S.K.Giam等提出的方法可供選用［11］，該方法通過有限元計(jì)算找出應(yīng)力水平較高點(diǎn)(記作P點(diǎn))，由P點(diǎn)出發(fā)在相鄰兩土條內(nèi)構(gòu)造局部滑裂面.圖1給出了由P點(diǎn)在第N+1個土條內(nèi)構(gòu)造局部滑裂面的過程，即以P點(diǎn)為基點(diǎn)呈放射狀按照1°傾角變化形成m個底邊長分別為ΔL的微小段，每一小段以其中點(diǎn)為控制點(diǎn)，控制點(diǎn)的應(yīng)力狀態(tài)由有限元計(jì)算結(jié)果得出，則由P點(diǎn)出發(fā)所得第i個條塊的安全系數(shù)為:

圖1 滑動面搜索示意Fig.1 Sketch of sliding surface search method

找出Fi中的最小值Fimin(i=1，2，…，m).Fimin即處于第N+1個土條所處局部滑裂面位置.再以Q點(diǎn)為基點(diǎn)找到第N+2個土條內(nèi)的局部滑裂面位置，重復(fù)上述過程，依次找到滑動面直至與坡面相交，即可形成整個邊坡的滑動面.對于復(fù)雜邊坡，建議多選擇幾個P點(diǎn)，有利于快速找到真正的滑動破壞面［12］.

1.3 有限元強(qiáng)度折減法

1.3.1 折減系數(shù) 有限元強(qiáng)度折減法是將巖土體材料的抗剪強(qiáng)度除以折減系數(shù)后用于有限元彈塑性計(jì)算，取邊坡達(dá)到破壞時的折減系數(shù)為安全系數(shù).假定巖土體的強(qiáng)度參數(shù)為c和φ，試算所取折減系數(shù)為Fr，則折減強(qiáng)度τr為:

則巖土體材料的折減強(qiáng)度參數(shù)cr，φr的計(jì)算表達(dá)式為:

1.3.2 屈服準(zhǔn)則

(1)Mohr-Coulomb屈服準(zhǔn)則 M-C準(zhǔn)則是邊坡穩(wěn)定性分析中廣泛使用的屈服準(zhǔn)則［12］，其表達(dá)式為:

式中:I1為應(yīng)力張量的第一不變量，滿足I1=σii;J2為偏應(yīng)力張量的第二不變量，滿足為應(yīng)力羅德角.

圖2 D-P準(zhǔn)則與M-C準(zhǔn)則的擬合關(guān)系(偏平面內(nèi))Fig.2 Fitting relationships of D-P criteria and M-C criterion(in deviatoric plane)

M-C準(zhǔn)則在應(yīng)力空間中是一個以靜水應(yīng)力軸為對稱軸的不等邊六角錐，存在棱角，給數(shù)值計(jì)算帶來不便［13］，故可用與M-C準(zhǔn)則相匹配的D-P準(zhǔn)則替代.

(2)與M-C相匹配的Drucker-Prager系列準(zhǔn)則 D-P準(zhǔn)則的通用表達(dá)式如下:

根據(jù)與M-C準(zhǔn)則不同的匹配條件，可以產(chǎn)生D-P系列準(zhǔn)則，記為D-P1～D-P5，偏平面內(nèi)的擬合關(guān)系見圖2，其匹配條件及系數(shù)見表1.

表1 與M-C相匹配的D-P系列準(zhǔn)則系數(shù)［13］Tab.1 Coefficients of the M-C criterion matching D-P series criteria

由式(10)以及表1可知:同一偏平面內(nèi)，系數(shù)α越大，屈服時所需的第二偏應(yīng)力不變量J2越大，故同一邊坡使用不同的D-P系列準(zhǔn)則所得安全系數(shù)滿足D-P1＞D-P2＞D-P5＞D-P4＞D-P3.

2 渠坡穩(wěn)定性比較計(jì)算

計(jì)算對象為南水北調(diào)工程中線河南境內(nèi)某典型復(fù)雜渠道邊坡，地層以相互間隔的泥灰?guī)r和弱膨脹性黏土巖為主，剖面示意圖見圖3.一級馬道寬4.66 m，以下設(shè)10 cm厚的混凝土襯砌板，坡比1∶2;二三級馬道均寬2.00 m，坡比1∶1.75，坡頂四級馬道相鄰上下段坡比均為1∶1.5，渠底寬9.66 m.渠道所處地層巖土材料的主要物理力學(xué)參數(shù)見表2.

取渠中心線一側(cè)渠坡作為計(jì)算對象，計(jì)算域范圍參見圖3.分別采用剛體極限平衡法中的Spencer法，彈塑性有限元極限平衡法以及有限元強(qiáng)度折減法進(jìn)行計(jì)算，屈服準(zhǔn)則取為M-C準(zhǔn)則.圖4為3種計(jì)算方法所得渠坡滑弧示意圖，其中強(qiáng)度折減法以等效塑性應(yīng)變貫通作為破壞標(biāo)準(zhǔn).表3給出了3種穩(wěn)定分析方法得到的渠坡抗滑穩(wěn)定安全系數(shù)，由表可知3種計(jì)算方法所得安全系數(shù)十分接近，安全系數(shù)相差在3%之內(nèi)，且所得滑弧位置基本一致.

圖3 渠坡典型剖面地層分布(單位:m)Fig.3 A section of a typical channel slope profile(unit:m)

表2 巖土體的物理力學(xué)參數(shù)Tab.2 Physical and mechanical parameters of the rock materials

圖4 3種計(jì)算方法所得渠坡滑弧示意Fig.4 Sketches of channel slope sliding surfaces obtained from three different methods

表3 不同屈服準(zhǔn)則所得安全系數(shù)Tab.3 Safety factors given by different yield criteria

此外，還根據(jù)不同形式的D-P系列屈服準(zhǔn)則計(jì)算了該渠坡的安全系數(shù).相對于M-C準(zhǔn)則，采用D-P準(zhǔn)則后抗剪強(qiáng)度參數(shù)折減過程中其計(jì)算結(jié)果具有很好的單調(diào)收斂性，不像采用M-C準(zhǔn)則時，對某些復(fù)雜邊坡在強(qiáng)度參數(shù)折減的過程中其計(jì)算結(jié)果有時會出現(xiàn)波動，不易趨于穩(wěn)定.這應(yīng)該與不等邊六角錐形的M-C屈服準(zhǔn)則存在棱角易引起數(shù)值計(jì)算奇異性有關(guān).限于篇幅，本文僅給出各類D-P系列屈服準(zhǔn)則所得安全系數(shù)(表3)，其中不同類型D-P系列屈服準(zhǔn)則下的屈服破壞面與圖4(c)很相近.由表3可知，不同形式的D-P準(zhǔn)則所得安全系數(shù)偏差較大，且D-P1所得安全系數(shù)最大;D-P3與D-P4所得安全系數(shù)與M-C準(zhǔn)則所得安全系數(shù)較為接近;若以M-C準(zhǔn)則為基準(zhǔn)，D-P1與M-C安全系數(shù)相差22.4%;不同的屈服準(zhǔn)則計(jì)算出的規(guī)律與1.3.2節(jié)中所述規(guī)律一致.事實(shí)上，表1給出的5種不同形式的D-P準(zhǔn)則之間有一定的相關(guān)性，各自屈服準(zhǔn)則下得到的不同安全系數(shù)之間可以相互轉(zhuǎn)換，其結(jié)果是一致的［14］.

3 結(jié)語

(1)基于有限元分析的極限平衡法以及強(qiáng)度折減法所得安全系數(shù)與剛體極限平衡法所得結(jié)果十分接近，且滑弧位置也基本一致，說明3種方法計(jì)算邊坡穩(wěn)定具有足夠的可靠性.

(2)鑒于M-C準(zhǔn)則的屈服面在應(yīng)力空間中存在棱角，為了使數(shù)值計(jì)算過程易于穩(wěn)定，建議考慮使用D-P準(zhǔn)則代替M-C準(zhǔn)則進(jìn)行計(jì)算.

(3)采用強(qiáng)度折減法計(jì)算邊坡穩(wěn)定時，屈服準(zhǔn)則的具體形式對計(jì)算結(jié)果影響很大.D-P系列準(zhǔn)則中D-P1(即M-C準(zhǔn)則的外接圓)與M-C準(zhǔn)則相差達(dá)到20%以上，不能直接作為安全系數(shù)計(jì)算成果(各準(zhǔn)則下的安全系數(shù)值可以進(jìn)行轉(zhuǎn)換).D-P3準(zhǔn)則以及D-P4準(zhǔn)則得到的結(jié)果較為精確，其強(qiáng)度折減系數(shù)可以直接作為邊坡穩(wěn)定的安全系數(shù)選用.

因此，在渠坡穩(wěn)定性計(jì)算時選擇何種計(jì)算方法對渠坡安全系數(shù)影響不大，使用有限元強(qiáng)度折減法計(jì)算時建議優(yōu)先采用D-P3或D-P4代替M-C強(qiáng)度準(zhǔn)則.

［1］DUNCAN J M.State of the art:Limit equilibrium and finite-element analysis of slopes［J］.Journal of Geotechnical Engineering，1996，122(7):577-596.

［2］曾亞武，田偉明.邊坡穩(wěn)定性分析的有限元法與極限平衡法的結(jié)合［J］.巖石力學(xué)與工程學(xué)報(bào)，2005，24(增刊2):5355-5359.(ZENG Ya-wu，TIAN Wei-ming.Slope stability analysis by combining FEM with limit equilibrium method［J］.Chinese Journal of Rock Mechanics and Engineering，2005，24(Suppl2):5355-5359.(in Chinese))

［3］劉真真，方朝陽，朱濤，等.有限元極限平衡法在邊坡穩(wěn)定性分析中的應(yīng)用［J］.水電能源科學(xué)，2010，28(12):93-94.(LIU Zhen-zhen，F(xiàn)ANG Chao-yang，ZHU Tao，et al.Application of FEM limit equilibrium method to slope stability analysis［J］.Water Resources and Power，2010，28(12):93-94.(in Chinese))

［4］ZIENKIEWICZ O C，HUMPHESON C，LEWIS R W.Associated and non-associated visco-plasticity and plasticity in soil mechanics［J］.Geotechnique，1975，25(4):671-689.

［5］GIFFITHS D V，LANE P A.Slope stability analysis by finite elements［J］.Geotechnique，1999，49(3):387-403.

［6］趙尚毅，鄭穎人，張玉芳.有限元強(qiáng)度折減法中邊坡失穩(wěn)的判據(jù)探討［J］.巖土力學(xué)，2005，26(2):332-336.(ZHAO Shang-yi，ZHENG Ying-ren，ZHANG Yu-fang.Study on slope failure criterion in strength reduction finite element method［J］.Rock and Soil Mechanics，2005，26(2):332-336.(in Chinese))

［7］鄭榮躍，梧松.基于Spencer法的邊坡穩(wěn)定性可靠度指標(biāo)分析［J］.巖土力學(xué)，2006，27(1):147-153.(ZHENG Rongyue，WU Song.Reliability index analysis of slope stability based on Spencer's method［J］.Rock and Soil Mechanics，2006，27(1):147-153.(in Chinese))

［8］殷宗澤.土工原理［M］.北京:中國水利水電出版社，2007:426-427.(YIN Zong-ze.Geotechnical principles［M］.Beijing:China Water＆ Power Press，2007:426-427.(in Chinese))

［9］王連捷，周春景，吳珍漢，等.用有限元極限平衡法分析邊坡的穩(wěn)定性［J］.中國地質(zhì)，2008，35(5):1031-1034.(WANG Lian-jie，ZHOU Chun-jing，WU Zhen-han，et al.Slope stability analysis by the limit equilibrium method based on the finite element method［J］.Geology in China，2008，35(5):1031-1034.(in Chinese))

［10］ZIENKIEWICZ O C，TAYLOR R L.The finite element method for solid and structural mechanics［M］.Oxford:Butterworth-Heinemann，2000:140-143.

［11］GIAM S K，DONALD I B.Determination of critical slip surfaces for slopes via stress-strain calculation［C］∥5thA N Z Conf Geomechanics，Sydney，1988:461-464.

［12］張學(xué)言，閆澍旺.巖土塑性力學(xué)基礎(chǔ)［M］.天津:天津大學(xué)出版社，2006:81-82.(ZHANG Xue-yan，YAN Shu-wang.Fundamentals of plastic geotechnics［M］.Tianjin:Tianjin University Press，2006:81-82.(in Chinese))

［13］鄧楚鍵，何國杰，鄭穎人.基于M-C準(zhǔn)則的D-P系列準(zhǔn)則在巖土工程中的應(yīng)用研究［J］.巖土工程學(xué)報(bào)，2006，28(6):735-738.(DENG Chu-jian，HE Guo-jie，ZHENG Ying-ren.Studies on Drucker-Prager yield criterions based on M-C yield criterion and application in geotechnical engineering［J］.Chinese Journal of Geotechnical Engineering，2006，28(6):735-738.(in Chinese))

［14］鄭穎人，趙尚毅，鄧楚鍵，等.有限元極限分析法發(fā)展及其在巖土工程中的應(yīng)用［J］.中國工程科學(xué)，2006，8(12):39-59.(ZHENG Ying-ren，ZHAO Shang-yi，DENG Chu-jian，et al.Development of finite element limit analysis method and its applications in geotechnical engineering［J］.Engineering Science，2006，8(12):39-59.(in Chinese))