# Dynamics of photo-induced ferromagnetism in oxides with orbital degeneracy

Jonathan B. Curtis<EMAIL_ADDRESS>College of Letters and Science,
University of California, Los Angeles, CA 90095, USA  Department of Physics,
Harvard University, Cambridge, MA 02138, USA  Ankit Disa  Applied and
Engineering Physics, Cornell University, Ithaca, NY 14853, USA  Max Planck
Institute for the Structure and Dynamics of Matter, Hamburg, DE  Michael
Fechner  Max Planck Institute for the Structure and Dynamics of Matter,
Hamburg, DE  Andrea Cavalleri  Max Planck Institute for the Structure and
Dynamics of Matter, Hamburg, DE  Clarendon Laboratory, Department of Physics,
Oxford University, Oxford, UK  Prineha Narang<EMAIL_ADDRESS>College of
Letters and Science, University of California, Los Angeles, CA 90095, USA

###### Abstract

By using intense coherent electromagnetic radiation, it may be possible to
manipulate the properties of quantum materials very quickly, or even induce
new and potentially useful phases that are absent in equilibrium. For
instance, ultrafast control of magnetic dynamics is crucial for a number of
proposed spintronic devices and can also shed light on the possible dynamics
of correlated phases out of equilibrium. Inspired by recent experiments on
spin-orbital ferromagnet YTiO3 we consider the nonequilibrium dynamics of
Heisenberg ferromagnetic insulator with low-lying orbital excitations. We
model the dynamics of the magnon excitations in this system following an
optical pulse which resonantly excites infrared-active phonon modes. As the
phonons ring down they can dynamically couple the orbitals with the low-lying
magnons, leading to a dramatically modified effective bath for the magnons. We
show this transient coupling can lead to a dynamical acceleration of the
magnetization dynamics, which is otherwise bottlenecked by small anisotropy.
Exploring the parameter space more we find that the magnon dynamics can also
even completely reverse, leading to a negative relaxation rate when the pump
is blue-detuned with respect to the orbital bath resonance. We therefore show
that by using specially targeted optical pulses, one can exert a much greater
degree of control over the magnetization dynamics, allowing one to optically
steer magnetic order in this system. We conclude by discussing interesting
parallels between the magnetization dynamics we find here and recent
experiments on photo-induced superconductivity, where it is similarly observed
that depending on the initial pump frequency, an apparent metastable
superconducting phase emerges.

##  I Introduction

The idea of using strong optical fields to gain control over phases of quantum
matter is a tantalizing one. Optical control of ferroelectric [1, 2],
structural [3, 4, 5, 6, 7, 8], superconducting [9, 10, 11, 12, 13, 14, 15],
charge density wave [16, 17, 18], and magnetic orders [19, 20, 21, 22, 23, 24,
25, 26, 27, 28, 29, 30, 31, 32, 33, 34] have all been proposed theoretically
or demonstrated experimentally. Essentially, by inducing strongly
nonequilibrium scenarios, one can explore an enlarged nonequilibrium phase
diagram which may allow for the access of novel phenomena and functionalities.

One way to realize such strongly nonequilibrium scenarios is by resonantly
driving the system, inducing coherent oscillations in the Hamiltonian, thereby
breaking time-translation symmetry and driving the system away from the
thermal regime. These coherent oscillations can potentially excite parametric
resonances [35, 36, 37, 38, 39], induce novel topology [40, 41, 42, 43, 44,
45, 46], produce tunable interactions [47, 48, 49, 30], generate non-thermal
correlations [50], and even lead to effective cooling mechanisms [51, 52, 53].
In particular, the dynamics of correlated Mott insulators hosting orbital
degrees of freedom [54] is known to be quite rich, exhibiting spin-orbital
separation [55, 56], tunable exchange [57, 58], and hidden phases [59]
Similarly, ferromagnets with high levels of spin-rotation symmetry can exhibit
many exotic phenomena away from equilibrium [60, 61, 62].

The interplay between spin and orbital fluctuations can lead to dramatic
effects including magnon softening [63], entangled spin-orbital phases [64,
65, 66, 67, 68, 69, 70, 71], and pronounced magnetic fluctuations [72, 73] and
subsequent phase transition [74, 75, 76]. Inspired by recent experiments on
ferromagnetic Mott-insulator YTiO3 (YTO) [19], we examine the nonequilibrium
dynamics of magnons in a model quasi-degenerate orbital system driven out of
equilibrium by coherently oscillating optical phonons.

In YTO the $3d^{1}$ conduction band is formed from the titanium $t_{2g}$ shell
which naively has a three-fold degeneracy enforced by a cubic lattice symmetry
[66]. However, in YTO and many other compounds, this cubic symmetry is broken
at low temperatures by a GdFeO3-type structural distortion, which then lifts
the resulting orbital degeneracy [77, 78, 79, 80, 81, 82, 83]. In this case
the system has non-degenerate, but potentially low-lying orbital excitations
[84, 85, 86, 87, 70, 88]. In equilibrium settings, the lifted orbital
degeneracy leads to a decoupling between spin and orbital excitations and for
most magnetic purposes their cross-coupling can be ignored.

Can these orbital excitations, which are essentially absent from equilibrium
processes, significantly modify the out-of-equilibrium dynamics? This question
is not purely academic; the ability to control magnetic order on ultrafast
timescales [28, 25, 49, 21] is crucial for many spintronic technologies, and
may also help design better ferromagnets which can operate at higher
temperatures more efficiently. We answer this question in the affirmative,
provided the orbitals are relatively low-lying and may come close in energy to
relevant optically driven degrees of freedom, such as infrared-active phonons,
which can reside in the 1-20 terahertz regime. In this case, by judiciously
choosing the parameters of the optical driving applied, one can speed-up,
slow-down, and even reverse the magnetization dynamics. This then paves the
way for novel control routes in quasi-degenerate spin-orbital systems [19, 89]
as well as potentially other systems [20, 34, 49].

The remainder of this paper is structured as follows. In Sec. II we outline
the model system considered, and motivate various parameter choices. Then, in
Sec. III we show how in equilibrium the orbitals essentially serve to provide
a bath for angular momentum for the magnons. In Sec. IV we examine the
dynamics of this system in equilibrium and estimate the equilibrium relaxation
time-scale. In Sec. V we then explore the nonequilibrium dynamics of this
system following a simulated impulsive drive of optical phonons, presenting
the main results of this work. Finally, we discuss the implications of our
results in Sec. VI, where we conclude by discussing interesting parallels
between this system and recent experiments on light-induced superconductivity.
In Appendix A we show how to map the $t_{2g}$ orbital levels into an effective
angular momentum. In Appendix B we present details on the nonequilibrium
Keldysh technique as applied to the Holstein-Primakoff spin-wave expansion,
and in Appendix D we show how to reduce these equations to a simple equation
of motion for the magnon occupation.

##  II Model

Here we introduce a simple model for ferromagnetic spins interacting with a
quasi-degenerate orbital bath. Though this is inspired by YTiO3 (YTO), we
emphasize we are considering a more abstract model, and we expect our results
to be relevant to other high-symmetry ferromagnetic insulators with low-lying
orbital excitations. In particular, we consider a single electron occupying a
low-lying $t_{2g}$ orbital manifold. In the cubic limit there is a large
threefold degeneracy which can lead to pronounced orbital fluctuations. In
reality this cubic degeneracy is lifted by the GdFeO3 structural distortion
which renders the crystal structure orthorhombic and induces a finite crystal
field splitting $\Delta$ between the lowest and next-lowest orbitals on each
Ti site. This leads to a model where each site has an orbital pseudospin-1/2
$\hat{\bm{\tau}}_{j}$ in addition to the actual electron spin
$\hat{\mathbf{S}}_{j}$ on each site.

Figure 1: (a) Splitting of Ti $t_{2g}$ shell into non-degenerate levels by the
GdFeO3 distortion, which are then occupied by a single electron giving
$S=1/2$. Interorbital coherences lead to angular momentum $\mathbf{L}$, which
is largely quenched in equilibrium. (b) Focusing on the lowest two-levels
$|0\rangle$ and $|1\rangle$ we find spin-flip $T_{1}$ processes in the orbital
ground-state obtained from virtual orbital transitions. Spin-orbit coupling
$\lambda$ can lead to a simultaneous orbital excitation along with a spin-
flip. This is then followed by a spin-independent orbital decay with rate
$\Gamma$, shown in the level diagram. Ultimately, the decay rate is governed
by the spectral overlap of the orbital bath (shown on the left schematically)
with the spin-transition, which is small leading to a long-lifetime.

We consider a three-dimensional isotropic ferromagnetic Heisenberg model along
side a local orbital degree of freedom with Hamiltonian

$\hat{H}=-J\sum_{j,\bm{\delta}}\mathbf{\hat{S}}_{j}\cdot\mathbf{\hat{S}}_{j+\bm{\delta}}-J^{z}_{0}\sum_{j,\bm{\delta}}\hat{S}^{z}_{j}\hat{S}^{z}_{j+\bm{\delta}}\\\
+\sum_{j}\frac{\Delta}{2}\hat{\tau}^{3}_{j}+\lambda\mathbf{\hat{L}}_{j}\cdot\mathbf{\hat{S}}_{j}.$
(1)

The first two terms are the isotropic Heisenberg exchange, with $J\sim 2.75$
meV [90] for the case of YTO, and an easy-axis exchange which is chosen to
counter the orbital bath-induced Lamb shift, leading to the renormalized spin-
wave gap which for the case of YTO was estimated to be 0.02 meV, though the
upper bound was qutie a bit larger, of order 0.3 meV [90]. We will consider a
modestly-sized renormalized gap of $\Omega_{0}=0.1$ meV in this work. Here $j$
labels the lattice sites $\mathbf{R}_{j}$ and
${\bm{\delta}}=\mathbf{e}_{x},\mathbf{e}_{y},\mathbf{e}_{z}$ labels the
nearest-neighbors along the three principle axes.

The third term, involving $\hat{\tau}^{3}_{j}$, corresponds to the local
crystal-field excitation gap. There is considerable uncertainty about the
value of this parameter, with theoretical estimates ranging from nearly zero
to over 300 meV. In Ref. [91] it was estimated by Raman scattering that
$2\Delta\sim$ 50 meV, while Ref. [92] found energies closer to $2\Delta\sim$
235 meV. Using resonant inelastic x-ray scattering (RIXS), Ref. [87] found
evidence for collective orbital excitations with a gap of order 120 meV. We
will consider $\Delta=90$ meV here, though more experiments with greater
resolution and sensitivity are probably needed in the case of YTO. The last
term describes the atomic $\mathbf{L}\cdot\mathbf{S}$ spin-orbit coupling,
which leads to a torque on the spin in the presence of an orbital angular
momentum $\mathbf{L}$.

The orbital angular momentum can be obtained by projecting the full three-
dimensional $t_{2g}$ angular momentum onto the lowest crystal-field levels,
leading to the expression

$\mathbf{\hat{L}}_{j}=\mathbf{n}_{j}\hat{\tau}^{2}_{j}.$ (2)

The unit vector $\mathbf{n}_{j}$ is orthogonal to the two participating
orbitals, and characterizes the “soft” axis for orbital angular momentum (see
Appendix A). The operator $\hat{\tau}^{2}_{j}$ characterizes the instantaneous
orbital many-body state and in particular $\hat{\tau}^{2}_{j}$ is odd under
time-reversal (which squares to $+1$ for the $L=1$ orbitals), satisfying the
selection rules. This is coupled to the spin angular momentum by the atomic
spin orbit interaction $\lambda$, which in fact need not be small. Reasonable
estimates place $\lambda\sim$ 15 meV for a light $3d$ transition metal such as
Ti [82]. Finally, we note that $\mathbf{n}_{j}$ will in general point in a
different direction on each of the four Ti sublattices depending on the local
crystal-field environment. For more details, we refer to Appendix B. Finally,
we will assume that due to, e.g. phonons or orbital interactions the orbital
excitation itself obtains a finite $T_{1}$ linewidth $\Gamma$. We estimate
$\Gamma\sim$ 15 meV as well, although is not known with great certainty and
may appear significantly more broad in, e.g a two-orbital spectral function,
which may appear in the Raman and RIXS measurements.

We can imagine that in the magnetically ordered phase, each site has a local
level scheme as illustrated in Fig. 1(a), where we show the crystal-field
splitting of the Ti $t_{2g}$ states and their corresponding spin and orbital
angular momentum. The main idea is that magnetization dynamics is often
intrinsically slow due to the bottleneck associated to transfer of spin
angular momentum in to a bath, such as the orbital. One such route is
illustrated in Fig. 1(b), which shows how in second order perturbation theory
this model can give rise to a finite longitudinal magnetization relaxation
rate. We argue that the phonon dynamics induced by the strong optical pulse
can lead to an acceleration of this relaxation time out of equilibrium,
leading to the possibility of pump-enhanced magnetization dynamics. Note that,
unlike Ref. [58] which considered the impact of dynamics on the superexchange
interactions, we are more concerned here with the impact on the spin-orbit
coupling.

Finally, we comment on the coupling to the drive. In the experiment [19], the
pump was performed using a mid-IR pulse which strongly couples to lattice
degrees of freedom, rather than e.g. an optical pulse which traverses the Mott
gap. This pump was tuned to be resonant with various different infrared active
phonon modes and used to strongly drive these vibrations. Based on ab initio
calculations, we argue that one of the dominant effects of this pump is a
strong modulation of the crystal field matrix and in particular, we find that
for relevant fluences this may lead to a sizeable change in the eigenvector
$\mathbf{n}_{j}$. This in turn leads to a dynamical modulation of the orbital
angular momentum $\mathbf{L}_{j}=\mathbf{n}_{j}(t)\hat{\tau}^{2}_{j}$, which
will now acquire sidebands at twice the phonon frequency.

We model this by writing

$\mathbf{n}_{j}(t)\sim\mathbf{n}_{j}+\delta\mathbf{n}_{j}Q^{2}_{\rm IR}(t).$
(3)

Here $Q_{\rm IR}$ is the generalized coordinate describing the infrared-active
phonon which is directly driven by light (in general there may be different or
multiple modes which are excited depending on the frequency and polarization
used in the pump and the absoprtion spectrum of the material). The coupling to
$Q^{2}_{\rm IR}(t)$ is due to the fact that the orbital angular momentum is a
Raman active transition whereas the infrared active phonon is polar. This is
in fact very important since this will induce oscillations at twice the phonon
frequency. For an $\Omega_{\rm ph}=$ 9 THz optical phonon mode, this leads to
sidebands for the spin-orbit coupling at a frequency of $2\Omega_{\rm
ph}\sim$80 meV. This comes close to the orbital resonance in this model at 90
meV. We will study in particular how the dynamics depends on the drive
frequency $\Omega_{d}$. We now proceed to determine the equilibrium structure
of this model before proceeding on to compute the nonequilibrium dynamics. In
particular, we show that the orbitals can act as a bath for angular momentum
for the spins even in equilibrium.

