朱忠华^1,2, 王李管¹, 涂小腾¹, 毕林¹
1. 中南大学 资源与安全工程学院， 湖南 长沙 410083
2. 南华大学 核资源工程学院， 湖南 衡阳 421001
收稿日期: 2015-08-12
基金项目: 国家高技术研究发展计划项目(2011AA060407).
作者简介: 朱忠华(1981-)，男，湖南郴州人，中南大学博士研究生;
王李管(1964－)，男，山西临汾人，中南大学教授。

摘要: 为准确掌握自然崩落法崩落矿岩流动特性，提出属性块体建模与随机介质理论相结合的方法，进行崩落矿岩流动模拟.建立属性块体模型拟合目标区域，实现目标区域离散化；建立适宜程序实现的颗粒体数据结构和随机介质空位传递模型；提出放出指数和块度指数概念.分别模拟和分析固定目标放矿区域、不同放矿高度以及不同块度条件下的矿岩流动，并在数字矿山软件平台上用C++编程语言予以实现.结果表明：该方法得到了崩落矿岩流动过程中品位变化及形态发育规律.
关键词：矿岩流动模型随机介质颗粒流自然崩落法
Flow Characteristics of Discrete Ore Rock of Block Caving Based on Stochastic Medium Theory
ZHU Zhong-hua^1,2, WANG Li-guan¹, TU Xiao-teng¹, BI Lin¹
1. 1.School of Resources and Safety Engineering, Central South University, Changsha 410083, China
2. 2. Nuclear Resources Engineering College, University of South China, Hengyang 421001, China.
Corresponding author: ZHU Zhong-hua， E-mail： zzhnihao545@126.com
Abstract: In order to accurately grasp the flow characteristics of ore and rock of block caving, a new method that combining attribute block model with stochastic medium theory was proposed. The procedure of simulation flow is as follows: attribute block model was established to fit and discretize target region. Block attributes were obtained by attribute interpolation technique. A proper flow model for programming was established. The concept of draw-out index and block index was put forward to analyze the flow characteristics of the rock under different oredraw conditions. Flow simulation was realized based on the proposed method in 3D visualization platform with C++ language. The results show that with the draw height increasing, the length of the long axis of the draw out body increases linearly, and the length of the short axis and the eccentricity increase gradually, showing a power exponential relationship. Under constant block index, with the increase of draw-out index, the ellipsoid of the ore body is gradually developed, the eccentricity ratio is increased, and the average depth of the depression increases. Under constant draw-out index, with the increase of block index, the eccentricity of the ellipsoid is decreased.
Key Words: attribute block modelore and rock flow modelstochastic mediumparticle flowblock caving
自然崩落法是一种高能力、低成本和高技术含量的采矿方法，应用较好可以大幅降低采矿成本，产生很好的经济效益^[1]；崩落矿石含有品位和岩性等信息，并且和覆岩直接接触，在覆岩下出矿，若采场结构参数设计不合理，或者放矿管理不当，将引起较大的贫化及损失，造成矿产资源浪费和矿山企业经济效益下降.放矿问题解决的好坏直接与崩落采矿法应用成败相关联.准确掌握崩落矿岩的流动特性，具有重要的理论和现实意义.
近年来，三维建模技术和数字矿山技术的快速发展使得地质建模和三维可视化变得方便^[2-5].为构建离散模型，黄刚海等基于MFC进行了球颗粒模型可视化建模^[6].常见的颗粒流软件PFC^2d及PFC^3d基于牛顿第二力学定律模拟颗粒体流动^[7-8]，但是不能表征颗粒体的属性信息.块段模型可以表征属性，相关学者对其有过研究和应用，取得了较好成果^[5].在三维地质体模型和块段模型的基础上可以进行统计分析，如岩性分析、品位估值等，但是这种统计分析是静态的，不能反映地质体及其属性的时空变化.
为掌握属性崩落矿岩的流动特性，本文将块体建模技术和属性模型相结合，建立属性块体模型，对放矿目标区域进行离散，然后将随机介质理论应用于属性块体模型，建立适宜编程实现的随机介质理论流动模型及其数据结构，模拟崩落矿岩块体及颗粒体流动，并考虑不同块度条件下的块体流动特性；在数字矿山软件平台上，使用C++语言编程实现.
1 属性块体建模技术三维块体建模的实质是对目标区域进行三维离散化.目标区域可以由起始点和三个方向的尺寸界定，也可以由三维表面模型界定.离散单元可以是四面体或六面体单元等，常用的八叉树构建块段模型使用的是六面体单元.
本文属性块体单元为属性六面体.属性六面体在几何六面体的基础上增加一个三维点，即六面体的几何中心点，作为块体属性信息的承载对象，记录块体属性.离散后的块体具有本文定义的属性结构，其属性值可赋值或由已知取样值使用属性插值得到，并存储在对应的结构体struc_Particle对象中.
2 基于随机介质理论的崩落矿岩块体流动模拟2.1 随机介质理论空位传递模型将随机介质理论用于崩落矿岩流动的计算机模拟，关键是确定离散化目标区域的随机流动概型.在放矿学或岩土工程中，不同的端部条件对应的概率模型不同，常见有六角模型、七块模型和九块模型，如常用的D. Jolley模型是九块模型^[9].
流动概率模型由目标区域离散方法和流动概型确定.离散后的区域块体流动概率模型体现为空位传递模型.本文选用九块模型、六块模型和四块模型分别考虑不同边界条件下的空位传递，如图 1所示.四块模型对应目标区域“角落”边界条件；六块模型适用端部受限边界条件；而九块模型适用无边界限制的情况；每次移动对应一随机过程，其概率由Monte Carlo法确定.综合考虑这三种模型，可以完整模拟离散目标区域的流动.
图 1(Fig. 1)

