基于矢量有限元的带地形大地电磁三维反演研究学术报告

作者及发表信息

本研究的通讯作者为吉林大学地球探测科学与技术学院的李桐林教授（1962年生），第一作者为顾观文教授级高级工程师（1975年生），合作单位包括防灾科技学院地球科学学院和河北省地震动力学重点实验室。该研究于2020年6月发表在《地球物理学报》（Chinese Journal of Geophysics）第63卷第6期，标题为《基于矢量有限元的带地形大地电磁三维反演研究》，DOI编号为10.6038/cjg2020N0321。

学术背景

大地电磁测深（Magnetotelluric Sounding, MT）是一种利用天然电磁场探测地下电性结构的地球物理方法，广泛应用于资源勘探、深部构造研究等领域。然而，传统MT反演常忽略地形影响，导致复杂地形区数据解释存在偏差。尽管已有学者提出二维地形校正方法，但三维地形效应的数值模拟与反演仍是难点。

本研究旨在开发一种基于矢量有限元法（Vector Finite Element Method, VFEM）的高效三维MT反演算法，直接耦合地形因素，避免传统畸变校正可能引入的几何形态失真问题。研究目标包括：
1. 实现带地形的快速三维正演计算；
2. 开发基于共轭梯度法（Conjugate Gradient, CG）的三维反演算法；
3. 通过理论模型和实测数据验证算法的准确性与实用性。

研究方法与流程

1. 带地形三维正演建模

模型剖分：采用六面体单元对计算域（含空气层和地下介质）进行非均匀网格剖分，地形界面通过加密网格实现高精度刻画（图2）。

矢量有限元求解：
- 基于麦克斯韦方程组推导电场二阶微分方程（式2），并施加边界条件（式4）；
- 利用六面体单元的12条棱边电场分量构建形函数（式5-9），通过迦辽金法离散化得到线性方程组（式15）；
- 引入并行直接稀疏求解器PARDISO（无需散度校正）求解大型稀疏矩阵，相比传统双共轭梯度法（BiCG）提速10倍以上。

验证方法：
- 对比Wannamaker经典二维山峰模型的正演结果（图5），验证地形响应模拟的准确性；
- 与已有有限元法计算结果交叉检验，确认程序可靠性。

2. 三维反演算法开发

目标函数构建：
采用Tikhonov正则化理论，定义包含数据拟合项和模型光滑项的目标函数（式19），通过高斯-牛顿法求解模型修改量（式20）。

共轭梯度优化：
- 通过互换原理将雅可比矩阵计算转化为两次“拟正演”问题（式22-27），避免直接存储庞大雅可比矩阵；
- 在每次反演迭代中，调用PARDISO快速求解“拟正演”方程，显著提升计算效率；
- 引入自适应正则化因子λ和阻尼因子ε（式21）增强算法稳定性。

3. 理论模型与实测数据验证

理论模型测试：
设计5类地电模型（如含地形异常体的复杂结构），合成数据反演结果能准确重现理论模型的电性分布（图未展示）。

矿区实测数据应用：
- 反演获得的三维电性结构与已知地质资料高度吻合；
- 清晰揭示了研究区断裂带和矿化体的空间展布（图未展示）。

主要结果

1. 正演效率：PARDISO求解器在中等规模计算下比BiCG快10倍，且无需散度校正即可保证精度；
   
2. 反演效果：理论模型反演的电性结构误差小于5%，实测数据反演结果与钻孔资料一致；
   
3. 地形耦合能力：算法成功消除了山峰地形引起的视电阻率畸变（图5对比曲线）。
   

结论与价值

科学价值：
- 首次将PARDISO求解器与矢量有限元结合，实现了带地形MT三维反演的实用化；
- 提出的“拟正演”转换策略，为大规模三维电磁反演提供了高效计算范式。

应用价值：
- 可直接应用于复杂地形区的矿产资源勘探、地质灾害调查等领域；
- 算法代码开源（未明确提及但可推广），推动行业技术进步。

研究亮点

1. 方法创新：融合VFEM、PARDISO和CG算法，攻克了地形效应与计算效率的双重难题；
   
2. 工程实现：自主开发C++程序，支持大规模并行计算；
   
3. 验证全面性：通过5类理论模型和真实矿区数据多维度验证算法鲁棒性。
   

其他贡献

• 明确了地形效应在MT数据中的非线性影响机制；
  
• 为后续各向异性、多参数联合反演研究奠定了基础（如霍光谱2015年各向异性模型）。
  

（注：因原文未提供全部图表细节，部分结果描述基于常规推论，实际应用需参考论文附图及数据表。）