On Caputo-Type Cable Equation: Analysis and Computation

: In this paper, a special case of nonlinear time fractional cable equation is studied. For the equation deﬁned on a bounded domain, the existence, uniqueness, and regularity of the solution are ﬁrstly studied. Furthermore, it is numerically studied via the weighted and shifted Grünwald difference (WSGD) methods/the local discontinuous Galerkin (LDG) ﬁnite element methods. The derived numerical scheme has been proved to be stable and convergent with order O (∆ t 2 + h k +1 ) , where ∆ t , h , k are the time stepsize, the spatial stepsize, and the degree of piecewise polynomials, respectively. Finally, a numerical experiment is presented to verify the theoretical analysis.


Introduction
In this paper, we consider a special case of nonlinear time fractional cable equation in the following form, ∂u(x, t) ∂t + C D α 0,t u(x, t) − u xx (x, t) + f (u) = g(x, t), x ∈ Ω, t > 0, with the initial value condition, and the boundary value condition, where 0 < α < 1, Ω = (a, b) is a bounded domain, g, u 0 are given smooth functions, C D α 0,t is the α-th order Caputo derivative operator defined by Podlubny [Podlubny (1999)] in which Γ(·) denotes the Gamma function. We always suppose that the nonlinear source term f (u) satisfies Lipschitz continuity condition with respect to u, that is, there exists a positive constant L such that for all u 1 , u 2 , Cable equations with fractional order temporal operators were introduced to model electrotonic properties of spiny neuronal dendrites by Henry et al. [Henry, Langlands and Wearne (2008)]. The time fractional cable equation (TFCE) is similar to the traditional cable equation except that the order of derivative with respect to the time is fractional. If there is a nonlinear source term, the equation reads as (with 0 < α, β < 1) [Henry, Langlands and Wearne (2008)] ∂u(x, t) ∂t which has been numerically treated by a number of authors. For example, Lin et al. [Lin, Li and Xu (2011)] constructed a finite difference/Legendre spectral scheme for discretization of TFCE. Hu et al. [Hu and Zhang (2012)] proposed two implicit compact difference schemes for TFCE. A fourth-order compact finite difference scheme for 2D TFCE was studied by Yu et al. [Yu and Jiang (2016)]. Zheng et al. [Zheng and Zhao (2017)] developed and analyzed a time LDG method (LDG method is applied in time direction) for solving TFCE. Al-Maskari et al. [Al-Maskari and Karaa (2018)] discussed the lumped mass Galerkin finite element method for TFCE. Recently, a scheme combining a finite difference approach in time direction and LDG finite element method in space direction for TFCE was proposed by Li et al. ]. They proved that the derived scheme could reach 2-nd order in time direction, which was higher than the classical L1 method. Liu et al. [Liu, Du, Li et al. (2019)] considered some second-order θ schemes combined with Galerkin finite element method for TFCE.
There seems no work on the mathematical analysis and LDG method for it. This motivates our interest in studying Eq. (1).
LDG method is a special class of discontinuous Galerkin method, proposed by Cockburn et al. [Cockburn and Shu (1998)]. The main technique of LDG method is to rewrite higher-order derivative equation into an equivalent system containing only the first derivative, and then discretize it by the standard discontinuous Galerkin method. For more information about this method, see the review paper by Xu et al. [Xu and Shu (2010)].
Here we propose the LDG finite element methods to numerically study Eq. (1). The main contributions of this paper are twofold: One is to provide a complete mathematical analysis for Eq. (1), including existence, uniqueness, and regularity of the solution; The other is to numerically studied Eq. (1) using WSGD method in time domain and using the LDG finite element method in space domain. The derived numerical scheme is stable and convergent with order O(∆t 2 + h k+1 ).
The rest of this paper is organized as follows. In Section 2, we introduce some preliminaries, which will be used in the following section. In Section 3, we discuss the existence, uniqueness, and regularity for the solution to Eq. (1). In Section 4, a fully discrete LDG scheme is proposed and the stability and convergence of the presented scheme is analyzed too. A numerical experiment is given in Section 5 to illustrate the effectiveness of the proposed numerical method. Finally, the last section concludes this paper.

Notations
We first recall some notations and preliminary facts, which are used throughout this paper. The L 2 -norm and inner product on Ω are given by Likewise, we define the L ∞ -norm on Ω by u ∞ = sup x∈Ω |u|.
The Laplace transform of a given function v(t) is defined as Li et al. [Li and Zeng (2015)] and the inverse Laplace transform is given by where c 0 lies in the right half plane of the absolute convergence of the Laplace transform (6).

