Next Article in Journal
Error Bounds of a Finite Difference/Spectral Method for the Generalized Time Fractional Cable Equation
Next Article in Special Issue
Non-Local Seismo-Dynamics: A Fractional Approach
Previous Article in Journal
On Sharp Estimate of Third Hankel Determinant for a Subclass of Starlike Functions
Previous Article in Special Issue
Strong Solution for Fractional Mean Field Games with Non-Separable Hamiltonians
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Fast Compact Difference Scheme for Solving the Two-Dimensional Time-Fractional Cattaneo Equation

1
College of Mathematics and Statistic, Guangxi Normal University, Guilin 541004, China
2
College of Science, Lanzhou University of Technology, Lanzhou 730050, China
3
College of Science, Guilin University of Technology, Guilin 541004, China
*
Author to whom correspondence should be addressed.
Fractal Fract. 2022, 6(8), 438; https://doi.org/10.3390/fractalfract6080438
Submission received: 26 June 2022 / Revised: 8 August 2022 / Accepted: 9 August 2022 / Published: 11 August 2022

Abstract

:
The time-fractional Cattaneo equation is an equation where the fractional order α ( 1 , 2 ) has the capacity to model the anomalous dynamics of physical diffusion processes. In this paper, we consider an efficient scheme for solving such an equation in two space dimensions. First, we obtain the space’s semi-discrete numerical scheme by using the compact difference operator in the spatial direction. Then, the semi-discrete scheme is converted to a low-order system by means of order reduction, and the fully discrete compact difference scheme is presented by applying the L2-1 σ formula. To improve the computational efficiency, we adopt the fast discrete Sine transform and sum-of-exponentials techniques for the compact difference operator and L2-1 σ difference operator, respectively, and derive the improved scheme with fast computations in both time and space. That aside, we also consider the graded meshes in the time direction to efficiently handle the weak singularity of the solution at the initial time. The stability and convergence of the numerical scheme under the uniform meshes are rigorously proven, and it is shown that the scheme has second-order and fourth-order accuracy in time and in space, respectively. Finally, numerical examples with high-dimensional problems are demonstrated to verify the accuracy and computational efficiency of the derived scheme.

1. Introduction

Let Ω R 2 be a bounded rectangle domain and T be the fixed time. In this paper, we study the efficient numerical scheme for solving the following time-fractional Cattaneo equation:
t u ( x , t ) + κ C D 0 , t α u ( x , t ) = μ Δ u ( x , t ) + f ( x , t ) , in Ω × ( 0 , T ] ,
where the boundary and initial conditions are u ( x , t ) = 0 on Ω × ( 0 , T ] and u ( x , 0 ) = u 0 ( x ) , t u ( x , 0 ) = u 1 ( x ) in Ω . Here, x = ( x 1 , x 2 ) .The parameters κ and μ are the positive constants related to the relaxation time and diffusion properties, respectively. The source term f and the initial conditions u 0 and u 1 are given suitable functions. The Laplacian Δ is defined by Δ = x 1 2 + x 2 2 . The Caputo derivative C D 0 , t α with α ( 1 , 2 ] is given by
C D 0 , t α v ( t ) = 1 Γ ( 2 α ) 0 t ( t s ) 1 α v ( s ) d s ,
where Γ ( · ) is the Gamma function.
The time-fractional Cattaneo Equation (1) provides a flexible and efficient way to model the anomalous dynamics of physical diffusion processes, especially the general dynamics crossover behaviors [1,2]. In recent years, many efforts have been made in the theoretical and numerical study of such a time-fractional equation or its variants. See [2,3] for the theoretical case and [1,4,5,6] for the numerical case, just to name a few. One may refer to these two review papers [7,8] or this book [9] for further details.
In [4], Zhao and Sun proposed a compact Crank–Nicolson scheme for solving the time-fractional Cattaneo Equation (1). Ren and Gao also constructed an alternating direction implicit (ADI) compact difference scheme by adding a small term to improve the computational efficiency for the two-dimensional problem [5]. The temporal convergence order of both of the above two methods is the order 3 α , where 1 < α < 2 . Later on, Chen and Li proposed a temporal second-order ADI Galerkin scheme by applying the Galerkin finite element method in space and the L2-1 σ method in time [6]. In [1], Chen and Nong developed two efficient, fully discrete schemes for solving Equation (1) based on the Galerkin finite method in space and the convolution quadrature in time. Error estimates and numerical experiments show that their schemes are efficient when dealing with different data regularity in Equation (1). It seems that the schemes mentioned above still have some limitations, especially in the high-dimensional and non-smooth solution problems. This motivated us to develop an efficient numerical scheme to cope with such situations.
For the high-dimensional time-fractional model, there are two ways to improve the computational efficiency of the numerical scheme. One involves the spatial direction, such as with the ADI technique, fast discrete Sine transforms (DSTs), and so on (see [6,10]). The other one involves the temporal direction, such as with the fast convolution quadrature and sum-of-exponentials (SOE) techniques (see [11,12,13,14]). In this paper, we focus on the DST and SOE techniques.
The DST technique has attracted the interest of many scholars in recent years, mainly because it avoids the direct inversion for the coefficient matrix of the linear discrete systems and, at the same time, ensures the convergence accuracy of the difference operators [15,16]. In [15], Wang, et al. proposed an efficient fourth-order compact finite difference scheme for solving the Poisson equation based on the fast discrete Sine transform and greatly reduced the computational cost by avoiding matrix inversion for the discretized system. There has been some work in recent years on the application of the DST technique in numerically solving high-dimensional time-fractional equations (see [10,17,18]).
The SOE technique is an efficient way to reduce the computational complexity caused by the non-locality of fractional derivatives [13,19]. In [13], Yan et al. presented the fast L2-1 σ (denoted by FL2-1 σ ) formula by employing the SOE technique with the kernel function in the Caputo derivative. Subsequently, some numerical studies of time-fractional models based on this approach have emerged (see [10,19] and the corresponding references therein). In [19], Liang et al. proposed a fast difference scheme for solving the one-dimensional time-fractional telegraph equation based on the FL2-1 σ formula. In [10], Li et al. derived a fast compact difference scheme for solving subdiffusion equations by applying the FL2-1 σ in time and the compact difference operator in space. It is worth mentioning that the compact difference operator is implemented by the DST via the fast Fourier transform (FFT) algorithm. Therefore, their scheme is computationally efficient in both the time direction and spatial dimensions.
Motivated by the work in [10], we aim to develop a fast compact difference scheme for solving Equation (1). To this end, we derive a compact difference scheme by applying the L2-1 σ formula in time and compact difference operator in space. Then, the SOE and DST techniques are both adopted to improve the computational performance of the derived scheme. The contributions of our paper are as follows. First, we construct a compact difference scheme with fast computation in both the temporal and spatial directions, which effectively handles the high computational cost caused by the high spatial dimension and the non-locality of the Caputo derivative (see the numerical scheme (10)). Second, we present the rigorous stability and error estimate for the uniform mesh-based fast compact difference scheme (see Theorems 1 and 2). Third, extensive numerical experiments are designed to verify the accuracy and efficiency of the fast compact difference scheme (see the Section 5, Table 1 and Table 2 for the accuracy; Table 5 and Figure 1 for efficiency).Besides, the fast compact difference scheme based on graded meshes in time is adopted to handle non-smooth solution problems (see Table 3 and Table 4 in Section 5).
This paper is organized as follows. In Section 2, we present the fully discrete difference scheme based on the compact difference operator for the Laplacian in space and the L2-1 σ formula for the Caputo derivative in time. In order to improve the computational performance of the derived scheme, we further adopt the DST and SOE techniques to obtain the fast compact difference scheme in Section 3. In Section 4, the stability and error estimate of the fast compact difference scheme based on the uniform meshes are provided rigorously. Numerical examples and the conclusions of this paper are given in Section 5 and Section 6, respectively.

