%~Mouliné par MaN_auto v.0.27.3 2024-03-18 10:16:55
\documentclass[CRPHYS, Unicode, XML, screen,biblatex]{cedram}
%
%\TopicEN{Unsolicited articles}
%\TopicFR{Articles spontanés}
\addbibresource{crphys20221130.bib}
\usepackage{graphicx}
\usepackage{dcolumn}
\usepackage{bm}
\usepackage{amssymb}
%\usepackage{epsfig}
\usepackage{verbatim}
\usepackage{cases}
\usepackage{mathrsfs}
\definecolor{darkgreen}{rgb}{0.0, 0.5, 0.0}
\newcommand{\bea}{\begin{eqnarray}}
\newcommand{\eea}{\end{eqnarray}}
\newcommand{\ee}{\end{equation}}
\newcommand{\bi}{\begin{itemize}}
\newcommand{\ei}{\end{itemize}}
\newcommand{\ds}{\displaystyle}
\newcommand{\rr}{\mathbf{r}}
\newcommand{\kk}{{\mathbf{k}}}
\newcommand{\KK}{{\mathbf{K}}}
\newcommand{\uu}{\mathbf{u}}
\newcommand{\RR}{\mathbf{R}}
\newcommand{\LL}{\mathbf{L}}
\newcommand{\CC}{\mathbf{C}}
\newcommand{\cc}{\mathbf{c}}
\newcommand{\JJ}{\mathbf{J}}
\newcommand{\XX}{\mathbf{X}}
\newcommand{\vn}{\mathbf{0}}
\newcommand{\ra}{\rangle}
\newcommand{\la}{\langle}
\newcommand{\ktyp}{k_{\rm typ}}
\newcommand{\gr}{\bm \nabla}
\newcommand{\UP}{\uparrow}
\newcommand{\down}{\downarrow}
\newcommand{\Ar}{\mathcal{A}}
\newcommand{\Vr}{\mathcal{V}}
\newcommand{\Cr}{\mathcal{C}}
\newcommand{\Fr}{\mathcal{F}}
\newcommand{\Rr}{\mathcal{R}}
\newcommand{\Wr}{\mathcal{W}}
\newcommand{\Er}{\mathcal{E}}
\newcommand{\Sr}{\mathcal{S}}
\newcommand{\Nr}{\mathcal{N}}
\newcommand{\Or}{\mathcal{O}}
\newcommand{\Mr}{\mathcal{M}}
\newcommand{\Dr}{\mathcal{D}}
\newcommand{\Kr}{\mathcal{K}}
\newcommand{\Oo}{{\bm \Omega}}
\newcommand{\rrho}{\bm \rho}
\newcommand{\md}{\mathsf{m}}
\newcommand{\dK}{d^3\!K}
\newcommand{\dr}{d^3\!r}
\newcommand{\dc}{d^3\!c}
\newcommand{\dC}{d^3\!C}
\newcommand{\bbN}{\mathbb{N}}
\newcommand{\bbC}{\mathbb{C}}
\newcommand{\bbR}{\mathbb{R}}
\newcommand{\bl}{}
\newcommand{\rmm}{\mathrm{m}}
\newcommand{\rMax}{\mathrm{Max}}
\newcommand{\ZRM}{\mathrm{ZRM}}
\newcommand{\rIm}{\mathrm{Im}}
\newcommand{\rMin}{\mathrm{Min}}
\newcommand{\rmho}{\mathrm{ho}}
\newcommand{\gre}{\bl}
\newcommand{\brown}{\bl}
\newcommand{\dd}{d_2}
\newcommand{\ddd}{{\bl d_3}}
%\newcommand{\bcbp}{\begin{comment}}
%\newcommand{\ecbp}{\end{comment}}
\makeatletter
\def\cdr@preabstracthook{%
\vspace*{-5pt}
\vskip2pt plus 3pt minus 1pt
}
\makeatother
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\graphicspath{{./figures/}}
\newcommand*{\mk}{\mkern -1mu}
\newcommand*{\Mk}{\mkern -2mu}
\newcommand*{\mK}{\mkern 1mu}
\newcommand*{\MK}{\mkern 2mu}
%\hypersetup{urlcolor=purple, linkcolor=blue, citecolor=red}
\newcommand*{\relabel}{\renewcommand{\labelenumi}{(\theenumi)}}
\newcommand*{\romanenumi}{\renewcommand*{\theenumi}{\roman{enumi}}\relabel}
\newcommand*{\Romanenumi}{\renewcommand*{\theenumi}{\Roman{enumi}}\relabel}
\newcommand*{\alphenumi}{\renewcommand*{\theenumi}{\alph{enumi}}\relabel}
\newcommand*{\Alphenumi}{\renewcommand*{\theenumi}{\Alph{enumi}}\relabel}
\let\oldtilde\tilde
\renewcommand*{\tilde}[1]{\mathchoice{\widetilde{#1}}{\widetilde{#1}}{\oldtilde{#1}}{\oldtilde{#1}}}
\let\oldhat\hat
\renewcommand*{\hat}[1]{\mathchoice{\widehat{#1}}{\widehat{#1}}{\oldhat{#1}}{\oldhat{#1}}}
\let\oldexists\exists
\renewcommand*{\exists}{\mathrel{\oldexists}}
\let\oldforall\forall
\renewcommand*{\forall}{\mathrel{\oldforall}}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\RubricEF{Rubric information missing}{Information rubrique manquante}
\alttitle{Contact à trois corps pour les fermions. \\I. Relations générales}
\keywords{\kwd{Unitary gas}
\kwd{fermions}
\kwd{cold atoms}}
\altkeywords{\kwd{Gaz unitaire}
\kwd{fermions}
\kwd{atomes froids}}
\title{Three-body contact for fermions. \\I. General relations}
\alttitle{Contact à trois corps pour les fermions. \\I. Relations générales}
\keywords{\kwd{Unitary gas}
\kwd{fermions}
\kwd{cold atoms}}
\altkeywords{\kwd{Gaz unitaire}
\kwd{fermions}
\kwd{atomes froids}}
\author{\firstname{F\'elix} \lastname{Werner} \CDRorcid{0000-0002-5631-9024}}
\address{Laboratoire Kastler Brossel, Ecole Normale Sup\'erieure - Universit\'e PSL, CNRS, Coll\`ege de France, Sorbonne Universit\'e, 75005 Paris, France}
\email[F. Werner]{werner@lkb.ens.fr}
\author{\firstname{Xavier} \lastname{Leyronas} \CDRorcid{0000-0002-9499-6800}}
\address{Laboratoire de Physique de l'Ecole Normale Sup\'erieure, ENS - Universit\'e PSL, Sorbonne Universit\'e, Universit\'e Paris Cité, CNRS, 75005 Paris, France}
\email[X. Leyronas]{xavier.leyronas@phys.ens.fr}
\thanks{F. W. acknowledges the hospitality of the Aspen Center for Physics, and support from ERC (project \emph{Critisup2}, H2020 Adv-743159) and ANR (project \emph{LODIS}, ANR-21-CE30-0033)}
\begin{abstract}
We consider the resonant Fermi gas, that is, two-component fermions in three dimensions interacting by a short-range potential of large scattering length. We introduce a quantity, the three-body contact, that determines several observables. Within the zero-range model, the number of nearby fermion triplets, the large-momentum tail of the center-of-mass momentum distribution of nearby fermion pairs, as well as the large-momentum tail of the two-particle momentum distribution, are expressed in terms of the three-body contact. For a small finite interaction range, the formation rate of deeply bound dimers by three-body recombination, as well as the three-body contribution to the finite-range correction to the energy, are expressed in terms of the three-body contact and of a three-body parameter. This three-body parameter, which vanishes in the zero-range limit, is defined through the asymptotic behavior of the zero-energy scattering state at distances intermediate between the range and the two-body scattering length. In~general, the three-body contact has different contributions labeled by spin and angular momentum indices, and the three-body parameter can depend on those indices. We also include the generalization to unequal masses for $\uparrow$ and $\downarrow$ particles. With respect to the relation between three-body loss rate and number of nearby triplets stated in [Petrov, Salomon and Shlyapnikov, \emph{PRL} {\bf 93}, 090404 (2004)], the present work adds a derivation, expresses the proportionality factor in terms of the three-body parameter, and includes the general case where there are several contributions to the three-body contact and several three-body parameters.
\end{abstract}
\begin{altabstract}
Nous considérons le gaz de Fermi résonnant, à savoir des fermions avec deux états internes à trois dimensions avec des interactions à courte portée de grande longueur de diffusion. Nous introduisons une quantité, le contact à trois corps, qui détermine plusieurs observables. Pour le modèle de portée nulle, le nombre de triplets de fermions proches, la queue de la distribution selon l'impulsion du centre de masse des paires de fermions proches, ainsi que la queue de la distribution en impulsion à deux particules, sont exprimées en termes du contact à trois corps. Pour une portée non nulle, le taux de formation de dimères fortement liés par recombinaison à trois corps, ainsi que la contribution à trois corps à la correction de portée finie à l'énergie, sont exprimées en termes du contact à trois corps et d'un paramètre à trois corps. Ce paramètre à trois corps, qui tend vers zéro dans la limite de portée nulle, est défini \emph{via} le comportement asymptotique de l'état de diffusion d'énergie nulle à des distances intermédiaires entre la portée et la longueur de diffusion à deux corps. En général, le contact à trois corps a différentes contributions repérées par des indices de spin et de moment cinétique, et le paramètre à trois corps peut dépendre de ces indices. Nous incluons aussi la généralisation à des masses différentes pour les particules $\uparrow$ et $\downarrow$. Par rapport à la relation donnée dans [Petrov, Salomon et Shlyapnikov, \emph{PRL} {\bf 93}, 090404 (2004)] entre taux de pertes à trois corps et nombre de triplets de fermions proches, le présent travail ajoute une dérivation, exprime le facteur de proportionnalité en termes du paramètre à trois corps, et inclut le cas général où il y a plusieurs contributions au contact à trois corps et plusieurs paramètres à trois corps.
\end{altabstract}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\dateposted{2024-04-22}
\begin{document}
\maketitle
\section{Introduction}
Over the last twenty years, the two-component Fermi gas with zero-range interactions in three dimensions has become one of the most extensively studied quantum many-body problems. One considers particles with two internal states (denoted $\UP$ and $\down$) and an interaction of vanishing range characterized by its $s$-wave scattering length $a_2$. When $1/a_2$ changes from $-\infty$ to $+\infty$, the interaction changes from weakly to strongly attractive, leading to the BCS to BEC crossover. The strongly correlated regime is reached in the central region of the crossover, around the unitary limit $1/a_2=0$. While the model was historically introduced as a theoretical abstraction~\cite{Leggett_Chap_1980_book,Haussmann_Z_Phys,Haussmann_PRB}, it accurately describes ultracold gases of fermionic atoms in two hyperfine states near a Feshbach resonance, which are the subject of numerous experimental studies, see e.g.~\cite{RevueTrentoFermions,ZwergerBook2,
ChapZwergerZwerger,ThomasStronglyInteractingGaz,
GrimmCrossover,JinPairCondensate,KetterlePairCondensate,
SalomonCrossover,GrimmModes,ShinPhaseDiag,SylEOS,
NirEOS,KuEOS,LENSRepulsivePolaron,ValeGoldstone,
EsslingerWiedemannFranz,Boiling_MIT,
ZwierleinSound,SagiPolaronMol,MoritzBragg3D,
NavonRecomb,HuletClosedChannel,AustraliensC,
AustraliensT,ValeC_precise,JinUnivRel,NirEOS,
Jin_C_homogeneous,SagiContactRF,ValeC_T_hom,
ZwierleinC,LaurentC}.
From a theoretical viewpoint, this resonant Fermi gas is a difficult problem. As for most strongly correlated many-body problems in dimension $>1$, numerical methods are generally the only option to make precise predictions. Furthermore the zero-range nature of the interactions typically constitutes an additional difficulty for numerical computations. But zero-range interactions also give rise to specific exact relations, called Tan relations, involving a ubiquitous quantity, the two-body contact $C_2$~\cite{TanEnergetics,TanLargeMomentum,ChapLeggettBref,LeChapitreIn2,ChapBraatenBref,RanderiaRF_arxiv,ZwergerRFLong,RanderiaRF,BaymRF,BraatenC,TanViriel,FelixViriel,ZwergerRF,BraatenLong,WernerTarruellCastin,ZhangLeggettUniv,CombescotC,WernerCastinRelationsArxiv,HofmannAnomaly,TanTrapFunctional,Moelmer,WernerCastinRelationsFermions}. In particular, $C_2$ determines the probability to find two particles close to each other. For the resonant Fermi gas, $C_2$ was measured and computed in numerous studies, see e.g.~\cite{WernerTarruellCastin,HuletClosedChannel,AustraliensC,
AustraliensT,ValeC_precise,JinUnivRel,NirEOS,
Jin_C_homogeneous,
SagiContactRF,ValeC_T_hom,ZwierleinC,LaurentC} and~\cite{Baym_C,Hu_C,SunVirial3,StrinatiC,
ZwergerViscosity,Carlson_C_relations,Drut_C,
ValeC_precise,WetterichC,Goulko_UFG_2016,
RossiContact,AlhassidC} respectively.
Similarly, universal relations involving a two-body contact $C_2$ hold for two-component fermions in 2D~\cite{FelixViriel,CombescotC,WernerCastinRelationsArxiv,
Moelmer,BraatenRFshift2D,HofmannAnomaly,
WernerCastinRelationsFermions,MolmerAnyD} and 1D~\cite{FelixViriel,ZwergerRelations1D}, and for single-component bosons in 1D~\cite{FelixViriel,Olshanii_nk} and 2D~\cite{FelixViriel,WernerCastinRelationsBosons}. For a 2D Bose gas, $C_2$ was recently measured interferometrically~\cite{beugnonC}. Several two-body contacts appear in the general relations for single-component fermions with $p$-wave~\cite{UeadeRelationsPwave,ZhangThywissenRelations_pwave,
ZhouRelationsLwave,PengLiuHuRelations3Dpwave,
ThywissenRelationspwave_exp} or higher partial-wave~\cite{ZhangRelations_dwave} short-range interactions in~3D, and dipolar plus short-range interactions in 2D~\cite{HofmannDipolarRelations}.
For bosons in 3D with resonant interactions, in addition to the two-body contact $C_2$, measured in~\cite{jin_contact_BEC,HazibabicC3}, a three-body contact $C_3$ appears in several exact relations~\cite{FelixViriel,WernerCastinRelationsArxiv,
BraatenBosons,CastinWerner_nk_trimer,
WernerCastinRelationsBosons}. This appearance of $C_3$ is linked to the Efimov effect. In particular, a three-body parameter has to be included in the definition of the zero-range model, and $C_3$ is proportional to the derivative of the energy w.r.t. the three-body parameter. $C_3$~was measured interferometrically in~\cite{HazibabicC3} after an interaction quench to unitarity. A three-body contact also appears for single-component fermions with higher-partial-wave resonant interactions, both in 2D~\cite{ZhangRelations2Dpwave_superefi} (where a super-Efimov effect occurs) and in 1D~\cite{NishidaRelations1D,NishidaRelations1D_2}. Two-body and three-body contacts were also found to be useful to describe short-distances or large-momenta properties in clusters of helium atoms~\cite{BazakHeliumContact} and in nuclei~\cite{BarneaNuclContactsPRL2015,
BarneaNuclContactsPRC2015,HenPLB2018,Nucl_C_exp_Nature2020,
HenContactsNatPhys2021,
CLAScollab_PLB2021,GandolfiC3_nucl}, although the corresponding relations are only approximate because the interaction range is not much smaller than the interparticle distance.
Here we show that a three-body contact $C_3$ plays an important role for two-component fermions with resonant interactions in 3D, although there is no three-body Efimov effect so that the zero-range model is parameterized by the scattering length $a_2$ without any three-body parameter~\cite{Efimov,Petrov3fermions}. We work within the zero-range model in Sec.~\ref{sec:ZR}, and we consider models with a small finite interaction range in Sec.~\ref{sec:finite_b}. Within the zero-range model, the number of triplets of particles separated by a small distance (Sec.~\ref{sec:N3}), the tail of the center-of-mass momentum distribution of pairs separated by a small distance (Sec.~\ref{sec:tail_NP}), and the tail of the two-particle momentum distribution (Sec.~\ref{sec:N_k1k2}) are expressed in terms of $C_3$, and $C_3$ is also related to the third order density correlation function (Sec.~\ref{sec:g3}) and to the behavior of the many-body wavefunction when three particles approach each other (Sec.~\ref{sec:psi}). When the interaction range $b$ is non-zero but still small compared to the other typical lengthscales, we consider two additional observables, the formation rate of deeply bound dimers by three-body recombination $\Gamma_3$ (Sec.~\ref{sec:loss}), and the three-body contribution to the energy correction induced by the finite interaction range $\delta\!E_3$ (Sec.~\ref{sec:dE}). We express $\Gamma_3$ and $\delta\!E_3$ in terms of $C_3$, and of a three-body parameter~$a_3$ (which is small in the zero-range regime). We define~$a_3$ through the asymptotic behavior of the three-body zero-energy scattering state at distances $\gg b$ and $\ll |a_2|$.
We consider $N_\UP$ fermions of spin $\UP$ and $N_\down$ fermions of spin $\down$, either confined by a smooth external trapping potential, or in a box with periodic boundary conditions. We consider equal masses for $\UP$ and $\down$ particles, and discuss the unequal-mass case in App.~\ref{app:imbal}. We consider a stationary state throughout the article, and discuss statistical mixtures and non-stationary states in App.~\ref{app:non_stat}. In parallel with presenting the relations involving the three-body contact, we will recall for comparison the known Tan relations involving the two-body contact.
\section{Relations for the zero-range model} \label{sec:ZR}
In this Section we work within the zero-range model, where interactions are characterized by a single parameter, the two-body scattering length $a_2$. The zero-range model is defined in Eqs.~\eqref{eq:schro}, \eqref{eq:BP}. The zero-range limit of finite-range models is expected to be universally described by the zero-range model.\footnote{The absence of $(N_\UP+N_\down)$-body Efimov effect (in the equal-mass case) was shown for $(N_\UP,N_\down)=(2,1)$~\cite{Efimov,Petrov3fermions}, $(2,2)$~\cite{Petrov4body2004,EndoCastin2+2,Michelangeli2plus2,Seiringer2plus2}, $(3,1)$~\cite{CMP}, $(4,1)$~\cite{PetrovBazak5body}, $(5,1)$~\cite{Bazak6body}, and $(N,1)$ for any $N$~\cite{SeiringerNplus1}. The zero-range model was proven to be self-adjoint in the $(2,2)$~\cite{Seiringer2plus2} and $(N,1)$~\cite{SeiringerNplus1} cases. In the latter case, some rigorous results about the Tan relations were also obtained~\cite{SeiringerNplus1}. The convergence of finite-range models towards the zero-range model in the zero-range limit was confirmed by various theoretical studies, see e.g.~\cite{BlumeRevue,LudoYvanBoite} for single-channel models and~\cite{BraatenLong,WernerTarruellCastin,
ZhangLeggettUniv,MoraCastinPricoupenkoCRAS,castin_tignone} for two-channel models, and by numerous theory-experiment comparisons for the many-body problem, e.g.~\cite{Giorgini,GrimmModes,ShinPhaseDiag,SylEOS,
NirEOS,GezerlisXi,CarlsonAFQMC,ValeC_precise,
VanHouckeEOS,KuEOS,RossiEOS,RossiContact,ValeC_T_hom,
ZwierleinC,AlhassidC}.}
\subsection{Number of nearby fermion triplets} \label{sec:N3}
If one measures the positions of all particles, the average number of pairs of particles whose separation is smaller than some $\epsilon$ is given by
\begin{equation}\label{eq:N2}
N_2(\epsilon) \underset{\epsilon\,\to\,0}{\sim} \,C_2\ \frac{\epsilon}{4\pi}
\end{equation}
where $C_2$ is the two-body contact~\cite{TanEnergetics,TanLargeMomentum}. Similarly, let us consider the number of triplets of fermions separated by small distances. For three particles 1,2,3, let us define the hyperradius
\begin{equation}\label{eq:def_R}
R=\sqrt{\frac{2}{3}
\left(r_{12}^{\phantom{aa}2}+r_{13}^{\phantom{aa}2}+r_{23}^{\phantom{aa}2}\right)}
\end{equation}
where $r_{ij}$ is the distance between particles $i$ and $j$. If one measures the positions of all particles, the average number of triplets of particles with hyperradius $R<\epsilon$ is given by
\begin{equation}\label{eq:N3}
\boxed{N_3(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \ C_3\ \epsilon^{2s+2}}
\end{equation}
where the prefactor $C_3$ is what we call the three-body contact, while the exponent $s=1.772724267\ldots$ is the lowest positive solution different from 1 of
\begin{equation}\label{eq:transc_s}
\left(1-s^2\right)\,\sin\left(\frac{s\pi}{2}\right) -\frac{4}{\sqrt{3}}\,s\,\cos\left(\frac{s\pi}{6}\right) + 4\,\sin\left(\frac{s\pi}{6}\right) =0.
\end{equation}
The scaling $N_3(\epsilon) \propto \epsilon^{2s+2}$ was already obtained in~\cite{Petrov4body2004,TanScalingREVTEX} (see Sec.~\ref{sec:psi} for a rederivation). The anomalous exponent $2s+2$ comes from the analytical solution of the unitary three-body problem~\cite{Efimov}, and is directly linked to a hidden dynamical symmetry and a separability of the three-body problem in hyperspherical coordinates~\cite{LeChapitreIn2,CRAS,WernerSym}, or in a field theory point of view, to non-relativistic conformal invariance, with $s+5/2$ the scaling dimension of a three-fermion operator~\cite{ChapNishidaSon,son_cft,Mehen_cft}.
In general there are two different contributions to the three-body contact, coming from $\UP\UP\down$ and~$\UP\down\down$ spin configurations: Denoting by $N_{2,1}$ (resp.~$N_{1,2}$) the contributions to $N_3(\epsilon)$ from triplets of particles of~spins~$\UP\UP\down$ (resp.~$\UP\down\down$), we have
\begin{equation}\label{eq:N21}
\boxed{N_{2,1}(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \,C_{2,1}\ \epsilon^{2s+2}}
\end{equation}
\begin{equation}\label{eq:N12}
\boxed{N_{1,2}(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \,C_{1,2}\ \epsilon^{2s+2}}
\end{equation}
where $C_{2,1}$ (resp.~$C_{1,2}$) is what we call the $\UP\UP\down$ (resp.~$\UP\down\down$) three-body contact. Clearly,
\begin{equation}\label{eq:C_3_21_12}
\boxed{C_3 = C_{2,1}+C_{1,2}.}
\end{equation}
%\vskip 1cm
%\noindent \underline{\it Remarks:}
\begin{remas}\leavevmode
\begin{itemize}
\item Due to the antibunching effect associated to the Pauli exclusion between fermions with identical spins, the contribution to Eq.~\eqref{eq:N2} coming from pairs of particles with identical spins ($\UP\UP$~or~$\down\down$) is negligible in the $\epsilon\to0$ limit (it scales as $\epsilon^5$), and $N_2(\epsilon)$ is dominated by the contribution from pairs of particles with opposite spins ($\UP\down$).
Similarly, the contribution to Eq.~\eqref{eq:N3} coming from triplets of particles with identical spins ($\UP\UP\UP$~or~$\down\down\down$) is negligible in the $\epsilon\to0$ limit (it scales as $\epsilon^{10}$), and $N_3(\epsilon)$ is dominated by the contributions from triplets of particles with non-identical spins, $\UP\UP\down$ or $\UP\down\down$, in agreement with Eq.~\eqref{eq:C_3_21_12}.
\item For comparison, in the non-interacting case, the number of nearby pairs and triplets scales as
\begin{align}
N_2^{(0)}(\epsilon) &\propto \epsilon^3\label{eq:N2_0}\\
N_3^{(0)}(\epsilon) &\propto \epsilon^8\label{eq:N3_0}
\end{align}
(more generally, these scalings also hold with a finite-range interaction that does not diverge too strongly at small distance, so that the wavefunction is bounded).
Eq.~\eqref{eq:N2_0} [resp. Eq.~\eqref{eq:N3_0}] is dominated by the contribution from pairs (resp. triplets) of particles with non-identical spins.
Eq.~\eqref{eq:N3_0} includes the antibunching effect due to the Pauli exclusion between the two identical-spin fermions: The wavefunction vanishes linearly with the distance between these fermions, hence an $\epsilon^2$ suppression factor compared to the completely uncorrelated case of non-interacting distinguishable particles
\begin{equation}\label{eq:N3_dist}
N^{(0)}_{3,\ \rm distinguishable}(\epsilon) \propto \epsilon^6.
\end{equation}
\item Compared to the non-interacting case Eqs.~\eqref{eq:N2_0},\eqref{eq:N3_0}, the exponents in Eqs.~\eqref{eq:N2},\eqref{eq:N3} are reduced, i.e. the probability to find particles near to each other is enhanced. This bunching effect is due to the attractive effect of the resonant zero-range interaction, which causes the wavefunction to diverge: When the distance $r$ between two opposite-spin particles vanishes, $\psi \propto 1/r$, which yields Eq.~\eqref{eq:N2}, while in the limit of vanishing hyperradius $R$ between three particles, $\psi \propto R^{s-2}$, which yields Eq.~\eqref{eq:N3}. Note that since $s<2$, the wavefunction indeed diverges for $R\to0$, and the exponent in Eq.~\eqref{eq:N3} is smaller than in the uncorrelated case Eq.~\eqref{eq:N3_dist}, which means that the bunching effect due to the zero-range interactions overcompensates the antibunching effect due to Pauli exclusion.\footnote{On the other hand, the Pauli exclusion effect is not overcompensated by too much: Due to the repulsive effective three-body potential $s^2/R^2$ (see App.~\ref{app:unit_hyp}), we still have $N_3(\epsilon) / \epsilon^2 \propto \epsilon^{2s} \underset{\epsilon\to0}{\to} 0$, which implies that the three-body loss rate divided by the thermalization rate vanishes in the zero-density limit, a crucial ingredient for the zero-range model to be an accurate description of ultracold-atom experiments (see Sec.~\ref{sec:loss}).}
\end{itemize}
\end{remas}
\subsection{Density correlation functions} \label{sec:g3}
The probability density of finding a spin-$\UP$ particle at~$\rr_1$ and a spin-$\down$ particle at $\rr_2$ is given by the pair correlation function $g_2(\rr_1,\rr_2) = \la \hat{\psi}^\dagger_\UP(\rr_1) \hat{\psi}^\dagger_\down(\rr_2) \hat{\psi}_\down(\rr_2) \hat{\psi}_\UP(\rr_1) \ra = \la \hat{n}_\UP(\rr_1) \, \hat{n}_\down(\rr_2) \ra$. Similarly, the probability density of finding a spin-$\UP$ particle at $\rr_1$, a second spin-$\UP$ particle at $\rr_2$, and a spin-$\down$ particle at $\rr_3$ is given by the triplet correlation function
\begin{align*}
g_{2,1}(\rr_1,\rr_2,\rr_3) &= \left\la \hat{\psi}^\dagger_\UP(\rr_1) \hat{\psi}^\dagger_\UP(\rr_2) \hat{\psi}^\dagger_\down(\rr_3) \hat{\psi}_\down(\rr_3) \hat{\psi}_\UP(\rr_2) \hat{\psi}_\UP(\rr_1) \right\ra
\\
&=\left\la \hat{n}_\UP(\rr_1) \, \hat{n}_\UP(\rr_2) \, \hat{n}_\down(\rr_3) \right\ra
\ -\ \delta^3(\rr_1-\rr_2)\ \left\la \hat{n}_\UP(\rr_1) \,\hat{n}_\down(\rr_3)\right\ra\,.
\end{align*}
The second-order density correlation function has the short-distance asymptotic behavior~\cite{TanEnergetics}
\begin{equation*}
\int d^3\!c\ \, g_2\!\!\left(\cc+\frac{\rr}{2},\cc-\frac{\rr}{2}\right) \ \underset{r\to0}{\sim}\ \, \frac{1}{(4\pi)^2}\ \, \frac{C_2}{r^2}\,.
\end{equation*}
For the third-order density correlation function, we find the short-distance asymptotic behavior
\begin{equation}\label{eq:g21_C3}
\boxed{\int d^3\!C\ d^5\Omega\ \, g_{2,1}(\rr_1,\rr_2,\rr_3)\underset{R\,\to\,0}{\sim}\ C_{2,1}\,R^{2s-4}
\ \frac{32(s+1)}{3\sqrt{3}}}
\end{equation}
where the limit $R\to0$ is taken for fixed $\CC$ and~$\Oo$. Here a change of coordinates is implied between $(\rr_1,\rr_2,\rr_3)$ and $(\CC,R,\Oo)$, with $R$ the hyperradius defined in Eq.~\eqref{eq:def_R}, $\CC = (\rr_1+\rr_2+\rr_3)/3$ the center-of-mass, and $\Oo$ the hyperangles which are five dimensionless coordinates that remain to determine the positions of the three particles, see~Eq.~\eqref{eq:def_Oo}. Similarly, the $\UP\down\down$ triplet correlation function $g_{1,2}(\rr_1,\rr_2,\rr_3) = \la \hat{\psi}^\dagger_\UP(\rr_1) \hat{\psi}^\dagger_\down(\rr_2) \hat{\psi}^\dagger_\down(\rr_3) \hat{\psi}_\down(\rr_3) \hat{\psi}_\down(\rr_2) \hat{\psi}_\UP(\rr_1) \ra$ satisfies
\begin{equation}\label{eq:g12_C3}
\boxed{\int d^3\!C\ d^5\Omega\ \, g_{1,2}(\rr_3, \rr_1,\rr_2) \underset{R\,\to\,0}{\sim}\ C_{1,2}\,R^{2s-4}
\ \frac{32(s+1)}{3\sqrt{3}}.}
\end{equation}
\subsection{Link with the many-body wavefunction}\label{sec:psi}
When three particles approach each other, the many-body wavefunction has a singular asymptotic behavior, with a prefactor related to three-body contact. Let $\psi(\rr_1,\,\ldots,\,\rr_N)$ be the (orbital) many-body wavefunction. Without loss of generality, we can assume $N_\UP\geq2$, and consider that particles $(1,2,3)$ have spins $(\UP,\UP,\down)$. We take this convention throughout the article. This means that $\psi$ is antisymmetric w.r.t. exchange of $\rr_1$ and $\rr_2$ ($\psi$~is also antisymmetric w.r.t. exchange of any other pair of same-spin particles). We take the normalization $\int \dr_1\ldots \dr_N\ |\psi(\rr_1,\,\ldots,\,\rr_N)|^2 = 1$.
Within the zero-range model, the stationary Schr\"odinger equation
\begin{equation}\label{eq:schro}
\sum_{i=1}^N \left[-\frac{\hbar^2}{2m}\Delta_{\rr_i} + U(\rr_i) \right] \psi = E\,\psi
\end{equation}
contains an external trapping potential $U$ but no interaction potential; instead, $\psi$ should satisfy a contact condition in the limit where two particles of opposite spin approach each other:
\\There exists $A$ such that
\begin{equation}\label{eq:BP}
\psi(\rr_1,\,\ldots,\,\rr_N)
\ \underset{r\to0}{=}\
\left(
\frac{1}{r} - \frac{1}{a_2}
\right)
\,A(\cc;\rr_2,\rr_4,\,\ldots,\,\rr_N) + O(r)
\end{equation}
where $a_2$ is the two-body scattering length, $r = \| \rr_1 - \rr_3 \|$ is the distance between the opposite-spin particles 1 and~3, and $\cc = (\rr_1+\rr_3)/2$ is their center-of-mass. The limit $r\to0$ is taken for fixed $\cc$ and fixed positions of the remaining particles $(\rr_2,\rr_4,\,\ldots,\,\rr_N)$. By antisymmetry, a similar contact condition automatically also holds for all other pairs of opposite-spin particles, and Eqs.~(\eqref{eq:schro},\eqref{eq:BP}) are sufficient to define the eigenstates $\psi$ and energies $E$ of the zero-range model.\footnote{
Configurations with a vanishing interparticle distance are implicitly excluded in~\eqref{eq:schro}. In an equivalent alternative formulation, these configurations are included and regularized delta pseudopotential terms are added~\cite{Petrov3fermions,CRAS}.}
When particles 1, 2 and 3 approach each other, the wavefunction of any stationary state has the asymptotic behavior~\cite{Petrov4body2004,TanScalingREVTEX,LeChapitreIn2,Efimov}
\begin{equation}\label{eq:R0}
\boxed{\psi(\rr_1,\,\ldots,\,\rr_N)
\underset{R\,\to\,0}{\sim} R^{s-2}\sum_{{\rmm}=-1}^{+1} \phi_{\rmm}(\Oo)\ B_\md(\CC;\rr_4,\,\ldots,\,\rr_N).}
\end{equation}
Here, as in Eq.~\eqref{eq:g21_C3}, $R$ is the hyperradius of particles ($1$, $2$, $3$) defined in Eq.~\eqref{eq:def_R}, $\CC$ is their center-of-mass, and $\Oo$ denotes their hyperangles. The limit $R\to0$ is taken for fixed $\Oo, \CC,\rr_4,\,\ldots,\,\rr_N$. The unitary hyperangular wavefunctions $\phi_\md(\Oo)$ are such that $R^{s-2}\,\phi_\md(\Oo)$ is a solution of the three-body problem at zero energy and infinite scattering length with total angular momentum quantum numbers $l=1$ and $\md\in\{-1,0,1\}$, see App.~\ref{app:unit_hyp} for more details.\footnote{
There is a similarity between the two-body and three-body short-distance asymptotic behaviors Eqs.~\eqref{eq:BP} and~\eqref{eq:R0}, given that $1/r - 1/a_2$ is a solution of the two-body problem at zero energy.
}$^{,}$\footnote{
An asymptotic behavior similar to~\eqref{eq:R0} holds when any three particles with spins $\UP\UP\down$ approach each other, the functions $B_\md$ corresponding to different triplets of particles being simply related to each other by antisymmtery.}
It is known~\cite{TanEnergetics} that $C_2$ is given by the norm of the function $A$ that appears in Eq.~\eqref{eq:BP},
\begin{equation}\label{eq:C2_AA}
C_2 = (4\pi)^2\, N_\UP\,N_\down
\int \left|A(\cc;\rr_2,\rr_4,\,\ldots,\,\rr_N)\right|^2 \,d^3\!c\,\dr_2\,\dr_4\ldots \dr_N.
\end{equation}
Similarly, the three-body contact is given by the norm of the function $B$ that appears in Eq.~\eqref{eq:R0},
\begin{equation}\label{eq:C3_BB}
\boxed{ C_{2,1} = N_\UP \left(N_\UP-1\right) N_\down
\ \frac{3\sqrt{3}}{32\,(s+1)}\
\\ \sum_{\md=-1}^{+1}
\int \left|B_\md(\CC;\rr_4,\,\ldots,\,\rr_N)\right|^2 d^3\!C\, \dr_4 \ldots \dr_N.}
\end{equation}
The expression~\eqref{eq:N21} for the number of nearby $\UP\UP\down$ fermion triplets, together with Eq.~\eqref{eq:C3_BB}, simply follow from Eq.~\eqref{eq:R0} by integrating $|\psi|^2$ over the $R<\epsilon$ region.\footnote{
Here we used the change of integration variables $(\rr_1,\rr_2,\rr_3)\longrightarrow(\rr,\rrho,\CC)$ with $\rrho$ defined in~\eqref{eq:def_jaco}, of Jacobian $|\partial(\rr_1,\rr_2,\rr_3)/\partial(\rr,\rrho,\CC)| = (\sqrt{3}/2)^3$. We also used the property $\int d^5\Omega\ \phi_\md(\Oo)^*\,\phi_{\md'}(\Oo) = \delta_{\md,\md'}$.} Similarly, the relation involving $g_{2,1}$, Eq.~\eqref{eq:g21_C3}, follows immediately from Eqs.~\eqref{eq:R0},\eqref{eq:C3_BB} and from the expression of $g_{2,1}$ in first quantization, $g_{2,1}(\rr_1,\rr_2,\rr_3) = N_\UP(N_\UP-1)N_\down
\int \dr_4 \ldots \dr_N\,|\psi(\rr_1,\,\ldots,\,\rr_N)|^2$.
There is a completely analogous relation between $C_{1,2}$ and the behavior of the many-body wavefunction when three particles of spins $\UP\down\down$ approach each other (provided $N_\down \geq 2$). Specifically, considering that particle 4 has spin $\down$ (while particles 1, 3 still have spins $\UP, \down$) and denoting by $\tilde{R}$, $\tilde{\Oo}$, $\tilde{\CC}$ the hyperradius, hyperangles and center-of-mass associated to particles $3,4,1$ [obtained by replacing $(\rr_1,\rr_2,\rr_3)$ with $(\rr_3,\rr_4,\rr_1)$ in Eqs.~\eqref{eq:def_jaco}, \eqref{eq:def_RR}, \eqref{eq:def_Oo}] we have
\begin{equation}\label{eq:R0t}
\boxed{\psi(\rr_1,\,\ldots,\,\rr_N)
\underset{\tilde{R}\to0}{\sim}
\tilde{R}^{s-2}\!\!\sum_{{\rmm}=-1}^{+1}\!\! \phi_{\rmm}(\tilde{\Oo})\ \tilde{B}_\md\left(\tilde{\CC};\rr_2,\rr_5,\,\ldots,\,\rr_N\right)}
\end{equation}
which yields the relations for $N_{1,2}$ and $g_{1,2}$, Eqs.~\eqref{eq:N12}, \eqref{eq:g12_C3}, with
\begin{equation}\label{eq:C3_BBt}
\boxed{C_{1,2} = N_\down (N_\down-1) N_\UP
\ \frac{3\sqrt{3}}{32\,(s+1)}
\ \, \sum_{\md=-1}^{+1}
\ \, \int \left|\tilde{B}_\md\left(\CC;\rr_2,\rr_5,\,\ldots,\,\rr_N\right)\right|^2 d^3\!C\, \dr_2\,\dr_5 \ldots \dr_N.}
\end{equation}
Here we assumed $N_\down\geq2$; if $N_\down=1$, then $N_{1,2}$ is obviously zero, and $C_{1,2}=0$.
%\vskip 0.1cm
%\noindent\underline{\it Remark:}
\begin{rema}
Higher-body contacts can be defined in the same way than the three-body contacts. When $j_\UP$ particles of spin $\UP$ and $j_\down$ particles of spin $\down$ approach each other, the $N$-body wavefunction factorizes into the product of {\it (i)} a function of the relative positions of the $j=j_\UP+j_\down$ nearby particles, given by a zero-energy solution of the $j_\UP+j_\down$ body problem at $a_2{=}\infty$, proportional to the $j$ body hyperradius to some power, and {\it (ii)} a function of the center-of-mass of the $j$ nearby particles and of the positions of the $(N{-}j)$ other particles~\cite{TanScalingREVTEX,LeChapitreIn2}. The $L^2$ norm of the latter function defines the $j_\UP+j_\down$ body contact $C_{j_\UP,j_\down}$ (up to a prefactor which is a matter of definition).
\end{rema}
\subsection{Large-momentum tail of the center-of-mass momentum distribution of nearby pairs}\label{sec:tail_NP}
Since $C_2$ and $C_3$ determine short-distance singularities, it is natural that they also determine large-momentum tails. $C_2$ determines the leading tail of the single-particle momentum distribution~\cite{TanLargeMomentum}
\begin{equation}\label{eq:nk_C2}
N_\sigma(\kk) \ \underset{k\,\to\,\infty}{\sim} \ \frac{C_2}{k^4}
\end{equation}
with the normalization $\int N_\sigma(\kk)\, d^3k/(2\pi)^3 = N_\sigma$ (in the case of periodic boundary conditions, momenta become discrete and momentum integrals should be replaced by sums).
$C_3$ also determines a large-momentum tail. Suppose that one measures, for a pair of particles with opposite spin, both their spatial separation $r$ and their center-of-mass momentum $\KK$ (this is allowed since the corresponding operators commute). Let $N_2(\epsilon, \KK)$ be the probability distribution over $\KK$ conditional to $r<\epsilon$, with the normalization $\int N_2(\epsilon,\KK)\, \dK/(2\pi)^3 = N_2(\epsilon)$. In other words, $N_2(\epsilon, \KK)$ is the center-of-mass momentum distribution of the pairs of particles separated by a distance $<\epsilon$, normalized to the total number of such pairs.\footnote{More formally, $N_2(\epsilon, \KK)$ is the expectation value of the operator
\[
\ds \sum_{i:\UP, j:\down} \theta\left(\epsilon-\hat{r}_{ij}\right)\ (2\pi)^3\,\delta^3\left(\hat{\KK}_{ij}-\KK\right),
\]
where $\theta$ is the Heaviside function, while $\hat{r}_{ij} = \| \hat{\rr}_j - \hat{\rr}_i \|$ and $\hat{\KK}_{ij}=\hat{\kk}_i + \hat{\kk}_j$ are the operators corresponding to the relative distance and the center-of-mass momentum of particles $i$ and $j$.} We have
\begin{equation}\label{eq:def_NP}
N_2(\epsilon, \KK) \ \underset{\epsilon\,\to\,0}{\sim}\ \, N_P(\KK)\ \frac{\epsilon}{4\pi}
\end{equation}
where $N_P(\KK)$ is what we call the center-of-mass momentum distribution of nearby fermion pairs. With this definition, one simply has the normalization
\begin{equation}\label{eq:int_NP}
\int N_P(\KK)\,\frac{d^3\!K}{(2\pi)^3} = C_2
\end{equation}
as a consequence of~\eqref{eq:N2}.\footnote{
$N_P(\KK)$ appears naturally in the diagrammatic formalism: In a homogeneous system, $N_P(\KK)$ divided by the volume equals $-(m^2/\hbar^2)\,\Gamma(\KK,\tau=0^-)$ where $\Gamma$ is the pair-propagator defined e.g. in~\cite{RossiContact}. This can be shown using a lattice model~\cite{WernerCastinRelationsFermions}, for which $1/(r\,r')$ in~\eqref{eq:rho2_BP} can be replaced by $\phi(\rr)\phi(\rr')$ with $\phi$ the zero-energy two-body scattering state, and setting $\rr=\rr'=\vn$.
}
The tail of $N_P(\KK)$ is determined by the three-body contact\footnote{The fact that $\Gamma(K,\tau=0^-)$ has a tail $\propto C_3 / K^{2s+4}$ was pointed out to us by Shina Tan (private communication, Aspen, 2011).}:
\begin{equation}\label{eq:NP_tail}
\boxed{
\bar{N}_P(K) \ \underset{K\to\infty}{\sim} \ \ \Mr\ \frac{C_3}{K^{2s+4}}}
\end{equation}
with the prefactor
\begin{equation}\label{eq:Mr}
\boxed{\Mr = 32\, \pi^3\,\frac{4^s}{3^{s+1/2}}\,(s+1)\,\Gamma(s+2)^2\,\sin^2(s\pi)\,\Nr^2}
\end{equation}
whose numerical value is
\begin{equation}\label{eq:Mr_val}
\boxed{\Mr \simeq 2272.}
\end{equation}
Here $\bar{N}_P(K)$ stands for the angular average $\int N_P(\KK)\,d{\bf \hat{K}}/(4\pi)$ (with ${\bf \hat{K}} := \KK / K$, and $d{\bf \hat{K}}$ the differential solid angle, so that $d^3\!K = d{\bf \hat{K}}\ K^2\,d\!K$).
To derive this result, we consider the two-body reduced density matrix $\rho_2(\rr_\UP, \rr_\down; \rr'_\UP, \rr'_\down) := \la \hat{\psi}^\dagger_\UP(\rr'_\UP) \hat{\psi}^\dagger_\down(\rr'_\down) \hat{\psi}_\down(\rr_\down) \hat{\psi}_\UP(\rr_\UP) \ra$ or equivalently in first quantization
\begin{equation}\label{eq:rho2_1q}
\rho_2\left(\rr_\UP, \rr_\down; \rr'_\UP, \rr'_\down\right) = N_\UP N_\down \int \dr_2\,\dr_4 \ldots \dr_N
\
\psi^*\!\left(\rr'_\UP, \rr_2, \rr'_\down,\rr_4,\, \ldots,\, \rr_N\right)
\ \psi\left(\rr_\UP, \rr_2, \rr_\down,\rr_4,\,\ldots,\,\rr_N\right).
\end{equation}
Inserting the two-body contact condition~\eqref{eq:BP} into~\eqref{eq:rho2_1q} yields
\begin{equation}\label{eq:rho2_BP}
\rho_2\left(\cc - \frac{\rr}{2}\,,\, \cc + \frac{\rr}{2} \,;\,
\cc'-\frac{\rr'}{2}\,,\, \cc'+\frac{\rr'}{2} \right)
\ \underset{r,r'\,\to\,0}{\sim}
\ \ \frac{g_P(\cc,\cc')}{(4\pi)^2\, r\, r'}
\end{equation}
where
\begin{equation}\label{eq:P_vs_A}
g_P(\cc,\cc') = (4\pi)^2\ N_\UP N_\down \int \dr_2\, \dr_4 \ldots \dr_N
\ A^*\!\left(\cc'; \rr_2, \rr_4,\, \ldots,\, \rr_N\right)
\ A(\cc; \rr_2, \rr_4, \,\ldots,\, \rr_N)\,.
\end{equation}
We see that $g_P$ can be physically interpreted as a coherence function for pairs of nearby fermions. Accordingly, $N_P$ is related to $g_P$ by Fourier transformation,
\begin{equation}\label{eq:NP_gP}
N_P(\KK) = \int \dc\,\dc'\ \,e^{-i\, \KK \cdot (\cc - \cc')} \ g_P(\cc,\cc')
\end{equation}
as can be formally shown using the definition~\eqref{eq:def_NP} of $N_P$. Hence
\begin{equation}\label{eq:N_P_AA}
N_P(\KK) = (4\pi)^2\,N_\UP\,N_\down\,\int \dr_2\,\dr_4\ldots \dr_N
\\ \left| \int \dc\ e^{-i\, \KK \cdot \cc}\ A(\cc; \rr_2, \rr_4, \,\ldots,\, \rr_N) \right|^2.
\end{equation}
We then follow a reasoning resembling the one used to derive Eq.~\eqref{eq:nk_C2} in~\cite[Sec.~IV.A]{WernerCastinRelationsFermions}.
\noindent In the large $K$ limit, the Fourier transform with respect to $\cc$ in~\eqref{eq:N_P_AA} is dominated by the contributions from the singularities of $A(\cc; \rr_2, \rr_4, \,\ldots,\,\rr_N)$, which occur when $\cc$ approaches one of the $\rr_i$ ($i=2, 4, \,\ldots,\, N$). This corresponds to particles 1, 3 and~$i$ being close to each other [since $A(\cc; \rr_2, \,\ldots)$ determines the wavefunction $\psi$ when particles 1 and 3 are close to $\cc$, according to Eq.~\eqref{eq:BP}]. For example, for $i=2$, the behavior of $A(\cc; \rr_2, \rr_4,\, \ldots,\, \rr_N)$ in the limit $\cc\to\rr_2$ is determined by the asymptotic behavior of $\psi$ when the three particles 1, 2, 3 are close, given by Eq.~\eqref{eq:R0}. Therefore, we just need to take the limit where particles 1 and 3 approach each other in~\eqref{eq:R0} to obtain
\begin{equation*}
A(\rr_2-\uu; \rr_2, \rr_4, \,\ldots,\, \rr_N)
\underset{u\to0}{\sim}\ \
\Nr\ \frac{s}{2}\ \cos\left(\frac{s\pi}{2}\right)
\left(\frac{2\,u}{\sqrt{3}} \right)^{s-1}
\, \sum_{\md=-1}^1
\,Y_1^\md({\bf \hat{u}})
\, B_\md(\rr_2; \rr_4, \,\ldots,\, \rr_N)
\end{equation*}
where we used the expression~(\eqref{eq:phi}, \eqref{eq:varphi_l=1}) of $\phi_\md(\Oo)$. There is a similar singularity of $A(\cc; \rr_2, \rr_4, \,\ldots,\, \rr_N)$ when $\cc$ approaches $\rr_i$ with $4\leq i\leq N$; when particle $i$ has spin $\down$, the function $\tilde{B}_\md$ introduced in~\eqref{eq:R0t} appears instead of $B_\md$. This gives
\begin{multline}\label{eq:TF_A}
\int \dc\ e^{-i\, \KK \cdot \cc}\ A(\cc; \rr_2, \rr_4, \,\ldots,\, \rr_N)
\underset{K\,\to\,\infty}{\sim}
\ \ \frac{2^{s-2}}{3^{(s-1)/2}}
\ \Nr\, \cos\left(\frac{s\pi}{2}\right)
\sum_{\md=-1}^1 I_\md(\KK)\\
\times \left[
\sum_{i:\UP, i\neq1} (1-2\,\delta_{i,2})\,e^{-i\KK\cdot\rr_i}
\, B_\md\left(P_{2i}(\rr_2;\rr_4,\,\ldots,\,\rr_N) \right)
\right.\\
\left. + \sum_{i:\down, i\neq3} (1-2\,\delta_{i,4})\,e^{-i\KK\cdot\rr_i}
\, \tilde{B}_\md\left(P_{4i}(\rr_2;\rr_4,\ldots,\rr_N)\right)
\right]
\end{multline}
where the first (resp. second) sum over $i$ is taken over particles with spin $\UP$ (resp. $\down$), $P_{ji}(\rr_2;\rr_4,\,\ldots,\,\rr_N)$ is obtained from $(\rr_2;\rr_4,\,\ldots,\,\rr_N)$ by exchanging $\rr_j$ with $\rr_i$ ($P_{ii}$ is the identity), and $I_\md(\KK) := s \int d^3\!u\, e^{i\KK\cdot\uu}\, u^{s-1}\, Y_1^\md({\bf \hat{u}})$. The latter integral can be evaluated analytically: Using $\int d{\bf \hat{u}}\ e^{i \KK\cdot\uu}\, Y_1^\md({\bf \hat{u}}) = 4\pi i \ Y_1^\md(\hat{K})\, j_1(K u)$ with $j_1(t) = (\sin t - t \cos t) / t^2$, and evaluating the remaining integral over $u$ by integrating along a closed contour including the positive real axis and negative imaginary axis, we get
\begin{equation}\label{eq:I_K}
I_\md(\KK) = 4\pi i\ \Gamma(s+2)\,\sin\left(\frac{s\pi}{2}\right)\ \frac{Y_1^\md({\bf \hat{K}})}{K^{s+2}}.
\end{equation}
Inserting~\eqref{eq:TF_A}, \eqref{eq:I_K} into~\eqref{eq:N_P_AA}, expanding the modulus squared, and neglecting in the large $K$ limit the cross terms coming from two different values of $i$,\footnote{By power counting, these cross terms give rise to a $\propto 1/K^{2 s_4 + 4}$ tail of $\bar{N}_P(K)$, where $s_4$ is the smallest scaling exponent of the unitary four-fermion problem; this is indeed negligible compared to the leading $1/K^{2 s+{\gre 4}}$ tail of $\bar{N}_P(K)$, given that $s_4 = 2.509(1)$ is larger than $s \equiv s_3$. This value of $s_4$ follows from the four-body ground-state energy $E_4$ in an isotropic harmonic trap computed in~\cite{DailyBlume_4body_spectrum} and the relation~\cite{TanScalingREVTEX,WernerSym} $E_N = (s_N + 5/2) \hbar\omega$.} we obtain the result~\eqref{eq:NP_tail}, \eqref{eq:Mr}, \eqref{eq:Mr_val}, where we used the value of $\Nr$ given in Eq.~\eqref{eq:Nr} of App.~\ref{app:unit_hyp}.
\subsection{Large-momentum tail of the two-particle momentum distribution}\label{sec:N_k1k2}
The three-body contact also determines the asymptotic behavior at large momenta of the two-particle opposite-spin momentum distribution function, defined by
\begin{equation}\label{eq:def_Nk1k2}
N(\kk_1,\kk_2) = \left\langle \hat{N}_\UP(\kk_1)\,\hat{N}_\down(\kk_2) \right\rangle
\end{equation}
where $\hat{N}_\sigma(\kk) := \hat{c}^\dagger_\sigma(\kk)\,\hat{c}_\sigma(\kk)$ with $\hat{c}_\sigma(\kk) = \int \dr \ e^{- i \kk \cdot \rr}\ \hat{\psi}_\sigma(\rr)$ the annihilation operator of a particle of spin $\sigma$ in the state $|\kk\ra$ defined by $\la \rr | \kk \ra = e^{i\kk\cdot\rr}$.
Since $\int \hat{N}_\sigma(\kk) \, \frac{d^3k}{(2\pi)^3} = N_\sigma$, we have the normalization $\int N(\kk_1,\kk_2) \, \frac{d^3k_1}{(2\pi)^3}\, \frac{d^3k_2}{(2\pi)^3} = N_\UP\,N_\down$.
Experimentally, the two-particle momentum distribution can be accessed from the statistics of time-of-flight images, as was recently demonstrated in 2D~\cite{Jochim_nk_noise}.\footnote{In~3D, early measurements in the BEC regime were reported in~\cite{Jin_nk_noise}; see also~\cite{Braun_n_k1_k2} for a recent numerical study. In optical lattices, detailed experimental studies were carried out in recent years using metastable helium atoms~\cite{Clement_HBT,Clement_pairs_opposite_k,Clement_FCS}.}
Taking the limit of a large relative momentum $k\to\infty$, and integrating over the center-of-mass momentum $\KK$, one obtains a tail proportional to $C_2$,
\begin{equation}\label{eq:Nkk_C2}
\int N\left(\frac{{\bf K}}{2}-{\bf k}\,,\,\frac{{\bf K}}{2}+{\bf k}\right)\ \frac{\dK}{(2\pi)^3}\; \underset{k\,\to\,\infty}{\sim}\; \frac{C_2}{k^4}
\end{equation}
as pointed out in~\cite{BarneaFewBodySys}. If we instead send $K$ to infinity, and average over the direction of $\KK$, we obtain a tail proportional to $C_3$,
\begin{equation}
\boxed{\lim_{K\to\infty}\ \lim_{k\,\to\,\infty} \ K^{2s+4}\ k^4\ \frac{1}{4\pi} \int d{\bf \hat{K}}\ \ N\left(\frac{{\bf K}}{2}-{\bf k}\,,\,\frac{{\bf K}}{2}+{\bf k}\right)\ \ = \ \Mr\ C_3}
\end{equation} where the constant $\Mr$ was given in~\eqref{eq:Mr}.
This result immediately follows from~\eqref{eq:NP_tail} and from the relation
\begin{equation}\label{eq:Nkk_NP}
N\left(\frac{{\bf K}}{2}-{\bf k}\,,\,\frac{{\bf K}}{2}+{\bf k}\right) \ \ \underset{k\,\to\,\infty}{\sim}\ \ \frac{N_P(\KK)}{k^4}\,,
\end{equation}
obtained in~\cite{BarneaFewBodySys} and rederived in the sequel.\footnote{Relation~\eqref{eq:Nkk_C2} also follows from~\eqref{eq:Nkk_NP}, given~\eqref{eq:int_NP}.}$^{,}$\footnote{
In Ref.~\cite{BarneaFewBodySys}, $N_P(\KK)$ was defined by
\begin{equation}\label{eq:NP_Atilde}
N_P(\KK) = (4\pi)^2\ \int \ \frac{d^3k_2}{(2\pi)^3}\ \frac{d^3k_4}{(2\pi)^3} \ldots \ \frac{d^3k_N}{(2\pi)^3}\ \left| \tilde{A}(\KK; \,\kk_2, \kk_4, \ \ldots \ ,\kk_N) \right|^2
\end{equation}
with $\tilde{A}$ the Fourier transform of $A$
\[
\tilde{A}(\KK; \kk_2, \kk_4, \ \ldots \ ,\kk_N) \ := \ \int \dc \ \dr_2\ \dr_4\,\ldots \, \dr_N\ e^{-i \left(\KK \cdot \cc + \kk_2 \cdot \rr_2 + \kk_4 \cdot \rr_4 + \ldots + \kk_N \cdot \rr_N\right) }\ A(\cc ; \rr_2, \rr_4, \, \ldots \,, \rr_N)\,,
\]
and~\eqref{eq:Nkk_NP} was deduced from the expression
\[
N(\kk_1,\kk_3) = N_\UP\,N_\down\ \int \ \frac{d^3k_2}{(2\pi)^3}\ \frac{d^3k_4}{(2\pi)^3} \ \ldots \ \frac{d^3k_N}{(2\pi)^3}\ \left| \tilde{\psi}(\kk_1, \, \ldots \,, \kk_N) \right|^2
\]
in terms of the momentum-space wavefunction
\[
\tilde{\psi}(\kk_1, \, \ldots \,, \kk_N) \, := \, \int\,\dr_1 \, \ldots \dr_N\ e^{-i \left(\kk_1 \cdot \rr_1 + \ldots + \kk_N \cdot \rr_N\right)}\,\psi(\rr_1, \, \ldots \,, \rr_N)
\]
which has the asymptotic behavior
\[
\tilde{\psi}\left(\frac{\KK}{2}-\kk, \, \kk_2, \, \frac{\KK}{2}+\kk, \, \kk_4, \, \ldots \, , \kk_N\right) \ \underset{k\,\to\,\infty}{\sim}\ \frac{4\pi}{k^2}\ \tilde{A}(\KK; \kk_2, \kk_4, \, \ldots \, , \kk_N)\,.
\]}
The definition~\eqref{eq:def_Nk1k2} yields
\begin{equation}\label{eq:N2_rho2}
N\left(\frac{{\bf K}}{2}-{\bf k}\,,\,\frac{{\bf K}}{2}+{\bf k}\right)=
\int \dc\, \dc' \,\dr\, \dr' e^{i{\bf K}\cdot({\bf c}'-{\bf c})}e^{i{\bf k}\cdot({\bf r}'-{\bf r})}
\ \rho_2\!\left({\bf c}-\frac{{\bf r}}{2},{\bf c}+\frac{{\bf r}}{2};{\bf c}'-\frac{{\bf r}'}{2},{\bf c}'+\frac{{\bf r}'}{2}\right)
\end{equation}
where $\rho_2$ is the two-body reduced density matrix, which has the diverging behavior~\eqref{eq:rho2_BP} when $r$ and $r'$ tend to zero. This short-distance divergence leads to a $k\to\infty$ tail of the Fourier transform Eq.~\eqref{eq:N2_rho2}, which can be computed by replacing $\rho_2$ with its asymptotic expression~\eqref{eq:rho2_BP}, and using the identity (in the sense of distributions) $\int \dr \, e^{-i\kk\cdot\rr}/r = 4\pi/k^2$. Using~\eqref{eq:NP_gP} then yields~\eqref{eq:Nkk_NP}.
\section{Relations for finite-range models} \label{sec:finite_b}
In this section, we go beyond the zero-range model, and consider interactions of small but non-zero range. We express two observable in terms of the three-body contact: the rate of three-body recombinations towards deeply bound dimers (Sec.~\ref{sec:loss}), and the three-body contribution to the energy difference between finite-range and zero-range models (Sec.~\ref{sec:dE}).
\subsection{Three-body loss rate} \label{sec:loss}
In ultracold atom experiments, three-body losses generically take place, being a manifestation of the fact that the true equilibrium state at such low temperatures is not gaseous (with the exception of polarized hydrogen). In this section we relate the rate of three-body losses to the three-body contact.
\subsubsection{Simple finite-range interaction} \label{subsec:simple}
To describe three-body losses, we need to go beyond the zero-range model. In this subsection we consider a simple model where fermions of different spin interact through a rotationally invariant potential $V_2(r)$, of finite range $b$. \footnote{Typically, $b$ is set by the true range $b_0$ of the potential (i.e. the length such that $V_2(r)$ decays quickly for $r\gg b_0$). More generally, $b = {\rMax}(b_0, |r_e|)$ where $r_e$ is the effective range.} We consider the resonant regime where the two-body scattering length $a_2$ is large,
\begin{equation}\label{eq:a>>b}
|a_2|\gg b.
\end{equation}
In this regime, there are two kinds of two-body bound states:
\begin{itemize}
\item the weakly bound dimer, of binding energy $\approx \hbar^2/(m a_2^{\phantom{2}2})$, which exists for $a_2>0$
\item deeply bound dimer(s), of binding energy ${\bl \gtrsim} \, \hbar^2/(m b^2)$, which exist if the interaction potential is deep enough (as in generic cold atom experiments).
\end{itemize}
We consider a stationary solution of the $N$-body Schr\"odinger equation
\begin{equation}\label{eq:schro_V2}
H\ \psi = E\ \psi
\end{equation}
\begin{equation}\label{eq:H_V2}
H = \sum_{i=1}^N \left[-\frac{\hbar^2}{2m}\Delta_{\rr_i} + U(\rr_i) \right]
\ +\ \sum_{\substack{i:\UP,\,j:\down}} V_2(r_{ij})
\end{equation}
in the zero-range regime
\begin{equation}\label{eq:ZRreg}
1/\ktyp \gg b
\end{equation}
where $1/\ktyp$ is defined as the smallest scale of variation of the stationary wavefunction $\psi$ in the region where all interparticle distances are $\gg b$.\footnote{For example, for the ground state of the homogeneous unpolarized gas, $\ktyp$ is $\sim k_F$ for $a_2<0$ and $\sim{\rMax}(k_F, 1/a_2)$ for $a_2>0$, where $k_F$ is the Fermi momentum; for the ground state of a few particles in an isotropic harmonic trap of frequency $\omega$, $1/\ktyp$ is $\sim a_{\rmho}$ for $a_2<0$ and $\sim {\rMin}(a_{\rmho}, a_2)$ for $a_2>0$, where $a_{\rmho}:=\sqrt{\hbar/(m\omega)}$ is the harmonic oscillator length.}
Let us first consider the case where there are no deeply bound states. For simplicity we also assume that the spectrum is discrete (which is the case in a trapping potential of infinite depth -- e.g. a harmonic trap-- or in a box with periodic boundary conditions). Then, in the zero-range regime~\eqref{eq:ZRreg}, the zero-range model is valid for any stationary state $\psi$, in the sense that standard observables tend to their respective values within the zero-range model. This includes observables such as the energy, as well as $N_2(\epsilon)$ and $N_3(\epsilon)$ provided $\epsilon\gg b$ [to reach the asymptotic regimes of Eqs.~\eqref{eq:N2},\eqref{eq:N3} one also needs $\epsilon \ll 1/\ktyp$].
We turn to the experimentally relevant case where deeply bound dimers exist. These deeply bound dimers can be formed through recombination processes between three atoms. Let us denote by $\Gamma_3$ the number of such events per unit of time. The recombination products (the deeply bound dimer and the third atom) escape from the trapped gas, provided the trapping potential $U(\rr)$ has a finite depth much smaller than the binding energy $(\sim \frac{\hbar^2}{m b^2})$ of deeply bound dimers. In typical experiments, this condition holds, and other loss processes are negligible, which allows one to measure $\Gamma_3$ from the decay of the number of trapped atoms, $\dot{N} = - 3\,\Gamma_3$~\cite{JinLifetime,SalomonCrossover,
PetrovJPhysB,ThomasLosses,LuoL3Narrow}.
In the zero-range regime, this decay is slow compared to the other timescales of the problem, as we will see. A standard way to describe such a slowly decaying state in quantum mechanics is to consider a quasi-stationary Gamow state, i.e., a solution of the Schr\"odinger equation with a complex energy and an outgoing-wave asymptotic behavior~\cite{Gamow,LandauLifschitzMecaQnote,Messiah_Vol1_gamow,
Taylor_Scattering_gamow,BlattWeisskopf,MichelGamowShellModel}. Accordingly, we will consider a solution $\psi$ of \eqref{eq:schro_V2}, \eqref{eq:H_V2} with a complex $E$, and an outgoing-wave asymptotic behavior corresponding to the recombination products (a deeply bound dimer + an atom) flying apart towards large distances.\footnote{An alternative approach would be to use the Lindblad equation. We expect that this would lead to the same result for the loss rate, as was checked for three-body losses for bosons in~\cite{BraatenLindblad}.} For such a quasi-stationary state, in the zero-range regime, standard observables again tend to their respective values within the zero-range model.\footnote{An appropriate normalization of the Gamow state will be given below in Eq.~\eqref{eq:norm}. Similarly, the expectation value of an observable $\hat{O}$ should be defined as $\int_\Rr d^{3N}\!X \ \int_\Rr d^{3N}\!X' \ \, \psi^*\!(\XX)\ \la\XX|\hat{O}|\XX'\ra\ \psi(\XX')$.}$^{,}$\footnote{Within the zero-range model, it is convenient to add steep infinite walls to the trapping potential at the boundary of $\Rr_{\rm trap}$, in order to have truly stationary states, thereby neglecting the exponentially suppressed evaporation process discussed in Fn.~\ref{fn:evap}.} The three-body loss rate $\Gamma_3$, however, is simply zero within the zero-range model. To compute $\Gamma_3$, one thus needs to go beyond the zero-range model. As we will see, one only needs to do so for the three-body problem, in order to define a three-body parameter $a_3$. We then find
\begin{equation}\label{eq:Gamma3}
\boxed{\Gamma_3 \simeq -\frac{\hbar}{m}\ 8 s (s+1) \ C_3\ {\rIm} \, a_3 }
\end{equation}
where $C_3$ can be evaluated within the zero-range model. Furthermore, breaking up $\Gamma_3$ into the sum of the two contributions $\Gamma_{2,1}$ and $\Gamma_{1,2}$ corresponding to $\UP\UP\down$ and $\UP\down\down$ loss processes, we have
\begin{equation}\label{eq:Gamma2,1}
\boxed{\Gamma_{2,1} \simeq -\frac{\hbar}{m}\ 8 s (s+1) \, C_{2,1}\ {\rIm} \, a_3}
\end{equation}
\begin{equation}\label{eq:Gamma1,2}
\boxed{\Gamma_{1,2} \simeq -\frac{\hbar}{m}\ 8 s (s+1) \, C_{1,2}\ {\rIm} \, a_3\,.}
\end{equation}
We expect these relations to be asymptotically exact in the resonant zero-range regime~\eqref{eq:a>>b},\eqref{eq:ZRreg}.
To define the three-body parameter $a_3$, we consider the zero-energy free-space solution of the Schr\"odinger equation~\eqref{eq:schro_V2}, \eqref{eq:H_V2} for three particles of spins $\UP\UP\down$ and angular-momentum quantum numbers $(l=1,\md)$ whose asymptotic behavior has the form
\begin{equation}\label{eq:def_a3}
\boxed{\Psi_\md(\RR)
\simeq
\left(R^s - \frac{a_3}{R^s}
\right)
\frac{1}{R^2}\
\phi_\md(\Oo)}
\end{equation}
in the region $\{ b \ll r_{ij} \ll |a_2|, \ \forall i0$). For evaporation of individual atoms, the trap depth can be defined as $[{\rMin}_{\rr\,\in\,\partial\Rr_{\rm trap}} U(\rr)]-[{\rMin}_{\rr\,\in\,\Rr_{\rm trap}} U(\rr)]$, the evaporation rate is $\Phi(\Sr_{\rm trap})$ with $\Sr_{\rm trap} = \partial((\Rr_{\rm trap})^N) \cap \Rr = \{ \XX\in \Rr \, | \, \exists i,\ \rr_i\in \partial\Rr_{\rm trap}\}$, and this rate is exponentially suppressed because $\JJ(\XX)$ with $\XX \in \Sr_{\rm trap}$ is exponentially suppressed.}, we have
\begin{equation}
\Gamma =
\Phi(\Sr_2) + \Phi(\Sr_3)
\end{equation}
where
\begin{equation*}
\Sr_2 = \partial\Rr_2 \cap \Rr\,,
\qquad
\Sr_3 = \partial\Rr_3 \cap \Rr\,.
\nonumber
\end{equation*}
\begin{figure}[!thbp]
\includegraphics[width=0.9\columnwidth]{loss.pdf}
\caption{Geometric illustration of the three-body loss process. The total decay rate is given by the probability flux exiting from the region~$\Rr$ (grey shaded area). Neglecting the flux through $\Sr_2$ (red straight lines) which corresponds to two-body losses, and the flux through $\Sr_{\rm trap}$ (dashed red lines) which corresponds to evaporation, the dominant contribution is the flux through $\Sr_3$ (green circular arcs) which corresponds to three-body losses. The blue arrows represent the \emph{deep-dimer + atom} outgoing wave, corresponding to a deeply bound dimer and an atom flying apart with a large relative momentum and escaping from the trap (this wave propagates in the region~$\Rr_2$). For the purpose of making this illustrative drawing two-dimensional, we considered $N=3$ particles in one space dimension, and fixed the center-of-mass coordinate to $\CC=\vn$. The positions of the three particles are then determined by the Jacobi coordinates $\rr$ and $\rrho$. The trapping region was simply assumed to be a symmetric interval around the origin.
\label{fig:illust}}
\end{figure}
A visual representation is shown in Fig.~\ref{fig:illust}. On physical grounds, we identify $\Phi(\Sr_2)$ as the two-body loss rate $\Gamma_2$, and $\Phi(\Sr_3)$ as the three-body loss rate $\Gamma_3$\,, which determine the average number of lost atoms per unit of time: $-\dot{N} = 2\,\Gamma_2 + 3\,\Gamma_3$. In the considered regime of small range and large trap depth, we expect that $\Gamma_2$ is given by the relation~\eqref{eq:Gamma2}, and that in the considered case where $a_2$ is real, $\Gamma_2$ is negligible compared to $\Gamma_3$ so that $\Gamma \simeq \Gamma_3$\,.\footnote{Indeed, in the zero-range limit, the reasoning of footnote~\ref{fn:Gamma2} yields the expression~\eqref{eq:Gamma2} for $\Phi(\Sr_2)=\Gamma_2$\,; moreover, in the case where $a_2 \in \bbR$, we expect $\Gamma_2 \ll \Gamma_3$, because there is no mechanism that would generate a non-negligible flux exiting $\Rr$ through $\Sr_2$ (hence entering into $\Rr_2$) and propagating in $\Rr_2$ with an initial wavevector high enough to climb the trapping potential barrier and escape from $(\Rr_{\rm trap})^N$.}
It remains to determine $\Gamma_3 = \Phi(\Sr_3)$. We will use the notation
\begin{equation*}
\Sr^{(ijk)} \ := \ \partial B_{ijk} \, \cap \, \Rr\, = \big\{ \, \XX \in \Rr \, \big| \, R_{ijk} = d_3 \big\}\,.
\end{equation*}
From~\eqref{eq:defR3}, we have $\ds \Sr_3 = \bigcup\nolimits_{i\,<\,j\,:\,\UP,\,k\,:\,\down \ {\rm or} \ i\,<\,j\,:\,\down,\,k\,:\,\UP} \Sr^{(ijk)}$. Hence $\Phi(\Sr_3) = \Gamma_{2,1} + \Gamma_{1,2}$ with
\begin{itemize}
\item $\ds \Gamma_{2,1} = \sum_{i\,<\,j\,:\,\UP, k\,:\,\down} \Phi\big(\Sr^{(ijk)}\big)$ the contribution from $\UP\UP\down$ triplets
\item $\ds \Gamma_{1,2} = \sum_{i\,<\,j:\,\down,\,k\,:\,\UP} \Phi\big(\Sr^{(ijk)}\big)$ the contribution from $\UP\down\down$ triplets.
\end{itemize}
\noindent By antisymmetry, each $\UP\UP\down$ triplet gives the same contribution, so that
\begin{equation}\label{eq:Gamma3_J2}
\Gamma_{2,1} =
\frac{N_\UP \, (N_\UP - 1)\,N_\down}{2}\,
\ \Phi\big(\Sr^{(123)}\big)
\,.
\end{equation}
This can be rewritten as \footnote{To justify this rewriting, we mostly follow the bosonic case treated in~\cite[App.~B]{WernerCastinRelationsBosons} (we will also correct in passing a minor error in an intermediate step in~\cite{WernerCastinRelationsBosons}). We need to evaluate $Q := \int_{\Sr^{(123)}} (\psi^* \, \gr_\XX \psi -
\psi \, \gr_\XX \psi^*) \cdot {\bf d^{3N-1}S}$. The constraint $R \equiv R_{123} > \ddd$ which defines the domain $B_{123}$ does not impose any constraint on $(\rr_4,\,\ldots,\,\rr_N)$. Therefore the differential surface vector ${\bf d^{3N-1}S}$, being normal to $\partial B_{123}$, only has its 9 first coordinates which are non-zero, while its coordinates 10 to $3N$ are vanishing. More precisely, denoting by ${\bf d^{3N-1}S}_t$ the 9-dimensional vector whose coordinates are equal to the 9 first coordinates of ${\bf d^{3N-1}S}$, we have ${\bf d^{3N-1}S}_t = {\bf d^{8}S}\ \dr_4 \ldots \dr_N$ where ${\bf d^{8}S}$ is the differential surface vector of the domain $\{(\rr_1,\rr_2,\rr_3) \ / \ R>\ddd \}$. Hence $Q = \int_{\Sr^{(123)}} \,\dr_4 \ldots \dr_N\ {\bf d^{8}S} \cdot (\psi^* \, \gr_{\tilde{\XX}} \psi -\psi \, \gr_{\tilde{\XX}} \psi^*)$ where ${\tilde{\XX}} := (\rr_1,\rr_2,\rr_3)$. Applying the divergence theorem backwards, this can be rewritten $Q = \int_{\Sr^{(123)}} \,\dr_4 \ldots \dr_N\ d^9\tilde{X}\ \ \gr_{\tilde{\XX}} \cdot (
\psi^* \, \gr_{\tilde{\XX}} \psi -
\psi \, \gr_{\tilde{\XX}} \psi^*
)$. We then perform the change of variables $\ \tilde{\XX} \longrightarrow (\CC,\RR)\ $, and rewrite the integrand as $ \psi^* \Delta_{\tilde{\XX}} \psi - \psi \, \Delta_{\tilde{\XX}} \psi^* =
\frac{1}{3} (\psi^* \Delta_{\CC} \psi - \psi \, \Delta_{\CC} \psi^*) + 2 (\psi^* \Delta_{\RR} \psi - \psi \, \Delta_{\RR} \psi^*) =
\frac{1}{3}\ \gr_\CC \cdot
(\psi^* \, \gr_{{\CC}} \psi -
\psi \, \gr_{{\CC}} \psi^*) + 2 \ \gr_\RR \cdot (
\psi^* \, \gr_{{\RR}} \psi -
\psi \, \gr_{{\RR}} \psi^*)$. By the divergence theorem, the integral over $\CC$ of the $\gr_\CC$ term vanishes, while the integral over $\RR$ of the $\gr_\RR$ term yields the result $Q = -\frac{3\sqrt{3}}{4}\ \ddd^5 \ \int_{\Sr^{(123)}} \dC\ \dr_4 \, \ldots \, \dr_N
\ d^5\Omega\
(\psi^* \, \partial_R
\psi -
\psi \, \partial_R
\psi^*)$\,.}
\begin{equation}\label{eq:Gamma3_partial}
\Gamma_{2,1} = - N_\UP \, (N_\UP - 1)\,N_\down\,
\frac{3\sqrt{3}}{8}\ \ddd^5
\int_{\Sr^{(123)}}\ \dC\ \dr_4 \, \ldots \, \dr_N \ d^5\Omega
\ \ \frac{\hbar}{m} \ {\rIm} \left[\psi^* \partial_R \psi \right]\,.
\end{equation}
We then simplify this expression in two steps:
\begin{enumerate}\romanenumi
\item\label{(i)} replace $\psi$ with its asymptotic behavior~\eqref{eq:R0_b}, with $B_\md$ evaluated within the zero-range model,
\item\label{(ii)} replace the integration domain $\Sr^{(123)}$ by the entire region $\partial B_{123}$.
\end{enumerate}
Step~\eqref{(i)} is justified since the conditions~\eqref{eq:R1/2$, we get $|\delta\!E_3/ \delta\!E_2| \ll 1$ in the zero-range regime, i.e. the~three-body correction to the energy is of higher order than the two-body correction.\footnote{\bl In the mass-imbalanced case discussed in App.~\ref{app:imbal}, the situation is reversed beyond a critical mass ratio, as was already noted in~\cite{GaoEndoCastin}.}$^,$\footnote{Apart from the leading-order term~\eqref{eq:dE2} which is of order $b$, there are also higher-order contributions to $\delta\!E_2$, which we expect to contain no term of order $b^{2s}$ (since the three-body physics does not enter in $\delta\!E_2$) but only integer powers of~$b$. Accordingly, the $\propto b^{2s}$ contribution to $\delta\!E$ should be entirely given by~\eqref{eq:dE_21_12} or~\eqref{eq:dE3}.}
\item Let us assume that~\eqref{eq:dE3} can be analytically continued to complex $a_3$, with $C_3$ still evaluated within the zero-range model. We then recover the expression~\eqref{eq:Gamma3} of the three-body loss rate, simply by substituting ${\rIm}\, E = {\rIm}\ \delta\!E_3$ into~\eqref{eq:ImE}, and using the fact that {\bl $\Gamma=\Gamma_3$ and} $C_3\in\bbR$. Similarly, applying this procedure to~\eqref{eq:dE_21_12} yields $\Gamma_3$ in agreement with the sum of Eqs.~(\eqref{eq:Gamma2,1_real}, \eqref{eq:Gamma1,2_real}).
\end{itemize}
\end{remas}
To derive~\eqref{eq:dE_21_12}, we consider a stationary state $\psi^{(0)}$ and the associated eigenenergy $E^{(0)}$ of the zero-range model, and the stationary state $\psi^{(1)}$ of the finite-range model whose energy $E^{(1)}$ is close to $E^{(0)}$ in the zero-range regime. We are interested in $\delta\!E = E^{(1)} - E^{(0)}$.
Using the shorthand notations
\begin{align*}
V(\rr_1,\,\ldots,\,\rr_N) \ &{:=} \
\sum_{\substack{i\,:\,\UP,\,j\,:\,\down}} V_2\left(\rr_{ij}\right)
\ +\ \sum_{\substack{i\,<\,j\,:\UP,\,k\,:\,\down}} V_{2,1}\left(\RR_{ijk}\right)
\ +\ \sum_{\substack{i\,<\,j\,:\,\down,\,k\,:\,\UP}} V_{1,2}\left(\RR_{ijk}\right)
\\
H_0 \ &{:=} \ \sum_{i=1}^N \left[-\frac{\hbar^2}{2m}\ \Delta_{\rr_i} + U(\rr_i) \right]
\end{align*}
we have $(H_0 + V) \psi^{(1)} = E^{(1)} \psi^{(1)}$, while $\psi^{(0)}$ satisfies $H_0 \psi^{(0)} = E^{(0)} \psi^{(0)}$ together with the two-body contact condition~\eqref{eq:BP}. Let us consider
\begin{equation*}
\Delta \ := \ \left\la \psi^{(0)}, (H_0+V)\,\psi^{(1)} \right\ra
\ -\
\left\la H_0\,\psi^{(0)}, \psi^{(1)} \right\ra.
\end{equation*}
We have $\Delta \ =\ \delta\!E\ \la \psi^{(0)} | \psi^{(1)} \ra$, and hence $\delta\!E / \Delta\to 1$ in the zero-range limit.
To evaluate $\Delta$, we write it as $\Delta = \Delta_0 + \ \la \psi^{(0)} | V | \psi^{(1)} \ra$ where $\Delta_0 \ := \ \la \psi^{(0)}, H_0\,\psi^{(1)} \ra -
\la H_0\,\psi^{(0)}, \psi^{(1)} \ra$. We have $\Delta_0=0$, assuming for simplicity that the interaction potentials ($V_2$, $V_{2,1}$ and $V_{1,2}$) are finite everywhere ({\it i.e.} excluding hard-wall potentials) and do not diverge too quickly at short distances.\footnote{\bl Specifically, by using the divergence theorem as in Sec.~\ref{app:imbal,subsec:dE}, we find the following sufficient conditions: $R^{2+s} \psi^{(1)}$ and $R^{3+s} \partial \psi^{(1)} / \partial R$ tend to zero when $R\to0$, and similarly, $\tilde{R}^{2+s} \psi^{(1)}$ and $\tilde{R}^{3+s} \partial \psi^{(1)} / \partial \tilde{R}$ tend to zero when $\tilde{R}\to0$.} Hence
\begin{equation}\label{eq:DeltaV}
\Delta
\ \simeq \ \int \dr_1 \ldots \dr_N\ \left(\psi^{(0) *} \ V\ \psi^{(1)} \right)\!(\rr_1,\,\ldots,\,\rr_N).
\end{equation}
Since the potentials are short-ranged, the integral is dominated by the contributions from three regions, outside of which $V$ becomes negligible:
\begin{itemize}
\item two particles of spins $\UP\down$ are nearby (at distance $\lesssim b_2$) while all other interparticle distances are $\gg b$
\item there is a triplet of particles of spins $\UP\UP\down$ which are nearby (their hyperradius is $\lesssim b$) while all interparticle distances other than the ones within that triplet are $\gg b$
\item there is a triplet of particles of spins $\UP\down\down$ which are nearby (their hyperradius is $\lesssim b$) while all interparticle distances other than the ones within that triplet are $\gg b$.
\end{itemize}
Denoting the contributions of these regions by $\Delta_2\,$, $\Delta_{2,1}$ and $\Delta_{1,2}$ respectively, we thus have $\Delta \simeq \Delta_2 + \Delta_{2,1} + \Delta_{1,2}$. The two-nearby-particle contribution $\Delta_2$ is given by the r.h.s. of~\eqref{eq:dE2}, as shown in App.~\ref{app:dE2}, in agreement with~\cite{WernerCastinRelationsFermions}. It remains to evaluate the contribution $\Delta_{2,1}$ coming from three nearby particles of spins $\UP\UP\down$. We introduce a length $d_3$ satisfying~\eqref{eq:<0$ is a normalization constant such that $(\phi|\phi)=1$. Note that $\Nr$ does not depend on $\md$, as follows from the relation $L_{\pm}\ \tilde{\phi}_\md = \hbar \sqrt{ l (l+1) - \md (\md+1) }\ \tilde{\phi}_{\md\pm1}$ where $L_{\pm} := L_x \pm i L_y$ and $\tilde{\phi}_\md := (1-\hat{P}) \ \frac{\varphi(\alpha)}{\sin(2\alpha)}\ Y_l^\md(\hat{\bm \rho})$.
\\For $l=1$,
\begin{equation}\label{eq:varphi_l=1}
\varphi(\alpha) = -s\,\cos\!\left[s\left(\frac{\pi}{2}-\alpha\right)\right] + \sin\!\left[s\left(\frac{\pi}{2}-\alpha\right)\right]\, \tan \alpha
\end{equation}
(for arbitrary $l$, see, e.g.,~\cite{LeChapitreIn2} and refs. therein). We refer to the $\{ \phi_\nu(\Oo) \}$ as unitary hyperangular wavefunctions.\footnote{Like the standard hyperspherical harmonics, the $\phi_\nu$ are eigenstates of the Laplacian on the hypersphere, Eq.~\eqref{eq:eigen_hyperang}; but they also satisfy the unitary-limit contact condition Eq.~\eqref{eq:cc_ang} (together with fermionic antisymmetry) which leads to non-integer eigenvalues $s^2$.} We use the shorthand notation $\phi_\md := \phi_{(l=1,\md,n=0)}$.
\end{itemize}
The value of the normalization constant, for $l=1$, is
\begin{equation}\label{eq:Nr}
\Nr = 0.31149\ldots
\end{equation}
We computed this value as follows. We have $\Nr=\sqrt{3/J(s)}$ with $J(s)=\sum_{\md=-1}^1(\tilde{\phi}_{\md}|\tilde{\phi}_{\md})$.
We evaluate the sum over $\md$ thanks to the identity: $\sum_{\md=-1}^{1}Y_1^{\md}({\bf \hat{u}})^*\,Y_1^{\md}({\bf \hat{u}'})=\frac{3}{4\,\pi}P_1({\bf \hat{u}}\cdot{\bf \hat{u}'})$ for any unit vectors ${\bf \hat{u}}$ and ${\bf \hat{u}'}$. This yields
\begin{equation}\label{eq:Js}
J(s)=\frac{3}{4\pi}\,
\int d^5\Omega\,\left[
\left(\frac{\varphi(\alpha)}{\sin(2\alpha)}\right)^2 +
\left(\frac{\varphi(\alpha')}{\sin(2\alpha')}\right)^2 -2\,P_1\!\left(\hat{\bm \rho}\cdot{\bf \hat{\bm \rho}'}
\right)\frac{\varphi(\alpha)}{\sin(2\alpha)}\,\frac{\varphi(\alpha')}{\sin(2\alpha')}
\right]
\end{equation}
where $\alpha'$ and ${\bf \hat{\bm \rho}'}$ are obtained from $\alpha$ and $\hat{\bm \rho}$ by permutation of particles 1 and~2. The integrand in~\eqref{eq:Js} only depends on $\alpha$ and $\alpha'$, since $\hat{\bm \rho}\cdot\hat{\bm \rho}'=\frac{\sin^2(\alpha)+\sin^2(\alpha')-5/4}{\cos(\alpha)\,\cos(\alpha')}$. Moreover $\alpha'$ is a function of $\alpha$ and of $u=\hat{\bf r}\cdot\hat{\bm \rho}$, since $u=\frac{2}{\sqrt{3}}\times\frac{3/4-\sin^2(\alpha)/2-\sin^2(\alpha')}{\cos(\alpha)\sin(\alpha)}$. Hence the integrand in~\eqref{eq:Js} is a function of $\alpha$ and $u$. To evaluate the integral we use the formula $\int d^5\Omega\ f(\Oo) =
\int_0^{\pi/2}d\!\alpha\ \sin^2(\alpha) \, \cos^2(\alpha)
\int d\!\hat{\bf r} \, d\!\hat{\bm \rho}\, f(\alpha,\hat{\bf r},\hat{\bm \rho})$. Since the integrand is independent of $\hat{\bf r}$ and of the azimuthal angle of $\hat{\bm \rho}$ w.r.t. $\hat{\bf r}$, we can integrate over them, which just gives a factor $4\pi\times2\pi$. We are left with the integral over $\alpha$ and $u$. The change of variable $u \longrightarrow \alpha'$ then yields
\begin{multline}\label{eq:Js2}
J(s)=8\sqrt{3}\pi\int_{\Dr}d\!\alpha\,d\!\alpha'\,
\cos(\alpha)\sin(\alpha)\cos(\alpha')\sin(\alpha')\\
\times\left[
\left(\frac{\varphi(\alpha)}{\sin(2\alpha)}\right)^2 +
\left(\frac{\varphi(\alpha')}{\sin(2\alpha')}\right)^2 -2\,P_1\left(
\frac{\sin(\alpha)^2+\sin(\alpha')^2-\frac{5}{4}}{\cos(\alpha)\cos(\alpha')}
\right)\frac{\varphi(\alpha)}{\sin(2\alpha)}\frac{\varphi(\alpha')}{\sin(2\alpha')}
\right]
\end{multline}
where $\Dr$ is the domain $\{
(\alpha,\alpha') \, |\,\frac{\pi}{3}\,\leq\alpha+\alpha'\leq\,\frac{2\pi}{3},\, |\alpha-\alpha'|\leq\,\frac{\pi}{3}\}$. We evaluate analytically the integrals over $\alpha$ and $\alpha'$ for the first and second term of~\eqref{eq:Js2}; for the third term, we evaluate analytically the integral over $\alpha'$, and perform numerically the remaining integration over $\alpha$.
\section{Homogeneous gas} \label{app:C3m}
Let us show that for the homogeneous gas at equilibrium, $C^{(\md)}_{2,1}$ is independent of $\md$. We will use the following lemma: Let $\hat{\Or}$ be an operator acting on functions of $\Oo$, such that
\begin{equation*}
(\hat{\Or}\,\phi_\md)(\Oo) = \lambda_\md\,\phi_\md(\Oo).
\end{equation*}
Then,
\begin{equation}
\frac{N_\UP\,(N_\UP-1)\,N_\down}{2}\ \left\la
\hat{\Or}
\ \theta(\epsilon-R) \right\ra \underset{\epsilon\to0}{\sim}\ \epsilon^{2s+2} \sum_{\md=-1}^1 \lambda_\md\ C^{(\md)}_{2,1}
\end{equation}
with $\theta$ the Heaviside function. This follows from Eq.~\eqref{eq:R0} and the fact that the $\phi_\md(\Oo)$ are orthonormal.
We then consider $Q := \langle \hat{L}_z\ \,\theta(\epsilon-R) \rangle$. In the absence of time-reversal symmetry breaking, we have $Q=0$. On the other hand, applying the above lemma to $\hat{\Or}=\hat{L}_z$ gives $Q \propto C^{(1)}_{2,1} - C^{(-1)}_{2,1}$. Hence $C^{(1)}_{2,1} = C^{(-1)}_{2,1}$.
Since the system is isotropic, the quantity $\langle \hat{L}_i^{\phantom{i}2}\,\theta(R-\epsilon) \rangle$ is the same for $i=x, y$ and $z$; hence $\langle \hat{L}_z^{\phantom{i}2}\,\theta(R-\epsilon) \rangle = \langle \hat{\bf L}^{2}\,\theta(R-\epsilon) \rangle / 3$. Applying the lemma to both sides then yields $C^{(1)}_{2,1} + C^{(-1)}_{2,1} = 2 \sum_{\md=-1}^1 C^{(\md)}_{2,1} / 3$, which gives $C^{(1)}_{2,1} = C^{(0)}_{2,1}$.
Finally we note that for the unpolarized gas ($N_\UP = N_\down$) at equilibrium, the two states $\UP$ and $\down$ play a symmetric role, so that $C^{(\md)}_{2,1} = C^{(\md)}_{1,2}$.
\section{Two-body contribution to the energy correction} \label{app:dE2}
In this Appendix we show that $\Delta_2$ is given by the r.h.s. of~\eqref{eq:dE2}, to leading order in the zero-range regime. In the mass-balanced case, to leading order, we have $\delta\!E \simeq \Delta_2$, so that we recover the result of~\cite{WernerCastinRelationsFermions}. We provide the present derivation because it is not identical (although similar) to the one in Ref.~\cite{WernerCastinRelationsFermions} (the quantity $\Delta_2$, as defined here, does not explicitly appear in~\cite{WernerCastinRelationsFermions}).
Let us introduce a length $d$ such that
\begin{equation}\label{eq:range_d2}
b \ll d \ll 1/\ktyp\,.
\end{equation}
$\Delta_2$ is given by the contribution to the integral~\eqref{eq:DeltaV} coming from the region of $(\rr_1,\,\ldots,\,\rr_N)$ such that there is a pair of particles $(i,j)$ of spins $(\UP,\down)$ separated by a distance $r_{ij} < d$ while all interparticle distances are $\gg d$ [the result will not depend on the value of $d$ within the range~\eqref{eq:range_d2}].
\\Let us denote by $\Rr_{13}$ the region where these conditions are met for the pair $(1,3)$,
\begin{equation*}
\Rr_{13} := \left\{ (\rr_1,\,\ldots,\,\rr_N) \ {\rm such\ that}\ \ r\equiv r_{13} m_\down$ for definiteness. Experimentally, this is realized in a mixture of two fermionic species, such as $^{40}$K-$^6$Li~\cite{GrimmUltrafast}, $^{161}$Dy-$^{40}$K~\cite{GrimmDyK}, or $^{53}$Cr-$^6$Li~\cite{ZaccantiLiCrFeshbach}. In Sec.~\ref{app:sec:gen}, we extend the relations obtained in the main text to the mass-imbalanced case. In Sec.~\ref{app:sec:simpl} we show that, when the mass ratio exceeds a critical value, a~conceptual simplification takes place: A generalized zero-range model can be introduced, within which relations involving the three-body parameters can be formulated and derived more directly. The results of this Appendix are valid if $m_\UP/m_\down$ is smaller than the threshold where the five-body Efimov effect appears~\cite{PetrovBazak5body}, which implies that the four-body Efimov effect~\cite{CMP} and three-body Efimov effect~\cite{Efimov73} do not take place either; moreover $m_\UP/m_\down$ should not be too close to the four-body and five-body Efimov thresholds, as discussed in Sec.~\ref{sec:max_ratio}.
Obviously, in the $N$-body Schr\"odinger equation, the mass $m$ is replaced by the mass $m_i$ of particle $i$ ($m_i = m_\UP$ or $m_\down$ depending on the spin of particle $i$):
\begin{equation}\label{eq:schro_imbal}
\sum_{i=1}^N \left[-\frac{\hbar^2}{2m_i}\Delta_{\rr_i} + U(\rr_i) \right] \psi = E\,\psi
\end{equation}
for the zero-range model, and
\begin{equation}\label{eq:schro_V2_imbal}
\sum_{i=1}^N \left[-\frac{\hbar^2}{2m_i}\Delta_{\rr_i} + U(\rr_i) \right] \psi
\ +\ \sum_{\substack{i:\UP,\,j:\down}} V_2(r_{ij}) \ \psi = E\,\psi,
\end{equation}
for the finite-range model.
It will prove convenient to introduce the angles $\varphi$ and $\tilde{\varphi}$ related to the mass ratio by
\begin{equation*}
{\rm sin}\ \varphi \ = \ \frac{m_\UP}{m_\UP + m_\down}\ ,
\qquad \qquad {\rm sin}\ \tilde{\varphi} \ = \ \frac{m_\down}{m_\UP + m_\down}.
\end{equation*}
The definitions of the center-of-mass and Jacobi coordinates should be generalized to the unequal-mass case. The center-of-mass of particles 1 and 3, which appears in the two-body contact condition~\eqref{eq:BP}, is now $\cc = (m_\UP \rr_1 + m_\down \rr_3)/(m_\UP+m_\down)$. The center-of-mass of particles 1,2,3 (assumed to have spins $\UP,\UP,\down$) is now $\CC = [m_\UP (\rr_1 + \rr_2) + m_\down \rr_3]/(2 m_\UP + m_\down)$, and the Jacobi coordinate $\rrho$ is now $(\rr_2 - \cc) / {\rm cos}\, \varphi$, while we still have $\rr := \rr_3 - \rr_1$ and $\RR := (\rr,\rrho)$. Similarly, the center-of-mass of particles 1,3,4 (of spins $\UP,\down,\down$) is now $\tilde{\CC} = [m_\UP \rr_1 + m_\down (\rr_3+\rr_4)]/(m_\UP + 2 m_\down)$, and the Jacobi coordinate $\tilde{\rrho}$ is now $(\rr_4 - \cc) / {\rm cos}\, \tilde{\varphi}\ $, while we still have $\tilde{\rr} := \rr_1 - \rr_2$ and $\tilde{\RR} := (\tilde{\rr},\tilde{\rrho})$. The three-body Schr\"odinger equation is then still given by~\eqref{eq:schro_3b} provided we define $m$ as twice the reduced mass,
\begin{equation*}
m \ :=\ 2\,\frac{m_\UP m_\down}{m_\UP+m_\down} \,.
\end{equation*}
The continuity equation and the probability current are still given by \eqref{eq:continu},\eqref{eq:current} provided we define
\begin{equation}\label{eq:Ximbal}
\XX \ :=\ \left(\sqrt{\frac{m_1}{m}}\ \rr_1, \ \ldots \ , \, \sqrt{\frac{m_N}{m}}\ \rr_N\right)\,.
\end{equation}
While in the mass-balanced case, there was a single exponent $s$ associated to the unitary three-body problem, in the mass-imbalanced case we have two exponents $s$ and $\tilde{s}$, associated respectively to the $\UP\UP\down$ and $\UP\down\down$ unitary three-body problems~\cite{Efimov73,Petrov3fermions}. As a function of $m_\uparrow/m_\downarrow$, $s$ is continuously decreasing and vanishes at the Efimov-effect threshold $m_\uparrow/m_\downarrow = 13.6069 \ldots\ $, while $\tilde{s}$ is continuously increasing and tends to 2 for $m_\uparrow/m_\downarrow \to \infty$.
\subsection{Extension of the relations from the mass-balanced case} \label{app:sec:gen}
The number of nearby $\UP\UP\down$ triplets is still given by
\begin{equation}\label{eq:N21_imb}
\boxed{N_{2,1}(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \,C_{2,1}\ \epsilon^{2s+2}}
\end{equation}
whereas the number of nearby $\UP\down\down$ triplets is now given by
\begin{equation}\label{eq:N12_imb}
\boxed{N_{1,2}(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \,C_{1,2}\ \epsilon^{2\tilde{s}+2}\,.}
\end{equation}
For our convention $m_\UP > m_\down$, we have $s < \tilde{s}$. Hence the total number of nearby triplets $N_3(\epsilon)$ is dominated by the $\UP\UP\down$ contribution, so that
\begin{equation}\label{eq:N3_imb}
\boxed{N_3(\epsilon)
\underset{\epsilon\,\to\,0}{\sim} \ C_{2,1}\ \epsilon^{2s+2}\,.}
\end{equation}
The remarks at the end of Sec.~\ref{sec:N3} remain valid. In particular, the bunching effect due to the zero-range interactions still overcompensates the antibunching effect due to Pauli exclusion, both for $N_{2,1}(\epsilon)$ and $N_{1,2}(\epsilon)$, because $s < \tilde{s} < 2$.
When three particles of spins $\UP, \UP, \down$ approach each other, the asymptotic behavior of the many-body wavefunction is still given by~\eqref{eq:R0}. This yields [using the Jacobian $|\partial(\rr_1,\rr_2,\rr_3) / \partial(\CC,\RR) | = {\rm cos}^3 \varphi\ $]
\begin{equation*}
\boxed{ C_{2,1} = \sum_{\md=-1}^{+1} C_{2,1}^{(\md)} }
\end{equation*}
with
\begin{equation}\label{eq:C21_BB_imbal}
\boxed{ C_{2,1}^{(\md)} = N_\UP (N_\UP-1) N_\down
\ \frac{\cos^3\!\varphi}{4\,(s+1)}\
\int \left|B_\md(\CC;\rr_4,\,\ldots,\,\rr_N)\right|^2 d^3\!C\, \dr_4 \ldots \dr_N.}
\end{equation}
When three particles of spins $\UP, \down, \down$ approach each other, the asymptotic behavior of the many-body wavefunction is given by~\eqref{eq:R0t} with $s$ replaced by $\tilde{s}$, i.e.
\begin{equation}\label{eq:R0t_imbal}
\psi(\rr_1,\,\ldots,\,\rr_N)
\underset{\tilde{R}\,\to\,0}{\sim}
\tilde{R}^{\tilde{s}-2}\!\!\sum_{{\rmm}=-1}^{+1}\!\! \phi_{\rmm}(\tilde{\Oo})\ \tilde{B}_\md\left(\tilde{\CC};\rr_2,\rr_5,\,\ldots,\,\rr_N\right)\,.
\end{equation}
This yields
\begin{equation*}
\boxed{ C_{1,2} = \sum_{\md=-1}^{+1} C_{1,2}^{(\md)} }
\end{equation*}
with
\begin{equation*}
\boxed{ C_{1,2}^{(\md)} = N_\down (N_\down-1) N_\UP
\ \frac{\cos^3\!\tilde{\varphi}}{4\,(\tilde{s}+1)}\
\int \left|\tilde{B}_\md(\tilde{\CC}; \rr_2, \rr_5,\,\ldots,\,\rr_N)\right|^2 d^3\!\tilde{C}\, \dr_2\,\dr_5 \ldots \dr_N.}
\end{equation*}
The triplet correlation functions satisfy
\begin{equation}\label{eq:g21_C3_imb}
\boxed{\int d^3\!C\ d^5\Omega\ \, g_{2,1}(\rr_1,\rr_2,\rr_3)\underset{R\,\to\,0}{\sim}\ C_{2,1}\,R^{2s-4}
\ \frac{4(s+1)}{\cos^3\!\varphi}}
\end{equation}
\begin{equation}\label{eq:g12_C3_imb}
\boxed{\int d^3\!\tilde{C}\ d^5\tilde{\Omega}\ \, g_{1,2}(\rr_1,\rr_3,\rr_4)\underset{\tilde{R}\,\to\,0}{\sim}\ C_{1,2}\,\tilde{R}^{2\tilde{s}-4}
\ \frac{4(\tilde{s}+1)}{\cos^3\!\tilde{\varphi}}}
\end{equation}
The leading large-momentum tail of $\bar{N}_P(K)$ is given by
\begin{equation}\label{eq:NP_tail_imb}
\boxed{ \bar{N}_P(K) \ \underset{K\to\infty}{\sim} \ \Mr\ \frac{C_{2,1}}{K^{2s+4}} }
\end{equation}
where $\Mr$ is still given by the expression~\eqref{eq:Mr} in terms of $\Nr$, and $\Nr$ now stands for the normalization constant of the unitary hyperangular wavefunction of the $\UP\UP\down$ problem.\footnote{The contribution from the $\UP\down\down$ three-body problem gives rise to a higher-order subleading tail of $\bar{N}_P(K)$, given by $\ds C_{1,2}\ K^{-2\tilde{s}-4}
\times 32\, \pi^3\,4^{\tilde{s}}\ 3^{-\tilde{s}-1/2}\,(\tilde{s}+1)\,\Gamma(\tilde{s}+2)^2\,\sin^2(\tilde{s}\pi)\,\tilde{\Nr}^2$ with $\tilde{\Nr}$ the normalization constant of the unitary hyperangular wavefunction of the $\UP\down\down$ problem.}
The two-particle momentum distribution has the tail
\begin{equation}\label{eq:N_k1k3_imb}
\boxed{\lim_{K\,\to\,\infty}\ \lim_{k\,\to\,\infty} \ K^{2s+4}\ k^4\ \frac{1}{4\pi} \int d\hat{\bf K}\ \ N\left(\frac{m_\UP}{m_\UP+m_\down}\,\KK-{\bf k}\,,\,\frac{m_\down}{m_\UP+m_\down}\,\KK+{\bf k}\right)\ \ = \ \Mr\ C_{2,1}\,.}
\end{equation}
The expression~\eqref{eq:Gamma2,1_real} for the $\UP\UP\down$ three-body loss rate remains valid. In the expression~\eqref{eq:Gamma1,2_real} for the $\UP\down\down$ three-body loss rate, $s$ is replaced by $\tilde{s}$, i.e. $\Gamma_{1,2} \simeq -\frac{\hbar}{m}\ 8 \tilde{s} (\tilde{s}+1) \, \sum_{\md=-1}^1 C^{(\md)}_{1,2} \ {\rIm} \,a^{(\md)}_{1,2}$. Since the $\UP\down\down$ three-body parameters $a^{(\md)}_{1,2}$ are now of order $b^{2\tilde{s}}$ whereas the $\UP\UP\down$ three-body parameters $a^{(\md)}_{2,1}$ are of order $b^{2s}$, the $\UP\down\down$ three-body loss rate $\Gamma_{1,2} \propto b^{2\tilde{s}}$ is negligible (in the zero-range regime) compared to the $\UP\UP\down$ three-body loss rate $\Gamma_{2,1} \propto b^{2s}$. Hence
\begin{equation}\label{eq:Gamma3_imbal}
\boxed{\Gamma_3 \simeq \Gamma_{2,1} \simeq -\frac{\hbar}{m}\ 8\, s\, (s+1) \, \sum_{\md=-1}^1 C^{(\md)}_{2,1} \ {\rIm} \,a^{(\md)}_{2,1}.}
\end{equation}
The expression~\eqref{eq:dE2} for the two-body contribution to the energy correction remains valid. For the three-body contribution to the energy correction, we obtain
\begin{equation}\label{eq:dE3_imbal}
\boxed{\delta\!E_3 \ \simeq\ \
\frac{\hbar^2}{m}\ 4\, s\ (s+1)\ \sum_{\md=-1}^1 C_{2,1}^{(\md)}\ a_{2,1}^{(\md)} }
\end{equation}
because the leading-order contribution again comes from $\UP\UP\down$ triplets.\footnote{Indeed, the contribution to $\delta\!E_3$ from $\UP\UP\down$ triplets of nearby particles, given by the r.h.s. of~\eqref{eq:dE3_imbal}, is $\propto b^{2s}$, whereas the $\UP\down\down$ contribution is $\simeq
\frac{\hbar^2}{m}\ \ 4\, \tilde{s}\ (\tilde{s}+1)\ \ \sum_{\md=-1}^1 C_{1,2}^{(\md)}\ a_{1,2}^{(\md)}\,\propto b^{2\tilde{s}}$.}
A peculiar situation takes place for $s<1/2$ (i.e. for $m_\UP/m_\down > 12.313\ldots\ $): $\ \delta\!E_3\propto~\!\!\!b^{2s}$ dominates over $\delta\!E_2 \propto\!\!b\,$; i.e., the finite-range correction mainly comes from configurations with three nearby particles, rather than two nearby particles. This was already pointed out in~\cite{GaoEndoCastin} at the level of the third virial coefficient. Few-body and many-body numerical computations are often performed with finite-range interactions and extrapolated to the zero-range limit, see e.g.~\cite{BlumeRelations,Blume3bodyResPRA} and~\cite{zhenyaPRL,zhenyaNJP,zhenyas_crossover,GezerlisXi,CarlsonAFQMC,Carlson_C_relations,ValeC_precise,GandolfiProceeding,KaplanTrap,KaplanXi,Goulko_UFG_2016,SchonenbergConduit,AlhassidPairing,AlhassidC,LeeCondFrac} respectively; to accurately perform such $b\to0$ extrapolations in the regime $s<1/2$, it may be important to include the $b^{2s}$ scaling.\footnote{Naturally, the $\propto b^{2s}$ scaling may be hard to distinguish from the regular $\propto b$ scaling if $s$ is close to $1/2$ (see e.g.~\cite[Fig.~9(b)]{Blume3bodyResPRA}).}$^{,}$\footnote{We take this opportunity to recall two other subtleties relevant to zero-range extrapolations (i.e. to continuum extrapolations in the case of lattice models, where the interaction range is set by the lattice spacing). Being due to breaking of Galilean invariance, they arise for any mass ratio. The first subtlety, reported in~\cite{zhenyaNJP,LeChapitreIn2,WernerCastinRelationsFermions}, is that for lattice models (where the interaction range $b$ is set by the lattice spacing), the $\propto r_e$ term given by the r.h.s. of~\eqref{eq:dE2} is not the only $\propto b$ contribution to the two-body finite-range energy-correction $\delta\!E_2$: There is a second contribution to $\delta\!E_2$, proportional to a parameter $R_e$ (that parameter $R_e$ quantifies the dependence of the two-body vacuum T-matrix on the center-of-mass momentum, which arises from lattice-induced breaking of Galilean invariance). Therefore, if one wishes to cancel the $\propto b$ term in the $b\to0$ continuum extrapolation, one needs to tune to zero not only $r_e$ (as done in~\cite{KaplanTrap,KaplanXi,LuuLuscher} and for one of the dispersion relations considered in~\cite{CarlsonAFQMC}) but also $R_e$ (and if the dispersion relation has a cusp at the edge of the Brillouin zone, there is a third $\propto b$ contribution to $\delta\!E_2$, in a finite box with periodic boundary conditions~\cite{WernerCastinRelationsFermions}). The second subtlety, reported in~\cite{WernerCastinRelationsFermions} and further evidenced in~\cite{AlhassidPairing}, arises if one restricts single-particle momenta to a ball of radius $\Lambda$ (with $\Lambda = \pi/b$ for lattice models, with $b$ the lattice spacing), i.e. if one takes a dispersion relation equal to $+\infty$ for momenta larger than $\Lambda$, as was done in~\cite{bulgacQMC,BulgacThermodynamics,BulgacCrossover,zhenyas_crossover,BulgacPG,BulgacPG2,BulgacPairing2013,BulgacViscosity}: The universal zero-range model is not obtained in the limit $\Lambda\to\infty$, because the hard cutoff induces a dependence of the two-body T-matrix on the center-of-mass momentum, and this dependence surprisingly survives for $\Lambda\to\infty$. We expect the same problem in~\cite{KaplanTrap,KaplanXi} where such a spherical cutoff was also used.}
\subsection{Relations within the generalized zero-range model}\label{app:sec:simpl}
In this Section, we assume that $m_\UP/m_\down$ is larger than $8.619\ldots\ $, so that $s < 1$. This ensures that a three-body wavefunction diverging as $R^{-s-2}$ at small $R$ remains square integrable. This allows one to define a \emph{generalized zero-range model} (GZRM) by supplementing the two-body contact condition~\eqref{eq:BP} (involving the two-body scattering length) by the following three-body contact condition (involving the three-body parameters):
There exists $B_\md$ such that
\begin{equation}\label{eq:cc_a21}
\boxed{\psi(\rr_1,\,\ldots,\,\rr_N)
\underset{R\,\to\,0}{\sim} \
\ \frac{1}{R^2}
\ \sum_{\md=-1}^{+1} \
\left(- \frac{a_3^{(\md)}}{R^{s}} + R^{s} + o(R^{s})
\right)
\ \phi_\md(\Oo)\ {B}_\md(\CC;\rr_4,\,\ldots,\,\rr_N)}
\end{equation}
in the limit $R\to0$ where particles 1,2,3 approach each other, while keeping fixed their hyperangles $\Oo$, their center-of-mass $\CC$, and the positions of the other particles $\rr_4, \ldots, \rr_N$. By antisymmetry, Eq.~\eqref{eq:cc_a21} imposes a similar condition when any triplet of particles with spins $\UP\UP\down$ approach each other. Apart from the two-body and three-body contact conditions~\eqref{eq:BP}, \eqref{eq:cc_a21}, the GZRM is defined (like the standard zero-range model) by the Schr\"odinger equation without any interaction potential~\eqref{eq:schro_imbal}. Note that for vanishing three-body parameters, the GZRM reduces to the standard zero-range model (ZRM).
In the regime where $\ktyp^{\phantom{t}2s}\ |a_3^{(\md)}| \ \gtrsim \ 1$ (for at least one value of $\md$), the ZRM becomes irrelevant, whereas the GZRM remains applicable (provided $\ktyp b \ll 1$). This regime can be reached
\begin{itemize}
\item near a ($\UP\UP\down$) three-body resonance~\cite{LeChapitreIn2,Petrov3fermions,WernerSym,NishidaSonTan3bRes,Blume3bodyResPRL,GandolfiCarlson3bRes,Blume3bodyResPRA,Sadeghpour_Blume_3bodyRes,Kartavtsev3bRes}
\item when the mass ratio $m_\UP/m_\down$ is only slightly smaller than the threshold $13.607\ldots$ where the three-body Efimov effect appears, so that $s$ is small~\cite{GaoEndoCastin}.
\end{itemize}
%\noindent \underline{\it Remarks:}
\begin{remas*}\leavevmode
\begin{itemize}
\item If the underlying microscopic interactions are rotationally-invariant, then $a_3^{(\md)}$ does not depend on~$\md$.
\item In the present GZRM, in the limit where three particles of spins $\UP\down\down$ approach each other, the asymptotic behavior of the wavefunction has the same form~\eqref{eq:R0t_imbal} than for the ZRM. In other words, the behavior $\psi \propto R^{-2-\tilde{s}}$ is forbidden. Such a wavefunction would not be a normalizable at small~$R$, since $\tilde{s}$ is always larger than one.\footnote{\bl Hence there are no $\UP\down\down$ three-body parameters within the GZRM, which is why we denoted the $\UP\UP\down$ three-body parameters by $a_3^{(\md)}$ instead of the more cumbersome notation $a_{2,1}^{(\md)}$.} This fact does not prevent one from rederiving the relations~(\eqref{eq:Gamma3_imbal},\eqref{eq:dE3_imbal}), which come from configurations with nearby $\UP\UP\down$ triplets.
\item Extensions of the GZRM to the regime $s\geq 1$, which were introduced recently~\cite{PricoupenkoNbodyRes}, are beyond the scope of this work.
\item For mathematically rigorous studies of the GZRM in the three-particle case, see~\cite{Minlos3body_res_2014,Teta_3body_res}.
\end{itemize}
\end{remas*}
\subsubsection{Derivative of the energy with respect to the three-body parameters}\label{app:imbal,subsec:dE}
Within the GZRM, the derivatives of the energy w.r.t. the three-body parameters are given by the three-body contacts,
\begin{equation}\label{eq:dE_da3}
\boxed{ \frac{\partial E}{\partial a_3^{(\md)}} \ =\ \frac{\hbar^2}{m}\ 4\,s\,(s + 1)\ C_{2,1}^{(\md)}\,.}
\end{equation} Here the derivative w.r.t. $a_3^{(\md)}$ is taken at fixed value of $a_2$ and of the other $a_3^{(\md')}$ with ${\md'\neq\md}$.
\\The three-body contacts $C_{2,1}^{(\md)}$ are still defined by~\eqref{eq:C21_BB_imbal}.
%{\it \underline{Remark:}}
\begin{rema}
Within the ZRM, the derivative of the energy w.r.t. $a_2$ is given by $C_2$, see~\eqref{eq:dE_da2}. Within the GZRM, this relation also holds, with $dE/d(-1/a_2)$ replaced by the partial derivative $\partial E/\partial(-1/a_2)$ taken at fixed $a_3^{(\md)}$. This can be justified by using the derivation presented
in~\cite[Sec.~IV.C]{WernerCastinRelationsFermions} (at~fixed three-body parameters, there is no additional contribution to the energy variation coming from nearby triplets of particles, as we will see below).
\end{rema}
To derive~\eqref{eq:dE_da3} we proceed as follows.\footnote{This derivation is similar to the case of three bosons treated in of~\cite[App.~A]{WernerCastinRelationsBosons}.} We consider two wavefunctions $\psi$ and $\tilde{\psi}$, which are stationary states of the zero-range model for different sets of three-body parameters $(a_3^{(\md)})$ and~$(\tilde{a}_3^{(\md)})$, but the same two-body scattering length $a_2$.
\\Denoting the corresponding eigenenergies by $E$ and $\tilde{E}$, we have
\begin{equation*}
\delta := \left\la \psi, H \tilde{\psi} \right\ra - \left\la H \psi, \tilde{\psi} \right\ra = \left(\tilde{E}-E\right)\,\left\la \psi| \tilde{\psi} \right\ra.
\end{equation*}
On the other hand, from the Schr\"odinger equations for $\psi$ and $\tilde{\psi}$,
\begin{equation}\label{eq:delta_sum_i}
\delta = - \sum_{i=1}^N \frac{\hbar^2}{2m_i} \int \dr_1 \ldots \dr_N
\, \left(\psi^* \Delta_{\rr_i} \tilde{\psi} - \tilde{\psi}\, \Delta_{\rr_i} \psi^* \right).
\end{equation}
As we will see below, there is a contribution $\delta_0$ to $\delta$ coming from the configurations where particles 1, 2 and 3 are close to each other. By symmetry, there is an identical contribution from the configurations where any set of three particles with spins $\UP\UP\down$ are close to each other. Hence $\delta$ equals $\delta_0$ times the number of such three-particle sets, $N_\UP (N_\UP-1) N_\down / 2$. To evaluate $\delta_0$, we only need to keep the terms $i=1,2,3$ in Eq.~\eqref{eq:delta_sum_i}, which gives, after the change of coordinates $(\rr_1,\rr_2,\rr_3) \longrightarrow (\CC,R,\Oo)$,
\begin{multline}\label{eq:delta_R}
\delta_0 = -\frac{\hbar^2}{m}
\cos^3\!\varphi\
\int_0^\infty dR\,R^5
\int d^5\Omega
\,\dC
\, \dr_4 \ldots \dr_N\\
\left\{
\Fr^*
\left(
\frac{\partial^2}{\partial R^2} + \frac{1}{R} \frac{\partial}{\partial R} + \frac{T_\Oo}{R^2} + \frac{\Delta_\CC}{3}
\right)
\tilde{\Fr} - \left[\Fr^* \leftrightarrow \tilde{\Fr}\right]
\right\}
\end{multline}
where $\Fr := R^2\,\psi$, $\tilde{\Fr} := R^2\,\tilde{\psi}$, and $T_\Oo$ is defined by~\eqref{eq:def_TOm}. Since the two-body scattering length is the same for $\psi$ and $\tilde{\psi}$, we only need to keep the terms involving derivatives with respect to $R$ in Eq.~\eqref{eq:delta_R}.\footnote{For a more detailed justification of this step, see the reasoning around~\cite[Eq.~(A7)]{WernerCastinRelationsBosons}.} Transforming the integral over $R$ into a boundary term at $R\to0$, we then get
\begin{equation*}
\delta_0 = \frac{\hbar^2}{m}
\cos^3\!\varphi\,
\int \,\dC
\, \dr_4 \ldots \dr_N
\,\int d^5\Omega
\ \lim_{R\,\to\,0} R \left(
\Fr^* \ \frac{\partial \tilde{\Fr}}{\partial R} -
\tilde{\Fr} \ \frac{\partial \Fr^*}{\partial R}
\right).
\end{equation*}
The result~\eqref{eq:dE_da3} follows by using the three-body contact conditions [given by Eq.~\eqref{eq:cc_a21} for $\psi$, and the same condition with $\tilde{a}_3^{(\md)}$ for $\tilde{\psi}$], taking the limit $\tilde{a}_3^{(\md)} \to {a}_3^{(\md)}$, and using the expression~\eqref{eq:C21_BB_imbal} of the three-body contacts.
\subsubsection{Three-body loss rate}
Let us consider the GZRM with complex three-body parameters, $a_3^{(\md)}\in\bbC$. The three-body loss rate is then given by
\begin{equation}\label{eq:Gamma3_GZRM}
\boxed{\Gamma_3 = -\frac{\hbar}{m}\ 8\, s\, (s+1) \, \sum_{\md=-1}^1 C^{(\md)}_{2,1} \ {\rIm} \,a^{(\md)}_{3}.}
\end{equation}
Within the GZRM, there are only $\UP\UP\down$ losses, and no $\UP\down\down$ losses, i.e. $\Gamma_3 = \Gamma_{2,1}$ and $\Gamma_{1,2}=0$.
To derive this relation, we proceed very similarly to the bosonic case of~\cite{WernerCastinRelationsBosons}. We consider a stationary state of the GZRM, i.e. a solution $\psi(\rr_1,\,\ldots,\,\rr_N)$ of the stationary Schr\"odinger equation~\eqref{eq:schro_imbal} together with the two-body and three-body contact conditions (\eqref{eq:BP},\eqref{eq:cc_a21}). The corresponding solution of the time-dependent Schr\"odinger equation is $\Psi(\rr_1,\,\ldots,\,\rr_N;t) = \psi(\rr_1,\,\ldots,\,\rr_N) e^{-i E t / \hbar}$, and $\Gamma_3 = -\frac{\partial}{\partial t}|_{t=0} \int \dr_1 \ldots \dr_N \, |\Psi(\rr_1,\,\ldots,\,\rr_N;t)|^2$, with $\psi$ normalized to unity. Excluding from the integration domain the regions where $R_{ijk}<\epsilon$ and taking the limit $\epsilon\to0$, the continuity equation and the three-body contact conditions leads to the final expression, where all extra mass-ratio dependent factors, arising e.g. from~\eqref{eq:Ximbal}, divide out.
\subsection{Validity conditions}\label{sec:max_ratio}
Let us denote by $s_{j_\UP,j_\down}$ the scaling exponent of the unitary $(j_\UP+j_\down)$ body problem with $j_\UP$ particles of spin $\UP$ and $j_\down$ particles of spin $\down$ (so that $s_{2,1} \equiv s$ and $s_{1,2} \equiv \tilde{s}$\,). We expect the relations~\eqref{eq:N21_imb},\eqref{eq:N3_imb},\eqref{eq:NP_tail_imb},\eqref{eq:N_k1k3_imb} to be valid under the condition ${\rMin}\{ s_{j_\UP,j_\down}; \ j_\UP \geq 2, \ j_\down \geq 1, \ j_\UP + j_\down > 3 \} > s$. Indeed, we expect a contribution from the $(j_\UP+j_\down)$ body problem to $N_{2,1}(\epsilon)$ at small $\epsilon$ [resp. to $\bar{N}_P(K)$ at large $K$] scaling as $\epsilon^{2s_{j_\UP,j_\down}+2}$ (resp.~$1/K^{2s_{j_\UP,j_\down}+4}$), which dominates over the contribution from the three-body problem when $s_{j_\UP,j_\down}~~ 3 \}>\tilde{s}$. The former condition breaks down when $m_\UP/m_\down$ is near the thresholds where the four-body~\cite{CMP} and five-body~\cite{PetrovBazak5body} Efimov effects appear, where $s_{3,1}$ and $s_{4,1}$ become small. Based on existing data, we expect that both of the above validity conditions are satisfied at least in the range $m_\UP/m_\down \leq 10$. Indeed, in this range, we have $s_{3,1}>s$ and $s_{4,1}>s$~\cite{PetrovBazak5body}, $s_{5,1}>s$~\cite{Bazak6body}, $s_{2,2}>\tilde{s}$~\cite{Blume3bodyResPRA}, $s_{1,3} > \tilde{s}$~\cite{PetrovBazak5body}, and the trends of available data suggest that the conditions will also hold for larger values of $j_\UP$ or $j_\down$. We conservatively restricted to $m_\UP/m_\down \leq 10$ because $s_{2,2}$ was not computed beyond this range, but the conditions $s_{3,1}>s$ and $s_{4,1}>s$ actually hold up to at least $m_\UP/m_\down = 13.2$~\cite{PetrovBazak5body}.
\section{Generalization to statistical mixtures and non-stationary states}\label{app:non_stat}
Many of the relations derived for stationary states in the main text and in App.~\ref{app:imbal} are directly generalizable to non-stationary states and statistical mixtures, similarly to the relations involving the two-body contact~\cite{TanEnergetics,TanLargeMomentum,WernerCastinRelationsFermions}. Indeed, Eqs.~\eqref{eq:R0}, \eqref{eq:C3_BB}, \eqref{eq:R0t}, \eqref{eq:C3_BBt} remain valid for any non-pathological linear combinations of stationary states (including solutions of time-dependent problems, where $a_2$ and the trapping potential $U$ can depend on time); and relations~\eqref{eq:N3}, \eqref{eq:N21}, \eqref{eq:N12}, \eqref{eq:g21_C3}, \eqref{eq:g12_C3}, \eqref{eq:NP_tail} remain true for arbitrary non-pathological statistical mixtures of such pure states (including the case of thermal equilibrium). Here, non-pathological means that the occupation probabilities of stationary states should decay sufficiently quickly (which includes the simple case where only a finite number of states are populated); more specifically, for a pure state, it is necessary that Eqs.~\eqref{eq:BP}, \eqref{eq:R0_b} still hold, while for a statistical mixture $\hat{\rho} = \sum_i c_i |\psi_i\rangle \langle \psi_i |$, it is necessary that the three-body contact of the mixture [defined by Eq.~\eqref{eq:N3}] equals $\sum_i c_i \,C_{3,i}$ with $C_{3,i}$ the three-body contact of state~$\psi_i$.
Moreover, at thermal equilibrium, the thermally averaged loss rates are given by~\eqref{eq:Gamma3},\eqref{eq:Gamma2,1}, \eqref{eq:Gamma1,2} for simple interactions; for more general interactions they are given by~\eqref{eq:Gamma2,1_real}, \eqref{eq:Gamma1,2_real}, and by~\eqref{eq:Gamma3_real} for a homogeneous unpolarized gas. Furthermore, \eqref{eq:dE_21_12}, \eqref{eq:dE3} remain valid in the canonical ensemble, with $\delta\!E_3$ the energy difference (between finite-range and zero-range models) taken at fixed entropy (which equals the free-energy difference at fixed temperature).
As for the wavevector $\ktyp$ that appears in the validity conditions, it should be defined by the same procedure as before and then taking the maximum over all significantly populated eigenstates; for example for the balanced unitary gas at thermal equilibrium, $\ktyp = {\max}(k_F, k_T)$ with $k_T$ the thermal wavevector, defined by $\hbar^2 k_T^2 / (2 m) = k_B T$.
%\bibliographystyle{crunsrt}
%\bibliography{Werner}
\printbibliography
\end{document}
~~