Solution (source code)

= Solution

The stationary <Schrödinger equation> is
$$
-\frac{\hbar^2}{2m}\Delta\psi+V\psi=E\psi.
$$
The <Laplacian> represents the kinetic-<energy> operator. For a <spherically symmetric function>, $\Delta\psi=\psi''+2\psi'/r$. The attractive <Coulomb potential> of the nucleus is $V=-q^2/(4\pi\epsilon_0r)$; $E$ is the <energy> <eigenvalue>. Multiplication by $-2m/\hbar^2$ produces the displayed radial equation with <Bohr radius>
$$
\boxed{a=\frac{4\pi\epsilon_0\hbar^2}{mq^2}.}
$$
For a normalizable bound state, write $\kappa^2=-2mE/\hbar^2$ and take $\lambda>0$. Substitution of $\psi=(A+Br)e^{-\lambda r}$ in
$$
\psi''+\frac2r\psi'+\left(\frac2{ar}-\kappa^2\right)\psi=0
$$
gives, after division by the exponential, the coefficients
$$
\begin{array}{c|c}
r&(\lambda^2-\kappa^2)B\\
1&(\lambda^2-\kappa^2)A-4\lambda B+2B/a\\
1/r&2B-2\lambda A+2A/a.
\end{array}
$$
If $B=0$ and $A\neq0$, these force $\lambda=1/a$ and $\kappa=\lambda$. If $B\neq0$, the first coefficient gives $\kappa^2=\lambda^2$, the second gives $\lambda=1/(2a)$, and the third gives $B=-A/(2a)$. Up to nonzero scalar factors, the two radial bound states and their energies are consequently
$$
\boxed{\begin{aligned}
\psi_1(r)&=e^{-r/a},& E_1&=-\frac{\hbar^2}{2ma^2},\\
\psi_2(r)&=(r-2a)e^{-r/(2a)},&E_2&=-\frac{\hbar^2}{8ma^2}.
\end{aligned}}
$$
They are linearly independent as functions and correspond to different <energy> <eigenvalue>[eigenvalues]. Strictly, they are not two solutions of one fixed-<energy> differential equation: each is a solution at its own <energy>.

For the higher-<energy> state, the spherical <Born rule> includes the radial factor $4\pi r^2$. Its normalization constant and the angular factor cancel in the ratio
$$
P(r<a)=\frac{\int_0^a r^2(r-2a)^2e^{-r/a}\,dr}
{\int_0^\infty r^2(r-2a)^2e^{-r/a}\,dr}.
$$
The given <integral> formula makes the denominator $8a^5$ and the numerator $8a^5-21a^5/e$. Thus
$$
\boxed{P(r<a)=1-\frac{21}{8e}.}
$$