图 1 空位传递模型Fig.1 Transfer model of null block (a)—九块模型; (b)—六块模型; (c)—四块模型.

2.2 随机介质理论和属性块体相结合将随机介质理论和属性块体相结合进行块体流动模拟，流动对象用块体中心点代替，将块体的流动简化为颗粒的流动，可以减少存储空间，又不产生主要信息的失真.
2.2.1 流动颗粒结构体设计设计流动颗粒结构体，进行属性块体的流动模拟，其数据结构设计如下：
struct struc_particle??//流动颗粒体结构体
{
??dmDPoint origin_location; //初始位置
??int ijkID^[3]; //初始定位数组
??int curijkID^[3]; //当前定位数组
??double grade; //品位等属性
??bool drawn;//状态，是否放出
??dmDPoint cur_location; //当前位置
??int oreFlag;//0--岩石，1--矿石
??int index; //索引号
??bool used; //状态，是否已使用
????……????//其他可扩展的属性
};
2.2.2 颗粒体流动模型及模拟流程流动过程中，流动对象是颗粒，而属性块体作为离散单元，作用是对颗粒进行定位及属性赋值.颗粒流动模型如图 2所示.
图 2(Fig. 2)

图 2 颗粒流动模型Fig.2 Flow model of particle

离散颗粒点具有本文2.2.1节定义的流动颗粒属性结构，存放于颗粒结构体数组m_vecParticleStruc之中；由Monte Carlo法确定出口定位数组，以此开始属性颗粒体的流动模拟，模拟流程如图 3所示.通过模拟得到放出体、残留体、接触面，并对其形态和属性进行分析.
图 3(Fig. 3)

图 3 属性块体流动模拟流程图Fig.3 Flow chart of attribute block flow simulation

颗粒流动模拟实现倚赖以下变量和函数：表示放出颗粒数量的整形变量nDrawOut，由放出指数κ确定，放出指数κ为放出颗粒数量与总颗粒数量的比值；表示空位定位数组的数组curNullIJK，preNullIJK；表示临时位置的dmDPoint类型变量TmpLocation,TmpLocation2,preLocation；pointInFun()函数，用以判断流动过程中在漏口范围内的颗粒，将满足条件的颗粒对应的定位数组存放于pointInFun数组中，对元素为n的pointInFun数组，第i个颗粒的定位数组为pointInFun[i].IJK ,0<i<n； FlowMode()函数，用以确定流动模型，由随机介质理论流动模型决定；而函数FromIJKtoVecScript()功能则是由定位数组得到下标值.
3 不同块度条件下块体流动特性矿岩块度对放矿的矿岩运动性能、放出强度和矿石损失贫化都有很大影响，决定着放出体形态和放矿漏斗尺寸，影响放矿结构参数设计、出矿设备选型及生产能力.
在放矿过程中，矿岩块度大小、形状及块度组成对放矿都有重要影响.实际生产过程中，随着矿岩流动和相互摩擦和碰撞，矿岩块度呈现动态的变化.由2.2节可知，基于随机介质放矿流动模拟，首先就是将矿岩堆进行离散化，然后再用流动模型进行放矿计算机模拟.如考虑矿岩块体形状及块度组成，将极大增加放矿模拟的困难，为此，本文仅考虑均匀矿岩块度和规则形状，以简化模型，增加计算机放矿模拟的可行性，同时又不失一般性和代表性.
考虑矿岩块度的矿岩颗粒流动模型与颗粒流动模型一致，见图 3.图中块度通过矿岩颗粒和放矿漏斗尺寸的比值来体现，称为块度指数，在矿岩流动模拟中，矿岩块度通过与漏斗尺寸相结合来考虑，由块度指数来体现.块度指数β定义为

(1)

