Shear viscosity and Strong-Coupling Corrections in the BCS-BEC crossover Regime of an Ultracold Fermi Gas

TL;DR

This paper theoretically studies shear viscosity in an ultracold Fermi gas across the BCS-BEC crossover, revealing how strong interactions influence transport properties and aligning with recent experimental findings.

Contribution

It provides a detailed analysis of shear viscosity behavior in the crossover regime using a strong-coupling self-consistent T-matrix approach, highlighting the minimum near the BEC side.

Findings

Shear viscosity diverges in weak and strong coupling limits.

Minimum shear viscosity occurs slightly on the BEC side, not at unitarity.

Anomalous temperature dependence of shear viscosity near T_c in the BCS regime.

Abstract

We theoretically investigate the shear viscosity $\eta$ in the BCS-BEC crossover regime of an ultracold Fermi gas with a Feshbach resonance. Within the framework of the strong-coupling self-consistent $T$-matrix approximation, we examine how a strong pairing interaction associated with a Feshbach resonance affects this transport coefficient, in the normal state above the superfluid phase transition temperature $T_{\rm c}$. We show that, while $\eta$ diverges in both the weak-coupling BCS and strong-coupling BEC limits, it becomes small in the unitary regime. The minimum of $\eta$ is obtained, not at the unitarity, but slightly in the strong-coupling BEC side. This deviation is consistent with the recent experiment on a $^6$Li Fermi gas. In the weak-coupling BCS regime, we also find that $\eta$ exhibits anomalous temperature dependence near $T_{\rm c}$, which is deeply related to the pseudogap phenomenon originating form strong pairing fluctuations.

Equations

$$\frac{\eta }{s}\ge \frac{\hslash }{4\pi k_{\mathrm{B}}}.$$

$$H=\sum _{\mathbf{p},\sigma }{\xi }_{\mathbf{p}}{c}_{\mathbf{p},\sigma }^{†}{c}_{\mathbf{p},\sigma }-U\sum _{\mathbf{p},{\mathbf{p}}^{\prime },\mathbf{q}}{c}_{\mathbf{p}+\mathbf{q},\uparrow }^{†}{c}_{{\mathbf{p}}^{\prime }-\mathbf{q},\downarrow }^{†}{c}_{{\mathbf{p}}^{\prime },\downarrow }{c}_{\mathbf{p},\uparrow },$$

$$\frac{4\pi a_s}{m}=\frac{-U}{1-U{\sum }_{\mathbf{p}}^{p_{\mathrm{c}}}\frac{1}{2{\epsilon }_{\mathbf{p}}}},$$

$$G(\mathbf{p},i{\omega }_n)=\frac{1}{i{\omega }_n-{\xi }_{\mathbf{p}}-\Sigma (\mathbf{p},i{\omega }_n)},$$

$$\Sigma (\mathbf{p},i{\omega }_n)=T\sum _{\mathbf{q},{\nu }_m}\Gamma (\mathbf{q},i{\nu }_m)G(\mathbf{q}-\mathbf{p},i{\nu }_m-i{\omega }_n).$$

$$\Gamma (\mathbf{q},i{\nu }_m)$$

$$=\frac{4\pi a_s}{m}\frac{1}{1+\frac{4\pi a_s}{m}\left[\Lambda (\mathbf{q},i{\nu }_m)-{\sum }_{\mathbf{p}}\frac{1}{2{\epsilon }_{\mathbf{p}}}\right]}$$

$$\Lambda (\mathbf{q},i{\nu }_m)=T\sum _{\mathbf{p},{\omega }_n}G(\mathbf{p},i{\omega }_n)G(\mathbf{q}-\mathbf{p},i{\nu }_m-i{\omega }_n),$$

$$1=-\frac{4\pi a_s}{m}\left[\Lambda (\mathbf{q}=0,i{\nu }_m=0)-\sum _{\mathbf{p}}\frac{1}{2{\epsilon }_{\mathbf{p}}}\right].$$

$$N=2T\sum _{\mathbf{p},{\omega }_n}G(\mathbf{p},i{\omega }_n),$$

$$\eta =-\underset{\omega \to 0}{\lim }\frac{1}{\omega }\mathrm{Im}[\Xi (\omega )].$$

$$\Xi (\omega )=-i\int d^3\mathbf{r}{\int }_0^{\infty }dt{e}^{i\omega t}⟨[{\hat{\Pi }}_{x,y}(\mathbf{r},t),{\hat{\Pi }}_{x,y}(0,0)]⟩$$

$${\hat{\Pi }}_{x,y}(\mathbf{r})$$

$$\Xi (i{\nu }_l)=-\int d^3\mathbf{r}{\int }_0^{\beta }{e}^{i{\nu }_l\tau }⟨T_{\tau }{\hat{\Pi }}_{x,y}(\mathbf{r},\tau ){\hat{\Pi }}_{x,y}(0,0)⟩,$$

$$\Xi (i{\nu }_l)$$

$$\times G(\mathbf{p},i{\omega }_n+i{\nu }_l).$$

$$T_{x,y}^{\mathrm{MT}}(\mathbf{p},i{\omega }_n,i{\omega }_n+i{\nu }_l)$$

$$=T\sum _{\mathbf{q},{\nu }_m}\Gamma (\mathbf{q},i{\nu }_m){\tilde{T}}_{x,y}(\mathbf{q}-\mathbf{p},i{\nu }_m-i{\omega }_n-i{\nu }_l,i{\nu }_m-i{\omega }_n),$$

$$T_{x,y}^{\mathrm{AL}}(\mathbf{p},i{\omega }_n,i{\omega }_n+i{\nu }_l)$$

$$=-2T\sum _{\mathbf{q},{\nu }_m}{\tilde{S}}_{x,y}(\mathbf{q},i{\nu }_m,i{\nu }_m+i{\nu }_l)G(\mathbf{q}-\mathbf{p},i{\nu }_m-i{\omega }_n).$$

$${\tilde{T}}_{x,y}(\mathbf{p},i{\omega }_n,i{\omega }_n+i{\nu }_l)$$

