軟巖巷道開挖支護的顆粒離散元模擬

廖九波 ，李夕兵 ，周子龍 ，鄒洋

(1.中南大學 資源與安全工程學院，湖南 長沙，410083；2.中南大學 深部金屬礦產開發(fā)與災害控制湖南省重點實驗室，湖南 長沙，410083)

隨著人類對礦產資源需求量的日益增大，已使淺部資源不斷減少和枯竭，國內外很多礦山相繼進入深部開采，深部高應力軟巖巷道開挖與支護問題也變得日益突出[1?2]。深部高應力軟巖巷道一般具有如下特征：(1) 巷道埋深大、構造應力較大或受采動影響；(2) 來壓時間快、持續(xù)時間長、初期變形量大；(3) 圍巖較松軟破碎、流變性和蠕變性較大；(4) 圍巖遇水極易崩解、強度降低明顯等特點。若相關措施采取不當，巷道圍巖變形將越演越烈，最終致使巷道失穩(wěn)破壞。因此，研究高應力軟巖巷道的破壞機制和支護方式成為礦業(yè)界亟待解決的重大問題[3?8]。PFC2D 既可解決靜態(tài)問題也可解決動態(tài)問題，既可用于參數(shù)預測，也可用于在原始資料詳細情況下的實際模擬。Hazard等[9]應用顆粒流分析軟件PFC研究了巖石材料的力學實驗，發(fā)現(xiàn)得到的應力應變曲線與伺服控制室內試驗所得結果能很好地吻合，表明離散元顆粒流軟件PFC具有模擬巖石基本力學特性的能力。Wang等[10]對巖石中的開挖問題運用PFC2D進行了模擬，探討了圓洞周圍的破壞情況。本文作者以貴州開磷集團馬路坪礦某巷道施工過程中出現(xiàn)的軟弱圍巖為背景，結合現(xiàn)場實際情況，對巷道圍巖的破裂區(qū)、塑性區(qū)和彈性區(qū)進行彈塑性力學分析，探討其應力和位移變化規(guī)律。同時，用顆粒流分析軟件PFC2D對巷道支護前后圍巖的位移變化進行了數(shù)值模擬，這有助于進一步研究軟巖巷道的破壞變形特性和評價支護結構的效果。

1 工程地質概況及軟巖特性

馬路坪礦750中段某軟巖巷道主要布置在紅頁巖中，巷道為三心拱斷面。巷道所在位置總體為一單斜構造，屬洋水背斜東翼，地層傾向 120°～130°，傾角30°～35°。圖1所示為具有代表性的此巷道獨頭掘進的紅頁巖斷面，紅頁巖呈明顯的傾斜層狀分布，下部為在潮濕空氣中風化崩解脫落的紅頁巖碎塊，從該處取回巖樣標本，在地質顯微鏡下進行了切片分析。

樣本的外觀特征描述為：紫紅色，層理構造，泥質結構，塊狀構造，礦物顆粒細小。經檢驗，其主要礦物成分見表1。

圖1 典型紅頁巖斷面Fig.1 Section of representative red shale

表1 紅頁巖的主要礦物成分(質量分數(shù))Table 1 Major mineral components of red shale %

顯微鏡下，該巖石樣本的結構組成如下：樣品具粉砂結構，泥質結構，層理構造。填隙物含量約為55%，碎屑含量約為45%。碎屑成分主要為鉀長石、白云母、絹云母、石英以及微量鋯石；膠結類型為孔隙膠結；填隙物主要為絹云母、泥質物、石英、方解石和鐵質物。

此類紅頁巖遇水極易崩解，即使在空氣中也會因吸水而逐漸解體。在這種軟弱巖層中施工巷道面臨很大難題：圍巖變形較大，支護異常困難，施工工期延長。

2 深埋軟巖巷道圍巖變形破壞的彈塑性力學分析

