The case of magnetic Ln atoms is more intricate than that of Cr because of
their electronic structure is not spherically symmetric. This results in an
anisotropy of the van der Waals interaction, which changes the short-range
physics Kotochigova and Petrov (2011); Petrov  _et al._ (2012); Kotochigova
(2014). Because both short- and long-range interaction potentials are
anisotropic in Lns, they both induce coupling between molecular channels,
which leads to FRs. Both anisotropic potentials were predicted to
substantially contribute to the character, distribution, and prevalence of FRs
in collisions among bosonic atoms Petrov  _et al._ (2012). A large number of
scattering channels contribute to each resonance, making it hard for coupled-
channel calculations based on partial-wave decomposition to converge. Indeed,
in Ref. Petrov  _et al._ (2012), $\ell$ up to 10 were considered, and later,
this was extended up to 20 Frisch  _et al._ (2014). An analytical model
estimates that $\ell$’s up to $\gtrsim 40$ must be considered to reproduce
experimental observations Frisch  _et al._ (2014). Needless to say, no
perturbative treatment can be safely applied.

First experimental observations of FRs in magnetic Lns were performed on 168Er
and revealed an unusually large number of FRs within a narrow magnetic field
range Aikawa  _et al._ (2012). Soon after, theory work of Ref. Petrov  _et
al._ (2012) predicted that the high Feshbach spectral density is a general
feature of magnetic Ln atoms. Systematic high-resolution trap-loss
spectroscopy later extended the observations of dense FRs in two bosonic Er
isotopes from 0 to 70 G Frisch  _et al._ (2014); see Fig. 3a. A statistical
analysis of the Feshbach spectra revealed correlations between the resonance
locations, with very similar characteristics for the two bosonic isotopes.
Based on the formalism of random matrix theory, these correlations were
quantified and found to be consistent with chaotic behaviour in the scattering
of these atoms; see Fig. 3b. The fermionic case was also studied and shows a
ten-fold larger density of FRs, increasing from 2.7 resonances per Gauss for
the bosons to 25.7 resonances per Gauss for the fermions. This was attributed
to the additional role of the hyperfine structure. A similarly high density of
FRs was concurrently reported for four Dy isotopes, three bosonic and one
fermionic, up to 6 G Baumann  _et al._ (2014). These measurements also
revealed a much higher density for the fermionic isotopes as well as an
intriguing temperature dependence of the Feshbach spectrum.

The chaotic behaviour of both Er and Dy was thoroughly analysed and compared
in the collaborative work of Ref. Maier  _et al._ (2015b). Interestingly, even
though Dy shows a higher density of FRs (4.3 resonances per Gauss), the degree
of correlation in the resonance locations is similar for the two species, and
in fact, slightly larger for Er. This work also introduced a new scheme for
coupled-channel calculations, employing a basis comprised of $B=0\,$G
Hamiltonian eigenstates. Such a basis is more amenable to non-central
scattering potentials, allowing more rapid computational convergence. This
theoretical analysis established the leading role played by the van der Waals
anisotropy in the chaotic scattering behaviour: the DDI alone cannot account
for the correlations. The slightly larger degree of correlation observed in Er
FRs may be attributed to a larger short-range anisotropy.

The temperature dependence was further analysed in Ref. Maier  _et al._
(2015b). At a statistical level, spectra at higher temperature exhibit a
larger density of resonances with the degree of correlation that is nearly
unchanged. The temperature dependence of individual resonances reveals
behaviour compatible with $d$-wave entrance channel predictions, similar to
the case in Cr Beaufils  _et al._ (2009); see Fig. 2b. The weakening and
disappearance of higher partial-wave entrance channels in the ultracold regime
were proposed as an explanation for the variation of the FR density with
temperature. A subsequent work using Ho atoms also revealed a dense FR
spectrum as well as a similar temperature dependence Khlebnikov  _et al._
(2019). In contrast to the Er system Maier  _et al._ (2015b), the chaotic
character of the Ho FR spectrum was observed to increase with $T$, ranging
from intermediate at $T~{}2\mu$K to almost fully chaotic at $T~{}12\mu$K. In
the high-temperature regime, the shifts of Ho $d$-wave resonances deviate from
a linear dependence on $T$, which is expected at low temperatures, as observed
in Er Maier  _et al._ (2015b).The authors speculate that this behaviour is
responsible for the observed change in FR statistics.

Because of the complexity and non-perturbative character of the coupling
between channels, the FRs cannot in general be assigned to a particular bound
state. A few assignments were, however, assigned in the low magnetic field
region of Er Frisch  _et al._ (2015). This was accomplished by comparing
measurements of the molecular state binding energy with coupled channel
calculations.

Broader FR resonance can be observed on top of the forest of narrow FR
features in the Ln atoms. Such features were already observable in the early
FR scans of Ref. Frisch  _et al._ (2014), where, e.g., a 3.5 G-wide FR is
found at 57 G for 168Er; see Fig. 3a. A detailed study of such features was
performed in 164Dy Maier  _et al._ (2015a) up to 600 G and revealed several
more broad features. Two of these broad features were extensively
characterised via spectroscopy of the molecular state binding energy.
Universal properties were found that are characteristic of open-channel-
dominated $s$-wave FRs. Particularly noticeable is the decoupling of the broad
FR from the chaotic background of the narrow FRs, allowing one to observe the
universal behaviour of the halo dimer binding energies over several Gauss to
the low-$B$ side of the FRs; see Fig. 4. A similar study was later reported in
162Dy Lucioni  _et al._ (2018): two broad FRs with $s$-wave character are
found at magnetic fields in the 20–30 G range. Broad FRs are also observed in
Ho Khlebnikov  _et al._ (2019).

##### Scattering length measurements

As is the case for other multi-electron atoms Pethick and Smith (2002), ab
initio calculations fail to predict the scattering lengths of Ln atoms.
Instead, these must be experimentally measured. A common technique involves
the observation of cross-rethermalisation: one observes the rethermalisation
dynamics of a thermal sample after exciting it along one of its the three
spatial directions Monroe  _et al._ (1993); the thermalization rate is
proportional to the elastic cross section, which in turn is proportional to
the square of the scattering length. The proportionality relies on a theory
prediction for the average number of collisions necessary for thermalisation
and depends on the particular partial waves and interaction involved in the
collision DeMarco  _et al._ (1999). Bohn et al.Bohn and Jin (2014) calculated
this for the anisotropic DDI, which was verified experimentally in Ref. Aikawa
_et al._ (2014b) using spin-polarised fermionic Er. This theory allowed Tang
et al. Tang  _et al._ (2015b, 2016a) to extract scattering lengths from cross-
rethermalisation measurements of two Dy isotopes Tang  _et al._ (2015b,
2016a). The theory-experiment comparison resulted in $a=122(10)\,a_{0}$ for
162Dy and $a=92(8)\,a_{0}$ for 164Dy at $B=1.58\,$G. This field is at least
$0.1\,$G away from any FR in either species. These $a$-values, as well as
their $B$-dependence, were also measured in Ref. Tang  _et al._ (2016b) using
a different method involving a detailed analysis of the anisotropic expansion
of thermal gases after release from their confining potentials., The analysis
relied on the competition of the DDI and the contact interaction; see Sec.
IV.1.2 for details. The results are consistent with the cross-rethermalisation
experiments: $a=154(22)\,a_{0}$ for 162Dy and $a=96(22)\,a_{0}$ for 164Dy, but
with larger error. A similar measurement near a 5 G FR in 162Dy yielded a
background scattering length of $a=157(24)\,a_{0}$, which is consistent with
the lower-field measurement, but larger than that obtained with cross-
rethermalization Tang  _et al._ (2015b, 2016a). A third scheme relies on the
universal scaling of the molecular binding energy close to broad FRs to
extract scattering lengths in Dy. Such an analysis used a relatively high-
field resonance in 164Dy to obtain $a_{\rm{bg}}=91(16)\,a_{0}$ Maier  _et al._
(2015a), which is in good agreement with the low-$B$ values of Refs. Tang  _et
al._ (2015b, 2016b). A similar study in 162Dy yielded
$a_{\rm{bg}}=220(50)\,a_{0}$ Lucioni  _et al._ (2018). This is consistent with
the $a=154(22)\,a_{0}$ measurement in Ref.Tang  _et al._ (2016b), but larger
than both the 5-G measurement in Ref.Tang  _et al._ (2016b) and the
$a=122(10)\,a_{0}$ measurement of Ref. Tang  _et al._ (2015b).

For Er, preliminary cross-rethermalisation measurements were first reported in
Ref. Frisch (2014) for the four bosonic isotopes in the low-field regime. A
simpler analysis of the rethermalisation rates than used in Ref. Tang  _et
al._ (2015b) yielded $a_{\rm{bg}}=81(10)\,a_{0}$ for 164Er,
$a_{\rm{bg}}=72(13)\,a_{0}$ for 166Er, $a_{\rm{bg}}=200(23)\,a_{0}$ for 168Er,
and $a_{\rm{bg}}=221(22)\,a_{0}$ for 170Er. In 170Er, the background
scattering length was first thought to be negative, but more recent
measurements have questioned this Trautmann  _et al._ (2018). These Er
$a_{\rm{bg}}$ measurements were taken by averaging data obtained at different
$B$-values between 0.2 G and 1 G, away from FRs. A more precise estimate was
later developed to extract $a$ and its $B$ dependence by loading a quantum gas
into a deep 3D lattice and performing lattice-modulation spectroscopy Baier
_et al._ (2016). This probes the Mott-insulator excitation gap, which is
related to the scattering length through the on-site interaction energy; see
also Sec. VII.2. This technique yielded $a=137(1)a_{0}$ at 0.4 G for 168Er
Baier  _et al._ (2016), and for 166Er, it was used to make a sequence of
measurements of $a$ versus $B$ from 0-to-2.5 G Chomaz  _et al._ (2016). In
this range, the data are well-described by overlapping FRs with a
$B$-dependent background scattering length $a_{\rm{bg}}(B)=62(4)+kB$, with
$k=5.8(1.2)\,a_{0}/{\rm G}$.

In the case of Ho, recent experiments reported scattering length values using
cross-rethermalisation. The first such measurement yielded a background value
of $a_{\rm{bg}}=144(38)\,a_{0}$ Khlebnikov  _et al._ (2019), while a
subsequent experiment reported $a_{\rm{bg}}=90(11)\,a_{0}$ via a more
extensive set of measurements Davletov  _et al._ (2020).

We note that a debate regarding the exact values of the scattering lengths in
bosonic Er and Dy in the quantum degenerate regime has arisen in light of
measurements of these parameters extracted from the many-body physics of the
droplet and supersolid states; see also Sec. V. These differ from the two-body
collision-based methods described above. Central to the debate is the
appropriateness of the description of the many-body physics based on
perturbative mean-field and local-density approximations, on the one hand, and
on the other, regarding the relevance, in the regime of interest, of
corrections beyond the Born approximation and the role of the momentum
dependence in the two-body scattering description. That is, the many-body
theory used to extract scattering lengths in such finite-size quantum systems
could be as problematic as the use of ultracold thermal gases in the former
collision-based measurements. Those questions are beyond the scope of the
current section, but additional discussions may be found in Secs. V.2 and
V.3.4.

Concerning the scattering length itself, the most extended discussions focus
on the 164Dy case. Studies of single quantum droplets provide estimates of the
background scattering length via a theory-experiment comparison of different
many-body properties, namely the critical atom number for the existence of a
self-bound droplet state and the frequencies of droplet’s elementary
excitations Schmitt  _et al._ (2016); Ferrier-Barbut  _et al._ (2018a);
Böttcher  _et al._ (2019a). These investigations yield values of
$a_{\rm{bg}}=62.5\,a_{0}$ and $a_{\rm{bg}}=69(4)\,a_{0}$ at $B\approx 6.6\,$G
and $B=0.8\,$G, respectively. While they are mutually compatible, they are
consistently lower than the measurement results in Refs. Tang  _et al._
(2015b); Maier  _et al._ (2015a); Tang  _et al._ (2016b) obtained with thermal
gases. Ferrier-Barbut et al. Ferrier-Barbut  _et al._ (2018a) speculate that
the discrepancy could come from the different temperatures of the samples used
in those two sets of experiments and related momentum dependence of the
scattering 777We note that a similar down-shift of the scattering length can
be observed in Er based on lattice modulation spectroscopy of a Mott insulator
versus that from cross-rethermalisation of a thermal sample. This might also
be interpreted as arising from the momentum dependence of the scattering..
Reference Böttcher  _et al._ (2019a) extends the measurements of the critical
atom number of a self-bound quantum droplet (see Sec. V.3.5) and was able to
extract the ratio between the background scattering lengths of the two
dysprosium isotopes to be $a_{\rm bg,162}/a_{\rm bg,164}=2.03(6)$. Using
$a_{\rm bg}=69(4)a_{0}$ for 164Dy Ferrier-Barbut  _et al._ (2018a) yields a
value of $a_{\rm bg}=140(7)a_{0}$ for the background scattering length of
162Dy. This agrees with most of the values obtained from thermal samples
within their uncertainty Tang  _et al._ (2015b); Lucioni  _et al._ (2018);
Tang  _et al._ (2016b). In the Er case, we note that sensitive probes based on
collective excitations of single macro-droplets Chomaz  _et al._ (2016),
supersolids Natale  _et al._ (2019), and roton excitations of dBECs Chomaz
_et al._ (2018); Petter  _et al._ (2019) show a good agreement between
experiment and theory using the lattice-calibrated scattering lengths values.

Finally, we note that little information is available on Ln spin-dependent
scattering, as opposed to the well-characterised Cr system. One exception is
the case of the two lowest energy weak-field seeking spin states of fermionic
167Er, whose scattering have been experimentally investigated in Ref. Baier
_et al._ (2018); see also Sec. VI. Extensive loss spectroscopy measurements
were performed in a limited $B$-field range (up to $2$ G) and for varying
mixture compositions. This enabled the identification of both intra-spin FRs
and inter-spin FRs; spin-polarised loss spectra have also been observed in
161Dy for its lowest energy weak-field seeking spin states Burdick  _et al._
(2016b). Lattice modulation spectroscopy also characterised the $B$-to-$a$
conversion close to a relatively broad inter-spin FR. The background
scattering for the interaction of the two lowest spin-states is
$a_{\rm{bg}}=91(8)\,a_{0}$. Theoretical investigations have also investigated
spin-dependant resonances Li  _et al._ (2018).

An important difference between Lns and Cr lies in the origin of spin-
dependent scattering. In contrast to Cr, the dominant effect in Lns arises
from their orbital anisotropy (i.e., $L\neq 0$; see also Sec. II.2). This
induces an interaction term that couples the angular momentum of each
individual atom to the orbital momentum of the pair and results in sets of
molecular potentials with very different character versus Cr atoms. More
systematic studies, especially experimental, are needed to fully elucidate the
spin-dependent collisional physics. A relevant open question relates to the
density of the Feshbach spectra in Ln spin mixtures and in fermionic Ln
isotopes. Regarding the latter, the density of FR spectra at high field in
fermionic Dy seem to be spin-dependant Burdick  _et al._ (2016b). Ericson
fluctuations, known from resonances in nuclear physics, may provide an
appropriate framework for explicating the field dependence of the overlapping
loss spectra Mayle  _et al._ (2012) observed in Ref. Burdick  _et al._
(2016b). Interspecies scattering involving Ln atoms is also of recent interest
Trautmann  _et al._ (2018); Ravensbergen  _et al._ (2018b, 2020); Zaremba-
Kopczyk  _et al._ (2018).

####  II.4.3 Feshbach molecules of magnetic atoms

Besides tuning of $a$, FRs provide the ability to associate a pairs of free
atoms into a loosely bound molecule Donley  _et al._ (2002); Chin  _et al._
(2010); Köhler  _et al._ (2006); Ferlaino  _et al._ (2009); Moses  _et al._
(2017). A FR couples two states, the state of the two free particles
scattering in the entrance channel and the bound state of the closed channel.
The relative energy of the two states is tunable via the $B$-field thanks to
their different magnetic moments. The coupling of the two states yields an
avoided crossing at $B=B_{0}$ where their bare energies coincide. By ramping
$B$ across the resonance (from $B>B_{0}$ to $B<B_{0}$), one can then drive the
population of the molecular state by (adiabatically) following the low-energy
state Donley  _et al._ (2002); Chin  _et al._ (2003); Regal and Jin (2003);
Cubizolles  _et al._ (2003); Herbig  _et al._ (2003); Hodby  _et al._ (2005).
Alternative approaches use small-amplitude modulation of $B$ at a frequency
resonant with the bound-state energy for $B\lesssim B_{0}$ Hanna  _et al._
(2007); Thompson  _et al._ (2005), or even further away from the FR via
radiofrequency coupling Weber  _et al._ (2008); Klempt  _et al._ (2008). Such
protocols are commonly used in alkali metal experiments, with typical
conversion efficiency of few tens–of–percents in the bosonic case and close to
unity in the fermionic case. Note that these molecules are not in their
internal ground-state, but in a highly excited (rovibrationnal) state. Close
to the FR, the properties of the weakly-bound molecule, referred to as a
_halo-dimer_ , follow a universal behaviour.

