Abstract A proposed paradigm for out-of-equilibrium quantum systems is that an analog of quantum phase transitions exists between parameter regimes of qualitatively distinct time-dependent behavior. Here, we present evidence of such a transition between dynamical phases in a cold-atom quantum simulator of the collective Heisenberg model. Our simulator encodes spin in the hyperfine states of ultracold fermionic potassium. Atoms are pinned in a network of single-particle modes, whose spatial extent emulates the long-range interactions of traditional quantum magnets. We find that below a critical interaction strength, magnetization of an initially polarized fermionic gas decays quickly, while above the transition point, the magnetization becomes long-lived because of an energy gap that protects against dephasing by the inhomogeneous axial field. Our quantum simulation reveals a nonequilibrium transition predicted to exist but not yet directly observed in quenched s-wave superconductors.

INTRODUCTION The challenge faced in understanding out-of-equilibrium systems is that the powerful formalism of statistical physics, which has allowed a classification of quantum phases of matter based on simple principles such as minimization of free energy, does not apply. A diverse range of nonequilibrium phenomena has been observed, including synchronization (1–5), self-organization (6–9), quantum chaos (10, 11), Loschmidt echo singularities (12), and time crystals (13, 14). A proposed organizing principle is that transitions, reminiscent of those found between thermodynamic ground states, can also be found between dynamical phases (15–19). In general terms, a nonequilibrium phase transition is characterized by the existence of a critical point separating phases with distinct dynamical properties in many-body systems. The analog of thermodynamic order parameters is found in long-time average observables, which have a nonanalytic dependence on system parameters. In driven open systems, energy and particle number are not necessarily conserved, and nonequilibrium transitions are typically signaled by different steady states that occur upon varying system parameters such as pump or loss rates (20–22), independently of initial conditions. In closed systems, dynamics are often initiated by quenching control parameters, with qualitatively distinct behaviors observed below, above, or at a critical point (15–18) that can depend on the initial state of the system. The label “dynamical phase transition” has been applied not only to the boundary between two dynamical phases but also to the nonanalytic behavior in real-time dynamics of the return probability amplitude (12, 23), which does not require an order parameter to be defined (24). The phenomenon under investigation in our work is the former case, which we will refer to as a transition between dynamical phases (TDP) to avoid confusion. The theoretical study of these transitions has encompassed a broad range of platforms including collective spin models (25–27), nonequilibrium phases of superconductors (28–30), interacting fermions and bosons on the lattice (15–17, 31, 32), and quantum field theories (18, 33); however, experimental investigations have so far been restricted to self-trapping transitions in bosonic systems (34–38) and the transverse-field Ising model realized with trapped-ion chains (19). Here, we report the observation of a transition between two dynamical phases of a quantum degenerate Fermi gas. The sample under investigation consists of neutral potassium atoms (40K) confined in a harmonic optical trap and cooled to nanokelvin temperatures. The controllable interactions of this closed quantum system enable a broad search for nonequilibrium phenomena that arise from the interplay of atomic contact interactions, quantum statistics, and motion. Using collective magnetization as an order parameter, system dynamics are observed directly and compared to theoretical models at various levels of approximation. We understand and analyze our system through a mapping of the single-particle eigenstates of the harmonic trap onto a lattice in mode space, as depicted in Fig. 1. By tuning the interaction strength to suppress collisions that would change the occupancy of the modes, the atoms become pinned on the conceptual lattice, enabling the description of our system with a spin model (39–41); here H ^ / ℏ = ∑ i h i s ^ i Z − ∑ i , j J ij s ^ i ⋅ s ^ j (1)where s ^ i = 1 / 2 { σ ^ i X , σ ^ i Y , σ ^ i Z } are spin-1/2 operators acting on the ith atom and X, Y, and Z denote orientations in Bloch space. This is the collective Heisenberg model, a canonical model for magnetism (42) in which the nonlocal spin-spin couplings J ij compete with an inhomogeneous axial field, h i (39–41). Similar treatments of fermionic systems using a spin model have been used successfully in optical lattice clocks (39, 43–46) in the microkelvin regime. However, in those experiments, undesirable inelastic collisions limited the number, N ≲ 50, and prevented the study of transitions between well-defined dynamical phases. The stability of low-lying hyperfine states and control over interactions in our experiment allow us to explore many-body dynamics in macroscopic ultracold samples of N ≃ 3 × 104 alkali atoms at nanokelvin temperatures and test the spin model in this new regime. Fig. 1 Simulation of the collective Heisenberg model with a local inhomogeneous axial field using weakly interacting fermions in a mode space lattice. Each site in mode space (left) has an occupancy of 0 or 1 and experiences a local field h i that causes single-spin precession at a rate that depends on the mode (top right). Atoms experience long-range spin-exchange interactions J ij (middle right). Mode-changing collisions would move atoms between sites in mode space (bottom right) but are not included in the spin model. We find that below a critical interaction strength, the magnetization of an initially polarized gas quickly decays, while above this critical point, the magnetization becomes long-lived and is protected by an energy gap against inhomogeneous field-induced dephasing. These observed dynamical phases are a manifestation of an emergent property in a many-body dynamical quantum system. The transition would be absent for small particle number since the emergence of a sufficiently strong gap from weak two-body collisions relies upon a collective nonlinearity. We then implement a many-body echo sequence (4, 47, 48), which is a direct test of reversibility. This allows us to identify the boundaries of the parameter regime in which the complex far-from-equilibrium dynamics of interacting fermions is quantitatively described by the collective Heisenberg model. The successful mapping allows us to implement quantum simulations of the nonequilibrium phases predicted to exist in quenched s-wave superconductors by Richardson-Gaudin models (28–30) but not yet directly observed, given the need for ultrafast probes (49).