Though we do not specifically consider YTO, for certain rough estimates of
parameter values and feasibility analysis we have used ab initio calculations
based on YTO. We estimate interactions between phonon and the crystal field
parameters by performing first-principles calculations in the framework of
Density Functional Theory (DFT). All technical details are listed in Appendix
C. Our approach is inspired by Ref. [79, 78], where we first compute the full
DFT bandstructure within the local density approximation (LDA). Using this we
construct localized $t_{2g}$ Wannier-functions using appropriate projectors
according to Ref. [93]. To estimate the modulation of the crystal field
parameters due to the phonon distortion, we performed frozen phonon
computations. Therefore, we recalculated the electronic structure and Wannier-
functions for crystal structures which have been modulated according
eigenvectors of polar eigenmodes ($Q_{\rm IR}$). Tthis allows us to estimate
changes in the crystal Hamiltonian for distinct polar distortions.

##  III Orbital Bath

We now analyze the spin-orbit coupling in equilibrium. To later accommodate
the nonequilibrium calculations, we will implement at the outset the
Schwinger-Keldysh formalism for describing this system. We begin by treating
the orbitals in a Gaussian approximation, valid for small orbital excitation
amplitudes. For details we again refer to Appendix B, though for a more
complete treatment we refer the reader to Ref. [94].

In the Keldysh formalism we have a doubling of the degrees of freedom, which
can be arranged into a “classical” part $\mathbf{S}_{j,cl}$ characterizing the
expectation value, and the “quantum” part $\mathbf{S}_{j,q}$ which
characterizes the fluctuations about the expectation value. Applying this
formalism to the spin-orbit interaction we find a Keldysh action of

$S_{\rm soc}=-\lambda\sum_{j}\int
dt\left[\mathbf{L}_{jq}(t)\cdot\mathbf{S}_{j,cl}(t)+\mathbf{L}_{jcl}(t)\cdot\mathbf{S}_{j,q}(t)\right].$
(4)

This then appears as a contribution to a path integral
$Z=\int\mathcal{D}[\mathbf{S},\mathbf{L}]e^{iS}$ which can be used to generate
nonequilibrium correlation functions, such as the magnetization
$\langle\mathbf{S}_{j}(t)\rangle$.

We remark to the reader that in this expression, $\mathbf{S}_{j\alpha}(t)$
($\alpha=cl,q$) should be understood as a stand-in for an appropriate
representation of the spin operator in terms of a canonical bosonic or
fermionic field. In this work we will focus on the dynamics in the ordered
phase, wherein the operators $\mathbf{S}_{j}(t)$ can be expanded in terms of
the Holstein-Primakoff bosons perturbatively in $1/S$ where $S$ is the spin
length. This is strictly valid only at low-temperatures with $T\ll T_{C}$ and
even then, it suffers from the fact that in YTO $S=\frac{1}{2}$ is small. It
therefore remains an important problem for future studies to extend this
treatment to include the fluctuation regime near $T_{C}$ by, e.g. an expansion
in terms of Schwinger bosons instead, which can better handle the dynamics in
the disordered phase. Nevertheless, we expect that for low magnon densities,
this ought to be at least qualitatively acceptable.

We proceed by integrating out the local orbital angular momentum, treating it
as a bath under a Gaussian approximation. This bath can be characterized by
the correlation functions

$\displaystyle\mathbb{\hat{D}}^{R}(t,t^{\prime})=-i\langle\mathbf{L}_{j,cl}(t)\mathbf{L}_{j,q}(t^{\prime})\rangle=\mathbf{n}_{j}(t)\mathbf{n}_{j}(t^{\prime})D^{R}(t,t^{\prime})$
(5a)
$\displaystyle\mathbb{\hat{D}}^{A}(t,t^{\prime})=-i\langle\mathbf{L}_{j,q}(t)\mathbf{L}_{j,cl}(t^{\prime})\rangle=\mathbf{n}_{j}(t)\mathbf{n}_{j}(t^{\prime})D^{A}(t,t^{\prime})$
(5b)
$\displaystyle\mathbb{\hat{D}}^{K}(t,t^{\prime})=-i\langle\mathbf{L}_{j,cl}(t)\mathbf{L}_{j,cl}(t^{\prime})\rangle=\mathbf{n}_{j}(t)\mathbf{n}_{j}(t^{\prime})D^{K}(t,t^{\prime}).$
(5c)

This is in turn expressed in terms of a scalar dynamical response function
$D(t,t^{\prime})$ and the unit vectors $\mathbf{n}_{j}(t)$, as emphasized in
the second equalities, and elaborated on in Appendix B. In equilibrium, these
are completely determined given knowledge of the orbital spectral function and
the thermal occupation function $\coth(\beta\omega/2)$.

In the Gaussian approximation, we can model the spectral function for
$\hat{\tau}^{2}_{j}$ as that of a damped harmonic oscillator. In particular,
we assume the orbital angular momentum has linear response equations of motion
of

$\displaystyle\frac{d\tau^{2}}{dt}=\Delta\tau^{1}-\Gamma\tau^{2}$ (6a)
$\displaystyle\frac{d\tau^{1}}{dt}=-\Delta\tau^{2}-F_{\rm ext}(t).$ (6b)

$F_{\rm ext}(t)$ is an external force which acts on the angular momentum
$\mathbf{L}\propto\hat{\tau}^{2}$, which is canonically conjugate to
$\hat{\tau}^{1}$ (the $x$ Pauli matrix whose expectation value corresponds to
an interorbital density rather than angular momentum). This leads to a
spectral function of

$\mathcal{A}(\omega)=-\frac{1}{\pi}\Im\frac{\Delta}{\omega^{2}+i\omega\Gamma-\Delta^{2}}=\frac{1}{\pi}\frac{\omega\Gamma\Delta}{(\omega^{2}-\Delta^{2})^{2}+\Gamma^{2}\omega^{2}}.$
(7)

In the limit of $\Gamma\to 0$, this reduces to the spectrum found from the
Hamiltonian (1). A realistic estimate for $\Gamma$, based on Raman data [91]
is that $\Gamma\sim 15$ meV, though it seems there is a great amount of
uncertainty about this parameter 111It also seems like there is a good amount
of uncertainty about the resonance frequency as well, and whether features are
one-orbital, two-orbital, or of other origin all together.. From the spectral
function, we find the frequency domain Green’s functions from Kramers-Kronig
relations as

$\displaystyle
D^{R}(\omega)=\frac{\Delta}{\omega^{2}+i\Gamma\omega-\Delta^{2}}$ (8a)
$\displaystyle
D^{A}(\omega)=\frac{\Delta}{\omega^{2}-i\Gamma\omega-\Delta^{2}}$ (8b)
$\displaystyle D^{K}(\omega)=-2\pi
i\coth\left(\frac{\beta\omega}{2}\right)\mathcal{A}(\omega).$ (8c)

This then generates an effective action for the spin after integrating out the
bath in the Gaussian approximation (for the orbitals) of

$S_{\rm eff}=-\frac{\lambda^{2}}{2}\int dt\int
dt^{\prime}\sum_{j}\left(\mathbf{S}_{j,cl}(t^{\prime})\cdot\mathbf{n}_{j}(t^{\prime}),\mathbf{S}_{j,q}(t^{\prime})\cdot\mathbf{n}_{j}(t^{\prime})\right)\begin{pmatrix}0&D^{A}(t^{\prime},t)\\\
D^{R}(t^{\prime},t)&D^{K}(t^{\prime},t)\\\
\end{pmatrix}\begin{pmatrix}\mathbf{n}_{j}(t)\cdot\mathbf{S}_{j,cl}(t)\\\
\mathbf{n}_{j}(t)\cdot\mathbf{S}_{j,q}(t)\\\ \end{pmatrix}.$ (9)

As a final step, we simplify by averaging over the four titanium sublattices.
For a long-wavelength spin-wave, it is reasonable to expect the magnon to only
be sensitive to the average of the four titanium sites. It is worth pointing
out that this possibly fails for a short-wavelength magnon, which is localized
to the order of one unit cell. In this case it is possible that the distinct
nature of the orbital bath on each site may be important and could be an
important source of quantum fluctuations, though we leave this for future
work.

Proceeding on, if we average over the four sites, we generate an effective
local action describing the orbital bath of

$S_{\rm eff}=-\frac{\lambda^{2}}{2}\int dt\int
dt^{\prime}\sum_{j}\left[\mathbf{S}_{j,cl}(t^{\prime})\cdot\mathbb{\hat{D}}^{A}(t^{\prime},t)\cdot\mathbf{S}_{j,q}(t)+\mathbf{S}_{j,q}(t^{\prime})\cdot\mathbb{\hat{D}}^{R}(t^{\prime},t)\cdot\mathbf{S}_{j,cl}(t)+\mathbf{S}_{j,q}(t^{\prime})\cdot\mathbb{\hat{D}}^{K}(t^{\prime},t)\cdot\mathbf{S}_{j,q}(t)\right].$
(10)

This involves the sublattice averaged anisotropy tensor
$\mathbb{N}(t^{\prime},t)=\overline{\mathbf{n}(t^{\prime})\otimes\mathbf{n}(t)}$
through

$\check{\mathbb{{D}}}(t^{\prime},t)=\mathbb{N}(t^{\prime},t)\check{D}(t^{\prime},t).$
(11)

This tensor is presented in Fig. 2, which illustrates the matrix elements
along the $x,y,z$ axes (note the $x,y,z$ axes may not necessarily align with
the $a,b,c$ axes of the crystal but rather are defined by the orientation of
the crystal-field levels).

Figure 2: Matrix elements of sublattice averaged anisotropy tensor
$\mathbb{N}=\overline{\mathbf{n}\otimes\mathbf{n}}$ in equilibrium. We see
that the off-diagonal elements are zero, indicating the eigenbasis is aligned
with the crystalline axes. We also see the eigenvalues are nondegenerate, due
to the orthorhombicity of the crystal. The matrix is further partitioned into
the $ab$-plane ($x,y$) components and the $c$ axis ($z$) components. To a good
approximation, the matrix projection along the $c$-axis is zero, while the
anisotropy in the $ab$-plane is appreciable but not extreme.

##  IV Equilibrium Spin-Orbit Coupling

We now analyze the spin-orbit coupling in equilibrium. At low temperatures, we
can expand around the fully-polarized $|\uparrow\uparrow\uparrow...\rangle$
ground state via the Holstein-Primakoff expansion. We then describe magnons in
terms of canonical bosons $\hat{b}_{j}$ via the formal mapping

$\displaystyle\hat{S}^{z}_{j}=S-\hat{b}^{\dagger}_{j}\hat{b}_{j}$ (12a)
$\displaystyle\hat{S}^{+}_{j}=\sqrt{2S-\hat{b}^{\dagger}_{j}\hat{b}_{j}}\hat{b}_{j}$
(12b)
$\displaystyle\hat{S}^{-}_{j}=\hat{b}_{j}^{\dagger}\sqrt{2S-\hat{b}^{\dagger}_{j}\hat{b}_{j}}.$
(12c)

We then expand in the large-$S$ limit up to order $1/S$ to find the linear
spin-wave Hamiltonian. The Heisenberg interaction (along with the external
field along $\mathbf{e}_{z}$) gives the standard form, which is diagonalized
in momentum space to give

$\hat{H}^{(2)}_{0}=\sum_{\bf p}\Omega_{\bf p}\hat{b}_{\bf
p}^{\dagger}\hat{b}_{\bf p},$ (13)

with dispersion relation (for spin $S=\frac{1}{2}$)

$\Omega_{\bf p}=6SJ\left[1-\frac{1}{3}(\cos p_{x}+\cos p_{y}+\cos
p_{z})\right]+6SJ^{z}_{0}.$ (14)

This has a gap set by the easy-plane anisotropy energy $6SJ^{z}_{0}$, and has
a bandwidth of order $12SJ\sim 18$\mathrm{meV}$$. At this point, we still need
to include the Lamb shift due to orbital fluctuations. This dispersion
relation is depicted in Fig. 3(a), along with the corresponding single-
particle density-of-states (DOS) in Fig. 3(b), computed using Monte Carlo
sampling.

Figure 3: (a) Linear spin-wave dispersion relation along $\Gamma-X-R-\Gamma$
cut in Brillouin zone using idealized cubic model with $J=2.75$ meV and magnon
gap of $\Omega_{0}=0.1$ meV (after renormalizing away the Lamb shift from the
orbital bath). (b) Magnon density of states (DOS) $\rho(E)=\int_{\bf
q}\delta(E-\Omega_{\bf q})$ obtained by Monte Carlo sampling dispersion
relation in (a). This is used later when we evaluate the integrals of the
kinetic equation. (c) Equilibrium occupation of magnons with dispersion in
(a), in terms of the magnetic moment $M_{z}=2\mu_{B}(S-n_{\rm eq})$. We expect
the Holstein-Primakoff expansion to underestimate the role of fluctuations
since $1/S$ is not small in reality.

We now include the orbital self-energy which can be written in the frequency
domain due to the time-translational invariance in equilibrium. We expand to
order $S$, neglecting the linear term which should vanish when expanding
around the ground state. We find a Gaussian action for the magnons of

$S^{(2)}=\int_{p}\overline{b}_{p}\left[\check{G}_{0}^{-1}(p)-\check{\Sigma}(p)\right]b_{p}.$
(15)

We expand the orbital bath term to $O(S)$ in Appendix B in order to find the
magnon self-energy $\check{\Sigma}(p)$.

One can calculate the anisotropy due to the orbital fluctuations and find that
due to the orthorhombic nature of YTO, it has three distinct eigenvectors.
Calculating the projections, we find essentially no projection along the $c$
axis, with approximately 75% and 25% along the two in-plane directions. This
has the result of inducing an easy and hard axis for fluctuations, which leads
to quantum fluctuations manifested by the anomalous correlation functions of
$\langle b_{p}b_{-p}\rangle$. In our simplistic treatment, we neglect these,
finding an approximately $U(1)$ system, with dominant eigenvalue of $.5$ and
an anisotropy of $.25$ splitting the two principal axes.

For simplicity, we will neglect the anisotropy so we may obtain a diagonal
self-energy, leading to a Green’s function of

