ABSTRACT

Recently many research works have been conducted and published regarding fractional order differential equations. There are several approaches available for numerical approximations of the solution of fractional order diffusion equations. Spectral collocation method based on Lagrange’s basis polynomials to approximate numerical solutions of one-dimensional (1D) space fractional diffusion equations are introduced in this research paper. The proposed form of approximate solution satisfies non-zero Dirichlet’s boundary conditions on both boundaries. Collocation scheme produce a system of first order Ordinary Differential Equations (ODE) from the fractional diffusion equation. We applied this method with four different sets of collocation points to compare their performance.

Keywords:

Fractional Diffusion Equation, Spectral Method, Collocation Method, Lagrange’s Basis Polynomial

1. Introduction

In recent times, a huge number of research articles have been published by researchers around the world regarding development of various methods for fractional differential equations. The sudden growth of attention around fractional differential equation is because of its use to describe complex physical phenomena like super diffusive process in diverse fields has also been explored. Such applications are cited in [1] . One-dimensional fractional partial differential equations can be put into two immediate categories: space fractional and time fractional differential equations. The order of the spatial derivative has a fractional value instead of integer value in space fractional partial differential equations. Among fractional partial differential equations diffusion and wave equations are most popular. 1D space fractional diffusion equation, contrary to classical diffusion equation spatial derivative will have fractional order rather than integral order; will be considered here in our research. Fractional derivative which is associated here is Caputo fractional derivative.

Several popular approaches have been introduced by researchers for numerical approximations of both space and time fractional diffusion and wave equations. Among them spectral collocation method is the most popular and effective. In spectral collocation method there are wide ranges of polynomials and set of collocation points which are available to choose from. Khader [2] used Chebyshev spectral collocation method to reduce space fractional diffusion equation into a system of ODE for time variable and then solved the system by finite difference method. He chose Chebyshev polynomials for space variable along with roots of shifted Chebyshev polynomials as collocation points. Azizi and Loghmani [3] also used Chebyshev spectral collocation method but they reduced the fractional diffusion equation into a set of algebraic equations using Chebyshev polynomials and Gauss-Lobatto nodes in both space and time domain. Xie et al. [4] also used Chebyshev polynomials to express the solution both in space and time but used Tau method to transform the fractional convection diffusion equation into a system of linear algebraic equations. Bhrawy [5] used shifted Legendre polynomials with Gauss-Lobatto nodes in space for a 2D space fractional diffusion equation and hence reduced it to a system of ODE which is then solved by fourth-order implicit Runge-Kutta method. Lin and Xu [6] introduced a method where they applied Legendre spectral scheme in space and finite difference in time to approximate the solution of time fractional diffusion equation. Bahsi and Yalcinbas [7] reduced the fractional order diffusion equation into a system of linear algebraic equations by expanding trial solution in terms of Fibonacci polynomials in both space and time and then using collocation technique with evenly spaced collocation points. Pirim and Ayaz [8] introduced Hermite collocation method with evenly spaced nodes for numerical approximations of fractional order system of differential equations. Huang and Zheng [9] presented a spectral method to calculate fractional derivative described in Riemann-Liouville sense using Jacobi orthogonal polynomials. They also discussed spectral collocation method based on Lagrange’s basis polynomials and Gauss-Lobatto nodes. With these methods they solve space fractional diffusion equation that have Dirichlet’s boundary conditions; zero in one boundary and non-zero in another. Spectral expansion with Lagrange interpolation polynomial is used for both space and time by Huang [10] for numerical approximations of time fractional differential equations. He used Jacobi-Gauss nodes in time domain and Gauss-Lobatto nodes in space as collocation points.

In this research paper, we will present a spectral collocation method where approximate solution will be expressed in terms of Lagrange’s basis polynomials in space and then a system of first order ODE for time variable is generated by collocation scheme from space fractional diffusion equation. In our proposed technique, approximate solution satisfies non-zero Dirichlet’s boundary conditions on both boundaries. We considered four different sets of collocation points to demonstrate their performance into proposed spectral collocation scheme. The four sets of collocation points are generated from Gauss-Lobatto nodes, roots of Chebyshev polynomials of first kind, roots of Legendre polynomials and equally spaced nodes over the space domain.

Remaining of this research paper is presented as follows: Preliminaries of Caputo fractional derivative and brief introduction of different polynomials are given in Section 2. Then detailed Spectral collocation scheme based on Lagrange’s basis polynomials along with error calculation are provided in Section 3. After that in Section 4, numerical solutions of two examples of space fractional diffusion equation are generated using four different sets of collocation points and absolute local error curves are given. Finally, Section 5 deals with the conclusion.

2. Preliminaries

Caputo Fractional Derivative: Caputo fractional derivative operator of order α is denoted by $D^{\alpha }$ and defined by:

$D^{\alpha }f\left(x\right)=\frac{1}{\text{Γ}\left(m-\alpha \right)}\underset{0}{\overset{x}{{\displaystyle \int }}}\frac{f^{\left(m\right)}\left(t\right)}{{\left(x-t\right)}^{\alpha -m+1}}\text{d}t,\alpha >0$ (1)

with $m-1<\alpha <m,m\in \mathbb{N},x>0,x\in ℝ$ .

Then for a constant $c$ , we have $D^{\alpha }c=0$ and

$D^{\alpha }{x}^n=\{\begin{array}{ll}0,\hfill & \text{for}n\in {\mathbb{N}}_0\text{and}n<\lceil \alpha \rceil \hfill \\ \frac{\text{Γ}\left(n+1\right)}{\text{Γ}\left(n+1-\alpha \right)}{x}^{n-\alpha },\hfill & \text{for}n\in {\mathbb{N}}_0\text{and}n\ge \lceil \alpha \rceil \hfill \end{array}$ (2)

where ${\mathbb{N}}_0=\left\{0,1,2,3,\cdots \right\}$ and $\mathbb{N}=\left\{1,2,3,\cdots \right\}$ .

Like classical integer order derivative, Caputo fractional order derivative is also a linear operator. Also it is evident from Equation (2) that for $\alpha \in \mathbb{N}$ , Caputo fractional order derivative coincides with the classical integer order derivative.

Lagrange Basis Polynomials: For $\left(p+1\right)$ points $x_1,x_2,\cdots ,x_p,x_{p+1}$ Lagrange basis polynomials $L_n\left(x\right);n=1,2,\cdots ,p+1$ is defined as follows:

$pp\left(x\right)={\displaystyle \prod }_{k=1}^{p+1}\left(x-x_k\right)$ (3)

$L_n\left(x\right)=\frac{pp\left(x\right)}{p{p}^{\prime }\left(x\right)\left(x-x_n\right)}={\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)x^r;n=1,2,\cdots ,p+1$ (4)

with the property $L_n\left(x_m\right)={\delta }_{mn}$ , where ${\delta }_{mn}$ is the kronecker delta function. Here $w\left(n,r\right)$ is the coefficient of $x^r$ in $L_n\left(x\right)$ and $p{p}^{\prime }\left(x\right)$ is the derivative of $pp\left(x\right)$ .

Legendre Polynomials: Legendre polynomials $P_n\left(w\right)$ are solutions of the Legendre differential equations and are orthogonal over the domain $\left[-1,1\right]$ . Explicit formula for $P_n\left(w\right)$ is

$P_n\left(w\right)={\displaystyle \sum }_{k=0}^n\left(\begin{array}{c}n\\ k\end{array}\right)\left(\begin{array}{c}-n-1\\ k\end{array}\right){\left(\frac{1-w}{2}\right)}^k$ (5)

Chebyshev Polynomials: Chebyshev polynomials $T_n\left(y\right)$ are solutions of the Chebyshev differential equations and are orthogonal over the domain $\left[-1,1\right]$ . Explicit formula for $T_n\left(y\right)$ is

$T_n\left(y\right)=y^n{\displaystyle \sum }_{k=0}^{\lfloor \frac{n}{2}\rfloor }\left(\begin{array}{c}n\\ 2k\end{array}\right){\left(1-y^{-2}\right)}^k$ (6)

Roots of Legendre and Chebyshev polynomials are within the interval $\left(-1,1\right)$ , later in Section 4, we shifted these roots to the required interval according to the problem.

3. Spectral Collocation with Lagrange’s Basis Polynomial

Here we present spectral collocation method which is based on Lagrange’s basis polynomials for numerical approximations of the solution of following 1D space fractional diffusion equation:

$\frac{\partial u\left(x,t\right)}{\partial t}=d\left(x,t\right)\frac{{\partial }^{\alpha }u\left(x,t\right)}{\partial x^{\alpha }}+s\left(x,t\right);a<x<b,0\le t\le T$ (7)

$u\left(x,0\right)=u^0\left(x\right),a<x<b$ (8)

$u\left(a,t\right)=v_1\left(t\right)\text{and}u\left(b,t\right)=v_{p+1}\left(t\right)$ (9)

Here parameter $\alpha$ represents the fractional order of spatial derivative where $1<\alpha \le 2$ and associated fractional derivative described in Caputo sense. Equation (8) is the initial condition; Equation (9) is the boundary condition, $d\left(x,t\right)$ is the diffusion coefficient and $s\left(x,t\right)$ is known as source function. With $\alpha =2$ , Equation (7) is the classical diffusion equation.

To approximate numerical solution of fractional diffusion equation given in Equation (7), we first divide the space domain $\left[a,b\right]$ into $p$ parts that results into following $\left(p+1\right)$ points along with boundary points:

$a=x_1<x_2<\cdots <x_p<x_{p+1}=b$

Later in this section we use these points as collocation points and these points can be chosen from anywhere within the domain in no specific pattern.

Now using Equations (3), (4) and above points we can form Lagrange’s basis polynomials $L_n\left(x\right)$ for $n=1,2,\cdots ,\left(p+1\right)$ .

We approximate the solution $u\left(x,t\right)$ of the fractional diffusion equation as $u_p\left(x,t\right)$ by finite sum of Lagrange’s basis polynomials $L_n\left(x\right)$ :

$u_p\left(x,t\right)={\displaystyle \sum }_{n=1}^{p+1}{v}_n\left(t\right)L_n\left(x\right)$ (10)

$=v_1\left(t\right)L_1\left(x\right)+{\displaystyle \sum }_{n=2}^p{v}_n\left(t\right)L_n\left(x\right)+v_{p+1}\left(t\right)L_{p+1}(x)$

$=v_1\left(t\right){\displaystyle \sum }_{r=0}^{p}w\left(1,r\right)x^r+{\displaystyle \sum }_{n=2}^p{\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)v_n\left(t\right)x^r+v_{p+1}\left(t\right){\displaystyle \sum }_{r=0}^{p}w\left(p+1,r\right)x^r$

$={\displaystyle \sum }_{r=0}^p\left(v_1\left(t\right)w\left(1,r\right)+v_{p+1}\left(t\right)w\left(p+1,r\right)\right)x^r+{\displaystyle \sum }_{n=2}^p{\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)v_n\left(t\right)x^r$ (11)

The unknowns $v_n\left(t\right),n=2,3,\cdots ,p$ in the trial solution needed to be determined. It is clear that this trial solution $u_p\left(x,t\right)$ automatically satisfies conditions on both boundaries: $u\left(a,t\right)=v_1\left(t\right)$ and $u\left(b,t\right)=v_{p+1}\left(t\right)$ . With the form of trial solution given in Equations (10), (11) and from the definition of fractional derivative given in Equation (2) along with its linear property we can write the required derivatives of trial solution as

$\frac{\partial u_p}{\partial t}={\displaystyle \sum }_{n=1}^{p+1}{L}_n\left(x\right)\frac{\text{d}{v}_n\left(t\right)}{\text{d}t}$ (12)

and

$\begin{array}{c}\frac{{\partial }^{\alpha }{u}_p}{\partial x^{\alpha }}={\displaystyle \sum }_{r=0}^p\left(v_1\left(t\right)w\left(1,r\right)+v_{p+1}\left(t\right)w\left(p+1,r\right)\right)D^{\alpha }\left(x^r\right)\\ +{\displaystyle \sum }_{n=2}^p{\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)v_n\left(t\right)D^{\alpha }\left(x^r\right)\end{array}$ (13)

Now using above two derivatives of trial solution from Equations (12) & (13) into Equation (7) we have

$\begin{array}{c}{\displaystyle \sum }_{n=1}^{p+1}{L}_n\left(x\right)\frac{\text{d}{v}_n\left(t\right)}{\text{d}t}=d\left(x,t\right){\displaystyle \sum }_{n=2}^p{\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)v_n\left(t\right)D^{\alpha }\left(x^r\right)+s\left(x,t\right)\\ +d\left(x,t\right){\displaystyle \sum }_{r=0}^p\left(v_1\left(t\right)w\left(1,r\right)+v_{p+1}\left(t\right)w\left(p+1,r\right)\right)D^{\alpha }\left(x^r\right)\end{array}$ (14)

Then using trial solution from Equation (10) into Equation (8) we have

$u_p\left(x,0\right)=u^0(x)$

$\Rightarrow v_1\left(0\right)L_1\left(x\right)+{\displaystyle \sum }_{n=2}^p{v}_n\left(0\right)L_n\left(x\right)+v_{p+1}\left(0\right)L_{p+1}\left(x\right)=u^0(x)$

$\Rightarrow {\displaystyle \sum }_{n=2}^p{v}_n\left(0\right)L_n\left(x\right)=u^0\left(x\right)-v_{p+1}\left(0\right)L_1\left(x\right)-v_{p+1}\left(0\right)L_{p+1}\left(x\right)$ (15)

In the trial solution there are $\left(p-1\right)$ unknowns $v_n\left(t\right);n=2,3,\cdots ,p$ and among $\left(p+1\right)$ points $x_i;i=1,2,3,\cdots ,p,\left(p+1\right)$ trial solution automatically satisfies the boundary conditions at $a=x_1$ and $b=x_{p+1}$ . So, in collocation method to determine the unknowns we will force Equations (14) & (15) to satisfy at each $x_i;i=2,3,\cdots ,p$ . That means from Equation (14) we write

$\begin{array}{l}{\displaystyle \sum }_{n=1}^{p+1}{L}_n\left(x_i\right)\frac{\text{d}{v}_n\left(t\right)}{\text{d}t}=d\left(x_i,t\right){\displaystyle \sum }_{n=2}^p{\displaystyle \sum }_{r=0}^{p}w\left(n,r\right)v_n\left(t\right)D^{\alpha }\left(x_i^r\right)+s\left(x_i,t\right)\\ +d\left(x_i,t\right){\displaystyle \sum }_{r=0}^p\left(v_1\left(t\right)w\left(1,r\right)+v_{p+1}\left(t\right)w\left(p+1,r\right)\right)D^{\alpha }\left(x_i^r\right);\text{for}i=2,3,\cdots ,p\end{array}$ (16)

where $D^{\alpha }\left(x_i^r\right)$ denotes the value of $D^{\alpha }\left(x^r\right)$ at $x_i$ .

Thus we have the following matrix equation

$A\left[\begin{array}{c}{v^{\prime }}_2\\ {v^{\prime }}_3\\ ⋮\\ {v^{\prime }}_p\end{array}\right]=B\left(t\right)\left[\begin{array}{c}{v}_2\\ v_3\\ ⋮\\ v_p\end{array}\right]+\left[\begin{array}{c}s\left(x_2,t\right)\\ s\left(x_3,t\right)\\ ⋮\\ s\left(x_p,t\right)\end{array}\right]+F\left(t\right)\left[\begin{array}{c}{v}_1\left(t\right)w\left(1,0\right)+v_{p+1}\left(t\right)w\left(p+1,0\right)\\ v_1\left(t\right)w\left(1,1\right)+v_{p+1}\left(t\right)w\left(p+1,1\right)\\ ⋮\\ v_1\left(t\right)w\left(1,p\right)+v_{p+1}\left(t\right)w\left(p+1,p\right)\end{array}\right]$

which in short we write as

$A{v}^{\prime }\left(t\right)=B\left(t\right)v\left(t\right)+s\left(t\right)+se\left(t\right)$ (17)

where

$A=\left[\begin{array}{ccc}{L}_2\left(x_2\right)& L_3\left(x_2\right)& \begin{array}{cc}\cdots & L_p\left(x_2\right)\end{array}\\ L_2\left(x_3\right)& L_3\left(x_3\right)& \begin{array}{cc}\cdots & L_p\left(x_3\right)\end{array}\\ \begin{array}{c}⋮\\ L_2\left(x_p\right)\end{array}& \begin{array}{c}⋮\\ L_3\left(x_p\right)\end{array}& \begin{array}{cc}\begin{array}{c}\ddots \\ \cdots \end{array}& \begin{array}{c}⋮\\ L_p\left(x_p\right)\end{array}\end{array}\end{array}\right]=\left[\begin{array}{ccc}1& 0& \begin{array}{cc}\cdots & 0\end{array}\\ 0& 1& \begin{array}{cc}\cdots & 0\end{array}\\ \begin{array}{c}⋮\\ 0\end{array}& \begin{array}{c}⋮\\ 0\end{array}& \begin{array}{cc}\begin{array}{c}\ddots \\ \cdots \end{array}& \begin{array}{c}⋮\\ 1\end{array}\end{array}\end{array}\right]=I$

$F\left(t\right)=\left[\begin{array}{cccc}d\left(x_2,t\right)D^{\alpha }\left(x_2^0\right)& d\left(x_2,t\right)D^{\alpha }\left(x_2^1\right)& \cdots & d\left(x_2,t\right)D^{\alpha }\left(x_2^p\right)\\ d\left(x_3,t\right)D^{\alpha }\left(x_3^0\right)& d\left(x_3,t\right)D^{\alpha }\left(x_2^1\right)& \cdots & d\left(x_3,t\right)D^{\alpha }\left(x_2^p\right)\\ ⋮& ⋮& \ddots & ⋮\\ d\left(x_p,t\right)D^{\alpha }\left(x_p^0\right)& d\left(x_p,t\right)D^{\alpha }\left(x_p^1\right)& \cdots & d\left(x_p,t\right)D^{\alpha }\left(x_p^p\right)\end{array}\right]$

$v\left(t\right)=\left[\begin{array}{c}{v}_2\\ v_3\\ \begin{array}{c}:\\ v_p\end{array}\end{array}\right],v^{\prime }\left(t\right)=\left[\begin{array}{c}{v^{\prime }}_2\\ {v^{\prime }}_3\\ \begin{array}{c}:\\ {v^{\prime }}_p\end{array}\end{array}\right],s\left(t\right)=\left[\begin{array}{c}s\left(x_2,t\right)\\ s\left(x_3,t\right)\\ \begin{array}{c}:\\ s\left(x_p,t\right)\end{array}\end{array}\right]$

and

$se\left(t\right)=F\left(t\right)\cdot \left[\begin{array}{c}{v}_1\left(t\right)w\left(1,0\right)+v_{p+1}\left(t\right)w\left(p+1,0\right)\\ v_1\left(t\right)w\left(1,1\right)+v_{p+1}\left(t\right)w\left(p+1,1\right)\\ \begin{array}{c}⋮\\ v_1\left(t\right)w\left(1,p\right)+v_{p+1}\left(t\right)w\left(p+1,p\right)\end{array}\end{array}\right]$

Since $A=I$ , Equation (17) immediately becomes

$v^{\prime }\left(t\right)=B\left(t\right)v\left(t\right)+s\left(t\right)+se\left(t\right)$ (18)

Thus with the help of collocation method we reduce Equation (7) which was a fractional diffusion equation into Equation (18), a system of ODE.

Now by forcing Equation (15) to satisfy at each $x_i;i=2,3,\cdots ,p$ we can write

${\displaystyle \sum }_{n=2}^p{v}_n\left(0\right)L_n\left(x_i\right)=u^0\left(x_i\right)-v_1\left(0\right)L_1\left(x_i\right)-v_{p+1}\left(0\right)L_{p+1}\left(x_i\right);i=2,3,\cdots ,p$

$\begin{array}{l}\Rightarrow v_2\left(0\right)L_2\left(x_i\right)+v_3\left(0\right)L_3\left(x_i\right)+\cdots +v_p\left(0\right)L_p\left(x_i\right)\\ =u^0\left(x_i\right)-v_1\left(0\right)L_1\left(x_i\right)-v_{p+1}\left(0\right)L_{p+1}\left(x_i\right);i=2,3,\cdots ,p\end{array}$

$\Rightarrow \left[\begin{array}{ccc}{L}_2\left(x_2\right)& L_3\left(x_2\right)& \begin{array}{cc}\cdots & L_p\left(x_2\right)\end{array}\\ L_2\left(x_3\right)& L_3\left(x_3\right)& \begin{array}{cc}\cdots & L_p\left(x_3\right)\end{array}\\ \begin{array}{c}⋮\\ L_2\left(x_p\right)\end{array}& \begin{array}{c}⋮\\ L_3\left(x_p\right)\end{array}& \begin{array}{cc}\begin{array}{c}\ddots \\ \cdots \end{array}& \begin{array}{c}⋮\\ L_p\left(x_p\right)\end{array}\end{array}\end{array}\right]\left[\begin{array}{c}{v}_2\left(0\right)\\ v_3\left(0\right)\\ \begin{array}{c}⋮\\ v_p\left(0\right)\end{array}\end{array}\right]=\left[\begin{array}{c}{u}^0\left(x_2\right)\\ u^0\left(x_3\right)\\ \begin{array}{c}⋮\\ u^0\left(x_p\right)\end{array}\end{array}\right]$

$\Rightarrow \left[\begin{array}{c}{v}_2\left(0\right)\\ v_3\left(0\right)\\ \begin{array}{c}⋮\\ v_p\left(0\right)\end{array}\end{array}\right]=v\left(0\right)=\left[\begin{array}{c}{u}^0\left(x_2\right)\\ u^0\left(x_3\right)\\ \begin{array}{c}⋮\\ u^0\left(x_p\right)\end{array}\end{array}\right]$ (19)

Solution of Equation (18) will give us the unknowns $v_n\left(t\right);n=2,3,\cdots ,p$ in trial solution of Equation (7). Approximate solution of system of ordinary differential equation in Equation (18) with its initial condition in Equation (19) can be obtained by very well-known Euler’s method. Instead of continuous approximation to the solution $v\left(t\right)$ , approximations will be generated at mesh points

$t_j>0$ . With step size $\Delta t=\frac{T}{n_t}>0;n_t\in \mathbb{N}$ , we define the mesh points as

$t_j=j\Delta t;j=0,1,2,3,4,\cdots ,n_t$

Then Euler’s method becomes

$v\left(t_0\right)=\left[\begin{array}{c}{u}^0\left(x_2\right)\\ u^0\left(x_3\right)\\ \begin{array}{c}⋮\\ u^0\left(x_p\right)\end{array}\end{array}\right]$

$v\left(t_{j+1}\right)=\left[I+hB\left(t_j\right)\right]v\left(t_j\right)+hs\left(t_j\right)+hse\left(t_j\right);j=0,1,2,3,\cdots ,\left(n_t-1\right)$ (20)

Equation (20) is the difference equation for the Euler’s method.

Finally substituting approximations of $v\left(t\right)$ at various mesh points $t_j$ into the trial solution at Equation (10) will produce the approximations of $u\left(x,t\right)$ at mesh points $t_j$ as $u_p\left(x,t_j\right)$ .

Error Calculation: Here we discuss only error calculation for the above method. Our main objective is to calculate local errors and global errors. At $t=T$ we define the absolute local error function as

$\epsilon \left(x\right)=\left|u_p\left(x,T\right)-u\left(x,T\right)\right|$

Then we define the global error at $t=T$ as

$\text{global\_error}=\underset{a}{\overset{b}{{\displaystyle \int }}}\epsilon \left(x\right)\text{d}x$

In this study, we focus only on the numerical results and the resultant errors. Details of error analysis caused by different types of collocation points and Euler methods for system of ODE are left for further research.

4. Numerical Comparisons and Discussions

Now, we apply spectral collocation method discussed in the previous section with $p=5$ to solve fractional diffusion equation with four different sets of collocation points and will compare the results obtained. The four sets of collocation points are generated from Gauss-Lobatto nodes, roots of Chebyshev polynomials of first kind, roots of Legendre polynomials and equally spaced nodes over the space domain. The four sets of points are generated by the following way:

From Gauss-Labatto nodes, to generate the points on the interval $\left[a,b\right]$ we consider:

$x_k=\frac{b-a}{2}+\frac{b-a}{2}\cos \frac{\left(k-1\right)\text{π}}{p};k=1,2,\cdots ,p+1$

Let $y_k$ be the roots of the Chebyshev polynomials of first kind $T_{p-1}\left(y\right)$ with $y_k<y_{k+1}$ for $k=1,2,\cdots ,p-2$ . Since roots of $T_{p-1}\left(y\right)$ are in the interval $\left(-1,1\right)$ we consider the following nodes along with $x_1=a$ and $x_{p+1}=b$

$x_{k+1}=\frac{b-a}{2}+\frac{b-a}{2}{y}_k;k=1,2,\cdots ,p-1$

Similarly, let $w_k$ be the roots of the Legendre polynomial $P_{p-1}\left(w\right)$ with $w_k<w_{k+1}$ for $k=1,2,\cdots ,p-2$ . Like Chebyshev polynomials, roots of $P_{p-1}\left(w\right)$ are in the interval $\left(-1,1\right)$ we consider the following nodes along with $x_1=a$ and $x_{p+1}=b$

$x_{k+1}=\frac{b-a}{2}+\frac{b-a}{2}{w}_k;k=1,2,\cdots ,p-1$

For equally spaced nodes over the space domain we consider the following points

$x_k=a+\left(k-1\right)\frac{b-a}{p};k=1,2,\cdots ,p+1$

Now we will consider the performance of these four sets of collocation points into proposed spectral collocation scheme with two examples of space fractional diffusion equations. Since in both examples $a=0$ and $b=1$ , with $p=5$ the above four sets of points are calculated as follows:

$\begin{array}{l}\text{gauss\_lobatto}=\left\{0,0.0954915,0.345492,0.654508,0.904508,1\right\}\\ \text{chebyshev}=\left\{0,0.0380602,0.308658,0.691342,0.96194,1\right\}\\ \text{legendre}=\left\{0,0.0694318,0.330009,0.669991,0.930568,1\right\}\\ \text{equally\_spaced}=\left\{0,0.2,0.4,0.6,0.8,1\right\}\end{array}$

The exact solution of both examples can be verified by using Equation (2).

Example 1: We consider the following fractional diffusion equation used by Bahsi and Yalcinbas [7] :

$\frac{\partial u\left(x,t\right)}{\partial t}=d\left(x,t\right)\frac{{\partial }^{1.5}u\left(x,t\right)}{\partial x^{1.5}}+s\left(x,t\right);0<x<1,0\le t\le 2$

$u\left(x,0\right)=\left(x^2+1\right)\sin 1,0<x<1$

$u\left(0,t\right)=\sin \left(t+1\right)\text{and}u\left(1,t\right)=2\sin \left(t+1\right)$

where,

$d\left(x,t\right)=\Gamma \left(1.5\right)x^{0.5}$

$s\left(x,t\right)=\left(x^2+1\right)\cos \left(t+1\right)-2x\sin \left(t+1\right)$

The exact solution of Example 1 is $u\left(x,t\right)=\left(x^2+1\right)\sin \left(t+1\right)$ .

Using proposed spectral collocation scheme, absolute local errors for the four different sets of collocation points along with exact solution at $T=2$ using $\Delta t=0.0025$ are given in Table 1.

The absolute local error curves for the four sets of points are given in Figure 1.

Four sets of collocation points give following global errors for Example 1:

$\begin{array}{l}{\text{global\_error}}_{\text{Gauss-Lobatto}}=0.0000360065\\ {\text{global\_error}}_{\text{Chebyshev}}=0.0000149475\\ {\text{global\_error}}_{\text{Legendre}}=0.0000312289\\ {\text{global\_error}}_{\text{Equally\_spaced}}=0.0000407057\end{array}$

Example 2: We consider the following fractional diffusion equation used by Huang and Zheng [9] :

$\frac{\partial u\left(x,t\right)}{\partial t}=d\left(x,t\right)\frac{{\partial }^{1.2}u\left(x,t\right)}{\partial x^{1.2}}+s\left(x,t\right);0<x<1,0\le t\le 2$

$u\left(x,0\right)=x^2,0<x<1$

$u\left(0,t\right)=0\text{and}u\left(1,t\right)={\text{e}}^{-t}$

where,

$d\left(x,t\right)=\frac{\Gamma \left(1.8\right)x^{2.2}}{2}$

$s\left(x,t\right)=x^2\left(x+1\right){\text{e}}^{-t}$

The exact solution of Example 2 is $u\left(x,t\right)=x^2{\text{e}}^{-t}$ .

Using proposed spectral collocation scheme absolute local errors for the four different sets of collocation points along with exact solution at $T=2$ using $\Delta t=0.0025$ are given in Table 2.

The absolute local error curves for the four sets of points are given in Figure 2.

Four sets of collocation points gives following global errors for Example 2:

$\begin{array}{l}{\text{global\_error}}_{\text{Gauss-Lobatto}}=0.000129786\\ {\text{global\_error}}_{\text{Chebyshev}}=0.000139473\\ {\text{global\_error}}_{\text{Legendre}}=0.000137378\\ {\text{global\_error}}_{\text{Equally\_spaced}}=0.000202811\end{array}$

From absolute local error curves and global errors of both examples we can say that spectral collocation method based on Lagrange’s basis polynomials give very satisfactory approximations to the solution of space fractional diffusion equation. Though we used lower order and very simple Euler’s method to solve resultant system of ODE but yet get very satisfactory approximations. The accuracy can be improved by using higher order method than Euler’s method to solve system of ODE. About performance of different collocation points it is evident that there is no way to declare which one is better since there are variations among absolute local errors over the space domain for different sets of

points. Even if we consider global error, performance of different sets of collocation points vary from one problem to another. We observe another intriguing feature, that is, in case of Example 1 errors due to Chebyshev nodes and in case of Example 2 errors due to equally spaced nodes fluctuates over the domain where other nodes does not show such fluctuations.

5. Conclusion

There are several spectral collocation methods available for different types of partial fractional differential equations. We discussed spectral collocation method based on Lagrange’s basis polynomials for 1D space fractional diffusion equation where our proposed form of trial solution can handle non-zero Dirichlet’s boundary conditions on both boundaries effectively. Properties of Lagrange’s basis polynomials reduce the volume of the calculations needed and result simpler equations. We implemented the method into two examples with four different sets of collocation points and found excellent match with exact solution in each case. We compared absolute local errors and global errors for each set of points. No clear conclusion can be drawn about which set of points give better approximations for space fractional diffusion equation in spectral collocation method which we discussed here.

Cite this paper

Nova, M.H., Molla, H.U. and Banu, S. (2017) Comparison of Numerical Approximations of One-Dimensional Space Fractional Diffusion Equation Using Different Types of Collocation Points in Spectral Method Based on Lagrange’s Basis Polynomials. American Journal of Computational Mathematics, 7, 469-480. https://doi.org/10.4236/ajcm.2017.74034

References