考慮裂紋閉合影響的巖石變形破壞過程模擬研究

彭 勇，陳文富

(1.中國電建集團貴州電力設計研究院有限責任公司，貴州 貴陽 550081；2.貴州鑫達建設工程質量檢測有限公司，貴州 貴陽 550081)

0 引 言

巖石本身具有非均一性、各向異性、非連續性，固其結構特性較為復雜多變，進而導致其多樣的變形破壞模式和錯綜復雜的力學性質。巖石的本構關系涉及到巖石的變形破壞機理，然而目前研究并沒有實質性的進展。損傷理論有助于揭示巖石破壞機理，結合損傷理論的本構模型研究也取得了一定的突破，但仍需向全方位發展，而不僅僅是基于小細節上的改動，特別是針對裂紋擴展情況，使其適用于實際巖石工程。

目前，關于巖石損傷本構模型的研究多數從微觀角度，假定巖石的微元強度符合某種概率分布函數，進而推導損傷演化關系，基于某種強度準則構建損傷本構模型，此類研究已取得了一定進展[1- 4]。較為常用的幾種概率分布函數，包括Weibull分布、正態分布、對數正態分布、冪函數分布和指數函數分布，基于這些分布參數建立的損傷本構模型多數不具備通用性，無法反應不同圍壓下的試驗情況[5]。如Weibull分布函數并不適用于準脆性材料，這是因為準脆性材料具有尺寸效應，若材料是一種準脆性巖石，那么建立損傷本構模型時依托Weibull分布函數就會存在紕漏。此外，正態分布具有兩頭低、中間高的特征，若使巖石的微元強度服從該分布，可能導致巖石強度參數值出現負數，因此依托正態分布函數構建損傷本構模型也有明顯缺陷?；趯嫡龖B分布、冪函數分布和指數函數分布建立的本構模型亦存在較為明顯缺陷，未考慮巖體內部初始裂紋，以及巖體變形過程中裂紋逐漸擴展這一動態現象。

本文為建立一種適用范圍更廣、考慮裂紋擴展情況的損傷本構模型，在現有研究基礎上，假設巖石微元強度服從Harris分布函數，并考慮裂紋擴展情況，進而引入裂紋閉合系數，建立考慮Harris分布和裂紋閉合效應的巖石損傷本構模型，提出了Harris分布參數的確定方法，并基于室內試驗得出的粉砂巖應力應變數據及參考相關文獻的巖石數據對本文模型進行了雙重驗證。

1 巖石損傷模型的建立

首先，基于J.Lemaitre應變等價理論，構建損傷本構模型，表示如下

[σ*]=[σ]/(1-hD)=[E][ε]/(1-hD)

(1)

式中，[E]為巖石材料彈性矩陣；[σ*]為有效應力矩陣；[σ]為名義應力矩陣；[ε]為應變矩陣；D為巖石的損傷變量；h為裂紋閉合系數，取值范圍0≤h≤1。

若巖石破壞前的微元強度表現為線彈性特征，且認定為符合虎克定律，則式(1)可改寫為[6-7]

[ε]=[E][σ*]=[E][σ]/(1-hD)

(2)

據文獻[8]可知，依托Harris分布函數，該函數具有衰減特性，公式如下

(3)

式中，F為自變量；S為因變量；a、b為大于0的模型參數。

假若將Harris概率分布函數作為巖石微元強度服從的一種準則，只需基于式(3)，求自變量F的導數，即可得出巖石微元強度的概率密度函數，用p(F)表示，即

(4)

將p(F)作為損傷變量D，可得巖石的損傷演化方程為

(5)

聯立式(4)、(5)可得

(6)

由式(6)可知，由于p(F)>0，損傷演化方程具有單調遞增的特性。2種極端情況：當F=0時，D=0，對應沒有損失的完整巖石；當F趨向于無窮大時，D=1，對應巖石完全破壞，不再有承載能力。

三軸壓縮情況下，只需聯立式(2)、(6)，即可獲得三軸情況下的損傷統計本構模型，即

(7)

式中，σ1為最大主應力；E為彈性模量；ε1為軸向應變；σ3為圍壓；μ為泊松比。

由于巖石屈服準則多樣性且形式復雜，如若將軸向應力作為微元強度，則不利于巖石損傷模型的建立及求解，其求解形式將十分復雜，推廣應用較難；如若將軸向應變作為微元強度，相比較而言，模型表達式更為簡易、方便求解。因此，微元強度選取軸向應變，即令F=ε1，則根據式(7)可得巖石損傷本構模型，即

(8)

2 分布參數的確定方法

(9)

由條件②可得

(10)

由式(9)、(10)可以求出參數a與b，即

(11)

(12)

聯立式(11)、(12)代入式(8)，構建了考慮裂紋閉合效應的巖石損傷本構模型，該模型包含了Harris分布函數的參數a與b。

3 巖石室內試驗

試驗材料主要為粉砂巖，試驗采用伺服控制巖石三軸壓力試驗儀器，圍壓σ3設定為0、5、10、20 MPa。巖石應力-應變曲線見圖1。經過觀察，巖石破壞形式有共軛剪切、張剪破裂并存、以剪切破壞為主3種。巖石各參數與圍壓的關系見圖2。從圖2可知，巖石峰值應力、峰值應變和殘余強度均隨圍壓的增加而增大，三者與圍壓的擬合均呈較好的線性關系。此外，低圍壓下粉砂巖表現為較強脆性，高圍壓下表現為較強塑性。高圍壓下巖樣的殘余承載能力越強，殘余強度與峰值應力之比越來越大。無圍壓時，比值約為0.17；圍壓為10 MPa時，比值約為0.44；圍壓為10 MPa時，比值約為0.63；而圍壓為20 MPa時，比值約為0.8。

圖1 不同圍壓下的應力-應變關系

圖2 巖石各參數與圍壓的關系

依據試驗曲線獲取本文提出模型曲線。圖3為粉砂巖的試驗曲線與模型曲線對比結果。從圖3可知，各圍壓下的試驗曲線與模型曲線吻合度均較好，峰前部分兩者相似度更甚，峰后匹配效果較峰前差，預測的殘余強度整體高于試驗獲得的殘余強度；此外，本文提出的模型能夠一定程度反映峰后的脆延轉換特性，表明本文提出的模型適用于描述粉砂巖的變形破壞特征。

4 工程實例與模型的進一步驗證

為說明本文提出模型是否適用其他巖性，本文采用參考文獻[9]的巖石三軸壓縮試驗數據，提出試驗曲線數據，獲取4種圍壓下的數據對比驗證。表1為4種圍壓下的模型參數。

表1 不同圍壓下本文模型的參數

裂紋閉合系數值較難確定，本文選定不同的h值，如當h=0.8、0.85、0.9和1時，得出不同裂紋閉合系數下的本構方程。很明顯，當h=1時，建立的本構方程即為前人研究的未考慮裂紋閉合效應的巖石損傷本構方程，故本文模型能夠選用不同的模型曲線進行擬合。

試驗曲線與理論模型曲線比較見圖4。從圖4可以看出，試驗曲線與理論曲線擬合效果較好。特別指出的是，峰前階段的試驗曲線部分與理論曲線部分基本重合，峰后階段也具有相當的匹配程度，說明本文模型是合理和可行的?；谔岢瞿Ｐ瞳@取的應力-應變曲線能夠完全呈現峰前巖石的應力-應變關系和一定程度反映峰后的破壞過程曲線，一定程度上也反映了巖石的應變軟化特性。h值可有效調節應力-應變曲線，峰值前，h值對曲線影響不大，峰值后影響較為明顯。巖石的殘余強度隨著h的減小而增大，然而當h值更大時，更能反映巖石的脆性特征，h=0.9時，試驗曲線與理論曲線的重合度最高。因此，選定不同的h值可以獲取更為廣泛的損傷本構模型曲線，進而取得更好的模擬效果。

圖3 粉砂巖的試驗曲線與模型曲線對比

圖4 試驗曲線與理論曲線的比較

圖5 不同圍壓下的損傷演化曲線

將表1中的Harris分布參數值代入式(6)獲取損傷演化曲線，見圖5。從圖5可以看出，損傷演化曲線大致呈現“S”形，是一個遞增的過程。損傷演化曲線初始階段對應巖石破壞過程的彈性階段，損傷變量很小，其值近乎為0，曲線趨于水平；曲線中間階段對應于屈服階段，該階段維持時間較短，曲線斜率迅速增加；末端曲線又趨于水平，對應于峰后破壞階段，曲線斜率迅速降低，曲線趨于水平。

5 討 論

準確確定概率分布函數中的分布參數十分重要。本文的模型參數是通過多元函數取極值的方法進行計算的，可以看出用該方法具有較好的模擬效果。但筆者認為，模型參數還可以通過曲線擬合的方法獲得。將式(8)變為

(13)

進一步化簡為

(14)

對式(14)進行轉換，公式兩端取對數，得出

(15)

基于式(15)可知，等式左右均為對數關系式，可將每個對數關系式作為一個線性關系式，故可將式(15)用線性關系擬合，可將參數b視作擬合直線的斜率，lna視作擬合直線的截距。進一步依據參考文獻的應力-應變數據進行擬合，獲取模型參數a、b。

本文引入了裂紋閉合系數h，并通過h的取值對試驗曲線進行模擬效果較好。關于裂紋閉合系數的具體求解公式，本文提供一種思路：假設殘余強度處對應的損傷變量為1，根據式(8)可得

(16)

式中，σ1r為殘余強度；σ3r為殘余圍壓；ε1r為最大主應變。

6 結 語

本文考慮裂紋閉合系數，采用Harris分布函數建立了巖石損傷本構模型，確定了模型參數的計算方法，基于粉砂巖的三軸壓縮室內試驗對提出的損傷本構模型進行了驗證，得出以下結論：

(1)巖石峰值應力、峰值應變和殘余強度均隨圍壓呈線性增大，本文模型曲線與試驗曲線較為一致，表明本文提出的模型適用于描述粉砂巖的變形破壞特征。

(3)模型能完全反映巖石的峰前階段應力-應變關系，并能很好地反映峰后階段的破壞過程。選定不同的裂紋閉合系數獲取不同的損傷本構模型，與試驗結果比較，本文模型曲線匹配度較高，符合巖石的變形破壞規律。

(4)探討了模型參數的不同求解方法，后期可通過對比不同方法求解的參數對模型的影響，使巖石損傷本構模型具有更廣的使用范圍。此外，還提供了裂紋閉合系數的一種求解方法。