The ageing effect in topological insulators: evolution of the surface electronic structure of Bi2Se3 upon K adsorption

Topological insulators (TIs) have attracted a lot of interest in recent years due to their topologically protected surface states, as well as exotic proximity-induced phenomena and device applications for TI heterostructures. Since the first experimental studies of TIs, angle-resolved photoemission spectra (ARPES) showed that the electronic structure of the topological surface states significantly changes as a function of time after cleavage. The origin and underlying mechanism of this ageing effect are still under debate, despite its importance. Here we investigate the evolution of the surface Dirac cone for Bi2Se3 films upon asymmetric potassium (K) adsorption, using density-functional theory and a tight-binding model. We find that the K adatoms induce short-ranged downward band bending within 2–3 nm from the surface, due to charge transfer from the adatoms to the TI. These findings are in contrast to earlier proposals in the literature, that propose a long-ranged downward band bending up to 15 nm from the surface. Furthermore, as the charge transfer increases, we find that a new Dirac cone, localized slightly deeper into the TI than the original one, appears at the K-adsorbed surface, originating from strong Rashba-split conduction-band states. Our results suggest possible reinterpretations of experiments because the new Dirac cone might have been observed in ARPES measurements instead of the original one that appears immediately after cleavage. Our findings are consistent with ARPES data and provide insight into building TI-heterostructure devices by varying the band-bending potential or film thickness.

reinterpretations of experiments because the new Dirac cone might have been observed in ARPES measurements instead of the original one that appears immediately after cleavage. Our findings are consistent with ARPES data and provide insight into building TI-heterostructure devices by varying the bandbending potential or film thickness.

Introduction
Topological insulators (TIs) have attracted great attention due to topologically protected surface states [1][2][3][4][5] and possible applications in spintronics and thermoelectrics. Exotic phenomena such as the quantum anomalous Hall effect, the topological magneto-electric effect and the formation of Majorana modes, were proposed and explored for TIs interfaced with other types of materials [3][4][5][6][7][8]. The topologically protected surface states bridge the bulk band gap, induced by a band inversion due to strong spin-orbit coupling (SOC). Moreover, the surface states have a helical Dirac dispersion near the Fermi level which is crossed only an odd number of times between time reversal invariant k points [1][2][3][4][5]. Consequently, elastic backscattering is strongly suppressed at the surface, because spin-flip processes are forbidden as long as time-reversal symmetry is preserved. These properties were observed in angle-resolved photoemission spectra (ARPES) on Bi-based alloys and crystals [4,5,9,10]. Among those materials, Bi 2 Se 3 has been studied the most because it has a single Dirac cone at the surface, with the largest bulk band gap.
Since the first experimental studies of TIs, it was reported that the electronic structure of the topological surface states significantly changes as a function of time after cleavage [11][12][13][14][15][16][17][18][19][20][21][22][23]. It was observed that the Dirac point shifts to higher binding energies when time passes, indicating more electron-rich surfaces. Furthermore, Rashba-split states appear in the conduction band, together with M-shaped valence-band states. This ageing effect was explained as a consequence of band bending near the interface and different proposals were formulated in the literature for the underlying mechanism. One explanation, motivated by the work of Urazdin et al [24], was an expansion of the van der Waals gap between the quintuple layers (QLs) of the structure after cleavage [25]. However, recent combined experimental and theoretical works [16,17,21,22] showed that the ageing effect is similar to the effect of adsorption, based on experimental studies of controlled atomic and/or molecular adsorption on TIs. The resultant band bending was explained in terms of an electrostatic potential variation near the surface. In these studies, a weak potential that penetrates down to 10-15 QLs (10-15 nm) was needed to explain the experimental results. Surprisingly, so far, an explicit ab initio study of the ageing effect, using the adsorption of atoms or molecules on TIs, has not been made. A deeper understanding of the underlying mechanism behind the band bending should prove useful for the study of the effect of interfaces on the topological surface states in various proposed systems [3][4][5][6][7][8].
Here we present a theoretical study of the evolution of the surface Dirac cone for Bi 2 Se 3 in the presence of asymmetric potassium (K) adsorption, using first-principles and tight-binding calculations. We have chosen K adsorption as a model system for general atomic and molecular adsorption, since a large charge transfer can be expected with K. In this way, small and large charge transfers to Bi 2 Se 3 can be studied by varying the concentration or the distance to the surface of the K adatoms, accordingly. Furthermore, K adsorption on Bi 2 Se 3 was experimentally realized [19,22]. We simulate K atoms adsorbed on the top surface of thin Bi 2 Se 3 films of one QLs using density-functional theory (DFT) including SOC, and study up to 20 QLs using a tight-binding model with an effective band-bending potential.
We show that the experimental findings, such as the appearance of strong Rashba-split conduction-band states and M-shaped valence-band states, can be explained by charge transfer due to K adsorption. Furthermore, our results show, in contrast to earlier proposals [16,17,21,22], that the band bending affects only the top two to three QLs (2-3 nm), instead of 10-15 QLs (10-15 nm). Moreover, the binding energy of the top Dirac cone increases until it drops below the M-shaped valence-band states, while a new Dirac cone appears close to the top surface, developed from the Rashba-split conduction-band states. These results suggest possible reinterpretations of experiments, because the new Dirac cone might have been observed in ARPES measurements, instead of the original one.
To anticipate a little further, our results show that even for a two-QL slab, the topologically protected surface states at the bottom surface can be recovered by K adsorption at the top surface for high K concentration. Our findings suggest the possibility to engineer hybrid TI structures for desirable properties by tuning parameters such as adsorption layers and slab thickness.

