A Discrete Method Based on the CE-SE Formulation for the Fractional Advection-Dispersion Equation

Jerez, Silvia; Dzib, Ivan

doi:https://doi.org/10.1155/2015/303857

Discrete Dynamics in Nature and Society

On this page

Abstract Introduction Conclusions Appendix References Copyright Related Articles

Research Article | Open Access

Volume 2015 | Article ID 303857 | https://doi.org/10.1155/2015/303857

A Discrete Method Based on the CE-SE Formulation for the Fractional Advection-Dispersion Equation

Silvia Jerez¹and Ivan Dzib¹

Academic Editor: Manuel De la Sen

Received16 Nov 2014

Revised15 Jan 2015

Accepted17 Jan 2015

Published08 Feb 2015

Abstract

We obtain a numerical algorithm by using the space-time conservation element and solution element (CE-SE) method for the fractional advection-dispersion equation. The fractional derivative is defined by the Riemann-Liouville formula. We prove that the CE-SE approximation is conditionally stable under mild requirements. A numerical simulation is performed for the one-dimensional case by considering a benchmark with a discontinuous initial condition in order to compare the results with the analytical solution.

1. Introduction

Recently, the analytical and numerical framework of fractional differential equations has attracted much attention because of its extended use in the modeling of nonlocal transport phenomena which appear in many fields like magnetic plasma, turbulence, and flows in porous media [1–5]. A nonlocal transport can lead an anomalous diffusion because of large tracer displacements (super-diffusion) or trapping structure vortices (sub-diffusion). The standard local diffusion description is linked to Gaussian processes and diffusion differential models are obtained. However, a nonlocal transport flow presents a non-Gaussian particle motion and in this case fractional diffusion models are derived [6–8]. One of the most considered proposals for the modeling of the mass flow in media with anomalous diffusion is the fractional advection-dispersion equation (FADE); see, for example, [7–10].

Different definitions for the fractional derivative can be used to describe the FADE, and the most common ones are the Caputo, Riemann-Liouville, and Grünwald-Letnikov derivatives [11–13]. Here, we consider the one-dimensional FADE given by and use the Riemann-Liouville formula to describe the fractional derivatives: where and represent the space and the time, respectively, is the solute concentration, is the fractional dispersion, is the skewness of the diffusive transport taking values in the interval , and is the fractional order. In our case we consider and by definition of the ceil function we get . We can rewrite the fractional equation (1) in a more compact form as follows: by defining the fractional operator as

The analytical solution of (3) is only known under particular initial and boundary conditions, so numerical solutions are calculated and studied for most of the problems. In the last years, a lot of works have been focused on the development of numerical algorithms for fractional differential equations. The first efforts are based on extensions of standard methods used for the nonfractional case like the finite difference method (FDM) [14–17], the finite volume method (FVM) [18, 19], and the finite element method (FEM) [20, 21]. More recently other methods have been developed, for example, the spectral method [22], the collocation method [23], or a combination of both [24]. Generally, the FDM is constructed using the shifted Grünwald formula for the discretization of the fractional derivative, the FVM under a conservative formulation is based on the discretization of the integral using the Riemann-Liouville definition or the Caputo definition of the fractional derivative, and the FEM and collocation methods are variational techniques based on the discretization of the integral on some interpolating polynomial that approximates the solution. Therefore, all these techniques generate algorithms with dense coefficient matrices for which it is complicated to analyze their stability and a high computational cost to solve them is required.

As an alternative to avoid the discretization of the integral associated with the fractional derivative, in this paper we propose a numerical algorithm based on the space-time conservation element and solution element (CE-SE) method developed by Chang and To [25]. For this method, the solution is approximated by a first-order Taylor expansion which must fulfill the integral form of the equation in the conservation elements and the differential form in the solution elements. Therefore, the space-time domain is discretized two times by: a rhomboid mesh and a rectangular mesh. The theoretical framework of the CE-SE method for the nonfractional case is well established and it has been applied to a large number of real problems with successful results [26–33]. An important property of the fractional differential operators is the nonlocality which is achieved with the discretization of the integral associated with the fractional derivative. In our case we seek to keep the nonlocality by providing information of the numerical solution at nearby points in both meshes and using nonlocal approximations of the standard derivatives. It is important to remark that we consider this work as a first step to open a venue for a different way to approximate the fractional differential equations like the FADE equation.

