Next Article in Journal
The Fixed Point Property of the Infinite M-Sphere
Previous Article in Journal
Efficiency of China’s Listed Securities Companies: Estimation through a DEA-Based Method
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A New Kind of Parallel Natural Difference Method for Multi-Term Time Fractional Diffusion Model

School of Mathematics and Physics, North China Electric Power University, Beijing 102206, China
*
Author to whom correspondence should be addressed.
Mathematics 2020, 8(4), 596; https://doi.org/10.3390/math8040596
Submission received: 10 March 2020 / Revised: 7 April 2020 / Accepted: 9 April 2020 / Published: 15 April 2020
(This article belongs to the Section Computational and Applied Mathematics)

Abstract

:
Multi-term time fractional diffusion model is not only an important physical subject, but also a practical problem commonly involved in engineering. In this paper, we apply the alternating segment technique to combine the classical explicit and implicit schemes, and propose a parallel nature difference method alternating segment pure explicit–implicit (PASE-I) and alternating segment pure implicit–explicit (PASI-E) difference schemes for multi-term time fractional order diffusion equations. The existence and uniqueness of the solutions are proved, and stability and convergence analysis of the two schemes are also given. Theoretical analyses and numerical experiments show that the PASE-I and PASI-E schemes are unconditionally stable and satisfy second-order accuracy in spatial precision and 2 α order in time precision. When the computational accuracy is equivalent, the CPU time of the two schemes are reduced by up to 2 / 3 compared with the classical implicit difference method. It indicates that the PASE-I and PASI-E parallel difference methods are efficient and feasible for solving multi-term time fractional diffusion equations.

1. Introduction

Fractional differential equations have been widely used in medicine, mechanics, control theory, environmental science, and finance [1,2,3,4]. In addition, fractional differential equations have unique advantages in describing the memory and genetic properties of matter. Therefore, the numerical solution of fractional differential equations has become one of the more active research fields in the world [5,6,7].
It is well known that the diffusion process of pollutants in complex structural soils or underground aquifers can span multiple scales in many cases, and needs to consider the impact of changes in time scales on the entire process. This paper considers multi-term time fractional diffusion equation [2,5].
P α , α 1 , , α m ( D t ) u ( x , t ) = 2 u ( x , t ) x 2 + f ( x , t ) , ( x , t ) ( 0 , L ) × ( 0 , T ] , u ( x , 0 ) = μ ( x ) , x [ 0 , L ] , u ( 0 , t ) = ϕ ( t ) , u ( L , t ) = v ( t ) , t ( 0 , T ] .
where μ ( x ) , v ( t ) , and ϕ ( t ) are three functions that are known to be properly smooth. P α , α 1 , , α m ( D t ) is defined as
P α , α 1 , , α m ( D t ) = D t α + i = 1 m l i D t α i , l i > 0 , m N + , 0 < α 1 < α 2 < α m < α < 1 .
D t α u ( x , t ) is a Caputo derivative operator for t, and is defined as
D t α u ( x , t ) = 1 Γ ( 1 α ) 0 t u ( x , ξ ) ξ d ξ ( t ξ ) α , 0 < α < 1 .
where Γ ( · ) is Gamma function.
For the study of multi-term time fractional diffusion equations, Luchko (2011) [8] applied the Fourier variable separation method and multi-term Mittag–Leffler (M-L) functions to construct a generalized solution. Li et al. (2015) [9] treated the low order fractional term as the perturbation of the high order fractional term, and used some important properties of M-L function and feature function expansion method to study the initial and boundary value problems of multi-time fractional differential equations. With the fractional Laplace operator, Luchko’s theorem, and multi-term M-L functions, Sin et al. (2017) [10] gave the analytical solutions of multi-term time-space Caputo–Riesz fractional diffusion equations in finite regions. However, the analytical solutions of fractional differential equations are difficult to be explicitly obtained in general. Even if the analytical solutions are given, most of them contain special functions, such as M-L functions, Wright function, H-function, hypergeometric function, etc. Therefore, it is particularly important to develop efficient numerical algorithms for solving fractional differential equations.
There are some related results on the numerical solution of multi-term time fractional differential equations. Ye et al. (2014) [11] proved a maximum principle for the multi-term time-space Riesz–Caputo fractional differential equations and applied a fractional predictor–corrector method combining the L1 and L2 discrete schemes to numerical solving the equation. The L1 and L2 formulas approximate the Caputo fractional derivative D t α u ( x , t ) of 0 < α < 1 and 1 < α < 2 by piecewise interpolation, respectively. For the multi-term time fractional diffusion equation in bounded convex polyhedron domain, Jin et al. (2015) [12] gave the numerical algorithm based on the standard Galerkin finite element method of space discretization and the finite difference method of time discretization, and discussed its stability and error estimation. Shiralashetti and Deshi (2016) [13] applied the Haar wavelet collocation method for solving multi-term fractional differential equations using the fractional order operational matrix of integration. Li et al. (2018) [14], based on the mixed finite-element method and finite difference method, gave the numerical algorithms of the multi-term time-fractional diffusion equations and diffusion-wave equations with Caputo fractional derivative. The unconditional stability and convergence results are proved. Wang et al. (2018) [15] applied the conforming triangular element method to numerically solve the two-dimensional multi-time fractional diffusion equation and carry out the accuracy analysis. Based on the bilinear finite element method in the spatial direction and L1 formula and Crank–Nicolson (C-N) formula in the temporal direction, Wei et al. (2018) [16] established unconditionally stable fully discrete approximation schemes for two-term mixed time fractional diffusion-wave equations and proved their unconditional stability. Furthermore, the superconvergence result is obtained, which improves the accuracy of numerical approximation without greatly increasing the computational complexity.
The numerical methods of fractional partial differential equations are still dominated by the finite difference method and the series method. The theoretical analysis tools mainly include Fourier method, energy estimation, matrix method, and mathematical induction. The research on the finite difference method for solving multi-term time fractional diffusion equations is as follows. Liu et al. (2013) [17] gave two kinds of implicit difference algorithms for multi-term time fractional wave equations with nonhomogeneous Dirichlet boundary conditions. Ren and Sun (2014) [18] proposed a high precision difference algorithm for one- and two-dimensional multi-term time fractional diffusion equations. For multi-term time fractional diffusion-wave equation, Dehghan et al. (2015) [19] combined with finite difference method and Galerkin spectral method to give a numerical algorithm with fourth-order precision. The energy method was applied to prove that the algorithm is unconditionally stable. For multi-term time fractional diffusion equation, Gao et al. (2017) [20] firstly constructed a numerical difference formula with second-order accuracy to approximate multiple Caputo type time fractional derivatives, and then applied fourth-order difference format in space and proposed the high precision difference scheme for time fractional diffusion equations. Yang et al. (2019) [21] applied explicit–implicit and implicit–explicit difference schemes for numerical solving double-term time fractional sub-diffusion equation. However, the existing numerical methods for solving multi-term time fractional diffusion equations are mostly serial algorithms. The computational efficiency of those algorithms is low and the computation time is long. It is difficult to simulate long-term or large computational domains, even with the application of high performance computers.
With the rapid development of multi-core and cluster technology, parallel algorithms have become one of the mainstream technologies to improve computational efficiency [22,23,24,25]. In recent years, some interesting results on the parallel algorithms of fractional partial differential equations have been obtained [26,27]. For one-dimensional space fractional diffusion equation, Wang et al. (2010) [28] presented a fast O(Nlog 2 N) algorithm for the difference scheme, based on the special structure of difference matrix in the constructed scheme. Diethelm (2011) [29] performed parallel computing on second-order Adams–Bashforth–Moulton method for fractional derivatives, and discussed the accuracy of the parallel method. Wang and Basu (2012) [30] further extended the fast algorithm to solve two-dimensional space fractional diffusion equation, which is an early attempt to apply the fast algorithm to the numerical simulation of fractional differential equation. Gong et al. (2013) [31] parallelized the explicit difference scheme of the space fractional reaction-diffusion equation. Sweilam et al. (2014) [32] constructed a class of parallel C-N difference schemes for time fractional parabolic equations. The core of the method is to use the precondition conjugate gradient method to solve discrete algebraic equations. Wang et al. (2016) [33] proposed an efficient parallel algorithm for Caputo fractional reaction–diffusion equation with implicit difference scheme. They developed a new tridiagonal reduced system with elimination method. Yang and Dang (2019) [34] constructed a class of improved alternating segment C-N difference scheme for time fractional reaction-diffusion equation. The parallel difference scheme has second-order spatial accuracy and 2 α -order temporal accuracy. Fu and Wang (2018) [35] developed a fast parareal finite difference method for space-time fractional partial differential equation. At each time step, they explored the structure of the stiffness matrix to develop a matrix-free preconditioned fast Krylov subspace iterative solver for the finite difference method without resorting to any lossy compression. Most of the algorithms for fractional partial differential equations are studied on the parallel algorithm of algebraic equations from the perspective of numerical algebra.
At present, the parallel algorithms for integer order differential equations are relatively mature [25,26]. However, the existing parallel algorithms cannot be directly applied to numerically solving fractional differential equations. To obtain parallel algorithms with higher precision and more relaxed stability conditions, we use the parallelization of traditional differential schemes, and hope to skip the difficulties of numerical algebra. As far as we know, research on parallel nature algorithms for multi-term time fractional diffusion equations has not been reported.
In this paper, we try to construct a class of alternating segment pure explicit–implicit (PASE-I) and pure implicit–explicit (PASI-E) parallel nature difference methods for solving multi-term time fractional diffusion equations. In Section 2, we alternately apply the classical explicit scheme and implicit scheme to segment the solution region, and obtain a new kind of parallel nature difference (PASE-I) scheme. Theoretical analysis of the PASE-I parallel scheme is given in Section 3. The construction and analysis of the PASI-E scheme are given in Section 4. Finally, numerical experiments are used to verify the correctness of theoretical analysis.

