基于ANSYS 的巖質邊坡軟弱結構面非線性接觸分析

蘇 向 震， 張 程 宏， 李 銳

(1.中國市政工程西南設計研究總院，四川成都 610081;2.長江巖土工程總公司，湖北武漢 430010)

1 引言

在巖質邊坡中存在各種形式的軟弱結構面，如裂隙、斷層、軟弱夾層等。巖質邊坡的穩定性往往受控于軟弱結構面的變形和受力特征，因此，軟弱結構面成為巖質邊坡穩定性分析的熱點。目前，處理軟弱結構面問題的有限元方法主要有兩種，即連續體法和接觸單元法。連續體法實質上是強度問題，其優點是方法簡單、力量成熟、計算耗時短，但是其不能很好的模擬軟弱結構面的張開和滑移特征，精度不高;接觸單元法是用界面的正應力和剪應力以及界面位移約束來描述界面的粘結開合狀態，其特點是接觸界面的區域大小和相互位置以及接觸狀態事先未知，隨時間變化需要在求解過程中確定，且接觸條件屬單邊性，表現出強烈的非線性，能夠很好的模擬結構面的變形和受力特性。

筆者采用增廣拉格朗日算法，借助有限元分析軟件ANSYS，建立了某巖質邊坡軟弱結構面接觸模型，考慮到軟弱結構面的工作性態對邊坡穩定的影響，利用非線性有限元技術為邊坡穩定性分析提供參考。

2 接觸分析的基本理論

接觸問題是一種高度的非線性問題，主要包括以下三種非線性:接觸邊界條件的非線性，即接觸區的范圍一般事先未知，并在接觸過程中可能還會變化;材料非線性，表現在應力集中產生的局部塑性變形;幾何非線性，在應力集中情況下往往產生局部大變形。因此，對未知接觸面和接觸作用力以及接觸體之間摩擦力的模擬，無論從數學理論上，還是計算實施中都比經典線性結構力學問題復雜。

求解接觸問題的數值方法主要有Lagrange法、罰函數法以及增廣Lagrange法。將接觸邊界條件用Lagrange乘子相乘，并與總勢能一起構成修正的勢能，再求駐值以得到最后的控制方程，該方法即為Lagrange法。Lagrange法可以很精確的滿足約束條件。但由于Lagrange乘子的引入，增加了自由度，系統的求解規模增大了，而且在系統的控制方程中出現了零主元，給計算帶來了麻煩。罰函數法是求解接觸問題的典型方法之一，它將接觸區域的非嵌入條件以及其他條件作為懲罰項引進接觸系統的能量泛函中，使原來的條件約束變分問題轉化為罰優化問題，實際上就是將接觸非線性問題轉化為材料非線性問題。罰函數法的最大優點在于引入接觸條件時并不增加系統的自由度，而且很容易從物理意義上解釋。缺點在于罰函數往往導致方程病態，隨著罰值增加，病態減弱，而約束條件只有在罰值很大時才能精確滿足。由于Lagrange法和罰函數法各有優缺點，于是將兩者聯合使用，產生了增廣Lagrange法，通過一個罰函數強迫滿足一個特定的關鍵約束。增廣Lagrange法已成功應用于不可壓縮彈性有限變形、無摩擦接觸問題及彈塑粘性問題中。

設兩個接觸物體的張開函數為gN，初始張開為g0≥0，其接觸條件為:

式中 u為接觸邊界的位移向量;n為接觸邊界的法向應力;σ為接觸體的應力;pN(u)為接觸邊界的法向應力。由coulomb摩擦定律及kuhntucker條件:

可以得到弱形式虛功方程:

在 S邊界上，如果 gn=0，則 δu·n≤0。

式中， PT(u)為接觸邊界的切向應力;uT為位移向量的切向分量;f，t分別為施加于接觸體上的體力和面力;υ，S，S0則分別表示積分域是整個接觸體、接觸面及面力作用邊界。

由于式(3)中的約束難以處理和求解，于是引入獨立的Lagrange乘子λN，λT和罰規劃，則式(3)可表示為:

式中 ξ為罰因子;Ф=‖λT‖－μλN，μ為摩擦系數;△N，△T為罰參數。

式(4)即為摩擦接觸問題的增廣Lagrange方程。

3 軟弱結構面——巖石相互作用非線性分析

相互作用問題存在兩種非線性:一是由于材料非彈性引起的材料非線性;一種是由于介質之間產生的局部滑移、脫離而造成的狀態非線性。

(1)軟弱結構面材料的模擬。

