Dipolar physics: a review of experiments with magnetic quantum gases

There is no permanent electric dipole in an atom or in a molecule in its non-degenerate rotational ground state due to their rotational symmetry. Yet when an external electric field E couples to the electric dipole moment operator, it mixes eigenstates of opposite parity. As the rotational symmetry is broken, an electric dipole moment is induced. The field to induce the electric dipole moment is lowest when the two states of opposite parity are closest in energy.

Some systems possess degenerate states of opposite parity, which allows the induced electric moment to arise at vanishingly small electric field [31]. Rydberg states in a hydrogen atom are an example of such a system: they possess electronically excited states of opposite parity that can be arbitrarily close. The electric dipole moment scales as n², with n the principal quantum number. Associating a Rydberg atom with a ground-state atom allows to form a Rydberg molecule with a permanent electric dipole moment [32]. Despite short lifetimes due to spontaneous emission, black-body radiation [33], and collisions, Rydberg gases have, over the last 10 years, been the centre of much experimental and theoretical activity. In particular, experiments are ongoing that investigate strongly correlated dipolar gases, lattice spin models, and Rydberg molecules [21, 22, 24, 34–38]. They are even the basis for a competitive quantum computing platform, which has been pushed forward by several newly founded companies (https://pasqal.io/, https://coldquanta.com/, www.quera.com/, www.atom-computing.com/).

A second, very productive field of research is the manipulation of heteronuclear molecules. In these, an electric field mixes two rotational states (for example N = 0 and N = 1) within the electronic molecular ground state. Ultracold molecular systems with a large electric dipole moment include: KRb [39–42], NaK [43–45], RbCs [46, 47], NaRb [48], KCs [49], LiCs [50], NaLi [51, 52], SrF [53, 54], H₂CO [55], CaF [56, 57], BaF [58]and YO [59, 60], HO [61, 62]. Due to their intrinsic complexity, cooling molecules has been an extremely challenging task. Recently, after many years of dedicated efforts [40, 43], the first quantum degenerate gas of polar molecules has been achieved with KRb [41, 42]. Many of the molecular systems have been shown to experience strong and rapid losses, which unfortunately presents an additional challenge for creating dense and ultracold samples. For the case of KRb, it is believed that the exo-energetic reaction KRb + KRb $\longrightarrow$ K$_2 +$ Rb₂ drives the decay [63]. For other molecular systems such as NaK and NaRb, for which the equivalent reactions are endo-energetic, the lifetime also appears rather short at large densities for reasons that are yet to be fully understood [64]. The impact of losses could be reduced thanks to an ingenious control of their spatial dependence: confining molecules in a quasi-two dimensional geometry enables one to take control off the stereodynamics of molecular reactions [65]. Producing molecules in three-dimensional optical lattices [40, 66] or optical tweezers [67, 68] prevents molecules from inelastically colliding due to their physical separation. Ultracold molecules now constitute a fast-expanding and promising field, especially for quantum simulation [19, 20].

In contrast to the situation with electric dipoles, elementary particles can have permanent magnetic dipoles even at zero field [69]. As a consequence, the effect of magnetic DDIs on quantum gases can be studied under full rotational symmetry at arbitrarily small magnetic fields. The magnetic dipole moment in atoms is primarily associated with the spin ($\hat {\boldsymbol S}$) and orbital ($\hat {\boldsymbol L}$) angular momentum of the electrons. The nucleus may also have a magnetic dipole moment, although it is three orders of magnitude smaller than the electron's. Nevertheless, the nuclear spin ($\hat {\boldsymbol I}$) couples to the electronic spin within the atom, giving rise to the hyperfine structure ($\hat {\boldsymbol F} = \hat {\boldsymbol L}+\hat {\boldsymbol S}+\hat {\boldsymbol I}$). Therefore, the sensitivity of a given Zeeman sublevel to magnetic fields indirectly depends on the nuclear spin. Only the fully stretched atomic state (i.e. maximal $F = L+S+I$ and $|m_F| = F$) reaches the full magnetic moment provided by the electrons. Here, and all along this review, X and m_X are the quantum numbers associated with the norm of the angular momentum operator $\hat {\boldsymbol X}$ and its projection along the quantization axis, respectively. Additionally, we use the dimensionless version of the vectors and operators of angular momenta and spins such that the eigenvalues associated to the norm and projection of $\hat {\boldsymbol X}$ are simply $\sqrt{X(X-1)}$ and m_X . In the above, $X = \{F,L,S,I\}$. Throughout the remainder of this review, we will usually denote by $\hat{\boldsymbol S}$ the total angular momentum of a magnetic particle.

It is possible to study dipolar physics with alkali atoms [70]. However, the energy scale associated with DDIs is rather small, typically in the Hz range. Therefore, significant focus has been on so-called highly magnetic atoms, such as chromium (Cr; with a dipole moment of 6 Bohr magnetons, µ_B ), erbium (Er; $7\,\mu_B$) and dysprosium (Dy; $10 \mu_B$). In principle, other highly magnetic atoms can be studied, such as holmium (Ho; $9\,\mu_B$), thulium (Tm; $4\,\mu_B$) or europium (Eu; $9\,\mu_B$) [71–74]. Moreover, it was demonstrated that one can use a Feshbach resonance (FR) to combine two Er atoms into a loosely bound molecule [75], which may possess up to twice the magnetic moment of the original atoms. Likewise, a nearly 20µ_B -large magnetic moment is accessible with Dy₂ molecules [76]. These systems are described in detail in section 2.

Whether an interaction, in particular the DDI, is long range depends on the exact system under study, its dimensionality, and on the exact physical question addressed. We discuss below a number of physical questions that lead to slightly different definitions of the long-range character of the interaction at hand, with particular focus on power law potentials $U(r) \propto 1/r^{\,n}$.

• Collisional point of view (in 3D). Physically, for short-range interactions, particles need to approach at small distances to interact. By decomposing the relative motion of the particles into the relative orbital angular momentum eigenstates (so-called partial waves), denoted by the quantum numbers (l, m) for the momentum's norm and projection eigenvalues, one finds that at low collision energy, the centrifugal barrier $\frac{2l(l+1) \hbar^2}{m r^2}$ prevents particles from approaching in higher partial waves l > 0. That is, the contributions from high partial waves vanish. For a $1/r^{\,n}$ interacting potential, the scattering phase shift $\delta_l(k)$ at low collision momentum k scales as $k^{2l+1}$ if $l\lt (n - 3)/2$ and as $k^{\,n-2}$ otherwise [78, 79]. Therefore, for $n\geqslant 4$, the interaction is purely s-wave at low energy and short ranged. In contrast, for n = 3, $\delta_l(k) \propto k$ for all partial waves. Therefore, all partial waves contribute to the scattering process even at low collision energy. The interaction is then long range and can be felt beyond the centrifugal barrier. The long-range character of the DDI is spectacularly manifest in the fact that polarised fermionic dipolar gases thermalize despite the absence of s-wave interactions (due to the Pauli exclusion principle). This is in contrast to nondipolar polarised Fermi gases; see section 3. Note: A thorough treatment of the above should account for the fact that the DDI is not a pure central potential $U(r) \propto 1/r^3$ due to its anisotropic character. This is of particular consequence for inelastic dipolar collisions, which necessarily involve the anisotropic character of the interaction—see section 1.3.2—and are actually short-range processes at large magnetic field despite the same $1/r^3$ scaling as their elastic counterparts. See section 3.3.3.We also remark that the scattering picture can be modified in the presence of strong confinement, in particular in reduced dimensions; see, e.g. section 7.1.
• Thermodynamic point of view. Short-range interactions lead to an energy that is thermodynamically extensive. This is true when $\int^{\infty}_0 U(r) d^D r$ converges, which only happens when n > D, where D is the spatial dimension. Thus, from this point of view, $1/r^3$ interactions are long range in 3D, but short range in 2D and 1D. Long-range-interacting systems possess peculiar thermodynamic properties, such as the non-equivalence of thermodynamic ensembles, the possibility for negative specific heat, and the spontaneous formation of structures. These arise because the hypothesis that the energy is additive in the thermodynamic limit—i.e. given the energy of two subsystems A and B, $E(A \cup B) = E(A)+E(B)$—breaks down when the interaction between the subsystems cannot be neglected [80]. The status of DDIs in 3D is therefore marginal since while there are distant couplings between sub-systems A and B, the integral $\int^{\infty} U(\boldsymbol{r}) d^3 \boldsymbol{r}$ does converge due to the peculiar $d-$wave shape of DDIs.
• Many-body physics perspective. In the context of many-body physics, DDIs may lead to qualitatively new behaviour, even when D < 3. For example, the Mermin–Wagner theorem precludes the possibility of spontaneous breaking of a continuous symmetry and of long-range order in (homogeneous) low-D systems. However, in 2D, this applies only for short-range interacting systems with n > 4 [17, 81]. Therefore, in this context, DDIs can be seen as long-ranged even in 2D. Indeed, it has been predicted that ferromagnetic ordering should be stable in 2D for DDIs [82]. The meaning of long range in 1D for the DDI will be addressed in section 7.1.
• Mathematical physics perspective. Let us for completeness also briefly mention the mathematical physics point of view of the meaning of 'long-range interaction.' For an interaction potential ${\sim}1/r^{\,n}$, the scattering wave is only described by an asymptotic outgoing spherical wave weighted by an angle-dependent scattering amplitude for n > 1. For $n\leqslant1$, e.g. the Coulomb potential, there are logarithmic corrections to the general form of the outgoing spherical wave, which defines another border between long and short-range potentials.

The anisotropic character of the DDI greatly impacts the properties of dipolar gases. It introduces profound differences from the point of view of two-body physics and scattering properties—see section 3—and also on the collective many-body properties of quantum degenerate dipolar gases. In particular, the stability diagram of dipolar condensates is affected by an interplay between the anisotropy of the trap and the anisotropy of the interactions; this will be described in sections 4 and 5. We now briefly describe a few basic consequences of this anisotropy.

The DDI is attractive in one direction and repulsive in the other two directions. The shape of the interaction follows a d-wave form mathematically described by the j components of the spherical harmonics ${Y_2^{j}}$, in particular with j = 0 in a fully polarized situation. This means that when integrated over all space in 3D, the DDI between polarised dipoles converges to zero for a 3D homogeneous gas. Consequently, the mean-field (MF) physics of dipolar gases is dominated by border and boundary effects: the average interaction between particles will strongly depend on the shape of the cloud. In particular, an elongated trap along the axis of the dipoles will favour the collapse of the gas due to the predominately attractive interaction. The stability of dipolar condensates as a function of geometry is described in section 4.

Another consequence of anisotropy is the existence of a special angle between the dipoles and the interatomic axis, $\theta_m = \arccos \sqrt{1/3} \approx 54.74^\circ$, at which DDIs vanish. More generally, controlling this angle can be used to tune the strength of DDIs, especially when performing experiments in reduced dimensions, as described in section 7.1. In a scheme inspired from NMR techniques, it has been suggested [83] and demonstrated [84] that by using time-varying magnetic fields, it is possible to time-average the DDI to reduce its amplitude or reverse its sign.

Finally, the anisotropy of the interaction has fundamental consequences from the point of view of collisions. The interaction potential is not central, and therefore the orbital angular momentum does not need to be conserved during a collision. The expansion in spherical harmonics yields the selection rule for angular momentum transitions $\Delta l = (0, \pm 2)$. Moreover, partial waves of differing l that contribute to the scattering become coupled. Finally, the angular momentum of the atoms' internal state may also change during the collisions, opening the possibility for inelastic processes.

In view of describing the physical processes at play when two dipolar particles collide, it is useful to rewrite the dipolar potential between atoms 1 and 2 in terms of quantum operators [85]:

$$\begin{align} V_{\mathrm{dd}}(\boldsymbol{r}) = \frac{d^{\,2}}{r^3} & \bigg[\bigg(S_1^z.S_2^z + \frac{1}{2} \left( S_1^+.S_2^- + S_1^-.S_2^+\right) \nonumber \\ & \;\; -\frac{3}{4} \left( 2 z S_1^z + r^- S_1^+ + r^+ S_1^- \right) \nonumber \\ & \;\; \times\left( 2 z S_2^z + r^- S_2^+ + r^+ S_2^- \right) \bigg], \end{align} \tag{ 2 }$$

where $(x,y,z)$ is the normalised unit vector connecting both atoms, $r^+ = (x + iy)$, $r^- = (x - iy)$, $S^+ = (S^x + iS^y)$, and $S^- = (S^x - iS^y)$. Note that both here and throughout this review, we use the dimensionless version of the vectors and operators of angular momenta and spins.

We describe three physical processes that arise from this expression:

• (a)  
  Elastic dipole–dipole interactions, where the spin of each atom is conserved in time:

  $$\begin{align} V_{\mathrm{dd}}^\mathrm{\,el}(\boldsymbol{r}) = \frac{d^{\,2}}{r^3} S_1^z.S_2^z \left(1-3 z^2 \right). \end{align} \tag{ 3 }$$

  The experimental manifestation of this anisotropic Ising term on quantum degenerate dipolar gases has been extensively studied. It is the main process at play for most of the results presented in this review article: see section 3 on scattering physics, sections 4 and 5 on the collective properties of dipolar gases, the stability diagram, the instability dynamics and the stabilisation of so-called dipolar droplets, section 7.1 on integrability breaking in 1D gases, and section 7.2 on the extended Bose–Hubbard model.
• (b)  
  Exchange interactions, where two atoms exchange one unit of spin (Zeeman) excitation, while the total magnetisation and energy is conserved (in the absence of quadratic Zeeman effects):

  $$\begin{align} V_{\mathrm{dd}}^\mathrm{ex}(\boldsymbol{r}) = -\frac{1}{4} \frac{d^{\,2}}{r^3} \left(S_1^+.S_2^- + S_1^-.S_2^+\right) \left(1-3 z^2 \right). \end{align} \tag{ 4 }$$

  This exchange term can drive spin dynamics at constant magnetisation as described in section 6.2.3 and dictates the physics of spinor dipolar gases in deep lattices, which is the topic of section 7.3. We note that the elastic and the exchange terms result in an anisotropic Heisenberg-like term (the so-called XXZ model).
• (c)  
  Relaxation terms describe the modification of the longitudinal magnetisation of the pair of atoms during the collision. There are two possible processes:

  $$\begin{align} V_{\mathrm{dd}}^{\,\mathrm{rel}_1}(\boldsymbol{r}) & = -\frac{3}{4} \frac{d^{\,2}}{r^3} (r^+)^2 S_1^-.S_2^-, \nonumber \\ V_{\mathrm{dd}}^{\,\mathrm{rel}_2}(\boldsymbol{r}) & = -\frac{3}{2} \frac{d^{\,2}}{r^3} z r^+ (S_1^z.S_2^-+S_2^z.S_1^-), \end{align} \tag{ 5 }$$

  plus the conjugate processes. Spin momentum and angular orbital momentum exchange while the magnetic energy is transferred into kinetic energy. The second process describes single spin flips, while the first describes double spin flips (i.e. both atoms flipping their spin). These terms underlie most of the results presented in sections 3.3 and 6, the latter describing spinor physics with free magnetisation.

The cross section classically describes the area, transverse to the relative motion, within which two particles must meet to scatter. In other words, the scattering cross section is related to the typical distance at which the wavefunction of the relative motion is distorted by the interaction. Employing the Heisenberg uncertainty principle, this distance $r_{\mathrm{d}}$ for the DDI is typically set by an interplay between the DDI strength $S^2 d^{\,2}/r_{\mathrm{d}}^3$ and the energy cost to bend the wave function by an amount $r_{\mathrm{d}}$. Setting $S^2 d^{\,2}/r_{\mathrm{d}}^3 = \hbar ^2 /m r_{\mathrm{d}}^2$ defines the dipolar length [86]:

$$\begin{align} {a_\mathrm{dd}} \equiv \frac{r_{\mathrm{d}}}{3}= \frac{S^2 d^{\,2} m}{3 \hbar^2}= \frac{C_\mathrm{dd} m}{3 \hbar^2}, \end{align} \tag{ 6 }$$

where m is the atomic mass. The order of magnitude of the scattering cross section is:

$$\begin{align} \sigma \approx r_{\mathrm{d}}^2 = \frac{S^4 d^4 m^2}{\hbar^4}. \end{align} \tag{ 7 }$$

Likewise, one defines the range of the van der Waals potential $V_\mathrm{VdW} = -C_6/r^6$ as $r_\mathrm{vdW} = (m C_6/\hbar^2)^{1/4}$, which sets the typical scattering cross section due to short-range interactions. The lengths r_VdW and ${a_\mathrm{dd}}$ are typically in the nm range; i.e. much larger than both the Bohr radius a₀ and the typical impact parameter at room temperature.

In section 3, we will describe the scattering theory for both dipolar and van der Waals interactions. The DDI cross sections are presented based on a first-order Born approximation, and the role of exchange statistics in these expressions is discussed.

In contrast to the van der Waals case, the dipolar cross section depends on only the mass of the atoms and their dipole moment. Because it is independent of the details of the molecular potentials, dipolar scattering assumes a universal character. One remarkable aspect is that the dipolar cross section follows the same universal scaling of equation (7) (up to numerical factors) regardless of particle exchange statistics. In particular, identical fermions have a finite dipolar cross section even at vanishingly small collision energy. This is a direct consequence of the long-range character of DDI, as discussed above. This topic will be discussed in section 3.2. Inelastic dipolar scattering will be discussed in section 3.3.

Finally, by integrating over all directions of the collision, the scattering theory outlined above obscures one of the central features of dipolar scattering, the anisotropic dependence on the colliding angle, which has been observed in both Er and Dy. This is the topic of section 3.4.

The form of the Fourier transform of the interaction potential often provides insight regarding the physics of interacting particles. For example, it facilitates the description of two-body scattering physics because scattering theory tends to formulate the wavefunction in terms of momentum states. It also proves convenient in discussing elementary excitations of a quantum gas, which, in a uniform system, are characterised by a well-defined momentum.

The Fourier transform of the elastic part of the DDI, equation (3), is:

$$\begin{align} \widetilde{V_{\mathrm{dd}}}({\boldsymbol{k}}) = \int e^{i {\boldsymbol{k}}\boldsymbol{r}} V_{\mathrm{dd}}^\mathrm{\,el}(\boldsymbol{r})d^3r = \frac{C_\mathrm{dd}}{3} \left[3 \cos(\theta_k)^2 -1 \right], \end{align} \tag{ 8 }$$

where θ_k is the angle between k and the polarisation of the dipoles. This form of the interaction is remarkable when compared to the contact interaction. While neither the Fourier transforms of the contact interaction nor the DDI depend on k, the Fourier transform of the DDI retains a nontrivial angular dependence. Such a feature can give rise to an anisotropic dispersion relation of excitations, as described in section 4.

That $\widetilde{V_{\mathrm{dd}}}({\boldsymbol{k}})$ does not depend on the modulus of k can be understood from dimensional analysis: $\int d^D r \exp(i k r) \frac{1}{r^3}$ is independent of k for D = 3. On the other hand, in a 2D system (D = 2), we expect a different behaviour with a linear dependence on k for small k. Quite generally, the fact that the Fourier transform of $V_{\mathrm{dd}}$ has a k dependence is an important feature for D ≠ 3 systems, in that the DDI introduces a tendency in these systems to develop structured excitations (rotons, solitons) and exotic phases (supersolid, crystals of quantum droplets, etc). These excitations are described in section 4. The quantum droplets occur when the gas spontaneously forms stable spatial arrangements of liquid-like droplets in dipolar Bose–Einstein condensates (dBECs) driven to the instability point of mechanical collapse. This is related to the emergence of new phases stabilised by beyond-mean-field (BMF) effects and is described in section 5; see also section 1.4.1.

1.4.1.1. Spin-polarised dipolar Bose quantum gases in the mean-field regime.

Many-body physics is often intractable. However, most of the first experiments on ultracold gases of magnetic atoms have been performed with weakly interacting BECs, and the associated theory is tractable because interatomic correlations are small and so MF theories apply.

Due to Bose stimulation, in which the population of bosonic atoms at low energy favours the occupation of a unique single-particle orbital, it is natural to propose a variational ansatz where the many-body wavefunction is assumed to be $\phi(r_1,\ldots,r_N;t) = \prod_{i = 1,\ldots,N} \psi(r_i,t)$.

The single-particle wave-function $\psi(r,t)$ is taken as a variational parameter to minimise the system's total energy. This approach leads to the well-known Gross–Pitaevskii equation (GPE), which has been found to describe most of the properties of dilute BECs [9, 87–89].

When all atoms are polarised (and the polarisation axis set to z), the GPE of a dBEC is:

$$\begin{align} i \hbar\frac{\partial}{\partial t} \psi = \left[\frac{- \hbar^2 \nabla ^2}{2 m}+V_\mathrm{tr}(r) +g \left| \psi \right|^2+\Phi_\mathrm{dd}(r,t)\right] \psi, \end{align} \tag{ 9 }$$

where $V_\mathrm{tr}(r)$ is the trap potential, $g = \frac{4 \pi \hbar^2}{m} a$ is the coupling constant describing contact interactions of s-wave scattering length a; see also sections 2.4 and 3 [9, 90]. The MF associated with the DDI is $\Phi_\mathrm{dd}(r,t)$ [91, 92]:

$$\begin{align} \Phi_\mathrm{dd}(r,t) & = \int dr^{^{\prime}} \left| \psi(\boldsymbol r^{^{\prime}},t) \right|^2 U_\mathrm{dd}(\boldsymbol r- \boldsymbol r^{^{\prime}}), \end{align} \tag{ 10 }$$

$$\begin{align} U_\mathrm{dd}(\boldsymbol r) & = C_\mathrm{dd} \frac{1-3 \cos^2(\theta)}{\left| r \right|^3}, \end{align} \tag{ 11 }$$

where θ is the angle between $\boldsymbol r$ and the polarisation axis z. Only the elastic part of the DDI, equation (3), contributes due to polarisation in a fully stretched Zeeman substate. This term is non-linear and nonlocal. To quantify the strength of the DDI with respect to the contact interactions within a BEC, it is useful to introduce the dimensionless parameter:

$$\begin{align} \varepsilon_{\mathrm{dd}}=\frac{C_\mathrm{dd} m}{3 \hbar^2 a}=\frac{{a_\mathrm{dd}}}{a}. \end{align} \tag{ 12 }$$

We note that writing the GPE of a dipolar BEC in the form of equation (9) is in fact non-trivial. Its validity, which relies on describing the total inter-particle interactions via an effective pseudo-potential that is the simple sum of the contact pseudo-potential and the DDI potential, has been long debated. The efforts to prove the validity of this treatment as well as identify its limitation will be reviewed in section 4.1.1. The applicability of the nonlocal GPE equation (9) for the case of weakly interacting trapped BECs of magnetic atoms in the stable regime, e.g. $\varepsilon_{\mathrm{dd}}\lt 1$ in a 3D isotropic trap, has been supported by numerous theory and experimental works. In this regime, the anisotropic and nonlocal character of equation (10) substantially modifies the static and dynamical properties of the BEC compared to contact-only BECs. This will be extensively discussed in section 4.

1.4.1.2. Spinor dipolar Bose quantum gases.

In the presence of spin degrees of freedom, the exact form of the GPE depends on the spin of the atoms and can be found in [93]. Taking, for example, the case of a spin-1 atom—i.e. the simplest example pertaining to bosonic physics—the GPE takes the form:

$$\begin{align} i \hbar \frac{d \psi_m}{dt}&= \left[\frac{- \hbar^2 \nabla ^2}{2 m}+V_\mathrm{tr}(r) -pm+q m^2 \right] \psi_m \nonumber \\ &\quad+c_0 n \psi_m + c_1 \sum_{m^{^{\prime}}=-1}^1 {\boldsymbol S}.{\boldsymbol s}_{m,m^{^{\prime}}} \psi_{m^{^{\prime}}}\nonumber\\ &\quad +C_\mathrm{dd} \sum_{m^{^{\prime}}=-1}^1 \textbf{b}. {\boldsymbol s}_{m,m^{^{\prime}}} \psi_{m^{^{\prime}}}, \end{align} \tag{ 13 }$$

where ψ_m denotes the macroscopic wave function associated with the spin state of projection quantum number m. The terms in p and q describe the linear and quadratic Zeeman energy shifts of the spin states, respectively. The trap is assumed to be spin-independent. The terms proportional to c₀ and c₁ are spin-independent and spin-dependent contact interactions, respectively. The spin density vector is ${\boldsymbol S}$, and ${\boldsymbol s} = \{s^x,s^y,s^z\}$ are the spin matrices. The DDIs are described by the term proportional to $C_\mathrm{dd}$, where the effective dipole field b is defined by:

$$\begin{align} b_{\nu} = \int dr^{^{\prime}} \sum_{\nu^{^{\prime}}} Q_{\nu,\nu^{^{\prime}}} (\boldsymbol{r}-\boldsymbol{r}^{^{\prime}}) S^{\nu^{^{\prime}}} (r^{^{\prime}}), \end{align} \tag{ 14 }$$

With:

$$\begin{align} Q_{\nu,\nu^{^{\prime}}}(\boldsymbol{r}) = \frac{\delta_{\nu,\nu^{^{\prime}}}-3 r_{\nu} r_{\nu^{^{\prime}}}}{r^3}, \end{align} \tag{ 15 }$$

and $\nu,\nu^{^{\prime}} = \{x,y,z\}$.

Equation (13) is central to the description of spinor dipolar physics, which is the subject of section 6. In general, the DDIs cannot be neglected when $C_\mathrm{dd}$ is comparable to either c₀ or c₁. If $c_1\ll c_0$, then the DDIs can be significant even when $\varepsilon_{\mathrm{dd}}\ll1$ [94]. Furthermore, magnetisation-changing processes effected by the term in equation (14) have no analogue in systems with only spherically-symmetric contact interactions. Such processes start to play a role when $C_\mathrm{dd} n \simeq \{p,q\}$.

Finally, it is useful to stress that in the MF regime, DDIs described by equations (10) and (14) correspond to the average magnetic field produced by all atoms within the condensate. This is due to the fact that correlations between atoms have been neglected, which is the essence of the MF approximation. In the case of spinor gases, the effect of quantum fluctuations and correlations can, however, be significant, even in the weakly interacting regime. This is due to entanglement or squeezing naturally arising in the spin degrees of freedom. The consequences of DDIs on the properties and behaviours of gases with spin degrees of freedom will be detailed in section 6.

1.4.1.3. Elementary excitations of a spin-polarised dipolar Bose quantum gas.

The elementary excitations of a BEC are usually well described within a Bogoliubov treatment, which matches a linearisation of the GPE around the ground-state wavefunction [9, 87]. The theory yields a simple dispersion relation for a uniform 3D gas ($V_\mathrm{tr} = 0$):

$$\begin{align} \epsilon({\boldsymbol{k}}) =\sqrt{\frac{\hbar^2k^2}{2m}\left(\frac{\hbar^2k^2}{2m}+2n\widetilde{V_\mathrm{int}}({\boldsymbol{k}})\right)}. \end{align} \tag{ 16 }$$

This describes the energy of the elementary excitation of momentum k in a BEC of density n. Here, $\widetilde{V_\mathrm{int}}({\boldsymbol{k}})$ is the Fourier transform of the total interaction potential. In contact-interacting gases, $\widetilde{V_\mathrm{int}}({\boldsymbol{k}}) = g$. For a dBEC, $\widetilde{V_\mathrm{int}}({\boldsymbol{k}}) = g+\widetilde{V_{\mathrm{dd}}}({\boldsymbol{k}})$. From the form of $\widetilde{V_{\mathrm{dd}}}({\boldsymbol{k}})$ (see equation (8)), one can infer that the energy of a collective excitation of a dipolar fluid depends not only on the magnitude of its wavevector but also on its propagation direction. The dispersion relation of elementary excitations of an isotropic, homogeneous dipolar fluid is anisotropic. For a dBEC, the dispersion retains the initial linear phonon character, but with an anisotropic speed of sound c; the dispersion relation remains monotonic.

This picture is modified in a constrained geometry, where one externally lifts the spatial symmetry along at least one dimension by, e.g. imposing anisotropic trapping confinement [95]. The trap along the constrained dimension yields a new length scale. Because of the anisotropy and long-range character of the DDI, this length scale also becomes relevant to the description of the physics of the otherwise translationally invariant directions. In particular, it affects their elementary excitations. That is, a BEC that is more tightly trapped along the dipoles' direction than transversely possesses a favoured wavelength in its dispersion relation at which the energy of the transverse excitations reaches a minimum [95–99]. This is referred to as a roton minimum in analogy to a similar minimum found in the dispersion relation of liquid helium [100–102]. These properties of dBECs are the topic of section 4.1.3.

1.4.1.4. Mean-field instability and collapse.

1.4.1.5. Quantum states stabilized by fluctuations: droplets and supersolids.

Fermionic dipolar atoms are also of great interest for exploring new physics. A remarkable property of dipolar Fermi gases lies in the fact that polarised samples remain interacting even in the ultracold regime. This is unlike nondipolar Fermi gases; see section 1.3.1. Yet, the MF theory developed above does not appropriately describe fermionic ensembles because the ansatz used to write the many-body wavefunction is incompatible with the Pauli exclusion principle: it must be antisymmetrized due to fermionic exchange statistics. Therefore, it is generally not possible to neglect correlations in a Fermi gas at low temperature, even for small interactions. This makes a theoretical treatment of fermionic gases challenging. The simplest treatment of MF theory that includes the antisymmetrization of the wavefunction replaces the product ansatz used in section 1.4.1 for $\psi(\boldsymbol{r}_1,\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N)$ by a Slater determinant. This procedure is known to be sufficient for a pure state without interactions. With interactions, it may still be sufficient, but with the single-particle wave-functions modified compared to the noninteracting case. This approach constitutes the Hartree–Fock theory [117, 118].

The mean DDI energy for an ensemble of N atoms in a state $\psi(\boldsymbol{r}_1,\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N)$ is generally written as:

$$\begin{align} E_\mathrm{dd} &= \int dr_1 \ldots dr_N \psi^*(\boldsymbol{r}_1,\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N)\nonumber\\ &\quad\times\sum_{i,j} V_{\mathrm{dd}}(\boldsymbol{r}_i-\boldsymbol{r}_j) \psi(\boldsymbol{r}_1,\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N)\nonumber\\ &= \frac{N(N-1)}{2} \int dr dr^{^{\prime}} \int dr_3 \ldots dr_N \nonumber\\ &\quad\times \psi^*(\boldsymbol{r},\boldsymbol{r}^{^{\prime}},\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N) V_{\mathrm{dd}}(\boldsymbol{r}-\boldsymbol{r}^{^{\prime}}) \psi(\boldsymbol{r},\boldsymbol{r}^{^{\prime}},\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N). \end{align} \tag{ 17 }$$

Because of antisymmetrization, the integral $N(N-1)\times \int dr_3 \ldots dr_N \psi^*(\boldsymbol{r},\boldsymbol{r}^{^{\prime}},\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N) \times \psi(\boldsymbol{r},\boldsymbol{r}^{^{\prime}},r_{3}, \ldots ,r_N)$ does not reduce to $n(\boldsymbol{r}) n(\boldsymbol{r}^{^{\prime}})$ as in the bosonic case. Though by using the Slater determinant ansatz, it can be simplified to $\rho_{1}(\boldsymbol{r},\boldsymbol{r}) \rho_{1}(\boldsymbol{r}^{^{\prime}},\boldsymbol{r}^{^{\prime}}) - \rho_{1}(\boldsymbol{r},\boldsymbol{r}^{^{\prime}}) \rho_{1}(\boldsymbol{r}^{^{\prime}},\boldsymbol{r})$, where

$$\begin{align} \rho_{1}(\boldsymbol{r},\boldsymbol{r}^{^{\prime}})&= N(N-1)\times \int dr_2dr_3 \ldots dr_N\; \psi^*(\boldsymbol{r},\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N)\nonumber\\&\quad \times \psi(\boldsymbol{r}^{^{\prime}},\boldsymbol{r}_2,\boldsymbol{r}_{3}, \ldots ,\boldsymbol{r}_N), \end{align} \tag{ 18 }$$

is the one-body density matrix. Therefore, the DDI MF energy for the fermionic gas consists of two parts: the usual term, also called the direct or Hartree term:

$$\begin{align} E_{dd}^\mathrm{\,dir} = \frac{1}{2} \int dr dr^{^{\prime}} V_{\mathrm{dd}} (\boldsymbol{r}-\boldsymbol{r}^{^{\prime}}) n(\boldsymbol{r}) n(\boldsymbol{r}^{^{\prime}}), \end{align} \tag{ 19 }$$

and an unusual term, called the Fock or exchange term, resulting from the requirement for an antisymmetric wavefunction upon particle exchange:

$$\begin{align} E_\mathrm{dd}^\mathrm{\,exc} =- \frac{1}{2} \int dr dr^{^{\prime}} V_{\mathrm{dd}} (\boldsymbol{r}-\boldsymbol{r}^{^{\prime}}) \rho_{1} (\boldsymbol{r}^{^{\prime}},\boldsymbol{r}) \rho_{1} (\boldsymbol{r},\boldsymbol{r}^{^{\prime}}). \end{align} \tag{ 20 }$$

This exchange term is zero in the case of a BEC. Based on equations (19) and (20), and by performing the variational minimisation of the total energy with respect to ψ, one can derive semiclassical (Hartree–Fock) equations for the degenerate Fermi gas (DFG). To describe a trapped Fermi gas, one can use a local-density approximation (LDA), which assumes that the atoms feel a local DDI [118–123]. This particular exchange interaction term, arising from the interplay of fermionic statistics and the nonlocal DDI, has several physical consequences—these will be the topic of section 4.2.

As previously discussed in section 1.4.1, BMF effects are typically weak for BECs in the weakly interacting regime (far from any instability). This is because the interaction energy is too small to create short-wavelength correlations in the gas. One way to reach strong correlations is to load the atoms into tight anisotropic traps or standing waves of light (so-called optical lattice). Confined trapping geometries effectively reduce the atoms' kinetic energy by restricting motion. By doing so, they allow the interaction and kinetic energies to play competing roles in the determination of how the system organises [7, 8, 81, 124–128]. In this review, we will discuss the experimental progress based on magnetic atoms in confined geometries; see section 7. We note that important advances have been made with systems of polar molecules [19, 20, 40, 129] as well as of Rydberg atoms [21, 23, 24, 130–132]. Focusing on magnetic atoms, we discuss three main areas: (a) the physics in low-dimensional spaces and in particular 1D, where the motion of the particles arise only in some directions of space and is frozen transversely; (b) the physics of spinless particles whose motion occurs, this time, along the specially confined directions in space, and in particular, along directions of a periodic external potential formed by an optical lattice. This realises extended Hubbard models for spinless dipolar particles; (c) the case of spinful dipolar particles in such periodic external potentials, leading to quantum magnetism and XYZ models.

1.4.3.1. Dipolar gases in lower dimensions.

1.4.3.2. Extended Hubbard model.

1.4.3.3. Spin physics in optical lattices.

The engineering of synthetic coupling involving the particle's spin, such as spin-spin or spin–orbit coupling (SOC), opens the door to the realisation of exotic states such as, on the many-body level, topological superfluids [147–154] as well as highly nonclassical or topological spin states [155–158], that can even be produced at the single-atom level for large-spin atoms. Spin-dependent light shifts or optical Raman dressing of internal spin states may be used to effect such coupling for the introduction of abelian (magnetic field-like) or non-abelian (SOC-like) gauge fields [159, 160]. Unfortunately, however, light fields also heat the atoms due to spontaneous emission, limiting the lifetime of systems. This may be circumvented by exploiting the level structure of magnetic atoms such as Dy and Er, while allowing dipolar interactions to play a role in the physics. Dipolar relaxation then sets new limits on the lifetime of these gases. Fortunately, there exist platforms in which the rate of dipolar relaxation is significantly reduced. One method uses a tightly confining potential in one or more spatial dimensions to suppress relaxation via phase-space restriction; see [161–163] and section 3.3.5. Another method employs fermions in a large magnetic field so that dipolar relaxation may be suppressed due to Fermi statistics; see [161, 164] and section 3.3.4. Exploiting such possibilities has yielded the realisation of long-lived SOC Fermi gases [165]. Related achievements and prospects are discussed in section 6.4.1.

The electronic configuration of magnetic Lns is $[Xe]4f^{m}6s^{2}$. It is characterised by a xenon-like core, an inner open 4 f shell with m valence electrons, and an outer closed 6 s shell. Due to the electron vacancies in the inner shell, magnetic Ln are often called submerged-shell atoms. The unfilled 4 f shell plays a particularly important role in their high magnetism and orbital anisotropy.

Figure 1 shows the atomic energy spectra for the cases of Er and Dy, up to a wavenumber of $25\,000 {\text{ cm}^{-1}}$. Both species have an even-parity ground state with a large J and many excited states of odd and even parity of various quantum numbers. A comprehensive set of spectroscopic data of all Ln elements can be found in references [216–218]. Note that theory results that are based on the Cowan suite of codes [219, 220] and used to estimate the polarizability of Lns predict atomic transitions that have not yet been observed [212, 213]; see also section 2.3. This reveals the still incomplete knowledge of Ln atomic spectra and their properties.

Zoom In Zoom Out Reset image size

Unusually for those accustomed to alkali atoms, many of the energy levels in the electronic spectrum of Ln atoms do not conform to the usual LS-coupling scheme. In the LS scheme, the total electron spin S and angular orbital momentum L couple to form $\boldsymbol{J} = \boldsymbol{L}+\boldsymbol{S}$, but the large SOC of the electrons in Lns renders this scheme sub-optimal for some of the electronic levels. In this case, a $J_{1}J_{2}$-coupling scheme is more appropriate [221]. In this $J_{1}J_{2}$ scheme, the electrons in each shell couple independently in a LS-coupling scheme, and then the total angular momentum quantum numbers from the different shells sum together. For instance, 4 f electrons give rise to J₁ and the 6 s ones to J₂, and $\boldsymbol{J} = \boldsymbol{J_{1}}+\boldsymbol{J_{2}}$ with quantum number J, which is denoted as $(J_{1},J_{2})_{J}$. For Lns, the LS scheme remains relevant for the ground state while the excited states (where SOC is stronger) typically follow a $J_{1}J_{2}$-coupling scheme. Finally, the bosonic isotopes of Dy and Er have nuclei with an even number of protons and neutrons. This results in a zero nuclear spin I and no hyperfine structure. In contrast, the fermionic isotopes have an even number of protons and an odd number of neutrons, resulting in a nuclear spin $\mathbf{I} = 5/2$ and $\mathbf{I} = 7/2$ for Dy and Er, respectively. We note that other magnetic Lns, such as Tm, Eu or Ho, have only bosonic isotopes and those isotopes have a hyperfine structure.

Lanthanides' atomic spectra offer a rich collection of $J\longrightarrow J+1$ optical lines, including broad-, narrow- and ultra-narrow-linewidth transitions. Many of such transitions can be used for optical manipulation and laser cooling, and they are readily accessible with common lasers. As a rule of thumb, the linewidth of the transitions is larger for short-wavelength (i. e.  high energy) transitions. The strongest line in Er (Dy) at 401 (421) nm, the intercombination line at 583 (626) nm, and the narrow line at 841 (741) nm are used in experiments for laser cooling, as discussed in the next section. The narrower resonances have also been used for spin manipulation or coupling, see section 6. Furthermore, the spectra of the $J\longrightarrow J$ and $J\longrightarrow J-1$ transitions are equally rich and relevant for optical manipulation schemes [222]. Among the rich spectra of Lns, an interest is also growing for the ultranarrow resonances of Hz-level linewidth, see e.g. [223–226]. Finally, the orbital anisotropy of Lns also yields a dependence of the matter-light interaction on the atomic spin. This presents both challenges (how to obtain equal trapping of all spin states?) and advantages (realisation of spin-dependent potentials, SOC, etc) for ultracold gas experiments; see also latter discussions.

Trapping and cooling of open-shell Ln atoms have been achieved using ZSs, MOTs, and ODTs, but without the use of magnetic trapping. The overall efficiency is similar to that of alkali-metals [192], though there are some key differences that we now describe.

2.2.2.1. Optical cooling and MOTs.

2.2.2.2. Optical dipole trapping and cooling to quantum degeneracy.

In presence of a (monochromatic) laser field, a semiclassical approximation may be used for $\hat{\boldsymbol E}$: $\hat{\boldsymbol E} = {\boldsymbol E}(\hat{\boldsymbol r},t) = {\mathcal E}_0(\hat{\boldsymbol r})\left({\hat \varepsilon}\exp(-i\omega_0 t+\phi_0(\hat{\boldsymbol r}))+c.c.\right)$ with ${\mathcal E}_0$, ω₀, φ₀ and $\hat \varepsilon$ the amplitude, frequency, phase and polarisation unit vector of the laser field. Within the semi-classical approximation, the atomic dipole operator is non-zero thanks to the charges' displacement induced by the electric field. This results in the simplified form:

$$\begin{align} \hat{\boldsymbol D}(\boldsymbol r,t)=\alpha\left(\omega_0 \right){\boldsymbol E}(\hat{\boldsymbol r},t), \end{align} \tag{ 22 }$$

here α is the complex dynamical atomic polarisability at frequency ω₀. It is a sum of all electronic transition contributions, which typically scale linearly with the transition linewidth (setting the electronic transition strength) and vary inversely to the detuning of ω₀ to the transition frequency [245, 252]. The real part of α relates to the (reactive) dipole force while its imaginary part relates to the (dissipative) radiation pressure force and to the photon-scattering rate. Practically, the former determines the AC Stark shifts induced by the light on the atomic levels' energies and, thus, to the induced dipole trap depth. The latter limits the trap lifetime.

Because of the internal structure of the atoms, the atom-light coupling depends on the total angular momentum of the electronic state $\boldsymbol F$ and on the light polarisation . This yields a (rank-2) tensor structure for the polarisability α. It is again a sum over all electronic transitions, but now accounting for the various electric dipole transition elements between the ground and excited states' sublevels as well as their possible spatial anisotropy. The polarisability tensor α can be decomposed into a scalar α_s , vector α_v , and tensor α_t contribution such that [220, 253, 254]:

$$\begin{align} \alpha(\omega_0) & = \alpha_s(\omega_0)\boldsymbol{\bar{\bar{1}}}+ \alpha_v(\omega_0)\frac{({\boldsymbol \epsilon}^* \times {\boldsymbol \epsilon}).{\boldsymbol F}}{2F}\nonumber\\ & \quad + \alpha_t(\omega_0)\left[\frac{3({\boldsymbol \epsilon}.{\boldsymbol F})({\boldsymbol \epsilon^*}.{\boldsymbol F})+3({\boldsymbol \epsilon^*}.{\boldsymbol F})({\boldsymbol \epsilon}.{\boldsymbol F})-2{\boldsymbol F}^2}{2F(2F-1)}\right]. \end{align} \tag{ 23 }$$

The scalar part is independent of the angular momentum operator, while the higher-rank contributions yield spin-dependent terms. In this respect, vector and tensor polarisabilities set the strength of Raman coupling between atomic Zeeman sublevels. They also provide tools for optical spin manipulation and spin-dependent trapping. Their contribution's dependence on the light polarisation additionally induces an anisotropy of the atom-light coupling, varying with the relative angle between the quantisation axis (set by the magnetic field) and the light-beam polarisation vector.

A major difference in the atom-light interactions involving magnetic atoms versus, e.g. alkali atoms, lies in their large vectorial and tensorial polarisabilities. This arises from the large electronic spin of the magnetic atoms, which yields strong SOC in electronic states with L ≠ 0, both in ground and excited states. This strong SOC results in differences between the electric dipole matrix elements coupling between the different fine-structure levels of the ground and the electronically excited states, which ultimately leads to spin-dependent and anisotropic behaviour of the atom-light coupling. In the Cr case, the anisotropy comes from SOC in the excited states [255–257], while in the Ln case, it arises from the combined effect of SOC in both the ground and excited states [176, 211–215, 247, 258]. The existence of fine structure in both the excited and the ground state electronic configurations of Lns—and the large energy splitting within these fine structure manifolds—results in a different scaling of the vector and tensor polarisabilities versus that of the alkalies [157, 259]. That is, both the vector and tensor polarisabilities scale as $\propto \Delta_a^{-1}$, while in alkali metals the scaling is $\propto\Delta_a^{-2}$ for all detunings $\Delta_a$ greater than the hyperfine splitting, where $\Delta_a$ is the (average) detuning of the light from the excited states.

The large vector and tensor polarisabilities of magnetic atoms have been investigated via different means. In Cr, spin-dependent optical trapping was indirectly revealed in the study of many-body spinor dynamics, both in bulk [257] and in lattices [260]. Spin-dependent quadratic light shifts, resulting from the tensor polarizability close to the 427.6 nm-transition (laser at 427.85 nm), were also used to produce controllable spin mixtures [256].

A detailed experimental investigation of the characteristics of atom-light coupling in Lns using the newly available ultracold samples was crucial for a proper understanding of these systems. Indeed, up to now, theoretical predictions for Ln atom-light interactions have been difficult and remain incomplete. This is due to the complexity of the many-body electronic structure of the atoms and the limited knowledge available. Much progress has recently been made in developing sophisticated numerical tools to analyse the electronic level structure of Dy [211, 213], Er [212], Ho [258], and Tm [261], see also section 2.2.1. Recent measurements of Er, Dy, and Tm dynamic (total, scalar) polarisabilities based on trap frequencies measurements have resulted in close agreement between this theory and the experiments [176, 215, 261, 262]. In addition, thanks to a direct comparison with an alkali atom with a precisely known polarizability (potassium), unprecedented accuracy and precision on the Dy polarizability have been obtained in experiment, allowing theory test [176].

The unusually large vector and tensor polarisabilities of Lns (as compared to alkali atoms like Rb) have also been experimentally investigated. The large vector and tensor polarisability of Dy was first studied using a laser field close to the 741 nm transition [214]. Near this transition, a so-called 'tune-out' wavelength was found, at which the total light shift vanishes. Kao et al observed a strong dependence of the tune-out wavelength on the light polarisation angle, which is a consequence of the large tensor polarisability of Dy. This introduces a novel experimental knob with which to tune the wavelength where the AC polarisability vanishes. Becher et al [215] probed the large tensor part of the polarisability of Er ground-state, by testing the dependence of the optical trap's frequencies (thus its depth) on the light's polarisation axis. A sizeable anisotropy effect (few percent) was observed at conventional wavelengths (1570 nm, 1064 nm, and 532 nm) far away from atomic transitions. The spin-dependent light shifts were later used to rapidly control the spin dynamics in an assembly of lattice-confined fermionic Er [263], see also section 7.2. A similar scheme employing the ellipticity of light polarisation was applied to gases of Tm atoms and enabled the measurement of both the tensor and vector parts of its ground-state polarisability at 532 nm [262].

The locations of other special wavelengths at which the ground and excited state polarizibilities match have been predicted in Dy [211]. These so-called 'magic' wavelengths are of great utility to research involving atomic clocks and optical lattices for quantum simulation [264]. By making use of Dy's large tensor polarisability at the conventional trapping wavelength of 1070 nm, the authors [247] developed an analogous 'magic polarisation' scheme wherein the AC Stark shifts of the ground and excited states of a transition (here at 626 nm) are tuned to be identical. The authors then apply their findings to develop a Doppler cooling scheme on a trapped sample with improved efficiency. Large tensor light shifts close to the 626 nm transition of Dy were also used to generate effective spin-coupling terms [180, 181]; see also section 6.4.2. Other important consequences of these large vector and tensor polarisabilities in Lns relate to increasing the strength of Raman coupling schemes and to the possibility of efficient SOC; see section 6.4.1.

The unconventional features of atom-light coupling in these magnetic atoms come into play in various ways in this review; see in particular sections 6 and 7.2. Experiments have only begun to explore the possibilities that are offered by the magnetic atoms' electronic structure for the purpose of manipulating the properties of dipolar quantum gases.

The bosonic isotopes of Cr have a spherical orbital wave function in its ground state ($^7S_3$) and no nuclear spin. Thus, the short-range van der Waals scattering of two ground-state bosonic Cr atoms is isotropic and no interaction of hyperfine origin between collision channels arises. Consequently, a single anisotropic DDI potential explains the emergence of FRs in this case. The DDI couples channels with orbital angular momentum $\Delta \ell = 0,2$, and the total spin can vary by $\Delta S = 0,2$; see also sections 1.3.2 and 1.3.3. Thanks to its relative simplicity, the scattering features in this case are very well accounted for by multichannel scattering calculations.

The first observation of FRs between spin-polarised ⁵²Cr atoms in their absolute ground state (leading to a total spin of the entrance channel S = 6, $m_S = -6$) revealed the existence of 14 FRs between 0 and 600 G [195, 295]. Coupled-channel calculations show that the relevant closed channels are those with $\ell = 2,4$ and with $S = 2,4$ and 6. This means that second and fourth-order mixing are relevant. In addition, it was shown that essentially a single closed-channel bound state contributes to each FR. In this way, all except one of the FRs could be assigned. The DDI provides the necessary coupling terms, but only slightly affects the FR positions through the molecular potential $V_\mathrm{ch}(r)$ themselves. Finally, the background scattering length can also be theoretically estimated from the comprehensive coupled-channel calculations, giving $a_{\mathrm{bg}} = 112(14)a_0$.

Remarkably, two FRs, the ones with smallest B₀, are explained by a d-wave entrance channel resonant with a s-wave closed channel [195]. Such a d-wave entrance channel is intrinsic to the DDI coupling mechanism. One of these two FRs was later experimentally characterised [293]. The corresponding atom loss feature is asymmetric in B and shows an unconventional dependence with decreasing temperature, shifting to lower B and with decreasing width and amplitude; see figure 2(a). The atom losses are shown to occur via three-body recombination, but with a recombination coefficient depending linearly on the density n, in contrast to the usual n² dependence of three-body processes. This unusual behaviour is attributed to the slow tunnelling through the centrifugal barrier in the collision channel, and the adiabatic elimination of the fast collision with a third colliding atom. A precise analysis of this resonance yields a more accurate estimate of the scattering length: $a_{\mathrm{bg}} = 102.5(4)a_0$ [293].

Zoom In Zoom Out Reset image size

Besides the scattering between atoms in their lowest spin states, the analysis of Feshbach spectra and their comparison to coupled-channel calculations also provides information about the strengths of the spin-dependent scattering. Generally speaking, scattering strengths depend on the total spin of the pair of atoms colliding, which is preserved during the collision in absence of the DDI. The total spin $\mathcal{S}$ of the pair defines a spin channel, corresponding to a given molecular potential and thus to a distinct scattering length $a_{\mathcal{S}}$. The spin-polarised case of Cr corresponds to S = 6 and the background scattering length mentioned above is actually that of $a_{S = 6}$. In [195], additional spin channels were characterised, with $a_{S = 4} = 56\,a_0$, $a_{S = 2} = -7\,a_0$ [195]. Additional analysis of bulk spin relaxation dynamics resulted in a determination of $a_{S = 0} = 13.5\,a_0$ [196]; see also section 3. The large differences between these values provide evidence for the importance of spin-dependent scattering in Cr. Knowing this is important for investigating spinor physics with Cr; see sections 6 and 7.

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 [205, 206, 209]. 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 [205]. 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 [205], $\ell$ up to 10 were considered, and later, this was extended up to 20 [296]. An analytical model estimates that $\ell$'s up to ${\gtrsim} 40$ must be considered to reproduce experimental observations [296]. Needless to say, no perturbative treatment can be safely applied.

First experimental observations of FRs in magnetic Lns were performed on ¹⁶⁸Er and revealed an unusually large number of FRs within a narrow magnetic field range [172]. Soon after, theory work of [205] 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 [296]; see figure 3(a). 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 figure 3(b). 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 [297]. These measurements also revealed a much higher density for the fermionic isotopes as well as an intriguing temperature dependence of the Feshbach spectrum.

Zoom In Zoom Out Reset image size

The chaotic behaviour of both Er and Dy was thoroughly analysed and compared in the collaborative work of [294]. 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 [294]. 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 [293]; see figure 2(b). 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 [298]. In contrast to the Er system [294], 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 [294]. 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 [75]. 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 [296], where, e.g. a 3.5 G-wide FR is found at 57 G for ¹⁶⁸Er; see figure 3(a). A detailed study of such features was performed in ¹⁶⁴Dy [76] 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 figure 4. A similar study was later reported in ¹⁶²Dy [186]: 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 [298].

Zoom In Zoom Out Reset image size

2.4.2.1. Scattering length measurements.

Besides tuning of a, FRs provide the ability to associate a pairs of free atoms into a loosely bound molecule [20, 268, 289–291]. 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\gt B_0$ to $B\lt B_0$), one can then drive the population of the molecular state by (adiabatically) following the low-energy state [289, 317–321]. Alternative approaches use small-amplitude modulation of B at a frequency resonant with the bound-state energy for $B \lesssim B_0$ [322, 323], or even further away from the FR via radiofrequency coupling [324, 325]. 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 µ is also nearly doubled, resulting in a dipolar length ${a_\mathrm{dd}} \propto m\mu^2$ eight times longer.

An experiment using Er has demonstrated and studied the production of an ultracold sample of Er₂ Feshbach molecules using B-field ramps [75]. Pure molecular samples were produced by the resonant removal of the remaining free atoms. The gas contained about $2\times 10^4$ Er₂ at 300 nK with typical densities of ${\sim}8\times 10^{11}$ cm⁻³ and lifetimes of ${\sim} 10\,$ms. Ultracold samples of molecules in four different molecular states, associated with four different low-field FRs of ¹⁶⁸Er, were produced. The properties of the Feshbach molecules have been further measured, in particular the magnetic moment, which were found to be $\mu/\mu_{\mathrm{B}}\gt 11$ for three of the four states. By contrast, the bare atomic value is $\mu/\mu_{\mathrm{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 [293, 326]. 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.

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 [76, 186], see also section 2.4.2. The corresponding magnetic moments range from $\mu \approx 18\mu_{\mathrm{B}}$ to $\mu \approx 20\mu_{\mathrm{B}}$; see figure 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.

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_{\mathrm{c},l}(r) = \hbar^2l(l+1)/m r^2$. For the isotropic case, asymptotically, the outgoing wave differs by only a phase shift δ_l from the incoming wave. Therefore, the general scattering amplitude f decomposes into:

$$\begin{align} f(k,\theta) = \frac{1}{2ik}\sum_{l=0}^{\infty}(2l+1)\left(e^{2i\delta_l(k)}-1\right)P_l(\cos\theta) = \sum_{l=0}^{\infty}f_l(k,\theta), \end{align} \tag{ 25 }$$

where θ describes the angle between k and r and varies between 0 and π. Here, P_l denote the ordinary Legendre polynomials. The scattering cross section then reads:

$$\begin{align} \sigma(k) &= \sum_{l=0}^{\infty}\sigma_l(k), \end{align} \tag{ 26 }$$

$$\begin{align} \sigma_l(k) &= \frac{4\pi}{k^2}(2l+1)\sin^2[\delta_l(k)]. \end{align} \tag{ 27 }$$

For an interaction of short range b, we can define the temperature below which the system is in the ultracold collisional regime as $k_BT\ll \hbar^2/m b^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.

3.1.2.1. Scattering length.

3.1.2.2. 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 section 1, 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\lt (n - 3)/2$ and as $k^{\,n-2}$ otherwise [78, 79]. 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_BT\sim \hbar^2/m b^2$. Associating b with the van der Waals length $r_\mathrm{vdW} = (m C_6/\hbar^2)^{1/4}$, we find that the ultracold limit is typically reached for temperatures of hundreds of ${\mu\mathrm{K}}$.

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

$$\begin{align} \sigma_l(E) \propto k^{4l}\propto E^{2l}, \end{align} \tag{ 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 (section 3.2) that all partial waves contribute.

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 equations (25) and (26). These expressions may be rewritten as:

$$\begin{align} f(k,\theta) = \frac{1}{ik}\sum_{l=0}^{\infty}\epsilon(l)(2l+1)\left(e^{2i\delta_l(k)}-1\right)P_l(\cos\theta), \end{align} \tag{ 29 }$$

where θ 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:

$$\begin{align} \sigma(k) = \sum_{l=0}^{\infty}\sigma_l(k), \end{align} \tag{ 30 }$$

$$\begin{align} \sigma_l(k) = \frac{8\pi}{k^2}\epsilon(l)(2l+1)\sin^2[\delta_l(k)]. \end{align} \tag{ 31 }$$

For van der Waals interactions, the Wigner threshold law, equation (28), implies that the cross section for the lowest partial wave accessible to fermions, l = 1, vanishes as E², while for bosons it tends towards 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.

The DDI in equation (2) 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,\pm2$ 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 [85, 161, 343]. The appropriate expression, neglecting all other interatomic potentials [161], but including direct and exchange terms, is:

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

where $\widetilde{V_{\mathrm{dd}}}(\boldsymbol{k})$ is the Fourier-transformed matrix elements of equation (2) and $\tilde{\epsilon} = [1,0,-1]$ for identical bosons, distinguishable particles, and identical fermions, respectively [16, 85]. The $\boldsymbol{k_\mathrm{i}}$ and $\boldsymbol{k_\mathrm{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 [85, 161]:

$$\begin{align} \sigma_{0} & = \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)], \end{align} \tag{ 35 }$$

$$\begin{align} \sigma_{1} & = \frac{8\pi}{15}F^{3}\left(\frac{\mu_{0}(g_{F}\mu_{B})^{2}m}{4\pi\hbar^{2}}\right)^{2}\left[1 + \tilde{\epsilon} h\left(\frac{k_\mathrm{f}}{k_\mathrm{i}}\right)\right]\frac{k_\mathrm{f}}{k_\mathrm{i}} , \end{align} \tag{ 36 }$$

$$\begin{align} \sigma_{2} & = \frac{8\pi}{15}F^{2}\left(\frac{\mu_{0}(g_{F}\mu_{B})^{2}m}{4\pi\hbar^{2}}\right)^{2}\left[1+ \tilde{\epsilon} h\left(\frac{k_\mathrm{f}}{k_\mathrm{i}}\right)\right]\frac{k_\mathrm{f}}{k_\mathrm{i}}. \end{align} \tag{ 37 }$$

In these expressions, σ₀ is the elastic cross section discussed in the previous section; σ₁ and σ₂ are the cross sections for single spin-flip and double spin-flip dipolar relaxation processes, respectively. The kinematic terms $\frac{k_\mathrm{f}}{k_\mathrm{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}$ [85, 161, 339].

The inelastic cross sections in equations (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 [161] 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 equations (36) and (37) are good approximations at low magnetic field for which $k_\mathrm{f}a\ll1$, the effect of the molecular potential becomes important when B is such that $k_\mathrm{f}a\gtrsim 1$. This leads to deviations from universal behaviour (see below).

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_\mathrm{rel} \rangle_{\textrm{th}}$, where $v_\mathrm{rel} = 2\hbar k_{i}/m$ describes the initial relative velocity between a pair of particle and $\langle . \rangle_{\textrm{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_\mathrm{rel} \rangle_{\textrm{th}}$.

Following earlier studies of weakly magnetic atoms [344–346], dipolar relaxation was quantitatively studied in bosonic Cr [85, 347]. In these studies, ⁵²Cr 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 [85]. Moreover, the data were in reasonable agreement with cross section predictions based on the first-order Born approximation given above, i.e. equations (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 [161] and Dy [164] at lower temperatures (100s of nK) were performed using optical trapping [348]. In optical traps, because of limited trap depth, the particles involved in dipolar relaxation are typically lost (except at very low magnetic field, see section 6 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 σ₁ for Dy at B = 0 can be up to 100×-larger than that of the alkali metal caesium. In addition, [164] pointed out that, while the ratio of $\sigma_1/\sigma_2 = F^{-1}\ll1$ 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, σ₂ for these states can be larger than σ₁ for the stretched states: one cannot avoid dipolar relaxation simply by employing spin states near $m_F = 0$.

While the elastic dipolar interaction is long-ranged, [161] 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 equation (34) for σ₀ involves incoming and outgoing waves oscillating at the same frequency ($k_\mathrm{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_\mathrm{dr}$, scales as $1/k_\mathrm{f}$. This is set by the release of Zeeman energy in the spin relaxation. Hence, the range of dipolar relaxation decreases with B following:

$$\begin{align} R_\mathrm{dr}(B) \approx \frac{\hbar}{\sqrt{m g_F \mu_{\mathrm{B}} B}}. \end{align} \tag{ 38 }$$

In terms of partial-wave decomposition, dipolar relaxation occurs at the interparticle distance $R_\mathrm{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 [161]. 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_\mathrm{dr}$ at which dipolar relaxation occurs. This provides a factor of 20 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 ⁵²Cr and it cannot be directly extended to the other highly magnetic atoms due in particular to the complexity of their FR spectrum.

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 section 3.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 equations (36) and (37). At high B and low T, the kinematic terms diverge as $\langle \frac{k_\mathrm{f}}{k_\mathrm{i}} \rangle_\mathrm{th} \propto \sqrt{B/T}$, where the $\mathrm{th}$ denotes a thermal average. In this limit, σ₁ 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: (a) Dipolar relaxation rates grow worse for identical bosons versus B, e.g. σ₁ for bosonic Dy is 10³ times larger than Cs's at only a 1 G field; (b) $\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 [164]. 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. [164] 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 [161]. 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 section 4.2) [13]. The stability of these states was also exploited for the long-lived implementation of artificial SOC in Dy gases; see section 6.4.1.

Zoom In Zoom Out Reset image size

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 equation (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 [161–163] for the case of 1D, 2D, and 3D optical lattices. In [161], 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 ⁵²Cr [162]. 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 figure 6. In the specific case of a deep 3D optical lattice [163], 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 SOC in the lattice [163, 349].

Zoom In Zoom Out Reset image size

Modifications to dipolar relaxation rates were also demonstrated in [161] 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.

Section 6 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_\mathrm{B}$100 nK. This induces spontaneous, incoherent demagnetisation from the absolute ground state. Further control of the field down to the 100 µ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 [350]. 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_\mathrm{dr} \gt \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 [351–353]. Strongly rotating and interacting Fermi gases approaching a Laughlin state may appear [354].

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 [164]. For example, long-lived mixtures of $|F,-(F-1);F,-F\rangle$ were achieved in fermionic Dy and Er [164, 313] (see section 6.2.4) and the same set of states was used to generate SOC in fermionic Dy [251]; see section 6.4.1. The study of out-of-equilibrium dynamics is also possible for time-scales up to the dipolar relaxation time scale [257, 355]; see section 6.2.3.

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 λ 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 [299]. One then observes the time for the temperature in the third direction to equilibrate. The time constant τ for this typically exponential relaxation process may be related to the cross section via $\tau = \alpha/n v \sigma$, where n is the average number density, v is the mean relative speed, and σ is the elastic cross section. The anisotropy of the differential cross section appears through the parameter α, which is the number of collisions, on average, for rethermalisation. For s-wave (p-wave) collisions, α = 2.5 (25/6 ≈ 4.17) [300]. However, for dipolar particles, α becomes angle dependent due to the anisotropic differential cross section [301]. Its value is plotted in figures 8(a) and (b) versus the angle λ between the axis along which energy is imparted and $\hat{\boldsymbol\varepsilon}$. We see that α can change by more than a factor two versus angle for fermions, while slightly less than two for bosons. The control of λ 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^∘).

Zoom In Zoom Out Reset image size

The large anisotropy in the rethermalisation of dipolar gases was experimentally demonstrated in ultracold gases of fermionic ¹⁶⁷Er [302]. For fermionic atoms, σ is simply given by the universal formula equation (33) so that a measurement of τ directly provides α. Figure 8(c) shows measurements of α in agreement with the dependence predicted in [301]. The authors also reported the reduction in rethermalisation rate due to Pauli-suppression of scattering once T was sufficiently lower than T_F (see section 3.2 and equations (32) and (33)). Unlike α, the suppression factor was observed to be independent of λ, 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 ¹⁶⁴Dy, ¹⁶²Dy [303] and of ¹⁶⁴Er, ¹⁶⁶Er, ¹⁶⁸Er and ¹⁷⁰Er [75]. These studies measured the unknown values of the scattering lengths a of these atoms by using either the theoretically predicted α or by using more extensive dynamics simulations to convert the measured τ into values of a; see equation (32) and section 2.4.2.1 for more details. In [303], the authors used a thorough theoretical analysis informed by data from different values of the angle λ to extract a. This study confirmed the anisotropic behaviour as well as provided evidence of hydrodynamical effects in the rethermalisation dynamics [356].

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 [356], 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 [303]. Nevertheless, cross sections were able to be extracted from the data through close comparison to Monte Carlo simulations. Section 2.4.2.1 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 [305]; see also section 4.1.2.2.

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. [165] describes such an experiment: one creates two dipolar BECs of ¹⁶²Dy counter-propagating at momentum $\pm2\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 figure 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.

Zoom In Zoom Out Reset image size

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 [357–360]—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 η, as expected from equation (39). The shapes of the halos in the data agree with Monte Carlo simulations based on equation (39) and described in [356].

4.1.1.1. Interaction potential in dilute gases.

As introduced in section 1.4.1, polarised dipolar bosons in an ultracold dilute ensemble can be considered to be interacting via a potential of the form:

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

Again, a is the s-wave scattering, while ${a_\mathrm{dd}} = \mu_0m\mu^2/12\pi\hbar^2$ is the dipolar scattering length. The interparticle separation vector is $\boldsymbol r$, with norm r and angle θ 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

$$\begin{align} V_{\mathrm{tr}}=\frac m2\left(\omega_\rho^2\,\rho^2+\omega_z^2\,z^2\right). \end{align} \tag{ 44 }$$

The trap aspect ratio is:

$$\begin{align} \lambda=\omega_z/\omega_\rho. \end{align} \tag{ 45 }$$

The validity of equation (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 [90]) and the DDI potential equation (11) has long been debated [9, 16]. 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 [91, 92, 339, 361–363]. While it is a good approximation for weak dipoles and away from scattering resonances, modifications may be necessary for systems with large $\varepsilon_{\mathrm{dd}}$. In particular, the s-wave scattering length may depend on d [91, 92, 364], requiring a renormalisation of the dipolar length in equation (43) [340, 365]. In the case of the most magnetic of the atoms, Dy, with $\varepsilon_{\mathrm{dd}} \sim 0.5$–2, corrections to ${a_\mathrm{dd}}$ were predicted to be of 10% at temperature of 100 nK and only 2% at 10 nK [366].

4.1.1.2. 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 [9, 87–89]. 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{n{a_\mathrm{dd}}\,^3}\,\ll1$. A proper treatment of fluctuations around MF values, outlined in section 5, leads to a more restrictive condition: $\sqrt{na^3}\ll1,\varepsilon_{\mathrm{dd}}\lesssim 1$. Note that the length ${a_\mathrm{dd}}$ is defined such that $a = {a_\mathrm{dd}}$ ($\varepsilon_{\mathrm{dd}} = 1$) marks the limit of the mechanical stability of a homogeneous isotropic dBEC, and other authors sometimes call $3{a_\mathrm{dd}}$ the dipolar length (which arises from a consideration of the two-body problem [339]); see section 1.3.4.

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

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

with $g = 4\pi\hbar^2a/m$ and $\Phi_{\mathrm{dd}}(\boldsymbol r) = \int d\boldsymbol{r^{^{\prime}}}U_{\mathrm{dd}}(\boldsymbol{r^{^{\prime}}}-\boldsymbol{r})|\psi(\boldsymbol{r^{^{\prime}}})|^2$; see also equation (10).

The corresponding energy functional is:

$$\begin{align} E[\psi]&=\int \left[-\frac{\hbar^2}{2m}|\nabla\psi|^2+V_{\mathrm{tr}}(\boldsymbol{r})|\psi|^2\right.\nonumber\\&\qquad\qquad\left. +\frac g2|\psi|^4+\frac12|\psi|^2\Phi_{\mathrm{dd}}(\boldsymbol r)\right]dr, \end{align} \tag{ 47 }$$

which can be minimised to obtain the equilibrium ground state of the BEC [13, 14, 16, 367]. In these equations, the DDI yields a nonlocal contribution $\Phi_{\mathrm{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 [365, 368].

4.1.2.1. Magnetostriction of dBECs.

When placed in a cylindrical harmonic trap, a nondipolar gas ($\varepsilon_{\mathrm{dd}} = 0$) assumes the geometry of the trap. A dipolar system behaves differently because of the presence of $\Phi_{\mathrm{dd}}$ in equation (9). Specifically, the dipole term $\Phi_{\mathrm{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 [369]. 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.

• Application of the Thomas-Fermi approximation

Rather than numerically solving equation (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 [9]. This allows one to neglect the kinetic energy term in equation (9). For $\varepsilon_{\mathrm{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$ [9]. In the dipolar case, [370, 371] showed that a dBEC retains the parabolic shape in this (TF) limit. However, its aspect ratio:

$$\begin{align} \kappa=R_\rho/R_z, \end{align} \tag{ 48 }$$

defined by the ratio of the TF radii perpendicular (R_ρ ) 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_{\mathrm{dd}}$ to the mean energy in equation (47), which one can calculate analytically:

$$\begin{align} E_{\mathrm{dd}}=-\frac47\varepsilon_{\mathrm{dd}}\,g\,\frac{n_0\,N}{2}f_{\mathrm{dip}}(\kappa), \end{align} \tag{ 49 }$$

where n₀ is the density at the trap centre. The function $f_{\mathrm{dip}}(\kappa)$ encompasses the anisotropy of the interaction energy, it is a decreasing function of κ [372]. The complete expression is [373]:

$$\begin{align} f_{\mathrm{dip}}(\kappa)=\frac{1+2\kappa^2}{1-\kappa^2}-\frac{3\kappa^2\mathrm{arctanh}\sqrt{1-\kappa^2}}{(1-\kappa^2)^{3/2}}. \end{align} \tag{ 50 }$$

Note that an angular integration of the DDI over an isotropic distribution is zero, leading to $f_{\mathrm{dip}}(1) = 0$. $f_{\mathrm{dip}}$ is bounded by two limits: fully collinear dipoles (attractive DDI) $\lim_{\kappa\to 0}f_{\mathrm{dip}} = 1$, while side-by-side dipoles (repulsive DDI) yield $\lim_{\kappa\to\infty}f_{\mathrm{dip}} = -2$. Thus, one readily sees that $E_{\mathrm{dd}}$ is reduced by lowering κ. In the TF approximation, minimising E leads to a transcendental equation for κ [92, 371, 373]:

$$\begin{align} \lambda=\kappa\left(\frac{1+2\varepsilon_{\mathrm{dd}}-\frac{3\varepsilon_{\mathrm{dd}}f_{\mathrm{dip}}(\kappa)}{1-\kappa^2}}{1-\varepsilon_{\mathrm{dd}}+\frac{\kappa^2}{2}\frac{3\varepsilon_{\mathrm{dd}}f_{\mathrm{dip}}(\kappa)}{1-\kappa^2}}\right)^{1/2}. \end{align} \tag{ 51 }$$

The solution of this equation for a given λ thus gives the degree of magnetostriction. A solution exists for any λ provided that $\varepsilon_{\mathrm{dd}}\leqslant1$. 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 [91, 92, 96, 371, 374]. This approximation will be later reviewed in section 4.1.4.1. Interestingly, and despite the difference in the approximations, such a Gaussian variational ansatz for ψ yields the same equation as equation (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 sections 4.1.4.3 and 5.4. In this paragraph, we focus on the shape of stable dBECs.

• 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 [184] and in-situ magnetostriction images first presented in [375]. An example is shown in figure 10.

• Magnetostriction in the TOF dynamics of BECs

Zoom In Zoom Out Reset image size

Observing the (hydrodynamic) expansion of the gas trap release is a common way to reveal BEC physics [376, 377], and in particular the impact of interactions on quantum gases. The hydrodynamic equations for the expansion of a (nondipolar) BEC were first presented in [378], showing inversion of ellipticity during TOF. Solutions to the hydrodynamic equations for dipolar gases can also be found [370, 371, 379]. 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 ⁵²Cr BEC [380]. 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_{\mathrm{dd}}$ for ⁵²Cr away from FRs was extracted [381]. The value agrees with most precise values of the background scattering length $a_{\mathrm{bg}}$ obtained more recently [161]; see section 2.4. Using a FR to lower a—and thus to enhance dipolar effects into the $\varepsilon_{\mathrm{dd}}\gt 1$ regime—magnetostriction in situ and in TOF was strong enough so that a complete suppression of ellipticity inversion was observed [382]. The variation of the aspect ratio of an expanding cloud with $\varepsilon_{\mathrm{dd}}$ can be seen in figure 11.

Zoom In Zoom Out Reset image size

4.1.2.2. 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_\mathrm{dip}/E_\mathrm{k}$. For a non-degenerate gas, $E_\mathrm{dip}\sim S^2d^{\,2}n$ and $E_\mathrm{k}\sim k_\mathrm{B}T$, such that $\eta\propto(n\lambda_\mathrm{th}^3)\times(k_\mathrm{th}{a_\mathrm{dd}})$, with $k_\mathrm{th} = 2\pi/\lambda_\mathrm{th} = \sqrt{mk_\mathrm{B}T}/\hbar$ the thermal wavenumber, and $\lambda_\mathrm{th}$ the thermal wavelength. To maximise η in a thermal gas, we may choose Dy, which provides a ${a_\mathrm{dd}} = 131\,a_0$; see section 2). At typical experimental densities for a gas with a temperature just above its BEC T_c , η remains $\ll1$. 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.

Tang et al [305] theoretically and experimentally studied the anisotropic expansion of a thermal dipolar Bose gases of ¹⁶⁴Dy and ¹⁶²Dy 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 [16]; see figure 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 section 2.4.2.1 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_\mathrm{th}^3)\times(k_\mathrm{th}{a_\mathrm{dd}})$, as expected, and arises primarily from the two-body collisions and the direct Hartree contribution to the MF energy.

Zoom In Zoom Out Reset image size

A similar effect arises in DFGs, 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 section 1.4.2 and the discussion in section 4.2.

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_{\mathrm{dd}}$ are possible. Besides, stringent modifications of the excitation spectrum at $\varepsilon_{\mathrm{dd}} \sim 1$ were also evidenced.

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

• Continuous Bogoliubov spectrum

Elementary excitations of BECs were introduced in section 1.4.1.3 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 equations (8) and (16), the spectrum of a uniform 3D dBEC of density n is [374]:

$$\begin{align} \hbar\,\omega(\boldsymbol k)=\sqrt{\frac{\hbar^2k^2}{2m}\left(\frac{\hbar^2k^2}{2m}+2\,gn\,\left(1+\varepsilon_{\mathrm{dd}}(3\cos^2\theta_k-1)\right)\right)}, \end{align} \tag{ 52 }$$

where θ_k is the angle between the direction of excitation propagation and the dipole orientation. The sound velocity for BECs with only contact interactions ($\varepsilon_{\mathrm{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_{\mathrm{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_{\mathrm{dd}}\geqslant1$; see section 4.1.4.

• Bragg spectroscopy

Using the standard method of two-photon Bragg spectroscopy (see, e.g. [383–387]) on a dBEC of ⁵²Cr, Bismut et al [388] 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 equation (52) (which holds in the homogeneous case). This limit comes from the fact that the BEC smallest size ($R_\mathrm{TF\,min}$) must be larger than the excitation wavelength: $k\,R_\mathrm{TF\,min}\gg1$. In this limit, the LDA 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 figure 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.

• Critical superfluid velocity measurement

Zoom In Zoom Out Reset image size

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_\mathrm{c} = \underset{\boldsymbol{k}}{\mathrm{min}}\,(\frac{\omega_{\boldsymbol{k}}}{\hat{\boldsymbol{v}}\cdot\boldsymbol{k}})$ [100, 389]. In the homogeneous case discussed above, equation (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 [389, 390]. 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 [389]. The critical velocity is thus systematically smaller than the speed of sound in the direction of motion. For anisotropic confinement, as discussed in section 4.1.3.3, the critical velocity should also be affected by the existence of low-energy high-momentum modes, such as the roton mode [390, 391].

Superfluid velocity measurements were performed on a ¹⁶²Dy dBEC using a local light defect driven linearly along one axis [375]. 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 figure 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.

Zoom In Zoom Out Reset image size

4.1.3.2. 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 [92, 96, 370, 392, 393]. MF methods to extract the collective modes nature and frequency rely on expanding the energy functional equation (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 sections 4.1.2 and 4.1.4.1.

Experimentally, Bismut et al [394] studied the influence of the DDI on such modes with a ⁵²Cr 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 figure 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.

Zoom In Zoom Out Reset image size

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

As highlighted in equation (16), the DDI contributes to the continuous excitation spectrum via an effective coupling strength $\widetilde{V_{\mathrm{dd}}}({\boldsymbol{k}})$, which adds to the constant contact contribution g. As observed in equation (52), in the 3D homogeneous case, the DDI contribution, simply given by its Fourier transform (equation (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_{\mathrm{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 equation (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 [100–102]. Here, roton excitation were first speculated to explain the exotic macroscopic properties of this superfluid [100, 101, 395], long before their observation [102]. 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) [396–398] (although one should note that the phase transition to solid helium is not due to roton softening).

• 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 [97] and, shortly after, in that of magnetic or electric dipolar system by Santos et al [98]. These seminal works consider semi-infinite trapping geometries, i.e. infinite along one ($\omega_z = 0$) [97] or two ($\omega_\rho = 0$) [98] 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_{\mathrm{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_{\mathrm{dd}}}({\boldsymbol{k}})$ becoming attractive at large k. In the quasi-infinite geometries of [97, 98], 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_{\mathrm{dd}}}(\boldsymbol k)\gt 0$, while, for $k\gtrsim 1/\ell_z$, $\widetilde{V_{\mathrm{dd}}}(\boldsymbol k)\lt 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^2k^2}{2m}$ which ultimately dominates at very large k, the effect of $\widetilde{V_{\mathrm{dd}}}(\boldsymbol k)\lt 0$ is the strongest at $k\sim 1/\ell_z$ and, for weak enough s-wave coupling strength g, a minimum arises at $k_{\mathrm{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 Δ can be lowered by increasing the density or increasing $\varepsilon_{\mathrm{dd}}$. When $\Delta = 0$, [98] finds $k_{\mathrm{rot}} = \sqrt{2}/\ell_z$, independent of the density and the interaction parameters $g,\,\varepsilon_{\mathrm{dd}}$.

• Roton excitation spectrum in anisotropic fully trapped dBECs

Numerous subsequent theoretical works describe the roton in dipolar gases confined in finite geometries (equation (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 [99, 112–114, 391, 399–409]. In the case where the trap is tightly confining along the dipoles and sufficiently anisotropic (e.g. $\lambda\gg1$), a roton-like spectrum arises for large enough $\varepsilon_{\mathrm{dd}}$. The relevant range of $\varepsilon_{\mathrm{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 ω [400]. 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 [401]. This confinement effect can be understood, within the LDA, as the roton energy decreasing with increasing density. This translates into a finger-like feature in the (discrete) dispersion relation [406]. 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 [402, 405, 407], strong correlation emerging from an interaction quench [409], and peculiar structures in the collapse dynamics [113, 114]; see also section 4.1.4.3.

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 [112, 410]; see also discussion in section 4.1.4.3. 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 [399].

• 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 section 5. This leads to a fast change of the BEC's properties, including its excitation spectrum. The roton mode was experimentally observed in a ¹⁶⁶Er BEC [311, 312] 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 [311], relying on instability dynamics, will be described in section 4.1.4.3. In [312], Petter et al [312] reported on roton spectrum measurements based on Bragg spectroscopy, similar to section 4.1.3.1. 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 figure 16(a). When increasing $\varepsilon_{\mathrm{dd}}$, a preferential softening of $\epsilon(k)$ at large k is observed, finally forming a minimum at $k = k_{\mathrm{rot}}\sim 1.3/\ell_z$; see figure 16(a). The minimum is observed only in a narrow range of scattering lengths and the roton gap shows a fast decrease with $\varepsilon_{\mathrm{dd}}$ towards instability ($\epsilon(k_{\mathrm{rot}}) \approx 0$), see figure 16(b). The Bragg measurements also proved, through the increase of the Bragg response, the enhancement of density-density correlations at $k_{\mathrm{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 [411], following early theoretical studies [402, 405, 407]. 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 [412]. Here 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 [406, 410] and reveal as such in the density fluctuations.

Zoom In Zoom Out Reset image size

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.

4.1.4.1. Global mechanical stability for a single BEC.

Three-dimensional, homogeneous BECs under contact interactions are mechanically stable for positive compressibility [413] $\chi = \frac{1}{nmc_0^2}>0$ (where c₀ is the speed of sound; see section 4.1.3.1), they are thus unstable for a < 0 [9, 109]. Finite-sized [414] and harmonically trapped [3, 105] 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_{\mathrm{dd}} = 1$ as the sound velocity perpendicular to the dipole orientation then cancels, $c^2|_{\theta_k = \pi/2} = 0$; see section 4.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.

• 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_\mathrm{MF}\sim-n^2$ and is thus minimised for infinite density. Using the TF approximation and neglecting the kinetic energy contribution, $E_\mathrm{MF}\lt 0$ leads to a singularity of the ground state density, thus giving an instability. From equation (12), the sign of the MF dipolar energy of a cylindrical dBEC depends on the cloud aspect ratio κ. Taking into account both interactions, one can show that the total MF interaction energy scales like:

$$\begin{align} E_\mathrm{MF}\sim gn_0N(1-\varepsilon_{\mathrm{dd}} f_{\mathrm{dip}}(\kappa)). \end{align} \tag{ 53 }$$

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

• The Gaussian ansatz for a BEC near collapse

Zoom In Zoom Out Reset image size

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 equation (47). The variational Gaussian ansatz provides a semi-analytical approach:

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

inserted into equation (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:

$$\begin{align} E_\mathrm{k} = N\,\frac{\hbar^2}{4\,m}\sum_i\frac{1}{\sigma_i^2}, \end{align} \tag{ 55 }$$

for the kinetic energy,

$$\begin{align} E_\mathrm{trap} = N\,\frac m4\sum_i\omega_i^2\sigma_i^2, \end{align} \tag{ 56 }$$

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

$$\begin{align} E_\mathrm{MF}=N^2\,\frac{g}{2(2\pi)^{3/2}\bar\sigma^3}(1-\varepsilon_{\mathrm{dd}}f_{\mathrm{dip}}(\kappa_x,\kappa_y)), \end{align} \tag{ 57 }$$

where we allowed non-axisymmetric traps. The expression for the generalised $f_{\mathrm{dip}}$ function in terms of the two aspect ratios $\kappa_{x,y} = \sigma_{x,y}/\sigma_z$ was first worked out in [416]. One can then minimise the total energy $E = E_\mathrm{k}+E_\mathrm{trap}+E_\mathrm{MF}$ to find the ground-state size $\sigma_{x,y,z}$. In general, for $\varepsilon_{\mathrm{dd}}>1$, there is always a global singularity with $\kappa_{x,y}\ll1$ 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. figure 17(a). At low trap aspect ratios $\lambda\ll1$, the kinetic energy acts to stabilise the gas against collapse, as in the contact case. However, for $\lambda\gg1$ the effect is opposite, kinetic energy increases the smallest size, which adds a little attraction via dipolar effects and destabilises the gas.

• First measurements of the instability threshold

The instability threshold was studied in ⁵²Cr BECs of about 25 000 atoms in tunable traps [415]: When quenching a downward using a magnetic FR, the authors observed an abrupt and quick disappearance below a critical value, $a_\mathrm{crit}$, of the BEC peak in TOF images. The experimental values of $a_\mathrm{crit}$ as a function of λ were found to agree well with the Gaussian ansatz predictions; see figure 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 [415].

4.1.4.2. 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 (NNI) can modify the stability diagram. Indeed, the first observation of nearest-neighbour dipolar effects was reported in [417] 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 [418] have observed this effect with a ⁵²Cr 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 figure 18.

Zoom In Zoom Out Reset image size

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 [420]. 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 section 5. Simulations of the GPE confirmed that no atomic losses were expected in trap, but were induced by the deconfinement and following collapse.

4.1.4.3. Modulational instabilities

• 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 [98] predicting a roton-type dispersion relation in anisotropic semi-infinite dBECs, see section 4.1.3.3 and [99]. Tuning the parameters of the dBEC (i.e. its scattering length and trapping geometry) can lead to a softening of this roton minimum (ω 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 section 4.1.3.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 [112, 421, 422]. In this case, an instability can be triggered when $\varepsilon_{\mathrm{dd}}>1$. It differs from a global instability in that several local attractors gather the atomic density leading to an ensemble of local collapses [113, 114, 410]. This kind of instability can be seen as the softening ($\omega^2\lt 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 [406, 410]. 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_{\perp}\lesssim1$ where $R_{\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 [114, 410, 412, 421, 422]. 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 [423] and investigated in this case [424, 425]. 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 ω). The modulational instability here arises when other modes are also unstable and have a larger growth rate. The mode with the highest growth rate has a finite k, which is then favoured rather than global collapse [426, 427].

• Experimental evidence for local collapse

In dBECs, the role played by local collapse was first noted by an extensive theoretical analysis of the data from [415] by Wilson et al [410], going beyond the Gaussian ansatz (equations (54)–(57)) and solving exactly the GPE (equation (9)). Their predictions, showing a better agreement with the measured stability threshold at large trap anisotropy with $\lambda\gtrsim 3$ (figure 17(b)), implies the occurrence of local collapse. Yet the instability mechanism was not experimentally resolvable and later studies of the collapse dynamics also let modulational instabilities remain elusive, see section 5.1 [428].

More recently, the impact of such instabilities was experimentally studied on finite-size anisotropically-trapped dBECs of ¹⁶⁴Dy, either in pancake traps [184] or in cigar shaped ones [429, 430] with a tight confinement along the dipoles. Here, after quenching a down, remarkable long-lived in-situ density structures were observed. The observed density structures have been attributed to the occurrence of a modulational instability following the quench [184]. Because of the limited system size of the original pancake-shaped geometry in [184], the instability is expected to be driven by an angular mode with $\langle k\rangle\times R_\mathrm{BEC}\sim 1$ and m > 1 [431], as also experimentally evidenced in a recent set of experiment [412]. Following their first observation, the authors of the above-cited experimental works characterised the absence or existence of a modulational instability as a function of the trap geometry via its signature on the final density distribution; i.e. the observation of single or multiple droplets long after the quench [432]. They identify a critical trap aspect ratio $\lambda = 1.87(14)$ (compared to the dipole orientation, z) above which the modulational instability exists; see also section 5.3.1. The instability itself was not experimentally investigated in these works as they focused on the long-time behaviour, which results from a subsequent intricate non-linear dynamics loosing track of the unstable mode driving the collapse. The authors found that these final density distributions were surprisingly stable, which revealed an unpredicted stabilisation mechanism. This discovery set the ground for a new paradigm of quantum fluids, which will be discussed in the next section 5.

• Experimental investigation of the roton instability

As introduced in section 4.1.3.3, Chomaz et al [311] first observed the roton excitation by probing the instability dynamics of a large, cigar-shaped dBEC of ¹⁶⁶Er with transverse magnetisation. Chomaz et al performed a fast interaction quench and studied the short-time evolution of dBEC. The authors reported the transient appearance of remarkable structure in the momentum distribution of the dBEC with a high-amplitude central peak framed along the cigar-long axis by two lower-amplitude symmetric side peaks, see figure 19(a). Based on a Bogoliubov picture relevant for the short-time dynamics, this was interpreted as the coherent population of finite-momentum excitation modes thanks to its privileged dynamical softening to imaginary energies, indicative of an unstable roton mode. By further studying the peak position and the time evolution of its population, the authors demonstrated the characteristic scalings $k_{\mathrm{rot}}\sim 1/\ell_z$ and $\omega_\mathrm{rot}\propto (a-a^*)^{1/2}$ for the unstable regime ($a\lt a^*$); see figures 19(b) and (c). These observations are in agreement with theory predictions based on an analytical model as well as on GPE simulations.

Zoom In Zoom Out Reset image size

The phenomenology of MF unstable, strongly dipolar BECs of Er and Dy is fundamentally modified by the BMF effects. These effects can affect not only the dynamics of the BEC when driven to an unstable regime but also modify its stability itself, as it introduces a stabilisation mechanism. Indeed, the LHY correction of equation (60) provides an additional conservative potential, which acts as an effective higher-order (in 3D, effective 5/2-body interaction) interaction term in the GPE and which remains repulsive within its whole validity range [461]. This higher-order repulsive term has a density dependence of a higher power than MF (2-body interaction). In the attractive MF regime, it thus could be thought to stabilise a high-density state for which BMF repulsion compensates MF attraction.

Let's first consider the simplest case of single-peaked ground state—standard BEC in the MF stable regime and later called 'droplet' in the MF unstable regime. Using a Gaussian ansatz for its wave function ψ (see also equations (54)–(57)), one can calculate the MF and BMF contributions to the ground-state energy per particle and single-out the key ingredients that lead to stability:

$$\begin{align} \frac{E_\mathrm{MF}}{N}&=\frac{1}{2^{3/2}}\,\frac g2\,n_\mathrm{c}\,(1-f_\mathrm{dip}(\kappa)\,\varepsilon_{\mathrm{dd}}), \end{align} \tag{ 61 }$$

$$\begin{align} \frac{E_\mathrm{BMF}}{N}&=\left(\frac25\right)^{3/2}\,\frac g2\,n_\mathrm{c}\,\frac{128}{15\sqrt \pi}\,\sqrt{n_\mathrm{c} a^3}\,Q_5(\varepsilon_{\mathrm{dd}}), \end{align} \tag{ 62 }$$

where $n_\mathrm{c}$ is the central density. For $f_\mathrm{dip}(\kappa)\varepsilon_{\mathrm{dd}}>1$, $E_\mathrm{MF}\lt 0$ while $E_\mathrm{BMF}>0$. $E_\mathrm{BMF}$ has a stronger dependency in the density by a power $n_\mathrm{c}^{1/2}$ compared to $E_\mathrm{MF}$, therefore, by increasing $n_\mathrm{c}$, the ground-state energy can be minimized at a finite value, corresponding to a finite density. This stabilisation mechanism, relying on the mere effect of quantum fluctuation, is what we coin dipolar quantum stabilisation.

Remarkably, in presence of MF attraction, the single-peaked state considered above can be stabilised even in the absence of trapping, by competing MF attraction and the BMF repulsion, the equilibrium between these two forces fixing the peak density $n_\mathrm{c}$ (MF attraction dominates at low density, and BMF repulsion at high density).

Such a stabilisation mechanism is reminiscent of a liquid phase of matter. In ordinary liquids, the weak attraction on large distance (e.g. of van der Waals or covalent types) is counterbalanced by a strongly repulsive core arising from the electromagnetic forces and the Pauli exclusion between the atoms' (or molecules') electron cloud. This stabilises the liquid at a large density at which the repulsion becomes effective $n_\mathrm{c}{a_0^3}\sim{1}$. The stabilisation resulting from competing MF and BMF effects in strongly dipolar gases can similarly be seen as resulting in a liquid state.

Considering for now such an untrapped system, and neglecting kinetic energy, the equilibrium peak density can be estimated by imposing $\frac{\partial E}{\partial n_\mathrm{c}} = 0$, with $E = E_\mathrm{MF}+E_\mathrm{BMF}$. One gets:

$$\begin{align} n_\mathrm{c}\sim\frac{1}{a^3}\left(\frac{f_\mathrm{dip}(\kappa)\,\varepsilon_{\mathrm{dd}}-1}{Q_5(\varepsilon_{\mathrm{dd}})}\right)^2. \end{align} \tag{ 63 }$$

We note that equation (63) also provides a good estimate of the peak density in a trapped stabilised state, as long as one can neglect contributions of the kinetic and external trapping energies.

Interestingly, from equation (63) one can read-off that a stabilisation of the density occurs before reaching the dense regime $n_\mathrm{c}{a^3}\sim{1}$, thus forming 'ultradilute' liquid states. Two main effects enable to maintain the diluteness: the equilibrium density is reduced by the fact that there are two competing MF interactions, resulting in an effective MF attraction of much smaller amplitude than each of the two contributions (numerator in equation (63)). Second there is an amplified BMF contribution that leads to a further reduction of the density (denominator in equation (63)) [462].

The relatively low stabilising density resulting from these two ingredients protects the sample against an immediate destruction by three-body recombination. The dipolar collapse observed with Cr is thus prevented. Instead the BEC is stabilised via the effect of its quantum fluctuations, resulting in a distinct phase of liquid-like properties, the clouds formed are thus named quantum (dipolar) droplets. These droplets are denser than the MF-stable gaseous BEC, yet much more dilute than ordinary liquid (by ≈8 orders of magnitude). They offer a distinct paradigm of quantum fluid where the BMF effects are predominant yet tractable, see e.g. sections 5.2 and 5.3.4.

We note that these exact same ingredients are present in another experimental system, namely mixtures of contact-interacting BECs with repulsive intra- and attractive inter-species interactions [463–465]. The stabilisation mechanism of bosonic mixtures was in fact proposed prior to the observations on dipolar BECs by Petrov [466]. For more information on these systems, see [459].

Finally, we note that, despite the fact that we have for now neglected trapping effects in this first description, they may also play a crucial role in the newly stabilised quantum states. In the case of dipolar atoms, this effect is not only quantitative (energy shifts) but also qualitative. Indeed, the interplay between trap anisotropy and DDI yields new features already in the MF-stable BEC, as reviewed in sections 4.1.3.3 and 4.1.4.3. These new features in the excitation spectrum, that affect the MF instability, when combined with quantum stabilisation, may yield new ground states. This effects will be the focus of section 5.4.

Following the experimental observation of quantum droplets that we will more thoroughly describe in section 5.3.2, several theory works [431, 456, 467] developed the framework based on an extension of the GPE (eGPE) to include the first order BMF effects. As introduced at the beginning of this section 5, this consists in extending equation (9) with the term equation (60) (and potentially equation (58)). This perturbative approach is justified because the sample remains dilute with $na_s\ll 1$. Here the BMF correction comes into play not because the gas parameter becomes large but rather because the overall MF terms are tuned small while both contact and dipolar interaction remains individually large. Then the BEC quantum depletion remains limited.

We note that the theory model relies on the use of the LDA for the LHY term. This approximation is justified if the dominant contribution to the LHY term has momentum larger than the inverse size of the ground state. Physically, the LHY correction, which corresponds to the zero point motion of the elementary excitations (see above), are dominated by contributions of the hard modes, i.e. the most energetic ones. Crucially, when decreasing a (as to drive the MF-instability), the modes whose excitation occurs perpendicular to the dipoles get softer (and lead to the instability) while the one along the dipole get harder. The latter will then dominate the LHY correction. Remarkably, the BEC close to instability and the quantum stabilised states get very extended in the direction of the dipole, under the effect of magnetostriction (see section 4.1.2). Therefore, the LDA may remain a surprisingly good approximation in this regime, and even beyond the instability threshold.

Quantitatively, Wächtler and Santos [431] show that, using an anisotropic momentum cutoff, with an anisotropy matching the anisotropy of the droplet itself, at 1% of the inverse of the droplet extent, one still recovers most of the BMF effects (80%). This proves that the LDA is a good qualitative, and even quantitative, approximation, as long as the droplet remain elongated enough. The contribution of the long-wavelength modes may lead to small corrections, yet the authors note that this would mainly modify the prefactor of the LHY correction while the scaling should remain that of equation (60). This domination of the LHY correction by hard modes which makes relevant the use of the LDA even in small sized quantum-stabilised states is also found in quantum droplets of bosonic mixtures as originally proposed by Petrov [466], see also [459, 463–465].

The relevance of eGPE framework was later quantitatively studied both in theory and experiments. Results from the eGPE and quantum Monte-Carlo (QMC) simulations of dipoles with hard-sphere repulsion were compared in [457] and found to be in good quantitative agreement. In experiment, the degree of agreement was found to vary depending on the exact settings, in particular on the gas geometry and on the density, see e.g. [189, 246, 249, 307, 311, 312, 328, 329] and later discussions. Typical discrepancies could be simply accounted by shifting the scattering length value by a few, and up to a few tens, of percents. More recently, a quantitative study compared eGPE theory, diffusive QMC simulations using finite-range interactions, and experiments [308]. A good agreement of the QMC results and the experiments was found, whereas the eGPE results are systematically shifted. Boettcher et al elaborate on the eGPE mismatch and show that it may be accounted by a more sophisticated description of the scattering, in particular by accounting for the effects of finite collision energy on the scattering properties, due to the finite temperature of the samples. Such effects are known to result in an effective renormalisation of the dipolar length [340], see also section 3.1. By including such corrections in the eGPE yields a better agreement with the experimental data. Note that corrected eGPE and QMC constitutes two complementary theories that are able to describe experimental observations by accounting for different effects. The corrected eGPE accounts for temperature effects on the scattering properties but only includes quantum fluctuations at a perturbative and approximate level. QMC fully accounts for the quantum fluctuations, yet neglects thermal effects. The two theories account for different effects, and suggest that different sources of corrections with respect to the standard eGPE theory may be relevant. A full theoretical modelling, accounting for the different effects at once, and revealing their respective role, is yet missing.

Following the seminal work establishing a theory description of the quantum-stabilisation mechanism, theoretical works tackle the question of the ground-state phase diagram of a dipolar gas and in the presence of the newly discovered stabilising term [456, 467].

The presence of an anisotropic harmonic trap with aspect ratio $\lambda = \omega_z/\omega_r$ (see also equation (44)) was considered, yet cylindrical symmetry around the dipoles direction z was assumed. Mean trap frequencies and atom numbers following the experimentally relevant values were considered. In these works, only single-peaked ground states were predicted (i.e. no droplet assemblies, see section 5.4.1). When decreasing a at fixed atom number and trap geometry, the ground state was found to change from a low-density phase, matching a standard MF-stabilised BEC, to a high-density phase, stabilised by the LHY term, forming a 'quantum droplet'. The low-density BEC has a geometry roughly following that of the trap, while the droplet state is not. In particular, the droplet is always elongated in the direction of the dipoles, which is why it is sometimes called 'filament'. For traps elongated along the dipoles ($\lambda\gtrsim1$), the ground state smoothly evolves from the low-density MF-repulsive BEC to the high-density MF-attractive quantum-droplet phase when decreasing a. This smooth crossover can be apprehended as the trap elongation along the dipole direction enables the two states to have similar geometries (elongated along the dipoles), so that one can continuously evolve into the other, while the BMF term suppresses the dipolar collapse described in section 5.1. In the opposite case of large λ, a discontinuous transition between the two states is found, with an intermediate region of bistability, where the two states form local energy minima. This bistability can be apprehended from the following argument: in such traps, a low-density state can be stabilised even for $\varepsilon_{\mathrm{dd}}>1$ as the trap asymmetry forces the dipoles to lie side-by-side enhancing the MF DDI repulsion. This is exactly how Cr BECs have been stabilised at low a, as described in section 4.1.4, and the addition of the BMF term only weakly modifies this behaviour. On the other hand, the BMF term stabilises an other solution which is elongated along the dipoles and where MF DDI is attractive. The bistability occurs as there is an intermediate regime of a where both solutions are local energy minima, while in the absence of the BMF term (as relevant for Cr) in this region the MF repulsive BEC is metastable.

The theoretical phase diagram is represented in figure 23. The exact position of the boundaries and critical point depend on the exact experimental parameters (number of atoms, mean trap frequency, etc). The crossover and the bi-stable region are separated at a critical aspect ratio $\lambda_\mathrm{c}$, typically $\lambda_\mathrm{c}\approx 1.5$–2. In the remainder of this section, we review the extensive experimental exploration of this phase diagram.

Zoom In Zoom Out Reset image size

We note that the description given above does not give a full picture of the dipolar gas phase diagram in presence of quantum stabilisation. This description encompasses the behaviour in a low-atom-number regime. For larger atom numbers (and/or tighter traps), distinct ground states may arise, and in particular spontaneous density modulation may occur, or, in other words, ground states bearing several droplets may be found. This regime, achieved more recently in experiments, will be the topic of section 5.4.

The first evidence for the absence of dipolar collapse and the existence of a stabilisation mechanism was reported by Kadau et al [184] on ¹⁶⁴Dy, before any of the theoretical development described in sections 5.2–5.3.1. Starting from a pancake-shaped BEC of ¹⁶⁴Dy atoms in a trap with aspect ratio $\lambda\simeq3$, a was lowered down to the MF unstable regime. In this geometry, we note that a modulational instability is expected, see section 4.1.4.3. Below the instability threshold, ordered droplet ensembles were observed in the in-situ radial density distribution of the cloud, as can be seen in figure 24. Contrary to the expectations at the time, these droplets were observed to have very long lifetimes of several hundreds of ms, simply limited by atom loss due to three-body recombination.

Zoom In Zoom Out Reset image size

As shown in the early theory works [431, 456], the observation of multiple droplets arises from a crossing of the boundary between the BEC region and the bi-stable region of the phase diagram. In this regime, a continuum of metastable states, with different numbers of droplets, is expected to exist. These works suggest that, following the fast change of the interaction strength, the system ends up in one of these metastable states, explaining the observed long lifetime of the crystal as well as the shot-to-shot variability of the structures. Experimentally, the system's bistability was evidenced via the hysteresis in the appearance and disappearance of the droplet ensembles when varying a down and back up, thus supporting the theoretical picture. As observed after a long hold time, the arrangement of the droplets was found to result from the repulsion between the corresponding macroscopic dipoles (see also [429, 430]).

A following set of experiments [429] showed that the droplets have an elongated shape along the dipoles with an axial length of a few micrometers and a radial size estimated to half a ${\mu\mathrm{m}}$ or less (this lies below the imaging resolution). Furthermore, matter-wave interference fringes observed from several expanding droplets proved them to be individually superfluid. They were observed to have extremely low expansion velocities when released into a waveguide. In [429], Ferrier-Barbut et al, inspired from the work of Petrov [466] on BEC mixtures, proposed that BMF effects act as a stabilisation mechanism. By comparing additional expansion and lifetime measurements to a model following the lines given in section 5.2.2 (neglecting kinetic and trapping effects), they validate this idea and invalidate a possible mechanism based on three-body forces proposed in [450–452], see also the earlier works of [448, 449]. In particular, the lifetime τ of the sample, being set by the three-body loss processes, gives information on the sample's density via $\tau = L_3\langle n^2\rangle$ (see equation (58)). The estimated density scaling of density as a function of a deduced from measuring τ as a function of a was found to reasonably agree with the simple density scaling of equation (63), and disagree with the scaling expected from three-body forces stabilisation.

In this setup, the bistability, yielding excited metastable states of multiple droplets after an interaction quench, prevents the formation of a large single droplet in the ground state. This feature has limited the study of the properties of the liquid-like state. Shortly after, it was found out, in particular thanks to the theory development described in section 5.3.1, that this issue can be resolved by taking advantage of the 'crossover region' of the phase diagram, i.e. changing the trap geometry, as we will now discussed.

Following the first droplet observations [184, 429, 430] and the resulting theoretical development [456, 467], Chomaz et al tested the universality of the stabilisation effect (which relies on the sole quantum-mechanical nature of the fluid, provided that the interactions are strong enough), by realising the first quantum droplet of a distinct chemical species, using ¹⁶⁶Er [249]. They also used a complementary geometry to the previous Dy observations, using a cigar-shaped trap with $\lambda\ll1$, and thus they observed a smooth crossover from a BEC to a single large droplet of Er atoms containing all the condensed atoms. Because the created droplet is isolated, the authors could study its properties such as its elementary excitation, and its expansion dynamics, see figure 25. The authors performed a systematic comparison of the measurements with the eGPE predictions (equation (9) including also equations (58)–(60), see section 5.2) using independently measured values of a and L₃. The quantitative agreement reached confirmed the stabilisation scenario and validated the eGPE framework in this setting. We note again, that quantum droplets have later been observed in mixtures of contact-interacting BECs [463, 465] following Petrov's initial proposal, establishing the quantum-stabilisation mechanism as universal even across different interaction types, provided that competing interactions lead to a balance of MF and BMF contributions.

Zoom In Zoom Out Reset image size

Finally, the two regimes of bi-stability and crossover where experimentally connected in [432], where similar experiments were performed using ¹⁶⁴Dy in traps of variable aspect ratio λ. The formation of a single droplet resulting from lowering a in the crossover was observed for small λ up to a critical aspect ratio $\lambda_\mathrm{c}$ which marked the onset of a modulational instability as in [184]. The experimental value of $\lambda_\mathrm{c}$ is in agreement with the expected critical aspect ratio marking the separation between the crossover and bi-stable regions.

Following the first experimental observations, several works investigated the specific properties of the droplet states properties.

Theoretically, it is interesting to highlight that, while the trap plays a role in the phase diagram and the transition between the different phases (see section 5.3.1), it plays a much lesser role on the properties of the single droplet state itself, at least for large enough atom number and $\varepsilon_{\mathrm{dd}}$. By studying full ground-state solution from the eGPE, [430, 456, 467] show for instance a very weak dependence of the droplet density on the trapping potential. At very low atom number (i.e. on the order of a few thousand), the one-body kinetic energy is non-negligible, and its interplay with the interaction energies leads to a dependence of the central density on N. However, at high atom number, the interaction energies fully dominate. Following the simple description of section 5.2.2, the density then reaches a saturation value independent of N, given only by the balance of MF and BMF energies. This marks the low compressibility typical of a liquid phase, and which stems here from the high energy cost of increasing density due to the BMF term. While early experiments have mostly explored the low atom regime, first measurements showing the onset of this density saturation have been recently reported in [459].

Besides its evidence via density saturation, the low compressibility of the dipolar quantum droplets has been revealed in different sets of measurements, in particular focusing on the collective modes. First investigations were performed by Chomaz et al, studying a particular collective mode, with compressional character, of a quantum gas of ¹⁶⁶Er in the BEC-droplet crossover [249]. As the scattering length is lowered through the crossover, a steep increase in the mode frequency was observed, changing by $50$% when varying a by less than $20$% while the trap frequency remained constant, see figure 25. This shows that the compressibility quickly drops when going from the BEC to the droplet phase.

In addition to being weakly compressible, a dipolar quantum droplet is anisotropic. This stems of course from the anisotropy of the DDI with respect to the external magnetic field axis. The consequence of this external breaking of rotational invariance is the existence of a collective mode, deemed the scissors mode. It corresponds to a rigid-body angular oscillation around the magnetic field direction. This scissors mode has been observed first in atomic nuclei [468–470] and in contact-interacting BECs in anisotropic traps [471, 472]. In the presence of the DDI, the rotational symmetry is broken even in isotropic traps. This was showed to result in a scissors mode for MF dipolar BECs in [473]. This mode is well-defined only for angular oscillations in the xz-plane with an amplitude lower than $\theta_\mathrm{max} = (\langle z^2\rangle-\langle x^2\rangle)/(\langle z^2\rangle+\langle x^2\rangle)$ where z is the field direction [471]. It is thus very challenging to observe with dipolar MF BECs, but was observed in dipolar quantum droplets of ¹⁶⁴Dy thanks to their considerable anisotropy [307]. In conclusion, studies of the collective modes of the liquid state are a sensitive probe for its properties [474] and will likely be pushed further.

The possibility of creating a single large quantum droplet by tuning the trap aspect ratio, as described in section 5.3.3, was key in the demonstration of the self-bound liquid nature of this phase. A quantum droplet living in equilibrium between the repulsive BMF and the attractive MF interactions, the question of the necessity of the trap for the existence of this phase then arises. Two theory papers demonstrated indeed the existence of a self-bound state in the absence of a trap, characterised by a non-vanishing peak density in infinite volume [467, 475]. This defines a liquid state, existing from the mere balance between high-density repulsion and low-density attraction. However, this liquid is in essence quantum, which provides it with macroscopic new properties. First, it remains phase-coherent within a droplet. Second, the self-binding potential resulting from MF attraction from other atoms within the droplet must counter-act kinetic energy. As a consequence if the volume is too small, kinetic energy prevents the existence of a bound state. This yields a phase diagram as a function of atom number and scattering length, with a single line separating the self-bound state from a gas at low atom number and high scattering length, see figure 26. The scattering length dependence of the atom number value on the separation line is sharp. Indeed, using scaling arguments, one can show that the minimal atom number for the existence of the self-bound state scales as $N_\mathrm{c}\sim(1-\varepsilon_{\mathrm{dd}}\,f_\mathrm{dip}(\kappa))^{-2/5}$. First indications of a self-bound behaviour were already evidenced in the early work [429] through the absence of expansion in a waveguide and record-low expansion velocities in free space. Yet, the self-bound character could not be observed because of the limited atom number. Because of the large loss rate occurring in Er droplet (see figure 21), the evidences of the self-bound behaviour of the macro-droplet of [249] was partial, limited by the state lifetime. Here, only a slowing down of the expansion dynamics was observed. An unambiguous self-bound behaviour was then observed on ¹⁶⁴Dy in [246]. Here a single droplet was created using an axially elongated dBEC [476] of 3000 Dy atoms and cruising the crossover to the single droplet state, similarly to [249]. The trap was then smoothly released while the atoms were levitated against gravity using a magnetic-field gradient (see figure 27). Unexpanded droplets were observed up to 90 ms of levitation. Using the atom loss as a probe, the liquid-to-gas transition was mapped out in the (a, N)-space. This phase-diagram was later expanded to larger atom number by using samples of ¹⁶²Dy atoms [308], see figure 26. The experimental points of the two isotopes match well together and agree with eGPE predictions, at the expenses of a substantial renormalisation of the background scattering length of the ¹⁶⁴Dy isotope, by a factor of 3/4. The need of such an important change of $a_{\mathrm{bg}}$ rose questions about the validity of the eGPE treatment in particular in the low-atom-number and large-$\varepsilon_{\mathrm{dd}}$ regimes (see also the discussion in section 5.2), or again of the methods employed to extract $a_{\mathrm{bg}}$ (see section 2.4), in particular in the Dy case.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To conclude, the self-bound liquid-gas phase diagram, being simply due to the quantum mechanical nature of the system, is not restricted to the dipolar case. It may in fact be generally expected in quantum gases where two interactions of different origin compete, and as such it has been observed in Bose-Bose mixtures [463, 465]. It is interesting to note that in lower dimensions this phase diagram is completely modified. This has been the topics of several theoretical works, both for dipolar and mixtures systems, see e.g. [477, 478].

Supersolidity is a paradoxical phase of matter in which the antithetical properties of crystal arrangement and of superfluid flow coexist. It has been suggested more than half a century ago as a paradigmatic manifestation of a state in which two continuous symmetries of distinct nature are simultaneously broken [479–481]. Originally predicted in quantum solids with mobile bosonic vacancies [482–484], the search for supersolidity has spread in many different fields of physics. Observation of supersolidity in helium was claimed [485] but the claim was withdrawn by the same authors a few years later [486] and the quest for supersolidity in helium is still open [481, 487]. The possibility of supersolid states in quantum gases were theoretically proposed long ago, see [488–498]; this possibility linking back to the seminal work from E. Gross [480] on assemblies of bosons with momentum-dependent interactions, see also [397, 499].

Besides the supersolids made from dipolar quantum fluids alternative approaches to supersolidity with ultracold atoms include miscible two-component BECs with SOC as well as optically pumped superfluids in a cavity that mediated interactions via the scattered light [500, 501]. All concepts have in common a momentum-dependent interaction that leads to a minimum at finite momentum in the dispersion relation on the superfluid side of the phase transition. This is called the roton minimum. In both the dipolar supersolid and the SOC systems, two branches of the excitation spectrum exist on the supersolid side of the phase transition. They correspond to the Goldstone modes of the two double symmetry-breaking processes, the breaking of the translational symmetry and the breaking of the U(1) phase invariance of the condensate. More recently, a continuous U(1) translational symmetry was created in system of multimode cavity light coupled to a BEC. The breaking of this symmetry resulted in a supersolid with transverse vibrations exhibiting a Goldstone dispersion of phonon excitations [502]. We will now focus the discussion on the simultaneous quest for a supersolid supporting crystal phonon excitations like a real solid in a dipolar gas only.

This authentic supersolid behaviour can be expected if the simultaneous breaking of the two symmetries arise from the intrinsic interactions between the particles. As already highlighted by E. Gross, the relevant platforms for observing supersolidity are quantum gases with inter-particle interactions yielding large-momentum attraction [397, 480]. Practical examples are Rydberg-dressed potentials, spontaneous or light-induced DDIs in confined geometries [489–498]. In these settings, a roton-type excitation is induced in the superfluid' excitation spectrum by the large-k attraction (see also section 4.1.3.3). Its full softening may indicate the transition to a density-modulated state as the roton signals an intrinsically favoured length scale and has been seen as a precursor of crystallisation [397, 499]. The observation of a roton mode in dBECs confined in cigar-shaped geometries [311, 312] was key to point promising grounds to observe supersolidity [498], see also section 4.1.3.3.

A recurrent hindrance to the interaction-driven formation of supersolid in quantum gases lies in the predicted MF collapse of the gas at the roton instability, see sections 4.1.4.3, 5.1 and [503, 504]. Protocols to stabilise the gases beyond its instability thanks to the engineering of higher-order (three-body) interaction potentials have been proposed [497]. However, these protocols have not been implemented in experiments so far. The discovered BMF stabilisation mechanism provides the key ingredient for an intrinsic stabilisation.

In 2017, theoretical works based on the eGPE [505, 506] demonstrated that not only a single droplet (as described in section 5.3.1) but also assemblies of multiple droplets, or in other words density-modulated states, could constitute the ground state of dipolar quantum gases. This was found for both pancake and cigar geometries, by simply using tighter trap and/or larger atoms numbers than in the previous works of [456, 467, 474]. Similar observations were also achieved via QMC simulations with large densities [507, 508]. Starting from the picture of a single droplet state, the physical argument to the formation of multiple-droplet ground-states can be formulated as follow: the liquid phase has a weak compressibility that yields a sharp increase in energy at high density, see also section 5.3.4. If one compresses an (isotropic) liquid in one or two directions, it then deforms to keep a constant density at the small cost of increasing its surface energy. Yet for the dipolar quantum liquid, the DDI dictates an anisotropy, and deformation costs a high amount of (dipolar) energy. Therefore, when a dipolar quantum droplet is compressed along its long axis, it might be energetically favourable to split it in several droplets (of smaller radial sizes), recovering thus an anisotropy closer to the one naturally imposed by the DDI. Ground states can thus form spontaneous density modulation. Calculations with periodic boundary conditions along one direction were also performed, see [505], proving that a continuous translation symmetry can indeed be broken (at the thermodynamic limit). These predictions indicate possibilities for quantum-stabilised dipolar 'crystallised' ground states and in particular open the doors for quantum-stabilised dipolar supersolids.

The work of [505] also experimentally produced spontaneously density-modulated states in cigar-shaped Dy gases that show close features to the expected ground states yet global phase coherence was found to be absent. The lack of phase coherence was interpreted based on a Josephson-Junction formalism, which describes the maintenance of a phase relation via tunnelling processes between the individual droplets, and its preclusion via quantum or thermal fluctuations. In this experimental realisation, the tunnelling rate between the individual droplets, or in other words, the wave-function overlap between them, was not enough to lock the droplets in phase. The authors of [505, 506] also showed that assemblies of droplets with sizeable wave-function overlap were theoretically possible by either increasing the trap frequencies or the atom numbers compared to the current experimental configurations.

At the end of 2018, the progressive understanding of the many key features of Bose gases of highly magnetic atoms combined in an acute picture. These features include the discovery of the BMF stabilisation in such gases (see section 5.3.2), the observation of the roton mode and its softening in cigar-shaped clouds (see section 4.1.3.3–4.1.4.3), the possibility of density-modulated ground states, matching an assembly of quantum-stabilised droplets (see section 5.4.1), and the need of making the droplets more extensively overlap in order to maintain the global phase coherence in such assemblies (see section 5.4.1). Building on this knowledge, a set of experiments [189, 328, 329] observed hallmarks of supersolid behaviours in cigar-shaped gases of ¹⁶²Dy, ¹⁶⁶Er, as well as ¹⁶⁴Dy atoms. Accompanying theoretical works [328, 329, 498, 510, 511] related these observations to an underlying supersolid ground state of the gas, see section 5.4.3.

The three experimental works proved two hallmarks of supersolidity in experiments, namely the simultaneous occurrence of a spontaneously formed density modulation and of a global phase coherence. Density modulation was revealed by the occurrence of modulated patterns in the gas absorption images, either in situ [329] or in time of flight. Here the density patterns form from the self-interference of the gas via free expansion of the in situ modulation. The global coherence was demonstrated by a statistical analysis of multiple repetitions of this self-interference TOF patterns. Global coherence is here marked by the stability of the interference patterns, differentiating supersolids from an incoherent array of droplets. Global coherence was mostly quantified through an analysis of the complex values of the Fourier transform of the TOF density profiles, see figure 28. We note that the coexistence of density modulation and phase coherence does not prove the superfluidity of the state, which can be demonstrated only by probing dynamic-related properties of the state, see section 5.4.4.

Zoom In Zoom Out Reset image size

All three experiments used relatively shallow cigar-shaped traps with transverse magnetisation, similar to [311]. The spontaneous density modulation occurred only along one axis, the long axis of trap. The coexistence of density modulation and phase coherence was observed in narrow ranges of scattering length values, of a few a₀ wide, and survives a few tens [189, 328, 329] up to a few hundreds of ms [328] in these first experiments, later extended up to few seconds [512]. Building on this long lifetime, [328] additionally established a different route than standard interaction tuning towards supersolid states, which is based on direct evaporative cooling starting from a thermal state, see also [512].

In all three works, the states with supersolid properties could be achieved by ramping down the scattering length starting from a stable BEC, using a slower ramp and a finer tuning of a than [311]. In the earlier works of [430, 505] observing droplet assemblies in cigar shaped traps, the rougher tuning of a as well as the lower initial BEC atom numbers are also thought to have prevented the observation of supersolid properties (see also sections 5.3.2 and 5.4.1). Finally, in the seminal work of [184], a distinct pancake-shaped geometry as well as smaller initial BEC atom numbers were used. Several following theoretical as well as experimental works indicates that the use of a cigar-shaped geometry was crucial for more easily achieving supersolid states in experiment. In this geometry, the MF instability is driven by the softening of a single (doubly-degenerate) roton mode dictating the dominant wavelength of the density fluctuations [311, 312, 411]. In contrast, the pancake case is more complex, with several radial and angular rotons, corresponding to different structures of density fluctuations, simultaneously softening [412, 513]. Furthermore, in the cigar-shaped case, it is expected that the transition to a supersolid state can occur continuously, with the supersolid modulation directly connecting to the softened roton mode, see e.g. [514].

Following the observations of the roton mode population in [311] and shortly preceding the works of [189, 328, 329], Rocuzzo and Ancilotto [498] theoretically explored the phase diagram of an Er quantum gas in an infinite cigar-shaped geometry with periodic boundary conditions. They relied on the eGPE framework developed in the context of the quantum droplet studies (see section 5.2) and calculated the ground state as a function of a. In this setting, they demonstrated the existence of a SSP in an intermediate range of a, separating an array of insulating droplet (ID) at low a and a regular BEC at large a. Both the SSP and the ID are density modulated ground states stabilised by quantum fluctuations, similar to the single-droplet phase described in section 5.3.1. The SSP distinguishes itself by bearing a density modulation of finite contrast. In addition, Rocuzzo and Ancilotto directly computed, via dynamical simulation, the superfluid density of the density-modulated states, thus rigorously establishing the connection between the occurrence of non-fully-contrasted density modulation and non-zero superfluid fraction. This work proved the relevance of the SSP in cigar geometries, beyond finite-size effects, see also [514].

Together with their experimental results, both [328, 329] reported on eGPE calculations of the ground-state phase diagrams in the finite experimental geometries, see figures 29(a) and (f). The ground states were calculated as a function of the atom number N and the scattering length a for a given trap and atomic species. At low N, a transition occurs from the regular BEC to the single-droplet phase when decreasing a, as described in section 5.3.1, see also [328, 456, 467]. In contrast, for large-enough atom numbers, three different regimes could be identified similar to the infinite case results of [498]: when decreasing a, the regular BEC state transitions first to a density-modulated state of finite contrast, identifying a SSP, and then to a density-modulated state of contrast almost unity, forming an ID array, see figures 29(b)–(e) and (g)–(i). We note that the atom number for the occurrence of density-modulated states versus single-droplet depends on the atomic species (being smaller for Dy than for Er), and on the trap, both on its overall tightness and on its shape. The exploration of the most favourable parameters for supersolidity as well as achieving such a phase in different settings are interesting directions that are currently under investigation, see e.g. section 5.4.5.

Zoom In Zoom Out Reset image size

Based on the eGPE framework, [189, 329, 515] also reported on simulations of the real-time evolution induced by the finite scattering-length ramp used to reach the various BEC/SSP/ID states in the experiments. These show that in the SSP regime, the dynamical state has strong similarities with the expected ground state and shows limited phase fluctuations. In contrast, when decreasing the scattering length lower, to an ID regime, the dynamical state deviates from the ground-state expectation and show large phase fluctuations.

Besides the study of the static properties of the state, as the hallmarks reported in the first studies (see section 5.4.2), a huge interest is drawn by the special dynamical properties associated to supersolidity. Indeed, such a dynamics intrinsically connects to the superfluid character of the supersolid and would characterise the related rigidity of its phase. The study of the dynamical properties covers a wide range of phenomena and concepts, spanning from the study of the spectrum of elementary excitations, to that of transport properties, or of the response to an external rotation for instance. Supersolidity brings exotic characteristics to these different features.

The first of these aspects to be explored in dipolar gas experiments concerns the properties of the spectrum of elementary excitations of a supersolid at low momenta. Here, not one but two branches of low-energy (gapless) 'phonons' modes [482, 498, 516–519] are expected, the different branches bearing excitations of different characters. Following Goldstone's idea [520], the number of gapless branch connects to the number of spontaneously broken symmetries—two in the case of the supersolid: one branch holds phonons of the crystal and relate to a broken translational symmetry, while the other branch holds phonon of the superfluid (phase phonons) and relate to a broken gauge symmetry. The phase branch corresponds to the branch of lower velocity, and, when evolving through the BEC-SSP-ID phase diagram, it is observed to soften and its weight in the response to decrease, up to vanishing in the ID case. On the contrary the crystal branch slightly harden. These properties of the supersolid excitation spectrum were theoretically investigated in various systems [516–519]. The excitation spectrum behaviour in dipolar supersolids was further theoretically investigated in various works [310, 498, 515, 521].

Several experiments show signatures of these two branches [190, 310, 411, 522]. Due to the finite system size, the phonon branches are here discretised (see also section 4.1.3), and the experiments probe the response of some specific low-lying modes of the trapped dipolar supersolid states. In [190, 310], the evolution of the BEC's lowest-lying quadrupole mode within the BEC-ID-SSP phase diagram is probed, see figure 30. In the BEC regime, the system responds at a single and roughly constant frequency, and when reaching the SSP regime, several frequencies are observed in the system's response. Furthermore, these frequencies are found to organise in two branches as a function of a. One branch is softening and one is hardening when decreasing a. These properties are indicative of the emergence of modes of dominant phase and crystal characters, respectively. Tanzi et al [190] identifies the different modes via their distinct signatures in the time-evolution of the gas's self-interference patterns. Natale et al [310] identifies resonant frequencies via an unbiased principal component analysis of the self-interference patterns and observes a mixed character in the associated principal component structures. The different modes' characters between [190, 310] may be attributed to the different numbers of density peaks in the underlying supersolid states. Distinctly, Guo et al [522] probe the lowest lying mode of the trapped supersolid, which has a frequency lying below the dipole mode and is a Goldstone mode of phase character, see figure 31. The frequency of the mode itself is not probed due to its low value, instead signatures of a spontaneous population of the low-energy excited mode are observed via a statistical analysis of the in situ density patterns. A peculiar feature of this mode is that the motion it induces preserves the centre of mass position thanks to superfluid flow. This implies a correlation between the array's displacement and the population imbalance between the density peaks, see figure 31. By repeatedly producing and imaging steady-state samples, Guo et al [522] evidenced such correlations in the limited a range where the state is supersolid. The presence of such correlations points to the existence of the phase Goldstone mode and implies an underlying superfluid flow. An interesting point to note is that the emergence of this Goldstone mode in the excitation spectrum connects to the softening of the antisymmetric roton mode from the BEC [522], see also section 4.1.3.3. In contrast, the symmetric roton mode connects to an amplitude mode, also called Higgs mode, which sharply hardens when moving away from the SSP-BEC transition in the SSP [521]. More recently, a detailed study of the in situ density fluctuations in the SSP also revealed the existence of density and crystal phonons [411], see also section 4.1.3.3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Besides the study of low-lying excitations, investigation of higher energy modes of supersolids has been recently undertaken. This includes studies of the response to rotational excitations [523], or to large-momentum excitations [515]. Rotation excitations in [523] probe a scissor-type excitation, similar to that described in section 5.3.4 for droplet states. Here the observed scissor-mode frequency decreases at the transition from BEC to modulated states. This indicates an increase of the state's momenta of inertia. It is however not straightforward to relate this behaviour to the superfluid fraction of the state. [515] probes the scattering of quasi-free particles via Bragg excitations (see also section 4.1.3) along the BEC-SSP-ID phase diagram. Here a decrease in the response amplitude is observed when crossing from BEC to modulated states. From impulse-approximation theory, this decrease probes a reduction of the population of the zero-momentum state. The observed decrease surpasses expectations obtained from the eGPE steady states. This was related to coherent phase evolution, induced by the dynamical crossing of the quantum phase transition [524].

Besides studies of specific excited modes, the supersolid's dynamical properties have also been unveiled in the time evolution induced by parameter quenches or ramps. This includes interaction quenches [509] and evaporative cooling ramps [512]. Ilzhöfer et al [509] probes the dephasing and rephasing dynamics of the gas' global coherence when dynamically crossing the ID-SSP transition. The observed dynamic was understood by comparing to a simple Josephson-Junction-Array model where the Josephson tunnelling amplitude was quenched, while dissipation is empirically introduced, using a Langevin formalism. This model shows a good qualitative agreement with experiments, while quantitative discrepancies are attributed to the inherently soft nature of the supersolid's crystal and the effects brought in by the possible excitation of the crystal's phonons. Studying the formation dynamics of supersolids via direct evaporative cooling (ramp of ODTs), Sohmen et al [512] shows that density modulation appears before global phase coherence, but after local phase coherence, is established. Sohmen et al [512] also studies the subsequent decay process of the supersolids, occurring spontaneously under the effect of three-body recombination. In this decay, in contrast to the formation dynamics, global coherence is observed to survive longer than density modulation, while the temperature remains roughly constant. Furthermore, in this process, the strength of the density modulation, for a fixed value of the number of locally coherent atoms, is observed to depend on the temperature, being weaker when colder, see also [525].

The various dynamical studies performed up to now provide valuable insights into the supersolid state special behaviours, including direct and indirect signatures of the existence of phase and crystal phonons (near k = 0) and their mixing, of superfluid flow and transport, or of coherent and incoherent phase dynamics in the experimentally produced states. Interesting prospects include direct measurements of the superfluid properties, and in particular of the superfluid density of supersolids, and the creation of vortices in the supersolid rotation, see e.g. [526, 527].

As discussed in section 5.4.2, the first observations of supersolidity in dipolar gases have been conducted in elongated and relatively shallow, cigar-shaped traps, resulting in supersolids in which (a) the breaking of the translational invariance occurs only along one dimension, (b) the number of density modulations found in the finite experimental system is small (typically a handful) and thus finite size effects might be important. Yet, dipolar supersolid states with larger and/or more complex crystalline structures, and in particular where crystallisation occurs in two directions of space have attracted intense interest. Theoretically, a series of works have focused on the phase-diagram and excitation spectra of two-dimensional supersolids in isotropic and anisotropic traps [528–532], on the possibility of creating vortex excitations [526, 527, 533], and on the emergence of exotic crystalline structures [513, 528, 529]. Very recently, two-dimensional supersolidity has been observed in experiments with Dy atoms, using anisotropic traps [530]. As a function of the anisotropy of the trap in the directions perpendicular to the atomic dipoles, experiments have demonstrated evaporative phase transitions to ground states of various supersolid patterns from a linear chain, to a zig-zag crystalline structures [530], and finally to a hexagonal structures in circular traps [531]. In a related effort, two-dimensional angular roton modes, analogous to the linear roton mode in elongated traps [98, 311, 312], has been observed in radially extended traps [412], see also section 4.1.3.3 for more details. By further tuning the interplay between trapping configurations, interaction parameters, and atom numbers, even more exotic ground states are expected to appear, such as honeycombs, ring and labyrinth-like SSPs [513, 528, 529]. With the considerable interest that the discovery of dipolar supersolids have generated, we expect a quick development of the field and a blossom of works probing the various intriguing properties of this phase as well as exploring its possibility in different settings.

Unlike electrons, atoms possess a composite spin that may be large compared to electrons' $S = 1/2$. This enriches spinor physics. For alkali atoms, the total spin arises from the coupling of the spin-1/2 electron to its orbital angular momentum L and the nuclear spin I. This large spin allows for the study of a new type of quantum fluid, involving the interplay between magnetism and superfluidity [94]. The study of degenerate quantum gases with a large spin degree of freedom $s>1/2$ initially focused on bosonic atoms, i.e. Bose spinor gases [93]. Experiments explored the spinor physics of F = 1 Na atoms and F = 1 and 2 Rb atoms [94]. More recently, these studies were extended to spinor fermions using the $F = 9/2$ state of fermionic K [534].

Strongly magnetic atoms such as Cr and Lns are multi-electron systems, which can possess an even larger total spin in the ground state than alkalis. For example, bosonic Cr atoms have an electronic spin S = 3 in the ground state (L = 0), while bosonic Dy and Er atoms have $L = 6,\,S = 2$ and $L = 5,\,S = 2$ respectively. Including the hyperfine structure for the Fermi isotopes, one can reach, for example, $F = 21/2$ using ¹⁶¹Dy atoms. Besides the richer physics at the many-body level, such large-spin atoms have also attracted attention for studying highly non-classical behaviours involving quantum coherence and entanglement, either at the many, few, or even single-particle levels [180, 181, 535]. These prospects may have important implications and applications for quantum information processing, sensing, or metrology purposes. In particular, we note that a single magnetic atom provides the realisation of large spin states (length F) in an Hilbert space of relatively moderate size ($2F+1$) and increasing linearly with the spin size. In contrast, a spin of similar length realised by a set of spin-$1/2$ particles lives in an Hilbert space whose size increases exponentially with the spin length. The reduction of the Hilbert space offered by magnetic atoms is beneficial to temper the effect of decoherence in the system.

We note that, in atomic gases, research has also explored small-spin systems, e.g. those involving only two spin states to form an (effective) spin-1/2 system. Such a system can be accessed either by considering a hyperfine level of $F = 1/2$ (as for fermionic Li and K [377, 536–539]) or by isolating two states of a larger-F level; see e.g. [540, 541]. For magnetic atoms, however, isolation is complicated by the spin-changing DDI; see also section 3.3. Nevertheless, as we have seen in section 3.3.4, thanks to the quantum statistics, some subspaces can be effectively preserved in spin mixtures of ultracold fermions on long time scales because decay rates are strongly suppressed [164]. This is in particular the case for the subspace formed by the two lowest spin states of a fermionic magnetic atom. This feature has been used in experiments to study effective spin-1/2 dipolar systems under the effect of an artificial SOC [251] and in bulk [313].

In order to study pristine spinor physics, i.e. to isolate the effects arising solely from spin-dependent interactions, it is convenient to first consider a spin-independent one-body Hamiltonian (typically the sum of a kinetic term and external potential; see equation (9)). Using a spin-space generated from the hyperfine level of an atomic species, the kinetic part of the Hamiltonian is spin independent. Then one need only consider the conservative trap confining the atoms. This potential should be independent of the magnetic sublevels. For alkali atoms, this is typically obtained by using ODTs created by far-off-resonant focused laser beams, in which case the AC-Stark shift is approximately independent of the magnetic sublevel [245]. As discussed in section 2.3, the situation of magnetic atoms differs. Both for Cr and for Lns, the large electronic spin/angular momentum results in large vector and tensor parts of the atomic polarisability. This yields a significant dependence of the light shifts on the magnetic sublevel, even when using light far away of optical transitions. Typically, this is not negligible in experiments exploring dipolar magnetism. It is then expected that spinor phases will be driven by an interplay between spin-dependent interactions, the Zeeman effect, and the tensor light-shift.

To study of spin physics with large-spin atoms it is important to recall that the interaction between the atoms can depend on the spin channel $\mathcal{S}$ of the collision. Here, $\mathcal{S}$ describes the total spin of two atoms and is preserved during the collision in the absence of a DDI, see sections 2.4 and 3.

Quantum statistics plays a major role in selecting which molecular potentials should be taken into account. Generally speaking, regardless of the statistics of the particles involved, $\mathcal{S}+l$ must be even, where l is the relative angular momentum; see also section 3.1. For short-range interactions, $s-$wave scattering dominates at low collision energy. Under these circumstances, only even $\mathcal{S}$ need to be considered, both for bosons and fermions.

The different molecular potentials have identical multipole expansion parameters at long distance where only electrostatic interactions contribute. However, at short distance, where electronic orbitals may overlap, quantum statistics also plays a pivotal role. This is because the molecular potentials associated with different spin channels are often quite different. Therefore, different spin channels $\mathcal{S}$ typically correspond to different scattering lengths $a_{\mathcal{S}}$. This is clearly true in the case of Cr because of the large electronic spin of the atom. This leads to seven very different molecular potentials, as calculated in [542] and confirmed by experiments [161, 195]. Scattering lengths are also (most likely) spin dependent in Lns. Because the outermost electronic shell has zero spin (identical to the case in Yb), all molecular potentials have a similar shape even at short distances. However, slight differences between the molecular potentials (due to, e.g. orbital anisotropy of unfilled submerged f-shell) are most likely sufficient to lead to spin-dependent contact interactions [206]. As highlighted in section 2.4, scattering properties of Ln are extremely difficult to predict, while, on the experimental side, spin-dependent scattering lengths have remained unexplored in these atoms; see also section 2.4.

In gases of large-spin magnetic atoms, there is, in addition to spin-dependent contact interactions, the DDI that must be taken into account. Dipolar interactions are not only spin-dependent, but they also induce direct SOC; see equations (3)–(5). Indeed, when including the spin degree-of-freedom, the DDI leads to two extra terms in the Hamiltonian. These are in addition to the bare (yet also spin-dependent) elastic term equation (3) that we have considered up till now. The first term is the so-called spin-exchange term of equation (4): it changes the spins of atoms at long distances while conserving the total longitudinal magnetisation of the pair of particles. The second term is the relaxation term of equation (5). This leads to a change in magnetisation and converts spin angular momentum into orbital angular momentum to conserve the total angular momentum in the system, thus providing a form of SOC. Due to this peculiar SOC, and to the long-range coupling between spins, large-spin systems thus allow the exploration of magnetism beyond paradigms inherited from solid-state physics in addition to beyond what could be achieved with ultracold atomic gases up till now. Experiments can be performed both with bosonic and fermionic isotopes (see section 2), greatly enhancing the scope for new magnetic behaviour.

Because of the spin-dependence of the contact interactions described above, the Hamiltonian of a spinor gas depends on the spin s of the interacting particles even in the absence of a DDI. The Hamiltonian becomes increasingly complicated for increasing s [93]. It is useful to examine the simplest case for what can be considered as a large-spin system, i.e. one beyond the spin-1/2 case. This is the case of s = 1 atoms. In this case, only two spin channels exist for the short-range interactions, $\mathcal{S} = 0$ and $\mathcal{S} = 2$, with the corresponding scattering lengths $a_{\mathcal{S} = 0}$ and $a_{\mathcal{S} = 2}$. In the Hamiltonian, given by equation (13) within the MF picture, the difference between the two scattering lengths is the key ingredient for spinor physics, giving rise to a so-called c₁-interaction term with $c_1 \propto \left(a_2 - a_0 \right)$.

This spin-dependent interaction term for example results in spin-exchange dynamics that is experimentally apparent by monitoring the evolution of the population in the different Zeeman states as a function of time. The spin-exchange processes correspond to the Forster-like exchange of a spin excitation $(m_s = 0,m_s = 0) \rightarrow (m_s = -1,m_s = 1)$ and are driven in the MF picture at a rate:

$$\begin{align} \Gamma_\mathrm{exc} \propto \frac{4 \pi \hbar^2}{m} n (a_{\mathcal{S}=2}-a_{\mathcal{S}=0}). \end{align} \tag{ 64 }$$

Spin dynamics has been widely studied in the community for both Bose and Fermi quantum degenerate and thermal gases [94, 534, 543–550]. More generally, the study of excitations within a spinor condensate is expected to be very rich, with a number of possible topological excitations observable, such as line defects, point defects, skyrmions [551], and knots [93].

Most importantly, spin-exchange processes associated with contact interactions conserve the total longitudinal magnetisation of the pair of particles. This conservation stems from the isotropy of contact interactions. A very important practical consequence is that the linear Zeeman effect is gauged-out and does not contribute to the dynamics or to the phase diagram. Spinor phases are rather determined by the interplay between spin-dependent contact interactions and the quadratic Zeeman effect [94], leading to quantum phase transitions [545–547]. At low magnetic fields, where the quadratic Zeeman effect may be neglected, the spinor ground state strongly depends on interactions. For s = 1, when $a_2\lt a_0$, the BEC is ferromagnetic (favouring collisions in the S = 2 channel), while the predicted ground state when $a_0\lt a_2$ is polar [552, 553]. On the other hand, a large positive quadratic Zeeman effect favours the polar state. As the spin increases, the wealth of possible phases makes spinor quantum gases a promising area of research [93], including for the study of non-Abelian spinors [156].

One of the important features of spinor gases is that quantum correlations can naturally occur, allowing quantum fluctuations to play an important role. For example, spin-exchange dynamics vanish at the MF level when a spinor condensate is prepared in the initial state $m_s = 0$. Quantum fluctuations then drive the onset of spin dynamics, with analogies to parametric amplification [550, 554, 555]. Atomic quantum optical effects lead to the generation of entanglement [556], which may have important applications for atom-based squeezing and atom interferometry [557].

Quantum fluctuations can also play a key role in determining the nature of the many-body ground state. A formal mapping to quantum optics was proposed by Law et al [558], who showed that the ground state differs from MF predictions [552, 553], at least for mesoscopic systems. Depending on the sign of the spin-dependent interactions, the ground state is strongly degenerate (with spontaneous symmetry breaking likely), or a condensate of spin-singlet pairs. So-called 'fragmented BECs' have also been studied by Ho and Yip [559], and were very recently experimentally produced in [560].

While most of the experiments have been performed with bosonic atoms, large-spin fermions are also a promising direction for research on large-spin systems. First studies have been performed in Hamburg [534, 549]. These large-spin Fermi systems possess increased spin fluctuations due to the large spin [561], new SU(N) symmetries for purely nuclear spins, as it is the case for Yb and Sr [562–566], BCS pairing with a non-singlet character [567], and greater than two-particle clustering [568].

The impact of dipolar interactions on spinor gases crucially depends on whether the dipolar field is larger or smaller than the external magnetic field. For a large magnetic field, magnetisation-changing collisions merely lead to atomic losses, and therefore, within the dipole-relaxation-limited lifetime, the experiment is essentially performed at constant magnetisation. That is, it evolves under terms that conserve the longitudinal magnetisation in the dipolar Hamiltonian.

Because they are spin-dependent, the DDIs modify the properties of spinor gases when they are not negligible compared to spin-dependent contact interactions. To judge a priori the impact of the DDI, one needs to compare the dipolar length ${a_\mathrm{dd}}$ to differences between scattering lengths [16]. As a consequence, dipolar effects can be prominent even for the case of alkali metal atoms, such as Na or Rb for which spin-dependent contact interactions are much weaker than spin-independent interactions. For example, the long-range and anisotropic nature of the DDI introduces a nonlocal nonlinearity that can modify spin textures, as experimentally investigated using alkali atoms [569, 570]. As we will see, dipolar interactions also impact spin dynamics.

The modifications introduced by the DDI on spinor physics are even more pronounced when the magnetic field is not large compared to the dipolar field. Then, the magnetisation-changing terms in the dipolar Hamiltonian cannot be neglected. Unlike contact interactions, the DDI, by allowing changes in magnetisation, yields an intrinsic SOC which breaks rotational symmetry; see section 6.1.3. It thus modifies the classification of the possible spinor phases. The exact nature of the dipolar spinor ground state phases at low magnetic field, which may depend on the trapping geometry and other spin-dependent forces, is still a matter of theoretical investigation; see for example [571] for the s = 1 dipolar system in the single mode approximation and [156, 157].

Spin–orbit coupling associated with the DDI corresponds to the exchange between orbital angular momentum and spin angular momentum. The DDI can, for example, trigger the demagnetisation of a sample, leading to spontaneous rotation, analogous to the Einstein–de-Haas effect [349, 351–353, 572, 573]. Furthermore, because of magnetisation-changing processes, the dipolar spinor gas is sensitive to magnetic field, which provides new ways to study phase transitions driven by the interplay between spin-dependent interactions and the linear Zeeman effect [350]. For example, it is expected that a ferromagnetic spinor dipolar condensate may display spontaneous circulation in the ground state at low magnetic fields [574].

In the following, we will describe the experiments that address the role the DDI has played in spinor physics, distinguishing between the cases where magnetisation is conserved from those where it is free to evolve.

One of the prominent features of the DDI, associated with its non-trivial momentum dependence shown in equation (8), is a tendency to develop spatial structures. This has already been discussed in the spin-polarised case in sections 4.1.3.3, 4.1.4.3 and 5. In the context of spinor dipolar gases, this tendency translates into the possibility of developing various spin-textures [93], i.e. inhomogeneous distributions of the spin vector.

The influence of the DDI on spin domains of Bose quantum gases was experimentally studied both in Berkeley [569] (for a Rb F = 1 BEC with ferromagnetic contact interactions) and in Tokyo [570] (for a Rb F = 2 BEC whose contact interactions disfavour ferromagnetism). In both experiments, spin domains were created using a well-defined magnetic field gradient. In [569], the spin textures decayed towards a spatially modulated structure of spin domains. The crucial role of the DDI was demonstrated by eliminating the DDI by a sequence of rf pulses, in which case the authors observed a suppression of the formation of the short-range domains. In [570], the spin textures observed after propagation within the magnetic field gradient are compared to numerical simulations of the GPE. This reveals that the observed spatial modulation of the longitudinal magnetisation is due to the spin precession in the effective magnetic field produced by the DDI. Both these results demonstrate that the DDI has considerable effect, even on spinor condensates of weakly dipolar alkali metal atoms.

Phonon excitations have been studied in the context of scalar BECs through, e.g. Bragg excitation spectra [384–388, 575]; see also section 4.1.3 for the dipolar case. Moreover, there exists, in spinor BECs, spinfull excitations that also behave like quasiparticles. Surprisingly, these excitations remain relatively unexplored experimentally. One of the first investigations was performed at Berkeley on a ferromagnetic Rb BEC [576]. Magnon excitations were precisely studied and shown to consist of a standing-wave of spin-excited atoms above a ferromagnetic BEC. While the Goldstone theorem predicts that a ferromagnetic state has gapless excitations (for symmetric interactions), a gap was measured and ascribed to the presence of a weak DDI that breaks rotational symmetry. Perhaps related to this feature, the authors also observed a larger effective mass than expected by MF and BMF theories involving only s-wave interactions. Although these effects remain small due to the weakness of the DDI in alkali atoms, they provide qualitative departures from the paradigms of magnetism in the presence of symmetric short-range interactions. This constitutes but one illustration of the possibilities dipolar spinor gases offer.

In the case of Cr, trapped magnon excitations were created by applying, on a polarised BEC, magnetic field gradients perpendicular to the magnetic field axis [260]. In that experiment, the wavelength of the excitation is not small compared to the cloud's size. The magnons that are created are quantised in energy, and manifest themselves as collective modes that couple the spin and the orbital degrees of freedom. The lowest energy mode, whose frequency is set by the zero-point energy of the BEC, consists of a sinusoidal oscillation of the local spin around its original axis, with an oscillation amplitude that linearly depends on the spatial coordinates. The observations are in excellent agreement with hydrodynamic equations. In the regime the experiment was performed, the observed spin mode has a universal character, independent of the atomic spin and spin-dependent contact interactions, and is therefore rather insensitive to the DDI. It was nevertheless predicted that dipolar interactions could alter such a mode, provided that the atomic sample has a size larger than a natural wavelength set by the DDI [260].

Strong dipolar effects have been observed in the spin dynamics of magnetic atoms even on time scales short compared to the dipolar relaxation processes. Studies where performed on S = 3 ⁵²Cr atoms, initially polarised in the lowest energy spin state $m_s = -3$, and after homogeneously quenching a BEC into a spin-excited state. This was first performed using an engineered tensor light shift to promote all atoms into a well-defined Zeeman state $m_s = -2$ [355]. In a latter experiment, the spin excitation was performed by simply rotating the spins using a rf pulse [257].

In both cases, out-of-equilibrium spin dynamics was monitored by means of a Stern–Gerlach procedure to measure the population of atoms in the different Zeeman states as a function of time. It was found that spin dynamics was driven by an interplay between spin-dependent contact interactions and DDI. Spin-dependent contact interaction played the dominant role in the spin dynamics because they are relatively strong in the case of ⁵²Cr atoms.

In the case where the atomic spins are tilted compared to the magnetic field by an rf pulse, it was nevertheless demonstrated that initial spin dynamics were entirely triggered by the DDI. This is a consequence of the SU(2)-symmetric nature of contact interactions, which cannot trigger dynamics after a simple rotation of the atomic spin initially in a stretched (ferromagnetic) state. Therefore, the onset of spin dynamics after rotation of the spins seen in [257] is a purely dipolar effect. Finally, in the specific case where the spins were tilted by $\pi/2$ compared to the magnetic field, it was found that the spin dynamics vanishes. The MF theory provides a natural explanation for this phenomenon, as the torque associated with the inhomogeneous magnetic field carried by the atoms vanishes when the magnetisation of each atom vanishes [577].

In the $\pi/2$ rotation case, spin-dynamics could be recovered by applying a magnetic field gradient providing the necessary SOC to trigger dynamics. It was then discovered that spin dynamics develops while preserving the local spin length of the condensate [257]. This protection of the initial ferromagnetic character of the gas was attributed to an energy gap provided by sufficiently large spin-dependent contact interactions. This dynamical protection of ferromagnetism is the reason why trapped magnon excitations could be observed in the Cr experiments and described using the hydrodynamic equation of a ferrofluid, indicating that the Cr BEC behaves like a genuine ferrofluid; see section 6.2.2 and [260].

It should be noted that the protection of ferromagnetism in the ⁵²Cr condensates is tied to the relatively small strength of the DDI compared to spin-dependent contact interactions; it is likely that similar experiments performed (for example with Ln atoms) in the regime where dipolar interactions overwhelm spin-dependent contact interactions may lead to qualitatively different behaviour, since the DDI does not preserve the spin length. As we will see in section 7, working in optical lattices is another way to investigate a purely dipolar spin system. In this case, the collective spin length is indeed reduced during dynamics [578].

Beyond the case of BEC described above, recent work has been dedicated to the study of spinor systems of dipolar fermions [164, 251, 313]. This system is of great interest for the study of how the DDI affects superfluid pairing and the celebrated crossover from delocalised Cooper pairs in the BCS regime to a BEC of molecules [13, 127]. By choosing specific spin mixtures, spinor gases of constant magnetisation can be studied over long time scales, even for highly magnetic species of Er and Dy. This is due to the quantum statistical suppression of dipolar relaxation [164]; see also section 3.3.4. This feature has been recently used to realise a long-lived spin-1/2 spin-orbit-coupled degenerate dipolar Fermi gas [251]. This experiment will be discussed in section 6.4.1.

A spin-1/2 degenerate Fermi mixture of strongly magnetic atoms ¹⁶⁷Er was created [313], and the collisional behaviour of this mixture experimentally investigated. This experiment used a deep 3D optical lattice as a tool during the preparation of the mixtures, which serves to inhibit collisional losses induced by large magnetic field sweeps and the crossing of numerous FRs: A quantum-degenerate spin-polarised sample is first loaded into a deep 3D optical lattice at low magnetic field, such that double occupancy are precluded and fermions are spatial separated. The magnetic field can there be ramped up to a value where quadratic Zeeman shifts are significant, allowing spin state preparation using using rf sweeps, and later be ramped down to the low-field region, where the physics of the spin mixture is studied. To do so, the 3D lattice is slowly ramped down and the Fermi mixture loaded back into a 3D trap. As already observed for the spin-polarized case—see section 2.4)—the collisional properties of the spin mixtures as a function of the bias magnetic field amplitude, B, is marked by a large number of intra- and inter-spin FRs. Baier et al [313] mapped the FR spectra for B in the $[0,2]$G range via loss spectroscopy. Many narrow and overlapping features are observed. A comparatively broad and isolated interspin FR is also identified. Baier et al used this resonance to reach the strongly-interacting regime. Here, a large collisional stability of the balanced spin-mixture at $T/T_F\approx 0.3$ is observed, in particular on the repulsive side of the FR. This paves the way towards study of BEC-BCS physics in such systems.

Thermodynamics of a Bose gas with a spin degree of freedom is derived by including the single-particle magnetic energy in the Bose occupation factor. For non-interacting particles,

$$\begin{align} f_{k,m_s}=\frac{1}{\exp \left[ \left(\epsilon_{k,m_s}- \mu\right)/ k_\mathrm{B} T\right] -1}, \end{align} \tag{ 65 }$$

where µ is the chemical potential, and $\epsilon_{k,m_s}$ is the single-particle energy of the (trapped) states labelled by the index k, in the Zeeman state m_s . When magnetisation is free (case of dipolar particles) $\epsilon_{k,m_s}$ includes the linear Zeeman effect, whereas when magnetisation is fixed (as for example of Rb and Na atoms), the linear Zeeman effect is gauged out and $\epsilon_{k,m_s}$ includes only the quadratic Zeeman effect.

For the generic case of a F = 1 spinor Bose gas interacting solely via contact interactions (and thus at constant magnetisation and without the quadratic Zeeman effect), the phase diagram has been worked out by Isoshima et al [579] as a function of the temperature T and the magnetisation M. Two phase transitions were predicted. The first transition, at a critical temperature $T_1(M)$ separates the normal phase, and a phase where a condensate forms in the (most populated) stretched state. Below a second critical temperature $T_2(M)$, all other spin components condense simultaneously (see figure 32).

Zoom In Zoom Out Reset image size

First investigations of spinor physics with free magnetisation were performed using cold Cr atoms in an ODT. It was observed that the spin degrees of freedom equilibrated to the gas mechanical temperature [580]. At high temperature, the population of the different Zeeman states followed a distribution close to the Boltzmann distribution. The spin temperature was equal to the mechanical temperature, and thermal equilibrium between spin and orbital degrees of freedom was ascribed to magnetisation-changing dipolar collisions. Below the critical temperature for BEC, a spontaneous accumulation into the Zeeman state of lowest energy was observed. This led to a bimodal spin distribution, arising from BEC and the spontaneous accumulation of the atoms in the lowest-energy spin state.

By varying the temperature and the magnetic field, it was possible to map out the corresponding critical line $T_1(M)$ (see figure 32). However, the second predicted phase transition at $T_2(M)$ could not be studied. Indeed, when magnetisation is free, Bose stimulation towards the single-particle lowest-energy state insures that the BEC is only produced in this polarised state [580]. The second transition was recently observed in the case of sodium atoms with negligible DDI by Frapolli et al [581]. The authors also studied the impact of spin-dependent interactions on the double-condensation scenario.

In a later study, it was possible to reach a multicomponent phase with Cr atoms by applying rapid, forced evaporative cooling to a depolarised gas, thus performing the experiment fast compared to the time scale for magnetisation-changing collisions. However, the spinor BECs were very small due an interplay between Bose condensation and spin dynamics [582]. This result pointed out the difficulty of fully thermalising the spin degrees of freedom. This is a prominent effect to be taken into account for very large spin systems such as Dy and Er where all spin states must be saturated for a stable multicomponent BEC to be produced.

The connection between magnetic order and superfluid order and the BEC transition is an intriguing question that remains largely unexplored. While these orders are intrinsically connected due to Bose stimulation, it was predicted that strong spin-dependent interactions induce spin ordering at a finite T above the BEC transition [583]. In low dimensions, the connection between magnetic order and superfluidity promises to be especially interesting. In particular, in the 2D polarised case, superfluidity does not arise from condensation but by the pairing of vortices with opposite circulation through the Berezinskii–Kosterlitz–Thouless mechanism. This gives rise to a topological order [81, 584–586]. The case of a depolarised gas promises new scenarios for superfluidity, due to the existence of topological excitations such as half-quantum vortices that include the spin degree of freedom [587].

We note that the case of fermions will drastically differ from the scenario described above due to the presence of the additional energy scale provided by the Fermi energy. One consequence is that spontaneous demagnetisation of a noninteracting Fermi sea is expected provided the Larmor energy is smaller than the Fermi energy (Pauli paramagnetism). How interactions modify this picture is an open question and connected to Stoner instability.

One of the fascinating avenues of research for dipolar spinor physics is the investigation of the zero-temperature phase diagram at low magnetic fields when the Larmor energy is comparable to spin-dependent interactions. For the experiments performed with alkalis, spin-dependent contact interactions are extremely small compared to the linear Zeeman effect. In addition, the gas magnetisation remains constant for all practical purposes because contact interactions are isotropic, and anisotropic DDIs between alkali atoms are small. Consequently, the true ground state of these system, with free magnetisation at extremely low magnetic fields, has never been experimentally investigated.

Quantum gases made of strongly magnetic atoms, on the other hand, offer the new possibility of free magnetisation introduced by the DDI. Moreover, spin-dependent contact interactions provide an energy scale that corresponds to a magnetic field in the few 100 µG range for Cr atoms (compared to 10 µG for alkalis). It therefore becomes possible to investigate phase transitions driven by the competition between the Zeeman effect and spin-dependent interactions. For example, the phase diagram of Cr atoms has been calculated (ignoring the DDI) by [352, 588], and is expected to display a number of transitions separating ferromagnetic, cyclic, polar, and biaxial nematic phases as a function of the magnetic field. The phase diagram of Ln atoms is unknown to the best of our knowledge, though some have speculated [156]. We note that the case of Lns may be quite different, because contact interactions themselves contain an anisotropy associated with the large electronic angular momentum in the ground state.

The exploration of the spinor phase diagram at very low magnetic field was started with a Cr BEC in the $B\sim 100\,\mu$G regime [350]. While the BEC remains polarised in the lowest energy single-particle state for large enough magnetic fields, it was then shown that a dipolar BEC spontaneously depolarises when the Zeeman energy is quenched below a (density-dependent) critical value (see figure 33). It was furthermore shown that, while the DDI was presumably too small to significantly modify the observed phases, the dynamics of depolarisation was driven by dipolar interactions.

Zoom In Zoom Out Reset image size

A number of open questions remain after the first results described in [350]. One interesting question concerns the possibility of characterising the excitations that may be produced when the phase transition is crossed. As magnetisation is dynamically modified, the sample should start rotating, in the spirit of the Einstein–de-Haas effect. Such rotation could not be observed. However, simulations indicate that the threshold for demagnetisation could originate in the resonant dynamics associated with this exchange between spin and orbital momentum [589].

Furthermore, the impact of the DDI on the phases at low field also remains unexplored. Note that it would be especially interesting to study the regime where the Zeeman energy is smaller than the MF energy associated with the dipolar field, which could provide a scenario for symmetry breaking and modify the possible quantum phases [571]. The requirements in terms of magnetic field, below 50 µG, are demanding, yet within reach of state-of-the-art experiments [590].

6.3.3.1. Demagnetisation cooling.

6.3.3.2. Purification of a BEC by spin-filtering.

Another closely related cooling scheme, which is efficient below quantum degeneracy, has been demonstrated by Naylor et al [594]. This scheme also relies on redistributing population between different spin states, with free magnetisation. The key idea is that only non-condensed atoms may populate spin-excited states since Bose thermodynamics enforce a fully polarised condensate; see section 6.3.1. Therefore, expelling spin-excited atoms from the trap provides a way to engineer losses specific to non-condensed atoms, thus cooling and purifying the condensate.

The scheme, experimentally demonstrated using ⁵²Cr atoms, starts with a partially condensed Bose gas polarised in the lowest energy spin state. Demagnetisation of the thermal component is triggered by lowering the magnetic field, such that the Larmor energy is comparable to the thermal energy. Then, the spin-excited thermal components produced by magnetisation-changing collisions are filtered out by a magnetic field gradient. It was found that when the initial BEC fraction is high enough, this scheme ends up with a polarised BEC in the lowest energy state with an increased BEC fraction. This provides a thermodynamic cycle that, in principle, decreases the entropy by a factor of up to $(2s+1)^{3/4}$ (where s is the atomic spin), and could also be in principle repeated. The obtained reduction in entropy was typically a factor of 2 for one cycle in the Cr experiment.

In the experiment, technical limitations arose from the difficulty to control the magnetic field at very low values. This limitation is directly related to the fact that dipolar gases are sensitive to the linear Zeeman effect. It was suggested in [594] that the purification of a BEC by spin-filtering could be extended to non-dipolar species by using spin-dependent contact interactions and spin-exchange interactions. Both the dipolar and the non-dipolar cooling schemes were theoretically explored in [595]. A related cooling scheme was demonstrated with Rb atoms [596].

Note that, despite their highly magnetic character, such cooling techniques could be difficult to apply to Lns, because the Feshbach spectrum of these atoms might be too dense at low field, below typically 10 mG. Up to now, these techniques haven not yet been reported using Ln species.

Zoom In Zoom Out Reset image size

Both approaches described here are ways to use the spin degrees of freedom to efficiently store and remove entropy from a gas. We point out that ideas using the spin degrees of freedom have already been developed in the context of lattice-based large-spin Fermi gases [565], mixtures of Bose gases [597], and Bose and Fermi gases [598]. Alternatively, it has also been suggested that the spin degrees of freedom, and/or spin-changing collisions, may be employed for thermometry purposes down to extremely low temperatures. Such temperatures are typically impossible to measure in BECs with usual thermometric techniques [596, 599].

6.4.1.1. Synthetic spin–orbit coupling and artificial gauge fields.

The coupling of a particle's momentum and spin underlies many important phenomena in quantum systems, and in particular in electronic, solid-state systems [600]. For example, topologically nontrivial materials often require the coupling of the electron spin to momentum [601, 602]. Here, SOC typically arises from the movement of electrons through the crystal electric field. The effect of a magnetic field on charged particles also yields intriguing phenomena, the most famous of them being the quantum Hall effect [603]. This effect has its roots in the special form of the Lorentz force, which, at a quantum level, yields a kinetic term of the form:

$$\begin{align} \hat{H}_k = \frac{[\hat{\mathbf{p}}-\hat{\mathbf{\mathcal{A}}}(\hat{\mathbf{r}})]^2}{2m}, \end{align} \tag{ 66 }$$

where $\hat{\mathbf{p}}$ is the canonical momentum and $\hat{\mathbf{\mathcal{A}}}(\hat{\mathbf{r}})$ is the vector potential or gauge potential. Even more exotic behaviour, such as the fractional quantum Hall states [604], Laughlin liquids [605], and topological superconductors harbouring Majorana excitations [151] can arise if $\hat{\mathbf{\mathcal{A}}}(\hat{\mathbf{r}})$ is not just a classical vector (i.e. three complex numbers), but a set of three operators $\hat{A}_x(\hat{\mathbf{r}})$, $\hat{A}_y(\hat{\mathbf{r}})$, and $\hat{A}_z(\hat{\mathbf{r}})$ that do not commute. This yields a non-Abelian field.

Interestingly, ultracold atoms enable the combination of SOC and gauge potentials by the engineering of a realisation of equation (66). This mimics the effect of a Lorentz force on neutral matter while, additionally, the set of operators $\hat{A}_i$ connect to the atomic spin operators. If realised in high dimensions (i.e. greater than one), non-Abelian gauge fields may then be synthesised. This technique relies on laser-field dressing [160]. In the following, we will review this SOC protocol and discuss its extension to magnetic atoms. We note that several other strategies, not connected to SOC, have been developed over time to realise artificial gauge fields on neutral atoms, e.g. by engineering trap rotations [606].

As introduced in section 1.3.3, and reviewed in detail in the above section 6.3, the DDI provides intrinsic SOC through the dipolar relaxation terms of equation (5). However, this is not generally the sort of SOC that results in the physics of artificial gauge fields mentioned above. It gives rise to an interaction term that does not conform to that of equation (66). (See [354] for a method that does.) In contrast, light-induced SOC schemes applied to ultracold quantum gases do realise a Hamiltonian term of the form of equation (66). The basic idea can be formulated as follows: Light fields are applied to couple Zeeman sublevels via a two-photon Raman transition. Recoil momentum is transferred as the optically coupled spin flips. The adiabatic evolution of these dressed states as the atom moves in the Raman field creates a synthetic gauge field [607]. This field provides SOC when these states form a (near) degenerate manifold. The general form of the SOC Hamiltonian under the 1D Raman coupling along the direction $\hat{x}$ of two Zeeman ground states is:

$$\begin{align} \hat{H}(\mathbf{p}) = \frac{[\mathbf{p}-\hbar\mathbf{k}\hat{\sigma}_z]^2}{2m} + \frac{\hbar\Omega}{2}\hat{\sigma}_x + \frac{\hbar\Delta}{2}\hat{\sigma}_z, \end{align} \tag{ 67 }$$

where $\hat{\sigma}_{x,y,z}$ are the Pauli spin matrices. Ω and Δ are the two-photon Raman Rabi frequency and detuning, respectively. The form of this coupling is the sum of Rashba ($k_z\sigma_z+k_y\sigma_y$) and Dresselhaus ($k_z\sigma_z-k_y\sigma_y$) SOC with equal weights. The resulting dispersion relation is of the form of a double well centred at k = 0, as is typical of SOC systems. That is, the momentum of spin up atoms is oriented in the opposite direction from that of spin down atoms. In the 1D case described by equation (67), $\hat{A}_z = \hbar\mathbf{k}\hat{\sigma}_z$ and $\hat{A}_{x,y} = 0$, which does not describe a non-Abelian coupling because the single component of $\hat{A}$ commutes with itself. Generalisations to 2D and 3D with a truly non-Abelian gauge potential have been proposed [160] and realised in the case of 2D SOC [608, 609].

6.4.1.2. Spin-orbit-coupling in quantum gases: limitations and prospects.

6.4.1.3. Spin-orbit-coupling in assemblies of magnetic atoms.

Burdick et al [251] first reported the realisation of SOC in DFGs of Ln atoms using Dy. The SOC was induced using Raman light near the 741 nm transition. Its 1.8 kHz width reduces spontaneous emission rates below the background lifetime of the gas with only modest ∼GHz detuning, which is large compared to the hyperfine splittings, but smaller than the fine-structure; see [214, 238]. The lifetime of SOC gases was then limited, not by spontaneous emission, but rather by dipolar relaxation.

Due to the suppression of inelastic dipolar scattering (see section 3.3.4 for details), relatively long-lived SOC Fermi gases could be realised while working at a bias field of a few tens of G. Such magnetic fields were required, not only because the suppression effect scales as $\sqrt{B}$, but also because the isolation of an effective pseudospin $1/2$ requires a large the quadratic Zeeman shift obtainable only at large B. However, as discussed in section 2.4, the extreme density of FRs in fermionic Lns, increasing with B, complicates the situation, and [251] even reported the overlapping of resonances at these fields. Nevertheless, a field region at 33.846(5) G was found that allowed fermionic spin mixtures to live enough to allow spin-orbit coupled ¹⁶¹Dy gases to be created with a lifetime of ∼400 ms at a Raman coupling strength of $\hbar\Omega = 1E_R$; see figure 35(b). At this field, quadratic Zeeman shifts ensures that the Raman fields couple only the $m_F = -21/2$ and $m_F = -19/2$ states, providing maximum fermionic suppression of dipolar relaxation; see figure 35(a). Note that by contrast, the lifetimes of SOC ⁴⁰K and ⁶Li were lower by ∼10 and 100, respectively [609, 611, 612]. Moreover, the SOC ¹⁶¹Dy lifetime is similar to that of free-space bosonic SOC alkali gases [610] and ∼10× longer than that achieved in a bosonic lattice system [623].

Zoom In Zoom Out Reset image size

Burdick et al [251] also reported that Raman-coupled bosonic Dy has a short lifetime of less than 10 ms at low B field. This is just as short as ⁶Li, demonstrating the importance of fermionic statistics in preventing fast relaxation. Finally, the effect of dipolar interactions on Rabi oscillations was observable, highlighting the interacting dipolar character of the fermionic SOC system.

6.4.1.4. Dipolar relaxation of Raman-dressed spins.

6.4.1.5. Realisation of a synthetic quantum Hall system using magnetic atoms.

In recent works using ultracold bosonic gases of ¹⁶²Dy and ¹⁶⁴Dy, a spin-dependent Hamiltonian was generated using laser light close to the 626 nm atomic line [180, 181]. The laser field has a linear polarisation along x while the magnetic field points along z. The resulting atom-light interaction yields an effective coupling term $\hbar\omega \hat{J}_x^2$, with $\hat{J}$ being the total atomic angular momentum. This scheme, similar to that previously investigated with a large ensemble of room temperature Cs atoms [535], realises a so-called one-axis twisting Hamiltonian [155], which has also been investigated in many other systems [630–633]. Typical coupling strengths ω are on the MHz scale and can be made to greatly exceed the kHz-scale Larmor precession frequency as well as the dipolar relaxation rate.

Highly non-classical spin states could be generated under the effect of this coupling term [155]. Chalopin et al [180] experimentally demonstrates the creation of Schrödinger cat states resulting from the coherent superposition of the two stretched spin states $|\pm J\rangle$. In contrast to ensembles of $s = 1/2$ particles, for which entanglement can only be generated between different atoms of the ensemble, here the non-classical states arise because of the entanglement between electrons within each individual atom [535]. This possibility is a key feature provided by large-spin systems. These states appear in the time evolution of the atomic spin (at the single particle level) after a quarter of coupling period, $T = 2\pi/\omega$. Furthermore, several revivals of the cat state spaced by $T/2$ as well as a repolarisation into the Fock states $|J\rangle$ and $|- J\rangle$ at intermediate times (half integer multiples of T) occur during the time evolution. Coherence is seen to decay under the effect of classical B-field noise, to which the cat state is more sensitive than a classical (coherent) state. Nevertheless, the cat state's lifetime reaches up $60\mu\mathrm{s}$. Focusing on a similar dynamics but shorter time scales (set by $\tau = 1/\left(\sqrt{2J}\omega\right)$), [181] probes the formation of non-Gaussian 'oversqueezed' states.

In a quasi-1D geometry, the external motion of the atoms is effectively frozen in two directions when $\hbar\omega_\perp \gg k_\mathrm{B}, \mu$. The effective form of the DDI in quasi-1D geometries can be derived from the 3D expression, equation (1), by integrating out the transverse degrees of freedom, as shown in [136–138]. Under a single-mode approximation, this yields:

$$\begin{align} U^\textrm{1D}(u) = V(\theta)\left[ {V}^\textrm{1D}(u) - \frac{8}{3} \delta (u) \right], \end{align} \tag{ 70 }$$

where

$$\begin{align} V(\theta) = \frac{\mu_0 \mu^2}{4\pi} \frac{1-3\cos^2{\theta}}{4 a_\perp^3}, \end{align} \tag{ 71 }$$

and

$$\begin{align} {V}^\textrm{1D}(u) = -2|u| + \sqrt{2\pi} (1+u^2) e^{u^2/2} \textrm{erfc}(|u|/\sqrt{2}). \end{align} \tag{ 72 }$$

Here, $u = x/a_\perp$, and $\textrm{erfc}(u)$ is the complementary error function. There are two contributions to the short-range part of the 1D DDI. The first is the δ-function term. It comes from the point limit of an extended dipole [668] and has an opposite sign to ${V}^\textrm{1D}(u)$. The second contribution is an effective delta-function term that arises from the fact that ${V}^\textrm{1D}(u)$ becomes more sharply peaked as $a_\perp$ shrinks [137, 138]. At distances $|x|\gg a_\perp$, ${V}^\textrm{1D}(u) \rightarrow 4/|u|^3$. This x⁻³ long-range potential is just like the DDI in 3D, but with diminished magnitude at ranges on the order of $a_\perp$. Within a distance set by $a_\perp$, the DDI ceases to seem 1D to atoms that approach near each other. The underlying 3D spatial nature is manifest and a small attractive (repulsive) contribution to the DDI emerges from the part of their wavefunctions that extend transversely by $a_\perp$, even if the long-range interaction along the 1D axis is repulsive (attractive) [137, 138].

The DDI in 1D thus has both long and short-range components [137, 138, 662]. The chemical potential remains intensive in 1D, which is indicative of a short-range interacting system. However, like a long-range interacting system, there is no asymptotic phase shift (scattering length) that can be defined for two-body collisions in 1D [662]. Away from collisional resonances, the short-range part of the DDI can add to the van der Waals contact pseudopotential to yield a total hard-core contact interaction strength $g_\textrm{1D}$, as considered in [140, 669]. In addition to the FR-tunability of the van der Waals contact interaction (discussed in previous sections), the 1D DDI provides wide tunability in the properties of 1D many-body systems because both the short and long-range parts of the DDI in 1D can be set to a positive or negative sign, or made to completely vanish at the magic angle $\theta_m = 54.7^\circ$.

Dipolar 1D gases have been created by loading BECs of either ⁵²Cr or ¹⁶²Dy into 2D optical lattices creating arrays of 1D tubes [140, 162]. In this geometry, the DDI affects different aspects of the interactions: (a) the dependence of $g_\textrm{1D}(\theta)$ on the short-range component of the 1D DDI; (b) the intratube long-range DDI; and (c) the mutual long-range DDI of atoms in nearby tubes. The intertube DDI does not change the dimensionality of the system, but rather is a multichannel effect that couples different flavors of 1D quasiparticles [140]; i.e. coupled 1D systems remain 1D because the phase space degrees of freedom are still 1D when there is no transport between tubes.

The first experimental work with 1D dipolar gases used Cr and explored magnetisation-changing collisions and their suppression in 2D optical lattices [162]; see figure 6 and section 3.3.5 for more details. The first experiments with strongly interacting gases entering the fermionised γ > 1 regime used Dy with ∼50 Dy atoms per tube [140].

Tang et al [140] investigated the thermalization of a near-integrable quantum system. The Lieb–Liniger model is an example of a quantum integrable Hamiltonian, in which there exists an extensive number of conserved quantities, resulting in regular, non-chaotic, non-thermalizing dynamics for any value of $g_\textrm{1D}$ [139, 640]. This is due to the 1D geometrical constraint imposed on interactions: two-body interactions permit non-diffractive scattering between nearest neighbours only, leaving the set of incoming momenta unchanged with respect to the set of outgoing momenta [670]. As a result, the scattering matrix obeys the Yang–Baxter equation, a sufficient condition for integrability [640, 670–672], since the binary collisions are unable to alter the distribution functions. Therefore, an empirical 'smoking gun' for integrability would be the persistence of an out-of-equilibrium momentum distribution beyond the intrinsic dynamical time scale.

The group of David Weiss at Penn State University observed this persistent non-equilibrium momentum distribution in 1D gases of contact-interacting Rb atoms, thereby experimentally establishing the integrability of the Lieb–Liniger model. Their experiment relied on quantum quenches using what they called a 'quantum Newton's cradle' [655]: atoms in 1D traps were set in motion using a Bragg diffraction pulse as a system quench. The atoms oscillated in counterpropagating packets and collided twice each period under the strong-coupling condition of $\gamma \gg 1$. This is akin to the desktop Newton's cradle toy, except that instead of the metal spheres reflecting upon each collision, the Rb atoms could also pass through one another as a manifestation of the quantum nature of the system. Rather than rapidly come to a steady state (i.e. a stationary Gaussian momentum distribution), the packets were observed to oscillate many times (far longer than what one would expect in a 3D gas) before the onset of heating from spontaneous emission. Thus, the Weiss group established that a strongly interacting integrable quantum system could be studied in the laboratory. Although both longitudinal confinement and transverse (virtual) motional excitation break integrability, these detrimental effects are suppressed in the $\gamma \gg 1$ regime [673–675].

Tang et al [140] used the magnetic DDI in a 1D Dy gas as the controllable interaction with which to break integrability in a dipolar version of the quantum Newton's cradle experiment. It is known from classical physics that magnetic spheres cause the motion of a toy Newton's cradle to be chaotic due to the long-range interaction; the long-range DDI should likewise break integrability in the quantum systems, since it allows for 'diffractive' collisions among atoms. Indeed, the DDI allows for 'diffractive' collisions among atoms: i.e. in addition to two colliding atoms swapping their momentum, as in the case of contact collisions (which are non-diffractive), their momentum may also be imparted to a third particle. This violates the Yang–Baxter condition, and hence breaks integrability. Since the DDI strength falls off slowly in space, the probability for a diffractive interaction is not suppressed by the usual need for three particles to be at the same place at once. Furthermore, by simply controlling θ, the strength of the integrability-breaking perturbation can be tuned. Figure 37 depicts the dipolar quantum Newton's cradle.

Zoom In Zoom Out Reset image size

The dipolar quantum Newton's cradle opens new avenues to explore how quantum thermalization arises upon the introduction of a perturbation that lifts integrability. For example, theoretical consensus is lacking regarding whether relaxation involves two distinct timescales or three for strongly interacting quantum near-integrable systems (i.e. first prethermalization, then a prethermal plateau of a near steady-state prethermal state, and finally a decay to the thermal state) [676–681]. More generally, thermalisation of a near-integrable system is an old question with a celebrated answer in the realm of classical physics [682]. In the 1950s, Kolmogorov et al developed what became known as KAM theory, establishing that chaos, and hence thermalization, sets in as a smooth crossover as the strength of a nonlinear, integrability-breaking perturbation is increased [683]. In the quantum realm, it has long been wondered whether any meaningful analogue to KAM theory exists, and a general theory for quantum thermalization in near-integrable systems is lacking despite much work [651, 676, 677, 679–681, 684–693].

In their experiment to investigate this issue, Tang et al [140] observed relaxation and thermalization dynamics in the time evolution of the momentum distribution of the kicked 1D gas. The distance to a thermal distribution is quantified and is observed to undergo two distinct exponential decay regimes. The first evolution is a prethermalization (dephasing) decay, and the second consists of the relaxation from the prethermal state towards a Gaussian thermal distribution determined by the Gibbs ensemble.

Tang et al found that the thermalization rate in this second step obeys a simple scaling expression that depends only on θ without any free parameters. The expression accounts for the DDI perturbation using Fermi's golden rule formula, multiplied by the cross section for short-range interactions to occur. This describes the dominant integrability-breaking diffractive interaction term: that of two atoms scattering via the contact interaction while simultaneously interacting with a third atom at long range via the DDI. It is remarkable that such a simple description appropriately describes the thermalization dynamics of such a strongly interacting system near an integrable point.

The experiment was supported by an exact diagonalization calculation of a two-rung XXZ quantum magnetism model [694] with long-range hopping serving as the perturbation. Similar two-stage thermalization dynamics was observed despite the microscopic dissimilarity of the systems.

Future work could attempt to understand how generally applicable is this simple, Fermi golden rule-like expression for a wider variety of quantum quench experiments. One could also explore how more sophisticated descriptions based on, e.g. adaptations of generalised hydrodynamics, may better describe the dipolar 1D near-integrable system. More generally, this work sets the stage for a wide array of inquiries into the physics of strongly interacting quantum systems near integrability. In the next section, we describe one such investigation, that of quantum many-body prethermal 'scar' states, enabled by the realisation of dipole-stabilised excited 1D quantum gases.

As mentioned above, the TG state is one in which 1D bosons with a divergently repulsive coupling strength behave like fermions. The two-body wavefunction vanishes when the atomic positions coincide, and so possesses exclusion correlations as if it were an ideal Fermi gas. A quantum quench of the TG gas prepares an eigenstate of the attractive Lieb–Liniger model, the so-called super-TG gas mentioned above in section 7.1. The TG and sTG state wavefunctions and energies are smoothly connected as a function of $g_\textrm{1D}$, so while the effective model exhibits a discontinuity, the system remains adiabatic throughout the parameter quench. The result is a highly excited sTG state wherein the 1D bosons that attractively interact at short range behave as if they were ground-state fermions repulsively interacting at long range. That is to say, the bosons develop even stronger exclusion correlations than the ideal Fermi gas. The two-body wavefunction node inherited from the quench from the TG state extends into a pair of nodes separated by the length $a_\textrm{1D}$, effectively introducing a rigid exclusion zone similar to the classical hard rod model. Moreover, its excited-state many-body wavefunction exhibits many more nodes than the TG wavefunction, which enhances the stiffness of the sTG gas [141, 695].

Such metastable attractive gases typically collapse into bound states (cf the 'Bose-nova' implosion of strongly attractive BECs in 3D [105]). By contrast, however, the strong antibunching correlations in the sTG gas wavefunction prevent atoms from approaching each other close enough to bind into cluster-like states. The sTG gas is expected to never collapse, despite the multitude of bound states of lower energy. The deeper reason for this surprising metastability lies in the aforementioned integrability of the Lieb–Liniger model, for any $g_{1D}$ value. That is, the sTG state is also a solution to the Bethe ansatz equations of the Lieb–Liniger model [644].

Experimentally, however, the sTG gas does collapse for attractive interactions weaker than those in the unitary $g_\textrm{1D}\rightarrow -\infty$ regime, as first observed in the nondipolar Cs system [142]. While the system remains effectively integrable in the unitary regime, since the interaction term dominates all others in the perturbed Lieb–Liniger Hamiltonian, the system can collapse when the interactions become sufficiently weak outside the unitary regime. That is, the breaking of the integrability of the Lieb–Liniger model becomes manifest once the system is no longer unitary near the resonance. Experimental imperfections such as the longitudinal harmonic trap and the virtual excitation of the higher motional bands of the transverse 2D optical lattice forming the array of 1D traps all break integrability, if only weakly [143, 673–675]. These yield nonzero matrix elements coupling the initial state and molecular/cluster-like bound states in the near-integrable system that enable the collapse of the sTG state at finite negative $g_\textrm{1D}$ coupling strength.

This may be understood in analogy to the classical 1D gas of hard rods of length $a_\textrm{1D}$. The rod length $a_\textrm{1D}$ extends an exclusion zone that prevents rods from overlapping. When more rods N are stuffed into a box of length L than can fit end-to-end (i.e. when $Na_\textrm{1D}>L$), then the system 'collapses' by kinking the rods out of 1D alignment. In the quantum system, $a_\textrm{1D}$ is the 1D scattering length definable through $g_\textrm{1D} = 2\hbar^2a/ma_\perp^2 = -2\hbar^2/ma_\textrm{1D}$. Both $g_\textrm{1D}$ and $a_\textrm{1D}$ are tunable using a confinement induced resonance (CIR) [641]. These arise in quasi-1D traps due to the open-channel-like role of higher-energy transverse motional states of the trap and were first observed in Cs [696]. The CIR provides tunability through the dependence of a on B-field near a 3D FR:

$$\begin{align} g_\textrm{1D}(B)=\frac{2\hbar^2 a(B)}{ma_\perp^2}\frac{1}{1-C a(B)/a_\perp}, \end{align} \tag{ 73 }$$

where $a_\perp$ is the transverse oscillator length and C ≈ 1.46 [641].