By producing Feshbach molecules of magnetic atoms, one creates ultracold gases
of dipolar particles with even larger permanent dipole moments. Indeed,
compared to the bare atoms, the mass $m$ of the molecule is doubled while the
magnetic moment $\mu$ is also nearly doubled, resulting in a dipolar length
$a_{\rm dd}\propto m\mu^{2}$ eight times longer.

An experiment using Er has demonstrated and studied the production of an
ultracold sample of Er2 Feshbach molecules using $B$-field ramps Frisch  _et
al._ (2015). Pure molecular samples were produced by the resonant removal of
the remaining free atoms. The gas contained about $2\times 10^{4}$ Er2 at 300
nK with typical densities of $\sim 8\times 10^{11}$ cm-3 and lifetimes of
$\sim 10\,$ms. Ultracold samples of molecules in four different molecular
states, associated with four different low-field FRs of 168Er, were produced.
The properties of the Feshbach molecules have been further measured, in
particular the magnetic moment, which were found to be $\mu/\mu_{\rm{B}}>11$
for three of the four states. By contrast, the bare atomic value is
$\mu/\mu_{\rm{B}}\approx 7$. The impact of the confinement geometry and
relative dipole orientation on the scattering properties of the Feshbach
molecules has also been investigated. In particular, a reduction of the
relaxation rate in quasi-two dimensional geometry with out-of-plane dipole
orientation was observed.

In the Cr case, the study of Feshbach molecules has focused on a single
peculiar resonance with an $d$-wave entrance channel Beaufils  _et al._ (2009,
2010). However, no ultracold sample has yet been produced. In this case, the
molecular state is very short lived and studies have focused on the remarkable
atom-loss behaviour associated with the special character of the resonance.

Figure 4: Atom-loss spectrum (top) and measured molecular binding energy
(bottom) of 164Dy near a broad FR. The atom loss spectrum exhibits a very high
density of narrow resonances and one broad feature. The binding energy of the
weakly bound dimer displays a universal quadratic behaviour (dashed line). The
very weak variation of the binding energy as a function of magnetic field
shows that the dimer has a magnetic moment very close to that of two free
atoms $\mu\approx 20\mu_{\rm{B}}$. Data originally published in Maier  _et
al._ (2015a).

In the Dy case, studies have focused on the properties of the halo dimers
close to the broad FRs found at relatively large $B$ in two bosonic isotopes
Maier  _et al._ (2015a); Lucioni  _et al._ (2018), see also Sec. II.4.2. The
corresponding magnetic moments range from $\mu\approx 18\mu_{\rm{B}}$ to
$\mu\approx 20\mu_{\rm{B}}$; see Fig. 4. A striking observation is that these
halo dimers decouple from the many bound states coming from the overlaid
forest of narrow FRs. This makes the broad FRs very promising for the
production of ultracold gases of highly magnetic molecules. They have not yet
been investigated, but this constitutes an interesting prospect for future
work.

Table 1: Summary of the main atomic properties, scattering properties, and available quantum gases of highly magnetic atoms made of Cr, Er or Dy. The electronic ground state is denoted “g.s." The exchange statistic—fermionic (f) or bosonic (b)—is listed under “stat." The natural abundance is under “abound. denotes.: $I$ is the nuclear spin quantum number. $\varepsilon_{\rm{dd}}$ is the background value $a_{\rm dd}/a_{\rm{bg}}$. “Degenerate gas" lists the properties of the most populous gas of that isotope produced, as reported in journal articles. For all these, the BEC’s are almost pure; i.e., the thermal fraction is small.  atom | g.s. | $\mu_{m}(\mu_{\rm{B}})$ | stat.:isotope | abund. | $I$ | $a_{\rm dd}(a_{0})$ |  $a_{\rm{bg}}(a_{0})$ |  $\varepsilon_{\rm{dd}}$ |  degenerate gas
---|---|---|---|---|---|---|---|---|---
Cr | ${}^{7}S_{3}$ | 6 | b:52Cr | 83.79% | 0 | 15.3 |  102.5(4) Pasquiou  _et al._ (2010) |  0.150(6) |  BEC, $5\times 10^{4}$ Griesmaier  _et al._ (2005); Beaufils  _et al._ (2008)
|  |  | f:53Cr | 9.5% | $3/2$ | 15.6 |  – |  – |  DFG, $10^{3}$, $\frac{T}{T_{\rm F}}\sim 0.6$ Naylor  _et al._ (2015a)
|  |  | b:50Cr | 4.35% | 0 | 14.7 |  40(15) Schmidt  _et al._ (2003b) |  0.37(14) |  none
|  |  | b:54Cr | 2.36% | 0 | 15.9 |  unknown |  |  none
Er | ${}^{3}H_{6}$ | 7 | b:168Er | 26.8% | 0 | 66.3 |  137(1) Baier  _et al._ (2016) |  0.484(3) |  BEC, $1.5\times 10^{5}$ Frisch (2014); Aikawa  _et al._ (2012); Phelps  _et al._ (2020)
|  |  | b:166Er | 33.6% | 0 | 65.5 |  68(5) at B=1G Chomaz  _et al._ (2016) |  0.96(7) |  BEC, $10^{5}$ Chomaz  _et al._ (2016)
|  |  | f:167Er | 23% | 7/2 | 65.9 |  – |  – |  DFG: $4\times 10^{4}$, $\frac{T}{T_{\rm F}}=0.1$ Frisch (2014); Aikawa  _et al._ (2014a); Phelps  _et al._ (2020)
|  |  | b:170Er | 15% | 0 | 67 |  preliminary estimate 221(22) Frisch (2014) |  0.30(3) |  BEC, $~{}10^{4}$ Trautmann  _et al._ (2018)
|  |  | b:164Er | 1.6% | 0 | 64.7 |  preliminary estimate 81(10) Frisch (2014) |  0.8(1) |  none
|  |  | b:162Er | 0.14% | 0 | 63.9 |  unknown |  |  none
Dy | ${}^{5}I_{8}$ | 10 | b:164Dy | 28.3% | 0 | 130.7 |  disagreeing estimates: 69(4) Böttcher  _et al._ (2019a); Ferrier-Barbut  _et al._ (2018a) or 92(8) Tang  _et al._ (2015b); Maier  _et al._ (2015a); Tang  _et al._ (2016b) |  1.9(1) or 1.4(1) |  BEC, $1.5\times 10^{4}$ Lu  _et al._ (2011a); Maier (2016); Kadau  _et al._ (2016), $3.5\times 10^{4}$ Chomaz  _et al._ (2019), with typically sizeable thermal fraction remaining.
|  |  | b:162Dy | 25.5% | 0 | 129.2 |  most precise estimates are 122(10) Tang  _et al._ (2015b)–157(4) Tang  _et al._ (2016b) |  0.92(5)-1.06(9) |  BEC, $10^{5}$ Tang  _et al._ (2015a),$5\times 10^{4}$Tanzi  _et al._ (2019a); Böttcher  _et al._ (2019b)
|  |  | f:163Dy | 24.9% | 5/2 | 130 |  – |  – |  none
|  |  | f:161Dy | 18.9% | 5/2 | 128.4 |  – |  – |  DFG: $4\times 10^{4}$, $\frac{T}{T_{\rm F}}=0.1$ Burdick  _et al._ (2016b), $8\times 10^{4}$, $\frac{T}{T_{\rm F}}=0.085(10)$Ravensbergen  _et al._ (2018b)
|  |  | b:160Dy | 2.3% | 0 | 127.6 |  unknown |  |  BEC, $10^{3}$ Tang  _et al._ (2015a)

###  II.5 Comparison between the different systems

Table 1 summarises the electronic, isotopic, scattering-related properties as
well as the characteristics of all the available degenerate quantum gases for
all isotopes, both bosonic and fermionic, of the highly magnetic elements, Cr,
Er and Dy. The background $\varepsilon_{\rm{dd}}$ varies among the elements
and isotopes; e.g., it is $\varepsilon_{\rm{dd}}\sim 0.16$ for Cr, varies
between $0.9\lesssim\varepsilon_{\rm{dd}}\lesssim 1.9$ for Dy isotopes, and is
between $0.3\lesssim\varepsilon_{\rm{dd}}\lesssim 1$ for Er isotopes. All
isotopes exhibit FRs, which provide opportunities for tuning
$\varepsilon_{\rm{dd}}$; see Secs. II.4, IV, and V.

With respect to degenerate gases, all three species offer relatively large
BECs of roughly similar sizes. The maximum population of the BEC depends on
the isotope, not only because of its natural abundance, but also because of
its scattering properties. In particular, isotopes with a background value of
$\varepsilon_{\rm{dd}}\gtrsim 1$, or $\varepsilon_{\rm{dd}}<0$ (i.e., with a
negative scattering length), are more challenging to condense. They typically
require a tuning of $a$ during the evaporation process by setting $B$ close to
a FR to stabilise the quantum gas at the mean-field level. This is in
particular the case for the most abundant isotope of Dy, 164Dy, for which the
BEC numbers are otherwise smaller than for the 162Dy isotope and Er’s most
populous BEC. In fact, we point out that, due to the extreme density of FRs in
the spectra of open-shell Lns, the choice of the $B$ value at which the
evaporative cooling scheme is performed is a crucial parameter to adjust for
all isotopes of Er and Dy, and even more so for their fermionic isotopes. The
population of the achieved quantum gases greatly depend on it. Large and cold
degenerate Fermi gases of Dy and Er are also possible due to the ability to
perform direct evaporative cooling. The quantum degeneracy of fermionic Cr has
been achieved, but so far, is smaller and not as cold with respect to the
Fermi temperature.

##  III Ultracold dipolar scattering

The notion of scattering in the ultracold regime has been briefly introduced
in Secs. I and II. We will now extend our overview on ultracold collisional
physics before describing the special features of the scattering induced by
the DDI. Topics include the universal properties of elastic dipolar scattering
(i.e., coming from spin-conserving collisions; see Sec. III.2); the special
features of inelastic dipolar scattering (i.e., dipolar relaxation), including
its local character and suppression via exploiting quantum statistics and
confinement (Sec. III.3); and finally the anisotropic character of (elastic)
dipolar scattering (Sec. III.4).

###  III.1 Relevant aspects of scattering theory

Collisions among atoms in dilute ultracold gases are described within quantum
scattering theory, which accounts for the long-distance behaviour of the
wavefunction encompassing the relative motion of two colliding atoms, the so-
called _collisional wavefunction_ Landau and Lifshitz (1965); Pethick and
Smith (2002); Dalibard (1999). In brief, any incident wave is decomposed in
plane waves of momentum $\boldsymbol{k}$, whose scattering yields outgoing
spherical waves of amplitude $f(\boldsymbol{k},{\boldsymbol{n}})$ in the
direction ${\boldsymbol{n}}=\boldsymbol{r}/r$. The scattering amplitude $f$
directly relates to the interaction potential $V(\boldsymbol{r})$ and, within
the first-order ’Born’ approximation, $f$ simply reduces to the Fourier
transform of $V(\boldsymbol{r})$. The scattering cross section $\sigma$ then
corresponds to the spherical integration of the square norm of $f$ over the
scattering direction
$\sigma(\boldsymbol{k})=\int|f(\boldsymbol{k},{\boldsymbol{n}})|^{2}d^{2}n$.

####  III.1.1 Partial wave expansion

In the case of an isotropic potential $V(\boldsymbol{r})=V(r)$, expanding the
collisional wavefunction in spherical harmonics is a very powerful tool as
each of the resulting radial waves are decoupled from one another. These are
indexed by $l$ (and $m_{l}$), corresponding to the quantum numbers for the
orbital angular momentum norm (and projection). For an isotropic potential,
only spherical harmonics with projection $m_{l}=0$ contribute, while in the
case of anisotropic interactions, $m_{l}\neq 0$ harmonics may also play a
role. Each spherical harmonic component of the expansion is an independent
solution of the 1D Schrödinger equation, where the potential $V(r)$ is
augmented by a centrifugal term $U_{{\rm c},l}(r)=\hbar^{2}l(l+1)/mr^{2}$. For
the isotropic case, asymptotically, the outgoing wave differs by only a phase
shift $\delta_{l}$ from the incoming wave. Therefore, the general scattering
amplitude $f$ decomposes into

$\displaystyle f(k,\theta)$ $\displaystyle=$
$\displaystyle\frac{1}{2ik}\sum_{l=0}^{\infty}(2l+1)(e^{2i\delta_{l}(k)}-1)P_{l}(\cos\theta)$
(25) $\displaystyle=$ $\displaystyle\sum_{l=0}^{\infty}f_{l}(k,\theta),$

where $\theta$ describes the angle between $\boldsymbol{k}$ and
$\boldsymbol{r}$ and varies between 0 and $\pi$. Here, $P_{l}$ denote the
ordinary Legendre polynomials. The scattering cross section then reads

$\displaystyle\sigma(k)$ $\displaystyle=$
$\displaystyle\sum_{l=0}^{\infty}\sigma_{l}(k),$ (26)
$\displaystyle\sigma_{l}(k)$ $\displaystyle=$
$\displaystyle\frac{4\pi}{k^{2}}(2l+1)\sin^{2}[\delta_{l}(k)].$ (27)

####  III.1.2 Ultracold limit and short-range interactions

For an interaction of short range $b$, we can define the temperature below
which the system is in the ultracold collisional regime as
$k_{B}T\ll\hbar^{2}/mb^{2}$. In this regime, the centrifugal barrier at $r=b$
is much larger than the typical kinetic energy of the colliding particles such
that they do not feel the interaction potential in $l>0$ harmonics (they are
reflected at $r>b$). Therefore, only $l=0$ contributes to the scattering.

##### Scattering length:

In the ultracold limit, the scattering amplitude and the cross section tend to
finite values given by the $k\rightarrow 0$ limit of their $l=0$
contributions: $f\approx f_{0}\rightarrow-a_{\rm s}$, and
$\sigma\approx\sigma_{0}\rightarrow 4\pi a_{\rm s}^{2}$. This is for the case
of distinguishable particles given by the sum over all $l$ in Eq. (26); see
Sec. III.1.3 for the role of quantum statistics in scattering. The $s$-wave
scattering length $a$ is then the only (non-universal) parameter describing
the ultracold scattering physics and it relates to the $s$-wave phase shift
via $a=\lim_{k\rightarrow 0}\delta_{0}(k)/k$. Hence, in the ultracold regime,
the details of the short-range interaction are all represented by the value of
a single parameter. Moreover, potentials yielding the same value of $a$ are
interchangeable. It is therefore useful to employ a particularly simple form
of the potential, the contact pseudo-potential $U_{s}(r)=g\delta(r)$, with
$g=\frac{4\pi\hbar^{2}a}{m}$ Pitaevskii and Stringari (2016); Pethick and
Smith (2002); Dalibard (1999); Fermi (1934); Huang and Yang (1957) 888This
simple formula is not regular everywhere. As described in Refs. Dalibard
(1999); Fermi (1934); Huang and Yang (1957), a more rigorous regular version
exists..

##### Power-law potentials and the case of van der Waals interactions:

The scaling of $\delta_{l}(k)$ can be deduced from the scattering potential
form by solving the associated 1D Schrödinger equation. As introduced in Sec.
I, in the case of a power-law interacting potential
$V(r)=V_{n}(r)=-C_{n}/r^{n}$, one finds that $\delta_{l}(k)$ scales as
$k^{2l+1}$ if $l<(n-3)/2$ and as $k^{n-2}$ otherwise Landau and Lifshitz
(1977); Dalibard (1999). Therefore, the ultracold regime described above only
holds for potentials decreasing sufficiently fast, as $n>3$.

Typically, this condition is satisfied by the dominant interaction between two
atoms in a dilute gas, the van der Waals interaction, which is an $n=6$ power
law. This scaling arises from the virtual exchange of photons between two
fluctuating electric dipoles: while on average the electric dipole of each
atom is zero (in the absence of an electric field), a dipolar fluctuation of
the electronic distribution of one atom can distort the other’s, leading to
the product $r^{-3}{\cdot}r^{-3}$ for the interaction dependence. The
ultracold limit is defined (in the absence of a direct DDI) by
$k_{B}T\sim\hbar^{2}/mb^{2}$. Associating $b$ with the van der Waals length
$r_{\rm vdW}=(mC_{6}/\hbar^{2})^{1/4}$, we find that the ultracold limit is
typically reached for temperatures of hundreds of $\mu{\rm K}$.

In the ultracold limit, the scattering cross section for partial waves
$l={0,1}$ at low collision energy $E=\hbar^{2}k^{2}/m$ follows the so-called
Wigner threshold law Wigner (1948):

$\sigma_{l}(E)\propto k^{4l}\propto E^{2l},$ (28)

while for $l>2$, $\sigma_{l}(E)\propto k^{6}\propto E^{3}$. Only the $l=0$
cross section does not vanish in the $E\rightarrow 0$ limit. This is in
contrast to the case of interacting polarised dipoles in which $n=3$. In this
case, the above scaling laws—based on the assumption of a short-range
potential—are not valid, and we will see (Sec. III.2) that all partial waves
contribute.

####  III.1.3 Role of quantum statistics

Quantum statistics plays a prominent role in the ultracold regime, as the
wavefunction should be appropriately symmetrised for collisions among
identical particles, i.e., atoms of the same isotope in the same spin state.
Indeed, the total collisional wavefunction of any pair of atoms—which is
comprised of the tensor product of wavefunctions describing both their spin
and orbital degrees of freedom—must be symmetric (antisymmetric) under
particle exchange for bosons (fermions). For example, if their spin state is
identical, then their motional state must be symmetric (antisymmetric), and
thus only even (odd) $l$ partial waves contribute to the partial wave
decomposition of Eqs. (25) and (26). These expressions may be rewritten as:

$f(k,\theta)=\frac{1}{ik}\sum_{l=0}^{\infty}\epsilon(l)(2l+1)(e^{2i\delta_{l}(k)}-1)P_{l}(\cos\theta),$
(29)

where $\theta$ now varies between 0 and $\pi/2$ and $\epsilon(l)$ depends on
the statistics of the particles. For bosons (fermions), $\epsilon(l)=0$ for
$l$ odd (even) and 1 in the diametric case $l$ even (odd). Similarly:

$\displaystyle\sigma(k)$ $\displaystyle=$
$\displaystyle\sum_{l=0}^{\infty}\sigma_{l}(k),$ (30)
$\displaystyle\sigma_{l}(k)$ $\displaystyle=$
$\displaystyle\frac{8\pi}{k^{2}}\epsilon(l)(2l+1)\sin^{2}[\delta_{l}(k)].$
(31)

For van der Waals interactions, the Wigner threshold law, Eq. (28), implies
that the cross section for the lowest partial wave accessible to fermions,
$l=1$, vanishes as $E^{2}$, while for bosons it tends toward a constant $8\pi
a^{2}$. One key consequence lies in the fact that while ultracold bosonic
gases can be evaporatively cooled to degeneracy via short-range elastic
thermalizing collisions, identical fermions cannot.

###  III.2 Universal elastic dipolar scattering

The above cross sections for van der Waals interactions may be said to be non-
universal, in that they depend on a parameter $a$ that cannot be calculated
from fundamental constants and fixed parameters of the system alone (like mass
or quantum numbers). Indeed, $a$ is calculable ab initio within state-of-the-
art methods only for light systems consisting of few electrons, such as
hydrogen and helium Przybytek  _et al._ (2010); Piszczatowski  _et al._
(2009). Accurately predicting the scattering length for more complex systems
is very challenging. In contrast, we will see that it is possible to derive
the dipolar cross section from first principles. This is due to the long-range
nature of $V_{\rm dd}$, which often allows one to neglect the short-range part
of the interaction potential. The longer-range character of the DDI ($n=3$)
also leads to a different Wigner threshold law; namely, one that exhibits a
constant, non-zero cross section even for identical fermions.

In contrast to the form of the $\delta_{l}(k)$ scaling presented above, for
$n=3$, $\delta_{l}(k)\sim k$ for all partial waves such that they all
contribute to the scattering cross section even in the low-energy limit.
Specifically, one can write $\delta_{l}(k)\approx Ak^{2l+1}+Bk$, where $A$ is
the strength of the short-range part of the potential and $B$ is the strength
of the long-range dipolar portion beyond the centrifugal barrier. $A$ and $B$
may be calculated under the first-order Born approximation 999This
approximation assumes that the interaction is sufficiently weak that there is
no rescattering during the collision. This means that the incident field can
be taken as the total field at the scatterer. Lagendijk  _et al._ (1986);
Shlyapnikov  _et al._ (1994); Fedichev  _et al._ (1996b); Bohn  _et al._
(2009b); Bohn and Jin (2014). This yields the total cross sections
$\sigma_{l\textrm{ even}}$ ($\sigma_{l\textrm{ odd}}$) for collision in even
(odd) partial waves:

$\displaystyle\sigma_{l\textrm{ even}}$ $\displaystyle=$ $\displaystyle 4\pi
a^{2}+\frac{16\pi}{45}D^{2},$ (32) $\displaystyle\sigma_{l\textrm{ odd}}$
$\displaystyle=$ $\displaystyle\frac{16\pi}{15}D^{2},$ (33)

with $D=S^{2}d^{2}m/\hbar^{2}=3a_{dd}$ and $a_{dd}$ as defined in Eq. (6). For
indistinguishable bosons (fermions) $\sigma=2\cdot\sigma_{l\textrm{ even}}$
($\sigma=2\sigma_{l\textrm{ odd}}$), while for distinguishable particles
$\sigma=\sigma_{l\textrm{ even}}+\sigma_{l\textrm{ odd}}$.

Regardless of the exchange statistics of the atoms or the value of $a$, the
cross sections of dipolar atoms are finite and energy independent as
$k\rightarrow 0$. Since $D$ contains only fundamental constants and known
quantum numbers, we may consider its contribution to the cross section
universal—i.e., not depending on the details of short-range physics Bohn  _et
al._ (2009b). However, this statement ignores the effects of molecular
potentials at certain values of magnetic field, as pointed out in Refs.
Lagendijk  _et al._ (1986); Pasquiou  _et al._ (2010) and later discussed in
Sec. III.3.1. In addition, as the dipole strength is increased, corrections
beyond the first-order Born approximation can cause $D$ to differ from its
universal value Ołdziejewski and Jachymski (2016).

An important consequence of the dipolar Wigner threshold law is the fact that
identical dipolar fermions can evaporatively cool themselves even at low
temperatures. This effect was used to cool spin-polarised fermionic polar
molecules of KRb, using forced evaporation in an optical dipole trap, to
$T\approx 2T_{F}$ Ni  _et al._ (2008)—most recently $T\approx 0.6T_{F}$ was
achieved in a similar system Valtolina  _et al._ (2020). Subsequently, a Fermi
gas of identical, spin-polarised 161Dy atoms was evaporated to $T/T_{F}\sim
0.7$ with ${\sim}$3$\times$103 atoms remaining Lu  _et al._ (2012). The deeply
degenerate regime was then reached with identical fermions of 167Er Aikawa
_et al._ (2014a), achieving $T/T_{F}\approx 0.1$ with several times $10^{4}$
atoms remaining. Similar temperatures and trap populations have also been
reached with 161Dy Burdick  _et al._ (2016b).

In Ref. Aikawa  _et al._ (2014a), the spin purity of the sample was
demonstrated and the elastic cross section was measured via cross-dimensional
rethermalisation measurements Monroe  _et al._ (1993); see also Sec. III.4.1.
The extracted elastic cross-section at small $T/T_{\rm F}$ was in good
agreement with the universal value of Eq. (33) without adjustable parameters.
Moreover, this study showed that the evaporation efficiency is among the
highest achieved in an ultracold gas experiment. This stems from the Pauli
suppression of inelastic processes for single-spin fermionic assemblies.
Together, these studies demonstrated the efficacy of universal elastic dipolar
scattering for degenerate Fermi gas production. This is a far simpler
technique than sympathetic cooling with spin or species mixtures Ketterle and
Zwierlein (2008).

###  III.3 Universal inelastic dipolar scattering: Dipolar relaxation

A large dipole moment enhances not only elastic dipolar scattering but also
the inelastic rate; see also Sec. I.3.3. While rapid elastic scattering is
useful for efficient evaporative cooling, inelastic scattering—dipolar
relaxation in this case—usually leads to gas heating and population loss.
Inelastic dipolar scattering arises from the anisotropy of the interaction and
causes atoms to change their spin state. In the presence of a finite external
magnetic field, the degeneracy between the sublevels of the ground state
hyperfine manifold is lifted via the Zeeman effect. In magnetic fields
inducing Zeeman splittings larger than $k_{B}T$, atoms in metastable spin
states relax to lower-energy, stronger-field-seeking Zeeman substates. The
Zeeman energy is released into the orbital motion of the collision pair,
conserving total orbital momentum by increasing the atoms’ kinetic energy.
Atoms in magnetic traps whose spins have flipped to strong-field seeking
states are lost. Lost too are atoms confined in traps—magnetic or
otherwise—that are shallower than the released Zeeman energy. Otherwise, the
relaxed atoms may remain confined but the change in their motional state leads
to an overall heating after rethermalisation via elastic collisions. We note
that for Zeeman splittings similar to or smaller than $k_{B}T$, the reverse
process called dipolar promotion is allowed. This was used for demagnetisation
cooling of a chromium gas; see Sec. VI.3.3 for details. For Zeeman splittings
similar to the harmonic trap excitation energy $\hbar\omega$, the change in
the motional state of the colliding atoms is spectroscopically resolvable.
This constitutes an atomic equivalent of the Einstein-de-Haas effect; see Sec.
VI.3 for details and Ref. Lahaye  _et al._ (2009).

Only atoms in their lowest-energy spin state—i.e., the maximally stretched
strong-field-seeking state—are immune to inelastic dipolar scattering (and
spin exchange) Pethick and Smith (2002). In this section, we will show that
dipolar relaxation can be seen as a short-range process for a large-enough
magnetic field. As a direct consequence, we will describe how confinement can
be used to control dipolar relaxation. We will also see that exchange
statistics plays a role in these inelastic scattering processes. In
particular, in Fermi gases, dipolar relaxation is suppressed for certain spin
configurations.

####  III.3.1 Inelastic cross sections

The DDI in Eq. (I.3.3) admits non-zero matrix elements between pairs of atoms
with differing spin and orbital angular momentum. The relevant selection rules
are $\Delta m^{i}_{F}=0,\pm 1$ for each of the two atoms $i={1,2}$ and $\Delta
l=0,\pm 2$ for their orbital motion. To conserve total angular momentum,
$\Delta m_{l}+\sum_{i=1,2}\Delta m^{i}_{F}=0$ must be satisfied. The cross
sections may be calculated in a first-order Born approximation treatment
Hensler  _et al._ (2003); Pasquiou  _et al._ (2010)101010Note that previous
theoretical treatments appeared in Refs. Lagendijk  _et al._ (1986);
Shlyapnikov  _et al._ (1994); Fedichev  _et al._ (1996b).. The appropriate
expression, neglecting all other interatomic potentials Pasquiou  _et al._
(2010), but including direct and exchange terms, is

$\displaystyle\sigma=\left(\frac{m}{4\pi\hbar^{2}}\right)^{2}\frac{1}{k_{\rm
i}k_{\rm f}}\left[\int|\widetilde{V_{\rm{dd}}}(\boldsymbol{k}_{\rm
i}-\boldsymbol{k^{\prime}})|^{2}\delta(|\boldsymbol{k^{\prime}}|-k_{\rm
f})\textrm{d}\boldsymbol{k^{\prime}}\right.$ (34)
$\displaystyle+\left.\tilde{\epsilon}\int\widetilde{V_{\rm{dd}}}(\boldsymbol{k}_{\rm
i}-\boldsymbol{k^{\prime}})\widetilde{V_{\rm{dd}}}^{*}(-\boldsymbol{k}_{\rm
i}-\boldsymbol{k^{\prime}})\delta(|\boldsymbol{k^{\prime}}|-k_{\rm
f})\textrm{d}\boldsymbol{k^{\prime}}\right],$

where $\widetilde{V_{\rm{dd}}}(\boldsymbol{k})$ is the Fourier-transformed
matrix elements of Eq. (I.3.3) and $\tilde{\epsilon}=[1,0,-1]$ for identical
bosons, distinguishable particles, and identical fermions, respectively
Hensler  _et al._ (2003); Lahaye  _et al._ (2009). The $\boldsymbol{k_{\rm
i}}$ and $\boldsymbol{k_{\rm f}}$ are the initial and final relative
wavevectors, respectively.

The cross sections for a collision channel starting with a two-body stretched
spin-state $|F,m_{F}=+F;F,m_{F}=+F\rangle$ and ending in either the same state
($\Delta m_{l}=0$) or a single ($\Delta m_{l}=1$) or double ($\Delta m_{l}=2$)
spin-flipped state form a hierarchy in powers of the spin $F$ Hensler  _et
al._ (2003); Pasquiou  _et al._ (2010):

$\displaystyle\sigma_{0}$ $\displaystyle=$
$\displaystyle\frac{16\pi}{45}F^{4}\left(\frac{\mu_{0}(g_{F}\mu_{B})^{2}m}{4\pi\hbar^{2}}\right)^{2}[1+\tilde{\epsilon}h\left(1\right)],$
(35) $\displaystyle\sigma_{1}$ $\displaystyle=$
$\displaystyle\frac{8\pi}{15}F^{3}\left(\frac{\mu_{0}(g_{F}\mu_{B})^{2}m}{4\pi\hbar^{2}}\right)^{2}[1+\tilde{\epsilon}h\left(\frac{k_{\rm
f}}{k_{\rm i}}\right)]\frac{k_{\rm f}}{k_{\rm i}},$ (36)
$\displaystyle\sigma_{2}$ $\displaystyle=$
$\displaystyle\frac{8\pi}{15}F^{2}\left(\frac{\mu_{0}(g_{F}\mu_{B})^{2}m}{4\pi\hbar^{2}}\right)^{2}[1+\tilde{\epsilon}h\left(\frac{k_{\rm
f}}{k_{\rm i}}\right)]\frac{k_{\rm f}}{k_{\rm i}}.$ (37)

In these expressions, $\sigma_{0}$ is the elastic cross section discussed in
the previous section; $\sigma_{1}$ and $\sigma_{2}$ are the cross sections for
single spin-flip and double spin-flip dipolar relaxation processes,
respectively. The kinematic terms $\frac{k_{\rm f}}{k_{\rm
i}}=\sqrt{1+\frac{m\Delta E_{\Delta m_{l}}}{\hbar^{2}k^{2}_{i}}}$ are dictated
by the energy conservation condition and depend on the number $\Delta m_{l}$
of spin flips involved in the collision. This is because $\Delta E_{\Delta
m_{l}}$ depends on the change in kinetic energy of the atomic pair. Typically
$\Delta E_{\Delta m_{l}}=\Delta m_{l}\Delta E$, where the Zeeman energy
splitting $\Delta E\propto B$. The value of $h\left(x\right)$, defined on
$(1,\infty)$, monotonically increases from $h\left(1\right)=-1/2$ to
$h(x\rightarrow\infty)=1-4/x^{2}$ Hensler  _et al._ (2003); Bohn  _et al._
(2009b); Pasquiou  _et al._ (2010).

The inelastic cross sections in Eqs. 36 and 37 contain only fundamental
constants, atomic quantum numbers, and kinematic variables. Therefore, just
like the elastic dipolar collisions, one may consider inelastic dipolar
scattering a universal process in the first-order Born approximation limit. We
note that Ref. Pasquiou  _et al._ (2010) extended the cross section
calculation to consider the effect of molecular potentials, and in particular,
that of the short-range van der Waals interaction. The authors demonstrate
theoretically and experimentally the effect of a finite scattering length on
the dipolar relaxation cross section. While Eqs. (36-37) are good
approximations at low magnetic field for which $k_{\rm f}a\ll 1$, the effect
of the molecular potential becomes important when $B$ is such that $k_{\rm
f}a\gtrsim 1$. This leads to deviations from universal behaviour (see below).

