張 超，杜翠鳳?，宋衛(wèi)東，付建新

1) 北京科技大學(xué)土木與資源工程學(xué)院，北京 100083 2) 金屬礦山高效開采與安全教育部重點(diǎn)實(shí)驗(yàn)室，北京 100083

在礦山開采不斷向深部進(jìn)行過(guò)程中，高應(yīng)力與高水壓的耦合作用愈發(fā)嚴(yán)重，復(fù)雜環(huán)境加劇了開采的難度；另一方面應(yīng)力場(chǎng)與滲流場(chǎng)的相互作用和相互影響對(duì)巖體產(chǎn)生變形破壞，影響巖體的整體穩(wěn)定性. 因此分析在應(yīng)力–滲流耦合作用下巖石的變化特征已成為當(dāng)前巖石力學(xué)領(lǐng)域的重要課題[1-3].

國(guó)內(nèi)外學(xué)者對(duì)此的研究主要集中在室內(nèi)試驗(yàn)、理論分析和數(shù)值模擬方面，并取得了較多卓有成效的成果. 在室內(nèi)試驗(yàn)方面，張培森等[4-6]、張俊文等[7]通過(guò)控制圍壓和孔隙水壓，研究了不同條件下砂巖的強(qiáng)度特性、滲透特性和不同損傷程度下的滲透率演化規(guī)律；李俊平等[8-9]借助單軸壓縮試驗(yàn)研究了4 種不同類型巖石在滲流與應(yīng)力耦合下的聲發(fā)射特征；胡少華等[10]采用氣體瞬態(tài)壓力脈沖法對(duì)不同圍壓下花崗巖試樣在三軸壓縮過(guò)程中的滲透率變化進(jìn)行了試驗(yàn)研究，得到花崗巖滲透率變化與損傷演化趨勢(shì)整體上具有一致性；王如賓等[11]設(shè)計(jì)三軸卸荷滲流力學(xué)試驗(yàn)，研究了初始圍壓和卸荷速率對(duì)巖石卸荷變形破壞過(guò)程中的滲透性演化規(guī)律. 在理論分析方面，Zhou 等[12]構(gòu)建了巖石裂隙–孔隙耦合滲流演化模型和裂隙–孔隙耦合損傷蠕變模型，定量表征了巖石的蠕變過(guò)程；Liu 等[13]引入損傷變量，建立了彈塑性損傷本構(gòu)模型，并通過(guò)滲流理論推導(dǎo)出巖石微觀和宏觀應(yīng)力滲流耦合特征的滲透率演化方程；Yang 等[14]提出了考慮滲流與應(yīng)力耦合效應(yīng)的滲流–應(yīng)力交叉耦合各向異性模型，初步研究了裂隙巖體的等效滲透率和損傷張量；王偉等[15]在損傷力學(xué)理論基礎(chǔ)上，構(gòu)建考慮孔隙水壓力的巖石統(tǒng)計(jì)損傷本構(gòu)模型，較好地反映巖石在應(yīng)力滲流耦合作用下應(yīng)力–應(yīng)變過(guò)程. 在數(shù)值模擬方面，顏丙乾等[16]在總結(jié)國(guó)內(nèi)外多場(chǎng)耦合數(shù)值模擬軟件及耦合計(jì)算程序的開發(fā)等方面的研究現(xiàn)狀的基礎(chǔ)上，展望多場(chǎng)多相耦合作用下巖石力學(xué)數(shù)值分析的研究方向；Zhang 等[17]利用FLAC3D–CFD 仿真器，采用非達(dá)西模型模擬了流體控制系統(tǒng)的演化過(guò)程，對(duì)地下礦山開采應(yīng)力–滲流–裂隙場(chǎng)耦合效應(yīng)做出準(zhǔn)確動(dòng)態(tài)分析；Liu 等[18]、Zeng 等[19]基于顆粒流方法建立了流體耦合DEM 模型，深入研究了不同圍壓下飽和巖石的三軸壓縮試驗(yàn). 由于巖石自身屬性差異，微裂隙的不同造成了巖石孔隙率的多樣性，進(jìn)而對(duì)巖石宏觀強(qiáng)度產(chǎn)生差異性影響，學(xué)者們?cè)谶M(jìn)行巖石應(yīng)力–滲流耦合試驗(yàn)和理論推導(dǎo)過(guò)程中孔隙率所產(chǎn)生的影響還有待于深入研究；對(duì)應(yīng)力–應(yīng)變曲線孔隙壓密階段的準(zhǔn)確研究有助于實(shí)現(xiàn)本構(gòu)模型的精確化，而在該部分還少有研究.

另外，在進(jìn)行巖石的損傷破壞分析中，采用能量演化特征也有利于從本質(zhì)上揭示巖石的破壞機(jī)制[20-22]. 本文在前人研究基礎(chǔ)上，設(shè)計(jì)花崗巖應(yīng)力–滲流耦合試驗(yàn)，改變圍壓和孔隙水壓，分析耦合條件下深部花崗巖的力學(xué)特征；引入孔隙率參數(shù)，推導(dǎo)巖石應(yīng)力–滲流耦合的本構(gòu)模型并進(jìn)行驗(yàn)證；結(jié)合試驗(yàn)結(jié)果對(duì)巖石破壞過(guò)程的能量耗散進(jìn)行分析，深入研究花崗巖耦合加載下的破壞規(guī)律. 研究結(jié)果可為地下礦山深部開采過(guò)程巖石的變形特征提供理論指導(dǎo).

