循环爆破累积损伤诱发似层状铝土矿采空区群失稳机制

姜立春，沈彬彬，陈敏

(1. 华南理工大学土木与交通学院，广东广州，510640；2. 华南理工大学安全科学与工程研究所，广东广州，510640；3. 广东省冶金建筑设计研究院有限公司，广东广州，510640)

山西沉积型铝土矿床整体上呈层状、似层状产出，矿体平均厚度薄，条带式空场法是常用的开采方法，回采结束后，地下将遗留大面积采空区群[1]。爆破作业是主要的破岩方式，随着工作面的不断推进，相邻采场爆破应力波将对采空区群稳定性产生一定的影响，应力波的破坏作用主要体现在扩展节理裂隙、劣化岩体物理力学性质、降低安全系数等方面[2-4]。长期循环累积损伤作用将破坏薄弱单元空区的矿柱和顶板稳定性，进而诱发整个空区群垮塌，引发地表塌陷[5-7]。因此，开展循环爆破累积损伤诱发似层状铝土矿采空区群失稳机制研究，有利于防控此类矿山灾害，提高矿山本质安全水平。

目前国内外学者开展了大量循环爆破损伤的研究[8-13]。胡刚等[8]依据监测实验获取的振速时程曲线，采用HHT 理论将其变换为瞬时能量谱，分析隧道围岩结构爆破累积损伤效应；曹峰等[9]基于HJC本构模型，引入损伤变量，模拟循环荷载作用下小净距隧道保留岩体的累积损伤演化过程，探讨了岩体损伤程度、损伤范围与爆破次数之间的关系；左双英等[10]采用数值模拟的方法探讨了圆形隧道的围岩损伤分布规律及影响范围；MALMGREN等[11-12]研究了爆破作业后矿山开采区的损伤分布范围及岩体的损伤；凌同华等[13]提出基于多元判别分析理论的岩体爆破震动损伤安全判据方法。总体来说，该领域的研究主要依靠数值模拟方法，从理论角度开展循环爆破损伤下诱发采空区群失稳的相关研究相对较少。

采空区群是由多个单元空区开挖形成的地下复杂结构体，是在地下原生矿岩上开挖而成，具有鲜明的结构动力学特征[14]。采空区垮塌实质是在重力和外部荷载作用下，内部结构发生破坏[15]。本文作者在现场调研基础上，构建铝土矿采空区群结构模型，引入矿岩结构损伤因子Dn，研究循环爆破作用下空区群内部薄弱单元空区顶板、间柱的应力及位移规律，研判在相同爆破当量和爆源条件下爆破频次与采空区群稳定性的关系，利用数值模拟计算和工程监测的结果验证分析理论计算结果的科学性，指导矿山安全生产。

1 结构动力学理论分析

1.1 结构模型构建

铝土矿地下采空区群由r个单元采空区组成，各单元空区由1条间柱、1个空区、1个顶板和1个底板组成，各空区的顶板和底板分别连成一个整体[15]。空区上部受覆岩重力P0作用，受相邻采场的爆破应力qd间歇作用。图1所示为山西某地下铝土矿矿山采空区群实景图，内含3个单元空区；图2所示为其结构受力分析示意图。

图2中，x轴沿空区群走向方向为正，y轴沿重力方向为正；l为空区顶板跨度；h为矿柱高度；b为宽度；ω，ν和θ分别为任意节点i的x轴方向位移、y轴方向位移和xy平面内转角。

为简化采空区群计算分析，进行如下假设：

1)采空区群岩体为连续、弹性介质，内部构造完整，仅考虑循环爆破对岩体物理力学性质的劣化作用，不考虑节理裂隙扩展损伤作用。

2)循环爆破的炸药当量相同、爆源位置相近，同次爆破应力波在同一岩体内部传播速度相等。

3)重力荷载做功既定，仅对比分析爆破荷载做功。

1.2 位移函数构建及求解

