# Convergence of an Eighth-Order Compact Difference Scheme for the Nonlinear Schrödinger Equation

Authors: Tingchun Wang
Journal: Advances in Numerical Analysis
 (2012)

---

## Abstract

A new compact difference scheme is proposed for solving the nonlinear Schrödinger equation. The scheme is proved to conserve the total mass and the total energy and the optimal convergent rate, without any restriction on the grid ratio, at the order ofO(h8+τ2) in the discrete L∞-norm with time step τ and mesh size h. In numerical analysis, beside the standard techniques of the energy method, a new technique named “regression of compactness” and some lemmas are proposed to prove the high-order convergence. For computing the nonlinear algebraical systems generated by the nonlinear compact scheme, an efficient iterative algorithm is constructed. Numerical examples are given to support the theoretical analysis.


---

## Body

## 1. Introduction

The time-dependent Schrödinger equation is one of the most important equations in quantum mechanics [1] such as the Bose-Einstein condensate (BEC), which is also used widely in many other different fields [2, 3] such as plasma physics, nonlinear optics, water waves, and bimolecular dynamics. There are many studies on numerical approaches, including finite difference [4–16], finite element [17–19], and polynomial approximation methods [20–24], of the initial or initial-boundary value problems of the Schrödinger equations.Recently, there has been growing interest in high-order compact methods for solving partial differential equations (PDEs) [25–35]. It was shown that the high-order difference methods play an important role in the simulation of high-frequency wave phenomena. However, there is few proof of unconditional H1 norm convergence of any compact difference scheme for nonlinear equations.In this paper, we consider the following cubic nonlinear Schrödinger (NLS) equation:(1.1)i∂tu(x,t)+12∂xxu(x,t)+V(x)u(x,t)-q|u(x,t)|2u(x,t)=0,x∈ℝ,t>0,
subject to periodic boundary condition
(1.2)u(x,t)=u(x+L,t),x∈ℝ,t>0,
and initial condition
(1.3)u(x,0)=φ(x),x∈ℝ,
where q is dimensionless constant characterizing the interaction (positive for repulsive interaction and negative for the attractive interaction) between particles in BEC. V(x) is a real function corresponding to the external trap potential and it is often chosen as a harmonic potential, that is, a quadratic polynomial, in most experiments. u:=u(x,t) is the unknown periodical complex-valued function whose initial value φ(x) is a given periodical complex-valued function.  L>0 is the period of the u(x,t).DenoteΩ=(0,L); it is easy to show that the problem (1.1)–(1.3) conserves the total mass(1.4)Q(u(·,t)):=∫Ω‍|u(x,t)|2dx≡Q(u(·,0))=Q(φ),t≥0,
and the global energy(1.5)E(u(·,t)):=14∫Ω‍[|∇u|2+2V(x)u+q|u|4]dx≡E(φ),t≥0.Fei et al. pointed out in [15] that the nonconservative schemes may easily show nonlinear blowup, and they presented a new conservative linear difference scheme for nonlinear Schrödinger equation. In [16], Li and Vu-Quoc also said, “… in some areas, the ability to preserve some invariant properties of the original differential equation is a criterion to judge the success of a numerical simulation.”In this paper, a high-order compact (HOC) difference scheme is proposed for solving the problem (1.1)–(1.3). On the proposed HOC scheme, we obtain the following novel results.(i)
The proposed scheme is eighth order in space and second order in time.(ii)
The new scheme conserves the discrete totalmass and the discrete global energy corresponding to (1.4) and (1.5), respectively.(iii)
A new technique named “regression of compactness” is introduced to analyze the convergence of the proposed.(iv)
An efficient iterative algorithm is constructed to compute the proposed scheme.The remainder of this paper is arranged as follows. In Section2, some notations are introduced, and a HOC difference scheme is constructed for the NLS equation. In Section 3, some useful lemmas are introduced or proved. In Section 4, discrete conservation laws of the proposed scheme are studied, then the convergence rate is proved at, without any restriction on the mesh ratio, the order of O(h8+τ2) in the discrete maximum norm. In Section 5, an efficient iterative algorithm is designed to solve the nonlinear scheme. Lastly, numerical experiments are presented in Section 6 and some remarks are given in the concluding section.


## 2. Notations and Compact Finite Difference Schemes

Before giving the compact difference scheme, some notations are firstly introduced. For a positive integerN, choose time-step τ=T/N and denote time steps tn=nτ, where 0<T<Tmax with Tmax the maximal existing time of the solution; choose mesh size h=L/J, with J a positive integer and denote grid points as
(2.1)xj=jh,j=0,1,…,J.
Denote the index sets
(2.2)𝒯J={j∣j=0,1,2,3,…,J-1},𝒯Je={j∣j=-2,-1,0,1,…,J+1}.
Denote Ujn to be the numerical approximation and, respectively, ujn be the exact value of u(x,t) at the point (xj,tn) for j∈𝒯Je and n=0,1,2,…,N, and denote Un∈ℂJ+3 to be the numerical solution and, respectively, un∈ℂJ+3 be the exact solution at time t=tn. For a grid function u={ujn∣j∈𝒯J,n=0,1,2,…,N}, introduce the following finite difference operators:
(2.3)δx+ujn=uj+1n-ujnh,δx-ujn=ujn-uj-1nh,δx〈1〉ujn=uj+1n-uj-1n2h,δx2ujn=δx+δx-ujn,δx〈2〉ujn=δx〈1〉δx〈1〉ujn,δt+ujn=ujn+1-ujnτ,We denote the space(2.4)XJ:={u=(uj∣j∈𝒯Je)∣uj=uj+Jwhenj=-2,-1,0,1}⊆ℂJ+3
and define seminorms and discrete inner product over XJ as
(2.5)∥u∥pp=h∑j=0J-1‍|uj|p,0<p<∞,∥u∥∞=maxj∈𝒯J|uj|,(u,v)=h∑j=0J-1‍ujvj¯.For simplicity, we write∥v∥2 as ∥v∥. Throughout the paper, we adopt the standard Sobolev spaces and their corresponding norms let C denote a generic constant which is independent of mesh size h and time step τ, and use the notation w≲v to present that there is a generic constant C which is independent of time step τ and mesh size h such that w⩽Cv. For the exact solution of the initial-boundary value problem (1.1)–(1.3), we assume that
(2.6)∥un∥⩽C,∥δx+un∥⩽C,∥δx〈1〉un∥⩽C,∥un∥∞⩽C.


### 2.1. Some Properties of Cyclic Matrix

A matrix in the form of(2.7)A=(a0a1a2⋯an-2an-1an-1a0a1⋯an-3an-2an-2an-1a0⋯an-4an-3⋮⋮⋮⋯⋮⋮a1a2a3⋯an-1a0)
is called a cyclic matrix. Because the matrix A is determined by the entries in the first row, the matrix can be denoted as
(2.8)A=C(a0,a1,…,an-1).
Denote G=C(0,1,0,…,0), which is called as a basic cyclic matrix. By a simple calculation, we can get that
(2.9)Gk=C(0,…,0︸k,1,0,…,0),k=0,1,…,n-1,Gn=C(1,0,…,0)=I,Iistheunitmatrix.
Then, any cyclic matrix A=C(a0,a1,…,an-1) can be written as a polynomial of G; that is,
(2.10)A=a0I+a1G+a2G2+⋯+an-1Gn-1.
This implies that, for any two cyclic matrices A and B, AB=BA. Now we list some useful lemmas as follows.Lemma 2.1.
If A and B are two cyclic matrices, AB is still a cyclic matrix.Proof.
For any two cyclic matrices(2.11)A=C(a0,a1,…,an-1),B=C(b0,b1,…,bn-1),
by using the multiplication of polynomials and the properties  we obtain
(2.12)AB=(a0I+a1G+a2G2+⋯+an-1Gn-1)(b0I+b1G+b2G2+⋯+bn-1Gn-1)=c0I+c1G+a2G2+⋯+cn-1Gn-1=C(c0,c1,…,cn-1),
where
(2.13)c0=a0b0+a1bn-1+⋯+an-1b1,c1=a0b1+a1b0+⋯+an-1b2,⋮cn-1=a0bn-1+a1bn-2+⋯+an-1b0.Lemma 2.2.
If a cyclic matrixA is invertible, then its inverse matrix A-1 is still a cyclic matrix.Proof.
IfA=C(a0,a1,…,an-1) is invertible, we assume that the inverse matrix is B=C(b0,b1,…,bn-1), then
(2.14)(a0I+a1G+a2G2+⋯+an-1Gn-1)B=I.
Multiplying b0I+b1G+b2G2+⋯+bn-1Gn-1 with (2.14) yields
(2.15)(c0I+c1G+c2G2+⋯+cn-1Gn-1)B=b0I+b1G+a2G2+⋯+bn-1Gn-1.
If (c0,c1,c2,…,cn-1)=(1,0,0,…,0), we obtain the linear system . Noting that |A|≠0, we obtain that the system  has a unique nonzero solution (b0,b1,…,bn-1). This together with (2.15) gives that B=b0I+b1G+b2G2+⋯+bn-1Gn-1 where (b0,b1,…,bn-1) is the unique nonzero solution of . It follows that the assumption is right; that is, B=C(b0,b1,…,bn-1) is a cyclic matrix.Lemma 2.3.
For any a real cyclic matrixA=C(a0,a1,…,an-1), all the eigenvalues are in the form of
(2.16)f(ɛk),k=0,1,2,…,n-1,
where f(x)=a0+a1x+a2x2+⋯+an-1xn-1 and ɛk=cos(2kπ/n)+isin(2kπ/n).Proof.
On the one hand, it follows from|ɛI-G|=0 that ɛn=1 which indicates that the eigenvalues of G are
(2.17)ɛk=ei(2kπ/n)=cos(2kπn)+isin(2kπn),k=0,1,2,…,n-1.
On the other hand, just as (2.10), the real cyclic matrix A can be written as
(2.18)A=f(G)=a0I+a1G+a2G2+⋯+an-1Gn-1,
where a0,a1,…,an-1 are n real numbers. This gives that the eigenvalues of A are in the form of
(2.19)λk=f(ɛk),k=0,1,2,…,n-1.


### 2.2. Compact Approximation of Second Derivative

The key point of the HOC method is to discretize the derivative with the fewest nodes to get the expected accuracy. In order to achieve higher accuracy, we adopt implicit compact scheme which equals a combination of the nodal derivatives to a combination of the nodal values of the function. Next, by Taylor’s formula, expanding all the nodal derivatives and nodal values of the function at the same node, we can determine the coefficients of the combinations under the constrains of accuracy order.The nodal derivatives are implicitly evaluated herein by solving linear algebraic equations in strip. The stencil and the weighted factors are restricted to be symmetric so that the resulting schemes are nondissipative. Supposef(x) is a L-periodic function, then for the discretization of the second-order derivatives fxx(x), we have the formula [31]:
(2.20)β(fxx)j-2+α(fxx)j-1+(fxx)j+α(fxx)j+1+β(fxx)j+2=bδx〈2〉fj+aδx2fj+O(hp),
where fj=f(xj), (fxx)j=fxx(xj), a, b, α, and β are undetermined parameters which depend on the accuracy order constraints. After tedious calculation, we can get these constraints:
(2.21)a+b=1+2α+2β,p=2,a+22b=4!2!(α+22β),p=4,a+24b=6!4!(α+24β),p=6,a+26b=8!6!(α+26β),p=8.
According to (2.20), an eighth-order scheme is obtained with
(2.22)α=3441179,β=232358,a=9601179,b=9301179,
which is the unique solution of the constraint equation .Thus, the approximation (2.20) to (fxx)j can be written in the form of finite-dimensional linear operator
(2.23)𝒜h(fxx)j=bδx^δx^fj+aδx+δx-fj+O(h8)⇒(fxx)j=b𝒜h-1δx^δx^fj+a𝒜h-1δx+δx-fj+O(h8),j∈𝒯J,
where the finite difference operator 𝒜h is defined as
(2.24)𝒜hfj:=βfj-2+αfj-1+fj+αfj+1+βfj+2.
The corresponding matrix form is
(2.25)Af~xx=bδx〈2〉f~+aδx2f~,
where
(2.26)f~xx=((fxx)0,(fxx)1,…,(fxx)J-1)⊤,f~=(f0,f1,…,fJ-1)⊤,A=C(1,α,β,0,…,0,β,α),
is a positive, symmetric, and cyclic pentadiagonal matrix.Lemma 2.4 (see [36]).
If a matrix H has real entries and is symmetric and positive definite, then(a)
H is invertible and its inverse matrix A-1 also has real entries and is symmetric and positive definite;(b)
H can be decomposed as
(2.27)H=D2,where D also has real entries and is symmetric and positive definite.


### 2.3. High-Order Compact Scheme

Imposing the approximations given in Section2.2 on the spatial direction and the Crank-Nicolson discretization on the temporal direction of the NLS equation (1.1) gives
(2.28)iδt+ujn+b4𝒜h-1δx〈2〉(ujn+ujn+1)+a4𝒜h-1δx2(ujn+ujn+1)+12Vj(ujn+ujn+1)+q4(|ujn|2+|ujn+1|2)(ujn+ujn+1)=ξjn,j∈𝒯J,0⩽n⩽N-1,(2.29)uj0=φ(xj),j∈𝒯Je,(2.30)un∈XJ,0⩽n⩽N,
where Vj=V(xj) and ξjn is the local truncation error. The corresponding matrix form of the (2.28)–(2.30) is
(2.31)iδtun+b4Mδx〈2〉(un+un+1)+a4Mδx2(un+un+1)+12V(un+un+1)+q4(|un|2+|un+1|2)(un+un+1)=ξn,0⩽n⩽N-1,uj0=φ(xj),j∈𝒯Je,un∈XJ,0⩽n⩽N,
where V=diag(V1,V2,…,VJ) and M=A-1.On the local truncation errorξjn, using Taylor's expansion, we obtain the follwing lemma.Lemma 2.5.
Suppose thatφ(x)∈H1(ℝ) and u(x,t)∈C10,4(ℝ×(0,T]), there is
(2.32)|ξjn|⩽C(h8+τ2),j∈𝒯J,0⩽n⩽N-1,|δtξjn|⩽C(h8+τ2),j∈𝒯J,0⩽n⩽N-2.Omitting the small termξjn in (2.28) yields the following difference scheme
(2.33)iδt+Ujn+b4𝒜h-1δx〈2〉(Ujn+Ujn+1)+a4𝒜h-1δx2(Ujn+Ujn+1)+12Vj(Ujn+Ujn+1)+q4(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)=0,j∈𝒯J,0⩽n⩽N-1,Uj0=φ(xj),j∈𝒯J,Un∈XJ,0⩽n⩽N;
that is,
(2.34)i𝒜hδtUjn+b4δx〈2〉(Ujn+Ujn+1)+a4δx2(Ujn+Ujn+1)+12𝒜h[Vj(Ujn+Ujn+1)]+q4𝒜h[(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)]=0,j∈𝒯J,0⩽n⩽N-1,Uj0=φ(xj),j∈𝒯J,Un∈XJ,0⩽n⩽N.Remark 2.6.
Though and  are two different forms of a same scheme, the latter is suitable for computing in implementation.


## 2.1. Some Properties of Cyclic Matrix

