考慮時(shí)效損傷的巖石黏彈塑性本構(gòu)模型

曹利軍，馬 超，王 偉

（水發(fā)規(guī)劃設(shè)計(jì)有限公司，濟(jì)南250014）

蠕變特性是指材料受外界荷載、化學(xué)作用、溫度、水等因素的影響，其應(yīng)變表現(xiàn)出顯著的時(shí)效特性［1］。巖石蠕變特性關(guān)系到巖體工程長期穩(wěn)定性，是實(shí)際工程中必須重視的內(nèi)容之一［2］。巖石在外界應(yīng)力條件下蠕變變形逐漸累積，外界應(yīng)力條件改變影響蠕變變形發(fā)展，達(dá)到一定應(yīng)力水準(zhǔn)時(shí)，巖石便表現(xiàn)出加速蠕變行為，短時(shí)間內(nèi)巖石破壞，巷道隧洞、高邊坡、核廢料儲(chǔ)存室、油田套管失穩(wěn)破壞常常與應(yīng)力條件改變導(dǎo)致的巖石加速蠕變行為有關(guān)，因此在巖體工程全壽命周期中，巖石蠕變特性不可忽視［3，4］。

本構(gòu)模型是描述巖石蠕變特性的核心手段，通過不同元件自由組合而成的元件模型由于其靈活性得到了廣泛應(yīng)用，但建立能精準(zhǔn)識(shí)別不同類型巖石蠕變行為（尤其是加速蠕變）的本構(gòu)模型仍是一個(gè)難點(diǎn)。Shibata等［5］通過搭建蠕變應(yīng)變率與流變參數(shù)之間的聯(lián)系，建立可較好描述凝灰?guī)r蠕變行為的本構(gòu)模型；Günther 等［6］開展不同應(yīng)力和溫度下的鹽巖蠕變?cè)囼?yàn)，通過建立穩(wěn)態(tài)蠕變速率和外界影響因素之間的關(guān)系，得到了較為精準(zhǔn)的鹽巖蠕變模型；唐佳等［7］以白云石英片巖和綠泥鈉長片巖為研究對(duì)象，構(gòu)建了能反映黏彈塑性偏量特征的改進(jìn)Burgers模型；盧軍振［8］針對(duì)硬巖的脆性蠕變破壞現(xiàn)象，探究硬巖脆性蠕變破壞時(shí)間和破壞后的應(yīng)力變形之間的聯(lián)系，從而建立能較好表征硬巖蠕變行為的本構(gòu)模型；黃海峰等［9］引入分?jǐn)?shù)階微積分來描述軟巖蠕變的黏彈性和黏塑性應(yīng)變，與損傷彈性體串聯(lián)，得到一個(gè)新的蠕變損傷模型，并證明能較好地反映軟巖蠕變變形特征；金俊超等［10］提出一種基于應(yīng)變軟化指標(biāo)的蠕變力學(xué)模型，能較好地模擬泥巖和砂巖蠕變現(xiàn)象。

目前已有大量蠕變本構(gòu)模型，但多數(shù)模型難以同時(shí)描述軟巖和硬巖的加速蠕變行為，鑒于此，本文以廣義Burgers 模型作為基礎(chǔ)模型，引入損傷力學(xué)理論，建立考慮時(shí)效損傷的黏塑性體，從而建立一維巖石黏彈塑性蠕變損傷模型，并推廣至三維應(yīng)力狀態(tài)。給出模型參數(shù)求取方法，并進(jìn)行參數(shù)敏感性分析。采用所建模型辨識(shí)5 種不同堅(jiān)硬程度巖石的蠕變?cè)囼?yàn)數(shù)據(jù)，驗(yàn)證所建模型的合理性和適用性。

1 蠕變損傷模型的建立

1.1 基礎(chǔ)模型的選擇

巖石作為一種復(fù)雜的非線性地質(zhì)材料，其內(nèi)部同時(shí)賦存著彈性、黏性、黏彈性和黏塑性特征，元件模型理論［1］中常用的有彈性元件、黏性元件和塑性元件，但單一的某個(gè)元件難以準(zhǔn)確描述材料蠕變性質(zhì)，因此需對(duì)不同元件進(jìn)行自由組合從而得到組合元件模型。傳統(tǒng)Burgers 模型是組合元件模型中的一類經(jīng)典模型，被眾多學(xué)者［1，7，8］證明能較好地描述巖石衰減、穩(wěn)定蠕變，但部分學(xué)者發(fā)現(xiàn)［11，12］傳統(tǒng)Burgers 模型描述硬巖蠕變特征時(shí)誤差較大，由此學(xué)者們?cè)趥鹘y(tǒng)Burgers 模型的基礎(chǔ)上，串聯(lián)了一個(gè)Kelvin 體，得到了適用范圍更廣的廣義Burgers 蠕變模型［11］，廣義Burgers蠕變模型示意圖如圖1所示。

