1 Introduction
There have been increasing interests in the description of the physical and chemical processes by means of equations involving fractional derivatives over the last decades. And, fractional derivatives have been successfully applied into many sciences, such as physics b1 , biology b2 , chemistry b3 , hydrology b4 ; b5 ; b6 ; b7 , and even finance b8 . In groundwater hydrology the fractional advectiondispersion equation (FADE) is utilized to model the transport of passive tracers carried by fluid flow in a porous medium b6 ; b9 .
Considerable numerical methods for solving the FADE have been proposed. Kilbas et al. Kilbas06 introduced the theory and applications of fractional differential equations. Meerschaert and Tadjeran b10 developed practical numerical methods to solve the onedimensional space FADE with variable coefficients on a finite domain. Liu et al. b6
transformed the space fractional FokkerPlanck equation into a system of ordinary differential equations (method of lines), which was then solved using backward differentiation formulas. Momani and Odibat
b9 developed two reliable algorithms, the Adomian decomposition method and variational iteration method, to construct numerical solutions of the spacetime FADE in the form of a rapidly convergent series with easily computable components. Zhuang et al. b11 discussed a variableorder fractional advectiondiffusion equation with a nonlinear source term on a finite domain. Liu et al. b12 proposed an approximation of the LévyFeller advectiondispersion process by employing a random walk and finite difference methods. In addition, other finite difference methods b13 , finite element method b23 , finite volume method b24 , homotopy perturbation method b25 and spectral method b26 ; b27 are also employed to approximate the FADE.In this paper, we consider the following RFADE:
(1) 
subject to the initial condition:
(2) 
and the zero Dirichlet boundary conditions:
(3) 
where , , and represent the average fluid velocity and the dispersion coefficient. The Riesz space fractional operators and on a finite domain are defined respectively as
where , , and
where represents the Euler gamma function.
The fractional kinetic equation (1) possesses a physical meaning (see b21 ; b22 for further details). Physical considerations of a fractional advectiondispersion transport model restrict , , and we assume and so that the flow is from left to right. In the case of and , Eq.(1) reduces to the classical advectiondispersion equation (ADE). In this paper, we only consider the fractional cases: when , Eq.(1) reduces to the Riesz fractional diffusion equation (RFDE) b21 and when , the Riesz fractional advectiondispersion equation (RFADE) is obtained b22 .
For the RFADE (1), Anh and Leonenko b18 presented a spectral representation of the meansquare solution without the nonhomogeneous part and for some range of values of the parameters. Later, Shen et al. b19 derived the fundamental solution of Eq.(1) and discussed the numerical approximation of Eq.(1) using finite difference method with first convergence order. Another method based on the spectral approach and the weak solution formulation was given in Leonenko and Phillips b20 . In addition, Zhang et al. b14 use the Galerkin finite element method to approximate the RFADE. Besides, Yang et al. b28 applied the approximation method, the standard/shifted Grünwald method, and the matrix transform method (MTM) to solve the RFADE. And, Ding et al. b29 also consider the numerical solution of the RFADE by using improved matrix transform method and the (2, 2) Pade approximation. Most of the numerical methods proposed by these authors are low order or lack stability analysis.
In this paper, based on the weighted and shifted Grünwald difference (WSGD) operators to approximate the Riesz space fractional derivative, we obtain the second order approximation of the RFADE. Furthermore, We propose the finite difference method for the RFADE and obtain the CrankNicolson scheme. Moreover, we prove that the CrankNicolson scheme is unconditionally stable and convergent with the accuracy of and improve the convergence order to by applying the Richardson extrapolation method.
The outline of the paper is as follows. In Section 2, the WSGD operators and some lemmas are given. In Section 3, we first present the finite difference method for the RFADE, and then derive the CrankNicolson scheme. We proceed with the proof of the stability and convergence of the CrankNicolson scheme in Section 4. Besides, we further improve the convergence order by applying the Richardson extrapolation method. In order to verify the effectiveness of our theoretical analysis, some numerical examples are carried out and the results are compared with the exact solution in Section 5. Finally, the conclusions are drawn.
2 The approximation for the RiemannLiouville fractional derivative
First, in the interval , we take the mesh points , , and , , where , , i.e., and are the uniform spatial step size and temporal step size. Now, we give the definition of the RiemannLiouville fractional derivative.
Definition 1 (b15 ).
The order left and right RiemannLiouville fractional derivatives of the function on , are given by