A matrix in the form of(2.7)A=(a0a1a2⋯an-2an-1an-1a0a1⋯an-3an-2an-2an-1a0⋯an-4an-3⋮⋮⋮⋯⋮⋮a1a2a3⋯an-1a0)
is called a cyclic matrix. Because the matrix A is determined by the entries in the first row, the matrix can be denoted as
(2.8)A=C(a0,a1,…,an-1).
Denote G=C(0,1,0,…,0), which is called as a basic cyclic matrix. By a simple calculation, we can get that
(2.9)Gk=C(0,…,0︸k,1,0,…,0),k=0,1,…,n-1,Gn=C(1,0,…,0)=I,Iistheunitmatrix.
Then, any cyclic matrix A=C(a0,a1,…,an-1) can be written as a polynomial of G; that is,
(2.10)A=a0I+a1G+a2G2+⋯+an-1Gn-1.
This implies that, for any two cyclic matrices A and B, AB=BA. Now we list some useful lemmas as follows.Lemma 2.1.
If A and B are two cyclic matrices, AB is still a cyclic matrix.Proof.
For any two cyclic matrices(2.11)A=C(a0,a1,…,an-1),B=C(b0,b1,…,bn-1),
by using the multiplication of polynomials and the properties  we obtain
(2.12)AB=(a0I+a1G+a2G2+⋯+an-1Gn-1)(b0I+b1G+b2G2+⋯+bn-1Gn-1)=c0I+c1G+a2G2+⋯+cn-1Gn-1=C(c0,c1,…,cn-1),
where
(2.13)c0=a0b0+a1bn-1+⋯+an-1b1,c1=a0b1+a1b0+⋯+an-1b2,⋮cn-1=a0bn-1+a1bn-2+⋯+an-1b0.Lemma 2.2.
If a cyclic matrixA is invertible, then its inverse matrix A-1 is still a cyclic matrix.Proof.
IfA=C(a0,a1,…,an-1) is invertible, we assume that the inverse matrix is B=C(b0,b1,…,bn-1), then
(2.14)(a0I+a1G+a2G2+⋯+an-1Gn-1)B=I.
Multiplying b0I+b1G+b2G2+⋯+bn-1Gn-1 with (2.14) yields
(2.15)(c0I+c1G+c2G2+⋯+cn-1Gn-1)B=b0I+b1G+a2G2+⋯+bn-1Gn-1.
If (c0,c1,c2,…,cn-1)=(1,0,0,…,0), we obtain the linear system . Noting that |A|≠0, we obtain that the system  has a unique nonzero solution (b0,b1,…,bn-1). This together with (2.15) gives that B=b0I+b1G+b2G2+⋯+bn-1Gn-1 where (b0,b1,…,bn-1) is the unique nonzero solution of . It follows that the assumption is right; that is, B=C(b0,b1,…,bn-1) is a cyclic matrix.Lemma 2.3.
For any a real cyclic matrixA=C(a0,a1,…,an-1), all the eigenvalues are in the form of
(2.16)f(ɛk),k=0,1,2,…,n-1,
where f(x)=a0+a1x+a2x2+⋯+an-1xn-1 and ɛk=cos(2kπ/n)+isin(2kπ/n).Proof.
On the one hand, it follows from|ɛI-G|=0 that ɛn=1 which indicates that the eigenvalues of G are
(2.17)ɛk=ei(2kπ/n)=cos(2kπn)+isin(2kπn),k=0,1,2,…,n-1.
On the other hand, just as (2.10), the real cyclic matrix A can be written as
(2.18)A=f(G)=a0I+a1G+a2G2+⋯+an-1Gn-1,
where a0,a1,…,an-1 are n real numbers. This gives that the eigenvalues of A are in the form of
(2.19)λk=f(ɛk),k=0,1,2,…,n-1.


## 2.2. Compact Approximation of Second Derivative

The key point of the HOC method is to discretize the derivative with the fewest nodes to get the expected accuracy. In order to achieve higher accuracy, we adopt implicit compact scheme which equals a combination of the nodal derivatives to a combination of the nodal values of the function. Next, by Taylor’s formula, expanding all the nodal derivatives and nodal values of the function at the same node, we can determine the coefficients of the combinations under the constrains of accuracy order.The nodal derivatives are implicitly evaluated herein by solving linear algebraic equations in strip. The stencil and the weighted factors are restricted to be symmetric so that the resulting schemes are nondissipative. Supposef(x) is a L-periodic function, then for the discretization of the second-order derivatives fxx(x), we have the formula [31]:
(2.20)β(fxx)j-2+α(fxx)j-1+(fxx)j+α(fxx)j+1+β(fxx)j+2=bδx〈2〉fj+aδx2fj+O(hp),
where fj=f(xj), (fxx)j=fxx(xj), a, b, α, and β are undetermined parameters which depend on the accuracy order constraints. After tedious calculation, we can get these constraints:
(2.21)a+b=1+2α+2β,p=2,a+22b=4!2!(α+22β),p=4,a+24b=6!4!(α+24β),p=6,a+26b=8!6!(α+26β),p=8.
According to (2.20), an eighth-order scheme is obtained with
(2.22)α=3441179,β=232358,a=9601179,b=9301179,
which is the unique solution of the constraint equation .Thus, the approximation (2.20) to (fxx)j can be written in the form of finite-dimensional linear operator
(2.23)𝒜h(fxx)j=bδx^δx^fj+aδx+δx-fj+O(h8)⇒(fxx)j=b𝒜h-1δx^δx^fj+a𝒜h-1δx+δx-fj+O(h8),j∈𝒯J,
where the finite difference operator 𝒜h is defined as
(2.24)𝒜hfj:=βfj-2+αfj-1+fj+αfj+1+βfj+2.
The corresponding matrix form is
(2.25)Af~xx=bδx〈2〉f~+aδx2f~,
where
(2.26)f~xx=((fxx)0,(fxx)1,…,(fxx)J-1)⊤,f~=(f0,f1,…,fJ-1)⊤,A=C(1,α,β,0,…,0,β,α),
is a positive, symmetric, and cyclic pentadiagonal matrix.Lemma 2.4 (see [36]).
If a matrix H has real entries and is symmetric and positive definite, then(a)
H is invertible and its inverse matrix A-1 also has real entries and is symmetric and positive definite;(b)
H can be decomposed as
(2.27)H=D2,where D also has real entries and is symmetric and positive definite.


## 2.3. High-Order Compact Scheme

Imposing the approximations given in Section2.2 on the spatial direction and the Crank-Nicolson discretization on the temporal direction of the NLS equation (1.1) gives
(2.28)iδt+ujn+b4𝒜h-1δx〈2〉(ujn+ujn+1)+a4𝒜h-1δx2(ujn+ujn+1)+12Vj(ujn+ujn+1)+q4(|ujn|2+|ujn+1|2)(ujn+ujn+1)=ξjn,j∈𝒯J,0⩽n⩽N-1,(2.29)uj0=φ(xj),j∈𝒯Je,(2.30)un∈XJ,0⩽n⩽N,
where Vj=V(xj) and ξjn is the local truncation error. The corresponding matrix form of the (2.28)–(2.30) is
(2.31)iδtun+b4Mδx〈2〉(un+un+1)+a4Mδx2(un+un+1)+12V(un+un+1)+q4(|un|2+|un+1|2)(un+un+1)=ξn,0⩽n⩽N-1,uj0=φ(xj),j∈𝒯Je,un∈XJ,0⩽n⩽N,
where V=diag(V1,V2,…,VJ) and M=A-1.On the local truncation errorξjn, using Taylor's expansion, we obtain the follwing lemma.Lemma 2.5.
Suppose thatφ(x)∈H1(ℝ) and u(x,t)∈C10,4(ℝ×(0,T]), there is
(2.32)|ξjn|⩽C(h8+τ2),j∈𝒯J,0⩽n⩽N-1,|δtξjn|⩽C(h8+τ2),j∈𝒯J,0⩽n⩽N-2.Omitting the small termξjn in (2.28) yields the following difference scheme
(2.33)iδt+Ujn+b4𝒜h-1δx〈2〉(Ujn+Ujn+1)+a4𝒜h-1δx2(Ujn+Ujn+1)+12Vj(Ujn+Ujn+1)+q4(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)=0,j∈𝒯J,0⩽n⩽N-1,Uj0=φ(xj),j∈𝒯J,Un∈XJ,0⩽n⩽N;
that is,
(2.34)i𝒜hδtUjn+b4δx〈2〉(Ujn+Ujn+1)+a4δx2(Ujn+Ujn+1)+12𝒜h[Vj(Ujn+Ujn+1)]+q4𝒜h[(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)]=0,j∈𝒯J,0⩽n⩽N-1,Uj0=φ(xj),j∈𝒯J,Un∈XJ,0⩽n⩽N.Remark 2.6.
Though and  are two different forms of a same scheme, the latter is suitable for computing in implementation.


