Nondestructive photon counting in waveguide QED

TL;DR

This paper introduces a method to perform nondestructive photon counting in waveguide QED by entangling photons with atomic arrays, enabling high-fidelity, scalable, and robust photon number measurements.

Contribution

The authors propose a novel scheme that maps photon counting to atom counting in waveguide QED, allowing nondestructive, high-fidelity measurements scalable with atom number.

Findings

Fidelity increases with atom number

Scheme is robust to disorder and finite Purcell factors

Can implement nondestructive photon-number measurement

Abstract

Number-resolving single-photon detectors represent a key technology for a host of quantum optics protocols, but despite significant efforts, state-of-the-art devices are limited to few photons. In contrast, state-dependent atom counting in arrays can be done with extremely high fidelity up to hundreds of atoms. We show that in waveguide QED, the problem of photon counting can be reduced to atom counting, by entangling the photonic state with an atomic array in the collective number basis. This is possible as the incoming photons couple to collective atomic states and can be achieved by engineering a second decay channel of an excited atom to a metastable state. Our scheme is robust to disorder and finite Purcell factors, and its fidelity increases with atom number. Analyzing the state of the re-emitted photons, we further show that if the initial atomic state is a symmetric Dicke state, dissipation engineering can be used to implement a nondestructive photon-number measurement, in which the incident state is scattered into the waveguide unchanged. Our results generalize to related platforms, including superconducting qubits.

Equations

$$H=$$

$$-$$

$${\dot{\sigma }}_{ge}^{(n)}=({\sigma }_{gg}^{(n)}-{\sigma }_{ee}^{(n)})$$

$$+{\sigma }_{gs}^{(n)}$$

$$\dot{\mathbf{b}}=(-i{\mathsf{H}}_{\mathrm{eff},1}-{\Gamma }^{\prime }/2)\mathbf{b}+{\mathsf{L}}_{\nu ,1}{\mathbf{a}}_{\mathrm{in},\nu ,1}.$$

$${\mathbf{a}}_{\mathrm{out}}(\omega )$$

$$\equiv \mathsf{S}(\omega ){\mathbf{a}}_{\mathrm{in}}(\omega ).$$

$${\mathsf{S}}^{-1}\mathsf{S}=1-{\mathsf{L}}^{†}\mathsf{M}(\omega )\left[\mathsf{L}{\mathsf{L}}^{†}-i({\mathsf{H}}_{\mathrm{eff}}-{\mathsf{H}}_{\mathrm{eff}}^{†})\right]{\mathsf{M}}^{†}(-\omega )\mathsf{L}=1,$$

$$\mathsf{M}(\omega )\equiv (-i\omega +i{\mathsf{H}}_{\mathrm{eff}}+{\mathsf{\Gamma }}^{\prime }/2{)}^{-1}.$$

$${\mathsf{H}}_{\mathrm{eff},1,mn}=-i\frac{{\Gamma }_g}{4}\left[e^{i{k}_0|x_m-x_n|}-e^{i{k}_0(x_m+x_n)}\right].$$

$$\dot{B}(t)=-\frac{{\Gamma }_{\mathrm{tot}}}{2}B(t)+\sqrt{N{\Gamma }_g}{a}_{\mathrm{in}}(t),$$

$$a_{\mathrm{out}}(t)=a_{\mathrm{in}}(t)-\sqrt{N{\Gamma }_g}B(t),$$

$$⟨a_{\mathrm{out}}^{†}(\omega )a_{\mathrm{out}}(\omega )⟩={\left|1-\frac{N{\Gamma }_g}{{\Gamma }_{\mathrm{tot}}/2-i\omega }\right|}^2⟨a_{\mathrm{in}}^{†}(\omega )a_{\mathrm{in}}(\omega )⟩.$$

$${\mathsf{H}}_{\mathrm{eff},1,mn}=-i{\Gamma }_g\exp (i{k}_0|x_m-x_n|),$$

$$\mathsf{S}(\omega )=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)-\left(\begin{array}{cc}{\mathbf{L}}_{-}\mathsf{M}(\omega ){\mathbf{L}}_{+}& {\mathbf{L}}_{-}\mathsf{M}(\omega ){\mathbf{L}}_{-}\\ {\mathbf{L}}_{+}\mathsf{M}(\omega ){\mathbf{L}}_{+}& {\mathbf{L}}_{+}\mathsf{M}(\omega ){\mathbf{L}}_{-}\end{array}\right),$$

$$\mathsf{S}(\omega )=\left(\begin{array}{cc}A(\omega )& B(\omega )e^{2i\varphi }\\ B(\omega )& A(\omega )\end{array}\right),$$

$${\mathsf{S}}_{\mathrm{AMC}}(\omega =0)=\frac{1}{2}\left(\begin{array}{cc}1& -1\\ -1& 1\end{array}\right),$$

$${\beta }_{\pm }=\frac{{\gamma }_{\pm }}{{\gamma }_{+}+{\gamma }_{-}+{\Gamma }^{\prime }},t_{\pm }=1-2{\beta }_{\pm },r_{\pm }=-2\sqrt{{\beta }_{+}{\beta }_{-}}.$$

$$A_{+}\equiv 1-|t_{+}{|}^2-|r_{+}{|}^2=4{\beta }_{+}\left(1-{\beta }_{+}\frac{{\gamma }_{-}}{{\gamma }_{+}}\right).$$

$${\eta }_{\mathrm{chiral}}=\left(1-\frac{{\gamma }_{-}}{{\gamma }_{+}+{\Gamma }^{\prime }}\right)\frac{{\Gamma }_s}{{\Gamma }^{\prime }}.$$

$$|N-m,m,0⟩\equiv \sqrt{\frac{(N-m)!}{N!m!}}[S_{sg}{]}^m|G⟩,$$

$${\dot{S}}_{ge}=(S_{gg}-S_{ee})\left(\sqrt{{\Gamma }_g}{a}_{\mathrm{in},1}-\frac{{\Gamma }_g}{2}{S}_{ge}\right)$$

$$+S_{gs}\left(\sqrt{{\Gamma }_s}{a}_{\mathrm{in},2}-\frac{{\Gamma }_s}{2}{S}_{se}\right)-\frac{{\Gamma }_{\mathrm{free}}}{2}{S}_{ge},$$

$$a_{\mathrm{out},1}=a_{\mathrm{in},1}-\sqrt{{\Gamma }_g}{S}_{ge},a_{\mathrm{out},2}=a_{\mathrm{in},2}-\sqrt{{\Gamma }_s}{S}_{se}.$$

$$S_{gs}{S}_{se}|\psi ⟩=(S_{ss}+1)S_{ge}|\psi ⟩.$$

$${\dot{S}}_{ge}=-$$

$${\Gamma }_g(N/2-m_0)={\Gamma }_s(N/2+m_0+1)$$

$${\mathcal{F}}_1=\frac{\kappa }{\kappa +{\kappa }_{\mathrm{int}}}\left\{1-\sqrt{2\pi }{e}^{2{\kappa }^2{\tau }^2}\kappa \tau [1-\mathrm{erf}(\sqrt{2}\kappa \tau )]\right\},$$

$${\mathcal{F}}_{m_{\mathrm{p}}}=1-\left(1-\frac{4m_{\mathrm{p}}^3-6m_{\mathrm{p}}^2+2m_{\mathrm{p}}}{3N^2}\right){\mathcal{F}}_1^{m_{\mathrm{p}}}.$$

$$H=\frac{2\pi g_k^2}{\Delta L}\sum _{ij}{\sigma }_{sg}^{(i)}{\sigma }_{gs}^{(j)}{e}^{-|x_i-x_j|/L}\simeq g_{\mathrm{eff}}{S}_{sg}{S}_{gs},$$

$$H=g_{\mathrm{eff}}[S_{sg}^{(1)}{S}_{gs}^{(2)}+S_{sg}^{(2)}{S}_{gs}^{(1)}].$$

$${\hat{L}}_{g,\mathrm{eff}}$$

$${\hat{L}}_{s_i,\mathrm{eff}}$$

Full text

I Introduction

Single-photon detectors have a long history Morton (1949), with a plethora of technologies available Hadfield (2009). Applications in quantum optics, such as quantum state preparation, quantum metrology Giovannetti et al. (2011), entanglement distribution Gisin et al. (2002), and quantum computing Knill et al. (2001); Kok et al. (2007) have placed a renewed focus on single-photon detectors capable of resolving the number of incoming photons. Perhaps most promising are superconducting transition-edge sensors, which have been demonstrated to achieve a (per-photon) detection efficiency (percentage of detected photons) of $\eta\approx 95\%$ Lita et al. (2008) and to distinguish up to seven photons, with a negligible dark count rate (clicks in the absence of incoming photons). They are based on the principle that near the critical temperature, the resistance of a superconductor is very sensitive to temperature changes, down to the level of single-photon energies. While very impressive, the device is limited to optical photons, destroys the photonic state, and is difficult to scale.

One strategy that in principle also allows for nondestructive measurements is based on quantum memory. An itinerant photon may be caught by an atom in a cavity Specht et al. (2011) or by a cavity with tunable coupling Yin et al. (2013), which however requires time-dependent tuning of the atom–cavity or cavity–waveguide coupling in accordance with the photon wavepacket shape. More generally, electromagnetically-induced transparency (EIT) can be used to slow down a light pulse such that it fits within an atomic cloud Fleischhauer and Lukin (2000); Fleischhauer et al. (2005). This allows the storage of arbitrary photon pulses within a length given by system parameters. The number of polaritons can in principle be read out by employing a cycling transition Imamoğlu (2002); James and Kwiat (2002). However, this scheme also requires knowledge about the arrival of the photon wavepacket. Furthermore, the combination of quantum memory and quantum nondemolition detection of one or few excitations in a 3D cloud of atoms makes this very challenging. EIT in a waveguide QED setting is one strategy to alleviate the imaging problem (see, e.g., Ref. Caneva et al. (2015)).

Here, we instead design a detector that does not require knowledge about the shape or arrival time of incoming wavepackets and therefore does not require implementation of a time-dependent Hamiltonian. This can be achieved either through continuous measurement, for example through dispersive coupling Chen et al. (2011); Govia et al. (2012); Poudel et al. (2012); Kyriienko and Sørensen (2016); Oelsner et al. (2017); Schöndorf et al. (2018), or if the detector permanently changes its state and is read out later, as in impedance-matched $\Lambda$-systems Pinotsi and Imamoglu (2008); Romero et al. (2009a, b); Peropadre et al. (2011); Koshino et al. (2013); Inomata et al. (2016); Bechler et al. (2018). The former class of detectors suffer from measurement backaction on the photonic state Helmer et al. (2009); Royer et al. (2018), whereas the latter does not produce ‘clicks’ to indicate detection event, and therefore does not localize photons in the waveguide, which avoids measurement backaction. The advantage of this strategy is that no knowledge about the photon wavepacket is needed, and that photons do not have to be destroyed to be detected, which opens the way toward nondestructive photon-number measurements.

