裂隙巖體的三維非線性有限元分析

摘要：討論天然裂隙巖體的機械力學行為及力學問題，建立了相應的有限元格式，介紹了裂隙巖體剛度方程的建立及彈塑性狀態(tài)下的應力——應變關(guān)系是非線性的表達式，闡述了非線性問題的分類及解法。

關(guān)鍵詞：裂隙巖體；接觸單元；非線性問題

工程巖體經(jīng)受過復雜的地質(zhì)作用后，其內(nèi)部存在不同規(guī)模、性質(zhì)各異的地質(zhì)結(jié)構(gòu)面（節(jié)理、裂隙、斷層等），使得巖體的結(jié)構(gòu)特征與力學性質(zhì)錯綜復雜。一方面，不連續(xù)面的存在使得巖體介質(zhì)類型不同于其他任何一種力學材料。在進行裂隙巖體力學問題分析時，準確選擇巖體力學模型成為必須解決的首要問題。另一方面，巖體中的不連續(xù)面使得巖體性質(zhì)表現(xiàn)出不連續(xù)、非均勻、各向異性及尺寸效應等特征。這就使得傳統(tǒng)彈性有限元法不能更好地模擬真實情況。

近年來，節(jié)理單元被廣泛應用于模擬巖體中的斷層、剪切帶和不連續(xù)區(qū)域。在節(jié)理單元中，由于引入了巖石燕子節(jié)理的法向和切向的剛度，使得節(jié)理單元的應用受到了一定的限制。這是因為剛度系數(shù)通常都是從野外實測得到的，測試存在困難和不確定性。Katona于1983年提出了一種不用剛度系數(shù)的簡單接觸摩擦單元，它能模擬兩物體之間的滑動摩擦、張開和閉合過程，適用于不計厚度影響的節(jié)理、斷層和位移不連續(xù)場的有限元分析。由于這種單元采用的是常接觸力二節(jié)點簡單單元，因而很難適應接觸面復雜的問題。同時，由于該單元選取節(jié)點接觸力作為基本未知量，接觸應力是由接觸力平均得到的，這不僅降低了計算接觸應力的精度，而且使得這種簡單求應力的方法很難推廣至三維問題。雷曉燕于1994年在此基礎上提出了一種新的接觸摩擦單元。在這種單元中，直接選取節(jié)點接觸應力作為基本未知量，同時采用六節(jié)點的等參單元模擬接觸面復雜的幾何形狀。接觸單元的幾何和靜力約束是作為補充方程包含在剛度方程中的。運用虛位移原理導出接觸摩擦問題的等效單元剛度——約束矩陣和等效載荷向量。

1、非線性問題的分類

引起結(jié)構(gòu)產(chǎn)生非線性變形的原因很多，大致可以分為三類：材料非線性問題、幾何非線性問題、狀態(tài)非線性問題。

1）材料的非線性問題：指材料的物理定律即其應力——應變關(guān)系是非線性的，但只考慮小位移和小應變的情況，也就是指結(jié)構(gòu)整個幾何形狀的變化及結(jié)構(gòu)材料內(nèi)部的應變與結(jié)構(gòu)尺寸相比，是無限小的，這樣可以忽略微元體的局部應變。例如，計算應力時可以采用原來的、未變形的微元體面積。至于應變——位移關(guān)系則采用線性的應變——位移關(guān)系式。各種小應變、小位移問題的結(jié)構(gòu)彈塑性分析及巖土工程中的問題一般都屬于這一類。

2）幾何非線性問題：雖然假設線性的應力——應變關(guān)系，但非線性的應變——位移關(guān)系和幾何形狀的有限變化卻引起幾何非線性問題。大多數(shù)的幾何非線性問題是小應變和大位移，當然，也包括大位移和大應變的情況。

3）狀態(tài)非線性問題：例如，一根只能拉伸的電纜可能是松弛的，也可能是張緊的；軸承套可能是接觸的，也可能是非接觸的；凍土可能是凍結(jié)的，也可能是非凍結(jié)的；巖土體中的斷層結(jié)構(gòu)面間相對運動。由于系統(tǒng)狀態(tài)的改變，從而引起系統(tǒng)的剛度在不同的值之間突然變化。狀態(tài)改變也許和載荷直接有關(guān)（如電纜情況），也可能由某種外部原因引起（如凍土情況）。

巖土體材料的一個重要特征是其應力——應變關(guān)系具有明顯的非線性性質(zhì)。因此，對于巖土體通常都按材料非線性問題來考慮。材料非線性問題又可以分為兩種情況：一種是非線性彈性問題；另一種就是非線性彈塑性問題，系材料超過屈服極限以后就呈現(xiàn)出非線性性質(zhì)，各種結(jié)構(gòu)和巖土介質(zhì)若采用彈塑性的本構(gòu)模型進行分析就是這類問題。在加載過程中，這兩種非線性問題在本質(zhì)上相同，但是卸載過程就會出現(xiàn)不同的現(xiàn)象，非線性彈性問題是可逆過程，卸載后結(jié)構(gòu)或介質(zhì)會恢復到加載前的位置，非線性彈塑性問題是不可逆的，它將會出現(xiàn)殘余變形。

