基于贝叶斯神经网络的冲击地压预测与影响因素权重分析

乔建永，王志强，罗健侨，武 超，李敬凯，林 陆

(1.中国矿业大学(北京)，北京 100083； 2.北京邮电大学，北京 100876)

煤炭是我国发展依赖的主要资源，随着经济的发展，煤炭的需求量不断增大，截至2019年底，原煤产量38.5亿t，在煤炭资源需求稳中有升的背景趋势下，随着浅部水平和缓倾斜煤层的不断开采，越来越多的煤矿已经进入千米深部开采。在此趋势下，冲击地压动力灾害发生的频次和能级随着采深和开采强度的持续增加更加剧烈，严重威胁到煤矿的生产和工人的生命安全[1]。冲击地压启动发生机理成为国内外专家学者的研究热点，同时也是对冲击地压预测预报与防治措施的重要理论指导条件。

近年来，我国在冲击地压监测预警技术方面取得了长足进展，一些先进的、具有自主知识产权的监测设备被用于冲击地压监测预警。目前，我国用于煤矿冲击地压监测预警的主要方法有直接接触式监测法和地球物理方法，能够根据微震事件和围岩变形量等特征对冲击地压进行监测预警[2]。但是综合分析多影响因素耦合作用关系对工作面冲击地压启动传递过程的影响机理仍需深入探索。本文以华丰煤矿1411工作面这一数据记录全面的典型冲击地压现场条件作为案例，重点分析其冲击地压发生力源条件的影响因素对冲击能量的影响关系，并通过引入深度学习方法对冲击地压的启动与显现状况进行预测，弥补了综合指数法预测冲击的部分缺陷。

此外，通过模拟实验和解析式计算发现，对基本顶岩层的水力压裂弱化[3]或水射流切割减小悬臂梁结构对实体煤载荷的方法是一种能够更加有效降低1411工作面前方煤体力源、减小冲击启动风险的防冲措施。综上考虑多种冲击地压影响因素参量变化之间的互相量变作用，并且最终对冲击地压的启动和显现产生质的影响作用。以往的打分和综合指数方法并不能很好地考虑各因素变化相互影响的特点，因此利用神经网络方法进行分析，并且为了解决以往用来实现预测的BP神经网络遇到局部最优、过拟合及梯度扩散等问题[4]，改进连接函数使网络结构能够实现贝叶斯概率分析以及预测功能[5]。

本文针对华丰煤矿进行神经网络建模，并分析其影响因素权值，拟进一步实现贝叶斯网络结构的模型和计算方法在其他工况条件下的应用，并且随着输入集向量样本的增加和深入学习，神经网络结构的预测分析效果也会更具理论指导意义[6]，研究成果可为类似地质生产条件的工作面安全生产作业提供理论指导和技术支持。

1 工程背景

1.1 工作面天然地质条件

新汶集团华丰煤矿4#煤层1411综放工作面位于井田-1 000 m水平一采区三区段，下为尚未开采的1412工作面，1411工作面已于2012年停产结束回采。华丰煤矿位于两个全球性大断裂——郯庐断裂和昆仑断裂交点处，根据区域动力规划方法[7]研究得出，1411工作面开采范围内存在的6条断裂也是采掘过程中导致冲击地压发生的重要地质因素。4#煤层平均埋深960 m，平均煤厚6.2 m左右，平均倾角32°，具有强烈冲击倾向性，其直接顶具有中等冲击倾向性。1411工作面范围内煤层节理裂隙较发育，表现为单轴抗拉强度与单轴抗压强度之比较小(约为0.05)，易在支承压力和覆岩弯曲回转作用下产生破坏。图1为1411工作面综合柱状图。

前人在华丰煤矿冲击地压机理研究中将冲击发生因素分为两大类：一类是自然地质类因素，另一类是回采技术类因素[7]。超过千米埋深、上覆直至地表的巨厚砾岩、较大煤层倾角、中等冲击倾向性顶板、强烈冲击倾向性煤层都是华丰煤矿4#煤层发生冲击地压的天然倾向条件，因此在掘巷与回采阶段选择合理的开采工艺和防冲措施尤为重要。

1.2 工作面回采技术条件

由于顶煤易于冒放，1411工作面采用长壁后退式综放开采技术，全部垮落法管理顶板，采2.2～2.4 m，放3.8～4.2 m，正常割煤倾斜长度143 m，在采用错层位巷道布置方式的1411工作面推进速度由1 m/d增加到2.4 m/d。上方1410工作面整层开采条件下所有动力现象几乎集中发生在上平巷位置，图2为1410工作面上平巷发生的一起震级2级，能量2.2×107J的冲击地压事故平面图，结合倾角条件说明冲击地压防治重点在工作面上巷[8]。

图1 1411工作面综合柱状图Fig.1 Comprehensive histogram of 1411 working face

图2 华丰煤矿“9.9”冲击地压事故平面图Fig.2 Plan of Huafeng Coal Mine “9.9” rock burst accident

针对华丰煤矿冲击地压防治，通过优化采煤方法与工艺可以实现安全生产，其最显著因素来源于改革原实体煤巷道布置在采空区下方，这一特征实际上是改善了巷道围岩扰动的影响，即改善了动静载荷的作用特征[9]，因此，分析静载荷导致的弹性应变能量来源、动载扰动的诱发位置成为优化工作面开采技术参数降低冲击风险的必要途径。

2 工作面冲击启动机理与传递方式

工作面前方较近处煤体受超前支承压力的作用发生塑性破坏后，所受应力降低至原岩应力以下，这部分煤体不会聚积弹性应变能。而煤壁深处的弹性区应力升高位置，容易积聚大量的体变和形变弹性能量，在此部分煤体强度处极限平衡临界状态时，煤体中聚集的弹性能量也处于极限平衡状态。这部分煤体三向受荷载，由于未超过煤体抗压强度，应力在加载方向做功全部转化为煤体的形状改变和体积压缩。这部分能量在受到上覆直接顶或基本顶岩块下沉、回转、滑落、搭接过程中对煤体的铅直载荷缓慢、突然的增大或水平约束力突然的减小都可能会导致其猛烈的释放。根据煤壁前方煤体的应力及储能状态将其分为冲击能量释放-传递-显现过程的三个区域，分别为冲击启动区、冲击阻力区和被抛出的破碎煤体区，如图3所示。

图3 冲击能量启动-传递-显现分区Fig.3 Burst energy initiate-transfer-display zone

根据最小能量理论中对三维应力状态下煤岩体发生破坏的条件描述[10]，当煤岩体破坏启动后岩体应力调整，应力状态迅速转变为二维，最终转变为单向应力状态，三维状态下储存在岩体中的大量弹性能仅需消耗一维压缩或剪切所需的部分能量，其余能量则用于破坏或大变形塑性破坏；冲击地压启动理论[11]中对集中动载荷与静载荷型冲击地压的启动理论分析，认为集中静载荷可以独立导致冲击启动，而集中动载荷必须通过静载荷集中区完成。此外，工作面煤壁前方位于冲击启动区的煤体由于三向承载，选取其中任一单元煤体A进行分析，其内部聚积储存大量弹性应变能，处于满足冲击启动能量准则和强度准则的临界状态，此时外力做功对煤体A输入的弹性应变能U计算式为式(1)。

(1)

式中：σ1、σ2、σ3分别为最大主应力、中间主应力和最小主应力；E和μ分别为储能煤体的弹性模量和泊松比。当其承受最大主应力超过三轴强度极限，则单元煤体将发生破坏，储存在煤体内部的一部分弹性应变能Ud=U-Ue将在破坏过程中耗散；另一部分可释放弹性能Ue将对相邻煤体B做功，能量以做机械功的形式传递，如图4所示。

此时相邻单元煤体B在原有储存弹性应变能量U’的基础上继续受到单元煤体A破坏释放的部分能量Ue，当输入能量Ue和做功方向上的应力满足煤体B冲击破坏的能量和强度准则时，相邻煤体B也将发生破坏，释放能量并瞬时传递至下一相邻单元煤体，则冲击启动-传递-显现过程从冲击启动区某一单元煤体A开始依次传递，经过冲击阻力区能量逐渐衰减，再到不能储存弹性应变能的破碎区煤体能量传递过程结束，最终的能量做功将破碎煤体猛烈抛向采出空间。由冲击启动区到冲击阻力区结束，整个能量传递过程像一个破坏释能-衰减传递-能量叠加继续传递的链式过程，如图5所示。

图4 煤体能量释放“链式”传递做功Fig.4 Coal energy release “chain-type” transfer work

图5 冲击能量“链式”传递做功力学模型Fig.5 The mechanical model of “chain-type” transfer of energy

因此，研究工作面前方实体煤中的主应力分布和弹性应变能量场随天然地质因素和回采技术因素变化的演化规律对工作面冲击的预测、防治和机理研究有重要意义[12]。

3 冲击启动区应力与应变能计算

3.1 主应力计算力学模型

为求出工作面前方冲击启动区域单元体实体煤受主应力情况和应变能分布情况，根据弹塑性理论和叠加原理[13]，将弹性状态占绝大部分的工作面前方煤层视为理想状态下均质、各向同性体，并由弹性力学方法给出该半无限体。图6为冲击启动区应力计算力学模型。

图6 冲击启动区应力计算力学模型Fig.6 Mechanical model of stress calculation in burst-initiation zone

利用(0，q1-γH)、(x0，KγH)和(L2，0)求解两段直线斜截式后得到直线与坐标轴围成的三部分区域下任一点M(x，z)的应力分量，并取无限多段长度为dξ、高度为斜截式y值的微元面积，求解M(x，z)处的三向应力增量[14]，计算见式(2)。

(2)

式中：k1和b1为破碎区和冲击阻力区直线斜率和截距；k2和b2为冲击启动区直线斜率和截距；L1和L2分别为破碎区与冲击阻力区宽度和弹性区影响范围，计算结果见式(3)。

(3)

式中：γ为体积力，kN/m3；H为埋深，m；m为采厚，m；ξ=(1+sinθ)/(1-sinθ)；p1为支护强度，MPa；C为内聚力，MPa；φ为内摩擦角，(°)；f为顶底板摩擦因数；K为应力集中系数；β为侧压系数；τ0cotφ煤体的自撑力[15]，τ0取5 MPa。

由式(2)和式(3)计算得到三向应力计算结果，见式(4)。

(4)

根据现场地应力测量结果将三向应力计算结果代入应力主轴偏转公式(式(5))。

σi=(σj+Δσjcosα)sinθ+(σk+Δσksinα)cosθ

(5)

式中：α为最大主应力与水平方向夹角，(°)；θ为煤层倾角，(°)；σi为最大主应力、中间主应力和最小主应力，MPa；σj和σk为单元煤体所受铅直方向、水平方向和切向方向的载荷，MPa。将式(5)计算得到的最大主应力、中间主应力和最小主应力σ1、σ2和σ3代入到式(1)得到工作面前方非塑性区煤体内任意一点的弹性应变能量计算公式(式(6))，利用Matlab软件将上述公式迭代得到，由于弹性应变能量最终解析表达式过于冗长复杂，限于篇幅在此不予列出。

(6)

分析上述表达式可以得知，工作面前方冲击阻力区和冲击启动区单元煤体所受应力和聚积弹性应变能量的数值解大小，除了与埋深、煤岩体物理力学性质、采厚、煤层倾角和峰值应力集中系数等经典因素有关外，还和计算选取位置与工作面水平距离、顶板垂直距离以及顶底板岩层对煤层的层间力学性质等因素有关。

3.2 冲击启动区最大主应力和能量分布规律

以1411工作面实际生产地质参数条件进行理论数值分析，埋深为900～1 100 m，采厚为6.2 m左右，4#煤层弹性模量为2.2 GPa，内聚力为1.88 MPa，内摩擦角为38°，采用控制变量法对一项影响因素进行研究时，其他影响因素参数值始终控制为中间组方案大小，并且为了验证解析公式是否具有针对华丰煤矿的冲击现象普遍适应性，参数调整时上下阈值将1407工作面、1409工作面和1410工作面的相关参数也包含在内进行分析，并且利用FLAC3D软件对应力和能量的分布演化规律辅以佐证。

3.2.1 埋深对主应力和应变能量的影响

根据经验公式估算破碎4#煤层物理力学性质的单元煤体需要的能量不小于0.142 MJ，华丰煤矿1411工作面地质条件下发生冲击的临界深度处需约1.57×105倍于0.142 MJ的能量来破坏双向受力状态下的单元煤体，显然与实际情况相差较大。将实际参数代入式(6)计算得出埋深逐渐增加时对工作面前方煤壁不同位置处最大主应力和弹性应变能量的影响规律，如图7所示。

随着计算选取位置埋深的增加，最大主应力和弹性应变能基本与埋深成一定比例的增大关系，在1 000 m以深的1411工作面前方冲击启动区实体煤处所受最大主应力基本为上覆直至地表的岩体重量[16]。

将相同参数代入FLAC3D模型中计算发现，随着埋深的增加，应力和应变能与理论数值解析解也具有类似的变化规律，计算结果数值上也较为适应(图8)。

最大主应力和主应力差值是影响弹性应变能聚积的主要因素之一，分析图8中不同埋深对应不同的最大主应力场分布方式以及阈值可以得出，随着埋深的增加煤壁前方同一位置的最大主应力值随之增大，从21.5 MPa增加至31.5 MPa，弹性应变能峰值也从4.8×105J增大到1.2×106J。最大主应力位置也随埋深的增加而靠近工作面，使得冲击启动区所受扰动载荷更容易达到临界值，且在冲击能量传递做功时阻力区长度更短，塑性破坏更严重，从而削弱对冲击能量的阻碍作用。

结合能量聚积-储存-链式传递做功力学模型，分析图7和图8中埋深对1411工作面前方实体煤中的应力和储存能量的影响规律发现，上覆岩层累重使冲击启动区煤体处于高地应力作用下，巨大的围压使连续的单元煤体屈服应力呈正二阶导数式增长，可以储存更多能够用于发生冲击破坏的弹性应变能量。随着埋深的增加，诱导冲击启动区煤体发生破坏释放冲击能量需要的临界扰动载荷也更小[17]。

图7 埋深对主应力和应变能量的影响Fig.7 Buried depth affects principal stress and strain energy

3.2.2 覆岩力学性质对主应力和应变能量的影响

图9展示了关键层物理力学性质与主应力和应变能的动态变化关系。由图9可知，随着关键层岩层的抗拉强度和内聚力的增大，岩层对上覆巨厚砾岩重量的承载能力变强，传递到下方煤层中的应力也随之减小，单元煤体内聚积的弹性应变能量也减小，而且可以看出，抗拉强度和内聚力对最大主应力和应变能量的影响效果在数值上显著小于覆岩重量的影响，结合前述研究成果分析表明，关键层对上覆500～960 m巨厚砾岩层的重量承载效果相对有限，巨厚砾岩层对工作面前方煤体施加的静载荷也是1411工作面冲击启动的重要影响因素之一。

图9 关键层物理力学性质与主应力和应变能的 动态变化关系Fig.9 The dynamic relationship between the physical and mechanical properties of the key layer and the principal stress and strain energy

覆岩关键层中细砂岩层厚度由40 m增大至70 m时，弹性应变能量从2.279×106J增大至2.663×106J；抗拉强度增大至25 MPa时，弹性应变能量从4.9×106J减小至2.3×105J；同时可以看出关键岩层的内聚力的变化对弹性应变能量的影响并不显著。在数值模拟模型中表现出对最大主应力和弹性应变能量的影响效果如图10所示。由图10可知，随着覆岩强度的增高对上覆直至地表岩层运动和重量的承载效果也越好，因此煤层承受的铅直应力也显著减小，最大主应力峰值从78.96 MPa减小至69.21 MPa。但当老顶发生周期性破断时，作用在煤层上的动载荷也会更加明显。因此应该从减小工作面前方煤体应力集中程度和增加动载扰动输入能量两个方面综合考虑顶板覆岩对冲击地压启动与传递的影响。

图10 随顶板岩石力学性质整体增强的 工作面前方最大主应力分布规律变化Fig.10 Changes in the distribution law of the maximum principal stress in front of the working face with the overall enhancement of the mechanical properties of the roof rock

3.2.3 工作面前方冲击启动区与工作面中心位置距离对主应力和应变能量的影响

针对1411工作面前方煤体的储能状态进行分区讨论以及解析式数值计算，发现冲击启动区中心位置与工作面的距离与选取计算位置M(x，z)的最大主应力和弹性应变能量之间有明显动态变化关系。

图11展示了冲击启动区位置与最大主应力和弹性应变能量关系。分析图11可发现，随着x值的增大计算选取位置向煤壁深处逐渐移动，最大主应力先增大至峰值71.2 MPa，后逐渐减小直至30 m左右后回到与原岩应力相近的40 MPa，与工作面前方铅直应力分布规律类似但峰值点位置和数值大小不同，弹性应变能量呈现先增大至最大主应力约4.5×106J后趋于稳定，这是因为峰值应力位置煤体处于极限平衡状态储能效果差，峰后曲线下的冲击启动区煤体处于弹性受载状态，随着应力趋于原岩应力状态，能量增量也逐渐放缓最终减小趋于稳定。

图12展示了冲击启动区距离对最大主应力和弹性应变能的影响。由图12可知，随冲击启动区与工作面距离的增加，最大主应力与弹性应变能的峰值大小和分布范围也逐渐减小，据此规律可以认为，通过增大冲击启动区中心位置与工作面之间的水平距离，即增大冲击阻力区的长度，可以降低工作面冲击风险，达到通过钻孔卸压、煤层注水或水压致裂切割顶板实现增大冲击阻力区范围、降低工作面冲击风险的效果。

