Kondo holes in topological Kondo insulators: Spectral properties and surface quasiparticle interference
A fascinating type of symmetry-protected topological states of matter are topological Kondo insulators, where insulating behavior arises from Kondo screening of localized moments via conduction electrons, and non-trivial topology emerges from the structure of the hybridization between the local-moment and conduction bands. Here we study the physics of Kondo holes, i.e., missing local moments, in three-dimensional topological Kondo insulators, using a self-consistent real-space mean-field theory. Such Kondo holes quite generically induce in-gap states which, for Kondo holes at or near the surface, hybridize with the topological surface state. In particular, we study the surface-state quasiparticle interference (QPI) induced by a dilute concentration of surface Kondo holes and compare this to QPI from conventional potential scatterers. We treat both strong and weak topological-insulator phases and, for the latter, specifically discuss the contributions to QPI from inter-Dirac-cone scattering.
In the exciting field of topological insulators,(1); (2) topological Kondo insulators (TKIs) play a particularly interesting role: in these strongly correlated electron systems, theoretically proposed in Refs. (3); (4), a topologically non-trivial bandstructure emerges at low energies and temperatures from the Kondo screening of -electron local moments due to the specific form of the hybridization between conduction and electrons. As with standard topological insulators, TKIs can exist in both two and three space dimensions, in the latter case as strong and weak topological insulators. TKIs display helical low-energy surface states which are expected to be heavy in the heavy-fermion sense, i.e., with a strong mass renormalization and a small quasiparticle weight.
The material SmB was proposed to be a three-dimensional (3D) TKI,(3); (4); (5) and a number of recent experiments appear to support this hypothesis: transport studies have been interpreted in terms of quantized surface transport,(6) quantum oscillation measurements indicate the presence of a two-dimensional Dirac state,(7) and results from photoemission measurements(8); (9) and scanning tunneling spectroscopy (STS)(10) appear consistent with this assertion. However, to date the topological nature of the surface states of SmB has not been unambiguously verified. Moreover, it has been suggested(11) that the observed surface metallicity is polarity-driven, raising questions about the proper interpretation of the experimental data. This calls for more detailed studies of the surface-state physics of SmB and other candidate TKI materials.
A powerful probe of the surface electronic structure is Fourier-transform scanning tunneling spectroscopy(12); (13) (FTSTS), applied in recent years, e.g., to both cuprate and iron-pnictide superconductors. Within FTSTS energy-dependent spatial variations of the local density of states (LDOS) are analyzed in terms of quasiparticle interference (QPI), i.e., elastic quasiparticle scattering processes due to impurities. Such experiments were performed on the topological insulators BiSb and BiTe, and the results were found to be consistent with a suppression of backscattering, , due to the spin-momentum locking of the helical surface state.(14); (15); (16); (17); (18)
Theoretically, impurity scattering and QPI on the surface of 3D topological insulators have been studied for lattice models,(19) and within effective surface theories for non-magnetic,(20); (21) magnetic,(20); (22) and Kondo (23); (24) impurities. In these works the surface electrons were assumed to be non-interacting, such that the interplay of impurity and strong-correlation effects, expected to be important for TKIs, has not been covered.
In this paper we aim at closing this gap, by studying the physics of local defects in Anderson lattice models of TKIs. In particular we will focus at local-moment vacancies, so-called Kondo holes.(25); (26) QPI from Kondo holes has been considered before(27); (28) for conventional heavy-fermion metals and has been found to be particularly revealing due to the interplay of defect and Kondo physics. Kondo holes on the surface of TKIs promise to be interesting also because they represent strong scatterers for which the simplest arguments of topological protection no longer apply.(29); (30) Here we will employ a fully self-consistent mean-field description of the Kondo insulator, taking into account the local modification of Kondo screening by defects. Applying this methodology to both the weak topological insulator (WTI) and strong topological insulator (STI) phases, we will calculate the electronic structure and the surface QPI patterns for dilute surface Kondo holes as well as for other types of impurities. We will also present selected results for a finite concentration of Kondo holes.
i.1 Summary of results
Our main results can be summarized as follows. Kondo holes on the surface of TKIs tend to create localized states, which hybridize with surface states. This gives rise to distinct features in the LDOS in the immediate vicinity of the hole, with an energy dependence mainly determined by the degree of particle–hole symmetry breaking: in the present model, the WTI phase occurs closer to the Kondo limit and is less particle–hole asymmetric, such that a strong in-gap resonance occurs. In the STI phase, particle–hole symmetry is strongly broken, and the hole-induced weight in the LDOS is shifted to elevated energies.
As expected, QPI patterns closely reflect the dispersion of surface states and are distinctly different for STI and WTI phases, with one and two surface Dirac cones, respectively. In the STI case the QPI signal close to the Dirac point is weak and weakly momentum-dependent, due to forbidden backscattering within a single Dirac cone. In contrast, the WTI displays a strong and strongly peaked QPI signal arising from intercone scattering.
To gain analytical insights into intercone scattering, we have extended the continuum Born-limit calculation of Ref. (20) to two Dirac cones for non-magnetic impurities, and we also sketch the extension for magnetic ones.
Our comparison of different types of impurities reveals surprisingly strong differences in the resulting QPI patterns, arising from (i) extended scattering regions (as compared to point-like defects) for Kondo holes due to a modification of the Kondo effect in the hole’s vicinity and (ii) real parts of Greens functions entering the QPI signal invalidating the naive joint-density-of-states picture. In turn, this implies that experimental QPI results, in connection with careful modelling, can be used to determine the nature of the underlying scatterers.
For a finite concentration of Kondo holes, we find the expected disorder-induced broadening of the surface states. In the WTI phase, the low-energy resonances hybridize to yield an impurity-induced band.
On a technical level, we note that the Kondo effect is strongly modified both at the surface and near vacancies as compared to the bulk of the system, rendering fully self-consistent calculations necessary for a reasonably accurate description of QPI.
The body of the paper is organized as follows. In Section II we briefly describe the model for TKIs and the type of impurities we studied. Section III summarizes the slave-boson mean-field treatment for the translationally invariant case; its modifications for systems with surfaces and/or impurities are described in Section IV. In particular, we discuss how to efficiently calculate propagators for the case of isolated Kondo holes with fully self-consistent mean-field parameters. Numerical results are shown in the remainder of the paper, starting with the clean system in Section V. The main body of results is given in Section VI for isolated impurities, covering the impurity-induced density of states and the QPI patterns. Finally, Section VII presents single-particle spectra for disordered systems with finite concentration of surface Kondo holes. In Section VIII we present the conclusions of our work.
ii.1 Anderson lattice model for topological Kondo insulator
Our work utilizes a tight-binding lattice model for a three-dimensional topological Kondo insulator. Following Refs. (3); (4), we consider a periodic Anderson lattice on a simple cubic (more precisely, tetragonal) lattice, with the Hamiltonian
in standard notation. The model entails two doubly degenerate orbitals per site, labelled for conduction electrons and for localized -shell electrons, respectively. The index denotes the spin of electrons, while the index corresponds to the pseudo-spin of the electrons. Both hopping and hybridization terms are assumed to be non-zero for pairs of nearest-neighbor sites only, with , , . Note that non-zero hopping is required to yield a finite band gap within the slave-boson approximation described below. We will employ as energy unit unless otherwise noted.
The operators describe the lowest Kramers doublet of the electrons once spin-orbit interaction and crystal-field splitting are taken into account, for details see Ref. (4). Here we choose a situation corresponding to tetragonal symmetry, with the lowest doublet being constituted by states with and ():
where , here denotes the spin degree of freedom for electrons, while , , is , the azimuthal quantum number for angular momentum . Experimentally, this may be realized, e.g., in tetragonal Ce compounds with dominant configuration, provided that the electron resides in the chosen doublet.
We note that other ground-state doublets may be considered; however, to our knowledge, this is the only choice compatible with tetragonal symmetry which grants bulk-insulating behavior in 3D. For example, the doublet used in Ref. (31) and in Ref. (32) (together with ) generates an insulator in 2D, but only a semimetal in 3D, as it provides no hybridization along the third direction. An appealing alternative, more tailored towards SmB, would be to consider a cubic environment.(33) However, this implies a degeneracy of the multiplet to be 4 rather than 2, thus significantly complicating the theoretical analysis. We expect that most of the features we find are generic, i.e., would also apply to the cubic case, but we leave a more detailed study of the latter for future work.
The non-trivial topological behavior of the model is encoded in the hybridization form factor(32) between a electron at site with spin , and an electron at site with pseudo-spin . It is defined (up to a constant factor) by the overlap between their wavefunctions:
Using the explicit expressions of the spherical harmonics
or, absorbing the factor in :
ii.2 Kondo holes and other impurities
Most generally, we will consider random potentials on both the and orbitals, described by a disorder Hamiltonian
The full Hamiltonian is then given by , where we can define site-dependent potentials according to and .
Kondo holes represent sites with missing -orbital degrees of freedom (i.e. non-magnetic ions); they are modelled by and (in practice we use ). We will also consider weak scatterers in either the or the band, described by small non-zero or , respectively.
A large part of the paper is devoted to the study of isolated impurities, where only a single site has non-vanishing or , but Section VII will also consider the case of a finite number (or finite concentration ) of defect sites with non-vanishing or .
Iii Mean-field theory: Translation-invariant case
When translation symmetry holds, the one-body part of Eq. (1) can be Fourier-transformed to yield:
where is a momentum from the first Brillouin zone (BZ) ,, further
and the hybridization is the Fourier transform of Eq. (11):
iii.1 Slave-boson approximation
The approximation is based on taking the limit , i.e., excluding doubly occupied orbitals. The remaining states of the local Hilbert space are represented by auxiliary particles, and , for empty and singly occupied orbitals, respectively, such that . The Hilbert space is constrained by . It is convenient to choose bosonic and fermionic, and to employ a saddle-point approximation . With fluctuations of frozen, the above constraint is imposed in a mean-field fashion using a Lagrange multiplier . takes the bilinear form:
with the total number of electrons and the number of lattice sites. We have introduced the chemical potential as the Lagrange multiplier enforcing the average electron number to be , with the Kondo insulator reached at .
which determine the parameters , , and .
The diagonalization of yields single-particle energies and corresponding eigenvectors , with quantum numbers and . The expectation value of a momentum-diagonal single-particle operator is then given by
where is the Fermi-Dirac distribution function and the temperature. Most of our calculations below are intended to be for ; practically we used to avoid discretization errors.
iii.2 Mean-field phases
Within the slave-boson approximation the model in Eqs. (1,11) has been shown (32) to have four different phases as a function of its parameters. For small , one encounters a decoupled phase with which may be classified as a fractionalized Fermi liquid(39) (FL), i.e., an orbital-selective Mott state. Upon increasing , transitions occur to a WTI phase with topological indexes , a STI phase with , and finally a trivial band-insulating (BI) phase with , see also Fig. 2 below. This is in contrast with the standard Doniach model (40) with on-site hybridisation , which only shows a transition from the decoupled to the trivial insulating phase. We note that more complicated mean-field phase diagrams can arise(41) when introducing second- and third-nearest neighbor hopping into , but no qualitatively new phases appear. Beyond the present mean-field theory, antiferromagnetism can be expected for small , but the phases at larger are likely robust.
As shown in Appendix A, the mean-field Hamiltonian is equivalent to the common cubic-lattice four-band model(42); (43) used in the topological-insulator literature. However, in the presence of boundaries and impurities, the physics of the Kondo-insulator model is richer due to the additional self-consistency conditions.
Iv Real-space mean-field theory
In situations without full translation symmetry the local mean-field parameters and become site-dependent, which requires to formulate the mean-field theory in real space. We shall assume the bare hopping matrix elements , , to be position-independent as in Eq. (1), but we treat the case with arbitrary on-site energies in . Then
The local mean-field equations for and read:
where denotes a nearest neighbor of site . Practically, Eqs. (23,24,25) need to be solved numerically for finite-size systems. For sites with Kondo holes, i.e. no degrees of freedom, we formally set which correctly excludes hopping to these sites.
The chemical potential remains a global parameter controlling the electron concentration . In the thermodynamic limit with fixed, will be insensitive to the existence of surfaces as well as to a finite number of impurities (). This is no longer true for finite systems. In our simulations we will fix to its value determined for the translation-invariant case, which implies that can differ slightly from 2 in the presence of surfaces. The advantage of this protocol is to avoid complications arising from a size-dependent .
To improve accuracy within a finite-size self-consistent calculation, we have employed supercells (equivalent to an average over twisted periodic boundary conditions). Unless noted otherwise, a supercell grid was used.
iv.1 Clean system in slab geometry
Surface states are efficiently modelled in slab systems of size , with open boundary conditions along and periodic boundary conditions along and directions. Then the in-plane momentum remains a good quantum number, and the mean-field parameters depend on only. After a (partial) Fourier transform in the plane, electron operators carry indices and , and Eq. (23) can be written as
where means , and with
The eigenstates of , , carry quantum numbers of in-plane momentum and band index .
iv.2 Single impurity: Embedding procedure
In the presence of impurities, the system becomes fully inhomogeneous such that the real-space equations (23,24,25) need to be solved. This restricts the numerical calculation to relatively small system sizes () which are insufficient to accurately study Friedel oscillation and quasiparticle interference.
Notably, the problem can be simplified for isolated impurities using scattering-matrix techniques. The basic observation is that mean-field parameters parameters and are locally perturbed by each impurity, but these perturbations decay on short length scales. Therefore, to a good approximation, electron scattering off each impurity can be described in terms of small-size scattering regions where and deviate from their bulk values.
Specifically, for a single impurity in the surface layer we employ the following procedure, Fig. 1. We determine and in a fully self-consistent inhomogeneous calculation for a small system of size in slab geometry. This scattering region is then embedded into a much larger system; for computational efficiency we further reduce the size of the scattering region to layers, because layers far from the impurity are only weakly perturbed.
Impurity-induced changes of electron propagators are then calculated using the T-matrix formalism:
where the scattering matrix is determined as
All matrices depend on indexes (, , , if , if ). The interaction matrix is given by
where is the mean-field Hamiltonian with defect and self-consistently determined and , and the mean-field Hamiltonian of the clean slab. , reflecting both the defect and its induced changes of mean-field parameters, is taken to be non-zero only in the small scattering region.
The Green’s function of the impurity-free slab, , is diagonal in in-plane momentum ,
with from Eq. (27), and is an artificial broadening parameter. Fast Fourier transform is used to obtain in real space. Finally, the LDOS is computed through
For weak scatterers, the lowest-order Born approximation is sufficient, .
For these scattering calculations we employed to achieve a high momentum resolution for QPI (see below), combined with a broadening for the WTI and for the STI, to have a smooth energy-dependent DOS. was found to be sufficient to obtain converged results for the surface LDOS.
iv.3 Surface quasiparticle interference
The QPI signal is obtained from the energy-dependent surface LDOS by Fourier transformation in the plane:
where only the impurity-induced change in the density of states is considered,
the homogeneous background would contribute a signal at only. We note that is in general a complex quantity. However, in the case of a single impurity is inversion-symmetric w.r.t. the impurity site, such that is real.
When relating the surface LDOS and to the signal in an actual STM measurement, complications arise from the fact that the differential tunneling current is not simply proportional to the LDOS: and signals are weighted differently, and an interference term is also present. (44); (45); (46); (47); (48) Such corrections can be taken into account, but the required ratio of the different tunneling matrix elements into and orbitals is usually not known. Therefore we refrain from doing so; we anticipate that no qualitative changes to our conclusions would arise, although the energy dependence of features in the tunneling spectra may be modified.
V Results: Clean system
We have investigated both the weak and strong topological-insulator phases of the model Eq. (1), with parameters chosen to obtain sufficiently large Kondo temperature and bulk gap, as to avoid finite-size effects, and a moderate surface-state Fermi velocity, as otherwise surface-state QPI is restricted to a tiny range in momentum space.
v.1 Phase diagram
Fig. 2 shows a zero-temperature mean-field phase diagram, obtained from Eqs. (18-21), as function of the hybridization for the choice , , . All phases except for FL have ; the WTI, STI, and BI phases have been detected by calculating the relevant topological invariants of the mean-field bandstructure;(4); (49); (50) the transitions WTISTI and STIBI can be detected via the closing of the bulk gap.
A general property of the model (1) is that the WTI phase is found deep in the Kondo regime, whereas the STI phase is realized in a regime of stronger valence fluctuations.
v.2 Band structure and surface states
Subjecting the system to open boundary conditions along , metallic states appear on the two (001) surfaces in both the WTI and STI phases. For the WTI we found two Dirac cones at the two inequivalent momenta and of the surface Brillouin zone; the STI has a single Dirac cone at .
We note that, depending on parameters, the surface states may disperse such that constant-energy cuts near the Dirac energy display multiple band crossings in addition to those arising from the Dirac cones, which in turn complicates the QPI analysis. In our choice of parameters we tried to avoid such situations.
Weak topological insulator
Explicit results for the WTI phase have been obtained using , , , . The resulting mean-field band structure is shown in Fig. 3 for both the periodic and slab cases. The surface Dirac cones at momenta and are clearly visible, with the Dirac point at ; for our choice the cones display a tiny finite-size gap due to the coupling between opposite surfaces. We note that the choice of a relatively large value of is dictated by the necessity to obtain a sizeable bulk gap (which is zero when ), in order to have a sufficient energy window in which Dirac cones and the associated QPI can be studied. With our parameters, the bulk gap evaluates to .
Fig. 4 displays the the layer-resolved spectral intensity in the slab case, defined as
with from Eq. (34), illustrating the weight distribution for bulk and surface states.
All quasiparticle states are mixtures of and electrons. This may be quantified by the band- and momentum-dependent peak weight in the -electron spectral function, usually dubbed quasiparticle weight, . We define an energy-dependent quasiparticle weight according to
where and are, respectively, the and contributions to the total density of states, and is a Lorentzian of width . As common for Kondo systems, is small near the Fermi level, as shown in Fig. 3(b). As a result the surface states are primarily of character, with . We note that does not directly correspond to the effective-mass ratio because the bare band is dispersive and hence the -electron self-energy momentum-dependent.
The bulk Kondo temperature is estimated as from the location of the mean-field phase transition where becomes non-zero upon cooling; this transition is well-known to become a crossover upon including corrections beyond mean fields. Let us point out that surface effect on the mean-field parameters are sizeable: In the bulk we have , whereas on the surface and , i.e., Kondo screening is suppressed at the surface. This strongly influences the Dirac-cone velocity: its value is much smaller than which would be obtained from a non-self-consistent slab calculation using bulk values of and . In other words, the interplay of Kondo and surface physics increases the mass of the surface quasiparticles – an effect only captured by fully self-consistent calculations.
Strong topological insulator
For the STI phase we used parameters , , , . This choice (positive , small ) stems from the necessity to have (i) a small Fermi velocity of the surface Dirac cone and (ii) a sizeable energy window around the Dirac point without additional surface states, in order to be able to study Dirac-cone QPI. As a result, is much larger than the bandwidth due to strong valence fluctuations. Despite this, the quasiparticle weight for surface states remains small, . Further we have , , bulk values of , , surface values , and a surface Fermi velocity of . Band structure and layer-resolved intensities are shown in Figs. 5 and 6, with a single Dirac cone at .
Dirac cone spin structure
For a full characterization of the surface states we have analyzed their spin–momentum locking, with details given in Appendix B. We find the STI Dirac cone to be described by the effective Hamiltonian
where is measured from the center of the cone at – this corresponds to the standard situation with spin perpendicular to momentum.(1) The WTI Dirac cones have somewhat different spin structures, described by
where momenta are measured from the centers of the cones at and . This unusual spin-momentum locking should be measurable by spin-polarized photoemission experiments and is illustrated in Figs. 12(a) and 13(a) below.
Vi Results: Dilute defects
Kondo holes in Kondo insulators are known (25); (26) to create a bound state in the gap, or close to the band edge. We have verified that this also applies to Kondo holes in the bulk of a topological Kondo insulator, essentially because a bound state emerges generically(51); (52) from strong scattering, provided that particle–hole symmetry is not too strongly broken, and is protected by the gap. The situation becomes more interesting for a Kondo hole at or near the surface of a topological Kondo insulator, which is metallic, such that the bound state turns into a resonance which can be in principle observed by STM.
Therefore we now consider Kondo holes in either the surface layer or the layer below. To this end, we perform fully self-consistent mean-field calculations for a system of size in slab geometry, typically with and . Sample results for the mean-field parameters in the case of a Kondo hole in the WTI phase are shown in Fig. 7. These inhomogeneous mean-field parameters are then used as input for the T-matrix calculation, as described in Section IV.2 above, to determine the LDOS and the QPI spectra from Kondo holes. These QPI spectra will be compared to those from weak impurities, as the latter allow us to gain some analytical insight.
vi.1 Local density of states (LDOS) for Kondo holes
The surface-layer LDOS, detectable in an STM experiment and obtained from the T-matrix calculation, is shown in Fig. 8 for the weak topological Kondo insulator with parameters as in Section V.2.1. For the case of a Kondo hole in the surface layer, Fig. 8(a), we see that a resonance appears in the bulk gap, and it hybridizes with surface states. It is mainly localized on the four sites surrounding the hole, with a rapid spatial decay. Similar resonances have been predicted in non-Kondo TIs.(53); (29); (54) For comparison, we have also performed a non-self-consistent calculation where the changes of mean-field parameters due to the impurity have been ignored, such that the T-matrix is non-zero on a single site only. The two results differ significantly concerning the energetic position of the resonance, underlining that full self-consistency is important.
For a hole in the second layer, Fig. 8(b), the resonance appears weaker and at higher binding energy, close to the van Hove singularities of surface states. In the surface-layer LDOS, it is visible essentially only on the site above the hole. We note that the energetic location of the resonance depends on microscopic details, i.e., we have also encountered cases with a sharp low-energy resonance for a hole in the second layer.
In contrast, for the strong topological Kondo insulator with parameters as in Section V.2.2 we have not found sharp low-energy resonances for any parameter set investigated. The reason is that the STI phase of model (1) only occurs in the mixed-valence regime, which in turn implies strong particle–hole asymmetry. Since, for scattering in Dirac systems, the resonance energy is a function of both scattering strength and particle–hole asymmetry,(51); (54) determined essentially by , increasing particle–hole asymmetry shifts the resonance of a strong scatterer away from the Dirac energy. For our parameters, we only found minor impurity signatures in the low-energy LDOS, Fig. 9. Impurity-induced changes are visible at higher energies, but spoiled by the influence of bulk states.
When analyzing the impurity-induced changes in the LDOS, , as a function of the distance from the hole, Friedel oscillations with a wavelength can be observed for energies close to , Fig. 10(a). As long as warping effects can be neglected, the decay is isotropic, and proportional to in the WTI phase and to in the STI phase, in agreement with earlier results for graphene (55); (56) and for STIs.(57); (19); (58) At higher energies, when warping effects cannot be neglected, the decay becomes anisotropic, and a strong focusing effect can be observed if a nesting of the Fermi surface can be achieved, Fig. 10(b).
vi.2 QPI from weak impurities
Before analyzing the quasiparticle interference signal caused by a Kondo hole, it is useful to analyze a few QPI properties of weak impurities, where we neglect any impurity-induced changes of mean-field parameters. A comparison of numerical QPI results for both Kondo holes and weak impurities can be found in Figs. 12 and 13 below.
In the following discussion we restrict our attention to energies within the bulk gap. This enables an analytical treatment close to the Dirac points, using the effective surface Hamiltonian Eqs. (40) or (41). This approach has been taken in the literature before, and we start with reviewing these results. In what follows we measure energies relative to the Dirac energy: .
Point-like effective impurity
We focus first on a strong topological insulator with the surface Hamiltonian (40). The unperturbed surface Green’s function is
For scatterers which act as point-like impurities in the effective theory – this assumption is in general not justified, see below – the T-matrix is momentum-independent, and the scattering matrix equation (31) gives, for intracone scattering,
Here we have used that is real for the single-impurity case considered here. For a non-magnetic impurity, both and are proportional to the identity in the spin space, so can be treated as scalar quantities, and we get:
and we have made explicit that describes intracone scattering . This result has been derived before in the context of graphene,(59) STIs,(20) high-temperature superconductors,(60) by going into Matsubara frequencies , applying the Schwinger-Feynman parametrization trick, then coming back to real frequencies. It only depends on the magnitude of the transferred momentum , and, up to a momentum-independent additive term that we omit, reads: