Quantum spin Hall insulators in centrosymmetric thin films composed from topologically trivial BiTeI trilayers

The quantum spin Hall insulators predicted ten years ago and now experimentally observed are instrumental for a break- through in nanoelectronics due to non-dissipative spin-polarized electron transport through their edges. For this transport to persist at normal conditions, the insulators should possess a sufficiently large band gap in a stable topological phase. Here, we theoretically show that quantum spin Hall insulators can be realized in ultra-thin films constructed from a trivial band insulator with strong spin-orbit coupling. The thinnest film with an inverted gap large enough for practical applications is a centrosymmetric sextuple layer built out of two inversely stacked non-centrosymmetric BiTeI trilayers. This nontrivial sextuple layer turns out to be the structure element of an artificially designed strong three-dimensional topological insulator Bi2Te2I2. We reveal general principles of how a topological insulator can be composed from the structure elements of the BiTeX family (X = I, Br, Cl), which opens new perspectives towards engineering of topological phases.

Two-dimensional (2D) topological insulators (TIs)-a new electronic phase also referred to as a quantum spin Hall (QSH) insulator-are characterized by an absolute band gap induced by spin-orbit coupling (SOC) and helical gapless edge states inside the gap [1].These states protected by timereversal symmetry provide perfectly conducting spin-filtered channels, meeting the demands of low-power nanoelectronics and spintronics.The existence of such states as the fingerprint of a topologically nontrivial 2D insulator was first predicted in Refs.[2,3].It was also suggested that the QSH effect can be observed in graphene, where SOC opens a gap at the two inequivalent Dirac points.This gap in graphene appears to be too small for practical use, so heavy-elements based analogs of graphene must be sought.Actually, the 2D materials with low-buckled honeycomb-lattice structures [4,5,6,7,8,9, 10]-silicene, germanene, and stanene-possess a significantly larger SOC-induced gap at the Dirac points (up to ∼ 0.1 eV in stanene), and the spin-polarized edge channels could be detected at easily accessible temperatures.However, the QSH effect in such systems has not been experimentally observed so far.
Further effective enhancement of the SOC to make the gap larger can be realized by chemical functionalization of the above 2D materials [11,12].Such a functionalization substantially enlarges the gap, in fact "destroying" the Dirac cones, and it may lead to a SOC-induced band inversion at the time reversal invariant momentum (TRIM) (normally at k=0) with an absolute gap of several hundred meV at this momentum.If the band inversion occurs, the resulting 2D system is a 2D TI that should support the QSH effect [12] similar to the inverted HgTe quantum wells predicted to be QSH insulators in Ref. [13].It is important that this prediction has found experimental confirmation: In the inverted HgTe/CdTe and InAs/GaSb quantum wells [14,16,15] the QSH effect was observed despite the very small gaps in these quantum wells, less than 10 meV.It has spurred a rising tide of theoretical propositions of different 2D TIs with honeycomb-or square-lattice structures and a large inverted gap enabling room-temperature operating [17].
Ab initio approaches to electronic structure, especially those based on the density functional theory (DFT), have become a powerful tool to search for new materials with unique properties.At the same time, the effective models that proved indispensable in predicting the QSH effect in graphene-like systems and quantum wells are currently widely used to analyze the effect of strain, quantum confinement, and external fields in 2D TIs, i.e., to solve the problems that presently are not accessible with ab initio methods.Thus, to efficiently model the nanoelectronics and spintronics devices, the microscopic methodology must be bridged with the effective Hamiltonian approach based on symmetry considerations and on the k•p perturbation theory.
With a few exceptions, none of the theoretically proposed 2D materials has been hitherto fabricated [17].Thus, the intensive search of robust and easily fabricated materials remains to be actual.In particular, it was suggested that 2D TIs can be produced from a thin film of layered 3D (three-dimensional) TIs of the Bi 2 Se 3 family, where the hybridization between the opposite surfaces of the film opens a gap at the Dirac point (DP).Depending on film thickness, the 2D system may "oscillate" between band insulator and QSH insulator as was predicted by the 4-band effective k•p model (see Refs. [18,19]).The thinnest known topologically nontrivial film consists of at least two structural elements-quintuple or septuple layers [20].
Besides the studies on the thin films of 3D TIs, recently it was heuristically suggested that a 3D TI can be constructed artificially via stacking 2D bilayers that are topologically trivial [21].It encourages our search for 2D TIs built out of trivial band-insulator constituents.These constituents should have strong spin-orbit coupling (SOC), and prospective candidates are bismuth tellurohalides BiTeX with X =I, Br, and Cl, among which the polar semiconductor BiTeI demonstrates the strongest spin-orbit coupling providing the biggest known Rashba spin-splitting of bulk and surface states [22,23].The structure element of BiTeI is a trilayer (TL) with the I-Bi-Te stacking.A single TL that possesses the Rashba spin-split band structure [24] can be grown epitaxially on a suitable substrate or be easily exfoliated from the bulk BiTeI, where the adjacent TLs couple through a weak van-der-Waals (vdW) interaction.The samples of BiTeI always contain a large number of randomly distributed bulk stacking faults, which leads to a mixture of terminations at the surface, as experimentally observed in Refs.[25,28,27,26,29,30].This implies that adjacent TLs may have different sequence order along the hexagonal z axis.
Here, based on DFT calculations we demonstrate that a centrosymmetric sextuple layer (SL) constructed from two BiTeI TLs with facing Te-layer sides and a typical vdW spacing is a 2D TI with the gap of 70 meV at Γ.The vdW interaction between these TLs is crucial to realize such a QSH insulator phase: The SL becomes topologically trivial with increasing the vdW spacing by 5% only.We consider the nontrivial SL as a structure element, a repetition of which along the z axis results in thin films that are found to "oscillate" between trivial and nontrivial phases with the number of SLs.The corresponding bulk system composed of SLs turns out to be a strong 3D TI (hereafter referred to as Bi 2 Te 2 I 2 ).It is energetically unfavourable by only 0.5 meV compared with the non-centrosymmetric BiTeI.This makes it plausible to suppose that crystals of BiTeI grown by the Bridgman method already contain the desired SLs, and that the alternative stacking can be experimentally observed and controllably manufactured.To describe the low-energy properties of Bi 2 Te 2 I 2 and its films, we derive four-band k•p Hamiltonians from the ab initio wave functions.They are similar to the Hamiltonians constructed for Bi 2 Se 3 -family 3D TIs and their thin films [31,18,32].For a more accurate description of the SL, we derive an eightband Hamiltonian that involves Rashba-split valence and conduction bands of the stand-alone TLs.We thus demonstrate that due to the bonding-antibonding splitting the inversion occurs between one of the Te-related valence bands and one of the conduction bands formed by Bi orbitals.The proposed materials illustrate the effectiveness of the new way to design 2D TIs from trivial band insulators with giant-Rashba spit bands for room-temperature operating.

Results
Figures 1(a)-1(c) show the band structure of 1 and 5 SLs and the bulk crystal of Bi 2 Te 2 I 2 obtained with the extended linearized augmented plane wave (ELAPW) method [33] within the local density approximation (LDA) for the exchange-correlation functional and with the use of the full potential scheme of Ref. [34].(Details on the equilibrium bulk atomic structure, the bulk-truncated slab geometry of the related thin films, and the calculations performed can be found in Supplementary Note 1.) The 1SL film is constructed from two BiTeI trilayers with facing Te-layer sides, Fig. 1(d).It is noteworthy that the band structure of this film with the gap of 56 meV (70 meV in the relaxed geometry, see Supplementary Fig. 1) differs substantially from that of its constituents (cf.Supplementary Fig. 2): there is no trace of Rashba-type split bands.The band structure of the 5SL film exhibits a gapless Dirac state residing in the band gap of 151 meV, see Fig. 1(b) and Supplementary Fig. 3.This is a signature of the topological character of the respective bulk band structure (Fig. 1(c)), which has an inverted gap of 234 meV at Γ and the fundamental gap of 169 meV in the Γ−A line close to the A point.As seen in Fig. 2(a), the Dirac surface state almost completely resides within the outer SL.Moreover, this state is localized stronger than the Dirac state of TIs like Bi 2 Te 3 , since 70% of its weight falls in the outermost half of the SL, i.e., in the surface BiTeI TL (see also Fig. 1(e)).
The spin texture of the Dirac state is illustrated in Figures 2(b) and 2(c), which show spin-resolved constant energy contours for the lower and upper cones of the Dirac surface state.Apart from the in-plane polarization -clockwise above the DP and counterclockwise below it, both contours also have an out-of-plane spin component, which is an intrinsic feature of the hexagonal surface.However, in this case S z is extremely small and varies in the range of ≈ ±0.01.
As has been shown for the Dirac state in Bi 2 Se 3 both experimentally and theoretically [35,36,37,38,39], the spin textures of p x , p y , and p z orbitals are remarkably different, which leads to the dependence of the spin polarization of photoelectrons on the polarization of light.The spin texture provided by p z orbitals has clockwise (counterclockwise) chirality for the upper (lower) cone, while the projections of the total spin on p x and p y orbitals are not chiral, and their spins are opposite to each other.Similar spin-orbital texture we find in the Bi 2 Te 2 I 2 , see Figures 2(d)-2(f) for the upper Dirac cone (for the lower Dirac cone, the coupling of spin and orbital textures is opposite, not shown).As can be seen, the spin orientations for p x and p y projections are antiparallel at each k point, whereas the spin orientation of p z projection coincides with the total spin.
It is noteworthy that similar spin-orbital texture has been observed for the spin-polarized Rashba state in BiTeI [40].For the Te-terminated surface, it was found that the outer Rashba branch demonstrates the same spin orientations for p x , p y , and p z projections as those for the upper Dirac cone in Bi 2 Te 2 I 2 , and the inner Rashba branch has opposite texture, i.e., the same as in the lower cone.Because the surface of the Bi 2 Te 2 I 2 slab has iodine termination and its spin-texture is reversed due to the opposite orientation of the z axis, the spin-orbit texture of the upper (lower) Dirac cone in Bi 2 Te 2 I 2 is the same as the texture of the inner (outer) Rashba branch in BiTeI.
To construct a simple effective k•p model for the centrosymmetric Bi 2 Te 2 I 2 we derive a model Hamiltonian of a desired dimension and accurate up to the second order in k from the LDA spinor wave functions Ψ n↑(↓) of the doubly degenerate bands E n found at k=0 (see Supplementary Note 2 and Ref. [41] for details.The subscripts ↑ or ↓ in Ψ n↑(↓) refer to the z-component of the total angular momentum J = L + S in the atomic sphere that has the largest weight in the n-th band, see Fig. 1(e).The Hamiltonian is constructed in terms of the matrix elements [42] where n and m run over the relativistic bands (from semi-core levels up to high-lying unoccupied bands).Here, σ is the vector of the Pauli matrices, and V (r) is the crystal potential.
For the bulk Bi 2 Te 2 I 2 , in the basis of the two valence bands Ψ bulk v↑(↓) and two conduction bands Ψ bulk c↑(↓) , our ab initio four-band Hamiltonian reads: where , and the direct matrix product of the Pauli matrices τ and σ is implied (the explicit matrix form of H bulk kp is presented in Supplementary Note 3).Note that this Hamiltonian is the same (to within a unitary transformation) as that constructed for Bi 2 Se 3 in Ref. [31] within the theory of invariants.The matrices τ and σ in Eq. ( 1) have different meaning: τ operates in the valence-conduction band space, while σ refers to the total angular momentum J.
The parameters in Eq. ( 1) obtained within the LDA are: C 0 = 0.03 eV, C z = 0.13 a.u., C = 4.19 a.u., M 0 = −0.12eV, M z = 1.35 a.u., M = 5.88 a.u., V = 0.52 a.u., and V z = 0.13 a.u.(we use Rydberg atomic units: = 2m 0 = e 2 /2 = 1).Since the basis functions explicitly refer to the valence and conduction bands rather than to atomic orbitals, the parameter M 0 that defines the band gap at k=0 is negative and does not change sign upon moving from the topologically non-trivial insulator to the trivial one.The eigenvalues E(k) of the Hamiltonian (1) with the above parameters are shown in Fig. 1(c) by red lines, nicely reproducing the LDA curves over a quite large k-region and providing an absolute gap in the k•p spectrum.Moreover, these parameters reflect the band inversion and meet the condition of the existence of topological surface states (see, e.g., Ref. [43]) in accord with the Z 2 topological invariant ν 3D = 1 obtained from the parities of the bulk LDA wave functions at the TRIM points [44].Actually, the diagonal dispersion term M z( ) is positive, and it is larger than the electron-hole asymmetry: For the Bi 2 Te 2 I 2 thin films, we derive the Hamiltonian in the basis Ψ slab v↑ , Ψ slab c↓ , Ψ slab c↑ , Ψ slab v↓ as where , and τ refers now to the two decoupled sets of massive Dirac fermions.The Hamiltonian (2) is similar to the one obtained for 3D TI thin films within the effective continuous model based on the substitution k z → −i∂ z in the Hamiltonian of Ref. [31] and on the imposition of the open boundary conditions (see, e.g., Refs.[18] and [43]).The crucial difference is that in our ab initio approach within the same formalism for 3D and 2D systems we obtain the Hamiltonian and its parameters from the original spinor wave functions.We do not a priori impose the form of the Hamiltonian based on symmetry arguments and do not resort to the fitting of ab initio band dispersion curves or to a solution of 1D Schrödinger equations derived by using the above substitution with special boundary conditions.All the considered Bi 2 Te 2 I 2 films are characterized by the velocity V = 0.45 ± 0.01 a.u. and the electron-hole asymmetry C = 4.15 ± 0.10 a.u., which are weakly sensitive to the number of SLs, where the ± ranges indicate the variations of V and C in moving from 1 to 5 SLs.On the contrary, as seen in Figures 1(g) and 1(h) the parameters M 0 and M depend strongly on the film thickness, approaching monotonically zero.
In order to explicitly indicate whether a given film is a QSH insulator, in Fig. 1(g) we also plot the gap parameter ∆ = 2M 0 (−1) 1+ν2D with ν 2D being the Z 2 invariant obtained from the parities of the wave functions at the TRIM points of the 2D Brillouin zone.This parameter is negative for a topologically non-trivial film and positive for a trivial one.As follows from the figure, ∆ "oscillates" with the period of 2 SLs within the examined thickness interval.(The parity of Ψ slab v↑(↓) and Ψ slab c↑(↓) is (+) and (−), respectively, for ∆ < 0, and it is (−) and (+) for ∆ > 0.) As in 3D TI films [45], the thickness dependence of ∆ may be sensitive to the quasi-particle approximation employed, and it may change if many-body corrections beyond DFT are introduced.However, even the simplest quasi-particle method, the GW approximation for the self-energy, is methodologically challenging and computationally too demanding to study a large series of complex systems.Thus, DFT remains the method of choice, and its good performance for a wide range of TIs justifies the use of the Kohn-Sham band structure as a reasonable starting point.
The diagonalization of the Hamiltonian (2) then leads to E(k) shown by red lines in Figures 1(a) and 1(b).The absence of the absolute gap in the resulting k•p spectrum is the general feature of all the films studied.It is caused by the rather big electron-hole asymmetry C compared with the diagonal dispersion parameter M , Fig. 1(h).It should be noted that the conclusion on whether the edge states exist in a TI film is often made based on the signs and relative values of the parameters M 0 , M , and C .On the contrary, we find that the asymmetry |C | is larger than |M | everywhere, breaking one of the conditions for the film to be a QSH insulator, see, e.g., Refs.[46] and [18].Focusing on the behaviour of the diagonal dispersion M (as, e.g., in the topology analysis of Ref. [15]), we note that it is positive for all the thicknesses, Fig. 1(h).Along with the negative M 0 , this should signify an inverted band gap for the respective films.However, it does not correlate with the oscillating ∆, Fig. 1(g).
Let us now analyze the behaviour of the diagonal dispersion parameter M together with the topological invariant under a continuously varying geometry.We choose the 1SL film-the thinnest film, for which the k•p prediction of the band inversion does not contradict the actual topological property-and gradually expand the van-der-Waals spacing d vdW .The evolution of the band structure with increasing d vdW is shown in Fig. 3.According to the gap parameter ∆, see Fig. 4(a), a topological phase transition occurs at d vdW that is just around the mentioned 5% larger than its bulk value, and the 1SL film becomes topologically trivial.Further expansion leads to a larger band gap at Γ, which is not inverted anymore.It is noteworthy that such a behaviour of ∆ as a function of d vdW with the topological phase transition around 5% is stable with respect both to the choice of the approximation to the DFT exchange-correlation functional (LDA, GGA, dispersion corrected GGA) and to the SL geometry (bulk truncated or relaxed).In the limit of very large d vdW , when the BiTeI trilayers composing the 1SL film are too far from each other, the band structure is identical to that of a free-standing BiTeI trilayer (see Supplementary Fig. 2).Similarly, artificial reduction of the spin-orbit interaction strength λ relative to its actual value λ 0 in the equilibrium SL leads to a decrease in the gap, which closes at λ/λ 0 = 0.95.A further decrease in λ causes a widening of the already uninverted gap of the trivial phase.In general, the dependence of the relative gap-width on the spin-orbit interaction strength is almost linear and can be approximated as ∆(λ)/|∆(λ 0 )| = −20.9(λ/λ0 ) + 19.9.
The 1SL parameters of the 4-band k•p Hamiltonian (2) strongly depend on d vdW (the respective eigenvalues E(k) of this Hamiltonian are shown by red lines in Figures 3(a)-3(c)).With the d vdW expansion ξ (given in percents of the bulk value d (0) vdW ) up to 50%, the velocity V decreases monotonically from 0.470 a.u. to 0.342 a.u., and the electron-hole asymmetry C becomes smaller as well, Fig. 4(c).At ξ = 33%, C is already smaller than M , ensuring an absolute gap in the 4-band k•p spectrum, see Fig. 3(c).With further increasing ξ it even becomes negative, but it remains |C | < M .A stepwise behaviour of the parameter M that changes sign at the small ξ indicates that M keeps following the actual ν 2D and, thus, predicts a gap without inversion.With increasing d vdW this parameter again goes through zero around ξ = 30%, telling us that the band gap becomes inverted again, and at ∼ 35% with the given C and M 0 meets the conditions of the existence of the edge states [46,18].However, as seen in Fig. 4(a), the 1SL film is too far from a topological phase transition at such ξ.With this example we illustrate the strong limitations of the predictive capabilities of the effective continuous model.
Let us now analyze the formation of the SL band structure with the inverted band gap.Starting from well-separated layers, Fig. 3(f), and going back to the bulk value d As seen in Fig. 4(b), at Γ in the large-d vdW limitthere are two doubly degenerate energy levels, Upon decreasing d vdW , the TLs start to interact primarily by their Te-layer sides to cause the bonding-antibonding splitting of the two degenerate levels: The Te-related energies as a function of d vdW disperse stronger than those of Bi.Near the bulk value d (0) vdW , the splitting is large enough to invert the order of the E v 2 and E c 1 levels, ensuring the topological phase transition.Thus, the stacking procedure that leads to the 3D TI is based on SL building blocks principally different from the Rashba bilayers used in Ref. [21].It is essential that in our case the two Rashba constituents of the block (the stand-alone TLs) bring not only the Rashba-split conduction band but also the valence band, see Supplementary Note 3. Then the gap in the SL (which may be inverted or not) is quite naturally the gap between the valence and conduction bands, in contrast to the scenario of Ref. [21], where the band gap in the bilayer block is achieved by a dispersive "finite quantum tunneling" between the two Rashba constituents -the 2D electron gases of the adjacent layers.
Figure 4(d) shows the behaviour of the inverse effective masses of the chosen bands over the d vdW interval considered.We find that the conduction-band inverse masses, which are equal in the large- vdW the parameter M c 2 becomes twice as large, while M c 1 falls below zero.On the contrary, the valence-band inverse masses (M v 1 = M v 2 in the large-d vdW limit) "diverge" because the band v 1 moves down and "goes through" the I-orbital dominated bands, and v 2 moves up and hybridizes with Te p x,y bands, see Fig. 3 and Supplementary Fig. 2. Finally, at d vdW the parameter M v 2 reaches its large d vdW limit, while M v 1 becomes negative.Thus, in the topologically non-trivial 1 SL we have M v 2 ≈ C − M and M c 1 ≈ C + M , where C and M are the 1SL parameters of the Hamiltonian (2).At that, the interband coupling of the bands v 2 and c 1 is equal to V of the 4-band k•p description.This reveals a close relation between the 4-band and 8-band Hamiltonians.However, already with 8 bands there is an absolute gap (see Fig. 3(a)), which is reasonably accurate and quite suitable for the theoretical research on linear response, Hall conductance, and motion of Dirac fermions in external fields.
of the van-der-Waals spacing, Fig.3(a), we retrace the valence bands (v 1 with the energy E v 1 and v 2 with E v 2 ) with the predominant contribution coming from the p z orbitals of Te and the conduction bands (c 1 with E c 1 and c 2 with E c 2 ) mainly formed by Bi p z orbitals, see Supplementary Fig.2.We derive an 8-band Hamiltonian H 1SL kp which is presented in Supplementary Note 3. Its eigenvalues are shown in Fig.3by blue lines, and the corresponding parameters as a function of the d vdW expansion are depicted in Figures 4(b) and 4(d), see also Supplementary Fig. 4.

Figure 1 :Figure 2 :
Figure 1: Ab initio calculations for Bi 2 Te 2 I 2 .LDA band structure (black curves) of 1SL (a) and 5SL (b) films as well as the bulk crystal (c) of Bi 2 Te 2 I 2 .Red curves correspond to the eigenvalues E(k) of the model four-band Hamiltonian.Atomic layers maximally contributing to the valence and the conduction band at the center of the 2D Brillouin zone in the case of the Bi 2 Te 2 I 2 films of different thickness are indicated in graph (e) by green and red bars, respectively.The dependence of the model Hamiltonian parameters on the number of SLs is presented in graphs (g) and (h).The atomic structure of the bulk Bi 2 Te 2 I 2 (unit cell) and its 1SL film (d) with the corresponding 3D and 2D Brillouin zone (f) are also shown.

Figure 3 :
Figure 3: Band structure of the 1SL Bi 2 Te 2 I 2 film with different van-der-Waals spacing.The expansion of the spacing is 10% (b), 40% (c), 100% (d), 200% (e), 1100% (f) as indicated in percents of the bulk value, (a).The spectra provided by the 4-band (red lines) and 8-band (blue lines) k•p models are presented in comparison with the LDA bands (black lines).Light green lines show the results of the 36-band Hamiltonian to demonstrate the convergence respect the size of the basis set.

Figure 4 :
Figure 4: The parameters of the k•p Hamiltonian for the 1SL Bi 2 Te 2 I 2 film.The parameters of the 4-band (a, c) and 8-band (b, d) Hamiltonian are shown as a function of the expansion of the van-der-Waals spacing, which is given in percents of the bulk value.In graph (b), the parity of the respective LDA wave functions is also shown.