巷道圍巖蝶形塑性區蝶葉方向性形成機制及工程應用

李 季，強旭博，馬念杰，張榮光，李 博

(1.西安科技大學 能源學院，陜西 西安 710054; 2.中國礦業大學(北京) 能源與礦業學院，北京 100083)

一直以來，巷道礦壓理論分析的基礎是經典的支承壓力理論，其核心是利用支承壓力集中系數來反應圍巖所處應力場的大小，但并未涉及圍巖所處應力場的方向。由于該理論具有簡單易行的優點，在巷道礦壓理論分析中被廣泛應用并發揮了重要的作用。然而，相關研究和實踐表明[1]，深部巷道在受到大埋深和強采掘的影響下，其所處應力場會變得相對復雜，其復雜性不僅體現在應力場大小的增減，而且其應力場方向也會發生改變[2-3]，進而導致深部巷道出現一些新的巷道礦壓現象，例如矩形巷道尖角處圍巖變形破壞并非最嚴重、巷道圍巖呈現非均勻變形破壞特征等[4-7]。眾所周知，地下巷道圍巖的變形破壞是由圍巖產生的塑性區引起的，塑性區的形態和范圍決定了圍巖變形破壞的模式和程度，而塑性區的產生是圍巖所處應力場和自身力學性質共同作用的結果。此時，采用經典支承壓力理論進行分析時，便會忽視應力場方向對深部巷道圍巖塑性區分布特征的影響，導致其無法很好地解釋深部巷道出現的一些新的礦壓現象。因此研究應力場方向對圍巖變形破壞特征的影響對于深部巷道礦壓理論分析具有重要的理論意義。

針對應力場方向對巷道圍巖變形破壞的影響，國內外學者先后進行了有益地探索。最早關于應力方向對巷道穩定性的研究即為澳大利亞學者GALE W J提出的最大水平應力理論[8]，該理論認為巷道圍巖的所受應力具有明顯的方向性，且巷道頂底板破壞主要受水平應力影響，即巷道與水平應力的夾角越小，頂底板的穩定性越好;此后，眾多學者以此為基礎進一步研究了圍巖應力方向與巷道軸向之間的關系對巷道圍巖穩定性的影響，孫玉福[9]利用數值模擬得到圍巖塑性區會隨著巷道軸向與最大水平主應力間夾角的不同而發生改變，且夾角與巷道圍巖穩定性有很大關系;朱志潔等[10]通過分析得出巷道圍巖穩定性會隨著最大水平主應力與巷道軸向間夾角的改變呈現復雜變化;高圣元等[11]通過研究發現背斜構造中巷道軸向布置角度會對圍巖穩定性產生重要影響，夾角在60°附近時可降低圍巖變形;趙亮等[12]研究發現巷道與最大水平主應力以一定角度斜交時，巷道頂底板的變形會偏向某一側;高富強[13]認為最大水平主應力方向與巷道軸向的夾角越大時，巷道穩定性越差;趙龍[14]通過研究發現巷道軸線與最大主應力方向夾角在20°～70°時，對巷道圍巖變形及破壞范圍影響最為嚴重;田成東和白海波[15]認為巷道變形量會隨著巷道軸向與最大水平應力間夾角的增大而增大。近年來，部分學者開始關注采動應力場方向對圍巖變形破壞的影響，黃炳香等[16]通過研究發現在高地應力和強烈采動應力共同作用下，圍巖行為迅速表現為復雜的非穩態、非線性特征;張文忠[17]利用數值模擬的方法研究了斷層附近采動主應力方向偏轉情況及斷層面剪應力及剪應力增量的分布規律;周鋼等[18]研究認為原巖地應力的最大主應力以鉛垂應力為主，在采動過程中最大主應力方向與水平面夾角逐漸減小，且會嚴重影響巷道圍巖穩定性;陳上元等[19]認為采動應力造成巷道圍巖應力場的大小和方向發生改變，從而導致了巷道的非對稱變形;王家臣等[20]采用理論分析、數值模擬及現場實測的方法研究了采動應力旋轉角度對圍巖穩定性的影響，采動應力旋轉造成圍巖承載能力降低，采動應力旋轉角度愈大，圍巖穩定性愈差。

以上研究分別從圍巖應力方向與巷道軸向關系，采動應力場方向變化等方面入手，研究了圍巖應力場方向與巷道布置的關系，并開始關注到采礦活動引起采場頂板巖層運動，會造成圍巖所處應力場方向改變[16-18]，但并未涉及圍巖應力場方向對塑性區形態影響的定量規律。因此，筆者以巷道蝶形破壞理論為基礎，采用理論計算和數值模擬等方法，重點研究主應力方向對圍巖蝶形塑性區形態的影響規律及蝶形塑性區蝶葉方向性的形成機制。

1 塑性區蝶葉形態表征

實際工程中，巷道所處的天然應力場并非均勻應力場，而是會因巷道埋深、地質構造作用等因素，水平應力與垂直應力存在大小差異(側壓系數≠1)，在受采掘擾動后，井下巷道圍巖所處非靜水壓力場的非均勻程度將更加明顯。因此，研究巷道處于這一非靜水壓力場中的圍巖塑性區具有重要工程意義。關于非靜水壓力場下圓形巷道圍巖塑性區問題，于學馥等[21]給出了雙向非等壓條件下圓孔周圍的塑性區邊界方程，此后由多位學者進行分析修正，后經馬念杰等[22-23]詳細推導后對圓形巷道圍巖塑性區形態進行了歸納總結，提出了蝶形塑性區的概念，并對其形成機制及工程意義進行了詳細闡述。雖然對非靜水壓力場中圓形巷道圍巖塑性區邊界方程的研究已經基本成熟，但對于塑性區蝶葉方向性及其工程應用的相關研究卻鮮有涉及。

在非靜水壓力場中，巷道圍巖塑性區形態將不再是經典的巷道圍巖塑性區理論(如卡斯特奈公式和芬納公式)所計算出的如圓形或橢圓形等均勻形態，而是呈軸對稱的“蝶形”非均勻形態分布[22-23]，如圖1所示。由圖1可以看出，蝶形塑性區的邊界輪廓分別在4個象限內凸出，將此凸出部分定義為塑性區蝶葉[24]。蝶葉最大半徑即為塑性區的最大半徑，其位置用極坐標(rmax，θmax)表示，a為圓形巷道的半徑;P1為最大主應力;P3為最小主應力。

圖1 蝶形塑性區形態示意Fig.1 Schematic diagram of butterfly plastic zone

同時根據筆者已有研究成果[25]表明，在受到圍巖巖性和應力方向等因素的影響下，井下巷道圍巖蝶形塑性區蝶葉形態并不會和圖1所示的蝶葉形態嚴格一致，而是會發生一定程度地改變，如圖2所示。

圖2 巷道圍巖塑性區形態變化示意Fig.2 Schematic diagram of plastic zone morphology change of roadway surrounding rock

根據蝶形破壞理論，導致此類現象出現的原因不僅是由于反映應力場大小的非均勻系數η對蝶形塑性區蝶葉分布范圍及尺寸產生影響，應力場方向也會對蝶形塑性區蝶葉位置產生重要影響。

2 蝶葉方向性形成機制

2.1 理論分析

于學馥[21]和馬念杰[26-27]等給出了非等壓條件下圓形巷道圍巖塑性區邊界隱性方程為

(1)

為分析非靜水壓力場條件下圍巖塑性區蝶葉方向性，在式(1)中引入圍巖塑性區主應力方向偏轉角度α并進行因式分解可得

(2)

cos22(θ-α)+sin22(θ-α)=1

(3)

將式(3)代入式(2)并提出公因式，變形可得

(4)

令:

則式(4)可化簡為圍巖塑性區主應力方向偏轉角度α與極坐標位置角度θ的余弦四次隱性方程:

K1cos4(θ-α)+K2cos2(θ-α)+K3=0

(5)

