We study the electronic structure and correlated phases of twisted bilayers of platinum diselenide using large-scale *ab initio* simulations combined with the functional renormalization group. PtSe_{2} is a group-X transition metal dichalcogenide, which hosts emergent flat bands at small twist angles in the twisted bilayer. Remarkably, we find that Moiré engineering can be used to tune the strength of Rashba spin–orbit interactions, altering the electronic behavior in a novel manner. We reveal that an effective triangular lattice with a twist-controlled ratio between kinetic and spin–orbit coupling (SOC) scales can be realized. Even dominant SOC can be accessed in this way and we discuss consequences for the interaction driven phase diagram, which features pronounced exotic superconducting and entangled spin-charge density waves.

## 1.Introduction

The advent of Moiré heterostructures and the demonstration of superconducting (SC), correlated insulating and topological phases of matter in these materials [1–16], has triggered a surge of theoretical and experimental studies. Common to these is the idea that a slight lattice constant mismatch or rotation between adjacent layers of two-dimensional van der Waals materials can significantly quench kinetic energy scales and alter the effective electronic band structure that dictates the low-energy behavior. Moiré heterostructures were envisioned to allow novel control of the ratio between kinetic and potential energies by superlattice engineering, allowing an exploration of strong electronic correlations in a tunable condensed matter setting [17]. Recent experimental findings suggest that superconductivity in twisted sheets of bilayer graphene is indeed of unconventional nature [18] and that twisted trilayer graphene may favor triplet pairing [19–22]. Substantial efforts have been made to unravel the nature of correlated states in related graphitic Moiré materials. Experiments report, among others, correlated states in twisted mono-bilayer graphene [23, 24], twisted double bilayer graphene [25–29] and rhombohedral graphene aligned with hexagonal boron-nitride [30–32]. A vast amount of theoretical work predicts correlation effects for an even larger subspace of the possible twisted graphitic Moiré systems [33–48]. In addition, in twisted bilayer graphene, control over topological properties has already been demonstrated [49–58] showing correlated Chern insulating phases.

Beyond graphene, experiments have studied twisted sheets of transition metal dichalcogenides (TMDs), concentrating primarily on group-VI homo- or hetero-bilayers of MoS_{2}/Se_{2} or WeS_{2}/Se_{2}, with fascinating observations of strongly correlated behavior [59–69] and excitonic physics [70–76]. Further proposals for TMD materials include exotic SC states with topological features [36, 60], possible spin-liquid phases [61] and engineering of multi-orbital systems in group-IV TMDs as a realization of the elusive Kagome lattice with strong and tunable spin–orbit coupling (SOC) which exhibits fractional quantum anomalous Hall and Chern insulating states [77].

However, going beyond graphene- and TMD-based systems, the profusion of available van der Waals materials allows for even more exotic quantum materials design. For instance, the reduced rotation symmetry of monochalcogenides permits engineering quasi-one dimensional structures [78]; alternatively, by departing from the realm of few-layer systems, Moiré induced control of three dimensional materials becomes possible [37]. By considering oxides as a basis for Moiré engineering, exotic *d*-wave superconductivity is supposed to emerge [79] with a potential connection to the fascinating high-*T*_{c} phase of the cuprates. All of this is to show that in the field of Moiré engineering much more is expected to be possible by exploiting the different chemical compositions offered by the choice of materials to consider. This general concept of identifying novel phenomena to be controlled by Moiré engineering might culminate into a versatile new solid state-based platform [17] to access quantum materials behavior with unprecedented level of tunability. The discovery and characterization of limits and opportunities in new van-der-Waals materials platforms hence remains an important avenue of pursuit.