####  III.3.2 Inelastic cross section measurements

Experiments typically do not give direct access to the cross section, but
rather to collisional rates. In a gas at thermal equilibrium, the dipolar
relaxation collisional rate is deduced from the above cross section via
$\beta_{dr}=2\langle(\sigma_{1}+\sigma_{2})v_{\rm
rel}\rangle_{\textrm{\small{th}}}$, where $v_{\rm rel}=2\hbar k_{i}/m$
describes the initial relative velocity between a pair of particle and
$\langle.\rangle_{\textrm{\small{th}}}$ defines a thermal average over the
gas. This collisional rate has two contributions, one from collisions
involving single spin-flips and one from double spin-flips respectively. The
latter is $\beta_{dr}^{(1,2)}=2\langle\sigma_{1,2}v_{\rm
rel}\rangle_{\textrm{\small{th}}}$.

Following earlier studies of weakly magnetic atoms Stoof  _et al._ (1988);
Moerdijk and Verhaar (1996); Gerton  _et al._ (1999), dipolar relaxation was
quantitatively studied in bosonic Cr Hensler  _et al._ (2003); Weinstein  _et
al._ (2002). In these studies, 52Cr was magnetically trapped in the highest
energy (i.e., low-field-seeking) Zeeman sublevel
$|S,m_{s}\rangle=|3,3\rangle$. Dipolar relaxation in Cr takes the atoms
initially in the state $|S,m_{s};S,m_{s}\rangle=|3,3;3,3\rangle$ to
$|3,3;3,2\rangle_{\mathcal{S}}$ ($|3,2;3,2\rangle$) via a single (double) spin
flip, where $\mathcal{S}$ denotes the symmetrised total wave function. The
spin-flipped atoms may be lost from the trap or, if not, heat the gas due to
their increased kinetic energy. The single and double spin-flip processes were
characterised by studying both the rate of population loss from
$|S,m_{s}\rangle=|3,3\rangle$ and the rate of temperature increase in the gas.
The overall dipolar relaxation rate was found to be four orders–of–magnitude
larger than in weakly magnetic alkali-metal atoms Hensler  _et al._ (2003).
Moreover, the data were in reasonable agreement with cross section predictions
based on the first-order Born approximation given above, i.e., Eqs.(36) and
(37).

Dipolar relaxation processes limit the temperature and density achievable in
magnetically trapped gases of highly magnetic atoms. In contrast, optical
traps enable the trapping of all Zeeman sublevels, and in particular of the
lowest of these states, thus providing immunity to dipolar relaxation and
opening the way to reaching higher gas phase-space densities. Additional
dipolar relaxation data for Cr Pasquiou  _et al._ (2010) and Dy Burdick  _et
al._ (2015) at lower temperatures (100’s of nK) were performed using optical
trapping 111111We note that the dipolar relaxation rates have also been
measured in bosonic Dy at a few hundred $\mu$K Lu  _et al._ (2010) and Er and
Tm at hundreds of mK Connolly  _et al._ (2010), as well as Ho and Dy at
similar temperatures Newman  _et al._ (2011).. In optical traps, because of
limited trap depth, the particles involved in dipolar relaxation are typically
lost (except at very low magnetic field, see Sec. VI for details) and loss
spectroscopy enables the characterisation of the collision rate. These
additional experimental results also agreed well with predictions using the
first-order Born approximation. Because the first-order Born approximation is
valid for only weakly interacting systems (stronger interactions cause higher-
order terms to be non-negligible), these findings imply that the atoms’ DDI
are not so strong as to invalidate this convenient approximation. To conclude,
the first-order Born approximation should hold for dipolar collisions among
any element (Dy being the most magnetic stable element of the periodic table)
as long as the atomic density remains on a similar scale as investigated in
these works. Such densities are typical of repulsive three-dimensional gases.

We note that the magnitude of $\sigma_{1}$ for Dy at $B=0$ can be up to
100$\times$-larger than that of the alkali metal caesium. In addition, Ref.
Burdick  _et al._ (2015) pointed out that, while the ratio of
$\sigma_{1}/\sigma_{2}=F^{-1}\ll 1$ for the stretched states of large-spin
atoms, i.e., those with $|m_{F}|\approx F$, these same atoms exhibit large
ratios of $\sigma_{2}/\sigma_{1}$ for states near $m_{F}=0$. Indeed,
$\sigma_{2}$ for these states can be larger than $\sigma_{1}$ for the
stretched states: One cannot avoid dipolar relaxation simply by employing spin
states near $m_{F}=0$.

####  III.3.3 Interaction range for dipolar relaxation

While the elastic dipolar interaction is long-ranged, Ref. Pasquiou  _et al._
(2010) shows that the dipolar relaxation processes are intrinsically short-
ranged, despite the same $1/r^{3}$ scaling. This effect stems from the reduced
overlap between the incoming and outgoing waves due to the increase in kinetic
energy. Specifically, while the integral overlap of Eq. (34) for $\sigma_{0}$
involves incoming and outgoing waves oscillating at the same frequency
($k_{\rm i}$), for inelastic collisions, the incoming and outgoing
wavefunctions oscillate at different frequencies and become spatially
mismatched at long distance, averaging to zero. Assuming a vanishingly small
collision energy (low $T$), the range for this effective cancellation of the
overlap, $R_{\rm dr}$, scales as $1/k_{\rm f}$ . This is set by the release of
Zeeman energy in the spin relaxation. Hence, the range of dipolar relaxation
decreases with $B$ following:

$R_{\rm dr}(B)\approx\frac{\hbar}{\sqrt{mg_{F}\mu_{\rm{B}}B}}.$ (38)

In terms of partial-wave decomposition, dipolar relaxation occurs at the
interparticle distance $R_{\rm dr}$ that matches the distance at which the
energy released in dipolar relaxation $\Delta E=g_{F}\mu_{B}B$ is compensated
for by the centrifugal energy of the output channel. This simple statement
assumes that the input channel is $l=0$, which is usually the case at low $T$.
(To satisfy the $\delta l=2$ selection rule, the output channel is $l=2$.)

The localised character of dipolar relaxation has been revealed in the
interplay between dipolar relaxation process and other interatomic (molecular)
potentials in Ref. Pasquiou  _et al._ (2010). This results in variations of
the relaxation rate with $B$ that deviate from that expected from the first-
order Born approximation when neglecting the effect of molecular potentials;
specifically, the rate can dip far below this expectation at finite $B$. This
strong reduction of dipolar relaxation is shown to arise from a node in the
initial wavefunction, preventing the presence of pairs of atoms at an inter-
particle distance corresponding to the distance $R_{\rm dr}$ at which dipolar
relaxation occurs. This provides a factor of twenty reduction in the
relaxation rate versus the rate predicted by the first-order Born
approximation. This reduction relies on a non-universal molecular potential
configuration of 52Cr and it cannot be directly extended to the other highly
magnetic atoms due in particular to the complexity of their Feshbach resonance
spectrum.

####  III.3.4 Quantum statistics and dipolar relaxation: Bose enhancement and
Fermi suppression.

Now that we have discussed the short-range nature of dipolar relaxation, we
pivot to a discussion of how this fact combines with quantum statistics to
affect relaxation rates. As with the short-ranged van der Waals interaction
(see Sec. III.1.3), one might expect a suppression of collisions between
identical fermions due to quantum statistics. Indeed, this is the case,
because identical fermions cannot surmount the centrifugal barrier of the
lowest quantum statistically allowed partial wave $l=1$: since identical
fermions cannot closely approach, the short-ranged inelastic DDI does not lead
to dipolar relaxation of their spin. Conversely, we will show below that the
relaxation rate is enhanced for identical bosons, since they can collide on
the barrierless $l=0$ partial wave.

The dependence of the relaxation rate on quantum statistics is manifest in the
form of the terms containing the ratios of incoming and outgoing momenta in
Eqs. (36) and (37). At high $B$ and low $T$, the kinematic terms diverge as
$\langle\frac{k_{\rm f}}{k_{\rm i}}\rangle_{\rm th}\propto\sqrt{B/T}$, where
the $\rm th$ denotes a thermal average. In this limit, $\sigma_{1}$ vanishes
as $4\sqrt{T/B}$ for indistinguishable fermions ($\tilde{\epsilon}=-1$), while
it increases as $2\sqrt{B/T}$ for indistinguishable bosons
($\tilde{\epsilon}=+1$) and is $\sqrt{B/T}$ for distinguishable particles
($\tilde{\epsilon}=0$). The relative suppression ratio in this limit becomes
$\sigma_{1}^{\textrm{fermions}}/\sigma^{\textrm{bosons}}_{1}\propto 2T/B$.
This scaling exhibits two important facts: i) Dipolar relaxation rates grow
worse for identical bosons versus $B$, e.g., $\sigma_{1}$ for bosonic Dy is
$10^{3}$ times larger than Cs’s at only a 1-G field; ii) $\sigma_{1,2}$ can be
suppressed for identical fermions, e.g., dipolar relaxation in fermionic Dy is
only $10{\times}$-worse than Cs’s at a few G, as experimentally demonstrated
in Ref. Burdick  _et al._ (2015). It is remarkable, and counterintuitive, that
the more exothermic the fermionic dipolar relaxation would be, the less likely
they are to occur. This constitutes another striking example of the role
quantum statistics plays in ultracold collisions.

Figure 5 shows measured dipolar relaxation rates for spin polarised bosonic
and fermionic Dy. Reference Burdick  _et al._ (2015) also shows the fermionic
suppression of dipolar relaxation of spin mixtures for processes in which the
output state of the fermions consists of identical spins Pasquiou  _et al._
(2010). This corresponds to the time-reversed, state-changing process of
indistinguishable fermions colliding. Importantly, the stability of a spin
mixture in the lowest-two Zeeman states opens interesting prospects for the
study of many-body spin physics. In particular, one could explore how the DDI
impacts the BEC–to–BCS crossover (see also Sec. IV.2) Baranov  _et al._
(2012). The stability of these states was also exploited for the long-lived
implementation of artificial spin-orbit coupling in Dy gases; see Sec. VI.4.1.

Figure 5: Dipolar relaxation for bosonic 162Dy ($m_{F}=+8$) and fermionic
161Dy ($m_{F}=+21/2$). Two-body collisional loss rates for 162Dy (squares) and
161Dy (triangles). Curves are calculated collisional loss rates using Eq. (36)
at $T\approx 400$ nK with no free parameters. Adapted from Ref. Burdick  _et
al._ (2015).

####  III.3.5 Control of dipolar relaxation via confinement

Confinement of atoms in lower dimensions can also suppress dipolar relaxation
rates through the phase-space reduction of outgoing scattering channels. That
is, the relaxation rate can be lowered via the interplay between the range of
the dipolar relaxation process and the length scales of the confinement
potential. Following Eq. (38), when the magnetic field $B$ is small, dipolar
relaxation takes place at large inter-particle distances. Dipolar relaxation
thus strongly depends on the trapping geometry, provided the size of the cloud
is comparable to $R_{DR}$. This phenomenon was demonstrated in Pasquiou  _et
al._ (2010, 2011a); de Paz  _et al._ (2013a) for the case of 1D, 2D, and 3D
optical lattices. In Ref. Pasquiou  _et al._ (2010), a 1D optical lattice was
shown to reduce dipolar relaxation in Cr by a factor of seven compared to a 3D
trap when $\Delta E$ was smaller than the gap to first transverse excited
state. The dipolar relaxation rate was further suppressed—by three
orders–of–magnitude—in a 2D optical lattice of 52Cr Pasquiou  _et al._
(2011a). This was achieved below the threshold $B$-field determined by this
gap condition for $\Delta E$. Interband transitions mediated by dipolar
relaxation were observed above this threshold, as illustrated in Fig. 6. In
the specific case of a deep 3D optical lattice de Paz  _et al._ (2013a), it
was observed that dipolar relaxation becomes a resonant process as a function
of $B$ due to the quantisation of kinetic energy in the 3D-band structure. We
see that by adjusting confinement and the magnetic field, dipolar relaxation
can be either eliminated (e.g., to study spinor physics at constant
magnetisation), or be used to couple different lattice bands. This then
realises an intrinsic nonlinear spin-orbit coupling in the lattice de Paz  _et
al._ (2013a); Sun and You (2007).

Figure 6: Control of dipolar relaxation in reduced dimensions. The false
colour absorption image is a band mapping of the atom population trapped in a
2D optical lattice. (a) Integrated population profile along $z$ yields
populations in the different Brillouin zones (BZ). Below the threshold
magnetic field for dipolar relaxation, only the first BZ is populated (grey
profile with the sharper peak). Above the threshold, the first transverse
excited band (second BZ) is populated (red profile with the wider peak). (b)
Integrated population profile along $y$ shows the (1D) velocity distribution
in the sites of the 2D lattice. Above threshold, vibrational de-excitation of
atoms from excited bands create a non-Gaussian distribution. (c) Fraction of
atoms detected in the first excited band as a function of the Larmor frequency
after a 25-ms duration allowed for dipolar relaxation. The lattice depth is
$25E_{r}$ which sets a threshold Larmor frequency of 120 kHz. Below threshold,
dipolar relaxation is reduced by three orders of magnitude. Adapted from
Pasquiou  _et al._ (2011a).

Modifications to dipolar relaxation rates were also demonstrated in Ref.
Pasquiou  _et al._ (2010) using RF fields. While only an increase in
relaxation rates was experimentally demonstrated, the work showed that RF-
dressing might serve as a valuable tool to control dipolar relaxation.

####  III.3.6 Dipolar relaxation and spinor physics

Section VI will describe in more detail why the large spin of the highly
magnetic atoms presents attractive prospects for the construction of exotic
spinor states and many-body phases. Dipolar relaxation can play a positive
role in this regard: magnetisation-changing collisions open the door to the
exploration of spinor physics under conditions of free magnetisation. For
example, controlling the field at the mG level makes the Zeeman energy $\Delta
E$ comparable to thermal excitations of $\sim k_{\rm B}$100 nK. This induces
spontaneous, incoherent demagnetisation from the absolute ground state.
Further control of the field down to the 100-$\mu$G level allows $\Delta E$ to
be of the order of spin-dependent interactions. This would enable the
observation of spinor phases driven by contact interactions under free
magnetisation. This degree of control has already been achieved in current
ultracold atoms experiments Pasquiou  _et al._ (2011b). With even better
control, possible experimentally, one can enter the regime where dipolar
relaxation is relevant from the many-body point of view. For example, if
$R_{\rm dr}>\bar{d}$, where $\bar{d}=n^{-1/3}$ is the average particle
distance, then many-body effects may arise. Interesting vortex structures may
appear due to the Einstein-de-Haas effect since the orbital angular momentum
is increased from spin relaxation Kawaguchi  _et al._ (2006a); Santos and Pfau
(2006); Gawryluk  _et al._ (2007). Strongly rotating and interacting Fermi
gases approaching a Laughlin state may appear Peter  _et al._ (2013).

In general, however, dipolar relaxation has unfortunate consequences for the
study of many-body spinor physics in fields larger than those considered above
or outside of 2D or 3D lattices where it can be strongly suppressed. High
densities are needed to enhance the interaction energy, but fast dipolar
relaxation can render systems with metastable spin states too fragile to
observe in (near-)equilibrium situations. Indeed, dipolar relaxation for all
metastable spin states in bosons is particularly severe, and only a few
metastable spin configurations for fermions are sufficiently long-lived
Burdick  _et al._ (2015). For example, long-lived mixtures of
$|F,-(F-1);F,-F\rangle$ were achieved in fermionic Dy and Er Burdick  _et al._
(2015); Baier  _et al._ (2018) (see Sec. VI.2.4) and the same set of states
was used to generate spin-orbit coupling in fermionic Dy Burdick  _et al._
(2016b); see Sec. VI.4.1. The study of out-of-equilibrium dynamics is also
possible for time-scales up to the dipolar relaxation time scale de Paz  _et
al._ (2016); Lepoutre  _et al._ (2018a); see Sec. VI.2.3.

###  III.4 Anisotropic scattering

Figure 7: a) Total dipolar cross section
$\sigma=\int\frac{d\sigma}{d\Omega}(\boldsymbol{k},\boldsymbol{k}^{\prime})d\Omega_{\boldsymbol{k}^{\prime}}$
in units of the dipolar length $a_{\rm d}=a_{\rm dd}/3$ versus the angle
$\eta$ between the incident scattering direction $\hat{\boldsymbol{k}}$ and
the polarisation direction $\hat{\epsilon}$. The curves are for
distinguishable particles (solid black), identical fermions (dashed blue), and
identical bosons (dotted red). b) 3D plots of $\sigma$ as a function of
$\hat{\boldsymbol{k}}$ with $\hat{\epsilon}$ set to the vertical axis.
Reproduced from Ref. Bohn and Jin (2014).

