Convergence of a high-order compact finite difference scheme for the Klein-Gordon-Schrödinger equations

Accepted Manuscript Convergence of a high-order compact finite difference scheme for the Klein-Gordon-Schrödinger equations

Juan Chen, Fangqi Chen

PII: DOI: Reference:

S0168-9274(19)30057-1 https://doi.org/10.1016/j.apnum.2019.03.004 APNUM 3516

To appear in:

Applied Numerical Mathematics

Received date: Revised date: Accepted date:

23 June 2018 17 February 2019 15 March 2019

Please cite this article in press as: J. Chen, F. Chen, Convergence of a high-order compact finite difference scheme for the Klein-Gordon-Schrödinger equations, Appl. Numer. Math. (2019), https://doi.org/10.1016/j.apnum.2019.03.004

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Convergence of a high-order compact finite difference scheme for the Klein-Gordon-Schr¨odinger equations Juan Chena,b,∗, Fangqi Chena a

Department of Mathematics, Nanjing university of Aeronautics and Astronautics, Nanjing 210016, PR China. b Department of Mathematics and Statistics, Changshu Institute of Technology, Changshu 215500, PR China.

Abstract In this paper, a conservative and linearly implicit finite difference scheme with second order temporal accuracy and eighth order spatial accuracy by means of the Richardson extrapolation method is proposed for approximating solution to the initial-boundary problem for the Klein-Gordon-Schr¨odinger equations. Furthermore, the difference solution is shown to be bounded by taking advantage of the fact that the proposed scheme possesses two conservation laws, and with the aid of the discrete energy method, the difference scheme is demonstrated to be convergent in the maximum norm. Finally, numerical experiments are assigned to support the theoretical results. Keywords: Klein-Gordon-Schr¨odinger equations, Compact finite difference scheme, Conservation, Convergence in the maximum norm 2010 MSC: 65M06, 65M12 1. Introduction The Klein-Gordon-Schr¨odinger (KGS) equations: √ 1 i = −1, iφt + φxx + uφ = 0, x ∈ R, t  0, 2 utt − uxx + u − |φ|2 = 0, x ∈ R, t  0,

(1.1) (1.2)

∗

Corresponding author at: Department of Mathematics and Statistics, Changshu Institute of Technology, Changshu 215500, PR China. Email address: [email protected] (Juan Chen) Preprint submitted to Applied Numerical Mathematics

April 5, 2019

are classical models which describe the interaction between conservative complex neutron field and neutral meson Yukawa in quantum field theory [1, 2], where φ(x, t) and u(x, t) represent a complex scalar nucleon field and a real scalar meson field, respectively. Many numerical studies have been carried out for the KGS equations, such as spectral method [3, 4], Chebyshev pseudospectral multidomain method [5], multi-symplectic method [6, 7] and finite difference method and so on. Obviously, each numerical method has its merits and demerits, in particular, the theoretical precision of the difference schemes (see [8, 9]) is O(τ 2 + h2 ). In [10], Pan et al. presented a three-level linear difference scheme which is second-order accuracy in time and fourthorder accuracy in space. Furthermore, Wang [11] developed a fourth-order implicit compact difference scheme, and proved that the scheme is stable and convergent without any restrictions. Following the discussion of Xia et al. [12], the analytic solution of the KGS equations (1.1)–(1.2) decays rapidly to zero as |x|  0 for a fixed t. Therefore, numercially we can solve (1.1)–(1.2) in a bounded domain Ω = [a, b], where −a  0 and b  0. Consider the numerical solution of (1.1)-(1.2) with the initial conditions ∂u φ(x, 0) = φ0 (x), u(x, 0) = u0 (x), (x, 0) = u1 (x), x ∈ [a, b], (1.3) ∂t and the boundary conditions φ(a, t) = φ(b, t) = 0,

u(a, t) = u(b, t) = 0,

t ∈ (0, T ],

(1.4)

where φ0 (x), u0 (x) and u1 (x) are given smooth functions. As is well known, Problems (1.1)–(1.4) possess the following mass and energy conserved quantities([8]):  b Q(t) := |φ(x, t)|2 dx = const, (1.5) a  1 b E(t) := (|φx (x, t)|2 + |u(x, t)|2 + |ux (x, t)|2 + |ut (x, t)|2 )dx − 2 a  b |φ(x, t)|2 u(x, t)dx = const. (1.6) a

In the past few years, there has been growing interest in developing, analysing and implementing high-order compact finite difference methods 2

