基于Flac3D的某金矿回采方案优化数值模拟研究

王磊 闫刚 付佳杰 钟健

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

基金项目：国家自然科学基金青年项目（52004053）；中央高校基本科研业务专项资金资助-国防重大培育项目（N2101040）

作者简介：王  磊（1980—），男，工程师，从事矿山安全管理工作；E-mail：hdjkahk@163.com

*通信作者：钟  健（1999—），男，硕士研究生，从事地表岩层移动机理研究工作；E-mail：1184729590@qq.com

摘要：为了探究矿床开采顺序对地表沉陷的影响，以某金矿的生产实际问题为工程背景，拟定不同开采方案，并采用Flac3D软件数值模拟的方法，将地表监测的变形量及塑性区体积作为评价指标，分析不同回采方案下的地表移动规律。研究结果表明：2个中段同时开采产生的扰动会对地表移动有较大影响，第七中段与第八中段同时回采时，地表沉陷量最小；塑性区破坏体积与开采顺序有较大的关系，同时开采的中段总体积越小，矿山整体开采结束后被破坏的塑性区体积越小。研究结果为矿山安全高效开采提供参考依据。

关键词：充填采矿；回采顺序；数值模拟；地表移动；地表沉陷

中图分类号：TD853.34          文章编号：1001-1277（2024）06-0011-05

文献标志码：Adoi：10.11792/hj20240603

引  言

矿产资源开采引起的地表沉陷给后续的生产与安全带来了一系列消极影响，在开采过程中合理的回采—充填过程能够有效减缓地表的沉陷，不仅能够降低矿石回采的成本，更能提高矿山采矿效率与安全水平。国内外学者利用数值模拟的方式对回采顺序优化进行了大量研究，取得了丰硕的成果［1-8］。杨八九等［9］基于数值模拟结果发现部分开采区域有应力集中现象，安全隐患极大，通过多方案的采场结构参数及回采顺序论证，优化了后续开采方案，提高了开采安全系数。杨海彬等［10］针对某金矿深部开采地压显现问题，对不同矿房、矿柱跨度与分层采高的方案进行模拟，分别从应力、位移、塑性区进行分析，确定了最优采场结构参数优化方案。郭进平等［11］为分析四方金矿深部矿体由空场采矿法转崩落采矿法回采时的地压活动规律，对不同回采方案的位移、应力、塑性区分布规律进行了研究，确定了地压活动最小的最优回采方案。张瑞明等［12］在不同采场跨度和不同回采顺序条件下对采场进行了稳定性模拟分析，确定了最优的采场结构参数及回采顺序。吴振坤等［13］通过正交试验结果建立了多指标综合评价模型，在综合满意度最优条件下确定了最优回采顺序与采场跨度。陈顺满等［14］建立了多指标的响应面回归模型，基于整体满意度准则利用非线性目标优化算法对采场结构参数进行优化。GUO等［15］运用理论分析法对复杂倾斜薄矿体的采场结构参数进行了优化。BAI等［16］提出了网络流的方法，提高了采场结构参数的优化效率。

某金矿主要有Ⅰ、Ⅱ、Ⅲ号矿体，估算矿石资源储量达30.387 6万t，平均金品位1.96 g/t。矿体总体走向NE50°，倾向NW，倾角约30°，矿块平均长度126 m，平均水平厚度18.5 m，矿体总体呈大透镜状，赋存标高为-420～-380 m。根据矿体赋存特征、开采技术条件和矿山开采现状，采用上向水平分层充填采矿法回采矿体。根据矿体赋存情况布置矿房，矿块长20～40 m，宽为矿体厚度，高为中段高30～40 m，分段高10 m，分层高2.5 m，底柱高5 m，不留顶柱和间柱。爆破通风后进行顶板和上盘排险工作，局部不稳固地段采用锚杆金属网支护，锚杆长度一般为1.5～2.0 m，安装网度为1.0 m×1.0 m～0.5 m×0.5 m，与锚杆同时安装。

该金矿Ⅰ、Ⅱ号矿体赋存高度差较大且矿体走向相似，在开采设计时可一同考虑设计开采第一至第七中段。Ⅲ号矿体距Ⅰ、Ⅱ号矿体较远，赋存高度差较小，在回采时需单独考虑设计开采第八中段，Ⅲ号矿体中段可与Ⅰ、Ⅱ号矿体的其中一个中段同时回采。若回采过程中采动叠加影响过大，会影响地表建（构）筑物等的安全，而按照中段顺序进行逐中段开采会导致矿区生产能力小、成本高等问题。为降低采动叠加效应对地表的影响，保障矿山安全、高效生产，需对中段回采顺序进行研究。本文为探究Ⅲ号矿体中段的最佳开采顺序，利用Flac3D软件数值模拟分析矿体不同回采方案的回采—充填过程对地表移动规律的影响，对矿山回采顺序进行优化，以提高回采顺序的合理性，为安全、高效回采提供保证。

