Pressure effects on crystal and electronic structure of bismuth tellurohalides

The bulk electronic structure calculations were carried out within DFT using the projector augmented-wave  method realized in the VASP code and the PBE exchange-correlation functional [29]. The SOC was accounted by a second-variation method. Crystal lattice parameters and atomic positions were optimized for pressures up to 8, 10, 18 GPa for BiTeI, BiTeBr and BiTeCl, respectively, using PBE exchange-correlation functional. Under the structural relaxation the crystal symmetry kept the same. The optimization was performed for all structures of BiTeX (X = I, Cl, Br) considered in the paper. At zero pressure, for BiTeI and BiTeCl compounds overestimation of a and c parameters is 2% and 7%, respectively. For BiTeBr compound—2% and 6%, respectively (see table 1). We also performed band structure calculations within PBEsol [30] functional and DFT-D3 [31] scheme (next we will denote it simply as PBEsol and DFT-D3).

For simulation surface under pressure, the large slab has been constructed on the basis relaxed parameters of bulk structure under the pressure. The large slab Hamiltonian derived from the bulk one is then used to calculate surface Green functions [34–36]. So, TB models were constructed using WANNIER90 code [32, 33]. The chosen basis consists of six spinor p-type orbitals for each atom: $| {p}_{x}^{\uparrow }\rangle$, $| {p}_{y}^{\uparrow }\rangle$, $| {p}_{z}^{\uparrow }\rangle$, $| {p}_{x}^{\downarrow }\rangle$, $| {p}_{y}^{\downarrow }\rangle$, $| {p}_{z}^{\downarrow }\rangle$. The low-lying s orbitals are not taken into consideration. To study the bulk electronic spectra near the point of TPT two TB Hamiltonians are constructed: one for the topologically trivial phase of BiTeI and the other for the TI phase, ${\hat{H}}_{\mathrm{triv}}$ and ${\hat{H}}_{\mathrm{top}}$, respectively. For each intermediate pressure the Hamiltonian is taken as a linear combination: $\hat{H}=\eta {\hat{H}}_{\mathrm{top}}+(1-\eta ){\hat{H}}_{\mathrm{triv}}$. Here $\eta =0$ and $\eta =1$ correspond to a pressure of 3.6 GPa and 4 GPa, respectively.

The band-bending effect was accounted by shifting the on-site matrix elements of surface atoms. The value of shifting is determined by the potential gradient obtained in the first-principles calculations for slabs of BiTeX (X = I, Cl, Br) at zero pressure [9, 15]. These values of band bending shift were applied for all surfaces of BiTeX (X = I, Cl, Br) under pressure.

The chosen method for calculating band structure under pressure is qualitatively similar to that of [25]. However, the main parameter in our calculations is pressure. As a result, we consider TPT in relation to measured quantity. While in [25], some parameter η has been used.

Dynamical properties of the bulk bismuth tellurohalides were calculated within the density-functional perturbation theory [37] in a mixed-basis pseudopotential approach [38–40]. The details of the calculation as well as the SOC implementation within the mixed-basis pseudopotential method can be found in [41, 42].

BiTeI crystalizes in a hexagonal structure, P3m1 (figure 1(a)). At zero pressure, the compound is a semiconductor. Using fully optimized crystal structure parameters (a, c and atomic coordinates) within PBE exchange-correlation functional, we obtained the band gap width of 408 meV which is in good agreement both with experimental measurements (380 meV [7]) and GW-calculations (400 meV in the LDA+GW scheme [13]). Such an agreement which is unusual for PBE is explained by some overestimation of structural parameters in the relaxation process (see table 1). Using of PBEsol or DFT-D3 leads to accurate description of atomic lattice parameters and strong underestimation of band gap width and CPTs pressure as well. Hence, next we describe bulk and surface electronic structure with PBE exchange-correlation functional only. However, we give all main results obtained with PBEsol, DFT-D3 and PBE in table 1.

Zoom In Zoom Out Reset image size