In recent years, quantum emitters coupled to waveguides have emerged as a powerful experimental paradigm Chang et al. (2018); Türschmann et al. (2019). Examples of such systems include cold atoms levitated near optical fibres Vetsch et al. (2010); Douglas et al. (2015) or photonic crystal waveguides Goban et al. (2015), but also solid-state realizations such as quantum dots Akimov et al. (2007); Lodahl et al. (2015), superconducting qubits van Loo et al. (2013); Sundaresan et al. (2019), or nitrogen-vacancy centres Huck et al. (2011); Sipahigil et al. (2016); Evans et al. (2018). Due to the strong confinement of light, even a single emitter can have a profound effect on light propagation, and many emitters show remarkable collective effects Solano et al. (2017). At the same time, impressive experimental breakthroughs have enabled the control Schlosser et al. (2001) and, in particular, readout of neutral atoms 111In an array of 160 atoms, a recent experiment achieved $99.94\%$ readout fidelity, albeit at the expense of a constant atom loss rate that dominates this theoretical fidelity Wu et al. (2019). Atoms trapped in tweezer arrays can circumvent this problem, offering similar fidelity without losing the atoms Covey et al. (2019)., superconducting qubits Jeffrey et al. (2014); Walter et al. (2017), and trapped ions Harty et al. (2014). Thus, quantum emitters coupled to waveguides appear to be an ideal platform to produce Goban et al. (2015); Paulisch et al. (2018, 2019); Perarnau-Llobet et al. (2020) and control Zheng et al. (2013); Zhang and Mølmer (2019) quantum light. Here we show that they also allow for detection of quantum light.

The key idea here is to engineer atoms such that for each incident photon in the waveguide, exactly one atom changes its internal state, such that a subsequent measurement of the atomic state yields the number of photons in the scattered wavepacket. We do this by identifying conditions such that all photons are absorbed in one atomic transition ($g\to e$, blue in Fig. 1) and dissipated in a different one ($e\to s$, green). This way, the atomic array keeps a memory of the number of scattered photons Pinotsi and Imamoglu (2008); Romero et al. (2009a, b); Peropadre et al. (2011); Koshino et al. (2013); Inomata et al. (2016).

In order to identify such conditions, we study a realistic model for photon absorption by atomic arrays in the presence of free-space decay (red) and disorder. In view of experimental realizations, we provide recipes for engineering an additional decay channel, and pay particular attention to spatial disorder, as this fundamentally modifies the eigenmodes of the system. Surprisingly, we find that it is possible to engineer dissipation to achieve full absorption independent of disorder. As a consequence of collective enhancement in the atom-waveguide coupling, detection efficiency and bandwidth grow with atom number, such that even with moderate Purcell factor (ratio of waveguide to free-space decay rate), detection efficiencies $\eta>99\%$ can be reached for intermediate numbers of atoms ($N>20$). We also consider emitters coupled to chiral waveguides, which have distinct advantages due to the natural suppression of backscattering. Our detection scheme is scalable, as the number of atoms required to detect a certain number of photons with a fixed error grows only polynomially.

In a second part, we study a natural modification, in which the dissipated photons are emitted back into the waveguide (cf. Fig. 1d). If the the resulting output wavepacket coincides with the input wavepacket, this realizes a quantum nondemolition (QND) measurement 222We use the terminology adopted by most authors Thorne et al. (1978); Unruh (1978); Braginsky et al. (1980); Caves et al. (1980); Nogues et al. (1999); Raimond et al. (2001); Guerlin et al. (2007); Wiseman and Milburn (2009); Johnson et al. (2010), but note that there is some controversy whether this is an appropriate name Monroe (2011), and sometimes nondestructive is used instead Reiserer et al. (2013), but this term lacks the precise definition given by, e.g. Wiseman and Milburn Wiseman and Milburn (2009). . In order to achieve high fidelities, this requires collective enhancement of both decay channels, which requires the atomic array to be prepared in a symmetric Dicke state between ground states. We show that such states can be prepared through coherent interaction of two arrays and subsequent conditioning. Ultimately, this can reach a scaling where for a given error probability, the number of atoms that can be prepared in this way scales exponentially with Purcell factor $N\sim\exp(P)$, which we verify with numerical simulations.

The rest of this paper is structured as follows. In Section II we present the scattering theory for atomic arrays on waveguides. Building on the principle of coherent perfect absorption, we show explicitly that perfect absorption can be obtained by tuning only dissipation, even in the presence of arbitrary disorder. We apply this principle to establish the feasibility of number-resolving detection by absorption in a mirror geometry (Section III.1) and in an infinite waveguide (Section III.2). For completeness, we also discuss how dissipation can be used in a chiral waveguide to obtain a number-resolving detector (Section III.3). We extend these concepts to a QND measurement in Section IV. In Section V, we discuss the experimental implementations of our proposal, including engineered dissipation, atomic species, readout, dector reset, the effect of non-idealities, and bandwidth. We conclude in Section VI.

II Absorption in atomic arrays

Below, we first formulate a generic input-output theory for photons scattering off atomic arrays (Section II.1), following previous work Kien and Hakuta (2008); Le Kien and Rauschenbeutel (2014); Caneva et al. (2015). This allows us to generically identify the conditions under which all photons are absorbed in one transition and emitted in the other (Section II.2).

II.1 Scattering theory

A generic Hamiltonian describing the system-waveguide interaction is given through Kien and Hakuta (2008); Le Kien and Rauschenbeutel (2014); Caneva et al. (2015)

[TABLE]

In Eq. 1, the label $\nu$ runs over different sets of waveguide modes, the index $i$ denotes whether the waveguide mode couples to the $g\leftrightarrow e$ transition ($i=1$) or the $s\leftrightarrow e$ transition ($i=2$), the index $n$ runs of the $N$ atoms, the individual waveguide modes are labelled by their wavevector $k$, and the coupling of them to the atoms is given by the rate $g_{\nu,k,i}^{(n)}$.

This Hamiltonian describes a situation where the two transitions are coupled to different fields, which can either be waveguide modes or free-space modes. In an infinite waveguide (Fig. 1b), $g_{\nu,k,1}^{(n)}=\sqrt{2c\Gamma_{g}}\exp(ik\nu x_{n})$ and $\nu\in\pm$, corresponding to left- and right-moving modes, whereas for the semi-infinite waveguide (Fig. 1a), there is only one set of waveguide modes with coupling $g_{k,1}^{(n)}=\sqrt{c\Gamma_{g}}\sin(kx_{n})$. In these expressions, $\Gamma_{g}$ is the decay rate to $|g\rangle$ of an individual atom into the waveguide, and $\omega_{k}=ck$. Integrating out the bath modes yields quantum Langevin equations for the spin operators (see Appendix A),

[TABLE]

where the sum over input fields $\nu$ and atoms $m$ is implied. The non-Hermitian Hamiltonians $\mathsf{H}_{\mathrm{eff},1}$ ($\mathsf{H}_{\mathrm{eff},2}$) describes both coherent interaction of the quantum emitters and decay into the waveguide induced by the coupling of the $g\leftrightarrow e$ ($s\leftrightarrow e$) transition to the waveguide, and $\mathbf{a}_{\mathrm{in},\nu,i}$ are waveguide fields coupling to these transitions. This coupling is given by the generically non-square matrix $\mathsf{L}_{i}$. If coupled via an infinite waveguide, the $i^{\mathrm{th}}$ transition couples to two input fields, $\mathbf{a}_{\mathrm{in},i}=(a_{\mathrm{in,}+,i},a_{\mathrm{in},-,i})$, corresponding to right- and left-moving photons, such that $\mathsf{L}_{i}$ is a $N\times 2$ matrix. This is in contrast to a semi-infinite waveguide, which only has on input field, such that $\mathsf{L}_{i}$ is a vector. Independent free-space decay at rate $\Gamma_{\mathrm{free},i}$, can also be captured by adding a term $\mathsf{H}_{\mathrm{free},i}=-i\Gamma_{\mathrm{free},i}\mathbbm{1}/2$ $\mathsf{L}_{3}=\sqrt{\Gamma_{\mathrm{free}}}\mathbbm{1}$, and introducing $N$ independent noise operators $\mathbf{a}_{\mathrm{in,free},i}$. Note that the assumption of independent decay into free-space does not necessarily hold and depends on how closely spaced the quantum emitters along the waveguide are. Such a situation may lead to a reduction in free-space decay Asenjo-Garcia et al. (2017), which would benefit the detector proposed here. Multiple baths coupling to the same transition appear as several terms taking either the form of the first or the second line on the right-hand side.

The Langevin equation for $\sigma_{se}^{(n)}$ can be obtained by exchanging both $1\leftrightarrow 2$ and $g\leftrightarrow s$ everywhere in Eq. 2. The explicit form of $\mathsf{H}_{\mathrm{eff}}$ for the mirror geometry and the infinite waveguide are given in Eqs. 7 and 10 below.

If there are only few excitations compared to the number of atoms, one can linearize the Langevin equation using a Holstein-Primakoff transformation that sends $\sigma_{ge}\to b$ Gardiner and Zoller (2004); Chang et al. (2012); Caneva et al. (2015); Chang et al. (2018). In this approximation, $\sigma_{gs}^{(n)}\sigma_{se}^{(m)}\to\delta_{mn}\sigma_{ge}^{(n)}$, such that only the diagonal terms in the second line of Eq. 2 survive $\mathsf{H}_{\mathrm{eff},2}\to-i\Gamma_{s}\mathbbm{1}/2$ (green in Fig. 1), independent of whether the decay $e\to s$ corresponds to guided or non-guided modes. Combining this with free-space decay (shown in red in Fig. 1), we obtain uniform incoherent decay at rate $\Gamma^{\prime}=\Gamma_{\mathrm{free},1}+\Gamma_{s}$. We thus arrive at

[TABLE]

Here, we have neglected all input fields except the ones pertaining to the waveguide, which is valid if they are in vacuum. In order to describe destructive photon measurements, it is then sufficient to show that all incoming photons are absorbed and re-emitted into the bath modes coupling to the $s\leftrightarrow e$ transition. Strictly speaking, decay to $|s\rangle$ eliminates the atom from the dynamics, but this effect is neglected in the linearization.

