Topological features of magnetically-ordered and correlated crystals

Department of Physics Topological features of magnetically-ordered and correlated crystals Mikel García Díez Supervised by Juan Luis Mañes Palacios Maia García Vergniory Dissertation submitted to the University of the Basque Country UPV/EHU as partial fulfillment of the requirements for the Ph.D degree in Physics May 16, 2025 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 (cc) 2025 Mikel García Díez (cc by-nc-sa 4.0) REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 Abstract The classification of condensed matter systems in terms of topological classes has witnessed significant advancements over the past decade. Although transitions between topological phases cannot be rigorously established through changes in the symmetry of a local order parameter, symmetry nonetheless plays a crucial role in protecting topologically non-trivial states. This is particularly relevant in crystalline media, which, in addition to internal symmetries such as time-reversal and charge conjugation, exhibit a high degree of spatial symmetry described in terms of space groups. Recent progress in the comprehensive tabulation of magnetic space groups and their representations has facilitated the extension of symmetry-based topological methods, notably Topological Quantum Chemistry, to systems exhibiting magnetic ordering, such as ferromagnets. In this work, we employ these methods extensively to address a diverse range of problems concerning the classification and characterization of topological phases in crystalline media. Following a brief introduction to group theory, density functional theory, and topology in crystalline systems, we present an update to IrRep, a computational tool written in Python designed for symmetry analysis of ab-initio calculations. In this iteration, IrRep has been enhanced to support magnetic space groups and expanded to perform a comprehensive symmetry-based topological characterization using elementary band representations and symmetry indicators. With this improved tool, we first analyze two hexagonal compounds, Co and Fe 3 GeTe 2 , which are grouped together due to their shared magnetic space group symmetry. In both materials, nodal lines are protected by mirror planes, and we conduct an exhaustive search for these features near the Fermi level. In the case of Co, our theoretical predictions are corroborated by experimental angle-resolved photoemission spectroscopy measurements. For Fe 3 GeTe 2 , we further investigate the origins of the anomalous Hall conductivity in this material, identifying three principal contributions: nodal lines lifted by magnetic order, spin-orbit coupling gaps between spin-up and spin-down states, and Weyl nodes. Subsequently, we revisit the topology of the FeSe iron based superconductor. First, we model the electronic band structure of the doped FeTe 0.55 Se 0.45 crystal through DFT and a symmetry-based tight binding approach and shed light into the mechanism for i REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 ii the transition into a non-trivial Z2 state supporting a Dirac surface state as evidenced by the ARPES measurements of our collaborators. Second, we predict that uniaxial strain in pristine FeSe is capable of driving transitions into strong and weak topological phases by modifying the lattice constant along the a and c directions of the original tetragonal cell. Lastly, we explore the application of topological quantum chemistry to strongly correlated systems. In this context, we analyze the honeycomb Kitaev model and its analytical solution in terms of Majorana fermions to characterize the symmetry of spin excitations in the flux-free sector and predict the existence of topological edge states. This approach presents an alternative pathway for the systematic study of interacting phases, complementing existing methods based on Green’s functions and topological Hamiltonians. REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 List of publications I Topological semimetals without quasiparticles. Hu, H., Chen, L., Setty, C., GarciaDiez, M., Grefe, S. E., Prokofiev, A., Kirchner, S., Vergniory, M. G., Paschen, S., Cano, J., & Si, Q.. arXiv:2110.06182 (2021). II Cubic 3D Chern photonic insulators with orientable large Chern vectors. Devescovi, C., García-Díez, M., Robredo, I., De Paz, M. B., Lasa-Alonso, J., Bradlyn, B., Mañes, J. L., Vergniory, M. G., & García-Etxarri, A.. Nature Communications, 12(1) (2021). III Vectorial Bulk-Boundary correspondence for 3D photonic Chern insulators. Devescovi, C., García-Díez, M., Bradlyn, B., Mañes, J. L., Vergniory, M. G., & García-Etxarri, A.. Advanced Optical Materials, 10(20) (2022) IV Transversality-Enforced Tight-Binding Model for 3D Photonic Crystals aided by Topological Quantum Chemistry. Morales-Pérez, A., Devescovi, C., Hwang, Y., García-Díez, M., Bradlyn, B., Mañes, J. L., Vergniory, M. G., & García-Etxarri, A.. arxiv:2305.18257 (2023) V Orbital ingredients and persistent Dirac surface state for the topological band structure in FeTe 0.55 Se 0.45 . Li, Y., Chen, S., García-Díez, M., Iraola, M. I., Pfau, H., Zhu, Y., Mao, Z., Chen, T., Yi, M., Dai, P., Sobota, J. A., Hashimoto, M., Vergniory, M. G., Lu, D., & Shen, Z.. Physical Review X, 14(2) (2024). VI Axion topology in photonic crystal domain walls. Devescovi, C., Morales-Pérez, A., Hwang, Y., García-Díez, M., Robredo, I., Mañes, J. L., Bradlyn, B., GarcíaEtxarri, A., & Vergniory, M. G.. Nature Communications, 15(1) (2024). VII Band representations in Strongly Correlated Settings: The Kitaev Honeycomb Model. Fünfhaus, A., García-Díez, M., Vergniory, M. G., Kopp, T., Winter, S. M., & Valentí, R.. arxiv:2501.11396 (2025). VIII Origins of the anomalous Hall conductivity in the symmetry enforced Fe3GeTe2 nodal-line ferromagnet.García-Díez, M., Beidenkopf, H., Robredo, I., &Vergniory, M. G.. arxiv:2502.07420 (2025) iii REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 iv IX Flat band driven itinerant magnetism in the Co-pnictides (La,Ca)Co 2 (As,P) 2 . Subires, D., García-Díez, M., Kar, A., Lim, C.-., Li, V. M., Yannello, V., Carbone, D., Gargiani, P., Yilmaz, T., Dai, J., Tallarida, M., Vescovo, E., Shatruk, M., Vergniory, M. G., & Blanco-Canosa, S.. arxiv:2503.02728 (2025) X Concurrent multifractality and anomalous hall response in the nodal line semimetal Fe 3 GeTe 2 near localization. Mathimalar, S., Gupta, A., Roet, Y., Galeski, S., Wawrzynczak, R., Garcia-Diez, M., Robredo, I., Vir, P., Kumar, N., Schnelle, W., Karin, V. A., Küspert, J., Wang, Q., Chang, J., Sassa, Y., Stern, A., Felix, V. O., Vergniory, M. G., Felser, C., . . . Beidenkopf, H.. arxiv:2503.04367 (2025) XI Topological phase transitions of FeSe under uniaxial strain.García-Díez M., Profe, J.B., Davignon, A., Backes, S., Valentí, R. & Vergniory, M.G. In preparation. XII Manifold of magnetic nodal lines in and elemental ferromagnet. Clark, O.J., García-Díez, M., Fink, J., Rader, O., Miranda, R., Vergniory, M.G. & SánchezBarriga, J. In preparation. XIII Magnetic IrRep: symmetry and topological analysis of magnetic crystals from Density Functional Theory.García-Díez M., Iraola, M., Robredo, I., Mañes, J.L. S. Tsirkin, Stepan & Vergniory, M.G.. In preparation. REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 List of abbreviations BCS Bilbao Crystallographic Server BZ (first) Brillouin Zone DFT Density Functional Theory EBR Elementary Band Representation FQHE Fractional Quantum Hall Effect IQHE Integer Quantum Hall Effect (M)SG (Magnetic) Space Group (M)TQC (Magnetic) Topological Quantum Chemistry QAHI Quantum Anomalous Hall Insulator QSHI Quantum Spin Hall Insulator SI Symmertry Indicator SOC Spin-Orbit Coupling SSG Site-Symmetry Group TB Tight-Binding TI Topological Insulator TRIM Time-Reversal-Invariant Momentum TRS Time Reversal Symmetry WL Wilson Loop WP Wyckoff Position v REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 Contents Contents vii List of Figures ix List of Tables xi 1 Introduction 1 2 Group theory of crystalline systems and their band structures 11 2.1 The translational symmetry of lattices . . . . . . . . . . . . . . . . . . 12 2.2 The crystallographic point groups . . . . . . . . . . . . . . . . . . . . 17 2.3 Spacegroups................................ 18 2.4 Introducing time reversal: magnetic space groups . . . . . . . . . . . 25 2.5 Electrons in periodic systems . . . . . . . . . . . . . . . . . . . . . . . 32 3 Density functional theory 37 3.1 Statement of the problem . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.2 Independent electrons and Hartree-Fock approximation . . . . . . . . 39 3.3 The Hohenberg-Kohn theorems and the Kohn-Sham method . . . . . 40 3.4 Computation of the DFT ground state . . . . . . . . . . . . . . . . . . 44 4 Topology in condensed matter systems and its geometrical interpretation 47 4.1 Derivation of the geometric phase . . . . . . . . . . . . . . . . . . . . 48 4.2 Gauge theory formulation . . . . . . . . . . . . . . . . . . . . . . . . 50 4.3 Wilson loops: generalization of the Berry phase . . . . . . . . . . . . 51 4.4 The Chern number: a topological invariant . . . . . . . . . . . . . . . 55 4.5 Application to crystalline systems of electrons . . . . . . . . . . . . . 62 5 Topology from real-space symmetry: Magnetic Topological Quantum Chemistry 75 5.1 Elementary band representations . . . . . . . . . . . . . . . . . . . . . 76 vii REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 1. Introduction As condensed matter Physics researchers, one of the first things we often do when discussing a given material is to classify it. The highest-level classification one usually considers is between solid, liquid and gas, focusing mainly on the first group. Of course, this separation into three broad categories is insufficient for virtually any description that aims to accurately model the behavior of a physical system. For example, one soon finds that solids can be split into insulators, which are not able to conduct electricity under an external voltage difference, and metals, which do have low-energy excitations that can transport charge. Materials can change between phases depending usually on the thermodynamic conditions of the environment, typically pressure and temperature. Up to 1980, the description of all phase transitions was framed in the theory developed by Lev Landau, which distinguished each phase by its symmetries [1]. In this framework, the usual way to proceed is to expand the free energy of the system in terms of a local order parameter, usually physical quantities like a local magnetization, and model parameters that encode the effect of pressure and temperature. The actual value the order parameter takes is the one that minimizes the free energy given the external conditions. A change of such conditions incurs a change on the free energy landscape such that it is energetically favorable for the local order parameter to acquire a non-zero equilibrium value which, in turn, breaks at least one symmetry of the system. The prototypical case is a rotation-invariant system whose energy we expand in terms of a local magnetization that takes a non-zero value when the temperature is lowered. Since the local magnetic field chooses a direction in an otherwise invariant system, the symmetry is lowered to only operations that leave this axial vector invariant: rotations about its axis, inversion and their combinations. Landau theory has proven to be a very powerful method not only to discuss ferromagnetic transitions like in this example but superfluids and liquid crystals for example. However, as successful as Landau theory might have been, the experiments conducted by Klaus von Klitzing in that year on low-temperature two-dimensional electron gasses under a strong external magnetic field showed that this theory could not be the complete picture. In particular, von Klitzing observed that the resistivity of the samples transverse to the applied voltage displayed quantized plateaus of values h/e2·1/n , where n is integer, which increased with the value of the magnetic field [2]. This pointed at a series of phase transitions, one per plateau, that shared the same symmetries and hence could not be described under Landau’s paradigm. As we know now, this was the first realization of a topological phase transition and the integer quantum Hall effect (IQHE). This phenomenon also displayed a hallmark of topological effects: its robustness. The IQHE effect was not a mere fluke, it could be consistently reproduced for different samples and, moreover, the quantization of the 2 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 Figure 1.1: Longitudinal ρxx and transverse ρxy resistivity as a function of the magnetic field in a low-temperature 2D electron gas showing plateaus of quantized Hall resistivity. Reproduced from Ref. [2]. Hall conductance has now been measured to be extremely precise [3]. The fact that the IQHE is a topological effect was later suggested by Thouless et al. in 1982, when they derived an analytical formula for the Hall conductivity in a 2D electron gas [4]. The expression, which ultimately arises from the single-valuedness of the electronic wave function under closed loops in the periodic Brillouin zone (BZ), showed that the transversal conductivity followed the formula σxy =e2 2πℏC, where C is an integer called the first Chern number or, in the context of that article, the TKNN invariant. The robustness and precise quantization of the transverse conductance in the IQHE could then be explained by the fact that C cannot change smoothly unless there is a severe change in the magnetic field strength that causes a phase transition. The derivation also explained that the Hall resistivity changed due a number of protected surface states confined to the border of the 2D sample which are the origin of the transverse conductance of the electron gas. What makes the phenomenon perhaps even more striking is that the existence of these chiral edge states was predicted from bulk properties only, as considering periodic boundary conditions (which makes the BZ description possible) neglects the existence of sur3 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 1. Introduction faces completely. This feature was later recognized to be a case of the bulk-boundary correspondence, another general characteristic of topological effects in condensed matter systems. The robustness of these effects and exotic characteristics like the boundary modes attracted the attention of the condensed matter Physics community, both for the inherent interest in understanding topological phenomena and the possible applications in spintronics, quantum computing or superconductivity. In 1988, D. Haldane proposed a model of electrons on a honeycomb lattice that displayed IQHE in the absence of an external magnetic field, showing that there can be materials which are inherently topological [5]. This was modeled by introducing a complex hopping term whose circulation is zero, therefore showing no net magnetic flux through each plaquette, but nonetheless breaking time reversal symmetry (TRS). Otherwise, chiral edge states, which propagate in only one direction and are incompatible with a time-reversed picture, could not exist in the model. This was the first model of a Chern insulator, due to the topological invariant that distinguishes this phase. A few years later, in 2005, C.L. Kane and E.J. Mele proposed a model with both TRS and topologically protected edge states, essentially by coupling two layers of the Haldane model and considering spinful fermions with components coupled by spin-orbit coupling (SOC) effect [6]. Naturally, single boundary modes were forbidden by TRS but, tuning the parameters of the model, one can obtain an insulator displaying pairs of boundary modes. To satisfy the TRS constraints, the edge modes propagate in opposite directions and carry opposite spin. Due to this, these kind of systems were called quantum spin Hall insulators (QSHIs). Contrary to Chern insulators, where the invariant can be increased with no limit, QSHIs are described by a Z2 number taking only values {0,1} . This is because the hybridization of two different pairs of helical edge states is not forbidden by TRS and generically gaps any phase that has an even number of boundary modes, destroying the effect. Topological insulators were not confined to two dimensions, and experimental realizations of 3D systems like Bi 2 Se 3 , Bi 2 Te 3 or Sb2Te3were achieved in the following years [7,8]. As the topological nature of these phenomena was revealed, we realized why Landau’s theory of phase transitions failed to describe them. While the central object in the theory is a local order parameter, the information about topological effects is smeared globally throughout the system. Tracking how the electronic wave function changes along closed loops in the BZ, as in the TKNN formula, requires knowledge of the wave function along the path, in particular how the phase changes from one point in reciprocal space to the others. This does not mean that symmetry is irrelevant in the study of topological materials. The interplay between internal symmetries and topology crystallized in what we know as the ten-fold classification [9 – 11], which describes 4 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 the possible topological phases for different dimensions according to their transformation under time-reversal, charge conjugation and their combination. Additional work showed that crystal symmetries could also be related to the topological invariants. For example, L.Fu and C.Kane showed that the Z2 invariant in inversion-symmetric crystals could be computed directly with the inversion eigenvalues of the occupied set of bands at especial TRS-invariant wave vectors [12]. Fang et al. also discovered a correspondence between rotation eigenvalues and the Chern number [13]. These kind of materials, where the topology is protected by the crystal symmetries, were later called topological crystalline insulators (TCIs) by L.Fu [14]. Up to then, the diagnosis of topology in realistic materials was a difficult computational problem. Most of the analysis relied on complex calculations involving, for example, line integrals over closed paths in the Brillouin zone. For model Hamiltonians like Haldane’s the computation was fairly simple but real insulators, often described via Density Functional Theory (DFT) calculations, were much harder to analyze. However, The efforts to link symmetry eigenvalues to topological invariants showed a way towards a simpler method to diagnose topological properties on systems that, for the most part, showed a high degree of symmetry, as crystals. In 2016, a framework unifying all these methods was developed, called Topological Quantum Chemistry (TQC) [15 – 18]. TQC translated the problem of symmetry-protected topology from reciprocal to real space by relating the eigenvalues of symmetry operations at highsymmetry wave vectors with the localization [19,20] of Wannier functions sitting at specific Wyckoff positions in the unit cell of the crystal. The discovery of topological materials was reduced to an algebraic problem: compute the symmetry properties of the occupied states at all high-symmetry points of the BZ and check if the bands can be obtained by placing orbitals of certain symmetry at Wyckoff positions in the real-space cell. If this is possible, then the material is trivial or the topology is not protected by crystal symmetries, otherwise it is topological. The method also confirmed that, unless there is a gap closing and reopening in the electronic band structure, topology was insensitive to smooth (adiabatic) deformations of the energy levels. The key insight is that trivial insulators are adiabatically equivalent to arrays of isolated atoms at Wyckoff positions or atomic limits, so any other material that is smoothly tunable into that limit must also be trivial. Moreover, from that same analysis, one can extract the so-called symmetry indicators (SIs), which are quantities computed from symmetry eigenvalues that are related to topological invariants [21,22]. With this addition, TQC can not only identify topological materials but also tell the specific kind of topology, as long as it is protected by symmetry. Such an streamlined method opened the door to the systematic search for topological materials from all experimentally-realized stoichiometric crystals. The automation of a workflow including DFT calculations and the subsequent TQC analysis allowed 5 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 1. Introduction Symmetry eigenvalues Input structure DFT Is atomic limit? Not symmetry indicated Topological Trivial Wilson loops, etc. Symmetry indicators Topological invariants No Yes If symmetry indicated Figure 1.2: The TQC analysis workflow. Starting from an experimental input structure, DFT calculations are performed from which the symmetry eigenvalues at highsymmetry points in the BZ can be extracted. Then, one checks if this is an atomic limit. If it is not, then the material is topological, otherwise the topology cannot be diagnosed only by symmetry (2D electron gas for example). In that case, more complex methods such as Wilson loops are required. If the topological phase is symmetry-indicated, then the topological invariants can be inferred from the SIs. to increase the number of predicted topological materials to the order of tens of thousands [23,24], including old examples whose topology had been overlooked due to a lack of a straightforward way to identify it. Furthermore, all the machinery necessary for a TQC analysis such as space group symmetry tables, atomic limits and SI formulas, was made available to the public, accelerating the field of TCIs. However, the analysis of TRS-broken crystals such as those ferromagnets and antiferromagnets remained outside the reach of TQC, which only considered space group symmetry with TRS included. The efforts into tabulating all the representation for the 1651 Shubnikov groups, which describe all possible magnetic arrangements of atoms in a lattice with SOC, led to the extension of the method to magnetic TQC (MQTC) in 2021 [25]. With this complete enumeration of all space group symmetries, the analysis of important cases, such as Chern insulators which intrinsically break TRS due to the non-zero local magnetic moment distribution, was unlocked. The direct extension of TQC to magnetic materials immediately allowed the mass search of topological magnetic materials [26], with a complete topological classification thanks to the a 6 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 theory of SIs that included all Shubnikov groups. All of the previous cases are primarily understood withing the framework of non-interacting or weakly interacting electron models, in which we can rely on band theory to compute topological invariants, whether via symmetry indicators or other methods such as Wilson loops. However, real materials can show strong electron-electron correlations that significantly modify these topological features. The prime example of these correlated topological phases is the fractional quantum Hall effect (FQHE), where interactions between electrons leads to charge fractionalization and quasiparticles with charge smaller than the electron’s and exotic exchange statistics [27 – 32]. Electron correlations usually arise in crystalline materials due to the Coulomb repulsion between electrons in localized orbitals. This is especially important in 4f and 3d of heavy atoms where the stronger pull from the nucleus and the smaller number of nodes in the wave function leads to an enhanced localization. From the point of view of the band structure, localized electrons give rise to very little dispersing bands, which implies that their quasiparticle velocity is small and hence can interact more efficiently. Among these correlated materials, we find Mott insulators [33,34], heavy-fermion systems [35,36], Kondo insulators [37,38], quantum spin liquids [39] and twisted and Moiré systems [40,41]. There have been a number of proposals to address the topology of interacting system in a systematic way. For example, a number of expressions for topological invariants in term of of Green’s functions [42,43], which capture the dynamics of interacting media, have proposed with their main drawback being the difficulty to apply them. Other research shows that it is possible to use group cohomology, cobordism theory and higher-form symmetry to study “symmetry-protected” [44] and “symmetry-enriched” [45] topological phases, requiring complex mathematical concepts. An alternative path has also been put forward based on mapping the interacting system into a non-interacting picture, constructing a so-called “topological Hamiltonian” [46,47]. Recently, the use of TQC applied to this alternative single-particle picture has had its first results [48 – 50], opening the door to extending this systematic analysis to interacting systems. Research in this direction is still ongoing and rather promising, since otherwise the diagnosis of topology in correlated systems still relies on ad-hoc methods, tailored to specific cases. In this thesis, we will address the following problems and objectives: Extension of the open-source Python package IrRep to magnetic materials and upgrade of enable the full TQC analysis via elementary band representations and SIs. Previous versions of IrRep [51] could not be applied to materials showing order of the local magnetic moments of the atoms, for example with 3d and 4f orbitals. 7 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 1. Introduction Thanks to the tabulation of all the symmetry information for the Shubnikov groups, we can now extend the analysis to magnetic materials still maintaining a simple, easy-to-use interface with the user with new quality-of-life improvements. Moreover, the package is now able to perform a full symmetry-indicated topology analysis, performing the decomposition of the selected set of bands in terms of elementary band representations and computing all the corresponding symmetry indicators for both single and double space groups. This enables a fully local and automatized workflow for the study of topological materials at the single-particle level. Using magnetic space groups to diagnose systems with topological nodal lines and their anomalous Hall conductivity. We perform an exhaustive symmetry analysis of Fe 3 GeTe 2 and hexagonal close-packed (hcp) Co that reveals how the ferromagnetic order can give rise and also destroy topological nodal lines, which we characterize via Wilson loop calculations. Additionally, on the one hand, Fe 3 GeTe 2 has been experimentally shown to display large anomalous Hall conductivity and our precise symmetry understanding combined with materials modeling via Wannier function models unveils the intrinsic sources of this effect. On the other hand, we present experimental results on the topological nodal lines in hcp-Co, confirming the predicted features via density functional theory and symmetry analysis. Study of a pressure-induced topological transition in the FeSe superconductor. FeSe is one of the most studied members of the family of iron-based superconductors. Despite its structural simplicity, a number of competing phenomena such as nematicity, magnetic order and orbital-selective correlations interact with superconductivity. The sensitivity of the electronic properties with respect to the Fe-Se and Fe-Fe distances has been used to tune its properties (such as the critical temperature) through external pressure. Additionally, chemical pressure through Te doping has been shown to achieve a topological transition into a strong topological insulator (TI). In this work, we perform density functional theory calculations of FeSe under uniaxial strain and a subsequent analysis in terms of TQC and SIs. Contrary to what previous research suggests, FeSe is a weak TI at ambient pressure. Moreover, we also find that uniaxial strain drives the system into an orthorhombic phase and prompts a topological phase transition into a strong TI, suggesting another route to tune the topological properties of FeSe. Topological transitions in the FeSe superconductor via doping and pressure. FeSe is one of the most studied members of the family of iron-based superconductors. Despite its structural simplicity, a number of competing phenomena such as nematicity, magnetic order and orbital-selective correlations interact with superconductivity. The sensitivity of the electronic properties with respect to the Fe-Se and Fe-Fe dis8 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 tances has been used to tune its properties (such as the critical temperature) through external pressure. Additionally, chemical pressure through Te doping has been shown to achieve a topological transition into a strong topological insulator (TI). In this work, we perform density functional theory calculations of the FeTe 0.55 Se 0.45 , a very well studied doped sibling compound, and develop a symmetry-based tight-binding model that help explain the electronic features observed in angle-resolved-photoemission spectroscopy measurements. Our work sheds light into the the mechanism for the transition into a strong TI phase upon doping that hosts surface Dirac states, as evidenced by the from our collaborators. We also propose uniaxial strain in pristine FeSe as a way to drive similar topological phase transitions. Based on an analysis in terms of TQC and SIs, we find that modifying the a and c lattice constants of the original tetragonal cell gives rise to strong and weak TI phases, depending on the direction and the strength of the external perturbation. Results on the application of TQC to highly-correlated systems in the honeycomb Kitaev model. We consider the extension of TQC to systems with topological order and strong correlations. The particular case of the honeycomb Kitaev model allows us to directly apply the method using the non-interacting picture obtained by transforming the spin degrees of freedom into Majorana fermions, without relying on the definition of a topological Hamiltonian. In this setting, we propose a concept of “spin orbitals” composed of more than one spin operator which give rise to a band structure of excitation energies. The analysis of this band structure with TQC and SIs reveals that, upon introducing a time-reversal-breaking term corresponding to an external magnetic field, the system can be tuned into a Chern insulator phase with chiral edge modes and we prove it by computing the spectrum of a finite stripe geometry. Thus, this work sheds light on the question of how TQC can be generalized to strongly correlated systems, adding to the previous results we described above. The thesis is structured as follows. In Chapters 2 and 3 we give an introduction to the group theory of Shubnikov space groups and density functional theory. Chapters 4 and 5 concentrate on the topology of condensed matter systems and its diagnosis through TQC and SIs. The results start with Chapter 6, where we describe the new functionalities of IrRep for the topological analysis of magnetic materials. In Chapter 7, this is used to unveil the nodal lines in Fe 3 GeTe 2 and hcp-Co. Chapter 8 studies the effect of doping and uniaxial strain in the FeSe superconductor. Finally, Chapter 10 explores the application of TQC in the honeycomb Kitaev model of strongly correlated spins. 9 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 CHAPTER 2 Group theory of crystalline systems and their band structures 11 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures Primitive Body-centered Face-centered Figure 2.2: The three cubic Bravais lattices. Gray nodes correspond to new points due to the centering with respect to the primitive lattice. 2.3 Space groups Space groups are the combination of the translational an rotational symmetries of real crystals. Their operations are usually denoted with Seitz symbols {R|t} , where R is an element of a point group and t is a translation. From the composition of operations, we can derive the inverse of any of them {R2|t2}{R1|t1}={R2R1|R2t1+t2}=⇒ {R|t}−1={R−1|−R−1t}.(2.12) The simplest way to obtain a space group from a Bravais lattice and a point group is by directly combining the two to give a symmorphic space group. In these, both groups are essentially separated: the group of rotations is by itself a point group and the translations t are always translations of the Bravais lattice. Only 73 of them can be constructed in this way. The remaining 157, for a total of 230, are called non-symmorphic. In these, the operations may have a translational part v which is not part of the Bravais lattice. This translational part must be a fraction of a lattice translation according to the order of the operation. For example, we may find elements 18 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.3. Space groups Cubic Tetragonal Orthorhombic Trigonal Hexagonal Monoclinic Triclinic 120º Figure 2.3: The seven crystal systems. Edges with different labels are assumed to be of different length in general. Unlabeled angles are also assumed to be 90◦. 19 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures of the form {2001|0,0,1/2}, where 2001 is a rotation about the third axis of the basis and the translation is given also in this basis. In this example, it is easy to see that we cannot separate rotations and translations since {2001|0,0,1/2}2 is a translation of the lattice with no rotational part and thus there is not a closed point group within the space group. Representations of space groups The process of obtaining the representations of a given space group starts by picking a basis in which the translations are diagonal, e.g the Bloch basis. First, notice that a state ψkat kis transformed by an operation {R|t}into a state at Rk {E|t}{R|v}ψk={R|v}{E|R−1t}ψk=e−ik·R−1t{R|v}ψk=e−iRk·t{R|v}ψk, (2.13) where we have used that R is unitary so k·R−1r=Rk·r . Let G be the space group of the crystal, then we define Definition 2.3. The subgroup Gk⊆G of operations that leave the wavevector k invariant up to a reciprocal lattice translation is called the little group of k Gk={g∈G:|gk=k′+K},(2.14) where Kis a reciprocal lattice vector. Notice that, by definition, Gk is also a space group, since translations have no effect on reciprocal lattice vectors. The little co-group ¯ Gk of k is obtained by taking the different rotational parts of the elements of Gk and is isomorphic to a point group. Gkis either Gor a strict subgroup of G, which motivates the following definition Definition 2.4. Given a wavevector k , the set of non-equivalent vectors obtained by acting with all the operations of Gon kis called the star of k Star(k) = {ki∈1BZ|∃g∈G:gk=ki∧k=ki+K}.(2.15) The wavevectors of the star are usually called arms of the star. The little groups of the vectors of the star are actually isomorphic and conjugate in the sense that gk=q=⇒Gq=gGkg−1.(2.16) The usual way to proceed is to induce a representation of G from the smaller subgroup Gk and those of the arms of the star. All the irreps of G are first labeled 20 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.3. Space groups by a wavevector k which indexes a star, where k must be in the irreducible part of the BZ. Once k is fixed, the space of states of the star of k is invariant under any operation of G . The full space group irreps are formed by blocks corresponding to an irrep of the little group Gk. In more detail, we can decompose Ginto cosets of Gk G=[ gα/∈Gk gαGk(2.17) such that given a representation σ of Gk , we can induce a representation Σ of G , which is denoted Σ=σ↑G . Let (α, β) index the blocks of Σ and (i, j) the elements of Σ, then for any g∈G [Σ(g)]iα,jβ = [˜σ(g−1 αggβ]]ij,(2.18) where [˜σ(g)]ij =([σ(g)]ij ifg∈Gk 0 otherwise .(2.19) If σ is d -dimensional irrep of Gk and assuming there are m arms, the full irrep Σ is (m×d) -dimensional. The block structure then arises from the that given a general {R|t}∈G, we can always find a relation with an element {h|v}∈Gkof the form {R|t}{Rβ|tβ}={Rα|tα}{h|v},(2.20) where Rαk=kα and Rβk=kβ are arms of the star. The only non-zero (α, β) blocks of Σ induced from σ are precisely those for which the above relation holds. Intuitively, focusing on the rotational part R=RαhR−1 β , the effect of {R|t} on the states at kβ is related to bringing them to k , applying the rotation h according to the irrep of Gk and then moving them to kα. To find all the irreps of G one must also find all the irreps of Gk in the first place. This poses a problem since Gk is infinite and therefore has an infinite number of irreps. One can actually find a finite number of them by realizing that Gk is a space group that can be decomposed into left cosets Gk=T∪[ gi/∈T giT, (2.21) where gi∈Gk . We first point out that there is a finite number of cosets. When Gk is symmorphic, the point group operations of its little co-group and the translations are essentially separated, i.e., Gk=¯ GkT . In this case, a finite set of irreps can be found from the irreps of ¯ Gk , which is finite. This is not the case for non-symmorphic space 21 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures groups. In that situation, two additional constraints can be imposed. First, we will only accept irreps σin which pure lattice translations have the form ∆σ(t) = e−ik·t 1 d,(2.22) where t∈T , d is the dimension of σ and ∆σ(t) is the matrix representation of t in the irrep σ . These are called small irreps of Gk and have the property that elements which differ by a lattice translation have the same matrix representation up to a phase ∆σ(tgi) = e−ik·t∆σ(gi),(2.23) with gi∈Gk and tgi∈giT . Second, we also impose that two small irreps σ and σ′ are equivalent if there exists a unitary N (the same for all elements of Gk ) that relates them ∆σ(g) = N∆σ′(g)N−1,∀g∈Gk⇐⇒ σ′≡σ. (2.24) Along with the restriction on the representation of translations, this is enough to find only a finite set of non-equivalent irreducible small representations of Gk and thus a finite set of full irreducible representations of G. Although being familiar with the fundamental properties, especially the block structure of the representations of space groups, is certainly beneficial, it is not necessary to know the details in practice. All the irreps of all the space groups for high symmetry points and lines of the BZ are available at the Bilbao Crystallographic Server website [52], which it is assumed that is used from this point onward. Induction of space group representations from real space orbitals. Consider a crystal with the symmetry of some space group G . The structure is composed of atoms in a unit cell which are coupled when brought close to one another due to electromagnetic interactions of the electrons and ions. The discrete single-atom levels get distorted in such a way that, upon taking the Fourier transform and for well-defined quasi-particle energies (such as in the case of the single-electron approximation), a description in terms of energy bands arises. Each state (or set of degenerate states), which carries a label of a reciprocal k vector, transforms according to a space group representation in a space formed by the star of k , as explained previously. Therefore, we see that atomic orbitals located in the crystalline structure give rise to a set of bands transforming as a space group representation. Each atomic orbital displays a particular symmetry due to the its environment described by a symmetry group. It turns out that the representations of this group can be promoted to a full representation of the space group, as we will see next. 22 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.3. Space groups To begin the discussion, given a position q in the unit cell, we define the Wyckoff position (WP) as the set of locations obtained by acting with all the symmetry elements of G on q that still lie in the same unit cell. This is sometimes called the orbit of q , with the restriction that only positions in the same unit cell are retained. WPs are labeled by a number, which is the number of unique sites or multiplicity, and a letter which varies for each space group (e.g. 2a). WPs can be viewed as the real-space analogue of Star(k) : each site in the WP is invariant under a subset of G called the site-symmetry group (SSG), which we will denote Gq: Gq={g∈G|gq=q}.(2.25) The difference is that Gq is finite, since it cannot contain lattice translations by definition, and isomorphic to a point group, whereas the little-group of k is infinite (and a space group in fact). The SSGs at different sites are isomorphic and conjugate, that is, if qα=gαqthen Gqα=gαGqg−1 α. In general, WPs can have fixed positions or can be parameterized by some coordinate, for example (1/2,0, z) where z∈(0,1/2) , such that it connects two different WPs at (1/2,0,0) and (1/2,0,1/2) . A Wyckoff position is called maximal when the SSG of any of its sites is not a proper subgroup of any of the sites of another WP that it is connected to. In our example, (1/2,0,0 ) may display more symmetry than the more general (1/2,0, z)such that the SSG of the latter is a subgroup of the former. Once in a crystalline environment, the orbitals of an atom at q of a given WP are distorted. The single-atom orbitals transform under the O(3) rotation group with inversion, whose representations can be subduced into representations of Gq , such that we can form linear combinations of spherical harmonics that form a basis of each irrep of the SSG. This basis is more natural because it is adapted to the crystal symmetry and it is the one we will use from this point. Let ρ label one irrep of the Gq which belongs to a WP with sites {q,q1,q2...qn} such that qα=gαq for some gα∈G\Gq . The process of inducing an space group representation Σρ from ρ is similar to the one used to induce it from little-group representations, relying again in a coset decomposition of the form: G= n [ α=1 gα(Gq⋉T),(2.26) where T is again the group of lattice translations and ⋉ denotes the semi-direct product5. The full space group representation Σρis the induced from ρsuch that [Σρ(g)]iαt,jβt′= [˜ρ(g−1 α{E|t}g{E|t′}gβ]]ij,(2.27) 5The semi-direct product applies here because Tis a normal subgroup of Gand Gq∩T=∅ 23 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures where [˜ρ(g)]ij =([ρ(g)]ij ifg∈Gq 0 otherwise .(2.28) Therefore, we see that the same intuition applies as for the induction from the little group of k : the real-space representation Σρ moves orbitals between sites of the WP and rotates them according to ρ of Gq . For states at qα , they are first shifted to q , where the rotation according to ρ is applied, and then moved to the final position qβ . If g∈Gqα then qα=qβ . The role of the lattice translations {E|t} and {E|t′} is to adjust for the translations mismatches due to gα and gβ , which may contain translations by themselves (even non-integer ones). The choice of coset representatives gα amounts to the selection of a specific basis. For example, if ρ is spanned by pz orbitals only and q1 is related to q by the inversion {I|0} such that it is chosen as representative, then the basis of Σρwill have a −pzat q1. Once this process is established, Σρ can be brought to reciprocal space by a Fourier transform. Let ϕi α be the set of states transforming as ρ of Gq (or its analogues in the conjugate Gqα). A change of basis to reciprocal space gives ϕi α(k,r) = 1 √NX T eik·Tϕi α(r−T).(2.29) Following the expression for Σρ , the action of g∈G on ϕi α(r−T) at some unit cell Tcan be shown to be Σρ(g)ϕi α(r−T) = 1 √N dim(ρ) X j=1 [ρ(h)]jiϕi β(r−RT−tαβ),(2.30) where tαβ =gqα−qβif gqα=qβand h∈Gqis the unique element that satisfies g={E|tαβ}gβh. (2.31) It follows that the action in the reciprocal-space basis is Σρ(g)ϕi α(k,r) = e−iRk·tαβ √N dim(ρ) X j=1 [ρ(h)]ijϕi β(Rk,r),(2.32) which we see changes a state at k to one at Rk , as it should. The space-group representation so induced is called a band representation, usually denoted as ρ@WP ↑G indicating that it was induced to G from the site-symmetry group representation ρ at the Wyckoff position "WP”. Since there are n sites in the unit cell in the WP, dim(ρ)×nbands are obtained in reciprocal space. 24 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.4. Introducing time reversal: magnetic space groups 2.4 Introducing time reversal: magnetic space groups Conventional space groups only describe the spatial symmetries of the atomic arrangement of crystals. However, there may be additional physical properties that must be taken into account in the symmetry analysis of a condensed matter system. Such is the case of the magnetization of atoms in ferromagnetic, anti-ferromagnetic and ferrimagnetic crystals. In general, in the absence of spin or magnetic fields, the most common situation is that physical systems show time-reversal symmetry (TRS). When magnetization is included in the picture, the symmetry group may not include TRS by itself, but a combination of TRS and some space group operation. The resultant space group is called a magnetic space group (MSG), which are part of the set of all space groups, called Shubnikov groups, which include those that we have discussed in the previous sections. Kramers’ degeneracy Time reversal acts trivially on real space but may have an effect on other physical quantities, which for the purpose of this work are: • It reverses magnetic fields and spin components (and other angular momenta). •Changes the wavevector kto −k It is implemented by an anti-unitary operator, which we will label Θ , and acts on a Bloch state ψk(r)as Θaψk(r) = a∗Θψk(r) = a∗˜ ψk(r),(2.33) where we note that ¯ ψk(r) is the time-reversed partner of ψk and Θ is anti-linear. It can be proven that Θ is always represented by Θ = UK , where U is a unitary matrix of the dimension of the space it acts on and K denotes complex conjugation of every term sitting on its right. For example, for a spin-1/2 Θ = ηeiπSy/ℏK=−iησzK, (2.34) where Sy is the y component of the spin operator and σy is the associated second Pauli matrix. In practice, U can always be found using the fact that Θ commutes with all the space group operations and that it is unitary. 25 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures Object Inversion Time-reversal r−r r k−k−k S S −S Table 2.1: The transformation properties under inversion and time reversal of realspace vectors r, reciprocal-space kand spin S. Another important property is that it behaves differently for fermions and bosons, in particular Θ2=(−1fermions +1 bosons .(2.35) This has far-reaching consequences in the electronic systems that concern this thesis. Consider, in the most general case, a state |ϕ⟩ and its time-reversed partner Θ|ϕ⟩ . We may ask whether these two are really the same state, up to a phase: Θ|ϕ⟩? =eiφ|ϕ⟩.(2.36) Crucially, in the case of electrons Θ2=−1so we have Θ2|ϕ⟩=e−iφeiφ|ϕ⟩=|ϕ⟩=−|ϕ⟩.(2.37) We see that we arrive at a contradiction because our assumption of the two states being proportional is false. Therefore, for electronic states, the time-reversed partner is in fact a different state and, if Θis indeed a symmetry, they have the same energy, which is doubly degenerate. This simple derivation for a non-degenerate energy can be extended to the degenerate case, so all the energies in a fermionic, timereversal-invariant system are at least doubly-degenerate. This is called Kramers’ degeneracy. In system with translational symmetry, the action of Θon Hamiltonian at kis ΘHkΘ−1=H−k(2.38) and the state of energy En at k , |ψnk⟩ is related to Θ|˜ ψnk⟩=|ψn−k⟩ , which has the same energy. Additionally, if k≡ −k , if follows that all energies at k are at least doubly degenerate. This special set of points are the Time-Reversal-Invariant Momenta (TRIMs), which are located at the corners of the BZ cube (in reduced coordinates and three dimensions) and Γ = (0,0,0) . Moreover, if inversion I is also a symmetry, Iθ is an anti-unitary symmetry as well, which means that all energies at all kpoints are at least doubly-degenerate. 26 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.4. Introducing time reversal: magnetic space groups Classification of magnetic space groups Once the TR operator is introduced, we may divide the resulting 1,651 Shubnikov space groups into four categories. The most basic one is groups that simply do not contain the TR operator, which is what we have been working with implicitly, and are named Fedorov or Type-I groups. Note that some configurations of magnetic moments may not present TRS in any way, not even combined with space group operations. Another possibility is that Θ by itself is a symmetry of the system, which is the case of Type-II or gray space groups. These have the structure MII =G∪ΘG, (2.39) where G is a unitary (Fedorov) space group. Consequently, MII has double the operations of G : half of them correspond to the original elements and the other half to the same ones combined with Θ . In particular, the identity can be combined with TR so Θ is an element in its own right. Consequently, MII cannot describe crystals with a non-zero magnetic moments, since the ground-state magnetization at every site would be reversed by Θ, leading to a distinctly different state. The third case occurs when Θ is not a symmetry but some of its combinations with unitary space-group operations are, giving as a result a Type-III MSG with the general structure MIII =H∪Θ(G\H),(2.40) where H⊂G is the so-called halving subgroup because it has index two in G , i.e., it has half the operations. As a consequence, and as it happens for all magnetic groups other than Type-I, half the elements are unitary and the other half are anti-unitary. The former are obtained by combining the remaining elements in the set difference G\H with Θ . Type-III MSGs can describe both ferromagnetic and antiferromagnetic configurations. Finally, Type-IV MSGs differ from Type-III in that there are pure translations combined with TR. Consequently, the general structure is MIV =H∪Θt0H, G =H∪t0H, (2.41) where H is isomorphic to a Type-I MSG and t0 is a centering translation whose length is half of either of the following a+b+c,a+b,a+c,b+c,a,b,c,(2.42) 27 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2. Group theory of crystalline systems and their band structures where Hij couples the subspaces of irreps ρiand ρj. From Eq.(2.67), we have Γ(g)HΓ†(g) = ρmδmiHijρ† jδjn =ρmHmnρ† n=Hmn.(2.70) Therefore, considering H as a G -linear map, only the blocks between equivalent irreps are non-zero according to the first lemma. Furthermore, the second lemma states that every block between equivalent irreps is proportional to the identity matrix. For the sake of clarity, consider a case where Γ=ρ1⊕ρ2⊕ρ2 . Then Wigner’s theorem is the statement that the matrix Hhas the form H=   H(1) 11 0 0 0H(2) 11 H(2) 12 0H(2)† 12 H(2) 22   ,(2.71) where we have already used that ˆ H is Hermitian and H(i) is a block that acts on the subspace of an irrep ρi . Moreover, the blocks for irreps of multiplicity greater than one, such as ρ2 , can be block-diagonalized by reordering the bases. Assume that each ρ2 has a basis {v1, v2, . . . , vd} and {w1, w2, . . . , wd} , where d is the dimension of ρ2 . Then, by taking a basis {v1, w1, v2, w2, . . . , vd, wd} for ρ2⊕ρ2 we achieve a fully block-diagonal form of the Hamiltonian matrix H=   H(1) 11 0 0 0˜ H(2) 11 0 0 0 ˜ H(2) 22   .(2.72) We are interested in using this result for a system obeying the symmetry of a space group. First, we consider the Bloch basis, which diagonalizes the translation operators such that each translation Tt by a vector t of the lattice is assigned a diagonal matrix proportional to e−ik·t . This implies that H can be decomposed into a diagonal form with blocks H(k) which are decoupled by symmetry. Furthermore, each H(k) is invariant under the little group of k , so Wigner’s theorem can be further applied to infer the form of each block. Considering a symmetry-adapted basis for the subspace corresponding to k , we know that we can diagonalize H(k) analogously to Eq.(2.72). Finally, given this structure, we can see that each d -dimensional irrep of the little group of k gives rise to a d -dimensional subspace of states that have the same degenerate energy and span a basis of the irrep. This exact degeneracy is not accidental and is protected as long as the symmetries of the little group are preserved. However, group theory does not predict the actual values of these energies and they depend on the details of the crystal. Operations that are not in the little group of k relate the states at one wave vector with the rest of the star, implying that every energy is repeated on every arm with the same degeneracy. 34 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 2.5. Electrons in periodic systems System with spin-1/2 particles Before concluding this chapter, a final note on the treatment of spin-1/2 particles is due. When spin-orbit coupling (SOC) is neglected, the spin up and down sectors can be effectively decoupled, treating electrons in each subsystem as scalar. However, when SOC is considered, we have to give the excitations a full spinorial treatment. There are two approaches to the same problem. The first one consists in taking into account that the electron wave function has orbital and spin parts. An element g in the space group Gacts on both degrees of freedom as g|ψjσ⟩=D(g)ijS(g)σ′σ|ψiσ′⟩,(2.73) where |ψiσ⟩ denotes a state of orbital part i and spin σ . Therefore, we are effectively considering that the spinful excitations transform under a double-valued irrep Γ=D×S , which is the product of the orbital and spin representations. The caveat is that the full rotation group acting on spin is SU(2) rather than SO(3) . Because the former is a double cover of the latter, every real-space rotation by an angle θ actually corresponds to two possible spin rotations that differ by a minus sign: one through θ and another one through θ+ 2π . To a rotation through an angle θ along an axis ˆn , we can assign the SU(2) matrix for spin S=e−iθ 2(ˆn·σ)= cos(θ 2)−isin(θ 2)(ˆn ·σ),(2.74) where σ is the vector of Pauli matrices. The mapping makes explicit that there is oneto-two correspondence due to the θ/2 factor. Instead of bothering with the product representation, the second and equivalent approach is to consider that the space group has double the elements, giving rise to double groups, where each new element is distinguished by a 2π additional rotation from the original one. The matrices for the spin rotation and double irreps differ by a minus sign between the two flavors of the same real-space operation and the inversion is assigned the identity. 35 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 CHAPTER 3 Density functional theory 37 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3. Density functional theory Even in small samples, the sheer number of particles involved in the physical behavior of a material renders the problem completely intractable if exact accuracy is to be preserved. Some methods, such as exact diagonalization, can indeed solve the exact physics of a model Hamiltonian without any simplification. However, they are usually limited by the dimension of the Hilbert space of the system, which grows exponentially with its size. For this reason, as usual, we require certain degree of approximation that allows to retain the most important aspects of the physics involved in macroscopic materials. In particular, since we are interested in the topology of the electronic eigenstates of crystalline systems, we wish to be able to compute the energy levels and eigenfunctions for a wide range of materials. This chapter is dedicated to explain one of the most extended methods, called Density Functional Theory or, by its acronym, DFT. 3.1 Statement of the problem The physics of the ionic and electronic degrees of freedom in condensed matter systems such as crystals, leaving aside spin for simplicity, are governed by a universal Hamiltonian ˆ H=−ℏ2 2mX i∇2 i−X i,I ZIe2 |ri−RI|+1 2X i=j e2 |ri−rj| −X I ℏ2 2MI∇2 I+1 2X I=J ZIZJe2 |RI−RJ|. (3.1) The corresponding terms are, in order of appearance: 1. Kinetic energy of the electrons of mass m. 2. Coulomb attraction of each electron of charge e at position ri to every other ion of charge ZIeat position RI. 3. Coulomb repulsion of each electron at riwith every other electron at rj. 4. Kinetic energy of the ions with mass MI. 5. Coulomb repulsion of each ion at RIwith every other ion at RJ. We will neglect the kinetic energy of the atoms, as it lies in a much lower energy scale due to the higher masses MI compared to the mass of the electrons. The problem is thus simplified to a picture where the electrons move in a fixed, periodic array 38 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3.2. Independent electrons and Hartree-Fock approximation of ions of positive charge. This is called the Born-Oppenheimer approximation. From here onward, we will use Hartree atomic units, which amounts to setting ℏ=m=4π ϵ0= 1 , where ϵ0 is the vacuum permittivity. Thus, in this approximation, the Hamiltonian reads ˆ H=−1 2X i∇2 i−X i.I ZI |ri−RI|+1 2X i=j 1 |ri−rj|+ +1 2X I=J ZIZJ |RI−RJ|=ˆ T+ˆ VIe +ˆ Vee +VII, (3.2) taking into account that the ion-ion energy VII is a constant since the inter-atomic distance remains fixed (apart from small displacements due to lattice vibrations). The complexity resides precisely in the ˆ Vee term, since it involves the interaction of the electrons by pairs and, as such, cannot be separated in one-particle operators. As a consequence, the true ground state wave function of the system depends at the same time on all the electron positions Ψ(r1,...,rN) , with a total of 3N degrees of freedom (three per particle in three dimensions with N electrons), which is huge already even neglecting spin. The problem is now how to find methods that can approximate the electron-electron interaction term. 3.2 Independent electrons and Hartree-Fock approximation The most basic approximation to the electronic problem is to assume that the electronic degrees of freedom can indeed be separated. The wave function Ψ thus is a product of single-electron functions, which only depend on one position. Since electrons are identical fermions, Ψ should also be anti-symmetric with respect to the exchange of two particles. The resulting form of the state of the system is known as a Slater determinant Ψ = 1 √N!  ϕ1(r1, σ1)ϕ1(r2, σ2)··· ϕ1(rN, σN) ϕ2(r1, σ1)··· ··· ··· . . .. . .. . .. . . ϕN(r1, σ1)ϕN(r2, σ2)··· ϕN(rN, σN)  ,(3.3) where each electron has the position ri and spin σi degrees of freedom. The antisymmetry of the determinant when exchanging two rows or columns ensures that Ψ is anti-symmetric when exchanging (ri, σi)↔(rj, σj) . The pre-factor is there to normalize the state to unity. 39 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3. Density functional theory This approximation is known as the Hartree-Fock method. Although it uses a function of independent electrons, some correlations are still included. In particular, forcing the anti-symmetry of Ψ includes the correlation due to the Pauli exclusion principle. That is, the only two-body correlation included is the one that prevents two electrons with the same spin being at the same position in space. Using this form of the solution in the Hamiltonian in Eq.(3.2) yields the expectation value ⟨Ψ|H|Ψ⟩=X i,σ Zϕσ∗ i(r)−1 2∇2+VIe(r)ϕσ i(r)d3r+ +1 2X i,j,σi,σjZϕσi∗ i(r)ϕσj∗ j(r′)1 |r−r′|ϕσi i(r)ϕσj j(r′)d3r d3r′− −1 2X i,j,σ Zϕσ∗ i(r)ϕσ∗ j(r′)1 |r−r′|ϕσ j(r)ϕσ i(r′)d3r d3r′. (3.4) The second term is known as the direct interaction while the third is the exchange interaction. The effect of the latter is always to reduce the energy. The interpretation in this approximation is that each electron creates a "hole” due to other particles tending to be away from it because of the exclusion principle. The interaction of the negatively charged electron with the hole gives an overall negative contribution to the total energy. In the next sections, we will see that this kind of independent-particle approximation can be used in the DFT framework. 3.3 The Hohenberg-Kohn theorems and the Kohn-Sham method To begin with the formulation of DFT, we must first introduce the theoretical grounds provided by two important theorems, which we will now formulate without proof. Theorem 3.1. The ground-state density n0(r) of the system determines uniquely, up to a constant, the external potential, such as VIe , of any system of interacting particles. This means that the electronic Hamiltonian in Eq.(3.2) has a universal part comprised of the kinetic and interacting energy of the electrons, and a specific term for the interaction of the ions with the electrons. Remarkably, knowing n0(r) implies full knowledge of the Hamiltonian and, therefore, all the physical properties of the system. 40 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3.3. The Hohenberg-Kohn theorems and the Kohn-Sham method Theorem 3.2. For any external potential, we can define a universal functional E[n] , which depends on the density, whose global minimum is achieved exactly for the ground-state density E[n0].The functional has the form E[n] = T[n]+Eee[n] + ZVIe(r)n(r)d3r+EII,(3.5) where T and Eee are the kinetic and interaction energy of electrons and EII is the ion-ion repulsion energy. The second theorem provides a way to determine the ground-state density. Provided we found the functional E[n] , then n0(r) is computed simply by minimization of the energy. Moreover, E[n] depends only on the density, which certainly has less information that the true wave function of the system, as different states may lead to the same electronic density. Unfortunately, given the electronic density of a material, in practice there is no known way to extract all the physical properties, although the first theorem proves that it is indeed possible. For example, even evaluation of the kinetic energy T[n] is not possible using only the density, without giving an explicit form for the wave function from where it originates (e.g. a Slater determinant). The practical solution given by Kohn-Sham approach is precisely to state the problem in terms of a wave function that is a product of one-particle states, which allows us to evaluate all the terms in the functional. This wave function is not the true many-body state of the system, but that of an auxiliary, independent-particle system. The key insight is that, if we can design the functional of the auxiliary system such that it has the same ground state density as the true system, then the computed electronic density is also the exact density of the original Hamiltonian. Kohn-Sham equations The state of the auxiliary independent-particle system can be described by a Slater determinant as in Eq.(3.3), whose density is given by n(r) = X σ Nσ X i|ϕσ i(r)|2,(3.6) where ϕσ i are one-particle states of spin σ , with a total of Nσ particles per spin. The kinetic term can be explicitly evaluated as TKS =−1 2X σ,i ⟨ϕσ i|∇2|ϕσ i⟩=1 2X σ,i Z|∇ϕσ i(r)|2d3r, (3.7) 41 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3. Density functional theory where it is understood that there is an i -index summation for each σ . Although TKS is evaluated in terms of the wave function, it is in fact a functional of the density, as stated by the second Hohenberg-Kohn theorem. From the Hartree-Fock theory, we know that there is also a direct interaction term of the density with itself EHartree =1 2Zn(r)n(r′) |r−r′|d3r d3r′(3.8) and and exchange-correlation term Exc . It is precisely this last, unspecified term that is used to design the system such that it has the real ground-state energy, defining the Kohn-Sham functional as EKS[n] = TKS[n] + ZVIe(r)n(r)d3r+EHartree[n]+EII +Exc[n].(3.9) The formal expression of the exchange-correlation part is Exc[n] = ⟨ˆ T⟩−TKS[n]+⟨ˆ Vee⟩−EHartree[n],(3.10) accounting exactly for the difference between the real electronic energy ⟨ˆ T+ˆ Vee⟩ and the independent-particle energy. This guarantees that the energy and ground-state density of the Kohn-Sham system are exactly the real quantities. The expression for TKS is only known in terms of the wave function. However, using the chain rule and Eq.(3.6), we can minimize EKS[n]subject to the orthonormality constraints ⟨ϕσ i|ϕσ′ j⟩=δσσ′δij.(3.11) As a result, we obtain the Kohn-Sham equations Hσ KSϕσ i=−1 2∇2+Vσ KSϕσ i=ϵσ iϕσ i,(3.12) where the Kohn-Sham potential is defined as Vσ KS(r)=VIe(r) + δEHartree δnσ(r)+δExc δnσ(r).(3.13) The exchange-correlation term Although the Kohn-Sham method is exact, the problem of treating the interaction terms has been translated to the computation of the exchange-correlation functional, which is essential to reproduce the true electronic density. As a consequence, much of the effort in developing DFT-based methods is put into modelling of the exchangecorrelation term accurately and efficiently. One important aspect is that, given that 42 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 3.3. The Hohenberg-Kohn theorems and the Kohn-Sham method the long-range interaction terms are already accounted for in the Hartree term, Exc must have a local expression of the form Exc[n] = Zεxc([n],r)n(r),(3.14) such that Vσ xc ≡δExc δnσ(r)=εxc([n],r)+n(r)δεxc([n],r) δnσ(r).(3.15) Thus, the modeling of Exc is ultimately a problem of modeling εxc . Some of the most common approximations are: • Local-spin-density approximation (LSDA) [53]: they assume the exchangecorrelation energy of an homogeneous electron gas ELSDA xc [n↑, n↓] = Zn(r)ϵhom x(n↑(r), n↓(r))+εhom c(n↑(r), n↓(r)),(3.16) separating linearly the exchange and correlation parts. The advantage is that the exchange energy is known analytically εhom x(nσ) = −3 46 πnσ1/3 ,(3.17) while the correlation part has been fitted to Monte Carlo simulations. Surprisingly, despite its simplicity, LSDA can be rather accurate, especially for systems close to homogeneous electron gases such as nearly-free-electrons metals. • Generalized-gradient approximation (GGA): the exchange-correlation energy is expanded also in terms of derivatives of the density, which helps account for its non-homogeneity away from the constant LSDA approximation. This corrects its tendency to underestimate the exchange and overestimate the correlation. One common implementation is the Perdew-Burke-Enzerhof [54] (PBE) form. Beyond, the basic GGA, there exists the meta-GGA functionals, which further expand the energy in terms of second derivates of the density and other relevant quantities. One example is the modified Becke-Johnson potential [55,56] (mBJ), which takes the LSDA correlation potential and adds an exchange term which also depends on the local kinetic energy density. • Hybrid functionals: they combine the exact Hartree-Fock exchange with emperically fitted exchange-correlation terms and provide the most accurate approximation. However, the computational cost is also greatly increased. Some examples are the B3LYP [57,58], PBE0 [59,60] and M06 [61] functionals. 43 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation phase for every path λ(t) . Only if ⟨ψ′(λ)|∂µ|ψ′(λ)⟩ is expressed as the gradient of a scalar function can the phase be adjusted to vanish for any path, since the integral does not depend on the specific curve followed in parameter space. If λ traces a closed loop in parameter space, such that λ(0) = λ(t) , the single-valuedness of the wave function requires that ∆α=α(t)−α(0) = 2πm , where m∈Z . This can only be made zero for any loop if we can cancel the phase for any open path, since a loop can be decomposed of a series of open paths. This fact was popularized by Michael Berry in a paper published in 1984 [62] and it means that the geometric phase cannot be removed in closed loops, hence it was named the Berry phase in this context. He realized that the true dependence of the the eigenstates is on the parameters λ rather than just t . Whereas for the latter one can always express the integrand as a derivative of a function since it is a one-dimensional parameter space, for the former it is not always possible. 4.2 Gauge theory formulation The reader familiar with the topic may realize that this is in fact a gauge theory, where in this case the gauge group is U(1) , i.e., a phase. However, this is simplified with respect to other cases such as electromagnetism since the field degrees of freedom are static. The previous discussion can be reformulated in this language by noticing that the Berry phase can be re-written γ(t) = IAµdλµ(4.11) with Aµ≡ −i⟨ψ(λ)|∂µψ(λ)⟩ being the gauge potential, which is usually called Berry potential or Berry connection. Under a gauge transformation such as in Eq.(4.9), the gauge potential transforms with the familiar expression Aµ→Aµ+∂µα. (4.12) We can also define a gauge-covariant derivative Dµ≡∂µ−iAµ(4.13) such that under a gauge transformation D′ µ|ψ′⟩= (∂µ−iA′ µ)eiα|ψ⟩=eiαDµ|ψ⟩.(4.14) The field strength or Berry curvature associated to the connection is defined as Ωµν =∂µAν−∂νAµ−i[Aµ, Aν] = i[Dµ, Dν],(4.15) 50 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.3. Wilson loops: generalization of the Berry phase where [·,·] denotes the commutator, which vanishes in general when the gauge theory is Abelian, as in the U(1) Berry phase described above. However, this is not necessary. If instead of a single state |ψ(λ)⟩ we had tracked a subspace of degenerate eigenstates, the Berry phase would in fact be a matrix transforming these states into each other instead of a simple phase. Then the full expression of the curvature is important when dealing with non-Abelian gauge theories. For such cases, the transformation in Eq.(4.12) is not quite correct and must be replaced by the more general expression Aµ→U†AµU+iU†∂µU, (4.16) were U , in this context, is a unitary matrix that expresses the change of the subspace states after traversing the loop. It is easy to verify that if U=e−iα , everything in the expression commutes and Eq.(4.12) is recovered. Remaining now in the simplified U(1) version for the sake of clarity, suppose we λ(t) is restricted to a manifold of d dimensions and that it traces a loop C which encloses an area Σ. Then the Berry phase is the line integral γ(t) = IC Aµdλµ=ZΣ Ωµν dsµ∧dsν(4.17) where Stokes’ theorem was used to transform the line to a surface integral, assuming that it holds, i.e., that A is defined on λ(t) and Ω in Σ (it could be worked out in patches otherwise). dsµ∧dsν is a differential of the surface Σ in some coordinates. This is precisely the Berry phase we found in Section 4.1. To reflect on this, we have found that there is a freedom to choose any linear combination of degenerate states at each point along the path λ(t) to represent the physical state of the system. This freedom comes at the cost of having to define additional structure that expresses how to transform the states from one point to an infinitesimally close one, which is the connection A . We have also found that, after traversing a loop, the final state may differ from the initial one by a unitary transformation. This difference is precisely the flux of the Berry curvature through the area enclosed by the curve. In fact, Ωgives the infinitesimal contribution to this difference over a plaquette of differential area such that the total change is the integral over the surface. 4.3 Wilson loops: generalization of the Berry phase The derivation in Section 4.1 is a simplified version of a more general concept called Wilson loop. In that procedure, we tracked the change of |ψ(λ)⟩ such that, at every 51 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation point in parameter space, it was both normalized and perpendicular to states with other energies. This can be viewed as a form of parallel transport in the sense we identify the orthonormality condition with the fact that the “angle” of the vector with itself along a line tangent to the curve λ(t) remains constant. The parallel transport condition can be expressed as the vanishing of the covariant derivative along the curve λ(t) λµDµ|ψ(λ)⟩=λµ(∂µ|ψ(λ)⟩−iAµ|ψ(λ)⟩) = 0 (4.18) and taking the inner product with ⟨ψ(λ)|we just find again that ⟨ψ(λ)|∂µψ(λ)⟩=iAµ.(4.19) Therefore, what we found in the derivation of the geometric phase is that, despite translating |ψ(λ)⟩ in a parallel fashion along the curve, the final and initial states may differ by a phase. This idea is generalized by considering a closed subspace of eigenstates of H(t){|ψn(λ)⟩} such that we expand the instantaneous state of the system as |χ(λ)⟩=cn(λ)|ψn(λ)⟩,(4.20) where we assume summation over repeated indices. The parallel transport condition reads Dµcm=∂µcm−i(Aµ)mncn= 0,(4.21) noting that Aµ is matrix-valued. The transformation from the initial to final coefficients is given by the formal expression cm(λ(t)) = hPeiRλ(t) λ(0) Aµdλµimn cn(λ(0)),(4.22) where we identify the Wilson line operator W(λ)mn =hPeiRλ(t) λ(0) Aµdλµimn .(4.23) The symbol P means that the integration is path-ordered in the direction of the curve. We emphasize that the spectrum of W(λ) depends both on the path λ(t) and the connection A , which means that it is not gauge invariant. It can be shown that, under a unitary transformation (a gauge) of the basis states |˜ ψn(λ)⟩=Umn|ψm(λ)⟩ , W(λ) transforms as ˜ Wmn(λ)=U†(λ(t))miWij(λ)U(λ(0))jn.(4.24) Thus, we see that in the special case that the line is closed such that λ(0) = λ(t) , Eq.(4.24) is actually a similarity transformation, which means that both the eigenvalues and the trace are gauge invariant. This specific case of the Wilson line operator 52 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.3. Wilson loops: generalization of the Berry phase is called a Wilson loop 1 and is a fundamental object to study the topology that will be used throughout this work. Wilson loops have some interesting properties that we enumerate here: • Considering all the loops that start from a point λ0 , Wilson loops along all these paths form a group called the holonomy group Φλ0(M, A) . The transformation corresponding to a compositions of loops is the product of matrices of each path. The inverse transformation can be found by traversing the loop in the opposite direction, while the identity matrix corresponds to the constant loop (no loop at all). • The choice of base point along the path is irrelevant for the spectrum. To see this, consider loops α,β with base points λα,λβ respectively as in Fig.4.1. Assuming the path-connectedness of the parameter space, the loop α can be composed by going along a path from λα to λβ , performing the β loop and the going back from λβto λα. The corresponding operator is Wα=Wλα→λβWβWλβ→λα=W† λβ→λαWβWλβ→λα(4.25) and since the matrices are unitary, this is a similarity transformation. • The Wilson loops over paths that are smoothly deformable into each other (homotopic) are in general different. However, they are equal if the connection is flat, so the curvature vanishes everywhere. In fact, the integral of the curvature over the area enclosed by two homotopic loops measures the difference in parallel transport along both paths. All in all, gauge freedom allows us to define a basis at each point in parameter space and the connection is the object that expresses how the basis changes from one point to another. When parallel-transported, the precise way in which the different spaces at each point in parameter space are glued together may result in states returning to a transformed version of themselves when tracing a loop. This change is precisely the holonomy and it is computed by the Wilson loop operator. Additionally, the set of such transformations found by considering all loops through a base point forms the holonomy group. Why Wilson loops probe the topology So far, all the concepts we have introduced such as parallel transport, connection and curvature, are purely geometrical. However, the holonomy introduces a way to relate 1 The term “Wilson loop” may sometimes be used to refer to the trace of the Wilson line operator over a closed path. 53 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation Figure 4.1: Sketch for the Wilson loop construction that shows its independence from the base point. The loop starting and ending at λα is equivalent to going from λα to λβ , performing the loop at λβ and going back to λα . This is mathematically equivalent to a similarity transformation and, therefore, the spectrum does not change. these ideas to the topology of the underlying manifold. As we will see in following sections, the final link is given by relating the holonomy to the curvature. In this section, we will first see how Wilson loops capture the topology of a system. First, we need to introduce the concept of homotopy. Consider a manifold 2M . A path on M starting at x0 and ending at x1 is just a map λ from the unit interval I≡[0,1] to points in M such that λ(0) = x0 and λ(1) = x1 . A loop is path whose start and end points are the same. A homotopy is a continuous map F:I×I→M that maps loops α(s), β(s)with the same base point x0to each other such that F(s, 0) = α(s), F(s, 1) = β(s)and F(0, t)=F(1, t)=x0.(4.26) Two loops are therefore homotopic if there exist such F relating them. With this notion in mind, the fundamental group π1(M, x0) of loops with base point x0 is formed by all loops that start and end at x0 where paths that are homotopic are considered the same. π1(M, x0) knows about the topology of M because there may be non-contractible loops on M that cannot be deformed to neither the identity or constant loop (because they are non contractible) nor to other loops, which may happen when there are holes in M . As a typical example, the fundamental group of the sphere S2 is trivial, since all closed paths can be deformed into each other without leaving S2 . However, π1 of the torus S1×S1 is not trivial and in fact π1(T2) = Z×Z , where we give one integer for the winding along one direction enclosing the “handle” of the torus and another one for the other direction that goes around the hole, which cannot be deformed into each other (see Fig.(4.2). Note that in these examples, the 2In fact, it need not be a differentiable manifold, but only a topological space 54 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.4. The Chern number: a topological invariant Figure 4.2: The two kinds of non-contractible loops in a 2-torus, with two examples of the homotopic families in red and blue. The fundamental group π1(T2) = Z×Z , corresponding to the number of times a loop goes around each handle. Because T2 is path-connected, π1(T2)is independent of the base point. base point of π1 is irrelevant, which happens for every manifold where any two points can be connected by a continuous path. We assume that this is the case from here onward. The connection of π1 to the holonomy group, comes from the fact that Wilson loops are computed precisely over closed paths such as those in the fundamental group. In fact, there is a map between the loops in π1(M) and the Wilson loops of the holonomy group Φ(M, A) such that to each loop in π1(M) we assign the Wilson loop operator along it. This map is in general not surjective because Φ(M, A) is larger than π1 because Wilson loops along homotopic paths are not equal unless there is no curvature. From the homotopies that identify equivalent loops, we also obtain a continuous map between Wilson loops whose paths are homotopic. It can be shown that Wilson loops along paths that can be deformed to the constant path are connected via this continuous deformation to the identity, while the ones that cannot be deformed in this way are completely disconnected. Apart from the relationship with the fundamental group, the following section shows that the Wilson loops can also be used to compute topological invariants, adding to the connection between geometry and topology. 4.4 The Chern number: a topological invariant The topology of manifolds can be distinguished by computation of topological indices 3 . In this section, we introduce perhaps the most widely known of them in the context of topological materials called the Chern number and show how it can be computed from the Berry curvature tensor. Moreover, we will also see that there 3 Manifolds with different indices are topologically inequivalent. However, two manifolds with the same set topological indices cannot be said to be equivalent, since there could be other unknown indices that distinguish them. 55 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation α λ A1 A2 Figure 4.3: Sketch of the two level system. Since |λ| is fixed, the parameter manifold is a sphere S2 . Since eigenstates have vortex singularities either at the north or the south pole, two gauges must be chosen corresponding to two different connections A1 and A2 . These are related by a gauge transformation given by α , whose winding number is the Chern number. is a link between Wilson loops and the curvature which completes the relationship between topology and geometry. First, suppose that the parameter space is a sphere M=S2 . This exercise will cover the most common example of topology, which is the two level system, whose Hamiltonian in some basis can be described by a two-dimensional matrix H(λ)=aσ0+λ·σ,(4.27) with energies E±=a±|λ|,(4.28) where λ is a set of parameters which may evolve adiabatically with time and σ are the Pauli matrices assembled into a vector, σ0 being the identity. For the sake of our discussion, the overall shift in the energy levels given by aσ0 , a being a real number, can be neglected. Furthermore, the modulus of λis fixed so the only real degrees of freedom if we take standard spherical coordinates are two angles. This means that λ is restricted to a sphere S2 as we have introduced. This general Hamiltonian describes the energy levels of a spin1/2 particle in a magnetic field B if we take λ≡B , for example. Consider now the integral over the sphere of the Berry curvature ZS2 Ωµν dsµ∧dsν(4.29) 56 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.4. The Chern number: a topological invariant where dsµ∧dsν is an infinitesimal area in some coordinates. We could conclude naively that this integral vanishes since, by application of Stokes’ theorem, the expression can be transformed into the integral of A over the boundary ∂S2 , but the sphere has no boundary (we should not be confused by the picture of embedding of a sphere in three-dimensional space). However, we will see that Stokes’ theorem cannot be blindly applied here. Notice first that if |λ| = 0 there are two eigenstates separated by a spectral gap. This means that both eigenstates are not mixed as long as the path traced by λ is restricted to S2 , so we can track, for example, the evolution of the lower energy state |ψ0(λ)⟩ . The connection A is computed from |ψ0(λ)⟩ and it is a well known fact that a smooth vector field cannot be defined over all S2 , which implies that we have to work by patches, for example dividing it into hemispheres. This forces us to use at least to different connections A1 and A2 , which are related in by a gauge transformation implemented on the states by eiα(λ) on the area where the two patches overlap A2 µ=A1 µ−∂µα(4.30) where it is most convenient in this case to work in spherical coordinates (i.e. µ refers to θ and ϕ angles). The actual form of the eigenstates for the Hamiltonian in Eq.(4.27) is |ψ0(λ)⟩=cos(θ/2) sin(θ/2)eiϕ,|ψ1(λ)⟩=sin(θ/2) −cos(θ/2)eiϕ,(4.31) where λ=|λ|(sin θcos ϕ, sin θsin ϕ, cos θ) and the eigenvalues are E0=−|λ| and E1=|λ| . We see that |ψ0(λ)⟩ has a phase singularity at the south pole where θ=π . However, this can be undone by a gauge e−iϕ such that the resulting vector is illdefined at the north pole. Proceeding now to compute the integral in Eq.(4.29), we find ZS2 Ωµν dsµ∧dsν=ZM1 Ωµν dsµ 1∧dsν 1+ZM2 Ωµν dsµ 2∧dsν 2(4.32) =Z∂M1 A1 µdλµ−Z∂M2 A2 µdλµ(4.33) =IL ∂µα dλµ.(4.34) where L denotes the curve that delimits the two hemispheres. In going from the first to the second line, we have to take into account that the two hemispheres M1 and M2 have the opposite orientation so the curves must be traversed in opposite directions, giving rise to the minus sign. Remembering that the single-valuedness of the wave function |ψ0(λ)⟩ requires that the gauge transformation be periodic, we finally arrive at ZL ∂µα dλµ= 2πC (4.35) 57 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation where C is known as the Chern number 4 . It is an invariant of the state |ψ0(λ)⟩ that has the very interesting property of being an integer. This means that, as its value cannot continuously change, extreme changes in the structure of the eigenstate as a function of λ must take place in order for it to vary, a feature which is general for topological invariants. In particular, unless the gap between |ψ0(λ⟩) and the higher-energy state is closed, the Chern number will remain unchanged, no matter the actual values of the energy levels over S2 . In passing, we have also seen that C is in fact related to the winding number of the gauge transformation over the sphere. When generalizing from the two-level system, we could track more than one eigenstate such that Ais non-Abelian. Then Eq.(4.29) is rewritten for a manifold M ZM Tr(Ωµν)dsµ∧dsν,(4.36) the gauge transformation between connections in two patches is A1 µ=U†A2 µU+iU†∂µU(4.37) and the Chern number is IL ∂µχ dλµ,(4.38) where we define χ≡Im log det U . Taking into account that U is unitary so that its eigenvalues are complex numbers of unit modulus i(U†)ij(∂µU)jk =ie−iϕiδij(i)(∂µϕj)eiϕjδjk =−∂µϕiδik (4.39) so χ is just the sum of the angles of the phase eigenvalues of U5 . Note also that taking the trace of Ωremoves the commutator term in Eq.(4.15). Going back to the two level system, we can compute the Berry connection from the eigenstate |ψ0(λ)⟩in Eq.(4.27) Aθ=−i⟨ψ0(λ)|∂θ|ψ0(λ)⟩= 0,(4.40) Aϕ==−i⟨ψ0(λ)|∂ϕ|ψ0(λ)⟩= sin2θ 2,(4.41) which is valid for any patch on the sphere that does not include θ=π . The only non-zero independent component of the Berry curvature is then Ωθϕ =−Ωϕθ =∂θAϕ−∂ϕAθ=∂θAϕ= sin θ 2cos θ 2.(4.42) 4 This is the first Chern number. There could be others related to integrals of powers of Ω . However, as a differential form, Ω is already of order two, which means that higher powers are already of dimension higher that S2and hence vanish. 5We have used a basis where Uis diagonal without loss of generality. 58 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.4. The Chern number: a topological invariant Performing the integral over the whole sphere Z2π 0 dϕ Zπ 0 dθ sin θ 2cos θ 2= 2π−1 2cos θϕ 0 = 2π, (4.43) which means that the Chern number is C= 1 . One can also show that, in Cartesian coordinates, the curvature has only a radial component of the form Ωr∝1 r,(4.44) which is similar to the expression for the magnetic field of a hypothetical magnetic monopole. Relating the Chern number to the holonomy In Section4.3 we saw how the holonomy expresses the change of the eigenstates when they are parallel-transported along closed paths. We have also defined the Chern number as the integral of the Berry curvature over the parameter space and we have seen how it probes the non-triviality of the Berry connection defined in that same space. We now show that the relation between the holonomy of a connection and the curvature relates the Chern number with the Wilson loop operators. The missing link is provided by the Ambrose-Singer theorem. In a broad sense, it tells us that the information of the holonomy group at a point λ0 can be found in the curvature at that same point. In other words, Ω gives the differential contribution to the holonomy of an infinitesimally small loop C at λ0 . Given the form of the Wilson line (or loop) in Eq.(4.23), by expanding the path-ordered exponential to second order [63], one can find that W(C) = Pexp iIC Aµdλµ≈1+iPIC Aµdλµ+ +i2 2PIC Aµdλµ2 +··· (4.45) The first order term can be expressed by Stokes’ theorem as PIC Aµdλµ=ZZΣ 1 2(∂µAν−∂νAµ)dλµ∧dλν,(4.46) where Σis the area delimited by C, and the second order term is 1 2PIC Aµdλµ2 =−1 2ZZΣ [Aµ, Aν]dλµ∧dλν+··· (4.47) 59 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation Therefore, the spectrum of the Wilson loop operator is related to that of the position operator 8 . Furthermore, the center of the n -th Wannier function ⟨wn0|x|wn0⟩=xn is the re-scaled logarithm of the n -th diagonal element of W , provided that {|˜ ψnk⟩} were chosen to be maximally localized. Naturally, the spectrum of the Wilson loop does change under gauge transformations and so do the Wannier centers. However, the trace is gauge-invariant so the sum of centers also remains unchanged, at least up to a lattice vector translation. This indeterminacy comes from the invariance of the Wilson loop operator when an additional 2π phase is included. Thus, there are some gauge transformations that shift the lattice of Wannier centers by some vector, but the lattice itself remains the same since it repeats indefinitely. Therefore, we can regard this sum of centers as the actual centers around which the charge of the electrons in the selected bands is centered. The above interpretation can be generalized to 3D cases by considering hybrid Wannier functions. These are a set of states which are localized in only one direction r⊥ obtained by Fourier transforming only in the associated reciprocal space direction k⊥: |wnR⊥,k∥⟩=a 2πZ2π/a 0 e−ik⊥R⊥|˜ ψnk⊥,k∥⟩dk⊥,(4.78) where the cell has length a in the r⊥ direction and k∥ are the two remaining momenta in the periodic directions. These states can be visualized as "sheets” that extend over r∥ but have a definite center in r⊥ . Just as in the 1D case, the Wilson loop computed by integrating over k⊥ provides the position of these positions. The set of hybrid Wannier functions is in fact a set of eigenstates of the position operator r⊥ projected onto the set of bands if they are maximally localized in r⊥ . However, in general it is not possible to find a set of functions that diagonalizes two position operators at the same time, since they do not commute. It is only possible when the Berry curvature vanishes at all points of the BZ. Hybrid Wannier functions and the Chern number We can finally provide a physical interpretation for the Chern number in this picture. Assume for simplicity a 2D system for which we have constructed hybrid Wannier functions {|wnR⊥,k∥⟩} which are smooth in k⊥ . We can compute Wilson loops W(k∥) by integrating over k⊥ which are dependent on the remaining k∥ and, by the previous explanations, these give the positions of the hybrid Wannier charge centers. 8 More precisely, it is related to the spectrum of the position operator projected onto the set of bands of interest. 66 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.5. Application to crystalline systems of electrons Furthermore, we have seen that tracking the evolution of W(k∥) over closed loops is equivalent to computing the Chern number. For the Chern number to be non-zero, the necessary but not sufficient condition is that time-reversal is broken. This is because upon time-reversal, the curvature transforms as Ω(k)TR −→ −Ω(−k)(4.79) which means that the integral over a closed manifold must vanish by cancellation of opposite wave vectors. In passing, we mention that inversion transforms Ω(k)inversion −−−−→ Ω(−k),(4.80) which implies that when both TR and inversion are present, the curvature must vanish everywhere. A non-zero Chern number C⊥ computed by integrating along k⊥ implies that the lattice of hybrid centers is shifted by C⊥ unit cells in the R⊥ direction upon returning to the starting point. In this sense, it is interpreted as charge being pumped in that direction. This also implies that there is an obstruction to computing exponentially-localized Wannier functions in all directions, since they are not periodic in the k∥ . We also saw in the two-level system that the Chern number is related to the winding of a gauge transformation, which is necessary when a global gauge cannot be found. Therefore, C⊥= 0 prevents us from finding a smooth gauge in all directions, so the Fourier transform gives badly localized functions as a result. Quantum Anomalous Hall insulators Insulators with non-zero Chern number are known as quantum anomalous Hall insulators (QAHI). They are the hallmark of non trivial topology in condensed matter systems, providing a manifestation of these features in the form a Hall current even in the absence of an external magnetic field. This can happen for example when the crystal presents an arrangement of magnetic moments obeying the symmetry of a Shubnikov space group without TRS, as the breaking of this symmetry is required. For simplicity, consider first the case of a 2D material. Suppose that we have constructed hybrid Wannier functions in the direction r1 and that the Wilson loop winds over a loop along k2 by 2π . This implies that the Chern number C= 1 and that the Wannier centers are shifted by a lattice vector R1 , which means that there is charge flowing in that direction. A consequence of this is that, when the system is cut perpendicular to R1 , there should be charge accumulating in the surface. However, we know that this is not possible since the charge in the surface does not change. After all, the Hamiltonian is the same at the base point of the loop so the contribution of the bands to the surface charge cannot vary. 67 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation The paradox is resolved by the presence of a pair of states, one at each surface, in the gap of the insulator. While one of the states injects one charge in the surface from the bulk, the other one must extract the same charge from the surface to the bulk. Moreover, these states are chiral: they move in opposite direction. This is required by the breaking of time-reversal, as the TR-partner state moving in the opposite direction cannot exist on the same surface. These states are responsible for the quantized Hall conductivity σQHAI =Ce2 h,(4.81) where h is the Planck constant. This is an example of the bulk-boundary correspondence: the properties of the bulk material remarkably require features on the surface. Systems in 3D are actually characterized by a triplet of Chern numbers (C1, C2, C3) , corresponding to the evolution of hybrid Wannier charge centers in each of the directions in real space. They can be regarded as weakly coupled layers of 2D QAHIs stacked in the direction of the Chern triplet. For example, a Chern insulator with (0,0, C) displays a Hall conductivity as in Eq.(4.81) when made finite in (R1, R2) . In this picture, the edge states of every 2D layer combine to form states that are extended over the surfaces of the full crystal. To conclude, we provide a heuristic explanation of the bulk-boundary theorem. Regarding vacuum (or air, for all intents and purposes) as a trivial insulator with vanishing Chern number, the interface between a Chern insulator in such environment must present some gap closing. This is because topological invariants can only vary upon closure and re-opening of an energy gap. Since when moving from the topological insulator to the vacuum the Chern number varies, there must be gapless edge states that are located precisely in the gap. Quantum Spin Hall insulators While Chern insulators cannot exist when TRS persists, it turns out that there another classification for TR-invariant insulators with a different topological index. Consider a Hamiltonian H of a 2D system of spinful electrons. H can be decomposed into spinup and down parts, with coupling terms between the two spin values corresponding to SOC H=H↑∆† ∆H↓.(4.82) Suppose now that the system is a combinations of two subsystems, one for each spin orientation, such that each subsystem of spinless fermions has a Chern number C↑,↓ 68 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.5. Application to crystalline systems of electrons E R R R a) b) c) Figure 4.5: Sketch of bands in a TR-invariant, one-dimensional insulator between two TRIMs Γ and R.. a) Quantum spin Hall insulator b) Trivial insulator without SOC c) Same trivial insulator with SOC. The blue and green regions represent the valence and conduction bands, respectively. The arrows over the bands indicate the spin component. when they are decoupled. As a whole, H , with coupled spins, must be TR-invariant since TRS interchanges H↓ and H↑ while also transforming spin-up electrons to spindown. Therefore, we must have that CT=C↑+C↓= 0 , such that C↑=−C↓ . This means that no net charge can flow on the boundary of the 2D system, since spin-up and down electrons move in opposite directions due to the opposite Chern number. However, we say that there is a spin current along the edge, corresponding to what we call the quantum spin Hall effect (QSHE), which is quantized and protected unless there is a gap closing because C↑,↓cannot change continuously. To define the topological invariant identifying the quantum spin Hall insulators (QSHIs), we need to take into account the action of TRS, denoted Θ , on Bloch states. Recall from Chapter 2 that the Bloch Hamiltonian satisfies ΘHkΘ−1=H−k(4.83) and that the action fo TRS on states is Θ|ψ↑ k⟩=|ψ↓ −k⟩, Θ|ψ↓ k⟩=−|ψ↑ −k⟩.(4.84) From this, we see that TR-partner states are degenerate at TRIMs, where k is equivalent to −k , even with SOC. This enforces two types of connectivity of the bands between two TRIMs, as sketched in Figure 4.5 for a 1D insulator. In any case, TRS requires that every state at any TRIM is doubly degenerate, forcing a connection of 69 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation two bands at those points. Figure 4.5a represents a QSHI, where a TR-preserving perturbation can move the pairs up and down into the conduction and valence bands (in green and blue), but a single in-gap state always remains. The remaining TR-partner edge state is along the path related to Γ−R by TRS. In Figure 4.5b, an insulator without SOC and C↓=−1, C↑= 1 is depicted. Notice that spin-up (red arrows) and down (blue) edge states cross without coupling because there is no SOC and Sz is a good quantum number. Finally, in Figure 4.5c, we activate SOC and the two spin orientations are coupled. However, this is a trivial insulator because perturbations that respect TRS, which keep states paired up at TRIMs, can push the bands to the conduction and valence regions such that there is no in-gap state left. QSHIs are therefore described by a Z2 index, since adding a second pair of edge states trivializes the topology due to hybridization of the spins via SOC, as opposed to Chern insulators that can display an unbounded number of modes and are described by a Z invariant. A visual way to compute the Z2 number is to trace an horizontal line from Γ to R and count the number of edge states it crosses. If there is any line that cuts an even number of bands, then the insulator is trivial. For inversion-symmetric crystals, Fu and Kane [12] found that the topological index can be computed by knowing the parity ξm of the Bloch states |um(k)⟩ for every band m at all the TRIMs Λaas (−1)ν=Y aY m ξm(Λa).(4.85) The same discussion we have used for the Hamiltonian can be done for the Wilson loop operators. In particular, the ν index can be computed by calculating the flow of Hybrid Wannier functions of the occupied bands of the 2D system over half the BZ and counting the number of WL lines crossed by an imaginary horizontal line. In 3D crystals, Fu, Kane and Mele [65] demonstrated that there are four invariants describing the topological properties of the QSHIs, which we arrange in the form of a four-tuple (ν0,;ν1, ν2, ν3) . The six ki= 0,1/2 planes in the BZ are effectively TR-invariant systems that can be assigned a Z2 invariant. It turns out that only four of them are independent due to Kramers’ degeneracy. The corresponding formulas of inversion eigenvalues are (−1)ν0=Y aY m ξm(Λm), (−1)νi=Y bY m ξm(Λb(ki= 1/2)),(4.86) where Λb(ki= 1/2) are the TRIM points on the plane ki= 1/2 . Equivalently, one can compute the flow of Wannier centers on the planes ki= 0,1/2 and do the horizontal 70 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.5. Application to crystalline systems of electrons line construction to obtain the indices ν′ i, νi , respectively. This results in six indices that are related as (−1)ν0= (−1)ν1+ν′ 1= (−1)ν2+ν′ 2= (−1)ν3+ν′ 3,(4.87) which forces νi, ν′ i to be different when ν0= 1 and equal if ν0= 0 . Any insulator with at least one non-zero index is considered non-trivial. When ν0= 1 , regardless of the values of νi , the QSHI is called a strong topological insulator, while crystals with at least one νi= 1 and ν0= 0 are called weak. Weak topological insulators can be considered as a series of 2D quantum spin Hall layers stacked in the direction of the (ν1, ν2, ν3) triplet. Correspondingly, these show surface states only on select facets. For example, a crystal with triplet (0,0,1) has boundary states on the x and y surfaces. There is no such layer construction for strong QSHIs, which show surface states on all the faces. Finding Wannier functions from DFT calculations The computation of a maximally localized Wannier basis from DFT can be carried out in a straightforward way thanks to the Wannier90 [66] software. Here, we will briefly describe the procedure followed by Wannier90, with a focus on the practical points when performing our own calculations. The role of the ab-initio is to compute the Bloch eigenstates |ψmk⟩ . However, these eigenstates will not be optimally smooth in k , for which we have to find a k dependent unitary transformation that will iron out the gauge throughout the BZ, obtaining a new set |˜ ψnk⟩=Uk mn|ψmk⟩.(4.88) This optimally smooth set is the Fourier transformed to obtain the maximally localized Wannier basis. The purpose of Wannier90 is to find the Uk transformation by minimizing the spread σof the WFs, defined as σ2=X n⟨wn0(r)|r2|wn0(r)⟩−|⟨wn0(r)|r|wn0(r)⟩|2.(4.89) The spread has a gauge-invariant component σ2 I and another component ˜σ2 which can be reduced by a suitable choice of Uk σ2=σ2 I+ ˜σ2=σ2 I+σ2 D+σ2 OD,(4.90) where σ2 D and σ2 OD are the “diagonal” and “off-diagonal” components. In a realistic DFT calculation, we may not find an isolated set of bands or, if it exists, this set may be formed by a large number of states. Wannier90 implements for this cases the so called 71 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4. Topology in condensed matter systems and its geometrical interpretation “disentanglement procedure”. The method requires that we set an energy window to select the states we want to include in the calculation for a total of Nk states at each k point. If we want to obtain N WFs, the program will find a rectangular Nk×N matrix Uk dis that minimizes σ2 Iinside the window, obtaining a set of optimal states |˜unk⟩= Nk X m=1 Uk dis|umk⟩.(4.91) Furthermore, we can also select an inner “frozen window” where the states are included as-is, which helps accurately reproduce the band structure inside that energy range and is usually set to be narrow and centered around the Fermi level. In practice, what we need to provide is a set of trial localized functions |ϕn⟩ . When no better chemical intuition is available (e.g. hybrid orbitals), the way to select these orbitals is to include the atomic orbitals with the largest contribution around the Fermi level. From this, Wannier90 will compute •Projection of the Bloch eigenstates onto the localized functions Ak mn =⟨ψmk|ϕn⟩.(4.92) •An overlap matrix of the cell-periodic part of the Bloch states Mk,b mn =⟨umk|unk+G⟩.(4.93) Apart from this, we will select the outer and frozen energy windows and a sufficient number of steps both in the disentanglement and gauge-smoothing procedures. It is usually beneficial to converge an additional set of bands above and below the selected set and completely disregard states which are very far isolated and very far in energy. Computing the Berry phase numerically The computation of the Wannier basis allows us to express the Hamiltonian in terms of transition amplitudes between localized orbitals, that is, a tight-binding (TB) model. In practical calculations of Berry phases and curvatures, we will use this kind of approximation since it is computationally more tractable. This is because we need a finer k -grid to obtain a version of the Berry phase using a discretized version in reciprocal space. For simplicity, consider first the case of an isolated band with cell-periodic eigenstates |unki⟩ along a closed path C discretized in N wave vectors ki . For two consecutive points which are sufficiently close Im ⟨umki|unki+1 ⟩(4.94) 72 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 4.5. Application to crystalline systems of electrons gives a notion of the phase difference between the states at ki and ki+1 . Therefore, we can define the Berry phase along Cbu the following product of overlaps ϕ=−Im ln ⟨unk0|unk1⟩⟨unk1|unk2⟩···⟨unkN−1|unk0⟩,(4.95) since Im ln reiϕ=ϕ . In the multi-band case, we have to calculate the overlap matrix Mki,ki+1 mn =⟨umki|unki+1 (4.96) and we obtain the following formula for the traced holonomy ϕ=−Im ln Y i det Mki,ki+1 .(4.97) Finally, to approximate the Berry curvature, we compute the eigenstates in a sufficiently finely-grained k -grid. Given nature of the curvature as the local contribution to the holonomy, we can calculate it by computing the Berry phase along small loops in the grid and dividing it by the area of the loop. 73 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 CHAPTER 5 Topology from real-space symmetry: Magnetic Topological Quantum Chemistry 75 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package In Chapter 2 we introduced the structure and representation theory of the Shubnikov space groups describing the symmetry of crystalline solids. Considering its application to the electronic band structure, one of the key takeaways was that eigenstates of the Hamiltonian at a given k point transform as an irreducible representation of the little group of k . Moreover, a d -dimensional irrep corresponds to a d -fold set of states with the same degenerate energy. In Chapter 5, we showed that an isolated set of bands can be represented by a vector containing the multiplicity of each irrep at every high-symmetry k point in the BZ from the point of view of an analysis of symmetry-indicated topology. This remarkable fact implies that, in many cases, we can characterize the topology of system by using algebraic methods instead of more demanding calculations such as Wilson loops (as in Chapter 4). Given the above and considering that numerical methods, in particular DFT, are necessary to compute realistic materials, any way to automatically identify the symmetry content and the subsequent topological analysis of electronic band structures obtained from ab-initio methods is extremely useful. Moreover, the problem is naturally well suited for a computer implementation since, as we will see, we only have to compute the trace of symmetry operations over subspaces corresponding to the same energy, which is in essence a tensor operation. Automation of the task, however, requires a big effort by itself. In particular, an application which is easy to use and general for any crystal has to take into account that the DFT calculation can be done in an arbitrary unit cell, while the symmetry tables are usually expressed in the conventional setting. In this chapter, we present a tool that streamlines this process in a user-friendly way and which we have enhanced to work with MSGs and analyze the topology in terms of EBR decomposition and symmetry indicators. 6.1 Computing the symmetry of eigenstates from DFT calculations irrep [51] is a Python package whose main purpose is to calculate the symmetry content of the electronic bands computed by a DFT plane-wave software in a simple and easy-to-use way. One of its strengths in comparison to similar programs is that the correct symmetry SG is automatically computed given the input structure for the calculation. Other solutions require the DFT calculation be performed in a particular setting which agrees with the standard unit cell in the some SG irrep database. However, irrep can be used in any setting and all the symmetry operation transformations are computed automatically to match them with the symmetry tables of the Bilbao Crystallographic Server (BCS) [52]. The consistent use of the BCS notation 82 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6.1. Computing the symmetry of eigenstates from DFT calculations also prevents any confusion regarding the labeling of the irreps, denomination of the space groups and unit-cell setting, for example. Another key feature is that irrep can be interfaced with some of the most popular plane-wave DFT suites, such as VASP, Abinit and QuantumEPRESSO, as well as Wannier90 for the analysis of interpolated tight-binding models in a Wannier basis calculated from ab-initio methods. However, irrep was previously unable to account for the symmetry properties of crystals displaying magnetic ordering. This has limited its application for the analysis of ferromagnetic, anti-ferromagnetic and ferrimagnetic systems, which have been shown to be a promising platform for the study of topological phases [67]. Moreover, many compounds with elements containing valence d and f orbitals display some kind of ground-state, spontaneous magnetization [68 – 70]. This leads to an more complex classification of the crystal symmetry in 1621 MSGs, compared to the 230 regular, TR-invariant SGs. One can also argue that the correct analysis of paramagnetic phases with SOC is best done using MSG theory, since TRS can be better understood in the context of gray or Type-I MSGs. irrep is packed with all the symmetry tables for high-symmetry points required for a MTQC analysis for all MSGs obtained from the BCS [25,26], enabling the local analysis without any external dependencies. Furthermore, previous iterations of the program where limited to only the computation of the little-group irreps at high-symmetry points and only limited topological analysis was available. This update introduces the capability to compute all the symmetry indicators for all MSGs, given in the so-called physical basis and thus with well known interpretation. Additionally, one can check the topological classification of a given set of electronic bands, computing its decomposition into EBRs [15] when the system displays trivial of fragile topology. Working with a plane-wave basis set Many DFT suites such as VASP, Abinit and QuantumESPRESSO use a plane wave basis to represent the electronic wave functions. In particular, the eigenstates 1|ψnk⟩ of H(k) at a given reciprocal space point are expanded in a series of plane waves |k+G⟩|ψnk⟩=X G Cnk(G)|k+G⟩(6.1) where we impose a G cutoff to a maximum plane-wave energy of Ecut such that ℏ2(k+G)2/2m < Ecut . The actual cutoff used by our program need not be the 1 When using PAW pseudo-potentials, we obtain the expansion for the pseudo wave functions representing the states outside the core. However, since they are related to the true wave function by a linear transformation that commutes with the symmetry operations, their transformation properties are the same. 83 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package calculation cutoff. In general, Ecut ∼50 eV is enough to calculate the traces of symmetry operations without errors. To express the transformation of the eigenstates using the plane-wave coefficients, we note that an operation g={R|t} acts on |k+G⟩as g|k+G⟩=e−iR(k+G)·t|R(k+G)⟩.(6.2) It follows that the trace can be computed as ⟨ψnk|g|ψnk⟩=X G,G′ C∗ nk(G′)Cnk(G)⟨k+G′|g|k+G⟩ =X G,G′ C∗ nk(G′)Cnk(G)e−iR(k+G)·t⟨k+G′|g|R(k+G)⟩ =X G,G′ C∗ nk(G′)Cnk(G)e−iR(k+G)·tδG′,Rk−k+G =X G,G′ C∗ nk(Rk−k+G)Cnk(G)e−iR(k+G)·t, (6.3) where we have used the orthonormality condition ⟨k′+G′|k+G⟩=δG′,k−k′+G.(6.4) When using spinor wave functions |ψσ nk⟩ where σ is the spin component, an additional transformation of the spin part of the eigenstate must be taken into account ⟨ψnk|g|ψnk⟩=X σ,σ′⟨ψσ nk|g|ψσ′ nk⟩.(6.5) The SU(2) matrix S(g) can be obtained by the mapping from the real-space operation S({R|t}) = 1 2cos(θ 2)−i(ˆn ·σ) sin(θ 2),(6.6) where 1 2 is the 2×2 identity matrix, σ is the vector of Pauli matrices, ˆn is the axis of R and θ its angle. Therefore, for every k point, the program groups eigenstates by identical (up to a threshold) energy and computes the trace of all unitary coset representatives in the little group of k and the corresponding irreps can be computed by using Equation 2.66. Handling changes of basis between different unit cells One of the strengths of irrep is that it allows the user to choose any unit cell for the calculation of electronic eigenstates. This is especially useful since information about the symmetry operations and representations of space groups is conventionally 84 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6.1. Computing the symmetry of eigenstates from DFT calculations a1 a2 a3 p x' x c1 c2 c3 Figure 6.1: Sketch of two axes settings related by a matrix P and shift of origin p . The same point is expressed in different coordinates x′in the cisetting and xin ai. defined with axes forming a non-primitive cell, while the calculations are usually done in the primitive unit cell, which shows the true periodicity of the lattice and requires a smaller basis due to the smaller number of atoms. This means that a user-friendly program must match the information (matrix representations, traces, k points, etc.) from the DFT calculation, which depends on how the user defined the calculation, with the symmetry tables, which is fixed to the conventional setting. Consider the calculation cell with vectors (a1,a2,a3) defining the unit cell and the basis vectors of the conventional setting of the symmetry tables (c1,c2,c3) . The two cells are related by a linear transformation Pas (c1,c2,c3) = (a1,a2,a3)·P, (6.7) where P is invertible. There can also be a shift of origin between the cells given by the vector p=Oc−Oa,(6.8) where Oc and Oa are the origins of the conventional and calculation cells, respectively. The expression of a real-space point rin both settings is r=x′ ici=x′ iPjiaj+pjaj=xjaj,(6.9) We can express it in matrix notation as   x1 x2 x3 =P·  x′ 1 x′ 2 x′ 3 +  p1 p2 p3 .(6.10) 85 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package The inverse transformation is (P, p)−1= (P−1,−P−1p)and, consequently   x′ 1 x′ 2 x′ 3 =P−1·  x1 x2 x3 −P−1  p1 p2 p3 .(6.11) Upon a cell transformation, the reciprocal lattice vectors are only affected by the transformation P . The vectors qi corresponding to the standard cell must be related to the kiof the calculation cell by some yet-to-be-determined transformation qi=Ajikj.(6.12) The key to identify Ain terms of Pis that both bases must satisfy ki·aj=qi·cj= 2πδij.(6.13) Using the real-space basis transformation qi·cj=X n,m Amikm·anPnj = 2πX m AmiPmj = 2πX m [At]imPmj = 2πδij.(6.14) This implies that A= [Pt]−1 . It is straightforward to see that, if the bases are related by a pure rotation (with a possible translation) then A=P . This is not the case for transformations between unit cells of different volume, which are indeed possible. The transformation of reciprocal space coordinates of a vector l=y′ iqi is therefore given by   y1 y2 y3 =A·  y′ 1 y′ 2 y′ 3 = [Pt]−1·  y′ 1 y′ 2 y′ 3 .(6.15) The representation of a symmetry operation {R|t} also changes from the standard to the calculation cell. The action on a real-space position express in the standard setting is ˜ r={R′|t′}r= (R′ ijx′ j+t′ i)ci= ˜x′ ici.(6.16) From the point of view of the calculation cell, the same transformation is expressed as ˜ r= (Rijxj+ti)ai= ˜xiai.(6.17) Denoting the vectors of coefficients in both cells as ˜x′and ˜x, we have ˜x′=R′·x′+t′,˜x =R·˜x +t.(6.18) Since ˜x′and ˜x must agree upon cell transformation, we have ˜x =P˜x′+p=P[R′·x′+t′]+p=P·R′·x′+P·t′+p =P·R′(P−1·x−P−1·p) + P·t′+p=P·R′·P−1·x −P·R′·P−1·p+P·t′+p=R·x+t, (6.19) 86 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6.1. Computing the symmetry of eigenstates from DFT calculations which means that {R|t}={PR′P−1|P·t′+p−PR′P−1·p}.(6.20) This expression implies that a symmetry with only fractional translation in one setting may be a combination of a non-symmorphic operation and a lattice translation in the other, which has meaningful effects for points other than Γ. Spin SU(2) rotations require special attention due to the one-to-two mapping from real-space operations to spinor space. By inspecting Equation 6.6, we see that there is freedom to choose the axis ˆn and angle θ of rotation. The problem is exacerbated by the change of sign by a shift in the angle θ→θ+ 2π , which also means that twofold rotations through π and −π are different. All in all, this means that there can be a minus sign between the matrix in the symmetry tables and the ones generated by finding ˆn and θ programatically and applying Equation 6.6. To circumvent this problem, we implemented the following method: 1. Given the known real-space rotation matrices for all Shubnikov groups, we compute a SU(2) matrix using Equation 6.6 axis ˆn computed from the realspace matrix and angle θ . In the BCS, the SU(2) matrix is expressed in the spin axes setting, which is different from the conventional cell setting. We found that all SGs in the trigonal and hexagonal crystal systems share the same spin axes, while the five remaining systems share another different set of axes (see Table 6.1). This means that we can generate the SU(2) matrix for any operation in the tables up to a sign, due to the arbitrary choice of ±ˆn and θ. 2. Transform the axis from the conventional cell to the calculation cell using the change of basis Mca =P−1 which is automatically found or given by the user ˆna=Mca ·ˆn.(6.21) If |Mca|>1 , we have to normalize the axis. We also impose the constraint |Mca|>0so the handedness of the axes does not change. 3. Change the axis ˆna from the calculation cell axes to Cartesian coordinates, since implemented DFT suites express spin in this convention. The transformation matrix Mae is simply (a1,a2,a3) and another normalization of the axis ˆne= Maeˆnamay be required. 4. Find the SU(2) matrix in Cartesian coordinates using ˆne and θ in Equation 6.6 and multiply by the sign found in step 1. This ensures that the exact same rotation as in the symmetry tables is performed in the calculation setting. 87 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package trigonal and hexagonal rest s1(0,−1,0) (1,0,0) s2(−√3 2,1 2,0) (0,1,0) s3(0,0,−1) (0,0,1) Table 6.1: Spin axes used in the BCS for the SU(2) matrix depending on the crystal system. 6.2 New functionalities of irrep Analysis of magnetic-group symmetry With these new efforts, irrep now supports crystal displaying magnetic ordering, allowing the user to automatically detect the correct Shubnikov SG and analyze the symmetry properties of the energy eigenstates. The user is required to create a plain-text file with the magnetic moments of all the atoms in the unit cell, a row with three numbers (in Cartesian coordinates) for each one. irrep can be prompted to use MSG symmetry as follows >> irrep -magmom path/to/magnetic/momenta/file [other options] With this new option, irrep will detect the MSG for the input structure (e.g. from a POSCAR file) and magnetic momenta by using the spglib [71,72] library, retrieving the symbol and number in the BNS convention. Although irrep identification only requires the unitary symmetries, irrep will also list al the anti-unitary operations in both the DFT calculation and standard settings, according to the BCS. This new capability requires the full list of little-group coreps at all high-symmetry points for all 1,651 Shubnikov space groups. This information is available at the BCS in the COREPRESENTATIONS 2 application. The new irrep tables are made available locally and packed in a companion package, which is automatically installed when irrep is downloaded. ---------- INFORMATION ABOUT THE SPACE GROUP ---------- Space group Cm'cm'(# 63.464) has 4 unitary symmetry operations ### 1 rotation : | 1 0 0 | rotation : | 1 0 0 | | 0 1 0 | (refUC) | 0 1 0 | 2https://www.cryst.ehu.es/cgi-bin/cryst/programs/ corepresentations.pl 88 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6.2. New functionalities of irrep | 0 0 1 | | 0 0 1 | gk = [kx, ky, kz] | refUC: gk = [kx, ky, kz] spinor rot. : | 1.000+0.000j 0.000+0.000j | | 0.000+0.000j 1.000+0.000j | spinor rot. (refUC) : | 1.000+0.000j 0.000+0.000j | | 0.000+0.000j 1.000+0.000j | translation : [ 0.0000 0.0000 0.0000 ] translation (refUC) : [ 0.0000 0.0000 0.0000 ] axis: [0. 0. 1.] ; angle = 0 , inversion : False Since DFT programs may not take into account the full MSG, we recommend adjusting the -degenThresh parameter, which controls the minimum difference between two energies to consider them degenerate. MTQC analysis and decomposition of DFT bands into EBRs irrep now includes all the symmetry vectors of all the EBRs at maximal k points for the 1,651 Shubnikov SGs, as available in MBANDREP 3 application of the BCS. Additionally, the Smith normal form of every EBR matrix is stored as part of the dataset. The EBR analysis can by activated using the –ebr-decomposition command line flag, requiring the identification of the irreps at all maximal k points in the same irrep call. The set of bands to analyze is restricted via the -IBstart and -IBend flags. The program will do the following: 1. Find the symmetry vector of the input bands and check if there are integer solutions to the EBR decomposition problem using the Smith normal form. This allows the identification of strong topological materials. 2. If there exists an integer solution, some possible linear combinations of EBRs will be computed. This covers both the trivial and fragile topology cases. The computation EBR decomposition is found using the openly-available OR-Tools [73] package, casting the task into a constrained programming Boolean satisfiability (CP-SAT) problem. It first looks for linear combinations with all-positive coefficients. If none is found, the subroutine will relax the constraints to allow negative values. Among the many viable solutions, we chose to show those with smaller coefficients. 3https://www.cryst.ehu.es/cgi-bin/cryst/programs/mbandrep.pl 89 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package Search for transformation to conventional cell Previous versions of irrep required that the DFT calculation (in the (a1,a2,a3) basis) be done directly in the conventional cell (with vectors (c1,c2,c3) ) or that users provide the correct cell transformation via the -refUC and -shiftUC command-line arguments such that   c1 c2 c3 = refUC ·  a1 a2 a3 .(6.22) The settings may also be related by a change of origin shiftUC = Oconv −ODF T . This new version of irrep can find both the rotation and translation relating the calculation and conventional cell, which is required for the irrep identification since the tables are written in the latter setting. This can be done with the argument -searchcell . If prompted to do so, irrep will use the rotation provided by spglib and try to match the symmetry operations against the irrep tables. If it does not work, the program can still find the correct transformation by considering the following two cases: • If the space group is centrosymmetric, the origin of the primitive cell is placed in the inversion center. • If it does not have inversion symmetry, different centering translations are tried for the corresponding crystal system. High-symmetry points in conventional cell Since the program allows the use of any computation cell, this also has an impact in the reciprocal space coordinates of the high-symmetry points. More specifically, the coordinates of a reciprocal space vector G written in the conventional cell reciprocal basis (q1,q2,q3) are usually not the same in the calculation reciprocal basis (k1,k2,k3) . This is important since the high-symmetry point coordinates of the calculation must match with those of the tables upon transformation, indicating that the correct kpoints where computed in the calculation cell. To account for the relation between reciprocal bases in Equation 6.14, irrep now includes the pre-processing argument –print-hs-kpoints to display the highsymmetry points in both the DFT and conventional cells: 90 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6.2. New functionalities of irrep ########## High-symmetry k points ########## Reading symmetries from tables for SG 63.464 Change of coordinates from conventional to DFT cell: | 0.50 -0.50 0.00 | | 0.50 0.50 0.00 | | 0.00 0.00 1.00 | Change of coordinates from DFT to conventional cell: | 1.00 1.00 0.00 | | -1.00 1.00 0.00 | | 0.00 0.00 1.00 | Coordinates in symmetry tables: GM : 0.000000 0.000000 0.000000 R : 0.500000 0.500000 0.500000 S : 0.500000 0.500000 0.000000 T : 1.000000 0.000000 0.500000 Y : 1.000000 0.000000 0.000000 Z : 0.000000 0.000000 0.500000 Coordinates for DFT calculation: GM : 0.000000 0.000000 0.000000 R : 0.000000 0.500000 0.500000 S : 0.000000 0.500000 0.000000 T : 0.500000 0.500000 0.500000 Y : 0.500000 0.500000 0.000000 Z : 0.000000 0.000000 0.500000 The user can then compute the eigenstates at the correct high-symmetry points in order to identify the irreps later. The example clearly shows that the coordinates in both reciprocal bases can be different. Calculation of symmetry indicators Although calculation of the SIs through the EBR matrix is possible as they are all available in irrep , we choose to implement the SI formulas in the physical basis for all single and double SG s and MSGs. To enable this functionality, one must add the flag –symmetry-indicators to the options used to identify the irreps at all maximal k points. We have stored all the SI formulas in terms of irrep multiplicities as obtained from the BCS, which can be used once the identification or irreps is done. We also report the cases where the (M)SG does not have a non-trivial SI group. One note is due regarding the set of SIs returned for a given (M)SG. In accordance to the information given in the BCS, we also chose to report every SI of the (M)SG 91 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 6. Diagnosing band topology from ab-initio calculations using the IrRep package Solution 1 1 x [ -E->G(4) @ 2a(4'3m',23) ] + 1 x [ -1Eu->G(4) @ 4b(3m',3) ] + 1 x [ -2Eu->G(4) @ 4b(3m',3) ] + 1 x [ -Eu->G(4) @ 4b(3m',3) ] + 1 x [ -1Eu->G(4) @ 4c(3m',3) ] + 1 x [ -2Eu->G(4) @ 4c(3m',3) ] + 1 x [ -Eu->G(4) @ 4c(3m',3) ] + 1 x [ -1E->G(12) @ 12f(2'2'2,2) ] + 2 x [ -2E->G(12) @ 12f(2'2'2,2) ] ... 98 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 CHAPTER 7 Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co 99 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co In Chapter 2, we described the tools to analyze the symmetry of a physical system, with an emphasis on crystalline structures, which show translational and rotational symmetry. We also showed how to treat the presence of a local magnetic ordering by considering the type-I, III and IV Shubnikov space groups. Chapter 4 also explained how topology can give rise to characteristic responses in the material, such as the anomalous Hall conductance and how it is related to features in the electronic band structure such as Weyl nodes. In this chapter, we will make extensive use of these tools to investigate two specific materials: Fe3GeTe2and hexagonal Co. Coincidentally, both Fe 3 GeTe 2 and this phase of Co share the same parent space group when magnetization is not considered, which is why the two cases of study are presented in this same chapter. Although this facilitated the analysis of both systems in the first place, each one shows unique features that should be studied in detail. 7.1 Mirror-protected nodal lines Nodal lines are closed, continuous lines in reciprocal space where two or more bands intersect. Although crossings between two bands can generically appear in a 3D crystal with no additional requirements (i.e. Weyl or Dirac nodes), there is no guarantee that they will assemble into a continuous set of nodal points, which is why nodal lines require some kind of symmetry protection. To see why symmetry may give rise to protected crossings, consider two Bloch states |ψ(k)⟩ and |ϕ(k)⟩ transforming as different coreps D and ˜ D of the little group Gk of k , which we consider to be one-dimensional for simplicity. A given g∈Gk acts on the states with D(g) and ˜ D(g) respectively. Because g is in the little group, the Bloch Hamiltonian H(k)must be invariant g−1H(k)g=H(gk) = H(k),(7.1) which follows because, by definition, gk=k+G≡k , where G is a vector of the reciprocal lattice. It follows that the following expectation value is invariant ⟨ψ(k)|g−1H(k)g|ϕ(k)⟩=⟨ψ(k)|H(k)|ϕ(k)⟩.(7.2) Equivalently, we could apply the operators to the states instead of transforming H(k) , which would yield D(g)−1˜ D(g)⟨ψ(k)|H(k)|ϕ(k)⟩=⟨ψ(k)|H(k)|ϕ(k)⟩.(7.3) This shows that this matrix element can only be non-zero if D(g)−1˜ D(g) is the identity, that is, |ψ(k)⟩ and |ϕ(k)⟩ transform as equivalent representations. It follows then 100 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.1. Mirror-protected nodal lines that hybridization and gapping of bands between states of the different symmetry is forbidden. For Fe 3 GeTe 2 and Co, we are interested in mirror symmetries due to the particular space group they belong to. For example, consider a mirror plane m001 perpendicular to the c axis. Any point on a plane of the form kP= (kx, ky, kz = 0,1/2) is invariant with respect to m001 , which means that m001 ∈GkP . In more detail, if GkP has two coreps whose trace for the m001 symmetry element differ, we say that the crossing between the states of different irreps is protected by the mirror symmetry. The hybridization is forbidden at all of the points kP on the plane, which means that two bands can intersect and give rise to a continuous line of degeneracies, a nodal line. This is precisely the configuration we have in both Fe3GeTe2and Co. How to detect mirror-protected nodal lines In this work, we characterized the symmetry-protected nodal lines by a method based on Wilson loops (WLs). Following our example of the m001 mirror symmetry, we computed the evolution along paths on the kP plane of the trace of WLs computed along the perpendicular direction kz . Due to the mirror symmetry, loops of this kind are quantized to either 0 or π , modulo 2π . A discontinuity of π identifies as symmetry-protected nodal line that cannot be removed by a small, local perturbation that preserves the m001 symmetry [79]. The reason for this is that the trace of the WLs of the form W(kx, ky) is related to the number of states of positive or negative eigenvalue which are occupied for a given point on the plane. Therefore, as sketched in Figure 7.1, when tracking the evolution of the bands just below a nodal crossing, the number of states with positive or negative eigenvalue changes by one as one crosses a nodal line. In practice, we model real materials from ab-initio calculations in the framework of DFT in this case. Although computing WLs and other quantities related to the Berry connection is possible from DFT, computing the eigenstates in a sufficiently dense grid of k points is out of the question, since it would make the calculations extremely costly. For this reason, we first find a basis of Wannier functions that accurately reproduce the electronic band structure in a window around the Fermi level using Wannier90. This procedure also provides an interpolated tight-binding Hamiltonian written in the Wannier basis, which for a reduced number of bands of interest is more manageable. Obtaining the eigenstates and energies for a k point is therefore reduced to diagonalization of the Bloch Hamiltonian, which is fast and can be very easily parallelized. To compute the Wannier functions, one specifies a set of starting orbitals to project the Bloch states onto, usually those with the larger contribution around the Fermi level. In essence, Wannier90 finds then a set of Bloch states from the initial projections that is optimally smooth in k and constructs the basis set. 101 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co k a) b) kx ky Energy Figure 7.1: Sketch of nodal line and identification via WLs. a) Energy dispersion where the color indicates the mirror eigenvalue of the states. The broken line separates the states tracked by the WLs (below) from the rest (above). b) Cut of BZ for a constant kz= 0,1 2 in fractional coordinates. The circle depicts the nodal line in (a) and the broken line is a directed path along which we compute the WLs integrating in the perpendicular direction out of the plane. The trace of the WL undergoes a π jump when it intersects the circle and, at the same time, the number of states with a given mirror eigenvalue (red band in (a)) changes by one. Upon computing the Wannier Hamiltonian, we use a post-processing tool called WannierTools [80], which enables two calculations in a simple way. First, we can compute all the crossing points between one band and next one throughout the BZ, allowing us to locate the nodal line candidates. Secondly, we can also perform WL calculations along any path in the BZ to confirm the nodal crossings. This is necessary in practice because the algorithm for searching crossing points may fail to find all of them, depending on the k grid and specific features of the band structure. Therefore, lines that may seem broken in the calculation may be symmetry-protected nodal lines, which are unequivocally detected with this method. Drum-head states The Wilson loop method has the advantage of being able to identify protected surface states that live inside the projection of the nodal line onto the 2D surface BZ, called drum-head states [81]. The name refers to these states extending all over the projection of the nodal line, resembling the batter head of a drum whose hoop is precisely the projection of the crossings (i.e. inside the circle in Figure 7.1b). Drum-head states are particularly interesting because their dispersion is usually fairly flat, which promotes the development of strong correlations. 102 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.1. Mirror-protected nodal lines a) b) Figure 7.2: Sketch of the two possible positions of the hybrid Wannier function centers of the nodal lines for a slab of three cells in the c -axis direction. The colored spheres represent Fe (orange-brown) Te (brass) and Ge (purple) atoms. The red square is the position of the hybrid Wannier function of the occupied bands of a Wilson loop calculation along kz .a) The trace of the Wilson loop is ϕ= 0 and the HWF center is at z= 0 .b) ϕ=π and the center is at z= 1/2 , which gives rise to surface states since the surface termination is exactly at the HWF center. We recall from Chapter 4 that the trace of the WL operator along kz (out of the invariant plane) is proportional to the average z position of the centers of the hybrid Wannier functions, extended in x, y and localized in z . The trace is quantized to 0 or π , meaning that the charge center of the bands just below de nodal crossing is either at the origin or at the unit cell boundary. If, when making the crystal finite in this direction, we cut the structure at the point where the charge center is, there will be a dangling charge at the surface, which is precisely the drum-head states. This is shown in Figure 7.2, where we have sketched the two situations for Fe3GeTe2. Very importantly, the actual value of the WL trace depends on the position of the mirror plane inside the unit cell, which is convention dependent. For example, a shift of half the c axis of all atomic sites would exchange the 0 and π (modulo 2π ) values of this kind of WLs. This is precisely because the position of the hybrid charge center inside the unit cell depends of the reference origin. With this in mind, Figure 7.2b is the case that shows surface drum-head states in our specific choice of origin for Fe 3 GeTe 2 , which corresponds to a WL trace inside the nodal-line projection of ϕ=π . This dependence on the choice of origin raises an important question, since the trace of the WL operator must be gauge invariant, as shown in Eq.4.24. A displacement 103 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co of the atoms by half the c axis corresponds to a gauge transformation of the Bloch states by a phase, which is a unitary transformation, leading to an apparent contradiction. However, gauge invariance holds for an isolated set of bands, which is not this case since the bands just below the nodal crossing are indeed not isolated from at least the first above, as sketched in Figure 7.1a. 7.2 Study of Fe3GeTe2 Fe 3 GeTe 2 is an itinerant ferromagnet with a Curie temperature Tc≈220K and a hexagonal layered structure [82]. The bulk material is formed by slabs of Fe 3 Ge lying between layers of Te bound by weak van der Waals forces to form a 3D structure, which makes it easily exfoliable [83] or grown by deposition [84]. The two nonequivalent Fe positions in the unit cell display a local magnetic moment of around 1.4µB per atom along the stacking c axis [85,86]. Previous studies have shown that Fe 3 GeTe 2 displays a large anomalous Hall conductivity [87], anomalous Nernst effect [88] and non-linear quantum Hall effect [89], which has drawn the attention for their applicability in the development of novel spintronic devices. Its magnetic properties have been studied in-depth [90 – 94], discovering the possibility to tune its magnetic configuration via strain [95] or applied voltage [96]. Moreover, other studies claim that Fe 3 GeTe 2 shows Kondo lattice behavior [97] and skyrmionic phases [98]. Full knowledge of the symmetry of Fe 3 GeTe 2 is necessary to understand its topological properties. It offers another platform to investigate the interplay between topology and symmetry in magnetically ordered materials, since its ground state is described by a MSG or a spin space group when SOC is neglected. The non-zero local magnetic moments break TRS intrinsically, which is necessary for some features to appear, such as Weyl nodes [99]. Previous research also shows that Fe 3 GeTe 2 hosts one nodal line located very near the Fermi level, when SOC is neglected [100]. In particular, the band structure shows a twofold-degenerate crossing at the K= (1/3,1/3,0) point. This degeneracy extends along the high-symmetry P line joining the K and H= (1/3,1/3,1/2) wave vectors. When SOC is taken into account, the nodal line was shown to be gapped, leading to an avoided crossing with a very small energy gap along which the Berry curvature is concentrated, which results in a large flux. This, in turn, induces a large anomalous Hall conductivity (AHC), which was forbidden by the TRS in the non-magnetic structure. This mechanism was put forward to explain the high AHC reported in experiments. 104 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.2. Study of Fe3GeTe2 Symmetry analysis, DFT and wannierization methods In order to study the electronic band structure of Fe 3 GeTe 2 , we performed DFT calculations as implemented in the Vienna Ab-initio Simulation Package (VASP) [74,75]. We used Projector Augmented Wave [101] (PAW) pseudo-potentials with the Perdew-Burke-Ernzhenhof [54] (PBE) implementation of the Generalized Gradient Approximation (GGA) for the exchange-correlation functional. Additionally, van der Waals forces were included using the DFT-D3 method with Becke-Johnson damping function [102]. The energy cut-off of the plane wave basis was set to 600eV and the k point grid was a Gamma-centered, regular Monkhorst-Pack point set of dimensions 11×11×7 for a uniform density of points given the unit cell taken from the reference literature [82]. To compute the maximally localized Wannier functions from DFT, we started with a basis of d orbitals for Fe and p for Te and Ge, with frozen energy window of 3 eV centered at the Fermi level. We therefore obtain a total of 96 basis functions, taking into account the doubling due to spin. Symmetry analysis of the band structure Bands without SOC Fe 3 GeTe 2 displays a ferromagnetic (FM) structure where the magnetic moments at Fe sites align parallel to the c axis. The space group symmetry corresponds to SG P63/mmc (No.194). The Fe atoms sit at tow different WPs labeled 2b and 4a. The calculations show that these two inequivalent positions host a magnetic moment per atom of approximately 1.5µB and 2µB , respectively, which is in good agreement with the reported experimental measurement of 1.4µB per atom. In the absence of SOC (Fig.7.3a), the correct symmetry description is given by a Type-1 collinear spin space group (L194.1.1). This kind of group arises when SOC is not taken into account. In that situation, there are no terms of the type L·S , which couple the spin and real-space rotations, meaning that atoms and their magnetic local magnetic moments can be transformed separately. In that case, symmetry operations are usually denoted by a symbol {S|R|t} , where S acts on spin and R and t are real-space rotation and translation, respectively. The FM order without SOC respects the continuous rotation symmetry of the spins about the alignment axis. Using the spin space group theory, one can show that this decouples the bands into spin-up and spin-down subsets, each one showing the symmetry of a discrete subgroup of L194.1.1 which takes the form G×SZT 2,(7.4) where G contains only spatial operations with identity spin part and SZT 2 is generated by θUx(π) . Here, θ denotes TRS while Ux(π) is a π rotation about the x spin axis. 105 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co K M spin spin 1 0.75 0.5 0.25 0 -0.25 -0.5 -0.75 -1 E - EF (eV) K M H |K |K H 0.5 0.375 0.25 0.125 0 -0.125 -0.25 -0.375 -0.5 Sz (ħ units) a) c) d) b) M L H K A ky kz kx Fe Ge Te Figure 7.3: Electronic band structure along high-symmetry paths a) without SOC and colored spin sectors b) with SOC and Sz spin component c) Crystal structure as seen from the a and c directions. The arrows on the atoms represent the magnetic moments. d) BZ for the given structure with the high-symmetry path of the band plots marked in red. θUx(π) is represented by iσzK , where K is the complex conjugation operation. The operator acts as an analogue of TRS but, since (ΘUx(π))2= +1 , each sector obeys the symmetry of a single-valued gray group. Therefore, each spin sector can be analyzed separately as spinless fermions of a non-magnetic crystal. The previously reported nodal line occurs along the high-symmetry line P joining the points K and H . Since along the nodal line two bands become degenerate, we must only consider 2D irreps at the endpoints, therefore leaving as possible coreps K5, K6 and H1, H2, H3 respectively. Because the subduction of the coreps onto the P line must match at both ends of the path, we can discard the H3 symmetry. This is due to to K5 , K6 , H1 and H2 all subducing to P3 whereas H3 subduces to P1⊕P2 , violating the compatibility relations. We also note that there are non-degenerate energies for a given spin-sector at K (e.g. spin-down bands at around −0.6 eV) since there are one-dimensional coreps of the little group at that point. On the contrary, there are only 2D coreps at Hand all the bands are doubly degenerate. 106 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.2. Study of Fe3GeTe2 Generators Invariant plane Point group Coreps {m001|0,0,1 2},Θ{¯ 1|0}(kx, ky, kz= 0,1/2) {m110|0},Θ{¯ 1|0}(−k, k, kz) {m100|0},Θ{¯ 1|0}(kx= 0,1/2, ky, kz)2’/m A′, A′′ {m1¯ 10|0,0,1 2},Θ{¯ 1|0}(k, k, kz) {m120|0,0,1 2},Θ{¯ 1|0}(−2k, k, kz) {m210|0,0,1 2},Θ{¯ 1|0}(−k, 2k, kz) Table 7.1: Generators of the little groups (with translations), invariant coordinates, little co-group label and coreps for the band structure without SOC that protect crossings that may lead to a nodal line. All invariant planes share the same point group and coreps. The line P= (1/3,1/3, w) is related to (1/3,1/3,−w) by symmetry, which implies that the aforementioned twofold degeneracy forms a nodal line that closes due to periodic boundary conditions of the BZ. This is markedly different from the mirrorprotected loops explained before, since they are confined to a single high-symmetry line. For this SOC-less case, since both spin sectors are decoupled, the crossing between up and down bands is never gapped and such accidental degeneracies can give rise to nodal lines [103]. Apart from this, the mirror-protection also applies to this material, even without SOC. In particular, there are six mirror planes which protect nodal line crossings, as show in Table 7.1 For all of the perpendicular planes, the little co-group is isomorphic to 2′/m , which only has two coreps A′ and A′′ which differ by their mirror eigenvalue. Bands with SOC When SOC is considered, spin and spatial operations cannot be decoupled and we can directly apply the theory of MSGs, as explained in Chapter 2. As a consequence, the symmetry group is reduced to the Type-III MSG P63/mm′c′ (No.194.270), which admits a coset decomposition M=G+θ{2110|0}G, (7.5) where G is the unitary subgroup generated by {6+ 001|0,0,1/2} , {2001|0,0,1/2} , {¯ 1|0} and translations {1|ai} by the vectors ai defining the standard hexagonal unit cell. This MSG shares the Bravais lattice and symmetries of the unitary subgroup P63/m (No,176), while also respecting additional space-group operations of No.194 combined with TRS. As shown in Figure 7.3b, the spin sectors remain well differentiated except at crossing points or points where a gap between spin up and down bands is opened, 107 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co -1 800 600 400 200 0 80 100 60 40 20 0 -0.75 -0.50 -0.25 0 0.25 0.5 0.75 1 E - EF (eV) Absolute AHC (S/cm) Weyl poiint count Figure 7.8: Computed absolute AHC σAHC xy (black) and the energy distribution of the non-zero-chirality Weyl nodes (red) in a window extending 1.5 eV above and below the Fermi energy. -1 -0.75 -0.50 -0.25 00.25 0.5 0.75 1 centered cube full BZ 800 600 400 200 0 Absolute AHC (S/cm) E - EF (eV) a) b) 1 0.75 0.5 0.25 0 -0.25 -0.5 -0.75 -1 1 0.75 0.5 0.25 0 -0.25 -0.5 -0.75 -1 (normalized) E - EF (eV) KM Figure 7.9: a) AHC from a cube enclosing the Γ point (blue), which estimates the contribution of the SOC-related gaps in (b). It shows that the signal in the -0.5 to 0.1 eV range around the Fermi level can be attributed to this mechanism. b) Electronic band structure along high-symmetry paths on the kz= 0 plane with SOC and normalized Ωxy component of the Berry curvature. 114 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.3. Nodal lines in hexagonal Co the prediction of topological features and to explain the observed anomalous Hall conductivity. The symmetry of the FM phase, with SOC considered, is described by the MSG No.194.270 whose {m001|0,0,1/2} mirror plane, in particular, protects nodal lines on the invariant kx= 0,1/2 planes. We have found that Fe 3 GeTe 2 hosts many of these features (Fig.7.4) and given a rigorous identification of the nodal crossings through the discontinuities of Wilson loops over paths that intersect them. We have used the trace of the Wilson loop operator to predict the appearance of drum-head surface states in the projection of the nodal lines and have detected one example in the computed surface DOS (Fig.7.5). Furthermore, we have computed the AHC for a range of chemical potentials and detected the three main mechanisms that give rise to the AHC response. Firstly, while the paramagnetic phase also protects nodal lines on kx= 0,1/2 and ky= 0,1/2 , the FM transition gaps them and they contribute to the anomalous transport due to being points where the Berry curvature concentrates (Fig.7.6). However, they are not the only source of AHC, since its contribution is far from the total response computed. Secondly, we have identified many Weyl points in the band structure and computed their chirality (Fig.7.7). Their contribution to the AHC has to be complemented by the two other sources, especially in the -0.6 to -0.1 eV range, where Weyl crossings are almost absent (Fig.7.8). Finally, with magnetic order but neglected SOC, the spin up and down bands are decoupled and free to cross, each set respecting the symmetries of the single-valued gray MSG No.194.264. These degeneracies are lifted when SOC is brought into the picture, giving rise to small gaps where the curvature is also prominent (Fig.7.3c and 7.9a)). The calculations also show that slightly shifting the chemical potential by approximately 0.3 eV would further increase the AHC in Fe 3 GeTe 2 from around 200 S/cm to approximately 800 S/cm, suggesting electron doping as a good mechanism to obtain an enhanced effect. This prediction is supported by DFT calculation with varying number of electrons that account for the effect of electron and hole doping. These show that the magnetic moments vary only around 5% with respect to the neutral configuration and that the electronic band structure is almost unchanged, implying that the effect of doping can be approximated by a simple shift in the chemical potential. Our computations show that a doping of two electrons per unit cell (one per formula unit) achieves an almost optimal shift in the chemical potential that leads to the fourfold AHC enhancement (Figure7.10). 7.3 Nodal lines in hexagonal Co Hexagonal close-packed Co (hcp-Co) shares the same parent paramagnetic group as Fe 3 GeTe 2 and hence the band structure is subjected to similar symmetry constraints. Experimentally realized Co films exhibit different in-plane (IP) and out-of-plane (OP) magnetic anisotropy A , where A(x, y)> A(z) . As a consequence, the resulting 115 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co E-EF (eV) E-EF (eV) -2 -1 0 +1 +2 Figure 7.10: Absolute AHC with respect to the chemical potential. The numbers correspond to the variation in electron number with respect to the neutral configuration and the corresponding location of the chemical potential in the AHC curve is shown in the inset with vertical bars of matching color. hysteresis loop is not perfectly square neither in the IP nor the OP directions and the samples develop a canted FM order when magnetized with an electromagnetic pulse. Moreover, there are three IP magnetization directions separated by 60 degrees and this IP axis is randomly chosen for every domain when the external pulse is applied. However, no matter the magnetized IP axis, all domains will show an OP component and, as a result, most band features can be explained by the symmetry of the OP configuration since this contribution is consistent across domains. The symmetry corresponding to the IP magnetization is thus “smeared out” by averaging the inequivalent IP axes over all domains. Symmetry analysis, DFT and wannierization methods To perform the DFT calculations, we employed a similar procedure to Fe 3 GeTe 2 in VASP. To account for the correlations of the localized d electrons in Co, we introduced a Hubbard-like repulsive term (DFT+U) in the spherically symmetric approximation [110]. The repulsive parameter was adjusted to U= 2 eV to reproduce de experimentally reported magnetic moment per Co atoms of around 1.8 µ B [111]. The energy cutoff of the plane wave basis was set to 600 eV and the k point grid was a Gamma-centered, regular Monkhorst-Pack point set of dimensions 11 ×11 ×7 . The crystal structure found in the literature [112] was slightly expanded according to the experimental measurements we performed to better reproduce the exact sample features. 116 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.3. Nodal lines in hexagonal Co b)a) E - EF (eV) 0.0 0.5 1.0 1.5 2.0 -0.5 -1.0 -1.5 -2.0 KAL H A MKAL H A M Figure 7.11: DFT electronic band structure and unit cell with a) OP and b) IP magnetization. The localized Wannier functions as were obtained starting from a basis of s, p, d orbitals per Co atom for a total of 36 Wannier functions given the doubling due to spin. The bands were frozen in an energy window of 2 eV centered around the Fermi level during the disentanglement procedure to exactly reproduce the bands in that energy range. Symmetry analysis and DFT band structure Due to the experimental conditions described above, we will analyze the IP and OP configurations separately. With OP magnetization (Fig. 7.11a), the symmetry of the system is exactly the same as in FM Fe 3 GeTe 2 (No.194.270). Therefore, the same symmetry-protection mechanisms of nodal lines are present and the analysis is straightforward. For IP magnetization (Fig. 7.11b) along the specific [1¯ 10] direction, the symmetry is reduced to MSG Cm′cm′ (No.63.464), where one has to take into account that using the hexagonal cell corresponds to a non-standard setting. Because {m001|0,0,1/2} is not a symmetry by itself but only combined with TRS, nodal lines are not protected and therefore are expected to be gapped with respect to the OP configuration. 117 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co E - EF (eV) -0.4 0.4 1.91 1.84 1.77 1.70 1.64 1.57 1.50 1.36 1.22 1.15 1.08 1.02 kx (Å-1) 0.0 -0.4 -0.8 kx (Å-1) kx (Å-1) UL LL 1.0 0.5 -0.5 -1.0 -2.0 -1.0 0.0 0.0 1/2 -1/2 kz Fermi surface Nodal line a) b) c) d) H LA M M K K E Figure 7.12: Experimental data of a nodal line along A−L ( kz= 1/2 ). a) Band dispersions along paths perpendicular to A−L showing a cone-like dispersion whose crossing point is identified as a nodal line. b) Fermi surface map on kz= 1/2 where the contribution of the upper legs (UL) and lower legs (LL) of the cone feature in a in the Fermi surface is marked. c) Sketch of the BZ and surface BZ. d) Sketch of the band dispersion giving rise to the nodal line along A−L. Experimental observations ARPES experiments were performed to study the electronic band structure of hcp-Co. Figure 7.12 shows the experimental data on kz= 1/2 . A series of band dispersions along cuts perpendicular to the A−L (or ¯ Γ−¯ K in the surface BZ, see Fig. 7.12c) are presented in Figure 7.12a showing a conic band dispersion whose crossing point raises above the Fermi level and goes back down. This feature can also be observed in the Fermi surface (FS) map in Figure 7.12b, since the upper and lower “legs” of the cone (UL and LL, respectively) contribute to the Fermi surface at kz= 1/2 as evidenced in the band dispersions. This structure is sketched in Figure 7.12d. In Figure 7.13, we show the experimental data for the kz= 0 plane. The FS map reveals a flower-like pattern of two interlocking bands, which are signatures of a nodal ring sketched in Figures 7.13c-d. The band dispersion below the Fermi level along the red path in Figure 7.13a is shown in Figure 7.13b, where we have overlain the OP and IP DFT calculations. The plot shows that the band structure and in particular this nodal line are reproduced, as expected, only for the OP configuration and the intersection with this path is marked with an orange dot. There is also an additional band with large spectral weight at around -1.2 Å −1 . which is not seen in 118 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.3. Nodal lines in hexagonal Co E E1 E2 IP OP ky (Å-1) kx (Å-1) ky (Å-1) 0.5 0.0 -0.5 E-EF (eV) 0.0 -0.4 -0.8 -0.8 -0.8 -0.8 -0.8 -0.80.0 0.0 0.0 0.0 0.0 0.8 -2.0 -1.0 0.0 1.0 2.0 a) b) c) d) Node Node Figure 7.13: Experimental observation of a flower-like nodal line on kz= 0 .a) Fermi surface map. b) ¯ Γ−¯ K band dispersions along the red path in afor different kz ranging from kz= 0 (left) to kz= 1/2 (right). c) Sketch of band crossing structure where formed by two bands in red and blue. d). Sketch of the isolated nodal line resulting from c. the theoretical calculations. However, this is identified as a surface state since its dispersion does not change with kz. Theoretical analysis of the nodal lines Given the previous symmetry analysis and experimental features, we turn to identify the nodal crossings in the DFT calculations and match them with the observational data. We will restrict to the OP configuration, since we know this is the only case where a protection mechanism on the kz= 0,1/2 holds. First, we identify the nodal line along A−L (R line) on Figure 7.14b, marked in red. MSG No.194.270 has a single irrep ¯ R3¯ R4 along this high-symmetry line which is subduced into ¯ E3⊕¯ E4 at generic points in the kz= 1/2 plane. This means that the twofold degeneracy along R splits into to bands away from the line and corresponds to a crossing. Moreover, R= (u, 0,1/2) is equivalent to (−u, 0,1/2) by symmetry therefore allowing nodal lines extending from (−1/2,0,1/2) to (1/2,0,1/2) that close due to periodic boundary conditions. Since we know that the conic feature observed in Figure 7.12a must cross the Fermi level, we can pin-point the exact nodal along R out of all the possible ones in our DFT 119 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7. Magnetic-symmetry-enforced topological nodal lines: the cases of Fe3GeTe2and Co K M kx (Å-1) ky (Å-1) -2 02 -4 -2 0 2 4 E - EF (eV) 0 1 2 -1 -2 ky (1/Ang) kx (1/Ang) -0.5 -1 00.5 E - EF (eV) 1 1 0.5 0 -0.5 -0.7 -0.75 -0.80 -0.85 -0.9 -0.95 -1 -1 0.5 0 k (fractional) 0 0.2 0.4 0.6 0.8 1 WL trace / 2π a) b) c) d) e) OP IP E - EF (eV) 0.0 0.5 1.0 1.5 2.0 -0.5 -1.0 -1.5 -2.0 KAL H A M Figure 7.14: Computed nodal lines in hcp-Co with OP magnetization. a) A selection of nodal lines on kz= 0 lying near the Fermi level. b) Electronic band structure with intersection with nodal rings in amarked in matching colors. c) Dispersion of a line perpendicular to A−L for OP (left) and IP (right) magnetization. The OP dispersion displays the cone corresponding to a nodal line observed in Figure 7.12a and marked with a red line in b.d) Energy dispersion of the green nodal line in a.c) WL analysis of some nodal lines on kz= 0 , including the one in dwhich is manifested in the discontinuity marked with a red arrow. The WLs are computed along kz for varying kyalong the red path of the inset. 120 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 7.3. Nodal lines in hexagonal Co calculation and compute the conic band dispersion (Fig. 7.14c left) 1 . This does not hold in the IP configuration since, although there is also a single 2D corep (Q)¯ GP2¯ GP2 along the A−L line 2 , this is the only 2D corep for all the kz= 1/2 plane, implying that the degeneracy along the high-symmetry line is not a crossing as shown in Figure 7.14c (right). The flower-like nodal line on kz= 0 can be found by identical methods as the ones used for Fe 3 GeTe 2 . In Figure 7.14a, we show a selection of nodal crossings on this plane that lie near the Fermi level, as can be seen by their intersections with matching colors in the high-symmetry path of Figure 7.14b. The nodal ring marked in green, whose computed dispersion with kx, ky is shown in Figure 7.14d, corresponds exactly to the feature observed in the experiment, shown in Figure 7.13d. In Figure 7.14e, we plot evolution with ky of WL computed along kz , whose third discontinuity is due to this protected nodal ring. Summary In conclusion, we have seen that hcp-Co displays nodal lines protected by the {m001|0,0,1/2} symmetry, although the magnetic anisotropy causes the magnetization to be canted with respect to the c axis of the hexagonal cell. We have identified through ARPES measurements two kinds of nodal line structures in this compound. First, a flower-like nodal degeneracy cause by two interlocking bands on the kz= 0 plane. Second, a nodal line along the A−L high-symmetry line on kz= 1/2 . Using a Wannier TB model that reproduces the low energy spectrum from DFT, we have identified these features and analyzed them using Wilson loop calculations and symmetry considerations to confirm their status as symmetry-protected nodal degeneracies. Our computations have also detected many other nodal lines that are present in a window of 1 eV around the Fermi level. 1 This figure shows explicitly why there cannot be any nodal lines on kz= 1/2 that cross the R line in both Fe3GeTe2and OP hcp-Co due to the symmetry constraints 2This is not a high-symmetry line in MSG No.63.464. 121 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 CHAPTER 8 Topological transitions in the FeSe superconductor via doping and pressure. 123 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8. Topological transitions in the FeSe superconductor via doping and pressure. In particular, we note that no operation exchanges irreps so, for example, aij will not transform to bkl under any operation g∈G . WPs positions are not mixed either since, by definition, they contain all the sites related by all symmetry operations (inside the same unit cell). Upon taking the Fourier transform αij(k) = 1 √NX R e−ik·(R+ri)αij(R),(8.5) where rα is the position of αij inside the unit cell and N is the total number of cells, the R degrees of freedom are changed by k and the Hamiltonian is decoupled into independent blocks H=X k H(k).(8.6) The transformation properties of αij(k) can be induced from real space as in Equation 2.32 {R|t}αij(k) = 1 √NX T e−ik·(T+ri){R|t}αij(T) =1 √NX T e−ik·(T+ri)Dα(g)jkαlk(T′), (8.7) where T′+rl=R(T+ri)+t. Changing the sum to T′ {R|t}αij(k) = 1 √NX T′ e−ik·R−1(T′+rl−t)Dα(g)jkαlk(T′) =eiRk·tDα(g)jkαlk(Rk) (8.8) When {R|t} is in the little group Gk of k , then we can write the representation as a 12 ×12 matrix. In that case, since we can read Dα(g) off the BCS, we can readily construct the representation in our basis, with matrices D(g) , for the generators of G D({2001|1/2,1/2,0})=ei 2(kx+ky) 1 6⊗−i0 0i(8.9) D({2010|0,1/2,1/2})=ei 2(ky+kz)  A0 0 0B0 0 0 B , A=    0 0 0 e−iπ/4 0 0 e−i3π/40 0e−iπ/40 0 e−i3π/40 0 0     B=    0 0 0 ei3π/4 0 0 eiπ/40 0ei3π/40 0 eiπ/40 0 0     (8.10) 130 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8.1. Surface states in FeTe0.55Se0.45 D({4+ 001|1/2,0,0})=ei 2kx  C0 0 0F0 0 0 J , C =    0 0 ei3π/40 000e−i3π/4 ei3π/4000 0e−i3π/40 0     F=    0 0 e−iπ/40 000eiπ/4 e−iπ/4000 0eiπ/40 0     , J =e−π/40 0eiπ/4⊗ 1 2 (8.11) D{¯ 1|0}) = 1 3⊗    0 0 1 0 0 0 0 1 1 0 0 0 0 1 0 0     (8.12) There are many strategies to obtain all the invariant terms in a Hamiltonian. In this work, we listed the coupling in order of increasing distance and sequentially introduced the terms in the expression. For each distance, we select a coupling and obtain al the symmetry related hopping terms by applying the generators of G until the term is invariant. As in the k·p model, one then has two additionally impose TRS and hermiticity. We present the final expression for the Hamiltonian in k space H(k) =   ϵa0 0 0ϵb0 0 0 ϵc ⊗ 1 3+ ∆(k) + ∆†(k),(8.13) where ∆(k) =                     0 0 f10000 0f+ 2(kx) 0 f+∗ 2(ky) 0 000f10 0 0 0 0 f− 2(kx) 0 f−∗ 2(ky) 0 0 0 0 0 0 0 0 f− 2(ky) 0 f−∗ 2(kx) 0 0 0 0 0 0 0 0 0 0 f+ 2(ky)′0f+∗ 2(kx) 0 0 0 0 0 0 f′ 10f3(kx) 0 f∗ 2(ky) 0 000 0000f′ 10f3(kx) 0 f∗ 3(ky) 0 0 0 0 0 0 0 0 f3(ky) 0 f∗ 3(kx) 0 0 0 0 0 0 0 0 0 0 f3(ky) 0 f∗ 3(kx) 0 0 0 0 0 0 0 0 0 0 f40 0 0 0 0 0 0 0 0 0 0 0 f4 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0                     . (8.14) 131 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8. Topological transitions in the FeSe superconductor via doping and pressure. The matrix elements depend on real parameters of the model as follows: ϵa=a+lcos kz, ϵb=b+mcos kz, ϵc=c+ (cos kx+ cos ky)(e−ikzξ+eikzξ′)+ncos kz, f1=δcos kx 2cos ky 2, f′ 1=δ′cos kx 2cos ky 2, f± 2(s)=±ite−ikz(zc−1) cos s 2+q′e−ikzzccos s 2, f3(s)=t′e−ikz(zc−1) cos s 2+q′e−ikzzccos s 2, f4=τe2ikz(zc−1) cos kx 2cos ky 2 (8.15) In the expressions above, zc and 1−zc are the chalcogen positions in the unit cell along the caxis and the rest of parameters correspond to: •a, b, c: onsite energies of aij, bij, cij orbitals •l, m, n : couplings between orbitals aij, bij, cij in adjacent cells in the z direction •ξ, ξ′ : hopping between cij orbitals along lattice vectors (−1,0,−1) , (0,−1,−1) , (0,1,−1) and (1,0,−1) for ξand their opposites for ξ′ •δ, δ′ : a1j - a2j and b1j - b2j couplings, respectively, for lattice vectors (0,−1,0) , (0,0,0),(1,−1,0) and (1,0,0) •τ : coupling of c1j and c2j for vectors (−1,−1,2) , (−1,0,2) , (0,−1,2) and (0,0,2). •t, t′ : nearest-neighbors hybridization terms among aij and bij , respectively, with c1j separated by (0,0,−1) and (1,0,−1) , and with c2j by (0,−1,1) and (0,0,1) The parameters can be fitted to reproduce the DFT band dispersion along Γ−Z (see Table 8.2) and the result is shown in Figure 8.4a. Adjustments to the DFT fit have to be made in order to correctly reproduce the highly-correlated phase detected by experiments. In going, to Figure 8.4b, one has to adjust the d band dispersion by a renormalization factor of 3 [159,160,177], which amounts to decreasing the l, m parameters. Furthermore, we increase the strength of the SOC coupling between dxz, dyz and pz orbitals [178], obtaining Figure 8.4c. A coupling term is considered to include SOC if either •Couples spin-up and down orbitals 132 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8.1. Surface states in FeTe0.55Se0.45 Parameters of fit a−2.33 δ5.30 n1.625 b−4.40 δ′9.25 ξ−1.0 c−6.13 τ−14.0ξ′0 t0.16 l−0.049 q0.4 t′0.19 m−0.049 q′−0.01 Table 8.2: Parameters of the model to reproduce the DFT band dispersion around the Fermi level along Γ−Z. EF -E (eV) 0.2 0.3 0.0 a) b) c) d) Figure 8.4: Fit of the TB model dispersion to DFT adjusted to experimental evidence. a) Fit of DFT dispersion b) Adjusted to account for d band renormalization factor of 3 [159, 160, 177] c) With SOC coupling between dxz and pz increased [178] d) Comparison with experimental measurements. • Couples spin-up with spin-up with a different amplitude than spin-down and spin-down Since the second option is forbidden by TRS, given the change of basis in Equation 8.4 and the form of the Hamiltonian, t, t′, q, q′ are the only four SOC-enabled hopping amplitudes. We also find that the size of the pz - d gap is mainly governed by q, q′ , since it couples aij, bij to cij , which have dxz, dyz and pz character, respectively. Figure 8.4d shows the comparison with the experimental result, showing that the pz band completely disappears. Since this situation cannot be obtained without substantial change in the fitted parameters, this is in line with correlations in FeSe being very strong and not being captured by a free-electron model correctly. Summary We have shown that, according to the DFT calculations, FeTe 0.55 Se 0.45 shows a band inversion that may be the cause of the topological transition into a strong Z2 states, which hosts Dirac surface states around gamma in a slab configuration finite along 133 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8. Topological transitions in the FeSe superconductor via doping and pressure. the c axis. We have performed a symmetry analysis to construct an EBR-based TB model that reproduces the low energy spectrum, especially along the Γ−Z line, where the inversion occurs. However, due to the strong renormalization, evidenced by the ARPES measurements, the DFT picture is not accurate. To match our ab-initio calculations with the observed band features, we show that one can empirically adjust the band dispersion by a renormalization factor of 3, in line with previous results. All in all, our work provided strong evidence of the non-trivial band topology in FeTe 0.55 Se 0.45 and the theoretical calculations helped identify the precise reason for this phase transition in the doped compound, originating from a band inversion along Γ−Zand a gap opening through dand porbital hybridization. 8.2 Engineering topological phases in FeSe via uniaxial strain Our work on FeTe 0.55 Se 0.45 helped us understand the importance of the chalcogen position and the sensibility of the low energy spectrum along Γ−Z for driving topological phase transitions in this particular system through band inversions. Inspired by these observations, it is natural to ask the question of whether a similar effect can be achieved by slightly modifying the lattice structure through external pressure. In the following sections, we show that uniaxial strain can indeed achieve topologically non-trivial phases. Moreover, the direction and strength of this external perturbation can select two different weak and strong TI phases. Methods Electronic structure calculations were performed within density functional theory (DFT) with the Vienna Ab-initio Simulation Package (VASP) software [74,75], choosing Projector Augmented Wave (PAW) pseudo-potentials [101] with the Perdew-BurkeErnzenhof (PBE) implementation of the Generalized Gradient Approximation (GGA) for the exchange-correlation functional [54]. The effect of spin-orbit coupling (SOC) was included in all calculations. The energy cutoff of the plane wave basis we set to 600 eV. As k point grid we employed a Gamma-centered, regular Monkhorst-Pack point set of dimensions 11 ×11 ×10 , providing a sampling of uniform density in all directions of the Brillouin zone (BZ). A gaussian broadening of 0.1 eV has been used throughout the calculations. Structural optimization of the ion positions was performed using the conjugate gradient method implemented in VASP until all forces exerted on the atoms had a modulus smaller than 10−3 eV/Å. To relax the strained structures we set the respective components of the stress tensor to zero before the gradient step. 134 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8.2. Engineering topological phases in FeSe via uniaxial strain P4/nmm No.129 Pmmn No.59 a1 a2 a3 Figure 8.5: Atomic structure of FeSe. Strain is applied along the a1 direction of the unit cell (red arrows) of the uncompressed tetragonal structure belonging to SG No.129. Once enough uniaxial pressure is applied, the systems adopts an orthorhombic symmetry described by SG No.59. The length of the a2 axis is exaggerated to better display the orthorhombic symmetry. The parameters of the tetragonal and orthorhombic phases are specified in in the appendix or supl. inf. In order to construct a basis set of maximally localized Wannier functions for any of the crystal structures, we employed Wannier90 [66] using a starting basis of d orbitals for Fe and p for Se. This results, after the wannierization process, in a total of 32 basis functions, taking into account the doubling due to spin. A frozen disentanglement window of approximately 2 eV centered at the corresponding Fermi level is chosen to reproduce the low energy spectrum and the gaps where surface states appear. Electronic structure at ambient pressure In regular conditions, FeSe crystallizes in a tetragonal structure whose symmetry is described by the space group (SG) P4/nmm , No.129 in the Belov-Neronova-Smirnova (BNS) notation [78]. The group is generated by the identity, lattice translations and the elements {4+ 001|1/2,0,0} , {2010|0,1/2,0} with inversion {¯ 1|0,0,0} , in Seitz notation. The conventional, primitive unit cell contains 2 Fe atoms in the 2a Wyckoff position (WP) and 2 Se atoms in the 2c WP. Previous studies [130,149] found that the band structure, especially along the Λ high-symmetry line joining Γ = (0,0,0) and Z= (0,0,1/2) , in fractional coordinates, is very sensitive to the z position of the Se atoms inside the unit cell, one at z and the other at −z . Therefore, doping FeSe with iso-electronic Te impurities located at the chalcogen positions modifies the band dispersion along Λ and it has been shown to result in a topological phase transition due to band inversion [149]. The introduction of Te modifies therefore the coupling between Te/Se atoms in adjacent cells in the a3 direction both by changing 135 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8. Topological transitions in the FeSe superconductor via doping and pressure. M X RZAZ 1.0 0.5 0.0 -0.5 -1.0 -1.5 -2.0 E - EF (eV) Figure 8.6: Electronic band structure of the experimentally determined tetragonal crystal structure of FeSe. The bands up to the green, red and blue sets correspond to fillings of 20, 24 and 28 electrons respectively. The states at Γ and Z are labeled by irreps of their little groups, according to the symmetry of SG 129. the z coordinate of chalcogen atoms in the unit cell and increasing the overlap due to more extended pzorbitals [149]. We first performed DFT calculations of the undoped FeSe crystal structure as measured experimentally, shown in Figure 8.6. Contrary to doped crystals, one can check that the structure at ambient pressure does not have a gap for a filling of 24 electrons, near the Fermi level because the top red band and the lowest blue one at Z have ¯ Z9 and ¯ Z7 , respectively, which subduce to ¯ Λ6 and ¯ Λ7 along the GM −Z line and thus cannot hybridize since they have different symmetry [179]. Consequently, one cannot define a notion of topology near the Fermi level. Note that throughout this work “filling” refers to the band index rather than the actual number of electrons obtained by integrating the density of states (DOS) up to the chemical potential. The reported sensitivity of the dispersion of the down-crossing pz -dominant band and the energy of the flatter d bands suggest that modifying the band structure through pressure may achieve a gap opening and a topological phase transition. Effect of strain To study the effect of uniaxial strain on the original FeSe structure, we mimicked uniaxial strain by fixing one of the cell directions to a reduced value and relaxing the remaining lattice parameters and atomic positions until the force on every atom is less than 1 meV / Å. The pressure corresponding to one strain configuration can be extracted from the component of stress tensor in the a1 direction, which provides an approximate value necessary for the experimental realization. 136 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8.2. Engineering topological phases in FeSe via uniaxial strain Structure SG 24 electrons SIs TI 28 electrons SIs TI acompression 1% 59 NLC z2w,3= 0 Strong NLC z2w,3= 1 Strong z4= 3 z4= 1 ccompression 3% 129 SEBR z2w,3= 1 Weak LCAO z2w,3= 0 Trivial z4= 2 z4= 0 cexpansion 1% 129 NLC z2w,3= 0 Strong - - - z4= 3 Table 8.3: Symmetry and topology for different FeSe structures under uniaxial strain. The first column indicates the ratio between the experimental and computed length of the a1 , a2 and a3 vectors. The third and sixth columns are the irrep decomposition: split linear combination of atomic orbitals (LCAO), split EBR (SEBR) and non-linear combination of EBRs (NLC). The columns 4-5 and 7-8 contain the SIs and the TI classification for fillings of 24 and 28, respectively. M X RZAZ 1.0 0.5 0.0 -0.5 -1.0 -1.5 -2.0 E - EF (eV) Figure 8.7: Electronic band structure for crystal under strain along the a1 direction which reduces its length by 1% and drives a transition to an orthorhombic phase (SG No.59). The bands up to the green, red and blue sets correspond to fillings of 20, 24 and 28 electrons respectively. The states at Γ and Z are labeled by their little group irreps. 137 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8. Topological transitions in the FeSe superconductor via doping and pressure. Intensity (a.u.) M M M M a) b) E-EF(eV) 0.22 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 Figure 8.8: Density of states for finite geometry in the a3 direction computed using the Wannier TB model. Surface spectrum along the direction M−Γ−M for the crystal a) compressed in the a1 direction b) expanded in the a3 direction. The insets show a zoom into the surface Dirac cone due to the strong Z2-odd topology. Compression along aaxis Following the above described method, we find that a contraction along a1 or a2 favors a topological phase transition by increasing the band dispersion of the Se pz dominant band, similar to the effect previously achieved by Te doping [149]. Figure 8.7 shows the band structure when the lattice parameter along the a1 direction is reduced by 1%, corresponding to a pressure of 0.89 GPa. The symmetry is lowered to an orthorhombic phase described by SG Pmmn (No.59), which is essential for the band structure to show a gap for a filling of 24 electrons, as shown in Figure8.7. In particular, the subduction of irreps from SG No.129 to No.59. ¯ Γ6,¯ Γ7→¯ Γ5¯ Γ8,¯ Γ9→¯ Γ6(8.16) and the fact that now there is only one irrep along the Λ line mean that the dispersive band, with more pronounced dispersion than in the unstrained structure, can now hybridize and gap out with the flatter one. The disconnected bands for a filling of 20 electrons is a non-linear combination (NLC) of Elementary Band Representations (EBRs) and the next four bands up to 24 electrons are trivial. Therefore, the gap at 24 electrons is topological since overall the 24 bands cannot be decomposed as a sum of EBRs with integer, positive coefficients, meaning there is no trivial atomic limit that reproduces this set of bands. The filling of 28 electrons (blue) is again a NLC because the four bands marked in blue correspond to a a branch of the EBR induced from irrep ¯ Ag¯ Agof the site-symmetry group ¯ 1at WP 4c. We can further refine the topological classification using the irrep decomposition to find the SIs [25,180], which can be obtained from the inversion eigenvalues at time138 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341 8.2. Engineering topological phases in FeSe via uniaxial strain M X RZAZ 1.0 0.5 0.0 -0.5 -1.0 -1.5 -2.0 E - EF (eV) M X RZAZ 1.0 0.5 0.0 -0.5 -1.0 -1.5 -2.0 E - EF (eV) a) b) Figure 8.9: Electronic band structure of FeSe under uniaxial strain along the c axis. a) Reducing a3 by 3%. b) Expanding a3 by 1%. The symmetry is still tetragonal and described by SG No.129. The bands up to the green, red and blue sets correspond to fillings of 20, 24 and 28 electrons respectively. The states at Γ and Z are labeled by their little group irreps. reversal-invariant momenta (TRIMs). For a filling of 24, the values z2w,3= 0 mod 2 and z4= 3 mod 4 for SG No.59 indicate that strain induces a transition to a strong TI protected by TRS. A strong TI cannot be deformed into a stack of 2D quantum spin Hall insulators and displays surface states on any of the surfaces. Likewise, the filling of 28 shows z2w,3= 1 and z4= 1 which, although the weak index is odd, also indicates a strong TI phase because the z4 index is also odd. Figures 8.8a shows the density of states of the orthorhombic phase along the M−Γ−M path computed from the Wannier TB model for the finite slab in the a3 direction, where we identify a Driac surface cone, expected given the classification as strong TI for that filling. Compression and expansion of caxis Our calculations also show that variation of the a3 length is also an effective method to induce topological transitions, while preserving the tetragonal symmetry of the undisturbed structure. In particular, we predict two ways of achieving two different topological phases. First, compression of the c axis by 3% (Figure 8.9a) drives the system into a weak TI phase at a filling of 24 (red bands), corresponding to 1.3 GPa of uniaxial pressure. This is because the first 16 bands are trivial while the next 8 correspond to a branch of a split EBR (SEBR) and thus cannot be an EBR by themselves. Computing the SIs yields z4= 2 , which implies that this is a weak TI, and z2w,3= 1 , so the phase is equivalent to a series of layers of Quantum Spin Hall insulators (QSHIs) stacked in the a3 direction. Thus, this structure should host topological surface states only on the faces perpendicular to a1 and a2 only. A second way to obtain a topological phase is by expanding the c axis by as little as 1% (Figure 8.9b), which 139 REGISTRO TELEMÁTICO Sarreren Erregistro Orokorra / Registro General de Entradas 22/05/2025 07:32 EHU2025E021341