由圖1可看出，廣義Burgers模型與傳統(tǒng)Burgers模型不同之處在于多串聯(lián)了一下Kelvin體，其蠕變狀態(tài)方程為

圖1 廣義Burgers蠕變模型Fig.1 Generalized Burgers creep model

式中：σ0為初始應(yīng)力；σ1、ε1、E1和η1分別為Maxwell體的應(yīng)力、應(yīng)變、彈性模量和黏滯系數(shù)；σ2、ε2、E2和η2分別為Kelvin-1體的應(yīng)力、應(yīng)變、彈性模量和黏滯系數(shù)；σ3、ε3、E3和η3分別為Kelvin-2體的應(yīng)力、應(yīng)變、彈性模量和黏滯系數(shù)；上標(biāo)圓點(diǎn)表示某變量對(duì)時(shí)間t的一階導(dǎo)數(shù)。

通過對(duì)式（1）進(jìn)行求解得到廣義Burgers 模型的一維蠕變方程：

引用文獻(xiàn)［13，14］中的大理巖和輝綠巖蠕變?cè)囼?yàn)數(shù)據(jù)，利用數(shù)學(xué)軟件Origin，基于Levenberg-Marquardt 算法，驗(yàn)證一維廣義Burgers 模型和一維傳統(tǒng)Burgers 模型描述硬巖穩(wěn)定蠕變行為的差異，對(duì)比驗(yàn)證結(jié)果如圖2所示。

由圖2 可知，傳統(tǒng)Burgers 模型擬合曲線形態(tài)相對(duì)固定，衰減、穩(wěn)定蠕變階段之間部分較為“外凸”，理論值遠(yuǎn)高于試驗(yàn)值，穩(wěn)定蠕變階段理論值曲線近似為直線性線段，該階段快結(jié)束部分的曲線明顯高于理論值，傳統(tǒng)Burgers模型對(duì)大理巖和花崗巖的擬合效果一般，R2分別為0.933 8 和0.957 7，廣義Burgers 模型對(duì)大理巖和花崗巖的擬合效果較好，R2分別為0.999 0 和0.990 5。鑒于廣義Burgers模型描述硬巖穩(wěn)定蠕變的優(yōu)越性，考慮到蠕變模型對(duì)不同種類、不同硬度巖石的適用性，擇取適用范圍更廣的廣義Burgers模型作為本文模型的基礎(chǔ)模型。

圖2 對(duì)比驗(yàn)證曲線Fig.2 Comparison validation curves

1.2 考慮時(shí)效損傷的黏塑性體

廣義Burgers 蠕變模型結(jié)構(gòu)為H-N-H｜N-H｜N-N｜S，具備一定的反映巖石彈性、黏性、黏彈性變形特征的能力，但無法描述巖石的黏塑性變形。由此引入一個(gè)黏塑性體，其微分本構(gòu)方程為

當(dāng)應(yīng)力超過損傷應(yīng)力閾值（即長期強(qiáng)度），巖石材料才能進(jìn)入加速蠕變階段，該階段黏塑性變形急劇增加，巖石材料內(nèi)部微裂紋不斷發(fā)育和擴(kuò)展，內(nèi)部損傷不斷累積，造成巖石力學(xué)性能的宏觀劣化［9］。由此在黏塑性體的基礎(chǔ)上，引入損傷力學(xué)理論，根據(jù)Lemaitre等效應(yīng)力原理［15］，將式（3）變形為：

式中：D（t）為損傷變量。

這里需注意的是，部分學(xué)者［16，17］按照Lemaitre 應(yīng)變等效原理進(jìn)行損傷演化時(shí)，直接將本構(gòu)方程中的應(yīng)力替換成等效應(yīng)力，譬如將直接替換為，這樣的做法是不對(duì)的。實(shí)際上損傷變量是關(guān)于時(shí)間的函數(shù)，積分后D（t）會(huì)發(fā)生變化，正確的方法應(yīng)該是將微分本構(gòu)方程中的應(yīng)力替換成等效應(yīng)力，例如本文式（3）演化為式（4）的方式［18］。

