深部地下硐室與應力場軸變關系及其圍巖損傷破裂分析

劉力源，張 樂，紀洪廣?

1) 北京科技大學土木與資源工程學院，北京 100083 2) 北京科技大學城市地下空間工程北京市重點實驗室，北京 100083

隨著我國經濟與社會的不斷發展，淺部資源日益枯竭，礦產資源和地下空間開發利用的需求日益增長[1-3].理論上，地球內部從地表到地層以下10000 m處均為成礦空間.因此，為滿足人類社會發展的需要，深部地層礦產資源開采是解決資源短缺問題的有效途徑[4-5].深部開采和地下空間工程成為未來發展的必然趨勢.

隨著地層深度增加，地層應力、地溫和孔隙水壓顯著提高，深部圍巖處于“三高一擾動”的復雜地質力學環境，其力學行為與變形機制不同于淺部地層[6-7]，在深部開采中，會有巖爆、片幫等特殊的破壞現象出現[8-10].國內外很多學者針對深部地下硐室開挖和支護進行了理論分析、數值模擬及原位試驗等工作，取得了一系列成果[11-15].值得特別關注的是，深部地下硐室的形狀、大小及其布置對圍巖損傷破裂和穩定性具有重要的影響.在“三高一擾動”的環境中，高地應力是引起地下硐室圍巖變形、損傷和破壞的根本作用力.巖體中地應力絕大多數是以水平應力為主的三向不等壓的空間應力場.于學馥和喬端[16]，以及于學馥[17]提出地下硐室軸比變化對圍巖變形和破壞起重要控制作用的理論，稱軸變論.

在實際工程中，深部硐室開挖形成的圍巖損傷區域（EDZ），損傷區域的產生顯著影響深部地下硐室的圍巖穩定性，該區域的巖石滲透率通常要比原巖高得多[18]，因此對硐室開挖損傷區域的評估是設計深部硐室的一項關鍵任務.Falls和Young[19]通過聲發射與超聲波手段，指出處于高應力條件下巖體，應力重分布與擾動造成的損傷要嚴重于低應力條件下.Hakami[20]、Read等[21]通過彈性分析與彈塑性分析的數值方法對損傷區域進行評估.Chang等[22]研究了一種基于室內試驗的損傷模型并使用Fish將其應用于開挖損傷區的模擬.張小波等[23]對巷道圍巖力學狀態進行極限平衡分析，推導出圍巖損傷破裂半徑及應力場分布的解析表達式.蔡德所等[24]基于三峽工程基巖開挖爆破建立了基巖爆破的損傷本構模型，將其嵌入SHALE程序中，實現了對基巖開挖爆破損傷范圍的數值模擬.楊棟等[25]通過巖體聲波測試，確定梅花山隧道的圍巖損傷范圍，并采用FLAC軟件模擬表明地應力大小對圍巖損傷分布有著顯著影響.

本文針對深部地下硐室軸比、地應力與圍巖穩定性關系，開展了硐室斷面形狀、側壓系數和構造應力場對圍巖損傷破裂的影響.隨著地層深度增大，硐室圍巖損傷破裂程度加劇，提出了臨界側壓系數概念用以指導深部地層硐室設計及穩定性分析.以三山島金礦西嶺礦區擬建的井深為2005 m的副井為工程依托，開展硐室斷面設計與損傷破裂關系數值模擬，以期助力深部地層巷道科學合理布置.

1 深部地下硐室軸比和地應力關系

地下工程中常采用橢圓形硐室斷面，圖1給出了不同角度（β）地應力作用下硐室開挖硐壁應力計算簡圖.根據彈性力學的解，將兩個主應力作用下求得的應力疊加，橢圓形硐室硐壁的應力為：

圖1 不同地應力場作用下橢圓形硐室應力計算簡圖Fig.1 Stress calculation diagram of an oval chamber under various in situ stress fields