Solution representation
Let w = u − u 0 , then (1) can be rewritten as the following equivalent system ∂w(x,t) ∂t Using Laplace transform, we obtain which further implieŝ Then by inverse Laplace transform and convolution rule, the solution of Eq. (9) can be represented by where the operators E (t), F (t) : X → X are defined by For fixed δ > 0 and θ ∈ ( π 2 , π), the contour of integration Γ θ,δ is defined by and with Ims increasing.
Now we obtain a representation of the solution of Eq. (1), 3 Regularity of the solution Before we present the main theorem of this section, we would like to give two useful lemmas here. Lemma 3.1. For the operators E (t) and F (t), the following estimates hold for any t ∈ (0, T ], m ∈ N 0 , and ν = 0, 1: The proof line of (i) is similar to that of Theorem 2.1 in Al-Maskari et al. [Al-Maskari and Karaa (2018)] (we refer to Theorem 2.1 for the case α 2 = 1 ). Since F (t) = E (t), (ii) immediately follows from (i).
Note that AE = E A : X → X is continuous with respect to t ∈ [0, T ]. Then taking t → 0 in (10) implies (iii). This ends the proof. Lemma 3.2 (Zeng, Cao and Li (2013), Gronwall's inequality). Let q(t) be continuous and nonnegative on [0, T ]. If where 0 ≤ µ < 1, c(t) is nonnegative monotonic increasing continuous function on [0, T ], and h is a positive constant, then Now we consider the existence, uniqueness, and regularity of the solution to Eq. (1). Theorem 3.1. For a given T > 0, suppose that u 0 ∈ D and g ∈ C [0, T ]; H 2 (Ω) . f : R → R is Lipschitz continuous. Then Eq. (1) has a unique solution u such that Proof.
Step 1: Existence and uniqueness. Following the idea in Li et al. ], we first define a map M : where the space Then we only have to prove that for some λ > 0, the map M has a unique fixed point. For where we have used Lemma 3.1 with ν = m = 0 in the first inequality and f is Lipschitz continuous in the second inequality.
By choosing a sufficiently large λ such thatC = C/λ < 1 and taking maximum of the left hand side of (15) with respect to t ∈ [0, T ], there holds Finally, applying the Banach fixed point theorem, we can obtain that Eq.
Step 2: C α [0, T ]; X regularity. Consider the following difference quotient for ∆t > 0 which will be estimated respectively as follows.
Applying Lemma 3.1, we arrive at where C α,T is a positive constant depending on α and T . By using Lemma 3.1 and the Lipschitz continuity of f again, we have and Similarly, the two terms I 4 and I 5 are shown to be bounded, respectively, by and Denoting W (t) = ∆t −α u(t + ∆t) − u(t) Ω and substituting the estimates of I i (i = 1, 2, . . . , 5) into (16), we obtain which together with Lemma 3.2 yields u ∈ C α [0, T ]; X . The assertion C D α 0,t u ∈ C [0, T ]; X is a direct result of the C α [0, T ]; X regularity and the mapping property of Caputo derivative.
Step 3: C [0, T ]; D regularity. By applying the operator A to both sides of (10), we arrive at Then by Lemma 3.1 and the regularity assumption of g, we have Step 4: Estimate of u (t). By differentiating (10) with respect to t, we have It follows from Lemma 3.1 that Using Lemma 3.2 again yields the assertion of (13). The proof is thus complete.

The LDG method and its convergence
In this section, we first present the semidiscrete scheme and fully discrete scheme for problem (1), where the time fractional derivative is discretized by WSGD method and the spatial derivative by the LDG method. Then we prove that the fully discrete LDG scheme is stable and convergent.
The usual notations of LDG method are introduced here. Assume that the mesh consisting of cells where P k (I j ) denotes the space of polynomials in I j of degree at most k ≥ 0.
As the usual treatment, we would like to introduce the Gauss-Radau projections P ± h [Castillo, Kanschat, Schotzau et al. (2002)]: for any scalar function q ∈ H 1 (Ω), the projection is the unique element in V h , satisfying for any j = 1, 2, . . . , N .
Suppose q ∈ H k+1 (Ω), then by a standard scaling argument [Ciarlet (1978)], there holds where C is a positive constant independent of h.

Semidiscrete scheme
On the space V h , the L 2 (Ω)-orthogonal projection P h : L 2 (Ω) → V h and the discrete Laplacian and respectively.
Replacing the exact solutions by the numerical solutions, then we can define the semidiscrete LDG scheme as follows: find u h (·, t) ∈ V h such that where g h = P h g. By a similar argument as (10), the solution to (24) can be represented by As proved in Theorem 3.1, we have the following similar results.
Theorem 4.1. Suppose that f , g and u 0 satisfy the conditions in Theorem 3.1. Then Eq.
(24) has a unique solution u h such that

Fully discrete LDG scheme
Let ∆t = T /M be the time mesh size, t n = n∆t, n = 0, 1, . . . , M be the mesh point, M ∈ Z + . For simplicity of notations, we denote u n+1 = u(x, t n+1 ) and δ n+1 , then the time fractional derivative (4) at time t n+1 can be approximated as Wang et al. [Wang and Vong (2014) where In what follows, we would like to introduce several lemmas which are very important in obtaining the error estimate.
given as above, the following be defined as above. Then for any positive integer k and real vector (v 1 , v 2 , . . . , v k ) ∈ R k , it holds that In order to get the LDG formulation, we firstly rewrite Eq. (1) into the following lower-order system of two equations by introducing an auxiliary variable p = ∂u/∂x Then we can get the weak form of Eq. (26) at t n+1 as where When n = 0, we take u −1 = 2u 0 − u 1 + O(∆t 2 ) by Taylor expansion.
Let u n+1 h , p n+1 h ∈ V h be the approximation of u n+1 and p n+1 , respectively. We get the fully discrete LDG scheme as follows: find u n+1 h , p n+1 where u −1 h = 2u 0 h − u 1 h . The "tilde" terms are the so-called "numerical fluxes", which are taken as the "alternating" numerical flux u n h = (u n h ) − , p n h = (p n h ) + .
Remark 4.1. The choice for the fluxes (29) is not unique. We can also take the numerical fluxes as u n h = (u n h ) + , p n h = (p n h ) − on each cell interface.

Stability analysis
In this subsection, we consider the stability of the LDG scheme (28). Let (U n h , P n h ) be the perturbed solution of (u n h , p n h ), i.e., (U n h , P n h ) and (u n h , p n h ) satisfy (27) with different initial conditions. Theorem 4.2. Suppose that f and u 0 satisfy the conditions in Theorem 3.1, then the fully discrete LDG scheme (28) with flux (29) is stable.
, we obtain the following perturbation equation Let v h = e n+1 u h and w h = e n+1 p h . Then (30) can be written as Firstly, we prove that the theorem holds true for n = 0. Taking n = 0 in (31) and adding the two equations together lead to Multiplying (32) by 2∆t and using Cauchy-Schwarz inequality, we have Noticing that b α (0) > 0 and b α (1) < 0, we get Now we are going to prove the case of n ≥ 1. Adding the two equations of (31) together results in Multiplying (34) by 4∆t and using Cauchy-Schwarz inequality and Young's inequality, we arrive at Summing n from 1 to k and using Lemma 4.2, we can get Then from discrete Gronwall's lemma (i.e., Lemma 4.3) and (33), it yields that Combining (33) with (35), we complete the proof of this theorem.

Error estimate
In this subsection, we will give the error estimate for the fully discrete LDG scheme (28). Theorem 4.3. Assume that f , g and u 0 satisfy the conditions in Theorem 3.1. Let u n+1 be the exact solution of (1) and u n+1 h be the numerical solution of the fully discrete LDG scheme (28) with fluxes (29). If we assume that u(x, t) ∈ C 2 [0, T ]; H k+1 (Ω) , then there holds where C is a positive constant independent of ∆t and h.
Proof. Denote According to (21), we have Thus in what follows, we will focus on the estimate for ξ n u .
In order to estimate ξ n u , we would like to set up the corresponding error equation first. Subtracting (28) with (27) and using the flux (29), we get Substituting (37) into (39), we can get the following equations.
Owing to the property (21), we get Multiplying (42) by 2∆t, it is easy to see that Then applying the Cauchy-Schwarz inequality, Young's inequality, (21), (43), as well as (44), we have As a consequence, if we let ∆t < 1 2L , we can get Next we are going to prove Case II. By taking (v h , w h ) = (ξ n+1 u , ξ n+1 p ) in (41), we can derive Multiplying (46)   where Cauchy-Schwarz inequality, Young's inequality, and (43) are used in the last step.
With the help of (21), we can further obtain By virtue of Lemma 4.2, we have Then summing (47) for n from 1 to K − 1 leads to where the property (21) and the result of (45) are used for the second inequality. Finally, it follows straightforwardly by using Lemma 4.3 that which combine the triangle inequality to complete the proof of this theorem. Remark 4.2.
(1) We must remark here that the above error estimate is optimal both in time and space.
(2) Compared with the classical L1 method with time convergence rate of (2 − α), our scheme can arrive at second order in time.
(3) Our discussions focus on Caputo-type partial differential equation, it may be interesting to extend the analysis to Riesz-type fractional differential equation ].

Numerical example
In this section, we present a numerical example to verify the theoretical results. Example 5.1. Consider the following Caputo-type cable equation with compactly supported boundary condition, on Ω = (0, 2π), t ∈ (0, 1]. The initial value condition is u(x, 0) = 0, x ∈ (0, 1), and the source term is The exact solution of (49) is given by u(x, t) = t 2 sin(x).
The results presented in this paper indicate that the proposed LDG scheme enjoys the same accuracy as the spectral schemes in Liu et al. [Liu and Lü (2019) ;Yang, Jiang and Zhang (2018)]. However, if the geometry and boundary conditions are complicated, LDG method may be more suitable and can achieve the uniformly high-order accuracy, which is what I will do next [Xu and Shu (2010)]. Besides, it is of much interest to investigate the blow-up phenomenon of the solution, see for example [Cao, Song, Wang et al. (2019)]. In the future work, I will consider using LDG method to deal with variable order fractional differential equations.