利用离散相似法[16]对图2所示的采空区群模型进行离散化处理，采用位移法进行变形位移求解。各梁和柱交界点的i和j的位移函数式[17]为

式中：{v}m，{ω}m及{θ}m分别为节点位移分量矩阵；Xm(y)，Ym(y)及Zm(y)为广义梁函数；a为单元空区数量；r为采空区群中单元空区总数量；[Ф]为离散函数[11]，[Ф]=[Ф0(x)，Ф1(x)，Ф2(x)，…，Ф(m)x]；m为沿水平采空区群的条带间柱数量。

将式(1)转换为矩阵形式：

式中：[N]为结构形函数；[ε]为节点广义位移向量。

在x-y坐标系下，采空区内的单元梁柱e的位移阵列由各连接节点i，j的位移列阵构成，即

采空区群内部单元空区的某子单元(梁或柱)e位移矩阵为

则某梁或某柱单元e的总势能泛函数Пe为

式中：[T]为坐标变换矩阵；[k]e为单元刚度矩阵；[F]e为单元荷载矩阵。

1.3 远场爆破荷载做功矩阵

爆破荷载qd取值主要是由矿山每次爆破炸药当量而定。由文献[18]可知，单位爆破等效荷载表达式为

将式(6)转换为矩阵形式为

式中：Δt为设爆破荷载振弦周期，Δt= 0.001 s，则在n次爆破作用下，爆炸扰动力做功矩阵为

根据最小势能原理，采空区群结构总势能泛函数П等于各梁和柱组成的全部单元总势能Пe之和，则

根据最小势能原理，得

1.4 远场爆破岩体损伤修正

研究表明[19]，应力波与岩体的速度和岩体物理力学性质有关，若岩体结构致密、形态完整、物理力学强度高，则应力波在介质中传播速度就大，反之则小。因此，理论上可以利用应力波在空区岩体中的三轴矢量传播速度变化，间接描述采空区岩体的损伤程度。

爆破作业后岩体的损伤因子D为

式中：v0为损伤前岩体应力波传播速度；v为损伤后岩体的应力波传播速度。

假设采空区某处岩体受到n次爆破作用，则采空区岩体的损伤因子Dn可表示为

式中：Dn为n次爆破后岩体的损伤因子；vn为n次爆破后应力波在岩体内部传播速度；ηn为n次爆破后应力波速度与第一次爆破后应力波相比的速度下降比。

考虑爆破损伤作用，修正采空区整体结构总势能泛函数，可得修正后空区群的总刚度矩阵[：

2 采空区群失稳判据及分析

顶板、间柱失稳容易引起单元采空区失稳。长期循环爆破作用，各单元采空区内部岩体结构容易产生不同程度损伤，损伤严重的单元空区将成为采空区群中的薄弱单元空区。因此，利用岩石强度理论的极限抗拉强度、极限抗压强度和容许极限位移等作为研判指标，分析薄弱单元空区的顶板、间柱的稳定性，可以推测薄弱单元采空区的稳定状态，进而分析整个采空区群的稳定性。

2.1 顶板

空区形成后，在上覆矿岩重力荷载P0作用下，顶板将发生不同程度的下向挠曲位移变形，当挠曲变形幅值超过一定阈值时，顶板会发生失稳垮塌。由于岩体的单轴抗压强度远大于抗拉强度，通常顶板的变形破坏主要是由拉伸破坏引起的。为此，下面分别从顶板极限抗拉强度、挠曲位移变形等角度进行分析。

1)抗拉强度。由结构力学理论知，采空区顶板最大拉应力σmax一般出现在其中心位置[20]：

式中：E为矿体弹性模量；I为顶板截面对中性轴的惯性矩；Wz为顶板抗弯截面系数；θi和θj分别为顶板两端角位移，可由单位位移阵列[ε]求得。

当顶板岩体的最大拉应力σmax超过其承受的拉应力阈值[σ]t时，岩体塑性破坏将导致顶板垮塌[21]。因此，σmax须满足

