A Trapezoidal-Like Integrator for the Numerical Solution of One-Dimensional Time Dependent Schrödinger Equation

In this paper, the one-dimensional time dependent Schrödinger equation is discretized by the method of lines using a second order finite difference approximation to replace the second order spatial derivative. The evolving system of stiff Ordinary Differential Equation (ODE) in time is solved numerically by an L-stable trapezoidal-like integrator. Results show accuracy of relative maximum error of order 10−4 in the interval of consideration. The performance of the method as compared to an existing scheme is considered favorable.


Introduction
Numerical methods have played an important role in researches in physics.Most often Partial Differential Equations (PDEs) are used to model physical problems.Proficiency in the numerical approach to the solution of PDEs is a great tool in the hand of a mathematical physicist.According to [1], numerical methods of PDEs have played an important role in various areas of mathematical physics.In numerical relativity a set of highly nonlinear systems of coupled differential equations of general relativity have to be solved under very general symmetry conditions in order to predict the profile of gravitational waves expected to be observed by ground-based gravitational wave detectors.In Euler's equation, the hydro-dynamical equations are solved by numerical me-thods in order to explain and predict astrophysical phenomena.In the evolution of a Bose condensate, the solution of the time dependent Schrodinger equation involved is sought for by numerical method.
Schrodinger equation has a wide range of application in physics varying from optics, propagation of electric field in optical fibers, self-focusing and self-modulation of trains in monochromatic wave in nonlinear optics, collapse of Langmuir waves in plasma physics as well as modeling of deep water freak waves in ocean [2]- [4].Reference [5], used the Schrodinger equation in modeling the 2004 Indonesian Tsunami while [6] used the Schrodinger equation in studying Bose-Einstein condensation using ultra cold neutral bosonic gases.Reference [7] used the Schrödinger wave equation in explaining the quantum behavior of a conventional wireless communication system.There exists therefore, the need for further researches seeking more accurate numerical approximations to the Schrödinger equation.
In this paper, we shall discretize the 1-dimensional Schrödinger equation by the Method of Lines (MOL).The numerical solution of the resulting system of coupled stiff Ordinary Differential Equations (ODEs) shall be sought for using an L-stable implicit trapezoidal-like integrator.This paper is organized as follows.In Section 2 we give a brief overview of the MOL.In Section 3 we discretize the 1-dimensional time dependent Schrödinger equation by the MOL under some assumptions.In Section 4 we discuss the stability property of the trapezoidallike integrator.In Section 5 we show how absolute and relative errors are obtained in the computation and in Section 6 we discuss the results obtained.

The Method of Lines
A well-known numerical approach for solving a time dependent PDE, whose solutions vary both in time and in space, is the Method of Lines.In the Method of Lines approach to the solution of PDEs, the solution of the PDE is sought for in a finite domain with a time coordinate t and a spatial coordinate, x.The spatial coordinate is defined as a discrete set of points given by j x j x = ∆ , and the boundaries correspond to the points 0 x and n x .Time, n t n t = ∆ is also defined only for certain values of the continuum time.Thus a function ( ) , u x t is defined only for the values of x and t that correspond to points ( ) in the mesh by ( ) , n j u x t .For a uniformly discretized domain, the resolution in the spatial and time coordinates respectively.In using the MOL to solve PDEs the spatial derivative is approximated by its finite difference approximation thereby, generating a coupled system of Ordinary Differential Equations (ODEs) in the time dependent variable t, which is associated with the initial value.Numerical approximations to the solutions of the PDE are obtained by marching forward in time on the grid.It is at this point that existing and generally well established numerical methods for initial value ODEs can be used to compute numerical solutions to the PDE.In [8], the applicability of the MOL in solving problems in Physics, fluid dynamics, rector models as well as automatic control is demonstrated.

Discretization of the Schrödinger Equation by the MOL
Given the Schrödinger equation where: ( ) According to [9], by letting ( ) and normalizing  and 2m to unity Equation (1), can be written in the form where q is arbitrary parameter.
In subscript notation Equation ( 2) can be written as Replacing the spatial derivative xx u by a second order finite difference approximation in Equation ( 3) and following [10], the function, ( ) , u x t , is replaced by the sum of its real part ( ) , v x t and its imaginary part ( ) Substituting Equation (4) into Equation ( 3), the two equations: are obtained.
Equating the real and imaginary parts of Equations ( 5a) and (5b) to zero, the following nonlinear system of PDE ( ) From Equation (6a) the equations: ( ) are obtained.
From Equation (6b), ( ) Replacing the second order partial derivative ( ) According to [11], the analytical solution of Equation ( 3) is To analyze ( ) , u x t as two real functions it is obtained that Equations ( 9a) and (9b) constitute the discretized form of the real and imaginary parts of Equation ( 3) by the MOL respectively.The real and the imaginary parts of the solution of the Schrodinger equation are independent, identically distributed standard normal random variables with mean 0 and variance 1 [12].Hence Equation (9a) or (9b) can be chosen for computational purposes.In this paper we shall use Equation (9a) and 0 q = , for the numerical integration.

L-Stable Implicit Trapezoidal-Like Integrator
The trapezoidal-like integrator shall be used in the numerical integration, where: ) 1 e e e 1 e e e e e l l l l and λ λ are the first and second eigenvalues of the discretization matrix respectively; l is the time step; and u′ denotes differentiation with respect to time.
The derivation of the method (12) and the analysis of the order of accuracy are as discussed in [13].

Stability Properties of the Method (12)
In this section, we state and prove one lemma and two theorems which relate to the stability property of the method (12) as follows: λ and 2 λ be the eigenvalues of the system of ODE evolved by discretizing the Schrödinger equation by the MOL; and let 1 0 λ < and 2 0 λ < .Suppose the following conditions hold:  (15) Therefore, ) Resolving inequality ( 16) into partial fractions, we have ) From Equations ( 17), (13a) and (13b), we have that 0 S R − > .This proves the Lemma 1.
Applying the test equation To the method in Equation ( 12) and from the general linear multistep method.The Schur's polynomial is thus given by ( ) ( ) ( ) But ( ) 0 S R − > as proved in Lemma 1, hence 0. q < The numerical integrator ( 12) is therefore A-stable.

Theorem 2
The numerical integrator is L-stable.
Substituting the scalar test Equation (18) into Equation ( 12) yields Hence ( ) Substituting for R and S in Equation ( 27), simplifying and taking the limit by use of L'Hospital's rule we obtain Hence, the proof of the Theorem 2. Thus the integrator ( 12) is L-stable.

Computation of Absolute and Relative Errors
In this section we explain how the absolute and relative errors of the methods shown in the Table 1 and Table 2 were obtained.

Absolute Errors
The absolute errors of the scheme were computed by the use of the formula: ( ) where the numerical solution at the grid point ( ) x t is ij u and the analytical solution at the same grid point is ( )

Relative Errors
Relative errors of the method were computed by use of the formula: where the numerical solution at the grid point ( ) x t is ij u and the analytical solution at the same grid point is ( ) , .
i j u x t

Results and Discussion
On the implementation of the L-stable implicit trapezoidal-like integrator for the solution of one-dimensional time dependent Schrödinger equation after discretizing with the MOL, the errors and relative errors of the scheme were computed as shown in Table 1 and Table 2 respectively.The L-stable property of the scheme is supported by the fact that the maximum relative errors at the various discretization levels are not significantly different as shown in Table 3.The performance of the trapezoidal-like integrator compares favorably with the      method of [10] as shown in Table 4.This indicates that the method is considered to give more accurate and stable results at the various spatial discretization levels considered.All computations and plotting of Figure 1 and Figure 2 in this paper were carried out using Maple 15.
7b) by the three point central difference ap-8b) where h is the space between the discretized line and n is the number of grids, the following equations are obtained.

∵
All the quantities on the L.H.S. of inequality (13) are positive.Expanding and rearranging the L.H.S. of (13), we have

Table 1 .
Absolute errors of the trapezoidal-like method at various spatial discretization levels.

Table 2 .
Relative errors of the trapezoidal-like method at various spatial discretization levels.

Table 3 .
Maximum relative errors of the trapezoidal-like integrator at various spatial discretization levels.

Table 4 .
[10]arison of the relative maximum errors of the trapezoidal-like integrator and[10]at various spatial discretization levels.