A numerical solution of the time-independent Schroedinger equation for hydrogenic atoms with reduced mass m@m_e.

The program Hatom1.txt is written in QuickBasic.  It demonstrates a method for solving the radial equation for a hydrogenic atom.  This is a boundary value problem.  We solve for the eigenvalues E of the Hamiltonian for which the wavefunction and its derivative are continuous and for which the wavefunction takes on given values at the boundaries of the integration volume.  The program plots the radial functions R_kl(r) corresponding to the calculated eigenvalues.  The nuclear charge Z and the angular momentum quantum number l can easily be varied.

Explanation:

The wavefunction y(r) is of the form We solve the radial equation for u_kl(r).

We integrate this equation using the Numerov method.  (See Program 3, QM1.)   For each l we find all eigenvalues E_kl<0.5Z²eV.  We denote the successive energy levels by k=1, 2, 3, ... .  We expect that E_kl depends only on n=k+l.

We assume that u_kl(r) is zero at r=0 and at r=r_max=6(l+1)Z².

As r goes to zero the effective potential varies rapidly, while for large r it is close to zero and varies slowly.  Integration in regions where k varies rapidly requires small steps, while larger step sizes can be used when k varies slowly.  Using a small step size everywhere increases CPU time and numerical inaccuracies, while using a large step size everywhere will not produce an acceptable result.  We divide the region between r=0 and r=r_max into two intervals to demonstrate how variable step sizes can be used.  In the first interval the step size is 100 times smaller than in the second interval.  Ideally we would like to change the step size continuously.  We will demonstrate how this can be done in  another program.

We start integrating at r=dr.  (At r=0  V_eff blows up resulting in a "divide by zero" error.)  We choose a starting energy lower than the expected ground state energy.  We integrate outward to r_max, counting the number of nodes in the wavefunction as a function of energy, as we increase the energy in small steps DE.  When the number of nodes changes from one energy to the next the wavefunction is close to zero at r=r_max and we are close to an eigenvalue.  We decrease DE once and then take as the eigenvalue the average value of the energy just before and just after the increase in the number of nodes.  This interpolation can easily be refined, for example by using a secant search.  (See Program 3, QM1.)

The program Hatom1.txt plots the radial functions R_kl(r) from r=0 to r=50 Å corresponding to the calculated eigenvalues.

Link:

• The Hydrogen Atom