# Phase-space views into dye-microcavity thermalised and condensed photons

Jakov Marelic  Benjamin T. Walker  Robert A. Nyman<EMAIL_ADDRESS>Quantum Optics and Laser Science group, Blackett Laboratory, Imperial College
London, Prince Consort Road, SW7 2BW, United Kingdom

###### Abstract

We have observed momentum- and position-resolved spectra and images of the
photoluminescence from thermalised and condensed dye-microcavity photons. The
spectra yield the dispersion relation and the potential energy landscape for
the photons. From this dispersion relation, we find that the effective mass is
that of a free photon not a polariton. We place an upper bound on the
dimensionless two-dimensional interaction strength of $\tilde{g}\lesssim
10^{-3}$, which is compatible with existing estimates. Both photon-photon and
photon-molecule interactions are weak. The temperature is found to be
independent of momentum, but dependent on pump spot size, indicating that the
system is ergodic but not perfectly at thermal equilibrium. Condensation
always happens first in the mode with lowest potential and lowest kinetic
energy, although at very high pump powers multimode condensation occurs into
other modes.

The concept of a phase-space distribution played an important role in the
development of statistical mechanics before quantum mechanics divided the
universe into discrete states Gibbs 111We would have cited Ludwig Boltzmann
here but our German language skills are insufficient to find which of his
articles is most appropriate.. The Boltzmann equation Huang  describes the
evolution of the phase-space distribution for a system made of particles which
interact among themselves, through a term known as the collision integral. For
indistinguishable quantum particles, the collision integral can be modified
and the Boltzmann equation is still validWu97 . There exist time-invariant
phase-space distributions which are functions only of a phase-space local
energy (called “ergodic”) for which entropy is maximised for some fixed
control parameters (called “thermal equilibrium”). At thermal equilibrium, the
statistical distribution can be shown on purely information-theoretic grounds
to be independent of the physical phenomena by which equilibrium is obtained
Jaynes57 .

True thermal equilibrium takes an infinite time to be achieved. In dynamical
thermalising systems with finite lifetime, thermal equilibrium cannot
therefore occur except approximately. Two important examples of dynamical-
equilibrium quantum systems are trapped atomic gases and microcavity
polaritons. Ultracold atoms can be modelled using the Boltzmann equation
Luiten96 , with the collision integral directly evaluated, for example using
Monte Carlo molecular-trajectory methods Wu96 . Loss is considered as a
perturbation, either a low-probability event or a high-energy truncation of a
near-equilibrium distribution. Microcavity polaritons, even when non-
resonantly driven, can be driven far from equilibrium, and may show very non-
thermal steady states Kasprzak08 .

Photons in a dye-filled microcavity Klaers10a ; Klaers10b  are unusual in that
the thermal equilibrium process is not among particles, but between particles
and a thermal bath deLeeuw13 ; Chiocchetta14  via the vibrational relaxation
of the dye molecules. Interparticle interactions are very weak so the
collision integral is effectively zero Nyman14 ; VanDerWurff14 , yet still
bosonic exchange statistics lead to condensation Klaers10b ; Snoke13 . A
classical rate equation for the populations of states can be derived from a
quantum master equation Kirton13 ; Kirton15 ; Keeling16 , but the Boltzmann
equation itself does not directly apply, because the collision integral would
be zero. Non-equilibrium distributions have been observed with time-dependent
Schmitt15  or steady-state Marelic15 ; Marelic16  pumping. Approximate thermal
equilibrium is obtained if the scattering rate for photons from dye molecules
is faster than the cavity loss rate, and under those circumstances standard
thermodynamic properties can be observed, such as the lambda transition
associated with Bose-Einstein condensation Damm16 .

One can think of the photoluminescence emitted from the microcavity as
representing the phase-space distribution of the intra-cavity photons. A
position-space image is the phase-space distribution marginalised over all
momenta. Similarly the observed spectrum can be thought of as the phase-space
distribution integrated over the phase space wherever the local-energy
function is matched. The spectrum is not always sufficient to determine
whether states are occupied thermally, since the density of states may not be
known. First- and second-order correlations are simply related in thermal
equilibrium, and dye-microcavity photons have been seen to be compatible with
thermal equilibrium Schmitt16 .

We begin this manuscript with a discussion of phase-space distributions for
non-interacting thermal-equilibrium distributions, and how they can be
observed. We then describe our experimental setup for dye-thermalised
microcavity photons, including the optics required for various phase-space
views. We present momentum-resolved spectra both below and above threshold
pump power for condensation. From the observed dispersion relation we infer
upper bounds for the strength of photon-molecule and photon-photon
interactions. We also show momentum-space images, noting that condensation
happens in the lowest-momentum, lowest-energy mode. Turning to position-
resolved spectra, we map the potential energy landscape.

##  I Views into phase-space at thermal equilibrium

The phase-space distribution, the population as a function of momentum for all
positions, of a non-interacting gas contains all of the relevant information
about the gas, when combined with the phase-space local energy. Non-
equilibrium, or non-ergodic systems can be fully described via the phase-space
distribution. We shall write a general phase-space distribution as $f({\bf
r,p})$, as function of position ${\bf r}$ and momentum ${\bf p}$. The
microcavity used for thermalising photons permits two free dimensions, but the
following discussion can easily generalised to other dimensionalities.

One simple view into phase space is a position-space image $n({\bf r})$ which
is simply $f({\bf r,p})$ marginalised over all momenta:

$\displaystyle n({\bf r})=\int{\rm d}^{2}{\bf p}\,f({\bf r,p})$ (1)

Similarly the observed spectrum $n(E)$ can be thought of as the phase-space
distribution integrated over the phase-space wherever the local-energy
function $\epsilon({\bf r,p})$ is matched 222Here we are making an
approximation that the phase-space distribution varies slowly compared to
discreteness of states, i.e. Planck’s constant is small.:

$\displaystyle n(E)=\iint{\rm d}^{2}{\bf r}\,{\rm d}^{2}{\bf p}\,f({\bf
r,p})\delta\left[{E-\epsilon({\bf r,p})}\right]$ (2)

The position-resolved spectrum $n({\bf r},E)$ is an alternative view into
phase space:

$\displaystyle n({\bf r},E)=\int{\rm d}^{2}{\bf p}\,f({\bf
r,p})\delta\left[{E-\epsilon({\bf r,p})}\right]$ (3)

The momentum-resolved spectrum, $n({\bf p},E)$ is determined similarly,
marginalising over positions.

Thermal equilibrium distributions are determined by the phase-space local
energy only. For light in a near-planar microcavity the local potential energy
$V({\bf r})$ is determined by the shape of the mirrors that form the cavity
Lugiato87 ; Chiao99 ; Chiao00 . The kinetic energy $Q({\bf p})$ is momentum-
space local. For non-interacting particles of mass $m$, $Q({\bf
p})=\frac{p^{2}}{2m}$. The total energy is then $\epsilon({\bf r,p})=Q({\bf
p})+V({\bf r})$.

For example the Boltzmann distribution of occupancies of states of energy $E$
is $f({\bf r,p})=\frac{1}{Z}{\rm e}^{\left[{-\frac{{\epsilon({\bf
r,p})}}{k_{B}T}}\right]}$, where $\iint{\rm d}^{2}{\bf r}\,{\rm d}^{2}{\bf
p}\,f({\bf r,p})=1$ defines the partition function $Z$, and $T$ is the
temperature. It is straightforward to show that, since the energy is a
quadratic, radially symmetric function of the magnitude of momentum in two
dimensions for free massive particles:

$\displaystyle n({\bf r},E)=\begin{cases}\frac{2\pi m}{Z}{\rm
e}^{\left[{-\frac{{E}}{k_{B}T}}\right]}&\text{, if }E>V({\bf r})\\\ 0&\text{
otherwise}\end{cases}$ (4)

The position-resolved spectrum then gives information about both the thermal
distribution, through the temperature, and the potential energy landscape,
through a position-dependent cutoff in the spectrum, below which no light is
detectable.

The momentum-resolved spectrum is defined similarly: $n({\bf p},E)=\int{\rm
d}^{2}{\bf r}\,f({\bf r,p})\delta\left[{E-\epsilon({\bf r,p})}\right].$ In the
case of a quadratic radially symmetric potential (the isotropic harmonic
oscillator, typically the case for dye-microcavity thermalised photons) with
$V({\bf r})=\frac{1}{2}m\Omega^{2}r^{2}$ for trap angular frequency $\Omega$,
the momentum-resolved spectrum looks very similar to the position-resolved
spectrum:

$\displaystyle n({\bf p},E)=\begin{cases}\frac{2\pi}{Zm\Omega^{2}}{\rm
e}^{\left[{-\frac{{E}}{k_{B}T}}\right]}&\text{, if }E>\frac{p^{2}}{2m}\\\
0&\text{ otherwise}\end{cases}$ (5)

and the dispersion relation can be detected through the momentum-dependent
cutoff. For distributions other than the Boltzmann, the exponential energy
factor is replaced by the appropriate function, e.g $\left[{\rm
e}^{(E-\mu)/kT}-1\right]^{-1}$ for a Bose-Einstein distribution with chemical
potential $\mu$. Phase space must be discretised carefully in units of the
Planck constant.

These principles apply to photons in a dye-filled microcavity. In the limit of
weak photon-photon interactions, as well as paraxial approximation, the
photons have an effective mass $m_{ph}$ given by the cavity length:
$m_{ph}=h\,n^{2}/(c\lambda_{0})$ where $\lambda_{0}$ is the wavelength (in
free space) associated with light of the energy of the lowest available cavity
mode Nyman14 , which we will call the cutoff wavelength. Here $n$ is the
refractive index of the medium, $c$ the speed of light in free space and $h$
Planck’s constant. The potential energy landscape is defined by the shape of
the cavity mirrors. It is straightforward to show that for a plano-spherical
cavity $V({\bf r})=\left(\frac{hc}{\lambda_{0}}\right)\frac{1}{2L_{0}R}r^{2}$
with $L_{0}$ being the length of the cavity at its greatest, and $R$ the
mirror radius of curvature, so that a harmonic trap is formed of
characteristic angular frequency $\Omega=(c/n)/\sqrt{L_{0}R}$. In a planar-
planar microcavity, the potential energy is a constant and $f({\bf
p},E)=\frac{1}{Z}{\rm
e}^{\left[{-\frac{{E}}{k_{B}T}}\right]}\delta\left[{E-Q({\bf p})}\right]$,
i.e. the photoluminescence is concentrated on a line which follows the
dispersion relation, weighted by the Boltzmann factor.

Dispersion relations have been measured in microcavities, for example showing
Bose-Einstein condensation of exciton-polaritons Kasprzak06  and organic
polaritons Daskalakis14 ; Plumhof14 , the Bogoliubov dispersion relation for
interacting Bose condensates Utsunomiya08  as well as the deviation from the
simple massive-particle dispersion induced by strong coupling between light
and matter Houdre94 ; Lidzey99 . It is worth noting that in ultracold atoms
Stenger99 ; Steinhauer03  and superfluid helium Donnelly81 , an equivalent
quantity is often referred to as the structure factor.

