Fourth-order Compact Iterative Scheme for the Two-dimensional Time Fractional Sub-diffusion Equations

Abstract The fractional diffusion equation is an important mathematical model for describing phenomena of anomalous diffusion in transport processes. A high-order compact iterative scheme is formulated in solving the two-dimensional time fractional sub-diffusion equation. The spatial derivative is evaluated using Crank-Nicolson scheme with a fourth-order compact approximation and the Caputo derivative is used for the time fractional derivative to obtain a discrete implicit scheme. The order of convergence for the proposed method will be shown to be of O(τ3−α + h). Numerical examples are provided to verify the high-order accuracy solutions of the proposed scheme.


Introduction
Fractional differential equations have gained more attention over the last few years because of their application in various fields of science and technology such as in physics, engineering, and biology [1][2][3][4]. Most of the fractional differential equations, however cannot be solved analytically; therefore researchers turn to numerical methods for alternative solvers. Numerous numerical methods have been formulated in solving various types of fractional differential equations [5][6][7][8][9][10][11][12][13][14][15]. In particular, two-dimensional time fractional diffusion equations have been solved by Zhuang and Liu [12], Balasim and Ali [13], Cui [14], and Abbaszadeh and Mohebbi [15] with promising results. However, formulations of high order ac-curacy solvers are still in its infancy particularly for two-dimensional fractional differential equations. In this paper, we present a high order method using Caputo time derivatives in hybrid with Crank-Nicolson scheme for the spatial derivatives in solving the two-dimensional time fractional diffusion equations. The method will be shown to be fourth-order accurate in space. The formulation of the scheme will be discussed in the next section. Numerical results are presented in section 3 and concluding remarks are given in section 4.

Formulation of the proposed scheme
The two-dimensional time fractional diffusion is described as Here C 0 D α t u (0 < α < 1) represents Caputo fractional derivative of defined by [12], i.e.
where Γ(.) is a Euler Gamma function. Using finite difference approximations to the time and space derivatives of (1), let h > 0 be the space step and k > 0 be the time step and the step size is taken be equal in both x and y directions. Define x = ih, y = jh, {i, j = 0, 1, 2, ..., n}, t k = τ k, k = 0, 1, 2, ..., l and steps sizes of spatial variables h x = h y = h = 1 n , where n is an arbitrary positive integer. Many approximation formulas could be obtained for (1) at the point (x i , y j , t k ). Consider a Taylor series expansion of each of function values of u(x i , y j , t k ) about the point (x i , y j , t k ): If the control operator is introduced as and then using above Taylor series expansions at points u k i+1,j and u k i,j+1 , we get (4) After re-arranging (3) and (4), the following is obtained For the time fractional derivative, we use the Crank-Nicolson approximation for Caputo derivative [13] ∂ α u(xi,yj ,t k+ 1 where b s = (s + 1) 1+α − (s) 1+α , s = 0, 1, 2, ..., n. Since we use Crank Nicolson scheme which is the average of implicit and explicit scheme, so replacing k by k + 1 2 in (5), (6) and substituting (5), (6) and (7) into (1), we have the following Multiplying both sides by 1 + 1 12 δ x 2 1 + 1 12 δ y 2 Γ(2−α) and after rearranging we get Since we know After simplifying to the point u k+1 i,j , we get where F k+ 1 Figure (1) shows nine grid points involved in the in the updates using (9). A compact scheme can be constructed by iterating on each point in the solution domain using Equation (9) until a certain convergence is achieved. Fig. 2 shows all the points involved at the different time levels in updating using Equation (9), where

Numerical experiments
Two examples are used to verify the effectiveness of the high-order compact scheme in solving the two-dimensional time fractional sub-diffusion equations. The experiments were conducted using the method formulated in Section II in  hybrid with a Successive Over Relaxation (SOR) technique for different mesh sizes ( n =8, 16, 24, 30,36) and for different time steps (τ = 0.12, 0.1, 0.05, 0.062, 0.041, 0.033, 0.025, 0.02, 0.015). Results were obtained using a PC with Core i7, 3.40 GHz, 4GB of RAM with Windows 7 Professional and Math-ematica software. The maximum absolute error (L ∞ ) with tolerance ε = 10 −5 was used for the conver-gence criteria. We calculated the computational orders of the proposed method in space variables with [15] where L ∞ is the maximum norm.

The analytical solution is
u(x, y, t) = e x+y t 1+α .
From Table 3 and Table 4, it is observed that with decreasing mesh size, maximum and average errors are also reduced; the effectiveness of the scheme is more pronounced. In Table 1 and Table 2, the C 2 −order of convergence is checked for the different values of in Example 2 and Example 1 respectively, which shows that the computational spatial accuracy of the scheme is in agreement with the theoretical spatial accuracy.

Conclusions
We have solved the two-dimensional time fractional subdiffusion equation with a higher-order compact Crank-Nicolson scheme using Caputo derivative. It is observed that by decreasing the grid size, the maximum norm (L ∞ ) error educe quite significantly, which shows that our proposed scheme is accurate and reliable. The theoretical order of convergence for the scheme was proven to be of O(τ 3−α + h 4 ). C 2 computational orders of the spatial accuracy of the numerical results have also been shown to be in agreement with the theoretical spatial accuracy.