Research Article  Open Access
Lijuan Su, Pei Cheng, "A HighAccuracy MOC/FD Method for Solving Fractional AdvectionDiffusion Equations", Journal of Applied Mathematics, vol. 2013, Article ID 648595, 8 pages, 2013. https://doi.org/10.1155/2013/648595
A HighAccuracy MOC/FD Method for Solving Fractional AdvectionDiffusion Equations
Abstract
Fractionalorder diffusion equations are viewed as generalizations of classical diffusion equations, treating superdiffusive flow processes. In this paper, in order to solve the fractional advectiondiffusion equation, the fractional characteristic finite difference method is presented, which is based on the method of characteristics (MOC) and fractional finite difference (FD) procedures. The stability, consistency, convergence, and error estimate of the method are obtained. An example is also given to illustrate the applicability of theoretical results.
1. Introduction
The history of fractional calculus is almost as long as integerorder calculus. However, because of lack of application background, fractional calculus was developed very slowly. Up to now, fractional calculus has been applied to almost every field of hydrology [1], physics [2, 3], engineering [4], mathematics [5], and science [6, 7], such as diffusion [8], oscillation, and viscoellastic dynamics. The fractional derivatives in space are often used to model anomalous diffusion because particles spread faster than the classical models predict. The main physical purpose for investigating diffusion equation of fractional order is to describe phenomena of anomalous diffusion in transport progresses through complex or disordered systems including fractal media, and fractional kinetic equations have been proved especially useful in the context of anomalous slow diffusion.
There are several different methods to solve differential equations of fractional order [9], for instance, the method of images, the Fourier transform methods, the method of separation of variables, the Mellin transform methods, and the Laplace transform methods. Wyss [10] considered the time fractional diffusion equation and the solution is given in closed form in terms of the Fox functions. In the last decade, different numerical methods have been developed to solve the fractional diffusion equations [11, 12]. Several explicit finite difference schemes were proposed and analyzed in [13]. Liu et al. [14] used L2 scheme form [5] to discretize the fractional FokkerPlanck equation. Deng [15] proposed finite element method for the space and time fractional FokkerPlanck equation. Meerschaert and Tadjeran [16] considered a GrünwaldLetnikov approximation for the twosided spacefractional partial differential equations.
Many difficult problems arise in the numerical simulation of fluid flow within porous media in petroleum reservoir simulation and in subsurface contaminant transport and remediation [17]. For example, their mathematical models are basically advection dominated. Fractional derivatives play an important role in modelling particle transport in anomalous diffusion [18]. The spacefractional advectiondominated diffusion equation [19] is where the function is the drift of the process, is the fractional order of the spatial derivative, physical considerations restrict , the functions may be interpreted as the coefficient of dispersion, and is a source/sink term. Moreover, when (1) is convection dominated diffusion equation, we have .
Equation (1) contains a RiemannLiouville fractional derivative of order of a function , for , , defined by [9] where is an integer such that , is the gamma function.
The standard finite difference methods (SFDM) for solving (1) were discussed in [16, 19, 20], and their stability was given in those papers, too. However, these methods often generate numerical solutions with severe nonphysical oscillations. Because of the nature of hyperbolic types [17, 21], the convectiondiffusion problems involve moving sharp fronts or boundary layers. While upwind finite difference methods eliminate nonphysical oscillations, they tend to generate numerical solutions with excessive numerical diffusion and smear the moving steep fronts. The efficient grid refinements are required if the standard methods are used in computations. But for high dimensional largescale problems, it may lead to very large linear systems and need long computation time.
In this paper, we will consider a new numerical difference scheme proposed by the characteristic method, which is named the fractional characteristic finite difference method (CFDM). It is proved that the new method is unconditionally stable, consistent, and convergent; therefore, there is no stability limitation on the size of . Moreover, the CFDM allows large time steps to be used without the loss of accuracy, which is much better than the SFDM and other numerical methods. The reason is that, for problems with significant convection, the solution changes much less rapidly in the characteristic direction than in the direction.
The outline of the paper is as follows. In Section 2, we describe the fractional characteristic finite difference schemes for solving (1). In Section 3, the stability, consistency, and convergence are proved, and the maximum error estimate is derived, too. The numerical experiments are given in Section 4. The conclusions of this paper are given in Section 5.
2. Fractional Characteristic Difference Method
In this section, we consider the combination of characteristic methods with finite difference techniques [22] and seek a new numerical method of (1), which is the convectiondominated problem and reflects the almost hyperbolic nature.
2.1. The Idea of Characteristic Methods
Let and denote the characteristic direction associated with the operator by (see [22]), where Then, the left part of (1) can be rewritten as
Next, in order to formulate the fractional CFDM of (1), we give some notations. Let be the temporal mesh or time step, , , , is the spatial mesh, the coordinates of the mesh points are , , , with and are positive integers. The values of the solution on these grid points are , and the numerical estimate of the exact value of at the point is denoted by . Finally, define , , and .
The characteristic derivative is approximated in the following manner: where .
Then, by (5) the differential operators at the point are discretized to be where , is calculated by the piecewiselinear interpolation.
Let , and let symbol be the integral part of . Let and , using linear interpolation to , we get [23] Generally, the quantity is called the Courant (or CFL) number.
2.2. The Discrete Versions of the Fractional Derivative
First, we develop the rightshifted Grünwald discrete approximation to the fractional derivatives by means of improving the standard form of the GrünwaldLetnikov definition [16]: where is a positive integer, , and , which can be calculated by
Following (5), we can discretize the RiemannLiouville operator [12] where is a positive integer, , and is the integral part of .
Here, according to the definition of , we can draw the following conclusion [9, 12]: when , ; when , .
According to the shifted GrünwaldLetnikov formulae (8) and (10), the definition (2) can be discretized by with .
Inserting (6), (7), and (11) into (1), neglecting the truncation error, after rearranging the terms, we finally get the fractional characteristic difference scheme where is associated with the diffusion coefficient, and the initial values are .
3. Stability and Convergence Analysis
In this section, before studying the stability and convergence of the CFDM, we give the following lemma, which was proved previously [16].
Lemma 1. According to the definition of in (9), when , the coefficients have the following properties:, and , , and ,, .
In the following section, we study the stability and convergence of the fractional CFDM for (1) with the given initial and boundary conditions.
3.1. Stability Analysis
Theorem 2. The fractional characteristic finite difference scheme (12) of (1) with , based on the shifted Grünwald approximation (8) to the space fractional derivative, is unconditionally stable.
Proof. Incorporating the boundary conditions, (12) is rewritten as , where , , and are two matrices.
When and , we have
where , so , too. Thus, holds.
For matrix , we can easily get
Therefore,
For matrix , we can easily see that for all , and
At the same time, according to Lemma 1, we have
for , we can obtain
Therefore, is diagonally dominant by rows. According to the matrix theory [24, 25],
Then, from (15) and (19), we can get
Thus, the characteristic finite difference scheme (12) is unconditionally stable.
In addition, taking into account (6) and (11), we have so this method is consistent with a local truncation error which is .
3.2. Convergence Analysis
Let ; be the exact solution of (1) at mesh point . Define , ; and . Using , substitution into (12) leads to where , and , ; .
From (6), (11), and (21), we have the following.
Lemma 3. There is a positive constant , such that where ; .
Let , because is finite, we can obtain the following result.
Theorem 4. Let be the numerical solution computed by the CFDM (12). Then there is a positive constant , such that
Proof. From (22) and , we can get
According to (19) and Lemma 3, using (25), we finally get
Because , from (26), we obtain the error estimate . Therefore, we get the result , where is a positive constant.
Finally, according to Theorem 4, we draw the following conclusion: when and , holds; that is, the characteristic finite difference scheme (12) is convergent. In summary, the fractional CFDM for (1) is unconditionally stable, consistent, and convergent.
4. Numerical Simulations
In this section, we give some numerical results for a special case of (1) to confirm our theoretical analysis discussed above. In addition, we also compare our new CFDM with the SFDM and point out the advantages of the new method.
For (1), if we assume that and are two constant functions, the function , and the initial condition is the dirac delta function. Then, under the zero Dirichlet boundary conditions, the exact solution of (1) could be found by the Fourier transform methods [9]:
In the numerical experiments, the data are chosen as follows: , , , and . For simplicity, we assume that the initial time , and the initial values are calculated from (27). The numerical solutions are obtained from the CFDM scheme (12) discussed above. The errors of numerical solutions are measured in and norms which are defined as follows: where and are the numerical solution and analytical solution at , respectively. is the space step, and is the time step. The ratios of convergence in time are calculated by where and are two time steps, and very small spatial step size is taken in computation. Similarly, the ratios of convergence in space can be calculated, too.
4.1. The Properties of the CFDM
In this subsection, we use the CFDM (12) to compute the numerical solutions to (1), and list the computational results in the figures and tables. From Figure 1, we can see that the numerical solutions of the CFDM coincide with the analytic solutions, with the diffusion coefficient , the space step , and the time step at , for and 1.6. It is clearly shown that the CFDM is stable for (1), which allows the large time step to be used and gains a very good numerical approximation to the exact solution.
Next, we examine the ratio of convergence in space for the CFDM. Table 1 shows and errors and ratios in space with the diffusion coefficient and the small time step at time , for two different fractional orders and 1.6, respectively. The second and fourth columns show the and errors, respectively. Column three shows the ratio of the CFDM in norm, and column five shows the ratio of the CFDM in norm. It is clearly shown that the CFDM method is of firstorder accuracy in space in both and norms.

