RBF-FD analysis of 2D time-domain acoustic wave propagation in heterogeneous media

Radial Basis Function-generated Finite Differences (RBF-FD) is a popular variant of local strong-form meshless methods that do not require a predefined connection between the nodes, making it easier to adapt node-distribution to the problem under consideration. This paper investigates an RBF-FD solution of time-domain acoustic wave propagation in the context of seismic modeling in the Earth’s subsurface. Through a number of numerical tests, ranging from homogeneous to highly-heterogeneous velocity models including non-smooth irregular topography, we demonstrate that the present approach can be further generalized to solve large-scale seismic modeling and full waveform inversion problems in arbitrarily complex models enabling more robust interpretations of geophysical observations

💡 Research Summary

The paper presents a comprehensive study of the Radial Basis Function‑Generated Finite Differences (RBF‑FD) method applied to two‑dimensional time‑domain acoustic wave propagation in heterogeneous Earth models. Unlike traditional grid‑based schemes such as Finite Difference Method (FDM), Finite Element Method (FEM), or Spectral Element Method (SEM), RBF‑FD does not require a predefined connectivity between computational nodes. Instead, it constructs local stencils (support domains) around each node, computes radial basis function (RBF) weights on‑the‑fly, and stores them for repeated use during time stepping.

The authors adopt Gaussian RBFs (Φ(r)=exp(−r²/σ_B²)) with a shape parameter σ_B=70 m and a stencil size of n=7 nearest neighbours. For each node i, the weights w_i are obtained by solving a small linear system that enforces exactness of the differential operator L (e.g., the Laplacian) on the set of RBFs centred at the stencil nodes. Because the Gaussian kernel yields a symmetric positive‑definite matrix when the stencil points are distinct, the system is guaranteed to be nonsingular. Once computed, the weights enable O(n) evaluation of spatial derivatives via a simple dot product, making the method computationally efficient despite the upfront cost of weight calculation.

The acoustic wave equation is formulated under the constant‑density assumption: (1/v_p²)∂²u/∂t² = ∇²u + δ(x−x_s, z−z_s)s(t), where v_p(x,z) denotes the spatially varying P‑wave velocity, u is the pressure field, and s(t) is a Ricker wavelet source. The source’s spatial delta function is regularized by a circular distribution of radius ε=4 m, ensuring adequate representation on any node layout. To suppress artificial reflections from the truncated computational domain, the authors employ Cerjan’s absorbing boundary condition (ABC) with an exponential damping factor G(d)=exp