$$=G(\mathbf{p},i{\omega }_n)T_{x,y}(\mathbf{p},i{\omega }_n,i{\omega }_n+i{\nu }_l)G(\mathbf{p},i{\omega }_n+i{\nu }_l),$$

$${\tilde{S}}_{x,y}(\mathbf{q},i{\nu }_m,i{\nu }_m+i{\nu }_l)=\Gamma (\mathbf{q},i{\nu }_m)\Gamma (\mathbf{q},i{\nu }_m+i{\nu }_l)$$

$$\times T\sum _{\mathbf{p},{\omega }_n}G(\mathbf{q}-\mathbf{p},i{\nu }_m-i{\omega }_n){\tilde{T}}_{x,y}(\mathbf{p},i{\omega }_n,i{\omega }_n+i{\nu }_l).$$

$$\rho (\omega )=-\frac{1}{\pi }\sum _{\mathbf{p}}\mathrm{Im}[G(\mathbf{p},i{\omega }_n\to \omega +i\delta )],$$

$$\Gamma (\mathbf{q},i{\nu }_m)\simeq \frac{8\pi }{m^2{a}_s}\frac{1}{i{\nu }_m-{\xi }_{\mathbf{q}}^{\mathrm{B}}}.$$

$$U_{\mathrm{B}}=\frac{4\pi a_{\mathrm{B}}}{m_{\mathrm{B}}},$$

$$\Sigma (\mathbf{p},i{\omega }_n)$$

$$\simeq G_0(-\mathbf{p},-i{\omega }_n)\times T\sum _{\mathbf{q},{\nu }_m}\Gamma (\mathbf{q},i{\nu }_m)$$

$$\equiv -{\Delta }_{\mathrm{pg}}^2{G}_0(-\mathbf{p},-i{\omega }_n).$$

$$G(\mathbf{p},i{\omega }_n)=-\frac{i{\omega }_n+{\xi }_{\mathbf{p}}}{{\omega }_n^2+{\xi }_{\mathbf{p}}^2+{\Delta }_{\mathrm{pg}}^2}.$$

$$\gamma (\mathbf{q},\omega \ge 0)$$

$$\times [1-f_{\mathrm{F}}(E_{\mathbf{q}/2+\mathbf{p}})-f_{\mathrm{F}}(E_{\mathbf{q}/2-\mathbf{p}})]\delta (\omega -[E_{\mathbf{q}/2+\mathbf{p}}+E_{\mathbf{q}/2-\mathbf{p}}])$$

$$-\frac{\pi }{2}\sum _{\mathbf{p}}\left[1+\frac{{\xi }_{\mathbf{q}/2+\mathbf{p}}}{E_{\mathbf{q}/2+\mathbf{p}}}\right]\left[1-\frac{{\xi }_{\mathbf{q}/2-\mathbf{p}}}{E_{\mathbf{q}/2-\mathbf{p}}}\right]$$

$$\times [f_{\mathrm{F}}(E_{\mathbf{q}/2+\mathbf{p}})-f_{\mathrm{F}}(E_{\mathbf{q}/2-\mathbf{p}})]\delta (\omega -[E_{\mathbf{q}/2+\mathbf{p}}-E_{\mathbf{q}/2-\mathbf{p}}]),$$

$$\gamma (\mathbf{q},\omega )$$

Full text

1 Introduction

Since the realization of the superfluid phase transition in 40K[1] and 6Li Fermi gases[4, 3, 2], thermodynamic properties of this system in the BCS (Bardeen-Cooper-Schrieffer)-BEC (Bose-Einstein condensation) crossover region[13, 6, 12, 11, 14, 5, 10, 9, 8, 7] have extensively been studied[16, 17, 15, 26, 27, 29, 25, 23, 24, 22, 21, 20, 30, 19, 18, 28], by maximally using the advantage that a paring interaction in this system can be tuned by adjusting the threshold energy of a Feshbach resonance[15]. Furthermore, as a new research direction, non-equilibrium transport properties in the BCS-BEC crossover region have recently been discussed[41, 39, 21, 32, 31, 38, 40, 36, 37, 35, 34, 33].

Regarding this new direction, the shear viscosity $\eta$ has particularly attracted much attention[42, 62, 59, 58, 66, 57, 65, 64, 56, 55, 54, 53, 61, 60, 52, 51, 50, 49, 48, 63, 47, 46, 45, 44, 43]. One reason is that transport coefficients are related to the particle mean-free-path originating from a pairing interaction, so that they involve useful information about collisional properties of the system in the BCS-BEC crossover region[42]. Another reason is the stimulation by the conjecture by Kovtun, Son, and Starinets (KSS)[67], stating that the ratio of $\eta$ to the entropy density $s$ has the lower bound as,

[TABLE]

This was originally proposed in relativistic quantum field theories at finite temperature and zero chemical potential; however, because Eq. (1) does not involve the speed of light, KSS speculated that it may be valid for the non-relativistic case, at least for a single-component gas with spin 0 or 1/2. Although some formal counter-examples are known for this KSS conjecture[72, 71, 70, 69, 68, 73], the ratio $\eta/s$ is still considered as a useful quantity for the study of interaction effects on fluid properties. Indeed, Ref. [74] evaluated $\eta/s$ for the some quantum fluids from experimental data: (i) liquid 4He: \eta/s\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 8.8, (ii) unitary 6Li Fermi gas: \eta/s\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 6.3, and (iii) quark-gluon plasma: \eta/s\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 5.0, in unit of $\hbar/(4\pi k_{\rm B})$. For comparison, it is useful to note that $\eta/s\simeq 380\times\hbar/(4\pi k_{\rm B})$ in water under the normal condition[67]. Thus, the results in the above-mentioned quantum fluids are very close to the so-called KSS bound $\eta/s=\hbar/(4\pi k_{\rm B})$.

From the viewpoint of BCS-BEC crossover physics, because the shear viscosity diverges in an ideal gas, the above result, \eta/s\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 6.3\times\hbar/(4\pi k_{\rm B}), obtained in a unitary 6Li Fermi gas[74] implies that the shear viscosity becomes the smallest around the unitary limit. (Note that the system is reduced to an ideal gas in both the BCS and BEC limits.) Regarding this, it has recently been observed in a 6Li Fermi gas[40] that, although the shear viscosity $\eta$ becomes small in the unitary regime, the minimum value is obtained, not at the unitarity, but slightly in the strong-coupling BEC side. This makes us expect that the minimum value of $\eta/s$ is also obtained (slightly) away from the unitarity limit.

Motivated by these, in this paper, we theoretically investigate the shear viscosity $\eta$ in an ultracold Fermi gas in the BCS-BEC crossover region. Including strong pairing fluctuations associated with a Feshbach-induced tunable pairing interaction within the framework of the self-consistent $T$-matrix approximation (SCTMA)[14, 5, 75, 64], we calculate this transport coefficient, by using the linear response theory. In the first step toward the assessment of the KSS conjecture, this paper focuses on the shear viscosity $\eta$. The ratio $\eta/s$ will separately be discussed in our subsequent paper. We clarify how $\eta$ behaves in the normal state above $T_{\rm c}$ in the BCS-BEC crossover region. We also examine where $\eta$ becomes minimum in the unitary regime, comparing the recent experiment on a 6Li Fermi gas[40].

This paper is organized as follows. In Sec. 2, we explain our formulation. The calculated shear viscosity in the BCS-BEC region is shown in Sec. 3. We also compare our results with the recent experiments on 6Li Fermi gases in this section. In the followings, we set $\hbar=k_{\rm B}=1$, and the system volume is taken to be unity (except in appendix A), for simplicity.

2 Formulation

We consider a two-component unpolarized Fermi gas, described by the ordinary BCS Hamiltonian,

[TABLE]

where $c_{\bm{p},\sigma}^{\dagger}$ is the creation operator of a Fermi atom with pseudospin $\sigma=\uparrow,\downarrow$, describing two atomic hyperfine states. $\xi_{\bm{p}}=\varepsilon_{\bm{p}}-\mu=p^{2}/(2m)-\mu$ is the kinetic energy of a Fermi atom, measured from the Fermi chemical potential $\mu$, where $m$ is an atomic mass. $-U~{}(<0)$ is a pairing interaction, which is assumed to be tunable by the Feshbach-resonance technique[15]. The interaction strength is conveniently measured in terms of the $s$-wave scattering length $a_{s}$, which is related to the bare interaction $-U$ as

[TABLE]

where $p_{\rm c}$ is a momentum cutoff. In this scale, the weak-coupling BCS regime and the strong-coupling BEC regime are, respectively, characterized as (k_{\rm F}a_{s})\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ -1, and $(k_{\rm F}a_{s})^{-1}\gtrsim 1$ (where $k_{\rm F}$ is the Fermi wave-length). The region, -1\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 1, is referred to as the BCS-BEC crossover regime.

Strong-coupling corrections to single-particle excitations can be described by the self-energy $\Sigma({\bm{p}},i\omega_{n})$ in the single-particle thermal Green’s function,

[TABLE]

where $\omega_{n}$ is the fermion Matsubara frequency. In SCTMA[5, 14, 75, 64], $\Sigma({\bm{p}},i\omega_{n})$ is obtained from the diagrams in Fig. 1, which gives

[TABLE]

Here, $\nu_{m}$ is the boson Matsubara frequency, and

[TABLE]

is the SCTMA particle-particle scattering matrix, describing fluctuations in the Cooper channel (see also Fig. 1(b)). Here,

[TABLE]

is the pair correlation function.

In this scheme, we determine $T_{\rm c}$ from the Thouless criterion[76], stating that the superfluid instability occurs, when the particle-particle scattering matrix in Eq. (6) has a pole at ${\bm{q}}=\nu_{m}=0$, which gives

[TABLE]

As in the ordinary BCS-BEC crossover theories[13, 6, 12, 11, 14, 5, 10, 9, 8, 7], we actually solve the $T_{\rm c}$-equation (8), together with the equation for the number $N$ of Fermi atoms,

[TABLE]

to self-consistently determine $T_{\rm c}$ and $\mu(T_{\rm c})$. Above $T_{\rm c}$, we only treat the number equation (9), to evaluate $\mu(T)$. We briefly show the self-consistent solution for $\mu(T\geq T_{\rm c})$ in Fig, 2, which is used in calculating $\eta(T)$.

In the linear response theory[77, 78], the shear viscosity $\eta$ is given by[80, 79, 66, 65],

[TABLE]

Here,

[TABLE]

is the shear-stress response function, where $\hat{\Pi}_{x,y}$ is the $xy$ component of the stress tensor operator. In the present BCS model in Eq. (2), it has the form,

[TABLE]

where the field operator $\psi_{\sigma}({\bm{r}})$ describes Fermi atoms with pseudospin $\sigma=\uparrow,\downarrow$. For the derivation of Eq (12), see Appendix A.

In this paper, we first evaluate the corresponding thermal response function,

[TABLE]

where $T_{\tau}$ is the imaginary-time-ordered-product and $\nu_{l}$ is the boson Matsubara frequency. $\Xi(\omega)$ is then obtained as $\Xi(\omega)=\Xi(i\nu_{l}\to\omega_{+}=\omega+i\delta)$, where $\delta$ is an infinitesimally small positive number. In this paper, we numerically carry out this analytic continuation by the Padé approximation[81].

In evaluating the thermal response function in Eq. (13), we need to choose diagrams so as to satisfy the Ward-Takahashi identity[82, 84, 85, 83], which is a required condition in any consistent theory. For the shear-stress response function $\Xi(i\nu_{l})$, this identity is derived from the momentum conservation law[64, 83]. In SCTMA, it is diagrammatically described in Fig. 3[64], giving

[TABLE]

Here, the three point vertex function $T_{x,y}=p_{x}p_{y}/m+T_{x,y}^{\rm MT}+T_{x,y}^{\rm AL}$ consists of the bare term $p_{x}p_{y}/m$, Maki-Thompson (MT) term,

