Discrete Dynamics in Nature and Society
Volume 2009 (2009), Article ID 283959, 27 pages
doi:10.1155/2009/283959
Research Article

Large Time-Stepping Spectral Methods for the Semiclassical Limit of the Defocusing Nonlinear Schrödinger Equation

Zongqi Liang

School of Science, Jimei University, Xiamen 361021, China

Received 28 August 2009; Accepted 20 October 2009

Academic Editor: Guang Zhang

Copyright © 2009 Zongqi Liang. 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

We analyze a class of large time-stepping Fourier spectral methods for the semiclassical limit of the defocusing Nonlinear Schrödinger equation and provide highly stable methods which allow much larger time step than for a standard implicit-explicit approach. An extra term, which is consistent with the order of the time discretization, is added to stabilize the numerical schemes. Meanwhile, the first-order and second-order semi-implicit schemes are constructed and analyzed. Finally the numerical experiments are performed to demonstrate the effectiveness of the large time-stepping approaches.

1. Introduction

Consider the initial-value problem for the semiclassical limit of defocusing nonlinear Schrödinger equation (DFNLS) in one space dimension 𝑖 𝜀 𝜕 𝑢 + 𝜀 𝜕 𝑡 2 2 𝑢 𝑥 𝑥 − 𝜆 | 𝑢 | 2 𝑢 𝑢 = 0 ( 𝑥 ∈ Ω , 𝑡 > 0 ) , 0 ( 𝑥 ) = 𝐴 0 ( 𝑥 ) 𝑒 𝑖 𝑆 0 ( 𝑥 ) / 𝜀 , ( 1 . 1 ) where 𝑖 2 = − 1 , the domain Ω is an open set in ℝ , 𝑢 ( 𝑥 , 𝑡 ) is the complex-valued wave function, 0 < 𝜀 ≪ 1 is the scaled Planck constant, and 𝜆 ≥ 0 is a real constant. 𝐴 0 ( 𝑥 ) is given with the initial amplitude and 𝑆 0 ( 𝑥 ) is the real initial phase function with 𝐴 0 ( 𝑥 ) and 𝑆 0 ( 𝑥 ) decaying rapidly as | 𝑥 | → ∞ . To avoid boundary effects and obtain enough information of statistics, we compute the solution in a sufficiently large interval [ − 𝑀 , 𝑀 ] .

The DFNLS (1.1) has infinitely many conservation laws (see, e.g., [1–3]), two of which are

𝜕 𝑡 𝜌 + 𝜕 𝑥 𝜕 𝐽 = 0 , 𝑡 𝐽 + 𝜕 𝑥 𝐽 2 𝜌 + 1 2 𝜕 𝑥 𝜌 2 = 𝜀 2 4 𝜕 𝑥 𝜌 𝜕 2 𝑥 , l o g 𝜌 ( 1 . 2 ) where 𝜌 (the position density) and 𝐽 (the current density) are primary physical quantities and can be computed from the wave function 𝑢 ( 𝑥 , 𝑡 ) :

| | | | 𝜌 ( 𝑥 , 𝑡 ) = 𝑢 ( 𝑥 , 𝑡 ) 2 , 𝐽 ( 𝑥 , 𝑡 ) = 𝜀 𝐼 𝑚 𝑢 ( 𝑥 , 𝑡 ) 𝑢 𝑥 , ( 𝑥 , 𝑡 ) ( 1 . 3 ) where “ − “ denotes complex conjugation. The two conservation laws can also be denoted by the 𝐿 2 norm and energy density conservation laws

( ‖ 𝑢 ⋅ , 𝑡 ) ‖ 2 = 𝑅 | | | | 𝑢 ( 𝑥 , 𝑡 ) 2 ‖ ‖ 𝑢 𝑑 𝑥 = 0 ‖ ‖ ( 𝑥 ) 2 , 𝜀 𝐸 ( ⋅ , 𝑡 ) = 2 2 𝑅 | | 𝑢 𝑥 | | ( 𝑥 , 𝑡 ) 2 𝜆 𝑑 𝑥 + 2 𝑅 | | | | 𝑢 ( 𝑥 , 𝑡 ) 4 𝑑 𝑥 , ( 1 . 4 ) where ‖ ⋅ ‖ is the standard 𝐿 2 -norm in ℝ .

The DFNLS is connected with many applications in science and technology. One of the most important applications of the semiclassical limit of (1.1) arises in nonlinear optics [4], in particular regarding the production and propagation of laser pulses in fiber optics. Here 𝜀 depends on the dimensions on the problem (amplitude of the initial pulse, dimensions of the fiber, etc.). Another application appears in the study of the so-called Bose-Einstein condensate [5, 6], where 𝜀 now does represent the Planck constant.

The limit 𝜀 → 0 is called the semiclassical limit and considerable attention has been given recently to the investigation of its existence and structure [7–11]. The dynamics of the limit is an open problem. Wright et al. have provided the exact solution of the geometric optics approximation of the defocusing nonlinear Schrödinger equation in [12]. Much progress has been made recently in understanding semiclassical limits of the linear Schrödinger equation. Particularly by the introduction of tools from microlocal analysis, such as defect measures [13], H-measures [14], and Wigner measures [15, 16].

Numerically this is a notoriously difficult problem. The oscillatory nature of the solutions of the nonlinear Schrödinger equation with small 𝜀 provides severe numerical burdens. Even for stable numerical approximations (or under mesh size and time step restrictions which guarantee stability), the oscillations may very well pollute the solution in such a way that the quadratic macroscopic quantities and other physical observables come out completely wrong unless the spatial-temporal oscillations are fully resolved numerically. For a long-time considerable consideration has been given to this problem by many researchers. In [17, 18], Markowich et al. utilized the Wigner measure in analyzing the semiclassical limit for the the linear Schrödinger equation with small 𝜀 , to study the finite difference approximation. In [19, 20], Bao et al. simulated the solution by the time-splitting spectral approximation for the linear and nonlinear Schrödinger equation. Hector D. Ceniceros presented the semiclassical limit of the focusing nonlinear Schrödinger equation by moving mesh method [21] and Ceniceros and Tian studied it by nonlinear Fourier filtering method in [22]. But there are few research on DFNLS.

The purpose of this paper is to provide an effective numerical method to solve the DFNLS with small parameter 𝜀 . In space direction, we use Fourier spectral method, and in time direction we use difference method. In the case that the space direction has spectral accuracy, we expect to improve the computational efficiency by increasing the time step, which is beneficial to the numerical simulation and the discussion on the long-time behavior of this problem. In numerical computation, for one thing, because of the strong nonlinear character of (1.1), the numerical stability cannot be completely ensured by the traditional numerical methods. For another thing, because of the oscillation of high frequency of the solution of this equation, an ideal discretization of space and time should decompose bordering surface precisely. Hence a high resolution ratio of space is needed, and very small time step must be used to satisfy the demands of computational stability if we take standard explicit or semi-implicit scheme. When 𝜀 is very small, the explicit treatment for nonlinear item 𝜆 | 𝑢 2 | 𝑢 usually causes strict restriction on time step, and it is not practical to have a long-time computation in this case. To relax restriction of this kind of problem, we provide a semi-implicit difference Fourier spectral method. That is to say, to improve the stability of algorithm, we add a routine item to the semi-implicit scheme where a larger time step is allowed in the numerical experiment [23–25]. First, for the treatment of time discretization, the semi-implicit scheme is used. Second, for the second-order item in (1.1), the implicit scheme is employed to reduce the restriction on time step. Third, for the nonlinear item, the explicit scheme is used to avoid solving systems of nonlinear equations at each time step as this explicit treatment can be computed by FFT. We expect that the computation will be easy, large time step is allowed in the solution procedure, and the numerical stability can be guaranteed at the same time.

The organization of this paper is as follows. In Section 2, a theoretical analysis for the semi-implicit method is provided, and the first-order time-stepping method stability properties are investigated. The second-order semi-implicit methods are mainly investigated in Section 3. Full discretization schemes and accuracy tests are presented in Section 4. Finally, some computational results are given in Section 5.

2. Stability Analysis for First-Order Time Discretization for (1.1)

In this section, we present Fourier spectral approximations of the problem (1.1) with periodic boundary conditions.

We define periodic Sobolev space as

𝐻 𝑠 𝑣 ( Ω ) = ( 𝑥 , 𝑡 ) ∈ 𝐻 𝑠 l o c ( Ω ) ∣ 𝑣 ( 𝑥 + 𝐿 , 𝑡 ) = 𝑣 ( 𝑥 , 𝑡 ) . ( 2 . 1 )

Define inner product and norm as

( 𝑢 , 𝑣 ) = Ω 𝑢 𝑣 𝑑 𝑥 , ‖ 𝑣 ‖ 2 = ( 𝑣 , 𝑣 ) . ( 2 . 2 )

For any 1 ≤ 𝑝 ≤ ∞ , let 𝐿 𝑝 ( Ω ) = { 𝑣 ∣ ‖ 𝑣 ‖ 𝐿 𝑝 < ∞ } , where

‖ 𝑣 ‖ 𝐿 𝑝 = ⎧ ⎪ ⎨ ⎪ ⎩ Ω | 𝑣 | 𝑞 𝑑 𝑥 1 / 𝑝 , 1 ≤ 𝑝 < ∞ , e s s s u p 𝑥 ∈ Ω | | | | 𝑣 ( 𝑥 ) , 𝑝 = ∞ . ( 2 . 3 )

First we consider the classical first-order semi-implicit method, the scheme is obtained by combining first-order backward difference (BD1) for 𝜕 𝑢 / 𝜕 𝑡 and a first-order extrapolation (EP1) for the nonlinear item at the same time 𝑢 𝑖 𝜀 𝑛 + 1 − 𝑢 𝑛 + 𝜀 ▵ 𝑡 2 2 𝑢 𝑛 + 1 𝑥 𝑥 | | 𝑢 − 𝜆 𝑛 | | 2 𝑢 𝑛 = 0 ( 𝑛 > 0 ) , ( 2 . 4 ) where ▵ 𝑡 is the time-step and 𝑡 𝑛 = 𝑛 ▵ 𝑡 , 𝑢 𝑛 is an approximation to 𝑢 ( 𝑥 , 𝑡 𝑛 ) .

We expect that the implicit treatment for the second-order item 𝑢 𝑥 𝑥 in (2.4) can relax the restriction on time step. And we notice that it will not increase the complexity of the algorithm through the spectral method because the second-order item here is complete linear. In Section 4, it is showed that large step is not allowed in scheme (2.4) when 𝜀 is small. That is to say, to guarantee the stability, a very small time step must be used to get the precise solution. To overcome this shortage, we add an item to the scheme and have the following modified scheme (BD1/EP1):

𝑢 𝑖 𝜀 𝑛 + 1 − 𝑢 𝑛 + 𝜀 ▵ 𝑡 2 2 𝑢 𝑛 + 1 𝑥 𝑥 𝑢 − 𝐴 𝑛 + 1 − 𝑢 𝑛 | | 𝑢 − 𝜆 𝑛 | | 2 𝑢 𝑛 = 0 ( 𝑛 > 0 ) , ( 2 . 5 ) where 𝐴 is an undetermined positive constant, which will be given in the following results. The purpose of adding an item is to improve the stability of the original scheme, and a larger time step can be allowed in the computation.

As we can see, (2.5) can be denoted as the following weak form:

𝑢 𝑖 𝜀 𝑛 + 1 − 𝑢 𝑛 − 𝜀 ▵ 𝑡 , 𝑣 2 2 𝑢 𝑥 𝑛 + 1 , 𝑣 𝑥 𝑢 − 𝐴 𝑛 + 1 − 𝑢 𝑛 | | 𝑢 , 𝑣 − 𝜆 𝑛 | | 2 𝑢 𝑛 , 𝑣 = 0 , ∀ 𝑣 ∈ 𝐻 1 ( Ω ) . ( 2 . 6 )

Theorem 2.1. If 𝐴 in (2.5) is sufficiently large, the following energy inequality holds: 𝐸 𝑢 𝑛 + 1 ≤ 𝐸 ( 𝑢 𝑛 ) , ( 2 . 7 ) where 𝐸 ( 𝑢 𝑛 𝜀 ) = 2 2 ‖ ‖ 𝑢 𝑛 𝑥 ‖ ‖ 2 + 𝜆 2 ‖ 𝑢 𝑛 ‖ 4 𝐿 4 . ( 2 . 8 ) More precisely, the positive constant 𝐴 satisfies 𝜆 𝐴 ≥ 2 | | 𝑢 𝑛 | | 2 + 𝜆 4 | | 𝑢 𝑛 + 1 + 𝑢 𝑛 | | 2 . ( 2 . 9 )

Proof. Let 𝑣 = 2 ( 𝑢 𝑛 + 1 − 𝑢 𝑛 ) in (2.6), take the real part of the result 𝜀 2 𝑢 R e 𝑥 𝑛 + 1 , 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 𝐴 𝑢 + 2 R e 𝑛 + 1 − 𝑢 𝑛 | | 𝑢 + 𝜆 𝑛 | | 2 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 = 0 . ( 2 . 1 0 )
Now we estimate the terms in (2.10), respectively,
𝑢 2 R e 𝐴 𝑛 + 1 − 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 𝑢 = 2 𝐴 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 , 𝜀 2 𝑢 R e 𝑥 𝑛 + 1 , 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 = 𝜀 2 𝑢 R e 𝑥 𝑛 + 1 + 𝑢 𝑛 𝑥 2 + 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 2 , 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 = 𝜀 2 2 ‖ ‖ 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 ‖ ‖ 2 + 𝜀 2 2 ‖ ‖ 𝑢 𝑥 𝑛 + 1 ‖ ‖ 2 − ‖ ‖ 𝑢 𝑛 𝑥 ‖ ‖ 2 , | | 𝑢 2 𝜆 R e 𝑛 | | 2 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 𝑢 = 2 𝜆 R e 𝑛 | | 2 , 𝑢 𝑛 + 1 + 𝑢 𝑛 2 − 𝑢 𝑛 + 1 − 𝑢 𝑛 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 𝑢 = 𝜆 𝑛 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 − | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 | | 𝑢 = 𝜆 𝑛 + 1 | | 2 + | | 𝑢 𝑛 | | 2 2 − | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 − 𝜆 Ω | | 𝑢 𝑛 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 = 𝜆 𝑑 𝑥 2 ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 4 𝐿 4 − ‖ 𝑢 𝑛 ‖ 4 𝐿 4 − 𝜆 2 Ω | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 𝑑 𝑥 − 𝜆 Ω | | 𝑢 𝑛 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 𝑑 𝑥 . ( 2 . 1 1 )
From the inequality
Ω | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 𝑑 𝑥 ≤ Ω | | 𝑢 𝑛 + 1 + 𝑢 𝑛 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 𝑑 𝑥 ( 2 . 1 2 ) and the above estimations, we obtain 𝐸 𝑢 𝑛 + 1 + 𝜀 2 2 ‖ ‖ 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 ‖ ‖ 2 + Ω 𝜆 2 𝐴 − 2 | | 𝑢 𝑛 + 1 + 𝑢 𝑛 | | 2 | | 𝑢 − 𝜆 𝑛 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 𝑑 𝑥 ≤ 𝐸 ( 𝑢 𝑛 ) . ( 2 . 1 3 )
The last term mentioned previously can be made nonnegative if the following inequality holds:
𝜆 𝐴 ≥ 2 | | 𝑢 𝑛 | | 2 + 𝜆 4 | | 𝑢 𝑛 + 1 + 𝑢 𝑛 | | 2 . ( 2 . 1 4 )

Remark 2.2. From the proof of Theorem 2.1, we know that the selection of 𝐴 in (2.14) depends on solution of the equation in the ( 𝑛 + 1 ) th step. But essentially it is not difficult to choose a suitable 𝐴 . Actually, we can determine it approximately from the following expression: 𝐴 ≥ m a x 𝑥 ∈ Ω 3 𝜆 2 | | 𝑢 𝑛 | | ( 𝑥 ) 2 . ( 2 . 1 5 )

Remark 2.3. In numerical computation, condition (2.14) can be guaranteed by a posterior method; that is, when computation of each step finishes, we check whether the condition (2.14) is satisfied to guarantee the decrement or the conservation of discrete energy.

3. Higher-Order Methods for Semi-Implicit Time Discretization

3.1. Second-Order Scheme: BD2/EP2

Similar to the first-order scheme, Higher-Order time discretization can be constructed. For example, by combining the second-order backward difference (BD2) for 𝜕 𝑢 / 𝜕 𝑡 , and the second-order extrapolation (EP2) for the explicit treatment of nonlinear item, we have the following second-order scheme:

𝑖 𝜀 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 + 𝜀 2 Δ 𝑡 2 2 𝑢 𝑛 + 1 𝑥 𝑥 | | − 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = 0 , ( 3 . 1 ) here we add an 𝐴 -item similar to the item added in the first-order scheme, and expect that this item can be consistent with the global order of the scheme and keep the scheme stable. Therefore an 𝑂 ( ▵ 𝑡 2 𝜕 𝑡 𝑡 𝑢 ) item is added in the above scheme, and we arrive at the following modified second-order scheme (BD2/EP2):

𝑖 𝜀 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 + 𝜀 2 Δ 𝑡 2 2 𝑢 𝑛 + 1 𝑥 𝑥 𝑢 − 𝐴 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | − 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = 0 . ( 3 . 2 ) As usual, to start the iteration, 𝑢 0 ( 𝑥 ) is given by the initial condition and 𝑢 1 ( 𝑥 ) is computed by the first-order scheme (2.5).

Equation (3.2) can also be denoted by the following weak form: 𝑖 𝜀 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 − 𝜀 2 Δ 𝑡 , 𝑣 2 2 𝑢 𝑥 𝑛 + 1 , 𝑣 𝑥 𝑢 − 𝐴 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | , 𝑣 − 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 − 𝑢 𝑛 − 1 , 𝑣 = 0 , ∀ 𝑣 ∈ 𝐻 1 ( Ω ) , ( 3 . 3 ) where 𝐴 is a positive constant.

Lemma 3.1. One has R e 2 𝑢 𝑛 + 1 , 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 = ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 2 − ‖ 𝑢 𝑛 ‖ 2 + ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 + ‖ ‖ 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 ‖ ‖ 2 . ( 3 . 4 )

Proof. We have R e 2 𝑢 𝑛 + 1 , 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 = R e 2 𝑢 𝑛 + 1 , 𝑢 𝑛 + 1 − 𝑢 𝑛 + 2 𝑢 𝑛 + 1 − 𝑢 𝑛 − 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = R e 2 𝑢 𝑛 + 1 , 𝑢 𝑛 + 1 − 𝑢 𝑛 + R e 2 𝑢 𝑛 + 1 , 2 𝑢 𝑛 + 1 − 𝑢 𝑛 − 2 𝑢 𝑛 − 𝑢 𝑛 − 1 , R e 2 𝑢 𝑛 + 1 , 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 = R e 𝑛 + 1 + 𝑢 𝑛 + 𝑢 𝑛 + 1 − 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 = ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 2 − ‖ 𝑢 𝑛 ‖ 2 𝑢 + R e 𝑛 + 1 − 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 , R e 2 𝑢 𝑛 + 1 , 2 𝑢 𝑛 + 1 − 𝑢 𝑛 − 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = R e 2 𝑢 𝑛 + 1 − 𝑢 𝑛 + 2 𝑢 𝑛 − 𝑢 𝑛 − 1 − 𝑢 𝑛 + 𝑢 𝑛 − 1 , 2 𝑢 𝑛 + 1 − 𝑢 𝑛 − 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 + R e − 𝑢 𝑛 + 𝑢 𝑛 − 1 , 2 𝑢 𝑛 + 1 − 3 𝑢 𝑛 + 𝑢 𝑛 − 1 = ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 + R e − 𝑢 𝑛 + 𝑢 𝑛 − 1 𝑢 , 2 𝑛 + 1 − 𝑢 𝑛 − 𝑢 𝑛 − 𝑢 𝑛 − 1 = ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 𝑢 − 2 R e 𝑛 − 𝑢 𝑛 − 1 , 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 + R e 𝑛 − 𝑢 𝑛 − 1 , 𝑢 𝑛 − 𝑢 𝑛 − 1 , ( 3 . 5 ) then R e 2 𝑢 𝑛 + 1 , 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 = ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 2 − ‖ 𝑢 𝑛 ‖ 2 + ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 𝑢 + R e 𝑛 + 1 − 𝑢 𝑛 , 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 − 2 R e 𝑛 − 𝑢 𝑛 − 1 , 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 + R e 𝑛 − 𝑢 𝑛 − 1 , 𝑢 𝑛 − 𝑢 𝑛 − 1 = ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 2 − ‖ 𝑢 𝑛 ‖ 2 + ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 + ‖ ‖ 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 ‖ ‖ 2 . ( 3 . 6 )

Lemma 3.2. One has 𝑢 R e 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 , 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 = ‖ ‖ 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 ‖ ‖ 𝑢 + 2 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 ‖ ‖ 2 . ( 3 . 7 )

Theorem 3.3. Let { 𝑢 𝑛 ( 𝑥 ) , 𝑛 = 1 , 2 , … } be the solution of the semidiscrete problem (3.2). The corresponding generalized energy 𝐸 𝑛 is defined as follows: 𝐸 𝑛 = 𝜀 2 2 ‖ ‖ 𝑢 𝑛 𝑥 ‖ ‖ 2 + ‖ ‖ 2 𝑢 𝑛 𝑥 − 𝑢 𝑥 𝑛 − 1 ‖ ‖ 2 ‖ ‖ 𝑢 + 2 𝐴 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 + 𝜆 2 ‖ 𝑢 𝑛 ‖ 4 𝐿 4 + ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 4 𝐿 4 . ( 3 . 8 ) Then, for large enough 𝐴 , the following discrete energy inequality holds 𝐸 𝑛 + 1 ≤ 𝐸 𝑛 . ( 3 . 9 )

Proof. Let 𝑣 = − 2 ( 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 ) in (3.3), taking the real part of the result we have 𝜀 2 2 R e 3 𝑢 𝑥 𝑛 + 1 − 4 𝑢 𝑛 𝑥 + 𝑢 𝑥 𝑛 − 1 , 2 𝑢 𝑥 𝑛 + 1 + 2 𝐴 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 + 2 𝜆 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = 0 . ( 3 . 1 0 ) Using Lemmas 3.1 and 3.2, the first term of the left of (3.10) is 𝜀 2 2 R e 3 𝑢 𝑥 𝑛 + 1 − 4 𝑢 𝑛 𝑥 + 𝑢 𝑥 𝑛 − 1 , 2 𝑢 𝑥 𝑛 + 1 = 𝜀 2 2 ‖ ‖ 𝑢 𝑥 𝑛 + 1 ‖ ‖ 2 − ‖ ‖ 𝑢 𝑛 𝑥 ‖ ‖ 2 + ‖ ‖ 2 𝑢 𝑥 𝑛 + 1 − 𝑢 𝑛 𝑥 ‖ ‖ 2 − ‖ ‖ 2 𝑢 𝑛 𝑥 − 𝑢 𝑥 𝑛 − 1 ‖ ‖ 2 + ‖ ‖ 𝑢 𝑥 𝑛 + 1 − 2 𝑢 𝑛 𝑥 + 𝑢 𝑥 𝑛 − 1 ‖ ‖ 2 , 2 𝐴 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 ‖ ‖ 𝑢 = 2 𝐴 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 2 − ‖ ‖ 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 2 ‖ ‖ 𝑢 + 2 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 ‖ ‖ 2 , ( 3 . 1 1 ) the third term of the left of (3.10) is 2 𝜆 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 − 𝑢 𝑛 − 1 = 𝜆 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑢 𝑛 + 1 − 2 𝜆 R e 3 𝑢 𝑛 + 1 − 4 𝑢 𝑛 + 𝑢 𝑛 − 1 , | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | = 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 + | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 + | | 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 | | − 2 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 + 2 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 | | = 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 + | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 | | 𝑢 − 3 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 | | − 2 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 = 𝐼 1 + 𝐼 2 + 𝐼 3 + 𝐼 4 , 𝐼 ( 3 . 1 2 ) 1 | | = 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 = 𝜆 2 | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 + | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − 𝜆 2 | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 = 𝜆 2 ‖ ‖ 2 𝑢 𝑛 + 1 − 𝑢 𝑛 ‖ ‖ 4 𝐿 4 − ‖ ‖ 2 𝑢 𝑛 − 𝑢 𝑛 − 1 ‖ ‖ 4 𝐿 4 − 𝜆 2 Ω | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝐼 𝑑 𝑥 , 2 | | = 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 | | 𝑢 = 𝜆 𝑛 + 1 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 + | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − | | 𝑢 𝑛 + 1 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 = 𝜆 2 | | 𝑢 𝑛 + 1 | | 2 + | | 𝑢 𝑛 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 + 𝜆 2 | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 | | + 𝜆 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − | | 𝑢 𝑛 + 1 | | 2 , | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 = 𝜆 2 ‖ ‖ 𝑢 𝑛 + 1 ‖ ‖ 4 𝐿 4 − ‖ 𝑢 𝑛 ‖ 4 𝐿 4 + 𝜆 2 Ω | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 𝑑 𝑥 + 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − | | 𝑢 𝑛 + 1 | | 2 | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 𝐼 𝑑 𝑥 , 3 = − 2 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝐼 𝑑 𝑥 , 4 = − 3 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 𝑑 𝑥 , ( 3 . 1 3 ) from R e − 2 | 𝑎 | 2 | | 𝑏 | | + 2 2 | | | | − 2 𝑎 − 𝑏 2 = − 4 R e 𝑎 ( 𝑎 − 𝑏 ) , ( 3 . 1 4 ) we have 𝐼 3 + 𝐼 4 = − 2 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝑑 𝑥 − 3 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 𝑑 𝑥 = − 4 𝜆 R e Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 𝑑 𝑥 − 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 𝑑 𝑥 . ( 3 . 1 5 ) Using the definition of 𝐸 𝑛 and from the estimation of (3.11), (3.12), then 𝐸 𝑛 + 1 + 𝜀 2 2 ‖ ‖ 𝑢 𝑥 𝑛 + 1 − 2 𝑢 𝑛 𝑥 + 𝑢 𝑥 𝑛 − 1 ‖ ‖ 2 𝐸 + 𝑟 = 𝑛 , ( 3 . 1 6 ) where ‖ ‖ 𝛿 𝑟 = 4 𝐴 𝑡 𝑡 𝑢 𝑛 ‖ ‖ 2 + 𝜆 2 Ω | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 𝜆 𝑑 𝑥 − 2 Ω | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 2 𝑑 𝑥 − 𝜆 Ω | | 𝑢 𝑛 + 1 | | 2 − | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 𝑑 𝑥 − 4 𝜆 R e Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 𝑑 𝑥 − 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 | | 2 𝑑 𝑥 . ( 3 . 1 7 ) Obviously, (3.9) holds provided that 𝑟 > 0 . In the following we will show that 𝑟 can be made nonnegative if 𝐴 is chosen large enough. To simplify we denote 𝑢 𝑛 + 1 − 2 𝑢 𝑛 + 𝑢 𝑛 − 1 by 𝛿 𝑡 𝑡 𝑢 𝑛 . By simple calculation and applying Young's inequality 𝑎 𝑏 ≤ 𝜖 𝑎 2 + 𝑏 2 / 4 𝜖 ( 𝜖 > 0 ) , to the suitable term, we get ‖ ‖ 𝛿 𝑟 ≥ 4 𝐴 𝑡 𝑡 𝑢 𝑛 ‖ ‖ 2 + 𝜆 2 Ω | | 𝑢 𝑛 + 1 | | 2 − | | 𝑢 𝑛 | | 2 2 − 𝜆 𝑑 𝑥 2 R e Ω | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 + 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 𝑑 𝑥 − 𝜆 R e Ω 𝑢 𝑛 + 1 + 2 𝑢 𝑛 − 𝑢 𝑛 − 1 𝑢 𝑛 + 1 + 𝑢 𝑛 𝑢 𝑛 + 1 − 𝑢 𝛿 𝑡 𝑡 𝑢 𝑛 𝑑 𝑥 − 4 𝜆 R e Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 𝑢 𝑛 + 1 − 𝑢 𝑛 𝛿 𝑡 𝑡 𝑢 𝑛 𝑑 𝑥 − 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 ‖ ‖ 𝛿 𝑑 𝑥 ≥ 4 𝐴 𝑡 𝑡 𝑢 𝑛 ‖ ‖ 2 + 𝜆 2 Ω | | | 𝑢 𝑛 + 1 + 𝑢 𝑛 𝑢 𝑛 + 1 − 𝑢 𝑛 | | | 2 − 𝜆 𝑑 𝑥 4 Ω | | | 𝑢 𝑛 + 1 + 𝑢 𝑛 𝑢 𝑛 + 1 − 𝑢 𝑛 | | | 2 𝑑 𝑥 − 𝜆 Ω | | 𝑢 𝑛 + 1 + 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 𝑑 𝑥 − 𝜆 𝜖 1 Ω | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − 𝜆 𝑑 𝑥 4 𝜖 1 Ω | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 𝑑 𝑥 − 𝜆 𝜖 2 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 − 𝑑 𝑥 4 𝜆 𝜖 2 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 𝑑 𝑥 − 𝜆 Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝛿 𝑡 𝑡 𝑢 𝑛 | | 2 𝑑 𝑥 , ( 3 . 1 8 ) that is, 𝑟 ≥ 4 𝐴 − 𝜆 m a x 𝑥 ∈ Ω | | 𝑢 𝑛 + 1 + 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − m a x 𝑥 ∈ Ω 𝜆 4 𝜖 1 | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − 4 𝜆 𝜖 2 m a x 𝑥 ∈ Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − 𝜆 m a x 𝑥 ∈ Ω | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 ‖ ‖ 𝛿 𝑡 𝑡 𝑢 𝑛 ‖ ‖ 2 + 𝜆 Ω 1 4 | | 𝑢 𝑛 + 1 + 𝑢 𝑛 | | 2 − 𝜖 1 | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 − 𝜖 2 | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 | | 𝑢 𝑛 + 1 − 𝑢 𝑛 | | 2 𝑑 𝑥 . ( 3 . 1 9 )
We choose 𝜖 1 , 𝜖 2 such that the integral term in the above estimate is nonnegative, then constant 𝐴 holds
4 𝐴 ≥ 𝜆 m a x 𝑥 ∈ Ω | | 𝑢 𝑛 + 1 + 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 + 1 4 𝜖 1 | | 2 𝑢 𝑛 + 1 + 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 + 4 𝜖 2 | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 + | | 2 𝑢 𝑛 − 𝑢 𝑛 − 1 | | 2 . ( 3 . 2 0 )
Under condition (3.20) on 𝐴 , we have
𝐸 𝑛 + 1 ≤ 𝐸 𝑛 . ( 3 . 2 1 ) This completes the proof of this theorem.

3.2. Third-Order Scheme: BD3/EP3

𝐴 third-order scheme for solving the semiclassical limit of the defocusing nonlinear Schrödinger equation can be constructed in a similar manner as used in the subsection. Specifically, we can obtain the BD3/EP3 scheme in the following form:

𝑖 𝜀 1 1 𝑢 𝑛 + 1 − 1 8 𝑢 𝑛 + 9 𝑢 𝑛 − 1 − 2 𝑢 𝑛 − 2 + 𝜀 6 Δ 𝑡 2 2 𝑢 𝑛 + 1 𝑥 𝑥 𝑢 − 𝐴 𝑛 + 1 − 3 𝑢 𝑛 + 3 𝑢 𝑛 − 1 − 𝑢 𝑛 − 2 | | − 𝜆 3 𝑢 𝑛 − 3 𝑢 𝑛 − 1 + 𝑢 𝑛 − 2 | | 2 3 𝑢 𝑛 − 3 𝑢 𝑛 − 1 + 𝑢 𝑛 − 2 = 0 ( 𝑛 ≥ 2 ) , ( 3 . 2 2 ) where, in order to start the iteration, 𝑢 1 , 𝑢 2 are calculated via a first- and second-order scheme, respectively. The stability analysis of the scheme (3.22) requires some very detailed energy estimates and will not be presented here. The numerical results obtained in the next section indicate that the third-order time discretization of type (3.22) is also stable as long as the constant 𝐴 is sufficiently large. But comparing with the lower-order scheme, the influence on the stability from 𝐴 -item in the third-order scheme is not very notable.

4. Full Discretization Scheme and Accuracy Tests

A complete numerical algorithm also requires a discretization strategy in space. Since Fourier spectral method is one of the most suitable spatial approximation methods for periodic problems [26–28], here we get the space discretization of semi-implicit scheme from different time iteration schemes by this method.

Let 𝑆 𝑁 𝜑 = s p a n 𝑘 𝑁 ( 𝑥 ) ∣ − 2 𝑁 ≤ 𝑘 ≤ 2 . − 1 ( 4 . 1 )

We use the following Fourier transformations:

1 ̂ 𝑢 ( 𝑘 , 𝑡 ) = 𝑁 𝑁 − 1 𝑗 = 0 𝑢 𝑥 𝑗 𝑒 , 𝑡 − 𝑖 𝑘 𝑥 𝑗 𝑁 , − 2 𝑁 ≤ 𝑘 ≤ 2 𝑢 − 1 , 𝑁 𝑥 𝑗 = , 𝑡 𝑁 / 2 − 1 𝑘 = 𝑁 / 2 ̂ 𝑢 ( 𝑘 , 𝑡 ) 𝑒 𝑖 𝑘 𝑥 𝑗 , 0 ≤ 𝑗 ≤ 𝑁 − 1 . ( 4 . 2 )

The expansion of (2.5) is required to satisfy the following weak form:

𝑢 𝑖 𝜀 𝑁 𝑛 + 1 − 𝑢 𝑛 𝑁 − 𝜀 ▵ 𝑡 , 𝑣 2 2 𝑢 𝑛 + 1 𝑁 𝑥 , 𝑣 𝑥 𝑢 − 𝐴 𝑁 𝑛 + 1 − 𝑢 𝑛 𝑁 | | 𝑢 , 𝑣 − 𝜆 𝑛 𝑁 | | 2 𝑢 𝑛 𝑁 , 𝑣 = 0 , ∀ 𝑣 ∈ 𝑆 𝑁 . ( 4 . 3 )

Similarly, we have the full discretization Fourier spectral approximation of problem (3.2):

𝑖 𝜀 3 𝑢 𝑁 𝑛 + 1 − 4 𝑢 𝑛 𝑁 + 𝑢 𝑁 𝑛 − 1 − 𝜀 2 Δ 𝑡 , 𝑣 2 2 𝑢 𝑛 + 1 𝑁 𝑥 , 𝑣 𝑥 𝑢 − 𝐴 𝑁 𝑛 + 1 − 2 𝑢 𝑛 𝑁 + 𝑢 𝑁 𝑛 − 1 | | , 𝑣 − 𝜆 2 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 | | 2 2 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 , 𝑣 = 0 , ∀ 𝑣 ∈ 𝑆 𝑁 . ( 4 . 4 )

The space discetization of BD3/EP3 (3.22) has a similar form.

It deserves to mention that the corresponding spectral coefficient ̂ 𝑢 𝑛 + 1 can be fully solved because all implicit items about 𝑢 𝑁 𝑛 + 1 in (4.3) and (4.4) are linear.

Choosing test function 𝑣 as the primary function 𝑒 − 𝑖 𝑘 𝑥 of 𝑆 𝑁 , we obtain a system of linear equations for each module 𝑘 in Fourier space:

𝑖 𝜀 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) − ̂ 𝑢 𝑛 ( 𝑘 ) − 𝜀 ▵ 𝑡 2 2 𝑘 2 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) − 𝐴 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) = 𝐴 ̂ 𝑢 𝑛 | | ( 𝑘 ) + 𝜆 ̂ 𝑢 𝑛 | | ( 𝑘 ) 2 ̂ 𝑢 𝑛 ( 𝑘 ) , ( 4 . 5 ) The corresponding second-order BD2/EP2 scheme is as follows:

𝑖 𝜀 3 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) − 4 ̂ 𝑢 𝑛 ( 𝑘 ) + ̂ 𝑢 𝑛 − 1 ( 𝑘 ) − 𝜀 2 Δ 𝑡 2 2 𝑘 2 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) − 𝐴 ̂ 𝑢 𝑛 + 1 ( 𝑘 ) = 𝐴 2 ̂ 𝑢 𝑛 ( 𝑘 ) − ̂ 𝑢 𝑛 − 1 | | ( 𝑘 ) + 𝜆 2 ̂ 𝑢 𝑛 ( 𝑘 ) − ̂ 𝑢 𝑛 − 1 | | ( 𝑘 ) 2 2 ̂ 𝑢 𝑛 ( 𝑘 ) − ̂ 𝑢 𝑛 − 1 ( 𝑘 ) = 0 . ( 4 . 6 ) Visibly solving (4.5) and (4.6) is efficient because we can solve the system of linear equations determined by (4.5) and (4.6) through obtaining the inverse of a simple diagonal matrix. For the problem of fully discretization (4.5) and (4.6) we derive two theorems similar to the energy inequalities in Theorems 2.1 and 3.3 (its proof will be omitted here).

Theorem 4.1. Consider the numerical scheme (4.3), if 𝜆 𝐴 ≥ 2 | | 𝑢 𝑛 𝑁 | | 2 + 𝜆 4 | | 𝑢 𝑁 𝑛 + 1 + 𝑢 𝑛 𝑁 | | 2 , ( 4 . 7 ) the solution of (4.3) satisfies 𝐸 𝑢 𝑁 𝑛 + 1 𝑢 ≤ 𝐸 𝑛 𝑁 , ∀ 𝑛 ≥ 0 , ( 4 . 8 ) where 𝐸 𝑢 𝑛 𝑁 = 𝜀 2 2 ‖ ‖ 𝑢 𝑛 𝑁 𝑥 ‖ ‖ 2 + 𝜆 2 ‖ ‖ 𝑢 𝑛 𝑁 ‖ ‖ 4 𝐿 4 . ( 4 . 9 )

Theorem 4.2. If 𝜆 𝐴 ≥ 4 m a x 𝑥 ∈ Ω | | 𝑢 𝑁 𝑛 + 1 + 2 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 | | 2 + 𝑐 1 | | 2 𝑢 𝑁 𝑛 + 1 + 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 | | 2 + 𝑐 2 | | 2 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 | | 2 , ( 4 . 1 0 ) then the solution of (4.4) satisfies 𝐸 𝑢 𝑁 𝑛 + 1 ≤ 𝐸 𝑢 𝑛 𝑁 , ∀ 𝑛 ≥ 0 , ( 4 . 1 1 ) where 𝑐 1 , 𝑐 2 is positive constant and 𝐸 𝑛 = 𝜀 2 2 ‖ ‖ 𝑢 𝑛 𝑁 𝑥 ‖ ‖ 2 + ‖ ‖ 2 𝑢 𝑛 𝑁 𝑥 − 𝑢 𝑛 − 1 𝑁 𝑥 ‖ ‖ 2 ‖ ‖ 𝑢 + 2 𝐴 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 ‖ ‖ 2 + 𝜆 2 ‖ ‖ 𝑢 𝑛 𝑁 ‖ ‖ 4 𝐿 4 + ‖ ‖ 2 𝑢 𝑛 𝑁 − 𝑢 𝑁 𝑛 − 1 ‖ ‖ 4 𝐿 4 . ( 4 . 1 2 )

Then we investigate the stability properties of different schemes by numerical experiments, where special attention will be paid to the improvements of stability caused by the 𝐴 -item. In the computation we choose Ω = [ − 4 , 4 ] or Ω = [ − 8 , 8 ] as the domain of computation, and choose random number as the initial condition. First, we define Δ 𝑡 𝑐 as the largest time step which makes the computation stable; that is, when the time step is larger than Δ 𝑡 𝑐 the numerical solution will “blow up." For different 𝐴 , Table 1 lines up the values of Δ 𝑡 𝑐 in the Fourier spectral schemes of (2.5), (3.2), and (3.22), respectively, and the Fourier module used in the computation is 𝑁 = 2 5 6 . From observing the data in Table 1, we draw the following conclusion.

tab1
Table 1: Stability comparison with different 𝐴 and 𝜀 when 𝜆 = 𝜀 .

When 𝜆 = 𝜀 is large, we see that the Δ 𝑡 𝑐 is also large. Since all standard semi-implicit schemes are stable, small 𝐴 -item is not necessarily required. When 𝜆 = 𝜀 is small, smaller time step is needed to satisfy the requirement of stability if we use a standard semi-implicit approach. This phenomenon is prominent especially for the case of higher-order or very small 𝜀 . For the first-order scheme BD1/EP1, just as proved in Theorem 2.1, when 𝐴 is large enough the unconditional stability of the computational process can be guaranteed. For the second-order scheme, 𝐴 -item improves the stability obviously. Comparing with the case of 𝐴 = 0 , about 10 times time-step is allowed when we use 𝐴 -item. For the third-order scheme, although the effect of adding 𝐴 -item is not very obvious, adding the 𝐴 -item scheme is still effective in improving the stability and enlarging the time step to a certain extent.

Now we turn to time accuracy comparison. Since the exact solution of (1.1) is unknown, we use numerical results of the BD2/EP2 with Δ 𝑡 = 1 0 − 5 , 𝜀 = 0 . 0 0 2 5 and 𝑁 = 1 0 2 4 as the “exact" solution. We discuss it in the following two cases.

Case 1 (weak 𝑂 ( 𝜀 ) cubic defocusing nonlinearity). We set 𝜆 = 𝜀 , the final time 𝑡 = 1 0 . Table 2 shows the 𝐿 2 -errors obtained, at the same time we derive the rate of convergence about time by the following formulation: 𝐸 l o g 𝑖 / 𝐸 𝑖 + 1 l o g Δ 𝑡 𝑖 / Δ 𝑡 𝑖 + 1 , ( 4 . 1 3 ) where 𝐸 𝑖 is the 𝐿 2 -error when using time-step Δ 𝑡 𝑖 , the rate may be verified experimentally by solving a problem on a sequence of finer and finer partitions and approximation. From Table 2, we can easily calculate 𝑟 ( 𝐵 𝐷 1 / 𝐸 𝑃 1 ) ≈ 1 , 𝑟 ( 𝐵 𝐷 2 / 𝐸 𝑃 2 ) ≈ 2 , 𝑟 ( 𝐵 𝐷 3 / 𝐸 𝑃 3 ) ≈ 3 , which is coincident with our theoretical conclusion.

tab2
Table 2: Numerical errors obtained by using the BD1/EP1, BD2/EP2, and BD3/EP3 with 𝜆 = 𝜀 = 0 . 0 0 2 5 and 𝑁 = 2 5 6 .

Case 2 (strong 𝑂 ( 1 ) cubic defocusing nonlinearity). Figure 1 shows the position density 𝜌 and the current density 𝐽 for 𝜀 = 0 . 0 0 2 5 , Ω = [ − 8 , 8 ] with different values of ( 𝐴 , Δ 𝑡 ) = ( 0 , 0 . 0 0 1 ) , ( 0 . 1 , 0 . 0 1 ) ; ( 0 , 0 . 0 0 2 ) , ( 0 . 1 , 0 . 0 2 ) and in different time 𝑡 = 0 . 5 , 1 . 5 , 1 0 , 2 0 . It is observed that there is a good agreement between the numerical results obtained by using standard semi-implicit time-stepping method (i.e., 𝐴 = 0 ) with small Δ 𝑡 and the modified method (3.2) with larger Δ 𝑡 .

fig1
Figure 1: Numerical solutions for the strong O(1) defocusing nonlinearity by using ED2/EP2. 𝜆 = 1 , 𝜀 = 0 . 0 0 2 5 , 𝑁 = 2 5 6 . “ − “: ( 𝐴 , Δ 𝑡 ) = ( 0 , 0 . 0 0 1 ) , ( 0 , 0 . 0 0 2 ) , “+ + +", “***": ( 𝐴 , Δ 𝑡 ) = ( 0 . 1 , 0 . 0 1 ) , ( 0 . 1 , 0 . 0 2 ) . (a) 𝑡 = 0 . 5 , (b) 𝑡 = 1 . 5 , (c) 𝑡 = 1 0 , (d) 𝑡 = 2 0 .

It is seen that once the methods are stable, the expected order of convergence (in time) is obtained, and larger time-steps can be used by adding an 𝐴 -term.

5. Numerical Examples

In this section, we present the numerical results by simulating the semiclassical limit of defocusing nonlinear Schrödinger equations. The simulations are carried out in the domain Ω = [ − 4 , 4 ] or Ω = [ − 8 , 8 ] , where periodic boundary conditions are used in the spatial directions. The second-order schemes (BD2/EP2), that is, (3.2) are used in our simulations. In the following computations we let 𝑁 = 2 5 6 . The initial condition (1.1) is always chosen in the classical WKB form: 𝑢 0 ( 𝑥 , 𝑡 = 0 ) = 𝑢 0 ( 𝑥 ) = 𝐴 0 ( 𝑥 ) 𝑒 𝑖 𝑆 0 ( 𝑥 ) / 𝜀 = √ 𝜌 0 ( 𝑥 ) 𝑒 𝑖 𝑆 0 ( 𝑥 ) / 𝜀 . ( 5 . 1 )

The analytic solutions of the semiclassical limit are available from [29] and are used for numerical simulations.

Example 5.1 (weak 𝑂 ( 𝜀 ) cubic defocusing nonlinearity). The initial condition is taken as 𝐴 0 ( 𝑥 ) = 𝑒 − 𝑥 2 , 𝑆 0 𝑥 ( 𝑥 ) = − 2 2 + 𝜀 𝑒 − 2 𝑥 2 1 𝑙 𝑛 𝜀 , 𝑥 ∈ ℝ . ( 5 . 2 )
We choose these initial data for this example such that the weak limits of the position density 𝜌 and current density 𝐽 can be obtained analytically [29]. The weak limits 𝜌 0 , 𝐽 0 of 𝜌 , 𝐽 , respectively, as 𝜀 → 0 , are given in [29], for example, before breaking
𝜌 0 1 ( 𝑥 , 𝑡 ) = 𝑒 1 − 𝑡 − 2 𝑥 2 / 1 − 𝑡 , 𝐽 0 𝑥 ( 𝑥 , 𝑡 ) = − ( 1 − 𝑡 ) 2 𝑒 − 2 𝑥 2 / 1 − 𝑡 , 0 ≤ 𝑡 < 1 . ( 5 . 3 ) When 𝑡 → 1 , they are singular distributions (“ 𝛿 -like"). Figures 2, 3, 4, and 5 show variations of the numerical results about the position density 𝜌 and current density 𝐽 for different 𝜀 = 0 . 0 4 , 0 . 0 1 , 0 . 0 0 2 5 , 0 . 0 0 0 6 2 5 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking) and at 𝑡 = 1 0 .

fig2
Figure 2: Numerical solutions of Example 5.1 by ED2/EP2 for the weak 𝑂 ( 𝜀 ) defocusing nonlinearity ( 𝜀 = 0 . 0 4 ). (a) Position density 𝜌 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (b) Current density 𝐽 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (c) Evolution of the position density 𝜌 (left) and the current density 𝐽 (right).
fig3
Figure 3: Numerical solutions of Example 5.1 by ED2/EP2 for the weak 𝑂 ( 𝜀 ) defocusing nonlinearity ( 𝜀 = 0 . 0 1 ). (a) Position density 𝜌 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (b) Current density 𝐽 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (c) Evolution of the position density 𝜌 (left) and the current density 𝐽 (right).
fig4
Figure 4: Numerical solutions of Example 5.1 by ED2/EP2 for the weak 𝑂 ( 𝜀 ) defocusing nonlinearity ( 𝜀 = 0 . 0 0 2 5 ).
fig5
Figure 5: Numerical solutions of Example 5.1 by ED2/EP2 for the weak 𝑂 ( 𝜀 ) defocusing nonlinearity ( 𝜀 = 0 . 0 0 0 6 2 5 ).

Example 5.2 (weak O ( 𝜀 3 / 2 ) cubic defocusing nonlinearity, with 𝜆 = 𝜀 3 / 2 ). The initial condition is taken as 𝐴 0 ( 𝑥 ) = 𝑒 − 𝑥 2 , 𝑆 0 𝑥 ( 𝑥 ) = − 2 2 , 𝑥 ∈ ℝ . ( 5 . 4 ) The semiclassical limits of the position density 𝜌 and current density 𝐽 are given analytically in [29]. For convenience, we only calculate the variations about the position density 𝜌 and current density 𝐽 for 𝜀 = 0 . 0 4 , 0 . 0 1 , 0 . 0 0 2 5 , 0 . 0 0 0 6 2 5 at 𝑡 = 0 . 5 (before breaking), 𝑡 = 1 . 5 (after breaking), and 𝑡 = 1 0 , respectively, (Figures 6, 7, 8, and 9). The computation scheme we used is BD2/EP2 with 𝑁 = 2 5 6 , Δ 𝑡 = 0 . 0 0 0 1 , 𝐴 = 0 . 0 5 .

fig6
Figure 6: Numerical solutions of Example 5.2 by ED2/EP2 for the weak 𝑂 ( 𝜀 3 / 2 ) defocusing nonlinearity ( 𝜀 = 0 . 0 4 ). (a) Position density 𝜌 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (b) Current density 𝐽 at 𝑡 = 0 . 5 (before breaking), at 𝑡 = 1 . 5 (after breaking), at 𝑡 = 1 0 . (c) Evolution of the position density 𝜌 (left) and the current density 𝐽 (right).
fig7
Figure 7: Numerical solutions of Example 5.2 by ED2/EP2 for the weak 𝑂 ( 𝜀 3 / 2 ) defocusing nonlinearity ( 𝜀 = 0 . 0 1 ).
fig8
Figure 8: Numerical solutions of Example 5.2 by ED2/EP2 for the weak 𝑂 ( 𝜀 3 / 2 ) defocusing nonlinearity ( 𝜀 = 0 . 0 0 2 5 ).
fig9
Figure 9: Numerical solutions of Example 5.2 by ED2/EP2 for the weak 𝑂 ( 𝜀 3 / 2 ) defocusing nonlinearity ( 𝜀 = 0 . 0 0 0 6 2 5 ).

Example 5.3 (strong O ( 1 ) cubic defocusing nonlinearity, 𝜆 = 1 ). The initial condition is taken as 𝐴 0 𝑆 ( 𝑥 ) = 1 − | 𝑥 | , | 𝑥 | < 1 , 0 , o t h e r w i s e , 0 ( 𝑥 ) = − 𝑙 𝑛 ( 𝑒 𝑥 + 𝑒 − 𝑥 ) , 𝑥 ∈ ℝ . ( 5 . 5 ) These initial data are not analytic at 𝑥 = 0 and 𝑥 = ± 1 . Resolving the solution for nonanalytic conditions is a challenging task. For numerical study of NLSs with cubic defocusing nonlinearity and analytic initial data, we refer to [17, 18, 20, 30]. To test the numerical method, for each fixed 𝜀 = 0 . 0 4 , 0 . 0 1 , 0 . 0 0 2 5 , 0 . 0 0 0 6 2 5 (Figures 10, 11, 12, and 13), we compute the numerical solution with 𝑁 = 5 1 2 , Δ 𝑡 = 0 . 0 0 0 1 , 𝐴 = 0 . 0 0 5 .

fig10
Figure 10: Numerical solutions of Example 5.3 by ED2/EP2 for the strong 𝑂 ( 1 ) cubic defocusing nonlinearity ( 𝜆 = 1 ) , 𝜀 = 0 . 0 4 .
fig11
Figure 11: Numerical solutions of Example 5.3 by ED2/EP2 for the strong 𝑂 ( 1 ) cubic defocusing nonlinearity ( 𝜆 = 1 ) , 𝜀 = 0 . 0 1 .
fig12
Figure 12: Numerical solutions of Example 5.3 by ED2/EP2 for the strong 𝑂 ( 1 ) cubic defocusing nonlinearity ( 𝜆 = 1 ) , 𝜀 = 0 . 0 0 2 5 .
fig13
Figure 13: Numerical solutions of Example 5.3 by ED2/EP2 for the strong 𝑂 ( 1 ) cubic defocusing nonlinearity ( 𝜆 = 1 ) , 𝜀 = 0 . 0 0 0 6 2 5 .

6. Conclusions

The semiclassical limit of the defocusing nonlinear Schrödinger equation presents a great computational challenge. Not only is the numerical method required to resolve the solution accurately, but also required to have high resolution ratio in temporal and spatial directions due to its high oscillation of solutions. In this work, a large time-stepping spectral approximation for the semiclassical limit of the defocusing nonlinear schrödinger equation is studied. It is demonstrated that the classical semi-implicit method can be further improved by simply adding a linear term consistent with the truncation errors in time direction. This treatment can be used to increase the time-step size and computational stability. Our numerical results show that: this method is very effective in capturing oscillatory solutions of the NLS for small 𝜀 . We believe that the method presented here is a tool of great value and can be used to learn more about the structure of the solutions of limiting behavior for the semiclassical limit of the defocusing nonlinear schrödinger equation and some other numerical simulations.

Acknowledgments

The research is supported by Fujian National Science Foundation of China Grant 2008J0198. The author would like to thank Professor Chuanju Xu for insightful discussions regarding the nonlinear Schrödinger equation and Professor Yinnian He for helpful suggestions and comments.

References

  1. S. Jin, C. D. Levermore, and D. W. McLaughlin, “The semiclassical limit of the defocusing NLS hierarchy,” Communications on Pure and Applied Mathematics, vol. 52, no. 5, pp. 613–654, 1999. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  2. A. Jüngel, “Nonlinear problems in quantum semiconductor modeling,” Nonlinear Analysis, vol. 41, pp. 669–688, 2000. View at Publisher · View at Google Scholar · View at MathSciNet
  3. A. Jüngel, Quasi-Hydrodynamic Semiconductor Equations, vol. 41 of Progress in Nonlinear Differential Equations and Their Applications, Birkhäuser, Basel, Switzerland, 2001. View at Zentralblatt MATH · View at MathSciNet
  4. L. Laudau and F. Lifschitz, Quantum Mechanics: Non-Relativistic Theory, Pergamon Press, New York, NY, USA, 1977.
  5. M. H. Anderson, J. R. Ensher, M. R. Matthews, C. E. Wieman, and E. A. Cornell, “Observation of Bose-Einstein condensation in a dilute atomic vapor,” Science, vol. 269, no. 5221, pp. 198–201, 1995. View at Publisher · View at Google Scholar · View at PubMed
  6. J. R. Anglin and W. Ketterle, “Bose-Einstein condensation of atomic gases,” Nature, vol. 416, no. 6877, pp. 211–218, 2002. View at Publisher · View at Google Scholar · View at PubMed
  7. P. D. Miller and S. Kamvissis, “On the semiclassical limit of the focusing nonlinear Schrödinger equation,” Physics Letters A, vol. 247, no. 1-2, pp. 75–86, 1998. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  8. J. C. Bronski and J. N. Kutz, “Numerical simulation of the semi-classical limit of the focusing nonlinear Schrödinger equation,” Physics Letters A, vol. 254, no. 6, pp. 325–336, 1999. View at Publisher · View at Google Scholar
  9. M. G. Forest and K. T.-R. McLaughlin, “Onset of oscillations in nonsoliton pulses in nonlinear dispersive fibers,” Journal of Nonlinear Science, vol. 8, no. 1, pp. 43–62, 1998. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  10. S. Kamvissis, “Critical and subcritical cases of the Toda shock problem,” in Singular Limits of Dispersive Waves (Lyon, 1991), vol. 320 of NATO Advanced Science Institutes Series B: Physics, pp. 257–271, Plenum, New York, NY, USA, 1994. View at Zentralblatt MATH · View at MathSciNet
  11. M. G. Forest and J. E. Lee, “Geometry and modulation theory for the periodic nonlinear Schrödinger equation,” in Oscillation Theory, Computation, and Methods of Compensated Compactness (Minneapolis, Minn., 1985), vol. 2 of IMA Vol. Math. Appl., pp. 35–69, Springer, New York, NY, USA, 1986. View at Zentralblatt MATH · View at MathSciNet
  12. O. C. Wright, M. G. Forest, and K. T.-R. McLaughlin, “On the exact solution of the geometric optics approximation of the defocusing nonlinear Schrödinger equation,” Physics Letters A, vol. 257, no. 3-4, pp. 170–174, 1999. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  13. P. G'rard, “Microlocal defect measures,” Communications in Partial Differential Equations, vol. 16, no. 11, pp. 1761–1794, 1991. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  14. L. Tartar, “H-measures, a new approach for studying homogenisation, oscillations and concentration effects in partial differential equations,” Proceedings of the Royal Society of Edinburgh A, vol. 115, no. 3-4, pp. 193–230, 1990. View at Zentralblatt MATH · View at MathSciNet
  15. I. Gasser and P. A. Markowich, “Quantum hydrodynamics, Wigner transforms and the classical limit,” Asymptotic Analysis, vol. 14, no. 2, pp. 97–116, 1997. View at Zentralblatt MATH · View at MathSciNet
  16. P. A. Markowich, N. J. Mauser, and F. Poupaud, “A Wigner-function approach to (semi)classical limits: electrons in a periodic potential,” Journal of Mathematical Physics, vol. 35, no. 3, pp. 1066–1094, 1994. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  17. P. A. Markowich, P. Pietra, and C. Pohl, “Numerical approximation of quadratic observables of Schrödinger-type equations in the semi-classical limit,” Numerische Mathematik, vol. 81, no. 4, pp. 595–630, 1999. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  18. P. A. Markowich, P. Pietra, C. Pohl, and H. P. Stimming, “A Wigner-measure analysis of the Dufort-Frankel scheme for the Schrödinger equation,” SIAM Journal on Numerical Analysis, vol. 40, no. 4, pp. 1281–1310, 2002. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  19. W. Bao, S. Jin, and P. A. Markowich, “Numerical study of time-splitting spectral discretizations of nonlinear Schrödinger equations in the semiclassical regimes,” SIAM Journal on Scientific Computing, vol. 25, no. 1, pp. 27–64, 2003. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  20. W. Bao, S. Jin, and P. A. Markowich, “On time-splitting spectral approximations for the Schrödinger equation in the semiclassical regime,” Journal of Computational Physics, vol. 175, no. 2, pp. 487–524, 2002. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  21. H. D. Ceniceros, “A semi-implicit moving mesh method for the focusing nonlinear Schrödinger equation,” Communications on Pure and Applied Analysis, vol. 1, no. 1, pp. 1–18, 2002. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  22. H. D. Ceniceros and F.-R. Tian, “A numerical study of the semi-classical limit of the focusing nonlinear Schrödinger equation,” Physics Letters A, vol. 306, no. 1, pp. 25–34, 2002. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  23. C. Xu and T. Tang, “Stability analysis of large time-stepping methods for epitaxial growth models,” SIAM Journal on Numerical Analysis, vol. 44, no. 4, pp. 1759–1779, 2006. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  24. Y. He, Y. Liu, and T. Tang, “On large time-stepping methods for the Cahn-Hilliard equation,” Applied Numerical Mathematics, vol. 57, no. 5–7, pp. 616–628, 2007. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  25. B. Li and J.-G. Liu, “Thin film epitaxy with or without slope selection,” European Journal of Applied Mathematics, vol. 14, no. 6, pp. 713–743, 2003. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  26. J. P. Boyd, Chebyshev and Fourier Spectral Methods, Dover, Mineola, NY, USA, 2nd edition, 2001. View at MathSciNet
  27. D. Gottlieb and S. A. Orszag, Numerical Analysis of Spectral Methods, SIAM, Philadelphia, Pa, USA, 1977.
  28. E. Tadmor, “Super-viscosity and spectral approximations of nonlinear conservation laws,” in Numerical Methods for Fluid Dynamics, 4 (Reading, 1992), pp. 69–81, Oxford University Press, New York, NY, USA, 1993. View at Zentralblatt MATH · View at MathSciNet
  29. R. Carles, “Remarques sur les mesures de Wigner,” Comptes Rendus de l'Académie des Sciences. Série I. Mathématique, vol. 332, no. 11, pp. 981–984, 2001. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  30. S. Jin, C. D. Levermore, and D. W. McLaughlin, “The behavior of solutions of the NLS equation in the semiclassical limit,” in Singular Limits of Dispersive Waves (Lyon, 1991), vol. 320 of NATO Advanced Science Institutes Series B: Physics, pp. 235–255, Plenum, New York, NY, USA, 1994. View at Zentralblatt MATH · View at MathSciNet