## 3. Some Useful Lemmas

Lemma2.3 gives that the eigenvalues of the cyclic matrix A=C(1,α,β,0,…,0,β,α) are in the form of
(3.1)λj=1+2αcos(2jπJ)+2βcos(4jπJ),j∈𝒯J.
This gives
(3.2)1-2α-2β⩽λj⩽1+2α+2β.
It is easy to verify that 1-2α-2β>0, which indicates that the matrix A is positive definite.For the real, positive definite, symmetric and cyclic matrixA=C(1,α,β,0,…,0,β,α), it follows from Lemma 2.2 that there are J real numbers m0,m1,…,mJ-1 such that M=C(m0,m1,…,mJ-1). Lemmas 2.3 and 2.4 give that M is real positive definite and symmetric. And then we can introduce the following discrete norm:
(3.3)|||u|||:=(Mu,u)1/2.
On the relation between |||·||| and ∥·∥, we have the following lemma.Lemma 3.1.
The discrete norms|||·||| and ∥·∥ are equivalent; in fact, there is
(3.4)c1∥u∥⩽|||u|||⩽c2∥u∥,foranyagridfunctionu∈XJ,
where c1=1/1+2α+2β, c2=1/1-2α-2β.Proof.
It follows from (3.2) that the eigenvalues of M satisfy
(3.5)11+2α+2β⩽μj⩽11-2α-2β,j∈𝒯J.
This gives the spectral radius ρ(M)⩽(1/(1-2α-2β)), and consequently
(3.6)∥M∥=ρ(M)⩽11-2α-2β.
Because
(3.7)11+2α+2β(u,u)⩽(Mu,u)⩽∥M∥·∥u∥2⩽11-2α-2β(u,u),
then (3.5) follows from (3.3) and (3.7).The following lemmas will be used frequently in this paper.Lemma 3.2.
For any a grid functionu∈XJ, there are
(3.8)δx+(Mu)=Mδx+u,(3.9)δx〈1〉(Mu)=Mδx〈1〉u.Proof.
For thekth entry of the vector Mu, there is
(3.10)(Mu)k=∑j=0J-1‍mjuj+k,
where the periodicity of u is used. This gives
(3.11)δx+(Mu)k=1h(∑j=0J-1‍mjuj+k+1-∑j=0J-1‍mjuj+k)=∑j=0J-1‍mj(uj+k+1-uj+k)h.
For the kth element of the vector Mδx+u, there is
(3.12)(Mδx+u)k=∑j=0J-1‍mj(uj+k-uj+k-1)h.
Then (3.8) is gotten from (3.11) and (3.12). Similarly, (3.9) can be proved.Lemma 3.3.
For any two grid functionsu,v∈XJ; there are
(3.13)(δx2u,v)=-(δx+u,δx+v),(δx〈2〉u,v)=-(δx〈1〉u,δx〈1〉v).Proof.
Direct calculation and noticing the periodicity of the two grid functions yield.Lemma 3.4.
For the approximationun,un+1∈XJ, there are
(3.14)ℑ(Mδx2(un+un+1),un+un+1)=0,(3.15)ℑ(Mδx〈2〉(un+un+1),un+un+1)=0,(3.16)ℜ(Mδx2(un+un+1),δt+un)=1τ(|||δx+un+1|||2-|||δx+un|||2),(3.17)ℜ(Mδx〈2〉(un+un+1),δt+un)=1τ(|||δx〈1〉un+1|||2-|||δx〈1〉un|||2),
where “ℑ(·)” and “ℜ(·)” mean taking the imaginary part and the real part, respectively.Proof.
Lemmas3.2, 3.3, and definition (3.3) together with the symmetry of M give
(3.18)ℑ(Mδx2(un+un+1),un+un+1)=ℑ(δx2(un+un+1),M(un+un+1))=ℑ(δx+(un+un+1),δx+M(un+un+1))=ℑ(δx+(un+un+1),Mδx+(un+un+1))=ℑ(Mδx+(un+un+1),δx+(un+un+1))=ℑ|||δx+(un+un+1)|||2=0.
Similarly, (3.15) can be proved. Also using Lemmas 3.2, and 3.3, definition (3.3), and the symmetry of M, we obtain
(3.19)ℜ(Mδx2(un+un+1),δt+un)=1τℜ(δx+un+1+δx+un,Mδx+un+1-Mδx+un)=1τ(|||δx+un+1|||2-|||δx+un|||2).
Similarly, (3.17) can be proved.Lemma 3.5 (see [33]).
For time sequencesw={w0,w1,…,wn,wn+1} and g={g0,g1,…,gn-1,gn}, there is
(3.20)|2τ∑l=0n‍glδtwl|⩽|w0|2+τ∑l=1n‍|wl|2+|wn+1|2+|g0|2+τ∑l=0n-1‍|δtgl|2+|gn|2.


## 4. Numerical Analysis

The corresponding matrix form of is
(4.1)iδtUn+b4Mδx〈2〉(Un+Un+1)+a4Mδx2(Un+Un+1)+12V(Un+Un+1)+q4(|Un|2+|Un+1|2)(Un+Un+1)=0,0⩽n⩽N-1,(4.2)Uj0=φ(xj),j∈𝒯Je,(4.3)Un∈XJ,0⩽n⩽N,
where
(4.4)A=C(1,α,β,0,…,0,β,α),((|Un|2+|Un+1|2)(Un+Un+1))j=(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1).


### 4.1. Solvability of the Nonlinear Scheme

For proving the solvability, we need the following Brouwer fixed point theorem.Lemma 4.1 (see [37]).
Let(H,〈·,·〉) be a finite-dimensional inner-product space, ∥·∥ the associated norm, and g:H→H continuous. Assume, moreover, that
(4.5)∃α>0,∀z∈H,∥z∥=α,ℜ〈g(z),z〉>0.
Then, there exists a z*∈H such that g(z*)=0 and ∥z*∥⩽α.On the solvability of the scheme, we have the following lemma.Lemma 4.2 (Solvability of difference scheme).
Under the assumption(A), for any initial data u0∈XJ, there exists a solution Un∈XJ of  for n≥1.Proof.
The HOC scheme is equivalent to the difference scheme  whose matrix form is (4.1)–(4.3), so we here just prove the unique solvability of the difference scheme (4.1)–(4.3).
In (4.1), for given Un∈XJ(n≥1), we first prove the existence. For j∈𝒯J, rewrite (4.1) as
(4.6)Un+(1/2)=Un+ib4τMδx〈2〉Un+(1/2)+ia4τMδx2Un+(1/2)+iτ2VUn+(1/2)+iq4τ(|Un|2+|2Un+(1/2)-Un|2)Un+(1/2).
We define a mapping ℱ: u*∈XJ→ℱ(u*)∈XJ as
(4.7)ℱ(U*)=U*-Un-ib4τMδx〈2〉U*-ia4τMδx2U*-iτ2VU*-iq4τ(|Un|2+|2U*-Un|2)U*,
which is obviously continue. Computing the inner product of (4.7) with u* and then taking the real part yield
(4.8)ℜ(ℱ(U*),U*)=∥U*∥2-(Un,U*)≥∥U*∥(∥U*∥-∥Un∥),
where Lemma 3.4 was used. Let ∥U*∥=∥Un∥+1, then it follows from (4.8) that ℜ(ℱ(U*),U*)≥0. This together with Lemma 4.1 indecaites that there exists a solution Un+1∈XJ satisfying ℱ(Un+1)=0. Then Un+1 satisfies the scheme .


