巖石長時效變形損傷過程的分數階模型研究

肖書文

(山西澤州天泰錦辰煤業有限公司, 山西 晉城 048000)

1 前言

隨著淺部煤礦資源開采殆盡,煤炭回采工作逐漸向深部邁進,深部巷道圍巖長期處于高地應力的環境中會發生蠕變變形[1-3],合理的表征巷道圍巖在長時間作用下的力學變形行為,量化變形規律,對于深部煤礦安全開采具有重要作用。

傳統的巖石蠕變本構理論模型主要包含2大類,分別是經驗模型、組合元件模型。其中經驗模型主要根據實驗數據進行多項式擬合,利用對數函數、指數函數、冪函數等相似函數進行組合,在形式上表征巖石應變與時間的關系,然而此種模型的參數意義不明確,普適性較低,不具備完整的理論力學基礎[4-5]。組合元件模型,主要通過彈性體(胡克體)、黏彈性體(牛頓黏性體)、黏塑性體(圣維南體)進行組合,建立蠕變本構模型。此種模型的物理含義更為明確,具有完備的理論依據,在此基礎上馮夏庭等采用Burgers模型建立了考慮體積蠕變的本構模型。徐衛亞等在黏性元件與塑性元件并聯組合的模型上引入了蠕變系數,再通過與五元件模型進行串聯,得到了可以反映第三蠕變階段的七元件模型[6]。夏才初等構建了包含黏性、黏彈性、黏塑性、黏彈塑性的四元件蠕變模型[7]。

為更好的表征巖石蠕變行為,現有的蠕變模型往往采取增加元件的方式進行優化,然而過多的元件將會導致本構方程的參數增加,大大降低了工程意義。同時,深部巖體通常處于三維應力狀態,傳統的單軸狀態下的蠕變本構方程難以充分考慮三維應力環境。因此,綜合考慮上述問題,采用分數階導數對蠕變元件模型進行改進,可以縮減本構方程的物理參數,同時構建三維的蠕變本構方程,并通過實驗進行有效驗證。

2 分數階導數元件模型

2.1 分數階導數定義

分數階導數相比于整數階導數的區別在于其導數階次為分數,分數階導數具有獨立的數學定義方式,其中Riemann-Liouville微積分的定義方式最為廣泛。設f在(0,+∞)上連續分布,并且在[0,+∞]的任何有限子區間可積,當t>0,Re(γ)>0時,稱式(1)為函數f(t)的γ階Riemann-Liouville分數階積分。

(1)

式中,Γ(γ)為Gamma函數。

通過對式(1)進行逆運算,則得到Riemann-Liouville分數階導數,如式(2)所示。

(2)

式中,m是大于γ的最小整數。

根據分數階導數的定義,將傳統的整數階導數蠕變元件模型進行優化,可得到分數階導數蠕變元件模型。

2.2 分數階導數蠕變元件模型

牛頓黏壺在組合元件模型中用來表示黏性變形,其本構方程如式(3)所示。

(3)

式中,σ為應力,MPa;η為黏性系數;ε為應變,%;t為時間,h;γ為導數的階次。

根據Riemann-Liouville分數階微積分算子理論可將牛頓黏壺變換為以時間為自變量的應變方程。

(4)

在巖石長時效蠕變的過程中,巖石在恒定應力的作用下內部微觀結構會逐漸發生損傷,隨著微觀損傷的逐漸積累,內部微裂隙會聚合擴展為宏觀裂縫,使巖石進入具有強時效特征的第三蠕變階段,并最終破壞。為表征巖石的損傷積累過程,考慮將損傷變量引入到巖石的力學參數中,用于表示巖石力學參數的弱化過程,引入的損傷變量如式(5)所示。

ηγ(D)=ηγ(1-D)

(5)

式中,D為損傷變量。考慮到損傷的積累是非線性的過程,因此,采用指數函數進行表示,如式(6)所示。

D=1-exp(-αt)

(6)

式中,α為損傷劣化系數。

將式(3)、(5)、(6)進行聯立,最終得到考慮損傷過程的分數階導數牛頓黏壺本構方程如式(7)所示。

(7)

3 分數階導數蠕變本構模型

3.1 Nishihara蠕變本構模型

在巖石力學蠕變理論中,蠕變主要包含三個主要階段,分別是第一蠕變階段(初始蠕變階段)、第二蠕變階段(穩定蠕變階段)、第三蠕變階段(加速蠕變階段),在蠕變過程中存在彈性變形、黏彈性變形、黏塑性變形以及損傷演化,可通過將元件模型進行組合得到相應的物理力學過程。日本學者提出的Nishihara蠕變模型被廣泛用于描述第一蠕變階段(初始蠕變階段)、第二蠕變階段(穩定蠕變階段)。Nishihara蠕變模型如圖1所示,主要包括彈性體、黏彈性體與黏塑性體,其總應變表達式為

圖1 Nishihara蠕變模型

ε=εe+εve+εvp

(8)

式中,ε為總應變,%;εe彈性應變,%;εve為黏彈性應變,%;εvp為黏塑性應變,%。

3.2 分數階蠕變本構方程

