A sparse spectral method for fractional differential equations in one-spatial dimension

TL;DR

The paper introduces a sparse spectral framework for 1D fractional PDEs on $\mathbb{R}$ that combines weighted Chebyshev functions of the second kind with their Hilbert transforms in a sum space. The key idea is that the operator $\mathcal{L}_{\lambda,\mu,\eta}$ decouples across different affine interval transforms, yielding independent sparse solves and enabling parallel computation with overall $\mathcal{O}(n)$ complexity per interval. A dual sum space is established to enable efficient RHS expansions, and a diagonal preconditioner guarantees conditioning independent of discretization size and interval count. Numerical experiments on manufactured, discontinuous, and nonsmooth data, as well as fractional heat and wave problems, demonstrate spectral convergence and strong performance, including substantial speedups when evolving time-dependent problems. The method offers a robust path toward extending sparse spectral techniques to broader fractional Laplacians and higher dimensions.

Abstract

We develop a sparse spectral method for a class of fractional differential equations, posed on $\mathbb{R}$, in one dimension. These equations can include sqrt-Laplacian, Hilbert, derivative and identity terms. The numerical method utilizes a basis consisting of weighted Chebyshev polynomials of the second kind in conjunction with their Hilbert transforms. The former functions are supported on $[-1,1]$ whereas the latter have global support. The global approximation space can contain different affine transformations of the basis, mapping $[-1,1]$ to other intervals. Remarkably, not only are the induced linear systems sparse, but the operator decouples across the different affine transformations. Hence, the solve reduces to solving $K$ independent sparse linear systems of size $\mathcal{O}(n)\times \mathcal{O}(n)$, with $\mathcal{O}(n)$ nonzero entries, where $K$ is the number of different intervals and $n$ is the highest polynomial degree contained in the sum space. This results in an $\mathcal{O}(n)$ complexity solve. Applications to fractional heat and wave equations are considered.

Figures (16)

• Figure 1: The domain and range of the operators $(-\Delta)^{1/2}$ and $\mathcal{I}$ for our approximating spaces.
• Figure 2: Plots of $\tilde{T}_1(x)$ and $W_0(x)$ centred at the three different intervals $[-3,-1]$, $[-1,1]$ and $[1,3]$. We note that $\tilde{T}_1(x)$ has global support and those centred on different intervals always overlap. Whereas, $W_0(x)$ has compact support on the interval it is centred on.
• Figure 3: Spy plots of the matrices in the linear systems after discretizing and decoupling \ref{['eq:ex2:case2']} interval-wise. Here $\lambda = \eta = \mu = 1$ and we use the additional functions $v_0$, $\tilde{u}_{-1}$, $v_{n+1}$, $\tilde{u}_{n+1}$. The matrices are sparse and almost banded.
• Figure 4: Error in the numerical solution $u(x)$ as measured by \ref{['eq:man-error']} with the right-hand sides \ref{['eq:ex2:case1']} and \ref{['eq:ex2:case2']} for increasing truncation degree $n$. The approximation space contains $10n +31$ functions for each value of $n$. The convergence is spectral.
• Figure 5: The $l^\infty$-norm of the coefficient vector of the sum space expansion of the right-hand side for increasing truncation degree $n$. The approximation space contains $10n +31$ functions for each value of $n$.

Theorems & Definitions (59)

• Remark 1.1
• Lemma 2.1: Anti-self adjointness of the Hilbert transform
• Theorem 2.1
• proof
• Theorem 2.2: Existence
• proof
• Remark 3.1
• Proposition 3.1: Recurrence relation
• proof
• Proposition 3.2: Hilbert transform