### 4.2. Discrete Conservation Laws

Corresponding to (1.4)-(1.5), the proposed scheme also possesses two counterpart discrete conservation laws.Lemma 4.3.
The difference solution of scheme conserves the discrete mass and discrete energy in the form of
(4.9)Qn=⋯=Q0,n=0,1,2,…,N,(4.10)En=⋯=E0,n=0,1,2,…,N,
where
(4.11)Qn≔∥Un∥2,En≔b4|||δx〈1〉Un|||2+a4|||δx+Un|||2+12h∑j=0J-1‍Vj|Un|2+q4∥Un∥44.Proof.
Computing the inner product of (4.1) with Un+Un+1 then taking the imaginary part yields
(4.12)1τ(∥Un+1∥2-∥Un∥2)+b4ℑ(Mδx〈2〉(Un+Un+1),(Un+Un+1))+a4ℑ(Mδx2(Un+Un+1),Un+Un+1)+12ℑ(h∑j=0J-1‍Vj(|Un+Un+1|2))+q4ℑ((|Un|2+|Un+1|2)(Un+Un+1),Un+Un+1)=0.
This together with Lemma 3.4 gives (4.9).
Computing the inner product of (4.1) with δt+Un then taking the real part yields
(4.13)ℜ(iδt+Un,δt+Un)+b4ℜ(Mδx〈2〉(Un+Un+1),δt+Un)+a4ℜ(Mδx2(Un+Un+1),δt+Un)+12ℜ(h∑j=0J-1‍Vj(|Un+1|2-|Un|2))+q4ℜ((|Un|2+|Un+1|2)(Un+Un+1),δt+Un)=0.
This together with
(4.14)((|Un|2+|Un+1|2)(Un+Un+1),δt+Un)=(∥Un+1∥44-∥Un∥44)
and Lemma 3.4 gives (4.10).


### 4.3. A Priori Estimate

In proof of the convergence, we need the followinga priori estimate.Lemma 4.4.
The difference solution of the scheme satisfies
(4.15)∥Un∥⩽C,∥δx〈1〉Un∥⩽C,∥δx+Un∥⩽C,∥Un∥∞⩽C.
This implies the stability of the scheme .Proof.
It follows from (4.9) that
(4.16)∥Un∥⩽C.
Lemma 3.1 gives
(4.17)|||δx+Un|||2≥c12∥δx+Un∥2.
Using the discrete Sobolev inequality and (4.16), we obtain
(4.18)∥Un∥44⩽c04(∥δx+Un∥1/4∥Un∥3/4+∥Un∥)4⩽8c04(∥δx+Un∥∥Un∥3+∥Un∥4)⩽ɛ∥δx+Un∥2+16c08ɛ∥Un∥6+8c04∥Un∥4⩽ɛ∥δx+Un∥2+C,
where c0 is a positive constant independent of the grid parameters. Taking ɛ=ac12/2|q|, then it follows from (4.16)–(4.18) and (4.10) that
(4.19)∥δx+Un∥⩽C,∥δx〈1〉Un∥⩽C.
Utilizing the discrete Sobolev estimate and Lemma 4.3, we obtain from (4.16) and (4.19) that
(4.20)∥Un∥∞⩽C.
This completes the proof.


### 4.4. Unconditional Convergence in Maximum Norm

Theorem 4.5.
Suppose thatu(x,t)∈C10,4(R×(0,T]) and τ is sufficiently small, then the solution ujn of the proposed scheme  converges to the solution u(xj,tn) of the problem (1.1)–(1.3) with order O(h8+τ2) in the discrete L∞-norm.Proof.
Denote the error function as(4.21)en=un-Un,n=0,1,2,…,N.
Subtracting (4.1)–(4.3) from , respectively, we obtain the equation of error function as follows:
(4.22)iδt+en+b4Mδx〈2〉(en+en+1)+a4Mδx2(en+en+1)+12Vj(en+en+1)+q4gn+1=ξn,n=0,1,…,N-1,(4.23)ej0=0,j∈𝒯J,(4.24)en∈XJ,n=0,1,2,…,N.
where
(4.25)gjn+1=(|ujn|2+|ujn+1|2)(ujn+ujn+1)-(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1),j∈𝒯J.
Noting that
(4.26)gjn+1=(|ujn|2+|ujn+1|2)(ujn+ujn+1)-(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)=(ejnUjn¯+ujnejn¯+ejn+1Ujn+1¯+ujn+1ejn+1¯)(ujn+ujn+1)+(|Ujn|2+|Ujn+1|2)(ejn+ejn+1)=(ejn(ujn¯-ejn¯)+ujnejn¯+ejn+1(ujn+1¯-ejn+1¯)+ujn+1ejn+1¯)×(ujn+ujn+1)+(|ujn-ejn|2+|ujn+1-ejn+1|2)(ejn+ejn+1),j∈𝒯J,
we obtain from (4.25), (2.6), and Lemma 4.4 that
(4.27)||gn+1||2⩽C(||en||2+||en+1||2).
Computing the inner product of (4.22) with en+en+1 and then taking the imaginary part, we obtain
(4.28)1τ(∥en+1∥2-∥en∥2)+q4ℑ(gn+1,en+en+1)=ℑ(ξn,en+en+1),
where Lemma 3.4 was used.
For the last two terms of (4.28), using Lemma 4.4 and (4.26) we have
(4.29)|q4ℑ(gn+1,en+en+1)|≲∥en∥2+∥en+1∥2,|ℑ(ξn,en+en+1)|≲∥ξn∥2+∥en∥2+∥en+1∥2.
Substituting  into (4.28) yields
(4.30)1τ(∥en+1∥2-∥en∥2)≲∥ξn∥2+∥en∥2+∥en+1∥2.
By virtue of Gronwall's inequality and Lemma 2.5, noting that ∥e0∥=0, we obtain from (4.30) that
(4.31)∥en∥⩽C(h8+τ2),n=1,2,…,N,
if τ is sufficiently small.
Computing the inner product of (4.22) with δt+en and then taking the real part, we obtain
(4.32)b4ℜ(Mδx〈2〉(en+en+1),δt+en)+a4ℜ(Mδx2(en+en+1),δt+en)+q4ℜ(gn+1,δt+en)=ℜ(ξn,δt+en),n=1,2,…,N-1.
This together with Lemma 3.4 gives
(4.33)-b4τ(|||δx〈1〉un+1|||2-|||δx〈1〉un|||2)-a4τ(|||δx+un+1|||2-|||δx+un|||2)+q4ℜ(gn+1,δt+en)=ℜ(ξn,δt+en).
It follows from (4.22) that
(4.34)δt+en=ib4Mδx〈2〉(en+en+1)+ia4Mδx2(en+en+1)+i2Vj(en+en+1)+iq4gn+1-iξn,n=1,2,…,N-1.
Thus, for the third term on the left side of (4.33), we have
(4.35)q4ℜ(gn+1,δt+en)=q4ℜ(gn+1,ib4Mδx〈2〉(en+en+1)+ia4Mδx2(en+en+1)+i2Vj(en+en+1)+iq4gn+1-iξn)=-q4ℑ(gn+1,b4Mδx〈2〉(en+en+1))-q4ℑ(gn+1,a4Mδx2(en+en+1))-q4ℑ(gn+1,12V(en+en+1))+q4ℑ(gn+1,ξn).
Noting (4.26), (4.31), and
(4.36)δx+(ejn(ujn¯-ejn¯))=δx+ejn(uj+1n¯-ej+1n¯)+ejn(δx+ujn¯-δx+ejn¯)=δx+ejnUj+1n¯+ejn(δx+ujn¯-δx+ejn¯),δx+((ujn-ejn)(ujn¯-ejn¯))=(δx+ujn-δx+ejn)(uj+1n¯-ej+1n¯)+(ujn-ejn)(δx+ujn¯-δx+ejn¯)=(δx+ujn-δx+ejn)Uj+1n¯+Ujn(δx+ujn¯-δx+ejn¯),|(δx+enen,en)|⩽∥δx+en∥2+∥en∥44⩽2∥δx+en∥2+C(∥en∥6+∥en∥4),
we obtain
(4.37)|q4ℑ(gn+1,a4Mδx2(en+en+1))|=|-q4ℑ(δx+gn+1,a4Mδx+(en+en+1))|≲∥M∥|(δx+gn+1,δx+(en+en+1))|≲|(δx+gn+1,δx+(en+en+1))|(4.38)≲∥δx+en∥2+∥δx+en+1∥2+(en+en+1)2.
Similarly we obtain
(4.39)|q4ℑ(gn+1,b4Mδx〈2〉(en+en+1))|≲∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+(h8+τ2)2.
For the last two terms of (4.35), using the Cauchy inequality, noting that V(x)∈C(Ω), we obtain
(4.40)|-q4ℑ(gn+1,12V(en+en+1))+q4ℑ(gn+1,ξn)|≲∥en∥2+∥en+1∥2+∥ξn∥2≲(h8+τ2)2.
Substituting (4.38)–(4.40) into (4.35) yields
(4.41)|q4ℜ(gn+1,δt+en)|≲∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+∥δx1en∥2+∥δx1en+1∥2+(h8+τ2)2.
Substituting (4.41) into (4.33) yields
(4.42)|||δx〈1〉un+1|||2-|||δx〈1〉un|||2+|||δx+un+1|||2-|||δx+un|||2≲τ[∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+∥δx1en∥2+∥δx1en+1∥2+(h8+τ2)2-ℜ(ξn,δt+en)].
Summing up for the superscript n from 0 to m and then replacing m by n and noting ∥|δx〈1〉u0|∥=∥|δx1u0|∥=0, we obtain
(4.43)|||δx〈1〉un+1|||2+|||δx+un+1|||2≲τ∑l=0n+1‍[∥δx〈1〉el∥2+∥δx1el∥2]+|τ∑l=0n+1‍(ξl,δt+el)|+(h8+τ2)2.
By virtue of Lemmas 3.5, 2.5, and (4.31), we obtain
(4.44)|τ∑l=0n+1‍(ξl,δt+el)|≲(h8+τ2)2.
Substituting (4.44) into (4.43) and then utilizing Gronwall's inequality, we obtain
(4.45)∥δx〈1〉en∥≲h8+τ2,∥δx+en∥≲h8+τ2,n=1,2,…,N,
if τ is sufficiently small. This together with (4.31) gives
(4.46)∥en∥∞≲h8+τ2,n=1,2,…,N.
This completes the proof.


