Logo
ScienceStack
Table of Contents
Fetching ...
Sign In

Numerical approximation of Caputo-type advection-diffusion equations in one and multiple spatial dimensions via shifted Chebyshev polynomials

Francisco de la Hoz, Peru Muniain

TL;DR

This work develops a stable, high-accuracy pseudospectral method for Caputo-type time-fractional advection–diffusion equations in one and multiple spatial dimensions by leveraging shifted Chebyshev polynomials. Central to the approach is the construction of fractional differentiation and integration matrices using arbitrary-precision arithmetic to counter instability at large $N$, which then enables transforming space–time discretizations into Sylvester (and Sylvester tensor) equations that can be solved efficiently. The authors provide complete Matlab implementations and demonstrate spectral-like accuracy for Riemann–Liouville integrals and robust, highly accurate handling of Caputo derivatives even for oscillatory time dependence. The methodology scales to higher dimensions via Sylvester tensor solvers and illustrates practical applicability through 1D and multi-dimensional PDE examples with rigorous error benchmarks. This work thus offers a rigorous, implementable framework for accurate fractional-time discretizations in complex PDE settings with potential impact on simulations requiring nonlocal temporal dynamics.

Abstract

In this paper, using a pseudospectral approach, we develop operational matrices based on the shifted Chebyshev polynomials to approximate numerically Caputo fractional derivatives and Riemann-Liouville fractional integrals. In order to make the generation of these matrices stable, we use variable precision arithmetic. Then, we apply the Caputo differentiation matrices to solve numerically Caputo-type advection-diffusion equations in one and multiple spatial dimensions, which involves transforming the discretization of the concerning equation into a Sylvester (tensor) equation. We provide complete Matlab codes, whose implementation is carefully explained. The numerical experiments involving highly oscillatory functions in time confirm the effectiveness of this approach.

Numerical approximation of Caputo-type advection-diffusion equations in one and multiple spatial dimensions via shifted Chebyshev polynomials

TL;DR

This work develops a stable, high-accuracy pseudospectral method for Caputo-type time-fractional advection–diffusion equations in one and multiple spatial dimensions by leveraging shifted Chebyshev polynomials. Central to the approach is the construction of fractional differentiation and integration matrices using arbitrary-precision arithmetic to counter instability at large , which then enables transforming space–time discretizations into Sylvester (and Sylvester tensor) equations that can be solved efficiently. The authors provide complete Matlab implementations and demonstrate spectral-like accuracy for Riemann–Liouville integrals and robust, highly accurate handling of Caputo derivatives even for oscillatory time dependence. The methodology scales to higher dimensions via Sylvester tensor solvers and illustrates practical applicability through 1D and multi-dimensional PDE examples with rigorous error benchmarks. This work thus offers a rigorous, implementable framework for accurate fractional-time discretizations in complex PDE settings with potential impact on simulations requiring nonlocal temporal dynamics.

Abstract

In this paper, using a pseudospectral approach, we develop operational matrices based on the shifted Chebyshev polynomials to approximate numerically Caputo fractional derivatives and Riemann-Liouville fractional integrals. In order to make the generation of these matrices stable, we use variable precision arithmetic. Then, we apply the Caputo differentiation matrices to solve numerically Caputo-type advection-diffusion equations in one and multiple spatial dimensions, which involves transforming the discretization of the concerning equation into a Sylvester (tensor) equation. We provide complete Matlab codes, whose implementation is carefully explained. The numerical experiments involving highly oscillatory functions in time confirm the effectiveness of this approach.

Paper Structure

This paper contains 17 sections, 80 equations, 3 figures.