We now return to the discussion of elastic dipolar scattering and explore
manifestations of the anisotropic nature of the DDI on scattering into
different directions in space. This 3D scattering is characterised by the
dipolar differential cross section $\frac{d\sigma}{d\Omega}$. From the
scattering theory reviewed in Sec. III.1, we see that the differential cross
section directly relates to the scattering amplitude
$\frac{d\sigma}{d\Omega}(\boldsymbol{k},{\boldsymbol{n}})=|f(\boldsymbol{k},{\boldsymbol{n}})|^{2}$.
The expressions $\frac{d\sigma_{\rm B/F}}{d\Omega}$, including both the
contact interaction and DDI, for identical bosons and fermions in the first-
order Born approximation are, respectively:

$\displaystyle\frac{d\sigma_{\rm
B}}{d\Omega}(\boldsymbol{k},\boldsymbol{k}^{\prime})=\frac{D^{2}}{8}\left[\frac{4}{3}-\frac{2a}{3a_{\rm
dd}}\right.$
$\displaystyle\left.-\frac{2(\hat{\boldsymbol{k}}.\hat{\boldsymbol{\varepsilon}})^{2}+2(\hat{\boldsymbol{k}}^{\prime}.\hat{\boldsymbol{\varepsilon}})^{2}-4(\hat{\boldsymbol{k}}.\hat{\boldsymbol{\varepsilon}})(\hat{\boldsymbol{k}}^{\prime}.\hat{\boldsymbol{\varepsilon}})(\hat{\boldsymbol{k}}.\hat{\boldsymbol{k}}^{\prime})}{1-(\hat{\boldsymbol{k}}.\hat{\boldsymbol{k}}^{\prime})^{2}}\right]^{2},$
(39) $\displaystyle\frac{d\sigma_{\rm
F}}{d\Omega}(\boldsymbol{k},\boldsymbol{k}^{\prime})=$
$\displaystyle\frac{D^{2}}{8}\left[\frac{4(\hat{\boldsymbol{k}}.\hat{\boldsymbol{\varepsilon}})(\hat{\boldsymbol{k}}^{\prime}.\hat{\boldsymbol{\varepsilon}})-2\left[(\hat{\boldsymbol{k}}.\hat{\boldsymbol{\varepsilon}})^{2}+(\hat{\boldsymbol{k}}^{\prime}.\hat{\boldsymbol{\varepsilon}})^{2}\right](\hat{\boldsymbol{k}}.\hat{\boldsymbol{k}}^{\prime})}{1-(\hat{\boldsymbol{k}}.\hat{\boldsymbol{k}}^{\prime})^{2}}\right]^{2},$
(40)

where $\boldsymbol{k}$ and $\boldsymbol{k}^{\prime}$ denote the relative
momenta before and after the collision, with directions
$\hat{\boldsymbol{k}}=\frac{\boldsymbol{k}}{|\boldsymbol{k}|}$ and
$\hat{\boldsymbol{k}}^{\prime}=\frac{\boldsymbol{k}^{\prime}}{|\boldsymbol{k}^{\prime}|}$
Bohn and Jin (2014); Sykes and Bohn (2015). The dipole moments are aligned
along the $B$-field direction $\hat{\boldsymbol{\varepsilon}}$. These are
nontrivial functions of the relative orientations of the three directions in
the problem, the incident and scattered wavevectors and the polarisation
direction of the dipoles. The total dipolar cross sections are found by
integrating these expressions over $d\Omega_{\hat{}\boldsymbol{k}^{\prime}}$
and are shown in Fig. 7. For a general angle $\eta$ between
$\hat{\boldsymbol{k}}$ and $\hat{\epsilon}$, these are

$\displaystyle\sigma_{B}(\hat{\bf p}_{\rm rel})$ $\displaystyle=$
$\displaystyle D^{2}\frac{\pi}{9}[72a^{2}-24a(1-3\cos^{2}\eta)$ (41)
$\displaystyle\,\,+11-30\cos^{2}\eta+27\cos^{4}\eta]$
$\displaystyle\sigma_{F}(\hat{\bf p}_{\rm rel})$ $\displaystyle=$
$\displaystyle
D^{2}\frac{\pi}{3}\left[3+18\cos^{2}\eta-13\cos^{4}\eta\right].$ (42)

Note that the results of Eqs. 32 and 33 correspond to angle-averaged versions
of the above expressions.

####  III.4.1 Anisotropy of rethermalisation

Figure 8:  (a,b) Upper panels: Number of collisions per rethermalisation,
$\alpha$, versus the angle between the quantization axis $\hat{\epsilon}$ and
the trap axis along which rethermalisation occurs, $\lambda$. Lower panels: 3D
plots of $d\sigma/d\Omega$ as a function of the scattering direction
$\hat{}\boldsymbol{k}^{\prime}$, for different values of $\eta$. $\eta$ is the
angle between $\hat{}\boldsymbol{k}$, set to the vertical $z$ axis, and the
polarisation direction $\epsilon$ in the $xz$-plane. Plots are for a)
identical bosons and b) identical fermions. Reproduced from Ref. Bohn and Jin
(2014). c) Comparison between theory and measurement of $\alpha$ for ultracold
dipolar fermions of 167Er. Reproduced from Ref. Aikawa  _et al._ (2014b).

The anisotropy of the differential cross section may reveal itself in the
thermalisation dynamics, and in particular, via the thermalisation rate as a
function of the angle $\lambda$ between $\hat{\boldsymbol{\varepsilon}}$ and
the axis along which energy is imparted. This can be studied in cross-
dimensional rethermalisation experiments. Cross-dimensional rethermalisation
(relaxation) is a method for measuring the total elastic cross section wherein
the gas is brought out of thermal equilibrium by heating it along one or two
of three directions Monroe  _et al._ (1993). One then observes the time for
the temperature in the third direction to equilibrate. The time constant
$\tau$ for this typically exponential relaxation process may be related to the
cross section via $\tau=\alpha/nv\sigma$, where $n$ is the average number
density, $v$ is the mean relative speed, and $\sigma$ is the elastic cross
section. The anisotropy of the differential cross section appears through the
parameter $\alpha$, which is the number of collisions, on average, for
rethermalisation. For $s$-wave ($p$-wave) collisions, $\alpha=2.5$ (25/6
$\approx$ 4.17) DeMarco  _et al._ (1999). However, for dipolar particles,
$\alpha$ becomes angle dependent due to the anisotropic differential cross
section Bohn and Jin (2014). Its value is plotted in Fig. 8(a,b) versus the
angle $\lambda$ between the axis along which energy is imparted and
$\hat{\boldsymbol{\varepsilon}}$. We see that $\alpha$ can change by more than
a factor two versus angle for fermions, while sightly less than two for
bosons. The control of $\lambda$ is therefore an important tool for
optimising, e.g., evaporative cooling efficiency: rethermalisation of fermions
(bosons) is far more efficient at 45∘ (90∘) than at 0∘ or 90∘ (0∘).

The large anisotropy in the rethermalisation of dipolar gases was
experimentally demonstrated in ultracold gases of fermionic 167Er Aikawa  _et
al._ (2014b). For fermionic atoms, $\sigma$ is simply given by the universal
formula Eq. (33) so that a measurement of $\tau$ directly provides $\alpha$.
Figure 8(c) shows measurements of $\alpha$ in agreement with the dependence
predicted in Ref. Bohn and Jin (2014). The authors also reported the reduction
in rethermalisation rate due to Pauli-suppression of scattering once $T$ was
sufficiently lower than $T_{F}$ (see Sec. III.2 and Eqs.32 and 33). Unlike
$\alpha$, the suppression factor was observed to be independent of $\lambda$,
implying that the occupation of the final density of states leading to Pauli-
blocking is mostly unaffected by the anisotropy of the DDI.

Cross-dimensional rethermalisation experiments in ultracold gases where also
performed in bosonic isotopes of 164Dy, 162Dy Tang  _et al._ (2015b) and of
164Er, 166Er, 168Er and 170Er Frisch  _et al._ (2015). These studies measured
the unknown values of the scattering lengths $a$ of these atoms by using
either the theoretically predicted $\alpha$ or by using more extensive
dynamics simulations to convert the measured $\tau$ into values of $a$; see
Eq. (32) and Sec. II.4.2 for more details. In Ref. Tang  _et al._ (2015b), the
authors used a thorough theoretical analysis informed by data from different
values of the angle $\lambda$ to extract $a$. This study confirmed the
anisotropic behaviour as well as provided evidence of hydrodynamical effects
in the rethermalisation dynamics Sykes and Bohn (2015).

These hydrodynamic effects were manifest as modifications of the mechanical
oscillations of the gas after the gas was kicked out of thermal equilibrium.
It was observed that the large magnetic dipole moment of Dy provides a
sufficiently large elastic collision cross section for the gas to lie near the
hydrodynamic collisional regime. These effects were predicted in Ref. Sykes
and Bohn (2015), where it was shown that the relaxation no longer follows a
simple exponential in the hydrodynamic regime when the collisional rate is
similar to the trapping frequencies, as observed in Ref. Tang  _et al._
(2015b). Nevertheless, cross sections were able to be extracted from the data
through close comparison to Monte Carlo simulations. Section II.4.2 describes
the scattering length measurements obtained as a result. Hydrodynamic
behaviour arising from the DDI has also been observed in the expansion aspect
ratio of Dy thermal gases Tang  _et al._ (2016b); see also Sec. IV.1.2.

####  III.4.2 Anisotropic scattering halos

The anisotropy of the differential cross sections may be directly visualised
by observing the shape of scattering “halos” of atoms from two gases that have
undergone a head-on collision. Reference Burdick  _et al._ (2016a) describes
such an experiment: One creates two dipolar BECs of 162Dy counter-propagating
at momentum $\pm 2\hbar k_{R}$ using Bragg diffraction. This energy is much
larger than the internal momentum distribution width of the individual BEC.
When the BECs propagate through one another, atoms scatter away from the
forward and backward directions of the BECs’ motion. A halo rapidly expands
from the centre-of-mass position, as shown in the data of Fig. 9. The BECs are
sufficiently dilute that atoms scatter at most one time. The shape of the halo
therefore reveals information about the two-body differential cross section.

Figure 9: 2D images of the scattering halos from the collision of two BECs of
162Dy that are counter-propagating along the horizontal direction of the
images. Left column are simulations, while the right column are experimental
data. Dipoles are aligned using a $B$-field oriented: a) along the vertical
direction, b) out of the page, and c) along the horizontal direction, as
indicated by the arrows. Adapted from Ref. Burdick  _et al._ (2016a).

The data visually reveal a novel regime of quantum scattering. Rather than
exhibiting halos with shapes dominated by a single partial wave—e.g., the
spherically symmetric halos from the $s$-wave scattering of ultracold
identical bosons Legere and Gibble (1998); Chikkatur  _et al._ (2000); Thomas
_et al._ (2004); Buggle  _et al._ (2004)—anisotropic shapes are manifest. This
is a consequence of the superposition of a large number of partial waves due
to high-order angular momentum coupling through the elastic DDI. Moreover, the
particular halo shape strongly depends on $\eta$, as expected from Eq. (39).
The shapes of the halos in the data agree with Monte Carlo simulations based
on Eq. (39) and described in Ref. Sykes and Bohn (2015).

##  IV Fully polarised repulsive bulk dipolar gases

In this chapter, we will consider dipolar quantum gases of atoms fully
polarised in the lowest Zeeman state. We will cover only harmonically trapped
gases in the absence of a three dimensional lattice; see Sec. VII for dipolar
lattice systems. Most of the experimental activity has been focused on dipolar
bosonic gases and in particular dBECs (Sec. IV.1). Nevertheless, we will also
include the relevant results on dipolar DFGs (Sec. IV.2) as well as results on
dipolar gases in lower dimensions, in particular in 1D (Sec. VII.1).

###  IV.1 Weakly interacting three-dimensional dipolar Bose gases

A majority of the phenomena arising from the DDI that have been observed in
bulk BECs, in particular in the regime where the overall interactions are
effectively repulsive, are well accounted for by a mean-field (MF) description
of the interactions. After a short description of the MF theory of dipolar
Bose gases in Sec. IV.1.1, this section will focus on experimental
observations of related phenomena. These observations can be understood as
hydrodynamic measurements, as they concern the linear response from the
equilibrium state of the quantum fluid; e.g., its time evolution in a trap or
in free space as well as its elementary collective excitations. These aspects
will be reviewed in Secs. IV.1.2 and IV.1.3. Beside modifying the speed of
sound and giving it a directional dependence, the DDI profoundly affects the
structure of the elementary excitation spectra of trapped BECs by the
interplay of long-range anisotropic interactions with the trap geometry. In
particular, a dispersion relation showing a local minimum in energy, i.e., a
roton mode, was predicted and observed, as discussed in Sec. IV.1.3.
Furthermore, the DDI influences the stability of BECs in traps. Due to the
anisotropy of the DDI, the effective sign of the MF interaction depends on the
density distribution, which, in turn, is not only determined by the trap
geometry but also by the DDI. This interplay results in a non-trivial
stability diagram as a function of $a$ and of the anisotropy of the trapping
potential, as well as instability behaviours which also connects to the
underlying modifications of the excitations spectrum of the BEC mentioned
above. This will be reviewed in Sec. IV.1.4. The stability analysis and
related stabilised states of dipolar Bose gases beyond a MF treatment is the
subject of Sec. V.

####  IV.1.1 The mean-field description

##### Interaction potential in dilute gases

As introduced in Sec. I.4.1, polarised dipolar bosons in an ultracold dilute
ensemble can be considered to be interacting via a potential of the form:

$\displaystyle
V_{\rm{int}}(\boldsymbol{r})=\frac{4\pi\hbar^{2}\,a}{m}\delta(\boldsymbol{r})+\frac{3\hbar^{2}\,a_{\rm
dd}}{m}\frac{1-3\cos^{2}\theta}{r^{3}}.$ (43)

Again, $a$ is the $s$-wave scattering, while $a_{\rm
dd}=\mu_{0}m\mu^{2}/12\pi\hbar^{2}$ is the dipolar scattering length. The
interparticle separation vector is $\boldsymbol{r}$, with norm $r$ and angle
$\theta$ with respect to the dipole orientation. In the following, we will
assume a polarising magnetic field aligned along the $z$-axis, and a harmonic
trap with frequencies $\omega_{x,y,z}$. Unless otherwise specified, we also
consider that the trap has a cylindrical symmetry along $z$ such that
$\omega_{x}=\omega_{y}=\omega_{\rho}$, and

$\displaystyle
V_{\rm{tr}}=\frac{m}{2}\left(\omega_{\rho}^{2}\,\rho^{2}+\omega_{z}^{2}\,z^{2}\right).$
(44)

The trap aspect ratio is

$\displaystyle\lambda=\omega_{z}/\omega_{\rho}.$ (45)

The validity of Eq. 43 as an effective pseudo-potential written as a simple
sum of a contact pseudo-potential $g\delta(\boldsymbol{r})$ (or its more
rigorous, regularised version Huang and Yang (1957)) and the DDI potential Eq.
(11) has long been debated Lahaye  _et al._ (2009); Pitaevskii and Stringari
(2016). This expression for the potential results from the first-order Born
approximation applied to a molecular potential of the dipole-dipole type
dominating at long distance plus a van der Waals potential dominating at short
range Yi and You (2000, 2001); Derevianko (2003, 2005a, 2005b); Bohn  _et al._
(2009b). While it is a good approximation for weak dipoles and away from
scattering resonances, modifications may be necessary for systems with large
$\varepsilon_{\rm{dd}}$. In particular, the $s$-wave scattering length may
depend on $d$ Yi and You (2000, 2001); Mazets and Kurizki (2007), requiring a
renormalisation of the dipolar length in Eq. 43 Bortolotti  _et al._ (2006);
Ołdziejewski and Jachymski (2016). In the case of the most magnetic of the
atoms, Dy, with $\varepsilon_{\rm{dd}}\sim 0.5-2$, corrections to $a_{\rm dd}$
were predicted to be of 10% at temperature of 100 nK and only 2% at 10 nK
Ołdziejewski and Jachymski (2017).

