Motion and Interaction of Envelope Solitons in Schrödinger Equation Simulated by Symplectic Algorithm

The expression of Gaussian envelope soliton in Schrödinger equations are given and proved in this paper. According to the characteristics of the Gauss envelope soliton, further proposed that the interaction between Gaussian envelope solitons exists in Schrödinger equation. The symplectic algorithm for solving Schrödinger equation is proposed after analysis characteristics of Schrödinger equation. First, the Schrödinger equation is transformed into the standard Hamiltonian canonical equation by separating the real and imaginary parts of wave function. Secondly, the symplectic algorithm is implemented by using the Euler center difference method for the canonical equation. The conserved quantity of symplectic algorithm is given, and the stability of symplectic algorithm is proved. The numerical simulation experiment was carried out on Schrödinger equation in Gauss envelope soliton motion and multi solitons interaction. The experimental results show that the proposed method is correct and the symplectic algorithm is effective.


Introduction
The research and application of solitons have gone deep into various fields, such as mathematics, fluid dynamics, nonlinear electromagnetism, optical fiber communication, solid state physics and so on. The definition of soliton is not the same in all research fields. Scholars of pure theory, such as mathematicians, define the soliton very strictly. From the mathematical point of view, the soliton is a stable, energy limited undispersed solution of some nonlinear partial differential equations, that is, it can always keep the wave and velocity unchanged. Scholars who value applications, such as some applied physicists, are more relaxed about the definition of soliton. Physicists believe that as long as the energy of the wave is limited and is distributed in a limited space or time range, even if the wave changes during the propagation process, such as the high order optical soliton in the optical fiber, they are also called soliton. There are many forms of solitons, such as bell shaped soliton, ring soliton, twisted soliton, enveloping soliton, anti soliton, sentry soliton, respiratory soliton, etc.
Soliton has experienced a historical process from discovery to establish mathematical model and application. In 1834, John Scott Russell accidentally discovered the soliton waves in the canal. In 1895, D J Korteweg and De Vries proposed the KdV equation of one-way wave propagation in fluids, theoretically proved the existence of solitons. In 1955, Enrico Fermi, John Pasta and Stan Ulam confirmed the existence of solitons using computer numerical experiments. In 1962, Perring and Skyrme used Sine-Gordon equation to study the basic particles. The calculation results show that the soliton is spread out, even if the two soliton collision also still maintain the original shape and speed. In 1965, M. D. Kruskal and N. Zabusky used numerical simulation method to study the nonlinear interaction process of the soliton collision in plasma. The results further confirmed that the soliton wave remains constant shape and speed after the collision. At present, the soliton research has entered the application stage, especially the application of optical soliton in optical fiber communication field has achieved a lot of results [1][2][3][4].
The form of soliton is mainly obtained by decoupling the equation. For the interaction process of solitons, it can be observed, verified and tested by numerical simulation [5][6].
There are many methods of numerical simulation, and the most commonly used method is the finite difference method. Symplectic algorithm is a numerical algorithm suitable for Hamiltonian system equations.
In 1984, Feng Kang [7] made a presentation report entitled "difference scheme and the symplectic geometry", first proposed the symplectic geometry algorithm for Hamiltonian system, in the Beijing conference of differential geometry and differential equation computation of partial differential equations. Since then, Feng Kang and his team system [8][9] developed this algorithm in. The symplectic algorithm has obvious advantages compared with then on symplectic algorithm when the numerical calculation for the Hamiltonian systems. Because of maintaining the symplectic structure, the symplectic algorithm makes the numerical calculation of each step be a symplectic transformation, eliminate the artificial dissipation introduced non symplectic algorithm in the calculation process. Therefore, it has the long time tracking ability and high stability for the evolutionary computation [10][11][12][13][14][15].
The Schrödinger equation belongs to Hamiltonian system, it's mathematical framework is the theory of Hilbert spaces. Calculate the object of symplectic algorithm is a Hamiltonian system, it's mathematical framework is symplectic geometry. Therefore, the transform is necessary from the Schrödinger's equations to the symplectic algorithm. The Schrödinger equation is complex. If the real part and the imaginary part of the wave function ware regarded as the generalized coordinates and generalized momentum. The Schrödinger equations can be written in standard form canonical equations of Hamilton system.
In the Schrödinger equation, the researchers are concerned with the mode square (wave packet) of the wave function instead of the direct solution of the Schrödinger equation.
Because the direct solution of the Schrödinger equation is a conjugate complex number in pairs which physical meaning is unknown, but the wave packet is a real number which physical meaning is clear.
The researching objects of this paper is the Gauss envelope soliton in Schrödinger equations, their motion and interaction process, using the method of combining the analytical method and symplectic algorithm. On the one hand, the correctness of symplectic algorithm will be verified, on the other hand the characteristics of solitons motion and interaction will be also verified by the numerical simulation using symplectic algorithm.

Analytical Solutions and Interaction of Envelope Solitons in the Schrödinger Equation
The one-dimensional stationary Schrödinger equation is shown in (1).
Where, i is the imaginary symbol; is the potential function. x is the position of quantum which is one-dimensional (1D), 2D or 3D.

The Schrödinger Oscillator Equation
In theory, any potential function can be approximation using the parabolic at the nearby of minima in the vicinity can use according to the Taylor series expansion principle. In practice, any vibration, as long as the amplitude is small enough, can be seen as simple harmonic vibration.
In the case of harmonic vibration, potential function is shown in (2).
Where, ω is the resonant frequency. The Schrödinger equation (1) can be written to the form such as (3) when potential function is the harmonic oscillator function.

The Envelope Soliton Solution of Schrödinger Harmonic Oscillator Equation
The Schrödinger harmonic oscillator equation has the envelope soliton solution. The envelope soliton solutions must satisfy two conditions: A. The wave function satisfies the Schrödinger harmonic oscillator equation as (3); B.The wave function does not change with time which always keeps the properties of soliton. The (4) is a envelope soliton solution of Schrödinger harmonic oscillator equation.
Where, a is a real number presenting the length scale.
The condition A is proved as the following: are calculated such as The Schrödinger harmonic oscillator equation will be established if (5) and (7) substitute into (3).
The condition B is proved as the following: The wave packet is shown in (8) calculated by wave function (4).
From equation (9), (10) and (11) can be seen: (1) the shape of x t ψ is invariant and Gaussian; (2) the center µ of wave packet is cosine with the amplitude a and the oscillating frequency ω ; (3) The wave packet width is determined by σ , and its size is inversely proportional to ω . σ is smaller, wave packet more narrow, kurtosis larger, the corresponding oscillation frequency ω larger, i.e. motion speed greater. In the Hamilton system, the generalized coordinates and the generalized momentum are a pair of variables. The other observations variables can be to transform from them. Therefore the calculation of these two variables is very important. The Schrödinger equation belongs to the Hamilton system. The generalized coordinate and the generalized momentum are the variables of the probability type. Their mathematical expectations are usually.
The mathematical expectations of the coordinate is calculated by (12) are shown.
Where, y is the intermediate variable, and cos y x a t ω = − . From (12) it can be seen: the mathematical expectation of coordinate is the cosine wave with amplitude a and frequency ω . The result is equal to the center position of the wave packet, which is the same (10). It reflected from a side that the wave packet maintains the symmetry of shape in movement. The mathematical expected is only equal to the value of the center if the wave packet is symmetric.
The mathematical expectation of momentum can be computed using (13).
From (13) it can be seen that the mathematical expectation of momentum is the sine wave with the amplitude maω and the frequency ω .
The relationship between coordinate and momentum mathematical expectation can be shown as (14) from (12) and (13).
According to the formula (14), the relationship between the coordinate and momentum mathematical expectation is elliptical. The trajectory is closed and elliptic in the phase space. In fact, it conclusion is consistent with the classical harmonic oscillator.

The Interaction Amongst Envelope Solitons in the Schrödinger Equation
According to the formula (4), the shape and motion of wave function are relevant with a and ω . In order to make the distinction of multi solitons, the wave function is expressed as ( ) , , , x a t ψ ω in this paper. When multiple solitons interact, their wave functions overlay as shown in (15).
Where, n is the number of solitons. The wave packet is the square of the wave function's models as shown in (16) when the multiple solitons interact.

Symplectic Algorithm for Schrödinger Equation
There are a variety of numerical solutions for the Schrödinger equation. Compared with other numerical methods, symplectic algorithm has outstanding advantages. It is shown that the algorithm can preserve the square of the x t ψ conservation of the wave function. This advantage is the symplectic properties of symplectic algorithm for solving Schrödinger equation numerically.

The Symplectic Algorithm for Schrödinger Equation
In (1), H is the Hamiltonian operator as shown in (17).
Therefore, (1) can be rewritten in matrix form as shown in (18).
Where,ψ is a complex matrix.
The (23) is available after separating the real and imaginary parts of (22).
The (23) can be abbreviated to (24) which is the Hamiltonian canonical form of Schrödinger equation.
Since H is a Hermite matrix, (27) is established.
According to (28), it can be proved that K is a symplectic matrix, that is (29).
Where, J is a symplectic unit matrix as shown in (30).
Where, E is a unit matrix. For (26), there are some numerical calculation methods.
The Euler center format is symplectic, as shown in (31).
Where, I is a unit matrix. Therefore, the iterative numerical solution of the symplectic scheme for (25) is shown in (34). Where,

The Symplectic Conservation Quantity in Schrödinger Equation
The general solution of the Schrödinger equation (18)  In the (37), ( ) is the time evolution matrix, and is called ( ) . Because G is a complex exponential matrix, its determinant is 1, that is (38).

( )
After the comparison (34) and the (37), it can be found that F is equivalent to the time evolution matrix G . The difference is that F is 2 2 complex exponential matrix. Because F is a symplectic matrix, its determinant is 1 according to the properties of symplectic matrix.  Where, n is the number of iterations. According to the formula (40), the form energy of the symplectic algorithm of (34) is the square of the modulus of wave function ( ) It is in conformity with the physical meaning of the Schrödinger equation which is that the cumulative probability is a constant.

The Stability of the Symplectic Algorithm of the Schrödinger Equation
According to the formula (35), the F is a complex matrix. It can been seen that K is an anti Hermite matrix from (29). So there is the (41).
According to (35), F is the Cayley transformation of 2 t ∆ K . According to the characteristics of the Cayley transformation, the available (42) is available.
So, F is a unitary matrix. According to the properties of unitary matrix, the modulus of all the eigenvalues of F are 1, that is: Since all the eigenvalues of F are not greater than 1, the (34) is always stable.

Simulation Experiment Results and Analysis
This paper simulated the soliton motion and multiple solitons interaction using symplectic algorithm to numerical calculating. The calculation results using symplectic algorithm were compared with that using theoretical analytical solution. The simulation parameters are: Planck constant

A Soliton Motion Process and Phase Trajectory
The wave function parameters are 1 a = , 5 ω = π . First, ( ) , 0 x ψ was calculated by (4) as the initial wave function. Then, it was substituted into (35) and iterative calculate. The snap shots of wave package at 0, 2,5 t = are shown in Figure 1. This results are consistent with the theoretical calculation equation (9).
In the case of a single envelope soliton motion, the energy n E calculated according to the (40) is shown in Figure 2. As shown in Figure 2, the conserved quantity n E of the symplectic algorithm has not changed in the numerical iterative calculation process, which is consistent with the theoretical analysis.
The coordinate mathematical expectation values x of a soliton motion change with time as shown in Figure 3. As you can see from Figure 3: The coordinate mathematical expectation value x is a cosine wave, which the amplitude and frequency are consistent with (12) of the theoretical calculation.
The relationship of single soliton motion between the coordinate mathematical expectation value x and the momentum expectation value p is shown in Figure 4. As you can see from Figure 4: (1) The relationship between x A and p is elliptical, which is consistent with the theoretical calculation (14). (2) After 2000 times iterations, four cycle calculation, x p − still maintains a good numerical relations, which the four cycle trajectories completely overlap. It embodies the advantage of the symplectic algorithm.

Two Solitons Interaction
The wave function parameters are 1    As shown in Figure 6, the conserved quantity n E of the symplectic algorithm has not changed in the numerical iterative calculation process, which is consistent with the theoretical analysis.

Three Solitons Interaction
The wave function parameters are 1    As shown in Figure 8, the conserved quantity n E of the symplectic algorithm has not changed in the numerical iterative calculation process, which is consistent with the theoretical analysis.

The Distribution of Eigenvalues of Matrix F
When 0.02 t ∆ = , x ∆ is 2, 0.20.02, the distribution of eigenvalues of F is shown in Figure 9. The all eigenvalues are distributed on the unit circle in the figure 9, that is =1 i λ , which is consistent with (46). This shows that the algorithm is stable.

Conclusions
When potential well is the harmonic oscillator, there are envelope solitons of the wave packet in Schrödinger equation. The shapes of solitons do not vary with time. The solitons can maintain the shapes respectively after two or more solitons colliding each other. Symplectic algorithm has the advantages in the numerical computation for the Hamiltonian system. The stability of numerical computation is very outstanding. Next, we will continue to study the shape and interaction of envelope solitons in 2D and 3D Schrödinger equation.