Spectral collocation approach with shifted Chebyshev sixth-kind series approximation for generalized space fractional partial di ﬀ erential equations

: In this paper, we propose a numerical scheme to solve generalized space fractional partial di ﬀ erential equations (GFPDEs). Besides, the proposed GFPDEs represent a great generalization of a signiﬁcant type of FPDEs and their applications, which contain many previous reports as a special case. Moreover, the proposed scheme uses shifted Chebyshev sixth-kind (SCSK) polynomials with spectral collocation approach. The fractional di ﬀ erential derivatives are expressed in terms of the Caputo’s deﬁnition. Furthermore, the Chebyshev collocation method together with the ﬁnite di ﬀ erence method is used to reduce these types of di ﬀ erential equations to a system of algebraic equations which can be solved numerically. In addition, the classical fourth-order Runge-Kotta method is also used to treat the di ﬀ erential system with the collocation method which obtains a great accuracy. Numerical approximations performed by the proposed method are presented and compared with the results obtained by other numerical methods. The introduced numerical experiments are fractional-order mathematical physics models, as advection-dispersion equation (FADE) and di ﬀ usion equation (FDE). The results reveal that our method is a simple and e ﬀ ective numerical method.


Introduction
Many phenomena such as biology, physics, and fluid mechanics can be modeled by certain fractional order partial differential equations (FPDEs). For example not exclusively, Magin in [1] presented fractional calculus models of complex dynamics in biological tissues, where, Kumar and Baleanu presented applications in physics related to fractional calculus [2], and in [3] their authors presented a real world applications of fractional calculus in science and engineering. So that, fractional calculus becomes a major branch of mathematical and numerical analysis. The importance of the numerical solution of FPDEs becomes a major, because of the difficulty of obtaining their analytical solutions. Vanani and Aminataei propose numerical solutions of FPDEs see [4]. In addition, the authors of [5] presented a numerical approach for FPDEs using Legendre functions, also [6] presents a high-precision numerical approach to solving space fractional Gray-Scott model. On the other hand, finite difference and finite element methods have received a great attention for treating the FPDEs see [7][8][9][10][11][12][13][14].
Spectral methods have been developed through the past few decades by a huge number of researchers, the most popular and widely used methods are: tau [15,16], collocation [17][18][19], Galerkin [20], kernel method and Legendre polynomial [21], Fourier spectral method [22] and many others. The principal feature of these methods lies in their ability to reach acceptably accurate results with substantially fewer degrees of freedom. In recent years, Chebyshev polynomials have become increasingly important again in numerical analysis, when a new two classes of polynomials appear, namely fifth and sixth kinds. In the Ph.D. thesis of Masjed-Jamei [23], 2006 he introduces a generalized polynomials using an extended Sturm-Liouville problem. These generalized polynomials generate Chebyshev polynomials of the first, second, third, and fourth kinds, in addition to the two new classes fifth and sixth kinds at special values of the given parameters, some additional details are also provided in [24] for fifth kind and in [25] for sixth kind. Therefore, many researchers recently began to apply numerical methods to solve mathematical models using fifth and sixth kinds. Abd-Elhameed et al, used the fifth-kind with spectral solution for convection-diffusion (CD) equation [26], while in [27] the authors used fifth Chebyshev polynomials to fractional partial integro-differential equations (FPIDEs). Additionally, very recently Sadri and Aminikhah presented an algorithm based on fifth-kind to treat multi-term variable-order time-fractional diffusion-wave (FDW) equation [28], however, the authors in [29] used the Galerkin approach with the fifth-kind for a kind of PDEs.
The objective of this research paper is to present a spectral scheme according to the collocation approach for the generalized space-fractional partial differential equations (GFPDEs) that we have introduced. The proposed GFPDEs are chosen to be linear and the fractional derivatives are expressed in terms of Caputo's definition. The method of solution aims to apply shifted sixth-kind Chebyshev polynomials using the collocation method to discretize the proposed equation, and then generate a linear system of ordinary differential equations (SODEs), which reduces the proposed problem.
Additionally, to treat the generated SODEs, the classical fourth-order Runge-Kotta method (RK4) and the finite difference method (FDM) as well, are used. The proposed equation is presented as: defined on a finite domains 0 < x H; 0 < t T and the parameters γ k refer to the fractional orders of a special derivative with k < γ k < (k + 1) n. The function f (x, t) is the source term, the functions Q k (x) are well defined and known, and P is a real constant. We also assume the initial condition (IC) as: 2) and the boundary conditions (BCs): In addition, the proposed GFPDEs (1.1) represent a great generalization of significant types of many applications. As special cases: at γ 1 0, γ k = 0, equation (1.1) reduces to a space-fractional order diffusion equation, and when γ 0 , γ 1 0, γ k = 0, then (1.1) becomes a space-fractional order advection-dispersion equation, which they will be studied in the application section, and more.
The rest of the paper is organized as: section two contains some notations of Chebyshev sixth kind and its properties; also, some properties of the Caputo's derivative are briefly listed. While, in section three the description of the solution process is presented. In section four, the numerical scheme based on the collocation method is obtained. Finally, section five contains the numerical applications and results, also comparisons with the previous works listed literately.

