The role of causality in tunable Fermi gas condensates

There has been a very fruitful exchange of ideas between particle/astro-physicists and condensed matter theorists over several decades. Specifically, black holes and other phenomena require an improved understanding of semiclassical gravity and the role of the even more inaccessible Planck scale. Causality plays a key role in this and, given that direct tests are experimentally impossible, analogies have been sought in condensed matter physics that seek to replicate this causality in the laboratory. Analogies with black holes are directly due to the pioneering work of Unruh [1], which showed the similarities between event horizons in fluids (most simply when the flow rate exceeds the local speed of sound) and the event horizons of black holes. The programme of 'analogue gravity', begun to exploit these similarities, has generated over 700 papers since its implementation (see the recent review article [2]).

One of the great achievements of the last several years has been the construction of condensates in which the speed of sound, c_s, is tunable. Systems whose causal behaviour is so simply controlled are potentially good candidates with which to explore analogies. In this paper we explore the way in which tunable gases of cold Fermionic atoms, which pair to form bosonic condensates, can produce and propagate phonons in mimickry of particle production in the early universe. Before doing so we need to address a basic problem in analogue gravity modelling. This is that there is an immediate mismatch; the early (and late) universe is relativistic, usually locally Lorentzian, whereas condensed matter systems and condensates are Galilean in their space–time symmetry. In practice we rely on the gapless phonon mode (with dispersion relation ω = c_sk) to mimic the relativistic behaviour of light in the geodesic equations that control its behaviour. However, at best the parallel breaks down at the Debye scale (the minimum wavelength of the phonon mode). This breakdown of speed of sound (SoS) Lorentzian behaviour has been considered [3, 4] a counterpart to the possible breakdown of true Lorentzian behaviour at the Planck scale.

For the purpose of drawing analogies with phonon production, and more generally, it is important to develop a formalism that interpolates between the SoS Lorentzian behaviour of the phonon, where appropriate, and the underlying Galilean invariance that provides the bedrock to all our analysis. The main result of this paper is to provide such a formalism for Fermi gases whose speed of sound is tunable in an external magnetic field by virtue of a narrow Feshbach resonance.

The organization of the paper is as follows. We first show how the speed of sound in di-fermion (or diatomic) condensates is tunable by the application of an external magnetic field. This is because the effective bosonic theory is not the phonon Goldstone model of the familiar Gross–Pitaevskii (GP) equation, but a Higgs–Goldstone model where the Higgs field is represented by the gapped diatom density fluctuations. We have discussed this elsewhere [5–8], to which the reader is referred for greater detail. Here we take a further step by tracing out the density fluctuations to give an effective low-energy phonon model for the phonon field θ alone. We then use this to discuss spontaneous phonon production in a rapid magnetic field quench.

There is a problem of familiarity (or lack of it) with Galilean invariance. We work with a Lagrangian formalism rather than the more familiar canonical Hamiltonian approach. Inevitably, we have to make approximations and it is crucial that these approximations preserve the underlying Galilean symmetry. Our results are straightforward to those readers familiar with explicit Galilean invariant approaches to condensed matter theory as in [9–11], for example. For those readers less familiar, suffice to say that the fundamental Galilean phonon scalar is the combination

$$\begin{eqnarray} &&G(\theta )=\dot {\theta }+(\nabla \theta )^{2}/2 M \end{eqnarray} \tag{ 1 }$$

for diatoms with mass M = 2m, where m is the atomic mass.

Our main result is that the phonon action takes the non-local form

$$\begin{eqnarray} {S}_{\mathrm{eff}}[\theta ]&=&\int {\mathrm{d}}^{4}x \hspace{0.167em} \left \{ \frac{{N}_{0}}{4}~{G}^{2}(\theta (x))-\frac{1}{2}{\rho }_{0}G(\theta (x))\right \} \nonumber\\ &&+~\frac{{\alpha }^{2}}{2}\int {\mathrm{d}}^{4}x \hspace{0.167em} \int {\mathrm{d}}^{4}{x}^{\prime}\hspace{0.167em} G(\theta (x))\nonumber\\ &&\times ~{D}_{\mathrm{ F}}(x,{x}^{\prime})\hspace{0.167em} G(\theta ({x}^{\prime})), \end{eqnarray} \tag{ 2 }$$

