基于Lorenz规范条件下磁矢势和标势耦合方程的频率域电磁法三维正演

陈辉,殷长春,邓居智

1 吉林大学地球探测科学与技术学院，长春　130026 2 东华理工大学放射性地质与勘探技术国防重点学科实验室，南昌　330013



基于Lorenz规范条件下磁矢势和标势耦合方程的频率域电磁法三维正演

陈辉1,2,殷长春1*,邓居智2

1 吉林大学地球探测科学与技术学院，长春130026 2 东华理工大学放射性地质与勘探技术国防重点学科实验室，南昌330013

为了克服空气层和地表耦合以及避免一次场计算，开发适合不同类型场源、不同应用范围的频率域三维正演模拟统一平台，本文从麦克斯韦基本方程出发，推导基于Lorenz规范条件的磁矢势和标势耦合方程；通过将不同类型场源分解成一系列短导线(电性)源组合，采用交错网格采样和有限体积技术对方程进行离散得到对称大型稀疏线性方程组，并采用Jacobi迭代预处理QMR(Quasi-Minimum-Residual，拟最小残差)算法进行求解，我们成功实现不同类型场源、不同应用范围的频率域电磁法三维正演模拟.通过层状模型下大地电磁法以及有限长接地导线和大回线磁性源激发下的电磁场响应模拟，并与一维解析解对比验证算法的有效性.进而，我们利用该算法平台的模拟结果对典型地电模型在不同场源激发下频率域电磁法响应特征进行对比分析.本文算法研究及实现为建立频率域电磁法三维正反演统一框架打下基础.

频率域电磁法；有限体积法；正演模拟；Lorenz规范；耦合势方程

1　引言

频率域电磁法(Frequency-Domain Electromagnetic,FDEM)通过观测和分析不同频率的人工场源或天然场源激发地下介质产生的二次场或总电磁场分布规律以探明地下电性分布特征.该方法种类繁多，按照场源性质(Nabighian，1988)可分为平面波场的大地电磁法(MT)、音频大地电磁法(AMT)、甚低频电磁法(VLF)；接地导线激发的可控源电磁法(包括地面可控源音频大地电磁CSAMT和海洋可控源电磁法MCSEM)；不接地回线激法的电磁剖面法(包括地面电磁和航空电磁法AEM)等.这些方法在深部地球结构探测 (Selway,2014)、地热勘探(Muoz,2014)、矿产勘查(Smith,2014)、油气勘探(Strack,2014)以及坏境与工程勘探(Sheard et al.,2005)等各领域起到关键作用.

频率域电磁法三维正、反演已成为目前电磁勘探研究的热点和难点问题.其中，MT三维正反演技术已基本成熟，并开始走向实用化(Siripunvaraporn，2012；Miensopust et al.,2013).Avdeev(2005)，Börner(2010)和Newman(2014)对电磁法三维正反演技术发展和挑战进行综述，指出在正演算法中积分方程法(Wannamaker,1991;Song et al.,2002;Farquharson et al.,2006;Zhdanov et al.,2006;Avdeev and Knizhnik,2009;Bakr and Mannseth,2009)通常利用Green函数求解关于二次电场的Fredholm积分式，方法已基本成熟，但对复杂模型的正演精度不高；有限单元法(Yoshimura and Oshiman,2002;Mukherjee and Everett,2011;Pardo et al.,2011;Yahya et al.,2012;Schwarzbach and Haber,2013；Cai et al.,2014)通常利用变分原理或加权余量法求解关于电场、磁场或耦合势的微分方程，能够模拟复杂地形和地下结构，但计算量巨大、求解速度慢；有限差分或有限体积法既能模拟复杂模型，又具有较快计算速度，已成为三维电磁正反演中主流算法.Mackie等(1993)使用有限差分法实现了三维大地电磁正演模拟并将计算结果与积分方程法进行了对比，验证了方法的正确性；Smith(1996a,1996b)对交错采样有限差分法的原理、误差分析以及加快计算速度的双共轭梯度法进行了详细阐述；谭捍东等(2003),Tan等(2006)系统论述了消去电场分量得到关于磁场的大地电磁三维交错采样有限差分数值模拟算法，并对实现过程中交错网格剖分、积分公式离散化、边界条件、方程组求解、三维张量阻抗的计算等内容做出研究，并实现了并行计算；陈辉等(2011)在此基础上提出利用磁场散度实校正方法(RRCM)提高三维正演精度和加快正演速度；沈金松(2003)利用交错网格有限差分法消去磁场分量得到关于电场的线性方程组，实现了三维频率域电磁响应正演模拟，并且在文中提出在迭代过程中施加了散度校正的方法；Weiss和Constable(2006)利用交错网格有限差分求解二次电场的双旋度方程实现海洋电磁法三维数值模拟；邓居智等(2011)利用交错网格有限差分求解二次磁场的双旋度方程实现可控源音频大地电磁模拟；殷长春等(2014)利用交错网格有限差分求解二次电场的双旋度方程实现海洋电磁法各向异性三维正演模拟；Haber等(2000)针对基于电或磁场的双旋度方程在空气电阻率无穷大引起非椭圆方程和地面强烈耦合作用造成求解不准问题，提出利用有限体积法求解Coulomb规范条件下矢势和标势耦合方程，以实现频率域电磁法三维正演，其离散线性方程组是正定的，但非对称.Hou等(2006)和张烨等(2012)利用有限体积法求解Coulomb规范条件下二次矢势和标势耦合方程，实现感应测井三维各向异性模拟；周建美等(2014)利用该方法实现海洋电磁法三维各向异性正演模拟.综上所述，为了解决人工场源源项处理和场源畸变问题，通常采用数值算法求解二次电场和磁场双旋度方程；考虑到空气层和地下介质耦合问题，人们近年提出求解Coulomb规范条件下二次场标势和矢势耦合方程，但该算法具有网格节点数巨大，求解一次场非常耗时及Coulomb规范条件产生的耦合方程非对称等缺点.

