Strain driven topological phase transitions in atomically thin films of group IV and V elements in the honeycomb structures

We have investigated topological electronic properties of freestanding bilayers of group IV (C, Si, Ge, Sn, and, Pb) and V (As, Sb, and, Bi) elements of the periodic table in the buckled and planar honeycomb structures under isotropic strain using first-principles calculations. Our focus is on mapping strain driven phase diagrams and identifying topological phase transitions therein as a pathway for guiding search for suitable substrates to grow two-dimensional (2D) topological insulators (TIs) films. Bilayers of group IV elements, excepting Pb, generally transform from trivial metal → ?> topological metal → ?> TI → ?> topological metal → ?> trivial metal phase with increasing strain from negative (compressive) to positive (tensile) values. Similarly, among the group V elements, As and Sb bilayers transform from trivial metal → ?> trivial insulator → ?> TI phase, while Bi transforms from a topological metal to TI phase. The band gap of 0.5 eV in the TI phase of Bi is the largest we found among all bilayers studied, with the band gap increasing further under tensile strain. Differences in the topological characteristics of bilayers of group V elements reflect associated differences in the strength of the spin–orbit coupling (SOC). We show, in particular, that the topological band structure of Sb bilayer becomes similar to that of a Bi bilayer when the strength of the SOC in Sb is artificially enhanced by a factor of 4. This study provides the first report that As can be a 2D TI under tensile strain. Notably, we found the existence of TI phases in all elemental bilayers we studied, except Pb.


Introduction
Two-dimensional (2D) topological insulators (TIs), also known as quantum spin Hall (QSH) insulators, are novel nonmagnetic insulators, which are currently attracting intense worldwide interest [1][2][3][4][5]. A unique characteristic of a QSH insulator is the presence of gapless edge state, which carry two counter-propagating spin-polarized currents. These gapless edge states hamper backscattering due to constraints of time-reversal symmetry, making them very desirable for spintronics and other applications. Moreover, 2D TIs are more robust against nonmagnetic impurities than three-dimensional (3D) TIs since the only available scattering channel, the backscattering channel in 2D, is forbidden. Although graphene was proposed to support a QSH state, the spin-orbit coupling (SOC) in graphene is very weak, and the associated gap is too small to be accessible experimentally [6].
Recent theoretical work has shown that group IV [7] and V [8] elements of the periodic table as well as III-V [9] compounds can also form atomically thin crystals with honeycomb structure. Freestanding silicene and its group IV and V elemental counterparts are 'beyond' graphene 2D materials with advantages over graphene of a stronger SOC and atomic bonds that are buckled and not flat. Very recent, first-principles calculations show many of these materials to be 2D TIs at the equilibrium structure, including silicene, germanene [10], stanene [11], bismuth [12], InBi, GaBi, and TlBi [13]. Furthermore, Sb [14], BBi, and AlBi [13] have been shown to become 2D TIs under tensile strain. Robustness of the topological phase in Sb [14] and Bi [15,16] against strain and electric field have also been studied. Furthermore, recent studies indicate that band topology in group IV and V honeycombs could be altered via chemical adsorptions [17,18]. Due to the buckled atomic bonds, the band gap in beyond graphene materials can be tuned with a perpendicular electric field, which breaks the inversion symmetry, and topological phase transitions can be realized by gating control [19]. For group IV 2D honeycomb thin films, the resulting field-tunable nondegenerate states possess nearly 100% spin-polarization, providing a basis for designing high efficiency spin filters, and other devices for manipulating spin currents [19][20][21].
In this study, we use first-principles calculations to investigate the atomic and electronic structures of bilayers of group IV (C, Si, Ge, Sb, and Pb) and group V (As, Sb, and Bi) elements in buckled and planar (PL) honeycomb structures. Our focus is on identifying ranges of strains where topological phases exist in these atomically thin films, providing a tangible guide for finding suitable substrates for practically realizing these 2D materials. We discuss similarities and differences between the strain driven phase diagrams of bilayers of group IV and group V elements. Although aspects of thin films discussed in this study have been reported in the literature, we are not aware of an earlier comprehensive study of the group IV and V bilayers as we have attempted in this article.

Computational methods
Computations were carried out within the generalized gradient approximation (GGA) [22] to the density-functional theory [23] using the projector-augmented-wave method [24] as implemented in the Vienna ab initio simulation package [25]. The cutoff kinetic energy was set at 500 eV. SOC was included in the band structure calculations. To model the single bilayer, a vacuum layer of 20 Å was included in the supercell. A 21 × 21 × 1 Monkhorst-Pack grid was used to sample the 2D Brillouin zone (BZ). We constructed both buckled and PL honeycomb structures for all bilayers investigated, and optimized the atomic positions within symmetry constraints at every lattice constant considered. Atomic positions were relaxed until the residual force on each atom was smaller than 0.0001 eV Topological characteristics of the band structures of films were assessed by computing the Z 2 topological invariant through an analysis of the parity of wave functions at the time-reversal invariant momentum points in the 2D-BZ. This procedure is applicable to crystals with inversion symmetry [3], as is the case here. The Z 2 invariant, ν, is then given in terms of products of parity eigenvalues of Kramerʼs doublets of occupied valence states at the four timereversal invariant k-points (K i ) as follows ( ) ( ) where ξ = ± is the parity eigenvalue and the number of occupied bands is 2N. The Z 2 invariant in all cases discussed in this paper was determined by this method. The value of Z 2 was computed at each value of strain and the evolution of the associated band gap was monitored as a function of strain to find possible topological phases transitions.
We also carried out a number of computations in which strength of the SOC was varied to ascertain whether or not the inversion of the band gap in a particular bilayer was driven by the SOC. This was done by reducing strength of the SOC in the band structure calculation by weighting the SOC contribution via a scaling factor WSOC, which varied between 0 and 1.

Results and discussion
3D bulk C, Si, Ge, and Sn assume a diamond structure, bulk Pb is face-centered cubic, and bulk phases of group V elements are in the rhombohedral lattice. With the exception of Pb, both the diamond and the rhombohedral lattice can be viewed as a stacking of 2D buckled honeycombs with two atoms/unit cell along a (111) direction. Therefore, we start our bilayer computations using a 2D honeycomb structure, with and without buckling, and increase/decrease lattice constant in small steps, optimizing structure in all cases.
Figures 1(a)-(h) summarize the results of our computed total energies and vertical distances between the two planes in the bilayers of various group IV and V elements as a function of lattice constant, or equivalently strain, for both PL and buckled structures. The lattice constant for zero strain is defined to correspond to the value obtained by projecting the 3D bulk phase on to a 2D hexagonal lattice. Three optimized equilibrium states can be identified from these results. A low-buckled (LB) honeycomb structure is formed at larger lattice constants, whereas the high-buckled (HB) honeycomb has lower energy at smaller lattice constants. As we increase the lattice constant, differences between the total energies of PL and buckled honeycomb become smaller and the PL honeycomb (d = 0) is seen to be stabilized in a number of cases (e.g., As in figure 1(f)). Our results are in accord with earlier studies, which have shown that there are three local minima (LB, HB, and PL) in the total energy landscape of bilayers of group IV and V elements [7]. Although the energy of the HB state for Ge, Sn, and Pb bilayers is lower, the structure becomes unstable due to the presence of negative frequency phonon modes [7]. Notably, the LB state is seen to be energetically favorable compared to the PL state in all cases with the exception of graphene in which the PL state has a lower energy over a wide range of lattice constants as a result of its sp 2 bonding. Lattice constants and other Since HB films are not stable, we will limit further discussion to LB and PL films [7,17,26]. In this connection, we examined the phonon spectra of the LB and PL films using the phonopy package [27] with a 7 × 7 supercell. No negative frequency was found for any of the LB films considered expecting carbon, while the PL films were found to be stable only for carbon. The PL thin films for other elements would require suitable substrates to stabilize the structure. Topological characteristics of the associated band structures were identified via a parity analysis.
In order to determine whether the band structure is metallic or insulating, it is useful to consider a 'system band gap', which is defined as the energy difference between the lowest level of the conduction bands and the highest level of the valence bands. A positive system gap implies an insulating state, while a negative value means that the system is metallic. Values of the system band gap are shown as a function of lattice constant in figure 2 for various bilayers. We monitored topological phase transitions by examining changes in the band structure as the  strain is varied gradually. At a topological transition point, the conduction and valence states meet each other at certain k-points. The system is a trivial with Z 2 = 0 on one side of the transition and topologically non-trivial with Z 2 = 1 on the other side of the transition point. Note that Z 2 is well-defined only when a gap between the conduction and valence bands exists throughout the BZ. A metal can have non-trivial band topology when the conduction and valence bands are separated at all the k-points, which is obviously not the case in a metal that the conduction and valence bands are connected.
In figure 2, the non-trivial and trivial insulator phases are identified with filled and open symbols, respectively. Among the group IV bilayers, the range of lattice constants over which the TI phase exists (filled symbols) is seen to decrease in figure 2(a) with increasing atomic number as we go from Si to Ge to Sn (Pb does not yield a TI phase). Within group IV, the largest non-trivial system band gap of about 72 meV is exhibited by a Sn bilayer. Among the group V bilayers in figure 2(b), on the other hand, Bi presents the widest range of lattice constants over which it remains topologically non-trivial, this range being the smallest in As, which is opposite to the trend in group IV bilayers with increasing atomic number. The nontrivial band gap supported by the Bi bilayer at optimal lattice constant is also the largest among the group V elements with a value of about 501 meV. Although the non-trivial band gap in Bi does not increase monotonically with lattice constant, it reaches a maximum value of 773 meV under 12% tensile strain. The aforementioned band gaps are well above the energy scale of thermal excitations at room temperature. Figure 3 shows total energies of the LB (black points and lines) and PL (red points and lines) bilayers as a function of the lattice constant or strain as in figure 1, with regions of various topological and non-topological phases separated by vertical lines. In films of group IV elements (C, Si, Ge, and Sn) in panels (a)-(d), phase transformations in going from compressive to tensile strains in the LB state generally follow the pattern (black dashed vertical lines): trivial metal (M) → topological metal (TM) → TI → TM → trivial metal (M). The sequence of phases for PL films shown by dot-dashed vertical lines is similar, although it is shifted and rescaled as a function of strain. These differences between the behavior of LB and PL films are driven by the size of the buckling of bilayers in the LB state [14,16]. In Ge and Sn, as the film is compressed, the HB state intervenes, and the trivial metal (M) phase gets skipped. As already pointed out in connection with figure 2 above, in films of group IV elements, the range of the TI phase in the LB state decreases progressively with increasing atomic numbers. Pb is an outlier, which is different from other group IV elements, displaying only the trivial metal (M) and insulator (I) phases in figure 3(e). Figure 4 gives further insight into the nature of topological phase transitions discussed above by presenting schematics of how the key band gaps at Γ and K and the parities of the associated energy levels evolve in various phases as strain is varied. Focusing first on the group IV bilayers in figures 4(a)-(e), we see that the size of the gap at the K-point is fairly insensitive to strain, while the band gap at Γ undergoes large variations. Figures 4(a)-(c) illustrate the transition from TM to TI to the TM phase without band closing, implying that the topology of bands remains unchanged through this series of transformations. In contrast, the closing and reopening of the band gap with strain through a critical point in figure 4(d) is the result of the system switching between trivial and non-trivial states.
Turning to group V bilayers, As and Sb bilayers undergo trivial metal (M) to trivial insulator transition with strain in going from the left hand to the right hand side in figures 3(f) and (g), followed by the realization of the TI phase before the system becomes PL. Notably, this is the first report that an As bilayer could become a 2D-TI under tensile strain. The schematic band structure in As and Sb films at Γ is shown in figures 4(f)-(h) along with parities of the key energy levels undergoing changes with strain. At zero strain, the parities at Γ are in the order ++−, figure 4(f), going from low to high energies. At strains larger than 6%, the band gap becomes inverted with parities of the two levels around the band gap interchanged, see figure 4(h). This change in topology yields the TI phase in As and Sb for tensile strains beyond 6%.
Bi bilayers are seen from figure 3(h) to follow a different path than As and sb 6 in that Bi films transform from being a TM to a TI phase and then again a TI [16] with increasing lattice HB, high-buckled. 6 We note that thin films of group V element P exhibit behavior similar to that of As and Sb films, but the SOC strength in P is much weaker so that the band gap in the non-trivial phase is smaller.
constant. The associated schematic band structures at Γ are compared in figures 4(i)-(k). The parities at Γ are in the order −++ (bottom to top) in each of the figures 4(i)-(k) at compressive or zero strain, and the bilayer assumes a TM or TI phase. With further increase in the strain, the band gap reduces to zero and then reopens but without swapping parities of eigenvalues, so that the system remains in the TI phase. As to PL bilayers, As, Sb, and Bi all undergo a transition from M → I phase with decreasing compressive strain, see figures 3(f)-(h). We have carried out computations based on the hybrid functional HSE06 [28] to examine the robustness of our results. Specifically, we performed calculations at the equilibrium lattice constant in the LB and the PL phase. These calculations show that our GGA based results correctly capture the evolution of topological phases, some differences in the exact values of the strain at which phase transitions occur notwithstanding.
The preceding discussion shows that the topological behavior of Bi bilayers is different from that of As and Sb films even though all three are group V elements. We expect, of course, these differences to be driven by those in the strength of the SOC in these elements. In this connection, we considered band structure of Sb bilayer with artificially enhanced SOC strength by scaling it with a weighting factor WSOC. Band structure of an Sb bilayer at the lattice constant of 4.5 Å at which the corresponding Bi bilayer is a TI phase is shown in figure 4(l) for WSOC = 1. The band structure of the Sb bilayer in figure 4(l) is similar to that in figure 4(f) and represents a trivial metal. However, when the SOC strength is increased by a factor of 4 by setting WSOC = 4 in the computations, the resulting band structure of a hypothetical Sb bilayer in figure 4(m) becomes similar to that of a Bi bilayer as seen by comparing the order of band parities at Γ in figure 4(m) with that in figure 4(k).

Conclusions
To summarize, the topological phases of freestanding bilayers of group-IV (C, Si, Ge, Sn, and, Pb) and V (As, Sb, and, Bi) elements in PL and buckled honeycomb structure are discussed under wide ranges of isotropic compressive and tensile strains using first-principles calculations. Our results for the group IV bilayers show that, as strain varies from negative (compressive) to positive (tensile) values, these bilayers generally transform from a trivial metal → TM → TI → TM → trivial metal phase with the exception of Pb. In contrast, among the group V bilayers, As and Sb transform from a trivial metal → insulator → a TI phase with strain, but Bi transforms from a TM to a TI. We delineate how topological characteristics evolve in various bilayer systems in terms of the evolution of underlying band orderings with strain at the high symmetry points, and how differences between bilayers of various elements are driven by differences in the strength of SOC of the elements. The bilayers of group IV elements remain insulating over fairly narrow ranges of strains due to the presence of relatively small band gaps at Γ. In group V bilayers, the strength of SOC can drive band inversions at Γ as is the case in Bi compared to As and Sb bilayers. The band gap of 0.5 eV in the TI phase of Bi is the largest we found among all bilayers studied, with the band gap increasing further under tensile strain. Furthermore, our study is the first, to our knowledge, to predict the possibility of realizing a TI phase in As bilayers under tensile strain. The comprehensive phase diagrams of bilayers of group IV and V elements presented in this study would provide a tangible basis for guiding search of viable substrates for growing and realizing topological phase transitions in 2D-TI systems.