基于Flac3D的某金礦回采方案優(yōu)化數(shù)值模擬研究

王磊 閆剛 付佳杰 鐘健

收稿日期：2023-12-26； 修回日期：2024-02-15

基金項(xiàng)目：國(guó)家自然科學(xué)基金青年項(xiàng)目（52004053）；中央高?；究蒲袠I(yè)務(wù)專項(xiàng)資金資助-國(guó)防重大培育項(xiàng)目（N2101040）

作者簡(jiǎn)介：王? 磊（1980—），男，工程師，從事礦山安全管理工作；E-mail：hdjkahk@163.com

*通信作者：鐘? ?。?999—），男，碩士研究生，從事地表巖層移動(dòng)機(jī)理研究工作；E-mail：1184729590@qq.com

摘要：為了探究礦床開采順序?qū)Φ乇沓料莸挠绊?，以某金礦的生產(chǎn)實(shí)際問題為工程背景，擬定不同開采方案，并采用Flac3D軟件數(shù)值模擬的方法，將地表監(jiān)測(cè)的變形量及塑性區(qū)體積作為評(píng)價(jià)指標(biāo)，分析不同回采方案下的地表移動(dòng)規(guī)律。研究結(jié)果表明：2個(gè)中段同時(shí)開采產(chǎn)生的擾動(dòng)會(huì)對(duì)地表移動(dòng)有較大影響，第七中段與第八中段同時(shí)回采時(shí)，地表沉陷量最??；塑性區(qū)破壞體積與開采順序有較大的關(guān)系，同時(shí)開采的中段總體積越小，礦山整體開采結(jié)束后被破壞的塑性區(qū)體積越小。研究結(jié)果為礦山安全高效開采提供參考依據(jù)。

關(guān)鍵詞：充填采礦；回采順序；數(shù)值模擬；地表移動(dòng)；地表沉陷

中圖分類號(hào)：TD853.34????????? 文章編號(hào)：1001-1277（2024）06-0011-05

文獻(xiàn)標(biāo)志碼：Adoi：10.11792/hj20240603

引? 言

礦產(chǎn)資源開采引起的地表沉陷給后續(xù)的生產(chǎn)與安全帶來了一系列消極影響，在開采過程中合理的回采—充填過程能夠有效減緩地表的沉陷，不僅能夠降低礦石回采的成本，更能提高礦山采礦效率與安全水平。國(guó)內(nèi)外學(xué)者利用數(shù)值模擬的方式對(duì)回采順序優(yōu)化進(jìn)行了大量研究，取得了豐碩的成果［1-8］。楊八九等［9］基于數(shù)值模擬結(jié)果發(fā)現(xiàn)部分開采區(qū)域有應(yīng)力集中現(xiàn)象，安全隱患極大，通過多方案的采場(chǎng)結(jié)構(gòu)參數(shù)及回采順序論證，優(yōu)化了后續(xù)開采方案，提高了開采安全系數(shù)。楊海彬等［10］針對(duì)某金礦深部開采地壓顯現(xiàn)問題，對(duì)不同礦房、礦柱跨度與分層采高的方案進(jìn)行模擬，分別從應(yīng)力、位移、塑性區(qū)進(jìn)行分析，確定了最優(yōu)采場(chǎng)結(jié)構(gòu)參數(shù)優(yōu)化方案。郭進(jìn)平等［11］為分析四方金礦深部礦體由空?qǐng)霾傻V法轉(zhuǎn)崩落采礦法回采時(shí)的地壓活動(dòng)規(guī)律，對(duì)不同回采方案的位移、應(yīng)力、塑性區(qū)分布規(guī)律進(jìn)行了研究，確定了地壓活動(dòng)最小的最優(yōu)回采方案。張瑞明等［12］在不同采場(chǎng)跨度和不同回采順序條件下對(duì)采場(chǎng)進(jìn)行了穩(wěn)定性模擬分析，確定了最優(yōu)的采場(chǎng)結(jié)構(gòu)參數(shù)及回采順序。吳振坤等［13］通過正交試驗(yàn)結(jié)果建立了多指標(biāo)綜合評(píng)價(jià)模型，在綜合滿意度最優(yōu)條件下確定了最優(yōu)回采順序與采場(chǎng)跨度。陳順滿等［14］建立了多指標(biāo)的響應(yīng)面回歸模型，基于整體滿意度準(zhǔn)則利用非線性目標(biāo)優(yōu)化算法對(duì)采場(chǎng)結(jié)構(gòu)參數(shù)進(jìn)行優(yōu)化。GUO等［15］運(yùn)用理論分析法對(duì)復(fù)雜傾斜薄礦體的采場(chǎng)結(jié)構(gòu)參數(shù)進(jìn)行了優(yōu)化。BAI等［16］提出了網(wǎng)絡(luò)流的方法，提高了采場(chǎng)結(jié)構(gòu)參數(shù)的優(yōu)化效率。

