考慮初始壓密階段的千枚巖變形破壞全過程統(tǒng)計損傷本構(gòu)模型

李占甫，莘子健 ，付聰，李玉豪，張家銘

(1.安徽省交控建設(shè)管理有限公司，安徽 合肥 230088；2.中國地質(zhì)大學(xué)（武漢）工程學(xué)院，湖北 武漢 430074)

1 引言

巖石是由具有一定結(jié)構(gòu)構(gòu)造的礦物集合體組成的非均質(zhì)材料，因其生成條件、礦物成分、賦存環(huán)境等不同，使得巖石內(nèi)部存在大量隨機分布的微觀裂隙。巖石的宏觀破壞可看作其內(nèi)部微觀裂隙的擴展和累積，該過程稱之為損傷[1]。損傷理論以材料內(nèi)部微觀裂隙的萌生、擴展和相互作用為基礎(chǔ)，研究試件由微觀直至發(fā)生宏觀破壞的全過程，是目前研究巖石等含有天然微裂隙材料的有效方法[2]。基于巖石內(nèi)部微觀裂隙的分布及演化具有隨機性，Krajcinovic等人[3]將統(tǒng)計強度理論與連續(xù)損傷理論相結(jié)合，開辟了一條能較好反映巖石損傷破壞過程和應(yīng)力—應(yīng)變關(guān)系的本構(gòu)模型建立途徑，即統(tǒng)計損傷本構(gòu)模型。

自統(tǒng)計損傷理論提出以來，國內(nèi)外學(xué)者利用其對巖石變形破壞過程的模擬方法進行了大量探索，并取得較為豐碩的研究成果。曹文貴[4]先基于Lemaitre應(yīng)變等價性假設(shè)[5]定義了傳統(tǒng)的巖石損傷變量，石崇[6]、曹瑞瑯[7]引入不同修正系數(shù)對其修正，但因該定義忽略了巖石破壞后的殘余強度，與實際不符，隨后曹文貴[8]將巖石骨架部分抽象為損傷與未損傷兩部分，對巖石的損傷重新定義，并建立可以反映巖石殘余強度的損傷模型，尹杰[2]在此基礎(chǔ)上基于不同準則進行計算和修正。周永強[9]和劉冬橋[10]分別考慮損傷閾值和微元強度的度量方式建立統(tǒng)計損傷模型。盡管當前理論模型能夠較好描述巖石部分階段的變形特征，但對于巖石變形破壞全過程的模擬仍存在一定不足，特別是在反映巖石初始壓密階段。究其原因，當前沿用較廣的統(tǒng)計損傷模型忽略了巖石內(nèi)部在初始狀態(tài)下存在的大量初始缺陷，并在模型推導(dǎo)過程中將巖石彈性模量視為常數(shù)。事實上巖石在受力過程中，因內(nèi)部微結(jié)構(gòu)面的存在其巖石宏觀變形與力學(xué)特性會產(chǎn)生較大變化，如千枚巖，因其特殊的千枚狀構(gòu)造，其微結(jié)構(gòu)面發(fā)育，該類巖石在初始壓密階段的應(yīng)力應(yīng)變關(guān)系會呈現(xiàn)較強的非線性變形特征。此外，大量試驗表明，巖石的彈性模量與初始缺陷的變形及閉合情況密切相關(guān)，并在一定范圍內(nèi)隨所受圍壓水平的變化而改變，將其視為常數(shù)并不恰當。

為此，曹文貴[11]將巖石抽象為空隙與骨架兩部分，采用宏觀與微觀相結(jié)合的方法，建立了可以反映脆性巖石初始宏觀非線性變形力學(xué)行為的模擬方法，但其模型在模擬煤巖在圍壓為8MPa以下時擬合效果不理想；張超[12]、李修磊[13]在文獻[11]的基礎(chǔ)上對模型進行了改進，雖較好地擬合了砂巖在圍壓40MPa下的應(yīng)力應(yīng)變關(guān)系，但該巖石變形破壞過程中初始壓密階段并不明顯，不能很好地驗證該類模型對于反映含大量微結(jié)構(gòu)面巖石初始壓密階段的變形特征。此外，前述三位學(xué)者均采用Mohr-Coulomb準則建立相關(guān)模型，該準則在描述巖石處于低應(yīng)力狀態(tài)下會存在一定偏差[2]。