##  II Experimental setup

Figure 1: Left: a schematic of the optics for momentum- and position-resolving
spectroscopy. Double solid-dashed lines are beamsplitters. Lenses are marked
with double-headed arrows. Right: an example momentum-resolved spectrum. The
only post-processing performed on the colour image is taking the square-root
for display purposes. Calibration of the axes is described later in this
manuscript. Dark lines are dirt stuck on the slit in momentum space.

The core of our experimental apparatus consists of the same optically-pumped,
dye-filled microcavity described in our previous articles Marelic15 ;
Marelic16 : see Fig. 1 333The planar mirror used is the same as before, but
the spherical mirror has a different specification for the reflectivity
coating. It now lets through 30 ppm of light, marginally limiting the cavity
finesse but giving increased signal strength.. The light is chopped into 500
ns pulses using an acousto-optical modulator (AOM), with the pulse repetition
rate lower for higher pump powers, ensuring that few dye molecules are shelved
in the unwanted triplet state. The cavity thermalised light is in one
longitudinal mode, typically the 8th mode, so that the cavity is 4 wavelengths
of light long. The photoluminescence is imaged to infinity by a objective
lens, and the light split on a beamsplitter (shown as a double solid-dashed
line), with half sent to be shared among an imaging camera and a commercial
spectrometer. The rest of the light goes to the resolving spectroscopy optics,
which can be configured in multiple ways for measurement and calibration. The
main resolving spectroscopy components are: a telescope which creates images
for position- or momentum-spatial filtering; a diffraction grating to separate
wavelengths of light; and cylindrical optics which focus the light onto the
spectroscopy camera. Intra-cavity in-plane momentum corresponds to angle of
emission, implying that the momentum image is the Fourier transform of the
position-space.

The light can also be sent directly to a diagnostic camera, which is initially
focussed on position-space features of the microcavity, usually one edge of
the 1 mm diameter planar mirror. A 1:1 telescope is then inserted into the
optical path, and focussed such that the image on the diagnostic camera is re-
created. The telescope uses a pair of 50 mm focal length achromatic doublets
(the same as the objective) to reduce aberrations, both chromatic and off-
axis.

A slit is placed in the telescope, in one of two positions, either the
position- or momentum-space plane. The 20 $\mu{\rm m}$ wide slit is vertical,
reducing a 2D image to one dimension. To find the position-space plane, we
send photoluminescence through the slit and move the slit until its image is
in focus on the diagnostic camera. We use collimated light, inserted just
before the objective as a single-momentum light source for placing the slit in
momentum-space. To focus the space (vertical) axis of the spectroscopy camera
we rotate the slit to horizontal, and focus transmitted photoluminescence to a
thin horizontal line, then rotate the slit back to vertical.

We focus the energy (horizontal) axis after the vertical axis, and use a
photon Bose-Einstein condensate, which is known to be narrow-band just above
threshold pump power Schmitt16  as a single-energy reference light source.
Calibration of energy uses the commercial spectrometer as a reference to
generate two known pixel-energy pairs. We then extend the cavity to about 100
wavelengths of light, and use the many longitudinal modes as a comb of
frequency references, to which we fit a cubic calibration of energy as a
function of pixel number. The calibration is nearly linear, with the gradient
varying by no more than 10% across the image. The wavelength resolution is
limited by vibrations of the microcavity, and can be as large as 0.5 nm. For
the very shortest cavities (not used in this article) the mirrors bump
together reducing vibrations, leaving the resolution, measured using a photon
condensate, as good as 0.05 nm, which is equivalent to just 2 pixels on
camera.

##  III Momentum space

###  III.1 Calibration of momentum

The back focal plane of the objective lens is a momentum-space image plane.
The correspondence between displacement from the optic axis in the back focal
plane $y_{bfp}$ and vertical intracavity momentum $p_{y}$ is:

$\displaystyle p_{y}=y_{bfp}\times\frac{h}{\lambda_{0}f^{\prime}}$ (6)

where $f^{\prime}$ is the focal length of the objective lens. In deriving this
relation, we are required to make use of the Minkowski formula for canonical
momentum ($h/\lambda$) in a medium, not the Abraham formula for kinetic
momentum Barnett10 . Near the back focal plane of the objective, we can insert
a calibrated graticule of black lines printed on acetate. The known size of
the graticule features provides the calibration, even though the magnification
of the imaging system from back focal plane to camera is not a priori known.
The measured magnification is very insensitive to the exact position of the
graticule, so we place it only approximately in the back focal plane of the
objective.

In Fig. 2 we show, with the graticule in place, an image in 2D momentum space
(left) and the momentum-resolved spectrum of the cavity photoluminescence
(right). The images are real colour, but the square-root of the raw value for
display. The graticule had a pitch of 800 $\mu{\rm m}$, the cutoff wavelength
we use is $\lambda_{0}\simeq 588$ nm and the objective has a nominal focal
length of 50 mm, so from one line to the next corresponds to a difference of
$1.8\times 10^{-29}$ N s. The momentum resolution in the centre is
approximately 5 times smaller than the spacing between lines: $4\times
10^{-30}$ N s.

Figure 2: Momentum calibration using a graticule with 800 $\mu{\rm m}$ pitch
in the back focal of the imaging objective. Each line corresponds to a
difference of in-plane $1.8\times 10^{-29}$ N s. Left: a momentum space image.
Right: momentum-resolved spectroscopy, with the horizontal axis being energy
increasing to the left.

Both images show distortions. The alignment of the light path, especially
through the telescope, is centred to avoid astigmatism. The central
undistorted portion shows a linear calibration between camera position and
momentum up to about $4\times 10^{-29}$ N s. For larger momenta, some
pincushion distortion which is typical of spherical aberration is noticeable.
The momentum-resolved spectrum shows a change in magnification with energy
relative to the cutoff energy, which persists when the cutoff energy is
changed by changing the cavity length. This distortion means that our
calibration is uncertain to within about 30%, and it probably comes from
position-space off-axis aberrations.

###  III.2 Momentum-resolved spectra: $n({\bf p},E)$

In Fig. 3 we show calibrated, false-colour momentum-resolved spectra, $n({\bf
p},E)$. Background light as well as camera hot pixels are eliminated by
subtracting a background image, taken under the same conditions but with the
pump light shut off by the AOM.

Figure 3: Momentum-resolved spectra $n({\bf p},E)$ of photoluminescence below
threshold (top) and BEC just above (middle) and far above threshold (bottom).
Pump spot size was $70\,\mu{\rm m}$, optimised to produce a thermal
distribution of photons at room temperature. The grey dashed lines are
dispersion relations for particles of mass $7.8\times 10^{-36}$ kg. The blue
dot-dashed lines on the lower two panels are the Bogoliubov dispersion
relations assuming a dimensionless two-dimensional interaction parameter
$\tilde{g}=10^{-4}$ (see text for more details).

The expected dispersion relation for microcavity photons can be thought of as
the requirement for light with in-plane momentum to satisfy the boundary
conditions at the mirror surfaces. Alternatively, in the Schrödinger equation
which describes the motion of the photons, the effective mass comes from the
coefficient of the spatial derivative term, i.e. the kinetic energy. The grey
parabola on the top panel of Fig. 3 represents the dispersion relation for a
free particle of mass $m_{ph}=7.8\times 10^{-36}$ kg, using $n=1.44$ for the
refractive index of the intracavity medium (mostly ethylene glycol).

At higher pump powers, Bose-Einstein condensation occurs. We see in Fig. 3
(lower panels) that the macroscopic occupation occurs in the lowest energy
state, for low momentum.

Below threshold, we have tested the validity of Eqn. (5). For each momentum,
we fitted the spectrum with an exponential decay with a cutoff. The cutoff is
directly related to the dispersion relation. From Fig. 3 (top panel), it is
clear that the dispersion relation follows that of a massive particle for
small momenta. For large momenta (large emission angles), the pincushion
distortion of our imaging system prevents us from drawing any conclusion. The
exponential decay constant is the inverse temperature, as plotted in Fig. 4.
The temperature is largely independent of momentum, compatible with Eqn. (5),
but the temperature does depend on the pump spot size. It is known that
equilibration with room temperature is not complete, but limited by
spontaneous emission from the dye molecules out of the cavity Keeling16 . A
similar conclusion can be drawn from the observation that the position-space
size of the photon cloud depends on the pump spot size Marelic15 . The
implication is that even when thermal equilibrium is broken, ergodicity is
still respected, i.e. a distribution in terms of energy alone is enough to
describe the system.

Figure 4: Temperature vs. momentum for several pump spot sizes, derived from
momentum-resolved spectra. Thermalisation is not complete, and but larger spot
sizes lead to higher temperatures. For comparison, the thermal-equilibrium
cloud size is 70 $\mu{\rm m}$. Temperature is independent of momentum,
implying that populations depend only on energy so ergodicity is respected.

####  III.2.1 Effective mass

If the photons are strongly coupled to molecules of mass $m_{molecule}$, they
form polaritons, whose effective masses can be written in terms of the
Hopfield coefficients Deng10 :
$|X_{\pm}|^{2}=\frac{1}{2}\left\\{1\pm\left[\frac{\Delta}{\sqrt{\Delta^{2}+4g_{0}^{2}}}\right]\right\\}$
where $\Delta$ is the cavity detuning from the molecular resonance, known as
the zero-photon line (ZPL) and $g_{0}$ the light-molecule coupling strength.
The effective masses of the two kinds of polaritons (usually known as “upper”
and “lower” polaritons) are

$\displaystyle\frac{1}{m_{\mp}}=\frac{|X_{\pm}|^{2}}{m_{molecule}}+\frac{|X_{\mp}|^{2}}{m_{ph}}.$
(7)

Within the range of validity of the calibration of momentum, the below-
threshold dispersion relation matches the kinetic energy of a microcavity
photon of the expected mass, $m_{ph}$. We have measured this effective mass
for cutoff wavelengths between 520 and 610 nm. Data were taken within 1 nm of
the ZPL at 545 nm. We found the effective mass to be compatible with the bare
photon mass to within the variation of wavelength calibration: about 30%
across the camera. From this, we conclude that the coupling strength $g_{0}<1$
THz, and that the particles involved are really un-dressed cavity photons and
not polaritons.

From the known spontaneous emission lifetime of 4 ns and the volume of the
smallest cavity mode, $4\times 10^{-17}$ m3, we expect the coupling of a
single molecule to the cavity mode to be of order 1 GHz Loudon . There are at
least $10^{7}$ molecules within a cavity mode. It seems that there is no
collective enhancement of the light-matter coupling, as expected here since
the same mechanism which helps photons thermalise also de-phases the
electronic states of the molecules.

####  III.2.2 Photon-photon Interactions

