Coulomb drag in topological materials

TL;DR

This paper investigates Coulomb drag phenomena in topological insulator films hosting Dirac fermions, exploring effects in both magnetic and non-magnetic cases, and analyzing the influence of the anomalous Hall effect on drag resistivity.

Contribution

It provides a theoretical analysis of Coulomb drag in topological insulator layers, including both massless and massive Dirac fermions, with new insights into the role of the anomalous Hall effect.

Findings

Coulomb drag resistivity depends on temperature and electron density.

Magnetic topological insulators show significant anomalous Hall effect influence.

Analytical and numerical models elucidate the dynamics of spin density and drag effects.

Abstract

Dirac fermions are at the forefront of modern condensed matter physics research. They are known to occur in materials as diverse as graphene, topological insulators, and transition metal dichalcogenides, while closely related Weyl fermions have been discovered in other materials. They have been predicted to lend themselves to a variety of technological applications, while the recent prediction and discovery of the quantized anomalous Hall effect of massive Dirac fermions is regarded as a potential gateway towards low-energy electronics. Some materials hosting Dirac fermions are natural platforms for interlayer coherence effects such as Coulomb drag and exciton condensation. The top and bottom surfaces of a thin topological insulator film provide such a prototype system. Here we describe recent insights into Coulomb drag between two layers of Dirac fermions relying primarily on topological insulator films as a minimal model. We consider both non-magnetic topological insulators, hosting massless Dirac fermions, and magnetic topological insulators, in which the fermions are massive. We discuss in general terms the dynamics of the thin-film spin density matrix, outlining numerical results and approximate analytical expressions where appropriate for the drag resistivity \r{ho}D at low temperatures and low electron densities. In magnetic topological insulators with out-of-plane magnetizations in both the active and passive layers we analyze the role of the anomalous Hall effect in Coulomb drag.

Tables

Table 1. Table 1: A s 𝒌 ′ ​ s 𝒌 1 ′ ( i ) subscript superscript 𝐴 𝑖 subscript 𝑠 superscript 𝒌 ′ subscript 𝑠 subscript superscript 𝒌 ′ 1 A^{(i)}_{s_{\bm{k}^{\prime}}s_{{\bm{k}}^{\prime}_{1}}} and P s 𝒌 ​ s 𝒌 1 ​ s 𝒌 ′ ( i ) subscript superscript 𝑃 𝑖 subscript 𝑠 𝒌 subscript 𝑠 subscript 𝒌 1 subscript superscript 𝑠 ′ 𝒌 P^{(i)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}s^{\prime}_{\bm{k}}}

$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(1\right)}=-⟨s_{𝒌}|s_{{𝒌}_1}⟩\left[1-f_{0{𝒌}_1,s_{{𝒌}_1}}^{\left(\text{p}\right)}\right]⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩f_{0𝒌,s_{𝒌}^{\prime }}^{\left(\text{p}\right)}$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(2\right)}=⟨s_{𝒌}|s_{{𝒌}_1}⟩f_{0{𝒌}_1,s_{{𝒌}_1}}^{\left(\text{p}\right)}⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩\left[1-f_{0𝒌,s_{𝒌}^{\prime }}^{\left(\text{p}\right)}\right]$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(3\right)}=P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(8\right)}=-⟨s_{𝒌}|s_{{𝒌}_1}⟩f_{0{𝒌}_1,s_{{𝒌}_1}}^{\left(\text{p}\right)}⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(4\right)}=⟨s_{𝒌}|s_{{𝒌}_1}⟩⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩f_{0𝒌,s_{𝒌}^{\prime }}^{\left(\text{p}\right)}$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(5\right)}=\left[1-f_{0𝒌,s_{𝒌}}^{\left(\text{p}\right)}\right]⟨s_{𝒌}|s_{{𝒌}_1}⟩f_{0{𝒌}_1,s_{{𝒌}_1}}^{\left(\text{p}\right)}⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(6\right)}=-f_{0𝒌,s_{𝒌}}^{\left(\text{p}\right)}⟨s_{𝒌}|s_{{𝒌}_1}⟩\left[1-f_{0{𝒌}_1,s_{{𝒌}_1}}^{\left(\text{p}\right)}\right]⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩$
$P_{s_{𝒌}{s}_{{𝒌}_1}{s}_{𝒌}^{\prime }}^{\left(7\right)}=f_{0𝒌,s_{𝒌}}^{\left(\text{p}\right)}⟨s_{𝒌}|s_{{𝒌}_1}⟩⟨s_{{𝒌}_1}|s_{𝒌}^{\prime }⟩$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(1\right)}={\sum }_{{\varsigma }_{{𝒌}^{\prime }}}=⟨s_{{𝒌}^{\prime }}|s_{{𝒌}_1^{\prime }}⟩⟨s_{{𝒌}_1^{\prime }}|{\varsigma }_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },{\varsigma }_{{𝒌}^{\prime }}{s}_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(2\right)}={\sum }_{{\varsigma }_{{𝒌}_1^{\prime }}}=⟨s_{{𝒌}^{\prime }}|s_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },s_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}⟨{\varsigma }_{{𝒌}_1^{\prime }}|s_{{𝒌}^{\prime }}⟩$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(3\right)}=⟨s_{{𝒌}^{\prime }}|s_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },s_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}⟨{\varsigma }_{{𝒌}_1^{\prime }}|{\varsigma }_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },{\varsigma }_{{𝒌}^{\prime }}{s}_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(4\right)}=⟨s_{{𝒌}^{\prime }}|s_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },s_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}⟨{\varsigma }_{{𝒌}_1^{\prime }}|{\varsigma }_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },{\varsigma }_{{𝒌}^{\prime }}{s}_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(5\right)}={\sum }_{{\varsigma }_{{𝒌}_1^{\prime }}}=⟨s_{{𝒌}^{\prime }}|{\varsigma }_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },{\varsigma }_{{𝒌}_1^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}⟨s_{{𝒌}_1^{\prime }}|s_{{𝒌}^{\prime }}⟩$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(6\right)}={\sum }_{{\varsigma }_{{𝒌}^{\prime }}}⟨{\varsigma }_{{𝒌}^{\prime }}|s_{{𝒌}_1^{\prime }}⟩⟨s_{{𝒌}_1^{\prime }}|s_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },s_{{𝒌}^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(7\right)}=⟨s_{{𝒌}_1^{\prime }}|s_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },s_{{𝒌}^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}⟨{\varsigma }_{{𝒌}^{\prime }}|{\varsigma }_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },{\varsigma }_{{𝒌}_1^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}$
$A_{s_{{𝒌}^{\prime }}{s}_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}_1^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(8\right)}=⟨s_{{𝒌}_1^{\prime }}|s_{{𝒌}^{\prime }}⟩f_{{𝒌}^{\prime },s_{{𝒌}^{\prime }}{\varsigma }_{{𝒌}^{\prime }}}^{\left(\text{a}\right)}⟨{\varsigma }_{{𝒌}^{\prime }}|{\varsigma }_{{𝒌}_1^{\prime }}⟩f_{{𝒌}_1^{\prime },{\varsigma }_{{𝒌}_1^{\prime }}{s}_{{𝒌}_1^{\prime }}}^{\left(\text{a}\right)}$

Equations

$${\rho }_{\text{D}}=-\frac{\hslash }{e^2}\frac{\zeta (3)}{16\pi }\frac{(k_{\text{B}}T{)}^2}{A^2{r}_s^2{n}_{\text{a}}^{\frac{3}{2}}{n}_{\text{p}}^{\frac{3}{2}}{d}^4},$$

$$H_{0\mathbf{k}}={\tau }_z\otimes h_{\mathbf{k}}+\text{diag}\{M_{\text{a}},-M_{\text{a}},M_{\text{p}},-M_{\text{p}}\},$$

$$V(\mathbf{q},\omega )=\frac{v_q{\text{e}}^{-qd}}{ϵ(\mathbf{q},\omega )}.$$

$$\begin{array}[]{rl}\epsilon({\bm{q}},\omega)&\displaystyle=[1-v_{q}\Pi_{\text{a}}({\bm{q}},\omega)][1-v_{q}\Pi_{\text{p}}({\bm{q}},\omega)]\\[12.91663pt] &\displaystyle-[v_{q}\text{e}^{-qd}]^{2}\Pi_{\text{a}}({\bm{q}},\omega)\Pi_{\text{p}}({\bm{q}},\omega),\end{array}$$

$${\Pi }_l(\mathbf{q},\omega )=-\frac{1}{L^2}\sum _{\mathbf{k}s{s}^{\prime }}\frac{(f_{0\mathbf{k},s}^{(l)}-f_{0{\mathbf{k}}^{\prime },s^{\prime }}^{(l)})}{{\epsilon }_{\mathbf{k},s}^{(l)}-{\epsilon }_{{\mathbf{k}}^{\prime },s^{\prime }}^{(l)}+\hslash \omega +i{0}^{+}}{F}_{\mathbf{k}{\mathbf{k}}^{\prime }}^{s{s}^{\prime }},$$

$$\frac{\text{d}\hat{F}}{\text{d}t}+\frac{i}{\hslash }[\hat{H},\hat{F}]=0.$$

$$\frac{\text{d}{f}_{\mathbf{k}}}{\text{d}t}+\frac{i}{\hslash }[H^{1e},f_{\mathbf{k}}]+{\hat{J}}_{ee}(f_{\mathbf{k}})=0.$$

$$\begin{array}[]{rl}&\displaystyle J_{ee}^{i}(f_{\bm{k}})_{s_{\bm{k}}s^{\prime}_{\bm{k}}}=\frac{\pi}{\hbar L^{4}}\sum_{{\bm{k}}_{1}{\bm{k}}^{\prime}{{\bm{k}}^{\prime}_{1}}}|v^{(\text{pa})}_{|{\bm{k}}-{\bm{k}}_{1}|}|^{2}\delta_{{\bm{k}}+{\bm{k}}^{\prime},{\bm{k}}_{1}+{\bm{k}}^{\prime}_{1}}\big{\{}\\[8.61108pt] &\displaystyle\sum_{i=1}^{4}\!P^{(i)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}s^{\prime}_{\bm{k}}}A^{(i)}_{s_{\bm{k}^{\prime}}s_{{\bm{k}}^{\prime}_{1}}}\delta[\varepsilon^{(\text{p})}_{k_{1},s_{{\bm{k}}_{1}}}\!-\varepsilon^{(\text{p})}_{k,s^{\prime}_{\bm{k}}}\!+\varepsilon^{(\text{a})}_{k^{\prime}_{1},s_{{\bm{k}}^{\prime}_{1}}}\!-\varepsilon^{(\text{a})}_{k^{\prime},s_{\bm{k}^{\prime}}}]+\\[8.61108pt] &\displaystyle\sum_{i=5}^{8}\!P^{(i)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}s^{\prime}_{\bm{k}}}A^{(i)}_{s_{\bm{k}^{\prime}}s_{{\bm{k}}^{\prime}_{1}}}\delta[\varepsilon^{(\text{p})}_{k_{1},s_{{\bm{k}}_{1}}}\!-\varepsilon^{(\text{p})}_{k,s_{\bm{k}}}\!+\varepsilon^{(\text{a})}_{k^{\prime}_{1},s_{{\bm{k}}^{\prime}_{1}}}\!-\varepsilon^{(\text{a})}_{k^{\prime},s_{\bm{k}^{\prime}}}]\big{\}},\end{array}$$

$$\text{vspace}-0.5cm\frac{\text{d}{f}_{\mathbf{k}}^{(\text{a})}}{\text{d}t}+\frac{i}{\hslash }[H_{0\mathbf{k}}^{(\text{a})},f_{\mathbf{k}}^{(\text{a})}]+{\hat{J}}_0(f_{\mathbf{k}}^{(\text{a})})=\frac{e\mathbf{E}}{\hslash }\cdot \frac{\partial f_{0\mathbf{k}}^{(\text{a})}}{\partial \mathbf{k}},$$

$$\frac{\text{d}{f}_{\mathbf{k}}^{(\text{p})}}{\text{d}t}+\frac{i}{\hslash }[H_{0\mathbf{k}}^{(\text{p})},f_{\mathbf{k}}^{(\text{p})}]+{\hat{J}}_0(f_{\mathbf{k}}^{(\text{p})})=J_{ee}^i(f_{\mathbf{k}}{)}_{s_{\mathbf{k}}{s}_{\mathbf{k}}^{\prime }},$$

$$\begin{array}[]{rl}J^{i}_{ee}(f_{\bm{k}})_{s_{\bm{k}}s^{\prime}_{\bm{k}}}&\displaystyle=\!\frac{J^{i}_{ee}(f_{\bm{k}})_{++}+J^{i}_{ee}(f_{\bm{k}})_{--}}{2}\sigma_{0}\\[12.91663pt] &\displaystyle+\frac{J^{i}_{ee}(f_{\bm{k}})_{++}-J^{i}_{ee}(f_{\bm{k}})_{--}}{2}\sigma_{z}\\[12.91663pt] &\displaystyle+\frac{J^{i}_{ee}(f_{\bm{k}})_{-+}+J^{i}_{ee}(f_{\bm{k}})_{+-}}{2}\sigma_{x}\\[12.91663pt] &\displaystyle+\frac{J^{i}_{ee}(f_{\bm{k}})_{-+}-J^{i}_{ee}(f_{\bm{k}})_{+-}}{2i}\sigma_{y},\end{array}$$

$$\text{T}=\left(\begin{array}[]{cc}\displaystyle-i\text{e}^{i\theta_{\bm{k}}}\sqrt{\frac{1+b_{k}}{2}}&\quad\displaystyle\frac{a_{k}}{\sqrt{2(1+b_{k})}}\\[12.91663pt] \displaystyle i\text{e}^{i\theta_{\bm{k}}}\sqrt{\frac{1-b_{k}}{2}}&\quad\displaystyle\frac{a_{k}}{\sqrt{2(1-b_{k})}}\end{array}\right),$$

$$\begin{array}[]{rl}J^{i}_{ee}(f_{\bm{k}})_{\text{ll,lh}}&\displaystyle=-\!\frac{J^{i}_{ee}(f_{\bm{k}})_{++}+J^{i}_{ee}(f_{\bm{k}})_{--}}{2}\sigma_{0}\\[12.91663pt] &\displaystyle-\frac{J^{i}_{ee}(f_{\bm{k}})_{++}-J^{i}_{ee}(f_{\bm{k}})_{--}}{2}{\bm{\sigma}}\cdot\hat{\bm{\Omega}}_{\bm{k}},\\[12.91663pt] J^{i}_{ee}(f_{\bm{k}})_{\text{hl,hh}}&\displaystyle=\frac{J^{i}_{ee}(f_{\bm{k}})_{-+}\!+\!J^{i}_{ee}(f_{\bm{k}})_{+-}}{2}{\bm{\sigma}}\cdot\hat{\bm{z}}_{\text{eff}}\\[12.91663pt] &\displaystyle-\frac{J^{i}_{ee}(f_{\bm{k}})_{-+}\!-\!J^{i}_{ee}(f_{\bm{k}})_{+-}}{2i}{\bm{\sigma}}\cdot\hat{\bm{k}}_{\text{eff}},\end{array}$$

$${\sigma }_{\text{D}}=\frac{e^2\beta }{16\pi }\sum _{\mathbf{q}}\int \text{d}\omega \frac{|v_q^{(\text{pa})}{|}^2\text{Im}{\chi }_{\text{a}}(\mathbf{q},\omega )\text{Im}{\chi }_{\text{p}}(\mathbf{q},\omega )}{{\sinh }^2\frac{\beta \hslash \omega }{2}},$$

$${\rho }_{\text{D}}\approx -\frac{{\sigma }_{\text{D}}}{{\sigma }_{\text{a}}{\sigma }_{\text{p}}}=-\frac{\hslash }{e^2}\frac{\zeta (3)}{16\pi }\frac{(k_{\text{B}}T{)}^2}{A^2{r}_s^2{n}_{\text{a}}^{\frac{3}{2}}{n}_{\text{p}}^{\frac{3}{2}}{d}^4},$$

$$\begin{array}[]{rl}P_{\parallel}\hat{J}(f^{(l)}_{{\bm{k}}\parallel})&\displaystyle=\frac{n_{i}}{\pi\hbar A}\int d{k^{\prime}_{x}}|\overline{U}_{{\bm{k}}{\bm{k}}^{\prime}}|^{2}(f^{(l)}_{{\bm{k}}\parallel}-f^{(l)}_{{\bm{k}}^{\prime}\parallel})\sigma_{z}\delta(k_{x}-k^{\prime}_{x}),\\[12.91663pt] \end{array}$$

$$\begin{array}[]{rl}J^{i}(f_{\bm{k}})_{s_{\bm{k}}}\displaystyle=&\displaystyle-\frac{2\pi}{\hbar}\!\sum_{{\bm{k}}^{\prime}}\!|v_{|{\bm{k}}-{\bm{k}}_{1}|}^{(\text{pa})}|^{2}F^{\text{(p)}}_{s_{\bm{k}}s_{{\bm{k}}_{1}}}F^{\text{(a)}}_{s_{{\bm{k}}^{\prime}}s_{{\bm{k}}^{\prime}_{1}}}\\[12.91663pt] &\displaystyle\delta[\varepsilon^{(\text{p})}_{k_{1},s_{{\bm{k}}_{1}}}\!\!-\!\varepsilon^{(\text{p})}_{k,s_{\bm{k}}}\!+\!\varepsilon^{(\text{a})}_{k^{\prime}_{1},s_{{\bm{k}}^{\prime}_{1}}}\!-\!\varepsilon^{(\text{a})}_{k^{\prime},s_{\bm{k}^{\prime}}}\!]\\[12.91663pt] &\displaystyle\times\big{\{}f^{(\text{p})}_{{\bm{k}},s_{\bm{k}}}[1\!-\!f^{(\text{p})}_{{\bm{k}}_{1},s_{{\bm{k}}_{1}}}]f^{(\text{a})}_{{\bm{k}}^{\prime},s_{\bm{k}^{\prime}}}[1\!-\!f^{(\text{a})}_{{\bm{k}}^{\prime}_{1},s_{{\bm{k}}^{\prime}_{1}}}]\\[12.91663pt] &\displaystyle-[1\!-\!f^{(\text{p})}_{{\bm{k}},s_{\bm{k}}}]f^{(\text{p})}_{{\bm{k}}_{1},s_{{\bm{k}}_{1}}}[1\!-\!f^{(\text{a})}_{{\bm{k}}^{\prime},s_{\bm{k}^{\prime}}}]f^{(\text{a})}_{{\bm{k}}^{\prime}_{1},s_{{\bm{k}}^{\prime}_{1}}}\!\big{\}}.\end{array}$$

$$H_{\mathbf{k}}=A{\tau }_z\otimes [\mathbf{\sigma }\cdot (\mathbf{k}\times \hat{\mathbf{z}})]+t{\tau }_x\otimes 1,$$

$$\begin{array}[]{rl}\displaystyle\hat{J}(f_{\bm{k}})=&\displaystyle\frac{n_{i}\epsilon_{k}}{4\pi A^{2}\hbar}\int_{0}^{2\pi}\,d\theta_{\bm{k}^{\prime}}\,|U_{{\bm{k}}{\bm{k}^{\prime}}}|^{2}(f_{\bm{k}}-f_{\bm{k}^{\prime}})\\[12.91663pt] &\displaystyle\times\bigg{(}1+\frac{t^{2}}{t^{2}+A^{2}k^{2}}+\frac{A^{2}k^{2}\cos\gamma}{t^{2}+A^{2}k^{2}}\bigg{)},\end{array}$$

Full text

††journal: Journal of Physics and Chemistry of Solids

1 Introduction

The energetic study of Dirac fermions in condensed matter physics is currently in full bloom and is driven primarily by the realization of their considerable potential for spin electronics, thermoelectricity, magnetoelectronics and topological quantum computing [1]. Exotic physical phenomena associated with Dirac fermions have been identified in materials ranging from graphene [2] to topological insulators, [3] transition metal dichalcogenides [4] and Weyl and Dirac semimetals. [5, 6, 7] The linear dispersion of massless Dirac fermions has aroused considerable interest experimentally [8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21] and theoretically [22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34]. Among Dirac fermion systems, three-dimensional topological insulators (TIs) stand out as a novel class of bulk insulating materials exhibiting conducting surface states with a chiral spin texture[2, 3, 35, 36, 37, 38, 39, 40, 41], with the surface carriers being massless Dirac fermions. Magnetic topological insulators, hosting massive Dirac fermions, have also been successfully manufactured [42, 12, 43]. Time-reversal symmetry breaking gives Dirac fermions a finite mass and results in a non-trivial Berry curvature [4, 44, 45, 46, 47], which is at the root of the anomalous [48, 28] and quantum anomalous Hall effects [49] detected in these materials [50, 51, 52]. In magnetic TI slabs a dissipationless quantized anomalous Hall effect has been discovered [53, 54, 52, 55], which has already been harnessed successfully [56], stimulating an intense search for device applications. In a related context hybrid TI/superconductor junctions have been fabricated [57, 58], in which topological superconductivity and Majorana fermions may be present [59, 60].

Three-dimensional topological insulator slabs are inherently two-layer systems naturally exhibiting interlayer coherence effects such as Coulomb drag [61, 62, 63], in which the charge current in one layer drags a charge current in the adjacent layer through the interlayer electron-electron interactions, as shown in Fig. 1. Physically, Coulomb drag is caused by the transfer of momentum between electrons in different layers due to the interlayer electron-electron scattering. It constitutes an experimental probe of electron-electron interactions [61, 64, 65] and has recently attracted considerable attention in massless Dirac fermion systems such as graphene [66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79]. A significant motivation for the intensive experimental effort dedicated to Coulomb drag in Dirac fermion systems is the search for exciton condensation in electron-hole bilayers [80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 64, 90, 69, 70, 72, 66, 67, 68, 71, 73, 74, 75, 76, 77, 78, 79, 91, 92, 93, 94]. At the same time Coulomb drag of massive Dirac fermions is associated with an important question: if topological terms are present in the drag current they could be exploited in longitudinal transport, potentially enabling a topological transistor.

In this article we discuss Coulomb drag of both massless and massive Dirac fermions. We take as our prototype model the surface states of 3D topological insulators and we consider two cases: (i) massless Dirac fermions, in which neither layer is magnetized; (ii) massive Dirac fermions, in which an out of plane magnetization exists in each layer, but no external magnetic field. The physics discussed is therefore that of 2D Dirac fermions: the effects are driven by the electrons in the surface conduction band(s) of a topological insulator. Case (i) provides a useful contrast to graphene. Unlike graphene, in TIs the spin and orbital degrees of freedom are coupled by the strong spin-orbit interaction, TIs have an odd number of valleys on a single surface, and the relative permittivity is different, while in known band TIs screening is qualitatively and quantitatively different, since it does not involve the interplay of the layer and valley degrees of freedom. All these features impact the drag current. We introduce a density matrix method to calculate the Coulomb drag current in topological insulator films, which fully takes into account the spin degree of freedom and interband coherence. The drag resistivity in case (i) takes the analytical form

[TABLE]

where $k_{\text{B}}$ is the Boltzmann constant, $A$ is the TI spin-orbit constant, $r_{s}$ is the Wigner-Seitz radius (effective fine structure constant) which represents the ratio of the electrons’ average Coulomb potential and kinetic energies, $d$ is the layer separation and $n_{\text{a},\text{p}}$ are the electron densities in the active and passive layers, respectively. For a single-valley system $r_{s}=e^{2}/(2\pi\epsilon_{0}\epsilon_{r}A)$, with $\epsilon_{r}$ the relative permittivity. Whereas the full numerical results are more complicated, the same trends are present as in the analytical formula.

Case (ii) is significantly more complicated. The drag current in the passive layer has a longitudinal component and a Hall component, which can be measured separately. The Hall drag current in principle has contributions from (i) the longitudinal current in the active layer via a transverse drag force, (ii) the Hall current in the active layer, which may be termed direct Hall drag (Fig. 2). Contribution (i) represents straightforward longitudinal charge transport in the presence of impurities and other scattering mechanisms. Contribution (ii) is topological: the anomalous Hall effect in Dirac fermion systems is predominantly driven by a topological term originating in the Berry curvature associated with the band structure. The Hall drag current stemming from the transverse drag force of (i) in effect comes from the conventional longitudinal drag current.

Interestingly, the (topological) anomalous Hall current in the active layer does not generate a drag current at all in the passive layer - both the longitudinal and the Hall components stemming from these topological terms vanish identically. Secondly, the remaining terms in the anomalous Hall drag current only depend on $M_{\text{p}}$, the magnetization of the passive layer. The dependence on $M_{\text{p}}$ is non-monotonic, with a peak at an intermediate value of $M_{\text{p}}$, which becomes pronounced at low densities. These results reveal the absence of a Hall drag per se in anomalous Hall systems, and the fact that the measured Hall drag current is simply the anomalous Hall response to the longitudinal drag force. On a deeper level, given that the Berry curvature in solid state systems represents interband coherence, the results discussed here imply that interband coherence cannot be transferred between layers, at least not by electron-electron scattering. This is in sharp contrast to conventional Coulomb drag in ordinary Hall systems [64, 65, 78, 79, 76]. There the Hall current in the active layer is caused by the Lorentz force rather than topology or interband coherence and makes a significant contribution to the longitudinal and Hall drag response, the latter of which depends on the applied magnetic field.

Although our discussion is general and applies to materials with Rashba spin-orbit interactions even beyond topological insulators, the question of drag has a particular significance to the latter. Following their initial observation [95, 96, 97, 98, 99, 100, 101, 8, 9], improvements in TI growth have made them suitable for fundamental research [102, 16, 103, 19]. A number of experiments have successfully identified signatures of the surface states in transport [9, 16, 103, 102, 13]. A current induced spin polarization also constitutes a signature of surface transport [104, 105] and was reported in recent experimental studies [106, 14, 20]. However, the reliable and reproducible identification of the surface states in transport, which remains the key to topological insulators becoming technologically important, has remained elusive. At the same time the interplay of strong spin-orbit coupling and electron-electron interactions in topological insulators is at present not completely understood [107, 108, 109, 30, 110, 111, 31, 33], with transport experiments and theoretical mostly focused on single-particle effects [29, 112, 113, 114, 15, 115, 51, 116, 117, 118, 119, 24, 18, 10, 120, 121, 122]. As this work will show, plenty of scope exists for investigating interacting setups such as Coulomb drag geometries.

The outline of this paper is as follows. In Sec. 2 we discuss a general theoretical framework for understanding Coulomb drag in systems in which spin/interband coherence is important, which requires the use of the spin density matrix. We discuss Coulomb drag for massless Dirac fermions in Sec. 3 and massive Dirac fermions in Sec. 4. Generalisations of our picture beyond the prototype system of a topological insulator film are presented in Sec. 5. We end with a summary and conclusions.

2 Transport theory

Quite generally, Dirac fermions in two dimensions are described by an effective Hamiltonian $H_{\text{D}}=A\,{\bm{\sigma}}\cdot({\bm{k}}\times\hat{{\bm{z}}})+M\sigma_{z}$, with ${\bm{\sigma}}=(\sigma_{x},\sigma_{y},\sigma_{z})$ the usual Pauli matrices, ${\bm{k}}=(k_{x},k_{y})$ the 2D wave vector, $A$ representing the Fermi velocity and $M$ a generic mass term. In the limit $M\rightarrow 0$ the quasi-particle dispersion is linear.

We consider a TI film with both the top (active) and bottom (passive) surfaces magnetized either through doping with magnetic impurities or proximity coupling to ferromagnets. The magnetizations of the two surfaces are allowed to differ. Without loss of generality we assume (i) the carrier number density and hence the Fermi energy is the same in each layer and (ii) sidewall states do not participate in transport, an assumption that recent work has shown to be justified [123]. The chemical potential lies in the surface conduction band of each layer. We require $\varepsilon_{\text{F}}\tau_{l}/\hbar\gg 1$ in each layer $l\in\{\text{a}\equiv\text{active},\text{p}\equiv\text{passive}\}$, with $\varepsilon_{\text{F}}$ the Fermi energy located in the surface conduction bands and $\tau_{l}$ the momentum scattering time. The two-layer effective band Hamiltonian

[TABLE]

where the (Rashba) Hamiltonian of a single layer $h_{{\bm{k}}}=A\,{\bm{\sigma}}\cdot({\bm{k}}\times\hat{{\bm{z}}})\equiv-Ak{\bm{\sigma}}\cdot\hat{\bm{\theta}}$ with $\hat{\bm{\theta}}$ the tangential unit vector corresponding to ${\bm{k}}$. The Pauli matrix $\tau_{z}$ represents the layer degree of freedom. The eigenvalues of Eq. (2) are $\varepsilon_{l\pm}={\pm}\sqrt{A^{2}k^{2}+M^{2}_{l}}\equiv\hbar\Omega_{k}^{(l)}$, the band index $s_{\bm{k}}=\pm$ with $+$ the conduction band and $-$ the valence band.

The single-particle Hamiltonian $\hat{H}^{1e}=\hat{H}_{0}+\hat{H}_{E}+\hat{U}$, where $\hat{H}_{0}$ is the band Hamiltonian defined in Eq. (2), $\hat{H}_{E}=e(\hat{\bm{E}}\otimes\sigma_{0})\cdot\hat{\bm{r}}$ is the electrostatic potential due to the driving electric field with $\hat{\bm{r}}$ the position operator and $\hat{U}$ the disorder potential, which is assumed to be a scalar in spin space. Adding the two-particle interaction term we write the total Hamiltonian as $\hat{H}=\hat{H}^{1e}+\hat{V}^{ee}$ with the single-particle term expressed generically as $\hat{H}^{1e}\!=\!\sum_{\alpha\beta}H_{\alpha\beta}c^{{\dagger}}_{\alpha}c_{\beta}$ and the Coulomb interaction term $\hat{V}^{ee}\!=\!\frac{1}{2}\sum_{\alpha\beta\gamma\delta}V^{ee}_{\alpha\beta\gamma\delta}c^{\dagger}_{\alpha}c^{\dagger}_{\beta}c_{\gamma}c_{\delta}$. The indices $\alpha\equiv{\bm{k}}s_{\bm{k}}l$ represent wave vector, band, and layer indices respectively. The matrix element $V^{ee}_{\alpha\beta\gamma\delta}$ in a generic basis $\{\phi_{\alpha}({\bm{r}})\}$ is given by $V^{ee}_{\alpha\beta\gamma\delta}=\int\!d{\bm{r}}\int\!d{\bm{r}}^{\prime}\ \phi^{*}_{\alpha}({\bm{r}})\phi^{*}_{\beta}({\bm{r}}^{\prime})V^{ee}_{{\bm{r}}-{\bm{r}}^{\prime}}\phi_{\delta}({\bm{r}})\phi_{\gamma}({\bm{r}}^{\prime})$, where $V^{ee}_{{\bm{r}}-{\bm{r}}^{\prime}}$ is the Coulomb interaction in real space.

The dynamical screening function is found by solving the Dyson equation for the two-layer system in the random phase approximation (RPA) as in Ref. [64]. The dynamically screened interlayer Coulomb interaction

[TABLE]

The dielectric function of the coupled layer system is

[TABLE]

in which the polarization function is obtained by summing the lowest bubble diagram and takes the form

[TABLE]

with $f^{(l)}_{0{\bm{k}},s}\equiv f^{(l)}_{0}(\varepsilon_{{\bm{k}}s})$ the equilibrium Fermi distribution function and $F^{ss^{\prime}}_{{\bm{k}}{\bm{k}}^{\prime}}=\langle{\bm{k}},s,l|{\bm{k}}^{\prime},s^{\prime},l^{\prime}\rangle\langle{\bm{k}}^{\prime},s^{\prime},l^{\prime}|{\bm{k}},s,l\rangle$ is the wave function overlap. In magnetic TIs, due to the typical smallness of the gap opened by the magnetization, $F^{++}_{{\bm{k}}{\bm{k}}^{\prime}}\approx\frac{1+\cos\gamma}{2}$, with $\gamma$ is the angle between the incident (${\bm{k}}$) and the outgoing wave vector (${\bm{k}}^{\prime}$).

Two parameters quantify the strength of interactions in topological insulator films and the limits of validity of the RPA. Firstly, we define an effective background dielectric constant $\epsilon_{r}$. The physics of films is determined by their thickness $d$ and the Fermi wave vector $k_{F}$ [123]. We take a Bi2Se3 film as an example, with $\epsilon_{r,Bi2Se3}\approx 100$, grown on a semiconductor substrate with $\epsilon_{r,s}\approx 11$. For $k_{F}d\gg 1$, the film is thick, and the two surfaces are independent [124]. For the top surface, where one side is exposed to air, $\epsilon_{r,top}=(\epsilon_{r,Bi2Se3}+1)/2\approx 50$. For the bottom surface $\epsilon_{r,btm}=(\epsilon_{r,Bi2Se3}+\epsilon_{r,s})/2\approx 55$. Both are independent of $d$. For $k_{F}d\ll 1$, the film is ultrathin and can be approximated as a pure 2D system, with $\epsilon_{r}=(1+\epsilon_{r,s})/2\approx 6$, also independent of $d$. However, since the TI bulk cannot be eliminated, $\epsilon_{r}\approx 6$ is an ideal lower bound. In films studied here, thick enough that there is no interlayer tunnelling, implying $d>5$ nm at the very least, $\epsilon_{r}$ has contributions from both the TI bulk and the semiconductor substrate, and for Bi2Se3 can range between 6 and 55. Two experiments have extracted $\epsilon_{r}\approx 30$ for relatively thick films of Bi2Se3 (10 nm $<d<$ 20 nm) [125]. Hence, $\epsilon_{r}$ is treated as a phenomenological parameter to be measured separately for each film. In general $r_{s}\ll 1$ making the RPA approach applicable. The approximation $k_{\text{F}}d\gg 1$ applies for TI film thicknesses up to $6\text{nm}$, with electron densities $n\thicksim 10^{12}\text{cm}^{-2}$. This is in contrast to graphene, where $k_{\text{F}}d\ll 1$ is usually satisfied. For other systems with small enough $k_{\text{F}}d$ and larger $r_{s}$, the interlayer Coulomb interaction is in general not small, and it remains to be determined whether it is necessary to take into account higher-order contributions in the Coulomb interaction.

The many-body density matrix $\hat{F}$ obeys the quantum Liouville equation [126]

[TABLE]

The one-particle reduced density matrix is the trace $\text{tr}(c^{\dagger}_{\eta}c_{\xi}\hat{F})\equiv\langle c^{\dagger}_{\eta}c_{\xi}\rangle$. Its ${\bm{k}}$-diagonal part can be written as a $4\times 4$ matrix in the joint spin/layer pseudospin subspace, and we refer to this matrix as $f_{\bm{k}}$. To second order in the electron-electron interaction, $f_{\bm{k}}$ satisfies [107]

[TABLE]

The term $\hat{J}_{ee}(f_{\bm{k}})$ in Eq. (7) represents intralayer and interlayer electron-electron scattering operator. Since the intralayer electron-electron scattering does not contribute to the drag current, we concentrate on the interlayer term which is the matrix form of operator $\hat{J}_{ee}(f_{\bm{k}})$,

[TABLE]

where $L^{2}$ is the area of the 2D system and the quantities $A^{(i)}_{s_{\bm{k}^{\prime}}s_{{\bm{k}}^{\prime}_{1}}}$ and $P^{(i)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}s^{\prime}_{\bm{k}}}$, given explicitly in the Table \colorblue 1, are functions of the occupations of the active and passive layers respectively. The interlayer momentum transfer ${\bm{q}}={\bm{k}}-{\bm{k}}_{1}={\bm{k}}^{\prime}_{1}-{\bm{k}}^{\prime}$, and $v^{(\text{pa})}_{|{\bm{k}}-{\bm{k}}_{1}|}$ is the interlayer Coulomb interaction. The interlayer electron-electron scattering processes that contribute to Coulomb drag are shown in Fig. 3.

At this stage one may include explicitly disorder and driving electric field in the one particle Hamiltonian and write $\hat{H}^{1e}=\hat{H}_{0}+\hat{H}_{E}+\hat{U}$, where $\hat{H}_{0}$ is the band Hamiltonian of TIs, $\hat{H}_{E}=e\hat{\bm{E}}\cdot\hat{\bm{r}}$ is the electrostatic potential due to the driving electric flied with $\hat{\bm{r}}$ is a position operator and $\hat{U}$ is the disorder potential. Referring to the layer index (l), we write $f^{(l)}_{\bm{k}}=n^{(l)}_{\bm{k}}\mathbb{1}+S^{(l)}_{\bm{k}}$, with $S^{(l)}_{\bm{k}}$ a $2\times 2$ Hermitian matrix which can be decomposed in terms of the Pauli spin matrices. The single-particle density matrix is diagonal in the layer index $f_{\bm{k}}={\rm{diag}}\,(f_{\bm{k}}^{(\text{a})},f_{\bm{k}}^{(\text{p})})$. The scattering term $J^{i}_{ee}(f_{\bm{k}})_{s_{\bm{k}}s^{\prime}_{\bm{k}}}$ represents interlayer coherence: in the active layer we solve for the non-equilibrium correction to $f_{\bm{k}}^{(\text{a})}$, feed it into $J^{i}_{ee}(f_{\bm{k}})_{s_{\bm{k}}s^{\prime}_{\bm{k}}}$, and the result is an effective driving term for the passive layer. This then enables us to solve for the non-equilibrium to correction $f_{\bm{k}}^{(\text{p})}$ in response to this effective driving term. The kinetic equation of the two-layer system

[TABLE]

with $H^{(l)}_{0{\bm{k}}}=h^{(l)}_{\bm{k}}+M_{l}\sigma_{z}$ the band Hamiltonian and $f^{(l)}_{0{\bm{k}}}$ the equilibrium density matrix of each magnetic layer. The electron-impurity scattering integral is \hat{J}_{0}(f^{(l)}_{\bm{k}})\leavevmode\nobreak\ =\leavevmode\nobreak\ \bigg{\langle}\int^{\infty}_{0}\frac{\text{d}t^{\prime}}{\hbar^{2}}[\hat{U},e^{-i\hat{H}t^{\prime}/\hbar}[\hat{U},\hat{f}]e^{i\hat{H}t^{\prime}/\hbar}]\bigg{\rangle}_{{\bm{k}}{\bm{k}}}, where $\langle\rangle$ denotes the average over impurity configurations. In EQ. (9b), $J^{i}_{ee}(f_{\bm{k}})_{s_{\bm{k}}s^{\prime}_{\bm{k}}}$ can be rewritten as

[TABLE]

with $\sigma_{0}$ is $2\times 2$ unit matrix and $\sigma_{x,y,z}$ are Pauli matrix. We should transform Eq. 10 from eigenstates representation (ER) to that of of $\sigma_{z}$ (ordinary representation, OR) with

[TABLE]

and for any $2\times 2$ matrix $\mathcal{M}$, $\mathcal{M}^{\text{OR}}=\text{T}^{{\dagger}}\cdot\mathcal{M}^{\text{ER}}\cdot\text{T}$. We have $\sigma_{z}\to{\bm{\sigma}}\cdot\hat{{\bm{\Omega}}}_{\bm{k}},\sigma_{x}\to-{\bm{\sigma}}\cdot\hat{{\bm{z}}}_{\text{eff}},\sigma_{y}\to{\bm{\sigma}}\cdot\hat{{\bm{k}}}_{\text{eff}}$. In $\sigma_{z}$ representation, according to the contributions to various drag currents, the electron-electron scattering term is divided into

[TABLE]

where ll,lh,hl,hh represent contributions to the four drag current components in Fig. 2.

The above encapsulates the philosophy or our approach: to begin with, we consider an external electric field applied to the active layer, we solve for the non-equilibrium density matrix in the active layer without any reference to electron-electron scattering, and we feed the solution for the active-layer density matrix into the interlayer electron-electron scattering integral. The result of this is an electric-field dependent term that acts as a new driving term (i.e. the drag force) for the passive layer. We then solve for the non-equilibrium density matrix in the passive layer with this driving term. The solution thus obtained represents the non-equilibrium density matrix in the passive layer, which may be used to calculate expectation values. Specifically, its trace with the current operator yields the drag current.

3 Massless Dirac fermions

Physically the drag current is a result of the rectification by the passive layer of the fluctuating electric field created by the active layer [64]. At low temperatures the predominant contribution to drag is due to intraband processes near the Fermi surface. When the Fermi level is finite, i.e. above the Dirac point, electrons take more energy to transition from the valence band to the conduction band than to transition within the conduction band or valence band, and with a small excitation energy the channel involving interband transitions becomes inaccessible [69]. The small interlayer momentum transfer and excitation energy is the dominant part of the spectrum contributing to the drag current.

The longitudinal drag conductivity $\sigma_{D}$ is

[TABLE]

where $\text{Im}\chi_{l}({\bm{q}},\omega)$ is the imaginary part of nonlinear susceptibility and $\beta\!=\!\frac{1}{k_{\text{B}}T}$. Using standard approximations we can express the drag resistivity analytically as

[TABLE]

where $\sigma_{\text{a,p}}$ is the longitudinal conductivity of each layer.

In Fig. 4(a) we present numerical results for the dependence of the Coulomb drag resistivity on the electron number density. The drag resistivity displays a $\frac{1}{n^{\alpha}_{\text{a}}n^{\alpha}_{\text{p}}}$ dependence where $\alpha$ is the power exponent of the density $n_{\text{a,p}}$, with $\alpha<1.5$ for $d=10,20,40\leavevmode\nobreak\ \text{nm}$. With increasing electron density $n$ the coefficient $\alpha$ approaches $1.5$. The fact that the exact numerical results shown in Fig. 4(a) disagree more strongly with the analytical results for smaller values of $n$ and $d$ is understandable, since Eq. (14) applies only in the $k_{\text{F}}d\gg 1$ limit with $k_{\text{F}}=\sqrt{4\pi n}$. This trend of an increasing quantitative failure of the asymptotic analytical drag formula for small $k_{\text{F}}d$ has also been noted in graphene [69, 68, 71] and 2DEG systems [64, 65]. The analytical result becomes more accurate with increasing $k_{\text{F}}d$.

In Fig. 4(b) we show the Coulomb drag resistivity as a function of temperature $T$ for three different thicknesses $d=10,20,40\leavevmode\nobreak\ \text{nm}$. The overall temperature dependence of the drag resistivity increases nearly quadratically and there is no logarithmic correction due to the absence of backscattering in TIs. The $T^{2}$ dependence stems from the allowed phase space where electron-electron scattering occurs at low temperature, and is expected for any interaction strength between the top and bottom layers of TIs as Fig. 4(b) shows, provided that the carriers can be described using a Fermi liquid picture. In addition, in TIs the acoustic phonon velocity is smaller than in graphene. These facts make the contribution of electron-phonon scattering processes to the resistivity much more important in the surface of 3DTIs than in graphene. For the surfaces of 3DTIs the effect of electron-phonon scattering events becomes important already for $T$ as low as $10\text{K}$. For this reason we consider temperatures up to $10\text{K}$ in our numerical calculations. In graphene this effect becomes relevant only beyond $T\gtrsim 200\text{K}$. It is also evident that Eq. (14) becomes increasingly accurate and approaches the numerical results as the layer separation $d$ increases.

The behavior of the drag resistivity $\rho_{\text{D}}$ as a function of layer separation $d$ is shown in Fig. 4(c) with $T=5\text{K}$ and $k_{\mathrm{F_{\mathrm{a}}}}=k_{\mathrm{F_{\mathrm{p}}}}=0.5\leavevmode\nobreak\ \mathrm{nm}$. The trend of the exact numerical results changes more slowly, a fact that is also embodied in Fig. 4(a). Interestingly, in drag experiments on graphene [75], the $d$ dependence of the drag resistivity is much slower than the $1/d^{4}$ expected in the weakly interacting regime, varying approximately as $1/d^{2}$ for $d>4\text{nm}$, which is comparable to TIs.

4 Massive Dirac fermions

In a magnetic system the current in the active layer has a longitudinal and a Hall component, the latter being a result of the anomalous Hall effect. The principal mechanisms behind the anomalous Hall effect are: [48, 28] (i) a topological mechanism related to spin precession under the combined action of spin-orbit coupling and the external electric field, which is strongly renormalized by scalar scattering (ii) three spin-dependent scattering mechanisms, which were shown to be of secondary importance in TIs, and will not be considered further. The drag current correspondingly has a longitudinal component and a Hall component, which can be measured separately. If we regard the current in the active layer as giving rise to an effective drag force, it is natural to ask which components of this drag force make the most sizable contribution to the drag current. Indeed, if we focus on the Hall component of the drag current we identify two potential contributions having physically distinct origins. Firstly, the longitudinal current in the active layer can give rise to a transverse drag force in the passive layer, resulting in a Hall current. Secondly, the Hall current in the active layer, stemming from a topological mechanism, can itself directly drag a Hall current in the passive layer, in what may be termed direct Hall drag.

We divide the drag current into four contributions, corresponding to the picture presented in Fig. 2 in the introduction. The electric field $\parallel\hat{\bm{x}}$ gives rise to a longitudinal (${\bm{j}}_{x}$) and an anomalous Hall (${\bm{j}}_{\text{AHE}}$) current in the active layer. In the passive layer there are four contributions to the drag currents: ${\bm{j}}^{\text{ll}}$ is the longitudinal current dragged by ${\bm{j}}_{x}$; ${\bm{j}}^{\text{lh}}$ is the transverse current dragged by ${\bm{j}}_{\text{AHE}}$; ${\bm{j}}^{\text{hl}}$ is the anomalous Hall current generated by ${\bm{j}}^{\text{ll}}$; while ${\bm{j}}^{\text{hh}}$ is the anomalous Hall current generated by ${\bm{j}}^{\text{lh}}$, and it flows longitudinally. In the regime $M_{\text{a,p}}\ll\varepsilon_{\text{F}}$, which is applicable to all samples studied experimentally, the longitudinal drag conductivity $\sigma^{xx}_{D}$ does not depend on the magnetization of either layer, and has the same form as in Eq. (14). Next, the transverse drag ${\bm{j}}^{\mathrm{hl}}$, which is the anomalous Hall current due to ${\bm{j}}^{\mathrm{ll}}$ yields $\rho^{yx}_{\text{D}}\!\approx\!\frac{\sigma^{yx}_{\text{pa}}}{\sigma_{\text{a}}\sigma_{\text{p}}}$. Finally ${\bm{j}}^{\text{lh}}={\bm{j}}^{\text{hh}}=0$. Referring to Fig. 2 this implies that the anomalous Hall current in the (doped) active layer does not generate a corresponding drag current in the passive layer.

The physics can be understood in two steps. Firstly, an external electric field drives the electrons in the active layer longitudinally, and these in turn exert a longitudinal drag force on the passive-layer electrons. The drag force acts as an effective longitudinal driving term for the passive-layer electrons. The anomalous Hall component of the drag current, which at low temperatures can be sizeable compared to the longitudinal component, represents the transverse response of the passive layer to this effective longitudinal driving force, and depends on the passive-layer magnetization $M_{\text{p}}$. The electric field also generates an anomalous Hall current in the active layer. This response is dominated by topological terms of the order of the conductivity quantum, which represents a re-arrangement of charge carriers among spin-momentum locked energy states. Physically this is because Coulomb drag occurs as a result of the interaction between the charge densities in different layers, whereas the anomalous Hall current flowing in the active layer is not associated with a change in the charge density: it does not arise from a shift in the Fermi surface, but from the Berry phase acquired through the rearrangement of carriers among spin-momentum locked states. This implies that the anomalous Hall drag current is quite generally independent of the active-layer magnetization $M_{\text{a}}$. This can be easily verified experimentally.

We analyze the parameter dependence of $\rho^{yx}_{\text{D}}$ in Fig. 5 and Fig. 6. The relationship between $\rho^{yx}_{\text{D}}$ and the magnetization of the passive layer is illustrated in Fig. 5. There is an upward trend at small magnetizations followed by a relatively slow decrease at larger values of $M_{\text{p}}$. The trend can be understood as follows. For $M_{\text{p}}\ll Ak_{\text{F}}$ one may expand ${\bm{j}}^{\text{hl}}$ in $M_{\text{p}}$, which reveals that the current increases nearly linearly with $M_{\text{p}}$. In the opposite limit in which $M_{\text{p}}\gg Ak_{\text{F}}$, [127] the anomalous Hall current vanishes, since the effect of spin-orbit coupling (chirality), which scales with $k_{\text{F}}$, becomes negligible. In fact at large $M_{\text{p}}$ one may expand the band energies in $Ak$, with the leading $k$-dependent term scaling as $k^{2}$, which suggests the system in this limit behaves as a regular, non-magnetic 2DEG. This explains the slow downward trend with increasing $M_{\text{p}}$, and hence the presence of the peak as a function of $M_{\text{p}}$. At larger electron densities the peak occurs at stronger values of the magnetization, since higher electron densities imply higher values of $Ak_{\text{F}}$, increasing the effect of chirality at the expense of the magnetization. TI magnetizations have been measured by superconducting quantum interference device magnetometers.[42, 12] Because the momentum scattering time is in principle known from the longitudinal conductivity of each layer, [128, 129] the trend exhibited by $\rho^{yx}_{\text{D}}$ as a function of the magnetization $M_{\text{p}}$ can be verified experimentally.

We examine the dependence of the anomalous Hall drag resistivity on additional experimentally measurable parameters. In Fig. 6(a) the electron density dependence of $\rho^{yx}_{\text{D}}$ for interlayer separations $d=10,20,40\leavevmode\nobreak\ \text{nm}$ is shown. Compared with the longitudinal drag resistivity, the anomalous Hall drag resistivity has a weaker dependence on electron density. As compared with the longitudinal drag resistivity, the group velocity appearing in the susceptibility of the passive layer is replaced by the Berry curvature, leading to a weaker density dependence, yet no topological contribution. Next, Fig. 6(b) illustrates the temperature dependence of $\rho^{yx}_{\text{D}}$ for separations $d=10,20,40\leavevmode\nobreak\ \text{nm}$. We find that, much like $\rho^{xx}_{\text{D}}$, $\rho^{yx}_{\text{D}}$ also increases nearly quadratically with temperature. The $T^{2}$ dependence stems from the allowed phase space for electron-electron scattering at low temperature, and is expected for any interaction strength between the top and bottom layers of TIs, provided the carriers can be described using a Fermi liquid picture. Moreover, due to the absence of backscattering in TIs there is no correction logarithmic in temperature. Fig. 6(c) presents the layer separation dependence of $\rho^{yx}_{\text{D}}$ for $k_{\text{F}}=0.4,0.5,0.6\text{nm}^{-1}$.

Experimentally, for TI films in the large-surface limit, non-topological contributions from the bulk and the side surface are negligible, [130] and we expect the effects described in this work to be observable. Aside from commonly used materials such as Bi2Se3 and Bi2Te3, a small-gap three-dimensional TI has also been identified in strained HgTe. The TI surface states in this material are, however, spatially extended and could be peaked quite far from the wide quantum well edges, reducing the effective 2D layer separation. [123]

Our theory predicts that in a sample in which the active layer is undoped, so that it produces no longitudinal current but only a quantized anomalous Hall current, there will be no drag current at all in the passive layer. Whereas this is in agreement with Ref. [[131]] its subtleties warrant a separate discussion. In this context, we note that the description of quantum Hall effects can be formulated in terms of a bulk model that accounts for the Berry phase, as well as in terms of an edge state model. The two must evidently agree. In the latter context Ref. [[131]] demonstrated that for chiral edge states in quantum spin-Hall systems this charge imbalance does not translate into a drag current due to the absence of backscattering. Quantum anomalous Hall systems, which have a single edge state, are a subset of quantum spin-Hall systems, which have two chiral edge states. The QAHE can also be understood as a charge imbalance between chiral edge states. Hence one expects the drag current to vanish quite generally at zero doping in quantum anomalous Hall systems. The model used in this work, which is also applicable to quantum spin-Hall systems, can be used to calculate explicitly the edge contribution (Supplement). In complete agreement with the bulk formulation described above, when our model is applied to the edge states it yields zero drag. Yet the edge states are poorly defined at high doping, $\varepsilon_{\text{F}}\tau_{l}/\hbar\gg 1$, which is the focus of this work. In addition, the Dirac fermions studied here live in TIs, which have only one shared edge between the top and bottom surfaces, making it tricky to describe thicker samples using an edge state model. In fact, at zero doping, the TI drag experiment is not well defined, since in that parameter regime leakage into the sidewall states is unavoidable.

Experimentally, for TI films in the large-surface limit, non-topological contributions from the bulk and the side surface are negligible [130], and we expect the effects described in this work to be observable. Aside from commonly used materials such as Bi2Se3 and Bi2Te3, a small-gap three-dimensional TI has also been identified in strained HgTe. The TI surface states in this material are, however, spatially extended and could be peaked quite far from the wide quantum well edges, reducing the effective 2D layer separation [123].

It is well known that the undoped limit of a magnetic topological insulator is a special case. [28] When the chemical potential lies in the middle of the gap opened by the magnetization between the surface conduction and valence bands the anomalous Hall conductivity is quantized: in three dimensions a single TI surface contributes exactly $e^{2}/(2h)$ to the Hall conductivity, referred to as the quantized anomalous Hall effect (QAHE). In a doped system the QAHE effectively yields an offset to the measured anomalous Hall conductivity. In the picture presented here the QAHE emerges naturally as the special case in which the carrier density in the conduction band is taken to zero. The effective two-dimensional model we have used throughout this work predicts that in the special case of the quantum anomalous Hall effect, when the system is undoped, the surface conduction band is empty, and the chemical potential lies in the magnetization gap, the drag current is identically zero. In the previous section we have discussed above the physical interpretation of this result that emerges from our effective two-dimensional picture. It is known, however, that the quantum anomalous Hall effect is associated with a set of chiral edge modes which are well-defined at zero doping (and only then). In this subsection we demonstrate that an effective one-dimensional model for the edge modes at zero doping yields the same result as the two-dimensional model, and gives additional physical insight.

To this end we consider the broader case of Coulomb drag between two identical quantum spin-Hall systems, each with one Kramers pair on its edge. A current $I_{1}$ is driven along the upper edge of the lower quantum spin-Hall system and through electron-electron interactions a voltage $V_{2}$ is induced in the lower edge of the upper quantum spin-Hall system. Each quantum spin-Hall edge state can be described by Hamiltonian $H_{k}=Ak_{x}\sigma_{z}$, and the edge dispersion is $\epsilon_{k}^{\pm}=\pm Ak_{x}$. Hence impurity scattering term $\hat{J}_{0}(f^{(l)}_{\bm{k}})$ becomes

[TABLE]

where ${\bm{k}}=k_{x}\hat{\bm{e}}_{x}$ for 1D-case and $P_{\parallel}\hat{J}(f^{(l)}_{{\bm{k}}\parallel})=0$. In the meantime, the direct inter-layer electron-electron scattering term of Eq. (8) becomes

[TABLE]

If the distribution function is the equilibrium one, $J^{i}(f_{\bm{k}})_{s_{\bm{k}}}=0$. In QSH 1D-case, the wave function overlap $F^{(l)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}}=\langle s_{\bm{k}}l|s_{{\bm{k}}_{1}}l\rangle\langle s_{{\bm{k}}_{1}}l|s_{\bm{k}}l\rangle=1$ with $s_{\bm{k}}=s_{{\bm{k}}_{1}}$ and $F^{(l)}_{s_{\bm{k}}s_{{\bm{k}}_{1}}}=\langle s_{\bm{k}}l|s_{{\bm{k}}_{1}}l\rangle\langle s_{{\bm{k}}_{1}}l|s_{\bm{k}}l\rangle=0$ with $s_{\bm{k}}\neq s_{{\bm{k}}_{1}}$. We found $J^{i}(f_{\bm{k}})_{s_{\bm{k}}}=0$. Because the spectrum is linear there is no contribution to the drag from forward scattering, and backscattering is forbidden by time-reversal symmetry. [131] Hence the Coulomb drag between two identical chiral quantum spin-Hall systems is identically zero. [131]

The argument that chiral edge states do not give rise to drag is general and applies just as well to the Hall conductivity as to the longitudinal conductivity. The key physics is the absence of backscattering. In chiral quantum spin-Hall systems the drag current is zero because backscattering is forbidden by the linear Dirac-like dispersion. In this sense quantum anomalous-Hall systems can be regarded as a special case of quantum spin-Hall systems: they have only one state on each edge, rather than a Kramers pair. So even if backscattering were allowed by the quasiparticle dispersion there would be no states to backscatter into. Hence the drag current in quantum anomalous-Hall systems is identically zero. This result applies also to helical Fermi liquids, such as HgTe quantum wells and thin Bi2Se3 in the 2D limit. Yet the edge states are poorly defined at high doping, $\varepsilon_{\text{F}}\tau_{l}/\hbar\gg 1$, which is the focus of this work. In addition, the Dirac fermions studied here live in TIs, which have only one shared edge between the top and bottom surfaces, making it tricky to describe thicker samples using an edge state model. In fact, at zero doping, the TI drag experiment is not well defined, since in that parameter regime leakage into the sidewall states is unavoidable.

5 Beyond topological insulator films

The focus of the paper up to now has been on thin films of topological band insulators in which both layers are doped with the same type of carrier, that is, either electron-electron or hole-hole layers. Our conclusions can straightforwardly be extended to structures beyond those considered thus far. We will now discuss the necessary modifications for treating the cases of multi-valley TIs and ultra-thin films, the possibility of exciton condensation, and massive Dirac fermions beyond magnetic topological insulators.

Certain materials, such as the topological Kondo insulator SmB6, have more than one valley. In this case, a valley degeneracy factor $g_{v}$ will need to be introduced. In the following, we will fix the electron density $n$. The polarisation $\Pi\propto\sqrt{g_{v}},$ while the susceptibility $\chi\propto g_{v}$, thus the screening function $\epsilon(q,\omega)$ is proportional to $g_{v}$ [132]. The drag current, and hence the drag conductivity, will remain the same as in the single-valley case. For the drag resistivity, based on Eq. (14) we will get the linear $g_{v}$ dependence because $\sigma_{a,p}\propto 1/\sqrt{g_{v}}.$ Note that there will be an intervalley impurity scattering contribution if short-range disorder is present in the system, but this will simply result in a renormalisation of the scattering time tau by an intervalley scattering term.

The case of ultra-thin films deserves special attention. When the TI film is ultra thin tunnelling is enabled between the top and bottom surfaces. The tunnelling between the surface states on the top and bottom surfaces may open an energy gap in the energy spectrum [32]. In this case, the massless Dirac Hamiltonian $H_{0{\bm{k}}}^{(l)}$ needs to be augmented by a series of tunnelling terms, and is generally written as:

[TABLE]

where $\tau$ matrices represent the layer pseudospin space, with $\tau_{z}=1$ symbolising the up surface and $\tau_{z}=-1$ the bottom surface. Here $\hat{\bm{z}}$ is the unit vector in the direction of $\bm{z}$, and the term $t$ represents the tunnelling matrix element between two opposite topological surfaces. After the direct diagonalization, the energy spectrum of the TI thin film is given by $\epsilon_{\bm{k}}=\pm\sqrt{t^{2}+A^{2}k^{2}}$, which has a gap of size $2t$. The disorder scattering term becomes:

[TABLE]

where $\gamma=\theta_{\bm{k}}-\theta_{\bm{k}^{\prime}}$ is the angle between the incident and scattered wave vectors. Note that the density matrix entering this term is the full density matrix $f_{\bm{k}}$ of the double-layer system, rather than its projection onto each individual layer. Two limiting cases can be identified: weak tunnelling $t\ll\varepsilon_{\text{F}}$ and strong tunnelling $t\gg\varepsilon_{\text{F}}$.

For weak tunnelling, $t\ll Ak_{F}$, the carrier wave functions are overwhelmingly located in one of the two layers (surfaces), and the notion of Coulomb drag can be retained to a good approximation. The effect of tunnelling can be taken into account perturbatively. For example, the momentum relaxation time becomes $(1/\tau)\rightarrow(1/\tau)\,[1+t^{2}/(A^{2}k_{\text{F}^{2}})]$. Interlayer tunnelling slightly decreases the momentum relaxation time. Screening is qualitatively different [32], but for weak $t$ it is a good approximation to retain the screening function defined in Eq. (4). We expect, rigorously, $\rho_{\text{D}}\rightarrow\rho_{\text{D}}[1-2t^{2}/(Ak_{F})^{2}]$. The numerics display a similar trend.

For strong interlayer tunnelling $t\sim Ak_{F}$ the notion of Coulomb drag is not applicable to the thin film system. In this case, the carrier wave functions spread over the two layers, hence the picture of the Coulomb interaction causing charges in one layer to drag charges in the other is no longer valid. The main effect of electron-electron interactions will be through the Coulomb renormalization of the conductivity [32]. We note, however, that the film needs to be extremely thin for the tunnelling gap to be noticeable [133], therefore we expect realistic samples to lie in the weak tunnelling limit.

When the active layer is doped with electrons and the passive layer is doped with holes, or vice-versa, exciton condensation may occur. This effect is driven by an exchange term in the interlayer Coulomb interaction. In principle, for interlayer exchange to be nonzero tunnelling also has to be nonzero, but one can think of a situation in which the tunnelling is negligible but the exchange is not. The interesting problem concerning the way exciton condensation impacts $\rho_{xx}$ has been considered in great detail in Ref. [134], where it was shown that the drag resistivity exhibits an upturn at low temperatures described by a logarithmic dependence on the temperature.

The theory of Ref. [134] has recently been shown to be a good description of experimental observations in graphene [135], where an exciton condensate phase has been identified, with the critical temperature estimated at 10 - 100mK. This estimate was for a sample grown on a GaAs substrate, in which the effective dielectric constant is expected to be $\epsilon_{r}\approx 6$, whereas in current TI films the lowest experimentally reported $\epsilon_{r}\approx 30$ [125], largely due to screening by the unavoidable bulk of the film. These observations suggest that the critical temperature in TI films could be at least an order of magnitude smaller than in graphene, placing it in the range 1 - 10 mK, which would make exciton condensation rather difficult to detect experimentally in currently available samples.

Massive Dirac fermions are also found in graphene with a staggered sublattice potential, and $\text{MoS}_{2}$ thin films where inversion symmetry is broken. The band Hamiltonian for a single layer is given by $H^{(l)}_{0{\bm{k}}}=at(\tau k_{x}\sigma_{x}+k_{y}\sigma_{y})+\frac{\Delta}{2}\sigma_{z}$ with $\tau=\pm 1$ the valley index, where $a$ is the lattice constant, $t$ is the effective hopping integral, and $\Delta$ is the energy gap. These Dresselhaus-like Hamiltonians can be directly mapped onto the Rashba Hamiltonian considered in this work [136]. The longitudinal drag current is identical for both TIs and other massive Dirac fermion systems, because the physical mechanism behind the longitudinal drag phenomenon is a result of rectification by the passive layer of the fluctuating electric field generated by the active layer. However, when inversion symmetry is broken in a 2D hexagonal lattice, a pair of valleys which are time-reversal of each other are distinguishable by their opposite values of magnetic moment and Berry curvature. Therefore, there will be no transverse drag current in graphene or $\text{MoS}_{2}$ with broken inversion symmetry because the Berry curvatures have opposite signs in the opposite valleys. In a Dirac semimetal, each Dirac point is four-fold degenerate and can be viewed as consisting of two Weyl nodes with opposite chiralities. Consequently, transverse drag currents also vanish in Dirac semimetals. At small magnetizations the longitudinal drag currents in these materials will be independent of the magnetizations of either layer.

The Hall drag results discussed in this work are also applicable to Rashba 2DEGs, though measuring a strong anomalous Hall effect [137] can be challenging. A sizable fraction of the conductivity quantum is obtained if the two Rashba sub-bands experience a large magnetization splitting and $\varepsilon_{\text{F}}$ lies in the bottom sub-band, but that is challenging experimentally. When $\varepsilon_{\text{F}}\gg M$ the effect vanishes altogether.

6 Conclusions

We have discussed a complete picture of Coulomb drag in Dirac fermion systems incorporating both massless and massive Dirac fermions. For massless Dirac fermions the longitudinal drag resistivity is proportional to $T^{2}d^{-4}n^{-3/2}_{\text{a}}n^{-3/2}_{\text{p}}$ at low temperature and electron densities. We have shown that the drag effect is expected to be weaker in topological insulators than in graphene, and that different regimes are accessible experimentally in these two massless Dirac fermion systems. For massive Dirac fermions we have shown that the drag due to the topological terms on the active surface vanishes and the only contribution to anomalous Hall drag comes from the anomalous Hall current generated by the longitudinal drag force experienced by the passive layer. In a system in which the active layer is undoped and only a quantized anomalous Hall effect exists in the active layer there is no drag current at all. In a doped system the transverse drag current has a non-monotonic dependence on the magnetization of the passive layer, with a peak at a value of the magnetization that becomes pronounced at lower densities.

7 acknowledgments

This research was supported by the Australian Research Council Centre of Excellence in Future Low-Energy Electronics Technologies (project number CE170100039) and funded by the Australian Government. The authors thank Di Xiao, Haizhou Lu, Yaroslav Tserkovnyak, Marco Polini, Greg Fiete, Vladimir Zyuzin, A. H. MacDonald, Oleg Sushkov, Tommy Li, Euyheon Hwang, Zhenyu Zhang, Wenguang Zhu, Zhenhua Qiao, Changgan Zeng, Shun-Qing Shen, and W. K. Tse for enlightening discussions.