\section{Introduction}

Since the pioneering work of Ashkin \citep{ashkin2006optical}, trapping
dielectric particles by optical forces has found a wide range of applications
\citep{grier2003revolution}, ranging from sensing ultra weak forces
\citep{hebestreit2018sensing,ranjit2015attonewton,ranjit2016zeptonewton,geraci2010short,geraci2015sensing,hempston2017force,blakemore2019three,monteiro2020force,monteiro2017optical}
and temperature gradients \citep{bykov2015flying,zeltner2018flyingparticle,zhang2014novel}
to manipulating biological systems \citep{ashkin1987cells,ashkin1987bacteria,lang2003resource,neuman2004optical,marago2013optical,grier2003revolution}.
Furthermore, dielectric particles levitated in vacuum are highly isolated
from environmental noise and decoherence, offering an ideal playground
to control the dynamics of micro- and nanoresonators \citep{gieseler2013thermal,gieseler2014dynamic,millen2014nanoscale,ricci2017optically,arita2013laserinduced,hebestreit2018measuring,hoang2016torsional,rashid2018precessionmotion}.
This has brought forth unprecedented performance in cavity and feedback
cooling of  levitated particle motion \citep{delic2019cavity,windey2019cavitybased,tabbenjohanns2019colddamping,conangla2019optimalfeedback,gieseler2012subkelvin,li2011millikelvin}
and the realization of fast rotation in high vacuum \citep{reinmann2018ghzrotation,ahn2018opticallylevitated}.
The high degree of isolation also makes levitated particles a key
component in many proposals for exploring quantum mechanics with mesoscopic
objects \citep{chang2010cavity,romero2011optically,romero2010toward,bateman2014near,scala2013matter,magrini2021real,tebbenjohanns2021quantum}.

Magnetically ordered dielectrics introduce magnetization as an additional
degree of freedom to levitated particles. For example, due to conservation of angular momentum, a change in the magnetization
of a freely moving body can
lead to its mechanical rotation and vice versa -- the so called Einstein-de
Haas effect \citep{einstein_de_haas} and its reciprocal, the Barnett
effect \citep{barnett1915magnetization}. Furthermore, the magnetization
couples to the particle orientation, defined by the specific magneto-crystalline
anisotropy of the material \citep{chikazumi2009physics}. These effects
are prominent in magnetically levitated systems and can be harnessed
for probing mesoscopic quantum mechanics including rotational degrees of freedom, as well as for sensing \citep{romero2012quantum,wang2019dynamics,timberlake2019acceleration,lewandowski2021high,slezak2018cooling,hsu2016cooling,prat2017ultrasensitive,gieseler2020single,o2019magneto,houlton2018axisymmetric,johnsson2016macroscopic,cirio2012quantum,vinante2020ultralow,perdriat2021spin}. 

The magnetization of dielectrics also can couple to light via magneto-optical
effects, the best known being the Faraday effect: the magnetization-induced
rotation of the polarization of light. Conversely, light can also
induce a small precession of the magnetization (the inverse Faraday
effect) \citep{landau1984electrodynamics, stancil2009spin}. This generally weak coupling can be enhanced by confinement
of light inside of the dielectric, which has given rise to the field of  (cavity optomagnonics) \citep{violakusminskiyCoupledSpinlightDynamics2016a,liu2016optomagnonics,sharma2017light,almpanis2018dielectric,almpanis2020spherical,osada2018orbital,kusminskiy2019cavity,PhysRevLett.117.133602,haigh2018selection,PhysRevLett.116.223601,osada2018brillouin,zhang2016optomagnonic, rameshti2021cavity}.
The combination of optomagnonic effects with levitation provides a
novel system within which light, magnetism and angular motion are intertwined.

In this work, we propose a new scheme for optomagnonically controlling and probing the coupled dynamics of levitated magnetic dielectric microparticles.
We consider the system depicted in Fig.$~$\ref{fig:schematic_model}:
a spherical dielectric magnet is levitated in the presence of an external static magnetic field. A pump laser drives two internal electromagnetic field modes
of the particle, which couple to the magnetization dynamics via the
cavity-enhanced optomagnonic coupling. In turn, the magnetization
of the particle couples to its angular motion via both magneto-crystalline
anisotropy and dissipation of magnons into the lattice. We show that
the optically driven magnetization dynamics induce angular oscillations,
which exhibit low damping due to the lack of direct dissipation of
the angular motion. The driven angular motion is limited in amplitude
by the optical quality factor of the particle and the driving
power that the system can support. The coupling between magnetization and angular motion can be
probed via the power spectrum of the cavity modes, which exhibits two
main resonance peaks that can be linked to the angular oscillations and to the
spin dynamics. Finally, we study the case of a driven ferromagnetic
resonance (FMR) by adding an oscillatory magnetic field perpendicular
to the static magnetic field. We show that this configuration can enhance the visibility of the characteristic peaks in the spectrum and can  induce, via the magnetization, fast rotations
of the particle around its anisotropy axis. 

The paper is structured as follows. In Sec. \ref{sec:Model} we present
the model and derive the corresponding Hamiltonian. In Sec. \ref{sec:Equations-of-motion}
we derive the coupled Heisenberg equations of motion with their conserved
quantities and discuss the torque generated by the applied laser via
the optomagnonic coupling. Sec. \ref{sec:SteadyState-1} presents
the steady state configurations of the system. In Sec. \ref{sec:Force-and-dynamical}
we investigate the nonlinear dynamics of the system numerically. Sec.
\ref{sec:Power-Spectrum} is devoted to the computation of the power
spectrum of the cavity field fluctuations and the resonance peaks.
In Sec.$~$\ref{sec:Driven-System} we show the effects of an additional
driving magnetic field and we discuss some generalizations to the
model in Sec.$~$\ref{sec:Generalisations}. We draw conclusions and
discuss further directions in Sec. \ref{sec:Conclusion}. 


\section{Model\label{sec:Model}}

We consider a Faraday-active dielectric microparticle that is levitated
in the presence of an external static magnetic field $\mathbf{B}_0$
as depicted in Fig.$~$\ref{fig:schematic_model}. The magnetization
of the particle is saturated along $\mathbf{B}_0$. The particle supports
internal electromagnetic modes. If the radius of the particle is much
larger than the driving laser wavelength, light is trapped by total
internal reflection, forming optical whispering gallery modes \citep{heebner2008optical}.
Otherwise, for radii comparable to the driving light wavelength, such
modes correspond to Mie resonances \citep{bohren2008absorption,kuznetsov2016optically}.
In case of an unmagnetized particle the underlying Mie mode frequencies are solely
determined by the polarization and the angular mode number, and the modes
are degenerate with respect to the azimuthal number. The magnetization
of the particle breaks this degeneracy and causes a splitting of the
Mie modes into modes with different azimuthal numbers \citep{almpanis2020spherical,ford1978scattering}.

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{schematic}
\caption{Schematic of the  model, showing the couplings between the relevant
degrees of freedom. (a) A levitated ferrimagnetic microsphere serves
as an optomagnonic cavity. The microsphere exhibits uniaxial anisotropy
and we assume that its magnetization is homogeneous, pointing initially along
the direction of an external magnetic field. (b) The intracavity field
couples via optomagnonic effects to the magnetization which in turn
couples to the rotation through the magnetic anisotropy. (c) Schematic
illustration of the relevant frequencies of the system and the pumping
scheme. Two optical modes have frequencies $\omega_{1}$ and $\omega_{2}$
separated by the frequency of the magnon mode, i.e. $\omega_{1}-\omega_{2}=\omega_{B}$.
One of the two modes is pumped in resonance with the laser frequency
$\omega_{\mathrm{L}}$.}
\label{fig:schematic_model}
\end{figure}

Driving the particle with light can induce a precession of the magnetization
around its equilibrium orientation. The associated quantized magnetic excitations
are denoted magnons. The mean number of magnons in a coherent state is directly related to the amplitude of the magnetization precession. In the following we focus on the homogeneous
magnon mode, known as the Kittel mode, consisting of a uniformly precessing
magnetization. Its frequency
can be tuned by the external magnetic field and usually lies in the GHz regime. For a driving electromagnetic field at optical frequencies, the resonant dipolar coupling between the magnetization and the magnetic field is suppressed due to the mismatch of frequencies.
The coupling to the electromagnetic field takes place instead through a  Raman scattering
type of process involving two photons, and it is dominated by the
dipolar coupling to the \emph{electric component} of the field \citep{fleury1968scattering}.

In this work, we consider that a pump laser of suitable wavelength can drive simultaneously two internal
optical modes of the levitated particle which satisfy the so-called \emph{triple-resonance condition} where the frequency
difference between the cavity modes equals the magnetization's precession
frequency \citep{PhysRevLett.117.133602,zhang2016optomagnonic}. In
other words, transitions between the modes can occur through emission
or absorption of a magnon. The triple-resonance condition can be fulfilled
in the visible/near-infrared part of the electromagnetic spectrum,
where the frequency splitting between neighboring Mie modes for micrometer-sized particles is in
the range of a few GHz, the typical frequency of magnons \citep{almpanis2020spherical}. We consider the only source of birefringence (and thus
of torque) to be magnetic. This is valid for example for an Yttrium Iron Garnet (YIG) sphere. For magnetic dielectric
particles which also exhibit optical or geometric anisotropy \citep{nieminen2001optical,friese1998optical,simpson2007polarization,donato2016light}, the model can be generalized accordingly. Our model takes optical losses and magnetic damping into consideration. The latter corresponds to damping of the magnetization precession,
i.e., the relaxation to its equilibrium configuration \citep{gilbert2004phenomenological}.
For the coupling of the confined light to the particle magnetization we
focus on the Faraday effect (circular birefringence), relegating consideration
of the Cotton-Mouton effect (linear birefringence) to Sec.$~$\ref{sec:Cotton-Mouton-Effect}.
The magnetization and the orientation of the particle are coupled
by a magnetocrystalline anisotropy, such that a mechanical torque
can be induced. We assume that the material exhibits uniaxial crystalline
anisotropy, such that the magnetization, in the absence of a magnetic
field, is aligned either parallel to the anisotropy axis, in the case
of an easy-axis magnet, or perpendicular to it, in the case of a hard-axis
magnet. We assume that the translational degrees of freedom are completely
decoupled from the rest, namely, we neglect effects related to the
center of mass motion of the particle. We therefore do not explicitly
include the trapping potential in our model. We discuss briefly deviations
from these assumptions in Sec.$~$\ref{sec:Generalisations}. 

A suitable material for an experimental realization of this system
is YIG, since it exhibits ferrimagnetic properties at room temperature
and shows high transparency in the near-infrared. Further, it has
a large Faraday rotation per unit length $\theta_{\mathrm{F}}=240^{\circ}/\mathrm{cm}$
with a dielectric index of $\varepsilon\approx5$ \citep{stancil2009spin,antonov1969optical}
at these frequencies, making it the material of choice in the state-of-the-art
optomagnonic experiments \citep{PhysRevLett.117.133602,haigh2018selection,PhysRevLett.116.223601,osada2018brillouin,zhang2016optomagnonic}.
In the following we investigate the proposed model for parameter values
corresponding to YIG, although the framework applies in general for
other materials too.


\subsection{Hamiltonian}

The Hamiltonian describing the system reads
\begin{align}
\hat{H} & =\hat{H}_{\mathrm{Sp}}+\hat{H}_{\mathrm{Cav}}+\hat{H}_{\mathrm{OM}}\label{eq:initialHamiltonian}
\end{align}
where $\hat{H}_{\mathrm{Sp}}$ is the Hamiltonian of the freely rotating
magnetic sphere in the presence of the magnetic field and $\hat{H}_{\mathrm{Cav}}$
that of the optical modes confined in the dielectric sphere. $\hat{H}_{\mathrm{OM}}$
describes the coupling between the optical modes and the magnetization
of the sphere. The operator character of the quantities is represented
by the ``hat'' notation. In the following we present each term in detail.

The first term in Eq.$~$\eqref{eq:initialHamiltonian} is given by
\begin{align}
\hat{H}_{\mathrm{Sp}} & =\dfrac{\hbar^{2}\hat{\mathbf{L}}^{2}}{2I}-\hbar\gamma\hat{\mathbf{S}}\cdot\mathbf{B}_0-\hbar^{2}D(\hat{\mathbf{S}}\cdot\mathbf{e}_{{\rm 3}}(\hat{\Omega}))^{2}\,,\label{eq:HSp}
\end{align}
where the first term is the kinetic energy of the rotational motion
with dimensionless angular momentum $\hat{\mathbf{L}}=(\hat{L}_{x},\hat{L}_{y},\hat{L}_{z})$
and $I=2mR^{2}/5$ the moment of inertia of a sphere of radius $R$
and mass $m$. In the following we will express the inertia in terms
of a rotation frequency of the rigid solid $\omega_{I}=\hbar/I$.

The second term in Eq.$~$\eqref{eq:HSp} corresponds to the Zeeman energy describing the interaction
between the particle's magnetization and the homogeneous magnetic field
$\mathbf{B}_0$. For the Kittel mode, all spins precess in phase and
can be described by the dimensionless macrospin $\hat{\mathbf{S}}=(\hat{S_{x}},\hat{S}_{y},\hat{S}_{z})$,
which relates to the magnetization as $\hat{\mathbf{M}}=\hat{\mathbf{S}}M_{\mathrm{S}}/S$
with $S=|\mathbf{S}|$ the total spin of the sample and $M_{\mathrm{S}}$ the saturation magnetization. $\gamma=1.76\cdot10^{11}\,\mathrm{s}^{-1}\mathrm{T}^{-1}$
denotes the absolute value of the gyromagnetic ratio. We consider
a magnetic field parallel to the $z$ axis, namely $\mathbf{B}_0=B_0\mathbf{e}_{z}$,
such that the Zeeman energy of the Hamiltonian reduces to $-\hbar\gamma\hat{S}_{z}B_0$
and we denote $\gamma B_0=\omega_{B}$ as the Larmor precession frequency.
For a sphere, this corresponds to the frequency of the Kittel mode. A magnon corresponding to this mode has energy $\hbar \omega_{B}$.

The third term in Eq.$~$\eqref{eq:HSp} accounts for uniaxial magnetocrystalline anisotropy (cubic anisotropy will be discussed in Sec.$~$\ref{subsec:Cubic-anisotropy}).
The anisotropy axis is fixed to the particle and labeled as $\mathbf{e}_{3}(\hat{\Omega})$.
As $\mathbf{e}_{3}(\hat{\Omega})$ depends on the Euler angles $\hat{\Omega}$,
this term couples the spin with the orientation of the sphere with
coupling constant $\omega_{D}=\hbar D$. Here $|D|=K_{u}V/(\hbar S)^{2}$
\citep{gatteschi2006molecular, rusconi2017quantum} denotes
the anisotropy strength with $K_{u}$ the anisotropy energy density
and $V$ the volume of the particle. For $D>0$, $\mathbf{e}_{3}(\hat{\Omega})$
is called an ``easy axis'', i.e. the interaction energy is minimized
if the magnetization and $\mathbf{e}_{3}$ are parallel. For $D<0$,
$\mathbf{e}_{3}(\hat{\Omega})$ is a ``hard axis'', and in this case
it is energetically favorable for the magnetization to lie in the
plane perpendicular to $\mathbf{e}_{3}$. 

The orientation of $\mathbf{e}_{3}(\hat{\Omega})$ is defined by an
operator specified by the Euler angles $\hat{\Omega}=\{\hat{\alpha},\hat{\beta},\hat{\gamma}\}$,
such that $\mathbf{e}_{3}(\hat{\Omega})$ is given in terms of elements
of the rotation matrix $R(\hat{\Omega})$ (Eq.$~$\eqref{eq:rotationMatrix-1}
and \eqref{eq:e3-1-1}) connecting the fixed laboratory frame $O\mathbf{e}_{x}\mathbf{e}_{y}\mathbf{e}_{z}$
to the rotating frame $O\mathbf{e}_{1}\mathbf{e}_{2}\mathbf{e}_{3}$
(which we call the body frame): $(\mathbf{e}_{1},\mathbf{e}_{2},\mathbf{e}_{3})^{T}=R(\hat{\Omega})(\mathbf{e}_{x},\mathbf{e}_{y},\mathbf{e}_{z})^{T}$
(see Fig.$~$\ref{Fig:euler_angles}). Note that from now on for simplicity
we write $\mathbf{e}_{3}$ without the explicit dependence on $\hat{\Omega}$.

The cavity Hamiltonian $\hat{H}_{\mathrm{Cav}}$ describes the electromagnetic
field inside the particle and is given by 
\begin{equation}
\hat{H}_{{\rm Cav}}=\dfrac{1}{2}\int\mathrm{d}^{3}\mathbf{r}\left(\varepsilon_{0}\mathbf{\hat{E}}^{\ast}(\mathbf{r},t)\bar{\mathbf{\boldsymbol{\varepsilon}}}\mathbf{\hat{E}}(\mathbf{r},t)+\dfrac{1}{\mu_{0}}\hat{\mathbf{B}}^{2}(\mathbf{r},t)\right)\label{eq:H_cav}
\end{equation}
with $\varepsilon_{0}$ and $\mu_{0}$ denoting, respectively, the
vacuum permittivity and permeability, and $\bar{\boldsymbol{\varepsilon}}$
the static part of the effective relative permittivity tensor that depends
on the magnetization $\mathbf{M}$. Considering only the Faraday effect and a saturation magnetization
along the $z$ axis, the permittivity tensor is given by 
\begin{equation}
\boldsymbol{\varepsilon}(\mathbf{M})=\bar{\boldsymbol{\varepsilon}}+\delta\boldsymbol{\varepsilon},\label{eq:permittivity_tensor_Faraday}
\end{equation}
with the static
\begin{equation}
\bar{\boldsymbol{\varepsilon}}=\begin{pmatrix}\varepsilon & -\mathrm{i}fM_{z} & 0\\
\mathrm{i}fM_{z} & \varepsilon & 0\\
0 & 0 & \varepsilon
\end{pmatrix},\label{eq:permittivity_static}
\end{equation}
and dynamical part
\begin{equation}
\delta\boldsymbol{\varepsilon}=\begin{pmatrix}0 & 0 & \mathrm{i}fM_{y}\\
0 & 0 & -\mathrm{i}fM_{x}\\
-\mathrm{i}fM_{y} & \mathrm{i}fM_{x} & 0
\end{pmatrix}.\label{eq:permittivity_dynamic}
\end{equation}
Here $\varepsilon$ denotes the relative isotropic permittivity of
the unmagnetized material, and $f$ a material-dependent constant.
It is related to the Faraday rotation per unit length $\theta_{\mathrm{F}}$, 
\begin{equation}
\theta_{\mathrm{F}}=\dfrac{\omega fM_{\mathrm{S}}}{2c\sqrt{\varepsilon}},\label{eq:faraday_rot}
\end{equation}
for optical frequency $\omega$ and vacuum speed of light $c$. The decomposition into a static and dynamical part as given in Eqs. \eqref{eq:permittivity_tensor_Faraday}-\eqref{eq:permittivity_dynamic} is valid for small deviations of the magnetization with respect to the $z$-axis, such that $M_z$ can be considered approximately constant $M_z\approx M_{\mathrm{S}}$ ($M_{x,y}/M_{\mathrm{S}}\ll 1$). The static magnetization causes a Zeeman-like splitting of the optical Mie modes \citep{almpanis2020spherical},
i.e. it lifts the $(2l+1)$- degeneracy with respect
to the azimuthal number $-l\leq m\leq l$, where $l$ is the angular mode number. This spliting is in the GHz range \citep{almpanis2020spherical}.

The electric and magnetic intracavity field, $\hat{\mathbf{E}}$ and
$\hat{\mathbf{B}}$, are quantized in the usual way such that $\hat{H}_{\mathrm{Cav}}$
is written in terms of the cavity operators $\hat{a}_{k}$ and the
corresponding frequencies $\omega_{k}$ (with $k$ labelling the mode
indices) as
\begin{equation}
\hat{H}_{\mathrm{Cav}}=\sum_{k}\hbar\omega_{k}\bigg(\hat{a}_{k}^{\dagger}\hat{a}_{k}+\dfrac{1}{2}\bigg).\label{eq:H_cav_quant}
\end{equation}

The optomagnonic Hamiltonian $\hat{H}_{\mathrm{OM}}$ is derived from
the time-averaged electromagnetic energy \citep{landau1984electrodynamics,violakusminskiyCoupledSpinlightDynamics2016a} 
\begin{equation}
\bar{U}=\dfrac{1}{4}\varepsilon_{0}\int\mathrm{d^{3}\mathbf{r}}\sum_{ij}E_{i}^{\ast}(\mathbf{r},t)\delta\varepsilon_{ij}(\mathbf{M})E_{j}(\mathbf{r},t).\label{eq:em_energy}
\end{equation}
We are interested in the interaction between the Kittel mode and a single pair of WGMs for
which the optomagnonic interaction Hamiltonian, after quantizing the
electric field $\mathbf{E}(\mathbf{r},t)\rightarrow\hat{\mathbf{E}}^{+}(\mathbf{r},t)=\sum_{k}\mathbf{E}_{k}(\mathbf{r})\hat{a}_{k}(t)$
(with $\mathbf{E}_{k}(\mathbf{r})$ indicating the electric field of the $k$-th eigenmode), can be expressed as \citep{PhysRevLett.117.133602}
\begin{equation}
\hat{H}_{\mathrm{OM}}=\hbar g(\hat{S}_{+}\hat{a}_{1}^{\dagger}\hat{a}_{2}+\hat{S}_{-}\hat{a}_{2}^{\dagger}\hat{a}_{1}),\label{eq:H_om_quant}
\end{equation}
where $\hat{S}_{\pm}=\hat{S}_{x}\pm\mathrm{i}\hat{S}_{y}$. The optomagnonic
coupling constant $g$ is defined as 
\begin{equation}
g=-\dfrac{\varepsilon_{0}fM_{\mathrm{S}}}{4\hbar S}\tilde{g}\label{eq:g}
\end{equation}
with $\tilde{g}$ (given in Appendix$~$\ref{sec:Transition-Amplitude})
denoting the transition amplitude between two Zeeman-split Mie modes. Here the selection rule for the initial and final azimuthal mode
numbers $m_{\mathrm{f}}-m_{\mathrm{i}}=1$ is fulfilled, and we assume
that the frequency splitting between two neighboring modes satisfies
the triple-resonance condition. The triple resonance can be achieved by tuning the frequency of the Kittel mode by the external magnetic field  $\mathbf{B}_0$. For a micron sized YIG sphere $g\sim0.03\,\mathrm{Hz}$.
This coupling is enhanced by driving one of the two optical modes, by a factor given by the square root of the number of circulating photons in the cavity \citep{kusminskiy2019cavity,kusminskiy2019quantum}. 

\begin{table}[t]
\caption{Definition and values of the relevant frequencies and couplings in
the system scaled by the total spin $S$. The values correspond to
the physical parameters given in Table$~$\ref{tab:parameters} and
are given in terms of the magnon frequency $\omega_{B}\equiv\gamma B_0=10^{10}\,\mathrm{s^{-1}}$,
corresponding to an applied magnetic field of $\sim60~\mathrm{mT}$.
Note that the full dependence of the optomagnonic coupling $g$ on $R$ cannot be obtained analytically, due to the non-trivial dependence of the transition amplitude $\tilde{g}$ on $R$ (see Fig.$~$\ref{fig:coupling_radius}). }
\begin{tabular}{>{\raggedright}p{0.5\columnwidth}>{\raggedright}p{0.5\columnwidth}}
\toprule 
Definition & Value (for $R$ in $\mu$m)\tabularnewline
\midrule
\multicolumn{1}{l}{$\omega_{I}S\equiv\hbar S/I=5\hbar S/(2mR^{2})$} & $2.15\cdot10^{-8}\omega_{B}/R^{2}$\tabularnewline
$\omega_{D}S\equiv\hbar DS=K_{u}V/(\hbar S)$ & $2.32\cdot10^{-1}\omega_{B}$\tabularnewline
$gS=\dfrac{\varepsilon_{0}fM_{\mathrm{S}}}{4\hbar}\tilde{g}$ & $5.596\cdot10^{-2}\omega_{B}/R$\tabularnewline
\bottomrule
\end{tabular}\label{tab:couplings}
\end{table}

\begin{table}[t]
\caption{Values of the physical parameters used in the calculations. The values
correspond to YIG \citep{stancil2009spin,wu2013recent}, however,
the framework holds generally. Note that we choose $m_{\mathrm{i}}=0$
such that the coupling strength is maximized (see Appendix \ref{sec:Transition-Amplitude}).}
\begin{tabular}{ll}
\toprule 
Parameter and values & Description \tabularnewline
\midrule
$\rho_{m}=5\cdot10^{3}\,\mathrm{kg}\,\mathrm{m^{-3}}$ & Mass density\tabularnewline
$S=\rho_{m}/738\,\mathrm{amu}$ & Total spin\tabularnewline
$K_{u}=10^{3}\:\mathrm{J/\mathrm{m^{3}}}$\footnote{For the cubic symmetry of YIG one finds two constants for the anisotropy
density in the literature. Here we approximate the uniaxial constant
with a value in the order of the first-order cubic anisotropy constant.
See Sec. \ref{subsec:Cubic-anisotropy} for details.} & Magnetic anisotropy constant \tabularnewline
$\varepsilon=5$ & Relative permittivity\tabularnewline
$n=2.2$ & Refractive index\tabularnewline
$\theta_{\mathrm{F}}=240^{\circ}/\mathrm{cm}$ & Faraday rotation\tabularnewline
$M_{\mathrm{S}}=140\,\mathrm{kA}/\mathrm{{m}}$ & Saturation magnetization\tabularnewline
$\lambda_{0}=1500\,\mu\mathrm{m}$ & Laser wavelength\tabularnewline
$m_{\mathrm{i}}=0$ & Initial azimuthal mode number\tabularnewline
$\eta_{\mathrm{G}}=10^{-4}$ & Gilbert damping parameter\tabularnewline
\bottomrule
\end{tabular}\label{tab:parameters}
\end{table}

The optical modes are driven externally by a laser, described by adding
a driving term $i\hbar\epsilon_{\mathrm{L}1,2}\left(\hat{a}_{1,2}^{\dagger}\mathrm{e}^{-\mathrm{i\omega_{Las}t}}-\hat{a}_{1,2}\mathrm{e}^{\mathrm{i}\omega_{\mathrm{Las}}t}\right)$
to the Hamiltonian. The parameter $\epsilon_{Li}=\sqrt{\dfrac{\kappa_{\mathrm{rad}}P_{i}}{\hbar\omega_{\mathrm{Las}}}}$
depends on the laser power $P_{i}$ and the radiative decay rate $\kappa_{\mathrm{rad}}$
(see Eq.$~$\eqref{eq:kappa}). Collecting all terms and working in
a frame rotating with the laser frequency $\omega_{\mathrm{Las}}$,
the total system Hamiltonian reads 
\begin{align}
\nonumber
\hat{H}= & \dfrac{1}{2}\hbar\omega_{I}\hat{\mathbf{L}}^{2}-\hbar\omega_{B}\hat{S}_{z}-\hbar\Delta_{1}\hat{a}_{1}^{\dagger}\hat{a}_{1}-\hbar\Delta_{2}\hat{a}_{2}^{\dagger}\hat{a}_{2}\\\nonumber
 & -\hbar\omega_{D}\left(\hat{\mathbf{S}}\cdot\hat{\mathbf{e}}_{3}\right)^{2}+\hbar g(\hat{S}_{+}\hat{a}_{1}^{\dagger}\hat{a}_{2}+\hat{S}_{-}\hat{a}_{2}^{\dagger}\hat{a}_{1})\\
 & +\mathrm{i}\hbar\epsilon_{\mathrm{L1}}\big(\hat{a}_{1}^{\dagger}-\hat{a}_{1}\big)+\mathrm{i}\hbar\epsilon_{\mathrm{L2}}\big(\hat{a}_{2}^{\dagger}-\hat{a}_{2}\big)
\label{eq:hamiltonian_reduced}
\end{align}
with detuning $\Delta_{i}=\omega_{\mathrm{Las}}-\omega_{i}$. The
assumption that the drive laser couples to both modes is a simplification
that can be refined \footnote{For instance, one could consider individual driving terms at different
frequencies for each mode, which would require additional unitary
transformations to a rotating frame in which $\hat{H}$ is time-independent.}.

In Table$~$\ref{tab:couplings} we summarize the relevant frequencies
and couplings and show their dependence on the size of the particle
in Fig.$~$\ref{fig:coupling_radius}(a). The optomagnonic coupling
$g$ and the anisotropy strength $\omega_{D}$ have a similar magnitude
while the inertia frequency $\omega_{I}$ is orders of magnitude smaller
and decreases faster with particle radius. All the relevant values
for the physical parameters appearing in the system are listed in
Table$~$\ref{tab:parameters}. 

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{eulerAngles+SteadyState}\caption{Definition of the Euler angles. (a) The rotation matrix $R(\hat{\Omega})$
aligns the fixed laboratory frame $O\mathbf{e}_{x}\mathbf{e}_{y}\mathbf{e}_{z}$
with the body frame $O\mathbf{e}_{1}\mathbf{e}_{2}\mathbf{e}_{3}$
where $\mathbf{e}_{3}(\hat{\Omega})$ denotes the anisotropy direction.
(b) The Euler angles correspond to three successive rotations that
are necessary for aligning the lab frame with the frame in which $\mathbf{e}_{3}$
is fixed. (c) Steady state configurations of the system. }
\label{Fig:euler_angles}
\end{figure}


\begin{figure}[t]
\includegraphics[width=0.95\columnwidth]{coupling2}\caption{Frequencies and couplings defined in Table \ref{tab:couplings} as
a function of the radius $R$. The values are given in units of the
Larmor frequency $\omega_{B}=\gamma B_0$ which is controlled by the
applied magnetic field and is typically in the GHz range. Note that
the vertical axis is in logarithmic scale.}
\label{fig:coupling_radius}
\end{figure}


\section{Equations of Motion\label{sec:Equations-of-motion}}