DISCUSSION In summary, we demonstrate the existence of a TDP in a neutral Fermi gas in a regime of reversible dynamics near the zero crossing of a Feshbach resonance. We outline a direct connection to nonequilibrium phases in the Richardson-Gaudin models for superconductivity, thereby extending experimental observations of TDP’s beyond prior manifestations in Josephson- and Ising-type systems (19, 34–38). Moreover, the excellent agreement between spin-model calculations and a two-axis exploration of the dynamical phase diagram with >104 spins provide experimental evidence of the scaling behavior and universal character of TDPs. The collective nature of the dynamics observed here could protect many-body states in other systems of interest for applied quantum technologies. For example, gap protection would increase the coherence time in optical lattice clocks operated in the quantum degenerate regime (61). Furthermore, using the effective time-reversal capability demonstrated here, together with technical improvements to magnetic field stability and homogeneity, our system could provide a fruitful platform to measure out-of-time order correlations and scrambling of quantum information (48) or test spectroscopic protocols that use time reversal to relax the detection resolution required for spectroscopy beyond the standard quantum limit (62).

MATERIALS AND METHODS The neutral atomic sample was prepared using laser cooling, magnetic and optical trapping, and evaporative and sympathetic cooling in an atom-chip apparatus described previously (63, 64). The hyperfine states used to encode spin information are the F = 9/2, m F = −9/2 (for ∣↓〉), and m F = −7/2 (for ∣↑〉) states. At the beginning of a simulation sequence, the fermionic ensemble had a typical temperature of several hundred nanokelvin, with data taken in the range of T ~ 0.3 to 0.5 E F /k B , where E F = ( 6 N ) 1 / 3 ℏ ω ¯ is the Fermi energy and ω ¯ is the geometric mean trap frequency. There was no optical lattice: Confinement was produced by a crossed-beam 1064-nm optical dipole trap, and the map in Fig. 1 is purely conceptual. The magnetic field and its gradients were controlled using a combination of microfabricated wires on the atom chip ≈200 μm from the atoms and macroscopic coils external to the vacuum system. The field was calibrated through rf spectroscopy of the ∣↓〉 to ∣↑〉 transition. During a typical experimental run, drifts are ≲1 μT. Taking advantage of the symmetry of the many-body dynamics with a, we determined the magnetic field corresponding to a = 0 by observing the magnetization after t = 100 ms at various magnetic fields and then fitting the resulting profile to find the center. The mean result of 10 such trials taken across 6 months results in B zc = 20.907(1) mT, which is an improvement in uncertainty by an order of magnitude over previous measurements (see the Supplementary Materials for data and literature comparison). Magnetic field gradients, which lead to periodic oscillation of the spin clouds (see the Supplementary Materials), were measured by displacing the trap center and repeating rf spectroscopy. In the optimal configuration, gradients are {13(1),12(1),2(5)} μT/m. The differential displacements Δ i resulting from these gradients are small compared to the harmonic oscillator length a i 2 = ℏ / m ω i , as dimensionless displacements Δ i /a i = {2.1,0.4,0.1} × 10−2. In this Δ i ≪ a i regime, the more general XXZ spin model (41) reduces to the Heisenberg model used here. The effective axial field h i in the Heisenberg Hamiltonian is the potential energy differential between the ∣↓〉 and ∣↑〉 states, as sampled by the occupied motional eigenstates. Since one can always subtract a constant potential term, dynamics depend only on the inhomogeneity in h i , quantified through its root mean square spread h ∼ . There are both magnetic and optical contributions to h i . The leading-order magnetic field contribution is curvature in real space. Direct spectroscopic data can bound this to ≲500 μT/m2, which is compatible with the 200 μT/m2 best-fit curvature in the model. The optical contribution to h ∼ is tuned using the polarization of the trapping light. For far-detuned light, the differential fractional vector light shift (VLS) experienced by the atoms is approximately Pg F (δ FS /3δ), where P is the polarization of the light, g F is the g-factor of the hyperfine sublevels, δ is the frequency detuning, and δ FS is the fine-structure splitting of the electronic excited state. For a full range of polarization P ∈ [ −1,1] of the optical dipole trapping beam propagating along the z axis, the resulting VLS should modify the trap frequency along x by roughly ±0.06 Hz or ±0.015%. This control over h ∼ allows the vertical exploration in the dynamical phase diagram shown in Fig. 2. The magnetization 2S/N was determined using repeated measurements and a maximum-likelihood estimator, as follows. Each experimental image measures N ↑ , N ↓ at the end of a Ramsey sequence. The fraction of spin-↑, f = N ↑ /N, relates to the collective spin as f = 1/2 + ( S Y /N)cosθ + ( S X /N)sinθ, where θ is the phase lag of the second π/2 pulse. However, the typical evolution times (100 ms) exceed the clock coherence time (1.5 ms) so that the accumulated phase ϕ randomizes the orientation of the transverse spin, S Y = Scosϕ, S X = Ssinϕ. This yields f = 1/2 + (S/N) cos (θ′), where θ′ = θ − φ is a random phase. To reconstruct the offset O and amplitude A of a Ramsey fringe F = O + Acos (θ′) from a set of fractions f i acquired in several experimental runs with the same conditions, we assumed a probability distribution p ( f ; A , O ) = ∫ 0 π d θ ′ 2 π σ π exp [ − ( A cosθ ′ + O − f ) 2 2 σ 2 ] which is a convolution of the noise-free probability distribution and Gaussian noise with width σ, calibrated to be σ = 0.01. We constructed a log-likelihood function for a set of n fraction measurements, ℓ ( A , O ; { f i } ) = 1 n ∑ i = 1 n log ( p ( f i ; A , O ) ) , from which we numerically computed the maximum-likelihood amplitude A and fringe offset O, as well as confidence intervals. Where O falls outside the range of 0.50 ± 0.05, we discarded the data; otherwise, the peak-to-peak amplitude 2A was taken as the best estimate of 2S/N.

SUPPLEMENTARY MATERIALS Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/5/8/eaax1568/DC1 Fig. S1. Order parameters predicted by Lax vector analysis in a 1D system. Fig. S2. Approximate form of the Lax vector in a 3D system. Fig. S3. Effective mean-field potential. Fig. S4. Magnetization dynamics for noninteracting particles. Fig. S5. Determination of the Feshbach zero crossing. Fig. S6. Spin-echo amplitude near the Feshbach zero crossing. Fig. S7. Fits to time series. Fig. S8. Equilibrium scattering rate versus temperature. Fig. S9. Nonequilibrium scattering rate. Table S1. Determination of 40K Feshbach resonance parameters. References (65–75)

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 A. Koller and C. Luciuk for early work on this project and V. Gurarie, B. Lev, J. Thompson, M. Foster, and D. Stamper-Kurn for discussions. Funding: This work is supported by NSERC; the Air Force Office of Scientific Research grants FA9550-13-1-0063 and FA9550-18-1-0319 and its Multidisciplinary University Research Initiative grant (MURI); the Army Research Office grants W911NF-15-1-0603, W911NF-16-1-0576, and W911NF-19-1-0210; the Defense Advanced Research Projects Agency (DARPA); the NSF grant PHY1820885; JILA-NSF grant PFC-173400; and the National Institute of Standards and Technology. Author contributions: The work was conceived by A.M.R., J.H.T., and S.T. Experiments were performed by S.S., B.A.O., H.S., K.G.J., and S.T. Data were analyzed by S.S., P.H., and B.A.O. Theoretical models and simulation were done by P.H., J.M., J.H.T., and A.M.R. All authors contributed to manuscript preparation. 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.The datasets generated and analyzed during the current study are available from the corresponding authors upon reasonable request.