1 試驗(yàn)研究條件及試驗(yàn)方案

1.1 試驗(yàn)研究條件

本次試驗(yàn)選取的巖塊均為山東黃金集團(tuán)某金礦-400 m 水平礦體下盤圍巖，巖性主要為碎裂花崗巖、鉀化花崗巖和絹英巖化花崗巖. 進(jìn)行試樣成分化驗(yàn)，主要成分包括斜長(zhǎng)石(質(zhì)量分?jǐn)?shù)45%～50%)、鉀長(zhǎng)石(15%～25%)、石英(15%～30%)、黑云母(15%～25%)等. 選取的巖塊對(duì)描述礦山巖性具有代表性，且表面無(wú)明顯節(jié)理. 巖樣加工為《水利水電工程巖石試驗(yàn)規(guī)程》(SL264—2020)要求的標(biāo)準(zhǔn)尺寸(直徑50 mm×高度100 mm)，上、下端的不平整度均小于0.02 mm.

為保證試驗(yàn)的準(zhǔn)確性，試樣首先進(jìn)行真空飽水，借助蘇州紐邁分析儀器股份有限公司的Meso-MR23-060H-I 核磁共振成像分析儀(圖1)進(jìn)行低場(chǎng)核磁共振(NMR)試驗(yàn)，測(cè)量試樣初始孔隙度，探究試驗(yàn)前試樣孔隙賦存規(guī)律. 采用北京科技大學(xué)Rock 600-50 VHT 型巖石高溫三軸流變?cè)囼?yàn)系統(tǒng)(圖2)進(jìn)行試件應(yīng)力–滲流耦合試驗(yàn)，該儀器可控制最大圍壓60 MPa，最大軸壓500 MPa，最大孔隙水壓50 MPa，精度可控制在±0.01 MPa，可實(shí)時(shí)記錄應(yīng)力、應(yīng)變及流量變形.

圖1 核磁共振成像分析儀Fig.1 Magnetic resonance imaging analyzer

圖2 巖石三軸流變?cè)囼?yàn)系統(tǒng)Fig.2 Triaxial rheological test system for rock

1.2 試驗(yàn)方案設(shè)計(jì)

共進(jìn)行8 組試驗(yàn)，分別改變圍壓、孔隙水壓，考慮在無(wú)孔隙水壓條件下、圍壓恒定改變孔隙水壓條件下、孔隙水壓恒定改變圍壓條件下花崗巖力學(xué)特性及破壞規(guī)律，每組選取3 個(gè)試件進(jìn)行相同試驗(yàn)，選取較好的試驗(yàn)結(jié)果進(jìn)行分析. 為盡量減少試樣之間的差異性，確保試驗(yàn)結(jié)果的可靠性，分別測(cè)試各試件物理性質(zhì)，采用NM-4B 非金屬超聲檢測(cè)分析儀篩選出縱波波速在3850±1% m·s-1范圍內(nèi)的試樣，并測(cè)定其密度為2655±1.5% kg·m-3.

試驗(yàn)研究環(huán)境為高地應(yīng)力、高滲透壓條件，根據(jù)礦區(qū)的實(shí)際應(yīng)力測(cè)試結(jié)果(式(1))，埋深400～1000 m 范圍內(nèi)巖體初始地應(yīng)力σz在10～30 MPa 之間，因此本文選取試驗(yàn)圍壓分別為10、20、30 MPa.

式中：H為深度，m；σhmax為水平最大主應(yīng)力，MPa；σhmin為水平最小主應(yīng)力，MPa；σz為垂直應(yīng)力，MPa.

圍壓為10 MPa 時(shí)，改變孔隙水壓為0、6 MPa；圍壓為20 MPa 時(shí)，改變孔隙水壓為0、3、6、9 MPa；圍壓為30 MPa 時(shí)，改變孔隙水壓為0、6 MPa. 本試驗(yàn)中首先控制進(jìn)、出水端孔隙水壓均為預(yù)設(shè)值，待穩(wěn)定后降低出水端孔隙水壓與大氣連通，形成最終滲透壓差. 圍壓加載過(guò)程中控制加載速率為2 MPa·min-1，軸壓加載過(guò)程中控制加載速率為3 MPa·min-1. 采用穩(wěn)態(tài)法進(jìn)行試驗(yàn)中滲透率的計(jì)算，在耦合過(guò)程中滲流符合達(dá)西定律，具體計(jì)算公式見式(2).

式中：Ki為試樣在Δti時(shí)間內(nèi)的平均滲透率，m2；μ為流體黏滯系數(shù)，Pa·s；L為試樣高度，m；ΔQi為Δti時(shí)間內(nèi)滲過(guò)試樣的水流體積，m3；A為試樣橫截面積，m2；ΔP為試樣上下游滲透壓差，ΔP=P1-P2，Pa；Δti為記錄間隔時(shí)間.

