Pairing mechanism in multiband superconductors

It has been a long-standing puzzle why electrons with repulsive interactions can form pairs in unconventional superconductors. Here we develop an analytic solution for renormalization group analysis in multiband superconductors, which agrees with the numerical results exceedingly well. The analytic solution allows us to construct soluble effective theory and answers the pairing puzzle: electrons form pairs resonating between different bands to compensate the energy penalty for bring them together, just like the resonating chemical bonds in benzene. The analytic solutions allow us to explain the peculiar features of critical temperatures, spin uctuations in unconventional superconductors and can be generalized to cuprates where the notion of multibands is replaced by multipatches in momentum space.


Wen-Min Huang 1✉ & Hsiu-Hau Lin 2✉
it has been a long-standing puzzle why electrons with repulsive interactions can form pairs in unconventional superconductors. Here we develop an analytic solution for renormalization group analysis in multiband superconductors, which agrees with the numerical results exceedingly well. The analytic solution allows us to construct soluble effective theory and answers the pairing puzzle: electrons form pairs resonating between different bands to compensate the energy penalty for bring them together, just like the resonating chemical bonds in benzene. the analytic solutions allow us to explain the peculiar features of critical temperatures, spin uctuations in unconventional superconductors and can be generalized to cuprates where the notion of multibands is replaced by multipatches in momentum space.
It has been over one hundred years since Heike Kamerlingh Onnes discovered the resistance of mercury suddenly drops to zero 1 when cooled down by liquid helium in 1911 and marked the birth of superconductivity. The exotic phenomena of superconductors remained mysterious until Bardeen, Cooper and Schrieffer (BCS) 2 came up with a complete theory in 1957. There are two key ingredients in BCS theory. Firstly, a negatively charged electron distorts the nearby lattice formed by positively charged ions and another electron thus feels an effective attraction toward the original electron. In professional jargons, one can say that the effective attraction is mediated by electron-phonon interactions. The second key ingredient is pair formation: an electron with specific momentum and spin pairs up with another electron with opposite momentum and spin via the effective attraction. These electron pairs, referred as Cooper pairs, form a Bose-Einstein condensate in the superconducting state 3,4 . Despite the celebrating success of BCS theory, there are other superconductors including cuprates 5-7 , heavy-fermion compounds [8][9][10] , organic superconductors 11 , iron-based superconductors 12,13 and recently found superconducting graphene superlattices 14 , remaining queer 15,16 and unexplained by the electron-phonon interactions alone. It is generally believed that the Coulomb interaction between electrons may be responsible for the emergent superconductivity 4 . However, there is a huge mismatch of the energy scales -the critical temperature is much smaller than the bare electronic interactions. Perhaps the most intriguing puzzle for these unconventional superconductors is the pairing mechanism: what is the glue to pair up electrons from mutual repulsive interactions? The experimental evidences seem to suggest that pairing in the unconventional superconductors is not due to electron-phonon interactions 13 . Due to strong magnetic correlations [17][18][19][20] in these materials, it is proposed that spin fluctuations [21][22][23][24][25][26][27] may play the role of glue to pair electrons up. The renormalization-group (RG) studies 3,[26][27][28][29] seems to indicate that the antiferromagnetic spin fluctuations may be the cause of unconventional superconductivity in iron-based materials.
In this Article, we investigate the pairing mechanism in multiband iron-based superconductors by RG analysis. As in previous RG studies 3,4,26,[30][31][32] , we integrate out quantum fluctuations at shorter length scales within the one-loop approximation and seek for the ground state in the low-energy limit. Because the RG equations are intrinsically non-linear, numerical analysis is already challenging, rendering simple understanding beneath the messy numerics inaccessible in most cases. We adapt the classification scheme of relevant couplings developed in previous works [33][34][35] to overcome the challenge and analyze how superconductivity emerges in RG transformations.
In addition, we make use of the scaling relation to bridge the results obtained by RG analysis and those in the conventional mean-field theory. While it is well known that the mean-field theory provides a self-consistent description for various physical quantities, it is less familiar to the scientific community that the mean-field theory also lays out the RG equations for its parameters. For example, for the generalized multiband BCS Hamiltonian, the vanilla version of the mean-field theory is characterized by two important parameters: g and g ⊥ , standing for intraband and interband pair hopping. The functional form of the gap function computed within the mean-field