with coefficients to be discussed later, where iD_F(x,x') describes the exchange of virtual density fluctuations in the condensate as diatoms associate and disassociate (with constant overall density). This density field can be understood as the Galilean scalar 'Higgs' boson (mentioned earlier) to complement the Goldstone phonon. Its properties, as well as its coupling strength α, depend on the external magnetic field with which we tune the binding energy of the Feshbach resonance. The form (2) is a straightforward generalization of the Galilean invariant action of [10] and a natural relation to that of [11], for example.

The SoS Lorentzian limit is obtained from the Galilean formulation by the substitution [9, 10]

$$\begin{eqnarray} &&G(\theta )\rightarrow (\nabla \theta )^{2}/2 M,\tqs {G}^{2}(\theta )\rightarrow {\dot {\theta }}^{2}. \end{eqnarray} \tag{ 3 }$$

The outcome is the action

$$\begin{eqnarray} {S}_{\mathrm{eff}}^{\mathrm{qu}}[\theta ]&=&\frac{1}{4}\int {\mathrm{d}}^{4}x \hspace{0.167em} \left \{ {N}_{0}{\dot {\theta }}^{2}(x)-\frac{{\rho }_{0}}{2 m}(\nabla \theta )^{2}\right \} \nonumber\\ &&+~\frac{{\alpha }^{2}}{2}\int {\mathrm{d}}^{4}x \hspace{0.167em} \int {\mathrm{d}}^{4}{x}^{\prime}\dot {\theta }(x)\hspace{0.167em} {D}_{\mathrm{ F}}(x-{x}^{\prime})\hspace{0.167em} \dot {\theta }({x}^{\prime}). \end{eqnarray} \tag{ 4 }$$

The details will be given later but it is already clear that this action essentially describes an SoS Lorentzian phonon, allowing for some 'fuzziness' in the sound-cone due to density fluctuations.

To bring this back to familiar territory, in those regimes for which the action (2) is approximately local with the fuzziness smoothed away, the resulting equations of motion can be rewritten in the Gross–Pitaevskii form

$$\begin{eqnarray} &&\mathrm{i}\hbar \dot {\psi }+\frac{{\hbar }^{2}}{2 m}{\nabla }^{2}\psi +2 m{c}_{\mathrm{ s}}^{2}\psi -\frac{2 m{c}_{\mathrm{s}}^{2}}{{\rho }_{0}}\psi \vert \psi \vert ^{2}=0, \end{eqnarray} \tag{ 5 }$$

where $\psi =\sqrt{\rho }\hspace{0.167em} \exp (\mathrm{i}\theta )$ for a condensate of density ρ, with θ the phonon field with speed of sound c_s. We have used this GP equation elsewhere [7] to examine possible spontaneous vortex production in field quenches, but we shall not pursue that here, and work from (2).

The action (2) permits different approximations. A phonon moving with variable speed c_s in space and time can be interpreted as moving along a geodesic path in a varying space–time metric. This is the basis of analogue gravity, for which the geodesic equation for the phonon following from (2) in the hydrodynamic limit is

$$\begin{eqnarray} &&\frac{1}{\sqrt{-g}}{\partial }_{\mu }(\sqrt{-g}\hspace{0.167em} {g}^{\mu \nu }{\partial }_{\nu }\theta )=0, \end{eqnarray} \tag{ 6 }$$

where g is the acoustic (or hydrodynamic) metric of the condensate. For example, for homogeneous condensates it is relatively simple experimentally to quench the condensate so that its acoustic metric takes Freedman–Robinson–Walker (FRW) form. Spontaneous phonon production in this metric is the counterpart to the familiar process of particle production in the early FRW universe and it is this analogy that we pursue, to demonstrate how simple the formalism is. Despite the underlying Galilean invariance, it has been argued [12–14] that there is still a window in which geodesic behaviour of (6) dominates. We will examine this below.

We stress that for us the default is the underlying Galilean invariance for which the SoS Lorentzian limit may sometimes be appropriate.

We now show how a cold Fermi gas can lead to a condensate described by the action (2). The forces between alkali Fermi atoms are such that, for weak coupling, they form Cooper pairs (the BCS regime), correlated in momenta, whereas, for strong coupling, they form molecules or diatoms (the BEC regime), correlated in position. In a two-channel model with a narrow Feshbach resonance the binding energy of this resonance can be changed by the application of an external magnetic field because of the different Zeeman effects between the channels. The outcome is to change the strength of the force and permit us to take the system smoothly from the BCS regime (with negative s-wave scattering length a_S) to the BEC regime (with positive a_S) through the unitary limit at which a_S diverges. The speed of sound also varies smoothly through this crossover from ${c}_{\mathrm{BCS}}={v}_{\mathrm{F}}/\sqrt{3}$ in the deep BCS regime, where v_F is the Fermi velocity, to c_s → 0 in the deep BEC regime. This remarkable control over the sound speed by the simple application of an external field gives us control over the causal properties of the gas.

Analytically this is most simply understood if we take the resonance to be very narrow, for which the mean-field equation, necessary for analytic approximations, can be trusted [15]. Further, for a narrow resonance, the condensate order parameter is the Feshbach resonance field itself, which, as we have shown [8], leads to a single-fluid model in the hydrodynamic limit.

Our starting point is the exemplary 'two-channel' microscopic action [15, 16] (in units in which ħ = 1)

$$\begin{eqnarray} S&=&\int \mathrm{d}t \hspace{0.167em} {\mathrm{d}}^{3}x \hspace{0.167em} \Biggl \{ \mathop{\sum }\limits _{\uparrow ,\downarrow }\psi _{\sigma }^{\ast }(x)~\left [\mathrm{i}~{\partial }_{t}+\frac{{\nabla }^{2}}{2 m}+\mu \right ]~{\psi }_{\sigma }(x)\nonumber\\ &&+~{\varphi }^{\ast }(x)~\left [\mathrm{i}~{\partial }_{t}+\frac{{\nabla }^{2}}{2 M}+2 \mu -\nu \right ]~\varphi (x)\nonumber\\ &&-~g\left [{\varphi }^{\ast }(x)~{\psi }_{\downarrow }(x)~{\psi }_{\uparrow }(x)+\varphi (x){\psi }_{\uparrow }^{\ast }(x)~{\psi }_{\downarrow }^{\ast }(x)\right ]\Biggr \} \end{eqnarray} \tag{ 7 }$$

for cold (T = 0) Fermi fields ψ_σ with spin label σ = (↑,↓). The diatomic field φ describes a narrow bound-state (Feshbach) resonance with tunable binding energy ν and mass M = 2m. The idealization of ignoring the self-interactions between fermions and diatomic molecules is well suited for describing the crossover phenomena, but breaks down in the deep BCS and deep BEC regimes [15], although not in such a way as to compromise our results.

We restrict ourselves to the mean-field approximation. S is quadratic in the Fermi fields. Integrating them out gives an effective action for the order parameter φ(x) =− |φ(x)| e^iθ(x) alone. This action possesses a U(1) invariance under θ → θ + const., which is spontaneously broken by space–time constant gap solutions |ϕ(x)| = |ϕ₀| ≠ 0. We expand in fluctuations about ϕ₀, but not with the decomposition ϕ = ϕ₀ + δϕ, since we need to preserve Galilean invariance at each step of the approximation.

The Galilean invariants of the theory are the density fluctuation $\delta \vert \phi \vert =\vert \phi \vert -\vert {\phi }_{0}\vert ,G(\theta )=\dot {\theta }+(\nabla \theta )^{2}/4 m$ as mentioned earlier, and ${D}_{t}(\delta \vert \phi \vert ,\theta )=\dot {(\delta \vert \phi \vert )}+\nabla \theta .\nabla (\delta \vert \phi \vert )/2 m$. D_t is the comoving time derivative of δ|ϕ| in a fluid with fluid velocity ∇θ/2m. θ(x) is not small. Let us rescale δ|ϕ| to δ|ϕ| = κ. The resulting Galilean invariant effective action for the long wavelength, low-frequency condensate takes the form [5, 6]

$$\begin{eqnarray} {S}_{\mathrm{eff}}[\theta ,\epsilon ]&=&\int {\mathrm{d}}^{4}x \hspace{0.167em} \Biggl [\frac{{N}_{0}}{4}~{G}^{2}(\theta ,\epsilon )-\frac{1}{2}{\rho }_{0}G(\theta ,\epsilon )\nonumber\\ &&-~\alpha \epsilon G(\theta ,\epsilon )+\frac{1}{4}\eta {D}_{t}^{2}(\epsilon ,\theta )-\frac{1}{4}{\bar {M}}^{2}{\epsilon }^{2}\Biggr ]. \end{eqnarray} \tag{ 8 }$$

The scale factor κ, defined in the appendix, is chosen so that on extending G(θ) to $G(\theta ,\epsilon )=\dot {\theta }+(\nabla \theta )^{2}/4 m+(\nabla \epsilon )^{2}/4 m,\epsilon$ has the same coefficients as θ in its spatial derivatives.

The quadratic part ${S}_{\mathrm{eff}}^{\mathrm{qu}}$ of S_eff is

$$\begin{eqnarray} {S}_{\mathrm{eff}}^{\mathrm{qu}}[\theta ,\epsilon ]&=&\int {\mathrm{d}}^{4}x \hspace{0.167em} \Biggl [\frac{{N}_{0}}{4}{\dot {\theta }}^{2}-\frac{{\rho }_{0}}{8 m}(\nabla \theta )^{2}\nonumber\\ &&-~\alpha \epsilon \dot {\theta }+\frac{1}{4}\eta {\dot {\epsilon }}^{2}-\frac{{\rho }_{0}}{8 m}(\nabla \epsilon )^{2}-\frac{1}{4}{\bar {M}}^{2}{\epsilon }^{2}\Biggr ]. \end{eqnarray} \tag{ 9 }$$

We see immediately from (9) that the action describes a gapless (i.e. SoS Lorentzian, or SoS relativistic) Goldstone mode θ, the phonon, and a gapped (but also relativistic) Higgs mode. The underlying Galilean invariance is only present in their time-derivative coupling with strength α.

As for the coefficients in (9), ${\rho }_{0}={\rho }_{0}^{\mathrm{F}}+{\rho }_{0}^{\mathrm{B}}$ is the total (fixed) fermion number density where ${\rho }_{0}^{\mathrm{F}}$ is the explicit fermion density and ${\rho }_{0}^{\mathrm{B}}=2\vert {\phi }_{0}\vert ^{2}$ is due to molecules (two fermions per molecule). For the evolving system the molecular or diatomic density is ${\rho }^{\mathrm{B}}=2\vert \phi \vert ^{2}={\rho }_{0}^{\mathrm{B}}+4 \delta \vert \phi \vert \hspace{0.167em} \vert {\phi }_{0}\vert$. This shows that  ∝ δρ^B is the scaled molecular density fluctuation, describing the repeated dissociation of molecules into atom pairs and their reconversion into molecules [16, 17]. Otherwise, the coefficients are somewhat opaque and relegated to the appendix. We give exemplary plots in figure 1 of the most important combinations.

Zoom In Zoom Out Reset image size

We observe that S_eff[θ,] of (8) is quadratic in . Integrating out the field gives the effective non-local action of (2) in terms of the phonon alone,

$$\begin{eqnarray} {S}_{\mathrm{eff}}[\theta ]&=&\int {\mathrm{d}}^{4}x \hspace{0.167em} \left \{ \frac{{N}_{0}}{4}~{G}^{2}(\theta (x))-\frac{1}{2}{\rho }_{0}G(\theta (x))\right \} \nonumber\\ &&+~\frac{1}{2}\int {\mathrm{d}}^{4}x \hspace{0.167em} \int {\mathrm{d}}^{4}{x}^{\prime}\hspace{0.167em} \alpha (x)G(\theta (x))\nonumber\\ &&\times ~{D}_{\mathrm{ F}}(x,{x}^{\prime})\hspace{0.167em} \alpha ({x}^{\prime})G(\theta ({x}^{\prime})), \end{eqnarray} \tag{ 10 }$$

where iD_F(x,x') = 〈T((x)(x'))〉₀, determined from the quadratic action S_eff[0,] of (9) when θ ≡ 0, and we now allow α to be space–time dependent.

As we said in section 1, this is the main formal result of the paper, showing intuitively how virtual (Higgs) density fluctuations mediate the phonon field in an explicitly Galilean invariant way. For constant parameters the quadratic part of S_eff[θ] was already given in (4) as

$$\begin{eqnarray} {S}_{\mathrm{eff}}^{\mathrm{qu}}[\theta ]&=&\frac{1}{4}\int {\mathrm{d}}^{4}x \hspace{0.167em} \left \{ {N}_{0}{\dot {\theta }}^{2}(x)-\frac{{\rho }_{0}}{2 m}(\nabla \theta )^{2}\right \} \nonumber\\ &&+~\frac{{\alpha }^{2}}{2}\int {\mathrm{d}}^{4}x \hspace{0.167em} \int {\mathrm{d}}^{4}{x}^{\prime}\hspace{0.167em} \dot {\theta }(x)\hspace{0.167em} {D}_{\mathrm{ F}}(x-{x}^{\prime})\hspace{0.167em} \dot {\theta }({x}^{\prime}) \end{eqnarray} \tag{ 11 }$$

but with the further information that

$$\begin{eqnarray} {D}_{\mathrm{F}}(x-{x}^{\prime})&=&-\frac{2}{(2 \pi )^{4}}\int \mathrm{d}\omega \hspace{0.167em} {\mathrm{d}}^{3}k\frac{{\mathrm{e}}^{\mathrm{i}\omega (t-{t}^{\prime})}{\mathrm{e}}^{-\mathrm{i}\mathbf{k}(\mathbf{x}-{\mathbf{x}}^{\prime})}}{\eta {\omega }^{2}-{\rho }_{0}{k}^{2}/2 m-{\bar {M}}^{2}}.\nonumber\\ \end{eqnarray} \tag{ 12 }$$

To understand the fuzziness of the sound-cone due to diatom density fluctuations we observe that, in the BEC regime, for which η ≈ 0, we can ignore the ηω² term in the denominator of the propagator for all but the shortest wavelengths. The resulting propagator is local in time, making the sound-cone smooth. However, as we move into the BCS regime both η and ${\bar {M}}^{2}$ increase and the two time derivatives in the second term of (11) are split. There is compensation in that α becomes smaller, vanishing in the deep BCS regime, but there is a characteristic frequency of density oscillations associated with the non-local smearing of the sound-cone. It is mainly in the BEC regime that the speed of sound changes sufficiently rapidly for causal effects to be important. What the above shows, and which we have shown elsewhere by different means [7], is that we cannot push the system much out of the BEC regime, if at all, before analogies break down because the oscillations are too slow to be ignored [16, 17].

From equation (11) we read off the dispersion relation

$$\begin{eqnarray} &&{\omega }^{2}\Biggl ({N}_{0}+\frac{4{\alpha }^{2}}{-\eta {\omega }^{2}+{\rho }_{0}{k}^{2}/2 m+{\bar {M}}^{2}}\Biggr )-\frac{{\rho }_{0}}{2 m}{\mathbf{k}}^{2}=0, \end{eqnarray} \tag{ 13 }$$

a result not achieved so transparently before. In general we allow all parameters to vary in space–time, although subsequently we restrict ourselves to homogeneous systems.

