The BCS-BEC crossover: From ultra-cold Fermi gases to nuclear systems

# The BCS-BEC crossover: From ultra-cold Fermi gases to nuclear systems

## Abstract

This report aims at covering the overlap between ultra-cold Fermi gases and nuclear matter in the context of the BCS-BEC crossover. By this crossover, the phenomena of Bardeen-Cooper-Schrieffer (BCS) superfluidity and Bose-Einstein condensation (BEC), which share the same kind of spontaneous symmetry breaking, are smoothly connected through the progressive reduction of the size of the fermion pairs involved as the fundamental entities in both phenomena. This size ranges, from large values when Cooper pairs are strongly overlapping in the BCS limit of a weak inter-particle attraction, to small values when composite bosons are non-overlapping in the BEC limit of a strong inter-particle attraction, across the intermediate unitarity limit where the size of the pairs is comparable with the average inter-particle distance.

The BCS-BEC crossover has recently been realized experimentally, and essentially in all of its aspects, with ultra-cold Fermi gases. This realization, in turn, has raised the interest of the nuclear physics community in the crossover problem, since it represents an unprecedented tool to test fundamental and unanswered questions of nuclear many-body theory. Here, we focus on the several aspects of the BCS-BEC crossover, which are of broad joint interest to both ultra-cold Fermi gases and nuclear matter, and which will likely help to solve in the future some open problems in nuclear physics (concerning, for instance, neutron stars). Similarities and differences occurring in ultra-cold Fermi gases and nuclear matter will then be emphasized, not only about the relative phenomenologies but also about the theoretical approaches to be used in the two contexts. Common to both contexts is the fact that at zero temperature the BCS-BEC crossover can be described at the mean-field level with reasonable accuracy. At finite temperature, on the other hand, inclusion of pairing fluctuations beyond mean field represents an essential ingredient of the theory, especially in the normal phase where they account for precursor pairing effects.

After an introduction to present the key concepts of the BCS-BEC crossover, this report discusses the mean-field treatment of the superfluid phase, both for homogeneous and inhomogeneous systems, as well as for symmetric (spin- or isospin-balanced) and asymmetric (spin- or isospin-imbalanced) matter. Pairing fluctuations in the normal phase are then considered, with their manifestations in thermodynamic and dynamic quantities. The last two Sections provide a more specialized discussion of the BCS-BEC crossover in ultra-cold Fermi gases and nuclear matter, respectively. The separate discussion in the two contexts aims at cross communicating to both communities topics and aspects which, albeit arising in one of the two fields, share a strong common interest.

## 1 Introduction

### 1.1 Historical background

The idea behind the BCS-BEC crossover dates back just after the birth of the BCS theory in 1957 Bardeen-1957 (). The authors of this theory made a point to emphasize the differences between their theory for superconductors based on strongly-overlapping Cooper pairs and the Schafroth-Butler-Blatt theory Schafroth-1957 () resting on non-overlapping composite bosons which undergo Bose-Einstein condensation at low temperature. Subsequently, the interest in these two different situations has been kept disjoint for some time, until theoretical interest arose for unifying them as two limiting (BCS and BEC) situations of a single theory where they share the same kind of broken symmetry. In this way, the physical system passes through an intermediate situation, where the size of the pairs is comparable to the average inter-particle spacing. (For a short-range (contact) interaction, this intermediate situation is nowadays referred to as the “unitary” regime.) Pioneering work in this sense was done by the Russian school, motivated by the exciton condensation in semiconductors Keldysh-1965 () or simply by intellectual curiosity Popov-1966 (). The theory of the BCS-BEC crossover took shape initially through the work by Eagles Eagles-1969 () with possible applications to superconducting semiconductors, and later through the works by Leggett Leggett-1980 () and Nozières and Schmitt-Rink Nozieres-1985 () where the formal aspects of the theory were developed at zero temperature and above the critical temperature, respectively.

The interest in the BCS-BEC crossover grew up with the advent of high-temperature (cuprate) superconductors in 1987, in which the size of the pairs appears to be comparable to the inter-particle spacing Randeria-1989 (); Randeria-1990 (); Micnas-1990 (); Randeria-1992 (); Drechsler-1992 (); Haussmann-1993 (); Pistolesi-1994 (); Casas-1994 (). Related interest in the BCS-BEC crossover soon spread to some problems in nuclear physics Baldo-1995 (), but a real explosion of this activity appeared starting from 2003 with the advent of the fully controlled experimental realization essentially of all aspects of the BCS-BEC crossover in ultra-cold Fermi gases (see Refs. Inguscio-2007 () and Zwerger-2012 () for an experimental and theoretical overview, respectively). This realization, in turn, has raised the interest in the crossover problem especially of the nuclear physics community, as representing an unprecedented tool to test fundamental and unanswered questions of nuclear many-body theory. The Fermi gas at the unitary limit (UL), where fermions of opposite spins interact via a contact interaction with infinite scattering length, was actually introduced as a simplified model of dilute neutron matter Baker-2000 (); Baker-1999 (), and the possibility to realize this limit with ultra-cold atoms was hence regarded as extremely important for this field of nuclear physics.

