Advances in Numerical Analysis

Volume 2012, Article ID 913429, 24 pages

http://dx.doi.org/10.1155/2012/913429

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

School of Mathematics and Statistics, Nanjing University of Information Science and Technology, Nanjing 210044, China

Received 24 May 2012; Accepted 7 August 2012

Academic Editor: Ting-Zhu Huang

Copyright © 2012 Tingchun Wang. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### 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 of in the discrete -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.

#### 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 norm convergence of any compact difference scheme for nonlinear equations.

In this paper, we consider the following cubic nonlinear Schrödinger (NLS) equation: subject to periodic boundary condition and initial condition where is dimensionless constant characterizing the interaction (positive for repulsive interaction and negative for the attractive interaction) between particles in BEC. 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. is the unknown periodical complex-valued function whose initial value is a given periodical complex-valued function. is the period of the .

Denote ; it is easy to show that the problem (1.1)–(1.3) conserves the total *mass*
and the global *energy*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 total *mass* 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 Section 2, 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 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 integer , choose time-step and denote time steps , where with the maximal existing time of the solution; choose mesh size , with a positive integer and denote grid points as Denote the index sets Denote to be the numerical approximation and, respectively, be the exact value of at the point for and , and denote to be the numerical solution and, respectively, be the exact solution at time . For a grid function , introduce the following finite difference operators:

We denote the space and define seminorms and discrete inner product over as

For simplicity, we write as . Throughout the paper, we adopt the standard Sobolev spaces and their corresponding norms let denote a generic constant which is independent of mesh size and time step , and use the notation to present that there is a generic constant which is independent of time step and mesh size such that . For the exact solution of the initial-boundary value problem (1.1)–(1.3), we assume that

##### 2.1. Some Properties of Cyclic Matrix

A matrix in the form of is called a cyclic matrix. Because the matrix is determined by the entries in the first row, the matrix can be denoted as Denote , which is called as a basic cyclic matrix. By a simple calculation, we can get that Then, any cyclic matrix can be written as a polynomial of ; that is, This implies that, for any two cyclic matrices and , . 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
by using the multiplication of polynomials and the properties we obtain
where

Lemma 2.2. *If a cyclic matrix is invertible, then its inverse matrix is still a cyclic matrix.*

*Proof. *If is invertible, we assume that the inverse matrix is , then
Multiplying with (2.14) yields
If , we obtain the linear system . Noting that , we obtain that the system has a unique nonzero solution . This together with (2.15) gives that where is the unique nonzero solution of . It follows that the assumption is right; that is, is a cyclic matrix.

Lemma 2.3. *For any a real cyclic matrix , all the eigenvalues are in the form of
**
where and .*

*Proof. *On the one hand, it follows from that which indicates that the eigenvalues of are

On the other hand, just as (2.10), the real cyclic matrix can be written as
where are real numbers. This gives that the eigenvalues of are in the form of

##### 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. Suppose is a -periodic function, then for the discretization of the second-order derivatives , we have the formula [31]: where , , , , , and are undetermined parameters which depend on the accuracy order constraints. After tedious calculation, we can get these constraints: According to (2.20), an eighth-order scheme is obtained with which is the unique solution of the constraint equation .

Thus, the approximation (2.20) to can be written in the form of finite-dimensional linear operator where the finite difference operator is defined as The corresponding matrix form is where 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)* is invertible and its inverse matrix also has real entries and is symmetric and positive definite;*(b)* can be decomposed as
**where also has real entries and is symmetric and positive definite.*

##### 2.3. High-Order Compact Scheme

Imposing the approximations given in Section 2.2 on the spatial direction and the Crank-Nicolson discretization on the temporal direction of the NLS equation (1.1) gives where and is the local truncation error. The corresponding matrix form of the (2.28)–(2.30) is where and .

On the local truncation error , using Taylor's expansion, we obtain the follwing lemma.

Lemma 2.5. *Suppose that and , there is
*

Omitting the small term in (2.28) yields the following difference scheme that is,

*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

Lemma 2.3 gives that the eigenvalues of the cyclic matrix are in the form of This gives It is easy to verify that , which indicates that the matrix is positive definite.

For the real, positive definite, symmetric and cyclic matrix , it follows from Lemma 2.2 that there are real numbers such that . Lemmas 2.3 and 2.4 give that is real positive definite and symmetric. And then we can introduce the following discrete norm: On the relation between and , we have the following lemma.

Lemma 3.1. *The discrete norms and are equivalent; in fact, there is
**
where , .*

*Proof. *It follows from (3.2) that the eigenvalues of satisfy
This gives the spectral radius ), and consequently
Because
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 function , there are
*

*Proof. *For the entry of the vector , there is
where the periodicity of is used. This gives
For the element of the vector , there is
Then (3.8) is gotten from (3.11) and (3.12). Similarly, (3.9) can be proved.

Lemma 3.3. *For any two grid functions ; there are
*

*Proof. *Direct calculation and noticing the periodicity of the two grid functions yield .

Lemma 3.4. *For the approximation , there are
**
where “” and “” mean taking the imaginary part and the real part, respectively.*

*Proof. *Lemmas 3.2, 3.3, and definition (3.3) together with the symmetry of give
Similarly, (3.15) can be proved. Also using Lemmas 3.2, and 3.3, definition (3.3), and the symmetry of , we obtain
Similarly, (3.17) can be proved.

Lemma 3.5 (see [33]). *For time sequences and , there is
*

#### 4. Numerical Analysis

The corresponding matrix form of is where

##### 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 be a finite-dimensional inner-product space, the associated norm, and continuous. Assume, moreover, that
**
Then, there exists a such that and .*

On the solvability of the scheme , we have the following lemma.

Lemma 4.2 (Solvability of difference scheme). *Under the assumption , for any initial data , there exists a solution of for .*

* 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 , we first prove the existence. For , rewrite (4.1) as
We define a mapping : as
which is obviously continue. Computing the inner product of (4.7) with and then taking the real part yield
where Lemma 3.4 was used. Let , then it follows from (4.8) that . This together with Lemma 4.1 indecaites that there exists a solution satisfying . Then 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
**
where
*

*Proof. *Computing the inner product of (4.1) with then taking the imaginary part yields
This together with Lemma 3.4 gives (4.9).

Computing the inner product of (4.1) with then taking the real part yields
This together with
and Lemma 3.4 gives (4.10).

##### 4.3. A Priori Estimate

In proof of the convergence, we need the following * a priori* estimate.

Lemma 4.4. *The difference solution of the scheme satisfies
**
This implies the stability of the scheme .*

*Proof. *It follows from (4.9) that
Lemma 3.1 gives
Using the discrete Sobolev inequality and (4.16), we obtain
where is a positive constant independent of the grid parameters. Taking , then it follows from (4.16)–(4.18) and (4.10) that
Utilizing the discrete Sobolev estimate and Lemma 4.3, we obtain from (4.16) and (4.19) that
This completes the proof.

##### 4.4. Unconditional Convergence in Maximum Norm

Theorem 4.5. * Suppose that and is sufficiently small, then the solution of the proposed scheme converges to the solution of the problem (1.1)–(1.3) with order in the discrete -norm.*

*Proof. *Denote the error function as
Subtracting (4.1)–(4.3) from , respectively, we obtain the equation of error function as follows:
where
Noting that
we obtain from (4.25), (2.6), and Lemma 4.4 that
Computing the inner product of (4.22) with and then taking the imaginary part, we obtain
where Lemma 3.4 was used.

For the last two terms of (4.28), using Lemma 4.4 and (4.26) we have
Substituting into (4.28) yields
By virtue of Gronwall's inequality and Lemma 2.5, noting that , we obtain from (4.30) that
if is sufficiently small.

Computing the inner product of (4.22) with and then taking the real part, we obtain
This together with Lemma 3.4 gives
It follows from (4.22) that
Thus, for the third term on the left side of (4.33), we have
Noting (4.26), (4.31), and
we obtain
Similarly we obtain
For the last two terms of (4.35), using the Cauchy inequality, noting that , we obtain
Substituting (4.38)–(4.40) into (4.35) yields
Substituting (4.41) into (4.33) yields
Summing up for the superscript from to and then replacing by and noting , we obtain
By virtue of Lemmas 3.5, 2.5, and (4.31), we obtain
Substituting (4.44) into (4.43) and then utilizing Gronwall's inequality, we obtain
if is sufficiently small. This together with (4.31) gives
This completes the proof.

#### 5. An Iterative Algorithm

Obviously, the scheme is nonlinear and implicit. In order to obtain the solution at the time level , an efficient iterative algorithm needs to be proposed. In this section, we construct an iterative algorithm to compute the scheme .

For a fixed , the difference scheme can be solved by the following iterative algorithm: with Denote then the iterative algorithm can be written as the following matrix form where is a cyclic pentadiagonal constant matrix.

Though the algebraical system generated by the proposed scheme is nonlinear and then need iteration in implementation, the coefficient matrix is cyclic pentadiagonal and invariable at different levels and then the matrix 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 to control the iteration; that is, .

*Example 6.1. * Periodic initial-boundary value problem.

Consider
This problem has the following plane wave solution:

Now, we compute Example 6.1 on the finite domain . In order to test the maximum norm convergence with order , we let in computation to test the eighth-order convergence in space and take a small 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:” Then we list the numerical results of Example 6.1 in Tables 1-2.

*Example 6.2. * Homogenous initial-boundary value problem.

Consider
The initial value problem has the one-soliton solution
Because decays to zero rapidly as , so for sufficiently large number and and a finite time , the initial problem is consistent with the following initial-boundary value problem:
where .

In computation of Example 6.2, we assume that the value of the exact solution is zero when and list the numerical results of Example 6.2 in Tables 3-4.

Tables 1–4 show that the new scheme is eighth-order convergent in space and second order in time for in the discrete -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 of in the discrete -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.

#### Acknowledgment

This work is supported by the National Natural Science Foundation of China, nos. 11126292 and 11201239.

#### References

- D. J. Griffiths,
*Introduction to Quantum Mechanics*, Prentice-Hall, Englewood Cliffs, NJ, USA, 1995. - C. R. Menyuk, “Stability of solitons in birefringent optical fibers,”
*Journal of the Optical Society of America B*, vol. 5, pp. 392–402, 1998. View at Google Scholar - M. Wadati, T. Iizuka, and M. Hisakado, “A coupled nonlinear Schrödinger equation and optical solitons,”
*Journal of the Physical Society of Japan*, vol. 61, no. 7, pp. 2241–2245, 1992. View at Publisher · View at Google Scholar - G. D. Akrivis, “Finite difference discretization of the cubic Schrödinger equation,”
*IMA Journal of Numerical Analysis*, vol. 13, no. 1, pp. 115–124, 1993. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - T. F. Chan and L. J. Shen, “Stability analysis of difference schemes for variable coefficient Schrödinger type equations,”
*SIAM Journal on Numerical Analysis*, vol. 24, no. 2, pp. 336–349, 1987. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - Q. Chang, E. Jia, and W. Sun, “Difference schemes for solving the generalized nonlinear Schrödinger equation,”
*Journal of Computational Physics*, vol. 148, no. 2, pp. 397–415, 1999. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - W. Z. Dai, “An unconditionally stable three-level explicit difference scheme for the Schrödinger equation with a variable coefficient,”
*SIAM Journal on Numerical Analysis*, vol. 29, no. 1, pp. 174–181, 1992. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - M. Dehghan and A. Taleei, “A compact split-step finite difference method for solving the nonlinear Schrödinger equations with constant and variable coefficients,”
*Computer Physics Communications*, vol. 181, no. 1, pp. 43–51, 2010. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - F. Ivanauskas and M. Radžiūnas, “On convergence and stability of the explicit difference method for solution of nonlinear Schrödinger equations,”
*SIAM Journal on Numerical Analysis*, vol. 36, no. 5, pp. 1466–1481, 1999. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - P. L. Nash and L. Y. Chen, “Efficient finite difference solutions to the time-dependent Schrödinger equation,”
*Journal of Computational Physics*, vol. 130, no. 2, pp. 266–268, 1997. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - J. M. Sanz-Serna, “Methods for the numerical solution of the nonlinear Schroedinger equation,”
*Mathematics of Computation*, vol. 43, no. 167, pp. 21–27, 1984. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - Z. Z. Sun and X. Wu, “The stability and convergence of a difference scheme for the Schrödinger equation on an infinite domain by using artificial boundary conditions,”
*Journal of Computational Physics*, vol. 214, no. 1, pp. 209–223, 2006. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - L. Wu, “Dufort-Frankel-type methods for linear and nonlinear Schrödinger equations,”
*SIAM Journal on Numerical Analysis*, vol. 33, no. 4, pp. 1526–1533, 1996. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - J. M. Sanz-Serna and J. G. Verwer, “Conservative and nonconservative schemes for the solution of the nonlinear Schrödinger equation,”
*IMA Journal of Numerical Analysis*, vol. 6, no. 1, pp. 25–42, 1986. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - Z. Fei, V. M. Pérez-García, and L. Vázquez, “Numerical simulation of nonlinear Schrödinger systems: a new conservative scheme,”
*Applied Mathematics and Computation*, vol. 71, no. 2-3, pp. 165–177, 1995. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - S. Li and L. Vu-Quoc, “Finite difference calculus invariant structure of a class of algorithms for the nonlinear Klein-Gordon equation,”
*SIAM Journal on Numerical Analysis*, vol. 32, no. 6, pp. 1839–1875, 1995. View at Publisher · View at Google Scholar · View at Zentralblatt MATH - G. D. Akrivis, V. A. Dougalis, and O. A. Karakashian, “On fully discrete Galerkin methods of second-order temporal accuracy for the nonlinear Schrödinger equation,”
*Numerische Mathematik*, vol. 59, no. 1, pp. 31–53, 1991. View at Publisher ·