关于疆锋铁矿矿体开采离散元数值模拟研究

杨国永

（景洪市龙鑫矿业有限责任公司，云南 景洪 666100）

自1970年Cundall第一次提出离散单元DEM（DistincElementMethod）模型后，这一技术已经在数值模拟、工程应用领域获得了很大的进步［1-2］。该技术在体现出岩土之间出现的位移的同时还反映他们之间的应力应变情况。所以，所有的岩体的本构模型都能应用到模型中来。此方法的第二个特点是他可以通过动态松弛法来解答动力平衡方程，这种方法可以很好的解答动力稳定和非线性大位移问题。所以自离散元方法出现至今，已经获得了非常广泛的应用和深入的研究。我国最早是由王泳嘉教授在节理岩体的数值分析中进行了放矿数值模拟研究［3］。

1 分析方法

离散元理论根据牛顿第二定律，建立起了力与速度以及位移之间的关系，针对不连续和分散的单元模拟仿真［4］。块体之间的相互作用主要为角与角、角与边及边与边的接触，如图1所示。

图1 块体单元接触关系图

块体之间的作用力取决于它们的叠合值U，我们将块体之间的法向、切向叠合值作为Un、Us，则对应的法向和切向的力Fn、△Fs分别为式（1）和式（2）：

式中：法向和切向刚度系数分别为Kn、Ks［5］。

式（1）和式（2）适用于弹性状态下。而对于塑性破坏，岩块在法向和切向的力会因为岩块的脱离而消失，因此，切向力Fs的大小是我们判断是否发生剪切破坏的依据。所以，根据库仑定律有式（3）成立：

式中：C为单元之间的粘聚力；φ为单元之间的内摩擦角。

根据岩块之间的相互关系以及形状，根据上述基本原理，便能得到作用于一个给定岩块上的Fx1、Fx2、Fy1、Fy2等力，并可计算出合力Fx、Fy和力矩M。

法向和切向的刚度系数Kn和Ks可以通过下式进行计算。

两个块体如图2所示，设其长为ɑ，宽为b，弹性模量为E，泊松比为μ。所以通过弹性力学理论可以得到式（4）～式（6）：

图2 块体计算接触模型

从而可得到法向刚度系数：

切向刚度系数可由法向刚度系数求得：

为了消除不稳定解的情况，一般通过阻尼的方法来得到稳定的解［6］。因此，可得沿X、Y方向以及转动（角）的速度为式（7）～式（9）：

式中：m为岩块质量；M为转矩；I为转动惯量；t为时步；ɑ为阻尼。

计算时，不断对（△t）迭代，通过上一步得出的结果作为下一步的起点，不断重复直到最终平衡。

2 计算模型

本次研究所采用的计算模型，综合考虑了疆锋铁矿Ⅱ号矿段矿体与围岩、断层之间的实际关系以及所要解决的技术问题，选取了具有代表性的24剖面建立了所需要的离散元模型，力求找出主矿体不同开采深度下对周边岩层的影响范围。

本研究所使用的计算模型，根据矿体与围岩、断层的实际相互作用和需要处理的问题，建立计算模型，模型如图3所示。矿体以及距矿体较近的围岩单元大小为3m×3m的结构单元，较远离矿体的围岩单元大小为6m×6m的结构单元，（因电脑计算能力的影响，单元体比较现场数据进行放大），结构单元倾角与现场实际调查数据的倾斜角相同。计算模型单元划分为47586个。

图3 24剖面主矿体与围岩、断层模型单元划分示意图

本模型将底部和左右两边采取位移边界的限制，不对形成空区后的顶部和内部施加约束，允许围岩在自然状况下产生的移动和破坏，从而研究矿体开采后采空区的冒落情况及影响范围。

3 计算方案

开挖方案采用分步开挖的方式，分四步先建立初始模型，不进行开挖，后依次将矿体开挖至460、340、220m。模拟不同开挖步骤上盘围岩的变形、冒落形式及范围。

建模时先固定模型的左右两边以及底部，这是为了可以研究单一重力因素下的岩石冒落的运动规律。离散元提供了剖面的位移图、速度图、主应力分布图等大量接触关系的信息。但鉴于各剖面的规律基本一致，本研究仅对24剖面进行了模拟。

4 模拟结果分析

4.1 原始状态

为了研究围岩在自然状况下的状态，不去约束形成采空区后的模型内部，通过岩体的弹性模量和泊松比作为法向刚度系数和切向刚度系数，用单元的几何尺寸进行换算得出。模拟需选用合适的模拟步长，经过大量模拟，结果分析如下：

模拟起始时，不开展开挖作业，将边界约束条件以及外部荷载添加到研究的区域当中，建立出初始应力场。模型由上到下，位移依次减小，最大位移值出现在顶端。这是因为上部施加的外部荷载和块体自身重力作用下的压实作用所造成的。这一模拟过程中，建立了矿石和岩体之中的初始应力场，其中，弹性形变是造成块体之间发生位移的主要原因。

4.2 矿体开采至460m

模拟第二步是将矿体在压实的前提下开采至460m，让上盘围岩依照自身的力学机制发生冒落。从模拟结果可以看出，矿体开采至460m后，上盘围岩发生了大面积冒落，几乎充填了整个采空区，采空区上盘放大图如图4所示，由图4可知，下盘围岩冒落范围以超出F4断层，而下盘围岩主要沿F6断层滑移。围岩移动范围图如图5所示，通过图5画出的放大图里岩石移动范围发现，24剖面块体开采至460m时下盘岩石移动角约为64°，上盘为78°。

图4 24剖面主矿体开采至460m上部放大图

图5 24剖矿体开采至460m围岩移动范围

4.3 矿体开采至340m

模拟第三步是在矿体开采至460m，上盘冒落围岩区域稳定的基础上将矿体开采至340m，让上盘围岩依照崩落法的特点，依照自身的力学机制发生冒落。从模拟结果图可以看出，主矿体开采至370m后，上盘围岩在原有基础上继续向下移动，从采空区上盘放大图如图6所示，由图6可知，由于保安矿柱的存在，上部围岩冒落至地表后继续移动的趋势较小。围岩移动范围如图7所示，通过图7画出的放大图里的岩石移动范围发现，当24剖面块体开采至340m时下盘岩石移动角约为62°，上盘为77°。

图6 24剖面矿体开采至340m上部放大图

图7 24剖矿体开采至340m围岩移动范围

5 结论

1.采用离散元对金属矿山地下采空区进行数值模拟，能够科学的评价对采矿区域顶板、附近围岩、矿柱等力学行为，能直观的体现所研究的结论，同时也适应了岩石力学的客观规律，是一个有效的方法。

2.数值模拟过程中，允许围岩依照自身的力学机制发生冒落，从模拟结果看，24剖面块体开采至460m时下盘岩石移动角约为64°，上盘为78°；开采至340m时下盘岩石移动角约为62°，上盘为77°。空场条件下，Ⅱ号矿段的开采会引发上部围岩大范围的冒落，下盘以延伸至F4断层，上盘延伸至F10断层。