ANSYS提供了一種以Drucker-Prager準則為屈服準則的材料，屈服面不隨材料的逐漸屈服而改變，即沒有強化準則;塑性行為假定為理想彈塑性，采用相關聯流動法則。文中所述的軟弱結構面材料模型采用DP材料，為確定破壞面，共需輸入3個參數，即粘聚力C，內摩擦角φ，膨脹角φf。為保守起見，文中取φf=0。

(2)軟弱結構面與巖石接觸面上的狀態非線性模擬。

大多數情況下，巖質邊坡的失穩破壞主要表現為結構面切割形成的巖塊沿軟弱結構面的滑移，巖塊與基巖內部應力基本上處于彈性和彈塑性狀態，因此，邊坡塊體穩定問題主要反映在軟弱結構面上。一般而言，軟弱結構面在承受一定外荷的情況下就有可能逐步發展成為滑動面。在邊坡非線性穩定分析中，對滑動面的模擬即為狀態非線性模擬。文中利用ANSYS程序的接觸單元實現接觸分析。軟弱結構面作為接觸面，巖石表面剛度相對要大一些，將其作為目標面，在軟弱結構面表面形成接觸單元，在巖石表面上形成目標單元，通過相同的實常數將相對應的接觸單元和目標單元定義為一個接觸對，通過選擇合理的參數，實現軟弱結構面與巖石界面上的粘結、滑移、張開、再閉合狀態的模擬。

4 算例分析

4.1 計算模型

某巖質邊坡坡高H=25 m，邊坡淺表層為強風化Ⅳ類巖體，深部為弱風化Ⅲ類巖體，淺表層巖體中發育有兩組結構面，結構面摩擦系數為0.4，邊坡截面示意圖見圖1，有限元計算模型如圖2所示，各巖層的計算參數如表1所示。

圖1 邊坡截面示意圖

圖2 有限元計算模型示意圖

表1 計算參數表

4.2 計算結果分析

(1)邊坡應力狀態及特征

邊坡大、小主應力等值線圖及主應力矢量圖見圖3、4和5。從圖中可以看出，在自重作用下，受斜坡地形、巖體質量類別以及斷層巖脈軟弱結構面的影響，邊坡天然應力場呈現明顯的分區、分帶特征，即淺表地層的低應力帶、中部的應力過渡區、深部巖體的應力平穩區以及斷層巖脈附近的低應力條帶。坡體大主應力等值線與坡面近于平行，且隨水平和垂直埋深的增大，量值呈現增大的趨勢;小主應力基本呈現順坡分布的特點，且隨埋深的增大而增大;同時，無論σ1還是σ3，由于受斷層的影響，在覆蓋層、斷層及巖脈附近的等值線均不平穩，即出現一定的拉應力區域，拉應力最大值約 0.14 MPa。

圖3 大主應力等值線圖(單位:Pa)

圖4 小主應力等值線圖(單位:Pa)

圖5 主應力矢量圖

(2)軟弱結構面物理力學特征。

邊坡在軟弱結構面上存在的滑動趨勢使得作為接觸面的結構面上產生摩擦應力，如圖6所示。從圖6中可以看出，兩條軟弱結構面交匯處的摩擦應力較大，破碎帶坡肩處摩擦應力幾乎為0，坡肩部位處于變形拉裂區，軟弱結構面上的有效應力減小，故摩擦應力較小。

圖6 軟弱結構面摩擦應力圖(單位:Pa)

軟弱結構面上的接觸壓力如圖7所示，接觸壓力分布圖與摩擦應力分布圖相似。一般接觸壓力為0的位置，摩擦應力也為0，坡肩部位接觸壓力亦為0，說明該部位巖體有脫離母巖滑動的趨勢，需對其進行加固處理。

圖7 軟弱結構面接觸壓力圖(單位:N)

5 結語

筆者采用增廣拉擱朗日算法，運用ANSYS有限元分析軟件及其自帶的接觸單元模擬巖質邊坡軟弱結構面的粘結、滑移及受力特性，通過對接觸部位的狀態和應力進行研究，得出了接觸部位的應力分布規律，合理反映出基巖與軟弱結構面之間相互作用過程中的變形協調和應力傳遞關系，為邊坡穩定性分析提供了可靠的依據。

[1]孫林松，王德信，謝能剛.接觸問題有限元分析方法綜述[J].水利水電科技進展，2001，21(3):18－20.

[2]許有飛，何江達，梁照江.混合式閘墻與基巖間的非線性接觸分析[J].紅水河，2004，23(4):71－75.