gRBF-FD: Mesh-Free Finite Difference Method

• gRBF-FD is a mesh-free framework that uses RBF interpolation combined with polynomial augmentation to discretize differential operators on arbitrary geometries.
• It integrates overlap parameters, least-squares stabilization, and hybrid finite differences to improve accuracy, stability, and computational efficiency.
• Adaptive node layouts and manifold-based extensions enable this method to effectively address PDEs and quadrature problems in complex, scattered data settings.

The generalized radial basis function-generated finite difference (gRBF-FD) method encompasses a family of mesh-free numerical algorithms for discretizing and solving differential equations, quadrature problems, and operator equations on arbitrary geometries and node sets. It unifies and extends the classical RBF-FD framework by incorporating overlap parameters, polynomial augmentation, least-squares stabilization, hybridization with classical finite differences, and adaptive or manifold-based stenciling, often yielding substantial improvements in stability, accuracy, and computational efficiency across diverse applications in scientific computing.

1. Conceptual Foundations and General Methodology

The gRBF-FD approach constructs local finite difference stencils using radial basis function interpolation, typically augmented with polynomials, on potentially scattered nodes $\{x_j\}_{j=1}^N$ in a domain $\Omega \subset \mathbb{R}^d$. For each center node, a local stencil $S(x)$ of $n$ neighbors is selected, and a local interpolant is defined: $$I_\phi f(x') = \sum_{j} \alpha_j \phi(\|x' - x_j\|) + \sum_{\ell} \beta_\ell p_\ell(x'),$$ where $\phi$ is an RBF (e.g., polyharmonic spline $r^q$, Gaussian, inverse multiquadric), and $p_\ell$ spans the space of polynomials up to total degree $m$. Imposing interpolation conditions and polynomial side-constraints yields a local block system whose solution provides local weights for discretizing differential operators or quadrature integrands. This methodology is consistent across standard RBF-FD (Shankar, 2016), surface PDE solvers (Shankar et al., 2014), adaptive strategies (Oanh et al., 2016), least-squares stabilization (Tominec et al., 2020), hybrid FD formulations (Rogan et al., 20 May 2025), high-dimensional manifolds (Li et al., 22 Nov 2025), and meshless quadrature (Reeger, 2019).

2. Overlapped RBF-FD: The $\delta$ Parameter and Local Lebesgue Stabilization

The overlapped RBF-FD method (Shankar, 2016) parameterizes stencil overlap via $\delta \in [0,1]$: $\delta=1$ recovers classic per-center stencils; $\delta=0$ yields full decoupling, i.e., weight computation at all stencil nodes. Given each center $x_k$ and stencil $P_k$, the retention ball $\mathbb{B}_k$ is defined with radius $r_k = (1-\delta)\rho_k$, where $\rho_k$ is the stencil radius. The ORBFD algorithm computes and retains weights at all nodes in $\mathbb{B}_k$ from a shared local system.

Automatic stabilization is achieved by evaluating $\mathcal{L}$-Lebesgue functions $\Lambda^k_{\mathcal{L}}(x)$, accepting computed weights only for nodes where $\Lambda^k_{\mathcal{L}}(x)$ does not exceed that at the center. This criterion, coupled with Gershgorin disk theory, promotes stability (eigenvalues in disks with non-positive real parts) and controls local error, as formalized in error bounds involving the Lebesgue constant. Quantitative speed-up estimates $\eta$ show $\sim16\times$ acceleration in 2D and $\sim60\times$ in 3D for differentiation matrix assembly with negligible loss in accuracy for moderate $\delta$.

3. Hybridization and Least Squares Generalizations

Hybrid gRBF-FD formulations (Rogan et al., 20 May 2025) combine RBF interpolation with classical finite difference stencils, e.g., 5-point or 9-point Laplacian formulas, interpolating the virtual FD stencil points via local RBF supports. This enables direct reproduction of polynomial degrees matched to the chosen FD stencil and exploits RBF smoothness on scattered nodes. Tuning the virtual stencil spacing $\delta \approx h$ (where $h$ is the fill distance) is crucial for optimal accuracy.

Least squares gRBF-FD (RBF-FD-LS) (Tominec et al., 2020) replaces collocation-based enforcement of PDE and boundary conditions with a rectangular, oversampled system. Stencil weights are reused across evaluation points, and the global PDE operator is solved in the discrete $\ell_2$ least squares sense. This yields enhanced stability (robust under Neumann BCs), maintains high-order convergence ($O(h^{p-1})$ for polynomial degree $p$), and is computationally competitive with collocation RBF-FD.

4. Manifold and Surface Extensions

gRBF-FD has been extended to embedded surfaces and manifolds (Shankar et al., 2014, Li et al., 22 Nov 2025), including closed surfaces in $\mathbb{R}^d$ and randomly sampled point clouds. Differential operators (e.g., Laplace–Beltrami) are approximated using extrinsic coordinates via projection matrices or via local Monge charts with tangent space coordinates. For manifold problems, coefficients are computed in two stages: first, local generalized moving least squares (GMLS) regression with polynomial weights, then correction by PHS interpolation of the residual (Li et al., 22 Nov 2025). Weight structures (diagonal spike criteria) and automatic tuning of stencil size $K$ further promote stability. Rigorous error estimates show consistency $O(D_\mathrm{max}^{l-1})$ as a function of local stencil diameter $D_\mathrm{max}$ and point count $N$.

5. Adaptive Strategies and Node Layouts

Adaptive gRBF-FD (Oanh et al., 2016) targets singularities and rapid variations by iterative refinement of node clouds and stencil supports, using Zienkiewicz–Zhu error indicators constructed from local least squares recoveries. The refinement algorithm adds nodes by marking edges with high discrepancy in meshless directional derivatives, updating the cloud until error equilibration is achieved. In all tested scenarios, the adaptive gRBF-FD matches or slightly outperforms adaptive finite element methods for elliptic PDEs with singularities.

Smoothly varying-density node layouts (Milovanović, 2018) are generated using a radius function $R(x)$ and advancing-front point placement, augmented by boundary repulsion iterations. This approach minimizes abrupt local changes in mesh spacing, keeping local interpolation matrix condition numbers moderate across large $N$ and ensuring stable stencil-weight computation. Polynomial augmentation ($p\geq 4$), odd-degree PHS kernels ($q=5,~7$), and large stencils ($n \approx 5m$) combine for high-order accuracy and stability—empirically, second-order global convergence is consistently observed in non-smooth financial PDE test cases.

6. Meshless Quadrature via gRBF-FD Principles

The gRBF-FD framework has also been applied to high-order quadrature on arbitrary scattered nodes in $\mathbb{R}^3$ (Reeger, 2019). Instead of operator weights, local RBF+polynomial interpolants are built on small stencils associated with tetrahedral integration cells (from a domain Delaunay tessellation), and basis-function integrals (moments) are enforced exactly to define quadrature weights. This yields a global quadrature rule with convergence rates $O(N^{-m/3})$ for polynomials up to degree $m$, drastically improving flexibility and accuracy compared to structured quadrature schemes.

7. Computational Complexity, Parameter Selection, and Practical Guidelines

The gRBF-FD methods exhibit computational complexity scaling as $O(N)$ or $O(N\log N)$, depending on nearest-neighbor searches, tessellation, and local system assembly. Augmented block matrices (RBF+poly) are small for moderate $m$; LU factorization and reuse are exploited for efficiency, particularly in hybrid and least-squares variants.

Key parameter guidelines from the literature:

Parameter	Recommendation	Source
Augmentation degree	$p\geq$ FD reproduction deg; $p=4$ typical	(Rogan et al., 20 May 2025, Milovanović, 2018)
PHS degree	Odd $q=5$ or $7$ for stability	(Milovanović, 2018)
Stencil size	$n \approx 5m$	(Milovanović, 2018)
Overlap $\delta$	$0.2 \leq \delta \leq 1$ for good accuracy/speed-up	(Shankar, 2016)
LS oversampling	$q=M/N \in [1,3]$	(Tominec et al., 2020)
Virtual stencil $\sigma$	$\sigma \approx 1$	(Rogan et al., 20 May 2025)

Implementation strategies such as adaptively tuning overlap or stencil size, enforcing Lebesgue constant control, and weighting schemes for stability are prevalent across both finite difference and quadrature contexts.

---

gRBF-FD methods represent a unifying, mesh-free framework for high-order approximation of differential operators and quadrature on scattered data and manifolds. Their generalization mechanisms—overlap parameters, polynomial augmentation, stabilization via Lebesgue functions and least-squares enforcement, hybridization with classical FD, adaptive node selection, and automatic parameter tuning—enable robust, scalable, and accurate discretizations for a wide spectrum of scientific and engineering problems (Shankar, 2016, Shankar et al., 2014, Li et al., 22 Nov 2025, Milovanović, 2018, Rogan et al., 20 May 2025, Tominec et al., 2020, Oanh et al., 2016, Reeger, 2019).