基于Maxwell分布的砂巖本構模型研究

童立紅，吳琳琳，徐長節(jié)

(1.華東交通大學 省部共建軌道交通基礎設施性能監(jiān)測與保障國家重點實驗室，江西 南昌 330013；2.華東交通大學 江西省巖土工程基礎設施安全與控制重點實驗室，江西 南昌 330013)

人類進行生產(chǎn)和活動的地球表面層主要由巖石圈構成，其中80%為沉積巖[1]，受地理環(huán)境和自然風貌的影響，各地巖石的礦物成分、顆粒結構和物理力學性質(zhì)等皆不相同，甚至同一巖石的各部分性質(zhì)也因其內(nèi)部不規(guī)則分布的微裂紋和微孔隙的不同而略有差異。通過研究巖石本構關系來了解巖石在受力過程中的變形和破壞特性一直是學者們的研究熱點之一。早在1982年，KRAJCI‐NOVIC等[2]基于統(tǒng)計強度理論思想對損傷進行了定義，并將之融入到巖石本構關系中，以此研究巖石力學破壞機制。國內(nèi)外學者在此基礎上展開了一系列研究，帶來許多創(chuàng)新成果。大量試驗表明，在傳統(tǒng)本構關系中引入數(shù)理統(tǒng)計方法來描述巖石應力-應變關系或損傷演化關系是合理的。CAO等[3]使用正態(tài)分布來描述微元強度分布，進一步建立了損傷本構模型，該模型僅包含常見的明確的物理參數(shù)，應用簡單，但其可能出現(xiàn)參數(shù)為負的情況，不符合客觀事實。蔣維等[4]選擇對數(shù)正態(tài)分布來構建本構模型，很好地解決了使用正態(tài)分布過程中出現(xiàn)的參數(shù)為負的現(xiàn)象。然而，更多學者[5-10]認為巖石微元強度破壞服從Weibull分布，以此構建的損傷本構模型能夠很好地描述巖石變形的全過程，尤其是低應變水平下的初始階段和屈服后的殘余變形階段。在此基礎上，通過運用不同的強度準則來描述巖石在不同環(huán)境下的破壞特征。例如：基于Mohr-Coulomb準則構建的巖石損傷演化模型，劉樹新等[11-12]引入第3個參數(shù)對模型進一步修正，更好地解釋了巖石峰后變形特征，王蘇生等[13]引入孔隙水壓力，研究了孔隙率對損傷閾值的影響；王輝等[14-15]利用Drucker-Prager強度準則建立巖石損傷演化模型，結果表明Weibull分布參數(shù)F0，m在不同受力情況下都具有合適的物理基礎，進一步證明了模型的適用性；高瑋等[16]假定應力水平滿足Hoek-Brown準則，從而建立泥質(zhì)砂巖破裂本構模型，并用數(shù)值模擬試驗再次驗證了模型的可靠性；劉之喜[17]認為裂紋數(shù)服從Weibull分布且裂紋擴展服從Griffith準則，建立了高圍壓下巖石循環(huán)加-卸載統(tǒng)計損傷本構模型。此外，周峙等[18]假設裂土微元強度服從Laplace分布，并考慮初始損傷的影響，建立了干濕循環(huán)作用下的損傷演化模型。然而，上述統(tǒng)計分布模型的使用均是基于經(jīng)驗而提出，缺乏相關的嚴格推導和理論基礎。尤其是使用最廣泛的Weibull分布，其分布參數(shù)F0，m的求解方法不一[7]，早期人們使用線性回歸法和反演分析法，前者計算簡便但參數(shù)本身無明確物理意義，后者相對精確但對試驗數(shù)據(jù)要求較高，后期研究發(fā)現(xiàn)峰值點法不僅處理簡便且物理意義較清晰。然而，使用峰值點法得到的分布參數(shù)的物理意義是對試驗結果的唯象描述，缺乏本質(zhì)的物理基礎。TONG等[19]基于統(tǒng)計力學基本理論，認為破壞微元強度需要的能量以最概然分布的形式存在，經(jīng)過嚴密的論證證明顆粒材料的微觀接觸狀態(tài)符合Maxwell分布。盡管TONG等[19]的推導是針對無黏性顆粒材料提出，但是經(jīng)過數(shù)學分析發(fā)現(xiàn)，其可以直接適用于黏性或者膠結顆粒組成的材料系統(tǒng)。因此，本文在其研究基礎上，假設砂巖微觀顆粒狀態(tài)服從Maxwell分布，基于Lemaitre應變等價原理構建砂巖本構模型，并通過室內(nèi)試驗對其進行驗證。

1 構建本構模型