In nuclear physics there is the paradigmatic example of a proton-neutron bound state (the deuteron) for which in symmetric or asymmetric nuclear matter, as a function of density, one may find a transition from BEC to BCS Alm-1993 (); Baldo-1995 (); Stein-1995 (); Lombardo-2001b () (see also Ref.Andrenacci-1999 () for the density-induced BCS-BEC crossover). Such a scenario may be realized in expanding nuclear matter generated from heavy ion collisions Baldo-1995 () or in proto-neutron stars Heckel-2009 (). More recently, this density driven crossover has been studied together with the competing liquid-gas phase transition Jin-2010 (). Deuteron condensation also heavily competes with alpha-particle condensation, which is somehow related to pairing of pairs discussed also in condensed-matter physics. Nuclear physics is a precursor of the theory for quartet condensation Roepke-1998 (), but theoretical studies of quartets came up later also in the area of ultra-cold atoms Capponi-2007 (). Alpha-particle condensation is presently very much discussed in finite nuclei (Hoyle state) Yamada-2012 (). It was also predicted that in the non-condensed phase, the deuterons give rise to a pseudo-gap formation Schnell-1999 (). With today’s ultra-cold atoms experiments, it has become possible to test such theories Stewart-2008 (); Gaebler-2010 (); Feld-2011 (). For neutrons, no bound state exists, but rather a virtual state at almost zero energy. As a consequence, a dilute gas of neutrons, as it exists in the inner crust of neutron stars, is almost in the unitary limit mentioned above. Recent studies of the dilute neutron gas, from (almost) unitarity at low density to the BCS limit at high density, were done within BCS theory Matsuo-2006 (), Quantum-Monte-Carlo calculations Abe-2009 (); Gezerlis-2010 (), and the Nozières-Schmitt-Rink approach Ramanan-2013 (). One can say that the equation of state and pairing properties of very dilute neutron matter, although inaccessible in experiments, are now known thanks to the analogy with ultra-cold atoms.

As these examples show, there are several aspects of the BCS-BEC crossover which are of broad joint interest to both ultra-cold atoms and nuclear communities. This paper is thus meant to provide a comprehensive review which focuses mainly on these common aspects of ultra-cold atoms and nuclear physics. Along these views, this paper provides also a pedagogical review of the main essential aspects of the BCS-BEC crossover. In the process of writing, we have mostly adopted the quantum many-body diagrammatic techniques in line with our own technical expertise, and we have mainly focused on the topics to which we have provided original contributions over the last several years. Accordingly, for readers interested in complementary theoretical approaches to the problem we refer to other reviews which cover various aspects of the BCS-BEC crossover. In particular, we can refer to reviews on the application to strongly-interacting Fermi gases (in particular, at unitarity) of Quantum Monte Carlo methods Bulgac-2012 (); Carlson-2012 (), functional-renormalization-group techniques Diehl-2010 (), epsilon Nishida-2012 () and virial Liu-2013 () expansions. For the application of functional-integral approaches (in particular, to the superfluid phase) we refer instead to the original research works of Refs. Pistolesi-1996 (); Hu-2006 (); Diener-2008 ().

In addition, it should be mentioned that further reviews cover a number of aspects on the BCS-BEC crossover which share a partial overlap with the material discussed here. Specifically, the unitary limit of Section 4.3 and the Tan contact of Section 4.5 have been of most interest in other reviews owing to the widespread recent interest in these topics, which have been treated in Refs. Zwerger-2012 (); Randeria-2014 (); Pitaevskii-2016 (); Zwerger-2016 (). The Fano-Feshbach resonances of Section 4.1, which are at the heart of the interaction-induced crossover, have been discussed in Refs. Chin-2010 (); Gubbels-2013 (); Randeria-2014 (); Pitaevskii-2016 (); Zwerger-2016 (). Polarized Fermi gases (considered here in Sections 2.6 and 3.7) have also been treated in Refs. Giorgini-2008 (); Radzihovsky-2010 (); Chevy-2010 (); Zwerger-2012 (); Gubbels-2013 (); Pitaevskii-2016 (). The topic of the single-particle spectral function and pseudo-gap of Section 3.3 can be found discussed also in Refs. Chen-2005 (); Randeria-2014 (); Zwerger-2016 (). Also some aspects of pairing in nuclear systems discussed in Section 5 can be found in the reviews of Refs. Brink-2005 (); Dean-2003 (); Gezerlis-2014 (). Finally, Refs. Pitaevskii-2016 () and Gezerlis-2014 () cover also some introductory material treated here in Sections 1.2, 2.1, 2.4, 3.4, 4.7, 4.8, and 5.1. However, it should be remarked that no other reference thus far has emphasized the aspects of the BCS-BEC crossover common to ultra-cold atoms and nuclear physics as we have done here.

