Spectral tau solution of the linearized time-fractional KdV-Type equations

1.   Introduction

Fractional differential equations (FDEs) are increasingly being used in several fields of research and engineering applications. Because of their extensive applications in the applied sciences, FDEs have sparked a lot of attention in recent years by several researchers. FDEs are well suited to represent a wide range of phenomena in electromagnetics, viscoelasticity, fluid mechanics, solid mechanics, electrochemistry, biological population models, and signal processing (see, for example ^[1,2]). Obtaining exact solutions for FDEs is not available in most cases, so the employment to the different numerical methods becomes necessary to propose effective numerical solutions to treat the different FDEs. In this regard, wavelets methods are used to treat some FDEs, see ^[3,4]. Operational matrix methods are also utilized in a variety of papers, see for example ^[5,6,7]. For some other numerical methods that used for treating different FDEs, one can be referred for example to ^[8,9,10].

Chebyshev polynomials have garnered a lot of attention in recent decades, from both theoretical and practical aspects. There are well-known four kinds of Chebyshev polynomials that are special ones of the Jacobi polynomials. The first and second kinds are symmetric polynomials, and they are used extensively to numerically treat the different types of differential equations. For example, Chebyshev polynomials of the first kind were employed along with the spectral Galerkin method in ^[11] to treat numerically the linear and non-linear hyperbolic telegraph type equations. In addition, numerical solutions of systems of certain FDEs were proposed in ^[12] based on utilizing shifted Chebyshev polynomials of the first kind. Regarding Chebyshev polynomials of the second kind, they also were used in a variety of papers. For example, the authors in ^[13] employed these polynomials together with the collocation method to solve the space fractional advection–dispersion equation. In ^[14], the authors followed an operational matrix technique based on the second-kind Chebyshev polynomials to solve the variable order fractional differential-integral equation. In addition, the authors in ^[15] used the shifted second-kind Chebyshev polynomials to obtain spectral solutions of the fractional Riccati differential equation. Moreover, Adomian decomposition method along with the Laguerre polynomials and the second kind Chebyshev polynomials were used in ^[16] to treat some types of differential equations. Regarding Chebyshev polynomials of the third and fourth kinds, they also were utilized in solving different types of differential equations. For some contributions regarding these polynomials, one can be referred to ^[17,18,19]. Other two kinds of Chebyshev polynomials namely, Chebyshev polynomials of the fifth and sixth kinds that were investigated in ^[20] were also utilized in a variety of papers, see for example ^[21,22].

Spectral methods are one of the most extensively used numerical strategies for solving various types of differential equations. There are some advantages of these methods if compared with other numerical methods. For example, these methods are global methods unlike the finite element methods. The main characteristic of the spectral method is based on selecting two families of basis functions, namely, trial functions and test functions. These two families of functions are often expressed in terms of suitable orthogonal polynomials or combinations of them. The choice of the trial and test functions depend on the method we choose. It is well-known that there are three main categories of spectral methods, namely Galerkin, collocation and tau methods. In Galerkin method, the test functions are the same as the trial functions and they are chosen such that each member of them satisfies the boundary/initial conditions governed by the given differential equation. Regarding the tau method, it differs from the Galerkin method that no restrictions on choosing the basis functions. This of course makes its application to different types of differential equations is easier than the application of the Galerkin method. So, as a result it is used for solving several types of differential equations. In this concern, the authors in the two papers ^[21,23] used respectively Chebyshev polynomials of the fifth- and sixth-kinds to treat multi-term fractional differential equations. Abd-Elhameed in ^[24] derived expressions for the high-order derivatives of Chebyshev polynomials of the sixth-kind and he employed some specific derivatives formulas along with the application of the tau method to treat numerically the non-linear one-dimensional Burgers' equation. The authors in ^[25] proposed a spectral tau algorithm for handling certain coupled system of FDEs based on selecting generalized Fibonacci polynomial as basis functions. The collocation method can be applied to any type of differential equation, see for example ^[26].

