Lagrange’s Spectral Collocation Method for Numerical Approximations of Two-Dimensional Space Fractional Diffusion Equation ()
1. Introduction
Theory of Fractional calculus is almost of same age as that of the classical calculus. The concept of a fractional derivative initiates the theory of fractional calculus when this concept first appeared in a letter of L’Hospital to Leibnitz in 1695 [1] . But in 1819, Lacroix first mentioned the derivative of arbitrary order in a text [2] . Later from 1832 Liouville [3] starts to provide the initial foundations of the fractional calculus. A brief history and fundamental theory of fractional calculus can be found in [4] . Though many important physical processes in almost all branches of sciences can be described using various tools of classical calculus, there are also many complex phenomena that cannot be modeled using classical calculus. Being one of the most important generalizations of classical calculus, fractional calculus has been used to describe many such complex processes recently. In last few decades, scientists and engineers have found fractional order derivative very powerful to describe problems in anomalous diffusion, signal processing, control processing, fractional stochastic system and viscoelasticity. In space fractional diffusion equation, spatial derivative of order
replaces the classical second order spatial derivative in classical diffusion equation and results into super diffusive flow model. Unlike classical order derivative, fractional derivative has several kinds of definitions such as Caputo, Riemann-Liouville, Grünwald-Letnikov and Riesz. Detailed discussion about various definitions of fractional derivative can be found in [5] [6] .
Spectral method is very effective and efficient and hence popular among scientists and engineers for numerical approximations of both classical and fractional differential and integral equations. In spectral method unknown solution is first approximated by the trial solution as a finite sum of a set of known basis functions. Then expansion coefficients in trial solution can be determined by several techniques like collocation, Galerkin, Petrov Galerkin, tau etc. In all of these techniques expansion coefficients are determined by minimizing the difference between exact and trial solutions. For spectral approximations of fractional differential equations collocation technique is the most efficient and widely used technique. In collocation technique expansion coefficients are determined by forcing the trial solution to satisfy the differential equation at some suitably chosen points from the domain known as collocation points. For a p parameter solution p collocation points are required within the domain and which results into p residual equations. In recent times, researchers showed immense interests to try out various types of polynomials as basis sets in spectral approximations and different kinds of collocation points. For numerical approximations of 1D space fractional diffusion equations Khader [7] used Chebyshev polynomials of space variable as the basis set in spectral approximations and roots of shifted Chebyshev polynomials as the collocation points to determine the expansion coefficients. Chebyshev polynomials are also used for both time and space domain as basis set in spectral approximations of 1D fractional diffusion equation by Azizi and Loghmani [8] and Xie et al. [9] . But Azizi and Loghmani applied collocation technique with Gauss-Lobatto nodes whereas Xie et al. used tau method to determine expansion coefficients. For numerical approximations of 1D fractional diffusion equation Bahsi and Yalcinbas [10] chosen Fibonacci polynomials to express the trial solution in both space and time domain and then used evenly spaced collocation points. Pirim and Ayaz [11] also chosen evenly spaced collocation points but expressed the trial solution in terms of Hermite polynomials for approximations of fractional order system of differential equations. Huang and Zheng [12] used Jacobi polynomials and Lagrange’s basis polynomials for spectral approximations of 1D space fractional diffusion equations and used Gauss-Lobatto nodes for collocation. Doha et al. [13] presented spectral tau method for 1D space fractional diffusion equation where they expand the solution by Jacobi polynomials. For spectral expansion of trial solution for time fractional diffusion equations Lagrange interpolation polynomials are used in both 1D space and time domain by Huang [14] . For collocation purpose he chose Jacobi-Gauss nodes for time domain and for space domain he chosen Gauss-Lobatto nodes. Zayernouri and Karniadakis [15] introduced fractional Lagrange interpolants and developed a new fractional spectral collocation method for 1D fractional differential equations. Krishnaveni et al. [16] proposed a hybrid method for 1D space fractional diffusion equation where trial solution is approximated by spectral expansion of fractional shifted Legendre polynomials. To approximate the solution of 2D space fractional diffusion equation by spectral collocation method Bhrawy [17] used shifted Legendre polynomials and Gauss-Lobatto nodes for 2D space domain. Other numerical methods such as ADI and Finite-Difference methods [18] [19] [20] [21] [22] are also very efficient for 2D problems.
In this research, we consider two dimensional space fractional diffusion equation with an initial condition and non-homogeneous Dirichlet’s boundary conditions where spatial fractional derivatives are described in Riemann-Liouville sense. Variety of spectral collocation methods are available for 1D fractional diffusion equations but spectral collocation method haven’t been tried with lots of variations in 2D problems. To the best of our knowledge Lagrange’s basis polynomials haven’t been used in spectral collocation method for 2D space fractional diffusion equations. Hence we use Lagrange’s basis polynomials for spectral expansion of trial solution in 2D space domain. To determine the expansion coefficients collocation technique is used which transforms the diffusion equation into a first order system of ODE of time variable. System of ODE is then solved by fourth order Runge-Kutta [23] method. Lagrange’s polynomials depends on the selected nodes. Also in collocation method choice of collocation points affect the resultant accuracy. So we consider four different types of nodes in Lagrange’s spectral collocation method to generate Lagrange’s basis polynomials and as collocation points to investigate their performances in the proposed method. Four different sets of nodes are derived from shifted roots of Chebyshev polynomials of 1st and 2nd kinds, Legendre polynomials and Jacobi polynomials. The performances in terms of resultant accuracy by four sets of nodes in proposed Lagrange’s spectral collocation method are illustrated and compared through numerical examples. Similar investigation on 1D space fractional diffusion equation is carried out by Nova et al. [24] recently.
Remaining parts of this paper are organized as follows: Section 2 contains the preliminaries of different polynomials and fractional derivative. Then in Section 3, detailed formulations of Lagrange’s spectral method for 2D space fractional diffusion equation are given. Next in Section 4, numerical examples are considered to demonstrate the accuracy of the proposed method and to compare the performance of four types of nodes. Finally, conclusion is drawn about this research in Section 5.
2. Preliminaries
In this section, we present a brief introduction of Legendre, Jacobi and Chebyshev polynomials of 1st and 2nd kinds. Also short description of Lagrange’s basis polynomials and Caputo’s fractional order derivative are given.
Legendre Polynomials: The nth degree Legendre polynomial
over the domain
are explicitly defined as
(1)
Chebyshev Polynomials of 1st Kind: The nth degree Chebyshev polynomials of 1st kind
over the domain
are explicitly defined as
(2)
Chebyshev Polynomials of 2nd Kind: The nth degree Chebyshev polynomials of 2nd kind
over the domain
are explicitly defined as
(3)
Jacobi Polynomials: The nth degree Jacobi polynomials
on the interval
are defined as
(4)
Lagrange Basis Polynomials: For
points
Lagrange basis polynomials
each of degree p are defined as follows:
(5)
with the property
, where
is the Kronecker delta function. Here
is the derivative of
.
Riemann-Liouville Fractional Derivative: Riemann-Liouville fractional derivative operator of order α is denoted by
and defined by:
(6)
with
.
Then for
and
(7)
where
.
Like classical integer order derivative, Riemann-Liouville fractional order derivative is also a linear operator.
3. Lagrange’s Spectral Collocation Method
Lagrange’s Spectral Collocation method for numerical approximations of 2D space fractional diffusion equation will be presented in this section. In this research, we consider following form of diffusion equation
(8)
subject to the initial condition
(9)
and non-homogeneous Dirichlet’s boundary conditions
(10)
(11)
Here
and
are the diffusion coefficients,
is the force function and
is the unknown function to be solved.
For numerical approximations of Equation (8) by Lagrange’s spectral collocation method we first divide the interval
of space domain x into p subintervals and interval
of y domain into q subinterval by means of following nodes:
(12)
(13)
These nodes can be chosen with no specific pattern and from anywhere in the domain. Then using these nodes
from x-domain we form Lagrange’s basis polynomials
each of degree p and using nodes
from y-domain we form Lagrange’s basis polynomials
each of degree q. Now at these nods Lagrange’s basis polynomials satisfies
and
. These nodes are also used as collocation points in the later part of this section.
Now as spectral approximation of
, the solution of Equation (8) we consider trial solution
in the form
(14)
(15)
Here Lagrange’s basis polynomials
and
are used as the basis
functions in spectral approximations and
are the expansion
coefficients needed to be determined. Now we write the residual
for Equation (8) as
(16)
where
(17)
(18)
(19)
Here
is the
order Riemann-Liouville fractional derivative of
with respect to x and
is the
order Riemann-Liouville fractional derivative of
with respect to y.
To determine the expansion coefficients we chose collocation technique and for this we force the residual equation at Equation (16) to be zero at
interior points
of the spatial domain and
force the trial solution to satisfy the boundary conditions Equation (10) & (11) at the boundary lines.
First by setting the residual at Equation (16) to be zero at the interior points of the spatial domain, we have the equation
(20)
for
and
with
(21)
(22)
(23)
where
(24)
(25)
Thus Equation (20) gives us the following system of ODE for time variable t
(26)
with
(27)
Now forcing the trial solution to satisfy the initial condition Equation (9) at the interior points gives us the initial condition for system of ODE in Equation (26) as
(28)
(29)
Now forcing the trial solution to agree with the boundary conditions for boundary lines
and
, that is, to satisfy the boundary conditions at Equation (10) gives
(30)
(31)
Then at the points
Equation (30) and (31) becomes
and
for
(32)
Again forcing the trial solution to agree with the boundary conditions, now for boundary lines
and
, that is, to satisfy the boundary conditions at Equation (11) gives
(33)
(34)
Now at the points
Equation (33) and (34) becomes
and
for
(35)
System of ODE in Equation (26) with its initial condition Equation (29) can
be solved for the value of expansion coefficients
at any
.
Rest of the values of expansion coefficients can be obtained from Equation (32) and (35). Replacing the values of expansion coefficients
for any particular
into trial solution in Equation (14) will gives us the approximate solution
of Equation (8) in terms of Lagrange’s basis polynomials
and
.
4. Numerical Results
In this section, we consider two examples to apply the proposed method with four different types of nodes to investigate the accuracy of the method and nodes. System of ordinary differential equations in Equation (26) with initial condition
Equation (29) is solved for the value of expansion coefficients
at
by 4th order Runge-Kutta method with time step size 0.01.
In both examples spatial domain was
and we consider
. Four different types of nodes are thus derived as follows:
Since all the roots of Legendre, Jacobi and both kinds of Chebyshev polynomials are in the interval
, to generate nodes
over any arbitrary interval
, first we choose
and
. Then, the remaining
nodes are generated by shifting the roots of one of the above
degree polynomials into the interval
.
Since in both examples 2D space domain is square shaped and we consider
, following four sets of nodes are thus used for both x and y domain.
To compare the approximations given by the proposed method with the exact solution we consider absolute local error for any given
given by
We also consider maximum absolute error given by
For numerical simulation we develop codes for the proposed method using MATLAB.
Example 1: Consider the diffusion equation [17] [18]
with boundary conditions
and with initial condition
where
The exact solution of this problem is given by
.
Applying the proposed Lagrange’s spectral collocation technique for
into Example 1 with four different types of nodes resulting surface graph of absolute local error
are given into Figures 1-4. Then Figure 5 contains the absolute local error curve
and Figure 6 contains the absolute local error curve
for example 1.
Example 2: Consider the diffusion equation [17]
with boundary conditions
and with initial condition
where
Figure 1. Error graph for Chebyshev 1st kind.
Figure 2. Error graph for Chebyshev 2nd kind.
Figure 5. Absolute local error curve for
.
Figure 6. Absolute local error curve for
.
The exact solution of this problem is given by
.
Now by proposed Lagrange’s spectral collocation technique with
, surface graph of absolute local error
by four different types of nodes for Example 2
are given into Figures 7-10. Then Figure 11 contains the absolute
Figure 7. Error graph for Chebyshev 1st kind.
Figure 8. Error graph for Chebyshev 2nd kind.
Figure 11. Absolute local error curve for
.
local error curve
and Figure 12 contains the absolute local error curve
for Example 2.
Table 1 contains maximum absolute error
by four different sets of nodes for both examples.
Observing the surface graphs and error curves for absolute local error in both examples it is evident that proposed Lagrange’s collocation technique gives very high accuracy approximations of 2D space fractional diffusion equations for all four types of collocation points. The maximum absolute errors in both examples for Jacobi nodes, Chebyshev 2nd, Legendre and Chebyshev 1st nodes are found in increasing order. The same performance of these nodes is clearly visible in Figure 5, Figure 11 and Figure 12. From Figure 5 we see that in the last quarter of the x-domain absolute local error for Chebyshev 1st nodes decreased rapidly and for Legendre nodes error decreased slowly whereas errors by both Jacobi and Chebyshev 2nd nodes keeps increasing. But in first three quarter of the
Figure 12. Absolute local error curve for
.
x-domain performance of Jacobi nodes, Chebyshev 2nd, Legendre and Chebyshev 1st nodes are found in decreased order.
5. Conclusion
Spectral collocation method with Lagrange’s basis polynomials for approximations of 2D space fractional diffusion equations is introduced in this research. We used four different types of nodes generated from roots of Legendre, Jacobi, 1st and 2nd kinds of Chebyshev polynomials. Nodes are used to form Lagrange’s basis polynomials for spectral approximations and then as collocation points to determine the expansion coefficients. Performances of these four types of nodes in Lagrange’s spectral collocation method are investigated through two examples of 2D space fractional diffusion equation. Four types of nodes in proposed method give very high accuracy approximations. We found that among four types of nodes Chebyshev 1st kind gives lowest accuracy and Jacobi nodes give highest accuracy. We compared approximate solutions with exact solution and considered absolute local error and maximum absolute error. Four different types of nodes showed consistent performance in both examples.