Results
To put our feet on the firm ground, we choose iron-based superconductor as a demonstrating example. In the following, we elaborate the numerical and analytical details how the above conclusions are achieved. It is worth emphasizing that most conclusions are not limited to iron-based superconductors. Generalization of our theory will be discussed at later paragraphs.
RG for Multiband Superconductor. The kinetic energy in iron-based superconductors is described by a five-orbital tight-binding model with appropriate hopping matrix elements 45,46 . The mutual repulsion between electrons is model by the simplest on-site interactions, containing intra-orbital repulsion U 1 , inter-orbital repulsion U 2 and the Hund's coupling J H between different orbitals. Detail descriptions of the model can be found in Methods.
The above Hamiltonian has been studied by functional renormalization group (fRG) method, which shows that the ground state is a superconductor with sign-reversed s ± pairing symmetry 30,31 . Due to the s-wave symmetry, it is enough to sample one point on each Fermi surfaces for qualitative understanding. It turns out that this simplification is equivalent to put the iron-based superconductor on the four-leg ladder with periodic boundary conditions along the x direction. It will become clear later that this geometry does not change the ground state but simplify the algebra significantly. Later, we will come back to two-dimensional calculations and reaffirm the validity of our major claims.
In weak coupling, the relevant degrees of freedom can be labeled by the quantized momentum k x and the chirality R/L as shown in Fig. 1. Following the standard field-theory description [47][48][49] , the allowed interactions are Cooper scattering σ ρ c c , ij ij and forward scattering σ ρ f f , ij ij between different bands, where σ, ρ denote the spin and charge channels respectively. The RG equations for all couplings are given explicitly in Methods. Though the RG equations can be written down explicitly, the solutions are still too complicated. With l = In(Λ 0 /Λ) being the logarithm of the ratio between bare energy cutoff Λ 0 and the running cutoff Λ, we integrate the coupled differential equations numerically, and elaborate the details in Methods. We find all couplings are captured by the scaling Ansatz [33][34][35] , between the hole pocket centered at (0,0) (band 1) and the electron pocket at π ± ( , 0) (band 2). The positive sign of the interband Cooper scattering c 12 leads to the s ± pairing symmetry. Note that our RG Ansatz predicts the dominant superconducting bands with the correct pairing symmetry as obtained by the fRG method. In fact, a detail analysis including all subdominant relevant couplings lead to correct signs for all gap functions in different bands as those obtained in the fRG method. Thus, the ground state in the ladder geometry is the same as that in the two-dimensional system.
However, due to the explicit form of the RG equations for ladder geometry, the pairing mechanism in the multiband superconductor can be studied analytically. Construct the couplings, ≡ + + + can be extracted from the numerical RG flows as the following. First of all, since the RG flows diverge at l = l d , it fixes that , ensuring the coupling strengths are of the same order. Because c(0), c ⊥ (l) are given and l d is found in numerics already, the values of g(0) and g ⊥ (0) are pretty much determined.
The above simplification is inspiring: it implies that the RG flows for the intraband and the interband pair hopping c(l) and c ⊥ (l) from integrating the whole set of nonlinear differential equations (about 50 equations for this case) can be captured by the analytic solutions given in Eq. (3). Comparison between the numerical solutions c(l), c ⊥ (l) and the analytic ones g(l), g ⊥ (l) are shown in Fig. 2, agreeing with each other exceedingly well in the whole RG regime. Thus, if we can understand the physics for the analytic solutions, we can answer the long-standing puzzle of pairing mechanism.
Look for effective theory. To dig out the secret message behind the analytic solutions, let us switch gear to review the well-known single-band BCS Hamiltonian first. The gap function takes the form Δ = Λ g g ( ) 2 exp(1/ ), where g < 0 arisen from effective attraction and Λ is the Debye energy. But, it is not widely known that the RG equation for g can be derived from the dependence of the gap function. In the RG method, the gap function is rescaled with rescaling energy of systems in RG steps. By integrating high-energy modes, the rescaling of the gap function is described by the RG equation dΔ/dl = Δ. It is said that gap functions at different energy scales can be related by , which is the so-called scaling argument in RG transformations 50 . Taking the derivative, d/dl, on the both sides, the RG equation reads, dg/d = −g 2 l. Therefore, if the RG flow for some system is described by the derived RG equation, its effective theory is just the single-band BCS Hamiltonian.
We can generalize the scaling argument to multiband BCS Hamiltonian as elaborated in Methods. In particular, when two bands reign interband pair hopping, the corresponding RG equations are where g and g ⊥ stand for intraband and interband pair hopping. The above non-linear coupled equations can be solved exactly, giving the analytic solutions we discussed previously. Therefore, the effective theory for iron-based superconductors is the multiband BCS Hamiltonian, proven by matching the RG flows together. RG flows for repulsive interactions are shown in Fig. 3. If the intraband pair hopping is larger than the interband one (g = g ⊥ < 0), the ground state flows toward the Fermi-liquid fixed point and no superconductivity occurs. However, if the interband pair hopping is larger (g = g ⊥ < 0), the ground state flows toward the superconducting phase with s ± pairing symmetry. It is worth mentioning that, if the initial couplings are close to the symmetric ray g = g ⊥ , it will flow toward the Fermi-liquid fixed point first and then turns around to the unconventional superconducting state. This means that, upon cooling down the system toward the critical temperature, it exhibits non-trivial crossover properties over a wide range of temperatures. This crossover is shown in the orange-shaded regime in Fig. 3(b) where both couplings g and g ⊥ are small.
It is interesting to compare with the pairing instability in the absence of interband pair hopping (g ⊥ = 0). As is clear in Fig. 3(a), a negative g(0), i.e. attractive interaction, is required to trigger the superconducting instability. www.nature.com/scientificreports www.nature.com/scientificreports/ For repulsive interaction g(0) > 0, the RG flow always bring the system to the trivial Fermi-liquid fixed point. The RG flows in the special case g ⊥ = 0 give rise to the commonly accepted criterion that an effective attraction is necessary for Cooper pair formation. However, the criterion is obviously wrong as the RG flows in the upper plane is quite different from those in the horizontal axis.
The key to Cooper pair formation in the presence of interband pair hopping is characterize by the parameter in contrast to the intraband pair hopping g in conventional superconductors. The parameter δ offers a vivid picture for electron pairing. The penalty to form a Cooper pair within the same band is g. However, through interband pair hopping, a Cooper pair gains − ⊥ g benefit through resonating between different bands. Unlike the conventional superconductors, the pairing mechanism in multiband superconductors is resonating pair hopping between different bands, just like the resonating chemical bonds in benzene. Our findings are consistent with the previous discoveries of superconductivity in the interband-dominated regime [38][39][40][41][42][43][44] .
This picture leads to several significant changes of our common intuitions. First of all, the attractive interaction is no longer a necessity. As long as δ < 0, the generalized BCS Hamiltonian delivers a superconducting ground state. The second feature is that the critical temperature is not locked up the the energy scale of the bare interactions g and g ⊥ anymore. Following the textbook calculations, the critical temperature is where Λ is of the same order of electronic band width. Setting Λ = = t eV 1 , reasonable estimate for the hopping amplitude t, we extract δ from the numerical RG flows for different on-site interaction strengths and plot the critical temperatures in Fig. 4. For the particular interaction profile we choose, , the parameter δ is one order of magnitude smaller that the bare coupling U t /t. This explains why the critical temperature can be so small while the bare interaction strength is large. If one chooses = . H , the predicted critical temperature is about 56 K -quite a reasonable estimate. Finally, the interband pair hopping g ⊥ brings in another interesting property in unconventional superconductors. Making use of the instability analysis developed by Wang and Lee 3 , the interband pair hopping also enhances spin-density-wave (SDW) instability if the momentum Q = K 1 − K 2 , connecting the two Fermi surfaces, is close to half of the reciprocal lattice vectors, i.e. Q = G/2. In the case of iron-based superconductor in ladder geometry, Q = (π, 0) satisfies the above condition. Therefore, it is expected that the antiferromagnetic spin fluctuations at the momentum Q The critical temperature is dictated by the small parameter δ, which can be extracted from numerical RG flows and is one order smaller than the bare interaction strength.  In the absence of interband pair hopping, we have checked that the above RG equations reduce to those derived by Shankar with just one circular Fermi surface. The general RG equations involving integration over internal angles can be simplified to the previous result. Assuming the density of states and the bare couplings are rotationally invariant, i.e. θ θ θ θ 2 , the RG equations can be decoupled by the partial-wave decomposition. By identifying → C g PP and → ⊥ C g PP , we go back to the same set of RG equations in Eq. (3) again. Therefore, our claims can be readily applied to two dimensional systems as well.
The two-pocket model presented here provides another perspective on the relation between superconductivity and antiferromagnetic spin fluctuations. In iron-based superconductors, the momentum connecting two Fermi surfaces is π π = Q ( , 0), (0, ) and thus the interband pair hopping enhances the antiferromagnetic spin fluctuations at these momenta. However, for the two pocket model shown in Fig. 5, both Fermi surfaces locate at the same momentum. As a consequence, though the superconductivity is still driven by the interband pair hopping, it does not come with any enhanced antiferromagnetic spin fluctuations.