基于上述分析，盡管現(xiàn)有巖石統(tǒng)計損傷模型對于巖石部分變形階段特征可以較好描述，但能夠反映考慮巖石初始壓密階段的全過程變形特征的本構(gòu)模型較少。此外，多數(shù)模型所針對的巖石為砂巖、花崗巖以及粉砂質(zhì)泥巖等，而對于含有大量初始微結(jié)構(gòu)面的千枚巖研究鮮有報道。針對上述不足，本文在前人的研究基礎(chǔ)上，基于文獻[11]提出的考慮巖石初始宏觀非線性變形力學(xué)行為的模擬方法，選用Drucker-Prager準則針對千枚巖在較低圍壓狀態(tài)下的全過程變形特征進行模擬，以期建立出能反映巖石全過程變形特征的統(tǒng)計損傷模型。

2 巖石變形過程

巖石受外部荷載作用下的變形可劃分為巖石骨架和內(nèi)部微結(jié)構(gòu)面兩部分，且?guī)r石骨架在屈服破壞前的變形為線彈性。低荷載條件下巖石的骨架變形與微裂隙閉合同時發(fā)生，應(yīng)力應(yīng)變曲線呈明顯非線性；當荷載增加到一定水平時，微裂隙閉合完成，此時隨荷載的增加至巖石屈服破壞前，巖石僅發(fā)生線彈性的骨架變形；荷載進一步增加達到巖石屈服強度后，巖石的骨架變形呈非線性，并依次出現(xiàn)屈服、應(yīng)變軟化和完全破壞的現(xiàn)象。因此，巖石在受荷變形直至破壞過程中具有階段性變形特征：初始壓密(OA)、線彈性變形(AB)、屈服(BC)、應(yīng)變軟化(CD)和完全破壞階段（D點之后），如圖1所示。

圖1 千枚巖變形破壞全過程

3 巖石變形分析模型

為分析巖石在初始壓密階段內(nèi)部微裂隙產(chǎn)生的變形，在巖石內(nèi)取一代表性柱體單元，該柱體單元由巖石骨架和微裂隙兩部分組成，如圖2所示。由前述巖石變形機理可知，在初始壓密階段微裂隙變形與骨架變形并不協(xié)調(diào)，故設(shè)柱體單元在加載前初始高度為h0，其中微裂隙部分高度為，骨架部分高度為該柱體單元在某荷載應(yīng)力σi下總變形量為Δh，巖石內(nèi)部微裂隙和骨架部分變形量為Δhw和Δhr，則有如下關(guān)系：

圖2 千枚巖變形分析模型

故巖石在荷載應(yīng)力σi的作用下的總應(yīng)變（宏觀應(yīng)變）εi、微裂隙部分應(yīng)變和骨架部分應(yīng)變可表示為：

為建立巖石總應(yīng)變與其組成部分應(yīng)變間的關(guān)系，特定義柱體單元微裂隙部分占比為b0,即利用式（1）～（5）可得：

分析式（6）可看出巖石在初始壓密階段（b0≠0）的宏觀應(yīng)變由微裂隙與骨架兩部分組成，且呈非線性，當微裂隙完全閉合（b0=0）后，巖石宏觀應(yīng)變僅由骨架部分提供。由巖石變形破壞過程知，在巖石未達到屈服強度前骨架部分材料的變形線彈性狀態(tài)，進一步證明在初始壓密階段的非線性變形是由微裂隙閉合導(dǎo)致的，顯然式（6）作為巖石變形分析模型是符合實際的。

3.1 微裂隙部分變形分析方法

按前述分析思路，若想獲得巖石微裂隙在荷載應(yīng)力σi作用下的應(yīng)變需先求得相應(yīng)微裂隙部分（hw）的變形量Δhw，故將荷載σi劃分為n個等級組成的應(yīng)力增量并逐級施加，故有：

巖石微裂隙部分的閉合變形具有不可恢復(fù)的特點，采用真應(yīng)變能夠更準確地描述其變形特點，故將變形量Δhw可看作在某一級應(yīng)力增量作用下變形量的累加，設(shè)微裂隙部分（hw）在應(yīng)力增量作用后的高度變?yōu)椴捎梦墨I[14]方法可得此時微裂隙部分的應(yīng)變可表示為：

由于微裂隙部分材料的累計變形高度無法通過試驗測量，則假定服從廣義胡克定律，微裂隙部分材料在第s級應(yīng)力增量作用下產(chǎn)生的應(yīng)變增量可表示為：

聯(lián)立式（8）和式（9）并變形有：

由式（11）可知，若想得到空隙部分在荷載應(yīng)力σi下的應(yīng)變得先求出在各級荷載作用下相應(yīng)的變形模量和泊松比，但常規(guī)試驗無法獲得，為便于研究，假設(shè)微裂隙部分材料的力學(xué)參數(shù)不隨巖石變形而改變，從而式（11）可改寫為：

式中E1和μ1為巖石微裂隙部分材料的彈性模量和泊松比，相關(guān)參數(shù)的確定方法詳見第4.2。

3.2 巖石骨架部分變形分析方法

巖石骨架部分的變形可根據(jù)荷載是否達到屈服強度分為兩部分，未達到屈服強度時呈線彈性變形，超過屈服強度巖石骨架發(fā)生損傷，致使變形呈非線性，并依次出現(xiàn)應(yīng)變硬化、應(yīng)變軟化和完全破壞的現(xiàn)象，對此可采用基于連續(xù)介質(zhì)的統(tǒng)計損傷理論進行分析。目前常用的統(tǒng)計損傷模型基于Lemaitre應(yīng)變等效假設(shè)，并不能準確反映巖石骨架部分損傷的力學(xué)本質(zhì)，尤其是認為巖石完全破壞后不具備任何承載能力，即忽略了巖石的殘余強度，為此，曹文貴等[8]對此類模型進行修正，以解決該類本構(gòu)模型的不足，即：

式中，σi為名義應(yīng)力為有效應(yīng)力為巖石的軸向殘余強度；D為損傷變量，是材料損傷程度的度量。

將巖石骨架部分抽象為未損傷和損傷兩部分，骨架整體受力為σi，其作用面積為A，則未損傷材料面積為A1，損傷部分面積為A2，兩部分所受微觀應(yīng)力分別為如圖3所示。則損傷變量D可定義為損傷部分面積與巖石骨架總面積之比，即D=A2/(A1+A2)。

圖3 千枚巖骨架部分變形分析模型

基于Lemaitre應(yīng)變等效假設(shè)，可認為巖石骨架部分的宏觀應(yīng)變與未損傷部分微觀應(yīng)變相等。為骨架未損傷部分變形的微觀應(yīng)力，且骨架未損傷部分為線彈性體，服從廣義胡克定律，則有：

式中，E2和μ2為巖石骨架部分彈性模量和泊松比為骨架未損傷部分所受有效應(yīng)力。假設(shè)骨架部分材料的損傷僅沿軸向發(fā)生，即在側(cè)向方向上名義應(yīng)力等于有效應(yīng)力，則有：

由式（13）、（14）及上述假設(shè)可得：

故巖石骨架部分變形可表示為：

4 巖石統(tǒng)計損傷本構(gòu)模型及參數(shù)確定方法

4.1 巖石統(tǒng)計損傷本構(gòu)模型

對于巖石骨架部分的變形破壞為一連續(xù)損傷過程，可采用統(tǒng)計損傷理論進行模擬，模擬關(guān)鍵在于對巖石骨架部分微元強度的合理度量，經(jīng)對比不同強度準則的優(yōu)劣，本文選用Drucker-Prager準則進行相關(guān)計算。該準則表達式為：

式中F為巖石骨架部分微元強度，取值分兩種情況：當骨架部分達到損傷閾值時，F(xiàn)≥0未達到時，變形呈線彈性，F(xiàn)＜0；I1和J2分別為應(yīng)力張量的第一不變量和應(yīng)力偏量的第二不變量；a和k分別為與巖石骨架部分粘聚力（c）、內(nèi)摩擦角（φ）相關(guān)的常數(shù)。上述參量可表示為：

由式（14）代入式（18）和（19），并在常規(guī)三軸壓縮試驗條件下（σ2=σ3），可得：

為得到損傷變量D，采用可以較好表征骨架材料微元強度變形破壞隨機性的Weibull概率密度函數(shù)，則損傷變量D與微元強度F間的關(guān)系可表示為：

式中m與F0為巖石骨架部分微元強度Weibull分布參數(shù)，當巖石處于完全破壞階段時D=1，巖石未發(fā)生損傷時D=0。將式（12）、（13）、（16）代入式（6）可得到考慮巖石內(nèi)部微裂隙及殘余強度的巖石變形破壞全過程模型，即：

4.2 模型參數(shù)確定方法

①b0的確定方法

由前述巖石變形破壞過程可知，巖石宏觀變形由微裂隙閉合與骨架變形兩部分組成，當微裂隙完全閉合時且未達到屈服強度前，骨架部分材料處于線彈性狀態(tài)，此時E1趨近于0，式（25）中第1式中可變?yōu)椋?/p>

若σ3=0時，即巖石處于單軸壓縮狀態(tài)，式（26）可變?yōu)椋?/p>

式（27）經(jīng)轉(zhuǎn)換可得：

由式（28）可知，b0為巖石單軸壓縮試驗曲線線彈性段的延長線在應(yīng)變軸上的截距。

②E2和μ2的確定方法

μ2為巖石骨架材料的泊松比，故可通過對已完成微裂隙壓密階段的巖石進行泊松比試驗直接測定。E2為巖石骨架部分彈性模量，可通過式（26）變換并在三軸試驗中彈性變形階段取值代入求得，即：

③E1和μ1的確定方法

在已知巖石在初始壓密階段滿足式（25）中第1式和已求得模型參數(shù)b0、E2和μ2的基礎(chǔ)上，在初始壓密階段任意取值，采用基于最小二乘法原理的曲線擬合方法可得E1和μ1。

④微元強度隨機分布參數(shù)的確定

微元強度隨機分布參數(shù)m與F0為模型建立的關(guān)鍵，由巖石具有應(yīng)變軟化特性可知，其應(yīng)力應(yīng)變關(guān)系曲線在峰值點（σ1=σsc，ε1=εsc）處具有極值特性，即：

式中σsc和εsc分別為巖石三軸試驗中峰值點的應(yīng)力和應(yīng)變。

同時式（30）滿足式（25）中第2式，聯(lián)立方程即可確定參數(shù)m與F0：

至此，已建立考慮巖石內(nèi)部微裂隙及殘余強度的巖石變形破壞全過程模擬方法，并給出相關(guān)參數(shù)的確定方法，現(xiàn)對其合理性進行驗證分析。

5 實例分析與驗證

引入文獻[15]的試驗曲線對前述模型進行驗證分析。利用該文獻針對千枚巖在結(jié)構(gòu)面夾角為0°時，圍壓分別為0、2、4、6和8MPa下的巖石三軸試驗曲線，如圖4所示。由該文獻可知黏聚力C=11.66MPa，內(nèi)摩擦角φ=50.35°，按照前述參數(shù)確定方法可知：b0=0.0017，不同圍壓下相關(guān)計算參數(shù)見表1所示。將相關(guān)計算參數(shù)代入式（31）～式（34）可求得微元強度隨機分布參數(shù)m與F0。將試驗數(shù)據(jù)及各計算參數(shù)代入本文所提出的本構(gòu)模型，得到不同圍壓作用下千枚巖全過程變形應(yīng)力應(yīng)變理論曲線，并將其與試驗曲線相比較，如圖5所示。

模型參數(shù) 表1

圖4 千枚巖三軸試驗曲線

圖5 本文理論曲線與試驗曲線的比較

為進一步驗證本文理論模型和方法的合理性與可行性，分別根據(jù)文獻[11]、文獻[12]和文獻[13]所提出的本構(gòu)模型，以文獻[15]圍壓4MPa試驗數(shù)據(jù)為基礎(chǔ)，分別獲取各模型理論曲線，并與本文模型所獲理論曲線相對比，如圖6所示。

圖6 圍壓4MPa下不同模型理論曲線與試驗曲線的比較

由圖5可以看出，本文所提出的統(tǒng)計損傷本構(gòu)模型，針對千枚巖在不同圍壓條件下的計算結(jié)果與試驗曲線均有較高的吻合程度，很好地反映千枚巖變形破壞全過程應(yīng)力應(yīng)變特征。以圍壓4MPa試驗數(shù)據(jù)為例，比較已有模型和本文模型對其模擬效果，可以看出，由于文獻[11]所建立模型并未考慮巖石破壞后的殘余強度，導(dǎo)致其殘余強度變形階段與試驗曲線差別較大。此外，文獻[11]、文獻[12]和文獻[13]所建立模型均基于Mohr-Coulomb準則，致使其在低圍壓狀態(tài)下與試驗曲線差別較大，特別是在巖石屈服和應(yīng)變軟化階段。由此可見，本文所構(gòu)建模型能夠較好地模擬含大量微結(jié)構(gòu)面的千枚巖在低圍壓狀態(tài)下的全過程變形破壞特征。

6 結(jié)語

本文利用統(tǒng)計損傷理論，對能反映巖石在初始壓密、屈服破壞及殘余變形破壞等階段的巖石變形全過程模擬方法進行了研究。得出如下結(jié)論：

①將巖石抽象為骨架和內(nèi)部空隙兩部分，分別建立了空隙和骨架部分變形分析方法，進而建立了可模擬空隙巖石變形破壞全過程的統(tǒng)計損傷本構(gòu)模型；

②由于Mohr-Coulomb準則，并不能很好描述巖石處于低應(yīng)力狀態(tài)下的變形破壞特征，本文選用Drucker-Prager準則，并以含有大量初始微結(jié)構(gòu)面的千枚巖為對象，對其在低圍壓狀態(tài)下的全過程變形特征進行模擬研究；

③通過本文模型、現(xiàn)有相關(guān)模型和試驗曲線的比較，驗證了本文所建巖石損傷模型的有效性和合理性。