An Efficient Parallel Approximate Algorithm for Solving Time Fractional Reaction-Diffusion Equations
In this paper, we construct pure alternative segment explicit-implicit (PASE-I) and implicit-explicit (PASI-E) difference algorithms for time fractional reaction-diffusion equations (FRDEs). They are a kind of difference schemes with intrinsic parallelism and based on classical explicit scheme and classical implicit scheme combined with alternating segment technology. The existence and uniqueness analysis of solutions of the parallel difference schemes are given. Both the theoretical proof and the numerical experiment show that PASE-I and PASI-E schemes are unconditionally stable and convergent with second-order spatial accuracy and order time accuracy. Compared with implicit scheme and E-I (I-E) scheme, the computational efficiency of PASE-I and PASI-E schemes is greatly improved. PASE-I and PASI-E schemes have obvious parallel computing properties, which shows that the difference schemes with intrinsic parallelism in this paper are feasible to solve the time FRDEs.
Fractional differential equations can be used to describe some physical phenomena more accurately than the classical integer order differential equations. The time fractional reaction-diffusion equations (FRDEs) play an important role in dynamical systems of physics [1, 2], bioinformatics [3, 4], image processing , and other research areas [6, 7]. With the further application of FRDEs, the numerical solutions of FRDEs have become an urgent research work. The analytical solutions of FRDEs are difficult to give explicitly, even the analytical solution of linear fractional differential equation contains special functions, such as Mittag-Leffler function. The series corresponding to these functions converge very slowly; therefore, it is quite difficult to calculate these special functions in practical application. It makes the development of efficient numerical algorithms for FRDEs particularly important [8–10].
In recent years, many scholars have studied the numerical algorithm of the fractional diffusion model [11–13]. The finite difference method is still dominant among the existing algorithms now. Tadjeran et al.  examined a second-order accurate numerical method in time and in space to solve a class of initial-boundary value fractional diffusion equations with variable coefficients. Gao and Sun  derived a compact finite difference scheme for the subdiffusion equation, which is fourth-order accuracy compact approximate for the second-order space derivative. For time fractional subdiffusion equation with Dirichelt boundary value conditions, Luo et al.  established a novel collocation method via taking quadratic spline polynomialsas basic functions with order accuracy, when the solution has four-order continual derivative with respects to and . Moreover, the accuracy of the proposed method does not depend on the order . Luo et al.  proposed a novel numerical approximate method for the Caputo fractional derivative by employing the piecewise linear and quadratic Lagrange interpolation functions, which convergence order is with . Yaseen et al.  proposed a finite difference scheme with a cubic trigonometric B-spline function for time fractional diffusion-wave equation with reaction term. The scheme is unconditional stabile and convergent. For a class of one-dimensional and two-dimensional time FRDE with variable coefficients and time drift term, Zhao and Gu  presented an implicit finite difference scheme based on the formula. The unconditional stability and convergence of this scheme are proved rigorously by the discrete energy method, and the optimal convergence order in the is . For semilinear parabolic singularly perturbed systems of reaction-diffusion type, Clavero and Jorge  constructed a numerical method, which combined the central finite difference schemes to discretize in space and a linearized fractional implicit Euler method together with a splitting by components technique to integrate in time. The computational cost of the algorithm is considerably less than that of classical schemes. For three kinds of typical Caputo-type PDEs (including Caputo-type reaction-diffusion equation, reaction-diffusion-wave equation, and cable equation), Li and Wang [21, 22] proposed the local discontinuous Galerkin finite element methods. They studied the existence, uniqueness, and regularity of the solutions to these equations.
However, the existing implicit serial difference schemes are of relatively high computational complexity and relatively low computation efficiency. The reason is the fractional differential operators have nonlocal property [4, 11]. So, its numerical methods have full coefficient matrices, which require the storage of and computational cost of . Even with high-performance computers, it is difficult to simulate long time or large computational domain. How to construct an effectively numerical algorithm with high precision and low computational complexity becomes an urgent problem.
In order to overcome the difficulty of numerical algorithms of fractional order differential equations which required large amount of computational cost and storage, there are two main types of solutions in the existing literature. One is the fast algorithm which needs to combine the characteristics of the algorithm for reducing computational complexity of the algorithms. The another is the parallel algorithm which makes full use of multicore and cluster devices. Parallel algorithms have been successfully applied in integer order differential equations [23, 24].
The research achievements of fast algorithms of fractional differential equations are as follows. Wang et al.  proposed a fast algorithm for solving the spatial fractional diffusion equation based on the special structure of difference scheme. Lu et al.  proposed a fast numerical difference algorithm for time fractional subdiffusion equation. For time space fractional convection-diffusion equations, Gu et al.  presented a new unconditionally stable implicit difference method, which converges with the second-order accuracy in both time and space variables. Then, they solved the Toeplitz-like systems of linear equations by Krylov subspace solvers with suitable circulant preconditioners. Each time level of these methods reduces the memory requirement of the proposed implicit difference scheme from to and the computational complexity from to in each iterative step. Jiang et al.  have developed a fast algorithm for the evaluation of the Caputo fractional derivative . The algorithm relies on an efficient sum-of-exponentials approximation for the convolution kernel with the absolute error over the interval . The algorithm has nearly optimal complexity in both CPU time and storage. Zhao et al.  discussed a block lower triangular Toeplitz system arising from the time space fractional diffusion equation and applied the preconditioned biconjugate gradient stabilized method and the flexible general minimal residual method to obtain the solutions of the linear system. They develpoed a block bidiagonal Toeplitz preconditioner for the block lower triangular Toeplitz system, whose storage is of with being the spatial grid number.
With the rapid development of multicore and cluster technology, parallel algorithms have become one of the mainstream technologies to improve computing efficiency. For integer order parabolic equation, Evans and Abdullab  put forward the idea of group explicit (GE), and designed the alternating group explicit (AGE) scheme, which not only guarantees the stability of numerical calculation, but also has good parallel property. Inspired by the AGE method, Zhang et al.  proposed the idea of constructing alternating segment implicit scheme by using Saul’yev asymmetric scheme, then established a variety of alternating segment explicit-implicit (ASE-I) and alternating segment pure implicit parallel difference formulas by appropriate use of alternating technique. The methods have been well applied to numerical solving integer partial differential equations (PDEs) [32, 33]. The existing numerical methods cannot be applied to numerically solve fractional order PDEs directly and even produce completely different numerical analysis. How to extend the existing parallel difference method of integer order PDEs to fractional order PDEs is a great challenge to the computational mathematics (physics) fields.
The research achievements of parallel algorithms for fractional differential equations are as follows. Diethelm  implements parallel computation for Adams-Bashforth-Moulton method of the second order of the fractional derivative and discusses the accuracy of the parallel algorithm. Gong et al. [35, 36] parallelizes the explicit difference scheme of the time fractional diffusion equation and implicit scheme of the Riesz type space FRDE. The core of parallelization is the product of matrix and vector, the addition vector and vector. Sweilam et al.  constructed a class of parallel Crank-Nicolson difference schemes for time fractional parabolic equations. The preconditioned conjugate gradient method was applied to solve discrete algebraic equations in parallel. Wang et al.  studies the parallel algorithm of the implicit difference scheme for Caputo-type FRDE. Compared with the classical schemes, the parallel algorithm has some improvements in computational efficiency. Under the principle of minimizing communication, reasonably allocating computing tasks and minimizing the change of the original serial difference format, the serial algorithm is parallelized. Biala and Khaliq  developed a time stepping scheme for solving nonlinear time space fractional PDEs. The schemes are implemented in parallel using the distributed (MPI), shared memory systems (OpenMP), and a combination of both (Hybrid). They discussed the advantages of the parallel algorithms over the sequential version. Fu and Wang  developed a parallel finite difference scheme for space-time fractional PDEs. By exploring the structure of the stiffness matrices, they developed a matrix-free preconditioned fast Krylov subspace iterative method with significantly reduced computational work and memory requirement. Most of the parallel algorithms of fractional PDEs are from the perspective of numerical algebra.
For a long time, most of the existing parallel schemes are conditionally stable or unconditionally stable but with only first order accuracy in space. In order to obtain more accurate and unconditionally stable parallel difference schemes, we do not study parallel algorithms from the point of numerical algebra, but based on the parallelization of traditional difference schemes. When faced with long time history problem, large computational domain problem or high-dimensional problem, the amount of computation and storage of serial algorithms will increase exponentially . At this point, considering the calculation accuracy and efficiency, the parallel algorithm with second-order spatial accuracy will be the good choice. In this paper, for timeFRDEs, the pure alternative segment explicit-implicit (PASE-I) and implicit-explicit (PASI-E) parallel difference schemes are obtained by processing the explicit and implicit difference schemes with alternating segment technique. Comparing with the existing difference schemes, this kind of schemes can greatly improve the efficiency of numerical simulation of fractional order models.
The rest of the paper is organized as follows. In Section 2, we structure the PASE-I parallel difference scheme of time FRDE. In Section 3, we prove the uniqueness, stability, and convergence of PASE-I scheme and PASI-E difference scheme. In Section 4, we present some numerical experiments to demonstrate the accuracy and efficiency of the PASE-I and PASI-E methods. Finally, some concluding remarks are presented in Section 5.
2. Parallel Intrinsic Difference Schemes for Time Fractional Reaction-Diffusion Equation
2.1. Time Fractional Reaction-Diffusion Equation
When , by using the finite sinusoidal transformation and Laplace transformation, the analytic solution of equation (1) can be obtained [1, 8]:where , is a special function of type Mittag-Leffler, andwhere is error function, .
2.2. Construction of PASE-I Parallel Difference Scheme
We will divide the computational domain into mesh by taking the space step and the time step . The mesh points are , here , ; , ; are positive integers. Let the value of be approximated to the value of .
The time fractional derivative term can be approximated by the following format:
Let , .
In order to construct PASE-I scheme of equation (1), the classical explicit and implicit schemes of equation (1) are given. Let , .(1)The classical explicit scheme Simplify,(2)The classical implicit scheme
The classical explicit scheme (7) has ideal parallelism and is very suitable for parallel computing. However, it is conditionally stable, especially in multidimensional problems its time step is restricted severely. The classical implicit scheme (9) is unconditionally stable, and usually discrete into tridiagonal systems. The tridiagonal systems are solved in complexity using the Thomas algorithm. It is inconvenient to obtain the numerical results directly and quickly. We will apply the classical explicit scheme (7) and implicit scheme (9) combining with alternating segmentation technology to construct PASE-I scheme.
We construct PASE-I scheme as follows. Let ( are integer numbers and is an odd number). On the same time level, the space grid points are divided into segments, which are recorded in the order of . On the odd time level, each segment is calculated according to the rule of “classical explicit-classical implicit-classical explicit” from left to right; on the even time level, each segment is calculated according to the order of “classical implicit-classical explicit-classical implicit” from left to right. For each segment, the explicit and implicit schemes are applied alternately in time. A piecewise sketch of grid points is shown in Figure 1. In Figure 1, denotes the classical explicit format at the grid point, and denotes classical implicit schemes, respectively. The internal boundary of each implicit segment depended on the value of the first or last point of the adjacent explicit segment.
Therefore, the PASE-I scheme of equation (1) is constructed as follows:where , is order zero matrix.
3. Numerical Analysis of PASE-I Difference Scheme
3.1. Uniqueness of PASE-I Difference Scheme’s Solution
Theorem 1. The pure alternative segment explicit-implicit (PASE-I) difference scheme (12) for the time fractional reaction-diffusion equation is uniquely solvable.
Proof. From the definition of and , we can know that and are strictly diagonally dominant matrices, and the main diagonal elements of the two matrices are positive real numbers. and are nonsingular matrices. From that is a nonsingular matrix, the second equation of PASE-I scheme (12) is uniquely solvable. From that is a nonsingular matrix, the first and third equations of PASE-I scheme (12) are uniquely solvable. Therefore, Theorem 1 can be obtained.
3.2. Stability of PASE-I Difference Scheme
Lemma 1. If the matrix is a nonnegative real matrix, then there is an estimation formula for any parameter .
Proof. Make a transformation , thenTherefore, . That is the proof of Lemma 1.
For stability, assume that is the numerical solution of the PASE-I difference scheme and is the numerical solution of the PASE-I difference scheme from another initial value. Error , The error is substituted into scheme (12) to obtain the following formula:Let is any eigenvalue of the matrix , and meet takes the maximum value and takes the minimum value. is any eigenvalue of the matrix , and meet takes the maximum value and takes the minimum value. Because the minimum eigenvalues of the matrices and are both 0, which means that . Therefore, . From Lemma 1, let , we can knowWhen , we have . When , let , here . From , we can get and . .
When ,When ,Suppose when , the equation is established. When , we haveWhen , we have . When , let , here ; from , we can get , ; and . Therefore, .
When ,In the same way, we can get .
To sum up, , and the following theorem can be obtained.
Theorem 2. The pure alternative segment explicit-implicit (PASE-I) difference scheme (12) for time fractional reaction-diffusion equation is unconditionally stable.
3.3. Convergence of PASE-I Difference Scheme
First, the accuracy of PASE-I scheme is analyzed. The PASE-I scheme is alternatively applied explicit and implicit segments, not only in the direction of space, but also in the direction of time. Examination of grid points and . The classical explicit scheme is
The classical implicit scheme is
The two terms of have the same form as two terms of , but the opposite of the symbol. When the explicit format and implicit format are applied alternately, the errors of the two terms will be offset. And . When applying explicit and implicit formats alternately, the accuracy is 2 order in space direction and in time direction. Similarly, the same results can be obtained when implicit and explicit formats are applied alternately. Therefore, the truncation of PASE-I scheme is .
When , the first layer is calculated applying the explicit format.
When , the second layer is calculated applying the implicit format.
Define two grid functions as follows:where . We do Fourier expansion on and :where
Define norm , .
We have the following formula:
Set , and substitute (33), we have
Because , (38) can be converted to
Lemma 2. Let be the solution of (39), makes