在原巖中開挖巷道，圍巖應力將發(fā)生2個顯著變化：一是巷道徑向應力下降，周邊處為零；二是切向應力增加，出現(xiàn)了應力集中現(xiàn)象。另一方面，巖體由三向應力平衡狀態(tài)變成了近似二向，巖石的強度下降了。如果集中應力小于巖體強度，那么圍巖將處于彈塑性穩(wěn)定狀態(tài)；當應力超過圍巖強度之后，巷道周邊圍巖將首先破壞，并逐漸向深部擴展，直至在一定深度取得三向應力平衡為止，此時圍巖已過渡到破碎狀態(tài)[11?12]。以圓形巷道為例，巷道圍巖將沿周邊至縱深依次形成破碎區(qū)、塑性區(qū)和彈性區(qū)，圖2所示為巷道圍巖變形分區(qū)。

2.1 力學模型

現(xiàn)對巷道進行彈塑性力學分析，其力學模型見圖3，建模做如下基本假設：

圖2 巷道圍巖變形分區(qū)Fig.2 Deformation layers within surrounding rock of tunnel

(1) 深埋無限長拱形平巷。

假設埋深Z大于或等于 20倍的巷道當量半徑R0(或其寬、高)。當埋深Z≥20R0時，可忽略此巷道影響范圍(3～5倍的R0)內的巖石自重。則在無限長的巷道長度里，圍巖的性質一致，能滿足彈塑性力學計算要求。所選模型可以看成平面應變問題，研究巷道任一截面的圍巖分區(qū)狀態(tài)。

(2) 假設圍巖為理想彈塑性體。

(3) 假設原巖應力為各向等壓的靜水壓力。

圖3 深埋圓形巷道力學模型Fig.3 Mechanical model of deep circular tunnel

2.2 巷道圍巖彈塑性位移解

彈塑性邊界的位移由彈性區(qū)的巖體變形引起，彈性區(qū)的變形可按外邊界趨于無窮、內邊界為Rp的厚壁圓筒來處理。

根據(jù)彈塑性力學中的平衡方程、幾何方程、本構方程和邊界條件，并利用巖石庫倫準則可計算得到在無支護時圓形巷道周邊的位移公式[13]為：

由于拱形巷道圍巖位移公式缺乏精確的理論解，故采用當量半徑折算法，任意形狀巷道尺寸折算成圓形巷道表示如下：

式中：R0為巷道當量半徑；k為斷面形狀修正系數(shù)，拱形取1.1；S為實際巷道斷面積。

圖4所示為拱形巷道力學模型。

圖4 拱形巷道力學模型Fig.4 Mechanical model of arched tunnel

由圖4可見：馬路坪紅頁巖拱形巷道拱寬為2 m，拱高0.5 m，直墻高1.5 m，可算得拱形巷道面積S為3.79 m2，進而巷道當量半徑：

現(xiàn)取紅頁巖力學參數(shù)：E=1 522 MPa，ν=0.33，p0=20.67 MPa，c=2.12MPa，φ=43°。利用式(2)，可算此拱形巷道圍巖塑性區(qū)半徑為：

巷道斷面的收縮量是0.22 m2，約占巷道總面積的5.7%?？梢妵鷰r只要處于完整狀態(tài)，其變形量相對較小。

從上述彈塑性區(qū)的變形求解公式分析可知：紅頁巖巷道開挖后變形量較大，需要立即對巷道進行支護，也就是臨時支護，其作用一是封閉圍巖，防止具有膨脹性圍巖吸水弱化，二是改變巷道圍巖的應力狀態(tài)，防止巷道圍巖破碎區(qū)的進一步擴展。馬路坪礦目前采取的支護措施是先掛網后噴漿，這種做法是不盡妥當?shù)?，要先噴漿再掛網打錨桿，在巷道變形前及時地噴漿能很好地起到封閉圍巖的作用。

3 巷道開挖與支護的顆粒流模擬

3.1 顆粒流軟件PFC2D簡介

PFC2D采用離散單元法能直接模擬球形顆粒間的運動和相互作用的物理問題，適用于研究固體材料的細觀/宏觀裂紋擴展、破壞累積并斷裂、破壞沖擊和微震響應等問題。同時，PFC2D內置接觸模型包括：簡單的黏彈性模型、簡單的塑性模型以及位移軟化模型。因此，PFC2D顆粒流程序也適用于模擬固體材料的彈塑性力學問題。

PFC2D模型中離散的顆粒單元不斷進行平移和旋轉運動，分離的顆粒在計算過程中又可重新形成新的接觸。顆粒流方法的理論基礎是力?位移定律和牛頓第二定律，采取顯式時步循環(huán)運算規(guī)則對圓形顆粒進行循環(huán)計算。根據(jù)力?位移定律確定 2個實體(顆粒?墻體或顆粒?顆粒)間的相對運動，并計算其間的接觸力。根據(jù)牛頓第二定律計算每個顆粒由于體積力或接觸力變化引起的顆粒運動(速度和位置)。

用PFC2D本身自帶的FISH語言，PFC2D程序在宏觀方面可以監(jiān)測模型的宏觀應力?應變關系、局部位移及表面裂隙的擴展情況，細觀方面可以觀察到模型內部各點的應力?應變關系、接觸力、速度、裂紋及應變能等物理量的變化[14?15]。

3.2 數(shù)值分析模型

為了使模擬結果更接近實際情況，需選取合理的顆粒體細觀參數(shù)，首先做紅頁巖試樣單軸抗壓試驗(見圖 5)，得到其單軸應力?應變曲線(見圖 6)，然后用PFC2D建立巖石單軸抗壓數(shù)值模型，分析模型中細觀參數(shù)對宏觀性質的影響，不斷試算得到與紅頁巖單軸抗壓試驗值基本一致的力學模型的細觀參數(shù)，當數(shù)值模型與室內實驗兩者所得的單軸應力?應變曲線能最大程度地吻合，表明此時的細觀參數(shù)可用來模擬該巷道開挖引起的圍巖破壞，其基本細觀參數(shù)取值見表2。

根據(jù)實際工程情況，取馬路坪礦750中段某軟巖巷道斷面作為模擬斷面，數(shù)值計算采用PFC2D建立二維數(shù)值模型，共定義上、下、左、右4道墻體，模型寬10 m，高10 m，顆粒半徑服從高斯分布。固定模型下部邊界，對上部邊界施加與覆巖厚度相應的垂直應力，對兩側邊界分別施加相應的水平應力。對模型施加的水平應力和垂直應力均為2.067×107N/m2，此時側壓力系數(shù)k為1，即靜水壓力狀態(tài)。

圖5 紅頁巖試樣單軸抗壓試驗Fig.5 Uniaxial compressive test of red shale sample

圖6 室內單軸應力?應變曲線Fig.6 Stress?strain curve by laboratory uniaxial test

表2 PFC2D數(shù)值模擬參數(shù)Table 2 PFC2D numerical simulation parameters

顆粒流模型重力平衡后，挖去中間拱形巷道的部分顆粒，模擬拱形巷道開挖。顆粒流模型的顆粒數(shù)目為70 028個，共挖去2 673個顆粒，剩余67 355個顆粒。巷道的數(shù)值模型見圖7。

拱形巷道參數(shù)如下：寬度2.0 m，直墻高1.5 m，拱高0.5 m，其中心點坐標x0=5.0，y0=5.0，即為模型正中心。

圖7 巷道的數(shù)值模型Fig.7 Numerical model of tunnel

3.3 計算結果分析

圖8所示為開挖后無支護巷道拉應力分布。先研究巷道未支護情況下的圍巖穩(wěn)定性，巷道圍巖在開挖前都是以壓應力為主，但在開挖后許多關鍵部位應力集中比較明顯，出現(xiàn)了較大的拉應力，巷道頂部、底部和兩幫數(shù)值較大(見圖 8)，其最大拉應力值分別為0.28，0.53和0.47 MPa，說明這些區(qū)域易受開挖的影響，由于巖體抗拉強度很小，易導致巷道圍巖受拉伸破壞。這些拉應力集中的關鍵部位，不利于巷道的穩(wěn)定，易出現(xiàn)較大破壞，應在最佳支護時間段內加強支護。

圖8 開挖后無支護巷道拉應力分布Fig.8 Tensile stresses distribution of unlined tunnel after excavation