某金礦主要有Ⅰ、Ⅱ、Ⅲ號(hào)礦體，估算礦石資源儲(chǔ)量達(dá)30.387 6萬t，平均金品位1.96 g/t。礦體總體走向NE50°，傾向NW，傾角約30°，礦塊平均長(zhǎng)度126 m，平均水平厚度18.5 m，礦體總體呈大透鏡狀，賦存標(biāo)高為-420～-380 m。根據(jù)礦體賦存特征、開采技術(shù)條件和礦山開采現(xiàn)狀，采用上向水平分層充填采礦法回采礦體。根據(jù)礦體賦存情況布置礦房，礦塊長(zhǎng)20～40 m，寬為礦體厚度，高為中段高30～40 m，分段高10 m，分層高2.5 m，底柱高5 m，不留頂柱和間柱。爆破通風(fēng)后進(jìn)行頂板和上盤排險(xiǎn)工作，局部不穩(wěn)固地段采用錨桿金屬網(wǎng)支護(hù)，錨桿長(zhǎng)度一般為1.5～2.0 m，安裝網(wǎng)度為1.0 m×1.0 m～0.5 m×0.5 m，與錨桿同時(shí)安裝。

該金礦Ⅰ、Ⅱ號(hào)礦體賦存高度差較大且礦體走向相似，在開采設(shè)計(jì)時(shí)可一同考慮設(shè)計(jì)開采第一至第七中段。Ⅲ號(hào)礦體距Ⅰ、Ⅱ號(hào)礦體較遠(yuǎn)，賦存高度差較小，在回采時(shí)需單獨(dú)考慮設(shè)計(jì)開采第八中段，Ⅲ號(hào)礦體中段可與Ⅰ、Ⅱ號(hào)礦體的其中一個(gè)中段同時(shí)回采。若回采過程中采動(dòng)疊加影響過大，會(huì)影響地表建（構(gòu)）筑物等的安全，而按照中段順序進(jìn)行逐中段開采會(huì)導(dǎo)致礦區(qū)生產(chǎn)能力小、成本高等問題。為降低采動(dòng)疊加效應(yīng)對(duì)地表的影響，保障礦山安全、高效生產(chǎn)，需對(duì)中段回采順序進(jìn)行研究。本文為探究Ⅲ號(hào)礦體中段的最佳開采順序，利用Flac3D軟件數(shù)值模擬分析礦體不同回采方案的回采—充填過程對(duì)地表移動(dòng)規(guī)律的影響，對(duì)礦山回采順序進(jìn)行優(yōu)化，以提高回采順序的合理性，為安全、高效回采提供保證。

1? 礦區(qū)模型建立及回采方案設(shè)計(jì)

1.1? 計(jì)算區(qū)域與網(wǎng)格劃分

根據(jù)巖體構(gòu)造特點(diǎn)和礦區(qū)地形地貌，建立三維數(shù)值模型（見圖1）。根據(jù)圣維南原理，圍巖區(qū)域應(yīng)取采空區(qū)范圍的3～5倍來消除局部效應(yīng)的影響［17］。本文三維數(shù)值模型研究的尺寸為長(zhǎng)×寬×高=1 200 m×1 100 m×700 m。模型共分1 183 137個(gè)四面體單元，984 844個(gè)節(jié)點(diǎn)。

1.2? 邊界條件和初始條件

本次模擬邊界條件采用位移約束，即模型的前后、左右邊界施加水平約束，邊界水平初始位移為零；模型的底面邊界施加水平和垂直約束，邊界水平和垂直初始位移為零；模型頂面邊界為自由面，不施加約束條件。計(jì)算初始條件考慮自重應(yīng)力場(chǎng)，重力加速度取9.81 m/s2，對(duì)模型進(jìn)行計(jì)算，平衡后將位移場(chǎng)、速度場(chǎng)、塑性區(qū)清零。

1.3? 巖體力學(xué)參數(shù)和本構(gòu)模型

巖石是一種脆性材料，當(dāng)荷載達(dá)到屈服強(qiáng)度后將發(fā)生破壞、弱化，屬于彈塑性體。本次計(jì)算模擬的礦體及圍巖本構(gòu)模型采用莫爾-庫侖彈塑性模型，其破壞準(zhǔn)則為［18］：

fs=σ1-σ31+sin φ1-sin φ-2C1+sin φ1-sin φ （1）

式中：σ1、σ3分別為最大和最小主應(yīng)力（MPa）；C、φ分別為內(nèi)聚力（MPa）和內(nèi)摩擦角（°）。

當(dāng)fs>0時(shí)，材料將發(fā)生剪切破壞。在通常應(yīng)力狀態(tài)下，巖體的抗拉強(qiáng)度很低。因此，可根據(jù)抗拉強(qiáng)度準(zhǔn)則判斷巖體是否產(chǎn)生拉破壞。

數(shù)值模擬模型的礦巖及充填體物理力學(xué)參數(shù)見表1。

1.4? 回采方案設(shè)計(jì)

為確定Ⅲ號(hào)礦體的最佳回采順序，根據(jù)礦區(qū)礦體實(shí)際情況，采用自上而下的下行式開采，模擬Ⅲ號(hào)礦體第八中段分別與前7個(gè)中段同時(shí)回采的7種不同回采方案，見表2。各中段位置關(guān)系見圖2。

2? 回采方案優(yōu)化分析

把地表監(jiān)測(cè)的變形量及塑性區(qū)體積作為評(píng)價(jià)指標(biāo)對(duì)各個(gè)方案進(jìn)行評(píng)價(jià)，對(duì)不同模擬結(jié)果進(jìn)行分析，計(jì)算結(jié)果見表3。

2.1? 地表位移分析

不同方案地表總位移云圖見圖3。由圖3可知：方案5的地表總位移最大，為44.847 mm，方案7的地表總位移最小，為43.350 mm。這說明，第五中段與第八中段同時(shí)開采對(duì)地表沉降的影響最大，第七中段與第八中段同時(shí)開采對(duì)地表沉降的影響最小。在實(shí)際生產(chǎn)中，每個(gè)中段開采后的采空區(qū)體積均不同，在模擬計(jì)算中，每個(gè)中段的橫截面積設(shè)置為相同尺寸，僅采場(chǎng)長(zhǎng)度不同。前7個(gè)中段所采礦體體積中第五中段體積最大，為8.4×105 m3；第七中段體積最小，為4.5×105 m3。在模擬方案中，2個(gè)中段同時(shí)開采產(chǎn)生的擾動(dòng)會(huì)對(duì)地表移動(dòng)有較大影響，其中體積最大的第五中段與體積最小的第七中段分別對(duì)應(yīng)了地表沉陷量最大值與最小值。結(jié)合實(shí)際中段采空區(qū)體積，對(duì)比不同方案下的地表位移情況，方案7中第八中段與采空區(qū)體積最小的第七中段同時(shí)開采，產(chǎn)生的擾動(dòng)影響最低，對(duì)地表的保護(hù)作用最好。

2.2? 塑性區(qū)分析

不同方案中段采空區(qū)與塑性區(qū)體積分布見圖4。由圖4可知：塑性區(qū)體積與開采順序有較大的關(guān)系。其中，方案5開采后的塑性區(qū)體積最大，為11.46×105 m3，方案7開采后的塑性區(qū)體積最小，為7.9×105 m3。同時(shí)，開采的中段體積越小，礦山整體開采結(jié)束后被破壞的塑性區(qū)體積越小。結(jié)合塑性區(qū)分布圖（見圖5），從整體上看，被破壞的塑性區(qū)體積較小，僅在第三中段產(chǎn)生了較多的塑性區(qū)域，沒有出現(xiàn)大面積的塑性貫通現(xiàn)象。因此，可以認(rèn)為礦山開采后的圍巖穩(wěn)定。為了減少第三中段被破壞的塑性區(qū)體積，在第三中段充填時(shí)，可以適當(dāng)提高充填體強(qiáng)度，提高第三中段的充填體承載能力，保障后續(xù)中段的安全生產(chǎn)。

2.3? 地表移動(dòng)安全性綜合評(píng)價(jià)

對(duì)于地表移動(dòng)的安全性評(píng)價(jià)，反映地表變形沉降和移動(dòng)特征的指標(biāo)除了地表垂直沉降位移和水平移動(dòng)位移，一般還有地表的傾斜變形、曲率、水平變形［19］。

1）地表傾斜變形。相鄰點(diǎn)在豎直方向的相對(duì)移動(dòng)量與相鄰兩點(diǎn)間水平距離的比值。它反映地表沿某一方向的坡度，通常以i表示。

iAB=WB-WAlAB=ΔWABlAB（2）

式中：iAB為傾斜值（mm/m）；WA、WB分別為地表觀測(cè)點(diǎn)A、B的下沉值（mm）；lAB為地表觀測(cè)點(diǎn)A、B點(diǎn)間的水平距離（m）。

2）地表曲率。兩相鄰線段的傾斜差和兩線段中點(diǎn)間水平距離的比值。它反映了觀測(cè)線斷面上的彎曲程度，以K表示，地表下沉曲線上凸為正、下凹為負(fù)。

KB=iBC-iAB12（lBC+lAB）=2ΔiA-B-ClBC+lAB（3）

式中：KB為曲率（10-3/m）；iAB、iBC分別為地表A、B點(diǎn)間和B、C點(diǎn)間的平均斜率（mm/m）；lAB、lBC分別為地表A、B點(diǎn)間和B、C點(diǎn)間的水平距離（m）。

3）地表水平變形。相鄰兩點(diǎn)的水平移動(dòng)差值與兩點(diǎn)間水平距離的比值。它反映相鄰兩觀測(cè)點(diǎn)間單位長(zhǎng)度的水平移動(dòng)差值，通常以ε表示。水平變形正值表示拉伸，分布在移動(dòng)盆地拐點(diǎn)和邊界；負(fù)值表示壓縮，分布在移動(dòng)盆地兩拐點(diǎn)之間。

εAB=UB-UAlAB=ΔUABlAB（4）

式中：εAB為水平變形值（mm/m）；UA、UB分別為地表觀測(cè)點(diǎn)A、B的水平移動(dòng)值（mm）。

對(duì)地表沉陷量最小及被破壞塑性區(qū)體積最小的方案7的地表監(jiān)測(cè)點(diǎn)位移規(guī)律進(jìn)行分析，結(jié)果見圖6。由圖6可知：開采充填后的地表最大傾斜變形i=-0.735 mm/m，最大水平變形ε=-0.316 mm/m，最大曲率K=0.037×10-3/m，均符合GB 50771—2012 《有色金屬采礦設(shè)計(jì)規(guī)范》允許的地表變形值。

3? 結(jié)? 論

1）在模擬方案中，兩個(gè)中段同時(shí)開采產(chǎn)生的擾動(dòng)會(huì)對(duì)地表移動(dòng)有較大影響，第七中段與第八中段同時(shí)回采時(shí)，地表沉陷量最小，對(duì)地表的保護(hù)作用最好。

2）塑性區(qū)體積與開采順序有較大關(guān)系，同時(shí)開采的中段體積越小，礦山整體開采結(jié)束后被破壞的塑性區(qū)體積越小。為了減少被破壞的塑性區(qū)體積，在充填時(shí)可以適當(dāng)提高充填體強(qiáng)度及承載能力。

3）基于最優(yōu)回采順序的模擬結(jié)果，得出地表最大傾斜變形i=-0.735 mm/m，最大水平變形ε=-0.316 mm/m，最大曲率K=0.037×10-3/m，均符合GB 50771—2012 《有色金屬采礦設(shè)計(jì)規(guī)范》允許的地表變形值，地下開采活動(dòng)對(duì)地表河流及建（構(gòu)）筑物影響在安全范圍內(nèi)。

4）礦山采用優(yōu)選出的第七、第八中段回采順序方案，極大地降低了礦山開采過程中的安全風(fēng)險(xiǎn)，并提高了回采率，為礦山企業(yè)帶來顯著的經(jīng)濟(jì)效益。

［參 考 文 獻(xiàn)］

［1］? 姜永恒.鎮(zhèn)沅金礦無底柱分段崩落采礦法采場(chǎng)結(jié)構(gòu)參數(shù)的確定［J］.黃金，2022，43（5）：37-42.

［2］? 劉濱，李振理，尹旭巖，等.上向水平分層充填采礦法在新城金礦的應(yīng)用［J］.黃金，2018，39（1）：48-51.

［3］? 崔曉林.三鑫金銅礦深部緩傾斜薄礦體回采工藝的優(yōu)化選擇［J］.黃金，2023，44（1）：33-37.

［4］? 國(guó)紹林，宋召法，劉濱，等.踏步式上向水平分層充填采礦法回采步距優(yōu)化及應(yīng)用［J］.黃金，2019，40（1）：26-29.

［5］? 周亞博.基于Flac3D數(shù)值模擬的深部采場(chǎng)結(jié)構(gòu)參數(shù)優(yōu)化及應(yīng)用［J］.黃金，2022，43（5）：51-55.

［6］? 高慶偉，孫嘉，郝顯福.基于3D-σ有限元法的采場(chǎng)結(jié)構(gòu)參數(shù)優(yōu)化［J］.黃金，2017，38（6）：29-32.

［7］? 王亞軍，李向東，盛佳.前河金礦Ⅳ S3礦體280 m中段回采順序優(yōu)化［J］.礦業(yè)研究與開發(fā)，2018，38（5）：13-15.

［8］? 關(guān)凱，朱萬成，張洪訓(xùn)，等.新城金礦中段頂柱盤區(qū)回采順序優(yōu)化［J］.金屬礦山，2013（9）：44-49.

［9］? 楊八九，侯克鵬，孔學(xué)偉，等.基于數(shù)值模擬分析的采場(chǎng)應(yīng)力阻隔及地壓控制研究［J］.有色金屬工程，2021，11（4）：118-124.

［10］? 楊海彬，李文超，秦秀合，等.某金礦深部采場(chǎng)結(jié)構(gòu)參數(shù)與回采順序優(yōu)化研究［J］.礦業(yè)研究與開發(fā)，2020，40（12）：18-22.

［11］? 郭進(jìn)平，強(qiáng)浩，劉少青，等.深部礦體無底柱分段崩落法回采順序優(yōu)化研究［J］.礦業(yè)研究與開發(fā)，2020，40（9）：11-17.

［12］? 張瑞明，魏丁一，杜翠鳳，等.階段空?qǐng)鏊煤蟪涮畈傻V法采場(chǎng)跨度及回采順序優(yōu)化研究［J］.礦業(yè)研究與開發(fā)，2022，42（3）：10-14.

［13］? 吳振坤，金愛兵，陳帥軍.基于多指標(biāo)綜合評(píng)價(jià)的采場(chǎng)結(jié)構(gòu)參數(shù)及回采順序優(yōu)選［J］.礦業(yè)研究與開發(fā)，2020，40（2）：1-7.

［14］? 陳順滿，吳愛祥，王貽明，等.基于響應(yīng)面法的破碎圍巖條件下采場(chǎng)結(jié)構(gòu)參數(shù)優(yōu)化研究［J］.巖石力學(xué)與工程學(xué)報(bào)，2017，36（增刊1）：3 499-3 508.

［15］? GUO Y，HOU K，LI W.Numerical optimization of stope structural parameters of complex inclined thin orebody in Mengnuo lead-zinc mine［J］.Electronic Journal of Geotechnical Engineering，2014，19：9 479-9 490.

［16］? BAI X，MARCOTTE D，SIMON R.Underground stope optimization with network flow method［J］.Computers & Geosciences，2013，52（3）：361-371.

［17］? 杜坤，李夕兵，劉科偉，等.采空區(qū)危險(xiǎn)性評(píng)價(jià)的綜合方法及工程應(yīng)用［J］.中南大學(xué)學(xué)報(bào)（自然科學(xué)版），2011，42（9）：2 802-2 811.

［18］? 韓斌，吳愛祥，劉同有，等.金川二礦區(qū)多中段機(jī)械化盤區(qū)回采順序的數(shù)值模擬優(yōu)化研究［J］.礦冶工程，2004（2）：4-7.

［19］? LI H，GUO G，ZHAI S C.Mining scheme design for super-high water backfill strip mining under buildings：A Chinese case study［J］.Environmental Earth Sciences，2016，75（12）：1 017.

Numerical experimental study on optimization

of mining scheme in a gold mine based on Flac3D

Wang Lei1，Yan Gang1，F(xiàn)u Jiajie1，Zhong Jian2

（1.Canzhuang Gold Mine，Zhaojin Mining Industry Co.，Ltd.;

2.School of Resources and Civil Engineering，Northeastern University）

Abstract：In order to explore the effect of mining sequence on surface subsidence，taking the practical production problems of a gold mine as the engineering background，different mining schemes were designed.The Flac3D software was used for numerical simulation，with surface deformation magnitude and plastic zone volume monitored as evaluation indicators to analyze the movement patterns of the surface under different mining schemes.The research results indicate that disturbances generated by simultaneous mining of two levels have a significant impact on surface movement，and the smallest surface subsidence occurs when the seventh and eighth levels are simultaneously mined.The volume of plastic zone damaged is closely related to the mining sequence，and the smaller the total volume of levels mined simul-taneously，the smaller the volume of plastic zones damaged after the completion of overall mining in the mine.The research findings provide a reference basis for safe and efficient mining in the mine.

Keywords：filling mining;mining sequence;numerical simulation;surface movement;surface subsidence