The paper is organized as follows. In Section 2, we give the discretization of the space-time domain for the CE-SE method, its extension to the fractional advection-dispersion equation (4), and the von Neumann stability analysis of the proposed algorithm. In Section 3, we present numerical experiments of the fractional CE-SE method for the FADE with an initial condition type shock and it is compared with the analytical solution of the problem. Finally, this work is completed with a section of conclusions and an appendix with the calculation of the fractional Riemann-Liouville of a linear polynomial of two variables.

2. CE-SE Method for the 1D FADE

2.1. Discretization of the Space-Time Domain

An important feature of the CE-SE method is how the space-time domain is discretized. Let us start with the usual discretization by taking a set of points such that for and for each we take such that with being the left end of the spatial interval and the spatial and time step sizes are defined by and with being a constant of proportionality. Next, the space-time domain is partitioned by two types of mesh: nonoverlapping opened rhombi centered in the solution node which are named solution elements (SEs) and nonoverlapping closed rectangles which are called conservation element (CEs); see Figure 1.

By definition, is the interior of the space-time region bounded by the opened rhombus in Figure 2. In each , the flow variables are assumed to be continuous with possible discontinuities at the interface between the two rhombuses. The second mesh is given by nonoverlapping rectangular regions, , which fill the computational domain, see Figure 2. In each , a local space-time flux balance is enforced so a global flux conservation relation is satisfied over the entire space-time domain. Note that each consists of three different s and if we compute the solution in a integer time step then the spatial index is running in the half-integer points and vice versa, see Figure 1.

(a)

(b)

2.2. The Fractional CE-SE Method

Now, we construct an explicit algorithm based on the CE-SE method [25, 34] for the one-dimension FADE. First, the solution of (3) is approximated in each by a first-order Taylor expansion of two variables which is denoted by and is given by Denoting by , , and for all and integers, we need to know these coefficients to obtain the approximation (5). Assuming that they are known in the time step , we are going to calculate the coefficients and imposing that the numerical approximation (5) verifies the integral form of (3) in each . We split the rectangle in two subrectangles and solve the integral form in each subrectangle; see Figure 3. Consider for . Applying the Green theorem to the integral of the left hand side of (6) we have that where is the parameterization of given by:

(a)

(b)

Substituting the values of and for each and replacing and its partial derivatives by the coefficients ’s, ’s, and ’s in (7), we obtain Now, we rewrite the integral of the right hand side of (6) as a function of the coefficients ’s and ’s. For that, it is necessary to calculate the Riemann-Liouville fractional derivative of the function defined in (5). This calculation is carry out in the Appendix. Let be a spatial interval, and for any we write and ( the right-end index of ) and substituting them in (A.5), we get where Since the approximation defined in (5) must fulfill the FADE (1) in each , then the coefficient can be obtained by for . Denoting by the Courant number and by the anomalous diffusion number, we can rewrite the previous formula for the coefficient as Using (15) for and substituting (9) and (11) into (6) for the case of and for substituting (10) and (12) into (6), we have the following discrete system: The system (16) can be written in matrix form as follows: where , , and are the coefficient matrices that are given by

The CE-SE method is a family of schemes [35], and the developed algorithm (16) is called the a-scheme which is the core algorithm of the CE-SE family. Next, we are going to prove the convergence of this fractional a-scheme verifying its consistency with the linear FADE (1) and obtaining its stability conditions. For the stability analysis, we use the von Neumann technique [36] which is a local analysis. Thus, the boundary effects are not considered which is a reasonable assumption since in most cases the numerical instability appears in the interior nodes of the discretization of the domain.

2.3. Stability and Convergence Analysis

We start the von Neumann stability assuming that the solution and its partial derivative with respect to can be expressed as a sum of eigenmodes at each node of the mesh: where is the spatial wave number and is a complex number. Now, we substitute (19) in the discrete system (17) and obtain and multiplying both sides by yields with Taking the square of Euclidean norm in both sides of (21), we have where is the amplification factor of each mode in the expansions (19). By calculating we have where , , , and . The stability of the method is achieved if the magnitude of the amplification factor is less than 1. By definition, the values of the positive coefficients and are close; thus the following condition is fulfilled By definition of the coefficient , we get ; thus if we take small enough.

Now, we analyze the local truncation error to obtain the consistency of the fractional CE-SE method (17) with the FADE (1). Given that we assume the flow variables , , and to be continuous in each and considering that the numerical solution is a first-order Taylor approximation (5), we conclude that the fractional a-scheme (17) is second-order accurate in space and time. Given the above results, we have the following theorem using the Lax Equivalence theorem [36] which states that For a uniformly solvable linear finite difference scheme which is consistent with a well-posed linear evolution problem, the stability is a necessary and sufficient condition for its convergence.

Theorem 1. The fractional CE-SE scheme (17) is convergent for a small enough.

