楚長青， 馮 旻， 李云濤

(1.鞍鋼集團礦業(yè)有限公司， 遼寧 鞍山 114009； 2.黑龍江科技大學(xué) 機械工程學(xué)院，哈爾濱， 150022； 3.鞍鋼集團礦業(yè)有限公司眼前山鐵礦， 遼寧 鞍山 114051)

0 引 言

礦山由露天轉(zhuǎn)入地下開采后，隨著采礦工程的不斷推進，開采對采場圍巖體產(chǎn)生的影響日趨復(fù)雜，主要原因是隨著地下開采的進行，地應(yīng)力的影響愈來愈明顯，受到開采擾動后，巖體內(nèi)部儲存能量大于其破壞前儲能能力，從而引發(fā)冒頂、垮落、塌方、巖爆等突發(fā)性地壓災(zāi)害[1-3]。地下開采中不同圍巖體的復(fù)雜性和斷面多樣性，使得利用經(jīng)典力學(xué)理論解析解的方法對地下工程穩(wěn)定性的分析和計算難以實現(xiàn)。

對于自然狀態(tài)下邊坡穩(wěn)定性，國內(nèi)學(xué)者常用極限平衡法進行分析，并取得了大量研究成果。最近十幾年來，由于計算機技術(shù)的發(fā)展，學(xué)者基于極限平衡法開發(fā)出了數(shù)值分析方法，如有限元、離散元等。楊宇江等[1]通過FLAC3D分析了不同開采深度條件下覆蓋層對露天礦殘余邊坡穩(wěn)定的影響。羅浪等[2]基于FLAC3D數(shù)值模擬，研究了云南某緩傾斜礦床露天終了坡高50 m和80 m對采場覆蓋巖層的影響程度，從而確定了保證此礦山安全穩(wěn)定的露天終了坡高范圍。李小雙等[3]采用數(shù)值模擬分析了35°、45°和55°坡角下露天轉(zhuǎn)地下采動時對上覆巖層的影響，得出45°左右對采場覆巖穩(wěn)定性較好。查文華等[4]結(jié)合現(xiàn)場勘測，利用FLAC3D詳細(xì)分析了逆傾分步開采對含斷層的高陡采礦邊坡穩(wěn)定性的影響。劉艷章等[5-6]采用彈性力學(xué)理論分析和數(shù)值模擬相結(jié)合的方法，研究了充填體的彈性模量和泊松比對舞陽鐵礦鐵古坑采區(qū)露天轉(zhuǎn)地下采動時邊坡穩(wěn)定性的關(guān)系。以密云鐵礦為例，利用有限元法分析了不同隔離層厚度對邊坡穩(wěn)定性的影響。魯宇等[7]運用相似材料模擬和PFC2D，研究了露天轉(zhuǎn)地下采動后30、35、40 cm厚度的上覆巖層的運動規(guī)律以及與邊坡穩(wěn)定性的關(guān)系。張鵬等[8]通過Midas Nx 有限元軟件分析了雙重應(yīng)力條件下采空區(qū)的規(guī)模和相對位置對露天轉(zhuǎn)地下開采邊坡穩(wěn)定性的影響。杜逢彬等[9]采用有限差分軟件分析了露天轉(zhuǎn)地下回采隔離頂柱的厚度與隔離頂柱穩(wěn)定性的關(guān)系。龍周彪等[10]通過理論分析和FLAC3D數(shù)值模擬，研究了姑山鐵礦基于進路充填采礦法境界頂柱厚度的取值范圍。

基于以上研究，筆者采用COMSOL Multiphysics軟件開展眼前山鐵礦大規(guī)模數(shù)值模擬計算工作，模擬了從左向右開采、從右向左開采、從中間向兩端開采、從兩端向中間開采四種開挖方案，并且根據(jù)計算的位移場、應(yīng)力場分析了四種開挖方案的圍巖穩(wěn)定性，對預(yù)測“露天轉(zhuǎn)地下”過渡期地下采場潛在地質(zhì)災(zāi)害，評估災(zāi)害風(fēng)險等級和區(qū)劃災(zāi)害潛在區(qū)域具有重要意義。

1 眼前山鐵礦工程概況

眼前山鐵礦原采用露天開采，后轉(zhuǎn)入地下開采。地采采礦方法采用無底柱分段崩落法，采場結(jié)構(gòu)參數(shù)為：分段高度18 m、進路間距20 m和崩礦步距2 m。隨著礦山由-177 m分段逐步完全轉(zhuǎn)入坑底以下礦體開采，作業(yè)面急劇增多，由-177 m分段的19條回采進路增加到-195 m分段的45條回采進路，多處回采進路進入首采分段開采，為準(zhǔn)確了解回采工作完成后巖體冒落沖擊形成機理進行數(shù)值模擬分析，建立研究區(qū)域的空區(qū)三維模型，頂部厚度在6.7 m到43.6 m之間，如圖1所示。

圖1 眼前山鐵礦-195首采分段排位設(shè)計Fig. 1 Sectional layout design of Yanqianshan Iron Mine -195 first mining

2 數(shù)值模型的建立

2.1 模型建立

將眼前山鐵礦三維地質(zhì)模型導(dǎo)入COMSOL Multiphysics軟件，通過軟件自帶網(wǎng)格生成工具，構(gòu)建全四面體網(wǎng)格的實體模型，生成眼前山三維數(shù)值計算網(wǎng)格模型。數(shù)值模型尺寸為2 690 m×1 840 m×946 m，地表標(biāo)高從-177 m到231 m起伏不等，整個數(shù)值模型共計3 527 845個單元，845 751個單元節(jié)點。邊界條件是固定模型底面，對四周采用法向位移約束。

2.2 參數(shù)選取

文中模擬采用的參數(shù)是基于室內(nèi)巖石力學(xué)試驗，結(jié)合現(xiàn)場結(jié)構(gòu)面賦存情況，通過Hoke-Brown強度進行計算，最終確定巖體基本力學(xué)參數(shù)如表1所示。

表1 眼前山巖體基本力學(xué)參數(shù)

2.3 本構(gòu)模型

圖2為巖石細(xì)觀單元在單軸和三軸應(yīng)力條件下的本構(gòu)關(guān)系。假設(shè)巖石試樣內(nèi)各個單元滿足線彈性本構(gòu)關(guān)系，即一般定義的胡克定律。當(dāng)單元的應(yīng)力條件滿足最大拉應(yīng)力準(zhǔn)則和Mohr-Coulomb準(zhǔn)則時，單元可能發(fā)生損傷。

圖2 巖石細(xì)觀單元損傷本構(gòu)關(guān)系 Fig. 2 Damage constitutive relation of rock meso-elements

由于巖石的抗拉強度遠(yuǎn)遠(yuǎn)小于其抗壓強度，最大拉應(yīng)力準(zhǔn)則通常是首先被用于判斷的，文中默認(rèn)以拉為正，最大拉應(yīng)力準(zhǔn)則和Mohr-Coulomb準(zhǔn)則表示為

(1)

式中：σ1——最大主應(yīng)力；

σ2——中間主應(yīng)力；

σ3——最小主應(yīng)力；

σt0、σc0——細(xì)觀巖石單元的單軸抗拉強度和單軸抗壓強度；

φ——巖石的內(nèi)摩擦角；

θ——洛德角；

J2——第二偏應(yīng)力張量不變量。

當(dāng)細(xì)觀單元滿足損傷閾值時，根據(jù)損傷力學(xué)調(diào)整巖石的力學(xué)性能，描述巖石的力學(xué)性能劣化，而峰后進一步的塑形變形描述材料由于裂隙存在的擴展、摩擦滑移引起的不可恢復(fù)變形。

將巖石單元的損傷演化本構(gòu)關(guān)系定義為

σij=(1-D)(2Gεij+λδijεkk)，

(2)

(3)

式中：D——損傷變量，其值為0～1；

G——剪切模量；

εij、εkk、ε1、ε3、εt0、εc0——各應(yīng)力狀態(tài)下對應(yīng)的應(yīng)變；

λ——拉梅系數(shù)；

δij——應(yīng)變增量。

將上述公式嵌入到COMSOL Multiphysics數(shù)值模擬程序中用來分析巖石漸進破壞過程。

3 模擬方案

眼前山-195分段回采過程會伴隨產(chǎn)生不同面積的暴露空間，造成圍巖以及頂板的應(yīng)力重新分布。隨著生產(chǎn)持續(xù)開展，地下暴露空間不斷增大，采場的應(yīng)力分布不斷遷移并重新分布，而采空區(qū)的冒落形式與頂板圍巖受力狀態(tài)息息相關(guān)，為了保證礦山生產(chǎn)的安全，需在最小的空區(qū)暴露面積下使頂板進行冒落，以減小冒落沖擊帶來的危害，因此，對不同回采工況下的應(yīng)力分布特征進行數(shù)值模擬計算，確定合理的開挖順序。

按照以往的礦房回采經(jīng)驗，將冒落研究區(qū)域劃分為從左向右開采(方案一)、從右向左開采(方案二)、從中間向兩端開采(方案三)、從兩端向中間開采(方案四)四種開挖方案。礦體區(qū)域劃分如圖3所示。

圖3 礦體區(qū)域劃分Fig. 3 Orebody regional division

4 采場圍巖穩(wěn)定性結(jié)果分析

4.1 不同回采工況下應(yīng)力場計算分析

為了研究不同方案開挖后露天采場的應(yīng)力場，模擬計算不同的開采階段，鑒于篇幅限制，以方案一為例說明。應(yīng)用方案一在不同開采階段最大主應(yīng)力云圖如圖4所示。

圖4 不同開采階段方案一最大主應(yīng)力云圖 Fig. 4 Maximum principal stress nephogram of scheme one in different mining stages

由圖4可知，在礦房開挖體積是398 250 m3(開挖回采工作進行68.7%)之前，開挖礦房上方最大主應(yīng)力等值曲線的最大、最小值分別為6 000、3 000 N/m2，最大主應(yīng)力由開挖礦房上方向礦房邊界轉(zhuǎn)移但并未發(fā)生較大的應(yīng)力集中，表明在此開挖階段頂板巖層的穩(wěn)定性相對較好，繼續(xù)開挖，采空區(qū)暴露面積持續(xù)增大，在采空區(qū)上邊界處發(fā)生較大的應(yīng)力集中，最大主應(yīng)力等值線達到9 000 N/m2，表明在最后開挖階段應(yīng)力集中發(fā)生明顯轉(zhuǎn)折，當(dāng)表層巖塊之間的聯(lián)系不足以克服自身重力時，塊體便會脫離母巖自然掉落，誘發(fā)大規(guī)模冒落現(xiàn)象從而影響礦山安全生產(chǎn)。

對比四種開挖方案在不同開挖階段最大主應(yīng)力σmax極值，如圖5所示。方案一在前期開挖階段未應(yīng)力集中，隨開挖繼續(xù)應(yīng)力迅速轉(zhuǎn)移并產(chǎn)生應(yīng)力集中，該回采方案造成的冒落危害相對較大；方案二的最大主應(yīng)力極值升高速率最快，隨礦房回采應(yīng)力集中現(xiàn)象最為迅速，右側(cè)冒落研究區(qū)域右側(cè)礦房頂板巖層穩(wěn)定性最差，該回采方案可在較小的暴露面積下發(fā)生冒落，冒落塊體及冒落的所帶來的能量沖擊較??；方案三開挖前期頂板巖層具有較好的穩(wěn)定性，隨地下暴露空間的持續(xù)增大，在采空區(qū)上邊界處發(fā)生較大的應(yīng)力集中，該方案前期回采冒落風(fēng)險較小后期回采的冒落風(fēng)險較大，冒落塊體及冒落的所帶來的能量沖擊較大；方案四回采前期頂板巖層的穩(wěn)定性較好，隨回采的持續(xù)進行，位于冒落研究區(qū)域中間的礦房對頂板的支撐作用不斷降低，當(dāng)?shù)V房的支撐強度不能約束圍巖微裂紋的產(chǎn)生與擴展時，未釋放的冒落能量得以積蓄。一旦失去礦房的支撐，這種冒落能量便會突然釋放，使圍巖微裂隙迅速貫通，從而有可能發(fā)生批量冒落甚至大規(guī)模冒落。

圖5 不同開采方案的最大主應(yīng)力最大值Fig. 5 Maximum value of maximum principal stress of different mining schemes

4.2 不同回采工況下位移場計算分析

通過模擬礦坑表面z方向位移分析空區(qū)穩(wěn)定性，鑒于篇幅限制，僅以方案一z方向位移云圖為例，結(jié)果如圖6所示。

圖6 不同開采階段方案一z向位移云圖Fig. 6 z-displacement cloud map of scheme one at different mining stages

從圖6可以看，出回采礦房上方沉降幅度較小且沉降等值線主要為-7.78×10-7、-8.2×10-7、-8.6×10-7m。隨著開采的進行，礦方上方沉降幅度不斷增加，由-7.78×10-7m增加到-8.6×10-7m，且沉降幅度較大的等值線范圍逐漸增大，這說明隨著開采進行，采空區(qū)冒落范圍持續(xù)增加。

不同開挖方案在不同開挖階段的沉降等值曲線變化，如圖7所示。不同開采方案沉降等值線占比變化規(guī)律相似：-7.78×10-7m等值線占比隨開采階段不斷增加；-8.2×10-7m等值線占比隨開采階段不斷降低；-8.6×10-7m等值線占比隨開采階段先增加后降低。回采方案三和回采方案二在開挖1/3礦房時沉降等值曲線占比變化最大，結(jié)合z方向位移場云圖，方案三的開挖礦房上方沉降變化更大，從沉降水平來看方案三在開挖前期頂板穩(wěn)定性最差。

圖7 沉降等值曲線占比Fig. 7 Proportion of settlement contour curve

5 結(jié) 論

(1)對比四種開挖方案在不同開挖階段最大主應(yīng)力極值可知：方案一在前期開挖階段未應(yīng)力集中，隨開挖繼續(xù)應(yīng)力迅速轉(zhuǎn)移并產(chǎn)生應(yīng)力集中，該回采方案造成的冒落危害相對較大；方案二的最大主應(yīng)力極值升高速率最快，該回采方案可在較小的暴露面積下發(fā)生冒落，冒落塊體及冒落的所帶來的能量沖擊較?。环桨溉捌诨夭擅奥滹L(fēng)險較小后期回采的冒落風(fēng)險較大，冒落塊體及冒落的所帶來的能量沖擊較大；方案四回采前期頂板巖層的穩(wěn)定性較好，隨回采的持續(xù)進行，有可能發(fā)生批量冒落甚至大規(guī)模冒落。

(2)回采礦房上方沉降幅度較小且沉降等值線主要為-7.78×10-7、-8.2×10-7、-8.6×10-7m，對比可以得知不同開挖方案在不同開挖階段的沉降等值曲線變化，回采方案三和回采方案二在開挖1/3礦房時沉降等值曲線占比變化最大，結(jié)合Z方向位移場云圖，方案三的開挖礦房上方沉降變化更大，在開挖前期頂板穩(wěn)定性最差。