采空区全张量微重力异常响应特征分析

贾豫葛，刘志新

(中国矿业大学 资源学院，江苏 徐州 221116)

1 引 言

煤矿采空区是指由于煤炭资源的开采而造成的地下空腔，这些空腔充满气体、水体或者淤泥等。经过多年地下煤炭资源的开采，地下采空区空间越来越大，随着开采过程的进行，有些巷道与峒室将被破坏，而另一部分则随着矿井的废弃而遗留在地层深处，成为永久性的废弃地下空间。这些废弃地下空间为后续的安全生产和生活都留下安全隐患，而且还会引发地面沉陷和地下水污染等一系列环境问题。

经过一段时间的失稳变形之后，煤矿采空区上部地质围岩在上覆荷载自重应力作用下发生沉陷和塌落，周边围岩破碎，出现应力撕裂裂缝。一般情况下，采空区的塌陷在垂直方向上可分为三带，如图1所示：①冒落带，煤层采空上部岩层出现坍落；②裂隙带，坍落带上方岩体因弯曲变形过大，在采空区上方产生较大的拉应力，两侧受到较大的剪应力，因而岩体出现大量裂隙，岩石的整体性受到破坏；③弯曲带，裂隙带以上直到地面，在自重应力作用下产生弯曲变形而不再破裂。

图1 煤层采空区塌陷垂直“三带”示意图

查明采空区的位置、范围和深度等,需要应用地球物理勘察方法。采空区形成后随着时间的推移，可能出现积水、充填等各种状态，使采空区周围的地质构造出现变化，从而引起采空区的地球物理特征的变化。刘菁华等[1]和程建远等[2]探讨了煤矿采空区及采空区塌陷的地球物理特征及地球物理探查方法。采空区及采空区塌陷的物理性质(如密度、磁性、电性、弹性等)与围岩相比具有较大差异，采空区具有进行地球物理探查的条件。煤矿采空区探查常用的方法有重力勘察、磁法勘察、地震勘探、电法勘探等，都能够有效地探查煤矿采空区的位置、范围和深度等，电法对赋水采空区响应尤其敏感。

本文基于煤矿采空区塌陷垂直“三带”模型，计算采空区或采空区塌陷在地表产生的全张量重力异常，为高精度微重力探查采空区或采空区塌陷提供可靠的理论依据。

2 地质模型

当煤层被采空时，短期内形成一定规模的充气空间，造成采空区相应地层的密度与围岩不同，采空区质量亏损，在地面会产生相应的负重力异常，根据重力异常的分布特征，可以推断采空区的大小、空间分布和深度等问题。经过一段时间后，采空区发生塌陷变形，塌陷垂直“三带”与围岩的密度差异，会在地表产生负重力异常。

弯曲带只产生弯曲沉降，并不产生破裂，带内的岩层仍为层状结构，因此引起的物性差异普遍比冒落带以及裂隙带微弱得多，可以忽略弯曲带产生的重力异常。本文对采空区塌陷的模拟主要是冒落带和裂隙带的重力异常响应。采空区塌陷不充水时，亏损的质量不变，其中裂隙带的形状类似马鞍形，密度小于围岩。冒落带脱离母岩，杂乱无章，填充采空区。采空区塌陷充水时，亏损的质量得到一定的补偿。

采空区或采空区塌陷，不同情况下，理论上都可以探测到负的布格剩余重力异常，采用高精度微重力测量，可以对采空区或采空区塌陷进行探测。

为了研究采空区及采空区塌陷的地球物理场响应，本文建立采空区和采空区塌陷的地质地球物理模型[3，4]。将采空区和采空区塌陷的简化地质背景模型设为三层，从浅到深分别为砂岩层、煤层和灰岩层。采空区简化地质模型剖面图和平面图如图2和图3所示，采空区塌陷简化地质模型剖面图和平面图如图4和图5所示。

图2 采空区模型剖面示意图

图3 采空区模型平面示意图

图4 采空区塌陷模型剖面示意图

图5 采空区塌陷模型平面示意图

采空区和采空区塌陷简化地质模型参数如表1所示。模拟采空区时，采空区大小为长200 m，宽200 m，高10 m，采空区剩余密度为-2 086 kg/m3。

表1 采空区和采空区塌陷模型参数

3 重力异常和全张量重力异常

3.1 引力加速度三分量

本文采用网格化模型空间进行数值计算，将模型空间进行三维网格划分，如图6所示，假定第q个网格质量单元dm质心的位置(x0,y0,z0)，它在地表任一计算点(x,y,z)对单位质量产生的引力加速度在x,y,z三个方向上的分量为[5，6]：

图6 网格质量引力加速度三分量示意图

(1)

(2)

(3)

3.2 全张量重力异常

图7 基于差分的全张量重力异常示意图

(4)

全张量重力异常梯度是一个对称张量，Gyx=Gxx，Gzx=Gxz，Gzy=Gyz。因位场V(x,y,z)满足Laplace′s方程，则有Gxx+Gyy+Gzz=0。

因此，只需计算或测量全张量重力异常梯度G的对角线上两个上角(或下角)的三个分量即可。全张量各元素反映了重力异常的空间变化率。

4 采空区和采空区塌陷重力异常和全张量重力异常模拟

4.1 采空区

采空区埋深380 m，剩余质量Δm=Δρ·ΔV=-8.344×107kg。重力异常即为异常体产生的引力加速度的垂直分量，即Δg=Gz，如图8(a)所示，采空区产生重力异常等值线似为同心圆，因采空区埋深较大，异常体的形态效应不明显[11，12]。重力异常的单位为μGal(微伽)，重力异常极大值约为Δgmax≈-36.98 μGal，其中1 μGal=10-8m/s2。

图8 重力异常平面图

图9 全张量重力异常梯度

张量元Gxx的0值大致对应x方向的边界，张量元Gyy的0值大致对应y方向的边界，张量元Gxy突出增强对角方向的信息。张量元Gxz的极大值和极小值大致对应x方向的边界，张量元Gyz的极大值和极小值大致对应y方向的边界，张量元Gzz突出增强聚焦异常中心的信号强度，与重力异常平面图相比，异常范围变小，便于识别异常体的边界，极大值约为Δgz≈-0.184 2 μGal/m。

4.2 冒落带均匀

模拟采空区塌陷时，不考虑地面沉降，冒落带大小约为长200 m，宽200 m，高50 m，冒落带均匀时，缺失质量不变，密度为1 993.1 kg/m3，则冒落带横向与煤层对应区域，剩余密度为-92.9 kg/m3，横向与砂岩对应区域，剩余密度为-296.9 kg/m3。裂隙带剩余密度为-90 kg/m3。

采空区塌陷，质量亏损不变，缺失的质量分布上移且分布范围变大，比塌陷前采空区产生更大的负的局部重力异常。采空区塌陷冒落带均匀时，重力异常如图8(b)所示，重力异常的单位为μGal，重力异常极大值约为Δgmax≈-48.6 μGal。

采空区塌陷模型产生的全张量重力异常，如图10所示，单位都为μGal/m。图10(a)、图10(b)和图10(c)分别为Gxx、Gxy和Gyy，图10(d)、图10(e)和图10(f)分别为Gxz、Gyz、和Gzz。与采空区全张量重力异常的形态也都相似，重力异常垂向一阶导数极大值约为Δgz≈-0.266 μGal/m。

图10 全张量重力异常梯度

4.3 冒落带密度随深度增大之一

当采空区塌陷的冒落带经过一段时间压实，密度分布不均匀，密度随深度连续增大，变化率为每10 m增加20 kg/m3，参见表2第3、4列所示。冒落带的顶部质量缺失，冒落带的下部质量增加，相当于质量分布下移，即顶部缺失质量距地表更近，影响更大，比冒落带均匀产生更大的负局部重力异常。

表2 冒落带不均匀时的参数

如图11(a)所示，采空区塌陷产生更大的负重力异常，重力异常的单位为μGal，重力异常极大值约为Δgmax≈-48.75 μGal，比图8(b)均匀冒落带情况下的重力异常极大值减小了约0.15 μGam，冒落带纵向密度变化不大时，对重力异常值的影响较小。

图11 重力异常平面

采空区塌陷模型产生的全张量重力异常，如图12所示，单位都为μGal/m。图12(a)、图12(b)和图12(c)分别为Gxx、Gxy、和Gyy，图12(d)、图12(e)和图12(f)分别为Gxz、Gyz、和Gzz。与采空区全张量重力异常的形态也都相似，重力异常垂向一阶导数极大值约为Δgz≈-0.267 2 μGal/m。

图12 全张量重力异常梯度

4.4 冒落带密度随深度增大之二

冒落带随深度密度增大，变化率为每10 m增50 kg/m3，参见表3所示5、6列，原理分析同4.2节。如图11(b)所示，此种情况下，采空区塌陷产生更大的负重力异常，重力异常的单位为μGal，重力异常极大值约为Δgmax≈-50.98 μGal，比图11(a)重力异常极大值减小了约2 μGal。总的来说，冒落带纵向上的密度变化对重力异常值影响不大。

采空区塌陷模型产生的全张量重力异常，如图13所示，单位都为μGal/m。图13(a)、图13(b)和图13(c)分别为Gxx、Gxy、和Gyy，图13(d)、图13(e)和图13(f)分别为Gxz、Gyz、和Gzz。与采空区全张量重力异常的形态也都相似，重力异常垂向一阶导数极大值约为Δgz≈-0.292 μGal/m。

图13 全张量重力异常梯度

5 结 论

通过模型试验，本文模拟计算几种采空区和采空区塌陷模型的全张量重力异常响应，比较分析结果，得出以下结论。

1)煤矿采空区和采空区塌陷造成的局部质量亏损，会在地表产生明显的重力异常和全张量异常响应，异常响应的强度和分布范围与采空区和采空区塌陷的埋深、规模、形状等有关。

2)本文根据前人的研究成果，建立采空区和采空区塌陷简化模型，据此模拟计算了采空区和采空区塌陷重力异常和全张量重力异常响应，并分析采空区和采空区塌陷重力异常和全张量重力异常分布特征。

3)本文数值模拟的结果表明，对煤矿采空区和采空区塌陷进行重力勘探排查可以为识别采空区和采空区塌陷提供可靠的依据，对老矿区地面的安全生产和生活也很有指导意义。

4)采空区塌陷冒落带的纵向密度变化，对重力异常的影响较小，如果采空区埋深较浅，冒落带的压实对重力异常的影响会增大。

5)如果考虑地面沉降，缺失的质量上移，产生的负重力异常，其幅值会更大。