The Korteweg–De Vries (KdV) type equations and their modified ones appear in many applications. There are many types of KdV equations. The third-order KdV equations is an equation that depicts the behavior of one-dimensional shallow water waves of small but finite amplitude, see for example ^[27,28]. In ^[29], a modified third-order KdV equation was treated based on the application of the Petrov-Galerkin finite element method. Fifth-order KdV equations were investigated in some articles. For example, the authors in ^[30] proposed a numerical approach for handling the fifth-order KdV equations. Another approach based on the decomposition method is followed in ^[31] to treat the fifth-order KdV equation. Regarding the seventh-order KdV equations, some explicit solutions of KdV-Burgers' and Lax's seventh-order KdV equations were presented in ^[32]. Other types of fractional KdV-equations were investigated by some authors. For example, a type of fractional KdV equation was treated in ^[33] using fractional natural decomposition method. A Green's function was employed in ^[34] for the fractional KdV equation. A Galerkin method is followed in ^[35] for the time fractional KdV equation with weak singularity solution.

The numerical solutions of the time FPDEs with a second-order partial derivative were proposed by many authors, see for example ^[36,37,38], however, the numerical investigations of the time FPDEs with a third-order partial derivative are not enough. This motivates our interest in investigating such problems. In ^[39], the authors proposed a Petrov-Galerkin spectral method for treating the linearized time fractional KdV equation. There are some methods for treating some fractional KdV equations, see for example ^[40,41].

The rest of the article can be organized as follows. The next section displays some properties of Chebyshev polynomials of the second kind and their shifted ones. Moreover, some fundamental definitions of the fractional calculus theory are stated. Section 3 is confined to proposing a numerical strategy to treat numerically the linearized KdV equation based on the application of the spectral tau method. Investigation of the convergence analysis of the double second kind Chebyshev approximation is discussed in Section 4. Moreover, some estimates concerning the truncation and global errors are also introduced in this section. Numerical discussions are introduced via presenting three illustrative examples accompanied with some comparisons in Section 5 to test the applicability and accuracy of the presented method. Finally, some conclusions are presented in Section 6.

2.   Some fundamentals and useful formulas

This section is confined to presenting some formulas concerned with Chebyshev polynomials of the second kind $\mathcal{U}_{m}(x), \, m\ge 0$, and their shifted ones. Moreover, some standard definitions of the fractional calculus theory will be presented. Additionally, some important formulas concerned with certain combinations of the second-kind Chebyshev polynomials that will be useful in our study in the upcoming sections are also presented in this section.

2.1.   Some fundamental properties of Chebyshev polynomials of the second kind

The second-kind Chebyshev polynomials $\mathcal{U}_{m}(z), \, z\in [-1, 1]$, (see ^[42]) are special ones of the Jacobi polynomials. They can be also defined by the following trigonometric representation:

In addition, they may be constructed by means of the following recursive formula:

The orthogonality relation of $\mathcal{U}_{m}(z), \, m\ge 0$ on the interval $[-1, 1]$ is given by

Several properties and relations concerned with $\mathcal{U}_{m}(x)$ can be found in ^[42].

Now, the shifted second-kind Chebyshev polynomials $\mathcal{U}^*_{m}(z)$ on $[0, 1]$ are defined as: $\mathcal{U}_{m}^{*}(z) = \mathcal{U}_{m}(2z-1)$. These shifted polynomials may be constructed with the aid of the following recursive formula:

It is well-known that the family of the shifted polynomials $\left\{\mathcal{U}^{*}_{m}(t)\right\}_{m\ge 0}$ is orthogonal on $[0, 1]$ in the sense that:

2.2.   Some definitions and properties of fractional calculus

In the following, we display some fundamental definitions and relations of the fractional calculus theory, (see ^[43]).

Definition 2.1. On the standard Lebesgue space $L_1[0, 1]$, the Riemann-Liouville fractional integral operator $I^{\gamma}$ of order $\gamma$ is defined as

The following identities are satisfied by this operator:

with $\gamma, \delta\ge 0,$ and $\nu > -1.$

Definition 2.2. The Riemann-Liouville fractional derivative of order $\gamma > 0$ is defined by

Definition 2.3. On the standard Lebesgue space $L_1[0, 1]$, the Riemann-Liouville fractional integral operator $_0I tnu$ of order $\nu$ is defined as: for all $t\in (0, 1)$

3.   Spectral solution of the linearized KdV equation

This section is devoted to developing spectral solutions for the linearized time fractional KdV equation. In fact, we will employ a double Chebyshev second-kind basis functions to solve the equation. Before we proceed in our proposed algorithm, the following formulas concerned with the Chebyshev polynomials of the second kind and their shifted ones are needed in the sequel.

3.1.   Useful formulas of Chebyshev polynomials and certain related polynomials

