2012, Vol. 55
地-井瞬变电磁法是井中瞬变电磁法中应用最广的一种装置形式，其发送回线通常采用铺设于地面的边长由几百米到一千多米的矩形单匝绝缘不接地大回线，供以脉冲电流产生激发电磁场，切断场源后通过布设于钻孔中的探头接收地层中导电地质体由一次场激发而感应产生的二次场(轴向分量)，通过分析二次场特征获取关于目标体的有用信息[1-2].
地-井TEM 装置下，作用于局部导体的一次场可近似看作指向某个方面的均匀场.在其激励下，局部导体中产生感应涡流.一次场消失后，涡流于早期分布于局部导体表面，后因欧姆损耗逐渐衰变(中期)，最终分布状况处于稳定且按指数规律衰减(晚期)[1-3].
通过对局部导体的井中TEM 响应进行数值及物理模拟研究，总结模型参数的特征关系曲线，可为实际工作中的定性定量解释提供依据，前人做了大量相关研究工作.物理模拟方面，D.V.Woods(1975)[4]采用cronePEM 地-井TEM 系统进行了比例模拟实验研究；G.Buselli和S.K.Lee(1996)[5]进行了均质半空间和导电覆盖条件下斜井的地-井TEM 模拟实验研究.数值模拟方面，前人对球体、圆柱体等简单规则形体的响应问题进行了详细阐述[1， 6]；针对板状导体井中TEM 响应的等效模型模拟，早期研究者采用求解薄板状体内本征电流场编制了模拟计算程序，此后研究中采用了“本征电流"等模型模拟板体，如Maxwell4.0 中TEMH 程序，进行板体组合的钻井中三分量正演模拟，“等效涡流"模型利用单个方形涡流环模拟晚期TEM 响应，在研究走向较长的板状体时可用一对平行反向线电流进行等效[1-2， 6]；P.A.Eaton和G.W.Hohmann(1984)[7]，R.C.West和S.H.Ward(1988)[8]研究了导电围岩对二维、三维导体井中TEM 响应的影响.资料解释方面，J.Macnae和G.Staltari(1987)[9]阐述了导电背景中地-井TEM 响应符号变化；A.V.Dyck(1984)[10]系统分析总结了地-井TEM 勘探中资料解释方法.
时域有限差分法(FDTD)是通过有限差分式代替Maxwell时域旋度方程中的微分式，并随时间逐步推进求解空间电磁场的一种数值模拟技术.自被提出以来[11]，经过几十年的发展和各方面的改进而日渐成熟.该方法被引入瞬变电磁场模拟研究后，以其简便直观的特点，得到了广泛的使用.
实际条件下，局部导体(如矿床)往往由多个导体组合在一起，并且受半空间介质和低阻覆盖层的影响，观测到的结果是多种响应因素共同影响的结果，其响应规律较为复杂，本文基于前人针对TEM数值模拟研究的方法成果[12-15]，对该类问题进行研究分析.
2 时域有限差分法基本原理时域有限差分对半空间瞬变电磁场的求解是将偏微分方程用差分方法离散计算出某时刻空间各点的场值，然后随时间推进计算此后各时间点的场值，可分为以下几个基本步骤:
(1) 对求解区域进行网格剖分，将连续问题离散化，用有限节点代替连续空间.
(2) 将微分方程离散，以差分近似代替微分，推导出离散的时间步进差分方程.
(3) 根据给定的初始条件和边界条件按时间推进求解差分方程.
2.1 基本关系式与差分方程在无源、均匀、非磁性且各向同性的线性介质中，忽略位移电流，准静态下近似无源Maxwell微分方程组为
|
(1) |
|
(2) |
|
(3) |
|
(4) |
为得到显式、稳定的步进解，利用矢量恒等式∇×∇×E=∇ ∇·E-∇2E，将(1)式两边取旋度后，将式(2)(3)带入其中并简化，可得y方向电场的扩散方程
|
(5) |
(1) -(5)式中，E为电场强度，H为磁场强度，σ 为电导率，μ 为磁导率.
在直角坐标系下，将求解区域以均匀平面网格进行剖分，使空间中连续电场值离散为各节点处电场值，以任意节点电场值Ei，j为中心与周边四个节点所构成矩形作为平面单元应用(5)式进行面积分.采用二维Green公式，并用中心差商近似替代电场分量对空间的偏导数的方法解出(5)式左端的离散方程，利用对空间电导率的网格剖分与DuFortFrankel法对(5)式右端电场分量对时间的偏导数进行近似替代的方法解出其差分格式，最终综合两端式子可得扩散方程的五点差分格式，详细推导过程见文献[13-14， 16]，本文仅列出最后结果:
|
(6) |
式中En为tn时刻的电场值，
时域有限差分算法中，源可利用均质半空间的解析式作为初始条件带入，前提是异常体要有一定的埋深，以满足均质半空间的限制，源的形式可根据需要选择[2].位于地表的无限长电流源可用作模拟矩形大回线的一个边[16]，本文选取铺设于地表的两条平行无限长线源(沿y轴走向且通电反向，由这样正负两个平行线源的场叠加)模拟等效半空间地电模型下的矩形回线作为激发源[13].线源均匀半空间介质的解析式推导见文献[16]，本文仅列出y走向单线源所激发的电场解析式的最后结果:
|
(7) |
式中
通过y轴走向上的这样一组平行线源即可得x-z剖面上二维回线激发产生的电场分量，相应的初始时刻t0 和t1 等参数可用网格扩散时间与近地表均质网格层数的试算得到，将该等效源转化为初始条件加入差分方程[13]，为保证采样足够和场的扩散性质，初始时刻可取值稍大一些，同时文中涉及到低阻覆盖层，故网格数应适当加密.
2.3 边界条件本文采用Mur吸收边界条件，针对二维模拟，应用其一阶近似式处理侧边和底边截断界的棱角节点较为简单方便，推导过程见文献[17]，以侧边X= Xmax 为例，列其表达式为:
|
(8) |
经推导可得侧边上的节点差分格式，仍以X=Xmax 界为例[17-18]:
|
(9) |
(8) -(9)式中δ为网格节点间距，
根据上述算法模拟半空间的二维矩形回线源瞬变场，以地-井观测方式计算其响应.本文模型参数选取尽量保持一致性:二维大地介质为800m×400m, 地层电导率为0.01S/m, 对整个介质采取均匀网格剖分，钻孔位于大地介质剖面中心点(x=0)处，接收线圈面积归一化，测量深度0~400 m, 测量点距10m；场源选取的两条线源分别位于(x= -100，z=0)和(x=100，z=0)，电流强度20A，以此等效x边走向长度为200m 的单匝矩形回线.
3.1 半空间TEM场的一般响应特征图 1为几种常见情况中地下瞬变电磁场Ey分量在x-z截面的半剖面等值线，图 1(a、c、d、e、f)约为同一时刻(约0.136ms).
|
图 1 瞬变电磁场(Ey分量)半剖面等值线图 (a)均质半空间；(b)均质半空间于0.266ms时刻；(c)含低阻覆盖层的半空间介质；(d)均质大地中含水平低阻板体；(e)均质大地中含垂直低阻板体；(f)均质大地中含垂直高阻板体. Fig. 1 Electric field (Ey) contour in hal--space of homogeneous medium (a) Half--space of uniform medium；(b) IIal--space of uniform medium at 0.266 ms； (c) Uniform earth with overburden； (d) Uniform earth with a horizontal conductor； (e) Uniform earth with a vertical conductor； (f) Uniform earth with a vertical resistor. |
图 1(a、b)为不同时刻均匀半空间介质中TEM场扩散情况:大地感应涡流产生的TEM 场随时间的延迟向下向外扩散传播，场强随时间推移而衰减；图 1c为均质半空间上部覆盖低阻层(虚线为分界面，层厚50m, 电导率0.1S/m)的TEM 场情况:扩散速度相对图 1a较慢，两层介质分界处电场等值线发生凹缩，显示出瞬变场在不同电性介质中传播速度的差异.图 1d为均匀介质中含一个水平低阻板状体(图中虚线所示，大小为200 m×10 m, 电导率1S/m, 深度200 m)情况，图 1e为均匀介质中含一个垂直低阻板状体(虚线所示，大小为10m×200m, 电导率1S/m, 顶部深200 m, x向位置200 m)情况:相比正常场，图 1d和图 1e在低阻导体及附近的电场等值线发生明显畸变，等值线密集，表明电磁波扩散到低阻导体时扩散速度降低，说明涡流在低阻导体中密度大、衰变慢，反应出瞬变场对低阻导体较敏感；图 1f为均匀介质中含一个垂直薄板状高阻体(虚线所示，电导率0.001S/m, 其他参数同于图 1e)情况:在高阻异常体及附近的电场等值线没有发生明显变化，反应出瞬变场对高阻异常体有较好的穿透性，实际计算结果与正常场仍有差别.
图 2为位于发射回线的中心位置的钻孔中接收的响应曲线；由于该情况下接收点水平位置始终与地层感应涡流的中心部位重合，则各测点在所选观测时段内接收响应同号(设为正)；纵坐标为响应幅值，横坐标为测点深度，参数与第3节的同.图 2a为0.136ms时刻均匀半空间介质情况下响应曲线:参见图 1a整个感应电场分布情况，曲线规律大致为浅部测量点响应较大，深部测量点响应较小，所接收的响应最大值出现在涡流场极大值所处深度附近；图 2b为均质半空间情况下不同时刻的接收响应按照最小值(一般是最深处测点所接收的响应)进行归一化所得曲线:随着时间延迟，接收到的响应幅度逐渐变小，说明瞬态场整体分布随时间推移而趋于均匀；同时注意到响应最大值位置在逐步变深，反映出涡流场极大值随时间延迟逐渐向下向外移动(参见图 1(a、b)).图 2c为均匀半空间介质上部覆盖低阻层(层厚50m, 电导率0.1S/m, 约0.136ms时刻)的响应曲线:参见图 1c感应电场分布情况，由于场扩散速度较慢，该时刻涡流场强值区域分布于覆盖层中，加之覆盖层本身电性特征，使得浅层位置接收的响应值较大，反映出地-井TEM 观测方式下“导电覆盖层，在理论上可视为一个水平导电薄层"[1]；图 2d为同样含低阻覆盖层的均质半空间情况下不同时刻的接收响应按照最小值(最深处测点所测响应)进行归一化所得曲线:对比图 2b，虽然随时间延迟其接收到的响应曲线幅度也变小，但其强度衰减速度较慢，故需要更大的延时才能减小其对整个观测值的影响.
|
图 2 均质半空间与含低阻覆盖层介质的地-井TEM 响应曲线(中心钻孔) (a)均质半空间响应；(b)归一化后不同时刻的响应(均质半空间)；(c)含低阻覆盖层半空间介质响应；(d)归一化后不同时刻的响应(含低阻覆盖层的半空间介质). Fig. 2 Surface-hole TEM responses for half-space of homogeneous medium or contain overburden in mid-point dill-hole (a) Response curve of uniform medium； (b) Response curves of normalized (uniform medium) ； (c) Response curve of uniform medium with overburden； (d) Response curves of normalized (uniform medium with overburden). |
图 3(a、b)为包含大小分别为200 m×10 m 和300m×10 m 的单个水平板状导体(电导率均为1S/m, 深度200m)的均质半空间中，0.084~0.657ms时刻(间隔约0.052 ms)共12 个时道沿钻孔(位于中心)所测量的响应曲线，单板体如图 3所示位于发射回线Tx正下方，Rx为接收线圈.
|
图 3 均匀半空间中板状导体的地-井TEM 响应曲线 (a、b)单个板状导体；(c、d)组合板状导体. Fig. 3 Surface-hole TEM responses for uniform medium that contain conductive plates in mid-point dilf-hole (a, b) A conductive plate； (c, d) Multiple conductive plates. |
由图 3(a、b)可知，在较早时道，所接收的响应曲线受导电半空间影响较大，如接收点约在0~150m处，随着时间的延迟，导电半空间对该地电模型总响应的影响逐渐变小，同时局部导体本身的异常响应值越大(与异常体大小和电性有关)，观测曲线受导电半空间的影响就越小；这些都反映出半空间响应对总响应的影响情况:早延时较强，随延时增加逐渐变弱，这种影响主要在涡流场极大值所处深度附近，离场极大值越远的接收位置(越深)受影响越小.需要注意，总响应特征还和激发源与异常体的耦合关系有关，在耦合情况较差或采样时道过早时，导电薄板体(目标体)的响应甚至会湮没在总响应中，难以显示出来[1].
图 3(c、d)为包含组合水平板状导体(大小和电性分别与(a、b)相同)的均匀半空间介质中，12个时道(0.334~0.911ms时刻，间隔约0.052 ms)沿中心钻孔所测量的响应曲线，两个板状导体深度分别为200m 和300m, 如图 3 所示位于发射回线正下方.可知，在较早时道，从总响应中较难分辨出两个单独板体(下部薄板体响应特征容易被湮没，且在更早时道更难分辨)，随着时间推移，逐渐可明显分辨出组合中的单体；在板体间距相同时，响应较大的板体组合(图 3d)需要更大的时间延迟才能分辨出组合中的单体.这是由于组合中板体与激发源的耦合差异以及半空间影响共同造成的.
图 4模型为图 3模型加盖一层厚50m、电导率0.1S/m 的低阻覆盖层，其他参数一致.图 4(a、b)为0.124~0.701 ms时刻(间隔约0.052ms)共12个时道沿钻孔(位于中心)所测量的响应曲线.可知总响应为覆盖层响应与局部导体响应及地层半空间响应的相互影响叠加的结果，较之图 3(a、b)，导电薄板体响应在总响应中更不明显，覆盖层对总响应的影响在所选整个观测时间段内都很明显，且多数测点接收的响应(包括响应值和衰减特征)均受其影响，这反映出由于覆盖层靠近激发源，其层内首先产生感应涡流，随时间扩散和衰减，这种涡流将作用于其下局部体并激发涡流且同样随时间扩散衰减[1， 5]，由于低阻覆盖层影响较之均质半空间影响更大，需要更大的延时才能减小其影响，故在相同条件下需更长的延迟时间才能在响应曲线中分辨出薄板体响应.图 4(c、d)为该覆有低阻层的均匀半空间介质中包含组合水平板状导体的情况，12个时道(0.649~1.226ms时刻，间隔约0.052ms)沿中心钻孔所测的响应曲线，两板体深度分别为200m和300m.可知，由于半空间、覆盖层和两板体响应差异的共同影响，相比图 3(c、d)，需要更大的时间延迟才能完全分辨组合中的两个单体.
|
图 4 含低阻覆盖层的半空间介质中板状导体的地-井TEM 响应曲线 (a、b)单个板状导体；(c、d)组合板状导体. Fig. 4 Surface-hole TEM responses for uniform medium with overburden and conductors in mid-point dil--hole (a, b) A conductive plate； (c, d) Multiple conductive plates. |
本文通过时域有限差分法，利用一对线源(等效二维回线源)为激发源，采用Mur吸收边界条件和均匀网格及时步长，模拟了地-井观测方式下均匀半空间介质和覆有低阻层的介质中导电薄板体的TEM 响应，取得了如下成果和结论:
模拟了地-井观测方式下均匀和含低阻覆盖层的半空间介质的TEM 响应，简述了响应特征规律和成因；基于此模拟了这两种地电模型中薄板导体和板状体组合的响应情况，分析了其特征规律，认识到选取较大的延时(晚延时)可有效的降低半空间响应产生的影响以显示出异常体响应特征，且能够更好分辨出组合板体中的单体响应特征；
低阻覆盖层对总响应的影响较均质半空间更大，且在相当大的观测时段内都存在，使得异常体响应特征在总响应中更不明显，在相同观测条件下需要相对更大的延时才能有效的分辨出异常体响应和导电薄板体组合中的各单体响应.
地-井TEM 响应特征规律复杂，还存在很多不足和问题:文中笔者仅考虑了最简单观测条件(垂直钻孔位于回线中心，异常体位于场源正下方，测量时道选取较早)下的响应特征和规律，参数选择较为单一，很多包括如倾斜钻孔观测、激发源和观测钻孔位置变化(如“定源动接收"、“动源定接收")而引起所接收异常体响应和覆盖层响应的“变号"情况及其特征规律、总响应的动态特征等问题都需要进一步研究；此外，笔者所选取的场源仅是二维模拟计算，为了得到更好的模拟效果，需选择更适当的场源及初始条件进行三维模拟，以讨论不同方位上的各场值分量的测量；同时为了提高模拟精度，要对边界条件进行改进，这些都是FDTD 模拟中的重要问题.
| [1] | 蒋邦远. 实用近区磁源瞬变电磁法勘探. 北京: 地质出版社, 1998 . Jiang B Y. Applied Near Zone Magnetic Source Transient Electromagnetic Exploration (in Chinese). Beijing: Geological Publishing House, 1998 . |
| [2] | 牛之琏. 时间域电磁法原理. 长沙: 中南大学出版社, 2007 . Niu Z L. ,The Theory of Time Domain Electromagnetic Methods (in Chinese). Changsha: Central South University of Technology Press, 2007 . |
| [3] | 李貅. 瞬变电磁测深的理论与应用. 西安: 陕西科学技术出版社, 2002 . Li X. The Theory and Application of Transient Electromagnetic Sounding (in Chinese). Xi'an: Shanxi Science Technology Press, 2002 . |
| [4] | Woods D V. A model study of the Crone borehole pulse electromagnetic (PEM) system. Unpublished M.Sc. thesis, Queen's University, Ontario, Canada , 1975. |
| [5] | Buselli G, Lee S K. Modelling of drill-hole TEM responses from multiple targets covered by a conductive overburden. Exploration Geophysics , 1996, 27(3): 141-153. DOI:10.1071/EG996141 |
| [6] | 布塞利, 等. 瞬变场法野外工作方法和数据解释手册. 北京: 地质出版社, 1992 . Buselli G, et al. Manual for Sirotem Field Procedures and Data Interpretation (in Chinese). Beijing: Geological Publishing House, 1992 . |
| [7] | Eaton P A, Hohmann G W. The influence of a conductive host on two-dimensional borehole transient electromagnetic responses. Geophysics , 1984, 49(7): 861-869. DOI:10.1190/1.1441732 |
| [8] | West R C, Ward S H. The borehole transient electromagnetic response of a three-dimensional fracture zone in a conductive half-space. Geophysics , 1988, 53(11): 1469-1478. DOI:10.1190/1.1442427 |
| [9] | Macnae J, Staltari G. Classification of sign changes in Borehole TEM decays. Exploration Geophysics , 1987, 18(3): 331-339. DOI:10.1071/EG987331 |
| [10] | Dyck A V, West G F. The role of simple computer models in interpretations of wide-band, drill-hole electromagnetic surveys in mineral exploration. Geophysics , 1984, 49(7): 957-980. DOI:10.1190/1.1441741 |
| [11] | Yee K S. Numerical solution of initial boundary value problems involving maxwell's equations in isotropic media. IEEE Trans. On Antennas and Propagation , 1966, 14(3): 302-307. DOI:10.1109/TAP.1966.1138693 |
| [12] | 宋维琪, 仝兆歧. 3D瞬变电磁场的有限差分正演计算. 石油地球物理勘探 , 2000, 35(6): 751–756. Song W Q, Tong Z Q. Forward finite differential calculation for 3D transient electromagnetic field. Oil Geophysical Prospecting (in Chinese) , 2000, 35(6): 751-756. |
| [13] | 闫述, 陈明生, 傅君眉. 瞬变电磁场的直接时域数值分析. 地球物理学报 , 2002, 45(2): 275–284. Yan S, Chen M S, Fu J M. Direct time-domain numerical analysis of transient electromagnetic fields. Chinese J. Geophys. (in Chinese) , 2002, 45(2): 275-284. |
| [14] | 杨海燕, 岳建华. 磁偶源2.5维瞬变电磁场全空间FDTD数值模拟. 物探与化探 , 2008, 32(3): 326–330. Yang H Y, Yue J H. 2.5-D numerical simulation on the whole-space TEM field with a magnetic dipole source. Geophysical & Geochemical Exploration (in Chinese) , 2008, 32(3): 326-330. |
| [15] | 岳建华, 杨海燕, 胡搏. 矿井瞬变电磁法三维时域有限差分数值模拟. 地球物理学进展 , 2007, 22(6): 1904–1909. Yue J H, Yang H Y, Hu B. 3D finite difference time domain numerical simulation for TEM in-mine. Progress in Geophysics (in Chinese) , 2007, 22(6): 1904-1909. |
| [16] | 纳比吉安. 勘查地球物理电磁法 (第1卷). 北京: 地质出版社, 1992 . Misac N N. Electromagnetic Methods in Applied Geophysics, Volume I, Theory Society of Exploration Geophysicists (in Chinese). Beijing: Geological Publishing House, 1992 . |
| [17] | 葛德彪. 电磁波时域有限差分方法. 西安: 西安电子科技大学出版社, 2002 . Ge D B. Finite-difference Time-domain Electromagnetic Method (in Chinese). Xi'an: Xidian University Press, 2002 . |
| [18] | 王长清. 电磁场计算中的时域有限差分方法. 北京: 北京大 学出版社, 1994 . Wang C Q. Calculation of Electromagnetic Fields Finite-Difference Time-domain Method (in Chinese). Beijing: Beijing University Press, 1994 . |
| [19] | Mur G. Absorbing boundary conditions for the finite-difference approximation of the time-domain electromagnetic-field equations. IEEE Trans. Electromagn. Compat. , 1981, 23(4): 377-382. |
| [20] | Liao Z P, Wong H L, Yang B P, et al. A transmitting boundary for transient wave analyses. Scientia Sinica (Series A) , 1984, 27(10): 1063-1076. |