Stabilizing Even-Parity Chiral Superconductivity in Sr$_2$RuO$_4$

Strontium ruthenate (Sr$_2$RuO$_4$) has long been thought to host a spin-triplet chiral $p$-wave superconducting state. However, the singlet-like response observed in recent spin-susceptibility measurements casts serious doubts on this pairing state. Together with the evidence for broken time-reversal symmetry and a jump in the shear modulus $c_{66}$ at the superconducting transition temperature, the available experiments point towards an even-parity chiral superconductor with $k_z(k_x\pm ik_y)$-like $E_g$ symmetry, which has consistently been dismissed based on the quasi-two-dimensional electronic structure of Sr$_2$RuO$_4$. Here, we show how the orbital degree of freedom can encode the two-component nature of the $E_g$ order parameter, allowing for an $s$-wave orbital-antisymmetric spin-triplet state that can be stabilized by on-site Hund's coupling. We find that this exotic $E_g$ state can be energetically stable once a complete, realistic three-dimensional model is considered, within which momentum-dependent spin-orbit coupling terms are key. This state naturally gives rise to Bogoliubov Fermi surfaces.

Introduction.-Based on early Knight shift [1], polarized-neutron-scattering [2], muon-spin-resonance [3], and polar Kerr measurements [4], Sr 2 RuO 4 has been widely thought to support a spin-triplet chiral p-wave superconducting state with E u symmetry [5][6][7].This proposed state has had difficulty reconciling other experimental results [7], including the absence of chiral edge currents [8], thermal transport consistent with a nodal state [9][10][11], apparent Pauli-limiting effects for in-plane fields [12], and the failure to observe a cusp-like behavior of the critical temperature under nematic strain [13,14].Plausible explanations for each of these inconsistencies have been presented [7,15,16].Recently, however, the Knight shift has been revisited [17,18] and, contrary to earlier results, a relatively large reduction of the Knight shift for in-plane fields in the superconducting state has been observed.This finding cannot be reconciled with the standard spin-triplet chiral p-wave state [6].
Although it now seems unlikely that Sr 2 RuO 4 is a spintriplet chiral p-wave superconductor, the observation of broken time-reversal symmetry [3,4] and a jump in the shear modulus c 66 [19,20] at the critical temperature still indicate a multi-component order parameter [21].The only other possible multi-component channel within D 4h symmetry belongs to the E g irreducible representation (irrep) [21].At the Fermi surface, a chiral order parameter in this channel resembles a spin-singlet d-wave state, which has horizontal line nodes.Such a state would appear to imply that the dominant pairing instability involves electrons in different RuO 2 layers, which is difficult to understand in view of the pronounced quasi-twodimensional nature of the normal state of Sr 2 RuO 4 .In-deed, no microscopic calculation for Sr 2 RuO 4 has found a leading weak-coupling E g instability [22][23][24].
In this Letter, we show that local interactions can lead to a weak-coupling instability in the E g channel, once we consider a complete three-dimensional (3D) model for the normal state.Physically, this E g state is an s-wave orbital-antisymmetric spin-triplet (OAST) state stabilized by on-site Hund's coupling.When the renormalized low-energy Hund's coupling J becomes larger than the inter-orbital Hubbard interaction U , this channel develops an attractive interaction [25][26][27][28][29][30].This pairing instability has been found in dynamical mean-field theory, which predicts it appears in the strong-coupling limit even when the unrenormalized high-energy J is less than U [31], and also in the presence of strong charge fluctuations [32].Pairing due to this type of interaction was considered for Sr 2 RuO 4 in Ref. [25] , where an A 1g pairing state was found to be stable.Motivated by the relevance of J for the normal state of Sr 2 RuO 4 [33], we revisit the local-pairing scenario.In the following, we show that an E g state can be stabilized over the A 1g state of Ref. [25] by including momentum-dependent spin-orbit coupling (SOC) corresponding to interlayer spin-dependent hopping with a hopping integral on the order of 10 meV.This small value leaves the quasi-two-dimensional nature of the band structure intact.Moreover, we use the concept of superconducting fitness [34,35] to understand the importance of this term in stabilizing the E g state.Finally, we show that this chiral multi-orbital E g state will display Bogoliubov Fermi surfaces [36,37], instead of line nodes.
Normal-state Hamiltonian.-Anaccurate description arXiv:1912.09525v1[cond-mat.supr-con]19 Dec 2019 of the normal-state Hamiltonian is crucial for understanding superconductivity in the weak-coupling limit.Our starting point is a tight-binding parametrization of the normal-state Hamiltonian that includes all terms allowed by symmetry [38].To determine the magnitude of each term, we carry out a fit to the densityfunctional-theory (DFT) results of Veenstra et al. [39].Details on the numerical procedures are provided in the Supplemental Material (SM) [40].However, angleresolved photoemission spectroscopy (ARPES) measurements [33,44] suggest that some DFT parameters differ appreciably from the measured values, in particular the SOC strengths [39].We therefore allow the SOC parameters to vary in order to understand how they affect the leading superconducting instability, under the constraint that the Fermi surfaces do not differ significantly from the DFT predictions and are hence still qualitatively in accordance with the ARPES results.
The relevant low-energy degrees of freedom (DOF) are the electrons in the t 2g -orbital manifold d yz , d xz , and d xy of Ru.Using the spinor operator Φ , where c † k,γσ creates an electron with momentum k and spin σ in orbital γ, we construct the most general three-orbital single-particle Hamiltonian as where the λ a are Gell-Mann matrices encoding the orbital DOF and the σ b are Pauli matrices encoding the spin (λ 0 and σ 0 are unit matrices), and h ab (k) are even functions of momentum.Time-reversal and inversion symmetries allow only for fifteen h ab (k) functions to be finite.The explicit form of the h ab (k) functions and the Gell-Mann matrices are given in the SM [40].
Interactions and superconductivity.-Weconsider onsite interactions of the Hubbard-Kanamori type [45], where c † iγσ (c iγσ ) creates (annihilates) an electron at site i in orbital γ with spin σ, and n iγσ = c † iγσ c iγσ .The first two terms describe repulsion (U, U > 0) between electrons in the same and in different orbitals, respectively.The third and fourth terms represent the Hund's exchange interaction and pair-hopping interactions respectively.We take J = J [45], where J > 0 is expected for Sr 2 RuO 4 .In the context of Sr 2 RuO 4 , H int is usually taken as the starting point for the calculation of the spin-and charge-fluctuation propagators which enter into the effective interaction [46,47].Here, we take a different approach [25,28] by directly decoupling the interaction in the Cooper channel, which, for U − J < 0, yields an attractive interaction for on-site pairing in an OAST state.This scenario has previously been applied to a two-dimensional model of Sr 2 RuO 4 , predicting an OAST A 1g state [25].Although a strong-coupling instability towards an OAST E g state in the absence of SOC has been predicted in Ref. [32], the superconductivity in Sr 2 RuO 4 is likely in the weak-coupling regime [7].It is therefore important to understand if an OAST E g state can be the leading instability in this limit.
In the spirit of Ref. [28], we treat H int as a renormalized low-energy effective interaction.We tabulate the allowed local gap functions, their symmetries, and the interactions in the respective pairing channels in Table I.Here, we adopt the common assumption of on-site rotational symmetry, which stipulates U = U +2J [45].This choice implies that all the OAST channels have the same attractive pairing interaction, which highlights the role of the normal-state Hamiltonian in selecting the most stable state.However, since the Ru sites have D 4h symmetry and not the assumed full rotational symmetry, the interaction strengths for the different pairing channels are generally different.Our results should therefore be in-terpreted as providing a guide to which superconducting states this form of attractive interaction can give rise to.
We write a free-energy expansion up to second order in the superconducting order parameter given by the gap matrices ∆i = ∆ i λ ai ⊗ σ bi (iσ 2 ), where i and j sum over all channels of a chosen irrep, g i are the corresponding interaction strengths from Table I, ω m = (2m + 1)πk B T are the fermionic Matsubara frequencies, and Ĝ = (iω m − Ĥ0 ) −1 and Ĝ = (iω m + ĤT 0 ) −1 are the normal-state Green's functions.Nontrivial solutions of the coupled linearized gap equations obtained from ∂F/∂∆ * i = 0 give the critical temperature T c and the linear combination of the ∆i corresponding to the leading instability.We include all channels in a chosen irrep, also repulsive ones, see Table I.In evaluating the last term in Eq. ( 3), we keep only intraband contributions; although the inclusion of interband terms will shift T c , this effect is negligible in the weak-coupling regime, as discussed in detail in the SM [40].
Results.-Weak-coupling OAST pairing states for attractive Hund's interaction require non-vanishing SOC [28,29,35].As described in the SM [40], SOC appears in five terms in the Hamiltonian Ĥ0 (k) in Eq. ( 1), representing a large parameter space to explore.We shall focus on the effects of the following terms: the zcomponent of the atomic SOC, h 43 = η z ; the in-plane atomic SOC, h 52 − h 61 = η ⊥ ; and the momentumdependent SOC associated with the interlayer hopping amplitude t SOC 56z between the d xy and the d xz and d yz orbitals, {h 53 , h 63 } = 8 t SOC 56z sin(k z c/2){cos(k x a/2) sin(k y a/2), − sin(k x a/2) cos(k y a/2)}.Here, we will ignore the anisotropy of the atomic SOC and set η z = η ⊥ = η.We have carried out a cursory exploration of the larger SOC parameter space and find that varying the other parameters within reasonable ranges such that the Fermi surfaces do not significantly deviate from the DFT predictions has little effect on the leading instability.
Figure 1(a) shows the phase diagram as a function of the atomic SOC η and the momentum-dependent SOC, parametrized by t SOC 56z .We find leading instabilities in the A 1g and E g channels.A 2g and B 2g states are not competitive anywhere in the phase diagram.A B 1g state is sometimes found as a sub-leading instability.The E g solution is dominated by the {[5, 3], [6,3]} channel and is stabilized for t SOC 56z η/4.Under the constraint of realistic Fermi surfaces, the E g state can be stabilized for t SOC 56z as small as about 5 meV, although this requires a rather small value of the on-site SOC.It is remarkable that such a small energy scale determines the relative stability of qualitatively different pairing states.As shown in Fig. 1(b), the Fermi surfaces for parameters stabilizing A 1g or E g states are indeed very similar.The SOC strength remains controversial [33,39,44], but here we have shown its importance for the determination of the most stable superconducting state.Our results are a proof of principle that an E g superconducting state can be realized in Sr 2 RuO 4 , even for purely local interactions, once one properly takes into account a complete and plausible 3D model for the normal state.Fig. 2 displays the projected gaps at the Fermi surfaces for representative A 1g and E g states.Note that in both cases the gap magnitude on the α sheet is very small, whereas the gaps on the β and γ sheets are comparable.This shows that we cannot simply identify the γ band [48] or the pair of almost one-dimensional α and β bands [46] as the dominant ones for superconductivity [47].
It is possible to understand why these SOC terms stabilize the respective ground states based on the notion of superconducting fitness [34,35].In particular, it has been shown for two-band superconductors that if the quantity FA (k) = ∆(k) H0 (k) + ∆(k) HT 0 (−k) is zero there is no intraband pairing and hence no weak-coupling instability [here, H0 (k) corresponds to Ĥ0 (k) with h 00 (k) set to zero].Hence, adding terms to the normal-state Hamiltonian such that FA (k) becomes nonzero for a particular gap function turns on a weak-coupling instability in this channel.The fitness analysis can be extended to our three-orbital model or, alternatively, we can construct an effective two-orbital model valid sufficiently far from the Brillouin-zone diagonals.Applying the fitness argument to the effective two-band model, we find that the on-site SOC η turns on both the A 1g and B 1g pairing channels, whereas the parameter t SOC 56z turns on the E g , { [5,3], [6,3]} channel, consistent with what we find numerically.Details of the fitness analysis are given in the SM [40].
In view of the Knight-shift experiments [17,18], it is important to comment on the spin susceptibility associated with the dominant E g , { [5,3], [6,3]} channel.Since it is a spin-triplet state with in-plane spin polarization of the Copper pairs, similar to the familiar chiral p-wave spin-triplet pairing with d-vector along the k z -direction, it might naively be expected to show a temperatureindependent spin susceptibility for in-plane fields.This is not the case, however, since the even parity of E g implies that the intraband pairing potential is a pseudo-spin singlet when expressed in the band basis and the low-energy response to a magnetic field is identical to a true spin singlet.This has been examined numerically for similar pairing states [30,49], where it was found that only a small fraction of the normal-state spin susceptibility persists at zero temperature in the superconducting state.
Bogoliubov Fermi Surfaces.-AnE g state is expected to have horizontal line nodes at k z = 0 and 2π/c [7,21].However, it has recently been shown that for an evenparity superconductor that spontaneously breaks timereversal symmetry, the excitation spectrum is either fully gapped or contains Bogoliubov Fermi surfaces (BFSs) [36,37].Indeed, the chiral E g state considered here has BFSs, which are shown in Fig. 3.These BFSs are very thin in the direction perpendicular to the normal-state Fermi surface, giving them a ribbon-like appearance that extends along the k z axis by about 0.2% of the Brillouin zone.This value is proportional to the gap amplitude, here set to 0.15 meV.While the total residual density of states from the BFSs is not large and may be difficult to observe [50], such a nodal structure implies that some experimental results require reinterpretation.In particular, given that the BFSs extend along the k z -axis, the argument that thermal conductivity measurements rule out the E g state because it has horizontal line nodes [10] no longer applies.The presence of BFSs may also require a reinterpretation of quasi-particle-interference experiments [51].We leave a detailed study of experimental consequences of the E g OAST state for future work.
Conclusions.-We have argued that an E g order parameter can be a realistic weak-coupling ground state for Sr 2 RuO 4 , once we consider a complete 3D model for the normal state and interactions of the Hubbard-Kanamori type.Key to our construction are the usually neglected even momentum-dependent SOC terms in the normal state.These terms can completely change the nature of the superconducting ground state, despite being so small that they do not significantly change the Fermi surfaces.Our theory reconciles the recent observation of a singlet-like spin susceptibility with measurements indicating a two-component order parameter and broken time-reversal symmetry.Han Gyeol Suh, Henri Menke, P. M. R. Brydon, Carsten Timm, Aline Ramires, and Daniel F. Agterberg