2. Construction of PASE-I Parallel Difference Scheme

Take two positive integers M and N and do equidistant rectangular meshing of the solution area { ( x , t ) | 0 x L , 0 t T } . Set h = L M for spatial mesh points x i = ( i 1 ) h , 1 i M + 1 and τ = T N for temporal mesh points t k = k τ , 0 k N , respectively.
For numerical approximation of the time Caputo derivative, the so-called L 1 formula is prepared below. The L 1 formula approximates the Caputo fractional derivative based on piecewise linear interpolation. Suppose f ( t ) C 2 [ t 0 , t k ] and 0 < α < 1 . Then,
D t α 0 C f ( t k ) = τ α Γ ( 2 α ) b 0 α f ( t k ) b k 1 α f ( t 0 ) j = 1 k 1 ( b j 1 α b j α ) f ( t k j ) + O ( τ 2 α ) ,
where b j α = ( j + 1 ) 1 α j 1 α , j 0 . τ is the time step [6,7].
Define the discrete operator L τ α u ( x i , t k ) : = τ α Γ ( 2 α ) [ b 0 α u ( x i , t k ) b k 1 α u ( x i , t 0 ) j = 1 k 1 ( b j 1 α b j α ) u ( x i , t k j ) ] . The discrete formula of P α , α 1 , , α m ( D t ) u ( x , t ) can be written as
P α , α 1 , , α m ( D t ) u ( x i , t k ) = L t α + i = 1 m l i L t α i u ( x i , t k ) .
Denote u i k = u ( x i , t k ) , f i k = f ( x i , t k ) . For 2 u x 2 in Equation (1), a second-order central difference formula is applied to approximate it. The classical explicit scheme and implicit scheme of Equation (1) are given, respectively.
The explicit difference scheme of Equation (1):
( L t α + i = 1 m l i L t α i ) u i k = r ( u i 1 k 1 2 u i k 1 + u i + 1 k 1 ) + f i k .
It can be written as
a 0 u i k = r u i 1 k 1 + ( a 0 a 1 2 r ) u i k 1 + r u j = 2 k 1 + j = 2 k 1 ( a j 1 a j ) u i k j + a k 1 u i 0 + f i k .
The implicit difference scheme of Equation (1):
( L τ α + i = 1 m l i L t α i ) u i k = r ( u i 1 k 2 u i k + u i + 1 k ) + f i k .
It can be written as
r u i 1 k + ( a 0 + 2 h 2 ) u i k r u i + 1 k = ( a 0 a 1 ) u i k 1 + j = 2 k 1 ( a j 1 a j ) u i k j + a k 1 u i 0 + f i k .
where r = 1 h 2 , a k = 1 Γ ( 2 α ) b k α + i = 1 m l i τ α α i Γ ( 2 α i ) b k α i , k = 1 , , N .
Before constructing the PASE-I parallel difference scheme, we firstly give the calculation format of the explicit segment and the implicit segment. For i 0 0 , consider the calculation of points ( i 0 + i , n + 1 ) , i = 1 , 2 , , l in an implicit (explicit) segment.
The implicit segment calculation scheme is
( a 0 I + A ) u i 0 + 1 n + 1 u i 0 + 2 n + 1 u i 0 + l 1 n + 1 u i 0 + l n + 1 = r u i 0 n + 1 0 0 r u i 0 + l + 1 n + 1 + j = 1 n w j u i 0 + 1 n + 1 j u i 0 + 2 n + 1 j u i 0 + l 1 n + 1 j u i 0 + l n + 1 j + a n u i 0 + 1 0 u i 0 + 2 0 u i 0 + l 1 0 u i 0 + l 0 + F n + 1 .
where w j = a j 1 a j , j = 1 , 2 , , n .
The explicit segment calculation scheme is
a 0 u i 0 + 1 n + 1 u i 0 + 2 n + 1 u i 0 + l 1 n + 1 u i 0 + l n + 1 = ( w 1 I A ) u i 0 + 1 n u i 0 + 2 n u i 0 + l 1 n u i 0 + l n + r u i 0 n 0 0 r u i 0 + l + 1 n + j = 2 n w j u i 0 + 1 n + 1 j u i 0 + 2 n + 1 j u i 0 + l 1 n + 1 j u i 0 + l n + 1 j + a n u i 0 + 1 0 u i 0 + 2 0 u i 0 + l 1 0 u i 0 + l 0 + F n .
where A = 2 r r r 2 r r r 2 r r r 2 r l × l .
Combining the classical explicit and implicit schemes and applied the alternating segment technique, the design of alternating segment pure explicit–implicit (PASE-I) for Equation (1) is as follows. Suppose M 1 = q l , where l N + and q is odd numbers ( l , q 3 ). The grid points to be calculated in the same even time layer are divided into q segments, which are sequentially calculated according to the rules of “(5)-(4)-(5)”. Similarly, the next odd layer is also divided into q segments, and the calculation rule becomes “(4)-(5)-(4)”. For example, the point diagram of the PASE-I scheme is shown in Figure 1, when q = 5 and l = 5 . Then, we get the PASE-I scheme for Equation (1) as follows.
( a 0 I + G 1 ) V n + 1 = ( w 1 I G 2 ) V n + w 2 V n 1 + w n V 1 + a n V 0 + b 1 n + F n + 1 , ( a 0 I + G 2 ) V n + 2 = ( w 1 I G 1 ) V n + 1 + w 2 V n + w n + 1 V 1 + a n + 1 V 0 + b 1 n + 2 + F n + 2 .
where b 1 n = ( r u 1 n , , r u M + 1 n ) , F n = ( f 2 n , f 3 n , , f M n ) , V n = ( u 2 n , u 3 n , , u M n ) , n = 0 , 2 , 4 , . I is M 1 order unit matrix, Q 1 is l 1 order zero matrix, and Q 2 is l 2 order zero matrix. G 1 and G 2 , M 1 order matrices, are defined as follows.
G 1 = Q 1 G 3 Q 2 G 3 Q 2 G 3 Q 1 ( M 1 ) × ( M 1 ) ,
G 2 = G 4 Q 2 G 3 Q 2 G 3 Q 2 G 5 ( M 1 ) × ( M 1 ) ,
where
G 3 = 0 0 r 2 r r r 2 r r 0 0 ( l + 2 ) × ( l + 2 ) ,
G 4 = 2 r r r 2 r r r 2 r r 0 0 ( l + 1 ) × ( l + 1 ) ,
G 5 = 0 0 r 2 r r r 2 r r r 2 r ( l + 1 ) × ( l + 1 ) .