图11 冲击启动区位置与最大主应力和 弹性应变能量关系Fig.11 The relationship between the burst-initiation zone position and the maximum principal stress and elastic strain energy

3.2.4 煤体力学性质对主应力和应变能量的影响

煤岩体从受载聚积储存能量到破坏释放能量的过程都与单元体的力学参数密切相关，选取解析式中的煤体内聚力和内摩擦角与最大主应力和应变能的关系进行研究，影响关系如图13所示。由图13可知，随着煤体内聚力的增大，单元煤体承受的最大主应力一定程度减小但并不显著；而单元煤体储存弹性应变能量的能力随煤体内聚力的增加显著提高，当内聚力达到2.5 MPa时，单元煤体能够储存约4.4×106J弹性应变能。

3.2.5 保护层开采对主应力和应变能量的影响

综合分析1411工作面所处顶底板赋存条件，应优先选择无冲击地压的下方6#煤层作为保护层开采。开采保护层后，在被保护层中确实受到保护的地区，可按无冲击地压煤层进行采掘工作。6#煤层开采对1411工作面前方煤体最大主应力和应变能分布规律的影响如图14所示。 由图14可知，开采保护层对于上方4#煤层应力集中程度和弹性应变能的减小影响十分明显，弹性应变能峰值从5.752×106J降低至3.22×106J。

3.2.6 回采速度对主应力和应变能量的影响

工作面回采速度对采场前方煤体能量转移和释放的速率，以及应变能量场的聚积程度均有影响[18]，不论是从顶板频繁回转滑落失稳带来的动载扰动，还是从采动影响的微震事件总量显著增加的现象来说，较高或者过低的回采速度对于冲击地压防控均属不利因素。保持合理的工作面回采速度需要结合采场前方应变能量和微震事件扰动两个方面综合考量(图15)。 由图15可知，随着工作面推进速度由0.5 m/d增加到3 m/d，工作面前方弹性应变能量峰值从2.14×106J增大到3.17×106J，峰值位置分布范围也增大，且根据微震事件发生的频次和能级分析，当工作面推进至见方和周期来压布局时，微震呈现成簇的规律性发生。

图12 冲击启动区距离对最大主应力和弹性应变能的影响Fig.12 The influence of impact start zone distance on the maximum principal stress and elastic strain energy

图13 煤体力学与最大主应力和弹性应变能的影响关系Fig.13 The relationship between coal mechanics and the influence of maximum principal stress and elastic strain energy

图14 保护层开采对弹性应变能量分布的影响Fig.14 The influence of protective layer mining on elastic strain energy distribution

图15 工作面推进速度对工作面前方弹性应变能量分布和微震事件的影响Fig.15 The influence of the advancing speed of the working face on the elastic strain energy distribution and microseismic events in front of the working face

3.2.7 开采厚度对主应力和应变能量的影响

随着总采厚从2 m增大到7 m，弹性应变能聚集量由3.45×105J显著增加至1.697×106J，最大主应力也显著增大，由于采厚的增加，支承应力峰值位置前移，应力集中程度增高。对照数值模拟得出的结果可以看出，随着开采厚度的增大，工作面前方应力集中程度和分布范围均显著增大，而且工作面后方覆岩能量集中程度也增大，聚集更大弹性应变能量的覆岩在发生滞后破断时会以机械波的形式将衰减后残存的能量传递至巷道围岩或工作面前方煤体，并在满足动载扰动最大能承应力时诱发冲击(图16和图17)[19]。

图16 开采厚度与最大主应力和弹性应变能的影响关系Fig.16 The relationship between mining thickness and maximum principal stress and elastic strain energy

4 贝叶斯神经网络分析与预测

在动静载叠加启动和分源防治力学结构下[20]进行冲击力源能量的分析预测过程中，冲击启动的弹性应变能量判据和动载扰动传递能量临界能量判据受距离和地质条件因素影响有量级上的差距，难以利用指数或打分方式统一判断冲击是否发生。因此建立贝叶斯前向神经网络结构结合前述应力与应变能量解析式，对1411工作面冲击地压的影响因素权重占比和任一位置的冲击启动区能量进行计算并预测冲击是否可能发生。

4.1 反馈贝叶斯神经网络影响因素权重分析

