This in an answer to the a p.SE question.
Supersymmetry
Well, for fixed \(l\), the degeneracy of \(m\) is because of SO(3) symmetry, we are just seeing a full representation of this group.
The big question is why all the radial hamiltonians \(H_l\) for different angular momenta have the same spectrum except a discrete number of eigenvalues.
Note that particularly the tower-spectrum for \(l\) and the tower for \(l+1\) only differ in one eigenvalue, the lowest energy one. This is the typical setup one can read in Witten’s Supersymmetric Quantum Mechanics: a pair of hamiltonians differing only in the vacuum eigenstate. So you should be able to build a supersymmetry generator Q such that \( H_1=QQ^+\) is the radial hamiltonian for angular momentum \(l\) and \(H_2=Q^+Q\) is the radial hamiltonian for angular momentum \(l+1\).
SUSY QM is simpler than QFT QM; it does not contemplate Spin; the state and the superpartner are just two levels in QM hamiltonians. It is just a little bit more advanced, mathematically, that the factorisation method; still it allows for some topological arguments on susy breaking that generalize to the QFT version, this was the idea of Witten when defining it.
Just now I am not sure if this the connection for every potential having $l,m$ degeneracy does exist, or only for Coulomb-Hydrogen case; to start with, it implies that the potential $V(r)$ must come from a superpotential, so surely it is not so trivial to do, not to classify all the families of radial potentials that allow to do this trick. But is is a twenty years old idea by now, so surely it is already done.
Ok, even there is an entry in the wikipedia. According it, the superpotential is
$$W = \frac{\sqrt{2m}}{h} \frac{\lambda}{2(l+1)} – \frac{h(l+1)}{r\sqrt{2m}}$$
So that the potentials
$$V_-=W^2-W’= -\lambda \frac{1}{r} + \frac{h^2 l (l+1)} {2m} \frac{1}{r^2}- \frac{\lambda^2 m}{2 h^2 (l+1)^2} $$
and
$$V_+=W^2+W’=-\lambda \frac{1}{r} + \frac{h^2 (l+1) (l+2)} {2m} \frac{1}{r^2} + \frac{\lambda^2 m}{2 h^2 (l+1)^2}$$
have the same spectrum except for the lowest energy eigenvalue of the first one, which is zero and can not appear in the second one (nice topological result).
The advantage of this explanation is that it can be extended to potentials without the full SO(4) symmetry and to more exotic cases where the pairing fails for other eigenvalues.
PS: It can be noticed that the superpotential for Coulomb problem is just a constant away from\(W(r)=1/r\). An interesting point is that this superpotential can be calibrated to pair with the free particle:\(V_+(r)=W^2+W’=0\); the superpotentials having this property generate the so-called “transparent potentials”, with special properties in the phase-shift. They can be thought as generalizing the radial equation to symmetric spaces, with\(1/r\) being the euclidean case.
connection with so(4) group representations (and Runge-Lenz vector?)
According last page of this lecture, the role of Runge-Lenz vector as a supercharge is analytically tricky. But at least we get some help from group theory, if we recall that so(4) ~ su(2) x su(2) and that so(3) ~ su(2). So for our purposes we could really write
$$so(4) \approx su(2) \oplus so(3)$$
The rotational part, so(3), gives us the $m$ degeneracy inside a representation of the group of rotations; this should exists for every central potential. The $su(2)$ part is the one whose ladder operator allows an interpretation as supersymmetry charge, where the degeneracy is not complete because of the difference in the lowest energy eigenstate; usually it disappears because $Q |\Omega>$ is not normalizable. (I am a bit fascinated that the susy ladder operator is related to an SU(2), because in susy QM this is not needed, or at least not explicit)
The supersymmetry generator Q (alt. $Q^+$) when applied to the eigenfunction of a hamiltonian produces the corresponding eigenfunction in the partner. This is the same role that the ladder generator used to produce states inside a representation of the symmetry group, but in this view it comes from the susy pairing: if
$$H_2 \psi = Q^+Q \psi = E\psi$$
then
$$H_1 (Q\psi) = (Q Q^+) (Q \psi) = Q H_2 \psi = E (Q\psi)$$
Note that the pairing fails if $Q\psi$ doesn’t exist; this is the case for the vaccum potential, but I remember that J Casahorran did some study for other eigenstates beyond the vacuum (it is tricky because of Witten’s results).