2. The Compact Difference Scheme

In this section, we first approximate the Laplacian Δ in Equation (1) by applying the compact difference operator and then find the low-order system for the space’s semi-discrete scheme with the aid of the reduced-order method so that the fully discrete scheme of the equation can be obtained by using the second-order L2-1 σ formula.
To this end, we introduce some useful notations in the following. Suppose the rectangle domain Ω = [ x 1 L , x 1 R ] × [ x 2 L , x 2 R ] . Let the spatial step size be h k = ( x k R x k L ) / M k and the grid points be x k , j k = x k L + j k h k ( j k = 0 , 1 , , M k ) , where M k ( k = 1 , 2 ) is a positive integer. The fully discrete grids on Ω are given by Ω ¯ h = { x h = ( x 1 , j 1 , x 2 , j 2 ) | 0 j k M k ( k = 1 , 2 ) } . We will use M without a subscript to denote M 1 (or M 2 ) when M 1 = M 2 , unless otherwise specified. The inner and boundary grids are denoted by Ω h = Ω ¯ h Ω and Ω h = Ω ¯ h Ω , respectively. We denote the space of the grid function as V h = { v | v = ( v h ) x h and v h = 0 for x h Ω h } .
For the grid function v h = v ( x h ) with the index vector h = ( i 1 , i 2 ) at the kth position ( k = 1 , 2 ) , the compact difference operator is given by Δ ¯ k v i k = δ k 2 H k v i k , where
H k v i k = I + h k 2 12 δ k 2 v i k .
Here, we have
δ k 2 v i k = δ k v i k + 1 2 δ k v i k 1 2 h k and δ k v i k + 1 2 = v i k + 1 v i k h k .
Therefore, for x h Ω h , we have Δ v ( x h ) = Δ ¯ h v h + O ( h 4 ) with Δ ¯ h v h = ( Δ ¯ 1 + Δ ¯ 2 ) v h and h = ( h 1 , h 2 ) .
When employing the compact difference operator Δ ¯ h for the discretization of the Laplacian Δ in (1), one has
t u ( x h , t ) + κ C D 0 , t α u ( x h , t ) = μ Δ ¯ h u ( x h , t ) + f ( x h , t ) + R x ( t ) ,
from which we obtain the following space semi-discrete scheme:
t u h ( t ) + κ C D 0 , t α u h ( t ) = μ Δ ¯ h u h ( t ) + f h ( t ) ,
where u h ( t ) u ( x h , t ) and the truncation error R x ( t ) = O ( h 4 ) .
Next, we introduce the L2-1 σ formula for the approximation of the Caputo derivative C D 0 , t β with β ( 0 , 1 ) . Let the temporal step size be τ = T / N t . Here, N t is a positive integer. Denote t n = n τ , where n 0 and t n + 1 2 = ( t n + t n + 1 ) / 2 . For a given function g ( t ) , denote
δ t g n + 1 2 = g ( t n + 1 ) g ( t n ) τ , g n + 1 2 = g ( t n + 1 ) + g ( t n ) 2 , δ ^ t g n = 1 2 τ ( ( 2 σ + 1 ) g ( t n + 1 ) 4 σ g ( t n ) + ( 2 σ 1 ) g ( t n 1 ) ) ,
and g n + σ = σ g ( t n + 1 ) + ( 1 σ ) g ( t n ) , where σ = 1 β / 2 . For more details about the difference operator δ ^ t , one may refer to Lemma 2.1 of [20].
We state the L2- 1 σ formula as follows [21]. If g ( t ) C 3 [ 0 , T ] , then
C D 0 , t β g ( t ) | t = t n + σ = D ¯ τ β g n + σ + O ( τ 3 β ) ,
where D ¯ τ β g n + σ = k = 0 n w n k ( n + 1 ) ( g ( t k + 1 ) g ( t k ) ) and for n = 0 , w 0 ( 1 ) = a 0 τ β / Γ ( 2 β ) , where n 1 , we have
w k ( n + 1 ) = τ β Γ ( 2 β ) a 0 + b 1 , k = 0 , a k + b k + 1 b k , k = 1 , 2 , , n 1 , a n b n , k = n .
Here, a 0 = σ 1 β , a k = ( k + σ ) 1 β ( k 1 + σ ) 1 β , k 1 , and
b k = ( k + σ ) 2 β ( k 1 + σ ) 2 β 2 β ( k + σ ) 1 β ( k 1 + σ ) 1 β 2 , k 1 .
By setting ϕ h ( t ) = t u h ( t ) for (2), we find the following lower-order system:
ϕ h ( t ) + κ C D 0 , t α 1 ϕ h ( t ) = μ Δ ¯ h u h ( t ) + f h ( t ) .
Therefore, by replacing C D 0 , t α 1 with the difference operator D ¯ τ α 1 of the L2-1 σ Formula (3) in the space semi-discrete scheme (2), one has
ϕ h ( t n + σ ) + κ D ¯ τ α 1 ϕ h n + σ = μ Δ ¯ h u h ( t n + σ ) + f h ( t n + σ ) + R x t n , δ t u h 1 / 2 = ϕ h ( t 1 / 2 ) + r 0 , δ ^ t u h n = ϕ h ( t n + σ ) + r n ,
where the truncation error R x t n = O ( τ 2 + h 4 ) and r 0 = r n = O ( τ 2 ) . By dropping the small terms R x t n , r 0 , and r n , we have the following compact difference scheme for solving Equation (1). We find the numerical solution ϕ h n of ϕ ( x h , t n ) for n 1 such that
ϕ h n + σ + κ D ¯ τ α 1 ϕ h n + σ = μ Δ ¯ h u h n + σ + f h n + σ , δ t u h 1 / 2 = ϕ h 1 / 2 , δ ^ t u h n = ϕ h n + σ ,
where u h 0 = u 0 ( x h ) , ϕ h 0 = u 1 ( x h ) , and u ( x h ) | x h Ω h = 0 . From (6), we obtain the numerical solution u h n of u ( x h , t n ) to (1): u h 1 = u h 0 + τ ( ϕ h 1 + ϕ h 0 ) / 2 and for n 1 :
u h n + 1 = 2 τ 2 σ + 1 ( σ ϕ h n + 1 + ( 1 σ ) ϕ h n ) + 1 2 σ + 1 ( 4 σ u h n ( 2 σ 1 ) u h n 1 ) .

3. The Derivation of the Fast Compact Difference Scheme

In this part, we use the DST and SOE techniques to improve the computational performance of the compact difference scheme (6).
By utilizing the discrete sine transform, one can derive that
v j k ^ v j k ^ 12 h k 2 · s j k 1 s j k + 5 : = v j k ^ λ j k , M k ,
where s j k = cos ( j k π M k ) and 1 j k M k 1 ( k = 1 , 2 ) (see [15] for more details about the derivation).
Therefore, the scheme in Equation (6) is equivalent to
ϕ ^ ν n + σ + κ D ¯ τ α 1 ϕ ^ ν n + σ = μ λ j 1 , M 1 + λ j 2 , M 2 u ^ ν n + σ + f ^ ν n + σ , δ t u ^ ν 1 / 2 = ϕ ^ ν 1 / 2 , δ ^ t u ^ ν n = ϕ ^ ν n + σ ,
where ν is the index set defined by ν = { ( j 1 , j 2 ) | 1 j k M k 1 ( k = 1 , 2 ) } . We can obtain the numerical solution to Equation (1) by the following steps:
(a)
Compute ϕ ^ ν 0 , u ^ ν 0 , and f ^ ν n from ϕ h 0 , u h 0 , and f h n by means of DST, where n 1 ;
(b)
Compute ϕ ^ ν n by solving the system (7), and then obtain u ^ ν n , where n 1 ;
(c)
Obtain the numerical solutions ϕ h n and u h n from ϕ ^ ν n and u ^ ν n , respectively, by means of the inverse of DST, where n 1 .
It is worth mentioning that if we solve numerical scheme (6) directly, the inverse of the coefficient matrix needs to be solved at each time layer, and this leads to a computational cost of O ( M 1 2 M 2 2 ) . If the DST technique based on FFT, that is the scheme (7) is used, then the computational cost will be reduced to O ( M 1 M 2 log ( M 1 M 2 ) ) .
Next, we consider the SOE technique to further reduce the computational complexity of the scheme (7) caused by the non-locality of the Caputo derivative.
Without loss of generality, we adopt the graded meshes t n = T ( n N t ) γ in time, where n = 0 , , N t and γ 1 is a proper chosen temporal mesh grading parameter. Denote the time step size τ n = t n t n 1 . Set t n + σ = t n + σ τ n + 1 and D τ g k = g k g k 1 τ k . The fast L2-1 σ formula based on the graded meshes is defined as
F D ¯ τ n β g n + σ = ϖ 0 D τ g n + 1 + k = 1 N exp ϖ k e s k σ τ n + 1 H k ( t n ) ,
where ϖ 0 = ( σ τ n + 1 ) 1 β Γ ( 2 β ) and H k ( t j ) ( j 0 , k 1 ) is obtained by the following recurrence:
H k ( t j ) = e s k τ j H k ( t j 1 ) + a k , j D τ g j + b k , j ( D τ g j + 1 D τ g j ) ,
with the initial value H k ( t 0 ) = 0 [10]. Here, the coefficients a k , j = 1 s k ( 1 e s k τ j ) and b k , j = 2 τ j + τ j + 1 t j 1 t j ( t t j 1 2 ) e s k ( t j t ) d t . The weights ϖ k ( k 1 ) and the points s k ( k 1 ) in Formula (8) are chosen to satisfy
t β Γ ( 1 β ) k = 1 N exp ϖ k e s k t ϵ , t [ t 1 , T ] ,
where ϵ is the absolute tolerance error and the number of exponentials N exp satisfies
N exp = O log 1 ϵ log log 1 ϵ + log T τ 1 + log 1 τ 1 log log 1 ϵ + log 1 τ 1 .
We remark that the letter F in the difference operator F D ¯ τ n β (see the left-hand side of the Formula (8)) represents the word “fast”, which refers to the fact that the corresponding formula uses the SOE technique to improve the computational performance. It is shown in [13] that the number of exponentials N exp = O ( log N t ) when T 1 and N exp = O log 2 N t when T 1 . The total computational cost for the fast L2-1 σ Formula (8) is O ( N t N exp ) with storage O ( N exp ) .
By replacing the classical L2-1 σ formula in the numerical scheme (7) with the above fast L2-1 σ formula based on graded meshes, we have the following fast compact difference scheme:
ϕ ^ ν n + σ + κ F D ¯ τ n α 1 ϕ ^ ν n + σ = μ λ j 1 , M 1 + λ j 2 , M 2 u ^ ν n + σ + f ^ ν n + σ , δ ˜ t u ^ ν 1 / 2 = ϕ ^ ν 1 / 2 , δ ^ ˜ t u ^ ν n = ϕ ^ ν n + σ ,
where δ ˜ t u ^ ν 1 / 2 = u ^ ν 1 u ^ ν 0 τ 1 and
δ ^ ˜ t u ^ ν n = τ n + 2 σ τ n + 1 τ n + τ n + 1 D τ u ^ ν n + 1 τ n + 1 ( 2 σ 1 ) τ n + τ n + 1 D τ u ^ ν n .
We will simply denote the fast compact difference scheme (10) as the FL2-1 σ ( γ ) scheme in what follows if no ambiguity arises.
One can observe that, compared with the numerical scheme (6), the FL2-1 σ ( γ ) scheme (10) reduces the overall computational cost from O ( N t 2 M 1 2 M 2 2 ) to O ( N t N exp M 1 M 2 log ( M 1 M 2 ) ) , which greatly improves the computational efficiency of the original scheme (6).
When γ = 1 , the graded meshes recover the uniform one. Thus, in this case, we shall use τ instead of τ n in (10) unless otherwise specified. The corresponding scheme based on uniform meshes is called the FL2-1 σ (1) scheme. When γ > 1 , the grid points in the time direction are concentrated near the origin, which is beneficial for the numerical solution to capture the weak singularity solution, thus improving the accuracy loss of the scheme caused by the non-smooth solution. This will be verified by numerical examples in Section 5.

4. Stability and Error Estimates

