Spherically symmetric solutions
Let’s first try to find some very simple function that satisfies the horrible equation in (19.7). Although the wave function $\psi$ will, in general, depend on the angles $\theta$ and $\phi$ as well as on the radius $r$, we can see whether there might be a special situation in which $\psi$ does not depend on the angles. For a wave function that doesn’t depend on the angles, none of the amplitudes will change in any way if you rotate the coordinate system. That means that all of the components of the angular momentum are zero. Such a $\psi$ must correspond to a state whose total angular momentum is zero. (Actually, it is only the orbital angular momentum which is zero because we still have the spin of the electron, but we are ignoring that part.) A state with zero orbital angular momentum is called by a special name. It is called an “$s$-state”—you can remember “$s$ for spherically symmetric.”2
 Now if $\psi$ is not going to depend on $\theta$ and $\phi$ then the entire Laplacian contains only the first term and Eq. (19.7) becomes much simpler: \begin{equation} \label{Eq:III:19:8} \frac{1}{r}\,\frac{d^2}{dr^2}\,(r\psi)= -\frac{2m}{\hbar^2}\biggl(E+\frac{e^2}{r}\biggr)\psi. \end{equation} Before you start to work on solving an equation like this, it’s a good idea to get rid of all excess constants like $e^2$, $m$, and $\hbar$, by making some scale changes. Then the algebra will be easier. If we make the following substitutions: \begin{equation} \label{Eq:III:19:9} r=\frac{\hbar^2}{me^2}\,\rho, \end{equation} and \begin{equation} \label{Eq:III:19:10} E=\frac{me^4}{2\hbar^2}\,\epsilon, \end{equation} then Eq. (19.8) becomes (after multiplying through by $\rho$) \begin{equation} \label{Eq:III:19:11} \frac{d^2(\rho\psi)}{d\rho^2}= -\biggl(\epsilon+\frac{2}{\rho}\biggr)\rho\psi. \end{equation} These scale changes mean that we are measuring the distance $r$ and energy $E$ as multiples of “natural” atomic units. That is, $\rho=r/r_B$, where $r_B=\hbar^2/me^2$, is called the “Bohr radius” and is about $0.528$ angstroms. Similarly, $\epsilon=E/E_R$, with $E_R=me^4/2\hbar^2$. This energy is called the “Rydberg” and is about $13.6$ electron volts.  Since the product $\rho\psi$ appears on both sides, it is convenient to work with it rather than with $\psi$ itself. Letting \begin{equation} \label{Eq:III:19:12} \rho\psi=f, \end{equation} we have the more simple-looking equation \begin{equation} \label{Eq:III:19:13} \frac{d^2f}{d\rho^2}=-\biggl(\epsilon+\frac{2}{\rho}\biggr)f. \end{equation}  Now we have to find some function $f$ which satisfies Eq. (19.13)—in other words, we just have to solve a differential equation. Unfortunately, there is no very useful, general method for solving any given differential equation. You just have to fiddle around. Our equation is not easy, but people have found that it can be solved by the following procedure. First, you replace $f$, which is some function of $\rho$, by a product of two functions \begin{equation} \label{Eq:III:19:14} f(\rho)=e^{-\alpha\rho}g(\rho). \end{equation} This just means that you are factoring $e^{-\alpha\rho}$ out of $f(\rho)$. You can certainly do that for any $f(\rho)$ at all. This just shifts our problem to finding the right function $g(\rho)$.  Sticking (19.14) into (19.13), we get the following equation for $g$: \begin{equation} \label{Eq:III:19:15} \frac{d^2g}{d\rho^2}-2\alpha\,\ddt{g}{\rho}+\biggl( \frac{2}{\rho}+\epsilon+\alpha^2\biggr)g=0. \end{equation} Since we are free to choose $\alpha$, let’s make \begin{equation} \label{Eq:III:19:16} \alpha^2=-\epsilon, \end{equation} and get \begin{equation} \label{Eq:III:19:17} \frac{d^2g}{d\rho^2}-2\alpha\,\ddt{g}{\rho}+ \frac{2}{\rho}\,g=0. \end{equation}  You may think we are no better off than we were at Eq. (19.13), but the happy thing about our new equation is that it can be solved easily in terms of a power series in $\rho$. (It is possible, in principle, to solve (19.13) that way too, but it is much harder.) We are saying that Eq. (19.17) can be satisfied by some $g(\rho)$ which can be written as a series, \begin{equation} \label{Eq:III:19:18} g(\rho)=\sum_{k=1}^\infty a_k\rho^k, \end{equation} in which the $a_k$ are constant coefficients. Now all we have to do is find a suitable infinite set of coefficients! Let’s check to see that such a solution will work. The first derivative of this $g(\rho)$ is \begin{equation*} \ddt{g}{\rho}=\sum_{k=1}^\infty a_kk\rho^{k-1}, \end{equation*} and the second derivative is \begin{equation*} \frac{d^2g}{d\rho^2}=\sum_{k=1}^\infty a_kk(k-1)\rho^{k-2}. \end{equation*} Using these expressions in (19.17) we have \begin{equation} \label{Eq:III:19:19} \sum_{k=1}^\infty \!k(k\!-\!1)a_k\rho^{k\!-\!2}\!-\!\! \sum_{k=1}^\infty \!2\alpha ka_k\rho^{k\!-\!1}\!+\!\! \sum_{k=1}^\infty \!2a_k\rho^{k\!-\!1}\!= 0. \end{equation} It’s not obvious that we have succeeded; but we forge onward. It will all look better if we replace the first sum by an equivalent. Since the first term of the sum is zero, we can replace each $k$ by $k+1$ without changing anything in the infinite series; with this change the first sum can equally well be written as \begin{equation*} \sum_{k=1}^\infty(k+1)ka_{k+1}\rho^{k-1}. \end{equation*} Now we can put all the sums together to get \begin{equation} \label{Eq:III:19:20} \sum_{k=1}^\infty[(k+1)ka_{k+1}-2\alpha ka_k+2a_k]\rho^{k-1}=0. \end{equation}  This power series must vanish for all possible values of $\rho$. It can do that only if the coefficient of each power of $\rho$ is separately zero. We will have a solution for the hydrogen atom if we can find a set $a_k$ for which \begin{equation} \label{Eq:III:19:21} (k+1)ka_{k+1}-2(\alpha k-1)a_k=0 \end{equation} for all $k\geq1$. That is certainly easy to arrange. Pick any $a_1$ you like. Then generate all of the other coefficients from \begin{equation} \label{Eq:III:19:22} a_{k+1}=\frac{2(\alpha k-1)}{k(k+1)}\,a_k. \end{equation} With this you will get $a_2$, $a_3$, $a_4$, and so on, and each pair will certainly satisfy (19.21). We get a series for $g(\rho)$ which satisfies (19.17). With it we can make a $\psi$, that satisfies Schrödinger’s equation. Notice that the solutions depend on the assumed energy (through $\alpha$), but for each value of $\epsilon$, there is a corresponding series.  We have a solution, but what does it represent physically? We can get an idea by seeing what happens far from the proton—for large values of $\rho$. Out there, the high-order terms of the series are the most important, so we should look at what happens for large $k$. When $k\gg1$, Eq. (19.22) is approximately the same as \begin{equation*} a_{k+1}=\frac{2\alpha}{k}\,a_k, \end{equation*} which means that \begin{equation} \label{Eq:III:19:23} a_{k+1}\approx\frac{(2\alpha)^k}{k!}. \end{equation} But these are just the coefficients of the series for $e^{+2\alpha\rho}$. The function of $g$ is a rapidly increasing exponential. Even coupled with $e^{-\alpha\rho}$ to produce $f(\rho)$—see Eq. (19.14)—it still gives a solution for $f(\rho)$ which goes like $e^{\alpha\rho}$ for large $\rho$. We have found a mathematical solution but not a physical one. It represents a situation in which the electron is least likely to be near the proton!  It is always more likely to be found at a very large radius $\rho$. A wave function for a bound electron must go to zero for large $\rho$.  We have to think whether there is some way to beat the game, and there is. Observe! If it just happened by luck that $\alpha$ were equal to $1/n$, where $n$ is any positive integer, then Eq. (19.22) would make $a_{n+1}=0$. All higher terms would also be zero. We wouldn’t have an infinite series but a finite polynomial. Any polynomial increases more slowly than $e^{\alpha\rho}$, so the term $e^{-\alpha\rho}$ will eventually beat it down, and the function $f$ will go to zero for large $\rho$. The only bound-state solutions are those for which $\alpha=1/n$, with $n=1$, $2$, $3$, $4$, and so on.  Looking back to Eq. (19.16), we see that the bound-state solutions to the spherically symmetric wave equation can exist only when \begin{equation*} -\epsilon=1,\,\frac{1}{4},\,\frac{1}{9},\,\frac{1}{16},\dotsc, \frac{1}{n^2},\dotsc \end{equation*} The allowed energies are just these fractions times the Rydberg, $E_R=me^4/2\hbar^2$, or the energy of the $n$th energy level is \begin{equation} \label{Eq:III:19:24} E_n=-E_R\,\frac{1}{n^2}. \end{equation} There is, incidentally, nothing mysterious about negative numbers for the energy. The energies are negative because when we chose to write $V=-e^2/r$, we picked our zero point as the energy of an electron located far from the proton. When it is close to the proton, its energy is less, so somewhat below zero. The energy is lowest (most negative) for $n=1$, and increases toward zero with increasing $n$.  Before the discovery of quantum mechanics, it was known from experimental studies of the spectrum of hydrogen that the energy levels could be described by Eq. (19.24), where $E_R$ was found from the observations to be about $13.6$ electron volts. Bohr then devised a model which gave the same equation and predicted that $E_R$ should be $me^4/2\hbar^2$. But it was the first great success of the Schrödinger theory that it could reproduce this result from a basic equation of motion for the electron.  Now that we have solved our first atom, let’s look at the nature of the solution we got. Pulling all the pieces together, each solution looks like this: \begin{equation} \label{Eq:III:19:25} \psi_n=\frac{f_n(\rho)}{\rho}=\frac{e^{-\rho/n}}{\rho}\,g_n(\rho), \end{equation} where \begin{equation} \label{Eq:III:19:26} g_n(\rho)=\sum_{k=1}^na_k\rho^k \end{equation} and \begin{equation} \label{Eq:III:19:27} a_{k+1}=\frac{2(k/n-1)}{k(k+1)}\,a_k. \end{equation} So long as we are mainly interested in the relative probabilities of finding the electron at various places we can pick any number we wish for $a_1$. We may as well set $a_1=1$. (People often choose $a_1$ so that the wave function is “normalized,” that is, so that the integrated probability of finding the electron anywhere in the atom is equal to $1$. We have no need to do that just now.)  For the lowest energy state, $n=1$, and \begin{equation} \label{Eq:III:19:28} \psi_1(\rho)=e^{-\rho}. \end{equation} For a hydrogen atom in its ground (lowest-energy) state, the amplitude to find the electron at any point drops off exponentially with the distance from the proton. It is most likely to be found right at the proton, and the characteristic spreading distance is about one unit in $\rho$, or about one Bohr radius, $r_B$.  Putting $n=2$ gives the next higher level. The wave function for this state will have two terms. It is \begin{equation} \label{Eq:III:19:29} \psi_2(\rho)=\biggl(1-\frac{\rho}{2}\biggr)e^{-\rho/2}. \end{equation} The wave function for the next level is \begin{equation} \label{Eq:III:19:30} \psi_3(\rho)=\biggl(1-\frac{2\rho}{3}+\frac{2}{27}\,\rho^2\biggr) e^{-\rho/3}. \end{equation} The wave functions for these first three levels are plotted in Fig. 19–2. You can see the general trend. All of the wave functions approach zero rapidly for large $\rho$ after oscillating a few times. In fact, the number of “bumps” is just equal to $n$—or, if you prefer, the number of zero-crossings of $\psi_n$ is $n-1$. 