Skyrmion-induced anomalous Hall conductivity on topological insulator surfaces

TL;DR

This paper theoretically explores how magnetic skyrmions on topological insulator surfaces induce an additional anomalous Hall conductivity, highlighting a unique phase effect in Dirac electrons distinct from conventional topological Hall phenomena.

Contribution

It introduces a novel mechanism for skyrmion-induced anomalous Hall conductivity on topological insulators, emphasizing the role of Dirac electron phase factors at skyrmion boundaries.

Findings

Skyrmions contribute to anomalous Hall conductivity without in-plane magnetic textures.

The phase factor acquired by Dirac electrons at skyrmion boundaries is crucial.

Distinct from conventional topological Hall effects in metals.

Abstract

Electron spin-momentum locking together with background magnetic textures can significantly alter the electron transport properties. We theoretically investigate the electron transport at the interface between a topological insulator and a magnetic insulator with magnetic skyrmions on the top. In contrast to the conventional topological Hall effect in normal metals, the skyrmions yield an additional contribution to the anomalous Hall conductivity even in the absence of in-plane magnetic texture, arising from the phase factor characteristic to Dirac electrons acquired at skyrmion boundary.

Equations

$$\mathcal{H}=v_{\mathrm{F}}(\hat{\mathbf{p}}\times \mathbf{\sigma }{)}_z-\Delta \mathbf{n}(\mathbf{r})\cdot \mathbf{\sigma },$$

$$\mathbf{n}(\mathbf{r})=\left(\begin{array}{c}-\sqrt{1-n_z^2(\rho )}\sin \varphi \\ \sqrt{1-n_z^2(\rho )}\cos \varphi \\ n_z(\rho )\end{array}\right),$$

$${\psi }_j(\mathbf{r})=\left(\begin{array}{c}{u}_j(\rho )e^{i(j-\frac{1}{2})\varphi }\\ v_j(\rho )e^{i(j+\frac{1}{2})\varphi }\end{array}\right).$$

$$\displaystyle\psi_{k_{\mathrm{F}},j}(\rho,\phi)\sim\frac{1}{\sqrt{\rho}}\left(\begin{array}[]{c}\sin\tfrac{\zeta}{2}\cos(\xi_{j}+\frac{\pi}{4}-\delta_{j})e^{i(j-\frac{1}{2})\phi}\\ -\cos\tfrac{\zeta}{2}\sin(\xi_{j}+\frac{\pi}{4}-\delta_{j})e^{i(j+\frac{1}{2})\phi}\end{array}\right),$$

$$\cot {\delta }_j=-\frac{1}{2d}\left[(1+d)T_{j+\frac{1}{2}}({\xi }_{\mathrm{F}})-(1-d)T_{j-\frac{1}{2}}({\xi }_{\mathrm{F}})\right]$$

$$\tan {\delta }_j=-\frac{2\pi d}{d+\mathrm{sgn}j}\frac{\tilde{j}}{(\tilde{j}!{)}^2}{\left(\frac{{\xi }_{\mathrm{F}}}{2}\right)}^{2\tilde{j}}+O\left({\xi }_{\mathrm{F}}^{2\tilde{j}+2}\right),$$

$$\tan {\delta }_j=-\frac{d\cos (2{\xi }_{\mathrm{F}})}{(-1{)}^{j-\frac{1}{2}}+d\sin (2{\xi }_{\mathrm{F}})}+O({\xi }_{\mathrm{F}}^{-1}),$$

$$\displaystyle\Psi(\boldsymbol{r})\sim e^{ik_{\mathrm{F}}x}\left(\begin{array}[]{c}-i\sin\frac{\zeta}{2}\\ \cos\frac{\zeta}{2}\end{array}\right)+\frac{e^{ik_{\mathrm{F}}\rho}}{\sqrt{\rho}}\left(\begin{array}[]{c}f_{\uparrow}(\phi)\\ f_{\downarrow}(\phi)\end{array}\right),$$

$$\displaystyle\left(\begin{array}[]{c}f_{\uparrow}(\phi)\\ f_{\downarrow}(\phi)\end{array}\right)=\frac{1}{\sqrt{2\pi k_{\mathrm{F}}}}e^{i\frac{\pi}{4}}\Phi(\phi)\left(\begin{array}[]{c}\sin\frac{\zeta}{2}e^{-i\phi/2}\\ -\cos\frac{\zeta}{2}e^{i\phi/2}\end{array}\right),$$

$$\Phi (\varphi )=\sum _j{e}^{ij\varphi }(1-e^{-2i{\delta }_j}),$$

$$\frac{d\sigma (\varphi )}{d\varphi }=\frac{j^{(\mathrm{out})}\rho }{j^{(\mathrm{in})}}=\frac{2E_{\mathrm{F}}}{v_{\mathrm{F}}{k}_{\mathrm{F}}}\mathrm{Im}\left[e^{-i\varphi }{f}_{\uparrow }^{\ast }{f}_{\downarrow }\right]=\frac{|\Phi (\varphi ){|}^2}{2\pi k_{\mathrm{F}}}.$$

$$\left(\begin{array}{c}F\\ F_{\parallel }\\ F_{\perp }\end{array}\right)$$

$$F\simeq 2\sum _{\pm }{\sin }^2{\delta }_{\pm 1/2}=\frac{16{\pi }^2{d}^2}{1-d^2}{\xi }_{\mathrm{F}}^4+O({\xi }_{\mathrm{F}}^6).$$

$$F\simeq 2\sum _{|j|<{\xi }_{\mathrm{F}}}{\sin }^2{\delta }_j=\frac{4{\xi }_{\mathrm{F}}{d}^2(1+d^2){\cos }^2(2{\xi }_{\mathrm{F}})}{(1+d^2{)}^2-4d^2{\sin }^2(2{\xi }_{\mathrm{F}})}+O({\xi }_{\mathrm{F}}^0).$$

$$F_{\parallel },|F_{\perp }|\lesssim 4{\xi }_{\mathrm{F}}{d}^2+O({\xi }_F^0,d^4).$$

$$-e\mathbf{E}\cdot {\nabla }_{\mathbf{k}}f(\mathbf{k})={\left(\frac{df(\mathbf{k})}{dt}\right)}_{\mathrm{coll}},$$

$${\left(\frac{df(\mathbf{k})}{dt}\right)}_{\mathrm{coll}}=-\frac{1}{{\tau }_0}\left[f(\mathbf{k})-f_0(\mathbf{k})\right]$$

$$-n_{\mathrm{S}}\int d\varphi |\mathbf{v}(\mathbf{k})|\left[\frac{d\sigma (\varphi )}{d\varphi }f(\mathbf{k})-\frac{d\sigma (-\varphi )}{d\varphi }f({\mathbf{k}}^{\prime })\right],$$

$$\delta f(\mathbf{k})$$

$$\left(\begin{array}{c}{C}_k\\ D_k\end{array}\right)=\frac{v_{\mathrm{F}}^2}{E_{\mathrm{F}}}\left(\begin{array}{c}{\tau }_0-{\tau }_0^2{\tau }_{\parallel }^{-1}\\ -{\tau }_0^2{\tau }_{\perp }^{-1}\end{array}\right).$$

$${\tau }_X^{-1}=n_{\mathrm{S}}\frac{v_{\mathrm{F}}^2}{E_{\mathrm{F}}}{F}_X.(X=\parallel ,\perp )$$

$$\left(\begin{array}{c}{\tilde{\sigma }}_{xx}\\ {\tilde{\sigma }}_{xy}\end{array}\right)=\frac{e^2}{4\pi }\frac{E_{\mathrm{F}}^2-{\Delta }^2}{E_{\mathrm{F}}}\left(\begin{array}{c}{\tau }_0-{\tau }_0^2{\tau }_{\parallel }^{-1}\\ {\tau }_0^2{\tau }_{\perp }^{-1}\end{array}\right),$$

$${\tau }_X=\frac{E_{\mathrm{F}}}{n_{\mathrm{S}}{v}_{\mathrm{F}}^2{F}_X}\sim \frac{10\mathrm{ps}}{F_X}.$$

$$\displaystyle\psi_{j}(\rho,\phi)=\left(\begin{array}[]{c}u_{j}(\rho)e^{il_{\uparrow}\phi}\\ v_{j}(\rho)e^{il_{\downarrow}\phi}\end{array}\right),$$

$$\displaystyle\left(\begin{array}[]{cc}-\Delta n_{z}(\rho)-E&v_{\mathrm{F}}(-\partial_{\rho}-l_{\downarrow}\rho^{-1})\\ v_{\mathrm{F}}(\partial_{\rho}-l_{\downarrow}\rho^{-1})&\Delta n_{z}(\rho)-E\end{array}\right)\left(\begin{array}[]{c}u_{j}(\rho)\\ v_{j}(\rho)\end{array}\right)=0,$$

$$\displaystyle\left(\begin{array}[]{cc}-\Delta-E&v_{\mathrm{F}}(-\partial_{\rho}-l_{\downarrow}\rho^{-1})\\ v_{\mathrm{F}}(\partial_{\rho}-l_{\downarrow}\rho^{-1})&\Delta-E\end{array}\right)\left(\begin{array}[]{c}u_{j}(\rho)\\ v_{j}(\rho)\end{array}\right)=0.$$

$$v_j(\rho )$$

$$\left[{\partial }_{\rho }^2+\frac{1}{\rho }{\partial }_{\rho }+\frac{E^2-{\Delta }^2}{v_{\mathrm{F}}^2}-\frac{l_{\uparrow }^2}{{\rho }^2}\right]u_j(\rho )=0.$$

$$\left[{\partial }_{\xi }^2-\frac{1}{\xi }{\partial }_{\xi }+1-\frac{l_{\uparrow }^2}{{\xi }^2}\right]u_j(\xi )=0,$$

$${\chi }_{k,j}^{(+)}(\rho ,\varphi )$$

$$\left({\partial }_{\xi }\mp \frac{l}{\xi }\right)J_l(\xi )=\mp J_{l\pm 1}(\xi ).$$

Full text

I Introduction

Physics of magnetic textures has been a significant topic in recent studies on magnetic materials, to make use of them as efficient carriers of information Zutic_Fabian_DasSarma ; Brataas_Bauer_Kelly ; Tatara_Kohno_Shibata . Skyrmion, in particular, is a particlelike magnetic excitation with a swirling texture in two dimensions (2D), in which the spin at the core and those at the perimeter point in the opposite directions Nagaosa_Tokura ; Fert . Such non-collinear magnetic textures are observed in non-centrosymmetric magnets, such as helimagnetic conductor $\mathrm{MnSi}$ (Refs. Muhlbauer, ; Tonomura_Tokura, ) and magnetic insulator (MI) $\mathrm{Cu}_{2}\mathrm{O}\mathrm{Se}\mathrm{O}_{3}$ (Refs. Seki_Tokura, ; Adams, ). Skyrmions form a periodic lattice, called skyrmion crystal, with the lattice constant around $10$-$100\mathrm{nm}$. One of the most important features of skyrmions is the idea of emergent electromagnetic fields for conduction electrons, arising from nonzero Berry curvature in the real space Nagaosa_Tokura_PhysScr ; Nagaosa_Yu_Tokura . A conduction electron traveling across the magnetic texture feels the emergent magnetic field, leading to the so-called topological Hall effect (THE) Ye_THE ; Tatara_Kawamura ; Lee_THE ; Neubauer_THE .

In the presence of strong spin-orbit coupling (SOC), on the other hand, the band topology in the momentum space significantly affects the electron transport as well TKNN . A typical example is the surface state of a topological insulator (TI), where the electrons show 2D linear dispersion with a single band-touching (Dirac) point Hasan_Kane ; Qi_Zhang . Under a finite out-of-plane magnetization, the Dirac point is gapped out and the system acquires a finite Chern number $C=\pm 1$ from the momentum-space Berry curvature, leading to the well-known quantum anomalous Hall effect (quantum AHE, QAHE) Qi_Wu_Zhang ; Qi_Hughes_Zhang ; Yu ; Nomura_QAHE ; Checkelsky ; Chang . Since the exchange coupling between the local in-plane magnetization and the TI surface electron can be regarded as the emergent gauge field (vector potential) for the conduction electrons Garate_Franz , topological magnetic textures on TI surfaces can give rise to even richer electronic properties: vortices and domain walls, for instance, are supposed to host zero modes, leading to electric charging of those textures Nomura_DW ; Tserkovnyak_Loss ; Tserkovnyak_Pesin_Loss . Recent theoretical studies predict that magnetic skyrmions on TI surface can be charged as well Hurst_Galitski ; Andrikopoulos . Transport measurements in magnetic TI heterostructures discovered the coexistence of the THE related to the real-space Berry curvature from skyrmions and the AHE related to momentum-space counterpart from SOC Yasuda , while the transport is dominated by the bulk states, leaving the surface transport ambiguous.

In this work, we theoretically examine the electron transport on TI surface in the presence of skyrmions, which possibly accounts for the interface between a TI and a non-centrosymmetric MI. Due to electron spin-momentum locking feature on TI surfaces, the electron spin does not follow the local magnetic texture even under a strong exchange splitting, which makes the conventional THE scenario arising from the real-space Berry curvature, based on the adiabatic approximation, unreliable. Hence we fully solve the electron scattering problem by a single skyrmion, and apply its result to the Boltzmann transport theory to estimate the longitudinal and Hall conductivities. We find that the skyrmions give a sizable additive contribution to the AHE, due to the skewness in the electron scattering at a skyrmion. The origin of this skewness is the phase factor acquired at the skyrmion boundary, which is characteristic to Dirac electrons and is absent in Schrödinger electron systems. It arises even in the absence of in-plane magnetic texture, which is totally different from the conventional THE in normal metals. In this paper, we take $\hbar=1$ and restore it in the final numerical results.

II Single-skyrmion problem

Let us set up a heterostructure of a TI (e.g. $\mathrm{Bi}_{2}\mathrm{Se}_{3}$) and a MI (e.g. $\mathrm{Cu}_{2}\mathrm{O}\mathrm{Se}\mathrm{O}_{3}$), with its interface taken parallel to the layers of the TI so that the coupling between the TI and the MI should be homogeneous. We first place a single skyrmion at the center of the 2D infinite space, to examine how an incoming electron plane wave gets scattered by the skyrmion. The electron at the TI-MI interface under the magnetic texture $\boldsymbol{n}(\bm{r})$ is described by the Hamiltonian

[TABLE]

where the coefficient $v_{\mathrm{F}}$ is the Fermi velocity, $\hat{\bm{p}}=-i\boldsymbol{\nabla}$ is the momentum operator, and $\Delta$ is the spin splitting energy from the exchange interaction between the electron and the magnetization. Since the interface is parallel to the layers of the TI, the spin of the interface electron is helical in the momentum space, i.e. spin-momentum locked Silvestrov ; Zhang_Kane_Mele . Although the TI surface Hamiltonian in realistic systems exhibit the terms beyond the linear order in the momentum $\hat{\boldsymbol{p}}$, their effect on the electron scattering process that we are interested in is just a slight modulation of electron spin-momentum locking. It may shift the scattering amplitude quantitatively, but would not qualitatively alter the angular profile of the scattering behavior (differential scattering cross section). Thus we will neglect those higher order terms in our present calculation. In the numerical results below, we fix the parameters $v_{\mathrm{F}}=0.5\times 10^{6}\;\mathrm{m/s}$ and $\Delta=10\;\mathrm{meV}$, which are typical values in TI-MI heterostructures, such as $\mathrm{Bi_{2}Se}_{3}$ and yttrium iron garnet (YIG) Liu_Wang ; Kubota .

We here fix the skyrmion texture unaffected by the interface electrons, assuming that it is stabilized by the Dzyaloshinskii-Moriya interaction (DMI) in the MI. Here a skyrmion with a vortical magnetic texture is energetically more stable than that with a hedgehog texture Nagaosa_Tokura ; Fert . A vortex skyrmion can be parametrized by the cylindrical coordinate $\bm{r}=(\rho,\phi)$ as

[TABLE]

where $n_{z}(\rho)$ is a scalar function taking the value between $-1$ at the center and $+1$ at the perimeter. Here we should note that the in-plane magnetic texture in this “vortex skyrmion” can be removed by the local U(1) gauge transformation $U(\rho)=\exp[i\frac{\Delta}{v_{\mathrm{F}}}\int_{0}^{\rho}d\rho^{\prime}\sqrt{1-n_{z}^{2}(\rho^{\prime})}]$, which states that only the sign flip in $n_{z}(\rho)$ is essential in our calculation. We thus employ the “hard-wall” approximation $n_{z}(\rho)=\mathrm{sgn}(\rho-R_{\mathrm{S}})$, with $R_{\mathrm{S}}$ the skyrmion radius Hurst_Galitski . Such a structure is likely to be realized under a strong out-of-plane magnetic anisotropy. It is not responsible for the THE in normal metals due to the absence of the real-space Berry curvature, while it still poses nontrivial effects to Dirac electrons.

It is inadequate to treat the skyrmion texture on the TI surface perturbatively, since it fails to incorporate the topological characteristics, namely the change in the Chern number inside the skyrmion. We therefore need to solve the scattering problem of a Dirac electron by a single skyrmion non-perturbatively Ferreira ; Ferreira_2014 . The scattering process is characterized by the “phase shifts” of the eigenstates, which are the angles quantifying how much the eigenstates are altered by the skyrmion R_G_Newton .

Since the Hamiltonian is invariant under the simultaneous rotation in the real and spin spaces around the $z$-axis, the total angular momentum $j=\pm 1/2,\pm 3/2,\cdots\in\mathbb{Z}+1/2$ serves as a good quantum number. The eigenstate wave function with the angular momentum $j$ can be parametrized as

[TABLE]

The radial equations for $u_{j}(\rho)$ and $v_{j}(\rho)$ can be analytically solved under the hard-wall approximation (see Appendix for detail). The radial solution behaves as a cylindrical wave, with the asymptotic wave number $k_{\mathrm{F}}=v_{\mathrm{F}}^{-1}\sqrt{E_{\mathrm{F}}^{2}-\Delta^{2}}$ away from the skyrmion center $(\rho\rightarrow\infty)$, for the conduction electron at the Fermi level $E_{\mathrm{F}}(>\Delta)$. Using the phase shift $\delta_{j}$, the asymptotic behavior can be expressed as

[TABLE]

where $\xi_{j}\equiv k_{\mathrm{F}}\rho-\frac{\pi}{2}(j+\frac{1}{2})$ and $\cos\zeta\equiv d\equiv\Delta/E_{\mathrm{F}}$. The skyrmion-free solution can be reproduced by setting $\delta_{j}=0$. The phase shift $\delta_{j}$ is determined to make the solution continuous at the skyrmion edge $R_{\mathrm{S}}$, given by

[TABLE]

for each angular momentum mode $j$, with $\xi_{\mathrm{F}}\equiv k_{\mathrm{F}}R_{\mathrm{S}}$. Here $T_{l}(\xi)\equiv Y_{l}(\xi)/J_{l}(\xi)$, with $J_{l}$ the Bessel function and $Y_{l}$ the Neumann function. As can be see from this relation, the electron scattering behavior by the skyrmion is characterized by two dimensionless parameters, $d\equiv\Delta/E_{\mathrm{F}}$ and $\xi_{\mathrm{F}}\equiv k_{\mathrm{F}}R_{\mathrm{S}}$.

The asymptotic behavior of the phase shift $\delta_{j}$ can be derived from that of the Bessel functions $J$ and $Y$, the details of which are shown in the Appendix. In the long-wavelength limit $\xi_{\mathrm{F}}\ll|j|$, $\delta_{j}$ is given as

[TABLE]

where $\tilde{j}\equiv|j|+1/2$. Therefore, for fixed $\xi_{\mathrm{F}}$, the scattering phase shift $\delta_{j}$ approaches zero as the angular momentum $|j|$ becomes larger than $\xi_{\mathrm{F}}$. It can be well understood in the semiclassical picture: if the impact parameter $|y|\sim|j|/k_{\mathrm{F}}$, namely the distance from the scattering center to the incident trajectory, is much larger than the skyrmion radius $R_{\mathrm{S}}$, the incident particle is almost unaffected by the skyrmion. In the short-wavelength limit $\xi_{\mathrm{F}}\gg|j|$, on the other hand, $\delta_{j}$ can be approximated as

[TABLE]

which implies that $\delta_{j}=\delta_{j+2}$ up to the difference of $O(\xi_{\mathrm{F}}^{-1})$.

If a plane wave with the incident wave vector $\boldsymbol{k}=(k_{\mathrm{F}},0)$ is scattered by the skyrmion placed at the origin, as shown in Fig. 1, we can set up an ansatz for the total wave function away from the skyrmion $(\rho\rightarrow\infty)$ as

[TABLE]

where the first term comes from the incoming plane wave in the $x$-direction. The second term corresponds to the outgoing scattered waves, which we denote $\psi^{\mathrm{(out)}}(\rho,\phi)$ here. Since the total wave function on the right hand side can be given as a linear combination of the eigenfunctions $\{\psi_{k_{\mathrm{F}},j}\}$, we can estimate the scattering amplitudes $f_{\uparrow/\downarrow}(\phi)$ by the partial wave decomposition: by decomposing both sides of Eq. (14) and comparing the radial parts for each partial wave component, we obtain the scattering amplitudes

[TABLE]

where the angular function $\Phi(\phi)$ is composed of the phase shifts as

[TABLE]

with $j\in\mathbb{Z}+1/2$. As shown in Eq. (8), since the phase shift approaches zero for sufficiently large $j$, namely $|j|\gg k_{\mathrm{F}}R_{\mathrm{S}}$, we shall limit the partial waves to $|j|<2\xi_{\mathrm{F}}+1$ in the calculations below.

Using the outgoing wave function $\psi^{\mathrm{(out)}}$ obtained above, the differential cross section, namely the ratio between the incoming current density $j^{\mathrm{(in)}}$ and the outgoing flux $j^{\mathrm{(out)}}(\phi)\rho d\phi$ toward the angle $\phi$, is given by

[TABLE]

The angular profile of the differential cross section $d\sigma(\phi)/d\phi$ is shown in Fig. 2, for several Fermi energies $E_{\mathrm{F}}$ and skyrmion radii $R_{\mathrm{S}}$. In contrast to normal scattering processes, such as the scattering of Schrödinger electrons by a symmetric Coulomb potential, the scattering process shows a skewness, asymmetric about $\phi=0$.

The strengths of the whole scattering, the back scattering, and the skew scattering, are measured by the integrals

[TABLE]

respectively, with $\bar{\delta}_{j}=\delta_{j+1}-\delta_{j}$. It is obvious that the back and skew scattering strengths do not exceed the whole scattering strength, i.e. $F_{\parallel},|F_{\perp}|<F$. If the skyrmion size is smaller than the electron wavelength ($\xi_{\mathrm{F}}\ll 1$), the asymptotic behavior for $\delta_{j}$ in Eq. (8) applies for all the partial wave modes. The lowest angular momentum modes $(j=\pm 1/2)$ give the dominant contribution to $F$, given by

[TABLE]

On the other hand, if the skyrmion is large enough $(\xi_{\mathrm{F}}\gg 1)$, the modes with $|j|\ll\xi_{\mathrm{F}}$ yield the dominant contribution to $F$. Using the asymptotic behavior in Eq. (9), the leading-order contribution to the total scattering amplitude is given as

[TABLE]

If the exchange gap $\Delta$ is small enough compared with the electron energy $E_{\mathrm{F}}$, i.e. $d\ll 1$, it further reduces to the simplified behavior, $F=4\xi_{\mathrm{F}}d^{2}\cos^{2}(2\xi_{\mathrm{F}})+O(\xi_{F}^{0},d^{4})$. Both the backscattering strength $F_{\parallel}$ and the skew scattering strength $|F_{\perp}|$ are limited by this upper bound, yielding

[TABLE]

This asymptotic relation shall be used for the qualitative estimation of the transport properties in the next section.

The scattering skewness shown here can be traced back to the geometric phase acquired at the skyrmion boundary, which is characteristic to Dirac electrons. When an electron is transmitted through the boundary between two regions with opposite magnetizations, the scattered wave acquires a phase factor $e^{i\gamma(\theta)}$ relative to the incident wave, depending on its incident angle $\theta$. In the electron scattering process treated here, the incident angle $\theta$ at the skyrmion boundary depends on the impact parameter $y$, at sufficiently long wavelength such that the uncertainty in $y$ can be neglected (see Fig. 1). Therefore, the phase factor $e^{i\gamma(\theta)}$ depends on $y$; since this phase can be expanded as $\gamma(\theta(y))=\pi+\tilde{k}y+O(y^{2})$ around the scattering center, the scattered wave bears the transverse component $e^{i\tilde{k}y}$, leading to the scattering skewness (see Appendix for detail).

III Semiclassical transport analysis

The electron scattering by skyrmions alters the electronic transport properties, the conductivity in particular. Here we set up a uniformly magnetized 2D system with an ensemble of skyrmions, with the number density $n_{\mathrm{S}}$. As long as the skyrmion distribution is random and dilute enough, quantum interference from multiple scattering processes is ruled out. Here the interplay between skyrmion scattering and other scattering process, such as skew or side-jump scattering at normal impurities in the presence of SOC, is quite small, so that scattering by skyrmions contributes to the conductivity additively. Hence we estimate this additive contribution by the semiclassical Boltzmann theory.

The electron transport driven by an electric field $\boldsymbol{E}=E_{x}\boldsymbol{e}_{x}$ is described by the Boltzmann equation,

[TABLE]

for the steady-state electron distribution $f(\boldsymbol{k})$. The scattering process is incorporated in the collision term

[TABLE]

with the velocity $\boldsymbol{v}(\bm{k})=\boldsymbol{\nabla}_{\bm{k}}E(\bm{k})$, where we omit the intrinsic contribution from the $\bm{k}$-space Berry curvature. The second term describes the skyrmion-induced scattering process from momentum $\bm{k}$ into $\bm{k}^{\prime}$ and vice versa, with $\phi$ the angle between them. We also introduce the backscattering effect by normal impurities phenomenologically by the first term, in terms of relaxation-time approximation, with the transport relaxation time $\tau_{0}$ and the equilibrium distribution $f_{0}(\bm{k})$. The effect of skew and side-jump scattering processes by normal impurities are so far neglected, whose contributions to the conductivity shall be incorporated in our discussion later. The momenta $\bm{k}$ and $\bm{k}^{\prime}$ are limited to the Fermi surface, i.e. $|\bm{k}|=|\bm{k}^{\prime}|=k_{\mathrm{F}}$, since the scattering here is elastic.

The steady-state distribution $f(\bm{k})$ can be obtained up to the linear response in the electric field $\boldsymbol{E}$. Taking the ansatz Ferreira_2014

[TABLE]

for the deviation of the distribution function $\delta f(\bm{k})=f(\bm{k})-f_{0}(\bm{k})$, we can straightforwardly solve the Botlzmann equation at zero temperature. Since we are neglecting the interplay between skyrmion and other scattering processes, we evaluate the solution up to the first order in the skyrmion concentration $n_{\mathrm{S}}$ (we shall check the validity of this approximation in the numerical calculation later), yielding

[TABLE]

Here $\tau_{\parallel}$ and $\tau_{\perp}$ are the “scattering times” characterizing the backscattering and skew scattering effects from the randomly distributed skyrmions, respectively, given by

[TABLE]

Once we obtain the distribution shift $\delta f(\bm{k})$ driven by the electric field, we can estimate its extrinsic contribution to the current by $\tilde{\boldsymbol{j}}=-(e/V)\sum_{\bm{k}}\boldsymbol{v}(\bm{k})\delta f(\bm{k})$, leading to the conductivity coefficients

[TABLE]

up to the first order in $n_{\mathrm{S}}$. The longitudinal part $\tilde{\sigma}_{xx}$ consists of the backscattering effect by the normal scatterers, proportional to $\tau_{0}$, and that by the skyrmions, proportional to $\tau_{0}^{2}\tau_{\perp}^{-1}$.

On the other hand, as we have mentioned, $\tilde{\sigma}_{xy}$ obtained here is the additive contribution to the Hall conductivity induced by the skyrmions; other contributions, such as the ordinary Hall effect contribution from the external magnetic field, the intrinsic contribution from the momentum-space Berry curvature Qi_Wu_Zhang , and the extrinsic (skew and side-jump) contribution from the normal scatterers Sinitsyn ; Sinitsyn_review ; Culcer ; Lu_2013 ; Titov , come additively to the Hall conductivity as well, which shall be discussed at the end of this section.

Now we calculate the skyrmion contribution to the Hall conductivity $\tilde{\sigma}_{xy}$ as functions of Fermi energy $E_{\mathrm{F}}$ and the skyrmion radius $R_{\mathrm{S}}$, which are shown in Fig. 3. We take screened Coulomb impurities as the normal scatterers, which give the transport relaxation time $\tau_{0}$ proportional to $n_{\mathrm{c}}^{-1/2}$, with $n_{\mathrm{c}}=k_{\mathrm{F}}^{2}/4\pi$ the carrier density Nomura_MacDonald_2006 ; Nomura_MacDonald_2007 ; Ando_JPSJ ; Hwang_DasSarma . Here we employ $\tau_{0}$ shown in Ref. Hwang_DasSarma, , which was calculated for graphene deposited on $\mathrm{SiO}_{2}$ substrate, yielding $\tau_{0}=0.6\mathrm{ps}\times(k_{\mathrm{F}}/0.1\mathrm{nm}^{-1})$ under the impurity density $n_{\mathrm{i}}=10^{11}\mathrm{cm}^{-2}$. The skyrmion number density $n_{\mathrm{S}}$ is fixed to $2.89\times 10^{9}\mathrm{cm}^{-2}$, corresponding to the triangular skyrmion lattice with the spacing $200\mathrm{nm}$. With those parameters, we can estimate the skyrmion scattering time $\tau_{X}\ (X=\parallel,\perp)$ around $E_{\mathrm{F}}\sim 50\mathrm{meV}$ as

[TABLE]

At the skyrmion size $R_{\mathrm{S}}$ around $100\mathrm{nm}$, $F_{X}$ can be roughly estimated as $F_{X}\lesssim 4\xi_{\mathrm{F}}d^{2}\sim 2$, leading to $\tau_{X}\gtrsim 5\mathrm{ps}$. Therefore, the skyrmion scattering effect $\tau_{X}^{-1}$ is weaker than the normal impurity effect $\tau_{0}^{-1}$, so that the linear approximation by $n_{\mathrm{S}}$ (or $\tau_{X}^{-1}$) will be reliable in this calculation.

Our calculation results show that the contribution to the Hall conductivity is suppressed at small $E_{\mathrm{F}}$ around the band bottom. Such a behavior can be traced back to two reasons: since spin-momentum locking becomes less significant around the band bottom, the phase factor effect at the skyrmion boundary becomes weaker, suppressing the skew scattering effect. Moreover, the Coulomb impurities give a larger scattering rate at lower electron concentration, which makes the skyrmion effect relatively small. On the other hand, if the chemical potential is set far away from the bandgap, the band structure asymptotically reaches the gapless Dirac dispersion so that the electron eigenstates will become insensitive to the magnetic texture, suppressing the skyrmion effect. Thus we need a moderate chemical potential to obtain a sizable skyrmion contribution to the Hall conductivity, which is around $50\mathrm{meV}$ in our calculation.

We also find a quantum oscillation-like behavior in $\tilde{\sigma}_{xy}$ under the modulation of the skyrmion size $R_{\mathrm{S}}$ and the Fermi energy $E_{\mathrm{F}}$. We expect that this oscillation may originate from the electron resonance states formed inside the skyrmion, in which an electron standing wave should be formed within the skyrmion diameter. Such resonance states become less stable under a smooth magnetic texture without the hard-wall approximation, which may relax the oscillation behavior.

From the results of Ref. Sinitsyn, , we can make an order estimation of the normal impurity contribution to the Hall conductivity to compare it with the skyrmion contribution, which is beyond the scope of our Boltzmann analysis. If the impurity distribution is Gaussian, its effect on the Hall conductivity is insensitive to the impurity potential strength and scales by $\lesssim(e^{2}/4\pi)(\Delta/E_{\mathrm{F}})$ for $E_{\mathrm{F}}\gg\Delta$. If the distribution is beyond Gaussian, it acquires an additional term from the third-order moment of the impurity distribution, corresponding to the skew scattering effect by the impurities. It is known to be inversely proportional to the scattering strength, scaling as $\sim(e^{2}/4\pi)(\Delta/n_{\mathrm{i}}V_{\mathrm{i}})$, where $V_{\mathrm{i}}$ parametrizes the impurity potential in terms of delta-function potential. At high $E_{\mathrm{F}}$, the impurity potential $V_{\mathrm{i}}$ gets reduced by strong screening, leading to enhancement of this third-order moment contribution. For $E_{\mathrm{F}}\sim 50\mathrm{meV}$ and $n_{\mathrm{i}}\sim 10^{11}\mathrm{cm}^{-2}$, the denominator $n_{\mathrm{i}}V_{\mathrm{i}}$ is around $10\mathrm{meV}$, estimated from the screened Coulomb potential Hwang_DasSarma . Therefore, the skyrmion contribution to the anomalous Hall conductivity is comparable or even superior to the normal impurity contribution at $E_{\mathrm{F}}\sim 50\mathrm{meV}$, which is possibly captured in transport measurements.

IV Discussion

In this paper, we have investigated the electronic transport properties on TI surfaces in the presence of magnetic skyrmions. We have demonstrated that the skyrmions give rise to the extrinsic AHE, originating from the interplay between the magnetic texture and the Dirac band structure. Variation of the Fermi energy $E_{\mathrm{F}}$ (or the Fermi momentum $k_{\mathrm{F}}$) can be experimentally realized by tuning the gate voltage or the carrier doping, while the skyrmion size can also be modulated by external magnetic fields.

It should be noted that, in realistic magnetic materials, an external magnetic field is essential to stabilize skyrmions. In $\mathrm{Cu}_{2}\mathrm{O}\mathrm{Se}\mathrm{O}_{3}$, for instance, a magnetic field of $\sim 0.05\mathrm{T}$ is required, which corresponds to the magnetic length $\sim 100\mathrm{nm}$. The skyrmion size should not exceed this scale; otherwise the Landau quantization inside the skyrmion takes place and the scattering property obtained in this paper is no longer reliable.

We have so far relied on the hard-wall approximation on the vortex-type skyrmion texture. Such an approximation is justified under strong out-of-plane magnetic anisotropy, which stabilizes out-of-plane magnetization pattern. On the other hand, our results do not apply to the hedgehog skyrmion, since it cannot be gauged away due to the emergent magnetic flux. Such type of skyrmion is expected to be induced by frustrated exchange interactions or four-spin exchange interactions Nagaosa_Tokura , whose effect on the electron transport on TI surfaces remains an open question.

Acknowledgements.

The authors acknowledge T. Chiba, C. Lee, Y. Ominato, and K. Sato for fruitful discussions. Y. A. is supported by JSPS KAKENHI Grant Number JP17K14316. K. N. is supported by JSPS KAKENHI Grant Numbers JP15H05854 and JP17K05485.

Appendix A Eigenstate wave functions under a single skyrmion

In this part, we show the detailed derivation of the eigenstate wave functions under a single skyrmion, which leads to the phase shift [Eq. (3)]. For the partial wave component with the total angular momentum $j(=\pm 1/2,\pm 3/2,\cdots)$, we can set up the ansatz

[TABLE]

with $l_{\uparrow/\downarrow}\equiv j\mp 1/2$ the orbital angular momentum for each component. With this ansatz, the Dirac equation $\mathcal{H}\psi=E\psi$ reduces to the radial differential equation,

[TABLE]

under the hard-wall approximation on the skyrmion texture.

Let us solve the equation without skyrmion, namely with $n_{z}(\rho)=1$ in the entire system, as the starting point. Here the equation reads

[TABLE]

Using the relation

[TABLE]

we obtain the second-order differential equation for $u$,

[TABLE]

With the radial wavenumber $k\equiv v_{\mathrm{F}}^{-1}\sqrt{E^{2}-\Delta^{2}}$ and the change of variable $\xi\equiv k\rho$, Eq. (47) reduces to the well-known Bessel’s differential equation,

[TABLE]

whose linearly independent solutions are given by the Bessel function $J_{l_{\uparrow}}(\xi)$ and the Neumann function $Y_{l_{\uparrow}}(\xi)$. Since $Y_{l_{\uparrow}}$ is ruled out for solution due to the singularity at $\rho=0$, we can set $u_{j}(\rho)=C_{k,j}^{(+)}J_{l_{\uparrow}}(k\rho)$ with $C$ the normalization constant. Substituting this solution to Eq. (46), we obtain the eigenstate wave function

[TABLE]

where we have used the recurrence relation

[TABLE]

The normalization condition

[TABLE]

fixes the normalization constant as

[TABLE]

where we have used the orthogonality relation

[TABLE]

Therefore, the normalized eigenstate wave function is given by

[TABLE]

Here the angle $\zeta$ is defined by $E=\Delta/\cos\zeta$, so that it should satisfy

[TABLE]

We should note here the relation

[TABLE]

for integer $l$. If the magnetization is flipped, i.e. $n_{z}(\rho)=-1$, we should substitute $\zeta$ by $\pi-\zeta$, leading to the eigenfunction

[TABLE]

In the presence of a hard-wall skyrmion, i.e.

[TABLE]

we should solve Eq. (45) inside and outside the skyrmion separately, and connect the solution at the boundary. Inside the skyrmion, the solution is given by Eq. (63). Outside the skyrmion, on the other hand, the solutions is not simply given by Eq. (58), since we do not need to respect the normalizability at the origin. The Neumann counterpart of $\chi_{k,j}^{(+)}$, given by

[TABLE]

also contributes to the solution. Thus the solution outside the skyrmion should be given by the linear combination,

[TABLE]

where the coefficients $A_{k,j}$ and $B_{k,j}$ are determined by the boundary condition

[TABLE]

Thus we obtain the coupled equations for those coefficients,

[TABLE]

The equations can be exactly solved, yielding

[TABLE]

with

[TABLE]

In order to meet the orthogonality relation

[TABLE]

the whole wave function $\psi_{k,j}(\bm{r})$ should be divided by the factor $C_{k,j}\equiv\sqrt{A_{k,j}^{2}+B_{k,j}^{2}}$, yielding the final solution

[TABLE]

The total wave function $\Psi(\bm{r})$ is given by the linear combination of the eigenstates, namely

[TABLE]

The asymptotic behavior of this solution in the limit $\rho\rightarrow\infty$ is important in deriving the scattering amplitude. Here we fix the incident wave number to the Fermi momentum $k_{\mathrm{F}}=v_{\mathrm{F}}^{-1}\sqrt{E_{\mathrm{F}}^{2}-\Delta^{2}}$. Using the asymptotic behavior of the Bessel functions for $\xi/|l|\gg 1$,

[TABLE]

their linear combination satisfies

[TABLE]

where $\delta$ serves as the phase shift of the cylindrical wave in comparison with the original wave $J_{l}$. Thus, by introducing the phase shift

[TABLE]

with $\xi_{\mathrm{F}}=k_{\mathrm{F}}R_{\mathrm{S}}$ and $d=\Delta/E_{\mathrm{F}}$, we obtain the asymptotic behavior

[TABLE]

For large angular momentum $(|j|\gg k_{\mathrm{F}}R_{\mathrm{S}}=\xi_{\mathrm{F}})$, we can use the asymptotic behavior of the Bessel functions

[TABLE]

for $|l|\gg\xi$. The asymptotic behavior of the phase shift is given by

[TABLE]

where $\tilde{j}=|j|+1/2$. Thus we can see that the phase shift for large angular momentum is suppressed, as mentioned in the main text.

On the other hand, for large momentum, or for a large skyrmion, i.e. $\xi_{\mathrm{F}}=k_{\mathrm{F}}R_{\mathrm{S}}\gg|j|$, we can rely on the asymptotic behavior of the Bessel functions shown in Eqs. (83) and (84). As a result, we obtain

[TABLE]

where we have taken the leading order terms and have neglected the subleading terms of $O(\xi_{\mathrm{F}}^{-1})$.

Appendix B Derivation of scattering amplitude

Here we show the derivation process of the electron scattering amplitude by a single skyrmion, whose final result is given in Eq. (6). As mentioned in the main text, we consider here a case where an incident plane wave with the wave vector $\bm{k}$ and the energy $E=\sqrt{v_{F}^{2}|\bm{k}|^{2}+\Delta^{2}}$ is scattered by the skyrmion and comes out as a cylindrical wave with the wave number $k=|\bm{k}|$. The total wave function away from the skyrmion can be defined as the linear combination of the incident and outgoing waves,

[TABLE]

with the incident plane wave

[TABLE]

Here the outgoing wave component is defined by

[TABLE]

with $f_{\uparrow/\downarrow}(\phi)$ the scattering amplitude for each spin component.

As in the scattering theory of Schrödinger electrons, the scattering amplitude can be calculated in terms of the partial-wave decomposition. Taking the $x$ axis to the direction of the incident wave vector $\bm{k}$, the plane wave can be decomposed into partial waves as

[TABLE]

which shows the asymptotic behavior

[TABLE]

For the scattering amplitude, we define the partial-wave decomposition

[TABLE]

By comparing the both sides of Eq. (107) for each partial wave component in the limit $\rho\rightarrow\infty$, we obtain a sequence of equations

[TABLE]

where we use shorthand notations $\xi\equiv k\rho$, $\xi_{j}\equiv k\rho-(\pi/2)j$, and $l_{\uparrow/\downarrow}=j\mp 1/2$. Since the second term in the right-hand side of each equation does not contribute to the incoming cylindrical wave $e^{-ik\rho}$, the incoming parts of Eqs. (125) and (126) give the same relation,

[TABLE]

which fixes the value of the coefficient $\alpha_{j}$ as

[TABLE]

Substituting this relation to the outgoing $(e^{ik\rho})$ parts, we obtain

[TABLE]

Thus the total scattering amplitude is given by

[TABLE]

which further reduces to Eq. (6).

Appendix C Phase factor on transmitted plane wave at skyrmion boundary

Here we give a detailed discussion on the geometric phase factor acquired at the skyrmion boundary, which leads to the scattering skewness, as briefly mentioned in the main text. We first start from the simplified geometry, where the region with $n_{z}=+1$ and that with $n_{z}=-1$ is separated by the straight boundary at $x=0$, as shown in Fig. 4(a). If one takes a plane wave with the wave number $k$ and the incident angle $\theta$ as the incident wave, its wave function is given by

[TABLE]

Due to the translational symmetry in $y$-direction, the $y$-component of the momentum is conserved through the transmission and reflection processes, yielding the transmitted and reflected wave functions

[TABLE]

respectively. The complex transmission coefficient $t(\theta)$ and the reflection coefficient $r(\theta)$ should satisfy the boundary condition at $x=0$,

[TABLE]

Solving this equation, we obtain

[TABLE]

One can easily check that these factors satisfy the conservation of flux,

[TABLE]

The scattered wave, namely the deviation of the transmitted wave from the incident wave, is characterized by the factor

[TABLE]

which implies that the scattered wave acquires the geometric phase factor

[TABLE]

dependent on the incident angle $\theta$.

Next we set up a curved boundary, with the curvature radius $R_{\mathrm{S}}$, as shown in Fig. 4(b). Here we focus on the behavior at a long wavelength so that the uncertainty in the position can be neglected. Introducing an incident plane wave with the wave number $k$ in the $x$-direction, its incident angle $\theta$ on the boundary depends on the lateral incident position $y$, as $\theta(y)=\arcsin(y/R_{\mathrm{S}})$, leading to the phase factor $e^{i\gamma(\theta(y))}$ for the scattered wave. The phase $\gamma(\theta)$ defined by Eq. (154) can be expanded around $\theta=0$ as

[TABLE]

where $\tilde{k}=\mu/R_{\mathrm{S}}$. Noting the inequality

[TABLE]

the factors $\mu$ and $\tilde{k}$ here are negative. Therefore, the phase factor can be approximated as

[TABLE]

which can be interpreted as the transverse wave component, accounting for the scattering skewness. Since the phase $\gamma(\theta)$ takes the opposite sign if the sign of $n_{z}(\bm{r})$ is flipped, the scattering skewness is sensitive to the skyrmion number.