##### Mean-field approximation

The MF approximation, which underlies common approaches to describe the
properties of bulk BECs, consists of replacing the field operator by its mean
value, $\psi(\boldsymbol{r},t)$. The approximation assumes that the BEC mode
is macroscopically occupied Bogoliubov (1947); Gross (1961); Pitaevskii
(1961); Stringari and Pitaevskii (2016). This effectively ignores fluctuations
of the bosonic field around its mean value. In the absence of the DDI, the
corrections to this approximation scale as $\sqrt{na^{3}}$. In the presence of
the DDI, one naively expects the MF approximation to be valid when
$\sqrt{na^{3}},\,\sqrt{na_{\rm dd}\,^{3}}\,\ll 1$. A proper treatment of
fluctuations around MF values, outlined in Sec. V, leads to a more restrictive
condition: $\sqrt{na^{3}}\ll 1,\varepsilon_{\rm{dd}}\lesssim 1$. Note that the
length $a_{\rm dd}$ is defined such that $a=a_{\rm dd}$
($\varepsilon_{\rm{dd}}=1$) marks the limit of the mechanical stability of a
homogeneous isotropic dBEC, and other authors sometimes call $3a_{\rm dd}$ the
dipolar length (which arises from a consideration of the two-body problem Bohn
_et al._ (2009b)); see Sec.I.3.4.

Considering a gas of $N$ atoms interacting via Eq. (43), the MF approximation
results (at zero temperature) in a nonlocal Gross-Pitaevskii equation (GPE)
given in Eq. (9) and reprinted here:

$\displaystyle i\hbar\frac{\partial\psi}{\partial
t}(\boldsymbol{r},t)=\left[-\frac{\hbar^{2}}{2m}\nabla^{2}+V_{\rm{tr}}(\boldsymbol{r})+g|\psi|^{2}+\Phi_{\rm{dd}}(\boldsymbol{r})\right]\psi,$
(46)

with $g=4\pi\hbar^{2}a/m$ and $\Phi_{\rm{dd}}(\boldsymbol{r})=\int
d\boldsymbol{r^{\prime}}U_{\rm{dd}}(\boldsymbol{r^{\prime}}-\boldsymbol{r})|\psi(\boldsymbol{r^{\prime}})|^{2}$;
see also Eq. (10). The corresponding energy functional is

$\displaystyle E[\psi]=$
$\displaystyle\int\left[-\frac{\hbar^{2}}{2m}|\nabla\psi|^{2}+V_{\rm{tr}}(\boldsymbol{r})|\psi|^{2}\right.$
(47)
$\displaystyle\left.+\frac{g}{2}|\psi|^{4}+\frac{1}{2}|\psi|^{2}\Phi_{\rm{dd}}(\boldsymbol{r})\right]dr,$

which can be minimised to obtain the equilibrium ground state of the BEC
Lahaye  _et al._ (2009); Baranov  _et al._ (2002a, 2008, 2012). In these
equations, the DDI yields a nonlocal contribution
$\Phi_{\rm{dd}}(\boldsymbol{r})$, whose sign depends on the shape of the
density distribution. We note that various analyses demonstrated the
applicability of the nonlocal GPE in accurately describing the properties of
trapped dipolar gases Ronen  _et al._ (2006b); Bortolotti  _et al._ (2006).

####  IV.1.2 Magnetostriction

##### Magnetostriction of dBECs

When placed in a cylindrical harmonic trap, a nondipolar gas
($\varepsilon_{\rm{dd}}=0$) assumes the geometry of the trap. A dipolar system
behaves differently because of the presence of $\Phi_{\rm{dd}}$ in Eq. (9).
Specifically, the dipole term $\Phi_{\rm{dd}}(\boldsymbol{r})$ has a saddle-
shape, tending to elongate the distribution along the magnetic field. That is,
to minimise its energy, the dBEC elongates along the dipolar axis due to the
attractive part of the DDI. This effect is called magnetostriction and is of
course not unique to quantum gases as it has been known since Joule discovered
the effect in iron Joule (1847). Note that the magneto- or electro-strictive
effects deform a system solely due to the presence of a homogeneous field
which breaks rotational symmetry. The forces responsible are due only to the
interactions between the particles because the external field has no gradient.

$\bullet$ Application of the Thomas-Fermi approximation:

Rather than numerically solving Eq. (9), magnetostriction can be analysed by
the use of the Thomas-Fermi (TF) approximation, valid for high atom numbers
and sufficiently large total interaction strength Stringari and Pitaevskii
(2016). This allows one to neglect the kinetic energy term in Eq. (9). For
$\varepsilon_{\rm{dd}}=0$ (i.e., a nondipolar BEC), one can easily show that
the density of a BEC acquires an inverted parabolic shape with a so-called TF
radius that scales as $R_{i}\propto 1/\omega_{i}$, $i=x,y,z$ Pitaevskii and
Stringari (2016). In the dipolar case, Refs. O’Dell  _et al._ (2004); Eberlein
_et al._ (2005) showed that a dBEC retains the parabolic shape in this (TF)
limit. However, its aspect ratio

$\displaystyle\kappa=R_{\rho}/R_{z},$ (48)

defined by the ratio of the TF radii perpendicular ($R_{\rho}$) to the radius
along ($R_{z}$) the dipole direction, no longer reflects the aspect ratio of
the trap. This is due to the contribution of $\Phi_{\rm{dd}}$ to the mean
energy in Eq. (47), which one can calculate analytically:

$\displaystyle
E_{\rm{dd}}=-\frac{4}{7}\varepsilon_{\rm{dd}}\,g\,\frac{n_{0}\,N}{2}f_{\rm{dip}}(\kappa),$
(49)

where $n_{0}$ is the density at the trap centre. The function
$f_{\rm{dip}}(\kappa)$ encompasses the anisotropy of the interaction energy,
it is a decreasing function of $\kappa$ 121212Expression (49) is valid for an
inverted parabola solution of the TF approximation. If one uses a Gaussian
variational ansatz for the wavefunction of the BEC, it becomes
$E_{\rm{dd}}=-\varepsilon_{\rm{dd}}\,g\,\frac{n_{0}\,N}{2}f_{\rm{dip}}(\kappa)$..
The complete expression is Giovanazzi  _et al._ (2003)

$f_{\rm{dip}}(\kappa)=\frac{1+2\kappa^{2}}{1-\kappa^{2}}-\frac{3\kappa^{2}\rm{arctanh}\sqrt{1-\kappa^{2}}}{(1-\kappa^{2})^{3/2}}.$
(50)

Note that an angular integration of the DDI over an isotropic distribution is
zero, leading to $f_{\rm{dip}}(1)=0$. $f_{\rm{dip}}$ is bounded by two limits:
fully collinear dipoles (attractive DDI) $\lim_{\kappa\to 0}f_{\rm{dip}}=1$,
while side-by-side dipoles (repulsive DDI) yield
$\lim_{\kappa\to\infty}f_{\rm{dip}}=-2$. Thus, one readily sees that
$E_{\rm{dd}}$ is reduced by lowering $\kappa$. In the TF approximation,
minimising $E$ leads to a transcendental equation for $\kappa$ Yi and You
(2001); Giovanazzi  _et al._ (2003); Eberlein  _et al._ (2005):

$\displaystyle\lambda=\kappa\left(\frac{1+2\varepsilon_{\rm{dd}}-\frac{3\varepsilon_{\rm{dd}}f_{\rm{dip}}(\kappa)}{1-\kappa^{2}}}{1-\varepsilon_{\rm{dd}}+\frac{\kappa^{2}}{2}\frac{3\varepsilon_{\rm{dd}}f_{\rm{dip}}(\kappa)}{1-\kappa^{2}}}\right)^{1/2}.$
(51)

The solution of this equation for a given $\lambda$ thus gives the degree of
magnetostriction. A solution exists for any $\lambda$ provided that
$\varepsilon_{\rm{dd}}\leqslant 1$. The absence of a MF solution will be
discussed at the end of this section, and for now we assume there exists one.

We note that assuming a 3D Gaussian shape for the dBEC is another common
approximation Eberlein  _et al._ (2005); Yi and You (2000, 2001); Santos  _et
al._ (2000); Góral and Santos (2002). This approximation will be later
reviewed in Sec. IV.1.4. Interestingly, and despite the difference in the
approximations, such a Gaussian variational ansatz for $\psi$ yields the same
equation as Eq. (51) for the aspect ratio. We note that both the TF and the
Gaussian approximations can describe a regular BEC wavefunction with only
maximal density at the trap centre. They fail in describing more complex
phenomena that can occur in dBECs, e.g., when operating close to the MF
instability threshold; see Secs. IV.1.4 and V.4. In this paragraph, we focus
on the shape of stable dBECs.

$\bullet$ In-situ measurements:

We now describe the distortion of the gas shape inside the trap. Quantitative
in-situ measurements of magnetostriction have only recently been possible: The
first in-situ images of dipolar BECs were reported in Ref.Kadau  _et al._
(2016) and in-situ magnetostriction images first presented in Ref. Wenzel  _et
al._ (2018a). An example is shown in Fig. 10.

Figure 10: In situ imaging of magnetostriction of a 162Dy dBEC. The atoms are
held in a pancake-shaped trap, nearly isotropic in the plane. The field is out
of plane in the left image, and in-plane for the right image, while the trap
is kept unchanged. Adapted from Ref. Wenzel  _et al._ (2018a).

$\bullet$ Magnetostriction in the TOF dynamics of BECs:

Observing the (hydrodynamic) expansion of the gas trap release is a common way
to reveal BEC physics Holland  _et al._ (1997); O’Hara  _et al._ (2002), and
in particular the impact of interactions on quantum gases. The hydrodynamic
equations for the expansion of a (nondipolar) BEC were first presented in Ref.
Castin and Dum (1996), showing inversion of ellipticity during TOF. Solutions
to the hydrodynamic equations for dipolar gases can also be found O’Dell  _et
al._ (2004); Eberlein  _et al._ (2005); Giovanazzi  _et al._ (2006a). The
effect of the DDI is the same as in trap: energetics favour a distortion in
the expansion that aligns the dipoles head-to-tail.

The first observation of magnetostriction in a gas was observed in the
hydrodynamic expansion of a 52Cr BEC Stuhler  _et al._ (2005). This was
manifest as a dependence of the inversion of ellipticity on the magnetic field
direction. Comparing the experimental dynamics with solutions to the
hydrodynamic equations, a value for $\varepsilon_{\rm{dd}}$ for 52Cr away from
FRs was extracted Griesmaier  _et al._ (2006). The value agrees with most
precise values of the background scattering length $a_{\rm{bg}}$ obtained more
recently Pasquiou  _et al._ (2010); see Sec. II.4. Using a FR to lower $a$—and
thus to enhance dipolar effects into the $\varepsilon_{\rm{dd}}>1$
regime—magnetostriction in situ and in TOF was strong enough so that a
complete suppression of ellipticity inversion was observed Lahaye  _et al._
(2007). The variation of the aspect ratio of an expanding cloud with
$\varepsilon_{\rm{dd}}$ can be seen in Fig. 11.

Figure 11: Aspect ratio of an expanding 52Cr BEC, measured 2 ms after its
release from a trap, as a function of $\varepsilon_{\rm{dd}}$.
$\varepsilon_{\rm{dd}}$ was controlled using a a FR centred about $B_{0}$.
Blue (red) dots: data with $B<B_{0}$ ($B>B_{0}$). Dashed line: background
$\varepsilon_{\rm{dd}}$. Solid line: prediction from Eq. (49) without
adjustable parameters. Shaded area: uncertainties in $\omega_{\rho,z}$. Inset,
same data as a function of $B-B_{0}$. From Lahaye  _et al._ (2007).

##### Magnetostriction and TOF dynamics of ultracold, non-condensed gases

Observation of magnetostriction in gases requires low temperatures and high
densities. Otherwise, thermal energies are orders-of-magnitude larger than the
DDI. As mentioned above, the use of dBECs enabled the observation of
magnetostriction both in TOF and in situ. Magnetostriction effects could also
be observed in ultracold, but non-degenerate thermal gases of sufficiently
strong dipolar atoms.

The dimensionless parameter relevant for DDI corrections to the ideal gas law
is obtained is obtained by comparing the mean dipolar energy to the mean
kinetic energy: $\eta=E_{\rm dip}/E_{\rm k}$. For a non-degenerate gas,
$E_{\rm dip}\sim S^{2}d^{2}n$ and $E_{\rm k}\sim k_{\rm B}T$, such that
$\eta\propto(n\lambda_{\rm th}^{3})\times(k_{\rm th}a_{\rm dd})$, with $k_{\rm
th}=2\pi/\lambda_{\rm th}=\sqrt{mk_{\rm B}T}/\hbar$ the thermal wavenumber,
and $\lambda_{\rm th}$ the thermal wavelength. To maximise $\eta$ in a thermal
gas, we may choose Dy, which provides a $a_{\rm dd}=131\,a_{0}$; see Sec. II).
At typical experimental densities for a gas with a temperature just above its
BEC $T_{c}$, $\eta$ remains $\ll 1$. Thus, at best, the DDI only weakly
modifies the TOF expansion dynamics. Nevertheless, DDI effects on time-of
flight dynamics have been observed, as we now discuss.

Figure 12:  The anisotropic aspect ratio of dipolar thermal gases after 16 ms
of TOF. Aspect ratio for (a) 162Dy and (b) 164Dy. Red (blue) is for magnetic
field along $\hat{z}$ and $\hat{y}$. The solid (dashed) curves are the full
(partial, MF-only) theory developed in Ref. Tang  _et al._ (2016b). The Bose-
condensation temperature is marked by the vertical grey line. From Ref. Tang
_et al._ (2016b).

Tang et al. Tang  _et al._ (2016b) theoretically and experimentally studied
the anisotropic expansion of a thermal dipolar Bose gases of 164Dy and 162Dy
just above their degeneracy temperature. Each gas, after TOF expansion,
exhibits an aspect ratio that depends on the polarisation angle of the
dipoles, as had already been noted in dBECs Lahaye  _et al._ (2009); see Fig.
12. The DDI strength of Dy is sufficient to cause the magnetostriction of even
a dilute thermal gas. To predict the experimental aspect ratios versus
temperature and dipolar angle, the authors developed a theory of the expansion
that accounts for the Hartree-Fock MF interactions, Bose-enhanced scattering,
and hydrodynamic effects that partially cross-dimensionally rethermalise the
gas during the expansion. By doing so, the authors were able to quantitatively
match the theory to the data. The theory fits provided a method to extract the
gas temperature and is a relatively simple method for determining the
scattering length of the gas, even near a FR; see Sec. II.4.2 for the results
of the $a$ measurements using this technique, which are compatible with
previous measurements. Moreover, the momentum distribution deformation scales
with the ratio $\eta=(n\lambda_{\rm th}^{3})\times(k_{\rm th}a_{\rm dd})$, as
expected, and arises primarily from the two-body collisions and the direct
Hartree contribution to the MF energy.

A similar effect arises in degenerate Fermi gases, although with a few
fundamental differences. The kinetic energy is dominated by the Fermi energy,
so the distortion in momentum space is smaller. This Fermi surface deformation
stems from the exchange contribution to the mean interaction energy, which
does not vanish in the degenerate Fermi many-body state due to the Pauli
exclusion principle; see Sec. I.4.2 and the discussion in Sec. IV.2.

####  IV.1.3 Elementary excitations of dBECs

In the previous section, we were concerned with the influence of the DDI at
equilibrium: i.e., the magnetostrictive effect on the gas wavefunction in real
space, and its consequences on the free-space expansion dynamics, which allows
observations of momentum-space magnetostriction. In contrast, elementary
excitations around equilibrium provide a window into dynamical behaviour of
the dipolar quantum gases. Due to the sensitivity of spectroscopic
measurements, experimental studies of DDI effects on elementary excitations
even at low $\varepsilon_{\rm{dd}}$ are possible. Besides, stringent
modifications of the excitation spectrum at $\varepsilon_{\rm{dd}}\sim 1$ were
also evidenced.

##### Elementary excitations in a (homogeneous) dBEC: from phonons to free-
particles.

$\bullet$ Continuous Bogoliubov spectrum:

Elementary excitations of BECs were introduced in Sec. I.4.1 in the uniform
BEC case. In presence of DDI, the well-known Bogoliubov spectrum, obtained by
linearising the GPE around the ground state, with its linear phonon branch
followed by the quadratic free-particle dependence, is modified. By combining
Eqs. (8) and (16), the spectrum of a uniform 3D dBEC of density $n$ is Santos
_et al._ (2000)

