Interesting relationship between the 2D Harmonic Oscillator and Pauli Spin matrices I have an isotropic 2D Harmonic Oscillator in cartesian coordinates
\begin{equation}
H = \frac{p_x^2}{2m} + \frac{p_y^2}{2m} + \frac{1}{2} m\omega^2 (x^2 + y^2)
\end{equation}
In terms of the usual creation and annihilation operators for the $x$ and $y$ modes, $n_x$ and $n_y$ this can be written as
\begin{equation}
H = \hbar\omega(n_x^\dagger n_x + n_y^\dagger n_y + 1)
\end{equation}
Now we can 'construct', three operators that commute with the Hamiltonian, apart from the (trivial) number operators $N_x$ and $N_y$:
\begin{align}
V_x &= a_x^\dagger a_y + a_y^\dagger a_x\\
V_y &= i(a_y^\dagger a_x - a_x^\dagger a_y)\\
V_z &= a_x^\dagger a_x - a_y^\dagger a_y
\end{align}
Now these operators can been proved to satisfy the commutation relations of angular momentum. In fact,
\begin{equation}
V_i = a^\dagger \sigma_i a
\end{equation}
Where $a = \begin{pmatrix} a_x & a_y \end{pmatrix}$ and $\sigma_i$ are the Pauli matrices. This is really surprising to me as I don't see why and how this must be true. Can someone shed some light on this? I think for $2s+1$ dimensions, we can 'construct' such operators using the matrix representations of spin $s$. However, I haven't seen a proof of this.
 A: It is well known the operators
$$
\hat C_{ij}=a_i^\dagger a_j\, ,i, j=1,\ldots, n
$$
span the Lie algebra $\mathfrak{u}(n)$.  If you choose $\hat C_{ij}, i\ne j$ and use $h_i= \hat C_{ii}-\hat C_{i+1,i+1}$, you get instead the Lie algebra $\mathfrak{su}(n)$.  Yours is just the special case with $n=2$.
The $n=3$ case is discussed in details in

Fradkin DM. Three-dimensional isotropic harmonic oscillator and
SU$_3$. American Journal of Physics. 1965 Mar;33(3):207-11,

and indeed the connection with harmonic oscillator for general $n$ goes back to

Jauch JM, Hill EL. On the problem of degeneracy in quantum mechanics. Physical Review. 1940 Apr 1;57(7):641.

Note that general representations (not of the fully symmetric type) of $\mathfrak{su}(n)$  having Young diagrams with $m<n$ rows usually start with $U(m\times n)\supset U(m)\otimes U(n)$ so the bosons $a^\dagger_{\alpha i}$, $\alpha=1,\ldots m$, $i=1,\ldots,n$,
etc and the operators
$$
\hat C_{ij}=\sum_\alpha a^\dagger_{\alpha i}a_{\alpha j}\, .
$$
A: The Hamiltonian in question has cylindrical symmetry, and can be transformed to cylindrical coordinates (i.e. $x,y\longrightarrow r, \phi$). The wave function then decomposes into radial and angular part, and the eigenstates obtained called Fock-Darwin states. Interestingly, the solution is simple also works with magnetic field, which affects only the angular momentum. You can google to find more information, e.g., here.
Perhaps another way to see it is in terms of Schwinger bosons (google again), which is a representation of angular momentum in terms of two bosons, i.e., two oscillators.