Kachanov［19］認(rèn)為在一維應(yīng)力狀態(tài)下，蠕變條件下的損傷發(fā)展方程為：

式中：A和ν為材料參數(shù)。

對(duì)式（5）進(jìn)行積分可得蠕變破壞時(shí)間te為：

聯(lián)立式（5）和式（6）可得D的演化規(guī)律：

將式（7）代入式（4）得：

對(duì)式（8）進(jìn)行積分可得：

式中：C為積分常數(shù)。

當(dāng)t=0時(shí)，黏塑性應(yīng)變?yōu)榱?，將該條件代入式（9）有：

將式（10）代入式（9）變形可得：

式（11）即為考慮時(shí)效損傷的黏塑性體。

1.3 一維條件下的蠕變損傷模型

為使蠕變本構(gòu)模型可同時(shí)反映巖石彈性、黏性、黏彈性和黏塑性特征，將考慮時(shí)效損傷的黏塑性體與廣義Burgers模型串聯(lián)，得到一個(gè)新的蠕變損傷模型，如圖3所示。

圖3 蠕變損傷模型Fig.3 Creep damage model

結(jié)合式（2）和（11）得到一維條件下的蠕變本構(gòu)方程：

式中：σ0為初始應(yīng)力；σs為長期強(qiáng)度。

式（12）即為本文所建一維條件下的蠕變損傷模型。

1.4 三維蠕變損傷本構(gòu)模型

在巖體工程中，巖石多處于三向應(yīng)力狀態(tài)，由此需構(gòu)建三維條件下的蠕變本構(gòu)模型。假設(shè)巖石各向同性，依據(jù)廣義Hooke定律有：

式中：σm為球應(yīng)力張量；Sij為偏應(yīng)力張量；K為體積模量；G為剪切模量。

分解巖石內(nèi)部的應(yīng)力張量，于是有：

式中：σij為應(yīng)力張量；δij為單位張量。

一般認(rèn)為，巖石在三維應(yīng)力下不發(fā)生體積蠕變［13］，故本文三維蠕變方程中不考慮體積蠕變，并引入屈服函數(shù)F來描述巖石的屈服。由此將式（12）推廣為三維應(yīng)力狀態(tài)：

式中：（Sij）0為初始偏應(yīng)力張量；（Sij）S為對(duì)應(yīng)長期強(qiáng)度的偏應(yīng)力張量；G1、G2和G3為對(duì)應(yīng)E1、E2和E3的剪切模量；H1、H2、H3和H4為三維黏滯系數(shù)；F0為屈服函數(shù)的初始值。

屈服函數(shù)可以取如下形式［13］：

式中：J2為應(yīng)力偏量第二不變量。

在等圍壓三軸壓縮應(yīng)力狀態(tài)中有：

將式（17）代入式（15）可得：

式（18）即為本文所建三維巖石蠕變損傷模型的軸向蠕變方程。

2 模型可行性驗(yàn)證及參數(shù)求解

2.1 巖石蠕變?cè)囼?yàn)

試樣取自某水電工程邊坡平硐內(nèi)頁巖巖層，運(yùn)回實(shí)驗(yàn)室后制備成直徑50 mm、高100 mm 的圓柱樣。蠕變?cè)囼?yàn)采用RLW-2000 型巖石三軸流變?cè)囼?yàn)系統(tǒng)，將圍壓分別設(shè)置為1、2 和3 MPa，共施加5級(jí)軸向荷載。通過試驗(yàn)數(shù)據(jù)采集系統(tǒng)，得到3種圍壓下的分級(jí)加載蠕變曲線（偏應(yīng)力標(biāo)示在曲線上），如圖4所示。

圖4 分級(jí)加載蠕變曲線Fig.4 Creep curves under step loading

由圖4可看出，巖石在應(yīng)力作用下首先出現(xiàn)瞬時(shí)彈性變形，接著進(jìn)入衰減蠕變階段，該階段蠕變速率遞減，當(dāng)蠕變速率衰減到某一恒定值附近時(shí)，便進(jìn)入了穩(wěn)定蠕變階段，當(dāng)應(yīng)力加載到某一級(jí)別時(shí)，巖石還表現(xiàn)出歷時(shí)較短的加速蠕變，隨后巖石破壞?？傮w上，隨著圍壓的提升，對(duì)應(yīng)級(jí)別的應(yīng)變量值遞增。由于式（18）針對(duì)每一級(jí)加載，故還需根據(jù)玻爾茲曼原理［7，8］對(duì)圖4進(jìn)行處理，得到分別加載蠕變曲線，從而對(duì)每一級(jí)加載的蠕變曲線進(jìn)行擬合辨識(shí)。

2.2 參數(shù)求解

（1）參數(shù)G1和K。根據(jù)廣義Hooke定律有。

式中：μ為泊松比。

根據(jù)等圍壓三軸壓縮試驗(yàn)可求得E1和μ，從而確定G1和K，限于篇幅，這里不再給出等圍壓三軸壓縮試驗(yàn)結(jié)果曲線。

（2）參數(shù)G2、G3、H1、H2和H3。巖石蠕變?cè)囼?yàn)前4 級(jí)加載中，未出現(xiàn)加速蠕變階段，通過式（18）的第一式求取參數(shù)，通過數(shù)學(xué)軟件Origin，基于Levenberg-Marquardt 算法，進(jìn)行非線性擬合便可確定參數(shù)G2、G3、H1、H2和H3。

（3）參數(shù)H4、te、ν和長期強(qiáng)度（σ1-σ3）s。當(dāng)巖石外界應(yīng)力水平超過長期強(qiáng)度時(shí)，還會(huì)表現(xiàn)出加速蠕變階段。采用式（18）的第二式求取參數(shù)H4和ν，根據(jù)試驗(yàn)總歷時(shí)確定蠕變破壞時(shí)間te，根據(jù)等時(shí)應(yīng)力-應(yīng)變曲線法［20］確定長期強(qiáng)度，三種圍壓下的長期強(qiáng)度分別為4.75、5.68 和6.64 MPa，限于篇幅，這里不再給出求取過程。

2.3 頁巖蠕變模擬

通過式（18）辨識(shí)頁巖分別加載蠕變曲線，同時(shí)利用傳統(tǒng)三維Burgers 模型作為對(duì)比驗(yàn)證，得到理論值與試驗(yàn)值對(duì)比曲線，如圖5所示。模型參數(shù)如表1所列。

表1 模型參數(shù)Tab.1 Parameters of model

圖5 理論值與試驗(yàn)值對(duì)比曲線Fig.5 Comparison curves between theoretical value and experimental value

由圖5 和表1 可看出，本文所建模型對(duì)頁巖蠕變?cè)囼?yàn)數(shù)據(jù)擬合效果較好，平均R2達(dá)到了0.986 8。傳統(tǒng)三維Burgers 模型的平均R2僅有0.926 4，對(duì)前4 級(jí)蠕變數(shù)據(jù)擬合誤差相對(duì)較小，但擬合值始終略低于試驗(yàn)值，無法辨識(shí)最后一級(jí)蠕變數(shù)據(jù)，難以描述巖石加速蠕變階段。

3 參數(shù)分析及適用性驗(yàn)證

3.1 參數(shù)敏感性分析

為了研究材料參數(shù)ν對(duì)巖石蠕變發(fā)展的影響，判斷其敏感性，將不同ν值和表1 中剩余已確定模型參數(shù)代入式（18），以頁巖圍壓1 MPa 下最后一級(jí)應(yīng)力水平為例，繪制不同ν值的蠕變曲線，如圖6（a）所示。以同樣的方法繪制不同H4值的蠕變曲線，如圖6（b）所示。

圖6（a）中ν=0.35 時(shí)的蠕變曲線即為圖5（a）中最后一級(jí)理論值曲線，當(dāng)ν值逐漸增加時(shí)，蠕變曲線逐漸遠(yuǎn)離時(shí)間軸，曲線斜率遞增。材料參數(shù)ν影響蠕變速率和極限蠕變值，ν值的變化明顯影響曲線形態(tài)。當(dāng)ν值從0.25 增長至1.00 時(shí)，ν值增長了300%，應(yīng)變值從0.032 4 增加至0.061 9，應(yīng)變值增加了91.05%，ν值平均每1%的變化引起了應(yīng)變值0.30%的變化。