Interactions are expected to cause the dispersion relation to deviate from
parabolic, following the Bogoliubov dispersion. In a local-density
approximation, that dispersion relation is:

$\displaystyle Q({\bf
p})=\frac{p^{2}}{2m_{ph}}\left(1+\frac{4\hbar^{2}\tilde{g}\rho_{ph}}{p^{2}}\right)^{1/2}$
(8)

where $\tilde{g}$ is the dimensionless two-dimensional interaction parameter
and $\rho_{ph}$ the number density of photons in the centre of the condensate.

The number of photons in the condensate $N_{C}$ was estimated from the
commercial spectrometer data (not shown), being about 12% and 70% of the non-
condensate photon number for middle (just above threshold) and bottom (far-
above threshold) panels respectively. We assume that the total number of non-
condensate photons is given by the equilibrium saturation value
$N_{th}=\frac{pi^{2}}{6}\left(\frac{k_{B}T}{\hbar\Omega}\right)^{2}\simeq
27\,000$, for mode number $q=8$ and mirror radius of curvature 250 mm. The
central two-dimensional density is given by the condensate number and harmonic
oscillator length $a_{ho}=\sqrt{\hbar/m_{ph}\Omega}\simeq 6\,\mu{\rm m}$ here:
$\rho_{ph}=N_{C}/2\pi a_{ho}^{2}$.

In Fig. 3 we have overlaid the Bogoliubov dispersion relation (blue dot-dashed
lines) onto the momentum-resolved spectra, for $\tilde{g}=10^{-4}$. Larger
values of $\tilde{g}$ are not compatible with the observed data in the bottom
panel. Given that the local-density approximation is not valid for such weak
interactions, that the condensate is probably multimode, and that
imperfections in the cavity mirrors distort the condensate, we conclude a more
conservative, loose upper bound of $\tilde{g}\lesssim 10^{-3}$. This value is
consistent with previously published estimates based on experimental data
which are in the range $7\times 10^{-4}$ Klaers10b  to $5\times 10^{-7}$
VanDerWurff14 . Interactions might alternatively be inferred by observing
frequency shifts in anisotropic potentials deLeeuw16 ; Vyas14 .

###  III.3 Momentum-space images: $n({\bf p})$

Using the diagnostic camera, we obtain momentum-space images $n({\bf p})$, as
shown in Fig. 5. A broad thermal cloud is seen. Above threshold, a condensate
peak appears on top of the thermal cloud. At very high pump powers, multimode
condensation occurs, which leads to a distorted, broadened condensate peak.
The same sorts of peaks are observed in position-space imaging, also
broadening when multimode condensation is seen. Since position-space and
momentum-space are mutually Fourier transforms, we conclude that the multiple
modes are incoherent, in agreement with Ref. Marelic16 .

Figure 5: Momentum-space ${\bf p}=(p_{x},p_{y})$ images, $n({\bf p})$. As pump
power is increased (top to bottom), to the thermal distribution (top panel) is
added a momentum-narrow condensate (middle), followed by a broadened multi-
mode condensate (bottom panel). The calibration is accurate for $|{\bf
p}|<4\times 10^{-29}$ Ns. Figure 6: Cuts through the same data as Fig. 5, top
two panels. The greyed out region is known to suffer from pincushion
distortion, rendering the calibration inaccurate. The dashed line is a thermal
distribution at room temperature. The dot-dashed line is the expected momentum
distribution for the non-interacting Bose-Einstein condensate.

The width of the condensate peak is approximately equal to the characteristic
momentum width of the harmonic oscillator $p_{ho}=\sqrt{\hbar\Omega m_{ph}}$,
since photon-photon interactions are negligible. For the experiment of Fig. 5,
$p_{ho}=1.1\times 10^{-29}$ Ns. The typical thermal momentum scale for
temperature $T$ is $\sqrt{m_{ph}k_{B}T}=1.8\times 10^{-28}$ Ns. In Fig. 6 we
see cuts through the momentum space images just below and above condensation
threshold (the same data as the top two panels of Fig. 5). The greyed out
sections show significant distortion, so the calibration therein of vertical
momentum $p_{y}$ is not trustworthy. Gaussian curves of widths appropriate to
the thermal (dashed line) and the condensate (dash-dotted line) momentum size
are also shown. In the calibrated region the width of the condensate is
compatible with the harmonic oscillator momentum width and the below-threshold
momentum width is compatible with the typical thermal momentum scale.

Repulsive interactions would increase the position-space size of the
condensate, which decreases the momentum-space size. We have seen no such
effect, which confirms our bound of $\tilde{g}\lesssim 10^{-3}$ based on the
dispersion relation.

##  IV Position space

We now turn to position-space resolved spectra. Where the cavity length
locally varies, the energy of light which matches the mirror boundary
conditions also varies locally. Thus, mirror shape is mapped to potential
energy landscape, $V({\bf r})$. Similar spectra have been demonstrated for
condensates of exciton-polaritons when a confining potential is applied
Balili07 .

###  IV.1 Calibration of position

One of the two cavity mirrors is planar, the other spherical with a radius of
curvature $R=0.25$ m. The light penetrates the dielectric mirrors, which have
a different refractive index to the dye-solvent mixture, which means that the
optical path length is not exactly the same as the true cavity length. We
write the cavity length as a function of longitudinal mode number $q$ so that

$\displaystyle L_{0}(q)=\frac{\lambda^{*}q}{2}-L_{off}$ (9)

where $L_{off}$ is the distance the light penetrates the mirrors and
$\lambda^{*}=\lambda_{0}/n$ is the wavelength of light in the intracavity
medium. In principle $L_{off}$ is also a function of wavelength, but in
practice it makes no difference if we take it to be a constant here.

We can use Eqn. (9) to calibrate the position-space magnification of the
position-resolved spectra. We expect the potential energy landscape to appear
as a parabolic cutoff energy as a function of position, as shown in Eqn. (4).
The parabolic coefficient in object space $V(r)=ar^{2}$ is given by:

$\displaystyle\frac{1}{a(q)}={R}\left(\frac{\lambda_{0}}{hc}\right)\left[{\lambda^{*}q}-2L_{off}\right]$
(10)

The energy calibration of the spectrometer is performed as for the momentum-
resolved spectroscopy. In Fig. 7 (top) we show the measured inverse of the
parabolic curvature $a$ in nm per (camera pixel)2 for various longitudinal
mode numbers. From the fit to Eqn. (10) both the mirror penetration and
magnification of the system are determined. An example position-resolved
spectrum with magnification calibrated is shown in Fig. 7 (bottom). The fitted
parabola is overlaid (dashed line).

Figure 7: Top: The measured parabolic curvature $a$ of the position dependence
of the confining potential for photons depends on the longitudinal mode number
(dots). Data were all taken with cutoff $\lambda_{0}=591\pm 0.7$ nm. A fit to
Eqn. (9) (solid line) allows us to calibrate the magnification in the
position-resolved spectra. Part of the light penetrates the cavity mirrors, to
an optical depth equivalent to $L_{off}$, also obtained from the fit. Bottom:
an example calibrated position-resolved spectrum $n({\bf r},E)$ with the
fitted parabolic potential.

###  IV.2 Position-resolved spectra $n({\bf r},E)$ and images $n({\bf r})$

We show in Fig. 8 (left) typical position-resolved spectra $n({\bf r},E)$ of a
our dye-microcavity photoluminescence far below threshold, just above
threshold for single-mode condensation and far above threshold for multimode
condensation. On the right are position-space images of the same conditions.

The photoluminescence is most intense near the potential energy, whereas Eqn.
(4) predicts uniform intensity for each energy, for all energies higher than
the minimum set by the potential energy. The explanation is that only
photoluminescence with low kinetic energy, i.e. small emission angle, is
admitted by our optics. The slit in the telescope, positioned at a real-space
image, acts in part as an anisotropic momentum-space filter. A side effect is
that the position-resolution is no better than about 15 $\mu{\rm m}$ in the
resolved spectrometer.

Figure 8:  Position-resolved spectra $n({\bf r},E)$ (left side) and images
$n({\bf r})$ with ${\bf r}=(x,y)$ (right side) of photoluminescence below
threshold (top) and BEC just above (middle) and far above threshold (bottom).
The parabolic curve of the harmonic potential energy is clearly visible.
Condensation occurs at the centre, where potential energy is lowest. At higher
pump powers, multimode condensation is observed.

##  V Discussion

###  V.1 Limitations of the optical setup

The experiment was conceived in an attempt to measure the Bogoliubov
dispersion relation which gives information about the photon-photon
interactions. In Ref. Nyman14  it was shown that wavelength energy resolution
0.04 nm and momentum resolution $1.4\times 10^{-30}$ N s would be sufficient
to measure a dimensionless interaction parameter $\tilde{g}=10^{-5}$, which is
in the expected range VanDerWurff14 . The measured wavelength resolution
depends on the mechanical stability of the cavity, which ranges from 0.05 nm
(for $q$ so small that the mirrors touch) to 0.5 nm (in more normal
operation). The minimum position and momentum resolution we have achieved here
are 15 $\mu{\rm m}$ and $4\times 10^{-30}$ N s respectively. For comparison,
the smallest cavity mode covers 6 $\mu{\rm m}$ and $1.1\times 10^{-29}$ N s:
we resolve well in momentum but not position. The wavelength resolution is
sufficient to resolve interactions at the $10^{-5}$ level, but not so
momentum. The multimode condensation phenomenon also makes interpretation in
terms of interactions tricky.

The resolution limits in position are set mainly by spherical aberrations in
the telescope. The achromatic doublets used give distortions for angles
greater than about 35 mrad BornWolf , equivalent to $4\times 10^{-29}$ N s,
exactly where our calibration is seen to become invalid. Resolution could be
improved by replacing the achromatic doublets with lenses which correct more
aberrations or by removing the telescope altogether. The latter option would
require a slit in the back focal plane of the objective, precluding
simultaneous momentum- or position-resolved spectroscopy and imaging as
presented in Fig. 8. In addition, momentum-space imaging optics could be moved
closer, occupying a large numerical aperture, improving the momentum-space
resolution.

The momentum resolution is in practice limited by the inverse of the size of
the lowest cavity mode. To improve effective momentum resolution we would
therefore use larger radius of curvature mirrors. More gently curved mirrors
also cause there to be more photons at threshold, enhancing the effects of
interactions, at the cost of decreasing the energy range over which
interaction effects occur. The energy resolution can be improved to compensate
by using a more finely spaced grating than the present 1800 lines/mm.

###  V.2 Conclusions

We have observed various projections of the phase-space distribution of dye-
microcavity thermalised and condensed photons: position- and momentum-
resolved, with and without energy resolution. The momentum-resolved spectrum
gives the dispersion relation, which results in upper bounds for both the
photon-photon and photon-molecule interaction strengths. Both are very weak,
indicating that the light is acting as free photons. Reasonable improvements
to the experimental setup might well make it possible to directly measure the
photon-photon interaction strength, which is still not well known for dye-
microcavity photons. Since we find that the effective temperature is not a
phase-space local quantity, it is clear that ergodicity is retained even when
the temperature is not room temperature because of incomplete thermalisation.