2 試驗(yàn)結(jié)果分析

2.1 NMR 試驗(yàn)結(jié)果分析

低場(chǎng)核磁共振是利用孔隙內(nèi)部液相水的橫向弛豫時(shí)間T2值反映該部位水分所處的孔隙結(jié)構(gòu).橫向弛豫時(shí)間T2是描述氫原子被磁化后，橫向磁化矢量的衰減快慢. 當(dāng)橫向磁化矢量由最大值衰減至最大值的37%時(shí)所需要的時(shí)間，被定義為T2.T2的計(jì)算公式為[23-24]，

式中：T2為橫向弛豫時(shí)間，ms；S為孔隙表面積，μm2；V為孔隙流體體積，μm3；ρ2為表面弛豫強(qiáng)度，μm·ms-1. 假設(shè)所有的孔隙均為單一形狀，則公式(3)可改寫為，

式中：Fs為孔隙形狀因子，球狀孔隙時(shí)，F(xiàn)s為3，柱狀孔隙時(shí)，F(xiàn)s為2；R為孔隙半徑，μm.

根據(jù)公式(4)可知，T2與孔隙半徑R成正比，T2越大則孔徑越大，因此可以根據(jù)T2來(lái)分析試件的孔隙結(jié)構(gòu)變化情況. 按照公式(4)對(duì)待測(cè)樣品NMR 試驗(yàn)結(jié)果進(jìn)行反算，得到試樣的孔隙率分布圖，如圖3 所示. 參考Lutz 等[23]的孔隙分類標(biāo)準(zhǔn)，依據(jù)孔徑大小將孔隙分為微孔、次級(jí)孔和主干孔，具體對(duì)應(yīng)關(guān)系見表1.

圖3 試樣孔隙率分布圖Fig.3 Porosity distribution of the sample

分別計(jì)算各試樣孔隙率，并考慮孔隙率曲線圖與邊界孔徑所圍成的面積作為各孔隙所占比例，分析結(jié)果如表2 所示. 試樣均為致密花崗巖，孔隙率均小于5%，孔隙以微孔和次級(jí)孔為主，極少出現(xiàn)使試樣產(chǎn)生較大劣化效應(yīng)的主干孔.

表2 孔隙率及孔隙比例計(jì)算Table 2 Calculation of porosity and pore proportion

考慮各孔隙所占比例，進(jìn)行孔隙率線性擬合，擬合公式如下，

式中：n為孔隙率，%；x1、x2、x3分別為微孔、次級(jí)孔、主干孔所占總孔隙比例，%；R2為相關(guān)系數(shù).

繪制試驗(yàn)值與擬合值之間的對(duì)比圖(圖4)，可以看出孔隙度試驗(yàn)值與擬合值之間的殘差大部分都控制在10%范圍內(nèi). 結(jié)合公式內(nèi)系數(shù)值可知，主干孔是造成孔隙率變化的主要原因，直接影響孔隙率大小和試樣強(qiáng)度；孔隙率的大小與次級(jí)孔、主干孔有關(guān)，微孔所占比例很難對(duì)孔隙率產(chǎn)生影響.

圖4 孔隙度結(jié)果的對(duì)比圖Fig.4 Comparison of porosity results

2.2 應(yīng)力–滲流耦合試驗(yàn)結(jié)果分析

依據(jù)試驗(yàn)方案進(jìn)行巖石應(yīng)力–滲流耦合試驗(yàn)，通過(guò)試驗(yàn)數(shù)據(jù)繪制不同條件下試件的全應(yīng)力–應(yīng)變曲線和滲透率曲線，如圖5 所示(圖中σ1為軸壓，MPa；σ3為圍壓，MPa；σp為孔隙水壓，MPa；ε1為軸向應(yīng)變；ε3為環(huán)向應(yīng)變). 可以看出，改變圍壓和孔隙水壓，巖石的應(yīng)力–應(yīng)變曲線均按照特定的形態(tài)進(jìn)行演化. 三軸壓縮試驗(yàn)中圍壓的施加使得試件在環(huán)向產(chǎn)生預(yù)緊力，抑制了試件環(huán)向的擴(kuò)散膨脹，峰值強(qiáng)度有了顯著提升；而隨著圍壓的升高峰值強(qiáng)度逐漸增大，主要原因是圍壓的增大使得環(huán)向的抑制作用增強(qiáng)，進(jìn)而使得軸向承載力逐漸加大.應(yīng)力–滲流耦合試驗(yàn)下，孔隙水壓一定時(shí)隨著圍壓的升高峰值強(qiáng)度逐漸增大，圍壓一定時(shí)隨著孔隙水壓的升高峰值強(qiáng)度逐漸減小，主要原因是孔隙水壓的存在會(huì)弱化巖石的力學(xué)性能，巖石內(nèi)部裂紋受孔隙水壓的影響擴(kuò)展逐漸加快，且隨著孔隙水壓增大弱化效果加劇.

