[Paper Review] An exponentially convergent discretization for space-time fractional parabolic equations using $hp$-FEM

Jens Markus Melenk, Alexander Rieder|arXiv (Cornell University)|Feb 4, 2022

Fractional Differential Equations Solutions31 references4 citations

TL;DR

This paper presents an $hp$-finite element method (hp-FEM) combined with sinc-quadrature for solving space-time fractional parabolic equations with Caputo time derivatives and spectral fractional Laplacian operators. The method achieves exponential convergence in both space and time, even in the presence of startup singularities due to data incompatibilities, by leveraging Riesz-Dunford functional calculus and geometric mesh refinement near corners.

ABSTRACT

We consider a space-time fractional parabolic problem. Combining a sinc-quadrature based method for discretizing the Riesz-Dunford integral with $hp$-FEM in space yields an exponentially convergent scheme for the initial boundary value problem with homogeneous right-hand side. For the inhomogeneous problem, an $hp$-quadrature scheme is implemented. We rigorously prove exponential convergence with focus on small times $t$, proving robustness with respect to startup singularities due to data incompatibilities.

Motivation & Objective

• To develop a robust, high-order numerical scheme for space-time fractional parabolic problems with Caputo time derivatives and fractional Laplacian operators.
• To address the challenge of startup singularities arising from incompatible initial data and boundary conditions in fractional diffusion problems.
• To achieve exponential convergence rates in both space and time using $hp$-FEM and sinc-quadrature, even for non-smooth initial data.
• To extend the applicability of $hp$-FEM to time-dependent fractional problems beyond the classical parabolic case, using functional calculus and quadrature techniques.
• To prove robust convergence in a weaker space-time energy norm under an abstract assumption on the initial condition, overcoming time-norm degeneracy at small times.

Proposed method

• Uses the Riesz-Dunford functional calculus to represent the fractional power operator $\mathcal{L}^\beta$ via a contour integral.
• Applies sinc-quadrature to discretize the contour integral, with parameters optimized for exponential convergence.
• Employs $hp$-finite elements in space with geometric meshes graded toward domain corners to resolve singularities.
• Implements an $hp$-quadrature scheme for the inhomogeneous problem, ensuring exponential convergence in time.
• Reorders computations to minimize the number of linear solves by reusing solutions of the shifted problems $(z_j - \mathcal{L})^{-1}$.
• Uses the NGSolve software package for finite element assembly and geometric mesh generation with $hp$-refinement.

Experimental results

Research questions

• RQ1 Can an $hp$-FEM approach combined with sinc-quadrature achieve exponential convergence for time-dependent space-time fractional parabolic problems with Caputo derivatives?
• RQ2 How can the method handle startup singularities caused by data incompatibilities in the initial condition and boundary data?
• RQ3 What is the impact of the time-norm degeneracy at small times $t>0$ on convergence estimates, and how can it be mitigated?
• RQ4 Can the method maintain exponential convergence in a weaker space-time energy norm under minimal assumptions on the initial data?
• RQ5 How does the choice of quadrature parameters $k$ and $\mathcal{N}_{\textrm{q}}$ affect the convergence rate and robustness?

Key findings

• The method achieves exponential convergence in the $L^2(\Omega)$-norm with respect to the number of degrees of freedom $\mathcal{N}_\Omega$, as confirmed by numerical experiments with smooth solutions.
• For incompatible data, exponential convergence is still observed with respect to $\mathcal{N}_\Omega^{1/4}$, demonstrating robustness to singularities.
• The sinc-quadrature error estimate is improved from $t^{-\gamma}$ to $t^{-\gamma/2}$, enhancing stability at small times.
• The method requires only $\mathcal{N}_{\textrm{q}}$ linear system solves due to efficient reordering, making it computationally efficient.
• Numerical results confirm exponential convergence for both compatible and incompatible data, even when the exact solution is unknown.
• The convergence rate is robust across different end times $t=0.001$, $t=0.1$, and $t=1$, with error decaying as $\mathcal{O}(e^{-0.7\mathcal{N}_\Omega^{1/4}})$ in the incompatible case.