We studied the stability and error estimates for the FL2-1 σ (1) scheme (10) in this section.
For any given grid function v V h , the discrete L 2 norm is given by v = ( v , v ) h with the discrete inner product ( u , v ) h = ( h 1 h 2 ) x h Ω h u h v h . The discrete H 1 semi-norm and H 1 norm are denoted as | v | 1 = h v h 2 = δ 1 v h 2 + δ 2 v h 2 and v 1 = v 2 + | v | 1 2 , respectively. Here, h = ( δ 1 , δ 2 ) . In light of the embedding theorem, we can conclude that for any given v V h , the discrete H 1 semi-norm | v | 1 and H 1 norm v 1 are equivalent.
Using the recursive relationship of the historical term H k ( t j ) in Equation (9), one may reformulate the expression of the difference operator F D ¯ τ β based on the uniform meshes as follows (or see [10,13] more details on the derivation):
F D ¯ τ β g n + σ = k = 0 n W n k ( n + 1 ) ( g ( t k + 1 ) g ( t k ) ) ,
where W n k ( n + 1 ) are defined as W 0 ( 1 ) = w 0 ( 1 ) when n = 0 , and for n 1 , we have
W k ( n + 1 ) = A 0 + B 1 , k = 0 , A k + B k + 1 B k , k = 1 , 2 , , n 1 , A n B n , k = n .
Here, A 0 = W 0 ( 1 ) , A n k = 1 τ t k t k + 1 j = 0 N exp ϖ j e s j ( t n + σ s ) d s , and
B n k = 1 τ 2 t k t k + 1 j = 0 N exp ϖ j e s j ( t n + σ s ) ( s t k + 1 2 ) d s .
We need the following lemma:
Lemma 1.
The weights W k ( n + 1 ) in Equation (12) have the following properties:
k = 1 m 1 ( W k 1 ( k + 1 ) W k ( k + 1 ) ) 9 ϵ ( m 1 ) 4 Γ ( 1 β ) + 4 3 β τ β Γ ( 3 β ) ( m 1 + σ ) 1 β ,
and
k = 1 m 1 W k ( k + 1 ) ϵ ( m 1 ) Γ ( 1 β ) + ( m 1 + σ ) 1 β τ β Γ ( 2 β ) .
Proof. 
In light of Lemma 3.3 in [19] and the definitions of W k ( n + 1 ) in Equation (12), one may readily obtain the desired results. □
Let the symbol c be a positive constant which is independent of the temporal and spatial step sizes. The stability of the scheme (10) is described as follows:
Theorem 1.
The FL2-1 σ ( γ ) compact difference scheme (10) with γ = 1 is stable in the sense that
u h m 1 2 c ( ϕ h 0 2 + h u h 0 2 + τ n = 1 m 1 h ε h n + σ 2 + τ h ε h 1 2 2 + τ n = 0 m 1 f h n + σ 2 ) ,
where ε h n + σ and ε h 1 2 are the perturbations of δ ^ t u h n = ϕ h n + σ + ε h n + σ and δ t u h 1 2 = ϕ h 1 2 + ε h 1 2 , respectively.
Proof. 
In light of the equivalence of schemes (6) and (7), one only needs to consider the following scheme instead of (10):
ϕ h n + σ + κ F D ¯ τ α 1 ϕ h n + σ = μ Δ ¯ h u h n + σ + f h n + σ , δ t u h 1 / 2 = ϕ h 1 / 2 , δ ^ t u h n = ϕ h n + σ .
We first consider the case n 1 for (13). By taking the discrete inner product on both sides of (13) with φ h n + σ , we have
( ϕ h n + σ , ϕ h n + σ ) h + ( κ F D ¯ τ α 1 ϕ h n + σ , ϕ h n + σ ) h = ( μ Δ ¯ h u h n + σ , ϕ h n + σ ) h + ( f h n + σ , ϕ h n + σ ) h .
We can use Lemma 1 and the estimate for the discrete operator Δ ¯ h presented in (3.14) of [10] such that
3 2 Δ h u h n + σ , ϕ h n + σ h < Δ ¯ h u h n + σ , ϕ h n + σ h < Δ h u h n + σ , ϕ h n + σ h ,
From this, we derive that
ϕ h n + σ 2 + κ 2 F D ¯ τ α 1 ϕ h n + σ 2 μ ( h u h n + σ , h ϕ h n + σ ) h + ( f h n + σ , ϕ h n + σ ) h .
Note that δ ^ t u h n = ϕ h n + σ + ε h n + σ , so it follows from Lemma 3.5 of [20] and the Cauchy–Schwarz inequality that
μ ( h u h n + σ , h ϕ h n + σ ) h = μ ( h u h n + σ , h δ ^ t u h n ) h + μ ( h u h n + σ , h ε h n + σ ) h μ 4 τ ( E n + 1 E n ) + μ h u h n + σ h ε h n + σ ,
where E n + 1 = ( 2 σ + 1 ) h u h n + 1 2 ( 2 σ 1 ) h u h n 2 + 2 σ 2 + σ 1 h u h n + 1 h u h n 2 .
By using the inequality a b ϵ a 2 + 1 4 ϵ b 2 with ϵ > 0 , we obtain
h u h n + σ h ε h n + σ 1 4 h u h n + σ 2 + h ε h n + σ 2 .
Therefore, the inequality (15) has the following estimate:
ϕ h n + σ 2 + κ 2 F D ¯ τ α 1 ϕ h n + σ 2 + μ 4 τ E n + 1 μ 4 τ E n + μ 4 h u h n + σ 2 + μ h ε h n + σ 2 + ( f h n + σ , ϕ h n + σ ) h ,
In other words, we have
ϕ h n + σ 2 + κ 2 W 0 ( n + 1 ) ϕ h n + 1 2 + κ 2 k = 1 n ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 + μ 4 τ E n + 1 κ 2 W n ( n + 1 ) ϕ h 0 2 + μ 4 τ E n + μ 4 h u h n + σ 2 + μ h ε h n + σ 2 + ( f h n + σ , ϕ h n + σ ) h .
We sum up n from n = 1 to n = m 1 for the above inequality and obtain
n = 1 m 1 ϕ h n + σ 2 + κ 2 n = 1 m 1 W 0 ( n + 1 ) ϕ h n + 1 2 + κ 2 n = 1 m 1 k = 1 n ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 + μ 4 τ E m κ 2 n = 1 m 1 W n ( n + 1 ) ϕ h 0 2 + μ 4 τ E 1 + μ 4 n = 1 m 1 h u h n + σ 2 + μ n = 1 m 1 h ε h n + σ 2 + n = 1 m 1 ( f h n + σ , ϕ h n + σ ) h .
Since
n = 1 m 1 k = 1 n ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 = n = 1 m 1 k = 2 n ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 + ( W n ( n + 1 ) W n 1 ( n + 1 ) ) ϕ h 1 2 = k = 2 m 1 n = k m 1 ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 + n = 1 m 1 ( W n ( n + 1 ) W n 1 ( n + 1 ) ) ϕ h 1 2 = k = 2 m 1 ( W m k ( m ) W 0 ( k + 1 ) ) ϕ h k 2 + n = 1 m 1 ( W n ( n + 1 ) W n 1 ( n + 1 ) ) ϕ h 1 2 ,
we have
n = 1 m 1 W 0 ( n + 1 ) ϕ h n + 1 2 + n = 1 m 1 k = 1 n ( W n k + 1 ( n + 1 ) W n k ( n + 1 ) ) ϕ h k 2 = k = 2 m W m k ( m ) ϕ h k 2 + n = 1 m 1 ( W n ( n + 1 ) W n 1 ( n + 1 ) ) ϕ h 1 2 C f k = 2 m ϕ h k 2 9 ϵ ( m 1 ) 4 Γ ( 1 β ) + 4 3 β τ β Γ ( 3 β ) ( m 1 + σ ) 1 β ϕ h 1 2 ,
where Lemma 1 and the boundness of W m k ( m ) presented in Lemma 4.1 of [13] are used. Here, C f = min { W 0 ( 1 ) , W n ( n + 1 ) } and β = α 1 .
The term E m has the estimate E m 1 σ h u h m 2 (see Lemma 3.5 of [20]). Thus, one has the further estimate for the inequality (16):
n = 1 m 1 ϕ h n + σ 2 + κ C f 2 k = 2 m ϕ h k 2 + μ σ 1 4 τ h u h m 2 9 ϵ ( m 1 ) 4 Γ ( 1 β ) + 4 3 β τ β Γ ( 3 β ) ( m 1 + σ ) 1 β ϕ h 1 2 + ϵ ( m 1 ) Γ ( 1 β ) + ( m 1 + σ ) 1 β τ β Γ ( 2 β ) ϕ h 0 2 + μ 4 τ E 1 + μ 4 n = 1 m 1 h u h n + σ 2 + μ n = 1 m 1 h ε h n + σ 2 + n = 1 m 1 ( f h n + σ , ϕ h n + σ ) h .
Since
n = 1 m 1 ( f h n + σ , ϕ h n + σ ) h 1 2 n = 1 m 1 1 4 ϵ 1 f h n + σ 2 + ϵ 1 ϕ n + σ 2 + 1 4 ϵ 2 f h n + σ 2 + ϵ 2 ϕ n + σ 2 c n = 1 m 1 1 ϵ 1 f h n + σ 2 + ϵ 1 n = 1 m 1 ϕ n + σ 2 + ϵ 2 n = 2 m ϕ n 2 + ϕ 1 2 ,
by choosing suitable ϵ 1 and ϵ 2 values, the inequality (17) yields
h u h m 2 c ( ϕ h 0 2 + ϕ h 1 2 + h u h 0 2 + h u h 1 2 + τ n = 1 m 1 h u h n 2 + τ n = 1 m 1 h ε h n + σ 2 + τ n = 1 m 1 f h n + σ 2 ) .
Next, we consider the bound of ϕ h 1 2 + h u h 1 2 on the right-hand side of (18). To this end, we operate the discrete inner product on both sides of the numerical scheme (6) with ϕ h 1 2 for n = 0 and find
( ϕ h σ , ϕ h 1 2 ) h + κ ( F D ¯ τ α 1 ϕ h σ , φ h 1 2 ) h = μ ( Δ ¯ h u h σ , ϕ h 1 2 ) h + ( f h σ , ϕ h 1 2 ) h .
Observing that F D ¯ τ α 1 ϕ h σ = W 0 ( 1 ) ( ϕ h 1 ϕ h 0 ) , and using the equality δ t u h 1 2 = ϕ h 1 2 + ε h 1 2 , we obtain
( ϕ h σ , ϕ h 1 2 ) h + κ W 0 1 2 ( ϕ h 1 2 ϕ h 0 2 ) = μ ( h u h σ , h δ t u h 1 2 ) h + μ ( h u h σ , h ε h 1 2 ) h + ( f h σ , ϕ h 1 2 ) h .
Following the idea in the proof of Theorem 4.1 of [6] for the case n = 0 , we further have
( ϕ h σ , ϕ h 1 2 ) h = 1 2 ( σ ϕ h 1 + ( 1 σ ) ϕ h 0 , ϕ h 1 + ϕ h 0 ) h = 1 2 σ ϕ h 1 2 + ( 1 σ ) ϕ h 0 2 + ( ϕ h 0 , ϕ h 1 ) h 1 2 σ ϕ h 1 2 + ( 1 σ ) ϕ h 0 2 1 2 ϵ 1 ϕ h 1 2 + 1 4 ϵ 1 ϕ h 0 2 ) = 1 4 ϕ h 1 2 + 1 2 1 σ 1 4 ( σ 0.5 ) ϕ h 0 2 ,
where we choose ϵ 1 = σ 1 2 > 0 .
For the first term on the right-hand side of Equation (19), we have
μ ( h u h σ , h δ t u h 1 2 ) h = μ τ ( σ h u h 1 + ( 1 σ ) h u h 0 , h u h 1 h u h 0 ) h = μ τ ( σ h u h 1 2 ( 1 σ ) h u h 0 2 ) + μ τ ( 2 σ 1 ) ( h u h 0 , h u h 1 ) h μ τ ( σ h u h 1 2 ( 1 σ ) h u h 0 2 ) + μ τ ( 2 σ 1 ) ϵ 2 h u h 1 2 + 1 4 ϵ 2 h u h 0 2 = 1 2 τ h u h 1 2 + μ τ 1 σ + ( 2 σ 1 ) 2 4 μ ( σ 0.5 ) h u h 0 2 ,
where we let ϵ 2 = σ 0.5 2 σ 1 > 0 .
The second term on the right-hand side of Equation (19) has the following estimate:
μ ( h u h σ , h ε h 1 2 ) h = μ ( σ h u h 1 + ( 1 σ ) h u h 0 , h ε h 1 2 ) h μ ( σ h u h 1 2 + ( 1 σ ) h u h 0 2 ) h ε h 1 2 2 μ σ ϵ 3 h u h 1 2 + μ σ 1 4 ϵ 3 h ε h 1 2 2 + μ ( 1 σ ) ϵ 4 h u h 0 2 + μ ( 1 σ ) 1 4 ϵ 4 h ε h 1 2 2 .
By suitably choosing the parameters ϵ 3 and ϵ 4 ( e.g., let ϵ 3 = 1 σ ( 1 2 τ 1 ) > 0 and ϵ 4 = 1 σ 4 > 0 ), and noting that W 0 ( 1 ) = w 0 ( 1 ) = t σ 1 β τ Γ ( 2 β ) , we can derive that
ϕ h 1 2 + h u h 0 2 c ϕ h 0 2 + h u h 0 2 + τ h u h 0 2 + τ f h σ 2 + τ h ε h 1 2 2 .
This together with (18) leads to
h u h m 2 c ( ϕ h 0 2 + h u h 0 2 + τ n = 0 m 1 h u h n 2 + τ n = 1 m 1 h ε h n + σ 2 + τ h ε h 1 2 2 + τ n = 0 m 1 f h n + σ 2 ) .
By applying the Gronwall inequality to the above estimate and the equivalence of the H 1 norm and H 1 semi-norm, one has the desired result. All of this ends the proof. □
Finally, we present the convergence result for the FL2-1 σ ( 1 ) scheme (10). Let e h n = u ( x h , t n ) u h n and e ¯ h n = ϕ h ( x h , t n ) ϕ h n . From Lemma 2.3 of [13], we have
C D 0 , t β g ( t ) | t = t n + σ = F D ¯ τ β g n + σ + O ( τ 3 β + ϵ ) .
Therefore, in light of (2), (6), and (13), one may find the following error equations:
e ¯ h n + σ + κ F D ¯ τ α 1 e ¯ h n + σ = μ Δ ¯ h e h n + σ + R ˜ h n + σ , δ t e h 1 / 2 = e ¯ h 1 / 2 + r ˜ h 1 / 2 , δ ^ t e h n = e ¯ h n + σ + r ˜ h n + σ ,
where R ˜ h n + σ = O ( τ 2 + h 4 + ϵ ) , r ˜ h 1 / 2 = O ( τ 2 ) , and r ˜ h n + σ = O ( τ 2 ) . By applying the stability result of Theorem 1, we obtain the error bound in the H 1 norm as follows:
e h n 1 2 c ( e ¯ h 0 2 + h e h 0 2 + τ n = 1 m 1 h r ˜ h n + σ 2 + τ h r ˜ h 1 2 2 + τ n = 0 m 1 R ˜ h n + σ 2 ) c ( τ 2 + h 4 + ϵ ) 2 .
The above estimate leads to the following error estimate of the FL2-1 σ ( 1 ) scheme (10). The corresponding proof is thus omitted:
Theorem 2.
Let u ( x h , t n ) and u h n be the solutions to Equation (1) and numerical scheme (10), respectively. Suppose the solution u ( x , t ) belongs to the function space C x , t 6 , 3 ( Ω × [ 0 , T ] ) . Then, for γ = 1 , we have
( u ( x h , t n ) u h n ) 1 2 c ( τ 2 + h 4 + ϵ ) .
Remark 1.
In practice, the absolute tolerance error ϵ in the scheme (10) is always set to a very small number so as not to affect the temporal or spatial convergence accuracy. In this sense, the convergence result in Theorem 2 can be understood as O ( τ 2 + h 4 ) .
Remark 2.
Although the nonuniform mesh-based FL2-1 σ ( γ ) scheme (i.e., γ > 1 ) can effectively handle non-smooth solution problems, its stability and error estimates are more difficult. In [10], Li et al. successfully gave a rigorous proof of the stability and convergence of the FL2-1 σ scheme based on the nonuniform meshes, but their derived scheme mainly deals with the subdiffusion problem. It seems that the application of their ideas here is not obvious, and further investigation is needed.

