,罗周全,秦亚光,孙杨中南大学 资源与安全工程学院, 湖南 长沙 410083
收稿日期: 2014-09-28
基金项目: "十二五"国家科技支撑计划项目(2012BAk09B02-05); 国家自然科学基金资助项目(51274250).
作者简介: 汪 伟(1990-), 男, 安徽安庆人, 中南大学博士研究生; 罗周全(1966-), 男, 湖南邵阳人, 中南大学教授, 博士生导师.
摘要: 为实现某铅锌矿山无底柱深孔后退式崩矿法安全高效开采，开展采场结构参数优化研究.实测-650 m中段原岩应力，运用Surpac-FLAC3D模型生成技术构建数值分析模型，进而采用多元应力回归方法反演初始地应力场.通过弹性力学小薄板理论分析得到顶板稳定性随采场结构参数的变化关系，结合实际采准条件，提出采场结构参数初选方案.基于反演地应力场，进行初选方案开采数值分析，获得各方案力学响应指标.引入多目标理想点法决策，考虑安全和经济指标，计算方案评价指标集与理想解的向量近似度，优化确定该矿无底柱深孔后退式崩矿法采场结构参数为矿房采场宽10 m，矿柱采场宽8 m.
关键词: 初始地应力场采场结构参数弹性力学数值分析理想点法决策
Stope Parameters Optimization of Non-pillar Longhole Retreat Caving
WANG Wei
,LUO Zhou-quan,QIN Ya-guang,SUN YangSchool of Resources and Safety Engineering, Central South University, Changsha 410083.
Corresponding author: WANG Wei, professor, E-mail: csu424@163.com
Abstract: In order to realize safe and efficient mining of non-pillar longhole retreat caving in a lead-zinc mine, stope parameters optimization was conducted. Field measurement of in-situ stress was carried out on-650 m level, and a numerical model was constructed by Surpac-FLAC3D modeling technique. On these grounds, in-situ stress field was inversed by multiple regression method. The change relation of roof stability with stope parameters was obtained by using the small thin elastic mechanics theory, and the primary stope parameters schemes were put forward by combining with actual mining condition. Based on the inversion of ground stress field, the mining numerical analysis of primary schemes was conducted to get mechanical response indexes of different schemes. The safety and economy indexes were comprehensively considered, and the multi-objective TOPSIS decision method was introduced, the superior order of parameter schemes was determined by calculating the vector similarity of evaluation indexes set to ideal solution. Optimal stope parameters of non-pillar longhole retreat caving were finally determined, i.e., chamber width 10 m, pillar width 8 m.
Key words: in-situ stress fieldstope parameterelastic mechanicsnumerical analysisTOPSIS decision 采场结构参数优化可使采场内部应力分布趋于均衡，避免局部应力集中，使采场处于安全有利状态.同时，采场结构参数的改变直接影响生产能力、采掘成本等经济指标，是决定开采效益的重要因素.因此，采场结构参数优化是实现矿山安全高效开采的重要手段.
诸多采场结构参数优化方法中，以数值方法应用最为普遍[1]，它既克服了传统经验类比法主观性大的缺点，又改善了纯数理计算效率低的不足，并可与试验手段或理论方法有效结合，适用性很强.Zhao等[2]通过Ansys有限元试算构建采场结构参数与稳定性指标对应的目标函数，利用混沌搜索获得最优参数方案，成功用于试验采场生产.罗周全等[3]采用MIDAS/GTS模拟开采进程，通过比较不同方案力学响应指标，从安全角度确定了某铜矿盘区隔离矿柱最优开采参数.吴爱祥等[4]基于大间距结构参数理论，开展相似放矿试验和放矿数值模拟，据此采用回归分析方法优化确定某铁矿采场结构参数.然而，采场结构参数优化是涉及多层次、多因素和多目标的决策过程[5]，当前的优化方法仍存在一些不足：未准确考虑初始地应力场影响，而这恰是合理数值分析的关键；确定采场结构参数优化区间或初选方案时，缺乏合理依据，主观性大；多是从安全角度优化决策，缺乏对其他因素，特别是经济指标的综合考虑.
为此，本文将针对现存采场结构参数优化方法的不足，以某铅锌矿山无底柱深孔后退式崩矿法采场为对象，开展初始地应力场反演、弹性力学分析和数值模拟等工作，考虑安全和经济指标，引入多目标理想点法决策，综合确定采场结构参数，以期为矿山安全高效开采提供合理依据.
1 初始地应力场反演1.1 原岩应力实测-650 m中段为该矿山主要开采中段，段内多为倾斜中厚矿体，为准确获取中段原岩应力分布特征，采用改进型LUT套孔应力测定系统[6]进行现场实测.经与矿山地质及采矿工程师讨论，决定于N3号穿脉、S4号穿脉与S5号穿脉之间各布置1个地应力测点.套取带有LUT探头的完整岩芯，采集岩芯解除应变，并进行双轴率测定试验.将岩芯解除应变和弹性参数输入LUT-str应力计算程序，得到原岩应力实测结果，见表 1.
表 1(Table 1)
 表 1 主应力实测结果 Table 1 Measured results of principal stress