##  VI Acknowledgements

We thank the UK Engineering and Physical Sciences Research Council for
supporting this work through fellowship EP/J017027/1 and the Controlled
Quantum Dynamics CDT EP/L016524/1. The data underlying this article, as well
as the scripts used for analysis, are available at
http://dx.doi.org/10.5281/zenodo.60901.

## References

  * (1) J. W. Gibbs, Elementary principles in statistical mechanics (Scribner, New York, 1902), available as e-book from Project Gutenberg.
  * (2) We would have cited Ludwig Boltzmann here but our German language skills are insufficient to find which of his articles is most appropriate.
  * (3) K. Huang, Statistical mechanics (Wiley, New Jersey, USA, 1987).
  * (4) H. Wu, E. Arimondo, and C. J. Foot, Phys. Rev. A 56, 560 (1997).
  * (5) E. T. Jaynes, Phys. Rev. 106, 620 (1957).
  * (6) O. J. Luiten, M. W. Reynolds, and J. T. M. Walraven, Phys. Rev. A 53, 381 (1996).
  * (7) H. Wu and C. J. Foot, Journal of Physics B: Atomic, Molecular and Optical Physics 29, L321 (1996).
  * (8) J. Kasprzak, M. Richard, A. Baas, B. Deveaud, R. André, J.-P. Poizat, and L. S. Dang, Phys. Rev. Lett. 100, 067402 (2008).
  * (9) J. Klaers, F. Vewinger, and M. Weitz, Nature Physics 6, 512 (2010).
  * (10) J. Klaers, J. Schmitt, F. Vewinger, and M. Weitz, Nature 468, 545 (2010).
  * (11) A.-W. de Leeuw, H. T. C. Stoof, and R. A. Duine, Phys. Rev. A 88, 033829 (2013).
  * (12) A. Chiocchetta and I. Carusotto, Phys. Rev. A 90, 023633 (2014).
  * (13) R. A. Nyman and M. H. Szymańska, Phys. Rev. A 89, 033844 (2014).
  * (14) E. C. I. van der Wurff, A.-W. de Leeuw, R. A. Duine, and H. T. C. Stoof, Phys. Rev. Lett. 113, 135301 (2014).
  * (15) D. Snoke and S. Girvin, J. Low Temp. Phys. 171, 1 (2013).
  * (16) P. Kirton and J. Keeling, Phys. Rev. Lett. 111, 100404 (2013).
  * (17) P. Kirton and J. Keeling, Phys. Rev. A 91, 033826 (2015).
  * (18) J. Keeling and P. Kirton, Phys. Rev. A 93, 013829 (2016).
  * (19) J. Schmitt, T. Damm, D. Dung, F. Vewinger, J. Klaers, and M. Weitz, Phys. Rev. A 92, 011602 (2015).
  * (20) J. Marelic and R. A. Nyman, Phys. Rev. A 91, 033813 (2015).
  * (21) J. Marelic, L. F. Zajiczek, H. J. Hesten, K. H. Leung, E. Y. X. Ong, F. Mintert, and R. A. Nyman, arXiv 1510.05562 (2015), to appear in New J. Phys.
  * (22) T. Damm, J. Schmitt, Q. Liang, D. Dung, F. Vewinger, M. Weitz, and J. Klaers, Nature communications 7, (2016).
  * (23) J. Schmitt, T. Damm, D. Dung, C. Wahl, F. Vewinger, J. Klaers, and M. Weitz, Phys. Rev. Lett. 116, 033604 (2016).
  * (24) Here we are making an approximation that the phase-space distribution varies slowly compared to discreteness of states, i.e. Planck’s constant is small.
  * (25) L. A. Lugiato and R. Lefever, Phys. Rev. Lett. 58, 2209 (1987).
  * (26) R. Y. Chiao and J. Boyce, Phys. Rev. A 60, 4114 (1999).
  * (27) R. Y. Chiao, Opt. Commun. 179, 157 (2000).
  * (28) J. Kasprzak, M. Richard, S. Kundermann, A. Baas, P. Jeambrun, J. Keeling, F. Marchetti, M. Szymańska, R. Andre, J. Staehli, et al., Nature 443, 409 (2006).
  * (29) K. Daskalakis, S. Maier, R. Murray, and S. Kéna-Cohen, Nat. Mater. 13, 271 (2014).
  * (30) J. D. Plumhof, T. Stöferle, L. Mai, U. Scherf, and R. F. Mahrt, Nat. Mater. 13, 247 (2014).
  * (31) S. Utsunomiya, L. Tian, G. Roumpos, C. Lai, N. Kumada, T. Fujisawa, M. Kuwata-Gonokami, A. Löffler, S. Höfling, A. Forchel, et al., Nature Physics 4, 700 (2008).
  * (32) R. Houdré, C. Weisbuch, R. P. Stanley, U. Oesterle, P. Pellandini, and M. Ilegems, Phys. Rev. Lett. 73, 2043 (1994).
  * (33) D. G. Lidzey, D. D. C. Bradley, T. Virgili, A. Armitage, M. S. Skolnick, and S. Walker, Phys. Rev. Lett. 82, 3316 (1999).
  * (34) J. Stenger, S. Inouye, A. P. Chikkatur, D. M. Stamper-Kurn, D. E. Pritchard, and W. Ketterle, Phys. Rev. Lett. 82, 4569 (1999).
  * (35) J. Steinhauer, N. Katz, R. Ozeri, N. Davidson, C. Tozzo, and F. Dalfovo, Phys. Rev. Lett. 90, 060404 (2003).
  * (36) R. J. Donnelly, J. A. Donnelly, and R. N. Hills, Journal of Low Temperature Physics 44, 471 .
  * (37) The planar mirror used is the same as before, but the spherical mirror has a different specification for the reflectivity coating. It now lets through 30 ppm of light, marginally limiting the cavity finesse but giving increased signal strength.
  * (38) S. M. Barnett, Phys. Rev. Lett. 104, 070401 (2010).
  * (39) H. Deng, H. Haug, and Y. Yamamoto, Rev. Mod. Phys. 82, 1489 (2010).
  * (40) R. Loudon, The quantum theory of light (OUP , ADDRESS, 2000).
  * (41) A.-W. de Leeuw, E. C. I. van der Wurff, R. A. Duine, D. van Oosten, and H. T. C. Stoof, Phys. Rev. A 94, 013615 (2016).
  * (42) V. M. Vyas, P. K. Panigrahi, and J. Banerji, Physics Letters A 378, 1434 (2014).
  * (43) R. Balili, V. Hartwell, D. Snoke, L. Pfeiffer, and K. West, Science 316, 1007 (2007).
  * (44) M. Born and E. Wolf, Principles of optics (Cambridge University Press, Cambridge, UK, 2000).


# Microwave control of Rydberg atom interactions

S Sevinçli1,3 and T Pohl2 1Department of Physics and Astronomy, Aarhus
University, 8000 Aarhus C, Denmark  2 Max Planck Institute for the Physics of
Complex Systems, Nöthnitzer Strasse 38, 01187 Dresden, Germany
<EMAIL_ADDRESS>$3$ Present address: Department of Physics, İzmir
Institute of Technology, IZTECH, TR35430, İzmir, Turkey.

###### Abstract

We investigate the interaction between Rydberg atoms, whose electronic states
are dressed by multiple microwave fields. Numerical calculations are used for
an exact description of the microwave induced interactions, and employed to
benchmark a perturbative treatment that yields simple insights into the
involved mechanisms. Based on this theory, we demonstrate that microwave
dressing provides a powerful approach to control dipolar as well as van der
Waals interactions and even permits to turn them off entirely. In addition,
the proposed scheme also opens up possibilities for engineering dominant
three-body interactions.

###### pacs:

32.80.Ee, 32.80.Rm

##  1 Introduction

In recent years the physics of cold Rydberg atoms has attracted widespread
interest due to diverse applications, ranging from quantum information science
[1, 2, 3, 4] and quantum simulations of magnetism [5, 6] over precision
measurements [7] to nonlinear quantum optics [8, 9]. Much of these prospects
stem from the exaggerated properties of Rydberg states and, in particular,
their strong mutual interactions. Due to the strong $C_{6}\sim n^{11}$ scaling
of the van der Waals coefficient, $C_{6}$, with the atoms principal quantum
number, $n$, Rydberg-Rydberg atom interactions can exceed those of ground
state atoms by many orders of magnitude [10, 11, 12, 13, 14]. Even at
comparatively large distances of several $\mu$m the resulting interactions can
exceed all remaining energy scales in the system, and, e.g., lead to a strong
excitation blockade that inhibits simultaneous excitation of Rydberg atoms
within a characteristic blockade radius $R_{\rm bl}$. In addition, the strong
$n$-dependence provides a convenient way to tune the interactions and study
the transition between the weak and strong interaction regime via moderate
changes of the principle quantum number.

On the other hand, the ability to directly tune interactions in a state-
selective manner would open up new possibilities for controlling the
properties and time evolution of strongly correlated Rydberg systems and,
thereby, add a valuable degree of flexibility to the aforementioned
applications. Here one can make use of the high susceptibility of Rydberg
states to electric and magnetic fields [15], which permits to manipulate their
interactions by applying moderate external fields. Experiments have explored
effects of static magnetic and electric [14, 16, 17, 18, 19] fields, where the
latter can be used to induce static dipole moments or dipole-coupled pair
resonances, both resulting in direct dipole-dipole interactions. Time-varying,
i.e. microwave fields, offer a more refined way of control, since they permit
to couple different Rydberg levels in a state selective manner. This has been
exploited in theoretical studies of polar molecules [20, 21, 22], showing how
the combination of a static electric and a microwave field can be used to
shape molecular interactions. For cold Rydberg atoms resonant coupling to a
single microwave field has been investigated theoretically [23] and
experimentally [24, 25, 26, 27, 28], demonstrating an enhanced excitation
blockade due to field-induced dipole-dipole interactions. Recent work
demonstrates that microwave control of Rydberg states provides a powerful
tool, e.g., for manipulating the atomic sensitivity to external fields [29],
implementing quantum gates [28] or to probe and control Rydberg state dynamics
in atomic beam experiments [30].

In this article, we study far off resonant Rydberg state coupling by a
combination of several microwave fields, which is shown to permit versatile
control of dipolar as well as van der Waals potentials and to realize
conditions with dominating three-body interactions. Numerical calculations for
two interacting atoms combined with a perturbative treatment of many-body
systems provide an exact description as well as an intuitive understanding of
the field induced interactions. We discuss several types of realizable
interactions and provide simple analytical formulae for the corresponding
microwave parameters, including conditions for which the long-distance tail of
the interactions can turned off entirely.

