Heuristic bounds on superconductivity and how to exceed them

What limits the value of the superconducting transition temperature ($T_c$) is a question of great fundamental and practical importance. Various heuristic upper bounds on $T_c$ have been proposed, expressed as fractions of the Fermi temperature, $T_F$, the zero-temperature superfluid stiffness, $\rho_s(0)$, or a characteristic Debye frequency, $\omega_0$. We show that while these bounds are physically motivated and are certainly useful in many relevant situations, none of them serve as a fundamental bound on $T_c$. To demonstrate this, we provide explicit models where $T_c/T_F$ (with an appropriately defined $T_F$), $T_c/\rho_s(0)$, and $T_c/\omega_0$ are unbounded.


I. Introduction
While superconducting transition temperatures are non-universal properties, and hence not generally amenable to a simple theoretical analysis, understanding what physics determines T c is of self-evident importance. One approach to this challenge is to focus on a key physical process that contributes to the development of superconducting order, and to formulate an upper bound -either rigorous or heuristic -on T c [1][2][3][4][5][6][7].
In this paper, we examine three proposed bounds on T c that are expressed as a fraction of a measurable physical quantity of a given system: an appropriately defined Fermi temperature, a characteristic phonon frequency, or the zero-temperature superfluid phase stiffness. While these putative bounds are physically motivated, and provide valuable intuition in many cases of practical importance, we show by explicit counter-examples that they can be violated by an arbitrary amount. In addition to the fundamental importance of these results, we hope they suggest routes to further optimize T c .

II. Summary of results
We briefly summarize our key results here: 1. The notion of an upper bound on T c in terms of an appropriately defined Fermi energy comes from the fact that, in many situations, as E F → 0, the electrons have no kinetic energy. Thus, in this limit, the superfluid stiffness must seemingly go to zero. What sets T c in the limit of small E F is pertinent to moiré superconductors [8][9][10][11][12][13], where the bands can be tuned to be narrow. To make this question precise, we must define E F in a strongly interacting system. We propose two such definitions of E F , in terms of (i) the difference in the chemical potential * debanjanchowdhury@cornell.edu † kivelson@stanford.edu ‡ erez.berg@weizmann.ac.il between a system with a given density of electrons and a system with a vanishing density, or (ii) the energy dispersion of a single electron added to the empty system. We show that there is no general bound on T c /E F by either definition, by studying two explicit models. In the first model, a flat band is separated by an energy gap from a broad band with pairhopping interaction between the two. The second model consists of a pair of perfectly flat bands with an on-site electron-electron attraction. We show explicitly that the first model violates any putative T c /E F bound when using the first definition of E F above, while the second model violates the bound using either definition of E F . Both models have been defined on a two-dimensional lattice for convenience 1 , but generalized versions of the same models in any D > 1 can be easily seen to exhibit qualitatively similar behaviors.
2. In two-dimensional systems where T c is limited by phase fluctuations, an intuitive bound on T c is given in terms of the zero-temperature phase stiffness, ρ s (0). This comes from the relation [14] T c = πρ s (T − c )/2, and the (often physically reasonable) assumption that ρ s (T ) is a decreasing function of T , and hence ρ s (0) ≥ ρ s (T ). We construct an explicit counter-example in a twoband model of bosons (or, equivalently, tightlybound Cooper pairs), where ρ s (0) can be made arbitrarily smaller than T c . In this model, ρ s (0) can even vanish while T c > 0, implying that there is a reentrant transition into the non-superconducting state below T c .
3. In electron-phonon superconductors, a heuristic bound on T c /ω 0 (ω 0 being the characteristic phonon frequency) was proposed [6,7,15]. The reasoning behind this bound is that, as the dimensionless electron-phonon coupling constant λ increases past an O(1) value, the system tends to become unstable, either towards the formation of localized bipolarons or towards a charge density wave state. At the same time (and relatedly), the phonon frequency is renormalized downward as λ increases, suppressing T c .
Here, we construct an explicit d−dimensional model where these strong-coupling instabilities are avoided, and T c increases without bound upon increasing λ. The model includes N electronic bands interacting with N 2 phonon modes. The model is solvable asymptotically in the large-N limit; then, the famous Allen-Dynes result [16] T c ∝ ω 0 √ λ is valid for large λ, so long as λ N , and hence T c /ω 0 is unbounded as N → ∞. Note that at the heuristic level, it is difficult to identify physical circumstances where more than a few phonons are comparably strongly coupled to the relevant electrons. Nevertheless, our analysis suggests that generically, the larger the number of phonon modes coupled to the electrons, the larger the λ at which the suppression of T c onsets.