I. MICROSCOPIC MODEL
In this section, we construct a 3D tight-binding model for Sr 2 RuO 4 .We take into account the full 3D Fermi surfaces (FSs) of Sr 2 RuO 4 , based on the DFT band structure obtained by Veenstra et al. [S1], who showed that despite the 2D shape of the FSs, the orbital and spin polarization vary along k z .To account for the presence of orbital mixing on the different FS sheets, we include the t 2g manifold of the Ru d yz , d xz , and d xy orbitals (we will assume this order throughout).
We parametrize the orbital space by the the Gell-Mann matrices, which are the generators of SU(3).We use the convention We write the normal-state Hamiltonian in terms of the spinor In the presence of inversion and time-reversal symmetries, only a subset of fifteen h ab (k) terms are allowed.Table S1 lists the symmetry-allowed terms, the associated irrep for the matrices λ a ⊗ σ b , the physical process to which these correspond, and their momentum dependence.Note that Table S1 has entries which are in accordance with previous literature [S2, S3, S4] but there are also new terms associated with hopping along the z-direction or momentum-dependent SOC, which are usually neglected.Here we take η z = η ⊥ = η as the parameter for the on-site atomic SOC.The intra-orbital hoppings ξ We now focus on terms corresponding to k-dependent SOC, usually not taken into account in the standard parametrization of the normal-state Hamiltonian.The first matrix in the list, λ 5 ⊗ σ 1 + λ 6 ⊗ σ 2 , which is of A 2g symmetry, will be ignored because the lowest-order polynomial basis function of this irrep is of order 4 (g-wave), which only appears at next-next-next-nearest-neighbor hopping and is therefore assumed to be negligible.We also take the other k-dependent SOC terms at the lowest order at which they appear.This concludes the construction of the microscopic model, which is characterized by a Hamiltonian with 26 free parameters.

