Wednesday, September 15, 2021

Solving the Schrödinger Equation with Numerov’s Algorithm

 

The Schrödinger equation describes the energy and time-evolution of a particle or system of particles, and is one of the fundamental building blocks of modern physics. In it’s general form, the (time-independent) Schrödinger equation looks like this:1

\(\frac{-\hbar^2}{2m} \frac{\partial^2}{\partial x^2}\psi(x) + V(x) \psi(x) = E\psi(x)\)

There are relatively few situations in which the Schrödinger equation can be solved analytically, and numerical methods and approximations are one way around that analytical limitation. To demonstrate how this is possible and how a numerical solution works, what better way than to solve a system which can be solved analytically and comparing the results.

In solving the Schrödinger equation, we will start with one of the simplest interesting quantum mechanical systems, the quantum mechanical harmonic oscillator.2 Let’s first define our quantum harmonic oscillator. The general form of the Schrödinger equation for a one-dimensional harmonic oscillator reads thus:

22m2z2ψ(z)+mz22ψ(z)=Eψ(z)(1)

We will make use of the Numerov algorithm which is particularly suited to solving second order differential equations of the form y(x)+k(x)y(x)=0

. You can read more about it elsewhere, including its derivation etc., but its form is quite simple and easy to code:

(1+112h2kn+1)yn+1=2(1512h2kn)yn(1+112h2kn1)yn1+O(h6)(2)

As you can see, it provides 6th order accuracy which is pretty impressive for such a simple algorithm. In the above equation, h

is the step size between each iteration, and n is the index of iteration; y and k

relate to those in the formula in the paragraph above.

Thus we need to manipulate (1)

into a (dimensionless) form which the Numerov algorithm can solve: using a substitution E=εω and z=xmω we can rearrange (1)

into the form:

ψ(x)+(2εx2)ψ(x)=0(3)

Now the Schrödinger equation is in the correct form, we can simply plug it into the Numerov algorithm and see what comes out.

Finding the Eigenvalues Numerically

To determine the eigenvalues, the program scans a range of energies, denoted by the Greek letter ε

in the above equations, and tests for when the tail of the graph changes from + to

or vice versa. When that happens, the tail must have crossed zero, and therefore it must have stepped over a solution.3 The program then goes backwards and so on with increased resolution, honing in until it finds all of the solutions we want.

Given the substitution above, we should expect the eigenvalues to be ε=n+12where n  

is an integer from zero (representing the ground state) upwards.4 Hit F12 to pull up the web console before you run the simulation to see what results you actually get.


blob:null/4730e6c2-37f3-4795-a528-64d20bd0306f


Note: On a smartphone? Possibly(?) wait till you’re on a more powerful machine (unless you’re hard; let me know in the comments how it worked out). It works fine! Likewise, using an old Internet Explorer? – well, don’t blame me if it crashes or just nothing happens. ;)

 

 https://mtdevans.com/2013/07/solving-the-schrodinger-equation-with-numerovs-algorithm/index.html

 


No comments:

Post a Comment