The Leapfrog Method

Last time we saw that the naive implementation of the finite difference approximation is unstable. We maintain that the central difference approximation yields better results than the forward difference approximation at the cost of additional memory.

The central difference approximates the time derivative using the difference between the previous and the next timestep.

t Ψ x t Ψ x t + Δ t - Ψ x t - Δ t 2 Δ t

We can reuse the same logic we applied to the Schrödinger equation the first time, to lead us to the an implementation of this modified approach. This typifies how understanding general principles and concepts is far more valuable than focusing on specific recipes. While some of the details are different, the concepts and thought processes that we employ now are very much the same as the first time. Also seeing the process again, with a slightly different application, will strengthen your understanding.

Substituting this approximation into the Schrödinger equation and isolating the t + Δ t term on the left gives our guiding equation for evolving the wave function forward in time.

i Ψ x t + Δ t i Ψ x t - Δ t - 2 x 2 Ψ x t Δ t + V x Ψ x t 2 Δ t

Remember that the terms on the right all became texture reads in the code. This equation shows that we will now need to read, and therefore keep, wavefunction values for t - Δ t as well as t . The wave function for t + Δ t becomes the wavefunction from t - Δ t plus adjustments from t . This gives raise to the name of the method, leapfrog.

The Leapfrog method adds Ψ x t - Δ t and hence another texture to our calculations.

Once again we have to cast the equation in a form that can be manipulated in JavaScript and GLSL shaders, which lack an explicit complex data type. Start by breaking the wave function into real and imaginary parts.

Ψ x t = Re ( Ψ x t ) + i Im ( Ψ x t )

We put this into the Schrödinger equation, then gather the real and imaginary parts into their own equations.

t Re ( Ψ x t ) = - 1 2 2 x 2 Im ( Ψ x t ) + V x Im ( Ψ x t ) t Im ( Ψ x t ) = 1 2 2 x 2 Re ( Ψ x t ) - V x Re ( Ψ x t )

We reuse our earlier substitution for 2 x 2 Ψ x t , which the clever reader will note uses elements of the wave function from x + Δ x and x - Δ x , and is already a central difference expression.

2 x 2 Ψ x t Ψ x + Δ x t - 2 Ψ x t + Ψ x - Δ x t Δ x 2

Now we have approximations of both the time and space derivatives in terms of values on our grid. We use these in the Schrödinger and write out equations for the evolution of the real and imaginary parts of the wave function.

For the real part of the wave function we have

Re ( Ψ x t + Δ t ) Re ( Ψ x t - Δ t ) - Im ( Ψ x + Δ x t ) - 2 Im ( Ψ x t ) + Im ( Ψ x - Δ x t ) Δ x 2 Δ t + 2 V x Im ( Ψ x t ) Δ t

and to evolve the imaginary component of the wave function in time

Im ( Ψ x t + Δ t ) Im ( Ψ x t - Δ t ) + Re ( Ψ x + Δ x t ) - 2 Re ( Ψ x t ) + Re ( Ψ x - Δ x t ) Δ x 2 Δ t - 2 V x Re ( Ψ x t ) Δ t

These initial conditions exploded using the forward difference in time. The leapfrog method is far more stable.

Creative Commons License
This work is licensed under a Creative Commons Attribution 4.0 International License.