應力–滲流耦合下深部花崗巖力學行為及破壞特征

張 超，杜翠鳳?，宋衛東，付建新

1) 北京科技大學土木與資源工程學院，北京 100083 2) 金屬礦山高效開采與安全教育部重點實驗室，北京 100083

在礦山開采不斷向深部進行過程中，高應力與高水壓的耦合作用愈發嚴重，復雜環境加劇了開采的難度；另一方面應力場與滲流場的相互作用和相互影響對巖體產生變形破壞，影響巖體的整體穩定性. 因此分析在應力–滲流耦合作用下巖石的變化特征已成為當前巖石力學領域的重要課題[1-3].

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

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

1 試驗研究條件及試驗方案

1.1 試驗研究條件

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

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

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

圖2 巖石三軸流變試驗系統Fig.2 Triaxial rheological test system for rock

1.2 試驗方案設計

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

試驗研究環境為高地應力、高滲透壓條件，根據礦區的實際應力測試結果(式(1))，埋深400～1000 m 范圍內巖體初始地應力σz在10～30 MPa 之間，因此本文選取試驗圍壓分別為10、20、30 MPa.

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

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

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

2 試驗結果分析

2.1 NMR 試驗結果分析

低場核磁共振是利用孔隙內部液相水的橫向弛豫時間T2值反映該部位水分所處的孔隙結構.橫向弛豫時間T2是描述氫原子被磁化后，橫向磁化矢量的衰減快慢. 當橫向磁化矢量由最大值衰減至最大值的37%時所需要的時間，被定義為T2.T2的計算公式為[23-24]，

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

式中：Fs為孔隙形狀因子，球狀孔隙時，Fs為3，柱狀孔隙時，Fs為2；R為孔隙半徑，μm.

根據公式(4)可知，T2與孔隙半徑R成正比，T2越大則孔徑越大，因此可以根據T2來分析試件的孔隙結構變化情況. 按照公式(4)對待測樣品NMR 試驗結果進行反算，得到試樣的孔隙率分布圖，如圖3 所示. 參考Lutz 等[23]的孔隙分類標準，依據孔徑大小將孔隙分為微孔、次級孔和主干孔，具體對應關系見表1.

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

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

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

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

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

繪制試驗值與擬合值之間的對比圖(圖4)，可以看出孔隙度試驗值與擬合值之間的殘差大部分都控制在10%范圍內. 結合公式內系數值可知，主干孔是造成孔隙率變化的主要原因，直接影響孔隙率大小和試樣強度；孔隙率的大小與次級孔、主干孔有關，微孔所占比例很難對孔隙率產生影響.

圖4 孔隙度結果的對比圖Fig.4 Comparison of porosity results

2.2 應力–滲流耦合試驗結果分析

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

圖5 全應力–應變曲線及滲透率曲線. (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

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

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

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

圖7 力學參數三維曲面關系圖. (a) 圍壓–孔隙水壓–峰值強度關系；(b) 圍壓–孔隙水壓–峰值應變關系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

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

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

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

3 討論

3.1 應力–滲流耦合本構模型

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

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

分別對固體骨架和孔隙進行應變分析，可表示為，

式中，n0為孔隙部分長度與總長度之比，即n0=/l0，本文認為n0值與巖石初始孔隙度n相等，其值為NMR 試驗結果.

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

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

式中，為孔隙部分在名義應力作用后總長度.由于采用室內試驗無法測得孔隙部分變形長度，故利用非線性變形特點，采用分級加載方法進行計算，將孔隙應變近似為由n級應力增量引起的應變增量的累加，可得，

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

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

以Lemaitre 應變等價原理為基礎，在圍巖三軸壓縮試驗中整理可得孔隙部分本構方程為，

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

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

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

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

由于巖石的抗拉強度遠小于抗壓強度，為簡化計算模型，則假設巖石固體骨架單元強度滿足最大拉伸應變屈服準則，即，

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

(3) 應力–滲流耦合本構模型分析.

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

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

韋伯分布參數w和k，可采用二次對數運算法計算求得. 需要通過若干試驗數據進行計算，為保持各曲線計算的統一，各試驗數據均選取峰值應力點和屈服應力點進行w、k的求解.

(4) 本構模型驗證.

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

表3 本構模型參數計算結果Table 3 Calculation results of constitutive model parameters

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

圖10 理論和試驗應力–應變曲線比較. (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 巖石破壞能量耗散分析

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

式中，U為總應變能；Ue為彈性應變能，Ud為耗散能. 應力–滲流耦合三軸試驗中，彈性能的計算公式為，

考慮三軸壓縮試驗中的做功轉化，軸向主要產生壓縮變形，σ1主要對試樣做正功；而環向主要產生膨脹變形，σ3主要對試樣做負功. 總做功轉化為應變能為，

式中，U1為σ1對試樣做正功轉化的應變能；U2、U3為σ2、σ3對試樣做負功轉化的應變能；U0為試樣本身所儲存的應變能，可根據彈性力學公式直接求出. 各應變能具體計算公式為，

綜合上式，可得到耗散能計算公式為，

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

圖11 不同加載條件下巖石能量參數關系曲線. (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

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

圖12 巖石破壞能量變化規律分析. (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

由于圍壓、孔隙水壓分別與應變能呈線性關系，因此以兩者相互關系為基礎分別構建總能量、彈性能、耗散能曲面擬合關系式，可以看出擬合關系式均具有較高的準確性.

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

圖13 巖石破壞耗能比變化規律分析. (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

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

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

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

4 結論

采用室內試驗及理論分析，系統的研究了深部花崗巖應力–滲流耦合下的力學行為及破壞特征，主要的研究結果有：

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

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

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

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