本文从麦克斯韦方程组出发，在Lorenz规范条件下建立磁矢势和标势耦合天然对称方程，采用交错采样网格和有限体积法离散磁矢势和标势方程，并将不同场源类型分解成一系列短导线(电性)源组合，对离散线性方程组采用Jacobi迭代的预处理和拟最小残差法(QMR)求解，从而实现任意场源激发下频率域三维电磁正演.最后，我们通过与各种不同类型一维理论解进行对比验证了该算法的有效性，并进一步对比分析不同场源激发下典型地电模型的电磁响应特征.

2　磁矢势和标势耦合方程及边界条件

在电磁勘探中，对大多数地下岩石我们可以忽略位移电流.假定时间变化因子为eiω t，则准静态条件下麦克斯韦方程式可写为

(1a)

(1b)

(1c)

(1d)

式中，E为电场，H为磁场，Je为电流密度，μ为磁导率，ε为介电常数，σ为电导率.

(2)

假设介质磁导率μ和介电常数ε均为常数，将式(2)代入式(1b)可得

(3)

为了保证磁矢势A和标势φ唯一及方程的对称性，我们引入洛伦兹(Lorenz)规范条件：

(4)

其中c=iωμ.将式(4)代入式(3)可得

(5)

对式(3)两边取散度,并代入Lorenz规范条件(4)可得

(7)

可以看出，该耦合方程具有天然对称性，对离散后形成的线性方程组求解精度和稳定性有一定提升.由于地下介质为有损介质，电磁场随传播距离增加呈指数衰减，在正演模拟中可选择足够大的求解区域Ω，在求解区域Ω外的电磁场值非常小.因此可以采用简单的截断边界条件.与其对应的Lorenz规范条件下磁矢势和标势的边界条件可表示为

(8)

3　频率域三维有限体积电磁正演模拟

为了求解磁矢势和标势耦合方程(7)，我们将系统阐述有限体积法离散过程、场源处理及线性方程组合成求解等.

3.1区域离散

首先，将求解区域Ω采用一系列平行于x、y、z坐标轴的三组平面剖分成若干个小的六面体网格单元.假设沿x轴方向被剖分成Nx段，含有l个节点(l=Nx+1)；节点编号沿x轴方向递增i=1,2,…,l，网格间距为dxi(i=1,…,Nx)；沿y轴方向被剖分成Ny段，含有m个节点(m=Ny+1)，节点编号沿y轴方向递增，j=1,2,…,m，网格间距为dyj(j=1,…,Ny)；沿z轴方向被剖分成Nz段，含有n个节点(n=Nz+1)，节点编号沿z轴方向递增，k=1,2,…,n，网格间距为dzk(k=1,…,Nz);共计有Nx×Ny×Nz个长方体剖分单元.