5. Numerical Examples

Two numerical examples will be presented to verify the convergence theoretical results and efficiency of the FL2-1 σ scheme (10). In the numerical examples, the computational domain is Ω = ( 0 , 1 ) 2 , and the parameters κ and μ are both one. Using the equivalence of the H 1 norm and H 1 semi-norm, we calculated the H 1 norm error at t = t n by e ( n , h ) = | u ( x h , t n ) u h n | 1 . The corresponding convergence orders in time and in space were obtained by log ( e ( n , h ) / e ( 2 n , h ) ) and log ( e ( n , h ) / e ( n , h / 2 ) ) , respectively. We always set the absolute tolerance error ϵ = 10 10 in (10) so that it did not contaminate the temporal or spatial convergence orders.
Example 1 (Accuracy).
Consider the following problem with zero Dirichlet boundary conditions:
t u ( x , y , t ) + C D 0 , t α u ( x , y , t ) = Δ u ( x , y , t ) + f ( x , y , t ) , ( x , y , t ) Ω × ( 0 , 1 ] , u ( x , y , 0 ) = sin ( π x ) sin ( π y ) , t u ( x , y , 0 ) = 0 ,
where
f ( x , y , t ) = sin ( π x ) sin ( π y ) η t η 1 + 2 π 2 ( 1 + t η ) + Γ ( η + 1 ) Γ ( η + 1 α ) t η α .
The exact solution is u = ( 1 + t η ) sin ( π x ) sin ( π y ) with the fixed parameter η > 1 .
We considered two cases, η = 3.5 and η = 1.5 , to verify the accuracy of the FL2-1 σ (γ) scheme (10). For these two cases, we applied the uniform mesh-based and nonuniform mesh-based FL2-1 σ (γ) scheme (10), respectively, to obtain the numerical results. See Table 1 and Table 2 for the first case and Table 3 and Table 4 for the second case. We remark that the H 1 norm errors here and below were obtained at the final time T.
Table 1. The H 1 -norm errors in time for the smooth case in Example 1 with h = 1 / 64 .
Table 1. The H 1 -norm errors in time for the smooth case in Example 1 with h = 1 / 64 .
N t α = 1.1 α = 1.5 α = 1.9
H 1 ErrorRate H 1 ErrorRate H 1 ErrorRate
40 2.66 × 10 4 - 6.73 × 10 4 - 1.29 × 10 3 -
80 6.69 × 10 5 1.99 1.69 × 10 4 2.00 3.23 × 10 4 2.00
160 1.68 × 10 5 1.99 4.22 × 10 5 2.00 8.07 × 10 5 2.00
320 4.27 × 10 6 1.98 1.05 × 10 5 2.01 2.01 × 10 5 2.00
Table 2. The H 1 -norm errors in space for the smooth case in Example 1 with τ = 1 / 4000 .
Table 2. The H 1 -norm errors in space for the smooth case in Example 1 with τ = 1 / 4000 .
M α = 1.1 α = 1.5 α = 1.9
H 1 ErrorRate H 1 ErrorRate H 1 ErrorRate
4 2.31 × 10 3 - 2.22 × 10 3 - 2.15 × 10 3 -
8 2.72 × 10 4 3.09 2.62 × 10 4 3.09 2.53 × 10 4 3.09
16 2.05 × 10 5 3.73 1.96 × 10 5 3.74 1.89 × 10 5 3.74
32 1.38 × 10 6 3.89 1.26 × 10 6 3.96 1.16 × 10 6 4.03
Table 3. The H 1 -norm errors in time for the non-smooth case in Example 1 with h = 1 / 64 and α = 1.5 .
Table 3. The H 1 -norm errors in time for the non-smooth case in Example 1 with h = 1 / 64 and α = 1.5 .
N t γ = 1 γ = 1.5 γ = 2
H 1 ErrorRate H 1 ErrorRate H 1 ErrorRate
40 5.86 × 10 5 - 1.85 × 10 4 - 3.70 × 10 4 -
80 3.98 × 10 5 0.56 4.42 × 10 5 2.06 9.28 × 10 5 2.00
160 2.22 × 10 5 0.84 1.03 × 10 5 2.10 2.31 × 10 5 2.00
320 1.15 × 10 5 0.94 2.27 × 10 6 2.18 5.52 × 10 6 2.07
Table 4. The H 1 -norm errors in space for the non-smooth case in Example 1 with τ = 1 / 1000 and α = 1.5 .
Table 4. The H 1 -norm errors in space for the non-smooth case in Example 1 with τ = 1 / 1000 and α = 1.5 .
M γ = 1 γ = 1.5 γ = 2
H 1 ErrorRate H 1 ErrorRate H 1 ErrorRate
4 2.52 × 10 3 - 2.52 × 10 3 - 2.52 × 10 3 -
8 3.00 × 10 4 3.07 2.97 × 10 4 3.09 2.96 × 10 4 3.09
16 2.56 × 10 5 3.55 2.22 × 10 5 3.74 2.06 × 10 5 3.84
32 5.05 × 10 6 2.34 1.33 × 10 6 4.07 3.41 × 10 7 5.92
When η = 3.5 , the solution of the equation satisfied the smoothness requirement of Theorem 2. From the numerical results in Table 1 and Table 2, it can be observed that the accuracy of the scheme (10) was O ( τ 2 + h 4 ) for different fractional orders α, which is consistent with the theoretical result of the convergence accuracy.
When η = 1.5 , the second-order partial derivative of u with respect to t was t 2 u ( x , y , t ) = η ( η 1 ) t η 2 sin ( π x ) sin ( π y ) , which was unbounded at t = 0 when η < 2 . Therefore, the solution was not sufficiently smooth in this case. In order to capture the weak singularity solution at the initial point, we let the temporal mesh grading parameter be γ > 1 . From the numerical results in Table 3 and Table 4, one can see that by selecting a different γ value (i.e., γ = 1.5 and γ = 2 ), the theoretical convergence of the scheme (10) was restored, which shows that the scheme (10) based on graded meshes indeed captured the weak singularity solution and could significantly improve the convergence accuracy.
Example 2 (Computational efficiency).
In this example, we studied the computational performance for the FL2-1 σ scheme (10). For the sake of simplicity, the case η = 3.5 in the above example was used. To better reveal that the scheme (10) had improved computational efficiency, we compared it with the two schemes (6) and (7). The numerical results are demonstrated in Table 5 and Figure 1.
Figure 1. Comparison of computational costs with fixed α = 1.5 . (The dashed lines below the numerical curves are the corresponding fitted curves, which are shifted down here for the convenience of observation.)
Figure 1. Comparison of computational costs with fixed α = 1.5 . (The dashed lines below the numerical curves are the corresponding fitted curves, which are shifted down here for the convenience of observation.)
Fractalfract 06 00438 g001
From Table 5, we can observe that when the spatial step size h was fixed, the accuracy of all three schemes increased as the number of temporal nodes increased. However, their computational cost also increased. With the fixed same number of spatial and temporal nodes, the accuracies of these threes scheme were almost the same, but the scheme (6) based on the direct solver took the most time, followed by the scheme (7), while the FL2-1 σ ( γ ) scheme (10) took the least time, and this time consumption gap increased with the increase in the number of temporal nodes.
For the case where α = 1.5 , when the time step size τ was fixed to 1 / 4 , and the spatial grid became denser, similar phenomena to those in Table 5 could be observed in the left subplot of Figure 1. In particular, we compared the computational efficiency of the two numerical schemes (7) and (10) by plotting the CPU time versus N t in a log-log style. In order to numerically quantify the complexity of these two schemes, we expressed the CPU time in the form of N t r using least squares fitting and reported the corresponding order r in the legend of the right subplot of Figure 1. One can see that the order r of the FL2-1 σ ( γ ) scheme (10) was about one, while that of the DST-based scheme (7) was near two. This indicates that the proposed scheme (10) based on the DST and SOE techniques is more competitive in solving high-dimensional time-fractional problems. For other cases, such as α = 1.1 and 1.9 , similar numerical results were also obtained, and they are not presented here for the sake of brevity.

6. Conclusions

In this paper, we proposed a fast compact difference scheme with O ( τ 2 + h 4 ) convergence accuracy to solve the two-dimensional time-fractional Cattaneo Equation (1). This scheme is based on the DST and SOE techniques and has the ability of fast calculation in both time and space. The stability and error estimate of the numerical scheme (10) based on uniform meshes are presented rigorously. Numerical examples show that this scheme is efficient and suitable for numerically solving high-dimensional time-fractional problems.
We remark that the results obtained in this paper can easily be generalized to other high-dimensional problems. Aside from that, one may notice that although numerical experiments have shown that the scheme (10) based on graded meshes can improve the accuracy of solving high-dimensional non-smooth solution problems, the analysis of its corresponding numerical theory is not an easy job. Such an issue deserves further study and will be one of our upcoming research directions.
We note that the scheme presented in this paper is only for linear problems with constant coefficients, while in practical problems, the equations are often accompanied by variable coefficients or nonlinear terms. In the following work, we will try to extend the methods of this paper to solve variable coefficients and nonlinear problems.

Author Contributions

Methodology, L.N., Q.Y., J.C. and A.C.; validation, L.N., Q.Y., J.C. and A.C.; formal analysis, L.N., Q.Y., J.C. and A.C.; investigation, L.N., Q.Y., J.C. and A.C.; writing—original draft preparation, L.N., Q.Y., J.C. and A.C.; writing—review and editing, L.N., Q.Y., J.C. and A.C.; funding acquisition, L.N., Q.Y., J.C. and A.C. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by Guangxi Natural Science Foundation grant numbers 2021GXNSFBA196027, 2020GXNSFBA297121, 2018GXNSFBA281020, and 2018GXNSFAA138121, National Natural Science Foundation of China grant number 11901266, and Gansu Natural Science Foundation grant number 21JR7RA253.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The data of the numerical simulation used to support the findings of this study are included within the paper.