### 1.2 The BCS wave function and its BEC limit

The starting point is the BCS ground-state wave function, of the form Bardeen-1957 (); Schrieffer-1964 ()

 |ΦBCS⟩=∏k(uk+\varvkc†k↑c†−k↓)|0⟩. (1)

In this expression, is the vacuum state, is a fermionic creation operator for wave vector and spin projection , and and are probability amplitudes given by:

 \varv2k=1−u2k=12(1−ξkEk) (2)

where ( being the fermion mass and the chemical potential) and where is the BCS order parameter (sometimes referred to as the superconducting gap) here taken at zero temperature.  For simplicity, we have assumed that is independent of the wave vector, as it is the case for a contact interaction. More generally, in the expression of will be replaced by the wave-vector (and temperature) dependent value . [Throughout this paper, we use units where the reduced Planck constant and the Boltzmann constant are set equal to unity.]

The BCS wave function (1) has the important property that it is the vacuum to the so-called quasi-particle operators (see Eq. (28) below and Ref. Ring-1980 ()), that is, . This relation facilitates considerably the evaluation of expectation values. One can readily show that is the occupation number which goes over to the Fermi step in the limit as tends to zero, and that is the so-called anomalous density (known also as the “pairing tensor” in nuclear physics) which characterizes the BCS wave function.

With reference to the BCS-BEC crossover, it has long been known that the BCS wave function (1) contains the Bose-Einstein condensation of composite bosons as a limit. This is because, upon setting , the expression (1) can be rewritten in the form (see, e.g., Ref. Ring-1980 ()):

 |ΦBCS⟩=(∏k′uk′)exp[∑kgk§c†k↑c†−k↓]|0⟩ (3)

since owing to Pauli principle. Here, the operator contains fermion pairs but it is not a truly bosonic operator, to the extent that the commutator is not a c-number but explicitly contains the fermionic operators . However, under some circumstances one may consider that for all practical purposes, provided for all of physical relevance. As a consequence, represents a bosonic coherent state (that is, a Bose-Einstein condensate) with a non-vanishing broken-symmetry average .

It is clear from Eq. (2) that the condition can be satisfied for all provided the fermionic chemical potential becomes large and negative. This condition can be achieved when a bound-state with binding energy occurs for the two-body problem in vacuum with positive scattering length , and the coupling parameter becomes large such that . In this limit, approaches the value , which amounts to saying that all fermions are paired up in tight (composite) bosons with a vanishing residual interaction among the bosons. This result for the BEC limit of the fermionic chemical potential can directly be obtained from the mean-field gap equation (at zero temperature), which in the case of a wave-vector dependent interaction reads:

 Δk=−∫dk′(2π)3Veff(k,k′)Δk′2Ek′ (4)

where now . (In the case of a contact potential, and reduce, respectively, to the coupling constant of Eq. (7) and to the constant gap given by the gap equation (8) below.) Using , , and Eq. (2), one can rewrite Eq. (4) in the form:

 2ξkϕk+(1−2nk)∫dk′(2π)3Veff(k,k′) ϕk′=0. (5)

Provided for all , Eq. (5) is just the Schrödinger equation for the relative motion of two particles of equal mass which are mutually interacting via the potential . The negative eigenvalue of this equation thus corresponds to (minus) the two-body binding energy as stated above.

More generally, as we shall see below, it is the fermionic chemical potential that takes the key role of the driving field which enables the system to pass from the BCS to the BEC limits of the BCS-BEC crossover. As an illustration, Fig. 1 shows the occupation number at zero temperature for different (from positive to negative) values of the chemical potential , which correspond to increasing values of the gap parameter . Note that only when the curves have an inflection point at , which highlights the presence of an underlying Fermi surface even for a system with attractive inter-particle interaction. As a consequence, when becomes negative, the Fermi sea gets completely dissolved and the occupation number becomes quite small for all .

### 1.3 Pairing correlations

The BCS wave function (1), or its equivalent form (3), treats all fermion pairs on the same footing to the extent that a single wave function is assigned to each pair. This mean-field type of approach is appropriate to describe a system of fermions with a mutual attractive interaction when the inter-particle correlations extend much beyond the average inter-particle distance, in such a way that many pairs are contained within the size of a given pair and different pairs strongly overlap with each other. This is definitely the case for the (BCS) weak-coupling limit, to which the BCS theory of superconductivity was originally meant to apply Bardeen-1957 (). When considering the BCS-BEC crossover, however, the range of the inter-particle correlations can decrease down to the size of the bound pair, which, in turn, can be much smaller than the inter-particle distance. Under these circumstances, the gas of composite bosons becomes dilute and one accordingly anticipates that pairing fluctuations beyond mean field can acquire a major role, especially at finite temperature when they are accompanied by thermal fluctuations.