對于材料非線性問題，由于應力——應變關(guān)系是非線性的，對于這類問題的微分方程的求解在數(shù)學上有一定的困難。然而有限元法處理非線性問題卻是十分有效的。用有限元法處理非線性問題的基本思想是用一系列線性問題的解來逐步逼近非線性問題的解。因而，非線性問題可以理解為一系列線性解進行迭代過程的結(jié)果。

2材料非線性問題的解法

用有限元解材料非線性問題，通常采用三種基本的方法來解析，即增量法、迭代法以及增量迭代法（或混合法）。在說明問題時，可以只考慮單個單元的非線性平衡方程：

[k （δ）]e{δ}={R} （1）

這是是非線性出現(xiàn)在單元剛度矩陣[k]e 中，k [δ]e是非線性材料性質(zhì)[D （ε）]的函數(shù)。可以把[k]e 寫成[k]e =[k{δ}，{R}]的形式，因而[k]e 中的材料參數(shù)不再是常數(shù)。圖1中的載荷{R}與節(jié)點位移{δ}相對的非線性應力——應變關(guān)系。根據(jù)這個應力——應變關(guān)系或本構(gòu)關(guān)系，進而確定非線性分析時變化的矩陣[D （ε）] 。本文主要介紹迭代法。

圖1 非線性的載荷——位移曲線圖

2.1常剛度迭代法

對于方程（1）：[k （δ）]e{δ}={R}，假設R= R0已知，則有限元解法的實質(zhì)就是在R-δ曲線與R= R0直線的交點的橫坐標。為此可用試探法，如果已經(jīng)找到一個試探值δn，但不滿足式（1），如圖1所示，它在曲線上的縱坐標值Rn低于R0值，用Pn（R n，δn）表示。現(xiàn)在可以通過點Pn（R n，δn）繪一條直線p，其斜率為kn ，則

kn （δn+1-δn ） = R0 - Rn （2）

這樣，直線p與R0 - Rn的交點的橫坐標δn+1就是一個改進了的試探解。如果按δn+1 由式（1）計算得的Rn+1仍小于R0，就可以從 Pn+1（R n+1，δn+1）出發(fā)，重復前面過程，雙可以得到進一步的改進試探解R n+2。這樣，循環(huán)迭代，直到前后兩次的位移值接近為止（圖2）。

所謂常剛度，是指取kn=k0為某個常數(shù)，即在迭代時直線p的斜率不變，故稱為常剛度法。

圖2 常剛度迭代法圖

2.2割線變剛度迭代法

用一系列的割線，即直線p與R= R0直線的交點P1S，P2S，P3S，……取逐次逼近R-δ曲線與R= R0直線的交點PS時，這些直線p的斜率并不一定相同。如對第n條直線，其斜率為kn，則迭代的方程為

R0= Rn+Kn（δn+1-δn） （3）

因而，可得第n+1次近似解為

δn+1=δn+ kn （δn+1-δn） （4）

這里的割線剛度kn=kns= Rn /δn。其中，kns為從原點出發(fā)到所迭代時與曲線交點割母的斜率。因而這個過程稱為割線迭代法（圖3）。

圖3 割線變剛度迭代法圖

2.3 切線變剛度迭代法

變切線剛度的表達式可表示為

kn=kn’= （5）

這里kn在物理上代表切線剛度，因而這個過程稱為切線變剛度迭代法（圖4）。

圖4 切線變剛度迭代法圖

2.4 New-Raphson 法

對應于非線性方程式（1），可以考慮一個非線性的初始試探解δ0，它與真實解δ*之差為一個小量，任何一個具有一階導數(shù)的連續(xù)函數(shù)Ψ（δ），在δ0附近光滑，則可按一階泰勒級數(shù)展開為

Ψ（δ）= Ψ（δ0）+ （δ-δ0）+…= 0 （6）

如僅取其線性項，則

Ψ（δ）= Ψ（δ0）+ （δ-δ0）= 0 （7）

所以，δ=δ0- （8）

這個迭代過程，直至在某個允許精度內(nèi) <ε。這種迭代方法就稱為New-Raphson 法，簡稱N-R法。

3結(jié)語

節(jié)理、裂隙、斷層及軟弱夾層是地質(zhì)中常見的現(xiàn)象，并且它們都對巖體的強度和穩(wěn)定性有極大的影響。通常，對于工程巖體而言，有以下兩種處理方法：①把節(jié)理巖體視作不連續(xù)介質(zhì)，由巖石結(jié)構(gòu)體和結(jié)構(gòu)面組成，因此分別研究巖石和結(jié)構(gòu)面的力學性質(zhì)及巖石和結(jié)構(gòu)面共同作用時的耦合原理，這種方法的代表有離散無法、塊體理論和剛彈性法；②把節(jié)理巖體視作宏觀上的連續(xù)體，建立巖體的等效本構(gòu)關(guān)系，這種方法的代表有當量體法、斷裂力學法和損傷力學法。本文根據(jù)巖體中含有節(jié)理、層理、斷層及軟弱夾層等不連續(xù)結(jié)構(gòu)面，研究了塑性有限元分析中常采用的迭代法的詳細過程。