随着我国经济的快速发展,浅部的煤炭资源已基本开采殆尽,因此,深部找矿是当前能源开发的重点研究方向。瞬变电磁法是当前寻找深部矿产资源中应用最广的物探方法之一,而地—井TEM作为其中一种装置形式,由于其接收探头在钻孔或井下巷道中,靠近目标异常体,具有电磁干扰小、有用信号强、纵横向分辨率高等优点,而使国内外学者对该方法进行了越来越深入的研究。其中P.A.Eaton和G.W.Hohmann^[1],S.H.Ward和R.C.West^[2]等人用积分方程法研究了长导线源激发的地—井瞬变电磁场异常响应特征;孟庆鑫^[3,4]运用有限差分方法研究了均匀半空间中板状导体的地—井TEM响应;徐正玉^[5,6,7]对覆盖层下含有低阻体的垂直断层面进行了地—井TEM正演数值模拟。

本文运用二维时域有限差分算法,建立地球物理模型,模拟出地面—井筒观测方式和地面—巷道观测方式下的板状导体的地—井瞬变电磁响应,并通过绘制多测道曲线,定量求得异常构造所处的位置,为地—井瞬变电磁深部找矿提供可靠的理论依据。

设大地为线性、各向同性、非磁性的导电媒质,在准静态近似下无源Maxwell 方程为^[8]

(1)

式中:E(r, t)、H(r, t)、μ、σ分别为电场强度、磁场强度、磁导率和电导率,r为空间矢量,t为电磁波扩散时间。对式(1)两边取旋度后利用矢量恒等式▽×▽×A=▽▽·A-▽²A,可导出电场的扩散方程:

(2)

求得电场后,根据电场和磁场的关系,可以得到磁场垂直分量:

(3)∂bz(t)∂t=∂ex(t)∂y-∂ey(t)∂x;

通过推导便可得到DuFort-Frankel差分方程^[9,10,11]:

(4)Ei,jn+1=1-4r̅i,j1+4r̅i,jEi,jn-1+2ri,jz1+4r̅i,j·ΔzjΔz̅jEi,j+1n+Δzj+1Δz̅jEi,j-1n+2ri,jx1+4r̅i,jΔxiΔx̅iEi+1,jn+Δxi+1Δx̅iEi-1,jn,Δz̅j=Δzj+1+Δzj2,Δx̅i=Δxi+1+Δxi2,ri,jz=Δtμσ̅i,jΔzjΔzj+1,ri,jx=Δtμσ̅i,jΔxiΔxi+1,r̅i,j=ri,jx+ri,jz2。

初始时间步由网格扩散时间Δt= 14μσ_minδmin2给出,以后随着迭代次数的增加逐步加大。式中σ_min是最小电导率,δ_min为最小网格边长。

具体网格差分格式见图1。

本文采用Mur吸收边界条件,推导过程见文献[12],以侧边x=x_max为例,其表达式为:

∂∂x-1c∂∂tE|x=max=0。

经推导可得侧边上的节点差分格式。

仍以x=x_max界为例:

(5)Emaxn+1=Emax-1,kn+cΔt-δcΔt+δ[Emax-1,kn+1-Emax,kn],

式中:δ为网格节点间距;c=1/μ0ε0;Δt为时间步进增量,一般取电磁波穿过1/2个网格间距的时间即可保持步进稳定,为保证精度,本文选取的Δt更小。

时域有限差分算法中,源可利用均质半空间的解析式作为初始条件代入,前提是异常体要有一定的埋深,以满足均质半空间的限制。源的形式可根据需要选择,此次研究选取回线源作为激发源。

回线源均匀半空间介质的解析式推导见文献[8],文中仅列出y走向所激发的电场解析式的最后结果:

(6)Ey=Iπσr2{2θ2z2e-θ2r2+x2-z2r2erfc(θz)-e-θ2r2-2πθze-θ2r21-2θx1+1θ2r2F(θx)};θ=μσ4t,r2=x2+z2。

式中F(θx)为Dawson积分。

经过matlab编程,得到如图2所示的不同模型在xOz平面观测到的垂直磁场对时间的导数∂b_z/∂t的等值线图,其中x、z分别为水平位置和埋深。图2a表示在纵向方向受到2个板状低阻体的影响,一个低阻体埋深在55~60 m,另一个低阻体埋深在100~105 m位置(黑线所示),磁场的等值线在低阻体处发生明显变化;图2b表示在横向方向受到2个低阻体影响,一个低阻体水平位置在350~400 m,埋深在80~100 m,另一个低阻体水平位置在600~650 m, 埋深在80~100 m位置(黑色方框所示), 磁场的等值线在低阻体处发生变化。说明瞬变电磁法在横向和纵向都对低阻体有较好的分辨力。

本文选取两种观测方式对井地瞬变电磁响应特征进行研究:地面—井筒观测方式和地面—巷道观测方式。

所谓地面—井筒观测方式指在地面发射,在井孔或钻孔内接收电磁场的瞬变响应。

图3是地面—井筒观测方式示意。发射线框400 m×400 m,发射电流1 A。大地电阻率为100 Ω·m。井孔位于大地介质剖面中心(x=0)处,测量深度0~500 m,点距10 m。

图4是地—井观测方式下的TEM响应曲线。选取6个不同时刻的响应组成多测道曲线。