Table of Contents

  1. Introduction
  2. A matrix-based pseudospectral method method to approximate numerically $D_t^{\alpha} f(t)$ and $I_t^{\alpha} f(t)$
  3. Obtention of $\{\hat{f}_k\}$
  4. Matrix-based obtention of $\{\hat{f}_k\}$ from $\{f_j\}$
  5. Construction of the fractional differentiation and integration matrices
  6. Implementation in Matlab
  7. Obtention of the coefficients $\{f_k\}$
  8. Computation of $\hat{\mathbf D}_t^\alpha$ and $\hat{\mathbf E}_t^\alpha$
  9. Computation of $\mathbf D_t^\alpha$ and $\mathbf E_t^\alpha$
  10. Matlab codes to generate $\hat{\mathbf D}_t^\alpha$, $\hat{\mathbf E}_t^\alpha$, $\mathbf D_t^\alpha$ and $\mathbf E_t^\alpha$
  11. Estimation of the number of digits
  12. Numerical tests
  13. A full numerical test in Matlab
  14. Numerical approximation of $D_t^{\alpha} e^{imt}$ and $I_t^{\alpha} e^{imt}$ for large natural values of $m$ and different values of $\alpha$ and $N$
  15. Numerical simulation of Caputo-type advection-diffusion equations
  16. ...and 2 more sections

Figures (3)

  • Figure 1: $\|\hat{\mathbf D}_t^\alpha(\mathtt{dig}) - \hat{\mathbf D}_t^\alpha(\mathtt{dig}+1)\|_{\max}$ (red), $\|\mathbf D_t^\alpha(\mathtt{dig}) - \mathbf D_t^\alpha(\mathtt{dig}+1)\|_{\max}$ (blue), $\|\hat{\mathbf E}_t^\alpha(\mathtt{dig}) - \hat{\mathbf E}_t^\alpha(\mathtt{dig}+1)\|_{\max}$ (green) and $\|\mathbf E_t^\alpha(\mathtt{dig}) - \mathbf E_t^\alpha(\mathtt{dig}+1)\|_{\max}$ (cyan), for $\mathtt{dig} = 2, 3, 4, \ldots$. In all cases, we have taken $N = 100$, $\alpha = 0.37$ and $T = 1.2$.
  • Figure 2: Left: minimum values of dig for which both $\|\hat{\mathbf D}_t^\alpha(\mathtt{dig}) - \hat{\mathbf D}_t^\alpha(\mathtt{dig}+1)\|_{\max} = 0$ and $\|\mathbf D_t^\alpha(\mathtt{dig}) - \mathbf D_t^\alpha(\mathtt{dig}+1)\|_{\max} = 0$ (red), and minumum values of dig for which both $\|\hat{\mathbf E}_t^\alpha(\mathtt{dig}) - \hat{\mathbf E}_t^\alpha(\mathtt{dig}+1)\|_{\max} = 0$ and $\|\mathbf E_t^\alpha(\mathtt{dig}) - \mathbf E_t^\alpha(\mathtt{dig}+1)\|_{\max} = 0$ (blue), together with $\max_{1\le i\le N}|c_{iN}|$, for $2\le N \le 500$. Right: elapsed time of the experiments on the left-hand side. In all cases, we have taken $\alpha = 0.37$ and $T = 1.2$.
  • Figure 3: Left: Errors $\operatorname{err}_{\hat{\mathbf D}_t^\alpha}$ (red) and $\operatorname{err}_{\mathbf D_t^\alpha}$ (blue) as defined in \ref{['e:errhatDta']} and \ref{['e:errDta']}, respectively, for $f(t) = e^{i110t}$, i.e., $m = 110$, and $\alpha\{0, 0.005, \ldots, 4\}$. For the sake of comparison, we also plot $\operatorname{err}_{\mathbf D^\alpha}$, as defined in \ref{['e:errhatDa']}, for $\alpha\in\{1, 2, 3, 4\}$ (thick black point). Right: Errors $\operatorname{err}_{\hat{\mathbf D}_t^\alpha}$ (red) and $\operatorname{err}_{\mathbf D_t^\alpha}$ (blue), for $f(t) = e^{i110t}$, and $N\in\{100, 105, \ldots, 1000\}$.