The valence and conduction gap edge bands are composed of tellurium and bismuth p states (figure 1(c)), respectively. A strong spin–orbit interaction leads to a large Bychkov–Rashba type spin splitting of the bulk states and thus two pairs of extrema along the H–A–L line are formed: (${E}_{\mathrm{gap}}^{{\rm{H}}-{\rm{A}}}$) and (${E}_{\mathrm{gap}}^{{\rm{A}}-{\rm{L}}}$). Under pressure the value of ${E}_{\mathrm{gap}}^{{\rm{H}}-{\rm{A}}}$ diminishes and at a certain pressure, P_c, shrinks to zero that indicates normal semiconductor–Weyl semimetal transition. These does not happen with gap along the A–L direction, ${E}_{\mathrm{gap}}^{{\rm{A}}-{\rm{L}}}$. Upon further increase of pressure a band gap appears with inverted edges in the vicinity of the A point: now the lowest conduction band is formed by tellurium p states while the top valence band consists of bismuth p orbitals (figure 1(c)). The calculated value of P_c equal to ∼3.73 GPa (figure 1(d)) agrees well with the experimental one, P_c = 3.5 GPa [22, 23].

Due to a small overestimation of the band gap compared to the experimental value (∼7%), the P_c parameter is also slightly overestimated. The calculated value of P_c = 3.5 GPa contradicts P_c = 10 GPa obtained in the GW calculation [27]. As a result according to the [27] the P_c = 10 GPa indicates the impossibility of TPT in BiTeI due to the CPT at ∼9 GPa. The discrepancy is explained by a strong overestimation of the band gap width at zero pressure in the GW calculations.

The behavior of band gap ${E}_{\mathrm{gap}}^{{\rm{H}}-{\rm{A}}}$ and the trajectory of crossing (Weyl) points near the TPT are plotted in figures 1(e) and (f),respectively. ${E}_{\mathrm{gap}}^{\min }$ is the smallest value of ${E}_{\mathrm{gap}}^{{\rm{H}}-{\rm{A}}}$ in the entire BZ. The position of Weyl points in the reciprocal space depends on the value of pressure. The Weyl points are formed at the BZ boundary, ${k}_{z}=\pi /c$, along each AH direction. At the phase transition pressure, P_c, each Weyl point splits up into pair of nodes with a clockwise and counterclockwise propagation on the (${k}_{x},{k}_{y}$) plane (figure 1(f)). Also, these two nodes shift in opposite directions along k_z. The Weyl points are observed up to ∼3.93 GPa and then the band degeneracy is lifted and BiTeI converts into a TI. Thus, the Weyl semimetal phase exists within the pressure interval from 3.73 to 3.93 GPa. The range is not small, so the phase can be experimentally observed.

The pairs of Weyl points appeared along the AH directions have opposite topological charges (chirality) which are defined by Chern numbers as the flux of Berry phase gauge field over a sphere around each Weyl point: ${C}_{n}=1/(2\pi )\oint {{\boldsymbol{\Omega }}}_{n}({\bf{k}}){\bf{n}}{\rm{d}}S$, where ${\bf{n}}$—surface normal vector and ${{\boldsymbol{\Omega }}}_{n}({\bf{k}})={{\boldsymbol{\nabla }}}_{k}\times {{\bf{A}}}_{n}({\bf{k}})$ is the Berry curvature. Here ${{\bf{A}}}_{n}({\bf{k}})$ is the Berry connection defined as ${{\bf{A}}}_{n}={\rm{i}}\langle n{\bf{k}}| {\boldsymbol{\nabla }}| n{\bf{k}}\rangle$ for nth Bloch state, $| n{\bf{k}}\rangle$, (in our case—the highest occupied band) with quantum number ${\bf{k}}$. The calculation reveals that the Weyl points which move clockwise (counterclockwise) have a positive (negative) Chern number, C_n = 1 (−1), that corresponds to drain (source) points of the Berry gauge field. The in-plane components of the Berry curvature at ${k}_{z}=\pi /c-0.012$ Å⁻¹ are demonstrated in figure 2.

Zoom In Zoom Out Reset image size

To illustrate surface electronic properties of the systems with Bychkov–Rashba-type spin splitting of bands in the topological and Weyl semimetal phases we considered within TB method a hypothetical model of the bismuth tellurohalides surface under pressure. As was mentioned above to construct the surface under pressure we used large slab with fully relaxed parameters of the bulk structure under the pressure. To do the model more realistic, we included the effect of band bending because of lack of inversion symmetry in the bismuth tellurohalides.

The TPT is accompanied by qualitative changes in the surface electronic structure. The BiTeI has two possible surface terminations: iodine or tellurium with positive and negative band bending, respectively [9, 10]. Let us first trace the evolution of surface electronic states with pressure without considering the effect of band bending (see figures 3(a)–(f)). In the absence of external pressure (figures 3(a) and (d)), the surface states on the Te(I)-terminated surface appear along the outer (inner) brunch of unoccupied spin–orbit split bulk bands and are mainly localized in the topmost trilayer.

Zoom In Zoom Out Reset image size

In the pressure-induced Weyl semimetal phase (figures 3(b) and (e)), the surface states on the Te-terminated surface remain along the outer edge of the bulk bands while in the case of iodine termination they shift into the band gap and exhibit a cone-like dispersion with crossing at the BZ center. The surface states touch the crossing (Weyl) points of the conduction and valence bulk bands. In the topological phase (figures 3(c) and (f)), the surface states become gapless linking the valence and conduction bands. In the case of I-terminated surface, the crossing of two surface states occurs in the band gap slightly above the bulk conduction band minima while on the Te-terminated surface the crossing (Dirac) point is inside off the valence bands due to the mixing of surface and bulk electronic bands.

The surface electronic structure obtained with taking into account the effect of band bending as well as its evolution under pressure is different. The overall band structure is strongly modified because a set of well-defined Bychkov–Rashba type spin split surface bands appearing due to the band-bending effect [20]. In the topologically trivial phase (figures 3(g) and (j)), the lowest surface band on the Te-terminated surface lies within the energy gap and is strongly localized in the topmost trilayer. The states with higher energies are confined within three upper trilayers. On the I-terminated surface, a set of surface bands splits off from the valence bulk band. The occupied surface electronic states appear mainly inside or at the edge of the bulk bands except the topmost state which is shifted upwards into the energy gap. Such a hierarchy of surface states was also observed in ARPES experiments [11, 14]. In the Weyl semimetal phase (figures 3(h) and (k)), when passing through the bulk Weyl nodes the surface states localized in the upper trilayer become discontinuous. In the topological phase (figures 3(i) and (l)), the surface states are breaking up into two parts near the bulk edge extrema, and gapless surface states with a linear dispersion in the band gap are formed. On both surface terminations, the crossing point of the gapless surface states is inside the bulk bands where they turn into surface resonances.

The surface electronic structure of Weyl semimetals is characterized by the presence of Fermi arcs connecting the Weyl points with opposite chirality. Figure 4 shows the electronic spectra at the Fermi level of the (001) surface for both terminations. The isoenergetic surface spectra were calculated with and without the account of the band-bending at a pressure of 3.84 GPa which corresponds to a Weyl semimetal phase (see figure 1(e)). Also shown are the positions of Weyl points on the (${k}_{x},{k}_{y}$) plane. The pairs of Weyl points lie on the trajectory connecting the positions of Weyl points at the P_c pressure (3.73 GPa) and at the annihilation point (3.93 GPa) (see figure 1(f)).

Zoom In Zoom Out Reset image size

Without taking band bending into account the Fermi arcs link up the pair of Weyl points induced originally along the A–H ($\overline{{\rm{\Gamma }}}$–$\overline{{\rm{K}}}$) direction of the bulk (surface) BZ. The form of Fermi arcs depends on the surface termination. On the Te-terminated surface, they form arc-like curves coinciding with the trajectory of a pair of nodes which move clockwise and counterclockwise over the BZ. In the case of I-termination, the Fermi arcs connect the pairs of Weyl points by a straight line. On further increase in pressure, the length of Fermi arcs increases and when the Weyl semimetal–TI phase transition occurs the arcs connect with each other forming a gapless topological surface state.

With the account of band bending, in addition to the Fermi arcs a circular surface state appears which originates from the outer branches of Bychkov–Rashba type spin split surface states (figures 3(h) and (k)). In the case of I-termination, the outer circular state is found farther from the surface BZ center as compared to the Te-terminated surface. However the isoenergetic curves for both terminations are qualitatively similar, unlike the previous case. In addition, with the account of band bending the Fermi arcs link up the Weyl points from neighboring pairs. Upon increasing pressure, the length of Fermi arcs decreases up to zero at ∼3.9 GPa when the transition to topological phase occurs. After the transition the isoenergetic electronic structure is presented by a single closed curve of the topological surface state.

The formation of Fermi arcs is closely related to the spin texture of the surface states. Without account of band bending, the single closed contour in the topological phase has a clockwise spin helicity [21]. In the Weyl semimetal phase, Fermi arcs inherit this behavior—the helicity also has a clockwise character. The spin texture 'connects' the source and drain Weyl points in the clockwise order, so the Fermi arcs link up the pair of points emerged along the A–H direction.

Subject to the band bending, the situation is opposite. The Fermi arcs get the spin texture from the inner closed contour of Bychkov–Rashba type spin split surface states in the trivial phase which have a counterclockwise helicity [10]. The spin texture 'connects' the source and drain Weyl points in counterclockwise order, so the Weyl points from neighboring pairs are connected.

Thus, in the calculations the modification of the surface state under pressure depends on whether the band-bending effect is taken into account or not. The modifications are closely related to the alteration of the Fermi arcs which characterize the isoenergetic spectra of Weyl semimetals. Our additional calculations revealed that this picture does not change qualitatively with a strong modification of band bending potential. Also, all prognosticated peculiarities of the BiTeI surface electronic structure in topological phase can be revealed in other noncentrosymmetric alloys where the appearance of the topological or trivial phases can depends, for instance, on the concentrations. Note, such systems were proposed in [25] and shown that in cases LaBi${}_{1-x}$Sb_xTe₃, LuBi${}_{1-x}$Sb_xTe₃ the topological phase can be induced by variation of the Sb concentration, x.

Since BiTeI undergoes CPTs under hydrostatic compression [22, 24, 26, 27] it is necessary to compare the pressures corresponding to the TPT and CPT. The x-ray diffraction data [22, 24] showed that BiTeI remained in the P3m1 (ambient pressure) structure up to 8–9 GPa. The stability of the hexagonal phase up to ∼9 GPa was also confirmed by high-pressure Raman spectra measurements [26, 27]. For the higher pressure phase the orthorhombic Pnma structure with 12 atoms per unit cell were suggested (see figure 5) [24]. As our calculation show, this structure is normal non-direct gap semiconductor. We calculated the difference in enthalpy between the hexagonal P3m1 and orthorhombic Pnma phases depending on pressure (figure 6(a)). As follows from the calculation the CPT should occur at ∼6.5 GPa for PBE calculatoins. This value is somewhat smaller than the experimental CPT pressure. Despite this discrepancy, the CPT takes place after the TPT and at pressures which are beyond the pressure range (3.73–3.93 GPa) where the Weyl semimetal phase exists. Thus the CPT is accompanied by a TPT from a TI to a normal semiconducting phase.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

It should be noted, that using of PBEsol and DFT-D3 leads to strong underestimation of CPT pressure. So, using of PBE functional is more preferable for theoretical evaluation both TPT and CPT.

To additionally verify the proposed structure we have also calculated Raman mode frequencies for the P3m1 phase up to 8 GPa and for the Pnma structure at pressures 10, 12, and 14 GPa (figure 6(b)). In the hexagonal P3m1 phase with 3 atoms per unit cell, there are two A₁ and two E (twofold degenerate) zone center optical modes, which are both Raman and infrared active because of the lack of inversion symmetry. At pressures below the CPT, the theoretical data are in excellent agreement with experimental results with the exception of mode E¹ between 4–8 GPa where experimental data are absent apparently due to decreasing of the mode intensity. In the experimental Raman spectra [26, 27] three modes are visible up to 8 GPa with an expected increase upwards with pressure. A sudden change in the number and frequency of Raman active modes which points out to a structural transition is observed at ∼9 GPa.

The Raman spectrum for the orthorhombic phase has a rather complicated character. All the modes are non degenerate. Some of them coincide with the experimental data. In general, variation of Raman mode frequencies with pressure (10–12–14 GPa) resembles the behavior of the experimental Raman modes. However, the number of the calculated Raman modes in the Pnma structure (13 are shown in the figure while the total number of Raman modes is equal to 18) is almost twice as large as the experimental data. Although the experimental Raman spectrum at high pressures bears some similarity to the calculated Pnma spectrum the measured modes apparently correspond to a more symmetric structure. The character of the measured spectrum indicates some orthorhombic phase with similar enthalpy.

Another candidate [27] can be the hexagonal ${\rm{P}}{6}_{3}{\rm{m}}{\rm{c}}$ structure with 6 atoms per unit cell which is the structure of another compound, BiTeCl, at ambient pressure [43]. The Raman active modes calculated for BiTeI in the hexagonal ${\rm{P}}{6}_{3}{\rm{m}}{\rm{c}}$ structure at 10 and 12 GPa are continuations of the corresponding modes in the low-pressure phase. The number of E modes is doubled and the two lower E modes are shifted down in energy compared to the E¹. So the Raman spectrum does not reproduce the high-pressure experimental data and the ${\rm{P}}{6}_{3}{\rm{m}}{\rm{c}}$ structure cannot be considered as a high-pressure phase of BiTeI.

For BiTeBr and BiTeCl there is no experimental information about the possibility of TPT under pressure. It was thought that BiTeBr crystallizes in the hexagonal structure, ${\rm{P}}\overline{3}$ m1, with Te and Br atoms statistically distributed within two layers adjacent to the Bi layer [43]. Later an ordered stacking of Te and Br sublattices, like in BiTeI, was confirmed by x-ray diffraction measurements [28] and, moreover, ARPES images revealed a well-defined Bychkov–Rashba type spin split states [14]. The structure of BiTeCl at ambient pressure is hexagonal, ${\rm{P}}{6}_{3}{\rm{m}}$c, with 6 atoms per unit cell [43]. We also considered BiTeCl in the P3m1 crystal phase which was found in the calculation to be energy preferable up to the CPT. To evaluate the probability of CPT we also examine both compounds in the orthorhombic Pnma structure which was originally proposed as a high-pressure structure for BiTeI [24]. This phase is considered for BiTeBr and BiTeCl compounds because most V-VI-VII semiconductors crystallize in this structure [44–47].

As follows from figure 7(a) and table 1, the CPT in BiTeBr takes place at pressure of ∼6.3 GPa within PBE calculations. The band gap in trivial phase (figure 7(b)) decreases with pressure until CPT occurs. The CPT pressure found with PBE is in good agreement with experimental results of [28]. The data exclude the possibility of TPT in the low-pressure phase of BiTeBr. Note, if another high-pressure structure which is energy preferable compared to the P3m1 and Pnma structures is realized then the pressure of CPT would be lower than 6.3 GPa and in any case the CPT occurs before a TPT.

Zoom In Zoom Out Reset image size

Due to relatively close pressures of TPT and CPT and lack of experimental investigations of both transitions, additional calculations have been performed within PBEsol and DFT-D3 (table 1). Using of PBEsol (DFT-D3) lead to shifting CPT and TPT to smaller pressures: ∼4 GPa (∼4 GPa) and ∼4 GPa (∼6–7 GPa), respectively. This agrees with BiTeI calculation results. However, the conclusion is preserved: the TPT is absent in BiTeBr.

We also calculated Raman modes in P3m1 crystal phase of BiTeBr under pressure which are shown in figure 6(c). The dispersion of the modes is in a good agreement with the experimental results. The Raman shift values in BiTeBr are similar to those in BiTeI (see figures 6(b) and (c)) except for ${A}_{1}^{1}$ mode which is ≈10% higher than in BiTeI. This peculiarity is observed both in the experimental and theoretical data.

The case of BiTeCl is more complicated from the crystallographic point of view. First of all, our results indicate the possibility of CPT from P6₃mc to the P3m1 structure at very small pressures because the calculated difference in enthalpy (figure 7(c)), red line) turns out to be ∼0.8 meV even at zero pressure. With increasing pressure, the P3m1 phase remains energy preferable up to ∼2–4 GPa (table 1) when the CPT to the Pnma structure takes place.

However, these results do not affect the general conclusion about a possibility of TPT in this compound. A TPT does not occur either in the P3m1 or in the P6₃mc crystal phases because although the band gap width decreases with pressure (figure 7(d)) this decreasing practically ceases at a value of ∼100 meV (∼80 meV) for the P3m1 (P6₃mc) structure. A relatively large band gap width makes impossible the TPT both in the P3m1 and in the P6₃mc hexagonal phase. This conclusion is preserved for all exchange-correlation approximations used: PBE, PBEsol and DFT-D3.

This conclusion contradicts the data reported in [48] where by APRES measurements a TI phase in BiTeCl with Dirac-type surface states was found at zero pressure. The observation of the TI phase might be explained by a possible crystallization of the material in another hexagonal structure what is confirmed by a very small bulk band gap width (∼220 meV) obtained in the experiment. The latter contradicts the electronic structure calculations within the GW-approximation that reveal the band gap width of 800–900 meV [13]. For additional verification of the trivial character of BiTeCl in the hexagonal phases, we have calculated the surface electronic structure at a pressure of 18 GPa with taking into account the band bending effect. Figures 7(e) and (f) show the electronic spectra of the Te-terminated surface for the P3m1 and P6₃mc crystal phases, respectively. It is obvious that there are no topological surface states in the bulk band gap and therefore the system is a trivial semiconductor. Both crystal phases are only characterized by the presence of well known Bychkov–Rashba type split surface states. It should be noted, that additional extrema of the conduction band appear, one at the $\overline{K}$ point and the other in the $\overline{{\rm{\Gamma }}}$–$\overline{M}$ direction which results in an indirect band gap in the case of the P6₃mc structure.