Acknowledgments

The authors wish to express their appreciation to the reviewers for their valuable suggestions, which greatly improved the presentation of this paper.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Chen, A.; Nong, L. Efficient Galerkin finite element methods for a time-fractional Cattaneo equation. Adv. Differ. Equ. 2020, 2020, 545. [Google Scholar] [CrossRef]
  2. Awad, E.; Metzler, R. Crossover dynamics from superdiffusion to subdiffusion: Models and solutions. Fract. Calc. Appl. Anal. 2020, 23, 55–102. [Google Scholar] [CrossRef]
  3. Compte, A.; Metzler, R. The generalized Cattaneo equation for the description of anomalous transport processes. J. Phys. Math. Gen. 1997, 30, 7277–7289. [Google Scholar] [CrossRef]
  4. Zhao, X.; Sun, Z. Compact Crank–Nicolson Schemes for a Class of Fractional Cattaneo Equation in Inhomogeneous Medium. J. Sci. Comput. 2015, 62, 747–771. [Google Scholar] [CrossRef]
  5. Ren, J.; Gao, G. Efficient and stable numerical methods for the two-dimensional fractional Cattaneo equation. Numer. Algorithms 2015, 69, 795–818. [Google Scholar] [CrossRef]
  6. Chen, A.; Li, C. An alternating direction Galerkin method for a time-fractional partial differential equation with damping in two space dimensions. Adv. Differ. Equ. 2017, 2017, 356. [Google Scholar] [CrossRef]
  7. Li, C.; Chen, A. Numerical methods for fractional partial differential equations. Int. J. Comput. Methods 2018, 95, 1048–1099. [Google Scholar] [CrossRef]
  8. Jin, B.; Lazarov, R.; Zhou, Z. Numerical methods for time-fractional evolution equations with nonsmooth data: A concise overview. Comput. Methods Appl. Mech. Eng. 2019, 346, 332–358. [Google Scholar] [CrossRef]
  9. Li, C.; Zeng, F. Numerical Methods for Fractional Calculus; Chapman and Hall/CRC: Boca Raton, FL, USA, 2015. [Google Scholar]
  10. Li, X.; Liao, H.; Zhang, L. A second-order fast compact scheme with unequal time-steps for subdiffusion problems. Numer. Algorithms 2021, 86, 1011–1039. [Google Scholar] [CrossRef]
  11. Zeng, F.; Turner, I.; Burrage, K. A stable fast time-stepping method for fractional integral and derivative operators. J. Sci. Comput. 2018, 77, 283–307. [Google Scholar] [CrossRef]
  12. Yin, B.; Liu, Y.; Li, H.; Zeng, F. A class of efficient time-stepping methods for multi-term time-fractional reaction-diffusion-wave equations. Appl. Numer. Math. 2021, 165, 56–82. [Google Scholar] [CrossRef]
  13. Yan, Y.; Sun, Z.; Zhang, J. Fast evaluation of the Caputo fractional derivative and its applications to fractional diffusion equations: A second-order scheme. Commun. Comput. Phys. 2017, 22, 1028–1048. [Google Scholar] [CrossRef]
  14. Liao, H.; Tang, T.; Zhou, T. A second-order and nonuniform time-stepping maximum-principle preserving scheme for time-fractional Allen-Cahn equations. J. Comput. Phys. 2020, 414, 109473. [Google Scholar] [CrossRef]
  15. Wang, H.; Zhang, Y.; Ma, X.; Qiu, J.; Liang, Y. An efficient implementation of fourth-order compact finite difference scheme for Poisson equation with Dirichlet boundary conditions. Comput. Math. Appl. 2016, 71, 1843–1860. [Google Scholar] [CrossRef]
  16. Ren, Y.; Feng, H.; Zhao, S. A FFT accelerated high order finite difference method for elliptic boundary value problems over irregular domains. J. Comput. Phys. 2022, 448, 110762. [Google Scholar] [CrossRef]
  17. Nong, L.; Chen, A.; Yi, Q.; Li, C. Fast Crank-Nicolson compact difference scheme for the two-dimensional time-fractional mobile/immobile transport equation. AIMS Math. 2021, 6, 6242–6254. [Google Scholar] [CrossRef]
  18. Nong, L.; Chen, A. Fast high-order difference scheme for the modified anomalous subdiffusion equation based on fast discrete Sine transform. J. Funct. Spaces 2021, 2021, 9918955. [Google Scholar] [CrossRef]
  19. Liang, Y.; Yao, Z.; Wang, Z. Fast high order difference schemes for the time fractional telegraph equation. Numer. Methods Partial Differ. Equ. 2020, 36, 154–172. [Google Scholar] [CrossRef]
  20. Sun, H.; Sun, Z.; Gao, G. Some temporal second order difference schemes for fractional wave equations. Numer. Methods Partial Differ. Equ. 2016, 32, 970–1001. [Google Scholar] [CrossRef]
  21. Alikhanov, A.A. A new difference scheme for the time fractional diffusion equation. J. Comput. Phys. 2015, 280, 424–438. [Google Scholar] [CrossRef]
