Moving grids method has important applications in a variety of physical and engineering areas such as solid and fluid dynamics, combustion, heat transfer, material science, etc. This method is more efficient than the fixed grids and does not increase computing cost. We usually apply the finite element methods to the spatial domain, but choose difference methods with respect to the time variable for solving evolution partial differential equations. At the same time, different meshes of domain are used at different time level.

Several moving grids techniques were studied. Such as [1] considered the moving grids finite element method; [2] and [3] constructed and analyzed this method for the oil-water two-phase displacement problem; [4-8] analyzed the parabolic, Stokes problems, parabolic integro-differential equations, generalized nerve conductive equations and fractional diffusion equations with moving grids nonconforming finite element scheme respectively. But the analysis in the above studies relies on the regular condition or quasi-uniform assumption for meshes.

The Klein-Gordon equation is the most basic equation used in relativistic quantum mechanics and quantum field theory to describe a spin-zero particle. The equation is closely related to the physical problem and plays an important role in the study of soliton. In [9], authors studied the existence of a unique global solution under the condition that the parameter is small enough. In [10], a display difference scheme was established for one-dimensional Klein-Gordon equation of unbounded region, and the results of stability and convergence of the scheme were obtained by the energy analysis method. In [11], the numerical solution of one-dimensional Klein-Gordon equation was studied. However, the finite element method for the Klein-Gordon equation is rare.

In this paper, we mainly focus on the convergence theory, the finite element method of moving grids is introduced, and the Crank-Nicolson discrete scheme of the nonlinear Klein-Gordon equation is analyzed without requiring the subdivision to satisfy the regular hypothesis, and the corresponding optimal error estimation of the moving grid approach is derived. It is worth mentioning that, in the usual finite element method of moving meshes, it is necessary to use the Riesz projection to approximate the solution of the original problem, and this paper makes use of the particularity of the element structure, that is, $ u-\Pi u $ and the elements in the finite element space are orthogonal in the sense of energy mode, and the Riesz projection is used to simplify the proof process of the previous documents.

For the sake of convenience, like [12], let $ \Omega \subset R^{2} $ be a rectangular domain with boundary $ \partial\Omega $ parallel to the $ x $-axis or $ y $-axis in the plane, Let $ \Gamma_{h} $ be a family of rectangular subdivisions, i.e., $ \overline{\Omega} = \bigcup\limits_{K\in \Gamma_{h}}K $. For the general element $ K $, we denote the lengths of edges parallel to $ x $-axis and $ y $-axis by $ 2h_{x} $ and $ 2h_{y} $, respectively, and the barycenter of element $ K $ by $ (x_{K}, y_{K}) $, the four vertices of $ K $ are $ a_{1} = (x_{K}-h_{x}, y_{K}-h_{y}) $, $ a_{2} = (x_{K}+h_{x}, y_{K}-h_{y}) $, $ a_{3} = (x_{K}+h_{x}, y_{K}+h_{y}) $ and $ a_{4} = (x_{K}-h_{x}, y_{K}+h_{y}) $, respectively, and the four edges are $ l_{i} = \overline{a_{i}a_{i+1}} $ $ (i = 1, 2, 3, 4 $ (mod 4)), $ h_{K} = \max\limits_{K\in\Gamma_{h}}\{h_{x}, h_{y}\}, \ h = \max\limits_{K\in\Gamma_{h}}h_{K}. $

Let $ \hat{K} = [-1, 1]\times[-1, 1] $ be reference element on the $ \xi-\eta $ plane, the four vertices of $ \hat{K} $ are $ \hat{a}_{1} = (-1, -1) $, $ \hat{a}_{2} = (1, -1) $, $ \hat{a}_{3} = (1, 1) $ and $ \hat{a}_{4} = (-1, 1) $, respectively, $ \hat{l}_{i} = \overline{\hat{a}_{i}\hat{a}_{i+1}} $ $ (i = 1, 2, 3, 4 $ (mod4)) be the four edges of $ \hat{K} $.

We define the finite element $ (\hat{K}, \hat{P}, \hat{\sum}) $ (see [13])

$ \hat{\sum} = \{\hat{v}_{1}, \hat{v}_{2}, \hat{v}_{3}, \hat{v}_{4}, \hat{v}_{5}\}, \ \ \ \hat{P} = \text{span}\{1, \xi, \eta, \varphi(\xi), \varphi(\eta)\}, $

where $ \hat{v}_{i} = \frac{1}{|\hat{l}_{i}|}\int_{\hat{l}_{i}}\hat{v}d\hat{s}, i = 1, 2, 3, 4, \hat{v}_{5} = \frac{1}{|\hat{K}|}\int_{\hat{K}}\hat{v}d\xi d\eta, \varphi(t) = \frac{1}{2}(3t^{2}-1). $

$ \forall\hat{v}\in H^{1}(\hat{K}), $ interpolation functions $ \hat{\Pi}\hat{v} $ can be expressed as

$ \begin{equation} \hat{\Pi}\hat{v} = \hat{v}_{5}+\frac{1}{2}(\hat{v}_{2}-\hat{v}_{4})\xi +\frac{1}{2}(\hat{v}_{3}-\hat{v}_{1})\eta+\frac{1}{2}(\hat{v}_{2}+\hat{v}_{4}-2\hat{v}_{5})\varphi(\xi) +\frac{1}{2}(\hat{v}_{3}+\hat{v}_{1}-2\hat{v}_{5})\varphi(\eta). \end{equation} $

(2.1)

Lemma 2.1 [14] The interpolation operator $ \hat{\Pi} $ defined as $ (2.1) $ has the anisotropic property, i.e., for all $ \hat{\phi}\in H^{2}(\hat{K}) $ and $ \alpha = (a_{1}, a_{2}) $, when $ |\alpha| = 1, $ there holds

$ \begin{equation} \|\widehat{D}^{\alpha}(\hat{\phi}-\hat{\Pi}^{c}\hat{\phi})\|_{0, \widehat{K}} \leq C|\widehat{D}^{\alpha}\hat{\phi}|_{1, \widehat{K}}. \end{equation} $

(2.2)

Here and later, we use $ C $ denotes a positive constant real number independent of $ h_{K} $ and $ \frac{h_{K}}{\rho_{K}}, $ possibly different at different occurrences.

Define the affine mapping $ F_{K}:\hat{K}\rightarrow K $ by

