For a one-particle central-force problem, the wave function is given by (6.16) as $\psi=R(r) Y_{l}^{m}(\theta, \phi)$ and the radial factor $R(r)$ is found by solving the radial equation (6.17). The Numerov method of Section 4.4 applies to differential equations of the form $\psi^{\prime \prime}=G(x) \psi(x)$ [Eq. (4.66)], so we need to eliminate the first derivative $R^{\prime}$ in (6.17). Let us define $F(r)$ by $F(r) \equiv r R(r)$, so
\(
\begin{equation}
R(r)=r^{-1} F(r) \tag{6.136}
\end{equation}
\)
Then $R^{\prime}=-r^{-2} F+r^{-1} F^{\prime}$ and $R^{\prime \prime}=2 r^{-3} F-2 r^{-2} F^{\prime}+r^{-1} F^{\prime \prime}$. Substitution in (6.17) transforms the radial equation to
\(
\begin{gather}
-\frac{\hbar^{2}}{2 m} F^{\prime \prime}(r)+\left[V(r)+\frac{l(l+1) \hbar^{2}}{2 m r^{2}}\right] F(r)=E F(r) \tag{6.137}\
F^{\prime \prime}(r)=G(r) F(r), \quad \text { where } G(r) \equiv \frac{m}{\hbar^{2}}(2 V-2 E)+\frac{l(l+1)}{r^{2}} \tag{6.138}
\end{gather}
\)
which has the form needed for the Numerov method. In solving (6.137) numerically, one deals separately with each value of $l$. Equation (6.137) resembles the one-dimensional Schrödinger equation $-\left(\hbar^{2} / 2 m\right) \psi^{\prime \prime}(x)+V(x) \psi(x)=E \psi(x)$, except that $r$ (whose range is 0 to $\infty$ ) replaces $x$ (whose range is $-\infty$ to $\infty$ ), $F(r) \equiv r R(r)$ replaces $\psi$, and $V(r)+l(l+1) \hbar^{2} / 2 m r^{2}$ replaces $V(x)$. We can expect that for each value of $l$, the lowest-energy solution will have 0 interior nodes (that is, nodes with $0<r<\infty$ ), the next lowest will have 1 interior node, and so on.
Recall from the discussion after (6.81) that if $R(r)$ behaves as $1 / r^{b}$ near the origin, then if $b>1, R(r)$ is not quadratically integrable; also, the value $b=1$ is not allowed, as noted after (6.83). Hence $F(r) \equiv r R(r)$ must be zero at $r=0$.
For $l \neq 0, G(r)$ in (6.138) is infinite at $r=0$, which upsets most computers. To avoid this problem, one starts the solution at an extremely small value of $r$ (for example, $10^{-15}$ for the dimensionless $r_{r}$ ) and approximates $F(r)$ as zero at this point.
As an example, we shall use the Numerov method to solve for the lowest boundstate H -atom energies. Here, $V=-e^{2} / 4 \pi \varepsilon{0} r=-e^{\prime 2} / r$, where $e^{\prime} \equiv e /\left(4 \pi \varepsilon{0}\right)^{1 / 2}$. The radial equation (6.62) contains the three constants $e^{\prime}, \mu$, and $\hbar$, where $e^{\prime} \equiv e /\left(4 \pi \varepsilon_{0}\right)^{1 / 2}$ has SI units of $\mathrm{m} \mathrm{N}^{1 / 2}$ (see Table A. 1 of the Appendix) and hence has the dimensions $\left[e^{\prime}\right]=\mathrm{L}^{3 / 2} \mathrm{M}^{1 / 2} \mathrm{~T}^{-1}$. Following the procedure used to derive Eq. (4.73), we find the H -atom reduced energy and reduced radial coordinate to be (Prob. 6.47)
\(
\begin{equation}
E{r}=E / \mu e^{\prime 4} \hbar^{-2}, \quad r{r}=r / B=r / \hbar^{2} \mu^{-1} e^{\prime-2} \tag{6.139}
\end{equation}
\)
Use of (6.139) and (4.76) and (4.77) with $\psi$ replaced by $F$ and $B=\hbar^{2} \mu^{-1} e^{\prime-2}$ transforms (6.137) for the H atom to (Prob. 6.47)
\(
\begin{equation}
F{r}^{\prime \prime}=G{r} F{r}, \quad \text { where } G{r}=l(l+1) / r{r}^{2}-2 / r{r}-2 E_{r} \tag{6.140}
\end{equation}
\)
and where $F{r}=F / B^{-1 / 2}$.
The bound-state H -atom energies are all less than zero. Suppose we want to find the H -atom bound-state eigenvalues with $E{r} \leq-0.04$. Equating this energy to $V{r}$, we have (Prob. 6.47) $-0.04=-1 / r{r}$ and the classically allowed region for this energy value extends from $r{r}=0$ to $r{r}=25$. Going two units into the classically forbidden region, we take $r{r, \max }=27$ and require that $F{r}(27)=0$. We shall take $s{r}=0.1$, giving 270 points from 0 to 27 (more precisely, from $10^{-15}$ to $27+10^{-15}$ ).
$G{r}$ in (6.140) contains the parameter $l$, so the program of Table 4.1 has to be modified to input the value of $l$. When setting up a spreadsheet, enter the $l$ value in some cell and refer to this cell when you type the formula for cell B7 (Fig. 4.9) that defines $G{r}$. Start column A at $r{r}=1 \times 10^{-15}$. Column C of the spreadsheet will contain $F{r}$ values instead of $\psi{r}$ values, and $F{r}$ will differ negligibly from zero at $r{r}=1 \times 10^{-15}$, and will be taken as zero at this point.
With these choices, we find (Prob. 6.48a) the lowest three H -atom eigenvalues for $l=0$ to be $E{r}=-0.4970,-0.1246$, and -0.05499 ; the lowest two $l=1$ eigenvalues found are -0.1250 and -0.05526 . The true values [Eqs. (6.94) and (6.139)] are $-0.5000,-0.1250$, and -0.05555 . The mediocre accuracy can be attributed mainly to the rapid variation of $G(r)$ near $r=0$. If $s{r}$ is taken as 0.025 instead of 0.1 (giving 1080 points), the $l=0$ eigenvalues are improved to $-0.4998,-0.12497$, and -0.05510 . See also Prob. 6.48b.