砂巖主要由砂粒膠結形成，具備明顯的顆粒結構特征。云南砂巖屬于泥質(zhì)砂巖，顆粒細膩，其由不同尺寸的微觀顆粒相互膠結而成，一般認為，不同顆粒之間的膠結強度不同，隨著應變的增加，低強度的膠結鍵先破壞，高強度的膠結鍵后破壞[20]。承受荷載時，當微觀顆粒膠結鍵破壞達到一定數(shù)量時，表現(xiàn)出宏觀上的整體破壞。從統(tǒng)計力學的角度來看，不同強度膠結鍵的破壞規(guī)律符合某種統(tǒng)計分布，由于微觀顆粒的數(shù)量足夠多，可以認為膠結鍵破壞強度是連續(xù)分布的。這里假設每個膠結鍵都對應一個能量狀態(tài)，由于巖石經(jīng)過足夠長時間的沉積后，其達到了最穩(wěn)定的狀態(tài)，則該狀態(tài)對應的是膠結鍵能量的最概然分布。當外部荷載足夠大時，或者輸入的破壞能量足夠高時，部分膠結鍵的能量不足以抵抗外部破壞能量，則這部分膠結鍵會產(chǎn)生破壞，宏觀上表現(xiàn)為巖石強度的下降。假設每個膠結鍵的破壞能量以變形勢能的形式存儲在顆粒系統(tǒng)中，根據(jù)TONG等[19]的工作，顆粒接觸能量對應的應變符合Maxwell分布，即巖石系統(tǒng)顆粒產(chǎn)生破壞的概率密度函數(shù)為：

式中：εc為顆粒接觸能量對應的應變；β2為巖石系統(tǒng)能量密度的參數(shù)。

當外部施加荷載后(即輸入了破壞能量，這里以三軸壓縮為例)，巖石產(chǎn)生軸向應變ε。根據(jù)上面的闡述，對于接觸應變εc∈[0,ε]的膠結顆粒將會產(chǎn)生破壞，因此巖石系統(tǒng)顆粒產(chǎn)生破壞的概率為：

式中：erf(ε)為誤差函數(shù)。

根據(jù)Lemaitre應變等價理論[21]，名義應力作用在受損材料上引起的應變與真實應力作用在相同的無損材料上引起的應變是相同的，即：

式中：σ*為真實應力；E為彈性模量；σ為名義應力；D為損傷變量。

若受損巖石單元即受損巖體不再為整體提供強度，則：

假設在受力過程中未損傷巖體仍然滿足廣義虎克定律，可得：

式中：ν為泊松比。

把式(4)代入式(5)，且兩邊同時乘以1-D可得：

在本試驗中，先施加圍壓至預定值，再施加軸壓直到巖樣破壞。其中σ2=σ3，因此，可以對式(6)進一步簡化得到式(7)。

由于所施加軸向應力為偏應力，即試驗數(shù)據(jù)中的軸向應力σ如式(8)所示。因此，在試驗應力-應變曲線記錄之前會產(chǎn)生一個初始應變ε0，其大小如式(9)所示。

式中：σ3為圍壓；E50為50%峰值強度處的割線模量。

綜合式(7)~(9)，可以得到應力-應變關系如下：

由式(2)可知微元破壞的概率為P(ε)，其認為微觀顆粒間膠結鍵的方向均為垂直或水平結構，而實際上膠結鍵方向是隨機分布的。通過對試驗數(shù)據(jù)進行分析可知，對P(ε)進行α次冪處理，可以很好地解決該問題，使其更符合客觀事實。因此，破壞微元數(shù)量占總微元數(shù)量的比例為Pα(ε)，則，可定義損傷變量D如式(11)所示。

把式(11)代入式(10)中，可得巖石本構模型如下：

2 試驗與模型驗證

2.1 試驗

為驗證本文所提出的本構模型的適用性，本文采用型號為ZTRE-210的微機控制巖石三軸測試系統(tǒng)展開了一系列不同圍壓下的三軸壓縮試驗。試驗所用微機控制巖石三軸壓縮系統(tǒng)如圖1所示，使用軸向和環(huán)向應變計實時測量記錄軸向應變和環(huán)向應變。該試驗設備最大試驗荷載為2 000 kN，最大圍壓為100 MPa。按照國際巖石力學學會試驗委員會的規(guī)定，將來自云南省楚雄彝族自治州武定縣的紅砂巖沿沉積方向取芯制成?50 mm×100 mm的標準圓柱體試樣，試樣表面光滑、色澤均勻、上下切面平整，如圖2所示。

圖1 微機控制巖石三軸壓縮系統(tǒng)示意圖Fig.1 Schematic diagram of microcomputer controlled rock triaxial compression system

圖2 砂巖試樣Fig.2 Sandstone sample

試驗時，首先施加1 kN的壓力使加載探頭與試樣頂面接觸預壓，保證試樣處于軸心受壓狀態(tài)；接著施加圍壓至設定值(分別為0，4，7和10 MPa)，穩(wěn)定后施加軸壓直至試樣破壞，其中加載速率為0.02 mm/min；最后，先卸載軸壓再卸載圍壓。試樣破壞形式如圖3所示，均為剪切破壞。如圖3(a)和3(b)所示，在圍壓為0和4 MPa時，巖樣破壞比較嚴重，存在多條裂縫且縫隙較大；在圍壓為7 MPa和10 MPa時，如圖3(c)和3(d)所示，巖樣主要存在一條主裂縫且縫隙較小，說明較大圍壓約束了巖石剪切面的發(fā)展。

圖3 試樣破壞圖Fig.3 Specimen failure diagram

2.2 模型驗證

在本構模型式(12)中，圍壓σ3為已知量，割線模量E50和泊松比ν通過應力-應變曲線獲得，初始應變ε0通過式(9)計算得到，未知量Pα(ε)取決于參數(shù)α和β。因此，通過對試驗曲線進行非線性擬合可得到參數(shù)α和β，不同圍壓下的曲線擬合結果如圖4所示。由圖4可知，在不同圍壓條件下試驗數(shù)據(jù)與理論曲線都吻合良好，且理論曲線能較好地反映應力-應變關系的各個階段。

圖4 應力-應變曲線及其理論曲線Fig.4 Stress-strain curve and its theoretical curve

3 討論與分析

3.1 模型參數(shù)的物理意義

為了明確參數(shù)α和β的物理意義，將其中一個參數(shù)固定，改變另一個參數(shù)，由式(12)得到2組應力-應變曲線，如圖5所示。由圖5(a)可知，保持α不變，每當β上升20%，峰值強度也隨之上升20%，且峰值右移，但其峰值附近的彎曲程度基本一致，即β的大小反映了巖石的宏觀強度和韌性，β越大巖石的強度越高、延展性越好。由圖5(b)可知，保持β不變，每當α上升20%，峰值強度隨之上升10%，峰值也發(fā)生右移，說明α的大小對巖石強度和延展性也存在影響，但其影響峰值強度的能力弱于參數(shù)β。因此，可以將參數(shù)α和β視為巖石強度和延展性的指標。

圖5 參數(shù)變化示意圖Fig.5 Schematic diagram of parameter change

3.2 圍壓對模型參數(shù)的影響

根據(jù)式(12)擬合得到的模型參數(shù)如表1所示，參數(shù)α隨著圍壓的增大而減小，參數(shù)β隨著圍壓的增大而增大。為了更直觀地描述圍壓與模型參數(shù)之間的關系，對表1中的數(shù)據(jù)進行數(shù)學擬合，擬合結果如圖6所示。由圖4可知，隨著圍壓的增大，巖樣峰值強度增大，屈服過程逐漸平緩，圍壓增強了巖石抵抗壓縮變形的能力。結合上一節(jié)內(nèi)容，如圖6(a)所示，參數(shù)α與圍壓存在良好的線性遞減關系；如圖6(b)所示，參數(shù)β與圍壓呈二次遞增關系。

表1 參數(shù)擬合結果Table 1 Parameter fitting results

3.3 不同圍壓下的損傷演化規(guī)律

利用式(2)可得到微元破壞概率與應變的關系，如圖7所示。在不同圍壓下，當應變較小時，微元破壞概率趨于0，且隨著圍壓的增大其趨于0的應變區(qū)間變大。這是因為較大的圍壓提升了接觸顆粒的變形能，即顆粒儲備更高的抵抗破壞的能量，因此，在相同的能量輸入條件下，其內(nèi)部顆粒接觸產(chǎn)生破壞的總量減少，對應圖7中概率P的減小，即圍壓遏制了試樣的破壞發(fā)展。

由式(11)可得損傷與應變的關系，進而得到不同圍壓下的D～ε關系曲線，如圖8所示。損傷在經(jīng)過初始為0的階段后隨著應變的增大而增大，破壞時的損傷值也隨著圍壓的增大而增大。但當圍壓為10 MPa時，損傷曲線并未按照0～7 MPa下?lián)p傷的發(fā)展規(guī)律演化，主要是因為在施加高圍壓的同時，會產(chǎn)生部分初始的顆粒破壞損傷，在圖8中表現(xiàn)為損傷曲線的左移，即對應10 MPa圍壓的情況。但是，從圖8還可以發(fā)現(xiàn)，曲線發(fā)展后期的斜率隨著圍壓的增大而減小，說明圍壓的增大依然可以有效地抑制損傷的發(fā)展程度。

圖8 損傷與應變的關系Fig.8 Relationship between damage and strain

4 結論

1)運用數(shù)理分布模型來構建巖石本構關系是一種常見的方法，在以往研究中，分布模型的選用通?；诮?jīng)驗主義，且其分布參數(shù)缺乏明確的物理意義。為了解決上述問題，本文基于統(tǒng)計力學基本理論，認為微元強度破壞所需能量以最概然分布的形式存在，采用經(jīng)過嚴格推導論證的Maxwell分布構建了砂巖本構理論。通過室內(nèi)試驗對比分析，驗證了模型的適用性，結果表明，該本構模型能很好地反映應力-應變曲線的各個階段。

2)模型參數(shù)α和β具有明確的物理意義。從微觀角度上看，α反映微觀顆粒膠結鍵的方向，β表示試樣壓縮性能的強弱；從宏觀角度上看，α和β反映巖石的宏觀強度和延展性，且β對峰值強度的影響較α更強。

3)參數(shù)α與圍壓存在良好的線性遞減關系，而參數(shù)β與圍壓呈二次遞增關系。