根据谷歌Deep Mind团队在深度学习领域研究进展[21]，通常利用贝叶斯公式作为激活函数推论神经网络所需要的训练量非常大，并且极难求解，此类神经网络并不适合进行精确推导和权值求解。取而代之的是使用精确贝叶斯的近似变分求解方法进行深入学习，通过改变激活函数的概率模型和权值计算方法，改进得到适用于多参量归一化分析预测学习的反向贝叶斯神经网络(BBN)代入前述建立的影响因素与冲击启动区弹性应变能量的量化联系，对各因素权重占比进行分析。BBN网络模型结构如图18所示。

图17 不同采厚工作面前方弹性应变能分布云图Fig.17 Counter of elastic strain energy distribution in front of different mining thicknesses

图18 反馈贝叶斯神经网络模型结构图Fig.18 Feedback Bayesian neural network model structure diagram

由图18可知，通过经典后向传递过程求得权重修正值的方式和改进后通过BNN方法求得权值分布的方法，每个连接函数的求解过程中的权值都有独立的分布方式相比直接赋予定值更加灵活，学习的精度更高，弥补了普通前向神经网络的缺陷使其更加适用于对多因素的冲击影响源耦合训练与分析。

根据前述解析公式计算得到的数据作为贝叶斯神经网络训练的训练集，数值模拟运算结果作为验证集将现场实测数据作为网络结构的测试集。验证集用于对网络结构训练过程的矫正、连接函数参数迭代修正赋值及计算得到权值的独立分布；测试集用于校验每次训练结果是否满足精度、噪声、迭代次数等要求，满足则结束迭代过程，见表1[22]。

表1 贝叶斯优化算法过程Table 1 Bayesian optimization algorithm process

现场实测值给出的后验预测密度见式(7)。

(7)

取预测分布的均值作为实测结果z的预测值则见式(8)。

(8)

反馈贝叶斯神经网络结构(BBN)进行前向深度学习训练的过程相当于利用连接函数中的权值计算调整理输入层与输出层间的量化关系，一定程度上弥补了解析公式和数值模拟方法的不足。后向反馈反复迭代的过程即是改变权值从属的分布方式以及连接函数的常系数以优化网络结构使训练结果精度满足要求(图19)。

图19 反馈贝叶斯神经网络模型计算结构示意图Fig.19 Schematic diagram of the calculation structure of the feedback Bayesian neural network model

由于建立的冲击影响因素权重的网络结构有11个父结点两个子节点，结点数量较多连接结构比较复杂，为了学习的最终结果足够精确，训练过程可能需要1.00e07左右次量级的迭代，因此学习过程中的数值计算可由Matlab软件实现。输入层训练集数据为11×200的矩阵数表，由于篇幅限制仅列出冲击启动能量的上下限阈值和一组较有代表性的中间组实验数据见表2。

表2 神经网络输入层训练集数据Table 2 Neural network input layer training set data

BBN计算后验概率选择利用贝叶斯正则化方法，则连接函数适于处理输入值在(0，1)或(-1，1)区间范围内的数据组，而11项影响因素参数大小存在量级上的较大差异，因此需要利用归一化公式处理数据，见式(9)和式(10)。

(9)

(10)

式中：θ与S分别为用于连接神经元的树突上的条件事件和概率；x为输入向量参数；y为输出向量参数；min为某种影响因素输入集的下限阈值；max为某种影响因素输入集的上限阈值。

图20展示了反馈贝叶斯神经网络训练精度与效果。由图20可知，经过6小时43秒的1 000次迭代过程，BNN神经网络输出数据的总体精度达到R=0.854 39，结果相关性置信程度较高，同时由图20(b)可知，输出集置信水平较高而且根据较集中的置信区间范围宽度认为输入集样本数量足够大。由图20(d)、图20(e)和图20(f)可知，贝叶斯神经网络的训练学习精确度较高，误差较小且总体分布在0.1附近，拟合程度达到0.754 64，且认为用于预测趋势的网络结构精确程度达到0.854 39。综合上述贝叶斯神经网络的学习训练效果认为，训练结果得到的连接函数结构、隐含层层数、中间层层数、函数的选择以及修正权值的独立分布均能满足对输入层数据进行分析和趋势预测的要求。根据Matlab软件NNtool功能，直接输出各输入层神经元的权重，计算结果见表3。

分析BNN神经网络训练结果得到的权重数值以及排序情况可以得到，1411工作面前方煤体中储存着可以释放做功的弹性应变能量力源主要来自于深埋积累着的上覆岩层体积重量，埋深占到弹性应变能量聚积来源的约22%。工作面超过千米的深埋条件同时也是硐室采场空间挖出后前方峰值应力集中系数K较大，顶底板高应力集中程度以及顶底板岩层和煤层冲击倾向性较强的先决条件之一。因此在理论解析式的计算结果中，埋深与冲击启动区弹性应变能量的直接关系虽然并没有达到22%，但由于埋深同时间接影响着近水平方向主应力、冲击启动区中心位置水平距及其他并未被解析计算公式考量在内的因素，因此通过贝叶斯神经网络方法发现埋深是1411工作面弹性应变能量聚积和影响冲击地压启动的权重占比第一的影响因素。

根据煤层柱状图分析，工作面上方关键层是否能良好承担上覆直至地表的巨厚砾岩层同样也决定着煤层中的应力分布和弹性应变能量聚积情况。当关键层岩梁力学模型足够长时，抗拉强度对描述抵抗上覆岩层载荷的层间错动以及铅直方向位移的力学性质比较直观，并且当关键层自身厚度较大时，传递到煤层的载荷也会增大，关键层的物理力学性质以及巨厚砾岩层重量力源对工作面前方冲击启动区煤体的主应力和应变能量分布的影响属于相互作用相互影响的关系，在神经网络训练分析后得到各因素的独立权重分布，其中上覆巨厚砾岩的厚度权重值约为14.88%，关键层抗拉强度权重值约为9.48%，而关键层内聚力占比6.2%。

结合数值模拟实验的结果分析，对6#煤层保护层是否开采的试验方案在全部影响因素中对弹性应变能量的影响达到了9.08%；工作面不同的推进速度对冲击启动区弹性能量的影响达到7.33%。另外通过理论公式计算和模拟结果发现，冲击启动区中心位置的水平以及纵向距离对冲击启动也有着6.70%左右的影响。将11项影响因素权重分别统计如图21所示。

其中，天然地质类影响因素占76.57%，回采技术类因素占23.43%。天然地质类因素的成因是力源来自于自然条件，并非人为无法调控，如煤体、覆岩的力学性质、冲击启动区中心位置的距离等影响因素可以通过防冲卸压钻孔、水射流切顶割缝预裂、对高应力顶底板及煤体岩层水压致裂的方式改变覆岩以及煤体的力学性质并且扩大冲击阻力区和破碎区的长度，对冲击启动区能量传递至巷道或工作面做功起到良好的阻碍效果。因此从防冲措施技术条件看，又将影响因素分为人为可调控和不可调控的因素，其中，人为可调控因素占48.76%，不可调控因素约占51.24%。根据权重分析结果得出的人为可控因素中对冲击启动区弹性应变能的控制手段优先级排序为：①开采6#煤层保护层对1411工作面进行卸压防冲；②对工作面前方冲击阻力区煤体采用钻孔卸压、水力压裂等防冲措施，改变煤体的力学性质以降低其储存、释放和传递弹性应变能的能力，同时也可以达到增加冲击阻力区宽度阻挡冲击启动区煤体猛烈释放弹性应变能的冲击作用；③降低回采速度，可以一定程度减小工作面前方煤体弹性应变能量聚积范围和能量密度，并且增大冲击阻力区长度和范围以减弱冲击能量链式传递的显现剧烈程度；④合理控制开采厚度，尤其是合理优化采放比，减少源于坚硬厚顶板覆岩的断裂带来的微震事件频次和平均能量，弱化动载扰动对冲击启动区煤体的破坏引起的能量释放；⑤通过水射流预裂顶板、直接钻孔注水压裂或预置裂缝定向高压水压裂等方法降低顶板覆岩的物理力学性质以降低冲击分区内煤体的峰值应力集中情况、弹性应变能密度集中程度和分布范围来减轻自发型冲击的静载力源。

4.2 贝叶斯神经网络冲击能量与启动预测

相比于对多因素影响的权重占比求值，BNN网络对于训练组数据处理的逻辑过程更为快速可靠，在前述建立成功的贝叶斯神经网络的基础上，将归一化后的训练组数据正向迭代拟合，神经网络逻辑流程、预测结果分别如图22和图23所示。

对采区有记录发生冲击能级的15处位置进行预测检验，结果发现BNN神经网络成功预测了14处冲击能量与显现程度，认为对华丰煤矿一采区相近地质条件与开采工艺下的冲击启动显现预测准确程度达到93.3%。但是对于其他煤层工作面或其他采区的预测，该训练完成的BNN并不具备准确预测的功能，需要根据上述过程重新计算输入集数据并训练新的贝叶斯神经网络。