3. Theoretical Analysis of PASE-I Difference Scheme

3.1. The Existence and Uniqueness of PASE-I Scheme’s Solution

Lemma 1.
The matrices a 0 I + G 1 and a 0 I + G 2 defined by the PASE-I scheme are nonsingular matrices.
Proof. 
It is known that a 0 I + G 1 is a strictly diagonally dominant matrix and the main diagonal elements are positive real numbers, from a 0 > 0 and the definition of G 1 . Thus, a 0 I + G 1 is a nonsingular matrix and ( a 0 I + G 1 ) 1 exists. a 0 I + G 2 is also a nonsingular matrix and ( a 0 I + G 2 ) 1 exists. Thus, there is the following theorem. □
Theorem 1.
The solution of the PASE-I scheme for solving multi-term time fractional diffusion Equation (1) is unique.

3.2. Stability of PASE-I Scheme

Lemma 2.
If the matrix C is a nonnegative definite matrix, and 0 σ 1 σ 2 , then there is an estimate for ρ 0
( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 2 1 .
Proof. 
( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 2 2 = max φ R n , φ 0 ( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 φ , ( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 φ ( φ , φ )
Make a transformation ψ = ( σ 2 I + ρ C ) 1 φ ; then,
( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 2 2 = max φ R n , φ 0 ( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 φ , ( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 φ ( φ , φ ) = max φ R n , φ 0 ( σ 1 I ρ C ) ψ , ( σ 1 I ρ C ) ψ ( ( σ 2 I + ρ C ) ψ , ( σ 2 I + ρ C ) ψ ) = max φ R n , φ 0 σ 1 2 ( ψ , ψ ) 2 σ 1 ρ ( C ψ , ψ ) + ρ 2 ( C ψ , C ψ ) σ 2 2 ( ψ , ψ ) + 2 σ 2 ρ ( C ψ , ψ ) + ρ 2 ( C ψ , C ψ ) .
From 0 σ 1 σ 2 , we can get
( σ 1 I ρ C ) ( σ 2 I + ρ C ) 1 2 1 .
 □
We suppose that u ¯ i n is the approximate solution of Equation (1) and the error ε i n = u i n u ¯ i n , E n = ( ε 2 n , ε 3 n , , ε M n ) , 1 n N + 1 . E n is introduced into the PASE-I scheme. We can get
( a 0 I + G 1 ) E n + 1 = ( w 1 I G 2 ) E n + w 2 E n 1 + w n E 1 + a n E 0 , ( a 0 I + G 2 ) E n + 2 = ( w 1 I G 1 ) E n + 1 + w 2 E n + w n + 1 E 1 + a n + 1 E 0 ,
where N = 0 , 2 , 4 , .
When n 2 ,
E n + 2 = ( a 0 I + G 2 ) 1 ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 1 I G 2 ) E n + ( a 0 I + G 2 ) 1 ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 2 E n 1 + + w n E 1 + a n E 0 ) + ( a 0 I + G 2 ) 1 ( w 2 E n + + w n + 1 E 1 + a n + 1 E 0 ) .
Taking the norm on both sides of the above equation, we can get
E n + 2 ( a 0 I + G 2 ) 1 ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 1 I G 2 ) E n + ( a 0 I + G 2 ) 1 ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 2 E n 1 + + w n E 1 + a n E 0 ) + ( a 0 I + G 2 ) 1 ( w 2 E n + + w n + 1 E 1 + a n + 1 E 0 ) .
The growth matrix of the PASE-I scheme is
T = ( a 0 I + G 2 ) 1 ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 1 I G 2 ) .
According to the definition of G 1 , G 2 , matrices G 1 and G 2 have the same eigenvalues. Suppose that the eigenvalues of G 1 and G 2 are λ . Let T ˜ = ( a 0 I + G 2 ) T ( a 0 I + G 2 ) 1 . Then, we have
T = T ˜ = ( w 1 I G 1 ) ( a 0 I + G 1 ) 1 ( w 1 I G 2 ) ( a 0 I + G 2 ) 1 = max ( w 1 λ a 0 + λ ) 2 1 .
The following inequality E n E 0 is proved by mathematical induction.
When n = 0 , for E 1 ,
( a 0 I + G 1 ) E 1 = ( I G 2 ) E 0 , E 1 = ( a 0 I + G 1 ) 1 ( I G 2 ) E 0 E 0 .
For E 2 , Case 1. max w 1 , λ = w 1 ,
E 2 a 0 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 I G 2 E 0 + a 0 I + G 2 1 a 1 E 0 max w 1 λ a 0 + λ E 0 + max a 1 a 0 + λ E 0 max a 0 λ a 0 + λ E 0 E 0 .
Case 2. max w 1 , λ = λ ,
E 3 a 0 2 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 I G 2 E 0 + a 0 I + G 2 1 a 1 E 0 max w 1 λ a 0 + λ E 0 + max a 1 a 0 + λ E 0 max λ + a 1 w 1 a 0 + λ E 0 E 0 .
Assume that, when n k + 1 , the inequality E n E 0 is true. When n = k + 2 ,
Case 1. max w 1 λ max w 1 , λ w 1 ,
E k + 2 a 0 I + G ¯ 2 1 w 1 I G ¯ 1 a 0 I + G ¯ 1 1 w 1 I G ¯ 2 E k + a 0 I + G ¯ 2 1 w 1 I G ¯ 1 a 0 I + G ¯ 1 1 w 2 E k 1 + + w k E 1 + a k E 0 + a 0 I + G ¯ 2 1 w 2 E k + + w k + 1 E 1 + a k + 1 E 0 w 1 a 0 + λ 2 E 0 + w 1 1 w 1 a 0 + λ 2 E 0 + 1 w 1 a 0 + λ E 0 w 1 E 0 + 1 w 1 E 0 E 0 .
Case 2. max w 1 λ max w 1 , λ λ ,
E k + 2 a 0 I + G ¯ 2 1 w 1 I G ¯ 1 a 0 I + G ¯ 1 1 w 1 I G ¯ 2 E k + a 0 I + G ¯ 2 1 w 1 I G ¯ 1 a 0 I + G ¯ 1 1 w 2 E k 1 + + w k E 1 + a k E 0 + a 0 I + G ¯ 2 1 w 2 E k + + w k + 1 E 1 + a k + 1 E 0 λ a 0 + λ 2 E 0 + λ 1 w 1 a 0 + λ 2 E 0 + 1 w 1 a 0 + λ E 0 = λ a 0 + λ λ a 0 + λ + 1 w 1 a 0 + λ E 0 + 1 w 1 a 0 + λ E 0 λ a 0 + λ E 0 + 1 w 1 a 0 + λ E 0 E 0 .
In summary, there is the following theorem.
Theorem 2.
The PASE-I scheme of the multi-term time fractional diffusion Equation (1) is unconditionally stable.

3.3. Convergence of PASE-I Scheme

Firstly, the accuracy analyses of the explicit and implicit schemes are performed separately. The truncation errors of explicit and implicit schemes are T 1 ( τ , h ) and T 2 ( τ , h ) , respectively. The Taylor expansion is performed at grid point ( x i , t n + 1 ) . It is known that the discreteness of the D t α 0 C u ( x i , t n + 1 ) formula has 2 α order numerical precision [6,7].
The truncation error of explicit scheme T 1 ( τ , h ) is
T 1 ( τ , h ) = ( L t α + i = 1 m l i L t α i ) u i k 1 h 2 ( u i 1 k 1 2 u i k 1 + u i + 1 k 1 ) f i k = D t α u + i = 1 m l i D t α i u u x x f i k + τ u x x t + τ h 2 12 u x x x x t τ 2 2 u x x t t h 2 12 u x x x x + O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) = τ u x x t + τ h 2 12 u x x x x t τ 2 2 u x x t t h 2 12 u x x x x + O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) .
The truncation error of implicit scheme T 2 ( τ , h ) is
T 2 ( τ , h ) = ( L τ α + i = 1 m l i L t α i ) u i k 1 h 2 ( u i 1 k 2 u i k + u i + 1 k ) f i k = D t α u + i = 1 m l i D t α i u u x x f i k τ u x x t τ h 2 12 u x x x x t τ 2 2 u x x t t h 2 12 u x x x x + O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) = τ u x x t τ h 2 12 u x x x x t τ 2 2 u x x t t h 2 12 u x x x x + O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) .
For the PASE-I scheme, the explicit scheme and implicit scheme are alternately applied for each grid point in spatial direction. The coefficients of the two terms u x x t and u x x x x t are opposite numbers in T 1 ( τ , h ) , T 2 ( τ , h ) . Therefore, the two terms u x x t and u x x x x t of the truncation errors are offset for the PASE-I scheme. The accuracy of PASE-I scheme is O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) .
Suppose u ( x i , t n ) is the solution of Equation (1) at the mesh point ( x i , t n ) . The defined e i n = u ( x i , t n ) u i n , e n = ( e 2 n , e 3 n , , e M n ) T , ( 1 i M + 1 , 1 n N + 1 ) . e n and e 0 = 0 are substituted into the PASE-I scheme,
( a 0 I + G 1 ) e n + 1 = ( w 1 I G 2 ) e n + w 2 e n 1 + w n e 1 + a n e 0 + R , ( a 0 I + G 2 ) e n + 2 = ( w 1 I G 1 ) e n + 1 + w 2 e n + w n + 1 e 1 + a n + 1 e 0 + R ,
where R = τ α O ( τ 2 α + i = 1 m τ 2 α i + h 2 ) C ( τ 2 + h 2 τ α ) and C is a constant.
When n = 1 , for e 1 ,
e 1 = I + G 1 1 I G 2 e 0 + I + G 1 1 R = I + G 1 1 R .
Take the norm on both sides of above equation.
e 2 = I + G 1 1 R R C τ 2 + h 2 τ α = a 0 1 τ α C τ 2 α + h 2 .
For e 2 , we have
( a 0 I + G 2 ) e 2 = ( w 1 I G 1 ) e 1 + R .
Take the norm on both sides of above equation.
e 2 = ( a 0 I + G 2 ) 1 [ ( w 1 I G 1 ) e 1 + R ] ( a 0 I + G 2 ) 1 [ w 1 I G 1 + a 1 ] a 1 1 R .
Case 1. max { w 1 λ } max { w 1 , λ } w 1 ,
e 2 ( a 0 I + G 2 ) 1 [ w 1 I G 1 + a 1 ] a 1 1 R 1 a 0 [ ( w 1 + a 1 ) a 1 1 R ] a 1 1 R a 1 1 τ α C τ 2 α + h 2 .
Case 2. max { w 1 λ } max { w 1 , λ } λ ,
e 2 = ( a 0 I + G 2 ) 1 [ ( w 1 I G 1 ) e 1 + R ] ( a 0 I + G 2 ) 1 [ w 1 I G 1 + a 1 ] a 1 1 R 1 a 0 + λ [ ( λ w 1 + a 1 ) a 1 1 R ] 2 a 1 a 0 + λ a 0 + λ a 1 1 R a 1 1 R a 1 1 τ α C τ 2 α + h 2 .
when n k + 1 , assume that the inequality e n a n 1 1 R is true. When n = k + 2 ,
Case 1. max { w 1 λ } max { w 1 , λ } w 1 ,
e k + 2 a 0 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 w 1 I G 2 e k + a 0 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 w 2 e k 1 + + w k e 1 + a k e 0 + a 0 I + G 2 1 w 2 e k + + w k + 1 e 1 + a k + 1 e 0 + a 0 I + G 1 R w 1 a 0 + λ 2 R + w 1 1 w 1 a 0 + λ 2 R + 1 a 0 + λ w 2 + + w k + a k + 1 a k + 1 1 R w 1 2 a 0 + λ 2 + w 1 1 w 1 a 0 + λ 2 + 1 w 1 a 0 + λ a k + 1 1 R w 1 a 0 + λ 2 + 1 w 1 a 0 + λ a k + 1 1 R a k + 1 1 R a k + 1 1 τ α C τ 2 α + h 2 .
Case 2. max { w 1 λ } max { w 1 , λ } λ ,
e k + 2 a 0 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 w 1 I G 2 e k + a 0 I + G 2 1 w 1 I G 1 a 0 I + G 1 1 w 2 e k 1 + + w k e 1 + a k e 0 + a 0 I + G 2 1 w 2 e k + + w k + 1 e 1 + a k + 1 e 0 + a 0 I + G 1 R λ a 0 + λ 2 R + λ 1 w 1 a 0 + λ 2 R + 1 a 0 + λ w 2 + + w k + a k + 1 a k + 1 1 R λ 2 a 0 + λ 2 + λ 1 w 1 a 0 + λ 2 + 1 w 1 a 0 + λ a k + 1 1 R λ a 0 + λ λ a 0 + λ + 1 w 1 a 0 + λ + 1 w 1 a 0 + λ a k + 1 1 R λ a 0 + λ + 1 w 1 a 0 + λ a k + 1 1 R a k + 1 1 R a k + 1 1 τ α C τ 2 α + h 2 .
From
a n 1 n α = ( b n α + i = 1 m b n α i ) 1 n α ( b n α ) 1 n α 1 + i = 1 m ( b n α i ) 1 n α 2 ,
we have
lim n ( b n α 1 ) 1 n α 1 = lim n n α 1 n 1 α 1 ( n 1 ) 1 α 1 = lim n n 1 1 ( 1 1 n ) 1 α 1 = 1 1 α 1 ,
lim n a n 1 n α = lim n ( b n α ) 1 n α + lim n i = 1 m ( b n α i ) 1 n α = 1 1 α + i = 1 m 1 1 α i .
From the Equations (11), (12), (14) and (15), we have
e n + 1 a n 1 τ α C τ 2 α + h 2 a k + 1 1 n α n α τ α C τ 2 α + h 2 1 1 α + i = 1 m 1 1 α i ( n τ ) α C ( τ 2 α + h 2 ) 1 1 α + i = 1 m 1 1 α i T α C ( τ 2 α + h 2 ) C 1 ( τ 2 α + h 2 ) ,
where C 1 = ( 1 1 α + i = 1 m 1 1 α i ) T α C . We can get the following inequality.
u ( x i , t n ) u i n C 1 ( τ 2 α + h 2 ) .
In summary, the following theorem is obtained.
Theorem 3.
The PASE-I scheme of the multi-term time fractional diffusion Equation (1) is convergent, u ( x i , t n ) u i n C ( τ 2 α + h 2 ) , and C is a positive number.

4. PASI-E Parallel Difference Scheme

By changing the calculation order of the explicit segment and implicit segment, the PASI-E scheme of the multi-term time fractional diffusion Equation (1) can be obtained. The calculation rule is calculated in the order of “(4)-(5)-(4)” in the even time layer, and in the order of “(5)-(4)-(5)” in the odd time layer. Then, we can get the PASI-E scheme for solving multi-term time fractional diffusion Equations (1) as follows.
( a 0 I + G 2 ) V n + 1 = ( w 1 I G 1 ) V n + w 2 V n 1 + w n V 1 + a n V 0 + b 1 n + F n + 1 , ( a 0 I + G 1 ) V n + 2 = ( w 1 I G 2 ) V n + 1 + w 2 V n + w n + 1 V 1 + a n + 1 V 0 + b 1 n + 2 + F n + 2 ,
where n = 1 , 3 , 5 , The definition of G 1 , G 2 , F n , b 1 n is the same as above.
By the same proof process, there is the following theorem.
Theorem 4.
The PASI-E scheme of the multi-term time fractional diffusion Equation (1) is unconditionally stable and convergent, u ( x i , t n ) u i n C ( τ 2 α + h 2 ) , and C is a positive number.
Since the PASE-I scheme and the PASI-E scheme differ only in the order of calculation of the explicit and implicit schemes, the amount of computation of the two parallel schemes should theoretically be equivalent.

