Eigenvalues & Eigenfunctions of the 1-d Schrödinger Equation

nsolve.f90 is a Fortran 90 module solving the 1-d Schrödinger equation (especially the radial equation for 3-d problems) using the Numerov method.

Introduction

nsolve.f90 finds the solutions \( (f, \varepsilon) \) of the 1-d Schrödinger equation, \[ -\frac{1}{2} \frac{d^2 f}{d x^2} + v(x) f(x) = \varepsilon f(x)\,, \] subject to the boundary conditions \[ \quad f(a) = f_a, \quad f(b) = f_b\,, \] where \( f_a \) and \( f_b \) are constants. The function \( f \) is assumed to be defined on the closed interval \( [a,b] \).

For example, the radial equation for a particle in an isotropoic harmonic potential is \[ -\frac{\hbar^2}{2m} \frac{d^2 u}{d r^2} + \left[ \frac{\hbar^2 l (l+1)}{2 m r^2} + \frac{1}{2} m \omega^2 r^2 \right] u(r) = E u(r) \,,\] where \( u(r) = r R(r) \) and \( R(r) \) is the radial part of the wavefunction. When lengths are expressed in units of the oscillator length \( d=\sqrt{\hbar/m\omega} \) and energies are expressed in units of \(\hbar \omega\), the equation to be solved becomes \[ -\frac{1}{2} \frac{d^2 u}{d r^2} + v(r) u(r) = \varepsilon u(r)\,, \] where \( v(r) = l(l+1)/2 r^2 + r^2/2 \). The domain of this equation is \( [0,\infty) \), where the point at 0 is determined by the condition \( u(0) = 0 \).

To solve, we truncate the domain to \( [0,b] \), where \( b \) is large relative to the oscillator length (e.g. \( b = 10 d \)), and then apply the Numerov method to integrate the Schrödinger equation, along with the shooting method to find eigenvalues and eigenfunctions.

Download

Notes

To Do

  1. Use compensated summation to reduce rounding error in the Numerov integration.

References

  1. S. E. Koonin and D. C. Meredith, Computational Physics, Fortran Version, Addison-Wesley (1990).

Note that the Fortran code linked from the website is not based on the code given in Ref. 1.