$G^{R}(\omega,\mathbf{p})=\left[\omega-\Omega_{\bf
p}-S\lambda^{2}N_{+-}D^{R}(\omega)\right]^{-1}.$ (16)

We have $N_{+-}\sim.48$, which is the isotropic projection of the in-plane
components of the anisotropy tensor. We plot the magnon spectral function in
the $(\omega,\mathbf{p})$ plane in Fig. 4, however it is worth briefly
examining the effects of the orbital bath perturbatively.

Figure 4: Magnon spectral function $A_{\rm
mag}(\omega,\mathbf{p})=-\frac{1}{\pi}\Im G^{R}(\omega,\mathbf{p})$ including
orbital bath for $N_{+-}=.5$, $\lambda=15$ meV, $\Delta=90$ meV, and
$\Gamma=10$ meV. Plotted along same dispersion contour as Fig. 3. Damping is
approximately proportional to frequency, so that the linewidth $\gamma_{\bf
p}\sim\Omega_{\bf p}$. The Lamb shift due to the orbital bath is renormalized
away so that the magnon gap is the physically measured gap of $\Omega_{0}=0.1$
meV.

Due to the large separation of scales between the orbital and spin degrees of
freedom, we can analyze the corrections to the magnon spectrum perturbatively.
We find a Lamb shift due to the coupling to the reservoir which shifts the
magnon band gap (it is essentially a source of single-ion anisotropy of the
easy-plane type). We then find renormalized magnon gap of

$\Omega_{0}=6SJ_{0}^{z}-\lambda^{2}SN_{+-}/\Delta.$ (17)

This is used to fix the counterterm $J_{0}^{z}$ by matching this to
experiment. It is empirically observed that the gap for magnons is quite small
[90], which is in and of itself an interesting fact though we won’t dwell on
this here. We also find a finite lifetime is generated for the magnons via
their interaction with the bath. This has a strong energy dependence and is
found to be

$\gamma(E_{\bf p})=\pi S\lambda^{2}N_{+-}\mathcal{A}(E_{\bf
p})=\frac{S\lambda^{2}N_{+-}\Gamma E_{\bf p}}{\Delta^{3}}.$ (18)

In particular, the imaginary part scales with $E_{\bf p}$, indicating it is
essentially a form of Ohmic Gilbert damping due to the orbital bath. Taking
estimates for YTO parameters of $\Delta\sim 90$ meV, $\Gamma\sim 15$ meV, and
$\lambda\sim 15$ meV, we find a lifetime in ns of

$\tau_{\bf p}=3.6\textrm{ ns}\frac{\textrm{meV}}{E_{\bf p}}.$ (19)

We note to the reader that the $\tau_{\mathbf{p}}$ (the lifetime for magnon
with momentum $\mathbf{p}$) is completely distinct from
$\hat{\tau}^{1},\hat{\tau}^{2}$ which correspond to the orbital operators, and
also from $\tau_{d}$ which corresponds to the lifetime of the phonon ring-
down. These all should occur in separate contexts, but emphasize this
distinction here to avoid confusion. At $T\sim$ 10 K, we have typical magnon
energies of $E_{\bf p}\sim 1$ meV and thus we have a typical lifetime for the
magnetization relaxation of order 3.6 ns according to this model, though other
channels for spin-flip processes may reduce this time according to
Matthiessen’s rule. We now proceed on to study the nonequilibrium dynamics of
this system.

##  V Nonequilibrium Dynamics

We now discuss the effect of the strong optical pulse. Focusing our attention
on to the most striking of the three pump frequencies from Ref. [19], which is
the pump at 9 THz, we start by describing how this pump effects the orbital
state.

Figure 5: (a) Time scales of pump-induced dynamics. Incident THz pulse (red)
resonantly excites a phonon mode which then exhibits coherent oscillations
(blue) for time scale $\tau_{\rm ph}\gg 2\pi/\Omega_{\rm ph}$. (b) These
dynamics also lead to an acceleration of the spin-orbit mediated magnetization
dynamics which leads to faster dynamics during the oscillations due to the
appearance of a new channel for spin-flip decay via the phonon-induced
sidebands. (c) After the oscillations decay, the dynamics returns to the
slower time scale present in equilibrium.

###  V.1 Pump Model

As per the estimates of the experiment [19], we consider a terahertz pulse
which resonantly drives an IR active phonon mode; in the experiment [19] these
were at frequencies of 4 THz, 9 THz, and 17 THz. Even though the pump itself
is quite short, the coherent oscillations it initiates in the phonon mode are
estimated to live much longer, with a ring-down time of order of 20-30 ps. We
therefore focus on the magnon dynamics which are induced by these coherent
ring-down dynamics rather than the initial pulse which is a quite short
duration. We use a ring-down model of the form

$Q_{IR}(t)=Q_{0}e^{-t/\tau_{d}}\sin(\Omega_{d}t)\theta(t),$ (20)

with initial (and maximal) excitation amplitude $Q_{0}$, central frequency
$\Omega_{d}$, and ring-down time $\tau_{d}$.

For our purposes, we will assume the pump has two main effects; first, it is
assumed to induce a transient change in the spin-exchange $J$ due to a
standard spin-phonon coupling mechanism. The origin of this mechanism is not
the main focus of this work, though it may also be interesting. We simply
model this as a coupling between $Q_{\rm IR}$ and
$\mathbf{S}_{j}\cdot\mathbf{S}_{j+\delta}$ of the form

$H_{\rm sp-
ph}=\sum_{j,\bm{\delta}}-\beta\mathbf{\hat{S}}_{j}\cdot\mathbf{\hat{S}}_{j+\bm{\delta}}Q_{\rm
IR}^{2}(t).$ (21)

This leads to a transient change in $J$ such that we have instantaneous value
of $J(t)=J+\beta Q_{\rm IR}^{2}(t)$. We focus on the rectified part of this,
which leads to a change in the exchange of

$\Delta J(t)=\frac{1}{2}\beta Q_{0}^{2}e^{-2t/\tau_{d}}=\Delta
J(0)e^{-2t/\tau_{d}}.$ (22)

We consider two cases —pump-induced enhancement of $\Delta J(0)=.5J$ and pump-
induced destruction $\Delta J(0)=-.5J$.

In addition to this, we also have argued that the pump induces substantial
changes to the excited crystal-field eigenvector, which in turn leads to a
dynamical modulation of the spin-orbit coupling between the magnons and
orbital bath. This is motivated by Fig. 6, which shows how the orbital angular
momentum associated to the first-excited crystal-field transition changes with
$Q_{\rm IR}$ in the case of YTO.

Figure 6: Change in in-plane projection of angular momentum unit vector
$\mathbf{n}_{j}$ for 9 THz $B_{2u}$ polarized phonon mode in YTO. We expect
for realistic fluence, the peak amplitude is of order $1$ in these units,
leading to an appreciable change in the eigenvalues of the anisotropy tensor
$\mathbb{N}$ which oscillate at frequencies $\pm 2\Omega_{d}$. This motivates
the amplitude parameter $A_{d}$ of order $A_{d}\sim\sqrt{.1}=.3$.

We see there is a quadratic coupling between the phonon mode $Q_{IR}^{2}$ and
the crystal field eigenvector $\mathbf{n}_{j}$ such that we write

$\mathbf{n}_{j}\sim\mathbf{n}_{j}+\delta\mathbf{n}Q^{2}_{\rm IR}(t).$ (23)

This also has a rectified part, which may lead to an interesting pump-induced
renormalization of the magnetic anisotropy [29, 30], however here we will
focus instead on the dynamic harmonics, which can dramatically change the
nature of magnetic relaxation in this system. This is modeled as a change in
the anisotropy tensor, which is obtained by averaging this over the four Ti
sublattices (we refer to Appendices B and D).

We write this as

$\mathbb{N}(t,t^{\prime})=\mathbb{N}_{\rm eq}\\\
\times\left[A_{0}(t)+A_{1}(t)e^{-i\Omega_{d}t}+A_{-1}(t)e^{i\Omega_{d}t}\right]\\\
\times\left[A_{0}(t^{\prime})+A_{1}^{*}(t^{\prime})e^{+i\Omega_{d}t^{\prime}}+A_{-1}^{*}(t^{\prime})e^{-i\Omega_{d}t^{\prime}}\right]$
(24)

where roughly, $A_{0}^{2}+|A_{1}|^{2}+|A_{-1}|^{2}=1$ models the rotation of
the excited state unit vector without changing the net length, such that the
projection onto the $c$-axis changes by $|A_{-1}|^{2}+|A_{-1}|^{2}$. We assume
the phase is not important, and take

$A_{1}(t)=A_{-1}(t)=\frac{1}{2}A_{d}e^{-2t/\tau_{d}},$ (25)

in line with the same pump-profile as the one which drives the change in
exchange.

In order to determine the effect of this Floquet-driven coupling to the
orbital bath we will work in the approximation that we may still separate the
time scales associated to (i) the magnon dynamics, (ii) the transience of the
drive, (iii) the orbital dynamics. It is very interesting, however
challenging, to relax this hierarchy and allow for a complete breakdown of
separation of time scales. This is left open for future works to handle.
Additionally, though this model captures the essential physics, it is still
only qualitatively motivated by YTO calculations, and in future work a more
detailed calculation of how exactly the changes in crystal-field evolve for
each phonon mode would be warranted. In particular, it may be the case that
different phonon modes are more or less effective at modulating various
components of this tensor and may allow for a more selective control over the
effects we describe here.

###  V.2 Quasiparticle Dynamics

We now examine the magnon dynamics in the presence of a hypothetical Floquet
modulation of the spin-orbit interaction. As argued in the previous section,
this is a reasonable model of the pumped phonon’s effect on the spin-orbit
coupling. To simplify matters, we assume that the pump doesn’t actually change
the orbital correlations or fluctuations, but rather changes the coupling of
the magnons to the orbital bath.

By using the Keldysh technique we are able to calculate the real-time
dynamical evolution of the magnon correlation functions, as detailed in
Appendix B. The key object of interest in this work is the magnon occupation
function, which is encoded in the Keldysh correlation function $G^{K}_{\bf
p}(t,t^{\prime})$, here taken to be diagonal in momentum space. From this, we
can then obtain the net magnetization as a function of time.

We further utilize the separation of time-scales between the evolution under
the pump profile and the internal frequency scales by taking the Wigner-
transform of $G^{K}$, which encodes the full two-time dependence in terms of a
“center-of-mass” time, which corresponds to the slow evolution, and the
frequency, which encodes the rapid oscillations in the relative time
difference. The Wigner transformed Keldysh function is

$G^{K}_{\bf p}(T;\omega)=\int d\tau G^{K}_{\bf
p}(T+\frac{\tau}{2},T-\frac{\tau}{2})e^{i\omega\tau}.$ (26)

From this, we can extract the total magnon density as a function of time as

$n(t)=\frac{1}{2}\int_{\bf p}\left(\int\frac{d\omega}{2\pi}iG^{K}_{\bf
p}(t;\omega)-1\right),$ (27)

and the corresponding magnetization is then found to be

$M^{z}(t)=2\mu_{B}(S-n(t)).$ (28)

By systematically expanding in terms of gradients of the slowly-varying pump
profile, we derive in Appendix D an effective relaxation-time approximation
for this, which to the very lowest order reads

$\frac{\partial G^{K}_{\bf p}(T;\omega)}{\partial T}\\\
=2i\Im\Sigma^{R}(T;\omega)\left[-iG^{K}_{\bf p}(T;\omega)+2\pi\mathcal{A}_{\rm
mag}(T;\omega,\mathbf{p})F_{\rm orb}(\omega)\right].$ (29)

Here $\mathcal{A}_{\rm mag}(T;\omega,\mathbf{p})=-1/\pi\Im G^{R}_{\bf
p}(T;\omega)$ is the instantaneous magnon spectral function, which depends on
time in the instance where the pump changes, e.g. the spin-exchange, as it
does in this system. We also see the appearance of the orbital occupation
function, which we assume remains in equilibrium at temperature $T_{\rm orb}$,
such that $F_{\rm orb}(\omega)=\coth\frac{\omega}{2T_{\rm orb}}$.

We now study the dynamics of this system under the quasiparticle
approximation, such that we can replace the frequency dependence by the
instantaneous on-shell frequency. This gives us a simple equation we can solve
for the quasiparticle occupation function $f_{\bf p}(T)$ of

$\frac{\partial f_{\bf p}(T)}{\partial T}=-\frac{1}{\tau_{\bf
p}(T)}\left(f_{\bf p}(T)-f^{({\rm bath})}_{\bf p}(T)\right),$ (30)

where $1/\tau_{\bf p}(T)$ is the instantaneous relaxation rate at time $T$,
derived from the magnon self-energy, and $f^{({\rm bath})}_{\bf p}(T)$ is the
instaneous equilibrium occupation set by the orbital bath occupation function
projected onto the magnon spectral density. For the details, we refer to
Appendix D.

If we only include the change in $J$, and therefore only include the
instantaneous change in the spectral function, we see a meager response to the
pump. This is shown in Fig. 7, which shows the change in instantaneous
magnetization following a transient increase in $J$ due to the coherent phonon
rind-down, schematically illustrated above the numerical plot. We plot the
change in magnetization $\Delta M_{z}(t)$ as a percent relative to the maximum
possible change, which would be $2\mu_{B}(S-n(0))$ so that if the initial
moment is $.9\mu_{B}$ and it increases to $.95\mu_{B}$ this would by 50% of
the maximum possible increase.

Figure 7:  Change in magnetization following a pump-induced change in the
ferromagnetic exchange due to the rectified spin-phonon coupling, modeled here
as a transient $\Delta J(t)\sim Q_{\rm IR}^{2}(t)$ which follows the impulsive
initial pulse, illustrated atop the frame. For a ring-down time of order
$\tau_{\rm ph}\sim 30$ ps and an initial change in the exchange of $\Delta
J(0)/J_{\rm eq}=50\%$ we find magnetization dynamics in the plot below, which
showns the change in magnetization $\Delta M_{z}(t)$ in terms a percent of the
maximum possible enhancement $\Delta M_{\rm max}$ (corresponding to a complete
saturation of the magnetization). Ring-down period is shaded red. This is not
including the resonant enhancement of the magnetization dynamics.

Though the magnetization does generally follow the pump-induced change $\Delta
J$, which here was set to $.5J(0)$, it is a relatively mediocre response since
the dynamics are still quite bottlenecked by the long-relaxation time,
$\tau_{\bf p}$ which is of order nanoseconds for a thermal magnon, whereas the
duration of the pump-induced oscillations are at most 50 ps.