II. FIT TO DFT RESULTS
We employ the tight-binding model presented in the supplemental material of [S1] to determine the free parameters.The tight-binding Hamiltonian is derived from an LDA band structure that is down-folded onto the O-2p and the Ru-4d orbitals and therefore has a total of 17 bands.The hopping integrals are truncated at 10 meV.We henceforth refer to the LDA-derived tight-binding Hamiltonian as the "DFT model".For the calculation of the linearized gap equation, the DFT model is much too large and most of the bands are irrelevant for superconductivity.The states at the Fermi surface are determined by the t 2g manifold of the Ru-4d orbitals (d yz , d xz , d xy ) and we fit Eq. (S2) to several quantities extracted from the DFT model projected into this subspace.
We extract the Fermi momenta kF of the DFT model and denote the eigenvalues by and the associated eigenvectors by V .We define the following measure where the sum is over momenta kF on the DFT Fermi surfaces formed by the bands n = α, β, γ, n (k) are the band energies, d n xy (k) is the d xy -orbtial content, p n SOC (k) is the spin polarization, and v n (k) the in-plane velocity.Quantities with a tilde are from the DFT model.The d xy -orbtial content is determined by the corresponding eigenvector components The spin polarization is determined from the expectation value of the atomic spin-orbit coupling Hamiltonian H SOC = λ 5 σ 2 − λ 6 σ 1 − λ 4 σ 3 : For the in-plane Fermi velocity we use a simple two-point central finite differences stencil where ε x,y are small We minimize the measure (S6) using the derivative-free optimization algorithm BOBYQA of dlib [S5].
The fit yields very good agreement with the DFT model close to the Fermi energy, including good reproduction of the d xy -orbital content and the spin polarization.In Fig. S1, we compare the result of our fit with the DFT model in the k z = 0 plane.In Fig. S2, we show the full 3D Fermi surface produced by our fit, together with the d xy -orbital content and the spin polarization.The corresponding fit parameters are listed in Tab.S2.
It is important to note that because the different sheets of the Fermi surface have varying orbital and spin content, it is not possible to isolate one dominant band for superconductivity.The pairing state will in general have contributions from all three sheets.

