基于彈塑性擴展有限元的裂隙邊坡穩定性分析

鄧 幫 李 旋

(武漢市建筑工程質量監督站,湖北 武漢 430061)

1 概述

在分析邊坡穩定的時候，邊坡頂端的張拉裂紋通常被視作邊坡失穩的第一個信號，坡頂有無張拉裂紋也被視為邊坡穩定與否的一個很重要的證據。所以對于頂部含有裂紋的邊坡的分析在邊坡穩定分析中占有很重要的地位[1]。目前，邊坡穩定分析最常用的方法是極限平衡法。極限平衡法的優點在于能夠考慮復雜的地層、滲水和不同的加載情況；不足之處在于極限平衡法的前提是土體可以劃分成條，而這必然要提到相鄰條之間的作用力且滑移體內部一般假設是剛性的，這與實際情況不符。而將有限元的理念應用到邊坡穩定性分析相較于傳統的極限平衡法有諸多優勢[2]，如能夠對具有復雜地貌的邊坡進行計算、能考慮土體的非線性彈塑性本構關系以及變形對應力的影響、能夠模擬土坡的失穩過程及滑移面的形狀、能夠模擬土體與支護的共同作用、求解安全系數的時候，無需假設滑動面的形狀，不必進行條分等等。但是對于含有張拉裂紋的邊坡穩定性分析而言，如何模擬張拉裂紋在有限元中是一個難題，首先是必須不停地進行網格重繪以保證裂紋走向和單元邊界一致，其次是在裂紋尖端由于存在奇異性，為了得到較為準確的計算結果，必須進行大量的網格細分操作。擴展有限元(Extended Finite Element Method,XFEM)在面對這類含有裂紋的問題時，比傳統有限元更有優勢[3,4]。

強度折減法(Shear Strength Reduction，SSR)[5,6]已被多位研究者證明不僅可以用來計算邊坡的安全系數，更可以可視化地展示邊坡破壞的發展過程。本文主要研究將彈塑性擴展有限元法和強度折減法結合考查含有張拉裂紋的邊坡穩定性。

2 彈塑性擴展有限元和強度折減法的結合

2.1 彈塑性擴展有限元法

擴展有限元作為一種基于單元分解的方法，在有限元的近似模擬空間中加入了描述域函數局部特征的特殊函數。

(1)

其中，uI為節點上的常規自由度對應的未知量；aI為Heaviside函數H(x)對應的富集自由度未知量；bI為漸進的裂尖函數Φα(x)對應的富集自由度未知量；N為所有節點組成的集合；Ncr為形函數支持域被裂紋穿過的節點組成的集合；Ntip為形函數支持域包括裂尖的節點組成的集合。

對于線彈性體裂紋尖端附近所有點的位移，可以用下面這個空間中的基線性表達：

(2)

可以看出，用上面這組基表達出的函數，在裂紋尖端附近r=0，對r求偏導則會產生奇異，所以上面這組基只適合表達脆性材料裂紋尖端的位移，對于非脆性材料，則并不適合。本文采用了一種與裂紋模型無關的富集基[7]:

Φ(x)=[r1.02θ,r1.20θ,r1.57θ,r1.95θ]

(3)

2.2 強度折減法

強度折減法利用抗剪強度折減系數的概念，認為在極限情況下，外荷載所產生的實際剪應力與抵抗外荷載所發揮的最低抗剪強度(即按照實際抗剪強度折減后所確定的、實際中得以發揮的抗剪強度)相等。可見對于Mohr-Coulomb模型而言，存在下列關系：

τ=C+σntanφ

(4)

τf=Cf+σntanφf

(5)

其中，τ為坡體材料的剪切強度；τf為滑動面上的剪切強度。粘聚力Cf和內摩擦角φf可用式(6)，式(7)求得：

(6)

(7)

傳統的邊坡穩定極限平衡法中，安全系數定義為沿滑動面的抗剪強度和滑動面上實際剪力的比：

(8)

將式(8)兩邊同除以折減系數ω，則式(8)變為：

(9)

可見當強度折減ω以后，坡體達到極限狀態。擴展有限元強度折減法與傳統的邊坡穩定極限平衡法在本質上是一致的。

本文將XFEM和SSR相結合對邊坡穩定性進行分析的步驟如圖1所示，簡述如下：

1)建立模型，利用水平集方法描述張拉裂紋，施加自重，將初始的安全系數設置為F0=1.0。

2)根據式(6)和式(7)降低剪切強度。

3)求出由于自重造成的應力場、應變場、位移場(這里不考慮外部荷載如堆載之類的影響)。

4)計算內力。如果不平衡力小于容許值，返回第3)步，再次進入Newton-Raphson迭代，否則進入第5)步。

