大地电磁三维矢量有限元正演模拟

石明+冯德山+李开鹏+王珣

摘 要：从Maxwell方程出发，开展了三维大地电磁场所满足的边值问题研究，利用加权余量法导出了三维大地电磁有限元方程.介绍了三维矢量有限元六面体网格剖分方式、插值基函数选取，推导了三维大地电磁矢量有限元正演的单元刚度系数矩阵及离散格式.编制了三维矢量有限元大地电磁正演的Matlab程序.三维COMMEMI 3D-1模型的视电阻率曲线与国际通用的标准测试数据能很好地拟合，验证了作者编写的矢量有限元正演程序的正确性.通过对高、低阻异常体的阻抗张量形态分析，说明张量阻抗等值线图能用以大致判断异常体特性，丰富了大地电磁响应特征的表达方式.

关键词：矢量有限元；大地电磁；正演模拟；张量阻抗

中图分类号：P631 文献标识码：A

文章编号：1674-2974（2016）10-0119-07

Abstract：Starting from the Maxwell equations， this article studied the boundary conditions of 3D MT. By using the weighted residual method， we derived the three-dimensional MT finite element equation. The three-dimensional vector finite element hexahedral meshing mode was introduced and the basis functions were selected. Then we derived the three-dimensional magnetotelluric vector finite element stiffness coefficient matrix and discrete format. A three-dimensional vector finite element magnetotelluric forward Matlab program was done. The apparent resistivity curve of the dimensional COMMEMI 3D-1 model matches the international standard test data， which proves the correctness of 3D magnetotelluric forward program. With the analysis of high and low resistivity anomalies， it shows that tensor impedance map can roughly determine the anomaly characteristics， which enriches the magnetotelluric response characteristics of expression.

Key words：vector finite element method；magnetotelluric； forward modeling； impedance tensor

大地电磁（MT）是以电离层激发的天然交变电磁场为场源，在地表观测相互正交的电场、磁场分量来获取地电构造信息的一种重要地球物理勘探方法[1].MT不需要庞大的发射源设备，只需采用比较轻便的接收设备，野外工作方便、成本低，被广泛应用于地壳和上地幔电性结构的研究，在石油天然气勘探、矿产资源勘探、工程与环境普查等领域，发挥着举足轻重的作用[2-9].可以预见，三维MT勘探技术是地球物理中深层领域的研究热点及今后MT的发展趋势，而三维MT正演是理解MT勘探物理现象并认识地质体电磁响应规律的有效手段，显然尤其重要.

尽管矢量FEM拥有诸多优点，但在地球物理的电磁法正演领域中，其应用并不多见，尚需要进一步完善.目前的研究主要包括：Yoshimura 等 [10]开展了矢量FEM的MT响应数值模拟，并将矢量FEM的计算结果与交错网格FDM的计算结果进行了对比；Mitsuhata等 [11]利用矢量FEM和节点FEM耦合的方法对三维MT数值模拟；Nam [12]采用不规则六面体矢量FEM直接计算电场，研究了起伏地形下MT的电阻率和相位的变化规律；刘长生等[13]将完全非结构化四面体单元引入到矢量有限元中，实现了三维大地电磁h-型自适应矢量有限元正演；王烨[14]开展了高频率大地电磁法矢量有限元正演，并采用改进的威尔金森方法求解大型病态方程组，提高了迭代速度；顾观文等[15]开展了矢量有限元法MT三维地形数值模拟，研究了地形起伏下三维阻抗张量的变化规律；杨军等[16]采用非结构四面体单元的三维矢量FEM实现了海洋可控源电磁数值模拟；苏晓波等[17]采用规则六面体单元的三维矢量FEM实现了大地电磁数值模拟，并对网格剖分的重要性进行了研究.

在前人基础上，作者推导了三维大地电磁矢量FEM正演的离散形式，应用矢量FEM算法计算了三维COMMEMI 3D-1国际模型[18]的MT视电阻率及模型张量阻抗，研究了高低阻异常体的电磁响应特性，有效地指导了MT的资料解释.

图5（a），（b）分别为f=0.1 Hz时XY模式与YX模式下矢量FEM正演视电阻率曲线.分析图5（a），（b）可知，两幅图中的矢量FEM曲线与COMMEMI所提供的数据都能够很好地吻合，说明无论是在低频还是高频部分，应用矢量FEM开展三维大地电磁正演，都具有较高的精度，同时也验证了矢量FEM算法及程序的正确性.