It is well known that the a-scheme of the CE-SE family is a nondissipative algorithm. We are interested in flow problems with discontinuities, so we need to introduce some numerical dissipation in the approximation. In the next section, we modify the a-scheme using another way to approximate the derivative , and the modified algorithm is known as a--scheme.

2.4. An Approximation of the Derivative for Discontinuities

In order to introduce information about the behavior of the solution in the calculation of the partial derivative of with respect to , we propose another approximation for the coefficient . By adding both equations of the system (16), we can rewrite the fractional CE-SE algorithm as follows: As a first approximation, we consider the coefficient given in [25, 28]: where the ratio , with a positive constant, can be considered as a sloper limiter coefficient which measures the smoothness of the function by using the forward difference and the backward difference: with . As the points are not solution nodes, we can obtain as follows: The estimation (27) for the partial spatial derivative has been used in nonfractional differential problems with nonsmooth solutions obtaining successful results [25, 26, 28]. In the boundary nodes, the coefficient is calculated from the forward difference and the coefficient is obtained with the backward difference .

Remark 2. Notice that we have used a weighted finite difference approximation to approximate the coefficient in the step , but other approximations can be used. We think that this can be a way to introduce nonlocal information in the calculation of the numerical solution.

3. Numerical Results

For the fractional case it is important to remark that there are few benchmark problems with known analytical solution. So in this section, we consider (3) with a jump condition of the shock type: with a positive constant function and analyze the performance of the fractional CE-SE method (26)-(27) by comparing with the exact solution, given by Sousa [37]. So, taking and the boundary conditions and , the exact solution is where and the cumulative probability function is defined by with For all and , the cumulative probability function is obtained using the identity . Next, we show the numerical results of the fractional a--scheme for the FADE (1) with initial condition (30) using the same test problem given in [37].

Example 3. We consider the solution of (1) with (30) in the interval for , , and and the boudary conditions are and . In Figures 4 and 5, we show the profile of the analytical solution (31) and the CE-SE approximation (26)-(27) at for several values of the parameter . Moreover, we compare the global error of the fractional a- scheme versus some finite difference methods that appear in [37] for a fixed time using the norm given by where is the difference between the exact and the numerical solutions; see Table 1. Notice that for the three schemes as is closer to one the numerical error is larger than in the cases as is closer to two, since for smaller the solution loses differentiability. This fact is visible at the extremes of the shock in Figure 4, where we display the CE-SE behavior for and 1.4. Although the three methods present similar numerical errors, the CE-SE method has the advantage of holding a low computational cost compared to the other two methods. Finally in Figure 5, we show the accuracy of the fractional CE-SE method as the step size is diminished.

(a)

(b)

(a)

(b)

4. Conclusions

In this work we have constructed a discrete method under the CE-SE formulation for the fractional advection-diffusion equation described by the Riemann-Liouville derivative. The new fractional method avoids the fact that the calculation at the new point depends on the function values at points of further regions at that point. We have also proved the convergence of the fractional CE-SE method and show its performance for a benchmark problem with a shock initial condition. Future work will focus on the analysis of different proposals to approximate the partial spatial derivative of the solution to include the nonlocal information of the flow. We also consider extending the fractional CE-SE method to two dimensions in order to cover a wider range of flow problems important.

Appendix

Fractional Derivative of a Two-Variable Linear Polynomial

Let be a first-degree polynomial function of two variables as follows: where and are the independent variables with as the indices of the -discretization and -discretization, respectively, and are constant coefficients. Using the linearity property, we rewrite the Riemann-Liouville fractional derivative (4) of the polynomial function (A.1) as follows: for . Next, we calculate the fractional derivatives that appear in (A.2) using the formulas (2) and obtain Substituting (A.3) in (A.2), we have Later on algebraic manipulations, the Riemann-Liouville fractional derivativative of the function defined in (A.1) is given by

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

References