Here, we add to the catalog of phenomena realizable by Moiré engineering by considering the group-X TMD PtSe_{2}, which is exfoliable down to monolayer [80]. It has raised lots of interest for its outstanding optical and electrical properties and high air-stability [81, 82]. In the context of Moiré engineering, this material is interesting due to the substantial SOC of heavy transition metal ions. We demonstrate via an *ab initio* characterization of large unit-cell systems at small twist angles that, when the kinetic energy scales are quenched by twisting two sheets of PtSe_{2} with respect to each other, a controlled twist-dependent tuning of Rashba SOC and kinetic energy scales can be achieved. Surprisingly, we find that relatively large twist angles of about 6° are sufficient to quench kinetic energy scales small enough to promote SOC to be the dominant energy scale. In contrast to the strong SOC of twisted bilayers of ZrS_{2} [77], the SOC interaction in PtSe_{2} is mainly of Rashba type and relies existentially on broken inversion symmetry in the Moiré superstructure, hence realizing a new regime. We discuss consequences for correlated phases of matter using a weak-coupling functional renormalization group (fRG) approach, which can be viewed as an unbiased renormalization-group-enhanced random phase approximation. Our results indicate a rich phase diagram of intertwined charge-spin density waves (DWs), which in the case of SOC cannot be disentangled, and exotic mixed-parity SC phases with topologically non-trivial properties.

The paper is structured as follows: we start from a full *ab initio* characterization of the twist angle dependence of the electronic band structure including the SOC for twisted bilayers of PtSe_{2}. We demonstrate the twist-dependent reduction of the effective electronic bandwidth, which coincides with the emergence of strong Rashba interactions. The resulting Moiré bands span a triangular lattice with few nearest neighbor hoppings plus Rashba SOC. We then treat this model by adding a Hubbard interaction and outline the emerging phase diagram. A discussion concludes the paper, with details of the methods used appended below.

## 2. *Ab initio* characterization

To provide a first-principle characterization of the electronic structure of twisted bilayers of PtSe_{2} we first employ a density functional theory (DFT) based approach to the material, which at small twist angles can exhibit a very large unit cell (see methods). As Pt is a heavy element, SOC is important and it is included in the calculations, which breaks the underlying *SU*(2) spin rotation symmetry (see methods for more details). PtSe_{2} is a group-X TMD, and the atomic structure of a monolayer is visualized in figure1(a) as a top and side view. The Pt and the Se atoms are shown as gray and green spheres, respectively. We concentrate on the energetically stable1T lattice structure of PtSe_{2}. In panels (b) and (c) we report the DFT bands structure of the monolayer (b) and the untwisted bilayer (c) for reference. As shown in the figures, the states at the valence band edge are dominated by Se *p*_{x }, *p*_{y } and *p*_{z } orbitals. Among them, we find that the Se *p*_{z } states are very sensitive to interlayer coupling, which strongly hybridize in the bilayer and form bonding and anti-bonding states with a large energy splitting (indicated by the orange arrow in figure1(c)). Consequently, the top of the valence bands shifts from the Se *p*_{x } and *p*_{y } states in the monlayer to the Se *p*_{z } states in the bilayer. This is different from the case of another 1T TMD, ZrS_{2} [77], in which the top of the valence bands is dominated by *p*_{x } and *p*_{y } states, in both monolayer and bilayer. Because of the strong interlayer coupling of the Se *p*_{z } states at the valence band edge, we expect that flat bands can be formed in twisted bilayer of PtSe_{2} at relatively large twist angles.

Next we turn to non-zero twist angles. Panel (d) shows the real space Moiré pattern emerging for twisted PtSe_{2}. We mark different regions as 'AA', 'AB' and 'BA' stacking, with the local stacking arrangement of the atoms given as insets to the side. In panels (g) and (h) we report the DFT analysis of the band structure in the twisted system excluding (g) and including (h) SOC at a twist angle of 6°. Panels (e) and (f) show the DFT band structure including SOC at twist angles 9.43° and 7.34°, respectively. As twist angles are approached, the electronic band near the Fermi energy becomes very flat with a width of (6) meV in the calculations without (with) SOC. This is significantly lower than the corresponding bandwidth for twisted bilayer graphene [1], especially at such relatively large twist angle. Comparing the calculations with and without SOC, the relevance of including the latter becomes strikingly clear. The degeneracy of the up- and down-spin electronic bands without SOC is lifted by including the SOC in the calculations, by virtue of broken inversion symmetry. To quantify this effect we will next analyze the relevance of SOC versus kinetic energy scales in dependence of the twist angle using a tight-binding approach. As the bandwidth of twisted bilayer PtSe_{2} at 6° is small enough for correlation effects to be relevant, we do not perform further DFT calculations for smaller twist angles, which are increasingly expensive. Nevertheless, we expect the bandwidth to be further reduced at smaller angles.