圖9所示為拱形巷道圍巖無支護時開挖引起的位移場示意圖。由圖9可見：高應力軟巖巷道受開挖和原巖應力場的影響，拱頂、拱底及兩側的位移較大，說明這些位置受開挖的影響較大，而左上、左下、右上和右下4個角點處位移較小。位移最大值發(fā)生在拱底的位置，其值達到33.43 mm，易出現(xiàn)地壓底鼓現(xiàn)象。在深井軟巖中，頂板下沉量、底鼓量和兩幫移近量是整個巷道變形量的主要部分，在進行支護設計時主要就是控制巷道圍巖的這3種變形。數(shù)值模擬結果表明，軟巖巷道開挖后，應盡早地進行支護以便有效地控制圍巖的變形[16]。數(shù)值模擬中開挖后不支護所得巷道圍巖表面位移平均為29.99 mm，而前文彈塑性理論分析解為28.46 mm，兩者誤差僅為5.11%，可見數(shù)值模擬和理論分析兩者結果比較一致，進一步證明基本模型正確，進而利用模型模擬巷道支護問題，找到支護效果隨支護時間的變化規(guī)律。圖10所示為拱形巷道圍巖噴射混凝土和錨桿支護時位移場示意圖。由圖 10可見：支護已經起到了明顯效果，減少了圍巖的位移量。

圖9 無支護時位移場Fig.9 Displacement field in case of no liner

下面分別研究無支護和有支護2種工況下的圍巖位移變化規(guī)律，對支護效果進行評價。若開挖不及時采取支護措施將導致巷道圍巖出現(xiàn)不穩(wěn)定現(xiàn)象，在現(xiàn)場施工中應在開挖后盡早地開展支護工作。在計算分析過程中，將巷道豎向中心線選作軸線，模型具有左右對稱特性?，F(xiàn)取巷道左側的圍巖分析，在其中布置5個測量點，分別對拱頂、左上角、左拱腰、左下角、拱底處的圍巖位移發(fā)展狀況進行監(jiān)測，以此得到對應的位移來說明支護前后軟巖巷道圍巖的變形[17]。因為無法在開挖時就對危巖進行支護，支護總是在開挖一段時間之后才進行。在數(shù)值模擬過程中，這點表現(xiàn)在支護之前，圍巖在自重作用下已經發(fā)生了徑向位移。因此，分別做了開挖后時間分別為1 000，2 000，3 000，4 000，5 000，6 000和7 000步才開始支護的模擬，得到了一些結論。在支護這項數(shù)值模擬中，利用了 2種支護措施：噴射混凝土層厚度0.4 m，采用接觸連結；布置5根長1.5 m的錨桿，分別與x軸正方向成30°，60°，90°，120°和 150°，本文噴射混凝土層采用小顆粒模擬，由于墻體尺寸有限，一次生成顆粒厚度較小，為此進行兩次顆粒生成以增加噴漿的厚度。錨桿體由定義為剛性重疊連接的圓形大顆粒組成，定義顆粒的黏結強度、密度、直徑等指標，采用Clump命令，使顆粒相互重疊連接，形成錨桿體。模型中襯砌支護結構的力學參數(shù)見表3。

圖10 支護時位移場Fig.10 Displacement field in case of liner

現(xiàn)選取巷道開挖后時間為3 000步才開始支護為代表來說明問題，圖11～13所示分別為拱頂?shù)呢Q向位移、左拱腰水平位移和拱頂?shù)呢Q向位移隨時步的變化規(guī)律。由圖11～13可以看到：2種工況下巷道圍巖在1萬步左右時基本穩(wěn)定，現(xiàn)場觀測馬路坪礦該軟巖巷道在開挖1月后變形基本穩(wěn)定，將模型中時間步數(shù)換算成實際時間，發(fā)現(xiàn)與現(xiàn)場情況基本相符。

從圖11～13可以發(fā)現(xiàn)：無論是底鼓量、兩幫移近量和頂板沉降量都有開挖初期變形量大，增長速率高，它們的增長速率會隨時間的增加而逐漸降低，最后曲線漸漸趨向于水平達到穩(wěn)定平衡。從圖11可以看出：軟巖巷道拱頂豎向位移隨時間的變化規(guī)律，開挖造成的位移有28.74 mm，支護之后位移只有22.12 mm，采用支護措施后位移減小了6.62 mm，說明支護在抑制圍巖變形方面起到了一定的作用。圖12中，支護前左拱腰水平位移為27.82 mm，支護后，位移只有20.95 mm，減少了 6.87 mm，表明支護起到了作用。從圖13可以看出：支護前拱底的豎向位移為33.43 mm，而采用支護等措施后，位移僅發(fā)展了26.40 mm，拱底的豎向位移得到了控制。若在最佳支護時間對圍巖施加最優(yōu)的支護措施，能進一步縮減圍巖的變形量。

