Research ArticlePHYSICS

Dissecting spin-phonon equilibration in ferrimagnetic insulators by ultrafast lattice excitation

See allHide authors and affiliations

Science Advances  13 Jul 2018:
Vol. 4, no. 7, eaar5164
DOI: 10.1126/sciadv.aar5164


To gain control over magnetic order on ultrafast time scales, a fundamental understanding of the way electron spins interact with the surrounding crystal lattice is required. However, measurement and analysis even of basic collective processes such as spin-phonon equilibration have remained challenging. Here, we directly probe the flow of energy and angular momentum in the model insulating ferrimagnet yttrium iron garnet. After ultrafast resonant lattice excitation, we observe that magnetic order reduces on distinct time scales of 1 ps and 100 ns. Temperature-dependent measurements, a spin-coupling analysis, and simulations show that the two dynamics directly reflect two stages of spin-lattice equilibration. On the 1-ps scale, spins and phonons reach quasi-equilibrium in terms of energy through phonon-induced modulation of the exchange interaction. This mechanism leads to identical demagnetization of the ferrimagnet’s two spin sublattices and to a ferrimagnetic state of increased temperature yet unchanged total magnetization. Finally, on the much slower, 100-ns scale, the excess of spin angular momentum is released to the crystal lattice, resulting in full equilibrium. Our findings are relevant for all insulating ferrimagnets and indicate that spin manipulation by phonons, including the spin Seebeck effect, can be extended to antiferromagnets and into the terahertz frequency range.


In solids, vibrations of the crystal lattice have a significant impact on the orbital dynamics of the electrons (Fig. 1A). They strongly modify properties such as electrical conductivity and may even cause insulator-to-metal transitions (1). Likewise, the interplay between phonons and electron spins (Fig. 1A) is relevant for equally drastic phenomena including colossal magnetoresistance, femtosecond magnetization control (26), and the spin Seebeck effect (7, 8). Recently, ultrafast optical techniques have provided new insights into the coherent coupling of individual phonon and spin modes (914).

Fig. 1 Ultrafast probing of spin-phonon interactions.

(A) Experimental principle. A terahertz pump pulse resonantly and exclusively excites optical phonons of a ferrimagnet. The impact on the sample magnetization is monitored by the Faraday rotation θ of a subsequent femtosecond probe pulse. By using an electric insulator, the electronic orbital degrees of freedom remain unexcited (see red crosses). (B) Part of the unit cell of ferrimagnetic YIG. Magnetic Fe3+ ions at tetrahedral d sites and octahedral a sites comprise, respectively, the majority and minority spin sublattice of the ferrimagnet. The pump pulse resonantly excites a TO(Γ) optical phonon associated with a Fe-O stretch vibration at the tetrahedral d site.

Despite this progress, the microscopic origins of spin-phonon interactions remain an intriguing problem. Even the equilibration of crystal lattice and electron spins, arguably the conceptually simplest collective process, is far from being understood. This notion is highlighted in the model system yttrium iron garnet Y3Fe5O12 (YIG), which is one of the best-studied magnetically ordered solids (15, 16) and ubiquitous in the field of magnonics (8). Only estimates of the time constant of spin-phonon equilibration exist, and they extend over as many as six orders of magnitude from ~1 μs (17) to ~250 ps (18, 19) and down to ~1 ps (20). Here, we introduce an approach based on resonant phonon excitation that allows us to directly probe the interactions between the crystal lattice and electron spins over multiple time scales.


A schematic of our terahertz-pump magneto-optic–probe experiment is shown in Fig. 1A and detailed in Materials and Methods and fig. S1A. An incident, intense, ultrashort terahertz pump pulse (21) (photon energy of ~0.1 eV, duration of ~250 fs; see fig. S1B) selectively excites the crystal lattice by resonantly driving infrared-active transverse optical (TO) phonons. The impact on the sample’s magnetic order is monitored from femtoseconds to milliseconds by measuring the magneto-optic Faraday rotation of a time-delayed probe pulse (see Fig. 1A and fig. S1A). In addition, we measure isotropic changes in the optical transmittance of the sample on ultrafast time scales.

As samples, we choose model systems for spin-wave dynamics in magnetic insulators (7, 8, 16, 22, 23): pure YIG and bismuth/gallium-substituted YIG (BiGa:YIG). In these ferrimagnets, magnetic Fe3+ ions at a and d sites in the unit cell (Fig. 1B) comprise two nonequivalent, ferromagnetic sublattices with magnetization Ma and Md, respectively, which couple antiferromagnetically. The 2:3 ratio of a to d sites results in a nonzero net magnetization Ma + Md below the Curie temperature TC. The Faraday rotation θ of the probe pulse is determined by (24)θ=aaMa+adMd(1)where aa and ad are the sublattice magneto-optic coefficients. Compared to YIG, the Faraday effect of BiGa:YIG is enhanced by about one order of magnitude and of opposite sign (see Materials and Methods and fig. S1C) (24).

In both materials, we expect that the analysis of the pump-induced dynamics is simplified due to the following features. First, the sublattice magnetizations Ma and Md almost entirely arise from the spin rather the orbital angular momentum of the electrons. At room temperature, the g factor of YIG differs by only 0.13% from that of the free electron (25). Second, the sizeable electronic band gap (2.85 eV for YIG) implies that the electronic orbital degrees of freedom remain in their ground state, thereby significantly reducing the complexity of the pump-induced processes.

Figure 2A displays the absorptance spectrum of a 15-μm-thick BiGa:YIG film from 10 to 45 THz, where absorption is known to be due solely to infrared-active phonons (26). Our ab initio calculations show that the pump pulse centered around 21 THz (red spectrum in Fig. 2A) predominantly excites long-wavelength TO(Γ) lattice normal modes characterized by an asymmetric Fe-O stretch vibration (see Fig. 1B, fig. S4, and Materials and Methods).

Fig. 2 Ultrafast phonon-induced dynamics of magnetic order.