式中：Kv为岩体完整性系数；[σ]t为岩体的极限抗拉强度；ka为安全系数，根据文献[22]，这里ka取2.2。

2)挠曲位移。顶板任意一点的位移表达式为

式中：x为顶板任意一点的位置横坐标，当x取l/2 时，可求得顶板中点的竖向挠曲位移。

当顶板的竖向位移变化超过容许极限位移时，单元采空区将进入临界失稳状态，根据地下硐室容许极限位移[ω]的经验公式[23]，顶板的竖向最大挠曲位移变化幅值ωx=l/2应满足：

式中：f为岩体的普氏系数。山西沉积型铝土矿体属中等坚固岩体，这里f取4[24]。

2.2 间柱

在上覆围岩重力和爆破扰动力的共同作用下，间柱将产生不同程度位移变形，x方向上，条带间柱变形主要是沿矿柱轴线方向的压缩变形；y方向上，沿矿柱垂线空区方向变形(变形量相对较小)。间柱严重变形将发生失稳坍塌，破坏方式主要有局部垮落、滑移破坏、整体坍塌等。因此，这里从抗压强度方面进行分析。根据文献[25]可知，矿柱轴向应力分量σr表达式为

式中：qd为爆破应力波扰动力；γ为矿柱容重。

当矿柱所受的轴向压应力σr超出了其承载的压应力阈值[σ]时，将发生塑性变形破坏。σr应满足

式中：σc为岩体的抗压强度参数；α为同规格岩柱所对应的尺寸系数，取值与岩柱的宽度高度比密切相关。当b/h＜5时，α=1.0；当b/h≥5时，α=1.5[26]。

2.3 综合判据

综上因素分析，这里可以给出似层状采空区群的失稳判定方法。

1)当薄弱单元空区顶板、间柱二者任一处于不稳定状态时，则采空区群处于不稳定状态。

2)当薄弱单元采空区的顶板、间柱均处于稳定状态时，顶板的竖向位移变化幅值小于容许竖向变形阈值[ω]时，则采空区处于稳定状态。

3)当顶板、间柱二者均处于稳定状态时，顶板竖向位移变化幅值大于或等于容许竖向变形阈值[ω]时，则采空区处于亚稳定状态。

综上所述，整体结构模型判断采空区群稳定性的计算过程如图3所示。

3 计算分析

相邻采场的循环爆破的累积损伤作用将导致空区构成要素(顶板、间柱)的位移、最大拉应力、最大压应力变化，可以借助结构理论计算和数值计算方法进行分析，对比计算结果，确定循环爆破次数与采空区群稳定性状态变化的定量关系，从而为似层状铝土矿采空区群稳定性分析提供理论依据[15]。

3.1 参数选择

根据矿山初步设计及实地勘察数据，以图2为例，选取图2采空区群及周边赋存围岩的物理力学参数，根据实测数据确定采空区2 为薄弱单元空区，现场实测初始损伤因子D1为0.86。

1)采空区群几何参数。空区群长度为42 m；单元空区数量为3个，跨度为5.0 m，高度为4.0 m；条带间柱宽度为4.0 m，长度为100.0 m(图1)；

2)铝土岩体物理力学参数见表1。

表1 铝土矿岩物理力学参数Table 1 Mechanical parameters of bauxite

3)重力荷载P0。采空区上覆黄土、砂岩层厚度为120.0 m，自重应力为2.7 MPa，作用位置施加在矿柱顶板上部岩体接触面(图2)。

4)爆源及应力波。爆源为相邻采场，距相邻采空区的距离为30 m(图2)。爆破炸药为2 号岩石乳化炸药，密度为1 000 kg/m3，爆破应力波初始爆速为3 200 m/s。单次爆破炸药当量为20 kg。

3.2 结构力学模型法计算

利用Matlab 工具包，结合Newmark-β自编程序，求解式(10)和(13)，得到采空区整体坐标下各节点的位移矩阵[ε]；联立式(14)，(16)和(18)，分别计算不同爆破频次条件下，顶板最大拉应力、位移及间柱最大压应力。

3.2.1 顶板位移

通过式(16)计算并记录采空区2 顶板中心位置监测点的位移。图4所示为爆破次数n与顶板位移ω响应关系。

由图4可知：在爆源距离不变的情况下，顶板位移ω与爆破次数n呈正相关关系。随着爆破次数n增加，顶板位移ω也不断增大，且位移增加速率增大。

由式(17)计算结果可知，该单元空区顶板的竖向容许极限位移[ω]=6.5 mm。当爆破次数达到13次时，顶板位移接近竖向容许极限位移，当爆破次数达到14 次以上时，顶板位移大于竖向容许极限位移，采空区处于不稳定状态。

这表明相邻采场的循环爆破应力波严重劣化采空区内部岩体的物理力学性质，岩体的抗爆破应力波冲击能力急速下降。

3.2.2 顶板最大拉应力

式(14)记录采空区2 顶板中心位置的监测点的最大拉应力。图5所示为爆破次数n与顶板最大拉应力σmax响应关系。

由图5可知：在爆源距离不变的情况下，顶板最大拉应力σmax与爆破次数n呈正相关关系。随着爆破次数n增加，顶板最大拉应力σmax也不断增大，且顶板最大拉应力增加速率增大。

由式(15)计算结果可知，该单元空区顶板的极限抗拉强度[σ]t=2.357 MPa，当爆破次数达到16 次时，顶板的最大拉应力接近顶板抗拉强度；当爆破次数达到17 次以上时，顶板最大拉应力大于顶板抗拉强度。

上述分析表明，在环爆破作用下，顶板存在失稳风险，在不采取其他保护措施下，容易发生采空区坍塌事故。

3.2.3 间柱最大压应力

式(18)记录采空区2 间柱中心位置的监测点的最大压应力。图6所示为爆破次数n与间柱最大压应力σr响应关系。

由图6可知：当爆破距离不变的情况下，矿柱最大压应力σr变化与爆破次数n呈正相关关系，随着爆破次数n增加，矿柱损伤产生径向裂纹，应力发生重新分布，间柱的最大压应力σr不断增大，且间柱最大压应力变化速率逐渐降低。

由式(19)计算结果可知，该单元空区间柱极限抗压应力[σ]=8.34 MPa，当爆破次数达到14 次时，间柱最大压应力接近其极限抗压强度，当爆破次数达到15 次以上时，间柱最大压应力大于极限抗压强度。

这表明循环爆破作用下，间柱稳定性受到极大威胁，处于不稳定状态。

3.3 数值模拟法计算

现有技术条件下，数值分析计算是分析复杂工程体稳定的一种有效手段。为了验证上述理论分析结果的可靠性，以图2为例，利用Midas-GTS/NX，构建内含3单元采空区的数值模型(图7)，同样设定采空区2为薄弱单元采空区，选取同一采空区群及周边赋存围岩的物理力学参数和爆破相关参数，利用损伤因子Dn对空区岩体物理性质折减，分析薄弱单元空区关键节点的位移、应力变化规律。

3.3.1 数值模型构建

1)边界设置及破坏准则。模型四边设置黏滞阻尼边界，顶部设置自由边界，底部设置固定边界，岩体破坏服从Mohr-Coulomb准则。

2)爆源位置及爆破次数。爆源处于距离采空区30 m远的附近区域，爆破次数为20次。

3)应力和位移监测点。考虑选定采空区群的对称性，在采空区2顶板和条带间柱中部位置分别设置监测点1 和2，记录顶板和条带矿柱的应力和位移变化情况。

3.3.2 结果分析

循环爆破作用下薄弱单元空区2的顶板挠曲位移、最大拉应力、矿柱最大压应力变化曲线如图4～6 所示。从图4～6 可以看出：数值模拟法与结构力学模型法的结果基本吻合，变化趋势相似。

3.4 结构力学模型法与数值模拟法对比分析

综合图4～6 分析可知：当爆破次数在12 次以下时，采空区顶板、间柱强度均小于其极限强度，顶板位移处于允许范围，采空区处于稳定状态；当爆破次数为13～15次时，采空区顶板、间柱二者存在处于稳定状态情况，顶板位移超过允许位移阈值内，采空区处于亚稳定状态，但需要一定的支护处理；当爆破次数大于16次时，采空区间柱、顶板均处于不稳定状态，且顶板位移远远超过允许位移阈值，采空区处于失稳状态。

4 工程验证

选取山西某地下铝土矿+1 160 m 水平生产中段3号采空区群为对象(图1)，设置2个监测点，监测相邻采场循环爆破作业对采空区稳定性的影响。经现场调查，该区域采矿方法为条带空场法，采空区位于铝土矿矿体内，距离地表高度152 m。空区群长度为122.8 m，宽度为86.4 m；单元空区跨度为4.2 m，高度为4.5 m；条带矿柱长度为86.4 m，宽度为4.0 m。

2 个监测点分别位于条带采空区顶板及间柱，位置如图8和图9所示。监测仪器为应力计GZY25和顶板沉降仪Yud300。监测内容为矿柱、顶板位移应力变化。监测时间为2019年4月份至9月份。监测初始，采空区稳定，周边围岩结构完好。

图10 所示为采空区围岩裂缝实景图。根据采空区群失稳判据及分析计算可得，当爆破次数小于8 次时，采空区处于稳定状态；当爆破次数为9～13 次时，采空区处于亚稳定状态；当采空区大于13 次时，采空区处于失稳状态。现场实地监测结果显示，当相邻采场进行10 次爆破作业后，位移监测仪器表明距离爆源21 m 处采空区顶板竖向位移增大，幅度为8.1 mm 左右，采空区顶板最大拉应力为2.83 MPa，采空区矿柱最大压应力为8.76 MPa，21 m 处采空区处于失稳状态，43 m 处采空区顶板最大拉应力为1.562 MPa，间柱最大压应力为7.25 MPa，顶板竖向位移增幅为7.3 mm 左右，43 m 处采空区处于亚稳定状态，现场监测结果与理论计算结果相符。

在监测时间窗口后期，21 m 处采空区区域多次发生顶板失稳坍塌，条带矿柱侧壁出现大量爆破裂纹，矿柱力学性质严重劣化，采空区群出现严重坍塌问题(图11)。工程实例说明采场10 次爆破后，采空区处于亚稳定状态，在一段时间窗口期后，采空区群坍塌，验证了理论计算结果的合理性。

5 结论

1)在现场调研的基础上，构建沉积型似层状铝土矿采空区群结构模型，首次引入岩体损伤因子Dn，研究循环爆破对采空区群顶板、间柱的累积损伤效应，依据顶板的位移和拉应力以及间柱的压应力临界阈值，研判采空区稳定状态，分析在相同地点、相同炸药当量条件下爆破频次与采空区群稳定性的关系。

2)当相邻采场爆源的距离和炸药当量一定时，采空区的间柱最大压应力、顶板的最大拉应力及顶板的位移均随着爆破次数增加而增加；当爆破次数小于12 次时，采空区顶板的最大拉应力、间柱的最大压应力均小于其极限值，顶板的位移处于允许范围内，采空区处于稳定状态；当爆破次数为13～15 次时，顶板位移超过允许位移阈值6.5 mm，采空区顶板、间柱二者均处于稳定状态情况，采空区处于亚稳定状态；当爆破次数大于16 次时，采空区间柱、顶板均处于不稳定状态，采空区处于失稳状态。

3)数值模拟法与结构力学模型法计算的位移、应力响应曲线变化趋势相似。现场工程验证表明，在+1 160 m 水平生产中段爆破作业后，采空区群严重失稳，采空区区域多次发生顶板失稳坍塌，且条带矿柱侧壁出现大量爆破裂纹，间柱力学性质大幅劣化，采空区群出现大面积坍塌，验证了结构力学模型法计算结果的合理性。