## 4.1. Solvability of the Nonlinear Scheme

For proving the solvability, we need the following Brouwer fixed point theorem.Lemma 4.1 (see [37]).
Let(H,〈·,·〉) be a finite-dimensional inner-product space, ∥·∥ the associated norm, and g:H→H continuous. Assume, moreover, that
(4.5)∃α>0,∀z∈H,∥z∥=α,ℜ〈g(z),z〉>0.
Then, there exists a z*∈H such that g(z*)=0 and ∥z*∥⩽α.On the solvability of the scheme, we have the following lemma.Lemma 4.2 (Solvability of difference scheme).
Under the assumption(A), for any initial data u0∈XJ, there exists a solution Un∈XJ of  for n≥1.Proof.
The HOC scheme is equivalent to the difference scheme  whose matrix form is (4.1)–(4.3), so we here just prove the unique solvability of the difference scheme (4.1)–(4.3).
In (4.1), for given Un∈XJ(n≥1), we first prove the existence. For j∈𝒯J, rewrite (4.1) as
(4.6)Un+(1/2)=Un+ib4τMδx〈2〉Un+(1/2)+ia4τMδx2Un+(1/2)+iτ2VUn+(1/2)+iq4τ(|Un|2+|2Un+(1/2)-Un|2)Un+(1/2).
We define a mapping ℱ: u*∈XJ→ℱ(u*)∈XJ as
(4.7)ℱ(U*)=U*-Un-ib4τMδx〈2〉U*-ia4τMδx2U*-iτ2VU*-iq4τ(|Un|2+|2U*-Un|2)U*,
which is obviously continue. Computing the inner product of (4.7) with u* and then taking the real part yield
(4.8)ℜ(ℱ(U*),U*)=∥U*∥2-(Un,U*)≥∥U*∥(∥U*∥-∥Un∥),
where Lemma 3.4 was used. Let ∥U*∥=∥Un∥+1, then it follows from (4.8) that ℜ(ℱ(U*),U*)≥0. This together with Lemma 4.1 indecaites that there exists a solution Un+1∈XJ satisfying ℱ(Un+1)=0. Then Un+1 satisfies the scheme .


## 4.2. Discrete Conservation Laws

Corresponding to (1.4)-(1.5), the proposed scheme also possesses two counterpart discrete conservation laws.Lemma 4.3.
The difference solution of scheme conserves the discrete mass and discrete energy in the form of
(4.9)Qn=⋯=Q0,n=0,1,2,…,N,(4.10)En=⋯=E0,n=0,1,2,…,N,
where
(4.11)Qn≔∥Un∥2,En≔b4|||δx〈1〉Un|||2+a4|||δx+Un|||2+12h∑j=0J-1‍Vj|Un|2+q4∥Un∥44.Proof.
Computing the inner product of (4.1) with Un+Un+1 then taking the imaginary part yields
(4.12)1τ(∥Un+1∥2-∥Un∥2)+b4ℑ(Mδx〈2〉(Un+Un+1),(Un+Un+1))+a4ℑ(Mδx2(Un+Un+1),Un+Un+1)+12ℑ(h∑j=0J-1‍Vj(|Un+Un+1|2))+q4ℑ((|Un|2+|Un+1|2)(Un+Un+1),Un+Un+1)=0.
This together with Lemma 3.4 gives (4.9).
Computing the inner product of (4.1) with δt+Un then taking the real part yields
(4.13)ℜ(iδt+Un,δt+Un)+b4ℜ(Mδx〈2〉(Un+Un+1),δt+Un)+a4ℜ(Mδx2(Un+Un+1),δt+Un)+12ℜ(h∑j=0J-1‍Vj(|Un+1|2-|Un|2))+q4ℜ((|Un|2+|Un+1|2)(Un+Un+1),δt+Un)=0.
This together with
(4.14)((|Un|2+|Un+1|2)(Un+Un+1),δt+Un)=(∥Un+1∥44-∥Un∥44)
and Lemma 3.4 gives (4.10).


## 4.3. A Priori Estimate