圖5 全應(yīng)力–應(yīng)變曲線及滲透率曲線. (a) σ3 = 10 MPa, σp = 0 MPa; (b) σ3 = 10 MPa, σp = 6 MPa; (c) σ3 = 20 MPa, σp = 0 MPa; (d) σ3 = 20 MPa, σp =3 MPa; (e) σ3 = 20 MPa, σp = 6 MPa; (f) σ3 = 20 MPa, σp = 9 MPa; (g) σ3 = 30 MPa, σp = 0 MPa; (h) σ3 = 30 MPa, σp = 6 MPaFig.5 Stress–strain and permeability curves: (a) σ3 = 10 MPa, σp = 0 MPa; (b) σ3 = 10 MPa, σp = 6 MPa; (c) σ3 = 20 MPa, σp = 0 MPa; (d) σ3 = 20 MPa,σp = 3 MPa; (e) σ3 = 20 MPa, σp = 6 MPa; (f) σ3 = 20 MPa, σp = 9 MPa; (g) σ3 = 30 MPa, σp = 0 MPa; (h) σ3 = 30 MPa, σp = 6 MPa

巖石的滲透率變化規(guī)律與巖石破壞規(guī)律類似，分別以裂紋閉合應(yīng)力σA、裂紋起裂應(yīng)力σB、峰值應(yīng)力σC為界點(diǎn)分為四階段：孔隙壓密階段(OA)，由于巖樣內(nèi)部存在原生的節(jié)理裂隙，巖樣具有一定的滲透性，在孔隙水壓的作用下，水流很快流入節(jié)理裂隙中，隨著軸力的增大，原生的節(jié)理裂隙被壓密，滲流喉道受阻，且隨著壓力的增大，滲流困難程度加大，導(dǎo)致滲透率降低，該階段與壓密階段出現(xiàn)的拐點(diǎn)即為裂紋閉合應(yīng)力σA，在三軸試驗(yàn)過(guò)程中，由于圍壓的作用使得應(yīng)力–應(yīng)變曲線在壓密階段并未出現(xiàn)明顯的上凹現(xiàn)象，因此對(duì)于A點(diǎn)位置的判斷主要是依據(jù)滲透性變化趨勢(shì)，當(dāng)滲透率下降趨勢(shì)變緩即認(rèn)為進(jìn)入彈性變形階段；在彈性變形階段(AB)，滲透率變化趨勢(shì)基本保持穩(wěn)定，產(chǎn)生較小范圍的波動(dòng)，并未出現(xiàn)大幅度的增加或降低現(xiàn)象；塑形強(qiáng)化變形階段(BC)，在該階段巖石發(fā)生塑性變形，由彈性變形過(guò)渡至塑性變形的拐點(diǎn)稱為裂紋起裂應(yīng)力σB，該階段促使試樣發(fā)生塑性破壞，微裂紋擴(kuò)張并匯聚，逐漸形成宏觀裂紋并產(chǎn)生非穩(wěn)定擴(kuò)展，使得滲流喉道增多，滲透率增加；破壞后階段(CD)，試樣達(dá)到其抗壓極限即峰值應(yīng)力σC而發(fā)生破壞，產(chǎn)生宏觀裂隙，在破壞階段滲透率迅速增大，并且破壞后滲透率持續(xù)增大，最大值出現(xiàn)在峰后；在峰值點(diǎn)后由于軸壓的繼續(xù)施加而使得裂紋壓密，導(dǎo)致滲透率在達(dá)到最大值后又急劇下降，后續(xù)達(dá)到新的穩(wěn)定狀態(tài). 整體上滲透率的變化過(guò)程為：緩慢下降(OA)→穩(wěn)定發(fā)展(AB)→緩慢上升(BC)→加速上升后下降(CD). 圖6 為試樣力學(xué)參數(shù)變化規(guī)律圖. 可以看出，孔隙水壓增大后試樣的峰值強(qiáng)度、峰值應(yīng)變(峰值強(qiáng)度所對(duì)應(yīng)的軸向應(yīng)變)均呈線性減小，且隨著圍壓的增大逐漸提升，下降趨勢(shì)愈發(fā)明顯；圍壓增大后試樣的峰值強(qiáng)度、峰值應(yīng)變均呈線性增加，且提升速率隨著孔隙水壓的增大而逐漸減小. 整體上圍壓的增大對(duì)花崗巖力學(xué)參數(shù)的顯著性影響逐漸增強(qiáng)，而孔隙水壓的增大對(duì)花崗巖力學(xué)參數(shù)的顯著性影響逐漸降低.

圖6 力學(xué)參數(shù)變化規(guī)律. (a)孔隙水壓–峰值強(qiáng)度變化規(guī)律；(b)孔隙水壓–峰值應(yīng)變變化規(guī)律Fig.6 Variation law of mechanical parameters: (a) pore pressure–peak strength variation law; (b) pore pressure–peak strain variation law

以線性擬合為基礎(chǔ)，考慮圍壓和孔隙水壓的交互作用，繪制峰值強(qiáng)度、峰值應(yīng)變曲面關(guān)系，如圖7 所示，并擬合其曲面關(guān)系式，見式(6)～(7). 相關(guān)系數(shù)分別可達(dá)0.996、0.970，表明公式具有較高的擬合性，因此可認(rèn)為峰值強(qiáng)度、峰值應(yīng)變采用線性表達(dá)式擬合具有較好的準(zhǔn)確性.

圖7 力學(xué)參數(shù)三維曲面關(guān)系圖. (a) 圍壓–孔隙水壓–峰值強(qiáng)度關(guān)系；(b) 圍壓–孔隙水壓–峰值應(yīng)變關(guān)系Fig.7 Three-dimensional surface relation diagram of mechanical parameters: (a) confining pressure–pore pressure–peak strength relationship;(b) confining pressure–pore pressure–peak strain relationship

分別考慮圍壓和孔隙水壓改變對(duì)峰值滲透率(峰值強(qiáng)度所對(duì)應(yīng)的滲透率)影響，分析結(jié)果如圖8所示. 隨著圍壓的增大，峰值滲透率呈線性降低，平均下降速率為0.09×10-18m2·MPa-1，主要原因是由于圍壓的增大使得巖樣本身所受的環(huán)向束縛增大，內(nèi)部的微裂紋擴(kuò)展受到抑制，滲流的喉道變窄，滲透率進(jìn)而降低. 隨著孔隙水壓的增大，峰值滲透率呈線性增大，平均上升速率為0.07×10-18m2·MPa-1，孔隙水壓的增大使得巖樣本身所受的軸向壓力增大，與軸壓共同作用下加快了巖樣的破壞，內(nèi)部滲流喉道擴(kuò)展加快，滲透率增大.

圖8 滲透率變化規(guī)律分析Fig.8 Analysis of the permeability variation law

根據(jù)圖8 分析可知花崗巖試樣的滲透率與圍壓和孔隙水壓密切相關(guān)，但兩種因素對(duì)試樣滲透率的影響各有差異，即滲透率對(duì)二者的敏感度不同. 對(duì)圍壓和孔隙水壓進(jìn)行敏感度分析，當(dāng)孔隙水壓為6 MPa 時(shí)，圍壓自10 MPa 增加至30 MPa 過(guò)程中，花崗巖的峰值滲透率分別較小約32.5%、59.2%，平均減小約41.6%；當(dāng)圍壓為20 MPa 時(shí)，孔隙水壓自3 MPa 增加至9 MPa 過(guò)程中，花崗巖的峰值滲透率分別增大約5.7%、18.5%，平均增大約12.1%.對(duì)比分析表明，圍壓對(duì)花崗巖試樣的峰值滲透率影響更大，即峰值滲透率對(duì)圍壓的敏感更高.

3 討論

3.1 應(yīng)力–滲流耦合本構(gòu)模型

考慮巖石內(nèi)部孔隙在初始?jí)好茈A段的非線性變形以及巖石的破壞特征，進(jìn)行巖石應(yīng)力–滲流耦合損傷本構(gòu)模型構(gòu)建. 選取柱狀巖石進(jìn)行變形分析，如圖9 所示. 參考巖石二相體模型及曹文貴[25-27]等人對(duì)巖石損傷模擬的研究，柱狀巖石整體由固體骨架和孔隙兩部分組成，假定巖石中孔隙全部集中于一邊，加載前長(zhǎng)度為；固體骨架集中于另一邊，加載前長(zhǎng)度為；柱狀巖石總長(zhǎng)度為l0. 圖中，E1、μ1為固體骨架部分的彈性模量和泊松比，E2、μ2為孔隙部分的平均彈性模量和泊松比；固體骨架壓縮后長(zhǎng)度為ls，產(chǎn)生變形量為 ?ls；孔隙壓縮后長(zhǎng)度為lp，產(chǎn)生變形量為 ?lp.

圖9 巖石變形分析示意圖Fig.9 Schematic of rock deformation analysis

分別對(duì)固體骨架和孔隙進(jìn)行應(yīng)變分析，可表示為，

式中，n0為孔隙部分長(zhǎng)度與總長(zhǎng)度之比，即n0=/l0，本文認(rèn)為n0值與巖石初始孔隙度n相等，其值為NMR 試驗(yàn)結(jié)果.

(1) 孔隙部分本構(gòu)方程構(gòu)建.

巖石初期的變形為孔隙壓密階段，主要是由于孔隙的壓縮而形成的非線性變形. 在巖石進(jìn)行三軸試驗(yàn)中，孔隙的變形具有不可恢復(fù)的特點(diǎn)，因此根據(jù)材料變形力學(xué)分析，采用真應(yīng)變方法進(jìn)行孔隙部分變形分析[26]，即，

式中，為孔隙部分在名義應(yīng)力作用后總長(zhǎng)度.由于采用室內(nèi)試驗(yàn)無(wú)法測(cè)得孔隙部分變形長(zhǎng)度，故利用非線性變形特點(diǎn)，采用分級(jí)加載方法進(jìn)行計(jì)算，將孔隙應(yīng)變近似為由n級(jí)應(yīng)力增量引起的應(yīng)變?cè)隽康睦奂樱傻茫?/p>

式中，和分別為孔隙部分在第m級(jí)應(yīng)力增量作用下的彈性模量和泊松比；、、為三個(gè)主應(yīng)力方向在第m級(jí)的應(yīng)力增量(i,j,k=1, 2,3). 參考張超等[26]的研究，假設(shè)孔隙部分各級(jí)增量作用下彈性模量和泊松比保持不變，則聯(lián)立式(9)、式(12)得，

式中， σij為施加應(yīng)力；為名義應(yīng)力； σp為孔隙水壓； δij為單位二階張量值，當(dāng)i=j時(shí)，δij=1、當(dāng)i≠j時(shí)，δij=0.

以Lemaitre 應(yīng)變等價(jià)原理為基礎(chǔ)，在圍巖三軸壓縮試驗(yàn)中整理可得孔隙部分本構(gòu)方程為，

(2) 固體骨架本構(gòu)方程構(gòu)建.

在固體骨架中引入損傷變量D，即已破壞微元體數(shù)目與總微元體數(shù)目之比. 考慮損傷變量進(jìn)行巖石有效應(yīng)力計(jì)算，主要為，

式中，w，k為韋伯分布參數(shù)；λ為巖石微裂強(qiáng)度. 在外荷載作用下，損傷變量即微元破壞概率，可認(rèn)為損傷變量為，

巖石服從廣義虎克定律，聯(lián)立式(16)、式(17)、式(18)，得到巖石固體骨架部分本構(gòu)方程為，

由于巖石的抗拉強(qiáng)度遠(yuǎn)小于抗壓強(qiáng)度，為簡(jiǎn)化計(jì)算模型，則假設(shè)巖石固體骨架單元強(qiáng)度滿足最大拉伸應(yīng)變屈服準(zhǔn)則，即，

式中， εw為固體骨架最大拉應(yīng)變. 將式(20)代入式(10)中得，

(3) 應(yīng)力–滲流耦合本構(gòu)模型分析.

將式(15)、式(19)代入式(10)中，得巖石整體的應(yīng)力–滲流耦合損傷本構(gòu)模型為，

模型中的待定參數(shù)主要有E1、μ1、E2、μ2、w、k. 其中，E1、μ1為應(yīng)力–應(yīng)變曲線線性變形階段的主要參數(shù)，可根據(jù)線性變形部分試驗(yàn)數(shù)據(jù)求得.E2、μ2為應(yīng)力–應(yīng)變曲線孔隙壓密階段的主要參數(shù)，選取壓密階段的試驗(yàn)數(shù)據(jù)，通過(guò)最小二乘法原理擬合曲線求得.

韋伯分布參數(shù)w和k，可采用二次對(duì)數(shù)運(yùn)算法計(jì)算求得. 需要通過(guò)若干試驗(yàn)數(shù)據(jù)進(jìn)行計(jì)算，為保持各曲線計(jì)算的統(tǒng)一，各試驗(yàn)數(shù)據(jù)均選取峰值應(yīng)力點(diǎn)和屈服應(yīng)力點(diǎn)進(jìn)行w、k的求解.

(4) 本構(gòu)模型驗(yàn)證.

基于上述本構(gòu)模型，分別進(jìn)行不同圍壓、孔隙水壓條件下的理論模型計(jì)算，最終確定模型的各參數(shù)如表3 所示.

表3 本構(gòu)模型參數(shù)計(jì)算結(jié)果Table 3 Calculation results of constitutive model parameters

綜合分析各理論應(yīng)力–應(yīng)變曲線與實(shí)驗(yàn)應(yīng)力–應(yīng)變曲線（圖10），對(duì)比分析可知，各理論曲線均與實(shí)驗(yàn)曲線有較高的擬合度，準(zhǔn)確把控壓密階段是提高理論曲線精確性的關(guān)鍵，本文采用巖石初始非線性變形更好的表征了巖石的壓密階段，并引入孔隙率和孔隙水壓作為損傷本構(gòu)模型的基本參數(shù)，準(zhǔn)確的結(jié)合了花崗巖試樣的基本特性，最終使得本構(gòu)模型具有較高的普適性.

圖10 理論和試驗(yàn)應(yīng)力–應(yīng)變曲線比較. (a) σp = 0 MPa; (b) σp = 3 MPa; (c) σp = 6 MPa; (d) σp = 9 MPaFig.10 Comparison of theoretical and experimental stress–strain curves: (a) σp = 0 MPa; (b) σp = 3 MPa; (c) σp = 6 MPa; (d) σp = 9 MPa

3.2 巖石破壞能量耗散分析

