Regularities of the quantum spin Hall phase formation in three-dimensional tetradymite-like topological insulator thin ﬁlms

We present ﬁrst-principles calculations of electronic structure and topological properties of thin ﬁlms of tetradymite-like binary, ternary, and quaternary compounds. Two-dimensional topological insulators were found, including two septuple layer ﬁlms of PbBi 2 Te 4 with large band gap. We demonstrate that the Pb-based thin ﬁlms inherit Z 2 dependence on the ﬁlm thickness from the “parent” quintuple-layer compounds. An effective way to induce the topological phase transition in the trivial thin ﬁlms is proposed that does not require substantial changes in the in–plane lattice parameter and surface crystal structure.


I. INTRODUCTION
A theoretical prediction and further experimental verification of quantum spin Hall (QSH) effect induced by spinorbit interaction in two-dimensional (2D) insulating systems lead to the identification of a state of quantum mattertwo-dimensional topological insulators (2D TIs) [1][2][3][4][5][6][7]. Such systems are characterized by topological invariant Z 2 = 1 and helical gapless edge states protected by time-reversal symmetry from backscattering. As a result, charge carriers at the edges provide dissipationless spin-polarized current [8], making 2D TIs promising for applications in low-powerconsumption electronic devices.
To date, QSH effect has been theoretically predicted in a number of compounds containing heavy elements like Te, Bi, Pb, Mo, Ti, and Sn . In particular, such compounds are by nature of different types of monolayers [9,16], bilayers [12,14], 2D inversion-asymmetric TIs [20], chemically functionalized systems [21,22], heterostructures [33], and others. Despite this successful progress in the theoretical prediction of the QSH systems, most proposed 2D materials have not been fabricated and experimentally verified.
One of the central topics in this field is the search of 2D TI materials simultaneously combining both sizable band gap and easy fabrication. Promising candidates for realization of the QSH state in 2D systems can be thin films of 3D TIs with layered tetradymite-like crystal structure [26][27][28][29][30][31][34][35][36][37]. Such films can be easily produced by using the advanced technology of nano-ribbon growth method or molecular beam epitaxy. Note that the high-quality samples of thin films of 3D TIs, for instance Bi 2 Se 3 [38,39] and Bi 2 Te 3 [40,41], have already been fabricated. Additionally, tetradymite-like layered compounds can form homologous series of ternary and quaternary compounds, including a family of natural minerals [42][43][44]. As a result, the common factors in the QSH phase formation in such thin films can be found, which have not been investigated yet. Also in such systems, the oscillations of the * menshikova_t@mail.ru Z 2 topological invariant with an increase of the film thickness can be expected [29][30][31].
In this paper, we present a detailed investigation of the topological properties of binary, ternary, and quaternarycompound thin films. For the binary-compound thin films, we show the dependence of the topological properties on crystal structure parameters and calculation details. We demonstrate the dependence of the thin film topological properties on the chemical composition, and predict the existence of the QSH phase among Bi 2 Te 2 S, Bi 2 Te 2 Se, and Bi 2 Te 2.4 Se 0.6 thin films. In the latter case, we find that the type of disordering strongly influences the topological properties. We predict existence of 2D TIs among thin films of Pb-based compounds. Moreover, Pb-based compound thin films are characterized by the same dependence of the Z 2 topological invariant on the thickness that is found for "parent" quintuple-layer (QL) compounds. We also demonstrate the effective way of tuning topological properties, which provides an opportunity to induce the QSH state in ultrathin films.
The paper is organized as follows. Section II provides computational details. In Sec. III, we present a study of topological properties depending on the film thickness in binary and ternary compounds. The conclusions are given in Sec. IV.

II. CALCULATION METHODS
Electronic structure calculations were carried out using the projector augmented-wave (PAW) [45,46] and full-potential linearized augmented plane wave (FLAPW) methods implemented in the ABINIT [47] and FLEUR [48] code, respectively. The PAW data sets in the ABINIT code were taken from Ref. [49]. The exchange-correlation energy was treated using the generalized gradient approximation (GGA) with PBE functional [50], and local density approximation (LDA) with PW functional [51] for ABINIT, and PZ functional for FLEUR [52]. The Hamiltonian contained scalar relativistic corrections and the spin-orbit coupling was taken into account by the second variation method [53]. The Z 2 invariant was obtained Bi Te 2 3 Bi Te S from the parity of occupied electronic states at the timereversal invariant points of the 2D Brillouin zone [54] within both methods. The electronic structures of thin films were calculated by using various structural parameter sets. The first one includes several sets of experimental crystallographic data. The second set is a structure geometry obtained from optimization of interatomic distances at fixed experimental lattice parameters a and c in the bulk (denoted as "EXPn-opt," n is a variant of experimental data). In the case of Bi 2 Se 3 , full optimization of the bulk parameters was also done, including cell parameters and interatomic distances (see Supplemental Material (SM) [55] for more details). The last configuration of thin films includes the surface relaxation: the atomic layers in the outermost structural blocks (from vacuum side) were relaxed with inclusion of van der Waals (vdW) correction [56].
To treat the disordered phase, we employed a virtual crystal approximation (VCA) as implemented in the ABINIT code where the averaged ionic pseudopotential of a virtual atom occupying a site in the A/B sublattice is defined as a mixture V VCA = xV A + (1 − x)V B of A (V A ) and B (V B ) pseudopotentials. In our paper, the role of A/B sublattice plays Bi/Sb or Te/Se sublattice. In the ABINIT calculations for this case, we used GGA-PBE Hartwigsen-Goedecker-Hutter relativistic norm-conserving pseudopotentials, which include the spin-orbit interaction [57].

A. Topological properties of quintuple-layer thin films
The binary compounds we study have a tetradymite-like structure with alternation of weakly bounded vdW forces, QL blocks along the c axis [ Fig. 1(a)]. Despite many electronic structure calculations of binary compound thin films performed to date, the topological properties of these thin films are still under debate. For instance, in Ref. [30] the 2D QSH state in some thin films of Bi 2 Se 3 was shown to exist by both analytical models and first-principles calculations. However, recent GW (GW approximation is expansion of the self-energy in terms of the single particle Green's function G and the screened Coulomb interaction W) [26] and full potential ELAPW (extended linearized augmented plane wave method) [31] calculations demonstrated that thin films of Bi 2 Se 3 are topologically trivial. For this reason, we start with a detailed investigation of the topological properties of thin films Bi 2 Se 3 , Bi 2 Te 3 , and Sb 2 Te 3 in dependence on different structure parameters and exchange-correlation energy within PAW and FLAPW methods implemented in the ABINIT and the FLEUR code, respectively (see S1 in the SM [55]).
Summarizing the obtained results, we conclude that the GGA approximation, in general, gives more consistent results compared to those obtained with the LDA. The optimization of the crystal structure can dramatically affect the topological properties. Also, strong sensitivity of these properties to the calculation method was found in the case of Bi 2 Te 3 (S1 in the SM [55]). The analysis of the obtained results showed that the Bi 2 Se 3 thin films are trivial. In the case of Bi 2 Te 3 , the thin film topological invariant changes a minimum of twice (Table I), and the Sb 2 Te 3 thin films are 2D TIs with a thickness of 3 QLs and 4 QLs. The topologically nontrivial properties of thin films are reflected in the inversion of gap edges. A careful analysis shows that the atomic layers near the vdW gaps [the closest-Te(Se), and the second-Bi(Sb)] are an essential factor responsible for the inversion of the gap edges. It assumes that changing or replacing these layers with other atomic layers can strongly influence the topological properties of thin films. In this context, a relevant modification is replacing of the fourth or fifth layer (near the vdW gap) in the outermost QL blocks by their isoelectronic analog. Such a modification does not cause substantial changes in the in-plane lattice parameter and the surface crystal structure. Note that the mentioned modification will be well pronounced in the case of ultrathin films formed by two or three QLs. Also, these films are more interesting for practical applications due to larger band gaps in comparison with those in the thicker films.
To illustrate this idea, we consider the cases of topologically trivial thin films of binary compounds: 2 QLs of Bi 2 Se 3 , 3 QLs of Bi 2 Te 3 , and 2 QLs of Sb 2 Te 3 . First, let us consider simple replacement of the fourth or fifth layer from the vacuum side (near the vdW gap) in the outermost QL blocks by isoelectronic atomic layer without any relaxation. In the case of 2 QLs of Bi 2 Se 3 , the replacement of Se layers forming vdW by Te leads to the topological phase transition. The electronic structure of such a system is characterized by the inverted energy gap and Z 2 = 1. The same effect was found for 3 QLs of Bi 2 Te 3 upon substitution of Bi layers closest to the vdW spacing by Sb in the outermost QL blocks. The modification of 2 QLs of Sb 2 Te 3 by replacement of the fourth Sb layer by Bi in the QL blocks leads to a topologically nontrivial system. Note, similar results were found upon substitution of both Bi or Sb layers in the outermost QLs of 3 QLs of Bi 2 Te 3 and 2 QLs of Sb 2 Te 3 , respectively.
Second, as was mentioned above, the optimization of the crystal structure parameters can significantly affect the topological properties. In the case of 2 QLs of Bi 2 Se 3 with substituted Se layers, the optimization induces a sizable increase of vdW spacing that is responsible for formation of the topologically trivial phase [32]. In contrast, for modified 3 QLs of Bi 2 Te 3 and 2 QLs of Sb 2 Te 3 , the optimization does not affect the topological properties, Z 2 = 1, as it happens in the unrelaxed cases.
Additionally, in thin films, the formation of alloy systems like Bi 2−x Sb x Te 3 and Sb 2−x Bi x Te 3 is also possible. To examine the evolution of topological characteristics in such films, let us consider the 3 QLs of Bi 2−x Sb x Te 3 and 2 QLs of Sb 2−x Bi x Te 3 with systematically varying the concentration x. In this context, the topological transition occurs at a certain concentration of x that is characterized by closing the band gap [58]. As seen in Fig. 2, the topological phase transition takes place in Bi 1.48 Sb 0.52 Te 3 and Sb 1.55 Bi 0.45 Te 3 for 3 QLs and 2 QLs, respectively. In general, the dependence of the band gap on concentration is linear for both cases.
Despite the central atomic layer in QL blocks, it makes only a small contribution to the inversion of gap edges, and replacing such a layer in all blocks by another atomic layer can also affect the topological properties due to the changing of the chemical interaction in the system. This situation is realized in the case of ternary tetradymite-like compounds like Sb 2 Te 2 S, Sb 2 Te 2 Se, and Bi 2 Te 2 S, Bi 2 Te 2 Se [ Fig. 1(b)]. The presence of S and Se in the center of QLs leads to the increase of chemical bonding within the QLs in comparison with parent binary compounds. The reinforcement of bonds is directly related to higher electronegativity of S and Se in comparison with Te. Moreover, such an increase of chemical bonding is followed by the decrease of the interaction between QL blocks, leading to significant separation of adjacent covalent blocks [59]. In the case of Sb 2 Te 2 S and Sb 2 Te 2 Se, the vdW gap expands by ∼9% and 7%, respectively, in comparison with Sb 2 Te 3 . As a consequence, all thin films of these compounds are trivial ( Table I). The expansion of vdW gap also takes place in ternary compounds Bi 2 Te 2 S and Bi 2 Te 2 Se. However, the presence of Bi layers characterized by stronger spin-orbit coupling than Sb leads to the formation of the QSH phase in thin films of Bi 2 Te 2 S and Bi 2 Te 2 Se (Fig. 3). As can be seen from Table I, 3-QL and 4-QL-thin films of Bi 2 Te 2 S are topologically nontrivial. In the case of Bi 2 Te 2 Se thin films, the dependence of Z 2 on the film thickness is different. The QSH phase was found for smaller thickness of films 2 QLs and 3 QLs. Note, for both compounds, Bi 2 Te 2 S and Bi 2 Te 2 Se, the electronic spectra of the 2 QL-4 QL films are semimetallic.
At the same time, intriguing topological properties can be expected in the case of thin films of recently studied 3D TI Bi 2 Te 2.4 Se 0.6 [60,61]. A distinct feature of this compound is the existence of two possible types of disordering [59,62] [ Fig. 1(c)]. The first type of disorder is a substitution of Te atoms in the middle layer of QLs by Se atoms (denoted as type I) [62]. The second one is random distribution of Se atoms in both Te sublattices formed by the central layer and outer Te atoms of QL blocks (denoted as type II) [59]. Since both types of disorder are possible, we consider the topological properties of the Bi 2 Te 2.4 Se 0.6 thin films for type I and type II cases.
Analysis of topological properties of Bi 2 Te 2.4 Se 0.6 thin films with admixture of Se atoms in the middle layer of QLs reveals the same dependence of Z 2 on the film thickness as in the case of Bi 2 Te 3 (Table I). In fact, such a concentration of Se atoms is not enough for tuning thin film topological properties to the Bi 2 Te 2 Se case. At the same time, unexpected results were obtained for the random distribution of Se atoms over the Te sublattices (Table I). As can be seen from the Table I, the topologically nontrivial films are 2 QLs and 3 QLs, as in the case of Bi 2 Te 2 Se. This can be understood from analysis of the crystal structure. The analysis shows the QL blocks of Bi 2 Te 2.4 Se 0.6 (type II) and Bi 2 Te 2 Se are chemically close: the value and ratio of the interplanar distances are almost the same for these compounds.

B. Pb-based septuple-layer thin films
The Pb-based ternary compounds (PbBi 2 Se 4 , PbBi 2 Te 4 , PbSb 2 Te 4 ) belong to the homologous series of systems nA I V B V I -mA V 2 B V I 3 (A I V = Ge, Sn, Pb; A V = Bi, Sb; B V I = Te, Se) with n,m = 1. Their crystal structure is hexagonal (R3m) and the unit cell contains septuple-layer (SL) blocks aligned along the hexagonal c axis [63][64][65] [Fig. 1(d)]. Quaternary Pb-based compounds PbBi 2 Te 2 S 2 (aleksite) and PbBi 2 Te 2 Se 2 (poubaite) are members of a natural minerals family and have the same crystal structure as the respective ternary compounds [42,43] [ Fig. 1(e)]. As in the case of QL compounds for the Pb-based thin films, we investigate the influence of different structure parameters and calculation methods on their topological properties (see S2 in the SM [55]).
Analysis of the thin film electronic structure shows that all the studied films of PbBi 2 Se 4 have a band gap, which  (Fig. 4). PbBi 2 Te 2 S 2 thin films are topologically nontrivial at 3 SLs and 4 SLs, and the PbBi 2 Te 2 Se 2 ones-at 2 SLs and 3 SLs (Table II and Fig. 4). Summarizing the above-mentioned results, one can note that in the case of Pb-based compounds, the behavior of Z 2 topological invariant is the same as for the parent QL compounds. In particular, thin films of PbBi 2 Se 4 are topologically trivial, as for Bi 2 Se 3 . In the case of PbBi 2 Te 4 , the Z 2 topological invariant changes a minimum of twice, as in the case of Bi 2 Te 3 . Sb 2 Te 3 and PbSb 2 Te 4 thin films have a similar Z 2 topological invariant, which changes only once. The topological invariant of the PbBi 2 Te 2 S 2 and PbBi 2 Te 2 Se 2 thin films is identical to the one obtained for Bi 2 Te 2 S and Bi 2 Te 2 Se thin films, respectively. Such inheritance of the topological properties in Pb-based compounds is directly connected to their crystal structure and peculiarities of the gap edges inversion. The crystal structure of Pb-based compounds can be represented as QLs of binary (Bi 2 Se 3 , Bi 2 Te 3 , Sb 2 Te 3 ) or ternary (Bi 2 Te 2 Se, Bi 2 Te 2 S) compounds with PbTe, PbSe, or PbS bilayer inserted into the center of the QL blocks. At the same time, the analysis of the thin film electronic structure shows that the main role in the inversion of the gap edges plays the states formed by atomic layers near the vdW gaps. The central layer states are not involved in the inversion as in the case of QL compounds. In such a way, incorporation of the PbTe (PbSe or PbS) bilayer into the QLs center does not change much the orbital composition of the gap edges in comparison with the case of the parent QL compounds. Moreover, the increase of the number of such bilayer does not change the situation either which is corroborated by the results of the Ref. [29] devoted to the study of topological properties of nine-layer Pb-containing ternary compound thin films.

IV. SUMMARY AND CONCLUSION
In summary, we have studied the electronic structure and topological properties of thin films of binary, ternary, and quaternary compounds. In the example of binary compounds, we have demonstrated the GGA approximation gives more consistent results for electronic structure in comparison with the LDA. Within both ABINIT and FLEUR calculations, we have found the dependence of topological properties on the crystal structure parameters. The full correspondence between these methods was achieved at optimized parameters, including surface relaxation. For the system under debate, Bi 2 Se 3 , all films (1 QL-4 QLs) were found to be topologically trivial within both ABINIT and FLEUR calculations.
We have shown that in ternary QL thin films, the topological properties substantially change in comparison with isostructural binary systems. In particular, all Sb 2 Te 2 S and Sb 2 Te 2 Se thin films are topologically trivial. Such a situation is directly related to huge vdW gaps caused by weakening of the interaction between QL blocks. This is caused by the increase of chemical bonding within the QLs due to the presence of the S and Se atoms, which possess higher electronegativity and smaller size in comparison with Te.
We have predicted 2D TIs among Pb-based thin films characterized by the same dependence of Z 2 on the film thickness as the parent QL compounds. Such inheritance of the topological properties is directly related to the crystal structure and the fact that the atomic layers located close to the vdW gap are responsible for topological properties of these films. None of the studied ultrathin PbBi 2 Se 4 (0001) films are 2D topological insulators. We predict the QSH phase for PbSb 2 Te 4 and PbBi 2 Te 4 films. In the latter case, the films formed by the odd number of SL blocks (one and three) are not 2D topological insulators, while the films formed by the even number of blocks (two and four) are topological insulators. The PbBi 2 Te 4 (0001) film with 2 SL thickness is most promising for the experimental observation of topological edge electronic states because it is characterized by the largest gap. In the case of natural mineral thin films, the dependence of the topological properties is identical to the parent ternary compounds.
Finally, we have proposed an effective way to induce the topological phase transition in the trivial thin films. This way, the replacement of the layers closest to the vdW gap in the outermost structural blocks by isoelectronic atomic analog does not require substantial changes in the in-plane lattice parameter and surface crystal structure. This method will be very efficient in the case of ultrathin films formed by two or three structural blocks.
It is necessary to notice that results may depend on the exchange-correlation functional used. Therefore, all results and conclusions in this paper are valid for the PBE exchangecorrelation functional we used in the calculations. We hope that our findings help to shed light on the topological properties of 2D systems and will be useful for experimental investigations.