式中:d_bell表示漏斗的最小尺寸;d_bloc表示崩落矿岩体的最大尺寸或平均尺寸.矿岩流动模拟时，由Monte Carlo法确定出口定位数组，不同的β决定了矿岩从某一确定漏斗放出的位置，也就确定了出口定位数组，对模拟程序具有重要的作用.
4 系统开发以某数字矿山软件系统为平台，采用平台+插件的方式进行二次开发.该平台开发框架包括数据层、操作层和算法层以及显示层，具有较好的实用性和可扩展性，开发模式如图 4所示.用C++面向对象编程语言实现属性块体流动模拟插件子系统的开发，数据层使用自定义的数据结构；操作层使用平台提供的交互机制，完成插件的定义、对话框及数据的交互；算法层则是将属性块体流动模拟功能封装为BlockFlowOperator()类，完成属性块体的流动模型；显示层调用平台的显示接口，系统界面如图 5所示.
图 4(Fig. 4)

图 4 平台+插件的开发模式Fig.4 Platform and plug-in development mode

图 5(Fig. 5)

图 5 系统界面Fig.5 System interface

5 算 例5.1 固定目标放矿区域的矿岩流动目标区域尺寸为20 m×20 m×30 m的矿岩堆，离散属性块体模型尺寸为单位尺寸，记录属性为Cu元素品位，漏口尺寸3.5 m×3.5 m×3.2 m.建立块体模型后，通过属性插值技术得到块体的品位，存储在属性块体中心点对应的颗粒结构体上.目标区域离散颗粒模型见图 6.
图 6(Fig. 6)

图 6 目标区域离散颗粒Fig.6 Discrete particle in target region

采用本子系统模拟属性颗粒的流动，模拟分为6组，分别对应放出指数为5%,10%,15%,20%,25%,30%的颗粒流动，选取的放出指数不大，是因为过大的放出指数将导致部分颗粒因脊部残留效应放不出来；此外对于固定目标区域的流动模拟，过大的放出指数还将导致明显的截头椭球体的产生，从而不利于描述放出体的形态.模拟结果见图 7和表 1.
图 7(Fig. 7)

图 7 放矿模拟结果(κ=10%)Fig.7 Results of oredraw simulation(κ=10%) (a)—放出体; (b)—放出体和矿岩接触面复合图; (c)—残留体; (d)—残留体和矿岩接触面复合图.

表 1(Table 1)

表 1 不同放出指数下的颗粒流动模拟结果Table 1 Result of particle flow simulation under fluctuant draw-out index 放出指数κ/% 5 10 15 20 25 30 放出体Cu品位/% 0.49 0.49 0.53 0.51 0.50 0.51 放出体长半轴/m 9.16 11.56 14.5 16.2 17.0 17.5 放出体短半轴/m 5.55 6.5 8.0 8.1 8.3 8.5 偏心率 0.79 0.83 0.84 0.87 0.87 0.88 残留体体积/m3 11 407 10 820 10 212 9 636 9 048 8 456 凹陷坑平均深度/m 1.48 2.95 4.47 5.91 7.38 8.86

表 1 不同放出指数下的颗粒流动模拟结果 Table 1 Result of particle flow simulation under fluctuant draw-out index

放出体指从漏口放出的颗粒体流动之前在原始目标区域中所形成的空间体，从图 7可以看出，放出体形态近似为椭球体.从表 1可知该方法得到了放出体的Cu品位信息，随着放出指数的增加，椭球体离心率逐渐增加且趋近1，残留体体积逐渐减小，凹陷坑深度逐渐增加.
算例模拟结果图表说明本文的算法及开发的系统是可行和实用的.
5.2 不同放矿高度下的矿岩流动在2.2节的基础上，采用本子系统分别模拟放矿口尺寸为3.5 m×4.0 m条件下放矿高度H=20，40,60,80,100,150 m的放出体长半轴、短半轴，分析放出体形态发育规律，为放矿结构参数设计提供参考和依据.模拟结果放出体形态发育如图 8所示.
图 8(Fig. 8)

图 8 放出体形态发育Fig.8 Results of oredraw simulation (a)—H=20 m; (b)—H=40 m; (c)—H=60 m; (d)—H=80 m; (e)—H=100 m; (f)—H=150 m.

根据图 8中不同放矿高度时放出体的具体形态，测量对应的放出体长半轴和短半轴长度，分析放出体形态的发育规律，见表 2.
表 2(Table 2)

表 2 放出体参数表Table 2 Parameters of draw-out body 放矿高度/m 20 40 60 80 100 150 200 300 长半轴长度/m 10.56 21.34 31.75 42.05 52.78 74.68 106.6 157.5 短半轴长度/m 6.3 8.7 11.7 12.6 14.1 16.4 19.5 20.4 偏心率 0.80 0.91 0.93 0.95 0.96 0.98 0.98 0.99 Cu品位/% 0.51 0.54 0.55 0.52 0.51 0.49 0.50 0.47

表 2 放出体参数表 Table 2 Parameters of draw-out body

在表 2的基础上，根据曲线拟合的最小二乘法，得出放出体长半轴长度、短半轴长度及放矿体偏心率与放矿高度关系:

(2)

(3)

(4)

式中：A为放出体长半轴长度，m；H为放矿高度，m；B为放出体短半轴长度，m；ε为放出体偏心率.
放出体长半轴长度与放矿高度近似成线性关系，随着放矿高度增加，放出体长半轴长度呈线性增加；而随着放矿高度的增加，放出体短半轴长度与放出体偏心率逐渐增加，但是增加的幅度逐渐平缓，基本呈现幂指数关系.利用上述得出的放出体参数与放矿高度的关系式，可在后续软件计算中定义放出体参数，也可以为放矿结构参数选取提供依据和参考.
5.3 不同块度条件下的目标放矿区域矿岩流动采用和5.1节相同的放矿区域，使用本文第3节方法进行不同块度条件下的矿岩流动模拟，模拟结果见图 9及表 3.从表 3可得出，在放出指数不变情况下，放矿椭球体随块度指数的增加，放矿椭球体偏心率呈递减趋势.
图 9(Fig. 9)

图 9 不同块度下的放出体形态(κ=10%)Fig.9 Results of oredraw simulation(κ=10%) (a)—块度指数β=1； (b)—块度指数β=2； (c)—块度指数β=3.

表 3(Table 3)

表 3 不同块度指数下的颗粒流动模拟结果(κ=10%)Table 3 Result of particle flow simulation under fluctuant block index(κ=10%) 块度指数 放出体Cu品位/% 放出体短半轴/m 放出体偏心率 凹陷坑平均深度/m 1 0.50 6.3 0.84 2.93 2 0.49 6.5 0.83 2.95 3 0.51 7.5 0.80 2.94

表 3 不同块度指数下的颗粒流动模拟结果(κ=10%) Table 3 Result of particle flow simulation under fluctuant block index(κ=10%)

6 结 语得到了流动放出体的品位信息，放矿过程中随着放矿高度增加，放出体长半轴长度呈线性增加，放出体短半轴长度与放出体偏心率逐渐增加，但是增加的幅度逐渐平缓，基本呈现幂指数关系；在块度指数不变的条件下，随放出指数的增加，放矿椭球体逐渐发育，偏心率逐渐增大，凹陷平均深度增加；在放出指数不变的条件下，随着块度指数的增加，放矿椭球体偏心率呈递减趋势.
参考文献

[1]	冯兴隆.自然崩落法矿岩工程质量数字化评价及模拟技术研究[D].长沙：中南大学,2010. (Feng Xing-long.Study on digital assessment and simulation technology at rock mass engineering quality of block caving[D].Changsha:Central South University,2010.)(0)
[2]	Itasca Consulting Group,Inc.Particle flow code 3D, version 4.0[R].Minneapolis:Itasca Consulting Group,Inc, 2007.(0)
[3]	Pourrahimian Y, Askari-Nasab H, Tannant D D. A multi-step approach for block-cave production scheduling optimization[J].International Journal of Mining Science and Technology, 2013, 23(3): 739–750.(0)
[4]	Shewchuk J R.Delaunay refinement mesh generation[D].Pittsburgh:Carnegie Mellon University,1997.(0)
[5]	毕林, 王李管, 陈建宏, 等. 基于八叉树的复杂地质体块段模型建模技术[J].中国矿业大学学报, 2008, 37(4): 532–537. ( Bi Lin, Wang Li-guan, Chen Jian-hong, et al. Study of octree-based block model of complex geological bodies[J].Journal of China University of Mining & Technology, 2008, 37(4): 532–537.)(0)
[6]	黄刚海, 焦玉勇, 张秀丽. 岩质工程地质体三维球颗粒计算模型可视化建模[J].岩土工程学报, 2014, 36(4): 782–786. ( Huang Gang-hai, Jiao Yu-yong, Zhang Xiu-li. Visual modeling for 3D spherical particle for geological body model of rock engineering[J].Chinese Journal of Geotechnical Engineering, 2014, 36(4): 782–786.)(0)
[7]	Han K, Feng Y T, Owen D R J. Sphere packing with a geometric based compression algorithm[J].Powder Technology, 2005, 155(1): 33–41.(0)
[8]	Liu Z, Zhou N, Zhang J X. Random gravel model and particle flow based numerical biaxial test of solid backfll materials[J].International Journal of Mining Science and Technology, 2013, 23(4): 463–467.(0)
[9]	任凤玉. 随机介质放矿理论及其应用[M].北京: 冶金工业出版社, 1994. ( Ren Feng-yu. The random medium theory of oredrawing and its application[M].Beijing: Metallurgical Industry Press, 1994.)(0)