Meanwhile, Table 2 shows and errors and ratios in time with and at time , for two different fractional orders and 1.6, respectively. The table shows that the CFDM method is of firstorder accuracy in time in both and norms, too. Furthermore, from the table, we can see that even with the very large time steps size of to 1/30, which corresponds to the Courant numbers of 14.14 to 84.85, the CFDM can still generate the accurate numerical solutions without noticeable numerical errors. In conclusion, according to Tables 1 and 2, the accuracy of the CFDM method is of order , which confirms the result of Theorem 4.

4.2. Comparison with the SFDM
To gain a better understanding of the CFDM, we use the welldeveloped and wellreceived SFDM to simulate the same model problem for comparison, which usually contains the implicit and explicit finite difference schemes. We only consider the case of in this subsection.
For the implicit SFDM, it is also unconditionally stable. Table 3 and Figure 2 both show that the CFDM is much better than the implicit SFDM, if they have the same assumptions. Firstly, Table 3 shows the comparison of errors calculated by the CFDM (12) and the implicit SFDM for two different diffusion coefficients and , with some different space steps and time steps for at time . From the table, we can easily see that the errors generated by the CFDM is obviously much smaller than the errors generated by the implicit SFDM, which has one even twoorder higher accuracy than the latter.

Moreover, Figure 2 shows that the errors between exact solutions and numerical solutions calculated by the CFDM and the SFDM at different time levels with , the diffusion coefficient , the space step , and the time step . It is easily to see that the CFDM can generate accurate numerical solution even if relatively large time step is used, and this significantly improves efficiency [17, 26]. However, the implicit SFDM does not have this kind of property.
Next, we compare the CFDM and the implicit SFDM numerical solutions with and at time in Figure 3. To keep similar error, for the CFDM, the grid sizes are chosen as and , and it generates the error of . However, for the implicit SFDM, we have to choose very small grid sizes, which are at least and . In this case, the implicit SFDM generates the error of , which is still larger than that of the CFDM. As a result, it is evidently shown that the CFDM can allow large time steps and space steps without the loss of accuracy.
In Section 4 discussed above, we only compare the numerical results of the CFDM with the implicit SFDM. In fact, the SFDM for (1) also contains several other explicit schemes which are conditionally stable. There are the upwind explicit scheme, central explicit scheme, and the LaxWendroff explicit scheme, and so forth. Of course, our new fractional CFDM is obviously better than all of them for the spacefractional convectiondominated diffusion problems because the explicit SFDM schemes are conditionally stable and they does not allow large time steps.
Finally, we compare the CFDM with the explicit SFDM schemes simply. Table 4 shows that the errors of three different explicit schemes with different and for and at time . As a result, the much refined time steps and space steps have to be used to improve the accuracy of the numerical solution in Table 4, in order to obtain a numerical solution with comparable accuracy to the numerical solution by the CFDM with the relatively coarse spatial grids and time steps in Tables 2 or 3. Moreover, under the same assumptions of Table 3, the three explicit SFDM schemes are all unstable.