## 3.Tight-binding description and strong spin–orbit coupling

We model the electronic flat bands with a tight-binding model on the triangular lattice, taking into account hopping parameters connecting up to third nearest neighbors. Due to the combination of broken inversion and *SU*(2) spin rotation symmetry, we capture strong spin orbit coupling of the flat bands using a Rashba term. Additionally, the non-*SU*(2) nature of the system can lead to intrinsic, spin-dependent electric field effects described by a complex spin-dependent phase accompanying the kinetic hopping parameters [67, 83]. The kinetic part of the Hamiltonian then reads

where *B*_{ij } is a set of symmetry related directed bonds ** b **

_{ij }in the triangular lattice with equal length, are the kinetic hopping parameters and

*ϕ*

_{ij }is the Kane–Mele phase factors. By construction, this Hamiltonian fulfills both time reversal and symmetry. The Rashba term is given by

with the SOC hopping parameters. Finally, we include a chemical potential (*μ*), such that the full tight binding Hamiltonian becomes

Note that we do not account for changes in tight-binding parameters when varying the chemical potential, such that a change in *μ* is directly reflected in a change in filling *ν*.

Figure2(a) shows both the DFT band structure (black crosses) and the tight-binding fit (line) for a twist angle of *θ* = 6°. We list the fit parameters for 6° as well as for two additional twist angles (7.34° and 9.43°) that we calculated using DFT (cf figures1(e) and (f)), in table1 in the methods section. Additionally, we show the spin expectation value as three-dimensional arrows in figure2(a): along the path Γ–*K*, the spin has finite expectation value in the *x*–*y* plane which then gradually shifts towards spin-*z* at the *K* point. From *K* to *M*, the expectation value of the spin-*z* component is nonzero with a slight tilt towards the *x*–*y* plane close to *M*. From *M* to Γ, the expectation value fully lies in spin-*x* direction. The possibility to generate a finite expectation value of a specific spin component aside from *S*_{z } arises from making a specific choice of the part Γ–*M* of the irreducible path. We label the lower band by *b*_{1} and the upper band by *b*_{2} with figure2(b) showing a two-dimensional false color plot of the dispersion for the same twist angle. The hexagonal BZ is indicated as gray lines. Strong spin orbit coupling and the very small band width of approximately 6meV are clearly visible. The only degeneracy points of the band structure lie on the BZ boundary at the *M* points and at Γ.

**Table 1.**Tight-binding fit parameters of the Moiré flat bands of twisted PtSe_{2} for three twist angles.

θ | μ(eV) | |||
---|---|---|---|---|

6.00° | 2.52 × 10^{−3} | 5.62 × 10^{−4} | −1.18 × 10^{−4} | −6.37 × 10^{−5} |

7.34° | 1.96 × 10^{−2} | 3.41 × 10^{−3} | −2.57 × 10^{−5} | −9.25 × 10^{−5} |

9.43° | 1.13 × 10^{−2} | 2.55 × 10^{−3} | −1.42 × 10^{−3} | −2.46 × 10^{−4} |

θ | ϕ_{1} | ϕ_{2} | ϕ_{3} | |

6.00° | 3.56 × 10^{−1} | −2.09 × 10^{−10} | 4.53 × 10^{−1} | |

