Deconfinement of Mott Localized Electrons into Topological and Spin-Orbit Coupled Dirac Fermions
José M. Pizarro, Severino Adler, Karim Zantout, Thomas Mertz, Paolo Barone, Roser Valentí, Giorgio Sangiovanni, Tim O. Wehling
DDeconfinement of Mott Localized Electrons into Topological and Spin-Orbit CoupledDirac Fermions
Jos´e M. Pizarro,
1, 2, ∗ Severino Adler,
3, 4
Karim Zantout, Thomas Mertz, Paolo Barone, Roser Valent´ı, Giorgio Sangiovanni, and Tim O. Wehling
1, 2, † Institute for Theoretical Physics, University of Bremen, Otto-Hahn-Allee 1, 28359 Bremen, Germany Bremen Center for Computational Material Sciences,University of Bremen, Am Fallturm 1a, 28359 Bremen, Germany Institut f¨ur Theoretische Physik und Astrophysik and W¨urzburg-Dresden Cluster of Excellence ct.qmat,Universit¨at W¨urzburg, 97074 W¨urzburg, Germany Institute of Solid State Physics, TU Wien, A-1040 Vienna, Austria Institute of Theoretical Physics, Goethe University Frankfurt am Main, D-60438 Frankfurt am Main, Germany Consiglio Nazionale delle Ricerche, Institute for Superconducting and Innovative Materials and devices (CNR-SPIN),c/o University ”G. D’Annunzio” , Chieti 66100, Italy (Dated: December 15, 2020)The interplay of electronic correlations, spin-orbit coupling and topology holds promise for therealization of exotic states of quantum matter. Models of strongly interacting electrons on honey-comb lattices have revealed rich phase diagrams featuring unconventional quantum states includingchiral superconductivity and correlated quantum spin Hall insulators intertwining with complexmagnetic order. Material realizations of these electronic states are however scarce or inexistent. Inthis work, we propose and show that stacking 1T-TaSe into bilayers can deconfine electrons froma deep Mott insulating state in the monolayer to a system of correlated Dirac fermions subject tosizable spin-orbit coupling in the bilayer. 1T-TaSe develops a Star-of-David charge density wavepattern in each layer. When the Star-of-David centers belonging to two adyacent layers are stackedin a honeycomb pattern, the system realizes a generalized Kane-Mele-Hubbard model in a regimewhere Dirac semimetallic states are subject to significant Mott-Hubbard interactions and spin-orbitcoupling. At charge neutrality, the system is close to a quantum phase transition between a quan-tum spin Hall and an antiferromagnetic insulator. We identify a perpendicular electric field andthe twisting angle as two knobs to control topology and spin-orbit coupling in the system. Theircombination can drive it across hitherto unexplored grounds of correlated electron physics includinga quantum tricritical point and an exotic first-order topological phase transition. I. INTRODUCTION
Prospects of quantum information technologies havemotivated an intense search for systems which intertwinetopology and electronic correlations [1–4]. Strongly in-teracting and spin-orbit coupled electrons on the hon-eycomb lattice, as theoretically described by the Kane-Mele-Hubbard model [5], feature quantum spin Hall(QSH), Mott-Hubbard, collective magnetic and chiral su-perconducting states. While several weakly interactingQSH systems are now well-studied experimentally [2, 6],material realizations of honeycomb Kane-Mele-Hubbardfermions with strong correlations and spin-orbit couplingare rare [7, 8]. Here, we introduce a new ’van der Waalsengineering’ platform to serve this purpose.We show that stacking of 1T-TaSe into bilayers candeconfine electrons from a deep Mott insulating state re-alized in the monolayer to a system of correlated Diracfermions subject to sizable spin-orbit coupling. Centralto this transition is the possibility of van der Waals ma-terials to stack in different configurations. For a spe- ∗ Electronic address: [email protected] † Electronic address: [email protected] cific honeycomb arrangement (Fig. 1), the kinetic en-ergy associated with the electronic hopping t turns outto be of the same order of magnitude as the effectivelocal Coulomb repulsion U . The system features there-fore electronic correlations, which turn out to put thesystem right on the verge between QSH and correlatedantiferromagnetic insulating states at charge neutralityand support chiral superconductivity under doping. Wefinally demonstrate that tuning the system via electricfields and twisting of the layers relative to each othersensitively affects the low-energy electronic structure interms of emerging Dirac mass and spin-orbit couplingterms and leads to completely unexplored regimes of cor-related electrons. II. RESULTS AND DISCUSSIONA. From a correlated insulator to emergent Diracfermions
Layered group-V transition metal dichalcogenides(TMDCs) such as 1T-TaSe or 1T-TaS feature alow-temperature commensurate charge density wave(CCDW) where Ta-atoms are displaced into Star-of-David (SoD) patterns (Fig. 1a) [9, 10]. In this phase, a r X i v : . [ c ond - m a t . s t r- e l ] D ec y x topbottom xy t = m e V U ≈ m e V Non-twisted E bottomtop z y x bca Fig. 1:
Crystal structures of 1T-TaSe mono- and bilayers in the CCDW phase. a , Monolayer 1T-TaSe in theCCDW phase. Only Ta atoms are shown. The Ta atoms are distorted into a SoD pattern, where the central Ta atoms (red)are surrounded by two rings of in total twelve Ta atoms (black). The SoDs are marked with red lines as guide to the eyes. b ,Top and three-dimensional side view of honeycomb stacked CCDW 1T-TaSe bilayer. Only the Ta atoms in the SoD centersare shown, with blue spheres and shaded regions marking the bottom layer atoms, and red spheres and shaded regions for thetop layer. The inset illustrates the leading terms of the effective low-energy Hamiltonian, i.e. the nearest-neighbor hopping t = −
34 meV and the local interaction U ≈
130 meV. c , Side view of 1T-TaSe bilayers embedded in field effect transistorstructures for the application of vertical electric fields. A non-twisted and a 180 ○ twisted bilayer are shown with Ta atoms (redand blue) and Se atoms (yellow). the SoDs form a triangular √ × √
13 superlattice in ev-ery layer and host correlated electrons: 1T-TaS shows ametal-to-insulator transition when entering the CCDWphase [10]. In 1T-TaSe , the bulk remains conductive tilllowest temperatures, while the surface exhibits a Motttransition around 250 K [11, 12]. Recently, 1T-TaSe has been fabricated down to monolayer thickness [13–15]and a pronounced thickness dependence of the electronicstructure has been reported [15]. As bonds between thelayers are mainly of van der Waals type, different stackingconfigurations are observed in experiments [16, 17] andhave a strong impact on the electronic structure [18–20].We compare the CCDW state in the monolayer to abilayer with honeycomb stacking, where the SoD centersform a buckled honeycomb lattice (Fig. 1b). The bot-tom layer SoD centers form sublattice A and the top layersublattice B. This stacking is one of many possible con-figurations, which can generally differ, both, by the localatomic stacking and by the stacking of the SoD centers.The local stacking considered here, has the Ta atomsof the bottom layer approximately beneath the lower Seatoms of the top layer (Fig. 1c). While this kind oflocal stacking is not the most commonly observed oneof the T-phase TMDCs [15, 17], a corresponding stack-ing sequence has been found in transmission electron mi-croscopy studies of 1T-TaS [17] and turns out to bemetastable in density functional theory (DFT) simula-tions of 1T-TaSe [21]. Regarding the stacking of theSoD centers further van der Waals DFT total energycalculations (see Supplementary Information) show thatthe particular honeycomb arrangement studied here, ismetastable and energetically on the order of 10 meV performula unit above the lowest energy configuration. Thisis similar to the configuration observed experimentallyin Ref. [15] and it is thus plausible that also the hon-eycomb configuration considered here, is within reach ofexperiments. Experimental approaches including tear-and-stack [22] and STM voltage pulse based manipula-tion schemes [23] can present possible avenues to con-trol and switch metastable stacking configurations andto reach the honeycomb configurations discussed in ourpaper. At small twist angles all different kinds of con-figurations can be realized locally in the moir´e includinglikely those honeycomb cases discussed here.To study the electronic structure of such engineeredstackings, we combine ab-initio calculations in the frame-works of DFT and the random phase approximation(RPA) with effective low-energy models, which we in-vestigate with dynamical mean-field theory (DMFT)and two-particle self-consistent (TPSC) many-body ap-proaches, see Methods section.In the CCDW phase, the band structure of monolayer1T-TaSe obtained with non-spin-polarized DFT is char-acterized by a single (Ta) flat band at the Fermi level[10, 18–20], which has a bandwidth of less than 20 meV(Fig. 2a, left). Hence, the CCDW formation largelyquenches in-plane hopping of the electrons. In the hon-eycomb stacked bilayer (with no twist), two dispersive bands with a bandwidth of the order of 200 meV emergefrom the low-energy flat band of the monolayer (Fig. 2a,right). Comparison of the mono- and bilayer bandwidthsshows that interlayer hopping effects must dominate overintralayer hoppings by approximately an order of magni-tude. In this sense, CCDW TaSe bilayers are the exactopposite of graphene bilayer systems, since in the latterout-of-plane coupling is an order of magnitude weakerthan in-plane hopping [22, 23]. For the 1T-TaSe hon-eycomb bilayer, the upper and lower low-energy bandstouch as Dirac points at the Brillouin zone corners Kand K’. In the undoped system, these Dirac points areexactly at the Fermi level.We next construct a Wannier Hamiltonian to describethe Dirac bands with one Wannier function for each SoDcenter, i.e. two Wannier orbitals per bilayer CCDW su-perlattice unit cell (Supplemental Figure 2). The result-ing nearest-neighbor hopping between a sublattice A sitein the bottom and a neighboring sublattice B site in thetop layer amounts to t = −
34 meV and is the leading termof the Wannier Hamiltonian. There are further terms inthe Wannier Hamiltonian, which are, however, at leastan order of magnitude smaller than t .The effective Hubbard interaction U for the SoD Wan-nier orbitals of CCDW TaSe calculated in RPA is U ≈
130 meV, which is also in line with the experimental es-timates in Ref. [15] and calculations for TaS [19]. Theratio of hopping to Coulomb interaction is decisive indetermining the strength and kind of electronic corre-lation phenomena taking place. Our calculations yield U /∣ t ∣ ≈ . bilayer in the framework ofDMFT and the TPSC approach [24, 25]. The quasi-particle weight Z shown in Fig. 2c is a measure of theelectronic correlation strength. In the temperature range T = −
230 K, both DMFT and TPSC consistently yieldessentially constant Z ≈ .
75. Our system is thus at in-termediate local correlation strength and far away fromthe paramagnetic Mott-Hubbard transition taking placeabove U Mott / t ≈ .
2. The DMFT and TPSC studies ofthe full Hubbard model for the non-twisted CCDW 1T-TaSe bilayer are, thus, in line with studies of the ideal-ized Hubbard model on the honeycomb lattice [26].The TPSC calculations give insight to spin-fluctuations taking place in the system (SupplementalFigure 4). The inverse intra- and inter-sublattice termsof the static magnetic susceptibility at wave number q = as well as the antiferromagnetic correlation length ξ AFM are shown in Fig. 2d. We observe that antiferro-magnetic fluctuations with alternating spin orientationbetween the two sublattices (A, B) are dominant andstrongly enhanced at temperatures T ≲
100 K. Thesefluctuations indicate that the non-twisted CCDW 1T-TaSe bilayer is close to a quantum phase transition froma Dirac semimetal to an antiferromagnetic insulator,which occurs for ideal Hubbard honeycomb systems S e m e n o ff m a ss Hubbard U quantum spin Hallphasetrivial band insulator e l e c t r i c f i e l d U c bandtouchingpoint TaSe bilayer antiferromagneticinsulator Δ M
60 100 150 200
T [K] Z TPSCDMFT -0.08-0.06-0.04-0.0200.020.04 / χ s p ( q = , ω = ) [ e V ] ξ a f m T [K] AB ξ afm AA K202
M = 0.00 meV K M = 0.45 meV K M = 0.82 meV K M = 1.19 meV AB S u b l a tt i c e c h a r a c t e r EE F [ m e V ] a bc d e Fig. 2:
Electronic structure and phase diagram of CCDW 1T-TaSe systems. a , Band structures of CCDW 1T-TaSe monolayer (left) in comparison to honeycomb stacked non-twisted CCDW 1T-TaSe bilayer (right). Bottom panels show azoom around the Fermi level E F . The red and dashed blue lines mark the DFT low-energy bands without and with SOCincluded, respectively. From the flat band near E F in the monolayer (red solid line, bandwidth ≈
14 meV) two dispersive bandswith a bandwidth ≈
200 meV emerge in the bilayer case. The bilayer bands exhibit Dirac points in the Brillouin zone corners,K and K’. b , Influence of extrinsic Semenoff mass terms ∆ M on the low-energy band structure. The sublattice character iscolor coded. The system changes from QSH to trivial band insulator at ∆ M = .
82 meV, which corresponds to a verticalelectric field of E z ≈ . − . c , Quasi-particle weight Z for non-twisted CCDW 1T-TaSe bilayer calculated with DMFTand TPSC. Both approaches place the system consistently in the moderately correlated regime Z ≈ .
75 at all calculatedtemperatures. d , Temperature-dependent antiferromagnetic correlation length ξ AFM and inverse static spin susceptibilities ofnon-twisted CCDW 1T-TaSe bilayer at wave vector q = as calculated with TPSC. The intra-sublattice (1 / χ AA = / χ BB ) andinter-sublattice (1 / χ AB ) elements of the inverse susceptibility at wave number q = are shown. e , Schematic phase diagram ofhoneycomb stacked non-twisted CCDW 1T-TaSe bilayer as a function of extrinsic Semenoff mass ∆ M and interaction strength U . The region accessible for non-twisted CCDW 1T-TaSe bilayer through tuning with external electric fields is highlighted.The transition from QSH to band insulator is a continuous transition at small U (dashed line) and a first order transition atlarger U (solid line). The red area in the quantum spin Hall region indicates the increasing many-body character of this phase. exactly in the range of U c / t ≈ . − . B. Spin-orbit coupling
The aforementioned magnetic correlation phenomenaare sensitive to details of the low-energy electronic struc-ture and spin-orbit coupling (SOC), which we discussin the following based on a symmetry analysis. Thespace group of non-twisted CCDW 1T-TaSe bilayer inthe honeycomb structure is P ¯3 ( C around anaxis perpendicular to the bilayer. Imposing also time-reversal symmetry, every band must be two-fold degen-erate. Therefore, SOC-induced qualitative changes of theband structures can occur near the Dirac points at K andK’. A corresponding k ⋅ p expansion (see SupplementaryInformation) reads H = ̵ hv F ( τ k x S x + k y S y ) + λ SOC τ σ z S z + α R2 ( k x σ y − k y σ x ) S z , (1) where the pseudospin S describes the sublattice degreeof freedom, σ acts on the electron spin, and τ = ± H : avalley-spin-sublattice coupling λ SOC = .
74 meV, whichis often called Kane-Mele spin-orbit term [29, 30], anda sublattice-staggered Rashba term α R2 , which belongsto the R2 class according to the classification form Ref.[31]. A finite Kane-Mele term λ SOC opens a gap andturns the system described by H into a QSH insulator[29, 30]. Importantly, the Kane-Mele term here is en-hanced in comparison to its counterpart in graphene bytwo orders of magnitude [22, 23, 32], and corresponds toa temperature T ≈
10 K, which is well accessible in exper-iments. Given that U /∣ t ∣ ≈ .
8, the non-twisted CCDW1T-TaSe bilayer implements a material realization of theKane-Mele-Hubbard model in a regime very close to thetopological quantum phase transition from a QSH insu-lator to an antiferromagnetic insulator (Fig. 2e).Dirac fermions and their topology are affected by dif-ferent kinds of mass fields. An energy difference betweenelectrons localized in sublattice A and B leads to a so-called Semenoff mass term M , which would enter theHamiltonian H of equation (1) in the form M S z . Thisterm breaks sublattice invariance and thereby inversionsymmetry and leads to a transition from a QSH to aband insulator at ∣ M ∣ = ∣ λ SOC ∣ [29, 33]. In the non-twisted CCDW 1T-TaSe bilayer, the intrinsic Semenoffmass M = M = .
14 meV, which issmall as compared to all other relevant terms in the sys-tem. Possible origins of this small symmetry breakingcan be the Wannier constructions and also faint asym-metries accumulated during self-consistency iterations inthe DFT calculations. Vertical electric fields E z , as re-alizable in field effect transistor geometries (Fig. 1c),break inversion symmetry, translate into staggered sub-lattice potentials, and therefore corresponding extrinsicSemenoff mass contributions ∆ M (Fig. 2b). DFT calcu-lations (Supplementary Figure 3) yield the approximaterelation ∆ M ≈ eE z d / (cid:15) ⊥ , where d = . e is the elementary charge, and (cid:15) ⊥ = .
72 playsthe role of an effective dielectric constant. The QSH toband insulator transition is reached for E z ≈ . − =
50 kV cm − (Fig. 2b), which is well within reach of ex-periments [34].Taken together, our calculations show that thehoneycomb-stacked non-twisted CCDW 1T-TaSe bi-layer is located in a region of the phase diagram (Fig. 2e)with three different phases (QSH insulator, band insu-lator and antiferromagnetic insulator) coming together.Electron correlation is known to change the order of theQSH-band insulator transition from second to first order[35]. Contrary to the standard non-interacting QSH toband insulator transition, where the gap closes and re-opens continuously with vanishing gap at the transitionpoint, the QSH gap remains finite and the system changesdiscontinuously to a band insulating state at the transi-tion point. Application of vertical electric fields in thesystem at hand represents hence a possibility to realizethis exotic interaction-induced first order transition. C. Twisted CCDW 1T-TaSe bilayers Layered van der Waals systems allow to realize dif-ferent stacking configurations via twisting, i.e. relativerotations between the layers. General twist angles θ leadto incommensurate moir´e patterns superimposed to theCCDW lattice. Arguably the simplest case of twisting isa rotation angle of θ = ○ , which leads to a system withidentical Bravais lattice but different symmetry of the su-percell basis (Fig. 1c): in case of honeycomb stacking ofthe CCDW, the space group of the 180 ○ twisted struc- ture of CCDW TaSe bilayer is reduced to P3, meaningthat inversion symmetry is lost with respect to the non-twisted case. The resulting band structure (Fig. 3a) isqualitatively similar to the non-twisted honeycomb caseregarding the overall shape and width of the low-energybands. Thus, also the 180 ○ twisted case will be far awayfrom the paramagnetic Mott transition. However, thelow-energy band structure is markedly different in the180 ○ twisted case. First, the conduction band is almostflat between K and M . Second, inversion symmetrybreaking lifts band degeneracies: our DFT calculationsreveal a staggered potential and an associated intrinsicSemenoff mass term of M = .
55 meV, which opens agap at the K and K’ points already without SOC andwithout external electric fields. Furthermore, additionalSOC terms are now allowed by symmetry (see Supple-mentary Information) and completely lift the remainingband degeneracies except for the time-reversal symmetricpoints Γ and M.Based on a symmetry analysis, we obtain the followinglow-energy model in the vicinity of K and K’: H ○ = H + M S z + Bτ σ z + λ R ( τ σ y S x − σ x S y ) + α R1 ( k x σ y − k y σ x )+ λ D ( τ S y k x − S x k y ) σ z . (2)The additional terms with respect to the non-twistedhoneycomb structure equation (1) are the intrinsic Se-menoff mass ( M ), the spin-valley coupling term ( B ),the Kane-Mele-Rashba interaction ( λ R ), and a Rashbainteraction belonging to the R1 class [31] giving rise topure Rashba spin-polarization patterns ( α R1 ). The lastterm (with coupling constant λ D ) can also be seen asan effective k − dependent magnetic field parallel to the z axis. Our DFT calculations yield M = .
55 meV, B = − .
85 meV, ∣ λ SOC ∣ ≲ .
05 meV and λ R = .
21 meV.These terms affect the dispersion and imprint an intricatesublattice and spin structure to the low-energy bands,which can be manipulated by vertical electric fields asshown in Fig. 3b.Except for E z = − . − (∆ M = − . M = M + ∆ M = . Z topological invariant for the non-interacting 180 ○ twisted bilayer in comparison to the non-twisted bilayer case as well as for two cases in betweenwhere the ratio B / λ SOC is varied (Fig. 3c, and Supple-mental Information). In the non-twisted bilayer, the sys-tem is in a QSH state unless an extrinsic sufficiently largeSemenoff mass term or an additional Rashba SOC term λ R are added. In the 180 ○ twisted case, the situationis very different regardless whether or not the intrinsicSemenoff mass term is compensated by an external elec-tric field and regardless of λ SOC . Indeed, many changesin the SOC terms suppress the QSH state in the 180 ○ twisted bilayer: The comparably large Rashba λ R andthe spin-valley coupling terms B and a strong reductionin λ SOC . Each of these alone is sufficient to suppress theQSH state. At vertical electric field E z = − . − -1.49 -0.74 0.00 0.74 1.49 M [meV] R [ m e V ] TaSe SOC =0.74 meV -0.74 0.00 0.74 1.49
M [meV]B/
SOC =0.9 -0.74 0.00 0.74 1.49
M [meV]B/
SOC =1.0 touchingbandsQSHtrivial -7.42 -6.68 -5.94 -5.19
M [meV]
TaSe E z = . m V / Å M = -5.34meV M = -6.68meV M = -8.02meV AB S ub l a tt i c e a E z bc EE F [ m e V ] -0.2-0.100.10.2 K M E - E F [ e V ] no SOCSOC EE F [ m e V ] Γ Γ
K K K R [ m e V ] Fig. 3:
Influence of 180 ○ twisting and electric fields on band structure and topology. a , Band structure for 180 ○ twisted CCDW 1T-TaSe bilayer without (solid red) and with (dashed blue) SOC included. The degeneracy of the bands islifted due to the absence of the inversion symmetry after twisting. b , Influence of vertical electric fields on the low-energy bandstructure for 180 ○ twisted CCDW 1T-TaSe bilayer. The gap at K closes and a band touching point occurs when the Semenoffmass term is ∣ M ∣ = ∣ M + ∆ M ∣ = .
87 meV at E z = − . − . c , Non-interacting λ R vs ∆ M topological phase diagrams fornon-twisted, 180 ○ twisted CCDW bilayer and two cases in between, where the spin valley coupling B of the non-twisted CCDW1T-TaSe bilayer is varied. For B = B increases,the QSH region shrinks and disappears when B = λ SOC . For B > λ SOC , the system behaves as a trivial band insulator. Theyellow arrows in the leftmost and rightmost panels indicate the path in phase space accessible by varying the vertical electricfield E z . the gaps at K and K’ close, and a band touching pointemerges.While it is clear that there will be tendencies towardsinteraction-induced (quasi)ordered phases as well, here,the kind of ordering is likely different from the non-twisted case but largely unexplored. The band touch-ing at E z = − . − implements a situation similarto saddle points in a two-dimensional dispersion, wherealready arbitrarily weak interactions would trigger differ-ent kinds of magnetic or excitonic instabilities [36]. Howthese instabilities translate into the intermediately corre-lated and strongly spin-orbit coupled case of 180 ○ twistedTaSe is a completely open question. D. Conclusions and outlook
The field of twistronics with materials like bilayergraphene is based on the idea that weak interlayer cou-pling can flatten highly dispersive bands and therebyboost electronic correlations [37–39]. The system intro-duced here takes the opposite route of deconfining for-merly Mott localized electrons. This approach shouldbe generally applicable to interfaces of Mott localizedelectrons under two conditions: the interlayer couplingshould substantially exceed the in-plane one and at the same time define a connected graph linking all sites ofthe system. Possible example systems range from stack-ing faults in the bulk of CCDW layered Mott materials[17] to molecular systems [40].Especially the twisting degree of freedom opens new di-rections to experiments. Since interlayer hopping is thedominant kinetic term in deconfined Mott systems likebilayers of CCDW 1T-TaSe , we expect incommensura-bility effects to be much more pronounced than in twistedgraphene systems [37]. θ = ○ twisted CCDW 1T-TaSe bilayer should realize a quasicrystal with twelvefold ro-tation symmetry and provide an experimental route tocorrelated electrons and emerging collective states in aquasicrystalline environment.The prototypical case of CCDW 1T-TaSe bilayerdemonstrates how deconfinement of Mott electrons leadsto exotic states of quantum matter: the non-twisted bi-layer approaches a quantum tricritical region of compet-ing QSH, trivial band insulating, and antiferromagneticinsulating states. At 180 ○ twist angle different kinds ofelectrically controllable band degeneracies with associ-ated many-body instabilities, hypothetically of excitonictype, emerge. Clearly, the phase space for manipulatingdeconfined Mott electrons is high dimensional. We hereidentified the combination of twist angle and perpendic-ular electric field as decisive for TaSe bilayers. Furthermeans to control emerging electronic states include di-electric engineering [41] and charge doping. Our calcu-lations showed that the non-twisted CCDW 1T-TaSe bilayer in honeycomb stacking approximates the (Kane-Mele) Hubbard model on the honeycomb lattice with U /∣ t ∣ ≈ . d + id -type. Methods
DFT calculations.
We perform DFT [46, 47] calcu-lations by using the Vienna ab-initio simulation package(VASP) [48, 49] with the generalized gradient approx-imation of Perdew, Burke, and Ernzerhof (GGA-PBE)for the exchange-correlation functional [50, 51]. We ob-tain the total energies, relaxed structures and electronicstructure of mono- and bilayer 1T-TaSe systems. Wecalculate the total energies for various possible stack-ings in undistorted and CCDW bilayer 1T-TaSe usingΓ-centered k -meshes of 15 × × × ×
1, respectivelyand taking into account van der Waals (vdW) correc-tions within DFT-D2 and cross-checking with DFT-D3[52, 53], see Supplemental Figure 1. The ionic relaxationis done using the conjugate gradient algorithm until allforce components are smaller than 0 .
02 eV ˚A − .Since the DFT-D2 corrections yield correct interlayerdistances but do not correctly capture the CCDW dis-tortions, we adopted the following relaxation procedureto calculate the commensurate √ ×√
13 CCDW bilayerstructures:1. We perform relaxations for a √ × √
13 supercellof the monolayer (without vdW corrections), usinga superlattice constant of a = .
63 ˚A accordingto Ref. [15]. We fix the vertical positions of theTa atoms, while allowing for Ta in-plane displace-ments. The Se atoms are allowed to freely relax inall three directions.2. We include then a second layer and optimize theinterlayer distance, d , while keeping all relative in-tralayer distances fixed. We find d = . .3. We relax the CCDW bilayer following the sameprocedure described for the monolayer, i.e. with-out vdW corrections and vertical positions of theTa atoms according to the optimized interlayer dis-tance d fixed, while allowing for in-plane displace-ments. The Se atoms are allowed to freely relax inall three directions. We cross check the results obtained by this procedureagainst calculations with vdW corrections according tothe DFT-D3 method. [53]. In contrast to DFT-D2, theCCDW distortions are well described in the DFT-D3.Thus, full relaxations of all atomic positions have beenperformed in the DFT-D3 framework for the mono- andbilayer. Both our step-by-step procedure using DFT-D2 method, and the full relaxation using DFT-D3 giveequivalent results for the total energies (see SupplementalFigure 1b), for the crystal and band structures.For the non-collinear magnetic calculation, i.e. whenSOC is included, we set the net magnetic moment to zeroin all atoms of the unit cell, and use a Γ-centered k -meshof 6 × × Estimation of the screened Hubbard interactionU via RPA.
We estimate the local Hubbard interac-tion U for the flat bands around the Fermi level in theCCDW 1T-TaSe bilayer from the ab-initio calculationof the screened Coulomb interaction, using RPA for theundistorted bilayer. We follow a similar procedure as inRef. [54], which we summarize below: • We initially calculate the
Wannier90 tight-binding model for the three low-energy Ta bands C ,whose orbital character is mostly { d z , d x − y , d xy } (see Fig. 2(a) in Ref. [54]) in the undistorted mono-layer 1T-TaSe . • The static RPA-screened Coulomb interaction ten-sor W αβγδ ( q , ω → ) is calculated for undistortedmonolayer 1T-TaSe , where q is a reciprocal wavevector on a Γ-centered mesh of 18 × ×
1, and α, β, γ, δ ∈ C . We neglect q = terms in our RPAanalysis in order to avoid unphysical effects. • In the CCDW 1T-TaSe bilayer, the d z orbitalsfrom Ta atoms in the SoD centers have the largestcontribution for the bands around the Fermi level.Thus, for each q , we consider only the tensor ele-ment W ( q ) ≡ W αααα ( q ) with α = d z . • Then, the local Hubbard interaction U in a singlestar of David is calculated by averaging over the d z orbital weight from each Ta atom (labeled by w d z ( R ) ) in the star of David: U = ∑ R , R ′ ∈ (cid:67) w d z ( R ) U ( R − R ′ ) w d z ( R ′ ) (3)where U ( R ) is the discrete Fourier transform of W ( q ) . DMFT and TPSC many-body calculations.
Fornon-twisted CCDW 1T-TaSe bilayer, we constructtight-binding-Hubbard Hamiltonians of the type H TBH = ∑ ⟨ i,j ⟩ ,σ,σ ′ t σσ ′ ij c † iσ c jσ ′ + U ∑ i n i ↑ n i ↓ (4)from a Wannierization of the ab-initio DFT band struc-ture (Supplementary Information), and estimate the on-site repulsion U to be about 130 meV by means of RPAcalculations. We study these effective low-energy mod-els from a many-body perspective to judge the typeof correlations in the system. In DMFT the latticeHamiltonian is mapped onto a self-consistently deter-mined single impurity problem, solved – in our case –within numerical exact quantum Monte Carlo in the hy-bridization expansion flavour (CT-HYB) (for a review,see [55]). The resulting sublattice-resolved self-energy,Σ, is local ( k -independent) but it contains frequency-dependent non-perturbative corrections beyond Hartree-Fock to all orders and can account for Mott-Hubbardmetal-to-insulator transitions. All DMFT calculationsare performed using w2dynamics [56]. The double count-ing is accounted for using the fully-localized limit. Fortwo-dimensional systems it is important to estimate non-local effects at the level of the self-energy, not includedin DMFT. To this goal, we apply the TPSC method[57], which produces accurate results in the weak-to-intermediate coupling regime, if compared to latticequantum Monte Carlo calculations in the single bandHubbard model. For non-twisted CCDW 1T-TaSe bi-layer, which is modelled by a multi-band system we usethe multi-site formulation of TPSC [24] while neglectingthe Hartree term to avoid double counting of correlationeffects already accounted for in DFT. Moreover, to beable to apply TPSC to this system we project out spinoff-diagonal terms and take only the diagonal spin-upcontributions from DFT while still assuming a paramag-netic state. The combination of TPSC accounting for the k -dependence of Σ and DMFT, in which we can includeall off-diagonal terms between spin-orbitals and accessantiferromagnetic ordering at strong coupling consitutesa powerful tool to determine the many-body nature of1T-TaSe . Data Availability
The data that support the findings of this study isavailable from the corresponding author upon reasonablerequest.
Acknowledgements
We thank D. Di Sante, P. Eck, E. van Loon, M.Sch¨uler, and C. Steinke for useful conversations. JMPand TW acknowledge funding from DFG-RTG 2247(QM ) and the European Graphene Flagship. SA andGS are supported by DFG-SFB 1170 Tocotronics, andfurther acknowledges financial support from the DFGthrough the W¨urzburg-Dresden Cluster of Excellenceon Complexity and Topology in Quantum Matter – ct.qmat Competing Interests
The authors declare no competing interests.
Author Contributions
J.M.P. performed the DFT and RPA calculations. S.A.calculated the topological invariant diagrams. S.A., K.Z.and T.M. performed the DMFT and TPSC calculations.P.B. derived the k ⋅ p model. R.V., G.S. and T.W. super-vised the project. J.M.P., S.A., G.S. and T.W. analyzedthe results. All authors contributed to write the paper. [1] Tsui, D. C., Stormer, H. L. & Gossard, A. C. Two-dimensional magnetotransport in the extreme quantumlimit. Phys. Rev. Lett. , 1559–1562 (1982).[2] Qi, X.-L. & Zhang, S.-C. Topological insulators and su-perconductors. Rev. Mod. Phys. , 1057–1110 (2011).[3] Sato, M. & Ando, Y. Topological superconductors: areview. Rep. Prog. Phys. , 076501 (2017).[4] Rachel, S. Interacting topological insulators: a review. Rep. Prog. Phys. , 116501 (2018).[5] Hohenadler, M. & Assaad, F. F. Correlation effects intwo-dimensional topological insulators. J. Phys.: Con-dens. Matter , 143201 (2013).[6] Wehling, T. O., Black-Schaffer, A. M. & Balatsky, A. V.Dirac materials. Adv. in Phys. , 1–76 (2014).[7] Marrazzo, A., Gibertini, M., Campi, D., Mounet, N. &Marzari, N. Prediction of a large-gap and switchableKane-Mele Quantum Spin Hall insulator. Phys. Rev.Lett. , 117701 (2018).[8] Wu, X., Fink, M., Hanke, W., Thomale, R. & Di Sante,D. Unconventional superconductivity in a doped quan-tum spin Hall insulator.
Phys. Rev. B , 041117 (2019).[9] Wilson, J. A., Salvo, F. J. D. & Mahajan, S. Charge-density waves and superlattices in the metallic layeredtransition metal dichalcogenides.
Adv. in Phys. , 117–201 (1975).[10] Rossnagel, K. On the origin of charge-density waves inselect layered transition-metal dichalcogenides. J. Phys.:Condens. Matter , 213001 (2011).[11] Perfetti, L. et al. Spectroscopic signatures of abandwidth-controlled Mott transition at the surface of1T-TaSe . Phys. Rev. Lett. , 166401 (2003).[12] Colonna, S. et al. Mott phase at the surface of 1T-TaSe observed by Scanning Tunneling Microscopy. Phys. Rev.Lett. , 036405 (2005).[13] Nakata, Y. et al. Selective fabrication of Mott-insulatingand metallic monolayer TaSe . ACS Appl. Nano Mater. , 1456–1460 (2018).[14] B¨orner, P. C. et al. Observation of charge densitywaves in free-standing 1T-TaSe monolayers by transmis-sion electron microscopy. Appl. Phys. Lett. , 173103(2018). [15] Chen, Y. et al. Strong correlations and orbital texture insingle-layer 1T-TaSe . Nature Phys. , 218–224 (2020).[16] M¨uller-Caspary, K. et al. Atomic-scale quantification ofcharge densities in two-dimensional materials. Phys. Rev.B , 121408 (2018).[17] Hovden, R. et al. Atomic lattice disorder in charge-density-wave phases of exfoliated dichalcogenides (1T-TaS ). PNAS , 11420–11424 (2016).[18] Freericks, J. K., Krishnamurthy, H. R., Ge, Y., Liu,A. Y. & Pruschke, T. Theoretical description of time-resolved pump/probe photoemission in TaS : a single-band DFT+DMFT(NRG) study within the quasiequilib-rium approximation. physica status solidi (b) , 948–954 (2009).[19] Darancet, P., Millis, A. J. & Marianetti, C. A. Three-dimensional metallic and two-dimensional insulating be-havior in octahedral tantalum dichalcogenides. Phys.Rev. B , 045134 (2014).[20] Ritschel, T., Berger, H. & Geck, J. Stacking-driven gapformation in layered 1T-TaS . Phys. Rev. B , 195134(2018).[21] Lazar, P., Martincov´a, J. & Otyepka, M. Structure, dy-namical stability, and electronic properties of phases inTaS from a high-level quantum mechanical calculation. Phys. Rev. B , 224104 (2015).[22] Weston, A. et al. Atomic reconstruction in twisted bilay-ers of transition metal dichalcogenides. Nat. Nanotech-nol. , 592–597 (2020).[23] Ma, L. et al. A metallic mosaic phase and the origin ofMott-insulating state in 1T-TaS 2. Nat. Comm. , 10956(2016).[24] Amaricci, A., Budich, J., Capone, M., Trauzettel, B. &Sangiovanni, G. First-order character and observable sig-natures of topological quantum phase transitions. Phys.Rev. Lett. , 185701 (2015).[25] Castro Neto, A. H., Guinea, F., Peres, N. M. R.,Novoselov, K. S. & Geim, A. K. The electronic prop-erties of graphene.
Rev. Mod. Phys. , 109–162 (2009).[26] Katsnelson, M. I. Graphene: Carbon in Two Dimensions ,(Cambridge University Press, 2012).[27] Zantout, K., Altmeyer, M., Backes, S. & Valent´ı, R.Superconductivity in correlated BEDT-TTF molecularconductors: Critical temperatures and gap symmetries.
Phys. Rev. B , 014530 (2018).[28] Zantout, K., Backes, S. & Valent´ı, R. Effect of nonlocalcorrelations on the electronic structure of LiFeAs. Phys.Rev. Lett. , 256401 (2019).[29] Tran, M.-T. & Kuroki, K. Finite-temperature semimetal-insulator transition on the honeycomb lattice.
Phys. Rev.B , 125125 (2009).[30] Arya, S., Sriluckshmy, P. V., Hassan, S. R. & Trem-blay, A.-M. S. Antiferromagnetism in the Hubbard modelon the honeycomb lattice: A two-particle self-consistentstudy. Phys. Rev. B , 045111 (2015).[31] Raczkowski, M. et al. The Hubbard model on the honey-comb lattice: from static and dynamical mean-field the-ories to lattice quantum Monte Carlo simulations. Phys.Rev. B , 125103 (2020).[32] Kane, C. L. & Mele, E. J. Z Topological order and theQuantum Spin Hall effect.
Phys. Rev. Lett. , 146802(2005).[33] Ezawa, M. Valley-polarized metals and QuantumAnomalous Hall effect in silicene. Phys. Rev. Lett. ,055502 (2012). [34] Zhang, X., Liu, Q., Luo, J.-W., Freeman, A. J. & Zunger,A. Hidden spin polarization in inversion-symmetric bulkcrystals.
Nature Phys. , 387–393 (2014).[35] Kochan, D., Irmer, S. & Fabian, J. Model spin-orbitHamiltonians for graphene systems. Phys. Rev. B ,165415 (2017).[36] Di Sante, D. et al. Towards topological quasifreestand-ing stanene via substrate engineering. Phys. Rev. B ,035145 (2019).[37] Klein, J. et al. Electric-field switchable second-harmonicgeneration in bilayer MoS by inversion symmetry break-ing. Nano Lett. , 392–398 (2017).[38] Kotov, V. N., Uchoa, B., Pereira, V. M., Guinea, F.& Castro Neto, A. H. Electron-electron interactions ingraphene: Current status and perspectives. Rev. Mod.Phys. , 1067–1125 (2012).[39] Bistritzer, R. & MacDonald, A. H. Moir´e bands intwisted double-layer graphene. Proc Natl Acad Sci USA , 12233 (2011).[40] Cao, Y. et al. Unconventional superconductivity inmagic-angle graphene superlattices.
Nature , 43–50(2018).[41] Cao, Y. et al. Correlated insulator behaviour at half-filling in magic-angle graphene superlattices.
Nature ,80–84 (2018).[42] Tsukahara, N. et al. Evolution of Kondo resonance froma single impurity molecule to the two-dimensional lattice.
Phys. Rev. Lett. , 187201 (2011).[43] Pizarro, J. M., R¨osner, M., Thomale, R., Valent´ı, R. &Wehling, T. O. Internal screening and dielectric engi-neering in magic-angle twisted bilayer graphene.
Phys.Rev. B , 161102 (2019).[44] Black-Schaffer, A. M. & Doniach, S. Resonating va-lence bonds and mean-field d -wave superconductivity ingraphite. Phys. Rev. B , 134512 (2007).[45] Nandkishore, R., Levitov, L. S. & Chubukov, A. V.Chiral superconductivity from repulsive interactions indoped graphene. Nature Phys , 158–163 (2012).[46] Kiesel, M. L., Platt, C., Hanke, W., Abanin, D. A. &Thomale, R. Competing many-body instabilities and un-conventional superconductivity in graphene. Phys. Rev.B , 020507 (2012).[47] Black-Schaffer, A. M. & Honerkamp, C. Chiral d -wavesuperconductivity in doped graphene. J. Phys.: Con-dens. Matter , 423201 (2014).[48] Hohenberg, P. & Kohn, W. Inhomogeneous electron gas. Phys. Rev. , B864–B871 (1964).[49] Kohn, W. & Sham, L. J. Self-consistent equations in-cluding exchange and correlation effects.
Phys. Rev. ,A1133–A1138 (1965).[50] Kresse, G. & Hafner, J. Norm-conserving and ultrasoftpseudopotentials for first-row and transition elements.
J.Phys.: Condens. Matter , 8245–8257 (1994).[51] Kresse, G. & Joubert, D. From ultrasoft pseudopotentialsto the projector augmented-wave method. Phys. Rev. B , 1758–1775 (1999).[52] Perdew, J. P., Burke, K. & Ernzerhof, M. GeneralizedGradient Approximation made simple. Phys. Rev. Lett. , 3865–3868 (1996).[53] Perdew, J. P., Burke, K. & Ernzerhof, M. GeneralizedGradient Approximation made simple [Phys. Rev. Lett.77, 3865 (1996)]. Phys. Rev. Lett. , 1396–1396 (1997).[54] Grimme, S. Semiempirical GGA-type density func-tional constructed with a long-range dispersion correc- tion. Journal of Computational Chemistry , 1787–1799(2006).[55] Grimme, S., Antony, J., Ehrlich, S. & Krieg, H. A con-sistent and accurate ab initio parametrization of densityfunctional dispersion correction (DFT-D) for the 94 el-ements H-Pu. The Journal of Chemical Physics ,154104 (2010).[56] Kamil, E. et al. Electronic structure of single layer1T-NbSe : interplay of lattice distortions, non-local ex-change, and Mott–Hubbard correlations. J. Phys.: Con-dens. Matter , 325601 (2018). [57] Gull, E. et al. Continuous-time Monte Carlo methods forquantum impurity models. Rev. Mod. Phys. , 349–404(2011).[58] Wallerberger, M. et al. w2dynamics: Local one- andtwo-particle quantities from dynamical mean field the-ory. Computer Physics Communications , 388–399(2019).[59] Vilk, Y. M. & Tremblay, A.-M. Non-perturbative many-body approach to the Hubbard model and single-particlepseudogap.
Journal de Physique I , 1309–1368 (1997). Supplementary Material
S1. STACKING POTENTIAL ENERGY LANDSCAPES
Here we study the dependence of the total energies of 1T-TaSe bilayers on the local atomic stacking (SupplementaryFigure S1a) and the stacking of the CCDW centers (Supplementary Figure S1b).Dependencies on the local atomic stacking can be inferred from the stacking potential energy landscapes derivedfrom our DFT calculations for undistorted bilayers of 1T-TaSe . The non-twisted and 180 ○ twisted cases are shown inSupplementary Figure S1a. For non-twisted bilayer, Ta from the top layer above Ta from the bottom layer (stackingI) is the most stable configuration. Top layer Ta above bottom layer Se (stacking III) results to be a local minimum.For 180 ○ twisted bilayer, the most stable configurations are given for Ta on top of Se (stackings III and V). Regardingthe local atomic arrangement, the honeycomb stacking discussed in the main text corresponds to configuration III.Energy differences between different stacking configurations of the CCDW centers in the distorted non-twisted1T-TaSe bilayers are given in Supplementary Figure S1b with the labeling of the configurations being defined inSupplementary Figure S1c: AA, AB and AC refers to CCDW centers of the top layer on top of a Ta atom fromthe bottom layer (local configuration I). Both, AtX and AbX refer to arrangements, where the Ta atom in theCCDW centers from one layer is on top of / beneath Se atom sites in the other layer (local configurations III and V,respectively). In this notation, At3 is the honeycomb stacking that is considered in the main text.It becomes clear that the undistorted bilayer (Supplementary Figure S1a) already gives a good overall estimate ofthe total stacking differences energies of the CCDW bilayer (Supplementary Figure S1b). Hence, also in the CCDWcase, the stacking potential energy landscape is dominated by contributions from the local atomic stacking.In the CCDW non-twisted bilayer, the total stacking energies using DFT-D2 and DFT-D3 calculations are similar,with a small deviation for the most unfavorable stackings AbX (local configuration V). The resulting crystal and bandstructures in both methods (not shown) are also equivalent. The energy differences between Ta on top of Ta (AA,AB, AC) and Ta on top of Se above the Ta plane (At1, At2, At3) are relatively small, i.e. on the order of 10 meVper f.u.. This means that the potential energy landscapes are very flat, which is in line with realizations of multiplestacking configurations in experiments.Assuming that the 180 ○ twisted bilayer is similar to non-twisted case in that the potential energy landscape isdominated by contributions from the local atomic stacking, we expect that stackings AtX and AbX are the moststable ones in the 180 ○ twisted case. This would support formation of the honeycomb pattern in the 180 ○ twistedcase, while again flat potential energy landscapes are expected. S2. WANNIER TIGHT-BINDING MODEL
The relevant subspace B for the low-energy bands of (distorted) CCDW 1T-TaSe monolayer and bilayer containsmainly d z orbitals from Ta atoms in the centers of the stars-of-David (SoD, see Fig. 1). We construct a correspondingminimal tight-binding model using Wannier90 code [58]. The Wannier basis is
B = { d topz , d bottomz } for each spin, wherethe superscripts refer to SoD center Ta atom from top and bottom layers. We use a Γ-centered k -mesh of 9 × × k -mesh of 6 × × bilayers. S3. k ⋅ pk ⋅ pk ⋅ p MODEL OF LOW ENERGY BANDS
We derive the effective model describing the low-energy band structure around the K and K’ points by using the k ⋅ p perturbation theory [59, 60] for the perturbed Hamiltonian H = H I + H II + H III , including Pauli-type spin-orbitinteraction, where: H I = ̵ hm k ⋅ p (S1) H II = ̵ h m c (∇ U × p ) ⋅ σ ≡ ̵ h m c A ⋅ σ (S2) H III = ̵ h m c (∇ U × k ) ⋅ σ ≡ ̵ h m c ( k × σ ) ⋅ ∇ U. (S3)2 I II III IV V VI E t o t [ m e V /f. u .] Non-twisted
IV VIVI II III
Se Ta
AAABAC Ab2Ab3At1 At2At3 Ab1 a b c
Fig. S1:
Stacking potential landscapes for 1T-TaSe bilayers. a, Potential energies for undistorted 1T-TaSe bilayer inthe non-twisted and 180 ○ twisted cases. The inset shows the nomenclature of the different stacking configurations I-VI withbottom layer Ta atoms in red and Se atoms above and below the bottom Ta plane in green and pink, respectively. The stackingof the top layer is marked with blue dots and latin numbers. Blue dots refer to the position of the Ta atom in the top layer,starting from the perfectly aligned bilayer in I. Non-twisted bilayer shows a total minimum for stacking I (perfect alignement),and a local metastable minimum for stacking III (Ta on top of the Se above the Ta plane). For 180 ○ twisted bilayer, stackings IIIand V are degenerate and are the most stable configurations. In the ideal honeycomb stacking, the CCDW 1T-TaSe bilayershave a local stacking with the Ta in the top layer (approximately) above Se in the bottom layer, i.e. configurations III and V. b, Total energies per formula unit and interlayer distances d versus stacking configurations in the CCDW non-twisted bilayer1T-TaSe as obtained within the DFT-D2 and DFT-D3 approaches. c, Nomenclature used for the different CCDW stackings,where blue, green and magenta dots describe the location of the Star-of-David (SoD) centers from the top layer. Thus, stackingAA corresponds to perfectly aligned CCDW bilayer. Stackings AtX stand for SoD central Ta atom from the top layer alignedto the Se atom above the Ta plane from the bottom layer, and AbX for SoD center from the top layer above the Se atom belowthe Ta plane from the bottom layer. DFT-D2 and DFT-D3 calculated energies are similar, with a discrepancy of approximately20 meV per f.u. for the AbX stackings. For Ta on top of Ta (AA, AB, AC) and Ta on top of Se above the lower Ta plane (At1,At2, At3) energy differences are ≲
10 meV per f.u.. Stacking At3 is the one considered in the main text. -0.2-0.100.10.2 E - E F [ e V ] Non-twisted -0.2-0.100.10.2 -0.2-0.100.1 Γ K M -0.2-0.100.1
K MK M
Γ Γ
Fig. S2:
Ab-initio
DFT band structure versus Wannierization for non-twisted and 180 ○ twisted CCDW 1T-TaSe bilayer . DFT bands (red solid lines) are compared with our Wannier tight-binding model (black dashed lines) for non-twisted(left panels) and 180 ○ twisted (right panels) bilayers. Top panels show the case without SOC included, the bottom panelswith SOC taken into account. Our Wannier tight-binding models fit the DFT band structures obtained for CCDW 1T-TaSe bilayers very well. The non-zero matrix elements of such perturbed Hamiltonian can be determined from group theory and exploiting thetransformation properties of the one-electron wave functions under the symmetry operations of the high-symmetrypoint K. The knowledge of the little group of K allows to define a set of linear equations [61]: ⟨ u µν ∣ O αβ ∣ u ij ⟩ = h ∑ R ∑ ν ′ β ′ j ′ µ D ∗ ν ′ ν ( R ) α D β ′ β ( R )× i D j ′ j ( R ) ⟨ u µν ′ ∣ O αβ ′ ∣ u ij ′ ⟩ (S4)3where h is the order of the group, u µν represents the ν -th component of the basis function of a representation µ and O αβ is an operator that transforms like the β -th component of the basis function of a representation α , while i D j ′ j ( R ) is the ( j ′ j ) element in the matrix representative of the group element R in the i -th representation. The non-twisted(180 ○ -twisted) CCDW bilayer belongs to the P ¯3 ( P
3) space group; both structures display a three-fold rotation C around an axis perpendicular to the bilayer, while the non-twisted stacked bilayer is also centrosymmetric, the twolayers being inversion partners. At point K only C symmetry is preserved, and non-relativistic bands can be groupedin a one-dimensional single-valued irreducible representation (IR) K and in two-fold degenerate K K IR for thespace group P ¯3; when the symmetry is lowered to P
3, bands belonging to K and K are not degenerate anymore.The single-valued IRs, the corresponding characters and the basis functions are listed in Table S1. By introducing ageneral operator π with components π = p z , π = √ ( p x + i p y ) , π = √ ( p x − i p y ) , (S5)each trasforming as K , K and K , respectively, one can use equation (S4) to identify its symmetry-allowed non-zeroexpectation values on the basis functions φ , φ spanning the K K IR: ⟨ φ ∣ π ∣ φ ⟩ , ⟨ φ ∣ π ∣ φ ⟩ , ⟨ φ ∣ π ∣ φ ⟩ , ⟨ φ ∣ π ∣ φ ⟩ . (S6)In this basis, therefore, the following k ⋅ p model is found to fulfill the point-group symmetries of the wave vectorat K: H = ( λ σ z + α ( k x σ y − k y σ x ) −̵ hv f k + + λ R ( iσ x − σ y ) + λ D k + σ z −̵ hv f k − − λ R ( iσ x + σ y ) + λ D k − σ z λ σ z + α ( k x σ y − k y σ x ) ) (S7)where k ± = k x ± ik y and we introduced the following parametrization: v f = m ⟨ φ ∣ p x − ip y ∣ φ ⟩ λ = ̵ h m c ⟨ φ ∣ A z ∣ φ ⟩ λ = ̵ h m c ⟨ φ ∣ A z ∣ φ ⟩ λ R = − i ̵ h m c ⟨ φ ∣ A x − iA y ∣ φ ⟩ λ D = − i ̵ h m c ⟨ φ ∣(∇ U ) x − i (∇ U ) y ∣ φ ⟩ α = ̵ h m c ⟨ φ ∣(∇ U ) z ∣ φ ⟩ α = ̵ h m c ⟨ φ ∣(∇ U ) z ∣ φ ⟩ . (S8)The effective model must also obey time-reversal symmetry Θ = ˆ T K , where K is the complex conjugation and ˆ T = iσ y P = S x , swapping the basis functions, also maps the Hamiltonianequation (S7) from point K to point K’. Here S represents the sublattice pseudospin, spanning the two-dimensionalspace defined by the basis functions φ , φ .Imposing both time-reversal and inversion invariance, one finds that λ R = λ D = λ = − λ ≡ − λ SOC α = − α ≡ α R (S9)recovering the effective k ⋅ p model of equation (1) in the main text, describing the low-energy band structure aroundK, K’ points of non-twisted honeycomb CCDW bilayer. Imposing only time-reversal symmetry, as relevant for thetwisted CCDW bilayer, implies that all terms appearing in equation (S7) are allowed by symmetry. Furthermore,the lack of inversion symmetry (coinciding here with a sublattice-symmetry breaking) removes all degeneracies of the4 gap K -10123450 0.5 1 1.5 2 P a r a m e t e r s [ m e V ] E z [mV/Å] M SOC B R -202468-6 -5 -4 -3 -2 -1 0 1 2 E z [mV/Å] a b Non-twisted λ λ
Fig. S3:
Parameters of the k ⋅ p expansion around Brillouin zone corners K, K ′ of the Wannier Hamiltonitans .Parameters for a, non-twisted and b, ○ twisted CCDW 1T-TaSe bilayer in vertical electric fields. The Semenoff mass M ,Kane-Mele SOC λ SOC , spin-valley coupling B , and Rashba SOC λ R are shown as function of vertical electric field E z . unperturbed non-relativistic bands, thus introducing an effective mass term M that acts as a staggered potential. Byintroducing the following parametrization: B = ( λ + λ ) λ SOC = ( λ − λ ) α R = ( α + α ) α R = ( α − α ) (S10)we recover the low-energy model equation (2) in the main text, describing the band structure in the vicinity of theK, K’ points.We conclude this appendix noticing that the low-energy Hamiltonian equation (S7) almost coincides with the onederived in Ref. [32] for a graphene-based system with C v symmetry, but for the additional spin-momentum couplingterm parametrized by λ D and allowed here by the lower C symmetry. In fact, a minimal tight-binding modelreproducing the low-energy band structure in the vicinity of K and K’ points can be derived following the generalscheme outlined in Ref. [32] and taking into account the reduced symmetry. Alongside the intrinsic SOC (next-nearestneighbor spin-conserving hopping, that is sublattice dependent due to lack of sublattice symmetry, parametrized hereby λ , λ ), the lack of any horizontal reflection allows for the so-called “pseudospin inversion asymmetric” SOC [62](next-nearest neighbor spin-flipping hopping, also sublattice dependent, parametrized here by α , α ) and, togetherwith inversion-symmetry breaking, for the Rashba SOC [29] (nearest-neighbor spin-flipping hopping, parametrizedby λ R ). The no-go arguments of Ref. [32] also implies that the lack of any vertical reflections in the C structuralsymmetry allows for a purely imaginary nearest-neighbor spin-conserving SOC hopping, whose coupling constant canbe parametrized by λ D .We obtain all effective parameters entering in the low-energy models around K and K’ from a first-order Taylorexpansion in k -space of the Wannier tight-binding models described above in terms of the external field E z , seeSupplementary Figure S3. E + − basis functions K P z , A z { K ω ω ∗ P x + iP y , A x + iA y K ω ∗ ω P x − iP y , A x − iA y TABLE S1: Character table for the little group at K. Here ω = e i π , and P ( A ) stands for polar (axial) vector. Notice thatpolar and axial vectors transform in the same way under the symmetry operations of C point group. Σ AA (k,i ω ) [meV] −0.8−0.7−0.6−0.5 Γ K K‘M a Σ AB (k,i ω ) [meV] −8−404 Γ K K‘M b Fig. S4:
Momentum dependent intra- and inter-sublattice parts of the self-energy as obtained from the TPSCcalculations . The a, intra- and b, inter-sublattice parts of the self-energy are obtained at T = .
005 eV =
58 K. The antiferro-magnetic fluctuations lead to a significant non-local inter-sublattice self-energy Σ AB . S4. TOPOLOGY
The topological properties of the twisted and non-twisted TaSe bilayers without the effects of correlation are studiedusing the Wannier charge center evolution, as described in [63]. Using the Wannier tight-binding Hamiltonians asinput, we determine the Z invariant to assess whether they describe a trivial or a quantum spin Hall insulator.For the non-twisted structure we find that the Z invariant equals 1 for all values of the electric field between 0 and0.43 mV ˚A − while Z =0 for all other E -field strengths. For the 180 ○ twisted bilayer we find Z =0 for all calculated E -field strengths instead.In order to understand why the 180 ○ twisted bilayer does not show any topologically non trivial phases we analyzethe influence of the terms entering the Hamiltonian H ○ from equation (2) of the main text. We proceed from theWannier Hamiltonian of the non-twisted bilayer and artificially add and vary parameters occurring in H ○ fromequation (2). Specifically, we consider a parameter space, where we vary the Semenoff mass M , the Rashba-spinorbit term λ R and the spin-valley coupling term / valley Zeeman B . To determine the topological properties ofnon-interacting Hamiltonians in this parameter space, it turns out to be numerically efficient to track the minimalgap in the band structure and to identify connected regions of the parameter space with strictly non-zero gap, whichare bounded by a submanifold of the parameter space with zero gap. If one point inside has Z ≠ =
0) the wholeregion is topologically non-trivial (trivial). The resulting topological phase diagrams are shown in Fig. 3c of the maintext.It is useful to compare our system to the well-known case of the ideal Kane-Mele model [29]. The characteristicshape of the phase-diagram with the Rashba coupling ( λ R ) on the y -axis and the staggered potential M (referred to asSemenoff mass) on the x -axis in the ideal Kane-Mele model resembles that of an onion [29]. We observe this behaviorfor the non-twisted case, see leftmost panel of Fig. 3c. The yellow horizontal line represents the trajectory of ourHamiltonian upon changing the electric field E z . The main effect of E z is to change M according to M ≈ M + eE z d / (cid:15) ⊥ (see Supplementary Figure S3) and hence E z drives the 0 ○ twisted system from inside the topological region to theoutside.Upon twisting, we reduce the symmetry from D h to C v consequently switching on various terms in H ○ , asdescribed in Ref. [32], among which the most important ones are the previously-mentioned Semenoff mass M , Rashba-spin orbit λ R and the spin-valley coupling B . The main effect of B is to shift two of the bands having the samesublattice but different spin character at the K-point towards each other. This yields a deformed onion as topologicalregion, where the topological region on the λ R axis is reduced. (See Fig. 3c of the main text for the case of B / λ SOC = . ∣ B ∣ = ∣ λ SOC ∣ the topological non-trivial region is completely suppressed and only avertical band touching line remains, which, however, does not separate a trivial from a non-trivial region. The lattercase is similar to the phase diagram of the 180 ○ twisted TaSe bilayer. The effects of B and λ R can, thus, qualitativelyexplain the differences observed for the topological phases of the non-twisted and 180 ○ twisted structure. In the 180 ○ twisted structure, λ SOC is additionally strongly reduced, which further contributes to suppressing the QSH phase.Following this line of argumentation the electric field in the 180 ○ twisted case induces a gap closing and a reopeningbut does not induce a topological phase transition.6In the non-twisted case, where we find QSH states in absence of interactions, we study the impact of correlationson the topological phase diagram within the TPSC approach. These calculations confirm the generic shape of theschematic phase diagram shown in (Fig. 2e).It is known [35] that sufficiently strong Coulomb repulsion U can change the nature of the topological phasetransition between QSH and band insulator from continuous with a band-touching point at the transition – as inthe non-interacting Kane-Mele or Bernevig-Hughes-Zhang models – to first-order with discontinuous jump of the gap(solid line in Fig. 2e). Since the electric field E z directly controls M , one can possibly tune the non-twisted systemthrough this exotic first-order transition in experiments by varying E z . S5. DMFT AND TPSC
As shown in the main text (Fig. 2), TPSC and DMFT consistently yield quasiparticle weights Z ≈ .
75 fortemperatures in the range 60-230 K for non-twisted CCDW 1T-TaSe bilayer. Also double occupancies ⟨ n ↓ n ↑ ⟩ ≈ . bilayer at intermediate local correlation strength and clearly far away from a paramagnetic Mott Hubbardtransition.The onset of non-local correlations can be inferred from the temperature dependent enhancement of antiferromag-netic susceptibilities and correlation lengths, which we obtained with TPSC and show in the main text (Fig. 2).These correlations manifest in sublattice off-diagonal contributions to the self-energy shown in Supplementary FigureS4. [1] Tsui, D. C., Stormer, H. L. & Gossard, A. C. Two-Dimensional Magnetotransport in the Extreme Quantum Limit. Phys.Rev. Lett. , 1559–1562 (1982).[2] Qi, X.-L. & Zhang, S.-C. Topological insulators and superconductors. Rev. Mod. Phys. (4), 1057–1110 (2011).[3] Sato, M. & Ando, Y. Topological superconductors: a review. Rep. Prog. Phys. (7), 076501 (2017).[4] Rachel, S. Interacting topological insulators: a review. Reports on Progress in Physics (11), 116501 (2018).[5] Hohenadler, M. & Assaad, F. F. Correlation effects in two-dimensional topological insulators. J. Phys.: Condens. Matter (14), 143201 (2013).[6] Wehling, T. O., Black-Schaffer, A. M. & Balatsky, A. V. Dirac materials. Advances in Physics (1), 1–76 (2014).[7] Marrazzo, A., Gibertini, M., Campi, D., Mounet, N. & Marzari, N. Prediction of a Large-Gap and Switchable Kane-MeleQuantum Spin Hall Insulator. Phys. Rev. Lett. (11), 117701 (2018).[8] Wu, X., Fink, M., Hanke, W., Thomale, R. & Di Sante, D. Unconventional superconductivity in a doped quantum spinHall insulator.
Phys. Rev. B (4), 041117 (2019).[9] Wilson, J. A., Salvo, F. J. D. & Mahajan, S. Charge-density waves and superlattices in the metallic layered transitionmetal dichalcogenides.
Adv. in Phys. (2), 117–201 (1975).[10] Rossnagel, K. On the origin of charge-density waves in select layered transition-metal dichalcogenides. J. Phys.: Condens.Matter (21), 213001 (2011).[11] Perfetti, L. et al. Spectroscopic Signatures of a Bandwidth-Controlled Mott Transition at the Surface of 1T-TaSe . Phys.Rev. Lett. (16), 166401 (2003).[12] Colonna, S. et al. Mott Phase at the Surface of 1T-TaSe Observed by Scanning Tunneling Microscopy.
Phys. Rev. Lett. (3), 036405 (2005).[13] Nakata, Y. et al. Selective Fabrication of Mott-Insulating and Metallic Monolayer TaSe . ACS Appl. Nano Mater. (4),1456–1460 (2018).[14] B¨orner, P. C. et al. Observation of charge density waves in free-standing 1T-TaSe monolayers by transmission electronmicroscopy. Appl. Phys. Lett. (17), 173103 (2018).[15] Chen, Y. et al. Strong correlations and orbital texture in single-layer 1T-TaSe . Nature Physics (2020).[16] M¨uller-Caspary, K. et al. Atomic-scale quantification of charge densities in two-dimensional materials.
Phys. Rev. B (12), 121408 (2018).[17] Hovden, R. et al. Atomic lattice disorder in charge-density-wave phases of exfoliated dichalcogenides (1T-TaS ). PNAS (41), 11420–11424 (2016).[18] Freericks, J. K., Krishnamurthy, H. R., Ge, Y., Liu, A. Y. & Pruschke, T. Theoretical description of time-resolvedpump/probe photoemission in TaS : a single-band DFT+DMFT(NRG) study within the quasiequilibrium approximation. physica status solidi (b) (5), 948–954 (2009).[19] Darancet, P., Millis, A. J. & Marianetti, C. A. Three-dimensional metallic and two-dimensional insulating behavior inoctahedral tantalum dichalcogenides. Phys. Rev. B (4), 045134 (2014).[20] Ritschel, T., Berger, H. & Geck, J. Stacking-driven gap formation in layered 1T-TaS . Phys. Rev. B (19), 195134(2018). [21] Lazar, P., Martincov´a, J. & Otyepka, M. Structure, dynamical stability, and electronic properties of phases in TaS froma high-level quantum mechanical calculation. Phys. Rev. B , 224104 (2015).[22] Castro Neto, A. H., Guinea, F., Peres, N. M. R., Novoselov, K. S. & Geim, A. K. The electronic properties of graphene. Rev. Mod. Phys. (1), 109–162 (2009).[23] Katsnelson, M. I. Graphene: Carbon in Two Dimensions , (Cambridge University Press2012).[24] Zantout, K., Altmeyer, M., Backes, S. & Valent´ı, R. Superconductivity in correlated BEDT-TTF molecular conductors:Critical temperatures and gap symmetries.
Phys. Rev. B , 014530 (2018).[25] Zantout, K., Backes, S. & Valent´ı, R. Effect of Nonlocal Correlations on the Electronic Structure of LiFeAs. Phys. Rev.Lett. (25), 256401 (2019).[26] Tran, M.-T. & Kuroki, K. Finite-temperature semimetal-insulator transition on the honeycomb lattice.
Phys. Rev. B (12), 125125 (2009).[27] Arya, S., Sriluckshmy, P. V., Hassan, S. R. & Tremblay, A.-M. S. Antiferromagnetism in the Hubbard model on thehoneycomb lattice: A two-particle self-consistent study. Phys. Rev. B (4), 045111 (2015).[28] Raczkowski, M. et al. The Hubbard model on the honeycomb lattice: from static and dynamical mean-field theories tolattice quantum Monte Carlo simulations. arXiv:1908.04307 [cond-mat] (2019). ArXiv: 1908.04307.[29] Kane, C. L. & Mele, E. J. Z Topological Order and the Quantum Spin Hall Effect.
Phys. Rev. Lett. , 146802 (2005).[30] Ezawa, M. Valley-Polarized Metals and Quantum Anomalous Hall Effect in Silicene. Phys. Rev. Lett. , 055502 (2012).[31] Zhang, X., Liu, Q., Luo, J.-W., Freeman, A. J. & Zunger, A. Hidden spin polarization in inversion-symmetric bulk crystals.
Nature Phys. (5), 387–393 (2014).[32] Kochan, D., Irmer, S. & Fabian, J. Model spin-orbit Hamiltonians for graphene systems. Phys. Rev. B , 165415 (2017).[33] Di Sante, D. et al. Towards topological quasifreestanding stanene via substrate engineering. Phys. Rev. B (3), 035145(2019).[34] Klein, J. et al. Electric-Field Switchable Second-Harmonic Generation in Bilayer MoS by Inversion Symmetry Breaking. Nano Lett. (1), 392–398 (2017).[35] Amaricci, A., Budich, J., Capone, M., Trauzettel, B. & Sangiovanni, G. First-Order Character and Observable Signaturesof Topological Quantum Phase Transitions. Phys. Rev. Lett. (18), 185701 (2015).[36] Kotov, V. N., Uchoa, B., Pereira, V. M., Guinea, F. & Castro Neto, A. H. Electron-Electron Interactions in Graphene:Current Status and Perspectives.
Rev. Mod. Phys. (3), 1067–1125 (2012).[37] Bistritzer, R. & MacDonald, A. H. Moir´e bands in twisted double-layer graphene. Proc Natl Acad Sci USA (30),12233 (2011).[38] Cao, Y. et al. Unconventional superconductivity in magic-angle graphene superlattices.
Nature , 43 (2018).[39] Cao, Y. et al. Correlated insulator behaviour at half-filling in magic-angle graphene superlattices.
Nature , 80 (2018).[40] Tsukahara, N. et al. Evolution of Kondo Resonance from a Single Impurity Molecule to the Two-Dimensional Lattice.
Phys. Rev. Lett. (18), 187201 (2011).[41] Pizarro, J. M., R¨osner, M., Thomale, R., Valent´ı, R. & Wehling, T. O. Internal screening and dielectric engineering inmagic-angle twisted bilayer graphene.
Phys. Rev. B (16), 161102 (2019).[42] Black-Schaffer, A. M. & Doniach, S. Resonating valence bonds and mean-field d -wave superconductivity in graphite. Phys.Rev. B (13), 134512 (2007).[43] Nandkishore, R., Levitov, L. S. & Chubukov, A. V. Chiral superconductivity from repulsive interactions in doped graphene. Nature Phys (2), 158–163 (2012).[44] Kiesel, M. L., Platt, C., Hanke, W., Abanin, D. A. & Thomale, R. Competing many-body instabilities and unconventionalsuperconductivity in graphene. Phys. Rev. B (2), 020507 (2012).[45] Black-Schaffer, A. M. & Honerkamp, C. Chiral d -wave superconductivity in doped graphene. J. Phys.: Condens. Matter (42), 423201 (2014).[46] Hohenberg, P. & Kohn, W. Inhomogeneous Electron Gas. Phys. Rev. , B864–B871 (1964).[47] Kohn, W. & Sham, L. J. Self-Consistent Equations Including Exchange and Correlation Effects.
Phys. Rev. , A1133–A1138 (1965).[48] Kresse, G. & Hafner, J. Norm-conserving and ultrasoft pseudopotentials for first-row and transition elements.
Journal ofPhysics: Condensed Matter (40), 8245–8257 (1994).[49] Kresse, G. & Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B ,1758–1775 (1999).[50] Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. ,3865–3868 (1996).[51] Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized Gradient Approximation Made Simple [Phys. Rev. Lett. 77, 3865(1996)]. Phys. Rev. Lett. , 1396–1396 (1997).[52] Grimme, S. Semiempirical GGA-type density functional constructed with a long-range dispersion correction. Journal ofComputational Chemistry (15), 1787–1799 (2006).[53] Grimme, S., Antony, J., Ehrlich, S. & Krieg, H. A consistent and accurate ab initio parametrization of density functionaldispersion correction (DFT-D) for the 94 elements H-Pu. The Journal of Chemical Physics (15), 154104 (2010).[54] Kamil, E. et al. Electronic structure of single layer 1T-NbSe : interplay of lattice distortions, non-local exchange, andMott–Hubbard correlations. Journal of Physics: Condensed Matter (32), 325601 (2018).[55] Gull, E. et al. Continuous-time Monte Carlo methods for quantum impurity models. Reviews of Modern Physics ,349–404 (2011).[56] Wallerberger, M. et al. w2dynamics: Local one- and two-particle quantities from dynamical mean field theory. Computer Physics Communications , 388–399 (2019).[57] Vilk, Y. M. & Tremblay, A.-M. Non-Perturbative Many-Body Approach to the Hubbard Model and Single-ParticlePseudogap.
Journal de Physique I (11), 1309–1368 (1997).[58] Mostofi, A. A. et al. wannier90: A tool for obtaining maximally-localised Wannier functions. Computer Physics Commu-nications (9), 685 – 699 (2008).[59] Yu, P. Y. & Cardona, M.
Fundamentals of Semiconductors , (Springer-Verlag Berlin Heidelberg2005).[60] Dresselhaus, M. D., Dresselhaus, G. & Jorio, A.
Group Theory - Application to the Physics of Condensed Matter , (Springer-Verlag Berlin Heidelberg2008).[61] Slonczewski, J. C. & Weiss, P. R. Band Structure of Graphite.
Phys. Rev. , 272–279 (1958).[62] Gmitra, M., Kochan, D. & Fabian, J. Spin-Orbit Coupling in Hydrogenated Graphene.
Phys. Rev. Lett. , 246602(2013).[63] Soluyanov, A. A. & Vanderbilt, D. Computing topological invariants without inversion symmetry.
Phys. Rev. B83