In proof of the convergence, we need the followinga priori estimate.Lemma 4.4.
The difference solution of the scheme satisfies
(4.15)∥Un∥⩽C,∥δx〈1〉Un∥⩽C,∥δx+Un∥⩽C,∥Un∥∞⩽C.
This implies the stability of the scheme .Proof.
It follows from (4.9) that
(4.16)∥Un∥⩽C.
Lemma 3.1 gives
(4.17)|||δx+Un|||2≥c12∥δx+Un∥2.
Using the discrete Sobolev inequality and (4.16), we obtain
(4.18)∥Un∥44⩽c04(∥δx+Un∥1/4∥Un∥3/4+∥Un∥)4⩽8c04(∥δx+Un∥∥Un∥3+∥Un∥4)⩽ɛ∥δx+Un∥2+16c08ɛ∥Un∥6+8c04∥Un∥4⩽ɛ∥δx+Un∥2+C,
where c0 is a positive constant independent of the grid parameters. Taking ɛ=ac12/2|q|, then it follows from (4.16)–(4.18) and (4.10) that
(4.19)∥δx+Un∥⩽C,∥δx〈1〉Un∥⩽C.
Utilizing the discrete Sobolev estimate and Lemma 4.3, we obtain from (4.16) and (4.19) that
(4.20)∥Un∥∞⩽C.
This completes the proof.


## 4.4. Unconditional Convergence in Maximum Norm

Theorem 4.5.
Suppose thatu(x,t)∈C10,4(R×(0,T]) and τ is sufficiently small, then the solution ujn of the proposed scheme  converges to the solution u(xj,tn) of the problem (1.1)–(1.3) with order O(h8+τ2) in the discrete L∞-norm.Proof.
Denote the error function as(4.21)en=un-Un,n=0,1,2,…,N.
Subtracting (4.1)–(4.3) from , respectively, we obtain the equation of error function as follows:
(4.22)iδt+en+b4Mδx〈2〉(en+en+1)+a4Mδx2(en+en+1)+12Vj(en+en+1)+q4gn+1=ξn,n=0,1,…,N-1,(4.23)ej0=0,j∈𝒯J,(4.24)en∈XJ,n=0,1,2,…,N.
where
(4.25)gjn+1=(|ujn|2+|ujn+1|2)(ujn+ujn+1)-(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1),j∈𝒯J.
Noting that
(4.26)gjn+1=(|ujn|2+|ujn+1|2)(ujn+ujn+1)-(|Ujn|2+|Ujn+1|2)(Ujn+Ujn+1)=(ejnUjn¯+ujnejn¯+ejn+1Ujn+1¯+ujn+1ejn+1¯)(ujn+ujn+1)+(|Ujn|2+|Ujn+1|2)(ejn+ejn+1)=(ejn(ujn¯-ejn¯)+ujnejn¯+ejn+1(ujn+1¯-ejn+1¯)+ujn+1ejn+1¯)×(ujn+ujn+1)+(|ujn-ejn|2+|ujn+1-ejn+1|2)(ejn+ejn+1),j∈𝒯J,
we obtain from (4.25), (2.6), and Lemma 4.4 that
(4.27)||gn+1||2⩽C(||en||2+||en+1||2).
Computing the inner product of (4.22) with en+en+1 and then taking the imaginary part, we obtain
(4.28)1τ(∥en+1∥2-∥en∥2)+q4ℑ(gn+1,en+en+1)=ℑ(ξn,en+en+1),
where Lemma 3.4 was used.
For the last two terms of (4.28), using Lemma 4.4 and (4.26) we have
(4.29)|q4ℑ(gn+1,en+en+1)|≲∥en∥2+∥en+1∥2,|ℑ(ξn,en+en+1)|≲∥ξn∥2+∥en∥2+∥en+1∥2.
Substituting  into (4.28) yields
(4.30)1τ(∥en+1∥2-∥en∥2)≲∥ξn∥2+∥en∥2+∥en+1∥2.
By virtue of Gronwall's inequality and Lemma 2.5, noting that ∥e0∥=0, we obtain from (4.30) that
(4.31)∥en∥⩽C(h8+τ2),n=1,2,…,N,
if τ is sufficiently small.
Computing the inner product of (4.22) with δt+en and then taking the real part, we obtain
(4.32)b4ℜ(Mδx〈2〉(en+en+1),δt+en)+a4ℜ(Mδx2(en+en+1),δt+en)+q4ℜ(gn+1,δt+en)=ℜ(ξn,δt+en),n=1,2,…,N-1.
This together with Lemma 3.4 gives
(4.33)-b4τ(|||δx〈1〉un+1|||2-|||δx〈1〉un|||2)-a4τ(|||δx+un+1|||2-|||δx+un|||2)+q4ℜ(gn+1,δt+en)=ℜ(ξn,δt+en).
It follows from (4.22) that
(4.34)δt+en=ib4Mδx〈2〉(en+en+1)+ia4Mδx2(en+en+1)+i2Vj(en+en+1)+iq4gn+1-iξn,n=1,2,…,N-1.
Thus, for the third term on the left side of (4.33), we have
(4.35)q4ℜ(gn+1,δt+en)=q4ℜ(gn+1,ib4Mδx〈2〉(en+en+1)+ia4Mδx2(en+en+1)+i2Vj(en+en+1)+iq4gn+1-iξn)=-q4ℑ(gn+1,b4Mδx〈2〉(en+en+1))-q4ℑ(gn+1,a4Mδx2(en+en+1))-q4ℑ(gn+1,12V(en+en+1))+q4ℑ(gn+1,ξn).
Noting (4.26), (4.31), and
(4.36)δx+(ejn(ujn¯-ejn¯))=δx+ejn(uj+1n¯-ej+1n¯)+ejn(δx+ujn¯-δx+ejn¯)=δx+ejnUj+1n¯+ejn(δx+ujn¯-δx+ejn¯),δx+((ujn-ejn)(ujn¯-ejn¯))=(δx+ujn-δx+ejn)(uj+1n¯-ej+1n¯)+(ujn-ejn)(δx+ujn¯-δx+ejn¯)=(δx+ujn-δx+ejn)Uj+1n¯+Ujn(δx+ujn¯-δx+ejn¯),|(δx+enen,en)|⩽∥δx+en∥2+∥en∥44⩽2∥δx+en∥2+C(∥en∥6+∥en∥4),
we obtain
(4.37)|q4ℑ(gn+1,a4Mδx2(en+en+1))|=|-q4ℑ(δx+gn+1,a4Mδx+(en+en+1))|≲∥M∥|(δx+gn+1,δx+(en+en+1))|≲|(δx+gn+1,δx+(en+en+1))|(4.38)≲∥δx+en∥2+∥δx+en+1∥2+(en+en+1)2.
Similarly we obtain
(4.39)|q4ℑ(gn+1,b4Mδx〈2〉(en+en+1))|≲∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+(h8+τ2)2.
For the last two terms of (4.35), using the Cauchy inequality, noting that V(x)∈C(Ω), we obtain
(4.40)|-q4ℑ(gn+1,12V(en+en+1))+q4ℑ(gn+1,ξn)|≲∥en∥2+∥en+1∥2+∥ξn∥2≲(h8+τ2)2.
Substituting (4.38)–(4.40) into (4.35) yields
(4.41)|q4ℜ(gn+1,δt+en)|≲∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+∥δx1en∥2+∥δx1en+1∥2+(h8+τ2)2.
Substituting (4.41) into (4.33) yields
(4.42)|||δx〈1〉un+1|||2-|||δx〈1〉un|||2+|||δx+un+1|||2-|||δx+un|||2≲τ[∥δx〈1〉en∥2+∥δx〈1〉en+1∥2+∥δx1en∥2+∥δx1en+1∥2+(h8+τ2)2-ℜ(ξn,δt+en)].
Summing up for the superscript n from 0 to m and then replacing m by n and noting ∥|δx〈1〉u0|∥=∥|δx1u0|∥=0, we obtain
(4.43)|||δx〈1〉un+1|||2+|||δx+un+1|||2≲τ∑l=0n+1‍[∥δx〈1〉el∥2+∥δx1el∥2]+|τ∑l=0n+1‍(ξl,δt+el)|+(h8+τ2)2.
By virtue of Lemmas 3.5, 2.5, and (4.31), we obtain
(4.44)|τ∑l=0n+1‍(ξl,δt+el)|≲(h8+τ2)2.
Substituting (4.44) into (4.43) and then utilizing Gronwall's inequality, we obtain
(4.45)∥δx〈1〉en∥≲h8+τ2,∥δx+en∥≲h8+τ2,n=1,2,…,N,
if τ is sufficiently small. This together with (4.31) gives
(4.46)∥en∥∞≲h8+τ2,n=1,2,…,N.
This completes the proof.