7.34° | 1.63 × 10^{−1} | 1.11 × 10^{0} | 2.43 × 10^{−2} | |

9.43° | 6.83 × 10^{−1} | −1.37 × 10^{−1} | −1.23 × 10^{0} | |

θ | ||||

6.00° | 1.94 × 10^{−4} | −8.43 × 10^{−6} | 4.45 × 10^{−6} | |

7.34° | 1.97 × 10^{−4} | 5.85 × 10^{−6} | −1.55 × 10^{−5} | |

9.43° | 4.89 × 10^{−4} | 2.13 × 10^{−4} | 2.12 × 10^{−4} |

After having established an accurate tight-binding representation of our DFT results, we can quantify the strength of SOC as a function of twist angle by carrying out the fitting procedure at two other (commensurate) twist angles. The resulting kinetic and SOC hopping parameters *t*_{1}, *t*_{2}, *t*_{3} are shown in figure2(c). For all three twist angles considered in the scope of this work, we see that SOC is extremely relevant. Since the Moiré potential becomes increasingly relevant at smaller twist angles, the overall kinetic energy scale given by the nearest neighbor hopping is drastically reduced for *θ* = 6°. As a consequence, the SOC hopping parameter is around 40% of the non-SOC *t*_{1}. Furthermore, the influence of longer range hoppings *t*_{2} and *t*_{3} becomes smaller when decreasing the twist angles.

With this in mind, we continue our analysis of the tight binding model and complement the quenched kinetic energy Hamiltonian *H*^{0} with onsite Coulomb interactions *H*^{U }:

In the following, we will study the effect of *H*^{U } on the non-interacting Moiré Hamiltonian *H*^{0}. We therefore focus on *θ* = 6° as twist angle for two reasons. First, kinetic energy scales are strongly suppressed and second, the quality of the tight-binding fit is the best due to long range hoppings being least relevant (among the cases studied within this work).

## 4.Interaction-driven phases of matter

We approach the interacting quantum many-electron problem using the unbiased fRG [84]. The broken *SU*(2) symmetry renders even this two-band problem a significant challenge and we truncate the infinite hierarchy of flow equations set up within the fRG approach at the four-point vertex Γ^{(4)}. Furthermore, focusing on static quantities, we neglect frequency dependencies of Γ^{(4)} and further set the two-point vertex (self-energy) to zero. The fRG flow then amounts to solving a differential equation (see methods) for Γ^{(4)} as a function of Λ, the parameter that smoothly interpolates from the free theory at Λ = **∞** to the full, interacting theory at Λ = 0. During the flow, we search for divergences in Γ^{(4)} indicating a tendency towards long-range order. The four-point vertex is then analyzed at the final scale Λ_{c} that roughly corresponds to a critical temperature of the phase transition associated with the divergence.

Within our approach we can distinguish between charge/spin-DW or SC instabilities. The primary indicator for the type of divergence is given by the divergent channel during the fRG flow, which can either be of particle–particle (SC) or particle–hole (DW) type. If the vertex remains finite up to Λ = 0, the fRG does not predict long range order and thus a metallic phase.

Figure3(a) shows the resulting phase diagram for the flat bands of twisted bilayer PtSe_{2} at *θ* = 6° as a function of Moiré–Hubbard interaction strength *U* and filling factor *ν* parametrized by the chemical potential *μ*. The upper panel displays the system's density of states (DOS) with two main van Hove singularities at *ν* ≈ 0.4 and *ν* ≈ 1.0. These regions of high DOS are responsible for the instabilities. However, the lower panel reveals that both SC and DW ordering can emerge away from points with divergent DOS. Moreover, there is a rich phase structure with various regions of SC and DW order at a broad range of critical scales (encoded in color). The DW instabilities predominantly occur at large *U*, whereas the SC instabilities are *dominant* for a broad range of fillings (*ν* = 0.25–0.5) and interactions (*U* = 3.5–5meV), and are only flanked by metallic regions. At larger interaction strengths, we observe a second DW instability with very low critical scale at *ν* = 0.2 and a second SC instability driven by the van Hove peak at *ν* = 1.0.

Figures3(b) and (c) illustrate the two types of SC instabilities found. First, for most of the central SC region and the remote region at *ν* = 1.0 in the phase diagram (figure3(a)), the leading instability is of *d*-wave and *p*-wave type (see figure3(b)). The *d* and *p* wave symmetries mandate that the SC gap is doubly degenerate (second instability not shown). We plot the SC gap amplitude in the BZ for both the singlet [*ψ*(** k **)] and triplet [

**(**

*d***)] channel. In the case of strong SOC and lack of inversion symmetry, the decoupling of a SC order parameter into independent singlet and triplet components is impossible and mixed-parity SC states form. Nevertheless, we can transform the SC gap Δ**

*k*_{σσ'}(

**) to singlet and triplet space (see methods) [85, 86], but with instabilities that have weight in both spaces at the same time. In momentum space,**

*k**ψ*(

**) and**

*k***(**

*d***) must fulfill the (anti-)symmetry relations of singlet (triplet) gaps. Second, for the SC instabilities at**

*k**ν*= 0.25, we find a different order parameter (see figure3(c)) with dominant

*g*-wave (singlet) and

*f*-wave (triplet,

*d*

_{z }) components and little weight in the

*d*

_{x }and

*d*

_{y }components of the triplet vector. This order parameter is not degenerate and, by its

*g*/

*f*-wave symmetry, leads to a nodal SC state.

To gain an intuitive understanding of the weight distribution in the ** d **-vector, it is helpful to consider the spin polarization of the Fermi contour. In the presence of strong SOC (equation(2)) and time-reversal symmetry an arbitrary single particle state |

*k*_{F}, ↑⟩ at Fermi momentum

*k*_{F}and with spin

*σ*= ↑ measured relative to the spin-polarization axis at

*k*_{F}can be transformed to . Since inversion symmetry is no longer conserved in the system, the states |

*k*_{F}, ↑⟩ and |−

*k*_{F}, ↑⟩ are no longer required to be degenerate [86] such that opposite spin Cooper pairs are favored in this situation [87]. Indeed, we observe that at filling

*ν*= 0.4, the bands at the Fermi energy are mostly

*z*-polarized leading to dominant weight in the singlet

*ψ*(

**) and triplet**

*k**d*

_{z }(

**) component, while the other components of the**

*k***-vector are substantially suppressed. As the system is filled with more electrons, the spin polarization axis changes from**

*d**z*to the

*x*–

*y*plane and consequently shifts the weight in the

**-vector.**

*d*At and close to the central van Hove singularity, the system is susceptible to an intertwined magnetic/charge density order (i.e. divergence in the particle–hole channel) presented in figure3(d). The physical spin- and density channels *χ*^{lm }(** q **) with

*l*,

*m*∈ {0,

*x*,

*y*,

*z*} are obtained from the four-point particle–hole susceptibility which is in turn calculated from the vertex Γ

^{(4)}at the critical scale Λ

_{c}(see methods) [88]. By virtue of the strong SOC, the density–density response

*χ*

^{00}(

**) (first panel) and the spin–spin responses are intrinsically coupled. Most weight in the spin sector is in**

*q**χ*

^{zz }(

**), i.e. the spin-**

*q**z*response with dominant ordering vectors

*K*and

*K*'. The in-plane responses are much weaker; albeit they are non-negligible.

*χ*

^{yz }(

**) is not shown in figure3(d), though its form can be inferred from**

*q**χ*

^{xz }(

**) by symmetry. The less dominant DW instability at**

*q**ν*= 0.2 is of qualitatively different type with dominant terms in the density–density and

*χ*

^{zz }sectors for transfer momentum

**= 0.**

