基于微米壓痕試驗(yàn)的分?jǐn)?shù)階蠕變模型

穆景福， 張帆, 丁祥*

(1.陜西延長(zhǎng)石油(集團(tuán))有限責(zé)任公司研究院， 西安 710065； 2. 湖北工業(yè)大學(xué)土木建筑與環(huán)境學(xué)院， 武漢 430068;3. 湖北工業(yè)大學(xué)巖土力學(xué)與混凝土材料中法聯(lián)合研究中心， 武漢 430068)

作為一種新型能源，干熱巖利用具有供能持續(xù)穩(wěn)定、高效循環(huán)利用、可再生的特點(diǎn)，可減少溫室氣體排放，改善生態(tài)環(huán)境，在未來(lái)清潔能源發(fā)展中占重要地位,有望成為能源結(jié)構(gòu)轉(zhuǎn)型的新方向，成為全球熱點(diǎn)研究問(wèn)題[1-2]。目前干熱巖的開(kāi)發(fā)主要集中于花崗巖熱儲(chǔ)層中，因此花崗巖的力學(xué)性質(zhì)對(duì)增強(qiáng)型地?zé)嵯到y(tǒng)的設(shè)計(jì)、施工及維護(hù)具有重要影響。由于增強(qiáng)型地?zé)嵯到y(tǒng)生命周期長(zhǎng)達(dá)幾十年，因此花崗巖的長(zhǎng)期力學(xué)性能對(duì)干熱巖的經(jīng)濟(jì)性評(píng)價(jià)及長(zhǎng)期產(chǎn)能預(yù)測(cè)具有十分重要的影響。