由圖6（b）可知，應(yīng)變值隨著H4值的減小而遞增，黏滯系數(shù)H4值影響蠕變速率和極限蠕變值，同時(shí)影響曲線形態(tài)。當(dāng)H4值從300 GPa·h 增長至700 GPa·h 時(shí)，H4值增長了133.33%，應(yīng)變值從0.045 5 減小至0.031 2，應(yīng)變值減少了31.43%。H4值平均每1%的變化引起了應(yīng)變值0.24%的變化，與ν值對(duì)比發(fā)現(xiàn)，材料參數(shù)ν的敏感性高于黏塑性體的黏滯系數(shù)H4。

圖6 不同參數(shù)值的蠕變曲線Fig.6 Creep curves with different parameters

3.2 模型適用性驗(yàn)證

巖石按照?qǐng)?jiān)硬程度將巖石劃分為堅(jiān)硬巖、較硬巖、較軟巖、軟巖和極軟巖［21］，本文蠕變?cè)囼?yàn)對(duì)象頁巖屬于軟巖。為了驗(yàn)證本文所建模型描述不同類型、不同硬度巖石蠕變力學(xué)行為的適用性，引用文獻(xiàn)［22-25］中的石英砂巖、板巖、凝灰?guī)r和全風(fēng)化花崗巖的蠕變?cè)囼?yàn)數(shù)據(jù)，通過本文所建模型進(jìn)行擬合。由于Burgers模型無法描述巖石加速蠕變階段，引用文獻(xiàn)［26］中可反映加速蠕變過程黏塑性變形的三維Cvisc 模型進(jìn)行對(duì)比驗(yàn)證，得到理論值與試驗(yàn)值對(duì)比曲線，如圖7所示。

圖7 理論值與試驗(yàn)值對(duì)比曲線Fig.7 Comparison curves between theoretical value and experimental value

分析圖7，利用三維Cvisc 模型辨識(shí)4 種巖石的R2分別為0.928 9、0.902 4、0.983 3 和0.881 7，本文所建模型辨識(shí)4 種巖石的R2分別為0.997 3、0.975 2、0.991 6 和0.953 8。石英砂巖、板巖、凝灰?guī)r和全風(fēng)化花崗巖分別屬于堅(jiān)硬巖、較硬巖、較軟巖和極軟巖［21］，由于不同類型、不同硬度巖石的蠕變曲線形態(tài)各異，盡管三維Cvisc 模型擬合凝灰?guī)r蠕變數(shù)據(jù)時(shí)具備較好的模擬能力，但對(duì)石英砂巖、板巖和全風(fēng)化花崗巖的擬合效果較差，適用范圍較小。本文所建模型總體上擬合精度更高，能較好地識(shí)別5 種堅(jiān)硬程度巖石的蠕變力學(xué)行為，由此證明本文所建蠕變損傷本構(gòu)模型反映不同類型巖石蠕變特性的適用性。

4 結(jié)論

（1）根據(jù)巖石蠕變黏彈塑性變形特性，選擇廣義Burgers 模型作為基礎(chǔ)模型，引入損傷力學(xué)理論，引入損傷力學(xué)理論，建立考慮時(shí)效損傷的黏塑性體，與基礎(chǔ)模型串聯(lián)，得到新的考慮時(shí)效損傷的巖石黏彈塑性蠕變本構(gòu)模型，并推廣至三維應(yīng)力狀態(tài)。

（2）利用所建模型辨識(shí)頁巖蠕變?cè)囼?yàn)數(shù)據(jù)，給出參數(shù)求解辦法，并進(jìn)行參數(shù)敏感性分析，材料參數(shù)ν的敏感性高于黏塑性體的黏滯系數(shù)H4。對(duì)比試驗(yàn)曲線和理論曲線，證明所建模型的可行性和合理性。

（3）引用相關(guān)文獻(xiàn)中不同堅(jiān)硬程度的巖石蠕變?cè)囼?yàn)數(shù)據(jù)，利用所建模型和三維Cvisc 模型進(jìn)行擬合對(duì)比，證明了本文所建模型反映不同類型、硬度巖石蠕變力學(xué)行為的優(yōu)越性，從而驗(yàn)證所建考慮時(shí)效損傷的巖石黏彈塑性模型的適用性。目前多數(shù)模型難以反映不同巖石加速蠕變特性，所建模型通過驗(yàn)證可描述不同類型、硬度巖石蠕變行為，對(duì)不同加速蠕變曲線仍有良好辨識(shí)能力，具有較廣闊的工程應(yīng)用前景?！?/p>