(A) Terahertz absorptance of the BiGa:YIG film (black solid line). Pump intensity spectra are either resonant (red) or nonresonant (blue) with the TO(Γ) phonon absorption band. Open circles show the pump-induced Faraday signal 10 ps after sample excitation as a function of the pump pulse center frequency. (B) Pump-induced change Δθ in Faraday rotation for resonant and off-resonant pumping on ultrafast and (C) microsecond time scales normalized to the equilibrium Faraday signal θ0 = θ(−2 ps). The incident fluence is 10 mJ cm−2. Panel (B) also shows the isotropic transient change in the sample transmittance for resonant pumping (thin black solid line). Dashed lines in (B) and (C) are single-exponential fits with time constants of τfast = 1.6 ps and τslow = 90 ns, respectively. The inset of (C) displays the ultrafast Faraday signal Δθ(10 ps) as a function of the incident pump-pulse fluence. Data are taken at a temperature of 296 K.


Ultrafast spin dynamics

Figure 2B shows the relative pump-induced change Δθ/θ0 = θ(t)/θ0 − 1 in the Faraday rotation as a function of the delay t since excitation. Here, θ0 = θ(−2 ps) refers to the equilibrium case. When pumping off the TO(Γ) phonon resonances, at around 38 THz (blue pump spectrum in Fig. 2A), a relatively small signal is found (blue curve in Fig. 2B). In marked contrast, we witness a response more than one order of magnitude stronger for resonant phonon excitation around 20 THz (red pump spectrum in Fig. 2A): an ultrafast single-exponential drop of the Faraday signal with a time constant as short as τfast = 1.6 ps (red curve in Fig. 2B). On much longer time scales, an additional exponential reduction of θ with a time constant of τslow = 90 ns is found (Fig. 2C). Recovery back to the initial state occurs over about 1 ms (fig. S2D).

Because the pump-probe signal grows linearly with the pump fluence (inset of Fig. 2C), excitation is dominated by one-photon absorption, whereas strong-field effects such as field or impact ionization are negligible (27). This notion is further supported by Fig. 2A, which demonstrates that the sample absorptance and transient Faraday rotation are found to depend on the pump frequency in a very similar way.

We emphasize that almost identical dynamics are observed when a pure YIG sample, instead of the BiGa:YIG film, is used (see fig. S2, B to D). Additional control experiments confirm that, as soon as the pump pulse has left the sample, the transient Faraday signal Δθ(t) reliably reflects the dynamics of the sublattice magnetizations Ma and Md with time-independent coefficients aa and ad in Eq. 1 (see Materials and Methods and fig. S3).

Figure 2B also shows the relative optical transmittance change of the BiGa:YIG sample for the case of resonant pumping. The signal starts with a peak-like feature whose width coincides with that of the terahertz pump pulse. Already for delays t > 1 ps, the signal is almost constant, changing by less than 10% in the subsequent evolution. This relaxation to a quasi-steady state is substantially faster than seen for the transient Faraday rotation, which still doubles its value at t = 1 ps in the following 5 ps (Fig. 2B). Furthermore, the transmittance signal is found to be independent of the sample magnetization (fig. S3C). These features suggest that the transmittance change predominantly monitors the redistribution dynamics of the pump-deposited energy in the crystal lattice.

Although the relative signal changes are small and still in the perturbative regime, our results in Fig. 2 show a proof of concept that resonant phonon excitation provides an ultrafast manipulation of magnetic order. It only involves the crystal lattice and electron spins, yet no electronic orbital degrees of freedom (Fig. 1A). The picosecond spin dynamics observed here (Fig. 2B) are unexpected because they are five orders of magnitude faster than the spin coherence lifetime (>0.1 μs) of YIG, which is known to be one of the longest among magnetically ordered materials (8, 15).

Temperature dependence

To characterize the transient states established on the fast (τfast = 1.6 ps) and slow (τslow = 90 ns) time scales, respectively, we increase the sample temperature T0 from 300 to 420 K. Simultaneously, we measure the equilibrium Faraday rotation θ0 as well as its pump-induced change Δθ at ultrashort (t = 10 ps; Fig. 2B) and long delays (t = 1 μs; Fig. 2C) after phonon excitation.

The resulting θ0 versus T0 has the typical shape (15) of a ferrimagnet’s static magnetization curve (Fig. 3A). Its slope ∂θ0/∂T0 steepens with rising T0 until the transition into the paramagnetic phase occurs at the Curie temperature TC = 398 K of BiGa:YIG. In contrast, the pump-induced Faraday signal at t = 1 μs (Fig. 3B) increases with T0 and reaches a maximum right below TC, reminiscent of the derivative of the static curve (Fig. 3A). Indeed, we find that Δθ(1 μs) versus T0 closely follows (∂θ0/∂T0T (Fig. 3B). Here, ΔT = 0.39 K is the increase in equilibrium temperature as calculated from the energy density deposited by the 1-μJ pump pulse (see Materials and Methods). The agreement of both curves shows that ~1 μs after pumping, the BiGa:YIG film is in full thermodynamic equilibrium characterized by temperature T0 + ΔT.

Fig. 3 Two regimes of spin-lattice equilibration.

(A) Equilibrium Faraday rotation θ0 = θ(−2 ps) versus ambient temperature T0 along with a fit to an analytical function (thin solid line). (B) Pump-induced change Δθ in the Faraday rotation at t = 10 ps (red symbols) and 1 μs (blue) after pump-pulse arrival. The black curve is the change (∂θ0/∂T0T in the Faraday rotation expected from the increase ΔT of the sample temperature due to heating by the pump pulse. θ0(T0) is taken from (A) (thin solid line), and ΔT = 0.39 K is calculated from the absorbed pump energy and the heat capacity of the excited volume.

Figure 3B reveals that the changes in the Faraday signal at 10 ps are systematically yet nonuniformly smaller than at 1 μs, in agreement with Fig. 2 (B and C). While for temperatures T0 < 380 K, Δθ(10 ps) amounts to roughly 80% of Δθ(1 μs), this ratio decreases strongly when T0 is increased. For example, slightly above 380 K (where both curves reach their maxima) and slightly below 400 K, Δθ(10 ps) has reduced to only 50 and 20% of Δθ(1 μs), respectively. Consequently, Δθ(10 ps) versus T0 in Fig. 3B cannot be made to agree with Δθ(1 μs) versus T0 by a simple global rescaling of the Δθ(10 ps) values. In particular, the shape of Δθ(1 μs) versus T0 agrees much better with a suitably scaled derivative ∂θ0/∂T0 of the static Faraday rotation (black curve in Fig. 3B) than Δθ(10 ps) versus T0. Therefore, at time t = 10 ps after pump excitation, the spin system is in a state that is significantly different from the equilibrium state found at t = 1 μs.

Analysis of spin couplings

To understand the microscopic mechanism driving the ultrafast change in magnetic order after resonant phonon excitation (Fig. 2B), we note that solids exhibit only three fundamental spin couplings. They can be understood as effective magnetic fields exerting torques on spins. In the following, we discuss all of them in terms of their capability to modify the sublattice magnetizations Ma and/or Md on ultrafast time scales.

First, spin-orbit (SO) coupling is often considered to dominate the ultrafast demagnetization of laser-excited ferromagnetic metals in the absence of transport (28). In our experiment, however, we find identical ultrafast dynamics for pure YIG and BiGa:YIG (see fig. S2B), although the Bi ions are known to increase SO coupling. This notion is manifested in the magneto-optic effects, which are enhanced by more than one order of magnitude (24). Therefore, SO coupling plays a negligible role in the ultrafast spin dynamics seen in Fig. 2B.

The same conclusion can be drawn for the second type of coupling, spin-spin magnetic dipole (SSMD) interaction, which has a comparable strength to SO coupling in YIG (15). We note that SO and SSMD coupling also account for magnetic anisotropy, which can undergo strong photo-induced changes in the garnet YIG:Co (29). Ultrafast laser-induced modification of the anisotropy field of a YIG:Co film was even shown to drive precessional magnetization switching (29). The reported time constant of this process (~20 ps) is, however, still more than one order of magnitude longer than the ultrafast dynamics (τfast = 1.6 ps) seen here. In addition, we do not observe any sign of coherent precession, which is typical of anisotropy changes in the perturbative regime. Consequently, neither SO nor SSMD coupling can explain the spin dynamics observed on the τfast scale here.

We finally consider the only remaining fundamental spin coupling, isotropic exchange, which is the strongest spin interaction in most magnets (2). In YIG, it is responsible for the ferrimagnetic order and the high magnon frequencies extending to more than 20 THz at room temperature (15, 30). Therefore, exchange interaction may well account for the picosecond spin dynamics observed here (Fig. 2B). Because it conserves the total spin angular momentum, the magnetization changes of a and d sublattices cancel each other. Despite ΔMa + ΔMd = 0, the resulting change in Faraday rotation is nonzero because the magneto-optic coefficients aa and ad of YIG differ significantly (see Eq. 1) (24). In summary, our analysis of the fundamental spin couplings shows that only isotropic exchange interaction is capable of explaining the picosecond spin dynamics seen in our experiment.


We suggest the following scenario of exchange-mediated spin-phonon coupling: The resonantly excited TO(Γ) phonons decay into other modes significantly faster than the time constant τfast = 1.6 ps of the measured ultrafast spin dynamics (26), thereby heating up the crystal lattice. The assumption of ultrafast phonon relaxation is consistent with our transient transmittance data (Fig. 2B), which indicate that the phonons reach an approximately stationary state within less than 1 ps. As the heat capacity of the crystal lattice of YIG is two orders of magnitude higher than that of the spin system (31), the lattice acts akin to a bath whose temperature is suddenly increased from T0 to T0 + ΔT.

As a consequence of the ultrafast lattice heating, the O2− ions, which are the lightest, will undergo additional random deflection Δu(t) (Fig. 4A). This perturbation modulates the superexchange (32) of adjacent a-Fe3+ and d-Fe3+ spins and the associated coupling constant Jad by (16, 33, 34)ΔJad(t)=(Jad/u)Δu(t)(2)We now put this model to test by conducting calculations of both dynamic and equilibrium states and by their comparison to our experimental observations on the ultrafast (τfast; Fig. 2B) and the slower time scale (τslow; Fig. 2C).

Fig. 4 Atomistic spin-dynamics simulations.

(A) Schematic of our model of ultrafast spin-phonon coupling. Thermal motion of the O2− ion modulates the superexchange constant Jad of the adjacent a-Fe3+ and d-Fe3+ spins, thereby enabling transfer of spin angular momentum between a- and d-spin sublattices. (B) Evolution of a- and d-spin sublattice magnetizations as obtained by atomistic spin-dynamics simulations. From 0 to 0.5 ps, thermal modulation of the exchange constant Jad is switched on (orange square). The variance of the Jad fluctuation is proportional to the difference ΔT between crystal lattice and spin temperature. To obtain agreement of the slope of ΔMa(t) with that found in the experiment directly after pump excitation (−0.1% ps−1; see Fig. 2B), where ΔT = 0.39 K (see Fig. 3B), an approximately three times smaller ∂Jad/∂u of ~10 Jad Å−1 than used here has to be chosen.

Dynamics on the τfast scale

To determine the ultrafast rate of change of Ma und Md, we perform atomistic spin-dynamics simulations based on ~106 coupled spin equations of motion (see Materials and Methods) (30). We include time-dependent exchange parameters through Eq. 2, where Δu(t) is assumed to be random with a variance given by the pump-induced temperature increase ΔT of the ionic lattice.

Simulation results are shown in Fig. 4B. At times t < 0, the magnetizations of both sublattices fluctuate around their constant means Ma0 and Md0. However, when fluctuations ΔJad(t) of the exchange parameter are switched on at t = 0, Ma decreases linearly with time until ΔJad(t) is switched off. We find the opposite behavior for ΔMd. The sum curve ΔMa + ΔMd has more than one order of magnitude smaller amplitude. It arises from the much weaker random fields of the bath due to spin couplings other than isotropic exchange (see Materials and Methods). Figure 4B demonstrates that thermal modulation of Jad induces demagnetization of the two spin sublattices by the same amount, and thus transfers energy into the spin system.

The slope of the simulated ΔMa(t)/Md0 (Fig. 4B) can be compared directly to the initial slope of the pump-probe signal Δθ(t)/θ0 (Fig. 2B), which amounts to approximately −0.1% ps−1. Agreement between experiment and theory is obtained by choosing ∂Jad/∂u~10 Jad Å−1 (see Materials and Methods). This value and the theory of superexchange (32) imply that the overlap integral of the ground-state wave functions of adjacent O2− and Fe3+ ions changes by 100% when their distance changes by ~0.4 Å. As the distance between the two ions is one order of magnitude smaller than the lattice constant a of YIG (15), we estimate ∂u/∂a~0.1 and obtain ∂Jad/∂a~(∂Jad/∂u)(∂u/∂a)~1 Jad Å−1. This result is in agreement with recent ab initio calculations (34) of the exchange coupling constants as a function of a, which yielded ∂Jad/∂a~0.7 Jad Å−1. Therefore, phonon-modulated exchange coupling successfully and consistently explains the rate (∂Δθ/∂t)/θ0 at which magnetic order of YIG is observed to reduce immediately after resonant phonon excitation.

According to our model, ∂Ma/∂t, ∂Md/∂t, and, thus, ∂Δθ/∂t (Eq. 1) are proportional to the difference between lattice and spin temperature (Eq. 10), in agreement with the linear pump fluence dependence of the transient Faraday rotation Δθ seen in the experiment (inset of Fig. 2C). Consequently, the measured exponential decay of Δθ to a constant value (Fig. 2B) indicates that the electron spins approach the lattice temperature T0 + ΔT with the time constant τfast = 1.6 ps.

Note that the subsequent slower dynamics evolve with a time constant of τslow = 90 ns. We conclude that for all pump-probe delays t ≪ τslow, the transient changes in Ma and Md arise solely from inter-sublattice exchange coupling. Because this interaction conserves the total spin angular momentum, the dynamics are under the constraint ΔMa + ΔMd = 0 throughout that time window.

Dynamics on the τslow scale

As indicated by our temperature-dependent measurements (Fig. 3), the garnet film approaches full equilibrium at temperature T0 + ΔT with the time constant τslow. The total magnetization Ma + Md of this state is lower than for t ≪ τslow, where the change ΔMa + ΔMd is still zero. Thus, transfer of angular momentum from the spin system to the lattice has certainly occurred for times t ~ τslow and above.

The only microscopic spin couplings capable of changing the total spin and, thus, magnetization Ma + Md are SO and SSMD coupling. If SO coupling were dominant, then one would expect a difference in the spin dynamics in pure YIG and in BiGa:YIG in which SO coupling is enhanced by partial Bi substitution. However, the signal-to-noise ratio of the Faraday signal seen on the τslow scale of pure YIG (fig. S2B) currently does not allow us to draw a conclusion about which of the two couplings is dominant.

Constrained spin state

As discussed above, the spin system is under the constraint of constant magnetization Ma + Md = Ma0 + Md0 for times t ≪ τslow, whereas it is in full, unconstrained equilibrium for t ≫ τslow (Fig. 2C). Assuming that for t ≫ τfast, the spins and crystal lattice are characterized by a single temperature T0 + ΔT, both the constrained and unconstrained spin state can be fully described by equilibrium statistical physics (see Materials and Methods).

Using a mean-field approximation, we obtain the change (∂Ml/∂T0T in the sublattice magnetizations versus T0 for heating with and without the constraint of conserved Ma + Md (Fig. 5A). The constrained ∂Ma/∂T0 of minority spin-sublattice a (green line in Fig. 5A) closely follows its unconstrained counterpart. In contrast, the constrained ∂Md/∂T0 of majority spin-sublattice d (blue line in Fig. 5A) is systematically smaller than the unconstrained ∂Md/∂T0 because sublattice d can only demagnetize as much as sublattice a (ΔMd = −ΔMa). In particular, the difference between constrained and unconstrained ∂Md/∂T0 increases considerably when the temperature T0 approaches the Curie point (Fig. 5A). We emphasize that this behavior is in excellent qualitative agreement with the measured Faraday signals Δθ(10 ps) and Δθ(1 μs) versus T0 (Fig. 3B).

Fig. 5 Constrained state and spin-phonon equilibration in YIG.

(A) Calculated change in sublattice magnetization Ma and Md per increase of temperature T0, without and with the constraint of constant spin angular momentum. The constrained and unconstrained ∂Md/∂T0 curves exhibit good qualitative agreement with the measured Faraday signals Δθ(10 ps) and Δθ(1 μs) versus T0 shown in Fig. 3B. (B) Schematic of spin-phonon equilibration in YIG. [1] Pump-excited TO(Γ) phonons [2] increase the population of other lattice modes. The increased thermal modulation of the a-d exchange by ΔJad(t) leads to [3a] transfer of angular momentum between a- and d-spin sublattices, accompanied by [3b] energy transfer from the phonon to the spin system on the time scale τfast = 1.6 ps. The resulting state is constrained by ΔMa + ΔMd = 0 and decays by [4] transfer of angular momentum and energy between crystal lattice and electron spins on the τslow = 90 ns scale. This process is mediated by SO and/or SSMD coupling.

It is interesting to note that in the equilibrium formalism used here, the constant total magnetization is reinforced by a Lagrange multiplier, which can be interpreted as a virtual homogeneous magnetic field (see Eq. 13). One could also term this field “spin pressure,” in analogy to the pressure that is required to keep a gas of particles in a bottle of constant volume while the gas temperature is increased. In our experiment, the spin pressure builds up on the time scale given by τfast and is released by angular momentum transfer to the crystal lattice on the time scale τslow.

Picture of spin-phonon equilibration

In summary, our dynamic and equilibrium calculations are fully consistent with the time scales, fluence dependence, temperature dependence, and magnitude of the transient Faraday rotation found in the experiment. This agreement leads us to the following picture of the flow of energy and angular momentum of phonon-pumped YIG (Fig. 5B): [1] the pump pulse excites zone-center TO(Γ) phonons that [2] decay within the pump duration, thereby increasing the crystal lattice temperature. The additional thermal modulation of a-d exchange induces [3a] transfer of angular momentum between a- and d-spin sublattices and implies [3b] energy transfer from the phonon to the spin system with the time constant τfast = 1.6 ps (Fig. 2B). The excess magnetization of this constrained state decays by [4] transfer of angular momentum and energy between crystal lattice and spins with the time constant τslow = 90 ns, resulting in full, unconstrained equilibrium (Fig. 2C). The angular momentum transfer is mediated by SO and/or SSMD interactions (22), which do not conserve Ma + Md.


We use a terahertz pump magneto-optic–probe experiment to investigate spin-lattice equilibration from the femtosecond to the microsecond time scale. Combined with a spin-coupling analysis and atomistic spin-dynamics simulations, our results reveal that the speed of spin-phonon relaxation in ferrimagnetic insulators critically depends on whether one refers to energy or angular momentum.

On a picosecond scale, spins and phonons reach a quasi-equilibrium through phonon-modulated exchange interaction. It induces redistribution of energy among phonons and spins, whereas angular momentum is only transferred between spins. On a much longer time scale in the nanosecond range, the additional transfer of angular momentum between spins and crystal lattice results in full equilibration. The quasi-equilibrium persists on the intermediate time scale and can be understood as a transient spin state with elevated temperature yet unchanged net magnetization. It could be considered as a realization of magnon populations with nonzero chemical potential, which were recently introduced in the theory of magnon transport by the spin Seebeck effect (19).

Transfer of angular momentum between spin sublattices is also a key process in magnetic switching of ferrimagnetic metallic alloys by optical femtosecond laser pulses (3538). This process is considered to follow the preferential ultrafast demagnetization of one of the two spin sublattices. Note that the flow of angular momentum can proceed through various mechanisms such as direct exchange torque (36) or nonlocal processes (38) including the transport of spin-polarized electrons. Our results reveal a new mechanism for manipulating the exchange interaction on ultrafast time scales and highlight the important role phonons can play in the direct exchange dynamics between spin sublattices.

In terms of applications, our results suggest that resonant phonon excitation is a new pathway to the preparation of spin states with increased temperature yet unchanged total magnetization. This route may be particularly interesting for the manipulation of antiferromagnetic order (2), where spin angular momentum is inherently conserved. Finally, the ultrafast spin-lattice coupling of YIG implies that the magnon temperature follows the phonon temperature with a delay on the order of only 1 ps, thereby shifting the cutoff frequency of the bulk spin Seebeck effect (7) to the terahertz range.



Our iron garnet films were grown by liquid-phase epitaxy on substrates of gadolinium gallium garnet Gd3Ga5O12 (GGG). To prevent pump absorption by the substrate, we transferred several iron-garnet films on diamond windows after mechanically removing the GGG. Test measurements confirmed that GGG pump absorption was irrelevant to the ultrafast dynamics.

We studied two types of garnet samples: pure yttrium iron garnet YIG and bismuth/gallium-substituted iron garnet BixY3−xGayFe5−yO12 (x = 1.53, y = 1.33, BiGa:YIG). They have, respectively, a Curie temperature of 545 and 398 K, and in-plane and out-of-plane magnetic anisotropy. Film thicknesses covered a range from 7 to 20 μm. The chemical composition of the garnet films was determined by x-ray fluorescence measurements. Substitution of Y by Bi in the BiGa:YIG sample enhanced the magneto-optic Faraday effect by one order of magnitude (24). Typical hysteresis loops obtained by measuring the static Faraday rotation of the probe beam are shown in fig. S1C.

Terahertz spectra of the sample absorptance (which equals one minus the reflected and transmitted power, both normalized to the incident power) were obtained by combined reflection and transmission measurements with a broadband terahertz time-domain spectrometer. Figure 2A shows that maximum absorptance occurs at 21 THz, which is 2 THz above the highest-frequency TO(Γ) phonon resonance (26). According to Hofmeister and Campbell (26), absorption in this frequency range solely arises from TO phonons, while crystal-field and charge-transfer gap transitions are located at much higher frequencies (24).

Experimental design

A terahertz pump pulse was used to resonantly excite long-wavelength TO phonons of the iron-garnet film under study. The instantaneous magnetic state of the sample was determined by a subsequent optical probe pulse measuring the magneto-optic Faraday effect. In this way, spin dynamics were monitored over a large range of pump-probe delays ranging from femtoseconds to milliseconds. Details of our setup are shown in fig. S1.

To generate intense terahertz pump pulses, we used an amplified Ti:sapphire laser system delivering pulses (energy of 15 mJ, center wavelength of 800 nm, duration of 40 fs, and repetition rate of 1 kHz) to drive two optical parametric amplifiers that, in turn, generated intense infrared pulses (pulse energy of 1.5 and 1.2 mJ, center wavelength of 1280 and 1400 nm, and duration of 50 and 50 fs, respectively). By difference-frequency mixing of the two infrared pulses in a nonlinear-optical GaSe crystal (thickness of 1 mm) (21), we obtained intense terahertz pulses (tunable from 15 to 60 THz, energy of typically 3 to 20 μJ, pulse duration around 250 fs) with a stable carrier-envelope phase. For sample excitation, the terahertz pulse was focused onto the sample surface to a spot with a diameter of 180 μm. Its transient electric field E(t) was measured by electro-optic sampling (see below). The pulse spectrum was obtained by Fourier transformation of E(t) and shown in Fig. 2A.

Low-noise probe pulses (8 fs, 800 nm, 0.5 nJ) were derived from the pulse train of the Ti:sapphire seed laser (repetition rate of frep = 80 MHz). They traversed the sample collinearly with the pump after a delay t. To measure the probe’s transient polarization rotation Δθ(t), we used a polarimeter consisting of a Wollaston-prism polarizer followed by two fast photodiodes. Measurement of the probe ellipticity Δη(t) was accomplished by putting a quarter-wave plate before the Wollaston prism. The resulting photodiode current was a train of temporally isolated electrical pulses (duration of 5 ns each, pulse-to-pulse distance of 1/frep = 12.5 ns), which was sampled using a computer-controlled fast analog-digital converter card (39). In this way, we were able to measure signals at delays of t = τ + j/frep (with integer j) between −100 ns and ~1 ms within a single pump shot and with a delay spacing of 1/frep = 12.5 ns. The ultrashort offset τ was set from −2 to 200 ps by a variable mechanical delay stage.

Pump-probe signal traces Δθ±(t) were taken for the sample magnetization saturated in opposite directions using an external magnetic field Bext of ±15 mT, respectively. To measure samples with in-plane or out-of-plane anisotropy, we oriented the magnetic field 45° with respect to the incident pump and probe beams. To vary the sample temperature, we mounted the sample on a resistive heater. The sample temperature was measured with two thermocouples and an imaging infrared thermometer, all of them yielding consistent temperature values. To measure the terahertz electric field, we substituted the sample by a GaSe crystal and measured the transient ellipticity Δη(t)∝E(t) induced by the electro-optic effect (21).

Optionally, our setup also permitted measurement of the transient isotropic transmittance of the sample. For this purpose, we removed all polarization-sensitive components such as the Wollaston prism and wave plate and detected the probe-pulse energy with one of two fast photodiodes of the polarimeter. Considering the small attenuation length (26) of ~1 μm of the resonant pump pulse into the iron-garnet sample, we estimate that the temporal blurring of the pump-probe signals due to pump-probe velocity mismatch is smaller than 10 fs and, therefore, negligible.

Signal analysis

Up to first order in magnetization, the Faraday rotation θ (and likewise the ellipticity η) of the probe’s outgoing polarization is a weighted average (24) of the magnetization of the two spin sublattices l = a, d, that isθ=alMld+b(3)Here, al are the local magneto-optic coupling coefficients (Verdet constants) of YIG, Ml = eprMl is the magnetization of sublattice l projected on the propagation direction unit vector epr of the probe, d is the sample thickness, and b is an offset arising from any magnetization-independent optical anisotropy. The angular brackets 〈.〉 denote averaging over the probed sample volume.

In equilibrium, the al are constant throughout the probed magnetic volume, 〈alMl〉 = alMl〉. The impact of the pump pulse makes all quantities in Eq. 3 time-dependent, resulting in relatively small changes θ(t) – θ0 of the signal θ0 = θ(−2 ps) measured 2 ps before arrival of the pump pulse. According to Eq. 3, the nonmagnetic offset b cancels by calculating the asymmetric partΔθ(t)=Δθ+Δθ2[al0ΔMl+Ml0Δal](4)The first term on the right-hand side scales with a weighted average of the sublattice magnetizations along epr, which is the quantity we are interested in. The second term, however, makes a contribution independent of the ΔMl and arises from a possible pump-induced modulation Δal of the Verdet constants. Figure S3 shows that the second term is negligible and that Δθ reflects the true magnetization dynamics, Δθ = ∑al0〈ΔMl〉, once the pump pulse has left the sample.

Excitation profile and heating

Because the absorption length of the terahertz pump pulse (26) (~1 μm) is smaller than the sample thickness d (~10 μm), pump excitation is inhomogeneous along the z axis normal to the film plane. Because the pump-induced signal Δθ(t) depends linearly on the absorbed pump fluence Fabs (inset of Fig. 2C), the pump-induced change in the local Faraday rotation Δθ(z,t) is proportional to the pump energy density wabs(z) absorbed in the vicinity of the plane at z. That is, one has Fabs=0ddzwabs(z) and Δθ(z,t) = C(t)wabs(z), where C(t) captures the temporal dynamics. As a consequence, the total pump-induced signalΔθ(t)=Δθ(z,t)d=0ddzΔθ(z,t)=C(t)Fabs(5)is independent of the shape of the absorption profile wabs (z). Therefore, without loss of generality, we can consider the absorption profile to be homogeneous with wabs = Fabs/d. Our argumentation is valid as long as transport of heat from the YIG film into the substrate is negligible. Figure S2D shows that this assumption is fulfilled for pump-probe delays up to at least 1 μs.

Assuming that the sample has reached (local) thermal equilibrium at t = 1 μs, we estimate the average temperature increase ΔT (1 μs) of the BiGa:YIG film by dividing Fabs by the pumped sample volume [d × (pump diameter)2 × π/4], by the mass density, and by the specific heat capacity (31) of YIG. We obtain a value of ΔT (1 μs) = 0.39 K for an incident pump energy of 1 μJ, which can be interpreted as the temperature increase of a homogeneously excited YIG film.

We checked the effects due to accumulative heating of the sample by reducing the repetition rate of the pump pulses from 1000 to 500 Hz by means of a mechanical chopper. The magnitude of pump-induced signal that persists after 1 ms was found to be smaller than 10% of the pump-induced signal at t ~ 1 μs. It corresponds to a static homogeneous temperature increase of the probing volume of only ~39 mK, which has a negligible influence on our results, including the temperature dependence seen in Fig. 3B.

YIG ab initio calculations

Our ab initio phonon calculations were performed using the finite-displacement method implemented in the open-source package phonopy (40). The electronic-structure calculation and equilibrium geometry search of YIG were performed within the density functional theory framework, using the Vienna Ab initio Simulation Package (41) and adopting the generalized gradient approximation (42) of the exchange-correlation functional. Our calculations and resulting phonon dispersion relations are detailed in the Supplementary Materials.

From the element-resolved vibrational density of states (fig. S4), we found that a continuum of modes from ~7 to 20 THz contributes to the displacement of the O2− ion with approximately frequency-independent weight. This result implies that the vibration of the O2− ion has a mean frequency of ΩO/2π ~ 14 THz and a correlation time of τO ~ 10 fs.

Dynamic spin model

Our consideration of all possible spin-coupling mechanisms (SO, SSMD, exchange; see main text) leads us to the view that the ultrafast regime of the spin dynamics of YIG after phonon excitation arises from isotropic exchange interaction. More precisely, the thermal vibrations of the ionic lattice modulate the a-d exchange interaction (Fig. 4A), thereby transferring spin angular momentum from the a- to d-spin sublattice and quenching each sublattice magnetization at the same rate.

To test this hypothesis, we calculated both the rate of spin angular momentum transfer between sublattices (based on numerical simulations) and the temperature dependence of the sublattice magnetizations (based on an analytical treatment). The starting point is the Heisenberg spin Hamiltonian of YIG (15)H^=12lljjJlljjS^ljS^ljgμBljBextS^lj(6)where S^lj is the spin angular momentum operator of the total electron spin of an Fe3+ ion (spin quantum number S = 5/2) located at site j in spin-sublattice l. Orbital contributions to the magnetic moments are negligibly small (25). In YIG, there are two different crystallographic sites l = a, d for the Fe3+ ions, which are, respectively, surrounded by six O2− ions (octahedral sites) and four O2− ions (tetrahedral sites). The exchange constants were determined as follows (15): If for a given l, l′, the j and j′ represent next neighbors, then Jlljj=Jll. Otherwise, Jlljj=0. Thus, the exchange part of the Hamiltonian is fully determined by the three constants Jaa, Jdd, and Jad = Jda.

Following Eq. 2, the phonon-induced modulation of the exchange interaction is modeled as variationΔJad(t)=uJadΔu(t)(7)of the exchange coupling constant between l = a and l′ = d spin sublattices. Here, ∂u = ∂/∂u, and Δu is the pump-heat-induced deflection of the O2− ion mediating the superexchange of adjacent a- and d-type Fe3+ ions, as depicted in Fig. 4A. This assumption is justified because Jad predominantly involves the a-Fe3+ and d-Fe3+ ions and the O2− ion in between. Because the O2− ion is, by far, the lightest ion in the YIG unit cell, its motion is expected to modulate the a-d superexchange most strongly.

Atomistic spin-dynamics simulations

To calculate the rate of spin angular momentum transfer between sublattices, we conducted numerical simulations based on Eqs. 6 and 7. For the implementation, we assumed an ensemble of classical spins and replaced each spin operator S^i by si(S+1)S. Here, si is a vector of length one, and the compound index (lj) is summarized in one index i. Likewise, the Hamilton operator Ĥ turns into the Hamilton function H. For the exchange constants Jaa, Jdd, and Jad, we used the values 0.33, 1.15, and −3.43 meV, respectively (15).

By applying Langevin theory, a stochastic Landau-Lifshitz-Gilbert–type equation is derived for each spin (30, 43, 44)tsi=|γi|(1+λi2)μi[si×Beff,i+λisi×(si×Beff,i)](8)where γi is the gyromagnetic constant of a spin with magnetic moment μisi. Damping due to coupling to nonspin degrees of freedom (the so-called bath) was quantified by λi = 2 × 10–5. The effective field Beff,i acting on si has three contributions: a deterministic part due to the exchange fields Bexch,i0=H0/si of the adjacent spins for temporally constant Jij0 and two stochastic components ξi(t) and ΔBexch,i(t) that arise from coupling of the spins to the bath. The first field ξi(t) transfers spin angular momentum and energy between the spin system and the bath having temperature T0. According to the fluctuation-dissipation theorem, in equilibrium, the fluctuations ξi(t) are balanced by a friction force, the second, Gilbert-type term in Eq. 8, thereby driving the spin system into an equilibrium state with temperature T0 (30, 43, 44). The time constant of this relaxation process is set by the damping parameters λi and found to be irrelevant on the picosecond scale considered here.

The second stochastic field ΔBexch,ijΔJijsj arises from the phonon-induced modulation of the exchange interaction (Eq. 7). It is distinctly different from ξi(t) as it causes energy transfer into the spin system but leaves the total spin angular momentum unchanged. Because the noise field at temperature T0 is already accounted for by ξi(t), the variance of ΔJij(t) scales with the temperature increase ΔT of the crystal lattice induced by the pump pulse. More precisely, owing to the equipartition theorem, the pump-induced deflection Δu of the O2− ion obeys mOΩO2Δu2=mOΩO2Δu2=kBΔT, where ΩO/2π~14 THz is the mean O2−-ion vibration frequency (see above).

To determine the initial rate of phonon-driven angular momentum transfer between the spin sublattices, we assume phonons and spins to be thermalized at temperatures T0 + ΔT and T0, respectively (see main text). Following our ab initio results (see above), Δu(t) is assumed to have a correlation function with a width of τO ~ 10 fs, which is accordingly modeled asΔu(t)Δu(t)=kBΔTmOΩO2τOδ(tt)(9)The equation of motion (Eq. 8) was integrated numerically for N = 643 = 2.6 × 105 unit cells at equilibrium temperature T0 = 300 K. At each time step t, we calculated the two sublattice magnetizations Ml=jslj(S+1)S/N with l = a, d. In our numerical implementation, the pump-induced thermal noise of the exchange constants (Eqs. 7 and 9) can be switched on or off over arbitrary time intervals and with variable strength (∂uJad)2ΔT.

The results of our atomistic simulations are shown in Fig. 4B and fig. S5. We found the demagnetization rate of each sublattice scales according to|tMl/Md0|ΔT=(uJad)2f(10)with constant f = 2 × 1016 m2 J–2 K–1 ps–1. On the other hand, our experiment shows that the Faraday signal Δθ and thus |∂tMl/Ml0| scale with the pump fluence (see inset of Fig. 2C) and thus the initial ΔT as well. In other words, we have|tMl/Ml0|ΔT|tΔθ/θ0|ΔT=g(11)where g = 1.9 × 10−3 K−1 ps−1 was determined from the initial slope of Δθ/θ0 ~ ΔMl/Ml0 (Fig. 2B) seen in our experimental results. Comparison of Eqs. 10 and 11 yielded the estimate |uJad|g/f10JadÅ−1.

According to the theory of superexchange (32), Jad is proportional to the fourth power of the overlap integral of the Fe3+ and O2− ions. Therefore, upon changing the distance of the two ions, the overlap integral undergoes relative changes on the order of 2.5 Å−1.

Mean-field sublattice magnetization

Our previous results strongly indicate that spins and phonons reach quasi-equilibrium for times t ≫ τfast after phonon excitation. As exchange interaction leaves the total spin angular momentum of the electrons unchanged, this thermal state is constrained by the boundary condition of conserved total spin angular momentum. Such constraints are ubiquitous in thermostatics (such as for a gas contained in a closed bottle) and can be treated by equilibrium statistical physics (45).

Our goal is to determine the magnetization of the d- and a-spin sublattices of YIG in thermal equilibrium for two cases: the unconstrained (UC) situation without boundary condition and the constrained (C) situation with the boundary condition of constant total spin angular momentumljS^lj= constant(12)In both cases, we need to calculate the statistical operator D^ of the system, from which the expectation value A^= Tr(D^A^) of an observable  (such as a sublattice magnetization) can be derived. In thermal equilibrium and for both the UC and C cases, the statistical operator is given by (45)D^=1Zexp(βH^βpljS^lj)(13)where β = 1/kBT, with kB being the Boltzmann constant. The partition function Z is determined by the normalization condition Tr D^=1. In the UC case, we set p = 0 in Eq. 13, thereby resulting in the standard canonical statistical operator (45). In the C case, the constraint of conserved spin angular momentum (Eq. 12) has to be fulfilled by proper choice of the Lagrange multiplier p in Eq. 13. We note that p can be interpreted as a virtual homogeneous magnetic field that is adjusted to reinforce the boundary condition of Eq. 12. It is analogous to the chemical potential that keeps the number of particles in a given macroscopic system constant.

The general Heisenberg spin Hamiltonian Ĥ of Eq. 6 is too complex for an analytical calculation of D^. Consequently, we apply the mean-field approximation (32), S^ljS^ljS^ljS^lj+S^ljS^ljS^ljS^lj, omit the number S^ljS^lj, and allow for alignment exclusively along the z-axis unit vector ez, S^lj=S^l,zjez=S^ljez, as set by the external magnetic field Bext = Bextez. In addition, we assume the spin order homogeneously on each sublattice l, that is, S^lj=S^l0 [sublattice approximation (32)]. By evaluating the Jlljj of YIG (see above), we finally obtain the mean-field (MF) HamiltonianH^MF=12llJllNllS^l0jS^ljgμBljBextS^lj(14)Here, Jll is defined as above, and Nll denotes the number of the next neighbors a site in sublattice l has in sublattice l′. For YIG, one has Naa = 8, Ndd = 4, Nad = 6, and Nda = 4 (15). We can now write the statistical operator of Eq. 13 in a very compact form asD^MF=1Zexp(βlFljS^lj)(15)where the effective fieldFl=12lJllNllS^l0+gμBBext+p(16)captures the mean field (first and second terms) and the scalar virtual field p = pez imposed by the constraint of Eq. 12. By taking advantage of the fact that all spin operators in Eq. 15 commutate with each other, we obtain the expectation valueS^l0=1βFl Tr exp(βFlS^l0)(17)for any Fe3+ spin of sublattice l. Evaluation of this equation leads to an implicit relationship for S^l0S^l0=SBS(βSFl)(18)where S = 5/2, and BS is the Brillouin function (32). The magnetization Ml of sublattice l is proportional to vlS^l0, where vl is the number of l-Fe3+ ions per unit cell with vd = 12 and va = 8.

For the UC case (p = 0), Eq. 18 is a system of two coupled equations and solved numerically, thereby yielding the unconstrained S^l0UC for each temperature T0. For the C case, the boundary condition of Eq. 12, rewritten asνdS^d0C+νaS^a0C= constant(19)adds a third equation to Eq. 18, which determines p. Note that in our experiment, we start from an UC equilibrium given by the S^l0UC at temperature T0. Subsequently, the pump pulse increases T0 by a small amount ΔT while keeping the total spin angular momentum constant. Because we are only interested in small relative changes of the S^l0, we linearized Eqs. 18 and 19 in terms of ΔT and solved analytically for the small changes ΔS^l0C.

To be consistent with our atomistic spin-dynamics model (see above), we chose the ratio of the coupling parameters Jad, Jaa, and Jdd to be identical to the values given in (15), but rescaled by a factor of 0.36 to fit the experimentally observed critical temperature of Tc = 398 K. For comparison with our experiment, we finally determined the Faraday rotation resulting from the calculated sublattice magnetizations Ml (see Eq. 1). The results of this procedure are shown and discussed in Fig. 5A and fig. S6.


Supplementary material for this article is available at

Supplementary Methods

1. Details of YIG ab initio calculations.

Fig. S1. Experimental details.

Fig. S2. Supporting pump-probe data.

Fig. S3. Nonmagnetic signal contributions.

Fig. S4. Ab initio phonon calculations.

Fig. S5. Supporting data of atomistic spin-dynamics simulations.

Fig. S6. Calculated equilibrium sublattice magnetizations.

References (4648)

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We thank G. E. W. Bauer and S. M. Rezende for stimulating discussions. Funding: T.K. acknowledges funding through the grants European Research Council H2020 CoG TERAMAG (grant no. 681917) and DFG SFB/TRR 227 “Spin Dynamics” (project A05). I.R. acknowledges funding through German Federal Ministry of Education and Research (BMBF) grant 05K16BCA “Femto-THz-X.” P.M. and P.M.O. acknowledge financial support from the Swedish Research Council (VR), the K. and A. Wallenberg Foundation (grant no. 2015.0060), and the Swedish National Infrastructure for Computing. M.G. thanks for funding through BMBF grant no. 05K10KEB. R.V.P. acknowledges financial support from the Russian Scientific Foundation (grant no. 16-12-10456). M.W. acknowledges funding by the Max Planck Society. Author contributions: I.R. conceived the original experimental idea. S.F.M., I.R., and T.K. designed the experiment. S.F.M. built the setup, performed measurements, and analyzed data. I.R., A.M.K., and R.V.P. prepared the samples, which were characterized by S.F.M., I.R., A.M.K., A.P., and M.G. The theoretical model was developed by T.K. and S.F.M., with contributions from J.B., P.M., A.P., and P.M.O. Atomistic spin-dynamics simulations were conducted by J.B. Equilibrium magnetization curves were calculated by S.F.M. and T.K. Ab initio calculations of electronic and vibrational properties were performed by P.M. and P.M.O. The manuscript was written by T.K., J.B., and S.F.M., with contributions from A.P., P.M.O., P.M., M.W., and I.R. All authors contributed to the discussions of the results and commented on the manuscript. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article