5. Numerical Experiments

The experiment platform was laptop with Intel(R) Core(TM) i5-2400 CPU, 4 GB main memory and Windows 7 operating system. The CPU clock frequency is 3.10 GHz. The code was developed with Matlab R2014b [36]. We consider the following multi-term time fractional diffusion equation [7,12],
D t α 1 u ( x , t ) + D t α 2 u ( x , t ) = 2 u ( x , t ) x 2 + f ( x , t ) , ( x , t ) ( 0 , 1 ) × ( 0 , 1 ] , u ( x , 0 ) = x ( 1 x ) , x [ 0 , 1 ] , u ( 0 , t ) = u ( 1 , t ) = 0 , t ( 0 , 1 ] .
where 0 < α 2 < α 1 < 1 , f ( x , t ) = 2 t 2 α 1 Γ ( 3 α 1 ) + 2 t 2 α 2 Γ ( 3 α 2 ) ( x 2 + x ) + 2 ( 1 + t 2 ) . The exact solution of the above equation is u ( x , t ) = ( 1 + t 2 ) ( x 2 + x ) .
Take M = 20 , N = 100 , α 1 = 0.9 , α 2 = 0.5 , the error surfaces of the PASE-I scheme and PASI-E scheme are shown in Figure 2. In Figure 2, the numerical solutions of the two parallel schemes are consistent with the exact solution.
Next, we verify the calculation precision and convergence order of PASE-I and PASI-E parallel difference schemes for solving multi-term time fractional diffusion equations. Define E ( h , τ ) = max 0 i M u ( x i , t n ) u i n , O r d e r x = log 2 E ( 2 h , τ ) E ( h , τ ) , and O r d e r t = log 2 E ( h , 2 τ ) E ( h , τ ) . Firstly, we choose the optimal step size τ 2 α h 2 and α 1 = α 2 for the calculation precision and convergence order in space [7]. Then, we choose the values of M as 25, 50, 100, 200, 400, and 800, separately, and the values of α 1 and α 2 as 0.4 , 0.5 , and 0.6 , separately. Table 1 shows that the PASE-I and PASI-E parallel difference methods have a convergence order of 2 in spatial direction. The accuracy of the two parallel difference schemes is almost the same.
For the convergence order of the time direction, we take three cases α 1 = α 2 = 0.35 , α 1 = 0.4 and α 2 = 0.2 , and α 1 = 0.8 and α 2 = 0.2 , respectively. Take M = 200 and the value of N as 200 , 400 , 800 , 1600 , 3200 , respectively. When α 1 = α 2 = 0.35 , the time convergence orders of PASE-I and PASI-E parallel difference schemes are about 1.65 . It is consistent with the theoretical analyses, and the time convergence order is 2 α . For α 1 = 0.4 and α 2 = 0.2 and α 1 = 0.8 and α 2 = 0.2 , the time convergence orders are 1.67 and 1.36 , respectively, which are slightly larger than the theoretical analysis order 2 α .
Suppose r 1 r 2 . We have
E ( ( 2 τ ) r 1 + ( 2 τ ) r 2 + h 2 ) E ( τ r 1 + τ r 2 + h 2 ) E ( ( 2 τ ) r 1 + ( 2 τ ) r 2 ) E ( τ r 1 + τ r 2 ) = 2 r 2 E ( ( 2 τ ) r 1 r 2 + 1 ) E ( τ r 1 r 2 + 1 ) 2 r 2 .
It can be seen from the above equation that the time convergence order is 2 α 1 for α 1 = α 2 . When α 1 α 2 , the time convergence order is slightly larger than 2 α , α = max { α 1 , α 2 } .
Next, we verify the stability and computational accuracy of two parallel difference schemes from the perspective of numerical experiments. Table 1 and Table 2 show that the accuracy of PASE-I and PASI-E schemes are similar, thus the following analysis is represented by PASE-I scheme. Let the numerical solution u i n of the difference scheme be the perturbation solution, and the exact solution u ( x i , t n ) is the control solution. Definition of the difference total energy (DTE) of the error is as follows: D T E ( i ) = 1 2 n = 1 N ( u ( x i , t n ) u i n ) 2 .
IN Figure 3, the DTE of PASE-I scheme is within 10 3 . With the encryption of the spatial grid, the DTE is gradually reduced. The PASE-I and PASI-E parallel difference schemes have good computational accuracy.
For the stability of the PASE-I and PASI-E schemes, we define the relative error (RE) as follows:
R E ( j ) = i = 1 M | u ( x i , t j ) u j i | u ( x i , t j ) .
Take the value of M as 200 and the value of N as 400 , 800 , 1600 , and 3200, respectively. Figure 4 shows that the RE of the PASE-I scheme is within a certain range. With the time step decreases, the RE becomes smaller and smaller. The speed of growth also slows down as the time step decreases. These demonstrate that the PASE-I parallel difference scheme is computationally stable and consistent with theoretical analysis.
Finally, we investigate the effect of increasing the number of spatial grid points on the computational complexity of serial and parallel difference schemes. We define the speedup as S p = T 1 / T p ( T 1 is the CPU time of implicit and T p is the CPU time of parallel scheme) [37]. When the number of spatial grid points is 100, 500, 1000, 2000, 3000, 4000, and 5000, respectively, the CPU time of the three schemes is shown in Figure 5 and Table 3.
It can be seen that, when the number of spatial grid points becomes larger, the parallel difference schemes in this paper show obvious superiority in computational efficiency in Figure 5. When M = 5000 , the CPU time of the two parallel schemes can be reduced by up to 2 / 3 compared with the serial (classical implicit) difference scheme. When the number of spatial grid points becomes smaller, the CPU time of the serial scheme and parallel scheme is similar. With the small number of spatial grid points, the influence of data communication on the program loop will greatly reduce the efficiency of parallel computing.
In Table 3, we can see that the speedup of PASE-I and PASI-E schemes will become more prominent with the increase of computational domain. When the number of spatial grids is small (100), the speedup of parallel difference scheme is near one, because the communication between modules consumes a lot of CPU time. When the number of grid points is 5000, the speedup of parallel difference schemes is optimal in this example. Therefore, data communication problems need to be considered in parallel programming. When the amount of data (number of spatial grid points) is large, the impact of program loop execution is much greater than the impact of data communication. In this case, parallel computing is more effective, and, as the number of spatial grid points increases, the efficiency of parallel computing becomes more obvious.

