Generation and manipulation of Schr\"odinger cat states in Rydberg atom arrays

TL;DR

This paper demonstrates the deterministic creation and manipulation of large-scale Schr"odinger cat states in Rydberg atom arrays, advancing quantum information processing and metrology.

Contribution

It introduces a method to generate and control GHZ-type entangled states of up to 20 qubits in neutral atom arrays using energy spectrum engineering and optimal control.

Findings

Generated up to 20-qubit GHZ states

Achieved entanglement distribution across the array

Demonstrated potential for quantum metrology

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\"odinger cat' states of the Greenberger-Horne-Zeilinger (GHZ) type with up to 20 qubits. Our approach is based on engineering the energy spectrum and using optimal control of the many-body 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.

Tables

Table 1. Table S1: Measured GHZ data for all system sizes. Errors denote 68 % percent 68 68\% confidence intervals.

System size $N$	4	8	12	16	20
Raw populations	0.893(6)	0.797(8)	0.695(9)	0.629(12)	0.585(14)
Inferred	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)
Inferred	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)
Inferred	0.852(7)	0.745(12)	0.643(14)	0.582(17)	0.542(18)

Equations

$$|{\mathrm{GHZ}}_N⟩=\frac{1}{\sqrt{2}}\left(|0101\cdots ⟩+|1010\cdots ⟩\right)$$

$$\frac{H}{\hslash }=\frac{\Omega (t)}{2}\sum _{i=1}^N{\sigma }_x^{(i)}-\sum _{i=1}^N{\Delta }_i(t)n_i+\sum _{i<j}\frac{V}{{\left|i-j\right|}^6}{n}_i{n}_j$$

$$\mathcal{F}=⟨{\mathrm{GHZ}}_N|\rho |{\mathrm{GHZ}}_N⟩=\frac{1}{2}\left(p_{A_N}+p_{{\overline{A}}_N}+c_N+c_N^{\ast }\right)$$

$$D={\left(\frac{ds}{dt}\right)}^2\sum _{n>0}\frac{|⟨E_n(s)|{\partial }_{s}H(s)|E_0(s)⟩{|}^2}{(E_n(s)-E_0(s){)}^2}$$

$$\rho =$$

$$+\left(\beta |A_N⟩⟨{\overline{A}}_N|+{\beta }^{\ast }|{\overline{A}}_N⟩⟨A_N|\right)+{\rho }^{\prime }$$

$$\mathcal{F}=⟨{\mathrm{GHZ}}_N|\rho |{\mathrm{GHZ}}_N⟩=\frac{{\alpha }_1+{\alpha }_2}{2}+\mathrm{Re}(\beta )$$

$$H_{\mathrm{st}}=\frac{\hslash \delta }{2}\sum _{i=1}^N(-1{)}^i{\sigma }_z^{(i)}$$

$$E(T)$$

$$=\sum _n⟨n|\rho U^{†}(T){\mathcal{U}}^{†}\mathcal{PU}U(T)|n⟩$$

$$H_{\mathrm{st}}|n⟩=\frac{\hslash \delta }{2}{M}_n|n⟩\Rightarrow U(T)|n⟩=e^{-i\delta T{M}_n/2}|n⟩$$

$$E(T)=$$

$$=$$

$$F(T)=$$

$$=$$

$$\left|\beta \right|\ge C/2;\arg (\beta )=\varphi$$

$$|{\mathrm{GHZ}}_N⟩=\frac{1}{\sqrt{2}}\left(|+++\cdots ⟩+|---\cdots ⟩\right)$$

$$\mathcal{P}=\prod _j\left(-|+⟩{⟨+|}_j-|-⟩{⟨-|}_j+|0⟩{⟨0|}_j\right)$$

$$\frac{H_{\mathrm{int}}}{\hslash }=$$

$$+V|-⟩{⟨-|}_j|+⟩{⟨+|}_{j+1}$$

Full text

Greenberger-Horne-Zeilinger (GHZ) states constitute an important class of entangled many-body states Greenberger1989 . Such states provide an important resource for applications ranging from quantum metrology Pezze2018 to quantum error correction Nielsen2011 . However, these are among the most fragile many-body states because a single error on any one of the $N$ 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 off-diagonal terms in the $N$-particle density matrix. In contrast to quantifying the degree of entanglement in general many-body states, which is extremely challenging Amico2008 ; Guhne2009 ; Islam2015 , the GHZ state fidelity ($\mathcal{F}>0.5$) constitutes an accessible witness for $N$-partite entanglement Sackett2000 . For these reasons, GHZ state creation can serve as an important benchmark for characterizing the quality of any given quantum hardware. Such states have been previously generated and characterized by using systems of nuclear spins Laflamme1998 ; Neumann2008 , individually controlled optical photons Bouwmeester1999 ; Pan2001 ; Wang2018a , trapped ions Sackett2000 ; Leibfried2005 ; Monz2011 ; Friis2018 , and superconducting quantum circuits DiCarlo2010 ; Song2017a . Large-scale superposition states have also been generated in systems of microwave photons Vlastakis2013 and atomic ensembles without individual particle addressing Pezze2018 .

