冲击矿压是井巷或工作面周围煤岩体由于弹性变形能的瞬时释放而产生突然剧烈破坏的动力现象，常伴有煤岩体抛出、巨响及气浪等现象(姜耀东等，2014).它不仅能造成井巷破坏、人员伤亡、地面建筑物破坏，而且还会引起瓦斯、煤尘爆炸、火灾及水灾，干扰通风系统等，是煤矿重大灾害之一(窦林名和何学秋，2001).随着煤炭开采深度和强度的增大，目前冲击矿压已成为煤矿普遍的安全问题.

当前冲击矿压危险评估和监测的主要方法有：早期的数值模拟法、综合指数法及多因素耦合法(窦林名和何学秋，2001)、微震监测法(姜福兴等，2006；蔡武等，2014)、钻屑法(Gu et al．，2012)、钻孔应力法(刘金海，2014；曲效成等，2011)、电磁辐射法(王恩元等，2009；He X Q et al．，2011；Wang et al．，2011)、声发射法(窦林名和何学秋，2000；贺虎等，2011)、电荷感应技术(潘一山等，2013)等.近年来，层析成像技术被广泛应用于工程和地质诊断，现已成为井下采矿过程中勘探断层、松动圈、隐藏裂隙、应力状态等方面的有力工具(Meglis et al.，2005；Zhao，2000).

层析成像技术(Tomography)起源于1895年Wilhelm Conrad Roentgen(Stanton，1896)对X射线的发现.后来，Radon(1917)首次从理论上概括了层析成像概念：物体的内部结构可以通过分析单一轴面激发的射线穿透被测物体至另一边界过程中不同部位能量的变化进行重构成像.随着科技的发展，第一台医学CT仪器于1972年诞生.到1979年，科学家们将该技术应用于地球物理领域，即震动波层析成像技术(seismic tomography)(Dines and Lytle，1979).根据反演利用的波形参数(到时和振 幅)不同，震动波层析成像分为速度层析成像(velocity tomography)和衰减层析成像(attenuation tomography)(Westman et al.，1996; Westman，2004; Zhao et al.，2000).另外，根据震源的来源不同，速度层析成像又分为主动和被动两种(Luxbacher，2008).在井下煤矿开采中，主动震源一般由人工激发，包括炸药爆破(Cai et al.，2014; Dou et al.，2012;He et al.，2011)、锤击(Friedel et al.，1995; Friedel et al.，1997)、连续采煤机割煤引起的震动(Luo et al.，2009)等；被动震源一般采用自然发生的矿震(Bańka and Jaworski，2010; Cai et al.，2014，2015a; Hosseini et al.，2012a，2012b，2013; Lurka，2008; Luxbacher et al.，2008; Westman et al.，2012)，因此，被动震动波层析成像又称之为矿震震动波层析成像.

主动层析成像技术中震源位置已知(如锤击、爆破等)，射线覆盖几何形状及其覆盖密度人为可控，从而保证了该技术反演精度的可靠性，而且该技术可真正近似为一种瞬态反演.然而，由于主动层析成像在现场操作中需要人工激发震源和记录震源位置，同时操作时需停止一切正常采煤工序，且爆破(或锤击)时可能诱发煤岩动力灾害给施工人员带来的安全威胁，因此增加了该技术实施的劳动、经济及安全成本，进而阻碍了该技术在现场中的推广应用，同时也不适于长期远程监测.矿震震动波速度层析成像技术中的震源为开采引起的自然矿震事件，不需要人工激发，可随井下开采活动持续发生，这就保证了该技术的长期连续远程监测.然而，矿震事件的定位精度以及不充分或不规则的射线覆盖直接影响着该技术反演结果的精度及可靠性.

查阅文献资料发现，Friedel等(1995，1997)、Dou等(2012)、王书文等(2012)、Luo等(2009)采用主动震动波层析成像研究了波速与应力及强矿震分布之间的关系，发现高波速区与应力集中区及未来强矿震分布区对应一致性较好.Luxbacher等(2008)、 Hosseini等(2012a，2012b；2013)利用矿震震动波层析成像技术，以开采诱发的自然矿震作为激发源对工作面回采过程中的P波速度进行了反演，发现支承压力区与高波速区吻合较好.窦林名等(2014)、Cai等(2014，2015a)、Banka and Jaworski(2010)、Lurka(2008)结合微震实时监测系统，采用矿震震动波层析成像技术对煤矿工作面开采过程中的冲击危险区进行了实时探测评估，发现冲击或强矿震往往发生在高波速区；同时对卸压解危措施的实施及效果进行了指导和检验.相比于传统监测手段，矿震震动波层析成像技术探测范围广、成本低、监测过程安全、近似实时监测，因此该技术将成为探测评估冲击危险的一个强力手段，是未来国内外发展的一个趋势.然而，进一步完善该技术，以达到直接反映煤岩体内部真实应力场及定量描述冲击危险区域，还需做大量工作(Cai et al.，2015a).

综上所述，当前矿震震动波速度层析成像技术研究仅局限于间接或定性描述应力分布及潜在冲击危险区，该技术与传统监测手段对比如何，未见有文献报道.正如Luxbacher(2008)在其博士论文的展望中提到：“今后在开展被动震动波层析成像技术应用时，应与传统应力测量手段进行联合探测，因为这样才能实现震动波层析成像技术量化探测应力的目的，而不再是仅仅推断相对应力值的变化.”

基于上述出发点，本文采用矿震震动波速度层析成像技术详细研究了波速分布与数值模拟结果、钻屑量、钻孔应力值以及电磁辐射值之间的对应关系，力求进一步验证该技术评估冲击危险的可行性，同时也为进一步研究波速与应力或冲击危险的直接及定量关系提供参考.

利用开采引起的矿震事件和井下或地面安置的震动台站，根据台站与震源之间的距离L和台站接收到的初至时间T来反演“台站-震源”空间包络区域的波速分布V(x，y，z)，探测示意如图 1所示.

图 1

Fig. 1

图 1 矿震震动波速度层析成像技术探测示意图Fig. 1 Detection sketch of passive seismic velocity tomography technology

震动波在走时成像情况下以射线的形式在探测区域内部介质中传播(见图 1)，反演前首先精确计算矿震事件的位置，然后选取探测目标区域的矿震事件及其对应目标区域附近台站上标记的初至时间(P波初至或S波初至)作为反演原始数据，随后将“台站-震源”包络的空间区域介质划分为一系列小矩形网格，最终通过一个高度近似进行反演，其公式为(窦林名，2014；Luxbacher et al.，2008；Nolet，2008)：

通常，式(4)是一个欠定或超定方程组，求解此类方程组的有效算法一般是迭代算法.目前，大多数 引用的迭代算法有ART算法和SIRT算法(Gilbert，1972). 其中，SIRT算法不会出现奇异解，并能平稳收敛，计算简便快捷(Hosseini et al.，2012b)，为本文所采用.

河南义马跃进煤矿25110工作面采深1000 m左右，为25采区东翼第一个综放工作面，平均采高11 m，主采2号煤层.该煤层平均厚度11.5 m，平均倾角12°，煤层上方依次为18 m泥岩直接顶、1.5 m厚1-2煤、4 m泥岩和190 m巨厚砂砾岩老顶；下方依次为4 m泥岩直接底和26 m砂岩老底(如图 2所示).井下四邻关系(图 3)为：东为23采区下山保护 煤柱，南为25区下部未采煤层，东南部接近F16逆冲断层，西为25采区下山保护煤柱，北为大采空的25采区.其中，25110上巷布置于25090采空区下方煤层中，下巷接近F16逆冲断层.F16逆冲断层是受巨大挤压应力作用，在逆冲推覆作用下形成(Cai et al.，2015b)，为区域性逆冲断层.进一步的水平挤压应力作用使得该断层上盘岩层以断层面为支点发生翻转，最后呈现出直立(或倒转)的形态.该断层延展长度约45 km，走向近东西，构成义马煤田南部自 然边界，浅部倾角75°，深部倾角15°~35°，落差50~450 m，水平错距120~1080 m.图 3中地质剖面图所示为途径跃进煤矿的勘探线显示的F16逆冲断层剖面图.此外，跃进煤矿于2011年4月22日引进并安装波兰16通道ARAMIS M/E微震监测系统，台网布置如图 3所示.

图 2

Fig. 2

图 2 工作面煤层顶底板柱状 (J1：下侏罗统； J2：中侏罗统； J3：上侏罗统) Fig. 2 Roof and floor rock layers of coal seam (J1: Lower Jurassic; J2: Middle Jurassic; J3: Upper Jurassic)

图 3

Fig. 3

图 3 跃进煤矿采掘工程平面图、地质剖面图 T3: 上三叠统; K：白垩纪；微震监测系统台网布置: 实心正方形表示临时台站，可随着工作面的回采推进向前移动，实心三角形表示固定台站.Fig. 3 Mining engineering plan map, geological section map T3: Upper Triassic; K: Cretaceous, and layout of the microseismic monitoring system (squares indicate temporary stations that can be moved as the coal face advances and triangles indicate permanent stations) installed in the Yima Yuejin coal mine, Henan Province, China.

图 4

Fig. 4

图 4 25110工作面煤层位置的震动波速度层析成像切片图及未来矿震事件分布 Monitoring section表示反演期间工作面开采的区域： (a) 反演期间为2012年2月10日至2月24日，符号表示2012年2月25日至3月25日发生的矿震事件； (b) 反演期间为2012年4月1日至4月20日，符号表示2012年4月21日至5月21日发生的矿震事件； (c) 反演期间为2012年4月16日至5月8日，符号表示2012年5月9日至6月9日发生的矿震事件； (d) 反演期间为2012年5月8日至 6月7日，符号表示2012年6月8日至7月8日发生的矿震事件.Fig. 4 Plan view of seismic velocity tomograms at coal seam on LW 25110 Symbols show positions of the future seismic events. The monitoring section indicates total area mined over the inversion period. (a) Plan view of velocity tomogram at coal seam obtained from seismic events between 10 February 2012 and 24 February 2012. Symbols show positions of seismic events that occurred between 25 February 2012 and 25 March 2012. (b) Plan view of velocity tomogram at coal seam obtained from seismic events between 1 April 2012 and 20 April 2012. Symbols show positions of seismic events that occurred between 21 April 2012 and 21 May 2012. (c) Plan view of velocity tomogram at coal seam obtained from seismic events between 16 April 2012 and 8 May 2012. Symbols show positions of seismic events that occurred between 9 May 2012 and 9 June 2012. (d) Plan view of velocity tomogram at coal seam obtained from seismic events between 8 May 2012 and 7 June 2012. Symbols show positions of seismic events that occurred between 8 June 2012 and 8 July 2012.

数据采集所用的设备为现场安装的微震监测系统，如图 3所示，反演使用的震源为微震监测系统监测到的矿震事件.反演分析采用的是课题组自行编制的MINESOSTOMO软件(Gong，2010).反演计算之前，需要先筛选矿震事件，一般选取波形清晰(P波到时容易识别)、激发台站个数较多(至少5个通道)的事件.同时，为了尽可能减小反演计算模型的尺寸，提高反演效率，以及避免不规则的“震源-台站”空间分布引起反演结果的不可靠性，筛选事件时，一般选取发生在研究区域25110工作面的矿震事件，标记有效矿震事件通道时，一般标记研究区域附近的台站通道(8#、9#、12#、13#、14#、15#、16#).

影响反演计算结果的主要参数有(Hosseini et al.，2012b; Luxbacher，2008)：初始速度模型，“震源-台网”几何分布或射线覆盖密度，以及模型网格尺寸.其中，“震源-台网”几何分布因素在前期的数据筛选过程中就已经确定.此外，为了增加射线的覆盖密度，我们往往是通过增加震源事件数量或台站数量来解决.由于本文在反演计算迭代过程中选取的是SIRT算法，因此参与第一次迭代过程计算的初始速度模型显得尤为重要.模型网格尺寸参数直接影响着反演结果的精度，研究表明(Friedel et al.，1997；Williamson，1991)，理论上的最大反演精度为一个波长长度的距离.

根据现场原位试验确定的定位P波速度，赋予初始的常值速度模型值为4 km·s^-1.根据Nyquist采样定理，当采样频率大于信号中最高频率的2倍时，采样之后的数字信号才能完整保留原始信号中的信息，一般实际应用中保证采样频率为信号最高频率的5~10倍.目前现场微震监测系统普遍使用的采样频率为500 Hz，即能完整监测到的震动信号频率约为0~150 Hz，以P波速度值4 km·s^-1计算，得出P波波长数值为大于26.7 m.同时，考虑到研究区域震源平面最大定位误差约为30 m，再者考虑到煤矿井下实际条件，微震台站无法在垂直方向上得到很好约束，进而导致震源垂向定位误差偏大，其最大误差一般为70 m以上.综上考虑，决定选用网格尺寸划分参数为30 m×30 m×100 m.最后，为减小反演解的不确定性和提高反演计算效率，对最终反演解的最大波速值进行了约束，即6 km·s^-1.

如图 4所示为25110工作面回采不同时段煤层水平位置的波速分布切片图，图中不同符号表示了反演阶段未来一个月不同能量级别的矿震事件，Monitoring section所指的区域表示反演期间工作面开采的区域.由图可知，几乎所有的矿震事件(尤其是能量大于10⁵J的矿震事件)均发生在高波速区和高波速变化梯度区.再次验证了矿震震动波速度层析成像技术评估冲击危险的可行性.

为了验证波速反演结果推断应力分布，进而获得冲击危险区域的可行性，采用FLAC^3D数值软件对该矿适当简化后的25大采区模型进行了模拟.如图 5所示为数值计算模型，其尺寸为长1365 m×宽1050 m×高350 m，共267936个单元.对于重点研究区域(煤层巷道开挖区域)的单元采取了细化处理.模型各岩层力学参数及厚度根据图 2所示的岩层柱状描述赋予.

图 5

Fig. 5

图 5 FLAC3D数值计算模型Fig. 5 Numerical calculation model of FLAC3D

Anderson断裂机制表明(Anderson，1905)，逆断层的最小主应力轴为垂直状态，最大主应力轴与中间主应力轴为水平状态，受力形式如图 5所示.进一步结合F16逆冲断层应力场理论分析结果(Cai et al.，2015b)，作为一般性研究，赋予模型边界条件为：底部固定，最大水平主应力σ₁=29 MPa，中间 水平主应力σ₂=24 MPa，最小主应力σ₃=20.5 MPa.

经计算得出25110工作面回采过程中煤层附近层位的垂直应力分布，如图 6所示.对比图 4和图 6分析发现，图 6中显示的高应力分布与图 4中的高波速分布大体一致，尤其是两者对工作面超前支承压力分布带的反映基本一致.然而，两者整体上的分布仍存在一些差异，虽然，层析成像反演结果分辨率不如数值模拟结果，如数值模拟能精确确定采空区遗留小煤柱上的应力，但是波速层析成像是根据现场实际发生的矿震事件反演得到，更应该符合现场真实情况.其理由是，数值模拟是对实际开采技术条件、生产地质条件及煤岩力学参数的一种简化，如无法有效模拟开采速率对工作面周围应力的影响、煤岩体中隐伏的断层构造、煤岩体的非均质性及各向异性等，而地下采矿工程极为复杂，从而导致模拟结果大都是半定量化或规律性的结论.如图 6所示的数值模拟结果，随着工作面推进，垂直应力分布规律一般变化不明显，与实际工作面周围的应力分布偏差较大，尤其是很难给出巷道不同位置的冲击危险性.而图 4所示的层析成像结果显示，不同时段内的波速分布规律均不同，获得的信息量比数值模拟结果也丰富，尤其是能够反映不同时段不同区域内的冲击危险状态.从这一角度考虑，矿震震动波速度层析成像技术对现场的冲击矿压治理具有更好的指导意义.

图 6

Fig. 6

图 6 数值模拟结果 (a) 25110工作面推进580 m； (b) 25110工作面推进460 m.Fig. 6 Numerical simulation result (a) LW 25110 has been advanced 580 m； (b) LW 25110 has been advanced 460 m.

钻屑法理论基础表明(窦林名等，2001)，钻出煤粉量与煤体应力状态具有定量的正相关关系，即其他条件相同的煤体，当应力状态不同时，其钻孔的煤粉量也不同.

选取2012年4月16日至5月8日时段(图 7a)和2012年5月8日至6月7日时段(图 7b)的层析成像结果与相应时段内的钻屑量进行对比分析发现(其中钻屑量数据来自反演时段末10天左右实施大直径卸压孔时排出的煤粉量，钻孔实施地点为巷道上帮，将相邻区域各钻孔煤粉量累加后，最终取每米煤粉量作为钻屑监测结果)：

图 7

Fig. 7

图 7 钻屑监测结果Fig. 7 Results of drilling pulverized coal parameters

(1)图 7a显示，不同波速区对应的钻屑监测结 果差别较大，如高波速区钻屑量高，为70.5 kg·m^-1和65.2 kg·m^-1，而低波速区钻屑量较低，为41.5 kg·m^-1和43.2 kg·m^-1.钻屑量高的区域应力水平较高，与波速分布对应关系较好.

(2)图 7b显示，不同波速区对应的钻屑监测 结果差别较小，距停采线80 m左右位置为高波速 区，对应的钻屑量为65.2 kg·m^-1，其余的低波速 区对应的钻屑量分别为44.8 kg·m^-1、52.9 kg·m^-1和57.8 kg·m^-1.由此可见，该时段内钻屑量与波速分布同样具有较好的对应关系.

如图 4所示的波速反演结果显示，工作面前方均存在明显的波速分区特征，各自超前工作面的距离分别为170.5、185.0、112.8 m和142.2 m，平均152.6 m，这就是熟知的超前支承压力带.

由于工作面超前支承压力带特征是预测冲击矿压的主要依据，为了验证矿震震动波速度层析成像技术反演推断支承压力带的可行性，该矿于2011年10月26日在25110工作面运输平巷上帮煤体中安装了20组应力传感器，如图 8所示，图中共布置20组应力计，每组由两个应力传感器组成，安装深度分别为12 m和18 m，前10组安装间距为15 m，后11组间距为25 m.如图 9所示为工作面回采过程中11#、12#、13#应力传感器的监测结果曲线，从图中可以看出，当工作面推进距离应力传感器120 m至130 m，平均125 m时，各传感器数值均急剧上升，由此可以推断作为超前支承压力的影响范围.

图 8

Fig. 8

图 8 25110工作面运输平巷钻孔应力监测布置方案 图中共布置20组应力计，每组由两个应力传感器组成，安装深度分别为12 m和18 m， 前10组安装间距为15 m，后11组安装 间距为25 m.Fig. 8 Layout of 20 groups of pressure sensors installed in the headgate entry of LW 25110, at spacing of 15 m and 25 m for the former ten groups and the latter eleven groups, respectively. Each group consists of two pressure sensors installed in 12-m-depth and 18-m-depth boreholes, respectively.

图 9

Fig. 9

图 9 工作面回采过程中典型压力传感器(11#, 12#和13#)的监测结果Fig. 9 Monitoring results of the typical pressure sensors (11#, 12# and 13#) with the panel of LW 25110 retreat

值得注意的是，层析成像反演结果与钻孔应力监测结果相差较大，造成该差异的两大可能原因是：

(1)反演精度不同，震动波层析成像反演精度为30 m，即最终反演得出的超前支承压力范围应为122.6~152.6 m之间的某个数值；而钻孔应力曲线中统计的离散点数值为每天传感器连续监测数值的平均值，25110工作面每天平均推进1.2 m，即钻孔应力法监测得出的超前支承压力范围应为123.8~125.0 m之间的某个数值.

(2)震动波速度层析成像反演选取的原始数据一般长达20天，甚至一个月左右，该时段内工作面周围煤岩体结构发生改变，即波速也发生了变化，而反演过程依然将其设定为不变量，这将引起一定的误差.

综上所述，矿震震动波速度层析成像结果与钻孔应力监测结果基本一致，其数值约为125 m，可以用来反演推断工作面超前支承压力的影响范围，进而达到评估冲击危险的目的.

研究表明(王恩元等，2009；He et al．，2011； Wang et al．，2011)，电磁辐射信息能够综合反映冲 击矿压、煤与瓦斯突出等煤岩灾害动力现象的主要影响因素，其中电磁辐射强度主要反映了煤岩体的受载程度及变形破裂强度，脉冲数主要反映了煤岩体变形及微破裂的频次.

为了验证矿震震动波层析成像评估冲击危险的可行性，选取2012年5月20至5月30日的电磁辐射数据进行分析，如图 10所示为每天工作面前方电磁辐射的监测方案及其监测原理示意，图 11所示为监测时段内各监测点处电磁辐射强度最大值之和曲线以及平均值之和曲线.

图 10

Fig. 10

图 10 KBD5电磁辐射监测原理与方案Fig. 10 Monitoring principle and scheme of KBD5 monitor

图 11

Fig. 11

图 11 电磁辐射监测结果Fig. 11 Monitoring results of KBD5 monitor

对比图 11和图 4d分析发现：测点1到测点12范围的电磁辐射强度最大值之和以及平均值之和均较大，之间出现波动，其中最大值之和较为明显，达到了1752~3215 mV；在测点13到测点15区域内的电磁辐射强度最大值之和及平均值之和均较小，其中最大值之和较为明显，最小值为1322 mV；测点16位置，电磁辐射值较大，最大值之和达到了2222 mV.这与图 4d所示的波速反演结果对应较好，即高波速区对应的电磁辐射强度值较高，低波速区对应的电磁辐射强度值较低.

(1)对比矿震事件发现，强矿震主要发生在高波速区或高波速变化梯度区，再次验证了矿震震动波速度层析成像评估冲击危险的可行性.

(2)对比传统监测手段发现：波速反演获得的超前支承压力影响范围与数值模拟及钻孔应力监测结果对应较好；波速获得的信息量比数值模拟结果更丰富，尤其是能够反映不同时段不同区域内的冲击危险状态；波速值与钻屑量及电磁辐射监测值之间存在正相关关系，且对应一致性较好，即波速越高，冲击危险性越强.从而验证了矿震震动波速度层析成像评估冲击危险的可行性，同时也为进一步研究波速与应力或冲击危险的直接及定量关系提供参考.

(3)矿震震动波速度层析成像技术探测范围广、成本低、监测过程安全、近似实时监测，该技术将成为评估冲击危险的一个强力手段，是未来国内外发展的一个趋势.在将来的研究中，可联合主动层析成像技术与衰减层析成像技术，进一步获得更高精度的层析成像结果，进而为煤矿井下开采提供更安全更高效的生产环境.