图20 反馈贝叶斯神经网络训练精度与效果Fig.20 Feedback Bayesian neural network training accuracy and effect

表3 冲击地压影响因素权值计算结果Table 3 Weight calculation results of influencing factors of rock burst

图21 影响因素权重占比Fig.21 Proportion of influencing factor weight

4.3 对于巷道冲击地压预测的研究展望

冲击地压事故90%以上发生于巷道之中[23]，近5年来有记录的冲击地压事故均发生在工作面超前回采巷道。由于文中华丰煤矿采用错层位负煤柱式巷道布置方法，使得巷道在回采期间“有震无灾”，因此并未对巷道围岩稳定性和冲击启动发生情况深入讨论。在此为全面考虑冲击地压的预测以及影响因素权重分析的方法，使本文在对工作面的冲击预测防治起到作用的同时也可以充分考虑巷道与采场其他条件的综合作用对采出空间远近场围岩冲击能量聚集、启动、传递和显现的影响，引入巷道蝶形冲击地压理论[24]对巷道掘进时与掘出后出现的蝶形塑性区及其增量进行描述(图24)。

巷道蝶形冲击地压理论是根据巷道横截面蝶型塑性区边界曲线为安全临界曲线给出巷道围岩稳定性复动力学预警量化解析方法。为量化巷道围岩应力应变场分布范围，乔建永等[25]将塑性区域边界成为位移变换的最小不变子集(Julia集)，并根据Kastner方程[26]刻画了巷道横截面岩土的弹塑区域的解析规律。

图22 BNN预测结构模型Fig.22 BNN predictive structure model

引入上述理论方法以期后续的研究工作中可以将蝶形冲击地压理论对巷道围岩冲击能量与分布范围的定量刻画与贝叶斯神经网络方法相结合，实现对巷道冲击发生的时空预测和能级预测。

5 结 论

1) 工作面开切眼后，采出空间导致工作面前方原岩应力重新分布形成支承压力的不同区域，并且升高位置的支承压力将煤体破坏形成不会聚集弹性应变能的破碎区；支承压力附近处于极限平衡区的煤体可以聚集少量弹性应变能量，总体上对冲击能量起到阻碍和间接传递的作用并不会直接释放大量弹性应变能，该范围内煤体属于冲击阻力区；而释放能量直接对阻力区煤体做功，并将能量通过机械波和位移做功形式传递至破碎区，将塑性煤体抛向工作面的弹性状态单元煤体组成的区域是冲击启动传递的力源位置，为冲击启动区。

图23 BNN预测结果Fig.23 BNN prediction results

图24 掘进巷道冲击地压发生机理模型Fig.24 Mechanism model of rock burst in driving roadway

2) 冲击启动区部分煤体处于弹性应变能量储存的极限状态，受到外部扰动应力达到一定输入能量时，即打破极限平衡态，能量释放对周围煤体“链式”传递做功。据此过程，建立工作面前方煤体应力计算模型，并求解处支承压力简化曲线下任意点位置的最大主应力、中间主应力和最小主应力，而后可以求得该位置聚集的弹性应变能量大小。发现随着如煤体、覆岩物理力学性质以及采厚、工作面回采速度等影响因素的变化，冲击启动区煤体弹性应变能量并非随之正比例变化，而是呈现出动态变化规律，这说明了影响因素之间也在随着研究参量的变化相互影响。因此，需对冲击地压启动发生显现起到影响作用的影响因素机理综合研究分析权重占比。

3) 结合弹性应变能量解析式计算和现场实测数据记录，形成神经网络数据输入集，并根据输入集的矩阵维度建立合适的神经网络结构选取最适于分析权重的激活函数后得到反馈贝叶斯神经网络结构(BNN)，对训练集数据进行处理后，得出1411工作面冲击地压发生的影响因素中，最重要的因素分别是：埋深→覆岩厚度→覆岩抗拉强度→保护层开采→煤体内聚力→回采速度→采厚→冲击启动区中心位置水平距→覆岩内聚力→煤体内摩擦角→冲击启动区中心位置纵向距。

4) 根据训练完成的BNN网络对新的输入向量进行拟合预测，得到的弹性应变能预测值与解析公式计算结果相近，并且输入向量数据处的现场实际记录结果与神经网络预测结果也基本一致。认为根据现场实测的天然地质参数和回采技术参数训练完成的贝叶斯神经网络对于煤矿冲击地压的启动预测具有较好的效果。