Generation and manipulation of Schrödinger cat states in Rydberg atom arrays
Abstract
Quantum entanglement involving coherent superpositions of macroscopically distinct states is among the most striking features of quantum theory, but its realization is challenging, since such states are extremely fragile. Using a programmable quantum simulator based on neutral atom arrays with interactions mediated by Rydberg states, we demonstrate the deterministic generation of “Schrödinger cat” states of the GreenbergerHorneZeilinger (GHZ) type with up to 20 qubits. Our approach is based on engineering the energy spectrum and using optimal control of the manybody system. We further demonstrate entanglement manipulation by using GHZ states to distribute entanglement to distant sites in the array, establishing important ingredients for quantum information processing and quantum metrology.
GHZ states constitute an important class of entangled manybody states Greenberger1989 (). On the one hand, such states provide an important resource for applications ranging from quantum metrology Pezze2018 () to quantum error correction Nielsen2011 (). On the other hand, these states are among the most fragile manybody states, since a single error on any one of the qubits collapses the superposition resulting in a statistical mixture. Remarkably, despite their highly entangled nature, GHZ states can be characterized by just two diagonal and two offdiagonal terms in the particle density matrix. In contrast to quantifying the degree of general entangled states in manybody systems, which is extremely challenging Amico2008 (); Guhne2009 (); Islam2015 (), the GHZ state fidelity () constitutes an accessible witness for partite entanglement Sackett2000 (). For these reasons, GHZ state creation can serve as an important benchmark to characterize the quality of any given quantum hardware. Such states have been previously generated using systems of nuclear spins Laflamme1998 (); Neumann2008 (), individually controlled optical photons Bouwmeester1999 (); Pan2001 (); Wang2018a (), trapped ions Sackett2000 (); Leibfried2005 (); Monz2011 (), and superconducting quantum circuits DiCarlo2010 (); Song2017a (); Song2019 (). Largescale superposition states have also been generated in systems of microwave photons Vlastakis2013 () and atomic ensembles without individual particle addressing Pezze2018 ().
In this Report we demonstrate the preparation of particle GHZ states
(1) 
in a one dimensional array of individually trapped neutral atoms, where the qubits are encoded in an atomic ground state and in a Rydberg state . Our entangling operation relies on the strong vanderWaals interaction between atoms in states and on engineering the energy spectrum of the quantum manybody system to allow for a robust quantum evolution from an initial product state to a GHZ state. The basic ingredients for the manipulation of atomic states, both for generating and characterizing GHZ states, are illustrated in Fig. 1. All the atoms are homogeneously coupled to a Rydberg state via a twophoton transition with an effective coupling strength and detuning Labuhn2016 (); Bernien2017 (). In addition, we use addressing beams to introduce local energy shifts on specific sites along the array (Fig. 1A). The resulting manybody Hamiltonian is
(2) 
where is the qubit flip operator, is the local effective detuning set by the Rydberg laser and the local light shift, is the number of Rydberg excitations on site , and is the interaction strength of two Rydberg atoms on neighboring sites. The separation between adjacent sites is chosen such that the nearestneighbor interaction MHz results in the Rydberg blockade Jaksch2000 (); Wilk2010 (); Isenhower2010a (), forbidding the simultaneous excitation of adjacent atoms into the state .
To prepare GHZ states, we utilize arrays with an even number of atoms. For large negative detuning of the Rydberg laser, the manybody ground state of the Hamiltonian (2) is . For large uniform positive detuning , the ground state manifold consists of nearly degenerate classical configurations with Rydberg excitations. These include in particular the two target antiferromagnetic configurations and Islam2013 (), as well as other states with nearly identical energy (up to a weak secondnearest neighbor interaction), with both edges excited, such as . To isolate a coherent superposition of states and , we introduce local light shifts using offresonant laser beams at nm, generated by an acoustooptic deflector (AOD), which energetically penalize the excitation of edge atoms (Fig. 1A), and effectively eliminate the contribution of undesired components. In this case, the ground state for positive detuning is given by the GHZ state (1) and there exists, in principle, an adiabatic pathway that transforms the state into by adiabatically increasing from negative to positive values (Fig. 1B).
In practice, the time necessary to adiabatically prepare such a GHZ state grows with system size and becomes prohibitively long for large owing to small energy gaps in the manybody spectrum. To address this limitation, we employ optimal control methods to find laser pulses that maximize the GHZ state preparation fidelity while minimizing the amount of time necessary. Our specific implementation, the Remote dressed ChoppedRAndom Basis algorithm (RedCRAB) Rach2015 (); Heck2018 (), yields optimal shapes of the laser intensity and detuning for the given experimental conditions SOM (). For atoms, we perform this optimization using MHz light shifts on the edge atoms. For larger systems of , we found the preparation to be more robust by increasing the edge light shifts to MHz and adding MHz light shifts on the third site from both edges.
Our experiments are based on the optical tweezer platform and experimental procedure that have been described previously Bernien2017 (). Following the initialization of a defectfree atom array, the traps are switched off while the atoms are illuminated with the Rydberg and local light shift beams. The internal state of the atoms is subsequently measured by imaging state atoms recaptured in the traps, while Rydberg atoms are repelled by the trapping light deLeseleuc2018 (). Fig. 2 demonstrates the result of such experiments for a atom array. After applying the optimized pulse shown in Fig. 2B, we measure the probability of observing different patterns in the computational basis, where is the density operator of the prepared state. Fig. 2A shows the measured probability to observe each of the possible patterns in a atom array. The states and clearly stand out (blue bars) with a combined probability of and almost equal probability of observing each one.
To characterize the experimentally prepared state , we evaluate the GHZ state fidelity
(3) 
where and are the populations in the target components and is the offdiagonal matrix element, which can be measured by utilizing the maximal sensitivity of the GHZ state to a staggered magnetic field. Specifically, evolving the systems with the Hamiltonian , the amplitude acquires a phase at a rate of . Measuring an observable that oscillates at this frequency provides a lower bound on the coherence through the oscillation contrast Garttner2017 (); SOM (). In our experiments, the staggered field is implemented by applying offresonant focused beams of equal intensity at nm, generated by another AOD, to every other site of the array (Fig. 1C), resulting in a local energy shift SOM (). Subsequently, we drive the atoms resonantly, applying a unitary operation in order to change the measurement basis SOM (), such that a measurement of the parity becomes sensitive to the phase of . Fig. 2C shows the measured parity as a function of the phase accumulated on each atom in the array, demonstrating the coherence of the created state.
To extract the entanglement fidelity for large atomic states, we carefully characterized our detection process used to identify atoms in and , since it has a small but finite error. We have independently determined the probability to misidentify the state of a particle to be , and SOM (). Using these numbers, we extract a corrected probability of preparing states and to be (orange bars in Fig. 2A) and a corrected amplitude of oscillation of (orange points in Fig. 2C). In the population measurement, we independently confirmed that the most commonly observed incorrect patterns are fully consistent with the correct target states with a single detection error (Fig. S5A). From these measurements we extract a lower bound for the atom GHZ state fidelity of , certifying genuine partite entanglement.
This protocol was applied for multiple system sizes of , using s control pulses optimized for each individually. Consistent with expected GHZ dynamics (Fig. 1C and Leibfried2005 ()), the frequency of the measured parity oscillations grows linearly with (Fig. 3A). Extracting the GHZ fidelity from these measurements shows that we surpass the threshold of for all studied system sizes (Fig. 3B). We further characterized the lifetime of the created GHZ state by measuring the parity signal after a variable delay (Fig. 3C). These observations are most consistent with Gaussian decay, while characteristic lifetimes are reduced relatively slowly for increasing system sizes, indicating the presence of a nonMarkovian environment Nielsen2011 (); Monz2011 ().
As an application of our entanglement manipulation technique, we demonstrate its use for entanglement distribution between distant atoms. Specifically, we consider the preparation of Bell states between atoms at the two opposite edges of the array. Our approach is based on first creating the GHZ state using the above procedure, followed by an operation that disentangles all but two target atoms. The latter is realized by shifting the transition frequencies of the two target edge atoms using two strong, bluedetuned addressing beams at nm. Subsequently, we perform a reverse detuning sweep of the Rydberg laser that effectively disentangles all atoms except those at the edges. The resulting state corresponds to a coherent superposition of two pinned excitations that can be converted into a Bell state by applying a resonant pulse on the edge atoms (Fig. 4A).
To demonstrate this protocol experimentally, we prepare a GHZstate of 8 atoms, and turn on the detuned nm addressing beams on the edge atoms resulting in a shift of MHz. We then use an optimized Rydberg laser pulse to distribute the entanglement, and observe the patterns and with a total probability of after accounting for detection errors (Fig. 4B). We verify the coherence of the remote Bell pair by applying an additional pulse with a variable laser phase, and observe parity oscillations with an amplitude of (Fig. 4C). Combining these results, we obtain the edge atom Bell state fidelity of .
Turning to the discussion of our experimental observations, we note that the optimal control provides a substantial improvement over naive analytic pulses (Fig. 3B), while bringing our protocol close to the speed set by a more conventional protocol of building up entanglement through a series of twoqubit operations SOM (). In contrast, a simple linear detuning sweep only allows for the creation of GHZ states for within a fixed s window (Fig. 3B), even under ideal conditions. Our analysis reveals that the reason for this improvement stems from diabatic excitations and deexcitations in the manybody spectrum, related to the recently proposed mechanisms for quantum optimization speedup Farhi2014 (); Zhou2018 (); SOM ().
The measured entanglement fidelity is partially limited by imperfect qubit rotations used for parity measurements. Specifically, the qubit rotation operation in our experiment is induced by an interacting Hamiltonian, which complicates this step. The resulting evolution can be understood in terms of quantum manybody scars Bernien2017 (); Turner2018a (), which gives rise to coherent qubit rotations, even in the presence of strong interactions. The deviations from an ideal parity measurement arises from the Rydberg blockade constraint and longrange interactions SOM (). These grow with the system size, resulting in finite fidelities even for a perfect initial GHZ state (gray shaded area in Fig. 3B). Our quoted fidelity values do not include the correction for this imperfection and represent the lower bound on the actual GHZ state fidelities.
The entanglement generation, manipulation and lifetime are further limited by several sources of decoherence. The finite temperature of the atoms leads to random Doppler shifts on every site as well as position fluctuations that influence interaction energies. These thermal dephasing mechanisms lead to a Gaussian decay of the GHZ state coherence, which decreases with the system size as , in good agreement with our observations (Fig. 3B). Additionally, offresonant laser scattering introduces a small rate of decoherence on each site in the array. We find that numerical simulations of the state preparation accounting for these imperfections predict a higher GHZ fidelity than that obtained experimentally (Fig. 3B and SOM ()). We can attribute this discrepancy to several additional sources of errors. Laser phase noise likely contributes to the finite fidelity of the state preparation. Drifts in the beam positions of the Rydberg lasers can lead to changing light shifts, giving rise to uncontrolled detunings, while drifts in the addressing beam positions can lead to an imbalance in the local energy shifts and thereby in the populations of the two GHZ components, limiting the maximum possible coherence. This analysis highlights the utility of GHZ states for uncovering sources of errors. We emphasize that all of these known error sources can be mitigated via technical improvements SOM ().
Our experiments demonstrate a new promising approach for the deterministic creation and manipulation of largescale entangled states, enabling the certification of genuine partite entanglement in system sizes up to , the largest GHZ state demonstrated to date. These results show the utility of this approach for benchmarking quantum hardware, demonstrating that Rydberg atom arrays constitute a competitive platform for quantum information science and engineering. Specifically, the entanglement generation and distribution could be potentially utilized for applications ranging from quantum metrology and quantum networking to quantum error correction and quantum computation. Our method can be extended by mapping the Rydberg qubit states used here to groundstate hyperfine sublevels, such that the entangled atoms can remain trapped and maintain their quantum coherence over very long times Picken2018 (). This could enable the sophisticated manipulation of entanglement and realization of deep quantum circuits for applications such as quantum optimization Farhi2014 (); Zhou2018 ().
Note added: During the completion of our manuscript we became aware of related work demonstrating GHZ states of 18 superconducting qubits Song2019 (); Wei2019 ().
Acknowledgments
We thank Dries Sels and Christian Reimer for helpful discussions and acknowledge financial support from the Center for Ultracold Atoms, the National Science Foundation, Vannevar Bush Faculty Fellowship, the US Department of Energy and the Office of Naval Research. H.L. acknowledges support from the National Defense Science and Engineering Graduate (NDSEG) fellowship. G.S. acknowledges support from a fellowship from the Max Planck/Harvard Research Center for Quantum Optics. J.C., S.M. and T.C. acknowledge funding from the EC H2020 grants 765267 (QuSCo) and 817482 (PASQuANS), the DFG SPP 1929 (GiRyd) and the IQST Alliance.
References
 (1) D. M. Greenberger, M. A. Horne, A. Zeilinger, Bell’s Theorem, Quantum Theory and Conceptions of the Universe, M. Kafatos, ed., Fundamental Theories of Physics (Springer Netherlands, Dordrecht, 1989), pp. 69–72.
 (2) L. Pezzè, A. Smerzi, M. K. Oberthaler, R. Schmied, P. Treutlein, Rev. Mod. Phys. 90, 035005 (2018).
 (3) M. A. Nielsen, I. L. Chuang, Quantum Computation and Quantum Information: 10th Anniversary Edition (Cambridge University Press, New York, NY, USA, 2011), 10th edn.
 (4) L. Amico, R. Fazio, A. Osterloh, V. Vedral, Rev. Mod. Phys. 80, 517 (2008).
 (5) O. Gühne, G. Tóth, Physics Reports 474, 1 (2009).
 (6) R. Islam, et al., Nature 528, 77 (2015).
 (7) C. A. Sackett, et al., Nature 404, 256 (2000).
 (8) R. Laflamme, E. Knill, W. H. Zurek, P. Catasti, S. V. S. Mariappan, Philosophical Transactions of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences 356, 1941 (1998).
 (9) P. Neumann, et al., Science 320, 1326 (2008).
 (10) D. Bouwmeester, J.W. Pan, M. Daniell, H. Weinfurter, A. Zeilinger, Phys. Rev. Lett. 82, 1345 (1999).
 (11) J.W. Pan, M. Daniell, S. Gasparoni, G. Weihs, A. Zeilinger, Phys. Rev. Lett. 86, 4435 (2001).
 (12) X.L. Wang, et al., Phys. Rev. Lett. 120, 260502 (2018).
 (13) D. Leibfried, et al., Nature 438, 639 (2005).
 (14) T. Monz, et al., Phys. Rev. Lett. 106, 130506 (2011).
 (15) L. DiCarlo, et al., Nature 467, 574 (2010).
 (16) C. Song, et al., Phys. Rev. Lett. 119, 180511 (2017).
 (17) B. Vlastakis, et al., Science 342, 607 (2013).
 (18) H. Labuhn, et al., Nature 534, 667 (2016).
 (19) H. Bernien, et al., Nature 551, 579 (2017).
 (20) D. Jaksch, et al., Phys. Rev. Lett. 85, 2208 (2000).
 (21) T. Wilk, et al., Phys. Rev. Lett. 104, 010502 (2010).
 (22) L. Isenhower, et al., Phys. Rev. Lett. 104, 010503 (2010).
 (23) R. Islam, et al., Science 340, 583 (2013).
 (24) N. Rach, M. M. Müller, T. Calarco, S. Montangero, Phys. Rev. A 92, 062343 (2015).
 (25) R. Heck, et al., PNAS 115, E11231 (2018).
 (26) See Supplementary Materials
 (27) S. de Léséleuc, D. Barredo, V. Lienhard, A. Browaeys, T. Lahaye, Phys. Rev. A 97, 053803 (2018).
 (28) M. Gärttner, et al., Nature Physics 13, 781 (2017).
 (29) E. Farhi, J. Goldstone, S. Gutmann, arXiv:1411.4028 (2014).
 (30) L. Zhou, S.T. Wang, S. Choi, H. Pichler, M. D. Lukin, arXiv:1812.01041 [condmat, physics:quantph] (2018).
 (31) C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn, Z. Papić, Nature Physics 14, 745 (2018).
 (32) C. J. Picken, R. Legaie, K. McDonnell, J. D. Pritchard, Quantum Sci. Technol. 4, 015011 (2018).
 (33) C. Song, et al., arXiv:1905.00320 (2019).
 (34) K. X. Wei, et al., arXiv:1905.05720 (2019).
Appendix A Supplementary Materials
a.1 Experimental setup
The Rydberg excitations are enabled by a twocolor laser system at nm and nm wavelength. The nm light is derived from a frequencydoubled titanium sapphire laser (M Squared SolsTiS 4000 PSX F) locked to an ultrastable reference cavity (by Stable Laser Systems).
The nm light is obtained from a highpower fiber amplifier (ALSIR101510ASP by Azur Light Systems). The seed light is derived from a FabryPérot laser diode injection locked to an external cavity diode laser (CEL002 by MOGLabs) stabilized to the same reference cavity and filtered by the cavity transmission Levine2018 ().
The detuning of both Rydberg lasers to the intermediate state is approximately GHz. The individual Rabi frequencies of the two Rydberg lasers are MHz and MHz. This gives a twophoton Rabi frequency of MHz.
To drive the optimal control pulses, we modulate the nm Rydberg laser with an acoustooptic modulator (AOM) driven by an arbitrary waveform generator (AWG, M4i.6631x8 by Spectrum). We correct the nonlinear response of the AOM to the drive amplitude by a feedforward approach to obtain the target output intensity pattern. Furthermore, the AOM efficiency changes with changing frequency, which we compensate by feeding forward onto the waveform amplitude to suppress the intensity variations with frequency. In addition, the light shift on the Rydberg transition from the nm laser can be as large as MHz. While the pulse intensity changes, this light shift changes, modifying the detuning profile. We therefore correct the frequency profile as a function of the pulse intensity to compensate this shift. These steps ensure that the experimentally applied pulse is a faithful representation of the desired profile.
The local addressing beam patterns are generated by two AODs (DTSX400800 by AA Optoelectronic), each driven by multiple frequencies obtained from an arbitrary waveform generator (M4i.6631x8 by Spectrum).
a.2 Optimal control
Optimal control was originally developed as a tool to harness chemical reactions to obtain the largest amount of desired products with given resources, and then introduced in quantum information processing as a standard way of designing quantum protocols and quantum devices Montangero2007 (); Brif2010 (); Cui2015 (); Koch2016 () as well as in manipulating quantum manybody systems to exploit complex phenomena Doria2011 (); Glaser2015 (); Brouzos2015 (); Lloyd2014 (); vanFrank2016 (); Muller2013 (); Caneva2009 (); Caneva2011a (); Cui2017 (); Heck2018 (). Quantum optimal control theory identifies the optimal shape of a timedependent control pulse to drive a quantum manybody system to accomplish given task, e.g. state preparation or quantum gate implementation. The quality of the transformation is certified by a Figure of Merit (FoM) that can be calculated or measured, e.g. the fidelity of the final state with respect to the target one, the final occupation, or the energy.
In this work, the optimization is achieved through RedCRAB, the remote version of the dressed Chopped RAndom Basis (dCRAB) optimal control via a cloud server Doria2011 (); Caneva2011a (); Heck2018 (). Within the optimization, control fields such as the Rabi coupling are adjusted as , where is an initial guess function obtained from physical intuition or existing suboptimal solutions. The correcting function is expanded by randomized basis functions. In this work, we chose a truncated Fourier basis. Thus, , where are randomized Fourier frequencies with , is the final time, and is a fixed scaling function to keep the values at initial and final times unchanged, i.e., . The optimization task is then translated into a search for the optimal combination of with a given to maximize the fidelity between the target state and the time evolved state at . It can be solved by iteratively updating using a standard NelderMead algorithm Nelder1965 (). In the basic version of the CRAB algorithm, all are fixed and the local control landscape is explored for all frequencies simultaneously. This leads to a restriction in the number of frequencies that can be efficiently optimised. Using the dressed CRAB (dCRAB) algorithm, only one Fourier frequency is optimised at a time. We then move on to after a certain number of iterations of the CRAB routine. This enables the method to include an arbitrarily large number of Fourier components and deriving the solutions without – whenever no other constraints are present – being trapped by local optima Rach2015 ().
In the RedCRAB optimization, the server generates and transmits a trial set of controls to the client user, who will then evaluate the corresponding FoM and communicates the feedback information to the server finishing one iteration loop (Fig. S1). The optimization continues iteratively and the optimal set of controls, as well as the corresponding FoM are derived. In the RedCRAB optimization, the user can either evaluate the FoM by numerical calculation, namely openloop optimization, or by experimental measurement, which is called closedloop optimization. In this work, openloop optimization was carried out only. The resulting controls could later serve as the initial guess for a future closedloop optimization. This last step would ensure that the resulting controls are robust, since all unknown or not modelled experimental defects and perturbations would automatically be corrected for.
For the openloop optimization of the pulse, we constrained the preparation time to s and allowed the detuning to vary between MHz and MHz, while could vary between MHz. The resulting pulses are shown in Fig. S2.
a.3 Optimal control dynamics
To gain insight into the timescales required to prepare a GHZ state in our setup, we can compare our optimal control protocol with a minimal quantum circuit consisting of a series of twoqubit gates that would achieve the same task. In this circuit, a Bell pair is created in the middle of the array using the Rydberg blockade, which for our maximal coupling strength of MHz takes ns. The entanglement can be spread to the two atoms adjacent to this Bell pair by simultaneously applying local pulses of ns to those sites. A sequence of nine of these pulses on the outgoing pairs of atoms leads to the same GHZ state we prepare. This gate sequence requires approximately s, which is within of the total evolution time required in our optimal control sequence, which builds up the entanglement in parallel.
In addition, it is interesting to compare this required evolution time with a naive (unoptimized) parameter ramp that adiabatically connects the initial state to the GHZ state, for example with a linear detuning sweep. Numerical simulations show that such a sweep requires a total evolution time of s to reach similar fidelities as the optimal control pulse for . To understand the origin of the speedup through optimal control, we numerically simulate the corresponding evolution and analyze the population of the instantaneous energy eigenstates (Fig. S3). The optimal control dynamics can be divided into three different regions: (I) A fast initial quench, (II) a slow quench, and (III) a fast final quench. Even though the change in the Hamiltonian parameters in region (I) is rather rapid, the system remains mostly in the instantaneous ground state, with negligible populations of the exited states, since the energy gap is large. In contrast, in region (II) the parameters change slows down, reflecting the fact that the energy gap becomes minimal. Unlike the adiabatic case however, one can observe nontrivial population dynamics, with a temporary population of excited states. Importantly, the optimal control finds a path in the parameter space such that the population is recaptured in the ground state at the end of region (II). This suggests that it actively uses diabatic transitions that go beyond the adiabatic principle. Finally, in region (III) the gap is large again, such that the system can follow the instantaneous ground state even for a fast change of the parameters. This mechanism is related to the recently discussed speedup in the context of the quantum approximate optimization algorithm (QAOA) Farhi2014 (); Zhou2018 ().
a.4 Quantifying detection errors
The manybody dynamics involving coherent excitation to Rydberg states occurs during a fewmicrosecond time window in which the optical tweezers are turned off. After the coherent dynamics, the tweezers are turned back on, and atoms in the ground state are recaptured. However, there is a small but finite chance of losing these atoms. To quantify this error, we perform the GHZ state preparation experiment while disabling the nm Rydberg pulse. This keeps all atoms in state , and we measure the loss probability to find a detection fidelity.
Atoms in state on the other hand have a small chance of being misidentified as being in state , as these atoms can decay prematurely from the Rydberg state to the ground state and get recaptured by the tweezers. This error probability can be measured by preparing atoms at sufficiently large distances as to be noninteracting and applying a calibrated pulse to transfer all atoms to and measure the probability of recapturing them. However, part of this signal is given by the pulse infidelity, i.e. a small fraction of atoms which did not get excited to in the first place.
To quantify the pulse fidelity, we note that a Rydberg atom that decays and is recaptured can decay either into the or ground states with branching ratios and , respectively (). Because our initial optical pumping of atoms into has high fidelity , the final population of atoms should be given only by Rydberg atom decay/recapture events. Following a pulse to excite all atoms to the Rydberg state, the final measured population in is , where is the total decay and recapture probability of a Rydberg atom. Meanwhile, the final measured population in is , which includes both decay events from Rydberg atoms as well as residual population left from an imperfect pulse. Experimentally, we separately measure the total recaptured ground state population , as well as the population only (by a resonant pushout of atoms). We additionally can vary the overall recapture probability by changing the depth of the tweezers that we recapture atoms in, which changes the repulsive force exerted by the optical tweezers on Rydberg atoms deLeseleuc2018 (). We measure and at four different total recapture probabilities to extract the pulse infidelity as (Fig. S4). From these measurements, we conclude a Rydberg detection fidelity of .
Detection errors of can be mitigated by implementing groundstate cooling in the tweezers Kaufman2012a (); Thompson2013a (), which reduces the probability of loss after releasing the atoms. The detection fidelity of can be improved by using Rydberg states with a longer radiative lifetime, actively ionizing the Rydberg atoms by electric or optical fields, or by pulling them away from the trapping region with electric field gradients.
a.5 Correcting for detection errors
The small imperfections in state detection of single qubits leads to a prominent effect on the analysis of large systems. The probability for a single detection error is sufficiently low that multiple errors per chain are very unlikely, and we observe that the reduction in probability of observing the correct GHZ pattern is dominated by these errors, as opposed to excitations of the system (Fig. S5A). This conclusion is further confirmed by noting that nearideal correlations extend across the entire system (Fig. S6).
For the datasets presented in which we correct for these detection errors, we use the following procedure:
Coherences: The coherences are extracted from the amplitude of parity oscillations. Each point in the parity oscillation is analyzed from the measured distribution of the number of excitations in the system. We encode this measured probability distribution in the vector , where is the probability to observe exactly excitations in the system (). The true probability distribution of excitation numbers, prior to the effect of detection errors, is denoted . Detection errors transform this distribution according to a matrix , where encodes the probability that a state with excitations will be detected as having excitations. Each matrix element is calculated using combinatoric arguments from the measured detection fidelities. We determine the true distribution as the one that minimizes the cost function . (Fig. S5B). This procedure is similar to applying the inverse matrix to the measured distribution , but is more robust in the presence of statistical noise on the measured distribution. Error bars on the corrected values are evaluated by random sampling of detection fidelities, given our measured values and uncertainties.
Populations: We carry out a similar procedure to correct the population data; however, we are interested in assessing the probability of two particular target states, which are defined not only by their number of excitations but also by their staggered magnetizations . Our procedure therefore operates by grouping all possible microstates according to their common staggered magnetization and number of excitations (Fig. S5C). For particles, there are in general such groups. As before, we denote the measured distribution with respect to these groups as . We construct a detection error matrix that redistributes populations between groups according to the measured detection error rates. We optimize over all possible true distributions to find the corrected distribution that minimizes the cost function . After this correction, we sum the populations in the two groups that uniquely define the two target GHZ components with a staggered magnetization of , and excitations.
a.6 Bounding the GHZ state coherence
We expand an experimental GHZlike density matrix in the following form
(1)  
where and are the target GHZ components, characterizes the diagonal populations in these states , characterizes the offdiagonal coherence between these states , and contains all other parts of the density matrix. The GHZ fidelity of state is given by:
(2) 
To measure the coherence , we implement a staggered magnetic field to which the target GHZ state is maximally sensitive:
(3) 
Applying to the system for time results in unitary phase accumulation . We then apply a unitary to the system and measure in the computational basis. From repeated measurements, we calculate the expectation value of the global parity operator as a function of the phase accumulation time . Denote the timedependent expectation value , where .
We show that if has a frequency component that oscillates at a frequency of , then the amplitude of this frequency component sets a lower bound for . Importantly, this holds for any unitary used to detect the phase accumulation.
Proof: The expectation value can be written explicitly as the expectation value of the timeevolved observable . In particular,
(4)  
where labels all computational basis states. Since the phase accumulation Hamiltonian is diagonal in the computational basis, the basis states are eigenvectors of with eigenvalues denoting the phase accumulation. Specifically,
(5) 
where is the staggered magnetization of state defined earlier. The staggered magnetization of the state is maximal: , and the staggered magnetization of is minimal: . Note that all other computational basis states have strictly smaller staggered magnetizations. Inserting an identity operator in Eq. (4):
(6)  
The highest frequency component comes from the states with maximally separated staggered magnetization, and . Separating out this frequency component as , we obtain:
(7)  
We note that the parity matrix element is bounded as . Furthermore, the matrix element is realvalued and positive for the unitary considered in the experiment. Fitting to an oscillation with amplitude and phase according to , we produce our lower bound for the offdiagonal coherence :
(8) 
a.7 Parity detection
The ideal observable to measure GHZ phase is the parity . However, the presence of Rydberg interactions and the Rydberg blockade prevents us from rotating all qubits such that we can measure in this basis. Instead, in this work we generate a unitary by resonantly driving all atoms in the presence of these interactions given by for a fixed, optimized time (Fig. S7), and subsequently measure the parity in the computational basis. While it is not obvious that this observable is suitable, we can understand the parity oscillations in the picture of weakly interacting spin particles defined on dimers of neighboring pairs of sites.
For two adjacent sites, we can define eigenstates of a spin1 operator as , , and . In this notation, the antiferromagnetic GHZ state we prepare is given by a ferromagnetic GHZ state in the spin1 basis:
(9) 
We must express all operations on the GHZ state in this new notation. In particular, the transverse field of the form applied to individual atoms gets transformed to an operation on all dimers. Furthermore, the staggered field we apply to individual atoms to rotate the GHZ phase is equivalent to an operation of the form acting on individual dimers.
The parity operator in the singlequbit basis can be transformed into the dimer basis as
(10) 
by noting that the three dimer states are eigenstates of , i.e. and .
Assuming we begin from a GHZ state, applying a rotation on all dimers for a duration given by saturates the difference in between GHZ states of opposite phase. This shows that such a protocol would be optimal if the dimer approximation were exact. However, interactions between dimers cannot be neglected. In particular, the Rydberg blockade suppresses configurations of the form owing to the strong nearestneighbor interaction , and neighboring dimers of the same type such as have a weak interaction given by the nexttonearest neighbor interaction strength . We can thus express the interactions in the system as
(11)  
An exact simulation of the dimer rotation under the interaction Hamiltonian (11) shows that both these interaction effects reduce the parity contrast by a small amount. In the recently discussed context of quantum manybody scars Bernien2017 (); Turner2018a (); Khemani2019 (); Ho2019a (), these effects of residual interactions lead of small deviations from a stable periodic trajectory through phase space.
a.8 Staggered field calibration
To apply the staggered field (3), we address each of the even sites in the array with a focused offresonant laser beam at nm. However, the unitary in question requires a staggered field with opposite sign on every site. We compensate for the missing acquired phase on the sites in between the addressed ones by shifting the phase of the Rydberg laser. The intensity of each addressing beam is measured by applying a spinecho sequence with an addressing pulse of variable duration to determine the light shift on the Rydberg transition. We correct for inhomogeneous intensities so that all atoms are subject to the same light shift.
We measure and calibrate the staggered field by measuring the effect of the field on each atom individually. To do so, we alternately rearrange the atoms to form different subsets of the atom system that are sufficiently far apart to avoid interactions between them. In this configuration, every atom is then subject to a rotation about the axis, followed by the staggered field for variable duration, then a rotation about the axis, to distinguish positive from negative phase evolution. With an additional rotation about the axis, we perform a spin echo to mitigate effects of dephasing. The outcome of this protocol is shown in Fig. S8 and demonstrates the implementation of the staggered magnetic field. By switching the local addressing beams to the opposite set of alternating sites, we switch the sign of the staggered field, enabling the measurement of both positive and negative phase accumulation.
a.9 Measured GHZ fidelities
For each system size , we measure the GHZ populations and the GHZ coherence by parity oscillations (Figs. 2, 3 of the main text). We correct the raw measurements to estimate the true GHZ fidelity before the effect of detection errors. All measured values are shown in table S1. Error bars on raw populations represent a confidence interval for the measured value. Error bars on the raw coherences are fit uncertainties from the parity oscillations. Error bars on the corrected values include propagation of the uncertainty in the estimation of the detection fidelities.
System size  4  8  12  16  20 

Raw populations  0.893(6)  0.797(8)  0.695(9)  0.629(12)  0.585(14) 
Corrected  0.946(10)  0.892(17)  0.824(21)  0.791(29)  0.782(32) 
Raw coherence  0.710(12)  0.516(11)  0.371(10)  0.282(11)  0.211(11) 
Corrected  0.759(11)  0.598(16)  0.462(19)  0.373(19)  0.301(18) 
Raw fidelity  0.801(7)  0.657(7)  0.533(7)  0.455(8)  0.398(9) 
Corrected  0.852(7)  0.745(12)  0.643(14)  0.582(17)  0.542(18) 
a.10 Experimental Imperfections
We identify a number of experimental imperfections that to varying degrees can limit the coherent control of our atomic system.

Atomic temperature: The atom temperature of K leads to fluctuating Doppler shifts in the addressing lasers of order kHz, as well as fluctuations in atomic position that leads to variation in Rydberg interactions strengths. These fluctuations are included in the simulations shown in the main text Figure 3. These effects can be dramatically reduced by improved atomic cooling, most notably by sideband cooling within the optical tweezers to the motional ground state Kaufman2012a (); Thompson2013a ().

Laser scattering: The twophoton excitation scheme to our chosen Rydberg state leads to offresonant scattering from the intermediate state, . This scattering rate has a timescale of s for the two laser fields, and can be reduced by higher laser powers and further detuning from the intermediate state.

Rydberg state lifetime: The Rydberg state has an estimated lifetime of s Beterov2009a (), limited both by radiative decay and blackbodystimulated transitions. This effect could be mitigated by selecting a higher Rydberg state with a longer lifetime or by cryogenic cooling of the blackbody environment.
Additional error sources that may limit our coherence properties include laser phase noise, which can be mitigated by better laser sources and stabilization schemes, and fluctuations in local addressing beam intensities and positions, which can be addressed by active feedback on the beam positions and improved thermal and mechanical stability of the setup.
References
 (35) H. Levine, et al., Phys. Rev. Lett. 121, 123603 (2018).
 (36) S. Montangero, T. Calarco, R. Fazio, Phys. Rev. Lett. 99, 170501 (2007).
 (37) C. Brif, R. Chakrabarti, H. Rabitz, New J. Phys. 12, 075008 (2010).
 (38) J. Cui, F. Mintert, New J. Phys. 17, 093014 (2015).
 (39) C. P. Koch, J. Phys.: Condens. Matter 28, 213001 (2016).
 (40) P. Doria, T. Calarco, S. Montangero, Phys. Rev. Lett. 106, 190501 (2011).
 (41) S. J. Glaser, et al., Eur. Phys. J. D 69, 279 (2015).
 (42) I. Brouzos, et al., Phys. Rev. A 92, 062110 (2015).
 (43) S. Lloyd, S. Montangero, Phys. Rev. Lett. 113, 010502 (2014).
 (44) S. van Frank, et al., Scientific Reports 6, 34187 (2016).
 (45) M. M. Müller, et al., Phys. Rev. A 87, 053412 (2013).
 (46) T. Caneva, et al., Phys. Rev. Lett. 103, 240501 (2009).
 (47) T. Caneva, T. Calarco, R. Fazio, G. E. Santoro, S. Montangero, Phys. Rev. A 84, 012312 (2011).
 (48) J. Cui, R. van Bijnen, T. Pohl, S. Montangero, T. Calarco, Quantum Sci. Technol. 2, 035006 (2017).
 (49) J. A. Nelder, R. Mead, Comput J 7, 308 (1965).
 (50) A. M. Kaufman, B. J. Lester, C. A. Regal, Phys. Rev. X 2, 041014 (2012).
 (51) J. D. Thompson, T. G. Tiecke, A. S. Zibrov, V. Vuletic, M. D. Lukin, Phys. Rev. Lett. 110, 133001 (2013).
 (52) V. Khemani, C. R. Laumann, A. Chandran, Phys. Rev. B 99, 161101 (2019).
 (53) W. W. Ho, S. Choi, H. Pichler, M. D. Lukin, Phys. Rev. Lett. 122, 040603 (2019).
 (54) I. I. Beterov, I. I. Ryabtsev, D. B. Tretyakov, V. M. Entin, Phys. Rev. A 79, 052504 (2009).