*q*We further investigate the physical consequences that arise when the system is in the intertwined *d*/*p*-wave SC phase at *ν* = 0.4. To examine which linear combination of the two degenerate instabilities is energetically favored, we calculate the free energy in the SC phase (see methods) as a function of all possible complex superpositions of the two *d*/*p*-wave order parameters:

Figure4(a) indicates that the free energy in the SC phase is minimized for *ϑ*_{0} ≈ *π*/4 and *φ*_{0} ≈ *π*/2. This particular choice of *ϑ*_{0} and *φ*_{0} leads to an order parameter that preserves *C*_{3} rotational symmetry in the *ψ*(** k **) and

*d*

_{z }(

**) components while breaking time reversal symmetry. We display the associated chiral SC order parameter in figure4(b), where we encode the complex phase as color and the magnitude as lightness. Next, we set an amplitude of |Δ|**

*k*_{max}= 2.35meV and study the Bogoljubov–de-Gennes (BdG) bandstructure

*E*

_{b }(

**) in figure4(c). To assess whether the system is topologically nontrivial, we first numerically diagonalize the BdG Hamiltonian (see methods) on a cylindrical geometry. We periodically continue the system in**

*k*

*a*_{2}direction and open the boundary in the

*a*_{1}direction. Further, we color-code the IPR as a function of band index

*b*and momentum in

*a*_{2}direction (

*k*

_{y }). Additionally, we determine the Chern number of the two upper BdG bands, where

*C*= +2 and the two lower BdG bands with

*C*= −2. This leads us to conclude that the SC order is topologically non-trivial.

## 5.Discussion

Our results elevate twisted bilayer PtSe_{2} as a novel platform for engineering strong Rashba SOC in a tunable setting. Importantly, the strong SOC regime can be accessed in a controllable fashion, allowing a novel inroad into this evasive physical regime. We discussed consequences of the strong SOC as well as the exotic form of the engineered low-energy effective Hamiltonian, which shows prominent effects of the breaking of the *SU*(2) symmetry even without interactions. The bands found within our approach have non-trivial spin polarization and an intriguing spin-momentum locking of potential interest for novel nano-devices and spintronics [89, 90].

The physics becomes even more rich upon the inclusion of electronic interactions. We focus the analysis on the twist angle *θ* = 6° motivated by the fact that the tight binding fit dictating the low-energy band structure has highest quality in this case because long range hopping parameters become less relevant. Furthermore, the flat-band bandwidth for *θ* = 6° is the smallest among the twist angles considered in this work which leads to substantially quenched kinetic energy scales and enhanced interaction effects. The system exhibits two separate van Hove singularities which trigger a series of unconventional weak-coupling instabilities. We scrutinized these instabilities using unbiased renormalization group enhanced diagrammatic techniques, which point to extended regions where DWs or superconductivity emerge. Since spin and charge are entangled in non-*SU*(2) symmetric systems without inversion symmetry, the phase diagram and the classification of the phases of matter expected becomes extremely intricate. Our analysis shows that of the two types of superconductivity present in the phase diagram, one is topologically non-trivial and the other is trivial. The topologically trivial SC phase, which occurs only at low filling fractions of the flat bands, is still interesting for its high-angular-momentum form factors, being of the *g*- and *f*-type. The topologically non-trivial SC phase occupies a large fraction of the phase diagram as a function of filling (around half-filling) and interaction strength, which suggests PtSe_{2} as an interesting material to search for topological superconductivity.

Our work highlights another exciting example enabled by flexible Moiré engineering, concentrating this time on the less explored tailoring of SOC. Engineering SOC is an important topic in the field of quantum materials as SOC can trigger many fascinating topological transitions which might find a materials based application in Moiré materials for the first time.

## 6.Methods

*Density functional theory*—in our characterization of the large unit-cell twisted bilayer PtSe_{2} material we used the Vienna *Ab initio* simulation package (VASP). VASP was employed to determine the ground state of the system within the DFT [91] with the basis chosen to be plane waves and energy cutoff of 400eV. The pseudo potentials are generated using the projector augmented wave method [92] and the exchange–correlation functions are treated within the Perdew, Burke, and Ernzerhof (PBE) approximation [93]. We calculate the equilibrium lattice parameters of PtSe_{2} in the bulk phase and found that the optB86b van der Waals (vdW) functionals [94] provide better agreement with the experimental values, within less than 2% errors [95]. The optB86b vdW functionals are then adopted for all calculations. For these very large unit-cell simulations a 1 **×** 1 **×** 1 momentum grid is used to characterize the ground state and the mechanical relaxation. We construct the supercell of the considered bilayer system by using the optimized lattice constants of a 1 **×** 1 unit cell. DFT is most conveniently set up using periodic boundary condition and therefore, along the *z*-direction an auxiliary vacuum region larger than 15 Å is added. This region is chosen large enough that artificial interaction between the periodic slabs can be neglected. Our calculations are fully relaxed (w.r.t. all the atoms), which is known to be important in other Moiré systems to avoid artificial effects stemming from unrelaxed structures [96–98]. The relaxation procedure ensures that the force on each atom converges to values smaller than 0.01eV Å^{−1}. For all calculations, due to the relativistic effect in heavy element Pt, the SOC effect is considered while results without SOC are also calculated for comparison to estimate the effect of the SOC on the Moiré flat bands. The twisting angles of 6°, 7.34° and 9.43° contain in total 546, 366 and 222 number of atoms to consider in the unit cell.

*Tight-binding parameters*—we perform fitting of the tight-binding Hamiltonian (equation(3)) to DFT band structures of the Moiré flat bands of twisted PtSe_{2}. By that, we obtain the ten parameters and *ϕ*_{1,2,3} for each of the three twist angles considered. We tabulate the parameters in table1.

*Functional renormalization group*—we treat the interacting two-band, non-*SU*(2) tight binding model on the triangular lattice using the fRG. This method smoothly interpolates the free action *S*^{Λ=∞} to the full, interacting action *S*^{Λ=0}. We employ a sharp frequency cutoff scheme in the fermionic propagator:

with . Numerical treatment is rendered possible by approximating the infinite hierarchy of flow equations [84, 99, 100] and discarding all vertices that describe more than four-fermion interactions as well as setting the four-point vertex Γ^{(4),Λ}(1, 2, 3, 4) constant for all incoming and outgoing frequencies. As we are interested in static properties, we further neglect frequency dependencies on the two-point vertex (self-energy) and thus set Γ^{(2),Λ}(1, 2) ≡ 0. By these approximations, we arrive at the following flow equations of the four-point vertex:

The channel-projections read

with the momenta transformed as

Note that the diagram contributions *P*^{Λ} and *D*^{Λ} are written in the respective channel-projected bases, whereas the relation from *P*^{Λ} to *D*^{Λ} (equation(10)) is given in the ordering where 1, 2 are ingoing and 3, 4 outgoing indices. The fermionic particle–particle and particle–hole loops have to be momentum-reordered as well with , and read

with *u*_{σb }(** k **) the Bloch functions of the non-interacting tight-binding Hamiltonian and

_{b }(

**) its dispersion.**

*k*We calculate the fRG flow on a regular 24 **×** 24 momentum mesh in the 2D primitive zone and take both the full spin- and momentum-structure of the four point vertex into account. The summation over ** k ** in equations(8) and (9) is carried out on a finer momentum mesh with 649 points per coarse momentum point. The fine points are constructed to equally space out the Wigner–Seitz cells defined by the regular, coarse momentum mesh. To integrate the differential equation for Γ

^{(4),Λ}, we employ an enhanced Euler scheme with adaptive step size chosen that the maximal step size can never be above 10% of the current Λ. We consider the vertex diverged if its absolute maximal entry reaches 20 times the system's bandwidth.