從能量的角度來(lái)看，巖石變形破壞的實(shí)質(zhì)上是能量積蓄、耗散和釋放的全過(guò)程. 研究巖石破壞過(guò)程中的能量演化特征并建立其與損傷發(fā)育之間的關(guān)系，有利于從本質(zhì)上揭示巖石的破壞機(jī)制. 在進(jìn)行試驗(yàn)過(guò)程中，試樣在到達(dá)峰值強(qiáng)度之前不斷吸收設(shè)備釋放的能量，而在峰后不斷釋放彈性能量. 巖石的耗散能是試樣變形失穩(wěn)破壞的本質(zhì)屬性，它反映了受載巖石內(nèi)部微裂隙的發(fā)展演化、巖石強(qiáng)度不斷弱化并最終喪失的過(guò)程，試樣的能量主要分為耗散能和彈性能，耗散能是不可逆的，而彈性能是雙向可逆的. 依據(jù)熱力學(xué)第一定理，

式中，U為總應(yīng)變能；Ue為彈性應(yīng)變能，Ud為耗散能. 應(yīng)力–滲流耦合三軸試驗(yàn)中，彈性能的計(jì)算公式為，

考慮三軸壓縮試驗(yàn)中的做功轉(zhuǎn)化，軸向主要產(chǎn)生壓縮變形，σ1主要對(duì)試樣做正功；而環(huán)向主要產(chǎn)生膨脹變形，σ3主要對(duì)試樣做負(fù)功. 總做功轉(zhuǎn)化為應(yīng)變能為，

式中，U1為σ1對(duì)試樣做正功轉(zhuǎn)化的應(yīng)變能；U2、U3為σ2、σ3對(duì)試樣做負(fù)功轉(zhuǎn)化的應(yīng)變能；U0為試樣本身所儲(chǔ)存的應(yīng)變能，可根據(jù)彈性力學(xué)公式直接求出. 各應(yīng)變能具體計(jì)算公式為，

綜合上式，可得到耗散能計(jì)算公式為，

根據(jù)三軸試驗(yàn)結(jié)果，依據(jù)上述應(yīng)變能計(jì)算方法，求得各試驗(yàn)方案下能量變形曲線，并繪制能量曲線如圖11 所示(由于文章篇幅有限，僅選取σ3=10 MPa 試驗(yàn)結(jié)果進(jìn)行分析). 與巖樣的滲透率變化規(guī)律相對(duì)應(yīng)，在壓密階段，能量較小且增長(zhǎng)率較緩，試樣內(nèi)初始裂隙壓縮，該階段能量的轉(zhuǎn)化效率較低；在彈性階段，巖石主要發(fā)生彈性變形，隨著荷載增大總能量和彈性能均在不斷升高且增長(zhǎng)趨勢(shì)一致，在巖石未達(dá)到破壞之前耗散能均呈現(xiàn)為比較低的平緩趨勢(shì)；在塑性階段，巖石內(nèi)部損傷和塑性變形增大，耗散能增長(zhǎng)速率大幅度提升，彈性能持續(xù)緩慢增加；在破壞及破壞后階段，巖石宏觀裂紋擴(kuò)展顯著，結(jié)構(gòu)發(fā)生較大改變，在峰值點(diǎn)為彈性能極值點(diǎn)，峰值點(diǎn)之后彈性能迅速轉(zhuǎn)化為用于巖石損傷的耗散能，在該階段巖石耗散能急劇增大，最終趨近于總能量，而彈性能也在該階段降低最終趨近于零.

圖11 不同加載條件下巖石能量參數(shù)關(guān)系曲線. (a) σ3 = 10 MPa, σp = 0 MPa; (b) σ3 = 10 MPa, σp = 6 MPaFig.11 Relationship curve of rock energy parameters under different loading conditions: (a) σ3 = 10 MPa, σp = 0 MPa; (b) σ3 = 10 MPa, σp = 6 MPa

定義巖石峰值強(qiáng)度位置所對(duì)應(yīng)的總能量、彈性能、耗散能分別為破壞總能量UP、極限彈性能UeP、破壞耗散能UdP. 進(jìn)行能量變化規(guī)律分析如圖12所示，分析各應(yīng)變能均表現(xiàn)出明顯的圍壓效應(yīng)和孔隙水壓效應(yīng). 極限彈性能、破壞耗散能、破壞總能量均隨圍壓的增大呈線性增長(zhǎng)，隨孔隙水壓的增大呈線性減小. 隨著圍壓增大，巖石的峰值強(qiáng)度增加，巖石發(fā)生破壞過(guò)程中可吸收能量增大，需要聚集更多的能量才能發(fā)生破壞，使得破壞時(shí)的極限彈性能增大；而增大孔隙水壓，加快了滲流喉道的形成，破環(huán)了巖石原有的完整性，對(duì)巖石的強(qiáng)度產(chǎn)生弱化效應(yīng)，使巖石發(fā)生破壞的可吸收能量降低，相比較于完整巖石在較低軸壓下就可使巖石發(fā)生破裂失穩(wěn).

圖12 巖石破壞能量變化規(guī)律分析. (a)破壞總能量; (b)極限彈性能; (c)破壞耗散能Fig.12 Analysis of change law of rock failure energy: (a) total failure energy; (b) ultimate elastic energy; (c) failure dissipation energy

由于圍壓、孔隙水壓分別與應(yīng)變能呈線性關(guān)系，因此以兩者相互關(guān)系為基礎(chǔ)分別構(gòu)建總能量、彈性能、耗散能曲面擬合關(guān)系式，可以看出擬合關(guān)系式均具有較高的準(zhǔn)確性.

定義耗能比為各具體變形條件下巖石耗散能與總應(yīng)變能的比值(式(30))，繪制在不同條件下耗能比的變化曲線，如圖13 所示.

圖13 巖石破壞耗能比變化規(guī)律分析. (a) σp = 0 MPa; (b) σp = 6 MPa; (c) σ3 = 20 MPaFig.13 Energy dissipation ratio change law analysis of rock failure: (a) σp = 0 MPa; (b) σp = 6 MPa; (c) σ3 = 20 MPa

整體上耗能比隨著應(yīng)變的增加均呈現(xiàn)增大→減小→增加的“S”型變化規(guī)律. 前期耗能比的增大階段對(duì)應(yīng)力–應(yīng)變曲線中的壓密階段，在該階段巖石內(nèi)部初始微裂紋閉合，耗散能增大，耗能比對(duì)應(yīng)增大；在中期的減小階段對(duì)應(yīng)巖石彈性階段，在該階段巖石的內(nèi)部能量主要轉(zhuǎn)化彈性能，耗能比逐漸減小；而在后期的增大階段主要是由于巖石屈服階段裂隙擴(kuò)張加快，耗散能聚積增大，耗能比提升.

在同一孔隙水壓條件下(圖13(a)、圖13(b))，耗能比隨著圍壓的增加呈現(xiàn)前期(孔隙壓密階段)增大、后期(塑性變形階段)減小的變化規(guī)律，前期耗能比增大主要是由于高圍壓作用下限制了裂隙的變形，裂隙要實(shí)現(xiàn)壓密過(guò)程就需要耗散更多的能量，使得耗散比增大；而后期耗能比主要是由于高圍壓作用下巖石破壞會(huì)儲(chǔ)存更大的彈性應(yīng)變能，經(jīng)過(guò)彈性階段后彈性應(yīng)變能持續(xù)增加，致使耗散比逐漸降低. 在同一圍壓條件下(圖13(c))，耗能比隨著孔隙水壓的增加呈現(xiàn)整體增大趨勢(shì). 主要是由于孔隙水壓的存在增大了內(nèi)部能量向耗散能的轉(zhuǎn)化，這也說(shuō)明了高孔隙水壓致使耗散能提升，從而降低了巖石自身的強(qiáng)度，與試驗(yàn)結(jié)果分析相一致.

通過(guò)分析可知，孔隙水壓的存在對(duì)花崗巖具有明顯的弱化效應(yīng)，在金屬礦深部區(qū)域，地下水的存在是影響巷道圍巖變形破壞的重要原因之一.在礦山深部巷道支護(hù)設(shè)計(jì)中，應(yīng)充分考慮地下水對(duì)圍巖的弱化作用，在耗散能快速增大之前進(jìn)行巷道的有效支護(hù)對(duì)減緩礦山災(zāi)害程度方面具有顯著影響.

4 結(jié)論

采用室內(nèi)試驗(yàn)及理論分析，系統(tǒng)的研究了深部花崗巖應(yīng)力–滲流耦合下的力學(xué)行為及破壞特征，主要的研究結(jié)果有：

(1) 巖石的滲透率變化和能量變化規(guī)律分別以裂紋閉合應(yīng)力、裂紋起裂應(yīng)力、峰值應(yīng)力為界點(diǎn)劃分為壓密階段、彈性階段、塑形階段、破壞后階段，在各階段中滲透率呈現(xiàn)緩慢下降→穩(wěn)定發(fā)展→緩慢上升→加速上升后下降的變化規(guī)律.

(2) 巖石的峰值強(qiáng)度、峰值應(yīng)變隨孔隙水壓的增大呈線性減小且減小速率逐漸提升，隨圍壓的增大呈線性增加且增大速率逐漸變緩. 峰值滲透率隨著孔隙水壓的增大呈線性增大，隨圍壓的增大呈線性減小. 巖石力學(xué)參數(shù)、峰值滲透率對(duì)圍壓的敏感更高.

(3) 引入花崗巖初始孔隙率，采用巖石初始非線性變形構(gòu)建應(yīng)力–滲流耦合本構(gòu)模型，準(zhǔn)確把控壓密階段變形曲線，采用二次對(duì)數(shù)運(yùn)算法進(jìn)行模型參數(shù)獲取，并與試驗(yàn)曲線進(jìn)行對(duì)比，具有較高的擬合度. 該模型準(zhǔn)確的結(jié)合了花崗巖的基本特性，使得模型具有較高的普適性.

(4) 應(yīng)變能均表現(xiàn)出明顯的圍壓效應(yīng)和孔隙水壓效應(yīng)，均隨圍壓的增大呈線性增長(zhǎng)，隨孔隙水壓的增大呈線性減小. 耗能比隨著應(yīng)變的增加均呈現(xiàn)增大→減小→增加的“S”型變化規(guī)律；耗能比隨著圍壓的增加呈現(xiàn)前期增大、后期減小的變化規(guī)律，隨著孔隙水壓的增加呈現(xiàn)整體增大趨勢(shì).