Numerical Treatment of the Time Fractional Diffusion Wave Problem Using Chebyshev Polynomials

Abstract

This paper introduces an efficient numerical method based on applying the typical Petrov–Galerkin approach ($\mathrm{PGA}$) to solve the time fractional diffusion wave equation ($\mathcal{TFDWE}$). The method utilises asymmetric polynomials, namely, shifted second-kind Chebyshev polynomials ($\mathrm{SSKCPs}$). New derivative formulas are derived and used for these polynomials to establish the operational matrices of their derivatives. The paper presents rigorous error bounds for the proposed method in Chebyshev-weighted Sobolev space and demonstrates its accuracy and efficiency through several illustrative numerical examples. The results reveal that the method achieves high accuracy with relatively low polynomial degrees.

1. Introduction

Fractional differential equations (FDEs) are crucial in numerous areas of the practical sciences. They provide an excellent explanation of why traditional differential equations (DEs) fail to capture some phenomena. This is due to their exceptional ability to model memory and hereditary properties. Numerical analysis is usually employed to solve these equations due to their analytical impossibility. Various FDEs have been solved numerically using several methods. Among these methods are the Chelyshkov wavelet scheme [1], the Ritz-Piecewise Gegenbauer approach [2], the homotopy analysis transform method [3], the Adomian spectral method [4], the Adomian decomposition method [5], the finite difference method [6], the finite element method [7], and Laplace optimized decomposition [8].

Spectral methods are powerful techniques for solving differential equations (DEs). FDEs and high-order ordinary DEs can be solved effectively using these methods. For smooth problems, their ability to achieve exponential or high-order convergence makes them incredibly precise, which constitutes their main advantage over traditional numerical methods. These methods extend the solution based on global basis functions, usually special functions or special polynomials, providing a precise approximation with a minimal number of degrees of freedom. Many studies have focused on using these methods to solve various types of DEs. For instance, the authors in [9] used $\mathrm{PGA}$ to approximate the Kudryashov–Sinelshchikov equation. In [10], the authors used the collocation method to approximate the nonlinear ordinary and fractional Newell–Whitehead–Segel equation. Also, in [11], the author used $\mathrm{PGA}$ to solve the $\mathcal{TFDWE}$. For more studies, see [12,13,14,15,16,17,18,19,20].

Orthogonal polynomials are widely used in science and engineering due to their excellent computational and approximation capabilities. They can be classified as either symmetric or asymmetric. The Chebyshev polynomials (CPs) are highly regarded for their significance and are among the essential orthogonal polynomials. The four well-known CPs are all special cases of Jacobi polynomials. The first and second types of polynomials are symmetric, whereas the third and fourth types are asymmetric. In numerical analysis, they are very important for spectral techniques for solving DEs because they give very accurate results with few retained modes. Some applications of orthogonal polynomials can be found in [21,22,23,24]. Of paramount relevance in theory and practice are the CPs. The CPs and their shifted forms have been used to solve various types of DEs. The authors in [25] employed CPs to treat high-order DEs. The authors in [26] used the fifth-kind CPs spectral method to solve the convection-diffusion equation. The authors in [27] presented a numerical technique based on CPs. For more studies, see [28,29,30,31].

Our main contributions, including the novelty of our work, are listed in the following items:

Proposing new basis functions in terms of $\mathrm{SSKCPs}$ to treat this type of FDE.

Developing some new theoretical results of $\mathrm{SSKCPs}$, such as their definite integral formulas.

Designing a new PGA for treating the $\mathcal{TFDWE}$ based on the theoretical background of these polynomials.

The paper is structured as follows: The following section describes some of the characteristics of fractional calculus. Moreover, some features of $\mathrm{SSKCPs}$ are given in this section. Section 3 proposes a numerical approach to solving $\mathcal{TFDWE}$ with homogeneous conditions using $\mathrm{PGA}$. Section 4 investigates the error bound of the developed double expansion. Some numerical experiments are described in Section 5 to illustrate the efficiency and precision of our numerical scheme. Finally, some conclusions are reported in Section 6.

2. Some Fundamentals

2.1. Caputo Fractional Derivative

Definition 1

([32]). The Caputo fractional derivative of order ν is defined as:

${\displaystyle D_{\rho }^{\nu }\mathrm{Z}\left(\rho \right)=\frac{1}{\Gamma (s-\nu )}{\int }_0^{\rho }{(\rho -y)}^{s-\nu -1}{\mathrm{Z}}^{\left(s\right)}\left(y\right)dy,\nu >0,\rho >0,}$

where $s-1⩽\nu <s,s\in \mathbb{N}$.

The following features are classified by the operator $D_{\rho }^{\nu }$ for $s-1⩽\nu <s,s\in \mathbb{N}$,

$D_{\rho }^{\nu }c=0,\left(cisaconstant\right)$

$D_{\rho }^{\nu }{\rho }^s=\left\{\begin{array}{cc}0,\hfill & if s\in {\mathbb{N}}_{0}ands<\lceil \nu \rceil ,\hfill \\ {\displaystyle \frac{s!}{\Gamma (s-\nu +1)}}{\rho }^{s-\nu },\hfill & if s\in {\mathbb{N}}_{0}andsp\ge \lceil \nu \rceil ,\hfill \end{array}\right.$

where ${\mathbb{N}}_0=\mathbb{N}\cup \left\{0\right\}$, $\mathbb{N}=\{1,2,3,\dots \}$ and the notation $\lceil \nu \rceil$ denotes the ceiling function.

2.2. An Overview of $\mathrm{SSKCPs}$

The $\mathrm{SSKCPs}$ defined in the interval $[0,1]$ by ${\mathbf{U}}_k^s\left(\rho \right)=U_k(2\rho -1)$ are denoted by ${\mathbf{U}}_k^s\left(\rho \right)$. The definition of these polynomials is [33]

(1)${\mathbf{U}}_k^s\left(\rho \right)=\sum _{r=0}^k{\chi }_{r,k}{\rho }^r,k\ge 0,$

where

${\chi }_{r,k}={\displaystyle \frac{2^{2r}{(-1)}^{k+r}(k+r+1)!}{(2r+1)!(k-r)!}},$

and fulfills the orthogonality relation [33]:

(2)${\int }_0^1\widehat{w}\left(\rho \right){\mathbf{U}}_m^s\left(\rho \right){\mathbf{U}}_n^s\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}{\delta }_{m,n},$

where $\widehat{w}\left(\rho \right)=\sqrt{\rho -{\rho }^2}$ and ${\delta }_{m,n}$ is the well-known Kronecker delta.

The recurrence relation of ${\mathbf{U}}_m^s\left(\rho \right)$ is

${\mathbf{U}}_m^s\left(\rho \right)=2\left(2\rho -1\right){\mathbf{U}}_{m-1}^s\left(\rho \right)-U_{m-2}^s\left(\rho \right),$

where ${\mathbf{U}}_0^s\left(\rho \right)=1,{\mathbf{U}}_1^s\left(\rho \right)=2\rho -1.$

Moreover, the inversion formula is [33]

(3)${\rho }^k=\sum _{p=0}^k{B}_{p,k}{\mathbf{U}}_p^s\left(\rho \right),j\ge 0,$

where

$B_{p,k}={\displaystyle \frac{4\Gamma \left(k+\frac{3}{2}\right)(p+1)k!}{\sqrt{\pi }(k-p)!(p+k+2)!}}.$

Lemma 1

([33]). Let i and m be any two integers that are not negative. The moment formula for the ${\mathbf{U}}_m^s\left(\rho \right)$ is provided by

${\rho }^i{\mathbf{U}}_m^s\left(\rho \right)=\sum _{k=-i-m}^{i+m}{F}_{k,i,m}{\mathbf{U}}_k^s\left(\rho \right),$

where

$F_{k,i,m}={\displaystyle \frac{1}{2^{2i}}}\left({\displaystyle \genfrac{}{}{0pt}{}{2i}{i-k+m}}\right).$

Corollary 1

([34]). The first derivative of ${\mathbf{U}}_j^s\left(t\right)$ is

${\displaystyle \frac{d{\mathbf{U}}_j^s\left(t\right)}{dt}}=4\sum _{\begin{array}{c}p=0\\ (p+j)\mathrm{odd}\end{array}}^{j-1}(p+1){\mathbf{U}}_p^s\left(t\right),j\ge 1.$

Remark 1.

For $m\in \mathbb{Z}$ and $n\in {\mathbb{N}}_0$, the following relation is satisfied

(4) ${\int }_0^1\widehat{w}\left(\rho \right){\mathbf{U}}_m^s\left(\rho \right){\mathbf{U}}_n^s\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}{\lambda }_{m,n},$

where

${\lambda }_{m,n}=\left\{\begin{array}{cc}1,\hfill & ifm\ge 0,n=m,\hfill \\ -1,\hfill & ifm<-1,n+m=-2,\hfill \\ 0,\hfill & otherwise.\hfill \end{array}\right.$

3. Treatment for the $\mathcal{TFDWE}$ with Homogeneous Conditions

This section focuses on examining $\mathrm{PGA}$ to address the following $\mathcal{TFDWE}$ [35]:

(5)${\displaystyle \frac{{\partial }^{\nu }\mathrm{Z}(\rho ,t)}{\partial t^{\nu }}}+a{\displaystyle \frac{\partial \mathrm{Z}(\rho ,t)}{\partial t}}+b\mathrm{Z}(\rho ,t)-{\displaystyle \frac{{\partial }^2\mathrm{Z}(\rho ,t)}{\partial {\rho }^2}}=g(\rho ,t),1<\nu <2,$

subject to homogeneous initial and boundary conditions ($\mathcal{HIBC}s$)

(6)$\begin{array}{c}\mathrm{Z}(\rho ,0)={\mathrm{Z}}_t(\rho ,0)=0,0\le \rho \le 1,\hfill \\ \mathrm{Z}(0,t)=\mathrm{Z}(1,t)=0,0\le t\le 1,\hfill \end{array}$

where $g(\rho ,t)$ is the source term and a and b are the coefficients of the damping and reaction terms, respectively.

3.1. Basis Functions

Consider the following basis functions

(7)$\begin{array}{c}{\zeta }_i\left(\rho \right)=\rho (1-\rho ){\mathbf{U}}_i^s\left(\rho \right),\hfill \\ {\mathcal{P}}_j\left(t\right)=t^2{\mathbf{U}}_j^s\left(t\right).\hfill \end{array}$

Remark 2.

The basis functions specified in (7) meet the conditions listed below

${\mathcal{P}}_k\left(0\right)={\displaystyle \frac{d{\mathcal{P}}_k\left(0\right)}{d\rho }}={\zeta }_k\left(0\right)={\zeta }_k\left(1\right)=0.$

Also, by virtue of the orthogonality relation in (2), we obtain the following relations

${\int }_0^1{\zeta }_i\left(\rho \right){\zeta }_j\left(\rho \right){\omega }_4\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}{\delta }_{i,j},$