## 5. An Iterative Algorithm

Obviously, the scheme is nonlinear and implicit. In order to obtain the solution Ujn+1 at the time level n+1, an efficient iterative algorithm needs to be proposed. In this section, we construct an iterative algorithm to compute the scheme .For a fixedn, the difference scheme  can be solved by the following iterative algorithm:
(5.1)iτ𝒜h(Ujn+1(s+1)-Ujn)+b4δx〈2〉(Ujn+Ujn+1(s+1))+a4δx2(Ujn+Ujn+1(s+1))+12Ah[Vj(Ujn+Ujn+1(s))]+q4𝒜h[(|Ujn|2+|Ujn+1(s)|2)(Ujn+Ujn+1(s))]=0,j∈𝒯J,0⩽n⩽N-1,Un+1(s+1)∈XJ,0⩽n⩽N,
with
(5.2)Ujn+1(0)={Ujn,ifn=0,2Ujn-Ujn-1,ifn≥1.
Denote
(5.3)Bjn+1(s)=q4[(|Ujn|2+|Ujn+1(s)|2)(Ujn+Ujn+1(s))],λ=τh2,a~=iβ+bλ16,b~=iα+aλ4,c~=i-bλ8-aλ2,d~jn+1(s)=(2iβ-a~-τβ2Vj-2)Uj-2n+(2iβ-a~-τβ2Vj+2)Uj+2n+(2iα-b~-τα2Vj-1)Uj-1n+(2iα-b~-τα2Vj+1)Uj+1n+(2i-c~-τβ2Vj)Ujn-τ[β(Bj-2n+1(s)+Bj+2n+1(s))+α(Bj-1n+1(s)+Bj+1n+1(s))+Bjn+1(s)]+τβ2Vj-2Uj-2n+1(s)+τα2Vj-1Uj-1n+1(s)+τ2VjUjn+1(s)+τα2Vj+1Uj+1n+1(s)+τβ2Vj+2Uj+2n+1(s),
then the iterative algorithm can be written as the following matrix form
(5.4)Hun+1(s+1)=d~n+1(s),
where
(5.5)H=C(c~,b~,a~,0,…,0,a~,b~)
is a cyclic pentadiagonal constant matrix.Though the algebraical systemHun+1(s+1)=d~n+1(s) generated by the proposed scheme is nonlinear and then need iteration in implementation, the coefficient matrix H is cyclic pentadiagonal and invariable at different levels and then the matrix H should not be computed repeatedly. This together with the proposed iterative algorithm and the Thomas algorithm can improve greatly the computational efficiency. To some extent, the proposed scheme is more efficient than most linearized scheme.


## 6. Numerical Experiments

In this section, numerical tests are given to support our theoretical analysis on convergence and stability. In implementation, we chose a toleranceɛ=10-10 to control the iteration; that is, ∥Un+1(s+1)-Un+1(s)∥∞⩽10-10.Example 6.1.
Periodic initial-boundary value problem.
Consider(6.1)i∂u∂t(x,t)+∂2u∂x2(x,t)+2|u(x,t)|2u(x,t)=0,(x,t)∈R×[0,T]u(x,0)=e2ix,x∈R,u(x,t)=u(x+π,t),x∈R.
This problem has the following plane wave solution:
(6.2)u(x,t)=exp[2i(x-t)].Now, we compute Example6.1 on the finite domain [0,π]. In order to test the maximum norm convergence with order O(h8+τ2), we let τ=O(h4) in computation to test the eighth-order convergence in space and take a small h in computation to test second-order convergence in time. In order to quantify the numerical results, we introduce the following “error” functions and “convergence rate:”
(6.3)e(h,τ)=∥UN(h,τ)-uN(h,τ)∥∞,sorder=log(e(h,τ)/e(h/2,τ/24))log(2),torder=log(e(h,τ)/e(h,τ/2))log(2).
Then we list the numerical results of Example 6.1 in Tables 1-2.Table 1
Errors of the numerical solution of Example6.1 computed by the proposed scheme at t=1.
h
τ
Error
Sorder
π
/
5
0.2
1.2112
E
-
002
—
π
/
10
0.2/16
6.7642
E
-
005
8.1724
π
/
20
0.2
/
16
2
4.1327
E
-
007
8.0007
π
/
40
0.2
/
16
3
1.6156
E
-
009
7.9989Table 2
Errors of the numerical solution of Example6.1 computed by the proposed scheme with h=π/40 at t=1.
τ
Error
Torder
0.2
3.0062
E
-
002
—
0.1
6.7267
E
-
003
2.1600
0.05
1.6655
E
-
003
2.0140
0.025
4.1653
E
-
004
1.9995Example 6.2.
Homogenous initial-boundary value problem.
Consider(6.4)i∂u∂t(x,t)+∂2u∂x2(x,t)+2|u(x,t)|2u(x,t)=0,(x,t)∈R×[0,T],u(x,0)=e2ix,x∈R,u(x,t)→0,as|x|→0.
The initial value problem has the one-soliton solution
(6.5)u(x,t)=sech(x-4t)exp(2ix-3it).
Because u(x,t) decays to zero rapidly as |x|→∞, so for sufficiently large number -xL and xR and a finite time T, the initial problem is consistent with the following initial-boundary value problem:
(6.6)i∂u∂t(x,t)+∂2u∂x2(x,t)+2|u(x,t)|2u(x,t)=0,(x,t)∈Ω×(0,T],u(x,t)=0,(x,t)∈∂Ω×[0,T],u(x,0)=e2ixsech(x),x∈Ω¯,
where Ω=(-40,40).In computation of Example6.2, we assume that the value of the exact solution is zero when x∈Ωc:=R∖Ω and list the numerical results of Example 6.2 in Tables 3-4.Table 3
Errors of the numerical solution of Example6.2 computed by the proposed scheme at t=1.
h
τ
Error
Sorder
0.5
0.2
4.3204
E
-
001
—
0.5/2
0.2/16
2.1276
E
-
003
7.6658
0.5
/
2
2
0.2
/
16
2
8.2863
E
-
006
8.0043
0.5
/
2
3
0.2
/
16
3
3.2391
E
-
008
7.9990Table 4
Errors of the numerical solution of Example6.2 computed by the proposed scheme with h=0.15 at t=1.
τ
Error
Torder
0.1
1.5320
E
-
001
—
0.05
3.4611
E
-
002
2.1461
0.025
8.4527
E
-
003
2.0338
0.0125
210386
E
-
004
2.0064Tables1–4 show that the new scheme is eighth-order convergent in space and second order in time for U in the discrete L∞-norm, which verifies Theorem 4.5. Tables 2 and 4 also show that the proposed scheme is very robust.


## 7. Conclusion

In this paper, we construct an eighth-order compact different scheme for numerically solving the NLS equation, use the discrete energy method together with a new technique named “regression of compactness” to prove that the proposed HOC scheme is convergent, without any restriction on the grid ratio, at the order ofO(h8+τ2) in the discrete L∞-norm. What's better, the HOC scheme conserves the discrete total mass and energy, which changed the viewpoint that the compact difference scheme of a nonlinear equation cannot preserve the discrete conservation laws. In order to implement the nonlinear compact difference scheme, we propose an iterative algorithm which is numerically proved to be very reliable.



---

*Source: 913429-2012-10-17.xml*