1  矿区模型建立及回采方案设计

1.1  计算区域与网格划分

根据岩体构造特点和矿区地形地貌，建立三维数值模型（见图1）。根据圣维南原理，围岩区域应取采空区范围的3～5倍来消除局部效应的影响［17］。本文三维数值模型研究的尺寸为长×宽×高=1 200 m×1 100 m×700 m。模型共分1 183 137个四面体单元，984 844个节点。

1.2  边界条件和初始条件

本次模拟边界条件采用位移约束，即模型的前后、左右边界施加水平约束，边界水平初始位移为零；模型的底面边界施加水平和垂直约束，边界水平和垂直初始位移为零；模型顶面边界为自由面，不施加约束条件。计算初始条件考虑自重应力场，重力加速度取9.81 m/s2，对模型进行计算，平衡后将位移场、速度场、塑性区清零。

1.3  岩体力学参数和本构模型

岩石是一种脆性材料，当荷载达到屈服强度后将发生破坏、弱化，属于弹塑性体。本次计算模拟的矿体及围岩本构模型采用莫尔-库仑弹塑性模型，其破坏准则为［18］：

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

式中：σ1、σ3分别为最大和最小主应力（MPa）；C、φ分别为内聚力（MPa）和内摩擦角（°）。

当fs>0时，材料将发生剪切破坏。在通常应力状态下，岩体的抗拉强度很低。因此，可根据抗拉强度准则判断岩体是否产生拉破坏。

数值模拟模型的矿岩及充填体物理力学参数见表1。

1.4  回采方案设计

为确定Ⅲ号矿体的最佳回采顺序，根据矿区矿体实际情况，采用自上而下的下行式开采，模拟Ⅲ号矿体第八中段分别与前7个中段同时回采的7种不同回采方案，见表2。各中段位置关系见图2。

2  回采方案优化分析

把地表监测的变形量及塑性区体积作为评价指标对各个方案进行评价，对不同模拟结果进行分析，计算结果见表3。

2.1  地表位移分析

不同方案地表总位移云图见图3。由图3可知：方案5的地表总位移最大，为44.847 mm，方案7的地表总位移最小，为43.350 mm。这说明，第五中段与第八中段同时开采对地表沉降的影响最大，第七中段与第八中段同时开采对地表沉降的影响最小。在实际生产中，每个中段开采后的采空区体积均不同，在模拟计算中，每个中段的横截面积设置为相同尺寸，仅采场长度不同。前7个中段所采矿体体积中第五中段体积最大，为8.4×105 m3；第七中段体积最小，为4.5×105 m3。在模拟方案中，2个中段同时开采产生的扰动会对地表移动有较大影响，其中体积最大的第五中段与体积最小的第七中段分别对应了地表沉陷量最大值与最小值。结合实际中段采空区体积，对比不同方案下的地表位移情况，方案7中第八中段与采空区体积最小的第七中段同时开采，产生的扰动影响最低，对地表的保护作用最好。

2.2  塑性区分析

不同方案中段采空区与塑性区体积分布见图4。由图4可知：塑性区体积与开采顺序有较大的关系。其中，方案5开采后的塑性区体积最大，为11.46×105 m3，方案7开采后的塑性区体积最小，为7.9×105 m3。同时，开采的中段体积越小，矿山整体开采结束后被破坏的塑性区体积越小。结合塑性区分布图（见图5），从整体上看，被破坏的塑性区体积较小，仅在第三中段产生了较多的塑性区域，没有出现大面积的塑性贯通现象。因此，可以认为矿山开采后的围岩稳定。为了减少第三中段被破坏的塑性区体积，在第三中段充填时，可以适当提高充填体强度，提高第三中段的充填体承载能力，保障后续中段的安全生产。

2.3  地表移动安全性综合评价

对于地表移动的安全性评价，反映地表变形沉降和移动特征的指标除了地表垂直沉降位移和水平移动位移，一般还有地表的倾斜变形、曲率、水平变形［19］。

1）地表倾斜变形。相邻点在竖直方向的相对移动量与相邻两点间水平距离的比值。它反映地表沿某一方向的坡度，通常以i表示。

iAB=WB-WAlAB=ΔWABlAB（2）

式中：iAB为倾斜值（mm/m）；WA、WB分别为地表观测点A、B的下沉值（mm）；lAB为地表观测点A、B点间的水平距离（m）。

2）地表曲率。两相邻线段的倾斜差和两线段中点间水平距离的比值。它反映了观测线断面上的弯曲程度，以K表示，地表下沉曲线上凸为正、下凹为负。

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

式中：KB为曲率（10-3/m）；iAB、iBC分别为地表A、B点间和B、C点间的平均斜率（mm/m）；lAB、lBC分别为地表A、B点间和B、C点间的水平距离（m）。

3）地表水平变形。相邻两点的水平移动差值与两点间水平距离的比值。它反映相邻两观测点间单位长度的水平移动差值，通常以ε表示。水平变形正值表示拉伸，分布在移动盆地拐点和边界；负值表示压缩，分布在移动盆地两拐点之间。

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

式中：εAB为水平变形值（mm/m）；UA、UB分别为地表观测点A、B的水平移动值（mm）。

对地表沉陷量最小及被破坏塑性区体积最小的方案7的地表监测点位移规律进行分析，结果见图6。由图6可知：开采充填后的地表最大倾斜变形i=-0.735 mm/m，最大水平变形ε=-0.316 mm/m，最大曲率K=0.037×10-3/m，均符合GB 50771—2012 《有色金属采矿设计规范》允许的地表变形值。

3  结  论

1）在模拟方案中，两个中段同时开采产生的扰动会对地表移动有较大影响，第七中段与第八中段同时回采时，地表沉陷量最小，对地表的保护作用最好。

2）塑性区体积与开采顺序有较大关系，同时开采的中段体积越小，矿山整体开采结束后被破坏的塑性区体积越小。为了减少被破坏的塑性区体积，在充填时可以适当提高充填体强度及承载能力。

3）基于最优回采顺序的模拟结果，得出地表最大倾斜变形i=-0.735 mm/m，最大水平变形ε=-0.316 mm/m，最大曲率K=0.037×10-3/m，均符合GB 50771—2012 《有色金属采矿设计规范》允许的地表变形值，地下开采活动对地表河流及建（构）筑物影响在安全范围内。

4）矿山采用优选出的第七、第八中段回采顺序方案，极大地降低了矿山开采过程中的安全风险，并提高了回采率，为矿山企业带来显著的经济效益。

［参 考 文 献］

［1］  姜永恒.镇沅金矿无底柱分段崩落采矿法采场结构参数的确定［J］.黄金，2022，43（5）：37-42.

［2］  刘滨，李振理，尹旭岩，等.上向水平分层充填采矿法在新城金矿的应用［J］.黄金，2018，39（1）：48-51.

［3］  崔晓林.三鑫金铜矿深部缓倾斜薄矿体回采工艺的优化选择［J］.黄金，2023，44（1）：33-37.

［4］  国绍林，宋召法，刘滨，等.踏步式上向水平分层充填采矿法回采步距优化及应用［J］.黄金，2019，40（1）：26-29.

［5］  周亚博.基于Flac3D数值模拟的深部采场结构参数优化及应用［J］.黄金，2022，43（5）：51-55.

［6］  高庆伟，孙嘉，郝显福.基于3D-σ有限元法的采场结构参数优化［J］.黄金，2017，38（6）：29-32.

［7］  王亚军，李向东，盛佳.前河金矿Ⅳ S3矿体280 m中段回采顺序优化［J］.矿业研究与开发，2018，38（5）：13-15.

［8］  关凯，朱万成，张洪训，等.新城金矿中段顶柱盘区回采顺序优化［J］.金属矿山，2013（9）：44-49.

［9］  杨八九，侯克鹏，孔学伟，等.基于数值模拟分析的采场应力阻隔及地压控制研究［J］.有色金属工程，2021，11（4）：118-124.

［10］  杨海彬，李文超，秦秀合，等.某金矿深部采场结构参数与回采顺序优化研究［J］.矿业研究与开发，2020，40（12）：18-22.

［11］  郭进平，强浩，刘少青，等.深部矿体无底柱分段崩落法回采顺序优化研究［J］.矿业研究与开发，2020，40（9）：11-17.

［12］  张瑞明，魏丁一，杜翠凤，等.阶段空场嗣后充填采矿法采场跨度及回采顺序优化研究［J］.矿业研究与开发，2022，42（3）：10-14.

［13］  吴振坤，金爱兵，陈帅军.基于多指标综合评价的采场结构参数及回采顺序优选［J］.矿业研究与开发，2020，40（2）：1-7.

［14］  陈顺满，吴爱祥，王贻明，等.基于响应面法的破碎围岩条件下采场结构参数优化研究［J］.岩石力学与工程学报，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］  杜坤，李夕兵，刘科伟，等.采空区危险性评价的综合方法及工程应用［J］.中南大学学报（自然科学版），2011，42（9）：2 802-2 811.

［18］  韩斌，吴爱祥，刘同有，等.金川二矿区多中段机械化盘区回采顺序的数值模拟优化研究［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，Fu 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