We analyze the instabilities in a two-fold procedure. First, we determine the divergent channel by inspecting whether particle–particle (*P*^{Λ}) or particle–hole (*D*^{Λ}) are the dominant contributions to make Γ^{(4),Λ} diverge.

In the particle–particle case, we further investigate the SC instability by solving a linearized gap equation for Γ^{(4),Λ}:

with the particle–particle (and particle–hole) fermi-loops given by

with the Fermi function . As the eigenproblem in equation(17) is non-Hermitian and thus numerically highly unstable, we instead solve for the singular values and vectors of the matrix composed of Γ^{P,Λ} and *L*^{f,P,Λ}:

where the right singular vectors are the Fermi surface projected gap functions and the left singular vectors lack the Fermi surface structure but instead show the gap's symmetry more clearly. We transform the gap function (i.e. leading singular vector) to singlet and triplet space using the following identity [85, 86]:

with the vector of Pauli matrices and the identity matrix .

The leading instability is doubly degenerate for a large part of the phase diagram. Thus, we compute the complex superposition of the two leading instabilities (cf equation(5)). For *ϑ* ∈ (0, *π*) and *φ* ∈ (0, *π*), we evaluate the free energy of the system in the SC state via

where *E*_{b }(** k **) is the dispersion of the Bogoljubov–de-Gennes Hamiltonian

and is the pseudoinverse of the SC vertex as a matrix in (** k **,

*σ*

_{1},

*σ*

_{2}) and (

**',**

*k**σ*

_{3},

*σ*

_{4}). We evaluate

*F*

^{Λ}for Λ = Λ

_{c}. Next, we find the angles

*ϑ*

_{0},

*φ*

_{0}at which the free energy is minimized. For the analysis of the topology in the SC state, we use the physically realized instability at

*ϑ*

_{0}and

*φ*

_{0}.

In the particle–hole case, we instead extract the spin/density susceptibility from the four-point vertex at the final scale given by

Subsequently, we transform the four-point susceptibility to physical channels [88]:

The density–density response is given by *χ*^{00}(** q **) and the spin response functions by

*χ*

^{xx }(

**),**

*q**χ*

^{yy }(

**),**

*q**χ*

^{zz }(

**),**

*q**χ*

^{xy }(

**),**

*q**χ*

^{xz }(

**) and**

*q**χ*

^{yz }(

**).**

*q*## Acknowledgments

We thank J Beyer and J Hauck for useful discussions on the generation and analysis of non-*SU*(2) fRG results. This work is supported by the European Research Council (ERC-2015-AdG-694097), Grupos Consolidados (IT1249-19), and SFB925. MC is supported by a startup Grant from the University of Pennsylvania. AR is supported by the Flatiron Institute, a division of the Simons Foundation. We acknowledge funding by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under RTG 1995 and RTG 2247, within the Priority Program SPP 2244 '2DMP', under Germany's Excellence Strategy—Cluster of Excellence Matter and Light for Quantum Computing (ML4Q) EXC 2004/1-390534769 and—Cluster of Excellence and Advanced Imaging of Matter (AIM) EXC 2056-390715994. LX acknowledges the support from Distinguished Junior Fellowship program by the South Bay Interdisciplinary Science Center in the Songshan Lake Materials Laboratory and the Key-Area Research and Development Program of Guangdong Province of China (Grant No. 2020B0101340001). We acknowledge computational resources provided by the Simons Foundation Flatiron Institute, the MaxPlanck Computing and Data Facility, RWTH Aachen University under Project Number rwth0716 and the Platform for Data-Driven Computational Materials Discovery of the Songshan Lake laboratory. This work was supported by the MaxPlanck-New York City Center for Nonequilibrium Quantum Phenomena.

## Data availability

The raw data sets used for the presented analysis within the current study are available from the corresponding authors on reasonable request.

## Code availability

The tailored developed codes used in this work can be provided from the corresponding author on reasonable request. *Ab initio* calculations are done with the code VASP (version 5.4.4).