The hydrodynamic (or acoustic) approximation is obtained by approximating the SoS relativistic Higgs propagator D_F(x,x') by its long wavelength, low-energy limit, the instantaneous contact term

$$\begin{eqnarray} &&{D}_{\mathrm{F}}(x-{x}^{\prime})\approx \frac{2}{{\bar {M}}^{2}}\delta (x-{x}^{\prime}), \end{eqnarray} \tag{ 14 }$$

(whether ${\bar {M}}^{2}$ is constant or not). In this ultra-local approximation the density fluctuations do not propagate. That is, in terms of the two-field action (9) we are neglecting the spatial and temporal variation of , in comparison to itself, equivalent to setting ω and k to zero in the denominators of (12) and (13).

This corresponds to making the approximation

$$\begin{eqnarray} &&\fl {S}_{\mathrm{hyd}}[\theta ]=\int {\mathrm{d}}^{4}x \hspace{0.167em} [\frac{1}{4}({N}_{0}+4{\alpha }^{2}/{\bar {M}}^{2})~{G}^{2}(\theta (x))\nonumber\\ &&-~\frac{1}{2}{\rho }_{0}G(\theta (x))] \end{eqnarray} \tag{ 15 }$$

for the action (10) in which we have taken the coefficients to be constant. In the BCS regime, when α ≈ 0 this is just the form of the action given in [10] in the same approximation. (However, in this latter case there was no explicit resonance.) It is also the (exact) long wavelength limit of (13), in which the phonon has the linear dispersion relation ${\omega }^{2}={c}_{\mathrm{s}}^{2}{\mathbf{k}}^{2}$, with speed of sound

$$\begin{eqnarray} &&{c}_{\mathrm{s}}^{2}=\frac{{\rho }_{0}/2 m}{{N}_{0}+4{\alpha }^{2}/{\bar {M}}^{2}}, \end{eqnarray} \tag{ 16 }$$

as can be read directly from (15) on expanding G. In the deep BCS regime ${c}_{\mathrm{s}}\rightarrow {c}_{\mathrm{BCS}}=\sqrt{{\rho }_{0}/2 m{N}_{0}}={v}_{\mathrm{F}}/\sqrt{3}$. On the other hand, in the deep BEC regime (where N₀ and ${\bar {M}}^{2}$ are small, and α = ρ₀), then c_s → 0. To justify the acoustic/hydrodynamic terminology, we note that the equation of motion following from (15) can be written as a continuity equation in which the density fluctuations satisfy the Bernoulli equation [10].

An exemplary graph of ${c}_{\mathrm{s}}^{2}$ is given in figure 2. To a fair approximation, ${c}_{\mathrm{s}}^{2}$ can be approximated as

$$\begin{eqnarray} &&{c}_{\mathrm{s}}^{2}\approx ({c}_{\mathrm{ BCS}}^{2}/2)[1+\tanh (d(g)-b(g)/{k}_{\mathrm{ F}}{a}_{\mathrm{S}})]. \end{eqnarray} \tag{ 17 }$$

This has the great virtue that calculations can be performed analytically. Empirically, b(g) is insensitive to the coupling constant. In fitting d(g) it happens that, as we reduce g, the best-fit to d(g) becomes more negative. Further, on going from the BCS towards the BEC regime, for a smaller value of g the sound speed changes earlier and more quickly.

Zoom In Zoom Out Reset image size

We conclude with an observation on the importance for the calculations above of the Feshbach resonance being narrow. If the resonance is not narrow direct atomic interactions occur, manifest through a four-Fermi interaction in (7). The order parameter now has two terms, a resonance contribution and a di-fermion field, characterizing this interaction [5, 18]. In the same acoustic/hydrodynamic limit as above the equations of motion can now be understood, in general, as representing two coupled fluids [5, 6]. Only when the contact terms overwhelm the resonance exchange do we again recover a single-fluid model [10] (in the hydrodynamic approximation). Otherwise, the single GP equation (5) is replaced by two coupled GP equations [6] or the geodesic equation (6) has two metrics, neither of which lends itself to simple analogy.

3.1. Tuning the condensate

For our narrow resonance, applying a homogeneous external magnetic field changes c_s by changing the s-wave scattering length a_S as [15]

$$\begin{eqnarray} &&\frac{1}{{k}_{\mathrm{F}}{a}_{\mathrm{S}}}\propto (\mathcal{H}-{\mathcal{H}}_{0}), \end{eqnarray} \tag{ 18 }$$

where ₀ is the field required to achieve the unitary limit (|a_S| → ∞). For our first case of interest we pass from the BEC to the BCS regime as increases through ₀.

If we adopt (17), the resulting equations can then be solved analytically, a good approximation is [7], even though the tanh-behaviour slightly overestimates the rate at which the speed of sound vanishes.

In the first instance, we consider simple quenches in which increases uniformly in time, for which, on going from the BCS to the BES regime, we take $\dot {\mathcal{H}}/\mathcal{H}\vert _{{\mathcal{H}}_{0}}=-{{\tau }_{H}}^{-1}$. The time dependence of c_s(t) as changes can now be written as

$$\begin{eqnarray} &&{c}_{\mathrm{s}}^{2}(t)=({c}_{\mathrm{ BCS}}^{2}/2)[1+\tanh (d(g)-t/{\tau }_{\mathrm{ Q}})], \end{eqnarray} \tag{ 19 }$$

where τ_Q ∝ τ_H. t = 0 is the time at which the system is at the unitary limit (i.e. t < 0 in the BCS regime, t > 0 in the BEC regime). If going from the BEC to the BCS regime, we reverse t → −t.

To be concrete, we consider the narrow resonance in ⁶Li at ₀ = 543.25 G, discussed in some detail in [19]. As our benchmark we take the achievable number density ρ₀ ≈ 3 × 10¹² cm⁻³, whence _F ≈ 7 × 10⁻¹¹ eV and γ₀ ≈ 0.2. In terms of the dimensionless coupling $\bar {g}$, where ${g}^{2}=(6 4{\epsilon }_{\mathrm{F}}^{2}/3{k}_{\mathrm{F}}^{3}){\bar {g}}^{2},{\text{}}^{6}\mathrm{Li}$ at the density above corresponds to ${\bar {g}}^{2}\lesssim 1$.

If τ₀ = ħ/_F, the inverse Fermi energy (in units of ħ), then

$$\begin{eqnarray} &&\frac{{\tau }_{\mathrm{Q}}}{{\tau }_{0}}\approx \frac{1}{\dot {\mathcal{H}}}, \end{eqnarray} \tag{ 20 }$$

where $\dot {\mathcal{H}}$ is measured in units of Gauss (ms)⁻¹. Experimentally, it is possible to achieve quench rates as fast as $\dot {\mathcal{H}}\approx 0.1~\mathrm{G}~{\mathrm{ms}}^{-1}$ [19]. As we shall see below, even they are not sufficiently fast to shake off Galilean invariance.

Rather than consider event horizons we continue with the relatively simple case of a homogeneous external field, looking for analogies with Friedmann–Robertson–Walker (FRW) universes in which particle creation can be observed, proposed in [12–14] in particular. Specifically, the analysis in [12–14] is predicated on the representation of the condensate by a Gross–Pitaevskii (GP) mean bosonic field in which the strength of the self-coupling, and hence the scattering length, is given an explicit time dependence.

We have already indicated how, for cold Fermi gases, the linear behaviour of the long wavelength condensate dispersion relation can be derived from a GP equation. We have also anticipated short wavelength non-linear behaviour, crucial for understanding the applicability of (6), that is quantitatively similar to that shown in [12–14]. However, the difference that lies in the detail is sufficient to make what seemed a difficult task in their case an impossible task here.

We stress that, while familiarity makes analogies simpler, there is nothing intrinsically important about mimicking FRW universes, and we begin more generally, with an inhomogeneous background.

4.1. The acoustic metric

On allowing for explicit time dependence in the Lagrangian (15), the Euler–Lagrange equation now takes the form

$$\begin{eqnarray} &&\fl \frac{\mathrm{d}}{\mathrm{d}t}\Biggl [\frac{{\rho }_{0}}{{c}_{\mathrm{s}}^{2}}G(\theta )\Biggr ]+\nabla \cdot \Biggl [\frac{{\rho }_{0}}{2 m{c}_{\mathrm{s}}^{2}}(G(\theta )\nabla \theta )\Biggr ]-{\rho }_{0}{\nabla }^{2}\theta =0 \end{eqnarray} \tag{ 21 }$$

in which c_s is the local speed of sound as before, but now where ρ₀,N₀, etc vary in space and time.

With zero background velocity (21) can be linearized as the geodesic equation

$$\begin{eqnarray} &&\frac{1}{\sqrt{-g}}{\partial }_{\mu }(\sqrt{-g}\hspace{0.167em} {g}^{\mu \nu }{\partial }_{\nu }\tilde {\theta })=0, \end{eqnarray} \tag{ 22 }$$

with which we began, where in d spatial dimensions

$$\begin{eqnarray} {g}_{\mu \nu }=\Biggl (\frac{{\rho }_{0}}{{c}_{\mathrm{s}}}\Biggr )^{\frac{2}{d-1}}\left (\begin{array}{@{}cc@{}} \displaystyle -{c}_{\mathrm{s}}^{2}&\displaystyle 0\\ \displaystyle 0&\displaystyle {\delta }_{i j} \end{array}\right ). \end{eqnarray} \tag{ 23 }$$

4.2. Mimicking spatially homogeneous (and FRW) universes

In this we follow the approach of [12–14]. Consider a homogeneous condensate in a homogeneous magnetic field (t) varying in time. The phonon field θ then satisfies (6) with the metric g(t) of (23), in which ${c}_{\mathrm{s}}^{2}(t)$ of (17) is controlled by of (18). To make the situation simpler from the viewpoint of FRW analogue gravity, we follow the authors of [12–14] in assuming that the system is essentially two-dimensional (a pancake condensate). This permits us to use the given speed of sound in a 2D setting with a direct correspondence with the FRW metric.

As in [12], on taking d = 2 in (23) we can write

$$\begin{eqnarray} {g}_{\mu \nu }=\Biggl (\frac{{\rho }_{0}}{{c}_{\mathrm{BCS}}}\Biggr )^{2}\left (\begin{array}{@{}cc@{}} \displaystyle -{c}_{\mathrm{ BCS}}^{2}&\displaystyle 0\\ \displaystyle 0&\displaystyle a(t)^{2}{\delta }_{i j} \end{array}\right ). \end{eqnarray} \tag{ 24 }$$

If ${c}_{\mathrm{s}}^{2}(t)={c}_{\mathrm{BCS}}^{2}B(t)$ then a²(t) = B(t)⁻¹, where

$$\begin{eqnarray} &&B=\frac{{N}_{0}^{\mathrm{BCS}}}{{N}_{0}+4{\bar {\alpha }}^{2}/{\bar {M}}^{2}}. \end{eqnarray} \tag{ 25 }$$

${N}_{0}^{\mathrm{BSC}}$ is the value of N₀ in the deep BCS regime.

We shall follow this 2D approximation but it is not crucial for any of the analysis. For the simplest uniform quench from the BCS to the BEC regime the time dependence of B(t) as changes is given from (18) as

$$\begin{eqnarray} &&B(t)=(1/2)[1+\tanh (d-t/{\tau }_{\mathrm{Q}})]. \end{eqnarray} \tag{ 26 }$$

The FRW scale factor a(t) is then

$$\begin{eqnarray} &&{a}^{2}(t)=B(t)^{-1}=1+{\mathrm{e}}^{-2 d}{\mathrm{e}}^{2 t/{\tau }_{\mathrm{Q}}}. \end{eqnarray} \tag{ 27 }$$

We have normalized a(t) to unity in the deep BCS regime. To be specific we adopt the ⁶Li parameters that we quoted earlier. Since d < 0 then, insofar as the hydrodynamic approximation is reliable, in the BEC regime (t > 0) we have

$$\begin{eqnarray} &&a(t)\approx {\mathrm{e}}^{-d}{\mathrm{e}}^{t/{\tau }_{\mathrm{Q}}}, \end{eqnarray} \tag{ 28 }$$

corresponding to an effective de Sitter universe with Hubble parameter H = 1/τ_Q. That is, the simplest experimental situation of a constant quench rate initially looks to give one of the most interesting analogue models!

The standard approach to phonon production in condensates is given in great detail in [12]. In driving the system from the BCS to the BEC regime, with our analogue de Sitter expansion, we move from one flat space–time at time t₀ to another at time t. The phonon-free initial ground state evolves into a multi-phonon state with amplitudes derived from the Bogoliubov transformations that relate the initial annihilation and creation operators at t₀ to their later counterparts.

We shall not attempt to perform any calculation of phonon production with respect to the metrics above since, in practice, this simple picture is never implemented due to the breakdown of the simple acoustic model with its SoS Lorentz invariance, as we shall see.

5.1. The breakdown of the simple acoustic model

The action (10) is manifestly Galilean invariant, showing the extent to which the speed of sound Lorentzian geodesic equations (22) are a long wavelength limit. The underlying Galilean invariance manifests itself through non-linearity in the phonon dispersion relation.

The simplest approximation that goes beyond the acoustic approximation is that of a 'rainbow' metric [13, 14], in which phonons with different wavenumbers k travel at different speeds c_k. In this we ignore the time derivatives of the density field in comparison to space derivatives. In our formalism this is very simple. We set η = 0 in (13), which corresponds to approximating D_F(x − x') of (12) by

$$\begin{eqnarray} &&{D}_{\mathrm{F}}(x-{x}^{\prime})=\frac{2}{(2 \pi )^{3}}\delta (t-{t}^{\prime})\int \mathbf{dk}\frac{{\mathrm{e}}^{-\mathrm{i}\mathbf{k}(\mathbf{x}-{\mathbf{x}}^{\prime})}}{{\rho }_{0}{k}^{2}/2 m+{\bar {M}}^{2}}, \end{eqnarray} \tag{ 29 }$$

valid for homogeneous systems in which ρ₀ and ${\bar {M}}^{2}$ depend only on time.

The outcome is that the phonon dispersion relation is

$$\begin{eqnarray} &&{\omega }^{2}={c}_{k}^{2}{k}^{2} \end{eqnarray} \tag{ 30 }$$

where (c_s ≡ c_k=0)

$$\begin{eqnarray} &&{c}_{k}^{2}={c}_{\mathrm{ s}}^{2}\Biggl [1+\frac{{k}^{2}}{{\bar {M}}^{2}}\frac{4{\alpha }^{2}{c}_{\mathrm{s}}^{2}}{{\bar {M}}^{2}}+\cdots \Biggr ]. \end{eqnarray} \tag{ 31 }$$

The restoration of Galilean invariance is more rapid in the BEC regime since ${\bar {M}}^{2}$ becomes vanishingly small in (29). We find

$$\begin{eqnarray} &&{c}_{k}^{2}\approx {c}_{\mathrm{ s}}^{2}\Biggl [1+\frac{{k}^{2}}{{K}^{2}}\Biggr ] \end{eqnarray} \tag{ 32 }$$

where

$$\begin{eqnarray} &&K\sim \frac{2 \alpha {c}_{\mathrm{s}}}{{\rho }_{0}/2 m}\approx 4 m{c}_{\mathrm{s}}. \end{eqnarray} \tag{ 33 }$$

This looks more familiar as the Bogoliubov dispersion relation

$$\begin{eqnarray} &&{\omega }^{2}\approx {c}_{\mathrm{ s}}^{2}{k}^{2}+({k}^{2}/2 M)^{2}. \end{eqnarray} \tag{ 34 }$$

That is, as c_s vanishes, we have ω ≈ k²/2M, describing free diatoms with mass M = 2m.

To see the effects of this on particle production we repeat the analysis of [13, 14]. For modes of wavenumber k the FRW scale a(t) is modified in a mode dependent way from (27) to

$$\begin{eqnarray} &&{a}_{k}^{2}(t)\approx [B(t)+{k}^{2}/{K}_{0}^{2}]^{-1} \end{eqnarray} \tag{ 35 }$$

where ${K}_{0}\approx 4 m{c}_{\mathrm{BCS}}=(4/\sqrt{3}){k}_{\mathrm{F}}$ with k_F the Fermi momentum.

In the BEC regime we find

$$\begin{eqnarray} &&{a}_{k}^{2}(t)\approx \Biggl [{\mathrm{e}}^{2 d}{\mathrm{e}}^{-2 t/{\tau }_{\mathrm{Q}}}+\frac{{k}^{2}}{{K}_{0}^{2}}\Biggr ]^{-1}. \end{eqnarray} \tag{ 36 }$$

From our earlier comments, we take t > 0. Some caution is required in that the approximation (36), derived from (31), breaks down when the second term is too large. The difference between our results and those of [13, 14] essentially lies in the prefactor e^2d < 10⁻¹ for d ≈− 1.4. That is just sufficient to make sure that Galilean invariance is the dominant feature. The detailed demonstration of this is tedious but we have included it for completeness.

The 'rainbow' Hubble parameters in the BEC regime can be obtained as

$$\begin{eqnarray} &&{H}_{k}(t)\equiv \frac{{\dot {a}}_{k}(t)}{{a}_{k}(t)}\approx \frac{1}{{\tau }_{\mathrm{Q}}}\frac{{\mathrm{e}}^{2 d}{\mathrm{e}}^{-2 t/{\tau }_{\mathrm{Q}}}}{[{\mathrm{e}}^{2 d}{\mathrm{e}}^{-2 t/{\tau }_{\mathrm{Q}}}+\frac{{k}^{2}}{{K}_{0}^{2}}]}. \end{eqnarray} \tag{ 37 }$$

The transition when the non-linear dispersion relation becomes important for wavenumber k happens at time t_k,

$$\begin{eqnarray} &&\frac{{t}_{k}}{{\tau }_{\mathrm{Q}}}=d+\ln \frac{{K}_{0}}{k}, \end{eqnarray} \tag{ 38 }$$

in terms of which

$$\begin{eqnarray} &&{H}_{k}(t)\approx \frac{1}{2{\tau }_{\mathrm{Q}}}\frac{{\mathrm{e}}^{-(t-{t}_{k})/{\tau }_{\mathrm{Q}}}}{\cosh ((t-{t}_{k})/{\tau }_{\mathrm{Q}})}. \end{eqnarray} \tag{ 39 }$$

From our comments above this approximation breaks down for t ≳ t_k.

As a guide to phonon production during the sweep from BCS to BEC we also need the modified dispersion relation

$$\begin{eqnarray} &&\fl {\omega }_{k}(t)\approx {\omega }_{0}\Biggl [{\mathrm{e}}^{2 d}{\mathrm{e}}^{-2 t/{\tau }_{\mathrm{Q}}}+\frac{{k}^{2}}{{K}_{0}^{2}}\Biggr ]^{1/2}\nonumber\\ &&\approx \sqrt{2}{\omega }_{0}\Biggl (\frac{k}{{K}_{0}}\Biggr ){\mathrm{e}}^{-(t-{t}_{k})/2{\tau }_{\mathrm{Q}}}\{ \cosh ((t-{t}_{k})/{\tau }_{\mathrm{Q}})\} ^{1/2}, \end{eqnarray} \tag{ 40 }$$

where we have adopted the notation of [13, 14], in which ω₀ = |k|c_BCS. There is an infrared bound on k. For the exemplary condensates with N ≈ 10⁵ atoms, their width is ξ₀ ≈ 10²/k_F [7]. This gives

$$\begin{eqnarray} &&0\lesssim \frac{{t}_{k}}{{\tau }_{\mathrm{Q}}}\lesssim d+\ln 1 0 0, \end{eqnarray} \tag{ 41 }$$

or perhaps a little larger. With d <− 1 this gives

$$\begin{eqnarray} &&0\lesssim {t}_{k}\lesssim 3{\tau }_{\mathrm{Q}}, \end{eqnarray} \tag{ 42 }$$

at best (without having to restrict ourselves to the BEC regime a priori). Taking the example of figure 2 with b ≈ 1.2 this translates into a transition when the non-linear dispersion becomes important for a value of t for which

$$\begin{eqnarray} &&0\lesssim 1/{a}_{\mathrm{S}}{k}_{\mathrm{F}}\lesssim 2.5, \end{eqnarray} \tag{ 43 }$$

as we go from the shortest to the longest wavelengths. We stress that t_k (or the corresponding 1/a_Sk_F) marks the boundary between the applicability of equation (6) with its 'Lorentzian' structure and the restoration of Galilean invariance.

The relevant quantity is the ratio

$$\begin{eqnarray} &&{\mathcal{R}}_{k}(t)=\frac{{\omega }_{k}(t)}{{H}_{k}(t)}. \end{eqnarray} \tag{ 44 }$$

A quantum mode with wavenumber k only experiences significant amplification (and hence phonon production) when _k(t) ≪ 1. As before, this approximation breaks down when t ≳ t_k.

In the vicinity of t_k, where H_k is small, _k is correspondingly large, as it is for t much greater than t_k. In between it achieves its minimum

$$\begin{eqnarray} &&{\mathcal{R}}_{k}({t}^{\ast })\approx \frac{{\tau }_{\mathrm{Q}}}{{\tau }_{0}}\Biggl (\frac{k}{{k}_{\mathrm{F}}}\Biggr )^{2}. \end{eqnarray} \tag{ 45 }$$

In order to have any phonon production we must have as fast a quench as possible, with a current lower bound of τ_Q/τ₀ ≈ 10 and a lower bound of k/k_F of 10⁻², say, for our typical condensate. Then, for the lowest momentum phonons,

$$\begin{eqnarray} &&{\mathcal{R}}_{k}({t}^{\ast })\approx 1{0}^{-3}, \end{eqnarray} \tag{ 46 }$$

this minimum increasing as momentum increases. Thus, from (44), there is a window in which _k(t) is sufficiently small to expect phonon production. However, even then the number is negligible. For such low momentum phonons to be produced within the timescale t_k, the corresponding number density for each k mode can be approximated by

$$\begin{eqnarray} &&{\mathcal{N}}_{k}\approx {\mathrm{e}}^{\int \nolimits _{0}^{{t}_{k}}{H}_{k}(t)\mathrm{d}t}\lt {\mathrm{e}}^{{H}_{k}({t}^{\ast }){t}_{k}}\approx {\mathrm{e}}^{0.6(d+\ln \frac{{K}_{0}}{k})}, \end{eqnarray} \tag{ 47 }$$

where we have used (38) and (39). For a lower bound of the momentum, say k = 10⁻²k_F, again with d ≈− 1.4 and ${K}_{0}=(4/\sqrt{3}){k}_{\mathrm{F}}$, the corresponding  ≈ 1. For higher momentum phonons and somewhat slower quenches _k(t) ≫ 1 throughout and there is no phonon production. Even when they occur, the effective temperature of the few phonons produced by changing the metric is O(ε_F) (in units in which k_B = 1), much lower than that of the gas. As a result such phonons will be swamped by thermal phonons.

As for reversing the direction of the quench, we have seen that the non-linear effects are greater as c_s becomes smaller, making a quench beginning in the deep BEC regime problematical.

Either way, a cold Fermi gas is not a helpful system for exploring parallels with early universe cosmology, if by that we mean that phonon production is determined from the geodesic equation (6), despite the optimism of [12–14] for general condensates. Simulations for coupled bosonic binary condensates, which also have a gapless and a gapped mode, suggest that they are better candidates for gravitational analogies [20, 21], when details are taken into account.

However, instead of requiring the geodesic equation to be relevant, we can ask the general question of whether we can predict and observe phonon production beyond the acoustic approximation. This may be possible by applying an oscillating external field, corresponding to the metric of a cyclic universe [12], were the acoustic approximation relevant. This is an interesting case in that, with the final and initial states identical, particle production occurs only as a result of parametric excitation.

A simple choice is to take

$$\begin{eqnarray} &&\frac{1}{{k}_{\mathrm{F}}{a}_{\mathrm{S}}}=\frac{1}{{k}_{\mathrm{F}}{a}_{\mathrm{S}}^{0}}(1-\sin \omega t), \end{eqnarray} \tag{ 48 }$$

corresponding from (18) to an oscillating field in which the oscillation of ${a}_{\mathrm{S}}^{-1}$ about $({a}_{\mathrm{S}}^{0})^{-1}$ extends to the unitary regime.

As a result

$$\begin{eqnarray} {a}^{2}(t)&=&1+{a}_{0}^{2}\hspace{0.167em} {\mathrm{e}}^{-2 A\sin \omega t}\nonumber\\ &\approx &{a}_{0}^{2}\hspace{0.167em} {\mathrm{e}}^{-2 A\sin \omega t}, \end{eqnarray} \tag{ 49 }$$

where $A=b/{k}_{\mathrm{F}}{a}_{\mathrm{S}}^{0}$ and ${a}_{0}^{2}=\exp -2(d-b/{k}_{\mathrm{F}}{a}_{\mathrm{S}}^{0})$.

The corresponding rainbow metric scale factor is

$$\begin{eqnarray} &&{a}_{k}^{2}(t)\approx \Biggl [{a}_{0}^{-2}\hspace{0.167em} {\mathrm{e}}^{2 A\sin \omega t}+\frac{{k}^{2}}{{K}_{0}^{2}}\Biggr ]^{-1}, \end{eqnarray} \tag{ 50 }$$

from which the rainbow Hubble parameters follow as

$$\begin{eqnarray} &&{H}_{k}(t)\approx -A \omega \cos \omega t \hspace{0.167em} \frac{[-{a}_{0}^{-2}\hspace{0.167em} {\mathrm{e}}^{2 A\sin \omega t}+\frac{{k}^{2}}{{K}_{0}^{2}}]}{[{a}_{0}^{-2}\hspace{0.167em} {\mathrm{e}}^{2 A\sin \omega t}+\frac{{k}^{2}}{{K}_{0}^{2}}]}. \end{eqnarray} \tag{ 51 }$$

The type of prediction we can make is not yet clear. This is under consideration and we shall not pursue it further here, beyond noting that for our parameters the bulk of any phonon production is controlled again by the Galilean group.

For a completely different test of analogies with cold Fermi gases, we observe that a further consequence of this splitting of time derivatives in the action, when it arises, is to induce fluctuations in the time of flight of phonons or sound waves. Whether or not this can be cast in the language of stochastic fluctuations as in [22, 23] is being analysed at the moment.

We have continued our exploration of the properties of cold Fermi gases which can be tuned through a narrow Feshbach resonance. We had shown in previous work [7] that the system can be represented by a coupled two-field model of a gapless Goldstone phonon and a gapped Higgs mode. Here we have gone a step further to trace out the Higgs to give an effective purely phononic theory. In particular, the non-linearity of the phonon dispersion relation which can be calculated explicitly.

In choosing to work with tunable condensates, we had hoped that they would provide a straightforward system for establishing parallels with the causal properties of the very early universe. In particular we are interested in the spontaneous creation of particles (phonons) because of the rapid variation of the speed of sound Lorentzian acoustic metric which defines their geodesics.

The first observation is that, for there to be any measurable effects, the speed of sound must change fast. This happens best in the BEC regime. However, the parallels with the early universe are complicated, for very different reasons. For phonon production our ability to construct analogue FRW universes, for example, is compromised by the restoration of Galilean invariance at the expense of the SoS Lorentzian behaviour of (6), a counterpart [3, 4] to the possible (but unobserved) Lorentzian breakdown of the early universe at the Planck length. However, whereas the Planck energy is sufficiently high that its effects can usually be ignored, this restoration of Galilean invariance in the BEC regime is sufficiently rapid to make the predictions from the acoustic metric unreliable. The effect is that any phonon production has little, if anything, to do with the geodesic equation (6) with which we began this paper. Nonetheless, if we treat phonon production as an interesting question in its own right, without looking for analogies, we may still be able to make useful predictions. At the moment these are largely questions of principle since, in our work, we have only considered homogeneous condensates. We should implement trapping, otherwise the condensate will split apart.

We conclude with a comment on the relation of our work, both here and in our papers [5–8], to that of the condensate mainstream, often canonically Hamiltonian based. Apart from ourselves several authors have used a path integral approach, both without an explicit resonance (e.g. [10, 24]) and with a resonance (e.g. [18]), but none as far as we know in the same detail when a resonance is present. A key ingredient in our analysis is our emphasis on step by step Galilean invariance, without which we could not create our fluid models. In fact, the importance of preserving the underlying symmetry in all approximations, which goes back at least to [9], is very clear in discussions of the unitary regime where explicit Galilean invariance is crucial [7, 11]. Unfortunately, Galilean invariance is not always manifest in condensed matter papers, particularly from a canonical viewpoint, which makes comparison difficult. However, whenever we have compared our results to those obtained by other methods, we are in agreement, e.g., the dependence of the sound speed on external field (see [5]), the nature of damped oscillations (see [7]), and the strength of condensate interactions (again see [5]).

We thank Dani Steer of APC, University of Paris Diderot, Marek Tylutki of Jagiellon University, Krakov, Piyush Jain of UBC, Vancouver and Silke Weinfurtner of SISSA, Trieste, for helpful discussions. The work of JTH, DSL and CYL was supported in part by the National Science Council, Taiwan.

The explicit Fermi density is

$$\begin{eqnarray*} &&{\rho }_{0}^{\mathrm{F}}=\int {\mathrm{d}}^{3}\mathbf{p}/(2 \pi )^{3}~\left [1-{\varepsilon }_{\mathrm{ p}}/{E}_{\mathrm{p}}\right ] \end{eqnarray*}$$

where, in conventional notation, ε_k = k²/2m and ${E}_{\mathrm{p}}=({\varepsilon }_{\mathrm{p}}^{2}+{g}^{2}\vert {\phi }_{0}\vert ^{2})^{1/2}$.

N₀ is the density of states at the Fermi surface, but the other relevant coefficients have no immediate interpretation. After renormalization [5] they take the form

$$\begin{eqnarray} {N}_{0}&=&{g}^{2}\vert {\phi }_{0}\vert ^{2}\int \frac{{\mathrm{d}}^{3}\mathbf{p}}{(2 \pi )^{3}}\frac{1}{2{E}_{\mathrm{p}}^{3}}, \end{eqnarray} \tag{ 52 }$$

$$\begin{eqnarray} \alpha &=&2\vert {\phi }_{0}\vert {\kappa }^{1/2}\Biggl [1+\frac{1}{2}{g}^{2}\int \frac{{\mathrm{d}}^{3}\mathbf{p}}{(2 \pi )^{3}}\frac{{\varepsilon }_{\mathrm{p}}}{2{E}_{\mathrm{p}}^{3}}\Biggr ], \end{eqnarray} \tag{ 53 }$$

$$\begin{eqnarray} {\bar {M}}^{2}&=&-\frac{{g}^{2}\kappa m}{\pi {a}_{\mathrm{S}}}-2{g}^{2}\kappa \int \frac{{\mathrm{d}}^{3}\mathbf{p}}{(2 \pi )^{3}}\Biggl [\frac{{\varepsilon }_{\mathrm{p}}^{2}}{{E}_{\mathrm{p}}^{3}}-\frac{1}{({\mathbf{p}}^{2}/2 m)}\Biggr ], \end{eqnarray} \tag{ 54 }$$

$$\begin{eqnarray} \eta &=&{g}^{2}\kappa \int \frac{{\mathrm{d}}^{3}\mathbf{p}}{(2 \pi )^{3}}\frac{{\varepsilon }_{\mathrm{p}}^{2}}{2{E}_{\mathrm{p}}^{5}}, \end{eqnarray} \tag{ 55 }$$

and

$$\begin{eqnarray} \zeta &=&\int \frac{{\mathrm{d}}^{3}\mathbf{p}}{(2 \pi )^{3}}\Biggl [\frac{1}{8{E}_{\mathrm{p}}^{3}}\Biggl [\Biggl (1-3\frac{{g}^{2}\vert {\phi }_{0}\vert ^{2}}{{E}_{\mathrm{p}}^{2}}\Biggr )\frac{{\varepsilon }_{\mathrm{p}}}{m}\nonumber\\ &&+~\Biggl (5\frac{{g}^{2}\vert {\phi }_{0}\vert ^{2}}{{E}_{\mathrm{p}}^{2}}\Biggl (1-\frac{{g}^{2}\vert {\phi }_{0}\vert ^{2}}{{E}_{\mathrm{p}}^{2}}\Biggr )\Biggr )\frac{\vert \mathbf{p}\vert ^{2}\hspace{0.167em} \hat {\mathbf{p}}\boldsymbol{ \cdot }\hat {\nabla }}{{m}^{2}}\Biggr ]\Biggr ]. \end{eqnarray} \tag{ 56 }$$

In (56) $\hat {\mathbf{p}}$ and $\hat {\nabla }$ are the unit vectors along the direction p and the direction of the spatial variation of the phase mode θ respectively.

In terms of ζ the scale factor is

$$\begin{eqnarray} &&\kappa =\frac{{\rho }_{0}}{4 m{g}^{2}\zeta +2}. \end{eqnarray} \tag{ 57 }$$

In our approximations η is not particularly relevant, but suffice to say that η = N₀ in the deep BCS regime (when α is small) and falls off in the same way as N₀ in the deep BEC regime (when α ≈ ρ₀).

Finally, in the definition of ${\bar {M}}^{2}$, we have used the relationship between the s-wave scattering length a_S and the binding energy,

$$\begin{eqnarray} &&2 \mu -\nu =\frac{{g}^{2}m}{4 \pi {a}_{\mathrm{S}}}. \end{eqnarray} \tag{ 58 }$$

N₀ and α depend indirectly on a_S through the constancy of ρ₀.