5. Conclusions
In actual applications, we often use the spacefractional derivatives to model convectiondiffusion problems, and in many diffusion processes arising in physical problems, convection essentially dominates diffusion. In this paper, we propose a new scheme called the fractional CFDM and prove that it is unconditionally stable, consistent, and convergent. For this new method, even though the time steps are very coarse, the resulting CFDM solutions are more accurate than the SFDM solutions with much finer temporal grid sizes. Thus, the CFDM is more efficient and superior to the SFDM, especially for the convectiondominated problems.
Acknowledgments
This research was supported by the National Natural Science Foundations of China Grant (no. 11126179 and no. 11226247), the 211 project of Anhui University (no. 02303319 and no. KJTD002B), the Scientific Research Award for Excellent MiddleAged and Young Scientists of Shandong Province no. BS2010HZ012, and the Nature Science Foundation of Anhui Provincial (no. 1308085QA15). The authors acknowledge the anonymous reviewers for their helpful comments.
References
 D. Benson, S. W. Wheatcraft, and M. M. Meerschaert, “Application of a fractional advectiondispersion equation,” Water Resources Research, vol. 36, pp. 1403–1413, 2000. View at: Google Scholar
 A. Ashyralyev and Z. Cakir, “On the numerical solution of fractional parabolic partial differential equations with the Dirichlet condition,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 696179, 15 pages, 2012. View at: Google Scholar  Zentralblatt MATH  MathSciNet
 J. P. Roop, “Computational aspects of FEM approximation of fractional advection dispersion equations on bounded domains in ${R}^{2}$,” Journal of Computational and Applied Mathematics, vol. 193, no. 1, pp. 243–268, 2006. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 R. L. Magin, Fractional Calculus in Bioengineering, Begell House Publishers, 2006.
 K. B. Oldham and J. Spanier, The Fractional Calculus, Academic Press, New York, NY, USA, 1974. View at: MathSciNet
 H. Yang, “Existence of mild solutions for a class of fractional evolution equations with compact analytic semigroup,” Abstract and Applied Analysis, vol. 2012, Article ID 903518, 15 pages, 2012. View at: Publisher Site  Google Scholar
 M. De la Sen, “Positivity and stability of the solutions of Caputo fractional linear timeinvariant systems of any order with internal point delays,” Abstract and Applied Analysis, vol. 2011, Article ID 161246, 25 pages, 2011. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 V. J. Ervin, N. Heuer, and J. P. Roop, “Numerical approximation of a time dependent, nonlinear, spacefractional diffusion equation,” SIAM Journal on Numerical Analysis, vol. 45, no. 2, pp. 572–591, 2007. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 I. Podlubny, Fractional Differential Equations, Academic Press, New York, NY, USA, 1999. View at: MathSciNet
 W. Wyss, “The fractional diffusion equation,” Journal of Mathematical Physics, vol. 27, no. 11, pp. 2782–2785, 1986. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 R. Gorenflo and E. A. AbdelRehim, “Convergence of the GrünwaldLetnikov scheme for timefractional diffusion,” Journal of Computational and Applied Mathematics, vol. 205, no. 2, pp. 871–881, 2007. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 C. Lubich, “Discretized fractional calculus,” SIAM Journal on Mathematical Analysis, vol. 17, no. 3, pp. 704–719, 1986. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 E. Sousa, “Finite difference approximations for a fractional advection diffusion problem,” Journal of Computational Physics, vol. 228, no. 11, pp. 4038–4054, 2009. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 F. Liu, V. Anh, and I. Turner, “Numerical solution of the space fractional FokkerPlanck equation,” in Proceedings of the International Conference on Boundary and Interior Layers—Computational and Asymptotic Methods (BAIL '02), vol. 166, no. 1, pp. 209–219, 2004. View at: Publisher Site  Google Scholar  MathSciNet
 W. Deng, “Finite element method for the space and time fractional FokkerPlanck equation,” SIAM Journal on Numerical Analysis, vol. 47, no. 1, pp. 204–226, 2008. View at: Publisher Site  Google Scholar  MathSciNet
 M. M. Meerschaert and C. Tadjeran, “Finite difference approximations for fractional advectiondispersion flow equations,” Journal of Computational and Applied Mathematics, vol. 172, no. 1, pp. 65–77, 2004. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 H. Wang, R. E. Ewing, G. Qin, S. L. Lyons, M. AlLawatia, and S. Man, “A family of EulerianLagrangian localized adjoint methods for multidimensional advectionreaction equations,” Journal of Computational Physics, vol. 152, no. 1, pp. 120–163, 1999. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 P. Zhuang, F. Liu, V. Anh, and I. Turner, “New solution and analytical techniques of the implicit numerical method for the anomalous subdiffusion equation,” SIAM Journal on Numerical Analysis, vol. 46, no. 2, pp. 1079–1095, 2008. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 R. Schumer, D. A. Benson, M. M. Meerschaert, and B. Baeumer, “Multiscaling fractional advectiondispersion equations and their solutions,” Water Resources Research, vol. 39, pp. 1022–1032, 2003. View at: Google Scholar
 L. Su, W. Wang, and Z. Yang, “Finite difference approximations for the fractional advectiondiffusion equation,” Physics Letters A, vol. 373, no. 48, pp. 4405–4408, 2009. View at: Google Scholar
 H. Wang, H. K. Dahle, R. E. Ewing, M. S. Espedal, R. C. Sharpley, and S. Man, “An ELLAM scheme for advectiondiffusion equations in two dimensions,” SIAM Journal on Scientific Computing, vol. 20, no. 6, pp. 2160–2194, 1999. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 J. Douglas, and T. F. Russell, “Numerical methods for convectiondominated diffusion problems based on combining the method of characteristics with finite element or finite difference procedures,” SIAM Journal on Numerical Analysis, vol. 19, no. 5, pp. 871–885, 1982. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
 K. W. Morton and D. F. Mayers, Numerical Solution of Partial Differential Equations, Cambridge University Press, Cambridge, UK, 1994. View at: MathSciNet
 J. M. Varah, “A lower bound for the smallest singular value of a matrix,” Linear Algebra and Its Applications, vol. 11, pp. 3–5, 1975. View at: Google Scholar  Zentralblatt MATH  MathSciNet
 R. S. Varga, Matrix Iterative Analysis, Springer, Berlin, Germany, 2000. View at: Publisher Site  MathSciNet
 H. Wang, “An optimalorder error estimate for a family of ELLAMMFEM approximations to porous medium flow,” SIAM Journal on Numerical Analysis, vol. 46, no. 4, pp. 2133–2152, 2008. View at: Publisher Site  Google Scholar  Zentralblatt MATH  MathSciNet
Copyright
Copyright © 2013 Lijuan Su and Pei Cheng. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.