II.2 Complete absorption

We now examine in general how tuning $\Gamma^{\prime}$ in Eq. 3 can lead to complete absorption of photons by an atomic array. For the rest of this section, we drop the indices from $\mathsf{H}_{\mathrm{eff}}$ and $\mathsf{L}$, for sake of generality and simplicity. In order to find the linear scattering properties of the array, we use the input-output equations that relate the output field operators to the input fields $\mathbf{a}_{\mathrm{out}}(t)=\mathbf{a}_{\mathrm{in}}(t)-\mathsf{L}^{\dagger}\mathbf{b}(t)$ Gardiner and Zoller (2004). Solving Eq. 3, we obtain the scattering matrix

[TABLE]

A detector that counts the number of photons in a specific input port (say, $\mathbf{a}_{\mathrm{in},+}$) needs to absorb all of them and dissipate them via the transition to $|s\rangle$. This is captured by our key figure of merit, the detection efficiency $\eta$, which is the product of the probability that a photon is not reflected $p_{\mathrm{abs}}=1-\sum_{\nu\neq 1}|\mathsf{S}_{\nu 1}(\omega)|^{2}$, and the probability that it is dissipated via the engineered channel (rate $\Gamma_{s}$) rather than into free space, $\eta=p_{\mathrm{abs}}\Gamma_{s}/(\Gamma_{s}+\Gamma_{\mathrm{free}})$. Thus, a high fidelity requires $\Gamma_{s}\gg\Gamma_{\mathrm{free}}$ and $p_{\mathrm{abs}}\approx 1$. In waveguide QED this regime can be reached through the collective enhancement of the emitter–waveguide coupling as compared to the free-space decay.

Unity absorption is attained if one of the eigenvalues of the scattering matrix $\mathsf{S}(\omega)$ is zero. This corresponds to a pole of the inverse scattering matrix $\mathsf{S}^{-1}=[\mathsf{1}+\mathsf{L}^{\dagger}(-i\omega+i\mathsf{H}_{\mathrm{eff}}^{\dagger}+\Gamma^{\prime}/2)^{-1}\mathsf{L}]$. A pole of $\mathsf{S}^{-1}$ arises whenever $\omega$ coincides with an eigenvalue of $\mathsf{H}_{\mathrm{eff}}^{\dagger}-i\Gamma^{\prime}/2$, which implies that the scattering matrix $\mathsf{S}$ has a zero if $\omega-i\Gamma^{\prime}/2$ coincides with an eigenvalue of $\mathsf{H}_{\mathrm{eff}}$. Thus, tuning $\omega$ and $\Gamma^{\prime}$, this can always be achieved, for any $\mathsf{H}_{\mathrm{eff}}$. Absorption based on this principle has been observed in a variety of systems Yan et al. (1989); Kishino et al. (1991); Cai et al. (2000), and been termed coherent perfect absorption Chong et al. (2010); Baranov et al. (2017). Note that to reach this conclusion we did not have to assume anything about the form of $\mathsf{H}_{\mathrm{eff}}$ (apart from linearity), which is the reason it works for arbitrary disordered systems. If these conditions are fulfilled, there exists an eigenvector $\mathbf{e}_{0}$, such that $\mathsf{S}(\omega_{0})\mathbf{e}_{0}=0$. Generically, $\mathbf{e}_{0}$ describes a linear combination of various input fields of the system, which implies that coherent perfect absorption is only a sufficient condition for unity absorption property when there is only one input field, such as in the mirror geometry. Nevertheless, in the infinite waveguide efficient detection is still attained for large atom numbers.

One can verify explicitly that $\mathsf{S}^{-1}$ is the inverse scattering matrix (here, $\mathsf{\Gamma}^{\prime}$ is a matrix, for generality)

[TABLE]

where

[TABLE]

The square brackets in Eq. 5 vanish, which can be checked explicitly for the two examples below. It can also be shown to hold generically, since if $\mathsf{\Gamma}^{\prime}=0$, the scattering matrix is unitary $\mathsf{S}^{-1}=\mathsf{S}^{\dagger}$, which is only true if the term in square brackets vanishes.

Physically, this can be interpreted as the fluctuation-dissipation theorem, as the anti-Hermitian part of $\mathsf{H}_{\mathrm{eff}}$ specifies the damping, whereas $\mathsf{L}$ captures how strongly the modes are coupled to the input noise operator.

III Destructive photon counting

III.1 Mirror geometry

We now turn to the specific model that best illustrates these ideas, namely an array of atoms coupled to a waveguide terminated by a mirror as sketched in Fig. 1a. In this geometry, the non-Hermitian Hamiltonian induced by integrating out the waveguide photons reads

[TABLE]

where $\Gamma_{g}$ is single-atom decay rate for the transition $e\to g$ into the waveguide, $x_{n}$ the position of the $n$th atom, and $k_{0}$ is the wavevector of the emitted light (wavelength $\lambda=2\pi/k_{0}$). The coupling of the atoms via the waveguide contains both a term due to photons travelling directly in between them, accumulating a phase $k_{0}|x_{m}-x_{n}|$, and one mediated by photons being reflected from the mirror, which incurs a minus sign and a phase $k_{0}(x_{m}+x_{n})$. Since there is only one input and output field (cf. Fig. 1a), the matrix $\mathsf{L}$ is now a vector $L_{n}=\sqrt{\Gamma_{g}}\sin(k_{0}x_{n})$.

It is instructive to see an example of how perfect absorption manifests in this setup. Placing the atoms in the atomic mirror configuration at positions $x_{n}=(1/4+n)\lambda$ (due to the infinite-range interactions, the lattice need not have unity filling), the photonic field only couples to the symmetric collective atomic excitation $B=\sum_{n}b_{n}/\sqrt{N}$, which also is an eigenmode of the atomic array. All other modes are dark and do not participate in the dynamics. In terms of this collective mode, the governing equations reduce to the input-output equations for a one-sided cavity Collett and Gardiner (1984) with internal dissipation

[TABLE]

where we have introduced the total decay rate $\Gamma_{\mathrm{tot}}=N\Gamma_{g}+\Gamma^{\prime}$. As in our discussion above, $\Gamma^{\prime}=\Gamma_{\mathrm{free}}+\Gamma_{s}$ comprises both free-space decay and the decay from $e\to s$, which is assumed to be tunable. Solving Eqs (8a, b) in frequency space, we calculate the number of photons in the output field

[TABLE]

If decay is tuned such that $\Gamma_{\mathrm{tot}}=N\Gamma_{g}+\Gamma^{\prime}=2N\Gamma_{g}$, there is perfect absorption on resonance ($p_{\mathrm{abs}}=1$), with a bandwidth of $2N\Gamma_{g}$, which in this single-mode picture is equivalent to impedance matching or critical coupling. In order to match the collective decay, the engineered decay $\Gamma_{s}\propto N$, such that as the atom number is increased, the detection efficiency $\eta=p_{\mathrm{abs}}\Gamma_{s}/(\Gamma_{s}+\Gamma_{\mathrm{free}})$ can become arbitrarily close to 1.

In Fig. 2a we include spatial disorder and show how the photon loss on resonance $p_{\mathrm{loss}}\equiv 1-\eta$ scales with atom number. Clearly, while the setup works very well for low spatial disorder ($\sigma/\lambda<1\%$), it suffers significantly from disorder. In the following we show how this is mitigated. Note that in Fig. 2 and indeed all plots in this paper we choose the Purcell factor $P=10$, which close to the current state-of-the-art Nayak et al. (2019). However, the collective enhancement of the atom-waveguide coupling means that the primary effect of having a lower Purcell factor is that more atoms have to be employed and therefore read out.

In the presence of disorder, the energies and decay rates of the eigenmodes of the atomic array are shifted, and many collective atomic modes couple to the input field. Thus, the picture presented above breaks down. Yet, as we have shown following Eq. 4, full absorption can be attained generically, independent of disorder, by tuning the engineered dissipation $\Gamma_{s}$ to one of the eigenmodes. However, the eigenmodes in the presence of disorder are not known a priori, so this approach appears infeasible. Surprisingly, one can still vastly improve over naively setting $\Gamma_{s}=N\Gamma_{g}$ as we did in Fig. 2a, if the standard deviation $\sigma\equiv\sqrt{\langle\delta x^{2}\rangle}$ of atomic positions is known. For a given $N,\sigma$, one can then calculate the average largest eigenvalue $\langle\mu\rangle_{N,\sigma}$ and tune the engineered dissipation to its imaginary part $\Gamma_{s}=\operatorname{Im}[\langle\mu\rangle_{\sigma,N}]$. Note that $\operatorname{Re}\langle\mu\rangle_{\sigma,N}=0$ if every configuration is as likely as its reflection. This restores the favourable scaling of detection efficiency with $N$, as illustrated by the blue ($\sigma=0.1\lambda$) and green ($\sigma=\lambda$) curves in Fig. 2b. Most strikingly, this works even in the presence of disorder equal to the lattice spacing (green), which is essentially equivalent to a fully random configuration. The reason it works lies in the fact that the largest eigenvalue and thus the absorption bandwidth grows faster than the fluctuations of the largest eigenvalue around its mean.

If the disorder is fixed as a result of fabrication, such as in solid-state implementations, one can further improve the scaling by measuring the largest decay rate. In this case, $\Gamma_{s}$ can be tuned exactly to the largest eigenvalue. This situation corresponds to the red curve in Fig. 2b, calculated for completely random configurations, where in each case $\Gamma_{s}$ was set to coincide with the imaginary part of the largest eigenvalue. This clearly leads to a better result than without exact tuning (green). The scaling is still worse than with low disorder (blue), because the largest eigenvalue grows slower in completely disordered configurations compared with mostly ordered ones.

So far, we have just discussed absorption and detection on resonance. Equally important is the detection bandwidth, given by the engineered decay rate. Since also the detection efficiency $\eta$ depends on the ratio between engineered dissipation and total decay rate, best detection is achieved when tuning to the most dissipative eigenmode (cf. Fig. 2c). This is our choice for all plots.

III.2 Infinite waveguide

Let us now turn to an atomic array coupled to an infinite waveguide, which has the simpler effective Hamiltonian

[TABLE]

since there is only one path for a photon to travel from one atom to the next. As illustrated in Fig. 1b, there are now two input and two output modes, a right-moving one ($+$) and a left-moving one ($-$) [cf. Eq. 3]. The atomic lowering operators couple to the input operators via the $N\times 2$ matrix $\mathsf{L}_{n,\nu}=\sqrt{\Gamma_{g}}\exp(ik_{0}\nu x_{n})$, where $\nu\in\{\pm 1\}$ labels right- and left-propagating modes. The scattering matrix reads

[TABLE]

where $\mathsf{M}$ is defined as above [Eq. 6]. Since $\mathsf{M}$ is symmetric, transmission of right- and left-moving waves (the diagonal elements of $\mathsf{S}$) are equal.

One can show that for atoms arranged in a periodic array, parity symmetry implies that reflection amplitude of right- and left-moving photons differs only by a phase. This is because $\mathsf{J}\mathsf{M}\mathsf{J}=\mathsf{M}$ under the action of the exchange matrix $\mathsf{J}$, which consists of ones on the anti-diagonal and otherwise zeros, and $\mathsf{J}\mathbf{L}_{+}=\exp(i\phi)\mathbf{L}_{-}$ for some $\phi$. In this case, we can parameterize the scattering matrix as

[TABLE]

with eigenvalues $\mu=A\pm B\exp(i\phi)$. In this system, coherent perfect absorption (one eigenvalue is zero) is equivalent to $B=\exp(i\theta)A$ and $\exp(i\phi+i\theta)=\pm 1$. In the end, full absorption of a uni-directional wavepacket may only be attained if $\mathsf{S}=0$. Thus, parity symmetry implies that perfect absorption may only occur if the scattering matrix is zero, corresponding to an exceptional point. Interestingly, as the atom number $N\to\infty$, the scattering matrix $\mathsf{S}\to 0$ for all arrays except the atomic-mirror configuration. In the latter, the scattering matrix can be shown to reduce to (for full absorption $\Gamma_{s}=N\Gamma_{g}$)

[TABLE]

which gives perfect absorption for wavepackets that are symmetric superpositions of left- and right-propagating modes, but not for wavepackets incident from one direction, an effect seen before Romero et al. (2009a); Fan et al. (2013). This can be overcome through atomic lenses Li et al. (2015), a mirror as above Peropadre et al. (2011), or by using any other lattice spacing Romero et al. (2009a). The latter can be deduced by estimating the magnitude of the elements of the scattering matrix through $|t|^{N}\sim\exp(-N^{1-\alpha})$, having assumed that the largest eigenvalue scales as $N^{\alpha}\Gamma_{g}$. Only in the atomic mirror configuration is $\alpha=1$, otherwise $\alpha<1$.

We numerically check the behaviour for a range of other spacings as well as fully disordered arrays and find similar behaviour as long as the atomic positions differ sufficiently from the atomic mirror configuration. In Fig. 3b,c we demonstrate this for a selection of array spacings (without including disorder). We choose $k_{0}=\pi/(2\lambda)$ for no particular reason, except that this appears to be a good choice. Including disorder, and employing the same technique of tuning to the average largest eigenvalue, we find similar scaling behaviour as in the mirror geometry, illustrated in Fig. 3a. The upshot is that arbitrary detection efficiencies can again be attained by increasing atom number, independent of disorder.

III.3 Chiral atom–waveguide coupling

Interestingly, the recently demonstrated platforms for chiral atom-waveguide coupling Lodahl et al. (2017) are another architecture in which robust photon detection may be achieved. Such a coupling is realized in a range of situation, for example when the light field is strongly confined Luxmoore et al. (2013); Junge et al. (2013); Shomroni et al. (2014); Söllner et al. (2015), when giant atoms are tuned to give a chiral coupling Kockum et al. (2018), or in topological systems Barik et al. (2018, 2020).

By design, (almost) no backscattering occurs in these systems, there are no collective effects, and the spacing of the atoms is immaterial, making the analysis straightforward. If an atom coupled to right- and left-moving modes at different rates, transmission and reflection on resonance are captured by Lodahl et al. (2017)

[TABLE]

As before, $\Gamma^{\prime}=\Gamma_{s}+\Gamma_{\mathrm{free}}$. This allows us to calculate the per-atom absorption probability (for the $+$-mode)

[TABLE]

The probability that the photon is dissipated in the right channel is $\Gamma_{s}/\Gamma^{\prime}$ as before. In the limit of many atoms, the transmission probability vanishes, as all photons are either reflected or absorbed. The detection efficiency is therefore the ratio of absorbed photons times the probability they are dissipated to $|s\rangle$, which to first order in $\gamma_{-}/\gamma_{+}$ is given by

[TABLE]

Note that even for moderate $\gamma_{-}/\gamma_{+}$, the second term can be reduced arbitrarily by increasing $\Gamma_{s}$, with the caveat that a larger number of emitters is needed before complete extinction is attained. In the absence of backscattering, this scheme is intrinsically robust against disorder. On top of that, the detection bandwidth depends only on the bandwidth of chirality and thus is—at least in principle—independent of $\Gamma_{s}$. This comes again with the caveat that photons far detuned from resonance on a scale of $\Gamma_{s}$ can only be absorbed with a large number of emitters.

IV Nondestructive photon counting

IV.1 Outline

In what we have discussed so far, we have disregarded the photons emitted via the engineered decay. After the state of the atomic array is measured, the only information obtained concerns the number of photons in the pulse. Since for each photon absorbed, one is emitted via the engineered channel, the question is pertinent whether a situation can arise in which the emitted photonic state coincides with the input state. This requires the outgoing photons to be disentangled from the atoms, save for in the collective number basis, and that the photonic state is not distorted in any other way. If these conditions are fulfilled, the setup realizes a quantum nondemolition (QND) photon-number measurement.

It is immediately obvious that if the photons are emitted into free space and thus scattered in all directions, the outgoing photonic state is a) useless, b) still entangled with the atoms, and c) distributed over many different modes. In principle, this is mitigated to some degree if the photons are emitted back into the waveguide. However, if each atom were to emit independently, the probability of any one photon getting lost is $1/(1+P)$, where $P$ is the Purcell factor, which severely limits the fidelity of the scattered wavepacket for realistic Purcell factors $P$. Furthermore, a photon has a different phase depending on the atom it is emitted from, causing residual entanglement of the photons with the atomic state, except in the atomic-mirror configuration. Thus, QND detection requires the mirror geometry (cf. Fig. 1d), which we study in detail below 333We note that even in the mirror geometry the phase of the emitted photon depends on the atom is emitted from if its frequency differs from the frequency of absorbed photons. This effect is negligible if $\delta\lambda/\lambda\ll\lambda/L$..

It turns out that the same trick that can make absorption robust against free-space decay and disorder—collective decay—can also be used to protect the re-emission into the waveguide. This is achieved if the atoms are in a superposition of $|g\rangle$ and $|s\rangle$ instead of all in $|g\rangle$. As we demonstrate mathematically below, this yields collective enhancement of the atom-waveguide coupling for both decay channels. This mode of operation therefore holds one further big advantage over the non-QND operation, namely that the bandwidth of the overall detector is not limited to the single-atom bandwidth. If the final readout of the atomic state is to give information about the number of photons, the ground state superposition must initially posses a known number of atoms in $|g\rangle$. These requirements are fulfilled by Dicke states, fully symmetric states with a definite number of excitations

[TABLE]

where $S_{sg}\equiv\sum_{i=1}^{N}\sigma_{sg}^{(i)}$ is a collective spin operator, and $|G\rangle$ is the state in which all atoms are in $|g\rangle$.

IV.2 QND measurement with Dicke states

We now aim to describe how a number-resolving QND measurement works in practice. We first focus on analytical results and approximations to motivate and illustrate the idea, and then corroborate the results with numerics. In the atomic-mirror configuration ($a=\lambda$) the coupling via a semi-infite waveguide (7) simplifies to $\mathsf{H}_{\mathrm{eff},1,mn}=-i\Gamma_{g}/2$. We also assume that the $e\to s$ transition couples to either the same waveguide at a slightly different frequency or to another waveguide field, such that $\mathsf{H}_{\mathrm{eff},2,mn}=-i\Gamma_{s}/2$. In this case, the Langevin equations Eq. 2 reduces to a single equation of motion for the collective spin operator

[TABLE]

and another one with $g\leftrightarrow s$ and $1\leftrightarrow 2$, where $S_{\alpha\beta}=\sum_{i}^{N}\sigma^{(i)}_{\alpha\beta}$ are collective spin operators. The input-output equations in this case read

[TABLE]

Note that neglecting the input fields corresponding to free-space decay in Eq. 18 is an approximation, as they take the system out of subspace of Dicke states. The reason is that when a photon decays into free-space, it destroys the coherence, which leaves the system in a mixed state. Therefore, Eq. 18 fails to account properly for the dynamics after the first photon has been dissipated to free space. The reason we can still use it to calculate the fidelity of QND measurements of Fock states is that once a photon is lost, the fidelity immediately drops to zero and remains zero, such that the subsequent time-evolution of the system is immaterial.

To understand how impedance matching can be attained here, consider the system being in the symmetric state with $m_{e}$ excitations $|\psi\rangle=|N/2-m_{0}-m_{\mathrm{e}},N/2+m_{0},m_{\mathrm{e}}\rangle$, where we have defined the imbalance $m_{0}=m-N/2$. Acting on this generic symmetric state,

[TABLE]

Since this is true independent of $m_{0}$ and $m_{\mathrm{e}}$, it is an operator identity, but only in the subspace of symmetric states. Another identity can be obtained by exchanging $s\leftrightarrow g$. This allows us to rewrite the equation of motion in the simplified form

[TABLE]

where $\xi_{\mathrm{in}}$ is the same input noise as in Eq. 18. In symmetric states with sufficiently large population in both ground states, we can replace the operators approximately by their expectation value in the initial state, $S_{gg}-S_{ee}\approx N/2-m_{0}$ and $S_{ss}\approx N/2+m_{0}$. Thus we conclude that

[TABLE]

ensures impedance matching, such that the decay rate of an excitation via the first and via the second channel are equal up to $\mathcal{O}(1/N)$. Under this condition, the equation of motion (21) on resonance can be solved approximately $S_{ge}\approx a_{\mathrm{in},1}/\sqrt{\Gamma_{g}}$. On the one hand, this implies that $a_{\mathrm{out},1}$ is independent of $a_{\mathrm{in},1}$, since the whole signal is absorbed. On the other hand, $a_{\mathrm{out},2}\approx-S_{sg}a_{\mathrm{in},1}\sqrt{2\Gamma_{s}/\Gamma_{g}N^{2}}$. Given this, $\langle a_{\mathrm{out},2}^{\dagger}a_{\mathrm{out},2}\rangle\approx\langle a_{\mathrm{in},1}^{\dagger}a_{\mathrm{in},1}\rangle$, i.e., every photon incident on port 1 is transmitted to port 2. The presence of $S_{sg}$ in this expression indicates that for each incoming photon, an atom is transferred from $g$ to $s$. This still holds for the second, third, …, $n^{\mathrm{th}}$ photon, but with an error of order $\mathcal{O}(n/N)$, which is the same as for the non-QND detector, up to constants of $\mathcal{O}(1)$.

To make this intuition quantitative, we follow a simple argument. First of all, let us denote the overlap of the output with the input wavepacket in case of the scattering of a single photon in a given state as $\mathcal{F}_{1}=|\langle\psi_{\mathrm{out}}|\psi_{\mathrm{in}}\rangle|^{2}$. For any finite-bandwidth wavepacket this differs from unity. For example, a single photon with a Gaussian wavefunction of width $\tau$ is transmitted by an impedance-matched cavity with outcoupling rate $\kappa$ and internal decay rate $\kappa_{\mathrm{int}}$ with fidelity

[TABLE]

where $\operatorname{erf}(x)=(2/\sqrt{\pi})\int_{0}^{x}\exp(-y^{2})dy$ is the error function. To make contact with our simulations, we thus set $\kappa=\Gamma_{g}(N/2-m_{0})$ and the internal cavity decay $\kappa_{\mathrm{int}}=\Gamma_{\mathrm{free}}$. In this limit of few excitations, it is therefore valid to think of free-space decay as limiting the fidelity by introducing a branching ratio between decay back into the waveguide and decay into free space, which is captured explicitly by the prefactor in Eq. 23.

Given some single-photon fidelity $\mathcal{F}_{1}$, a linear device would transmit $m_{\mathrm{p}}$ photons in the same mode with a fidelity of $\mathcal{F}_{1}^{m_{\mathrm{p}}}$. However, for each atom that transitions from $|g\rangle$ to $|s\rangle$, the single-photon collective decay rate from $e\to g$ is reduced by $\Gamma_{g}$, whereas the collective decay rate from $e\to s$ is increased by $\Gamma_{s}$, which leads to imperfect impedance matching. Using Eq. 9 to estimate the resulting additional reflection, we find the probability for absorbing the $k^{\mathrm{th}}$ photon is reduced by $1-2k^{2}/3N^{2}$. Thus, the probability that $m_{\mathrm{p}}$ photons are absorbed is reduced by a factor of $1-2(2m_{\mathrm{p}}^{3}-3m_{\mathrm{p}}^{2}+m_{\mathrm{p}})/3N^{2}$ to third order in $1/N$. Combining this with the single-photon fidelity $\mathcal{F}_{1}$, we can estimate the $m_{\mathrm{p}}$-photon QND fidelity as

[TABLE]

In the following, we compare these predictions with numerical simulation and find they agree well.

IV.3 Numerical simulation

In order to verify these conclusions numerically, we study the scattering of a multi-photon Fock states with a Gaussian wavepacket of width $\Gamma_{g}$ from an array of atoms in the atomic mirror configuration in the mirror geometry using a recently proposed technique Kiilerich and Mølmer (2019). We present details of the simulation in Appendix D.

In Fig. 4a we show the fidelity for a single photon scattering off an array starting from the initial state $|N/2,N/2,0\rangle$. Clearly, Eq. 23 is a good approximation to the transmission fidelity of the atomic array. This is still true for any other symmetric starting state, as we illustrate through Fig. 4b, which shows the fidelity of single-photon scattering when starting with an initial state $|N/2-m_{0},N/2+m_{0},0\rangle$, with the imbalance $m_{0}$ ranging from $-N/2$ to $N/2$. We assume there is a maximally achievable decay rate $\Gamma_{\mathrm{max}}\geq\Gamma_{g},\Gamma_{s}$ and consequently lower either $\Gamma_{g}$ or $\Gamma_{s}$ to fulfil the above condition (22). Together, these results show that single-photon transmission is captured very well by the above equations.

Turning to multi-photon scattering, we simulate several photons in the same Gaussian wavepacket scattering off the atomic array and again calculate the overlap of the output wavepacket with the input wavepacket. The results are shown in Fig. 4c. We find that our simple argument captures the fidelity well, and that our proposal can in principle reach very high fidelities for modest atom numbers. Finally, in Fig. 4d we show an example of the time evolution of the system.

IV.4 Dicke state preparation

Following similar arguments as in the other sections, the QND detector is robust against spatial disorder. However, the suppression of free-space decay crucially relies on the preparation of a Dicke state between the two ground states. One way to obtain such a state is to start with all atoms in the ground state $|g\rangle$, apply a $\pi/2$-pulse on the ground states $\left\{g,s\right\}$, and finally perform a projective measurement of $S_{gg}$ (or $S_{ss}$ or $S_{gg}-S_{ss}$). This heralds a fully symmetric state with a binomial distribution of imbalances around zero and standard deviation $\sqrt{N/4}$. However, such measurements are difficult. In principle, one can apply an off-resonant probe in the waveguide and recording the phase shift of the reflected light, which is proportional to the number of atoms Béguin et al. (2014). In Appendix B we briefly analyze this kind of measurement and find it is fundamentally limited to atom numbers of the order of the Purcell factor $N\lesssim P$. It has been proposed to produce atomic states by manipulating the dark-state manifold González-Tudela et al. (2015), which however is limited in fidelity by $1-\mathcal{F}\propto N/(2\sqrt{P})$, where $N$ is the number of atoms and $P$ the Purcell factor. Neither method scales well to many atoms. Thus, in the following we find a fast preparation method that imposes much less stringent requirements on the Purcell factor.

We propose to produce Dicke states in a way that is in keeping with the core idea of this article: measuring atoms, not photons. This requires two arrays (halves of one array) that are individually addressable with external driving. We further require that the atomic transition frequency be in the bandgap of the waveguide, which allows the atoms to be coupled coherently without dissipation. As shown elsewhere Douglas et al. (2015), such a setup readily gives rise to a Hamiltonian that couples all spins

[TABLE]

where the combination of coupling to individual waveguide modes $g_{k}\approx\text{const}$, bound state decay length $L=\sqrt{\alpha/\Delta}$, detuning of the impurity from the band edge $\Delta=\omega_{0}-\omega_{b}$, and band curvature $\omega_{k}=\omega_{b}+\alpha k^{2}$ yield the effective coupling strength $g_{\mathrm{eff}}=2\pi g_{k}^{2}/(\Delta L)$.

Importantly, when using a Raman transition ($s\to e\to g$) to couple the atoms to the waveguide modes, both $\Delta$ and $g_{k}$ in Eq. 25 are tunable. Careful analysis (cf. Appendix C) reveals that the effective Purcell factor $P_{\mathrm{eff}}=g_{\mathrm{eff}}/\Gamma_{\mathrm{eff,free}}$ scales as the ratio of the effective waveguide density of state ($1/\sqrt{\Delta\alpha}$) and the constant free-space density of state ($\rho_{0}$), viz., $P_{\mathrm{eff}}\propto 1/(\rho_{0}\sqrt{\alpha\Delta})$. The detuning from the band edge $\Delta$ can in principle be made arbitrarily small without violating the adiabatic condition or the Markov approximation (cf. Appendix C). Note also that disorder in the positions of the atoms is not an issue here, as there is no position-dependent phase 444 As discussed in Ref. Douglas et al. (2015), a dependence on position may arise, if the modes at the band edge, through which the atoms are coupled, have spatial structure, for example a finite wavevector. If the atoms are placed below the band edge of a waveguide, this does not happen as the lowest energy modes are extended modes at zero wavevector. . Instead this scheme is likely ultimately limited by disorder in the coupling strengths or energy, as they destroy the symmetry of the effective Hamiltonian.

The protocol to prepare an approximately half-excited state of one of the arrays is as follows. One array is fully excited (by applying a $\pi$-pulse) and the other is left in the ground state. The time evolution under the Hamiltonian (25), shown in Fig. 5a, transfers excitations from one chain to the other, while leaving them in their individual symmetric subspaces. Notably, after a short time $g_{\mathrm{eff}}t_{0}\sim 1/\sqrt{N}$, corresponding to the first zero crossing in Fig. 5a, the average number of excitations in each array is equal. However, this comes with the caveat that while on average the two arrays hold $N/2$ excitations each, in fact their imbalance has a very wide probability distribution, as we illustrate in Fig. 5c for $N=100$ emitters. As we have shown above, a large imbalance does not invalidate our scheme, but it does mean that the usable atom number on average is halved. Since this is a constant penalty, it does not change the overall scaling with atom number $N$. Ultimately, this scheme requires $P_{\mathrm{eff}}\gg\sqrt{N}$, where $P_{\mathrm{eff}}$ is the Purcell factor enhanced through the proximity to the band edge.

There is another, intrinsically faster, way to prepare Dicke states, if the system is governed by the Hamiltonian

[TABLE]

As we detail in Section C.2, this Hamiltonian may be engineered through the interference of interactions via two different waveguides (or bands in the same waveguide), and constitutes a waveguide QED version of spin flip-flops recently realized in cavity QED Davis et al. (2019). While certainly more challenging to implement experimentally, this Hamiltonian has the advantage of equilibrating the number of excitations in each array on an asymptotic time scale of $g_{\mathrm{eff}}t_{0}\sim\log(N)/(2\sqrt{2}N)$ (cf. Section C.3). This is the fastest time that can be achieved for a given $g_{\mathrm{eff}}$, essentially saturating the time scale obtained by adding all average transition times $\sum_{n}^{N}||S_{+}^{(2)}S_{-}^{(1)}|N,n\rangle\otimes|N,N-n\rangle||^{-1}\sim\log(N)/2N$. As a result, using the Hamiltonian Eq. 26 reduces the requirement on the effective Purcell factor to $P_{\mathrm{eff}}\gg\log(N)$. Another advantage of this Hamiltonian is illustrated by Fig. 6c, which shows the probability distribution of Dicke states after a time $t_{0}$. The overall probability for the state to be close to $|N/2,N/2\rangle$ is larger as with the other Hamiltonian (cf. Fig. 5c).

V Experimental considerations

A tacit assumption in the preceding sections has been that the decay rates from $|e\rangle$ to $|g\rangle$ and/or $|s\rangle$ are tunable. In circuit quantum electrodynamics, decay rates can be tuned by changing the detuning of an intermediate resonator Pechal et al. (2014). In atomic systems, this can be done with Raman transitions, which we analyze in the following for the D2 line of 87Rb. For concreteness, we assume here that the waveguide efficiently couples to $\pi$ transitions.

V.1 Engineered decay: Destructive photon measurement

The level scheme we consider specifically is drawn in Fig. 7a. In it, we make the choice $|g\rangle\equiv|F=2,m_{F}=2\rangle$, $|e\rangle\equiv|1,1\rangle$, and $\{|s\rangle_{i}\}=\{|2,1\rangle,|2,-1\rangle,|1,0\rangle,|1,-1\rangle\}$ (green area), which are coupled via excited states (in $5^{2}$P3/2) $|f_{1}\rangle\equiv|2,2\rangle$, $|f_{2}\rangle\equiv|2,0\rangle$. While this is just one choice among many it has the advantage that the applied lasers do not couple to any transition of the large number of atoms in the ground state $|g\rangle$. Since in the destructive scheme we only need to count how many atoms have scattered photons, the final state is irrelevant, provided it is not $|g\rangle$. The Raman scheme allows for large tuneability of the relative decay rate, which is required to obtain $\Gamma_{e\to s}(\Omega_{2})\approx N\Gamma_{e\to g}(\Omega_{1})$.

To show that the additional decays drawn in red do not spoil the scheme, we derive the effective quantum master equation governing the double-$\Lambda$ model in their presence. Neglecting the energy shifts due to the pumps, the dynamics are purely dissipative, given by the jump operators

[TABLE]

where in the last expression $\Gamma_{1,g_{1}}=\Gamma_{1,g}$ and $\Gamma_{2,g_{2}}=\Gamma_{2,s}$. We denote the (sum of the) rates corresponding to these jump operators $\Gamma_{e\to g}(\Omega_{1}),\Gamma_{e\to s}(\Omega_{2})$, and $\Gamma_{ee,\mathrm{eff}}$, respectively.

We note that the deleterious decays $\Gamma_{i,e}$ induce dephasing described by $\hat{L}_{ee,\mathrm{ee}}$ that is not negligible, primarily due to the $i=2$ term, as $\Omega_{2}\gg\Omega_{1}$. However, the dephasing removes the coherence and thus the superradiant decay to $|g\rangle$. Since $\Gamma_{e\to g}/\Gamma_{e\to s}=\mathcal{O}(1/N)$, the corresponding excitation decays to $|s\rangle$, with an error $\mathcal{O}(1/N)$. We conclude that all potential errors analyzed here are suppressed with increasing $N$.

V.2 Engineered decay: Nondestructive photon measurement

Similar considerations apply to the operation of the QND detector. A possible choice, shown in Fig. 7b contains ground states $|s\rangle=|2,-1\rangle$, $|g\rangle=|2,1\rangle$, and excited states $|f_{1}\rangle=|2,1\rangle$ and $|f_{2}\rangle=|2,-1\rangle$. Assuming we are able to prepare a Dicke state $|N_{g}=m,N_{s}=N-m,N_{e}=0\rangle$, we require the two collective decay rates to be similar [cf. Eq. 22].

Unlike the destructive scheme, now all other decays are deleterious. However, both decays $e\to s$ and $e\to g$ are enhanced by a factor of $N/2$. This implies that photon loss, either into free space or through decay into another hyperfine state scales as $1/N$. All deleterious decay rates can be taken together to form $\Gamma_{\mathrm{free}}$ in the calculations of Section IV.

V.3 Dicke state preparation

In order to prepare Dicke states, we propose to start with one array of the atoms in one hyperfine ground state $|s\rangle$ and the other array in another hyperfine ground state $|g\rangle$. In this protocol, the requirements for the Purcell factor (ratio of coupling strength to free-space decay) is most stringent. The most favourable level scheme is therefore a closed $\Lambda$-system, where $|g\rangle=|F=2,m_{F}=2\rangle$ and $|s\rangle=|2,1\rangle$ are hyperfine ground states, but $|f\rangle=|3,3\rangle$ is an excited state in $5^{2}$P3/2. In this case, $f\to g$ is a cycling transition and decays outside this subspace are suppressed. Coherent driving between $|s\rangle$ and $|f\rangle$ can be implemented using a two-photon transition Porras and Cirac (2008). This leaves free-space decay as error source, which has been discussed in Section IV.4.

V.4 Other sources of disorder

We have neglected inhomogeneous broadening and atomic motion, which could be taken into account in the same way as positional disorder and do not modify our conclusions. Furthermore, the relative effect of disorder decreases as the number of atoms increases Romero et al. (2009b). On the other hand, fast atomic motion, analyzed in Appendix F, essentially only renormalizes the coupling strength of the atoms to the waveguide (thereby decreasing the Purcell factor), but is not a fundamental obstruction.

In typical experiments today, the atom number might fluctuate in unknown ways from one experiment to the next. Allowing some error in $N$ is like increasing the detuning error in $\Gamma_{\mathrm{eng}}$. As long as the relative error decreases with $N$, as it should, this does not affect the scaling with atom number.

V.5 Photon loss from coupling into the waveguide

It is challenging to couple photonic wavepackets travelling in free space into fibre. If this is done with $80\%$ fidelity, say, then the detector becomes already unreliable for three photons, as there is already a 50% chance of losing one in the first stage.

However, the detectors described here are also useful for entirely waveguide-based setups, which obviate the need for coupling free-space wavepackets into a fibre. Recent proposals show that atomic arrays make excellent sources of quantum light Paulisch et al. (2018, 2019); Perarnau-Llobet et al. (2020) especially in view towards quantum metrology. There also exist proposals for two-photon gates in this platform Zheng et al. (2013), and it has been found that such arrays coupled to waveguides support long-lived subradiant states that can be used for storage Zhang and Mølmer (2019); Albrecht et al. (2019) and also manipulation Paulisch et al. (2016) of quantum states of light.

V.6 Atom readout

After the photons have scattered, one needs to measure the number of atoms in one of the states, or ideally both. Readout of superconducting qubits is well studied due to the advent of quantum computers. Fidelities now reach $>99\%$ in realistic devices Jeffrey et al. (2014); Walter et al. (2017).

For neutral atoms, a range of techniques have been developed. Among the earliest was to employ a cycling transition between a ground state and some excited state, which in recent experiments has yielded fidelities of 98-99% Kwon et al. (2017); Martinez-Dorantes et al. (2017). A similar idea is to push one type of atoms out of the trap with a resonant drive and use a quantum gas microscope to detect occupied sites Nelson et al. (2007). Instead, one could “boil” one type of atom away through inelastic light scattering. This can in principle achieve extremely high fidelities (99.97%), but is limited by atom loss and has the drawback that the lattice has to be refilled. We note that detection of atom numbers has already been used in experiment Goban et al. (2015); Prasad et al. (2019). Likewise, a strong magnetic field Boll et al. (2016) or a state-dependent trap Wu et al. (2019) can be used to separate the states and image them afterwards. In these setups, imaging has been performed with 99.94% fidelity, but background gas collisions still cause atom loss Wu et al. (2019). One option to prevent this would be tweezer arrays Covey et al. (2019). While such capabilities have not yet been demonstrated for atomic arrays coupled to waveguides, and scattering from the fibre might make photon collection harder, a solution could be to move the optical lattice slowly away from the fibre before performing the imaging.

After successful detection, the detector can be reset by pumping the $s\to e$ transition, such that eventually all atoms decay to $|g\rangle$.

VI Conclusion

We have explored the use of arrays of quantum emitters coupled to waveguides for number-resolving photon detection. Paying particular heed to experimental limitations such as disorder and free-space decay, we have found that both can be overcome, leaving no fundamental limitation to the achievable detection efficiency. Moreover, we have shown that the same platform can also be used to perform QND number-resolving photon detection. To this end, we also propose a novel way to prepare Dicke states based on the interaction of two arrays and subsequent heralding by measuring the number of excitations in one of them.

In a nutshell, our proposal builds on four facts that together enable highly efficient detectors: (1), few-level systems allow for strong, projective measurements of their state due to their intrinsic nonlinearity, (2), nevertheless, sufficiently large ensembles of atoms are linear, (3), collective decay mitigates errors due to non-idealities, and (4), in linear systems one can always engineer dissipation to obtain complete absorption. We hope that the ideas outlined here will mark a step towards high-fidelity number-resolving photon detectors, both destructive and nondestructive.

Acknowledgements.

D.M. would like to thank Adam Smith, Clara Wanjura, and Petr Zapletal for enlightening discussions. D.M. and J.I.C. acknowledge funding from ERC Advanced Grant QENOCOBA under the EU Horizon 2020 program (Grant Agreement No. 742102).

Appendix A Derivation of Langevin equations

A.1 Two-level systems coupled to one semi-infinite waveguide field

We study quantum emitters coupled to a semi-infinite waveguide terminated at $x=0$ by a mirror. First assuming that the waveguide is also terminated on the other side after a length $L$, the bath eigenmodes have the wavefunction $\phi_{n}(x)=\sin(k_{n}x)\sqrt{2/L}$, where $k_{n}=\pi n/L$ for all natural $n$. Taking the length of the waveguide to infinity, we recover the Hamiltonian given in the main text [Eq. 1]

[TABLE]

Defining the sine transform and its inverse through

[TABLE]

we can write the field in the waveguide as

[TABLE]

Defined this way, the commutation relation $[\phi(x,t),\phi(x^{\prime},t)]=\delta(x-x^{\prime})$ (for positive $x$ only) implies canonical commutation relations $[a_{k}(t),a_{q}^{\dagger}(t)]=2\pi\delta(k-q)$. In terms of complex amplitudes, we have

[TABLE]

Solving the operator equations of motion

[TABLE]

yields

[TABLE]

This solution can be plugged into the equation of motion for $\sigma_{ge}^{(n)}$. We need the following integral

[TABLE]

The equation of motion for $\sigma_{ge}^{(n)}$ becomes

[TABLE]

We define the input field $a_{\mathrm{in}}(t)$ as the portion of the waveguide field that was at a position $x=ct$ at time $t=0$ and has since travelled all the way to the atoms. Thus, $a_{\mathrm{in}}(t)=-\sqrt{c/2}e^{i\omega_{0}t}a(ct,0)$, where the pre-factor is fixed by the commutation relations of $a_{\mathrm{in}}$, up to an arbitrary phase. This yields the Langevin equation

[TABLE]

As defined, $a_{\mathrm{in}}(t)$ is a slow variable, so if the dynamics of the system and the bandwidth of the input state around $\omega_{0}$ are slow compared to the time it takes for light to travel a distance $2x_{n}/c$, we can neglect the retardation, rendering our description Markovian. The same applies to the atomic lowering operators. Finally, we arrive at a time-local equation

[TABLE]

To calculate the output field, we take the integrated equation of motion for the light field and apply a sine transform. This is essentially the same as the right-hand side of the equation of motion for $\sigma_{ge}^{(n)}$, except evaluated at a different point in space. Choosing this point to be $x_{R}+\varepsilon$, i.e., a small distance to the right of the rightmost atom, and again neglecting retardation, we find

[TABLE]

Further choosing $\varepsilon$ such that $\sin[k_{0}(x_{R}+\varepsilon)]=1$, and defining $a_{\mathrm{out}}(t)=-\sqrt{c/2}e^{i\omega_{0}t}a(x_{R}+\varepsilon,t)$, we have

[TABLE]

Finally, let us define the decay rate $\Gamma_{g}=g^{2}/2c$.

[TABLE]

A.2 Three-level systems coupled to two semi-infinite waveguide fields

In the main text we consider (effective) three-level systems with state $|g\rangle,|e\rangle,|s\rangle$ coupled to two waveguide fields, which respectively couple to the transition $g\leftrightarrow e$ and $s\leftrightarrow e$. Starting from the Hamiltonian

[TABLE]

the derivation follows through as above, and yields

[TABLE]

and two more equations if $1\leftrightarrow 2$ and $g\leftrightarrow s$ are exchanged simultaneously. These equations reduce to Eq. 18 for atoms in the atomic-mirror configuration $k_{0}x_{n}=2\pi(n+1/4)$.

A.3 Two-level and three-level systems coupled to one and two infinite waveguide fields

The derivation for the infinite waveguide proceeds in much the same way and can be found elsewhere Caneva et al. (2015). A Hamiltonian that combines both bath couplings reads

[TABLE]

Here, $k_{L,i}$ are phases imparted by the laser. They have little effect on the detection of incoming light. There are two waveguide fields, $\alpha=1,2$, distinguished either in frequency, polarization, or by being in a different waveguide. As before, $\nu\in\{\pm 1\}$ labels right- and left-moving modes.

From Eq. 43, we can derive the bath equations of motion, integrate them up and Fourier transform them Caneva et al. (2015)

[TABLE]

For the other field, we have to exchange $1\leftrightarrow 2$ and $g\leftrightarrow s$. As above, we next derive the atomic equations of motion

[TABLE]

and replace the photon field [Eq. 44]

[TABLE]

Using the description above, the Langevin equations for two-level systems on an infinite waveguide read

[TABLE]

whereas the governing equations for three-level systems are

[TABLE]

and two more equations if $1\leftrightarrow 2$ and $g\leftrightarrow s$ are exchanged simultaneously.

To linearize these equations, we substitute $\sigma_{gg}\to 1$, $\sigma_{ee}\to 0$ and $\sigma_{ge}\to b$, where $[b,b^{\dagger}]=1$.

Appendix B Preparing a Dicke state with collective measurements

There is another way to prepare a Dicke state in our setup. First applying a $\pi/2$-pulse to all atoms, and then measuring $S_{gg}$, i.e., the number of atoms in $|g\rangle$, which projects the state onto a symmetric state with a definite number of excitations. For large $N$, $m$ has variance $\sqrt{N}$ around $N/2$, such that $|N/2-m|\ll N$.

A $\pi/2$ pulse can for example be produced by driving the two transitions $g\leftrightarrow e$ and $s\leftrightarrow e$ with lasers, each detuned by a large amount $\Delta$, or by applying a microwave tone. After the pulse, the atomic state becomes $|X\rangle=\prod_{i}(1/\sqrt{2})(1+\sigma_{sg}^{(i)})|G\rangle$, which could optionally be written as $|X\rangle=2^{-N/2}\sum_{m}\begin{pmatrix}N\\ m\end{pmatrix}|N-m,m,0\rangle$. Measuring $\hat{S}_{gg}$ thus probabilistically returns the state $|N-m,m,0\rangle$ where $m$ follows a binomial distribution.

A known method to measure the number of atoms is by measuring the phase shift of an off-resonant probe tone Béguin et al. (2014). In order to prevent excited atoms from decaying via emission of a free-space photon (essentially removing the photon from the symmetric state) or via emission of a photon in the $e\to s$ channel (thereby inducing a transition to another Dicke state), the probe has to be far detuned. However, the further detuned the probe is, the lower is the induced phase shift, thus requiring longer averaging times, which leads to a trade-off. The analysis below shows that this sort of measurement is limits the number of atoms to the Purcell factor.

If it is possible to turn off $\Gamma_{s}$, we only have to consider free-space decay. Using the input-output equations, we find the output field for a weak coherent probe with amplitude $\alpha_{\mathrm{in}}$ detuned by $\Delta\gg N_{g}\Gamma_{g}$

[TABLE]

where $N_{g}$ is the number of atoms in state $|g\rangle$. Thus, the per-photon phase shift is $\varphi\simeq-\Gamma_{g}/2\Delta$. The number of photons required to resolve single excitations therefore is $N_{\mathrm{p}}\propto 4\Delta^{2}/\Gamma_{g}^{2}$. On the other hand, for a given probe amplitude $\alpha_{\mathrm{in}}$ the probability of an atom to be excited is $\langle S_{ee}\rangle\simeq\Gamma_{g}N_{g}|\alpha_{\mathrm{in}}/2\Delta|^{2}$. Thus, during the time it takes to scatter $N_{\mathrm{p}}$ photons, $\Gamma_{\mathrm{free}}N_{g}/\Gamma_{g}=N_{g}/P$ photons are lost. Ultimately, this way of preparing Dicke states is thus limited to atom numbers lower than the Purcell factor. If $\Gamma_{s}$ cannot be set to zero, but instead is comparable to $\Gamma_{g}$, then this constitutes the dominant decay channel and the requirement on the Purcell factor becomes even more stringent $P>N_{g}N_{s}$.

Appendix C All-to-all interaction in the bandgap

C.1 Strong coupling

Here we give details on how using a Raman transition allows strong coupling of the atoms when brought close to the bandgap of a waveguide and discuss briefly some limitations.

Following Douglas et al. (2015), two-level systems interacting via the modes in the bandgap couple at a rate $g_{\mathrm{eff}}=2\pi g_{k}^{2}/\sqrt{\Delta\alpha}$, where $g_{k}$ is the coupling to the individual waveguide modes (assumed constant), $\Delta$ is the effective detuning of the renormalized transition frequency from the band edge and $\alpha$ is a parameter to characterize the band curvature at the band edge, via $\omega_{k}=\omega_{b}+\alpha k^{2}$. If $g_{k}$ is the effective coupling rate obtained from adiabatically eliminating an excited state $|e\rangle$ in a Raman transition, it takes the form $g_{k}=g_{k,0}(\Omega/\delta)$, where $\Omega$ is the pump Rabi frequency and $\delta$ is the detuning of the pump from the transition. Thus,

[TABLE]

In order for the Markov approximation and the adiabatic elimination to be valid, we require $\delta\gg\Omega$ and $\sqrt{\Delta\alpha}\gg g_{k}$, so $g_{\mathrm{eff}}$ will be slow. Fortunately, reducing $\Omega/\delta$ also reduces free-space decay in the same way, $\Gamma_{\mathrm{free,eff}}=\Gamma_{\mathrm{free}}(\Omega/\delta)^{2}$.

Since the atom-atom coupling additionally scales with the effective waveguide density of states $g_{\mathrm{eff}}\propto(\alpha\Delta)^{-1/2}$, the Purcell factor $P_{\mathrm{eff}}=g_{\mathrm{eff}}/\Gamma_{\mathrm{eff,free}}\propto(\alpha\Delta)^{-1/2}$ can be made very large, independent of the original Purcell factor, while preserving the Markov condition $\Delta\gg g$. As this analysis shows, physically this relies on increasing the effective waveguide density of states. At the same time, this has the effect of making the bound state extent and therefore the decay length of the induced interaction $L$ much larger than the extent of the atomic array, such that the all-to-all interaction on the right-hand side of Eq. 25 becomes a very good approximation.

C.2 Engineering $S_{+}^{(1)}S_{-}^{(2)}+\text{H.c.}$

In order to instead realize a Hamiltonian with the interaction $S_{+}^{(1)}S_{-}^{(2)}+\text{H.c.}$ between two atomic arrays, one needs to combine two bands of the same waveguide (or two waveguides), and additionally make use of the spatial profile of the induced atom-atom interaction. The full expression for the coupling induced between the atoms through a waveguide band is given through Douglas et al. (2015)

[TABLE]

which differs from Eq. 25 through the addition of the spatial profile $E_{k_{0}}(x)$ of the bound state induced by coupling an atom at position $x$, where $k_{0}$ is the wavevector at the bandgap. A simple model of a waveguide exhibits bands with a dispersion $\omega_{k}\sim\cos(ka)$. As such, there are two band edges at wavevectors $k_{0}=0$ and $k_{0}=\pi$. In more complex models the bandgaps might occur at different values of $k_{0}$, but these are generically not expected to all coincide, and certainly not in between two different waveguides. With an otherwise constant mode profile, we can approximate $E_{k_{0}}(x_{i})\approx e^{ik_{0}x_{i}}$. For simplicity we will take $k_{0}=\pi$, but this is by no means required.

Combining two waveguides (waveguide bands), one with positive detuning $\Delta$, the other with the opposite detuning $-\Delta$, and placing both arrays of atoms in atomic mirror configurations, but spaced by an odd half-integer-multiple of wavelengths from each other, the two waveguides induce the interaction

[TABLE]

as required. Note that we have taken $L$ to be larger than the whole configuration, and choosing $g_{\mathrm{eff}}$ to be the same for both waveguides, which is not unreasonable as they can be tuned with the Raman transition. Furthermore, we would like to note that the requirement of placing the atoms in two atomic mirror configurations spaced by an odd multiple of half a wavelength is not any more stringent that producing a single atomic mirror configuration, although it certainly is experimentally more challenging than the disordered configurations we have considered in other parts of the main text.

C.3 Interaction time scales

While the full Hilbert space of two arrays of $N$ spins each is $2^{2N}$ dimensional, the Hamiltonian Eq. 52 connects the fully symmetric state $|N,0\rangle=|N\rangle_{1}\otimes|0\rangle_{2}$, to only $N$ other states, $\{|N-m,m\rangle\}$. As in the main text, our convention here is that $|m\rangle$ denotes a fully symmetric spins state in which $m$ out of $N$ atoms are excited, whereas $|m_{1},m_{2}\rangle$ denotes two arrays of $N$ atoms each, with $m_{1}$ and $m_{2}$ excitations, respectively.

In this space of $N+1$ states, the Hamiltonian Eq. 52 is a matrix of the form

[TABLE]

where

[TABLE]

As an aside, we note that this is similar to the model studied in Ref. Christandl et al. (2004), except with the elements of the matrix squared. The connection arises since the model with $\sqrt{a_{m}}$ on the diagonal maps to the $N$-excitations subspace of two coupled harmonic oscillators $H=a^{\dagger}b+\mathrm{H.c.}$ (leading to perfect state transfer), while here we instead study two coupled spins $H_{2}=2g_{\mathrm{eff}}S_{+}^{(2)}S_{-}^{(1)}+\mathrm{H.c.}$

If the index $m=x$ is understood as a spatial coordinate, the Hamiltonian $H_{2}$ can be rewritten as

[TABLE]

where $[\hat{x},\hat{p}]=i$, such that $e^{i\hat{p}}$ generates a translation by $-1$. In order to get some insight into the dynamics of $H_{2}$, we consider the classical long-wavelength limit neglecting terms $\hat{p}^{n}$ with $n\geq 3$ and replacing $\hat{x}\to x,\hat{p}\to p$, which yields

[TABLE]

where

[TABLE]

This Hamiltonian gives rise to periodic trajectories. A more physical Hamiltonian would be obtained by sending $p^{2}\to-p^{2}$ and $H\to-H$, but the dynamics are the same.

Starting at $x(t=0)=p(t=0)=0$, the conserved energy of the particle is $H(0,0)=E_{0}=2g_{\mathrm{eff}}(3N-1)$. This allows us to solve for the momentum as a function of position

[TABLE]

which can be used to calculate the period of the orbit

[TABLE]

The upper limit of the integral is the point at which the particle turns around, which can be found from $p(x)=0$. For large $N$, the integral may be approximated through

[TABLE]

For completeness, we mention another way to arrive at this result. Instead of Eq. 52, consider the Hamiltonian

[TABLE]

where $a,b$ are the annihilation operators for two harmonic oscillators. In the subspace spanned by the states $|2N-2m\rangle_{1}\otimes|2m\rangle_{2}$, where now $\{|m\rangle\}$ are now Fock states of an oscillator, the Hamiltonian is again given by a matrix of the same form as Eq. 53, except that now

[TABLE]

For large $N,m$ these are essentially the same matrix elements as before, and one can check numerically that the dynamics is very similar for large enough $N$. In the classical limit, we replace the annihilation operators by the amplitudes of the corresponding coherent states. The classical mean-field equations of motion read

[TABLE]

with initial conditions $\beta(0)=1$, $\alpha(0)=\sqrt{2N}$.

Amplitude and phase degrees of freedom can be separated by change of variables to $\alpha=a\exp(ix)$, $\beta=b\exp(iy)$, which follow equations of motion

[TABLE]

where $\phi=2x-2y$ and we have redefined time to include the factor $2g_{\mathrm{eff}}$. Further defining the constant $K=a^{2}+b^{2}$, and the variable $r=a^{2}-b^{2}$, we obtain the equations of motion

[TABLE]

The shape of the resulting periodic orbits can be found by integrating $dr/d\phi$.

However, since the phase $\phi$ is ill-defined in the initial state, we can choose it to be $\phi=\pi/2$. In this case, the system is governed only by

[TABLE]

which is readily integrated to yield

[TABLE]

Since initially, $r/K\simeq 1-1/N$ and $K\simeq N$, we conclude that the natural timescale for the switch is

[TABLE]

the same result as before. We note that if another phase $\phi$ is chosen initially, $\phi$ first rapidly evolves to a value close to $\pi/2$, after which the system evolution is very similar.

Appendix D Numerical simulation of multi-photon scattering

In this Appendix we briefly outline the theory behind our simulation of multi-photon scattering, following Ref. Kiilerich and Mølmer (2019). The key idea is to the time-evolution of the system and bath for a specific input field using the Langevin equation (18) in combination with the input-output equations (19).

In general, an input wavepacket can be described through one or more modes of the waveguide in which photons are created on top of the vacuum $|0\rangle$,

[TABLE]

As a function of time, the field travels through the waveguide, and eventually the photons interact with the atoms locally. Instead of this spatiotemporal description, in which a quantum system couples with constant rate to the waveguide, but the wavefunction of the photons changes in space, it has been shown that the dynamics can be equivalently captured by modulating the coupling between the system and the input mode according to the mode shape Kiilerich and Mølmer (2019).

Specifically, here we assume that all $m_{p}$ photons of the input field are in one mode, such that the waveguide field can be written (at $t=0$, before the photons interact with the system)

[TABLE]

with

[TABLE]

Since this is assumed to be the wavepacket before the interaction, it has support only for negative $x$ (we assume the wavepacket is incident from the right in Fig. 1). After a time $t_{0}$, the centre of the wavepacket (travelling at speed $c$) has arrived at the system.

Replacing the spatiotemporal evolution of wavepackets in the waveguide with a time-dependent coupling constant, the input-output equations can be replaced by a quantum master equation Kiilerich and Mølmer (2019)

[TABLE]

where the Hamiltonian

[TABLE]

describes coupling to the input mode, and the jump operator

[TABLE]

describes the dissipation induced by the coupled waveguides. We determine the time-evolution given by Eq. 72 using QuTiP Johansson et al. (2013).

Note a subtlety when comparing with the equation corresponding to Eq. 73 derived by Kiilerich and Mølmer Kiilerich and Mølmer (2019), which features an additional term $g_{\mathrm{out}}^{*}(t)g_{\mathrm{in}}(t)a_{i}^{\dagger}a_{o}$. This term arises when calculating the scattering from one input mode to the corresponding output mode, which in our case might be $a_{\mathrm{in},1}\to a_{\mathrm{out},1}$ and can be thought of arising from the first term in the input-output equation Eq. 19. Here, we are interested in the scattering $a_{\mathrm{in},1}\to a_{\mathrm{out},2}$, where this term is absent.

In the above expressions, the variable couplings strengths are given by

[TABLE]

and

[TABLE]

Here, the initial time to start the simulation has been chosen to be $t=0$, and since photons are assumed to travel at speed $c$, the temporal wavefunction is related to the spatial wavefunction through $\tilde{\psi}_{\mathrm{in}}(x)=\psi_{\mathrm{in}}(-ct)$. For consistency, $g_{\mathrm{in}}(t)=g_{\mathrm{out}}(t)$ are taken to be zero before $t=0$, and $\psi_{\mathrm{in}}$ should be square-integrable for positive times $\int_{0}^{\infty}dt\,|\psi_{\mathrm{in}}(t)|^{2}=1$. This is true for our choice Eq. 71 as long as $\Gamma_{g}^{-1}\ll t_{0}$. Since $t_{0}$ is an arbitrary offset it can always be chosen to fulfil this condition.

If the system and the input field are initially in states of known excitation number, which we assume to be the case throughout, then it is clear that for $m_{p}$ photons in the input wavepacket, we need to consider at most $m_{p}$ excitations in the bosonic modes $a_{o}$ and $a_{i}$ as well as the system. Recall that we assume a symmetric Dicke state as starting point, $|\psi_{\mathrm{sys}}(0)\rangle=|N/2-m_{0},N/2,0\rangle$. As a result, the required Hilbert space has a dimension of $(m_{p}+1)^{3}$.

We note that this formalism does not specify what shape the output wavepacket $\psi_{\mathrm{out}}$ has. Indeed, it is generically not even true that it can be described through a single mode and neither is it true in general that all photons are emitted in the second channel rather than the first. The advantage of this formalism in our case is that the system dynamics are independent of the choice of $\psi_{\mathrm{out}}(t)$. Instead, if the set of output modes considered does not comprise all modes of the physical output field, the corresponding photons are lost. In the extreme case, when $g_{\mathrm{out}}(t)=0$, the QME Eq. 72 describes the system dynamics for a given input, but does not yield any information about the output photonic state.

Here, we use this property by setting the output mode equal to the input mode $\psi_{\mathrm{out}}(t)=-\psi_{\mathrm{in}}(t)$. This is equivalent to asking the question how many photons are transmitted without changing the shape of the wavepacket – i.e., performing a QND measurement. If the final number of excitations in mode $a_{o}$ is equal to the number of input photons $m_{p}$, all photons have been transmitted faithfully, and all photons have been dissipated via the second channel. This defines the fidelity $\mathcal{F}$, which we plot in Fig. 4. Mathematically, it is defined as

[TABLE]

where $|m_{p}\rangle$ is the $m_{p}^{\mathrm{th}}$-Fock state of the output mode $a_{o}$. This is the probability that the all photons have been transmitted in the specified output mode.

Appendix E Effect of finite Purcell factor on performance of QND detector

In this Appendix, we offer a side-to-side comparison of the fidelity with a finite Purcell factor, as shown in Fig. 4 in the main text, and infinite Purcell factor. This allows to discern which features come from free-space decay, and which from nonlinearities and finite bandwidth. It furthermore offers further verification of the predicted analytical approximate fidelity Eqs. 23 and 24 through simulations. The results are shown in Fig. 8 and are commented on the the caption.

Appendix F Fast thermal motion of atoms

In the main text we have only considered static disorder, which is valid for slowly moving atoms. If the thermal motion of atoms is fast, some of the effect of disorder will be averaged out. This can be modelled by instead averaging the atomic coupling over a distribution of atomic positions Porras and Cirac (2008). Assuming a Gaussian distributions of positions around the atomic mirror configuration, captured by the random variable $y_{m}$

[TABLE]

we can calculate the off-diagonal coupling

[TABLE]

If $m=n$, there should only be one integral, giving $\bar{g}_{nn}=-[1+\exp(-k_{0}^{2}\sigma^{2})]\Gamma_{g}/(4g)$. In the case $m\neq n$, we can straightforwardly evaluate the second term, which yields $-\Gamma_{g}e^{-k_{0}^{2}\sigma^{2}}/4$ overall. For the first term, we first shift $y_{m}\to y_{m}+y_{n}$, in which case the $y_{n}$ integral becomes a straightforward Gaussian giving a factor of $\sqrt{\pi\sigma^{2}}e^{y_{m}^{2}/4\sigma^{2}}$. The leftover integral reads

[TABLE]

Taken together, we get

[TABLE]

where $F(x)$ is the purely real Dawson integral. It is peaked at $x=1$, is odd, and obeys $|F(x)|<0.6$, and $F(x)\to 0^{\pm}$ as $x\to\pm\infty$. Equation 81 predicts that the coupling decreases exponentially in $(k_{0}\sigma)^{2}$. For low to moderate $k_{0}\sigma$, the effect of fast thermal motion can simply be accounted for by re-scaling the couplings, without affecting the conclusions in the main text.