Computational methods
We calculate the electronic structure for bulk Bi 2 Se 3 and Bi 2 Se 3 (111) slabs with K adatoms on the top surface using DFT including SOC self-consistently. We use the VASP code [26,27] with the projector-augmented-wave method [28], within the generalized-gradient approximation [29]. We find that the bulk band gap is 0.28 eV, in agreement with experimental values of 0.3-0.35 eV [30]. The DFT calculations are limited to six QLs, because the hybridization between top and bottom surfaces is already very small for five QLs. This is clear from the hybridization gap at for pristine slabs of five to six QLs, which is found to be of the order of meV, as shown in table 1 as a function of slab thickness. The optimal K position for one to three QLs with a given surface coverage , is found by relaxing the topmost QL while keeping the other QLs fixed at the experimental value [31]. For slabs with four to six QLs, the optimized K position of the three-QL slab is used, because of the high computational cost. We find that the optimal vertical distance between the K adatoms and the slab is 2.17 Å for all slabs, within the numerical accuracy. Compared to hcp, bridge, and on-top adsorption sites, we find that fcc adsorption sites are the most stable. For thicker slabs, we rely on a tight-binding model. In this model, we consider up to third nearest-neighbor hopping between the s and p orbitals of Bi and Se, and include SOC as an on-site contribution. We use the lattice constants and tight-binding parameters from [32] which were obtained from fitting the bulk bands to DFT results. The hopping matrices are expressed using the method of Slater and Koster [33]. More computational details for both the DFT and tight-binding model can be found in the supporting information 4 .

Results and discussion
The calculated band structures for K-adsorbed slabs of one to six QLs with a full surface coverage( =1)areshowninfigure1.Resultsforlowersurfacecoverage( =1/4)aregiven inthesupportinginformation(seefootnote4).Firstnotethatthedoubledegeneracypresent for pristine Bi 2 Se 3 slabs is lifted due to broken inversion symmetry, except at the time reversal invariant points such as and M.
This section is further organized as follows. First, we discuss the identification of the surface states (section 3.1) and the overall features of the calculated band structures (section 3.2). Next, we investigate the effect of charge transfer between the K adsorption layer and the Bi 2 Se 3 slab (section 3.3). Then we propose a band-bending mechanism (section 3.4) to explain the evolution of the surface and quantum-well (QW) states as a function of the amount of charge transfer (section 3.5). Finally, we compare our calculated results with reported ARPES data, and propose a possible reinterpretation of the experimental results (section 3.6).

Identification of the surface states
The surface states are identified by projecting the electron densities on the individual atomic layers. In this way, we can clearly distinguish surface states from bulk-like states [34]. This is explicitly shown for the fully K-adsorbed five-QL slab in figures 2(c) and (d) for ten states at the point. These states are labeled a1-a10, and are marked in the corresponding band structure, shown in figure 2(a). We proceed by discussing the characteristics of these states. We assign states a4 and a10 as surface states from the bottom (bss) and top (tss) surfaces, respectively, because they are mostly localized in the bottommost and topmost QL, similarly to the surface states of the bare slab. State a5 clearly represents a K-dominated or impurity band very close to the Dirac point of the bottom surface. The states a1-a3 and a6-a9 in figures 2(a)-(d) are identified as bulk-like or QW states. Although the band crossings a6 and a9 show Dirac-like dispersion near , they are not true surface states because they are mostly localized in the second topmost QL. The corresponding states for the pristine five-QL slab are labeled b1-b10 in figure 2(b). The same analysis has been made for two to six QLs upon K adsorption, for which the band structures are shown in figures 1(b)-(f). The projected electron densities are also shown for fully K-adsorbed two to three QLs in figure S4 of the supporting information (available from stacks.iop.org/NJP/15/113031/mmedia).

Overall features of the band structure
Several observations can be made regarding the band structures shown in figure 1. First, we focus on the similar band structures for four, five and six QLs in figures 1(d)-(f). Interestingly, we find that the surface states at both top and bottom surfaces have a gapless Dirac dispersion near . The top-surface Dirac point is shifted downward relative to the bottom-surface Dirac point, showing an increased binding energy of about 0.7 eV. Note that the conduction-band QW states have a large Rashba spin-splitting, which was observed in ARPES measurements on Bi 2 Se 3 with adsorption, surface impurities or aged bare Bi 2 Se 3 [16][17][18][19]. It was also discussed in a DFT study of Bi 2 Se 3 slabs with substitutional impurities [35]. For the thinnest system comprising one QL, a large band gap still persists (figure 1(a)). However, even for a very thin Bi 2 Se 3 slab of two QLs with full K coverage, two distinct Dirac cones are present ( figure 1(b)). This is in sharp contrast to a pristine two-QL slab where hybridization between both surfaces opens a gap. The gapless Dirac cones reappear because the cone that is localized at the top surface is shifted in energy because of the asymmetric K adsorption. This behavior can be explained by a minimal model for the surface states of a TI slab with structural inversion asymmetry [35,36]. Additionally, note that the band gap of the two-QL slab ( figure 1(b)) has increased to about 0.7 eV due to confinement. These results show a possibility of using confinement to increase the band gap, as requested for applications, while preserving the topologically protected surface states by separating them energetically.
Interestingly, for a fully K-adsorbed three-QL slab (figure 1(c)), the surface states at the top surface show a gapless Dirac dispersion, while the surface states at the bottom surface are not well defined at . In this case, strong coupling to bulk-like and K-dominated states causes an apparent energy gap of about 0.1 eV at in the bottom-surface Dirac cone. However, the gap does not open at itself and time-reversal symmetry is still preserved. The coupling strength can be adjusted by varying the coverage, slab thickness or the type of adsorption layer. Away from , the coupling reduces, which is shown in figures S4 and S5 of the supporting information.

Charge transfer: origin of local band bending
The origin of the local band bending is the charge transfer between the K-adsorbed layer and the Bi 2 Se 3 slab, which leads to the creation of a surface dipole. The direction and the amount of charge transfer are obtained from comparing the layer-projected electron density ρ(z) of a fully K-adsorbed Bi 2 Se 3 slab, a bare Bi 2 Se 3 slab and an isolated K layer, using the same unit cell, parameter values and k-point sampling. For this calculation, the geometries of the bare slab and isolated K layer are taken from the fully K-adsorbed slab. This is shown in figure 3(a) for a four-QL slab, as a function of the distance from the bottom surface z. Next, we calculate the electron density difference ρ(z) = ρ(adsorbed slab) − ρ(bare slab) − ρ(isolated K), as shown in figure 3(b). In order to find the total amount of charge transfer, we integrate ρ(z) from the bottommost Se layer, to the midpoint between the topmost Se layer and the K adatom indicated as the vertical dashed line in figure 3. We find that 0.112 electrons per unit cell surface area (14.86 Å 2 ) or 0.0075 e Å −2 are transferred from the K adatoms to the four-QL slab. For two For this calculation, we use the group velocity of the surface states close to which is about 5.0 × 10 5 m s −1 [4,5]. Going up to 0.7 eV, 0.06 electrons per unit cell surface area are needed to fill the cone. This is still smaller than the calculated charge transfer of the order of 0.1 electrons per unit cell surface area. However, while the top-surface Dirac cone has moved downwards, conduction-band states were also pulled down. Some of the conduction-band states fall below the Fermi level and become populated with a small fraction of electrons. This shows that a small charge transfer can indeed lead to a strong increase of the binding energy of the top-surface Dirac point, and an occupation of conduction-band states.

Band-bending mechanism
Based on the above results, we propose our band-bending scheme. In figure 4(a), we show the band diagram for a semiconductor-metal junction (Schottky barrier) before charge transfer, where the metal is given here by the adsorbed K layer. Note that the chemical potential of Here the metal is given by the adsorbed K layer. The chemical potential or charge neutrality point of the semiconductor, the K layer, and the surface, are given by µ CN , µ K and µ S , respectively. The charge-transfer range is given by δ and L is the system size. (c) The combination of the upward band bending caused by the n-type doping with a wide depletion region of depth z 0 , and the downward band bending due to the K layer over a narrow charge-transfer region of depth δ. (d) Band diagram of the downward band-bending region with the QW states and the bottom (blue) and top (red) surface Dirac points for a five-QL slab with a full K coverage. the K layer µ K is not necessarily aligned with the surface chemical potential µ S or the charge neutrality point µ CN of the slab. Indeed, the calculated work function for a free-standing K layer equals 2.39 eV, while the work function for a bare four-QL slab is 5.56 eV. In our case, this means that µ K is well above the bulk conduction band of the slab.
If µ K > µ S , electrons from the K layer are transferred to vacant surface states of the slab. After the charge transfer, a dipole layer is formed at the interface and µ K is aligned with µ S , as illustrated in figure 4(b). The dipole points toward the bottom surface and is strongly localized at the interface within the topmost two to three QLs denoted by δ in figure 4. This result is consistent with the observation that the energy difference E d between the Dirac points of the bottom and top surfaces does not depend on slab thickness for slabs thicker than three QLs (see table 1).
Experimentally, it has been known that Bi 2 Se 3 is naturally n-type doped. In the Schottky model, electron doping in semiconductors leads to upward band bending at the surface, as the surface states become occupied. This creates a depletion region into the bulk, with a characteristic depth z 0 of several tens of nanometers. On the other hand, as shown in figure 4(b), the charge transfer from the K adatoms to the surface leads to downward band bending in a charge region near the surface of depth δ, which is only a few nanometers wide. The combination of these two effects is illustrated in figure 4(c). In the Schottky model, this short-ranged downward band bending is ignored because z 0 δ for typical semiconductors. In the proposed band-bending schemes from [16,17,21,22], the downward band-bending region is much wider up to 15 nm. Therefore, we conclude that these schemes are not relevant to explain the ageing effect of Bi 2 Se 3 slabs. Figure 4(d) shows the final situation for the fully K-adsorbed five-QL slab. The small downward band-bending region leads to the formation of strongly localized QW levels as demonstrated by the layer-projected electron densities in figures 2(c) and (d). In this case, the lowest-energy conduction-band state (QWS1) has a lower energy than the bottom-surface Dirac point (bss) due to the strong downward band bending. Interestingly, the top-surface Dirac point (tss) appears below the lowest-energy valence-band level QWS3'. Notice that in all previous studies [16][17][18]21], it was always assumed that the top-surface Dirac point lies higher in energy than the valence-band QW states, originating from a much shallower band bending but penetrating much deeper than our scheme. Moreover, in ARPES experiments, another unidentified surface state was found below the the valence-band QW states [17,18]. We suggest that this surface state is the original Dirac point of the top surface. As pointed out in the next section, the hybridization with valence-band QW states causes the top-surface Dirac cone to become less localized at the surface, making it harder to track with a surface technique such as ARPES.
Lastly, we discuss the effect of the downward band-bending in the thinner slabs of one to three QLs, whose band structures are shown in figures 1(a)-(c), respectively. The persistence of a band gap for one QL, is understood from the fact that the slab is so thin that the charge transfer affects both surfaces equally enough, still resulting in strong overlap between them. For two to three QLs, the binding energy of the bottom-surface Dirac cone increases by around 0.4 and 0.1 eV, respectively. This result is consistent with the extent of charge transfer and the band-bending depth. For the three-QL slab, the slab thickness is comparable to the size of the charge-transfer region, which explains the hybridization between the bottom surface state and QWS1, discussed earlier.

Ageing effect: evolution of surface and quantum-well states as a function of charge transfer
The motivation to study K adsorption is partly based on the large charge transfer that is expected for K. A smaller charge transfer can be realized by increasing the distance of the K atoms to the surface. This also simulates the experimental ageing effect, where the number of adsorbed atoms, and thus the amount of charge transfer, increases over a time scale of hours. Therefore, by examining the evolution of the electronic band structure near as a function of the vertical distance between the K layer and the slab, we can investigate the ageing effect. We verified that results obtained in this way are similar to those for lower K concentrations, realized by taking a larger supercell. More details on this subject are included in the supporting information(availablefromstacks.iop.org/NJP/15/113031/mmedia)(seefootnote4).Duetothe high computational cost, we present only the result for the K-adsorbed four QLs using DFT, and we use a tight-binding model with an effective band-bending potential for thicker slabs.
First we discuss the DFT calculations for the four-QL slabs. We consider an optimized K-adsorbed slab, a bare slab and slabs with K adatoms positioned at 4, 6 and 10 Å above the top surface, giving different amounts of charge transfer. In these calculations, we use the experimental lattice constants, i.e. we do not relax the structure, except for the optimized slab. We find that 0.087, 0.054 and 0.025 electrons per unit cell surface area are transferred to the slab for K adatoms at 4, 6 and 10 Å above the surface, respectively. For the optimized K-adsorbed slab, we find a charge transfer of 0.112 electrons per unit cell surface area. Here we use the same procedure to find the charge transfer as in section 3.3.
We also investigate the band bending for a deeper but weaker potential with the tightbinding model, as proposed in the literature. The best fit of the band-bending profile in figure 1(c) of [17] to the effective potential, is obtained for w = 8 nm. The evolution of the band structure for V B from 0 to 0.5 eV is shown in movie 3 of the supplementary information (available from stacks.iop.org/NJP/15/113031/mmedia). Now there is no mixing between the top-surface Dirac cone and the valence-band states, and the Rashba splitting of the conductionband states is very weak. Additionally, the lowest-energy conduction-band states do not form a second Dirac cone as V B increases. This is also clear from the evolution at , which is shown in figure S3 of the supporting information. These observations are not in agreement with experimental results [16][17][18][19].

Comparison to experimental angle-resolved photoemission spectra data
Finally, we compare our results with ARPES measurements [17,19] on Bi 2 Se 3 with nonmagnetic adsorption. ARPES is a surface technique, and therefore we project the electron density onto the top three QLs for the K-adsorbed slabs thicker than three QLs. In this section, we first discuss our DFT results for slabs of four to six QLs (figure 6) and then those of the tight-binding model for the 20-QL slab ( figure 5(b)). For four and five QLs, the Dirac dispersion is gapless, while for six QLs there is an apparent energy gap in the projected band structure shown in figure 6(f). This gap is seen more clearly in figure 6(i) where we plot the projected band structure together with the original band structure, close to the bottom-surface Dirac cone. The entire band structure of the optimized K-adsorbed six-QL slab is shown in figure 7(a), together with the layer-projected densities in figures 7(b) and (c). The gap between states c5 and c6 in figures 6(i) and 7(a) appears due to coupling between the QWS1 and QWS1' at . This coupling causes the states c5 and c6 to become delocalized over the slab at ( figure 7(b)). Away from the coupling is reduced, as shown by the projected electron densities for states d1 and d2 in figure 7(c).
We demonstrate that the short-ranged band bending explains most of the important features observed in the ARPES experiments [17,19]. For example, the situation in figure 6(c) corresponds to figure 1(c) of [17]. The situation presented in figure 6(d) was not reached in thisexperimentalstudy.Thiscorrespondenceisalsosupportedbyourtight-bindingmodel calculationsforthe20-QLslab.Inmovie2ofthesupplementaryinformation(availablefrom stacks.iop.org/NJP/15/113031/mmedia)(seefootnote4),thesamebandstructureevolutionasin movie1isshown,butnowprojectedonthetopthreeQLs.AsthepotentialdepthV B increases, thetop-surfaceDiracconehybridizeswiththeM-shapedvalence-bandstates,andthelowest Rashba-splitstatesmovedownwardwhilebecomingincreasinglymorelocalizedinthetop three QLs. When the Rashba-split states reach the top of the valence band, they give rise to a new topological Dirac cone. From the DFT results, we find that the new Dirac cone does not have strict surface-state characteristics because the electron density is mostly localized in the second QL from the top, instead of the topmost QL. However, the spin structure of the new Dirac cone is similar to that of the original top-surface Dirac cone which is opposite to the spin structure of the bottom-surface Dirac cone, corresponding to the scheme shown in figure 6(h). Several ARPES experiments measured the evolution of the energy states at as a function of time after cleavage. Examples are figure 1 in [17] and figure 8 in [22]. We plot the same evolution as a function of the band bending depth V B for the 20-QL slab in figure 5(b). Only states of which 50% or more are localized in the top three QLs are shown. The hybridization between the top-surface Dirac cone and the valence-band QW levels is clearly visible. Similar features can be recognized in the experimental data mentioned above.

Conclusion
We examined the evolution of surface states of Bi 2 Se 3 films with asymmetric K adsorption as a function of film thickness and the amount of charge transfer, using DFT and a tight-binding model, in order to investigate the experimentally observed ageing effect. The asymmetric K adsorption creates well-separated top and bottom surface states with gapless Dirac dispersion even for a very thin slab of two QLs, as a consequence of charge transfer from the K adsorption layer to the Bi 2 Se 3 slab. Based on our calculations, we showed that the band-bending potential caused by the charge transfer has an unusually short range of 2-3 nm, in contrast to earlier proposals in the literature. Applying this short-ranged band-bending potential, we found that, as charge transfer increases, the top-surface Dirac cone moves downward in energy and can be strongly coupled to QW states. In addition, we were able to show the experimental observations of strongly Rashba-split conduction-band and M-shaped valence-band states. Furthermore, for large charge transfer, we found a new topological Dirac cone at the top surface, originating from the Rashba-split conduction-band states that moved downwards. This result suggests possible reinterpretations of experiments because this new Dirac cone might have been observed by ARPES measurements instead of the original one.
One caveat of our DFT calculations is that we considered periodic arrangements of K adatoms, while periodicity in K adatoms may be difficult to achieve experimentally. However, major qualitative features of our study still remain relevant to thin TI slabs with different types of adsorption layers. A good agreement between the DFT and the tight-binding calculations, supports the validity of our approach. The key elements in engineering surface states of TIs with asymmetric surfaces, are the amount of charge transfer at the adsorbed surface, the hybridization between the top and bottom surface, and the bulk band gap value.