2018, Vol. 61
回线源瞬变电磁法作为一种重要的电磁勘探方法，以其显著的优点(经济、无损、快速、信息丰富、分辨率强等)广泛应用于金属矿产勘探、煤田灾害探测、地下断层划分、地下水探测、隧道超前地质预报等领域(李貅和薛国强，2013).三维正演能够为勘探最优化设计提供有力帮助，是三维反演和成像技术发展的前提.薛国强等(2008)，Börner(2010)，李建慧和朱自强(2012)综述了瞬变电磁三维正演算法的研究进展.
按照空间离散的不同处理方法，模拟复杂地层模型的瞬变电磁三维正演算法主要包括有限差分法(Wang and Hohmann, 1993; 孙怀凤和李貅，2013；李展辉和黄清华，2014)，有限体积法(Haber and Oldenburg, 2002)，有限元法(Um, 2010).其中有限体积法是一种基于守恒特性的空间离散方法.该方法将计算区域离散为一个个控制体积，在控制体积的离散网格上对原始方程进行积分，利用斯托克斯定理得到控制体积边界上的环量积分形式，通过对边界上的环量进行离散得到相应的离散方程.有限差分法是直接根据微分方程导出，不涉及积分过程，导数离散借助于Taylor展开，离散精度主要取决于对导数处理的精度，因此不能保证物理量的守恒特性(Mackie et al., 1994).与有限差分法相比，有限体积离散是根据积分方程导出的(对每一个控制体积进行积分)，其精度不仅取决于对导数处理的精度，还取决于积分的精度，而且积分过程能够保证物理量的守恒特性.与有限元法一样，有限体积法也适用于结构化(Haber et al., 2002)或非结构化网格(Jahandari and Farquharson, 2015).基于这些特点，有限体积法在电磁场三维正演模拟中得到了较广泛的应用(Haber et al., 2002; Weiss and Constable, 2006; Haber and Heldmann, 2007; 杨波等, 2012; 张烨等, 2012; Marchant et al., 2014; Jahandari and Farquharson, 2015).
近年来一种改进的有限体积空间离散方法——基于模拟离散的有限体积法(MFVM)得到了大力的发展.模拟有限体积算法除了保持有限体积法的守恒特性和适用于非结构化网格的特性外，该算法还能够更有效的处理复杂构造，比如各向异性和强非均匀性区域的离散，同时避免产生各种非物理伪解(Castillo and Miranda, 2013).Hyman和Shashkov(1997a, 1997b)系统的给出了MFVM的数学基础，基于支撑算子理论构建了一次算子和导出算子(梯度GRAD，旋度CRUL和散度DIV)在逻辑正交网格下的离散，实现了任意二阶算子(比如DIV GRAD，CURL CURL)的离散，并给出了精确的内积离散方法.MFVM已经应用于多孔介质(Hyman et al., 2001)、流体动力学(Barbosa and Daube, 2005)、图像处理(Yuan and Sheng, 2007)、广义相对论(Di Bartolo, 2004)以及电磁场(Hyman and Shashkov, 1999; Haber and Ruthotto, 2014；彭荣华和胡祥云, 2016)的数值模拟.Lipnikov等(2014)总结了MFVM的最新进展.
随着电磁勘探方法的发展，越来越多的学者认识到电导率各向异性在地下构造中分布的普遍性以及对电磁响应数据解释的影响(Everett and Constable, 1999; Collins, 2006; Li and Dai, 2011; Dennis, 2012; Mattsson, 2013; Martí, 2014; Yin, 2016).本文研究采用模拟有限体积算法实现双轴各向异性地层(Everett and Constable, 1999)回线源瞬变电磁三维正演.首先引入内积定义，采用简单的自然边界条件，将瞬变电磁法的控制方程转化为弱形式表示.将计算区域划分为一系列的控制体积单元，采用交错网格对弱形式的控制方程进行有限体积空间离散，包括旋度算子离散和空间内积离散(Hyman and Shashkov, 1999).旋度算子离散基于斯托克斯定理的旋度积分定义公式.采用中点原则实现电导率双轴各向异性的空间内积离散，从而得到离散化的控制方程.时间步迭代采用无条件稳定的欧拉后向差分格式(Haber et al., 2002).并通过均匀全空间中稳定电流回线源的磁场解析表达式，给出回线源初始时刻的电磁场分布.为了同时保证计算精度和效率，本文采用分段等间隔的时间步迭代，利用直接法实现其快速求解.最后通过各向异性半空间模型、层状模型和三维模型的正演计算结果对比，验证本文算法的计算精度，通过典型三维各向异性模型的计算分析各向异性对瞬变电磁响应的影响.
1 正演算法瞬变电磁三维正演算法分为空间离散和时间迭代两部分内容.本文采用模拟有限体积算法实现空间离散，采用隐式时间步迭代算法实现时间迭代.
1.1 空间离散的模拟有限体积算法 1.1.1 控制方程瞬变电磁法对应的时间域Maxwell方程，忽略位移电流：
|
(1a) |
|
(1b) |
其中E是电场强度矢量，B是磁感应强度矢量，t是时间，σ=diag(σx, σy, σz)是各向异性的电导率张量，μ是磁导率，s是外加源项.
为了求解方程(1)，还需要补充边界条件和初始条件.本文采用简单的自然边界条件(Haber and Ruthotto, 2014)：
|
(2a) |
|
(2b) |
对于回线源瞬变电磁法，在t=0时刻空间中只存在稳定的磁场分布，即初始条件为
|
(3) |
其中B0为t=0时刻空间中的磁场分布.
1.1.2 弱形式表示与有限元方法一样，MFVM处理弱形式的控制方程.按照泛函分析理论，E和B位于不同的Sobolev空间(Volakis等, 2006)，即E∈H(Curl; Ω)，B∈H(Div; Ω).定义内积为：
|
(4) |
其中A和G为空间H(Curl; Ω)或H(Div; Ω)中的任意参数.
为了得到弱形式的控制方程，引入与E位于相同的Sobolev空间的参数W与B位于相同的Sobolev空间的参数F.将方程(1a)与F做内积，方程(1b)与W做内积，得到：
|
(5a) |
|
(5b) |
利用分部积分公式:
|
(6) |
根据边界条件B×n = 0|∂Ω，利用公式(6)，得到弱形式的控制方程：
|
(7a) |
|
(7b) |
由方程(7)可知，由于E和W位于相同的Sobolev空间H(Curl; Ω)，因此空间离散只需要求取空间H(Curl; Ω)参数的旋度即可.对比方程(1)，弱形式控制方程不需要求取B的旋度，由于求取旋度是一种微分运算，因此弱形式控制方程弱化了对电磁场的可微性.
1.1.3 空间离散设定整个网格剖分区域的6个面均为规则矩形，采用正交规则矩形网格，在x、y和z三个方向上离散网格单元数为nx、ny和nz.典型的控制体积网格单元如图 1所示，定义网格中心点为(i, j, k).根据方程(1a)，E的旋度对应B，定义E在网格棱边中心，三个方向的场点分别为Ei, j±1/2, k±1/2x、Ei±1/2, j, k±1/2y和Ei±1/2, j±1/2, kz，定义B在网格面中心，三个方向的场点分别为Bi±1/2, j, kx、Bi, j±1/2, ky和Bi, j, k±1/2z.
|
图 1 网格单元(i, j, k) Fig. 1 Grid cell (i, j, k) |
由公式(7)可知，空间离散主要包含两部分内容：旋度算子离散和空间内积离散.MFV采用积分形式的斯托克斯定理处理电场旋度的离散.对于网格单元(i, j, k)，在x、y和z三个方向上的单元网格长度记为hxi、hyj和hzk.单元网格6个表面分别记为Si±1/2, j, k、Si, j±1/2, k和Si, j, k±1/2.以表面Si+1/2, j, k为例，计算电场的旋度x方向的分量，即关于表面Si+1/2, j, k的投影为
|
(8) |
根据中点积分规则可以得到方程(8)右端项线积分的离散表示为
|
(9) |
同理计算电场的旋度y方向分量为
|
(10) |
计算电场的旋度z方向分量为
|
(11) |
根据公式(8)—(11)，可以将旋度算子整理为矩阵形式：
|
(12) |
其中，P为包含剖分网格所有表面面积的对角矩阵，L为包含剖分网格所有棱边长度的对角矩阵，e是网格中电场E的矩阵表示形式，C为包含0和±1的矩阵，其形式表达为
|
其中，nb为剖分网格所有表面数，ne为剖分网格所有棱边数，Dij(i, j=x, y, z)为各方向上的差分矩阵，如
|
nx, ny, nz分别表示x, y, z方向网格数.
由公式(7)可知，空间内积离散包含两种不同类型：位于面中心点的参数内积(B, F)和位于棱边中心的参数内积(σE, W).设定网格单元内部的电导率均一，内积计算可以采用简单的中心点算术平均求得.
对于网格单元(i, j, k)中，面中心点参数内积，即：
|
(13) |
由于B和F都定义在网格面中心，对于一个网格单元，在一个方向上有2个面中心，因此：
|
(14a) |
|
(14b) |
|
(14c) |
其中vi, j, k是网格单元(i, j, k)的体积.综合公式(13)和(14)，采用矩阵形式表示，即为：
|
(15) |
其中，f和b为F和B的矩阵表示形式，Mf的具体形式如下：
|
(16) |
其中，v为包含所有剖分网格单元体积的对角矩阵，Afr(r=x, y, z)分别对应Bx、By和Bz的平均，具体形式为
|
(17) |
其中，nc为剖分网格单元数，nr为剖分网格r方向表面数.
棱边中心的参数内积(σE, W)离散，同样是将参数平均到网格单元中心位置，与面中心点参数内积离散方法类似，主要的不同点在于：对于一个网格单元，在每个方向上有4条棱边，因此每条棱边上的场值比重为1/4(如图 2所示).定义在网格体中心处的电导率是各向异性的，即沿不同方向的电导率是不同的.电导率各向异性的这一特性正好与内积计算相结合(如图 2所示)，不同方向的内积(σE, W)计算采用不同的电导率，从而得到棱边中心的参数内积(σE, W)离散的矩阵表示为
|
图 2 网格单元中E与σ位置分布图 Fig. 2 Location of E and σ in grid cell |
|
(18) |
其中w为W的矩阵表示形式，Mσe的具体形式为
|
(19) |
vσx，vσy和vσz分别为包含所有网格单元体积与该网格单元不同方向电导率的乘积的对角矩阵，Aer(r=x, y, z)代表定义在网格棱边中心的Er平均到网格中心的转换矩阵.
1.1.4 矩阵表示在得到了旋度离散和内积离散的具体形式后，我们可以给出控制方程(7)的完整离散，用矩阵形式表示，即为
|
(20a) |
|
(20b) |
由于f和w是任意引入的参数F和W对应的矩阵形式，因此公式(20)两边可以消去f和w，得到控制方程空间离散的矩阵表示为
|
(21a) |
|
(21b) |
时间步离散采用无条件稳定的欧拉后向差分格式(Haber and Oldenburg, 2002)，即：
|
(22) |
将公式(22)代人到公式(21)，得到离散控制方程为
|
(23a) |
|
(23b) |
由公式(21a)可以得到btn=-CURLen.
1.2.1 初始场对于回线源瞬变电磁装置，在电流关断之前，空间中只存在稳定电流产生的静态磁场，该静态磁场分布与模型电导率无关，该静态磁场即为t=0时刻空间中的磁场分布b0=b(0).如果不考虑模型中的磁导率变化，即假定地下模型的磁导率与空气磁导率相同，均为真空磁导率μ0，则柱坐标系中全空间的磁矢势可以解析的表示为(Bladel, 2007)：
|
(24) |
其中回线源中心点位置(x0, y0, z0)，回线源半径为a，接收点位置(xr, yr, zr)，x=xr－x0，y=yr－y0，z=zr－z0，
|
(25) |
即可得到初始时刻的磁场.对比公式(25)和公式(1a)，可知A与E位于相同的Sobolev空间，即A与E在离散网格中位于同一位置.利用公式(24)，计算离散网格棱边中点位置的磁矢势，利用公式(25)，即可以得到t=0时刻空间中的磁场分布，表示为矩阵形式：
|
(26) |
采用磁矢势求得磁场而不是直接采用解析形式的磁场表示式，主要是考虑到公式(26)能够保证初始时刻的磁场b0是无散的(Hyman and Shashkov, 1999)，即∇·B(0)≡0，根据模拟有限体积法的特点，则能够保证之后的磁场bn总是无散的.
1.2.2 线性方程组的求解在得到初始场b0之后，通过求解时间步迭代的线性方程组(23)(sn=0)即可得到不同时刻的电磁场响应.该方程组的求解可以采用迭代法(Haber et al., 2002)或直接法(Um, 2010).直接法对矩阵条件数不敏感，能有效处理多右端项问题(Um, 2010; Yin, 2016)，本文采用直接法求解器PARDISO (Schenk and Gartner, 2004)求解线性方程组(23).为了同时保证计算精度和效率，本文选取分段等间隔的时间步长(Um, 2010).
2 数值实例采用层状地层模型和各向异性均匀半空间模型，通过与解析解的比较，验证本文算法(简称为MFVTD)的计算精度.通过对典型三维各向异性模型的计算，分析各向异性对瞬变电磁响应的影响.本文计算设备为32 G内存、四核主频3.6 GHz的Intel i7 CPU的台式电脑.
2.1 层状地层模型首先采用层状地层模型，通过与一维解析解结果(Ward and Hohmann, 1998)对比，验证本文算法的计算精度.层状地层模型参数如图 3所示，空气层电导率设置为10－6S·m-1，在电导率为0.01 S·m-1半空间地层中存在一个顶部埋深50 m、层厚50 m、电导率0.1 S·m-1的低阻层.发射源为地表圆形回线框，半径为10 m，发射电流为1 A，接收回线中心点的dBz/dt和Bz.采用非均匀网格剖分，最小网格长度为5m，网格放大系数为1.3，总的网格单元数为37×37×52.最小时间步长为1×10－7s，每间隔200步时间步长增大一倍，总的时间步迭代次数为1800次.采用直接法求解器PARDISO求解，需9次系数矩阵分解，1800次方程求解，总的计算时间为350 s.与1D解析解的对比结果如图 4所示，其中图 4a为接收回线中心点的dBz/dt响应，图 4b为dBz/dt三维响应与一维响应的相对误差；图 4c为接收回线中心点的Bz响应，图 4d为Bz三维响应与一维响应的相对误差.由图 4可知，Bz三维响应与一维响应的相对误差都在3%以内，dBz/dt三维响应与一维响应的相对误差稍大，但除了早期(0.02 ms以内)的相对误差略大，其他时间区域的相对误差都3%以内.采用本文MFVM算法计算的层状模型响应与解析解吻合的很好，说明本文算法是有效的.
|
图 3 层状地层模型 Fig. 3 Layered earth model |
|
图 4 层状模型3D解与1D解的对比 (a) dBz/dt；(b) dBz/dt相对误差；(c) Bz；(d) Bz相对误差. Fig. 4 Comparison of 3D results with 1D solutions for layered earth model (a) dBz/dt; (b) Relative errors for dBz/dt; (c) Bz; (d) Relative errors for Bz. |
采用Yin(2016)的各向异性检验模型，验证本文算法计算各向异性地层响应的有效性.模型如图 5所示，发射源为边长20 m的方形回线源，发射电流为1 A，距离地面高度为30 m，接收回线中心点的dBz/dt和Bz.分别计算地层电阻率为σ=
|
图 5 各向异性地层模型 Fig. 5 Anisotropic earth model |
|
图 6 MFVM解与Yin(2016)解的对比 (a) dBz/dt；(b) Bz. Fig. 6 Comparison of results by MFVM with Yin (2016) |
采用Li(2017)的三维垂直接触带模型，通过与三维矢量有限元计算结果比较，进一步验证本文算法的计算精度.该模型最早是由Commer和Newman(2004)在计算电性源瞬变电磁响应时设计的.模型如图 7所示，地表下方是厚度50 m、电导率为0.1 S·m-1的覆盖层，覆盖层下由两部分的垂直接触带构成，电导率分别为0.01 S·m-1和0.0033333 S·m-1.垂直接触带中间存在一个电导率为1 S·m-1的三维复杂形状的低阻体，沿走向的长度为400 m，宽度为100 m，厚度近似为500 m，具体的形状如图 7a所示.空气层电导率设置为10-6 S·m-1.发射线圈为100 m×100 m的方形回线框，发射线框的中心点坐标为(0，50，0)m，4个观测点分别位于(0，50，0)m，(0，150，0)m，(0，450，0)m和(0，1050，0)m.采用非均匀网格剖分，最小网格长度为10 m，网格放大系数为1.4，总的网格单元数为58×40×53，具体的剖分如图 7b所示，总的计算时间为1060 s.本文计算结果与Li(2017)采用的矢量有限元算法的对比结果如图 8所示，由图 8可知，两种方法计算的不同观测点处的响应重合的非常好，进一步说明了本文算法计算的结果是可靠的.
2.4 各向异性分析
计算典型的三维各向异性地层模型的回线源瞬变电磁响应，分析各向异性对瞬变电磁响应的影响.模型参数如图 9a所示，在电阻率为σ1半空间地层中存在一个顶部埋深50 m，长宽高尺寸为100 m×100 m×50 m、电阻率为σ2的低阻异常体.发射源为地表圆形回线框，半径为10 m，发射电流为1 A.计算2个观测点的响应：(1)发射回线中心对应低阻矩形块中心，即图 9a中S1处；(2)发射回线中心偏离低阻矩形块中心70 m，即图 9a中S2处.接收回线中心点的dBz/dt和Bz.分别分析σ1和σ2为各向异性时对瞬变电磁响应的影响.
|
图 9 三维地层模型 (a)模型示意图；(b)网格剖分图. Fig. 9 Layered earth model (a) Model; (b) Grid for model. |
首先计算发射回线中心位于S1处的瞬变电磁响应.计算σ2为各向异性对瞬变电磁响应的影响.即设定σ1=0.01 S·m-1，分别讨论电导率各向同性σ2=
|
图 10 S1处不同各向异性σ2模型的瞬变电磁响应 (a) dBz/dt；(b) Bz. Fig. 10 TEM response in S1 for different anisotropic σ2. |
计算σ1为各向异性对瞬变电磁响应的影响.即设定σ2=0.1S·m-1，分别讨论电导率各向同性σ1=
|
图 11 S1处不同各向异性σ1模型的瞬变电磁响应 (a) dBz/dt；(b) Bz. Fig. 11 TEM response in S1 for different anisotropic σ1 |
接下来计算发射回线中心位于S2处的瞬变电磁响应.分别计算σ2为向异性和σ1为各向异性对瞬变电磁响应的影响.模型参数设置和网格剖分与发射回线中心位于S1处时保持相同，计算结果如图 12和图 13所示，Bz和dBz/dt响应与发射源位于S1处的响应特征基本一致.综合图 10—13，可知主要是水平方向电导率各向异性(即σx=σz≠σy)对电磁响应产生影响，垂直方向的电导率各向异性(即σx=σy≠σz)对Bz和dBz/dt响应几乎没有影响.产生这一现象的原因是回线源产生的感应电流主要是水平方向的，因此主要受到水平方向电导率的影响，而垂直方向的电导率影响很小.
|
图 12 S2处不同各向异性σ2模型的瞬变电磁响应 (a) dBz/dt；(b) Bz. Fig. 12 TEM response in S2 for different anisotropic σ2 |
|
图 13 S2处不同各向异性σ1模型的瞬变电磁响应 (a) dBz/dt；(b)Bz. Fig. 13 TEM response in S2 for different anisotropic σ1 |
本文采用模拟有限体积算法实现了双轴各向异性地层回线源瞬变电磁三维正演.从弱形式表示的控制方程出发，采用交错网格对控制方程进行模拟有限体积空间离散，包括旋度算子离散和空间内积离散.基于斯托克斯定理的旋度积分定义公式实现旋度算子离散.中点平均实现电导率双轴各向异性的空间内积离散.时间步迭代采用无条件稳定的欧拉后向差分格式.并通过均匀全空间中稳定电流回线源的磁场解析表达式得到初始时刻的电磁场分布.为了同时保证计算精度和效率，本文采用分段等间隔的时间步迭代，利用直接法求解器PARDISO实现其快速求解.
通过层状模型和各向异性半空间模型的正演计算结果对比，验证了本文算法的计算精度和计算效率.计算三维双轴各向异性模型的正演响应可知，水平方向电导率变化对电磁响应产生显著影响，而垂直方向的电导率变化对电磁响应几乎没有影响.产生这一现象的主要原因是回线源产生的感应电流主要是水平方向的，因此响应主要受到水平方向电导率的影响，垂直方向的电导率影响很小.
模拟有限体积算法能够有效处理复杂模型的三维正演计算.未来的工作包括引入更复杂的内积形式、采用逻辑正交网格和非结构化网格，实现完全各向异性、包含不规则形状的复杂模型的三维快速正演响应.
Barbosa E, Daube O.
2005. A finite difference method for 3D incompressible flows in cylindrical coordinates. Computers & Fluids, 34(8): 950-971.
|
|
Bladel J G V. 2007.
Electromagnetic Fields. (2nd ed). Wiley IEEE Press.
|
|
Börner R U.
2010. Numerical modelling in geo-electromagnetics:Advances and challenges. Surveys in Geophysics, 31(2): 225-245.
DOI:10.1007/s10712-009-9087-x |
|
Castillo J E, Miranda G F. 2013.
Mimetic discretization methods. Boca Raton: CRC Press.
|
|
Collins J L, Everett M E, Johnson B.
2006. Detection of near-surface horizontal anisotropy in a weathered metamorphic schist at Llano Uplift (Texas) by transient electromagnetic induction. Physics of the Earth and Planetary Interiors, 158(2-4): 159-173.
DOI:10.1016/j.pepi.2006.05.008 |
|
Commer M, Newman G.
2004. A parallel finite difference approach for 3D transient electromagnetic modeling with galvanic sources. Geophysics, 69(5): 1192-1202.
DOI:10.1190/1.1801936 |
|
Dennis Z R, Cull J P.
2012. Transient electromagnetic surveys for the measurement of near-surface electrical anisotropy. Journal of Applied Geophysics, 76: 64-73.
DOI:10.1016/j.jappgeo.2011.10.014 |
|
Di Bartolo C, Gambini R, Pullin J.
2004. Consistent and mimetic discretizations in general relativity. Journal of Mathematical Physics, 46(3): 032501.
|
|
Everett M E, Constable S.
1999. Electric dipole fields over an anisotropic seafloor:theory and application to the structure of 40Ma Pacific Ocean lithosphere. Geophysical Journal International, 136(1): 41-56.
DOI:10.1046/j.1365-246X.1999.00725.x |
|
Haber E, Ascher U, Oldenburg D W. 2002. 3D forward modelling of time domain electromagnetic data.//72th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts.
|
|
Haber E, Heldmann S.
2007. An octree multigrid method for quasi-static Maxwell's equations with highly discontinuous coefficients. Journal of Computational Physics, 223(2): 783-796.
DOI:10.1016/j.jcp.2006.10.012 |
|
Haber E, Ruthotto L.
2014. A multiscale finite volume method for Maxwell's equations at low frequencies. Geophysical Journal International, 199(2): 1268-1277.
DOI:10.1093/gji/ggu268 |
|
Hyman J, Morel J, Shashkov M, et al.
2001. Mimetic finite difference methods for diffusion equations. Computational Geosciences, 6(3-4): 333-352.
|
|
Hyman J M, Shashkov M.
1997a. Natural discretizations for the divergence, gradient, and curl on logically rectangular grids. Computers & Mathematics with Applications, 33(4): 81-104.
|
|
Hyman J M, Shashkov M.
1997b. Adjoint operators for the natural discretizations of the divergence, gradient and curl on logically rectangular grids. Applied Numerical Mathematics, 25(4): 413-442.
DOI:10.1016/S0168-9274(97)00097-4 |
|
Hyman J M, Shashkov M.
1999. Mimetic discretizations for Maxwell's equations. Journal of Computational Physics, 151(2): 881-909.
DOI:10.1006/jcph.1999.6225 |
|
Jahandari H, Farquharson C G.
2015. Finite-volume modelling of geophysical electromagnetic data on unstructured grids using potentials. Geophysical Journal International, 202(3): 1859-1876.
DOI:10.1093/gji/ggv257 |
|
Li J H, Zhu Z Q, Zeng S H, et al.
2012. Progress of forward computation in transient electromagnetic method. Progress in Geophysics, 27(4): 1393-1400.
|
|
Li J H, Farquharson C G, Hu X Y.
2017. 3D vector finite-element electromagnetic forward modeling for large loop sources using a total-field algorithm and unstructured tetrahedral grids. Geophysics, 82(1): E1-E16.
DOI:10.1190/geo2016-0004.1 |
|
Li X, Xue G Q. 2013.
Study on Pseudo-seismic Migration Imaging of Transient Electromagnetic Method. Beijing: Science Press.
|
|
Li Y G, Dai S K.
2011. Finite element modelling of marine controlled-source electromagnetic responses in two-dimensional dipping anisotropic conductivity structures. Geophysical Journal International, 185(2): 622-636.
DOI:10.1111/gji.2011.185.issue-2 |
|
Li Z H, Huang Q H.
2014. Application of the complex frequency shifted perfectly matched layer absorbing boundary conditions in transient electromagnetic method modeling. Chinese Journal of Geophysics, 57(4): 1292-1299.
DOI:10.6038/cjg20140426 |
|
Lipnikov K, Manzini G, Shashkov M.
2014. Mimetic finite difference method. Journal of Computational Physics, 2014, 257: 1163-1227.
|
|
Mackie R L, Smith J T, Madden T R.
1994. Three-dimensional electromagnetic modeling using finite difference equations:The magnetotelluric example. Radio Science, 29(4): 923-935.
DOI:10.1029/94RS00326 |
|
Marchant D, Haber E, Oldenburg D W.
2014. Three-dimensional modeling of IP effects in time-domain electromagnetic data. Geophysics, 79(6): E303-E314.
DOI:10.1190/geo2014-0060.1 |
|
Martí A.
2014. The role of electrical anisotropy in magnetotelluric responses:From modelling and dimensionality analysis to inversion and interpretation. Surveys in Geophysics, 35(1): 179-218.
DOI:10.1007/s10712-013-9233-3 |
|
Mattsson J, Engelmark F, Anderson C.
2013. Towed streamer EM:The challenges of sensitivity and anisotropy. First Break, 31(6).
|
|
Peng R H, Hu X Y, Han B, et al.
2016. 3D frequency-domain CSEM forward modeling based on the mimetic finite-volume method. Chinese Journal of Geophysics, 59(10): 3927-3939.
DOI:10.6038/cjg20161036 |
|
Sun H F, Li X, Li S C, et al.
2013. Three-dimensional FDTD modeling of TEM excited by a loop source considering ramp time. Chinese Journal of Geophysics, 2013, 56(3): 1049-1064.
DOI:10.6038/cjg20130333 |
|
Schenk O, Gärtner K.
2004. Solving unsymmetric sparse systems of linear equations with PARDISO. Future Generation Computer Systems, 20(3): 475-487.
DOI:10.1016/j.future.2003.07.011 |
|
Um E S, Harris J M, Alumbaugh D L.
2010. 3D time-domain simulation of electromagnetic diffusion phenomena:A finite-element electric-field approach. Geophysics, 75(4): F115-F126.
DOI:10.1190/1.3473694 |
|
Volakis, J L, Sertel K, Usner B C, 2006. Frequency Domain Hybrid Finite Element Methods in Electromagnetics.//Synthesis Lectures on Computational Electromagnetics. Morgan & Claypool.
|
|
Wang T, Hohmann G W.
1993. A finite-difference, time-domain solution for three-dimensional electromagnetic modeling. Geophysics, 58(6): 797-809.
DOI:10.1190/1.1443465 |
|
Ward S, Hohmann G W. 1998. Electromagnetic theory for geophysical applications.//68th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts.
|
|
Weiss C J, Constable S.
2006. Mapping thin resistors and hydrocarbons with marine EM methods, Part Ⅱ-Modeling and analysis in 3D. Geophysics, 71(6): G321-G332.
DOI:10.1190/1.2356908 |
|
Xue G Q, Li X, Di Q Y.
2008. Research progress in TEM forward modeling and inversion calculation. Progress in Geophysics, 23(4): 1165-1172.
|
|
Yang B, Xu Y X, He Z X, et al.
2012. 3D frequency-domain modeling of marine controlled source electromagnetic responses with topography using finite volume method. Chinese Journal of Geophysics, 55(4): 1390-1399.
|
|
Yin C C, Qi Y F, Liu Y H.
2016. 3D time-domain airborne EM modeling for an arbitrarily anisotropic earth. Journal of Applied Geophysics, 131: 163-178.
DOI:10.1016/j.jappgeo.2016.05.013 |
|
Yuan G W, Sheng Z Q.
2007. Analysis of accuracy of a finite volume scheme for diffusion equations on distorted meshes. Journal of Computational Physics, 224(2): 1170-1189.
DOI:10.1016/j.jcp.2006.11.011 |
|
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, 55(6): 2141-2152.
DOI:10.6038/j.issn.0001-5733.2012.06.036 |
|
李建慧, 朱自强, 曾思红, 等.
2012. 瞬变电磁法正演计算进展. 地球物理学进展, 27(4): 1393–1400.
DOI:10.6038/j.issn.1004-2903.2012.04.013 |
|
李貅, 薛国强. 2013.
瞬变电磁法拟地震偏移成像研究. 北京: 科学出版社.
|
|
李展辉, 黄清华.
2014. 复频率参数完全匹配层吸收边界在瞬变电磁法正演中的应用. 地球物理学报, 57(4): 1292–1299.
DOI:10.6038/cjg20140426 |
|
彭荣华, 胡祥云, 韩波, 等.
2016. 基于拟态有限体积法的频率域可控源三维正演计算. 地球物理学报, 59(10): 3927–3939.
DOI:10.6038/cjg20161036 |
|
孙怀凤, 李貅, 李术才, 等.
2013. 考虑关断时间的回线源激发TEM三维时域有限差分正演. 地球物理学报, 56(3): 1049–1064.
DOI:10.6038/cjg20130333 |
|
薛国强, 李貅, 底青云.
2008. 瞬变电磁法正反演问题研究进展. 地球物理学进展, 23(4): 1165–1172.
|
|
杨波, 徐义贤, 何展翔, 等.
2012. 考虑海底地形的三维频率域可控源电磁响应有限体积法模拟. 地球物理学报, 55(4): 1390–1399.
|
|
张烨, 汪宏年, 陶宏根, 等.
2012. 基于耦合标势与矢势的有限体积法模拟非均匀各向异性地层中多分量感应测井三维响应. 地球物理学报, 55(6): 2141–2152.
DOI:10.6038/j.issn.0001-5733.2012.06.036 |
|