Solving Schrodinger equation using numerical methods
The code for this article can be found in the Colab notebook.
In this article, we will solve the time-independent Schrodinger equation using Numerov’s algorithm. Let’s get into it.
In 1D, the time-independent Schrodinger equation looks like
which can be rewritten as
where
and we ignored Plank’s constant and mass by assuming they are both equal to 1. (This only changes the scale of the solution).
We will discretize the coordinate x using a uniform grid; x_k = k∆x. Now, using Taylor expansion, the second equation above can be rewritten as the following recurrence relation (This is known as the Numerov’s method. See this doc for derivation, section 2.1)
Infinite well
We will try to solve this for the infinite well where potential is
We have 2 unknowns in the recurrence relation above: energy and wavefunction. We define a function that will solve for the wavefunction given a value of the energy:
Let’s find the wavefunction for some values of E (energy):
Let’s plot psi:
Do the same for E=1.5
We know that the wavefunction must vanish at x=1.0, but psi(1.0)>0 for E=1.0 and psi(1.0)<0 for E=1.5. So the actual value for E must be between 1.0 and 1.5. Let’s use the Shooting method to find the exact value of E. The method is very similar to the Binary search algorithm — we narrow down the interval around the actual energy value by constantly dividing up the interval in half. (See this and this for the description of the Shooting method)
We can also plot the intermediate wavefunctions which were used during the binary search
This is just the first (ground) state of the infinite well. There are other energy states which solve the Schrodinger equation:
We saw that wavefunction satisfies the right boundary condition (psi(1.0)=0) only for certain values of the energy. This is exactly where the discretization comes from — Schrodinger equation is solved only for certain energies and wavefunctions. There is an interactive version of the above graphs in the Colab notebook if you want to see the discretization in action.
Conclusion
We can get some intuition about why wavefunction looks like sine waves by rearranging the recurrence relation. In the region where V=0,
And the recurrence relation above reduces to
Take the ground state for example, where E=1.2337. Plugging it in,
We see that the middle value is slightly greater than the average of the neighboring values. You can see how that can lead to the negative second derivative. The figure below shows a few values of the wavefunction. The actual values are in blue and the average value of the neighboring values is in red
Unfortunately, this method breaks for other “soft” potentials such as the finite well and the harmonic oscillator. A few improvements to the method such as using the asymptotic form (described here, section 2.2) of the potential can help.
