Many-body dynamics of holes in a driven, dissipative spin chain of Rydberg superatoms

TL;DR

This paper investigates the many-body dynamics of holes in a driven, dissipative Rydberg superatom spin chain, revealing a nearly incompressible liquid phase with sub-Poissonian statistics and complex excitation behavior.

Contribution

It introduces an effective many-body model capturing hole mobility and interactions, providing new insights into the steady-state dynamics of Rydberg superatom arrays.

Findings

Steady state Rydberg excitation density is approximately 2/3.

Holes form a nearly incompressible liquid with sub-Poissonian statistics.

The model describes continuous creation, annihilation, and mobility of holes.

Abstract

Strong dipole-dipole interactions between atoms in high-lying Rydberg states can suppress multiple Rydberg excitations within a micron-sized trapping volume and yield sizable Rydberg level shifts at larger distances. Ensembles of atoms in optical microtraps then form Rydberg superatoms with collectively enhanced transition rates to the singly excited state. These superatoms can represent mesoscopic, strongly-interacting spins. We study a regular array of such effective spins driven by a laser field tuned to compensate the interaction-induced level shifts between neighboring superatoms. During the initial transient, a few excited superatoms seed a cascade of resonantly facilitated excitation of large clusters of superatoms. Due to spontaneous decay, the system then relaxes to the steady state having nearly universal Rydberg excitation density $\rho_{\mathrm{R}} = 2/3$. This state is characterized by highly-nontrivial equilibrium dynamics of quasi-particles -- excitation holes in the lattice of Rydberg excited superatoms. We derive an effective many-body model that accounts for hole mobility as well as continuous creation and annihilation of holes upon collisions with each other. We find that holes exhibit a nearly incompressible liquid phase with highly sub-Poissonian number statistics and finite-range density-density correlations.

Equations

$$\mathcal{H}=\sum _k\left[\Omega \left({\hat{\sigma }}_{eg}^k+{\hat{\sigma }}_{ge}^k\right)-\Delta {\hat{\sigma }}_{ee}^k\right]+\sum _{k\ne k^{\prime }}V({\vec{r}}_k,{\vec{r}}_{k^{\prime }}){\hat{\sigma }}_{ee}^k\otimes {\hat{\sigma }}_{ee}^{k^{\prime }},$$

$$\partial_{t}\hat{\rho}=-i[\mathcal{H},\hat{\rho}]+\sum_{k}\sum_{\mu=s,d}\left(L_{\mu}^{k}\hat{\rho}L^{k\dagger}_{\mu}-\mbox{$\frac{1}{2}$}\{L^{k\dagger}_{\mu}L_{\mu}^{k},\hat{\rho}\}\right),$$

$$\mathcal{H}\simeq\sum_{j}\left[\sqrt{N}\Omega\Bigl{(}\ket{G^{j}}\bra{E_{s}^{j}}+\ket{E_{s}^{j}}\bra{G^{j}}\Bigr{)}\right]-\Delta\hat{P}_{EE}^{j}\\ +\sum_{j>i}\frac{C_{p}}{|\vec{r}_{j}-\vec{r}_{i}|^{p}}\hat{P}_{EE}^{j}\otimes\hat{P}_{EE}^{i},$$

$$\frac{\partial }{\partial t}{\rho }_{EE}^j$$

$${\Gamma }_{\mathrm{ex}}({\Delta }_j)$$

$${\Gamma }_{\mathrm{de}}({\Delta }_j)$$

$${\Delta }_j=\Delta -\sum _{i\ne j}\frac{C_p}{|{\vec{r}}_i-{\vec{r}}_j{|}^p}{\rho }_{EE}^i$$

$${\Gamma }_{\uparrow }$$

$${\Gamma }_{\downarrow }$$

$${\Gamma }_{\mathrm{seed}}={\Gamma }_{\mathrm{ex}}({\Delta }_{\mathrm{eff}}=\Delta ).$$

$$g^{(2)}(d_k)\equiv \frac{⟨{\hat{n}}_h^{(j)}{\hat{n}}_h^{(j+k)}⟩}{⟨{\hat{n}}_h^{(j)}{⟩}^2},$$

$$L_c^{(j)}=\sqrt{{\Gamma }_c}{\hat{\sigma }}_{+}^{(j)}[1-{\hat{n}}_h^{(j+1)}][1-{\hat{n}}_h^{(j-1)}],$$

$$L_a^{(j)}$$

$$L_{a_{\pm }}^{(j)}$$

$$L_{t_{\pm }}^{(j)}=\sqrt{{\Gamma }_t}{\hat{\sigma }}_{+}^{(j\pm 1)}{\hat{\sigma }}_{-}^{(j)}[1-{\hat{n}}_h^{(j\pm 2)}].$$

$$Q\equiv \frac{⟨{\hat{n}}_h^2⟩-⟨{\hat{n}}_h{⟩}^2}{⟨{\hat{n}}_h⟩}-1,$$

$${\partial }_t{p}_{\mathrm{g}}$$

$${\partial }_t{p}_{\mathrm{e}}$$

$${\partial }_t{p}_{\mathrm{s}}$$

$${\partial }_t\vec{\rho }=\Lambda \vec{\rho },$$

$$\Lambda =\left(\begin{array}{ccccc}0& {\Gamma }_s& i\sqrt{N}\Omega & -i\sqrt{N}\Omega & {\Gamma }_s\\ 0& -{\Gamma }_s-\frac{N-1}{N}{\Gamma }_d& -i\sqrt{N}\Omega & i\sqrt{N}\Omega & \frac{1}{N}{\Gamma }_d\\ i\sqrt{N}\Omega & -i\sqrt{N}\Omega & -i\Delta -\gamma & 0& 0\\ -i\sqrt{N}\Omega & i\sqrt{N}\Omega & 0& i\Delta -\gamma & 0\\ 0& \frac{N-1}{N}{\Gamma }_d& 0& 0& -{\Gamma }_s-\frac{1}{N}{\Gamma }_d\end{array}\right).$$

$$\begin{gathered} {\rho }_{EE}\equiv {\rho }_{E_s{E}_s}+{\rho }_{E_{ns}{E}_{ns}} \\ =\frac{2N{\Omega }^2\gamma }{{\Omega }^2(N+1)\left({\Gamma }_d+\frac{2N}{N+1}{\Gamma }_s\right)+{\Gamma }_s({\gamma }^2+{\Delta }^2)}. \end{gathered}$$

$${\rho }_{ee}=\frac{{\Omega }^2}{2{\Omega }^2+\frac{{\Gamma }_s}{2\gamma }({\gamma }^2+{\Delta }^2)},$$

$${\rho }_{EE}\simeq \frac{N}{N+1},$$

$$w_{\mathrm{SA}}\simeq \sqrt{2N{\Omega }^2\gamma /{\Gamma }_s+{\gamma }^2}.$$

$${\Lambda }_3=\left(\begin{array}{ccc}-N\chi & N\chi +{\Gamma }_s& {\Gamma }_s\\ N\chi & -N\chi -{\Gamma }_s-\frac{N-1}{N}{\Gamma }_d& \frac{1}{N}{\Gamma }_d\\ 0& \frac{N-1}{N}{\Gamma }_d& -{\Gamma }_s-\frac{1}{N}{\Gamma }_d\end{array}\right),$$

$$\chi \equiv \frac{2{\Omega }^2\gamma }{{\gamma }^2+{\Delta }^2}.$$

$${\Lambda }_2=\left(\begin{array}{cc}-{\Gamma }_{\mathrm{ex}}& {\Gamma }_{\mathrm{de}}\\ {\Gamma }_{\mathrm{ex}}& -{\Gamma }_{\mathrm{de}}\end{array}\right),$$

$${\Gamma }_{\mathrm{ex}}$$

$${\Gamma }_{\mathrm{de}}$$

$$\Delta r\lesssim \frac{1}{p}\frac{r_{\mathrm{fac}}^{p+1}}{C_p/w_{\mathrm{SA}}}=\frac{1}{p}\frac{w_{\mathrm{SA}}}{\Delta }a.$$

$$p\log (2)\ge \log (\Delta /w_{\mathrm{SA}}).$$

$${\partial }_t{p}_{111}=$$

Full text

I Introduction

Strongly-interacting many body systems subject to external driving and coupled to (possibly tailored) reservoirs offer a new route to create and stabilize interesting states of matter. As a simple example, a quantum state can be made immune to particle losses if it is the stationary state of an open system coupled to a particle reservoir. Furthermore, the competition between coherent driving and dissipation can lead to exotic steady states Kapit et al. (2014); Vermersch et al. (2016) and phase transitions in open many body systems (Kessler et al., 2012; Lang and Büchler, 2015; Mendoza-Arenas et al., 2016; Wilson et al., 2016; Marcuzzi et al., 2016; Casteels et al., 2017).

Rydberg atoms are well suited to study the interplay between strong interaction and coupling to coherent laser fields and dissipative environments. They are thus prime candidates to investigate many-body physics of driven, dissipative spin models. A prominent and well studied consequence of the strong, long-range interaction between atoms in Rydberg states is the so called blockade phenomenon, whereby a Rydberg excited atom suppresses further excitations within a certain blockade distance Lukin et al. (2001). Rydberg blockade in a dilute gas or in a regular array of single atoms leads to short-range spatial ordering of excitations, as was studied theoretically Ates and Lesanovsky (2012); Höning et al. (2013); Petrosyan et al. (2013); Petrosyan (2013) and demonstrated experimentally Schauß et al. (2012); Labuhn et al. (2016). In the so-called anti-blockade regime, successive excitation of atoms at a certain distance from each other can be resonantly enhanced Ates et al. (2012); Lesanovsky and Garrahan (2014); Schempp et al. (2014); Urvoy et al. (2015); Marcuzzi et al. (2017); Carr et al. (2013); Malossi et al. (2014), which lead to the lively debate on the possibility of attaining bistable steady states Carr et al. (2013); Malossi et al. (2014); Weimer (2015); Šibalić et al. (2016); Letscher et al. (2017).

When many atoms are confined within the blockade distance from each other, they form an effective two-level system – Rydberg superatom – that can accommodate at most one Rydberg excitation Lukin et al. (2001); Honer et al. (2011); Petrosyan et al. (2013). The coupling of a superatom to the laser radiation is collectively enhanced, while the steady state probability of a single Rydberg excitation can approach unity. This permits the level of control of single collective spins represented by superatoms far exceeding that for individual atoms. Moreover, being composed of many atoms, superatoms are relatively insensitive to atom number fluctuations and losses. Regular arrays of spins represented by superatoms can then be prepared with less experimental effort, which should be contrasted with the sophisticated dynamical preparation techniques used to realize defect-free arrays of individual Rydberg atoms Barredo et al. (2016); Endres et al. (2016). This, together with the strong, long-range interactions between the Rydberg excitations, makes superatoms ideal building blocks for realizing dissipative many-body spin models and analyzing their dynamics.

Single Rydberg superatoms have been observed in several experiments Heidemann et al. (2007); Dudin et al. (2012); Ebert et al. (2015); Weber et al. (2015); Zeiher et al. (2015); Labuhn et al. (2016). A two-dimensional square lattice of superatoms with nearest-neighbor excitation blockade was studied in Hoening et al. (2014), demonstrating the possibility of a phase transition to an anti-ferromagnetic steady state with spontaneously broken lattice symmetry. In the complementary interaction regime of the Rydberg anti-blockade, little is known about the many-body dynamics and the steady state of a lattice of superatoms. Here we study a one-dimensional lattice of Rydberg superatoms [see Fig. 1(a)], in which an already excited superatom facilitates resonant excitation of its neighbor, but the presence of two excited neighbors suppresses the excitation. This systems exhibits interesting excitation dynamics and a highly-nontrivial steady state characterized by an almost universal density $\rho_{\mathrm{R}}=2/3$ of Rydberg excitations with strongly suppressed number fluctuations. We show that this behavior can be explained in terms of mobile excitation holes on the background of Rydberg excited lattice. The holes behave as a nearly incompressible liquid of hard rods with characteristic two-particle correlations [see Fig. 1(b)]. We derive and verify an effective many-body model for holes. Varying the parameters of the effective model, we find a cross-over between a liquid of holes with density-density correlations peaked at the distance of two lattice periods, $2a$, and the onset of crystalline order with period $3a$. In both cases the density of holes is $\rho_{h}=1/3$ with highly suppressed number fluctuations.

The paper is organized as follows. In Sec. II, we formulate the model for a regular array of superatoms and derive the formalism for the efficient treatment of the system. In Sec. III we present the results of numerical simulations of the dynamics of the chain of superatoms (driven spin chain) and introduce the effective hole model that leads to an intuitive physical picture for the equilibrium phase of the system. Section IV summarizes our results. Technical derivations are deferred to the Appendices A, B, and C.

II Chain of Rydberg superatoms

In this Section, and in Appendices A and B, starting from the fully quantum many-body master equation, we derive rate equations for the chain of laser-driven and mutually interacting Rydberg superatoms. These rate equations will then be used in Sec. III for numerical simulations of the dynamics and steady state of the many-body system.

II.1 The microscopic model

We consider an ensemble of cold atoms trapped in a regular array of microtraps Dumke et al. (2002); Würtz et al. (2009); Weber et al. (2015); Ebert et al. (2015); Nogrette et al. (2014) or a long-wavelength optical lattice with the period $a$ of a few microns. Each lattice site $j$ contains on average $N$ atoms, see Fig. 1(a). A laser field of carrier frequency $\omega$ drives the atoms on the transition from the ground state $\ket{g}$ to the excited Rydberg state $\ket{e}$ with the Rabi frequency $\Omega$ and detuning $\Delta=\omega-\omega_{eg}$. In the frame rotating with frequency $\omega$, the coherent excitation dynamics of the atoms is described by the Hamiltonian ($\hbar=1$)

[TABLE]

where $\hat{\sigma}_{\mu\nu}^{k}\equiv\ket{\mu}_{k}\!\bra{\nu}$ is the transition ($\mu\neq\nu$) or projection ($\mu=\nu$) operator for the $k$th atom, and $V(\vec{r}_{k},\vec{r}_{k^{\prime}})=\frac{C_{p}}{|\vec{r}_{k}-\vec{r}_{k^{\prime}}|^{p}}$ is the interaction potential between pairs of atoms at positions $\vec{r}_{k},\vec{r}_{k^{\prime}}$ excited to the Rydberg state $\ket{e}$. The usual van der Waals interaction corresponds to $p=6$.

Atoms excited to the Rydberg state $\ket{e}$ spontaneously decay to the ground state with the rate $\Gamma_{s}$, and are dephased with a typically much larger rate $\Gamma_{d}$ due to collisions of the Rydberg electron with the ground state atoms Balewski et al. (2013); Gaj et al. (2014), atomic motion in the inhomogeneous trapping potential, or intermediate state decay if $\ket{g}\to\ket{e}$ is a two-photon transition Schauß et al. (2012); Schauss et al. (2015). The dissipative dynamics is described by the master equation for the density matrix $\hat{\rho}$ of the system,

[TABLE]

where $L_{\mu}^{k}$ are the Lindblad jump operators for the two relaxation processes assumed acting independently on each atom $k$ as $L_{s}^{k}=\sqrt{\Gamma_{s}}\hat{\sigma}_{ge}^{k}$ and $L_{d}^{k}=\sqrt{\Gamma_{d}}\hat{\sigma}_{ee}^{k}$.

II.2 Rate equations for Rydberg superatoms

In Eq. (1) we can split the sum over all the atoms into two parts: the sum over the lattice sites $j$, and the sum over the atoms $k_{j}$ in each lattice site. We assume that all $N$ atoms within each lattice site are confined within a small spatial interval $\Delta r\ll a$, such that the interatomic interaction energy $C_{p}/(\Delta r)^{p}$ exceeds all the relevant energy scales pertaining to the atoms, namely, the laser Rabi frequency $\Omega$ and detuning $\Delta$, as well as the atomic spontaneous decay $\Gamma_{s}$ and dephasing $\Gamma_{d}$ rates and the resulting Rydberg-state excitation linewidth is $w\simeq 2\Omega\sqrt{\gamma/\Gamma_{s}}$ with $\gamma\equiv\frac{1}{2}(\Gamma_{d}+\Gamma_{s})$ Petrosyan et al. (2013). This allows us to neglect all the multi-atom states containing more than one Rydberg excitation per lattice site Lukin et al. (2001); Petrosyan et al. (2013). If on the spatial scale $\Delta r\lesssim 1\>\mu$m of such a Rydberg superatom the laser field can be assumed uniform, it would couple the collective ground state $\ket{G}\equiv\ket{g_{1}g_{2}\ldots g_{N}}$ only to the symmetric single excitation state $\ket{E_{s}}\equiv\frac{1}{\sqrt{N}}\sum_{k}\ket{g_{1}g_{2}\ldots e_{k}\ldots g_{N}}$ with the collectively enhanced Rabi frequency $\sqrt{N}\Omega$, see Fig. 2(a). There are, in addition, $(N-1)$ nonsymmetric single excitation states $\ket{E_{ns}}_{m}$, labeled by index $m$, decoupled from the laser field. We can then recast the Hamiltonian as

[TABLE]

where $\hat{P}_{EE}^{j}\equiv\ket{E_{s}^{j}}\bra{E_{s}^{j}}+\sum_{m}\ket{E_{ns}^{j}}_{m}\!\bra{E_{ns}^{j}}$ is the projector onto the manifold of $N$ single Rydberg excitation states of superatom at site $j$, and $\vec{r}_{j}=\frac{1}{N}\sum_{k_{j}}\vec{r}_{k_{j}}$ is its center of mass coordinate.

Within each Rydberg superatom, the dephasing couples incoherently the symmetric single excitation state $\ket{E_{s}}$ to all $(N-1)$ nonsymmetric states $\ket{E_{ns}}$ with the rate $\Gamma_{d}(N-1)/N$ which approaches $\Gamma_{d}$ for $N\gg 1$. In turn, the reverse coupling of the nonsymmetric states to the symmetric state has the rate $\Gamma_{d}/N$. All single excitation states decay back to the ground state $\ket{G}$ with rate $\Gamma_{s}$, see Fig. 2(a). In Appendix A we derive a simple rate equation model describing the excitation dynamics of a superatom in the limit of strong dephasing $\Gamma_{d}\gtrsim\Omega$. Starting from the density matrix equations for a single Rydberg superatom, we adiabatically eliminate all coherences and the population of the symmetric excited state, which scales as $\sim 1/N$. The superatom then reduces to an effective two level system, see Fig. 2(b), and its dynamics is governed by the rate equations for the populations of the ground $\rho_{GG}^{j}$ and excited $\rho_{EE}^{j}$ states,

[TABLE]

where the excitation and de-excitation rates are given by

[TABLE]

with $\chi(\Delta_{j})\equiv\frac{2\Omega^{2}\gamma}{\gamma^{2}+(\Delta_{j})^{2}}$. Here

[TABLE]

is the effective detuning of superatom $j$ which includes the Rydberg level shift due to the interaction with all the other superatoms in the Rydberg state. In Fig. 2(c), and in more detail in Appendix A, we compare the dynamics of a single superatom as obtained from the exact solution of the complete set of the density matrix equations and the solution of the rate equations. We observe that the rate equation model approximates well the relaxation timescale of a superatom and the steady state population of the excited state, $\rho_{EE}=\frac{\Gamma_{\mathrm{ex}}}{\Gamma_{\mathrm{ex}}+\Gamma_{\mathrm{de}}}$. Remarkably, unlike for a single two-level atom, the excitation probability of the superatom under continuous (near-)resonant driving and in the presence of strong dephasing $\gamma\gtrsim\Omega$ can approach unity, $\rho_{EE}\simeq\frac{N}{N+1}\to 1$, with increasing $N$.

II.3 Facilitated excitation of superatoms

The laser irradiates continuously and uniformly the 1D chain of superatoms. We set the detuning $\Delta$ of the laser to be equal to the interaction strength $V(a)=C_{p}/a^{p}\equiv V_{\mathrm{N}}$ between neighboring superatoms, $\Delta=V_{\mathrm{N}}$. We neglect the interaction between the next to nearest neighbors. The conditions for the validity of our treatment of the chain of superatoms are discussed in Appendix B. In this so-called facilitation regime Lesanovsky and Garrahan (2014), an already excited superatom shifts the frequency of its nearest neighbor into resonance with the laser. The excitation $\Gamma_{\uparrow}$ and de-excitation $\Gamma_{\downarrow}$ rates for the facilitated superatom, having one and only one, excited neighbor, are

[TABLE]

On the other hand, a superatom surrounded by non-excited neighbors is non-resonant with the laser and has a much lower excitation rate

[TABLE]

But once excited, it will play the role of a seed for a rapid growth of a cluster of excited superatoms. Naively, such clusters will grow until they collide and nearly all of the superatoms in a lattice will be excited. However, spontaneous decay of superatoms with rate $\Gamma_{s}$ will produce excitation holes – ground state superatoms surrounded by excited superatoms. A hole cannot be resonantly excited as its Rydberg state is shifted by the interaction with two excited neighbors by $2V_{\mathrm{N}}$, leading to the effective detuning $\Delta_{\mathrm{eff}}=\Delta-2V_{\mathrm{N}}=-\Delta$. Hence, the rate to refill the hole turns out to be the same highly-suppressed seed rate $\Gamma_{\mathrm{seed}}$.

In the facilitation regime that we consider, the typical hierarchy of the relevant rates is $\Gamma_{\mathrm{seed}}\ll\Gamma_{s},\Gamma_{\mathrm{\downarrow}}\ll\Gamma_{\mathrm{\uparrow}}$.

III The many-body dynamics

Upon turning on the excitation laser, the chain of superatoms under the facilitation conditions described above will exhibit transient dynamics on the timescale of $t\sim\Gamma_{s}^{-1}$ and then settle in a steady state. The excitation transients will be described later in this Section. First we discuss the steady-state characterized by the average Rydberg excitation density of $\rho_{\mathrm{R}}\simeq 2/3$ and highly nontrivial equilibrium dynamics of the excitation holes.

III.1 Steady-state distribution of holes

By definition, a hole is a ground state superatom surrounded by two excited superatoms. Holes originate from collisions of growing clusters of Rydberg excitations during the transient and spontaneous decay of excited superatoms inside the cluster. In Fig. 3 we show the density of holes in a long lattice of superatoms. After switching on the excitation laser, within a few lifetimes $\Gamma_{s}^{-1}$ of Rydberg excitations, the density of holes reaches an equilibrium. For large enough values of the interaction strength $V_{\mathrm{N}}\gg w$, and thereby the laser detuning $(\Delta=V_{\mathrm{N}})$, the steady-state density of holes approaches the value of $\rho_{h}=1/3$.

In Fig. 4 we show the spatial correlation function $g^{(2)}(d)$ for holes in the steady state. The two-particle correlation function is defined via

[TABLE]

where $\hat{n}_{h}^{(j)}$ is the hole number operator for site $j$ of the lattice, $d_{k}=ak$ with $k\in\mathbb{N}$ is the distance, and we assume that the average hole density $\langle\hat{n}_{h}^{(j)}\rangle=\rho_{h}$ is spatially uniform. We observe that, to a good approximation, holes behave as hard rods of length $2a$, with the average density of rods being $1/3$.

III.2 Effective model for holes

To understand the results of the numerical simulations for the hole density and correlation function, we have derived an effective model for the equilibrium dynamics of the holes. The derivation, details of which are given in Appendix C, is based on adiabatic elimination of short-lived configurations involving two or more neighboring superatoms in the ground state. Such configurations appear when an excited superatom next to a hole decays to the ground state with rate $\Gamma_{s}$, but the lifetime of these configurations is very short, $t\sim\Gamma_{\uparrow}^{-1}$, due to the excitation facilitation with the fast rate $\Gamma_{\uparrow}$. Hence, in the effective model, holes cannot be located on the neighboring sites, which gives the physically intuitive picture as to why they behave as hard rods of length $2a$.

There are three fundamental processes affecting the many-body dynamics of the holes on the lattice, as illustrated in the inset of Fig. 3. We describe these processes in terms of the Lindblad jump operators acting in the corresponding subspace for the holes.

(i) Holes are created with rate $\Gamma_{c}=\Gamma_{s}$. The corresponding Lindblad operator is given by

[TABLE]

where $\hat{\sigma}_{\pm}^{(j)}$ is the hole creation/annihilation operator for site $j$ of the lattice, and $\hat{n}_{h}^{(j)}\equiv\hat{\sigma}_{+}^{(j)}\hat{\sigma}_{-}^{(j)}$ is the number operator. The last two terms on the right hand side ensure that a hole cannot be created next to an existing one.

(ii) When two holes are separated by one excitation, one of the holes can be annihilated in one of the following three ways described by,

[TABLE]

The total annihilation probability is $\Gamma_{a}=\Gamma_{s}$. The remaining hole can then occupy either the middle site with half of the total probability, or one of the side sites, each with quarter of the total probability.

(iii) Holes can hop from site to site with transport rate $\Gamma_{t}=\Gamma_{\mathrm{\downarrow}}/2$,

[TABLE]

Here the last term ensured that the hole cannot hop to a site next to an existing hole.

In Fig. 4 we compare the spatial correlation function $g^{(2)}(d)$ for ground state superatoms obtained from the numerical simulations of the full superatom model and the effective hole model. We observe very good agreement between the full and effective models, including $g^{(2)}(1)\simeq 0$.

III.3 Liquid-crystal crossover for lattice holes

Although in the present setup the hole hopping rate cannot be made arbitrary small, $\Gamma_{t}\gtrsim\Gamma_{c,a}$, it is instructive to analyze how varying $\Gamma_{t}$ would affect the many-body steady state.

Consider first the hypothetical case of no hole transport, $\Gamma_{t}\to 0$. We can also neglect for now the probability of refilling the hole, due to smallness of the corresponding rate $\Gamma_{\mathrm{seed}}$. Then the only stable configuration corresponds to holes on every third site of the lattice, since neither the hole creation nor annihilation processes of Eqs. (9) and (10) affect the system. Due to translational invariance of the lattice, this configuration is triple degenerate. The hole density-density correlation function $g^{(2)}(d)$ is peaked at $d=3a$, and the steady state approaches a pure state with long-range crystalline order. Since the density of holes is exactly $\rho_{h}=1/3$, their number does not fluctuate. We can characterize the number fluctuations of the holes by the Mandel $Q$ parameter

[TABLE]

where $\hat{n}_{h}\equiv\sum_{j}\hat{n}_{h}^{(j)}$ is the total number of holes. $Q<0$ signifies sub-Poissonian number distribution, with $Q=-1$ corresponding to a precise number of holes with no fluctuations.

As we now increase $\Gamma_{t}$, the holes become mobile and the peak of the correlation function $g^{(2)}(d)$ at $d=3a$ is progressively reduced, see Fig. 5. The holes can now approach each other and annihilate, followed by hole creation on the allowed sites, which causes finite fluctuation of the hole number. Yet, their statistics remains highly sub-Poissonian, $Q\simeq-0.81$, even for large $\Gamma_{t}\gg\Gamma_{c,a}$. Since the hole creation and annihilation rates are the same, $\Gamma_{c}=\Gamma_{a}=\Gamma_{s}$, their mean density stays close to $\rho_{h}=1/3$. Interestingly, the crystal with periodicity $d=3a$ does not simply melt into a liquid with the same period. Rather, for a large hopping rate, the correlation function $g^{(2)}(d)$ exhibit period $d=2a$ with short correlation length $\xi\lesssim 1a$, see Fig. 4. Since holes cannot come closer than two lattice sites, they start to behave as mobile hard rods of length $2a$.

We finally note that when the interaction strength $V_{\mathrm{N}}$, and thereby the laser detuning $\Delta=V_{\mathrm{N}}$, are not sufficiently larger than the excitation linewidth, the seed rate $\Gamma_{\mathrm{seed}}$ is not negligible and holes can refill, leading to $\rho_{h}<1/3$, as can be seen in Fig. 3.

III.4 Transient excitation dynamics of the system

As promised above, we now discuss the dynamics of the system of superatoms initially in the ground state upon switching on the excitation laser. In Fig. 6 we show the density $\rho_{\mathrm{R}}$ of Rydberg excitations as a function of time, for different values of the interaction strength $V_{\mathrm{N}}$ (and laser detuning $\Delta=V_{\mathrm{N}}$). In the long time limit $t\gg\Gamma_{s}^{-1}$, for large enough $V_{\mathrm{N}}\gg w$, the system reaches the steady state with the Rydberg excitation density $\rho_{\mathrm{R}}\simeq 2/3$, consistent with the hole density $\rho_{h}\simeq 1/3$ analyzed above. For smaller values of $V_{\mathrm{N}}$, we have in the steady state $\rho_{\mathrm{R}}>2/3$, since holes can be refilled with appreciable seed rate $\Gamma_{\mathrm{seed}}$.

Consider now the transient regime $0<t\lesssim\Gamma_{s}^{-1}$. As seen in Fig. 6, after switching on the laser, the density of Rydberg excitations can peak at a large value and then relax to the lower steady-state density. At first sight, this observation is surprising, since such a behavior is reminiscent to partially coherent dynamics of a quantum system, such as, e.g., damped Rabi oscillations, while here our system is completely governed by rate equations and no coherences are involved. We now outline a macroscopic model that will explain the nature of the peak and the associated peak time $t_{0}$.

We consider three basic states of the system: The ground state (g), the fully excited state (e) and the final steady state (s). The corresponding probabilities are denoted by $p_{\mathrm{g}}$, $p_{\mathrm{e}}$, and $p_{\mathrm{s}}$. Initially the system is in the ground state, $p_{\mathrm{g}}=1$, and to start the dynamics we need a seed excitation. The probability of the seed $p_{\mathrm{seed}}$ is governed by the equation $\partial_{t}p_{\mathrm{seed}}(t)=p_{\mathrm{g}}\Gamma_{\mathrm{seed}}$. For short times, we can assume $p_{\mathrm{g}}\simeq 1$ and obtain linear growth of the seed probability, $p_{\mathrm{seed}}(t)\simeq p_{\mathrm{seed}}(0)+\Gamma_{\mathrm{seed}}t$. For longer times, this approximation breaks down, but once $p_{\mathrm{g}}$ is depleted, the role of the seed becomes unimportant, as will become clear shortly.

The equations of motion for the probabilities of the three basic states are

[TABLE]

Each seed excitation in the lattice triggers fast growth of facilitated excitations in both lattice directions. The ground state is then being depleted with the rate $2\Gamma_{\uparrow}$ and its population is transferred into the fully excited state. In turn, the Rydberg excitations in the fully excited state decay with rate $\Gamma_{s}$. Since the steady state corresponds to configurations with, on average, every third site non-excited, the decay rate of $p_{\mathrm{e}}$ to $p_{\mathrm{s}}$ is $3\Gamma_{s}$, which we have also verified via numerical simulations of the full system with all the superatoms initially excited. We also include in the above equations the process of refilling the holes with the corresponding rate $\Gamma_{\mathrm{seed}}$. Then the solution of these equations approximates remarkably well the mean density of Rydberg excitations \rho_{\mathrm{R}}(t)=p_{\mathrm{e}}(t)+\mbox{\frac{2}{3}}p_{\mathrm{s}}(t), even in the regime of sizable $\Gamma_{\mathrm{seed}}$, as can be seen in the main panel of Fig. 6.

The peak in the excitation density $\rho_{\mathrm{R}}$ is reached when the ground state probability $p_{\mathrm{g}}$ is depleted and the majority of superatoms are transferred to the excited state, $p_{\mathrm{e}}\simeq 1$. Integration of $\partial_{t}p_{\mathrm{g}}$ suggests the peak time scaling as $t_{0}\propto 1/\sqrt{\Gamma_{\mathrm{seed}}\Gamma_{\uparrow}}\propto V_{\mathrm{N}}$. The scaling of the peak time with the interaction strength $V_{\mathrm{N}}/w$ is verified in the inset of Fig. 6. Note finally that our macrocsopic model neglects cluster collisions, also producing holes, and thereby slightly overestimates the excitation density $\rho_{R}$ obtained from the simulations of the complete microscopic model.

IV Conclusions and outlook

To summarize, we have studied the excitation dynamics and steady state of a lattice of Rydberg superatoms driven by a laser in the facilitation regime. We have shown that the steady state of the system has nearly universal Rydberg excitation density of $\rho_{\mathrm{R}}=2/3$. More interestingly, it corresponds to an equilibrium dynamics of mobile quasi-particles – excitation holes. We have derived an effective hole model which involves hole creation and pair annihilation of holes separated by two lattice sites. We have found that the number fluctuations of the holes are characterized by the Mandel $Q$ parameter $Q\simeq-0.81$, and their spatial correlations decay on a distance of $\xi\lesssim 1a$ comparable to the lattice constant $a$.

That negative values of $Q$, which signify pronounced sub-Poissonian number distribution, have their origin in the hard rod constraint, has been pointed out already in Ates and Lesanovsky (2012) considering a Gibbs ensemble of Rydberg excitations. In our system, however, the pair annihilation of holes leads to a state with much stronger suppressed fluctuations. In the Gibbs state of Ref. Ates and Lesanovsky (2012), we would have to choose a rod length of about $d=0.611$ to obtain the average hole density $\rho_{h}=1/3$, but this would then lead to $Q\simeq-0.63$ which is significantly larger than what we obtain from our simulations.

Our model corresponds to an experimentally amenable regime with a large hole transport rate $\Gamma_{t}\sim\Gamma_{c,a}$ comparable to the hole creation and annihilation rates $\Gamma_{c,a}=\Gamma_{s}$, which, in turn, are determined by the spontaneous decay rate of Rydberg excitations. The resulting spatial correlations of the holes have a period of $d=2a$ and short correlation length $\xi\simeq a$. If we freeze the hole motion, $\Gamma_{t}\ll\Gamma_{s}$, we would obtain long-range order, $\xi\gg a$, with $d=3a$ periodicity. Such a regime can in principle be achieved, but at the expense of more complicated experimental setup involving additional lasers coupling the ground state of superatoms to a different long lived Rydberg or metastable state which would make the holes immobile.

Another interesting direction of research is to consider different lattice geometries, e.g., a two-dimensional square or triangular lattice. The latter might simulate dissipative frustrated spin models. Finally, the rate equation approach, amenable to large scale numerical simulations, is applicable in the regime of strong dephasing. Coherence effects might lead to interesting dynamics and yield long range correlations and entanglement. However, fully quantum many-body simulations are limited to small system sizes.

Acknowledgements.

This work was supported by DFG through SFB/TR185 (M.F. and F.L.) and the H2020 FET Proactive project RySQ (D.P.). F.L. is supported by a fellowship through the Excellence Initiative MAINZ (DFG/GSC 266). D.P. is grateful to the University of Kaiserslautern for hospitality and to the Alexander von Humboldt Foundation for support during his stay in Germany.

Appendix A Rate equations model for a superatom

Here we consider in some detail the dissipative dynamics of a single Rydberg superatom. The superatom consists of $N$ two-level atoms within a volume of linear dimension $\Delta r$ smaller than the Rydberg blockade distance $a_{\mathrm{B}}$. We define $a_{\mathrm{B}}$ as the distance below which the interatomic interaction strength $V(r\leq a_{\mathrm{B}})\geq w$ starts to exceed the steady-state excitation linewidth of an atom, $w=\sqrt{4\Omega^{2}\gamma/\Gamma_{s}+\gamma^{2}}\simeq 2\Omega\sqrt{\gamma/\Gamma_{s}}$ Petrosyan et al. (2013); Petrosyan (2013). Here $\gamma\equiv\frac{1}{2}(\Gamma_{d}+\Gamma_{s})$, and we assume $\Omega^{2}\gg\Gamma_{s}\Gamma_{d}$. Due to the strong Rydberg blockade, the superatom can accommodate at most one Rydberg excitation.

The collective ground state of a superatom $\ket{G}\equiv\ket{g_{1}g_{2}\ldots g_{N}}$ is coherently coupled by the laser to the symmetric single excitation state $\ket{E_{s}}\equiv\frac{1}{\sqrt{N}}\sum_{k}\ket{g_{1}g_{2}\ldots e_{k}\ldots g_{N}}$ with Rabi frequency $\sqrt{N}\Omega$, see Fig. 2(a). In addition, there are $N-1$ non-symmetric single excitation states $\{\ket{E_{ns}}_{m}\}$ which are not directly coupled to the ground state by the laser. All the excited states $\{\ket{E}\}$ spontaneously decay to the ground state $\ket{G}$ with rate $\Gamma_{s}$. The dephasing $\Gamma_{d}$ of the atomic Rydberg state $\ket{e}$, with respect to the ground state $\ket{g}$, leads to incoherent coupling of any single excitation state $\ket{E}$ of the superatom to any other such state $\ket{E^{\prime}}$ with rate $\Gamma_{d}/N$. We may replace the manifold $\{\ket{E_{ns}}_{m}\}$ with a single aggregate non-symmetric state $\ket{E_{ns}}$, obtaining for the superatom an effective three-level system, $\{\ket{G},\ket{E_{s}},\ket{E_{ns}}\}$. We describe this dissipative system with a “vector” of density matrix elements $\vec{\rho}=(\rho_{GG},\rho_{E_{s}E_{s}},\rho_{GE_{s}},\rho_{E_{s}G},\rho_{E_{ns}E_{ns}})^{T}$ which obeys the equation of motion

[TABLE]

where

[TABLE]

For $N=1$, we retrieve the well known optical Bloch equations for a single two level atom, with the decay rate $\Gamma_{s}$ of the excited state population $\rho_{ee}$ and the relaxation rate $\gamma=\frac{1}{2}(\Gamma_{s}+\Gamma_{d})$ of the coherence $\rho_{ge}$. For a superatom with $N>1$, dephasing induces population $\rho_{E_{ns}E_{ns}}>0$ of non-symmetric state(s), but their coherences remain decoupled, $\rho_{E_{ns}G}=\rho_{E_{ns}E_{s}}=0$. By definition, the superatom can contain at most one Rydberg excitation, $\rho_{GG}+\rho_{E_{s}E_{s}}+\rho_{E_{ns}E_{ns}}=1$.

Solving Eq. (13) in the steady state, $\partial_{t}\vec{\rho}=0$, we obtain the following expression for the total excited state population of the superatom,

[TABLE]

For $N=1$, this reduces to the excited state population of a single two level atom Petrosyan et al. (2013),

[TABLE]

which is bounded by $\rho_{ee}\leq 0.5$ even for resonant excitation $|\Delta|\ll\gamma$. On the other hand, for large $N$ and small decay rate $\Gamma_{s}\ll\Gamma_{d},\Omega^{2}/\Gamma_{d}$, the excited state population of the superatom is approximately given by

[TABLE]

which quickly approaches $\rho_{EE}\simeq 1$ with increasing $N$. We emphasize that the steady-state population inversion of the superatom is brought about by strong driving and dephasing, which tend to equalize populations of all $N+1$ states, i.e., $N$ single excitation states and the ground state. In contrast, without dephasing ($\Gamma_{d}\ll\Gamma_{s}$), the superatom reduces to a two-level system, with the ground state $\ket{G}$ coherently coupled to the symmetric excited state $\ket{E_{s}}$ which decays back to the ground state with rate $\Gamma_{s}$, and the resulting steady-state Rydberg excitation probability is $\rho_{EE}\simeq\rho_{E_{s}E_{s}}\leq 1/2$.

From Eq. (15), with strong dephasing $\Gamma_{d}\gg\Gamma_{s}$ and $N\gg 1$, we obtain the excitation linewidth of the superatom

[TABLE]

Comparing it with the single-atom excitation linewidth $w$ which follows from Eq. (16), we notice the analogy upon replacement $\Omega^{2}\to\frac{1}{2}N\Omega^{2}$, rather than $\Omega^{2}\to N\Omega^{2}$ as one would naively expect from the collective enhancement of the $N$-atom Rabi frequency $\sqrt{N}\Omega$. This factor of $2$ difference stems from the fact that the Rydberg excitation probability of a saturated superatom approaches unity, rather than $\frac{1}{2}$.

For strong dephasing $\Gamma_{d}\gtrsim\Omega$, the dynamics of a superatom can be described, to a good approximation, by rate equations which we now derive. In the equations for the density matrix elements, we adiabatically eliminate the coherences $\rho_{E_{s}G}$, $\rho_{GE_{s}}$, obtaining rate equations $\partial_{t}\vec{\rho}=\Lambda_{3}\vec{\rho}$ for populations $\vec{\rho}=(\rho_{GG},\rho_{E_{s}E_{s}},\rho_{E_{ns}E_{ns}})^{T}$ with

[TABLE]

where

[TABLE]

The laser tends to equalize populations $\rho_{GG}$ and $\rho_{E_{s}E_{s}}$ of the ground and symmetric excited states. In turn, strong dephasing quickly transfers the population of the symmetric state to the non-symmetric states. For large $N$, the reverse transfer from $\rho_{E_{ns}E_{ns}}$ to $\rho_{E_{s}E_{s}}$ is suppressed by a factor of $1/N$. The symmetric state then plays the role of an intermediate state having small population $\rho_{E_{s}E_{s}}\sim\rho_{E_{ns}E_{ns}}/N$. Upon adiabatic elimination of $\rho_{E_{s}E_{s}}$ we finally obtain rate equations $\partial_{t}\vec{\rho}=\Lambda_{2}\vec{\rho}$ for the populations $\vec{\rho}=(\rho_{GG},\rho_{EE})^{T}$ of the ground and excited states ($\rho_{GG}+\rho_{EE}\simeq 1$),

[TABLE]

where the excitation and de-excitation rates are given by

[TABLE]

This is illustrated in Fig. 2(b).

In Fig. 7(a) we compare the excitation dynamics of a resonantly driven superatom as obtained from the solution of the rate equations and the exact solution of the master equation for the density operator using the Monte Carlo wavefunction approach in the truncated single excitation Hilbert space. We observe that the rate equations provide accurate description of the dynamics of the system when $\gamma\gtrsim\sqrt{N}\Omega$, while they cannot account for the (damped) Rabi oscillations between the ground and the excited states when $\sqrt{N}\Omega\gg\gamma$. Nevertheless, the rate equations model approximates well the relaxation timescale and the steady state population of the superatom.

In Fig. 7(b) we show the spectrum of eigenvalues $\lambda_{k}$ of the matrix $\Lambda$ in Eq. (14). The eigenvalue $\lambda_{0}=0$ corresponds to the steady state of the system, while the (negative) real parts of the other eigenvalues, $\mathrm{Re}[\lambda_{k}]$ for $k=1,2,3,4$, characterize the relaxation rates of the superatom towards the steady state. In that figure, we also show the total relaxation rate of a superatom towards the steady state as given by the rate equations model, $\Gamma_{\mathrm{tot}}=\Gamma_{\mathrm{ex}}+\Gamma_{\mathrm{de}}$, which compares favorably with the exact relaxation rate for a wide range of parameters.

Appendix B Excitation facilitation conditions for a chain of superatoms

With the interatomic potential $V(r)=C_{p}/r^{p}$ (assuming repulsive interaction $C_{p}>0$) and the single atom excitation linewidth $w$, the Rydberg blockade distance is defined as $a_{\mathrm{B}}\equiv(C_{p}/w)^{1/p}$. Our model assumes that each superatom can accommodate at most one Rydberg excitation. We therefore require the spatial extent of a superatom to be small compared to blockade distance, $\Delta r\ll a_{\mathrm{B}}$.

Clearly, the superatom excitation probability is maximal at resonant driving and is suppressed when the laser is detuned by more than its excitation linewidth $w_{\mathrm{SA}}$. We consider a one-dimensional chain of superatoms driven by spatially uniform laser with detuning $\Delta$. We set the laser detuning to be equal to the interaction strength between neighboring superatoms separated by the lattice constant $a$, $\Delta=V(a)$. Then an already excited superatom will facilitate the excitation of the nearest neighbors by shifting the Rydberg energy level into resonance with the laser field. On the other hand, we require that $\Delta>w_{\mathrm{SA}}$ in order to suppress excitation of superatoms that have either two non-excited neighbors or two excited neighbors. Since $w_{\mathrm{SA}}>w$, the facilitation distance $r_{\mathrm{fac}}=(C_{p}/\Delta)^{1/p}$, and thereby the lattice constant $a=r_{\mathrm{fac}}$, are smaller than the blockade distance $a_{\mathrm{B}}$.

The interaction potential $V(r)$ is a convex monotonic function. Our approach is valid when the interaction induced Rydberg level shifts of all the atoms of the facilitated superatom are within the superatom linewidth $w_{\mathrm{SA}}$. Linearizing the interaction potential around the facilitation distance $r_{\mathrm{fac}}$, we obtain the following condition on the spatial extent $\Delta r$ of the superatom

[TABLE]

Finally, in our treatment of the 1D lattice of superatoms, we neglect the interaction shifts of the next neighbors, $V(2a)\lesssim w_{\mathrm{SA}}$. This can be rewritten as a condition on the power-law scaling $p$ of the interaction potential,

[TABLE]

Appendix C Derivation of the effective dynamical model for holes

Consider a small chain of superatoms shown in Fig. 8. We focus on a subset of the system containing three superatoms that can be in the ground state or excited to the Rydberg state. We denote the probability of each configuration $\{a,b,c\}$ by $p_{abc}$, where $a,b,c=0,1$ for non-excited (0) and excited (1) superatom on the corresponding site. Starting with all the superatoms excited, $\{1,1,1\}$, we can create a hole on any one site, e.g., $\{0,1,1\}$, with the spontaneous decay rate $\Gamma_{s}$. The reverse process of exciting a hole can be neglected due to smallness of the transition rate $\Gamma_{\mathrm{seed}}\ll\Gamma_{\uparrow,\downarrow},\Gamma_{s}$. Next, an excited superatom neighboring a non-excited one can be de-excited, e.g., $\{0,1,1\}\to\{0,0,1\}$, with rate $\Gamma_{\downarrow}$. An excited superatom between two non-excited ones can be de-excited, e.g., $\{0,1,0\}\to\{0,0,0\}$, with the spontaneous decay rate $\Gamma_{s}$. Finally, the probability of configuration with two or more neighboring non-excited superatoms will quickly decay, e.g., $\{0,0,1\}\to\{1,0,1\}$ or $\{0,1,1\}$, with the the facilitated excitation rate $\Gamma_{\uparrow}$. We thus obtain the following set of equations for the probabilities of various configurations (see Fig. 8 upper panel),

[TABLE]

The probability of full excitation $p_{111}$ decays to zero on the timescale $t>\Gamma_{s}^{-1}$. The probabilities in the last three equations relax with the very fast rates $\sim 2\Gamma_{\uparrow}$ and can be adiabatically eliminated. Setting there $\partial_{t}p_{abc}=0$, and substituting the resulting solutions in the other equations, we obtain

[TABLE]

where we assumed $\frac{\Gamma_{\downarrow}}{\Gamma_{s}}p_{001},\frac{\Gamma_{\downarrow}}{\Gamma_{s}}p_{100}\ll p_{010}$ and $\Gamma_{\downarrow}\ll 2\Gamma_{\uparrow}$.

By adiabatic elimination of the states with neighboring ground state superatoms, we projected the system onto a reduced configuration space. From the above rate equations, we can deduce three fundamental processes affecting the holes (see Fig. 8 lower panel):

• (i)

creation of holes with the rate $\Gamma_{c}=\Gamma_{s}$,

• (ii)

annihilation of one of the two holes separated by one excitation with the rate $\Gamma_{a}=\Gamma_{s}$,

• (iii)

transport of holes between neighboring sites with rate \Gamma_{t}=\mbox{\frac{1}{2}}\Gamma_{\downarrow}.

The corresponding jump operators, involving the constraints that no two neighboring holes are allowed, are given in Sec. III.1.