A novel and simple spectral method for nonlocal PDEs with the fractional Laplacian

TL;DR

The paper develops a novel spectral method for the nonlocal fractional Laplacian $(-\Delta)^{\alpha/2}$ based on a semi-discrete Fourier transform, yielding an exact discrete symbol $|\xi|^\alpha$ and producing a multilevel Toeplitz stiffness matrix suitable for fast FFT-based solutions. It establishes rigorous error estimates showing spectral accuracy for smooth data, along with stability and convergence results for fractional Poisson problems. The approach extends to high dimensions via two generalization pathways, one leading to analytically expressed $_1F_2$-based weights and the other to a radial-Hankel formulation that preserves dimension-independent coefficient evaluation. Numerically, the method demonstrates superior accuracy with far fewer grid points than finite-difference schemes, and scales efficiently in higher dimensions, offering a practical tool for simulating fractional PDEs and anomalous diffusion phenomena.

Abstract

We propose a novel and simple spectral method based on the semi-discrete Fourier transforms to discretize the fractional Laplacian $(-Δ)^\fracα{2}$. Numerical analysis and experiments are provided to study its performance. Our method has the same symbol $|\boldsymbolξ|^α$ as the fractional Laplacian $(-Δ)^\fracα{2}$ at the discrete level, and thus it can be viewed as the exact discrete analogue of the fractional Laplacian. This {\it unique feature} distinguishes our method from other existing methods for the fractional Laplacian. Note that our method is different from the Fourier pseudospectral methods in the literature which are usually limited to periodic boundary conditions (see Remark \ref{remark0}). Numerical analysis shows that our method can achieve a spectral accuracy. The stability and convergence of our method in solving the fractional Poisson equations were analyzed. Our scheme yields a multilevel Toeplitz stiffness matrix, and thus fast algorithms can be developed for efficient matrix-vector multiplications. The computational complexity is ${\mathcal O}(2N\log(2N))$, and the memory storage is ${\mathcal O}(N)$ with $N$ the total number of points. Extensive numerical experiments verify our analytical results and demonstrate the effectiveness of our method in solving various problems.

Figures (8)

• Figure 1: Comparison of our method and the existing Fourier pseudo-spectral methods in Kirkpatrick2016Duo2016 for approximating $(-\Delta)^\frac{\alpha}{2}u(x)$ on bounded domain.
• Figure 2: Numerical errors of our method in solving the benchmark Poisson problem (\ref{['elliptic-ex1']}) with $s = \alpha/2$, where the exact solution $u\in C^{0,\frac{\alpha}{2}}(\bar{\Omega})$.
• Figure 3: Comparison of our method and others (including Huang in Huang2014, Duo in Duo2018, Minden in Minden2020, and Hao in Hao2021) in approximating $(-\Delta)^{\frac{\alpha}{2}}u$ on $\Omega = (-1, 1)$, where $u$ is defined in (\ref{['ex2']}) with $a = 1$.
• Figure 4: Numerical results of $(-\Delta)^\frac{1}{2} u(x,y)$, where $u$ is from \ref{['ex3']} with $s = 2$ (a) or $s = 4$ (b).
• Figure 5: Numerical errors of our method in solving the Poisson problem (\ref{['elliptic-ex1']}).

Theorems & Definitions (29)

• Remark 1.1
• Definition 2.1: Semi-discrete Fourier transform
• Definition 2.2: Generalized hypergeometric function
• Lemma 2.1
• proof
• Remark 2.1
• Lemma 3.1
• Lemma 3.2
• Lemma 3.3
• proof