深井大规模进路充填开采扰动效应影响分析

吴昊燕,张爱民,王志远(中国恩菲工程技术有限公司,北京 100038)

1 前言

随着计算机技术的迅速发展,数值模拟凭借高效、安全、低成本的优势在矿山开采设计、工程研究中的地位愈加突出。目前广泛应用的采矿岩石力学的数值分析方法主要包括有限元法、边界元法、有限差分法、离散元法等。基于有限差分法FLAC3D的矿山开采数值模拟[1-9]为解决地表沉降、采场结构参数优化、采场稳定性分析、矿柱设置等工程和科学问题提供了重要的研究手段。这些研究一般针对某个问题建立局部区域模型,或简化矿体模型、回采顺序进行模拟,与实际工程存在一定的出入。本文通过建立深井大规模进路充填开采矿山大模型,对比了三种不同开采模拟方式的扰动效应。

2 FLAC3D数值计算

FLAC3D(Fast Lagrangian Analysis of Continua,连续介质快速拉格朗日分析)是由Itasca国际集团公司开发出的有限差分数值计算程序。软件可以计算岩土体在各种外荷载作用下产生的变形、应力、稳定性,尤其擅长计算岩土体破坏后的大变形和峰后特性等问题。

利用FLAC3D进行矿山开采岩石力学数值分析的主要步骤包括前处理、过程计算和后处理:

(1)建立三维岩体模型。根据地质模型及设计开采范围构建矿体、围岩、断层等岩体模型。

(2)网格划分。结合采场结构参数及回采顺序对关键区域进行网格加密。

(3)材料赋予。对不同分组的矿岩体进行物理力学参数赋予。

(4)条件输入。包括初始地应力条件参数、边界条件设置、本构模型选择等。

(5)模型计算。主要包括开挖和充填的设置。

(6)结果输出。

众所周知,前处理过程中材料参数、边界条件、力学条件、网格模型等输入条件决定着数值模拟计算的准确度。受到计算速度和存储能力的限制,不得不对模型进行简化和网格粗化处理,导致的直接后果是计算精度的牺牲。

然而在计算过程中,矿山开采模拟与实际矿山生产开采匹配度同样对计算精度产生较大影响。本文针对模型计算这一环节中的开采仿真度进行研究,在保证同等前处理条件下设置中段开采、盘区开采和进路开采三种不同维度的开采模拟方式。

3 多维度开采模拟研究

3.1 工程背景

金川镍矿是世界著名的多金属共生的大型硫化铜镍矿床之一,发现于1958年,集中分布在龙首山下长6.5 km、宽500 m的范围内。金川二矿区是承担着金川矿区主力生产任务的矿区,通过机械化生产使得年产量达400万t以上。在金川矿体中的1#矿体是二矿区主要进行开采的矿体,金川二矿区是较早利用下向充填采矿法实现大规模开采的矿山,并通过各种支护手段进行矿山安全的维护管理。

金川二矿区按50 m分成一个回风中段,2015年已掘中段为1 350 m、1 300 m、1 250 m、1 200 m、1 150 m、1 000 m、850 m。中段之间的分段高度为20 m,分为5个分层,每层高4 m。盘区方向垂直于矿体,宽100 m,盘区间全部开采。盘区上、下分层宽4 m、高4 m,分层进路宽5 m、高4 m,进路长不超过50 m。盘区开采顺序为先上盘后下盘,先两翼后中间,后退式开采。开采矿体范围位于6行至25行之间,采用机械化下向分层水平进路胶结充填采矿法。根据矿体走向和形态变化,灵活划分工区和盘区,全区目前共有6个工区,其中一、二、三工区位于1 000 m中段,且各工区均包含有3个盘区;四、五工区位于850 m中段,且四工区包含4个盘区;五工区包含3个盘区;六工区位于1 150 m中段,其主要包含Ⅳ、Ⅴ、Ⅶ这3个盘区。二矿区整体开采纵投影如图1所示。

图1 金川二矿区盘区分布图

3.2 模型概况

将金川二矿区1#矿体6～24行采区的地质模型导入FLAC3D,并划分围岩、超基性岩、贫矿、富矿分组。模型地表简化为平面考虑,模型尺寸为垂直矿体走向方向x轴 2 880 m,沿矿体走向y轴4 160 m,垂直方向z轴1 728 m,其中矿体最小单元尺寸细化至4 m,单元逐渐向外放大,围岩单元尺寸32 m,共生成约334余万个网格单元。导入开拓系统模型标定地下主要工程与采区位置关系,模型情况如图2所示。

图2 金川二矿区数值模型

综合现有研究报告,富矿、贫矿、围岩、充填体的岩体物理力学参数取值见表1。

表1 模型矿岩体物理力学参数表