III. LINEARIZED GAP EQUATION
In this section, we outline our solution of the linearized BCS gap equation.For convenience, we repeat the secondorder expansion of the free energy as given in Eq. ( 3) of the main text, where the gap functions are ∆i = ∆ i λ ai ⊗ σ bi (iσ 2 ) and the indices a i , b i , and interaction energies g i are given in Table 1 of the main text.We introduce an interaction scaling parameter s, and for concreteness choose the interaction energies to be given by U = 5/s, U = 1/s, and J = 2/s.Since we are interested in the weak-coupling limit we will later assume s to be large.The Green's functions and gap function are expressed in the energy eigenbasis by where U is a unitary matrix that diagonalizes the normal-state Hamiltonian Ĥ0 , U = (iσ 2 ) † U , and a are band energies.We define new gap matrices by The frequency summation yields where β = 1/k B T .The linearized gap equations are obtained by differentiating the free energy with respect to the gap amplitudes, ∂F/∂∆ * i = 0, written explicitly as where i and j run over all gap-structure indices of a given irrep, a and b run over band indices, [Λ i ] ab is a matrix element of Λ i , and gi is the value of g i when s = 1.First, consider a = b intraband terms in Eq. (S16), the k-integration is written as with Differentiating Eq. (S17) with β gives where F abij ( ) = F abij (0) + F abij (0) + . . .has been used.Note that when this is integrated with respect to β it yields the log β divergence in Eq. (S17).Next, consider the a = b interband terms, as β → ∞, S ab converges to θ( a b )/| a + b | , where θ is the Heaviside step function.Because this is a bounded function, there is no divergence in the interband contributions.In Eq. (S16) a non-trivial solution for the gap amplitudes ∆ j is found by considering i and j as matrix and taking the corresponding 6 × 6 matrix to be singular.Including both the intraband and interband contributions, the critical where ] ab S ab that remains after removing the log β divergent term.By definition, C ij (β) is convergent as β → ∞, so the last term in the determinant can be ignored when s is sufficiently large.More explicitly, in the weak-coupling limit, T c is given in the form where m is the smallest log β c /s solution found when C ij = 0 and δ(s) is a function that approaches a constant as s → ∞.Different channels (irreps) have different values of m, and the channel with the smallest m is the leading instability in the weak-coupling limit.Note that the definition of m does not depend on C ij and all the interband contributions go into C ij .Thus we can drop the interband terms in Eq. (S16).This changes δ(s) but does not change m.The resultant expression is which is the equation we solve numerically.The log β divergence originates from momenta near the Fermi surface, so we carry out the k-integration on adaptive meshes with finer resolution near the Fermi surface.We obtain log β c for several values of s and use linear regression to get the slope, which determines m.If the values of log β c at the sampling points are not linear in s we sample larger s values until we encounter linear behavior.In our calculation, an equidistant set of four sampling points is used for a linear regression and their R 2 measures are always greater than 0.999.Using this procedure, we get the slope m for each pairing channel and determine the leading instability at each point in the phase diagram displayed in Fig. 1(a) in the main text.While this procedure may seem more elaborate than a direct solution of Eq. (S20) with C ij = 0, it allows us to verify Eq. (S21) showing that our solution is in the weak coupling limit.