[TABLE]

as well as the Aslamazov-Larkin (AL) term,

[TABLE]

Here, ${\tilde{T}}_{x,y}$ and ${\tilde{S}}_{x,y}$ are, respectively, given as follows:

[TABLE]

[TABLE]

Before ending this section, we comment on our numerical calculations to obtain $\eta(T)$. In computing the shear viscosity, we have sometimes met the difficulty that the Padé approximation unphysically gives negative $\eta$ or positive but abnormally large/small $\eta$. We have also found that this problem depends on the detailed choice of momentum cutoff in numerically evaluating $\Xi(i\nu_{l})$ in Eq. (14), as well as on the number of Matsubara frequencies in executing the numerical analytic continuation by the Padé approximation. At this stage, we have no idea to completely overcome this problem. Thus, we have employed the following prescription in this paper: (1) We first calculate $\Xi(i\nu_{l})$ by introducing various values of the momentum cutoff $k_{\rm c}$ (10k_{\rm F}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ k_{\rm c}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 60k_{\rm F}) to the momentum summation in Eq. (14) (which is nothing to do with $p_{\rm c}$ in Eq. (3)). (2) For each result, we next numerically execute the analytic continuation by the Padé approximation[81], retaining $50\sim 100$ Matsubara frequencies $\nu_{l}$. (3) For these data set, we remove clearly unphysical negative data, and then remove the highest and lowest 10% of data to avoid influence of abnormal results. (4) For the remaining data, we evaluate the averaged value $\bar{\eta}$, as well as the standard derivation ${\bar{\sigma}}$. When $|{\bar{\sigma}}/\bar{\eta}|$ is less than 0.1, we plot it in Fig. 4. Otherwise, we judge that the result is not reliable, not to plot it in Fig. 4. The latter situation occurs near $T_{\rm c}$ in the BEC regime, so that $\eta$ is not shown there in this figure.

3 Shear viscosity $\eta$ in the BCS-BEC crossover regime of an ultracold Fermi gas

