Solving the Heisenberg equation of motion for the spin operators

This question is motivated by the first problem in the second chapter of Modern Quantum Mechanics (by J. J. Sakurai).

I want to find the time evolution of the spin operators ($$S_x,S_y,S_z$$) in the Heisenberg picture of an electron ($$e<0$$) subject to the following Hamiltonian: $$H = -\frac{eB}{mc}S_z=\omega S_z, \quad \omega \equiv -\frac{eB}{mc}S_z$$

I want to do it by solving the Heisenberg equation of motion.

On the notation used: $$A^H$$ is an operator in the Heisenberg formalism and $$A$$ in the Schrodinger formalism.

The Heisenberg equation of motion is:

$$\frac{dA^H}{dt} = \frac{1}{i\hbar} [A^H,U^\dagger H U]$$

For any $$S_i$$: $$\frac{dA^H}{dt} = (1/i\hbar)[U^\dagger S_i U,U^\dagger H U]$$ $$= (1/i\hbar) \omega [U^\dagger S_i U,U^\dagger S_z U]$$ $$=(\omega/i\hbar) (U^\dagger S_i U U^\dagger S_z U - U^\dagger S_i U U^\dagger S_z U )$$ usung $$UU^\dagger =1$$ $$= (\omega/i\hbar) U^\dagger (S_iS_z - S_zS_i) U = (\omega/i\hbar) U^\dagger [S_i,S_z] U$$ $$= \omega U^\dagger \epsilon_{i,z,j}S_j U$$

For the case $$S_z$$, we clearly have $$d_tS_z^H = 0$$. Therefore, $$S_z^H (t) =S_z^H (0) =S_z^S$$.

For the other cases, let's consider the time-evolution operator: $$U=exp(-iHt/\hbar) = exp(-i\omega S_zt/\hbar)$$. Then, we have:

$$\frac{dS_i^H}{dt} = \omega \epsilon_{i,z,j} e^{-i \omega S_z t/\hbar} S_j e^{i \omega S_z t/\hbar}= \omega \epsilon_{i,z,j} S_j^H$$

which gives a system of two coupled equations:

$$\frac{dS_x^H}{dt} = -\omega S_y^H$$ $$\frac{dS_y^H}{dt} = \omega S_x^H$$

which leads to: $$\frac{d^2S_x^H}{dt^2} = -\omega^2 S^H_x$$

and

$$\frac{d^2S_y^H}{dt^2}= -\omega^2 S^H_y$$

I am not sure how to proceed, since $$S^H_{x,y}$$ are operators, hence a solution of the type $$S^H_{x,y} = Ae^{i\omega t}+Be^{-i\omega t}$$ doesn't make sense.

Is this approach correct?

If not (most likely), where is the flaw?

If it makes sense, what is the best way to go forward?

Each operator is represented as a 2x2 matrix:

$$S_x^H(t)=\begin{pmatrix}a(t)&b(t)\\c(t)&f(t)\end{pmatrix}$$

Then the equation

$$\frac{d^2S_x^H}{dt^2}=-\omega^2S_x^H(t)$$

becomes

$$\begin{pmatrix}\frac{d^2a}{dt^2}&\frac{d^2b}{dt^2}\\\frac{d^2c}{dt^2}&\frac{d^2f}{dt^2}\end{pmatrix}=-\omega^2\begin{pmatrix}a(t)&b(t)\\c(t)&f(t)\end{pmatrix}$$

Solve component-by-component, and insert initial conditions.