根据金川矿区最近研究成果,二矿区最大主应力为水平构造应力,最大主应力方向为北北东向,与矿体走向垂直。根据地应力测试结果,在模型内部施加阶梯地应力条件。本构模型选择摩尔库伦模型。

3.3 中段回采模拟

按50 m分一个中段,各中段开采完成后进行充填,主要分单中段、双中段、三中段三种采充形式,共分6个采充阶段。中段回采模拟示意如图3所示,回采顺序见表2。

图3 中段回采示意图

表2 中段回采模拟顺序表

3.4 盘区回采模拟

沿矿体走向按100 m划分一个盘区,盘区开采模拟示意图如图4所示。盘区回采顺序参照中段回采表2,自上而下,盘区隔一采一。

图4 盘区回采模拟示意图

3.5 进路回采模拟

沿矿体走向每5 m划分一条进路,进路回采模拟示意如图5所示,回采顺序参照中段回采表2,自上而下每4 m一采,进路隔一采一。

图5 进路开采模拟示意图

4 结果分析

为便于展示开采充填过程对围岩及重要水平、中段的影响,现取14行剖面及1 000 m水平剖面做结果分析。

4.1 最大主应力

当开采至812 m水平,此时1 000 m水平矿柱回收完毕,850 m水平矿柱高56 m。图6至图8所示分别为14行剖面中段、盘区和进路三个不同层级回采模拟的最大主应力分布。中段回采模拟时,850 m水平矿柱表现为卸压状态,应力集中于812 m以下,呈月牙形分布,最高达94.1 MPa;盘区回采模拟时,850 m水平矿柱卸压状态不显著,应力水平与周边围岩相近,应力同样集中于812 m以下,呈月牙形分布,但应力集中程度小于中段回采模拟,最高达93.8 MPa;进路回采模拟时,850 m水平矿柱表现为应力集中,另有812 m水平近上盘矿岩交界处局部小范围应力集中,最大主应力高达94.1 MPa。

图6 中段开采模拟最大主应力

图8 进路开采模拟最大主应力

综合对比三种仿真层级的开采模拟,发现最大主应力数值相差无几,但是应力集中程度和分布区域相差较大。对于上下盘围岩应力的影响也不同,中段和盘区开采模拟时,上盘围岩卸压范围较大,进路开采模拟时,上盘围岩卸压范围较小。

4.2 塑性区

当开采至812 m水平,此时1 000 m水平矿柱回收完毕,850 m水平矿柱高56 m。图9至图11所示分别为中段、盘区和进路三个不同层级回采模拟的塑性区分布。中段回采模拟时,富矿体的开采对周边贫矿和围岩影响较大,塑性区范围最大;盘区回采模拟时,富矿体的开采对周边贫困和围岩影响略小,塑性区范围较中段回采模拟时小;进路回采模拟时,富矿体的开采对周边贫矿及围岩几乎无影响,塑性区范围最小。

图9 中段开采模拟塑性区分布

图10 盘区开采模拟塑性区分布

图11 进路开采模拟塑性区分布

4.3 竖向位移

当开采至756 m水平,此时850 m水平矿柱回收完毕,图12至图14所示分别为14行剖面中段、盘区和进路三个不同层级回采模拟的竖向位移分布。中段回采模拟时,最大沉降是1.03 m,最大隆起是0.58 m;盘区回采模拟时,最大沉降是0.85 m,最大隆起是0.52 m;进路回采模拟时,最大沉降是0.32 m,最大隆起是0.19 m。

图12 中段回采竖向位移

图13 盘区回采竖向位移

图14 进路回采竖向位移

三种回采模拟方式均表现出上盘围岩沉降,下盘围岩隆起的基本趋势,但是中段回采模拟时的竖向位移值最大,进路回采模拟时的竖向位移值最小。取标志性测点14行风井1 000 m水平处竖向位移进行对比。中段回采时,14行风井底部竖向隆起约0.18 m;盘区回采时,14行风井底部竖向隆起约0.1 m;进路回采时,14行风井底部竖向隆起约0.08 m。14行剖面竖向位移对比见表3。

表3 14行剖面竖向位移对比

5 结论

本文结合实际矿山设计过程中遇到的岩石力学模拟问题,通过FLAC3D按照中段、盘区和进路开采三种仿真维度模拟了大规模进路充填开采金川二矿区的开采过程。应力集中、塑性区及竖向位移结果表明,不同的开采仿真维度对应不同的力学响应。开采模拟越细,开采扰动越小,数值模拟与实际开采的匹配程度越高,结果表现为应力集成程度、塑性区范围、竖向位移越小。因此除前处理阶段网格划分、参数设置外,在计算阶段仍需注重提高开采过程的模拟精细程度。