The paper is organized as follows. First (section 2), we describe the
considered system consisting of an ensemble of Rydberg atoms and introduce the
different interaction terms arising from microwave driving and dipole-dipole
coupling of different Rydberg states. In section 3 we describe our numerical
procedure for calculating the dressing-induced interaction between two atoms.
Finally we present a perturbative many-body treatment of the resulting
interaction, which is compared to the exact calculations and used to identify
the required microwave parameters for obtaining different types of
interactions in section 4 .

Figure 1: (a)Two-photon excitation scheme for three-level system. States are
resonantly coupled by Rabi frequencies $\Omega_{p}$ and $\Omega_{c}$ where
$\gamma$ stands for spontaneous decay from the intermediate state. (b)
Microwave dressing of Rydberg states. Two linearly polarized and two
circularly polarized fields couple $|ns\rangle$ state to
$|p_{m}\rangle=|npm\rangle$ and $|p^{\prime}_{m}\rangle=|(n-1)pm\rangle$
states with Rabi frequencies $\Omega_{p},\Omega_{\sigma}$ and
$\Omega_{p}^{\prime},\Omega_{\sigma}^{\prime}$; detunings $\Delta$ and
$\Delta^{\prime}$, respectively.

##  2 Microwave dressing of Rydberg states

We consider an ensemble of cold atoms at positions $\mathbf{R}_{i}$,
$i=1,\ldots,N$. We focus on $ns$-Rydberg states ($|r\rangle$) of alkaline
atoms, which can be excited from the atomic ground state ($|g\rangle$) by two-
photon coupling via an intermediate $p$-state ($|e\rangle$) with Rabi
frequencies $\Omega_{p}$ and $\Omega_{c}$ (see figure 1a). For alkaline atoms
the van der Waals interaction between two $nS_{1/2}$-Rydberg states is
repulsive and predominantly originates from off-resonant dipole-dipole
coupling to ($np,(n-1)p$) pair states [31, 32]. Four different microwave
fields couple the atomic $ns$ state to these two $p$-states. As shown in
figure 1b, we consider two pairs with linear and circular polarization fields,
which have a frequency of $\omega$ and $\omega^{\prime}$, but different
amplitudes $F_{\pi}$, $F_{\pi}^{\prime}$ and $F_{\sigma}$,
$F_{\sigma}^{\prime}$, respectively. The two microwave frequencies for each
pair of fields are chosen near-resonant with the transition between the
$|ns\rangle$-Rydberg state and the two adjacent $|p\rangle$-states. In
particular, $\omega$ is assumed to be close to the
$|ns\rangle\rightarrow|np\rangle$ transition energy, while $\omega^{\prime}$
is assumed to be close to the $|ns\rangle\rightarrow|(n-1)p\rangle$ transition
energy, such that we can define two detunings $\Delta$ and $\Delta^{\prime}$
for the upper and lower transitions, respectively. With the radial matrix
elements $\mu$ and $\mu^{\prime}$ of the $ns\leftrightarrow np$ and
$ns\leftrightarrow(n-1)p$ transitions, the field amplitudes define the
corresponding the Rabi frequencies $\Omega_{\pi}=F_{\pi}\mu/2$,
$\Omega_{\pi}^{\prime}=F_{\pi}^{\prime}\mu^{\prime}/2$,
$\Omega_{\sigma}=F_{\sigma}\mu/2$ and
$\Omega_{\sigma}^{\prime}=F_{\sigma}^{\prime}\mu^{\prime}/2$, as indicated in
figure 1b. In what follows, we will assume the detunings to be larger than the
respective Rabi frequencies, i.e. that the atoms are weakly dressed instead of
driven resonantly by the microwave fields.

The total Hamiltonian of the system can be split into three parts

$\hat{H}=\hat{H}_{\rm A}+\hat{H}_{\rm AF}+\hat{H}_{\rm AA}\;.$ (1)

The first term denotes the sum of single-atom Hamiltonians, and yields the
unperturbed energy spectrum for each of the $N$ atoms. Following the preceding
discussion the atom-microwave coupling is given by

$\displaystyle\hat{H}_{\rm AF}=$ $\displaystyle\sum_{i}$ $\displaystyle
F_{\pi}z_{i}\cos(\omega
t)+F_{\pi}^{\prime}z_{i}\cos(\omega^{\prime}t)+F_{\sigma}\left[x_{i}\cos(\omega
t)+y_{i}\sin(\omega t)\right]$ (2) $\displaystyle+$ $\displaystyle
F_{\sigma}^{\prime}\left[x_{i}\cos(\omega^{\prime}t)+y_{i}\sin(\omega^{\prime}t)\right]\;,$

where the field amplitudes and frequencies have been defined above. Finally,
the dipole-dipole interaction between the Rydberg atoms is given by

$\hat{H}_{\rm AA}=\sum_{i<j}\frac{{\bf r}_{i}{\bf
r}_{j}}{R_{ij}^{3}}-3\frac{({\bf r}_{i}{\bf R}_{ij})({\bf r}_{j}{\bf
R}_{ij})}{R_{ij}^{5}}\>,$ (3)

where ${\bf R}_{ij}={\bf R}_{j}-{\bf R}_{i}$ denotes the distance vector
between two Rydberg atoms at positions ${\bf R}_{i}$ and ${\bf R}_{j}$.

For a description of $ns-ns$ interactions we may approximately neglect spin-
orbit coupling, e.g. the fine structure of the Rydberg states [33, 34], such
that each atomic state is characterized by the three quantum numbers, n, l and
m. Denoting the relevant Rydberg states by $|s\rangle\equiv|n,l=0,m=0\rangle$,
$|p_{0,\pm}\rangle\equiv|n,l=1,m=0,\pm 1\rangle$ and
$|p^{\prime}_{0,\pm}\rangle\equiv|n-1,l=1,m=0,\pm 1\rangle$, we can write the
isolated-atom Hamiltonian in our reduced basis as

$\hat{H}_{\rm
A}=E_{p}\sum_{i}(\hat{\sigma}_{p_{-}p_{-}}^{(i)}+\hat{\sigma}_{p_{0}p_{0}}^{(i)}+\hat{\sigma}_{p_{+}p_{+}}^{(i)})+E_{p^{\prime}}\sum_{i}(\hat{\sigma}_{p_{-}^{\prime}p_{-}^{\prime}}^{(i)}+\hat{\sigma}_{p_{0}^{\prime}p_{0}^{\prime}}^{(i)}+\hat{\sigma}_{p_{+}^{\prime}p_{+}^{\prime}}^{(i)})\>,$
(4)

where
$\hat{\sigma}_{\alpha\alpha}^{(i)}=|\alpha^{(i)}\rangle\langle\alpha^{(i)}|$
denote single-atom projection operators of the $i$th atom onto a given state
$|\alpha\rangle=|l,m,n\rangle$. The energies $E_{p}$ and $E_{p}^{\prime}$ are
the energies of the $np$-states and $(n-1)p$ states relative to the
$ns$-Rydberg state, respectively, and define the microwave detunings
$\Delta=\omega-|E_{p}|$ and
$\Delta^{\prime}=\omega^{\prime}-|E_{p^{\prime}}|$.

##  3 Numerical description

In addition to an analytical perturbative treatment of the microwave-induced
interactions, we also performed numerical calculations that account for all
possible Rydberg state couplings due to microwave driving and dipole
interactions. To this end, we use a Floquet representation [35] of the
Hamiltonian, which permits to eliminate the time-dependence in (2) and yields
a time-independent Floquet matrix eigenvalue problem. For two atoms the
corresponding Hilbert space is spanned by the basis states
$|\alpha^{(1)},\alpha^{(2)},\nu_{\pi},\nu_{\sigma},\nu_{\pi}^{\prime},\nu_{\sigma}^{\prime}\rangle$
where $\alpha^{(i)}$ refers to the atomic states of each of the two atoms and
the quantum numbers $\nu_{\pi}$, $\nu_{\pi}^{\prime}$, $\nu_{\sigma}$ and
$\nu_{\sigma}^{\prime}$ correspond to the number of microwave photons with the
respective polarization and frequency. Since the chosen Rabi frequencies are
considerably smaller than the respective transition energies it suffices to
consider three adjacent photon numbers $\nu=-1,0,1$ for each microwave mode,
which makes for a total of $3969$ pair states. Upon diagonalizing the
resulting $3969\times 3969$ Floquet matrix, we obtain the microwave-dressed
interaction potentials and corresponding molecular eigenstates, among which we
are interested in the ($ns-ns$)-asymptote, i.e. the molecular potential curves
that asymptotically converges to $|s,s,0,0,0,0\rangle$ upon decreasing the
microwave amplitudes.

##  4 Perturbative approach

In order to develop an analytical description we first transform into a frame
of reference that rotates with the applied fields at the frequencies $\omega$
and $\omega^{\prime}$, by transforming the $N$-atom Hamiltonian according to
$\hat{H}\rightarrow\hat{U}\hat{H}\hat{U}^{\dagger}+i\frac{\partial\hat{U}}{\partial
t}\hat{U}^{\dagger}$ where

$U=\exp\left[\sum_{i}i\omega
t(\hat{\sigma}_{p_{-}p_{-}}^{(i)}+\hat{\sigma}_{p_{0}p_{0}}^{(i)}+\hat{\sigma}_{p_{+}p_{+}}^{(i)})-i\omega^{\prime}t(\hat{\sigma}_{p_{-}^{\prime}p_{-}^{\prime}}^{(i)}+\hat{\sigma}_{p_{0}^{\prime}p_{0}^{\prime}}^{(i)}+\hat{\sigma}_{p_{+}^{\prime}p_{+}^{\prime}}^{(i)})\right]\>.$
(5)

The transformed Hamiltonian contains time dependent terms that vary
proportional to $e^{2i\omega t}$, $e^{2i\omega^{\prime}t}$ and $e^{\pm
i(\omega-\omega^{\prime})t}$. If the difference between $\omega$ and
$\omega^{\prime}$ is sufficiently large we can neglect all of these fast
oscillating terms to obtain a time-independent Hamiltonian within the rotating
wave approximation. The single atom contribution is of the familiar form

$\hat{H}_{A}=-\Delta\sum_{i}(\hat{\sigma}_{p_{-}p_{-}}^{(i)}+\hat{\sigma}_{p_{0}p_{0}}^{(i)}+\hat{\sigma}_{p_{+}p_{+}}^{(i)})+\Delta^{\prime}\sum_{i}(\hat{\sigma}_{p_{-}^{\prime}p_{-}^{\prime}}^{(i)}+\hat{\sigma}_{p_{0}^{\prime}p_{0}^{\prime}}^{(i)}+\hat{\sigma}_{p_{+}^{\prime}p_{+}^{\prime}}^{(i)})\>,$
(6)

and the atom field coupling is given by

$\displaystyle\hat{H}_{AF}$ $\displaystyle=$
$\displaystyle\Omega_{\pi}\sum_{i}(\hat{\sigma}_{sp_{0}}^{(i)}+\hat{\sigma}_{p_{0}s}^{(i)})+\sqrt{2}\Omega_{\sigma}\sum_{i}(\hat{\sigma}_{sp_{+}}^{(i)}+\hat{\sigma}_{p_{+}s}^{(i)})$
(7) $\displaystyle+$
$\displaystyle\Omega_{\pi}^{\prime}\sum_{i}(\hat{\sigma}_{sp_{0}^{\prime}}^{(i)}+\hat{\sigma}_{p_{0}^{\prime}s}^{(i)})+\sqrt{2}\Omega_{\sigma}^{\prime}\sum_{i}(\hat{\sigma}_{sp_{+}^{\prime}}^{(i)}+\hat{\sigma}_{p_{+}^{\prime}s}^{(i)})\>.$

The resulting expression for the Rydberg-Rydberg atom interaction Hamiltonian
(3) is given in A.

In the limit where both the Rabi frequencies as well as the dipole-dipole
couplings are smaller than the microwave detunings we can apply perturbation
theory to calculate the long-range part of the dressing-induced interactions.
The first order contribution vanishes. In the laboratory frame, the second
order energy correction yields the standard van der Waals interaction [36, 37]

$E_{\rm lab}^{(2)}=\sum_{i<j}\frac{C_{6}}{R_{ij}^{6}}\;,$ (8)

with $C_{6}=-12\mu^{2}\mu^{\prime 2}/(E_{p}+E_{p^{\prime}})$. In the rotating
frame of reference, the second order energy correction

$E^{(2)}=N\left(\frac{\Omega_{\pi}^{2}}{\Delta}-\frac{\Omega_{\pi}^{\prime
2}}{\Delta^{\prime}}+\frac{2\Omega_{\sigma}^{2}}{\Delta}-\frac{2\Omega_{\sigma}^{\prime
2}}{\Delta^{\prime}}\right)$ (9)

just yields the combined light shift of the applied microwave fields. Atomic
interactions arise from the microwave dressing in third order

$\displaystyle E^{(3)}$ $\displaystyle=$
$\displaystyle\sum_{i<j}\frac{2\mu^{2}}{R_{ij}^{3}}\left[\left(\frac{\Omega_{\pi}}{\Delta}\right)^{2}\left(1-3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)+\left(\frac{\Omega_{\sigma}}{\Delta}\right)^{2}\left(2-3\frac{X_{ij}^{2}+Y_{ij}^{2}}{R_{ij}^{2}}\right)-6\frac{\Omega_{\pi}\Omega_{\sigma}}{\Delta^{2}}\frac{X_{ij}Z_{ij}}{R_{ij}^{2}}\right]$
(10) $\displaystyle+\sum_{i<j}\frac{2\mu^{\prime
2}}{R_{ij}^{3}}\left[\left(\frac{\Omega_{\pi}^{\prime}}{\Delta^{\prime}}\right)^{2}\left(1-3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)+\left(\frac{\Omega_{\sigma}^{\prime}}{\Delta^{\prime}}\right)^{2}\left(2-3\frac{X_{ij}^{2}+Y_{ij}^{2}}{R_{ij}^{2}}\right)-6\frac{\Omega_{\pi}^{\prime}\Omega_{\sigma}^{\prime}}{\Delta^{\prime
2}}\frac{X_{ij}Z_{ij}}{R_{ij}^{2}}\right]\>,$

and fourth order

$\displaystyle E^{(4)}$ $\displaystyle=$
$\displaystyle\sum_{i<j}\frac{2\mu^{4}}{\Delta
R_{ij}^{6}}\left[\left(\frac{\Omega_{\pi}}{\Delta}\right)^{2}\left(1+3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)+\left(\frac{\Omega_{\sigma}}{\Delta}\right)^{2}\left(2+3\frac{X_{ij}^{2}+Y_{ij}^{2}}{R_{ij}^{2}}\right)+6\frac{\Omega_{\pi}\Omega_{\sigma}}{\Delta^{2}}\frac{X_{ij}Z_{ij}}{R_{ij}^{2}}\right]$
(11) $\displaystyle-\sum_{i<j}\frac{2\mu^{\prime
4}}{\Delta^{\prime}R_{ij}^{6}}\left[\left(\frac{\Omega_{\pi}^{\prime}}{\Delta^{\prime}}\right)^{2}\left(1+3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)+\left(\frac{\Omega_{\sigma}^{\prime}}{\Delta^{\prime}}\right)^{2}\left(2+3\frac{X_{ij}^{2}+Y_{ij}^{2}}{R_{ij}^{2}}\right)+6\frac{\Omega_{\pi}^{\prime}\Omega_{\sigma}^{\prime}}{\Delta^{\prime
2}}\frac{X_{ij}Z_{ij}}{R_{ij}^{2}}\right]$ $\displaystyle+$
$\displaystyle\sum_{i\neq
j,k}^{N}U_{3b}(\mathbf{R}_{i},\mathbf{R}_{j},\mathbf{R}_{k})\>,$

perturbation theory. The last term in (11) corresponds to microwave-induced
three-body interactions (see B), that can not be written in terms of binary
potentials. The interactions of the third order contribution have dipolar
character, which arise from a single excitation exchange between ($s,p$)-pairs
that are virtually populated by the off-resonant microwave fields (see figure
2a). As illustrated in figure 2b the next order binary terms correspond to two
consecutive exchange processes between an atom pair and therefore results in a
van der Waals type interaction $\sim R_{ij}^{6}$ (see 11). In a similar
fashion, the three-body terms also involve two pair-exchange processes but
between three atoms, as illustrated in figure 2c. Specific choices of the
microwave parameters significantly simplify the derived expressions (10) and
(11), and permit to realize qualitatively different types of interactions,
which will be discussed in the following.

Figure 2: Schematic explanation of various exchange mechanisms between the
pairs, which is controlled by microwave fields, result in different types of
interactions. A single excitation exchange between the pairs gives rise to (a)
dipolar interactions which is manifested in the third-order perturbation
theory. Similarly, two exchange processes yield (b) van der Waals and (c)
three-body interactions in the fourth order contribution between two and three
atoms, respectively.

###  4.1 Dipole-dipole interaction

We start by considering only the two linearly polarized fields, i.e.
$\Omega_{\sigma}=\Omega_{\sigma}^{\prime}=0$. In this most simple case, the
interactions (10) and (11) greatly simplify and depend on the molecular
orientation only through the $Z$-component of the inter-particle distance. If
we further require for the microwave parameters

$\Omega^{\prime}=\frac{|\mu|\Delta^{\prime}\Omega}{|\mu^{\prime}|\Delta}\>\>,\quad\Delta^{\prime}=\frac{\mu^{\prime
2}\Delta}{\mu^{2}},$ (12)

the fourth order contribution vanishes, $E^{(4)}=0$, and the third order term
yields pure dipole-dipole interactions

$E^{(3)}=\sum_{i<j}\frac{\Omega^{2}}{\Delta^{2}}\frac{4\mu^{4}}{R_{ij}^{3}}\left(1-3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)\>.$
(13)

The only remaining free parameter that controls the strength of the microwave
induced interactions is the ratio of the Rabi frequency and the detuning,
$\Omega/\Delta$.

To obtain the total atomic interaction we add both contributions in the
laboratory frame (8) and rotating frame of reference (13). Figure 3 shows the
resulting interaction potential between two Rubidium atoms in $55s$ Rydberg
states, microwave dressed according to (12) for $|\Omega/\Delta|=0.08$ and
$\Delta=-150$ MHz. A comparison to the numerical result demonstrates that this
procedure indeed yields a proper description of the total interactions. Here
we assumed that the internuclear axis is aligned with the polarization axis of
the microwave, such that (13) is attractive. A crossover between van der Waals
interactions and dipole-dipole interactions due to the microwave appears at a
critical distance $\sim(C_{6}\Delta^{2}/(8\mu^{4}\Omega^{2}))^{1/3}$. Hence,
the combination of the attractive dipole-dipole interaction with the repulsive
van der Waals term gives rise to a pronounced potential well around this
distance, whose depth and position can be widely tuned by the microwave
parameters. Since the total potential turns isotropically repulsive for
perpendicular orientation of the internuclear axis this type of potential can
be used to create long-range bound molecules [38, 40] which can be aligned
with the polarization direction and whose properties can be well controlled by
the external microwave fields and tuned dynamically. Similar interaction
potentials have also been considered for microwave-dressed polar molecules
[21, 39] and where shown to enable superfluid p-wave pairing in ultracold
fermionic ensembles [39].

Figure 3: Perturbative potential (dashed red line) compared to exact Floquet
calculations (blue line) where the assumptions in (12) are made for $n=55$ and
$|\Omega/\Delta|=0.08$ and $\Delta=-150$ MHz. The dotted black line shows the
bare vdW interaction.

###  4.2 van der Waals interaction

Figure 4: Perturbative potential (dashed red line) compared to the bare vdW
interaction (dotted black line) and the nonperturbative Floquet calculations
(blue line). We have chosen $n=55$ and $|\Omega/\Delta|=0.15$ (a), $0.2$ (b)
for $\Delta=-150$MHz, while the remaining parameters are chosen according to
(15) to cancel the long-range interaction.

Another important case can be realized by applying both linearly and
circularly polarized fields. In order to simplify the discussion we will set
the Rabi frequencies to

$\displaystyle\Omega_{\sigma}$ $\displaystyle=$
$\displaystyle\Omega_{\pi}=\Omega$ $\displaystyle\Omega_{\sigma}^{\prime}$
$\displaystyle=$ $\displaystyle\Omega_{\pi}^{\prime}=\Omega^{\prime}\>.$ (14)

If one further requires a slightly different relation between the parameters
for the driving fields of the upper and lower transition

$\Omega^{\prime}=\frac{|\mu|\Delta^{\prime}\Omega}{|\mu^{\prime}|\Delta}\>\>,\quad\Delta^{\prime}=-\frac{\mu^{\prime
2}\Delta}{\mu^{2}},$ (15)

one can cancel the third order contribution $E^{(3)}=0$, i.e. realize
vanishing dipole-dipole interactions in the presence of microwave driving.
However, (11) remains finite and yields isotropic van der Waals type
interactions. Hence the total interaction in the rotating and laboratory frame
can be expressed in terms of an effective van der Waals interaction

$E=\sum_{i<j}\frac{C_{6}^{\prime}}{R_{ij}^{6}}\>,$ (16)

with a van der Waals coefficient

$C_{6}^{\prime}=C_{6}+24\kappa\frac{\mu^{4}}{\Delta}\>,$ (17)

that can be tuned by the external fields. Again the dressing-induced
contribution is suppressed by the ratio $\kappa=\Omega^{2}/\Delta^{2}$, which
corresponds to the small fraction of admixed $p$-states. However, the
microwave detuning, $\Delta$, which enters in the denominator of (17), can be
made much smaller than the energy mismatch between the different pair states
in the laboratory frame. The latter determines the strength of the bare van
der Waals interaction, such the magnitude of the microwave-induced van der
Waals coefficient can be on the same order than the bare coefficient. However,
the additional potential can be made attractive or repulsive depending on the
sign of the microwave detuning. In particular, one can, thus, use the
microwave fields to cancel the binary van der Waals interactions entirely by
setting

$\Delta=-24\kappa\mu^{4}/C_{6}.$ (18)

Figure 4 shows the resulting interaction potentials for microwave dressed
Rb($55s$) atoms with $|\Omega/\Delta|=0.15,0.2$ and $\Delta=-150$MHz. A
comparison to the numerical results and the bare van der Waals interaction
demonstrates that the long-range tail of the latter can indeed be suppressed
to vanishingly small values. However, below a critical distance the dipole-
dipole coupling between the Rydberg states exceeds the microwave detuning,
causing a break down of the perturbation theory. As result, the interaction
starts to increase again for distances within this critical radius (see Fig.
4b).

Figure 5: Interaction energies relative to bare vdW interaction
$V=C_{6}/R^{6}$ at different interatomic distances and for different ratios of
$|\Omega/\Delta|=0.1$ (a), $0.15$ (b), $0.2$ (c) and for $n=55$. The solid
line shows the perturbative result.

To further analyze the possibility to cancel interactions, we show in Fig.5
the $\Delta$-dependence of the numerically obtained interaction energies for
different interatomic distances relative to the bare interaction energy (8).
With a decreasing ratio $\Omega/\Delta$ the data points tend to collapse onto
a single curve that approaches the analytical prediction (16) and (17) for
large distances. Despite showing significant deviations from the perturbative
results, comparatively large ratios still yield detunings for which the
interactions are strongly suppressed down to comparably small distances of a
few micrometers. For given microwave intensities one can, hence, find an
optimal detuning for which interactions are strongly suppressed over a broad
range of distances. One such case is shown in figure 6, corresponding to
$|\Omega/\Delta|=0.2$ and $\Delta=-170$ MHz.

The found modification of the interaction potential can have significant
consequences for the Rydberg blockade effect. For example, for Rydberg
excitation with a typical Rabi frequency $\sim 200$kHz and the interaction
curves of figure 6, microwave dressing decreases the corresponding blockade
radius by about a factor two. This implies a significant change of the
excitation dynamics, which for the bare van der Waals interaction would
require a $\sim 4000$-fold increase of the intensity of the laser that drives
the Rydberg transition. Such a reduction of the blockade radius implies a
eightfold decrease of the fraction of blockaded atoms, a quantity which has
been thoroughly studied in a number of previous measurements and, thereby,
provides a sensitive and well established experimental approach to study the
effects of the microwave dressing and to identify optimal field parameters for
maximum interaction control.

Figure 6: Perturbative potential (dashed red line) compared to the bare vdW
interaction (dotted black line) and the nonperturbative Floquet calculations
(blue line) for optimized detuning $\Delta=-170$ MHz, $|\Omega/\Delta|=0.2$
and for $n=55$.

###  4.3 Three-body interaction

Three-body interactions that appear in the fourth order energy correction (see
(11)) become significant under conditions where binary interactions are
minimized. The general expression for $U_{\rm 3b}$ is rather lengthy and given
in B. However, if we choose the parameters such that the long-range tail of
the binary interactions vanishes ( conditions (4.2), (15) and (18)) the
potential takes on a particularly simple form

$U_{\rm
3b}(\mathbf{R}_{i},\mathbf{R}_{j},\mathbf{R}_{k})=\sum_{i<j<k}-\frac{6\Omega^{2}}{\Delta^{3}}\mu^{4}\left\\{\frac{1}{R_{ij}^{3}R_{jk}^{3}}\left[1-3\frac{(\mathbf{R}_{ij}\cdot\mathbf{R}_{jk})^{2}}{R_{ij}^{2}R_{jk}^{2}}\right]+i\leftrightarrow
j+j\leftrightarrow k\right\\}\>.$ (19)

In figure 7, we show the interaction potential scaled by the bare van der
Waals potential $C_{6}/r_{12}^{6}$ for three particles at positions ${\bf
r}_{i}$ ($i=1,2,3$), fixing ${\bf r}_{1}$ and ${\bf r}_{2}$ and varying the
position ${\bf r}_{3}$ of the third atom. The interaction potential exhibits a
dipolar pattern, independent of the orientation of the three particle
compound. Such genuine three-body interactions have attracted great
theoretical interest, as they give rise to a wealth of exotic phenomena in
condensed matter systems [41, 42]. While systems with dominating multi-body
interactions are typically scarce in nature, lattices of ultracold polar
molecules have been proposed for realizing such conditions [20]. Recent
experiments have demonstrated the implementation of elementary spin models via
Rydberg excitation of atomic lattices [43], such that the combination with the
present microwave control scheme may offer a viable approach to artificial
quantum magnets with dominating three-body interactions.

Figure 7: Three-body interaction scaled with $C_{6}/r_{12}^{6}$ for a three-
body compound. The positions of two particles are fixed and are shown by red
circles at $r_{1}$ and $r_{2}$ in the figure.

##  5 Conclusion

In summary, we have shown that off resonant dressing of Rydberg states by
multiple microwave fields provides a versatile approach to control Rydberg-
Rydberg atom interactions. In addition to inducing familiar dipole-dipole
interactions, the presented scheme enables to induce isotropic interactions
and effective van der Waals potentials, allowing to turn off the long-range
tail of the field-free interactions entirely, for proper parameters realizes
conditions with dominating three-body interactions. For simplicity, we have
neglected the fine structure splitting of the Rydberg states, which is well
justified for describing the van der Waals interaction between ($ns,ns$) pairs
[36, 37]. Including fine structure for the microwave dressing one obtains
identical microwave-induced potentials, upon adding two dressing fields with
circular polarisation opposite to that of the $\sigma$\- and
$\sigma^{\prime}$-fields (cf. Fig.1) in order to assure symmetric dressing for
the $m_{J}=\pm 1/2$ states of the $nS_{1/2}$ Rydberg manifold.

Such control capabilities may find direct applications in ongoing experiments
towards quantum simulations of magnetism using ultracold Rydberg ensembles.
For example, the demonstrated reduction of the Rydberg blockade radius
suggests a promising approach to significantly increase the number of
excitable atoms to high lying, and thus long lived states in finite atomic
lattices, as studied in recent experiments. The ability, to tune the long-
distances behaviour of atomic interactions may also be important for
experimental studies of nonequilibrium phase transitions in driven,
dissipative Rydberg gases and lattices [44] for which the power-law tail of
the van der Waals potential where recently shown to play an important role
[45].

Precise measurements of van der Waals Rydberg-Rydberg atom interactions have
been recently demonstrated in experiments with two or three atoms confined in
separate dipole traps at well controllable distances [46]. Such settings
appear to be ideally suited to explore the effects of microwave dressing,
optimize the control parameters and study the emergence of exotic genuine
multi-body interactions due to microwave dressing.

Another important application of Rydberg-state interactions is the realization
of large optical nonlinearities in cold atomic gases under electromagnetically
induced transparency (EIT) conditions, which has been explored in recent
theoretical and experimental studies [27, 28, 47, 48, 49, 50, 51]. It is
proposed that interactions between atoms in the Rydberg states diminish the
EIT effect and thereby leads to highly nonlinear absorption and refraction,
corresponding to dissipative and unitary effective photon interactions[47, 48,
52, 53, 54] . Microwave control of the atomic interactions can therefore be
employed to control photonic interactions in such Rydberg-EIT media as shown
in a recent experiment [55].

We acknowledge the financial support from the EU through the FET-Open Grant
MALICIA No. 265522, through the FET grant HAIRS No. 612862 and through the
Marie Curie ITN ”COHERENCE”, and from IARPA through the project grant AQuA-64.

##  Appendix A Rydberg-Rydberg atom interaction in rotating frame

The Rydberg-Rydberg atom interaction Hamiltonian in the rotating frame of
reference can be conveniently written as a sum of two terms

$\hat{H}_{\rm AA}=\sum_{i<j}\hat{A}_{ij}+\hat{B}_{ij}$ (20)

where
$\hat{A}_{ij}=\hat{A}_{x_{i}x_{j}}+\hat{A}_{y_{i}y_{j}}+\hat{A}_{z_{i}z_{j}}+\hat{A}^{\prime}_{x_{i}x_{j}}+\hat{A}^{\prime}_{y_{i}y_{j}}+\hat{A}^{\prime}_{z_{i}z_{j}}$,
$\hat{B}_{ij}=\hat{B}_{x_{i}y_{j}}+\hat{B}_{x_{i}z_{j}}+\hat{B}_{y_{i}z_{j}}+\hat{B}^{\prime}_{x_{i}y_{j}}+\hat{B}^{\prime}_{x_{i}z_{j}}+\hat{B}^{\prime}_{y_{i}z_{j}}$
and

$\displaystyle\hat{A}_{x_{i}x_{j}}$ $\displaystyle=$
$\displaystyle\frac{\mu^{2}}{2R_{ij}^{3}}\left(1-3\frac{X_{ij}^{2}}{R_{ij}^{2}}\right)\left(\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}-\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}+\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}-\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}+{\rm
h.c.}\right)$ $\displaystyle\hat{A}_{y_{i}y_{j}}$ $\displaystyle=$
$\displaystyle\frac{\mu^{2}}{2R_{ij}^{3}}\left(1-3\frac{Y_{ij}^{2}}{R_{ij}^{2}}\right)\left(\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}+\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}+\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}+\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}+{\rm
h.c.}\right)$ $\displaystyle\hat{A}_{z_{i}z_{j}}$ $\displaystyle=$
$\displaystyle\frac{\mu^{2}}{R_{ij}^{3}}\left(1-3\frac{Z_{ij}^{2}}{R_{ij}^{2}}\right)\left(\hat{\sigma}_{sp_{0}}^{(i)}\hat{\sigma}_{p_{0}s}^{(j)}+{\rm
h.c.}\right)$ (21) $\displaystyle\hat{B}_{x_{i}y_{j}}$ $\displaystyle=$
$\displaystyle
3i\mu^{2}\frac{X_{ij}Y_{ij}}{R_{ij}^{5}}\left(\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}-\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}-{\rm
h.c.}\right)$ $\displaystyle\hat{B}_{x_{i}z_{j}}$ $\displaystyle=$
$\displaystyle-3\mu^{2}\frac{X_{ij}Z_{ij}}{\sqrt{2}R_{ij}^{5}}\left(\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{0}s}^{(j)}-\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{0}s}^{(j)}+\hat{\sigma}_{sp_{0}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}-\hat{\sigma}_{sp_{0}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}+{\rm
h.c.}\right)$ $\displaystyle\hat{B}_{y_{i}z_{j}}$ $\displaystyle=$
$\displaystyle
3i\mu^{2}\frac{Y_{ij}Z_{ij}}{\sqrt{2}R_{ij}^{5}}\left(\hat{\sigma}_{sp_{-}}^{(i)}\hat{\sigma}_{p_{0}s}^{(j)}+\hat{\sigma}_{sp_{+}}^{(i)}\hat{\sigma}_{p_{0}s}^{(j)}-\hat{\sigma}_{sp_{0}}^{(i)}\hat{\sigma}_{p_{-}s}^{(j)}-\hat{\sigma}_{sp_{0}}^{(i)}\hat{\sigma}_{p_{+}s}^{(j)}-{\rm
h.c.}\right)$ (22)

The operators $\hat{A}^{\prime}$ and $\hat{B}^{\prime}$ are obtained by
replacing $\mu\rightarrow\mu^{\prime}$ and $p_{0,\pm}\rightarrow
p_{0,\pm}^{\prime}$.

##  Appendix B Three-Body Interaction terms

Similarly we write the three-body interaction in (11) as $U_{\rm 3b}({\bf
R}_{i},{\bf R}_{j},{\bf R}_{k})=E_{ijk}+E_{kij}+E_{jki}$, where

$E_{ijk}=G_{ijk}-G^{\prime}_{ijk},$ (23)

and

$\displaystyle G_{ijk}$ $\displaystyle=$ $\displaystyle\frac{\mu^{4}}{\Delta
R_{ij}^{3}R_{jk}^{3}}\left[\left(\frac{\Omega_{\pi}}{\Delta}\right)^{2}(\alpha_{ij}\alpha_{jk}+\gamma_{ij}\gamma_{jk}^{*}+\gamma_{ij}^{*}\gamma_{jk})+2\left(\frac{\Omega_{\sigma}}{\Delta}\right)^{2}(\epsilon_{ij}\epsilon_{jk}+\beta_{ij}\beta_{jk}^{*}+\gamma_{ij}^{*}\gamma_{jk})\right.$
(24) $\displaystyle+$
$\displaystyle\left.\sqrt{2}\frac{\Omega_{\pi}\Omega_{\sigma}}{\Delta^{2}}(\alpha_{ij}\gamma_{jk}+\alpha_{ij}\gamma_{jk}^{*}+\epsilon_{ij}\gamma_{jk}+\epsilon_{ij}\gamma_{jk}^{*}-\beta_{ij}\gamma_{jk}^{*}-\beta_{ij}^{*}\gamma_{jk})\right],$

with

$\displaystyle\alpha_{ij}$ $\displaystyle=$
$\displaystyle\left(1-\frac{3Z_{ij}^{2}}{R_{ij}^{2}}\right),$
$\displaystyle\beta_{ij}$ $\displaystyle=$
$\displaystyle\frac{3}{2R_{ij}^{2}}(X_{ij}^{2}-Y_{ij}^{2}-2iX_{ij}Y_{ij}),$
$\displaystyle\gamma_{ij}$ $\displaystyle=$
$\displaystyle-\frac{3}{2R_{ij}^{2}}(X_{ij}Z_{ij}-iY_{ij}Z_{ij}),$
$\displaystyle\epsilon_{ij}$ $\displaystyle=$
$\displaystyle\left(1-\frac{3(X_{ij}^{2}+Y_{ij}^{2})}{2R_{ij}^{2}}\right).$
(25)

The contribution from the lower transition $G^{\prime}_{ijk}$ is obtained by
replacing $\mu\rightarrow\mu^{\prime}$,
$\Omega_{\pi}\rightarrow\Omega^{\prime}_{\pi}$ and
$\Omega_{\sigma}\rightarrow\Omega^{\prime}_{\sigma}$. $E_{kij}$ and $E_{jki}$
are written in the same form.

## References

## References

  * [1] Lukin M D, Fleischhauer M, Cote R, Duan L M, Jaksch D, Cirac J I and Zoller P 2001 Phys. Rev. Lett.  87 037901
  * [2] Jaksch D, Cirac J I, Zoller P, Rolston S L, Cote R and Lukin M D 2000 Phys. Rev. Lett. 85 2208
  * [3] Saffman M and Walker T G 2002 Phys. Rev. A 66 065403
  * [4] Dudin Y O and Kuzmich A 2012 Science 336 887
  * [5] Weimer H, Müller M, Lesanovsky I, Zoller P and Büchler H P 2010 Nat. Phys. 6 382
  * [6] Pohl T, Demler E and Lukin M D 2010 Phys. Rev. Lett. 104 043002
  * [7] Bouchoule I and Mølmer K 2002 Phys. Rev. A 65 041803(R)
  * [8] Hofmann C S, Günter G, Schempp H, Robert-de-Saint-Vincent M, G rttner M, Evers J, Whitlock S and Weidemüller M 2013 Phys. Rev. Lett. 110 203601
  * [9] Peyronel T, Firstenberg O, Liang Q -Y, Hofferberth S, Gorshkov A V, Pohl T, Lukin M D and Vuletić V 2012 Nature 488 57
  * [10] Gallagher T F 1994 Rydberg atoms (Cambridge Univ. Press)
  * [11] Li W, Tanner P J and Gallagher T F 2005 Phys. Rev. Lett. 94 173001
  * [12] Afrousheh K, Bohlouli-Zanjani P, Vagale D, Mugford A, Fedorov M and Martin J D D 2004 Phys. Rev. Lett. 93 233001
  * [13] Mudrich M, Zahzam N, Vogt T, Comparat D and Pillet P 2005 Phys. Rev. Lett. 95 233002
  * [14] Carroll T J, Claringbould K, Goodsell A, Lim M J and Noel M W 2004 Phys. Rev. Lett. 93 153001
  * [15] Pohl T, Sadeghpour H R, Schmelcher P 2009 Phys. Rep. 484 181
  * [16] Bohlouli-Zanjani P, Petrus J A and Martin J D D 2007 Phys. Rev. Lett. 98 203005
  * [17] Petrus J A, Bohlouli-Zanjani P and Martin J D D 2008 J. Phys. B 41 245001
  * [18] Han J and Gallagher T F 2009 Phys. Rev. A 79 053409
  * [19] Park H, Tanner P J, Claessens B J, Shuman E S and Gallagher T F 2011 Phys. Rev. A 84 022704
  * [20] Buechler H P, Micheli A and Zoller P 2007 Nat. Phys. 2007 3 726
  * [21] Buechler H P, Demler E, Lukin M, Micheli A, Prokof’ev N, Pupillo G and Zoller P 2007 Phys. Rev. Lett. 98, 060404
  * [22] Gorshkov A V, Rabl P, Pupillo G, Micheli A, Zoller P, Lukin M D and Büchler H P, 2008 Phys. Rev. Lett. 101 ,073201
  * [23] Müller M, Liang L, Lesanovsky I and Zoller P 2008 N. J. Phys. 10 093009
  * [24] Tanasittikosol M, Pritchard J D, Maxwell D, Gauguet A, Weatherill K J, Potvliege R M and Adams C S 2011 J. Phys. B 44 184020
  * [25] Mohapatra A K, Jackson T R and Adams C S 2007 Phys. Rev. Lett. 98 113003
  * [26] Bason M G, Mohapatra A K, Weatherill KJ and Adams C S 2008 Phys. Rev. A 77 032305
  * [27] Maxwell D, Szwer D J, Paredes-Barato D, Busche H, Pritchard J D, Gauguet A, Weatherill K J, Jones M P A and Adams C S 2013 Phys. Rev. Lett. 110 103001
  * [28] Paredes-Barato D and Adams C S 2014 Phys. Rev. Lett. 112 040501
  * [29] Jones L A, Carter J D and Martin J D D, 2013 Phys. Rev. A 87 023423
  * [30] Colombo A P, Zhou Y, Prozument K, Coy S L and Field R W 2013 J. Chem. Phys. 138 014301
  * [31] Saffman M, Walker T G and Mølmer K 2010 Rev. Mod. Phys. 82 2313
  * [32] Comparat D and Pillet P 2010 J. Opt. Soc. Am. 27 A208
  * [33] Zimmerman M L, Littman M G, Kash M M and Kleppner D 1979 Phys. Rev. A 20 2251
  * [34] Li W, Mourachko I, Noel M W and Gallagher T F 2003 Phys. Rev. A 67 052502
  * [35] Chu S-I and Telnov D A 2004 Phys. Rep. 390 1
  * [36] Singer K, Stanojevic J, Weidemüller M and Cote R 2005 J. Phys. B 38 S295
  * [37] Walker T G and Saffman M 2008 Phys. Rev. A 77 032723
  * [38] Boisseau C, Simbotin I and Côté R 2002 Phys. Rev. Lett. 88 133004
  * [39] Cooper N R and Shlyapnikov G V 2009 Phys. Rev. Lett. 103 155302
  * [40] Kiffner M, Park H, Li W and Gallagher T F 2012 Phys. Rev. A 86 031401(R)
  * [41] Moessner R and Sondhi S L 2001 Phys. Rev. Lett. 86 1881
  * [42] Motrunich O I and Senthil T 2002 Phys. Rev. Lett. 89 277004
  * [43] Bettelli S, Maxwell D, Fernholz T, Adams C S, Lesanovsky I and Ates C 2013 Phys. Rev. A 88 043436
  * [44] Schauss P, Cheneau M, Endres M, Fukuhara T, Hild S, Omran A, Pohl T, Gross C, Kuhr S and Bloch I 2012 Nature 491 87
  * [45] Zoubi H, Eisfeld A and Wüster S 2014 Phys. Rev. A 89 053426
  * [46] Béguin L, Vernier A, Chicireanu R, Lahaye T and Browaeys A 2013 Phys. Rev. Lett. 110 263201
  * [47] Sevinçli S, Henkel N, Ates C and Pohl T 2011 Phys. Rev. Lett. 107 153001
  * [48] Petrosyan D, Otterbach J and Fleischhauer M 2011 Phys. Rev. Lett. 107 213601
  * [49] Mohapatra A K, Bason M G, Butscher B, Weatherill K J and Adams C S 2008 Nat. Phys. 4 890
  * [50] Pritchard J D, Maxwell D, Gauguet A, Weatherill K J, Jones M P A and Adams C S 2010 Phys. Rev. Lett. 105 193603
  * [51] Ates C, Sevinçli S and Pohl T 2011 Phys. Rev. A 83 041802(R)
  * [52] Gorshkov A V, Otterbach J, Fleischhauer M, Pohl T and Lukin M D 2011 Phys. Rev. Lett. 107 133602
  * [53] Parigi V, Bimbard E, Stanojevic J, Hilliard A J, Nogrette F, Tualle-Brouri R, Ourjoumtsev A and Grangier P 2012 Phys. Rev. Lett. 109 233602
  * [54] Gorshkov A V, Nath R and Pohl T 2013 Phys. Rev. Lett. 110 153601
  * [55] Maxwell D, Szwer D J, Paredes-Barato D, Busche H, Pritchard J D, Gauguet A, Jones M P A and Adams C S 2014 Phys. Rev. A 89 043827