for solving partial differential equations [13–19]. However, a prior estimate in L∞ norm is often hard to be obtained because of the difficulty caused by compact operator. Therefore, the stability and convergence in L∞ norm of any compact difference scheme are difficult to demonstrate. Recently, by the Richardson extrapolation method, Wang [20] has presented a conservative high-order compact difference scheme for the nonlinear Schr¨odinger equation with periodic boundary conditions. He also proved that the numerical solution of the proposed scheme converges to the analytic solution in L∞ norm by means of a technique named ”regression of compactness”. Although the technique is invalid for equations with aperiodic conditions, [20] provides us with a new way of thinking on developing higher-precision schemes for computing partial differential equations. This article aims to construct a new difference scheme which has the following three advantages: (i) The proposed scheme is compact, linearized and decoupled. (ii) The scheme simulates the two conservative laws (1.5) and (1.6) well. (ii) The scheme is convergent for φ in the discrete L∞ , and respectively, for u in the discrete L∞ , at the order of O(τ 2 + h8 ) with mesh size h and time step τ . The rest of this paper is organized as follows: In Section 2, a linearly implicit and conservative difference scheme is proposed by the Richardson extrapolation method, and in Section 3, the difference solution is demonstrated to be bounded by means of the conservations of difference scheme. What is more, Section 4 is devoted to the unconditional convergence and stability of the proposed scheme. Finally, in Section 5, we verify our theoretical analysis by numerical experiments. 2. Finite Difference Scheme and Its Conservative Laws Before proposing the new difference scheme, we will introduce some notations that will be used in this article. Let (Φnj , Ujn ) and (φnj , unj ) be the numerical approximation and the exact solution of (φ(x, t), u(x, t)) at the point (xj , tn ), respectively. Denote J = [ b−a ], Ωh = {xj |xj = a + jh, j = −2, −1, 0, 1, · · · , J, J + 1, J + 2}, and h Ωτ = {tn |tn = nτ, n = 0, 1, 2, · · · , N = [ Tτ ]}, where h and τ are the spatial and temporal step size, respectively. In addition, let Vh0 = {v = (vj )|v−2 = v−1 = v0 = vJ = vJ+1 = vJ+2 = 0, j = −2, −1, 0, 1, · · · , J, J + 1, J + 2} be the space of complex grid functions on Ωh . We define the following difference 3

quotients, inner product and norms: n n − wjn wj+1 wjn − wj−1 , (wjn )x = , (wjn )x = h h wjn+1 − wjn wjn − wjn−1 , (wjn )t = , (wjn )t = τ τ   J−1 J−1    p (u, v) = h uj v j , vp =  h |vj |p , j=1

n n − wj−1 wj+1 , 2h wjn+1 − wjn−1 (wjn )tˆ = , 2τ

(wjn )xˆ =

v∞ = max |vj |, 1≤j≤J−1

j=1

where u, v ∈ Vh0 . For convenience, we write ·2 as ·, and let C be a general positive constant which may have different values on different occasions. To develop a high-order difference scheme for solving Problems (1.1)– (1.4), we first discretize the second-order derivative φxx in Eq.(1.1) by β(φxx )j−2 + α(φxx )j−1 + (φxx )j + α(φxx )j+1 + β(φxx )j+2 = γ(φj )xˆxˆ + ξ(φj )x¯x + O(h8 ), where α = Denote

344 ,β 1179

=

23 ,γ 2358

=

930 ,ξ 1179

=

960 1179

(2.1)

([17]).

Ah φj = βφj−2 + αφj−1 + φj + αφj+1 + βφj+2 .

(2.2)

Then, (2.1) can be rewritten as Ah (φxx )j = γ(φj )xˆxˆ + ξ(φj )x¯x + O(h8 ).

(2.3)

In the same way, for the derivative uxx in Eq.(1.2), we have Ah (uxx )j = γ(uj )xˆxˆ + ξ(uj )x¯x + O(h8 ).

(2.4)

Thus, applying (2.3), (2.4) to (1.1) and (1.2), and using the C-N method, we obtain a new compact finite difference scheme for KGS equations (1.1)– (1.4) as follows: 1 iAh (Φnj )tˆ + [γ(Φn+1 + Φn−1 )xˆxˆ + ξ(Φn+1 + Φn−1 )x¯x ] j j j j 4 1 + Ah Ujn (Φn+1 + Φn−1 ) = 0, 1 ≤ j ≤ J − 1, 1 ≤ n ≤ T − 1, (2.5) j j 2 1 1 Ah (Ujn )tt¯ − [γ(Ujn+1 + Ujn−1 )xˆxˆ + ξ(Ujn+1 + Ujn−1 )x¯x ] + Ah (Ujn+1 2 2 n−1 n 2 +Uj ) − Ah |Φj | = 0, 1 ≤ j ≤ J − 1, 1 ≤ n ≤ T − 1, (2.6) Φ0j = φ0 (xj ), Uj0 = u0 (xj ), (Uj0 )tˆ = u1 (xj ), Φn0 = ΦnJ = 0, U0n = UJn = 0, 0 ≤ n ≤ T. 4

1 ≤ j ≤ J − 1, (2.7) (2.8)

Obviously, the scheme (2.5)–(2.8) is suitable for parallel computing and can be easily solved without iterative technique. Noticing the scheme (2.5) cannot be started by itself, we need to propose another suitable two-level eighth-order scheme to calculate Φ1 , such as: 1 iAh (Φ0j )t + [γ(Φ1j + Φ0j )xˆxˆ + ξ(Φ1j + Φ0j )x¯x ] 4 1 1 + Ah (Uj + Uj0 )(Φ1j + Φ0j ) = 0, 1 ≤ j ≤ J − 1, 4 By making use of the following definitions

(2.9)