$\hbar\,\omega(\boldsymbol{k})=\sqrt{\frac{\hbar^{2}k^{2}}{2m}\left(\frac{\hbar^{2}k^{2}}{2m}+2\,gn\,\left(1+\varepsilon_{\rm{dd}}(3\cos^{2}\theta_{k}-1)\right)\right)},$
(52)

where $\theta_{k}$ is the angle between the direction of excitation
propagation and the dipole orientation. The sound velocity for BECs with only
contact interactions ($\varepsilon_{\rm{dd}}=0$) is $c_{0}=\sqrt{gn/m}$ and
acquires an angular dependence in presence of the DDI,
$c(\theta_{k})=c_{0}\sqrt{1+\varepsilon_{\rm{dd}}(3\cos^{2}\theta_{k}-1)}$.
Excitations propagating along the dipoles’ direction ($\theta_{k}=0$) have the
highest sound velocity (i.e., they are stiff modes), while those propagating
in the perpendicular plane ($\theta_{k}=\pi/2$) have the lowest velocity
(i.e., they are softer modes). This anisotropy of the dispersion relation is
one of the main feature of dBECs with respect to collective excitations and
has important physical consequences; namely, the anisotropic superfluid
critical velocity (see below) and a mechanical instability with anisotropic
nature manifest for gases with $\varepsilon_{\rm{dd}}\geq 1$; see Sec.IV.1.4.

$\bullet$ Bragg spectroscopy:

Using the standard method of two-photon Bragg spectroscopy (see, e.g., Refs
Stenger  _et al._ (1999a); Stamper-Kurn  _et al._ (1999); Steinhauer  _et al._
(2002); Vogels  _et al._ (2002); Ozeri  _et al._ (2002)) on a dBEC of 52Cr,
Bismut et al. Bismut  _et al._ (2012) probed the dispersion relation of such a
trapped dBEC away from FRs. In practice, the finite size of the trapped
samples sets a lower limit on the momentum at which one can probe and compare
against Eq. (52) (which holds in the homogeneous case). This limit comes from
the fact that the BEC smallest size ($R_{\rm TF\,min}$) must be larger than
the excitation wavelength: $k\,R_{\rm TF\,min}\gg 1$. In this limit, the local
density approximation can be used, allowing one to take into account the BEC’s
inhomogeneity. Bismut et al. observed a Bragg dispersion peak that depends on
the relative orientation between the Bragg momentum and the magnetic field
angle, in very good agreement with theory; see Fig. 13. At low momentum, a
departure from the homogeneous expectation was observed. Numerical simulations
taking finite-size effects into account reproduce the evolution of the
measured dispersion relation and its anisotropy.

Figure 13: Bragg excitation spectra of a 52Cr BEC, with the magnetic field in
two perpendicular directions. Lines are fit to the data. From Bismut  _et al._
(2012).

$\bullet$Critical superfluid velocity measurement:

The DDI effect on the dBEC’s excitation spectrum implies other changes to the
BEC’s physical behaviour. A prime example lies in the change to the superfluid
critical velocity. The famous Landau’s criterion relates the critical velocity
for an impurity moving in the direction $\hat{\boldsymbol{v}}$ in an (infinite
and uniform) superfluid to its dispersion relation: $v_{\rm
c}=\underset{\boldsymbol{k}}{{\rm
min}}\,(\frac{\omega_{\boldsymbol{k}}}{\hat{\boldsymbol{v}}\cdot\boldsymbol{k}})$
Landau (1941); Yu (2017). In the homogeneous case discussed above, Eq. (52)
implies that the critical velocity becomes anisotropic in a dBEC, depending on
whether the excitation is applied along or perpendicular to the dipole
orientation Ticknor  _et al._ (2011); Yu (2017). We note that the critical
velocity does not generally match the speed of sound, even in the homogeneous
case. This is because the dissipating excitations can occur in a different
direction than the impurity’s motion Yu (2017). The critical velocity is thus
systematically smaller than the speed of sound in the direction of motion. For
anisotropic confinement, as discussed in Sec. IV.1.3, the critical velocity
should also be affected by the existence of low-energy high-momentum modes,
such as the roton mode Wilson  _et al._ (2010); Ticknor  _et al._ (2011).

Superfluid velocity measurements were performed on a 162Dy dBEC using a local
light defect driven linearly along one axis Wenzel  _et al._ (2018b). These
showed that the anisotropy of the critical velocity as well as of the heating
rate above the critical velocity are in excellent agreement with dynamical
simulations based on the GPE; see Fig. 14. Large corrections compared to the
homogeneous predictions are observed and mainly attributed to inhomogeneous
density effects, as corroborated by numerical simulations of the GPE.

Figure 14: Temperature of the 162Dy dBEC after applying a stirring protocol.
The stirring beam is moved along $x$ (red squares) or $y$ (blue circles). The
dipoles are polarised along $x$. The trap is cylindrically symmetric in $x$
and either (c) elongated or (b) made narrow in this direction; see insets with
example of in situ images. A piecewise linear fit (dashed lines) quantifies
the anisotropy of the critical velocities (arrows) with respect to the dipole
orientation. For both trap geometries, the velocity, $v_{\perp}$, for a
perpendicular excitation is small than for a parallel one, $v_{||}$, with
$v_{\perp}=0.16(2)$mm/s ($0.12(3)$mm/s) and $v_{||}=0.36(3)$mm/s
($0.26(4)$mm/s) in b (c). Numerical simulations of the GPE (solid lines) show
excellent agreement with the experiment. Adapted from Wenzel  _et al._
(2018a).

##### Low-lying collective modes of trapped dBECs and oscillation measurements

For a trapped gas, if one lowers the excitation momentum down to the regime
where the corresponding wavelength is on the order of the cloud size, then one
reaches the regime of low-lying excitations where the spectrum is discrete and
momentum is not a good quantum number. The low-lying modes are typically
surface modes, implying an (out-of-phase) oscillation of radii in different
directions. These modes typically have a compressional character, yielding a
density (and thus interaction) dependence. Signatures of the DDI in these
systems have been theoretically studied Yi and You (2001, 2002); Góral and
Santos (2002); O’Dell  _et al._ (2004); Giovanazzi  _et al._ (2007). MF
methods to extract the collective modes nature and frequency rely on expanding
the energy functional Eq. (47), or corresponding hydrodynamic equations,
around the stationary solution. This can be done semi-analytically by either
applying the TF approximation or using the Gaussian ansatz; see Secs. IV.1.2
and IV.1.4.

Experimentally, Bismut et al. Bismut  _et al._ (2010) studied the influence of
the DDI on such modes with a 52Cr BEC. They investigated the second-lowest-
lying mode, a quadrupole mode in a non-axisymmetric trap. The experimental
results demonstrate that the collective mode frequency is dependent on the
relative orientation of the dipoles with the trap axes, and that the frequency
shift is dependent on the trap geometry, see Fig. 15. The experiment is in
good agreement with a TF approximation theory, which neglects the kinetic
energy of the atoms. When lowering the atom number below a few thousand, the
frequency shift is clearly reduced, demonstrating the importance of quantum
pressure for very small samples. Quantum pressure can be taken into account by
either performing full numerical simulations of the GPE, or by using a
Gaussian variational ansatz for the BEC wavefunction.

Figure 15: Time evolution of the aspect ratio of a dBEC oscillating due to the
excitation of a low-lying surface collective mode (taken after TOF expansion).
The two different colours show the time evolution for two different angles of
the $B$ field, showing the dependence of the frequency of the collective modes
on the dipoles direction due to the DDI. From Bismut  _et al._ (2010).

##### Spectrum of elementary excitations in trapped dBEC, roton mode.

As highlighted in Eq. (16), the DDI contributes to the continuous excitation
spectrum via an effective coupling strength
$\widetilde{V_{\rm{dd}}}(\boldsymbol{k})$, which adds to the constant contact
contribution $g$. As observed in Eq. (52), in the 3D homogeneous case, the DDI
contribution, simply given by its Fourier transform (Eq. (8)), yields an
orientation-dependence, but no $k$-norm dependence. When combining the DDI
with an anisotropic confinement, even more striking behaviours of the dBEC’s
excitation spectrum arise. Here, a $k$-norm dependence may arise from the
interplay of the DDI with the trap’s natural length scales. This $k$
dependence yields qualitative differences in the dispersion relation compared
to that of a non-dipolar BEC and of a uniform dBEC. In this section, we keep
the notation of $\widetilde{V_{\rm{dd}}}(\boldsymbol{k})$ for the DDI
effective coupling strength, which determines the excitation dispersion
relation, even if the gas in non-uniform and this differs from the Fourier
transform of Eq. (8).

A particularly interesting case is when the dispersion relation $\epsilon(k)$
of a dBEC becomes non-monotonic, presenting a local maximum (maxon) followed
by a local minimum (roton). Such a qualitative change in the dispersion
relation underlies important new physical behaviour. To get a first insight
into such changes, it is interesting to highlight that such dispersion
relation resembles the celebrated dispersion relation of superfluid helium
Landau (1941, 1947); Henshaw and Woods (1961). Here, roton excitation were
first speculated to explain the exotic macroscopic properties of this
superfluid Landau (1941, 1947); Feynman (1955), long before their observation
Henshaw and Woods (1961). Thanks to its low energy, the roton excitation
strongly influences the response of the superfluid to small excitations.
Furthermore, because of its large momentum, the roton underlies the tendency
of the fluid to crystallise (at the wavelength corresponding to its inverse
momentum) Nozières (2004); Kirzhnits and Nepomnyashchii (1971); Schneider and
Enz (1971) (although one should note that the phase transition to solid helium
is not due to Roton softening).

$\bullet$ Roton excitation spectrum in anisotropic semi-infinite dBECs:

In 2003, a dispersion relation presented a roton minimum was predicted to
occur in anisotropically trapped weakly interacting dipolar gases, first in
the context of light-induced DDI by O’Dell et al. O’Dell  _et al._ (2003) and,
shortly after, in that of magnetic or electric dipolar system by Santos et al.
Santos  _et al._ (2003). These seminal works consider semi-infinite trapping
geometries, i.e., infinite along one ($\omega_{z}=0$) O’Dell  _et al._ (2003)
or two ($\omega_{\rho}=0$) Santos  _et al._ (2003) directions of space, and
harmonically confined along the others. This treatment allows to account for
anisotropy effects while facilitating theoretical treatment, yielding semi-
analytical expressions of $\widetilde{V_{\rm{dd}}}(\boldsymbol{k})$ within the
TF approximation, and providing an intuitive picture of the effect.

The occurrence of a roton minimum arises from
$\widetilde{V_{\rm{dd}}}(\boldsymbol{k})$ becoming attractive at large $k$. In
the quasi-infinite geometries of Refs. Santos  _et al._ (2003); O’Dell  _et
al._ (2003), the confinement acts to limit the attractive contribution of the
DDI: The attraction dominates over the repulsive contribution only if the
momenta (along an unconfined direction) have a norm larger than the inverse
characteristic confinement length $\ell_{z}$. In this way, for excitations
along the unconfined directions, $k\ll 1/\ell_{z}$ yields
$\widetilde{V_{\rm{dd}}}(\boldsymbol{k})>0$, while, for $k\gtrsim 1/\ell_{z}$,
$\widetilde{V_{\rm{dd}}}(\boldsymbol{k})<0$. Therefore, the DDI stiffens the
dispersion relation in the phononic regime while it bends it down for large
$k$. Because of the additional contribution of kinetic term
$\frac{\hbar^{2}k^{2}}{2m}$ which ultimately dominates at very large $k$, the
effect of $\widetilde{V_{\rm{dd}}}(\boldsymbol{k})<0$ is the strongest at
$k\sim 1/\ell_{z}$ and, for weak enough $s$-wave coupling strength $g$, a
minimum arises at $k_{\rm{rot}}\sim 1/\ell_{z}$, matching a roton mode. Here
the roton wavelength is typically set by the confinement along the direction
of attractive DDI. Furthermore, the roton energy $\Delta$ can be lowered by
increasing the density or increasing $\varepsilon_{\rm{dd}}$. When $\Delta=0$,
Ref. Santos  _et al._ (2003) finds $k_{\rm{rot}}=\sqrt{2}/\ell_{z}$,
independent of the density and the interaction parameters
$g,\,\varepsilon_{\rm{dd}}$.

$\bullet$ Roton excitation spectrum in anisotropic fully trapped dBECs:

Numerous subsequent theoretical works describe the roton in dipolar gases
confined in finite geometries (Eq. (44) with $\lambda\gg 1$). They study how
the roton spectrum of the infinitely elongated geometries described above
survive in the fully trapped case and how the steady state and dynamical
behaviour of dBECs are affected Ronen  _et al._ (2007); Martin and Blakie
(2012); Blakie  _et al._ (2012); Jona-Lasinio  _et al._ (2013a); Bisset and
Blakie (2013); Bohn  _et al._ (2009a); Parker  _et al._ (2009); Wilson  _et
al._ (2008); Nath and Santos (2010); Wilson  _et al._ (2010); Klawunn  _et
al._ (2011); Bisset  _et al._ (2013); Blakie  _et al._ (2013); Jona-Lasinio
_et al._ (2013b); Natu  _et al._ (2014). In the case where the trap is tightly
confining along the dipoles and sufficiently anisotropic (e.g., $\lambda\gg
1$), a roton-like spectrum arises for large enough $\varepsilon_{\rm{dd}}$.
The relevant range of $\varepsilon_{\rm{dd}}$ depends on the exact trap
geometry. While, in the inhomogeneous case, the momentum is not longer a good
quantum number, the rotonic properties of the spectrum are most directly
revealed by the behaviour of the dynamic structure factor. This factor
$S(k,\omega)$ conveniently describes the system’s response to a Bragg
excitation at momentum $k$ and frequency $\omega$ Blakie  _et al._ (2012).
Furthermore, in harmonically trapped geometries, the elementary mode
contributing to the roton minimum in the dynamic structure factor has an
amplitude that vanishes away from the dBEC centre Jona-Lasinio  _et al._
(2013a). This confinement effect can be understood, within the local density
approximation, as the roton energy decreasing with increasing density. This
translates into a finger-like feature in the (discrete) dispersion relation
Bisset  _et al._ (2013). Other signatures of the existence of the roton mode
in finite dBECs include anomalously large density fluctuations at the roton
wavelength in an equilibrium gas Klawunn  _et al._ (2011); Bisset and Blakie
(2013); Blakie  _et al._ (2013), strong correlation emerging from an
interaction quench Natu  _et al._ (2014), and peculiar structures in the
collapse dynamics Parker  _et al._ (2009); Bohn  _et al._ (2009a); see also
Sec. IV.1.4.

When the trap anisotropy is reduced (e.g., $\lambda\gtrsim 1$) and
$R_{\perp}\gtrsim R_{z}$, the existence of a roton feature, occurring at large
$k$ compared to $1/R_{\perp}$, becomes progressively lost. Indeed, the
argument of a finite $k$ scale for the activation of the attractions starts to
break down. Yet the interplay between DDI and geometry persists, now including
also finite-size effects. This yields other interesting features in the
excitation spectrum, in particular the so-called angular roton mode that has a
momentum $k\sim 1/R_{\perp}$, but a non-zero angular momentum Ronen  _et al._
(2007); Wilson  _et al._ (2009); see also discussion in Sec. IV.1.4. Finally,
we note that while most of the theoretical work has focused on cylindrically
symmetric geometries about the dipole axis, extensions to transversely
anisotropic geometries have also been discussed Martin and Blakie (2012).

Figure 16: Roton excitation spectrum measured with a 166Er dBEC of cigar shape
(cigar along $y$, dipoles along $z$). (a) Measured dispersion relation (dots)
compared to predicted response for the Bragg measurement (colorcode) at
$a=52.5a_{0}$ [theory at $51.6a_{0}$]. (b) Measured $\epsilon(k_{\rm{rot}})$
vs. scattering length (circles) and corresponding theory prediction. The
shadings show the calculations over the prediction interval of $a$. The dashed
line shows the power-law fit to the experiment. Inset: zoom-in around the
instability threshold comparing $\epsilon(k_{\rm{rot}})$ (circles) to
$\epsilon(k_{\rm max})$ of the maxon [local maximum] (triangles), respectively
with the power-law fit and a guide to the eye (dashed lines). The region where
$0\leq\epsilon(k_{\rm{rot}})\leq\epsilon(k_{\rm max})$ is highlighted with a
white background. Adapted from Petter  _et al._ (2019).