式中，r為塑性區邊界上任一點的半徑;C為黏聚力;φ為內摩擦角;η為圍壓非均勻系數(定義其為最大主應力與最小主應力的比值);θ為塑性區邊界上任一點的極坐標角度;α為主應力方向偏轉角度(定義其為P3與豎直方向的夾角，順時針為正)。

分析式(5)隱形方程可知，Ki是關于巷道半徑、圍巖性質參數、圍壓以及塑性區半徑的函數，而當圍巖性質及圍壓大小一定時，巷道圍巖塑性區半徑也會為定值，則Ki可以看作常數。此時，塑性區邊界上對應點的極坐標角度與主應力方向偏轉角度存在如隱形方程所示的函數關系，2者之間的變化會同步發生。當塑性區邊界上每一點的極坐標角度因主應力方向偏轉而發生改變時，則圍巖塑性區整體會發生偏轉。

對于井下巷道而言，其原巖應力場中，最大主應力主要以水平應力為主，而采礦活動加載效應引起的采動支承壓力集中系數會達到2～3倍，導致巷道圍巖所處的應力場中最大主應力不再以水平應力為主。根據式(5)，通過改變主應力方向的偏轉角度α，固定其他計算參數不變，得到不同主應力方向下的塑性區邊界分布圖，并將其與初始主應力方向(最小主應力為豎直方向，最大主應力為水平方向)的塑性區邊界分布圖進行對比，如圖3所示。

圖3 圓形巷道圍巖塑性區邊界分布對比(理論計算)Fig.3 Comparison of the distribution of the plastic zone boundary of circular roadway rock surrounded (Theoretical calculation)

通過圖3可以看出，在非靜水壓力場環境下，均質圓形巷道圍巖會形成蝶形塑性區，且蝶形塑性區的4個蝶葉會隨著主應力方向的偏轉而隨之發生偏轉。當主應力方向未發生偏轉時，將第1,2,3,4象限的蝶葉分別標記為A，B，C，D，表1為主應力方向偏轉不同角度時塑性區4個蝶葉位置的極坐標。由表1可知，圍巖塑性區的4個蝶葉位置會隨著主應力方向的偏轉而改變，且蝶葉位置的偏轉角度與主應力方向的偏轉角度相等，如圖4所示。

圖4 主應力偏轉角度與蝶葉偏轉角度關系Fig.4 Diagram of relationship between deflection angle of principal stress and deflection angle of butterfly leaf

表1 主應力偏轉不同角度時塑性區蝶葉位置(理論計算)Table 1 Butterfly leaf position of plastic zone when principal stress deflects different angles(Theoretical calculation)

2.2 數值模擬驗證

為了驗證理論計算解的可靠性和普適性，利用FLAC3D數值模擬驗證圓形巷道圍巖塑性區蝶葉的方向性。數值模擬圍巖本構模型采用摩爾—庫倫模型，模型邊界采用應力約束，施加最大主應力60 MPa，最小主應力為26 MPa;當主應力方向發生偏轉時，根據彈塑性力學中的應力轉換公式進行主應力換算。數值模型巷道半徑a=3 m;其余數值模擬參數見表2。

表2 模型力學參數Table 2 Mechanical parameters of model

數值模擬得到的圓形巷道圍巖塑性區分布形態如圖5所示，以左上角的蝶葉為研究代表，隨著主應力方向發生偏轉，其與豎直方向的夾角隨之改變，蝶葉的偏轉角度與主應力方向的偏轉角度近似相等(2者差值不超過3°)，見表3。

表3 圓形巷道塑性區蝶葉偏轉角度(數值模擬)Table 3 Butterfly leaf deflection angle in plastic zone of circular roadway(Numerical simulation)

圖5 圓形巷道圍巖塑性區分布形態示意(數值模擬)Fig.5 Distribution pattern of plastic zone in circular roadway rock surrounded(Numerical simulation)

綜合理論分析和數值模擬可以看出，在相同參數條件下，理論分析和數值模擬得到的結果基本一致。當主應力方向發生偏轉時，圍巖蝶形塑性區的蝶葉會隨之發生偏轉，且偏轉的角度與主應力方向的偏轉角度基本相同，蝶形塑性區蝶葉具有顯著的方向性。

3 蝶葉方向性的適用性分析

3.1 巷道斷面形狀對蝶葉方向性的影響

目前，煤礦井下常見的巷道斷面形狀有圓形、矩形和拱形。因此，本節將分別通過理論計算和數值模擬的方法研究斷面形狀對蝶形塑性區蝶葉方向性的影響。受制于現有數學力學方法的限制，矩形巷道和拱形巷道圍巖的塑性區無法通過理論計算得到精確解。因此，對拱形巷道與矩形巷道塑性區邊界進行理論計算時，采取外接圓的方式并利用式(5)計算得到其主應力方向發生偏轉前后的塑性區形態對比圖，如圖6，7所示，并利用數值模擬進行驗證分析。矩形巷道尺寸為4 m×3 m;拱形巷道尺寸為上部半徑3 m，下部尺寸6 m×3 m。由圖6，7可以看出，隨著主應力方向的偏轉，拱形巷道與矩形巷道的蝶形塑性區蝶葉同樣會隨之發生偏轉，且兩者的偏轉角度相等(表4)。

圖6 拱形斷面巷道塑性區理論計算結果Fig.6 Theoretical calculation result of plastic zone in arched roadway

圖7 矩形斷面巷道塑性區理論計算結果Fig.7 Theoretical calculation result of plastic zone in rectangular roadway

表4 不同巷道斷面形狀塑性區蝶葉偏轉角度(理論計算)Table 4 Deflection angle of butterfly leaf in plastic zone of different roadway shapes (Theoretical calculation)

為了驗證上述理論計算的可靠性，采用數值模擬對矩形和拱形巷道的塑性區分布進行驗證分析，數值模擬參數與理論計算保持一致，見表5，其他模擬初始條件參數與圓形巷道相同。

表5 數值模擬力學參數Table 5 Numerical simulation of mechanical parameters

數值模擬得到的不同巷道斷面形狀的塑性區分布形態如圖8所示，并以左上角的蝶葉為例，得到了蝶葉隨主應力方向發生偏轉時，其與豎直方向的夾角。由圖8可以看出，當主應力方向發生偏轉時，不同巷道斷面形狀條件下的圍巖塑性區蝶葉位置均會發生偏轉，且蝶葉位置的偏轉角度與主應力方向的偏轉角度基本一致，誤差不超過3°，見表6。因此，圍巖蝶形塑性區蝶葉的方向性與巷道斷面形狀無關。

圖8 不同巷道斷面形狀塑性區分布形態Fig.8 Distribution pattern of plastic zone in different shape of roadway

表6 不同巷道斷面形狀塑性區蝶葉偏轉角度(數值模擬)Table 6 Deflection angle of butterfly leaf in plastic zone of different roadway shapes(Numerical simulation)

3.2 圍巖性質對蝶葉方向性的影響

由于煤系地層多為沉積形成的層狀巖體[28]，井下巷道圍巖基本為沉積巖，造成不同巖層之間的圍巖力學性質差異明顯。因此，為了研究圍巖性質對蝶形塑性區蝶葉方向性的影響，本節先對不同力學參數條件下圍巖蝶形塑性區蝶葉方向性進行研究，再通過數值模擬研究巷道復合頂板條件對蝶形塑性區蝶葉方向性的影響。

3.2.1力學參數

通過理論計算得到不同力學參數條件下圓形巷道圍巖蝶形塑性區分布圖，如圖9所示。由圖9可以看出，不同力學參數條件下，蝶形塑性區蝶葉均會隨著主應力方向的偏轉發生偏轉，且偏轉角度與主應力方向偏轉角度一致，見表7。因此，蝶形塑性區蝶葉的方向性與圍巖力學參數條件無關。

表7 不同力學參數下塑性區蝶葉偏轉角度(理論計算)Table 7 Deflection angle of butterfly leaf in plastic zone under different mechanical parameters (Theoretical calculation)

圖9 不同力學參數下塑性區邊界理論計算分布Fig.9 Theoretical calculation distribution of plastic zone boundary under different mechanical parameters

同時，通過數值模擬得到了不同圍巖力學參數條件下巷道圍巖塑性區分布圖，并以左上角的蝶葉為代表，得到其與豎直方向的夾角，如圖10所示。可以看出，不同圍巖力學參數條件下的圍巖塑性區的形態均呈蝶形形態，雖然蝶葉的形狀會因圍巖力學參數的不同存在差異，但當主應力方向發生偏轉時，圍巖塑性區都會隨之偏轉，蝶葉位置也會隨之發生改變，且偏轉角度都與主應力方向偏轉角度基本相同，誤差不超過3°，見表8。

圖10 不同力學參數下塑性區分布形態數值模擬Fig.10 Numerical simulation of plastic zone distribution under different mechanical parameters

表8 不同力學參數下塑性區蝶葉偏轉角度(數值模擬)Table 8 Deflection angle of butterfly leaf in plastic zone of different mechanical parameters (Numerical simulation)

3.2.2層狀巖層

由于層狀巖層條件下的塑性區無法利用理論計算得到，故利用數值模擬對4種類型的層狀巖層頂板結構的圍巖塑性區分布形態進行研究，分別為下軟上硬型、下硬上軟型、軟弱夾層型及堅硬夾層型頂板結構，各結構組成見表9，巖層及其巖石力學參數見表10，模擬結果如圖11所示。

表9 層狀結構頂板巖性組成Table 9 Lithologic composition of layered roof

表10 模型巖層力學參數Table 10 Mechanical parameters of model rock stratum

圖11 層狀巖層頂板條件巷道圍巖塑性區形態Fig.11 Plastic zone of surrounding rock of roadway under condition of layered rock roof

由圖11可知，由于頂板層狀巖層的圍巖性質存在差異，因此圍巖中較為堅硬巖層的存在會限制蝶葉塑性區的擴展，甚至使塑性區出現徹底隔斷[29]，同時，受到圍巖巖性的影響，塑性區蝶葉會發生不同程度的退化等變異現象[25]。但蝶葉方向性不會因發生變異或塑性區隔斷現象而消失，塑性區蝶葉同樣會隨著主應力方向的偏轉而偏轉，且兩者偏轉角度基本保持一致，誤差同樣不超過3°，見表11。

表11 層狀巖層頂板條件巷道圍巖塑性區蝶葉偏轉角度 (數值模擬)Table 11 Butterfly leaf deflection angle in plastic zone of roadway surrounding rock under roof condition of layered strata (numerical simulation)

綜上所述，不同巷道斷面形狀和圍巖性質條件下，蝶形塑性區蝶葉均會隨主應力方向偏轉而發生偏轉，且兩者的偏轉角度基本相等(誤差不超過3°)，蝶形塑性區蝶葉方向性與斷面形狀和圍巖性質無關，具有普適性。

4 蝶葉方向性對圍巖穩定性的影響

為了分析巷道圍巖變形量與圍巖蝶形塑性區形態的聯系，利用FLAC3D數值模擬得到了圍巖主應力方向偏轉40°時(圖12)的圍巖塑性區分布形態及巷道變形量等值線云圖。

由圖12可以看出，巷道頂板左側上方巖層、巷道左幫上部圍巖、巷道右幫下部圍巖和頂板右側下部圍巖由于形成了蝶葉塑性區，其塑性區尺寸較大，且與之對應的左側上方巖層的下沉量、巷道左幫上部圍巖和巷道右幫下部圍巖的變形量、頂板右側下部巖層的底臌量也比巷道其他位置圍巖變形量大，巷道圍巖變形量和塑性區尺寸呈正相關性。

圖12 主應力方向偏轉40°時塑性區與位移關系對比Fig.12 Comparison of the relationship between plastic zone and displacement when the principal stress direction is deflected by 40°