图6为应用矢量FEM正演计算COMMEMI3D-1模型得到的张量阻抗.由图可见，10 Hz与0.1 Hz两个频率下的张量阻抗形态基本一致，10 Hz的数值较0.1 Hz要大.对比图中4个不同的张量阻抗，可以发现，图6（a），（d），（e）和（h）中两个频率下的Zxx与Zyy分为四瓣，且阻抗值较小，四瓣的中心反映了异常体的边界，而图6（b），（c），（f）和（g）中的Zxy与Zyx阻抗值较大，反映了入射场的特性.根据张量阻抗理论可知，当构造为二维构造时，Zxx和Zyy为零，即当异常体走向方向越长，Zxy与Zyx越小，Zxy与Zyx差异也越大.由此，张量阻抗分解后，无需做反演即可以判断出异常体的简单特性.

为了进一步认识大地电磁的响应特性，对比高、低阻异常体张量阻抗的不同，在图3中COMMEMI3D-1测试模型的基础上，仅将低阻异常体改为1 000 Ω·m高阻异常体.其他参数均与国际模型相同.应用三维矢量FEM开展三维高阻异常体模型的张量阻抗研究.

图7为应用三维矢量FEM正演的10 Hz大地电磁张量阻抗图.分析图7（a）与图7（d）可知，高阻异常体张量阻抗中的Zxx与Zyy同样分为四瓣，且阻抗值较小，其四瓣的中心反映了异常体的边界.由于异常体x方向与y方向的比值为1∶2，图7（b）中的Zxy与图7（c）中的Zyx差异较大.对比高低阻异常10 Hz时的阻抗相位Zxx，虽然两者都为四瓣，但是阻抗值正负值的分布正好相反，低阻异常体四瓣的中心向外辐射，幅值变小趋于0；而高阻异常体四瓣的中心向外辐射，幅值变小趋于0之后会发生反转之后再次趋于0.Zyy具有相同的规律.对比Zxy和Zyx，高阻异常体中心仅出现一个闭合异常形态，而低阻异常体则形态更为复杂.

4 结 论

1） 介绍了三维矢量有限元区域剖分方式，对矢量FEM插值基函数以及单元插值方式进行了阐述，应用Galerkin算法，推导了三维矢量FEM大地电磁方程离散格式，编制了矢量FEM三维MT的Matlab模拟程序.

2） 设置三维COMMEMI 3D-1模型进行矢量FEM的计算，模拟结果与COMMEMI提供的数据拟合效果很好，验证了矢量有限元程序的正确性.通过对比高低阻异常体的张量阻抗，分析了不同异常下张量阻抗的特点，进一步认识了MT的响应特性.

参考文献

[1] 柳建新，童孝忠，郭荣文，等. 大地电磁测深勘探：资料处理反演与解释[M]. 北京：科学出版社，2012：1-12.

LIU jian-xin， TONG Xiao-zhong， GUO Rong-wen， et al. Magnetotelluric sounding exploration： data processing inversion and interpretation [M].Beijing： Science Press， 2012：1-12.（In Chinese）

[2] 底青云，王若. 可控源音频大地电磁数据正反演及方法应用[M]. 北京：科学出版社，2008：1-8.

DI Qing-yun， WANG Ruo. Controlled source audio magnetotelluric data inversion and methods application [M].Beijing： Science Press， 2008： 1-8.（In Chinese）

[3] 谭捍东，余钦范，JOHN B，等. 大地电磁法三维交错采样有限差分数值模拟[J]. 地球物理学报，2003，46（5）：705-711.

TAN Han-dong， YU Qin-fan， JOHN B， et al. Magnetotelluric three-dimensional modeling using the staggered-grid finite difference method[J]. Chinese Journal of Geophysics， 2003， 46（5）：705-711.（In Chinese）

[4] 陈辉，邓居智，谭捍东，等. 大地电磁三维交错网格有限差分数值模拟中的散度校正方法研究[J]. 地球物理学报，2011，54（6）：1649-1659.

CHEN Hui， DENG Ju-zhi， TAN Han-dong， et，al. Study on divergence correction method in three-dimensionalmagnetotelluric modeling with staggered-grid finite difference method [J]. Chinese Journal of Geophysics， 2011， 54（6）： 1649-1659. （In Chinese）

[5] 李焱，胡祥云，杨文采，等. 大地电磁三维交错网格有限差分数值模拟的并行计算研究[J]. 地球物理学报，2012，55（12）：4036-4043.

LI Yan， HU Xiang-yun， YANG Wen-cai， et al. A study on parallel computation for 3D magnetotelluric modeling using the staggered-grid finite difference method [J]. Chinese Journal of Geophysics， 2012， 55（12）： 4036-4043. （In Chinese）

[6] 徐凌华，童孝忠，柳建新，等. 基于有限单元法的二维/三维大地电磁正演模拟策略[J]. 物探化探计算技术，2009，31（5）：421-425.

XU Ling-hua， TONG Xiao-zhong， LIU Jian-xin， et al. Solution strategies for 2D and 3D magnetotelluric forward modeling based on the finite element method [J]. Computing Techniques for Geophysical and Geochemical Exploration， 2009， 31 （5）： 421-425. （In Chinese）

[7] MOGI T. Three-dimensional modeling of magnetotelluric data using finite- element method [J]. Journal of Applied Geophysics， 1996， 35（2）： 185-189.

[8] 梁生贤，张胜业，吾守艾力，等. 复杂三维介质的大地电磁正演模拟[J]. 地球物理学进展，2012，27（5）：1981-1988.

LIANG Sheng-xian， ZHANG Sheng-ye， WU Shouaili， et al. Magnetotelluric forward modeling in complex three- dimensional media [J]. Progress in Geophys， 2012， 27（5）： 1981-1988. （In Chinese）

[9] 金建铭.电磁场有限元方法[M]. 西安：西安电子科技大学出版社， 1998：164-177.

JIN Jian-ming. Finite element method of electromagnetic field [M]. Xian： Xian Electronic University Press， 1998：164-177.（In Chinese）

[10]YOSHIMURA R，OSHIMAN N. Edge-based finite element approach to the simulation of geoelectromagnetic induction in a 3-D sphere [J]. Geophysical Research Letters， 2002， 29（3）： 1039-1042.

[11]MITSUHATA Y， UCHIDA T. 3D magnetotelluric modeling using the T-Ω finite-element method [J].Geophysics， 2004， 69（1）：108-119.

[12]NAM M J， KIM H J， SONG Y， et al. 3D magnetotelluric modelling including surface topography [J]. Geophysical Prospecting， 2007， 55（2）： 277-287.

[13]刘长生，汤井田，任政勇，等. 基于非结构化网格的三维大地电磁自适应矢量有限元模拟[J]. 中南大学学报：自然科学版，2010，41（5）：1855-1860.

LIU Chang-sheng， TANG Jing-tian， REN Zheng-yong， et al. Three-dimension magnetotellurics modeling by adaptive edgefinite-element using unstructured meshes [J]. Journal of Central South University： Science and Technology，2010， 41（5）： 1855-1860. （In Chinese）

[14]王烨. 基于矢量有限元的高频大地电磁法三维数值模拟[D]. 长沙：中南大学地球科学与信息物理学院，2008：35-62.

WANG Ye. A study of 3D high frequency magnetotelluric modeling by edge-based finite element method [D]. Changsha： Central South University. School of Geosciences and Info-Physics， 2008：35-62. （In Chinese）

[15]顾观文，吴文鹂，李桐林. 大地电磁场三维地形影响的矢量有限元数值模拟[J]. 吉林大学学报：地球科学版，2014，44（5）：1678-1686.

GU Guan-wen， WU Wen-li， LI Tong-lin. Modeling for the effect of magnetotelluric 3D topography based on the vector finite-element method [J]. Journal of Jilin University： Earth Science Edition， 2014，44（5）：1678-1686. （In Chinese）

[16]杨军，刘颖，吴小平. 海洋可控源电磁三维非结构矢量有限元数值模拟[J]. 地球物理学报，2015，58（8）：2827-2838.

YANG Jun， LIU Ying， WU Xiao-ping. 3D simulation of marine CSEM using vector finite element method on unstructured grids [J]. Chinese Journal of Geophysics， 2015， 58（8）： 2827-2838. （In Chinese）

[17]苏晓波，李桐林，朱成，等.大地电磁三维矢量有限元正演研究[J].地球物理学进展，2015，30（4）：1772-1778.

SU Xiao-bo， LI Tong-lin， ZHU Cheng， et al. Study of three-dimensional MT forward modeling using vector finite element method [J].Progress in Geophysics， 2015，30（4）：1772-1778. （In Chinese）

[18]ZHDANOV M S， VARENTSOV I M， WEAVER J T， et al. Methods for modelling electromagnetic fields Results from COMMEMI—the international project on the comparison of modelling methods for electromagnetic induction [J]. Journal of Applied Geophysics， 1997， 37（3/4）：133-271.