6. Conclusions

In this paper, the PASE-I and PASI-E schemes of multi-term time fractional diffusion equation are constructed. The theoretical analysis of the two parallel schemes shows that both schemes are unconditionally stable and convergent. The methods are simple and feasible, and keep high precision of calculation. Numerical experiments verify the theoretical analysis, indicating that the advantages of PASE-I and PASI-E parallel methods are more and more obvious compared with classical implicit difference scheme with the increase of grid points. It is feasible to solve multi-term time fractional diffusion equations by the PASE-I and PASI-E parallel natural difference methods. At the same time, the methods can be easily extended to solve two-dimensional problems, especially suitable for MIMD computers.

Author Contributions

Conceptualization, X.Y. and L.W.; methodology, X.Y.; software, L.W.; validation, L.W.; formal analysis, L.W.; data curation, X.Y. and L.W.; writing-original draft preparation, L.W.; writing-review and editing, X.Y.; visualization, L.W.; supervision, X.Y.; project administration, X.Y.; funding acquisition, X.Y. and L.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Subproject of Major Science and Technology Program of China (2017ZX07101001-1) and the Fundamental Research Funds of the Central Universities (2018MS168).

Acknowledgments

The authors would like to deeply thank Hai Huang (School of Mathematics and Physics, North China Electric Power University) for his valuable suggestion and constructive comments during the preparation.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Uchaikin, V.V. Fractional Derivatives for Physicists and Engineers: Volume II: Applications; Higher Education Press: Beijing, China, 2013. [Google Scholar]
  2. Chen, W.; Sun, H.G. Fractional Differential Equations and Statistical Models of Anomalous Diffusion; Science Press: Beijing, China, 2017. (In Chinese) [Google Scholar]
  3. Sabatier, J.; Agrawal, O.P.; Machado, J.A.T. (Eds.) Advances in Fractional Calculus: Theoretical Developments and Applications in Physics and Engineering; Beijing World Publishing Corporation: Beijing, China, 2014. [Google Scholar]
  4. Tejado, I.; Emiliano, P.; Duarte, V. Fractional Derivatives for Economic Growth Modelling of the Group of Twenty: Application to Prediction. Mathematics 2020, 8, 50. [Google Scholar] [CrossRef] [Green Version]
  5. Guo, B.L.; Pu, X.K.; Huang, F.H. Fractional Partial Differential Equations and Their Numerical Solutions; Science Press: Beijing, China, 2015. [Google Scholar]
  6. Liu, F.W.; Zhuang, P.H.; Liu, Q.X. Numerical Solutions of Fractional Partial Differential Equations and Their Application; Science Press: Beijing, China, 2015. (In Chinese) [Google Scholar]
  7. Sun, Z.Z.; Gao, G.H. Finite Difference Methods for Fractional Differential Equations; Science Press: Beijing, China, 2015. (In Chinese) [Google Scholar]
  8. Luchko, Y. Initial-boundary-value problems for the generalized multi-term time-fractional diffusion equation. Math. Anal. Appl. 2011, 374, 538–548. [Google Scholar] [CrossRef] [Green Version]
  9. Li, Z.Y.; Liu, Y.K.; Yamamoto, M. Initial-boundary value problems for multi-term time-fractional diffusion equations with positive constant coefficients. Appl. Math. Comput. 2015, 257, 381–397. [Google Scholar] [CrossRef] [Green Version]
  10. Sin, C.S.; Ri, G.I.; Kim, M.C. Analytical solutions to multi-term time-space Caputo-Riesz fractional diffusion equations on an infinite domain. Adv. Differ. Equ. 2017, 306. [Google Scholar] [CrossRef] [Green Version]
  11. Ye, H.; Liu, F.W.; Anh, V.; Turner, I. Maximum principle and numerical method for the multi-term time-space Riesz-Caputo fractional differential equations. Appl. Math. Comput. 2014, 227, 531–540. [Google Scholar] [CrossRef] [Green Version]
  12. Jin, B.; Lazarov, R.; Liu, Y.; Zhou, Z. The Galerkin finite element method for a multi-term time-fractional diffusion equation. J. Comput. Phys. 2015, 281, 825–843. [Google Scholar] [CrossRef] [Green Version]
  13. Shiralashetti, S.C.; Deshi, A.B. An efficient Haar wavelet collocation method for the numerical solution of multi-term fractional differential equations. Nonlinear Dynam. 2016, 83, 293–303. [Google Scholar] [CrossRef]
  14. Li, M.; Huang, C.M.; Ming, W.Y. Mixed finite-element method for multi-term time-fractional diffusion and diffusion-wave equations. Comp. Appl. Math. 2018, 37, 2309–2334. [Google Scholar] [CrossRef]
  15. Wang, F.L.; Fan, M.Z.; Zhan, Y.M.; Shi, Z.G.; Shi, D.Y. High accuracy analysis of anisotropic linear triangular element for multi-term time fractional diffusion equations. Math. Numer. Sin. 2018, 40, 299–312. (In Chinese) [Google Scholar]
  16. Wei, Y.B.; Zhao, Y.M.; Tang, Y.F.; Wang, F.; Shi, Z.; Li, K. High-accuracy analysis of finite-element method for two-term mixed time fractional diffusion-wave equations. Sci. Scinica Inform. 2018, 48, 871–887. (In Chinese) [Google Scholar] [CrossRef]
  17. Liu, F.W.; Meerschaert, M.M.; Mcgough, R.J.; Zhuang, P.; Liu, Q. Numerical methods for solving the multi-term time-fractional wave-diffusion equation. Fract. Calcul. Appl. Anal. 2013, 16, 9–25. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  18. Ren, J.H.; Sun, Z.Z. Efficient and stable numerical methods for multi-term time fractional sub-diffusion equations. East Asian J. Appl. Math. 2014, 4, 242–266. [Google Scholar] [CrossRef] [Green Version]
  19. Dehghan, M.; Safarpoor, M.; Abbaszadeh, M. Two high-order numerical algorithms for solving the multi-term time fractional diffusion-wave equations. Comp. Appl. Math. 2015, 290, 174–195. [Google Scholar] [CrossRef]
  20. Gao, G.H.; Alikhanov, A.A.; Sun, Z.Z. The temporal second order difference schemes based on the interpolation approximation for solving the time multi-term and distributed-order fractional sub-diffusion equations. J. Sci. Comput. 2017, 73, 93–121. [Google Scholar] [CrossRef]
  21. Yang, X.Z.; Shao, J.; Sun, S.Z. A class of efficient difference methods for the double-term time fractional sub-diffusion equation. Acta Math. Appl. Sin. 2019, 42, 492–505. (In Chinese) [Google Scholar]
  22. Zhou, Y.L. Difference schemes with intrinsic parallelism for quasi-linear parabolic systems. Sci. China (Ser. A) 1997, 40, 270–278. [Google Scholar] [CrossRef]
  23. Zhang, B.L.; Yuan, G.X.; Liu, X.P.; Chen, J. Parallel Finite Difference Methods for Partial Differential Equations; Science Press: Beijing, China, 1994. (In Chinese) [Google Scholar]
  24. Petter, B.; Mitchell, L. Parallel Solution of Partial Differential Equations; Springer: New York, NY, USA, 2000. [Google Scholar]
  25. Zhang, B.L.; Gu, T.X.; Mo, Z.Y. Principles and Methods of Numerical Parallel Computation; National Defense Industry Press: Beijing, China, 1999. (In Chinese) [Google Scholar]
  26. Yuan, G.W.; Sheng, Z.Q.; Hang, X.D. Calculation Methods of Diffusion Equations; Science Press: Beijing, China, 2015. (In Chinese) [Google Scholar]
  27. Jiang, Y.L. New Methods for Engineering Mathematics; Higher Education Press: Beijing, China, 2013. (In Chinese) [Google Scholar]
  28. Wang, H.; Wang, K.; Sircar, T. A direct O(Nlog2N) finite difference method for fractional diffusion equations. J. Comput. Phys. 2010, 229, 8095–8104. [Google Scholar] [CrossRef]
  29. Diethelm, K. An efficient parallel algorithm for the numerical solution of fractional differential equations. Fract. Calcul. Appl. Anal. 2011, 14, 475–490. [Google Scholar] [CrossRef]
  30. Wang, H.; Basu, T.S. A fast finite difference method for two-dimensional space-fractional diffusion equations. SIAM J. Sci. Comput. 2012, 34, 2444–2458. [Google Scholar] [CrossRef]
  31. Gong, C.Y.; Bao, W.M.; Tang, G.J. A parallel algorithm for the Riesz fraction reaction-diffusion equation with explicit finite difference method. Fract. Calcul. Appl. Anal. 2013, 16, 654–669. [Google Scholar]
  32. Sweilam, N.H.; Moharram, H.; Moniem, N.K.A.; Ahmed, S. A parallel Crank-Nicolson finite difference method for time-fractional parabolic equation. J. Numer. Math. 2014, 22, 363–382. [Google Scholar] [CrossRef]
  33. Wang, Q.L.; Liu, J.; Gong, C.Y.; Tang, X.; Fu, G.; Xing, Z. An efficient parallel algorithm for Caputo fractional reaction-diffusion equation with implicit finite-difference method. Adv. Differ. Equ. 2016, 207. [Google Scholar] [CrossRef] [Green Version]
  34. Yang, X.Z.; Dang, X. A new parallel difference algorithm based on improved alternating segment Crank-Nicolson scheme for time fractional reaction-diffusion equation. Adv. Differ. Equ. 2019, 417. [Google Scholar] [CrossRef] [Green Version]
  35. Fu, H.F.; Wang, H. A preconditioned fast parareal finite difference method for space-time fractional partial differential equation. J. Sci. Comput. 2018, 78, 1724–1743. [Google Scholar] [CrossRef]
  36. Liu, W. Actual Combat Matlab Parallel Programming; Beihang University Press: Beijing, China, 2012. (In Chinese) [Google Scholar]
  37. Chi, X.B.; Wang, Y.W.; Wang, Y.; Liu, F. Parallel Computing and Implementation Technology; Science Press: Beijing, China, 2015. (In Chinese) [Google Scholar]
Figure 1. Point diagram of the PASE-I scheme.
Figure 1. Point diagram of the PASE-I scheme.
Mathematics 08 00596 g001
Figure 2. The error surfaces of numerical solutions of PASE-I and PASI-E schemes.
Figure 2. The error surfaces of numerical solutions of PASE-I and PASI-E schemes.
Mathematics 08 00596 g002
Figure 3. Distribution of the difference total energy in spatial grid points (M takes 50 , 100 , 200 , 400 ).
Figure 3. Distribution of the difference total energy in spatial grid points (M takes 50 , 100 , 200 , 400 ).
Mathematics 08 00596 g003
Figure 4. Changes in the relative error over time steps (N takes 400 , 800 , 1600 , 3200 ).
Figure 4. Changes in the relative error over time steps (N takes 400 , 800 , 1600 , 3200 ).
Mathematics 08 00596 g004
Figure 5. Comparison of computational efficiency of three schemes.
Figure 5. Comparison of computational efficiency of three schemes.
Mathematics 08 00596 g005
Table 1. E and O r d e r x of PASE-I and PASI-E schemes ( τ α h 2 ) .
Table 1. E and O r d e r x of PASE-I and PASI-E schemes ( τ α h 2 ) .
α 1 α 2 MNPASE-I SchemePASI-E Scheme
E Order x E Order x
0.40.4501321.960988 × 10 3 ——1.933065 × 10 3 ——
1003164.750145 × 10 4 2.0455374.715343 × 10 4 2.035455
2007521.175639 × 10 4 2.0145261.171268 × 10 4 2.009290
40017882.927383 × 10 5 2.0057612.921901 × 10 5 2.003092
80042547.293185 × 10 6 2.0049907.286329 × 10 6 2.003643
0.50.5501841.808959 × 10 3 ——1.783177 × 10 3 ——
1004644.441176 × 10 4 2.0261464.408679 × 10 4 2.016032
20011691.093836 × 10 4 2.0215441.097908 × 10 4 2.005589
40029472.716023 × 10 5 2.0098292.721104 × 10 5 2.012493
80074266.763548 × 10 6 2.0056436.757212 × 10 6 2.009691
0.60.6502671.344353 × 10 4 ——1.345261 × 10 4 ——
1007193.270205 × 10 5 2.0394593.271047 × 10 5 2.040061
20019378.019243 × 10 6 2.0278438.020099 × 10 6 2.028060
40052141.979933 × 10 6 2.0180141.979839 × 10 6 2.018236
800140364.901814 × 10 7 2.0140644.901704 × 10 7 2.014027
Table 2. E and O r d e r t of PASE-I and PASI-E schemes ( M = 200 ).
Table 2. E and O r d e r t of PASE-I and PASI-E schemes ( M = 200 ).
α 1 α 2 NPASE-I SchemePASI-E Scheme
E Order t E Order t
0.350.352001.540090 × 10 4 ——1.541745 × 10 4 ——
4004.981818 × 10 5 1.6282704.986837 × 10 5 1.628367
8001.579056 × 10 5 1.6576091.580609 × 10 5 1.657643
16004.971665 × 10 6 1.6672624.976502 × 10 6 1.667277
32001.561734 × 10 6 1.6705791.563244 × 10 6 1.670588
0.40.22001.641599 × 10 4 ——1.643375 × 10 4 ——
4005.320684 × 10 5 1.6254185.326029 × 10 5 1.625530
8001.681644 × 10 5 1.6617391.683288 × 10 5 1.661777
16005.266880 × 10 6 1.6748515.271969 × 10 6 1.674868
32001.643123 × 10 6 1.6805081.644700 × 10 6 1.680517
0.80.22006.070178 × 10 4 ——6.075953 × 10 4 ——
4002.349702 × 10 4 1.3692602.351861 × 10 4 1.369307
8009.139187 × 10 5 1.3623409.147204 × 10 5 1.362400
16003.598744 × 10 5 1.3445723.601678 × 10 5 1.344661
32001.445759 × 10 5 1.3156661.446810 × 10 5 1.315793
Table 3. Comparison of the three difference schemes’s CPU time (Unit: second).
Table 3. Comparison of the three difference schemes’s CPU time (Unit: second).
10050010002000300040005000
Implicit9.9043253.6988139.582359.223499.958900.6281462.45
PASE-I9.1451644.590991.2816171.389272.927368.622493.003
PASI-E10.567244.705291.4469171.517262.139368.364491.059
S P of PASE-I1.083011.204251.529142.095941.831832.443222.96641
S P of PASI-E0.937261.201171.526372.094371.907222.444942.97815

Share and Cite

MDPI and ACS Style

Yang, X.; Wu, L. A New Kind of Parallel Natural Difference Method for Multi-Term Time Fractional Diffusion Model. Mathematics 2020, 8, 596. https://doi.org/10.3390/math8040596

AMA Style

Yang X, Wu L. A New Kind of Parallel Natural Difference Method for Multi-Term Time Fractional Diffusion Model. Mathematics. 2020; 8(4):596. https://doi.org/10.3390/math8040596

Chicago/Turabian Style

Yang, Xiaozhong, and Lifei Wu. 2020. "A New Kind of Parallel Natural Difference Method for Multi-Term Time Fractional Diffusion Model" Mathematics 8, no. 4: 596. https://doi.org/10.3390/math8040596

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop