# Bulk and surface spin conductivity in topological insulators with hexagonal warping

###### Abstract

We investigate the spin conductivity of topological insulators taking into account both the surface and quasi-two-dimensional bulk states. We apply a low-energy expansion of the Hamiltonian up to the third order in momentum and take into account the vertex corrections arising due to the short range disorder. Hexagonal warping gives rise to the additional anisotropic components in the spin conductivity tensor. Typically, isotropic part of the spin conductivity is larger than anisotropic one. The helical regime for the bulk states, in which the electrons in the Fermi level have the same projection of the spin on the direction of momentum, have been studied in a more detail. In this regime, a substantial increase of the spin conductivity contribution from the bulk states at the Fermi level is observed. We find that the bulk spin conductivity is insensitive to disorder if Rashba spin-orbit coupling is larger than disorder strength, otherwise, it is strongly suppressed. The contribution to the spin conductivity from the surface states is almost independent of the chemical potential, robust to disorder and its value is comparable to the spin conductivity contribution from the bulk states per layer. The obtained results are in agreement with experimental data.

###### pacs:

03.67.Lx, 74.90.+n## I Introduction

Topologically-protected surface states form a Dirac cone in the electronic spectrum of the topological insulators (TI) Hasan and Kane (2010). The electron dispersion near Dirac points is linear. However, a hexagonal warping of the Dirac cone arises when we take into account the next-order terms in the momentum expansion of corresponding Hamiltonian of the TIs with the hexagonal lattices, such as BiTe Fu (2009) and BiSe Kuroda et al. (2010). The hexagonal warping influences not only the surface states but also quasi-two-dimensional bulk states in these systems. Effects of the hexagonal warping on the electronic properties of the TI have been studied extensively Wang and Yu (2011); Pal et al. (2012); Siu et al. (2014); Repin and Burmistrov (2015). In our recent paper Akzyanov and Rakhmanov (2018), we find that the presence of the hexagonal warping significantly affects the charge conductivity of the TI. In particular, it gives rise to the anisotropic anomalous in-plane magnetoresistance. Hexagonal warping also affects the quantum anomalous Hall effect and anomalous out-of-plane magnetoresistance.

A remarkable feature of the TIs is the existence of high spin conductivity in the absence of magnetic field, which is associated with an intrinsic spin Hall effect Sinova et al. (2015). This effect has been first predicted in Rashba spin-orbit coupled materials, such as GaAs Murakami et al. (2003); Sinova et al. (2004). However, the intrinsic spin Hall effect in such materials is weak due to short-range disorder (from a theoretical point of view, due to vertex corrections caused by this disorder) Inoue et al. (2004); Raimondi and Schwab (2005).

A change of direction of the magnetization in the magnetic material by a spin current is referred to as spin-transfer torque (STT) Ralph and Stiles (2008). The STT is closely related to the spin conductivity Chen et al. (2015). This effect can be used for the design of the fast and low dissipative magnetic memory Wang et al. (2013). Recent experiments reveal that STT in the TIs is by orders of magnitude larger than for any other material, which is a sign of a substantial spin conductivity in TIs Mellnik et al. (2014); Fan et al. (2014, 2016); Han et al. (2017). Experimental study of the STT in the TIs demonstrates some intriguing features. Both the in-plane and out-of-plane STT exist in the system, and the value of these effects is of the same order, which is unexpected from the spin-momentum locking argument Mellnik et al. (2014). Moreover, the sign of the spin conductivity may be different in different samples of the same material Yang et al. (2016). Spin conductivity in the TI is tuned by chemical potential and obeys a particle-hole asymmetry Kondou et al. (2016); Fan et al. (2016). Also, spin conductivity is suppressed in the bulk-insulating regime Kondou et al. (2016). It has been speculated that the large spin currents arise in the TIs due to the existence of the topologically-protected surface states Shiomi et al. (2014); Wang et al. (2015, 2016). However, in the other papers it is complained that the spin conductivity in the TI mainly comes from the bulk states Jamali et al. (2015); Han et al. (2017).

In general, the spin conductivity includes both contributions from the states at the Fermi surface and from all filled states Yang and Chang (2006); Kodderitzsch et al. (2015). While the contribution to the spin conductivity from the filled states can be calculated in a clean limit Murakami et al. (2003); Sinova et al. (2004), it is vital to treat the disorder correctly to describe the contribution from the states at Fermi level Inoue et al. (2004); Raimondi and Schwab (2005).

Unexpectedly small number of theoretical works are devoted to the spin conductivity in TIs. Recent DFT calculations of the contribution to the spin conductivity from the filled states show that quite large spin currents can exist in BiSb and the value of the spin conductivity can be tuned by the chemical potential variation Sahin and Flatte (2015). In Ref. Peng et al., 2016, the spin conductivity of the surface states in a thin film of TI with a cubic lattice has been studied neglecting the vertex corrections. The authors concluded that the dependence of the surface spin conductivity on the disorder and chemical potential is small. The spin conductivity in another Dirac material, graphene, attracted much more attention Sinitsyn et al. (2006); Liu et al. (2015); Garcia and Rappoport (2016). Provided the spin-orbit interaction is induced in graphene, quite reasonable spin currents can be obtained in it. Recent calculations also show that large spin currents can be induced in Weyl semimetal, another Dirac material with large spin-orbit interaction Sun et al. (2016).

We study the spin conductivity of the surface and bulk states in the TI in a low energy approximation with taking into account the hexagonal warping. Both contributions from the filled states and from the states at the Fermi surface are considered. We apply the Kubo formalism accounting the vertex corrections to the velocity operators arising due to the short-range disorder. We show that the presence of the hexagonal warping leads to the additional anisotropic terms in the spin conductivity. We get that the spin conductivity is robust against disorder. The spin conductivity of the surface states is comparable with the spin conductivity of the bulk states per layer. The obtained results are consistent with the experimental data.

The paper is organized as follows. In Section II we analyze the Hamiltonian describing the surface and bulk states in the TI. In Section III we introduce disorder and in Section IV calculate the vertex corrections to the velocity operator. In Section V we study the contribution to the bulk and surface spin conductivity from the states at the Fermi level. In Section VI we consider the contribution to the spin conductivity from the filled states. We estimate the values of the characteristic for TIs parameters in Section VII. In Sec. VIII we discuss the obtained results and compare them with the experiments and numerical calculations.

## Ii Model

Low energy surface and quasi-two-dimensional bulk states in the TI can be described by the Hamiltonian Liu et al. (2010) ()

where are the Pauli matrices acting in spin space, is the chemical potential, is the value of Rashba coupling (equal to the Fermi velocity for the surface states), is the inverse mass term, characterizes the next order correction in momentum to , and are the in-plane momentum components, is the hexagonal warping coefficient. The term in the Hamiltonian responsible for the hexagonal warping can be rewritten as and the Hamiltonian is invariant under rotation on the angle . The spectrum of the Hamiltonian (II) is given by

(2) |

If we measure the energy in terms of , then, the chemical potential, the next order correction to the spin-orbit coupling, and the hexagonal warping are conveniently characterized by the dimensionless values , , and , respectively.

Energy spectrum (2) is shown in Fig. 1 for different set of parameters characteristic of the bulk, (a) and (c), and surface, (e), states. A key feature of the surface states in the TI is the existence of a robust Dirac cone, which is the case if is sufficiently large, Fig. 1(e). Corresponding Fermi surface has a characteristic form of a snow-flake, Fig. 1(f). The bulk states corresponds to smaller values of , Figs. 1(a) and (c). In the latter case, the spectrum has an appearance characteristic of a two-dimensional electron gas with bands splitted due to the Rashba spin-orbit interaction. Corresponding Fermi surfaces with two pockets are shown in Figs. 1(b) and (d) for different values of the chemical potential. Note, that this model describes well ARPES data for the surface and bulk states in the TIs Fu (2009); Kuroda et al. (2010); Liu et al. (2010).

We obtain from Eq. (2) that the robust Dirac cone exists when . For , two spin split bands emerge in the system as it is expected for the Rashba spin-coupled electron gas. Therefore, we can formally to write down that

(3) | |||||

In the case of the bulk states, the function decreases at large and a proper momentum cut-off must be introduced to avoid arising fake Fermi surface pockets at large momentum. We define cut-off momentum as

We can calculate average spin projection of electrons as , where is the spin operator and are eigenfunctions corresponding to the bands . The in-plane spin polarization component is schematically shown in Figs. 1(b), (d), and (f). The calculated spin polarization lies in the plane if we neglect the hexagonal warping. We see that each band can be characterized by helicity, that is, the sign of the projection of the spin on the direction of momentum. The z-component of the spin polarization arises if we take into account that . If the bulk states have two splitted Fermi surfaces with different helicity, see panel (b). In the case of two Fermi surfaces have the same helicity, see panel (d), so, we call this regime for the bulk states as helical. Surface states are helical in any case, see panel (f).

In general, the spin conductivity can be presented as a sum of three terms Yang and Chang (2006); Kodderitzsch et al. (2015)

(4) |

where the first two items correspond to a contribution from the states at the Fermi surface and the third one from the filled states. Here and denote the in-plane coordinates and , respectively, and denotes the spin projection.

At zero temperature and can be written in the form Inoue et al. (2004); Yang and Chang (2006)

(5) | |||

(6) |

Here , is the is the velocity operator, is the velocity operator with vertex corrections, means the anticommutator, and are the retarded and advanced Green functions, which will be specified in the next section.

The contribution to the spin conductivity from the filled states is Sinova et al. (2004); Sinitsyn et al. (2004)

(7) |

Here is the energy of an electron in the -th band with the momentum , is the corresponding Bloch vector, , is the Fermi distribution function corresponding to (which is the Heaviside step function in the considered case of zero temperature), means impurity averaged, and is the disorder parameter or scattering rate. The latter will be also specified in the next section.

## Iii Disorder

We will describe disorder by a potential , where is the Dirac delta function, are positions of the randomly distributed point-like impurities with the local potential and concentration . We assume that the disorder is Gaussian, that is, and .

In the self-consistent Born approximation (SCBA), the impurity-averaged Green’s functions can be calculated as

(8) |

where are bare Green’s functions of the Hamiltonian (II)

(9) |

and is the self-energy, which is defined as

(10) |

In the case under consideration, we can calculate the self-energy using an expression similar to that derived in Ref. Shon and Ando, 1998

(11) |

The function under integral in Eq. (11) decays as when . Thus, the value of this integral is determined by zeros of the denominator.

The value is usually referred to as a disorder parameter or scattering rate. It determines the analytical properties of the Green’s functions , while is only a small correction to the chemical potential since we consider here only the case of small disorder. Thus, we can neglect the real part of the self-energy with the exception of some singular point, which will be specified below. If we put and have in mind that is small in the limit of small disorder, we derive from Eq. (11) an explicit formula for the scattering rate

(12) |

### iii.1 Bulk states

First, we consider the bulk states. In the simplest case, when and tends to zero, we can derive an explicit formula for the scattering rate in two opposite limits, and . If the chemical potential is negative, we obtain from Eq. (12) following Ref. Kato et al., 2007,

(13) |

This value is independent of the chemical potential and the strength of the spin-orbit interaction. It is convenient to introduce dimensionless disorder parameter . In these notations .

Effects of the spin-orbit coupling are not smeared by the disorder if or, equivalently, . We will call further the spin-orbit coupling strong if condition is satisfied. Otherwise, , the spin-orbit coupling is weak.

In the helical regime, , the behavior of depends on the system parameters. If the spin-orbit coupling is weak, , we get that and rapidly decays to zero with an increase of . In the opposite limit of strong strong spin-orbit coupling, , we found from Eq. (12) that the scattering rate increases if :

(14) |

When the chemical potential attains the singularity point, , the Fermi level crosses the bottom of the energy bands if and [see Eq. (2)]. If , the Fermi level occurs in the energy gap. If we apply self-consistent Eq. (11) we get that at the real part of the self-energy vanishes and

(15) |

Thus, in the helical regime the scattering rate increases significantly if the spin-orbit coupling is strong, .

In a more general case, the scattering rate for the bulk states was calculated numerically using Eq. (11). The results are shown in Fig. 2 for the case of the strong spin-orbit coupling characteristic of the TIs. As we can see from the figure, the higher order corrections to the spin-orbit coupling and the hexagonal warping has a little impact on the value of the scattering rate in the case of the bulk states. In particular, a characteristic peak in arises near the point .

### iii.2 Surface states

For the surface states, we neglect a correction to the value of due to the real part of the self-energy in the limit of weak disorder, similar to the case of the bulk states. Thus, we can use Eq. (12) to calculate .

In the simplest case we obtain a well-known result Hu et al. (2008); Chiba et al. (2017),

(16) |

When the chemical potential crosses the Dirac point, , we apply Eq. (11) and find that the real part of the self-energy vanishes while imaginary part is exponentially small Hu et al. (2008) , where is the momentum cut-off. That is, the scattering rate at the Dirac point is exponentially suppressed in the case of the strong spin-orbit coupling, .

In a more general case, the scattering rate for the surface states was calculated numerically with the help of Eq. (11). The dependence of is shown in Fig. 3. We see that the scattering rate is almost particle-hole symmetric since the spectrum of the surface states close to such a symmetry [see Fig. 1(e)].

The self-energy is proportional to the identity matrix. This allows to obtain an explicit expression for the impurity averaged Green function. We can rewrite Eq. (8) as or

(17) |

Therefore, the expression for is given by an equation similar to Eq. (9) for , in which is replaced by . We characterize disorder by a single value neglecting renormalization of the chemical potential.

## Iv Vertex corrections

In the SCBA, following the approach described in Ref. Shon and Ando, 1998, we can derive an equation for the vertex corrected velocity operator Chiba et al. (2017)

(18) |

We present here the derivation of ; results for can be obtianed just by the substitution .

It is easy to show that and , where and are scalars. In these notations we obtain from Eq. (18) that

(19) |

We begin our consideration with the bulk states and derive some analytical results in the simplest case, when and are zero. Under such conditions and if , the vertex corrected spin-orbit coupling is small, , either at weak, , or strong, , spin-orbit coupling in accordance with the results of previous works (see, e.g., Refs. Inoue et al., 2004; Raimondi and Schwab, 2005). However, in the case of strong spin-orbit coupling, when the chemical potential is positive and lies in the interval , the vertex corrected is comparable to its bare value

(20) |

Numerically calculated dependence of the vertex correction to the spin-orbit coupling for the bulk states on the chemical potential, , is shown in Fig. 4 for a more general situation. As we can see, the higher order corrections to the spin-orbit coupling and the hexagonal warping significantly enhances in the region . Nevertheless, its value is still much smaller than the bare value . When (in the helical regime), is of the order of and almost independent of and .

For the surface states, in the case and , we get that away from the Dirac point, , the vertex correction is , while vanishes at , as it have been obtained in Refs. Adroguer et al., 2012; Chiba et al., 2017. Numerically calculated dependence of the vertex correction on the chemical potential is shown in Fig 5 for the parameters characteristic of the surface states. The presence of the hexagonal warping slightly increases , while the existence of the finite mass term leads to the particle-hole asymmetry. Taking into account correction to the spin-orbit coupling results in a significant increase of if the chemical potential is away from the Dirac point.

## V Spin conductivity from the states at the Fermi surface

Now we use the results obtained in the previous sections and Eqs. (5) and (6) to calculate the contribution to the spin conductivity due to the states at the Fermi surface. On this way, we obtained that in the considered approach the term vanishes exactly. Thus, we should to compute only the term .

Isotropic tensor component is the only term that persists in the system in the case of zero hexagonal warping. All other components are anisotropic and they are non-zero only if . The measured value of the spin conductivity depends on the mutual orientation of the current and the crystal axes. So, it is convenient to relate the spin conductivity tensor components in the crystal axes with that related to the current direction, . New coordinates are obtained by anticlockwise rotation by the angle along the crystal axes. We assume that the current is directed along -axis. In this coordinates we have

(21) | |||||

Therefore, it is sufficient to calculate and .

We derive from Eq. (5)

(22) |

where is the spin conductivity quanta and is given by Eq. (2). As we can see, the isotropic spin conductivity component is proportional to if we neglect the higher-order correction to the spin-orbit coupling. This value increases significantly in the helical state . The anisotropic component is proportional to the hexagonal warping strength , the vertex corrections does not affect it.

First, we calculate the contribution to the spin conductivity from the bulk states. The results are shown in Fig. 6. As we can see from the top panel in Fig. 6, the isotropic spin conductivity component is suppressed when the chemical potential is negative. It occurs since the vertex correction is small in this region of . However, when , the value of is comparable to and the value of increases significantly. Note that and produces a weak effect on in this range of . If , the spin conductivity vanishes since the density of states on the Fermi disappears.

The results for the anisotropic component of the spin conductivity are presented in the bottom panel of Fig. 6. This value decreases almost linearly with an increase of the chemical potential if and, when , it has a small peak near . The value also demonstrate almost a linear growth with an increase of the coefficients and .

The dependencies of the isotropic and anisotropic components of the spin conductivity on for the surface states are shown in Fig. 7. Both these values have minima at the Dirac point and particle-hole asymmetry that is smaller for larger . Note that and for the surface states decreases with an increase of the next order correction to the spin-orbit coupling coefficient .

The effect of disorder on the spin conductivity is illustrated in Fig. 8. Both surface and bulk conductivities are robust against disorder in the weak scattering limit, . Moreover, the topologically protected surface terms are robust even in the case of higher disorder, , while the components of the bulk conductivity decrease significantly in this limit. However, the SCBA is not correct for a strong disorder and more advanced techniques are required to study the robustness of the spin conductivity of the surface states in such regime.

## Vi Spin conductivity from the filled states

Here we calculate the contribution to the spin conductivity from the filled states using Eq. (II) and the obtained above results for the disorder parameter . Note that the vertex corrections do not affect this part of the spin conductivity. Similar to the spin conductivity from the states at the Fermi surface, the spin conductivity from the filled states has an isotropic component and anisotropic one . The isotropic component, , is the only term that persists in the system in the absence of the hexagonal warping. The anisotropic components are non-zero only if . In the rotated coordinates , the components of tensor transform similar to , see Eq. (V).

We start with the bulk spin conductivity. In the clean limit, , and zero third order corrections, and , we get that if and if . The latter relation is a well-known result for a spin hall conductivity Sinova et al. (2004); Sinitsyn et al. (2004). In a more general case, the results were obtained numerically and presented in Fig. 9. As we can see from upper panel of this figure, the hexagonal warping has a little effect on the value of the isotropic spin conductivity, while the correction to the spin-orbit coupling enhances it. According to bottom panel of Fig. 9, the anisotropic part of the spin conductivity, , decays monotonically with the increase of the chemical potential. It increases almost linearly with the increase of the hexagonal warping strength . The bulk spin conductivity becomes zero when the chemical potential crosses the bottom of the conduction band and occurs in the energy gap.

The results for the contribution to the spin conductivity from the surface states are shown in Fig. 10. We see that both and have maxima at the Dirac point and decreases with the increase of (in contrast to the contribution from the states at the Fermi level, Fig. 7). These functions are more or less particle-hole symmetric. From Fig. 10, we see that the value of isotropic spin conductivity decreases with the increase of the higher order momentum corrections and . The value of anisotropic spin conductivity increases with the increase of hexagonal warping strength and decreases with the increase of .

The dependence of on the disorder parameter is shown in Fig. 11. We obtain that the spin conductivity from the filled states is robust against disorder if . If the disorder is stronger, , both the bulk and surface conductivities are suppressed.

## Vii Evaluation of characteristic parameters

In this section, we demonstrate that the values of the parameters used above for the calculation of the spin conductivity are reasonable.

We can extract information on the disorder strength from Ref. Chen et al., 2013. The imaginary part of the self-energy for the surface states in BiTe can be estimated from the ARPES data presented in Ref. Chen et al., 2013 as a half-width of the quasiparticle peak: meV and peak position corresponds to meV. Thus, we get . This is an upper limit for the disorder strength, since, for example, electron-phonon and electron-electron interactions also contribute to the blurring of the quasiparticle peak. The alternative indirect estimate we obtain as follows. The STM data from Ref. Cheng et al., 2010 shows that for a clean surface of BiTe there exists one defect approximately per . We suppose that a typical impurity potential is of the order of the chemical potential (which was about 200 meV). This assumption is true, e.g., for vacations. The Fermi velocity for the surface states was evaluated in Ref. Zhang et al., 2009 as eV. Then, we get . The value of the Rashba spin-orbit coupling is comparable for the bulk and surface states King et al. (2011). Therefore, the value of for the bulk states can be almost the same as for the surface states. Also, the electron mass for the bulk states Piot et al. (2016) is close to that for the surface states Nomura et al. (2014). Thus it is reasonable to expect that for the bulk states would of the same order as for the surface states.

In Table 1 we put estimated values of the dimensionless parameters for BiSe. These values were obtained by fitting the ARPES data presented in Refs. Nomura et al., 2014 and Rakyta et al., 2012. Here subscripts and stands for the bulk and surface states, respectively. In general, the positions of the Dirac cone for the surface states is different from the position of zero for the bulk states (as it was defined in Fig. 1). Thus, . Naturally, we can extract reliable values of parameters , , , and only for the surface states. We assume that the characteristics , , and are the same for the surface states and bulk states, while . We believe that such a choice does not affect the results within an order of magnitude.

We calculate the components of the spin conductivity for the set of parameters from Table 1 and for the dimensionless disorder strength estimated above. The results are presented in Table 2. We see that typically the isotropic, , and anisotropic, , components of the spin conductivity has the same order of magnitude. The contribution from the states at the Fermi level is comparable to the contribution from the filled states. The spin conductivity of the surface states is of the same order as the conductivity from the bulk states per layer.

BiSe | 0.1 | 0.1 | -0.5 | 0.07 | 0.7 |
---|

Bulk states per layer | 1.67 | 1.47 | 0.5 | 0.22 |

Surface states | 0.57 | 0.07 | 0.8 | 1.02 |

Using the estimated above values of parameters, we calculate the dependence of the total spin conductivity, , on the chemical potential. The results for the bulk and surface states are shown in Figs. 12. As we can see from the top panel in Fig. 12, the isotropic bulk spin conductivity slowly decreases from the plateau with an increase of , then, has a peak, and finally drops to zero. It is rather high in the helical state. As for the anisotropic component, , its value is considerably smaller than . It monotonously decreases to zero with the growth of the chemical potential. The spin conductivity due to the surface states does not depend crucially on the chemical potential, Fig. 12. The anisotropic component of the surface spin conductivity is much smaller than the isotropic one.

## Viii Discussion

The spin conductivity quanta can be expressed in the dimensional units as . In the previous section we estimate that the spin conductivity per conducting layer (either bulk or surface) is of the order of . The distance between the layers in the TIs is of the order of 1 nm. Then, the specific (volume) spin conductivity can be estimated as m. These values are close to that measured in BiSe Han et al. (2017).

In the present study, however, we can not explain a colossal spin conductivity in BiSb Huynh Duy Khang et al. (2017) and BiSbTe Fan et al. (2014), which is about two orders of magnitude higher than our estimation. This discrepancy may occur for the following reasons. First, it has been argued in Ref. Han et al., 2017 that different values of the spin conductivity can be a result of different fitting procedures of the experimental data. Second, in our consideration we do not take into account the effects of a magnetic field. The external magnetization (which is typically presented in the experiments Fan et al. (2014)) could drastically enhance the spin conductivity.

In Refs. Inoue et al. (2004); Raimondi and Schwab (2005) it has been shown that the vertex corrections in Rashba spin-orbit materials are small and, consequently, the (bulk in the TIs) spin conductivity from the states at the Fermi level is damped. However, in these papers the helical state, , has not been considered. According to our analysis, in such phase, the vertex correction to the Rashba spin-orbit coupling is comparable to its bare value . So, a quite significant contribution to the spin current can be observed even from the bulk states. The vertex correction also increases the contribution from the surface states to the spin conductivity.

It has been speculated that the surface states can generate large spin currents observed in the experiment Shiomi et al. (2014); Wang et al. (2015, 2016). According to our study, the surface states cannot produce very large spin current and the spin conductivity of the surface states typically has the same value as the spin conductivity of the bulk states per layer. So, our work confirms the experiments that show that spin conductivity mainly arises from the bulk states for a multilayer TI Jamali et al. (2015); Han et al. (2017) Also, it can explain the experiment, where the spin conductivity is small when the bulk of the TI sample is insulating, and the spin conductivity increases, when the bulk is conducting Kondou et al. (2016).

The bulk spin conductivity is robust against disorder if the spin-orbit coupling is large in comparison with a disorder, . Otherwise, these contributions to the spin conductivity is suppressed, see Figs. 8 and 11. The surface spin conductivity is robust against disorder even if a disorder is not weak, . The nature of robustness of the surface spin conductivity is similar to the robustness of the surface charge conductivity against disorder and arises due to suppression of the back-scattering. However, the study of the spin conductivity of the surface states in case of strong disorder deserves future studies.

The surface spin conductivity in a thin layer of TI with a cubic lattice have been studied in Ref. Peng et al., 2016 without taking into account the vertex corrections. The authors of the latter paper argue that the dependence of the surface spin conductivity on the disorder and chemical potential is weak. Our analysis confirms these results, see Figs. 8, 11, and 12.

The bulk spin conductivity can be tuned by a changing the chemical potential. Adjusting the chemical potential to the vicinity of the bottom of electron band, , we can attain the largest spin currents, see Fig. 12. The tuning of the spin conductivity contribution from the filled states by changing the chemical potential has been demonstrated numerically in Ref. Sahin and Flatte, 2015. The dependence of the spin-hall angle on the chemical potential has been measured in the experiments Fan et al. (2016). However, the observed result can be explained not only by the present analysis but also by the particle-hole asymmetry of the charge current.

In the experiment, the components of the spin conductivity tensor are measure in the coordinate axes related to the current. In Ref. Yang et al., 2016 the spin conductivity component was measured in BiSbTeSe. The authors observed different sign of this value in different samples. According to the conductivity tensor transformation presented in Section V, Eq. (V), the measured conductivity should be anisotropic and depends on the angle between the current and crystallographic -axis as . Thus, different sign of the measured spin conductivity may be due to the different orientation of the current leads with respect to the crystallographic axes in different samples.

## Acknowledgements

We acknowledge support from the Russian Scientific Foundation, Grant No 17-12-01544. RSA acknowledge the partial support by the Basis Foundation and ICFPM (MMK) of Education and Science of the Russian Federation, Grant No. 14Y26.31.0007.

## References

- Hasan and Kane (2010) M. Z. Hasan and C. L. Kane, “Colloquium: Topological insulators,” Rev. Mod. Phys. 82, 3045 (2010).
- Fu (2009) L. Fu, “Hexagonal Warping Effects in the Surface States of the Topological Insulator ,” Phys. Rev. Lett. 103, 266801 (2009).
- Kuroda et al. (2010) K. Kuroda, M. Arita, K. Miyamoto, M. Ye, J. Jiang, A. Kimura, E. E. Krasovskii, E. V. Chulkov, H. Iwasawa, T. Okuda, et al., “Hexagonally Deformed Fermi Surface of the 3D Topological Insulator ,” Phys. Rev. Lett. 105, 076802 (2010).
- Wang and Yu (2011) C. M. Wang and F. J. Yu, “Effects of hexagonal warping on surface transport in topological insulators,” Phys. Rev. B 84, 155440 (2011).
- Pal et al. (2012) H. K. Pal, V. I. Yudson, and D. L. Maslov, “Effect of electron-electron interaction on surface transport in the BiTe family of three-dimensional topological insulators,” Phys. Rev. B 85, 085439 (2012).
- Siu et al. (2014) Z. B. Siu, M. B. A. Jalil, and S. G. Tan, “Topological state transport in topological insulators under the influence of hexagonal warping and exchange coupling to in-plane magnetizations,” Scientific Reports 4, 5062 (2014).
- Repin and Burmistrov (2015) E. V. Repin and I. S. Burmistrov, “Surface states in a 3D topological insulator: The role of hexagonal warping and curvature,” Journal of Experimental and Theoretical Physics 121, 509 (2015).
- Akzyanov and Rakhmanov (2018) R. S. Akzyanov and A. L. Rakhmanov, “Surface charge conductivity of a topological insulator in a magnetic field: The effect of hexagonal warping,” Phys. Rev. B 97, 075421 (2018).
- Sinova et al. (2015) J. Sinova, S. O. Valenzuela, J. Wunderlich, C. Back, and T. Jungwirth, “Spin Hall effects,” Rev. Mod. Phys. 87, 1213 (2015).
- Murakami et al. (2003) S. Murakami, N. Nagaosa, and S.-C. Zhang, “Dissipationless Quantum Spin Current at Room Temperature,” Science 301, 1348 (2003).
- Sinova et al. (2004) J. Sinova, D. Culcer, Q. Niu, N. A. Sinitsyn, T. Jungwirth, and A. H. MacDonald, “Universal Intrinsic Spin Hall Effect,” Phys. Rev. Lett. 92, 126603 (2004).
- Inoue et al. (2004) J.-i. Inoue, G. E. W. Bauer, and L. W. Molenkamp, “Suppression of the persistent spin Hall current by defect scattering,” Phys. Rev. B 70, 041303 (2004).
- Raimondi and Schwab (2005) R. Raimondi and P. Schwab, “Spin-Hall effect in a disordered two-dimensional electron system,” Phys. Rev. B 71, 033311 (2005).
- Ralph and Stiles (2008) D. Ralph and M. Stiles, “Spin transfer torques,” Journal of Magnetism and Magnetic Materials 320, 1190 (2008).
- Chen et al. (2015) W. Chen, M. Sigrist, J. Sinova, and D. Manske, “Minimal Model of Spin-Transfer Torque and Spin Pumping Caused by the Spin Hall Effect,” Phys. Rev. Lett. 115, 217203 (2015).
- Wang et al. (2013) K. L. Wang, J. G. Alzate, and P. K. Amiri, “Low-power non-volatile spintronic memory: STT-RAM and beyond,” Journal of Physics D: Applied Physics 46, 074003 (2013).
- Mellnik et al. (2014) A. R. Mellnik, J. S. Lee, A. Richardella, J. L. Grab, P. J. Mintun, M. H. Fischer, A. Vaezi, A. Manchon, E.-A. Kim, N. Samarth, et al., “Spin-transfer torque generated by a topological insulator,” Nature 511, 449 (2014).
- Fan et al. (2014) Y. Fan, P. Upadhyaya, X. Kou, M. Lang, S. Takei, Z. Wang, J. Tang, L. He, L.-T. Chang, M. Montazeri, et al., “Magnetization switching through giant spinâorbit torque in a magnetically doped topological insulator heterostructure,” Nat Mater 13, 699 (2014).
- Fan et al. (2016) Y. Fan, X. Kou, P. Upadhyaya, Q. Shao, L. Pan, M. Lang, X. Che, J. Tang, M. Montazeri, K. Murata, et al., “Electric-field control of spinâorbit torque in a magnetically doped topological insulator,” Nat Nano 11, 352 (2016).
- Han et al. (2017) J. Han, A. Richardella, S. A. Siddiqui, J. Finley, N. Samarth, and L. Liu, “Room-Temperature Spin-Orbit Torque Switching Induced by a Topological Insulator,” Phys. Rev. Lett. 119, 077702 (2017).
- Yang et al. (2016) F. Yang, S. Ghatak, A. A. Taskin, K. Segawa, Y. Ando, M. Shiraishi, Y. Kanai, K. Matsumoto, A. Rosch, and Y. Ando, “Switching of charge-current-induced spin polarization in the topological insulator ,” Phys. Rev. B 94, 075304 (2016).
- Kondou et al. (2016) K. Kondou, R. Yoshimi, A. Tsukazaki, Y. Fukuma, J. Matsuno, K. S. Takahashi, M. Kawasaki, Y. Tokura, and Y. Otani, “Fermi-level-dependent charge-to-spin current conversion by Dirac surface states of topological insulators,” Nat Phys 12, 1027 (2016).
- Shiomi et al. (2014) Y. Shiomi, K. Nomura, Y. Kajiwara, K. Eto, M. Novak, K. Segawa, Y. Ando, and E. Saitoh, “Spin-Electricity Conversion Induced by Spin Injection into Topological Insulators,” Phys. Rev. Lett. 113, 196601 (2014).
- Wang et al. (2015) Y. Wang, P. Deorani, K. Banerjee, N. Koirala, M. Brahlek, S. Oh, and H. Yang, “Topological Surface States Originated Spin-Orbit Torques in ,” Phys. Rev. Lett. 114, 257202 (2015).
- Wang et al. (2016) H. Wang, J. Kally, J. S. Lee, T. Liu, H. Chang, D. R. Hickey, K. A. Mkhoyan, M. Wu, A. Richardella, and N. Samarth, “Surface-State-Dominated Spin-Charge Current Conversion in Topological-Insulator–Ferromagnetic-Insulator Heterostructures,” Phys. Rev. Lett. 117, 076601 (2016).
- Jamali et al. (2015) M. Jamali, J. S. Lee, J. S. Jeong, F. Mahfouzi, Y. Lv, Z. Zhao, B. K. NikoliÄ, K. A. Mkhoyan, N. Samarth, and J.-P. Wang, “Giant Spin Pumping and Inverse Spin Hall Effect in the Presence of Surface and Bulk SpinâOrbit Coupling of Topological Insulator Bi2Se3,” Nano Lett. 15, 7126 (2015).
- Yang and Chang (2006) M.-F. Yang and M.-C. Chang, “Středa-like formula in the spin Hall effect,” Phys. Rev. B 73, 073304 (2006).
- Kodderitzsch et al. (2015) D. Kodderitzsch, K. Chadova, and H. Ebert, “Linear response Kubo-Bastin formalism with application to the anomalous and spin Hall effects: A first-principles approach,” Phys. Rev. B 92, 184415 (2015).
- Sahin and Flatte (2015) C. Sahin and M. E. Flatte, “Tunable Giant Spin Hall Conductivities in a Strong Spin-Orbit Semimetal: ,” Phys. Rev. Lett. 114, 107201 (2015).
- Peng et al. (2016) X. Peng, Y. Yang, R. R. Singh, S. Y. Savrasov, and D. Yu, “Spin generation via bulk spin current in three-dimensional topological insulators,” Nature Communications 7, 10878 (2016).
- Sinitsyn et al. (2006) N. A. Sinitsyn, J. E. Hill, H. Min, J. Sinova, and A. H. MacDonald, “Charge and Spin Hall Conductivity in Metallic Graphene,” Phys. Rev. Lett. 97, 106804 (2006).
- Liu et al. (2015) Z. Liu, M. Zhu, and Y. Zheng, “Quantum transport properties of graphene in the presence of randomly distributed spin-orbit coupling impurities,” Phys. Rev. B 92, 245438 (2015).
- Garcia and Rappoport (2016) J. H. Garcia and T. G. Rappoport, “Kuboâ-Bastin approach for the spin Hall conductivity of decorated graphene,” 2D Materials 3, 024007 (2016).
- Sun et al. (2016) Y. Sun, Y. Zhang, C. Felser, and B. Yan, “Strong Intrinsic Spin Hall Effect in the TaAs Family of Weyl Semimetals,” Phys. Rev. Lett. 117, 146403 (2016).
- Liu et al. (2010) C.-X. Liu, X.-L. Qi, H. Zhang, X. Dai, Z. Fang, and S.-C. Zhang, “Model Hamiltonian for topological insulators,” Phys. Rev. B 82, 045122 (2010).
- Sinitsyn et al. (2004) N. A. Sinitsyn, E. M. Hankiewicz, W. Teizer, and J. Sinova, “Spin Hall and spin-diagonal conductivity in the presence of Rashba and Dresselhaus spin-orbit coupling,” Phys. Rev. B 70, 081312 (2004).
- Shon and Ando (1998) N. Shon and T. Ando, “Quantum Transport in Two-Dimensional Graphite System,” J. Phys. Soc. Jpn. 67, 2421 (1998).
- Kato et al. (2007) T. Kato, Y. Ishikawa, H. Itoh, and J. ichiro Inoue, “Anomalous Hall effect in spin-polarized two-dimensional electron gases with Rashba spinâorbit interaction,” New Journal of Physics 9, 350 (2007).
- Hu et al. (2008) B. Y.-K. Hu, E. H. Hwang, and S. Das Sarma, “Density of states of disordered graphene,” Phys. Rev. B 78, 165411 (2008).
- Chiba et al. (2017) T. Chiba, S. Takahashi, and G. E. W. Bauer, “Magnetic-proximity-induced magnetoresistance on topological insulators,” Phys. Rev. B 95, 094428 (2017).
- Adroguer et al. (2012) P. Adroguer, D. Carpentier, J. Cayssol, and E. Orignac, “Diffusion at the surface of topological insulators,” New Journal of Physics 14, 103027 (2012).
- Chen et al. (2013) C. Chen, Z. Xie, Y. Feng, H. Yi, A. Liang, S. He, D. Mou, J. He, Y. Peng, X. Liu, et al., “Tunable Dirac Fermion Dynamics in Topological Insulators,” Scientific Reports 3, 2411 (2013).
- Cheng et al. (2010) P. Cheng, C. Song, T. Zhang, Y. Zhang, Y. Wang, J.-F. Jia, J. Wang, Y. Wang, B.-F. Zhu, X. Chen, et al., “Landau Quantization of Topological Surface States in ,” Phys. Rev. Lett. 105, 076801 (2010).
- Zhang et al. (2009) H. Zhang, C.-X. Liu, X.-L. Qi, X. Dai, Z. Fang, and S.-C. Zhang, “Topological insulators in Bi2Se3, Bi2Te3 and Sb2Te3 with a single Dirac cone on the surface,” Nat Phys 5, 438 (2009).
- King et al. (2011) P. D. C. King, R. C. Hatch, M. Bianchi, R. Ovsyannikov, C. Lupulescu, G. Landolt, B. Slomski, J. H. Dil, D. Guan, J. L. Mi, et al., “Large Tunable Rashba Spin Splitting of a Two-Dimensional Electron Gas in ,” Phys. Rev. Lett. 107, 096802 (2011).
- Piot et al. (2016) B. A. Piot, W. Desrat, D. K. Maude, M. Orlita, M. Potemski, G. Martinez, and Y. S. Hor, “Hole Fermi surface in probed by quantum oscillations,” Phys. Rev. B 93, 155206 (2016).
- Nomura et al. (2014) M. Nomura, S. Souma, A. Takayama, T. Sato, T. Takahashi, K. Eto, K. Segawa, and Y. Ando, “Relationship between Fermi surface warping and out-of-plane spin polarization in topological insulators: A view from spin- and angle-resolved photoemission,” Phys. Rev. B 89, 045134 (2014).
- Rakyta et al. (2012) P. Rakyta, A. PÃ¡lyi, and J. Cserti, “Electronic standing waves on the surface of the topological insulator BiTe,” Phys. Rev. B 86, 085456 (2012).
- Huynh Duy Khang et al. (2017) N. Huynh Duy Khang, Y. Ueda, and P. N. Hai, “A conductive topological insulator with colossal spin Hall effect for ultra-low power spin-orbit-torque switching,” ArXiv e-prints (2017), eprint 1709.07684.