Discussion
Our previous calculations concentrate on the dominant bands in the iron-based superconductor to reveal the pairing mechanism in multiband superconductors on firm ground. But, it is reasonable to include all active bands in the generalized BCS Hamiltonian when quantitative accuracy is required. In addition, it is intriguing to check whether the fRG flows can be captured by the analytic solutions in Eq. (2) as well.
For single-band materials but with significant variations in density of states at different momenta, one should divide the Fermi surface into appropriate patches and applies the method developed here. We have carried out primitive RG analysis on cuprates by cutting the single-band Fermi surface into 16 patches (minimum to differentiate nodal and antinodal regimes with four-fold dihedral symmetry). The dominant patches locate at (π, 0) and (0, π) with sign-reversed gap functions, agreeing with the d-wave symmetry. It implies that pairing mechanism in cuprates arises from resonating pair hopping between antinodal regimes 51 . The interpatch pair hopping between antinodal regimes also explains the enhance spin fluctuations at (π, π). The pairing instability in the nodal regimes is only triggered when approaching close to the critical temperature. How to rigorously generalize the analysis presented here from multiband superconductors to multipatch single-band system remains an open question at this point.
We would like to emphasize that the weak-coupling RG analysis demonstrated here may work well in strong coupling if no quantum phase transition is encountered. That is to say, the generalized BCS Hamiltonian with multiple bands (or different momentum patches) may be the mother Hamiltonian for most unconventional superconductors in the low-energy limit. Even the pairing mechanism is revealed, non-trivial crossovers in the unconventional superconductors are expected at wide temperature range and a complete phenomenological description is still lacking. Finally, the pursue for higher critical temperatures should not be on seeking for stronger glues. In stead, it is time to search for optimal interaction profile which gives the largest interband pair hopping g ⊥ and suppresses the intraband one g. It is not yet clear whether the optimal electronic interaction for different materials share a generic profile, or the answer may vary from one material to the other.

Methods
RG equations in ladder geometry. The band structure of the iron-based superconductor can be described by the five-orbital tight-binding model, † 1, 2, , 5 label the five orbitals of iron, represented as , is the spin index. The 5 × 5 kinetic matrix K vw in the momentum space has been constructed in previous studies 45 . The Fermi surface with electron doping = .
To investigate the correlation effects, the most general on-site interactions 30 a a ( ) iv iv iv and i runs over all lattice sites. Note that, due to the multiple orbitals, even the simplest on-site interactions contain three terms: intra-orbital repulsion U 1 , inter-orbital repulsion U 2 and Hund's coupling J H between different orbitals. By putting the iron-based superconductor on the four-leg ladder with periodic boundary conditions along the x direction, we quantize the momentum k x and reduce the five active bands to five pairs of Fermi points, as show as Fig. 1a. For the ladder geometry, the RG equations for all couplings have been studied in previous literatures [33][34][35][47][48][49] , and can be written down explicitly, ij ij ij ij ij ij ij ij 2 2 The notation ≡  g dg dl / with = Λ Λ l ln( / ) 0 being the logarithm of the ratio between bare energy cutoff Λ 0 and the running cutoff Λ. The coefficient tensor depends on the Fermi velocity v i for each band. We intentionally separate the intraband and interband RG equations for clarity. Various symmetries ensure that c ij = c ji and f ij = f ji . To avoid double counting, we choose f ii = 0 here.
To compute the initial conditions of the RG transformation, we rewrite the field operators in the five-orbital tight-binding Hamiltonian by the chiral field decomposition 33,[47][48][49] Here u R/Lia , a = 1, 2, …, 5, are the five orbital components of the Bloch wave functions associated with the R/L (or +y/−y direction) members of the i-th Fermi pair, and ψ L R / are chiral-field operators. After some calculation, the initial conditions of c σ and c σ terms can be represented as  With the similar computation, we also obtain the initial conditions of the forward couplings,  We notice that the initial conditions of g and g ⊥ of our effective theory in the main text are = + + + Within the initial conditions for U 1 = 2.8, U 2 = 1.4 and J H = 0.7, we numerically solve the RG equations to obtain the RG flows with different filling factors. We find that the coupling constants diverge at the same length scale l d . By showing the RG flows versus 1/(l d − l) in the log-log plot, as illustrated in Fig. 6, we find that the hierarchy order of coupling constants can be differentiated by RG exponents in the scaling ansatz 34 . The RG exponents for different filling factors are showed in Fig. 1.
where α, β are spin indices, ω and K = (K x , K y ) are labelled as the frequency and the two-dimensional momentum respectively. The derivation is aided by dividing the Fermi surfaces into two bands and N patches for each band in the vicinity of the two Fermi surfaces, determined by the energy cut-off Λ 0 , as shown in Fig. 7. Since the Fermi momenta are different, each patches are not of the same size and will give rise to a factor related to the curvatures of Fermi surfaces. Under the patching approach, we rewrite the momentum integral as with k FP and v P respectively being the Fermi momentum and Fermi velocity of the band P. We also approximate that the dispersion of the two Fermi pockets are independent of the angular direction. As a consequence, the tree-level RG analysis, demonstrated by Shankar 50 , shows that the only marginal couplings are the forward and Cooper scattering in this con- 1, 2, , being the patching indices. We note that the patching index = + i i N/2 is represented as the time-reversal partner of the i-th patch. With the convention of the band indices = 1 2 = (2 1), the vertex function in the patching approach can be written as The L pp and L ph functions, depicting particle-particle and particle-hole bubble diagrams respectively, are denoted as and ω ω ω ω χ χ  . We note that the summation of patches is constrained by the integral, which is also shown in Shankar's single-pocket model 50 .
With the same spirit, we compute the full set of RG equations. Furthermore, to simplify the coefficients in the RG equations, we rescale the couplings as PQ P Q P Q PQ Thus, the full set of non-trivial RG equations is represented as, PP PP , we choose = F i i ( , ) 0 PP to avoid double counting. Take N → ∞ and the two dimensional RG emerges. Terms with 1/N factor but without summation vanish in the two dimensional limit. This simplification is the same as the elegant "phase-space argument" invented by Shankar 50 before. Terms with summations smoothly evolve back into angular integral,