In this post, we will derive the energy levels of the Hydrogen atom using only operator algebra, i.e., without dealing with wave functions. Our work will mainly consist in defining special operators which will assist us in determining the eigenvalues of the Hamiltonian $H$ directly. Our model of the Hydrogen atom is an electron of mass $m$ under the influence of a Coulomb potential; the associated Hamiltonian should be known by any physics student who has taken a quantum mechanics course: $$ \displaystyle H = \frac{{\bf p}^2}{2m} - \frac{e^2}{r} \label{post_d379f7e34709563115ebd4c41241ed5e_hamiltonian} $$ The first term on the right-hand side above represents the kinetic energy of the electron and the second term represents its electrostatic potential energy in the presence of a proton fixed at the origin. On the equation above, ${\bf p}$ and ${\bf x}$ are the momentum and position operators respectively, and $r = \|{\bf x}\|$.

Most of our derivation will rely on the properties of the self-adjoint operator ${\bf A}$ defined below (called the Laplace-Runge-Lenz vector operator): $$ \displaystyle{\bf A} = \frac{1}{2}({\bf L}\times{\bf p} - {\bf p}\times{\bf L}) + me^2\frac{\bf x}{r} $$ where ${\bf L} = {\bf x}\times{\bf p}$ is the angular momentum operator. The self-adjoint property of ${\bf A}$ comes from the fact that ${\bf x}$ and $r$ are self-adjoint and also because $({\bf L}\times{\bf p})^{\dagger} = - {\bf p}\times{\bf L}$. In what follows, we will extensively use the Einstein summation notation, i.e., a product $a_ib_i$ represents the sum $\sum_{i=1}^3 a_ib_i$ for two vector quantities ${\bf a} = (a_1, a_2, a_3)$ and ${\bf b} = (b_1, b_2, b_3)$ respectively.

Basic commutators

First, let us list some well-known commutators and compute some others which we will need later: $$ \begin{eqnarray} [x_i, x_j] &=& 0 \\[5pt] [p_i, p_j] &=& 0 \\[5pt] [x_i, p_j] &=& i\hbar\delta_{ij} \\[5pt] [x_i, L_j] &=& [x_i, \epsilon_{jkl}x_k p_l] = \epsilon_{jkl} x_k [x_i, p_l] = \epsilon_{jkl} x_k i\hbar\delta_{il} = i\hbar \epsilon_{jki} x_k = i\hbar \epsilon_{ijk} x_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_xi_Lj} \\[5pt] [p_i, L_j] &=& [p_i, \epsilon_{jkl}x_k p_l] = \epsilon_{jkl} [p_i, x_k] p_l = \epsilon_{jkl} (-i\hbar)\delta_{ik} p_l = -i\hbar\epsilon_{jil} p_l = i\hbar\epsilon_{ijk} p_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_pi_Lj} \\[5pt] \end{eqnarray} $$ where above we used multiple properties of the Levi-Civita symbol $\epsilon_{ijk}$. Using the following identity: $$ [A, BC] = ABC - BCA + BAC - BAC = [A,B]C + B[A,C] \label{post_d379f7e34709563115ebd4c41241ed5e_comm_A_BC} $$ we have that: $$ \begin{eqnarray} [L_i, L_j] &=& [\epsilon_{imn}x_m p_n, \epsilon_{jkl}x_kp_l] \nonumber\\[5pt] &=& \epsilon_{imn}\epsilon_{jkl} [x_m p_n, x_k p_l] \nonumber\\[5pt] &=& \epsilon_{imn}\epsilon_{jkl} \left( [x_m p_n, x_k ] p_l + x_k [x_m p_n, p_l] \right) \nonumber\\[5pt] &=& \epsilon_{imn}\epsilon_{jkl} \left(x_m [p_n, x_k] p_l + x_k [x_m, p_l] p_n\right) \nonumber\\[5pt] &=& \epsilon_{imn}\epsilon_{jkl} \left( x_m (-i\hbar)\delta_{kn} p_l + x_k i\hbar\delta_{ml} p_n\right) \nonumber\\[5pt] &=& i\hbar \left(-\epsilon_{imk}\epsilon_{jkl} x_m p_l + \epsilon_{iln}\epsilon_{jkl} x_k p_n\right) \nonumber\\[5pt] &=& i\hbar \left(\epsilon_{imk}\epsilon_{jlk} x_m p_l - \epsilon_{inl}\epsilon_{jkl} x_k p_n\right) \nonumber\\[5pt] &=& i\hbar \big((\delta_{ij}\delta_{ml} - \delta_{il}\delta_{mj}) x_m p_l - (\delta_{ij}\delta_{nk} - \delta_{ik}\delta_{nj}) x_k p_n\big) \nonumber\\[5pt] &=& i\hbar \left(\delta_{ij}x_m p_m - \delta_{il}\delta_{mj} x_m p_l - \delta_{ij}x_n p_n + \delta_{ik}\delta_{nj} x_k p_n\right) \nonumber\\[5pt] &=& i\hbar \left( - \delta_{il}\delta_{mj} x_m p_l + \delta_{ik}\delta_{nj} x_k p_n\right) \nonumber\\[5pt] &=& i\hbar \left( - \delta_{il}\delta_{jm} x_m p_l + \delta_{im}\delta_{jl} x_m p_l\right) \nonumber\\[5pt] &=& i\hbar \left( \delta_{im}\delta_{jl} - \delta_{il}\delta_{jm} \right)x_m p_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{mlk} x_m p_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{kml} x_m p_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk} L_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Lj} \end{eqnarray} $$ Additionally, since $r = \sqrt{x_ix_i}$: $$ \displaystyle [p_i, f(r)] = -i\hbar \frac{\partial f(r)}{\partial x_i} = -i\hbar \frac{d f(r)}{d r}\frac{\partial r}{\partial x_i} = -i\hbar \frac{d f(r)}{d r}\frac{x_i}{r} $$ we have that: $$ \displaystyle \left[p_i, \frac{1}{r}\right] = -i\hbar \frac{d}{d r}\left(\frac{1}{r}\right)\frac{x_i}{r} = -i\hbar \left(-\frac{1}{r^2}\right)\frac{x_i}{r} = i\hbar \frac{x_i}{r^3} \label{post_d379f7e34709563115ebd4c41241ed5e_eq_p1r} $$ Using equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_eq_p1r} and the fact that $[x_i, f(r)] = 0$, we can then show that: $$ \left[L_i, \frac{1}{r}\right] = \epsilon_{ijk}\left[x_j p_k, \frac{1}{r}\right] = \epsilon_{ijk}x_j \left[p_k, \frac{1}{r}\right] = \epsilon_{ijk} x_j i\hbar \frac{x_k}{r^3} = \frac{i\hbar}{r^3} \epsilon_{ijk} x_j x_k = 0 \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_1r} $$ since $\epsilon_{ijk}$ is an anti-symmetric tensor. Equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_1r} is directly related to the fact that $r$ is invariant under rotations.

Computation of $[L_i, A_j]$

Using the commutators we computed above, we can compute $[L_i, A_j]$ directly: $$ \begin{eqnarray} \displaystyle [L_i, A_j] &=& \left[L_i, \frac{1}{2}\epsilon_{jkl}(L_k p_l - p_k L_l) + me^2 \frac{x_j}{r}\right] \nonumber\\[5pt] &=& \frac{1}{2}\epsilon_{jkl}\left([L_i, L_k p_l] - [L_i,p_k L_l]\right) + me^2 \left[L_i,\frac{x_j}{r}\right] \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Aj_tmp} \end{eqnarray} $$ We will compute each of the terms above separately. From equations \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_pi_Lj}, \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_A_BC} and \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Lj}, we have that: $$ \begin{eqnarray} \epsilon_{jkl}[L_i, L_k p_l] &=& \epsilon_{jkl}([L_i, L_k] p_l + L_k [L_i, p_l]) \nonumber\\[5pt] &=& \epsilon_{jkl}i\hbar\epsilon_{ikm} L_m p_l + \epsilon_{jkl}L_k(-i\hbar)\epsilon_{lim} p_m \nonumber\\[5pt] &=& i\hbar\epsilon_{jlk}\epsilon_{imk} L_m p_l - i\hbar \epsilon_{jkl}\epsilon_{iml} L_k p_m \nonumber\\[5pt] &=& i\hbar\big((\delta_{ji}\delta_{lm} - \delta_{jm}\delta_{li}) L_m p_l - (\delta_{ji}\delta_{km} - \delta_{jm}\delta_{ki}) L_k p_m \big)\nonumber\\[5pt] &=& i\hbar\left(\delta_{ji} L_m p_m - \delta_{jm}\delta_{li} L_m p_l - \delta_{ji}L_m p_m + \delta_{jm}\delta_{ki} L_k p_m \right)\nonumber\\[5pt] &=& i\hbar\left(- \delta_{jm}\delta_{li} L_m p_l + \delta_{jm}\delta_{ki} L_k p_m \right)\nonumber\\[5pt] &=& i\hbar\left(- \delta_{jm}\delta_{li} L_m p_l + \delta_{jl}\delta_{mi} L_m p_l \right)\nonumber\\[5pt] &=& i\hbar\left(\delta_{im}\delta_{jl} - \delta_{il}\delta_{jm} \right) L_m p_l\nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{mlk} L_m p_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{kml} L_m p_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk} ({\bf L}\times{\bf p})_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Lkpl} \end{eqnarray} $$ $$ \begin{eqnarray} \epsilon_{jkl}[L_i, p_k L_l] &=& \epsilon_{jkl}([L_i, p_k] L_l + p_k [L_i, L_l]) \nonumber\\[5pt] &=& \epsilon_{jkl}(-i\hbar)\epsilon_{kim} p_m L_l + \epsilon_{jkl} p_k i\hbar\epsilon_{ilm} L_m \nonumber\\[5pt] &=& i\hbar\epsilon_{jlk}\epsilon_{imk} p_m L_l - i\hbar \epsilon_{jkl}\epsilon_{iml} p_k L_m \nonumber\\[5pt] &=& i\hbar\big((\delta_{ji}\delta_{lm} - \delta_{jm}\delta_{li}) p_m L_l - (\delta_{ji}\delta_{km} - \delta_{jm}\delta_{ki}) p_k L_m \big)\nonumber\\[5pt] &=& i\hbar\left(\delta_{ji} p_m L_m - \delta_{jm}\delta_{li} p_m L_l - \delta_{ji} p_m L_m + \delta_{jm}\delta_{ki} p_k L_m \right)\nonumber\\[5pt] &=& i\hbar\left(- \delta_{jm}\delta_{li} p_m L_l + \delta_{jm}\delta_{ki} p_k L_m \right)\nonumber\\[5pt] &=& i\hbar\left(- \delta_{jm}\delta_{li} p_m L_l + \delta_{jl}\delta_{mi} p_m L_l \right)\nonumber\\[5pt] &=& i\hbar\left(\delta_{im}\delta_{jl} - \delta_{il}\delta_{jm} \right) p_m L_l\nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{mlk} p_m L_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk}\epsilon_{kml} p_m L_l \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk} ({\bf p}\times{\bf L})_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_pkLl} \end{eqnarray} $$ Finally, using equations \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_xi_Lj} and \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_1r}, we can show that: $$ \displaystyle \left[L_i,\frac{x_j}{r}\right] = \frac{1}{r} [L_i, x_j] = \frac{1}{r} (-i\hbar) \epsilon_{jik} x_k = i\hbar \epsilon_{ijk} \frac{x_k}{r} = i\hbar \epsilon_{ijk}\left(\frac{{\bf x}}{r}\right)_k \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_xjr} $$ Using equations \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Lkpl}, \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_pkLl} and \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_xjr} on equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Aj_tmp} then yields: $$ \begin{eqnarray} \displaystyle [L_i, A_j] &=& \frac{1}{2}\epsilon_{jkl}\left([L_i, L_k p_l] - [L_i,p_k L_l]\right) + me^2 \left[L_i,\frac{x_j}{r}\right] \nonumber\\[5pt] &=& \frac{1}{2}\big(i\hbar \epsilon_{ijk} ({\bf L}\times{\bf p})_k - i\hbar \epsilon_{ijk} ({\bf p}\times{\bf L})_k\big) + me^2 i\hbar \epsilon_{ijk}\left(\frac{{\bf x}}{r}\right)_k \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk} \left( \frac{1}{2}\big( ({\bf L}\times{\bf p})_k - ({\bf p}\times{\bf L})_k\big) + me^2 \left(\frac{{\bf x}}{r}\right)_k \right) \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Aj_tmp2} \end{eqnarray} $$ Therefore: $$ \boxed{ [L_i, A_j] = i\hbar \epsilon_{ijk} A_k } \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Li_Aj} $$

Computation of $[A_i, A_j]$

The computation of $[A_i, A_j]$ is very long and tedious. Since it does not introduce any new techniques beyond what has been presented above, we will merely state the result of this commutator, but the reader is strongly encouraged to try to prove the equation below as its derivation is an excellent way to deeply understand everything we have done so far (a lot of helpful material for that can be found here): $$ \boxed{ [A_i, A_j] = -i\hbar 2m \epsilon_{ijk} L_k H } $$ where $H$ is the Hamiltonian given in equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_hamiltonian}.

The ${\bf K}^+$ and ${\bf K}^-$ operators

Suppose we have an eigenstate of $H$ representing a bound state with energy eigenvalue $E = -\kappa^2 / 2m$ for some $\kappa \gt 0$ (the assumption $E \lt 0$ is correct for a bound state since $E \geq 0$ means the electron has enough kinetic energy to escape the proton's electric potential). Let then ${\bf K}^{\pm}$ be defined as below: $$ \displaystyle {\bf K}^{\pm} = \frac{1}{2}{\bf L} \pm \frac{1}{2\kappa}{\bf A} $$ We are interested in computing $[K^+_i, K^-_j]$ and $[K^\pm_i, K^\pm_j]$: $$ \begin{eqnarray} \displaystyle [K^+_i, K^-_j] &=& \frac{1}{4}[L_i, L_j] - \frac{1}{4\kappa}[L_i,A_j] + \frac{1}{4\kappa}[A_i,L_j] - \frac{1}{4\kappa^2}[A_i,A_j] \nonumber\\[5pt] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k - \frac{1}{\kappa}\epsilon_{ijk}A_k - \frac{1}{\kappa} \epsilon_{jik}A_k - \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k H\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k - \frac{1}{\kappa}\epsilon_{ijk}A_k + \frac{1}{\kappa} \epsilon_{ijk}A_k - \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k H\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4} \left(\epsilon_{ijk}L_k - \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k H\right) \end{eqnarray} $$ Applying this result to the energy eigenstate mentioned above, we have $H \rightarrow E = -\kappa^2/2m$, so we get: $$ \begin{eqnarray} \displaystyle [K^+_i, K^-_j] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k - \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k \frac{(-\kappa^2)}{2m}\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k - \epsilon_{ijk}L_k\right) \end{eqnarray} $$ Therefore: $$ \boxed{ \displaystyle [K^{+}_i, K^{-}_j] = 0 } $$ Also: $$ \begin{eqnarray} \displaystyle [K^\pm_i, K^\pm_j] &=& \frac{1}{4} [L_i, L_j] \pm \frac{1}{4\kappa}[L_i, A_j] \pm \frac{1}{4\kappa}[A_i, L_j] + \frac{1}{4\kappa^2}[A_i,A_j] \nonumber\\[5pt] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k \pm \frac{1}{\kappa}\epsilon_{ijk}A_k \mp \frac{1}{\kappa} \epsilon_{jik}A_k + \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k H\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4}\left(\epsilon_{ijk}L_k \pm \frac{1}{\kappa}\epsilon_{ijk}A_k \pm \frac{1}{\kappa} \epsilon_{ijk}A_k + \frac{1}{\kappa^2} (-2m) \epsilon_{ijk}L_k H\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4} \epsilon_{ijk} \left( L_k \pm \frac{2}{\kappa} A_k + \frac{1}{\kappa^2}(-2m)L_k \frac{(-\kappa^2)}{2m}\right) \nonumber\\[5pt] &=& \frac{i\hbar}{4} \epsilon_{ijk} \left( 2L_k \pm \frac{2}{\kappa}A_k\right) \nonumber\\[5pt] &=& i\hbar \epsilon_{ijk} \left( \frac{1}{2} L_k \pm \frac{1}{2\kappa}A_k\right) \end{eqnarray} $$ Therefore: $$ \boxed{ [K^\pm_i, K^\pm_j] = i\hbar \epsilon_{ijk} K^{\pm}_k } \label{post_d379f7e34709563115ebd4c41241ed5e_comm_Ki_Kj} $$

Computation of ${\bf A}^2$

$$ \begin{eqnarray} {\bf A}^2 &=& {\bf A}\cdot{\bf A} \nonumber\\[5pt] &=& \frac{1}{4}\big(({\bf L}\times{\bf p} )\cdot( {\bf L}\times{\bf p}) - ({\bf L}\times{\bf p})\cdot( {\bf p}\times{\bf L})\big) \nonumber\\[5pt] &+& \frac{1}{4}\big(-({\bf p}\times{\bf L})\cdot({\bf L}\times{\bf p}) + ({\bf p}\times{\bf L})\cdot( {\bf p}\times{\bf L})\big) \nonumber\\[5pt] &+& \frac{me^2}{2}({\bf L}\times{\bf p} - {\bf p}\times{\bf L})\cdot \frac{\bf x}{r} + \frac{me^2}{2}\frac{\bf x}{r} \cdot ({\bf L}\times{\bf p} - {\bf p}\times{\bf L}) \nonumber\\[5pt] &+& m^2e^4\frac{\bf x}{r}\cdot \frac{\bf x}{r} \label{post_d379f7e34709563115ebd4c41241ed5e_AA_tmp} \end{eqnarray} $$ Each of the terms on the equation above will be computed separately: $$ \begin{eqnarray} ({\bf L}\times{\bf p} )\cdot( {\bf L}\times{\bf p}) &=& ({\bf L}\times{\bf p} )_i( {\bf L}\times{\bf p})_i \nonumber\\[5pt] &=&(\epsilon_{ijk}L_j p_k)(\epsilon_{ilm}L_l p_m) \nonumber\\[5pt] &=&\epsilon_{jki}\epsilon_{lmi}L_j p_k L_l p_m \nonumber\\[5pt] &=&(\delta_{jl}\delta_{km} - \delta_{jm}\delta_{kl})L_j p_k L_l p_m \nonumber\\[5pt] &=&L_j p_k L_j p_k - L_j p_k L_k p_j \end{eqnarray} $$ Since $\epsilon_{ijk}$ is anti-symmetric, we have that: $$ p_k L_k = {\bf p}\cdot{\bf L} = \epsilon_{ijk} p_i x_j p_k = 0 $$ The equation above then becomes: $$ \begin{eqnarray} ({\bf L}\times{\bf p} )\cdot( {\bf L}\times{\bf p}) &=& L_j p_k L_j p_k \nonumber\\[5pt] &=&L_j (L_j p_k + [p_k, L_j]) p_k \nonumber\\[5pt] &=&L_j L_j p_k p_k + i\hbar\epsilon_{kjm}p_m p_k \nonumber\\[5pt] &=&{\bf L}^2{\bf p}^2 \end{eqnarray} $$ Using a nearly identical sequence of derivation steps, we can also show that: $$ \begin{eqnarray} ({\bf L}\times{\bf p})\cdot( {\bf p}\times{\bf L}) &=& -{\bf L}^2{\bf p}^2 \\[5pt] ({\bf p}\times{\bf L})\cdot( {\bf L}\times{\bf p}) &=& -({\bf L}^2 + 4\hbar^2){\bf p}^2 \\[5pt] ({\bf p}\times{\bf L})\cdot( {\bf p}\times{\bf L}) &=& {\bf L}^2{\bf p}^2 \end{eqnarray} $$ Let us now proceed to the other terms: $$ \displaystyle({\bf L}\times{\bf p}) \cdot \frac{\bf x}{r} = \epsilon_{ijk} L_j p_k x_i\frac{1}{r} = - L_j \epsilon_{jik} p_k x_i\frac{1}{r} = - L_j L_j^\dagger\frac{1}{r} = -\frac{\bf L^2}{r} $$ $$ \displaystyle \frac{\bf x}{r} \cdot ({\bf p}\times{\bf L}) = \left( ({\bf p}\times{\bf L})^\dagger \cdot \left(\frac{\bf x}{r}\right)^\dagger \right)^\dagger = \left( -({\bf L}\times{\bf p}) \cdot \frac{\bf x}{r} \right)^\dagger = \left( \frac{\bf L^2}{r} \right)^\dagger = \frac{\bf L^2}{r} $$ $$ \begin{eqnarray} \displaystyle \frac{\bf x}{r} \cdot ({\bf L}\times{\bf p}) &=& \frac{1}{r} \epsilon_{ijk} x_i L_j p_k \nonumber\\[5pt] &=& \frac{1}{r} \epsilon_{ijk} (L_j x_i + [x_i, L_j]) p_k \nonumber\\[5pt] &=& \frac{1}{r} (\epsilon_{ijk} L_j x_i p_k + \epsilon_{ijk} i\hbar \epsilon_{ijm} x_m p_k) \nonumber\\[5pt] &=& \frac{1}{r} (-\epsilon_{jik} L_j x_i p_k + i\hbar \epsilon_{ikj} \epsilon_{imj} x_m p_k) \nonumber\\[5pt] &=& \frac{1}{r} (- L_j L_j + i\hbar (\delta_{ii}\delta_{km} - \delta_{im}\delta_{ki}) x_m p_k) \nonumber\\[5pt] &=& \frac{1}{r} (- {\bf L}^2 + i\hbar (3x_k p_k - x_i p_i)) \nonumber\\[5pt] &=& \frac{1}{r} (- {\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p}) \end{eqnarray} $$ $$ \begin{eqnarray} \displaystyle ({\bf p}\times{\bf L}) \cdot \frac{\bf x}{r} &=& \left( \left( \frac{\bf x}{r} \right)^\dagger \cdot ({\bf p}\times{\bf L})^\dagger \right)^\dagger \nonumber\\[5pt] &=& \left( \frac{\bf x}{r} \cdot (-{\bf L}\times{\bf p}) \right)^\dagger \nonumber\\[5pt] &=& -\left( \frac{1}{r} (- {\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p}) \right)^\dagger \nonumber\\[5pt] &=& -(-{\bf L}^2 - 2 i\hbar {\bf p}\cdot{\bf x}) \frac{1}{r} \nonumber\\[5pt] &=& \left({\bf L}^2 + 2 i\hbar p_i x_i\right) \frac{1}{r} \nonumber\\[5pt] &=& \left({\bf L}^2 + 2 i\hbar (x_i p_i + [p_i, x_i]) \right) \frac{1}{r} \nonumber\\[5pt] &=& \left({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 2 i\hbar (-i\hbar)\delta_{ii}\right) \frac{1}{r} \nonumber\\[5pt] &=& ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 6\hbar^2) \frac{1}{r} \nonumber\\[5pt] &=& \frac{1}{r} ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 6\hbar^2) + \left[ 2 i\hbar {\bf x}\cdot{\bf p}, \frac{1}{r}\right] \nonumber\\[5pt] &=& \frac{1}{r} ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 6\hbar^2) + 2i\hbar x_i \left[ p_i, \frac{1}{r}\right] \nonumber\\[5pt] &=& \frac{1}{r} ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 6\hbar^2) + 2i\hbar x_i i\hbar \frac{x_i}{r^3} \nonumber\\[5pt] &=& \frac{1}{r} ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 6\hbar^2) - 2\hbar^2 \frac{1}{r} \nonumber\\[5pt] &=& \frac{1}{r} ({\bf L}^2 + 2 i\hbar {\bf x}\cdot{\bf p} + 4\hbar^2) \end{eqnarray} $$ Putting all the results above on equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_AA_tmp} gives us: $$ \begin{eqnarray} {\bf A}^2 &=& \frac{1}{4}({\bf L}^2{\bf p}^2 + {\bf L}^2{\bf p}^2 + {\bf L}^2{\bf p}^2 + 4\hbar^2 {\bf p}^2 + {\bf L}^2{\bf p}^2) \nonumber\\[5pt] &+& \frac{me^2}{2r}\left(-{\bf L}^2 -{\bf L}^2 - 2i\hbar{\bf x}\cdot{\bf p} - 4\hbar^2 -{\bf L}^2 + 2i\hbar{\bf x}\cdot{\bf p} - {\bf L}^2\right) \nonumber\\[5pt] &+& m^2e^4 \nonumber\\[5pt] &=& {\bf L}^2{\bf p}^2 + \hbar^2{\bf p}^2 - \frac{2me^2}{r}{\bf L}^2 - \frac{2me^2\hbar^2}{r} + m^2e^4 \nonumber\\[5pt] &=& ({\bf L}^2 + \hbar^2)\left({\bf p}^2 - \frac{2me^2}{r}\right) + m^2e^4 \nonumber\\[5pt] &=& ({\bf L}^2 + \hbar^2)2mH + m^2e^4 \end{eqnarray} $$ Applying this result to our energy eigenstate, we have $H \rightarrow E = -\kappa^2/2m$, so we get: $$ \boxed{ {\bf A}^2 = -\kappa^2({\bf L}^2 + \hbar^2) + m^2e^4 } \label{post_d379f7e34709563115ebd4c41241ed5e_AA} $$

Determination of $E$ in terms of ${\bf K}^\pm$

First, since ${\bf L}$ and ${\bf A}$ are self-adjoint operators: $$ {\bf L}\cdot{\bf A} + {\bf A}\cdot{\bf L} = {\bf L}\cdot{\bf A} + ({\bf L}^\dagger\cdot{\bf A}^\dagger)^\dagger = {\bf L}\cdot{\bf A} + ({\bf L}\cdot{\bf A})^\dagger $$ Given that: $$ \displaystyle {\bf L}\cdot{\bf A} = {\bf L}\cdot\frac{1}{2}({\bf L}\times{\bf x} - {\bf x}\times{\bf L}) + me^2{\bf L}\cdot\frac{\bf x}{r} = me^2({\bf x}\times{\bf p})\cdot\frac{\bf x}{r} = 0 $$ then: $$ {\bf L}\cdot{\bf A} + {\bf A}\cdot{\bf L} = 0 \label{post_d379f7e34709563115ebd4c41241ed5e_LA_AL} $$ From equations \eqref{post_d379f7e34709563115ebd4c41241ed5e_AA} and \eqref{post_d379f7e34709563115ebd4c41241ed5e_LA_AL}, we have that: $$ \begin{eqnarray} \displaystyle ({\bf K}^\pm)^2 &=& {\bf K}^\pm\cdot{\bf K}^\pm \nonumber\\[5pt] &=& \frac{1}{4}{\bf L}\cdot{\bf L} \pm \frac{1}{4\kappa}({\bf L}\cdot{\bf A} + {\bf A}\cdot{\bf L}) + \frac{1}{4\kappa^2}{\bf A}\cdot{\bf A} \nonumber\\[5pt] &=& \frac{1}{4}{\bf L}^2 + \frac{1}{4\kappa^2}( -\kappa^2({\bf L}^2 + \hbar^2) + m^2e^4 ) \nonumber\\[5pt] &=& \frac{1}{4}{\bf L}^2 - \frac{1}{4}{\bf L}^2 - \frac{\hbar^2}{4} + \frac{m^2e^4}{4\kappa^2} \end{eqnarray} $$ Therefore: $$ ({\bf K}^+)^2 = ({\bf K}^-)^2 = \frac{1}{4}\left(-\hbar^2 + \frac{m^2e^4}{\kappa^2}\right) \label{post_d379f7e34709563115ebd4c41241ed5e_Kp2} $$ With this result, we then have that: $$ \boxed{ \displaystyle E = -\frac{\kappa^2}{2m} = -\frac{m^2e^4}{2m}\left(\hbar^2 + 4({\bf K}^\pm)^2\right)^{-1} } \label{%INDEX_E_tmp} $$

A final expression for $E$

As shown on equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_comm_Ki_Kj}, the operators $K^+_i$ and $K^-_i$ satisfy the commutation relations for angular momentum operators, i.e.: $$ \begin{eqnarray} [K^+_i, K^+_j] &=& i\hbar \epsilon_{ijk} K^+_k \\[5pt] [K^-_i, K^-_j] &=& i\hbar \epsilon_{ijk} K^-_k \end{eqnarray} $$ It then follows directly that: $$ \begin{eqnarray} [({\bf K}^+)^2, K^+_i] &=& 0 \\[5pt] [({\bf K}^-)^2, K^-_i] &=& 0 \end{eqnarray} $$ Using ladder operators, it is straightforward to prove that the eigenvalues of $({\bf K}^+)^2$ and $({\bf K}^-)^2$ have the form $\hbar^2 j(j+1)$ for values of $j$ which are either integers of half-integers, i.e., $j = (n-1)/2$ for integer values $n \geq 1$. Since we are assuming we have an energy eigenstate with energy $E = -\kappa^2/2m$, and since equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_Kp2} implies that such a state is also an eigenstate of both $({\bf K}^+)^2$ and $({\bf K}^-)^2$, then from equation \eqref{%INDEX_E_tmp} we have that: $$ \begin{eqnarray} \displaystyle E &=& -\frac{m^2e^4}{2m}\left(\hbar^2 + 4({\bf K}^\pm)^2\right)^{-1} \nonumber\\[5pt] &=& -\frac{m^2e^4}{2m}\left(\hbar^2 + 4\hbar^2j(j+1)\right)^{-1} \nonumber\\[5pt] &=& -\frac{m^2e^4}{2m}\left(\hbar^2 + 4\hbar^2\frac{n-1}{2}\frac{n+1}{2}\right)^{-1} \nonumber\\[5pt] &=& -\frac{m^2e^4}{2m}\left(\hbar^2 + \hbar^2(n^2 - 1)\right)^{-1} \end{eqnarray} $$ Finally, we have our desired result, which we express as $E_n$ instead of $E$ to make the dependence on $n$ explicit: $$ \displaystyle E_n = -\frac{me^4}{2\hbar^2n^2} $$ Since we used $4\pi\epsilon_0 = 1$ on equation \eqref{post_d379f7e34709563115ebd4c41241ed5e_hamiltonian}, we can replace $e^2$ with $e^2 / 4\pi\epsilon_0$ to get the energy levels in SI units: $$ \boxed{ \displaystyle E_n = -\frac{me^4}{2\hbar^2(4\pi\epsilon_0)^2n^2} \quad \textrm{for} \quad n = 1, 2, 3, \ldots } $$