Figure 4 shows the shear viscosity $\eta(T)$ in the BCS-BEC crossover regime of an ultracold Fermi gas above $T_{\rm c}$. When (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ -0.5 in the weak-coupling BCS side, we find that, with decreasing the temperature from the Fermi temperature $T_{\rm F}$, $\eta$ exhibits a dip structure. Deep inside the BCS regime ((k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ -2), Fig. 5(a), as well as the inset in this panel, show that this non-monotonic temperature dependence is well reproduced by the Boltzmann equation. As summarized in Appendix B, in the classical regime ($T\gg T_{\rm F}$), the Boltzmann-equation approach gives $\eta(T)\propto T^{3/2}$ ($T\gg 1/(ma_{s}^{2})$) and $\eta(T)\propto T^{1/2}$ ($T\ll 1/(ma_{s}^{2})$)[59, 62], both of which decreases with decreasing the temperature. In the Fermi degenerate regime ($T\ll T_{\rm F}$), this approach gives $\eta(T)\propto T^{-2}$, which increases as the temperature decreases. These explain the behavior of $\eta(T)$ above and below the dip temperature ($\equiv T^{\rm BCS}_{\rm dip}$), respectively. Thus, this dip structure is nothing to do with the BCS-BEC crossover phenomenon, but is a phenomenon which can be explained within the standard Boltzmann equation ignoring pairing fluctuations[86].

In addition to the dip, we also see in Fig. 4 a peak structure near $T_{\rm c}$ when (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ -0.5. As shown in Fig. 5(a), the Boltzmann equation cannot explain this anomaly, indicating that it is a many-body phenomenon originating from pairing fluctuations enhanced near $T_{\rm c}$ (that are ignored in the Boltzmann equation). Indeed, Fig. 6 shows that the peak temperature $T_{\rm peak}$ is close to the so-called pseudogap temperature $T^{*}$[87], below which the single-particle density of states,

[TABLE]

has a BCS-state like dip structure around the Fermi level $\omega=0$, due to strong pairing fluctuations (or the formation of preformed Cooper pairs).

To simply understand the role of pairing fluctuations around $T_{\rm peak}$, we conveniently employ the expression for the shear viscosity $\eta\sim nl_{\rm mfp}{\bar{p}}$[88], where $n$ is the number density of Fermi atoms, $l_{\rm mfp}$ the mean free path, and ${\bar{p}}$ is the averaged particle momentum. When the effective interaction between Fermi atoms described by the particle-particle scattering matrix $\Gamma$ in Eq. (6) is temperature-independent far above $T_{\rm c}$, the quasi-particle lifetime $\tau$ in the mean free path $l_{\rm mfp}\sim(k_{\rm F}/m)\tau$ behaves as, symbolically, $\tau^{-1}\sim|\Gamma|^{2}T^{2}\propto T^{2}$ when $T\ll T_{\rm F}$, as in the ordinary Fermi liquid theory. The resulting $\eta(T)\propto T^{-2}$ explains the temperature dependence of $\eta(T)$ in the region $T_{\rm peak}\leq T\leq T^{\rm BCS}_{\rm dip}$. Near $T_{\rm c}$, on the other hand, the particle-particle scattering matrix $\Gamma(0,0)$ is enhanced, to eventually diverge at $T_{\rm c}$ (Thouless criterion). This enhancement of $\Gamma(0,0)$ with decreasing the temperature shortens the quasi-particle lifetime $\tau$, which leads to the decrease of $\eta(T)\propto\tau\sim[|\Gamma|T]^{-2}$, giving the peak structure in Fig. 4.

We briefly note that low-energy properties of a repulsively interacting Fermi gas at low temperatures are known to be well described by the Fermi liquid theory with the quasi-particle lifetime $\tau\propto T^{-2}$. Regarding this, although the present interaction is attractive, one still sees the “Fermi liquid like behavior”, $\eta(T)\propto\tau\propto T^{-2}$, in the intermediate temperature region, T_{\rm peak}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ T\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ T_{\rm dip}^{\rm BCS}.

Figures 4 and 5 show that the above-mentioned dip-peak structure disappears in the unitary regime (-0.5\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 0.5), because of the smearing of the Fermi surface, as well as the enhancement of fluctuations in the Cooper channel, by the strong pairing interaction between Fermi atoms. The resulting shear viscosity $\eta(T)$ monotonically decreases with decreasing the temperature down to $T_{\rm c}$. We point out that such behavior has recently been observed in a 6Li unitary Fermi gas[37]; our result semi-quantitatively agrees with this experiment, as shown in Fig. 7. We briefly note that this behavior of $\eta(T)$ has also theoretically been obtained in Ref. [64].

As one passes through the unitarity regime to enter the strong-coupling BEC regime ((k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 0.7), the shear viscosity $\eta(T)$ again becomes large, as shown in Figs. 4 and 5(b). In addition, although we cannot completely obtain $\eta$ near $T_{\rm c}$ in the BEC regime because of the computational problem mentioned in Sec. 2, we still see in Fig. 5(b) that $\eta(T)$ again exhibits a non-monotonic behavior, when (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 0.7. Plotting the dip temperature ($\equiv T_{\rm dip}^{\rm BEC}$) in this regime, we find that it is close to the binding energy $E_{\rm bind}=1/(ma_{s}^{2})$ of a two-body bound state (see Fig. 6). This implies that the dip structure in the BEC regime is associated with the start of the formation of two-body bound molecules below $T_{\rm dip}^{\rm BEC}$, overwhelming thermal dissociation.

To explain how Bose-molecules affect the shear viscosity $\eta$ in the BEC regime, it is convenient to represent the shear-stress response function $\Xi(i\nu_{l})$ with respect to the Bose degrees of freedom: In the extreme BEC limit ($(k_{\rm F}a_{s})^{-1}\gg 1$), the particle-particle scattering matrix $\Gamma({\bm{q}},i\nu_{l})$ is known to describe tightly bound molecular bosons. To simply see this, approximately evaluating the pair correlation function $\Lambda({\bm{q}},i\nu_{m})$ in Eq. (7) by replacing the SCTMA dressed Green’s function with the bare one, $G_{0}({\bm{p}},i\omega_{n})=[i\omega_{n}-\xi_{\bm{p}}]^{-1}$, one finds[14], in the BEC limit,

[TABLE]

Here, $\xi_{\bm{q}}^{\rm B}={\bm{q}}^{2}/(4m)-\mu_{\rm B}$ is the molecular kinetic energy, measured from the Bose chemical potential $\mu_{\rm B}=2\mu+E_{\rm bind}$. Then, (1) regarding the particle-particle scattering matrix $\Gamma$ appearing in the diagrams of $\Xi(i\nu_{l})$ in Fig. 3 as the Bose Green’s function (wavy line in Fig. 8(a)), and (2) further introducing an effective Bose-Bose interaction $U_{\rm B}$ as Fig. 8(b), as well as an bare and dressed molecular three-point vertex functions $S_{0;x,y}$ and $S_{x,y}$ as Figs. 8(c) and (d), respectively, we find that the AL-type diagrams in Fig. 3 involves the contribution that can be written in the form of the Bose shear-stress response function, as shown in Fig. 8(e).

In SCTMA, the effective molecular interaction $U_{\rm B}$, that are mediated by four unpaired Fermi atoms (see Fig. 8(b)), is repulsive and is given by[14]

[TABLE]

where $m_{\rm B}=2m$ is a molecular mass and $a_{\rm B}=2a_{s}>0$ is the $s$-wave molecular scattering. In the BEC regime, because most Fermi atoms form tightly bound molecules described by the wavy line in Fig. 8, the shear viscosity is considered to be dominated by the molecular shear-stress response function in Fig. 8(e). This contribution becomes large with decreasing the molecular interaction strength $U_{\rm B}\propto 2a_{s}\to 0$ with approaching the BEC limit $(k_{\rm F}a_{s})^{-1}\to+\infty$. (Note that the shear viscosity diverges in an ideal Bose gas.)

To support the above diagrammatic discussion, we show in Fig. 9 the result ($\equiv\eta^{\prime}$) in the case when we only retain the first term in Fig. 3 and ignore all the MT and AL vertex corrections. In this figure, we see that $\eta^{\prime}$ simply decreases as one approaches the BEC regime, because it does not involve molecular contribution to the shear viscosity. In addition, as shown in Fig. 5(b), $\eta^{\prime}$ monotonically decreases with decreasing the temperature when ($(k_{\rm F}a_{s})^{-1}\geq 0.7$), being in contrast to the up-turn behavior of $\eta(T)$ including vertex corrections.

When we apply the simple expression for the shear viscosity $\eta\sim nl_{\rm mfp}{\bar{p}}$ to the molecular Bose gas in the BEC regime, the mean-free path $l_{\rm mfp}$ is expected to be closely related to the molecular lifetime $\tau_{\rm B}$. Regarding this, we note that, although SCTMA includes the molecular interaction $U_{\rm B}$ in Eq. (21), this strong-coupling theory involves effects of this repulsive interaction only in the mean-field level[14]. That is, the molecular lifetime $\tau_{\rm B}$ coming from $U_{\rm B}$ is completely ignored in SCTMA. However, because $\eta\sim nl_{\rm mfp}{\bar{p}}$ diverges when $l_{\rm mfp}\propto\tau_{\rm B}=\infty$, the converging $\eta(T)$ in Fig. 5(b) implies the presence of another scattering process contributing to the molecular lifetime $\tau_{\rm B}$.

As the origin of $\tau_{\rm B}$ in SCTMA, we note that, because the particle-particle scattering matrix $\Gamma({\bm{q}},i\nu_{m})$ in Eq. (6) is reduced to the single-particle Bose Green’s function in the BEC regime, the lifetime of the boson should come from the pair correlation function, ${\rm Im}[\Lambda({\bm{q}},i\nu_{m}\to\omega_{+})]~{}(\equiv\gamma({\bm{q}},\omega))$. To confirm this, we evaluate this molecular damping $\gamma({\bm{q}},\omega)$ by simply approximating the self-energy $\Sigma({\bm{p}},i\omega_{n})$ in Eq. (5) involved in the dressed Green’s function $G$ to[90]

[TABLE]

Here, we have approximately set ${\bm{q}}=\nu_{m}=0$ in the bare Green’s function $G_{0}$, using the fact that $\Gamma(0,0)$ is enhanced near $T_{\rm c}$. $\Delta_{\rm pg}$ is the so-called pseudogap parameter[91, 90], physically describing effects of pairing fluctuations on single-particle excitations. The resulting dressed Green’s function $G$ formally has the same form as the diagonal component of the BCS Green’s function[77, 78] as,

[TABLE]

Substituting Eq. (23) into Eq. (7), as well as executing the $\omega_{n}$-summation, we obtain, after the analytic continuation $i\nu_{m}\to\omega_{+}$,

[TABLE]

where $E_{\bm{p}}=\sqrt{\xi_{\bm{p}}^{2}+\Delta_{\rm pg}^{2}}$ and $f_{\rm F}(x)=1/[\exp(x/T)+1]$ is the equilibrium Fermi distribution function. In Eq. (24), the first term comes from the so-called inter-band excitation[92]. Since this excitation is accompanied by pair-breaking, it has the excitation threshold, $\omega_{\rm th}={\rm Min}[E_{{\bm{q}}/2+{\bm{p}}}+E_{{\bm{q}}/2-{\bm{p}}}]=2|\Delta_{\rm pg}|$. As a result, this term does not contribute to the damping $\gamma({\bm{q}},\omega)$, as far as we consider low-energy excitations $\omega\leq\omega_{\rm th}$.

The second term in Eq. (24) is associated with the intra-band excitation[92]. Since this is an excitation of an unpaired Fermi atom which has already been excited thermally, it has no excitation threshold. Thus, the molecular damping $\gamma({\bm{q}},\omega)$ is dominated by this process. Noting that $\mu\to-E_{\rm bind}/2=-1/(2ma_{s}^{2})\ll-\varepsilon_{\rm F}$[6, 12, 11, 14, 5] and $E_{\bm{p}}\simeq\varepsilon_{\bm{p}}+E_{\rm bind}/2$ in the strong-coupling BEC regime, one obtains, in the low-energy and low-momentum regime,

[TABLE]

Even in the BEC regime, the thermal dissociation of bound molecules gives unpaired Fermi atoms to some extent, which is reflected by the factor $e^{-E_{\rm bind}/(2T)}~{}(\ll 1)$ in Eq. (25). This small factor describes the longevity of molecular bosons in this regime:

[TABLE]

In the strong-coupling BEC regime where most Fermi atoms form tightly bound molecules, $\eta$ is dominated by the molecular contribution in Fig. 8(e). To see how the long-lived molecules affect $\eta(T)$ in this regime, we consider a model Bose gas with $\tau_{\rm B}=\lambda e^{E_{\rm bind}/(2T)}$, ignoring the energy- and momentum-dependence that Eq. (26) possesses. In this simple model, we obtain (For the derivation, see Appendix C.),

[TABLE]

where $n_{\rm B}(x)=1/[\exp(x/T)-1]$ is the Bose distribution function. We find from Eq. (27), as well as the result $T_{\rm dip}^{\rm BEC}\sim E_{\rm bind}$ shown in the inset in Fig. 6, that the enhancement of $\eta(T)$ below $T_{\rm dip}^{\rm BEC}$ in Fig. 5(b) originates from the appearance of very long-lived molecular bosons.

As mentioned previously, although SCTMA ignores molecular scatterings by the effective interaction $U_{\rm B}$ in Eq. (26), these should actually contribute the molecular lifetime $\tau_{\rm B}$. Although inclusion of this effect remains as our future problem, we speculate that it would decrease the magnitude of $\eta(T)$ in the BEC regime, especially below $T_{\rm dip}^{\rm BEC}$, compared to the present SCTMA result.

Figure 10(a) shows $\eta$ as a function of the interaction strength. In this figure, we find that the minimum of $\eta$ is obtained, not in the unitarity limit, but slightly in the BEC side, being consistent with the recent experiment on a 6Li Fermi gas[40]. However, when we directly compare our result with this experiment, we should note that the experiment was done in a trap potential, under the condition,

[TABLE]

where $N$ is the total particle number, and $p(\bm{r})$ is the local pressure in a trap. Although our uniform calculation cannot completely describe this trapped geometry, noting that the left hand side of Eq. (28) represents three times the grand potential per particle, we plot in Fig. 10(b) the calculated shear viscosity in SCTMA, imposing the condition $3P/N=\varepsilon_{\mathrm{F}}$ at each interaction strength (where $P$ is the pressure). This figure shows that our result semi-quantitatively explain the observed $\eta(T)$, which takes the minimum at $(k_{\rm F}a_{s})^{-1}\sim 0.25$. This makes us expect that the minimum of the ratio $\eta/s$ may also exist away from the unitarity limit, although we need to calculate the entropy density $s$ to confirm this, which remains as our future problem.

4 Summary

To summarize, we have discussed the shear viscosity $\eta(T)$ in the BCS-BEC crossover regime of an ultracold Fermi gas. Including strong-coupling effects within the framework of the self-consistent $T$-matrix approximation (SCTMA), we have evaluated this transport coefficient in the normal state above the superfluid phase transition temperature $T_{\rm c}$.

In the weak-coupling BCS regime ((k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ -0.5), we found that $\eta(T)$ exhibits a dip structure, which is followed by a peak, as one decreases the temperature from $T_{\rm F}$. As the background physics of this, we pointed out that, while the former dip can be understood as the crossover from the classical regime to the Fermi degenerate regime, the latter peak is associated with the enhancement of pairing fluctuations near $T_{\rm c}$. Indeed, the latter peak temperature $T_{\rm peak}$ is closed to the pseudogap temperature $T^{*}$, which is defined as the temperature below which strong pairing fluctuations bring about a BCS-state like dip structure around the Fermi level in the single-particle density of states, in spite of the vanishing superfluid order parameter above $T_{\rm c}$.

This non-monotonic temperature dependence once disappears in the unitary regime (-0.5\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ (k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{<}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 0.7), where $\eta(T)$ monotonically decreases with decreasing the temperature down to $T_{\rm c}$. This behavior was shown to be consistent with the recent experiment on a 6Li unitary Fermi gas.

As one enters the strong-coupling BEC regime ((k_{\rm F}a_{s})^{-1}\ \raise 1.29167pt\hbox{>}\kern-8.00003pt\lower 3.01385pt\hbox{\sim}\ 0.7), $\eta(T)$ again exhibits a monotonic dip structure as a function of the temperature. In this case, the dip temperature $T_{\rm dip}^{\rm BEC}$ is close to the binding energy $E_{\rm bind}$ of a two-body bound molecules. This indicates that the rapid increase of $\eta(T)$ below $T_{\rm dip}^{\rm BEC}$ is due to the appearance of long-lived molecular bosons overwhelming thermal dissociations.

The shear viscosity diverges in both the BCS and BEC limits, so that $\eta$ is expected to become minimum in the unitary regime. We have numerically confirmed this expectation; however, the minimum of $\eta(T)$ is obtained, not at the unitarity, but slightly in the BEC side. This deviation from the unitarity limit is consistent with the observed interaction dependence of this quantity in a 6Li Fermi gas.

We finally note that the present work still has room for improvement. First, we could not calculate $\eta$ in the BEC regime near $T_{\rm c}$ because the numerical analytic continuation by the Padé approximation did not work well there. To avoid this difficulty, an alternative approach which does not need the analytic continuation might be useful. Regarding this, the real-time formalism based on the Keldysh Green’s function is promising. Second, SCTMA only includes effect of the molecular interaction in the BEC regime within the mean-field level, so that the molecular lifetime coming from inter-molecular scatterings is ignored. To improve this, at least, the second-order self-energy diagram in terms of this effective interaction must be taken into account. It is an interesting future challenge to explore how to incorporate such molecule-molecule scatterings into the present approach. Since we have only treated the normal state, the extension to the superfluid phase also remains as our future challenge. Since the shear viscosity has recently attracted much attention in various fields in connection to the KSS conjecture, our results would contribute to the further understanding of strongly interacting fermions, from the viewpoint of transport properties.

Acknowledgements.

We thank H.Tajima and R.Hanai for discussions. This work was supported by KiPAS project in Keio University. D.K. was supported by KLL Ph. D. Program Research Grant from Keio University. Y.O. was supported by a Grant-in-aid for Scientific Research from MEXT and JSPS in Japan (No.JP18K11345, No.JP18H05406, and No.JP19K03689).

Appendix A Expression for $\hat{\Pi}_{x,y}$ in the BCS model

The stress tensor operator $\hat{\Pi}_{i,j}$ ($i,j=x,y,z$) is defined from the conservation law in terms of the momentum flux,

[TABLE]

where ${\hat{\bm{J}}}=\sum_{\sigma}[\psi_{\sigma}^{{\dagger}}(\nabla\psi_{\sigma})-(\nabla\psi_{\sigma}^{{\dagger}})\psi_{\sigma}]/(2i)$ is the mass current operator, with $\psi_{\sigma}$ being the Fermi field operator. For the ordinary Hamiltonian consisting of the kinetic term and a short-range interaction potential $V(r_{12})\equiv V(|{\bf r}_{1}-{\bf r}_{2}|)$, one has[93],

[TABLE]

To obtain the detailed expression for $r_{12}(\partial V(r_{12})/\partial r_{12})$ in the present BCS model in Eq. (2), we evaluate the pressure $P$ by using the Hellmann-Feynman theorem,

[TABLE]

We assume that a gas is confined in a box with the size $V=L\times L\times L$, and the periodic boundary condition is imposed. In this case, the momentum ${\bm{p}}$ is quantized as ${\bm{p}}=2\pi{\bm{n}}/L$ (where $\bm{n}=(n_{x},n_{y},n_{z})$ is the quantum number), which depends on the system size $L$. In addition, the cutoff momentum $p_{\rm c}=2\pi|n_{\rm c}|/L$ also depends on $L$ (where $n_{\rm c}$ is a large quantum number). Keeping these in mind, varying the system size as $L\to L(1+\epsilon)$, we obtain the pressure from Eq. (31) as,

[TABLE]

We briefly note that, when we use Eq. (3) to remove the cutoff $p_{\rm c}$, Eq. (32) is found to reproduce the Tan’s pressure relation [95, 96, 94], $PV=2E/3+C/(12\pi ma_{s})$, where $E$ is the internal energy, and $C=(m^{2}U^{2}/V)\sum_{\bm{p},\bm{p}^{\prime},\bm{q}}\left\langle c^{{\dagger}}_{\bm{p}+\bm{q},\uparrow}c^{{\dagger}}_{\bm{p}^{\prime}-\bm{q},\downarrow}c_{\bm{p}^{\prime},\downarrow}c_{\bm{p},\uparrow}\right\rangle$ is the Tan’s contact[97].

Comparing Eq. (32) with the relation between the pressure and the stress tensor,

[TABLE]

one finds,

[TABLE]

Substituting this into Eq. (30), we reach

[TABLE]

We note that, although $\hat{\Pi}_{i,j}^{\prime}\equiv\hat{\Pi}_{i,j}+\delta\hat{\Pi}_{i,j}$ with any symmetric tensor $\delta\hat{\Pi}_{i,j}$ satisfying $\sum_{j}\partial_{j}\delta\hat{\Pi}_{i,j}=0$ also satisfies Eq. (29), this non-uniqueness doesn’t affect the shear viscosity[65]. We also note that the interaction term in Eq. (35) is for the case when Eq. (3) is used. When one employs another relation between $-U$ and a renormalized quantity, e.g., the dimensional regularization, the interaction term in Eq. (35) is replaced by the corresponding form. The expression for th shear stress tensor in the case of the dimensional regularization is given in Ref. [98].

Appendix B Shear viscosity obtained from the Boltzmann equation

This appendix summarizes how to calculate $\eta$ in the Boltzmann-equation approach. For a two-component Fermi gas, the Boltzmann equation is given by

[TABLE]

where $f(\bm{r},\bm{p},t)$ is a quasi-classical distribution function, and

[TABLE]

is the collision integral, with $d\sigma/d\Omega$ being the differential cross section, $f_{1}=f(\bm{r},\bm{p}_{1},t)$, $f^{\prime}=f(\bm{r},\bm{p}^{\prime},t)$, and $f^{\prime}_{1}=f(\bm{r},\bm{p}^{\prime}_{1},t)$.

To evaluate the shear viscosity $\eta$, we consider the current-flowing steady state, where the fluid velocity $\bm{u}(\bm{r})=(u_{x}(y),0,0)$, $u_{x}$, as well as the gradient $\partial u_{x}/\partial y$, are assumed to be small. In this case, the shear viscosity $\eta$ is related to the shear stress $\Pi_{x,y}$ as,

[TABLE]

$\Pi_{x,y}$ is also related to the quasi-classical distribution function as,

[TABLE]

To calculate Eq. (39), we divide the distribution function $f=f_{\rm loc}+\delta f$ into the sum of the local equilibrium distribution function,

[TABLE]

and the deviation $\delta f$ from it. Following the standard technique, we parametrize the latter as

[TABLE]

Substituting $f=f_{\rm loc}+\delta f$ into Eq. (39), and using Eq. (38), we obtain[59, 62, 66, 60, 61, 63]

[TABLE]

The function $\Phi(\bm{p})$ is obtained from the linearized Boltzmann equation, which is derived by substituting $f=f_{\rm loc}+\delta f$ into Eq. (36) and Eq. (37), and retaining terms with $O(\partial u_{x}/\partial y)$. The result is

[TABLE]

Here, $f^{0}=f_{\rm F}(\xi_{\bm{p}})$, $f_{1}^{0}=f_{\rm F}(\xi_{\bm{p}_{1}})$, ${f^{0}}^{\prime}=f_{\rm F}(\xi_{\bm{p}^{\prime}})$, ${f^{0}_{1}}^{\prime}=f_{\rm F}(\xi_{\bm{p}^{\prime}_{1}})$.

Following the previous works[59, 62, 66, 60, 61, 63], we set $\Phi(\bm{p})=Ap_{x}p_{y}/m$, where $A$ is a constant. Substituting this into the linearized Boltzmann equation (43), and multiplying $\sum_{\bm{p}}p_{x}p_{y}/m$ in both the sides of this equation, we have

[TABLE]

Once $A$ is determined, $\eta$ is immediately obtained as

[TABLE]

In Fig. 5, ‘BE’ shows the result obtained from Eq. (45), where the SCTMA chemical potential $\mu(T)$ in Fig. 2 is used.

In 40K and 6Li Fermi gases, the effective range of the Feshbach-induced tunable interaction is negligibly small. In this case, the differential cross section has the resonance form[59, 62],

[TABLE]

where $k=|\bm{p}-\bm{p}_{1}|/2$. In the classical regime ($T\gg T_{\rm F}$), one obtains

[TABLE]

where

[TABLE]

Equation (47) is reduced to[59, 62]

[TABLE]

In the Fermi degenerate regime ($T\ll T_{\rm F}$), in the weak coupling BCS regime, Eq. (45) gives,[59, 62]

[TABLE]

We briefly explain the background physics of these results in the classical and quantum regime. In the simple expression $\eta\sim nl_{\rm mfp}\bar{p}$[88], the mean free path $l_{\rm mfp}\sim 1/(n\sigma)$[88] is related to the collisional cross section $\sigma$ (where $n$ is the number density), and the averaged momentum ${\bar{p}}$ is estimated as ${\bar{p}}\sim\sqrt{mT}$. When the thermal de-Broglie wave length $\lambda_{T}=1/\sqrt{mT}$ is much larger than the scattering length $|a_{s}|$, the area where collisions occur is estimated as $\sim a_{s}^{2}$, giving $\sigma\propto a_{s}^{2}$. In this case, one obtains $\eta\propto T^{1/2}$. When $\lambda_{T}\ll|a_{s}|$, on the other hand, the area where collisions occur is estimated as $\sim\lambda_{T}^{2}$, because collision cannot occur when wave-functions of particles don’t overlap. In this case, one has $\sigma\propto 1/(mT)$, leading to $\eta\propto T^{1.5}$.

When we still use the expression $\eta\sim nl_{\rm mfp}\bar{p}$ in the Fermi degenerate regime, one may take ${\bar{p}}\sim k_{\rm F}$, and $l_{\rm mfp}\sim v_{\rm F}\tau\propto T^{-2}$, where $v_{\rm F}=k_{\rm F}/m$ is the Fermi velocity, and $\tau\propto T^{2}$ is the lifetime of quasi-particles originating from a particle-particle interaction. These give $\eta(T)\propto T^{-2}$, as expected.

Appendix C Derivation of Eq. (27)

The bare molecular three-point vertex function $S_{0;x,y}$ in Fig. 8(c) consists of three Fermi Green’s functions as

[TABLE]

We evaluate Eq. (52) in the BEC regime within the replacement of the SCTMA dressed Green’s function $G$ by the bare one $G_{0}$. The result is, at $\nu_{m}=\nu_{l}=0$,

[TABLE]

Approximately using Eq. (53) for $S_{0;x,y}$ in the shear-stress response function $\Xi(i\nu_{l})$ in Fig. 8(e)[89], we obtain

[TABLE]

where the molecular Bose Green’s function has the form,

[TABLE]

Carrying out the $\nu_{m}$-summation, we execute the analytic continuation $i\nu_{l}\to\omega_{+}$ in Eq. (54). Then, we obtain from Eq. (10),

[TABLE]

Using the fact that $\tau_{\rm B}=\lambda e^{E_{\rm bind}/(2T)}$ is very large because of the exponential factor, we further employ the approximation,

[TABLE]

Substitution of this into Eq. (56) gives Eq. (27).