A high-accuracy conservative numerical scheme for the generalized nonlinear Schrödinger equation with wave operator

1.   Introduction

In ^[1], Matsunchi proposed the generalized nonlinear Schrödinger equation (NLSE) with wave operator, which reads:

where $\gamma$, $\alpha$, $\theta$, $\lambda$, and $\beta$ are real constants, $i^2 = -1$, $x\in \mathcal{R}$ and $0 < t < T$, which describes the nonlinear interaction between two quasi-monochromatic waves. The nonlinear Schrödinger equation has many important applications in different fields as a vital mathematics model, such as Langmuir wave envelope approximation in plasma physics ^[2], water waves, and bimolecular dynamics ^[3,4], nonlinear topics ^[5,6], and references therein.

To solve the Eq (1.1), we set it up in a compact subset $[x_l, x_r]$ ^[7]. Then, the initial and periodic boundary conditions are added as follows:

There is the energy-conserved property of (1.1)–(1.3) ^[8]:

In terms of numerical study, Wang and Kong et al. ^[8] considered the multi-symplectic preserving integrator for the Schrödinger equation with wave operator. Moreover, in the paper, the authors discussed mainly the conservative properties, without the necessary convergence analysis of the scheme. In the case of $\theta = \gamma = \lambda = 0$, Guo and Liang^[9] developed a nonconservative implicit difference scheme to solve the NLSE with a wave operator; Zhang ^[10] proposed an explicit conservative difference scheme that was conditionally stable for it. In ^[11], several unconditionally stable conservative schemes were shown. However, the above-mentioned schemes were all of second-order accuracy in space. Furthermore, Li and Zhang ^[12] designed a conservative scheme of (1.1). But the scheme is nonlinear implicit, so it is not suitable for parallel computation because it needs heavy iterative calculations.

Recently, high-accuracy computational methods have been attracted by many researchers. In recent works, the high-order accuracy approximation methods were proposed to study the Klein–Gordon equation^[13], the Schrödinger equation^[14,15], the Klein–Gordon–Schrödinger equation^[16,17], nonlinear wave equations^[18], respectively. In addition, for wide and interesting topics covered, the numerical studies should also be recalled in the literatures. Hu ^[19] presented compact conservative schemes for the coupled nonlinear Schrödinger system. Dehghan et al. ^[20,21] studied the high-order solution for Sine–Gordon equation, heat and advection–diffusion equations of one dimension, respectively. Later, this team ^[22,23] solved the NLSE with constant and variable coefficients, and 2D Rayleigh–Stokes problem by compact finite difference method. In ^[24], an efficient and compact finite difference scheme is developed for the Klein–Gordon–Zakharov equation. In ^[25], the NLSE was solved by Fourier pseudo-spectral method. Especially in ^[26,27,28], some efficient linearly implicit and high-order energy-preserving schemes were proposed for Hamiltonian systems, and monotonicity-preserving ones for wave equations. These useful methods inspire us to establish an efficient numerical method for the generalized NLSE with wave operator. In this article, a novel energy-preserving approximation scheme is designed to solve (1.1)–(1.3) with the following advantages: The scheme is high accuracy, unconditionally stable, and convergent, whose theoretical accuracy is $O(\tau^2+h^4)$; the scheme preserves the physical conservative property of the original system; and the proposed method is linearized, which significantly reduces the computational cost compared with the nonlinear one.

The outline of this article is as follows: In Section 2, a linearized high-accuracy energy-preserving scheme for (1.1) is described. The simulation of conservative property and error estimates of the scheme are shown in Section 3. In Section 4, we prove the convergence of the scheme. In Section 5, several useful numerical examples are given to test the theoretical results.

2.   High-accuracy numerical scheme

In this section, we define the solution domain $\bar{\Omega} = \{(x, t)|x\in[x_l, x_r], t\in[0, T]\}$, which is covered by the uniform grid $\bar{\Omega}_{h\times\tau} = \{(x_j, t^n)|x_j = x_l+jh, t^n = n\tau, 0\leq j\leq M, 0\leq n\leq N\}$, where $h = \frac{x_r-x_l}{M}$ is spatial step, $\tau = \frac{T}{N}$ is temporal step. Denote $U^n_{j}\approx u(x_j, t_n)$. Denote discrete grid function $\omega = \{\omega^n_j; j = 0, 1, 2, \cdots, M, n = 0, 1, 2, \cdots, N\}$ on $\bar{\Omega}_{h\times\tau}$. Define:

In the article, the constant $C$ is general positive and independent of mesh parameters $h$ and $\tau$ at different circumstances.

According to the operators above, the high-accuracy linearized energy-preserving scheme for (1.1)–(1.3) is derived:

where $v^{n+\frac{1}{2}}_j = \frac{v^{n+1}_j+v^{n}_j}{2}$.

Assume that $n = 0$ is valid for (2.1). Applying (2.2), we have

3.   Conservative property and error estimate

Let $Z^0_h = \{\mathbf{V}|\mathbf{V} = (V^n_0, V^n_1, \cdots, V^n_{M-1})^{T}, V^n_0 = V^n_M, V^n_{-1} = v^n_{M-1}, V^n_{-2} = V^n_{M-2}\}$. For $\forall \mathbf{\phi}, \mathbf{\varphi}\in Z^0_h$, define:

Next, we discuss the conservative property and error estimate of (2.1)–(2.4).

Lemma 3.1. ^[29] $\forall \; V, W\in Z^0_h$, we obtain

Then one has

Lemma 3.2. ^[30] For any grid function $V\in Z^0_h$, there is

Lemma 3.3.^[31] For $\forall \; V\in Z^0_h$, we obtain

Theorem 3.1. The difference scheme (2.1) inherits the property of energy conservation of the original system (1.1)–(1.3):

Proof. By Lemma 3.1, do the inner product of (2.1) with $\delta_tU^{n}+\delta_tU^{n-1}$. Then take the real part:

Noting that

By using (2.3), we obtain

Similarly, we obtain

Substituting (3.3)–(3.5) into (3.2). Let

which implies (3.1).

Lemma 3.4. Assume that $u_0\in H^1_{per}[x_l, x_r]$, then the solution of (1.1)–(1.3) is estimated:

Proof. It follows from (1.4) that

For the parameters $\lambda$ and $\theta$ that satisfy $1-\frac{1}{2}|\theta| > 0$ and $\lambda-\frac{1}{2}|\theta| > 0$, we obtain

Thus, $||u||_{L_{\infty}}\le C$ follows by Sobolev inequality.

Theorem 3.2. For the scheme of (2.1), its solution satisfies the following estimation: $||\delta_tU^{n}||\le C, ||U^n||\le C, ||\delta_x U^n||\le C$, $||U^n||_{\infty}\le C$.

Proof. From (3.1), we obtain

For $\mathbf{U^n}\in Z^0_h$, we have

Thus

Similarly, we obtain

This, together with (3.10), (3.8), and Lemma 3.2, gives the following:

For $\lambda$ and $\theta$ that satisfy $\frac{1}{2}-\frac{5}{12}|\theta| > 0$ and $\frac{1}{2}\lambda-\frac{5}{12}|\theta| > 0$, we have

It follows from Lemma 3.3 that

4.   Convergence

Let $\omega^n_j = u(x_j, t_n)$. Together with (2.4), the truncation error of (2.1)–(2.3) is defined:

Applying Taylor expansion, there is $|r^n|+|\sigma^0| = O(\tau^2+h^4)$.

Next, we shall carry out the convergence analysis of the present scheme.

Lemma 4.1. ^[32] Suppose that the mesh function $\{v^n, n = 1, 2, \cdots, N; N\tau = T\}$ satisfies the inequality

the constants $B$, $C$, and $A_n$ are nonnegative. Thus

where, $\tau$ small enough stafies $(B+C)\tau \le \frac{N-1}{2N} (N > 1)$.

Theorem 4.1. Suppose that $u(x, t)\in C^{6, 3}_{x, t}$. The numerical solution $U^n$ of the finite difference scheme (2.1)–(2.3) is convergent to the solution of (1.1)–(1.3) in the $||\cdot||_{\infty}$ norm with the convergent rate $O(\tau^2+h^4)$.

Proof. Let $e^n_j = \omega^n_j-U^n_j$. From (4.1)–(4.4) and (2.1)–(2.4), we have

Multiply both sides of (4.5) with $\delta_{t}e^n+\delta_{\bar{t}}e^n$. Take the real parts:

where

Computing the fifth, sixth, and last terms on the right-hand side of (4.9), we have

In addition, we have

and

It follows from (4.9)–(4.16) that

Let $\Phi^n = ||\delta_te^{n}||^2+\frac{1}{2}(||\delta_xe^{n+1}||^2+||\delta_xe^{n}||^2)+\frac{\lambda}{2}(||e^{n+1}||^2+||e^{n}||^2)$. By Lemma 3.3, (4.17) can be rewritten as follows:

Using Lemma 4.1 yields

In (4.6), multiplying $e^1$ yields

According to Lemma 3.2, (4.20) together with $|\sigma^0| = O(\tau^2+h^4)$, $(\sigma^0, e^1)\leq \frac{1}{2}(||\sigma^0||^2+||e^1||^2)$, and $\tau$ small enough gives

Consequently

It follows from (4.19) that

Thus

From Lemma 3.3, we obtain

5.   Numerical experiments

In this section, numerical tests are given to verify the theoretical analysis of the different solutions. To test the space accuracy of the scheme, we denote

5.1.   Plane wave solution

Consider the initial and periodic boundary problems:

The exact periodic solution of (5.1)–(5.3) is known as:

In computations, we take $x_l = 0, x_r = 2\pi$, and $T = 10$. Setting $t = 0$, $u_0(x)$ and $u_1(x)$ are derived from (5.4). Wang and Kong proposed a second-order accuracy MI scheme to study (5.1)–(5.3) ^[8]. The MI scheme and the presented scheme are denoted as Schemes I and II, respectively. In view of all the schemes having 2-order accuracy in time, we mainly consider the accuracy in space in the tests. The comparison between Scheme I and II is provided in Table 1 and Figure 1. Table 1 includes the errors, convergence order, and CPU time of both schemes. From Table 1 and Figure 1, it is clear that our method is much better than the other in ^[8]. Table 1 also demonstrates the present scheme (2.1)–(2.4) is fourth-accuracy in space. In view that the discrete energy is complex, the imaginary and real parts of the conservative property of $E^n$ for Scheme II have been shown in Figure 2 under the temporal $\tau = \frac{1}{640}$ and the spatial $h = \frac{pi}{40}$, respectively. Figure 2 shows that the scheme (2.1) inherits the energy-conservative property very well, which the original problem possesses.

5.2.   Solitary wave solution

The initial value problem of (1.1)–(1.3):

In experiments, $g_0(x) = A\operatorname{sech}(Kx)$, $g_1(x) = ivA\operatorname{sech}(Kx)$ are chosen with the parameters $A = |K|$, and $v = \frac{1}{2}(-1\pm \sqrt{1-4K^2})$ ^[33]. We take $K = \frac{1}{3}$, $v = \frac{1}{2}(-1-\sqrt{1-4K^2})$. The comparisons of convergence order, CPU time, and the errors in the $\|\cdot\|_{\infty}$ norm are shown in Table 2 under different spatial $h$ and temporal $\tau$. The values of $E^n$ are shown in Table 3. Both Tables 2 and 3 all demonstrate the accuracy and effectiveness of the present scheme in this article.

6.   Conclusions

In this paper, an attempt was made to design a novel energy-conserved numerical scheme to solve the initial and periodic boundary problem of the generalized nonlinear Schrödinger equation with a wave operator. The proposed scheme possesses the following merits: Coupling with the Richardson extrapolation, the scheme is linear, of high accuracy $O(h^4+\tau^2)$, and without any restrictions of mesh steps; the presented scheme is energy-conserved and inherits the conservative property that the original system possesses. The convergence analysis of the scheme is discussed in detail. Numerical tests further illustrate the effectiveness of the scheme.

Acknowledgments

This work is supported by the Natural Science Foundation of China (No.11401438), A Project of Shandong Province Higher Educational Science and Technology Program (No. J15LI56), The Doctoral Research Foundation of WFU (2019BS01).

Conflict of interest

The author declares no conflict of interest.