Mathematical Problems in Engineering

Volume 2013 (2013), Article ID 246025, 11 pages

http://dx.doi.org/10.1155/2013/246025

Research Article

New High-Order Compact ADI Algorithms for 3D Nonlinear Time-Fractional Convection-Diffusion Equation

Shuying Zhai,¹ Xinlong Feng,¹ and Zhifeng Weng²

¹College of Mathematics and Systems Science, Xinjiang University, Urumqi 830046, China

²School of Mathematics and Statistics, Wuhan University, Wuhan 430072, China

Received 14 January 2013; Revised 8 April 2013; Accepted 15 April 2013

Academic Editor: Jun Jiang

Copyright © 2013 Shuying Zhai et al. 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.

Abstract.

Numerical approximations of the three-dimensional (3D) nonlinear time-fractional convection-diffusion equation is studied, which is firstly transformed to a time-fractional diffusion equation and then is solved by linearization method combined with alternating direction implicit (ADI) method. By using fourth-order Padé approximation for spatial derivatives and classical backward differentiation method for time derivative, two new high-order compact ADI algorithms with orders and are presented. The resulting schemes in each ADI solution step corresponding to a tridiagonal matrix equation can be solved by the Thomas algorithm which makes the computation cost effective. Numerical experiments are shown to demonstrate the high accuracy and robustness of two new schemes.

1. Introduction

In this paper, we consider numerical approximations of the following 3D nonlinear time-fractional convection-diffusion equation:

with the initial and boundary conditions

where , , , , and are constants, and is some reasonable nonlinear function of . is a continuous convex domain in 3D space with suitable boundary conditions prescribed on its boundary . The solution and the source terms and are assumed to be sufficiently smooth and have the necessary continuous partial derivatives up to certain orders. The fractional derivative is Grünward-Letnikov fractional derivative defined by

It should be noted that Grünward-Letnikov fractional derivative will become Caputo fractional derivative when .

Fractional differential equations (FDEs) have become popular in recent years both in mathematics and in applications. They have been increasingly used in physical and chemical processes [1–3], relaxation in polymer systems, dynamics of protein molecules, and the diffusion of contaminants in complex geological formations [4]. Because of the absence or appropriate mathematical methods, most of the analytical solutions for FDEs are difficult to obtain. However, several numerical methods have been proposed to solve them.

For the one-dimensional (1D) time FDEs, the schemes in [5–10] are convergent with order not higher than two for the space variable, and it is interesting to seek high-order numerical methods for FDEs. Recently, such problem was solved by the compact finite difference scheme with convergence order in [11], and a higher order one can be found in [12, 13]. Lin and Xu [14], Jiang and Ma [15], and Wei et al. [16–18] proposed high order method based on the finite spectral method, the finite element method (FEM), and local discontinuous Galerkin (LDG) method for space discretization, respectively.

For the two-dimensional (2D) FDEs, Zhuang and Liu [19] considered a time-fractional diffusion equation on a finite domain, and a first-order implicit difference approximation was proposed to solve the equation. Additional references [20–22] treat explicit or/and implicit finite difference methods (FDMs) for 2D time FDEs, as well as FEM [23], LDG method [24], and spectral method [25] for such problems. Recently, Zhang and Sun [26] presented new numerical techniques for time FDEs and ADI schemes were considered, but both methods have only second-order accuracy in spatial direction. Then, Cui [27] proposed a high-order ADI compact finite difference scheme with order for time-fractional diffusion equation. References [28–34] address FDM and FEM for space FDEs, respectively.

In this paper, we focus on the 3D nonlinear time-fractional convection-diffusion equation. To the authors’ knowledge, the research for the 3D case is very lacking. The reason, of necessity, is that 3D problems present many computational challenges because of the prohibitive memory and time requirement often necessary to solve them on grids of sufficient resolution. Traditional numerical schemes have low accuracy and thus require fine discretization. Moreover, being opposite to differential equations of integer order, in which derivatives depend only on the local behaviour of the function, fractional differential equations accumulate the whole information of the function in a weighted form. This makes realistic 3D simulations of FDEs prohibitively expensive. Therefore, it is interesting to discuss high-order numerical methods for 3D time-fractional convection-diffusion equation.