$ \begin{equation} \left\{\begin{array}{ll} x = x_{K}+h_{x}\xi, \\ y = y_{K}+h_{y}\eta. \end{array} \right. \end{equation} $

(2.3)

Then $ \forall K\in\Gamma_{h} $, the associated finite element spaces $ V_{h} $ can be defined as

$ \begin{equation} V_{h} = \{v_{h}|\hat{v} = v_{h}|_{K}\circ F_{K}\in \hat{P}, \int_{F}[v_{h}]ds = 0, F\subset \partial K\}, \end{equation} $

(2.4)

where $ [v_{h}] $ stands for the jump of $ v_{h} $ across the edge $ F $ if $ F $ is an internal edge, and it is equal to $ v_{h} $ itself if $ F $ belongs to $ \partial \Omega $.

We define the following interpolations $ \Pi $ over spaces $ V_{h} $ as

$ \Pi:H^{2}(\Omega)\rightarrow V_{h}, \Pi|_{K} = \Pi_{K} = \hat{\Pi}\circ F_{K}^{-1}. $

It is easy to show, $ \|\cdot\|_{h} = (\sum\limits_{K\in {\mathcal{T}}_{h}}|\cdot|^{2}_{1, K})^{\frac{1}{2}} $ is the norm over $ V_{h} $.

We consider the nonlinear Klein-Gordon equation

$ \begin{equation} \left\{\begin{array}{ll} u_{tt}+\alpha u_{t}-\gamma\triangle u+g(u) = f(X, t), &(X, t)\in \Omega\times(0, T], \\ u(X, 0) = u_{0}(X),& X\in \Omega, \\ u_{t}(X, 0) = \varphi_{0}(X),& X\in \Omega, \\ u(X, t) = u_{t}(X, t) = 0, &(X, t)\in \partial\Omega\times(0, T], \end{array} \right. \end{equation} $

(3.1)

where $ X = (x, y) $, $ \alpha>0 $, $ \gamma>0 $, $ g(u) $ satisfies the Lipschitz continuous conditionon on the variable $ u $, and has the second order bounded partial derivative.

Let $ u_{t} = Q $, $ (3.1) $ is equivalent to the following question

$ \begin{equation} \left\{\begin{array}{ll} \frac{\partial Q}{\partial t}+\alpha Q-\gamma\triangle u+g(u) = f(X, t),& (X, t)\in \Omega\times(0, T], \\ \frac{\partial u}{\partial t} = Q,& (X, t)\in \Omega\times[0, T], \\ u(X, 0) = u_{0}(X), & X\in \Omega, \\ Q(X, 0) = \varphi_{0}(X),& X\in \Omega, \\ u(X, t) = Q(X, t) = 0,& (X, t)\in \partial\Omega\times(0, T]. \end{array} \right. \end{equation} $

(3.2)

The variational formulation for problem $ (3.2) $ is written as: $ \forall v\in H_{0}^{1}(\Omega) $,

$ \begin{equation} \left\{\begin{array}{ll} ( \frac{\partial Q}{\partial t}, v)+(\alpha Q, v)+a(u, v)+(g(u), v) = (f, v), \\ ( \frac{\partial u}{\partial t}, v) = (Q, v), \\ (u(0)-u_{0}, v) = 0, \\ (Q(0)-\varphi_{0}, v) = 0, \end{array} \right. \end{equation} $

(3.3)

where $ u(0) = u(X, 0), Q(0) = Q(X, 0), a(u, v) = \gamma\int_{\Omega}\nabla u\nabla vdX $.

Then the approximation problem corresponding to $ (3.3) $ reads as: find $ u_{h}, Q_{h}\in V_{h}, \forall v_{h}\in V_h $, such that

$ \begin{equation} \left\{\begin{array}{ll} ( \frac{\partial Q_{h}}{\partial t}, v_{h})+(\alpha Q_{h}, v_{h})+a_{h}(u_{h}, v_{h})+(g(u_{h}), v_{h}) = (f, v_{h}), \\ ( \frac{\partial u_{h}}{\partial t}, v_{h}) = (Q_{h}, v_{h}), \\ (u_{h}(0)-\Pi u_{0}, v_{h}) = 0, \\ (Q_{h}(0)-\Pi\varphi_{0}, v_{h}) = 0, \end{array} \right. \end{equation} $

(3.4)

where $ a_{h}(u_{h}, v_{h}) = \gamma\sum\limits_{K}\int_{K}\nabla u_{h}\cdot\nabla v_{h}dX. $

In this section we apply the idea of moving grids to problem $ (3.4) $ and develop the Crank-Nicolson discretization scheme for anisotropic finite element. Let $ 0 = t_{0}<t_{1}<\cdot \cdot \cdot< t_{N} = T $ be a partition of the time interval [0, T], $ \triangle t = \frac{T}{N}, \ t_{n} = n\triangle t \ (n = 0, 1, 2, \cdots, N) $, $ V^{h}_{n} = \{v (X, t_{n});v\in V_{h}\} $. We choose the approximating space $ S $ of $ u (X, t) $ in the following way: the approximation solution $ u^{h} (X, t)\in S $ is the piecewise linear function with respect to the time subdivisions $ 0 = t_{0}<t_{1}<\cdot \cdot \cdot< t_{N} = T $ based on $ u^{h} (x, t_{n})\in V^{h}_{n}\; (n = 0, 1, 2, \cdot\cdot\cdot, N) $.

Now, we introduce the Crank-Nicolson discretization scheme of anisotropic finite element to determine the function values $ u^{h}_{n} = u^{h}(X, t_{n}) $ as follows

$ \begin{eqnarray} (\widehat{u}^{h}_{0}-u_{0}, v_{h})& = &0, \forall v_{h}\in V^{h}_{1}\; \; \; \; \; \text{for} \ n = 0, \end{eqnarray} $

(3.5)

$ \begin{eqnarray} (\widehat{Q}^{h}_{0}-Q_{0}, v_{h})& = &0, \forall v_{h}\in V^{h}_{1}\; \; \; \; \; \text{for} \ n = 0, \end{eqnarray} $

(3.6)

$ \begin{eqnarray} (\widehat{u}^{h}_{n}-u^{h}_{n}, v_{h})& = &0, \forall v_{h}\in V^{h}_{n+1}\; \; \text{for} \ n>0, \end{eqnarray} $

(3.7)

$ \begin{eqnarray} (\widehat{Q}^{h}_{n}-Q^{h}_{n}, v_{h})& = &0, \forall v_{h}\in V^{h}_{n+1}\; \; \text{for} \ n>0, \end{eqnarray} $

(3.8)

$ \begin{equation} (Q^{h}_{n+1}-\widehat{Q}^{h}_{n}, v_{h})+(\alpha Q^{h}_{n+\frac{1}{2}}, v_{h})\triangle t+a_{h}(u^{h}_{n+\frac{1}{2}}, v_{h})\triangle t +(g(u_{n+\frac{1}{2}}), v_{h})\triangle t = (f_{n+\frac{1}{2}}, v_{h})\triangle t, \; \forall v_{h}\in V^{h}_{n+1}, \end{equation} $

(3.9)

$ \begin{equation} (\frac{u^{h}_{n+1}-\widehat{u}^{h}_{n}}{\triangle t}, v_{h}) = (\frac{Q^{h}_{n+1}+\widehat{Q}^{h}_{n}}{2}, v_{h}), \end{equation} $

(3.10)

where $ u^{h}_{n} = u^{h}(X, t_{n}), \ \ u^{h}_{n+\frac{1}{2}} = \frac{1}{2}(\widehat{u}^{h}_{n}+u^{h}_{n+1}), \ \ Q^{h}_{n+\frac{1}{2}} = \frac{1}{2}(\widehat{Q}^{h}_{n}+Q^{h}_{n+1}), \ \ f_{n+\frac{1}{2}} = \frac{1}{2}(f(X, t_{n})+f(X, t_{n+1})). $ $ (3.7), (3.8) $ implies that $ \widehat{u}^{h}_{n} = u^{h}_{n}, \widehat{Q}^{h}_{n} = Q^{h}_{n} $ when $ V^{h}_{n} = V^{h}_{n+1}. $ Furthermore, $ (3.7), (3.8) $ is a $ L^{2} $-projection modification scheme for the former space when the two spaces $ V^{h}_{n} $ and $ V^{h}_{n+1} $ have different meshes or different interpolation functions. $ (3.9) $ is the general trapezoid difference scheme. We get $ \widehat{u}^{h}_{n}, \widehat{Q}^{h}_{n} $ from $ u^{h}_{n}, Q^{h}_{n} $ in $ (3.7), (3.8) $ and get $ u_{n+1}^{h}, Q_{n+1}^{h} $ from $ \widehat{u}^{h}_{n}, \widehat{Q}^{h}_{n} $ in $ (3.9), (3.10) $. So by partial differential equation theory $ u^{h}_{n} $ and $ Q^{h}_{n} $ can be determined uniquely through $ (3.5) $–$ (3.10) $.

The main error between the solution $ u(X, t) $ and the approximation solution $ u^{h}(X, t) $ consists of three parts: the interpolation error with respect to the finite element method, the difference error with respect to the time, and the error of moving grids.

Lemma 4.1 [13] For anisotropic meshes, $ \forall u\in H^{2}(\Omega), $ then we have

$ \begin{equation} \Big|\sum\limits_K \int\nolimits_{\partial K}\frac{\partial u}{\partial n}\cdot v_h ds\Big| \leq Ch|u|_2\|v_h\|_h, \ \ \ \forall\ v_h\in V_h, \end{equation} $

(4.1)

here and later $ n = (n_{x}, n_{y}) $ denotes the outward unit normal vector to $ \partial K $.

Lemma 4.2 [13] $ \forall u\in H^{1}(\Omega), $ there hold $ a_{h}(u-\Pi u, v_{h}) = 0, \ \ \forall v_{h}\in V_h. $

Lemma 4.3 [13] $ \|v_{h}\|_0\leq C\|v_{h}\|_h, \ \ \ \forall\ v_{h}\in V_{h}. $

Lemma 4.4 Let $ u_{n} = u(X, t_{n}) $, $ \forall v_{h}\in V^{h}_{n+1}, $then there holds

$ \begin{equation} (Q_{n+1}-Q_{n}, v_{h})+(\alpha Q_{n+\frac{1}{2}}, v_{h})\triangle t+a_{h}(u_{n+\frac{1}{2}}, v_{h})\triangle t+(g(u_{n+\frac{1}{2}}), v_{h})\triangle t = (f_{n+\frac{1}{2}}, v_{h})\triangle t+E_{n}(v_{h}), \end{equation} $

(4.2)

where

$ \begin{eqnarray} |E_{n}(v_{h})|&\leq& C [ ( \int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}f}{\partial t^{2}} \|^{2}_{0}d\tau )^{\frac{1}{2}} + (\int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}Q}{\partial t^{2}} \|^{2}_{0}d\tau )^{\frac{1}{2}}\\ &&+ ( \int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{1}d\tau )^{\frac{1}{2}} + ( \int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{0}d\tau )^{\frac{1}{2}} ](\triangle t_{n})^{\frac{5}{2}}\|v_{h}\|_{h} +Ch\triangle t\|v_{h}\|_{h}. \end{eqnarray} $

(4.3)

Proof From $ (3.2) $, for all $ v\in V^{h}_{n+1} $,

$ \begin{equation} \begin{array}{ll} &(Q_{n+1}-Q_{n}, v_{h})+( \int_{t_{n}}^{t_{n+1}}\alpha Qd\tau, v_{h})+a_{h}( \int_{t_{n}}^{t_{n+1}}ud\tau, v_{h})+( \int_{t_{n}}^{t_{n+1}}g(u)d\tau, v_{h})\\ = &( \int_{t_{n}}^{t_{n+1}}fd\tau, v_{h})+ \int_{t_{n}}^{t_{n+1}}\Gamma_{h}(u, v_{h})d\tau. \end{array} \end{equation} $

(4.4)

Combining $ (4.2) $ and $ (4.4) $, we get

$ \begin{equation} \begin{array}{ll} E_{n}(v_{h}) = &( \int_{t_{n}}^{t_{n+1}}(f-f_{n+\frac{1}{2}})d\tau, v_{h}) -( \int_{t_{n}}^{t_{n+1}}\alpha (Q-Q_{n+\frac{1}{2}})d\tau, v_{h}) -a_{h}( \int_{t_{n}}^{t_{n+1}}(u- u_{n+\frac{1}{2}})d\tau, v_{h})\\ &-( \int_{t_{n}}^{t_{n+1}}(g(u)- g(u_{n+\frac{1}{2}}))d\tau, v_{h}) + \int_{t_{n}}^{t_{n+1}}\Gamma_{h}(u, v_{h})d\tau, \end{array} \end{equation} $

(4.5)

where $ \Gamma_{h}(u, v_{h}) = \sum\limits_{K}\int_{\partial K}\frac{\partial u}{\partial n}v_{h}ds $.

By Cauchy-Schwartz inequality and one-dimensional linear interpolation theory, we deduce that

$ \begin{eqnarray} | (\int^{t_{n+1}}_{t_{n}}(f-f_{n+\frac{1}{2}})d\tau, v_{h} ) | \leq C | (\int^{t_{n+1}}_{t_{n}}\frac{\partial^{2}f}{\partial t^{2}}(\triangle t)^{2}d\tau, v_{h} ) | \leq C (\int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}f}{\partial t^{2}} \|_{0}^{2}d\tau )^{\frac{1}{2}}(\triangle t)^{\frac{5}{2}}\|v_{h}\|_{h}. \end{eqnarray} $