同時，通過對比分析不同主應力方向下矩形巷道蝶形塑性區分布圖(圖8(a))可以看出，主應力方向未發生偏轉時，巷道塑性區蝶葉主要位于巷道4個尖角處的深部圍巖中，這與普遍認為矩形巷道尖角處圍巖的變形破壞較嚴重是相符的;而當主應力發生偏轉時，巷道塑性區蝶葉不再位于巷道尖角處位置圍巖中，其會隨著主應力方向偏轉而位于巷道不同位置的圍巖中，這就解釋了在某些深部巷道中其變形破壞嚴重的位置并非一定出現在巷道尖角處，蝶葉方向性會導致圍巖的變形破壞具有非對稱性。

綜上所述，蝶形塑性區的蝶葉即為圍巖塑性破壞范圍，而蝶形塑性區尺寸與圍巖變形量呈正相關性。當圍巖塑性區的形態不再呈現均勻形態而呈現為蝶形分布形態時，蝶形塑性區蝶葉方向性會使蝶葉位于巷道圍巖的不同位置，導致巷道圍巖不同位置的變形破壞呈現差異化分布特征。尤其當塑性區蝶葉位于巷道頂板位置時，蝶葉塑性破壞范圍會成為巷道頂板的冒頂高風險區域。在實際工程中，由受圍巖低抗拉特性及節理化等特征的影響，巷道的實際冒頂形態與頂板塑性區蝶葉形態無法完全吻合，但產生蝶葉塑性區的頂板位置仍是冒頂的高發區域。

5 工程應用

趙固二礦11030工作面運輸巷與11011采空區留設8 m窄煤柱掘進(圖13)，斷面尺寸高3.3 m×寬4.8 m，埋深700 m左右，巷道頂板巖層主要為泥巖、砂質泥巖和大占砂巖，底板巖層為砂質泥巖和泥巖。該巷道掘進期間，巷道斷面收斂嚴重，礦壓顯現劇烈。現場礦壓監測數據顯示，巷道圍巖變形呈現明顯的非均勻大變形特征，如圖14所示。

圖13 巷道位置與工作面布置Fig.13 Roadway location and working face layout

圖14 11030運輸巷非均勻變形特征Fig.14 Non-uniform deformation characteristics of 11030 headgate

11303運輸巷圍巖破壞現場探測如圖15所示，根據文獻[30]研究，當煤柱尺寸為8 m時，圍巖應力場主應力方向大致會逆時針偏轉20°左右，導致其圍巖塑性區分布特征如圖16所示，塑性區蝶葉分別位于巷道頂板煤柱幫側、底板煤壁幫側、煤柱幫煤體及煤壁幫煤體，與圍巖現場探測破壞形態基本吻合。

圖15 11030運輸巷圍巖破壞現場探測Fig.15 Field detection diagram of surrounding rock failure of 11030 headgate

通過對比11030運輸巷巷道變形特征和11030運輸巷由數值模擬得到的圍巖塑性區分布形態可以看出(圖16)，巷道圍巖塑性區分布形態與圍巖的非均勻大變形特征高度一致。圍巖應力場主應力方向偏轉導致塑性區蝶葉分別位于巷道圍巖不同位置，塑性區深度較大的蝶葉位置的圍巖變形量相對較大。因此，蝶葉方向性是導致巷道圍巖出現非均勻大變形的主要原因。

圖16 11030工作面運輸巷圍巖非均勻變形破壞機理Fig.16 Non-uniform deformation failure mechanism diagram of 11030 headgate rock surrounded

6 結 論

(1)當主應力方向發生偏轉時，蝶形塑性區蝶葉會隨之發生一定的偏轉，且主應力方向偏轉角度與蝶葉偏轉角度基本相等，蝶形塑性區蝶葉具有顯著的方向性。

(2)不同斷面形狀和圍巖性質條件下的巷道蝶形塑性區蝶葉，均會隨主應力方向偏轉而隨之發生偏轉，即蝶葉方向性不受斷面形狀及圍巖性質的影響，具有普適性。

(3)蝶形塑性區蝶葉方向性決定了蝶葉產生的位置，導致蝶葉位置的圍巖變形破壞程度相對嚴重，且巷道圍巖的變形破壞特征與蝶形塑性區分布特征高度一致。