中外學(xué)者對(duì)花崗巖的長(zhǎng)期力學(xué)性能開(kāi)展了廣泛的研究。劉凱等[3]基于20～300 ℃下的花崗巖蠕變?cè)囼?yàn),分析溫度對(duì)瞬時(shí)應(yīng)變和蠕變應(yīng)變速率的影響。Chen等[4]基于不同溫度下的花崗巖蠕變實(shí)驗(yàn)，提出了損傷蠕變模型。Dwivedi等[5]對(duì)印度粉花崗巖進(jìn)行了低圍壓30～160 ℃條件下長(zhǎng)達(dá)4個(gè)月蠕變?cè)囼?yàn)。李修磊等[6]在傳統(tǒng)西原模型和Burgers模型的基礎(chǔ)上，將指數(shù)形式的損傷變量、臨界損傷應(yīng)力以及與其有關(guān)的非線性黏性分量引入到流變微分方程，通過(guò)疊加原理推導(dǎo)了考慮溫度效應(yīng)的單軸和三軸壓縮條件下巖石的流變本構(gòu)關(guān)系，建立了巖石的熱-力耦合損傷蠕變本構(gòu)模型。王春萍等[7]針對(duì)北山花崗巖，開(kāi)展了應(yīng)力、溫度分級(jí)加載條件下的蠕變特性試驗(yàn)，研究了花崗巖蠕變破壞過(guò)程中的裂紋擴(kuò)展規(guī)律。王睿[8]開(kāi)展了不同圍壓下的花崗巖分級(jí)加載蠕變實(shí)驗(yàn),獲得了花崗巖三軸蠕變破壞特征，分析了穩(wěn)定蠕變、非穩(wěn)定蠕變階段下的變形和變形速率，建立了考慮損傷演化的花崗巖蠕變本構(gòu)模型。熊良宵等[9]總結(jié)了中國(guó)近20年在高溫作用下和高溫作用后巖石力學(xué)特性研究所取得的成果，發(fā)現(xiàn)開(kāi)展高溫作用下巖石的拉伸斷裂試驗(yàn)和蠕變?cè)囼?yàn)等相對(duì)偏少，中國(guó)學(xué)者主要基于統(tǒng)計(jì)損傷力學(xué)建立了高溫下和高溫后巖石的力學(xué)模型，但有些模型還存在不足，并未通過(guò)實(shí)際工程進(jìn)行驗(yàn)證模型是否適用。郤保平等[10]將巖石現(xiàn)象流變學(xué)與物理流變學(xué)結(jié)合起來(lái)，提出熱力耦合作用下巖石熱黏彈塑性流變?cè)W(xué)元件，在廣義西原模型的基礎(chǔ)上建立熱力耦合作用下花崗巖流變模型,推導(dǎo)出可描述150 MPa及600 ℃以內(nèi)花崗巖的流變本構(gòu)方程。左建平等[11]對(duì)基本元件在溫度和應(yīng)力作用的變形特性做了相應(yīng)的假設(shè)，然后基于西原流變模型，得出熱力耦合作用下的西原模型，當(dāng)把溫度影響因素去除，方程可退化為只受到載荷作用下西原模型的本構(gòu)方程。通過(guò)這些本構(gòu)關(guān)系可大致預(yù)測(cè)特定溫度變化條件下巖石的流變破壞時(shí)間，并得出溫度的變化會(huì)縮短巖石的流變破壞時(shí)間，溫度變化率越大，巖石的破壞時(shí)間越短，即溫度升高變化率每增加一個(gè)數(shù)量級(jí)，破壞時(shí)間就縮短一個(gè)數(shù)量級(jí)。周廣磊等[12]基于巖石變形與熱力學(xué)基本理論，建立了溫度-應(yīng)力耦合作用下脆性巖石時(shí)效蠕變損傷模型。李志強(qiáng)等[13]以片麻狀花崗巖為例，通過(guò)不同溫度下的室內(nèi)分級(jí)加載三軸蠕變?cè)囼?yàn)，獲得了熱力耦合作用下巖石穩(wěn)態(tài)蠕變率和強(qiáng)度的變化規(guī)律，揭示了考慮溫度影響的蠕變破壞模式。張強(qiáng)勇等[14]開(kāi)展了不同溫度和不同加載應(yīng)力路徑條件下片麻狀花崗巖的三軸蠕變?cè)囼?yàn)，系統(tǒng)分析了溫度、圍壓、軸壓對(duì)片麻狀花崗巖蠕變變形特征、蠕變強(qiáng)度和蠕變破壞模式的影響。王建美等[15]考慮溫度-應(yīng)力的耦合效應(yīng)，給出了模型參數(shù)隨溫度及應(yīng)力變化的關(guān)系式，用廣義開(kāi)爾文模型來(lái)表達(dá)高溫高壓環(huán)境中的花崗巖穩(wěn)態(tài)蠕變變形特性。王春萍等[16]通過(guò)將高溫?fù)p傷流變?cè)娼?jīng)典西原模型中Newton元件的方法，構(gòu)建了能夠描述不同溫度條件下花崗巖蠕變?nèi)^(guò)程的本構(gòu)模型。

然而上述對(duì)花崗巖長(zhǎng)期力學(xué)性能的研究基本上都是基于標(biāo)準(zhǔn)尺寸的巖心，一般取心較為困難，同時(shí)試驗(yàn)周期較長(zhǎng)、成本較高、個(gè)體差異性較大及對(duì)目的層精細(xì)研究效果不理想[17]。尋求高效、簡(jiǎn)單、可靠的力學(xué)參數(shù)測(cè)試方法一直是巖石力學(xué)領(lǐng)域的熱點(diǎn)關(guān)注問(wèn)題[18]。最近，一種新型的基于微米壓痕測(cè)試技術(shù)的巖石細(xì)觀力學(xué)特性定量評(píng)價(jià)方法在巖石力學(xué)性能評(píng)價(jià)上表現(xiàn)出良好的適用性[19-21]。

