Dirac node engineering and flat bands in doped Dirac materials

We suggest the tried approach of impurity band engineering to produce flat bands and additional nodes in Dirac materials. We show that surface impurities give rise to nearly flat impurity bands close to the Dirac point. The hybridization of the Dirac nodal state induces the splitting of the surface Dirac nodes and the appearance of new nodes at high-symmetry points of the Brillouin zone. The results are robust and not model dependent: the tight-binding calculations are supported by a low-energy effective model of a topological insulator surface state hybridized with an impurity band. Finally, we address the effects of electron-electron interactions between localized electrons on the impurity site. We confirm that the correlation effects, while producing band hybridization and Kondo effect, keep the hybridized band flat. Our findings open up prospects for impurity band engineering of nodal structures and flat-band correlated phases in doped Dirac materials.


I. INTRODUCTION
Impurity band engineering is at the core of modern semiconducting industry where impurity bands enable functionality of a semiconductor. Similarly, electronic structure and topology of Dirac materials 1 (DMs) can be manipulated by impurity doping. A well known example of such manipulation is the quantum anomalous Hall effect (QAHE) 2,3 , which is a new quantum state of matter observed in magnetically-doped three-dimensional (3D) topological insulators (TIs). It occurs as a result of a gap opening at the Dirac node of TIs due to broken timereversal symmetry. It is known that impurities give rise to low-energy resonant states near Dirac nodes [4][5][6] . For a magnetically-doped 3D TI, the magnetic energy gap is filled with impurity resonant states. Hence, disorder effects have significant implications for QAHE 7,8 . In this work, we consider another example of impurity band engineering in DMs. Specifically, we propose to use impurity bands to introduce flat bands and produce additional Dirac nodes in DMs.
Flat bands can occur in a variety of systems, including electronic materials such as superconducting wire networks and engineered 2D atomic lattices, in optical lattices of cold atoms, and in photonic systems such as waveguide arrays and exciton-polariton condensates 9 . Due to the quenching of kinetic energy, flat bands are highly susceptible to interactions. In particular, flat electronic bands are expected to give rise to interactiondriven quantum phases, such as superconductivity and Bose-Einstein condensation. Nearly flat electronic bands can be found in "slow" DMs, in which Dirac states exhibit extremely small velocity, resulting in a large coupling constant 10 . Recently, there has been a renewed interest in flat bands due to the discovery of superconductivity in twisted bilayer graphene near the so-called magic angles which host "slow" Dirac fermions close to the charge-neutrality point [11][12][13] .
Here we consider impurity-induced flat bands that emerge in doped DMs, e.g. in 3D TIs, with periodicallyarranged impurities. In a certain range of impurity potentials, these impurity bands appear near the Dirac nodes and hybridize with the Dirac states. We demonstrate, by using an effective model and tightbinding calculations for a typical 3D TI with a single Dirac node at the Γ point, that when the impurity resonance state is energetically close to the Dirac node, such hybridization results in the splitting of the original Dirac node at Γ and in the appearance of additional nodes at other high-symmetry points in the Brillouin zone.
We also show, using a slave boson approach, how a Dirac conduction band hybridizing with a strongly interacting localized band via a Kondo coupling J also leads to this physics, with a strongly renormalized chemical potential for the local orbitals, provided that J is sufficiently large (on the order of the conduction electron bandwidth).

II. MODEL
Our theoretical modelling is based on the sp 3 Slater-Koster tight-binding (TB) model for Bi 2 Se 3 3D TI, with parameters fitted to ab initio calculations obtained with Wien2k package 14,15 . The TB model for pristine Bi 2 Se 3 includes s and p orbitals and Slater-Koster hopping elements between atoms in the same atomic layer and between atoms in first and second nearest-neighbor layers. The spin-orbit interaction is incorporated in the intraatomic matrix elements.
For surface calculations with impurity doping on the (111) surface, we consider a slab consisting of five quintuple layers (QLs) and a 3 × 3 surface supercell 15 (unless stated otherwise). An impurity substitutes an atom in the topmost Bi layer and is described by a local onsite potential 4-6 . A point-like impurity potential U acts as a uniform shift to the on-site energy of the impurity atom. The impurity potential introduces localized impu-rity states that can affect electronic states in the vicinity of the Dirac node 5 . Microscopic tight-binding model calculations are accompanied by low-energy continuum model calculations for a TI surface state hybridized with a doubly degenerate impurity band.

A. Tight-binding calculation for 3D TI
The results of a representative TB calculation with U = 6 eV are shown in Fig. 1. The supercell structure with a substitutional impurity near the top surface is shown in Fig. 1(a). The spectral function corresponding to atomic layers along the growth direction of the TI slab is calculated by averaging the atom-resolved spectral function over the atoms in a specific layer. The spectral functions of the top and bottom surfaces are shown in Fig. 1(b) and (c), respectively. The presence of surface impurities leads to impurity resonance states that appear as nearly flat bands in the spectral function of the doped surface [ Fig. 1(b)], while the undoped surface exhibits unperturbed Dirac states. Figure 2 shows the calculated band structures and the density of states (DOS) for a range of impurity potentials U . For U = 0, the Dirac states of the top and bottom surfaces of the slab are degenerate. With increasing U , a nearly dispersionless band emerges in the valence band, and the degeneracy between the top and bottom Dirac states is lifted. The Dirac node of the undoped surface remains pinned at the position of the Dirac node of the pristine system, while the Dirac node of the doped surface is shifted in energy. With increasing the impurity potential further, the impurity band crosses the Dirac node of the top surface and shifts further up in energy. For large U , the impurity band merges with the conduction band and the degeneracy of the top and bottom Dirac states is restored.
The impurity states also appear as peaks in the DOS as shown in the right hand-side panels in Fig. 2. Their position is controlled by the non-magnetic impurity potential U and coincides with the position of the flat impurity band.
When the impurity resonance state falls in the vicinity of the node, which corresponds to U = 6 eV in this model [ Fig. 2(c)], the Dirac spectrum is drastically modified. The Dirac node of the doped surface is split into two nodes, displaced vertically in energy, and the surface states hybridize with the flat impurity band. Figure 3 shows the spatial distribution of the wavefunctions of the three nodes at Γ, which are marked as ε i in Fig. 2(c) and correspond to the Dirac nodes of the bottom surface and the split nodes of the top surface. The two new nodes at Γ are predominantly localized at the top surface. However, they are not pure surface states and are hybridized with the impurity band localized on the impurity site.
This feature was noticed in previous theoretical work on 3D TIs, which used the scattering matrix approach with and without disorder to study impurity resonances 7,16 . In contrast to our tight-binding calculations, this approach does not reply on the use of a supercell. The splitting of the Dirac node due to coherent impurity scattering and the appearance of dispersionless band at zero energy was shown in graphene in the presence of vacancies 17 . A similar effect was found in ab initio calculations of TI/normal semiconductor heterostructures, where a semiconductor valence band hybridizes with TI surface states 18 . This result was qualitatively explained by an effective model, based on the low-energy surface state Hamiltonian hybridized with a trivial band from a proximal semiconductor layer. More generally, the possibility of reshaping the topologically protected surface states with localized impurity resonances was demonstrated in 19 using numerical simulations and scanning tunneling microscopy (STM).
Although the splitting of the Dirac nodes in the presence of disorder was noticed in previous work 7, 16 , it has not been studied in the context of topology and impuritycontrolled nodal structure of DMs. To further illustrate the details of this effect, we show in Fig. 4 the calculated band structures for U = 6 eV and for different surface impurity concentration. As one can see from Fig. 4(a), the Dirac node of the top surface is split vertically into two nodes (red arrows) and is hybridized with the impurity band, while the bottom surface states remain unaffected by doping. The splitting between the two nodes at Γ decreases with increasing the doping concentration.
An additional doubly degenerate state appears at the Brillouin zone corner (M points) and is marked by a blue arrow in Fig. 4. Here, we consider non-magnetic doping, hence time-reversal symmetry guarantees Kramers degeneracy at the time-reversal invariant momenta Γ and M . The surface of Bi 2 Se 3 3D TI is a triangular lattice. Due to the symmetry of the corresponding hexagonal Brillouin zone [see the inset in Fig. 4(a)], there are three non-equivalent M points each hosting a doubly degenerate state.
The splitting and generation of new nodes leads to a natural question of whether the topology of the surface states is somehow affected by non-magnetic impurities. As confirmed by the calculations, the nodes remain gapless; however, their number and position change. We verified that the topological properties are preserved in the presence of non-magnetic doping despite the modified nodal structure. The Dirac states in a 3D TI are characterized by spin-momentum locking in the vicinity of the node. We define the helicity h as the eigenvalue of the helicity operator σ ·k, wherek = k/|k|. This can be visualized as the direction of rotation of the spin of an energy eigenstate as the momentum changes clockwise from +k y to −k y . We found numerically that the helicity of the doped surface, calculated by summing the helicities of the nodes at Γ and M , is the same as the helicity of the pristine surface. Thus, the peculiar splitting of the surface state node by impurity resonance necessitates the We will illustrate the splitting of the Dirac node by hybridization with an impurity band by using an effective low-energy model. We consider the following Hamiltonian where is the Hamiltonian of a 3D TI surface (or graphene); v F is the Fermi velocity. Here V = γÎ is the coupling matrix, γ is the coupling strength,Î is a 2 × 2 identity matrix. E imp = ε imp (k)Î is the impurity band Hamiltonian, where ε imp (k) is the impurity band dispersion. We consider ε imp (k) = ak 2 − µ imp ; the case a = 0 corresponds to a quadratic band, while a = 0 gives a flat band. The resulting bandstructures are shown in Fig. 5. The eigenvalues of the Hamiltonian in Eq. (1) can be found from the following equation The pair of doubly degenerate states at Γ are given by the solution at k x = k y = 0 For µ imp γ (γ = 0), we have E 1 = 0 and E 2 = µ imp .