以上只對開挖后時間為3 000步進行支護的情況開展了討論，為了能充分說明支護快慢影響位移值大小，有必要對巷道開挖后若干不同時步才支護的圍巖位移變化情況分別予以模擬。為此，整理出各監(jiān)測圈最終位移隨不同支護時步的變化規(guī)律如表4所示，最終位移隨不同支護時步的曲線如圖14所示?？梢姡洪_挖至支護的時間步數(shù)越少，即時間越短，拱頂、左拱腰和拱底的最終位移值將會越小，在高應力軟巖巷道開挖時，應盡早進行支護加固工作[18]。圍巖表面位移最大值發(fā)生在拱底的位置，此處易出現(xiàn)地壓底鼓現(xiàn)象，應加強支護。

圖11 拱頂豎向位移曲線Fig.11 Curve of vertical displacement for vault

圖12 左拱腰水平位移曲線Fig.12 Curve of horizontal displacement for left hance

表3 襯砌支護結構的力學性能Table 3 Mechanical properties for liner and timbering

表4 各監(jiān)測圈最終位移隨不同支護時步的變化Table 4 Variation of final displacement of every monitored circles versus different steps of lining mm

圖13 拱底豎向位移曲線Fig.13 Curve of vertical displacement for arch foot

圖14 最終位移隨不同支護時步的曲線Fig.14 Curve of final displacement versus different steps of lining

4 結論

(1) 根據(jù)理想彈塑性模型中的變形協(xié)調條件、基本方程和相應的邊界條件，并利用巖石強度理論—庫侖準則，將巷道圍巖分成彈性區(qū)、塑性區(qū)和破裂區(qū)，并建立了在特定原巖應力條件下軟巖拱形巷道開挖后的圍巖位移解析表達式。

(2) 用顆粒流分析軟件PFC2D對軟巖巷道開挖后圍巖變形破壞規(guī)律進行了數(shù)值模擬，可以清晰地看到軟巖巷道在開挖后的應力分布情況，詳細探討了圍巖各個關鍵部位的位移分布特征，所獲巷道表面位移與彈塑性解比較一致，證明基本模型正確，進而利用模型模擬巷道支護，找到支護效果隨支護時間的變化規(guī)律，結果表明開挖后支護時間越早，越能控制圍巖變形量。所獲研究結論，可為復雜應力條件下類似巷道的變形研究和支護設計提供理論依據(jù)。

(3) 分別運用彈塑性力學理論和PFC2D數(shù)值模擬得到了馬路坪礦特定紅頁巖巷道圍巖表面的位移變形量，位移最大值發(fā)生在拱底的位置，其值達到 33.43 mm，易出現(xiàn)地壓底鼓現(xiàn)象，應盡早地進行支護以便有效地控制圍巖變形。

[1] 古德生, 李夕兵.有色金屬深井采礦研究現(xiàn)狀與科學前沿[J].礦業(yè)研究與開發(fā), 2003, 23(2): 1?5.GU Desheng, LI Xibing.Science problems and research state of deep mining in metal and nonferrous mines[J].Mining research and development, 2003, 23(2): 1?5.

[2] 李夕兵, 古德生.深井堅硬礦巖開采中高應力的災害控制與破碎誘變[C]//科學前沿與未來(第6集).北京: 中國環(huán)境科學出版社, 2002: 101?108.LI Xibing, GU Desheng.The hazard control and cataclastic mutagenesis induced by high stress in hard rock mining at depth[C]//Science Foreland and Future(Volume 6).Beijing:China Environmental Science Press, 2002: 101?108.

[3] 何滿潮, 景海河, 孫曉明.軟巖工程力學[M].北京: 科學出版社, 2002: 38?45.HE Manchao, JING Haihe, SUN Xiaoming.Engineering mechanics of soft rock[M].Beijing: Science Press, 2002: 38?45.

[4] Sellers E J, Klerck E.Modeling of the effect of discontinuities on the extent of the fracture zone surrounding deep tunnels[J].Tunneling and Underground Space Technology, 2000, 15(4):463?469.

[5] Russo-Bello F, Murphy S K.Longwalling at great depth in a geologically disturbed environment.the Way forward[J].The Journal of South African Institute of Mining And Metallurgy,2000, 100(2): 91?100.

[6] 王連國, 李明遠, 王學知.深部高應力極軟巖巷道錨注支護技術研究[J].巖石力學與工程學報, 2005, 24(16): 2889?2893.WANG Lianguo, LI Mingyuan, WANG Xuezhi.Study on mechanisms and technology for bolting and grouting in special soft rock roadways under high stress[J].Chinese Journal of Rock Mechanics and Engineering, 2005, 24(16): 2889?2893.

[7] 王其勝, 李夕兵, 李地元.深井軟巖巷道圍巖變形特征及支護參數(shù)的確定[J].煤炭學報, 2008, 33(4): 364?367.WANG Qisheng, LI Xibing, LI Diyuan.Surrounding rock deformation properties and determination of support parameters of soft rock roadway in deep mine[J].Journal of China Coal Socitey, 2008, 33(4): 364?367.

[8] 楊新安, 陸士良, 葛家良.軟巖巷道錨注支護技術及其工程實踐[J].巖石力學與工程學報, 1997, 16(2): 171?176.YANG Xinan, LU Shiliang.GE Jiahang.The technology and practice of bolting and grouting support method at soft rock roadway[J].Chinese Journal of Rock Mechanics and Engineering, 1997, 16(2): 171?176.

[9] Hazard J F, Young R P.Simulating acoustic emissions in bonded-panicles models of rock[J].International Journal of Rock Mechanics and Mining Sciences, 2000, 37(5): 867?872.

[10] Wang C, Tannant D D, Lilly P A.Numerical analysis of the stability of heavily jointed rock slopes using PFC2D[J].International Journal of Rock Mechanics and Mining Science,2003, 40(3): 415?424.

[11] 蔣斌松, 張強, 賀永年, 等.深部圓形巷道破裂圍巖的彈塑性分析[J].巖石力學與工程學報, 2007, 26(5): 982?986.JIANG Binsong, ZHANG Qiang, HE Yongnian, et a1.Elastoplastic analysis of cracked surrounding rocks in deep circular openings[J].Chinese Journal of Rock Mechanics and Engineering, 2007, 26(5): 982?986.

[12] Sharan S K.Exact and approximate solutions for displacements around circular openings in elastic-brittle-plastic Hock-Brown rock[J].Int J Rock Mech Min Sci, 2005, 42(4): 542?549.

[13] Kyung H P, Yong J K.Analytical solution for a circular opening in an elastic-brittle-plastic rock[J].Int J Rock Mech Min Sci,2006, 43(4): 616?622.

[14] ITASCA Consulting Group lnc.PFC2D theory and background[M].Minnesota: Itasca Consulting Group Inc, 2008:25?86.

[15] 徐泳, 孫其減, 張凌, 等.顆粒離散元法研究進展[J].力學進展, 2003, 33(2): 25l?260.XU Yong, SUN Qicheng, ZHANG Ling, et a1.Advances in particle discrete element methods[J].Advances in Mechanics,2003, 33(2): 25l?260.

[16] Funatsu T, Hoshino T, Sawae H, et al.Numerical analysis to better understand the mechanism of the effects of ground supports and reinforcements on the stability of tunnels using the distinct element method[J].Tunneling and Underground Space Technology, 2008, 23(3): 561?573.

[17] 陶連金, 蒯本秋, 張波.松散軟巖巷道破壞的顆粒離散元模擬分析[J].地下空間與工程學報, 2010, 6(2): 318?322.TAO Lianjin, KUAI Benqiu.ZHANG Bo.PFC the fracture mechanism analysis of loose soft roadway by discrete element method[J].Chinese Journal of Underground Space and Engineering, 2010, 6(2): 318?322.

[18] 孟云偉, 肖世洪, 柴賀軍, 等.隧道開挖中破碎帶支護的顆粒離散元模擬研究[J].地下空間與工程學報, 2007, 3(4):673?677.MENG Yunwei, XIAO Shihong, CAI Hejun, et al.PFC Simulation of timbering in fragment strip during tunnel excavation[J].Chinese Journal of Underground Space and Engineering, 2007, 3(4): 673?677.