Figure 2 shows the evolution throughout the BCS-BEC crossover of the zero-temperature pair coherence length (dashed line) calculated for the BCS ground state (1). This quantity provides information about the intra-pair correlation established between fermions of opposite spins Pistolesi-1994 () and, at the mean-field level, corresponds to the spatial extension of the order parameter . Since can, in turn, be interpreted as a wave function as shown by Eq. (5), the r.m.s. radius of the pairs can be obtained from the expression:

 ξ2pair=∫drr2|ϕ(r)|2∫dr|ϕ(r)|2 (6)

where is the Fourier transform of . [For later convenience, Fig. 2 shows also the corresponding evolution of the healing length (full line), that instead provides information about the inter-pair correlation Pistolesi-1996 () and requires the inclusion of pairing fluctuations beyond mean field - see Section 3.2]. In Fig. 2, the average inter-particle distance is used as the unit for , where is the Fermi wave vector related to the particle density via the relation . In Fig. 2 and throughout this paper, the quantity plays the role of the coupling parameter of the theory. Depending on the sign of , this parameter ranges from characteristic of the weak-coupling BCS regime when , to characteristic of the strong-coupling BEC regime when , across the value at unitarity when diverges. In practice, the “crossover region” of most interest, to be discussed throughout this paper, turns out to be approximately limited to the interval .

Pairing fluctuations play a particularly important role in the normal phase above the critical temperature where the order parameter vanishes (the generalization of the BCS approach to finite temperature will be considered in Section 2). This is because, a “local” order is expected to survive above if the system is fluctuating, even though the long-range order characteristic of the superfluid phase is lost. These considerations have led people to associate pairing fluctuations with the occurrence of a pseudo-gap above (as seen, for instance, by a depression of the single-particle density of states about the chemical potential), in a similar way to what occurs in the superfluid phase below when the order parameter is instead nonvanishing. Such an analogy (between the order parameter below and the pseudo-gap above ) has been carried to the point that even the pseudo-gap state occurring in copper-oxide superconductors has been interpreted from the point of view of a BCS-BEC crossover scenario Chen-2005 ().

In this context, it can be useful to consider the similarity between the pseudo-gap physics resulting from pairing fluctuations above and the persistence of (damped) spin waves, which are present over regions of limited extent in the normal phase of ferromagnetic (or antiferromagnetic) materials when a strict long-range order is absent Mook-1973 (). Through this analogy, pseudo-gap phenomena in a Fermi system with an attractive inter-particle interaction are attributed to the persistence of a “local pairing order” above the superfluid temperature , which occurs even though the (off-diagonal) long-range order is absent. This local order, which is built up by pairing fluctuations above , makes the single-particle excitation spectrum to resemble that below , although with a frequency broadening due to the decay of the local excitations which cannot propagate over long distances in the absence of long-range order.

In addition, it is relevant to point out that, when raising the temperature and approaching from below, pairing correlations turn out to persist over a finite distance even within mean field, without the need of considering pairing fluctuations beyond mean field. Indeed, explicit calculations of the temperature dependence of the pair correlation function within mean field show that maintains a finite value at Marsiglio-1990 (), a result which remains true irrespective of coupling across the BCS-BEC crossover Palestini-2014 ()). This remark, too, points out the importance of including pairing fluctuations in the normal phase above , in order to obtain a meaningful finite value of even when approaching from above. Otherwise, would be discontinuous when passing from to , since above fermions become non-interacting within the BCS (mean-field) approximation. We shall return to this point more extensively in Section 3.2.

## 2 BCS mean field

### 2.1 The homogeneous infinite matter case for a contact interaction

The BCS theory of superconductivity considers an effective attractive interaction between otherwise free fermions of two different species (conventionally referred to as spin and ), which are embedded in a continuous medium at a distance Bardeen-1957 (); Schrieffer-1964 (). This attractive interaction, albeit weak, is responsible for the formation of Cooper pairs in the medium Cooper-1956 (). In condensed matter (like in metallic superconductors), detailed knowledge of the form of the attractive interaction is not required for most purposes and one may accordingly consider the simple form of a “contact” (zero-range) potential , where is a negative constant. This model interaction fully applies to ultra-cold Fermi gases (at least in the presence of a broad Fano-Feshbach resonance as discussed in Section 4.1), but does not apply to nuclear systems discussed in Section 5 for which a finite-range interaction should in principle be retained. In order to reduce the numerical difficulties, however, even in the treatment of nuclear superfluids the finite-range interaction is often approximated by a zero-range one, which requires one to introduce a density-dependent coupling constant (see, e.g., Ref. Garrido-2001 ()). This approach has further been elaborated (for instance, in Refs. Bulgac-2002 (); Yu-2003 ()) with the introduction of an efficient regularization scheme for finite systems (as discussed in Section 2.4). Here, like in most part of the paper, we consider the case of equal populations , while the new features arising under the more general condition will be explicitly considered in Section 2.6.