and

${\int }_0^1{\mathcal{P}}_i\left(t\right){\mathcal{P}}_j\left(t\right){\omega }_3\left(t\right)dt={\displaystyle \frac{\pi }{8}}{\delta }_{i,j},$

where ${\omega }_4\left(\rho \right)={\displaystyle \frac{\widehat{w}\left(\rho \right)}{{\rho }^2{(1-\rho )}^2}}$ and ${\omega }_3\left(t\right)={\displaystyle \frac{\widehat{w}\left(t\right)}{t^4}}.$

Theorem 1

([34]). It is possible to clearly express the second-derivative of ${\zeta }_i\left(\rho \right)$ in terms of ${\mathbf{U}}_j^s\left(\rho \right)$ as

${\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}=\sum _{j=0}^i{\mu }_{j,i}{\mathbf{U}}_j^s\left(\rho \right),$

where

${\mu }_{j,i}=\left\{\begin{array}{cc}-2(j+1),\hfill & if(i+j)even\mathrm{and}i>j,\hfill \\ (i+1)(i+2),\hfill & if i=j,\hfill \\ 0,\hfill & otherwise.\hfill \end{array}\right.$

3.2. $\mathrm{PGA}$ for the $\mathcal{TFDWE}$ with Homogeneous Conditions

Assume that the $\mathcal{TFDWE}$ (5), governed by the $\mathcal{HIBC}s$ (6).

Now, consider

$\begin{array}{c}{\mathbb{T}}_{\mathcal{L}}(\Omega )=\mathrm{span}\{{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right):i,j=0,1,\dots ,\mathcal{L}\},\hfill \\ {\mathbb{L}}_{\mathcal{L}}(\Omega )=\{\mathrm{Z}(\rho ,t)\in {\mathbb{T}}_{\mathcal{L}}(\Omega ):\mathrm{Z}(\rho ,0)={\mathrm{Z}}_t(\rho ,0)=\mathrm{Z}(0,t)=\mathrm{Z}(1,t)=0\},\hfill \end{array}$

where $\Omega ={[0,1]}^2.$ Then, any ${\mathrm{Z}}^{\mathcal{L}}(\rho ,t)\in {\mathbb{L}}_{\mathcal{L}}(\Omega )$ can be expressed as

${\mathrm{Z}}^{\mathcal{L}}(\rho ,t)=\sum _{i=0}^{\mathcal{L}}\sum _{j=0}^{\mathcal{L}}{c}_{ij}{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right)=\mathbf{\zeta }\mathit{C}{\mathsf{P}}^T,$

where

$\mathbf{\zeta }=[{\zeta }_0\left(\rho \right),{\zeta }_1\left(\rho \right),\dots ,{\zeta }_{\mathcal{L}}\left(\rho \right)],{\mathsf{P}}^T={[{\mathcal{P}}_0\left(t\right),{\mathcal{P}}_1\left(t\right),\dots ,{\mathcal{P}}_{\mathcal{L}}\left(t\right)]}^T,$

and the matrix of unknowns with order ${(\mathcal{L}+1)}^2$ is $\mathit{C}={\left(c_{ij}\right)}_{0\le i,j\le \mathcal{L}}$.

Now, the residual $\mathrm{R}(\rho ,t)$ of Equation (5) can be calculated

(8)$\begin{array}{cc}\hfill \mathrm{R}(\rho ,t)& ={\displaystyle \frac{{\partial }^{\nu }{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial t^{\nu }}}+a{\displaystyle \frac{\partial {\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial t}}+b{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-{\displaystyle \frac{{\partial }^2{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial {\rho }^2}}-g(\rho ,t).\hfill \end{array}$

The application of $\mathrm{PGA}$ leads to

(9)${\int }_0^1{\int }_0^1\mathrm{R}(\rho ,t){\mathbf{U}}_r^s\left(\rho \right){\mathbf{U}}_s^s\left(t\right)w(\rho ,t)d\rho dt=0,0\le r,s\le \mathcal{L},$

where $w(\rho ,t)=\widehat{w}\left(\rho \right)\widehat{w}\left(t\right)$.

Assume that

$\begin{array}{cc}\hfill \mathit{G}& ={\left({\mathrm{g}}_{r,s}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{g}}_{r,s}={\int }_0^1{\int }_0^{1}g(\rho ,t){\mathbf{U}}_r^s\left(\rho \right){\mathbf{U}}_s^s\left(t\right)w(\rho ,t)d\rho dt,\hfill \\ \hfill \mathsf{M}& ={\left({\mathrm{m}}_{i,r}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{m}}_{i,r}={\int }_0^1{\zeta }_i\left(\rho \right){\mathbf{U}}_r^s\left(\rho \right){\omega }_1\left(\rho \right)d\rho ,\hfill \\ \hfill \mathsf{P}& ={\left({\mathrm{p}}_{i,r}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{p}}_{i,r}={\int }_0^1{\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}{\mathbf{U}}_r^s\left(\rho \right){\omega }_1\left(\rho \right)d\rho ,\hfill \\ \hfill \mathsf{F}& ={\left({\mathrm{f}}_{j,s}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{f}}_{j,s}={\int }_0^1{\mathsf{P}}_j\left(t\right){\mathbf{U}}_s^s\left(t\right){\omega }_1\left(t\right)dt,\hfill \\ \hfill \mathsf{K}& ={\left({\mathrm{k}}_{j,s}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{k}}_{j,s}={\int }_0^1{\displaystyle \frac{d{\mathsf{P}}_j\left(t\right)}{dt}}{\mathbf{U}}_s^s\left(t\right){\omega }_1\left(t\right)dt,\hfill \\ \hfill \mathsf{Q}& ={\left({\mathrm{q}}_{j,s}\right)}_{(\mathcal{L}+1)\times (\mathcal{L}+1)},{\mathrm{q}}_{j,s}={\int }_0^1\left[D_t^{\nu }{\mathsf{P}}_j\left(t\right)\right]{\mathbf{U}}_s^s\left(t\right){\omega }_1\left(t\right)dt.\hfill \end{array}$

Therefore, Equation (9) can be rewritten in matrix form as

(10)${\mathsf{M}}^T\mathit{C}\mathsf{Q}+a{\mathsf{M}}^T\mathit{C}\mathsf{K}+b{\mathsf{M}}^T\mathit{C}\mathsf{F}-{\mathsf{P}}^T\mathit{C}\mathsf{F}=\mathit{G}.$

Finally, the system of Equations (10) of order ${(\mathcal{L}+1)}^2$ can be solved utilizing the Gauss elimination method.

Theorem 2.

The elements of matrices $\mathsf{M},\mathsf{P},\mathsf{F},\mathsf{K}$, and $\mathsf{Q}$ are given as

$\begin{array}{c}\left(1\right){\mathrm{m}}_{i,r}={\int }_0^1{\zeta }_i\left(\rho \right){\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}\left(\sum _{k=-1-i}^{1+i}{F}_{k,1,i}{\lambda }_{k,r}-\sum _{k=-2-i}^{2+i}{F}_{k,2,i}{\lambda }_{k,r}\right),\hfill \\ \left(2\right){\mathrm{p}}_{i,r}={\int }_0^1{\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}{\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}\sum _{j=0}^i{\mu }_{j,i}{\delta }_{j,r},\hfill \\ \left(3\right){\mathrm{f}}_{j,s}={\int }_0^1{\mathcal{P}}_j\left(t\right){\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt={\displaystyle \frac{\pi }{8}}\sum _{k=-2-j}^{2+j}{F}_{k,2,j}{\lambda }_{k,s},\hfill \\ \left(4\right){\mathrm{k}}_{j,s}={\int }_0^1{\displaystyle \frac{d{\mathcal{P}}_j\left(t\right)}{dt}}{\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt={\displaystyle \frac{\pi }{2}}\sum _{k=-2-j}^{2+j}\sum _{\begin{array}{c}p=0\\ (p+k)odd\end{array}}^{k-1}{F}_{k,2,j}(p+1){\lambda }_{p,s},\hfill \\ \left(5\right){\mathrm{q}}_{j,s}={\int }_0^1\left[D_t^{\nu }{\mathcal{P}}_j\left(t\right)\right]{\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt,\hfill \end{array}$

where

$\begin{array}{cc}\hfill q_{j,s}=& \sum _{k=0}^j{\displaystyle \frac{\pi 4^{k-1}(s+1)\Gamma (k+3){(-1)}^{j+k+s}\Gamma (j+k+2)\Gamma \left(k-\nu +\frac{7}{2}\right)}{\Gamma (2k+2)(j-k)!\Gamma (k-\nu +3)}}\hfill \\ \hfill & \times {}_3\tilde{F}_2\left.\left(\begin{array}{c}-s,s+2,-\nu +k+{\displaystyle \frac{7}{2}}\\ {\displaystyle \frac{3}{2}},-\nu +k+5\end{array}\right|1\right).\hfill \end{array}$

Proof. 

To prove the part (1), the application of Lemma 1 along with the definition of the basis function ${\zeta }_i\left(\rho \right)=\rho (1-\rho ){\mathbf{U}}_i^s\left(\rho \right),$ enables us to write

$\begin{array}{cc}\hfill {\zeta }_i\left(\rho \right)& =\rho (1-\rho ){\mathbf{U}}_i^s\left(\rho \right)\hfill \\ & =\sum _{k=-1-i}^{1+i}{F}_{k,1,i}{U}_k^{*}\left(\rho \right)-\sum _{k=-2-i}^{2+i}{F}_{k,2,i}{U}_k^{*}\left(\rho \right),\hfill \end{array}$

Now, the integration ${\int }_0^1{\zeta }_i\left(\rho \right){\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho$ can be written after using the last relation and the relation (4) as

${\int }_0^1{\zeta }_i\left(\rho \right){\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}\left(\sum _{k=-1-i}^{1+i}{F}_{k,1,i}{\lambda }_{k,r}-\sum _{k=-2-i}^{2+i}{F}_{k,2,i}{\lambda }_{k,r}\right).$

To prove the part (2), the application of Theorem 1 enables us to write ${\int }_0^1{\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}{\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ,$ as

${\int }_0^1{\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}{\mathbf{U}}_r^s\left(\rho \right){\omega }_1\left(\rho \right)d\rho =\sum _{j=0}^i{\mu }_{j,i}{\int }_0^1{\mathbf{U}}_j^s\left(\rho \right){\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ,$

which can be written after using the orthogonality relation (2) as

${\int }_0^1{\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}{\mathbf{U}}_r^s\left(\rho \right)\widehat{w}\left(\rho \right)d\rho ={\displaystyle \frac{\pi }{8}}\sum _{j=0}^i{\mu }_{j,i}{\delta }_{j,r}.$

To prove the part (3), the application of Lemma 1 along with the definition of the basis function ${\mathcal{P}}_j\left(t\right)=t^2{\mathbf{U}}_j^s\left(t\right),$ enables us to write

(11)$\begin{array}{cc}\hfill {\mathcal{P}}_j\left(t\right)& =t^2{\mathbf{U}}_j^s\left(t\right)\hfill \\ & =\sum _{k=-2-j}^{2+j}{F}_{k,2,j}{U}_k^{*}\left(t\right),\hfill \end{array}$

Now, using the last relation and relation (4), we get

${\int }_0^1{\mathcal{P}}_j\left(t\right){\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt={\displaystyle \frac{\pi }{8}}\sum _{k=-2-j}^{2+j}{F}_{k,2,j}{\lambda }_{k,s}.$

To prove part (4), the application of Corollary 1 and Equation (11), one has

${\displaystyle \frac{d{\mathcal{P}}_j\left(t\right)}{dt}}=4\sum _{k=-2-j}^{2+j}\sum _{\begin{array}{c}p=0\\ (p+k)\mathrm{odd}\end{array}}^{k-1}{F}_{k,2,j}(p+1){\mathbf{U}}_p^s\left(t\right).$

Now, using the last relation and relation (4), we get

${\int }_0^1{\displaystyle \frac{d{\mathcal{P}}_j\left(t\right)}{dt}}{\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt={\displaystyle \frac{\pi }{2}}\sum _{k=-2-j}^{2+j}\sum _{\begin{array}{c}p=0\\ (p+k)\mathrm{odd}\end{array}}^{k-1}{F}_{k,2,j}(p+1){\lambda }_{p,s}.$

To find the part (5): Using Equation (1) along with the orthogonality relation (2), one can write

$\begin{array}{cc}\hfill {\int }_0^1\left[D_t^{\nu }{\mathcal{P}}_j\left(t\right)\right]{\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt& =\sum _{k=0}^j{\displaystyle \frac{2^{2k}(k+2)!{(-1)}^{j+k}(j+k+1)!}{(2k+1)!(j-k)!\Gamma (k-\nu +3)}}{\int }_0^1{\mathbf{U}}_s^s\left(t\right)t^{k+2-\nu }\widehat{w}\left(t\right)dt\hfill \\ & =\sum _{k=0}^j{\displaystyle \frac{2^{2k}(k+2)!{(-1)}^{j+k}(j+k+1)!}{(2k+1)!(j-k)!\Gamma (k-\nu +3)}}\sum _{n=0}^s{\displaystyle \frac{2^{2n}{(-1)}^{n+s}(n+s+1)!}{(2n+1)!(s-n)!}}\hfill \\ & \times {\int }_0^1\sqrt{t(1-t)}{t}^{-\nu +k+n+2}dt,\hfill \end{array}$

the last equation may be written after using the relation given below

${\int }_0^1\sqrt{t(1-t)}{t}^{-\nu +k+n+2}dt={\displaystyle \frac{\sqrt{\pi }\Gamma \left(k+n-\nu +\frac{7}{2}\right)}{2\Gamma (k+n-\nu +5)}},$

as

(12)$\begin{array}{cc}\hfill q_{j,s}={\int }_0^1\left[D_t^{\nu }{\mathcal{P}}_j\left(t\right)\right]{\mathbf{U}}_s^s\left(t\right)\widehat{w}\left(t\right)dt& =\sum _{k=0}^j{\displaystyle \frac{2^{2k}(k+2)!{(-1)}^{j+k}(j+k+1)!}{(2k+1)!(j-k)!(-\nu +k+2)!}}\hfill \\ & \times \sum _{n=0}^s{\displaystyle \frac{\sqrt{\pi }{2}^{2n-1}{(-1)}^{n+s}\Gamma (n+s+2)\Gamma \left(k+n-\nu +\frac{7}{2}\right)}{\Gamma (2n+2)(s-n)!\Gamma (k+n-\nu +5)}}.\hfill \end{array}$

If we see the identity listed below:

(13)$\begin{array}{cc}\hfill & \sum _{n=0}^s{\displaystyle \frac{\sqrt{\pi }{2}^{2n-1}{(-1)}^{s+n}(s+n+1)!\Gamma \left(n-\nu +k+\frac{7}{2}\right)}{\Gamma (2n+2)(s-n)!\Gamma (n-\nu +k+5)}}\hfill \\ \hfill =& {\displaystyle \frac{\pi (s+1){(-1)}^s\Gamma \left(k-\nu +\frac{7}{2}\right)}{4}}{\times }_3{\tilde{F}}_2\left.\left(\begin{array}{c}-s,s+2,-\nu +k+{\displaystyle \frac{7}{2}}\\ {\displaystyle \frac{3}{2}},-\nu +k+5\end{array}\right|1\right),\hfill \end{array}$

Now, inserting Equation (13) into Equation (12), we get

$\begin{array}{cc}\hfill q_{j,s}=& \sum _{k=0}^j{\displaystyle \frac{\pi 4^{k-1}(s+1)\Gamma (k+3){(-1)}^{j+k+s}\Gamma (j+k+2)\Gamma \left(k-\nu +\frac{7}{2}\right)}{\Gamma (2k+2)(j-k)!\Gamma (k-\nu +3)}}\hfill \\ \hfill & \times {}_3\tilde{F}_2\left.\left(\begin{array}{c}-s,s+2,-\nu +k+{\displaystyle \frac{7}{2}}\\ {\displaystyle \frac{3}{2}},-\nu +k+5\end{array}\right|1\right).\hfill \end{array}$

□

Remark 3.

The expansion of ${\displaystyle \frac{d^2{\zeta }_i\left(\rho \right)}{d{\rho }^2}}$ in terms of ${\mathbf{U}}_j^s\left(\rho \right)$ is spectrally accurate since ${\mathbf{U}}_j^s\left(\rho \right)$ is complete in the weighted $L^2$ space. Additionally, the sparsity of the coefficients ${\mu }_{j,i}$ results in a sparse matrix $\mathsf{P}$, which significantly increases the method’s computational efficiency.

3.3. Transformation to the $\mathcal{HIBC}s$

Assuming the following $\mathcal{TFDWE}$ [35]:

(14)${\displaystyle \frac{{\partial }^{\nu }\mathcal{Y}(\rho ,t)}{\partial t^{\nu }}}+a{\displaystyle \frac{\partial \mathcal{Y}(\rho ,t)}{\partial t}}+b\mathcal{Y}(\rho ,t)-{\displaystyle \frac{{\partial }^2\mathcal{Y}(\rho ,t)}{\partial {\rho }^2}}=f(\rho ,t),1<\nu <2,$

governed by the following constraints:

(15)$\begin{array}{cc}\hfill & \mathcal{Y}(\rho ,0)=u_0\left(\rho \right),{\mathcal{Y}}_t(\rho ,0)=u_1\left(\rho \right),0\le \rho \le 1,\hfill \end{array}$

(16)$\begin{array}{cc}\hfill & \mathcal{Y}(0,t)=u_3\left(t\right),\mathcal{Y}(1,t)=u_4\left(t\right),0\le t\le 1,\hfill \end{array}$

In virtue of the following transformation:

$\mathrm{Z}(\rho ,t):=\mathcal{Y}(\rho ,t)+\widehat{\mathcal{Y}}(\rho ,t),$

where

$\begin{array}{cc}\hfill \widehat{\mathcal{Y}}(\rho ,t)=& -t\left((\rho -1){\mathcal{Y}}_t(0,0)-\rho {\mathcal{Y}}_t(1,0)+{\mathcal{Y}}_t(\rho ,0)\right)+(\rho -1)\mathcal{Y}(0,t)\hfill \\ & -\rho \mathcal{Y}(1,t)-(\rho -1)\mathcal{Y}(0,0)+\rho \mathcal{Y}(1,0)-\mathcal{Y}(\rho ,0).\hfill \end{array}$

The $\mathcal{TFDWE}$ (14) governed by (15) and (16) is converted into the adapted Equation (5) governed by (6), where

$g(\rho ,t)=f(\rho ,t)+{\displaystyle \frac{{\partial }^{\nu }\widehat{\mathcal{Y}}(\rho ,t)}{\partial t^{\nu }}}+a{\displaystyle \frac{\partial \widehat{\mathcal{Y}}(\rho ,t)}{\partial t}}+b\widehat{\mathcal{Y}}(\rho ,t)-{\displaystyle \frac{{\partial }^2\widehat{\mathcal{Y}}(\rho ,t)}{\partial {\rho }^2}}.$

4. Error Bound

The following Chebyshev-weighted Sobolev space is assumed to conform to

${\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}\left(I_1\right)=\{\eta :\eta \left(0\right)={\eta }^{\prime }\left(0\right)=0\mathrm{and}{D}_t^{\nu +k}\eta \in L_{{\omega }_3\left(t\right)}^2\left(I_1\right),0\le k\le m\},$

${\Psi }_{{\omega }_4\left(\rho \right)}^m\left(I_2\right)=\{\eta :\eta \left(0\right)=\eta \left(1\right)=0\mathrm{and}{D}_{\rho }^k\eta \in L_{{\omega }_4\left(x\right)}^2\left(I_2\right),0\le k\le m\},$

where $I_1=(0,1)$ and $I_2=(0,1)$, equipped with the norm, semi-norm, and inner product

$\begin{array}{c}{(\eta ,v)}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}=\sum _{k=0}^m{(D_t^{\nu +k}\eta ,D_t^{\nu +k}v)}_{L_{{\omega }_3\left(t\right)}^2}{,\left|\right|\eta \left|\right|}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}^2={(\eta ,\eta )}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}{,\left|\eta \right|}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}=\left|\right|D_t^{\nu +m}\eta {\left|\right|}_{L_{{\omega }_3\left(t\right)}^2},\hfill \\ {(\eta ,v)}_{{\Psi }_{{\omega }_4\left(\rho \right)}^m}=\sum _{k=0}^m{(D_{\rho }^k\eta ,D_{\rho }^{k}v)}_{L_{{\omega }_4\left(\rho \right)}^2}{,\left|\right|\eta \left|\right|}_{{\Psi }_{{\omega }_4\left(\rho \right)}^m}^2={(\eta ,\eta )}_{{\Psi }_{{\omega }_4\left(\rho \right)}^m}{,\left|\eta \right|}_{{\Psi }_{{\omega }_4\left(\rho \right)}^m}=\left|\right|D_{\rho }^m\eta {\left|\right|}_{L_{{\omega }_4\left(\rho \right)}^2},\hfill \end{array}$

where $1<\nu <2$ and $m\in \mathbb{N}.$

Remark 4.

The orthogonality relation (4), which includes ${\lambda }_{m,n}$, is very important for finding the elements of the system matrices $\mathsf{M},\mathsf{F},$ and $\mathsf{K}$ that come from the $\mathrm{PGA}$. It makes sure that the projection fits with ${\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}$, which makes the presented numerical technique more accurate and stable.

Now, suppose that the two-dimensional Chebyshev-weighted Sobolev space is represented by

$\begin{array}{cc}\hfill {\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{\ell ,\tau }(I_1\times I_2)=\{\mathrm{Z}:& \mathrm{Z}(\rho ,0)={\mathrm{Z}}_t(\rho ,0)=\mathrm{Z}(0,t)=\mathrm{Z}(1,t)=0\hfill \\ & \mathrm{and}{\displaystyle \frac{{\partial }^{\nu +p+q}\mathrm{Z}}{\partial {\rho }^p\partial t^{\nu +q}}}\in L_{\overline{\omega }(\rho ,t)}^2(I_1\times I_2),\ell \ge p\ge 0,\tau \ge q\ge 0\},\hfill \end{array}$

along with the norm and semi-norm

${\left|\right|\mathrm{Z}\left|\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{\ell ,\tau }}={\left(\sum _{p=0}^{\ell }\sum _{q=0}^{\tau }{\left|\left|{\displaystyle \frac{{\partial }^{\nu +p+q}\mathrm{Z}}{\partial {\rho }^p\partial t^{\nu +q}}}\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}^2\right)}^{{\displaystyle \frac{1}{2}}},{\left|\mathrm{Z}\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{\ell ,\tau }}={\left|\left|{\displaystyle \frac{{\partial }^{\nu +\ell +\tau }\mathrm{Z}}{\partial {\rho }^{\ell }\partial t^{\nu +\tau }}}\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2},$

where $\overline{\omega }(\rho ,t)={\omega }_3\left(t\right){\omega }_4\left(\rho \right),$ and $\ell ,\tau \in \mathbb{N}.$

Lemma 2

([36]). Let $m\ge 1,$ $m+\ell >1$ and $m+\tau >1$, where $\ell ,\tau ,$ are constants. One has

${\displaystyle \frac{\Gamma (m+\ell )}{\Gamma (m+\tau )}}\le {\mathbf{o}}_m^{\ell ,\tau }{m}^{\ell -\tau },$

where

${\mathbf{o}}_m^{\ell ,\tau }=exp\left({\displaystyle \frac{\ell -\tau }{2(m+\tau -1)}}+{\displaystyle \frac{1}{12(m+\ell -1)}}+{\displaystyle \frac{{(\ell -\tau )}^2}{m}}\right).$

Remark 5.

For the given $\ell ,\tau$, ${\mathbf{o}}_m^{\ell ,\tau }$ can be written as follows:

${\mathbf{o}}_m^{\ell ,\tau }=1+O\left(m^{-1}\right).$

Corollary 2.

Suppose $u^{\mathcal{L}}\left(\rho \right)={\sum }_{i=0}^{\mathcal{L}}{\widehat{u}}_i{\zeta }_i\left(\rho \right)$ is the approximate solution of $u\left(\rho \right)\in {\Psi }_{{\omega }_4\left(\rho \right)}^m\left(I_2\right).$ Then, for $0\le k\le m\le \mathcal{L}+1,$ we obtain

$\left|\right|D_{\rho }^k(u\left(\rho \right)-u^{\mathcal{L}}\left(\rho \right)){\left|\right|}_{L_{\widehat{{\omega }_4}\left(\rho \right)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-1}{4}}(m-k)}{\left|u\left(\rho \right)\right|}_{{\Psi }_{{\omega }_4\left(\rho \right)}^m}^2,$

where $a_1\lesssim a_2$ certifies that a constant ν exists such that $a_1\le \nu a_2.$

Proof. 

To derive the proof of this corollary, insert $\ell =1$ in Theorem 4.2 in Ref. [34]. □

Theorem 3.

Suppose $1<\nu <2,$ and ${\displaystyle v^{\mathcal{L}}\left(t\right)=\sum _{j=0}^{\mathcal{L}}{\widehat{\eta }}_j{\mathcal{P}}_j\left(t\right)}$ is the approximate solution of $v\left(t\right)\in {\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}\left(I_1\right).$ Then, for $0\le k\le m\le \mathcal{L}+1,$ we get

(17) $\left|\right|D_t^{\nu +k}(v\left(t\right)-v^{\mathcal{L}}\left(t\right)){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-3}{4}}(m-k)}{\left|v\left(t\right)\right|}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}^2,$

Proof. 

The definitions of $v\left(t\right)$ and $v^{\mathcal{L}}\left(t\right)$ enable us to derive

(18)$\begin{array}{cc}\hfill \left|\right|D_t^{\nu +k}(v\left(t\right)-v^{\mathcal{L}}\left(t\right)){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2& =\sum _{n=\mathcal{L}+1}^{\infty }|{\widehat{\eta }}_n{|}^2\left|\right|D_t^{\nu +k}{\mathcal{P}}_n\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2\hfill \\ & =\sum _{n=\mathcal{L}+1}^{\infty }|{\widehat{\eta }}_n{|}^2{\displaystyle \frac{\left|\right|D_t^{\nu +k}{\mathcal{P}}_n\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}{\left|\right|D_t^{\nu +m}{\mathcal{P}}_n\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}}\left|\right|D_t^{\nu +m}{\mathcal{P}}_n\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2\hfill \\ & \le {\displaystyle \frac{\left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}{\left|\right|D_t^{\nu +m}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}}{\left|v\left(t\right)\right|}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}^2.\hfill \end{array}$

To estimate the factor ${\displaystyle \frac{\left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}{\left|\right|D_t^{\nu +m}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}},$ we find $\left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2.$

$\begin{array}{cc}\hfill \left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2& ={\int }_0^1{D}_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right)D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\omega }_3\left(t\right)dt.\hfill \end{array}$

Equation (1) along with (2) allows us to write

$\begin{array}{c}\hfill D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right)=\sum _{r=k}^{\mathcal{L}+1}{\chi }_{r,\mathcal{L}+1}{\displaystyle \frac{(r+2)!}{\Gamma (r-\nu -k+3)}}{t}^{r+2-k-\nu },\end{array}$

and accordingly, we have

$\begin{array}{cc}\hfill \left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2& =\sum _{r=k}^{\mathcal{L}+1}{\chi }_{r,\mathcal{L}+1}^2{\displaystyle \frac{{(\Gamma (r+3))}^2}{{\Gamma }^2(r-k-\nu +3)}}{\int }_0^1{t}^{2(r-k-\nu )+{\displaystyle \frac{9}{2}}}{(1-t)}^{{\displaystyle \frac{1}{2}}}dt\hfill \\ & =\sum _{r=k}^{\mathcal{L}+1}{\chi }_{r,\mathcal{L}+1}^2{\displaystyle \frac{\sqrt{\pi }{(\Gamma (r+3))}^2\Gamma (2(r-k-\nu )+\frac{11}{2})}{2{\Gamma }^2(r-k-\nu +3)\Gamma (2(r-k-\nu )+5)}}.\hfill \end{array}$

The following inequality can be obtained after using Lemma 2:

${\displaystyle \frac{{(\Gamma (r+3))}^2\Gamma (2(r-\nu -k)+\frac{11}{2})}{{\Gamma }^2(r-k-\nu +3)\Gamma (2(r-k-\nu )+5)}}\lesssim r^{2(k+\nu )}{(r-k)}^{-{\displaystyle \frac{1}{2}}}.$

By virtue of Lemma 2, $\left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2$ can be written as

$\begin{array}{cc}\hfill \left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2& \lesssim {\lambda }^{*}{(1+\mathcal{L})}^{2(k+\nu )}{(1-k+\mathcal{L})}^{-{\displaystyle \frac{1}{2}}}\sum _{r=k}^{\mathcal{L}+1}1\hfill \\ & ={\lambda }^{*}{(1+\mathcal{L})}^{2(k+\nu )}{(1-k+\mathcal{L})}^{-{\displaystyle \frac{1}{2}}}(2-k+\mathcal{L})\hfill \\ & ={\lambda }^{*}{\left({\displaystyle \frac{\Gamma (2+\mathcal{L})}{\Gamma (1+\mathcal{L})}}\right)}^{2(k+\nu )}{\left({\displaystyle \frac{\Gamma (2-k+\mathcal{L})}{\Gamma (1-k+\mathcal{L})}}\right)}^{-{\displaystyle \frac{1}{2}}}\left({\displaystyle \frac{\Gamma (\mathcal{L}-k+3)}{\Gamma (\mathcal{L}-k+2)}}\right)\hfill \\ & \lesssim {\mathcal{L}}^{2(k+\nu )}{(\mathcal{L}-k)}^{{\displaystyle \frac{1}{2}}},\hfill \end{array}$

where ${\chi }^{*}=\underset{0\le r\le 1+\mathcal{L}}{max}\left\{{\displaystyle \frac{{\chi }_{r,\mathcal{L}+1}^2\sqrt{\pi }}{2}}\right\}$.

Similarly, we have

$\left|\right|D_t^{\nu +m}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2\lesssim {\mathcal{L}}^{2(m+\nu )}{(\mathcal{L}-m)}^{{\displaystyle \frac{1}{2}}},$

and accordingly, we have

(19)$\begin{array}{cc}\hfill {\displaystyle \frac{\left|\right|D_t^{\nu +k}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2}{\left|\right|D_t^{\nu +m}{\mathcal{P}}_{\mathcal{L}+1}\left(t\right){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}}}& \lesssim {\mathcal{L}}^{2(k-m)}{\left({\displaystyle \frac{\mathcal{L}-k}{\mathcal{L}-m}}\right)}^{-{\displaystyle \frac{1}{2}}}\hfill \\ & ={\mathcal{L}}^{-2(m-k)}{\left({\displaystyle \frac{\Gamma (\mathcal{L}-k+1)}{\Gamma (\mathcal{L}-m+1)}}\right)}^{{\displaystyle \frac{1}{2}}}\hfill \\ & \lesssim {\mathcal{L}}^{{\displaystyle \frac{-3}{2}}(m-k)}.\hfill \end{array}$

Inserting Equation (19) into Equation (18), one gets

$\left|\right|D_t^{\nu +k}(\eta \left(t\right)-\widehat{\eta }\left(t\right)){\left|\right|}_{L_{{\omega }_3\left(t\right)}^2}^2\lesssim {\mathcal{L}}^{{\displaystyle \frac{-3}{2}}(m-k)}{\left|v\left(t\right)\right|}_{{\Phi }_{{\omega }_3\left(t\right)}^{\nu ,m}}^2.$

Therefore, we obtain the desired result. □

Remark 6.

The bound in inequality (17) remains valid and precise as $\nu \to 1$ or $\nu \to 2$. The spectral coefficients ${\widehat{\eta }}_j$ exhibit a rapid decay, and the solution becomes smoother over time for such values. Therefore, the convergence of the approximation is enhanced in the weighted norm $L_{{\omega }_3\left(t\right)}^2$.

Corollary 3.

Let ${\mathrm{Z}}^{\mathcal{L}}(\rho ,t)$ be the approximate solution of $\mathrm{Z}(\rho ,t)\in {\mathbf{CH}}_{\ddot{\omega }}^{r,s}(\Omega )$; then, for $0\le p\le r\le \mathcal{L}+1,$ the following estimation holds:

${\left|\left|{\displaystyle \frac{{\partial }^p}{\partial {\rho }^p}}(\mathrm{Z}(\rho ,t)-{\mathrm{Z}}^{\mathcal{L}}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-1}{4}}(r-p)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{CH}}_{\overline{\omega }(\rho ,t)}^{r,0}}.$

Proof. 

By virtue of $\mathrm{Z}(\rho ,t)$ and ${\mathrm{Z}}^{\mathcal{L}}(\rho ,t),$ one gets

$\begin{array}{cc}\hfill \mathrm{Z}(\rho ,t)-{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)& =\sum _{i=0}^{\mathcal{L}}\sum _{j=\mathcal{L}+1}^{\infty }{c}_{ij}{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right)+\sum _{i=\mathcal{L}+1}^{\infty }\sum _{j=0}^{\infty }{c}_{ij}{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right)\hfill \\ & \le \sum _{i=0}^{\mathcal{L}}\sum _{j=0}^{\infty }{c}_{ij}{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right)+\sum _{i=\mathcal{L}+1}^{\infty }\sum _{j=0}^{\infty }{c}_{ij}{\zeta }_i\left(\rho \right){\mathcal{P}}_j\left(t\right).\hfill \end{array}$

Now, imitating the same steps as in Theorem 3, one has

${\left|\left|{\displaystyle \frac{{\partial }^p}{\partial {\rho }^p}}(\mathrm{Z}(\rho ,t)-{\mathrm{Z}}^{\mathcal{L}}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-1}{4}}(r-p)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{CH}}_{\overline{\omega }(\rho ,t)}^{r,0}}.$

□

Theorem 4.

Let ${\mathrm{Z}}^{\mathcal{L}}(\rho ,t)$ be the approximate solution of $\mathrm{Z}(\rho ,t)\in {\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{r,s}(\Omega )$. Then, for $0\le q\le s\le \mathcal{L}+1,$ one gets

${\left|\left|{\displaystyle \frac{{\partial }^q}{\partial t^q}}(\mathrm{Z}(\rho ,t)-{\mathrm{Z}}^{\mathcal{L}}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-3}{4}}(s-q)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}},$

Theorem 5.

Let ${\mathrm{Z}}^{\mathcal{L}}(\rho ,t)$ be the approximate solution of $\mathrm{Z}(\rho ,t)\in {\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{r,s}(\Omega )$, and assume that $1<\nu <2$. Consequently, for $0\le q\le s\le \mathcal{L}+1,$ one has

${\left|\left|{\displaystyle \frac{{\partial }^{\nu +q}}{\partial t^{\nu +q}}}({\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-\mathrm{Z}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\lesssim {\mathcal{L}}^{{\displaystyle \frac{-3}{4}}(s-q)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}},$

Proof. 

The proof of Theorems 4 and 5 are in the same manner as the proof of Theorem 3. □

Theorem 6.

Let $\mathrm{R}(\rho ,t)$ be the residual of Equation (8), then the following relation holds

$\begin{array}{cc}\hfill {\left|\left|\mathrm{R}(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}& \lesssim {\mathcal{L}}^{{\displaystyle \frac{-3s}{4}}}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}+a{\mathcal{L}}^{{\displaystyle \frac{-3}{4}}(s-1)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}\hfill \\ & +b{\mathcal{L}}^{{\displaystyle \frac{-3s}{4}}}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}+{\mathcal{L}}^{{\displaystyle \frac{-1}{4}}(r-2)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{r,0}}.\hfill \end{array}$

Proof. 

${\left|\left|\mathrm{R}(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\sigma ,t)}^2}$ of Equation (5) can be written as

$\begin{array}{cc}\hfill {\left|\left|\mathrm{R}(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}& ={\left|\left|{\displaystyle \frac{{\partial }^{\nu }{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial t^{\nu }}}+a{\displaystyle \frac{\partial {\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial t}}+b{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-{\displaystyle \frac{{\partial }^2{\mathrm{Z}}^{\mathcal{L}}(\rho ,t)}{\partial {\rho }^2}}-g(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\hfill \\ & \le {\left|\left|{\displaystyle \frac{{\partial }^{\nu }}{\partial t^{\nu }}}({\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-\mathrm{Z}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}+a{\left|\left|{\displaystyle \frac{\partial }{\partial t}}({\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-\mathrm{Z}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\hfill \\ & +b{\left|\left|({\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-\mathrm{Z}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}+{\left|\left|{\displaystyle \frac{{\partial }^2}{\partial {\rho }^2}}({\mathrm{Z}}^{\mathcal{L}}(\rho ,t)-\mathrm{Z}(\rho ,t))\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}.\hfill \end{array}$

Using Corollary 3, Theorems (4), and 5 lead to

$\begin{array}{cc}\hfill {\left|\left|\mathrm{R}(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}& \lesssim {\mathcal{L}}^{{\displaystyle \frac{-3s}{4}}}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}+a{\mathcal{L}}^{{\displaystyle \frac{-3}{4}}(s-1)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}\hfill \\ & +b{\mathcal{L}}^{{\displaystyle \frac{-3s}{4}}}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{0,s}}+{\mathcal{L}}^{{\displaystyle \frac{-1}{4}}(r-2)}{\left|\mathrm{Z}(\rho ,t)\right|}_{{\mathbf{H}}_{\overline{\omega }(\rho ,t)}^{r,0}}.\hfill \end{array}$

At the end, it is obvious that ${\left|\left|\mathrm{R}(\rho ,t)\right|\right|}_{L_{\overline{\omega }(\rho ,t)}^2}\to 0$ as $\mathcal{L}\to \infty .$ □

5. Examples

Example 1

([35]). Consider the following equation:

$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^{\nu }\mathcal{Y}(\rho ,t)}{\partial t^{\nu }}}-{\displaystyle \frac{{\partial }^2\mathcal{Y}(\rho ,t)}{\partial {\rho }^2}}=f(\rho ,t),\end{array}$

governed by

$\begin{array}{cc}\hfill & \mathcal{Y}(\rho ,0)=0,{\mathcal{Y}}_t(\rho ,0)=-sin\left(\pi \rho \right),0\le \rho \le 1,\hfill \\ & \mathcal{Y}(0,t)=\mathcal{Y}(1,t)=0,0\le t\le 1,\hfill \end{array}$

and $f(\rho ,t)$ is selected such that the analytical solution is $\mathcal{Y}(\rho ,t)=\left(t^2-t\right)sin\left(\pi \rho \right)$.

The absolute errors ($\mathsf{AEs}$) and the approximate solution for $\nu =1.2$ when $\mathcal{L}=6$ are shown in Figure 1. The $\mathsf{AEs}$ of our technique at $\mathcal{L}=7$ and the method in [35] at $\nu =1.8$ when $t=0.2$ are compared in Table 1. Additionally, a comparison of the error norm $L_2$ and $L_{\infty }$ at $\nu =1.9$ is presented in Table 2. These comparisons show that our method is superior to the method in [35]. The precision of the proposed method is demonstrated by the significant decrease in $\mathsf{AEs}$ as $\mathcal{L}$ increases from 4 to 7 when $\nu =1.5$, as illustrated in Figure 2. The proposed procedure is accurate for small choices of $\mathcal{L}$, as illustrated by these results.

Remark 7.

Table 3 illustrates the agreement between the theoretical and numerical results of the error norm $L_2$ presented in Table 2 for $\mathcal{L}=7$. For example, assume $s=8$, $q=0$, and the generic constant $\nu =0.005$ in Theorem 4.

Example 2.

Consider the following equation:

$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^{\nu }\mathcal{Y}(\rho ,t)}{\partial t^{\nu }}}+{\displaystyle \frac{\partial \mathcal{Y}(\rho ,t)}{\partial t}}+\mathcal{Y}(\rho ,t)-{\displaystyle \frac{{\partial }^2\mathcal{Y}(\rho ,t)}{\partial {\rho }^2}}=f(\rho ,t),\end{array}$

governed by

$\begin{array}{cc}\hfill & \mathcal{Y}(\rho ,0)={\mathcal{Y}}_t(\rho ,0)=0,0\le \rho \le 1,\hfill \\ & \mathcal{Y}(0,t)=t^2,\mathcal{Y}(1,t)=t^{2}e,0\le t\le 1,\hfill \end{array}$

and $f(\rho ,t)$ is selected such that the analytical solution is $\mathcal{Y}(\rho ,t)=t^2{e}^{\rho }$.

The $\mathsf{AEs}$ at $\nu =1.5$ and $\mathcal{L}=6$ are displayed in Table 4. The results of the $\mathrm{PGA}$ are extremely near to the precise solution, as illustrated by this table. The precision of the proposed method is shown in Figure 3, where the $\mathsf{AEs}$ decrease significantly as $\mathcal{L}$ increases from 3 to 6 when $\nu =1.9$. The error norms $L_{\infty }$ at several ν and $\mathcal{L}$ values are displayed in Table 5. These results demonstrate our technique’s outstanding efficiency. These findings show that the suggested strategy works well for small choices of $\mathcal{L}$.

Example 3

([35]). Consider the following equation:

$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^{\nu }\mathcal{Y}(\rho ,t)}{\partial t^{\nu }}}+{\displaystyle \frac{\partial \mathcal{Y}(\rho ,t)}{\partial t}}-{\displaystyle \frac{{\partial }^2\mathcal{Y}(\rho ,t)}{\partial {\rho }^2}}=f(\rho ,t),\end{array}$

governed by

$\begin{array}{cc}\hfill & \mathcal{Y}(\rho ,0)={\mathcal{Y}}_t(\rho ,0)=0,0\le \rho \le 1,\hfill \\ & \mathcal{Y}(0,t)=\mathcal{Y}(1,t)=0,0\le t\le 1,\hfill \end{array}$

and $f(\rho ,t)$ is selected such that the analytical solution is $\mathcal{Y}(\rho ,t)=t^2\rho (1-\rho )$.

Table 6 compares the $\mathsf{AEs}$ between our method at $\mathcal{L}=1$ and the approach in [35] at various values of ν when $t=0.1$. Furthermore, Table 7 compares the error norm $L_2$ and $L_{\infty }$ at $\nu =1.9$. These comparisons show that our method outperforms the method in [35]. The $\mathsf{AEs}$ and approximate solution at $\nu =1.3$ when $\mathcal{L}=1$ are finally depicted in Figure 4. This figure demonstrates that the approximate solution is extremely close to the exact solution. These findings show that the suggested method works well for small choices of $\mathcal{L}$.

Example 4

([35]). Consider the following equation:

$\begin{array}{c}\hfill {\displaystyle \frac{{\partial }^{\nu }\mathcal{Y}(\rho ,t)}{\partial t^{\nu }}}+\mathcal{Y}(\rho ,t)-{\displaystyle \frac{{\partial }^2\mathcal{Y}(\rho ,t)}{\partial {\rho }^2}}=f(\rho ,t),\end{array}$

governed by

$\begin{array}{cc}\hfill & \mathcal{Y}(\rho ,0)={\mathcal{Y}}_t(\rho ,0)=0,0\le \rho \le 1,\hfill \\ & \mathcal{Y}(0,t)=0,\mathcal{Y}(1,t)=t^2\sinh \left(1\right),0\le t\le 1,\hfill \end{array}$

and $f(\rho ,t)$ is selected such that the analytical solution is $\mathcal{Y}(\rho ,t)=t^{2}sinh\left(\rho \right)$.

Table 8 compares the $\mathsf{AEs}$ between our method at $\mathcal{L}=6$ and the method in [35] at various values of ν when $t=1$. The $\mathsf{AEs}$ and approximate solution at $\nu =1.4$ when $\mathcal{L}=6$ are depicted in Figure 5. This figure demonstrates that the approximate solution is extremely close to the exact solution. Finally, Table 9 also provides a comparison of the error norm $L_2$ and $L_{\infty }$ for $\nu =1.5$ at various values of t. The findings indicate that the proposed technique is precise for limited selections of $\mathcal{L}$.

Remark 8.

Table 1, Table 2, Table 6, Table 7, Table 8 and Table 9 demonstrate that the findings are precise for small selections of $\mathcal{L}$. Furthermore, these comparisons demonstrate the superior performance of our strategy compared to the method presented in [35].

6. Concluding Remarks

In this work, we present certain $\mathrm{SSKCPs}$-based spectral methods for numerically solving $\mathcal{TFDWE}$. We built an efficient spectral framework to handle this problem by proving operational identities for $\mathrm{SSKCPs}$. The spectral $\mathrm{PGA}$ was utilized for this purpose. The numerical results confirmed the excellent accuracy of our approach. We believe that this basis can be utilized to handle different forms of DEs. In the future, we hope to apply the theoretical findings developed in this study, as well as suitable spectral approaches, to address other FDEs. All codes were written and debugged using Mathematica 11 on an HP Z420 Workstation with an Intel (R) Xeon(R) CPU E5-1620 v2-3.70 GHz, 16 GB RAM DDR3, and 512 GB storage.

Author Contributions

Conceptualization, A.G.A.; Methodology, S.S.A., A.A.A. and A.G.A.; Software, A.G.A.; Validation, S.S.A., A.A.A. and A.G.A.; Formal analysis, S.S.A., A.A.A. and A.G.A.; Investigation, A.G.A.; Data curation, A.G.A.; Writing—original draft, A.A.A. and A.G.A.; Writing—review & editing, S.S.A. and A.G.A.; Visualization, A.G.A.; Supervision, S.S.A. and A.G.A. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Data are contained within the article.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1 The $\mathsf{AEs}$ (left) and approximate solution (right) of Example 1 at $\nu =1.2$ when $\mathcal{L}=6$.

Figure 2 The $\mathsf{AEs}$ of Example 1 at different values of $\mathcal{L}$ when $\nu =1.5$.

Figure 3 The $\mathsf{AEs}$ of Example 2 at different values of $\mathcal{L}$ when $\nu =1.9$.

Figure 4 The $\mathsf{AEs}$ (left) and approximate solution (right) of Example 3 at $\nu =1.3$ when $\mathcal{L}=1$.

Figure 5 The $\mathsf{AEs}$ (left) and approximate solution (right) of Example 4 at $\nu =1.4$ when $\mathcal{L}=6$.

Table 1

Comparison of $\mathsf{AEs}$ of Example 1 at $t=0.2$.

x	Method in [35]	Our Method
0.1	$4.673\times {10}^{-4}$	$4.99358\times {10}^{-10}$
0.2	$6.451\times {10}^{-4}$	$2.27053\times {10}^{-9}$
0.3	$6.474\times {10}^{-4}$	$3.4633\times {10}^{-9}$
0.4	$5.484\times {10}^{-4}$	$2.8527\times {10}^{-10}$
0.5	$3.978\times {10}^{-4}$	$4.1095\times {10}^{-9}$
0.6	$2.229\times {10}^{-4}$	$1.6001\times {10}^{-9}$
0.7	$5.111\times {10}^{-5}$	$3.17794\times {10}^{-9}$
0.8	$1.619\times {10}^{-4}$	$3.39455\times {10}^{-9}$
0.9	$3.037\times {10}^{-4}$	$1.51603\times {10}^{-9}$

Table 2

Comparison of the error norm $L_2$ and $L_{\infty }$ for Example 1 at $\nu =1.9$.

	Method in [35]		Our Method at $\mathcal{L}=7$
$\mathit{t}$	${\mathit{L}}_{\mathbf{2}}$ Error at $\mathit{N}=\mathbf{500}$	${\mathit{L}}_{\infty }$ Error at $\mathit{N}=\mathbf{500}$	${\mathit{L}}_{\mathbf{2}}$ Error	${\mathit{L}}_{\infty }$ Error
0.2	$4.02886\times {10}^{-8}$	$5.69766\times {10}^{-8}$	$2.42083\times {10}^{-9}$	$3.92399\times {10}^{-9}$
0.4	$1.58537\times {10}^{-7}$	$2.24204\times {10}^{-7}$	$9.60276\times {10}^{-9}$	$1.63421\times {10}^{-8}$
0.6	$3.63372\times {10}^{-7}$	$5.13884\times {10}^{-7}$	$2.18266\times {10}^{-8}$	$3.68398\times {10}^{-8}$
0.8	$6.43038\times {10}^{-7}$	$9.09391\times {10}^{-7}$	$3.8955\times {10}^{-8}$	$6.39571\times {10}^{-8}$
1	$9.30169\times {10}^{-7}$	$1.31546\times {10}^{-6}$	$6.10998\times {10}^{-8}$	$9.76926\times {10}^{-8}$

Table 3

Theoretical error of Example 1.

$\mathcal{L}$	7
Error in Theorem 4	${10}^{-8}$

Table 4

$\mathsf{AEs}$ of Example 2 at $\nu =1.5$.

$\mathit{\rho }$	$\mathit{t}=0.2$	$\mathit{t}=0.5$	$\mathit{t}=0.8$
0.1	$2.47968\times {10}^{-11}$	$5.19403\times {10}^{-12}$	$4.69373\times {10}^{-11}$
0.2	$2.25018\times {10}^{-13}$	$1.1339\times {10}^{-11}$	$4.2167\times {10}^{-11}$
0.3	$2.42502\times {10}^{-11}$	$4.687\times {10}^{-12}$	$7.35052\times {10}^{-11}$
0.4	$8.16285\times {10}^{-13}$	$1.85602\times {10}^{-11}$	$6.33866\times {10}^{-11}$
0.5	$2.08427\times {10}^{-11}$	$7.23566\times {10}^{-12}$	$6.52431\times {10}^{-11}$
0.6	$5.34696\times {10}^{-12}$	$2.42187\times {10}^{-11}$	$7.77238\times {10}^{-11}$
0.7	$1.63961\times {10}^{-11}$	$4.10921\times {10}^{-12}$	$5.54976\times {10}^{-11}$
0.8	$4.58489\times {10}^{-12}$	$1.78044\times {10}^{-11}$	$5.76983\times {10}^{-11}$
0.9	$1.18844\times {10}^{-11}$	$1.40076\times {10}^{-11}$	$7.75564\times {10}^{-11}$
CPU time	6.907	6.907	6.907

Table 5

The error norm $L_{\infty }$ of Example 2.

$\mathcal{L}$	2	3	4	5	6
$\nu =1.3$	$3.40553\times {10}^{-5}$	$1.72651\times {10}^{-6}$	$5.03061\times {10}^{-8}$	$1.77561\times {10}^{-9}$	$4.97715\times {10}^{-11}$
CPU time	1.922	2.203	2.829	4.	6.986
$\nu =1.7$	$3.36601\times {10}^{-5}$	$1.70784\times {10}^{-6}$	$5.04312\times {10}^{-8}$	$1.76958\times {10}^{-9}$	$4.31206\times {10}^{-11}$
CPU time	1.845	2.173	2.829	4.094	6.907

Table 6

Comparison of $\mathsf{AEs}$ at $t=0.1$ of Example 3.

	Method in [35] at $\mathit{N}=10$, $\mathbf{\Delta }\mathit{t}=0.0001$		Our Method
$\mathit{\rho }$	$\mathit{\nu }=\mathbf{1}.\mathbf{5}$	$\mathit{\nu }=\mathbf{1}.\mathbf{9}$	$\mathit{\nu }=\mathbf{1}.\mathbf{5}$	$\mathit{\nu }=\mathbf{1}.\mathbf{9}$
0.1	$1.46277\times {10}^{-7}$	$4.23821\times {10}^{-8}$	$2.1684\times {10}^{-19}$	$1.95156\times {10}^{-18}$
0.2	$2.73755\times {10}^{-7}$	$7.78422\times {10}^{-8}$	$4.33681\times {10}^{-19}$	$2.60209\times {10}^{-18}$
0.3	$3.70496\times {10}^{-7}$	$1.03582\times {10}^{-7}$	0	$2.1684\times {10}^{-18}$
0.4	$4.30339\times {10}^{-7}$	$1.19048\times {10}^{-7}$	0	$1.30104\times {10}^{-18}$
0.5	$4.50539\times {10}^{-7}$	$1.24200\times {10}^{-7}$	0	0
0.6	$4.30339\times {10}^{-7}$	$1.19048\times {10}^{-7}$	0	$1.73472\times {10}^{-18}$
0.7	$3.70496\times {10}^{-7}$	$1.03582\times {10}^{-7}$	0	$2.60209\times {10}^{-18}$
0.8	$2.73755\times {10}^{-7}$	$7.78422\times {10}^{-8}$	$4.33681\times {10}^{-19}$	$2.60209\times {10}^{-18}$
0.9	$1.46277\times {10}^{-7}$	$4.23821\times {10}^{-8}$	$2.1684\times {10}^{-19}$	$1.95156\times {10}^{-18}$

Table 7

Comparison of the error norm $L_2$ and $L_{\infty }$ for Example 3 at $\nu =1.9$.

	Method in [35] at $\mathbf{\Delta }\mathit{t}=0.001$, $\mathit{N}=10$		Our Method at $\mathit{N}=1$
$\mathit{t}$	${\mathit{L}}_{\mathbf{2}}$ Errors	${\mathit{L}}_{\infty }$ Errors	${\mathit{L}}_{\mathbf{2}}$ Errors	${\mathit{L}}_{\infty }$ Errors
0.2	$2.08301\times {10}^{-6}$	$2.89129\times {10}^{-6}$	$2.95081\times {10}^{-17}$	$1.04083\times {10}^{-17}$
0.4	$6.27643\times {10}^{-6}$	$8.76628\times {10}^{-6}$	$9.04516\times {10}^{-17}$	$3.46945\times {10}^{-17}$
0.6	$1.37256\times {10}^{-5}$	$1.92499\times {10}^{-5}$	$1.38309\times {10}^{-16}$	$6.93889\times {10}^{-17}$
0.8	$2.42064\times {10}^{-5}$	$3.40618\times {10}^{-5}$	$1.84334\times {10}^{-16}$	$1.11022\times {10}^{-16}$
1	$3.60397\times {10}^{-5}$	$5.08449\times {10}^{-5}$	$2.78356\times {10}^{-16}$	$1.11022\times {10}^{-16}$

Table 8

Comparison of $\mathsf{AEs}$ at $t=1$ of Example 4.

	Method in [35]		Our Method at $\mathcal{L}=6$
$\mathit{\rho }$	$\mathit{\nu }=\mathbf{1.25}$	$\mathit{\nu }=\mathbf{1.75}$	$\mathit{\nu }=\mathbf{1.25}$	$\mathit{\nu }=\mathbf{1.75}$
0.1	$7.07383\times {10}^{-9}$	$5.20813\times {10}^{-9}$	$3.94213\times {10}^{-11}$	$1.55187\times {10}^{-11}$
0.2	$1.34525\times {10}^{-8}$	$1.02913\times {10}^{-8}$	$5.0199\times {10}^{-11}$	$2.13433\times {10}^{-11}$
0.3	$1.91598\times {10}^{-8}$	$1.49610\times {10}^{-8}$	$1.20562\times {10}^{-11}$	$4.07002\times {10}^{-12}$
0.4	$2.38155\times {10}^{-8}$	$1.89463\times {10}^{-8}$	$6.83864\times {10}^{-11}$	$5.17046\times {10}^{-11}$
0.5	$2.69901\times {10}^{-8}$	$2.18956\times {10}^{-8}$	$1.67253\times {10}^{-11}$	$4.46376\times {10}^{-12}$
0.6	$2.81898\times {10}^{-8}$	$2.33512\times {10}^{-8}$	$6.33505\times {10}^{-11}$	$5.90379\times {10}^{-11}$
0.7	$2.68376\times {10}^{-8}$	$2.27226\times {10}^{-8}$	$1.84351\times {10}^{-11}$	$5.44031\times {10}^{-12}$
0.8	$2.22522\times {10}^{-8}$	$1.92641\times {10}^{-8}$	$4.33142\times {10}^{-11}$	$5.01216\times {10}^{-11}$
0.9	$1.36261\times {10}^{-8}$	$1.20571\times {10}^{-8}$	$1.49967\times {10}^{-11}$	$3.31435\times {10}^{-11}$

Table 9

Comparison of the error norm $L_2$ and $L_{\infty }$ for Example 4 at $\nu =1.5$.

	Method in [35]		Our Method at $\mathcal{L}=6$
$\mathit{t}$	${\mathit{L}}_{\mathbf{2}}$ Error	${\mathit{L}}_{\infty }$ Error	${\mathit{L}}_{\mathbf{2}}$ Error	CPU Time	${\mathit{L}}_{\infty }$ Error	CPU Time
0.2	$1.88811\times {10}^{-10}$	$2.72263\times {10}^{-10}$	$1.9977\times {10}^{-12}$	9.986	$2.84641\times {10}^{-12}$	8.017
0.4	$1.20289\times {10}^{-9}$	$1.69968\times {10}^{-9}$	$5.52307\times {10}^{-12}$	11.954	$9.22304\times {10}^{-12}$	8.017
0.6	$3.44298\times {10}^{-9}$	$4.82802\times {10}^{-9}$	$1.27564\times {10}^{-11}$	14.376	$2.20679\times {10}^{-11}$	8.017
0.8	$7.09972\times {10}^{-9}$	$9.93347\times {10}^{-9}$	$2.09841\times {10}^{-11}$	16.704	$3.39929\times {10}^{-11}$	8.032
1	$1.22658\times {10}^{-8}$	$1.71278\times {10}^{-8}$	$4.00609\times {10}^{-11}$	16.704	$6.98624\times {10}^{-11}$	8.032

References

1. Ordokhani, Y.; Sabermahani, S.; Rahimkhani, P. Application of Chelyshkov Wavelets and Least Squares Support Vector Regression to Solve Fractional Differential Equations Arising in Optics and Engineering. Math. Methods Appl. Sci.; 2025; 48, pp. 1996-2010. [DOI: https://dx.doi.org/10.1002/mma.10420]

2. Dehestani, H. Performance of Ritz-Piecewise Gegenbauer Approach for Two Types of Fractional Pantograph Equations Including Piecewise Fractional Derivative. Math. Methods Appl. Sci.; 2025; 48, pp. 6889-6903. [DOI: https://dx.doi.org/10.1002/mma.10724]

3. Jassim, H.; Hussein, G. Analytical Solutions to Fractional Differential Equations: A Comparative Study with the Natural Homotopy Perturbation Method. Proceedings of the 3rd International Conference on Mathematics, AI, Information and Communication Technologies (ICMAICT 2023); Erbil, Iraq, 27–28 April 2025; AIP Conference Proceedings AIP Publishing: Melville, NY, USA, 2025; Volume 3264, 050048.

4. Sayed, S.; Mohamed, A.; Abo El-Dahab, E.; Youssri, Y. Alleviated Shifted Gegenbauer Spectral Method for Ordinary and Fractional Differential Equations. Contemp. Math.; 2024; pp. 1344-1370. [DOI: https://dx.doi.org/10.37256/cm.5220244559]

5. Ramadan, M.; Mansour, M.; El-Shazly, N.; Osheba, H. The Double Ramadan Group Accelerated Adomian Decomposition Method for Solving Nonlinear Partial Differential Equations. Comput. Methods Differ. Equ.; 2025; [DOI: https://dx.doi.org/10.22034/cmde.2025.65368.3000]

6. Vargas, A. Finite Difference Method for Solving Fractional Differential Equations at Irregular Meshes. Math. Comput. Simul.; 2022; 193, pp. 204-216. [DOI: https://dx.doi.org/10.1016/j.matcom.2021.10.010]

7. Nedaiasl, K.; Dehbozorgi, R. Galerkin Finite Element Method for Nonlinear Fractional Differential Equations. Numer. Algorithms; 2021; 88, pp. 113-141. [DOI: https://dx.doi.org/10.1007/s11075-020-01032-2]

8. Maayaha, B.; Bushnaq, S.; Moussaoui, A. Numerical solution of fractional order SIR model of dengue fever disease via Laplace optimized decomposition method. J. Math. Comput. Sci.; 2024; 32, pp. 86-93. [DOI: https://dx.doi.org/10.22436/jmcs.032.01.08]

9. Samy, H.; Adel, W.; Hanafy, I.; Ramadan, M. A Petrov–Galerkin approach for the numerical analysis of soliton and multi-soliton solutions of the Kudryashov–Sinelshchikov equation. Iran. J. Numer. Anal. Optim.; 2024; 14, pp. 1309-1334.

10. Yassin, N.; Atta, A.; Aly, E. Numerical Solutions for Nonlinear Ordinary and Fractional Newell–Whitehead–Segel Equation Using Shifted Schröder Polynomials. Bound. Value Probl.; 2025; 2025, 57. [DOI: https://dx.doi.org/10.1186/s13661-025-02041-7]

11. Atta, A. Approximate Petrov–Galerkin Solution for the Time Fractional Diffusion Wave Equation. Math. Methods Appl. Sci.; 2025; 48, pp. 11670-11685. [DOI: https://dx.doi.org/10.1002/mma.10984]

12. Hafez, R.; Ahmed, H.; Alqubori, O.; Amin, A.; Abd-Elhameed, W. Efficient Spectral Galerkin and Collocation Approaches Using Telephone Polynomials for Solving Some Models of Differential Equations with Convergence Analysis. Mathematics; 2025; 13, 918. [DOI: https://dx.doi.org/10.3390/math13060918]

13. Taema, M.; Dagher, M.; Youssri, Y. Spectral Collocation Method via Fermat Polynomials for Fredholm–Volterra Integral Equations with Singular Kernels and Fractional Differential Equations. J. Math.; 2025; 14, pp. 481-492.

14. Salamaa, M.H.; Zedan, H.A.; Abd-Elhameed, W.M.; Youssri, Y.H. Galerkin Method with Modified Shifted Lucas Polynomials for Solving the 2D Poisson Equation. J. Comput. Appl. Mech.; 2025; 56, pp. 737-775.

15. Abdelhakem, M.; Abdelhamied, D.; El-Kady, M.; Youssri, Y. Two Modified Shifted Chebyshev–Galerkin Operational Matrix Methods for Even-Order Partial Boundary Value Problems. Bound. Value Probl.; 2025; 2025, 34. [DOI: https://dx.doi.org/10.1186/s13661-025-02021-x]

16. Abd-Elhameed, W.M.; Alkhamisi, S.O.; Amin, A.K.; Youssri, Y.H. Numerical contrivance for Kawahara-type differential equations based on fifth-kind Chebyshev polynomials. Symmetry; 2023; 15, 138. [DOI: https://dx.doi.org/10.3390/sym15010138]

17. Yağmurlu, N.M.; Karakaş, A.S. A novel perspective for simulations of the MEW equation by trigonometric cubic B-spline collocation method based on Rubin-Graves type linearization. Comput. Methods Differ. Equ.; 2022; 10, pp. 1046-1058.

18. Kutluay, S.; Yağmurlu, N.M.; Karakaş, A.S. A novel perspective for simulations of the Modified Equal-Width Wave equation by cubic Hermite B-spline collocation method. Wave Motion; 2024; 129, 103342. [DOI: https://dx.doi.org/10.1016/j.wavemoti.2024.103342]

19. Kutluay, S.; Yağmurlu, N.M.; Karakaş, A.S. A robust septic Hermite collocation technique for Dirichlet boundary condition heat conduction equation. Int. J. Math. Comput. Eng.; 2025; 3, pp. 253-266. [DOI: https://dx.doi.org/10.2478/ijmce-2025-0019]

20. Luo, M.; Qiu, W.; Nikan, O.; Avazzadeh, Z. Second-order accurate, robust and efficient ADI Galerkin technique for the three-dimensional nonlocal heat model arising in viscoelasticity. Appl. Math. Comput.; 2023; 440, 127655. [DOI: https://dx.doi.org/10.1016/j.amc.2022.127655]

21. Gottlieb, D.; Orszag, S. Numerical Analysis of Spectral Methods: Theory and Applications; SIAM: Philadelphia, PA, USA, 1977.

22. Brezinski, C. Padé-Type Approximation and General Orthogonal Polynomials; Springer: Berlin/Heidelberg, Germany, 1980; Volume 50.

23. Funaro, D. Polynomial Approximation of Differential Equations; Springer: Berlin/Heidelberg, Germany, 2008; Volume 8.

24. Gautschi, W. Orthogonal Polynomials: Computation and Approximation; Oxford University Press: Oxford, UK, 2004.

25. Gamal, M.; Zaky, M.; El-Kady, M.; Abdelhakem, M. Chebyshev Polynomial Derivative-Based Spectral Tau Approach for Solving High-Order Differential Equations. Comput. Appl. Math.; 2024; 43, 412. [DOI: https://dx.doi.org/10.1007/s40314-024-02908-y]

26. Abd-Elhameed, W.; Youssri, Y. New Formulas of the High-Order Derivatives of Fifth-Kind Chebyshev Polynomials: Spectral Solution of the Convection–Diffusion Equation. Numer. Methods Partial Differ. Equ.; 2024; 40, e22756. [DOI: https://dx.doi.org/10.1002/num.22756]

27. Sadiq, S.; ur Rehman, M. Numerical Technique Based on Generalized Laguerre and Shifted Chebyshev Polynomials. J. Appl. Anal. Comput.; 2024; 14, pp. 1977-2001. [DOI: https://dx.doi.org/10.11948/20220504]

28. Brahim, M.S.T.; Youssri, Y.H.; Alburaikan, A.; Khalifa, H.; Radwn, T. A Refined Galerkin Approach for Solving Higher-Order Differential Equations via Bernoulli Polynomials. Fractals; 2025; 2540183, 14.

29. Adebisi, A.; Uwaheren, O.; Oseni, W.; Ishola, C.; Peter, O. Numerical Solution of Fractional Order Differential Equations by Chebyshev Least Squares Approximation Method. Euler J. Ilm. Mat. Sains Teknol.; 2025; 13, pp. 38-44. [DOI: https://dx.doi.org/10.37905/euler.v13i1.31034]

30. Abd-Elhameed, W.; Ahmed, H. Spectral Solutions for the Time-Fractional Heat Differential Equation Through a Novel Unified Sequence of Chebyshev Polynomials. AIMS Math.; 2024; 9, pp. 2137-2166. [DOI: https://dx.doi.org/10.3934/math.2024107]

31. Ahmed, H.; Hafez, R.; Abd-Elhameed, W. A Computational Strategy for Nonlinear Time-Fractional Generalized Kawahara Equation Using New Eighth-Kind Chebyshev Operational Matrices. Phys. Scr.; 2024; 99, 045250. [DOI: https://dx.doi.org/10.1088/1402-4896/ad3482]

32. Podlubny, I. Fractional Differential Equations: An Introduction to Fractional Derivatives, Fractional Differential Equations, to Methods of Their Solution and Some of Their Applications; Elsevier: Amsterdam, The Netherlands, 1998; Volume 198.

33. Abd-Elhameed, W.; Youssri, Y. Explicit Shifted Second-Kind Chebyshev Spectral Treatment for Fractional Riccati Differential Equation. Comput. Model. Eng. Sci.; 2019; 121, pp. 1029-1049. [DOI: https://dx.doi.org/10.32604/cmes.2019.08378]

34. Abd-Elhameed, W.; Al-Sady, A.; Alqubori, O.; Atta, A. Numerical Treatment of the Fractional Rayleigh-Stokes Problem Using Some Orthogonal Combinations of Chebyshev Polynomials. AIMS Math.; 2024; 9, pp. 25457-25481. [DOI: https://dx.doi.org/10.3934/math.20241243]

35. Shafiq, M.; Abbas, M.; El-Shewy, E.; Abdelrahman, M.; Abdo, N.; El-Rahman, A. Numerical Investigation of the Fractional Diffusion Wave Equation with the Mittag–Leffler Function. Fractal Fract.; 2023; 8, 18. [DOI: https://dx.doi.org/10.3390/fractalfract8010018]

36. Zhao, X.; Wang, L.L.; Xie, Z. Sharp Error Bounds for Jacobi Expansions and Gegenbauer–Gauss Quadrature of Analytic Functions. SIAM J. Numer. Anal.; 2013; 51, pp. 1443-1469. [DOI: https://dx.doi.org/10.1137/12089421X]