| 表 1 主应力实测结果 Table 1 Measured results of principal stress | ||||||||||||||||||||||||||||||||||||||||||
1.2 数值建模采用Surpac构建涵盖地层、断层、矿体等信息的三维实体模型，通过Surpac块体约束功能将实体模型剖分为块体网格模型，经Surpac-FLAC3D格式转换技术导出研究区域数值分析模型(图 1)，考虑边界效应，建模范围定为(x，y，z)/m：(2 440~2 890，7 830~8 730，-750~-450)，模型单元数共562 500个，节点数共585 738个.
图 1(Figure 1)
 | 图 1 数值分析模型 Fig. 1 Numerical analysis model |
计算区域包含①地层：泥盘系中统东岗岭阶上亚组(D2db)，泥盘系中统东岗岭阶下亚组(D2da)，泥盘系上统天子岭下亚组(D3ta).②断层.③铅锌矿和黄铁矿.岩体参数来自课题组研究结果[7].
1.3 初始地应力场反演采用多元回归方法反演初始地应力场[8]，研究中考虑自重应力和构造应力，共6种载荷工况.其中，构造应力工况为①x方向和y方向边界水平挤压构造运动；②xy平面内剪切变形构造运动；③xz平面和yz平面内竖向剪切变形构造运动.依次按照上述6种载荷工况进行初始地应力场试算，获得6组独立初始地应力场，则应力回归值与应力实测值最小二乘逼近的函数表达式为
式中：k为测点序号；j为应力分量编号；σjk*为测点应力分量实测值；Li为多元回归系数；σjki为不同工况下测点应力分量试算值；n为工况数.
基于最小二乘原理，采用SPSS软件求解S残满足极小的方程组[9]，解得回归系数L1=2.333，L2=0.124，L3=1.685，L4=-4.014，L5=2.824，L6=-1.731，复相关系数r=0.854.研究区域内任一点P的应力回归值σjp可由该点在各独立工况下的应力试算值按下式迭加而得：
根据式(2)，计算应力回归值，即得到反演初始地应力场.
2 采场结构参数方案初选2.1 优化目标确定该矿山无底柱深孔后退式崩矿法的开采工艺：采场垂直于矿体走向布置，采场长度为矿体厚度，高度为40 m，宽度为8 m.上部凿岩硐室高为3 m，顶板采用锚杆和铁丝网联合支护，下部不留底柱，采用平底硐室出矿，出矿硐室高为3 m.采切工程完成后，由上部硐室向下钻小孔径深孔进行拉槽爆破，以爆破空场为自由面分次侧向崩矿.
采场垂直于矿体布置，其长度为矿体厚度，高度为中段高，采场长度与高度均由工程实际条件决定，可变的主要结构参数为采场宽度.目前，设计部门参照矿山上向水平分层充填采矿法的开拓采准参数，暂将无底柱深孔后退式崩矿法矿房及矿柱采场的宽度均定为8 m，考虑到矿柱回采时，其两帮均为强度较低的充填体，失稳概率较大，因此，需要研究解决矿房和矿柱采场宽度的优化取值问题.
2.2 力学分析以空场落矿为基本特征的采场，其稳定性状态多由顶板要素决定[10].将顶板简化为受上覆岩体自重Q作用的板状结构(图 2)，采场长、宽和顶板厚度依次为a，b和H，薄板四周为固定支撑.
图 2(Figure 2)
 | 图 2 采场顶板受力示意图 Fig. 2 Force diagram of stope roof |
基于弹性力学小变形薄板理论，顶板受自重载荷Q作用，运用瑞利-李兹法[11]解出其挠曲函数：
式中，i=j= 3,5,7….进一步推导可得
式中：v为泊松比；E为弹性模量；H为顶板厚.
由弹性力学理论解得顶板各向应力分布：
由式(5)求出顶板最大主应力为
式中k为安全系数.根据H.Tresca屈服准则[12]，顶板某点产生剪切屈服时的主应力满足
由式(6),式(7)解得顶板最大剪应力为
比较τmax与顶板抗剪强度τ的大小关系，即可判断顶板稳定性状态.
2.3 初选方案确定现以-650 m中段N6#S，N6#，N6#N采场为具体研究对象，N6#为设计矿柱采场.采场上覆围岩为D2db，覆岩容重γ=2 720 kg/m3，埋深h=730 m，Q=19.5 MPa，τ=21.6 MPa，a=48 m(上部硐室长)，H=1 m，v=0.23，k=2.根据式(8)计算获得采场顶板受力特征随宽度变化关系见图 3.
图 3(Figure 3)
 | 图 3 顶板最大剪应力与采场宽度关系 Fig. 3 Relationship between maximum shear stress and stope width |
由图 3可知：当采场宽度布置为12.5 m时，顶板最大剪应力为21.5 MPa，小于岩体抗剪强度(21.6 MPa)，当采场宽度大于12.5 m后，顶板最大剪应力将超过其抗剪强度.实际开采过程中，应力突变、工程扰动等不可预见因素对顶板稳定性影响很大，故采场宽度不宜大于12.5 m.结合当前采场布置宽度8 m及炮孔排距2 m等实际开采条件，拟取试验采场宽度为8，10和12 m，全部试验共9组方案.矿柱回采时两帮均为充填体，其宽度不宜大于矿房采场，故确定采场结构参数初选方案：①矿房宽8 m，矿柱宽8 m；②矿房宽10 m，矿柱宽8 m；③矿房宽12 m，矿柱宽8 m；④矿房宽10 m，矿柱宽10 m；⑤矿房宽12 m，矿柱宽10 m；⑥矿房宽12 m，矿柱宽12 m.
3 最优结构参数方案确定3.1 开采数值分析在确定初始地应力场和结构参数初选方案的基础上，进行开采数值分析，为免赘述，以方案②(矿房宽10 m，矿柱宽8 m)为例进行分析说明.
开采扰动影响为采场规模的2~3倍(图 4)，由于采场南北两帮及上部中段邻近区域均已充填，采场周边形成了较大规模的充填体贯通区域.充填体内形成了小范围的拉应力区，最大拉应力值为0.15 MPa.开采活动引起了较为明显的应力重分布，采场顶板及边界处均产生了应力集中，最大主应力值达到42.78 MPa.采场及邻近区域产生了小规模的拉伸和剪切破坏(图 5)，拉伸破坏分布于采场顶板，而剪切破坏是采场围岩破坏的主要形式，分布于采场顶底部及边界处.
图 4(Figure 4)
 | 图 4 最大主应力分布(Y=8 498 m) Fig. 4 Maximum principal stress distribution |
图 5(Figure 5)
 | 图 5 塑性区分布(Y=8 498 m) Fig. 5 Plastic zone distribution |
3.2 最优参数方案确定基于开采数值分析结果，获取各方案力学响应指标，结合矿山生产科提供的方案对应经济指标，汇总采场结构参数方案评价指标见表 2.
表 2(Table 2)
表 2 采场结构参数方案评价指标 Table 2 Evaluation indexes of stope parameter schemes
| 表 2 采场结构参数方案评价指标 Table 2 Evaluation indexes of stope parameter schemes | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
1)构建评价矩阵.评价矩阵综合考虑安全与经济因素，安全指标包括最大压应力(ai1)、最大拉应力(ai2)、塑性区体积(ai3)、顶板最大竖向位移(ai4)，经济指标为采切比(ai5)、采场生产能力(ai6).其中，采场生产能力(ai6)为效益型指标，取值越大越好；其他均为成本型指标，取值越小越好.为便于运算，将ai6的倒数定义为其新的指标取值，进行归一化处理，构建评价矩阵A：
式中Ai为初选方案的评价指标集，i=6.
2)确定正负理想解.以待选方案指标的最小值集合为正理想解(A+)，最大值集合为负理想解(A-)，得
3)向量相似度计算.计算指标向量集Ai与正负理想解A+，A-的距离为Di+，Di-：
Ai与评价指标集正负理想解的向量相似度为
对于成本型指标集，Ti值越大，其解越优[13].计算6组初选方案评价指标集与理想解的向量相似度Ti为0.789 3，0.828 9，0.788 9，0.822 5，0.732 0，0.190 8，故采场结构参数初选方案的优劣度排序依次为②，④，①，③，⑤，⑥.因而，确定方案②(矿房采场宽10 m，矿柱采场宽8 m)为综合考虑安全与经济因素的最优采场结构参数.
4 结论1) 采用弹性力学小薄板理论分析采场顶板受力特征随宽度的变化关系，并考虑应力突变、工程扰动等不可预见因素的影响，从力学角度确定采场安全开采的理论最大宽度为12.5 m.
2) 在原位实测进而反演初始地应力场的基础上，开展采场结构参数初选方案数值分析，综合考虑安全与经济指标，引入多目标理想点法决策，优化确定采场结构参数方案为：矿房采场宽10 m，矿柱采场宽8 m.
参考文献
| [1] | Henry B,Jung S J,Lim H W.Room and pillar stope design in highly fractured area[J].International Journal of Rock Mechanics and Mining Sciences,1997,34(3):674-678.( 1) |
| [2] | Zhao X D,Li L C,Tang C A,et al.Stability of boundary pillars in transition from open pit to underground mining[J].Journal of Central South University,2012,19(11):3256-3265.(1) |
| [3] | 罗周全,管佳林,冯富康,等.盘区隔离矿柱采场结构参数数值优化[J].采矿与安全工程学报,2012,29(2):261-264.(Luo Zhou-quan,Guan Jia-lin,Feng Fu-kang,et al.Stope structural parameters of panel isolation pillar numerical optimization[J].Journal of Mining & Safety Engineering,2012,29(2):261-264.)(1) |
| [4] | 吴爱祥,武力聪,刘晓辉,等.无底柱分段崩落法结构参数研究[J].中南大学学报(自然科学版),2012,43(5):1845-1850.(Wu Ai-xiang,Wu Li-cong,Liu Xiao-hui,et al.Study on structural parameters of sublevel caving[J].Journal of Central South University(Science and Technology),2012,43(5):1845-1850.)(1) |
| [5] | Luo Y Z,Wu A X,Liu X P,et al.Stability and reliability of pit slopes in surface mining combined with underground mining in Tonglushan mine[J].Journal of Central South University,2004,11 (4):434-439.(1) |
| [6] | Chen F,Xu G,Ma C D,et al.Thermal residual stress of polycrystalline diamond compacts[J].Transactions of Nonferrous Metals Society of China,2010,20(2):227-232.(1) |
| [7] | Luo Z Q,Xie C Y,Jia N,et al.Safe roof thickness and span of stope under complex filling body[J].Journal of Central South University,2013,20(12):3641-3647.(1) |
| [8] | Li K,Wang Y Y,Huang X C.DDM regression analysis of the in-situ stress field in a non-linear fault zone[J].International Journal of Minerals,Metallurgy and Materials,2012,19(7):567-573.(1) |
| [9] | Jiang Y,Wang Z L,Chu J X,et al.Study on the application of the fundamental statistic methods of SPSS in efficiency analysis of quercetin and resveratrol mixture[J].Advances in Information Sciences and Service Sciences,2012,17(4):380-387.(1) |
| [10] | Grasselli G.3D Behavior of bolted rock joints:experimental and numerical study[J].International Journal of Rock Mechanics & Mining Sciences,2005,42(1):13-24.(1) |
| [11] | Zheng T S,Hasebe N.Calculation of equilibrium position and dynamic coefficients of a journal bearing using free boundary theory[J].ASME Journal of Tribology,2000,122(3):616-621.(1) |
| [12] | Torkamani M A M.Elastic-plastic analysis for cyclic loading and Tresca yield condition[J].Computational Mechanics,1989,12(5):407-422.(1) |
| [13] | Lin M C,Wang C C,Chen M S,et al.Using AHP and TOPSIS approaches in customer-driven product design process[J].Computers in Industry,2008,59(1):17-31.(1) |
本文献在全文中的定位：
... 诸多采场结构参数优化方法中，以数值方法应用最为普遍[1]，它既克服了传统经验类比法主观性大的缺点，又改善了纯数理计算效 ...
1
本文献在全文中的定位：
... 率低的不足，并可与试验手段或理论方法有效结合，适用性很强.Zhao等[2]通过Ansys有限元试算构建采场结构参数与稳定性指标对应的目标函数， ...
1
本文献在全文中的定位：
... 数，利用混沌搜索获得最优参数方案，成功用于试验采场生产.罗周全等[3]采用MIDAS/GTS模拟开采进程，通过比较不同方案力学响应指标，从安全角 ...
1
本文献在全文中的定位：
... 应指标，从安全角度确定了某铜矿盘区隔离矿柱最优开采参数.吴爱祥等[4]基于大间距结构参数理论，开展相似放矿试验和放矿数值模拟，据此采 ...
1
本文献在全文中的定位：
... 数.然而，采场结构参数优化是涉及多层次、多因素和多目标的决策过程[5]，当前的优化方法仍存在一些不足：未准确考虑初始地应力场影响，而 ...
1
本文献在全文中的定位：
... 为准确获取中段原岩应力分布特征，采用改进型LUT套孔应力测定系统[6]进行现场实测.经与矿山地质及采矿工程师讨论，决定于N3号穿脉、S4号 ...
1
本文献在全文中的定位：
... ②断层.③铅锌矿和黄铁矿.岩体参数来自课题组研究结果[7]. ...
1
本文献在全文中的定位：
... 采用多元回归方法反演初始地应力场[8]，研究中考虑自重应力和构造应力，共6种载荷工况.其中，构造应力工 ...
1
本文献在全文中的定位：
... 于最小二乘原理，采用SPSS软件求解S残满足极小的方程组[9]，解得回归系数L1=2.333，L2=0.124， ...
1
本文献在全文中的定位：
... 以空场落矿为基本特征的采场，其稳定性状态多由顶板要素决定[10].将顶板简化为受上覆岩体自重Q作用的板状结构 ...
1
本文献在全文中的定位：
... 性力学小变形薄板理论，顶板受自重载荷Q作用，运用瑞利-李兹法[11]解出其挠曲函数： ...
1
本文献在全文中的定位：
... 为安全系数.根据H.Tresca屈服准则[12]，顶板某点产生剪切屈服时的主应力满足 ...
1
本文献在全文中的定位：
... 对于成本型指标集，Ti值越大，其解越优[13].计算6组初选方案评价指标集与理想解的向量相似度Ti为 ...
