A linearized second-order scheme for nonlinear time fractional Klein-Gordon type equations

TL;DR

This paper introduces a linearized second-order numerical scheme for nonlinear time fractional Klein-Gordon equations, avoiding iterative methods and ensuring convergence with second-order temporal accuracy.

Contribution

A novel linearized difference scheme for nonlinear time fractional Klein-Gordon equations that achieves second-order accuracy without iterative solvers.

Findings

The scheme converges with second-order accuracy in time.

Standard energy estimates are adapted using a new grid function.

The method simplifies computation by eliminating the need for iterative methods.

Abstract

We consider difference schemes for nonlinear time fractional Klein-Gordon type equations in this paper. A linearized scheme is proposed to solve the problem. As a result, iterative method need not be employed. One of the main difficulties for the analysis is that certain weight averages of the approximated solutions are considered in the discretization and standard energy estimates cannot be applied directly. By introducing a new grid function, which further approximates the solution, and using ideas in some recent studies, we show that the method converges with second-order accuracy in time.

Equations

$${}_0^C{D}_t^{\alpha }u(x,t)=\frac{{\partial }^{2}u(x,t)}{\partial x^2}-f(u)+p(x,t),x\in (a,b),t\in (0,T],$$

$$u(a,t)=u(b,t)=0,t\in (0,T],$$

$$u(x,0)=\phi (x),u_t(x,0)=\psi (x),x\in [a,b],$$

$${}_0^C{D}_t^{\alpha }u(x,t)=\frac{1}{\Gamma (2-\alpha )}{\int }_0^t\frac{{\partial }^{2}u(x,s)}{\partial s^2}\frac{ds}{{(t-s)}^{\alpha -1}},$$

$$|f({\varphi }_1)-f({\varphi }_2)|\le L|{\varphi }_1-{\varphi }_2|,\forall {\varphi }_1,{\varphi }_2\in \Omega .$$

$${\delta }_x{u}_{i-\frac{1}{2}}=\frac{u_i-u_{i-1}}{h},{\delta }_x^2{u}_i=\frac{{\delta }_x{u}_{i+\frac{1}{2}}-{\delta }_x{u}_{i-\frac{1}{2}}}{h}=\frac{u_{i+1}-2u_i+u_{i-1}}{h^2},$$

$$⟨u,v⟩=h\sum _{i=1}^{M-1}{u}_i{v}_i,\Vert u\Vert =\sqrt{⟨u,u⟩},|u{|}_1=\sqrt{h\sum _{i=1}^M\mid {\delta }_x{u}_{i-\frac{1}{2}}{\mid }^2},\Vert u{\Vert }_{\infty }=\underset{1\le i\le M-1}{\max }|u_i|.$$

$${\delta }_t{u}^{n+\frac{1}{2}}=\frac{u^{n+1}-u^n}{\tau },{\delta }_{\hat{t}}{u}^n=\frac{u^{n+1}-u^{n-1}}{2\tau }.$$

$$\frac{1}{\Gamma (2-\alpha )}{\int }_0^{t_{n+\theta }}\frac{v^{\prime }(s)}{(t_{n+\theta }-s{)}^{\alpha -1}}ds-{\Delta }_t^{\alpha }v(t_{n+\theta })=\mathcal{O}({\tau }^2),$$

$$c_{m}^{(n+1)}=\left\{\begin{array}[]{ll}a_{0}+b_{1},&m=0,\\ a_{m}+b_{m+1}-b_{m},&1\leq m\leq n-1,\\ a_{n}-b_{n},&m=n;\end{array}\right.$$

$$(a)0<b_k<\frac{\alpha -1}{2(3-\alpha )}{a}_k<\frac{1}{2}{a}_k,k\ge 1,(b)\sum _{k=1}^n{b}_k<\frac{1}{2}\sum _{k=0}^n{a}_k=\frac{1}{2}(n+\theta {)}^{2-\alpha },n\ge 1,$$

$$(c)c_n^{(n+1)}>\frac{2-\alpha }{2(n+\theta {)}^{\alpha -1}},(d)c_0^{(n+1)}>c_1^{(n+1)}>\cdots >c_{n-1}^{(n+1)}>c_n^{(n+1)},$$

$$(e)(2\theta -1)c_0^{(n+1)}-\theta c_1^{(n+1)}>0,(f)\sum _{k=0}^n{c}_k^{(n+1)}=\sum _{k=0}^n{c}_k^{(k+1)}+\sum _{k=1}^n{b}_k=(n+\theta {)}^{2-\alpha },$$

$$(g)c_{n-1}^{(n+1)}=c_{n-1}^{(n)}+b_n,n\ge 1;c_k^{(n+1)}=c_k^{(n)},0\le k\le n-2,n\ge 2,$$

$$(h)\sum _{k=0}^n\frac{1}{c_k^{(n+1)}}<\frac{2{(n+1)}^{\alpha }}{2-\alpha },n\ge 1.$$

$$\sum _{k=0}^n\frac{1}{c_k^{(n+1)}}<\sum _{k=0}^n\frac{1}{c_n^{(n+1)}}<\frac{2(n+1){(n+\theta )}^{\alpha -1}}{2-\alpha }<\frac{2{(n+1)}^{\alpha }}{2-\alpha },$$

$$d_{0}^{(1)}=\frac{c_{0}^{(1)}}{\mu}=\frac{\theta^{2-\alpha}\tau^{1-\alpha}}{\Gamma(3-\alpha)};\quad\quad d_{k}^{(n+1)}=\left\{\begin{array}[]{ll}\frac{c_{k}^{(n+1)}}{\mu},&0\leq k\leq n-1,\\ \frac{c_{n}^{(n+1)}-\frac{\theta}{1-\theta}b_{n}}{\mu},&k=n,\end{array}\right.\quad n\geq 1.$$

$$(a)\frac{(2-\alpha {)}^2{t}_{n+\theta }^{1-\alpha }}{\Gamma (4-\alpha )}<d_n^{(n+1)}<\frac{t_{n-1+\theta }^{1-\alpha }}{\Gamma (2-\alpha )},(b)d_0^{(n+1)}>d_1^{(n+1)}>\cdots >d_{n-1}^{(n+1)}>d_n^{(n+1)},$$

$$(c)(2\theta -1)d_0^{(n+1)}-\theta d_1^{(n+1)}>0,(d)\tau \sum _{k=0}^n{d}_k^{(n+1)}<\frac{t_{n+\theta }^{2-\alpha }}{\Gamma (3-\alpha )},$$

$$(e)\tau \sum _{k=0}^n{d}_k^{(k+1)}<\frac{t_{n+\theta }^{2-\alpha }}{\Gamma (3-\alpha )},(f)\tau \sum _{k=0}^n\frac{1}{d_k^{(k+1)}}<\frac{\Gamma (4-\alpha )T^{\alpha }}{(2-\alpha {)}^2}.$$

$$d_n^{(n+1)}=\frac{a_n-\frac{1}{1-\theta }{b}_n}{\mu }<\frac{a_n}{\mu }=\frac{2-\alpha }{\mu }{\int }_0^1\frac{ds}{(n+\theta -s{)}^{\alpha -1}}<\frac{t_{n-1+\theta }^{1-\alpha }}{\Gamma (2-\alpha )},$$

$$d_n^{(n+1)}=\frac{a_n-\frac{1}{1-\theta }{b}_n}{\mu }>\frac{2-\alpha }{(3-\alpha )\mu }{a}_n>\frac{(2-\alpha {)}^2}{(3-\alpha )\mu }{\int }_0^1\frac{ds}{(n+\theta -s{)}^{\alpha -1}}>\frac{(2-\alpha {)}^2{t}_{n+\theta }^{1-\alpha }}{\Gamma (4-\alpha )}.$$

$$(2\theta-1)d_{0}^{(n+1)}-\theta d_{1}^{(n+1)}=\big{[}(2\theta-1)c_{0}^{(n+1)}-\theta c_{1}^{(n+1)}\big{]}/\mu>0,~{}n\geq 2;$$

$$\tau \sum _{k=0}^n\frac{1}{d_k^{(k+1)}}<\tau \sum _{k=0}^n\frac{\Gamma (4-\alpha )t_{n+\theta }^{\alpha -1}}{(2-\alpha {)}^2}<\frac{\Gamma (4-\alpha )t_{n+1}^{\alpha }}{(2-\alpha {)}^2}\le \frac{\Gamma (4-\alpha )T^{\alpha }}{(2-\alpha {)}^2},$$

$$g(t_n)=(1-\theta )g(t_{n+\theta })+\theta g(t_{n-1+\theta })+\mathcal{O}({\tau }^2).$$

$$g(t_{n+\theta })=g(t_n)+\theta \tau g^{\prime }(t_n)+{\theta }^2{\tau }^2{\int }_0^1{g}^{\prime \prime }(t_n+\rho \theta \tau )(1-\rho )d\rho ,$$

$$\displaystyle g(t_{n-1+\theta})=g(t_{n})-(1-\theta)\tau g^{\prime}(t_{n})+(1-\theta)^{2}\tau^{2}\int_{0}^{1}g^{\prime\prime}\big{(}t_{n}-\rho(1-\theta)\tau\big{)}(1-\rho)d\rho.$$

$$U_i^n=u(x_i,t_n),\text{mbox}and{V}_i^{n+\sigma }=v(x_i,t_{n+\sigma })\text{mbox}for0\le \sigma \le 1,0\le i\le M,0\le n\le N.$$

$$\displaystyle^{C}_{0}D_{t}^{\alpha}u(x_{i},t_{n})=\frac{\partial^{2}u(x_{i},t_{n})}{\partial x_{i}^{2}}-f\big{(}u(x_{i},t_{n})\big{)}+p(x_{i},t_{n}),\quad 1\leq i\leq M-1,~{}1\leq n\leq N.$$

$${}_0^C{D}_t^{\alpha }u(x_i,t_n)=$$

$$=$$

$$(1-\theta ){\Delta }_t^{\alpha }{V}_i^{n+\theta }+\theta {\Delta }_t^{\alpha }{V}_i^{n-1+\theta }=$$

$$=$$

$$-\frac{\theta }{\mu }{b}_n(V_i^1-V_i^0)$$

$$=$$

$$+\frac{(1-\theta )}{\mu }{c}_n^{(n+1)}(V_i^1-V_i^0)-\frac{\theta }{\mu }{b}_n(V_i^1-V_i^0)$$

$$=$$

$$\displaystyle+\frac{1}{\mu}\big{[}(1-\theta)c_{n}^{(n+1)}-\theta b_{n}\big{]}(V_{i}^{1}-V_{i}^{0})$$

$$=$$

Full text

Key words: Linearized scheme; Time fractional differential equations; Nonlinear Klein-Gordon equations; Convergence

1 Introduction

This paper studies finite difference schemes for nonlinear fractional order Klein-Gordon type equations. Fractional differential equations have applications in physics, biology and petroleum industry. Interested reader can refer to [3, 1, 4, 5, 6, 2] for more details. One of the key features of fractional derivatives are their nonlocal dependence which make fractional differential equations suitable to model some phenomena. However, the nonlocal dependence causes difficulty in the study of these equations. In the past decade, many works have been done on the study of effective numerical method for time fractional differential equations, the most popular approach are finite difference, spectral and finite element, see [18, 16, 10, 20, 15, 14, 11, 12, 9, 8, 7, 17, 19, 13] and the references therein. Klein-Gordon equation is a basic equation to describe many phenomena. Solving it in numerically is an interesting topic. Many efficient methods have been employed to solve the linear and nonlinear Klein-Gordon, or sin-Gordon equations successfully. Such as the Adomian’s decomposition method [22, 21], the variational iteration method (VIM) [23], the He’s variational iteration method [24] and the Homotopy analysis method (HAM) [25], and so on. When studying this kind of equations with fractional order derivative, it would be more challenging.

In this paper, we consider finite difference schemes for nonlinear time fractional Klein-Gordon type equations with the following form:

[TABLE]

where $1<\alpha<2$, ${}^{C}_{0}D_{t}^{\alpha}$ denotes the Caputo’s derivative which is defined by

[TABLE]

and $f$ is a continuous function and satisfies the Lipschitz condition:

[TABLE]

Here $\Omega$ is a suitable domain, and $L$ is a positive constant only depends on the domain $\Omega$.

Lots of literatures have devoted to the study of time fractional Klein-Gordon (or sin-Gordon) type equations, see also [26, 28, 29, 30, 31, 32, 27] and references therein. The authors in [29, 30] proposed space compact schemes to solve the one and two dimensional time fractional Klein-Gordon type equations, respectively, and stability and convergence were analyzed by energy method. In [31], a fully spectral scheme with finite difference discretization in time and Legendre spectral approximation in space was derived. Moreover, the meshless method based on radial basis function was used in [32] to obtain an unconditionally stable discrete scheme for this kind of equation. We note that the finite difference schemes (or finite difference discretization in time) proposed in [29, 30, 31, 32] are nonlinear with temporal convergence order ${\cal O}(\tau^{3-\alpha})~{}(1<\alpha<2)$, or linear with convergence order ${\cal O}(\tau)$. These motivate us to investigate linearized and higher temporal convergence order scheme to solve the nonlinear equations. We remark that linearized scheme was shown to be very efficiently for dealing with nonlinear problems [37, 34, 36, 35]. The main advantage is that the nonlinear term is evaluated at previous time level so that iterative method is not needed to solve the solution at the current time level, which would be more convenience and save much computational costs. However, to our knowledge, the idea has not been applied to construct second temporal convergence order scheme for time-fractional differential equation. Our scheme in this paper is second-order in time. The proposed method is based on the descretization given in [33] and the idea of linearized scheme. We further note that the discretization formula for fractional derivatives developed in [33] is not given at grid points. This induce some technical difficulties for shifting the evaluation of nonlinear term to previous time level. Inspired by some estimates in the recent works [38, 39], we show that our proposed scheme converges with second-order in time.

The rest of the paper is organized as follows. In section 2, we first give some estimates of the discretization coefficients on the fractional derivative, and using the weighted approach we derive a linearized implicit scheme for the problem (1.1)–(1.3). The scheme is shown to be convergent with ${\cal O}(\tau^{2}+h^{2})$ and stable by discrete energy method in section 3. Spatial fourth-order compact scheme is proposed in section 4. In section 5, we test some numerical examples to confirm the theoretical results. A brief conclusion is followed in the last section.

2 Derivation of the difference scheme

2.1 Preliminary notations and lemmas

The following notations are needed to present our scheme. Let $\tau=\frac{T}{N}$ and $h=\frac{b-a}{M}$ be the temporal and spatial step sizes respectively, where $N$ and $M$ are some given integers. For $n=0,1,...,N$, and $i=0,1,...,M$, denote $t_{n}=n\tau$, $x_{i}=ih$, $t_{n+\theta}=(n+\theta)\tau$ for a constant $\theta\in[0,1]$, $\varphi_{i}=\varphi(x_{i})$ and $\psi_{i}=\psi(x_{i})$. We next introduce the grid function spaces \mathcal{V}_{h}=\big{\{}u|u=\{u_{i}|0\leq i\leq M\}\mbox{ and }u_{0}=u_{M}=0\big{\}} and $\mathcal{W}_{\tau}=\{w^{n}|0\leq n\leq N\}$. For any $u,v\in\mathcal{V}_{h}$, we denote

[TABLE]

and the inner product and norms

[TABLE]

For any $u^{n}\in\mathcal{W}_{\tau}$, we further consider

[TABLE]

Discretization on the fractional derivative of our scheme based on the following lemma, which is obtained straightly by replacing the parameter $\sigma=\frac{2-\alpha}{2}$ $(0<\alpha<1)$ in Lemma 2.3 of [38] with $\theta=\frac{3-\alpha}{2}$ $(1<\alpha<2)$ and, in fact, it also can be found in Lemma 2.3 of [39].

Lemma 2.1**.**

Suppose $1<\alpha<2,~{}\theta=\frac{3-\alpha}{2},~{}v(t)\in{\cal C}^{2}[0,T]\cap{\cal C}^{3}(0,T]$, and there exists a positive constant $C$ such that $v^{\prime\prime\prime}(t)\leq Ct^{\alpha-2}$ in $[0,T]$. Then

[TABLE]

where $\Delta_{t}^{\alpha}v(t_{n+\theta})=\frac{\tau^{1-\alpha}}{\Gamma(3-\alpha)}\sum_{s=0}^{n}c_{n-s}^{(n+1)}[v(t_{s+1})-v(t_{s})]$, and $c_{0}^{(1)}=a_{0}$ for $n=0$,

[TABLE]

for $n\geq 1$, in which $a_{0}=\theta^{2-\alpha}$, $a_{l}=(l+\theta)^{2-\alpha}-(l-1+\theta)^{2-\alpha}$, for $l\geq 1$; and $b_{l}=\frac{1}{3-\alpha}[(l+\theta)^{3-\alpha}-(l-1+\theta)^{3-\alpha}]-\frac{1}{2}[(l+\theta)^{2-\alpha}+(l-1+\theta)^{2-\alpha}].$

Lemma 2.2**.**

The coefficients $a_{k}$, $b_{k}$, and $c_{k}^{(n+1)}$ defined in Lemma 2.1 satisfy

[TABLE]

Proof.

The inequalities $(a)$–$(g)$ are obtained directly by replacing $\sigma=\frac{2-\alpha}{2}$ $(0<\alpha<1)$ in Lemma 2.1 and Lemma 2.2 of [38] with $\theta=\frac{3-\alpha}{2}$ $(1<\alpha<2)$, which also can be found in Lemma 2.1 and Lemma 2.2 of [39]. Using $(c)$–$(d)$, we have

[TABLE]

thus $(h)$ is verified. ∎

Denote $\mu=\tau^{\alpha-1}\Gamma(3-\alpha)$ and

[TABLE]

We further have the following lemma.

Lemma 2.3**.**

The above coefficients $d_{k}^{(n+1)}$ $(0\leq k\leq n\leq N-1,~{}n\geq 1)$ satisfy

[TABLE]

Proof.

Applying Lemma 2.2(a), we get

[TABLE]

So (a) is verified. Since $b_{n}>0$, it follows from Lemma 2.2(d) that $c_{0}^{(n+1)}>c_{1}^{(n+1)}>\cdots>c_{n-1}^{(n+1)}>c_{n}^{(n+1)}-b_{n}$, then the definition (2.1) yields the inequality (b). Similarly, Lemma 2.2(e) implies that

[TABLE]

and (2\theta-1)d_{0}^{(2)}-\theta d_{1}^{(2)}=\big{[}(2\theta-1)c_{0}^{(2)}-\theta c_{1}^{(2)}\big{]}/\mu>0. Thus (c) is verified. The definition (2.1) also implies $\tau\sum_{k=0}^{n}d_{k}^{(n+1)}<\frac{\tau}{\mu}\sum_{k=0}^{n}c_{k}^{(n+1)}$ such that the inequality (d) is obtained by using Lemma 2.2(f). The proof of (a) shows that $d_{k}^{(k+1)}<\frac{a_{k}}{\mu}$ for $k\geq 1$. Then Lemma 2.2(b) yields the inequality (e). It is easy to check that $\frac{\theta}{\Gamma(3-\alpha)}>\frac{(2-\alpha)^{2}}{\Gamma(4-\alpha)}$ for $1<\alpha<2$, then $d_{0}^{(1)}=\frac{\theta t_{\theta}^{1-\alpha}}{\Gamma(3-\alpha)}>\frac{(2-\alpha)^{2}t_{\theta}^{1-\alpha}}{\Gamma(4-\alpha)}$, and it follows by combining (a) that

[TABLE]

so inequality (f) is verified. ∎

We have the following lemma relating solution values at different points.

Lemma 2.4**.**

For any $g(t)\in{\cal C}^{2}[t_{n-1+\theta},t_{n+\theta}]$, it holds that

[TABLE]

Proof.

Using Taylor formula with integral remainder, we have

[TABLE]

Then the desired result can be obtained by a direct calculation. ∎

2.2 Weighted approximation to time fractional derivative

Denote $v=u_{t}$, then ${}^{C}_{0}D_{t}^{\alpha}u=\frac{1}{\Gamma(2-\alpha)}\int_{0}^{t}\frac{v_{t}(s)}{(t-s)^{\alpha-1}}ds$. We define the grid functions

[TABLE]

Consider equation (1.1) at the point $(x_{i},t_{n})$, we have

[TABLE]

Utilizing Lemma 2.4 and Lemma 2.1, we introduce a weighted approximation for the time fractional derivative. Let ${}^{C}_{0}D_{t}^{\alpha}u(x_{i},t)\in{\cal C}^{2}[t_{n-1+\theta},t_{n+\theta}]$ and $v(x_{i},t)\in{\cal C}^{2}[0,T]\cap{\cal C}^{3}(0,T]$, $1\leq n\leq N-1$, it follows that

[TABLE]

where $({\hat{R}}_{t})_{i}^{n}={\cal O}(\tau^{2})$.

Denote $V_{i}^{k+1-\theta}=(1-\theta)V_{i}^{k+1}+\theta V_{i}^{k}$ ($k\geq 0$). Then, for $n\geq 1$, using Lemma 2.2(g) and (2.1), we have

[TABLE]

If $n=0$,

[TABLE]

For $k\geq 1$, Taylor expansion gives

[TABLE]

and

[TABLE]

where

[TABLE]

[TABLE]

and

[TABLE]

So, by inserting (2.6) and (2.7) into (2.4), we can get the following approximation for the fractional derivative on grid function $U_{i}^{n}$

[TABLE]

where

[TABLE]

Utilizing the Taylor expansion, it is observed that

[TABLE]

These three equations above yield

[TABLE]

where $C_{1}$ is a positive constant.

Then, suppose $\frac{\partial^{4}u}{\partial t^{4}}$ is continuous on $[a,b]\times[0,T]$, we can obtain

[TABLE]

where Lemma 2.3(a),(d) have been used, and $C_{2},C_{3}$ are positive constants.

2.3 Approximation on first time level and space discretization

Note that the discretization (2.8) is devoted to solving the numerical solutions of $U^{n+1}$ ($n\geq 1$). For the approximation on the first time level, further construction is required. Inserting (2.7) into (2.5), we have

[TABLE]

in which

[TABLE]

hence we have |({\tilde{R}}_{t})_{i}^{\theta}|\leq\frac{2\theta^{2-\alpha}\tau^{1-\alpha}}{\Gamma(3-\alpha)}\big{(}|(R_{t})_{i}^{\frac{1}{2}}|+|(R_{v})_{i}^{\frac{1}{2}}|\big{)}\leq C_{3}\tau^{3-\alpha}.

We can use the Taylor formula to get

[TABLE]

Consider equation (1.1) at the point $(x_{i},t_{\theta})$, we have

[TABLE]

Then combining Lemma 2.1, (2.10) and (2.11), we obtain the approximation of first time level

[TABLE]

where $({\hat{R}}_{t})_{i}^{\theta}={\cal O}(\tau^{2})$.

For space discretization at each of grid points, suppose $\frac{\partial^{4}u}{\partial x^{4}}$ is continuous on $[a,b]\times[0,T]$. The Taylor expansion gives

[TABLE]

where $(R_{x})_{i}^{n}={\cal O}(h^{2})$.

2.4 The linearized scheme

To construct a stable implicit difference scheme, we need one more Lemma on the approximation of $U_{i}^{n}$, which is obtained by using Taylor expansion. This lemma plays an important role in analysis in the next section, see the equality (3.10). It reads as:

Lemma 2.5**.**

For $n\geq 1$, it holds that

[TABLE]

where

[TABLE]

and $(R_{w})_{i}^{n}={\cal O}(\tau^{2})$.

Therefore, by (2.3), (2.8), (2.9), (2.12), (2.13) and Lemma 2.5, one can approximate the equation (1.1)–(1.3) as following:

[TABLE]

in which

[TABLE]

and

[TABLE]

Denote $u_{i}^{n}$, $v_{i}^{n}$ and $w_{i}^{n}$ the numerical approximations of $U_{i}^{n}$, $V_{i}^{n}$ and $W_{i}^{n}$, respectively. Omitting the small terms $(R_{x})_{i}^{n}$, $({\hat{R}}_{t})_{i}^{n}$, $(R_{w})_{i}^{n}$ in (2.14), $({\hat{R}}_{t})_{i}^{\theta}$ in (2.15), $({\tilde{R}}_{t})_{i}^{n+1}$ in (2.18) and $({\tilde{R}}_{t})_{i}^{\theta}$ in (2.19), we obtain the following linearized difference scheme:

[TABLE]

3 Analysis of the proposed scheme

Before carrying out the convergence and stability of difference scheme (2.20)–(2.23), we first list some preliminary lemmas.

Lemma 3.1**.**

(Gronwall’s inequality [40]) Let $\{G_{n}\}$ and $\{k_{n}\}$ be nonnegative sequences satisfying

[TABLE]

where $K\geq 0$. Then

[TABLE]

Lemma 3.2**.**

([41]) Let $u\in\mathcal{V}_{h}$, it holds that

[TABLE]

Lemma 3.3**.**

([33]) If the positive sequence $\{g_{k}^{(n+1)}|0\leq k\leq n,n\geq 1\}$ is strictly decreasing and satisfies $(2\sigma-1)g_{0}^{(n+1)}-\sigma g_{1}^{(n+1)}>0$ for a constant $\sigma\in(0,1)$ . Then

[TABLE]

We need a special form of Lemma 3.3, which is stated as the following:

Lemma 3.4**.**

For any real sequence $F^{n}$, the following inequality holds:

[TABLE]

Proof.

Taking $g_{k}^{(n+1)}=d_{k}^{(n+1)}$, $y^{n}=v^{n-\theta}$ $(n\geq 1)$ and $y^{0}=v^{0}+\frac{1}{d_{n}^{(n+1)}}F^{n}$ in Lemma 3.3, and using Lemma 2.3(b),(c), we get

[TABLE]

∎

3.1 Convergence

Now we denote the error $e_{i}^{n}=U_{i}^{n}-u_{i}^{n}$, $0\leq i\leq M$, $0\leq n\leq N$, and

[TABLE]

in which

[TABLE]

Taking

[TABLE]

and subtracting (2.20)–(2.23) from (2.14)–(2.17), we obtain the following error system:

[TABLE]

where $R_{i}^{1}=({\hat{R}}_{t})_{i}^{\theta}+({\tilde{R}}_{t})_{i}^{\theta}$ and $R_{i}^{n+1}=(R_{x})_{i}^{n}-({\hat{R}}_{t})_{i}^{n}+(R_{w})_{i}^{n}+({\tilde{R}}_{t})_{i}^{n+1}$, $1\leq n\leq N-1$. Then there exists a positive constant $C_{4}$ such that

[TABLE]

Then we conclude the convergence of proposed scheme (2.20)–(2.23) as following.

Theorem 3.5**.**

*Let $u(x,t)$ be the solution of the problem (1.1)–(1.3) and smooth enough, and let $\{u_{i}^{n},0\leq i\leq M,0\leq n\leq N\}$ be the solution of the scheme (2.20)–(2.23). It holds that *

[TABLE]

*where {\bar{C}}=\frac{\exp(C_{7})}{2(1-\theta)}\Big{[}\frac{1}{2}\Big{(}\frac{(\theta-\frac{1}{2})C_{4}\Gamma(3-\alpha)}{\theta^{2-\alpha}}\Big{)}^{2}+4\Gamma(2-\alpha)T^{\alpha}\Big{(}C_{5}+C_{6}+\frac{2\Gamma(4-\alpha)T^{\alpha}}{(2-\alpha)^{2}}C_{4}^{2}\Big{)}\Big{]}^{\frac{1}{2}}, with

C_{5}=\Big{[}\frac{(3\theta-2\theta^{2}-\frac{1}{2})C_{3}\Gamma(3-\alpha)}{2\theta^{2-\alpha}}\Big{]}^{2} , $C_{6}=\frac{(1-\theta)\Gamma(3-\alpha)(3-2\theta)T^{2-\alpha}}{2\theta^{4-2\alpha}}$, and $C_{7}=\frac{\Gamma(2-\alpha)\Gamma(4-\alpha)T^{2\alpha}L^{2}}{\left[(2-\alpha)(1-\theta)\right]^{2}}$.*

Proof.

We have $\|e^{0}\|=0$ from (3.5). Now we use mathematical induction to prove

[TABLE]

It follows from (3.1), (3.3) and (3.6) that

[TABLE]

Hence (3.8) holds for $n=1$.

Suppose (3.8) is valid for $1\leq n\leq m$ ($1\leq m\leq N-1$), we then prove it is also valid for $n=m+1$.

Taking the inner product of (3.2) with

[TABLE]

we have

[TABLE]

where (R_{f})_{i}^{n+1}=-\big{[}f(U_{i}^{n})-f(u_{i}^{n})\big{]}+R_{i}^{n+1}.

It is easy to verify that

[TABLE]

Noticing ${\hat{v}}^{0}=0$ and utilizing Lemma 3.4, we get

[TABLE]

Substituting (3.12) and (3.13) into (3.11), we obtain

[TABLE]

where

[TABLE]

Summing up (3.14) for $n$ from $1$ to $m$ yield

[TABLE]

then we can deduce the following inequality

[TABLE]

It can be verified by using Cauchy-Schwarz inequality and Lemma 2.2(h) that

[TABLE]

Furthermore, the inequality $(y+z)^{2}\leq 2(y^{2}+z^{2})$ gives

[TABLE]

Consequently, it follows from (3.15)–(3.17) that

[TABLE]

where $t_{m+1}^{\alpha}\leq T^{\alpha}$ has been used, and

[TABLE]

We note that if $\|e^{m+1}\|\leq\|e^{m}\|$, (3.8) follows directly. Therefore, we only consider the situation that

[TABLE]

Then the triangular property of $L_{2}$ norm yields

[TABLE]

which implies that

[TABLE]

Combining (3.18) and (3.19), we get

[TABLE]

Then we estimate $B^{m}$ term by term. Recalling the definition of ${\hat{w}}_{i}^{1}$, a straightforward calculation shows

[TABLE]

By using Lemma 2.1, Lemma 2.2(b), (3.3) and (3.6), we have

[TABLE]

For the nonlinear term, assuming that the global Lipschitz condition (1.4) hold, we have

[TABLE]

Then utilizing Lemma 2.3(f), (3.6) and (3.23), we can conclude that

[TABLE]

Thus, (3.6), (3.9), (3.20)–(3.22) and (3.24) yield

[TABLE]

which shows that (3.8) is proved.

Consequently, we can apply Lemma 3.1 and Lemma 2.3(f) on (3.8) to conclude

[TABLE]

∎

Remark 3.6**.**

For functions which are not globally Lipschitz continuous, for example $f(u)=[u(x,t)]^{r}$, $r$ is positive integer and $r\geq 3$. Inspired by the approach for dealing with the nonlinear term in Theorem 4.2 of [37], one can also obtain convergence of the scheme by assuming $\tau=\nu h^{\frac{1}{2}+\epsilon}$, where $\nu$, $\epsilon$ are positive numbers. In fact, based on the smoothness assumption of the exact solution, there exists a positive constant $C_{0}$ such that

[TABLE]

Note that (3.7) is valid for $1\leq n\leq m$ by applying Lemma 3.1 and Lemma 2.3(f) on the inductive assumption of Theorem 3.5. If ${\bar{C}}(\nu^{2}h^{2\epsilon}+h)\leq{\tilde{C}}_{0}$, in which ${\tilde{C}}_{0}$ is a positive constant independent of $\tau$ and $h$, then it follows that

[TABLE]

and we get

[TABLE]

which implies that $u^{n}$ ($1\leq n\leq m$) is uniformly bounded. The bound (3.23) can still be obtained so long as $f$ is Lipschitz continuous on $[-{\tilde{C}}_{0}-C_{0},{\tilde{C}}_{0}+C_{0}]$. One can then follow the other parts of the proof to conclude convergence of the scheme.

3.2 Stability

Now we show the stability of propose scheme (2.20)–(2.23). Suppose that $\{{\tilde{u}}_{i}^{n},0\leq i\leq M,0\leq n\leq N\}$ is the solution of the following difference scheme:

[TABLE]

where

[TABLE]

with

[TABLE]

and

[TABLE]

Denoting the perturbation term

[TABLE]

and taking \xi_{i}=\big{[}({\tilde{\varphi}}_{xx}-\varphi_{xx})_{i}+\theta\tau({\tilde{\psi}}_{xx}-\psi_{xx})_{i}\big{]}-\big{[}f({\tilde{\varphi}}_{i}+\theta\tau{\tilde{\psi}}_{i})-f(\varphi_{i}+\theta\tau\psi_{i})\big{]}.

Then we can get the following perturbation system:

[TABLE]

where

[TABLE]

with

[TABLE]

and

[TABLE]

We have the following theorem to describe the stability of proposed scheme.

Theorem 3.7**.**

Let $\{\eta_{i}^{n},0\leq i\leq M,0\leq n\leq N\}$ be the solution of the perturbation system (3.30)–(3.33). It holds that

[TABLE]

*where ${\tilde{C}}=\exp(C_{7})$, C_{8}=\frac{1}{4(1-\theta)^{2}}\Big{\{}\Big{(}\frac{(\theta-\frac{1}{2})(b-a)\Gamma(3-\alpha)}{\theta^{2-\alpha}}\Big{)}^{2}+\\ 12\Gamma(2-\alpha)T^{\alpha}\Big{[}\Big{(}\frac{(3\theta-2\theta^{2}-\frac{1}{2})\Gamma(3-\alpha)}{\theta^{2-\alpha}}\Big{)}^{2}+\frac{2(b-a)^{2}C_{5}}{9}\Big{]}\Big{\}},

C_{9}=\frac{1}{4(1-\theta)^{2}}\Big{\{}{[(2\theta-1)^{2}(L+1)^{2}+(\frac{3}{2}-\theta)^{2}](b-a)^{2}}+\\ 12\Gamma(2-\alpha)T^{\alpha}\Big{[}{(3\theta-2\theta^{2}-\frac{1}{2})^{2}(L^{2}+1)+(\frac{3}{2}-3\theta+2\theta^{2})^{2}}+\frac{2(b-a)^{2}C_{5}L^{2}}{9}\Big{]}\Big{\}} and C_{10}=\frac{1}{4(1-\theta)^{2}}\Big{\{}(2\theta-1)^{2}\Big{[}6\Big{(}\frac{(2\theta-1)(b-a)}{2\sqrt{6}(1-\theta)}+L\Big{)}^{2}+(b-a)^{2}\Big{]}+12\Gamma(2-\alpha)T^{\alpha}\Big{[}\big{(}\frac{(2\theta-1)(3\theta-2\theta^{2}-\frac{1}{2})}{2-2\theta}\big{)}^{2}(L^{2}+1)+(3\theta-2\theta^{2}-1)^{2}+\frac{2(b-a)^{2}C_{5}L^{2}}{9}+\frac{(b-a)^{2}T^{2-\alpha}}{9\Gamma(3-\alpha)}\Big{]}\Big{\}}.*

Proof.

Obviously, (3.36) is valid for $n=0$. We use mathematical induction once again to prove

[TABLE]

It follows from (3.31), (3.34), and Lemma 3.2 that

[TABLE]

Note that

[TABLE]

and then (for sufficiently small $\tau$)

[TABLE]

Combining (3.38) and (3.39), we get

[TABLE]

Since $\sqrt{C_{9}}>\frac{\sqrt{6}(b-a)}{3}(L+1)$, then (3.37) is valid for $n=1$.

Suppose (3.37) hold for $1\leq n\leq m$ ($1\leq m\leq N-1$). Now we show (3.37) also hold for $n=m+1$.

Taking the inner product of (3.30) by

[TABLE]

we get

[TABLE]

where ({\tilde{R}}_{f})_{i}^{n+1}=-\big{[}f({\tilde{u}}_{i}^{n})-f(u_{i}^{n})\big{]}.

Then following the similar methodology in the proof of Theorem 3.5, we can obtain

[TABLE]

where

[TABLE]

Applying (3.35), (3.34), (3.39) and Cauchy-Schwarz inequality, we have

[TABLE]

Note that

[TABLE]

and the globally Lipschitz continuity (1.4) yields

[TABLE]

and then

[TABLE]

So we can derive form combining (3.39)–(3.45) and Lemma 3.2 that

[TABLE]

which shows that (3.37) is proved.

Therefore, applying Lemma 3.1 and Lemma 2.3(f) on (3.37), we finally get

[TABLE]

∎

Remark 3.8**.**

For nonlinear terms which are locally Lipschitz continuous as discussed in Remark 3.6, we assume that

[TABLE]

where $C_{11}$ and $\delta$ are positive constants. Note that (3.36) is valid for $1\leq n\leq m$ by applying Lemma 3.1 and Lemma 2.3(f) on (3.37). If ${\tilde{C}}C_{11}\sqrt{C_{8}+C_{9}+C_{10}}h^{\delta}\leq{\tilde{C}}_{0}$, it follows that

[TABLE]

Noticing (3.25), then

[TABLE]

which implies that ${\tilde{u}}^{n}$ is uniformly bounded. If $f$ is Lipschitz continuous on an interval containing $[-2{\tilde{C}}_{0}-C_{0},2{\tilde{C}}_{0}+C_{0}]$, one can still obtain the desired conclusion.

4 Compact scheme

In this section, we propose the spacial fourth-order scheme for the problem (1.1)–(1.3). For any $u\in\mathcal{V}_{h}$, we define the spatial high order operator ${\cal A}$ as follow:

[TABLE]

and we define the corresponding norm:

[TABLE]

Then it is easy to check that

[TABLE]

By Taylor expansion, if $\frac{\partial^{6}u}{\partial x^{6}}$ is continuous on $[x_{i-1},x_{i+1}]$, we have

[TABLE]

Performing the compact operator ${\cal A}$ on both sides of (2.2) and following the derivation of the difference scheme (2.20)–(2.23) in section 2, we obtain the compact difference linearized scheme for (1.1)–(1.3):

[TABLE]

Using similar approach with the proof of Theorem 3.7 in [42] or Lemma 4.2 in [39], we have the following Lemma.

Lemma 4.1**.**

For any real sequence $F^{n}$, the following estimate holds:

[TABLE]

Then, with the help of (4.1), (4.2) and Lemma 4.1, and under the same assumptions in Theorem 3.5 and Theorem 3.7, the theoretical results can be obtained following similar arguments in the proof of Theorem 3.5 and Theorem 3.7. We present the convergence conclusion in the following.

Theorem 4.2**.**

Let $u(x,t)$ be the solution of the problem (1.1)–(1.3) and smooth enough, and let $\{u_{i}^{n},0\leq i\leq M,0\leq n\leq N\}$ be the solution of the scheme (4.3)–(4.6). If (1.4) holds globally, we then have

[TABLE]

where {\hat{C}}=\frac{3\exp(\frac{9}{4}C_{7})}{4(1-\theta)}\Big{[}\frac{1}{2}\Big{(}\frac{(\theta-\frac{1}{2})C_{4}\Gamma(3-\alpha)}{\theta^{2-\alpha}}\Big{)}^{2}+4\Gamma(2-\alpha)T^{\alpha}\Big{(}C_{5}+C_{6}+\frac{2\Gamma(4-\alpha)T^{\alpha}}{(2-\alpha)^{2}}C_{4}^{2}\Big{)}\Big{]}^{\frac{1}{2}}.

Remark 4.3**.**

As in Remark 3.6, we can show that (4.7) holds for locally Lipschitz continuous nonlinear term, if we assume ${\hat{C}}(\nu^{2}h^{2\epsilon}+h^{3})\leq{\tilde{C}}_{0}$.

5 Numerical experiments

In this section, we carry out numerical experiments for the proposed finite difference schemes (2.20)–(2.23) and (4.3)–(4.6) to illustrate our theoretical statements. All our tests were done in MATLAB R2014a with a desktop computer (Dell optiplex 7020, configuration: Intel(R) Core(TM) i7-4790 CPU 3.60GHz and 16.00G RAM).

The $L_{2}$ norm errors between the exact and the numerical solutions

[TABLE]

are shown in the following tables. In the tables,

[TABLE]

is used to denote the temporal convergence order for sufficiently small $h$, and

[TABLE]

is the spatial convergence order for sufficiently small $\tau$.

5.1 Accuracy verification

We consider the problem (1.1)–(1.3) for $x\in[0,1]$, $T=1$ and the forcing term

[TABLE]

is chosen to such that the exact solution is $u(x,t)=\sin(\pi x)(t^{4}+1)$, where

[TABLE]

The numerical results for the above three cases by applying difference scheme (2.20)–(2.23) were recorded in Table 1 and Table 2, while the results for the three cases by applying compact scheme (4.3)–(4.6) were presented in Table 3 and Table 4.

Table 1 reports the numerical results in time direction of the proposed scheme (2.20)–(2.23) with fixed $h=\frac{1}{1000}$ and different choices of $\alpha$ are taken. Meanwhile, the numerical results in space direction with fixed $\tau=\frac{1}{1000}$ and different $\alpha$ are listed in Table 2. From these two tables, one can see that the convergence rate for the three cases are both 2 in time and space, which are in accordance with the theoretical statements.

Table 3 lists the numerical results of the compact scheme (4.3)–(4.6) in time direction with fixed $h=\frac{1}{100}$, and Table 4 shows the numerical results in space direction with fixed $\tau=\frac{1}{5000}$ and different $\alpha$. The second-order accuracy in time and fourth-order accuracy in space are apparent in these two tables, respectively. We note that, in our implementation to the the non-global Lipschitz continuous term $f(u)=u^{3}$, the condition $\tau=\nu h^{\frac{1}{2}+\epsilon}$ theoretically imposed in Remark 3.6 is not necessary.

5.2 Comparison with other numerical schemes

In this subsection, we give some comparisons between our proposed method and some standard methods for solving the equation (1.1)–(1.3). It will be shown that our linearized schemes have advantages in both theoretical analysis and numerical computation.

One may construct the numerical schemes where time fractional derivative was approximate by the widely used classical $L_{1}$ formula [43] to solve this kind of equations. The scheme will take the form as

[TABLE]

where $a_{k}=(k+1)^{2-\alpha}-k^{2-\alpha}$.

There are two commonly used approaches to approximate the nonlinear term $f(u_{i}^{n+\frac{1}{2}})$. The first is the central approximation, i.e. f(u_{i}^{n+\frac{1}{2}})=\big{[}f(u_{i}^{n+1})+f(u_{i}^{n})\big{]}/2+{\cal O}(\tau^{2}), and the resulting scheme will be accurate of temporal order ${\cal O}(\tau^{3-\alpha})$, e.g. [29, 30, 31], but this leads to nonlinear treatment so that iterative methods are required. However, the iterative methods would cause additional computational costs. Here we give a comparison by a numerical example, we apply our proposed scheme (2.20)–(2.23) and the scheme (5.1) with the above nonlinear approximation to solve the example in subsection 5.1, in which the nonlinear scheme is dealt with a fixed-point method. For simplicity of presentation, we only list results of Case 2 in Table 5, remarking that similar results can be obtained for other cases and different parameters. From this table, one can clearly see that our scheme works more efficiently and spend less CPU (seconds) time than the iterative method. Moreover, when applying iterative methods for solving nonlinear schemes, the convergence of the iterative methods is usually not easy to be established. The second way to approximate nonlinear term is a linear approach, e.g. $f(u_{i}^{n+\frac{1}{2}})=f(u_{i}^{n})+{\cal O}(\tau)$. The stability and convergence can be verified like [32], but its convergence order is only ${\cal O}(\tau)$ or no more than ${\cal O}(\tau^{3-\alpha})$ in time.

One may also use the formula in Lemma 2.1 to replace the $L_{1}$ formula and then achieve a ${\cal O}(\tau^{2})$ in fractional approximation, but this gives rise to some difficulties in theoretical analysis and, to our knowledge, the linear method like that in [32] may not be applied. Consequently, iterative methods are required.

All the discussions above illustrate that our proposed linearized schemes would be favorable to some usual numerical schemes.

6 Concluding remarks

In this paper, we consider a nonlinear fractional order Klein-Gordon type equation. We proposed a linearized finite difference scheme to solve the problem numerically. The main advantage is that the nonlinear term is evaluated on previous time level. As a result, iterative method is not needed for implementation. However, shifting the evaluation to previous level causes difficulties in theoretical analysis. Inspired by some recent studies, we show that the scheme converges with second-order in time. The results are justified by some numerical tests.

Acknowledgment

The authors would like to thank the referees for their comments which improve the paper significantly. We also want to thank Prof. Honglin Liao for helpful discussion on the estimates of Lemmas 2.2 and 2.3.