It is well known that Padé approximation has the advantages of the fourth-order accuracy to approximate the second-order derivatives and of keeping the desirable tridiagonal nature of the finite difference equations. So, the main purpose of this paper is to solve the nonlinear time-fractional convection-diffusion problem using the Padé approximation combined with ADI method and the one-dimensional tridiagonal Thomas algorithm. However, the convection terms will impede implementation of tridiagonal Thomas algorithm. So, we propose a transformation based on the compact finite difference method [35]. In this method, the original time-fractional convection-diffusion equation is transformed to a time-fractional diffusion equation; the resulting time-fractional diffusion problem is then solved by an efficient high-order method. Another advantage of our method is that the nonlinear term is transformed into a linear term which makes the computation cost effective and the approximation solutions more accurate.

This paper is organized as follows. In Section 2, two new high-order ADI algorithms are presented for solving 3D nonlinear time-fractional convection-diffusion equation. In Section 3, numerical examples are carried out to verify the high efficiency of our method. Finally, conclusions are drawn in Section 4.

2. Two New High-Order Compact ADI Algorithms

For simplicity of the presentation, we consider a rectangular domain , which is divided into an mesh with the spatial step size , while the time step size is chosen as and is the number of time coordinate direction.

To obtain an efficient scheme, it is important to approximate the nonlinear term in (1) explicitly. Otherwise, one would need to perform a nonlinear iteration in each time step which is quite time consuming. Let and . The following linear method is used to transform into linear term:

Next we discuss the linear equation

for .

In order to keep the fourth-order accuracy in space and the tridiagonal nature of the scheme by utilizing the Padé approximation, we first introduce a transformation that is similar to that of Liao [35] to eliminate the convection terms in (1).

Let

where and , and are functions to be determined. Substituting (6) into the first formula of (1), we have

Setting the coefficients of , and to zero leads to

Then, we have

Set

Combining (1) and (5)–(10), we obtain the following time-fractional diffusion equation satisfied by :

In the following, we will introduce a high-order compact difference method for (11).

The time-fractional derivative at is estimated by

where and the truncation error has the following form:

If we expend in Taylor series about at in (13) and use the analysis in [14], we will get

where is a positive constant only depending on .

As for the spatial derivatives, we can use the fourth-order Padé scheme. That is, the spatial derivatives can be approximated by

where , , and are the second-order central difference operators along - and - and -directions, respectively.

Combining with (12) and (15), then the implicit compact finite difference scheme for (11) is given as follows:

with the initial discretization

and the boundary discretization

where denotes . In the following, will be written in short for if there is no confusion about the notation.

By reformulating (16), we obtain an equivalent form

where

Due to the nature of three-time-level scheme, we need an additional method to compute the numerical solution at the first time step.

For the special case , the scheme simply reads

with

It is clear that (21) must be solved using some relaxation procedure because is nonlinear.

For simplicity of the presentation, we write (19) and (21) in a unified form

Multiplying on both sides of (23), we have

In the following, two new high-order algorithms are proposed for (24) that will coincide with the ADI methods employed. First, some finite difference operations are introduced that will be used hereafter as follows:

Algorithm 1. (1) For , add splitting term

to the left-hand side of (24), and rewrite it into an equivalent form

(2) Introducing two intermediate variables and , scheme can be solved in three steps as

Note that the intermediate values of and on the boundary are easily obtained from (28c) and (28b).

Algorithm 2. (1) For the first time step, can be obtained from Algorithm 1.

(2) For , add splitting term

to the left-hand side of (24), and rewrite it into an equivalent form

(3) Use the same way as the second step of Algorithm 1.

Remark 3. It is clear that the splitting error of II is higher order in than the local error for the time-fractional derivative approximation, so, for the time accuracy, scheme SII can ensure the local truncation error to be order. However, as for the splitting term I, it should be in order to maintain the time accuracy. That is, the order of local error for scheme SI is determined by the following equation:

Based on the previous analysis, we propose two linear ADI schemes, and the truncation errors are and , respectively. The resulting ADI schemes in each ADI solution step corresponding to a strictly diagonally dominant matrix equation which can be solved using the one-dimensional tridiagonal Thomas algorithm with a considerable saving in computing time.

Remark 4. The previous analysis shows that splitting term II is better than I when . However, when , the splitting error is essentially negligible if the time step size is small enough, and so both schemes can obtain the same accuracy in this case. Nevertheless, numerical experiments in Section 3 show that scheme SII is much better than schem SI, regardless of the value of . We would like to comment the fact that if the splitting error is much larger than the truncation error, as pointed out by Douglas and Kim [36], a good splitting term will play an important role in the accuracy of the solution, especially for the analytical solutions that decrease exponentially in time. So, we think that scheme SII is a better choice even with an additional computational cost than scheme SI.

3. Numerical Experiments

In this section, two numerical examples are presented to demonstrate the high efficiency and accuracy of the new method. The numerical results got from Algorithm 2 will be compared with Algorithm 1. The -norm error (denoted as -error) and the -norm error (denoted as -error) of the numerical solutions are shown. Note that the treatment of space is the same for two schemes all of them can achieve the same results when is small enough to eliminate the influence of the temporal approximation. So, in the following experiments, we only choose Algorithm 2 to investigate the spatial convergence rate.

Moreover, we compute the convergent order of accuracy for both schemes. Let us consider two mesh sizes and on and , respectively. The -norm errors of these two grids are denoted as and . If we set the convergent order of accuracy to be Rate, then we have the following form:

So, the convergent order of accuracy Rate can be computed as

The order of accuracy is formally defined when the mesh size approaches to zero. Therefore, when the mesh size is relatively large, the discretization scheme may not achieve its formal order of accuracy.

Problem 5. We firstly consider a pure diffusion problem. The exact analytical solution and the corresponding forcing term and initial condition are given by

Since it is difficult to accurately solve based on except for , in this problem, we only discuss the case with and the corresponding

We first investigate the time accuracy with for both algorithms. To this end, is chosen such that the errors stemming from the spatial approximation are negligible. Numerical results with are presented in Table 1. It is clearly show that Algorithm 2 is much better than Algorithm 1. For instance, in Table 1, the maximum absolute error of the computed solution by Algorithm 2 is around when . Whereas for Algorithm 1 with the same mesh size, the maximum absolute error is around . Similar comparison can be made to reach similar conclusions.

Meanwhile, we can see that the convergence rate for Algorithm 1 is around 1.5, which agrees with the theoretical estimates min(, ). Howevre, for Algorithm 2, the -rate is approximately 2.6, which is superconvergent. For better comparison of two algorithms, we plot the errors (as shown in Figure 1) at and with and for . Comparing both subplots, we can see that Algorithm 2 has substantially higher accuracy than Algorithm 1.

For the spatial convergence rate, we choose and to eliminate the influence of the temporal approximation. Numerical results using Algorithm 2 are presented in Table 2. The convergence order is indeed 4.

Figure 2 shows the contour plot of the exact solution in the - plane and numerical solutions using Algorithm 2 for different at and with and . We can see that when , the numerical solution is a good approximation of the exact solution.

Table 1: , and with different for Problem 5.

10
20
30
40
50

Algorithm 1
-error

-rate
—

-error

-rate
—

Algorithm 2
-error

-rate
—

-error

-rate
—

Table 2: ,  , and with different using Algorithm 2 for Problem 5.

5
10
20
30
40

-error

-rate
—

-error

-rate
—

(a) The -error using Algorithm 1

(b) The -error using Algorithm 2

Figure 1: The plot of the -error using both algorithms at and with , and for Problem 5.

(a)

(b)

(c)

(d)

Figure 2: The contour plot of the exact solution and numerical solutions for different at and   with and for Problem 5. (a) Exact solution; (b) ; (c) ; (d) .

Problem 6. In order to illustrate the effectiveness of our method, the following time-fractional convection-diffusion problem is considered:

The right-hand side function is specified to satisfy the given equation and the exact solution.

Numerical results for time with are presented in Tables 3–5, respectively. Again, as expected, it is shown that Algorithm 2 has higher accuracy in comparison with Algorithm 1. For instance, in Table 5, the maximum absolute error of the computed solution by Algorithm 2 is around when . Whereas for Algorithm 1 with the same mesh size, the maximum absolute error is around .

However, we note that the convergence orders for Algorithm 1 with and 0.5 are unsatisfying, as shown in Tables 3 and 4, respectively. We would like to comment on the fact that the order of accuracy is formally defined when the mesh size approaches to zero. Therefore, when the mesh size is relatively large, the discretization scheme may not achieve its formal order of accuracy. Moreover, the data show that the order is increasing now; hence, it will achieve its formal order with the decrease of the time step .

Figure 3 shows the errors at and with and for different . Again, as expected, it indicates that Algorithm 2 has substantially higher accuracy than Algorithm 1 regardless of .

Table 6 presents space results with and for different and , from which we can find that Algorithm 2 is stable and convergent and it can obtain fourth-order accuracy when mesh becomes finer, regardless of the value of .

Figure 4 shows the contour plot of the exact solution in the - plane and numerical solutions using Algorithm 2 for different at and with and . We can see that when , the numerical solution is a good approximation of the exact solution.

Table 3: ,   , and with different for Problem 6.

10
20
30
40
50

Algorithm 1
-error

-rate
—

-error

-rate
—

Algorithm 2
-error

-rate
—

-error

-rate
—

Table 4: , and with different for Problem 6.

10
20
30
40
50

Algorithm 1
-error

-rate
—

-error

-rate
—

Algorithm 2
-error

-rate
—

-error

-rate
—

Table 5: ,  , and with different for Problem 6.

10
20
30
40
50

Algorithm 1
-error

-rate

-error

-rate

Algorithm 2
-error

-rate

-error

-rate

Table 6: and with different and using Algorithm 2 for Problem 6.

5
10
20
30
40

-error

-rate
—

-error

-rate
—

-error

-rate
—

-error

-rate
—

-error

-rate
—

-error

-rate
—

(a) Algorithm 1,

(b) Algorithm 1,

(c) Algorithm 1,

(d) Algorithm 2,

(e) Algorithm 2,

(f) Algorithm 2,

Figure 3: The plot of the -error using both algorithms at and   with , and different for Problem 6.

(a)

(b)

(c)

(d)

Figure 4: The contour plot of the exact solution and numerical solutions for different at and with and   for Problem 6. (a) Exact solution; (b) ; (c) ; (d) .

4. Conclusions

In this paper, we have proposed two high-order compact ADI algorithms with order and , respectively, for solving the 3D nonlinear time-fractional convection-diffusion equation. Since it only involves three-point stencil for each one-dimensional operator, it can be solved by applying the one-dimensional tridiagonal Thomas algorithm with a considerable saving in computing time. Numerical experiments have shown that Algorithm 2 is very significantly reduced in the splitting perturbation error and is superconvergent, which is much better than Algorithm 1.

Compared with existing deducing methods, our approach is advantageous because

(1)the original time-fractional convection-diffusion equation is first transformed to a time-fractional diffusion equation, and the resulting time-fractional diffusion problem is then solved by ADI method, this is novel and much more simpler than other traditional methods;(2)the nonlinear term is transformed into a linear term which makes the computation cost effective and the approximation solutions more accurate;(3)we show that the critical modification in the ADI procedures can virtually eliminate the perturbation errors and produce identical approximations of the solution of the differential problem.

Acknowledgments

The authors would like to thank the editor and referees for their criticism, valuable comments, and suggestions which helped to improve the results of this paper. This work is in part supported by the NSF of China (nos. 11271313 and 61163027), the China Postdoctoral Science Foundation (nos. 201104702 and 2012M512056), the Key Project of Chinese Ministry of Education (no. 212197), the NSF of Xinjiang Province (no. 2013211B01), and the Excellent Doctor Innovation Program of Xinjiang University (no. XJUBSCX-2012003).

References

R. Gorenflo and F. Mainardi, “Fractional calculus: integral and differential equations of fractional order,” in Fractals and Fractional Caluculus in Continuun Mechanics, A. Carpinteri and F. Mainardi, Eds., vol. 378, pp. 223–276, Springer, Vienna, Austria, 1997.

F. Mainardi, “Fractional relaxation-oscillation and fractional diffusion-wave phenomena,” Chaos, Solitons and Fractals, vol. 7, no. 9, pp. 1461–1477, 1996.

I. Podlubny, Fractional Differential Equations, Academic Press, New York, NY, USA, 1999.

I. M. Sokolov, J. Klafter, and A. Blumen, “Fractional kinetics,” Physics Today, vol. 55, no. 11, pp. 48–54, 2002.

C.-M. Chen, F. Liu, I. Turner, and V. Anh, “A Fourier method for the fractional diffusion equation describing sub-diffusion,” Journal of Computational Physics, vol. 227, no. 2, pp. 886–897, 2007.

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.

C.-M. Chen, F. Liu, and K. Burrage, “Finite difference methods and a Fourier analysis for the fractional reaction-subdiffusion equation,” Applied Mathematics and Computation, vol. 198, no. 2, pp. 754–769, 2008.

R. Scherer, S. L. Kalla, L. Boyadjiev, and B. Al-Saqabi, “Numerical treatment of fractional heat equations,” Applied Numerical Mathematics, vol. 58, no. 8, pp. 1212–1223, 2008.

F. Liu, C. Yang, and K. Burrage, “Numerical method and analytic technique of the modified anomalous subdiffusion equation with a nonlinear source term,” Journal of Computational and Applied Mathematics, vol. 231, no. 1, pp. 160–176, 2009.

X. Zhao and Z.-z. Sun, “A box-type scheme for fractional sub-diffusion equation with Neumann boundary conditions,” Journal of Computational Physics, vol. 230, no. 15, pp. 6061–6074, 2011.

M. Cui, “Compact finite difference method for the fractional diffusion equation,” Journal of Computational Physics, vol. 228, no. 20, pp. 7792–7804, 2009.

G.-H. Gao and Z.-Z. Sun, “A compact finite difference scheme for the fractional sub-diffusion equations,” Journal of Computational Physics, vol. 230, no. 3, pp. 586–595, 2011.

J. C. Ren, Z.-Z. Sun, and X. Zhao, “Compact difference scheme for the fractional sub-diffusion equation with Neumann boundary conditions,” Journal of Computational Physics, vol. 232, pp. 456–467, 2013.

Y. M. Lin and C. J. Xu, “Finite difference/spectral approximations for the time-fractional diffusion equation,” Journal of Computational Physics, vol. 225, no. 2, pp. 1533–1552, 2007.

Y. J. Jiang and J. T. Ma, “High-order finite element methods for time-fractional partial differential equations,” Journal of Computational and Applied Mathematics, vol. 235, no. 11, pp. 3285–3290, 2011.

L. L. Wei, Y. N. He, and B. Tang, “Analysis of the fractional Kawahara equation using an implicit fully discrete local discontinuous Galerkin method,” Numerical Methods for Partial Differential Equations, 2012.

L. L. Wei, X. D. Zhang, S. Kumar, and A. Yildirim, “A numerical study based on an implicit fully discrete local discontinuous Galerkin method for the time-fractional coupled Schrödinger system,” Computers & Mathematics with Applications, vol. 64, no. 8, pp. 2603–2615, 2012.

L. L. Wei, Y. N. He, and Y. Zhang, “Numerical analysis of the fractional seventhorder KdV equation using an implicit fully discrete local discontinuous Galerkin method,” International Journal of Numerical Analysis & Modeling, vol. 10, pp. 430–444, 2013.

P. Zhuang and F. Liu, “Implicit difference approximation for the time fractional diffusion equation,” Journal of Applied Mathematics & Computing, vol. 22, no. 3, pp. 87–99, 2006.

S. B. Yuste and L. Acedo, “An explicit finite difference method and a new von Neumann-type stability analysis for fractional diffusion equations,” SIAM Journal on Numerical Analysis, vol. 42, no. 5, pp. 1862–1874, 2005.

C.-M. Chen, F. Liu, I. Turner, and V. Anh, “Numerical schemes and multivariate extrapolation of a two-dimensional anomalous sub-diffusion equation,” Numerical Algorithms, vol. 54, no. 1, pp. 1–21, 2010.

C.-M. Chen, F. Liu, V. Anh, and I. Turner, “Numerical methods for solving a two-dimensional variable-order anomalous subdiffusion equation,” Mathematics of Computation, vol. 81, no. 277, pp. 345–366, 2012.

X. D. Zhang, J. Liu, L. L. Wei, and C. X. Ma, “Finite element method for Grunwald-Letnikov time-fractional partial differential equation,” Journal of Applied Analysis, vol. 1, pp. 1–12, 2012.

L. L. Wei, X.D. Zhang, and Y. N. He, “Analysis of a local discontinuous Galerkin method for time-fractional advection-diffusion equations,” International Journal of Numerical Methods for Heat & Fluid Flow, vol. 23, no. 4, pp. 634–648, 2013.

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.

Y.-N. Zhang and Z.-Z. Sun, “Alternating direction implicit schemes for the two-dimensional fractional sub-diffusion equation,” Journal of Computational Physics, vol. 230, no. 24, pp. 8713–8728, 2011.

M. R. Cui, “Compact alternating direction implicit method for two-dimensional time fractional diffusion equation,” Journal of Computational Physics, vol. 231, no. 6, pp. 2621–2633, 2012.

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.

M. M. Meerschaert, H.-P. Scheffler, and C. Tadjeran, “Finite difference methods for two-dimensional fractional dispersion equation,” Journal of Computational Physics, vol. 211, no. 1, pp. 249–261, 2006.

C. Tadjeran and M. M. Meerschaert, “A second-order accurate numerical method for the two-dimensional fractional diffusion equation,” Journal of Computational Physics, vol. 220, no. 2, pp. 813–823, 2007.

S. Chen, F. Liu, I. Turner, and V. Anh, “An implicit numerical method for the two-dimensional fractional percolation equation,” Applied Mathematics and Computation, vol. 219, no. 9, pp. 4322–4331, 2013.

V. J. Ervin, N. Heuer, and J. P. Roop, “Numerical approximation of a time dependent, nonlinear, space-fractional diffusion equation,” SIAM Journal on Numerical Analysis, vol. 45, no. 2, pp. 572–591, 2007.

H. Zhang, F. Liu, and V. Anh, “Galerkin finite element approximation of symmetric space-fractional partial differential equations,” Applied Mathematics and Computation, vol. 217, no. 6, pp. 2534–2545, 2010.

Y. Y. Zheng, C. P. Li, and Z. G. 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.

W. Y. Liao, “A compact high-order finite difference method for unsteady convection-diffusion equation,” International Journal for Computational Methods in Engineering Science and Mechanics, vol. 13, no. 3, pp. 135–145, 2012.

J. Douglas Jr. and S. Kim, “Improved accuracy for locally one-dimensional methods for parabolic equations,” Mathematical Models & Methods in Applied Sciences, vol. 11, no. 9, pp. 1563–1579, 2001.