Probing many-body dynamics on a 51-atom quantum simulator.

TL;DR: This work demonstrates a method for creating controlled many-body quantum matter that combines deterministically prepared, reconfigurable arrays of individually trapped cold atoms with strong, coherent interactions enabled by excitation to Rydberg states, and realizes a programmable Ising-type quantum spin model with tunable interactions and system sizes of up to 51 qubits.

Abstract: Controllable, coherent many-body systems can provide insights into the fundamental properties of quantum matter, enable the realization of new quantum phases and could ultimately lead to computational systems that outperform existing computers based on classical approaches. Here we demonstrate a method for creating controlled many-body quantum matter that combines deterministically prepared, reconfigurable arrays of individually trapped cold atoms with strong, coherent interactions enabled by excitation to Rydberg states. We realize a programmable Ising-type quantum spin model with tunable interactions and system sizes of up to 51 qubits. Within this model, we observe phase transitions into spatially ordered states that break various discrete symmetries, verify the high-fidelity preparation of these states and investigate the dynamics across the phase transition in large arrays of atoms. In particular, we observe robust many-body dynamics corresponding to persistent oscillations of the order after a rapid quantum quench that results from a sudden transition across the phase boundary. Our method provides a way of exploring many-body phenomena on a programmable quantum simulator and could enable realizations of new quantum algorithms.

Summary

Introduction

• Here, the authors demonstrate a novel platform for the creation of controlled many-body quantum matter.
• The realization of fully controlled, coherent many-body quantum systems is an outstanding challenge in modern science and engineering.
• The authors approach makes use of atom-by-atom assembly to deterministically prepare arrays of individually trapped cold neutral 87Rb atoms in optical tweezers [16, 18, S1].
• Controlled, coherent interactions between these atoms are introduced by coupling them to Rydberg states (Fig. 1a).
• In general, within this platform, one can program the control parameters Ωi,∆i by changing laser intensities and detunings in time.

PROGRAMMABLE QUANTUM SIMULATOR

• In the case of homogeneous coherent coupling considered here, Hamiltonian (1) closely resembles the paradigmatic Ising model for effective spin-1/2 particles with variable interaction range.
• The ground state corresponds to a Rydberg crystal breaking Z2 translational symmetry that is analogous to antiferromagnetic order in magnetic systems.
• The authors then switch on the laser fields and sweep the two-photon detuning from negative to positive values using a functional form shown in Fig. 3a.
• As shown in Fig. 3, this fully coherent simulation without free parameters yields excellent agreement with the observed data when the finite detection fidelity is accounted for.
• With perfect Z2 order is by far the most commonly observed many-body state (Fig. 4b).

QUANTUM DYNAMICS ACROSS A PHASE TRANSITION

• In single instances of the experiment the authors observe long ordered chains where the atomic states alternate between Rydberg and ground state.
• As the system enters the Z2 phase, ordered domains grow in size, leading to a substantial reduction in the domain wall density (blue points in Fig. 5b).
• The red dots are the measured values and the blue dots are corrected for finite detection fidelity (Supplementary Information).
• The measured domain wall number distribution allows us to directly infer the statistics of excitations created when crossing the phase transition.
• After such a quench, the authors observe oscillations of many-body states between the initial crystal and a complementary crystal where each internal atomic state is inverted (Fig. 6a).

DISCUSSION

• Several important features of these experimental observations should be noted.
• More specifically, if an effective temperature is estimated based on the measured domain wall density, the corresponding thermal ensemble predicts a correlation length ξth = 4.48(3), which is significantly longer than the measured value ξ = 3.03(6).
• Indeed, the exact numerics predict that this simplified model exhibits crystal oscillations with the observed frequency, while the entanglement entropy grows at a rate much smaller than Ω, indicating that the oscillation persists over many cycles (Fig. 6d and Supplementary Information).
• This relatively slow thermalization is rather unexpected, since their Hamiltonian, with or without long-range interactions, is far from any known integrable systems [29], and features neither strong disorder [37] nor explicitly conserved quantities [38].
• Instead, their observations are associated with constrained dynamics due to Rydberg blockade and large separations of timescales Vi,i+1 Ω Vi,i+2 [S8], which gives rise to so-called quantum dimer models, with the Hilbert space dimension determined by the golden ratio ∼ (1 + √ 5)N/2N , that are known to possess non-trivial dynamics [40, 41].

OUTLOOK

• The authors observations demonstrate that Rydberg excitation of neutral atom arrays constitutes an exceptionally promising platform for studying quantum dynamics and quantum simulations in large systems.
• Further improvement in coherence and controllability can be obtained by encoding qubits into hyperfine sublevels of the electronic ground state and using state-selective Rydberg excitation [26].
• Implementing two-dimensional (2d) arrays could provide a path towards realizing thousands of traps.
• While their current observations already allow us to gain unprecedented insights into the physics associated with transitions into ordered phases and to explore novel many-body phenomena in quantum dynamics, they can be directly extended along several directions [46].
• Finally, the authors note that their approach is exceptionally well suited for the realization and testing of quantum optimization algorithms with systems that are well beyond the reach of modern classical machines [53, 54].

ACKNOWLEGEMENTS

• This work was supported by NSF, CUA, ARO, MURI, AFOSR, and Vannevar Bush Faculty Fellowship.
• H.B. acknowledges support by a Rubicon Grant of the Netherlands Organization for Scientific Research (NWO).
• A.O. acknowledges support by a research fellowship from the German Research Foundation (DFG).
• S.S. acknowledges funding from the European Union under the Marie Sk lodowska Curie Individual Fellowship Programme H2020-MSCA-IF-2014 (project number 658253).
• H.P. acknowledges support by the National Science Foundation (NSF) through a grant at the Institute of Theoretical Atomic Molecular and Optical Physics at Harvard University and the Smithsonian Astrophysical Observatory.

1.1 Trapping setup and experimental sequence

• The tweezers are generated by feeding a multi-tone RF signal into an acousto-optic deflector (AA Opto-Electronic model DTSX-400-800.850), generating multiple deflections in the first diffraction order, and focusing them into the vacuum chamber using a 0.5 NA objective (Mitutoyo G Plan Apo 50X).
• The traps are loaded from a magneto-optical trap (MOT), leading to individual tweezer single-atom loading probabilities of ∼ 0.6.
• The authors then turn off the traps, pulse the Rydberg lasers on a timescale of a few microseconds, and then turn the traps back on to recapture the atoms that are in the ground state |g〉 while pushing away the atoms in the Rydberg state |r〉, and finally take a third image.
• This provides a convenient way to detect them as missing atoms on the third image (with finite detection fidelity discussed in section 1.3).
• The entire experimental sequence, from MOT formation to the third image, takes ∼ 250 ms.

1.2 Rydberg lasers setup

• The van der Waals interaction between two 87Rb 71S atoms follows a 1/R6 power law and is on the order of 1 MHz at 10µm [S2], making it the dominant energy scale in their system for up to several lattice sites.
• The reflected beam from the cavity is sent on a fast photodetector (Thorlabs PDA8A), whose signal is demodulated and low-pass filtered to create an error signal which is fed into a high-bandwidth servo box (Vescent D2-125).
• The other part of the blue laser beam goes through an acousto-optic modulator (IntraAction ATM1002DA23), whose first diffraction order is used to excite atoms, providing frequency and amplitude control for the Rydberg pulses.
• After the traps are turned back on, a third EMCCD image is taken to detect Rydberg excitations with single-site resolution.

1.3 Detection fidelity

• Detection infidelity arises from accidental loss of atoms in |g〉 or accidental recapture of atoms in |r〉.
• In particular, for the 7-atom data shown in Figure 3 in the main text and the 51-atom data shown in Figure 4 and 5, the authors measured ground state detection fidelities of 98% and 99%, respectively.
• For an atom in state |r〉, the optical tweezer yields an anti-trapping potential, but there is a finite probability that the atom will decay back to the ground state and be recaptured by the tweezer before it can escape the trapping region.
• The authors quantify this probability by measuring Rabi oscillations between |g〉 and |r〉 (Fig S2) and extracting the maximum amplitude of the oscillation signal.
• Furthermore, the authors observe a reduced detection fidelity at lowerlying Rydberg states, which is consistent with the dependence of the Rydberg lifetime on the principal quantum number [S5].

2.1 Pulse optimization

• The detuning ∆ is set to truncate at minimum 11 and maximum values ∆min and ∆max, respectively.
• All parameters in (3) or (4) are iteratively optimized as to minimize the domain wall number, i.e. maximize the crystal preparation fidelity.
• After passing through the AOM, the 420 nm light is coupled into a fiber.
• The power throughout all frequency sweeps is ≥ 75% of the power at fopt.

2.2 Limitations

• When sweeping into the crystalline phase, the control parameter ∆(t) must be varied slowly enough that the adiabaticity criterion is sufficiently met.
• Preparation fidelity is therefore given by the probability that each atom in the array is still present for the Rydberg pulse, and that it is prepared in the correct magnetic sublevel:∣∣5S1/2, F = 2,mF = −2〉.
• The longitudinal position fluctuations add in quadrature, so they contribute less to fluctuations in distance.
• Typical Rabi oscillation, homogeneity and coherence for non-interacting atoms (a = 24µm, Ω Vi,i+1 ' 5 kHz), also known as S2.

3 CORRECTING FOR FINITE DETECTION FIDELITY

• The number of domain walls is a metric for the quality of preparing the desired crystal state.
• Boundary conditions make it favorable to excite the atoms at the edges.
• The appearance of domain walls can arise from nonadiabaticity across the phase transition, as well as scattering from the intermediate 6P state, imperfect optical pumping, atom loss, and other processes (see section 2.2).
• For this reason, the authors use a maximum-likelihood routine to estimate the parent distribution, which is the distribution of domain walls in the prepared state that best predicts the measured distribution.
• The authors use two methods to correct for detection infidelity, depending on whether they are interested only in the probability to generate the many- body ground state, or in the full probability distribution of the number of domain walls.

3.1 Many-body ground state preparation

• Having prepared the many-body ground state, the probability to correctly observe it depends on the measurement fidelity for atoms in the electronic ground state fg, the measurement fidelity for atoms in the Rydberg state fr, and the size of the system N .
• Therefore, if the authors observe the ground state with probability pexp, the probability of actually preparing this state is inferred to be pexp/pm.
• The blue data points in Fig. 4 a in the main text are calculated this way.

3.2 Maximum likelihood state reconstruction

• For this purpose, the authors assume that the density of domain walls is low, such that the probability of preparing two overlapping domain walls, meaning three consecutive atoms in the same state, is negligibly small.
• Under this assumption, misidentifying an atom within a domain wall shifts its location, but does not change the total number.
• B, Comparison of measured domain wall distribution (red) and predicted observation given the parent distribution in a (blue).
• The authors can find the most likely parent distribution, ~Wi, by minimizing the cost function over the different possible ~W ′i , under the constraint that that every element is between 0 and 1, and the sum of the elements is 1.
• For this purpose, the authors use a Sequential Least Square Programming routine.

4 COMPARISON WITH A CLASSICAL THERMAL STATE

• To gain some insight into the states obtained from their preparation protocol (Fig. 3a in the main text), the authors provide a quantitative comparison between experimentally measured quantities and those computed from a thermal ensemble.
• Also, the authors may consider the interactions only up to next-nearest neighbors as the coupling strengths for longer distances are weak compared to the maximum timescale accessible in their experiments.
• The eigenstates of this Hamiltonian are simply 2N classical configurations, where each atom is in either |g〉 or |r〉.
• Since β characterizes the thermal state completely, the authors can extract the corresponding domain wall distribution (Fig. S4c) and the correlation function (Fig. S4d) as described above.
• The authors find that the correlation length in the corresponding thermal state is 15 ξth = 4.48(3), which is significantly longer than the measured correlation length ξ = 3.03(6), from which they deduce that the experimentally prepared state is not thermal.

5.1 Matrix Product State Ansatz

• In addition, the authors replace the nearest neighbor interactions with hard constraints that two neighboring atoms cannot be excited at the same time; such an approximation is well controlled in the limit of Vi,i+1 Ω, as in the case of their experiments, for a time exponentially long in Vi,i+1/Ω [S8].
• In the simplest approximation, one can treat the array of atoms as a collection of independent dimers, |Ψ(t)〉 = ⊗ i |φ(t)〉2i−1,2i, where for each pair of atoms only three states are allowed due to the blockade constraint, |r, g〉, |g, g〉 and |g, r〉.
• This dimer model predicts that each atom flips its state with respect to its initial configuration after a time τ = √ 2π/Ω.
• The blue line shows the evolution of the domain wall density obtained from integrating the variational equation of motion eq. (9) with initial conditions θa = π/2, θb = 0, i.e. the crystalline initial state.
• Red lines correspond to the initial state |g〉⊗N , while blue lines correspond to crystalline initial states.

5.2 Decay of the oscillations and growth of entanglement

• In order to obtain more insight into the dynamics of their system beyond these variational models, the authors use exact numerical simulations to integrate the many-body Schrödinger equation.
• For the disordered initial state, the domain wall density quickly relaxes to a steady state value.
• Points and error bars represent measured values.
• Numerically, the authors treat the strong nearest neighbor interactions perturbatively – by adiabatic eliminations of simultaneous excitation of neighboring Rydberg states – while the weak interactions beyond nearest neighbors are treated exactly.
• From the growth of the entanglement entropy, the authors see that the crystalline initial state still thermalizes slower than the disordered initial state.

5.3 Time evolution via matrix product state algorithm

• The numerical data presented in Fig. 5b and Fig. 6b in the main text are obtained by simulating the evolution of the 51 atom array during the sweep across the phase transition as well as the subsequent sudden quench using a matrix product state algorithm with bond dimension D = 256.
• The authors simulate the entire preparation protocol to generate the Rydberg crystal [Fig. 5 b in the main text], and use the resulting state as an initial state for the time evolution after the sudden quench.
• The authors take into account only nearest neighbor and next-nearest neighbor interactions and neglect small interactions for atoms that are separated by 3 or more sites (as discussed also in Sec. 4).
• The authors account for finite detection fidelities that are determined independently, but otherwise do not include any incoherent mechanisms.
• Remarkably, for local quantities, such as the domain wall density, this fully coherent simulation agrees well with the experimentally measured values.

Figures

FIG. 6: Emergent oscillations in many-body dynamics. a, Schematic sequence (top, showing ∆(t)) involves adiabatic preparation and then a sudden quench to single-atom resonance. The heat map shows the single atom trajectories for a 9 atom cluster. We observe that the initial (left inset) crystal with a Rydberg excitation at every odd trap site collapses after the quench and a crystal with an excitation at every even site builds up (middle inset). At a later time the initial crystal revives (right inset). Error bars denote 68% CI. b, Density of domain walls after the quench. The dynamics decay slowly on a timescale of 0.88 µs. Shaded region represents the standard error of the mean. Solid blue line is a fully coherent MPS simulation [36] with bond dimension D = 256, taking into account measurement fidelity. c, Toy model of non-interacting dimers (see main text). d, Numerical calculations of the dynamics after a quench starting from an ideal 25 atom crystal, obtained from exact diagonalization. Domain wall density as a function of time (red), and growth of entanglement entropy of the half chain (13 atoms) (blue). Dashed lines take into account only nearest neighbor blockade constraint. Solid lines correspond to the full 1/R6 interaction potential.

FIG. 2: Phase diagram and buildup of crystalline phases. a, The schematic ground-state phase diagram of Hamiltonian (1) displays phases with various broken symmetries depending on the interaction range Rb/a (Rb blockade radius, a trap spacing) and detuning ∆ (see main text). Shaded areas indicate potential incommensurate phases (diagram adapted from [29]). b, The buildup of Rydberg crystals on a 13 atom array is observed by slowly changing the laser parameters as indicated by the red arrows in a (see also Fig. 3a). The bottom panel shows a configuration where the atoms are a = 5.9µm apart which results in a nearest neighbor interaction of Vi,i+1 = 2π × 24 MHz and leads to a Z2 order where every other atom is excited to the Rydberg state |r〉. The right bar plot displays the final, position-dependent Rydberg probability (error bars denote 68% CI). The configuration in the middle panel (a = 3.67µm, Vi,i+1 = 2π × 414.3 MHz) results in Z3 order and the top panel (a = 2.95µm, Vi,i+1 = 2π× 1536 MHz) in a Z4 ordered phase. For each configuration, we show a single-shot fluorescence image before (left) and after (right) the pulse. Red circles highlight missing atoms, which are attributed to Rydberg excitations.

FIG. 5: Quantifying Z2 order in a N = 51 atom array. a, Single-shot fluorescence images of a 51 atom array before applying the adiabatic pulse (top row) and after the pulse (bottom three rows correspond to three separate instances). Red circles mark missing atoms, which are attributed to Rydberg excitations. Domain walls are identified as either two neighboring atoms in the same state or a ground state atom at the edge of the array (Supplementary Information), and are indicated with ellipses. Long Z2 ordered chains between domain walls can be observed. b, Blue points show the mean of the domain wall density as a function of detuning during the sweep. Error bars show the standard error of the mean, and are smaller than the marker size. The red points are the corresponding variances, where the error bars represent one standard deviation. The onset of the phase transition is witnessed by a decrease in the domain wall density and a peak in the variance (see main text for details). Each point is obtained from ∼ 1000 realizations. The solid blue curve is a fully coherent MPS simulation without free parameters (bond dimension D = 256), taking measurement fidelities into account. c, Domain wall number distribution for ∆ = 2π× 14 MHz, obtained from 18439 experimental realizations (blue bars, top plot). Error bars indicate 68% CI. Owing to the boundary conditions, only even numbers of domain walls can appear (Supplementary Information). Green bars in the bottom plot show the distribution obtained by correcting for finite detection fidelity using a maximum likelihood method (Supplementary Information), which results in an average number of 5.4 domain walls. Red bars show the distribution of a thermal state with the same mean domain wall density (Supplementary Information). d, Measured correlation function (2) in the Z2 phase.

Full text

Probing many-body dynamics on a 51-atom quantum simulator
Hannes Bernien,
1
Sylvain Schwartz,
1, 2
Alexander Keesling,
1
Harry Levine,
1
Ahmed Omran,
1
Hannes Pichler,
3, 1
Soonwon Choi,
1
Alexander S. Zibrov,
1
Manuel Endres,
4
Markus Greiner,
1
Vladan Vuleti´c,
2
and Mikhail D. Lukin
1
1
Department of Physics, Harvard University, Cambridge, MA 02138, USA
2
Department of Physics and Research Laboratory of Electronics,
Massachusetts Institute of Technology, Cambridge, MA 02139, USA
3
ITAMP, Harvard-Smithsonian Center for Astrophysics, Cambridge, MA 02138, USA
4
Division of Physics, Mathematics and Astronomy, California Institute of Technology, Pasadena, CA 91125, USA
Controllable, coherent many-body systems provide unique insights into fundamental properties
of quantum matter, allow for the realization of novel quantum phases, and may ultimately lead
to computational systems that are exponentially superior to existing classical approaches. Here,
we demonstrate a novel platform for the creation of controlled many-body quantum matter. Our
approach makes use of deterministically prepared, reconfigurable arrays of individually controlled,
cold atoms. Strong, coherent interactions are enabled by coupling to atomic Rydberg states. We
realize a programmable Ising-type quantum spin model with tunable interactions and system sizes of
up to 51 qubits. Within this model we observe transitions into ordered states (Rydberg crystals) that
break various discrete symmetries, verify high-fidelity preparation of ordered states, and investigate
dynamics across the phase transition in large arrays of atoms. In particular, we observe a novel type
of robust many-body dynamics corresponding to persistent oscillations of crystalline order after
a sudden quantum quench. These observations enable new approaches for exploring many-body
phenomena and open the door for realizations of novel quantum algorithms.
The realization of fully controlled, coherent many-body
quantum systems is an outstanding challenge in modern
science and engineering. As quantum simulators, they
can provide unique insights into strongly correlated quan-
tum systems and the role of quantum entanglement [1],
and enable realizations and studies of new states of mat-
ter, even away from equilibrium. These systems also form
the basis for the realization of quantum information pro-
cessors [2]. While basic building blocks of such proces-
sors have been demonstrated in systems of a few coupled
qubits [3–5], the current challenge is to increase the num-
ber of coherently coupled qubits to potentially perform
tasks that are beyond the reach of modern classical ma-
chines.
A number of physical platforms are currently being
explored to reach these challenging goals. Systems com-
posed of about 10-20 individually controlled atomic ions
have been used to create entangled states and explore
quantum simulations of Ising spin models [6, 7]. Sim-
ilarly sized systems of programmable superconducting
qubits have been recently implemented [8, 9]. Quan-
tum simulations have been carried out in larger sys-
tems of over 100 trapped ions without individual ad-
dressing and control [10]. Strongly interacting quantum
dynamics has been explored using optical lattice simula-
tors [11]. These systems are already addressing computa-
tionally difficult problems in quantum dynamics [12] and
the fermionic Hubbard model [13]. Larger-scale Ising-
like machines have been realized in superconducting [14]
and optical [15] systems but these realizations lack either
coherence or quantum nonlinearity that are essential for
achieving full quantum speedup.
STRONGLY INTERACTING ATOM ARRAYS
Our approach makes use of atom-by-atom assembly to
deterministically prepare arrays of individually trapped
cold neutral
87
Rb atoms in optical tweezers [16, 18, S1].
Controlled, coherent interactions between these atoms
are introduced by coupling them to Rydberg states
(Fig. 1a). This results in repulsive van der Waals in-
teractions (V
ij
=
C
/R
6
ij
, C > 0) between Rydberg atom
pairs at a distance R
ij
[19]. Such interactions have al-
ready been used for realizing quantum gates [20–22], im-
plementing strong photon-photon interactions [23] and
studying many-body physics [24–26]. The quantum dy-
namics of this system is governed by the Hamiltonian
H
~
=
X
i
Ω
i
2
σ
i
x
−
X
i
∆
i
n
i
+
X
i<j
V
ij
n
i
n
j
, (1)
where ∆
i
are the detunings of the driving lasers from
the Rydberg state (Fig. 1b), σ
i
x
= |g
i
ihr
i
| + |r
i
ihg
i
| de-
scribes the coupling between the ground state |gi and
the Rydberg state |ri of an atom at position i, driven
at Rabi frequency Ω
i
, and n
i
= |r
i
ihr
i
|. In general,
within this platform, one can program the control param-
eters Ω
i
, ∆
i
by changing laser intensities and detunings
in time. Here, we focus on homogeneous coherent cou-
pling (|Ω
i
|= Ω, ∆
i
= ∆). The interaction strength V
ij
is
tuned by either varying the distance between the atoms
or coupling them to a different Rydberg state.
The experimental protocol that we implement is de-
picted in Fig. 1c. First, atoms are loaded from a
magneto-optical trap into a tweezer array created by an
acousto-optic deflector (AOD). We then use a measure-
ment and feedback procedure that eliminates the entropy
arXiv:1707.04344v1 [quant-ph] 13 Jul 2017

2
1013 nm 420 nm
V
ij
a
c
1. Load
2. Arrange
3. Evolve
4. Detect
b
d
Time (¹s)
0
0.5
1
0
0.5
Single Rydberg probability
0 0.5 1 1.5
0
0.5
FIG. 1: Experimental platform. a, Individual
87
Rb atoms
are trapped using optical tweezers and arranged into defect-
free arrays. Coherent interactions V
ij
between the atoms are
enabled by exciting them to a Rydberg state, with strength Ω
and detuning ∆. b, A two photon process is used to couple the
ground state |gi =


5S
1
/2
, F = 2, m
F
= −2

to the Rydberg
state |ri =


71S
1
/2
, J =
1
/2, m
J
= −
1
/2

via an intermediate
state |ei =


6P
3
/2
, F = 3, m
F
= −3

using circularly polarized
420 nm and 1013 nm lasers (typically δ ∼ 2π × 560MHz 
Ω
B
, Ω
R
∼ 2π × 60, 36 MHz). c, The experimental protocol
consists of loading the atoms into a tweezer array (1) and
rearranging them into a preprogrammed configuration (2).
After this, the system evolves under U (t) with tunable pa-
rameters ∆(t), Ω(t), V
ij
. This can be implemented in parallel
on several non-interacting sub-systems (3). We then detect
the final state by fluorescence imaging (4). d, For resonant
driving (∆ = 0), isolated atoms (blue points) display Rabi
oscillations between |gi and |ri. Arranging the atoms into
fully blockaded clusters of N = 2 (green) and N = 3 (red)
atoms results in only one excitation being shared between the
atoms in the cluster, while the Rabi frequency is enhanced
by
√
N. The probability to detect more than one excitation
in the cluster is ≤ 5%. Error bars indicate 68% confidence
intervals (CI) and are smaller than the marker size.
associated with the probabilistic trap loading and results
in the rapid production of defect-free arrays with over 50
laser cooled atoms as described previously [S1]. These
atoms are prepared in a preprogrammed spatial configu-
ration in a well-defined internal ground state |gi (Supple-
mentary Information). We then turn off the traps and
let the system evolve under the unitary time evolution
U(Ω, ∆, t), which is realized by coupling the atoms to
the Rydberg state |ri =


71S
1/2

with laser light along
the array axis (Fig. 1a). The final states of individual
atoms are detected by turning the traps back on, and
imaging the recaptured ground state atoms via atomic
fluorescence, while the anti-trapped Rydberg atoms are
ejected (Supplementary Information).
The strong, coherent interactions between Rydberg
atoms provide an effective coherent constraint that pre-
vents simultaneous excitation of nearby atoms into Ryd-
berg states. This is the essence of the so-called Rydberg
blockade [19], demonstrated in Fig. 1d. When two atoms
are sufficiently close so that their Rydberg-Rydberg inter-
actions V
ij
exceed the effective Rabi frequency Ω, then
multiple Rydberg excitations are suppressed. This de-
fines the Rydberg blockade radius, R
b
, for which V
ij
= Ω
(R
b
= 9 µm for |ri =


71S
1/2

and Ω = 2π × 2 MHz as
used here). In the case of resonant driving of atoms sep-
arated by a = 24 µm, we observe Rabi oscillations associ-
ated with non-interacting atoms (blue curve on Fig. 1d).
However, the dynamics change significantly as we bring
multiple atoms close to each other (a = 2.95 µm < R
b
).
In this case, we observe Rabi oscillations between the
ground state and a collective W-state with exactly one
excitation ∼
P
i
Ω
i
|g
1
...r
i
...g
N
i with the characteristic
√
N-scaling of the collective Rabi frequency [25, 27, 28].
These observations allow us to quantify the coherence
properties of our system (see Supplementary Information
for details). In particular, the contrast of Rabi oscilla-
tions in Fig. 1d is mostly limited by the state detection
fidelity (93% for |ri and ∼ 98% for |gi, Supplementary
Information). The individual Rabi frequencies are con-
trolled to better than 3% across the array, while the co-
herence time is ultimately limited by the small probabil-
ity of spontaneous emission from the intermediate state
|ei during the laser pulse (scattering rate 0.022/µs, Sup-
plementary Information).
PROGRAMMABLE QUANTUM SIMULATOR
In the case of homogeneous coherent coupling consid-
ered here, Hamiltonian (1) closely resembles the paradig-
matic Ising model for effective spin-1/2 particles with
variable interaction range. Its ground state exhibits a
rich variety of many-body phases that break distinct spa-
tial symmetries (Fig. 2a). Specifically, at large, negative
values of
∆
/Ω, its ground state corresponds to all atoms
in the state |gi, corresponding to paramagnetic or disor-
dered phase. As
∆
/Ω is increased towards large positive
values, the number of atoms in |ri rises and interactions
between them become significant. This gives rise to spa-
tially ordered phases where Rydberg atoms are regularly
arranged across the array, resulting in ‘Rydberg crys-
tals’ with different spatial symmetries [29–31], as illus-
trated in Fig. 2a. The origin of these correlated states

3
Z
4
ordered
Detuning ( )
disordered
Interaction
range ( )
Z
3
ordered
Z
2
ordered
a
-4 0 4 8
Detuning (MHz)
1
5
9
13
0 0.5 1
1
5
9
13
0 0.5 1
1
5
9
13
Position in cluster
0 0.5 1
0 0.5 1
Rydberg probability
b
FIG. 2: Phase diagram and buildup of crystalline phases. a, The schematic ground-state phase diagram of Hamilto-
nian (1) displays phases with various broken symmetries depending on the interaction range R
b
/a (R
b
blockade radius, a trap
spacing) and detuning ∆ (see main text). Shaded areas indicate potential incommensurate phases (diagram adapted from [29]).
b, The buildup of Rydberg crystals on a 13 atom array is observed by slowly changing the laser parameters as indicated by
the red arrows in a (see also Fig. 3a). The bottom panel shows a configuration where the atoms are a = 5.9 µm apart which
results in a nearest neighbor interaction of V
i,i+1
= 2π × 24 MHz and leads to a Z
2
order where every other atom is excited
to the Rydberg state |ri. The right bar plot displays the final, position-dependent Rydberg probability (error bars denote
68% CI). The configuration in the middle panel (a = 3.67 µm, V
i,i+1
= 2π × 414.3 MHz) results in Z
3
order and the top panel
(a = 2.95 µm, V
i,i+1
= 2π ×1536 MHz) in a Z
4
ordered phase. For each configuration, we show a single-shot fluorescence image
before (left) and after (right) the pulse. Red circles highlight missing atoms, which are attributed to Rydberg excitations.
can be understood intuitively by first considering the sit-
uation when V
i,i+1
 ∆  Ω  V
i,i+2
, i.e. blockade
for neighboring atoms but negligible interaction between
next-nearest neighbors. In this case, the ground state
corresponds to a Rydberg crystal breaking Z
2
transla-
tional symmetry that is analogous to antiferromagnetic
order in magnetic systems. Moreover, by tuning the pa-
rameters such that V
i,i+1
, V
i,i+2
 ∆  Ω  V
i,i+3
and
V
i,i+1
, V
i,i+2
, V
i,i+3
 ∆  Ω  V
i,i+4
, we obtain arrays
with broken Z
3
and Z
4
symmetries, respectively (Fig. 2).
To prepare the system in these phases, we dynamically
control the detuning ∆(t) of the driving lasers to adia-
batically transform the ground state of the Hamiltonian
from a product state of all atoms in |gi into crystalline
states [31, 32]. In the experiment, we first prepare
all atoms in state |gi =


5S
1
/2
, F = 2, m
F
= −2

by
optical pumping. We then switch on the laser fields
and sweep the two-photon detuning from negative
to positive values using a functional form shown in
Fig. 3a. Fig. 2b displays the resulting single atom
trajectories in a group of 13 atoms for three different
interaction strengths as we vary the detuning ∆. In
each of these instances, we observe a clear transition
from the initial state |g
1
, ..., g
13
i to an ordered state of
different broken symmetry. The distance between the
atoms determines the interaction strength which leads
to different crystalline orders for a given final detuning.
To achieve a Z
2
order, we arrange the atoms with a
spacing of 5.9 µm, which results in a nearest neighbor
interaction of V
i,i+1
= 2π × 24 MHz  Ω = 2π × 2 MHz,
while the next-nearest neighbor interaction is small
(2π × 0.38 MHz). This results in a buildup of antiferro-
magnetic order where every other trap site is occupied
by a Rydberg atom (Z
2
order). By reducing the spacing
between the atoms to 3.67 µm and 2.95 µm, Z
3
- and Z
4
-
orders are respectively observed (Fig. 2b).
We benchmark the performance of the quantum simu-
lator by comparing the measured Z
2
order buildup with
theoretical predictions for a N = 7 atom system, ob-
tained via exact numerical simulations. As shown in
Fig. 3, this fully coherent simulation without free pa-
rameters yields excellent agreement with the observed
data when the finite detection fidelity is accounted for.
The evolution of the many-body states in Fig. 3c shows
that we measure the perfect antiferromagnetic state with
54(4)% probability. When corrected for the known detec-
tion infidelity, we find that the desired many-body state
is reached with a probability of p = 77(6)%.
To investigate how the preparation fidelity depends on
system size, we perform detuning sweeps on arrays of
various sizes (Fig. 4a). We find that the probability of
observing the system in the many-body ground state at
the end of the sweep decreases as the the system size
is increased. However, even at system sizes as large as
51 atoms, the perfectly ordered crystalline many-body
state is obtained with p = 0.11(2)% (p = 0.9(2)% when
corrected for detection fidelity), which is remarkable in
view of the exponentially large 2
51
-dimensional Hilbert
space of the system. Furthermore, we find that this state

4
-20
0
20
¢ (MHz)
0
0.25
0.5
0.75
1
Rydberg
probability
1 2 3 4 5 6 7
0 1 2 3
Time t
stop
(¹s)
0
0.25
0.5
0.75
1
State probability
jggggggg
®
jrgrgrgg
®
jrgggggr
®
jrgrgrgr
®
jrgrgggr
®
jrggrggr
®
jrgggrgr
®
0
1
2
 (MHz)
t
stop
a
b
c
FIG. 3: Comparison with a fully coherent simulation.
a, The laser driving consists of a square shaped pulse Ω(t)
with a detuning ∆(t) that is chirped from negative to positive
values. b, Time evolution of Rydberg excitation probability
for each atom in a N = 7 atom cluster (colored points), ob-
tained by varying the duration of laser excitation pulse Ω(t).
The corresponding curves are theoretical single atom trajec-
tories obtained by an exact simulation of quantum dynamics
with (1), the functional form of ∆(t) and Ω(t) used in the
experiment, and finite detection fidelity. c, Evolution of the
seven most probable many-body states. The target state is
reached with 54(4)% probability (77(6)% when corrected for
finite detection fidelity). Error bars denote 68% CI.
with perfect Z
2
order is by far the most commonly ob-
served many-body state (Fig. 4b).
QUANTUM DYNAMICS ACROSS A PHASE
TRANSITION
We next present a detailed study of the transition into
the Z
2
phase in an array of 51 atoms (Fig. 5). In single in-
stances of the experiment we observe long ordered chains
where the atomic states alternate between Rydberg and
ground state. These ordered domains can be separated
by domain walls that consist of two neighboring atoms
in the same electronic state (Fig. 5a) [33].
The domain wall density can be used to quantify the
transition from the disordered phase into the ordered Z
2
phase as a function of detuning ∆. As the system enters
the Z
2
phase, ordered domains grow in size, leading to
a substantial reduction in the domain wall density (blue
points in Fig. 5b). Consistent with expectations for an
Ising-type second-order quantum phase transition [33],
we observe domains of fluctuating lengths close to the
transition point between the two phases, which is re-
flected by a pronounced peak in the variance of the den-
0 20 40
System size
10
-3
10
-2
10
-1
1
Ground state probability
Measured
Corrected
0 10 20
Number of occurences
1
10
10
2
10
3
10
4
Number of states
jr
1
g
2
r
3
: : : : r
49
g
50
r
51
®
a b
FIG. 4: Scaling behavior. a, Preparation fidelity of the
crystalline ground state as a function of cluster size. The red
dots are the measured values and the blue dots are corrected
for finite detection fidelity (Supplementary Information). Er-
ror bars denote 68% CI. b, Number of observed many-body
states per number of occurrences out of 18439 experimental
realizations in a 51-atom cluster. The most occurring state is
the ground state of the many-body Hamiltonian.
sity of domain walls. Consistent with predictions from
finite-size scaling analysis [29, 34], this peak is shifted
towards positive values of
∆
/Ω. The measured position
of the peak is ∆ ' 0.5Ω. The observed domain wall
density is in excellent agreement with fully coherent sim-
ulations of the quantum dynamics based on 51-atom ma-
trix product states (blue line); however, these simulations
underestimate the variance at the phase transition (Sup-
plementary Information).
At the end of the sweep, deep in the Z
2
phase (
∆
/Ω 
1) we can neglect Ω such that the Hamiltonian (1) be-
comes essentially classical. In this regime, the mea-
sured domain wall number distribution allows us to di-
rectly infer the statistics of excitations created when
crossing the phase transition. From 18439 experimen-
tal realizations we obtain the distribution depicted in
Fig. 5c with an average of 9.01(2) domain walls. From a
maximum-likelihood estimation we obtain the distribu-
tion corrected for detection fidelity (Supplementary In-
formation), which corresponds to a state that has on av-
erage 5.4 domain walls. These domain walls are most
likely created due to non-adiabatic transitions from the
ground state when crossing the phase transition, where
the energy gap becomes minimal [35]. In addition, the
preparation fidelity is also limited by spontaneous emis-
sion during the laser pulse (an average number of 1.1
photons is scattered per µs for the entire array, see Sup-
plementary Information).
To further characterize the created Z
2
ordered state,
we evaluate the correlation function
g
(2)
ij
= hn
i
n
j
i − hn
i
ihn
j
i (2)
where the average h···i is taken over experimental repe-
titions. We find that the correlations decay exponentially
over distance with a decay length of ξ = 3.03(6) sites (see
Fig. 5d and Supplementary Information).

5
0 8 16 24
0
0.1
0.2
0.3
Measured
0 8 16 24
Domain wall number
0
0.1
0.2
0.3
Probability
Corrected
Thermal
1 11 21 31 41 51
Position i
1
11
21
31
41
51
Position j
-0.15
0
0.15
g
(2)
i; j
-15 0 15
Detuning (MHz)
0
0.5
1
Domain wall density
Mean
Variance
MPS
a
b c d
FIG. 5: Quantifying Z
2
order in a N = 51 atom array. a, Single-shot fluorescence images of a 51 atom array before
applying the adiabatic pulse (top row) and after the pulse (bottom three rows correspond to three separate instances). Red
circles mark missing atoms, which are attributed to Rydberg excitations. Domain walls are identified as either two neighboring
atoms in the same state or a ground state atom at the edge of the array (Supplementary Information), and are indicated with
ellipses. Long Z
2
ordered chains between domain walls can be observed. b, Blue points show the mean of the domain wall
density as a function of detuning during the sweep. Error bars show the standard error of the mean, and are smaller than
the marker size. The red points are the corresponding variances, where the error bars represent one standard deviation. The
onset of the phase transition is witnessed by a decrease in the domain wall density and a peak in the variance (see main text
for details). Each point is obtained from ∼ 1000 realizations. The solid blue curve is a fully coherent MPS simulation without
free parameters (bond dimension D = 256), taking measurement fidelities into account. c, Domain wall number distribution
for ∆ = 2π ×14 MHz, obtained from 18439 experimental realizations (blue bars, top plot). Error bars indicate 68% CI. Owing
to the boundary conditions, only even numbers of domain walls can appear (Supplementary Information). Green bars in
the bottom plot show the distribution obtained by correcting for finite detection fidelity using a maximum likelihood method
(Supplementary Information), which results in an average number of 5.4 domain walls. Red bars show the distribution of a
thermal state with the same mean domain wall density (Supplementary Information). d, Measured correlation function (2) in
the Z
2
phase.
Finally, Fig. 6 demonstrates that our approach also en-
ables the study of coherent dynamics of many-body sys-
tems far from equilibrium. Specifically, we focus on the
quench dynamics of Rydberg crystals initially prepared
deep in the Z
2
ordered phase, as we suddenly change
the detuning ∆(t) to the single-atom resonance ∆ = 0
(Fig. 6a). After such a quench, we observe oscillations
of many-body states between the initial crystal and a
complementary crystal where each internal atomic state
is inverted (Fig. 6a). We find that these oscillations are
remarkably robust, persisting over several periods with
a frequency that is largely independent of the system
size. This is confirmed by measuring the dynamics of the
domain wall density, signaling the appearance and dis-
appearance of the crystalline states, shown in Fig. 6b for
arrays of 9 and 51 atoms. We find that the initial crys-
tal repeatedly revives with a period that is slower by a
factor ∼ 1.4 compared to the Rabi oscillation period for
independent, non-interacting atoms.
DISCUSSION
Several important features of these experimental ob-
servations should be noted. First of all, our Z
2
ordered
state cannot be characterized by a simple thermal en-
semble. More specifically, if an effective temperature is
estimated based on the measured domain wall density,
the corresponding thermal ensemble predicts a correla-
tion length ξ
th
= 4.48(3), which is significantly longer
than the measured value ξ = 3.03(6). Such a discrep-
ancy is also reflected in distinct probability distributions
for the number of domain walls (Fig. 5c). These observa-
tions suggest that the system does not thermalize within
the timescale of the Z
2
state preparation.
Even more striking is the coherent and persistent oscil-
lation of the crystalline order after the quantum quench.
With respect to the quenched Hamiltonian (∆ = 0), the
energy density of our Z
2
ordered state corresponds to
that of an infinite-temperature ensemble within the man-
ifold constrained by Rydberg blockade. Also, our Hamil-
tonian does not have any explicit conserved quantities
other than total energy. Nevertheless, the oscillations
persist well beyond the natural timescale of local relax-
ation ∼ 1/Ω as well as the fastest timescale, 1/V
i,i+1
.
To understand these observations, we consider a sim-
plified model where the effect of long-range interactions is
neglected, and nearest-neighbor interactions are replaced
by hard constraints on neighboring excitations of Ry-
dberg states [29]. In this limit, the qualitative behav-

Frequently asked questions

Q1. What are the contributions in "Probing many-body dynamics on a 51-atom quantum simulator" ?

Bernien, Sylvain Schwartz, 2 Alexander Keesling, Harry Levine, Ahmed Omran, Hannes Pichler, 1 Soonwon Choi, Alexander S. Zibrov, Manuel Endres, Markus Greiner, Vladan Vuletić, and Mikhail D. Lukin Department of Physics, Harvard University, Cambridge, MA 02138, USA Department of physics and Research Laboratory of Electronics, Massachusetts Institute of Technology, Cambridge 

Q2. What are the future works in "Probing many-body dynamics on a 51-atom quantum simulator" ?

While these considerations provide important insights into the origin of robust emergent dynamics, the authors emphasize that their results challenge conventional theoretical concepts and warrant further studies. 

Q3. What is the simplest way to quantify the transition from the disordered phase into the ordered?

The domain wall density can be used to quantify the transition from the disordered phase into the ordered Z2 phase as a function of detuning ∆. 

Q4. How do the authors propagate the state vector of up to 25 spins?

The authors make use of the constrained size of the Hilbert space (blockade of nearest neighboring excitations of Rydberg states), and propagate the state vector of up to 25 spins using a Krylov subspace projection method. 

Q5. How long does the initial crystal revive?

The authors find that the initial crystal repeatedly revives with a period that is slower by a factor ∼ 1.4 compared to the Rabi oscillation period for independent, non-interacting atoms. 

Q6. How do you tune the frequency of the red laser?

tune the frequency of the red laser over a full free-spectral range of the reference cavity (1.5 GHz) by tuning the driving frequency of the high-bandwidth EOM. 

Q7. Why is the tangent adiabatic sweep used for the main text?

The tangent adiabatic sweep has been used for datasets with 51 atoms shown on Figures 4 and 5 of the main text due to improved performance, whereas the cubic form has been used for all smaller system sizes and for the data on crystal dynamics shown on Figure 6 of the main text. 

Q8. What is the probability to measure exactly n domain walls in a thermal ensemble?

In particular, the probability to measure exactly n domain walls pn = tr {Pnρ} can be computed from a Fourier transform of the Kronecker delta function Pn ≡ δD,n = 1N+2 ∑N+1 k=0 exp[i 2π N+2k(n − D)] with n = 0, 1, 2, . . . 

Q9. How fast does the entanglement entropy grow?

The ad-dition of long-range interactions leads to a faster decay of the oscillations, with a timescale that is determined by ∼ 1/Vi,i+2, in good agreement with experimental observations (Fig. 6b), while the entanglement entropy also grows on this time scale (Fig. 6d). 

Q10. What is the optimum frequency for a tangent adiabatic sweep?

When sweeping into the crystalline phase, the control parameter ∆(t) must be varied slowly enough that the adiabaticity criterion is sufficiently met. 

Q11. How does the origin of the correlated states of the Hamiltonian be understood?

The origin of these correlated statescan be understood intuitively by first considering the situation when Vi,i+1 ∆ Ω Vi,i+2, i.e. blockade for neighboring atoms but negligible interaction between next-nearest neighbors. 

Q12. How many atoms are allowed to be a pair of independent dimers?

In the simplest approximation, one can treat the array of atoms as a collection of independent dimers, |Ψ(t)〉 = ⊗ i |φ(t)〉2i−1,2i, where for each pair of atoms only three states are allowed due to the blockade constraint, |r, g〉, |g, g〉 and |g, r〉. 

Q13. How can the authors evaluate the number of domain walls in a thermal ensemble?

using this approach, the authors can evaluate all measurable quantities for the thermal ensemble such as the average number of domain walls 〈D〉 = tr {Dρ}, where D is an operator counting the number of domain walls, i.e. D = ∑N−1 i=1 (nini+1 + (1− ni)(1− ni+1)) + (1− n1) + (1− nN ), the correlation function g(2)(d) = 1/(N − d) 

Q14. How long is the correlation length in the corresponding thermal state?

The authors find that the correlation length in the corresponding thermal state is15ξth = 4.48(3), which is significantly longer than the measured correlation length ξ = 3.03(6), from which the authors deduce that the experimentally prepared state is not thermal. 

Q15. How do the authors determine the contribution to the laser linewidth of the noise?

Following [S3] and [S4], the authors estimate that the contribution to the laser linewidth of the noise within the servo loop relative to the cavity is less than 500 Hz (see section 2.2 for further discussions on laser frequency noise). 

Q16. What is the fidelity of the atom array?

Preparation fidelity is therefore given by the probability that each atom in the array is still present for the Rydberg pulse, and that it is prepared in the correct magnetic sublevel:∣∣5S1/2, F = 2,mF = −2〉.