利用极化率参数监测煤矿滞后突水的可行性

刘志新， 薛国强， 张小楷,4

1 中国矿业大学资源与地球科学学院， 徐州 221116 2 中国科学院地质与地球物理研究所， 中国科学院矿产资源研究重点实验室， 北京 100029 3 中国科学院大学地球与行星科学学院， 北京 100049 4 贵州省地矿局第二工程勘察院， 遵义 563000

0 引言

在煤炭开采过程中，煤层底板的完整性会遭到一定程度的破坏，在裂隙发育区域，往往导致工作面底板突水事故的发生.特别是采空区，由于煤炭开采活动便底板产生一定深度的破坏区域，同时也会使原有隐伏构造的裂隙逐步扩大，如果隐伏构造与底板下承压含水层连通，承压水会沿着隐伏构造裂隙涌入工作面采空区，形成采空区底板滞后突水灾害(刘树才，2008；胡彦博，2020)，如图1所示.如果工作面底板存在隐伏断层、陷落柱等构造时，煤层开采活动可能诱发的滞后突水将更加严重.

图1 采空区滞后突水示意图Fig.1 Schematic diagram of lagging water inrush

滞后突水是矿井生产中经常遇到的特殊水害类型，因具有一定的隐蔽性、突发性，对矿井生产的危害极大.近几年来，煤矿工作面滞后出水发生较为频繁，华北石炭-二叠系煤层随着煤炭资源的高效产出，矿井逐渐向深部延伸，地层中广泛发育的岩溶及高压奥灰岩溶水，均严重威胁安全开采，其中以滞后性突水事故造成的人员伤亡和经济损失最为巨大(虎维岳，2010；虎维岳和尹尚先，2010).工作面滞后突水严重威胁着矿井安全生产，轻者使突水工作面关闭停产，严重影响矿井正常生产；重者会造成淹井事故，危及财产和工人生命安全.

地球物理方法在与环境领域勘查中发挥着重要的作用(Xue et al., 2007；底青云等，2020).通常采用地成电磁方法对煤矿含水结构体进行的探测或监测(Xue et al., 2013; Chen et al., 2015, 2019a,b).目前，在采煤工作面回采之前需要进行矿井地球物理勘探，采明工作面煤层内部、两侧及顶、底板采动影范围内是否存在构造异常或导含水通道，常用的矿井电磁法技术主要有：矿井瞬变电磁，无线电波透视，矿井直流电法，矿井音频电透视技术等(薛国强等,2007,2008,2021; Xue et al., 2019).但是，对于规模尺度较小的隐伏构造，由于探测技术与资料解释精度限制，传统的矿井地球物理方法不能对其进行有效识别，可能导致后期采空区滞后突水事故的发生.

在 20 世纪 90 年代中期就研发了煤层底板突水征兆监测设备，实现底板位移、钻孔水压和弹性波等的监测，开展了水害监测预警相关理论、硬件和软件等方面的研究(胡耀青等，2000；姜福兴等,2008),先后研发了松散层突水监测预警系统、底板突水灾害预警系统等(乔伟等,2020).针对华北型煤田突水类型及特点，提出了监测预警指标体系和矿井突水监测预警模式(刘德民,2015；刘德民等，2019)，即根据监测指标将水害预警划分为物理预警模式和化学预警模式两种主要类型，并制定了隔水层电阻率参数动态全面监测和关键部位单点多参数监测的现场实时监测方式.另外，以底板“下三带”理论为基础，提出了集多频连续电法充水水源监测、“井-地-孔”联合微震采动底板破坏带监测以及监测大数据智能预警为一体的煤层底板突水三维监测与智能预警技术思路(靳德武等，2020).目前，基于电磁法的矿井突水预测(杨天鸿等,2007)，主要是在已有巷道空间内进行电位或电阻率参数连续观测，实现矿井突水过程的实时监测，同时可实现对水质进行快速检测识别，判断水源类型.本论文提出的基于极化率参数的突水过程监测系统，将接收电极或电磁传感器通过钻孔布置在采空区范围内，也可利用下抽巷空间向两侧钻孔，采用行列式方式提前布置电极或电磁传感器，实现对采空区滞后突水过程的时实动态监测.

论文主要采用三维激发极化法装置形式，利用理论分析、数值模拟和物理模型实验的方法，研究了采空区底板下陷落柱构造在不同突水阶段极化率参数的变化规律，论证了利用极化率参数监测煤矿滞后突水的可行性.

1 极化率参数实时监测技术

1.1 激发极化效应

激发极化法是利用地壳中岩、矿石不同激电效应差异为基础，外部人为供直流电(时间域)或交流电(频率域)，通过研究和分析其激电场分布规律进行找矿和解决地质问题的一个电法勘探方法，本文正是基于时间域激电法，在矿井地质条件下进行相关研究分析.

激发极化效应(激电效应)是在外电流场作用下，岩、矿石与周围的水溶液产生物理-电化学反应，引起地质介质中出现电荷分离，产生附加“过电位”的一种物理化学现象.见图2，A和B为直流供电点，可通过MN测量总电位差ΔU(t)(一次电位ΔU1与二次电位差ΔU2(t)之和)和二次电位差ΔU2.极化体被极化是一个先快后慢的过程，经过Ts秒后趋于饱和状态.断电后，ΔU1已不存在，但极化引起的电位差ΔU2(t)先快后慢的逐渐减小，直至消失.

图2 激发极化现象示意图Fig.2 Schematic diagram of induced polarization phenomenon

一般为了得到激发极化的全过程，需要充电的时间T足够长(一般几分钟)，然后断电，应测量供电时和断电后电位差变化的全过程.但在实际中，这样费时费力，成本不划算，所以一般只取若干个时间点测量的电位差间接求得某些参数来描述激发极化特性.常用的参数是视极化率ηs：

(1)

1.2 极化率参数正演模拟算法

底板突水激发极化法的三维有限元正演模拟是利用“等效电阻率法”进行替换求解模拟.第一步，利用软件数值模拟得到一次场ΔU1(不含激电效应)；第二步，根据“等效电阻率法”替换相应地质体的电阻率，用软件再次数值模拟计算得到总场ΔU(含有激电效应)；第三步，根据软件两次模拟计算的得到的结果，用总场减去一次场得到相应的二次场ΔU2，然后由视极化率公式得到模拟视极化率.

激发极化法的视极化率参数正演计算是利用电阻率来求解的等效电阻率法.由于极化率η是二次场电位差ΔU2与总场电位差ΔU的比值：

(2)

因此，设大地电阻率为ρ，极化率为η，未产生极化效应的一次场电位差为：

(3)

设K为装置系数：

(4)

电阻率为：

(5)

当供电时间使岩矿石达到饱和时，测量电极MN之间的总场电位为：

(6)

其中ρ*为等效电阻率，由式(6)可得：

(7)

根据极化率η的计算式，可得η、ρ和ρ*三者间的关系：

(8)

于是有：

(9)

利用ρ*代替ρ来求含激电效应总场极化总场电位差为：

(10)

将ρ*代替ρ来求含激电效应总场的方法叫“等效电阻率法”，由式(3)、(10)得二次场电位差为：

(11)

两式相除，得极化率：

(12)

上面通过等效电阻率法计算得到的极化率η即为激发极化法三维正演模拟的视极化率.

1.3 监测系统设计

采空区滞后突水观测系统是采用三维激发极化法工作方式，观测系统采用单极-单极工作方式，即供电点A放置于采煤切眼中间位置，接收电极(电磁传感器)可通过有线或无线方式传输采集的电位数据，记录点位置为接收电极正下方.如果主采煤层下方存在瓦斯抽放巷道，可在下方巷道内通过钻孔向两侧按照行列式方式布置电磁传感器，具体位置如图中红点所示.如果是单一主采煤层，可利用采煤切眼两端巷道空间向采空区布置2～4个水平钻孔，将电磁传感器通过钻孔放置在指定位置.本文中要分析两个水平钻孔Z1、Z2内的电磁传感器对应极化率电性参数变化规律，如图3所示.通过两个钻孔中电磁传感器处极化率参数的变化规律及对应水平位置，实现对底板下隐伏突水构造进行实时动态监测的目的.

图3 底板监测系统示意图Fig.3 Schematic diagram of monitoring system for bottom plate

2 极化率参数实时监测数值模拟分析

2.1 模型建立

(1)整个地电模型为长方体，长宽高分别为400 m、300 m、300 m；根据煤矿采空区顶底板岩层分布，分别设置三层，底部为底板，厚度190 m，中部为煤层和采空区，厚度为10 m；上部为煤层顶板，厚度为100 m(图4)，各参数设置见表1.

表1 正演模型参数表Table 1 Forward model parameters table

(2)陷落柱形态设置为圆柱体，半径为30 m，极化率为8%，其顶部距离采空区底板高度h(图4).在模拟计算过程中h值分别为100 m、80 m、60 m、40 m、20 m 5个高度变化模拟水位上升.陷落柱中心点到煤层与采空区边界(采煤工作面)横向距离S，计算是分别取值130 m、100 m、70 m、40 m，模拟含税水陷落柱靠近采煤工作面；陷落柱中心偏离采空区中心线的纵向距离为L，取值0 m、40 m、80 m、110 m，模拟纵向含水陷落柱不同位置；以S和L的变化模拟陷落柱的不同位置对电位和极化率影响及其规律研究.

图4 模型剖面示意图Fig.4 Model profile diagram

(3)如图3所示，供电采用单点电流源A，电流大小10 A.理论上可分析采空区范围内任一测点的电性参数变化规律.考虑工作的可行性，主要研究Z1和Z2两个钻孔内各测点的参数变化情况.Z1、Z2是在工作面两侧巷道端向采空区内布置两个水平交叉钻孔.

2.2 底板突水过程监测模拟

2.2.1S=100 m，L=0 m，不同高度h引起的视极化率分布

图5中，A为供电电流源，假设水位逐渐上升，水位到煤层底板的距离h值分别设定为100 m、80 m、60 m、40 m、20 m，不同水位情况下的极化率动态参数变化分布如图5a—e所示.从图5中可以看出，陷落柱异常体投影在纵向和横向上位于采空区中心部位，陷落柱异常随着h降低(水位上升)逐渐出现并越加明显，极化率异常与陷落柱投影位置靠近.当h=100 m时极化率异常基本不显示，当h小于80 m时极化率随着水位上升，异常响应越来越明显.

图5 A3供电h变化底板极化率监测等值线图Fig.5 Monitoring isoline map of goaf floor polarizability of A3 power supply with h changed

分别提取不同水位情况下(水位到煤层底板的距离h)的异常极化率极值，获得极化率极大值随水位h变化曲线图(图6)，h值从100 m减小到80 m时(水位升高20 m)，视极化率参数值增加了6.0%，曲线斜率为0.268；h值从80 m减小到60 m时，视极化率参数值增加了11.04%，曲线斜率为0.522；h值从60 m减小到40 m时，视极化率参数值增加了19.5%，曲线斜率为1.024；h值从40 m减小到20 m时，视极化率参数值增加了34.1%，曲线斜率为2.142.极大值随水位高度增加而变化，水位上升，极大值增长率和曲线斜率越来越大，底板极化率变化愈加明显.

图6 极化率极大值随水位h变化曲线图Fig.6 Curve of maximum polarizability variation with water level h

各模型Z1测线和Z2测线上异常极化率变化图7所示，曲线总体呈现中间高两侧低的开口状，极化率值随着h减小而升高，在两条钻孔方向异常极值位置与陷落柱异常体位置相同，可通过两条钻孔测线不同方向的异常极值位置圈定含水陷落柱异常体的位置.

图7 模型Z1测线和Z2侧线上异常极化率变化Fig.7 Abnormal polarizability changes on Z1 and Z2 lateral lines of the model

2.2.2h=40 m，L=0 m，不同横向距离S引起的视极化率分布

图8为A3为供电电流源，异常陷落柱逐渐靠近采煤工作面，横向距离S分别取130 m、100 m、70 m、40 m变化的动态极化率参数变化分布图.随着陷落柱逐渐靠近采煤工作面，极化率异常响应越加明显，极化率异常中心逐步向采煤工作面移动，位置偏向异常体在底板投影位置附近偏离供电电极一侧，且随着横向距离S变小，极化率异常越来越明显.

图8 A供电横向距离S变化采空区底板极化率监测等值线图Fig.8 Monitoring isoline map of goaf floor polarizability of A power supply with S changed

提取图中极化率极大值，以不同横向距离S为横坐标，对应的异常极化率极大值以S，视极化率极大值为纵坐标得到图9.当S从130 m到100 m时，极化率极大值增长为6.2%，两点曲线斜率为0.771；当S从100 m到70 m时，极化率极大值增长5.6%，曲线斜率为0.743；当S从70 m到40 m时，极化率极大值增长1.3%，曲线斜率为0.186.随着异常体靠近煤层，当S从130 m逐渐靠近到S=40 m时，以30 m为步长，视极化率极大值增长率和曲线斜率由大变小.

图9 极化率极大值随横向距离S变化曲线图Fig.9 Curve of maximum polarizability variation with transverse distances S

随着横向距离增大(从40～130 m)，由于异常体极化效应较强，异常体对采空区底板极化率影响越来越小，导致底板总体极化率减小.当陷落柱靠近采空区与煤层边界时，靠近供电电源，极化率值相应增大.各模型Z1侧线和Z2侧线上异常极化率变化如图10所示，两条钻孔侧线方向曲线规律相同，极化率曲线总体呈现中间高两侧低的开口状，极化率的值随着S增大而逐渐减小，极大值对应钻孔坐标位置远离巷道方向移动，位置接近异常体中心.

2.2.3h=40 m，S=100 m，不同纵向距离L引起的视极化率分布

图11为异常体投影在横向上处于采空区底板中部，纵向上距离L取0 m、30 m、60 m、90 m变化的动态极化率参数变化分布图.图中显示极化率异常中心区域位于异常体投影远离供电电源一侧，随着陷落柱纵向距离L增大，异常位置向采空区上边界移动，异常体位置和极化率异常位置同步移动，随着L的增大，异常极值小于6%的幅度在减小.

图11 A3供电纵向距离L变化采空区底板极化率监测等值线图Fig.11 Monitoring isoline map of goaf floor polarizability of A3 power supply with L changed

从等值线图中看不出明显参数值变化，提取图中极化率极大值，以纵向距离L为横坐标，异常极化率极大值为纵坐标，得到图12.当L从0 m到30 m时，异常极化率极大值增长-0.12%，曲线斜率为-0.033；当L从30 m到60 m时，异常极化率极大值增长-1.0%，曲线斜率为-0.26；当L从60 m到90 m时，异常极化率极大值增长-5.1%，曲线斜率为-1.327.极化率异常极大值保持减小，增长率越来越小，总体值得变化幅度很小.

图12 极化率极大值随纵向距离L变化曲线图Fig.12 Curve of maximum polarizability variation with longitunal distances L

A供电随着纵向距离L增大(从0～90 m)，由于异常体极化率高，极化效应影响采空区底板极化率分布，水位没有变化的情况下，纵向距离L变化对底板极化率影响较小.所以当陷落柱靠近采空区边界时，极化率值变化不大.各模型Z1侧线和Z2侧线上异常极化率变化如图13所示，Z1方向上异常极值位置位于异常体在采空区中的投影位置远离巷道一侧，极大值随着L增大呈现先增大后减小的趋势，位置逐渐向右侧移动.Z2方向上极化率曲线随着L增大而减小，极值位置向着巷道移动.根据两条钻孔方向上极值位置做平行于采空区上边界的直线能大致圈定异常体的位置.

图13 模型Z1测线和Z2侧线上异常电位变化Fig.13 Abnormal potential changes on Z1 and Z2 lateral lines of the models

3 物理模型实验

3.1 实验监测方案设计

实验模型长宽高大小分别为380 cm×300 cm×180 cm，分为三部分：上部是覆盖层(黏土)30 cm；中部是弱含水层(细砂)130 cm，在其内部存在导水陷落柱(充填砾石)大小为100 cm×100 cm×130 cm，导水陷落柱位于模型中部；下部为含水层(粗砂)10 cm和底部隔水层(黏土)10 cm，图14是物理模型的设计图，通过从下部含水层和供水管道陷落柱异常体中持续供水，右侧安装水位监测管监测水位变化.

图14 物理模型实验平台设计图Fig.14 Physical model experiment design

实验仪器选用法国IRIS公司生产的SYSCAL电法仪，其中SYSCAL Pro单元是专门针对高效电法勘探所设计研发生产的电阻率&IP测量计.SYSCAL Pro单元能够测量一次场电压和衰减电压曲线值，从而可以计算出电阻率和极化率.

电极布设：根据观测系统的设计方案，本次实验电极布置方式采用三维激发极化法观测布极方式，观测系统采用单-单极方式，供电电极A位于36号点，另一供电电极B装置于室外20 m以外.供电电流为1 A.电极在模型上表面以间隔40 cm距离均匀布置，根据实验模型大小设置为7排9列，供63个电极.其中边界电极与模型边界距离都为30 cm.在模型的对角位置上布置两条斜向测线Z1和Z2，主要分析两条斜向测线Z1和Z2位置上的测点极化率曲线变化规律，如图15所示.

实验步骤如下：

(1)根据电极布设方案在模型顶面布设电极，如图15所示.

图15 实验电极布设图Fig.15 Experimental electrode layout

(2)向模型底部含水层充水，从水位监测管观察水位的高度.

(3)通过观察水位监测管，水位从80 cm开始第一次监测数据采集，上升过程中由于含水层水的渗流，降低水的注速度，保持水位不变直到完成数据采集.以水位上升20 cm为一个周期，重复上述采集过程完成80 cm、100 cm、120 cm、140 cm水位数据采集工作.

(4)完成数据采集工作之后，结束数据采集工作，整理仪器设备，进行数据反演处理与解释.

3.2 实验数据分析

首先对实验数据采进行三维成像，可以得出随着陷落柱构造内容水位上升，极化率逐渐增大.由于模型顶面由干燥的黏土组成，当水位较低的时候，水位对顶面的影响较小，随着水位上升，模型含水层水面上升，模型整体极化率逐渐升高.为了更清晰反应极化率参数随水位变化规律，提取了三维立体数据中深度为0.7 m处的水平切片数据，得到-0.7 m平面极化率等值线分布情况，如图16所示.在图16a是水位80 cm时-0.7 m深度截面极化率分布等值线图，从图中可以看出中部陷落柱位置出现一个低极化率异常区域，左侧出现一个高极化率异常区域和右侧靠近边缘出现一个低极化率异常区域.由于前期实验模型进行了第一次物理模型实验，排水时间不够导致模型弱含水层孔隙水含量高，反之陷落柱内部孔隙直径较大，残留的水含量相对含水层低.加之此时模型底板供水水位为80 cm，-0.7 m界面上水位未达到，因此陷落柱异常体极化率异常显示为低极化率异常.左侧高极化率异常是由于模型建好之后先进行了污染物模拟实验所残留的污染物吸附于细砂颗粒中引起.右侧设置了排水管，相对左侧排水效果更好，孔隙水含量低，导致了右侧低极化率异常.

图16 物理模型-0.7 m平面极化率分布等值线图Fig.16 Physical model -0.7 m plane isoline map of polarizability distribution

图16b、c分别为水位100 cm和120 cm时-0.7 m深度截面极化率分布等值线图，相对于图16a，图16b、c陷落柱位置及右侧低极化率异常区域减小，左侧高极化率区域增大，等值线密度增大.此时水位到达截面-0.6 m深度，超过了-0.7 m深度截面，由于陷落柱水位上升，其与四周含水层保持接触，水由于水压较高向四周扩散，导致深度为-0.7 m的截面极化率截面显示为高低极化率异常均存在的情况.图16d为水位为140 cm时-0.7 m深度截面极化率分布等值线图，与图16a、b、c相比，图16d低极化率异常区域几乎不存在，只有右侧存在很小区域低极化率异常，极化率等值线图中心陷落柱异常体位置极化率由于水位上升而极化率升高.物理模型实验表明极化率参数在监测底板陷落柱突水过程的变化规律与数值模拟结果相似，能较好的反应底板隐伏构造内水位上升情况，证明该方法对监测采空区底板下隐伏构造滞后突水的有效性.

图17为深度-0.7 m时，与数值模拟相同的两条钻孔测线Z1和Z2在不同供水水位情况下极化率曲线.从图中可以看出，在Z1和Z2中心坐标处，两条测线都显示中部低极化率异常，Z1左侧高极化率异常，Z2右侧的高极化率异常.随着水位上升，各水位对应的极化率曲线极大值和极小值差值逐渐减小.这是由于含水层本省含水，当陷落柱从不含水到含水时，极化率升高，同时水的极化率值不会无限升高，所以高极化率异常和低极化率异常区别越来越不明显.与数值模拟钻孔测线所显示规律一致，随着水位升高，陷落柱异常体极化率随之升高.

图17 斜向钻孔测线Z1和Z2极化率曲线Fig.17 Polarizability curves of oblique borehole survey lines Z1 and Z2

4 结论

(1) 当水位上升，异常体响应越来越明显，极化率参数的极大值增大，通过极大值点的位置能确定陷落柱异常体在采空区得平面位置.两条钻孔内极化率曲线异常极大值位置能圈定异常体位置，同时能反应水位的变化情况.当利用极化率参数对底板突水情况进行监测的时候，根据不同的测线所观测异常叠加的方法来圈定异常体的位置，并依据异常的值的大小变化监测其水位变化.

(2)在物理模型实验平台采用三维直流电法观测方式，以三维极化率立体图作为对比.物理模型监测结果与数值模拟结果具有很高的一致性，结果表明极化率参数对采空区底板突水监测技术的有效行和可行性.