Pseudo-spectral Method for Space Fractional Diffusion Equation

This paper presents a numerical scheme for space fractional diffusion equations (SFDEs) based on pseudo-spectral method. In this approach, using the Guass-Lobatto nodes, the unknown function is approximated by orthogonal poly-nomials or interpolation polynomials. Then, by using pseudo-spectral method, the SFDE is reduced to a system of ordinary differential equations for time variable t. The high order Runge-Kutta scheme can be used to solve the system. So, a high order numerical scheme is derived. Numerical examples illustrate that the results obtained by this method agree well with the analytical solutions.


Introduction
Anormalous diffusion model, where a particle spreads at a rate inconsistent with the classical Brownian motion model, has been applied in many fields, such as in fractured and porous media, in chaotic or turbulent processes [1][2][3][4][5][6].It is known that anomalous diffusion processes can be described by fractional partial differential equations where the fractional derivatives are defined as the Riemann-Liouville's representation.The former is referred as time fractional diffusion equation and the latter as space fractional diffusion equation.The difference between two cases can be shown using the interpolation of operator.It is well-known that the unknown   , u x t denotes the probability density function, which is the probability of finding the particle at position x and at time .For time fractional diffusion Equation (1.1), it can be rewritten as (under suitable conditions) where f  depends only on f .Hence, from the view point of mathematics, the Equation (1.3) can be thought as the interpolation of the equations in the case 0   and the case 1   .In a word, the particles described by the Equation (1.3) are diffused slowly in comparison with the classic situations, namely sub-diffusion.
On the other hand, the fractional derivative of Equation (1.2) can similarly be thought as the interpolation of operators x   and 2 2 x   .Note that u x   describes the free transport of particles and 2 u x   2 , the diffusion, due to the collisions of particles.So, in the model of space fractional diffusion the velocity of particle is larger than one of ordinary diffusion, namely super-diffusion.
Several methods have been developed for numerical solving the fractional diffusion equation.Langlands and Henry studied sub-diffusion equation based on L1 scheme [7].The authors proposed an implicit difference scheme which is unconditional stable.Based on Grünwald-Letnikov formula, an explicit difference scheme has been presented for sub-diffusion equation [8,9].
In contrast to sub-diffusion equation, numerical solu-tion of supper-diffusion equations have been studied and some methods are also developed.It is interesting to note that the explicit and implicit difference schemes based on Grünwald-Letnikov formula are all unconditional unstable for super-diffusion Equation (1.2) and advectiondiffusion equation [10,11].However, the authors proved that a shifted Grünwald-Letnikov formula can produce stable difference scheme [10,11].
The finite difference method is the most classic method for fractional differential equation [12][13][14][15][16].The recent works can see the references [17][18][19].However, high order accuracy schemes are seldom derived by finite difference method.In general, extrapolation method is applied in order to obtain a high accuracy [10,16].It is well-known that spectral methods are superior to finite different methods in many instances for partial differential equations [20][21][22][23][24].In the recent paper [25], the authors presented a spectral method to calculate the fractional derivative and integral, and studied the numerical solution of differential equations by the spectral method.
In this paper, we shall deal with the numerical solution of the one-dimensional variable coefficients space fractional diffusion equation with the initial condition for and boundary conditions and Here, is variable coefficient and source term.

  d x
Chebyshev spectral collocation method (often called pseudo-spectral) is used in this paper.By employing the Chebyshev polynomials and Gauss-Lobbato nodes, the unknown is approximated by using the orthogonal projection and interpolation.So, the spatial fractional derivatives are easily computed and a system of ordinary differential equations in time can be derived.Then, high order Runge-Kutta method can be utilized and a high order accuracy scheme is obtained.In a comparison with the finite difference method, the Chebyshev spectral method used here shows remarkably superiority in terms of accuracy and the number of grid points required.
The paper is arranged as follows.The Section 2 introduces some basic concepts of fractional derivative and Jacobi orthogonal polynomials and their properties.The third section presents the spectral method for calculating the Riemann-Liouville fractional derivative.The pseudospectral method and its implement are proposed in Section 4. In Section 5, several numerical examples are provided.These numerical examples illustrate the high accuracy and efficiency of our method.Finally, we give the conclusion in Section 6.

Fractional Derivative
The Riemann-Liouville fractional integral of order α(α > 0) for casual function   f x is defined by where   .


is the Gamma function.Riemann-Liouville fractional integral is an analogue of the well-known Cauchy formula, which reduces the n-fold integration     The Riemann-Liouville representation of fractional derivative of order   By the part integration formula, one can easily derive the following properties.

Jacobi Orthogonal Polynomials
Given a weight function   x  , the orthogonal polynomials sequence denotes weighted inner product of a Hilbert space.
For , the Jacobi polynomials is a orthogonal polynomials sequence , 1 where The following is some useful properties of Jacobi polynomials that will be used in present paper [26] (also refer to [25] and references therein).
Proposition 2.3. 1) High order derivative for Jacobi polynomials where 2) Expression by the derivative of Jacobi polynomials where By the previous Jacobi polynomials some special orthogonal polynomials can be obtained.Legendre and Chebyshev polynomials are two important polynomials of the special cases of the Jacobi polynomials.For the case of , the corresponding polynomials is said to Legendre polynomials; and the case of

Spectral Approximation to Fractional Derivative
Let N H denote the set of polynomials of degree not exceeding N. It is clear that Denote the weighted inner product for weight function and the weighted Lebesgue space Define the orthogonal projection  can be approximated by the orthogonal projection The expression coefficients is determined by Therefore, the fractional derivative can be approximated by can be the computed by recurrence formula [25].Clearly, By the three-term recurrence relation (2.5), In light of (2.10) and notice that , , as (2.11).Therefore, we obtain by (2.3) and (2.9) that Remark 1.By the standard theory of orthogonal projection, the spectral accuracy to approximate the fractional derivative

  D u x
 can be obtained (refer to the references [25,26]).

Collocation Method
This section presents the Chebyshev collocation method (often called pseudo-spectral method) for solving the space fractional diffusion Equation (1.4).Firstly, we shall consider the Gauss-Lobatto nodes in order to obtain the collocation equation.
Consider the Chebyshev polynomials   n T x of the first kind, which are the special case of Jacobi polynomials for 1 1 x    .Chebyshev polynomials are the weighted orthogonal polynomials with weight function One can easily obtain the three-term recurrence formula , here p, q is determined by satisfying .Then the roots of is named as the Chebyshev Gauss-Lobatto points (in this paper is also called Gauss-Lobatto points).Gauss-Lobatto points can be explicitly written as [26,27] where k  are the associated weights of Chebyshev Gauss-Lobatto integration 2 for 0, , , 1 for 1, , 1.
By orthogonality it can easily be derived that Next, let us consider the following Gauss-Lobatto points cos , 0,1, 2, , .
Assume that the approximate solution of the space fractional diffusion Equation (1.4) has separable formation with initial condition , where In addition, the boundary conditions for x and N x can be expressed by Equations (4.5) and (4.6) is called the collocation equation for fractional diffusion Equation (1.4), and k x defined by (4.3) is called the collocation points.
The Equation (4.5) is a system of ordinary differential equations.The calculation of can be implemented by  can be conveniently calculated by making use of (3.4).Remark 2. 1) Equations (4.5) and (4.6) is also said to be the strong form of the collocation method.
2) In order to get high order schemes, high order Runge-kutta methods can be used to solve the system of ordinary differential Equation (4.5).Now, we employ the discrete inner product to deal with the initial condition.Notice that Therefore, making use of (4.2) it gives In order to conveniently deal with the boundary conditions and initial condition, let us consider another form of the collocation method based on interpolation.Based on the Gauss-Lobatto nodes (4.3) above, the Lagrangian interpolation basis function Assume that the approximation be the interpolation where In order to compute conveniently the coefficient matrix M, we expanse Therefore, making use of Prop.2.1 and Prop.2.2 one can derive.

Numerical Examples
In this section, we consider the space fractional diffusion equations for different source terms and the values of  utilizing the collocation method.Here, the fourth order method is used for all the examples with time step . The first example is a fractional equation with 1

Conclusion
The collocation method, namely called pseudo-spectral method, is proposed in present paper.This kind of method can be efficiently applied to fractional partial differential equations.The remarkably superiorities, efficiency and high accuracy have been found through the numerical examples presented in this paper.The high accurate approximation only by a few grids can be derived.However, the space and time steps must satisfy certain condition in order to guarantee the stability for finite difference method.It must be stressed that the higher order collocation method based on Lagrange interpolation maybe result in numerical instability because of the ill-condition of the matrix .It is hoped that the error estimated on the spectral method for fractional differential equations will be studied in future work.

N  Example 1 .
and time variable belongs to   0,  .The second example is associated to the problem of time belonging to finite interval with 1.5 2    .For the first example, we use the second kind of collocation method, namely based on the Lagrange interpolation, with 3 N  .The second example is solve by the Chebyshev polynomials approximation with .3 Consider the following asymptotic problem with 1.2