The price to pay for the use of a contact interaction is that, when dealing with a homogeneous system, integrals over the wave vector may diverge in the ultraviolet since the Fourier transform is a constant (in what follows, we shall consider a three-dimensional system). This difficulty can be simply overcome by introducing a cut-off, given for metallic superconductors by the Debye frequency and for nuclei by a phenomenologically adjusted cut-off energy that reflects the range of the effective force. For ultra-cold gases, on the other hand, one can exploit the fact that a similar divergence affects also the two-body problem in vacuum, whereby the fermionic scattering length is obtained from the relation Sa-de-Melo-1993 ():

 m4πaF=1\varv0+∫|k|≤k0dk(2π)3mk2 (7)

with an ultraviolet cutoff . This regularization procedure entails the limits and to be taken simultaneously, in such a way that is kept at the desired value. Replacing the parameter with the physical quantity is especially relevant in the context of ultra-cold gases, where the natural cut-off given by the inverse of the interaction range would be orders of magnitude larger than all wave vectors of physical interest and where can be experimentally controlled Regal-2003a (). As we shall see in Section 3.1, the regularization (7) is also of help when dealing with the many-body problem based on this two-body interaction, to the extent that it greatly reduces the number of many-body diagrams that survive in the limit .

In this way, the BCS mean-field equations for the order (gap) parameter at the temperature and for the density become Schrieffer-1964 ():

 −m4πaF=∫dk(2π)3(1−2f(Ek)2Ek−mk2) (8) n=∫dk(2π)3(1−ξkEk(1−2f(Ek))) (9)

where is the Fermi function and . Note that the quantity that enters Eqs. (8) and (9) plays the dual role of the order parameter (which is non-vanishing between and ) and of the minimum value of the single-particle excitation energy at when . [The extension to finite temperature of the gap equation (4) with a finite-range interaction will be considered in Section 5.1 - cf. Eq. (115).]

For a weak attractive interaction, is small and negative such that . This limit characterizes what is called the “conventional” BCS theory Bardeen-1957 (). In this limit, the Fermi surface of non-interacting fermions is only slightly perturbed by the presence of the interaction, and the chemical potential at coincides with the Fermi energy . This is actually the limit where the mean-field approximation is expected to work best, to the extent that the size of a Cooper pair is much larger than the average inter-particle distance and a large number of pairs is contained within the size of Schrieffer-1964 ().

### 2.2 Solution at T=0 for a contact potential

For a generic value of the coupling that ranges between the BCS and BEC limits, the values of and at are determined by solving the coupled equations (8) and (9), where one sets since is always positive. Although a numerical solution of the ensuing equations is possible Carter-1995 (), it was found independently that the three-dimensional integrals over the wave vector occurring in these equations can be expressed either in terms of the complete elliptic integrals Marini-1998 () or in terms of the Legendre function Papenbrock-1999 (), in such a way that closed-form expressions for and are obtained. In Ref. Marini-1998 (), the results were expressed in terms of the parameter that ranges from (BEC limit) to (BCS limit), while in Ref. Papenbrock-1999 () the inverse of this parameter was used.

In Ref. Marini-1998 (), the following results were obtained:

 Δ0EF=1(x0J(x0)+I(x0))2/3, (10) μ0EF=x0Δ0EF, (11) 1kFaF=−4πx0I(x0)−J(x0)(x0J(x0)+I(x0))1/3, (12)

where

 I(x0)=12(1+x20)1/4F(π2,κ) (13)

and

 J(x0)=(1+x20)1/4E(π2,κ)−14y0(1+x20)1/4F(π2,κ). (14)

In these expressions, and are the complete elliptic integrals of the first and second kind, respectively, and

 y0=√1+x20+x02,κ2=y0√1+x20. (15)

In Ref. Papenbrock-1999 (), on the other hand, Eq. (12) was equivalently expressed in the form:

 1kFaF=(1+x20x20)1/4P1/2(−|x0|(1+x20)−1/2) (16)

where denotes the Legendre function. In addition, analytic expressions for the pair coherence length and the condensate fraction were obtained in Refs. Marini-1998 () and Salasnich-2005 (), respectively.

Analytic approximations can be readily obtained in the BCS and BEC limits. In the BCS limit where , one obtains in Eq.(15). One can thus approximate , , and , in such a way that:

 Δ0=8EFe2exp(π2kFaF). (17)

In the BEC limit, on the other hand, and such that . One thus finds:

 1kFaF≃(163π)1/3|x0|2/3, (18) Δ0EF≃(163π)2/3|x0|1/3≃(163π)1/2(1kFaF)1/2, (19) μ0EF≃−(163π)2/3|x0|4/3≃−12ε0EF. (20)

The full expressions (10), (11), and (12) interpolate smoothly between the above (BCS and BEC) limits through the crossover region . Plots of and as functions of are readily obtained from these equations, as shown in Figs. 3(a) and 3(b). Note, in particular, that the chemical potential: (i) In the BCS limit equals the Fermi energy ; (ii) In the crossover region gradually decreases and eventually changes its sign; (iii) In the BEC limit reaches the asymptotic value . It is just this behaviour of the chemical potential that drives the BCS-BEC crossover, as it was anticipated in Section 1.2.