Φn = (Φn1 , Φn2 , · · ·, ΦnJ−1 )T , u1 = (u1 (x1 ), u1 (x2 ), · · · , u1 (xJ−1 ))T , n )T , |Φn |2 = (|Φn1 |2 , |Φn2 |2 , · · ·, |ΦnJ−1 |2 )T , U n = (U1n , U2n , · · ·, UJ−1 n−1 T n + Φn−1 ), · · ·, UJ−1 (Φn+1 U n (Φn+1 + Φn−1 ) = (U1n (Φn+1 1 1 J−1 + ΦJ−1 )) , we rewrite the derived scheme (2.5)–(2.9) as 1 iH(Φn )tˆ + [γ(Φn+1 + Φn−1 )xˆxˆ + ξ(Φn+1 + Φn−1 )x¯x ] 4 1 + HU n (Φn+1 + Φn−1 ) = 0, 1 ≤ n ≤ T − 1, 2 1 H(U n )tt¯ − [γ(U n+1 + U n−1 )xˆxˆ + ξ(U n+1 + U n−1 )x¯x ] 2 1 1 ≤ n ≤ T − 1, + H(U n+1 + U n−1 ) − H|Φn |2 = 0, 2 1 iH(Φ0 )t + [γ(Φ1 + Φ0j )xˆxˆ + ξ(Φ1 + Φ0 )x¯x ] 4 1 + H(U 1 + U 0 )(Φ1 + Φ0 ) = 0, 4 0 Φj = φ0 (xj ), Uj0 = u0 (xj ), Utˆ0 = u1 , 1 ≤ j ≤ J − 1, U0n = UJn = 0, 0 ≤ n ≤ T, Φn0 = ΦnJ = 0, where

⎡

1 α β

⎢ ⎢ ⎢ ⎢ H=⎢ ⎢ ⎢ ⎣ 0 0

α β 0 0 0 ··· 1 α β 0 0 ··· α 1 α β 0 ··· .. .. .. .. .. .. . . . . . . 0 ··· ··· ··· ··· ··· 0 ··· ··· ··· ··· ··· 5

0 0 0 0 0 0

0 0 0

⎤

⎥ ⎥ ⎥ ⎥ . ⎥ ⎥ ⎥ α 1 α ⎦ β α 1 (J−1)×(J−1)

(2.10)

(2.11)

(2.12) (2.13) (2.14)

It is obvious that H is a real-value, symmetric and positive definite pentadiagonal matrix, therefore H is invertible and its inverse matrix M = H −1 also is a symmetric and positive definite matrix. Then, for the convenience of theoretical analysis, the scheme (2.10)–(2.14) is further rewritten as 1 i(Φn )tˆ + M [γ(Φn+1 + Φn−1 )xˆxˆ + ξ(Φn+1 + Φn−1 )x¯x ] 4 1 n n+1 + Φn−1 ) = 0, 1 ≤ n ≤ T − 1, + U (Φ 2 1 (U n )tt¯ − M [γ(U n+1 + U n−1 )xˆxˆ + ξ(U n+1 + U n−1 )x¯x ] 2 1 n+1 + (U + U n−1 ) − |Φn |2 = 0, 1 ≤ n ≤ T − 1, 2 1 i(Φ0 )t + M [γ(Φ1 + Φ0j )xˆxˆ + ξ(Φ1 + Φ0 )x¯x ] 4 1 1 + (U + U 0 )(Φ1 + Φ0 ) = 0, 4 Φ0j = φ0 (xj ), Uj0 = u0 (xj ), Utˆ0 = u1 , 1 ≤ j ≤ J − 1, U0n = UJn = 0, 0 ≤ n ≤ T. Φn0 = ΦnJ = 0,

(2.15)

(2.16)

(2.17) (2.18) (2.19)

Depending on the lemmas below, we will prove that the difference solution of the scheme (2.15)–(2.19) conforms to two conservation laws. Lemma 1. [21] For any two mesh functions u, v ∈ Vh0 , there are h

h

J−1 

(uj )x¯x v j = −h

J−1 

j=1

j=0

J−1 

J 

(uj )xˆxˆ v j = −h

j=1

(uj )x (v j )x ,

(2.20)

(uj )xˆ (v j )xˆ .

(2.21)

j=0

Lemma 2. [22] For any a real value symmetric positive definite matrix M , we have Im(M unx¯x , un ) = 0, Re(M (un+1 + un−1 )x¯x , un+1 − un−1 ) 2 n−1 2 = −(Run+1 x  − Rux  ). where R is obtained from M by the Cholesky decomposition. 6

(2.22) (2.23)

We can prove the following lemma in a similar way of Lemma 2. Lemma 3. For any a real value symmetric positive definite matrix M , there are Im(M unxˆxˆ , un ) = 0, Re(M (un+1 + un−1 )xˆxˆ , un+1 − un−1 ) n−1 2 2 = −(Run+1 x ˆ  − Rux ˆ  ).

(2.24) (2.25)

where R is obtained from M by the Cholesky decomposition. Theorem 1. Suppose that φ0 (x) ∈ H 1 ∩ L4 , u0 (x) ∈ H 1 , u1 (x) ∈ L2 , ϕ(x, t) ∈ C 10,4 and u(x, t) ∈ C 10,4 , then the scheme (2.15)–(2.19) is conservative for the discrete mass and the discrete energy, that is Qn = Qn−1 = · · · = Q0 , E n = E n−1 = · · · = E 0 ,

(2.26) (2.27)

where 1 Qn = (Φn+1 2 + Φn 2 ), 2

(2.28)

and 1 γ ξ E n = Utn 2 + (RUxˆn+1 2 + RUxˆn 2 ) + (RUxn+1 2 + RUxn 2 ) 2 4 4 1 γ ξ 2 n 2 n+1 2 + (U n+1 2 + U n 2 ) + (RΦn+1 ˆ  ) + (RΦx  x ˆ  + RΦx 4 4 4 J−1  1 +RΦnx 2 ) − h (U n |Φn+1 |2 + Ujn+1 |Φnj |2 ). (2.29) 2 j=1 j j proof. By computing an inner product of (2.15) with Φn+1 + Φn−1 and taking the imaginary part, we have 1 (Φn+1 2 − Φn−1 2 ) = 0, 2τ

(2.30)

and from this we deduce that 1 Qn = (Φn+1 2 + Φn 2 ) = Qn−1 = · · · = Q0 . 2 7

(2.31)

Computing an inner product of (2.15) with Φn+1 − Φn−1 , then taking the real part yields ξ γ n−1 2 2 n+1 2 n−1 2 − (RΦn+1 x ˆ  − RΦx ˆ  ) − (RΦx  − RΦx  ) 4 4 J−1  1 U n (|Φn+1 |2 − |Φn−1 |2 ) = 0. (2.32) + h j j 2 j=1 j Next, taking an inner product of (2.16) with U n+1 − U n−1 implies γ Utn 2 − Utn−1 2 + (RUxˆn+1 2 − RUxˆn−1 2 ) 2 ξ 1 + (RUxn+1 2 − RUxn−1 2 ) + (U n+1 2 − U n−1 2 ) 2 2 J−1  |Φnj |2 (Ujn+1 − Ujn−1 ) = 0. −h

(2.33)

j=1

Then, it follows from (2.32) and (2.33) that E n = E n−1 = · · · = E 0 .

(2.34)

This completes the proof. 3. Prior Estimate Lemma 4. [23] For any a real value symmetric positive definite matrix M , there exist two positive constants C0 and C1 , such that C0 un 2 ≤ (M un , un ) = Run 2 ≤ C1 un 2 ,

(3.1)

where R is obtained from M by the Cholesky decomposition. Lemma 5. (Discrete Sobolev’s inequality [24]) For any discrete function uh = {uj |j = 0, 1, · · · , J} on the finite interval [0, l] and for any given ε > 0, there exists a constant C, depending only on ε, such that δ k uh p ≤ εδ n uh 2 + Cuh 2 , where 2 ≤ p ≤ ∞, 0 ≤ k < n. 8

(3.2)

In order to study the unconditional convergence and stability of the discrete scheme (2.15)–(2.19), we need the following prior estimate of the difference solution. Theorem 2. Under the conditions listed in Theorem 1, we have Φn ∞ ≤ C,

U n ∞ ≤ C.

(3.3)

proof. From Theorem 1, we obtain Φn  ≤ C,

(3.4)

and γ ξ 2 n 2 n+1 2 Utn 2 + (RUn+1  + RUxn 2 ) ˆ  ) + (RUx x ˆ  + RUx 2 2 1 γ ξ 2 n 2 n+1 2 + (U n+1 2 + U n 2 ) + (RΦn+1 ˆ  ) + (RΦx  x ˆ  + RΦx 2 2 2 J−1  n 2 +RΦx  ) ≤ C + h (Ujn |Φn+1 |2 + Ujn+1 |Φnj |2 ). (3.5) j j=1

An application of Young’s inequality yields h

J−1 

(Ujn |Φn+1 |2 + Ujn+1 |Φnj |2 ) j

j=1

1 ≤ (U n 2 + U n+1 2 ) + (Φn+1 44 + Φn 44 ). 4 Furthermore, utilizing Lemma 5 and the estimate (3.4), we have Φn+1 44 ≤ Φn+1 2∞ Φn+1 2 ≤ CΦn+1 2∞ 2 n+1 2 n+1 2 ≤ εΦn+1 x  + C(ε)Φx  ≤ εΦx  + C,

(3.6)

(3.7)

and Φn 44 ≤ Φn 2∞ Φn 2 ≤ εΦnx 2 + C,

(3.8)

where ε is a positive number which can be arbitrarily small. Then, substituting (3.6), (3.7) and (3.8) into (3.5) and using Lemma 4, we have U n  ≤ C,

Φnx  ≤ C,

Uxn  ≤ C,

Φnxˆ  ≤ C,

Uxˆn  ≤ C.

(3.9)

Therefore, according to Lemma 5, we get Φn ∞ ≤ C,

U n ∞ ≤ C.

This completes the proof. 9

(3.10)

4. Convergence and Stability We are going to prove the convergence and stability of the difference scheme (2.15)–(2.19) in this section. First of all, we introduce some auxiliary lemmas. Lemma 6. (Gronwall’s inequality [24]) Suppose that the nonnegative mesh functions {w(n), ρ(n), n = 1, 2, · · · , N, N τ = T } satisfy the inequality: w(n) ≤ ρ(n) + τ

n 

Bl w(l),

l=1

where Bl (l = 1, 2, · · · , N ) are nonnegative constants. Then for any 0 ≤ n ≤ N , there is w(n) ≤ ρ(n) exp(nτ

n 

Bl ).

(4.1)

l=1

Lemma 7. (Gronwall’s inequality [21]) Assume {Gn |n ≥ 0} is a nonnegative sequence, and satisfies n+1

G

≤ φ + cτ

n 

Gl

l=1

where c and φ are nonnegative constants. Then G satisfies Gn ≤ φ exp(cnτ )

(4.2)

Lemma 8. [23] For time sequences w = {w0 , w1 , · · · , wn+1 } and g = {g 1 , g 2 , ·· ·, g n }, there is |2τ

n  l=1

g

l

wtˆl |

0 2

1 2

≤ |w | + |w | + τ

n−1 

|wl |2 + |wn |2 + |wn+1 |2

l=2

+|g 1 |2 + |g 2 |2 + τ

n−1 

|gtl |2 + |g n−1 |2 + |g n |2

(4.3)

l=2

In addition, we consider the truncation errors of the proposed scheme (2.15)–(2.19). 10

Define 1 r n = i(φn )tˆ + M [γ(φn+1 + φn−1 )xˆxˆ + ξ(φn+1 + φn−1 )x¯x ] 4 1 n n+1 + (u (φ + φn−1 )), 1 ≤ n ≤ T − 1, (4.4) 2 1 σ n = (un )tt¯ − M [γ(un+1 + un−1 )xˆxˆ + ξ(un+1 + un−1 )x¯x ] 2 1 n+1 + un−1 ) − |φn |2 , 1 ≤ n ≤ T − 1, (4.5) + (u 2 1 r 0 = i(φ0 )t + M [γ(φ1 + φ0 )xˆxˆ + ξ(φ1 + φ0 )x¯x ] 4 1 1 (4.6) + (u + u0 )(φ1 + φ0 ), 4 φ0j = φ0 (xj ), u0j = u0 (xj ), u0tˆ = u1 , 1 ≤ j ≤ J − 1, (4.7) n n n n φ0 = φJ = 0, u0 = uJ = 0, 0 ≤ n ≤ T, (4.8) n n where r n = (r1n , r2n , · · · , rJ−1 )T and σ n = (σ1n , σ2n , · · · , σJ−1 )T . Making use of Taylor expansion theorem, we obtain the following lemma.

Lemma 9. Under the conditions of Theorem 1, the trucation errors of the difference scheme (2.15)–(2.19) satisfy: |rjn | + |σjn | = O(τ 2 + h8 ) and |(rjn )t | + |(σjn )t | = O(τ 2 +h8 ) as τ → 0, h → 0, where 1 ≤ j ≤ J −1 and 0 ≤ n ≤ T −1. Finally, we analyze the convergence and stability of the scheme (2.15)– (2.19). Theorem 3. Suppose the conditions listed in Theorem 1 hold, then the difference solution of the scheme (2.15)–(2.19) is convergent to the solution of Problems (1.1)–(1.4) with the order O(τ 2 + h8 ) in L∞ norm. proof. Let en = φn − Φn ,

η n = un − U n .

(4.9)

Subtracting (2.15)–(2.19) from (4.4)–(4.8), respectively, we have 1 1 r n = i(en )tˆ + M [γ(en+1 + en−1 )xˆxˆ + ξ(en+1 + en−1 )x¯x ] + Gn1 , 4 2 1 ≤ n ≤ T − 1, (4.10) 11

1 σ n = (η n )tt¯ − M [γ(η n+1 + η n−1 )xˆxˆ + ξ(η n+1 + η n−1 )x¯x ] 2 1 n+1 + (η + η n−1 ) − Gn2 , 1 ≤ n ≤ T − 1, 2 1 1 r 0 = i(e0 )t + M [γ(e1 + e0 )xˆxˆ + ξ(e1 + e0 )x¯x ] + G01 , 4 4 e0j = 0, ηj0 = 0, ηtˆ0 = 0, 1 ≤ j ≤ j − 1,

en0 = e0J = 0,

η0n = ηJn = 0,

0 ≤ n ≤ T,

(4.11) (4.12) (4.13) (4.14)

where Gn1 = un (φn+1 + φn−1 ) − U n (Φn+1 + Φn−1 ) = η n (φn+1 + φn−1 ) + U n (en+1 + en−1 ), ¯ n, Gn2 = |φn |2 − |Φn |2 = φn e¯n + en Φ G01 = (u1 + u0 )(φ1 + φ0 ) − (U 1 + U 0 )(Φ1 + Φ0 ) = η 1 (φ1 + φ0 ) + e1 (U 1 + U 0 ). get

(4.15) (4.16) (4.17)

Noticing Gn10 = Gn1J = Gn20 = Gn2J = 0 and according to Theorem 2, we Gn1 2 ≤ C(η n 2 + en+1 2 + en−1 2 ), Gn2 2 ≤ Cen 2 ,

(4.18) (4.19)

Similarly, we obtain (Gn1 )x 2 ≤ C(η n 2 + en+1 2 + en−1 2 + ηxn 2 2 n−1 2 +en+1 (4.20) x  + ex  ), n+1 2 n 2 n 2 n+1 2 n−1 2 n 2 (G1 )xˆ  ≤ C(η  + e  + e  + ηxˆ  + exˆ  2 n+1 2 n−1 2 +en−1 (4.21) x ˆ  + ex  + ex  ). Computing an inner product of (4.10) with en+1 + en−1 , taking the imaginary part and using Lemma 2 and Lemma 3, one can get 1 (en+1 2 − en−1 2 ) + Im(Gn1 , en+1 + en−1 ) 2τ (4.22) = Im(r n , en+1 + en−1 ). From Cauchy-Schwarz inequality and (4.18), we obtain 1 |Im(Gn1 , en+1 + en−1 )| ≤ (Gn1 2 + en+1 + en−1 2 ) 2 n 2 ≤ C(η  + en+1 2 + en−1 2 ), 12

(4.23)

and |Im(r n , en+1 + en−1 )| ≤ C(r n 2 + en+1 2 + en−1 2 ).

(4.24)

Therefore, it is obvious that 1 (en+1 2 − en−1 2 ) ≤ Cτ (r n 2 + η n 2 + en+1 2 + en−1 2 ). (4.25) 2 Taking an inner product of (4.11) with η n+1 − η n−1 , then using Lemma 2 and Lemma 3, we obtain γ ξ 2 − Rηxn−1 2 ) + (Rηxn+1 2 ηtn 2 − ηtn−1 2 + (Rηxn+1 ˆ ˆ 2 2 1 −Rηxn−1 2 ) + (η n+1 2 − η n−1 2 ) − (Gn2 , η n+1 − η n−1 ) 2 n n+1 − η n−1 ). (4.26) = (σ , η It follows from (4.19) and Cauchy-Schwarz inequality that |(Gn2 , η n+1 − η n−1 )| ≤ Cτ (en 2 + ηtn 2 + ηtn−1 2 ).

(4.27)

Similarly, we have |(σ n , η n+1 − η n−1 )| ≤ Cτ (σ n 2 + ηtn 2 + ηtn−1 2 ).

(4.28)

Substituting (4.27), (4.28) into (4.26) yields γ 2 − Rηxn−1 2 ) ηtn 2 − ηtn−1 2 + (Rηxn+1 ˆ ˆ 2 1 ξ + (Rηxn+1 2 − Rηxn−1 2 ) + (η n+1 2 − η n−1 2 ) 2 2 n 2 n 2 n 2 ≤ Cτ (σ  + e  + ηt  + ηtn−1 2 ).

(4.29)

Then, adding (4.29) to (4.25), one can get 1 γ 2 (en+1 2 − en−1 2 ) + ηtn 2 − ηtn−1 2 + (Rηxn+1 ˆ 2 2 ξ 1 −Rηxn−1 2 ) + (Rηxn+1 2 − Rηxn−1 2 ) + (η n+1 2 − η n−1 2 ) ˆ 2 2 n 2 n 2 n 2 n 2 n 2 ≤ Cτ (r  + σ  + η  + e  + ηt  + ηtn−1 2 ). (4.30) 13

Denote 1 γ B n = (en+1 2 + en 2 ) + ηtn 2 + (Rηxn+1 2 + Rηxnˆ 2 ) ˆ 2 2 ξ 1 + (Rηxn+1 2 + Rηxn 2 ) + (η n+1 2 + η n 2 ). (4.31) 2 2 According to Lemma 9, we can rewrite (4.30) as B n − B n−1 ≤ Cτ ((r n 2 + σ n 2 ) + Cτ (B n + B n−1 ) ≤ Cτ ((τ 2 + h8 )2 ) + Cτ (B n + B n−1 ).

(4.32)

Thus, by Lemma 6, it can be immediately obtained that B n ≤ C(τ 2 + h8 )2 .

(4.33)

This together with Lemma 4 gives en  ≤ C(τ 2 + h8 ), η n  ≤ C(τ 2 + h8 ),

ηxn  ≤ C(τ 2 + h8 ), ηxnˆ  ≤ C(τ 2 + h8 ).

(4.34)

Hence, we have η n ∞ ≤ C(τ 2 + h8 ).

(4.35)

Next, computing an inner product of (4.10) with en+1 − en−1 and taking the real part, we obtain ξ γ n−1 2 2 n+1 2 n−1 2 − (Ren+1 x ˆ  − Rex ˆ  ) − (Rex  − Rex  ) 4 4 1 (4.36) + Re(Gn1 , en+1 − en−1 ) = Re(r n , en+1 − en−1 ). 2 From (4.10), we find 1 en+1 − en−1 = τ i{ M [γ(en+1 + en−1 )xˆxˆ + ξ(en+1 + en−1 )x¯x ] 2 n +G1 − 2r n }, (4.37) then, for the third term on the left side of (4.36), we have Re(Gn1 , en+1 − en−1 ) 14

1 = Re(Gn1 , τ i{ M [γ(en+1 + en−1 )xˆxˆ + ξ(en+1 + en−1 )x¯x ] + Gn1 − 2r n }) 2 τγ τξ = − Im(Gn1 , M (en+1 + en−1 )xˆxˆ ) − Im(Gn1 , M (en+1 + en−1 )x¯x ) 2 2 +2τ Im(Gn1 , r n ) τγ τξ = Im(R(Gn1 )xˆ , R(en+1 + en−1 )xˆ ) + Im(R(Gn1 )x , R(en+1 + en−1 )x ) 2 2 (4.38) +2τ Im(Gn1 , r n ). According to (4.20), (4.21), (4.34) and Lemma 4, one can get τξ Im(R(Gn1 )x , R(en+1 + en−1 )x )| 2 τξ ≤ (R(Gn1 )x 2 + R(en+1 + en−1 )x 2 ) 4 2 n−1 2 ≤ Cτ (η n 2 + en+1 2 + en−1 2 + ηxn 2 + en+1 x  + ex  ) 2 n−1 2 2 8 2 (4.39) ≤ Cτ (en+1 x  + ex  ) + Cτ (τ + h ) , τγ | Im(R(Gn1 )xˆ , R(en+1 + en−1 )xˆ )| 2 τγ (R(Gn1 )xˆ 2 + R(en+1 + en−1 )xˆ 2 ) ≤ 4 n−1 2 2 ≤ Cτ (η n 2 + en+1 2 + en−1 2 + ηxnˆ 2 + en+1 x ˆ  + ex ˆ  n+1 2 n−1 2 +ex  + ex  ) n−1 2 2 n+1 2 n−1 2 ≤ Cτ (en+1 x ˆ  + ex ˆ  + ex  + ex  ) (4.40) +Cτ (τ 2 + h8 )2 , n n n 2 n+1 2 n−1 2 n 2 |2τ Im(G1 , r )| ≤ Cτ (η  + e  + e  + r  ) (4.41) ≤ Cτ (τ 2 + h8 )2 .

|

Therefore, substituting (4.39), (4.40) and (4.41) into (4.38) implies n−1 2 2 |Re(Gn1 , en+1 − en−1 )| ≤ Cτ (en+1 x ˆ  + ex ˆ  n+1 2 n−1 2 2 +ex  + ex  ) + Cτ (τ + h8 )2 .

(4.42)

This together with (4.36) gives γ ξ n−1 2 2 n+1 2 n−1 2 (Ren+1 x ˆ  − Rex ˆ  ) + (Rex  − Rex  ) 4 4 n−1 2 2 n+1 2 n−1 2 ≤ Cτ (en+1  + e  + e x  + ex  ) x ˆ x ˆ +Cτ (τ 2 + h8 )2 − 2τ Re(r n , entˆ ). 15

(4.43)

Summing up (4.43) for the superscript n for 1 to m and then replacing m by n yields γ ξ 2 n 2 n+1 2 n 2 (Ren+1 ˆ  ) + (Rex  + Rex  ) x ˆ  + Rex 4 4 n  2 l 2 l+1 2 l 2 ≤ Cτ (el+1 ˆ  + ex  + ex  ) x ˆ  + ex l=0 2

8 2

+CT (τ + h ) − 2τ Re

n  l=1

(r l , eltˆ).

(4.44)

From Lemma 8 and Lemma 9, we have |2τ Re

n  l=1

(r l , eltˆ)| ≤ C(τ 2 + h8 )2 .

(4.45)

Consequently, substituting (4.45) into (4.44) and by making use of Lemma 4 and Lemma 7, we obtain enxˆ  ≤ C(τ 2 + h8 ),

enx  ≤ C(τ 2 + h8 ).

(4.46)

This together with (4.34) gives en ∞ ≤ C(τ 2 + h8 ).

(4.47)

This completes the proof. According to Theorem 3, we have Theorem 4. Assume the conditions in Theorem 1 hold, then the solution of the difference scheme (2.15)–(2.19) is stable in the sense of norm  · ∞ . 5. Numerical experiments In this section, we will present numerical results to demonstrate the efficiency of our difference scheme. In our test, the initial conditions √ 3 2 1 φ0 (x) = √ sech2 ( √ (x − x0 )) exp(ivx), (5.1) 2 4 1−v 2 1 − v2 3 1 u0 (x) = √ sech2 ( √ (x − x0 )), (5.2) 4 1 − v2 2 1 − v2 1 sech( 2√1−v 3v 2 (x − x0 )) , (5.3) u1 (x) = 3 1 4 (1 − v 2 )3 cosh ( 2√1−v2 (x − x0 )) 16

were used, where x0 is the initial phase and v is the propagating velocity of the solitary wave. Following the discussion of Xia et al. [12], the analytic solution of (1.1)–(1.2) is √ 3 2 1 φ(x, t) = √ sech2 ( √ (x − vt − x0 )) · 2 4 1−v 2 1 − v2 1 − v2 + v4 t)), (5.4) exp(i(vx + 2(1 − v 2 ) 3 1 u(x, t) = √ sech2 ( √ (x − vt − x0 )). (5.5) 2 4 1−v 2 1 − v2 For a fixed n, the numerical solution of the KGS equations(1.1)-(1.4) was calculated by the scheme (2.10)–(2.14) in practical computation. It can be seen that the scheme (2.10)–(2.14) results in two pentadiagonal linear systems, which can be easily solved by using the Thomas method. In order to process numerical investigations, we choose v = 0.8, x0 = 0 and set the spatial domain [a, b] = [−10, 30]. In addition, errors between analytic solution and numerical solution in the sense of L∞ norm are defined by m1 (h, τ ) = φn − Φn ∞ ,

m2 (h, τ ) = un − U n ∞ ,

(5.6)

and convergence orders of the proposed scheme can be calculated as m1 (h, τ ) ), m1 (h/2, τ /24 ) m2 (h, τ ) rate3 = log2 ( ), m2 (h/2, τ /24 )

rate1 = log2 (

m1 (h, τ ) ), (5.7) m1 (h, τ /2) m2 (h, τ ) rate4 = log2 ( ). (5.8) m2 (h, τ /2) rate2 = log2 (

Table 1 and 2 show L∞ norm errors and convergence orders for the derived scheme at T = 20. It is clearly seen that the numerical results are consistent with the theoretical ones that the numerical solution Φ, U are of order O(τ 2 +h8 ). Table 3 lists some values of discrete mass Qn and energy E n for Scheme (2.10)–(2.14). It is easily seen that this scheme possesses satisfactory conservative properties. Finally, Figure 1 and 2 plot the modulus of the analytic solution and numerical solution computed by Scheme (2.10)–(2.14) with h = 0.1, τ = 0.05, from which we can also conclude that the proposed scheme is accurate and suitable for long-term computation. 17

Table1.Errors and convergence rates of Φn and U n when T = 20. (h, τ ) m1 (h, τ ) rate1 m2 (h, τ ) rate3 (0.5, 0.1) 2.7187e-001 – 1.2056e-002 – (0.5/2, 0.1/16) 1.0797e-003 7.9761 4.5420e-005 8.0522 4.2375e-006 7.9932 1.7512e-007 8.0188 (0.5/22 , 0.1/162 )

Table2. Errors and convergence rates of Φn and U n with h = 0.2 when T = 20. τ m1 (h, τ ) rate2 m2 (h, τ ) rate4 0.01 2.7803e-003 – 1.0875e-004 – 0.005 6.9436e-004 2.0015 2.7555e-005 1.9806 0.0025 1.7275e-004 2.0070 7.2423e-006 1.9278

Table3. Discrete mass and energy computed by Scheme (2.10)–(2.14) with h = 0.1, τ = 0.005. tn Qn En 0 4.999999999999958 1.291407667737759 2 4.999999999999957 1.291445002215217 4 4.999999999999957 1.291444563749225 6 4.999999999999966 1.291408854486048 8 4.999999999999966 1.291443299744085

2.5 |φ| |Φ|

|φ| and |Φ|

2

1.5

t=0 t=5

1

t=10 t=15

0.5

0 -10

-5

0

5

10

15

20

25

30

x

Figure 1: Modulus at different levels computed by Scheme (2.10)–(2.14) with h = 0.1, τ = 0.05.

18

2.5 |u| |U|

2

|u| and |U|

t=0 1.5

t=5 t=10

1

t=15

0.5

0 -10

-5

0

5

10

15

20

25

30

x

Figure 2: Modulus at different levels computed by Scheme (2.10)–(2.14) with h = 0.1, τ = 0.05.

6. Conclusion In this paper, we propose a high-order linear compact finite difference scheme for the KGS equations and prove the conservations, convergence and stability of the new scheme. Numerical results demonstrate that the difference scheme is accurate and efficient. It is worth mentioning that our method can be directly extended to two dimensions and/or three dimensions, but some new techniques are required to be used to deal with the prior estimates and convergence of the difference solution.

Acknowledgment This work was supported by the National Natural Science Foundation of China (grant No.11872201, 11772148,). Declarations of interest: none.

Reference [1] I. Fukuda, M. Tsutsumi, On coupled Klein-Gordon-Schr¨ odinger eqautions, II, J. Math. Anal. Appl. 66 (2) (1978) 358–378. [2] Makhankov VG, Dynamics of classical solitons (in non-integrable systems), Phys. Rep.-Rev. Sec. Phys. Lett. 35 (1) (1978) 1–128. [3] X. Xiang, Spectral method for solving the system of equations of Schr¨odingerKlein-Gordon field, J. Comput. Appl. Math. 21 (2) (1988) 161–171.

19

[4] W. Bao, L. Yang, Efficient and accurate numerical methods for the KleinGordon-Schr¨ odinger equations, J. Comput. Phys. 225 (2) (2007) 1863–1893. [5] M. Dehghan, A. Taleei, Numerical solution of the Yukawa-coupled KleinGordon-Schr¨ odinger equations via a Chebyshev pseudospectral multidomain method, Appl. Math. Model. 36 (6) (2012) 2340–2349. [6] L. Kong, R. Liu, Z. Xu, Numerical simulation of interaction between Schr¨odinger field and Klein-Gordon field by multisymplectic method, Appl. Math. Comput. 181 (1) (2006) 342–350. [7] J. Hong, S. Jiang, C. Li, Explicit multi-symplectic methods for Klein-GordonSchr¨odinger equations, J. Comput. Phys. 228 (9) (2009) 3517–3532. [8] L. Zhang, Convergence of a conservative difference scheme for a class of KleinGordon-Schr¨ odinger equations in one space dimension, Appl. Math. Comput. 163 (1) (2005) 343–355. [9] T. Wang, Y. Jiang, Point-wise errors of two conservative difference schemes for the Klein-Gordon-Schr¨odinger equation, Commun. Nonlinear Sci. Numer. Simul. 17 (12) (2012) 4565–4575. [10] X. Pan, L. Zhang, High-order linear compact conservative method for the nonlinear Schr¨ odinger equation coupled with the nonlinear Klein-Gordon equation, Nonlinear Anal.-Theory Methods Appl. 92 (2013) 108–118. [11] T. Wang, Optimal point-wise error estimate of a compact difference scheme for the Klein-Gordon-Schr¨odinger equation, J. Math. Anal. Appl. 412 (1) (2014) 155–167. [12] J. Xia, M. Wang, Exact solitary solution of coupled Klein-Gordon-Schr¨odinger equations, Appl. Math. Mech. 23 (2002) 52–57. [13] G. Berikelashvili, M. M. Gupta, M. Mirianashvili, Convergence of fourth order compact difference schemes for three-dimensional convection-diffusion equations, SIAM J. Numer. Anal. 45 (1) (2007) 443–455. [14] G. Cohen, High-order numerical methods for Transient Wave equations, Springer, New York, 2002. [15] A. Gopaul, M. Bhuruth, Analysis of a fourth-order scheme for a threedimensional convection-diffusion model problem, SIAM J. Sci. Comput. 28 (6) (2006) 2075–2094.

20

[16] B. Gustafsson, E. Mossberg, Time compact high order difference methods for wave propagation, SIAM J. Sci. Comput. 26 (1) (2004) 259–271. [17] S. K. Lele, Compact finite difference schemes with spectral-like resolution, J. Comput. Phys. 103 (1) (1992) 16–42. [18] H. Liao, Z. Sun, H. Shi, Error estimate of fourth-order compact scheme for linear Schr¨ odinger equations, SIAM J. Numer. Anal. 47 (6) (2010) 4381–4401. [19] D. W. Zingg, Comparison of high-accuracy finite-difference methods for linear wave propagation, SIAM J. Sci. Comput. 22 (2) (2000) 476–502. [20] T. Wang, Convergence of an eighth-order compact difference scheme for the nonlinear Schr¨odinger equation, Adv. Numer. Anal. 2012 (2012) PP24. [21] Z. Sun, Numerical methods of the partial differential equations, Beiging: Science Press, 2005. [22] X. Hu, L. Zhang, Conservative compact difference schemes for the coupled nonlinear Schr¨odinger equation system, Numer. Meth. Part Differ. Equ. 30 (2014) 749–772. [23] T. Wang, B. Guo, Unconditional convergence of two conservative compact difference schemes for non-linear Schr¨odinger equation in one dimension (in Chinese), Sci. Sin. Math. 41 (2011) 1–27. [24] Y. Zhou, Application of discrete functional analysis to the finite difference method, International Academic Publishers, 1990.

21