Hydrogen Atom Wave Function...
Standard approaches to the Hydrogen atom often leap directly to Laguerre polynomials. Here, we maintain mathematical rigor by utilizing the Whittaker differential equation to solve the radial component.
We begin with the time-independent Schrödinger equation:
\begin{equation} \left( -\frac{\hbar^2}{2m_e}\Delta - \frac{e^2}{4\pi\varepsilon r} \right)\psi = E\psi \label{eq:a1_1} \end{equation}Where $m_e$ is the mass of the electron, which can be found mostly around the nucleous in a cloud of prabibility, since the hydrgen atom problem has a central potential $1/r$, it means it has a rotationally invariance symmetry $SO(3)$, which means we're going to adapt a Laplacien $\Delta$ adequately:
\begin{equation} \Delta f = \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 \frac{\partial f}{\partial r} \right) + \frac{1}{r^2 \sin \theta} \frac{\partial}{\partial \theta} \left( \sin \theta \frac{\partial f}{\partial \theta} \right) + \frac{1}{r^2 \sin^2 \theta} \frac{\partial^2 f}{\partial \phi^2} \label{eq:a1_2} \end{equation}We define the orbital angular momentum operator $L^2$ "Casimir invariant in $\mathfrak{so}(3)$" as follows:
\begin{equation} \hat{L}^2 f = \frac{1}{\sin \theta} \frac{\partial}{\partial \theta} \left( \sin \theta \frac{\partial f}{\partial \theta} \right) + \frac{1}{\sin^2 \theta} \frac{\partial^2 f}{\partial \phi^2} \label{eq:a1_3} \end{equation}Which makes our Laplacien takes the following form:
\begin{equation} \Delta f = \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 \frac{\partial f}{\partial r} \right) + \frac{\hat{L}^2 f}{r^2} \label{eq:a1_4} \end{equation}Upon applying the orbital angular momentum operator on a random state we get:
\begin{equation} \hat{L}^2 f(\theta,\phi) = \ell(\ell+1) f(\theta,\phi) \label{eq:a1_5} \end{equation}1. Radial Separation
Isolating the radial part $\psi_r$ and including the centrifugal term:
\begin{equation} \frac{\partial^2\psi_r}{\partial r^2} + \frac{2}{r}\frac{\partial\psi_r}{\partial r} + \left[ \frac{2m_eE}{\hbar^2} - \frac{\ell(\ell+1)}{r^2} + \frac{2me^2}{4\pi\varepsilon\hbar^2 r} \right]\psi_r = 0 \label{eq:a1_6} \end{equation}Next we try to compactify the form of the differential equation by defining some parameters, bohr radius $a$ and the wavevctor $k$ as follows:
\begin{equation} a=\frac{m_ee^2}{4\pi\varepsilon\hbar^2} \ \ \ \ \ \ k^2=\frac{2m_eE}{\hbar^2} \label{eq:a1_7} \end{equation}Which makes our differential equation look much cleaner
\begin{equation} \frac{\partial^2\psi_r}{\partial r^2} + \frac{2}{r}\frac{\partial\psi_r}{\partial r} + \left[ k^2 - \frac{\ell(\ell+1)}{r^2} + \frac{2}{a r} \right]\psi_r = 0 \label{eq:a1_8} \end{equation}2. Transformation to Whittaker Form
In order to find the Whittaker differential equation form, we need to perform a smal change of variable:
\begin{equation} \phi(r)= \frac{\psi_r(r)}{r} \label{eq:a1_9} \end{equation}This will end up removing this term $\frac{2}{r}\frac{\partial\psi_r}{\partial r}$, then we multiply our resultant differential equation with $1/\beta^2$
\begin{equation} \frac{\partial^2\phi}{\partial \beta^2 r^2} + \left[ \frac{k^2}{\beta^2} - \frac{\ell(\ell+1)}{\beta^2 r^2} + \frac{2}{a \beta^2 r} \right]\phi = 0 \label{eq:a1_10} \end{equation}We define a new variable $\rho$
\begin{equation} \rho = \beta r \label{eq:a1_11} \end{equation}Which allows us to write the differential equation as follows:
\begin{equation} \frac{\partial^2\phi}{\partial \rho^2} + \left[ \frac{k^2}{\beta^2} + \frac{2/a\beta}{\rho} - \frac{\ell(\ell+1)}{\rho^2} \right]\phi = 0 \label{eq:a1_12} \end{equation}By setting the constants to be equal to the followig:
\begin{equation} \frac{k^2}{\beta^2}=-\frac{1}{4} \ \ \ \ \ \ \frac{2}{a\beta}=n \label{eq:a1_13} \end{equation}The equation transforms into the Whittaker form[1]:
\begin{equation} \frac{\partial^2\phi}{\partial \rho^2} + \left[ -\frac{1}{4} + \frac{n}{\rho} - \frac{\ell(\ell+1)}{\rho^2} \right]\phi = 0 \label{eq:a1_14} \end{equation}This means that we uncovered what $\beta$ equals, meaning that the energy is pretty much found
\begin{equation} \beta=\frac{2}{na} \implies k^2=-\frac{1}{4}\beta^2 \implies \frac{2m_eE}{\hbar^2}=-\frac{1}{n^2a^2} \implies E_n=-\frac{\hbar^2}{2m_ea^2n^2} \label{eq:a1_15} \end{equation}Which is the correct answer for the energy eigenvalues of the hydrgen atom
3. The Total Normalized Wavefunction
Using the Whittaker function $M_{n, \ell+1/2}$, we define the radial normalization constant $A$:
\begin{equation} A = \sqrt{ \frac{2}{na \int_0^{+\infty} M_{n, \ell+1/2}^2(\rho) d\rho} } \label{eq:a1_16} \end{equation}The total wavefunction $\psi_{n,l,m,\alpha}$ including the spherical harmonics and the pauli spinor $\alpha$ is then:
\begin{equation} \psi_{n,\ell,m,\alpha}(\mathbf{r}) = \frac{1}{r} A \cdot M_{n, l+1/2}\left( \frac{2r}{na} \right) Y_{lm}(\theta, \phi) \otimes \alpha \label{eq:a1_17} \end{equation}Where the angular momentum coupling is defined as:
\begin{equation} |\ell m_\ell \otimes s m_s\rangle = \sum_{\tau=|\ell-s|}^{\ell+s} C_{\ell m_\ell s m_s}^{\tau m_\tau} |\tau m_\ell+m_s\rangle \label{eq:a1_18} \end{equation}