C. Hybridization via Kondo coupling
For completeness we also analyze the case of correlated impurity bands. While the presence of the flat bands due to impurities is natural, the question can be asked about stability of our results in the presence of correlations. To address this, we now consider an interacting model of a Kondo lattice. We use an analysis based on a localized, SU(2) × SU(N )-degenerate f -band of electrons coupled to a conduction electron band with Dirac cones. The

Hamiltonian is
where σ is an SU (2)  presumed to be strongly interacting, satisfying the constraint m,σ f † Rmσ f Rmσ = qN f , where q is the fixed filling fraction of the flavor orbitals in each unit cell, which is conserved by H. This constraint is enforced by Lagrange multipliers λ R at each site. Invoking the well-established mean field treatment 21 , the quartic term is decoupled via a Hubbard-Stratonovich transformation with complex local fields V R , and the Kondo term becomes Assuming a mean field solution where λ R = λ and V R = V are spatially uniform, the mean field Hamiltonian becomes At this point, we can work in the diagonal basis of h σσ (k), whose eigenvalues E α (k) are upper and lower bands which touch at certain k values where there are Dirac cones. We also define ε ≡ ε 0 + λ, which is the renormalized f -level energy. For simplicity of calculation, we adopt a model of the c-electron density of states g(E), with g(E) = 2 |E| W −2 Θ W − |E| . At T = 0, the dimensionless free energy per site per flavor, in units of the c-electron half-bandwidth W , is where is the non-fermionic contribution to the energy, are the energies of the hybridized c-f bands (upper and lower), and E * is defined by ξ − (E * ) = 0, where we position our Fermi level.
The mean field solution is obtained by extremizing ϕ with respect to the parameters |V | and ε. We defer the description of the complete solution to a future publication and describe here some limiting results. Writing ε ≡ xW , |V | ≡ √ r ε, we obtain a solution to the mean field equations only when J > J c = 1 2 W . Defining ζ ≡ 1 − W 2J , in the limit where 0 < ζ 1 we obtain the equations where x(ζ) is a solution of the first of these mean field equations. That a critical value of J on the order of the bandwidth is necessary in order to obtain a solution is expected from the work of Fradkin and others on magnetic impurities in Dirac systems 22,23 . In the present context, we find then that for J > J c , strong interaction physics within a localized f -band nominally located below the Fermi level results in a renormalization of the f -electron energy, pushing it up to just above the chemical potential, where the f -band can effectively hybridize with the Dirac c-band and yield the nearly flat bands discussed previously.

IV. DISCUSSION AND CONCLUSIONS
Although it can be problematic to achieve a regular lattice of impurities on a TI surface, there is a strong evidence for impurity-induced states in typical doped TI samples 24 . Moreover, recent ARPES studies suggested that overlapping of impurity resonances with the Dirac node, a situation depicted in Fig. 1(b), could be achieved in the experiment as the binding energy of the node changes with increasing the film thickness 25 . It is also possible to resolve, at least partially, the impurity bands by adjusting the photon energy of the pulse. Another promising platform that can be used to study impurity flat bands are artificially grown impurity superlattices on graphene. Such lattices can be realized by self-assembly of organic molecules on graphene deposited on the substrate 26 .
In Fig. 1(b), we considered a situation where material parameters are such that the flat bands occur exactly at the Dirac point of the pristine material. Assuming that the chemical potential is at the Dirac node, this is the most favorable situation which allows access to flat bands. In the new-generation of 3D TIs with composition (Bi,Sb) 2 B 3 (B=Se,Te), chemical potential can be positioned within 10 meV from the Dirac point. The location of the impurity resonances varies greatly with the type of material and dopant 24 . Recent experiments indicate that favourable conditions can be achieved at least for some samples and dopants, e.g. for Cr in Bi 2 Se 3 above the magnetic ordering temperature 25 .
However, in the majority of materials, the above conditions will not be satisfied. Further complications arise from the fact that impurity doping itself may change the position of the Fermi level. For such situations, optical pumping can be used to populate the flat bands. Superconducting or excitonic pairing may occur between carriers residing in the flat bands and will have a transient nature 27,28 . These ideas can be generalized to pumping of flat, or nearly-flat bands, which do not necessarily originate from impurity-induced states. One fascinating system is a magnetic Weyl semimetal such as the recently discovered Co 3 Sn 2 S 2 , in which the band connecting the Weyl nodes is flattened due to correlations of 3d Co electrons 29 .
In summary, we showed that flat impurity bands arise in impurity-doped Dirac materials. We demonstrate this impurity assisted band structure engineering for a specific case of a three-dimensional topological insulator with non-magnetic impurities on the surface. The impurity flat bands hybridize with the Dirac states and modify their nodal structure. Peculiar features, such as the splitting of the Dirac nodes and generation of additional nodes, are predicted and are explained by a low-energy effective model and topological band theory. We argue that these features are not model specific and are not the artefacts of the supercell approach, and are also present in disordered systems. To test the robustness of these predictions in the presence of electron correlation effects we applied the large-N SU(2) × SU(N ) hybridization model to probe the mean field bands and Kondo effect. We find that, while the bands are renormalized, the "extra" Dirac nodes and flat bands persist in interacting model provided the coupling is sufficiently large. These results reinforce our predictions of flat bands engineered via impurity bands and render them experimentally feasible. Predicted impurity-engineered flat bands present a natural platform for exploring interaction-induced phases, beyond the paradigm of twisted bilayer graphene.