矿井高密度三维电法成像数值模拟及应用

万豪豪，朱 鲁，陈金昊，李文喜

（山东科技大学 地球科学与工程学院，山东 青岛 266590）

直流高密度电阻率法在矿井中的应用已渐为成熟[1-2]，但随着煤矿开采深度的逐渐增加，地质条件变得尤为复杂，其中岩体突水已成为危害矿井安全的重大灾害之一，为确保井下开采工作的安全高效进行，科研工作者进行了深入研究。2003年，底青云分别运用佐迪法和2D积分法进行了高密度电阻率成像方面的研究[3]。2005年，为解决工作面采空区滞后突水问题，翟培合开发了底板突水动态监测系统，实现了工作面底板突水监测过程[4]；2008年，施龙青等运用三维高密度电法探测了煤矿巷道底板富水性[5]；刘志新、徐志刚和王大庆等[6-7]介绍了电阻率成像技术的特殊问题并对电阻率成像技术进行了展望；2016年，高卫富等使用ANSYS进行全空间直流电法异常体正演模拟，精准地探测出了异常体的位置[8]。

1 地电场三维有限元理论

有限元理论的数学基础是变分原理和剖分插值，将研究区剖分成有限个小单元，在每个单元上，将函数简化成线性函数、二次函数或高次函数。它适用于物理性质及参数复杂区域，其求解过程规范化，适于解决复杂数学模型问题[9]。

在直流电法数值模拟时，由于点电源周围电场会因存在诸多干扰因素使得无法精确计算，采用二次场方法解决数据的误差问题，使用正方形的网格对其进行划分，每个网格节点的潜在值被设置为1个未知数，未知数被代入拉普拉斯方程进行计算，即将模型的电导率σ细分为背景电导率σp和异常电导率σs，其分别产生一次场up和二次场us，则总场u为二者叠加：

三维电场电位的边值问题与下列变分问题等价：

式中：σ为地下介质的电导率；r为点电源到边界点的距离；Г∞为区域Ω的下边界；Г为2种介质的分界面；n为边界外法线方向的向量坐标；u0为正常电位。

2 模型建立与完善

在矿井中进行三维直流电法勘探时，根据矿井环境的特点，井下巷道之间大多数都是通过联络巷连接，但是每条联络巷并不都是在1条直线上，无法沿用地面勘测的传统布线方法，且巷道与联络巷两帮的围岩均可导电，使得顶板与底板之间有着电力联系，故并没有采用传统模拟实验中，将巷道、联络巷设定为类空气的超高阻体，而采用井下全空间地电模型[10]，将巷道与联络巷相结合形成测量通道。

巷道设置于地表深度600 m处，巷道长度设置为300 m，位于x坐标-100～200 m处，联络巷长度设置为100 m，位于y坐标-50～50 m处，总长为600 m，模型所要测量的每条边长250 m，宽100 m，无穷远边界长2 400 m，宽2 400 m，深1 500 m，地质模型位置平面图如图1。正方形GHKM表示无穷远边界，ABCD表示模型所要测量的范围，ABEF表示巷道和联络巷的范围。

对所建立的模型进行网格剖分[11]，沿着x轴剖分的平面数110个，沿着y轴剖分的平面数为110个，沿着z轴剖分的平面数为100个，三维立体空间总共剖分区块数约为1.2×106个。为了提高正演计算的精度，在网格剖分过程中对所研究的目标区域进行加密剖分，巷道，联络巷网格单元剖分大小为5 m，研究区域到无穷远边界处网格剖分逐渐变大，网格剖分大小不均匀。

图1 地质模型位置平面图

供电电极A与测量电极M分布在ABCD的巷道上，电极之间的距离选择10 m，电极总数为61个，电极编排的顺序为1～61，电极测量的顺序为CBAD，采用固定断面连续测量扫描方式，测量深度为A到M距离的一半，测量点O位于A到M距离的中点，二级装置AM[12]在三维空间中测量时，测量点O的位置如图2。二极装置A-M三维测量散点图如图2。

图2 二极装置A-M三维测量散点图

根据实际矿井突水隐患建立在层状介质中的地电模型，巷道和联络巷的电阻率为10 000 Ω·m，围岩的电阻率为2 000 Ω·m。模型1设置块状低阻异常体尺寸为20 m×20 m×20 m，电阻率值为10 Ω·m，位于2条巷道之间的煤层下方，顶部距离煤层底板20 m；模型2设置低阻陷落柱异常体模型尺寸为80 m×80 m×20 m，电阻率值为 10 Ω·m，位于 2条巷道之间的煤层下方，顶部距离煤层底板30 m；模型3在模型2的基础上改变陷落柱在底板的位置，其余设置同模型1；模型4设置正方体低阻异常体尺寸为 40 m×40 m×40 m，电阻率值为 1 Ω·m，位于 2条巷道之间的煤层下方，顶部距离煤层底板20 m。

3 数值模拟结果与分析

3.1 块状低阻异常体地电模型

块状低阻异常体地电模型经过有限元数值模拟计算后在经过反演得到xz切面图（图3）。相比于正常地层明显出现了低阻异常区，视电阻率值在10～30 Ω·m，异常体的位置和大小已经在图中用正方形线圈标出，对比地电模型中异常体的位置发生了偏移，大小形态发生了改变，但总体来说二极装置（A M）对于本次的块状低阻异常体的捕捉能力较强，分辨低阻异常体的能力较高。

图3 块状低阻异常体地电模型xz切面图（y=0）

块状低阻异常体地电模型经模拟计算后的三维切片图如图4。其中y=50 m和y=-50 m的xz切片在位于巷道底板处，有2条相互平行的条带，且条带所显示的视电阻率较高，推测是由于巷道高阻区的影响所致，无异常体区域模拟结果良好，巷道底板下方绿色区域的视电阻率值在100 Ω·m左右，与地电模型相符。但异常体x的坐标范围与实际存在偏差，因此，二极装置（AM）在对异常体测量过程中，在y方向上测量效果好，分辨率高。

图4 块状低阻异常体地电模型xz方向三维切片图

3.2 陷落柱低阻异常体地电模型

通过陷落柱低阻异常体地电模型xz切面图（图5），二极装置对本次分辨陷落柱低阻异常体的能力较强，断面图可以明显地反映低阻异常区的大小以及埋深，而且也能够判定异常体的形态，异常体的位置和大小已经在图中用长方形线圈标出，视电阻率值在10～30 Ω·m，其顶底埋深与实际模型基本吻合。

图5 陷落柱低阻异常体地电模型xz切面图（y=0）

由陷落柱低阻异常体地电模型xz方向三维切片图（图6）可知，异常体视电阻率值在100 Ω·m左右。从低阻区的视电阻率范围可推测异常体的y坐标范围为0～15 m，与实际地电模型中y坐标的范围0～20 m基本一致，异常体x的坐标范围与实际模型一致。因此，二极装置（A-M）在对异常体测量过程中，在y方向上测量效果好，分辨率高。中间切片的地层，除了相对低阻区外，与正常地层和块状异常体的视电阻率基本相同。

图6 陷落柱低阻异常体地电模型xz方向三维切片图

3.3 柱状低阻异常体地电模型

柱状低阻异常区的大小和埋深与块状低阻异常体地电模型保持一致，异常体的长度基本和地电模型中的一样，其顶底埋深与实际模型基本吻合。柱状低阻异常体地电模型xz切面图如图7。但低阻体的左边有明显的上移，视电阻率值在20～50 Ω·m，其相较于陷落柱低阻异常体视电阻率值偏大，出现这种现象的原因有可能是异常体的中心坐标发生偏移，离巷道较近，受到巷道高阻体的影响。

从柱状低阻异常体地电模型xz方向三维切片图（图8）观察可知，低阻区的视电阻率范围可推测异常体的y坐标范围为-25～15 m，与实际地电模型中y坐标的范围0～20 m基本一致。

图7 柱状低阻异常体地电模型xz切面图（y=0）

图8 柱状低阻异常体地电模型xz方向三维切片图

3.4 正方体低阻异常体地电模型

正方体低阻异常体地电模型xz切面图如图9。正方体低阻体模拟结果可以反映低阻异常区的大小和埋深，视电阻率值在1～30 Ω·m，图中的空白区域视电阻率值小于10 Ω·m。受低阻异常体的影响，主要异常区域集中在异常体的上部，地面下660 m处是异常区域的分界面，但660 m下方异常反应在实际成像中并没有出现，主要是受探测深度的影响，探测的越深，数据量越少，分辨率越低，因此，在700 m以下的位置视电阻率值并不准确。

图9 正方体低阻异常体地电模型xz切面图（y=0）

从正方体低阻异常体地电模型xz方向三维切片图（图10）可知，从低阻区的视电阻率范围可推测异常体的y坐标范围为-25～25 m，与实际地电模型中y坐标的范围-20～20 m基本一致，在y=0方向切片图中低阻异常区的范围大于实际地电模型，这可能是由于低阻异常体的背景电阻率值非常小，电阻率值为1 Ω·m所致。

图10 正方体低阻异常体地电模型xz方向三维切片图

4 工程应用

某矿1302上工作面（3煤）位于山西组中下部，3煤层开采过程中深受顶板砂岩水害的威胁，为保证1302上工作面开采的顺利、安全进行，对该矿的1302上工作面进行井下高密度三维电法勘探。要求顶板最大探测高度为顶板上150 m，构建顶板砂岩地层电阻率三维数据体，应用切片技术圈定顶板砂岩地层的富水区域，结合具体钻探情况指导1302上工作面开采过程中的防治水工作。

本次井下高密度三维电法勘探数据采集使用的是WDJD-4型高密度电阻率系统，采用A-M二极装置形式的三维电法采集方法，采集参数中仪器所使用的电极道数是60道，测量的剖面数为30个，采用滚动的测量方式，2个相邻电极之间的道间距为10 m，供电电压选择180 V的直流电压。

把RES3DINV反演软件所得出的三维立体成像图数据输出，然后在用Slicer Dicer和Voxler对三维立体成像图进行多方向的三维切片，即为最终的成像成果图。1302上工作面顶板砂岩地层电阻率三维数据体原点位于运输巷与1联络巷向交汇处，沿运输巷向里为y、沿1联络巷向下为正方向x，高度方向为垂直地层向上，本次勘探建立的1302上工作面外侧部分顶板砂岩地层电阻率三维数据体长度750 m、宽度140 m、高度280 m。

1302上工作面顶板地层电阻率三维数据体如图11，1302上工作面顶板地层电阻率三维数据体切片技术如图12，1302上工作面顶板地层电阻率三维数据体垂层切片如图13。

根据数据整体（图11）和数据切片（图12）可知，低阻区域明显，对比图13，可圈定3个低阻区进行解释。第1区域位于横向起始测点到280 m，顶板向上至高度180 m；第2区域位于横向320～400 m之间，顶板向上至高度90 m处；第3区域位于横向400～750 m之间，顶板向上至高度120 m处，以上3个区域，低阻情况明显，推测为含水区域。

图11 1302上工作面顶板地层电阻率三维数据体

图12 1302上工作面顶板地层电阻率三维数据体切片技术

图13 1302上工作面顶板地层电阻率三维数据体垂层切片（平行切眼、50 m间距）

5 结论

1）使用Visual Fortran编程软件，对一定深度范围内所建立的三维全空间地电模型进行正演数值模拟，反演后的结果能够准确反应出所建地电模型异常体区域。Surfer软件绘成的ρs视电阻率断面图，在一定的深度范围内，能大体上反映异常体的形态、位置和大小，分辨率高，基本上与地电模型吻合。

2）三维直流电法在煤矿中的应用实例，圈定了3煤层顶板砂岩地层富水区，后期经过钻探验证了富水区的存在，针对性强，极大地提高了钻探的效率，从生产实际角度论证了高密度电法的可行性，为指导矿井水防治补充了新的方法。