left RiemannLiouville fractional derivative:

right RiemannLiouville fractional derivative:
Generally, the standard GrünwaldLetnikov difference formula is applied to approximate the RiemannLiouville fractional derivative. Meerschaert and Tadjeran b10 showed that the standard GrünwaldLetnikov difference formula was often unstable for time dependent problems and they proposed the shifted Grünwald difference operators
whose accuracies are first order, i.e.,
where are integers and . In fact, the coefficients are the coefficients of the power series of the function ,
for all , and they can be evaluated recursively
Lemma 1 (b13 ).
Suppose that , then the coefficients satisfy
Inspired by the shifted Grünwald difference operators and multistep method, Tian et al. b16 derive the WSGD operators:
Lemma 3 (b16 ).
Supposing that , let , and
and their Fourier transforms belong to
, then the WSGD operators satisfyuniformly for , where , are integers and .
Tian et al. prove Lemma 3 under the additional conditions that . In fact, their proof also holds for the case. The proof proceeds the same as b16 , hence we will not repeat it here.
Remark 1.
Considering a well defined function on the bounded interval , if or , the function can be zero extended for or . And then the order left and right RiemannLiouville fractional derivatives of at each point can be approximated by the WSGD operators with second order accuracy
where and .
When , and , the discrete approximations for the RiemannLiouville fractional derivatives on the domain are
(4)  
(5)  
(6)  
(7) 
where
(8) 
(9) 
Now, we discuss the properties of the coefficients and .
Lemma 4.
Suppose that , then the coefficients satisfy
Proof.
Combining the definition of and the property of , it is easy to derive the value of , and . When , by the definition of
we have as for and . Moreover,
Recalling Lemma 1, when , we obtain when , i.e., . For the sum , we have
Since when , for . ∎
Lemma 5 (b16 ).
Suppose that , then the coefficients satisfy
3 The finite difference method for the RFADE
In this section, we utilize the Eqs.(4)(7) to approximate the Riesz space fractional derivative and derive the CrankNicolson scheme of the equation. We define , , let be a finite domain, setting be a uniform partition of , which is given by for , where and are the time and space steps, respectively. First, we present the semidiscrete form of Eq.(1),
(10) 
Let be the approximation solution of . Substituting (4)(7) into (10), we obtain
(11) 
where and . Denote
(12) 
(13) 
and
(14) 
Thus, Eq. (3) can be simplified as
(15) 
The boundary and initial conditions are discretized as
where .
4 Theoretical analysis of the finite difference method
4.1 Stability
Before giving the proof, we start with some useful lemmas.
Lemma 6 (b17 ).
Let be an order positive define matrix, then for any parameter , the following two inequalities
hold.
Now, we discuss the property of matrix .
Theorem 1.
Proof.
It is easy to obtain
where and . First, we consider the signs of . According to Lemma 4, when , , thus . According to Lemma 5, when , , thus . Therefore, when or . For the items and , we have
Since and , then
According to Lemmas 4 and Lemma 5, we have and , thus
as and . Now, for a given , we consider the sum
i.e.,
Thus, the proof is completed. ∎
According to the theorem, it is easy to conclude the following corollaries.
Corollary 1.
The matrix is strictly diagonally dominant as well. Therefore, is invertible and Eq. (15) is solvable.
Corollary 2.
The matrix is symmetric positive definite.
Proof.
In view of (14), the symmetry of is evident. Let
Comments
There are no comments yet.