Main notations
In this section, some definitions and properties for the sixth kind Chebyshev polynomials [30] and fractional derivative [31] are listed.

Sixth-kind Chebyshev polynomials
The basis polynomials used in this work are the Chebyshev polynomials of the sixth-kind Y n (x) and they are defined as: an orthonormonal polynomials in x of degree n defined on the closed interval [−1, 1]. The polynomials Y k (x), k = 0, 1, ..., n form an orthogonal system and the orthogonality relation is: if i=j, and i, j even, The shifted Chebyshev polynomials of the sixth-kind Y * n (x) are orthogonal on the closed interval [0, 1], and are generated by using the following recurrence relation From (2.1), it is not difficult to note that Y * n (x), n ≥ 0, form an orthonormal system on [0, 1], and they have an orthogonality relation as: if i=j, and i, j even, (2.4) Proposition 1. The shifted polynomials Y * n (x) are defined through the shifted second kind U * n (x) by the following formula: if n and k even, , if n and k odd, 0, other , .
Proof. The complete proof is given in [30].
According to Proposition 1, the following corollary is easy to prove.
Corollary 2.1. Shifted Chebyshev polynomials of the sixth-kind Y * n (x) are explicitly expressed in terms of U * n (x) in the following form: Corollary 2.2. shifted Chebyshev polynomials of the sixth-kind Y * n (x) are explicitly expressed in terms of x n , or the analytic form in the following form: where . is the floor function.
According to relations (2.3), (2.5), (2.6), (2.7) and (2.8), the first four terms of Y * n (x) are: where Y * 0 (x) and Y * 1 (x) are defined before in the recurrence relation (2.3). Proposition 1 gives the connection formulae of the sixth-kind Chebyshev polynomials and the second Chebyshev polynomials. In Ref [25], the authors drive the connection formulae, that can be express the sixth-kind Chebyshev polynomials in terms of ultraspherical polynomials functions. Moreover, the ultraspherical polynomials are generalized polynomials which may give Chebyshev polynomials of the first, second-kinds and Legendre polynomials as special ones of them. The connection formulas between the sixth-kind polynomials and these polynomials can be deduced as special cases of ultraspherical polynomials functions and by extension, the sixth-kind inherits from them its ability and convergence.
Lemma 2.1. Shifted Chebyshev polynomials of the sixth-kind Y * n (x) are bounded according to the following form: The full proof is in [30], and it directly given from the connection relation (2.5).

The Caputo's fractional derivative
The Caputo's fractional derivative operator D γ t (instead of C 0 D γ t for short) of order γ is characterized in the following form: where λ i and γ are constants.
• Property 2 The Caputo fractional differentiation of a constant is zero. Such that: • Property 3 where γ denotes to the smallest integer greater than or equal to γ. For more properties about the Caputo's derivative see [32,33] Remark 1. In this work, we write the fractional Caputo's operator symbol D γ instead of C 0 D γ x for short. Remark 2. In this work, the shifted Chebyshev polynomials of the sixth-kind Y * n (x) are considered, which are defined on [0, 1] then upper limit of the space argument H in (1.1) becomes 1.

Process the solution
In the spectral method, in contrast, the function Ψ(x) may be expanded by shifted Chebyshev polynomials of the sixth-kind series, which Ψ(x) is a square-integrable in [0, 1], [29,30,34]: Lemma 3.1. The series in (3.1) uniformly converges to Ψ(x), where the following relation holds: and L is some positive constant provided from: Lemma 3.2. The global error e N (x) for the function Ψ(x) defined in (3.1), such that: e N (x) = ∞ n=N+1 a n Y * n (x), is bounded and the following relation is valid: The proof of Lemma 3.2 is found in [30], and it refers to that the error almost tends to zero in the case of a large N. By truncate series (3.1) to N < ∞, then the approximation of Ψ(x) is given by a finite sum of (N + 1)−terms and expressed in the following form: The coefficients a n in relation (3.5) are given by the following relation: and ε n is given from: According to the definition of Caputo's fractional derivative (2.10), property 1 and the analytic form (2.8) the following theorem is introduced.
The fractional derivative of order γ for the polynomials Y * n (x) is given by: where, k,n is defined in Corollary 2.2.
Proof. According to (2.10) (the Caputo's operator) and the relation given in corollary 2.2 it is easy to obtain the result, for more details see [24,30,35].
is an approximated function of Ψ(x) in terms of shifted Chebyshev polynomials of the fifth kind as (3.5), then the Caputo fractional derivative of order γ when operating Ψ N (x) is given by: where, k,n is defined in Corollary 2.2.
Proof. According to Theorem .1 and relation (3.5) one writes: then the result (3.10) is easily obtained.

Numerical scheme
Consider the generalized space fractional partial differential equations of the type given in equation (1.1) with their given conditions. In order to use the Chebyshev sixth-kind collocation method, let us approximate u(x, t) as follows [36,37]: with the help of Theorem 1, then Now, we turn to collocate equation (4.3) at (N + 1) points, the collocation points are defined in the following form: By substituting the collocation points (4.4) in (4.3), we get Also, two additional equations may generate from the boundary conditions using relation (4.1) in (1.3) as: The collocated Eq (4.5), together with the generated equations of the boundary conditions (4.6), give us an ordinary system of differential equations with φ k (t) as the unknowns, which can be solved by a suitable technique. Using the initial conditions (1.2) and by the help of relation (4.1) and the orthogonality (2.4), we can generate initial conditions for the proposed system of differential equations, the IC may take the form: and by expanding h(x) in terms of Y * k (x) and comparing the coefficients, then we get the constants φ k in the initial case at t = 0, (φ k (0)). The produced system of ordinary differential equations according to (4.5) is linear and generally has the following matrix form: andQ is a square constant matrix representing the coefficients of the unknowns φ k (t), which is featured by the first column is null. Additionally, (4.8) may be written as: therefore, the system (4.9) is ready to be solved with a suitable solver technique, under the subjected initial conditions (4.7).

Numerical applications
In this section, several numerical applications (physical models) have been given to illustrate the accuracy and effectiveness of the method.

Example 1:
Consider the following space fractional order PDE: The IC is: and the BCs: By the same process, as Eqs (4.2)-(4.9), we have and by expanding h(x) = x 2 (the IC) in terms of Y * k (x) according to (2.4) and comparing the coefficients, then we get the initial conditions of the differential system as: Two additional equations may generate from the boundary conditions (5.3) using relation (4.1) in (5.3), then System (4.9) with matrices (5.5) and the initial conditions (5.6) is a system of differential equations, (equations (5.7) may be replaced with the last two equations in (4.9)) the Runge-Cotta method of the fourth-order (RK4) is used here with h step size equal to 0.01 with 100 iterations means that 0 ≤ t ≤ 1, (the regular algorithm for RK4 is coded by the authors using Mathematica.10. package) the numerical results obtained as:  According to (5.4), one obtains the approximate solution u 3 (x, 1) (at t = 1 ) using the last row in (5.8) as: As references [36][37][38], their numerical results were obtained using a finite difference method (FDM) for the differential system. We turn to solve the system (4.9) with matrices (5.5) using FDM. Then, Therefore, the system in Eq (4.9) with matrices (5.5), is discretized in the time and has the following form: where Hence, a sample of the numerical results for FDM is obtained as:  In Table 1, the comparison of the absolute errors for the present method with both RK4 and FDM at N = 3, ∆t = 0.01 where γ 0 = 0.9, also, shows the numerical values of the approximate solution using the proposed method (using both RK4 and FDM) with the exact solution. Also, Table 2 shows the L 2 error norm [39] at N = 3 at different values of T . In Figure 1 and Figure 2, the comparison of the exact and the approximate solutions with both RK4 and FD methods for example.1 with N = 3 and T = 1, 2. Additionally, the CPU time used for getting the approximate solution by using the present schemes at N = 3 is listed in Table 3. Table 3 shows that the time used for calculating the approximate solution by RK4 is less than that given by FDM. We note that: all of the examples are implemented by the help of Mathematica 7.1 package and handled on a usual machine (Intel(R)-core(TM)-i3, CPU-3.43 GHz).

Example 2:
Consider the following generalized space fractional order diffusion equation of the following type: then equation (5.13) has the exact solution of the form u(x, t) = x 2 (1 − x)e −t at γ 1 = 1.8, which is mentioned in [36][37][38]. The IC is: (5.14) and the BCs: u(0, t) = u(1, t) = 0, 0 < t ≤ T. (5.15) At N = 3, according to (4.1), (using same process (4.2)-(4.9)), we have and C not changed for N = 3 as example 1. In addition, by expanding h(x) = x 2 (1− x) in terms of C k (x) according to (2.4) and comparing the coefficients, then we get the initial conditions of the differential system as: The generated equations from the homogeneous boundary conditions (5.15) using relation (4.1) are: System (4.9) with matrices (5.16) and the initial conditions (5.17) is a system of differential equations, by replacing Eqs (5.18) with the last two equations in (4.9) the RK4 method used as example.1 with 0 ≤ t ≤ 2. The RK4 method's numeric results at t = 1, t = 2, N = 3 are obtained as: (5.20) In Table 4 the comparison of the absolute errors for the present two schemes at N = 3, T = 2 with the methods mentioned in [36][37][38]. Also, the numerical absolute errors are represented in Table 4 for the collocation method with Chebyshev first [38] second [36] and third [37] kinds. These values show that the sixth kind gives a more accurate approximate solution using the proposed method with RK4, but less accuracy is given when using regular FDM with the present method. Additionally, Table 5 shows the L 2 error norm at N = 3 at two values of T . The CPU time used for getting the approximate solution by using the present schemes at N = 3 is listed in Table 6. In Figure 3 and Figure 4 the comparison of the exact and the approximate solutions with both RK4 and FD methods for example.2 with N = 3 and T = 1, 2.
The generated equations from the homogenous boundary conditions (5.28) are the same as given in example 2 and example 3. The system (4.9) with matrices (5.29) has zero ICs, by replacing the generated equations from the homogenous boundary conditions with the last two equations in (4.9), the RK4 used as examples 2, 3 with 0 ≤ t ≤ 2. As ref [40] the numerical results were obtained using the non-standard FDM for the differential system with the aid of Vieta-Lucas polynomials. Therefore, as example 3 we turn to solve the system (4.9) with matrices (5.24) using FDM. Then we use same elements as examples 2, 3, for system (5.10), (5.11) but using matrices (5.29). The numerical comparisons will hold only with [40] because the results in [42,43] (collocation method with fourth and second Chebyshev kinds) are less than 10 −5 , it is much less accurate than indicated in our results. In Table 10 the comparison of the absolute errors for the present two schemes (PM with RK4 and FDM) at N = 3, where γ 1 = 1.5, γ 0 = 1, T = 0.5 with [40], while same comparison given in Table 11 but T = 0.5. Also, it shows the numerical values of the proposed method gives a highly accurate approximate solution with RK4, and [40] which uses a modified technique gives accuracy near PM with FDM. Table 12 gives the L 2 error norm along the interval [0, 1] at N = 3 with three values of T . The CPU time used for getting the approximate solution by using the present schemes at N = 3 is listed in Table 13. In Figure 7, Figure 8 and Figure 9 the comparison of the exact and the approximate solutions with both RK4 and FD methods for example 1 with N = 3 and T = 0.3, 0.5, 0.9. In the end, we conclude that the Chebyshev sixth-kind series approximation gives a great accuracy when using high appropriate accurate methods, and the Runge-Kota method remains one of the best methods in dealing with linear systems, as was shown in the last two examples. Table 10. Comparing absolute errors for present technique at N = 3, T = 0.5 with different methods for example 4.  Table 11. Comparing absolute errors for present technique at N = 3, T = 0.9 with different methods for example 4.

Conclusions
A numerical study for a generalized form of linear space-fractional partial differential equations (GFPDEs) is introduced using the Chebyshev sixth-kind series. The suggested general form represents many fractional-order mathematical physics models, as advection-dispersion equation (FADE) and diffusion equation (FDE). Additionally, the proposed scheme uses the shifted Chebyshev polynomials of the sixth-kind, where the fractional derivatives are expressed in terms of Caputo's definition. Therefore, the collocation method is used to reduce the GFPDE to a system of ordinary differential equations which can be solved numerically. Moreover, the classical fourth-order Runge-Kotta method is used to treat the differential system as well as the finite difference method which obtains a great accuracy. We have presented many numerical examples, where represent mathematical physical models, that greatly illustrate the accuracy of the presented study to the proposed GFPDE, and also show how that the sixth-kind polynomials are very competitive than others.

Data availability
The data used to support the findings of this study are available from the corresponding author upon request.