# The effects of a passive Bi-Polar Grid (BPG) on
Ion Back-Flow (IBF) and Resolution

V. Zakharov  E. Shulga  P. Garg  T. Hemmick  A. Milov

###### Abstract

Time Projection Chambers (TPC)s are excellent tracking detectors for high
multiplicity events and can intrinsically be high-rate, but are limited by the
ions created in their avalanche stage. GEMs and Micromegas can reduce IBF
through their geometry and $\vec{E}$-field ratios, but these can lead to gain
fluctuations and still leave IBF as the dominant source of space charge. An
active BPG can block all IBF ions, but their slow drift speed creates too much
dead time. A passive BPG will overcome this limitation by using an external
$\vec{B}$-field to allow the electrons to pass through while still blocking
all ions. Since the grid changes the electron’s trajectory, a loss of
resolution will occur. The trajectory is shifted symmetrically along the wires
so the wire alignment with respect to the detection pads is a specific
question not studied before. We present completed IBF analysis from data
collected at Weizmann Institute of Science (WIS), along with an intro to our
test on wire resolution.

##  1 Gating Grids

A passive BPG with $(V_{g}\pm\Delta g)$ can provide $\vec{E}$-field ratios
without creating fluctuations and would eliminate dead time by constantly
being opaque to ions while exploiting an external magnetic $\vec{B}$-field to
remain transparent to electrons.

$m\frac{d\vec{v}}{dt}=q\vec{E}+q(\vec{v}\times\vec{B})-\kappa\vec{v}$ (1.1)

Charges follow the Langevin equation and lighter electrons will experience a
continuous $(\vec{v}\times\vec{B})$ “push" allowing them to pass between the
wires. Electrons coming from above the wire will get their 1st push along the
length of the wire: $\hat{F_{1}}=\pm\hat{x}\times\hat{z}=\mp\hat{y}$. Then
they get a 2nd weaker push which allows them to get past the wires:
$\hat{F_{2}}=\mp\hat{y}\times\hat{z}=\mp\hat{x}$. Fig. 1 shows these effects
along with the distortion both pinched and shifted in their y-position.
Acceptable electron transparency with 0% ion transparency has been
demonstrated by previous groups [1].

Figure 1: Left: Illustration of electrons drifting from opposite sides of a
positive BPG wire. Right: [2] Simulation close-up showing electrons in orange
lines making it to the green detection plane past the BPG. Blue lines show ion
paths. The black oval indicates the electron distortions.

The general equation for the charge velocity in a gas as a function of
$\vec{E}$ and $\vec{B}$ fields is given by

$\vec{v}=\frac{\mu}{1+\omega^{2}\tau^{2}}\left(\vec{E}+\frac{\vec{E}\times\vec{B}}{|\vec{B}|}\omega\tau+\dfrac{(\vec{E}\cdot\vec{B})\vec{B}}{B^{2}}\omega^{2}\tau^{2}\right)$
(1.2)

where $\mu=e\tau/m$ is the charge mobility with $\tau$ as the average time
between collisions. $\omega=eB/m$ is the cyclotron frequency. Defining angle
$\alpha$ between $\vec{E}$ and the positive x-direction gives
$\vec{E}=(E\cos{\alpha},0,E\sin{\alpha})$. Together with $\vec{B}=(0,0,B)$
leads to

$v_{x}=c(E\cos{\alpha})\qquad v_{y}=c(E\omega\tau\cos{\alpha})\qquad v_{z}=\mu
E\sin{\alpha}$ (1.3)

If the BPG wires aren’t uniformly spaced along their tensioned direction an
$E_{y}$ component will be present and the charge velocity equations become

$v_{x}=c(E_{y}\omega\tau+E_{x})\qquad v_{y}=c(-E_{x}\omega\tau+E_{y})\qquad
v_{z}=c(E_{y}\omega^{2}\tau^{2}+E_{z})$ (1.4)

where $c=\dfrac{\mu}{1+\omega^{2}\tau^{2}}$. In both cases, the velocities of
main interest for ion blocking and electron transparency $v_{x}$ and $v_{z}$
are strongly dependent on $\omega\tau$. For electrons $\omega\tau\gg 1$ and
the velocities will follow the $\vec{B}$-field and ignore the local
$\vec{E}$-field of the BPG. For ions the converse is true.

##  2 Weizmann IBF Measurements

Ion and electron transparencies, as a function of $\vec{B}$-field and $\Delta
V$ volts on a linear BPG, were measured at Weizmann. A detector with MWPC gain
wires was constructed. A strong 55Fe X-Ray source was used and picoAmmeters
were implemented to measure the small changes in current as the grid was
slowly closed to ions. A mesh was used to help create a uniform drift field,
as seen in Fig. 2, and to create a region where the $\vec{E}$-field above-to-
below ratios can be better exploited.

Figure 2: Illustration showing all elements of the WIS detector. Orange and
red lines are electrons and ions respectively.

Field ratios were varied measuring the electron and ion transparencies from
the BPG and mesh $T^{g,m}_{e,i}$. Increased grid ratios meant decreased mesh
ratios pulling more ions through it. $T_{e}^{g}<100\%$ requires an increased
gain, which creates more IBF. A Figure of Merit (FoM) was constructed to
quantify the effectiveness of the BPG [3]. 1st term comes from $\vec{E}$
ratios, where higher ratios can extract more ions from the gain region. 2nd
term compensates for loss of primaries creating more IBF ions. Fig. 3 shows an
example for transparencies and FoM.

$FoM(\omega,\Delta
V)=\frac{T^{m}_{i}(\omega,0)}{T^{m}_{i}(1,0)}\cdot\frac{T^{g}_{i}(\omega,\Delta
V)}{T^{g}_{e}(\omega,\Delta V}\qquad\omega=E_{t}/E_{d}$ (2.1)

(a) Ion and electron grid transparencies

(b) FoM demonstrating BPG effectiveness

Figure 3: Some IBF results from Weizmann. Full results: [3]

##  3 BPG Resolution Distortion

From the equations of motion, and as seen in Fig. 1, the electrons fall "out
of plane", and alternating wires "pinch" them together, hurting resolution.
The error is cyclical with the same periodicity as the wires. These recurring
shifts are known as Differential Non-Linearity (DNL), and can be corrected.
Segmented pads are used to spread out the electron cloud charge collection
with finer detail in order to get a more precise track position measurement.

Zig-zag pads [4], which themselves have DNL-corrected error due to their
shape, are used. From the grid’s symmetric shift, perhaps BPG’s distortions
will average out over the full particle track, regardless of the wire’s
position with respect to the radial zig-zags. However, if the grid is aligned
with the pads then their independent DNLs might combine for better resolution.
A BPG with 2 linear and 2 radial wire configurations is placed above the pad-
plane and goes inside a prototype TPC. We’ve taken it to Argonne National Lab
(ANL) to test this hypothesis (Fig. 4).

(a) Wires of four different configurations placed above a quad-GEM stack

(b) Prototype TPC about to enter $\vec{B}=3T$ at ANL

Figure 4: Recording position resolution data using cosmic rays passing through
a TPC.

##  4 Conclusion

The passive BPG source test demonstrates ion blocking with high electron
transparency. The ANL test showed excellent resolution results and a paper
with the full test description and analysis will soon be published.

## References

  * [1] S.R. Amendolia et al., _Influence of a magnetic Field on the gating of a TPC_ , _NIM-A_ 234 (1985) 47.
  * [2] V. Zakharov, P. Garg, T. Hemmick and K. Dehmelt, _Study of a Passive Gating Grid for Ion Back Flow Suppression_ , _J. Phys. Conf. Ser._ 1498 (2020) 012 [1912.05005].
  * [3] E. Shulga, V. Zakharov, P. Garg, T.K. Hemmick and A. Milov, _Measurement of the ion blocking by the passive bipolar grid_ , _IEEE Transactions on Nuclear Science_ 68 (2021) 59.
  * [4] B. Azmoun et al., _Design studies for a tpc readout plane using zigzag patterns with multistage gem detectors_ , _IEEE Transactions on Nuclear Science_ 65 (2018) 1416.


# Renormalization group analysis of near-field induced dephasing of optical
spin waves in an atomic medium

S. Grava1,2, Y. He3, S. Wu3 and D. E. Chang1,2 1ICFO-Institut de Ciencies
Fotoniques, The Barcelona Institute of Science and Technology, 08860
Castelldefels, Barcelona, Spain.  2ICREA-Institució Catalana de Recerca i
Estudis Avançats, 08015 Barcelona, Spain.  3Department of Physics, State Key
Laboratory of Surface Physics and Key Laboratory of Micro and Nano Photonic
Structures (Ministry of Education), Fudan University, Shanghai 200433, China.
<EMAIL_ADDRESS>
###### Abstract

While typical theories of atom-light interactions treat the atomic medium as
being smooth, it is well-known that microscopic optical effects driven by
atomic granularity, dipole-dipole interactions, and multiple scattering can
lead to important effects. Recently, for example, it was experimentally
observed that these ingredients can lead to a fundamental, density-dependent
dephasing of optical spin waves in a disordered atomic medium. Here, we go
beyond the short-time and dilute limits considered previously, to develop a
comprehensive theory of dephasing dynamics for arbitrary times and atomic
densities. In particular, we develop a novel, non-perturbative theory based on
strong disorder renormalization group, in order to quantitatively predict the
dominant role that near-field optical interactions between nearby neighbors
has in driving the dephasing process. This theory also enables one to capture
the key features of the many-atom dephasing dynamics in terms of an effective
single-atom model. These results should shed light on the limits imposed by
near-field interactions on quantum optical phenomena in dense atomic media,
and illustrate the promise of strong disorder renormalization group as a
method of dealing with complex microscopic optical phenomena in such systems.

##  1 Introduction

The interaction of light with atomic ensembles provides the basis for numerous
potential applications, such as quantum memories for light [1, 2], quantum
nonlinear optics with strong photon-photon interactions [3, 4, 5], and quantum
metrology [6, 7, 8, 9, 10, 11]. In order to avoid the complexity associated
with the large microscopic number of degrees of freedom, such as the large
atom number and their positions, our standard theories for such systems
typically favor a macroscopic approach. For example, for atom-light
interactions in free space, the Maxwell-Bloch equations (MBE) [12, 13, 14, 15,
16] treat the atoms as a smooth polarizable quantum medium. The MBE have
yielded many important insights into the physics that enables the applications
above, as well as elucidating performance limitations [12, 17, 18].

Beyond macroscopic phenomena, many microscopic optical effects driven by
granularity, dipole-dipole interactions and multiple scattering have been
predicted, such as modifications of refractive indices and scattering rates
[19, 20, 21, 22, 23, 24], subradiance [25, 26], and coherent back-scattering
[27, 28, 29]. Besides being of foundational interest, such microscopic effects
could also have practical consequences on applications. For example, recently
it was experimentally shown in Ref. [30] that such effects lead to a
fundamental inhomogeneous broadening of optical transitions in an ensemble.
This manifests itself as an additional dephasing on top of spontaneous
emission decay for optical spin waves, with a rate that is exponential at
early times and is directly proportional to atomic density. It was argued that
this initial dephasing arises from the strong near-field interaction of a
small fraction of particularly close nearest neighbors, quantitatively
reproducing the experimental results. Separately, though, one might wonder
what governs the apparently non-exponential behavior at later times, or what
occurs at very high densities, when many atoms sit within a wavelength of each
other and experience strong near-field interactions. We also note that near-
field interactions have been recognized to play key roles in other collective
behavior, ranging from the modification of superradiance in small systems [31]
to late-time subradiance [32] in extended systems. Beyond exact numerics,
however, development of effective theories generally remains a challenge in
many-atom disordered systems.

Here, we provide a comprehensive theoretical picture of the spin-wave
dephasing phenomenon by applying a non-perturbative technique based on strong
disorder renormalization group (RG), which is a powerful method to elucidate
the physics in diverse disordered condensed matter systems [33, 34, 35, 36,
37, 38, 39, 40] and has also recently been applied to atom-light interactions
[41]. As in the short-time theory of dephasing, one key idea underlying this
approach is that for highly disordered atomic media, strong near-field
interactions between particularly close nearest neighbors allow such pairs to
be approximately diagonalized first. The resulting dynamics is equivalent to
replacing the pair with two, new effective atoms with renormalized
frequencies, as illustrated in Fig. 1 (a). The RG theory goes significantly
beyond this, however, by realizing that nearby, strongly interacting pairs
(including atoms previously renormalized) can continue to be identified and
diagonalized, i.e. the many-atom system interacting via the near field can be
thought of and diagonalized in terms of an extended hierarchy of strongly
interacting pairs (Fig. 1 (b)). The final result is that the original system
is optically equivalent to an inhomogeneously broadened medium with a well-
defined distribution of resonance frequencies $P(\omega)$, and with the strong
near-field interactions effectively removed. This approach was recently used
to predict that a disordered atomic medium has a limiting value of maximum
refractive index, regardless of its physical density [41]. Here, we show that
RG not only works to capture the stationary optical response of a dense gas,
but also to capture the above mentioned time-dependent dephasing dynamics of
spin waves, in a simple and non-perturbative way. The validity of the RG
approach is quantitatively verified by comparison with full, microscopic
coupled-dipole simulations of large ($N\sim 10^{4}$) atomic ensembles.

Figure 1: (a) A pair-wise approach to the many-atom optical dynamics. As
derived in Ref. [30], in a dilute atomic medium, a small fraction of pairs of
atoms (red circles) are separated by a distance much smaller than a
wavelength, and thus interact strongly via their near fields. These atoms can
be replaced with a non-interacting, dynamically equivalent pair with new
effective frequencies (indicated by different colors), which evolve with the
phase $e^{-i\omega_{\pm}t}$. Single isolated atoms instead will not
significantly contribute to time evolution. (b) Representation of the RG
scheme [41]. Each step is characterized by identifying the most strongly
interacting pairs and replacing them with two new atoms with different
frequencies, which do not interact anymore through the near field. Unlike in
Fig. 1 (a), one can continue this process (including the renormalization of
atoms previously renormalized) until all near-field interactions have been
eliminated. The overall system at the end is equivalent to an inhomogeneously
broadened ensemble of spectral distribution $P(\omega)$, plotted in Fig. 3.

The remainder of the paper is structured as follows. In Sec. 2 we briefly
review the microscopic theoretical description of photon-mediated dipole-
dipole interactions, which accounts for atomic positions, near-field
interactions, and multiple scattering of light, and which serves as the basis
for the microscopic simulations of spin-wave dynamics. In Sec. 3 we describe
the RG approach, which enables one to predict a universal inhomogeneous
broadening function for a disordered medium. From here, we then formulate a
simple, approximate, single-atom model for the spin-wave dephasing dynamics.
In Sec. 4, we present detailed numerical simulations of the spin wave dynamics
from dilute to high-density media, which show both the initial exponential
dephasing and non-exponential behavior at later times. We also compare these
results with the RG approach, which exhibits good quantitative agreement in
all regimes. We conclude and provide an outlook in Sec. 5.

##  2 Microscopic model of atom-light interaction dynamics

We consider a minimal model consisting of $N$ identical 2-level atoms at
fixed, random positions $\left\\{\mathbf{r}_{i}\right\\}_{i=1,\dots,N}$ that
are uniformly distributed within a spherical cloud. The ground and excited
states $|g_{j}\rangle$ and $|e_{j}\rangle$ have an electric dipole transition
characterized by resonance frequency $\omega_{0}=ck_{0}$ and wavelength
$\lambda_{0}=2\pi/k_{0}$, and a single-atom excited-state spontaneous emission
rate given by $\Gamma_{0}$. We also define a dimensionless density in terms of
the number of atoms per cubic wavelength $\eta=\lambda_{0}^{3}N/V$, where
$V=\frac{4}{3}\pi R^{3}$ is the volume of the ensemble and $R$ its radius.

Within the standard assumptions in quantum optics (dipole, rotating wave and
Markov approximations), the effects of photon-mediated dipole-dipole
interactions, multiple scattering, and wave interference in spontaneous
emission are captured by an effective atomic Hamiltonian [42, 31, 43]:

$H=-\frac{3\pi\Gamma_{0}}{k_{0}}\sum_{j\ell}\mathbf{d}^{*}\cdot\mathbf{G}_{0}(\mathbf{r}_{j}-\mathbf{r}_{\ell},\omega_{0})\cdot\mathbf{d}\
\sigma_{eg}^{j}\sigma_{ge}^{\ell}.$ (1)

Physically, the photon-mediated interactions between atoms are described by
the free-space Green’s tensor $\mathbf{G}_{0}(\mathbf{r},\omega_{0})$, the
fundamental solution of the wave equation [44]. The Green’s tensor
characterizes how a photon emitted by an atom at $\mathbf{r}_{j}$, via the
action of the atomic lowering operator $\sigma_{ge}^{j}=|g_{j}\rangle\langle
e_{j}|$, propagates to a second atom at $\mathbf{r}_{i}$, travelling a
distance $r_{j\ell}=|\textbf{r}_{j\ell}|=|\mathbf{r}_{j}-\mathbf{r}_{\ell}|$.
Moreover, as $\mathbf{d}$ here represents the orientation of the atomic dipole
matrix element, which we fix to be along $\hat{\mathbf{x}}$, we conveniently
define $\theta_{j\ell}$ as its angle with respect to $\mathbf{r}_{j\ell}$.
Doing so, the Hamiltonian
$H=\sum_{j\ell}H_{j\ell}\sigma_{eg}^{j}\sigma_{ge}^{\ell}$ is characterized by
the coefficients

$H_{j\ell}=-\frac{3\Gamma_{0}}{4}e^{ik_{0}r_{j\ell}}\left[\frac{3\cos^{2}\theta_{j\ell}-1}{(k_{0}r_{j\ell})^{3}}-i\frac{3\cos^{2}\theta_{j\ell}-1}{(k_{0}r_{j\ell})^{2}}-\frac{\cos^{2}\theta_{j\ell}-1}{k_{0}r_{j\ell}}\right].$
(2)

The Hamiltonian is symmetric and non-Hermitian as a consequence of
$\mathbf{G}_{0}$ enforcing reciprocity and being a complex quantity, including
both coherent and dissipative interactions.

Important to later discussions, the real part of the Green’s function
(describing coherent interactions) contains a $\sim 1/r^{3}$ near-field
component, which dominates at small inter-atomic distances ($k_{0}r\\!<\\!1$).
This near-field interaction explicitly reads:

$H_{j\ell}^{\text{near}}=-\frac{3\Gamma_{0}}{4(k_{0}r_{j\ell})^{3}}\left(3\cos^{2}\theta_{j\ell}-1\right).$
(3)

The dissipative part instead describes collective spontaneous emission as
arising from wave interference of the emitted light, while in the limit of a
single atom, predicts the known spontaneous emission rate of
$\Gamma_{0}=\omega_{0}^{3}d_{eg}^{2}/3\pi\hbar\varepsilon_{0}c^{3}$, with
$d_{ge}$ being the amplitude of the dipole matrix element for the atomic
transition. Generally, the presence of dissipation requires a master equation
treatment [31, 45, 46, 47], but the non-Hermitian Hamiltonian (1) is
sufficient to describe the single-excitation regime of interest in our work,
which generally reduces to solving a set of classical coupled dipole equations
of motion [19, 44, 16, 48, 43, 49, 48, 50, 51, 52].

We will specifically be interested in applying the Hamiltonian (2) above to
investigate the dynamics of a single-excitation “timed Dicke state” or spin
wave, defined as

$|\mathbf{k}\rangle=\frac{1}{\sqrt{N}}\sum_{j}e^{i\mathbf{k}\cdot\mathbf{r}_{j}}|e_{j}\rangle.$
(4)

These collective states with well-defined wavevector constitute a natural
basis to describe light-matter excitations. For example, phase-matched spin
waves with $|\mathbf{k}|\approx k_{0}=\omega_{0}/c$ are naturally and easily
excited by an incoming resonant short pulse. By reciprocity, it is well known
that they also efficiently emit into a narrow, well-defined direction centered
around $\mathbf{k}$ [53], with a collectively enhanced rate [54, 55, 56, 16,
50, 48, 57, 52, 58], $\Gamma_{|\mathbf{k}|\sim
k_{0}}^{\text{coll}}/\Gamma_{0}=1+\overline{\text{OD}}/4$, which linearly
scales with the average optical depth of the medium. This narrow emission
occurs due to constructive interference of the emitting atoms along the
$\mathbf{k}$ direction, and forms the basis of collective enhancement at the
heart of efficient atom-light interfaces [12, 59] and the applications
mentioned in the introduction. This behavior can be equally derived by
microscopic theories [58] or by the macroscopic MBE [54, 55, 16].

It should be noted that phase-matched spin wave excitations undergo non-
trivial macroscopic spatio-temporal propagation dynamics [60, 61]. This makes
it challenging to quantify the magnitude and effects of microscopic dephasing,
due to the difficulty in defining an ideal time-evolving reference state to
compare to, if dephasing could hypothetically be eliminated. An elegant,
robust solution to this problem was proposed in Ref. [30], with experimental
realization based on a series of time-domain spin-wave control techniques [58,
62]. In particular, one can excite mismatched spin-wave excitations
$|\mathbf{k}\rangle$, characterized by $|\mathbf{k}|\neq k_{0}$, which neither
couple to light efficiently nor have any preferred emission or propagation
direction. Due to the lack of any direction in which the emitted field
(averaged over random configurations) interferes, the average initial
spontaneous emission rate of this state reduces to the single-atom value
$\Gamma_{|\mathbf{k}|\neq k_{0}}^{\text{coll}}=\Gamma_{0}$ [25]. Since this
state does not exhibit any background macroscopic dynamics,
$|\mathbf{k}\rangle$ itself serves as a reference to compare against the
actual time-evolved state.

Figure 2: (a) Snaphots of the time evolution of a mismatched spin-wave. We
initially prepare a mismatched spin-wave $|\mathbf{k}\rangle$
($|\mathbf{k}|=6.0k_{0}$) in a particular configuration $\\{\mathbf{r}_{j}\\}$
of a disordered gas at density $\eta=10$, $N=10^{4}$ atoms and radius
$R/\lambda_{0}\sim 6$. The state is then let to evolve under the dipole-dipole
Hamiltonian (1) and the projection over the initial state is computed,
$\langle\mathbf{k}|\mathbf{k}(t)\rangle=\sum_{j}|f_{j}|e^{i\varphi_{j}}$,
extracting the time dependent amplitudes and phases in the single atom basis.
To create the snapshots above we consider all the atoms contained in a slice
of size $\Delta y=R/35$, plot their position along the $xz$ plane, and color
them according to their accumulated phase in time evolution $\varphi_{j}$, as
a consequence of the interaction, for different times. The global effect of
this dephasing is represented in (b) where we plot the global overlap
$O_{\mathbf{k}}=|\langle\mathbf{k}|\mathbf{k}(t)\rangle|^{2}$ (blue curve) and
population $P_{\mathbf{k}}=\langle\mathbf{k}(t)|\mathbf{k}(t)\rangle$ (purple)
of the same mismatched spin-wave. Quantifying the deviation of the overlap
with respect to the predicted decay $e^{-\Gamma_{0}t}$ (red dotted line) is
the main purpose of this work.

Specifically, we will be interested in the time evolution under the dipole-
dipole Hamiltonian (1) as given by
$|\mathbf{k}(t)\rangle=e^{-iHt}|\mathbf{k}\rangle$. Since the initial state
contains only a single excitation, as argued above, the dynamics can be
efficiently evaluated numerically, with the resulting equations of motion
equivalent to classical coupled dipole equations. From the time-dependent
state, we can construct two quantities of interest,

$\displaystyle P_{\mathbf{k}}(t)$
$\displaystyle=\langle\mathbf{k}(t)|\mathbf{k}(t)\rangle$ (5) $\displaystyle
O_{\mathbf{k}}(t)$
$\displaystyle=|\langle\mathbf{k}|\mathbf{k}(t)\rangle|^{2}.$

The first quantity $P_{\mathbf{k}}(t)\leq 1$ monotonically decreases and gives
the total remaining excited state population at any time $t$, with the rest
having been irreversibly lost due to (collective) spontaneous emission. The
second quantity, $O_{\mathbf{k}}(t)$, on the other hand, quantifies the
overlap with the initial spin wave, and thus describes the survival of the
spin-wave order. Importantly, it was shown to be a measurable quantity in the
experiments of Ref. [30].

While we will present a more systematic analysis in Sec. 4, we provide a
visual example of the physics encoded in $P_{\mathbf{k}}(t)$ and
$O_{\mathbf{k}}(t)$ in Fig. 2. In particular, we simulate the dynamics of an
initial spin wave $|\mathbf{k}\rangle$, with $|\mathbf{k}|=6k_{0}$, in a
particular configuration $\\{\mathbf{r}_{j}\\}$ of a disordered gas of
$N=10^{4}$ atoms with density $\eta=10$, and in a spherical volume of radius
$R/\lambda_{0}\approx 6$. The evolved state $|\mathbf{k}(t)\rangle$ under Eq.
(1) is calculated for several specific times $t$. In Fig. 2 (a) we consider
all the atoms contained within a slice $\Delta y=R/35$ of the center of the
cloud, and plot their positions in the $x$-$z$ plane. The colors represent the
accumulated phase of each atom relative to its initial value
$e^{i\mathbf{k}\cdot\mathbf{r}_{j}}$, with a strong dephasing evident at time
$\Gamma_{0}t=1.0$. In Fig. 2 (b), we then plot $P_{\mathbf{k}}(t)$ (purple)
and $O_{\mathbf{k}}(t)$ (blue) as calculated for the entire ensemble, along
with the single-atom spontaneous emission $e^{-\Gamma_{0}t}$ for reference
(dashed red). It can be seen that while the initial decay of the total excited
population $P_{\mathbf{k}}(t)$ occurs at a rate $\sim\Gamma_{0}$ (confirming
the absence of collective enhancement) before slowing down, the spin wave
survival ratio $O_{\mathbf{k}}(t)$ decays significantly faster than
$\Gamma_{0}$, due to the dephasing illustrated in Fig. 2 (a).

##  3 A Renormalization group approach

While the microscopic model (Eq. (1)) can be numerically solved for moderate
atom number, its complexity scales directly with the number of atoms $N$, and
with the number of configurations needed to obtain disorder-averaged results.
Furthermore, the exact numerics does not directly elucidate the underlying
physics. Motivated by that, here we introduce a simple model, based on a
strong-disorder RG approach, which clearly identifies the role that near-field
interactions have on the dynamics, and which allows the effects of such
interactions to be captured by a simple effective, single-atom theory.

In a previous work [30], we derived a theory to understand the short-time
dephasing rate of a spin wave for dilute ensembles with densities
$\eta\lesssim 1$. We then found that the dephasing is primarily attributable
to a small fraction of atomic pairs separated by a distance much smaller than
the optical wavelength, $k_{0}r\\!<\\!1$ (highlighted in Fig. 1 (a)), which
strongly interact via their near fields (Eq. (3)). The $\sim 1/r^{3}$ scaling
of the near field implies (in three dimensions) that the presence of other
atoms is just a weak perturbation on top of the strong pairwise interaction,
such that the pair can be separately and approximately diagonalized. In the
single-excitation manifold, diagonalizing the near-field interaction (3) of a
pair yields symmetric and anti-symmetric eigenstates,
$|\pm\rangle=(|eg\rangle\pm|ge\rangle)/\sqrt{2}$, which experience opposite
frequency shifts
$\omega_{\pm}=\pm\frac{3\Gamma_{0}}{4k_{0}^{3}r^{3}}(1-3\cos^{2}\theta)$
relative to the bare atomic transition frequency.

The time evolution for the two-body problem can now be studied in terms of its
normal modes. Concretely, in the single atom rotating frame
$e^{-i\omega_{0}t}$, an initially prepared two-body spin-wave
$|\mathbf{k}\rangle$, will evolve as
$\langle\mathbf{k}|\mathbf{k}(t)\rangle=e^{-i\omega_{+}t}|c^{+}_{\mathbf{k}}|^{2}+e^{-i\omega_{-}t}|c^{-}_{\mathbf{k}}|^{2}$,
having defined the projections
$c^{\pm}_{\mathbf{k}}=\langle\mathbf{k}|\pm\rangle$. Although the magnitude of
$\mathbf{k}$ might be constrained in an experiment, we can take the conceptual
limit where $k\rightarrow\infty$, or infinite mismatch. This implies that the
phase $e^{i\mathbf{k}\cdot\mathbf{r}_{j}}$ of each excited atom is effectively
random (being infinitely sensitive to the specific atomic position), which
implies that the actual spin wave should have on average equal overlap with
the $\pm$ eigenstates, i.e. that $|c^{\pm}_{\mathbf{k}}|^{2}\rightarrow 1/2$.
The dynamics can therefore be equivalently modeled by replacing the two
original atoms with two new atoms of new resonance frequencies $\omega_{\pm}$,
that now do not interact anymore through the near field, but evolve with their
“free” inhomogeneous phases $e^{-i\omega_{\pm}t}$ (see Fig. 1 (a)).

In Ref. [30], the short-time dynamics of a many-atom, dilute ensemble were
approximated by identifying a small fraction of close-by atomic pairs that
evolved with the phases $e^{-i\omega_{\pm}t}$ as argued above, while the
remaining atoms were assumed to undergo no evolution. By taking the position
dependence of the function $\omega_{\pm}(\mathbf{r})$ and combining with the
known distribution function of separations $\mathbf{r}$ of nearest neighbors
in a random ensemble, an initial decay of the spin wave survival order was
predicted and measured to be exponential, $|O_{\mathbf{k}}(\delta
t)|^{2}=e^{-\Gamma_{\mathbf{k}}\delta t}$, with a density-dependent rate
$\Gamma_{\mathbf{k}}/\Gamma_{0}=1+\xi\eta$, where $\xi=1/6\pi\sqrt{3}$ for
two-level atoms.

Now, in a dense ensemble (Fig. 1 (b)) or at longer evolution times, the basic
picture of replacing close-by atomic pairs with new effective atoms of
renormalized frequencies does not change. However, a key realization is that
after an atom has been renormalized, it can still see another atom close by
with which it can strongly interact (again highlighted by red circles in Fig.
1 (b)). This allows yet another renormalization step to take place, which now
will involve the diagonalization of a pair of atoms with the effective
frequencies previously obtained. Whereas only a small fraction of atoms
dictates the initial decay in a dilute ensemble, here, we must specify a
general procedure valid for any density, to repeatedly and hierarchically
identify the single most strongly interacting pair (including the possibility
that the pair contains already renormalized atoms) and replace them with two
new effective atoms.

Having anticipated that RG can involve the renormalization of atoms that have
already been renormalized in previous steps, here, we consider the more
general case of two inhomogeneous atoms of general resonance frequencies of
$\omega_{i}$ and $\omega_{j}$, interacting via the near field, as described by
the two-body Hamiltonian in the single excitation sector,

$H_{j\ell}^{2b}=\langle\omega_{j\ell}\rangle\mathbb{1}+\begin{pmatrix}\delta\omega_{j\ell}&H_{j\ell}^{\text{near}}\\\
H_{j\ell}^{\text{near}}&-\delta\omega_{j\ell}\end{pmatrix}.$ (6)

For convenience, we have defined
$\langle\omega_{j\ell}\rangle=(\omega_{j}+\omega_{\ell})/2$,
$\delta\omega_{j\ell}=(\omega_{j}-\omega_{\ell})/2$. To quantify the strength
of the interaction, we define the ratio between the off-diagonal and the
diagonal elements,
$\mathcal{K}_{j\ell}=\mathcal{L}_{j\ell}|H_{j\ell}^{\text{near}}|/(|\delta\omega_{j\ell}|+1)$,
where the matrix $\mathcal{L}$ keeps the information of whether a pair of
atoms has already been renormalized ($\mathcal{L}_{j\ell}=0$) or not
($\mathcal{L}_{j\ell}=1$). This prevents a pair of renormalized atoms from
being renormalized between themselves multiple times (although each atom from
the pair can be renormalized with other atoms). Intuitively, a large value of
$\mathcal{K}_{j\ell}$ (which requires $\mathcal{L}_{j\ell}=1$) means that the
strength of the interaction is able to further split the original frequency
difference $\delta\omega_{j\ell}$. Thus, the most strongly interacting pair is
identified as that with the largest value of $\mathcal{K}_{j\ell}$ (red
circles in Fig. 1). Once identified, the full diagonalization of (6) gives two
eigenvalues,
$\omega_{\pm}=\langle\omega_{j\ell}\rangle+\sqrt{\delta\omega_{j\ell}^{2}+\left(H_{j\ell}^{\text{near}}\right)^{2}}$.
The pair can be therefore replaced by an approximately equivalent one, made of
two atoms with the new resonance frequencies $\omega_{\pm}$ and that do not
interact anymore through the near field (setting $\mathcal{L}_{j\ell}=0$).

Repeated application of this algorithm, which constitutes the renormalization
group (RG) flow, continues until all near-field interactions have been removed
and the atoms have been assigned the new effective frequencies
$\left\\{\omega_{i}\right\\}$. When the RG scheme is applied to multiple
realizations of randomly distributed ensembles of atoms, we can build up the
probability distribution $P(\omega)$ of the effective frequencies, which we
illustrate in Fig. 3.

As the near-field interaction only depends on distance through the
dimensionless parameter $(k_{0}r)^{-3}$, the distribution when rescaled by
density, $P(\omega/\eta)$, should be a universal function for a sufficiently
large number of atoms and sufficiently large geometry, where boundary effects
are negligible. We can directly confirm this numerically in Fig. 3, where we
plot $P(\omega/\eta)$ obtained from RG for various densities.

Figure 3: Universal probability distribution of normalized frequencies
$\omega/\eta$. We apply the RG approach to spherical samples of $N=2500$ atoms
at dimensionless densities $\eta=500,1000,1500$ to extract the new effective
frequencies $\left\\{\omega_{i}\right\\}$, sampling approximately $\sim
10^{3}$ different configurations at each density to build the probability
distribution $P(\omega/\eta)$. The purple solid line instead corresponds to
the exact numerical spectrum of the near-field interaction Hamiltonian. This
is obtained considering an ensemble of atoms (same parameters as before) and
diagonalizing $H^{\text{near}}$ (defined in Eq. (3)) to get the eigenfrequency
statistics.

The high-frequency tails of $P(\omega)$ correspond to the most strongly
interacting pairs, which are renormalized at the beginning of the flow. This
perfectly matches the simpler theory presented in [30] for dilute atoms, based
on the probability distribution of nearest neighbors, where it was found that
asymptotically

$P(\omega)\sim_{\pm\infty}\xi\eta\frac{1}{2\pi\omega^{2}}$ (7)

The central part of the distribution $P(\omega)$ instead consists of atoms
that have been renormalized multiple times. In this sense, the presented RG
scheme and the resulting probability distribution is the correct way to
capture the near field induced inhomogeneous broadening of the medium and the
induced dephasing rate of spin-waves in dense media, as we are going to
discuss.

Formally, the RG procedure amounts to approximately diagonalizing the near-
field part of the Hamiltonian (1) (an $N\\!\times\\!N$ matrix), by repeatedly
identifying and diagonalizing a dominant interacting pair of atoms (a $2\times
2$ block). We can quantify the error by comparing the resulting frequency
distribution $P(\omega)$ obtained by RG, with the eigenvalue distribution
obtained by exact, numerical diagonalization of the real, symmetric $N\times
N$ Hamiltonian $H^{\text{near}}$ defined by only considering the near field
interaction (Eq.(3)). It can be seen in Fig. 3 that the two are essentially
indistinguishable, which validates the RG scheme.

We now discuss how the RG results can be applied to predict the time evolution
of a spin wave in an effective single-atom picture, focusing on the spin-wave
survival ratio $O_{\bf k}(t)$. In particular, while the coherence of a single,
isolated atom (without the rotating frame) is expected to evolve as
$\langle\sigma_{ge}(t)\rangle=\langle\sigma_{ge}(0)\rangle
e^{-i(\omega_{0}-i\Gamma_{0}/2)t}$, the distribution in resonance frequencies
of the new effective medium will result in an uncertainty of accumulated phase
in time

$O_{\textsc{rg}}(t)=\left|\int d\omega P(\omega)e^{i\omega
t}\right|^{2}e^{-\Gamma_{0}t},$ (8)

thus introducing microscopically-driven dephasing due to inhomogeneous
broadening. We note that the simplicity of Eq. (8) arises in part due to our
choice of a phase-mismatched spin wave as the initial state. In particular,
although RG removes the strong near-field interactions, the atoms can still
interact and radiate via the far field, which as we described earlier can in
general result in non-trivial spatio-temporal dynamics even within the MBE.
The use of a phase-mismatched spin wave suppresses the macroscopic dynamics,
leaving the predicted spin wave survival ratio to depend on the microscopic
quantity $P(\omega)$ alone.

##  4 Spin-wave dephasing

Figure 4:  Time evolution of an initially prepared ideal mismatched spin-wave
($|\mathbf{k}|\neq k_{0}$), at different values of the dimensionless density,
from small, $\eta=0.1$ to high $\eta=1000$. The blue line is the average time
evolution of the overlap $O_{\mathbf{k}}$ (see Eq.(5)) over different
realizations of the disordered gas, while the blue shaded region corresponds
to the standard deviation. The red dotted line shows an exponential decay with
a rate $\Gamma_{0}$, as predicted by treating the atomic medium as smooth
(MBE). The purple dashed line is the population of the time evolved state
divided by the atom number, $P_{\mathbf{k}}/N$. Finally, the orange line is
the overlap $O_{RG}$ as predicted by RG theory (see Sec. 3). We simulate the
time evolution of $N=10^{4}$ atoms, to guarantee that, at the maximum density
$\eta=1000$, the radius of the uniformly distributed spherical cloud is
$R/\lambda_{eg}=1.34>1$, such that the cloud is not subwavelength. All the
quantities are averaged over $N_{s}\sim 500$ different atomic samples. At
density $\eta=100$, the inset shows the short-time exponential dynamics. The
interval over which we fit for an exponential $e^{-(\Gamma_{0}+\gamma)t}$
(dashed black curve) is highlighted in green and corresponds to
$\Gamma_{\mathbf{k}}t<0.1$, while at longer times the overlap deviates from
this simple behavior, as predicted by our RG approach.

We now present the exact numerical simulations and analysis of the time
evolution of a mismatched spin wave for densities ranging from dilute
($\eta\ll 1$) to dense ($\eta\gg 1$), which we will then compare with the
simple RG prediction of Eq. (8). To be concrete, we take an initial state
consisting of a highly mismatched spin wave (Eq. (4)) with momentum
$|\mathbf{k}|\\!=\\!6k_{0}$, $\hat{x}-$polarized and directed along $\hat{z}$,
in an ensemble of $N\\!=\\!10^{4}$ atoms. Then, the time evolution of the
total excited-state population $P_{\mathbf{k}}(t)$ and of the overlap with the
initially prepared spin-wave order $O_{\mathbf{k}}(t)$ is calculated for
$N_{s}\\!=\\!500$ realizations of the disordered gas and averaged, as
indicated by $\langle P_{\mathbf{k}}\rangle$ for example. Numerical results
are represented in Fig. 4, where we also plot our prediction for time
evolution of the overlap, $O_{RG}$, made in Eq.(8) (orange dash-dotted lines),
based on the effective single-atom theory described in the previous section.

We first focus on the short-time dynamics. As introduced in the previous
section, the short-time decay of a spin-wave is predicted to be exponential
[30], $\langle O_{\mathbf{k}}(\delta t)\rangle\sim
e^{-\Gamma_{\mathbf{k}}\delta t}$, within a time interval
$\Gamma_{\mathbf{k}}\delta t\ll 1$. The rate
$\Gamma_{\mathbf{k}}=\Gamma_{0}+\gamma$ is given by the sum of the single-atom
emission rate, as predicted by the macroscopic MBE, and an additional density-
dependent dephasing rate, $\gamma=\Gamma_{0}\xi\eta$, with
$\xi=1/6\pi\sqrt{3}$. As shown in the inset of $\eta=100$ (Fig. 4), for
example, the microscopic dynamics reveal a short-time decay of the spin wave
order that is distinctly faster than the MBE prediction. The decay becomes
even more evident at higher densities ($\eta=10^{2},10^{3}$). We can confirm
the density dependence in the short-time decay rate $\Gamma_{\mathbf{k}}$ by
fitting the curve $\langle O_{\mathbf{k}}\rangle$ for each density in the
short-time window, defined by $\Gamma_{\mathbf{k}}\delta t<0.1$, to an
exponential of the form $e^{-(\Gamma_{0}+\gamma)t}$, and plotting the
dependence of the fit parameter $\gamma$ versus density parameter $\eta$ in
Fig. 5 (a). We also plot the prediction $\gamma=\Gamma_{0}\xi\eta$ in red. It
is important to observe that within our defined “short time” interval, already
$\sim 10\%$ of the initial spin-wave order is lost. An excellent agreement is
observed over a large range of densities, both changing size and number of
atoms (respectively blue squares and circles) of a random gas, confirming that
this effect does only depend on the dimensionless density parameter $\eta$.
This dephasing for moderate $\eta\lesssim 1$ was also observed experimentally
in Ref. [30].

Figure 5: (a) Short-time dephasing rate vs. dimensionless density $\eta$. Blue
dots are extracted from an exponential fit at short times ($\langle
O_{\mathbf{k}}\rangle\sim e^{-\Gamma_{\mathbf{k}}\delta t}$), within the time
interval $\Gamma_{\mathbf{k}}\delta t=0.1$, of the full time evolution of a
mismatched spin-wave ($|\mathbf{k}|=6k_{0}$) at different densities
$0.1\leq\eta\leq 10^{3}$. The red line is instead the theoretical prediction
for the short-time dephasing rate $\gamma=\Gamma_{e}\xi\eta$, discussed in the
main text. Simulations are performed with a cloud of atoms of fixed size
($R/\lambda_{eg}=1.34$) and averaging over $N_{s}\\!\sim\\!500$ realizations.
(inset) The same short time dephasing rate is evaluated (square dots), but now
to explore the low-density regime we simulate the spin-wave dynamics in an
ensemble with a fixed number of atoms, $N=10^{4}$ (again averaged $N_{s}\sim
500$ times) and varying the radius of the cloud. (b) Numerical simulations of
the time evolution of the averaged overlap $\langle O_{\mathbf{k}}\rangle$
(solid blue lines) for a mismatched spin-wave ($|\mathbf{k}|=6.0k_{0}$)
compared with the RG prediction $O_{RG}$ (orange dash-dotted line), and the
total excited population divided by atom number, $\langle P_{k}\rangle/N$
(purple dashed line). The simulations consider a system of fixed density
$\eta=100$, but different atom number $N=10^{2},10^{3},10^{4}$.

At longer times and at higher densities, the decay of spin wave order
noticeably deviates from exponential. Despite its simplicity, our single
effective atom model based on RG (Eq. (8)), displays excellent agreement
beyond the short-time interval. Viewed from the RG perspective, this non-
exponential contribution comes from the frequency components near the center
of the inhomogeneous broadening probability distribution $P(\omega/\eta)$
(Fig. 3), corresponding to atoms that are renormalized multiple times.
Interestingly, at even longer times, it can be seen that for each density, the
average $\langle O_{\mathbf{k}}\rangle$ deviates from our prediction and
saturates to a value that barely decreases over the range of times plotted. We
furthermore observe numerically that this value closely coincides with the
total excited state population remaining divided by the atom number, $P_{\bf
k}(t)/N$ (dashed purple curves).

The slow decay of population $P_{\mathbf{k}}(t)$ at long times is an effect
that has been studied extensively in recent years, and is known as late-time
subradiance [26, 63, 64, 65, 66]. While a microscopic derivation is difficult,
a heuristic argument can be made that the remaining population should be
roughly equally distributed throughout the ensemble, given a smooth initial
distribution. Furthermore, given the randomness of the dynamics, this
population will be statistically evenly distributed among any N extended modes
that can be defined for the system, such as our spin wave mode of interest.

As far as we can numerically check (e.g., up to $N=10^{4}$ atoms for a density
of $\eta=100$ in Fig. 5 (b)), we see that the RG prediction $O_{RG}$ follows
the actual spin wave survival ratio $O_{\mathbf{k}}$ for increasingly long
times as $N$ is increased, due to the decrease in the saturation value
$P_{\mathbf{k}}/N$. This strongly suggests that the single-atom RG prediction
should be interpreted as the correct description of the dephasing dynamics in
the thermodynamic limit, when the late-time population in any one mode $\sim
1/N$ becomes negligible.

##  5 Conclusions

In summary, we have developed an effective single-atom theory that describes
well the non-exponential dephasing dynamics of optical spin waves in
disordered atomic media, including at high densities and at long times. This
theory is based upon the technique of strong disorder renormalization group,
which treats the potentially strong near field interactions in such a medium
in a non-perturbative way.

We envision that our predictions, particularly in the high density regime,
could be immediately explored using solid-state emitter ensembles such as rare
earth doped crystals [67, 68, 69, 70], where many atoms per cubic wavelength
are typical. Separately, the remarkable accuracy by which RG is found to
reproduce the dephasing dynamics suggests that it can be a powerful tool to
quantitatively investigate and understand other microscopic optical phenomena
in disordered systems.

##  6 Acknowledgements

We acknowledge F. Andreoli for stimulating discussions. DEC acknowledges
support from the European Union’s Horizon 2020 research and innovation
programme, under European Research Council grant agreement No 101002107
(NEWSPIN), FET-Open grant agreement No 899275 (DAALI), and Quantum Flagship
project 820445 (QIA); the Government of Spain (Europa Excelencia program
EUR2020-112155, Severo Ochoa program CEX2019-000910-S, and MICINN Plan
Nacional Grant PGC2018-096844-B-I00); Generalitat de Catalunya through the
CERCA program, Fundació Privada Cellex, AGAUR Project No. 2017-SGR-1334,
Fundació Mir-Puig, and Secretaria d’Universitats i Recerca del Departament
d’Empresa i Coneixement de la Generalitat de Catalunya, co-funded by the
European Union Regional Development Fund within the ERDF Operational Program
of Catalunya (project QuantumCat, ref. 001-P-001644). YH and SW acknowledge
support from National Key Research Program of China under Grant No.
2017YFA0304204, and from NSFC under Grant No. 12074083.

## References

## References

  * [1] P. Vernaz-Gris, K. Huang, M. Cao, A. S. Sheremet, and J. Laurat,  Nat. Commun. 9 (2018).
  * [2] Y. Wang et al.,  Nat. Photonics 13, 346 (2019).
  * [3] H. Gorniaczyk, C. Tresp, J. Schmidt, H. Fedder, and S. Hofferberth,  Phys. Rev. Lett. 113, 053601 (2014).
  * [4] D. Tiarks, S. Baur, K. Schneider, S. Dürr, and G. Rempe,  Phys. Rev. Lett. 113, 053602 (2014).
  * [5] D. Tiarks, S. Schmidt-Eberle, T. Stolz, G. Rempe, and S. Dürr,  Nat. Phys. 15, 124 (2019).
  * [6] W. Wasilewski et al.,  Phys. Rev. Lett. 104, 133601 (2010).
  * [7] R. J. Sewell et al.,  Phys. Rev. Lett. 109, 253605 (2012).
  * [8] O. Hosten, N. J. Engelsen, R. Krishnakumar, and M. A. Kasevich,  Nature 529, 505 (2016).
  * [9] Z. Chen, J. G. Bohnet, J. M. Weiner, K. C. Cox, and J. K. Thompson,  Phys. Rev. A - At. Mol. Opt. Phys. 89, 043837 (2014).
  * [10] K. C. Cox, G. P. Greve, J. M. Weiner, and J. K. Thompson,  Phys. Rev. Lett. 116, 093602 (2016).
  * [11] R. J. Lewis-Swan, M. A. Norcia, J. R. Cline, J. K. Thompson, and A. M. Rey,  Phys. Rev. Lett. 121, 070403 (2018).
  * [12] K. Hammerer, A. S. Sørensen, and E. S. Polzik,  Rev. Mod. Phys. 82, 1041 (2010).
  * [13] C. M. Bowden and J. P. Dowling,  Phys. Rev. A 47, 1247 (1993).
  * [14] Y. Castin and K. Molmer,  Phys. Rev. A 51, R3426 (1995).
  * [15] M. Fleischhauer and S. F. Yelin,  Phys. Rev. A - At. Mol. Opt. Phys. 59, 2427 (1999).
  * [16] A. A. Svidzinsky, X. Zhang, and M. O. Scully,  Phys. Rev. A - At. Mol. Opt. Phys. 92, 013801 (2015).
  * [17] A. V. Gorshkov, A. André, M. Fleischhauer, A. S. Sørensen, and M. D. Lukin,  Phys. Rev. Lett. 98, 123601 (2007).
  * [18] A. V. Gorshkov, J. Otterbach, M. Fleischhauer, T. Pohl, and M. D. Lukin,  Phys. Rev. Lett. 107, 133602 (2011).
  * [19] O. Morice, Y. Castin, and J. Dalibard,  Phys. Rev. A 51, 3896 (1995).
  * [20] J. Ruostekoski and J. Javanainen,  Phys. Rev. Lett. 82, 4741 (1999).
  * [21] S. Jennewein et al.,  Phys. Rev. Lett. 116, 233601 (2016).
  * [22] S. Jennewein et al.,  Phys. Rev. A 97, 053816 (2018).
  * [23] S. D. Jenkins et al.,  Phys. Rev. A 94, 023842 (2016).
  * [24] N. J. Schilder, C. Sauvan, Y. R. Sortais, A. Browaeys, and J. J. Greffet,  Phys. Rev. Lett. 124, 073403 (2020).
  * [25] M. O. Scully,  Phys. Rev. Lett. 115, 243602 (2015).
  * [26] W. Guerin, M. O. Araújo, and R. Kaiser,  Phys. Rev. Lett. 116, 083601 (2016).
  * [27] G. Labeyrie et al.,  Phys. Rev. Lett. 83, 5266 (1999).
  * [28] Y. Bidel et al.,  Phys. Rev. Lett. 88, 2039021 (2002).
  * [29] C. M. Aegerter and G. Maret,  Prog. Opt. 52, 1 (2009).
  * [30] Y. He et al.,  ArXiv (2021).
  * [31] M. Gross and S. Haroche,  Phys. Rep. 93, 301 (1982).
  * [32] A. Cipris, R. Bachelard, R. Kaiser, and W. Guerin,  Phys. Rev. A 103 (2021).
  * [33] L. S. Levitov,  Phys. Rev. Lett. 64, 547 (1990).
  * [34] D. S. Fisher,  Phys. Rev. B 50, 3799 (1994).
  * [35] K. Damle, O. Motrunich, and D. A. Huse,  Phys. Rev. Lett. 84, 3434 (2000).
  * [36] O. Motrunich, S. C. Mau, D. A. Huse, and D. S. Fisher,  Phys. Rev. B - Condens. Matter Mater. Phys. 61, 1160 (2000).
  * [37] G. Refael and J. E. Moore,  Phys. Rev. Lett. 93, 260602 (2004).
  * [38] F. Iglói and C. Monthus,  Phys. Rep. 412, 277 (2005).
  * [39] R. Vosk and E. Altman,  Phys. Rev. Lett. 110, 067204 (2013).
  * [40] G. Refael and E. Altman,  Comptes Rendus Phys. 14, 725 (2013).
  * [41] F. Andreoli, M. J. Gullans, A. A. High, A. Browaeys, and D. E. Chang,  Phys. Rev. X 11, 11026 (2021).
  * [42] G. S. Agarwal,  Phys. Rev. A 2, 2038 (1970).
  * [43] A. Asenjo-Garcia, J. D. Hood, D. E. Chang, and H. J. Kimble,  Phys. Rev. A 95, 33818 (2017).
  * [44] L. Novotny and B. Hecht,  Principles of nano-optics,  Cambridge University Press, 2006.
  * [45] H. T. Dung, L. Knöll, and D. G. Welsch,  Phys. Rev. A - At. Mol. Opt. Phys. 66, 16 (2002).
  * [46] S. Y. Buhmann, L. Knöll, D. G. Welsch, and H. T. Dung,  Phys. Rev. A - At. Mol. Opt. Phys. 70, 52117 (2004).
  * [47] S. Y. Buhmann,  Dispersion Forces I,  Springer-Verlag, 2012.
  * [48] B. Zhu, J. Cooper, J. Ye, and A. M. Rey,  Phys. Rev. A 94, 023612 (2016).
  * [49] S. E. Skipetrov and I. M. Sokolov,  Phys. Rev. Lett. 112, 023905 (2014).
  * [50] M. O. Araújo, I. Krešić, R. Kaiser, and W. Guerin,  Phys. Rev. Lett. 117, 073002 (2016).
  * [51] R. J. Bettles, S. A. Gardiner, and C. S. Adams,  Phys. Rev. Lett. 116, 103602 (2016).
  * [52] S. L. Bromley et al.,  Nat. Commun. 7, 1 (2016).
  * [53] M. O. Scully, E. S. Fry, C. H. Ooi, and K. Wódkiewicz,  Phys. Rev. Lett. 96, 010501 (2006).
  * [54] M. O. Scully and A. A. Svidzinsky,  Science (80-. ). 325, 1510 (2009).
  * [55] R. Friedberg and J. T. Manassah,  Phys. Lett. Sect. A Gen. At. Solid State Phys. 374, 1648 (2010).
  * [56] T. Bienaimé, M. Petruzzo, D. Bigerni, N. Piovella, and R. Kaiser,  J. Mod. Opt. 58, 1942 (2011).
  * [57] W. Guerin, M. T. Rouabah, and R. Kaiser,  J. Mod. Opt. 64, 895 (2017).
  * [58] Y. He et al.,  Phys. Rev. Res. 2, 043418 (2020).
  * [59] M. Fleischhauer, A. Imamoglu, and P. J. Marangos,  Rev. Mod. Phys. 77, 633 (2005).
  * [60] F. Cottier, R. Kaiser, and R. Bachelard,  Phys. Rev. A 98, 013622 (2018).
  * [61] T. S. Do Espirito Santo et al.,  Phys. Rev. A 101, 013617 (2020).
  * [62] Y. He et al.,  Phys. Rev. Lett. 125, 213602 (2020).
  * [63] M. O. Araújo, W. Guerin, and R. Kaiser,  J. Mod. Opt. 65, 1345 (2018).
  * [64] P. Weiss, M. O. Araújo, R. Kaiser, and W. Guerin,  New J. Phys. 20, 063024 (2018).
  * [65] Y. A. Fofanov, I. M. Sokolov, R. Kaiser, and W. Guerin,  ArXiv (2020).
  * [66] G. Ferioli, A. Glicenstein, L. Henriet, I. Ferrier-Barbut, and A. Browaeys,  Phys. Rev. X 11, 021031 (2021).
  * [67] H. De Riedmatten, M. Afzelius, M. U. Staudt, C. Simon, and N. Gisin,  Nature 456, 773 (2008).
  * [68] M. P. Hedges, J. J. Longdell, Y. Li, and M. J. Sellars,  Nature 465, 1052 (2010).
  * [69] T. Zhong, J. M. Kindem, E. Miyazono, and A. Faraon,  Nat. Commun. 6, 8206 (2015).
  * [70] M. Zhong et al.,  Nature 517, 177 (2015).