A High-order nodal spectral element method for 3D magnetotelluric forward modeling

Summary We present a high-order tetrahedral spectral element (SE) method for the computation of three-dimensional (3D) magnetotelluric (MT) forward responses, designed to overcome the limitations of conventional SE methods that rely on hexahedral grids. Our approach utilizes tetrahedral grids, enabling the accurate simulations of large-scale, geophysically complex models, including intricate subsurface anomalies and irregular topography. Starting from Maxwell’s equations, we derive the governing SE equations using a magnetic vector potential A and an electric scalar potential Φ, incorporating the Coulomb gauge to suppress spurious solutions. The computational domain is discretized using the weighted residual Galerkin method, with Proriol-Koornwinder-Dubiner (PKD) polynomials serving as the weighting and shape functions within each tetrahedral element. Two coordinate transformations-affine and collapse transformations are applied during the solution process. To better leverage the properties of the basis functions, both the interpolation and integration nodes are chosen from the same Warp & Blend point set, rather than using two separate sets, which simplifies the computation of the coefficient matrix terms. The resulting global sparse linear system is solved efficiently using the PARDISO direct solver. We assess the accuracy and computational performance of our method through validation against well-established MT community models. Our evaluation, based on misfit (relative error), degrees of freedom (DOFs), computational time, and memory usage for various polynomial orders, demonstrates that the proposed SE method on tetrahedral grids offers a robust and efficient solution for high-precision forward modeling in MT applications.