Here, we demonstrate the preparation of $N$-particle GHZ states

[TABLE]

in a one dimensional array of individually trapped neutral ${}^{87}{\rm Rb}$ atoms, in which the qubits are encoded in an atomic ground state $\left|0\right\rangle$ and a Rydberg state $\left|1\right\rangle$ (phase convention is provided in SOM ). Our entangling operation relies on the strong van-der-Waals interaction between atoms in states $\left|1\right\rangle$ and on engineering the energy spectrum of the quantum many-body system to allow for a robust quantum evolution from an initial product state to a GHZ state. For both generating and characterizing GHZ states (Fig. 1), all the atoms were homogeneously coupled to the Rydberg state $\left|1\right\rangle$ by means of a two-photon transition with an effective coupling strength $\Omega(t)$ and detuning $\Delta(t)$ Labuhn2016 ; Bernien2017 . In addition, we used addressing beams to introduce local energy shifts $\delta_{i}$ on specific sites $i$ along the array (Fig. 1A). The resulting many-body Hamiltonian is

[TABLE]

where $\sigma_{x}^{(i)}=\left|0\right\rangle\!\left\langle 1\right|_{i}+\left|1\right\rangle\!\left\langle 0\right|_{i}$ is the qubit flip operator, $\Delta_{i}(t)=\Delta(t)+\delta_{i}$ is the local effective detuning set by the Rydberg laser and the local light shift, $n_{i}=\left|1\right\rangle\!\left\langle 1\right|_{i}$ is the number of Rydberg excitations on site $i$, and $V$ is the interaction strength of two Rydberg atoms on neighboring sites. The separation between adjacent sites was chosen so that the nearest-neighbor interaction $V=2\pi\cdot 24\,$MHz $\gg\Omega$ results in the Rydberg blockade Jaksch2000 ; Wilk2010 ; Isenhower2010a , forbidding the simultaneous excitation of adjacent atoms into the state $\left|1\right\rangle$.

To prepare GHZ states, we uesd arrays with an even number $N$ of atoms. For large negative detuning $\Delta$ of the Rydberg laser, the many-body ground state of the Hamiltonian (2) is $\left|G_{N}\right\rangle=\left|0000\cdots\right\rangle$. For large uniform positive detuning $\Delta_{i}=\Delta$, the ground-state manifold consists of $N/2+1$ nearly degenerate classical configurations with $N/2$ Rydberg excitations. These include in particular the two target antiferromagnetic configurations $\left|A_{N}\right\rangle=\left|0101\cdots 01\right\rangle$ and $\left|\overline{A}_{N}\right\rangle=\left|1010\cdots 10\right\rangle$ Islam2013 , as well as other states with nearly identical energy (up to a weak second-nearest neighbor interaction), with both edges excited, such as $\left|10010\cdots 01\right\rangle$. To isolate a coherent superposition of states $\left|A_{N}\right\rangle$ and $\left|\overline{A}_{N}\right\rangle$, we introduced local light shifts $\delta_{e}$ using off-resonant laser beams at $840\,$nm, generated with an acousto-optic 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 $\left|G_{N}\right\rangle$ into $\left|{\rm GHZ}_{N}\right\rangle$ by adiabatically increasing $\Delta(t)$ 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 $N$, owing to small energy gaps in the many-body spectrum. To address this limitation, we used 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 chopped-random basis algorithm (RedCRAB) Rach2015 ; Heck2018 , yields optimal shapes of the laser intensity and detuning for the given experimental conditions SOM . For $N\leq 8$ atoms, we performed this optimization using $\delta_{e}/(2\pi)\approx-4.5\,$MHz light shifts on the edge atoms. For larger systems of $N>8$, the preparation was found to be more robust by increasing the edge light shifts to $\delta_{e}/(2\pi)\approx-6\,$MHz and adding $\delta_{4,N-3}/(2\pi)\approx-1.5\,$MHz light shifts on the third site from both edges.

Our experiments are based on the optical tweezer platform and experimental procedure described previously Bernien2017 . After the initialization of a defect-free $N$-atom array, the traps were switched off while the atoms were illuminated with the Rydberg and local light shift beams. The internal state of the atoms is subsequently measured by imaging state $\left|0\right\rangle$ atoms recaptured in the traps, while Rydberg atoms are repelled by the trapping light deLeseleuc2018 . The results of such experiments for a $20$-atom array are demonstrated in Fig. 2. After applying the optimized pulse shown in Fig. 2B, we measured the probability of observing different patterns $p_{n}=\left\langle n\right|\rho\left|n\right\rangle$ in the computational basis, where $\rho$ is the density operator of the prepared state. The measured probability to observe each one of the $2^{20}$ possible patterns in a $20$-atom array is shown in Fig. 2A. The states $\left|A_{20}\right\rangle$ and $\left|\overline{A}_{20}\right\rangle$ clearly stand out (Fig. 2A, blue bars) with a combined probability of $0.585(14)$ and almost equal probability of observing each one.

To characterize the experimentally prepared state $\rho$, we evaluated the GHZ state fidelity

[TABLE]

where $p_{A_{N}}$ and $p_{\overline{A}_{N}}$ are the populations in the target components and $c_{N}=\left\langle\overline{A}_{N}\right|\rho\left|A_{N}\right\rangle$ is the off-diagonal matrix element, which can be measured by using the maximal sensitivity of the GHZ state to a staggered magnetic field. Specifically, evolving the system with the Hamiltonian $H_{p}=\hbar\delta_{p}/2\sum_{i=1}^{N}(-1)^{i}\sigma_{z}^{(i)}$, the amplitude $c_{N}$ acquires a phase $\phi$ at a rate of $\dot{\phi}=N\delta_{p}$. Measuring an observable that oscillates at this frequency provides a lower bound on the coherence $\left|c_{N}\right|$ through the oscillation contrast Garttner2017 ; SOM . In our experiments, the staggered field was implemented by applying off-resonant focused beams of equal intensity at $420\,$nm, generated by another AOD, to every other site of the array (Fig. 1C), resulting in a local energy shift $\delta_{p}$ SOM . Subsequently, we drove the atoms resonantly, applying a unitary operation $\mathcal{U}_{x}$ in order to change the measurement basis SOM , so that a measurement of the parity $\mathcal{P}=\prod_{i}\sigma_{z}^{(i)}$ becomes sensitive to the phase of $c_{N}$. The measured parity is shown in Fig. 2C as a function of the phase accumulated on each atom, 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 $\left|0\right\rangle$ and $\left|1\right\rangle$ because it has a small but finite error. We have independently determined the probability to misidentify the state of a particle to be $p(1|0)=0.0063(1)$, and $p(0|1)=0.0227(42)$ SOM . Subsequently, we use a maximum-likelihood estimation procedure to infer the properties of created states on the basis of the raw measurement results. Using this procedure, we infer a probability of preparing states $\left|A_{20}\right\rangle$ and $\left|\overline{A}_{20}\right\rangle$ to be $0.782(32)$ (Fig. 2A, orange bars) and an amplitude of oscillation of $0.301(18)$ (Fig. 2C, orange points). From these measurements, we extracted a lower bound for the $20$-atom GHZ state fidelity of $\mathcal{F}\geq 0.542(18)$.

This protocol was applied for multiple system sizes of $4\leq N\leq 20$, using $1.1\,\mu$s control pulses optimized for each $N$ individually. Consistent with expected GHZ dynamics (Fig. 1C) Leibfried2005 , the frequency of the measured parity oscillations grows linearly with $N$ (Fig. 3A). Extracting the GHZ fidelity from these measurements shows that we surpass the threshold of $\mathcal{F}=0.5$ for all system sizes studied (Fig. 3B and table S1). 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 non-Markovian 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 was based on first creating the GHZ state by 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 by using two strong, blue-detuned addressing beams at $420\,$nm. Subsequently, we performed 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 $\left|\Phi^{+}\right\rangle=(\left|00\right\rangle+\left|11\right\rangle)/\sqrt{2}$ by applying a resonant $\pi/2$ pulse on the edge atoms (Fig. 4A).

To demonstrate this protocol experimentally, we prepare a GHZ-state of 8 atoms, and turn on the detuned $420\,$nm addressing beams on the edge atoms, resulting in a shift of $\delta_{1,8}/(2\pi)=6\,$MHz. We then used an optimized Rydberg laser pulse to distribute the entanglement and observed the patterns $\left|00000000\right\rangle$ and $\left|10000001\right\rangle$ with a total probability of $0.729(9)$ after accounting for detection errors (Fig. 4B). We verified the coherence of the remote Bell pair by applying an additional $\pi/2$ pulse with a variable laser phase, and observed parity oscillations with an amplitude of $0.481(24)$ (Fig. 4C). Combining these results, we obtained the edge atom Bell state fidelity of $0.605(13)$.

Regarding our experimental observations, the optimal control provides a substantial improvement over naïve 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 two-qubit operations SOM . By contrast, a simple linear detuning sweep only allows for the creation of GHZ states for $N\leq 16$ within a fixed $1.1\,\mu$s window (Fig. 3B), even under ideal conditions. Our analysis reveals that the reason for this improvement stems from diabatic excitations and de-excitations in the many-body 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 $\mathcal{U}_{x}$ in our experiment is induced by an interacting Hamiltonian, which complicates this step SOM . The resulting evolution can be understood in terms of quantum many-body scars Bernien2017 ; Turner2018a , which gives rise to coherent qubit rotations, even in the presence of strong interactions. The deviations from an ideal parity measurement arise from the Rydberg blockade constraint and long-range interactions SOM . These grow with the system size, resulting in finite fidelities even for a perfect initial GHZ state (Fig. 3B, gray shaded area). Our quoted fidelity values do not include the correction for this imperfection and represent the lower bound on the actual GHZ state fidelities.

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, whose time scale decreases with the system size as $1/\sqrt{N}$, which is in good agreement with our observations (Fig. 3C). Additionally, off-resonant 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 higher GHZ fidelities than those obtained experimentally (Fig. 3B) 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, and 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 through technical improvements SOM .

Our experiments demonstrate a new promising approach for the deterministic creation and manipulation of large-scale entangled states, enabling the realization of GHZ-type entanglement in system sizes of up to $N=20$ atoms. 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 potentially be used for applications that range 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 ground-state hyperfine sublevels, so that the entangled atoms can remain trapped and maintain their quantum coherence over very long times Wilk2010 ; Isenhower2010a ; Picken2018 ; SOM . This could enable the sophisticated manipulation of entanglement and realization of deep quantum circuits for applications such as quantum optimization Farhi2014 ; Zhou2018 .

During the completion of this work, we became aware of related results demonstrating large GHZ state preparation using superconducting quantum circuits Song2019 ; Wei2019 .

Acknowledgments

We thank Dries Sels and Christian Reimer for helpful discussions. The authors 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), 817482 (PASQuANS), and QuantERA QTFLAG; the DFG SPP 1929 (GiRyd) and TWITTER; the IQST Alliance; and the Italian PRIN 2017.

I Supplementary Materials

I.1 Experimental setup

The Rydberg excitations are enabled by a two-color laser system at $420\,$nm and $1013\,$nm wavelength. The $420\,$nm light is derived from a frequency-doubled titanium sapphire laser (M Squared SolsTiS 4000 PSX F) locked to an ultrastable reference cavity (by Stable Laser Systems). The $1013\,$nm light is obtained from a high-power fiber amplifier (ALS-IR-1015-10-A-SP by Azur Light Systems). The seed light is derived from a Fabry-Pé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 $\left|6P_{3/2},F=3,m_{F}=-3\right\rangle$ is approximately $2\pi\times 2\,$GHz. The individual Rabi frequencies of the two Rydberg lasers are $\Omega_{\rm 420}/(2\pi)\approx 174\,$MHz and $\Omega_{1013}/(2\pi)\approx 115\,$MHz. This gives a two-photon Rabi frequency of $\Omega=\Omega_{420}\Omega_{1013}/(2\Delta)\approx 2\pi\times 5\,$MHz. We define the local phases of each atom’s states $\left|0\right\rangle$ and $\left|1\right\rangle$ in the reference frame associated with the local phases of Rydberg excitation lasers, such that the two GHZ components have a relative phase $\phi=0$ after state preparation.

To drive the optimal control pulses, we modulate the $420\,$nm Rydberg laser with an acousto-optic modulator (AOM) driven by an arbitrary waveform generator (AWG, M4i.6631-x8 by Spectrum). We correct the nonlinear response of the AOM to the drive amplitude by a feed-forward 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 $420\,$nm laser can be as large as $2\pi\times 4\,$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 (DTSX400-800 by AA Opto-electronic), each driven by multiple frequencies obtained from an arbitrary waveform generator (M4i.6631-x8 by Spectrum).

I.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 many-body 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 time-dependent control pulse to drive a quantum many-body 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 $\Omega(t)$ are adjusted as $\Omega(t)=\Omega_{0}(t)+f(t)$, where $\Omega_{0}(t)$ is an initial guess function obtained from physical intuition or existing suboptimal solutions. The correcting function $f(t)$ is expanded by randomized basis functions. In this work, we chose a truncated Fourier basis. Thus, $f(t)=\Gamma(t)\sum_{k=1}^{n_{c}}[A_{k}\sin(\omega_{k}t)+B_{k}\cos(\omega_{k}t)]$, where $\omega_{k}=2\pi(k+r_{k})/\tau$ are randomized Fourier frequencies with ${r_{k}}\in[-0.5,0.5]$, $\tau$ is the final time, and $\Gamma(t)$ is a fixed scaling function to keep the values at initial and final times unchanged, i.e., $\Gamma(0)=\Gamma(\tau)=0$. The optimization task is then translated into a search for the optimal combination of $\{A_{k},B_{k}\}$ with a given $r_{k}$ to maximize the fidelity between the target state and the time evolved state at $\tau$. It can be solved by iteratively updating $\{A_{k},B_{k}\}$ using a standard Nelder-Mead algorithm Nelder1965 . In the basic version of the CRAB algorithm, all $r_{k}$ are fixed and the local control landscape is explored for all $n_{c}$ 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 $\omega_{k}$ is optimised at a time. We then move on to $\omega_{k+1}$ 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 open-loop optimization, or by experimental measurement, which is called closed-loop optimization. In this work, open-loop optimization was carried out only. The resulting controls could later serve as the initial guess for a future closed-loop 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 open-loop optimization of the pulse, we constrained the preparation time to $1.1\,\mu$s and allowed the detuning $\Delta/(2\pi)$ to vary between $-20\,$MHz and $20\,$MHz, while $\Omega/(2\pi)$ could vary between $0-5\,$MHz. The resulting pulses are shown in Fig. S2. While shorter pulses can work sufficiently well for smaller system sizes, we use an equal pulse duration for all $N$ for better comparability. We find that the optimized pulses for larger systems appear smoother than for smaller system sizes, where the pulses bear less resemblance to an adiabatic protocol. However, the adiabaticity does not improve for larger system sizes, owing to the shrinking energy gaps.

I.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 two-qubit gates that would achieve the same task. In this circuit, a Bell pair is created in the first layer $p=1$ in the middle of the array using the Rydberg blockade, which for our maximal coupling strength of $\Omega/(2\pi)=5\,$MHz takes $100\,$ns$/\sqrt{2}$. The entanglement can be spread to the two atoms adjacent to this Bell pair by simultaneously applying a pair of local $\pi$ pulses of $100\,$ns to those sites, corresponding to controlled rotations. A sequence of such gate layers $p=2,...,10$, including operations on qubit pairs and the free evolution of other qubits, leads to the same GHZ state we prepare. This gate sequence requires approximately $1\,\mu$s, which is within $10\%$ of the total evolution time required in our optimal control sequence, which builds up the entanglement in parallel. Furthermore, the fidelity of each layer of such a circuit effectively acting on all $N=20$ qubits needs to be higher than $0.94$ to achieve the 20-qubit GHZ fidelity demonstrated in this work.

It is interesting to compare this required evolution time with a parameter ramp that tries to adiabatically connect the initial state to the GHZ state. To this end, we parametrize the detuning and Rabi frequency as $\Delta(s)=(1-s)\Delta_{0}+s\Delta_{1}$ and $\Omega(s)=\Omega_{\rm max}[1-\cos^{12}(\pi s)]$ respectively. A naïve (unoptimized) linear ramp of the detuning corresponds to choosing $s=t/T$. Alternatively, one can adjust the local ramp speed to minimize diabatic transitions, for example by choosing $s(t)$ minimizing

[TABLE]

during a ramp of duration $T$. Here $|E_{n}(s)\rangle$ are the instantaneous eigenstates of the Hamiltonian $H(s)$ specified by the parameters $\Omega(s)$ and $\Delta(s)$, with $|E_{0}(s)\rangle$ denoting the instantaneous ground state. In Fig. S3, we show the results of numerical simulations using both the linear sweep and a sweep that minimizes the strength of diabatic processes quantified by $D$. Both sweep profiles require larger total evolution time $T$ than the optimal control pulse to reach similar fidelities.

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. S4). 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 mostly recaptured in the ground state at the end of region (II). Finally, in region (III) the gap is large again and the system parameters are quickly changed to correct also for higher order contributions. This suggests that it actively uses diabatic transitions that go beyond the adiabatic principle. This mechanism is related to the recently discussed speedup in the context of the quantum approximate optimization algorithm (QAOA) Farhi2014 ; Zhou2018 .

II Quantifying detection

The many-body dynamics involving coherent excitation to Rydberg states occurs during a few-microsecond 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 $\left|0\right\rangle$ 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 $420\,$nm Rydberg pulse. This keeps all atoms in state $\left|0\right\rangle$, and we measure the loss probability to find a $0.9937(1)$ detection fidelity.

Atoms in state $\left|1\right\rangle$ on the other hand have a small chance of being misidentified as being in state $\left|0\right\rangle$, 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 non-interacting and applying a calibrated $\pi$ pulse to transfer all atoms to $\left|1\right\rangle$ and measure the probability of recapturing them. However, part of this signal is given by the $\pi$ pulse infidelity, i.e. a small fraction of atoms which did not get excited to $\left|1\right\rangle$ in the first place.

To quantify the $\pi$ pulse fidelity, we note that a Rydberg atom that decays and is recaptured can decay either into the $F=2$ or $F=1$ ground states with branching ratios $\alpha$ and $\beta$, respectively ($\alpha+\beta=1$). Our initial optical pumping of atoms into $\left|0\right\rangle$ has high fidelity $>0.998$, measured using microwave spectroscopy on different sublevels of the $F=2$ manifold. Thus, the final population of $F=1$ atoms should be given only by Rydberg atom decay/recapture events. Following a $\pi$ pulse to excite all atoms to the Rydberg state, the final measured population in $F=1$ is $p_{1}=p\times\beta$, where $p$ is the total decay and recapture probability of a Rydberg atom. Meanwhile, the final measured population in $F=2$ is $p_{2}=p\times\alpha+\epsilon$, which includes both decay events from Rydberg atoms as well as residual population $\epsilon$ left from an imperfect $\pi$ pulse. Experimentally, we separately measure the total recaptured ground state population $(p_{1}+p_{2})$, as well as the $F=1$ population $p_{1}$ only (by a resonant push-out of $F=2$ atoms). We additionally can vary the overall recapture probability $p$ 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 $p_{1}$ and $(p_{1}+p_{2})$ at four different total recapture probabilities to extract the $\pi$ pulse infidelity as $\epsilon=0.006(3)$ (Fig. S5). From these measurements, we conclude a Rydberg detection fidelity of $0.9773(42)$.

Detection errors of $\left|0\right\rangle$ can be mitigated by implementing ground-state cooling in the tweezers Kaufman2012a ; Thompson2013a , which reduces the probability of loss after releasing the atoms. The detection fidelity of $\left|1\right\rangle$ 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.

III Accounting for detection imperfections

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. S6A). This conclusion is further confirmed by noting that near-ideal correlations extend across the entire system (Fig. S7).

To properly infer the obtained fidelities, we account for these imperfections using 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 $\mathbf{W}$, where $W_{n}$ is the probability to observe exactly $n$ excitations in the system ($0\leq n\leq N$). The true probability distribution of excitation numbers, prior to the effect of detection errors, is denoted $\mathbf{V}$. Detection errors transform this distribution according to a matrix $M$, where $M_{mn}$ encodes the probability that a state with $n$ excitations will be detected as having $m$ excitations. Each matrix element is calculated using combinatoric arguments from the measured detection fidelities. We determine the true distribution $\mathbf{V}$ as the one that minimizes the cost function $|M\mathbf{V}-\mathbf{W}|^{2}$. (Fig. S6B). This procedure is similar to applying the inverse matrix $M^{-1}$ to the measured distribution $\mathbf{W}$, but is more robust in the presence of statistical noise on the measured distribution. Error bars on the inferred values are evaluated by random sampling of detection fidelities, given our measured values and uncertainties.

Populations: We carry out a similar procedure for 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 $M_{n}=\sum_{i=1}^{N}(-1)^{i}\langle\sigma_{z}^{(i)}\rangle$. Our procedure therefore operates by grouping all possible microstates according to their common staggered magnetization and number of excitations (Fig. S6C). For $N$ particles, there are in general $(N/2+1)^{2}$ such groups. As before, we denote the raw measured distribution with respect to these groups as $\mathbf{W}$. We construct a detection error matrix $M$ that redistributes populations between groups according to the measured detection error rates. We optimize over all possible true distributions to find the inferred distribution $\mathbf{V}$ that minimizes the cost function $|M\mathbf{V}-\mathbf{W}|^{2}$. Following this procedure, we sum the populations in the two groups that uniquely define the two target GHZ components with a staggered magnetization of $\pm N$, and $N/2$ excitations.

III.1 Bounding the GHZ state coherence

We expand an experimental GHZ-like density matrix in the following form

[TABLE]

where $\left|A_{N}\right\rangle=\left|0101\cdots\right\rangle$ and $\left|\overline{A}_{N}\right\rangle=\left|1010\cdots\right\rangle$ are the target GHZ components, $\alpha_{i}$ characterizes the diagonal populations in these states $(0\leq\alpha_{i}\leq 1)$, $\beta$ characterizes the off-diagonal coherence between these states $(0\leq|\beta|\leq 1/2)$, and $\rho^{\prime}$ contains all other parts of the density matrix. The GHZ fidelity of state $\rho$ is given by:

[TABLE]

To measure the coherence $\left|\beta\right|$, we implement a staggered magnetic field to which the target GHZ state is maximally sensitive:

[TABLE]

Applying $H_{\rm st}$ to the system for time $T$ results in unitary phase accumulation $U(T)=\exp\left(-iH_{\rm st}T/\hbar\right)$. We then apply a unitary $\mathcal{U}$ to the system and measure in the computational basis. From repeated measurements, we calculate the expectation value of the global parity operator $\mathcal{P}=\prod_{i}\sigma_{z}^{(i)}$ as a function of the phase accumulation time $T$. Denote the time-dependent expectation value $E(T)$, where $-1\leq E(T)\leq 1$.

We show that if $E(T)$ has a frequency component that oscillates at a frequency of $N\delta$, then the amplitude of this frequency component sets a lower bound for $|\beta|$. Importantly, this holds for any unitary $\mathcal{U}$ used to detect the phase accumulation.

Proof: The expectation value $E(T)$ can be written explicitly as the expectation value of the time-evolved observable $\mathcal{P}\to U^{\dagger}(T)\mathcal{U^{\dagger}PU}U(T)$. In particular,

[TABLE]

where $\left|n\right\rangle$ labels all computational basis states. Since the phase accumulation Hamiltonian $H_{\rm st}$ is diagonal in the computational basis, the basis states $\left|n\right\rangle$ are eigenvectors of $U(T)$ with eigenvalues denoting the phase accumulation. Specifically,

[TABLE]

where $M_{n}$ is the staggered magnetization of state $\left|n\right\rangle$ defined earlier. The staggered magnetization of the state $\left|A_{N}\right\rangle$ is maximal: $M_{A_{N}}=N$, and the staggered magnetization of $\left|\overline{A}_{N}\right\rangle$ is minimal: $M_{\overline{A}_{N}}=-N$. Note that all other computational basis states have strictly smaller staggered magnetizations. Inserting an identity operator in Eq. (4):

[TABLE]

The highest frequency component comes from the states with maximally separated staggered magnetization, $\left|n\right\rangle=\left|A_{n}\right\rangle$ and $\left|m\right\rangle=\left|\overline{A}_{n}\right\rangle$. Separating out this frequency component as $F(T)$, we obtain:

[TABLE]

We note that the parity matrix element is bounded as $0\leq\left|\langle\overline{A}_{N}|\mathcal{U^{\dagger}PU}|A_{N}\rangle\right|\leq 1$. Furthermore, the matrix element is real-valued and positive for the unitary $\mathcal{U}$ considered in the experiment. Fitting $F(T)$ to an oscillation with amplitude $C\geq 0$ and phase $\phi$ according to $F(T)=C\cos(N\delta T-\phi)$, we produce our lower bound for the off-diagonal coherence $\beta$:

[TABLE]

III.2 Parity detection

The ideal observable to measure GHZ phase is the parity $\mathcal{P}_{x}=\prod_{i}\sigma_{x}^{(i)}$. 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 $\mathcal{U}_{x}=\exp\left(-i\Omega t/2\sum_{i}\sigma_{x}^{(i)}-iH_{\rm int}t/\hbar\right)$ by resonantly driving all atoms in the presence of these interactions given by $H_{\rm int}$ for a fixed, optimized time (Fig. S8), and subsequently measure the parity $\mathcal{P}=\prod_{i}\sigma_{z}^{(i)}$ in the computational basis. The finite duration of the unitary $\mathcal{U}_{x}$ incurs a small amount of additional infidelity, owing both to dephasing and an additional laser scattering. However, we estimate that this effect should only lead to small losses in fidelity on the percent level.

While it is not obvious that the parity observable used here is suitable, we can understand the parity oscillations in the picture of weakly interacting spin-$1$ particles defined on dimers of neighboring pairs of sites. For two adjacent sites, we can define eigenstates of a spin-1 $S_{z}$ operator as $\left|\circ\bullet\right\rangle=\left|-\right\rangle$, $\left|\circ\circ\right\rangle=\left|0\right\rangle$, and $\left|\bullet\circ\right\rangle=\left|+\right\rangle$. In this notation, the antiferromagnetic GHZ state we prepare is given by a ferromagnetic GHZ state in the spin-1 basis:

[TABLE]

We must express all operations on the GHZ state in this new notation. In particular, the transverse field of the form $\hbar\Omega/2\sum_{i}\sigma_{x}^{(i)}$ applied to individual atoms gets transformed to an operation $\hbar\Omega/\sqrt{2}\sum_{j}S_{x}^{(j)}$ on all dimers. Furthermore, the staggered field $\hbar\delta/2\sum_{i}(-1)^{i}\sigma_{z}^{(i)}$ we apply to individual atoms to rotate the GHZ phase is equivalent to an operation of the form $\hbar\delta\sum_{j}S_{z}^{(j)}$ acting on individual dimers.

The parity operator in the single-qubit basis $\mathcal{P}=\prod_{i}\sigma_{z}^{(i)}$ can be transformed into the dimer basis as

[TABLE]

by noting that the three dimer states are eigenstates of $\mathcal{P}$, i.e. $\mathcal{P}\left|\pm\right\rangle=-\left|\pm\right\rangle$ and $\mathcal{P}\left|0\right\rangle=\left|0\right\rangle$.

Assuming we begin from a GHZ state, applying a rotation on all dimers for a duration given by $\Omega t=\pi/\sqrt{2}$ saturates the difference in $\mathcal{P}$ 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 $\left|\cdots-+\cdots\right\rangle$ owing to the strong nearest-neighbor interaction $V$, and neighboring dimers of the same type such as $\left|\cdots\pm\pm\cdots\right\rangle$ have a weak interaction given by the next-to-nearest neighbor interaction strength $V_{2}=V/2^{6}$. We can thus express the interactions in the system as

[TABLE]

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 many-body scars Bernien2017 ; Turner2018a ; Khemani2019 ; Ho2019a , these effects of residual interactions lead of small deviations from a stable periodic trajectory through phase space.

III.3 Staggered field calibration

To apply the staggered field (3), we address each of the even sites in the array with a focused off-resonant laser beam at $420\,$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, through a change in phase of the radio-frequency drive of the AOM. The intensity of each addressing beam is measured by applying a spin-echo 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 $20$-atom system that are sufficiently far apart to avoid interactions between them. In this configuration, every atom is then subject to a $\pi/2$ rotation about the $x$-axis, followed by the staggered field for variable duration, then a $\pi/2$ rotation about the $y$-axis, to distinguish positive from negative phase evolution. With an additional $\pi$ rotation about the $y$-axis, we perform a spin echo to mitigate effects of dephasing. The outcome of this protocol is shown in Fig. S9 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.

III.4 Measured GHZ fidelities

For each system size $N$, we measure the GHZ populations and the GHZ coherence by parity oscillations (Figs. 2, 3 of the main text). From the raw measurements, we infer the true GHZ fidelity using the maximum likelihood procedure discussed in Section “Accounting for detection imperfections”. All measured values are shown in table S1. Error bars on raw populations represent a $68\%$ confidence interval for the measured value. Error bars on the raw coherences are fit uncertainties from the parity oscillations. Error bars on the inferred values include propagation of the uncertainty in the estimation of the detection fidelities.

III.5 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 $\sim 10\,\mu$K leads to fluctuating Doppler shifts in the addressing lasers of order $\sim 2\pi\times 43~{}$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 Fig. 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 . 2. 2.

Laser scattering: The two-photon excitation scheme to our chosen Rydberg state leads to off-resonant scattering from the intermediate state, $6P_{3/2}$. This scattering rate has a timescale of $50-100\mu$s for the two laser fields, and can be reduced by higher laser powers and further detuning from the intermediate state. 3. 3.

Rydberg state lifetime: The $70S$ Rydberg state has an estimated lifetime of $150~{}\mu$s Beterov2009a , limited both by radiative decay and blackbody-stimulated 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. Simulations predict that we could go beyond the system sizes studied here. While GHZ states of $N=24$ could be within reach with current parameters, generation of even larger GHZ states should be feasible with the additional technical improvements discussed above.

IV Ground-state qubit encoding

The GHZ state parity could be more easily detected and manipulated if the qubits were encoded in a basis of hyperfine sublevels of the electronic ground state. In particular, one can consider two alternative qubit states $|\tilde{0}\rangle=\left|5S_{1/2},F=1,m_{F}=-1\right\rangle$ and $|\tilde{1}\rangle=\left|5S_{1/2},F=2,m_{F}=-2\right\rangle$. Rotations between these states are possible through stimulated Raman transitions or microwave driving, and the interactions are introduced by coupling $|\tilde{1}\rangle$ to a Rydberg level $\left|r\right\rangle$. This type of hyperfine encoding has been employed in multiple experiments with cold Rydberg atoms Wilk2010 ; Isenhower2010a ; Picken2018 . To prepare GHZ states in this basis, all atoms can be initialized in $|\tilde{1}\rangle$ and the system transferred to the state $|\tilde{1}r\tilde{1}r\cdots\rangle+|r\tilde{1}r\tilde{1}\cdots\rangle$ using the method described in this work. A ground-state qubit $\pi$ pulse followed by a $\pi$ pulse on the Rydberg transition transforms the state into $|\tilde{0}\tilde{1}\tilde{0}\tilde{1}\cdots\rangle+|\tilde{1}\tilde{0}\tilde{1}\tilde{0}\cdots\rangle$, enabling the long-lived storage of entanglement. Additionally, local qubit rotations can flip the state of every other site to prepare the canonical form of the GHZ state, $|\tilde{0}\tilde{0}\tilde{0}\cdots\rangle+|\tilde{1}\tilde{1}\tilde{1}\cdots\rangle$, which can achieve entanglement-enhanced metrological sensitivity to homogeneous external fields Pezze2018 . Incorporating this type of hyperfine qubit encoding with Rydberg qubit control will be important for realizing quantum gates and deeper quantum circuits in future experiments.