Atomic-Waveguide Quantum Electrodynamics

Atom arrays are an exciting new type of quantum light-matter interface. Here, we propose to employ one-dimensional ordered arrays as atomic waveguides. These arrays support guided modes that do not decay into free space. We show that these modes can be harnessed to mediate tunable, long-range interactions between additional"qubit"atoms coupled to the chain, without need for photonic structures. The qubits are strongly coupled to the waveguide. Moreover, owing to the two-level nature of atoms, these waveguides are intrinsically quantum. In contrast to classical waveguides, where photons do not interact with each other, atomic waveguides display strong non-linearities. We demonstrate that counter-propagating photons collide, leading to dissipation and scattering out of the chain. This tunable non-linearity opens the door to the exploration of many-body physics between guided photons and to the realization of switches and transistors at the single-photon level. This physics can also be observed in arrays of molecules or solid-state emitters.

The realization of efficient interactions between photons and atoms is a central challenge in quantum optics. Besides enabling the exploration of exotic many-body physics [1], they are also a critical resource to develop practical implementations of quantum information protocols, such as quantum memories and repeaters [2,3]. Deterministic light-matter interactions also form the underpinnings of quantum non-linear optics at the single photon level [4,5], as well as of metrology and sensing applications [6,7].
To control and enhance the interactions between light and atoms, it is generally believed they must be interfaced with nanophotonic structures. This has propelled the development of the field of cavity quantum electrodynamics (QED) and, more recently, of waveguide QED, where atoms are coupled to one-dimensional (1D) photonic reservoirs, such as fibers [8][9][10][11] and photonic crystal waveguides [12][13][14][15]. These platforms offer efficient lightmatter coupling. Moreover, by engineering the dispersion relation of the optical modes, they allow for the realization of almost arbitrary interactions between atoms [16]. However, interfacing quantum emitters with nanophotonics in a scalable manner has proven to be technically difficult.
Here, we suggest an alternative approach: to employ atomic arrays as one-dimensional waveguides that mediate interactions between distant "qubits". In ordered arrays, interference in photon emission leads to the emergence of subradiant states, which do not decay into free space. These states can be understood as guided modes of the atomic chain [17], and can be used to mediate both coherent and dissipative interactions between qubits that are coupled to the atomic waveguide. The coupling between these qubits and the waveguide can be remarkably efficient.
Atomic waveguides also provide a platform to observe many body physics between guided photons. We demonstrate that excitations display strong non-linear behavior, as these waveguides are intrinsically quantum (a) Schematic of the setup under consideration: a "qubit" atom (red) of resonace frequency ωq is located in the vicinity of a 1D chain of atoms of resonance frequency ω0 and spontaneous emission rate Γ0 (blue), at a distance ρq from the axis of the array. If the distance d between the array atoms is smaller than λ0/2, the chain behaves as a waveguide and supports guided modes that do not decay into free space. The qubit emission rates into the atomic-waveguide mode and into free space are Γ q 1D and Γ q , respectively. (b) Dispersion relation for the single-excitation mode of an infinite, 1D chain of atoms polarized parallel to the chain axis with spacing d = 0.1λ0. In the region enclosed within the light line (shaded), the chain does not guide light and the mode decays into free space. Beyond the light line (kz > k0 = ω0/c), the mode is guided.
(a single atom cannot be excited twice). In particular, two photons collide with each other whenever they spatially overlap, due to dissipative non-linearities. We note that this nonlinearity is different from that realized with Rydberg states [4,5]. Together, our results arXiv:1912.06234v1 [quant-ph] 12 Dec 2019 provide a comprehensive picture of linear and non-linear quantum optical phenomena that can be observed in 1D atomic waveguides.

Atom array as a waveguide
Ordered arrays of atoms support guided photons (in the form of spin waves) that do not scatter light into free space. These guided modes emerge due to destructive interference in photon emission. To understand this phenomenon, we write down a "spin model" that describes dipolar interactions between atoms. We consider an array of N two-level atoms of resonance frequency ω 0 separated by a distance d, as shown in Fig. 1(a). We describe the atoms' dynamics employing a stochastic wavefunction approach [18,19], where the atomic state evolves under the non-Hermitian Hamiltonian interrupted by quantum jumps that lower the number of excitations and occur at random times. The coherent and dissipative interaction rates between atoms i and j read [17] where ℘ is the dipole matrix element associated with the atomic transition. The Green's tensor G 0 (r i , r j , ω 0 ) is the propagator of the electromagnetic field between atoms i and j in vacuum. It admits the closed expression where k 0 = ω 0 /c. For a single atom, the spontaneous emission rate is simply , and the local frequency shift simply renormalizes the resonance frequency ω 0 . The Hamiltonian of Eq. (1) only contains spin degrees of freedom (i.e., the atomic coherence operatorŝ σ i ge = |g i e i | between the ground and excited states, and the population operatorσ i ee = |e i e i |). In the presence of a driving field of frequency ω, the equations are identical, but with the prescription ω 0 → ω.
The non-Hermitian Hamiltonian of Eq. (1) is that of an open, long-range XY model, and is derived within the Born-Markov approximation. This approximation allows for integrating out the electromagnetic degrees of freedom and requires two conditions. First, that the spectral response of the reservoir is flat compared to that of the atoms (such that the Green's function is evaluated at the atomic resonance frequency). Second, that retardation can be ignored (such that the Hamiltonian is local in time). These conditions are fulfilled in vacuum unless the separation between atoms is extremely large (of the order of a meter for typical optical transitions) [20,21].
In the single-excitation regime, guided modes emerge for inter-atomic distances d < λ 0 /2, with λ 0 = 2πc/ω 0 [17]. To demonstrate their waveguiding behavior, we analyze the physics of an infinite chain that extends along the z direction. The eigenstates of the Hamiltonian are spin waves with well defined wave-vector k z ∈ [−π/d, π/d], and are generated by the action of the collective spin operatorŜ on the ground state |g ⊗N . These Bloch modes satisfy is a complex frequency whose imaginary part describes the decay rate of the spin wave, and its real part accounts for a frequency shift with respect to the bare atomic resonance. In the above equation,G 0 (k z ) = j e −ikzzj G 0 (r j ) is the discrete Fourier transform of the free-space Green's tensor [see Supplementary Information (SI) [22] for details]. Figure 1(b) shows the dispersion relation (i.e., the real part of ω kz vs k z ) of a chain with lattice constant d = 0.1λ 0 , for atoms polarized along z, the direction of the array [23]. For |k z | > k 0 , the spin waves have a zero decay rate and are guided modes of the array [17,24]. For |k z | ≤ k 0 , the spin waves have a finite lifetime and decay due to photon emission. Guided modes also exist in a finite chain. Their decay rate is non-zero and scales as ∼ 1/N 3 [17].
The quantum character of the waveguide is revealed only in the presence of multiple excitations. The emergence of guided modes in the single excitation manifold is not a uniquely quantum property. Waveguiding behavior is also displayed by arrays of classical dipoles, such as subwavelength grating waveguides [25] and chains of dielectric [26] and metallic nanoparticles [27].

Coupled qubit decay rates
The decay rate of a "qubit" atom in the vicinity of the array is altered by the presence of the waveguide [see Fig. 1(a)]. The qubit can decay into free space (whose modes are modified by the presence of the waveguide) and into the guided mode of the array, exciting spin waves that propagate away from the qubit without scattering. We calculate the decay rates into free space and the guided mode by computing the Green's tensor of the surrounding environment, i.e., the vacuum and the waveguide. Exploiting the cylindrical symmetry of the infinite chain, we find an expression of the Green's tensor in terms of an integral over reciprocal space (see [22] for details). For a drive frequency ω, the Green's tensor For all plots, the detuning between qubit and waveguide atoms is chosen such that the guided wave-vector is k 1D = 0.7π/d (as shown by the purple line in Fig. 1). reads G(r, r , ω) = G 0 (r, r , ω) where k = ω/c and the field eigenmodes are given by (6b) In the above expressions,℘ = ℘/|℘|, H (1) 0 is a Hankel function of the first kind, ρ is the radial distance to the chain axis, and k ⊥ = k 2 − (k z + g) 2 is the transversal wavevector. The sum is performed over reciprocal-lattice vectors g = 2πn/d, with n ∈ Z, and accounts for Umklapp processes.
The efficiency of the coupling to the waveguide is given by the ratio between the guided (Γ q 1D ) and the freespace (Γ q ) decay rates. The analytical expression for the Green's function provides an elegant way to compute these rates separately. For atoms in free space, the decay rate is given by Eq. (2b). Similarly, we postulate that the decay rate of the qubits is proportional to the imaginary part of the generalized Green's tensor of Eq. (5). We thus trace out the waveguide atoms, and treat the chain as a bath for the qubits [as the photons were integrated-out to derive Eq. (1)]. This procedure is only exact within the single excitation subspace (the atoms are spins, not bosons) and under the Born-Markov approximation. This implies that retardation can be ignored (i.e. that the group velocity of the guided mode is not too small) and that the decay rate of the qubit is much smaller than the bandwidth of the waveguide. We discuss how to achieve this in the section "physical implementations".
The qubit decay into free space is given by radiative wave-vectors (i.e. |k z | < k). It reads where Γ q 0 , ℘ q , and r q are the qubit's vacuum spontaneous emission rate, dipole matrix element, and position, respectively.
The decay into the guided mode arises from the pole of the Green's function, and is readily found to be where |k 1D | > k is the guided mode wave-vector (i.e., the wave-vector k z for which ω kz ≡ ω), and v g = ∂ω kz /∂k z | kz=k 1D is the group velocity. The decay rate into the guided mode increases for low group velocities, as the mode has more time to interact with the qubit. Close to the band-edge, the group velocity is low, and the decay into the waveguide becomes large (the waveguide behaves more like a photonic crystal than a fiber in this region). Note that k k 0 as ω ω 0 , ω q except for deviations of the order of Γ 0 ω 0 , ω q . The qubit interacts efficiently with the atomic waveguide mode, as shown in Fig. 2. The ratio between guided and free-space decay rates -so-called optical depth D = Γ q 1D /Γ q -can be larger than 1. The optical depth displays a modulation along z related to the Bloch periodicity, with a contrast that decreases with the distance to the array. We find simple scaling laws for the optical depth when the qubit is placed exactly on top of a waveguide atom [i.e., at z q = 0, see Fig 2( for constant ρ q /d, and Γ q 1D /Γ q ∼ 1/ρ 6 q for fixed d and ρ q 0.4d, below which the coupling rates plateau. Remarkably, a waveguide that is one-atom thick provides an optical depth ∼ 30 times larger than that of a fiber (see SI for details [22]). We corroborate the analytic calculations with numerical simulations, by evolving a finite chain and coupled qubit under the non-Hermitian Hamiltonian of Eq. (1). The numerical and analytical results fully agree with each other for qubits in the central part of the chain, where finite size effects are negligible.
Surprisingly, we find a magic point where the emission into free space is strongly suppressed. At ρ q 0.4d and z q = d/2, the optical depth is extremely large while the total linewidth of the qubit remains small, as shown in Fig. 2(c). The existence of the magic point is solely due to interference. We study change in the location of the magic point as a function of lattice constant in the SI [22].

Tunable-range interactions between qubits
The atomic chain mediates long-range interactions between distant qubits, without the need of photonic structures. The character of these interactions (coherent or dissipative) is modified by tuning the qubit resonance frequency ω q . If ω q lies within the band, the interactions are dissipative. If, instead, ω q lies outside the band, the interactions are coherent.
Dissipative qubit-qubit interactions lead to superradiant decay. We analyze photon transport through an atomic waveguide with either one or five periodically spaced qubits coupled to it, as shown in Fig. 3(a). To calculate transmittance and reflectance spectra, we launch a spin wave of the form where ζ = 300d is the spatial spread,z i are the atomic positions relative to the center of the spin wave, and k 1D is chosen to determine the relative detuning between qubit and spin wave, i.e., such that ∆ ≡ (ω k 1D −ω 0 )−ω q . We discuss how to prepare such a state in the section "physical implementations". The evolution is performed under the Hamiltonian in Eq. (1), conditioning the results on no jumps. Interaction with the qubit -positioned 1000 sites from the initial spin wave center-leads to strong reflection of the spin wave. We calculate the transmission (reflection) from the population of the array atoms located past (before) the qubit, while the lost norm provides the scattering into free space. The spectra display the traditional Lorentzian lineshape with a width that scales as Γ q 1D + Γ q [20,28]. We perform a similar calculation for five qubits separated by a distance such that k 1D d q = 14π. In conventional waveguide QED, this corresponds to the mirror configuration, where the qubits behave collectively as a single qubit with a larger dipole moment and superradiantly decay (a) Transmission (blue) and reflection spectra (orange) of a guided mode that interacts with either a single qubit or with five qubits evenly spaced 20 lattice sites apart. In both cases, the qubits' vacuum spontaneous emission rate is Γ q 0 = 0.02Γ0, sit at ρq = d and are detuned from the waveguide atoms such that light with wavevector k 1D = 0.7π/d is near-resonant with the qubits. The chain consists of N = 4000 atoms, and the lattice constant is d = 0.1λ0. (b) Evolution of the excitedstate population of qubit 1 (red) and 2 (green) after fully inverting qubit 1 at the initial time. The Rabi oscillations reveal strong coupling mediated by photonic bound states. The resonance frequency of the qubits lies in the band gap of the atomic waveguide (at a frequency 4.5Γ0 from the band edge), with the qubits (of spontaneous emission rate Γ q 0 = 0.001Γ0) placed in between two array atoms at ρq = 0.4d and separated by a distance of 8d, for a chain of N = 99 atoms of lattice constant d = 0.05λ0. at a rate 5Γ q 1D + Γ q [20,28]. Our transmission spectra, shown in Fig. 3(a), deviates slightly from a Lorentzian profile due to non-Markovian effects associated with retardation [22]. The group velocity of the spin chain is remarkably slow compared to the speed of light in free space, scaling as v  Coherent qubit-qubit interactions give rise to spin exchange. They occur if the qubits' frequency sits beyond the band edge [29][30][31]. In this regime, spin waves cannot propagate and form bound states that are spatially localized around the qubit position. Mathematically, it is easy to see that a resonant excitation cannot propagate through the array as there is no pole in the integral of Eq. (5). Bound states mediate coherent, finite-range interactions between qubits, as shown in Fig. 3(b). Since the spin-exchange rate is small, we place the qubits at magic points, where the free space decay rate is strongly suppressed. This allows us to observe several oscillations before the dynamics is damped. Without the waveguide, given the qubits' separation, they would simply decay.

Quantum nonlinearity and photon collisions
The physics of two-photon transport is qualitatively different from that observed in a classical waveguide. For one excitation, the system behaves as a chain of classical dipoles, and photon transport obeys the laws of linear optics. However, a single atom cannot support two excitations, which leads to striking consequences for the statistics of many-body subradiant states. In particular, (Ŝ † kz ) 2 |g ⊗N is not an eigenstate of the Hamiltonian. Instead, the true eigenstates behave as fermions (or hardcore bosons), obeying a Pauli exclusion principle in space [17,32]. As we demonstrate below, this leads to strong dissipative interactions between photons.
The spatial overlap of two photons (spin waves) produces scattering into free space. If the two spin waves propagate in opposite directions, this can be interpreted as a photon-photon collision, as shown in Fig. 4(a). To observe such interaction, we initialize a two-excitation state where ζ = 15d andz i,j are the atom positions relative to centers 60 sites apart.
We evolve the wavefunction under the Hamiltonian in Eq. (1). In Fig. 4(b), we show the decrease of population ( σ ) as a function of time. The population loss ranges from less than 1% (for d = 0.1λ 0 ) to almost 60% (for d = 0.3λ 0 ) for |k 1D | = 0.7π/d; the inset in Fig. 4(b) shows that the loss grows with the distance as a power law. For large lattice constants, jumps lead to the emission of one and (most probably) two photons. For small distances, the probability of photon emission is negligible and the effect of the collision reduces to the acquisition of a global phase.
Dissipation during the photon collision is due to the spatial distortion of the guided modes. The collision generates population in radiative modes (those with |k z | < k 0 ), as shown in Fig. 4(c). To corroborate that the loss is due to photon emission, we calculate the radiated field at different times using an input output equation [17]. Before and after the collision there is no light emission as the spin excitations are guided. At the collision time the field leaks out of the waveguide, as shown in Fig. 4(d).
The degree of interaction between two excitations is controlled by tuning a few experimentally accesible parameters (frequency and temporal duration of the photons, and lattice constant). The decay probability is a function of the overlap between the distorted state and radiative modes, the rate of scattering of those modes, and the interaction time. Smaller lattice constants lead to larger group velocities, shorter interaction times, and smaller light cones.

Physical implementations
The implementation of an atomic waveguide involves two main challenges: to trap atoms at short distances and to efficiently excite guided modes. Coupling qubits to the waveguide comes with an extra set of difficulties. We discuss strategies for overcoming these challenges below. For most of the discussion, the experimental setup we have in mind consists in neutral atoms trapped in optical lattices [33][34][35] or optical tweezers [36][37][38][39][40], which have recently been suggested as quantum metasurfaces [41].
First, we require small inter-atomic separations (d < λ 0 /2). The diffraction limit can be overcome using two different atomic transitions: one to trap and one to drive the optical excitation. As an example, Strontium can be trapped at a magic wavelength with d = λ 0 /16.3 [42,43], driven on the λ 0 = 2.6 µm 3 P 0 → 3 D 1 transition. The bosonic species lacks hyperfine structure, which prevents additional difficulties [44,45]. Another possibility would be to use Ytterbium's telecom transition [46]. Quantum and classical disorder may affect the guiding properties of the waveguide. Nevertheless, guided modes have been shown to be robust against spatial disorder [24,45]. The finite spread of the atomic wavefunction adds an independent decay channel for each atom (Γ trap ∼ Γ 0 η 2 , where η is the Lamb-Dicke parameter [47]), but can be reduced using tight traps.
Second, we need to excite guided modes efficiently. The frequency of the external field selects the wave-vector of the spin wave that propagates in the array, and the temporal duration of the laser pulse sets its spatial width. Coupling is possible by focusing the external light into the array edge, either with a lens with high numerical aperture or with a spatial light modulator. One can also employ coupled qubits to inject spin waves into the waveguide. The coupling loss can be alleviated by using a near-field probe, such as a fiber tip close to the ar-ray. Finally, a phase could be imprinted, via magnetic or optical fields [48,49], into easily-accessible superradiant states.
We require the frequency of the qubits to be distinct from that of the waveguide atoms. The qubit frequency can be tuned with Stark shifts through optical tweezer beams. To realize Markovian interactions, the waveguide bandwidth has to be broad compared to the qubit linewidth. One option is to rely on compact chains, as the bandwidth decreases with d. Another is to use different atomic isotopes (e.g., 87 Sr and 88 Sr), as they have similar transition frequencies [50] but different linewidths due to hyperfine structure. For distances d/λ 0 = 16.3, the ratio between waveguide bandwidth and qubit linewidth is ∼ 400. Cold molecules are also interesting candidates, as they have dense frequency spectra and have been recently trapped in tweezer arrays [51].
This physics can also be observed in arrays of solid state qubits, such as color centers [52], rare-earth ions [53,54], and localized excitonic quantum dots or strain-generated defects in 2D materials [55,56]. While deteministic placement of solid state emitters is becoming a reality, these emitters have their own set of issues, mostly related to inhomogeneous broadening and non-radiative decay.

Outlook
We have demonstrated that atomic waveguides are versatile quantum light-matter interfaces. They support strong non-linearities that may enable a number of quantum applications, including single-photon transistors and switches. They also mediate waveguide-less long-range interactions between qubits, without the need for interfacing atoms with traditional nanophotonic structures. They are uniquely programmable and can be created "on demand". For instance, the propagation of excitations in the array can be dynamically altered by using electromagnetic fields to tune the resonance frequency of the atoms. They also provide opportunities to explore the rich physics of self-organization [57,58], as optomechanical degrees of freedom might play an important role in determining the geometry of the array.
Atomic waveguides can be harnessed to explore less traditional QED paradigms, such as chiral quantum optics and time-delayed interactions. The near field of the waveguide has longitudinal and radial components. Employing a set of qubits with the right hyperfine structure would result in cascaded propagation of excitations, just as if the qubits were coupled to a fiber [59]. Time-delayed interactions allow for the study of retardation effects in interacting quantum systems. As we discuss in the SI, atomic waveguides present multiple opportunities to explore non-Markovian physics, allowing for engineering coherent feedback loops for qubits [60]. All these examples serve to illustrate the wealth of quantum phenomena that can be explored in atomic waveguides.

Derivation of the decay rates
Here we derive the expressions for the qubit decay rates into free space and into the guided mode of the atomic array. The decay rate of an emitter is related to the imaginary part of the Green's tensor evaluated at the emitter's position. We thus begin by finding an expression for the propagator of the electromagnetic field in the presence of the chain. Following the main text, we consider that the atomic array is pumped by a weak coherent field (such that saturation is negligible and the dynamics is confined to the single excitation manifold). The equation of motion for the expectation value of the atomic coherence operator of atom j is readily found to bė where G 0 (r i , r j , ω) is the vaccuum Green's tensor, σ j eg ≡ σ j ge , E + p = Ê + p is the expectation value of the positivefrequency component of the driving field, ∆ = ω − ω 0 is the detuning between the driving and the atomic resonance frequencies, and r j = (R j , z j ) is the position of atom j. For an infinite chain, we can define a spin-wave operator, eg that creates an excitation of well-defined longitudinal momentum k z . In the steady state (i.e.σ j eg = 0), the expectation value of the spin-wave annihilation operator reads where is the spatial Fourier transform of the field. Employing the inputoutput equation [17]Ê we can write the expectation value of the field in any point in space as From this equation, we can obtain an expression for the Green's tensor of the medium consisting of vacuum modified by the presence of the atomic chain. To do so, we assume that the pump field is generated by a dipole-like source p at r p , which generates the current j(r, ω) = −iωpδ(r − r p ), such that Substituting the above expression into Eq. (14), and transforming the sum over k z into an integral over the Brillouin zone, i.e., we find that E + p (r) = µ 0 ω 2 G (r, r p , ω) · p, with a Green's tensor that now accounts for the presence of the chain and can be written as G(r, r , ω) = G 0 (r, r , ω) − µ 0 ω 2 d 2π where we have defined We now express the vacuum Green's tensor in cylindrical coordinates, by making use of the integral representation of spherical waves [61] e ik|r−rj | |r − r j | = i 2 where Li s (z) = ∞ =1 z −s is a polylogarithm function of order n. The decay rate of a "qubit" atom placed in the vicinity of the chain is directly related to the imaginary part of the Green's tensor through where r q and ℘ q are the qubit position and dipole matrix element, respectively. The integration path is shown in Fig. 5. The integrand displays several branch cuts (associated with the square root and polylogarithm functions), as well as simple poles for k z such that ω kz = ω. We can clearly separate two different contributions to the decay: emission into free space (in the region such that k z ∈ [−k, k]), and emission into the atomic waveguide mode (due to the pole at k z = k 1D ).

Free-space decay rate
The presence of the chain alters the vacuum modes and thus leads to a modified decay rate of the qubit, which is now calculated not only from the vaccuum's Green's tensor, G 0 , but also taking into account a contribution to the integral arising from wavevectors within the light cone, i.e., k z ∈ [−k, k]. The free-space decay rate is thus readily found to be There is also a frequency shift that arises from the real part of the Green's function, which can be calculated numerically by taking the real part of the integrals along I r , I c , and I cc , as shown in Fig. 5.

Guided-mode decay rate
For an infinite chain, we can perfectly isolate the decay into the guided mode of the atomic-waveguide as it appears as a pole in the integral. Beyond the light line, v kz = −u † kz and ω kz is real (as the guided mode has infinite lifetime, i.e., it is not "leaky"). This means that the imaginary part of the integral is zero everywhere outside the light cone, except for the poles where ω kz = ω. Note that if the driving frequency is detuned from the guided mode band, there is no pole contribution and thus no decay. For longitudinal polarization there are two poles (corresponding to forward and backward propagating guided waves at ±k 1D ). We close the integral around one of these poles, as shown in Fig. 5, and find Approximating ω kz ω 0 +v g k z , where v g is the group velocity at frequency ω, and applying Cauchy's residue theorem, we find a simple expression for the guided decay rate:

Functional form of the field modes
For waveguide atoms polarized along the direction of the chain, the different polarization components of u kz (r) readρ φ · u kz (r) = 0, The components of v kz (r) admit similar expressions.

Group velocity of the guided mode
For an infinite array, we can calculate the group velocity of the guided modes as the derivative of the dispersion relation of Eq. (24), i.e., v g = ∂ω kz /∂k z | kz=k 1D . We show the group velocity scalings with d and k z in Fig. 6. Note that k k 0 as ω ω 0 , ω q except for deviations of the order of Γ 0 ω 0 , ω q . Beyond the light line [corresponding to the peaks in Fig 6(b)] the mode is guided, and the group velocity tends to zero as k z approaches the edge of the Brillouin zone, though the dependence on k z is not trivial for small distances. For fixed k z d, v g ∼ d −1.7 where the exponent is approximate and varies slightly for different k z .

Spatial dependence of the decay rates
The decay rate into free space is not simply that of a qubit in vacuum (Γ q 0 ), but is modified by the presence of the atomic waveguide, which alters the vacuum modes. This decay rate displays a non-trivial dependence with the position of the qubit, as shown in Fig. 7. A similar scaling is followed by the decay rate into the waveguide mode. Generically, both decay rates are enhanced for short radial and longitudinal distances to the waveguide atoms. However, there are magic points -manifested as dark lines in the figures -where decay rates are strongly suppressed due to interference effects. For the free space scattering rate, these lines appear as a narrow band at ρ q ≈ 0.4d in between two array atoms (at z q = ±0.5d) and then move towards the central atom as ρ q decreases. For the guided-mode scattering rate, these positions draw virtually straight lines that appear at ρ q ≈ 0.6d in the middle of two array atoms. This translates into a ratio between guided-mode and free-space scattering that is strongly enhanced at z q = ±0.5d, as shown in Fig. 2 in the main text. This pattern is robust against changes in inter-atomic distance d and displays only minor changes when altering k 1D .
As we discuss in the main text, a waveguide that is one-atom thick provides an optical depth ∼ 30 times larger than that of a fiber. To estimate these numbers, we have considered a waveguide with lattice constant d = 0.2λ 0 and a qubit frequency such that k 1D = 0.7π, placed at z q = 0, ρ q = 0.1λ 0 [red line in Fig. 2(b)]. The fiber has radius k 0 r = 1.2 and is made of silicon nitride, with dielectric constant = 4. The qubit is located at ρ q from the surface of the fiber (this leads to Γ q 1D /Γ q 0.3 [17]). Note that the coupling to the atomic waveguide can be increased by placing the qubit frequency closer to the band-edge.

Non-Markovian behavior of coupled qubits
In the main text, we propose the use of an atomic array as a waveguide. By treating the array as a bath for a coupled qubit, we obtained the qubit decay rates into free space and guided mode of the array. For this to be an accurate description, we require a Markovian interaction between the qubit and the bath, i.e., the linewidth of the qubit, Γ q 1D + Γ q , must be small with respect to the bandwidth of the waveguide. To achieve this, we can either reduce the linewidth of the qubit with respect to the array atoms, or use small inter-atomic spacing to increase the bandwidth of the bath. For multiple qubits, we also require negligible retardation in the propagation of the atomic spin wave. Here, we discuss what happens outside of these assumptions.
• Non-Markovianity due to strong coupling. This occurs for qubits with large linewidth, so that they spectrally sample different regions of the waveguide dispersion relation. In Fig. 8(a), we plot the decay rate of an excited qubit into the array by evolving an initially excited qubit with the non-Hermitian Hamiltonian of Eq. (1). In the Markovian case, there is an exponential decay of the qubit into guided modes of the array and free space.
Here, the population scattered into the array is reabsorbed by the qubit and then reemitted. This results in oscillations between the qubit and a bound state of the array. This bound state is different to that appearing when the qubit frequency is in the band gap, where the spin excitation is bound because it cannot be guided.
Here, the bound state appears because the interaction with the qubit happens on a much faster timescale than the transport of the excitation along the array. The bound state is trapped by its strong interactions with the qubit. The qubit scatters into free space cyclically, while the array does not.
• Non-Markovianity due to retardation. The slow group velocity produces retardation in multi-qubit interactions. In Fig. 8(b), we show the evolution of the populations of five evenly spaced qubits coupled to the array. As described in the main text, we launch a spin wave and evolve the system under the non-Hermitian Hamiltonian of Eq. (1). The spin wave is detuned by 14.5Γ q 0 from the resonance frequency of the qubits, corresponding to a point on the shoulder of Fig. 3(a). Each of the five qubits is excited at slightly different times. The delay in the excitation is small, but enough to break the Markovianity of the waveguide [62][63][64][65]. The analysis of the ensemble as a single large qubit in the mirror configuration is not exactly correct, which leads to the minor differences from the expected Lorentzian shape discussed in the manuscript. We can implement tuneable time-delayed interactions by using this phenomenon. The length of the delay can be tuned by altering the distance between qubits or the group velocity, which is slowest for large k 1D and d. In Fig. 8(c), we launch a spin wave with low group velocity, and plot the evolution of the populations of two qubits spaced 800 sites apart coupled to the array. We see a long delay between each qubit being excited, such that the first qubit is almost completely de-excited before the second atom interacts with the pulse. There are small oscillations in the excitation, as the slow group velocity and high coupling produces some of the trapping effects described above.  . (b,c) Evolution of the populations of (b) five qubits evenly spaced 20 lattice sites apart and (c) two qubits spaced 800 sites apart, coupled to an array of 4000 atoms. In both plots, the qubits' vacuum spontaneous emission rate is Γ q 0 = 0.02Γ0, they sit at ρq = d and are detuned from the waveguide atoms such that light with wavevector (b) k 1D = 0.7π/d, or (c) k 1D = 0.95π/d is near-resonant with the qubits, and the initial state of the array is a spin wave with frequency shift detuned 14.5Γ q 0 below the coupled qubits. In all cases, d = 0.1λ0. Figure 9 replicates Fig. 4 in the main text for a chain with larger lattice constant (d = 0.3λ 0 ), showing significant population loss (around 60%). In this case, a large region of the Brillouin zone falls inside the light cone (as depicted by the dashed lines), and thus a larger proportion of the wave-vector components generated in the collision are radiative. Moreover, the slower group velocity increases the interaction time, leading to the population remaining inside the light cone for much longer time, with t f being more than ten times larger for d = 0.3λ 0 than for d = 0.1λ 0 . At t = ti, two counter-propagating spin waves are initialized with wave-vector k 1D = ±0.7π/d. They collide, leading to population loss (the maximum overlap occurs at t = tc). After the collision has occurred they propagate without loss until they hit the ends of the chain. Before this occurs, they first return to their original positions, at t = t f ≡ 2tc. (b) Probability amplitude in real space (top) and in reciprocal space (bottom) of the two-excitation wavefunction for the three different times shown in (a). The light lines are shown by dotted lines (the area delimited by the lines is the light cone and corresponds to radiative modes). At the collision time t = tc there is significant population inside the light cone, which leads to loss. (c) Emitted field intensity at t = tc (at the plane y = 5d). The atoms are depicted by white circles.