地球物理学报  2017, Vol. 60 Issue (5): 1937-1945   PDF    
无限长源二维各向异性地层中CSAMT有限元模拟
熊治涛1,2, 唐新功1,2     
1. “油气资源与勘探技术”教育部重点实验室 (长江大学), 武汉 430100;
2. 长江大学非常规油气湖北省协同创新中心, 武汉 430100
摘要: 可控源音频大地电磁测深法(CSAMT)基本上是基于电性各向同性理论的假设,而实际地层普遍表现为电性各向异性现象,因此研究各向异性地层条件下CSAMT方法的电磁响应具有重要的理论意义和应用价值.本文从有源Maxwell方程出发,推导出了二维各向异性地层条件下无限长线源频率域的电磁场响应,得到了一组平行于线源方向的电场分量Ex和磁场分量Hx的偏微分方程,并采用有限元法进行求解.模型的计算结果表明:不论是对均匀半空间还是一维层状介质中的各向异性二维体模型的计算,均能从结果中明显地看出各向同性背景场和各向异性异常体的存在;在计算中通过改变各向异性系数或各向异性旋转角,其结果在视电阻率和相位曲线上均存在显著区别.计算结果表明了在各向异性地层中,仍然采用各向同性的假设进行CSAMT资料的处理与解释,可能会带来比较严重的误差,这种情况下必须采用基于各向异性的地质模型对CSAMT资料进行处理和解释.研究结果对于提高CSAMT方法的勘探效果和应用水平具有重要的指导意义.
关键词: 各向异性      CSAMT      有限元      正演      线源     
2D CSAMT finite element simulation of anisotropic formation with an infinite source
XIONG Zhi-Tao1,2, TANG Xin-Gong1,2     
1. Key Laboratory of Exploration Technologies for Oil and Gas Resources of MOE, Yangtze University, Wuhan 430100, China;
2. Hubei Cooperative Innovation Center of Unconventional Oil and Gas, Wuhan 430100, China
Abstract: The controlled source audio-frequency magnetotellurics (CSAMT) is traditionally based on the electrical isotropic theory. However, it has been well recognized that the electrical anisotropic is widely present in the earth interior. So it is significant theoretically and practically to calculate EM responses with anisotropic earth model. Based on Maxwell equations with source terms, EM responses to an infinite linear source in 2D anisotropic formation in the frequency domain are derived. A set of partial differential equations of the electric field component Ex and the magnetic field component Hx parallel to the linear source direction are obtained and solved by the finite element method. The results indicate that the isotropic background and anisotropic anomalous body can be clearly figured out by calculating the 2D anisotropic anomalous body in either homogeneous half-space or one-dimensional layered media. When changing the anisotropic coefficients and rotation angle from main electrical axis of the anisotropic model, the results are obviously different from the anisotropic case in both apparent resistivity and impedance phase curves. The results indicate that errors will be very large if we still use electrical isotropic theory when dealing with CSAMT data in anisotropic media. In such situation, the calculation should be based on anisotropic theory instead of electrical anisotropic one. The result is significant for improvement of the application of the CSAMT method to exploration.
Key words: Anisotropy      CSAMT      Finite element      Forward modeling      Linear source     
1 引言

人工源电磁法相对于天然源方法而言优势突出,因此有不少的学者对其进行研究.在数值模拟方面 (唐新功等,2000陈小斌和胡文宝,2002闫述等,2002底青云等, 2004, 2009王若等,2006薛国强等,2007李展辉等,2009Huang et al., 2010) 等做了大量的工作,但这些多主要针对的是各向同性介质的情况,而实际地层在一定尺度范围内总是表现为各向异性的 (Løseth and Ursin, 2007).对于可控源音频大地电磁法 (CSAMT) 而言,当地层为各向同性情况时,CSAMT方法可以较好地提供地下地电结构的信息;但对于各向异性地层,若仍然按照各向同性的假设来处理,则结果通常不可靠甚至有时会是错误的反演解释结果 (Ramananjaona et al., 2011).因此研究电性各向异性问题对于CSAMT方法具有重要的理论与应用价值.

