Strong pairing in two dimensions: pseudogaps, domes, and other implications

This paper addresses the transition from the normal to the superfluid state in strongly correlated two dimensional fermionic superconductors and Fermi gases. We arrive at the Berezinskii–Kosterlitz–Thouless (BKT) temperature TBKT as a function of attractive pairing strength by associating it with the onset of ‘quasi-condensation’ in the normal phase. Our approach builds on a criterion for determining the BKT transition temperature for atomic gases which is based on a well established quantum Monte Carlo analysis of the phase space density. This latter quantity, when derived from BCS–BEC crossover theory for fermions, leads to non-monotonic behavior for TBKT as a function of the attractive interaction or inverse scattering length. In Fermi gases, this implies a robust superconducting dome followed by a long tail from the flat BEC asymptote, rather similar to what is observed experimentally. For lattice systems we find that TBKT has an absolute maximum of the order of 0.1EF. We discuss how our results compare with those derived from the Nelson–Kosterlitz criterion based on the mean field superfluid density and the approach to the transition from below. While there is agreement in the strict mean-field BCS regime at weak coupling, we find that at moderate pairing strength bosonic excitations cause a substantial increase in TBKT followed by an often dramatic decrease before the system enters the BEC regime.


Introduction
Recently there has been a resurgence of interest in superconductivity in (quasi-)2D materials. This has been driven by exciting discoveries of novel superconductors such as magic-angle twisted bilayer graphene [1], FeSe monolayers [2,3] and transition metal dichalcogenides [4][5][6]. Many of these and other interesting superconductors [1,[7][8][9] appear to belong to the more strongly correlated class which is distinct from BCS-Eliashberg superconductors and can be argued [10] to be intermediate between the BCS and Bose-Einstein condensation (BEC) limits. The challenge then is to develop an understanding of strongly correlated superconductivity in two dimensions where the long-range superconducting instability is replaced by a Berezinskii-Kosterlitz-Thouless (BKT) transition [11,12]. Meeting this challenge is essential: an in-depth understanding of these quasi-2D superconductors, requires that we abandon the predictions of BCS theory. At issue, also is whether strict BCS theory is appropriate for computing even the superfluid stiffness; one might expect that this should be obtained by including contributions of preformed pairs, not present in BCS theory, at the BKT transition temperature.
Arriving at this formalism is the goal of this paper which addresses BKT superconductivity in the presence of strong pairing correlations. Our attention is on the approach and calculation of T BKT from the normal state, following the extensive body of work on BKT in atomic Bose systems [13]. This is complementary to the research which addresses T BKT from the superfluid side [14][15][16][17][18][19]. In a seminal work [20], Halperin and Nelson have used a fluctuation approach to address the physics of approaching the transition from above. We argue in this paper that, in line with their thinking and with reference [21], the normal state in question should reflect stronger pairing correlations, particularly those that lead to a stable, observable 'pseudogap'.
We stress here that understanding BKT in fermionic systems is not as straightforward as in their bosonic counterparts. Indeed the experimental realization of the BKT model was established in superfluid helium films [22] many years ago. There is also a convincing case for the observation of BKT in atomic Bose gases [13]. The nature of the transition and whether or not it is present in superconducting films has been a subject of debate [23][24][25][26]. For this reason it is important to pursue a number of different approaches which address fermionic BKT. This provides the underlying motivation for our paper and leads us to study the transition when approached from the normal phase. We do so following the methodology introduced for atomic gases [13,27,28].
In determining physical variables and consequences, a notable complication is that plots of T BKT as a function of the attractive pairing interaction strength g are non-monotonic, so that knowing T BKT does not uniquely determine other fundamental properties. Indeed, quantum Monte Carlo (QMC) simulations [29,30] and other more analytic calculations for the case of a lattice dispersion [15], show plots of T BKT vs g which exhibit a superconducting 'dome' shape. It is generally argued [29] that this dome lies just beneath the intersection of two curves: an increasing trend on the BCS side and a decreasing contribution at larger g representing the BEC asymptotics, as is shown schematically in figure 1(a). Similar arguments are presented for the case of a Fermi gas, except that the BEC limit, rather than decreasing with increasing coupling constant g, reaches a constant asymptote, as shown in figure 1(b).
A central result of this paper is that when the instability is approached from the normal state, we, too, find robust domes for the lattice dispersion, and in addition we find they are present as well for the case of a gas dispersion. Importantly, these non-monotonicities appear in the intermediate coupling regime, away from the BEC regime. Here a Fermi surface is still present, as one would expect in any physical 2D superconducting system. The dome arises from a competition between a rising trend of BCS pairing on the BCS side and the strong suppression of T c due to formation of pairs and the concurrent onset of a pseudogap, well before the BEC regime is accessed.
The approach in this paper is to combine a pairing fluctuation theory [31][32][33] for BCS-BEC crossover with a description [13,34,35] of BKT for bosons, and, thereby, establish a 2D BCS-BEC crossover theory with a finite transition temperature. Because the BKT criterion approaches the instability from the normal state [13,34,35] it reflects the phenomenon of quasi-condensation or 'presuperfluidity' [28,36,37]. From a theoretical perspective quasi-condensation appears when the bosonic chemical potential becomes sufficiently small, but non-zero. This leads to a large number of non-condensed pairs having very small (but not strictly zero) momentum. When approaching the transition from above, it is found [13,34,35] that the BKT temperature depends on the ratio of the effective pair density, n B (T), representing the areal number density of bosons, to their effective mass, M B (T). The transition occurs when this ratio (which is proportional to the bosonic phase space density) reaches a critical value established from previous, (lattice) QMC calculations [34].
We emphasize that these bosons are a composite made up of fermions and the fundamental bosonic variables n B and M B must depend on the fermionic excitation gap Δ(T) (which is non-vanishing even at T BKT ). In a true Bose system n B and M B are fixed in temperature. In fermionic superfluids both depend on T and interaction strength. The non-monotonicity we observe depends on a competition between n B (T), which increases, and 1/M B (T), which decreases with increasing g. That the effective number of fermion pairs increases as the attraction becomes stronger should be clear 5 . The pair mass, on the other hand, can be understood as reflecting the inverse square of the pair size; the mass is light in strict BCS and becomes heavier with increased attraction while still remaining within the fermionic regime.
There is a large body of work which addresses the BKT transition as approached from below using the Nelson-Kosterlitz [38] condition. This criterion depends on a presumed form for the superfluid phase stiffness ρ s , which is given by the ratio of superfluid density n s over the fermion mass m. In earlier calculations it was often assumed that this quantity can be calculated at the mean field level generalized to include stronger g, via a mean field treatment of the crossover from BCS through BEC. This theory presumes that the destruction of the superfluid stiffness derives entirely from fermionic degrees of freedom. Application of this mean field picture for the case of a lattice dispersion [15] or a Fermi gas [14,17] generally leads to plots similar to figures 1(a) and (b).
Clearly, both theoretical approaches (using either the superfluid density or quasi-condensation) need to be simultaneously pursued by the community if progress is to be made. It should be cautioned, however, that neither of these two schemes explicitly accommodates the important vortex-antivortex excitations which presumably affect the size of T BKT , although how much is not precisely known. Importantly, in this paper we discuss the relation between the two, and demonstrate agreement at the weak coupling, BCS level. However, as the attractive coupling constant g is increased in magnitude, bosonic excitations become more significant. Through our comparison we are able to characterize these bosonic contributions; these also turn out to be non-monotonic, causing an increase in T BKT in the near 'unitary' regime and a decrease very close to the onset of the BEC regime.
We stress that our approach which is based on the onset of quasi-condensation is more directly connected to those experiments where BKT is most clearly observed as in 2D Bose superfluids [13,34,35] and in 2D Fermi superfluids [27,39] as well. A quasi-condensation approach presumes that the inter-boson interactions are sufficiently weak. In fermionic systems, while there may be strong inter-fermion interactions, the inter-boson correlations inherent in a BCS-like ground state are not presumed to be large. Indeed, in the BEC regime, the inter-pair interaction becomes progressively weaker as the inter-fermion attraction becomes stronger.
Finally, we end this section by noting that other consequences of strong pairing correlations should be a central feature of the normal state. Indeed, the foundation for using phase only (XY) models in 2D systems depends on having a substantial pairing at the transition temperature [21]. Thus, one should characterize a given superconductor by the pair of temperatures T BKT and the pseudogap onset temperature T * , which then removes the ambiguity associated with the non-monotonicity in the transition temperature. Pseudogap effects are enhanced in 2D systems and have been clearly observed in 2D atomic Fermi gases [39,40]. For a transition temperature of, for example 0.08E F , which is rather strong coupling, we find that the pseudogap onset temperature is about twice this temperature.

Outline
We now present an outline of the remaining sections of this paper. Section 2 of the paper presents a brief review of our BCS-BEC crossover theory based on a self-consistent T-matrix approximation. The goal of this discussion is to show how to obtain the important bosonic quantities n B and M B for general g from their fermionic counterparts.
In section 3 we discuss the case of 2D superfluids and present an expression for the BKT transition temperature which is widely used in the bosonic literature [13,35]. As in references [27,28], we show how to apply it to fermionic superfluids (with both lattice and continuum dispersions). Contrasting with this 'quasi-condensation' approach to the BKT temperature, is the more widely used criterion based on the superfluid phase stiffness, ρ s , as discussed in section 4. We present comparison plots of T BKT in the two approaches from above and below the transition. These are indistinguishable in the weak coupling BCS regime. However, at moderate or strong coupling, bosonic contributions, which are absent in the mean-field ρ s approach, become increasingly more important. By comparing these two schemes, we are able to characterize and quantify these bosonic contributions which are interestingly non-monotonic as a function of increasing g.
Reasonable quantitative comparisons with Fermi gas experiments are presented in section 5, along with predictions for the behavior in the lattice case. Comparisons with QMC results on the attractive Hubbard model indicate some deviation (roughly within a factor of 2). We show how to associate the measured transition temperature with other attendant properties such as the size of the pseudogap. Following a discussion in section 6, our conclusions are presented in section 7.

Theory of BCS-BEC crossover
The incorporation of stronger than BCS pairing attraction is an important component in addressing BKT in superconductors. In 2D the s-wave pairing instability is always accompanied by two-body bound states [41] suggesting that pairing effects are amplified. We thus invoke an extended form of BCS theory associated with a BCS-BEC crossover. This represents an analytic approach to the attractive Hubbard model. Just beyond the BCS regime, we show that the transition temperature for three dimensional systems is to be associated with a BEC-like condensation of preformed pairs. In two dimensions, T c ≡ 0, but, in a similar way, we argue quasi-condensation of preformed pairs depends on a very similar set of parameters which appear in both 2D and 3D: n B and M B .
Our approach is intimately tied to the observation [42] that the BCS ground state wavefunction has a greater applicability than had been appreciated at the time of its original proposal. Here u k and v k are the BCS coherence factors. This BCS ground state shows that the bosons do not interact directly, but because they are composed of fermions, there is a non-zero scattering length associated with indirect boson-boson interactions, via the fermions.
As the strength of the attractive pairing interaction g between fermions is increased, this wavefunction is capable of describing a continuous evolution from BCS-like behavior to a form of Bose-Einstein condensation. The BEC regime sets in when the fermionic chemical potential μ becomes negative and at this point, the underlying Fermi surface disappears.
We extend this ground state to finite temperature following an approach to BCS theory proposed by Kadanoff and Martin [31,43]. While in strict BCS theory all bosonic degrees of freedom (fermion pairs) appear only at (and below) the transition as condensed pairs, with stronger attraction non-condensed pairs are present. They are accompanied by the existence of a normal state pairing gap or pseudogap, which is particularly pronounced for a 2D system. For the purposes of our BKT analysis what will be important is to quantify the number of these preformed pairs n B and their effective mass M B .
To do so we can associate the underlying structure of BCS theory with a two particle propagator, called t pg , which is given by where G 0 is the bare fermion Green's function, and the dressed Green's function G assumes the BCS form with the Bogoliubov quasiparticle dispersion E k = ξ 2 k + Δ 2 . We define ξ k = k − μ. A detailed derivation of equation (2) can be found in reference [44]. Analytically continuing equation (2) [for four-vector Q ≡ (q, Ω) = 0] leads to Then, as a generalization of the usual Thouless condition, we recognize that the statement effectively leads to the usual BCS temperature dependent gap equation for Δ(T). We can think of this familiar gap equation as a BEC condition that the pairs associated with t pg (which are necessarily non-condensed) have zero chemical potential μ pair = 0 for all T T c . For the most part our interest will be on the long wavelength and low frequency limit, where the pair propagator can be approximated as Here a −1 0 characterizes the pair fluctuation strength, γ is the decay rate due to the two-fermion continuum, and M B is the effective pair mass. Both a 0 and M B can be determined via Taylor expansion of t −1 pg (Q) in Ω and q. In particular, a 0 Δ 2 = n/2 − k f (ξ k ) 6 . In the presence of a quasiparticle excitation gap, the pair decay rate at low frequencies vanishes. It is small compared to Ω q = q 2 /(2M B ) when finite momentum pairs make a significant contribution to the self energy (away from the BCS limit). From now on we omit iγ in the pair propagator. We should note that the pair mass M B is now accessible through equation (4). Next we focus on these non-condensed pairs in the normal state [31,32], where the pairs have non-zero chemical potential μ pair which smoothly vanishes at the transition into the ordered phase. Here we identify the pairing gap Δ with the pseudogap so that Δ ≡ Δ pg . This excitation gap is to be distinguished from the order parameter. The self consistency condition can be written as t −1 pg (0) = a 0 μ pair . In two dimensions μ pair (T) will be shown to assume small values, but never reach zero, except at T = 0. By contrast, in three dimensions μ pair (T) vanishes at and below a finite T c .
To obtain n B we note that the self energy associated with the dressed Green's function is more completely given by where in this last step we have assumed that the system is near an instability where t pg (Q) is strongly peaked at Q = 0. Equation (5) is a standard approximation in the cuprate literature for the pseudogap-related self-energy [33,45]. We stress that this second line in equation (5) is the only approximation used here, aside from the overarching assumption implicit in equation (1) that we are dealing with a BCS-like gap equation and ground state, importantly extended to BCS-BEC crossover. Note that this approximation effectively ignores Hartree as well as incoherent contributions to the fermionic self-energy, which may, as well, introduce particle-hole asymmetry effects. There is an additional complication (for the lattice case) near half filling associated with competing charge density wave order in the particle-hole channel [29]. For simplicity, we ignore this here. It is, however, advantageous to adopt the approximation in equation (5) in the vicinity of small μ pair (which is appropriate near T BKT ) for analytical tractability.
Combined with the parametrization in equation (4), we derive the following self-consistent equations for a fixed-density system consisting of fermionic and bosonic quasi-particles: where b(x) is the Bose-Einstein distribution function. Here a 0 and M B depend on the three parameters μ, Δ, T, and can be deduced through Taylor expansions. For a 2D system, these equations can be solved self-consistently for (Δ, μ, μ pair ) at low T and for given interaction strength g. The zero T solution can be taken as the limit of T → 0 so that μ pair remains finite in equation (7) 7 . We emphasize that both n B and M B are a function of temperature, and should be determined self-consistently via equations (6)-(8) when solving for T BKT .
What should be clear from this analysis is that in BCS-BEC crossover theory the normal state consists effectively of an admixture of fermions (with number n − 2n B and chemical potential μ) and bosons (with number n B and chemical potential μ pair ). Equation (7) appears physically reasonable in establishing the direct correspondence between the number of bosons and the energy scale Δ for binding fermions.

Behavior of the 3D transition temperature: hints about 2D
It is useful to present a few analytic results from this formalism. The 3D transition temperature for a gas dispersion is associated with the condition μ pair = 0 at T c . This enters in the boson number equation (7), 6 Expressions for M B in various situations can be found in reference [45]. 7 Alternatively, at very low T, where μ pair is very small, one can set it to zero in equation (6) so that equations (6) and (8) reduce to the BCS-Leggett mean-field equations [43], which can be solved for (μ, Δ). This then also fixes the value of a 0 and M B and thus n B . Finally, one determines μ pair as a function of (low) T via equation (7). While this alternative procedure is an approximation at nonzero T, it becomes exact in the T → 0 limit, where μ pair necessarily vanishes. Therefore, we conclude that at T = 0, the pair density n B is completely determined by the mean-field solution of the ground state, and so is (M B ). and after some algebra, leads to where both n B and M B are temperature dependent and calculated at T c . In a very compact way this equation encapsulates the behavior of BCS-BEC crossover theory, beyond the strict weak coupling limit. It should be viewed as reflecting the condensation temperature of preformed pairs. Importantly, these represent the emergent bosons which are central to a treatment of BCS-BEC crossover. We note that the dependence on n B is similar to what is found in an ideal Bose gas, but it should be stressed here that inter-boson interactions are present, as is reflected in the superfluidity [46] and in the collective modes [47,48] of BCS-BEC systems. Inter-boson interactions are associated with both the pairing interaction and the Pauli repulsion of the underlying fermionic constituents in the pairs. If, instead, one considers a 2D system, by analogy the associated number to mass ratio which determines the transition, the transition temperature might be expected to be n B /M B where n B (T), represents now the areal number density of bosons, and M B (T), their effective mass. We show in the next section that this same ratio (known as the phase space density) appears in the BKT criterion applied by the atomic Bose gas community [13]. Here, however, the bosonic variables are temperature dependent and depend on the fermionic excitation gap Δ(T) and chemical potential μ.
This ratio n B /M B and its 3D analogue determine the shape of the transition curves as a function of g. Indeed, the fractional power n 2/3 B is not very different quantitatively from n B , away from the n B → 0 limit. In this way, we will see that the shape of the curves in 2D BKT are not too dissimilar from their 3D counterparts.
To elucidate the physics, it is useful to present in figure 2 an anticipatory plot of T BKT as a function of coupling constant in a way which serves to identify the boson and fermion constituents. What we indicate in figure 2 is the relative admixture of broken pairs (fermions) and pairs (bosons) as the interaction strength is continuously varied. The small boxes in figure 2 should be viewed as representative 'cartoons' which characterize this pseudogap phase. A very small transition temperature is expected when the boson number is almost zero, as shown in the low g regime. The largest T BKT is found in an intermediate state consisting of bosonic and fermionic quasi-particles. At the highest value for g, all signs of the fermionic constituents are gone and the transition begins to approach zero as t 2 /g. The pseudogap is present whenever there are a finite number of pairs at the transition temperature; it becomes progressively larger, the larger the number of pairs.

BKT criterion as approached from the normal state: BKT in atomic gases
In the next two sections we discuss two types of criteria which have been used to establish T BKT and follow this with a comparison. It is useful to consider equation (7) next for the strictly two dimensional case where there is no true condensate, away from the ground state. This equation can be inverted exactly to give the pair chemical potential: The size of |μ pair |, which measures how close the normal fluid is to a long range-ordered superfluid phase, reflects the bosonic phase-space density: D B (T) ≡ n B (T)λ 2 B , where λ B ≡ 2π/M B T is the de Broglie thermal wavelength for the pairs. In this notation k B = 1 and = 1. Importantly, in two dimensions, D B determines the pair chemical potential μ pair , so that there is quasi-condensation [28,36] when D B is sufficiently large or |μ pair | is sufficiently small.
When approached from the high temperature side [13,49], the bosonic BKT transition is known to occur [34] when the temperature dependent phase space density reaches a critical value where the dimensionless coupling constantg reflects the size of the 3D inter-boson scattering length a B , along with the 2D localization length. The constant C ≈ 380 has been established by QMC [34], based on a tight binding lattice, but quite generally argued to be universal. If one parameterizes the 2D confinement by a trap of frequency ω 0 , it follows that [35]g = a B 8πM B ω 0 / . Estimates of D crit B for fermionic superfluids are available in the literature [27,28]. The values range from around 4.9 to 6.45. This can be compared with the counterparts in atomic Bose gases which are typically [37] around 8. Here, our approach is not to do fine tuning, but rather to associate D crit B (independent of the pairing strength) with the value which has been shown to best fit the data on Fermi gases [28]. Notably, this lowest value, D crit B = 4.9, in the range, is closest to the factor 4.0 in the usual BKT relation. These are not order of magnitude variations and the uncertainty does not significantly affect the shape of the curves for T BKT vs g; however, it does affect somewhat their position on the vertical axis. While one would imagineg varying with inter-fermion interaction strength, its effect on D crit B is necessarily weak through a logarithmic dependence.
Thus, based on the atomic Bose [13,35,37] and Fermi gas literature [28] we apply the BKT criterion at T = T BKT . In the above equation n B and M B reflect the fermionic degrees of freedom, through the pairing gap Δ(T).

BKT criterion derived from superfluid density: Nelson-Kosterlitz condition
When approached from the low temperature side, the BKT transition [16,50] occurs at a universal value of the bosonic superfluid density. The transition temperature can be defined [38] in terms of the superfluid component of D B such that D s B (T BKT ) = 4.
To tie the two approaches together, we can, however, extract an inequality This reflects the fact that the total bosonic phase space density must exceed its superfluid counterpart, given in equation (12) 8 . If we use BCS-BEC mean field theory [15] to evaluate the superfluid phase stiffness ρ MF s , equation (12) is equivalent to the condition that at T BKT , the superfluid density satisfies

Comparison of the two BKT criteria
Of central importance is to compare these two schemes for the BKT transition temperature obtained when approached from the normal state using equation (11) or alternatively using equation (14). The detailed numerical results in this section and the next are based on equations (6)-(8) along with equation (11). This comparison is presented in figure 3 for the case of a lattice dispersion and in figure 4 for the Fermi gas case.  1 4 ρ MF s calculation. The difference between the two curves reflects bosonic contributions to the destruction of superfluid phase stiffness which are not present in the mean field approach. Panel (b) shows the two components in the (more bosonic) quasi-condensation picture n B and M B . This figure indicates that it is a suppression of the pair mass at moderate interaction strength which leads to an enhanced maximum in the ratios, plotted above. Here E F is taken to be the non-interacting Fermi energy, with E F ≈ 0.604t, and we take the lattice constant to be unity. . We see that the agreement is very good in the strict BCS regime, for small coupling g. This largely derives from the fact that here T BKT is close to the pairing onset temperature T * . Both theories yield the same T * . We can refer to figure 2 to see that in this regime the two BKT criteria yield equivalent results, as the only quasi-particles in the normal state are fermionic and both are associated with the same pairing onset T * .
A central difference is that the quasi-condensation approach leads to a higher maximum at intermediate coupling and a more dramatic plummet in T BKT beyond the maximum. A slight kink appears exactly when the fermionic chemical potential reaches zero and one might expect this feature as n B has to have a discontinuity in slope. Here the boson number density is precisely half the fermion density and all fermions are paired. Asymptotically, in the BEC regime, the two curves also coincide as expected 9 .
Understanding the physical origin of these two principal differences is particularly important for arriving at a more complete physical picture of the BKT transition in a fermionic system. Both of these effects arise from the bosonic contributions to the phase stiffness which are missing in the mean field approximation to ρ s . We refer to the lower panels of figures 3 and 4 to help understand this behavior.
Plotted in these lower panels are the two components in the quasi-condensation picture n B and M B , with a rescaling for better visibility. This rescaling will not affect the deduced ratio plot (except for an overall normalization). The origin of the important non-monotonic effects in the ratio can now be seen. More specifically, we see that n B and M B rise just beyond the BCS regime where pseudogap effects associated with meta-stable non-condensed pairs begin to emerge. Notably, the pair mass increases more slowly (in the plots) giving rise to the maximum in the ratio, which overshoots the ρ s -based plot 10 . It should be noted that the region of this overshoot is just beyond the BCS regime, where fluctuations start to  1 4 ρ MF s in the two approaches for the gas case. The temperature used throughout this figure is taken to be the critical BKT temperature for the 1 4 ρ MF s calculation. The difference between the two curves reflects bosonic contributions to the destruction of stiffness, absent in the mean field approach. Panel (b) shows the two components n B and M B . This figure indicates that it is a suppression of the pair mass at moderate interaction strength which leads to an enhanced maximum in the ratios. be important. Here, however, we argue that these fluctuations are not strictly phase fluctuations but contain amplitude contributions as well. Thus, the physics is not solely characterized by ρ s , as would be case if the amplitude variations were frozen.
More specifically, the origin of this slower rise in M B appears to derive from an increased stability of non-condensed pairs which is associated with the onset of the pseudogap. Stabilization arises because the presence of a pseudogap means that there is an energy cost, inhibiting the dissociation of pairs (into fermions). We can think of M B as very roughly representing the inverse square of the pair size. Hence a smaller pair mass reflects an increased coherence length of pairs. In this way the transition temperature exhibits a higher maximum T BKT .
At increasingly stronger coupling, the bosons contribute a second structural feature in the T BKT plots which appears as a downturn after the maximum, but before the BEC regime is reached. This result has been anticipated [19,51]: bosonic quasi-particles are expected to provide alternative mechanisms for exciting the condensate. Hence they lead to a reduction in the phase stiffness and related transition temperature.
We end by summarizing the essential points from this comparison, which apply to both the lattice and gas dispersion cases. The mean field ρ s approach is missing bosonic contributions 11 . In the strict BCS regime these can be neglected and in that regime the two calculations of T BKT are equivalent. (This equivalence is insured by the particular T matrix used here.) The most important consequence of including non-condensed pairs in the gas case is that they lead to a maximum at intermediate coupling. This derives from the extended stabilization of pairs and concomitant reduction in their mass. In the lattice case, for the same general reason, these pairs enhance an existing (weaker) maximum.
We end this section by noting that correlation functions have also been addressed [28,36] within this quasi-condensation approach. Theoretically we find that a screened algebraic decay best fits our numerically obtained results.

Numerical results
We can compare our theoretical framework directly to Fermi gas experiments [27,28] on trapped superfluids (although, in contrast to reference [36] trap effects have not been included). We should stress that in addition to this experimental work which addresses the temperature dependent phase diagram, there is a substantial literature on other aspects of 2D Fermi gases, both from an experimental [52][53][54] and theoretical [55][56][57] Unlike the curve shown figure 4 for the mean field scheme, the experiments at intermediate coupling exhibit a non-monotonic behavior. In particular when ln(k F a 2D ) ≈ 1, there is an enhancement of the critical temperature. While the extreme BCS limit is not apparent in these experiments, T BKT must ultimately reach zero at weak coupling, so there is a dome like feature [27] followed by a nearly constant BEC asymptote.
The left panel of figure 5 presents a direct comparison between our calculated T BKT and the experimental data in units of E F versus −ln(k F a 2D ), with E F being the non-interacting Fermi energy and k F the Fermi wave-vector. From right to left on the horizontal axis represents the transition from BCS-like to BEC. The theory curve (black solid) is overlaid on top of the data showing colored contours of the quasi-condensate fraction, N q /N. While the data points are incomplete, with large error bars, in the phase diagram, the overall agreement between theory and data is reasonably good. The dome structure in the data for both T BKT and N q /N is most apparent in the edge of the red and green contours, for ln(k F a 2D ) > −1. We find a kink near ln(k F a 2D ) = 0 where the fermionic chemical potential μ = 0. Both theory and experiment have to exhibit a decrease in the transition temperature toward the BCS limit. Beyond the dip which establishes the BEC regime, T BKT ≈ 0.1E F . (The calculated asymptotic value is slightly different from 1 8 by a factor of 4/4.9, since the critical value D crit B is slightly larger than 4.) This figure is consistent with the expected asymptotic values for n B = n/2 and M B = 2m.
We emphasize that by presenting this figure we are not claiming absolute agreement with experiment. (Although, perhaps surprisingly, within error bars, our theory curve passes through all but one data point with no adjustable parameters.) The experimental figure should be viewed as a relevant benchmark to help the community arrive at an understanding of BKT in fermionic superfluids, which is a rather unique case where there is rather systematic data. Notably, here we are dealing with greater complications than, for example, in a prototypical BKT system such as helium-4.
We replot in the right panel of figure 5 from reference [36] the theoretically calculated T BKT curve and the contours showing the quasi-condensate fraction, when trap effects are included. Evidently, these trap effects do not qualitatively affect the general behavior we report above 12 .
We now focus exclusively on the lattice case. Figure 6 provides a summary of our results at two representative electron densities. Panels (a)-(c) are characteristic of low electron density n = 0.3. As shown in (a), at weak to intermediate couplings, T BKT has a dome shape followed by a long slow tail. Each dome we find is accompanied by a dip where the chemical potential μ changes sign. The downturn of T BKT on the stronger coupling side of the dome is caused by the increasing contributions of pairing fluctuations due to increasing pairing strength. The increasing pairing gap reduces the Fermi level. In addition, these fluctuations lead to a growing pseudogap at and above T BKT , which depletes the density of states and thus suppresses T BKT . These two combined effects are so strong that T BKT starts to decrease in the intermediate pairing strength regime, before the Fermi surface shrinks to zero when μ = 0. Beyond this point, the Fermi surface is gone, so that all fermions are paired up.
Panel (b) shows how the above picture can be regarded as driven by a competition between an increase in the density of Cooper pairs n B (which saturates to n/2 above g c ) and an even stronger increase in the mass M B . Here the critical coupling g c is associated with the point where μ changes sign, as depicted in panel (c). For strong coupling g > g c , the normal state essentially consists purely of bosonic pairs without unbound fermions (except at the highest T). Note that the pair mass scales linearly with g. This gives rise to the expected asymptotic tail in T BKT ∝ t 2 /g. We emphasize here that the dome at intermediate couplings is not determined by the t 2 /g asymptotics seen in strong coupling. For completeness, in figure 6(a) we also present the temperature T * where the pseudogap sets in. Over most of the BKT dome, the magnitude of the gap Δ(T BKT ) at the transition temperature is essentially unchanged from its zero-temperature value. Panels (d)-(f) are representative results for high electron densities (here we use n = 0.7 for illustrative purposes). Just as in the previous case with n = 0.3, there is also a superconducting dome in the range of g/E F 8. In addition, the maximal transition temperature T BKT ∼ 0.1E F in both cases. However, a notable difference is that we do not find the long asymptotic tail as it is not possible to achieve a purely bosonic regime where all electrons bind into Cooper pairs. This is reflected in the fact that the fermionic chemical potential [panel (f)] never changes sign before T BKT reaches zero. This occurs concurrently with the vanishing of 1 M B , corresponding to Cooper pair localization 13 . The fact that the fermionic regime is so robust at high densities is intimately connected to the (near-) particle-hole symmetry of the underlying lattice Hamiltonian. In a bipartite lattice at exactly half-filling, the fermionic chemical potential is pinned at μ = 2dt (where d is the dimension), regardless of the interaction strength. As a result a purely bosonic regime can never be achieved.
Interestingly, within our approach, we observe re-entrant superconductivity in a narrow range of intermediate electron densities around n = 0.55. Here in addition to the dome for g < g c , there is a strong coupling tail with t 2 /g asymptotic behavior that sets in at a slightly larger g. Similar re-entrant behavior has been observed elsewhere [58].
We can compare to earlier QMC data on the attractive Hubbard model at n = 0.7 [30]. There it was found that the BKT transition temperature reaches a maximum of about 0.175t which occurs at g = 5t, as compared with the maximum we find of 0.33t which occurs at g ≈ 6.8t. (The T BKT calculated using the mean-field superfluid density yields a maximum of 0.24t around 4t, also larger than the QMC result.) The QMC data do not extend beyond g = 8t. It is likely that the self-energy based approximation [33,45] we Table 1. Estimates of physical quantities for the case T c /E F ≈ 0.08 based on our calculations for n = 0.3 (in top line) and n = 0.7 in bottom line. To convert to units of temperature, we assume T c = 1.5 K. Here Δ is the pairing gap at the BKT transition. We find g/E F = 1.87 and 1.06 for the low and high densities respectively. make as shown in equation (5) leads to an over estimate of particle-hole symmetry and may be in part responsible for the differences from the QMC data. Additionally, the absence of particle-hole fluctuations, as in generic T-matrix approaches, may lead to over estimates of the transition temperature and pairing gap [59][60][61]. Also important may be short-ranged charge density wave fluctuations which are neglected in the present study.

Discussion
We turn to table 1 for a more quantitative summary of the various energy and length scales in the intermediate coupling regime; here for a given ratio of T BKT /E F , there are two possible values of the coupling strength g/E F . For concreteness we choose the ratio to be 0.08, motivated by estimates made for twisted bilayer graphene (TBG) [1]. We want to firmly stress that this paper does not incorporate the band structure or other complexities of this material. (Also note that the maximum transition temperature of the Monte Carlo calculations [30] does not appear to be sufficiently large to reach this value.) 14 . Recent studies suggest that superconductivity in TBG is associated with at least two coupled bands and, moreover, inter-band contributions arising from fermionic excitation terms in the superfluid density [62] can be very important. Currently underway is an extension of the present work, to treat two band models and to determine how these additional terms inter-relate with bosonic contributions to the superfluid density. Both of these become progressively more important at strong coupling. Nonetheless, as in more conventional BCS theory, once one knows the transition temperature a number of additional properties can be quantified regardless of the underlying microscopic details. Of particular interest are the size of the pseudogap Δ at the transition in comparison to T BKT and the pairing onset temperature. The lower of the two g values appears most reasonable physically when compared to estimates in TBG [63]. In both cases the amplitude of Δ is relatively the same at T BKT and T = 0; notably, for the smaller g, the chemical potential is close to E F , so that the system is far from BEC. For this more likely situation, we note that the pairing onset temperature T * is roughly twice T BKT . When it differs significantly from T BKT , this is a crucially important parameter as it suggests (from figure 3) that this particular material is outside of the regime where the mean field ρ s approach is applicable. Rather bosonic excitations must be included.
This emphasizes that there are two important temperature scales: T * and T BKT . In general, it is the pair of temperatures [21] which provides full characterization of a given BKT system. If it is known that T c /E F ≈ 0.08 with T c ≈ 1.5 K then one can read off from the phase diagram we present, the size of the pairing gap (around 4-5 K) and the size of the pairing onset temperature: (T * ≈ 3 K).

Conclusions
We have stressed that understanding BKT in fermionic systems is not as straightforward as its bosonic counterpart. Indeed the experimental realization of the BKT model was established in superfluid helium films [22] many years ago. There is also a pretty convincing case for the observation of BKT in atomic Bose gases [13]. Whether or not this model applies to superconducting films has been a subject of debate [23][24][25][26]. For this reason it is important to pursue a number of different approaches for addressing fermionic BKT. We argue that this provides the underlying motivation for our paper. Here we study the transition when approached from the normal phase, following the methodology introduced for atomic gases [13,27,28].
An additional motivation for this paper is based on the excitement behind the recent discoveries of novel 2D superconductors which appears to be largely based on the hope that these (often engineered) systems can produce new forms of high temperature superconductivity. Also exciting is the possibility that they will serve to teach us about mechanisms for existing high T c (say, cuprate) systems.
Our paper argues for a somewhat more modest perspective. Independent of the specifics of the attractive interaction mechanism, in these 2D systems, there is an absolute maximum to the transition temperature T BKT . It can be rather high, say of the order of 0.1E F as found here, or somewhat lower (0.05E F as found in Monte Carlo [30]), but it does ultimately set an important limit.