Similarly,

$ \begin{eqnarray*} | ( \int^{t_{n+1}}_{t_{n}}\alpha(Q-Q_{n+\frac{1}{2}})d\tau, v_{h} ) | &\leq& C (\int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}Q}{\partial t^{2}} \|_{0}^{2}d\tau )^{\frac{1}{2}}(\triangle t)^{\frac{5}{2}}\|v_{h}\|_{h}, \\ |a_{h} ( \int^{t_{n+1}}_{t_{n}}(u-u_{n+\frac{1}{2}})d\tau, v_{h} ) | &\leq& C (\int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}u}{\partial t^{2}} \|_{1}^{2}d\tau )^{\frac{1}{2}}(\triangle t)^{\frac{5}{2}}\|v_{h}\|_{h}, \\ | ( \int^{t_{n+1}}_{t_{n}}(g(u)-g(u_{n+\frac{1}{2}}))d\tau, v_{h} ) | &\leq& C (\int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2}u}{\partial t^{2}} \|_{0}^{2}d\tau )^{\frac{1}{2}}(\triangle t)^{\frac{5}{2}}\|v_{h}\|_{h}. \end{eqnarray*} $

Applying Lemma 4.1, we obtain

$ \begin{eqnarray} | \int_{t_{n}}^{t_{n+1}}\Gamma_{h}(u, v_{h})d\tau| \leq C(u)h\triangle t\|v_{h}\|_{h}. \end{eqnarray} $

Then (4.3) comes from the combination of above inequalities. Let

$ \begin{array}{ll} \nabla\xi_{0} = 0, \qquad \nabla\eta_{0} = 0, \qquad \lambda_{0} = 0, \qquad \theta_{0} = 0\\ \xi_{n} = u^{h}_{n}-\Pi_{n}u_{n}, \qquad \eta_{n} = u_{n}-\Pi_{n}u_{n}, \qquad n = 1, 2, \cdots, N, \\ \widehat{\xi}_{n} = \widehat{u}^{h}_{n}-\Pi_{n+1}u_{n}, \quad \widehat{\eta}_{n} = u_{n}-\Pi_{n+1}u_{n}, \quad n = 0, 1, 2, \cdots, N, \\ \lambda_{n} = Q^{h}_{n}-\Pi_{n}Q_{n}, \quad\, \, \theta_{n} = Q_{n}-\Pi_{n}Q_{n}, \quad\, \, \, \, n = 1, 2, \cdots, N, \\ \widehat{\lambda}_{n} = \widehat{Q}^{h}_{n}-\Pi_{n+1}Q_{n}, \, \, \widehat{\theta}_{n} = Q_{n}-\Pi_{n+1}Q_{n}, \, \, \, \, n = 0, 1, 2, \cdots, N. \end{array} $

Theorem 4.1 For the approximate solution $ u^h_{n}, Q^h_{n} $ and the solution $ u_{n}, Q_{n} $ of $ (3.5)-(3.10) $, the following optimal error estimate holds

$ \begin{equation} \max\limits_{1\leq n\leq N}(\|Q^{h}_{n}-Q_{n}\|^{2}_{0}+\|u^{h}_{n}-u_{n}\|^{2}_{h}) \leq C(h^{2}+(\triangle t)^{4}). \end{equation} $

(4.6)

Proof By definitions of $ \xi_{n}, \eta_{n}, \lambda_{n}, \theta_{n} $, combining $ (3.9) $ and $ (4.2) $, we have for all $ v_{h}\in V_{n+1}^{h}, $

$ \begin{equation} \begin{array}{ll} &(\lambda_{n+1}-\widehat{\lambda}_{n}, v_{h})+ (\alpha\frac{\lambda_{n+1}+\widehat{\lambda}_{n}}{2}, v_{h} )\triangle t + a_{h} (\frac{\xi_{n+1}+\widehat{\xi}_{n}}{2}, v_{h} )\triangle t + (g(u^{h}_{n+\frac{1}{2}})-g(u_{n+\frac{1}{2}}), v_{h} )\triangle t\\ = &(\theta_{n+1}-\widehat{\theta}_{n}, v_{h})+ (\alpha\frac{\theta_{n+1}+\widehat{\theta}_{n}}{2}, v_{h} )\triangle t + a_{h} (\frac{\eta_{n+1}+\widehat{\eta}_{n}}{2}, v_{h} )\triangle t-E_{n}(v_{h}). \end{array} \end{equation} $

