Xushan Lu和Colin G. Farquharson来自加拿大纽芬兰纪念大学地球科学系的研究团队，于2020年7月在期刊《Geophysics》发表了题为《3D finite-volume time-domain modeling of geophysical electromagnetic data on unstructured grids using potentials》的原创研究论文。该研究针对地球物理电磁勘探中的三维时域正演建模问题，提出了一种基于非结构化网格的有限体积法（finite-volume, FV）解决方案，重点探讨了利用磁矢量势（A）和电标量势（φ）的建模方法及其物理意义。

学术背景

该研究属于电磁法勘探（electromagnetic methods）领域，旨在解决复杂地质模型（如真实地形和矿体形态）的电磁响应模拟问题。传统结构化网格（如矩形或八叉树网格）难以精确表征复杂几何结构，而非结构化网格（unstructured grids）可通过局部加密灵活适应地质界面和发射源/接收器区域。尽管有限元法（finite-element, FE）在非结构化网格中应用广泛，但有限体积法具有矩阵稀疏性更高、计算效率更优的特点（Jahandari et al., 2017）。此外，时域电磁（TDEM）数据建模通常采用电场（E-field）方法或势函数（A-φ）方法，后者能直接计算磁场（B-field）且更易通过迭代法求解，但需规范势函数的唯一性（如库仑规范条件）。本研究的目标是开发基于非结构化Delaunay-Voronoi对偶网格的FV时域算法，并分析A-φ分解对电磁场“感应”（inductive）和“ galvanic”（galvanic）效应的物理阐释价值。

研究流程

1. 控制方程构建
   研究从准静态近似下的麦克斯韦方程组出发，推导了电场亥姆霍兹方程（Helmholtz equation）和基于A-φ的耦合方程（式6-7）。其中，A-φ方法需引入库仑规范条件（∇·A=0）以保证解的唯一性，并通过隐式（式8）或显式（式10a-c）两种方式实现规范约束。隐式方法通过添加−∇∇·A项修改控制方程，而显式方法将规范条件作为额外方程引入系统。

2. 空间离散化
   采用Delaunay-Voronoi对偶网格，利用TetGen生成四面体网格。关键操作包括：

   • 离散算子设计：通过面积分和体积分将连续算子（如∇×、∇·）转化为离散矩阵（如curl、div）。例如，磁矢量势A定义于Delaunay边，其旋度算子通过Voronoi面的面积分近似（式13-15）。
     
   • 电导率平均：提出体积加权算法（volume-weighted method，式26）将四面体单元电导率映射至Voronoi面，避免Delaunay退化情况下的计算困难。
     

3. 时间离散化
   采用后向欧拉法（backward Euler）进行时间步进，时间步长按分段策略动态调整：早期采用小步长（如10⁻⁸ s）捕捉快速变化的瞬态场，后期逐步增大步长以提高效率。每个时间步的线性系统通过直接求解器MUMPS处理，利用矩阵稀疏性优化计算。

4. 初始条件计算
   对于回线源（loop source），通过解析解（式29）求解静磁问题获得初始A；对于接地线源（grounded-wire source），需额外求解直流电阻率（DCR）方程（式31）获得初始φ。

5. 验证与对比

   • 均质半空间模型：对比A-φ方法与解析解的B-field响应，验证隐式/显式规范条件的精度（图5）。结果显示，A-φ方法直接计算的磁场与解析解吻合，而E-field方法通过数值积分得到的磁场在晚期存在显著误差（因未考虑导体响应的外推）。
     
   • 三维导电块体模型：对比FV与FE、积分方程（IE）法的dB/dt响应（图6），三者一致性验证了算法的正确性。计算效率分析表明，E-field FV方法耗时最短，但A-φ方法因能直接计算磁场而更适合B-field传感器数据解释。

6. 电磁场分解分析
   通过水平/垂直薄板导体模型（图7,12），可视化A-φ分解的电场分量：

   • 感应部分（−∂A/∂t）：主导背景场中的涡流行为，表现为环绕导体的环形电流（图8e, 14d）。
     
   • Galvanic部分（−∇φ）：源于导体-围岩界面的电荷积累，在导体边缘形成类直流场模式（图8f, 14c）。研究表明，即使以感应源为主的场景，galvanic效应仍通过边界电荷影响总场分布（图10g, 16g）。

主要结果

1. 算法性能

   • 隐式规范方法内存消耗为E-field方法的3倍，但计算时间仅增加80%（表1）。显式规范在精度相近的前提下内存略优。
     
   • 对于包含193,512边的网格，A-φ FV与FE方法的db/dt响应相对误差%（图6b），验证了非结构化网格FV法的可靠性。

2. 物理机制

   • 水平薄板模型中，早期（t=0.018 ms）感应与galvanic部分在导体内近乎抵消，导致弱电场（图8a）；晚期（t=31.623 ms）感应部分主导，但galvanic效应在边界显著（图8i）。
     
   • 非对称源配置下（图10），galvanic部分通过电荷分布调控总场，形成导体内局部涡流（图16g），表明“感应”响应实为两者协同作用的结果。

结论与价值

本研究提出的非结构化网格FV时域算法，为复杂地质模型的电磁响应模拟提供了高效工具。其科学价值在于：
1. 通过A-φ分解量化了galvanic效应对TDEM数据的影响，修正了传统“纯感应”假设的局限性；
2. 直接计算B-field的能力避免了数值积分的外推误差，适用于深部导体勘探。应用层面，该方法可优化矿产（如块状硫化物、石墨）和油气资源的电磁探测方案设计。

研究亮点

1. 方法创新：首次将FV法扩展至非结构化网格的时域A-φ建模，支持隐式/显式规范约束。
   
2. 物理洞察：通过矢量场可视化揭示了galvanic效应在瞬态电磁中的关键作用，为数据解释提供新视角。
   
3. 工程优化：动态时间步长策略与MUMPS求解器的结合，平衡了计算精度与效率。
   

此外，研究还开源了网格生成（TetGen）和线性求解（MUMPS）的集成方案，为后续开发奠定了基础。