However, as we argued before, the nonequilibrium dynamics induced by the pump
can potentially have exhibit accelerated timescales, as illustrated in Fig. 5.
Due to a combination of high-frequency oscillations at $2\Omega_{d}\sim 80$
meV and low-lying orbital excitations with $\Delta\sim 90$ meV or so, we can
find a transient acceleration of the relaxation rate, quantified by
$1/\tau_{\bf p}=-2/\pi\Im\Sigma^{R}(T;\omega)$, making the system essentially
relax faster than in equilibrium during the driving period. This is confirmed
by calculating the effective magnon lifetime in the presence of steady-state
coherent oscillations. In Fig. 8 we plot the magnon lifetime $\tau_{\bf p}$ as
a function of the magnon kinetic energy $E_{\bf p}$ for different pump
frequencies $\omega_{d}$ and amplitudes $A_{d}$ 222In fact, $A_{d}\sim Q_{\rm
IR}^{2}$ is itself scaling linearly with fluence, so that scaling
quadratically with $A_{d}^{2}$ implies quadratic fluence dependence..

Figure 8: (a) Plot of magnon lifetime $\tau_{\bf p}=1/\gamma_{\bf p}$ in the
presence of the coherent phonon enhancement as a function of magnon kinetic
energy $E_{\bf p}$ for different drive frequencies. We fix fluence parameter
$A_{d}=.3$ and fix orbital parameters to $\Delta=90$ meV and $\Gamma=10$ meV,
with $\lambda=15$ meV. The lifetime is reduced by the magnon appearance of
sidebands at $\Delta\pm 2\Omega_{d}$. Near $\Omega_{d}=40$ meV this process
nears resonance and the decay rate is maximally enhanced by nearly two orderes
of magnitude. (b) We study for varying drive fluence parameter $A_{d}$ at
fixed $\Omega_{d}=35$ meV for the same orbital parameters. The dependence in
this model is monotonic, though in a more refinded model we would expect some
saturation as $A_{d}\to 1$.

To see whether the increased relaxation rate has any effect in practice, we
carry out the simulations of Eq. (30) now including both the pump-induced
change in $J(t)$ as well as the pump-induced change in relaxation rate. This
is presented in Fig. 9 which shows the equivalent $\Delta J$ as in Fig. 7 but
now including the pump-accelerated relaxation rate for different frequencies
$\Omega_{d}$ at fixed fluence $A_{d}=.3$. We see that when the pump approaches
resonance with the orbital excitation, the dynamics greatly accelerates and as
a result, the magnetization can grow much more over the same $\sim 30$ ps
window of growth time.

Figure 9: Fractional change in magnon occupation following pump-induced change
in exchange $\Delta J$ while also including the enhancement of the relaxation
rate due to the phonon ring-down. For different pump frequencies (here we only
model the pump frequency as changing the spin-flip time) we see a dramatic
increase in the maximum change in magnetization upon approaching the resonance
condition around $\Omega_{d}\sim 40$ meV. For pump frequencies above this, the
effect quickly reverses and by $\Omega_{d}\sim$ 45 meV we see the dynamics has
actually slowed substantially.

Curiously, we see that around $\Omega_{d}=45$ meV, the effect seems to
completely dissappear, and the resulting magnetization growth is almost
completely stunted. In fact, this is a manifestation of the pump actually
passing through the orbital resonance and changing from red-detuning to blue-
detuning. If we continue to increase the drive frequency further, we find that
the relaxation rate actually becomes negative—an effect which is clearly
impossible in equilibrium. This negative relaxation rate essentially indicates
that in the rotating frame the orbital bath is population-inverted with
respect to the magnon system. Therefore, the bath actually acts as a gain
medium rather than a retarder. The resulting dynamics are shown in Fig. 10
where we simulate both an initial increase in exchange, as in Fig. 9, as well
a pump-induced reduction in $J$ of $\Delta J/J(0)=-50\%$. We see that the
response is most pronounced when $\Omega_{d}$ is around $\pm 5$ meV detuned
from the $\Delta/2=45$ meV point. We also see that the negative relaxation
rate essentially leads to an effectively reversed sign of $\Delta J$, leading
to growth in magnon number when it should become less ferromagnetic, and vice
versa.

Figure 10: Magnetization dynamics for frequencies below and above resonance.
(a) For a transient increase in $J$ of 50% ferromagnetism should increase in
equilibrium, however in the presence of a high-frequency drive this can
amplify, diminish, or even reverse as the frequency passes through resonance
with the bath. For $\Omega_{d}=40$ meV the relaxation rate reaches near
maximal enhancement and the magnetization is most responsive to the pump-
induced increase in $J$, while for $\Omega_{d}=50$ meV it has already passed
to the other side of the resonance. The bath now acts to induce “gain” rather
than loss and drives the magnetization opposite to the naive result, quite
dramatically. (b) If we consider instead a pump-induced reduction in $J$ of
=50% the same features qualitatively persist, with opposite directions. In
this case, driving above the resonance leads to a substantial enhancement of
magnetization.

Therefore, we see that not only can one try to accelerate magnetic dynamics
away from equilibrium by modulating the coupling to the orbital bath, but one
may even potentially slow the dynamics down (in our example, by tuning
$\Omega_{d}\sim 45$ meV) or reverse them altogether by changing from red- to
blue-detuning. This is a genuinely nonequilibrium process and may potentially
explain the apparent opposite trend between the equilibrium spin-phonon
coupling and pump-induced response in YTO in the recent experiment [19].

We can more systematically map this effect out by plotting the most extreme
value of the time-traces as a function of $\Omega_{d}$ and temperature, shown
in Fig. 11(a) as a density plot, and in Fig. 11(b) for two line-cuts at fixed
temperature $T$. We see quite clearly that the dynamics are most dramatically
affected near the resonance of $2\Omega_{d}=\Delta$, and upon passing through
the resonance the sign of the effect changes.

Figure 11: (a) Plot of maximum change in magnetization as a function of
initial temperature and pump frequency. (a) Color map in $(\Omega_{d},T)$
plane. We see that crossing through the resonance at $\Omega_{d}=\Delta/2$
there is a dramatic change in the sign of the effect, and that the greatest
change occurs in this region. (b) Line cuts at low temperature ($T=1$ K) and
high temperature ($T=10$ K). We see that the effect is slightly more efficient
at increasing the magnetization when temperature is low, while it is more
effective at reducing the magnetization at higher temperatures.

Thus, we see that going beyond the “quasistatic” picture and actually
considering how the coupling to the orbital bath changes in the presence of
nonequilibrium dynamics can lead to striking, and potentially useful changes
to the magnetization dynamics. Crucially, the effect we outline here relies on
a relatively low-lying orbital excitation which couples to the spins and also
the pumped phonon modes. If the orbitals are too low-lying then they will
exhibit strong fluctuations and cannot be treated as a bath, as we have here.
On the other hand, if they are too high in excitation energy, they cannot be
effectively coupled to by phonon oscillations and therefore are cannot
realistically participate in the dynamics. Thus, quasi-degenerate magnetic
insulators present a special opportunity for this type of “bath-control,”
although as we will discuss next, this type of physics may be able to be
extended to more general systems such as antiferromagnets, superconductors, or
potentially other correlated phases.

##  VI Discussion

We now summarize our findings. We considered a simple model for the
nonequilibrium dynamics of magnons in a Heisenberg ferromagnetic insulator
with low-lying “quasi-degenerate” orbitals, as may be be realized in
orthorhombic titanates $R$TiO3 ($R=$ Y, Sm, Gd), and possibly other compounds.
By using a powerful terahertz pulse to resonantly excited optical phonons we
argued that relatively long-lasting nonequilibrium dynamics can be induced by
the coherently oscillating phonon modes, which may have lifetimes lasting up
to 30 ps. These phonon oscillations may lead to transient modifications to the
superexchange though, e.g. the rectified part of the spin-phonon coupling
$\sim Q_{\rm IR}^{2}\mathbf{S}_{j}\cdot\mathbf{S}_{k}$; this may then lead to
dynamic changes in the magnetic free-energy landscape which can potentially be
used to optical drive the magnetization and control the phase diagram.
However, this dynamics is often plagued by a bottleneck due to small spin-
orbit coupling which leads to an approximate conservation law for
magnetization, leading to slow diffusive dynamics on the relevant time scales.

This bottleneck can be circumvented in a nonequilibrium setting, as we showed
in Sec. V. In particular, in the presence of low-lying orbital excitations,
the coupling between magnons and the orbital angular momentum can become
unquenched in the presence of phonon dynamics, which may lead to “stimulated
emission” type processes into the orbital bath. This can in principle lead to
a significant acceleration in time scale for the magnetization dynamics,
allowing for more effective optical control on relevant time scales.
Furthermore, we found that in principle it is even possible to reverse the
nature of the coupling to the bath by changing from red- to blue-detuning with
respect to the orbital bath, allowing for an even greater degree of control
over the magnetization dynamics.

More generally, our results should be able to be applied to other systems of
interest including antiferromagnetic insulators, spin liquids, and other
correlated insulators. The key component is the ability to induce a dynamical
coupling between the degrees of freedom of interest (such as spins) and the
bath degrees of freedom. In addition to controlling the bath decay rates this
may also allow to control the bath-induced Lamb shift, which in the case we
consider here enters as an effective single-ion anisotropy. Thus, it may also
be possible to control the anisotropy dynamically, as proposed in the recent
experiment [29] and theory [30]. Control over the isotropic superexchange
interaction may also be possible through the mechanism we outline here as well
as through similar mechanisms [49, 20, 57, 24, 31].

Our results may also be relevant to recent experiments on nonequilibrium
light-induced superconductivity [97] in fullerides [9, 12, 11], organics salts
[98, 13], and cuprates [99, 14, 15]. In this case, we argue that there are a
number of parallels which make it even more interesting to understand this
physics. Chief among these are the observations of pump-induced signatures of
the ordered phase above the equilibrium transition temperature, long-lived
resilience of this long-range order, and enhancement of “coherence” below the
ordering temperature. In the case of superconductors, the effects of pump-
induced order are seen most clearly in systems which are strongly coupled and
don’t exhibit a simple mean-field transition [100, 101] (e.g. cuprates,
fullerides), and this is also the case for the magnetic order in the recent
experiment on YTO in Ref. [19], which appears to exhibit a “magnetic
pseudogap.” It is also possible that a similar equilibrium slowing-down of
reaction pathways occurs in these systems, which reside near to a
metal=insulator transition [102, 103].

Although in the current work we don’t address the “pseudogap regime,” it
should be possible to extend our results to include strong magnetic
fluctuations via, e.g. the Schwinger boson technique or various slave-particle
mappings, which can be extended to nonequilibrium settings [104, 61, 105, 106,
107]. It may turn out that the equivalent problem in the superconducting case
will actually be more tractable since in this case the theory for a
fluctuating superconductor is more amenable to nonequilibrium diagrammatic
approaches [108, 101, 109].

We also comment that similar ideas have recently been discussed in the context
of “pump-induced sideband cooling” for various solid-state systems by various
groups [51, 52, 53]. In particular, it was proposed that recent experiments on
light-induced superconductivity [9] could be understood by a Floquet sideband
cooling utilizing an intermediate bath state provided by an internal excitonic
resonance of the fullerene molecules [51]. This was later extended to the case
of a quantum spin-system with a dynamical coupling induced to a complementary
bath system [52]. In this respect, this is very similar to the system we are
proposing here, where the orbitals serve as an analogue to the excitonic bath
of Ref. [51]. However, our results should still be present even if the true
sideband cooling does not materialize. In particular, it is likely that both
processes will be happening in a true driven system.

To conclude, we have examined the nonequilibrium spin-orbital dynamics in a
ferromagnetic insulator and found that away from equilibrium there is a rich
variety of dynamical processes which can happen even in a relatively simple
quasiparticle description. Experimentally, this is possibly relevant to
various ferromagnetic insulators realized in ferromagnetic rare-earth
titanates $R$TiO3, and may be more generally applicable to strongly-correlated
spin-orbital systems such as NiPS3 [29, 30], CuSb2O6 [89], other titanates
[25, 78, 64, 24], manganites [110, 111], vandates [64, 26, 112, 113], and a
number of other compounds [114]. There may also be connections to charge-
density wave physics, which can also be manipulated by light [16, 115].

We also argued that our results may be analogous to recent experiments on
photoinduced superconductivity. In future works it will be important to
consider extending our results to the strongly fluctuating regime near, and
above $T_{C}$ as well as to incorporate the truly dynamical terms which break
time-translational symmetry [116, 43, 37]. In addition, considering systems
which do have degenerate orbitals (or exhibit genunie spontaneous orbital
ordering) would be of great interest, with many exotic phenomena already known
to occur [56]. It will also be necessary to develop closer connection to
specific materials in order to make contact with current and future
experiments. Experiments using ultrafast x-ray scattering may be able to
directly confirm these nonequilibrium dynamics [55, 117], though this is
likely to be quite challenging theoretically.

###### Acknowledgements.

The authors. would like to acknowledge crucial discussions with Pavel
Dolgirev, Eugene Demler, Andrey Grankin, Andy Millis, David Hsieh, Mohammad
Maghrebi, Benedetta Flebus, Aaron Müller, and Zhiyuan Sun. This work is
primarily supported by the Quantum Science Center (QSC), a National Quantum
Information Science Research Center of the U.S. Department of Energy (DOE).
P.N. acknowledges support as a Moore Inventor Fellow through Grant No.
GBMF8048 and gratefully acknowledges support from the Gordon and Betty Moore
Foundation as well as support from a Max Planck Sabbatical Award that enabled
this collaborative project.

