Probing quench dynamics across a quantum phase transition into a 2D Ising antiferromagnet

Simulating the real-time evolution of quantum spin systems far out of equilibrium poses a major theoretical challenge, especially in more than one dimension. We experimentally explore the dynamics of a two-dimensional Ising spin system with transverse and longitudinal fields as we quench it across a quantum phase transition from a paramagnet to an antiferromagnet. We realize the system with a near unit-occupancy atomic array of over 200 atoms obtained by loading a spin-polarized band insulator of fermionic lithium into an optical lattice and induce short-range interactions by direct excitation to a low-lying Rydberg state. Using site-resolved microscopy, we probe the correlations in the system after a sudden quench from the paramagnetic state and compare our measurements to exact calculations in the regime where it is possible. We achieve many-body states with longer-range antiferromagnetic correlations by implementing a near-adiabatic quench and study the buildup of correlations as we cross the quantum phase transition at different rates.

In this work, we explore the dynamics of a two dimensional quantum Ising model using a nearly defect-free array of neutral atoms which are coupled with laser light to a low-lying Rydberg state in an optical lattice [38]. The spin coupling in the model arises due to a van der Waals interaction between atoms in the Rydberg state. If one atom is in a Rydberg state, the excitation of another atom to a Rydberg state is strongly suppressed within a blockade radius R b [39][40][41][42][43]. This is because the interaction between the Rydberg atoms within this radius is much larger than the laser coupling strength. Previous experiments in 2D arrays have studied the regime R b a l , where a l is the lattice spacing [25,28]. In this regime, the Rydberg blockade makes it difficult to access many-body states with a large Rydberg fraction. This significantly reduces the size of the relevant Hilbert space of the system from the maximum possible size of 2 N , where N is the number of sites, rendering simulation of the quantum dynamics feasible for the experimentally realized system sizes. Here we focus on the regime R b ∼ a l , where there is no such reduction of the Hilbert space size. This regime allows us to study quench dynamics across a quantum phase transition between a paramagnet and an antiferromagnet with broken Z 2 symmetry. Recently, Rydberg atoms in rearrangeable optical tweezers have explored this regime in 1D chains [29] and rings [28].
We realize a quantum Ising spin system with an array of 6 Li atoms in an optical lattice with near unitoccupancy. The lattice is deep enough to suppress tunneling over the timescale of the experiments. We prepare all the atoms in the same hyperfine ground state |↓ . Interactions are introduced by globally coupling the atoms with a single laser field to a Rydberg state |↑ . The van der Waals interaction between atoms in the Rydberg state is isotropic and takes the form V ij = C 6 /|r i − r j | 6 . The Hamiltonian of the system is given by Here S α i are the spin 1/2 operators for the ith lattice site and α = x, y, z. The first two terms of this Hamiltonian describe transverse and longitudinal magnetic fields that couple to the pseudospin. The Rabi frequency Ω that drives a transition between the ground and the Rydberg state for an isolated atom determines the transverse field, while the detuning ∆ of the laser frequency from atomic resonance determines the longitudinal field.
Vij 2 can be taken as a site independent detuning in a large system as ours. We work with an attractively interacting (V ij < 0) Rydberg state [44]. In the absence of the fields, the Hamiltonian's most excited state is a classical antiferromagnet, which is the ground state of the HamiltonianH = −H. For R b = (C 6 /Ω) 1/6 a l , the ground state phase diagram ofH in Ω/∆ parameter space contains multiple Rydberg crystalline phases with arXiv:1711.00887v1 [quant-ph] 2 Nov 2017 different Rydberg atom fractions [29,[45][46][47]. However for R b ∼ a l , the regime we study in this experiment,H can be approximated by a nearest-neighbor Ising Hamiltonian with coupling J = C 6 /a 6 l . A phase diagram for this model is shown in Fig. 1(a) and has only one ordered phase, the antiferromagnet [48,49]. The initial state in the experiment is the paramagnetic ground state ofH for positive detuning ∆ J Ω. In this work, we quench the system from this initial state to the antiferromagnet with varying degrees of adiabaticity and study the ensuing dynamics of the spin correlations.
We prepare nearly defect-free 2D arrays of atoms by taking advantage of Pauli blocking in a highly spinimbalanced degenerate Fermi gas loaded into a square optical lattice (for details see [37]). The spin mixture consists of the first |1 and third |3 lowest hyperfine ground states of 6 Li, with |1 as the majority. The minority atoms, needed to thermalize the gas while loading into the lattice, are subsequently removed with a pulse of resonant light. We focus our analysis on an annular region with outer (inner) radius of 9 (4) sites where the average occupancy of the remaining state |1 atoms, measured from repeated preparations of the system, is maximal and corresponds to 95.7(4) %.
We couple the state |1 atoms to the 23P (m l = 0, m s = −1/2, m I = 1) Rydberg state using singlephoton excitation with an ultraviolet (UV) laser at 230 nm. The experiments are performed at a bias magnetic field of 595 G pointing orthogonal to the 2D layer, allowing us to address a single |m l , m s , m I Rydberg state. Up to 60 mW of UV light is available from a frequency-quadrupled diode-laser system. The light is π-polarized and focused to a waist of 70 µm. The intensity and the frequency of the light can be changed rapidly to control the time dependence of the transverse and longitudinal fields in the Hamiltonian [44].
The atoms are located at the focus of a high resolution objective that can resolve individual sites of the optical lattice ( Fig. 1(b)). The Rydberg dynamics take place in a lattice of depth 55E R , where E R = (π ) 2 /2ma 2 l is the recoil energy and a l = 1064 nm/ √ 2. We image the distribution of ground state atoms after removing Rydberg atoms with an efficiency of 90(3) % by increasing the lattice depth to 2500E R , leading to rapid photoionization or expulsion of the anti-trapped Rydberg atoms ( Fig. 1(c)). We obtain site-resolved fluorescence images of the ground state atoms by collecting ∼ 1000 photons per atom scattered from laser beams in a Raman cooling configuration [37].
We calibrate the transverse and longitudinal fields of the Hamiltonian using sparse clouds where the average spacing between atoms is much larger than R b . The location of the Rydberg resonance (∆ = 0) is determined by finding the laser frequency which maximizes atom loss during a long exposure to the UV light, since atoms in the Rydberg state experience an anti-trapping optical potential. The Rabi frequency Ω is determined by measuring single atom Rabi oscillations, and we attain a maximum Rabi frequency Ω max = h × 5.4(1) MHz [44]. Ω varies 4.9(3) % over the region of interest due to the Gaussian intensity profile of the UV beam. The C 6 coefficient, which determines the strength of the van der Waals interaction, depends strongly on the principal quantum number. We obtain a theoretical C 6 /h = −1.92(6) MHz µm 6 = −10.6(3) MHz a 6 l for the 23P, m l = 0 state at an offset field of 595 G [44]. The angular dependence of the interaction potential in the P -state is unimportant in our experiments since the magnetic quantization axis is orthogonal to the plane of the lattice, leading to an isotropic interaction for atoms in the 2D plane. For these parameters, Ω max , J h/τ , where τ ∼ 20 µs is the lifetime of the Rydberg state [50], leading to negligible decay over the relevant timescales.
We first study dynamics in the Ising system after a sudden quench, where the transverse field is switched on quickly compared to h/Ω.
The system is initially in a product state, with all spins in |↓ , and we image the atoms after an evolution time T . From the images, we extract the spin correlators We observe a change in the sign of the nearest neighbor correlations as the detuning ∆ is varied.
For such short times, the correlations remain shortrange and the dynamics can be calculated. We implement a dynamical version [44,51] of the numerical linked cluster expansion (NLCE) to compare with our results [16,52]. The dynamics is computed on clusters of increasing size (the "order" of the expansion) and the re-sults are expected to converge if the correlation length is smaller than the cluster size. The 11th order NLCE results for the on-site and nearest-neighbor correlations are fit to the measured correlations after the quench with two free parameters: the van der Waals interaction coefficient C 6 and a scaling factor α corresponding to the Rydberg imaging efficiency. The NLCE dynamics calculations take into account interactions up to next-nearest neighbors and experimental imperfections including the finite rise and fall time of Ω and 2.8 % anisotropy of the lattice spacing [37], which translates to an 18 % anisotropy of the interactions on the nearest neighbor sites. We also compare the data to exact diagonalization results on a 4 × 4 lattice. From these fits, we obtain an experimental l and a scaling factor α = 0.89(1), which agrees with the expected detection efficiency [44]. The fitted value of C 6 is about 40% lower than the theoretically calculated C 6 , which has possible systematic errors due to uncertainties in the matrix elements in lithium, in particular at high magnetic fields, and finite wavefunction size of the atoms on the lattice sites [44].
To go beyond the regime where the dynamics can be calculated, we perform a longer quench with ΩT /h = 3π/2. The extracted correlators are shown in Fig. 2(e)-(h). In this case, even the next-nearest neighbor correlations exhibit a zero crossing as a function of detuning, showing that the system is building up longer range correlations. The different NLCE orders already stop converging at much earlier times. Thus, this data presents a challenge to state-of-the-art numerical methods for calculating dynamics.  To prepare many-body states with longer antiferromagnetic correlations, we investigate a more adiabatic quench scheme [46,53,54], illustrated in Fig. 3(a). We start from the same initial state but use a soft switch on and off of the Rabi frequency and a linear ramp of the detuning from ∆ i = 3.3J to a varying ∆ f . During the detuning ramp, the Rabi frequency is fixed at Ω 0 = 0.9(1)J. We explore a variety of detuning ramp rates∆ ranging from 1.6J 2 /h to 8.9J 2 /h. For each∆ we measure correlations at different times in the ramp. Fig. 3(b) shows the buildup of nearest neighbor and next-nearest neighbor antiferromagnetic correlations as the quantum phase transition is crossed with different speeds. The buildup of antiferromagnetic correlations starts approximately at the time the detuning ramp crosses ∆ = 0. For all quench speeds studied, we observe that the correlations reach a maximum at ∆/J ∼ −2, as expected from the phase diagram in Fig. 1(a). The peak value of the correlations initially increases as the quench speed is reduced as one might expect, but then decreases for slower ramps. This is likely due to decoherence playing a more important role in the slower quenches. We studied single spin coherence in our system with a Ramsey sequence in a sparse cloud and did not observe any significant decoherence over the relevant timescales [44]. However, strong attractive interactions between atoms in the Rydberg state can lead to many-body decoherence effects resulting from atomic motion. These effects are particularly strong in our system compared to previous optical tweezer experiments [28,29] due to the light mass of lithium and the relatively small lattice spacing. Displacement of the atoms, estimated in [44], leads to intrinsic decoherence by changing the spin coupling J and to lesser extent the detuning ∆ because of differential Stark shifts between the ground and Rydberg states induced by the optical lattice potential.
Near the end of the ramps, where significant antiferro- magnetic correlations have built up, we find that we can fit the decay of the correlations with distance to an exponential ( Fig. 4(a)), even though the spin system may not have reached thermal equilibrium. The fitted correlation lengths range from ξ = 0.74(6)a l to ξ = 1.9(2)a l depending on∆. Another way to characterize the short range antiferromagnets created by these slow quenches is by extracting the probabilities for observing a particular spin configuration in a sub-system. In Fig. 4(b), we show the probability of observing different spin configurations in 3 × 3 sub-systems, not correcting for detection fidelities.
The two antiferromagnetic states are the most probable states near the end of the ramp, with an enhancement of a factor of 16(2) over a uniform distribution in the Hilbert space.
In conclusion, we studied quench dynamics in a 2D Ising model realized with ultracold atoms coupled to a Rydberg state in an optical lattice. The use of a light fermionic atom, 6 Li, allows us to use Pauli blocking in a relatively large spacing lattice to create 2D atomic arrays with high-filling (∼ 96%), comparable to what is achieved in atom-by-atom assembler experiments [55,56]. Combining the large spacing with the use of a low-lying Rydberg state, we reached the strong correlation regime with R b ∼ a l and prepared states exhibiting strong shortrange antiferromagnetic correlations. We found good agreement of our data with state-of-the-art numerics for short-time quench dynamics. Additionally, we studied the dynamics in a regime outside the reach of exact theoretical techniques, providing test data for approximate techniques to calculate dynamics. Our new ultracold 6 Li Rydberg platform opens many interesting directions for future work. Rydberg excitation in a Fermi gas may allow the exploration of impurity dynamics in the presence of Pauli blocking effects [57,58]. Another possible direction is the use of Rydberg dressing techniques to realize a dipolar Fermi gas, taking advantage of the fast tunneling of lithium in an optical lattice to go beyond the frozen gas regime [27,[59][60][61].

SUPPLEMENTAL INFORMATION
The quantum gas microscopy apparatus used for preparing the atomic arrays used in this work and siteresolved imaging of these arrays is described in detail in the supplements of ref. [37].

Laser system for coupling to Rydberg states
For single-photon excitation to the Rydberg state, we use a deep UV laser system at 230 nm based on frequency-quadrupling light from a diode laser source (Laser & Electro-Optic Solutions). We start with 1.8 W of 920 nm light from a tapered amplifier seeded by an external grating diode laser. Next, a frequency-doubling cavity generates 750 mW at 460 nm. Finally, a second doubling cavity yields ≈50 mW of 230 nm light. A constant flow of oxygen was required to maintain intensity stability and increase the lifetime of the second harmonic generation crystal in the second cavity. The frequency of the 920 nm laser is stabilized to a Stable Laser Systems ultralow expansion glass (ULE) cavity with a measured finesse of ≈ 4000 and linewidth of 370 kHz. We use a fiber-EOM before the cavity to generate sidebands at two frequencies for Pound-Drever-Hall locking with offset frequency control of the laser with respect to the reference cavity [62]. To obtain fast amplitude and frequency control of the UV light, we used an AOM (IntraAction ASM-1501LA61) driven by a fast VCO (MiniCircuits ZX95-200A+).
For initial identification of the Rydberg lines, we performed a "V-scheme" spectroscopy [63] with the 230 nm beam and a diode laser driving the D 2 transition of lithium at 671 nm. The locking point on the cavity was recorded for a range of Rydberg states including the 23P state used in this experiment (Fig. S1). With the 920 nm laser locked, the final UV output had an intensity stability better than 10 %. The linewidth at 230 nm was measured to be ≈100 kHz.

Single atom Rabi oscillations and coherence
To measure the strength of the laser coupling to the Rydberg state, we measure Rabi oscillations in a sparse cloud where the interactions between the Rydberg atoms are negligible. A typical single atom Rabi oscillation is shown in Fig. S2(a). The decay of the Rabi oscillation is mainly due to shot-to-shot fluctuations of the laser intensity. In addition we measure the coherence of the atoms in a sparse cloud using a Ramsey echo sequence: π/2 -τ -π -τ -π/2 pulse, where τ is a delay time. The ground state fraction is measured at the end of the sequence (Fig. S2(b)). Even for τ = 0, corresponding to a 2π pulse, the measured ground state fraction is reduced to ∼ 0.8, because of laser intensity fluctuations. However, for a total delay 2τ = 1 µs, corresponding to 6h/J, we do not observe any decay of the ground state fraction, indicating that there is no loss of coherence.

Movement of the atoms during the quenches
In our calculations, we ignore the movement of atoms during the quenches. In reality, the atoms can move because of strong van der Waals forces between atoms in Rydberg states, photon recoil kicks and anti-trapping forces on atoms in the Rydberg state. From the experimentally fitted value for C 6 , we estimate Rydberg atoms accelerate towards each other with a(r) = 6C6 r 7 m where a ≈ 3.2 × 10 6 m/s 2 for neighboring Rydberg atoms and a( √ 2a lat ) ≈ 2.8 × 10 5 m/s 2 for next-neighbor Rydberg atoms. This leads to a displacement of 0.09a l (0.01a l ) for a typical quench time of 200 ns for nearest neighbor (next-nearest neighbor) Rydberg atoms. We note that the excitation of two Rydberg atoms on neighboring sites is largely suppressed due to the blockade. Recoil and anti-trapping forces lead to much smaller displacements of the atoms.

Theoretical calculation of C6
We calculate the value of C 6 for the state 23P, m l = 0, m s = 1/2 at 595 G with two different techniques. First we use a perturbation theory calculation in the m l basis which yields C 6 = −1.915 MHz µm 6 . As a second approach we use a pair state exact diagonalization in the m j basis at 595 G using [64] and perform a basis transformation to the m l -basis. A fit to the potential curve then yields C 6 . Depending on the inner cutoff r 0 for the fit we obtain values between C 6 = −1.974 MHz µm 6 for r 0 = 0.7 µm and C 6 = −1.864 MHz µm 6 at r 0 = 1 µm which leads us to the error estimate in the main text.

NLCE dynamics
The numerical linked cluster expansion (NLCE) algorithm used to calculate the dynamics in this work is an extension of the NLCE technique for thermodynamic quantities, reviewed in refs [16,52]. Our NLCE calculations take into account next-nearest neighbor (diagonal) interactions, lattice anisotropy, and the finite time for turning on and off the Rabi frequencies during "sudden" quenches. We discuss some of the major modifications that are made to the algorithm, assuming the reader's familiarity with the standard algorithm.
For systems with the symmetry of the square lattice and only nearest neighbor interactions, each embedding of a graph on the lattice is dependent only on the topology of the graph, which allows for a significant reduction in the number of clusters that need to be diagonalized. Taking into account next-nearest neighbor interactions and lattice anisotropy means that this is no longer true, as two topologically identical graphs may have different graph Hamiltonians. Thus, we break down topologically identical graphs further into classes of graphs with the same Hamiltonian up to graph symmetries, each of which we only need to solve once.
For each of these, we then perform a time evolution starting from the initial state using sparse representations of the Hamiltonian. The initial and final ramp is simulated by breaking down the overall ramp time into five time steps, and applying the time evolution between these steps with the appropriate time-averaged Hamiltonian. Then, the appropriate correlation functions can be extracted from the final state. The subgraph subtraction then proceeds as usual, except that each embedding should be treated independently, as the contribution of a graph to a particular correlator depends on its embedding in the lattice. Finally, we perform an Euler resummation starting from the 3rd order to reduce odd-even order fluctuation. Note that the graphical expansion done here is site-based, rather than link-based.
Finally, we have checked the effect of including beyondnext-nearest neighbor interactions, and the number of time steps, and found that our results are well converged with respect to them.

Fitting sudden quench dynamics to NLCE results
The correlator dynamics are computed using NLCE for a grid of ∆ and C 6 values at 9th order in the expansion, taking into account the independently calibrated Rabi frequency Ω = h × 4.05(2) MHz. The correlators C(0, 0), C(1, 0), and C(0, 1) are simultaneously fit to the results using two fit parameters: C 6 and a scaling factor α corresponding to a Rydberg atom detection efficiency. The scaling factor reduces the nearest neighbor correlators C(1, 0) and C(0, 1) as α 2 . Since C(0, 0) = n 2 − n 2 and n 2 = n since n is either 0 or 1, we obtain C(0, 0) = n − n 2 . This leads to a corrected correlator C * (0, 0) = α n − α 2 n 2 .