將Nishihara蠕變模型改進為分數階Nishihara蠕變模型的方法是將模型中原有的牛頓黏壺替換為分數階黏壺。

對圖2所示的分數階Nishihara蠕變模型的組合元件進行分析,推導其本構方程。對彈性體元件進行分析,其本構方程為

圖2 分數階Nishihara蠕變模型

(9)

式中,E0為彈性體的彈性模量,GPa。

對黏彈性體元件進行分析,其本構方程為

(10)

在分數階導數理論中,為便于計算將Riemann-Liouville導數轉化為Caputo導數,兩種導數具有的數學關系為:

(11)

當t=0時,εve=0,利用Laplace雙重變換得到黏彈性體分數階本構方程。

(12)

對黏塑性體元件進行分析,黏塑性體中包含摩擦滑塊,其主要物理意義在于表征應力閾值點,當應力大于閾值應力時,此元件啟動,即代表當應力達到閾值應力時,巖石會發生具有長時效特征的第三蠕變階段(加速蠕變階段),摩擦滑塊的公式為

(13)

式中,σp為摩擦滑塊應力,MPa;σs為閾值應力,MPa。

考慮到蠕變第三階段(加速蠕變階段)的損傷演化趨勢最為明顯,因此,將損傷變量帶入到黏塑性體元件中,經過Laplace雙重變換得到黏塑性體的本構方程。

(14)

將上述彈性體元件、黏彈性體元件、黏塑性體元件進行疊加,得到分數階蠕變本構方程方程。

(15)

為簡化表示式(15),在式中引入雙參數Mittag-Leffler函數(式(16)),將式(15)變換為式(17)。

(16)

式中,α、β為函數參數;z為函數自變量;k為函數階次。

(17)

4 三維分數階蠕變本構模型

深部巖石通常處于三維應力環境下,因此,忽略圍壓的本構方程往往不能有效地表征深部巖石地力學行為,三維應力環境下,蠕變的總應變為

(18)

根據廣義胡克定律,將蠕變變形分解為偏應變與球應變,偏應變狀態下的彈性體三維本構方程為

(19)

式中,eij為應變偏張量,%;sij為應力偏張量,MPa;G0為彈性體的剪切模量,GPa。

相應的,可得到偏應變狀態下的黏彈性體三維本構方程為

(20)

式中,G1為黏彈性體的剪切模量,GPa。

黏塑性體中包含塑性變形部分,因此,在建立三維狀態時,需要引入屈服函數,采用關聯流動法則,黏塑性體三維本構方程為:

(21)

(22)

(23)

式中,F為屈服函數;F0為屈服函數參考值,取1;φ為冪函數,取1;J2為應力偏量第二不變量。

假設兩向圍壓σ2=σ3,則得到以下條件:

(24)

將式(19)～(24)聯立,最終得到偏應變狀態下的三維分數階蠕變本構方程:

(25)

同理,按照相同方式,可得到球應力狀態下的三維分數階蠕變本構方程:

(26)

將式(25)與式(26)合并最終可得到三維全應力狀態下的分數階蠕變本構模型如式(27)所示:

(27)

式中,δij為狄利克雷函數。

5 巖石蠕變實驗驗證

5.1 巖石蠕變實驗方案

為驗證本文推導的三維分數階本構方程,采用晉城市錦辰煤業3采區巷道砂巖進行常規三軸加載實驗。實驗采用Rock Triaxial V4三軸流變儀如圖3所示。

圖3 三軸流變儀

實驗根據砂巖的賦存狀態,同步加載圍壓與軸壓,加載速度為2.0 MPa/min,當加載到20 MPa時,保持圍壓不變,繼續加載軸壓直至達到30 MPa,隨后保持三向應力恒定,直至砂巖破壞。

5.2 巖石蠕變曲線擬合

蠕變實驗共選取3個砂巖試樣開展,將實驗曲線與三維分數階本構模型進行擬合,得到圖4所示的擬合曲線。根據圖中所示,砂巖的蠕變曲線表現出明顯的三階段蠕變變形。對比結果顯示,推導的本構模型可以很好的描述深部砂巖的三軸蠕變實驗數據,擬合參數見表1。

表1 擬合參數表

圖4 三維分數階本構模型與實驗數據擬合曲線

6 結論

(1)根據分數階導數理論,將牛頓黏壺改進為分數階黏壺,考慮到巖石內部的損傷演化,引入指數形式的損傷變量描述蠕變的第三階段。通過改進Nishihara蠕變本構模型建立分數階Nishihara蠕變本構模型。

(2)考慮到深部巖石的三維應力狀態,根據廣義胡克定律,將分數階Nishihara蠕變本構模型分解為偏應力三維分數階蠕變本構方程與球應力三維分數階蠕變本構方程,最終建立全應力狀態下的三維分數階蠕變本構方程。

(3)采用晉城市錦辰煤業3采區巷道砂巖進行常規三軸加載實驗,實驗數據具有明顯的蠕變三階段變形特征。推導的本構模型可以很好的描述深部砂巖的三軸蠕變實驗數據,驗證了模型的可靠性。