P. J. Torvik and R. L. Bagley, “On the appearance of the fractional derivative in the behavior of real materials,” Journal of Applied Mechanics, vol. 51, no. 2, pp. 294–298, 1984.
View at: Publisher Site | Google Scholar
T. R. Ellsworth, P. J. Shouse, T. H. Skaggs, J. A. Jobes, and J. Fargerlund, “Solute transport in unsaturated soil: experimental design, parameter estimation, and model discrimination,” Soil Science Society of America Journal, vol. 60, no. 2, pp. 397–407, 1996.
View at: Publisher Site | Google Scholar
Y. Pachepsky, D. Benson, and W. Rawls, “Simulating scale-dependent solute transport in soils with the fractional advective-dispersive equation,” Soil Science Society of America Journal, vol. 64, no. 4, pp. 1234–1243, 2000.
View at: Publisher Site | Google Scholar
D. del-Castillo-Negrete, B. A. Carreras, and V. E. Lynch, “Fractional diffusion in plasma turbulence,” Physics of Plasmas, vol. 11, no. 8, pp. 3854–3864, 2004.
View at: Publisher Site | Google Scholar
B. M. Schulz and M. Schulz, “Numerical investigations of anomalous diffusion effects in glasses,” Journal of Non-Crystalline Solids, vol. 352, no. 42–49, pp. 4884–4887, 2006.
View at: Publisher Site | Google Scholar
A. A. Greenenko, A. V. Chechkin, and N. F. Shul’ga, “Anomalous diffusion and Lévy flights in channeling,” Physics Letters, A, vol. 324, no. 1, pp. 82–85, 2004.
View at: Publisher Site | Google Scholar | MathSciNet
D. del-Castillo-Negrete, “Fractional diffusion models of nonlocal transport,” Physics of Plasmas, vol. 13, no. 8, Article ID 082308, 2006.
View at: Publisher Site | Google Scholar | MathSciNet
W. Chen, H. Sun, X. Zhang, and D. Korošak, “Anomalous diffusion modeling by fractal and fractional derivatives,” Computers & Mathematics with Applications, vol. 59, no. 5, pp. 1754–1758, 2010.
View at: Publisher Site | Google Scholar
D. A. Benson, S. W. Wheatcraft, and M. M. Meerschaert, “The fractional-order governing equation of Lévy motion,” Water Resources Research, vol. 36, no. 6, pp. 1413–1423, 2000.
View at: Publisher Site | Google Scholar
Y. Zhang, D. A. Benson, M. M. Meerschaert, and E. M. LaBolle, “Space-fractional advection-dispersion equations with variable parameters: diverse formulas, numerical solutions, and application to the Macrodispersion Experiment site data,” Water Resources Research, vol. 43, no. 5, Article ID W05439, 2007.
View at: Publisher Site | Google Scholar
M. Weilbeer, Efficient Numerical Methods for Fractional Differential Equations and Their Analytical Background, Papierflieger, 2005.
D. Baleanu, K. Diethelm, E. Scalas, and J. J. Trujillo, Fractional Calculus Models Andnumerical Methods, vol. 3, World Scientific, Singapore, 2009.
View at: Publisher Site | MathSciNet
C. Li, D. Qian, and Y. Chen, “On Riemann-Liouville and Caputo derivatives,” Discrete Dynamics in Nature and Society, vol. 2011, Article ID 562494, 15 pages, 2011.
View at: Publisher Site | Google Scholar | MathSciNet
M. M. Meerschaert and C. Tadjeran, “Finite difference approximations for fractional advection-dispersion flow equations,” Journal of Computational and Applied Mathematics, vol. 172, no. 1, pp. 65–77, 2004.
View at: Publisher Site | Google Scholar | MathSciNet
C. Li and F. Zeng, “Finite difference methods for fractional differential equations,” International Journal of Bifurcation and Chaos in Applied Sciences and Engineering, vol. 22, no. 4, Article ID 1230014, 2012.
View at: Publisher Site | Google Scholar | MathSciNet
F. Liu, P. Zhuang, V. Anh, I. Turner, and K. Burrage, “Stability and convergence of the difference methods for the space-time fractional advection-diffusion equation,” Applied Mathematics and Computation, vol. 191, no. 1, pp. 12–20, 2007.
View at: Publisher Site | Google Scholar | MathSciNet
A. Ashyralyev and F. Dal, “Finite difference and iteration methods for fractional hyperbolic partial differential equations with the Neumann condition,” Discrete Dynamics in Nature and Society, vol. 2012, Article ID 434976, 15 pages, 2012.
View at: Publisher Site | Google Scholar | MathSciNet
X. Zhang, J. W. Crawford, L. K. Deeks, M. I. Stutter, A. G. Bengough, and I. M. Young, “A mass balance based numerical method for the fractional advection-dispersion equation: theory and application,” Water Resources Research, vol. 41, no. 7, Article ID W07029, 2005.
View at: Publisher Site | Google Scholar
H. Hejazi, T. Moroney, and F. Liu, “Stability and convergence of a finite volume method for the space fractional advection-dispersion equation,” Journal of Computational and Applied Mathematics, vol. 255, pp. 684–697, 2014.
View at: Publisher Site | Google Scholar | MathSciNet
Q. Huang, G. Huang, and H. Zhan, “A finite element solution for the fractional advection-dispersion equation,” Advances in Water Resources, vol. 31, no. 12, pp. 1578–1589, 2008.
View at: Publisher Site | Google Scholar
Y. Zheng, C. Li, and Z. Zhao, “A note on the finite element method for the space-fractional advection diffusion equation,” Computers & Mathematics with Applications, vol. 59, no. 5, pp. 1718–1726, 2010.
View at: Publisher Site | Google Scholar
X. Li and C. Xu, “A space-time spectral method for the time fractional diffusion equation,” SIAM Journal on Numerical Analysis, vol. 47, no. 3, pp. 2108–2131, 2009.
View at: Publisher Site | Google Scholar | MathSciNet
A. Saadatmandi, M. Dehghan, and M.-R. Azizi, “The Sinc-Legendre collocation method for a class of fractional convection-diffusion equations with variable coefficients,” Communications in Nonlinear Science and Numerical Simulation, vol. 17, no. 11, pp. 4125–4136, 2012.
View at: Publisher Site | Google Scholar | MathSciNet
M. Zayernouri and G. E. Karniadakis, “Fractional spectral collocation method,” SIAM Journal on Scientific Computing, vol. 36, no. 1, pp. A40–A62, 2014.
View at: Publisher Site | Google Scholar | MathSciNet
S. C. Chang and W. M. To, “A new numerical framework for solving conservation laws—the method of space-time conservation element and solution element,” Tech. Rep. 104495, NASA, Cleveland, Ohio, USA, 1991.
View at: Google Scholar
X. Y. Wang, S. C. Chang, and P. C. E. Jorgenson, “Prediction of sound waves propagatingthrough a nozzle withoutlwilh a shock wave using the space-time CE/SE method,” NASA 209937, 2000.
View at: Google Scholar
Z.-C. Zhang, S. T. J. Yu, and S.-C. Chang, “A space-time conservation element and solution element method for solving the two- and three-dimensional unsteady Euler equations using quadrilateral and hexahedral meshes,” Journal of Computational Physics, vol. 175, no. 1, pp. 168–199, 2002.
View at: Publisher Site | Google Scholar
S. Jerez, J. V. Romero, M. D. Rosell{\'o}, and F. J. Arnau, “A semi-implicit space-time CESE method to improve mass conservation through tapered ducts in internal combustion engines,” Mathematical and Computer Modelling, vol. 40, no. 9-10, pp. 941–951, 2004.
View at: Publisher Site | Google Scholar | MathSciNet
J. Wang, K. Liu, and D. Zhang, “An improved CE/SE scheme for multi-material elastic-plastic flows and its applications,” Computers & Fluids, vol. 38, no. 3, pp. 544–551, 2009.
View at: Publisher Site | Google Scholar
S. Qamar and S. Mudasser, “On the application of a variant CE/SE method for solving two-dimensional ideal MHD equations,” Applied Numerical Mathematics, vol. 60, no. 6, pp. 587–606, 2010.
View at: Publisher Site | Google Scholar | MathSciNet
S. C. Chang, A New Approach for Constructing Highly Stable High Order CESE Schemes, National Aeronautics and Space Administration, Glenn Research Center, 2010.
S. Sadighi, A. Ahmad, and M. Shirvani, “Dynamic simulation of a pilot scale vacuum gas oilhydrocracking unit by the spacetime CE/SE method,” Chemical Engineering and Technology, vol. 35, no. 5, pp. 919–928, 2012.
View at: Publisher Site | Google Scholar
S. Qamar and M. Yousaf, “The space-time CESE method for solving special relativistic hydrodynamic equations,” Journal of Computational Physics, vol. 231, no. 10, pp. 3928–3945, 2012.
View at: Publisher Site | Google Scholar | MathSciNet
D. Yang, S. Yu, and J. Zhao, “Convergence and error bound analysis for the space-time CESE method,” Numerical Methods for Partial Differential Equations. An International Journal, vol. 17, no. 1, pp. 64–78, 2001.
View at: Google Scholar | MathSciNet
S.-C. Chang, “Courant number insensitive CE/SE schemes,” in Proceedings of the 38th AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit, Indianapolis, Ind, USA, July 2002.
View at: Google Scholar
L. N. Trefethen, Finite Difference and Spectral Methods for Ordinary and Partial Differential Equations, Cornell University, Ithaca, NY, USA, 1996.
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 | MathSciNet

Copyright

Copyright © 2015 Silvia Jerez and Ivan Dzib. 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.

PDF Download Citation

Download other formats

Order printed copies

Views

900

Downloads

839

Citations