A Massively Parallel Three-Grid Preconditioner for the High-Frequency Helmholtz Equation

Abstract

Accurate simulation of three-dimensional time-harmonic wave propagation over many wavelengths requires control of phase error and efficient solution of large indefinite systems. We develop a three-grid solver based on the compact 27-point interpolated optimized finite-difference (IOFD) discretization. Its wavenumber-dependent stencil supports a fine-grid resolution of six points per shortest wavelength and an unshifted physical correction on the \(2h\) grid at only three points per shortest wavelength. The method retains unshifted IOFD operators on the \(h\) and \(2h\) grids, while a complex-shifted \(2h\)--\(4h\) auxiliary cycle preconditions a factorization-free iterative approximation of the coarse inverse. Restricting the shift to this auxiliary cycle preserves the propagative character of the coarse correction. Comparison with the outgoing Green function confirms phase and relative-amplitude accuracy on a sequence of meshes up to \revision{\(10240^3\)}, \revision{exceeding one trillion unknowns,} with the largest problem spanning approximately \revision{1704} wavelengths per coordinate. The same fixed solver configuration retains robust convergence across smooth, discontinuous, high-contrast, and geophysical velocity models and exhibits scalable parallel performance. In particular, a \revision{\(2560^3\)} problem spanning approximately \revision{425} wavelengths in each coordinate direction is solved in \revision{43.8 seconds} on just 64 NVIDIA A100 GPUs.