## Simulating magnetism out of equilibrium

A major goal of quantum simulation is to help us understand problems that are difficult to describe analytically or solve with conventional computers. This goal has been very challenging to reach experimentally, requiring, for example, extremely low temperatures to determine the effects of quantum magnetism in equilibrium. Brown *et al.* studied a nonequilibrium system of ultracold ^{87}Rb atoms in a two-dimensional optical lattice. They monitored the atoms' dynamics after a sudden change in the lattice parameters and were able to reach a regime where the magnetic interactions dominated the dynamics.

*Science*, this issue p. 540

## Abstract

The interplay of magnetic exchange interactions and tunneling underlies many complex quantum phenomena observed in real materials. We study nonequilibrium magnetization dynamics in an extended two-dimensional (2D) system by loading effective spin-1/2 bosons into a spin-dependent optical lattice and use the lattice to separately control the resonance conditions for tunneling and superexchange. After preparing a nonequilibrium antiferromagnetically ordered state, we observe relaxation dynamics governed by two well-separated rates, which scale with the parameters associated with superexchange and tunneling. With tunneling off-resonantly suppressed, we observe superexchange-dominated dynamics over two orders of magnitude in magnetic coupling strength. Our experiment will serve as a benchmark for future theoretical work as the detailed dynamics of this 2D, strongly correlated, and far-from-equilibrium quantum system remain out of reach of current computational techniques.

The interplay of spin and motion underlies some of the most intriguing and poorly understood behaviors in many-body quantum systems (*1*). A well-known example is the onset of superconductivity in cuprate compounds when mobile holes are introduced into an otherwise insulating two-dimensional (2D) quantum magnet (*2*). Understanding this behavior is particularly challenging because the insulating state has strong quantum spin fluctuations, and the introduction of holes impedes the use of numerically exact techniques for quantum-correlated systems in more than one dimension. Ultracold atoms in optical lattices realize tunable, idealized models of such behavior and can naturally operate in a regime where the quantum motion (tunneling) of particles and magnetic interactions (superexchange) appear simultaneously (*3*, *4*).

For ultracold atoms in equilibrium, the extremely small energy scale associated with superexchange interactions makes the observation of magnetism challenging, and short-range antiferromagnetic correlations resulting from superexchange have only recently been observed (*5*, *6*). Out-of-equilibrium superexchange-dominated dynamics has been demonstrated in isolated pairs of atoms (*7*), in 1D systems with single-atom spin impurities (*8*), and recently in the decay of spin-density waves (*9*). However, the perturbative origin of superexchange in these systems implies that it must be weak compared with tunneling, and thus the manifestation of superexchange requires extremely low motional entropy. Dipolar gases (*10*) and ultracold polar molecules (*11*) in lattices provide a promising route toward achieving large (nonperturbative) magnetic interactions (*12*), but technical limitations in these systems currently complicate the simultaneous observation of motional and spin-exchange effects.

Here, we study the magnetization dynamics of pseudospin-1/2 bosons in a 2D optical lattice after a global quench from an initially antiferromagnetically ordered state (*13*). The dynamics we observe are governed by a bosonic t-J model (*14*, *15*) and occur in a parameter regime where quantitatively accurate calculations are not currently possible. Using a checkerboard optical lattice, we continuously tune the magnetization dynamics from a tunneling-dominated regime into a regime where superexchange is dominant, even at relatively high motional entropies. We observe decay of the antiferromagnetic order with clearly identifiable time scales that correspond to the underlying Hamiltonian parameters for tunneling and superexchange. In contrast to the oscillatory dynamics observed in isolated pairs of atoms (*7*), we observe exponential decay of the magnetization. Some damping of the dynamics can be attributed to the extended, many-body nature of the experiment, but the degree of the expected damping is difficult to quantify because a full theoretical description of nonequilibrium 2D systems is lacking. This experiment bridges the gap between studies of nonequilibrium behavior in systems with exclusively motional (*16*–*20*) or spin (*21*, *22*) degrees of freedom, demonstrating the requisite control to explore the intriguing intermediate territory in which they coexist. In addition, the techniques that we demonstrate lay the groundwork for adiabatic preparation of low-entropy spin states relevant for studies of equilibrium quantum magnetism (*23*, *24*).

Our experiment uses two hyperfine states (denoted by and ) of ultracold ^{87}Rb atoms trapped in a dynamically controlled, 2D checkerboard optical lattice (*25*) composed of two sublattices, and (Fig. 1A). For most experimental conditions presented here, our system is well described by a Bose-Hubbard model (*3*) characterized by a nearest-neighbor tunneling energy , and an onsite interaction energy . In addition, we use the lattice to apply an energy offset between the and sublattices, consisting of a spin-independent part and a spin-dependent part acting as a staggered magnetic field, (*26*). All of these parameters can be dynamically controlled, which we exploit to prepare initial states with 2D antiferromagnetic order and to observe the resulting dynamics after a quench to different values of , , , and .