From the Hamiltonian in Eq.$~$\eqref{eq:hamiltonian_reduced} we
obtain the equations of motion for the set of operators $\hat{\boldsymbol{\xi}}=(\hat{\mathbf{e}}_{3},\hat{\mathbf{L}},\hat{\mathbf{S}},\hat{a}_{1},\hat{a}_{1}^{\dagger},\hat{a}_{2},\hat{a}_{2}^{\dagger})$
via the Heisenberg equation 
\begin{equation}
\dfrac{\mathrm{d}}{\mathrm{d}t}\hat{\boldsymbol{\xi}}=\dfrac{\mathrm{i}}{\hbar}[\hat{H},\hat{\boldsymbol{\xi}}],\label{eq:Heisenberg_EoM}
\end{equation}
and focus on the classical limit by taking the average value of the
operators over coherent states and disregarding any correlations and noise \citep{gardiner1985input} (see
Appendix$~$\ref{sec:Commutation-Relations}). In the following we
denote $\langle\boldsymbol{\hat{\xi}}\rangle=\boldsymbol{\xi}$. We
also add dissipative terms, which we discuss in detail below. The
coupled equations of motion read 
\begin{align}
\dot{\mathbf{e}}_{3} & =-\omega_{I}\mathbf{e}_{3}\times\mathbf{L},\nonumber \\
\dot{\mathbf{L}} & =2\omega_{D}\left(\mathbf{e}_{3}\cdot\mathbf{S}\right)\left(\mathbf{e}_{3}\times\mathbf{S}\right)-\boldsymbol{\tau}_{\mathrm{G}},\nonumber \\
\dot{\mathbf{S}} & =\mathbf{B}_{{\rm eff}}\times\mathbf{S}-2\omega_{D}\left(\mathbf{e}_{3}\cdot\mathbf{S}\right)\left(\mathbf{e}_{3}\times\mathbf{S}\right)+\boldsymbol{\tau}_{\mathrm{G}},\nonumber \\
\dot{a}_{1} & =\mathrm{i}\Delta_{1}a_{1}-\mathrm{i}gS_{+}a_{2}+\epsilon_{L1}-\dfrac{1}{2}\kappa a_{1},\nonumber \\
\dot{a}_{2} & =\mathrm{i}\Delta_{2}a_{2}-\mathrm{i}gS_{-}a_{1}+\epsilon_{L2}-\dfrac{1}{2}\kappa a_{2},\label{eq:eom}
\end{align}
where $\mathbf{B}_{\mathrm{eff}}=-\omega_{B}\mathbf{e}_{z}+g(a_{1}^{\dagger}a_{2}+a_{2}^{\dagger}a_{1})\mathbf{e}_{x}+\text{i}g(a_{1}^{\dagger}a_{2}-a_{2}^{\dagger}a_{1})\mathbf{e}_{y}$.
The evolution of the total angular momentum $\mathbf{J}=\mathbf{L}+\mathbf{S}$
is given by
\begin{equation}
\dot{\mathbf{J}}=\mathbf{B}_{{\rm eff}}\times\mathbf{S},\label{eq:eomJ}
\end{equation}
and the dynamics is given either in terms of $(\mathbf{S},\mathbf{L})$
or $(\mathbf{S},\mathbf{J})$.

The equations of motion include cavity damping terms $\kappa$, which
we consider to include only radiative losses \citep{buck2003optimal}, i.e. $\kappa=\kappa_\mathrm{rad}$,
although other sources, such as surface roughness, can also be relevant
depending on the size of the particle. Furthermore, we introduced
Gilbert damping to account for intrinsic damping of the spin by adding
the torque $\boldsymbol{\tau}_{\mathrm{G}}$ \citep{gilbert2004phenomenological,keshtgar2017magnetomechanical},
\begin{equation}
\boldsymbol{\tau}_{\mathrm{G}}=\dfrac{\eta_{\mathrm{G}}}{S}\dot{\mathbf{S}}\times\mathbf{S}-\frac{\eta_{\mathrm{G}}}{S}\omega_{I}\mathbf{S}\times(\mathbf{S}\times\mathbf{L}),\label{eq:gilbert_torque}
\end{equation}
where $\eta_{\mathrm{G}}$ is a dimensionless damping parameter (in
the order of $10^{-4}$ for YIG). The second term stems from the fact
that the Gilbert damping is defined in the rotating (body) frame,
describing the relative motion of the magnetization with respect to
the lattice. Thus a coordinate transformation is required to include
damping in the laboratory frame. This additional term gives rise to
the Barnett effect, i.e. magnetization dynamics induced by rotation
\citep{barnett1915magnetization,keshtgar2017magnetomechanical}. Since
Gilbert damping is solely due to internal forces, we subtracted $\boldsymbol{\tau}_{\mathrm{G}}$
in the equation for $\dot{\mathbf{L}}$, such that the total angular
momentum is only affected by the external torque (see Eq.$~$\eqref{eq:eomJ})
\citep{band2018dynamics}. This approximation neglects,
for example, heating effects or damping into other magnon modes.

The dynamics described by Eq.$~$\eqref{eq:eom} conserves two quantities:
\begin{equation}
\begin{aligned}|\mathbf{e}_{3}|= & \sqrt{R_{31}^{2}+R_{32}^{2}+R_{33}^{2}}=1,\\
|\mathbf{S}|= & \sqrt{S_{x}{}^{2}+S_{y}^{2}+S_{z}^{2}}=S.
\end{aligned}
\label{eq:conservation_rules}
\end{equation}
Furthermore, for negligible Gilbert damping (i.e. $\eta_{\mathrm{G}}/S\rightarrow0$)
 \footnote{For YIG particles of $\sim(1\,\mu{\rm m})^{3}$, this approximation
is justified as $\eta_{\mathrm{G}}/S<10^{-14}$.}, the equations of motion yield a third conserved quantity
\begin{equation}
\mathbf{e}_{3}\cdot\mathbf{L}=\text{const}.\label{eq:conservtion_Le3}
\end{equation}
This is obtained from the first two equations of motion for $\boldsymbol{\tau}_{G}=\mathbf{0}$
as it imposes that $\dot{\mathbf{e}}_{3}\cdot\mathbf{L}+\dot{\mathbf{L}}\cdot\mathbf{e}_{3}=\dfrac{\mathrm{d}}{\mathrm{d}t}(\mathbf{e}_{3}\cdot\mathbf{L})=0$.
In this case, if $\mathbf{L}(t=0)=\bm{0}$, at any instant of time
$\mathbf{e}_{3}(t)\cdot\mathbf{L}(t)=0$, and thus, the sphere does
not rotate around its anisotropy axis. The only possible angular motions
are precession and libration, with angular speeds given by the components
of $\mathbf{L}$ perpendicular to $\mathbf{e}_{3}$. 

The equation of motion for $\mathbf{L}$ can be recast in a more familiar
form
\begin{equation}
\hbar\dot{\mathbf{L}}=\bm{\tau}=\hbar\mathbf{e}_{3}\times2\omega_{D}\left(\mathbf{e}_{3}\cdot\mathbf{S}\right)\mathbf{S}=R\mathbf{e}_{3}\times\mathbf{F}.\label{eq:torque}
\end{equation}
The angular torque $\bm{\tau}$ is therefore generated by the force
$\mathbf{F}=\frac{2\hbar\omega_{D}}{R}\left(\mathbf{e}_{3}\cdot\mathbf{S}\right)\mathbf{S}$.
The torque $\bm{\tau}$ can be driven by a laser, which induces, via
optomagnonic coupling, the spin dynamics. In particular, if the system's
initial configuration is the undriven equilibrium point $\bm{S}\parallel\bm{e}_{3}=\bm{e}_{z}$,
the only source of initial motion is the laser drive. Any subsequent
time evolution of the spin, and consequently any angular motion of
the sphere, is generated by a combination of the optically induced
torque and the torque by the external magnetic field. The latter is
responsible for a Larmor precession of the spin around the direction
of the magnetic field. 

From the equations of motion Eq.$~$\eqref{eq:eom} it can be further seen
that in the absence of crystalline anisotropy, namely $\omega_{D}=0$,
the angular motion is completely decoupled from the optomagnonic part
for $\eta_{G}=0$, i.e. in this case one recovers the dynamics for
a pure optomagnonic system as treated for one optical mode in Ref.$~$\citep{violakusminskiyCoupledSpinlightDynamics2016a}. If Gilbert
damping is considered $\eta_{G}\neq0$, then there is a coupling of
the angular motion to the spin due to conservation of angular momentum.

We will now proceed to characterize the dynamics of this driven dissipative
system. We will first obtain the possible steady state configurations
given by Eq.$~$\eqref{eq:eom}, and then study the light-spin-rotation
dynamics for small perturbations around such steady states. 


\section{Steady State \label{sec:SteadyState-1}}

We find the
possible steady states of the system by setting $\dfrac{\mathrm{d}\xi}{\mathrm{d}t}=0$
$\forall\,\xi$ in Eq.$~$\eqref{eq:eom}. Denoting by $\bar{\xi}$
the steady state value of the quantity $\xi(t)$, we find
the following possible steady states, as depicted in Fig. \ref{Fig:euler_angles}(c): 
\begin{align}
(\mathrm{I)\;} & \bar{\mathbf{L}}=\mathbf{0},\bar{\mathbf{e}}_{3}\parallel\pm\mathbf{\bar{S}}\nonumber \\
\mathrm{(II)\;} & \bar{\mathbf{L}}=\mathbf{0},\bar{\mathbf{e}}_{3}\perp\bar{\mathbf{S}}\nonumber \\
\mathrm{(III)\;} & \mathbf{\bar{L}\neq0},\bar{\mathbf{e}}_{3}\parallel\pm\bar{\mathbf{L}}\parallel\pm\bar{\mathbf{S}}\label{eq:possible_steadyState}
\end{align}
with $\mathbf{\bar{S}\parallel\pm\bar{\mathbf{B}}_{\mathrm{eff}}}$.
In the case of a non-rotating sphere ($\bar{\mathbf{L}}=\mathbf{0}$), (I)
minimizes the angular-spin interaction energy for an easy-axis anisotropy,
whereas (II) minimizes the interaction energy for a hard-axis anisotropy.
If $\mathbf{\bar{L}\neq0}$, then the sphere rotates around $\mathbf{e}_{3}$
and $\bar{\mathbf{L}}$ must be parallel to $\bar{\mathbf{S}}$, which
implies that $\bar{\mathbf{e}}_{3}\parallel\mathbf{\pm\bar{S}}$.

For the steady states of the spin and the optical fields, we find a solution with a finite average number of circulating photons in both modes ($\vert\bar{a}_{i}\vert^2\neq0$) when both modes are simultaneously driven, $\epsilon_{Li}\neq0$ (Supplement 1, Sec.$~$4).
The light creates a deflection of the magnetization, such that in
the steady-state
\begin{align}
\bar{S}_{x}= & S\cos\phi\sin\theta,\nonumber\\
\bar{S}_{y}= & S\sin\phi\sin\theta,\nonumber\\
\bar{S}_{z}= & S\cos\theta.\label{eq:SpinSteadyState-1}
\end{align}
with the azimuthal angle $\phi$
\begin{equation}
\phi=\arctan\left(\mathrm{i}\dfrac{\bar{a}_{2}^{\ast}\bar{a}_{1}-\bar{a}_{1}^{\ast}\bar{a}_{2}}{\bar{a}_{2}^{\ast}\bar{a}_{1}+\bar{a}_{1}^{\ast}\bar{a}_{2}}\right),\label{eq:azimuthalDeflection}
\end{equation}
and the deflection $\theta$ from the $z$-axis
\begin{equation}
\theta=\arctan\left(\dfrac{2g}{\omega_{B}}|\bar{a}_1||\bar{a}_2|\right),\label{eq:deflectionAngle}
\end{equation}
which depends on the square root of the circulating photon number in each optical mode.

In the limit that the optically induced magnetic field is smaller
than the external magnetic field, i.e. $\omega_{B}\gg4g\vert\bar{a}_{1}\vert\vert\bar{a}_{2}\vert$,
the optical steady state reduces to
\begin{align}
\bar{a}_{1}= & -\frac{\epsilon_{L1}}{\mathrm{i}\left(\Delta_{1}+\frac{2g^{2}\vert\bar{a}_{2}\vert^{2}S}{\omega_{B}}\right)-\frac{1}{2}\kappa},\label{eq:a1_steadyState}
\end{align}
where $\bar{a}_{2}$ is given by 
\begin{equation}
\begin{aligned}\mathrm{i}\Delta_{2}\bar{a}_{2}+\mathrm{i}\frac{2g^{2}S}{\omega_{B}}\frac{\epsilon_{L1}^{2}}{\left(\Delta_{1}+\frac{2g^{2}\vert\bar{a}_{2}\vert^{2}S}{\omega_{B}}\right)^{2}+\frac{1}{4}\kappa^{2}}\bar{a}_{2}\\
-\frac{1}{2}\kappa\bar{a}_{2}+\epsilon_{L2}=0.
\end{aligned}
\label{eq:a2_steadyState}
\end{equation}
Note that it is equivalent to first express $\bar{a}_2$ as a function of $\bar{a}_1$ and then solve the corresponding equation for $\bar{a}_1$. By driving only one optical mode, the steady state of the undriven mode equals zero. As a consequence, the steady state value of the spin is unaffected by the light, as can be seen from Eqs.$~$ \eqref{eq:azimuthalDeflection} and \eqref{eq:deflectionAngle}, and the effect of the driving is purely dynamical. In the following we focus on a pumping scheme where both modes are driven (Fig. \ref{fig:schematic_model}(c)).

For the range of laser powers considered (up to $1\,\mathrm{W}$)
the deflection angle scales linearly with the power and we obtain
$\theta\approx0.2\sqrt{P_{1}P_{2}}{}^{\circ}/W$ for a sphere of
size $R=1\,\mu\mathrm{m}$ (taking $\Delta_{1}=-\omega_{B}$, $\Delta_{2}=0$).
Such small deflection angle justifies the static Mie splitting determined by $M_\mathrm{S}$
along the $z$-axis (see Sec. \ref{sec:Model}) also in the
presence of light. 

The steady states given in Eq.$~$\eqref{eq:possible_steadyState} are not necessarily
stable. In particular for $\mathbf{e}_{3}$ an easy-axis (hard-axis),
the perpendicular (parallel) steady state configuration is always
unstable. Further, we observe in the following that the dynamics do
not fully converge to the steady states in the considered timescales due
to the low dissipation in the system. 

\section{Dynamics \label{sec:Force-and-dynamical}}

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{dynamics} 
\caption{Dynamics of a sphere of size $R=1.5\,\mu\mathrm{m}$ with initial
configuration $\mathbf{S}\parallel\mathbf{e_{3}}\parallel\mathbf{e}_{z}$
and $\mathbf{L}=\bm{0}$. The system is red detuned with $P_{2}=0.02\,\mathrm{W}$
and $P_{1}=P_{2}/2$. (a) Amplitude fluctuation of the cavity mode
1. The inset shows the initial time evolution of the system. After
a transient period, the cavity dynamics is modulated by (b) the dynamics
of the anisotropy axis $\mathbf{e}_{3}$. Angular momentum (c) parallel
$\mathbf{L}_{\parallel}$ and (d) perpendicular $\mathbf{L}_{\perp}$
to $\mathbf{e}_{3}$. As schematically depicted in (e), $\mathbf{L}_{\parallel}$
gives the angular speed of rotation of the sphere around the anisotropy
axis via $\hbar$$\vert\mathbf{L}_{\parallel}\vert/I$ and $\mathbf{L}_{\perp}$
is related to the libration and precession of $\mathbf{e}_{3}$, shown
in the $xy$ plane (f). The time is given in terms of the magnon frequency
$\omega_{B}=10^{10}\,\mathrm{s^{-1}}$. }
\label{fig:a+e3}
\end{figure}

Having characterized the steady states of the system we now turn to
study numerically the dynamics of the system under an optical drive. For concreteness we consider an easy-axis magnet, i.e. $D>0$, and
an initial state in which the spin and the anisotropy axis are aligned
parallel along the external magnetic field, namely $\mathbf{e}_{3}(0)=\mathbf{S}(0)/S=\mathbf{e}_{z}$. For a hard-axis magnet we would retrieve the same results by choosing
an initial configuration in which the anisotropy axis is perpendicular
to the spin and therefore do not explicitly discuss them. For the
optical fields we assume that the initial conditions are given by
Eqs.$~$\eqref{eq:a1_steadyState} and \eqref{eq:a2_steadyState}. We
consider driving laser powers such that the deviation of the steady state with respect to the chosen
initial configuration is sufficiently small.
Thus, we expect the system to oscillate around this configuration
with an amplitude that depends on the driving power. Since we consider a Mie splitting given by the static magnetization along the $z$-axis, our model requires modifications for initial configurations
where the spin is not approximately along the $z$ direction, and it breaks down for pumping powers that exceed a critical value such that
the $z$ component of the spin deviates too far from the $z$-axis. In this case the so-called Suhl instabilities can further occur, opening new dissipation channels \citep{suhl1957theory}. 

We consider two pumping schemes depicted in Fig.$~$\ref{fig:schematic_model}(c):
First, the optical cavity-mode with higher frequency is pumped in resonance, such
that $\Delta_{1}=0$ and $\Delta_{2}=\omega_{B}$, which we refer
to as blue detuning and describes magnon emission; second, the mode
with lower frequency is pumped in resonance, such that $\Delta_{1}=-\omega_{B}$
and $\Delta_{2}=0$, which is related to a magnon absorption and which
we call red detuning. These processes can be interpreted as cavity-enhanced
Stokes and anti-Stokes Brillouin scattering. The Gilbert damping included
in the equations of motion Eq.$~$\eqref{eq:eom} renders the north pole
a stable equilibrium point in the absence of driving. With driving,
for the red detuning scheme, the north pole is stable whereas the
south pole results in runaway solutions which can lead to limit cycles
or chaotic dynamics. For blue detuning, the south pole can become
a stable point depending on the driving power and on the detuning.
For a more detailed discussion see Ref.$~$\citep{violakusminskiyCoupledSpinlightDynamics2016a}.
Note that we choose $P_{1}=P_{2}/2$ for our simulations.  This is
not a restrictive assumption and changing the ratio between
$P_{1}$ and $P_{2}$ does not change the results qualitatively. 

We show the time evolution of $\vert\delta a_{1}(t)\vert=\vert\bar{a}_{1}-a_{1}(t)\vert$,
i.e. the deviations of the amplitude of the higher frequency optical
mode from its steady state value, and of $\mathbf{e}_{3}$ in Figs.$~$\ref{fig:a+e3}(a)
and (b) for the red detuned pumping scheme. The initial configuration
is such that the particle has no initial rotational motion, i.e. $\mathbf{L}(t=0)=\mathbf{0}$.
We shortly discuss the dynamics of an initially rotating particle,
$\mathbf{L}(t=0)\neq\mathbf{0}$, in Appendix \ref{sec:Further-Dynamics}.
This latter case is especially relevant if one starts with an unmagnetized
particle at rest and then magnetizes it, which induces an initial
rotation due to the Einstein-de Haas effect \citep{einstein_de_haas}.

The cavity field dynamics are modulated by the angular motion and,
at short times, the dominant contributions are due to the optomagnonic
coupling and the fast spin dynamics. To gain further insight into
the driven angular motion, we show in Fig.$~$\ref{fig:a+e3}(c) the
dynamics of the absolute value of the components of $\mathbf{L}$
parallel and perpendicular to $\mathbf{e}_{3}$. The parallel component
of the angular momentum $\mathbf{L}_{\parallel}$ gives the angular
speed of rotation of the sphere around its magnetic anisotropy axis via $\hbar$$\vert\mathbf{L}_{\parallel}\vert/I$
in the lab frame, while the perpendicular component $\mathbf{L}_{\perp}$
is related to libration and precession of $\mathbf{e}_{3}$. We notice
that $\vert\mathbf{L}_{\parallel}\vert$ quickly increases but then
stays stationary. This is due to the quick spin dynamics during the
initial evolution from the starting configuration to close-to-steady-state
oscillation. After the transient dynamics the effects of spin damping,
as described by the Gilbert damping term, become negligible and $\vert\mathbf{L}_{\parallel}\vert$
is conserved.

\section{Power Spectrum Analysis\label{sec:Power-Spectrum}}

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{cavity_red}
\caption{Logarithmic power spectrum of the cavity optical modes (a) $\delta a_{1}$
and (b) $\delta a_{2}$ for a red detuned system with $P_{2}=0.02\,\mathrm{W}$
and $P_{1}=P_{2}/2$. The sphere has a radius of $R=1.5\,\mu\mathrm{m}$,
the initial state is $\mathbf{S}\parallel\mathbf{e_{3}}\parallel\mathbf{e}_{z}$,
and $\mathbf{L}=\mathbf{0}$, corresponding to the results presented in Fig.$~$\ref{fig:a+e3}(a).
Both modes yield a peak in the spectrum at a frequency of $\sim90\,\mathrm{kHz}$
(highlighted by the blue arrow) which corresponds to the induced angular
motion. At a frequency of $\sim2.33\,\mathrm{GHz}$ a second peak
can be observed which corresponds to the fast oscillations at small
times (orange arrow), related to the magntization's transient dynamics. Note that the horizontal axis in the right panel is in linear scale.
The frequency is given in terms of the magnon frequency $\omega_{B}=10^{10}\,\mathrm{s^{-1}}$.}
\label{fig:cavity_red}
\end{figure}

As we observed in the previous section, the dynamics of the cavity optical modes exhibit a transitory behavior followed by an evolution modulated
by the anisotropy axis dynamics. Namely, at short times the deviations
$\vert\delta a_{i}(t)\vert$ oscillate with the same frequency as
the spin oscillates and then adapts to the oscillations of the magnetic anisotropy
axis. The frequencies associated with each of these dynamics are imprinted
in the power spectrum of the light, as we show now. 

First, we write the cavity fields as
\begin{equation}
a_{i}(t)=\bar{a}_{i}+\delta a_{i}(t),\label{eq:cavity_field}
\end{equation}
such that $\delta a_{i}(t)$ represents the deviation of $a_{i}$
from its steady state value $\bar{a}_{i}$. We then consider the Fourier
transform of $\delta a_{i}(t)$
\begin{equation}
\mathcal{F}[\delta a_{i}(t)]=\int_{0}^{+\infty}\delta a_{i}(t)\mathrm{e}^{\mathrm{i}\omega t}\mathrm{d}t=\delta a_{i}(\omega),\label{eq:Fourier_trafo}
\end{equation}
from which we compute the power spectrum $|\delta a_{i}(\omega)|^{2}$. 

We show in Fig.$~$\ref{fig:cavity_red} the absolute value of the
fluctuations and the resulting spectrum for an initially non-rotating
sphere ($R=1.5\,\mu\mathrm{m}$) under a red detuned drive. The spectra
for both modes show a clear peak at a frequency of $\sim90\,\mathrm{kHz}$.
This peak corresponds to the oscillation frequency of the cavity field
amplitude after the initial transitory behavior. In fact, this frequency
is the same as the oscillation frequency of the magnetic anisotropy axis  and
is a fingerprint of the modulation of the cavity dynamics by the sphere's
angular motion. Further, the spectra yield a second peak at a frequency
$\sim2.33\,\mathrm{GHz}$, that is related to optomagnonic induced
transient spin oscillations including the effects of its coupling to the
angular motion (a frequency shift). The finite linewidth of the Lorentzian shaped peaks is governed by the low dissipation and the finite simulation time. The dip around $0.4\,\mathrm{GHz}$ in the power spectrum
of mode 1 can be attributed to a Fano resonance \citep{fano1961effects}
emerging from the coupling between spin and angular momentum which,
close to the steady state, act like two coupled harmonic oscillators.
Such Fano resonances have also been observed in optomechanical systems
\citep{wang2020optomechanically,qu2013fano}.

To better understand the origin of the observed resonance frequencies,
we consider now the linearized version of the equations of motion Eq.
\eqref{eq:eom} and derive the susceptibility of the cavity-field
fluctuations. For that, we consider small displacements from
the steady state 
\begin{align*}
\boldsymbol{\xi}(t) & =\bar{\boldsymbol{\xi}}+\delta\boldsymbol{\xi}(t),
\end{align*}
and retain in the equations of motion only linear terms in the fluctuations.
Since the considered optomagnonic effects are weak, to a first approximation
we assume that the steady state of the spin and the anisotropy axis
is sufficiently close to the north pole, i.e. $\bar{S}_{x,y}/S=\bar{e}_{3x,y}\approx0$
and $\bar{S}_{z}/S=\bar{e}_{3z}\approx1$, and take a non-rotating
sphere, $\mathbf{\bar{L}}=\mathbf{0}$. Relegating the detailed calculations
to the Appendix$~$\ref{sec:Resonance-Peak}, we obtain the following
susceptibility for the fluctuations $\delta a_1$ in frequency space: 
\begin{equation}
\chi_{1}(\omega)=\left[-\mathrm{i}\left(\omega+\Delta_{1}-2g^{2}S\vert\bar{a}_{2}\vert^{2}\chi_{+}(\omega)\right)+\frac{\kappa}{2}\right]^{-1},\label{eq:a1(omega)}
\end{equation}
where 
\begin{align*}
\chi_{+}^{-1}(\omega) & =\omega-\tilde{\omega}_{B}-\dfrac{|\bar{B}_{+}|^{2}}{2\omega}+\dfrac{4\omega_{D}^{2}\omega_{I}S^{3}}{2\omega_{D}\omega_{I}S^{2}-\omega^{2}}\\
 & +\dfrac{2g^{2}S|\bar{a}_{1}|^{2}\left(\omega-\Delta_{2}-\mathrm{i}\kappa/2\right)}{\kappa^{2}/4+(\omega-\Delta_{2})^{2}}-\dfrac{\bar{B}_{+}^{2}\bar{B}_{+}^{\ast2}}{4\omega^{2}}\chi_{-}(\omega),\\
\chi_{-}^{-1}(\omega) & =\omega+\tilde{\omega}_{B}-\dfrac{|\bar{B}_{+}|^{2}}{2\omega}-\dfrac{4\omega_{D}^{2}\omega_{I}S^{3}}{2\omega_{D}\omega_{I}S^{2}-\omega^{2}}\\
 & +\dfrac{2g^{2}S|\bar{a}_{1}|^{2}\left(\omega+\Delta_{2}-\mathrm{i}\kappa/2\right)}{\kappa^{2}/4+(\omega+\Delta_{2})^{2}}\\
 & -\dfrac{2g^{2}S|\bar{a}_{2}|^{2}\left(\omega-\Delta_{1}-\mathrm{i}\kappa/2\right)}{\kappa^{2}/4+(\omega-\Delta_{1})^{2}},
\end{align*}
with $\tilde{\omega}_{B}=\omega_{B}+2\omega_{D}S$ and $\bar{B}_{+}=2g\bar{a}_{1}\bar{a}_{2}^{\ast}$.
The peaks of the power spectrum are obtained by the poles of the imaginary part of $\chi_{1}[\omega]$, and are approximately given by the equation 
\begin{align}
(\omega+\Delta_{1})\left(\omega-\tilde{\omega}_{B}-\frac{|\bar{B}_{+}|^{2}}{2\omega}+\dfrac{4\omega_{D}^{2}\omega_{I}S^{3}}{2\omega_{D}\omega_{I}S^{2}-\omega^{2}}\right)\nonumber \\
-2Sg^{2}\vert\bar{a}_{2}\vert^{2}=0.\label{eq:peak_frequencies}
\end{align}
For the values used in Fig.$~$\ref{fig:cavity_red}, we obtain from
Eq.$~$\eqref{eq:peak_frequencies} two resonance frequencies: $\omega_{L}/(2\pi)\approx88.7\,\mathrm{kHz}$
and $\omega_{H}/(2\pi)\approx2.33\,\mathrm{GHz}$, in agreement with
the numerical results shown in Fig.$~$\ref{fig:cavity_red}.

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{phase_diagrams2}
\caption{Phase diagrams for a red detuned system. Values of the (a) low $\omega_{L}$
and (b) high $\omega_{H}$ resonance frequency depending on the radius
$R$ and the magnetic field $B_0$ for $P_{2}=0.02\,\mathrm{W}$ and
$P_{1}=P_{2}/2$. The inset shows the dependence on $B_0$ and the laser
power $P_{2}$ for an increased optomagnonic coupling $g'=10g$ and
$R=1.5\,\mu\mathrm{{m}}$. The resonance frequencies are obtained
from the imaginary part of the eigenvalues of the linearized system Eq.$~$\eqref{eq:eom_linearized}. }
\label{fig:phase_diagrams}
\end{figure}

The high-frequency peak $\omega_{H}$ can be obtained by discarding the terms $\dfrac{4\omega_{D}^{2}\omega_{I}S^{3}}{2\omega_{D}\omega_{I}S^{2}-\omega^{2}}$
and $|\bar{B}_{+}|^{2}/(2\omega)$ in Eq.$~$\eqref{eq:peak_frequencies}.
This corresponds to neglecting the dynamics of $\mathbf{e}_{3}$ (a
good approximation for describing the dynamics of the system for periods
$\ll1/\omega_{I}$), and to considering the optically induced magnetic
field negligible compared to the applied external magnetic field.
Under these approximations, the high-frequency peak is given by
$\omega_{H}\approx\tilde{\omega}_{B}$. The coupling between spin
and angular motion changes the spin precession frequency by a factor
$2\omega_{D}S$ compared to the uncoupled case \citep{kittel1948theory}. Note that this approximation for the shift in the spin precession frequency is volume independent and only depends on the material properties
like the anisotropy constant or the density. This can also be seen
in Fig.$~$\ref{fig:phase_diagrams} (for radii up to $R\approx2.5\,\mu\mathrm{{m}}$)
where we show the dependence of the resonance frequencies on the particle
radius and the external magnetic field. Regarding its dependence on
the laser power, we observe that in principle $\omega_{H}$ is shifted
to higher frequencies for larger powers. However, for the considered
values of the optomagnonic coupling, this effect is negligible. Nonetheless,
the dependence is shown in the inset of Fig.$~$\ref{fig:phase_diagrams}(b)
for an increased coupling $g'=10g.$ In this case we further find
that an additional peak close to the lower frequency peak appears
in the spectrum for power values $P_{2}>0.1\,\mathrm{{W}}$. The same
can be observed for the unchanged optomagnonic coupling strength for
smaller magnon frequencies. However, this peak is weak in amplitude
and thus we do not further investigate it. 

An initial rotation does not change the power spectrum qualitatively
but decreases the frequency of the lower peak $\omega_L$. Such effects
are typically small for $L\sim S$, which corresponds to frequencies
of $\sim100/(2\pi)\,\mathrm{Hz}$. For $L\gg S$, rotations
yield a more pronounced frequency shift, for instance, for $L=10^{3}S$
($L=10^{4}S$) the position of the lower peak is at $\sim80\,\mathrm{kHz}$
($\sim40\,\mathrm{kHz}$). For large $L$ gyroscopic effects become
more dominant such that the torque due to the fast angular rotation
stabilizes the anisotropy axis and thus decreases its oscillation
frequency. As mentioned before, the discussion presented here is valid
as long as the initial state of the system is close to its steady
state. %
\begin{comment}
The induced torque could be probed by measuring the polarization state
of the transmitted laser beam. The angular momentum transfer renders
the balance of left ($+\hbar$) and right ($-\hbar$) circular components
of the light such that the transmitted beam has a net angular momentum.
The induced torque is determined by difference of the power of left
and right circular components that is measurable with a torque detector,
as demonstrated in \citep{la2004optical}. 
\end{comment}

 

\section{RF Driven System \label{sec:Driven-System}}

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{rotation_driven3}
\caption{(a) Illustration of the driven system. An oscillating magnetic field
is applied perpendicular to the bias field $B_0$. (b) Angular speed of rotation
of the sphere around the anisotropy axis. The system is red detuned
with $R=1.5\,\mu\mathrm{m}$, $P_{2}=0.02\,\mathrm{W}$ and $P_{1}=P_{2}/2$
and initial configuration $\mathbf{S}\parallel\mathbf{e_{3}}\parallel\mathbf{e}_{z}$
and $\mathbf{L}=\bm{0}$. The amplitude of the driving magnetic field
is given by $a_{\mathrm{rf}}=0.01$. The angular speed reached after
$\sim20\,\mathrm{ms}$ is 12 orders of magnitude higher than the corresponding
rotation of the undriven system shown in Fig.$\,$\ref{fig:a+e3}(c).
The time is given in terms of the magnon frequency $\omega_{B}=10^{10}\,\mathrm{s^{-1}}$.}
\label{fig:rotationSpeed_driven}
\end{figure}
To enhance the peak in the power spectrum related to the spin oscillations
as well as the angular speed of rotation of the sphere, we take inspiration
from ferromagnetic resonance experiments where a periodic magnetic
field is applied perpendicular to the bias field, as schematically
depicted in Fig. \ref{fig:rotationSpeed_driven}(a). A resonant drive
corresponds to an oscillating driving field with frequency matching
the magnetization's precession frequency. We thus include in our model
a periodic magnetic field $\mathbf{B}_{\mathrm{rf}}=a_{\mathrm{rf}}B_0\cos(\omega_{\mathrm{rf}}t)\mathbf{e}_{y}$
where $a_{\mathrm{rf}}$ is a dimensionless factor describing the
ratio between the FMR drive amplitude and the bias field strength
(given by $B$). The Hamiltonian Eq.$~$\eqref{eq:hamiltonian_reduced} is
then augmented by the term
\begin{align}
\hat{H}_{\mathrm{\mathrm{rf}}}= & -\hbar a_{\mathrm{rf}}\omega_{B}\hat{S}_{y}\cos(\omega_{\mathrm{rf}}t).\label{eq:Hamiltonian_driven}
\end{align}
The amplitude of the driving field can be related to the excitation
power $P_{\text{mw}}$ by $a_{\mathrm{rf}}\omega_{\mathrm{rf}}=2\sqrt{\kappa_{\text{mw}}P_{\text{mw }}/(\hbar\omega_{\text{mw}})}$
with the microwave frequency $\omega_{\mathrm{mw}}$ and coupling to the antenna
$\kappa_{\mathrm{mw}}$. The value of the latter can be kept close
to the intrinsic magnon dissipation $\kappa_{m}\sim1\,\mathrm{MHz}$
\citep{tabuchi2016quantum,rameshti2021cavity} to enable a maximal
power flow into the magnet. Thus, for $\omega_{\mathrm{mw}}/(2\pi)=10\,\mathrm{GHz}$
this corresponds to $P_{\mathrm{mw}}\sim a_{\mathrm{rf}}^{2}\,\mathrm{pW}$
which is comparable to excitation powers used in Ref.$~$\citep{lachance2017resolving}.

This external drive not only enhances the visibility of the peaks
in the cavity spectrum but also introduces new features depending
on the frequency of the drive. Fig.$~$\ref{fig:powerSpectrum_driven}(a)-(c)
shows the power spectrum of cavity mode $1$ for different FMR frequencies $\omega_{\mathrm{rf}}$ which are visible by a corresponding peak.
If we choose the FMR frequency $\omega_{\mathrm{rf}}$ to match the
shifted spin precession frequency $\omega_{H}\approx\tilde{\omega}_{B}$,
the additional field enhances the amplitude of the higher frequency
peak $\omega_H$ and yields a third peak with intermediate frequency between $\omega_L$ and $\omega_H$. This extra resonance made visible by the rf driving field is due to a mismatch between the
optical modes frequency difference and the shifted magnetization precession
frequency due to the coupling to angular motion. This is a consequence of the nonlinearity of the system and is not captured by the analytically calculated susceptibilities, so we have to resort to numerical analysis. The position of the appearing peak depends on the driving
field amplitude (Fig.$~$\ref{fig:powerSpectrum_driven}(d)), namely,
for larger amplitudes this peak shifts to higher frequencies.

\begin{figure}[t]
\includegraphics[width=1\columnwidth]{powerSpectrum_driven}
\caption{Logarithmic power spectrum of the driven system of $\delta a_{1}$
for $a_{\mathrm{rf}}=10^{-3}$ and FMR frequency (a) $\omega_{\mathrm{rf}}=0.01\tilde{\omega}_{B}$,
(b) $\omega_{\mathrm{rf}}=0.1\tilde{\omega}_{B}$, (c) $\omega_{\mathrm{rf}}=\tilde{\omega}_{B}$.
The system is red detuned with $P_{2}=0.02\,\mathrm{W}$, $P_{1}=P_{2}/2$
and $R=1.5\,\mu\mathrm{m}$. The initial state is $\mathbf{S}\parallel\mathbf{e_{3}}\parallel\mathbf{e}_{z}$,
and $\mathbf{L}=\mathbf{0}$. The FMR frequency is imprinted in the spectrum
by a corresponding peak (marked by the green, dotted arrow). For $\omega_{\mathrm{rf}}=\tilde{\omega}_{B}$,
additionally to the peaks at $\sim100\,\mathrm{kHz}$ and $\sim2.33\,\mathrm{GHz}$,
the power spectrum yields a third peak (red, dashed arrow) that depends
on the amplitude of the driving field (d). The driving field amplitude
values are $a_{\mathrm{rf}}=10^{-4},10^{-3},10^{-2},10^{-1}$ from
bottom to top. The frequency is given in terms of the magnon frequency
$\omega_{B}=10^{10}\,\mathrm{s^{-1}}$.}
\label{fig:powerSpectrum_driven}
\end{figure}

The driving magnetic field also affects the angular motion of the
particle. In fact, it increases the angular speed of rotation of the
sphere around the anisotropy axis. This can be seen in Fig. \ref{fig:rotationSpeed_driven}(b)
where we show $\hbar \vert\mathbf{L}_{\parallel}\vert/I$ in the
red detuned scheme for the initial configuration $\mathbf{e}_{3}(0)=\mathbf{S}(0)/S=\mathbf{e}_{z}$,
and $\mathbf{L}(t=0)=\mathbf{0}$. We observe that the induced rotation is
always counterclockwise with respect to $\mathbf{\mathbf{e}_{3}}$.
Thus, an initially clockwise rotating particle will be decelerated.
The effects of rotation are only imprinted in the optical power
spectrum for high rotation speeds as discussed in the previous section. 

\section{Generalizations \label{sec:Generalisations}}

In this section we briefly discuss some possible generalizations
to the Hamiltonian Eq.$~$\eqref{eq:hamiltonian_reduced} of our model. 

\subsection{Cotton-Mouton Effect\label{sec:Cotton-Mouton-Effect}}

So far we considered that only the Faraday effect plays a role in
optomagnonic coupling. Another effect relevant in several optomagnonic
systems is the Cotton-Mouton effect \citep{stancil2009spin,osada2018orbital},
which we now include, discussing its possible impacts on our results. 

The Cotton-Mouton (CM) effect is quadratic in the magnetization and
is taken into account in the full effective permittivity tensor \citep{wettling1975relation} 
\begin{equation}
\varepsilon_{ij}(\mathbf{M})=\varepsilon_{0}\left(\varepsilon\delta_{ij}-\mathrm{i}f\sum_{k}\epsilon_{ijk}M_{k}+\sum_{kl}G_{ijkl}M_{k}M_{l}\right).\label{eq:permettivity_tensor}
\end{equation}
Notice that this includes a term $\propto M_{k}M_{l}$ which was not
included in our previous formalism {[}cf. Eq. \eqref{eq:permittivity_tensor_Faraday}{]}.
For cubic crystals, such as YIG, most of the tensor components $G_{ijkl}$
are zero by symmetry and the non-vanishing ones are \citep{pisarev1971magnetic,wettling1975relation}
$G_{iiii}=G_{11}$, $G_{iijj}=G_{12}$, $G_{ijij}=G_{44}$ for $i\neq j$
and $i,j\in\{1,2,3\}$.

The Hamiltonian obtained by quantizing the CM contribution to the
energy density is given by (for details see Appendix \ref{sec:CM-effect-app}) 
\begin{equation}
\begin{aligned}\hat{H}_{\mathrm{int}}^{\mathrm{CM}}= & \hbar g_{\mathrm{CM}}^{44}\left(\hat{S}_{z}\hat{S}_{+}\hat{a}_{1}^{\dagger}\hat{a}_{2}+\hat{S}_{-}\hat{S}_{z}\hat{a}_{2}^{\dagger}\hat{a}_{1}\right)\\
 & +\hbar g_{\mathrm{CM}}^{12}\hat{S}_{z}^{2}(\hat{a}_{1}^{\dagger}\hat{a}_{1}-\hat{a}_{2}^{\dagger}\hat{a}_{2}),
\end{aligned}
\label{eq:H_Int_CM_final-1}
\end{equation}
where we neglected the constant shift in the photon field. The coupling
constants $g_{\mathrm{CM}}^{44}$ and $g_{\mathrm{CM}}^{12}$ are
given by 
\begin{align}
g_{\mathrm{CM}}^{44} & =\dfrac{\varepsilon_{0}}{2\hbar}\left(\dfrac{M_{\mathrm{S}}}{S}\right)^{2}G_{44}\tilde{g},\nonumber \\
g_{\mathrm{CM}}^{12} & =\dfrac{\omega_{r}}{8\epsilon}\left(\dfrac{M_{\mathrm{S}}}{S}\right)^{2}(G_{11}-G_{12}),\label{eq:gCM}
\end{align}
under the approximation $\int\mathrm{d}\mathbf{r}|E_{h}|^{2}\approx\int\mathrm{d}\mathbf{r}|E_{v}|^{2}\approx\dfrac{\hbar\omega_{r}}{2\varepsilon_{0}\varepsilon}$. 

The equations of motion for the spin and the cavity fields including
the additional interacting term Eq.$~$\eqref{eq:H_Int_CM_final-1} read
\begin{align}
\dot{\mathbf{S}}= & \mathbf{\tilde{B}}_{{\rm eff}}\times\mathbf{S}-2\omega_{D}\left(\mathbf{e}_{3}\cdot\mathbf{S}\right)\left(\mathbf{e}_{3}\times\mathbf{S}\right)+\tau_{\mathrm{G}},\nonumber \\
\dot{a}_{1}= & \mathrm{i}(\Delta_{1}-g_{\mathrm{CM}}^{12}S_{z})a_{1}-\mathrm{i}(g+g_{\mathrm{CM}}^{44}S_{z})S_{+}a_{2}+\epsilon_{L1}-\dfrac{1}{2}\kappa a_{1},\nonumber \\
\dot{a}_{2}= & \mathrm{i}(\Delta_{2}+g_{\mathrm{CM}}^{12}S_{z})a_{2}-\mathrm{i}(g+g_{\mathrm{CM}}^{44}S_{z})S_{-}a_{1}+\epsilon_{L2}-\dfrac{1}{2}\kappa a_{2},\label{eq:eom_CM}
\end{align}
where the effective field is now given by $\tilde{\mathbf{B}}_{\mathrm{eff}}=\mathbf{B}_{\mathrm{eff}}+\mathbf{B}_{\mathrm{opt}}^{\mathrm{CM}}$
with
\begin{align}
\mathbf{B}_{\mathrm{opt}}^{\mathrm{CM}}= & g_{\mathrm{CM}}^{44}\left(\begin{array}{c}
(a_{1}^{\dagger}a_{2}+a_{2}^{\dagger}a_{1})S_{z}\\
\mathrm{i}(a_{1}^{\dagger}a_{2}-a_{2}^{\dagger}a_{1})S_{z}\\
(a_{1}^{\dagger}a_{2}+a_{2}^{\dagger}a_{1})S_{x}+\mathrm{i}(a_{1}^{\dagger}a_{2}-a_{2}^{\dagger}a_{1})S_{y}
\end{array}\right)\nonumber \\
 & -2g_{\mathrm{CM}}^{12}\left(\begin{array}{c}
(a_{1}^{\dagger}a_{1}-a_{2}^{\dagger}a_{2})S_{x}\\
(a_{1}^{\dagger}a_{1}-a_{2}^{\dagger}a_{2})S_{y}\\
0
\end{array}\right).\label{eq:Bopt_CM}
\end{align}

For YIG, $G_{44}M_{\mathrm{S}}^{2}=-1.14\cdot10^{-4}$ and $(G_{11}-G_{12}-2G_{44})M_{\mathrm{S}}^{2}=5.73\cdot10^{-5}$
\citep{stancil2009spin}. For small deflections of the spin from the
$z$ axis, i.e. $S_{z}\sim S$ and $S_{x,y}\sim0$, the additional
contribution to the optically induced magnetic field is of the same
structure and magnitude as the one caused by the Faraday effect, considering
that the value of $g_{\mathrm{CM}}^{44}$ mainly differs from $g$
given by Eq.$~$\eqref{eq:g} by a factor of $S^{-1}$. Therefore, the
inclusion of Cotton-Mouton effect terms do not change qualitatively
the results presented in the main text. 

\subsection{Anisotropy \label{subsec:Cubic-anisotropy}}

Some crystal systems, including YIG, exhibit cubic symmetry for which
the magnetocrystalline anisotropy energy can be written as (for a
homogeneous material) \citep{chikazumi2009physics,cullity2008introduction,dionne2009anisotropy,landau1984electrodynamics}
\begin{align}
E_{c}= & K_{1}V\left[(\mathbf{s}\cdot\mathbf{e}_{1})^{2}(\mathbf{s}\cdot\mathbf{e}_{2})^{2}+(\mathbf{s}\cdot\mathbf{e}_{1})^{2}(\mathbf{s}\cdot\mathbf{e}_{3})^{2}\right.\nonumber \\
+ & \left.(\mathbf{s}\cdot\mathbf{e}_{2})^{2}(\mathbf{s}\cdot\mathbf{e}_{3})^{2}\right]+K_{2}V(\mathbf{s}\cdot\mathbf{e}_{1})^{2}(s\cdot\mathbf{e}_{2})^{2}(s\cdot\mathbf{e}_{3})^{2},\label{eq:cubic_anisotropy}
\end{align}
where $K_{1}$ and $K_{2}$ denote the respective first- and second-order
cubic anisotropy constants and $\mathbf{s}=\mathbf{S}/S$. As the
term $\propto K_{2}$ is of 6th order it can be safely neglected.
If we assume only small deviations of the spin $\ensuremath{\mathbf{s}}$
from one of its equilibrium directions, e.g. $\ensuremath{\mathbf{s}\cdot\mathbf{e}_{3}}\sim1$,
and keep the other terms up to second order, Eq.$~$\eqref{eq:cubic_anisotropy}
can be approximated as
\begin{equation}
\begin{aligned}E_{c}\approx & K_{1}V[(\mathbf{s}\cdot\mathbf{e}_{1})^{2}+(\mathbf{s}\cdot\mathbf{e}_{2})^{2}]\\
= & K_{1}V\vert\mathbf{s}\times\mathbf{e}_{3}\vert^{2},
\end{aligned}
\label{eq:cubic_anisotropy_approx}
\end{equation}
where we have used that $\mathbf{e}_{3}=\mathbf{e}_{1}\times\mathbf{e}_{2}$
and the vector triple product expansion. By expressing $\sin\phi=\vert\mathbf{s}\times\mathbf{e}_{3}\vert$
and $\cos\phi=\mathbf{s}\cdot\mathbf{e}_{3}$, $\phi$ being the angle
between $\ensuremath{\mathbf{s}}$ and $\mathbf{e}_{3}$, the cubic
anisotropy term Eq.$~$\eqref{eq:cubic_anisotropy_approx} can be rewritten
as 
\begin{equation}
E_{c}\approx K_{1}V\left[1-(\mathbf{s}\cdot\mathbf{e}_{3})^{2}\right],\label{eq:cubic_anisotropy_approx_final}
\end{equation}
which has the same form as the one for uniaxial symmetry apart from
a constant shift. Therefore, for small oscillations of the magnetization
around one of the directions which minimizes the anisotropy energy
density, cubic crystals can be treated to a good approximation as
crystals with uniaxial anisotropy, as is often considered for YIG
\citep{streib2018damping,elyasi2020resources,pacewicz2019rigorous,rijnierse1975optical}.
Furthermore, other magnetic oxides, such as those incorporating bismuth,
exhibit large uniaxial anisotropies \citep{hansen1985magnetic}.

In case of a non-spherical particle there is an additional source
of magnetic anisotropy due to its shape. The shape of the sample generates
a demagnetizing field which is not equal in all directions \citep{cullity2008introduction}.
The energy density of such geometric anisotropy is given by $\mu_{0}\bm{M}\tilde{N}\bm{M}/2$,
where $\tilde{N}$ is the demagnetization. For example, for a prolate
spheroid with semi-major axis along the $\mathbf{e}_{z}$ direction
one has $\tilde{N}={\rm diag}\left[N_{T},N_{T},1-2N_{T}\right]$ with
$N_{T}$ the demagnetization coefficient. A non-spherical shape entails
changes not only in the magnetic energy but also in the rotational
part of the Hamiltonian. On the one hand, it gives rise to a more
complex inertia tensor and on the other hand the electromagnetic fields
inside the rotating particle cannot be assumed to be independent of
the particle orientation. These effects can be considered to be small
if the shape does not deviate much from a sphere. 

\section{Conclusion\label{sec:Conclusion}}

We proposed an optical method to probe the coupled spin-mechanics
of levitated magnetic microparticles via optomagnonic effects. We
showed that the spin dynamics can be also driven by the optical drive,
inducing almost lossless angular oscillations through the spin-mechanics
coupling. This coupling between magnetization and angular motion can
be probed via the power spectrum of the cavity modes, which exhibit
two main resonance peaks at $\sim100\,$kHz and $\sim2\,$GHz for micrometer-sized particles and an applied magnetic field of $\sim60\,$mT. These are attributed to angular oscillations and spin dynamics
that are shifted by the magnetocrystalline anisotropy. We found that
due to the damping of the magnetization (Gilbert damping) a rotation
of the microparticle around the magnetic anisotropy direction can
be optically induced for an initially non-rotating particle. However,
the induced angular frequency is low ($\sim0.2\,\mu\mathrm{Hz}$). We showed
that this angular frequency can be increased by adding an oscillatory
rf magnetic field perpendicular to the bias magnetic field. In particular,
if the drive is at resonance with the shifted spin precession frequency,
the aforementioned effects are enhanced and the setup can be used
to induce fast rotations of the particle around its magnetic anisotropy
axis. This resonant condition can be identified by the appearance
of an additional peak in the spectrum. 

In general, the optical quality
factor of the particle and the supported driving power limit the driven
angular motion. In experiments the quality factor will be restricted
by additional sources, like surface roughness, and not only by the
considered radiative losses. In order to improve the quality factor
the system can be modified by placing the particle in an external
optical cavity. 

We have studied the classical dynamics of the levitated
magnetic particle coupled to light, in particular neglecting thermal
and shot noise. These, together with the effect of the trapping potential,
need to be included in order to study the dynamics in the quantum regime. Experimentally, levitated
magnetic particles can provide a unique platform to probe transduction
of angular momentum in the quantum regime, and eventually serve as
ultrasensitive torque sensors \citep{la2004optical,ahn2020ultrasensitive}. 

\section*{Acknowledgements}
We acknowledge funding from the Max Planck Society and from the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) through Project-ID 429529648-TRR 306 QuCoLiMa ("Quantum Cooperativity of Light and Matter"). V.W. thanks K. Kustura, A. E. Rubio L\'opez and C. C. Rusconi for fruitful discussions.

\input{OpticalSignatures.bbl}

