劉曉明 ，羅周全，徐紀(jì)成

(1. 中南大學(xué) 資源與安全工程學(xué)院，湖南 長沙，410083；2. 中國礦業(yè)大學(xué) 煤炭資源與安全開采國家重點實驗室，江蘇 徐州，221116；3. 中南大學(xué) 現(xiàn)代分析測試中心，湖南 長沙，410083)

隱患空區(qū)動力失穩(wěn)過程是一個極其復(fù)雜的力學(xué)系統(tǒng)，但從局部而言，都是巖石的變形與破裂問題。隨著巖石力學(xué)的發(fā)展，各種數(shù)值計算方法如有限元法、有限差分法、邊界元法、離散元法、拉格朗日元法和流形元法等廣泛應(yīng)用于各種巖石力學(xué)工程問題中[1?4]。而巖層移動與破裂研究的數(shù)值計算方法主要包括 2類：一類是建立在連續(xù)介質(zhì)力學(xué)基礎(chǔ)上的數(shù)值計算方法，如有限元法、邊界單元法等；另一類則是非連續(xù)介質(zhì)力學(xué)分析方法，如離散元法。近年來，細(xì)觀力學(xué)分析方法為巖石破裂過程的數(shù)值模擬方法的發(fā)展提供了很好的思路。巖石破裂過程分析系統(tǒng) RFPA(rock failure process analysis)是基于彈性損傷統(tǒng)計模型的數(shù)值模擬工具，能較好地模擬巖石破裂失穩(wěn)過程，目前已廣泛應(yīng)用于巖石破裂機(jī)制以及巖石力學(xué)與工程問題的研究，取得了許多成果[5?6]。隱患空區(qū)失穩(wěn)根據(jù)失穩(wěn)類型可以分為靜力失穩(wěn)和動力失穩(wěn)。靜力失穩(wěn)即在原巖應(yīng)力場作用下進(jìn)行空區(qū)開挖，由于空區(qū)周圍應(yīng)力集中致使空區(qū)頂板或兩幫出現(xiàn)拉伸破壞或剪切破壞；動力失穩(wěn)指在靜力作用的基礎(chǔ)上，受到周圍爆破動荷載的作用，致使空區(qū)失穩(wěn)。與靜力失穩(wěn)相比，空區(qū)動力失穩(wěn)模式雖然更復(fù)雜，仿真難度更大，但是，更能反映實際空區(qū)破壞過程。此外，隱患空區(qū)形態(tài)邊界復(fù)雜，如采用矩形、梯形或圓等規(guī)則形狀進(jìn)行簡化，再用簡化的模型進(jìn)行空區(qū)動力失穩(wěn)仿真，空區(qū)的破壞過程與規(guī)律往往與實際不符，不能如實反映實際復(fù)雜隱患空區(qū)動力失穩(wěn)過程，因此，有必要繼承隱患空區(qū)三維信息獲取的結(jié)果，在此基礎(chǔ)上進(jìn)行動力失穩(wěn)仿真模型構(gòu)建。華錫集團(tuán)銅坑礦92號礦體的開采技術(shù)條件復(fù)雜，井下隱患空區(qū)空間分布復(fù)雜，而爆破地震波作用下空區(qū)發(fā)生動力失穩(wěn)的趨勢日益明顯。為保障92號礦體的安全高效開采，揭示空區(qū)動力失穩(wěn)的內(nèi)在機(jī)制，有效預(yù)防空區(qū)塌陷事故的發(fā)生，在銅坑礦隱患空區(qū) CMS(cavity monitoring system)探測并獲取準(zhǔn)確的三維信息的基礎(chǔ)上，借助RFPA提供的巖石破裂過程動力仿真分析功能，同時考慮開挖過程和爆破震動作用對空區(qū)穩(wěn)定性的影響，對空區(qū)動力失穩(wěn)過程進(jìn)行仿真，以便為開展空區(qū)災(zāi)害預(yù)防與控制、空區(qū)處理方案的制定、應(yīng)力監(jiān)測點布置、資源開采布局等工作提供基礎(chǔ)依據(jù)。

1 隱患空區(qū)三維信息獲取

1.1 工程概況

柳州華錫集團(tuán)銅坑礦主要開采對象為細(xì)脈帶、91號和92號礦體，這三大礦體在空間位置上重疊產(chǎn)出。目前細(xì)脈帶、91號等礦體開采即將結(jié)束，正在大量開采92號礦體及細(xì)脈帶及91號礦體邊緣薄貧礦體。在長期開采中形成了大量的采空區(qū)，局部地段空區(qū)的冒落或坍塌已構(gòu)成威脅礦山生產(chǎn)的重大安全隱患。為能采取有效措施控制因空區(qū)引發(fā)的地質(zhì)災(zāi)害，急需對現(xiàn)有空區(qū)進(jìn)行徹底探查和分析，以掌握空區(qū)空間分布、實際邊界、體積等基本情況，并在此基礎(chǔ)上對空區(qū)群危險性進(jìn)行分析，確定不穩(wěn)定或隨著周圍礦體的開采將不穩(wěn)定的空區(qū)，以便為空區(qū)綜合治理與控制、資源有效回收提供技術(shù)支撐[7?8]。

1.2 隱患空區(qū)三維信息獲取

空區(qū)激光探測系統(tǒng)(CMS)因其能準(zhǔn)確測定井下危險采空區(qū)及其他硐室空間的三維形態(tài)而被廣泛應(yīng)用于國內(nèi)外金屬地下礦山[9?12]。因此，研究采用CMS對銅坑礦92號礦體隱患空區(qū)進(jìn)行探測，輸出的點云數(shù)據(jù)經(jīng)過濾、拼接及坐標(biāo)、格式轉(zhuǎn)換，采用三維礦業(yè)軟件Surpac生成空區(qū)三維實體模型，并以此為基礎(chǔ)開展模型計算、三維可視化及剖面輸出等應(yīng)用。隱患空區(qū)三維信息獲取及可視化過程如圖1所示[13]。

圖1 隱患空區(qū)三維信息獲取及可視化過程Fig.1 Process of 3D information acquisition of disaster cavity and its visualization

隱患空區(qū)三維空間信息主要包括空區(qū)位置形態(tài)、體積、頂板面積及邊界等，這些信息將隨時間發(fā)生變化，因此，準(zhǔn)確、及時掌握隱患空區(qū)信息可為礦山進(jìn)行采礦貧損指標(biāo)控制、開采設(shè)計、空區(qū)管理、安全評價及災(zāi)害預(yù)測等提供可靠的基礎(chǔ)性依據(jù)。

1.3 計算區(qū)域的選取及其特征分析

通過分析銅坑礦 92號礦體探測空區(qū)群的三維空間信息，發(fā)現(xiàn)92號礦體西南采區(qū)內(nèi)的空區(qū)之間互相重疊、貫通，穩(wěn)定性較差，鄰近頻繁爆破生產(chǎn)作業(yè)誘發(fā)采空區(qū)失穩(wěn)破壞的趨勢日益明顯。為此，結(jié)合目前銅坑礦的開采現(xiàn)狀，選擇92號礦體西南采區(qū)內(nèi)的采空區(qū)及其周邊礦巖作為本次動力失穩(wěn)仿真技術(shù)研究的計算區(qū)域。

因空區(qū)三維形態(tài)復(fù)雜，構(gòu)建復(fù)雜空區(qū)三維數(shù)值計算模型難度大，為此，選用二維 RFPA，選擇典型空區(qū)剖面進(jìn)行分析，相比三維數(shù)值分析，二維計算結(jié)果安全性更高。因空區(qū)主要沿勘探線分布，為選取典型空區(qū)進(jìn)行動力失穩(wěn)過程模擬仿真，分別沿 205-2，203-2，201-2和202-2這4條勘探線截取剖面，分析剖面上空區(qū)動力失穩(wěn)的過程，其中剖面位置如圖2所示。不同剖面上隱患空區(qū)的分布特征見表1。

圖2 隱患空區(qū)及計算剖面位置Fig.2 Disaster cavity and section position

2 隱患空區(qū)RFPA數(shù)值計算模型

2.1 計算模型

目前，RFPA-dynamic尚未開發(fā)與其他三維建模軟件或數(shù)值軟件數(shù)據(jù)轉(zhuǎn)換接口，只能在軟件內(nèi)建立計算模型，軟件提供了完整的用于建立各種直線、曲線、不規(guī)則多邊形、矩形和圓等工具，方便建立各種邊界不規(guī)則的隱患空區(qū)。因此，充分結(jié)合205-2，203-2，201-2和202-2這4條勘探線剖面上空區(qū)形態(tài)邊界，在RFPA內(nèi)利用不規(guī)則多邊形完成復(fù)雜邊界隱患空區(qū)數(shù)值計算模型的建立，綜合考慮模型尺寸及軟件計算能力，針對不同模型選取不同的尺寸及單元尺寸，具體參數(shù)見表2，最后建立的RFPA計算模型如圖3所示。

表1 不同剖面上隱患空區(qū)的分布特征Table 1 Distribute characteristic in different sections

表2 RFPA數(shù)值計算模型參數(shù)Table 2 RFPA numerical calculation model parameters

2.2 基元賦值

RFPA通過 Monte-Carlo方法和統(tǒng)計描述相結(jié)合對基元進(jìn)行初始化賦值[14?15]。設(shè)網(wǎng)格中所有基元的彈性模量最大閾值為E0，φ(E)為具有某彈性模量E基元的分布值，則彈性模量Weibull分布函數(shù)的積分為

其中：φ(E)為具有彈性模量E的基元的統(tǒng)計量。由上式統(tǒng)計分布構(gòu)成的基元組成一個樣本空間，在閾值E0不變的情況下，由于m的差別，積分空間分布不完全一樣。這些基元構(gòu)成的巖石類介質(zhì)在宏觀上性質(zhì)可能大體一致(E0相同)，但是，由于細(xì)觀結(jié)構(gòu)的無序性，使得基元的空間排列方式有顯著不同。這種細(xì)觀上的無序性正好體現(xiàn)了巖石類介質(zhì)獨特的離散性特征。

圖3 銅坑礦隱患空區(qū)RFPA數(shù)值計算模型Fig.3 RFPA numerical calculation models of Tongkeng Mine disaster cavity

一般物理空間隨機(jī)分布的無序性可以通過Monte-Carlo方法來實現(xiàn)，其產(chǎn)生方法是：產(chǎn)生一組在(0, 1)區(qū)間上均勻分布的隨機(jī)數(shù)序列{iγ≤1|i=1, 2, …,n}；對于任何iγ，存在一個與{iγ≤1|i=1, 2, …,n}相對應(yīng)的隨機(jī)數(shù)序列{Ei|i=1, 2, …,n}。由隨機(jī)數(shù)序列{iγ}映射一組彈性模量參數(shù)序列{iγ}。這一組基元彈性模量參數(shù)隨機(jī)序列逐一賦予網(wǎng)絡(luò)中的每一個基元，其他力學(xué)參數(shù)(強(qiáng)度、泊松比、容重)同樣用此種方法賦值。這種既具有統(tǒng)計性、又具有隨機(jī)性(無序性)的基元力學(xué)參數(shù)賦值方法滿足了網(wǎng)格中基元力學(xué)參數(shù)非均勻性和隨機(jī)性的要求。

根據(jù)銅坑礦92號礦體的礦巖的基本條件，采用的基元參數(shù)見表3。

表3 巖體物理力學(xué)參數(shù)Table 3 Physical and mechanical parameters of rock mass

2.3 初始應(yīng)力與邊界條件

靜力計算邊界條件為：根據(jù)銅坑礦的原巖應(yīng)力場分布，采用上部受壓、側(cè)向受原巖應(yīng)力場作用的平面應(yīng)變模型，原巖應(yīng)力場的水平方向應(yīng)力為25.4 MPa，垂直方向應(yīng)力根據(jù)模型邊界以上巖石自重求取。

動力計算邊界條件為：在靜力計算邊界條件的基礎(chǔ)上，在空區(qū)底板設(shè)置爆破地震波。

2.4 爆破震動荷載輸入

將現(xiàn)場實測爆破地震波時間?速度數(shù)據(jù)進(jìn)行簡化，并轉(zhuǎn)換為RFPA能夠接受的應(yīng)力波應(yīng)力，如圖4所示。

圖4 爆破震動應(yīng)力時程曲線Fig.4 Schedule curve of blasting vibration stress

在計算軟件中選擇靜力和動力計算模塊，先設(shè)置靜力加載條件，再設(shè)置動力加載條件。其中動力加載可以通過應(yīng)力、位移、波速和波加速度這4種方式實現(xiàn)。本文選擇應(yīng)力加載方式，在空區(qū)底部施加動力荷載。

3 隱患空區(qū)動力失穩(wěn)過程模擬仿真

隱患空區(qū)動力失穩(wěn)仿真的整個計算過程分為靜力計算和動力計算。動力分析以靜力計算為前提，即在完成靜力計算的基礎(chǔ)上，再施加動力時程進(jìn)行計算。靜力分析時，首先進(jìn)行初始應(yīng)力平衡計算，然后，根據(jù)空區(qū)形成次序分步模擬空區(qū)形成過程，每步形成一個空區(qū)并計算至平衡狀態(tài)。動力計算在靜力計算結(jié)束后進(jìn)行，通過在空區(qū)底板施加爆破震動波，再進(jìn)行動力計算，直至空區(qū)破壞，最終獲得空區(qū)從形成到失穩(wěn)破壞的過程模擬仿真。

3.1 勘探線205-2上的空區(qū)

圖5所示為勘探線205-2空區(qū)動力擾動下失穩(wěn)破壞過程的剪應(yīng)力分布。從圖5可見：隨著T506空區(qū)的開挖，空區(qū)圍巖應(yīng)力發(fā)生重分布，T505空區(qū)形成后應(yīng)力進(jìn)一步增大，尤其表現(xiàn)在間柱上；施加爆破動載后，空區(qū)圍巖受力狀態(tài)有所改變，首先在 T506空區(qū)右?guī)统霈F(xiàn)裂紋；隨著應(yīng)力波的逐步向上傳播，裂紋不斷往 T505空區(qū)方向擴(kuò)展。由于空區(qū)間柱最小寬度僅為5.5 m，裂紋最先在該位置貫通，致使空區(qū)間柱完全破壞。

圖5 勘探線205-2剖面空區(qū)動力失穩(wěn)過程仿真Fig.5 Dynamic destabilization numerical simulations of cavity in section 205-2

3.2 勘探線201-2上的空區(qū)

圖6所示為勘探線201-2空區(qū)動力擾動下失穩(wěn)破壞過程的剪應(yīng)力分布。從圖6可見：

(1) 由于S1空區(qū)與其他空區(qū)間距高達(dá)82 m，因此，其他空區(qū)的開挖對其影響較小。動荷載作用后，空區(qū)左右兩幫周圍圍巖出現(xiàn)局部裂紋，但并未擴(kuò)展。

(2) T111與T112空間上呈現(xiàn)上下斜交，間柱最小寬度僅為6 m，靜力作用下隨著空區(qū)的逐步開挖，空區(qū)圍巖應(yīng)力發(fā)生重分布，間柱出現(xiàn)應(yīng)力集中。

(3) 動荷載作用后，T112空區(qū)右?guī)褪紫瘸霈F(xiàn)局部片幫跨落現(xiàn)象，隨著爆破應(yīng)力波的進(jìn)一步作用，片幫范圍逐步加大，并有往下擴(kuò)展的趨勢。

(4) 動荷載進(jìn)一步作用后，T110-111頂板靠T112處出現(xiàn)裂紋，擴(kuò)展與 T112片幫區(qū)域連通，最后，整個間柱出現(xiàn)斷裂破壞。

圖6 勘探線201-2剖面空區(qū)動力失穩(wěn)過程仿真Fig.6 Dynamic destabilization numerical simulation of cavity in section 201-2

3.3 勘探線202-2上的空區(qū)

圖7所示為勘探線202-2空區(qū)動力擾動下失穩(wěn)破壞過程的剪應(yīng)力分布。從圖7可見：

(1) 隨著空區(qū)的逐步開挖，空區(qū)圍巖應(yīng)力發(fā)生重分布，尤其在空區(qū)開挖第3步后，由于空區(qū)分布密集，空區(qū)間柱跨度小，使得空區(qū)間柱應(yīng)力集中。

(2) 空區(qū)開挖第 4步后，空區(qū)間柱應(yīng)力集中進(jìn)一步加強(qiáng)，并出現(xiàn)局部破壞，如T214與T213，T215與T216的間柱均不同程度出現(xiàn)裂紋，此時，裂紋并未擴(kuò)散。

(3) 動荷載作用后，應(yīng)力集中程度得到加強(qiáng)，間柱裂紋出現(xiàn)擴(kuò)展，并沿著間柱薄弱地帶發(fā)展，致使T214與T213間柱首先發(fā)生破壞。

(4) 動荷載進(jìn)一步作用后，T215與T216以及T211與 T212的間柱依次出現(xiàn)破壞，致使整個空區(qū)向下垮塌。

圖7 勘探線202-2剖面空區(qū)動力失穩(wěn)過程仿真Fig.7 Dynamic destabilization numerical simulation of cavity in section 202-2

3.4 勘探線203-2上的空區(qū)

圖8所示為勘探線203-2空區(qū)動力擾動下失穩(wěn)破壞過程的剪應(yīng)力分布。從圖8可見：

圖8 勘探線203-2剖面空區(qū)動力失穩(wěn)過程仿真Fig.8 Dynamic destabilization numerical simulation of cavity in section 203-2

(1) 模型中T214-316與T312呈上下疊加式，即T314，T315和T316空區(qū)的頂板為T312空區(qū)的底板，相交處為最危險處，容易出現(xiàn)斷裂。

(2) 空區(qū)開挖第 3步后，空區(qū)間柱出現(xiàn)明顯的應(yīng)力集中現(xiàn)象，T314，T315和T316空區(qū)頂板出現(xiàn)局部破壞。

(3) 動荷載作用后，空區(qū)圍巖應(yīng)力出現(xiàn)二次分布。T312空區(qū)靠T314底板出現(xiàn)裂紋，隨著動荷載的作用，T312空區(qū)右?guī)统霈F(xiàn)局部破壞。

(4) 動荷載進(jìn)一步作用后，原有裂紋繼續(xù)擴(kuò)展，最后空區(qū)間柱全部垮塌破壞。

4 結(jié)論

(1) 在研究 RFPA基本思想和巖石破裂過程仿真原理的基礎(chǔ)上，提出了運用RFPA開展隱患空區(qū)動力失穩(wěn)模擬仿真的方法。

(2) 在對銅坑礦隱患空區(qū)現(xiàn)場實測并準(zhǔn)確獲取相關(guān)信息的基礎(chǔ)上，結(jié)合隱患空區(qū)危險性可視化分級結(jié)果，對處于205-2，203-2，201-2和202-2這4條勘探線上的空區(qū)空間分布特征進(jìn)行了分析，為空區(qū)動力失穩(wěn)過程模擬仿真提供便利。

(3) 模擬了空區(qū)從形成到動力擾動作用下的破壞過程，獲得了空區(qū)破壞過程的裂紋擴(kuò)展規(guī)律，形成了金屬礦隱患空區(qū)動力失穩(wěn)過程模擬方法。

[1] 何滿潮, 李學(xué)元, 劉斌, 等. 深部開采巖體力學(xué)研究[J]. 巖石力學(xué)與工程學(xué)報, 2005, 24(16): 2803?2813.HE Manchao, LI Xueyuan, LIU Bin, et al. Study on rock mechanics in deep mining engineering[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(16): 2803?2813.

[2] 林杭, 曹平, 李江騰, 等. 基于Surpac的Flac3D三維模型自動構(gòu)建[J]. 中國礦業(yè)大學(xué)學(xué)報, 2008, 37(3): 339?342.LIN Hang, CAO Ping, LI Jiangteng, et al. Automatic generation of Flac3Dmodel based on Surpac[J]. Journal of China University of Mining & Technology, 2008, 37(3): 339?342.

[3] 王濤, 陳曉玲, 楊建. 基于3DGIS和3DEC的地下洞室圍巖穩(wěn)定性研究[J]. 巖石力學(xué)與工程學(xué)報, 2005, 24(19): 3476?3481.WANG Tao, CHEN Xiaoling, YANG Jian. Study on stability of underground cavern based on 3DGIS and 3DEC[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(19):3476?3481.

[4] Tang C A, Tang S B. Applications of rock failure process analysis (RFPA) method[J]. Journal of Rock Mechanics and Geotechnical Engineering, 2011, 3(4): 352?372.

[5] 朱萬成, 唐春安, 黃志平, 等. 靜態(tài)和動態(tài)載荷作用下巖石劈裂破壞模式的數(shù)值模擬[J]. 巖石力學(xué)與工程學(xué)報, 2005, 24(1):1?7.ZHU Wancheng, TANG Chunan, HUANG Zhiping, et al.Numerical simulation on splitting failure mode of rock under static and dynamic loadings[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(1): 1?7.

[6] ZHU Wancheng, TANG Chunan. Numerical simulation on shear fracture process of concrete using mesoscopic mechanical model[J]. Constructions and Building Materials, 2002, 16(8):453?463.

[7] 高峰. 頂板誘導(dǎo)崩落機(jī)理及次生災(zāi)變鏈?zhǔn)叫?yīng)控制研究[D].長沙: 中南大學(xué)資源與安全工程學(xué)院, 2009: 20?24.GAO Feng. Study on mechanism of induced caving roof and control of its secondary catastrophe chain effect[D]. Changsha:Central South University. School of Resources and Safety Engineering, 2009: 20?24.

[8] 蘇家紅. 大范圍隱患區(qū)資源開采安全控制技術(shù)研究[D]. 長沙:中南大學(xué)資源與安全工程學(xué)院, 2008: 1?2.SU Jiahong. Study on safety control to mining under large scale hidden trouble zone[D]. Changsha: Central South University.School of Resources and Safety Engineering, 2008: 1?2.

[9] LUO Zhouquan, LIU Xiaoming, ZHANG Bao, et al. Cavity 3D modeling and correlative techniques based on cavity monitoring[J]. Journal of Central South University of Technology, 2008, 15(5): 639?644.

[10] LIU Xiling, LI Xibing, LI Faben, et al. 3D cavity detection technique and its application based on cavity auto scanning laser system[J]. Journal of Central South University of Technology,2008, 15(4): 285?288.

[11] LIU Xiaoming, LUO Zhouquan, YUAN Wenni, et al. CMS applications in metal mines of China[C]// Proceedings of Pacific-Asia Workshop on Computational Intelligence and Industrial Application. Wuhan, China, 2009: 353?356.

[12] 劉曉明, 羅周全, 孟穩(wěn)權(quán), 等. 深井采場大規(guī)模垮塌三維探測及可視化計算[J]. 中南大學(xué)學(xué)報: 自然科學(xué)版, 2011, 42(1):158?163.LIU Xiaoming, LUO Zhouquan, MENG Wenquan, et al. 3D monitoring and visible calculation of large-scale collapse area in deep mine[J]. Journal of Central South University: Science and Technology, 2011, 42(1): 158?163.

[13] 劉曉明, 羅周全, 袁雯妮, 等. 基于 CMS的隱患空區(qū)三維特征信息獲取[J]. 科技導(dǎo)報, 2011, 29(5): 32?36.LIU Xiaoming, LUO Zhouquan, YUAN Wenni, et al. 3D information acquisition of disaster cavity based on CMS[J].Science & Technology Review, 2011, 29(5): 32?36.

[14] 謝林茂, 朱萬成, 王述紅, 等. 含孔洞巖石試樣三維破裂過程的并行計算分析[J]. 巖土工程學(xué)報, 2011, 33(9): 1447?1455.XIE Linmao, ZHU Wancheng, WANG Shuhong, et al.Three-dimensional parallel computing on failure process of rock specimen with a pre-existing circular opening[J]. Chinese Journal of Geotechnical Engineering, 2011, 33(9): 1447?1455.

[15] 黃志平, 唐春安, 馬天輝, 等. 卸載巖爆過程數(shù)值試驗研究[J].巖石力學(xué)與工程學(xué)報, 2011, 30 (S1): 3120?3128.HUANG Zhiping, TANG Chunan, MA Tianhui, et al. Numerical test investigation on unloading rockburst processes[J]. Chinese Journal of Rock Mechanics and Engineering, 2011, 30(S1):3120?3128.