### 2.3 Extension to finite temperature

At finite temperature, the presence of the Fermi function in the gap equation (8) and density equation (9) makes it impossible to find an analytic solution and one has then to resort to a numerical solution of these equations. An exception occurs, however, upon approaching the critical temperature from below, where analytic results can be obtained both in BCS and BEC limits (see, e.g., Ref. Strinati-2000 ()).

At the mean-field level, the critical temperature is defined by the vanishing of the BCS gap parameter . Equations (8) and (9) are thus replaced by:

 −m4πaF=∫dk(2π)3(tanh(ξk/2Tc)2ξk−mk2) (21) n=2∫dk(2π)31exp(ξk/Tc)+1. (22)

In the weak-coupling limit where , Eq. (22) yields the value of the chemical potential of a Fermi gas, which is only slightly smaller than the zero-temperature result provided . Correspondingly, Eq. (21) gives the following expression for :

 Tc=8eγEFπe2exp(π2kFaF) (23)

where is Euler’s constant (such that with reference to the result (17) that holds in the same limit). Since in the weak-coupling regime is negative and , Eq. (23) yields consistently with our assumptions.

In the strong-coupling limit where , the role of the two equations (21) and (22) is reversed. If we assume that , we may set in Eq. (21), making it to reduce to the bound-state equation and yielding at the leading order. The same kind of approximation, however, cannot be used in Eq. (22) since has to remain finite. In this equation, we set instead and obtain

 (24)

a result which coincides with the expression of the classical chemical potential at temperature . Using at this point the value as determined from Eq. (21), we can solve the expression (24) iteratively for yielding Sa-de-Melo-1993 ():

 Tc≃ε02ln(ε0EF)3/2 (25)

at the leading order in . Although consistently with our assumptions, the expression of given by Eq. (25) diverges in the BEC limit at fixed density, instead of approaching as expected the finite value at which the Bose-Einstein condensation of an ideal gas of composite bosons with mass and density occurs ( being the Riemann function of argument ). Physically, the mean-field temperature (25) corresponds to the dissociation of composite bosons, with the factor in the denominator originating from entropy effects. It is thus appropriate to associate the mean-field result (25) with the “pair dissociation temperature” of the composite bosons, which is completely unrelated to the BEC temperature at which quantum coherence is established among the composite bosons. The reason for this failure is that only the internal degrees of freedom of the composite bosons are taken into account by the mean-field expressions (21) and (22), thereby leaving aside the translational degrees of freedom of the composite bosons. To include these, pairing fluctuations beyond mean field need to be considered. The corresponding analysis will be discussed in Section 3.1.

### 2.4 Bogoliubov-de Gennes and Hartree-Fock-Bogoliubov equations for inhomogeneous Fermi systems

The BCS wave function (1), for the ground state of a homogeneous Fermi gas with an attractive inter-particle interaction, is a variational wave function with an associated gap parameter at zero temperature, which extends to the superfluid phase the Hartree-Fock approximation for the normal phase by retaining also operator averages of the type that do not conserve the particle number.

At finite temperature, the “mean-field” character of the BCS gap equation (8) for and of the associated coherence factors and given by the expressions (2), is better captured by an alternative procedure, whereby a mean-field (or Hartree-Fock-type) decoupling is directly performed at the level of the Hamiltonian by including particle non-conserving averages. Accordingly, one replaces the grand-canonical Hamiltonian by its mean-field approximation as follows

 K=∑kσξkc†kσckσ+\varv0∑kk′qc†k+q/2↑c†−k+q/2↓c−k′+q/2↓ck′+q/2↑→∑kσξkc†kσckσ−∑k(Δ∗c−k↓ck↑+Δc†k↑c†−k↓)−Δ2\varv0 (26)

where, for simplicity, we have considered the case of a zero-range interaction with strength (the generalization to finite-range interactions will be considered below). The parameter in Eq.(26) thus satisfies the condition

 Δ=−\varv0∑k⟨c−k↓ck↑⟩ (27)

where the thermal average is self-consistently determined through the mean-field Hamiltonian itself. The quadratic form (26) can be readily diagonalized by the canonical transformation

 ck↑=u∗kαk+\varvkβ†−k,ck↓=−\varvkα†−k+u∗kβk, (28)

where and are constrained by at any given . In terms of the new operators and , the mean-field Hamiltonian (26) becomes (see, e.g., Chapt. 6 of Ref. Ring-1980 ())

 K→∑k(α†kαk+β†kβk)[ξk(u2k−\varv2k)+2Δuk\varvk]+∑k(α†kβ†−k+β−kαk)[2ξkuk\varvk−Δ(u2k−\varv2k)]+2∑k(ξk\varv2k−Δuk\varvk)−Δ2\varv0 (29)

(here, for simplicity, all quantities are assumed to be real). Provided that and are solutions to the equation

 (30)