5)如果Newton-Raphson迭代收斂，用Fk=Fk-1+ΔFi更新安全系數，然后進入第2)步，否則得到最終的安全系數，程序結束。

3 實際工程裂隙邊坡穩定分析

3.1 工程概況

谷竹高速公路(谷城—竹溪)是湖北省規劃的“552(五縱五橫二環)”骨架公路網的重要組成路段，所經過的ZK172+950～ZK173+420東川滑坡群位于湖北省竹山縣溢水鎮東川村境內，路線從東川1號、2號滑坡前緣通過。滑坡前緣為河流階地。坡面局部為梯田和池塘，東川滑坡沿竹山斷裂方向一字展開，滑坡體的后緣可見較明顯的斷裂陡坎發育，相對高差1 m～3.5 m；滑坡體前部地形稍緩，坡度15°～20°。

以東川1號滑坡為例，滑坡形態：1號滑坡平面上呈梯形，縱長約180 m，前緣寬約170 m，后緣寬約80 m。后緣高程407.5 m左右，前緣高程370.9 m左右，相對高差約36 m；2號滑坡空間上呈圈椅狀，縱長約205 m，前緣寬約290 m，平均寬220 m。后緣高程420.5 m左右，前緣高程359 m左右，相對高差約61 m。根據現場地質調繪結合區域地質資料，東川滑坡群大地構造部位處于秦嶺褶皺系南沿。東川滑坡的滑體為殘坡積的粉質粘土夾碎片石及上部呈砂礫狀的強風化含炭鈣質板巖。從揭露的幾個鉆孔可見滑動帶位于強風化層偏上部,成分為強風化含炭鈣質板巖。滑床為下伏的基巖，巖性為碎塊狀的強風化含炭鈣質板巖及中風化含炭鈣質板巖。

3.2 擴展有限元分析

1)初始情況下，認為坡體不含裂紋且不開挖，用常規有限元法結合強度折減進行分析，通過計算得到邊坡安全系數為1.18。坡體達到極限狀態時的塑性區變形情況、水平位移分別如圖2,圖3所示。從圖中可以看出，當坡角處未開挖時，坡體較為穩定，強度折減后的潛在滑動面的后緣已經位于坡腳后方60 m左右。如果坡腳處沒有挖方工程，此邊坡無需支護，可以自穩。

2)根據實地勘測，當坡腳處高速公路施工產生挖方的時候，在標高26.5 m處產生張拉裂紋，張拉裂紋深度由式(10)確定：

(10)

其中，c,φ,Y分別為上覆粉質粘土的內聚力、內摩擦角和容重。計算得安全系數0.98。此時塑性區變形情況、剪應力分布、水平位移分別如圖4，圖5所示。從圖中可以看出，坡腳處挖方后，將在坡腳和坡體張拉裂紋之間形成滑動面，并在坡腳處產生一定程度的應力應變集中。將數值模擬結果與實地監測結果比較會發現，模擬結果比實際監測位移小，這是由于本文程序暫時沒有考慮孔隙水壓力的作用。

4 結語

通過對于擴展有限元方法中能夠完全描述跳躍富集函數和尖端富集函數兩種富集函數的富集基的選擇給出了一種與裂紋體材料特性無關的富集基選擇方式，以此為基礎在擴展有限元方法中實現了Mohr-Coulomb模型。然后結合強度折減法，實現了對含裂隙邊坡的擴展有限元分析，并結合谷竹高速公路東川滑坡段進行了分析。

結果表明，將擴展有限元法和強度折減法相結合，通過對于富集基的選擇來計入Mohr-Coulomb模型的影響，以實現對于邊坡的穩定性分析，是可行的。

[1] 趙學勐,陳運理.考慮垂直裂隙影響的均質黃土挖方邊坡穩定性[J].土木工程學報,1981(1):43-49.

[2] 鄭穎人,趙尚毅,時衛民,等.邊坡穩定分析的一些進展[J].地下空間,2001,21(4):262-271.

[3] 常建梅,宋思紋.基于擴展有限元法的邊坡穩定分析[J].圖學學報,2017,38(1):128-131.

[4] Sukumar N,Chopp D L,Moran B.Extended finite element method and fast marching method for three-dimensional fatigue crack propagation[J].Engineering Fracture Mechanics,2003(70):29-48.

[5] 方修君,金 峰,王進廷.用擴展有限元方法模擬混凝土的復合型開裂過程[J].工程力學,2007,24(1):46-52.

[6] Griffiths D V,Lane P A.Slope Stability Analysis by Finite Elements[J].Geotechnique,1999,49(3):387-403.

[7] Abbas S,Peter F T.Model-independent approaches for the XFEM in fracture mechanics.IOP Conference Series:Materials Science and Engineering，2010(10):12-45.