$\bullet$ Experimental measurements of the roton spectrum:

The roton mode has remained elusive in Cr experiments due to the weak dipolar
character of this species. That weakness makes the $a$-range of existence of a
potential roton minimum very narrow, as well introduces strong three-body
losses when tuning within this small $a$ – see also Sec. V. This leads to a
fast change of the BEC’s properties, including its excitation spectrum. The
roton mode was experimentally observed in a 166Er BEC Chomaz  _et al._ (2018);
Petter  _et al._ (2019) using a cigar-shaped trap geometry with the dipoles
aligned along a tightly confining direction. This geometry, simplifying
experimental observations, differs from most of the early theoretical works on
roton excitations, which rely on cylindrical symmetry around the dipole axis.
The observation of Ref. Chomaz  _et al._ (2018), relying on instability
dynamics, will be described in Sec. IV.1.4. In Ref. Petter  _et al._ (2019),
Petter et al. Petter  _et al._ (2019) reported on roton spectrum measurements
based on Bragg spectroscopy, similar to Sec. IV.1.3. The scattering length was
tuned close to instability and the excited momentum, along the trap axis,
varied from $k\ll 1/\ell_{z}$ to $k\sim 2/\ell_{z}$, thus probing the full
phonon-maxon-roton dispersion relation; see Fig. 16(a). When increasing
$\varepsilon_{\rm{dd}}$, a preferential softening of $\epsilon(k)$ at large
$k$ is observed, finally forming a minimum at $k=k_{\rm{rot}}\sim
1.3/\ell_{z}$; see Fig. 16(a). The minimum is observed only in a narrow range
of scattering lengths and the roton gap shows a fast decrease with
$\varepsilon_{\rm{dd}}$ toward instability ($\epsilon(k_{\rm{rot}})\approx
0$), see Fig. 16(b). The Bragg measurements also proved, through the increase
of the Bragg response, the enhancement of density-density correlations at
$k_{\rm{rot}}$ as the roton instability is approached, indicative of the
system tendency to crystallise. More recently, experiments performed on Dy
BECs revealed roton excitations by a direct analysis of the in-situ density
fluctuations Hertkorn  _et al._ (2021a), following early theoretical studies
Klawunn  _et al._ (2011); Bisset and Blakie (2013); Blakie  _et al._ (2013).
Fluctuation analysis provides access to the static structure factor. It also
allowed to identify the two degenerate roton modes, which correspond to
symmetric and antisymmetric density patterns to the trap centre. The study of
roton excitations via the induced density fluctuations was extended beyond the
cigar-shaped case, to oblate systems, in an independent set of experiments
Schmidt  _et al._ (2021). Here the the distinct softenings of several radial
and angular rotons were revealed. As introduced in the previous paragraph,
radial rotons correspond to the standard situation of a mode of large radial
momentum $\langle k\rangle$ and radial symmetry, angular rotons distinctly
have low $\langle k\rangle$ but large angular momentum and show azimuthal
patterns Wilson  _et al._ (2009); Bisset  _et al._ (2013) and reveal as such
in the density fluctuations.

####  IV.1.4 Mean-field stability of dBECs

Beyond the properties of the ground state and its collective excitations in a
stable regime, the DDI also impacts the very stability of the state, in
particular making it geometry dependent. In this section, we review in detail
the effect of the DDI on the stability of the quantum gases: first, how the
DDI makes stability depends on the trap geometry; second, how it affects the
stability of a cloud via long-range interactions between neighbouring clouds;
and finally how the DDI introduces, beyond the global stability condition, a
distinct (local) collapse mechanism in some special geometries. We will see
how the later case relates to the softening of excited modes distinct from the
lowest lying ones, arising from the interplay of DDI and either anisotropic
geometries or, additionally, finite-size effects.

##### Global mechanical stability for a single BEC

Three-dimensional, homogeneous BECs under contact interactions are
mechanically stable for positive compressibility
131313$\chi=-\frac{1}{V}\frac{\partial V}{\partial
P}=\frac{1}{n}\frac{\partial n}{\partial P}=\frac{1}{n^{2}}\frac{\partial
n}{\partial\mu}$ $\chi=\frac{1}{nmc_{0}^{2}}>0$ (where $c_{0}$ is the speed of
sound; see Sec. IV.1.3), they are thus unstable for $a<0$ Stringari and
Pitaevskii (2016); Kagan  _et al._ (1998). Finite-sized Eigen  _et al._ (2016)
and harmonically trapped Bradley  _et al._ (1997); Roberts  _et al._ (2001)
BECs can be stabilised by quantum pressure for small negative scattering
length. As dBECs experience competing interactions, their mechanical stability
criterion is more complex. Even at positive scattering length, this mechanical
instability occurs for $\varepsilon_{\rm{dd}}=1$ as the sound velocity
perpendicular to the dipole orientation then cancels,
$c^{2}|_{\theta_{k}=\pi/2}=0$; see Sec. IV.1.3. However, this holds for only
infinite, homogeneous, and isotropic BECs. Because of the anisotropic
character of the DDI, the anisotropy of real samples must be considered to
understand their stability.

$\bullet$ Stability within the TF approximation:

Global mechanical collapse can be alternatively understood from an energy
argument. The energy density of MF interactions, when attractive, scales like
$E_{\rm MF}\sim-n^{2}$ and is thus minimised for infinite density. Using the
TF approximation and neglecting the kinetic energy contribution, $E_{\rm
MF}<0$ leads to a singularity of the ground state density, thus giving an
instabibility. From Eq. (49), the sign of the MF dipolar energy of a
cylindrical dBEC depends on the cloud aspect ratio $\kappa$. Taking into
account both interactions, one can show that the total MF interaction energy
scales like

$E_{\rm MF}\sim gn_{0}N(1-\varepsilon_{\rm{dd}}f_{\rm{dip}}(\kappa)).$ (53)

Imposing an attractive MF interaction, $E_{\rm MF}<0$, leads to the
instability condition $a<a_{\rm dd}f_{\rm{dip}}(\kappa)$. Solving for
$\kappa(\lambda)$ through (48) leads to a stability diagram as a function of
$a$ and $\lambda$. Intuitively, for very prolate traps, $\lambda\ll 1$ and the
dipoles are aligned head-to-tails $\kappa\ll 1$ while the DDI is mostly
attractive ($f_{\rm{dip}}(\kappa)\approx 1$). Thus, $a$ must surpass $a_{\rm
dd}$ for stability. For the opposite case of very oblate traps ($\lambda\gg
1$), dipoles repel each other twice more strongly
$f_{\rm{dip}}(\kappa)\approx-2$, and $a$ must be large and negative to reach
an effective attractive MF energy, $a<-2\,a_{\rm dd}$. This simple reasoning
that considers only the MF interaction leads to a prediction that well
describes experimental observations Koch  _et al._ (2008); see Fig. 17(a).

Figure 17: Stability diagram of a 52Cr dBEC in a cylindrical trap of aspect
ratio $\lambda$ versus $a$. (a) Taken from Koch  _et al._ (2008). The thin
grey line is the result of the condition $a<a_{\rm dd}f_{\rm{dip}}(\kappa)$,
while the green line is the result of the Gaussian ansatz analysis. The red
line is for $\varepsilon_{\rm{dd}}=0$. In the large N limit this diagram is
universal, i.e., it also holds for electric dipoles. (b) Same data as (a), now
compared to either a Gaussian ansatz (dashed line) or a numerical solution,
minimising Eq. (49) (solid line). From Wilson  _et al._ (2009)

$\bullet$ The Gaussian ansatz for a BEC near collapse:

However, near collapse, kinetic energy effects also play a role, especially
for low atom numbers. One must depart from the TF approximation, keeping all
terms in Eq. (47). The variational Gaussian ansatz provides a semi-analytical
approach

$\psi(\boldsymbol{r})=\frac{N}{\pi^{3/4}\bar{\sigma}^{3}}e^{-\sum_{j=(x,y,z)}\frac{i^{2}}{2\sigma_{j}^{2}}},$
(54)

inserted into Eq. (47), the energy is minimised with respect to the
variational parameters $\sigma_{x,y,z}$, describing the BEC sizes along the
trap axes $x,y,z$. Here,
$\bar{\sigma}=(\sigma_{x}\sigma_{y}\sigma_{z})^{1/3}$. Once the integral is
carried out, the different contributions read:

$E_{\rm k}=N\,\frac{\hbar^{2}}{4\,m}\sum_{i}\frac{1}{\sigma_{i}^{2}}$ (55)

for the kinetic energy,

$E_{\rm trap}=N\,\frac{m}{4}\sum_{i}\omega_{i}^{2}\sigma_{i}^{2}$ (56)

for the external trapping energy, and the MF interaction energy reads

$E_{\rm
MF}=N^{2}\,\frac{g}{2(2\pi)^{3/2}\bar{\sigma}^{3}}(1-\varepsilon_{\rm{dd}}f_{\rm{dip}}(\kappa_{x},\kappa_{y}))$
(57)

where we allowed non-axisymmetric traps. The expression for the generalised
$f_{\rm{dip}}$ function in terms of the two aspect ratios
$\kappa_{x,y}=\sigma_{x,y}/\sigma_{z}$ was first worked out in Giovanazzi  _et
al._ (2006b). One can then minimise the total energy $E=E_{\rm k}+E_{\rm
trap}+E_{\rm MF}$ to find the ground-state size $\sigma_{x,y,z}$. In general,
for $\varepsilon_{\rm{dd}}>1$, there is always a global singularity with
$\kappa_{x,y}\ll 1$ and $E\rightarrow-\infty$, but a local metastable minimum
exists as well and sets the stability condition. While they coincide for very
large $N$, at finite $N$ the stability diagram from the Gaussian ansatz
differs from the one derived using a TF approximation, see, e.g., Fig. 17(a).
At low trap aspect ratios $\lambda\ll 1$, the kinetic energy acts to stabilise
the gas against collapse, as in the contact case. However, for $\lambda\gg 1$
the effect is opposite, kinetic energy increases the smallest size, which adds
a little attraction via dipolar effects and destabilises the gas.

$\bullet$ First measurements of the instability threshold:

The instability threshold was studied in 52Cr BECs of about 25,000 atoms in
tunable traps Koch  _et al._ (2008): When quenching $a$ downward using a
magnetic FR, the authors observed an abrupt and quick disappearance below a
critical value, $a_{\rm crit}$, of the BEC peak in TOF images. The
experimental values of $a_{\rm crit}$ as a function of $\lambda$ were found to
agree well with the Gaussian ansatz predictions; see Fig. 17. For the first
time, the stability and thus the very existence of a BEC was ensured solely by
dipolar interactions, and a purely dipolar BEC with $a=0$ was obtained Koch
_et al._ (2008).

##### Stability of dBEC assemblies

Consider a stack of pancake-shaped BECs, realised with a one-dimensional
optical lattice. The DDI being long-range, the total energy of a given layer
contains contributions from interactions with neighbouring layers. Since the
stability of the system is ensured by the existence of an energy minimum,
nearest-neighbour interaction can modify the stability diagram. Indeed, the
first observation of nearest-neighbour dipolar effects was reported in Fattori
_et al._ (2008) via dephasing of Bloch oscillations in a potassium optical
lattice interferometer. Fattori et al. then calculated the inter-site
interaction using a Gaussian ansatz for an individual layer. The contribution
is attractive if the dipoles are aligned with the lattice direction. This
negative contribution modifies the energy landscape and can suppress the local
minimum, destabilising the gas. This can be understood in the following way:
The neighbouring layers attract the atoms resulting in an effective repulsive
potential along this axis. The pull is stronger at the radial centre where the
neighbouring density is highest leading to a stronger effective radial
trapping. The total effect is destabilising the BEC. The instability threshold
in scattering length is thus higher. Müller et al. Müller  _et al._ (2011)
have observed this effect with a 52Cr gas in a one-dimensional optical
lattice, obtained with a retro-reflected 1064 nm beam. The difference in
critical scattering length with respect to a theory neglecting the
neighbouring layers was as high as $8\,a_{0}$, in agreement with expectations;
see Fig. 18.

Figure 18: Stability diagram of a dBEC in an optical lattice as a function of
lattice depth and scattering length. The experimental data is compared with
theory including (excluding) the influence of neighbouring traps as full (dot-
dashed) lines. Red lines correspond to a full numerical solution of the GPE
(Eq. (9)), while blue lines correspond to a Gaussian variational ansatz for
the on-site density distribution. From S. Müller, PhD thesis, Stuttgart
University.

All arguments given above consider the stability of a dBEC at equilibrium
inside a trapping potential (harmonic or optical lattice). The in-situ density
distribution determines the MF stability. However, even for a stable in-trap
density, it is possible that the release from an optical lattice into free-
space modifies the distribution in a way that makes it mostly attractive and
then induces collapse. This deconfinement-induced collapse has been observed
experimentally in Ref. Billy  _et al._ (2012). Experiments were performed in
an optical-lattice along the dipoles axis. They showed that, for large lattice
depths, no losses were observable in-trap. However the release from the trap
leads to a collapse visible in the d-wave shape of the density distribution,
see also Sec. V. Simulations of the GPE confirmed that no atomic losses were
expected in trap, but were induced by the deconfinement and following
collapse.

##### Modulational instabilities

$\bullet$ Modulational instabilities predicted from DDI:

The stability criterion developed above, based on the existence of a local
minimum in the energy landscape as a function of the BEC widths
($\sigma_{x,y,z}$), is only partial. It allows only shape-conserving
perturbations to the BEC profile. However, other perturbations containing
local density modulations are not taken into account by this criterion. The
existence of instabilities driven by higher lying modes in dBECs was first
noted in the seminal paper Santos  _et al._ (2003) predicting a roton-type
dispersion relation in anisotropic semi-infinite dBECs, see Sec. IV.1.3.
Tuning the parameters of the dBEC (i.e., its scattering length and trapping
geometry) can lead to a softening of this roton minimum ($\omega$ reaching
zero and becoming imaginary). Then, the dBEC becomes unstable and, in the
early dynamics, its population gets transferred to the roton mode leading to
the onset of density waves and, hence, a type of modulational instability.

The combined effect of axial trapping and DDI is also present in finite-size
dBECs, as long as a tighter confinement is applied along the dipoles; see also
Sec. IV.1.3. This mechanism was also shown to favour local density modulations
in the equilibrium profile of the dBEC itself, i.e., as a result of finite-
size effects Ronen  _et al._ (2007); Ronen and Bohn (2007); Dutta and Meystre
(2007). In this case, an instability can be triggered when
$\varepsilon_{\rm{dd}}>1$. It differs from a global instability in that
several local attractors gather the atomic density leading to an ensemble of
local collapses Parker  _et al._ (2009); Wilson  _et al._ (2009); Bohn  _et
al._ (2009a). This kind of instability can be seen as the softening
($\omega^{2}<0$) of a mode that is not one of the lowest-lying surface and
monopole modes: one quantum number is not the lowest possible one, for
instance the momentum $\langle k\rangle$ or the angular momentum $m$ Wilson
_et al._ (2009); Bisset  _et al._ (2013). Yet, due to confinement, this mode
still belongs to the discrete part of the spectrum where the wavenumber is not
a good quantum number : $\langle k\rangle\times R_{\rm\perp}\lesssim 1$ where
$R_{\rm\perp}$ is the typical size of the BEC in the plane perpendicular to
the magnetic field and $\langle k\rangle$ is calculated for the given mode. In
particular, in narrow regions of the $(\lambda,N)$-space, BECs with biconcave
shapes have been predicted by minimising the energy functional of the GPE
(i.e., in the MF regime) and their collapse are driven by angular roton Ronen
and Bohn (2007); Dutta and Meystre (2007); Parker  _et al._ (2009); Wilson
_et al._ (2009); Schmidt  _et al._ (2021). The softening of this intermediate-
low-lying modes also yields a type of modulational instability.

We note that, in ultracold atom experiments, modulational types of instability
were first predicted for contact interacting gases KEVREKIDIS and
FRANTZESKAKIS (2004) and investigated in this case Nguyen  _et al._ (2017);
Everitt  _et al._ (2017). These modulational instabilities are in fact of a
different nature to the one described above for dBECs. Indeed, in the contact
case, the lowest lying mode is always soft at the instability (imaginary
$\omega$). The modulational instability here arises when other modes are also
unstable and have a larger growth rate. The mode with the highest growth rate