图1　磁矢势和标势交错采样网格Fig.1　Staggered grids of magnetic vectors and scalar potentials

3.2方程离散

为实现有限体积法离散，需要定义磁矢势三个分量和标势体积单元大小.我们以Ax分量为例.参见图2a，体积单元沿x方向边长为dxi，沿y和z方向的边长分别为(dyj-1+dyj)/2和(dzk-1+dzk)/2.因此，该单元体积为

Vi+1/2,j,k=dxi(dzk-1+dzk)(dyj-1+dyj)/4.

(9)

按照等效电阻率原理，该体积单元电导率可以表示为

(10)

图2　磁矢势和标势体积单元示意图Fig.2　Sketch showing cells of magnetic vectors and scalar potentials

同理可得Ay,Az分量体积单元表达式.对于标势φ体积单元(参见图2b)，沿x，y和z方向的边长分别为(dxi-1+dxi)/2，(dyj-1+dyj)/2，(dzk-1+dzk)/2.由此，单元体积为

(11)

而其体积单元等效电导率为

(14)

(15)

(16)

对于磁矢势和标势耦合方程组(2)第二式，在某个标势φ体积单元上积分可得

(18)

仿照前面离散过程，(18)式中左端四项的离散公式为

(19)

(20)

(21)

(22)

3.3场源设置

3.4离散线性方程组合成及求解

将(14)—(17)式代入方程组(7)可以建立磁矢势Ax分量和标势φ的离散表达式(见图4a)，同理将(19)—(22)式代入方程组(7)可建立标势φ和磁矢势Ax、Ay、Az分量的离散表达式(见图4b).将方程组(7)离散后Ax、Ay、Az、φ四个等式，并且结合边界条件(8)式可以合成总体线性方程组：

(23)

式中Cx，Cy，Cz，Cφ及Dx，Dy，Dz为系数项，而rx，ry，rz，rφ为源项，待求未知数个数为

N=(nx-1)×(ny-2)×(nz-2)

+(nx-2)×(ny-1)×(nz-2)

+(nx-2)×(ny-2)×(nz-1)

+(nx-2)×(ny-2)×(nz-2).

(24)

对于线性方程组求解目前主要采用直接求解和迭代求解.迭代求解通常在Krylov子空间内进行.首先采用不完全LU分解、不完全Cholesky分解等预处理方法降低系数矩阵的条件数，然后采用共轭梯度法(CG),双共轭梯度法(BiCG)，广义最小残量法(GMRES),拟最小残差法(QMR)等迭代方法进行求解.该算法需要的内存小，计算速度快，适合求解场源较少的正反演问题.直接求解方法一般直接对系数矩阵进行LU分解，然后利用各种成熟的求解库MUMPS、PARDISO等进行求解.该方法所需内存大、计算速度慢，但特别适合多源问题.本文离散线性方程组(23)的系数矩阵为大型对称稀疏矩阵，由于仅考虑单一发射源问题，我们选用Jacobi迭代的QMR算法进行求解.

图3　频率域电磁法不同类型场源类型Fig.3　Transmitting sources of different types for FDEM

图4　标势和磁矢势离散关系图Fig.4　Template for discretizing magnetic vectors and scalar potentials

4　数值模拟与结果分析

为了验证本文的频率域三维电磁正演算法的可行性和有效性，我们利用本文算法结果和一维正演结果进行对比，研究不同场源激发条件下频率域视电阻率和相位响应特征.

4.1算法验证

选取水平层状三层H模型(如图5a).第一层厚度500 m，电阻率100 Ωm；第二层厚度1000 m，电阻率10 Ωm；第三层电阻率为1000 Ωm.求解区域大小为10000 m×10000 m×10000 m，网格剖分为100×100×100，分别向外延拓16层，延拓指数为1.5.

首先进行平面波场X极化和Y极化模式下的电磁场三维模拟，然后计算大地电磁阻抗视电阻率和相位.计算频率为1～1000 Hz.图5b为水平层状介质三维和一维正演结果响应对比.由图可见三维模拟的视电阻率和相位与一维结果完全吻合，视电阻率相对误差最大0.3%，相位最大相对误差2.9%，说明本文算法对大地电磁法三维模拟具有较高精度.

为验证本文算法的广谱有效性，我们进一步对有限长接地线源的CSAMT进行三维正演模拟.网格参数设置和上面模型一致.假设接地长导线发射源沿x方向置于地面，长度为2000 m，中心点相对投影(0,0,0)，发射频率为10 Hz.接收点设在位于发射源中垂线的正向y轴上，点距间隔为100 m.图6为层状介质CSAMT三维和一维模拟结果对比.从图可以看出，Ex分量的实部和虚部三维模拟结果几乎与一维结果重合，最大相对误差不超过3%；Hz分量实部和虚部也基本与一维正演结果吻合，除了Hz分量实部受场源的畸变影响在发射源附近两个测点相对误差达到10%，其他测点相对误差均小于5%.考虑到CSAMT通常在远区或过渡区进行观测，在场源附近发生的畸变并不会影响数据解释精度.

图5　三层水平层状模型及大地电磁响应结果对比Fig.5　Comparison of 1D and 3D MT responses for a 3-layer model

图6　接地长导线水平层状模型CSAMT一维和三维电磁响应对比Fig.6　Comparison of 1D and 3D CSAMT responses for the 3-layer model of Fig.4

图7　不接地回线三层水平层状模型一维和三维电磁响应对比Fig.7　Comparison of 1D and 3D EM responses for a wire-loop transmitter over the 3-layer model in Fig.4

最后对大回线源激发条件下的三维电磁场进行正演模拟.网格参数设置和上面模型一致.设置发射源为3000 m×3000 m方形回线，置于地面，中心点相对投影(0,0,0)，发射频率为10 Hz.接收点位于发射源中垂线的正向y轴上，点距100 m.图6为层状介质大回线激发下的电磁场三维和一维模拟结果对比.从图中可以看出，Ex分量的实部和虚部三维模拟和一维模拟结果重合，最大相对误差不超过3%；Hz分量实部和虚部也基本与一维正演结果吻合，除了Hz分量实部受场源的畸变影响在靠近发射源1.5 km附近相对误差可达10%，但其余各测点相对误差均小于5%.因此，本文算法对大回线激发下的电磁场三维模拟具有较高精度.

4.2典型地电模型

通过水平三层模型不同类型场源的数值模拟结果与解析解的对比，证明了本文算法适合各类场源的频域电磁法三维正演，并且具有较高精度.下面，我们研究不同类型场源复杂模型的频率域电磁响应特征，特别是进行视电阻率和相位响应结果的对比.如图8，我们设计一个低阻异常体模型大小为800 m×800 m×500 m，顶部埋深为500 m，异常体电阻率为10 Ωm，围岩电阻率为100 Ωm.将10000 m×10000 m×10000 m计算区域剖分为100×100×100个单元，沿x、y、z三个方向的扩展层数为16，空气扩展层为12.将异常体剖分成8×8×5个单元,计算频率为1～1000 Hz.人工源频率域电磁法采用赤道观测装置，有线长接地导线长度为2000 m，磁偶源为600 m×600 m回线，收发距都为6 km.

图8　低阻异常模型Fig.8　Low-resistivity anomaly model in a homogeneous half-space

图9为不同场源激发下低阻模型中心点视电阻率和相位结果对比.从图可以看出，在100 Hz以上四种不同类型场源(XY模式、YX模式、接地长导线、大回线磁偶源)的高频段视电阻率和相位曲线基本重合；随着频率降低XY模式、YX模式大地电磁视电阻率响应受三维异常体的影响出现分离，而相位曲线基本重合；其中XY模式受三维异常体影响要大于YX模式.对于接地长导线、大回线源，随频率降低观测场由远区逐渐进入近区，电阻率呈现快速上升，相位快速下降趋势.其中大回线磁偶源进入近区频率要高于接地长导线源，而且幅度大于长导线源.因此，实际电磁勘探工作中，为便于在远区观测通常采用接地长导线源作为发射源.

图10为低阻异常体模型不同类型场源激发条件下视电阻率和相位拟断面图.由图可以看出，当频率高于100 Hz时，即处于远区条件下，CSAMT、回线源赤道装置的视电阻率响应基本相同，在异常体上方近地表位置有一个不太明显的小范围高阻区，下方为低阻异常区；低阻异常区的分布范围与模型中低阻异常体的实际宽度基本相同；相位特征也具有类似特征，呈现高值异常.随着频率降低，XY模式视电阻率呈现低阻视电阻率异常和高值相位异常；YX模式呈现低阻加两个对称高阻视电阻率异常和高值相位外加两个对称低值异常；然而接地长导线和大回线源的赤道装置受近场影响，呈现高阻假异常和低值相位假异常.其中大回线源进入近场频率要高于接地长导线，而且异常响应范围和幅值都要小.因此开发不同场源的频率域电磁三维反演算法可克服三维异常体和场源引起的虚假异常，有助于提高数据解释质量.

图9　不同场源激发下低阻模型中心点视电阻率和相位结果对比Fig.9　Comparison of apparent resistivities and phases at central point in low-resistivity model for different transmitting sources

图10　不同类型场源低阻异常模型的视电阻率和相位响应拟断面图横坐标为算术坐标，纵坐标频率为对数坐标.Fig.10　Pseudo-cross section of apparent resistivity and phase for the model in Fig.8 for different transmitting sources

5　结论

本文从麦克斯韦方程组出发，在Lorenz规范条件下建立磁矢势和标势耦合对称方程.采用交错网格和有限体积法离散磁矢势和标势耦合方程，并将不同场源类型分解成一系列短导线电性源的组合；对离散线性方程组采用Jacobi迭代的预处理拟最小残差法(QMR)进行求解，我们成功实现任意场源下频率域三维电磁正演模拟.通过对水平层状模型不同类型场源的三维电磁模拟可以发现，本文开发的算法能够实现不同类型场源各种不同应用范围频率域三维电磁正演模拟，避免耗时的一次场计算，计算精度较高.由于采用特殊的源处理技术，使得频率域三维电磁正反演统一数据处理和解释平台成为可能.对不同场源激发下典型地电模型模拟发现，不同类型场源的视电阻率和相位特征和规律不同，在实际勘探中需要考虑场源的选择和设置，同时在电磁数据解释中应考虑场源的畸变效应.

References

Avdeev D B.2005.Three-dimensional electromagnetic modelling and inversion from theory to application.Surveys in Geophysics,26(6):767-799.Avdeev D,Knizhnik S.2009.3D integral equation modeling with a linear dependence on dimensions.Geophysics,74(5):F89-F94.

Bakr S A,Mannseth T.2009.Feasibility of simplified integral equation modeling of low-frequency marine CSEM with a resistive target.Geophysics,74(5):F107-F117.

Börner R U.2010.Numerical modelling in geo-electromagnetics:Advances and challenges.Surveys in Geophysics,31(2):225-245.

Cai H Z,Xiong B,Han M R,et al.2014.3D controlled-source electromagnetic modeling in anisotropic medium using edge-based finite element method.Computers &Geosciences,73:164-176.

Chen H,Deng J Z,Tan H D,et al.2011.Study on divergence correction method in three-dimensional magnetotelluric modeling with staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),54(6):1649-1659,doi:10.3969/j.issn.0001-5733.2011.06.025.

Deng J Z,Tan H D,Chen H,et al.2011.CSAMT 3D modeling using staggered-grid finite difference method.Progress in Geophysics (in Chinese),26(6):2026-2032,doi:10.3969/j.issn.1004-2903.2011.06.017.Farquharson C G,Duckworth K,Oldenburg D W.2006.Comparison of integral equation and physical scale modeling of the electromagnetic responses of models with large conductivity contrasts.Geophysics,71(4):G169-G177.Haber E,Ascher U M,Aruliah D A,et al.2000.Fast simulation of 3D electromagnetic problems using potentials.Journal of Computational Physics,163(1):150-171.

Hou J S,Mallan R K,Torres-Verdín C.2006.Finite-difference simulation of borehole EM measurements in 3D anisotropic media using coupled scalar-vector potentials.Geophysics,71(5):G225-G233.Mackie R L,Madden T R,Wannamaker P E.1993.Three-dimensional magnetotelluric modeling using difference equations—Theory and comparisons to integral equation solutions.Geophysics,58(2):215-226.Miensopust M P,Queralt P,Jones A G.2013.Magnetotelluric 3-D inversion—a review of two successful workshops on forward and inversion code testing and comparison.Geophysical Journal International,193(3):1216-1238.Mukherjee S,Everett M.2011.3D controlled-source electromagnetic edge-based finite element modeling of conductive and permeable heterogeneities.Geophysics,76(4):F215-F226.Muoz G.2014.Exploring for geothermal resources with electromagnetic methods.Surveys in Geophysics,35(1):101-122.Nabighian M N.1988.Electromagnetic Methods in Applied Geophysics.Tulsa:Society of Exploration Geophysicists.Newman G A.2014.A review of high-performance computational strategies for modeling and imaging of electromagnetic induction data.Surveys in Geophysics,35(1):85-100.

Pardo D,Nam M J,Torres-Verdín C,et al.2011.Simulation of marine controlled source electromagnetic measurements using a parallel Fourier hp-finite element method.Computational Geosciences,15(1):53-67.Schwarzbach C,Haber E.2013.Finite element based inversion for time-harmonic electromagnetic problems.Geophysical Journal International,193(2):615-634.

Selway K.2014.On the causes of electrical conductivity anomalies in tectonically stable lithosphere.Surveys in Geophysics,35(1):219-257.

Sheard S N,Ritchie T J,Christopherson K R,et al.2005.Mining,environmental,petroleum,and engineering industry applications of electromagnetic techniques in geophysics.Surveys in Geophysics,26(5):653-669.

Shen J S.2003.Modelling of 3-D electromagnetic responses in frequency domain by using staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),46(2):281-288.Siripunvaraporn W.2012.Three-dimensional magnetotelluric inversion:an introductory guide for developers and users.Surveys in Geophysics,33(1):5-27.

Smith J T.1996a.Conservative modeling of 3-D electromagnetic fields,Part I:Properties and error analysis.Geophysics,61(5):1308-1318.

Smith J T.1996b.Conservative modeling of 3-D electromagnetic fields,Part II:Biconjugate gradient solution and an accelerator.Geophysics,61(5):1319-1324.Smith R.2014.Electromagnetic Induction Methods in Mining Geophysics from 2008 to 2012.Surveys in Geophysics,35(1):123-156.Song Y,Kim H J,Lee K H.2002.An integral equation representation of wide-band electromagnetic scattering by thin sheets.Geophysics,67(3):746-754.Strack K M.2014.Future directions of electromagnetic methods for hydrocarbon applications.Surveys in Geophysics,35(1):157-177.

Tan H D,Tong T,Lin C H.2006.The parallel 3D magnetotelluric forward modeling algorithm.Applied Geophysics,3(4):197-202.

Tan H D,Yu Q F,Booker J,et al.2003.Magnetotelluric three-dimensional modeling using the staggered-grid finite difference method.Chinese Journal of Geophysics (in Chinese),46(5):705-711.Wannamaker P E.1991.Advances in three-dimensional magnetotelluric modeling using integral equations.Geophysics,56(11):1716-1728.

Weiss C J,Constable S.2006.Mapping thin resistors and hydrocarbons with marine EM methods,Part II—Modeling and analysis in 3D.Geophysics,71(6):G321-G332.

Yahya N,Akhtar M N,Nasir N,et al.2012.Guided and direct wave evaluation of controlled source electromagnetic survey using finite element method.Journal of Electromagnetic Analysis and Applications,4(3):135-146.

Yin C C,Ben F,Liu Y H,et al.2014.MCSEM 3D modeling for arbitrarily anisotropic media.Chinese Journal of Geophysics (in Chinese),57(12):4110-4122,doi:10.6038/cjg20141222.Yoshimura R,Oshiman N.2002.Edge-based finite element approach to the simulation of geoelectromagnetic induction in a 3-D sphere.Geophysical Research Letters,29(3):9-1-9-4.Zhang Y,Wang H N,Tao H G,et al.2012.Finite volume algorithm to simulate 3D responses of multi-component induction tools in inhomogeneous anisotropic formation based on coupled scalar-vector potentials.Chinese Journal of Geophysics (in Chinese),55(6):2141-2152,doi:10.6038/j.issn.0001-5733.2012.06.036.

Zhdanov M S,Lee S K,Yoshioka K.2006.Integral equation method for 3D modeling of electromagnetic fields in complex structures with inhomogeneous background conductivity.Geophysics,71(6):G333-G345.

Zhou J M,Zhang Y,Wang H N,et al.2014.Efficient simulation of three-dimensional marine controlled-source electromagnetic response in anisotropic formation by means of coupled potential finite volume method.Acta Physica Sinica (in Chinese),63(15):159101.

附中文参考文献

陈辉,邓居智,谭捍东等.2011.大地电磁三维交错网格有限差分数值模拟中的散度校正方法研究.地球物理学报,54(6):1649-1659,doi:10.3969/j.issn.0001-5733.2011.06.025.

邓居智,谭捍东,陈辉等.2011.CSAMT三维交错采样有限差分数值模拟.地球物理学进展.26(6):2026-2032,doi:10.3969/j.issn.1004-2903.2011.06.017.

沈金松.2003.用交错网格有限差分法计算三维频率域电磁响应.地球物理学报,46(2):281-288.

谭捍东,余钦范,Booker J等.2003.大地电磁法三维交错采样有限差分数值模拟.地球物理学报,46(5):705-711.

殷长春,贲放,刘云鹤等.2014.三维任意各向异性介质中海洋可控源电磁法正演研究.地球物理学报,57(12):4110-4122,doi:10.6038/cjg20141222.

张烨,汪宏年,陶宏根等.2012.基于耦合标势与矢势的有限体积法模拟非均匀各向异性地层中多分量感应测井三维响应.地球物理学报,55(6):2141-2152,doi:10.6038/j.issn.0001-5733.2012.06.036.

周建美,张烨,汪宏年等.2014.耦合势有限体积法高效模拟各向异性地层中海洋可控源的三维电磁响应.物理学报,63(15):159101.

(本文编辑何燕)

A finite-volume solution to 3D frequency-domain electromagnetic modelling using Lorenz-gauged magnetic vector and scalar potentials

CHEN Hui1,2,YIN Chang-Chun1*,DENG Ju-Zhi2

1 Geo-Exploration Science and Technology Institute,Jilin University,Changchun 130026,China 2 Key Laboratory of Radioactive Geology and Exploration Technology Fundamental Science for National Defense， East China Institute of Technology，Nanchang 330013,China

To eliminate the coupling between air and the earth and to avoid calculation of primary field,we develop a universal algorithm on forward modelling of the frequency-domain electromagnetic (EM)method for different source types and applications.We first present a magnetic vector and a scalar potential with Lorenz gauge based on Maxwell equation,then divide different sources into a lot of short wires.Next,we use staggered grids and the finite volume method to discrete the potential equations and implement the quasi-minimum-residual (QMR)iteration with Jacob to solve the large sparse and symmetrical linear equations system.Finally,we accomplish frequency-domain EM modelling for different source types and application areas.Through analyses and comparison of 1D and 3D MT,CSAMT and loop source responses,we demonstrate the efficiency and accuracy of the algorithm proposed in this study.Based on this,we analyze EM responses for different source types.The algorithm and numerical results presented in this paper build a framework for 3D frequency-domain EM forward modelling and inversion.KeywordsFDEM;Finite-volume method;Forward modelling;Lorenz-gauged;Coupled potentials

陈辉,殷长春,邓居智.2016.基于Lorenz规范条件下磁矢势和标势耦合方程的频率域电磁法三维正演.地球物理学报，59(8):3087-3097,

10.6038/cjg20160831.

Chen H,Yin C C,Deng J Z.2016.A finite-volume solution to 3D frequency-domain electromagnetic modelling using Lorenz-gauged magnetic vector and scalar potentials.Chinese J.Geophys.(in Chinese),59(8):3087-3097,doi:10.6038/cjg20160831.

国家青年基金项目(41404057),国家自然科学基金项目(41164003),国家重大科研装备研究项目(ZDYZ2012-1-03和20130523MTEM05)联合资助.

陈辉，男，1985年8月生，博士生，主要从事电磁勘查技术正反演研究.E-mail:schoolhui@163.com

殷长春，男，1965年生，教授，国家“千人计划”特聘专家，主要从事电磁勘探理论，特别是航空和海洋电磁方面的研究.

E-mail:yinchangchun@jlu.edu.cn

10.6038/cjg20160831

P631

2015-01-21，2015-10-08收修定稿