Spectral tau technique via Lucas polynomials for the time-fractional diffusion equation

1.   Introduction

Both integrals and fractional derivatives are non-local operators because they include integration, which is non-local. Because of this property, these operators help describe certain physical phenomena, such as hereditary features, asymptotic scaling, and long-term memory effects. Fractional calculus attracted many mathematicians, who worked hard to develop the subject further and devised distinctive methods of defining fractional-order integrals and derivatives. Many researchers are starting to pay more attention to the idea of fractional calculus because of all the places it might be useful, including fluid dynamics, electrochemistry of corrosion, biology, optics, signal processing, fluid dynamics, regular variation in thermodynamics, aerodynamics, and many more ^[1,2,3]. Researchers conduct extensive research using various numerical algorithms for fractional differential equations (FDEs). For example, in ^[4], the Adomian decomposition method was used. In ^[5], a collocation algorithm was used for a specific type of Korteweg–De Vries equation (KdV) equation. In ^[6], the authors applied a physics-informed neural network-based scheme for treating some FDEs. The authors of ^[7] used a collocation scheme and certain unified Chebyshev polynomials to solve the fractional heat equations. In ^[8], a matrix algorithm was applied to solve certain FDEs. A Haar wavelets method was followed to solve some pantograph FDEs in ^[9]. An approach based on an inverse Laplace transform was followed in ^[10] to handle some FDEs. The authors of ^[11] used meshless analysis and the generalized finite difference method to treat the fractional diffusion-wave equation. The same authors performed another study in ^[12].

Lucas polynomials and their generalized sequences are important in many fields, such as number theory, statistics, and computer science; see, for example, ^[13]. In numerical analysis, the roles of Lucas polynomials and their generalized polynomials have increased. They are utilized in solving differential equations (DEs) of all types. For example, the authors of ^[14] treated the time-fractional Burgers equation using a finite difference approach. In ^[15], a numerical scheme was proposed to treat the advection-diffusion equations using mixed Fibonacci Lucas polynomials. In ^[16], a spectral method was utilized to solve an electro-hydrodynamics flow model. In ^[17], the authors employed Lucas polynomials to solve certain multi-dimensional equations. The authors of ^[18] treated Cauchy integral equations based on Lucas polynomials. In ^[19], the authors used modified Lucas polynomials to treat some FDEs. Other multi-dimensional problems were treated via Lucas polynomials in ^[20]. High-order boundary value problems were handled using Lucas polynomials in ^[21].

To explain out-of-the-ordinary diffusion behaviors seen in a wide range of biological, financial, and physical systems, the TFDE adds a fractional time derivative to the classical diffusion equation. The traditional diffusion equation fails to effectively explain these behaviors because of nonlocal dynamics and memory effects. The TFDE plays a crucial role in capturing these complicated processes by providing a more precise and all-encompassing modeling framework. Types of TFDE were the focus of numerous contributions. The authors of ^[22] proposed a numerical procedure using splines to solve a class of TFDE. The authors of ^[23] used the Petrov-Galerkin method for treating the TFDE. In ^[24], the authors applied a linearized numerical scheme for a class of nonlinear TFDE. In ^[25], the authors utilized an implicit difference scheme for the TFDE. In ^[26], another difference scheme for the generalized TFDE was followed. The authors of ^[27] used an exponential-sum algorithm for variable-order TFDE.

For numerical analysis and scientific computing, spectral methods provide a collection of approaches for solving DEs. Their distinguishing feature is the ability to describe the solution throughout the entire domain of the problem using global basis functions, unlike approaches like finite difference and finite element. These techniques can provide great precision with limited degrees of freedom, particularly useful for problems with smooth solutions; see ^[28,29,30]. They exhibit high convergence for the approximate solution. These advantages make their application a target for many authors seeking solutions for various DEs. We assume that the solution of a given differential equation is a combination of some basis functions, which could be orthogonal or non-orthogonal polynomials. We extensively employ the various versions of spectral methods to solve different types of DEs. Various papers have utilized the Galerkin method to solve some DEs, see ^[31,32,33]. The tau method is more popular than the Galerkin method. The main reason for this is that there are no restrictions in choosing the trial and test functions if the tau method is utilized, unlike in the Galerkin method; see, for example, ^[34,35]. Various contributions extensively use the typical collocation method because it can treat all types of DEs. For example, the authors of ^[36] followed a meshless superconvergent stabilized collocation method to treat some linear and nonlinear elliptic problems. The collocation method can also be applied to treat FDEs ^[37,38,39]. In addition, matrix methods together with the collocation method were used in ^[40,41]. For some other applications regarding the spectral methods and their applications, one can refer to ^[42,43,44].

This article aims to analyze a new Lucas tau algorithm for handling the TFDE method. To tackle this method, we need formulas concerning the Lucas polynomials, such as their product and specific, definite integral formulas.

The structure of the rest of the paper is as follows: Section 2 presents some fundamentals and useful preliminaries. A tau approach to solving the TFDE is presented in Section 3. The error bound is given in Section 4. Four illustrative examples are provided in Section 5. Finally, Section 6 reports a few conclusions.

2.   Essentials and useful relations

This section introduces some basic definitions and fundamental formulas. The definition of Caputo's fractional derivative is given. An account of Lucas polynomials is given. Certain formulas of Lucas polynomials that will be beneficial are also presented.

2.1.   Caputo's fractional derivative

Definition 2.1. Caputo's fractional derivative is defined as ^[3]

We also have

where $\lceil \zeta \rceil$ is the ceiling function.

2.2.   An overview and some formulas of Lucas polynomials

This recursive formula is utilized to generate Lucas polynomials ^[23]:

and they can be expressed as

that can be alternatively written as

where

and

In addition, $L_{i}(\theta)$ can be expressed as

In the following, we give some important formulas of Lucas polynomials that will be useful in deriving our proposed algorithm.

Theorem 2.1. ^[45] Consider $i, q\in \mathbb{Z^{+}}$ with $i\ge q$. In terms of Lucas polynomials, we can write

where

Theorem 2.2. ^[46] Consider $r$ and $i$ to be two non-negative integers. The following linearization formula holds:

Lemma 2.1. ^[47] For $r\in \mathbf{Z}$, $\int_{0}^{1}L_{r}(\theta)\, d\theta$ is given explicitly by the following formula:

where

Remark 2.1. It is worth mentioning here that many generalizations are established for the standard Lucas polynomials. Among these generalizations, the generalized Lucas polynomials are generated by the following recursive formula:

It is clear that $L_{i}(\theta) = \psi^{1, 1}_{i}(\theta)$.

3.   Tau approach for the TFDE

We confine this section to analyzing a numerical algorithm for solving the TFDE in one dimension. We also account for another extended algorithm that treats a two-dimensional model effectively. We apply the spectral tau approach to obtain the desired approximate solutions in one and two dimensions.

3.1.   tau approximate solution for one-dimensional TFDE

Consider the following TFDE ^[22,23]:

governed by the following conditions:

where $\beta$ is a positive constant, $g(\theta), \, h_{1}(t),$ and $h_{2}(t)$ are given continuous functions, and $f(\theta, t)$ is the source term.

Now, define

and thus, we can assume that any function $\xi^{\mathcal{N}}(\theta, t)\in \mathcal{P}^{\mathcal{N}}$ may be expressed as

where $\mathit{\boldsymbol{L}}(\theta) = [L_1(\theta), L_2(\theta), \ldots, L_{\mathcal{N}+1}(\theta)], \quad \mathit{\boldsymbol{L}}(t)^{T} = [L_1(t), L_2(t), \ldots, L_{\mathcal{N}+1}(t)]^{T},$ and $\mathit{\boldsymbol{C}} = (c_{ij})_{1\leq i, j\leq \mathcal{N}+1}$ is the matrix of unknowns of dimension $(\mathcal{N}+1)^2$.

The residual $\mathbf{Res}(\theta, t)$ of Eq (3.1) has the following form:

As a result of applying tau method, we get

Now, consider the following matrices:

where

Therefore, Eq. (3.6) can be rewritten as

or in the following matrix form:

In addition, the conditions in (3.2) and (3.3) lead to the following equations:

Now, a suitable approach may be used to solve the resulting algebraic system of equations of order $(\mathcal{N}+1)^{2}$, which includes Eqs (3.13) and (3.14).

Theorem 3.1. The elements $a_{i, r}, \, {h}_{i, r}$, and $k_{j, s}$ that appear in the system (3.12) can be expressed in the following closed forms:

where $B_{r, j}$ and $M_{r}$ are given by (2.7) and (2.12), respectively.

Proof. To find the elements $a_{i, r}$, we make use of (2.11) to write

The last formula together with the integral in (2.12) leads to (3.15).

Now, to find the elements $h_{i, r}$, we make use of formula (2.10), after putting $q = 2,$ to get

which, after applying the product formula (2.11), may be transformed into

This shows formula (3.16). To obtain $k_{j, s}$, we use formula (2.1) to get

which immediately gives the following result:

This shows formula (3.17), and thus the proof of Theorem 3.1 is now complete. □

Remark 3.1. In Algorithm 1, we outline the methodology for addressing the TFDE with our spectral tau approach.

3.2.   An extension of our algorithm

Our algorithm can be extended to solve 2D models. Similar steps can be followed to obtain the proposed numerical solution. As an extension to our one-dimensional model, we can solve the following two-dimensional TFDE:

governed by the following conditions:

and

where $\beta_{1} > 0, \beta_{2} > 0$, and $\mathcal{G}(\theta, y, t)$ is the source term.

In this case, we may assume the approximate solution of the form

We follow similar steps to those followed in Section 3.1 to get a linear system of algebraic equations of dimension $(\mathcal{N}+1)^3$ in the unknown expansion coefficients $c_{ijk}$, which can be solved using the Gaussian elimination procedure.

4.   Error bound

This section discusses the error analysis of the proposed double Lucas expansion.

Theorem 4.1. ^[48] Let $u(\theta)$ be a function defined on $[0, 1]$, such that $|u^{(i)}(0)| \leq f^i$ for $i \geq 0$, where $f > 0$, and let it have the expansion $u(\theta) = \sum_{i = 0}^{\infty} \hat{u}_{i} \, \psi^{a, b}_{i}(\theta)$, where $\psi^{a, b}_{i}(\theta)$ denotes the generalized Lucas polynomials generated by (2.14). Then, we have

● $|\hat{u}_{i}| \leq \frac{|a|^{-i}\, f^{i}\, \cosh(2\, |a|^{-1}\, |b|^{\frac{1}{2}}\, f)}{i!}$.

● The series converges absolutely.

Corollary 4.1. Let $u(\theta)$ be a function defined on $[0, 1]$, with $|u^{(i)}(0)|\leq f^i$, $i\geq0$, and $f > 0$, and let it have the expansion $u(\theta) = \sum_{i = 0}^{\infty}\hat{u}_{i}\, L_i(\theta).$ Then, we have

● $|\hat{u}_{i}| \leq \frac{f^{i}\, \cosh(2\, f)}{i!}$.

● The series converges absolutely.

Proof. It is a special case of Theorem 4.1 only by setting $a = b = 1$. □

Theorem 4.2. Consider a function $\xi(\theta, t) = \xi_{1}(\theta)\, \xi_{2}(t)$ defined on $[0, 1] \times [0, 1]$, where $|\xi_{k}^{(i)}(0)| \leq \phi_{k}^i$ for $k = 1, 2$, $i \geq 0$, and $\phi_{k} > 0$. Assume that $\xi(\theta, t)$ has the expansion $\xi(\theta, t) = \sum\limits_{i = 0}^{\infty} \sum\limits_{j = 0}^{\infty}c_{ij}\, L_i(\theta)\, L_j(t)$. Then, we have

● $|c_{ij}| \leq \frac{\phi_{1}^{i}\, \phi_{2}^{j}\, \cosh(2\, \phi_{1})\, \cosh(2\, \phi_{2})}{i!\, j!}$.

● $\sum\limits_{i = 0}^{\infty} \sum\limits_{j = 0}^{\infty}c_{ij}\, L_i(\theta)\, L_j(t)$ converges absolutely.

Proof. Based on Lemma 1 of ^[48], and setting $a = b = 1$, the following relation can be obtained:

According to the assumption $\xi(\theta, t) = \xi_{1}(\theta)\, \xi_{2}(t)$, we have

Now, the expansion coefficients $c_{ij}$ can be written as

Taking into consideration that $|\xi_{k}^{(i)}(0)|\leq \phi_{k}^i$, it is possible to rewrite the final equation in the following form:

Now, based on Corollary 4.1 and after using similar steps as in ^[49], we get the following estimation:

This finalizes the proof. □

Theorem 4.3. Assuming that $\xi(\theta, t)$ meets the conditions of Theorem 4.2, we have

where $\lambda = \cosh \left(2 \phi_1\right)\, \cosh \left(2 \phi_2\right)\, e^{\sqrt{3}\, (\phi_1+\phi_2)}$.

Proof. From the definitions of $\xi(\theta, t)$ and $\xi^{\mathcal{N}}(\theta, t),$ we can write

Theorem 4.2 along with the inequality $|L_{i}(\theta)|\leq 2\, (3)^{\frac{i}{2}}$ yields ^[49]

Using the following inequalities:

we get the following estimation:

where $\lambda = 8\, \cosh \left(2 \phi_1\right)\, \cosh \left(2 \phi_2\right)\, e^{\sqrt{3}\, (\phi_1+\phi_2)}$, $\Gamma(.)$ and $\Gamma(., .)$ denote, respectively, gamma and upper incomplete gamma functions ^[50]. This completes the proof. □

5.   Some numerical tests

To show the usefulness and effectiveness of our suggested numerical algorithms described in Sections 3.1 and 3.2, we offer four test problems here.

Example 5.1. Consider the following equation:

controlled by

with the exact solution $\xi(\theta, t) = t^2\, \sin(\pi\, \theta)$.

Table 1 displays the absolute errors (AEs) at $\zeta = 0.5$ and $\mathcal{N} = 7$. Also, Table 2 presents the AE at $\zeta = 0.7$ and $\mathcal{N} = 6$. Figure 1 shows the AEs (left) and the numerical solution (right) at $\zeta = 0.9$ and $\mathcal{N} = 6$. Finally, Table 3 presents the maximum absolute error (MAE) at $\mathcal{N} = 6$ and various values of $\zeta$ when $0 < t < 1$.

Remark 5.1. We can demonstrate that the theoretical results of the error bound given in Section 4 agree with the numerical results presented in Section 5. As an example, if we set $\phi_1 = \phi_2 = 0.5$ in Eq. (4.5) of Theorem 4.3, then we can see from the results of Tables 1–3 that the values of the MAEs do not exceed those from the theoretical bound given in Table 4.

Example 5.2. ^[22] Consider the following equation:

controlled by

with the exact solution $\xi(\theta, t) = \sin (\pi\, t)\, \sin (\pi\, \theta)$.

In Table 5, we compare the MAE of our method and the method in ^[22] at various values of $\zeta$. Also, Figure 2 shows the AEs (left) and the numerical solution (right) at $\zeta = 0.2$ and $\mathcal{N} = 6$. Finally, Table 6 presents the AEs (left) and the numerical solution (right) at $\zeta = 0.6$ and $\mathcal{N} = 6$.

Example 5.3. ^[51] Consider the following equation:

controlled by

with the exact solution $\xi(\theta, t) = \left(t^2+t+1\right) e^\theta$.

In Table 7, we compare the MAE between our method and the method in ^[51] at different values of $\zeta$. Also, Figure 3 shows the AEs (left) and the numerical solution (right) at $\zeta = 0.4$ and $\mathcal{N} = 6$. Table 8 presents the MAEs at $\mathcal{N} = 6$ and different values of $\zeta$ when $0 < t < 1$. Finally, Figure 4 shows the AEs when $\theta = t$ at $\zeta = 0.25,$ and different values of $\mathcal{N}$.

Example 5.4. Consider the following equation:

governed by

and

where

and the exact solution of this problem is $\xi(\theta, y, t) = t^2\, e^{\theta+y}.$

Figure 5 shows the AEs when $\zeta = 0.3$ and $\mathcal{N} = 4$ at different values of $y$. Table 9 presents the AEs at $\mathcal{N} = 4$ and different values of $t$ when $\zeta = 0.9$.

6.   Conclusions

In this article, we offered two numerical algorithms for treating the TFDE in one and two dimensions. The tau method was used to propose the desired numerical solutions. The suggested numerical solution was an expansion of the double basis of the Lucas polynomials. Some specific formulas of the Lucas polynomials were the backbone for transforming the TFDE with their conditions into a solvable system of equations. We believe that the method may be useful in treating other DEs. The accuracy of the double expansion of Lucas polynomials was tested from a theoretical point of view. In addition, some examples were presented to test the algorithm from a numerical point of view. As an expected future work, we aim to employ this paper's developed theoretical results together with suitable spectral methods to treat some other problems. All codes were written and debugged by Mathematica 11 on HP Z420 Workstation, Processor: Intel(R) Xeon(R) CPU E5-1620 v2 - 3.70GHz, 16 GB Ram DDR3, and 512 GB storage.

Author contributions

WMA contributed to Conceptualization, Methodology, Validation, Formal analysis, Investigation, Project administration, Supervision, Writing – Original draft, and Writing - review & editing. AFAS contributed to Methodology and validation. MHA contributed to Validation and Supervision. AGA contributed to Conceptualization, Methodology, Validation, Formal analysis, Investigation, Software, Writing – Original draft, and Writing - review & editing. All authors have read and agreed to the published version of the manuscript.

Acknowledgments

The authors would like to express their sincere gratitude to the reviewers for their valuable comments and suggestions, which have improved the paper in its present form.

Conflict of interest

The authors declare that they have no competing interests.