Table 5. Comparison of CPU time of three schemes (6), (7), and (10) for Example 1 with h = 1 / 64 and τ = 0.004 / 2 k .
Table 5. Comparison of CPU time of three schemes (6), (7), and (10) for Example 1 with h = 1 / 64 and τ = 0.004 / 2 k .
α Scheme k = 1 k = 2 k = 3
H 1 ErrorCPU (s) H 1 ErrorCPU (s) H 1 ErrorCPU (s)
(6) 1.66 × 10 6 3.86 3.69 × 10 7 9.85 4.78 × 10 8 28.03
1.1(7) 1.80 × 10 6 1.59 5.18 × 10 7 5.50 1.96 × 10 7 19.56
(10) 1.81 × 10 6 1.03 5.04 × 10 7 1.99 2.08 × 10 7 3.98
(6) 4.39 × 10 6 3.67 1.14 × 10 6 9.59 3.28 × 10 7 27.68
1.5(7) 4.25 × 10 6 1.43 9.97 × 10 7 5.23 1.85 × 10 7 19.22
(10) 4.25 × 10 6 0.94 9.97 × 10 7 1.96 1.85 × 10 7 3.98
(6) 8.33 × 10 6 3.77 2.12 × 10 6 9.77 5.72 × 10 7 28.13
1.9(7) 8.19 × 10 6 1.56 1.98 × 10 6 5.38 4.34 × 10 7 19.66
(10) 8.19 × 10 6 0.97 1.98 × 10 6 1.98 4.34 × 10 7 4.07
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Nong, L.; Yi, Q.; Cao, J.; Chen, A. Fast Compact Difference Scheme for Solving the Two-Dimensional Time-Fractional Cattaneo Equation. Fractal Fract. 2022, 6, 438. https://doi.org/10.3390/fractalfract6080438

AMA Style

Nong L, Yi Q, Cao J, Chen A. Fast Compact Difference Scheme for Solving the Two-Dimensional Time-Fractional Cattaneo Equation. Fractal and Fractional. 2022; 6(8):438. https://doi.org/10.3390/fractalfract6080438

Chicago/Turabian Style

Nong, Lijuan, Qian Yi, Jianxiong Cao, and An Chen. 2022. "Fast Compact Difference Scheme for Solving the Two-Dimensional Time-Fractional Cattaneo Equation" Fractal and Fractional 6, no. 8: 438. https://doi.org/10.3390/fractalfract6080438

APA Style

Nong, L., Yi, Q., Cao, J., & Chen, A. (2022). Fast Compact Difference Scheme for Solving the Two-Dimensional Time-Fractional Cattaneo Equation. Fractal and Fractional, 6(8), 438. https://doi.org/10.3390/fractalfract6080438

Article Metrics

Back to TopTop