对于电磁法勘探资料中所呈现出的电性各向异性特征,一直是国内外学者关注的焦点.早在1845年Faraday就进行了电参数各向异性的研究,Schlumberger于1920年在进行电法勘探时发现了沉积地层的电导率各向异性现象,然而直到20世纪60年代对这个问题的研究才变得活跃起来,之后在理论上不断发展和完善.Parkhomenko (1967)Keller (1968)Hill (1972)分别讨论了岩石微观各向异性问题.O’Brien和Morrison (1967)研究了一维各向异性地层中电磁场的递推公式;Reddy和Rankin (1971, 1975) 分别研究了倾斜各向异性下的大地电磁场问题以及使用有限元方法对二维电性各向异性特征进行模拟;Pek等 (Pek and Verner, 1997; Pek and Santos, 2002) 在前人研究的基础上将有限差分法引入到对MT二维电性各向异性问题的求解中,并于2002年推导了一种联合计算层状各向异性介质中大地电磁阻抗张量和参数灵敏度矩阵的算法,较为全面、完善地阐述了1D各向异性的正演问题;Li (Li, 2002; Li and Pek, 2008) 详述了用有限元方法计算二维各向异性介质的电磁场问题,并于2008年采用自适应网格剖分法,取得了较为显著的应用效果;胡祥云等 (2013)霍光谱等 (2015)基于Pek (Pek and Verner, 1997; Pek and Santos, 2002) 的研究理论,采用有限差分法分别研究了MT二维各向异性以及带地形的二维各向异性正演问题,并将其应用于对实测资料的处理解释中.

近年来,对有源电磁法各向异性问题的研究也很活跃.Xiong (1989)给出了垂直各向异性介质中水平电偶极源激发的结果.Li等 (Li and Pedersen, 1991, 1992; Li et al., 2000) 讨论了层状方位角各向异性地层中CSAMT的正、反演问题.岳瑞永和徐义贤 (2004)讨论了电方位各向异性条件下CSAMT的近场效应问题,并以相对偏差作为划分近、远区场的判断条件.Løseth和Ursin (2007)提出了波数域中水平电偶极子激发下层状各向异性介质数值计算的理论,该算法使用基于矩阵传播算子来计算电磁场.罗鸣和李予国 (2015)采用二维傅里叶变化将根据Løseth和Ursin (2007)算法得到的波数域的场转换至空间域,并讨论了一维层状电阻率任意各向异性介质中的海洋可控源电磁 (MCSEM) 响应.陈桂波等 (陈桂波, 2009; 陈桂波等, 2009) 通过对三维积分方程算法的研究,并将其应用于各向异性地层中,实现了三维MT、CSAMT以及MCSEM的正演模拟.周建美等 (2012)实现了一维垂直各向异性地层介质中频率测深资料的全参数反演.Puzyrev等 (2013)以各向同性的全空间海水为背景并基于Coulomb位,采用有限元方法实现了海洋电磁三维各向异性模型数据的并行计算.蔡红柱等 (2015)在Badea和Puzyrev等人的基础上,以更接近真实情况的层状各向异性介质为背景,使用有限单元法实现了海洋电磁三维正演计算.殷长春等 (2014)采用交错网格有限差分,通过对张量电导率实行体积和空间电流密度平均,成功地离散了二次散射电场,实现了对三维任意各向异性介质下MCSEM方法的正演计算.尽管国内外的学者们一直在努力研究各向异性情况下各种电磁勘探方法的电磁场响应及其数值计算,然而基于电性各向异性的可控源电磁法不但还远没有走到实用化的程度,而且对于无限长源激发下电性各向异性地层的CSAMT的电磁场响应的计算也还没有见到过报道.本文采用有限元法计算了二维电各向异性介质中无限长线源的频率域电磁响应,并与各向同性的结果进行了对比,这对于深化电性各向异性情况下可控源电磁法的研究和认识,提高CSAMT方法的勘探效果和应用水平具有重要的指导意义.

2 CSAMT二维电性各向异性正演理论及有限元法求解

在笛卡尔坐标系中,设无限长线源和二维各向异性电性异常体的走向均平行于x轴,z轴垂直于xy平面向下.围岩背景模型为各向同性介质,二维体的各向异性参数设置为电性主轴的电导率值σx′=σz′>σy′,且只考虑测量坐标轴x轴与电性主轴x′轴间的旋转夹角θ(如图 1所示,假设测量坐标xy平面与电性主轴x′y′平面均平行于水平面),因此张量电导率σ可表示为

(1)

图 1 旋转夹角示意图 (x′y′平面绕逆时针方向旋转) Fig. 1 Schematic diagram of the included angle of rotation (x′y′ plane rotates counterclockwise)

带源的谐变电磁场的Maxwell方程组为

(2)

其中,场量取负时谐因子e-iωtI为单位矩阵,E为电场强度 (单位为V·m-1),B为磁感应强度 (单位为Wb·m-2),D为电位移矢量 (单位为C·m-2),H为磁场强度 (单位为A·m-1),JS为源电流密度 (单位为A·m-2),ω为角频率,μ为磁导率,ε为介电常数.

由 (1)、(2) 式可得 (考虑到构造走向与x轴平行,故∂/∂x=0):

(3)

(4)

其中,σ11=σx′cos2(θ)+σy′sin2(θ)-iωε

   σ22=σy′cos2(θ)+σx′sin2(θ)-iωε

   σ12=σ21=sin (θ) cos (θ)(σy′-σx′).

由 (3) 和 (4) 两式可知,若ExHx为已知,则EyEzHyHz均可求得,故可将其简化为只含ExHx的方程,公式为

(5)

(6)

这里采用有限元方法对方程进行求解.在处理边界条件时,采用第三类边界条件且外边界取得足够远,并利用伽辽金方法来推导以上两式的有限元方程,为此将整个计算区域划分成若干矩形小单元,在每个网格中进行双线性插值,并对单元的加权积分求和,公式为

(7)

(8)

以上两式即为各向异性条件下二维CSAMT的最终计算公式,其中Ω为单元区域,Ne为单元总数,Nm为边界单元总数,Nj(j=1, 2, 3, 4) 为单元插值函数,Nc(c=1, 2) 为边界单元插值函数,β为源到边界上任一点的矢径与边界法线方向的夹角.

为求解 (7) 和 (8) 两式,将其进行单元分析 (未考虑边界条件) 并得到如下线性方程为

(9)

其中,

Q是全为零的列向量,P为含源的列向量.对整个区域里每个单元,都可得到形如 (9) 式的方程,把这些方程组合成一个大型的矩阵方程,最后在将边界条件引入,即可得到最终求解的线性方程.

在场源的处理上,本文通过伪δ函数对源进行加载,无限长线源的表达式为

(10)

其中,I为电流,δi(i)(i=y, z) 为伪δ函数.

通过求解线性方程组,得到ExHx的节点值后,将其代入下面两式 (11) 和 (12) 式即可解得EyHy的节点值,公式为

(11)

(12)

其中求偏导时采用徐世浙 (1994)的方法,最后用求得的场值算出视电阻率和阻抗相位.区域网格剖分时,采用非等间距网格剖分方式,在源处和观测区域网格加密,其他计算区域 (含空气、边界) 网格适当稀疏.

3 模型计算

对于二维无限长线源CSAMT,计算发现TE模式的结果是不可靠的,因此下面只对TM模式进行讨论,这等同于标量CSAMT的结果.结果曲线均为异常体正上方同一测点的计算结果.

3.1 均匀各向同性半空间中各向异性二维异常体模型的计算结果 3.1.1 相同各向异性系数、不同旋转角时,顶面埋深400 m的二维异常体的计算结果

模型设置为:如图 2所示,在均匀半空间中嵌入一方形柱体,其顶面埋深400 m,yz方向的规模分别为500 m×200 m且沿x轴无限延伸,电性主轴方向电阻率的取值为ρx′=ρz′=10 Ωm和ρy′=20 Ωm,旋转角取值分别为θ=0°、30°、45°、60°、90°.均匀半空间的电阻率为200 Ωm.图 2为模型示意图,图 3为埋深400 m时的相同各向异性系数、不同旋转角二维各向异性介质模型的视电阻率和阻抗相位曲线.图中横坐标为频率 (单位为Hz),纵坐标为视电阻率或阻抗相位 (单位为°).同时图中还与等效各向同性 (ρ=14.4 Ωm) 条件下的结果进行了比较.其中,各向异性系数定义为y′方向电阻率值与x′方向电阻率值之比的平方根.

图 2 模型示意图 Fig. 2 Schematic diagram of the model
图 3 相同各向异性系数、不同旋转角二维各向异性介质顶面埋深400m的视电阻率 (a) 与阻抗相位 (b) 曲线 Fig. 3 Apparent resistivity (a) and phase (b) curves of a 2D anisotropic model berried 400m beneath the surface with the same anisotropy coefficients and different rotation angles

图 3中可以看出各向异性的影响是较大的.图 3中的视电阻率曲线 (图 3a) 在高频段的值接近于200 Ωm,而随着频率的降低其值逐渐偏离,最后以接近45°角直线上升且最终趋于重合,相当于进入到了CSAMT的近区.从相位曲线 (图 3b) 上看,在高频段各条曲线的值都接近-45°,并随着频率的降低而逐渐分开,在15 Hz左右各条曲线相交,过交点后随着频率的降低曲线略微分开且角度越小相位越大,并最终趋于重合.在探测到异常体的频段内,角度为0°时的值最小 (此时等同于电阻率ρ=ρx′的各向同性介质),而随着角度的增大,其结果逐渐增大,当角度为90°时值最大 (此时等同于电阻率ρ=ρy′的各向同性介质).等效各向同性的结果介于30°与45°角的结果之间.

3.1.2 相同旋转角、不同各向异性系数时,顶面埋深400 m的二维异常体的计算结果

模型设置为:如图 2所示,且仅将图 2中异常体处电性主轴方向的电阻率值改为ρx′=ρz′=10 Ωm、ρy′=10、20、30、40 Ωm,旋转角取为θ=30°,其余参数均不变.图 4为该模型的视电阻率与阻抗相位曲线.

图 4 相同旋转角、不同各向异性系数二维各向异性介质顶面埋深400 m的视电阻率 (a) 与阻抗相位 (b) 曲线 Fig. 4 Apparent resistivity (a) and phase (b) curves of a 2D anisotropic model berried 400 m beneath the surface with the same rotation angle and different anisotropy coefficients

图 4a可知,在高频时视电阻率值都在200 Ωm附近,后随着频率的降低逐渐分开,且在探测到异常体的频段内各向异性系数越大视电阻率值越大,最后以接近于45°角直线上升并趋于重合.从相位曲线 (图 4b) 看,在高频段先都趋于-45°,后随频率降低而逐渐分开且各向异性系数越大相位越大,在15 Hz左右曲线相交,后随频率降低而分开且各向异性系数越大其值越小,最终趋于重合.可见各向异性的影响是很明显的,不能忽略.

3.2 各向同性层状介质中各向异性二维异常体模型的计算结果 3.2.1 相同各向异性系数、不同旋转角时,顶面埋深600 m的二维异常体的计算结果

模型设置为:如图 5所示,在三层层状介质 (K型) 中嵌入一方形柱体,其中第一层厚400 m、电阻率为80 Ωm,第二层厚600 m、电阻率为600 Ωm,第三层电阻率为80 Ωm,二维体的顶面埋深为600 m,yz方向的规模分别为500 m×200 m且沿x轴无限延伸.异常体在电性主轴方向电阻率的取值为ρx′=ρz′=5 Ωm、ρy′=15 Ωm,旋转角取值分别为θ=0°、30°、45°、60°和90°.图 5为模型示意图,图 6为该模型的视电阻率与阻抗相位曲线.同时图中还与等效各向同性 (ρ=8.7 Ωm) 条件下的结果进行了比较.

图 5 模型示意图 Fig. 5 Schematic diagram of the model
图 6 相同各向异性系数、不同旋转角二维各向异性介质顶面埋深600 m的视电阻率 (a) 与阻抗相位 (b) 曲线 Fig. 6 Apparent resistivity (a) and phase (b) curves of a 2D anisotropic model berried 600 m beneath the surface with the same anisotropy coefficient and different rotation angles

图 6中可明显的看出结果曲线的走势与设置的地电模型基本相同.图 6a中视电阻率在异常体处随角度的增大,视电阻率值越大,且对各向异性的影响表现得很显著,曲线在高频段基本重合其值接近于80 Ωm,在进入近区后以接近45°角直线上升.而相位曲线也表现得很明显.等效各向同性的结果介于30°与45°角的结果之间.

3.2.2 相同旋转角、不同各向异性系数时,顶面埋深600m的二维异常体的计算结果

模型设置如图 5所示,且仅将图 5中二维异常体处电性主轴方向的电阻率值取为ρx′=ρz′=5 Ωm、ρy′=5、10、15、20 Ωm,旋转角取为θ=30°,其余参数均不变.图 7为该模型的视电阻率与阻抗相位曲线.

图 7 相同旋转角、不同各向异性系数二维各向异性介质顶面埋深600 m的视电阻率 (a) 与阻抗相位 (b) 曲线 Fig. 7 Apparent resistivity (a) and phase (b) curves of a 2D anisotropic model berried 600 m beneath the surface with the same rotation angle and different anisotropy coefficients
图 8 在一定各向异性系数情况下、改变旋转角θ时对各向异性的影响分析 Fig. 8 Sketch of effects of changing the rotation angle θ on anisotropy for certain anisotropic coefficients

图 7可见,结果曲线的走势与设置的地电模型基本相同.视电阻率曲线 (图 7a) 在异常体处的计算频段内对各向异性异常体表现明显且随各向异性系数的增大相应的计算结果也增大;在高频段其值都在80 Ωm附近,进入近区后曲线以接近45°角直线上升并趋于重合.相位曲线 (图 7b) 在高频段基本重合其值接近-45°,在进入异常体处的计算频段后各条曲线逐渐降低且各向异性系数越大相位越大,在6 Hz左右曲线相交,此后略微分开且各向异性系数越大相位越小,最终趋于重合.计算结果表明了各向异性的影响是不可忽略的.

3.3 各向异性的探测精度

在本文的计算中发现旋转角越大,各向异性对TM模式的影响越大.为了考察各向异性情况对模型计算结果的影响,这里分别把以上两类模型在不同各向异性系数下、旋转角θ从小到大改变 (0°≤θ≤90°) 时地表处电场的计算结果与各向同性结果之差归纳到表 1表 2中,并以加拿大V8仪器可探测的精度 (其电场强度的探测精度为10-9V·m-1) 进行了可探测性分析 (本文计算中未考虑误差因素).

表 1 均匀半空间背景的二维模型在不同各向异性系数下受到各向异性影响的角度分析 Table 1 Angles affected by anisotropy with different anisotropic coefficients for two-dimensional model embedded in uniform half space
表 2 一维层状介质背景的二维模型在不同各向异性系数下受到各向异性影响的角度分析 Table 2 Angles affected by anisotropy under different anisotropic coefficients for two-dimensional model embedded in a layering earth medium

表 1表 2中可见,在一定各向异性系数下,各向异性角度达到一定数值后,才能被仪器探测到,此时各向异性的影响开始出现.在能探测到各向异性的角度范围内,在临界角附近时能探测到的频率数较少,且随着角度的增大可探测频率数增加、各向异性的影响程度增大,当角度增大到一定值后可探测频率数增加不明显、但各向异性的影响程度继续增大.

4 结论

本文推导了二维各向异性地层条件下无限长线源CSAMT方法的电磁场表达式,并用有限元数值模拟方法进行了实现.在推导出的线源二维各向异性电磁场理论的基础上,设计了一系列地电模型,分别研究了不同各向异性角度和不同各向异性系数情况下的二维各向异性地层的电磁场响应,并与相应的各向同性的情况进行了比对,得到了以下几点认识:

(1) 在区域剖分时,根据发射源和异常体位置的不同采用了非均匀的网格剖分方法,即在发射源和异常体附近网格加密、远处逐渐稀疏的剖分方式,计算表明比单一的均匀网格剖分方式在精度不变的前提下,既节省了内存又提高了效率.

(2) 对于二维电性各向异性模型,当异常体表现为各向异性情况时,无论在视电阻率还是阻抗相位上,均与各向同性的情况存在差异.

(3) 计算结果表明,在各向异性情况下,除非各向异性参数很小,一般不能用各向同性的情形来近似,各向异性的影响程度通常不可忽略.

参考文献
Cai H Z, Xiong B, Zhdanov M. 2015. Three-dimensional marine controlled-source electromagnetic modelling in anisotropic medium using finite element method. Chinese J. Geophys., 58(8): 2839-2850. DOI:10.6038/cjg20150818
Chen G B. 2009. Integral equation method for 3-D electromagnetic modeling in layered anisotropic earth and its applications (in Chinese)[Ph. D]. Changchun: Jilin University.
Chen G B, Wang H N, Yao J J, et al. 2009. Modeling of electromagnetic responses of 3-D electrical anomalous body in a layered anisotropic earth using integral equations. Chinese J. Geophys., 52(8): 2174-2181. DOI:10.3969/j.issn.0001-5733.2009.08.028
Chen X B, Hu W B. 2002. Direct iterative finite element (DIFE) algorithm and its application to electromagnetic response modeling of line current source. Chinese J. Geophys., 45(1): 119-130. DOI:10.3321/j.issn:0001-5733.2002.01.015
Di Q Y, Unsworth M, Wang M Y. 2004. 2.5 D CSAMT modeling with finite element method. Progress in Geophysics, 19(2): 317-324. DOI:10.3969/j.issn.1004-2903.2004.02.018
Di Q Y, Wang G J, Wang M Y, et al. 2009. Frequency response characteristics of target electric property with long bipole large power control source. Chinese J. Geophys., 52(1): 275-280.
Hill D G. 1972. A laboratory investigation of electrical anisotropy in Precambrian rocks. Geophysics, 37(6): 1022-1038. DOI:10.1190/1.1440311
Hu X Y, Huo G P, Gao R, et al. 2013. The magnetotelluric anisotropic two-dimensional simulation and case analysis. Chinese J. Geophys., 56(12): 4268-4277. DOI:10.6038/cjg20131229
Huang Q H, Li Z H, Wang Y B. 2010. A parallel 3-D staggered grid pseudospectral time domain method for ground-penetrating radar wave simulation. J. Geophys. Res., 115(B12): B12101. DOI:10.1029/2010JB007711
Huo G P, Hu X Y, Huang Y F, et al. 2015. MT modeling for two-dimensional anisotropic conductivity structure with topography and examples of comparative analyses. Chinese J. Geophys., 58(12): 4696-4708. DOI:10.6038/cjg20151230
Keller G V. 1968. Electrical Prospecting for Oil, (Quarterly of the Colorado School of Mines, v.63, no.2). Colorado: Colorado School of Mines.
Li X B, Oskooi B, Pedersen L B. 2000. Inversion of controlled-source tensor magnetotelluric data for a layered earth with azimuthal anisotropy. Geophysics, 65(2): 452-464. DOI:10.1190/1.1444739
Li X B, Pedersen L B. 1991. The electromagnetic response of an azimuthally anisotropic half-space. Geophysics, 56(9): 1462-1473. DOI:10.1190/1.1443166
Li X B, Pedersen L B. 1992. Controlled-source tensor magnetotelluric responses of a layered earth with azimuthal anisotropy. Geophys. J. Int., 111(1): 91-103. DOI:10.1111/gji.1992.111.issue-1
Li Y G. 2002. A finite-element algorithm for electromagnetic induction in two-dimensional anisotropic conductivity structures. Geophys. J. Int., 148(3): 389-401. DOI:10.1046/j.1365-246x.2002.01570.x
Li Y G, Pek J. 2008. Adaptive finite element modeling of two-dimensional magnetotelluric fields in general anisotropic media. Geophys. J. Int., 175(3): 942-954. DOI:10.1111/gji.2008.175.issue-3
Li Z H, Huang Q H, Wang Y B. 2009. A 3-D staggered grid PSTD method for borehole radar simulations in dispersive media. Chinese J. Geophys., 52(7): 1915-1922. DOI:10.3969/j.issn.0001-5733.2009.07.027
Løseth L O, Ursin B. 2007. Electromagnetic fields in planarly layered anisotropic media. Geophys. J. Int., 170(1): 44-80. DOI:10.1111/gji.2007.170.issue-1
Luo M, Li Y G. 2015. Effects of the electric anisotropy on marine controlled-source electromagnetic responses. Chinese J. Geophys., 58(8): 2851-2861. DOI:10.6038/cjg20150819
O'Brien D P, Morrison H F. 1967. Electromagnetic fields in an N-layer anisotropic half-space. Geophysics, 32(4): 668-677. DOI:10.1190/1.1439882
Parkhomenko E I. Electrical Properties of Rocks. New York: Plenum Press, 1967.
Pek J, Santos F A M. 2002. Magnetotelluric impedances and parametric sensitivities for 1-D anisotropic layered media. Computers & Geosciences, 28(8): 939-950.
Pek J, Verner T. 1997. Finite-difference modelling of magnetotelluric fields in two-dimensional anisotropic media. Geophys. J. Int., 128(3): 505-521. DOI:10.1111/gji.1997.128.issue-3
Puzyrev V, Koldan J, de la Puente J, et al. 2013. A parallel finite-element method for three-dimensional controlled-source electromagnetic forward modelling. Geophys. J. Int., 193(2): 678-693. DOI:10.1093/gji/ggt027
Ramananjaona C, MacGregor L, Andréis D. 2011. Sensitivity and inversion of marine electromagnetic data in a vertically anisotropic stratified earth. Geophysical Prospecting, 59(2): 341-360. DOI:10.1111/gpr.2011.59.issue-2
Reddy I K, Rankin D. 1971. Magnetotelluric effect of dipping anisotropies. Geophysical Prospecting, 19(1): 84-97. DOI:10.1111/gpr.1971.19.issue-1
Reddy I K, Rankin D. 1975. Magnetotelluric response of laterally inhomogeneous and anisotropic media. Geophysics, 40(6): 1035-1045. DOI:10.1190/1.1440579
Tang X G, Hu W B, Yan L J. 2000. Transient electromagnetic response to a 3-D thin plate in a layered earth. Oil Geophysical Prospecting, 35(5): 628-633,650.
Wang R, Wang M Y, Di Q Y. 2006. Electromagnetic modeling due to line source in frequency domain using finite element method. Chinese J. Geophys., 49(6): 1858-1866. DOI:10.3321/j.issn:0001-5733.2006.06.035
Xiong Z H. 1989. Electromagnetic fields of electric dipoles embedded in a stratified anisotropic earth. Geophysics, 54(12): 1643-1646. DOI:10.1190/1.1442633
Xu S Z. FEM in Geophysics. Beijing: Science Press, 1994.
Xue G Q, Li X, Guo W B, et al. 2007. Characters of response of large-loop transient electro-magnetic field. Oil Geophysical Prospecting, 42(5): 586-590.
Yan S, Chen M S, Fu J M. 2002. Direct time-domain numerical analysis of transient electromagnetic fields. Chinese J. Geophys., 45(2): 275-284. DOI:10.3321/j.issn:0001-5733.2002.02.014
Yin C C, Ben F, Liu Y H, et al. 2014. MCSEM 3D modeling for arbitrarily anisotropic media. Chinese J. Geophys., 57(12): 4110-4122. DOI:10.6038/cjg20141222
Yue R Y, Xu Y X. 2004. Study on electric azimuth anisotropy of controlled source in frequency domain and near-field effects. Oil Geophysical Prospecting, 39(3): 342-347.
Zhou J M, Wang H N, Yao J J, et al. 2012. Parameterized inversion of frequency sounding data in horizontally stratified transverse isotropic medium. Acta Physica Sinica, 61(8): 089101.
蔡红柱, 熊彬, ZhdanovM. 2015. 电导率各向异性的海洋电磁三维有限单元法正演. 地球物理学报, 58(8): 2839–2850. DOI:10.6038/cjg20150818
陈桂波. 2009. 各向异性地层中电磁场三维数值模拟的积分方程算法及其应用[博士论文]. 长春: 吉林大学.
陈桂波, 汪宏年, 姚敬金, 等. 2009. 用积分方程法模拟各向异性地层中三维电性异常体的电磁响应. 地球物理学报, 52(8): 2174–2181. DOI:10.3969/j.issn.0001-5733.2009.08.028
陈小斌, 胡文宝. 2002. 有限元直接迭代算法及其在线源频率域电磁响应计算中的应用. 地球物理学报, 45(1): 119–130. DOI:10.3321/j.issn:0001-5733.2002.01.015
底青云, UnsworthM, 王妙月. 2004. 有限元法2.5维CSAMT数值模拟. 地球理学进展, 19(2): 317–324. DOI:10.3969/j.issn.1004-2903.2004.02.018
底青云, 王光杰, 王妙月, 等. 2009. 长偶极大功率可控源激励下目标体电性参数的频率响应. 地球物理学报, 52(1): 275–280.
胡祥云, 霍光谱, 高锐, 等. 2013. 大地电磁各向异性二维模拟及实例分析. 地球物理学报, 56(12): 4268–4277. DOI:10.6038/cjg20131229
霍光谱, 胡祥云, 黄一凡, 等. 2015. 带地形的大地电磁各向异性二维模拟及实例对比分析. 地球物理学报, 58(12): 4696–4708. DOI:10.6038/cjg20151230
李展辉, 黄清华, 王彦宾. 2009. 三维错格时域伪谱法在频散介质井中雷达模拟中的应用. 地球物理学报, 52(7): 1915–1922. DOI:10.3969/j.issn.0001-5733.2009.07.027
罗鸣, 李予国. 2015. 一维电阻率各向异性对海洋可控源电磁响应的影响研究. 地球物理学报, 58(8): 2851–2861. DOI:10.6038/cjg20150819
唐新功, 胡文宝, 严良俊. 2000. 层状地层中三维薄板的瞬变电磁响应. 石油地球物理勘探, 35(5): 628–633,650.
王若, 王妙月, 底青云. 2006. 频率域线源大地电磁法有限元正演模拟. 地球物理学报, 49(6): 1858–1866. DOI:10.3321/j.issn:0001-5733.2006.06.035
徐世浙. 地球物理中的有限单元法. 北京: 科学出版社, 1994.
薛国强, 李貅, 郭文波, 等. 2007. 大回线源瞬变电磁场响应特性. 石油地球物理勘探, 42(5): 586–590.
闫述, 陈明生, 傅君眉. 2002. 瞬变电磁场的直接时域数值分析. 地球物理学报, 45(2): 275–284. DOI:10.3321/j.issn:0001-5733.2002.02.014
殷长春, 贲放, 刘云鹤, 等. 2014. 三维任意各向异性介质中海洋可控源电磁法正演研究. 地球物理学报, 57(12): 4110–4122. DOI:10.6038/cjg20141222
岳瑞永, 徐义贤. 2004. 可控源频率域电方位各向异性与近场效应研究. 石油地球物理勘探, 39(3): 342–347.
周建美, 汪宏年, 姚敬金, 等. 2012. 水平层状横向同性地层中频率测深资料全参数快速迭代反演. 物理学报, 61(8): 089101.