III. Flat band superconductivity:
bound on Tc/TF ?
In most conventional superconductors, T c is determined by the energy scale associated with electron pairing. On the other hand, across numerous unconventional superconductors, T c is more strongly sensitive to the "phase ordering scale" [4]. In this context, an important recent advance is the result by Hazra, Verma, and Randeria (HVR) [5] of a rigorous upper bound on ρ s (T ), the temperature-dependent superfluid phase stiffness, in terms of the integral of the optical conductivity over frequency (the optical sum rule). However, since this integral includes all the bands, this upper bound is often of the order of several electron-volts in electronic systems of interest.
At the heuristic level, this bound has been interpreted [17] as implying a bound on T c /T F , where T F = E F /k B is the Fermi energy in units of temperature. 2 At the outset, it is important to define sharply what we mean by E F . A particular protocol that is often adopted in experiments to estimate E F is to use an effective mass, m * , obtained from quantum oscillations along with an estimate of the Fermi momentum (k F ) from a measurement of the the carrier density, and then defining E F = k 2 F /(2m * ). This procedure is only possible when there is a nearby Fermi liquid-like state that displays quantum oscillations. Below, we propose two different definitions of E F , that we can use in settings that do not rely on any underlying assumptions (e.g., that there is a nearby Fermi liquid) and are also amenable to an experimental interpretation. We consider the case in which we add a given density, n, of electrons to an insulating reference state. We can define E F as follows: • Starting from our reference state, we set the temperature to zero and consider the change in the chemical potential, µ(n, T = 0), as we fill in n electrons, as an effective definition of E F . Here µ 0 = lim n→0 µ(n, T = 0). Note that the above definition of E F includes all many-body corrections, which can furthermore be dependent on the density itself, and does not make any reference to any noninteracting limit.
• Alternatively, we can define the Fermi energy through the "non-interacting" density of states, ρ(ε), for adding a single electron to the insulating system. In this case, the Fermi energy (E where ε min is the energy of the ground state with one electron added on top to the insulating reference state. E F is accessible directly in e.g. STM measurements [18].
Below, we provide model Hamiltonians of interacting electrons in flat bands where the superconducting T c exceeds the Fermi energy by one or both of the above def-initions. Thus, these models exemplify flat band superconductivity, where T c is determined entirely by the interaction scale [19][20][21]. Analogous phenomena may also occur in semi-metals [22,23].

A. Flat band superconductivity induced by a nearby dispersive band
We consider a model consisting of a nearly-flat band and a dispersive band. 3 The single-particle part of the Hamiltonian is given by where c † k, ,σ creates an electron with quasi-momentum k in band = 1, 2 and spin polarization σ. We consider the lower band ( = 1) to be a flat band with bandwidth, W 1 , that we will ultimately take to be parametrically small (i.e. W 1 → 0). The upper band ( = 2) is separated from the flat band by an energy gap, ∆ gap , and has a large bandwidth, W 2 W 1 . The bands are topologically trivial and the Wannier functions are tightly localized on the lattice sites.
We now introduce an on-site interaction which scatters a pair of electrons between the flat band and the dispersive band: where R labels a lattice site. Let us focus on the case where the flat band is half-filled, such that the number of particles per unit cell is n = 1. We consider the case where V ∆ gap W 2 . Within mean-field theory, the superconducting transition temperature is given by (see Appendix A 1), where we have assumed a constant density of states per unit cell, ν 2 = 1/W 2 , in the dispersive band. The zerotemperature phase stiffness is given by (Appendix A 1) Hence ρ s (0) T MF , which implies that phase fluctuations are unimportant in determining T c [4], i.e. T c ≈ T MF .
We now examine the Fermi energies E (1,2) F defined in Eqs. (1) and (2) and compare them to T c . Adding a single particle to the empty system, we find that E can be made arbitrarily large. E (1) F is computed in Appendix A 1 by calculating the chemical potentials at n = 1 and n → 0. The result is E F is unbounded is presented in the next section.

B. Flat band superconductivity induced by spatial extent of Wannier functions
We now introduce a different model for superconductivity in a narrow band. The model is defined on a twodimensional square lattice with two electronic orbitals per unit cell. The Hamiltonian is given by Here, where the operator c † k, ,σ creates an electron with momentum k in orbital = 1, 2 with spin σ =↑, ↓. The Pauli matrices τ x,y,z and σ x,y,z act on the orbital and spin degrees of freedom, respectively, and δn r, ≡ σ c † r, ,σ c r, ,σ − 1 is the number of particles at site r and orbital , relative to quarter filling (n = 1). r, r denotes nearest-neighbor sites. The single-particle Hamiltonian H 0 leads to perfectly flat bands with energies ε = ±t. The function α k ≡ ζ(cos k x + cos k y ) controls the spatial extent of the Wannier functions in each band, tuned by the dimensionless parameter ζ. 4 Note that there is no obstruction towards constructing exponentially localized Wannier functions for the above model since the bands are topologically trivial. 5 The strength of on-site attractive Hubbard interaction is denoted U > 0, whereas V > 0 is a nearest-neighbor repulsion.
We are interested in the limit where T U, V ∆ gap (= 2t) and ζ 1. 6 In this limit, we project H U to the lower eigenband. The projected Hamiltonian is expanded in powers of ζ. The average density is set to n = 1 particles per unit cell.
For ζ = 0, the problem effectively reduces to a set of decoupled sites with a strong attractive interaction; the resulting ground state manifold is highly degenerate with local "Cooper pairs" but a vanishing phase stiffness. The linear corrections in ζ vanish due to a chiral symmetry and the orbital l independent interaction strength U . At second order in ζ, the projected interaction, H U , contains nearest-neighbor pair-hopping and density-density interactions: Here, we have also introduced the pseudospin operators ,c † r,σ creates an electron with spin σ in a Wannier-orbital localized around site r in the lower band of H 0 (with τ y = −σ for ζ = 0), and η a are Pauli matrices. The J ⊥ and J z terms correspond to hopping and a nearest-neighbor interaction of the Cooper pairs, respectively, and their strengths are For V = 0, the system has an emergent SU (2) symmetry that relates the pairing and charge order parameters [25]. This symmetry is weakly broken by terms of order (U/t) 2 , not included in Eq. (8). For 0 < V U , the problem is equivalent to a two-dimensional pseudospin ferromagnet with a weak easy-plane anisotropy. Parameterizing the anisotropy by ∆J = J ⊥ −J z , we can estimate the critical temperature of the BKT transition as [26] Note that in the limit ∆J → 0 we get T c → 0, as required by the Mermin-Wagner-Hohenberg theorem.
We now turn to estimating E (1,2) F . Due to the particlehole symmetry of the effective Hamiltonian in Eq. (8), the chemical potential at n = 1 (i.e., η z r = 0) is µ(n = 1) = 0. In the limit n → 0, the system consists of dilute Cooper pairs. In this limit, the interactions between the Cooper pairs can be neglected, and the chemical potential is equal to half the energy per Cooper pair: F . Clearly, since the lower band is completely dispersionless, E We conclude that T c can be made arbitrarily larger than the Fermi energy by either of the two definitions of Eqs. (1) and (2).
It is worth noting that, in the parameter regime we are considering, ρ s (0) ∼ J ⊥ ∼ ζ 2 U . Hence, the delocalization of the Cooper pairs is entirely due to the interactions and the spatial overlap between the Wannier function of the active band, as in Refs. [27][28][29][30][31][32]. 8 Here, however, for U V , we get ∆J IV. Bound on Tc/ρs(0)?
In this section, we turn to the question of whether the zero temperature phase-stiffness, ρ s (0), sets an upper bound on T c in two dimensions. ρ s (0) can be extracted directly from a measurement of the London penetration depth (λ 2 L (0) ∝ 1/ρ s (0)), or from the imaginary part of the low-frequency optical conductivity.
As is well known, in two spatial dimensions, the phase stiffness right below T c is related to T c by the inequality 9 ρ s (T → T − c ) ≥ 2T c /π. However, T c is not directly related to ρ s (0). On physical grounds, it often makes sense to identify ρ s (0) with a "phase ordering scale" that sets an upper limit on T c [4]. This is justified by the fact that ρ s (T ) is usually a monotonically decreasing function of temperature, i.e. ρ s (0) ≥ ρ s (T c ), and therefore T c can be bounded from above by ρ s (0). In conventional superconductors, ρ s (0) T c , and T c is almost entirely determined by the pairing scale. In contrast, in underdoped cuprates, ρ s (0) is close to T c , as illustrated by the famous Uemura plot [34]. This suggests that in these systems, phase fluctuations play an important role in limiting T c [4].
While this line of reasoning is likely correct in most circumstances, we will show here that -as a matter of principle -there is no bound on T c /ρ s (0). We outline a concrete model where ρ s (0) can be made arbitrarily smaller than ρ s (T c ).
Let us begin with a two-dimensional lattice model of two species of (complex) bosons, b 1 , b 2 , where ε 2 (k) is assumed to have a large bandwidth, W 2 , and ε 1 = ε 2 (0) − ε 0 forms a completely flat band at an 9 At a continuous BKT transition, ρs(T → T − c ) = 2Tc/π. However, if the transition is first order, ρs right below Tc can be larger than the universal BKT value. See, e.g., Ref. [33]. energy ε 0 below the bottom of the ε 2 band, i.e. the b 1 bosons are completely localized on individual sites. The dispersion of the two species of bosons is illustrated in Fig. 2. For the purpose of our discussion here, we can approximate ε 2 (k) ≈ k 2 /2m b near the bottom of the broad band. H b int includes an on-site (Hubbard) interaction of strength U 1,2 for the b 1,2 bosons. We take U 2 W 2 whereas U 1 → ∞, implying that the number of b 1 bosons on each localized site can only be 0 or 1. The total average number of bosons per unit cell is chosen to be n b > 1.
At temperatures near T c , the chemical potential is slightly above the bottom of the broad band. Then, assuming that ε 0 T , the average occupation of the localized sites is close to 1/2 (since the b 1 bosons are essentially hard-core bosons at effectively 'infinite' temperature), so there are approximately n b −1/2 bosons per unit cell left to occupy the broad band. The critical tem- [35][36][37]. 10 The momentum distribution of particles at T T c is shown schematically in Fig. 2a.
On the other hand, at T = 0, all the localized sites are filled with one boson. The density of bosons in the broad band is thus n b − 1, and the superfluid stiffness is The boson distribution function is illustrated in Fig. 2b. Clearly, the ratio T c /ρ s (0) can be made arbitrary large by letting n b → 1 + . If 1/2 < n b < 1, the ground state is not a superfluid, and there is a reentrant transition into a superconducting state with increasing T .
Note that in our simple model (Eq. (10a)), the numbers of the two boson species are separately conserved. However, we do not expect the key results to be changed qualitatively by the addition of a weak hybridization between the two species, that breaks this separate conservation of the two boson numbers. In particular, a small hybridization generically produces a perturbative correction to T c and ρ s (0). Indeed, a mild version of this sort of breakdown of the heuristic bound on T c /ρ s (0) has been documented experimentally in Zn-doped cuprates [38]. Here, the pristine material comes close to saturating the heuristic bound; light Zn doping suppresses T c but apparently suppresses ρ s (0) more rapidly, leading to a ratio that slightly exceeds the value proposed in Ref. [4]. This was explained -likely correctly -by the authors of Ref. [38] as being due to Zn-induced inhomogeneity of the superfluid response. This explanation is spiritually close to the model discussed above: Each Zn impurity destroys the superconductor in its vicinity, possibly due to pinning of local spin-density-wave order [39]. In some sense, this can be thought of as a state with localized d-wave pairs near the impurities. This effect depletes the condensate at low temperature, causing a decrease in the superfluid 10 Since we are in two spatial dimensions, in the absence of interactions, Tc = 0.
density. However, near T c , this effect weakens, as the spin-density wave order partially melts. From this perspective, it would be interesting to explore whether this violation can be made parametrically large with increasing Zn concentration -approaching the point at which T c → 0.
Recently, it has been proposed that T c in an electronphonon superconductor can never exceed a certain fraction of the characteristic phonon frequency, ω 0 [6,7]. This putative bound implies that Migdal-Eliashberg (ME) theory [40][41][42] must break down when the dimensionless electron-phonon coupling λ is of order unity [6,15], since according to ME theory, T c grows without limit with increasing λ [16]. In general, the failure of ME theory at λ = O(1) is a result of strong-coupling effects: (i) The lattice tends to become unstable for large λ, resulting in a charge density wave (CDW) transition; (ii) electrons become tightly bound into bipolarons, whose kinetic energy is strongly quenched in the strong coupling limit; and (iii) as λ increases, the phonon frequency renormalizes downward by an appreciable amount, ∆ω, suppressing T c [43,44].
These strong-coupling effects certainly play an important role in limiting T c in real systems, where it is typically found never to exceed about 0.1 ω 0 across numerous conventional superconductors [7]. Determinant Monte Carlo simulations of the paradigmatic Holstein model reveal that ME theory indeed fails for λ = O(1), and the maximal T c is significantly below ω 0 [6]. As we shall now show, however, this is not a rigorous bound on T c .
To demonstrate this, we consider a particular large−N variant of the electron-phonon problem [45,46] with N −component electrons and N × N −component ('matrix') phonons, defined on a d−dimensional hypercubic lattice. The Hamiltonian is given by Here, c † r,σ,a creates an electron at position r with spin σ in "orbital" a. The hopping parameters t rr and chemical potential µ are assumed to be identical for all orbitals. We have introduced a real, symmetric matrix of phonon displacements, X r and their canonically conjugate momenta, P r , with frequency ω 0 = K/M , assumed to be much smaller than the Fermi energy. The phonons are taken to be dispersionless for simplicity. The purely on-site electron-phonon coupling constant is denoted α, with a N −dependent normalization factor that ensures that the model has a finite energy density in the N → ∞ limit. The dimensionless electron-phonon coupling constant is defined as λ = α 2 ν(0)/K, where ν(0) is the electronic density of states at the Fermi level per orbital.
We are interested in the large−N limit of the model defined in Eq. (11a). Since the number of phonon degrees of freedom is much larger than the number of electron orbitals, the phonon dynamics are essentially unaffected by the coupling to the electrons, even when the electrons are strongly perturbed. This implies that the strongcoupling effects mentioned above are suppressed, even for λ 1. In particular, as we show in Appendix C 1, there is no lattice instability or polaron formation as long as λ N , and the softening of the phonon frequency is only of the order of ∆ω ∼ λω 0 /N .
To zeroth order in 1/N , the equations for the electron self-energy and the pairing vertex are exactly those given by Eliashberg theory, whereas the phonon self-energy is of order 1/N (see Fig. 3). Thus, the self-consistent equations for the pairing vertex are identical to those of ME theory neglecting the renormalization of the phonons, and hence the result is the same. In particular, for 1 λ N , T c ≈ 0.1827 ω 0 √ λ [16]. Implicit in the fact that Migdal-Eliashberg theory is exact at N → ∞ is an assurance that there is no suppression of T c by phase fluctuations. 11 More specifically, this follows from the 11 In the d = 2 version of our model, the BKT temperature dif-observation that the superfluid stiffness is O (N ). Thus, T c is unbounded.
The key ingredient in our model that allows us to take λ > 1 without suffering from lattice instabilities is that the different phonon modes couple to electron bilinears that do not commute with each other [see Eq. (11d)]. This limits the energy gain from distorting a given set of phonon modes when forming a CDW or a polaron bound state, since the resulting perturbations to the electronic Hamiltonian cannot be diagonalized simultaneously. In contrast, the contributions of the individual phonon modes add algebraically in the total dimensionless coupling λ that enters the equation for the pairing vertex (the same dimensionless coupling also determines the resistivity in the normal state of this model [45]).
It is worth noting that while these considerations may be of some use in searching for systems with ever higher T c , as a practical matter it may be difficult to significantly violate the proposed heuristic bound. To achieve T c ≈ ω 0 requires the extremely large value of λ ≈ 25. At the same time, to avoid polaron formation requires N λ, which means the number of distinctly coupled phonon modes would have to be N 2 λ 2 ∼ 625!

VI. Outlook
The notion of a fundamental upper bound on T c for models of interacting electrons is an attractive concept. In this paper, we have demonstrated that while there are numerous physical settings where such bounds can be formulated at a heuristic level, there exists no fundamental, universal upper bound on T c in terms of the characteristic energy scales of interest to us, which include an appropriately defined T F , ρ s (0) and ω 0 . We have constructed explicit counter-examples which violate these heuristic bounds by an arbitrary amount.
On the experimental front, it would be fruitful to look for candidate materials where the heuristic bounds are violated by a large amount. The fact that these bounds are usually satisfied is to be expected, since although the bounds are not rigorous, the physical reasoning behind them is quite robust. As our theoretical discussion illustrates, whenever such bounds are violated, there is an interesting underlying physical reason behind the violation; moreover, the mechanisms behind the violation of the heuristic bounds may suggest ways to optimize T c . Our work provides two such examples. Flat band systems with a large spatial extent of the Wannier functions are a promising platform for increasing T c . In electron-phonon systems, the instabilities that limit T c at large electronphonon interaction strength can be partially mitigated if the coupling is shared between several active phonon fers from the mean-field transition temperature only by a 1/N correction.
modes that couple to non-commuting electronic operators.
A. Additional details for flat band superconductivity induced by nearby dispersive band

Mean-field gap equations
Here, we derive the expressions for T MF , ρ s (0), and E The mean-field equations for the pairing potentials ∆ =1,2 are given by Here, we have set the width of the first band to zero, such that whereas, for simplicity, we have taken the dispersion of the second band to be parabolic with an effective mass m: The total width of this band is W 2 . The density of particles per unit cell is where E (k) = ξ 2 (k) + ∆ 2 . We set the density to n = 1. We first solve for the critical temperature T MF . We assume that µ ≈ 0, to be checked self-consistently. We then linearize Eqs. (A2,A3) in ∆ 1,2 at β = 1/T MF . Eq. (A3) gives ∆ 2 = V ∆ 1 /(4T MF ); inserting this in Eq. (A2) gives T MF = V 2 log (W 2 /∆ gap ) /(8W 2 ). Under our assumptions, we get that T MF ∆ gap , and therefore the number of particles in the second band is exponentially small (setting ∆ 2 = 0 at T = T MF ). Hence, Eq. (A6) indeed gives that µ ≈ 0 corresponds to a filling of n = 1, consistently with our starting assumption.
Next, we compute E F . The chemical potential at n → 0 can be obtained by solving the linearized gap equations at T = 0: which give the chemical potential µ at which superconductivity onsets. This gives µ(n → 0) = −2T MF . To find the chemical potential at n = 1 and T = 0, we solve the gap equations at T = 0. Assuming that µ(n = 1) ≈ 0, we obtain ∆ 2 = V /2, Inserting this back into Eq. (A6) gives µ(n = 1) ∝ , justifying our assumption that µ is small. Thus we find that E (1) F = µ(n = 1) − µ(n → 0) = 2T MF , the result quoted in the text.

Phase stiffness
Let us compute the phase stiffness at T = 0 within mean-field theory. We couple the system to an external gauge field, and differentiate the ground state energy twice with respect to vector potential. At T = 0, the entire contribution to the stiffness arises from the diamagnetic term. Using ∆ 2 = V /2 ∆ gap , we can further expand the result to leading order in |∆ 2 | such that, Note that, within mean-field theory, we can write the average particle number in the = 2 band as which can be similarly expanded to leading order in a small |∆ 2 | as, Thus, in the particular flat-band limit for = 1, we find that the phase stiffness is determined by the usual Galilean invariant result, ρ s (0) = n s /m, where n s now corresponds to the average particle density in the dispersive = 2 band whose mass is m.
As mentioned in the main text, we are interested in the limit where T U, V ∆ gap (= 2t) and ζ 1. For V = 0 and ζ = 0, we recover the Hubbard interaction with a renormalized coupling strength U l |ψ ri (l, r i )| 4 = U/2. Hence, each site r i is either empty or double occupied at low temperatures, T U .
For finite ζ, one might expect correlated hopping terms at linear order. However, those terms vanish identically. Due to the absence of τ z in Eq. (7) the system has a chiral symmetry, which guaranties that the wavefunctions of the low-energy band have equal weights on both orbitals. This, combined with the orbital-independent coupling strength U , implies that all terms linear in ζ originating from orbital l = 1 are canceled by l = 2 contributions. Note that such a hopping term, e.g., d † r,σ d r,σ d † r,σ d r+ x,σ , necessarily generates two single-occupied sites and hence creates an excitation with an energy of order U . Accordingly, terms linear in ζ cannot show up in the projected Hubbard interaction as leading order perturbation even if the chiral symmetry is weakly broken or if the coupling strength U was weakly orbital dependent.
Hence, the leading order corrections are quadratic in ζ. Nearest neighbor density interaction and pair hopping terms preserve the local parity and thus are relevant perturbations within the ζ = 0 ground-state sector. Their effective coupling strength are given by l |ψ ri (l, r i )| 2 |ψ ri (l, r j )| 2 and l ψ * ri (l, r i ) 2 ψ ri (l, r j ) 2 , respectively, up to combinatorial factors.
C. Additional details of the electron-phonon model

Stability in the large-N limit
We study the stability of the system towards CDW order and polaron formation at T = 0 in the limit M → ∞ keeping K fixed, i.e., ω 0 = K/M → 0. We assume that setting M finite or T > 0 will only make the system more stable, as these promote quantum and thermal fluctuations of the phonons, respectively.
To second order in the phonon displacements, the ground state energy of the system can be computed by perturbation theory. The result is where χ q is the static charge susceptibility: where L is the linear dimension of the system. The condition for stability with respect to CDW formation is that the term in parenthesis in Eq. (C1) be non-negative for all q. Note the 1/N prefactor of the second term in Eqn. C1. Hence, in the large-N limit, the system does not possess a CDW instability. By the same token, polaron formation is also suppressed, since a local deformation of the phonons accompanied by an electronic bound state is never favored. These instabilities only appear when λ becomes of the order of N .

Eliashberg equation
The Eliashberg equations for electron-phonon superconductivity and their solution at large λ are quite standard [15,16,42]. For completeness, we will nevertheless outline the derivation below.
In the large-N limit, the electron self-energy at T > T c obeys the following self-consistent equation, depicted in Fig. 3b: where D 0 (iΩ n ) is the bare phonon Green's function The phonon self-energy, shown in Fig. 3a, is of order 1/N in our model, and is thus neglected. In the second line of Eq. (C3), we have performed the k summation, assuming that at the frequencies of interest, W |ω n − Σ(iω n )|, where W is the electronic bandwidth. At temperatures large compared to ω 0 , the largest contribution to the self-energy comes from ω n = ω n , giving Σ(iω n ) ≈ −iπλT sgn(ω n ).
The linearized Eliashberg equation for the pairing vertex Φ(iω n ), represented diagrammatically in Fig. 3c, is given by It is convenient to define the Eliashberg Z factor: (C7) In the limit of large λ, T c is much larger than ω 0 . Importantly, the ω n = ω n term in the sum above contribution vanishes, which allows us to take the limit ω 0 /T → 0 in the deminator of the right-hand side. This results in the following equation: