In this paper, we investigate interface condition for linear Schrödinger equation with non-uniform grid. Schrödinger equation is the fundamental equation of quantum mechanics, it describes the transformation of system's quantum state. Schrödinger equation is a wave function which describe a probability amplitude, the amplitude square corresponds to the measured probability distribution of electrons. The Cauchy problem of Schrödinger equation is

$\begin{eqnarray}\label{eq_Schr} && i\hbar\frac{\partial \Psi}{\partial t}=-\frac{\hbar^{2}}{(2m)}\frac{\partial^2 \Psi}{\partial \mathbf{x^2}}+V(\mathbf{x})\Psi, \mathbf{x}\in R^n, t >0, \\ && \Psi^{0}(\mathbf{x})=\phi(\mathbf{x}), \mathbf{x}\in R^n, \nonumber \end{eqnarray}$

(1.1)

where $\hbar$ is Planck constant, $m$ is quality of the particle, $\Psi(\mathbf{x}, t)$ is the wave function, $V(\mathbf{x})$ is potential energy.

Equation (1.1) has many important applications, such as in semiconductors [1].

The finite difference method and finite element method can be employed for solving Schrödinger equation. In this article, we use finite difference method. The numerical methods for Schrödinger equation simulation have numerous studies. Feit et al. studied the spectral method for Schrödinger equation [2]. Taha et al. made several comparison of numerical methods for nonlinear Schrödinger equation [3]. Bao et al. studied time-splitting spectral approximations for the linear Schrödinger equation [4].

To numerically solve the equation, we need to discrete the space and time domain. But Schrödinger wave equation in a local region maybe is singular or nearly singular, in this case it has great advantage to use non-uniform mesh. Because the limitation of numerical CFL condition, the maximum time step is limited by the minimum space step. If we use non-uniform mesh to solve the equation, this limitation will waste computational resources. To accurate the computation and speed up the algorithm, we use the local time step method, i.e. use two different time steps corresponding the fine space mesh and the coarse mesh. The problem is that how to deal with the interface between the fine mesh and coarse mesh which generated by the local time step method. In this article, for simplicity and better illustrating our basic idea, we just consider the non-unform mesh with two space step scales.

Our idea is enlightened by a multi-scale coupling method [5], base on the localization property of depended region of Schrödinger wave equation. We treat the fine mesh part and coarse mesh part respectively, the coarse mesh part is easy to calculate. The interface arise in the fine mesh part, so our main technique means focus on the fine mesh part. Because the linearity of the Schrödinger equation, we split the fine mesh part equation into two equations, the initial value of one equation is the constant value of nearest coarse mesh point. The analysis solution can be derived from D'Alembert's principle. The initial value of other equation is the original equation minus the constant, at the short time scale, we can consider the equation is no source outside, then we treat this equation with artificial boundary condition.

The rest of this paper is organized as follows. In Section 2 we briefly recall several artificial boundary conditions for Schrödinger equation. The interface condition for local time step and non-uniform mesh is given in Section 3. In Section 4, several numerical examples are given to present the efficiency of the interface condition we proposed.

The artificial boundary conditions of Schrödinger equation have fruitful results which were collected in the review article wrote by Antoine [6] et al.

We concern the one-dimensional linear Schrödinger equation with zero potential function

$\begin{eqnarray}\label{eq:Schrsimple} && i\frac{\partial \Psi}{\partial t}=\frac{\partial^2 \Psi}{\partial x^2}, x\in R, t >0, \\ && \Psi^{0}(x)=\phi(x), x\in R.\nonumber \end{eqnarray}$

(2.1)

There are two basic kinds of artificial boundary conditions, transparent boundary conditions and perfect matched layer. But transparent boundary conditions are non-local in time and space, it need to do some local approximation, next we will introduce some absorbing boundary conditions for the linear Schrödinger equation which are local both in space and time.

Szeftel [7] developed the strategy of pseudo-differential calculus for the linear Schrödinger equation and constructed a family of absorbing boundary conditions for the linear Schrödinger equation on curved boundaries which are local both in space and time as below:

$\begin{eqnarray} && \partial_{x}\Psi|_{x=x_{N}}=i\beta \Psi_{N-1}+i\sum\limits_{k=1}^{m}a_{k}(\Psi_{N-1}-d_{k}\varphi_{k}), \\ && i\partial_t\varphi_{k}+d_{k}\varphi_{k}=\Psi(x_N, t), k=1, \cdots, m, \nonumber \end{eqnarray}$

(2.2)

where $x_{N}$ is the value of x in boundary and $\beta, a_1\cdots a_k, d_1\cdots d_k$ are the reflection coefficients.%

Next, we introduce the PMLs for Schrödinger equation, also consider the time-dependent Schrödinger equation in one dimension:

$\begin{equation*} i*\frac{\partial \Psi}{\partial t}=-\frac{\partial^2 \Psi}{\partial x^2}, (x, t)\in [0, x_{b}]\times R^+. \end{equation*}$

We want to construct a PMLs with layer width $d$ at $ x = x_{b}$ in the $x$ direction, then the new computation domain is $[0, x_{b} + d]$. Next, expand $\Psi$ in a Fourier series donate by $\psi_{j}$. Then Laplace transform the equation in time, then after some technical means, we can obtain decaying waves, we solve (2.3) instead of Schrödinger equation:

$\begin{equation}\label{pml} i\frac{\partial \Psi}{\partial t}=-\frac{1}{1+e^{i\gamma}\sigma(x)}\frac{\partial}{\partial x}\frac{1}{1+e^{i\gamma}\sigma(x)}\frac{\partial \Psi}{\partial x}, (x, t)\in [0, x_b+d]\times R^+, \end{equation}$

(2.3)

where

$\begin{equation*} \sigma(x)=\left\{\begin{array}{ll} 0, &x\in [0, x_b], \\ {\rm constant},&x\in[x_b, x_b+d]. \end{array}\right. \end{equation*}$

We can see that in the interval $x\in [0, x_b]$, equation (2.3) is same as Schrödinger equation. In the interval $[x_b, x_b+d]$, equation (2.3) have the exponential decay property. So the PMLs have this good property.

The numerical simulation of the Schrödinger equation need discrete the space and time domain. In most researches, they use uniform mesh in space discretion, and use proper time step to match the space mesh. It has a problem that how to choose the proper time step when using non-uniform space mesh to discrete space domain, just for the numerical method of ordinary wave equation, the non-uniform mesh will cause many problems, such as spurious reflection during simulation [8, 9].

To avoid the spurious reflection, we use local time step for the non-uniform mesh, this will generate an interface between fine grid and coarse grid. So we need introduce a proper interface condition for the numerical scheme. For ordinary wave equation, nonlinear conservation laws are approximated with an explicit time difference method [10], Tan et al. [11] use the moving mesh methods with local time stepping for solving time dependent PDEs. The interface condition we proposed is utilizing the property of wave equation's local dependance and the artificial boundary condition of wave equation. We use 1 dimensional Schrödinger equation as example to present the interface condition:

$\begin{eqnarray}\label{eq:1DSchrsimple} && i\frac{\partial \Psi}{\partial t}=\frac{\partial^2 \Psi}{\partial x^2}, (x, t)\in[a, b]\times R^+, \\ && \Psi^{0}(x)=\phi(x), x\in[a, b].\nonumber \end{eqnarray}$

(3.1)

The space domain we concerned is $x\in[a, b]$, and we set the interface between the fine grid and coarse grid is $x=x_{i}$, then the coarse mesh is $[a, x_{i}]$, the fine mesh is $[x_{i}, b]$. Assume that we have gotten the numerical solution of the problem at $t=t^{n}$, we need to numerically solve the equation at $t=t^{n+1}$.

Here, for simplicity, we consider the right side as fine mesh and coarse mesh on left side. We use a larger time step in coarse mesh, $\vartriangle t_{n}^{c}=t^{n+1}-t^{n}$, and smaller time step in fine mesh, $\vartriangle t_{n}^{f}=\vartriangle t_{n}^{c}/(m+1)$, $x=x_{i}$ is between fine mesh and coarse mesh, see Figure 1.

We deal with the coarse mesh and fine mesh respectively. First, for the coarse mesh, we solve the problem with the matched time step, then use the finite difference method to solve equation (3.2)

$\begin{eqnarray}\label{eq_coarse} && i*\frac{\partial \Psi}{\partial t}= -\frac{\partial^2 \Psi}{\partial x^2}, \ x\in[a, x_{i}], \\ && \Psi^{0}(x)= \Psi(x, t_{n}), \ x\in[a, x_{i}].\nonumber \end{eqnarray}$

(3.2)

And in fine mesh domain, we decompose the equation into two equations. The initial condition of first equation (3.3) is the original equation minus the constant value at interface point $x_{i}$,

$\begin{eqnarray}\label{eq_artifiacial} && i*\frac{\partial \Psi}{\partial t}= -\frac{\partial^2 \Psi}{\partial x^2}, \ x\in[x_{i}, b], \\ && \Psi^{0}(x)= \Psi(x, t_{n})-\Psi(x_{i}, t_{n}), \ x\in[x_{i}, b].\nonumber \end{eqnarray}$

(3.3)

The initial condition of second equation (3.4) is constant value at interface point $x_{i}$.

$\begin{eqnarray}\label{eq_constant} && i*\frac{\partial \Psi}{\partial t}=-\frac{\partial^2 \Psi}{\partial x^2}, \ x\in[x_{i}, b], \\ && \Psi^{0}(x)=\Psi(x_{i}, t_{n}), x\in[x_{i}, b].\nonumber \end{eqnarray}$

(3.4)

For equation (3.3), because of the local dependence of wave equation, the equation can be treated as a no source equation, so we can use the artificial boundary condition as boundary condition. Equation (3.4) have the constant initial value, it is easy to get its solution. At last, in consequence of the linearity of Schrödinger equation, we add up two equation solutions to get the solution on the fine mesh.

In our scheme, the complexity of scheme is not changed too much. Then we can adopt many artificial boundary conditions to our interface condition.

In this section, we will give several numerical examples to illustrate the efficiency of the interface condition. The numerical example settings are as follow, the space domain is $x\in$[-9, 5], the interface is at $x=-2$, The interval [-9, -2] is discrete as coarse mesh, the fine mesh is in the interval [-2, 5]. We give two different initial conditions corresponding to the wave propagate from the coarse grid to fine grid and fine grid to coarse grid.

$\begin{align*} \Psi_1(x, 0)&=\frac{\exp(-i\pi/4)}{\sqrt -i}\exp(\frac{ix^2-6x}{-i}), \\[1.5ex] \Psi_2(x, 0)&=\exp(\frac{-i(x+5)^2-5(x+5)}{i}). \end{align*}$

Accordingly, the real solutions are

$\begin{align*} \Psi_1(x, t)&=\frac{\exp(-i\pi/4)}{\sqrt{4t-i}}\exp(\frac{ix^2-6x-36t}{4t-i}), \\[1.5ex] \Psi_2(x, t)&=\sqrt{\frac{i}{-4t+i}}\exp(\frac{-i(x+5)^2-5(x+5)+25t}{-4t+i}). \end{align*}$

In numerical simulation, we use the artificial boundary condition as the boundary condition at boundary point $x=-9, x=5$. The coarse mesh size is $h_c$=0.04 and fine mesh size is $h_f$=0.02. First, we will give the exact solutions in Figures 2 and 3.

Then we will give an example to present the spurious reflection caused by the non-uniform space mesh.

Example 4.1  In this example, we use the non-uniform space mesh to discrete the space domain. We set time step corresponding to the minimal space mesh step $h_f$.

The numerical result of Example 4.1 are present in Figures 4 and 5. The $L^2$ error show in Figure 6. There is a huge spurious reflection after the wave propagated the interface between the fine and coarse mesh. The spurious reflection is harmful to the convergence and stability of the numerical simulation. We can use our interface condition to decrease this spurious reflection.

The next example we will use local time step and interface condition, the example shows the interface condition significantly decrease the spurious reflection when the wave propagate cross the interface, and the error of numerical solution is much smaller than uniform time step.

Example 4.2 In this example, we use Szeftel boundary condition as the interface condition. We set the local time step corresponding to fine mesh and coarse mesh respectively.

The results of Example 4.2 present in Figures 7 and 8, and the error norm is in Figure 9, we can see from the figures, the interface condition decrease the spurious reflection significantly.

The interface condition can also use other kind of artificial boundary condition, such as PMLs. We use the PMLs as the boundary condition and interface condition.

Example 4.3 In this example, we use PMLs boundary condition as the interface condition. We set the local time step corresponding to fine mesh and coarse mesh respectively.

The numerical simulation of Example 4.3 is showing in Figures 10 and 11. The interface condition also can decrease the reflection when the wave propagate cross the interface.

We have proposed a new interface condition for local time step simulation of linear Schrödinger's equations. Our method effectively reduce the reflections in interface, and it gives a good numerical accuracy with a low numerical cost. We can use all kinds of boundary conditions to match this method and the accuracy is decided by the ability of absorbing reflection of boundary conditions. In two dimensional case, we could still reduce the reflection on interface by similar method. It is possible to extend the interface condition to the non-linear Schrödinger's equations.