Shaking the lattice uncovers universality
Most of our knowledge of quantum phase transitions (QPTs)—which occur as a result of quantum, rather than thermal, fluctuations—comes from experiments performed in equilibrium conditions. Less is known about the dynamics of a system going through a QPT, which have been hypothesized to depend on a single time and length scale. Clark et al. confirmed this hypothesis in a gas of cesium atoms in an optical lattice, which was shaken progressively faster to drive the gas through a QPT.
Science, this issue p. 606
Abstract
The dynamics of many-body systems spanning condensed matter, cosmology, and beyond are hypothesized to be universal when the systems cross continuous phase transitions. The universal dynamics are expected to satisfy a scaling symmetry of space and time with the crossing rate, inspired by the Kibble-Zurek mechanism. We test this symmetry based on Bose condensates in a shaken optical lattice. Shaking the lattice drives condensates across an effectively ferromagnetic quantum phase transition. After crossing the critical point, the condensates manifest delayed growth of spin fluctuations and develop antiferromagnetic spatial correlations resulting from the sub-Poisson distribution of the spacing between topological defects. The fluctuations and correlations are invariant in scaled space-time coordinates, in support of the scaling symmetry of quantum critical dynamics.
Critical phenomena near a continuous phase transition reveal fascinating connections between seemingly disparate systems that can be described via the same universal principles. Such systems can be found in the contexts of superfluid helium (1), liquid crystals (2), biological cell membranes (3), the early universe (4), and cold atoms (5, 6). An important universal prediction is the power-law scaling of the topological defect density with the rate of crossing a critical point, as first discussed by T. Kibble in cosmology (4) and extended by W. Zurek in the context of condensed matter (1). Their theory, known as the Kibble-Zurek mechanism, has been the subject of intense experimental study that has largely supported the scaling laws (7). Recent theoretical works further propose the so-called universality hypothesis, according to which the collective dynamics across a critical point should be invariant in the space and time coordinates that scale with the Kibble-Zurek power law (8–10).
Atomic quantum gases provide a clean, well-characterized, and controlled platform for studying critical dynamics (6, 11, 12). They have enabled experiments on the formation of topological defects across the Bose-Einstein condensation transition (13–16) as well as critical dynamics across quantum phase transitions (17–23). Recent experiments using cold atoms in shaken optical lattices (24–26) have provided a vehicle for exploring phase transitions in spin models (27–29).
Here we study the critical dynamics of Bose condensates in a shaken optical lattice near an effectively ferromagnetic quantum phase transition. The transition occurs when we ramp the shaking amplitude across a critical value, causing the atomic population to bifurcate into two pseudo-spinor ground states (28). We measure the growth of spin fluctuations and the spatial spin correlations for ramping rates varied over two orders of magnitude. Beyond the critical point, we observe delayed development of ferromagnetic spin domains with long-range antiferromagnetic correlations due to the bunching of the domain sizes, which is not expected in a thermal distribution of ferromagnets. The times and lengths characterizing the critical dynamics agree with the scaling predicted by the Kibble-Zurek mechanism. The measured fluctuations and correlations collapse onto single curves in scaled space and time coordinates, supporting the universality hypothesis.
Our experiments use elongated three-dimensional (3D) Bose-Einstein condensates (BECs) of cesium atoms. We optically confine the condensates with trap frequencies of 
Hz, where the long (
) and short (
) axes are oriented at 
with respect to the 
and 
coordinates (Fig. 1A). The tight confinement along the vertical 
axis suppresses nontrivial dynamics in that direction (see the discussion on the dynamics in the 
direction below), which is also the optical axis of our imaging system. We adiabatically load the condensates into a 1D optical lattice (11) along the 
axis with a lattice spacing of 
nm and a depth of 8.86 
, where 
kHz is the recoil energy and 
is Planck’s constant.
(A) A BEC of cesium atoms (spheres) in a 1D optical lattice (pink surface) shaking with peak-to-peak amplitude 
can form ferromagnetic domains (blue and red regions). The elliptical harmonic confinement has principal axes rotated 45° from the lattice. (B) The transition occurs when the ground band evolves from quadratic for 
[paramagnetic (PM) phase], through quartic at the quantum critical point 
, to a double well for 
[ferromagnetic (FM) phase] with two minima at 
(28). (C) Kibble-Zurek picture. Evolution of the condensate crossing the phase transition becomes diabatic in the frozen regime (cyan) when the time 
remaining to reach the critical point is less than the relaxation time. Faster ramps cause freezing farther from the critical point, limiting the system to smaller domains. Sample domain images are shown for slow, medium, and fast ramps.
To induce the ferromagnetic quantum phase transition, we modulate the phase of the lattice beam to periodically translate the lattice potential by 
, where 
is the shaking amplitude and the modulation frequency 
is tuned to mix the ground and first excited lattice bands (fig. S1) (28, 30). The hybridized single-particle ground band energy 
can be modeled for small quasimomentum 
by
(1)where 
is the atomic mass, and the coefficients of its quadratic (
) and quartic (
) terms depend on the shaking amplitude (Fig. 1B). For shaking amplitudes below the critical value, the coefficient 
is positive and the BEC occupies the lone ground state at momentum 
. The quantum phase transition occurs when the quadratic term crosses zero at 
, where 
and 
. At this point, the speed of sound for superfluid excitations, formally studied in (31), drops to zero along 
but remains nearly constant along 
and 
. Even stronger shaking converts the ground band into a double well with 
, yielding two degenerate ground states with 
. Repulsively interacting bosons with this double-well ground band are effectively ferromagnetic, having two degenerate many-body ground states with all atoms either pseudo-spin up (
) or down (
) (28). Notably, transitioning to one of these two ground states requires the system to spontaneously break the symmetry of its Hamiltonian. Describing the dynamics across the critical point presents a major challenge because of the divergence of the correlation length of quasimomentum and the relaxation time (critical slowing).
The Kibble-Zurek mechanism provides a powerful insight into quantum critical dynamics. According to this theory, when the time remaining to reach the critical point inevitably becomes shorter than the relaxation time, the system becomes effectively frozen (Fig. 1C). The system only unfreezes at a delay time 
after passing the critical point, when relaxation becomes faster than the ramp. At this time, topological defects become visible, and the typical distance 
between neighboring defects is proportional to the equilibrium correlation length. The Kibble-Zurek mechanism predicts that 
and 
depend on the quench rate 
as
(2)
(3)where 
and 
are the equilibrium dynamical and correlation length exponents given by the universality class of the phase transition. Although the details of this picture may not apply to every phase transition, the general scaling arguments are very robust, and similar predictions hold for a variety of quench types across the transition (12) and for phase transitions that break either continuous or discrete symmetries (7).
For slow ramps, 
and 
diverge and become separated from other scales in the system, making them the dominant scales for characterizing the collective critical dynamics (8–10). This idea motivates the universality hypothesis, which can be expressed as
(4)indicating that the critical dynamics of any collective observable 
obeys the scaling symmetry and can be described by a universal function 
of the scaled coordinates 
and 
. The only effect of the quench rate is to modify the length and time scales.
We test the scaling symmetry of time by monitoring the emergence of quasimomentum fluctuations at different quench rates. Here, fluctuations refer to deviations of quasimomentum from zero, which vary across space and between individual samples; fluctuations should saturate to a large value when domains having 
are fully formed. After loading the condensates into the lattice, we ramp the shaking amplitude linearly from 
to values well above the critical amplitude 
nm (32) and interrupt the ramps at various times to perform a brief time-of-flight (TOF) before detection. After TOF, we measure the density deviation 
(32), which is nearly proportional to the quasimomentum distribution (fig. S2), where 
is the density profile of the 
th Bragg peak and the angle brackets denote averaging over multiple images. This detection method is particularly sensitive near the critical point when the quasimomentum just starts deviating from zero, indicating the emergence of fluctuations in the ferromagnetic phase where the ground states have nonzero quasimomentum. The spin density measurement used later to study spatial correlations is viable only when atoms settle to 
.
Over a wide range of quench rates, the evolution of quasimomentum fluctuation can be described in three phases (Fig. 2A). First, below the critical point, quasimomentum fluctuation does not exceed its baseline level. Second, just after passing the critical point, critical slowing keeps the system “frozen,” and fluctuation remains low. Finally, the system unfreezes and quasimomentum fluctuation quickly increases and saturates, indicating the emergence of ferromagnetic domains. We quantify this progression by investigating the fluctuation of contrast 
(Fig. 2B) that tracks quasimomentum fluctuation in our condensates via the fluctuation of 
(fig. S2), where 
is the total density and the angle brackets denote averaging over space and over multiple images. For comparison between different quench rates, we calculate the normalized fluctuation 
, where subscripts i and f indicate the fluctuation at early and late times, respectively (32). We find empirically that the growth of normalized fluctuations is well fit by the function
(5)where the time 
is defined relative to when the system crosses the critical point at 
, 
characterizes the delay time when the system unfreezes, and 
is the formation time over which the fluctuation grows. The measurement of fluctuation over time provides a critical test for both the Kibble-Zurek scaling and the universality hypothesis. First, both 
and 
exhibit clear power-law scaling with the quench rate 
varied over more than two orders of magnitude (Fig. 2C). Power-law fits yield the exponents of 
and 
, respectively. The nearly equal exponents are suggestive of the universality hypothesis, which requires all times to scale identically. Indeed, the growth of contrast fluctuation 
follows a universal curve when time is scaled by 
(Fig. 2D), strongly supporting the universality hypothesis (Eq. 4). Note that any observable time characterizing the collective dynamics can be chosen as 
in Eq. 4, including 
and 
.
(A) Sample images show the emergence of nonzero local quasimomentum via the density deviation between Bragg peaks 
as the system is linearly ramped across the ferromagnetic phase transition; four ramps with increasing quench rates (bottom to top) are shown. Each ramp exhibits three regimes: a subcritical regime before the transition; a frozen regime beyond the critical point where the fluctuation remains low; and a growth regime in which fluctuation increases and saturates, indicating domain formation. Time 
corresponds to the moment when the system reaches the critical point. (B) Quasimomentum fluctuation for ramp rates 
3.6 (triangles), 0.91 (squares), 0.23 (diamonds), and 0.06 nm/ms (circles) arises at a delay time 
over a formation time 
. Fluctuation is normalized for each ramp rate to aid comparison (32). The solid curves show fits based on Eq. 5. (C) The dependence of 
(circles) and 
(squares) on the quench rate is well fit by power laws (solid curves) with scaling exponents of 
and 
, respectively. The inset shows 
on a linear scale. (D) Fluctuations measured for 16 ramp rates from 0.06 to 10.3 nm/ms collapse to a single curve when time is scaled by 
based on the power-law fit. The solid curve shows the best fit based on the empirical function (Eq. 5), and the gray shaded region covers one standard deviation. Error bars in (B) and (C) indicate one standard error.
We next test the spatial scaling symmetry based on the structures of pseudo-spin domains that emerge after the system unfreezes. Here, we cross the critical point with two different protocols: The first is a linear ramp starting from 
, whereas the second begins with a jump to 
, followed by a linear ramp. We detect domains near the time 
in the spin density distribution 
based on the density 
of atoms with spin up/down (fig. S3). At this time, the spin domains are fully formed and clearly separated by topological defects (domain walls), as shown in Fig. 3A. Furthermore, choosing this time just after domain formation minimizes the time available for nonuniversal relaxation processes. We characterize the domain distribution with the spin correlation function (17, 28)
(6)averaged over multiple images (Fig. 3B). Both ramping protocols lead to similar correlation functions, suggesting that the domain distribution is insensitive to increases in the quench rate below the critical point.
(A) Two sample images at each quench rate exemplify spin domains measured near the time 
after crossing the phase transition. These images correspond to linear ramps starting at 
(
to 
nm/ms) and 
(
to 
nm/ms). (B) Spin correlation functions 
(Eq. 6) are calculated from ensembles of 110 to 200 images. (C) Cuts across the density-weighted correlation functions 
are shown for quench rates 
1.28 (triangles), 0.45 (squares), 0.16 (diamonds), and 0.056 nm/ms (circles). Solid curves interpolate the data to guide the eye. The typical domain size 
and the correlation length 
are illustrated for 0.056 nm/ms by the arrow and dashed envelope, respectively (32). (D) The dependence of 
(green) and 
(black) on the quench rate is well fit by power laws (Eq. 3) with spatial scaling exponents of 
and 
, respectively. Marker shape indicates linear ramps starting at 
(squares) or at 
(circles). The inset shows the results on a linear scale. Error bars indicate 1 SE. (E) Correlation functions for 
to 
nm/ms collapse to a single curve when distance is scaled by the domain size extracted from the power-law fit in (D). The solid curve shows the fit based on Eq. 7; the gray shaded area covers 1 SD. (F) The temporal scaling exponents 
and 
from Fig. 2C (magenta) and the spatial scaling exponents 
and 
from (D) (green) constrain the critical exponents 
and 
according to Eqs. 2 and 3 with 
(dark) and 
(light) confidence intervals. The cross marks the best values with contours of 
and 
overall confidence (32).
The spin correlations reveal rich domain structure that strongly depends on the quench rate. For slower ramps 
nm/ms, the structures are predominantly one-dimensional and the density of topological defects increases with the quench rate. The tighter confinement and finite speed of sound near the critical point along the 
and 
axes allow spin correlations to span the gas in those directions. The dynamics thus appear one-dimensional. When the quench rate exceeds 
nm/ms, defects start appearing along the 
axis, and the domain structures become multidimensional. We attribute this dimensional crossover to the unfreezing time becoming too short to establish correlation along the 
axis. For the remainder of this work, we focus on the slower quenches and investigate the spin correlations along the 
axis.
We examine the 1D correlations using line cuts of the density-weighted correlation functions 
(17, 28). The results exhibit prominent decaying oscillation (Fig. 3C). We extract two essential length scales from the correlation functions: the average domain size 
, or equivalently, the distance between neighboring topological defects, and the correlation length 
, indicating the width of the envelope function. These two scales are determined from the position and width of the peak in the Fourier transform of 
(32).
These length scales enable us to test the spatial scaling symmetry. The lengths 
and 
both display power-law scaling consistent with the Kibble-Zurek mechanism (Fig. 3D), with fits yielding exponents 
for the domain size and 
for the correlation length. Similarly, the correlations, measured at the same scaled time, collapse to a single curve in spatial coordinates scaled by the domain size 
(Fig. 3E). This result strongly supports the spatiotemporal scaling from the universality hypothesis (Eq. 4). An empirical curve
(7)provides a good fit to the universal correlation function, yielding 
= 1.01(1) and 
= 1.04(1), indicating that the width of the envelope is close to the typical domain size.
The most striking feature of the universal correlation function is the emergence of oscillatory, antiferromagnetic order in the ferromagnetic phase. In thermal equilibrium, ferromagnets are expected to have a finite correlation length but no anticorrelation. The appearance of strong anticorrelation at 
suggests that domains of size 
form preferentially during the quantum critical dynamics. A statistical analysis of the topological defect distribution reveals that the domain sizes are bunched with their standard deviation 
well below their mean, indicating that the topological defects are created by a sub-Poisson process (fig. S4).
Finally, the combined scaling exponents of space and time allow us to extract the equilibrium critical exponents based on the Kibble-Zurek mechanism (33) (Fig. 3F). Solving Eqs. 2 and 3, we obtain the dynamical exponent 
1.9(2) and correlation length exponent 
0.52(5), which are close to the mean-field values 
and 
up to our experimental uncertainty. Note that the dynamical critical exponent 
results from the unique quartic kinetic energy 
of our system at the critical point (32).
Direct identification of domain walls presents intriguing possibilities for future studies of the topological defects generated during critical dynamics. These opportunities would be particularly interesting if the shaking technique were extended to higher dimensions in such a way that the transition breaks a continuous symmetry. In addition, the scaling of the correlation functions suggests that the antiferromagnetic order may be a shared feature of quantum critical dynamics for phase transitions in the same universality class, meriting future experiments.