In this section, we will present some interesting formulas of the Chebyshev polynomials of the second-kind, their shifted ones, and some of their related polynomials. These formulas will be essential in the derivation of our proposed algorithm.

The following theorem gives an approximation to the fractional derivatives of the shifted second-kind Chebyshev polynomials in terms of their original ones.

Theorem 3.1. The following formula can be used to approximate the fractional derivatives of the polynomials $\mathcal{U}^{*}_{k}(t)$ as (^[15]):

where $n = \lceil \alpha \rceil.$ is the well-known ceiling notation, and $M_{\alpha} > > n$ is a sufficiently large positive integer.

Remark 3.1. Note that the generalized hypergeometric function that appears in the approximating formula (3.1) is defined as (^[44])

where $c_{1}, \, c_{2}, \dots, c_{r}, \, d_{1}, \, d_{2}, \dots, d_{s},$ are complex or real parameters, where the constants $d_{i}$ are all neither zeros nor negative integers for all $1\le i\le s$.

Now, the following lemma exhibits an expression for the the third-order derivative of the polynomials $\mathcal{U}_n(x)$ in terms of their original ones.

Lemma 3.1. For every integer $n$ with $n\ge 3$, one has the following third-order derivatives formula:

where the coefficients $\bar{E}_{n, i}$ can be computed explicitly by the following formula:

Proof. Since the Chebyshev polynomials of the second kind $\mathcal{U}_{n}(x)$ are special ones of the ultraspherical polynomials $C^{(\lambda)}_{n}(x)$ (see ^[45]), so their derivatives can be obtained as special cases of the derivatives of $C^{(\lambda)}_{n}(x)$. The $qth$-derivative of $C^{(\lambda)}_{n}(x), \ q\ge 1$ is given by (see ^[45]):

thus noting the well-known relation: $\mathcal{U}_{k}(x) = (k+1)\, C^{(1)}_{k}(x)$, it is easy to compute the $qth$-derivative of $\mathcal{U}_{k}(x)$ as

where

Setting $q = 3$ in (3.4) yields the following formula

where the coefficients $\bar{E}_{n, i}$ are given by:

Lemma 3.1 is now proved.

Lemma 3.2. If we let $\phi_{i}(x) = (1+x)(1-x)^2\ \mathcal{U}_{i}(x), \ i$ is any non-negative integer, then the following formula holds:

Proof. Based on the recurrence relation (2.1), it is easy to express the three terms $x\, \mathcal{U}_{i}(x), \ x^2\, \mathcal{U}_{i}(x)$ and $x^3\, \mathcal{U}_{i}(x)$ in the following forms:

and therefore, it is easy to see that formula (3.6) holds.

Lemma 3.3. Let $i$ and $j$ be two non-negative integers. The following integral formula holds:

with $w_{1}(x) = w_{1} = \sqrt{1-x^2}$, and

Proof. The integral in (3.7) can be easily computed by making use of formula (3.6) along with the orthogonality relation (2.2).

The following theorem gives an explicit expression for the third-order derivative of the polynomials $\phi_{i}(x)$ in terms of the second kind Chebyshev polynomials $\mathcal{U}_{i}(x)$.

Theorem 3.2. For every non-negative integer $i$, one has the following identity:

where the coefficients $d_{k, i}$ are given by the following formula:

Proof. First, we make use of Lemma 3.2 to express $\phi_{i}(x)$ in terms of $U_{i}(x)$. We will differentiate both sides of (3.6) to obtain the following relation:

Now, in virtue of Lemma 3.1, and after performing some algebraic computations, it can be shown that the following identity holds:

Some computations enable us to write the last formula in the following alternative form:

where the coefficients $d_{k, i}$ are given by the following formula

3.2.   Analyzing the proposed Galerkin algorithm

In this section, we consider the following linearized time fractional KdV equation of the form ^[39]:

subject to the following initial and boundary conditions:

where ${ }_{0}^{C} D_{t}^{\alpha} u(x, t), 0 < \alpha < 1$, is defined as the Caputo fractional derivatives of order $\alpha$ with respect to time $t$, given by ^[39]

Now, in order to propose our numerical algorithm, we choose the following two kinds of basis functions:

We suggest an approximate solution to Eq (3.13) as:

Moreover, assume the following double series approximation for $f(x, t)$