图4a为在z=200 m的深度,含有一个电阻率为10 Ω·m、厚0.2 m的单个板状低阻体的地电模型。在早期,响应幅值较大,随着时间的推移,响应幅值逐渐减弱。极大值出现在涡流场感应最大的地方,本模型在z=200 m深度有一板状低阻体,故在z=200 m处出现极大值。并且,局部导体本身的异常越大(包括异常体大小与电性),电磁场响应的异常幅值就越大。

图4b为在z=200 m和z=300 m的深度,含两个电阻率都为10 Ω·m、厚0.2 m的组合板状低阻体的地电模型。多测道曲线与图4a具有相同的特征,并且出现两个极值,极值点的位置分别在z=200 m和z=300 m附近,这与模型中两个低阻薄板的深度相符。故认为地—井观测方式对低阻体的埋深能很好地分辨出来。

所谓地面—巷道观测方式,指在地面发射,在地下巷道中接收电磁场的瞬变响应。为了便于对地下深部的小信号瞬变响应进行研究,笔者借鉴Greenfield等^[13]在对构造扰动煤层内电磁信号的衰减特性进行研究的过程中,引入参数E_a(单位:dB)来标定构造扰动煤层和完整煤层中电磁信号的幅值差。沿用参数E_a的标定形式,但将E_a重新定义为

Ea(x)=Hz(x)-Hz0(x),

式中:H_z(x)表示存在低阻异常体时所测得的垂直磁场值; Hz0(x)表示在其他参数一样的情况下,去除低阻异常体影响所测得的垂直磁场值,即E_a可看作低阻异常体的“纯响应”。

图5是地面—巷道观测方式示意。发射线框位于地表,大小为400 m×400 m,发射电流 1 A。大地电阻率为100 Ω·m。在z=400 m的深度,x=100 m处,有一个0.5 m×0.5 m大小的低阻异常体,电阻率为10 Ω·m。布置如5所示的3条测线,测线z方向分别在z=0 m、z=300 m和z=500 m的位置,测线x方向都是从-200~200 m。测点点距10 m。

图6是地面—巷道观测方式下的TEM多测道响应曲线。x轴是测点在x方向的位置。从磁场的多测道响应曲线上几乎看不出有低阻异常体的反应。

图7为测线1(z=0 m)的垂直磁场E_a多测道响应曲线,在x=100 m处相交,正好与低阻异常体所在的横向位置相符。图8为测线2(z=300 m)的垂直磁场E_a多测道响应曲线,在x=100 m处相交,正好与低阻异常体所在的横向位置相符。图9为测线3(z=500 m)的垂直磁场E_a多测道响应曲线,曲线在x=100 m处相交,正好与低阻异常体所在的横向位置相符。

与图6相比,图7、8、9都能明显地区分出低阻异常的横向位置。究其原因,低阻体对电磁场的“吸引”作用使得响应在此处不断增大,并随时间的推移,后时刻的“纯响应”超过前时刻的“纯响应”,故E_a曲线在此处发生交叉,正好位于异常体附近。

与图7相比,图8、图9的E_a值提高了2个数量级,说明越靠近低阻异常体,异常响应越大;同时也说明地—井TEM相对于常规的地面TEM,由于接收探头靠近目标异常体,而接收信号具有更强的异常响应。

为了进一步研究地面—巷道观测方式的横向分辨率,重新建立如图10所示的地电模型。模型参数与图5一致,发射线框位于地表,大小为400 m×400 m,发射电流1 A。大地电阻率为100 Ω·m。并在x方向x=-100 m处,增加一个与x=+100 m处对称的低阻异常体,异常体大小0.5 m×0.5 m,电阻率10 Ω·m。同样布置三条测线,测线深度分别在z=0 m、z=300 m和z=500 m的位置,x方向从-200~200 m,测点点距为10 m。同样选取t₁=0.15 ms、t₂=0.2 ms、t₃=0.25 ms、t₄=0.3 ms、t₅=0.35 ms、t₆=0.4 ms 6个时刻的数据进行分析。

图11为测线1(z=0 m)的垂直磁场EA多测道响应曲线,曲线有的不相交,有的相交,且交叉点位置不明确。故难以区分出2个低阻异常体的具体横向位置。

图12、图13分别为测线2(z=300 m)、测线3(z=500 m)的垂直磁场E_a多测道响应曲线,图中曲线均分别在x=-100 m和x=+100 m处相交,相交点清晰明确,并且正好是两个低阻异常体所在的横向位置。

与图11相比,图12、图13不仅E_a响应值较大,并且交叉点清晰明确,对异常体横向位置能准确区分。相对于常规的地面TEM,地—井TEM异常响应大,大大提高了对异常体的分辨率。

1) FDTD法是模拟分析瞬变电磁场特征规律的有效方法,利用该方法可以准确地模拟出不同地电模型的瞬变电磁响应特征。

2) 由于地—井瞬变电磁接收探头更接近目标体,相比于常规的地面瞬变电磁,其响应信号更明显,探测分辨率更高。

3) 当有多个异常体同时存在时,地—井瞬变电磁能够区分出不同的异常体存在,并且能够对不同异常体的埋深和横向位置准确定位。该研究结果为地—井瞬变电磁法定量解释异常体埋深和位置提供了参考依据。

4) 文中对于地—井瞬变电磁大深度小异常的定位问题,需要借助E_a多测道曲线来确定,而对异常体定位的直接算法还不明确,今后作者将对此作进一步研究。