At unit filling, for and , double occupation at each site is allowed only virtually, and the Bose-Hubbard model can be mapped onto a ferromagnetic Heisenberg model (*4*, *27*) with a nearest-neighbor magnetic interaction strength that is second-order in the tunneling energy. In the presence of hole impurities, first-order tunneling (with the much larger energy scale ) must be included, which substantially modifies the dynamics even at low hole concentrations (*9*). The offset provides the flexibility to tune the relative importance of first-order tunneling and second-order superexchange processes. For example, below unit filling, if and , the Bose-Hubbard model can be mapped onto a bosonic t-J model with a staggered energy offset.

The local spin operators are defined as , where a_{iσ} () annihilates (creates) a hardcore boson of spin on site , and is a vector of Pauli matrices. The notation indicates that the sum over and is restricted to nearest neighbors, and indicates that the sum is restricted to sites , such that are both nearest neighbors of . The superexchange energy (Fig. 1B) can be either ferromagnetic () or antiferromagnetic () (*7*). The last term in Eq. 1 describes hole-mediated exchange, where an atom on site interacts via superexchange with an atom on site , while simultaneously hopping to site (Fig. 1B). Here, , where −(+) applies when *j* is in the lower (higher) energy sublattice. In writing Eq. 1, we have ignored second-order processes (*28*, *29*) that conserve sublattice magnetization (*30*). When , first-order tunneling is resonant and dominates the magnetization dynamics except at extremely low hole densities. For , however, first-order tunneling is effectively suppressed, in which case superexchange dominates the magnetization dynamics. (The frequently ignored term should be included for a quantitative description at nonzero hole density.) Superexchange is resonant when but is suppressed when . The values of , , and are determined from an experimentally calibrated model of the lattice (*30*). Inhomogeneity in the system—arising, for example, from trap curvature—primarily enters the model via inhomogeneities in the parameters and .

The experiments begin with ^{87}Rb atoms loaded into a square 3D optical lattice with one atom per site (*30*), initially spin-polarized in the state . We use the hyperfine states and to represent the pseudospin-1/2 system. ( and are the quantum numbers labeling the total angular momentum and its projection along the quantization axis, respectively.) The 3D lattice is composed of a vertical lattice along , which confines the atoms to an array of independent 2D planes, along with the dynamic 2D checkerboard lattice in the *x-y* plane. The vertical lattice depth is typically , held constant throughout the experiment, and the 2D lattice depth is initially with no staggered offset, (the recoil energy , kHz, where is the mass of ^{87}RbRb and nm). The atoms occupy roughly 13 to 15 2D planes, with the central plane containing 800 to 1100 atoms. The ratio of surface lattice sites to total lattice sites of the trapped cloud is 15% and sets a zero-temperature lower bound for the number of sites with neighboring holes. Based on spectroscopic measurements and assuming a thermal distribution (*30*), we estimate the hole density at the center of the cloud to be about 8%.

To measure the spin population independently on each sublattice, we map the four spin/sublattice states , , and onto four distinct Zeeman states and determine their populations with absorption imaging (Fig. 2A). To perform the experiment, we start with a spin-polarized configuration (Fig. 2B) and construct an initial state with staggered magnetization by applying the addressing offset and transferring the -site atoms to (Fig. 2C). After returning to zero, we initiate dynamics by quenching to a given configuration with lattice depth and offsets and (Fig. 2D). The ramp time for the quench of 200 μs was chosen to be fast with respect to subsequent dynamics but slow enough to avoid band excitation. After a variable hold time , we freeze the dynamics by raising to and read out the populations , from which we determine the staggered magnetization and the sublattice population difference .

(2)The exchange terms in Eq. 1 conserve , whereas the first-order tunneling does not, allowing for population transport between sublattices. We also monitor the total spin imbalance and the population to quantify unwanted spin-changing processes that drive the atoms out of the pseudospin-1/2 manifold containing and . The measured time for depopulation of the pseudospin-1/2 states is greater than 6 s, and the atom number lifetime in the lattice is greater than 3 s.

For the lattice parameters studied in this paper, the magnetization dynamics is well described by exponential decay, with decay time scales ranging between 0.5 ms and 500 ms. Example decay curves are shown in Fig. 3A for a lattice depth and different offsets . For some and , the exponential decay clearly occurs on two well-separated time scales (Fig. 3A, inset): a fast time scale, , which dominates the behavior in shallow lattices when or , and a slow time scale, , which dominates the behavior in deep lattices with larger offset, .

To investigate the faster time scale, we measure and at a fixed decay time for different , as shown in Fig. 3B for 15(0.5) . The 5-ms decay time (vertical line in Fig. 3A, inset) was chosen so that nearly all of the fast decay but little of the slow decay occurred. The fast magnetization decay reveals resonant features at and , with the decay rate at twice as fast as at . (Near the resonance at , the condition is not satisfied and the system must be described with the full Bose-Hubbard Hamiltonian, rather than Eq. 1.) In addition, shows sublattice transport from to sites at , indicative of resonant first-order tunneling. At the demagnetizing sublattice transport and are balanced. We theoretically estimate the expected width of the resonance in to be Hz (*30*), which is narrower than the 530(60) Hz width that we observe experimentally, suggesting inhomogeneous broadening. We note, however, that the observed broadening is beyond what is expected from the measured trap curvature and is inconsistent with estimates of light-shift inhomogeneity from spectroscopic measurements (*30*). A residual could account for the width.

The measured decay times for a range of lattice depths are plotted against the calculated tunneling time in Fig. 3C, showing a decay rate linear in , with . This slope is comparable to a simple theoretical estimate taking only resonant tunneling into account, which predicts , with the lattice coordination number. Given the relatively large hole density near the surface of the cloud, the agreement with a noninteracting estimate is not surprising, although we would expect interactions to reduce the decay rate.

To investigate the slow dynamics, we measure the magnetization decay time for , where first-order tunneling is negligible and superexchange should dominate the dynamics. To determine the dependence of on the spin-dependent staggered offset , we measure the remaining staggered magnetization and population difference after a fixed wait time , as shown in Fig. 4A for a lattice depth 9.4(3) and offset kHz. As expected for superexchange dynamics, the magnetization decay is resonant in . The full width at half maximum of the Lorentzian fit to the resonance is 126(14) Hz, considerably narrower than the observed tunneling resonances shown in Fig. 3B, and is most likely dominated by inhomogeneous broadening (second-order exchange processes are sensitive to inhomogeneity in both and ). Figure 4A also shows that there is negligible sublattice transport associated with the demagnetization resonance. We note that at these values of , the ground state of the system would have all atoms on the lower sublattice, and the conservation of indicates that the spin dynamics occurs within a metastable manifold with respect to population.

Figure 4B shows the measured resonant decay times versus calculated for different and , with and chosen to be larger than but considerably less than the next excited band. The decay time collapses to a single curve over two orders of magnitude in . The solid gray line through the data is a fit to , where is needed to capture the apparent saturation of at large ; the fitted values are and s. A quantitative calculation of the decay rate in 2D, including the effects of holes, is extremely challenging. However, a short-time perturbative calculation supports the experimental observation that scales with the single energy scale , which is not a priori expected given the existence of other () energy scales in the Hamiltonian (*30*). Surprisingly, the perturbative calculations suggest that the decay rate is nearly independent of hole density, which can be attributed to the approximate cancellation of two competing effects of holes: they decrease the rate of superexchange dynamics but simultaneously open new demagnetization channels through the hole-mediated spin-exchange term in Eq. 1. We therefore compare the data to an analytical estimate, strictly valid at unit filling, of (black line in Fig. 4B). The empirically determined time scale (beyond which deviations from linear scaling become apparent in Fig. 4B) is shorter than the measured times for depopulation of the pseudospin manifold or loss of atoms. For times substantially larger than , interplane tunneling may also play a role in the short-time magnetization dynamics. This background decay rate may also be related to the nonzero relaxation processes observed outside the resonance in Fig. 4A. The mechanism for this off-resonant decay is not clear, but because the initial and final states differ in energy by substantially more than , it must arise from energy-nonconserving processes such as noise-assisted relaxation or doublon production (*20*). Corrections to caused by excited-band virtual processes (*8*), which we estimate to be on the order of 10 to 20% at the largest and smallest shown in Fig. 4B, may partially explain the observed saturation.

The scaling and resonant behavior of the fast and slow relaxation processes clearly reveal their origin as first-order tunneling and superexchange, respectively. For , our experiment realizes an unusual situation where tunneling, which is only active within a given sublattice, is comparable in strength to the superexchange coupling. This feature—which is crucial to our ability to observe superexchange-dominated dynamics—may have interesting implications for the equilibration of a doped antiferromagnetic state, because it determines the extent to which entropy (initially introduced in the motional degrees of freedom) is shared by the spin degrees of freedom. For smaller but nonzero , the ability to observe both tunneling and superexchange, often simultaneously and at experimentally accessible entropies, opens exciting opportunities to explore the nonequilibrium interplay of spin exchange and motion. Understanding the detailed dynamics of this strongly correlated, 2D quantum system is a formidable challenge, which may require the development of new theoretical techniques.

## Supplementary Materials

## References and Notes

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵Strictly speaking, the validity of Eq. 1 also requires and to avoid narrower second-order resonances that are associated with double and triple occupancy, respectively.
- ↵
- ↵
- ↵Materials and methods are available as supplementary materials on
*Science*Online. - ↵
- ↵
- ↵
**Acknowledgments:**This work was partially supported by the Army Research Office’s atomtronics Multidisciplinary University Research Initiative and NIST. M.F.F. and E.A.G. acknowledge support from the National Research Council Research Associateship program. We thank B. Grinkemeyer for his contributions to the data taking effort, E. Tiesinga and S. Paul for discussions about tight-binding models, and A. V. Gorshkov and S. Sugawa for a critical reading of the manuscript.