IV. SUPERCONDUCTING-FITNESS ANALYSIS
In this section, we present details of the superconducting-fitness analysis.We start with the more realistic threeorbital model and then consider an effective two-orbital model, which dramatically simplifies the analysis but gives consistent results.

A. Complete 3D three-orbital model
In previous works [S6, S7], a proof of the direct relation between the superconducting-fitness measures FC (k) and FA (k) (defined below) and the superconducting critical temperature was provided for the one-and two-orbital scenario.The first measure, quantifies how incompatible a given gap structure is for a specific normal state, namely, how much inter-band pairing there is.Here, H0 (k) = Ĥ0 (k) − h 00 (k) λ 0 ⊗ σ 0 and we have defined a normalized gap matrix ∆(k) = ∆(k)/| ∆(k)| such that average over the normal-state Fermi surface is ∆(k) ∆ † (k) FS = 1.The second measure, quantifies how much intra-band pairing there is, or what fraction of the gap survives upon projection onto the Fermi surface.For the two-orbital scenario, these measures satisfy Tr F † A (k) FA (k) + Tr F † C (k) FC (k) FS = 1, up to normalization of the normal-state Hamiltonian, which highlights their complementarity.The proof of this relation relies on the fact that the matrices associated with the orbital DOF are Pauli matrices for the two-orbital scenario and therefore form a totally anticommuting set, which greatly simplifies the calculations.On the other hand, for n > 2 orbitals, the basis matrices are the generators of SU(n), which do not form a totally anticommuting set and therefore do not allow a direct generalization of this relation for models with more than two orbitals.However, the physical meaning of FC (k) and FA (k) is preserved within some approximations, as discussed below.
For the three-orbital situation, the corresponding superconducting-fitness functions can be identified as Given that [A 2 , B] ± = A[A, B] ∓ [A, B]A, the core of the analysis still depends on the original form of the superconducting-fitness functions.Therefore, we use the form linear in Ĥ0 (k) to get some insight.Below, we will see that a simplified two-orbital model, for which the linear version of the fitness functions is valid rigorously, corroborates our analysis.We summarize the results for the complete three-orbital problem in Tables S3 and S4.The first row gives the irrep of each term in the normal-state Hamiltonian displayed in the second row as h ab .The first column gives the irrep of each order parameter displayed in the second column following the notation [a, b] corresponding to ∆ = λ a ⊗ σ b (iσ 2 ).The third column gives the interaction stemming from the Hubbard-Kanamori Hamiltonian for each channel.Finally, the numerical entries in the tables correspond to Tr F † A,C FA,C = cd (table entry) |h cd | 2 .Note that the order parameters with a potentially attractive interaction U − J are [4,3] and [5,2] 6,2] in B 2g , and finally { [4,2], −[4, 1]} and { [5,3], [6,3]} in E g .All these order parameters are associated with spin-triplet states.If we focus first on the largest terms in the normal-state Hamiltonian, namely h 80 and h 70 (intra-orbital hopping), h 10 (inter-orbital hopping), and h 43 and h 52−61 (atomic SOC), we conclude from Tables S3 and S4 that, among the one-dimensional irreps, the most stable state should be in the A 1g channel since these states are associated with larger entries for FA and smaller entries for FC .
Considering now the two-dimensional order parameters, for { [4,2], −[4, 1]}, we find that the terms stabilizing it, i.e., the ones with the largest contribution to FA , are h 51+62 , h 52+61 , and {h 42 , −h 41 }, all associated with momentumdependent SOC.However, these terms contribute with the same value to the detrimental fitness measure FC , suggesting that they overall do not favor this pairing state.On the other hand, the two-dimensional order parameter { [5,3], [6,3]} is stabilized by h 51+62 , h 52+61 , and {h 53 , h 63 }.Again, the terms h 51+62 and h 52+61 contribute with a large value to FC .On the other hand, the terms {h 53 , h 63 } contribute only moderately.This analysis suggests that the {[5, 3], [6,3]} channel should be the one driving the superconducting instability in the E g channel and can be stabilized by large terms {h 53 , h 63 }.
From this analysis, we can understand the tendencies observed in the numerical results as follows: the order parameters in A 1g , in particular [5, 2] − [6, 1], are strongly stabilized by atomic SOC, in particular by the term h 52−61 in the normal-state Hamiltonian, such that reducing the magnitude of this coupling is expected to weaken the superconducting instability in this channel.Moreover, the terms {h 53 , h 63 } primarily suppress the order parameter in this channel since their contribution to FC is larger than the one to FA .In contrast, the E g order parameters, in  S4.Superconducting-fitness measure FC for the 3D three-orbital model for Sr2RuO4.The same notation as in Table S3 has been used.
particular for { [5,3], [6,3]}, are primarily stabilized by {h 53 , h 63 } since for these terms the contribution to FA is larger than the one to FC , while atomic SOC is clearly detrimental.This analysis suggests that by reducing the atomic SOC and enhancing the terms {h 53 , h 63 } associated with nonlocal SOC even in momentum, the ground state should change from A 1g to E g .

B. Effective two-orbital model in the kxkz-plane
Sufficiently far from the Brillouin-zone diagonals k y = ±k x , the bands close to the Fermi energy are dominated by only two of the Ru d-orbitals.For concreteness, here we consider the k x k z -plane, but our conclusions remain qualitatively valid for general k, except close to k y = ±k x .
In the k x k z -plane, the dominant orbitals at the Fermi energy are d xz and d xy .Projecting into this subspace, we obtain an effective two-orbital Hamiltonian which is parametrized by where the hab (k) are real functions of momentum, τ a and σ b are Pauli matrices for a, b = 1, 2, 3 and the 2 × 2 identity matrix for a, b = 0, encoding the orbital and the spin DOF, respectively.There are, in principle, 16 parameters hab (k) but in the presence of time-reversal and inversion symmetries these are constrained to only six, including the term proportional to the identity.The symmetry-allowed terms are listed in Table S5; we classify them in terms of the irreps of D 2h , which is the little group for D 4h in the k x k z -plane.Analogously, we can parametrize the s-wave gap matrices in the orbital basis as The irreps associated with each [a, b] combination are the same as for the normal-state Hamiltonian, given in the first two columns of Table S5.
The superconducting-fitness analysis, which is summarized in Table S6, is very much simplified in the two-orbital scenario since the symmetry-allowed matrices form a totally anticommuting set.From the table, one can see that the results concerning FA (k) and FC (k) are complementary.Note that the trivial order parameter, [0, 0], is stabilized by all the terms in the Hamiltonian while the remaining order parameters of the form [a, b] need the associated term hab in the Hamiltonian to develop a weak-coupling instability.There is an attractive interaction in the orbital-singlet spin-triplet channels [2, b].The order parameter [2,1] in A g is stabilized by the atomic SOC term h21 .The other

FIG. 1 .
FIG. 1.(a) Phase diagram showing the stability of A1g and Eg pairing states as a function of the SOC parameters η and t SOC 56z .The vertical dashed lines indicate the minimum distance between two Fermi surfaces.Percentages are defined as fractions of 2π/a.For small η, the separation between the β and γ bands becomes too small, in view of the ARPES data [33].The thin solid lines indicate the maximum variation of the Fermi surface along the kz-direction.For large t SOC 56z , the Fermi surfaces become too dispersive.The blue and magenta dots denote the parameter choices for Eg and A1g stable solutions used in panel (b).(b) Fermi surface shapes, projected onto the kxky-plane, for representative points in the A1g (red) and Eg (blue) regions in (a).For A1g, η = 57 meV and t SOC 56z = 10 meV, while for Eg, η = 40 meV and t SOC 56z = 12 meV.

FIG. 2 .
FIG. 2. Projected gaps at the Fermi surfaces for a representative (a) A1g and (b) Eg state in the first Brillouin zone.Parameters are the same as in Fig. 1 (b).The color code is normalized to the maximum value of the A1g gap.

FIG. 3 .
FIG.3.BFSs for the chiral Eg state.The Fermi surfaces in red, green, and blue correspond to inflated nodes stemming from the α, β, and γ band, respectively.
where we introduce the numbering of the orbitals 1 = d xz , 2 = d yz , 3 = d xy .In terms of the Gell-Mann and Pauli matrices, we write the Hamiltonian H 0

FIG. S1 .
FIG. S1.Comparison of the DFT model (red dashed lines) with our fit (blue solid lines) in the kz = 0 plane.(a) Fermi surface in the first quadrant of the Brillouin zone.(b) dxy-orbital content, (c) spin polarization, and (d) in-plane velocity as functions of the angle θ = arctan(ky/kx) in the first quadrant.The three columns pertain to the α, γ, and β band.
FIG. S2.Full 3D Fermi surface obtained from the best fit to the DFT results of Ref. [S1], with color indicating (a) the dxy-orbital content and (b) the spin polarization.The three columns pertain to the α, γ, and β band.

o
Desenvolvimento da UNESP (FUNDUNESP) process 2338/2014-CCP.C. T. acknowledges financial support by the Deutsche Forschungsgemeinschaft through the Collaborative Research Center SFB 1143, Project A04, the Research Training Group GRK 1621, and the Cluster of Excellence on Complexity and Topology in Quantum Matter ct.qmat (EXC 2147).H. G. S. and H. M. contributed equally to this work.Material for Stabilizing Even-Parity Chiral Superconductivity in Sr 2 RuO 4

TABLE S1
. List of the fifteen symmetry-allowed terms in the normal-state Hamiltonian Ĥ0(k) in Eq. (S2).For each (a, b), the basis function h ab (k) must belong to the same irrep of D 4h as the matrix λa ⊗ σ b .The table gives the irrep, the associated physical process ("Type"), where "k-SOC" means momentum-dependent (nonlocal) SOC, and the momentum dependence of h ab (k).For the two-dimensional irrep Eg, the entries are organized such that the first transforms as x and the second as y.
11,22,33 (k) are included up to next-next-nearest neighbors in plane and next-nearest neighbors out of plane.The inter-orbital hopping λ(k) between the d xz and the d yz orbitals is kept up to next-nearest neighbors in plane and nearest neighbors out of plane.For the inter-orbital hopping {(3, 0), −(2, 0)} between the d xz and d xy (d yz and d xy ) orbitals, we only keep the nearest-neighbor component out of plane.The explicit form of the functions not given explicitly in Table S1 is x a cos k y a + cos 2k y a cos k x a) + 2t33zz (cosk z c − 1) − µ 1 ,(S4)λ(k) = 4t 12 z sin(k x a/2) sin(k y a/2) cos(k z c/2) − 4t 12 xy sin k x a sin k y a − 4t 12xxy (sin 2k x a sin k y a + sin 2k y a sin k x a).

TABLE S2 .
) Parameters of the Hamiltonian (S2) determined from the fit to the DFT model.All values are in meV.

TABLE S3 .
) Superconducting-fitness measure FA for the 3D three-orbital model for Sr2RuO4.The first column gives the irreps of the order parameters associated with matrix form [a, b] (second column) and the third column displays the local interaction in the respective channel, where the potentially attractive channels are highlighted in boldface.Columns 4-17 give the results for the fitness function such that Tr F † A FA = cd (table entry) |h cd | 2 , for each term h cd in the normal-state Hamiltonian, indicated in the second row with the associated irrep given in the first row.We highlight in boldface the h cd terms which are present in the standard 2D models for Sr2RuO4, while the terms in normal typeface are either momentum-dependent SOC or interlayer couplings.