Algorithm for Acoustic Wavefield in Space-Wavenumber Domain of Vertically Heterogeneous Media Using NUFFT

Balancing efficiency and accuracy is often challenging in the numerical solution of three-dimensional (3D) point source acoustic wave equations for layered media. To overcome this, an efficient solution method in the spatial-wavenumber domain is proposed, utilizing the Non-Uniform Fast Fourier Transform (NUFFT) to achieve arbitrary non-uniform sampling. By performing a two-dimensional (2D) Fourier transform on the 3D acoustic wave equation in the horizontal direction, the 3D equation is transformed into a one-dimensional (1D) space-wavenumber-domain ordinary differential equation, effectively simplifying significant 3D problems into one-dimensional problems and significantly reducing the demand for memory. The one-dimensional finite-element method is applied to solve the boundary value problem, resulting in a pentadiagonal system of equations. The Thomas algorithm then efficiently solves the system, yielding the layered wavefield distribution in the space-wavenumber domain. Finally, the wavefield distribution in the spatial domain is reconstructed through a 2D inverse Fourier transform. The correctness of the algorithm was verified by comparing it with the finite-element method. The analysis of the half-space model shows that this method can accurately calculate the wavefield distribution in the air layer considering the air layer while exhibiting high efficiency and computational stability in ultra-large-scale models. The three-layer medium model test further verified the adaptability and accuracy of the algorithm in calculating the distribution of acoustic waves in layered media. Through a sensitivity analysis, it is shown that the denser the mesh node partitioning, the higher the medium velocity, and the lower the point source frequency, the higher the accuracy of the algorithm. An algorithm efficiency analysis shows that this method has extremely low memory usage and high computational efficiency and can quickly solve large-scale models even on personal computers. Compared with traditional FEM, the algorithm has much higher advantages in terms of memory usage and efficiency. This method provides a new approach to the numerical solution of partial differential equations. It lays an essential foundation for background field calculation in the scattering seismic numerical simulation and full-waveform inversion of acoustic waves, with strong theoretical significance and practical application value.