## References

  * Nova  _et al._ [2019] T. F. Nova, A. S. Disa, M. Fechner, and A. Cavalleri, Metastable ferroelectricity in optically strained SrTiO3, Science 364, 1075 (2019), 1812.10560 .
  * Li  _et al._ [2019] X. Li, T. Qiu, J. Zhang, E. Baldini, J. Lu, A. M. Rappe, and K. A. Nelson, Terahertz field–induced ferroelectricity in quantum paraelectric SrTiO3, Science 364, 1079 (2019), 1812.10785 .
  * Wall  _et al._ [2018] S. Wall, S. Yang, L. Vidas, M. Chollet, J. M. Glownia, M. Kozina, T. Katayama, T. Henighan, M. Jiang, T. A. Miller, D. A. Reis, L. A. Boatner, O. Delaire, and M. Trigo, Ultrafast disordering of vanadium dimers in photoexcited VO2, Science 362, 572 (2018).
  * Sood  _et al._ [2021] A. Sood, X. Shen, Y. Shi, S. Kumar, S. J. Park, M. Zajac, Y. Sun, L.-Q. Chen, S. Ramanathan, X. Wang, W. C. Chueh, and A. M. Lindenberg, Universal phase dynamics in VO2 switches revealed by ultrafast operando diffraction, Science 373, 352 (2021), 2102.06013 .
  * Cavalleri  _et al._ [2001] A. Cavalleri, C. Tóth, C. W. Siders, J. A. Squier, F. Ráksi, P. Forget, and J. C. Kieffer, Femtosecond Structural Dynamics in VO2 during an Ultrafast Solid-Solid Phase Transition, Physical Review Letters 87, 237401 (2001).
  * Juraschek  _et al._ [2020] D. M. Juraschek, Q. N. Meier, and P. Narang, Parametric Excitation of an Optically Silent Goldstone-Like Phonon Mode, Physical Review Letters 124, 117401 (2020), 1912.07349 .
  * McLeod  _et al._ [2020] A. S. McLeod, J. Zhang, M. Q. Gu, F. Jin, G. Zhang, K. W. Post, X. G. Zhao, A. J. Millis, W. B. Wu, J. M. Rondinelli, R. D. Averitt, and D. N. Basov, Multi-messenger nanoprobes of hidden magnetism in a strained manganite, Nature Materials 19, 397 (2020), 1910.10361 .
  * Tobey  _et al._ [2008] R. I. Tobey, D. Prabhakaran, A. T. Boothroyd, and A. Cavalleri, Ultrafast Electronic Phase Transition in La1/2Sr3/2MnO4 by Coherent Vibrational Excitation: Evidence for Nonthermal Melting of Orbital Order, Physical Review Letters 101, 197404 (2008).
  * Mitrano  _et al._ [2016] M. Mitrano, A. Cantaluppi, D. Nicoletti, S. Kaiser, A. Perucchi, S. Lupi, P. D. Pietro, D. Pontiroli, M. Riccò, A. Subedi, S. R. Clark, D. Jaksch, and A. Cavalleri, Possible light-induced superconductivity in K3C60 at high temperature, Nature 530, 461 (2016), 1505.04529 .
  * Rajasekaran  _et al._ [2016] S. Rajasekaran, E. Casandruc, Y. Laplace, D. Nicoletti, G. D. Gu, S. R. Clark, D. Jaksch, and A. Cavalleri, Parametric amplification of a superconducting plasma wave, Nature Physics 12, 1012 (2016), 1511.08378 .
  * Buzzi  _et al._ [2021a] M. Buzzi, G. Jotzu, A. Cavalleri, J. I. Cirac, E. A. Demler, B. I. Halperin, M. D. Lukin, T. Shi, Y. Wang, and D. Podolsky, Higgs-Mediated Optical Amplification in a Nonequilibrium Superconductor, Physical Review X 11, 011055 (2021a), 1908.10879 .
  * Budden  _et al._ [2021] M. Budden, T. Gebert, M. Buzzi, G. Jotzu, E. Wang, T. Matsuyama, G. Meier, Y. Laplace, D. Pontiroli, M. Riccò, F. Schlawin, D. Jaksch, and A. Cavalleri, Evidence for metastable photo-induced superconductivity in K3C60, Nature Physics 17, 611 (2021), 2002.12835 .
  * Buzzi  _et al._ [2021b] M. Buzzi, D. Nicoletti, S. Fava, G. Jotzu, K. Miyagawa, K. Kanoda, A. Henderson, T. Siegrist, J. A. Schlueter, M.-S. Nam, A. Ardavan, and A. Cavalleri, Phase Diagram for Light-Induced Superconductivity in $\kappa$-(ET)2-X, Physical Review Letters 127, 197002 (2021b), 2106.14244 .
  * Cremin  _et al._ [2019] K. A. Cremin, J. Zhang, C. C. Homes, G. D. Gu, Z. Sun, M. M. Fogler, A. J. Millis, D. N. Basov, and R. D. Averitt, Photoenhanced metastable c-axis electrodynamics in stripe-ordered cuprate La1.885Ba0.115CuO4, Proceedings of the National Academy of Sciences 116, 19875 (2019), 1901.10037 .
  * Liu  _et al._ [2020] B. Liu, M. Först, M. Fechner, D. Nicoletti, J. Porras, T. Loew, B. Keimer, and A. Cavalleri, Pump Frequency Resonances for Light-Induced Incipient Superconductivity in YBa2Cu3O6.5, Physical Review X 10, 011053 (2020), 1905.08356 .
  * Kogar  _et al._ [2020] A. Kogar, A. Zong, P. E. Dolgirev, X. Shen, J. Straquadine, Y.-Q. Bie, X. Wang, T. Rohwer, I.-C. Tung, Y. Yang, R. Li, J. Yang, S. Weathersby, S. Park, M. E. Kozina, E. J. Sie, H. Wen, P. Jarillo-Herrero, I. R. Fisher, X. Wang, and N. Gedik, Light-induced charge density wave in LaTe3, Nat. Phys. 16, 159 (2020), 1904.07472 .
  * Lee  _et al._ [2019] M.-C. Lee, C. H. Kim, I. Kwak, C. W. Seo, C. Sohn, F. Nakamura, C. Sow, Y. Maeno, E.-A. Kim, T. W. Noh, and K. W. Kim, Strong spin-phonon coupling unveiled by coherent phonon oscillations in Ca2RuO4, Physical Review B 99, 144306 (2019), 1712.03028 .
  * Rohwer  _et al._ [2011] T. Rohwer, S. Hellmann, M. Wiesenmayer, C. Sohrt, A. Stange, B. Slomski, A. Carr, Y. Liu, L. M. Avila, M. Kalläne, S. Mathias, L. Kipp, K. Rossnagel, and M. Bauer, Collapse of long-range charge order tracked by time-resolved photoemission at high momenta, Nature 471, 490 (2011).
  * Disa  _et al._ [2021] A. S. Disa, J. Curtis, M. Fechner, A. Liu, A. v. Hoegen, M. Först, T. F. Nova, P. Narang, A. Maljuk, A. V. Boris, B. Keimer, and A. Cavalleri, Optical Stabilization of Fluctuating High Temperature Ferromagnetism in YTiO3, arXiv (2021), 2111.13622 .
  * Afanasiev  _et al._ [2021a] D. Afanasiev, J. R. Hortensius, B. A. Ivanov, A. Sasani, E. Bousquet, Y. M. Blanter, R. V. Mikhaylovskiy, A. V. Kimel, and A. D. Caviglia, Ultrafast control of magnetic interactions via light-driven phonons, Nature Materials 20, 607 (2021a), 1912.01938 .
  * Kimel  _et al._ [2009] A. V. Kimel, B. A. Ivanov, R. V. Pisarev, P. A. Usachev, A. Kirilyuk, and T. Rasing, Inertia-driven spin switching in antiferromagnets, Nature Physics 5, 727 (2009).
  * Liu  _et al._ [2019] B. Liu, W. Niu, Y. Chen, X. Ruan, Z. Tang, X. Wang, W. Liu, L. He, Y. Li, J. Wu, S. Tang, J. Du, R. Zhang, and Y. Xu, Ultrafast Orbital‐Oriented Control of Magnetization in Half‐Metallic La0.7Sr0.3MnO3 Films, Advanced Materials 31, 1806443 (2019).
  * Seifert and Balents [2019] U. F. P. Seifert and L. Balents, Optical excitation of magnons in an easy-plane antiferromagnet: Application to Sr2IrO4, Physical Review B 100, 125161 (2019), 1905.01313 .
  * Gu and Rondinelli [2018] M. Gu and J. M. Rondinelli, Nonlinear phononic control and emergent magnetism in Mott insulating titanates, Physical Review B 98, 024102 (2018), 1710.00993 .
  * Lovinger  _et al._ [2020a] D. J. Lovinger, E. Zoghlin, P. Kissin, G. Ahn, K. Ahadi, P. Kim, M. Poore, S. Stemmer, S. J. Moon, S. D. Wilson, and R. D. Averitt, Magnetoelastic coupling to coherent acoustic phonon modes in the ferrimagnetic insulator GdTiO3, Physical Review B 102, 085138 (2020a), 2009.10222 .
  * Lovinger  _et al._ [2020b] D. J. Lovinger, M. Brahlek, P. Kissin, D. M. Kennes, A. J. Millis, R. Engel-Herbert, and R. D. Averitt, Influence of spin and orbital fluctuations on Mott-Hubbard exciton dynamics in LaVO3 thin films, Physical Review B 102, 115143 (2020b), 2009.10219 .
  * Ron  _et al._ [2020] A. Ron, S. Chaudhary, G. Zhang, H. Ning, E. Zoghlin, S. D. Wilson, R. D. Averitt, G. Refael, and D. Hsieh, Ultrafast Enhancement of Ferromagnetic Spin Exchange Induced by Ligand-to-Metal Charge Transfer, Physical Review Letters 125, 197203 (2020), 1910.06376 .
  * Torre  _et al._ [2022] A. d. l. Torre, K. L. Seyler, M. Buchhold, Y. Baum, G. Zhang, N. J. Laurita, J. W. Harter, L. Zhao, I. Phinney, X. Chen, S. D. Wilson, G. Cao, R. D. Averitt, G. Refael, and D. Hsieh, Decoupling of static and dynamic criticality in a driven Mott insulator, Communications Physics 5, 35 (2022), 2112.08397 .
  * Afanasiev  _et al._ [2021b] D. Afanasiev, J. R. Hortensius, M. Matthiesen, S. Mañas-Valero, M. Šiškins, M. Lee, E. Lesne, H. S. J. v. d. Zant, P. G. Steeneken, B. A. Ivanov, E. Coronado, and A. D. Caviglia, Controlling the anisotropy of a van der Waals antiferromagnet with light, Science Advances 7, eabf3096 (2021b), 2010.05062 .
  * Seifert  _et al._ [2022] U. F. P. Seifert, M. Ye, and L. Balents, Ultrafast optical excitation of magnetic dynamics in van der Waals magnets: Coherent magnons and BKT dynamics in NiPS3, Physical Review B 105, 155138 (2022), 2111.03674 .
  * Maehrlein  _et al._ [2018] S. F. Maehrlein, I. Radu, P. Maldonado, A. Paarmann, M. Gensch, A. M. Kalashnikova, R. V. Pisarev, M. Wolf, P. M. Oppeneer, J. Barker, and T. Kampfrath, Dissecting spin-phonon equilibration in ferrimagnetic insulators by ultrafast lattice excitation, Science Advances 4, eaar5164 (2018), 1710.02700 .
  * Kimel  _et al._ [2002] A. V. Kimel, R. V. Pisarev, J. Hohlfeld, and T. Rasing, Ultrafast Quenching of the Antiferromagnetic Order in FeBO3: Direct Optical Probing of the Phonon-Magnon Coupling, Physical Review Letters 89, 287401 (2002).
  * Titov  _et al._ [2021] S. V. Titov, W. T. Coffey, Y. P. Kalmykov, M. Zarifakis, and A. S. Titov, Inertial magnetization dynamics of ferromagnetic nanoparticles including thermal agitation, Physical Review B 103, 144433 (2021).
  * Disa  _et al._ [2020] A. S. Disa, M. Fechner, T. F. Nova, B. Liu, M. Först, D. Prabhakaran, P. G. Radaelli, and A. Cavalleri, Polarizing an antiferromagnet by optical engineering of the crystal field, Nature Physics 16, 937 (2020), 2001.00540 .
  * Michael  _et al._ [2020] M. H. Michael, A. v. Hoegen, M. Fechner, M. Först, A. Cavalleri, and E. Demler, Parametric resonance of Josephson plasma waves: A theory for optically amplified interlayer superconductivity in YBa2Cu3O6+x, Physical Review B 102, 174505 (2020), 2004.13049 .
  * Knap  _et al._ [2016] M. Knap, M. Babadi, G. Refael, I. Martin, and E. Demler, Dynamical Cooper pairing in nonequilibrium electron-phonon systems, Physical Review B 94, 214504 (2016), 1511.07874 .
  * Babadi  _et al._ [2017] M. Babadi, M. Knap, I. Martin, G. Refael, and E. Demler, Theory of parametrically amplified electron-phonon superconductivity, Physical Review B 96, 014512 (2017), 1702.02531 .
  * Dolgirev  _et al._ [2021] P. E. Dolgirev, A. Zong, M. H. Michael, J. B. Curtis, D. Podolsky, A. Cavalleri, and E. Demler, Periodic dynamics in superconductors induced by an impulsive optical quench, arXiv (2021), 2104.07181 .
  * Bukov  _et al._ [2015] M. Bukov, S. Gopalakrishnan, M. Knap, and E. Demler, Prethermal floquet steady states and instabilities in the periodically driven, weakly interacting bose-hubbard model, Phys. Rev. Lett. 115, 10.1103/physrevlett.115.205301 (2015), 1507.01946 .
  * Oka and Aoki [2009] T. Oka and H. Aoki, Photovoltaic Hall effect in graphene, Phys. Rev. B 79, 081406 (2009), 0807.4767 .
  * Kitagawa  _et al._ [2010] T. Kitagawa, E. Berg, M. Rudner, and E. Demler, Topological characterization of periodically driven quantum systems, Phys. Rev. B 82, 235114 (2010), 1010.6126 .
  * Lindner  _et al._ [2011] N. H. Lindner, G. Refael, and V. Galitski, Floquet topological insulator in semiconductor quantum wells, Nat. Phys. 7, 490 (2011), 1008.1792 .
  * Seetharam  _et al._ [2015] K. I. Seetharam, C.-E. Bardyn, N. H. Lindner, M. S. Rudner, and G. Refael, Controlled Population of Floquet-Bloch States via Coupling to Bose and Fermi Baths, Physical Review X 5, 041050 (2015).
  * Tonielli  _et al._ [2020] F. Tonielli, J. C. Budich, A. Altland, and S. Diehl, Topological Field Theory Far from Equilibrium, Physical Review Letters 124, 240404 (2020), 1911.07834 .
  * Vogl  _et al._ [2021] M. Vogl, M. Rodriguez-Vega, B. Flebus, A. H. MacDonald, and G. A. Fiete, Floquet engineering of topological transitions in a twisted transition metal dichalcogenide homobilayer, Physical Review B 103, 014310 (2021), 2011.10948 .
  * Nuske  _et al._ [2020] M. Nuske, L. Broers, B. Schulte, G. Jotzu, S. A. Sato, A. Cavalleri, A. Rubio, J. W. McIver, and L. Mathey, Floquet dynamics in light-driven solids, Physical Review Research 2, 043408 (2020), 2005.10824 .
  * Shan  _et al._ [2021] J.-Y. Shan, M. Ye, H. Chu, S. Lee, J.-G. Park, L. Balents, and D. Hsieh, Giant modulation of optical nonlinearity by Floquet engineering, Nature 600, 235 (2021).
  * Chaudhary  _et al._ [2020] S. Chaudhary, A. Ron, D. Hsieh, and G. Refael, Controlling ligand-mediated exchange interactions in periodically driven magnetic materials, arXiv (2020), 2009.00813 .
  * Mikhaylovskiy  _et al._ [2020] R. V. Mikhaylovskiy, T. J. Huisman, V. A. Gavrichkov, S. I. Polukeev, S. G. Ovchinnikov, D. Afanasiev, R. V. Pisarev, T. Rasing, and A. V. Kimel, Resonant Pumping of $d$-$d$ Crystal Field Electronic Transitions as a Mechanism of Ultrafast Optical Control of the Exchange Interactions in Iron Oxides, Physical Review Letters 125, 157201 (2020).
  * Dolgirev  _et al._ [2020] P. E. Dolgirev, M. H. Michael, A. Zong, N. Gedik, and E. Demler, Self-similar dynamics of order parameter fluctuations in pump-probe experiments, Physical Review B 101, 174306 (2020), 1910.02518 .
  * Nava  _et al._ [2018] A. Nava, C. Giannetti, A. Georges, E. Tosatti, and M. Fabrizio, Cooling quasiparticles in A3C60 fullerides by excitonic mid-infrared absorption, Nature Physics 14, 154 (2018), 1704.05613 .
  * Fabrizio [2018] M. Fabrizio, Selective Transient Cooling by Impulse Perturbations in a Simple Toy Model, Physical Review Letters 120, 220601 (2018), 1801.06622 .
  * Werner  _et al._ [2019] P. Werner, M. Eckstein, M. Müller, and G. Refael, Light-induced evaporative cooling of holes in the Hubbard model, Nature Communications 10, 5556 (2019), 1904.00822 .
  * Kugel’ and Khomskii [1982] K. Kugel’ and D. Khomskii, The Jahn-Teller effect and magnetism: transition metal compounds, Sov. Phys. Usp. 25, 231 (1982).
  * Müller  _et al._ [2021] A. Müller, F. Grandi, and M. Eckstein, Ultrafast control of spin-orbital separation probed with time-resolved RIXS, arXiv (2021), 2111.12514 .
  * Wohlfeld  _et al._ [2011] K. Wohlfeld, M. Daghofer, S. Nishimoto, G. Khaliullin, and J. v. d. Brink, Intrinsic Coupling of Orbital Excitations to Spin Fluctuations in Mott Insulators, Physical Review Letters 107, 147201 (2011).
  * Barbeau  _et al._ [2019] M. Barbeau, M. Eckstein, M. Katsnelson, and J. Mentink, Optical control of competing exchange interactions and coherent spin-charge coupling in two-orbital Mott insulators, SciPost Physics 6, 027 (2019), 1803.03796 .
  * Liu  _et al._ [2018] J. Liu, K. Hejazi, and L. Balents, Floquet Engineering of Multiorbital Mott Insulators: Applications to Orthorhombic Titanates, Physical Review Letters 121, 107201 (2018), 1801.00401 .
  * Li  _et al._ [2018] J. Li, H. U. R. Strand, P. Werner, and M. Eckstein, Theory of photoinduced ultrafast switching to a spin-orbital ordered hidden phase, Nature Communications 9, 4581 (2018), 1806.01567 .
  * Rodriguez-Nieva  _et al._ [2022a] J. F. Rodriguez-Nieva, A. Schuckert, D. Sels, M. Knap, and E. Demler, Transverse instability and universal decay of spin spiral order in the Heisenberg model, Physical Review B 105, L060302 (2022a).
  * Rodriguez-Nieva  _et al._ [2022b] J. F. Rodriguez-Nieva, D. Podolsky, and E. Demler, Probing hydrodynamic sound modes in magnon fluids using spin magnetometers, Phys. Rev. B 105, 174412 (2022b).
  * Bhattacharyya  _et al._ [2020] S. Bhattacharyya, J. F. Rodriguez-Nieva, and E. Demler, Universal prethermal dynamics in Heisenberg ferromagnets, Phys. Rev. Lett. 125, 230601 (2020).
  * Singh  _et al._ [2010] D. K. Singh, B. Kamble, and A. Singh, Orbital fluctuations, spin-orbital coupling, and anomalous magnon softening in an orbitally degenerate ferromagnet, Physical Review B 81, 064430 (2010), 0910.5321 .
  * Khaliullin and Maekawa [2000] G. Khaliullin and S. Maekawa, Orbital Liquid in Three-Dimensional Mott Insulator: LaTiO3, Physical Review Letters 85, 3950 (2000), cond-mat/0008032 .
  * Khaliullin and Okamoto [2002] G. Khaliullin and S. Okamoto, Quantum Behavior of Orbitals in Ferromagnetic Titanates: Novel Orderings and Excitations, Physical Review Letters 89, 167201 (2002), cond-mat/0209434 .
  * Khaliullin and Okamoto [2003] G. Khaliullin and S. Okamoto, Theory of orbital state and spin interactions in ferromagnetic titanates, Physical Review B 68, 205109 (2003), cond-mat/0510175 .
  * Oleś  _et al._ [2005] A. M. Oleś, G. Khaliullin, P. Horsch, and L. F. Feiner, Fingerprints of spin-orbital physics in cubic Mott insulators: Magnetic exchange interactions and optical spectral weights, Physical Review B 72, 214431 (2005), cond-mat/0601381 .
  * Oleś  _et al._ [2006] A. M. Oleś, P. Horsch, L. F. Feiner, and G. Khaliullin, Spin-Orbital Entanglement and Violation of the Goodenough-Kanamori Rules, Physical Review Letters 96, 147205 (2006), cond-mat/0604417 .
  * Brzezicki  _et al._ [2015] W. Brzezicki, A. M. Oleś, and M. Cuoco, Spin-Orbital Order Modified by Orbital Dilution in Transition-Metal Oxides: From Spin Defects to Frustrated Spins Polarizing Host Orbitals, Physical Review X 5, 011037 (2015), 1408.1838 .
  * Brink  _et al._ [1999] J. v. d. Brink, P. Horsch, F. Mack, and A. M. Oleś, Orbital dynamics in ferromagnetic transition-metal oxides, Physical Review B 59, 6795 (1999), cond-mat/9812123 .
  * Khaliullin [2001] G. Khaliullin, Order from disorder: Quantum spin gap in magnon spectra of LaTiO3, Physical Review B 64, 212405 (2001), cond-mat/0108016 .
  * Feiner  _et al._ [1997] L. F. Feiner, A. M. Oleś, and J. Zaanen, Quantum Melting of Magnetic Order due to Orbital Fluctuations, Physical Review Letters 78, 2799 (1997), cond-mat/9701098 .
  * Mochizuki and Imada [2001] M. Mochizuki and M. Imada, Magnetic and Orbital States and Their Phase Transition of the Perovskite-Type Ti Oxides: Strong Coupling Approach, Journal of the Physical Society of Japan 70, 1777 (2001), cond-mat/0102223 .
  * Pal  _et al._ [2018] P. Pal, X. Liu, M. Kareev, D. Choudhury, and J. Chakhalian, Emergent antiferromagnetism of YTiO3 in YTiO3-CaTiO3 superlattices, Physical Review B 98, 045420 (2018), 1804.10333 .
  * Wang  _et al._ [2021] Z. Wang, D. Gautreau, T. Birol, and R. M. Fernandes, Strain-tunable metamagnetic critical end-point in Mott insulating rare-earth titanates, arXiv (2021), 2105.01559 .
  * Katsufuji  _et al._ [1997] T. Katsufuji, Y. Taguchi, and Y. Tokura, Transport and magnetic properties of a Mott-Hubbard system whose bandwidth and band filling are both controllable: R1-xCaxTiO3+y/2, Physical Review B 56, 10145 (1997).
  * Pavarini  _et al._ [2004] E. Pavarini, S. Biermann, A. Poteryaev, A. I. Lichtenstein, A. Georges, and O. K. Andersen, Mott Transition and Suppression of Orbital Fluctuations in Orthorhombic 3d1 Perovskites, Physical Review Letters 92, 4 (2004), cond-mat/0309102 .
  * Pavarini  _et al._ [2005] E. Pavarini, A. Yamasaki, J. Nuss, and O. K. Andersen, How chemistry controls electron localization in 3d1 perovskites: a Wannier-function study, New Journal of Physics 7, 188 (2005), cond-mat/0504034 .
  * Zhang  _et al._ [2020] X.-J. Zhang, E. Koch, and E. Pavarini, Origin of orbital ordering in YTiO3 and LaTiO3, Physical Review B 102, 035113 (2020).
  * Mochizuki and Imada [2004] M. Mochizuki and M. Imada, Orbital physics in the perovskite Ti oxides, New Journal of Physics 6, 154 (2004), cond-mat/0410650 .
  * Solovyev [2006] I. V. Solovyev, Lattice distortion and magnetism of 3d-t2g perovskite oxides, Phys. Rev. B 74, 054412 (2006), cond-mat/0602035 .
  * Solovyev [2009] I. Solovyev, Superexchange interactions in orthorhombically distorted titanates $R$TiO3 ($R$ = Y, Gd, Sm, La), N.J. Phys. 11, 093003 (2009), 0906.2035 .
  * Mizokawa  _et al._ [1999] T. Mizokawa, D. I. Khomskii, and G. A. Sawatzky, Interplay between orbital ordering and lattice distortions in LaMnO3, YVO3, and YTiO3, Physical Review B 60, 7309 (1999), cond-mat/9906010 .
  * Ishihara  _et al._ [2002] S. Ishihara, T. Hatakeyama, and S. Maekawa, Magnetic ordering, orbital ordering, and resonant x-ray scattering in perovskite titanates, Physical Review B 65, 064442 (2002), cond-mat/0111249 .
  * Ishihara [2004] S. Ishihara, Orbital wave and its observation in orbital-ordered titanates and vanadates, Physical Review B 69, 075118 (2004).
  * Ulrich  _et al._ [2008] C. Ulrich, G. Ghiringhelli, A. Piazzalunga, L. Braicovich, N. B. Brookes, H. Roth, T. Lorenz, and B. Keimer, Orbital excitations in YTiO3 and LaTiO3 probed by resonant inelastic soft x-ray scattering, Phys. Rev. B 77, 113102 (2008).
  * Ulrich  _et al._ [2009] C. Ulrich, L. J. P. Ament, G. Ghiringhelli, L. Braicovich, M. M. Sala, N. Pezzotta, T. Schmitt, G. Khaliullin, J. v. d. Brink, H. Roth, T. Lorenz, and B. Keimer, Momentum Dependence of Orbital Excitations in Mott-Insulating Titanates, Physical Review Letters 103, 107205 (2009), 0908.2497 .
  * Krivenko [2012] S. Krivenko, Interplay between crystal-field splitting and superexchange interaction in the t2g orbital Mott insulator, Physical Review B 85, 064406 (2012).
  * Maimone  _et al._ [2018] D. T. Maimone, A. B. Christian, J. J. Neumeier, and E. Granado, Coupling of phonons with orbital dynamics and magnetism in CuSb2O6, Physical Review B 97, 174415 (2018), 1803.00559 .
  * Ulrich  _et al._ [2002] C. Ulrich, G. Khaliullin, S. Okamoto, M. Reehuis, A. Ivanov, H. He, Y. Taguchi, Y. Tokura, and B. Keimer, Magnetic Order and Dynamics in an Orbitally Degenerate Ferromagnetic Insulator, Physical Review Letters 89, 167202/1 (2002), cond-mat/0209400 .
  * Sugai  _et al._ [2006] S. Sugai, A. Kikuchi, and Y. Mori, Raman scattering of orbital waves in YTiO3, Physical Review B 73, 22 (2006).
  * Ulrich  _et al._ [2006] C. Ulrich, A. Gössling, M. Grüninger, M. Guennou, H. Roth, M. Cwik, T. Lorenz, G. Khaliullin, and B. Keimer, Raman Scattering in the Mott Insulators LaTiO3 and YTiO3: Evidence for Orbital Excitations, Physical Review Letters 97, 157401 (2006).
  * Marzari and Vanderbilt [1997] N. Marzari and D. Vanderbilt, Maximally Localized Generalized Wannier Functions for Composite Energy Bands, Phys. Rev. B 56, 12847 (1997).
  * Kamenev [2011] A. Kamenev,  _Field Theory of Non-Equilibrium Systems_ (Cambridge University Press, Cambridge, 2011).
  * Note [1] It also seems like there is a good amount of uncertainty about the resonance frequency as well, and whether features are one-orbital, two-orbital, or of other origin all together.
  * Note [2] In fact, $A_{d}\sim Q_{\rm IR}^{2}$ is itself scaling linearly with fluence, so that scaling quadratically with $A_{d}^{2}$ implies quadratic fluence dependence.
  * Cavalleri [2018] A. Cavalleri, Photo-induced superconductivity, Contemp. Phys. 59, 31 (2018).
  * Buzzi  _et al._ [2020] M. Buzzi, D. Nicoletti, M. Fechner, N. Tancogne-Dejean, M. A. Sentef, A. Georges, T. Biesner, E. Uykur, M. Dressel, A. Henderson, T. Siegrist, J. A. Schlueter, K. Miyagawa, K. Kanoda, M.-S. Nam, A. Ardavan, J. Coulthard, J. Tindall, F. Schlawin, D. Jaksch, and A. Cavalleri, Photomolecular High-Temperature Superconductivity, Physical Review X 10, 031028 (2020), 2001.05389 .
  * Hoegen  _et al._ [2022] A. v. Hoegen, M. Fechner, M. Först, N. Taherian, E. Rowe, A. Ribak, J. Porras, B. Keimer, M. Michael, E. Demler, and A. Cavalleri, Amplification of Superconducting Fluctuations in Driven YBa2Cu3O6+x, Physical Review X 12, 031008 (2022).
  * Jotzu  _et al._ [2021] G. Jotzu, G. Meier, A. Cantaluppi, A. Cavalleri, D. Pontiroli, M. Riccò, A. Ardavan, and M.-S. Nam, Superconducting fluctuations observed far above Tc in the isotropic superconductor K3C60, arXiv (2021), 2109.08679 .
  * Larkin and Varlamov [2005] A. Larkin and A. Varlamov,  _Theory of Fluctuations in Superconductors_ (Oxford University Press, 2005).
  * Imada  _et al._ [1998] M. Imada, A. Fujimori, and Y. Tokura, Metal-insulator transitions, Reviews of Modern Physics 70, 1039 (1998).
  * Sayyad and Eckstein [2016] S. Sayyad and M. Eckstein, Slowdown of the Electronic Relaxation Close to the Mott Transition, Physical Review Letters 117, 096403 (2016), 1601.02994 .
  * Schuckert  _et al._ [2018] A. Schuckert, A. P. Orioli, and J. Berges, Nonequilibrium quantum spin dynamics from two-particle irreducible functional integral techniques in the Schwinger boson representation, Physical Review B 98, 224304 (2018), 1806.02347 .
  * Babadi  _et al._ [2015] M. Babadi, E. Demler, and M. Knap, Far-from-Equilibrium Field Theory of Many-Body Quantum Spin Systems: Prethermalization and Relaxation of Spin Spiral States in Three Dimensions, Phys. Rev. X 5, 041005 (2015), 1504.05956 .
  * Kiselev and Oppermann [2000] M. N. Kiselev and R. Oppermann, Schwinger-Keldysh Semionic Approach for Quantum Spin Systems, Phys. Rev. Lett. 85, 5631 (2000), cond-mat/0008060 .
  * Shnirman and Makhlin [2003] A. Shnirman and Y. Makhlin, Spin-Spin Correlators in the Majorana Representation, Physical Review Letters 91, 207204 (2003), cond-mat/0305064 .
  * Lemonik and Mitra [2019] Y. Lemonik and A. Mitra, Transport and spectral signatures of transient fluctuating superfluids in the absence of long-range order, Phys. Rev. B 100, 10.1103/physrevb.100.094503 (2019), 1905.08906 .
  * Lemonik and Mitra [2018] Y. Lemonik and A. Mitra, Quench dynamics of superconducting fluctuations and optical conductivity in a disordered system, Phys. Rev. B 98, 214514 (2018), 1804.09280 .
  * Brink  _et al._ [2002] J. v. d. Brink, G. Khaliullin, and D. Khomskii, Orbital effects in manganites (Springer, Dordrecht, 2002) pp. 263–301.
  * Feiner and Oleś [1999] L. F. Feiner and A. M. Oleś, Electronic origin of magnetic and orbital ordering in insulating LaMnO3, Physical Review B 59, 3295 (1999), cond-mat/9805011 .
  * Fujioka  _et al._ [2010] J. Fujioka, T. Yasue, S. Miyasaka, Y. Yamasaki, T. Arima, H. Sagayama, T. Inami, K. Ishii, and Y. Tokura, Critical competition between two distinct orbital-spin ordered states in perovskite vanadates, Physical Review B 82, 144425 (2010).
  * Khaliullin  _et al._ [2000] G. Khaliullin, P. Horsch, and A. M. Oleś, Spin Order due to Orbital Fluctuations: Cubic Vanadates, Physical Review Letters 86, 3879 (2000), cond-mat/0103300 .
  * Rückamp  _et al._ [2005] R. Rückamp, E. Benckiser, M. W. Haverkort, H. Roth, T. Lorenz, A. Freimuth, L. Jongen, A. Möller, G. Meyer, P. Reutler, B. Büchner, A. Revcolevschi, S.-W. Cheong, C. Sekar, G. Krabbes, and M. Grüninger, Optical study of orbital excitations in transition-metal oxides, New Journal of Physics 7, 144 (2005), cond-mat/0503405 .
  * Shi  _et al._ [2019] X. Shi, W. You, Y. Zhang, Z. Tao, P. M. Oppeneer, X. Wu, R. Thomale, K. Rossnagel, M. Bauer, H. Kapteyn, and M. Murnane, Ultrafast electron calorimetry uncovers a new long-lived metastable state in 1T-TaSe2 mediated by mode-selective electron-phonon coupling, Science Advances 5, eaav4449 (2019), 1901.08214 .
  * Genske and Rosch [2015] M. Genske and A. Rosch, Floquet-Boltzmann equation for periodically driven Fermi systems, Physical Review A 92, 062108 (2015), 1508.04551 .
  * Mitrano and Wang [2020] M. Mitrano and Y. Wang, Probing light-driven quantum materials with ultrafast resonant inelastic X-ray scattering, Communications Physics 3, 184 (2020), 2009.11315 .
  * Georges  _et al._ [2013] A. Georges, L. d. Medici, and J. Mravlje, Strong Correlations from Hund’s Coupling, Annual Review of Condensed Matter Physics 4, 137 (2013), 1207.3033 .
  * Kresse and Furthmüller [1996] G. Kresse and J. Furthmüller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B 54, 11169 (1996).
  * A. and Tanaka [2015] T. A. and I. Tanaka, First principles phonon calculations in materials science, Scripta Mat. 108, 1 (2015).
  * Kresse and Joubert [1999] G. Kresse and D. Joubert, First principles phonon calculations in materials science, Phys. Rev. B 59, 1758 (1999).
  * Monkhorst and Pack [1976] H. Monkhorst and J. Pack, Special points for Brillouin-zone integrations, Phys. Rev. B 13, 5188 (1976).

##  Appendix A Orbital Angular Momentum

Here we present details of the orbital angular momentum projected onto the
$t_{2g}$ states. The orbital angular momentum of the full $d$ shell is in
general characterized by a total $L=2$ operator. In the presence of a cubic
crystal field splitting, this is further split into two $e_{g}$ levels with
quenched angular momentum, and three $t_{2g}$ levels which in general may have
an unquenched effective angular momentum of $L_{\rm eff}=1$ roughly
corresponding to the vector representation of the three orbitals [118].

Explicitly, we have the representation of the effective angular momentum
operators given in terms of the $t_{2g}$ states
$|a\rangle=|yz\rangle,|b\rangle=|zx\rangle,|c\rangle=|xy\rangle$ on site $j$
as

$\hat{L}_{j}^{l}=-i\epsilon_{lmn}\left(|n;j\rangle\langle
m;j|-|m;j\rangle\langle n;j|\right).$ (31)

Note that this is odd under time-reversal symmetry and has purely off-diagonal
matrix elements written in terms of the Cartesian orbitals.

In the presence of further splitting induced by the GdFeO3 lattice distortion,
the orbital angular momentum becomes quenched. We can then project it onto the
lowest two states of the intra-$t_{2g}$ distortion, which we call $|0\rangle$
and $|1\rangle$. If we express the angular orbital momentum in the basis of
the crystal-field matrix we find the relevant operator is

$\mathbf{\hat{L}}_{j}=\left(-i|1;j\rangle\langle 0;j|-|0;j\rangle\langle
1;j|\right)\mathbf{e}_{2}=\mathbf{{e}}_{2}\hat{\tau}_{j}^{2},$ (32)

where the second equality expresses this in a vector representation in terms
of the unit vector which points along the direction $\mathbf{e}_{z}$ for the
$xy$ orbital, and so on in a cyclic way. The Pauli matrix $\hat{\tau}_{j}^{2}$
is the relevant second-quantized operator resulting from the projection of the
full orbital state operator onto the two lowest orbitals. More generally, we
may have the crystal-field matrix evolving with a parameter (such as a phonon
coordinate), in which case we express this in terms of the relevant crystal-
field wavefunctions as

$\mathbf{\hat{L}}_{j}={\mathbf{n}}\hat{\tau}^{2}_{j},$ (33)

where the matrix elements of the orbital angular momentum vector are obtained
via the Levi-Civita symbol as

$\mathbf{n}_{c}=\mathbf{e}_{a}\times\mathbf{e}_{b}$ (34)

where $\mathbf{e}_{a}$ is the $a^{\rm th}$ (real, by time-reversal symmetry)
eigenvector of the crystal-field.

##  Appendix B Keldysh Holstein-Primakoff

We introduce the spin fields on the forward and backward contours, expanded up
to $O(S^{0})$ as

$\mathbf{S}_{j,\pm}=S\mathbf{e}_{3}+\sqrt{S}\left[\mathbf{e}_{-}b_{j\pm}+\mathbf{e}_{+}\overline{b}_{j\pm}\right]-\mathbf{e}_{3}\overline{b}_{j\pm}b_{j\pm}$
(35)

with

$\mathbf{e}_{\pm}=\frac{1}{\sqrt{2}}(\mathbf{e}_{1}\pm i\mathbf{e}_{2}).$ (36)

We therefore find the expansion for the classical/quantum fields of

$\displaystyle\mathbf{S}_{j,cl}=S\frac{1}{\sqrt{2}}\mathbf{e}_{3}+\sqrt{S}\left[\mathbf{e}_{-}b_{jcl}+\mathbf{e}_{+}\overline{b}_{jcl}\right]-\mathbf{e}_{3}\frac{1}{\sqrt{2}}(\overline{b}_{jccl}b_{jcl}+\overline{b}_{jq}b_{jq})$
(37a)
$\displaystyle\mathbf{S}_{j,q}=\sqrt{S}\left[\mathbf{e}_{-}b_{jq}+\mathbf{e}_{+}\overline{b}_{jq}\right]-\mathbf{e}_{3}\frac{1}{\sqrt{2}}(\overline{b}_{jcl}b_{jcl}+\overline{b}_{jq}b_{jcl}).$
(37b)

We require the products
$\mathbf{S}^{\alpha}_{jcl}(t)\mathbf{S}^{\beta}_{jq}(t^{\prime})$ and
$\mathbf{S}^{\alpha}_{jq}(t)\mathbf{S}^{\beta}_{jq}(t^{\prime})$ up to order
$O(S)$. We find at quadratic order the contributions

$\mathbf{S}_{jq}(t)\mathbf{S}_{jcl}(t^{\prime})=S\left[\mathbb{T}_{++}\overline{b}_{jq}(t)\overline{b}_{jcl}(t^{\prime})+\mathbb{T}_{+-}\overline{b}_{jq}(t)b_{jcl}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jcl}(t^{\prime})+\mathbb{T}_{--}{b}_{jq}(t)b_{jcl}(t^{\prime})\right]\\\
-\frac{S}{2}\mathbb{T}_{33}\left(\overline{b}_{jcl}(t)b_{jq}(t)+\overline{b}_{jq}(t)b_{jcl}(t)\right),$
(38)

and

$\mathbf{S}_{jq}(t)\mathbf{S}_{jq}(t^{\prime})=S\left[\mathbb{T}_{++}\overline{b}_{jq}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{+-}\overline{b}_{jq}(t)b_{jq}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{--}{b}_{jq}(t)b_{jq}(t^{\prime})\right].$
(39)

Here we have introduced the tensors

$\displaystyle\mathbb{T}_{++}=\mathbf{e}_{+}\otimes\mathbf{e}_{+}$ (40a)
$\displaystyle\mathbb{T}_{+-}=\mathbf{e}_{+}\otimes\mathbf{e}_{-}$ (40b)
$\displaystyle\mathbb{T}_{-+}=\mathbf{e}_{-}\otimes\mathbf{e}_{+}$ (40c)
$\displaystyle\mathbb{T}_{--}=\mathbf{e}_{-}\otimes\mathbf{e}_{-}$ (40d)
$\displaystyle\mathbb{T}_{33}=\mathbf{e}_{3}\otimes\mathbf{e}_{3}.$ (40e)

We have an effective action due to the orbital bath of

$S_{\rm
eff}=-S\frac{\lambda^{2}}{2}\sum_{j}\int_{t,t^{\prime}}\mathbf{tr}\bigg{\\{}\mathbf{n}_{j}(t^{\prime})\otimes\mathbf{n}_{j}(t)\cdot\bigg{[}D^{K}(t,t^{\prime})\left(\mathbb{T}_{++}\overline{b}_{jq}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{+-}\overline{b}_{jq}(t){b}_{jq}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{--}{b}_{jq}(t){b}_{jq}(t^{\prime})\right)\\\
+D^{R}(t,t^{\prime})\left(\mathbb{T}_{++}\overline{b}_{jq}(t)\overline{b}_{jcl}(t^{\prime})+\mathbb{T}_{+-}\overline{b}_{jq}(t){b}_{jcl}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jcl}(t^{\prime})+\mathbb{T}_{--}{b}_{jq}(t){b}_{jcl}(t^{\prime})-\frac{1}{2}\mathbb{T}_{33}\left(\overline{b}_{jcl}(t)b_{jq}(t)+\overline{b}_{jq}(t)b_{jcl}(t)\right)\right)\\\
+D^{A}(t,t^{\prime})\left(\mathbb{T}_{++}\overline{b}_{jcl}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{+-}\overline{b}_{jcl}(t){b}_{jq}(t^{\prime})+\mathbb{T}_{-+}{b}_{jcl}(t)\overline{b}_{jq}(t^{\prime})+\mathbb{T}_{--}{b}_{jcl}(t){b}_{jq}(t^{\prime})-\frac{1}{2}\mathbb{T}_{33}\left(\overline{b}_{jcl}(t^{\prime})b_{jq}(t^{\prime})+\overline{b}_{jq}(t^{\prime})b_{jcl}(t^{\prime})\right)\right)\bigg{]}\bigg{\\}}.$
(41)

Here the trace is taken over the spin tensor indices. We will herein replace
the sublattice specific angular momentum vectors $\mathbf{n}$ with the
sublattice averaged matrix

$\mathbb{N}(t^{\prime},t)=\overline{\mathbf{n}_{j}(t^{\prime})\otimes\mathbf{n}_{j}(t)}.$
(42)

We also will for the most part throw away the anomalous correlations, assuming
they are small, though this may be an interesting direction for the future. We
then find an effective $U(1)$ symmetry for the magnons, getting

$S_{\rm
eff}=-S\frac{\lambda^{2}}{2}\sum_{j}\int_{t,t^{\prime}}\mathbf{tr}\bigg{\\{}\mathbb{N}(t^{\prime},t)\cdot\bigg{[}D^{K}(t,t^{\prime})\left(\mathbb{T}_{+-}\overline{b}_{jq}(t){b}_{jq}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jq}(t^{\prime})\right)\\\
+D^{R}(t,t^{\prime})\left(\mathbb{T}_{+-}\overline{b}_{jq}(t){b}_{jcl}(t^{\prime})+\mathbb{T}_{-+}{b}_{jq}(t)\overline{b}_{jcl}(t^{\prime})-\frac{1}{2}\mathbb{T}_{33}\left(\overline{b}_{jcl}(t)b_{jq}(t)+\overline{b}_{jq}(t)b_{jcl}(t)\right)\right)\\\
+D^{A}(t,t^{\prime})\left(\mathbb{T}_{+-}\overline{b}_{jcl}(t){b}_{jq}(t^{\prime})+\mathbb{T}_{-+}{b}_{jcl}(t)\overline{b}_{jq}(t^{\prime})-\frac{1}{2}\mathbb{T}_{33}\left(\overline{b}_{jcl}(t^{\prime})b_{jq}(t^{\prime})+\overline{b}_{jq}(t^{\prime})b_{jcl}(t^{\prime})\right)\right)\bigg{]}\bigg{\\}}.$
(43)

At this point we can read out the retarded self-energy

$\Sigma^{R}(t,t^{\prime})=\frac{S\lambda^{2}}{2}\left[N_{+-}(t^{\prime},t)D^{R}(t,t^{\prime})+N_{-+}(t,t^{\prime})D^{A}(t^{\prime},t)\right]\\\
-\frac{S\lambda^{2}}{4}\left(\int
dt^{\prime\prime}N_{33}(t^{\prime\prime},t)D^{R}(t,t^{\prime\prime})\delta(t-t^{\prime})-\frac{1}{2}\int
dt^{\prime\prime}N_{33}(t^{\prime},t^{\prime\prime})D^{A}(t^{\prime\prime},t^{\prime})\delta(t-t^{\prime})\right).$
(44)

and the Keldysh self energy as

$\Sigma^{K}(t,t^{\prime})=c\frac{S\lambda^{2}}{2}\left[N_{+-}(t^{\prime},t)D^{K}(t,t^{\prime})+N_{-+}(t,t^{\prime})D^{K}(t^{\prime},t)\right].$
(45)

In the retarded self-energy, the last two terms describe the drive-induced
dephasing ($T_{2}$ process), which only enters when the effective magnon gap
due to the orbital fluctuations is time-dependent. We will leave this study to
future works, and ignore it in this case as we assume the projection of the
angular momentum matrix elements is small along the $c$ axis.

To summarize, once we discard the anomalous terms and the pump-induced
dephasing we are left with the magnon self-energies of

$\displaystyle\Sigma^{R}(t,t^{\prime})=\frac{S\lambda^{2}}{2}\left[N_{+-}(t^{\prime},t)D^{R}(t,t^{\prime})+N_{-+}(t,t^{\prime})D^{A}(t^{\prime},t)\right]$
(46a)
$\displaystyle\Sigma^{A}(t,t^{\prime})=\frac{S\lambda^{2}}{2}\left[N_{+-}(t^{\prime},t)D^{A}(t,t^{\prime})+N_{-+}(t,t^{\prime})D^{R}(t^{\prime},t)\right]$
(46b)
$\displaystyle\Sigma^{K}(t,t^{\prime})=\frac{S\lambda^{2}}{2}\left[N_{+-}(t^{\prime},t)D^{K}(t,t^{\prime})+N_{-+}(t,t^{\prime})D^{K}(t^{\prime},t)\right].$
(46c)

We are now tasked with using these to solve the equations of motion in the
driven case. Note that

$N_{+-}(t^{\prime},t)=\overline{\mathbf{e}_{+}\cdot\mathbf{n}(t^{\prime})\mathbf{n}(t)\cdot\mathbf{e}_{-}}=N_{-+}(t,t^{\prime}).$
(47)

##  Appendix C Density Functional Theory Calculations

We performed our computations with the Vienna ab-initio simulation package
VASP.6.2 [119]. For the phonon calculations we used the Phonopy software
package [120] and the Wannier90 package for Wannierization [93]. Our
computations further utilized pseudopotentials generated within the Projected
Augmented Wave (PAW) [121] method. Specifically, we take the following
configurations for default potentials: Ti $3p^{6}4s^{1}3d^{3}$, Y
$4s^{2}4p^{6}5s^{2}4d^{1}$, and O $2s^{2}2p^{4}$. We applied the Local Spin
Density Approximation (LsDA) approximation for the exchange-correlation
potential, which we augment with the Hubbard $U-J$ parameter to account for
the localized nature of the d-states of Ti. We use $U =4$ eV and $J =0.0$ eV.
As a numerical setting, we used a $9\times 9\times 7$ Monkhorst [122]
generated $k$-point-mesh sampling of the Brillouin zone and a plane-wave
energy cutoff of 600 eV. We iterate self-consistent calculations until the
change in total energy has converged up to $10^{-8}$ eV.

##  Appendix D Two-Time Equations

Here we elaborate on the details associated to computing the various non-
equilibrium Green’s function which enter into the magnon kinetic equation. We
focus on the time-frequency domain transforms, assuming the space and momentum
dependencies are trivial.

To begin with, we invoke the formula for the relation of the Wigner transform
of two products. We consider two correlation functions with known Wigner
transforms $A(T_{1};\omega_{1})$ and $B(T_{2};\omega_{2})$. We want the Wigner
transform of their convolution, expressed in terms of the two-time functions
$A,B$ as

$C(T;\Omega)=\int d\tau e^{i\Omega\tau}\int dtA(T+\tau/2,t)B(t,T-\tau/2).$
(48)

This expression can be found in [94] and is formally given as an exponential
derivative operation as

$C(T;\Omega)=A(T;\Omega)\exp\left(-\frac{i}{2}\left[\overleftarrow{\partial_{T}}\overrightarrow{\partial_{\Omega}}-\overleftarrow{\partial_{\Omega}}\overrightarrow{\partial_{T}}\right]\right)B(T;\Omega).$
(49)

This is only useful if one can expand the relevant functions in terms of
slowly-varying in both time and frequency, which in turn relies on a
separation of scales between the dynamics and frequencies.

We also use the related formula, relevant for the Wigner transform of the
point-wise product,

$D(T;\omega)=\int d\tau
e^{i\Omega\tau}A(T+\frac{\tau}{2},T-\frac{\tau}{2})B(T+\frac{\tau}{2},T-\frac{\tau}{2}),$
(50)

which yields

$D(T;\Omega)=\int\frac{d\omega}{2\pi}A(T;\Omega-\omega)B(T;\omega).$ (51)

We now apply this to the Green’s functions. First, we consider the magnon
retarded Green’s function, which obeys the integral equation

$\left(i\partial_{t}-\Omega_{\bf p}\right)G^{R}(t,t^{\prime})-\int
dt^{\prime\prime}\Sigma^{R}(t,t^{\prime\prime})G^{R}(t^{\prime\prime},t^{\prime})=\delta(t-t^{\prime}).$
(52)

In the absence of the drive, this is solved in the frequency domain, and we
obtain the standard result which in particular amounts to a form of Gilbert
damping at low frequencies.

In this work we will still retain the separation between the evolution times,
which are of order of 20-2000 ps, and the time-scales of the internal degrees-
of-freedom which are from 20-800 fs or so. This allows us to efficiently
employ the equations of motion using the Wigner transformations and the Moyal
expansions.

The lowest order in the Moyal expansion is simply the product. We retain
expansion up to first order, giving equation of motion for the retarded
Green’s function

$\frac{i}{2}\left[1-\partial_{\omega}\Sigma^{R}(T;\omega)\right]\partial_{T}G^{R}_{\bf
p}(T;\omega)+\left[\omega-\Omega_{\bf
p}-\Sigma^{R}(T:\omega)+\frac{i}{2}\partial_{T}\Sigma^{R}(T;\omega)\partial_{\omega}\right]G^{R}_{\bf
p}(T;\omega)=\mathds{1}.$ (53)

This yields a first order differential equation for the Green’s function,
though it remains non-local in frequency space due to the the changing self-
energy. When solving, we also supplement with the initial condition that

$G^{R}_{\bf p}(-\infty.\omega)=\frac{1}{\omega-\Omega_{\bf
p}-\Sigma^{R}(\omega)}.$ (54)

To complete this, we need to express the self-energy as a Wigner transform as
well. We use the product formula to find Wigner-transform (applied to $R,A,K$
self-energies)

$\check{\Sigma}(T;\Omega)=\frac{1}{N^{\rm
eq}_{+-}}\int\frac{d\omega}{2\pi}N_{+-}(T;\omega)\check{\Sigma}_{\rm
eq}(\Omega-\omega).$ (55)

Here $\Sigma_{\rm eq}(\omega)$ is the equilibrium self-energy and depends only
on frequency. $N_{+-}^{\rm eq}$ is the equilibrium angular momentum
projection, while $N_{+-}(T;\omega)$ is the Wigner transform of the modulated
angular momentum tensor. We model the modulation via

$N_{+-}(t,t^{\prime})=N_{+-}^{\rm
eq}(A_{0}(t)+A_{1}(t)e^{-i\Omega_{d}t}+A_{-1}(t)e^{i\Omega_{d}t})\\\
\times(A_{0}^{*}(t^{\prime})+A_{1}^{*}(t^{\prime})e^{i\Omega_{d}t^{\prime}}+A_{-1}^{*}(t^{\prime})e^{-i\Omega_{d}t^{\prime}}),$
(56)

where we have expressed this in terms of a Floquet expansion in the drive-
frequency $\Omega_{d}$, along with slowly varying envelope functions
$A_{0},A_{\pm 1}$, which vary over times of order
$\tau_{d}\gg\Omega_{d}^{-1}$. This gives, in the slowly varying envelope
approximation for $A$’s of

$N_{+-}(t,t^{\prime})/N_{+-}^{\rm
eq}=\left[|A_{0}(T)|^{2}+A_{1}(T)A_{-1}^{*}(T)e^{-2i\Omega_{d}T}+A_{-1}(T)A_{1}^{*}(T)e^{2i\Omega_{d}T}\right]2\pi\delta(\omega)\\\
+|A_{1}(T)|^{2}2\pi\delta(\omega-\Omega_{d})+|A_{-1}(T)|^{2}2\pi\delta(\omega+\Omega_{d})\\\
+\left[A_{0}(T)A_{1}^{*}(T)e^{i\Omega_{D}T}+A_{0}^{*}(T)A_{1}(T)e^{-i\Omega_{D}T}\right]2\pi\delta(\omega-\Omega_{d}/2)\\\
+\left[A_{0}(T)A_{-1}^{*}(T)e^{-i\Omega_{D}T}+A_{0}^{*}(T)A_{-1}(T)e^{i\Omega_{D}T}\right]2\pi\delta(\omega+\Omega_{d}/2).$
(57)

This involves a number of terms, including some which couple the slow-dynamics
to the fast degrees of freedom. These terms involve oscillatory couplings like
$e^{i\Omega_{d}T}$. While these are important close to parametric resonance,
or in the steady-state Floquet system, where the separation of time scales
completely disintegrates, or must be treated non-perturbatively, we limit
ourselves to the regime where the dynamics are still able to be disentangled.
We therefore only keep in this expansion those terms which don’t average out
over long times $T$. This leaves only the terms

$\check{\Sigma}(T;\omega)=|A_{0}(T)|^{2}\check{\Sigma}(\omega)\\\
+|A_{1}(T)|^{2}\check{\Sigma}(\omega-\Omega_{d})+|A_{-1}(T)|^{2}\check{\Sigma}(\omega+\Omega_{d}).$
(58)

In fact, the object we are interested in is the magnon Keldysh occupation
function, whose equal-time value reflects the time-dependence of the total
number of magnons. At the Gaussian level, one can find that the this Green’s
function is given by

$\mathbb{G}^{K}=\mathbb{G}^{R}\circ\Sigma^{K}\circ\mathbb{G}^{A}.$ (59)

In order to proceed, we manipulate this to obtain an equation of motion of the
form

$(\mathbb{G}^{R})^{-1}\circ\mathbb{G}^{K}-\mathbb{G}^{K}\circ(\mathbb{G}^{A})^{-1}=-(\mathbb{G}^{R}\circ\Sigma^{K}-\Sigma^{K}\circ\mathbb{G}^{A}).$
(60)

We now utilize the fact that this is diagonal in momentum space and take the
Wigner transform of this equation. In general, this will not yield a closed
form since the Wigner transform is over convolutions of the functions. In the
very lowest-order limit of a slowly-varying change in the self-energy, we get

$\left[i\frac{\partial}{\partial
T}-(\Sigma^{R}(T;\omega)-\Sigma^{A}(T;\omega))\right]\mathbb{G}^{K}_{\bf
p}(T;\omega)=-\left(\mathbb{G}^{R}_{\bf p}(T;\omega)-\mathbb{G}^{A}_{\bf
p}(T;\omega)\right)\Sigma^{K}(T;\omega).$ (61)

This is formulated in terms of the occupation and spectral functions as

$\frac{\partial\mathbb{G}^{K}_{\bf p}(T;\omega)}{\partial
T}=\left[\Sigma^{R}(T;\omega)-\Sigma^{A}(T;\omega)\right]\left(-i\mathbb{G}^{K}_{\bf
p}(T;\omega)+2\pi\mathcal{A}_{\rm mag}(T;\omega,{\bf p})F_{\rm
orb}(T;\omega)\right).$ (62)

This is the simple frequency-dependent relaxation-time approximation. We find
a relaxation of the instantaneous magnon occupation towards the bath
temperature with the relaxation rate given by the bath coupling.

To conclude, we implement the quasiparticle approximation, which assumes the
linewidth of the magnon is much smaller than its central frequency. In this
case we can derive a simple equation solely for the total magnon occupation
function

$f_{\bf p}(T)=i\int\frac{d\omega}{2\pi}G^{K}_{\bf p}(T;\omega)$ (63)

as

$\partial_{t}f_{\bf p}(t)=-\frac{1}{\tau_{\bf p}(t)}\left[f_{\bf p}(t)-f_{\bf
p}^{(0)}(t)\right],$ (64)

where the instantaneous relaxation rate is given by

$\frac{1}{\tau_{\bf
p}(t)}=-2\int\frac{d\omega}{2\pi}\Im\Sigma^{R}(\omega;t)\mathcal{A}_{\rm
mag}(\omega,\mathbf{p};t),$ (65)

and the instantaneous occupation function is

$f_{\bf
p}^{(0)}(t)=\frac{\int\frac{d\omega}{2\pi}\Im\Sigma^{R}(\omega;t)\mathcal{A}_{\rm
mag}(\omega,\mathbf{p};t)F_{\rm
orb}(\omega)}{\int\frac{d\omega}{2\pi}\Im\Sigma^{R}(\omega;t)\mathcal{A}_{\rm
mag}(\omega,\mathbf{p};t)}.$ (66)

We approximate the spectral function as

$\mathcal{A}_{\rm mag}(\omega,\mathbf{p};t)=\delta(\omega-\Omega_{\bf p}(t)),$
(67)

since the quasiparticle decay rate is expected to be small.