A Second-order method on graded meshes for fractional Laplacian via Riesz fractional derivative with a singular source term

Abstract: The high-order numerical analysis for fractional Laplacian via the Riesz fractional derivative, under the low regularity solution, has presented significant challenges in the past decades. To fill in this gap, we design a grid mapping function on graded meshes to analyse the local truncation errors, which are far less than second-order convergence at the boundary layer. To restore the second-order global errors, we construct an appropriate right-preconditioner for the resulting matrix algebraic equation. We prove that the proposed scheme achieves second-order convergence on graded meshes even if the source term is singular or hypersingular. Numerical experiments illustrate the theoretical results. The proposed approach is applicable for multidimensional fractional diffusion equations, gradient flows and nonlinear equations.

• The paper presents a novel second-order scheme based on difference-quadrature on graded meshes to approximate the fractional Laplacian with singular sources.
• It combines finite difference and quadrature techniques with tailored mesh grading to resolve boundary-induced low regularity and achieve optimal convergence.
• The use of a right preconditioner ensures strict diagonal dominance, enabling robust error control even for hypersingular source terms.

Second-Order Methods for Fractional Laplacian with Riesz Derivatives and Singular Sources on Graded Meshes

Introduction and Problem Setting

The paper "A Second-order method on graded meshes for fractional Laplacian via Riesz fractional derivative with a singular source term" (2502.11117) addresses the design and analysis of second-order accurate numerical schemes for the fractional Laplacian operator, formulated via the Riesz fractional derivative, particularly in scenarios where the solution exhibits low regularity due to singular or hypersingular source terms. The fractional Laplacian is a central object in modeling anomalous diffusion, appearing in mathematical physics, finance, and image processing. However, its nonlocal and boundary-singular nature presents considerable challenges for the development of high-order discretization methods, especially on local grids.

Existing approaches in the literature—including finite difference, finite element, and spectral methods—typically experience severe order reduction when the solution lacks regularity near the boundary or the source term is singular, with global convergence guarantees deteriorating to $O(h^a)$ or worse for $a \in (1,2)$. The key technical hurdle is controlling the local and global truncation errors induced by the boundary singularity inherent to both the solution and the data.

Methodology: Difference-Quadrature on Graded Meshes

The authors propose a difference-quadrature (DQ) scheme that synergistically combines finite difference approximations for the local (differential) part of the operator with quadrature rules targeting the nonlocal (integral) components. The domain is partitioned using a graded mesh characterized by a grading parameter $r \geq 1$, which produces non-uniform node distributions dense near the boundary. This mesh grading is essential: it enables better resolution of boundary layers where the solution exhibits reduced regularity, thus mitigating error inflation inherent to uniform meshes.

The spatial discretization employs continuous piecewise linear finite elements (hat functions) as a basis. Collocation is performed at the nodal points, resulting in an algebraic system $A U = F$, where $A$ encapsulates both the differential and fractional terms, and $F$ incorporates the (potentially singular) source data.

To address the truncation error at the mesh boundary, which dominates the global error, the analysis leverages a carefully constructed grid mapping function and a "natural-skew" ordering of difference quotients. This framework enables precise partitioning and estimation of local truncation errors. A critical theoretical contribution is the construction of a right preconditioner for the discretized system; this not only ensures strict diagonal dominance of the system matrix but is also crucial for recovering second-order convergence globally, even as the local error at the boundary can be substantially less regular than in the interior.

Theoretical Results: Regularity and Error Estimates

The regularity theory revisits the boundary singularity of the solution to the fractional Laplacian: for $f \in L^\infty(\Omega)$, it is known that solutions have the leading-order behavior $u(x)\sim C \delta(x)^{a/2}$ as $x\rightarrow \partial\Omega$, where $\delta(x)$ is the distance to the boundary. The authors develop and employ weighted Hölder norms adapted to this setting.

Main Theorems Highlight:

• Under $f \in L^\infty \cap C^{B}(\Omega)$, with the mesh grading parameter $r \geq 1$, the scheme achieves a global convergence rate of $O(h^{\min\{r a/2,2\}})$.
• With a special choice of grading parameter (e.g., $r=4/a$ for the regular case, $r=2/o$ for source $f(x)=x^{-o}$), optimal quadratic convergence in the maximum norm is achieved even for source terms that are singular/hypersingular at the boundary.
• The proposed right preconditioner is explicitly shown to guarantee strict diagonal dominance of the resulting system. Without it, a significant reduction in global order is theoretically and numerically observed.

Contrary to much of the prior literature, second-order global accuracy is not lost at the boundary layer, provided the mesh grading and preconditioner are carefully chosen—even when the source term itself is singular.

Numerical Experiments

The effectiveness of the proposed method is robustly demonstrated through numerical experiments on prototypical test cases, including both regular and singular source terms. Errors are systematically measured in the maximum norm, with convergence rates given by $\log_2 (E_{N/2}/E_N)$ for mesh refinements.

For a regular source ($f=1$), the scheme delivers rates matching $O(h^2)$ for both uniform ($r=1$) and optimal graded meshes ($r=4/a$). For a singular source $f(x)=x^{-o}$, the method still achieves $O(h^2)$ on appropriately graded meshes ($r=2/o$), while in the uniform mesh case, the convergence stagnates at $O(h^o)$, reflecting the theoretical analysis.

Implications and Further Directions

The proposed DQ method on graded meshes, equipped with theoretically sound preconditioning, resolves a long-standing issue in the robust high-order numerical solution of nonlocal fractional elliptic problems with non-smooth data. This result carries several important practical and theoretical implications:

• General Applicability: The analysis and algorithm are extensible to multi-dimensional domains, time-dependent fractional diffusion equations, and gradient flows governed by similar nonlocal operators.
• Handling Hypersingular Terms: By rigorously quantifying errors under minimal regularity, the scheme is directly applicable to models with measure-valued or even hypersingular sources, closing a significant gap left by previous research.
• Preconditioning Strategy: The explicit preconditioning approach for the stiffness matrix in the presence of strong heterogeneity or ill-conditioning is relevant beyond this specific context, potentially benefiting a broad range of nonlocal PDE discretizations.

Future work may focus on extending the analysis to nonlinear equations, adaptivity with dynamic mesh refinement, and coupling to time-fractional or variable-coefficient models prevalent in realistic anomalous diffusion scenarios.

Conclusion

This work establishes a rigorous, mesh-sensitive numerical analysis for the fractional Laplacian based on the Riesz derivative under low-regularity conditions, overcoming the boundary-induced order reduction that plagues standard schemes. The theoretical framework, blending graded meshes, difference-quadrature discretization, and matrix preconditioning, is substantiated by sharp convergence proofs and confirmed by numerical experiments. The results decisively clarify the conditions under which second-order accuracy is attainable for fractional Laplacian problems with singular sources, setting a new standard for numerical analysis in this domain.