(4.7)

Applying Lemma 4.2, we obtain $ a_{h} (\frac{\eta_{n+1}+\widehat{\eta}_{n}}{2}, v_{h} )\triangle t = 0. $ By $ (3.10) $, we have

$ \begin{equation} \frac{\lambda_{n+1}+\widehat{\lambda}_{n}}{2} = \frac{\xi_{n+1}-\widehat{\xi}_{n}}{\triangle t} +\Pi_{n+1} (\frac{u_{n+1}-u_{n}}{\triangle t}-\frac{Q_{n+1}+Q_{n}}{2} ). \end{equation} $

(4.8)

Substituting $ (4.8) $ into $ (4.7) $ with $ v_{h} = \lambda_{n+1}+\widehat{\lambda}_{n} $, we deduce that

$ \begin{equation} \begin{array}{ll} &\|\lambda_{n+1}\|_{0}^{2}-\|\widehat{\lambda}_{n}\|_{0}^{2} +\frac{\alpha}{2}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t +\|\xi_{n+1}\|_{h}^{2}-\|\widehat{\xi}_{n}\|_{h}^{2}\\ = &(\theta_{n+1}-\widehat{\theta}_{n}, \lambda_{n+1}+\widehat{\lambda}_{n}) + (\alpha\frac{\theta_{n+1}+\widehat{\theta}_{n}}{2}, \lambda_{n+1}+\widehat{\lambda}_{n} )\triangle t - (g(u^{h}_{n+\frac{1}{2}})-g(u_{n+\frac{1}{2}}), \lambda_{n+1}+\widehat{\lambda}_{n} )\triangle t\\ &+ a_{h} (\xi_{n+1}+\widehat{\xi}_{n}, \Pi_{n+1} (\frac{u_{n+1}-u_{n}}{\triangle t} -\frac{Q_{n+1}+Q_{n}}{2} ) )\triangle t-E_{n}(\lambda_{n+1}+\widehat{\lambda}_{n}). \end{array} \end{equation} $

(4.9)

Applying Cauchy-Schwartz inequality and Young inequality to yield

$ \begin{equation} \begin{array}{ll} & |(\theta_{n+1}-\widehat{\theta}_{n}, \lambda_{n+1}+\widehat{\lambda}_{n}) | \leq\|\theta_{n+1}-\widehat{\theta}_{n}\|_{0}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}\\ \leq& C\frac{\|\theta_{n+1}-\widehat{\theta}_{n}\|_{0}^{2}}{\triangle t}+C\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t\\ \leq& Ch^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d \tau +C(\|\lambda_{n+1}\|_{0}^{2}+\|\widehat{\lambda}_{n}\|_{0}^{2})\triangle t, \end{array} \end{equation} $

(4.10)

$ \begin{equation} \begin{array}{ll} | (\alpha\frac{\theta_{n+1}+\widehat{\theta}_{n}}{2}, \lambda_{n+1}+\widehat{\lambda}_{n} )\triangle t | &\leq C\|\theta_{n+1}+\widehat{\theta}_{n}\|_{0}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}\triangle t\\ &\leq Ch^{4}\triangle t+\frac{\alpha}{8}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t. \end{array} \end{equation} $

(4.11)

Using Lemmas 4.3, we deduce that

$ \begin{eqnarray} && | (g(u^{h}_{n+\frac{1}{2}})-g(u_{n+\frac{1}{2}}), \lambda_{n+1}+\widehat{\lambda}_{n} )\triangle t | \leq C\|u^{h}_{n+\frac{1}{2}}-u_{n+\frac{1}{2}}\|_{0}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}\triangle t\\ &\leq& C(\|\xi_{n+\frac{1}{2}}\|+\|\eta_{n+\frac{1}{2}}\|_{0}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0})\triangle t\\ &\leq& C(\|\nabla\xi_{n+1}\|_{0}^{2}+\|\nabla\xi_{n}\|_{0}^{2})\triangle t+Ch^{4}\triangle t+\frac{\alpha}{8}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t. \end{eqnarray} $

(4.12)

$ \begin{eqnarray} && |a_{h} (\xi_{n+1}+\widehat{\xi}_{n}, \Pi_{n+1} (\frac{u_{n+1}-u_{n}}{\triangle t} -\frac{Q_{n+1}+Q_{n}}{2} ) )\triangle t |\\ &\leq& \|\nabla (\xi_{n+1}+\widehat{\xi}_{n})\|_{0} \|\nabla \Pi_{n+1} (\frac{u_{n+1}-u_{n}}{\triangle t} -\frac{Q_{n+1}+Q_{n}}{2} ) \|_{0}\triangle t\\ &\leq& C\|\nabla(\xi_{n+1}+\widehat{\xi}_{n})\|_{0} \|\nabla \int^{t_{n+1}}_{t_{n}} (\frac{\partial u}{\partial t}-Q_{n+\frac{1}{2}} )d\tau \|_{0}\\ &\leq& C\, (\|\nabla \xi_{n+1}\|_{0}^{2}+\|\nabla \widehat{\xi}_{n}\|_{0}^{2} )\triangle t +C ( \int^{t_{n+1}}_{t_{n}} \|\frac{\partial^{2} Q}{\partial t^{2}} \|_{1}^{2}d\tau )(\triangle t)^{4}. \end{eqnarray} $

(4.13)

By $ (4.3) $ and Young inequality, we have

$ \begin{equation} \begin{array}{ll} |E_{n}(\lambda_{n+1}+\widehat{\lambda}_{n})| \leq& C ( \int^{t_{n+1}}_{t_{n}} ( \| \frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}Q}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}f}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{1} )d\tau ) (\triangle t)^{4}\\ &+ \frac{\alpha}{8}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t +Ch^{2}\triangle t+\frac{\alpha}{8}\|\lambda_{n+1}+\widehat{\lambda}_{n}\|_{0}^{2}\triangle t. \end{array} \end{equation} $

(4.14)

Substituting $ (4.10)-(4.14) $ into $ (4.9) $, we get

$ \begin{equation} \begin{array}{ll} (1-C\triangle t) (\|\lambda_{n+1}\|^{2}_{0}+\|\xi_{n+1}\|^{2}_{h} ) -(1+C\triangle t)(\|\widehat{\lambda}_{n}\|^{2}_{0}+\|\widehat{\xi}_{n}\|^{2}_{h}) \leq Ch^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d \tau +\phi_{n}, \end{array} \end{equation} $

(4.15)

where

$ \phi_{n} = C ( \int^{t_{n+1}}_{t_{n}} ( \| \frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}Q}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}f}{\partial t^{2}} \|^{2}_{0}+ \|\frac{\partial^{2}u}{\partial t^{2}} \|^{2}_{1} )d\tau ) (\triangle t)^{4} +Ch^{2}\triangle t+Ch^{4}\triangle t. $

On the other hand, by $ (3.8) $, we obtain

$ \begin{equation} (\widehat{\lambda}_{n}-\lambda_{n}, v_{h}) = (\widehat{\theta}_{n}-\theta_{n}, v_{h}), \; \forall \ v_{h}\in V_{n+1}^{h}. \end{equation} $

(4.16)

Choosing $ v_{h} = \widehat{\lambda}_{n} $ into $ (4.16) $, applying $ \varepsilon $-Cauchy inequality to yield

$ \begin{equation} \|\widehat{\lambda}_{n}\|^{2}_{0}\leq\frac{1}{1-\varepsilon}\|\lambda_{n}\|^{2}_{0} +\frac{1}{(1-\varepsilon)\varepsilon}\|\widehat{\theta}_{n}-\theta_{n}\|^{2}_{0}. \end{equation} $

(4.17)

Similarly, by $ (3.7) $, we obtain

$ \begin{equation} \|\widehat{\xi}_{n}\|^{2}_{0}\leq\frac{1}{1-\varepsilon}\|\xi_{n}\|^{2}_{0} +\frac{1}{(1-\varepsilon)\varepsilon}\|\widehat{\eta}_{n}-\eta_{n}\|^{2}_{0}. \end{equation} $

(4.18)

Substituting $ (4.17)-(4.18) $ into $ (4.15) $, we get

$ \begin{equation} \begin{array}{ll} & \frac{(1-\varepsilon)(1-C\triangle t)}{1+C\triangle t}(\|\lambda_{n+1}\|^{2}_{0}+\|\xi_{n+1}\|^{2}_{h}) -(\|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h})\\ \leq& Ch^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d\tau +\phi_{n}+\frac{1}{\varepsilon}\|\widehat{\theta}_{n}-\theta_{n}\|^{2}_{0}+\frac{1}{\varepsilon}\|\widehat{\eta}_{n}-\eta_{n}\|^{2}_{0}. \end{array} \end{equation} $

(4.19)

When there is no change between the two layers of the meshes,

$ \begin{equation} \begin{array}{ll} \frac{1-C\triangle t}{1+C\triangle t}\, (\|\lambda_{n+1}\|^{2}_{0}+\|\xi_{n+1}\|^{2}_{h}) -(\|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h} ) \leq Ch^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d\tau+\phi_{n}. \end{array} \end{equation} $

(4.20)

We use $ M $ to represent the number of mesh changes and assume that $ M $ is a bounded number independent of $ h $ and $ \triangle t $. Assume that the former $ M $-level mesh changes, and the $ N-M $ layer mesh does not change, we have

$ \begin{eqnarray} && \frac{(1-\varepsilon) (1-C\triangle t)}{1+C\triangle t}\, (\|\lambda_{n+1}\|^{2}_{0}+\|\xi_{n+1}\|^{2}_{h} ) - (\|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h} ) \end{eqnarray} $

(4.21)

$ \begin{eqnarray} &\leq& Ch^{2} \int^{t_{n+1}}_{t_{n}}\|\frac{\partial Q}{\partial t}\|_{1}^{2}d\tau+\phi_{n} +\frac{1}{\varepsilon}\|\widehat{\theta}_{n}-\theta_{n}\|^{2}_{0}+\frac{1}{\varepsilon}\|\widehat{\eta}_{n}-\eta_{n}\|^{2}_{0}, \qquad n = 1, 2, \cdots, M-1, \\ && \frac{1-C\triangle t}{1+C\triangle t}\, (\|\lambda_{n+1}\|^{2}_{0}+\|\xi_{n+1}\|^{2}_{h} ) - (\|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h} )\\&\leq& Ch^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d\tau +\phi_{n}, \qquad n = M, M+1, M+2, \cdots, N-1.\nonumber \end{eqnarray} $

(4.22)

Take $ n = 1, $ notice the selection of $ u_{1} $ and $ Q_{1} $, $ \|\lambda_{1}\|_{0} = \|\xi_{1}\|_{0} = 0 $, by (4.21), we obtain

$ \begin{equation} \begin{array}{ll} \|\lambda_{2}\|^{2}_{0}+\|\xi_{2}\|^{2}_{h} \leq \frac{1+C\triangle t}{(1-\varepsilon)(1-C\triangle t)} (Ch^{2} \int^{t_{n+1}}_{t_{n}}\|\frac{\partial Q}{\partial t}\|_{1}^{2}d\tau+\phi_{1} +\frac{1}{\varepsilon}\|\widehat{\theta}_{1}-\theta_{1}\|^{2}_{0}+\frac{1}{\varepsilon}\|\widehat{\eta}_{1}-\eta_{1}\|^{2}_{0} ), \end{array} \end{equation} $

(4.23)

when $ \triangle t $ is full small

$ \begin{equation} (\frac{1+C\triangle t}{1-C\triangle t} )^{n}\leq (\frac{1+C\triangle t}{1-C\triangle t} )^{N} = (1+\frac{2C\triangle t}{1-C\triangle t} )^{\frac{T}{\triangle t}}\leq e^{2CT}. \end{equation} $

(4.24)

Again by the $ \phi_{n} $, we obtain

$ \|\lambda_{2}\|^{2}_{0}+\|\xi_{2}\|^{2}_{h}\leq C_{1}(h^{2}+(\triangle t)^{4}). $

Furthermore, we have

$ \|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h}\leq C_{n-1}(h^{2}+(\triangle t)^{4}), \qquad 3\leq n\leq M+1, $

where $ C_{n-1} \ (n = 2, \cdots, M+1) $ are bounded number.

On both sides of (4.22) is multiplied by $ (\frac{1-C\triangle t}{1+C\triangle t})^{n} $, we have

$ \begin{equation} \begin{array}{ll} & (\frac{1-C\triangle t}{1+C\triangle t} )^{N}(\|\lambda_{N}\|^{2}_{0}+\|\xi_{N}\|^{2}_{h}) - (\frac{1-C\triangle t}{1+C\triangle t} )^{M+1}(\|\lambda_{M+1}\|^{2}_{0}+\|\xi_{M+1}\|^{2}_{h})\\ \leq &C\sum\limits_{n = M+1}^{N-1} (h^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d\tau +\phi_{n}) (\frac{1-C\triangle t}{1+C\triangle t} )^{n}, \end{array} \end{equation} $

(4.25)

that is

$ \begin{equation} \begin{array}{ll} \|\lambda_{N}\|^{2}_{0}+\|\xi_{N}\|^{2}_{h} \leq &Ce^{2CT} (\|\lambda_{M+1}\|^{2}_{0}+\|\xi_{M+1}\|^{2}_{h} + T\max\limits_{0\leq n\leq N} (\phi_{n}+h^{2} \int^{t_{n+1}}_{t_{n}} \|\frac{\partial Q}{\partial t} \|_{1}^{2}d\tau ) ), \end{array} \end{equation} $

(4.26)

take $ \triangle t $ full small, by Gronwall inequalities, we have

$ \|\lambda_{N}\|^{2}_{0}+\|\xi_{N}\|^{2}_{h} \leq C(h^{2}+(\triangle t)^{4}). $

Furthermore

$ \max\limits_{M+1\leq n\leq N} (\|\lambda_{n}\|^{2}_{0}+\|\xi_{n}\|^{2}_{h} ) \leq C(h^{2}+(\triangle t)^{4}). $

Again

$ \|Q^{h}_{n}-Q_{n}\|_{0}^{2}+\|u^{h}_{n}-u_{n}\|^{2}_{h} \leq \| \lambda_{n}\|^{2}_{0}+\| \theta_{n}\|^{2}_{0}+\| \xi_{n}\|^{2}_{h}+\| \eta_{n}\|^{2}_{h}, $

by Lemma 2.1 and the finite element interpolation theory, we obtain (4.6), the proof is completed.