Now, in order to be able to apply a suitable spectral method, we first compute the residual of Eq (3.13). If we make use of the proposed approximation of $u_N(x, t)$ in (3.18) together with Theorems 3.1 and 3.2, and formula (3.19), then the residual of Eq (3.13) is given by the following formula:

Let us agree on the following weight functions,

Also, let

Now, to get the expansion coefficients $\hat{u}_{ij}$, we apply the spectral tau method, by letting the inner product between $\mathcal{R}(x, t)$ and $\Phi_{rs}$ with respect to the wight $w(x, t)$ be zero, that is, we get

The residual $\mathcal{R}(x, t)$ in (3.20) enables one to convert Eq (3.21) into the following form:

Now, based on the two orthogonality relations of $\mathcal{U}_{n}(x)$ and $\mathcal{U}^*_{n}(x)$, (3.22) can be converted into

where $\delta_{ir}$ denotes the well-known Kronecker Delta function.

Now, we choose $M_{\alpha}\geq N$, an accordingly, we will get

that can be reduced into the following relation

Now, since the choice of the basis functions $\phi_{i}(x)$ guarantees that the boundary conditions (3.15) are satisfied. This means that the initial condition (3.14) should be set as an additional constraint.

Now assume that $u_{0}(x)$ can be approximated in terms of $\mathcal{U}_i(x)$, that is we can write

and therefore, we have

Taking the inner product for both sides of (3.24) with $\mathcal{U}_{r}$ with respect to the wight $w_1(x)$, we get

Merging the last relation together with (3.23), the following system of equations in the unknown expansion coefficients $\left(\hat{u}_{ij}\right)_{0\le i, j\le N}$ can be obtained:

In virtue of Lemma 3.3, the last system can be turned into the following one

where $\epsilon_{ij}$ are given by (3.8). Equation (3.25) reduces the model problems (3.13–3.15) into a neat explicit system of algebraic equations in the $(N+1)^2$ unknown expansion coefficients $\{\hat{u}_{r, s}: (r, s)\in \{0, 1, \dots, N\}\times\{0, 1, \dots, N\}\}$ of dimension $N(N+1)+N+1 = (N+1)^2$, this system can be efficiently inverted via the Gaussian-elimination criteria to get the expansion coefficients $\hat{u}_{r, s}$ with ease.

Remark 3.2. It is worth mentioning here that the significance of the numerical algorithm that is described in Subsection 3.2 is due to the following two points:

● The proposed method is novel. It depends on using new expressions for the fractional and integer derivatives for certain basis functions.

● The method reduces the linearized time-fractional KdV equation (3.13) governed by the underlying conditions (3.14) and (3.15) into a system of Eqs (3.25) that can be efficiently solved.

4.   Investigation of the convergence analysis

This section is confined to investigating the convergence and error analysis of the proposed double expansion, Additionally, some estimates regarding the truncation and global errors are also presented. Let us first agree on the following notation; if $A_n, B_n$ are two real-valued sequences, by writing $A_n \lesssim B_n$, we mean that there exists a generic positive constant $K$ such that $A_n\leq K\, B_n$. The following results are needed.

Lemma 4.1. The following estimates are valid:

(a) ^[42] For all $i\geq0$, we have the following estimate

(b) ^[15] For all $j\geq0$, we have the following estimate

(c) ^[15] If $y(z)\in C^{q}[-1, 1]$ for some $q > 3$ is approximated as $y_L(z) = \sum\limits_{\ell = 0}^L\hat{y}_{\ell}\, \mathcal{U}_{\ell}(z)$, then we have the following estimate:

(d) ^[15] Under the same assumptions in (c), the following estimate is valid

Based on the above lemma, we are ready to prove the following important theorem, in which, we find and estimate for the following:

● The unknown expansion coefficients, to ensure the convergence of the approximate solution.

● The truncation error to ascertain the convergence of the approximate solution to the exact one.

● An estimate for the Caputo temporal derivative of the truncated series.

● An estimate for the third space derivative of the truncated series.

● The global error estimate.

Theorem 4.1. If $u$ and $u_N$ are respectively the exact and spectral-approximate solutions of (3.13)–(3.15), and if in addition $u\in C^4(\Lambda)$, where $\Lambda = (-1, 1)\times(0, T)$ is separable, i.e., $u = u_1(x)\, u_2(t)$, then the following estimates are justifiable:

(A) $|\hat{u}_{i, j}|\lesssim i^{-n_1}\, j^{-n_2}$, for some $n_1, n_2 > 4$.

(B) $|u-u_{N}|\lesssim N^{1-min(n_1, n_2)}$, for the same $n_1, n_2$ in (A).

(C) $|{ }_{0}^{C} D_{t}^{\alpha}\, (u-u_{N})|\lesssim N^{1+\alpha-min(n_1, n_2)}$, \ $0 < \alpha\leq1$, for the same $n_1, n_2$ in (A).

(D) $|D_{xxx}\, (u-u_{N})|\lesssim N^{4-min(n_1, n_2)}$, for the same $n_1, n_2$ in (A).

(E) $|{ }_{0}^{C} D_{t}^{\alpha} u_N+u^{(3, 0)}_N-f|\lesssim N^{1+\alpha-min(n_1, n_2)}.$

Proof.

Part (A)

If we note that

that can expressed in the following form:

from the separability hypothesis, we can write

By Lemma 4.1(c), part (A) can be proved.

Part (B)

Now, noting the simple inequality: $(1+x)(1-x)^2\leq\frac{32}{27}$, by repeated use of Lemma 4.1(a), and by the result of part (A), part (B) is proved.

Part (C)

Starting with,

By the application of Lemma 4.1(b), followed by (a), and the use of part (A), part (C) can be proved.

Part (D)

The proof is similar to (C), but we take into consideration that $|D_{xxx}\, \mathcal{U}_i(x)|\lesssim i^4$.

Part (E)

Noting that ${ }_{0}^{C} D_{t}^{\alpha} u+u^{(3, 0)} = f$, we have

By the triangle inequality and the application of (C), and (D), part (E) is proved.

5.   Numerical results

In this section, we present three numerical examples accompanied with some comparisons to show the accuracy and applicability of our proposed algorithm.

Example 5.1. Consider the following linearized KdV equation (^[39]):

subject to the following homogeneous initial and boundary conditions:

with the exact solution

and, the corresponding compatible forcing term

In Figures 1 and 2, respectively, we depict the LogErrors as well as the approximate solution for different values of $\alpha$. In Table 1, we report the maximum point-wise error of Example 1 for some values of $\alpha, t$ and $N$.

Example 5.2. Consider the following linearized KdV equation (^[39]):

subject to the following homogeneous initial and boundary conditions:

with the following exact solution with a limited regularity

and the corresponding compatible forcing term

In Figures 3 and 4, respectively, we depict the LogErrors as well as the approximate solution for different values of $\alpha$.

Example 5.3. Consider the following linearized KdV equation (^[39]):

subject to the following homogeneous initial and boundary conditions:

with the exact solution

and the corresponding compatible forcing term

In Figures 5 and 6, respectively, we depict the LogErrors as well as the approximate solution for different values of $\alpha$.

Remark 5.1. For the sake of comparisons for our numerical results obtained from the application of our proposed algorithm with those obtained from the application of the algorithm presented in Ref. [39], Table 2 displays a comparison between the orders of the best errors resulted from the two methods in the presented three examples. The results of this table show that our method is more accurate than the method developed in ^[39].

6.   Conclusions

Herein, we have developed a numerical method for the numerical treatment of a type of FDEs, namely, the linearized time-fractional KdV equation. An approximate solution was expressed as a double expansion in terms of Chebyshev polynomial of the second-kind. The tau method is applied to obtain the desired numerical solution. Some theoretical results served in the derivation of our desired approximate solutions. A formula that approximates the fractional derivatives of the shifted Chebyshev polynomials of the second kind was utilized. Moreover, a third-order derivatives formula of a certain basis functions that are expressed as a certain combination of the Chebyshev polynomials of the second kind was given in terms of the second kind Chebyshev polynomials themselves. The algorithm was tested through presenting some examples to ensure the the high accuracy and efficiency of the presented numerical algorithm. We do believe that our numerical algorithm can be applied to treat some other similar FDEs involving non-linear terms by performing some modifications. All codes in this paper were written and debugged by Mathematica 11 on HP Z420 Workstation, Processor: Intel (R) Xeon(R) CPU E5-1620-3.6 GHz, 16 GB Ram DDR3, and 512 GB storage.

Acknowledgements

This work was funded by the University of Jeddah, Jeddah, Saudi Arabia, under grant No. (UJ-21-DR-141). The authors, therefore, acknowledge with thanks the University of Jeddah technical and financial support.

Conflict of interest

The authors declare that they have no competing interests.