whereby and , the coefficient of the second term on the right-hand side of Eq. (29) vanishes identically, while the coefficient of the first term equals . Correspondingly, the self-consistent condition (27) for reduces to:

 Δ=−\varv0∑kuk\varvk[1−2f(Ek)] (31)

with . Note further that the last term of Eq. (29), which can be written as with the help of Eq. (30), corresponds to the (grand-canonical) ground-state energy associated with the BCS state.

What is relevant here is that the above procedure can be readily generalized to situations in which depends on the spatial position , owing, for instance, to the presence of an external potential (like a trapping potential or a barrier), or of a magnetic field, or of the self-consistent mean field in finite nuclei, which introduce spatial inhomogeneities in the system. In this case, the algebraic matrix equation (30) is replaced by a pair of coupled Schrödinger-like equations for a two-component single-particle fermionic eigenfunction, of the form de-Gennes-1966 ():

 (H(r)Δ(r)Δ∗(r)−H(r))(uν(r)\varvν(r))=Eν(uν(r)\varvν(r)) (32)

where contains the external potential (in addition, in the presence of an external vector potential the gradient operator in the kinetic energy is replaced as usual by ). Note that the Hartree term is absent in the diagonal elements of the matrix (32), owing to the use of a contact inter-particle interaction which makes this term to vanish in the limit . [In weak coupling, however, one often introduces in a Hartree-like-term of the form Bruun-1999 (); Grasso-2003 ().] These equations, known as the Bogoliubov-de Gennes (BdG) equations, have to be solved up to self-consistency for the local gap parameter , which acts like an off-diagonal pair potential and satisfies the local condition:

 Δ(r)=−\varv0∑νuν(r)\varv∗ν(r)[1−2f(Eν)]. (33)

In addition, the eigenfunctions obey the orthonormality condition:

 ∫dr[u∗ν(r)uν′(r)+\varv∗ν(r)\varvν′(r)]=δνν′ (34)

where the Kronecker delta on the right-hand side is readily generalized to account for continuous eigenvalues. Physical quantities like the local number density, the current density, and the energy density can all be expressed in terms of the eigenfunctions and eigenvalues of Eq. (32) de-Gennes-1966 ().

Note that the complex character of (as well as of the eigenfunctions ) has been restored in Eqs. (32)–(34), to allow for the presence of a particle current. Note also that it is the particle-hole mixing characteristic of the BCS pairing to be responsible for the coupling of the two components of the eigenfunctions of Eq. (32), in contrast to the normal Fermi gas where particle- and hole-excitations are separately good quasi-particles.

When solving the BdG equations for a contact inter-particle interaction with coupling constant (which is the case of ultra-cold Fermi gases), technical problems arise from the need to regularize the self-consistency condition (33) for the gap parameter since the sum over therein diverges in the ultraviolet. In the homogeneous case with a uniform gap parameter this regularization can be readily implemented with the help of Eq.(7), by expressing the bare coupling constant that enters the gap equation in terms of the scattering length of the two-body problem [cf. Eq.(8)]. This procedure, however, cannot be utilized when the gap parameter has a spatial dependence, and a new strategy is required. A first generalization of the regularization procedure for fermions on the BCS side of the crossover trapped in a harmonic potential was given in Ref. Bruun-1999 () in terms of a pseudo-potential method in real space, but it turns out that this method is not easily implemented in numerical calculations. A numerically more efficient regularization procedure was described in Ref. Grasso-2003 () for the solution of the BdG equations (32) in a harmonic oscillator basis, which for practical purposes has to be truncated at some energy cutoff . Following Refs. Bulgac-2002 (); Yu-2003 (), the contribution of states above the cutoff is included within the Thomas-Fermi approximation. In this way, one obtains the following regularized version of the gap equation:

 Δ(r)=−geff(r)Eν

where the numerical solution of the BdG equations (32) is explicitly performed only for eigenvalues up to that appear on the right-hand side of Eq.(35). The pre-factor plays the role of a (cutoff dependent) effective coupling constant and is given by

 1geff(r)=m4πaF−R(r) (36)

where

 R(r)=m2π2[kc(r)−kμ(r)2lnkc(r)+kμ(r)kc(r)−kμ(r)]. (37)

The position-dependent wave vectors and are related to the chemical potential and energy cutoff by and . As noticed in Refs. Bulgac-2002 (); Grasso-2003 (), Eq. (37) can also be used for negative values of by allowing for imaginary values of .

This approach was extended in Ref. Simonucci-2013 () to the whole BCS-BEC crossover and generic inhomogeneous situations, by combining the introduction of the cutoff in the quasi-particle energies such that with the derivation of the Gross-Pitaevskii equation for composite bosons in the BEC limit that was done in Ref. Pieri-2003 () (cf. Section 2.5). Under the assumption that the energy cutoff is the largest energy scale in the problem (such that , , and , where is the Fermi energy associated with the mean density - conditions that can be somewhat relaxed in the weak-coupling limit), the regularized version of the gap equation reads eventually Simonucci-2013 ():

 {−m4πaF+R(r)−[12I02(r)−13I13(r)]∇24m+2Vext(r)14I02(r)+14I03(r)|Δ(r)|2}Δ(r)=Eν

