[Paper Review] A Tutorial and Open Source Software for the Efficient Evaluation of Gravity and Magnetic Kernels

Jarom D. Hogue, Rosemary A. Renaut+1 more|arXiv (Cornell University)|Dec 15, 2019

Geophysical and Geoelectrical Methods16 references29 citations

TL;DR

This paper presents an open-source MATLAB framework for efficiently computing gravity and magnetic forward modeling using 2D Fast Fourier Transforms (FFT) on block-Toeplitz Toeplitz-block (BTTB) sensitivity matrices. By exploiting the BTTB structure—symmetric for gravity and non-symmetric for magnetic fields—the method enables fast matrix-vector multiplication without explicitly forming the full kernel matrices, achieving significant speedups and reduced memory usage, especially for large 3D volumes with variable padding.

ABSTRACT

Fast computation of three-dimensional gravity and magnetic forward models is considered. Measurement data is assumed to be obtained on a uniform grid which is staggered with respect to the discretization of the parameter volume. Then, the resulting kernel sensitivity matrices exhibit block-Toeplitz Toeplitz-block (BTTB) structure. These matrices are symmetric for the gravity problem but non-symmetric for the magnetic problem. In each case, the structure facilitates fast forward computation using two-dimensional fast Fourier transforms. The construction of the kernel matrices and the application of the transform for fast forward multiplication, for each problem, is carefully described. But, for purposes of comparison with the transform approach, the generation of the unique entries that define a given kernel matrix is also explained. It is also demonstrated how the matrices, and hence transforms, are adjusted when padding around the volume domain is introduced. The transform algorithms for fast forward matrix multiplication with the sensitivity matrix and its transpose, without the direct construction of the relevant matrices, are presented. Numerical experiments demonstrate the significant reduction in computation time that is achieved using the transform implementation. Moreover, it becomes feasible, both in terms of reduced memory requirements and computational time, to implement the transform algorithms for large three-dimensional volumes. All presented algorithms, including with variable padding, are coded for optimal memory, storage and computation as an open source MATLAB code which can be adapted for any convolution kernel which generates a BTTB matrix. This work, therefore, provides a general tool for the efficient simulation of gravity and magnetic field data, as well as any formulation which admits a sensitivity matrix with the required structure.

Motivation & Objective

• To develop an efficient, open-source implementation for forward modeling of 3D gravity and magnetic fields using BTTB-structured sensitivity matrices.
• To extend existing FFT-based BTTB methods to handle non-symmetric magnetic kernels and arbitrary padding in x and y directions.
• To provide optimized, memory-efficient algorithms for computing matrix-vector products and their transposes without explicitly constructing large kernel matrices.
• To enable scalable, high-performance inversion workflows by integrating fast forward and transpose operations into geophysical modeling pipelines.

Proposed method

• Leverages the BTTB structure of gravity and magnetic kernel matrices arising from staggered uniform grids in x and y.
• Uses 2D FFTs on block-circulant-circulant-block (BCCB) embeddings to accelerate matrix-vector multiplication.
• Applies symmetric and non-symmetric BTTB algorithms separately for gravity (symmetric) and magnetic (non-symmetric) kernels.
• Implements variable padding in x and y directions via asymmetric boundary extensions, preserving FFT efficiency.
• Derives explicit formulas for computing the forward and transpose operations using FFTs of embedded BCCB matrices.
• Provides a modular, open-source MATLAB implementation with safety checks for memory usage and extensibility to any convolution kernel generating BTTB structure.

Experimental results

Research questions

• RQ1 How can the BTTB structure of gravity and magnetic kernel matrices be exploited to accelerate forward modeling without explicit matrix formation?
• RQ2 What is the performance gain of using 2D FFTs on BCCB-embedded BTTB matrices compared to direct computation for large 3D models?
• RQ3 How can domain padding be incorporated into the BTTB-FFT framework while maintaining computational efficiency and numerical stability?
• RQ4 What are the computational and memory savings when computing matrix-vector products and their transposes using the FFT-based approach versus direct methods?
• RQ5 To what extent can the same framework be generalized to other convolution-type kernels beyond gravity and magnetic fields?

Key findings

• The FFT-based BTTB method reduces computation time by orders of magnitude compared to direct matrix multiplication, especially for large 3D models.
• The method achieves significant memory savings by avoiding explicit construction of large kernel matrices, making large-scale modeling feasible.
• The inclusion of asymmetric padding in x and y directions is supported without loss of FFT efficiency, enabling flexible domain discretization.
• The algorithm for non-symmetric magnetic kernels is successfully implemented with only moderate increases in memory and computation compared to the symmetric gravity case.
• Numerical experiments confirm that the relative mean errors in forward and transpose operations remain below 1e-10 for both gravity and magnetic cases, even with padding.
• The open-source MATLAB code is production-ready, with built-in memory safety checks and modular design for easy adaptation to other BTTB-generating kernels.