其中，σθ為極坐標下圍巖切向應力；σr為極坐標下圍巖徑向應力；τrθ為極坐標下圍巖剪應力；Z為x軸正半軸a與y軸上半軸b的比值，Z=a/b；θ為硐壁上任意一點到橢圓形硐室中心點連線與x軸的夾角；β為側向主應力作用線與x軸的夾角，β≠0時，橢圓形硐室為構造應力狀態；λ為地應力側壓系數；σv為豎向主應力，β=0時，σv為上覆巖層自重應力.

針對β=0時，橢圓形硐室豎向軸端點（θ=π/2和3π/2）切向應力σθt和水平軸端點（θ=0和π）切向應力σθ1的變化相對應，其中一對為最大值，另一對為最小值.

需要特別指出的是，存在圍巖應力均勻分布的硐室軸比，于學馥先生稱之為“等應力軸比”，也稱“諧硐”.等應力軸比具有以下特點，圍巖應力均勻分布，不出現拉應力，出現的最大應力值最小，即由于開挖引起的應力集中最小.等應力軸比與地應力測壓系數λ的關系，經計算求得Z=λ.

將β=0和Z=λ代入式（1），即得：

由式（4）可知，等應力軸比條件下硐壁的切向應力與θ無關，并且在λ≠1的條件下為均勻壓應力，且其應力值小于圓形硐室λ=1時硐壁切向應力值.因此，等應力軸比規律確定的橢圓形硐室斷面是最穩定的幾何形狀.圖2給出了不同側壓系數下等應力軸比硐室圍巖應力分布情況.

圖2 等應力軸比硐室應力集中系數變化Fig.2 Stress concentration factor for a chamber where the axial ratio equals the stress ratio

2 圍巖損傷破裂數值模擬

采用彈性損傷力學理論對不同工況下的深部地下硐室圍巖損傷破裂進行數值模擬，以探究硐室應力場軸變與圍巖損傷破裂的對應關系.硐室埋深2000 m，豎向主應力為50 MPa，側壓系數λ為1.5，即水平主應力為75 MPa.地層所處的花崗巖單軸抗壓強度為132.3 MPa，拉伸強度為12.4 MPa，彈性模量為42.7 GPa.借助于Weibull分布函數生成深部地層非均質巖體力學參數，圍巖非均質系數取8，如圖3所示，硐室開挖數值模型邊長為30 m，開挖斷面面積為20 m2，模型上邊界為豎向應力σv，模型右邊界為最小水平主應力σh，模型的下邊界和左邊界為滾動約束.采用最大拉應力準則和Mohr-Coulomb準則進行損傷判別，并基于COMSOL Multiphysics with MATLAB平臺進行硐室圍巖損傷破裂數值模擬分析[26-29].

圖3 硐室開挖數值模型Fig.3 Numerical model for chamber excavation

2.1 硐室斷面軸比對圍巖損傷破裂影響

軸變論指出地應力條件不變時，橢圓形硐室的長短軸之比是影響巷道穩定性的重要因素，并提出了“等應力軸比”、“零應力軸比”、“壓應力軸比”與“拉應力軸比”的不同硐室形狀.由第1節分析可知，等應力軸比的圍巖應力均勻分布，不出現拉應力且最大應力值為最小，等應力軸比條件下巷道的圍巖是最穩定的.

為探究不同軸比對圍巖損傷破裂影響情況，構建側壓系數λ=1.5的地應力條件下不同軸比硐室，分別設定Z=1/3、1/2、2/3、1時的拉應力軸比，Z=1.5時的等應力軸比及Z=2、3時的壓應力軸比；且硐室斷面面積均為20 m2.各軸比條件下損傷破裂區及彈性模量演化如圖4和圖5所示，硐室頂部和中線兩幫應力分布如圖6所示.

圖4 硐室斷面軸比對圍巖損傷破裂影響.（a）Z=1/2；（b）Z=2/3；（c）Z=1；（d）Z=3/2；（e）Z=2Fig.4 Effect of the chamber axial ratio on the damage of country rock: (a) Z=1/2; (b) Z=2/3; (c) Z=1; (d) Z=3/2; (e) Z=2

圖5 硐室斷面軸比對圍巖彈性模量的影響.（a）Z=1/2；（b）Z=2/3；（c）Z=1；（d）Z=3/2；（e）Z=2Fig.5 Effect of the chamber axial ratio on the elastic modulus of country rock: (a) Z=1/2; (b) Z= 2/3; (c) Z=1; (d) Z=3/2; (e) Z=2

圖6 不同軸比硐室X方向應力分布情況.（a）拱頂處；（b）水平中線處Fig.6 Stress distribution in the X direction of the chamber with various axial ratios: (a) chamber roof; (b) middle route of the chamber

數值模擬結果表明，側壓系數λ= 1.5的地應力條件下，拉應力軸比硐室的應力集中區域主要表現在拱頂與拱底位置，損傷破裂區亦分布在拱頂和拱底.壓應力軸比圍巖的應力集中情況與損傷區域均與拉應力軸比相反.等應力軸比時，應力圍繞硐室均勻分布，且最大應力為130 MPa，均小于其余軸比硐室最大應力值，等應力軸比條件下，損傷區呈圍繞硐室零星分布.

不同的軸比條件下，硐室圍巖損傷破裂區大小不同.如圖7所示，等應力軸比條件下，圍巖損傷區域面積最小，為2.21 m2；壓應力軸比損傷面積次之；拉應力軸比相對損傷面積最大.因此，在實際工程中硐室斷面設計應避免拉應力軸比情況.此外，地應力條件不變的情況下，硐室軸比與等應力軸比差別越大，損傷破裂區面積越大.

圖7 圍巖損傷破裂區與硐室軸比條件關系Fig.7 Relationship between the damaged zone and chamber axial ratio

2.2 硐室斷面形狀對圍巖損傷破裂影響

根據地下工程的實際用途，往往會采用不同形狀的硐室斷面設計.不同的硐室斷面形狀對應于不同程度的應力集中，對硐室圍巖損傷破裂具有重要影響.本節構建硐室面積為20 m2的圓形、橢圓形（軸比Z=1.5）、馬蹄形、矩形（長短軸比為1.5）和正方形硐室進行圍巖損傷破裂過程的數值模擬.圖8和圖9分別給出了橢圓形、矩形、圓形、馬蹄形和正方形硐室圍巖損傷破裂區和彈性模量分布云圖.

圖8 硐室斷面形狀對圍巖損傷破裂影響.（a）橢圓形；（b）矩形；（c）圓形；（d）馬蹄形；（e）正方形Fig.8 Effect of the chamber shape on the extent of damage of country rock: (a) ellipse; (b) rectangle; (c) circle; (d) horseshoe; (e) square

圖9 硐室斷面形狀對圍巖彈性模量影響.（a）橢圓形；（b）矩形；（c）圓形；（d）馬蹄形；（e）正方形Fig.9 Effect of the chamber shape on the elastic modulus of country rock: (a) ellipse; (b) rectangle; (c) circle; (d) horseshoe; (e) square

如圖10所示，相較于圓形、橢圓形硐室，其余形狀硐室在直角處均產生較大程度應力集中，應力最大值高達531 MPa，誘使該區域圍巖發生損傷破裂現象.矩形、正方形硐室圍巖損傷破裂演化過程表明，損傷破裂區輪廓向近似于橢圓形發展，證實了于學馥先生提出的橢圓形硐室是各種形狀巷道破壞后的最終形態[16-17].需要指出的是，損傷力學計算結果表明，盡管在矩形、馬蹄形和正方形硐室邊角處產生了較大應力集中，其損傷破裂區范圍不會大幅增大或無限擴展.圍巖的損傷破裂將會大幅降低邊角處應力集中區域的應力集中程度，并調整該區域力學平衡.

圖10 不同形狀硐室X方向應力分布情況.（a）拱頂處；（b）水平中線處Fig.10 Stress distribution in the X direction of a chamber with various shapes: (a) chamber roof; (b) middle route of the chamber

統計各斷面形狀硐室圍巖損傷破裂區面積可知，等應力軸比橢圓形硐室圍巖損傷破裂區面積最小2.21 m2，其次為與側壓系數相等軸比的矩形硐室3.52 m2，正方形硐室損傷破裂區面積最大5.58 m2，如圖11所示.盡管矩形硐室和正方形硐室斷面面積相同，然而矩形硐室損傷破裂區面積遠小于正方形硐室.由此可見，相同地應力條件下，硐室斷面設計應先以地應力條件為主要參考，依據硐室使用功能設計最符合地應力條件的硐室斷面，從而有效控制圍巖損傷破裂區范圍.

圖11 圍巖損傷破裂區與硐室形狀條件關系Fig.11 Relationship between the damaged zone and chamber shape

2.3 地應力場對圍巖損傷破裂影響

地應力是引起地下硐室圍巖變形破壞的根本原因[30]，是巖體力學中的重要組成部分.深部地應力分布規律多表現為水平主應力較大，豎向主應力為最小主應力；水平主應力與豎向主應力的比值即側壓系數，一般情況下為0.5～5.5，多數情況大于1.地質構造對地應力場具有重要影響，地應力量值和方向都有較大的變化；同樣地，由于地下硐室群開挖或采礦等活動，區域應力場往往也會發生較大變化.

2.3.1 側壓系數對圍巖損傷破裂影響

隨著地層深度的增加，硐室周圍原巖應力水平不斷增大，尤其是在地質構造強烈地區，水平主應力遠大于豎向主應力.側壓系數的增大，加劇了深部硐室開挖后圍巖損傷破裂程度，也使得開挖后的支護變得更為復雜[31-32].為探究側壓系數對深部硐室變形與破壞影響，模擬了豎向主應力為50 MPa，側壓系數分別為1、1.5、1.75、2和2.2，橢圓形硐室圍巖損傷破裂情況.數值模擬結果表明，低側壓系數條件下，硐室周邊應力表現為壓應力；側壓系數大于硐室軸比后，硐室的拱頂拱底部位出現較高的壓應力集中區，損傷情況變化加劇.側壓系數越大，應力越集中于頂、底板兩側，應力集中區域多以拉伸損傷為主.圖12和13給出了各側壓系數下橢圓形硐室損傷破裂區及彈性模量分布情況.如圖14（a）所示，橢圓形硐室頂部Mises應力分布隨側壓系數增大而不斷增大.圖14（b）表明，橢圓形硐室中線處Mises應力分布隨側壓系數增大而不斷減小.圖15圍巖損傷破裂區面積統計結果表明，保持橢圓形硐室斷面軸比不變，圍巖損傷破裂區面積隨著側壓系數增大不斷增大，損傷破裂區面積與側壓系數呈指數關系增長；側壓系數大于硐室軸比后，損傷破裂面積快速增大.需要特別關注的是，當側壓系數增大到一定數值時，硐室圍巖將會發生非穩定破裂；將此時的側壓系數定義為臨界側壓系數.圖16給出了隨著豎向應力增大，臨界側壓系數變化規律.隨著埋深的增大，硐室的軸比越趨近于1.

圖12 側壓系數對圍巖損傷破裂影響.（a）λ=1；（b）λ=1.5；（c）λ=1.75；（d）λ=2；（e）λ=2.2Fig.12 Effect of the lateral pressure coefficient on the damage of country rock: (a) λ=1; (b) λ=1.5; (c) λ=1.75; (d) λ=2; (e) λ=2.2

圖13 側壓系數對圍巖彈性模量影響.（a）λ=1；（b）λ=1.5；（c）λ=1.75；（d）λ=2；（e）λ=2.2Fig.13 Effect of the lateral pressure coefficient on the elastic modulus of country rock: (a) λ=1; (b) λ=1.5;(c) λ=1.75; (d) λ=2; (e) λ=2.2

圖14 不同側壓系數硐室X方向應力分布情況.（a）拱頂處；（b）水平中線處Fig.14 Stress distribution in the X direction of a chamber with various lateral pressure coefficients: (a) chamber roof; (b) middle route of the chamber

圖15 圍巖損傷破裂區面積與地應力側壓系數條件關系Fig.15 Relationship between the damaged zone and in situ stress lateral pressure coefficient

圖16 不同地應力條件下對應的臨界側壓系數Fig.16 Critical lateral pressure coefficient under various in situ stress conditions

保持地層條件不變，分析不同豎向應力作用下臨界側壓系數變化規律.數值模擬結果表明，隨著埋深增大，硐室臨界側壓系數不斷減小.因此，對于深部地層條件而言，側壓系數的變化對硐室圍巖損傷破裂影響顯著，側壓系數的微小改變可能造成硐室圍巖大范圍損傷破裂.

2.3.2 構造應力場對圍巖損傷破裂影響

構造應力由地殼運動及各板塊間的拉伸擠壓構成[33].在實際工程中，一些硐室和巷道受構造應力的影響，可能會發生大變形、坍塌甚至破壞的風險，影響工程的安全高效運行.因此，構造應力場作用下硐室的布置與支護設計應與一般情況不同[34-35].

本節對側壓系數和硐室軸比不變條件下，不同β（β= 0°, 15°, 30°, 45°）構造應力場作用下硐室圍巖損傷破裂過程進行了數值模擬.當β= 0°時，硐室圍巖應力均勻分布且相等；β增大后，硐室圍巖垂直于最大主應力處應力增大，產生明顯的應力集中現象.數值模擬結果表明，硐室圍巖損傷破裂區面積隨著β增大呈增大趨勢；β= 45°時，損傷破裂區面積最大.如圖17和18所示，構造應力的方向對硐室圍巖損傷破壞區位置具有重要影響.隨著β增大，損傷破裂區由圍繞硐室均勻分布，逐漸向垂直于最大主應力方向的應力集中處轉移.圖19（a）和（b）分別給出了橢圓形硐室頂部和中線處Mises應力分布，構造應力場對Mises應力分布具有重要影響.因此，深部地下硐室設計時應根據實際地應力場進行動態調整，從而有效降低構造應力對硐室圍巖穩定性的不利影響，減小圍巖損傷破裂區范圍.圖20給出了圍巖損傷破裂區面積與構造應力角度間的關系.

圖17 構造應力角對圍巖損傷破裂影響.（a）β=0°；（b）β=15°；（c）β=30°；（d）β=45°；（e）β=60°Fig.17 Effect of the tectonic stress dip on the damage of country rock: (a) β=0°; (b) β=15°; (c) β=30°; (d) β=45°; (e) β=60°

圖18 構造應力角對圍巖彈性模量影響.（a）β=0°；（b）β=15°；（c）β=30°；（d）β=45°；（e）β=60°Fig.18 Effect of the tectonic stress dip on the elastic modulus of country rock: (a) β=0°;(b) β=15°;(c) β=30°;(d) β=45°;(e) β=60°

圖19 不同構造應力角硐室X方向應力分布情況.(a)拱頂處;(b)水平中線處Fig.19 Stress distribution in the X direction of a chamber with different tectonic stress dips: (a) chamber roof; (b) middle route of the chamber

圖20 圍巖損傷破裂區面積與構造應力角度關系Fig.20 Relationship between the damaged zone area and tectonic stress dip

3 地下硐室斷面設計與地壓調控討論

地壓是引起圍巖變形破壞的根本作用力，因此工程所在處的實測地應力是力學分析的前提條件.不同的地應力場，硐室斷面設計及其力學分析和支護設計不盡相同.根據實測地應力，進行硐室斷面設計，降低硐室圍巖應力集中系數，從而科學論證硐室垮落和自穩演變過程.不同的地壓具有不同的變形破壞特征和規律，應采用不同的力學理論進行機理分析，確定合理的圍巖穩定性控制策略.

三山島金礦西嶺礦區規劃設計了國內擬建的最深豎井，井深為2005 m.三山島金礦西嶺礦區深部地層地應力測試結果表明，埋深1900 m處地層最大水平主應力為72 MPa，最小水平主應力為48 MPa，豎向應力為50 MPa.根據軸變論的思想，巷道走向應沿最大水平主應力方向布置，而后根據最小水平主應力與豎向應力比值得到橢圓形硐室設計軸比Z= 0.96.然而，實際工程中金屬礦巷道一般采用拱形巷道，本文建立了寬4 m，高4 m的三心拱巷道，并對其進行了損傷破裂分析.如圖21～23所示，等應力軸比硐室應力均勻分布在巷道周邊圍巖，且損傷區域也圍繞巷道周邊零星分布，根據數值模擬結果對損傷單元進行統計可得出損傷區域面積為4.36 m2；三心拱巷道則在拱腳處產生較強的應力集中，損傷區域集中在拱腳與拱頂處，損傷區域面積為5.58 m2，略大于圓形巷道.相同巷道斷面面積條件下，三心拱巷道可利用有效面積及實用性均優于圓形硐室.通過對三心拱巷道拱頂和拱腳進行適當錨噴支護可保障巷道安全服役.

圖21 三山島金礦設計巷道圍巖損傷破裂.（a）圓形；（b）三心拱形Fig.21 Damage to a roadway in the Sanshandao gold mine: (a) circle;(b) three-centered arch

圖22 三山島金礦設計巷道圍巖彈性模量.（a）圓形；（b）三心拱形Fig.22 Elastic modulus of the roadway in the Sanshandao gold mine:(a) circle; (b) three-centered arch

圖23 三山島金礦設計巷道 Mises應力分布情況.（a）圓形；（b）三心拱形Fig.23 The von mises stress distribution in a roadway in the Sanshandao gold mine: (a) circle; (b) three-centered arch

4 結論

（1）地下硐室開挖圍巖損傷破裂程度與地層地應力和硐室斷面息息相關.硐室形狀與軸比是影響圍巖損傷破裂的重要因素，相同開挖面積下等應力軸比的橢圓形硐室應力分布均勻，應力集中效應最小，損傷破裂區零星分布在硐室周圍，損傷破裂區面積小，是硐室設計的最優解.

（2）地應力場是造成地下硐室圍巖變形破壞的根本原因，側壓系數越大，硐室頂、底板處應力峰值越大，損傷破裂越嚴重，應力集中區圍巖以拉伸破裂為主.側壓系數大于硐室軸比后，圍巖損傷區域面積隨側壓系數增大呈指數性增大.數值模擬結果表明，隨著地層深度增大，圓形和橢圓形硐室的臨界側壓系數逐漸減小并趨于1，即硐室開挖軸比的微小變化將會引起圍巖嚴重損傷破裂.

（3）構造應力場作用下地下硐室圍巖損傷破裂程度和區域發生顯著變化.構造應力場角度越大，圍巖損傷破壞程度越嚴重，損傷破裂區由均勻隨機分布逐漸向垂直于最大主應力方向的圍巖應力集中處轉移.構造應力場作用下地下硐室開挖誘發冒頂和巖爆的風險大幅升高.

（4）實際工程中，地下硐室的設計與布置應結合實測原位地應力條件，確保硐室形狀和走向符合地應力條件.硐室軸向設計應沿最大水平主應力方向，從而有效減小硐室圍巖應力集中程度，最大程度降低地應力場對圍巖損傷破裂及穩定性的不利影響.

致謝

謹以此文緬懷我國著名巖石力學與采礦工程專家于學馥先生.