The quantities that enter the left-hand side of this equation are defined as follows:

 Iij(r)=∫|k|>kc(r)dk(2π)3(k22m)i(k22m+Vext(r)−μ)j. (39)

The right-hand side of Eq.(38), which results from an explicit numerical integration of the BdG equations over the reduced energy range , acts as a source term on the non-linear differential equation for associated with the left-hand side. By implementing the numerical calculation of the BdG equations for an isolated vortex embedded in a uniform superfluid, it was shown in Ref. Simonucci-2013 () that the inclusion of the various terms on the left-hand side of Eq.(38) [from the linear () term, to the linear plus cubic () terms, and finally to the linear plus cubic plus Laplacian () terms] plays an increasingly important role in reducing the total computational time and memory space at any coupling throughout the BCS-BEC crossover, by decreasing the value of the cutoff up to which the eigenfunctions of the BdG equations have to be explicitly calculated Simonucci-2013 ().

In practice, the BdG equations have been solved numerically up to self-consistency only for relatively simple geometries since their numerical solution becomes rapidly too demanding due to computational time and memory space. For ultra-cold Fermi gases, the BdG equations have been solved to account for the spatial dependence of the order parameter in balanced Bruun-1999 (); Grasso-2003 (); Ohashi-2005 () and imbalanced Castorina-2005 (); Jensen-2007 () trapped gases, as well as for the microscopic structure of a single vortex at zero temperature on the BCS side Nygaard-2003 () and across the BCS-BEC crossover Sensarma-2006 (), and also at finite temperature across the BCS-BEC crossover Simonucci-2013 (). The microscopic structure of a single vortex for the unitary Fermi gas has been analyzed in Ref. Bulgac-2003 () with BdG-like equations resulting from the superfluid local-density approximation of Ref. Bulgac-2002b (), which has later been used also to calculate the profile of the order parameter in balanced trapped gases Bulgac-2007 () or in the FFLO phase of imbalanced Fermi gases throughout the BCS-BEC crossover Bulgac-2008b () (see Section 2.6). In addition, the BdG equations have been used to determine the occurrence of vortex lattices in a rotating trap on the BCS side of the crossover Feder-2004 (); Tonini-2006 (). A study of the Josephson effect with a one-dimensional barrier throughout the BCS-BEC crossover at zero temperature Spuntarelli-2010 () will be reported in Section 4.6.

In nuclear physics, where in the most sophisticated cases finite-range (instead of zero-range) forces are used in the mean-field and gap equations (a famous example in this context being the Gogny-force Decharge-1980 ()), both matrix elements in (32) become non-local. The corresponding equations, usually referred to as the Hartree-Fock-Bogoliubov (HFB) equations, are thus more general than the BdG equations and read (for simplicity, we do not write down the spin and isospin structure):

 ∫dr′(H(r,r′)Δ(r,r′)Δ∗(r,r′)−H∗(r,r′))(uν(r′)\varvν(r′))=Eν(uν(r)\varvν(r)) (40)

where

 H(r,r′)=[−∇2r2m+ΓH(r)]δ(r−r′)+ΓF(r,r′) (41)

and

 (42)

In the above expressions, and are the Hartree and Fock fields, respectively, with the density matrix given by (such that is the number density). In nuclear physics, the phenomenological effective forces are generally density dependent, also in the pairing channel, as indicated in the argument of in Eq.(42).

The non-local HFB equations (40), that hold for finite nuclei, are numerically much harder to solve than the local BdG equations (32). Yet, they are almost routinely solved even for deformed and rotating nuclei, in general the results of the calculations being in excellent agreement with experiments (see Refs. Bender-2003 (); Vretenar-2005 () for reviews). Owing to the high numerical cost to solve these equations, however, approximate local approximations have been developed for the mean-field terms of Eq.(41) (such as the Skyrme Energy Density Functionals (EDF) Erler-2011 ()) and for the pairing interaction of Eq.(42).

The need to regularize the gap equation arises also for the nuclear problem. In this context, it is worth mentioning that the renormalization scheme of Eq. (35) discussed above was originally developed in Refs. Bulgac-2002 (); Yu-2003 () for the nuclear problem, in the case when the gap equation reduces to a local form. In this case, Eqs. (36) and (37) are replaced by an effective coupling constant which takes the form:

 1geff(r)=1g(r)−m∗(r)2π2[kc(r)−kμ(r)2lnkc(r)+kμ(r)kc(r)−kμ(r)]. (43)

Here, in contrast to Eqs. (36) and (37), the effective mass and the bare coupling constant are now position dependent through the value of the local density, which can be fitted such that as a function of reproduces the gap obtained with the Gogny force Garrido-2001 () in infinite matter. In addition, the wave vectors and