上述基于花崗巖蠕變實(shí)驗(yàn)多基于常規(guī)元件模型的組合。近年來(lái)，分?jǐn)?shù)階導(dǎo)數(shù)在描述各類(lèi)復(fù)雜力學(xué)與物理行為方面取得較為廣泛的應(yīng)用，特別是在記憶性較強(qiáng)過(guò)程建模過(guò)程中，如各類(lèi)流變模型的建模中[22]。國(guó)內(nèi)外學(xué)者在分?jǐn)?shù)階流變模型方面也開(kāi)展了廣泛研究。陳文[23]探討了軟物質(zhì)的時(shí)間-空間分?jǐn)?shù)階導(dǎo)數(shù)力學(xué)建模問(wèn)題。徐明瑜等[24]介紹了分?jǐn)?shù)階算子理論、方法的最新進(jìn)展以及在現(xiàn)代力學(xué)各個(gè)領(lǐng)域中的應(yīng)用。郭霄怡[25]研究了分?jǐn)?shù)階微積分在黏彈性流體力學(xué)中的某些應(yīng)用。殷德順等[26]用分?jǐn)?shù)階微積分理論提出等應(yīng)變率加載情況下的軟土應(yīng)力-應(yīng)變關(guān)系。宋勇軍等[27]提出了一種基于分?jǐn)?shù)階微積分的巖石非線性蠕變損傷力學(xué)模型。唐皓[28]建立了基于分?jǐn)?shù)階微積分構(gòu)建的大理巖黏彈塑性流變模型。丁靖洋等[29]通過(guò)引入分?jǐn)?shù)階微積分理論及損傷理論，構(gòu)建基于損傷演化的鹽巖分?jǐn)?shù)階三元件模型。何明明等[30]構(gòu)建了巖石在單三軸壓縮、松弛及蠕變條件下的分?jǐn)?shù)階體積應(yīng)變-軸向應(yīng)變關(guān)系。閆云明等[31]根據(jù)泥巖在不同圍壓下分級(jí)加載蠕變?cè)囼?yàn)結(jié)果，給出了一種帶應(yīng)變觸發(fā)的分?jǐn)?shù)階蠕變本構(gòu)模型。吳池等[32]對(duì)兩種不同雜質(zhì)鹽巖進(jìn)行分級(jí)加載下的三軸蠕變?cè)囼?yàn)，建立了耦合雜質(zhì)的鹽巖分?jǐn)?shù)階非定常蠕變模型。周宏偉等[33]通過(guò)將分?jǐn)?shù)階Abel黏壺代替經(jīng)典西原模型中Newton黏壺的方法, 構(gòu)建了基于分?jǐn)?shù)階導(dǎo)數(shù)的鹽巖流變本構(gòu)模型。李任杰等[34]利用水泥砂漿制作非貫通硬性結(jié)構(gòu)面試樣，進(jìn)行剪切蠕變?cè)囼?yàn)，并對(duì)蠕變特性進(jìn)行分析，采用元件組合方式，構(gòu)建基于分?jǐn)?shù)階微積分的巖石硬性結(jié)構(gòu)面蠕變本構(gòu)模型。

可見(jiàn)，在花崗巖長(zhǎng)期力學(xué)性能及分?jǐn)?shù)階建模上取得了較好的進(jìn)展，尚未見(jiàn)到基于細(xì)觀壓痕試驗(yàn)對(duì)花崗巖長(zhǎng)期力學(xué)性能的研究。因此，現(xiàn)利用花崗巖微米壓痕蠕變?cè)囼?yàn)，建立花崗巖分?jǐn)?shù)階蠕變模型，研究花崗巖長(zhǎng)期力學(xué)性能，以期為干熱巖設(shè)計(jì)施工提供有力支撐。

1 微米壓痕蠕變?cè)囼?yàn)

1.1 試驗(yàn)設(shè)備

采用湖北工業(yè)大學(xué)CB500型微納米力學(xué)綜合測(cè)試系統(tǒng)進(jìn)行，如圖1所示。CB500是美國(guó)NANOVEA公司推出的一款微納米力學(xué)綜合測(cè)試系統(tǒng)，該系統(tǒng)采用模塊化設(shè)計(jì)，可選配納米壓痕儀、納米劃痕儀、微米壓痕儀、微米壓痕儀4個(gè)功能模塊，可實(shí)現(xiàn)納米、微米、大載荷3個(gè)尺度下的壓痕，劃痕與摩擦磨損測(cè)試，進(jìn)而可得到硬度、彈性模量、蠕變信息、彈塑性、斷裂韌度、應(yīng)力-應(yīng)變曲線、膜基結(jié)合力、劃痕硬度、摩擦因數(shù)、磨損率等微觀力學(xué)數(shù)據(jù)。

圖1 CB500型微納米力學(xué)綜合測(cè)試系統(tǒng)Fig.1 CB500 micro/nano mechanical testing system

該系統(tǒng)壓痕試驗(yàn)符合國(guó)際ISO14577[35]與美國(guó)ASTM E3-01[36]標(biāo)準(zhǔn)。采用閉環(huán)載荷加載垂直加載，準(zhǔn)確性優(yōu)于傳統(tǒng)的開(kāi)環(huán)載荷加載技術(shù)及懸臂加載技術(shù)；同時(shí)采用高精度壓電陶瓷驅(qū)動(dòng)，精度優(yōu)于電磁力驅(qū)動(dòng)。

1.2 試驗(yàn)樣品

花崗巖試樣取自湖南岳陽(yáng)，外觀完整，無(wú)可見(jiàn)裂紋。根據(jù)X射線衍射儀測(cè)試顯示，該花崗巖樣品含石英29.52%，鈉長(zhǎng)石20.24%，鉀長(zhǎng)石22.18%，黑云母11.95%，綠泥石9.59%，伊利石6.52%。

將試件加工成直徑50 mm、高度30 mm的圓餅狀試件，如圖2所示。用專(zhuān)業(yè)砂紙對(duì)試件表面進(jìn)行機(jī)械拋光，然后進(jìn)行離子二次拋光，以保證花崗巖試件表面粗糙度滿足要求。

圖2 花崗巖試件Fig.2 Granite specimen

1.3 試驗(yàn)方法

采用網(wǎng)格化壓痕技術(shù)[37]，在滿足試驗(yàn)壓痕間距要求條件下盡可能多地布置壓入點(diǎn)位，如圖3所示。試驗(yàn)過(guò)程中，壓痕點(diǎn)位從左至右，從上至下依次進(jìn)行壓入試驗(yàn)。試驗(yàn)過(guò)程中實(shí)際壓入深度不超過(guò)20微米，按照間距大于20倍最大壓入深度原則，需保證壓入點(diǎn)間距400 μm(即0.4 mm)以上。

每一個(gè)壓痕點(diǎn)試驗(yàn)加卸載過(guò)程如圖4所示，通過(guò)預(yù)壓試驗(yàn)確定持荷階段最大載荷為15 N，采用線性加載，以0.5 N/s的加載速率加載至15 N，維持該荷載180 s，然后采用0.5 N/s的速率卸載，加載曲線如圖4中藍(lán)色線條所示；試驗(yàn)過(guò)程中記錄壓頭嵌入深度，如圖4中紅色線條所示。

圖3 網(wǎng)格化壓痕點(diǎn)布置Fig.3 Layout of grid-indentation

圖4 加卸載過(guò)程示意圖Fig.4 Illustration of loading and unloading

2 分?jǐn)?shù)階流變模型

2.1 分?jǐn)?shù)階微積分

分?jǐn)?shù)階微積分有多種定義方式，常見(jiàn)的有Riemann-Liouville定義、Caputo定義、Grtinwald-Letnikov定義等，研究采用Riemann-Liouville定義[38]。

假定函數(shù)f(t)在區(qū)間(0,+∞)上連續(xù)，且在[0,+∞)上任意有限子集上可積，則對(duì)于t>0，Re(α)>0，則

(1)

設(shè)f∈C，α>0，m為大于μ的最小整數(shù)，α=m-μ，則

Dμf(t)=Dm[D-αf(t)]

(2)

式(2)為函數(shù)f(t)的μ階Riemann-Liouville分?jǐn)?shù)階微分。

2.2 Abel黏壺

Abel黏壺本構(gòu)方程可以寫(xiě)為

(3)

式(3)中：當(dāng)α=0時(shí)，描述胡克固體的本構(gòu)關(guān)系；當(dāng)α=1時(shí)，描述牛頓流體的本構(gòu)關(guān)系；當(dāng)0≤α≤1時(shí)，描述介于中間狀態(tài)的黏彈性體的本構(gòu)關(guān)系。

令式中σ(t)=σ0H(t)，H(t)為Heaviside函數(shù)，當(dāng)t<0時(shí)，H(t)=0；當(dāng)t≥0時(shí)，H(t)=1，此時(shí)可模擬Abel黏壺在恒定應(yīng)力下的蠕變行為。則式(3)寫(xiě)成

(4)

將式(4)進(jìn)行拉普拉斯變換可得

(5)

可解得

(6)

將式(6)進(jìn)行拉普拉斯逆變換可得

(7)

故Abel黏壺在常應(yīng)力下的蠕變模量為

(8)

根據(jù)Boltzmann疊加原理可得

(9)

(10)

將式(9)和式(10)進(jìn)行拉普拉斯變換可得

ε(s)=J(s)sσ(s)

(11)

σ(s)=G(s)sε(s)

(12)

通過(guò)式(11)和式(12)可得

(13)

從式(6)可得

(14)

結(jié)合式(13)和式(14)可得

(15)

將式(15)進(jìn)行拉普拉斯逆變換可得常應(yīng)變下的松弛模量為

(16)

根據(jù)式(8)和式(16)繪制不同α?xí)r蠕變曲線[圖5(a)]和松弛曲線[圖5(b)]。

圖5 Abel黏壺不同α取值Fig.5 Abel dashpot with varying α

圖5(a)中，α=1描述了牛頓流體的特征，牛頓流體在剪切力作用下應(yīng)變隨時(shí)間線性增長(zhǎng)；α<1時(shí)，應(yīng)變隨時(shí)間非線性增長(zhǎng)，并且隨著α的減小，應(yīng)變?cè)鲩L(zhǎng)速率也在減小，描述了蠕變應(yīng)變的非線性增長(zhǎng)。圖5(b)中，α<1其描述的應(yīng)力松弛過(guò)程也較牛頓流體的松弛緩慢，應(yīng)力有緩慢減小的過(guò)程。α取值不同，Abel具有不同的特性。當(dāng)α=0時(shí)，Abel黏壺可以描述彈簧元件的行為，當(dāng)α=1時(shí)，可以描述牛頓黏壺的行為。因此，Abel黏壺可以視為一個(gè)同時(shí)包含彈簧及牛頓黏壺的組合元件。

2.3 分?jǐn)?shù)階Maxwell模型

傳統(tǒng)的Maxwell模型由一個(gè)彈簧元件和一個(gè)牛頓黏壺串聯(lián)組成[42]，如圖6(a)所示。Maxwell模型的總體變形由胡克彈簧的瞬時(shí)彈性應(yīng)變和牛頓黏壺的黏性應(yīng)變組成。圖6(b)為分?jǐn)?shù)階Maxwell模型[43-45]，它與傳統(tǒng)Maxwell模型的不同在于用Abel黏壺取代了傳統(tǒng)Maxwell模型中的牛頓黏壺。

圖6 Maxwell模型Fig.6 Maxwell model

根據(jù)2.2節(jié)對(duì)Abel黏壺的研究，結(jié)合分?jǐn)?shù)階Maxwell模型中各元件的應(yīng)力及應(yīng)變與總體變形及應(yīng)力的關(guān)系，給出分?jǐn)?shù)階Maxwell模型的本構(gòu)方程為

(17)

式(17)中：ε、εe、εv分別為總體應(yīng)變、彈性應(yīng)變和黏性應(yīng)變；σ、σe、σv分別為總應(yīng)力、彈性應(yīng)力和黏性應(yīng)力；E為彈簧元件的彈性模量；ηα為黏性系數(shù)。將式(17)整理后得到

(18)

將式(18)進(jìn)行拉普拉斯變換可得

(19)

又因?yàn)?/p>

(20)

結(jié)合式(19)和式(20)可得

(21)

將式(21)進(jìn)行拉普拉斯逆變換可得

(22)

因此，恒定載荷下分?jǐn)?shù)階Maxwell模型預(yù)測(cè)的應(yīng)變?yōu)?/p>

(23)

3 結(jié)果分析

(24)

(25)

式中：β為壓頭形狀參數(shù)；F為蠕變?cè)囼?yàn)時(shí)施加的恒定荷載；h(t)為壓入深度。對(duì)式(25)變換可得應(yīng)變?yōu)?/p>

(26)

式(26)中：h(0)為蠕變起始時(shí)刻的壓入深度。

共進(jìn)行了188次壓入試驗(yàn)，其中13次為無(wú)效數(shù)據(jù)，剩余有效數(shù)據(jù)175組。隨機(jī)挑選部分壓入點(diǎn)位，利用式(24)及式(26)分別計(jì)算微米壓痕蠕變?cè)囼?yàn)中應(yīng)力及應(yīng)變，并通過(guò)分?jǐn)?shù)階Maxwell模型擬合蠕變?cè)囼?yàn)曲線。擬合參數(shù)如表1所示，擬合結(jié)果如圖7所示。同時(shí)，對(duì)所有175個(gè)有效壓入點(diǎn)數(shù)據(jù)擬合后，統(tǒng)計(jì)了參數(shù)E、ηα、α的分布規(guī)律，如圖8所示。

僅關(guān)注研究對(duì)象花崗巖的時(shí)間相關(guān)變形，其他參數(shù)如硬度、斷裂韌性等參數(shù)未重點(diǎn)關(guān)注。圖7中紅色點(diǎn)位為測(cè)試過(guò)程中波動(dòng)，采樣間隔約為0.06 s。藍(lán)色曲線為分?jǐn)?shù)階Maxwell模型擬合的壓入深度隨時(shí)間變形關(guān)系，較好地描述了恒定載荷下壓頭壓入深度隨時(shí)間的增長(zhǎng)，因此分?jǐn)?shù)階Maxwell模型是一種良好的描述花崗巖時(shí)間相關(guān)細(xì)觀變形的工具。

表1 擬合參數(shù)Table 1 Fitting parameters

圖7 蠕變階段壓入深度隨時(shí)間變化Fig.7 Indentation depth with time during creep stage

圖8 擬合參數(shù)分布Fig.8 Distribution of fitting parameters

從條形統(tǒng)計(jì)圖[圖8(a)～圖8(c)]可見(jiàn)，擬合參數(shù)E、ηα、α大致服從正態(tài)分布；從盒須圖[圖8(d)～圖8(f)]可見(jiàn)，E中位數(shù)為74.96 GPa，ηα中位數(shù)為20.39 GPa/sα，α中位數(shù)為0.282。

值得注意的是，擬合參數(shù)的E、ηα、α與花崗巖組成成分之間關(guān)系尚不明確，有待下一步工作中進(jìn)一步開(kāi)展。

4 結(jié)論

基于網(wǎng)格化微米壓痕蠕變?cè)囼?yàn)，研究了花崗巖細(xì)觀力學(xué)參數(shù)，并用分?jǐn)?shù)階Maxwell模型擬合了花崗巖蠕變變形，基于上述研究，主要得到以下結(jié)論。

(1)不同于傳統(tǒng)標(biāo)準(zhǔn)尺寸蠕變?cè)囼?yàn)，微米壓痕試驗(yàn)可以利用少量巖心開(kāi)展大量試驗(yàn)，試驗(yàn)周期相對(duì)較短，是一種具有良好前景的巖石力學(xué)參數(shù)測(cè)試方法。

(2)Abel黏壺是一種介于線彈性固體和線性流體之間的黏彈性體，通過(guò)調(diào)整分?jǐn)?shù)階階數(shù)α描述介于牛頓流體及彈簧之間的黏彈性體變形特征。

(3)分?jǐn)?shù)階Maxwell模型較好地描述了花崗巖恒定載荷下的細(xì)觀變形特征，是一種描述花崗巖時(shí)間相關(guān)變形的良好工具。

(4)基于分?jǐn)?shù)階Maxwell模型的花崗巖細(xì)觀力學(xué)參數(shù)與花崗巖成分及宏觀力學(xué)參數(shù)間的關(guān)系尚不明確，有待進(jìn)一步研究。