Quantum Kinetic Theory of the Chiral Anomaly

TL;DR

This paper develops a comprehensive quantum kinetic theory for magnetotransport in weakly disordered crystals, capturing Berry phase effects and applying it to Weyl semimetals to analyze the chiral anomaly.

Contribution

It introduces a quantum kinetic equation that fully incorporates Berry phase effects and provides a systematic way to analyze magnetotransport phenomena.

Findings

Chiral anomaly persists only with very weak intervalley scattering.

The theory captures Berry phase corrections to the density of states.

Application to Weyl semimetals demonstrates the impact of scattering on magnetoconductivity.

Abstract

We present a general quantum kinetic theory of low-field magnetotransport in weakly disordered crystals that accounts fully for the interplay between electric-field induced interband coherence, Bloch-state scattering, and an external magnetic field. The quantum kinetic equation we derive for the Bloch-state density matrix naturally incorporates the momentum-space Berry phase effects whose influence on Bloch-state wavepacket dynamics is normally incorporated into transport theory in an ad hoc manner. The Berry phase correction to the momentum-space density of states in the presence of an external magnetic field implied by semiclassical wavepacket dynamics is captured by our theory as an intrinsic density-matrix response to a magnetic field. We propose a simple and general procedure for expanding the linear response of the Bloch-state density matrix to an electric field in powers of magnetic field. As an illustration, we apply our theory to magnetotransport in Weyl semimetals. We show that the chiral anomaly (positive magnetoconductivity quadratic in magnetic field) that appears when separate Fermi surface pockets surround distinct Weyl points survives only when intervalley scattering is very weak compared to intravalley scattering.

Tables

Table 1. Table 1: Schematic comparison of our formalism with semiclassical wavepacket dynamics in the presence of a magnetic field 𝑩 𝑩 \bm{B} . Semiclassical wavepacket dynamics implies the Berry phase correction to the density of states, ( 2 ​ π ) − D ​ [ 1 + ( e / ℏ ) ​ 𝑩 ⋅ 𝛀 𝒌 m ] superscript 2 𝜋 𝐷 delimited-[] 1 ⋅ 𝑒 Planck-constant-over-2-pi 𝑩 subscript superscript 𝛀 𝑚 𝒌 (2\pi)^{-D}[1+(e/\hbar)\bm{B}\cdot\bm{\Omega}^{m}_{\bm{k}}] , while the electron distribution function ℱ 𝒌 m superscript subscript ℱ 𝒌 𝑚 \mathcal{F}_{\bm{k}}^{m} remains unchanged Xiao2005 . In our quantum kinetic formalism, however, it is the equilibrium electron density matrix that is modified by a magnetic field, while the momentum-space density of states remains unchanged [see Eqs. ( 41 ) and ( 43 )]. Arrows in this table indicate the changes in physical quantities due to the presence of a magnetic field.

Formalism	Momentum-space density of states	Momentum-space distribution function
Semiclassical wavepacket dynamics Xiao2005	${\left(2\pi \right)}^{-D}\to {\left(2\pi \right)}^{-D}\left[1+\left(e/\hslash \right)𝑩\cdot {𝛀}_{𝒌}^m\right]$	${ℱ}_{𝒌}^m\to {ℱ}_{𝒌}^m$
Our quantum kinetic formalism	${\left(2\pi \right)}^{-D}\to {\left(2\pi \right)}^{-D}$	${ℱ}_{𝒌}^m\to {ℱ}_{𝒌}^m\left[1+\left(e/\hslash \right)𝑩\cdot {𝛀}_{𝒌}^m\right]$

Equations

$$S_{\theta }=\frac{e^2}{4{\pi }^2\hslash c}\int dt{d}^{3}r\theta (\mathbf{r},t)\mathbf{E}\cdot \mathbf{B},$$

$$\frac{\partial N_i}{\partial t}=Q_i\frac{e^2}{4{\pi }^2{\hslash }^{2}c}\mathbf{E}\cdot \mathbf{B},$$

$$Q_i=\int \frac{d^{3}k}{2\pi \hslash }\frac{\partial f_0({\epsilon }_{\mathbf{k}}^m)}{\partial {\epsilon }_{\mathbf{k}}^m}{\mathbf{v}}_{\mathbf{k}}^m\cdot {\mathbf{\Omega }}_{\mathbf{k}}^m$$

$${\sigma }_{zz}^{\mathrm{CA}}(B_z^2)=\mathcal{C}\frac{e^2}{4{\pi }^2\hslash c^2}\frac{(e{B}_z{)}^2{v}_F^3}{{\mu }^2}\tau ,$$

$$H_0|m,\mathbf{k}⟩={\epsilon }_{\mathbf{k}}^m|m,\mathbf{k}⟩,$$

$$H_0=\sum _{\mathbf{k},m}{\epsilon }_{\mathbf{k}}^m|m,\mathbf{k}⟩⟨m,\mathbf{k}|=\sum _{\mathbf{L}{\mathbf{L}}^{\prime }\alpha {\alpha }^{\prime }}{H}_{\mathbf{L}{\mathbf{L}}^{\prime }}^{\alpha {\alpha }^{\prime }}|\alpha ,\mathbf{L}⟩⟨{\alpha }^{\prime },{\mathbf{L}}^{\prime }|$$

$$|\alpha ,\mathbf{L}⟩$$

$$|\alpha ,\mathbf{k}⟩$$

$$⟨\mathbf{r}|\alpha ,\mathbf{k}⟩=e^{i\mathbf{k}\cdot \mathbf{r}}{u}_{\mathbf{k}}^{\alpha }(\mathbf{r}),$$

$$\rho =\sum _{\mathbf{k}{\mathbf{k}}^{\prime }m{m}^{\prime }}{\rho }_{\mathbf{k}{\mathbf{k}}^{\prime }}^{m{m}^{\prime }}|m,\mathbf{k}⟩⟨m^{\prime },{\mathbf{k}}^{\prime }|=\sum _{\mathbf{L}{\mathbf{L}}^{\prime }\alpha {\alpha }^{\prime }}{\rho }_{\mathbf{L}{\mathbf{L}}^{\prime }}^{\alpha {\alpha }^{\prime }}|\alpha ,\mathbf{L}⟩⟨{\alpha }^{\prime },{\mathbf{L}}^{\prime }|,$$

$$\frac{\partial \rho }{\partial t}+\frac{i}{\hslash }[H,\rho ]=0,$$

$$\frac{\partial ⟨\rho ⟩}{\partial t}+\frac{i}{\hslash }[H_0,⟨\rho ⟩]+\frac{i}{\hslash }⟨[U,g]⟩=0,$$

$$\frac{\partial g}{\partial t}+\frac{i}{\hslash }[H_0,g]+\frac{i}{\hslash }[U,g]-\frac{i}{\hslash }⟨[U,g]⟩$$

$$\frac{\partial ⟨\rho ⟩}{\partial t}+\frac{i}{\hslash }[H_0,⟨\rho ⟩]+K(⟨\rho ⟩)=0,$$

$$K(⟨\rho ⟩)=\frac{1}{{\hslash }^2}{\int }_0^{\infty }d{t}^{\prime }⟨\left[U,[e^{-i{H}_0{t}^{\prime }/\hslash }U{e}^{i{H}_0{t}^{\prime }/\hslash },⟨\rho (t)⟩]\right]⟩,$$

$$[I(⟨n⟩){]}_{\mathbf{k}}^{mm}=\frac{2\pi }{\hslash }\sum _{m^{\prime }{\mathbf{k}}^{\prime }}$$

$$[J(⟨n⟩){]}_{\mathbf{k}}^{m{m}^{\prime \prime }}=\frac{\pi }{\hslash }\sum _{m^{\prime }{\mathbf{k}}^{\prime }}⟨U_{\mathbf{k}{\mathbf{k}}^{\prime }}^{m{m}^{\prime }}{U}_{{\mathbf{k}}^{\prime }\mathbf{k}}^{m^{\prime }{m}^{\prime \prime }}⟩[(n_{\mathbf{k}}^m-n_{{\mathbf{k}}^{\prime }}^{m^{\prime }})$$

$$\times \delta ({\epsilon }_{\mathbf{k}}^m-{\epsilon }_{{\mathbf{k}}^{\prime }}^{m^{\prime }})+(n_{\mathbf{k}}^{m^{\prime \prime }}-n_{{\mathbf{k}}^{\prime }}^{m^{\prime }})\delta ({\epsilon }_{\mathbf{k}}^{m^{\prime \prime }}-{\epsilon }_{{\mathbf{k}}^{\prime }}^{m^{\prime }})],$$

$$⟨\rho {⟩}_{\mathbf{p}}^{mn}(\mathbf{r})=\int d^D\mathbf{R}{e}^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}}⟨m,{\mathbf{r}}_{+}|\rho |n,{\mathbf{r}}_{-}⟩,$$

$$⟨m,{\mathbf{r}}_{+}|[H_0,\rho ]|n,{\mathbf{r}}_{-}⟩$$

$$=\int d^3{r}^{\prime }\sum _{m^{\prime }}[⟨m,{\mathbf{r}}_{+}|H_0|m^{\prime },{\mathbf{r}}^{\prime }⟩⟨m^{\prime },{\mathbf{r}}^{\prime }|\rho |n,{\mathbf{r}}_{-}⟩$$

$$-⟨m,{\mathbf{r}}_{+}|\rho |m^{\prime },{\mathbf{r}}^{\prime }⟩⟨m^{\prime },{\mathbf{r}}^{\prime }|H_0|n,{\mathbf{r}}_{-}⟩].$$

$$⟨m,\mathbf{q}|H_0({\mathbf{p}}_{+}+{\mathbf{A}}_{+})|m^{\prime },{\mathbf{q}}^{\prime }⟩$$

$$=⟨m,\mathbf{q}|H_0(\mathbf{q})|m^{\prime },{\mathbf{q}}^{\prime }⟩-\frac{1}{2}i\hslash ⟨m,\mathbf{q}|{\nabla }_{\mathbf{q}}{H}_0(\mathbf{q})\cdot \nabla |m^{\prime },{\mathbf{q}}^{\prime }⟩$$

$$+\frac{1}{2}e⟨m,\mathbf{q}|{\nabla }_{\mathbf{q}}{H}_0(\mathbf{q})\cdot \left[\left(\mathbf{R}\cdot \nabla \right)\mathbf{A}\right]|m^{\prime },{\mathbf{q}}^{\prime }⟩$$

$$\int d^D\mathbf{R}{e}^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}}⟨m,{\mathbf{r}}_{+}|\left\{\frac{\partial \rho }{\partial t}+\frac{i}{\hslash }[H_0,\rho ]\right\}|n,{\mathbf{r}}_{-}⟩=0.$$

$$e^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}}\nabla =\left\{\nabla -(ie/\hslash )[\nabla (\mathbf{A}\cdot \mathbf{R})]\right\}{e}^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}},$$

$$e^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}}\partial /\partial t=\left[\partial /\partial t+(ie/\hslash )\mathbf{E}\cdot \mathbf{R}\right]e^{-(i/\hslash )\mathbf{P}\cdot \mathbf{R}},$$

$$\nabla (\mathbf{A}\cdot \mathbf{R})-\left(\mathbf{R}\cdot \nabla \right)\mathbf{A}=\mathbf{R}\times \mathbf{B},$$

$$⟨m,\mathbf{q}|{\nabla }_{\mathbf{q}}{H}_0(\mathbf{q})|m^{\prime },{\mathbf{q}}^{\prime }⟩=\delta (\mathbf{q}-{\mathbf{q}}^{\prime })\times$$

$$\left[{\nabla }_{\mathbf{q}}{\mathcal{H}}_{0\mathbf{q}}^{m{m}^{\prime }}+⟨u_{\mathbf{q}}^m|{\nabla }_{\mathbf{q}}{u}_{\mathbf{q}}^{n^{\prime }}⟩{\mathcal{H}}_{0\mathbf{q}}^{n^{\prime }{m}^{\prime }}-{\mathcal{H}}_{0\mathbf{q}}^{m{n}^{\prime }}⟨u_{\mathbf{q}}^{n^{\prime }}|{\nabla }_{\mathbf{q}}{u}_{\mathbf{q}}^{m^{\prime }}⟩\right],$$

$$⟨n^{\prime },{\mathbf{q}}^{\prime \prime }|\mathbf{R}|m^{\prime },{\mathbf{q}}^{\prime }⟩⟨\rho {⟩}_{\mathbf{q}}^{m^{\prime }n}=i\hslash \delta ({\mathbf{q}}^{\prime \prime }-{\mathbf{q}}^{\prime })\delta ({\mathbf{q}}^{\prime }-\mathbf{q})\times$$

$$\left[{\nabla }_{\mathbf{q}}⟨\rho {⟩}_{\mathbf{q}}^{n^{\prime }n}+⟨u_{\mathbf{q}}^{n^{\prime }}|{\nabla }_{\mathbf{q}}{u}_{\mathbf{q}}^{m^{\prime }}⟩⟨\rho {⟩}_{\mathbf{q}}^{m^{\prime }n}-⟨\rho {⟩}_{\mathbf{q}}^{n^{\prime }{m}^{\prime }}⟨u_{\mathbf{q}}^{m^{\prime }}|{\nabla }_{\mathbf{q}}{u}_{\mathbf{q}}^n⟩\right],$$

$$\frac{\partial ⟨\rho ⟩}{\partial t}+\frac{i}{\hslash }[{\mathcal{H}}_0,⟨\rho ⟩]+\frac{1}{2\hslash }\left\{\frac{D{\mathcal{H}}_0}{D\mathbf{k}}\cdot \nabla ⟨\rho ⟩\right\}+K(⟨\rho ⟩)$$

$$=D_E(⟨\rho ⟩)+D_B(⟨\rho ⟩),$$

Full text

I Introduction

The electrical transport properties of weakly disordered quantum degenerate crystalline conductors are normally controlled by occupation probabilities of states close to the Fermi surface that change in response to electric fields or temperature gradients. Transport properties then depend on the shapes of all Fermi surfaces, on the distributions of group velocities, and the details of the processes that scatter electrons between Bloch states on those surfaces. It has recently become more widely appreciated that intrinsic effects, independent of disorder and related more to Bloch-state wave functions than energies, are sometimes important. In special cases these effects are nearly completely characterized by momentum-space Berry curvatures Xiao2010 , which can be nonzero only in systems that have broken time-reversal symmetry, or broken inversion symmetry, or both. One widely recognized and important example of a transport effect that is often dominated by Berry curvature physics is the anomalous Hall effect Nagaosa2010 , i.e., the Hall effect in the absence of a magnetic field in a crystal in which time-reversal symmetry is broken by magnetic order. This paper is motivated by the need for a practical theory suitable for application to real crystals that can provide a general description of magnetotransport in systems in which Bloch-state wave function properties play an important role.

Large momentum-space Berry curvatures are often associated with nontrivial band topology. For example the anomalous Hall conductivity of a two-dimensional (2D) insulator Thouless1982 ; Haldane1988 is equal to the quantum unit of conductance, $e^{2}/h$, times the sum of the Chern numbers of all occupied bands. The Chern number of a band is a topological index that is simply the integral of its Berry curvature over the full 2D Brillouin zone divided by $2\pi$, and must be an integer. Indeed, recent interest in the topological classification of crystalline matter started with theories of the quantum Hall effect, and was later extended to time-reversal invariant topological insulators Kane2005 ; Kane2005a ; Bernevig2006 ; Hasan2010 ; Qi2011 ; Ando2013 . The classification of topological phases has now been extended to 3D gapless metallic systems referred to as Weyl semimetals Volovik2003 ; Murakami2007 ; Wan2011 or Dirac semimetals Wang2012 ; Young2012 ; Wang2013 ; Yang2014 . Although Weyl semimetals cannot be characterized by bulk topological invariants, nondegenerate band touching points (Weyl points) are topologically stable and can be regarded as monopole momentum-space sources of phase flux Fang2003 ; AHM2004 . Dirac points, which can be viewed as a superposition of two degenerate Weyl points, are not generically stable, but can be stabilized by crystalline symmetries Yang2014 . The experimental identification of Dirac semimetals Liu2014 ; Liu2014a ; Neupane2014 ; Borisenko2014 and Weyl semimetals Lv2015 ; Yang2015 ; Xu2015 ; Lu2015 has motivated a growing effort, with both theoretical and experimental components, aimed at identifying and exploring novel phenomena in these materials. Three dimensional topological semimetals host a variety of transport and magnetotransport properties that are dependent on their unusual Bloch-state wave function properties, some of which we highlight in the following paragraphs, and therefore require the type of transport theory discussed in this paper for a proper theoretical description.

The topologically nontrivial electronic band structures of Weyl and Dirac semimetals can lead to an approximate realization of the chiral anomaly in condensed matter physics. In quantum field theory the chiral anomaly, which occurs only in systems with odd space dimensions, refers to the violation of axial current conservation, $\partial_{\mu}J_{5}^{\mu}\neq 0$ as a combined response to electric and magnetic fields. In the case of Weyl semimetals these response properties are conveniently described in terms of a $\theta$ term added to the system’s electromagnetic action Wilczek1987 :

[TABLE]

where $\theta(\bm{r},t)$ is a coupling coefficient, and $\bm{E}$ and $\bm{B}$ are external electric and magnetic fields. The presence of a $\theta$ term implies a charge current $\bm{j}(\bm{r},t)=\delta S_{\theta}/\delta\bm{A}=(e^{2}/4\pi^{2}\hbar c)[\nabla\theta(\bm{r},t)\times\bm{E}+\dot{\theta}(\bm{r},t)\bm{B}]$ Wilczek1987 . Here, $\bm{A}$ is the electromagnetic field’s vector potential. The electric-field induced current is the anomalous Hall effect, since it is perpendicular to the electric field, and the magnetic-field induced current is referred to as the chiral magnetic effect Fukushima2008 ; Zyuzin2012 ; Son2012 which we discuss below.

Weyl semimetals with two Weyl points are characterized by the following approximate expression $\theta(\bm{r},t)=2(\bm{b}\cdot\bm{r}-b_{0}t)$ Zyuzin2012 ; Son2012 ; Grushin2012 ; Wang2013a ; Goswami2013 ; Burkov2015 . Here, $\bm{b}$ is the vector connecting the distinct Weyl points in momentum space, and $b_{0}$ is the difference between distinct local chemical potentials, if one is somehow established in different regions of momentum space. Such a chemical potential difference is of course absent in equilibrium; the possible existence of an equilibrium chiral magnetic effect has been considered theoretically Vazifeh2013 ; Chang2015 ; Buividovich2015 ; Yamamoto2015 ; Ma2015 ; Zhong2016 ; Zubkov2016 , but can be ruled out in crystalline solids as discussed in Ref. Vazifeh2013 . A related dynamical realization of the chiral magnetic effect was recently proposed for axionic insulators Sekine2016 ; Sekine2016a . A chemical potential difference can be generated by the combined influence of electric and magnetic fields as explained below, and when present is responsible for a related negative magnetoresistance (positive magnetoconductance) that is quadratic in magnetic field Son2013 ; Burkov2014 ; Burkov2015a ; Spivak2016 . A positive magnetoconductance stands in striking contrast to the familiar negative magnetoconductance due to the Lorentz force. Remarkably, negative magnetoresistance has recently been observed in the low magnetic field regime in the Dirac semimetals Na3Bi Xiong2015 and Cd3As2 Li2015 ; Li2016 , and ZrTe5 Li2016a , and in the Weyl semimetals TaAs Huang2015 and TaP Arnold2016 . The peculiar positive magnetoconductance behavior occurs only for parallel electric and magnetic fields, which suggests that it is related in some way to an $\bm{E}\cdot\bm{B}$ contribution to the electromagnetic action. This very specific effect is a particularly attractive target of the magnetotransport theory developed here, which is able to realistically address its partial realization in real materials.

In this paper, we develop a general quantum kinetic theory of low-field magnetotransport in weakly disordered crystals that accounts fully for the electric-field induced interband coherence responsible for Berry phase contributions to wavepacket dynamics and to the momentum-space density-of-states, and at the same time to account for the interplay between Bloch-state scattering and the presence of an external magnetic field. We take the effect of magnetic fields into account using a semiclassical approximation that we expect to be accurate when the magnetic field is weak enough that Landau quantization can be neglected. Our main result is a quantum kinetic equation which includes driving terms associated with both external electric and magnetic fields, and which we expect to be valuable for understanding the interesting magneto-transport anomalies present in any conductor that has band crossings at energies close to the Fermi energy.

Our paper is organized as follows. In Sec. II we provide a brief overview and physical explanation of our main results, neglecting all technical details. In Sec. III we derive a general quantum kinetic equation that accounts fully for momentum-dependent Bloch-state wave functions and the associated interband coherence response in the presence of external electric and magnetic fields. In Sec. IV we describe a general scheme to apply our theory to magnetotransport phenomena, by performing a systematic low-field expansion. In Sec. V we explain in detail when the intervalley scattering rate emerges as an important parameter in transport properties of many-valley electronic systems, and when it does not. In Sec. VI, as an application of our formalism, we explain how our formalism quite generally captures the pumping between valleys that occurs in systems with Fermi surface pockets that enclose Weyl points with a net chirality. In Sec. VII we apply our quantum kinetic equation to a commonly employed toy model of a two-node Weyl semimetal, highlighting the complicated interplay between disorder, free-evolution, and driving terms in the kinetic equation. In Sec. VIII we focus on the positive magnetoconductance quadratic in magnetic field induced by the chiral anomaly in Weyl and Dirac metals. Finally in Sec. IX we discuss our results and comment on some other potentially interesting applications of our transport theory.

II Chiral Magnetotransport Anomaly

The anomalous positive magnetoconductance of Weyl semimetals was predicted theoretically by observing that when Berry-phase density-of-states and anomalous-velocity corrections are included, the theory of Bloch state wavepacket dynamics predicts that electrons will be pumped between Fermi surface pockets when $\bm{E}\cdot\bm{B}\neq 0$. Under disorder-free semiclassical dynamics, the number of electrons in a given valley changes at the rate Son2012 ; Son2013

[TABLE]

where $i$ is a valley label and

[TABLE]

is the chirality of the associated band. Here, $\bm{\Omega}^{m}_{\bm{k}}$ is the momentum-space Berry curvature of band $m$ (a band that possesses a Fermi surface), $\bm{v}^{m}_{\bm{k}}$ is its Bloch-state group velocity, and $f_{0}(\varepsilon_{\bm{k}}^{m})$ is the Fermi-Dirac distribution function whose derivative with respect to energy $\varepsilon^{m}_{\bm{k}}$ provides a $\delta$-function at the Fermi energy $E_{F}$. It is easy to show that the chirality of each valley is an integer equal to the sum of the chiralities AHM2004 of all the Weyl (band crossing) points enclosed by its Fermi surface, and that the sum of chirality over all Fermi surfaces vanishes.

The goal of our quantum kinetic theory is to capture all effects that are captured as Berry phase corrections in semiclassical wave-packet dynamics, and at the same time to account for disorder in a systematic fashion. The chiral-anomaly-related magnetotransport anomalies of Weyl semi-metals are therefore an excellent test case for the theory. We find that the chiral magnetotransport anomaly survives disorder in the limit that intervalley scattering is very weak compared to intravalley scattering. In the following we take the direction of the parallel electric and magnetic fields to be the $z$-direction.

In a solid, magnetic and electric fields perturb electrons in qualitatively different ways since equilibrium can be established in the presence of the former, but not in the presence of the later. When disorder-free semiclassical wavepacket dynamics is used to derive Eq. (2), the origin of the term on the right hand side proportional to $\Omega_{z}$ is the change in the equilibrium electronic density of states induced by the magnetic field, whereas the terms proportional to $\Omega_{x}$ and $\Omega_{y}$ originate from the anomalous velocity induced by the electric field. We show below that our quantum kinetic theory is able to completely account for both effects.

The robustness of the intervalley pumping effect in the presence of disorder is most easily understood using the Landau level picture illustrated in Fig. 1, in which the influence of the magnetic field is accounted for by solving the equilibrium quantum problem. Note that at nonzero magnetic field each valley carries a nonzero current in equilibrium because the number of Landau levels with positive and negative velocities along the field direction in an individual valley differ, but this does not lead to important observable effects. The intervalley pumping effect central to the chiral anomaly is present only when the electric field is added. In the semiclassical theory, the $\Omega_{z}$ contribution to intervalley pumping accounts for the influence of Landau quantization on the equilibrium density of states, whereas the $\Omega_{x}$ and $\Omega_{y}$ corrections account for the changing character of Landau level wave functions as they are swept through momentum space by the electric field and appear as anomalous velocity contribution to the Lorentz force. Although it is not possible to establish equilibrium in the presence of an electric field, in the presence of disorder scattering it is normally possible to establish a steady state with time-independent values of currents and other observables that can be measured experimentally. The goal of the transport theory we present is to describe such a transport steady state and thereby to make connection with experimental observables while allowing for the possibility that momentum space Berry curvature effects survive disorder. For the chiral anomaly, for example, it is important to establish the connection between free Bloch-state dynamics and the transport steady state, and to distinguish the roles of intervalley and intravalley scattering. We will show that the chiral magnetotransport anomaly is observable only when intervalley scattering at the Fermi energy is very weak compared to intravalley scattering. For a microscopic toy model of Weyl semimetals we compute the positive quadratic magnetoconductivity induced by the chiral anomaly $\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})$ in a fully systematic way. In the low temperature limit we find that the contribution to the magnetoconductivity from a single isotropic Weyl cone is

[TABLE]

where $\mu$ is the chemical potential, $v_{F}$ is the Fermi velocity, $\tau$ is the intervalley scattering time, and $\mathcal{C}\sim 1$ is a coefficient that depends on electronic structure and disorder details. In the case of isotropic Weyl points we find that $\mathcal{C}=1$ in the limit of extremely weak intervalley scattering, independent of the details of the disorder scattering.

III Quantum Kinetic Equation

In this section, we derive a quantum kinetic equation for Bloch electrons in the presence of electric and magnetic fields. To this end, we start with some general considerations related to Bloch Hamiltonians and basis function choices. Throughout this paper, we work in the basis of the disorder-free Hamiltonian eigenstates, which we refer to as the eigenstate basis:

[TABLE]

where $H_{0}$ is the crystal Hamiltonian, $\varepsilon^{m}_{\bm{k}}$ is an eigenvalue of $H_{0}$, $\bm{k}$ is a momentum in the crystal’s Brillouin-zone, and $m$ is a band index. The eigenstates of a crystal Hamiltonian are Bloch states, products of plane-waves with wave-vector $\bm{k}$ and periodic functions that are eigenfunctions of the $\bm{k}$-dependent $\bm{k}\cdot\bm{p}$ Hamiltonian. A Bloch state does not, of course, have to be an eigenstate of the crystal Hamiltonian. Our formalism assumes that a good approximation to the Hamiltonian matrix of the perfect crystal is known in a representation of $\bm{k}$-independent periodic functions. In the following we assume that the Hamiltonian is known in a representation of Wannier functions, but the formalism could be applied with little change if we used a representation of $\bm{k}\cdot\bm{p}$ eigenstates at a particular Brillouin-zone point of interest.

Given this starting point, we write

[TABLE]

with $H^{\alpha\alpha^{\prime}}_{\bm{L}\bm{L}^{\prime}}=\langle\alpha,\bm{L}|H_{0}|\alpha^{\prime},\bm{L}^{\prime}\rangle$. Here $|\alpha,\bm{L}\rangle$ is a Wannier function associated with orbital $\alpha$ and real-space lattice vector $\bm{L}$. For a particular orbital the Bloch and Wannier functions are related by

[TABLE]

where $N$ is the number of Bravais lattice sites in the crystal. $|\alpha,\bm{k}\rangle$ is a Bloch state with wave function

[TABLE]

where $\bm{r}$ is position and $u^{\alpha}_{\bm{k}}(\bm{r})$ is the cell-periodic part of the Bloch wave function. Our transport theory is formulated in terms of the single-particle density matrix. The density-matrix operator can be expressed in either the eigenstate or Wannier representations using

[TABLE]

where $\rho^{mm^{\prime}}_{\bm{k}\bm{k}^{\prime}}=\langle m,\bm{k}|\rho|m^{\prime},\bm{k}^{\prime}\rangle$ and $\rho^{\alpha\alpha^{\prime}}_{\bm{L}\bm{L}^{\prime}}=\langle\alpha,\bm{L}|\rho|\alpha^{\prime},\bm{L}^{\prime}\rangle$. Since we focus on circumstances in which translational symmetry is not broken, the nonequilibrium expectation value of the density matrix operator $\rho$ will always be diagonal in Bloch state wavevector. The eigenstate basis and the Wannier basis are related to each other through the unitary transformation related to band eigenvectors in the $\alpha$-representation: $|m,\bm{k}\rangle=\sum_{\alpha}z_{\alpha}^{m}|\alpha,\bm{k}\rangle$ with $z_{\alpha}^{m}=\langle\alpha,\bm{k}|m,\bm{k}\rangle$.

III.1 Introducing disorder

The main result of this paper is a generic quantum kinetic equation that accounts for disorder and electric and magnetic fields. Throughout this paper we take disorder into account within the Born approximation, implicitly assuming therefore that disorder is weak. We follow a procedure developed in an earlier paper Culcer2016 , generalizing it to allow for magnetic fields.

To establish some notation we first consider the case without external electric and magnetic fields. We start with the quantum Liouville equation

[TABLE]

where $\rho$ and $H$ are the density-matrix operator and the Hamiltonian of the system, respectively. In the absence of external fields, the total Hamiltonian of the system is $H=H_{0}+U$ where $U$ is the disorder potential. We decompose the density matrix $\rho$ into two parts, writing $\rho=\langle\rho\rangle+g$, where $\langle\rho\rangle$ is the density matrix averaged over disorder configurations, and $g$ is the deviation from this average. The quantum Liouville equation can then be decomposed into coupled equations for $\langle\rho\rangle$ and $g$:

[TABLE]

[TABLE]

In the Born approximation we can ignore the last two terms on the left hand side of Eq. (12). The equation for $g$, Eq. (12), can then be solved straightforwardly. By substituting the solution for $g$ into Eq. (11), we obtain

[TABLE]

where the scattering term $K(\langle\rho\rangle)$ is given by Culcer2016

[TABLE]

where $\langle\rho(t)\rangle=e^{-iH_{0}t/\hbar}\langle\rho\rangle e^{iH_{0}t/\hbar}$.

In the following, we do not exhibit the explicit time dependence of $\langle\rho(t)\rangle$ in order to simplify the notation. We separate the density matrix into band-diagonal and band-off-diagonal parts because these two types of components behave quite differently using the notation $\langle\rho\rangle=\langle n\rangle+\langle S\rangle$, where $\langle n\rangle$ is diagonal in band index and $\langle S\rangle$ is off-diagonal. The scattering kernel can be separated into four parts which map $\langle n\rangle$ and $\langle S\rangle$ to band-diagonal and band off-diagonal contributions to $\partial\langle\rho(t)\rangle/\partial t$. This separation is discussed at length later in connection with magnetotransport. Focusing on the elastic-scattering case and using the notation of Ref. Culcer2016 we find that the band-diagonal to band-diagonal part of $K(\langle\rho\rangle)$ from $\langle n\rangle$ is

[TABLE]

with $m$ and $m^{\prime}$ being band indices. This is exactly Fermi’s golden rule. Similarly, the band-diagonal to band-off-diagonal part of $K(\langle\rho\rangle)$ from $\langle n\rangle$ is

[TABLE]

where $m\neq m^{\prime\prime}$. The full expression for $K(\langle\rho\rangle)$, including the contributions from the off-diagonal density matrix $\langle S\rangle$, can be found in Ref. Culcer2016 .

III.2 Introducing electric and magnetic fields

Next, we take the effects of magnetic fields into account using a semiclassical approximation that we expect to be accurate when the weak magnetic fields condition $\omega_{c}\tau\ll 1$ is satisfied and Landau quantization can be neglected. Here, $\omega_{c}$ is the cyclotron frequency and $\tau$ is the transport relaxation time discussed at greater length below. In order to derive the kinetic equation for systems in magnetic fields, we apply the Wigner transformation to the quantum Liouville equation (10). We consider a generic single-particle $D$-dimensional Bloch Hamiltonian $H_{0}(\bm{p})$ with momentum operator $\bm{p}=-i\hbar\nabla$. In the presence of a vector potential $\bm{A}(\bm{r},t)$, minimal coupling results in $\bm{p}\rightarrow\bm{p}+e\bm{A}$ and adopt the notation that $e>0$ is the magnitude of the electron charge.

The Wigner distribution function in the presence of a vector potential $\bm{A}$ is defined by Vasko-book

[TABLE]

where $\bm{P}=\bm{p}-e\bm{A}$, $\bm{r}_{\pm}=\bm{r}\pm\bm{R}/2$, and $|n,\bm{r}\rangle=\sum_{m,\bm{k}}|m,\bm{k}\rangle\langle m,\bm{k}|n,\bm{r}\rangle=\sum_{\bm{k}}e^{i\bm{k}\cdot\bm{r}}|n,\bm{k}\rangle$ is the Fourier transform of $|n,\bm{k}\rangle$. Note that the sign in front of $e\bm{A}$ in $\bm{P}$ is different from the one in minimal coupling Vasko-book . We project the quantum Liouville equation (10) onto the $\{|m,\bm{r}\rangle\}$ space. Let us consider the term

[TABLE]

Since the Hamiltonian is diagonal in real space, we have $\langle m,\bm{r}_{1}|H_{0}|m^{\prime},\bm{r}_{2}\rangle=(\langle m^{\prime},\bm{r}_{2}|H_{0}|m,\bm{r}_{1}\rangle)^{*}=\sum_{\bm{q},\bm{q}^{\prime}}e^{-i\bm{q}\cdot\bm{r}_{1}}e^{i\bm{q}^{\prime}\cdot\bm{r}_{2}}\langle m,\bm{q}|H_{0}(\bm{p}_{1}+e\bm{A}_{1})|m^{\prime},\bm{q}^{\prime}\rangle$ with $\bm{p}_{1}=-i\hbar\partial_{\bm{r}_{1}}$. We expand $\bm{p}_{\pm}$ and $\bm{A}_{\pm}$ up to linear order in $\bm{R}$ as $\bm{p}_{\pm}=-i\hbar\left(\frac{1}{2}\nabla\pm\nabla_{\bm{R}}\right)$ and $\bm{A}_{\pm}=\bm{A}\pm\frac{1}{2}\left(\bm{R}\cdot\nabla\right)\bm{A}$, where $\nabla_{\bm{R}}\equiv\partial/\partial\bm{R}$. Finally we obtain

[TABLE]

with $\bm{q}=-i\hbar\nabla_{\bm{R}}+e\bm{A}=\bm{p}+e\bm{A}$. In reaching Eq. (19) we have observed that it follows from the definition of the Wigner distribution function (17) that $\bm{R}=i\hbar\nabla_{\bm{p}}$ and $\bm{p}=-i\hbar\nabla_{\bm{R}}$ .

We perform the Wigner transformation on the quantum Liouville equation (10) as

[TABLE]

Note that the scattering term $K(\langle\rho\rangle)$ is not changed after the Wigner transformation. We use the following identities

[TABLE]

where $\bm{E}=-\partial\bm{A}/\partial t$ is the electric field. Note that we can obtain the same electric-field dependent term as in Eq. (21) by taking the electric field into account using a scalar potential $\phi(\bm{r})=-\bm{E}\cdot\bm{r}$, i.e., by adding $H_{E}=-e\phi(\bm{r})=e\bm{E}\cdot\bm{r}$ to the Hamiltonian as $H_{0}\to H_{0}+H_{E}$ in Eq. (18) Culcer2016 . By multiplying Eq. (19) by $e^{-(i/\hbar)\bm{P}\cdot\bm{R}}$ from the left we see that the second term in the right hand side of Eq. (19) generates an additional term proportional to $\nabla(\bm{A}\cdot\bm{R})$. Finally, the vector potential can be replaced by the magnetic field using the identity

[TABLE]

where $\bm{B}=\nabla\times\bm{A}$ is the magnetic field. Notice that the matrix element of $\nabla_{\bm{q}}H_{0}(\bm{q})$ in Eq. (19) can be written as

[TABLE]

where $\mathcal{H}_{0\bm{q}}^{mm^{\prime}}=\langle m,\bm{q}|H_{0}(\bm{q})|m^{\prime},\bm{q}\rangle=\delta_{mm^{\prime}}\varepsilon^{m}_{\bm{q}}$ and $|u^{m}_{\bm{q}}\rangle$ is the periodic part of the Bloch function. The magnetic-field dependent term in Eq. (19) becomes $-\frac{1}{2}e\sum_{n^{\prime},\bm{q}^{\prime\prime}}\langle m,\bm{q}|\nabla_{\bm{q}}H_{0}(\bm{q})|n^{\prime},\bm{q}^{\prime\prime}\rangle\cdot\langle n^{\prime},\bm{q}^{\prime\prime}|\bm{R}|m^{\prime},\bm{q}^{\prime}\rangle\times\bm{B}$. Here, the term $\langle n^{\prime},\bm{q}^{\prime\prime}|\bm{R}|m^{\prime},\bm{q}^{\prime}\rangle$ (with $\bm{R}=i\hbar\nabla_{\bm{q}}$) acts on $\langle m^{\prime},\bm{r}^{\prime}|\rho|n,\bm{r}_{-}\rangle$ in Eq. (18), which finally results in

[TABLE]

where the third term on the right hand side comes from the Hermiticity of the equation and consistency with the single-band limit.

Combining Eq. (13) with the analysis above, we arrive at the final form for the quantum kinetic equation in the presence of disorder, an electric field $\bm{E}$, and a magnetic field $\bm{B}$:

[TABLE]

where $\bm{k}=\bm{q}/\hbar$ is the crystal wavevector, and the angle brackets in $\langle\rho\rangle$ imply an impurity average. Here and below $\{\bm{a}\cdot\bm{b}\}\equiv\bm{a}\cdot\bm{b}+\bm{b}\cdot\bm{a}$ (with $\bm{a}$ and $\bm{b}$ being vectors) denotes a symmetrized operator product. In Eq. (25) $D_{E}(\langle\rho\rangle)$ and $D_{B}(\langle\rho\rangle)$ are the electric and magnetic driving terms:

[TABLE]

and

[TABLE]

Eqs. (26) and (27) have been simplified by introducing the covariant derivative notation defined by

[TABLE]

where $\bm{\mathcal{R}}_{\bm{k}}=\sum_{a=x,y,z}\mathcal{R}_{\bm{k}}^{a}\bm{e}_{a}$ with $[\mathcal{R}_{\bm{k}}^{a}]^{mn}=i\langle u^{m}_{\bm{k}}|\partial_{k_{a}}u^{n}_{\bm{k}}\rangle$ being the generalized Berry connection. From another point of view, the covariant derivative is simply the derivative evaluated in the $\bm{k}$-independent Wannier state representation.

Note that the covariant $\bm{k}$-derivative acting on the density matrix accounts for changes in the density matrix due to the fact that its band-eigenstate representation elements are $\bm{k}$-dependent and also due to the fact that the band-eigenstates themselves are $\bm{k}$-dependent. The latter contributions capture the momentum-space Berry curvature contributions to semi-classical wavepacket dynamics, but now in a formalism that accounts consistently for the scattering terms that are necessary to establish the transport steady state. The covariant derivative involving $\mathcal{H}_{0}$ in Eq. (27) is defined in the same way, simply replacing $\langle\rho\rangle$ by $\mathcal{H}_{0}$ in Eq. (28). For our formulation of magnetotransport, it will be important that $D_{B}$ is linear in the density matrix. As we see below the $\bm{k}$-dependence of the eigenstates in this derivative play the essential role in capturing the chiral anomaly.

The covariant derivatives reduce to simple derivatives in a spin-independent single-band system, for example in a parabolic band system with $\mathcal{H}_{0}({\bm{k}})=\hbar^{2}\bm{k}^{2}/2m$. In the same limit, Eq. (25) reduces to the usual semiclassical Boltzmann equation:

[TABLE]

where we have defined the velocity $\bm{v}_{\bm{k}}=(1/\hbar)\nabla_{\bm{k}}\mathcal{H}_{0}(\bm{k})$ and $I(f_{\bm{k}})$ is the single-band version of Eq. (15). The quantum kinetic equation (25) we have derived can therefore be understood as a generalization of the simple Boltzmann equation (29) in which the velocity and distribution function scalars are replaced by matrices, the simple derivatives $\nabla_{\bm{k}}$ are replaced by covariant derivatives $D/D\bm{k}$, and scalar products are replaced by symmetrized matrix products $\frac{1}{2}\{\ \ \cdot\ \ \}$. Equation (25) is the principal result of this paper.

IV Magnetotransport Theory

In this section, we describe a general scheme to apply the quantum kinetic equation (25) to phenomena in systems in which the momentum-space Berry connection plays an important role. We start by obtaining the expression for the density matrix in linear response to an electric field in the absence of a magnetic field. Then we incorporate the effect of a magnetic field by performing a systematic low-magnetic-field expansion. Throughout this paper, we assume that the system is uniform, i.e., that $\nabla\langle\rho\rangle=0$.

IV.1 Transport at zero magnetic field

We first briefly summarize the $\bm{B}=0$ limit of our transport theory. In linear response, we write the electron density matrix $\langle\rho\rangle$ as $\langle\rho\rangle=\langle\rho_{0}\rangle+\langle\rho_{E}\rangle$, where $\langle\rho_{0}\rangle$ is the equilibrium density matrix and $\langle\rho_{E}\rangle$ is the correction to $\langle\rho_{0}\rangle$ in linear order in an electric field $\bm{E}$. With this notation we need to solve the kinetic equation in the form

[TABLE]

We divide the electron density matrix response $\langle\rho_{E}\rangle$ into a diagonal part $\langle n_{E}\rangle$ and an off-diagonal part $\langle S_{E}\rangle$, writing $\langle\rho_{E}\rangle=\langle n_{E}\rangle+\langle S_{E}\rangle$. Note that the equilibrium density matrix $\langle\rho_{0}\rangle$ is diagonal in the band index.

When only band-diagonal to band-diagonal terms are included in the scattering kernel, it is easy to solve for the steady-state value of $\langle n_{E}\rangle$. The kinetic equation (30) in this limit is

[TABLE]

where $\bm{v}^{m}_{\bm{k}}=(1/\hbar)\nabla_{\bm{k}}\varepsilon^{m}_{\bm{k}}$ and $\langle\rho_{0}\rangle^{mm}=f_{0}(\varepsilon^{m}_{\bm{k}})$ is the Fermi-Dirac distribution function. The equation for $n_{E}$ is therefore a familiar linear integral equation which is discussed at greater length below and yields $\langle n_{E}\rangle^{m}_{\bm{k}}=e\tau^{m}_{\mathrm{tr}\bm{k}}\bm{E}\cdot\bm{v}^{m}_{\bm{k}}\partial f_{0}(\varepsilon^{m}_{\bm{k}})/\partial\varepsilon^{m}_{\bm{k}}$, where $\tau^{m}_{\mathrm{tr}\bm{k}}$ is the transport lifetime which is often nearly constant across the Fermi surface.

Next we consider the solution for the off-diagonal part of the density matrix $\langle S_{E}\rangle$, which is independent of weak disorder. From Eq. (30) the kinetic equation for $\langle S_{E}\rangle$ is given by

[TABLE]

where $D^{\prime}_{E}(\langle\rho_{0}\rangle)$ is the off-diagonal part of the intrinsic driving term:

[TABLE]

The off-diagonal part of the driving term is responsible for the Berry phase contribution to the Hall conductivity of systems with broken time-reversal symmetry in the absence of a magnetic field. As we will see, it also plays an essential role in the anomalous magnetoconductivity response. The solution to this equation is Culcer2016

[TABLE]

where we have not explicitly exhibited the time dependences of $\langle\rho_{0}(t-t^{\prime})\rangle$ and $\langle n_{E}(t-t^{\prime})\rangle$. It can be further expanded in the eigenstate basis by inserting an infinitesimal $e^{-\eta t^{\prime}}$ and taking the limit $\eta\rightarrow 0$ to obtain

[TABLE]

Here we have written only the principal value part and omitted $\delta$-function terms which can be important when bands touch, giving rise for example to the Zitterbewegung contribution to the minimum conductivity Culcer2016 in graphene. We note that, as shown in Ref. Culcer2016 , the contribution from $J(\langle n_{E}\rangle)$ corresponds to the vertex correction in the ladder-diagram approximation of perturbation theory.

IV.2 Equilibrium response to a magnetic field

Next we look closely at the magnetic driving term (27) by examining linear response to magnetic field in the absence of an electric field. In this $\bm{E}=0$ limit, we expect that the steady state is actually an equilibrium state, and not one in which energy transfer from the external field is balanced by dissipation. We write the electron density matrix as $\langle\rho\rangle=\langle\rho_{0}\rangle+\langle\rho_{B}\rangle$, where $\langle\rho_{0}\rangle$ is the $\bm{B}=0$ density matrix, and $\langle\rho_{B}\rangle$ is the correction to $\langle\rho_{0}\rangle$ at linear order in the magnetic field $\bm{B}$. This equilibrium density matrix $\langle\rho_{B}\rangle$ satisfies

[TABLE]

To first order in magnetic field the right-hand side can be replaced by $D_{B}(\langle\rho_{0}\rangle)$, where $\langle\rho_{0}\rangle$ is the Fermi-Dirac distribution function of Bloch states in the absence of a magnetic field. For $\bm{B}=(0,0,B_{z})$ we have

[TABLE]

where $\{\ \ ,\ \ \}$ indicates a matrix anticommutator. We will show in Sec. VII.3 that $D_{B}\left(\langle\rho_{0}\rangle\right)$ describes valley-dependent electrical equilibrium currents induced by a magnetic field (i.e., the chiral magnetic effect) that are apparent in the Landau level representation illustrated in Fig. 1. These currents cancel when summed over valleys, but will not cancel when the chemical potentials of the two valleys differ.

We separate the electron density matrix response $\langle\rho_{B}\rangle$ into its diagonal $\langle\xi_{B}\rangle$ and off-diagonal $\langle S_{B}\rangle$ parts: $\langle\rho_{B}\rangle=\langle\xi_{B}\rangle+\langle S_{B}\rangle$. By writing the Hamiltonian and the $\bm{B}=0$ density matrix operators as $\mathcal{H}_{0}=\sum_{m,\bm{k}}\varepsilon_{\bm{k}}^{m}|m,\bm{k}\rangle\langle m,\bm{k}|$ and $\langle\rho_{0}\rangle=\sum_{m,\bm{k}}f_{0}(\varepsilon_{\bm{k}}^{m})|m,\bm{k}\rangle\langle m,\bm{k}|$, where $\varepsilon_{\bm{k}}^{m}$ is an energy eigenvalue of band $m$ with momentum $\bm{k}$ and $f_{0}(\varepsilon_{\bm{k}}^{m})$ is the Fermi-Dirac distribution function, it can be shown quite generally that only band-off-diagonal components in $D_{B}\left(\langle\rho_{0}\rangle\right)$ are nonzero, and hence that for weak disorder that $K(\langle\rho_{B}\rangle)=0$. (See Appendix A for a detailed discussion.) This means that unlike an electric field [Eq. (31)], a magnetic field does not generate a dissipative current when applied to a conductor, as we know. However both fields do alter the band off-diagonal part of the density matrix. We find that the solution for the kinetic equation (36) is

[TABLE]

where $m\neq m^{\prime}$. See Eq. (35) for the corresponding electric field equation.

There is another contribution to the equilibrium density matrix linear in magnetic field in the absence of an electric field. To calculate this, it is convenient to introduce an operator $P$ for an arbitrary density matrix $\langle\rho\rangle$ as

[TABLE]

Note that, in the eigenstate representation, the matrix $P$ is purely diagonal in both wavevector and in density-matrix element at a given wavevector, and that it is nonzero only for off-diagonal density-matrix elements. The notation $P$ for the Bloch-state evolution term $\partial\langle\rho\rangle/\partial t$ in the absence of fields and disorder is intended to suggest its role as the many-band generalization of the spin-precession terms in two-band models which have only a spin-$1/2$ degree of freedom for each momentum. Using Eq. (39) we can rewrite the kinetic equation (36) as

[TABLE]

which can be viewed as a quantum operator.

As mentioned above, $[D_{B}(\langle\rho_{0}\rangle)]^{mm}_{\bm{k}}=0$. However, it can also be shown quite generally that the band-diagonal components of $P^{-1}D_{B}(\langle\rho_{0}\rangle)$ are not zero, i.e., that $[P^{-1}D_{B}(\langle\rho_{0}\rangle)]^{mm}_{\bm{k}}\neq 0$. Here, the operator $P^{-1}$ acting on a driving term $D$ is defined by $[P^{-1}D]^{mm^{\prime}}=-i\hbar D^{mm^{\prime}}_{\bm{k}}/(\varepsilon^{m}_{\bm{k}}-\varepsilon^{m^{\prime}}_{\bm{k}})$. See Appendix A for a detailed derivation. We find that in the weak disorder ($W\tau_{\mathrm{tr}}\gg\hbar$ where $W$ is the scale of the typical energetic separation between Bloch bands)

[TABLE]

where $\Omega_{\bm{k},a}^{m}=\epsilon^{abc}\,i\langle\partial_{k_{b}}u_{\bm{k}}^{m}|\partial_{k_{c}}u_{\bm{k}}^{m}\rangle$ is the Berry curvature of band $m$. Note that the magnetic-field induced change in the diagonal density matrix $\langle\xi_{B}\rangle$ is quite distinct in character from the changes induced by an electric field $\langle n_{E}\rangle$ (31) because $\langle\xi_{B}\rangle$ is intrinsic (independent of disorder effects) and not related to dissipation. Equation (41) captures the same physics as the Berry phase correction to the density of states implied by semiclassical wavepacket dynamics, $(2\pi)^{-D}[1+(e/\hbar)\bm{B}\cdot\bm{\Omega}^{m}_{\bm{k}}]$ Xiao2005 , and can be viewed as an alternate derivation of that effect. In our quantum kinetic formalism, however, it is the equilibrium electron density matrix that is modified by a magnetic field, while the number of states per momentum space volume is unchanged. From Eq. (41) we can immediately obtain the Středa formula for the quantum Hall effect with a magnetic field along the $z$ direction Streda1982 :

[TABLE]

Finally we note that the same result is obtained when $D_{B}$ acts on any band diagonal density matrix $\langle n\rangle=\sum_{m,\bm{k}}\mathcal{F}_{\bm{k}}^{m}|m,\bm{k}\rangle\langle m,\bm{k}|$. The calculation details do not appeal to special properties of the equilibrium density matrix. In the general case we find that

[TABLE]

This property will play an important role in evaluating contributions to the magnetoconductivity (48). A schematic comparison between our formalism and semiclassical wavepacket dynamics is given in Table 1.

IV.3 Magnetotransport

Now we consider a density matrix in the presence of electric and magnetic fields. We write the electron density matrix as $\langle\rho\rangle=\langle\rho_{0}\rangle+\langle\rho\rangle_{F}$, where $\langle\rho_{0}\rangle$ is the density matrix in the absence of fields, and $\langle\rho\rangle_{F}$ is the field-induced density matrix. Then we can rewrite the steady-state uniform system limit of Eq. (25) in the form

[TABLE]

where we have defined an operator $\mathcal{L}\equiv P+K$ and have used the fact that $\mathcal{L}\langle\rho_{0}\rangle=0$. It follows that

[TABLE]

We then view the four terms $P$, $K$, $D_{E}$, and $D_{B}$ as matrices that act on vectors formed by all eigenstate-representation density-matrix components at all wavevectors, both band diagonal and band off-diagonal, ordered in any convenient way. Note also that the matrix $P$ is purely diagonal in both wavevector and in density-matrix element at a given wavevector, and that it is nonzero only for off-diagonal density-matrix elements.

Equation (45) can be used to derive a low-magnetic-field expansion for the linear response of the density-matrix, and hence of any single-particle observable, to an electric field $\bm{E}$. From Eq. (45) we have

[TABLE]

Here, the $N=N^{\prime}=0$ term is given by Eqs. (31) and (35). It is convenient to define a density matrix induced solely by the magnetic field as $\langle\rho_{B}\rangle\equiv\sum_{N\geq 1}(\mathcal{L}^{-1}D_{B})^{N}\langle\rho_{0}\rangle$. Note that this expression is the generalized solution for Eq. (36), i.e., for $(\mathcal{L}-D_{B})\langle\rho_{B}\rangle=D_{B}(\langle\rho_{0}\rangle)$. With this expression for $\langle\rho_{B}\rangle$ we have a simple expression for the low-magnetic-field expansion of the density matrix as

[TABLE]

At each order in $B$, contributions to $\langle\rho\rangle_{F}$ can quite generally be organized by their order in an expansion in powers of scattering strength $\lambda$ by letting $K\to\lambda K$ and identifying terms with a particular power of $\lambda$. The various low-field expansion terms are generated by repeated action of $D_{B}$ and $\mathcal{L}^{-1}$. In Sec. V we will take a close look at the properties of $\mathcal{L}^{-1}$, pointing out that even in systems with weak disorder, the terms that have the highest power of $\lambda^{-1}$ do not necessarily dominate. For states near the Fermi level, the low-field expansion is an expansion in powers of $\omega_{c}\tau$, where $\omega_{c}$ is the cyclotron frequency and $\tau$ is the intravalley scattering time discussed below.

We are now in a position to obtain the magnetoconductivity of a system, i.e., the electrical conductivity in the presence of a magnetic field. Since we are assuming that the magnetic field $\bm{B}$ is very weak, we may set $\langle\rho_{B}\rangle=\langle\xi_{B}\rangle$ in Eq. (47), which means that the correction to the Fermi-Dirac distribution function due to $\bm{B}$ is given by the Berry phase correction $\langle\xi_{B}\rangle$ (41). Then with the use of Eq. (47) the total conductivity $\sigma_{\mu\nu}$ is obtained from the definition $\sigma_{\mu\nu}=\mathrm{Tr}[(-e)v_{\mu}\langle\rho\rangle_{F}]/E_{\nu}=\sum_{N\geq 0}\sigma_{\mu\nu}(B^{N})$ ($\mu,\nu=x,y,z$). Here, $\sigma_{\mu\nu}(B^{N})$ ($N\geq 1$) is the magnetoconductivity contribution proportional to $B^{N}$, which is given by

[TABLE]

Here, $v_{\mu}$ is the $\mu$ component of the velocity operator matrix in the eigenstate representation (i.e., $v_{\mu}=D\mathcal{H}_{0}/Dk_{\mu}$), and $\mathrm{Tr}$ indicates both a matrix trace and an integration with respect to $\bm{k}$ over the Brillouin zone. Note that in Eq. (48) we have allowed for an arbitrary angle between the directions of the electric and magnetic fields. In Sec. VIII we will apply Eq. (48) to the positive magnetoconductivity proportional to $B^{2}$, which arises as a consequence of the chiral anomaly, in Dirac and Weyl semimetals.

V Fermi Surface Response in Multi-valley Systems

Anomalous transport is often related to band crossings at energies that are close to the Fermi energy and give rise to large Berry curvatures. These are often required by symmetry to have replicas at different Brillouin-zone points, for example at points related by time-reversal operations. For that reason, anomalous transport properties of the type that the present formalism is intended to describe often occur in systems with more than one Fermi surface pocket. In what follows, a pocket refers to a closed Fermi surface segment. We will follow the common practice in the literature on Weyl and Dirac semimetals by borrowing terminology from semiconductor physics and referring to these clearly separated regions of momentum space that contribute importantly to the physical properties of interest as valleys. It is common that scattering between valleys related to each other by time-reversal, or by some approximate crystal symmetry, is often weaker than scattering within valleys, giving rise to long intervalley relaxation times. In this section, we explain when these long relaxation times are physically relevant.

V.1 Properties of the operator $\mathcal{L}^{-1}$

The scattering term $K(\langle\rho\rangle)$ [Eq. (14)] may be viewed as a linear operator that acts on the density matrix and is the sum of four terms categorized by the following separation: processes that map diagonal $\langle\rho\rangle$ to diagonal $\partial\langle\rho\rangle/\partial t$ (defined as $K^{dd}$), processes that map diagonal $\langle\rho\rangle$ to off-diagonal $\partial\langle\rho\rangle/\partial t$ ($K^{od}$), processes that map off-diagonal $\langle\rho\rangle$ to diagonal $\partial\langle\rho\rangle/\partial t$ ($K^{do}$) and processes that map off-diagonal $\langle\rho\rangle$ to off-diagonal $\partial\langle\rho\rangle/\partial t$ ($K^{oo}$). For example, with this notation we can rewrite $I(\langle n\rangle)$ and $J(\langle n\rangle)$ as

[TABLE]

with $m\neq m^{\prime}$. Explicit expressions for $K^{dd}$ and $K^{od}$ can be read off from Eqs. (15) and (16). Those for $K^{do}$ and $K^{oo}$ can be found in Ref. Culcer2016 .

When block-diagonalized to separate band-diagonal and band-off-diagonal components of the density matrix, with this notation the Liouville operator

[TABLE]

The disorder ($K$) contributions to $\mathcal{L}$ are proportional to the disorder strength parameter $\lambda$, whereas $P$ is intrinsic and independent of disorder strength. Our expressions for magnetoconductivity are formulated in terms of repeated action of the operator

[TABLE]

The leading terms in the disorder strength expansion of $\mathcal{L}$ are $\sim\lambda^{-1}$. As discussed in Ref. Culcer2016 it follows that up to the order of $\lambda^{0}$

[TABLE]

Note that the ${\it dd}$ block in $\mathcal{L}^{-1}$ is proportional to disorder strength $\lambda^{-1}$ and diverges in the limit of weak disorder, while all other blocks are proportional to $\lambda^{0}$. We see from Eqs. (35) and (38) that the matrix $P^{-1}$ is diagonal in the density-matrix-component wavevector space. When it acts on a driving term $D$, it simply multiplies $D$ by a numerical factor:

[TABLE]

A similar expression applies when $P^{-1}$ acts on any Liouville-space density-matrix vector with nonzero off-diagonal components.

We conclude from these considerations that the leading order term in the disorder strength expansion of the contribution to the response of the diagonal part of the density matrix $\langle n_{N}\rangle$ at $N$-th order in magnetic field strength $B$ is obtained by setting $\mathcal{L}^{-1}$ in Eq. (47) to

[TABLE]

Then it follows that

[TABLE]

where $D^{d}_{B}$ is the diagonal part of the magnetic driving term $D_{B}$, and $\langle\rho_{N-1}\rangle=(\mathcal{L}^{-1}D_{B})^{N-1}\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)$ is the density matrix at $(N-1)$-th order in magnetic field. On the other hand, the leading order term in the disorder strength expansion of the contribution to the response of the off-diagonal part of the density matrix $\langle S_{N}\rangle$ is obtained by setting

[TABLE]

Then it follows that

[TABLE]

where $D^{o}_{B}$ is the off-diagonal part of the magnetic driving term $D_{B}$. Notice from Eq. (49) that $K^{od}\langle n_{N}\rangle=J(\langle n_{N}\rangle)$. Therefore, we see that Eq. (57) is similar to the $B=0$ equation $\langle S_{E}\rangle=P^{-1}[D^{o}_{E}(\langle\rho_{0}\rangle)-K^{od}\langle n_{E}\rangle]$ [see Eq. (35)]. Here, recall that, as shown in Ref. Culcer2016 , the contribution from $J(\langle n_{E}\rangle)=K^{od}\langle n_{E}\rangle$ corresponds to the vertex correction in the ladder-diagram approximation of perturbation theory. The same argument is applied to the case of magnetotransport. Thus, we conclude that the contribution from $K^{od}\langle n_{N}\rangle$ to $\langle S_{N}\rangle$ corresponds to a vertex correction in the ladder-diagram approximation of perturbative calculation. Finally, we note that $\langle S_{N}\rangle$ is accompanied by a diagonal density matrix

[TABLE]

This is the Berry phase correction to the diagonal part of the density matrix, Eq. (43). $\langle n_{N}\rangle$ will be dominant in the case where both $\langle n_{N}\rangle$ and $\langle\xi_{N}\rangle$ are present, since $(K^{dd})^{-1}$ is of the order of $\lambda^{-1}$ while $P^{-1}$ is of the order of $\lambda^{0}$. However, as will be shown in Sec. VIII.1, $\langle\xi_{N}\rangle$ can make a large contribution to the magnetoconductivity in systems with large momentum-space Berry curvatures.

V.2 Scattering times and the total number of electrons in a given valley

We have delayed to this point a discussion of some important properties of $K^{dd}$ and $(K^{dd})^{-1}$ that are related to particle number conservation and play an important role in anomalous magnetotransport. The components of these matrices are labelled by band and wavevector labels. Comparing with Eq. (15) we see that

[TABLE]

It follows that $\sum_{m^{\prime}\bm{k}^{\prime}}K^{dd}[m{\bm{k}};m^{\prime}{\bm{k}}^{\prime}]=0$, i.e.,that the response vector corresponding to total particle number is an eigenvector of $K^{dd}$ with eigenvalue [math]. Scattering does not relax total particle number. Since the leading diagonal response $\langle n_{E}\rangle$ is given as $\langle n_{E}\rangle=(K^{dd})^{-1}D_{E}^{d}$ with $D_{E}^{d}$ the diagonal part of $D_{E}$ [see Eq. (31)], it might seem that the response is divergent. However it is easy to see that the total particle number component of the driving term, i.e.,the driving term $D_{E}^{d}$ summed over all band and wavevector labels is also zero. In fact a stronger statement is valid, namely that the sum of the driving term vanishes when summed over all wavevectors associated with each individual Fermi surface pocket surrounding each relevant valley:

[TABLE]

where $N_{m}$ is the total particle number in band $m$ of a given valley. The electric field moves electrons through momentum space; it does not create or destroy particles. In the matrix equations we have been using to express response to electric fields, the total particle number component of the band-diagonal response has implicitly been projected out.

The eigenvalues of $(K^{dd})^{-1}$ have units of time and characterize the various time scales on which contributions to nonequilibrium populations on the Fermi surface relax. For example, the time scale for relaxation of the nonequilibrium distribution induced by an electric field in a single-valley system, the transport lifetime $\tau_{\mathrm{tr}}$, is given approximately by

[TABLE]

where $D^{m}(\mu)=\sum_{\bm{k}}\delta(\mu-\varepsilon^{m}_{\bm{k}})$ is the density of states at the Fermi energy $\mu$ for valley $m$, and $\theta_{\bm{k}\bm{k}^{\prime}}$ is the angle between the Bloch state group velocities in valley $m$ at $\bm{k}$ and $\bm{k}^{\prime}$ (i.e., $\cos\theta_{\bm{k}\bm{k}^{\prime}}=\bm{k}\cdot\bm{k}^{\prime}/|\bm{k}||\bm{k}^{\prime}|$). In multi-valley systems with weak intervalley scattering it is clear that deviations from equilibrium in the total number of electrons in a valley will relax especially slowly. As we have explained above, these are not driven directly by an electric field, but we will show later that they can be driven by the combined action of electric and magnetic fields.

In order to precisely define intervalley relaxation times explicitly it is necessary to separate the Bloch state Hilbert space into valleys in some precise way. The best way to do this depends on the system being studied. For the sake of definiteness, below we use band labels $m$ to distinguish valleys. When $K^{dd}$ is expressed in the representation in which it is diagonal in the absence of intervalley scattering, and intervalley scattering is parametrically weaker than intravalley scattering, its smallest eigenvalues (longest relaxation times) are associated with an $M\times M$ block with entries

[TABLE]

where $(m,m^{\prime})=1,\ldots,M$ and

[TABLE]

Here, $M$ is the number of Fermi surface pockets, and for simplicity we have assumed that all pockets have the same Fermi-level density of states because they are related by symmetry. Note that this total valley population block of the scattering kernel still has one zero eigenvalue, corresponding to total population summed over all valleys. In the $M=2$ case the nonzero eigenvalue has value $2/\tau^{1,2}_{\mathrm{inter}}$, which is the time scale for relaxing differences in population between the two valleys.

We are now in a position to explain when the magnetoresponse of semimetals with atomically smooth disorder is anomalous, and when it is not, by making the following observations: i) When the Fermi surface of a conductor consists of a few small valleys that are well separated in momentum space and disorder is smooth, scattering between valleys is parametrically weaker because it requires large momentum transfers; ii) When intervalley scattering is negligible, the diagonal-component response function $(K^{dd})^{-1}$ is divergent not only in the total particle number channel, but also in the valley-projected particle number channel; iii) because the electric driving term vanishes separately for each valley, very weak intervalley scattering does not lead to anomalous response at $\bm{B}=0$; but iv) in Weyl and Dirac semimetals $D_{B}\mathcal{L}^{-1}D_{E}(\rho_{0})$ does drive the valley-projected total particle number. The total particle number in a given valley is not conserved when the magnetic driving term summed over a given valley is not zero when it acts on the electric-field disturbed density matrix:

[TABLE]

Here, $\langle\rho_{E}\rangle=\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)$ represents a density matrix linear in electric field. As a consequence of the total particle number nonconservation in a given valley, the intervalley scattering time $\tau_{\mathrm{inter}}$ appears as an eigenvalue of $(K^{dd})^{-1}$.

So far we have focused on systems where valleys are not degenerate, i.e., are separated in momentum space, such as Weyl semimetals. Similar considerations apply to systems where valleys are degenerate, such as Dirac semimetals. 3D Dirac semimetals such as Cd3As2 and Na3Bi have two Dirac points which are protected by crystalline symmetry. The effective Hamiltonian around a Dirac point in such Dirac semimetals can be written in the block-diagonal form Yang2014

[TABLE]

where $A$ and $B$ denote two states after a unitary transformation, and $\mathcal{H}_{\pm}(\bm{q})$ is a $2\times 2$ Weyl Hamiltonian with chirality $\pm 1$. As is seen from Eq. (65), two Weyl Hamiltonians that belong to different states $A$ and $B$ are degenerate around the Dirac point. In other words, two different valleys are degenerate. In general, atomically smooth disorder does not allow scattering processes between two different states $A$ and $B$. Namely, intervalley scattering is much weaker than intravalley scattering, i.e., we have $\tau_{\mathrm{inter}}\gg\tau_{\mathrm{intra}}$. Thus we can apply the same argument as in the case of nondegenerate valleys to the case of degenerate valleys: the total particle number in a given valley is not conserved if the magnetic driving term summed over a given valley is not zero [see Eq. (64)].

We will demonstrate in Sec. VIII that the particle number nonconservation in a given valley characterized by Eq. (64) occurs in the process of calculating the magnetoconductivity proportional to $B^{2}$ in Dirac and Weyl semimetals, which is identified as a consequence of the chiral anomaly.

VI Chiral Anomaly in Three-Dimensional Semimetals

In Weyl semimetals the term “chiral anomaly” is sometimes used to refer to the property that the total number of electrons in a given valley is not conserved in the presence of simultaneous electric and magnetic fields, and sometimes to the enhanced magnetoconductivity that this lack of separate particle number conservation produces. In this section we focus on pumping of charge between valleys without referring to a specific microscopic model, which leads to observable effects only when intervalley scattering times are much longer than intravalley scattering times. For long intervalley scattering times, anomalous pumping leads to a difference between the effective chemical potentials of different valleys and thus to a enhanced current via the magnetoelectric effect, as explained in great deal for a Weyl metal toy model in the following sections. We set $\hbar=1$ in the rest of this paper.

We study a general model with the Hamiltonian $\mathcal{H}_{0}=\sum_{m,\bm{k}}\varepsilon_{\bm{k}}^{m}|m,\bm{k}\rangle\langle m,\bm{k}|$ and the equilibrium density matrix $\langle\rho_{0}\rangle=\sum_{m,\bm{k}}f_{0}(\varepsilon_{\bm{k}}^{m})|m,\bm{k}\rangle\langle m,\bm{k}|$, where $\varepsilon_{\bm{k}}^{m}$ is an energy eigenvalue and $f_{0}(\varepsilon_{\bm{k}}^{m})$ is the Fermi-Dirac distribution function. For concreteness, we choose $\bm{E}=(0,0,E_{z})$ and $\bm{B}=(0,0,B_{z})$. Let us consider the following quantity that is linear in both electric and magnetic fields,

[TABLE]

evaluated for a particular Fermi surface pocket associated with a particular valley. $\partial N/\partial t$ is the rate of change of the density matrix that is balanced by the scattering and free evolution contributions to $\mathcal{L}\langle\rho\rangle$ to yield the part of the density matrix that is linear in both electric and magnetic fields; the steady state density matrix is therefore obtained by acting on $\partial N/\partial t$ with $\mathcal{L}^{-1}$.

Let us consider the first term in the right-hand side of Eq. (66), i.e., $\mathrm{Tr}[D_{B}\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)]$. As we have explained earlier, $\langle\rho_{E}\rangle=\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)$, the linear-response density matrix in the absence of a magnetic field, contains both band-diagonal $\langle n_{E}\rangle$ and band off-diagonal $\langle S_{E}\rangle$ contributions: $\langle\rho_{E}\rangle=\langle n_{E}\rangle+\langle S_{E}\rangle$. It follows that $\mathrm{Tr}\,[D_{B}\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)]=\mathrm{Tr}\,[D_{B}(\langle n_{E}\rangle)]+\mathrm{Tr}\,[D_{B}(\langle S_{E}\rangle)]$, since $D_{B}$ is linear in the density matrix. First we evaluate the diagonal element $[D_{B}(\langle n_{E}\rangle)]^{mm}_{\bm{k}}$ with $\langle n_{E}\rangle=eE_{z}\sum_{m,\bm{k}}[\partial_{k_{z}}f_{0}(\varepsilon_{\bm{k}}^{m})]|m,\bm{k}\rangle\langle m,\bm{k}|$. In this case we can use the result of $[D_{B}(\langle\rho_{0}\rangle)]^{mm}_{\bm{k}}$ by replacing $f_{0}$ by $eE_{z}\partial_{k_{z}}f_{0}$. Namely, from Eq. (141) we get

[TABLE]

It is obvious that this is an odd function of $k_{x}$, $k_{y}$, and $k_{z}$, which means that $\mathrm{Tr}\,[D_{B}(\langle n_{E}\rangle)]=0$. Then we see that one of the two contributions to the nonconservation of the total electron number in a given valley induced by the chiral anomaly is associated with the action of the $D_{B}$ operator on the off-diagonal density matrix $\langle S_{E}\rangle$ as

[TABLE]

where FS represents the integration on the Fermi surface of the valley. After a calculation we obtain the following general expression for the Fermi surface contribution which can be expressed in terms of band velocities and Berry curvatures alone (see Appendix B for a detailed derivation):

[TABLE]

where $\Omega_{\bm{k},a}^{m}=\epsilon^{abc}\,i\langle\partial_{k_{b}}u_{\bm{k}}^{m}|\partial_{k_{c}}u_{\bm{k}}^{m}\rangle$ is the Berry curvature of band $m$. Thus the right-hand side of Eq. (68) is written as

[TABLE]

where $\bm{v}^{m}_{\bm{k}}$ is the Bloch state group velocity, and we have assumed that only the band $m$ intersects the Fermi surface, i.e., $\partial f_{0}(\varepsilon^{n}_{\bm{k}})/\partial\varepsilon^{n}_{\bm{k}}=\delta_{mn}\partial f_{0}(\varepsilon^{m}_{\bm{k}})/\partial\varepsilon^{m}_{\bm{k}}$ , which can in general apply to multi-valley systems.

Next let us consider the second term in the right-hand side of Eq. (66), i.e., $\mathrm{Tr}[D_{E}\mathcal{L}^{-1}D_{B}(\langle\rho_{0}\rangle)]$. As we have explained earlier, $\langle\rho_{B}\rangle=\mathcal{L}^{-1}D_{B}(\langle\rho_{0}\rangle)$, the linear-response density matrix in the absence of an electric field, contains both band-diagonal $\langle\xi_{B}\rangle$ and band off-diagonal $\langle S_{B}\rangle$ contributions: $\langle\rho_{B}\rangle=\langle\xi_{B}\rangle+\langle S_{B}\rangle$. In this case the calculation is easier than that of $\mathrm{Tr}[D_{B}(\langle\rho_{E}\rangle)]$. After a calculation we find that only the contribution from the Berry phase correction $\langle\xi_{B}\rangle$ survives:

[TABLE]

Combining Eqs. (70) and (71) we arrive at the final expression for the rate of pumping of electrons between valleys due to the chiral anomaly

[TABLE]

which is in complete agreement with the expression obtained by disorder-free semiclassical dynamics [Eq. (2)] Son2012 ; Son2013 . When disorder-free semiclassical wavepacket dynamics is used to derive Eq. (2), the origin of the term on the right hand side proportional to $\Omega_{z}$ is quite different from the terms proportional to $\Omega_{x}$ and $\Omega_{y}$. The later are due to the well known anomalous velocity of electrons in a magnetic field, which contributes to the Lorentz force and hence to the rate of motion of states through momentum space. The anomalous velocity effect is an interband wave-function polarization effect, akin to the dielectric screening response of insulators. The $\Omega_{z}$ contribution, on the other hand, is due to the correction of the density of states in momentum space in the presence of a magnetic field. The density-of-states correction is the semiclassical manifestation of the relationship between Berry curvatures (Berry phases) and Landau quantization. We have shown in Eq. (72) that our quantum kinetic theory is able to completely account for both effects.

Here, let us consider the case of Weyl semimetals. As depicted in Fig. 1, a magnetic field induces an anomalous ($N=0$) Landau level branch that has only one sign of velocity $v_{z}$ in a given valley. It follows that in each valley the density of states is increased for states with one sign of velocity and decreased for states with the other sign of velocity, and that the total current summed over a valley is already nonzero in equilibrium. When an electric field drives states through momentum space, the total number of states in a valley varies. Scattering within valleys can relax the current in each valley to its equilibrium value, but cannot establish a steady state because the number of states in each valley still changes at a constant rate. A steady state can be established only when scattering processes between valleys are present. When Bloch-state scattering between valleys is much weaker than scattering within valleys the anomalous $\bm{E}\cdot\bm{B}$ pumping leads to differences between the chemical potentials of different valleys.

VII Chiral Anomaly in Weyl Semimetals

So far we have i) formulated a general theory of low-field that accounts for the wavevector-dependence of Bloch states in multiband systems, ii) shown that the Berry-phase related density of states and anomalous velocities of semiclassical wavepacket dynamics are fully captured by our theory, and iii) demonstrated that Bloch state scattering, which is consistently described in our transport theory, can alter theoretical predictions made by introducing Berry phase effects in semiclassical wavepacket dynamics into the transport theory in an ad hoc manner. In this section, we present an analysis of the magnetotransport properties of a specific toy model of a Weyl semimetal which illustrates in more detail the interplay between intrinsic dynamics, scattering, and external field driving terms in the kinetic equation.

VII.1 Theoretical model

We start with the $4\times 4$ continuum toy model Hamiltonian for two-node Weyl semimetals with broken time-reversal symmetry Burkov2011 ; Vazifeh2013 ; Sekine2014 ; Burkov2015

[TABLE]

where the Pauli matrices $\tau_{i}$ and $\sigma_{i}$ are Weyl-node and spin degrees of freedom, respectively, $v_{F}$ is the Dirac velocity, and $\Delta$ is the mass of 3D Dirac fermions. The term $b\sigma_{z}$ can be regarded as a magnetic interaction of the $z$ component of spin such as $s$-$d$ coupling and Zeeman coupling. Therefore, it can be said that the above Hamiltonian represents a 3D (topological or normal) insulator doped with magnetic impurities. As will be shown just below, a Weyl semimetal is realized when $|b/\Delta|>1$. The time-reversal $\mathcal{T}$ and spatial inversion (parity) $\mathcal{P}$ operators are given by $\mathcal{T}=-i\sigma_{y}\mathcal{C}$ with $\mathcal{C}$ the complex conjugation operator and $\mathcal{P}=\tau_{x}$. Here, we have the following relations: $\mathcal{T}^{-1}=i\sigma_{y}\mathcal{C}$, $\mathcal{P}^{-1}=\tau_{x}$, $\mathcal{T}^{2}=-1$, and $\mathcal{P}^{2}=1$. Then it is easily seen that time-reversal symmetry of the system is broken such that $\mathcal{T}\mathcal{H}_{0}(\bm{k})\mathcal{T}^{-1}\neq\mathcal{H}_{0}(-\bm{k})$, but inversion symmetry of the system is preserved such that $\mathcal{P}\mathcal{H}_{0}(\bm{k})\mathcal{P}^{-1}=\mathcal{H}_{0}(-\bm{k})$.

Performing a canonical transformation such that $\sigma_{x,y}\rightarrow\tau_{z}\sigma_{x,y}$ and $\tau_{x,y}\rightarrow\sigma_{z}\tau_{x,y}$, Eq. (73) can be rewritten as $\tilde{\mathcal{H}}_{0}(\bm{k})=v_{F}(k_{x}\sigma_{x}+k_{y}\sigma_{y})+(b+v_{F}\tau_{z}k_{z}+\Delta\tau_{x})\sigma_{z}$. In the representation of eigenstates of $b+v_{F}\tau_{z}k_{z}+\Delta\tau_{x}$, Eq. (73) is block-diagonal with two $2\times 2$ Hamiltonians given by Burkov2015

[TABLE]

with $m_{\pm}(k_{z})=b\pm\sqrt{v_{F}^{2}k_{z}^{2}+\Delta^{2}}$. The two Weyl nodes are located at $W_{\pm}=(0,0,\pm k_{0})$ with $k_{0}=\sqrt{b^{2}-\Delta^{2}}/v_{F}$. Near the Weyl nodes, we have $m_{-}(q_{z})\approx\mp v_{F}^{2}(k_{0}/b)q_{z}$ with momentum $\bm{q}=\bm{k}-W_{\pm}$ ($\bm{q}^{2}\ll 1$). The eigenvectors of $\mathcal{H}_{t}(\bm{k})$ $(t=\pm)$ with eigenvalues $\varepsilon^{\pm}_{t\bm{k}}=\pm\varepsilon_{t\bm{k}}=\pm\sqrt{v_{F}^{2}(k_{x}^{2}+k_{y}^{2})+m_{t}^{2}}$ are given by

[TABLE]

where $e^{i\theta}=(k_{x}+ik_{y})/k_{\perp}$ with $k_{\perp}=\sqrt{k_{x}^{2}+k_{y}^{2}}$. We see that $\mathcal{H}_{-}(\bm{k})$ describes a subsystem with two Weyl nodes, while $\mathcal{H}_{+}(\bm{k})$ has a fully gapped spectrum. For convenience, we omit the subscripts $t=\pm$ at intermediate steps of the calculations below and restore them in the final results.

The generalized Berry connection in the eigenstate representation is given by $[\mathcal{R}_{\bm{k},\alpha}]^{mn}=i\langle u_{\bm{k}}^{m}|\partial_{k_{\alpha}}u_{\bm{k}}^{n}\rangle$ with $\alpha=x,y,z$ and $m,n=\pm$. The individual components are given explicitly by

[TABLE]

where $\tilde{\sigma}_{\alpha}$ are the Pauli matrices in the eigenstate basis. Also, the individual components of the Berry curvature, $\Omega^{\pm}_{\bm{k},a}=\epsilon^{abc}\,i\langle\partial_{k_{b}}u_{\bm{k}}^{\pm}|\partial_{k_{c}}u_{\bm{k}}^{\pm}\rangle$, are given by

[TABLE]

Here, we briefly review the derivation of the $\theta$ term from the microscopic four-band model (73). In order to describe a more generic Weyl semimetal, we add the term $\mu_{5}\tau_{z}$ to the Hamiltonian, which generates a chemical potential difference $2\mu_{5}$ between the two Weyl nodes, and replace $b\sigma_{z}$ by $\bm{b}\cdot\bm{\sigma}$. Note that this term breaks inversion symmetry. Figure 2 shows a schematic diagram of a two-node Weyl semimetal. We set $\Delta=0$ for simplicity. The action of the system in the presence of an external electromagnetic potential $A_{\mu}=(A_{0},-\bm{A})$ is given by

[TABLE]

where $e>0$, $\psi$ is a four-component spinor, $\bar{\gamma}=\psi^{\dagger}\gamma^{0}$, $\gamma^{0}=\tau_{x}$, $\gamma^{j}=\tau_{x}\tau_{z}\sigma_{j}=-i\tau_{y}\sigma_{j}$, $\gamma^{5}=i\gamma^{0}\gamma^{1}\gamma^{2}\gamma^{3}=\tau_{z}$, and $b_{\mu}=(\mu_{5},-\bm{b})$. After applying the Fujikawa method Fujikawa , i.e., applying an infinitesimal gauge transformation such that $\psi\rightarrow e^{-id\phi\theta(\bm{r},t)\gamma^{5}/2}\psi$ with $\theta(\bm{r},t)=-2x^{\mu}b_{\mu}=2(\bm{b}\cdot\bm{r}-\mu_{5}t)$ and $\phi\in[0,1]$ for infinite times, the action of the system becomes Zyuzin2012 ; Burkov2015

[TABLE]

where $S_{\theta}$ is the $\theta$ term [Eq. (1)] with $\theta(\bm{r},t)=2(\bm{b}\cdot\bm{r}-\mu_{5}t)$. Note that the first term in Eq. (79) represents the (trivial) action of massless Dirac fermions.

In general, regardless of whether the system is gapless or gapped, the four-current density $j^{\nu}$ can be obtained from the variation of the $\theta$ term with respect to the four-potential $A_{\nu}$ as $j^{\nu}=\delta S_{\theta}/\delta A_{\nu}=-(e^{2}/4\pi^{2})\epsilon^{\mu\nu\rho\lambda}[\partial_{\mu}\theta(\bm{r},t)]\partial_{\rho}A_{\lambda}$ Wilczek1987 . The induced current density in the presence of external electric and magnetic fields is given by

[TABLE]

where $\dot{\theta}=\partial\theta(\bm{r},t)/\partial t$, $\bm{E}=-\nabla A_{0}-\partial\bm{A}/\partial t$, and $\bm{B}=\nabla\times\bm{A}$. In the present case of Weyl semimetals with $\theta(\bm{r},t)=2(\bm{b}\cdot\bm{r}-\mu_{5}t)$, we obtain in the ground state a static current of the form

[TABLE]

where the electric-field induced and magnetic-field induced terms are the anomalous Hall effect and chiral magnetic effect, respectively Zyuzin2012 ; Son2012 ; Grushin2012 ; Wang2013a ; Goswami2013 ; Burkov2015 ; Fukushima2008 ; Vazifeh2013 ; Sekine2016 .

The anomalous Hall effect in Weyl semimetals can be understood straightforwardly. It is well known that the Hall conductivity of 2D massive Dirac fermions (such as those on the surface of 3D topological insulators) of the form (74) is given by $\sigma_{xy}^{\pm}(k_{z})=\mathrm{sgn}[m_{\pm}(k_{z})]e^{2}/2h$, when the chemical potential lies in the energy gaps. This (half-)quantized value holds valid even in the presence of disorder Nomura2011 ; Culcer2011a . We see that the 3D Hall conductivity is nonzero only in the region $-k_{0}\leq k_{z}\leq k_{0}$, which gives $\sigma^{\mathrm{3D}}_{xy}=\int_{-k_{0}}^{k_{0}}\,(k_{z}/2\pi)\sigma_{xy}^{\pm}(k_{z})=k_{0}e^{2}/\pi h$. This value is exactly the same as the first term in Eq. (81). On the other hand, the chiral magnetic effect in Weyl semimetals looks like a peculiar phenomenon. The chiral magnetic effect indicates a direct current generation along a static magnetic field (without electric fields) in the presence of a chemical potential difference between the two nodes, as is seen from the second term in Eq. (81). If the static chiral magnetic effect exists in real materials, there will be substantial possible applications. Discussions on the existence of the static chiral magnetic effect in Weyl semimetals have continued theoretically Vazifeh2013 ; Chang2015 ; Buividovich2015 ; Yamamoto2015 ; Ma2015 ; Zhong2016 ; Zubkov2016 . However, the existence of the static chiral magnetic effect would be ruled out in crystalline solids (i.e., lattice systems) as discussed in Ref. Vazifeh2013 , which is consistent with our understanding that static magnetic fields do not generate equilibrium currents.

Lastly, we refer to these two phenomena from the symmetry viewpoint. The term $\bm{b}\cdot\bm{\sigma}$ breaks time-reversal symmetry but preserves inversion symmetry: $\mathcal{T}\sigma_{i}\mathcal{T}^{-1}=-\sigma_{i}$ and $\mathcal{P}\sigma_{i}\mathcal{P}^{-1}=\sigma_{i}$ ($i=x,y,z$). This indicates that the occurrence of the anomalous Hall effect in Weyl semimetals requires the breaking of time-reversal symmetry, which is consistent with the generic requirement for the realization of the anomalous Hall effect. On the other hand, the term $\mu_{5}\tau_{z}$ breaks inversion symmetry but preserves time-reversal symmetry: $\mathcal{T}\tau_{z}\mathcal{T}^{-1}=\tau_{z}$ and $\mathcal{P}\tau_{z}\mathcal{P}^{-1}=-\tau_{z}$. This indicates that the occurrence of the chiral magnetic effect in Weyl semimetals requires the breaking of inversion symmetry. Note that both the anomalous and chiral magnetic effects are possible (at least theoretically) in systems with broken time-reversal and inversion symmetries.

VII.2 Anomalous Hall effect in Weyl metals

Let us see that the electric driving term (26) properly describes the anomalous Hall effect in Weyl metals. As described in Sec. IV, the diagonal part of the electric driving term $D_{E}$ results in usual Drude conductivity. Hence, we need the explicit matrix expression for the off-diagonal part of the electron density $\langle S_{E}\rangle$ to obtain the anomalous Hall conductivity. As is seen from Eq. (35), there are contributions from the electric driving term $D_{E}(\langle\rho_{0}\rangle)$ and the anomalous driving term $J(\langle n_{E}\rangle)$ to the off-diagonal part of the density matrix $\langle S_{E}\rangle$. Let us denote the resultant total anomalous Hall conductivity of the system as

[TABLE]

where $\sigma^{\rm I}_{xy}$ is the Hall conductivity from $D_{E}(\langle\rho_{0}\rangle)$ and $\sigma^{\rm II}_{xy}$ is that from $J(\langle n_{E}\rangle)$. In the following, we consider the case where an electric field is applied along the $y$ direction as $\bm{E}=E_{y}\bm{e}_{y}$, and set $\mu_{5}=0$. Also, we add a positive small chemical potential $\mu$ which lies sufficiently close to the Weyl nodes. Note that the sign of $\mu$ is not essential, since the Hamiltonian (74) has particle-hole symmetry.

First, we evaluate the intrinsic contribution $\sigma^{\rm I}_{xy}$, which originates from the Fermi sea of the electronic band. Using the expression for $D_{E}(\langle\rho_{0}\rangle)$ (33) and the Berry connection (76), the off-diagonal part of the electron density matrix $\langle S_{E}\rangle$ is given by

[TABLE]

where $f_{0}(\varepsilon^{\pm}_{\bm{k}})=1/[e^{(\varepsilon^{\pm}_{\bm{k}}-\mu)/T}+1]$ is the Fermi-Dirac distribution function. The velocity operator is given by $v_{x}=\partial\mathcal{H}_{t}/\partial k_{x}=v_{F}\sigma_{x}$. We work in the eigenstate basis, in which the matrix element $[v_{x}]^{mn}=v_{F}\langle u_{\bm{k}}^{m}|\sigma_{x}|u_{\bm{k}}^{n}\rangle$ ($m,n=\pm$) is given by

[TABLE]

The anomalous Hall conductivity $\sigma^{\rm I}_{xy}$ is calculated from the definition $\sigma^{\rm I}_{xy}=\mathrm{Tr}\left[(-e)v_{x}\langle S_{E}\rangle\right]/E_{y}$. Combining the above ingredients, we obtain

[TABLE]

where we have used the fact that $m_{-}(k_{z})\approx\mp v_{F}^{2}(k_{0}/b)q_{z}$ around the Weyl nodes with $\bm{q}=\bm{k}-W_{\pm}$ ($\bm{q}^{2}\ll 1$). Then we see that the contribution from $\varepsilon^{+}_{t\bm{k}}$ always vanishes as long as the chemical potential $\mu$ lies sufficiently close to the Weyl nodes. Also, it should be noted that the quantity $\Omega_{t\bm{k},z}=v_{F}^{2}m_{t}/2\varepsilon_{t\bm{k}}^{3}$ is the $z$ component of the Berry curvature of the two-band Hamiltonian (74). The anomalous Hall conductivity (85) is exactly the same as the value obtained from a field-theoretical approach, the first term in Eq. (81).

Next, we need to evaluate the extrinsic contribution $\sigma^{\rm II}_{xy}$, which originates from the Fermi surface and the presence of disorder. As described in Sec. III, we treat disorder within the Born approximation. We consider a short-range (on-site) disorder potential of the form $U(\bm{r})=U_{0}\sum_{i}\delta(\bm{r}-\bm{r}_{i})$, and assume that the correlation function satisfies $\langle U(\bm{r})U(\bm{r}^{\prime})\rangle=n_{\mathrm{imp}}U_{0}^{2}\,\delta(\bm{r}-\bm{r}^{\prime})$ with $n_{\mathrm{imp}}$ the impurity density. Then we obtain

[TABLE]

where $\gamma=\theta^{\prime}-\theta$ with $e^{i\theta^{\prime}}=(k^{\prime}_{x}+ik^{\prime}_{y})/k^{\prime}_{\perp}$. We also have the relations $\langle U^{+-}_{\bm{k}\bm{k}^{\prime}}U^{--}_{\bm{k}^{\prime}\bm{k}}\rangle=-\langle U^{++}_{\bm{k}\bm{k}^{\prime}}U^{+-}_{\bm{k}^{\prime}\bm{k}}\rangle$, $\langle U^{--}_{\bm{k}\bm{k}^{\prime}}U^{-+}_{\bm{k}^{\prime}\bm{k}}\rangle=-\langle U^{++}_{\bm{k}\bm{k}^{\prime}}U^{+-}_{\bm{k}^{\prime}\bm{k}}\rangle^{*}$, and $\langle U^{-+}_{\bm{k}\bm{k}^{\prime}}U^{++}_{\bm{k}^{\prime}\bm{k}}\rangle=\langle U^{++}_{\bm{k}\bm{k}^{\prime}}U^{+-}_{\bm{k}^{\prime}\bm{k}}\rangle^{*}$. From the expression for $J(\langle n_{E}\rangle)$ [Eq. (16)] with $\langle n_{E}\rangle=\mathrm{diag}[n^{+}_{E\bm{k}},n^{-}_{E\bm{k}}]$, we get

[TABLE]

where $n^{+}_{E\bm{k}}=e\tau_{+}E_{y}\partial f_{0}(\varepsilon^{+}_{\bm{k}})/\partial k_{y}$ with $\tau_{+}\propto 1/n_{\mathrm{imp}}U_{0}^{2}$ the scattering rime for the band $\varepsilon^{+}_{\bm{k}}$. Here, we have used $\delta(\varepsilon^{+}_{\bm{k}}-\varepsilon^{-}_{\bm{k}^{\prime}})=\delta(\varepsilon^{-}_{\bm{k}}-\varepsilon^{+}_{\bm{k}^{\prime}})=0$ and the fact that $n^{-}_{E\bm{k}}=n^{-}_{E\bm{k}^{\prime}}=0$, because we have considered the case of $\mu>0$. Similarly, we obtain $[J(\langle n_{E}\rangle)]^{-+}_{\bm{k}}=\{[J(\langle n_{E}\rangle)]^{+-}_{\bm{k}}\}^{*}$. Then the off-diagonal density matrix resulting from $J(\langle n_{E}\rangle)$ is given by

[TABLE]

where $\mathrm{Re}$ and $\mathrm{Im}$ represent the real and imaginary part, respectively. The calculation of the conductivity from the definition $\sigma^{\rm II}_{xy}=\mathrm{Tr}[(-e)v_{x}\langle S^{\prime}_{E}\rangle]/E_{y}$ is somewhat long. See Appendix C for a detailed calculation. Finally, we find that $\sigma^{\rm II}_{xy}=0$ when the chemical potential $\mu$ lies sufficiently close to the Weyl nodes.

In the end, the total anomalous Hall conductivity of the Weyl metal reads

[TABLE]

which means that the anomalous Hall conductivity of Weyl metals is purely intrinsic, i.e., is determined by the Berry curvature of the filled band, as long as the chemical potential $\mu$ lies sufficiently close to the Weyl nodes. This plateau-like behavior of the anomalous Hall conductivity is consistent with a calculation by Burkov Burkov2014 . As mentioned in Sec. IV.1 the contribution from $J(\langle n\rangle)$, i.e., $\sigma^{\rm II}_{xy}$, corresponds to the vertex correction in the ladder approximation. Therefore it can be said that the vertex correction to the anomalous Hall conductivity of Weyl metals described by the Hamiltonian (73) is absent, as long as the chemical potential $\mu$ lies sufficiently close to the Weyl nodes.

VII.3 Chiral magnetic effect in Weyl metals

Let us see that the magnetic driving term (27) properly describes the chiral magnetic effect in Weyl metals. The chiral magnetic effect, the second term in Eq. (81), is written as

[TABLE]

where $Q_{\nu}$ and $\mu_{\nu}$ are the chirality and chiral chemical potential of a Weyl node, respectively. In two-node Weyl semimetals, we have $Q_{\pm}=\pm 1$ and $2\mu_{5}=\mu_{+}-\mu_{-}$ (see Fig. 2). As is obvious, the chiral magnetic effect does not occur when there is no chemical potential difference between the nodes as $\mu_{+}=\mu_{-}$.

VII.3.1 The case of $\mu_{+}=\mu_{-}$

First, it is informative to check the absence of the chiral magnetic effect in our formulation when $\mu_{+}=\mu_{-}$. Since there is no electric fields, we start from the diagonal density matrix $\langle\rho_{0}\rangle=\mathrm{diag}[f_{0}(\varepsilon^{+}_{\bm{k}}),f_{0}(\varepsilon^{-}_{\bm{k}})]$ with $f_{0}(\varepsilon^{\pm}_{\bm{k}})=1/[e^{(\varepsilon^{\pm}_{\bm{k}}-\mu)/T}+1]$. Without loss of generality we can consider the case of $\bm{B}=(0,0,B_{z})$ and can set $\mu>0$. Explicit matrix expression for $D_{B}(\langle\rho_{0}\rangle)$ can be found in Appendix D. Using Eqs. (38) and (166), the linear-response off-diagonal electron density matrix induced by the magnetic field is obtained as

[TABLE]

By substituting the Berry connections (76) into this equation, we get a concrete expression for $\langle S_{B}\rangle$. Also, the intrinsic diagonal density matrix induced by the magnetic field (41) is given by

[TABLE]

The velocity operator is written in the eigenstate basis as

[TABLE]

We are in a position to calculate the electric current. Since the chemical potential $\mu$ lies sufficiently close to the Weyl nodes, we do not need to take into account the contribution from the $t=+$ band. Finally, we find that the electric current along the magnetic field $\bm{B}$ vanishes as expected:

[TABLE]

where we have used the fact that the integrand is an odd function of $k_{z}$. This shows that the chiral magnetic effect does not occur when there is no chemical potential difference between the Weyl nodes, i.e., when $\mu_{+}=\mu_{-}$.

VII.3.2 The case of $\mu_{+}\neq\mu_{-}$

Next, let us consider the case where there is a chemical potential difference between the two nodes as $2\mu_{5}=\mu_{+}-\mu_{-}$. As is shown above, the expression for the chiral magnetic effect (90) is derived analytically from the four-band model (73) with nonzero $\mu_{5}\tau_{z}$ term. However, it is difficult to obtain analytically the eigenvalues and eigenstates in the presence of the $\mu_{5}\tau_{z}$ term, unlike in the case of the $b\sigma_{z}$ term where we can obtain the eigenstates analytically as Eq. (75). To determine whether the magnetic driving term (27) can describe the chiral magnetic effect, we consider two-band Hamiltonians around each Weyl node separately and sum each contribution to the electric current. The Hamiltonian around each Weyl node reads

[TABLE]

which gives the eigenvalues $\varepsilon^{\pm}_{\bm{q}}=\pm v_{F}\sqrt{q_{x}^{2}+q_{y}^{2}+q_{z}^{2}}$. Here, $\nu=\pm$ is the node index and $Q_{\pm}=\pm 1$ is the chirality of each Weyl node. The Fermi-Dirac distribution function around each node is given by $f_{0}^{\nu}(\varepsilon^{\pm}_{\bm{q}})=1/[e^{(\varepsilon^{\pm}_{\bm{q}}-\tilde{\mu}_{\nu})/T}+1]$ with $\tilde{\mu}_{\nu}=\mu+\mu_{\nu}$, where $\mu$ is a uniform chemical potential and $\mu_{\nu}$ is a chiral chemical potential such that $\mu_{+}\neq\mu_{-}$. Without loss of generality, we can set $\tilde{\mu}_{\nu}>0$. In linear response, the diagonal part $\langle\xi_{B,\nu}\rangle$ and the off-diagonal part $\langle S_{B,\nu}\rangle$ of the density matrix induced by the magnetic field are given respectively by Eqs. (92) and (91), with $f_{0}(\varepsilon^{\pm}_{\bm{q}})$ replaced by $f_{0}^{\nu}(\varepsilon^{\pm}_{\bm{q}})$. The velocity operator is given by

[TABLE]

In the present case, we may approximate the total electric current by the sum of the contributions from each node. Then the total electric current along the magnetic field $\bm{B}$ is calculated as $j_{z}=\sum_{\nu=\pm}(j_{z,\nu}^{S}+j_{z,\nu}^{\xi})$, where $j_{z,\nu}^{X}\equiv\mathrm{Tr}[(-e)v_{z,\nu}\langle X_{B,\nu}\rangle]$ ($X=S,\xi$). The contribution from $\langle S_{B,\nu}\rangle$ is calculated to be

[TABLE]

where $\delta\Lambda_{\nu}$ is a small momentum space around each node, $q=\sqrt{q_{x}^{2}+q_{y}^{2}+q_{z}^{2}}$, and we have used $\partial f_{0}^{\nu}(\varepsilon^{+}_{\bm{q}})/\partial\varepsilon^{+}_{\bm{q}}=-\delta(\varepsilon_{\bm{q}}-\tilde{\mu}_{\nu})$. The contribution from $\langle\xi_{B,\nu}\rangle$ is calculated to be

[TABLE]

where $\Lambda$ is the cutoff of the Fermi sea of Weyl cones.

Combining Eqs. (97) and (98), we arrive at the final expression for the chiral magnetic effect:

[TABLE]

where we have used $\sum_{\nu=\pm}Q_{\nu}\tilde{\mu}_{\nu}=\sum_{\nu=\pm}Q_{\nu}\mu_{\nu}$ and $\sum_{\nu=\pm}Q_{\nu}\Lambda=0$. This is in full agreement with the result by a field-theoretical approach (90). Here, we comment on an important observation in this study. The portion 2/3 comes from the Fermi surface contribution, as seen from Eq. (97). The portion 1/3 comes from the Fermi sea contribution, as see from Eq. (98). This observation is in contrast to a derivation by semiclassical wavepacket dynamics Son2013 ; Spivak2016 in which the chiral magnetic effect comes entirely from the Fermi surface contribution.

VIII Negative magnetoresistance in Weyl and Dirac Metals

So far we have shown that the electric driving term $D_{E}$ (26) and magnetic driving term $D_{B}$ (27) properly describe the chiral-anomaly induced phenomena in Weyl metals, the anomalous Hall effect and chiral magnetic effect, respectively. Another phenomenon manifested by the chiral anomaly is a negative magnetoresistance (or equivalently positive magnetoconductance) quadratic in magnetic field for parallel electric and magnetic fields in Weyl and Dirac metals Son2013 ; Burkov2014 ; Burkov2015a ; Spivak2016 . Such an unusual negative magnetoresistance has recently been experimentally observed in Dirac semimetals Na3Bi Xiong2015 , Cd3As2 Li2015 ; Li2016 , and ZrTe5 Li2016a , and in the Weyl semimetals TaAs Huang2015 and TaP Arnold2016 . Here, note that the usual magnetoresistance due to Lorentz force is positive. As denoted in the introduction, the chiral anomaly is described by the $\theta$ term $S_{\theta}=\int dt\,d^{3}r\,(e^{2}/4\pi^{2})\theta(\bm{r},t)\bm{E}\cdot\bm{B}$. As we can see from this expression for the $\theta$ term, the $\bm{E}\cdot\bm{B}$ contribution becomes largest in the case of parallel electric and magnetic fields. The positive quadratic magnetoconductivity arising from the chiral anomaly first derived by Son and Spivak reads Son2013

[TABLE]

where $\mu$ is the chiral potential and $\tau$ is the intervalley scattering time. Since their derivation is based on semiclassical wavepacket dynamics alone, it is difficult to take disorder scattering into account, and the microscopic origin of the intervalley scattering time in Eq. (100) is not precisely clear.

In this section, we apply the theory we have developed so far to the positive quadratic magnetoconductivity, starting from the microscopic continuum model of Weyl semimetals (73). It should be mentioned that the positive quadratic magnetoconductivity has been experimentally observed in the low magnetic field regime. Thus, our semiclassical treatment of magnetic fields can be justified. In our theory, the formal expression for the $\mu\nu$-component of the quadratic magnetoconductivity is written from Eq. (48) in the form

[TABLE]

Here, note that the angle between the electric and magnetic fields is arbitrary in this formalism.

Our evaluation of Eq. (101) involves a number of subtleties. It is informed by noting (i) that the expression of $\sigma_{zz}(B_{z}^{2})$ (100) has a linear dependence on $\tau$, and (ii) that this expression is derived via the chiral magnetic effect Son2013 ; Burkov2014 ; Burkov2015a ; Spivak2016 ; Li2016a . We use the fact that a magnetic driving term obtained from an off-diagonal density matrix is purely diagonal, and that a magnetic driving term obtained from a diagonal density matrix has diagonal and off-diagonal components that lead to the chiral magnetic effect as we have seen in Sec. VII.3 (see Appendix D for explicit expressions for the magnetic driving term). Indeed, our explicit calculations summarized schematically in Fig. 3 confirmed that among all the possible contributions that are second order in magnetic field, only magnetoconductances originating from the chiral anomaly and the usual Lorentz force survive.

To obtain the quadratic magnetoconductivity, we divide the calculation of the matrices $(\mathcal{L}^{-1}D_{B})^{2}\mathcal{L}^{-1}D_{E}(\langle\rho_{0}\rangle)$ and $\mathcal{L}^{-1}D_{B}\mathcal{L}^{-1}D_{E}(\langle\xi_{B}\rangle)$ in Eq. (101) into four steps as follows:

[TABLE]

where $\langle\rho_{0}\rangle=\mathrm{diag}[f_{0}(\varepsilon^{+}_{\bm{k}}),f_{0}(\varepsilon^{-}_{\bm{k}})]$ with $f_{0}(\varepsilon^{\pm}_{\bm{k}})=1/[e^{(\pm\varepsilon_{\bm{k}}-\mu)/T}+1]$ being the Fermi-Dirac distribution function, $\langle n\rangle$ and $\langle\xi\rangle$ indicate diagonal density matrices, $\langle S\rangle$ indicates an off-diagonal density matrix, and $\tau$ is the intervalley scattering time. See Eqs. (55), (57) and (58) for the formal expressions for $\langle n\rangle$, $\langle S\rangle$ and $\langle\xi\rangle$ respectively in a magnetic field. As we shall show below, the appearance of the intervalley scattering time $\tau$ is due to that the driving terms $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$ have nonzero values when integrated over a given valley (Weyl cone). In the following, we consider the low-temperature case where $T\ll\mu$, with $T$ and $\mu>0$ being the temperature and chemical potential of the system, respectively.

VIII.1 Quadratic magnetoconductivity for $\vec{E}\parallel\vec{B}$

VIII.1.1 Contribution from the chiral anomaly

Let us consider the case of parallel electric and magnetic fields $\bm{E}=(0,0,E_{z})$ and $\bm{B}=(0,0,B_{z})$. We start by obtaining the off-diagonal part of the density matrix induced solely by the electric field, $\langle S_{E}\rangle$ in Eq. (102). As we have seen in Sec. VII.2, this off-diagonal density matrix results in an intrinsic effect, i.e., the anomalous Hall effect. Using the expression for the electric driving term $D_{E}(\langle\rho_{0}\rangle)$ (33) and the Berry connection (76), we get

[TABLE]

At this stage we recall from Eq. (35) that there exists an extrinsic contribution to the off-diagonal density matrix $\langle S_{E}\rangle$ from the anomalous driving term $J(\langle n_{E}\rangle)$. Here, $\langle n_{E}\rangle=\mathrm{diag}[n^{+}_{E\bm{k}},n^{-}_{E\bm{k}}]$ with $n^{\pm}_{E\bm{k}}=e\tau_{\mathrm{intra}}E_{z}\partial f_{0}(\varepsilon^{\pm}_{\bm{k}})/\partial k_{z}$ ($\tau_{\mathrm{intra}}$ is the intravalley scattering time). Let us consider the same setup as in the evaluation of the anomalous Hall conductivity in Sec. VII.2. Namely, we consider the case of a short-range (on-site) disorder potential of the form $U(\bm{r})=U_{0}\sum_{i}\delta(\bm{r}-\bm{r}_{i})$, and assume that the correlation function satisfies $\langle U(\bm{r})U(\bm{r}^{\prime})\rangle=n_{\mathrm{imp}}U_{0}^{2}\,\delta(\bm{r}-\bm{r}^{\prime})$ with $n_{\mathrm{imp}}$ the impurity density. In this case, we can show explicitly from Eqs. (155) and (156) that

[TABLE]

and hence $[J(\langle n_{E}\rangle)]_{\bm{k}}^{+-}=[J(\langle n_{E}\rangle)]_{\bm{k}}^{-+}=0$, as long as the chemical potential $\mu$ lies sufficiently close to the Weyl nodes. Hence, the expression for $\langle S_{E}\rangle$, Eq. (103), remains valid even after the anomalous driving term is taken into account in the case of sufficiently small $\mu$.

Second, we compute the diagonal density matrix $\langle n_{EB}\rangle$ proportional to $E_{z}B_{z}$ in Eq. (102). This $\langle n_{EB}\rangle$ is the most important quantity in our derivation of the magnetoconductivity induced by the chiral anomaly, as is understood from its form similar to $\bm{E}\cdot\bm{B}$. As described in Appendix D, the magnetic driving term obtained from an off-diagonal matrix is purely diagonal. Using the expressions for the magnetic driving term $D_{B}(\langle S_{E}\rangle)$ (172) and (173) and the Berry connection (76), we have

[TABLE]

with

[TABLE]

where $\bm{1}$ is the $2\times 2$ identity matrix and $c_{\bm{k}}=[f_{0}(\varepsilon^{+}_{\bm{k}})-f_{0}(\varepsilon^{-}_{\bm{k}})](v_{F}k_{\perp}/4\varepsilon_{\bm{k}}^{3})\partial m/\partial k_{z}$. Also, using the expressions for the electric driving term (31) and (33) and the Berry phase contribution (41), we have

[TABLE]

where $\langle S_{EB}\rangle$ is a purely off-diagonal matrix obtained by replacing $\langle\rho_{0}\rangle$ by $\langle\xi_{B}\rangle$ in Eq. (103).

Now we show that these $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$ have a special property which is not present when the system is driven by an electric field alone. Around a Weyl node with momentum $\bm{q}=\bm{k}-W_{\pm}$ ($\bm{q}^{2}\ll 1$), $m_{-}(k_{z})$ is approximated as $m_{-}(q_{z})\approx\mp v_{F}^{2}(k_{0}/b)q_{z}$. Then we find that $\mathcal{F}_{\bm{q}}$ and $\partial f_{0}(\varepsilon^{\pm}_{\bm{q}})\Omega^{\pm}_{\bm{q},z}/\partial q_{z}$ are both even functions of $q_{x}$, $q_{y}$, and $q_{z}$. Accordingly, we find that the integral of $D_{B}(\langle S_{E}\rangle)+D_{E}(\langle\xi_{B}\rangle)$ over the Fermi surface of a given valley has a nonzero value: $\int_{\mathrm{FS}}\frac{d^{3}q}{(2\pi)^{3}}\,\sum_{m}[D_{B}(\langle S_{E}\rangle)+D_{E}(\langle\xi_{B}\rangle)]_{\bm{q}}^{mm}\neq 0$. As discussed in Sec. V, this is a consequence of the total particle number nonconservation in a given Weyl cone. Namely we obtain the rate of pumping

[TABLE]

In the case of isotropic Weyl fermions [i.e., $m_{-}(q_{z})\equiv v_{F}q_{z}$], Eq. (106) can be simplified to be

[TABLE]

where $\Omega_{\bm{q},a}^{+}=-q_{a}/(2q^{3})$ is the Berry curvature and $q=\sqrt{q_{x}^{2}+q_{y}^{2}+q_{z}^{2}}$. Integrating this over the Fermi surface we have

[TABLE]

which is indeed consistent with the general expression (72).

It follows from Eq. (55) that the intervalley scattering time $\tau$ appears when $\mathcal{L}^{-1}$ acts on $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$. When intravalley scattering is much stronger than intervalley scattering, $\mathcal{L}^{-1}$ acting on $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$ will also change the way in which the pumped charge is distributed across the Fermi surface, altering it so that it is equally distributed across the Fermi surface and corresponds simply to a change in chemical potential. In that limit the remaining steps in the calculation correspond precisely to the magnetoelectric effect calculation outlined in the previous section. In general though some anisotropy remains when $\mathcal{L}^{-1}$ acts on the pumped charge. To illustrate its potential role, we continue the calculation here considering the opposite limit in which only intervalley scattering is present. Then the diagonal density matrix $\langle n_{EB}\rangle$ is obtained as

[TABLE]

where

[TABLE]

is the component which represents the Fermi surface response in Eqs. (105) and (107). Note that we have neglected the Fermi sea response in Eqs. (105) and (107).

Third, we compute the off-diagonal density matrix $\langle S_{EB^{2}}\rangle$ proportional to $E_{z}B_{z}^{2}$ in Eq. (102). Before doing this, notice from Eq. (164) that there also are the diagonal components in $D_{B}(\langle n_{EB}\rangle)$ as $[D_{B}(\langle n_{EB}\rangle)]^{++}_{\bm{k}}=e^{3}E_{z}B_{z}^{2}\tau(\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{y}}\frac{\partial\tilde{\mathcal{F}}^{++}_{\bm{k}}}{\partial k_{x}}-\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{x}}\frac{\partial\tilde{\mathcal{F}}^{++}_{\bm{k}}}{\partial k_{y}})$ and $[D_{B}(\langle n_{EB}\rangle)]^{--}_{\bm{k}}=-e^{3}E_{z}B_{z}^{2}\tau(\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{y}}\frac{\partial\tilde{\mathcal{F}}^{--}_{\bm{k}}}{\partial k_{x}}-\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{x}}\frac{\partial\tilde{\mathcal{F}}^{--}_{\bm{k}}}{\partial k_{y}})$, which represent the Lorentz force contribution. Since these components are odd functions of $k_{x}$ and $k_{y}$ around a Weyl node, $D_{B}(\langle n_{EB}\rangle)$ does not drive the total particle number in a given Weyl cone. Hence the diagonal density matrix $\langle n^{\prime}_{EB^{2}}\rangle$ obtained from these $[D_{B}(\langle n_{EB}\rangle)]^{mm}_{\bm{k}}$ contains the intravalley scattering time as $\langle n^{\prime}_{EB^{2}}\rangle^{mm}_{\bm{k}}=\tau_{\mathrm{intra}}[D_{B}(\langle n_{EB}\rangle)]^{mm}_{\bm{k}}=me^{3}E_{z}B_{z}^{2}\tau\tau_{\mathrm{intra}}(\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{y}}\frac{\partial\tilde{\mathcal{F}}^{mm}_{\bm{k}}}{\partial k_{x}}-\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{x}}\frac{\partial\tilde{\mathcal{F}}^{mm}_{\bm{k}}}{\partial k_{y}})$. However, we see that $\langle n^{\prime}_{EB^{2}}\rangle$ does not contribute to the magnetoconductivity proportional to $B_{z}^{2}$, since it is an odd function of $k_{x}$ and $k_{y}$ and thus vanishes when integrated over momentum space. From Eqs. (57) and (166), we obtain the relevant off-diagonal density matrix as

[TABLE]

Here, let us recall that there exists an extrinsic contribution to the off-diagonal density matrix $\langle S_{EB^{2}}\rangle$ from the anomalous driving term $J(\langle n^{\prime}_{EB^{2}}\rangle)$, as in the case of $\langle S_{E}\rangle$. We can show explicitly by replacing $\langle n_{E}\rangle$ in Eqs. (155) and (156) by $\langle n^{\prime}_{EB^{2}}\rangle$ that

[TABLE]

and hence $[J(\langle n^{\prime}_{EB^{2}}\rangle)]_{\bm{k}}^{+-}=[J(\langle n^{\prime}_{EB^{2}}\rangle)]_{\bm{k}}^{-+}=0$, as long as the chemical potential $\mu$ lies sufficiently close to the Weyl nodes. Hence, the expression for $\langle S_{EB^{2}}\rangle$, Eq. (113), remains valid even after the anomalous driving term is taken into account in the case of sufficiently small $\mu$.

Notice from Eqs. (43) and (58) that $\langle S_{EB^{2}}\rangle$ is accompanied by an intrinsic Berry phase contribution to the diagonal part of the density matrix proportional to $E_{z}B_{z}^{2}$. From Eq. (43) we readily obtain

[TABLE]

We are now in a position to evaluate the $zz$-component of the magnetoconductivity proportional to $B_{z}^{2}$, $\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})=\mathrm{Tr}\{(-e)v_{z}[\langle S_{EB^{2}}\rangle+\langle\xi_{EB^{2}}\rangle]\}/E_{z}$. Since we have considered the scattering of electrons at the Fermi surface, we do not need to take into account the contribution from the $t=+$ band. From Eqs. (93), (113) and (115), an explicit expression for $\sigma_{zz}(B_{z}^{2})$ at low temperatures such that $T\ll\mu$ is obtained as

[TABLE]

which takes a positive value. We numerically find that $\sigma^{\mathrm{CA}}_{zz}(B_{z}^{2})$ is proportional to $1/\mu^{2}$ and $v_{F}^{3}$. This is consistent with the result by Son and Spivak, Eq. (100). We also find that, as in the case of the chiral magnetic effect (99), the contributions from $\langle S_{EB^{2}}\rangle$ and $\langle\xi_{EB^{2}}\rangle$ to $\sigma^{\mathrm{CA}}_{zz}(B_{z}^{2})$ are respectively 2/3 and 1/3 of the total value (116).

$\mathcal{C}(b,\Delta,T)$ in Eq. (116) is a dimensionless parameter of order $1$. In Fig. 4(a) we show the $k_{0}$ dependence of $\mathcal{C}(b,\Delta,T)$ with $\Delta=0$, where $2k_{0}=2\sqrt{b^{2}-\Delta^{2}}/v_{F}$ is the distance between the two Weyl nodes in momentum space. We find that $\mathcal{C}$ is an increasing function of $k_{0}$ and that is saturates when $k_{0}$ is sufficiently large. As shown in Figs. 4(b), (c), and (d), a single Weyl cone splits into two Weyl cones when $k_{0}\neq 0$. The two Weyl cones partly overlap at the Fermi level and therefore there is only one Fermi surface in the region where $\mathcal{C}$ is increasing as $k_{0}$ becomes larger. The separate contributions from the two Weyl cones are therefore incomplete in this region. When two Fermi surfaces are present, the value of $\mathcal{C}$ saturates to the contributions from the two independent Weyl cones. In the present case, the value of $k_{0}$ at which the saturation begins is given by $k_{0}=\mu/v_{F}$. Here, note that we have not considered the $k_{0}$ dependence of the intervalley scattering time $\tau$. Although it is not easy to take into account such a $k_{0}$ dependence in $\tau$, it is expected that the value of $\tau$ will become larger as the value of $k_{0}$ becomes larger, i.e., as the distance between the two Weyl nodes becomes longer. This is because the strength of intervalley scattering processes will decrease as the value of $k_{0}$ becomes larger, since intervalley scattering processes with large $k_{0}$ require large momentum transfer. Thus, based on these viewpoints, we conclude that the value of $\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})$ in Weyl metals becomes larger as the distance between the two Weyl nodes becomes longer.

Finally, we note that the contribution to $\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})$ from a single isotropic Weyl cone in the low temperature limit $T/\mu\to 0$ is obtained by setting $\mathcal{C}(b,\Delta,T)=1/2$ and $b=\Delta=0$ in Eq. (116) (see Fig. 4), since our Hamiltonian (74) with $b=\Delta=0$ (i.e, $k_{0}=0$) possesses a single isotropic Weyl cone. Also, the saturated value $\mathcal{C}=1$ obtained at large $k_{0}$ with $T/\mu\to 0$ is the contribution from two isotropic Weyl cones. The value from a single isotropic Weyl cone, $\mathcal{C}=1/2$, is smaller than the value obtained by invoking semiclassical wavepacket dynamics and neglecting scattering [Eq. (100)] that corresponds to $\mathcal{C}=1$. This result can be interpreted as follows. The rate of pumping into a valley (110), $\partial N/\partial t=e^{2}E_{z}B_{z}/4\pi^{2}$, is the same as that obtained by invoking semiclassical wavepacket dynamics and neglecting scattering [Eq. (2)]. When intravalley scattering is dominant (i.e., in the limit of extremely weak intervalley scattering), by following the chiral magnetic effect we will find that $\mathcal{C}=1$. In general, however, the electron distribution on the Fermi surface $\tilde{\mathcal{F}}^{mm}_{\bm{k}}$ (130) can be anisotropic due to intravalley scattering. The numerical results in Fig. 4, obtained neglecting intravalley scattering, demonstrate that the precise value of $\mathcal{C}$ can be depend on the details of disorder.

VIII.1.2 The case of Dirac semimetals

Here, we briefly discuss the case of Dirac semimetals. To this end let us consider the contribution from a single Dirac cone, i.e., two degenerate Weyl cones. Time reversal symmetry is preserved in Dirac metals, and therefore we have to set $b=0$. Also, for the existence of a Dirac point we have to set $\Delta=0$. As mentioned above, Eq. (116) with $b=\Delta=0$ represents the contribution from a single Weyl cone. Hence, the value of $\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})$ coming from a single Dirac cone in the low temparature limit $T/\mu\to 0$ is obtained by setting $b=\Delta=0$ (i.e., $\mathcal{C}=1/2$) and multiplying a factor $2$ in Eq. (116):

[TABLE]

where $\tau_{\mathrm{Dirac}}$ is the intervalley scattering time between the two Weyl cones that are degenerate in momentum space. As discussed in Sec. V the effective Hamiltonian around a Dirac point in Dirac semimetals which are protected by crystalline symmetry can be written in the block-diagonal form Yang2014

[TABLE]

where $A$ and $B$ denote two states after a unitary transformation, and $\mathcal{H}_{\pm}(\bm{q})$ is a $2\times 2$ Weyl Hamiltonian with chirality $\pm 1$. In general, atomically smooth disorder does not allow scattering processes between two different states $A$ and $B$. Therefore, the intervalley scattering time $\tau_{\mathrm{Dirac}}$ is much larger than intravalley scattering time.

VIII.1.3 Contribution from the Lorentz force

Next, we compute the magnetoconductivity proportional to $B_{z}^{2}$ due to the Lorentz force. Such a magnetoconductivity does not originate from the intrinsic (i.e., the Berry curvature) contribution but originates solely from the Fermi surface contribution. Thus, only the diagonal components in the electric and magnetic driving terms are relevant to obtaining it. We start by obtaining the diagonal part of the density matrix induced by the electric field. It is easily seen that the integral of $D_{E}(\langle\rho_{0}\rangle)$ over the Fermi surface of a given Weyl cone is zero. Then, from the discussion in Sec. V, the diagonal density matrix obtained from $D_{E}(\langle\rho_{0}\rangle)$ contains the intravalley scattering time as $\langle n_{E}\rangle=\tau_{\mathrm{intra}}D_{E}(\langle\rho_{0}\rangle)$, where $n_{E\bm{k}}^{+}=eE_{z}\tau_{\mathrm{intra}}[\partial f_{0}(\varepsilon^{+}_{\bm{k}})/\partial k_{z}]$ and $n_{E\bm{k}}^{-}=0$ (since we consider the case of $\mu>0$). From Eq. (164), the diagonal components in $D_{B}(\langle n_{E}\rangle)$ are given by

[TABLE]

and $[D_{B}(\langle n_{E}\rangle)]^{--}_{\bm{k}}=0$. We see that Eq. (120) indeed represents the Lorentz force term, as $[D_{B}(\langle n_{E}\rangle)]^{++}_{\bm{k}}=e(\bm{v}_{\bm{k}}\times\bm{B})\cdot\nabla_{\bm{k}}n^{+}_{E\bm{k}}$ with $\bm{v}_{\bm{k}}=\nabla_{\bm{k}}\varepsilon_{\bm{k}}$. Here, note that Eq. (120) can be rewritten as

[TABLE]

from which we see clearly that $[D_{B}(\langle n_{E}\rangle)]^{++}_{\bm{k}}=0$ for the present case of the isotropic energy band $\varepsilon_{\bm{k}}$ in the $(k_{x},k_{y})$ plane, $\varepsilon_{\bm{k}}=\sqrt{v_{F}^{2}k_{x}^{2}+v_{F}^{2}k_{y}^{2}+m_{\pm}(k_{z})^{2}}$. Namely, there are no contributions from the Lorentz force in the present case, as is also understood from the fact that the Lorentz force does not act when the velocity is along the magnetic-field direction. However, when the energy band $\varepsilon_{\bm{k}}$ is anisotropic in the $(k_{x},k_{y})$ plane, we see that in general $[D_{B}(\langle n_{E}\rangle)]^{++}_{\bm{k}}\neq 0$. See, for example, Ref. Pal2010 for a detailed discussion on the conditions for nonzero contributions from the Lorentz force in the case of $\bm{E}\parallel\bm{B}$.

The Weyl cones in real materials will be anisotropic. In the following, let us assume that there is an anisotropy in $\varepsilon_{\bm{k}}$, for example, by replacing $m_{\pm}(k_{z})$ by $m_{\pm}(k_{x},k_{z})=b\pm\sqrt{v_{F}^{2}k_{z}^{2}+\alpha k_{x}^{2}+\Delta^{2}}$. In such as case we have $[D_{B}(\langle n_{E}\rangle)]^{++}_{\bm{k}}\neq 0$. Again we find that the integral of $D_{B}(\langle n_{E}\rangle)$ over the Fermi surface of a given Weyl cone is zero. Then we see from Eq. (55) the diagonal density matrix obtained from $D_{B}(\langle n_{E}\rangle)$ contains the intravalley scattering time as $\langle n^{\prime}_{EB}\rangle=\tau_{\mathrm{intra}}D_{B}(\langle n_{E}\rangle)$. The same procedure can be applied to obtain the diagonal density matrix proportional to $E_{z}B_{z}^{2}$, and we obtain $\langle n_{EB^{2}}\rangle=\tau_{\mathrm{intra}}D_{B}(\langle n_{EB}\rangle)=\mathrm{diag}[(eB_{z})^{2}\tau_{\mathrm{intra}}^{2}(\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{y}}\frac{\partial}{\partial k_{x}}-\frac{\partial\varepsilon_{\bm{k}}}{\partial k_{x}}\frac{\partial}{\partial k_{y}})^{2}n^{+}_{E\bm{k}},0]$. Finally, noting that only the $t=-$ band is relevant, the magnetoconductivity proportional to $B_{z}^{2}$ induced by the Lorentz force at low temperatures such that $T\ll\mu$ is calculated to be

[TABLE]

which takes a negative value. Here, $\sigma_{zz}^{0}\approx e^{2}(\mu^{2}/v_{F})\tau_{\mathrm{intra}}$ is the Drude conductivity, and $\omega_{c}=eB_{z}v_{F}^{2}/\mu$ is the cyclotron frequency.

VIII.1.4 Total magnetoconductivity of the system

We have seen that the contribution from the Lorentz force is absent when the Weyl cones are isotropic. In this case the total quadratic magnetoconductivity for parallel electric and magnetic fields is written in the low-field limit as

[TABLE]

On the other hand, when the Weyl cones are anisotropic, the contribution from the Lorentz force is present. In this case, from Eqs. (116) and (121), the total quadratic magnetoconductivity for parallel electric and magnetic fields is written in the low-field limit as $\sigma_{zz}(B_{z}^{2})=\sigma_{zz}^{\mathrm{CA}}(B_{z}^{2})+\sigma_{zz}^{\mathrm{LF}}(B_{z}^{2})$. The ratio of these two contributions to the magnetoconductivity is estimated as

[TABLE]

In general, the intervalley scattering time $\tau$ is much larger than the intravalley scattering time $\tau_{\mathrm{intra}}$, i.e., $\tau_{\mathrm{intra}}/\tau\ll 1$, since the intervalley scatterings require large momentum transfers, i.e., the number of intervalley scattering processes is much smaller than that of intravalley scattering processes. On the other hand, in usual high-mobility semiconductors the condition $(\mu\tau_{\mathrm{intra}})^{2}\gg 1$ is satisfied. When the contribution from the chiral anomaly dominates that from the Lorentz force, we have the total magnetoconductivity in the low-field limit as

[TABLE]

We have confirmed that the magnetoconductivity linear in $B$ vanishes in our formalism. Namely, the magnetoconductivity that has the lowest $B$-dependence is the one proportional to $B^{2}$, which is in agreement with experimental results on the magnetoconductivity in the low-field limit in the Dirac semimetal Cd3As2 Li2016 . On other hand, note that in this study we have not considered the contribution from the weak antilocalization such that $\sim a\sqrt{B}$ with $a<0$ in the low-field limit which is observed in the Weyl semimetal TaAs Huang2015 .

VIII.2 Quadratic magnetoconductivity for $\vec{E}\perp\vec{B}$

Let us consider the case of perpendicular electric and magnetic fields as $\bm{E}=(0,0,E_{z})$ and $\bm{B}=(B_{x},0,0)$. Even in this case, the intrinsic contribution to the magnetoconductivity proportional to $B_{x}^{2}$, which resembles that from the chiral anomaly (116), takes nonzero value in our formalism. However, it is not characterized by the intervalley scattering time $\tau$ which appeared in the case of parallel fields, but characterized by the intravalley scattering time $\tau_{\mathrm{intra}}$. We show that as the result the contribution from the Lorentz force dominates the intrinsic contribution, which is consistent with experimental results.

We compute the intrinsic contribution to the quadratic magnetoconductivity. The exactly same calculation as in the case of parallel electric and magnetic fields can be applied here. After a calculation, we obtain the magnetic driving term which is proportional to $E_{z}B_{x}$ as

[TABLE]

with

[TABLE]

We also have the electric driving term which is proportional to $E_{z}B_{x}$ as

[TABLE]

where $\langle S_{EB}\rangle$ is a purely off-diagonal matrix. Here, we show that these $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$ have the opposite property of the case of parallel electric and magnetic fields. Around a Weyl node with momentum $\bm{q}=\bm{k}-W_{\pm}$ ($\bm{q}^{2}\ll 1$), we find that $\mathcal{G}_{\bm{q}}$ and $\partial f_{0}(\varepsilon^{\pm}_{\bm{q}})\Omega^{\pm}_{\bm{q},x}/\partial q_{z}$ are odd functions of $q_{x}$ and $q_{z}$. Accordingly, we find that the integral of $D_{B}(\langle S_{E}\rangle)+D_{E}(\langle\xi_{B}\rangle)$ over the Fermi surface of a given valley is zero: $\int_{\mathrm{FS}}\frac{d^{3}q}{(2\pi)^{3}}\,\sum_{m}[D_{B}(\langle S_{E}\rangle)+D_{E}(\langle\xi_{B}\rangle)]^{mm}_{\bm{q}}=0$. As discussed in Sec. V, this is a consequence of the total particle number conservation in a given Weyl cone. Namely we obtain

[TABLE]

Then it follows from Eq. (55) that the intravalley scattering time $\tau_{\mathrm{intra}}$ appears when $\mathcal{L}^{-1}$ acts on $D_{B}(\langle S_{E}\rangle)$ and $D_{E}(\langle\xi_{B}\rangle)$. In the end, the diagonal density matrix $\langle n_{EB}\rangle$ is obtained as

[TABLE]

where

[TABLE]

is the component which represents the Fermi surface response in Eqs. (125) and (127). The off-diagonal part $\langle S_{EB^{2}}\rangle$ and diagonal part $\langle\xi_{EB^{2}}\rangle$ of the density matrix proportional to $E_{z}B_{x}^{2}$ can be obtained in a similar way as Eqs. (113) and (115), respectively.

Finally, noting that only the $t=-$ band is relevant, the $zz$-component of the quadratic magnetoconductivity from the intrinsic contribution at low temperatures such that $T\ll\mu$ is obtained as

[TABLE]

which takes a positive value with $\mathcal{C}^{\prime}(b,\Delta,T)\sim\mathcal{O}(0.1)$. As mentioned above, this form resembles that from the chiral anomaly (116). However, the appearance of the intravalley scattering time indicates that the total particle number in a given valley is conserved, while the chiral anomaly requires nonconservation of the total particle number in a given valley. Thus, the magnetoconductivity (131) does not originate from the chiral anomaly.

The quadratic magnetoconductivity due to the Lorentz force in the case of perpendicular electric and magnetic fields can be evaluated in the same way as Eq. (121). At low temperatures such that $T\ll\mu$ we have

[TABLE]

which takes a negative value. Here, $\sigma_{zz}^{0}\approx e^{2}(\mu^{2}/v_{F})\tau_{\mathrm{intra}}$ is the Drude conductivity, and $\omega_{c}=eB_{x}v_{F}^{2}/\mu$ is the cyclotron frequency. Note that, in contrast to the case of parallel electric and magnetic fields, $\sigma_{zz}^{\mathrm{LF}}(B_{x}^{2})$ is nonzero even when the Weyl cones are isotropic.

Therefore from Eqs. (131) and (132), the total quadratic magnetoconductivity for perpendicular electric and magnetic fields is written in the low-field limit as $\sigma_{zz}(B_{x}^{2})=\sigma_{zz}^{\mathrm{In}}(B_{x}^{2})+\sigma_{zz}^{\mathrm{LF}}(B_{x}^{2})$. The ratio of these two contributions to the magnetoconductivity is estimated as

[TABLE]

where we have used that in usual high-mobility semiconductors the condition $(\mu\tau_{\mathrm{intra}})^{2}\gg 1$ is satisfied. Thus, we have the total magnetoconductivity in the low-field limit as

[TABLE]

which is in agreement with experimental results in Dirac and Weyl semimetals Xiong2015 ; Li2015 ; Li2016 ; Li2016a ; Huang2015 ; Arnold2016 . As in the case of parallel electric and magnetic fields, we have confirmed that the magnetoconductivity linear in $B$ vanishes in our formalism.

IX Discussions

We have derived a multi-band quantum kinetic equation for the Bloch-state density matrix of a crystal that includes driving terms related to both electric and magnetic fields. The magnetic driving term $D_{B}$ is new, as far as we are aware, and has a simple form that elegantly generalizes the well known Lorentz force driving term of the scalar single-band Boltzmann quantum kinetic equation, $-e(\nabla_{\bm{k}}\varepsilon_{\bm{k}}\times\bm{B})\cdot\nabla_{\bm{k}}f_{\bm{k}}$, to a multi-band matrix form. We have also shown that our quantum kinetic equation captures two effects that have been previously identified by examining Bloch-state wavepacket dynamics in multi-band systems, namely an anomalous contribution to velocity that is proportional to $\dot{\bm{k}}\times\bm{\Omega}$ and a relative change in the momentum-space density of states associated with a particular band that is proportional to $\bm{B}\cdot\bm{\Omega}$, where $\bm{B}$ is an external magnetic field and $\bm{\Omega}$ is the momentum-space Berry phase curvature associated with the band. In transport the former effect is responsible for the intrinsic anomalous Hall effect of magnetic crystals. When disorder scattering is neglected Son2013 , the two effects combine to yield a remarkable condensed matter realization of the chiral anomaly in which charge is pumped between Fermi surfaces surrounding different Weyl points. As an illustration of the physics that can be captured using our quantum kinetic equation we have examined how the chiral anomaly is observably manifested in the magnetoconductance of Weyl semimetals, discovering a complex interplay between electric and magnetic field driving terms, free-particle dynamics, and scattering contributions to the quantum kinetic equations. We have found that the charge pumping effect survives only when intervalley scattering is very weak compared to intravalley scattering.

The free-evolution ($P$), scattering ($K$), and electric and magnetic driving term ($D_{E}$ and $D_{B}$) contributions to the quantum kinetic equation can all be viewed as operators that act on the Bloch-state density matrices that characterize the transport steady state. Only band-off-diagonal elements of this density matrix have time dependence in the absence of external fields and disorder, i.e., the projection of $P$ onto the band-diagonal part of the Liouville density-matrix vector space is zero. For this reason, the roles of diagonal and off-diagonal density-matrix components in the quantum kinetic equations are quite distinct. The magnetic driving term maps the equilibrium density matrix $\langle\rho_{0}\rangle$ (i.e., the Fermi-Dirac distribution function), which is diagonal in a band-eigenstate representation, to a purely band off-diagonal density matrix. In other words, $D_{B}\left(\langle\rho_{0}\rangle\right)$ is a purely band off-diagonal matrix. On the other hand, the matrix $P^{-1}D_{B}\left(\langle\rho_{0}\rangle\right)$ has diagonal components $\langle\xi_{B}\rangle$. Similar diagonal components are obtained whenever $P^{-1}D_{B}$ acts on any band-diagonal density matrix components. These magnetic-field induced changes in the occupation probabilities of Bloch states capture the Berry-phase related density-of-states changes in magnetic field of semiclassical wavepacket dynamics. As is shown in Sec. VII.3, $\langle S_{B}\rangle$ and $\langle\xi_{B}\rangle$ combine to produce a magnetic-field driven current when the local chemical potentials of Fermi surface pockets near different Weyl points are different, thus capturing what is referred to as the chiral magnetic effect. As shown in Sec. VIII.1, the density matrices $\langle S_{EB^{2}}\rangle$ and $\langle\xi_{EB^{2}}\rangle$, which result in the positive quadratic magnetoconductivity in Weyl and Dirac metals, are closely related $\langle S_{B}\rangle$ and $\langle\xi_{B}\rangle$. In the limit of extremely weak intervalley scattering, the chemical potential difference of the chiral magnetic effect is established by balancing charge pumping between valleys and intervalley relaxation. Although our calculations show that the precise numerical value of the positive quadratic magnetoconductivity can be quite complex and depend subtly on the details of specific systems, it is still true, as explained in earlier theoretical work Son2013 ; Burkov2014 ; Burkov2015a ; Spivak2016 ; Li2016a , that the positive quadratic magnetoconductivity is basically a simple consequence of the chiral magnetic effect. It is important to emphasize that the anomalous Hall effect and the chiral magnetic effect arise essentially from band off-diagonal parts of the electric and magnetic driving terms, respectively. This indicates that transport phenomena induced by the chiral anomaly are an important example of interband coherence response in conductors.

This study is an extension of Ref. Culcer2016 which developed a quantum kinetic theory for crystals that is sensitive to recent developments in electronic structure theory that allow accurate Wannier representations of the Bloch-state Hamiltonian which are convenient for transport theory to be constructed, and also sensitive to recent awareness of the importance of momentum-space Berry phase effects in a number of different contexts. The aim of the theory is to be able to account for intrinsic effect related to Bloch-state wave function properties and extrinsic effects related to disorder scattering to be accounted for on a consistent footing in real materials with complicated electronic structure. The quantum kinetic theory in Ref. Culcer2016 , which is the $\bm{B}=0$ limit of the theory presented in this paper, conveniently captures effects that appear in Kubo formula formulations of transport theory with disorder treated as a perturbation, as ladder-diagram vertex corrections, for example the well known absence of a spin Hall conductivity Inoue2004 and anomalous Hall conductivity Inoue2006 in certain Rashba models. An important aspect of the detailed calculations we have presented for the simplified Weyl semimetal toy model, is its demonstration that the corresponding vertex corrections are absent in the anomalous Hall conductivity (89) and positive quadratic magnetoconductivity (116) of Weyl metals, as long as the chemical potential lies sufficiently close to the Weyl points.

We briefly discuss other possible applications of our theory. In this connection, we emphasize that although our theory is formulated as an independent particle theory, we really intend it as a theory of independent quasiparticles in a mean-field approximation to an interacting electron theory. For example the mean-field quasiparticles can be viewed as the Kohn-Sham quasiparticles of density functional theory. Because our theory calculates the density-matrix response to electric and magnetic fields, it can be used to determine the response of any crystal observables, including observables like the spin density that can feed back into the steady state crystal Hamiltonian. All single-particle observables $\mathcal{O}$ maintain their crystal periodicity when they respond to spatially constant electric and magnetic fields and therefore have expectation values of the form

[TABLE]

where $\langle\rho\rangle$ is a density matrix we have calculated. The evaluation of field-induced spin currents and spin densities, which are related to the current-induced torques of spintronics is one practically important problem to which our quantum kinetic theory can be applied. Our theory can be flexibly applied to metals, insulators, nodal-line semimetals, and to systems with many other types of electronic structure. It can be applied to toy models that capture the essence of different phenomena or to realistic models of specific materials. We anticipate, for example, that it will have interesting implications for the properties of 2D multivalley systems such as graphene and transition metal dichalcogenides. The kinetic equation approach can also be applicable to interacting systems that are described at a mean-field theory level Culcer2011 ; Culcer2013 ; Culcer2016a . It will be interesting to investigate from the quantum kinetic theory viewpoint the interplay between nontrivial band topology and interactions in the presence of electric and magnetic fields.

X Summary

In summary, we have developed a general quantum kinetic theory of low-field magnetotransport in weakly disordered multiband systems. Our theory naturally incorporates momentum-space Berry phase effects, which are often discussed in the context of semiclassical wavepacket dynamics, into transport theory. By applying the Wigner transformation to the quantum Liouville equation, we have derived a quantum kinetic equation (25), which is the principle result of this study. From this equation we have obtained a generic expression for the magnetoconductivity that is applicable for arbitrary angle between electric and magnetic fields. We note that the purely band-diagonal contributions in our theory describe all the regular Fermi-surface dominated effects that are familiar from textbooks. Our theory is able to simply explain and predict when large momentum-space Berry curvatures yield important corrections to the standard Boltzmann theory of Fermi-surface dominated magnetotransport.

We have applied our theory to study transport phenomena induced by the chiral anomaly in a toy model of two-node Weyl metals (73). We have shown that the anomalous Hall effect in the model is purely intrinsic, which means that the vertex correction in the ladder-diagram approximation is absent. We have also shown, by constructing the linear response theory to a magnetic field in the absence of electric fields, that the magnetic driving term we have newly introduced describes properly the chiral magnetic effect in a continuum model of Weyl metals. We have obtained an explicit expression for the positive quadratic magnetoconductivity that includes the intervalley scattering time in parallel electric and magnetic fields. In the process of obtaining this expression, we have shown that the vertex correction in the ladder-diagram approximation is absent. Our study indicates that the positive quadratic magnetoconductivity is a consequence of the chiral magnetic effect. On the other hand, in the case of perpendicular electric and magnetic fields, we have obtained a negative quadratic magnetoconductivity due to the Lorentz force that includes the intravalley scattering time, which is in agreement with experimental results. We have clarified that the chiral anomaly is observable only when intervalley scattering at the Fermi energy is very weak compared to intravalley scattering.

Acknowledgements.

Work at Austin was supported by the Department of Energy, Office of Basic Energy Sciences under Contract No. DE-FG02-ER45958 and by the Welch foundation under Grant No. TBF1473. A.S. is supported by the JSPS Overseas Research Fellowship. D.C. is supported by the Australian Research Council Centre of Excellence in Future Low-Energy Electronics Technologies (Project No. CE170100039). A.H.M. acknowledges stimulating conversations with Anton Burkov and Boris Spivak.

Appendix A General properties of the magnetic driving term $D_{B}$

In this Appendix, we consider general properties of the magnetic driving term (27). We write the Hamiltonian and the density matrix of a system as $\mathcal{H}_{0}=\sum_{m}\varepsilon_{m}|m\rangle\langle m|$ and $\langle\rho_{0}\rangle=\sum_{m}f_{0m}|m\rangle\langle m|$ with $m$ a band index and $f_{0m}$ a Fermi-Dirac distribution function, where we have omitted the wavevector dependences in these equations to simplify the notation. For concreteness, we consider the case of $\bm{B}=(0,0,B_{z})$. Then the driving term (27) is written as

[TABLE]

First, we consider a diagonal component of $D_{B}(\langle\rho_{0}\rangle)$, i.e., $\langle m|D_{B}(\langle\rho_{0}\rangle)|m\rangle$. We immediately get

[TABLE]

where $\partial_{a}=\partial/\partial k_{a}$. Then we have

[TABLE]

from which a diagonal component is obtained as

[TABLE]

where we have used $\langle m^{\prime}|\partial_{a}n^{\prime}\rangle+\langle\partial_{a}m^{\prime}|n^{\prime}\rangle=\partial_{a}(\delta_{m^{\prime}n^{\prime}})=0$. Similarly we have

[TABLE]

Then we see that

[TABLE]

Next, we consider a diagonal component of $P^{-1}D_{B}(\langle\rho_{0}\rangle)$, i.e., $\langle m|P^{-1}D_{B}(\langle\rho_{0}\rangle)|m\rangle$. Here, the matrix $P$ is defined by $P\langle\rho\rangle\equiv\frac{i}{\hbar}[\mathcal{H}_{0},\langle\rho\rangle]$, i.e., $P$ is the precession term that accounts for the time dependence of the density matrix in the absence of fields and disorder. Recall from Eq. (35) that when the matrix $P^{-1}$ acts on a driving term $D$, it simply multiplies $D$ by a numerical factor:

[TABLE]

Let us take a closer look at the second term in the third line of Eq. (138). Noticing the fact that $\partial_{a}(\sum_{m^{\prime}}|m^{\prime}\rangle\langle m^{\prime}|)=0$, we can replace $\varepsilon_{m^{\prime}}$ by $(\varepsilon_{m^{\prime}}-\varepsilon_{n^{\prime}})$ in Eq. (137). Then we can rewrite it as

[TABLE]

from which we have

[TABLE]

Then a diagonal element of this matrix is obtained as $\langle m|P^{-1}\frac{D\mathcal{H}_{0}}{Dk_{y}}\frac{D\langle\rho_{0}\rangle}{Dk_{x}}|m\rangle=-i\hbar f_{0m}\langle\partial_{y}m|\partial_{x}m\rangle$. Similarly we have $\langle m|P^{-1}\frac{D\langle\rho_{0}\rangle}{Dk_{x}}\frac{D\mathcal{H}_{0}}{Dk_{y}}|m\rangle=+i\hbar f_{0m}\langle\partial_{x}m|\partial_{y}m\rangle$. In the end, we see that

[TABLE]

where $\Omega^{m}_{z}=i\langle\partial_{x}m|\partial_{y}m\rangle-i\langle\partial_{y}m|\partial_{x}m\rangle$ is the $z$ component of the Berry curvature of band $m$. Note that only the terms of the form of (143) give rise to nonzero contributions to the diagonal component of $P^{-1}D_{B}(\langle\rho_{0}\rangle)$.

Appendix B General expression for $\mathrm{Tr}\,[D_{B}(\langle S_{E}\rangle)]$ in parallel electric and magnetic fields

In this Appendix, we derive an explicit expression for the diagonal matrix element $[D_{B}(\langle S_{E}\rangle)]^{mm}$ in a general system with $\mathcal{H}_{0}=\sum_{m}\varepsilon_{m}|m\rangle\langle m|$ and $\langle\rho_{0}\rangle=\sum_{m}f_{0m}|m\rangle\langle m|$ with $m$ a band index and $f_{0m}$ a Fermi-Dirac distribution function, where we have omitted the wavevector $\bm{k}$ dependences in these equations to simplify the notation. For concreteness, we consider the case of $\bm{E}=(0,0,E_{z})$ and $\bm{B}=(0,0,B_{z})$. An off-diagonal component of the electric driving term (26) reads

[TABLE]

with $\partial_{a}=\partial/\partial k_{a}$, from which we obtain the off-diagonal part of the density matrix induced by the electric field

[TABLE]

where $n\neq n^{\prime}$. On the other hand, the magnetic driving term (27) is written as

[TABLE]

Since we are focusing on the Fermi surface response, we consider only the terms proportional to $\partial_{x}f_{0}$ in $D\langle S_{E}\rangle/Dk_{x}$:

[TABLE]

We also have the relevant term

[TABLE]

where we have used $\partial_{a}(\sum_{m^{\prime}}|m^{\prime}\rangle\langle m^{\prime}|)=0$. Note that the terms proportional to $\partial_{y}\varepsilon$ in $D\mathcal{H}_{0}/Dk_{y}$ do not contribute to the diagonal part of $D_{B}(\langle S_{E}\rangle)$. Then we obtain

[TABLE]

Similarly we have

[TABLE]

Finally we arrive at a general expression for the rate of pumping from the Fermi surface contribution:

[TABLE]

where $\Omega_{a}^{m}=\epsilon^{abc}\,i\langle\partial_{b}m|\partial_{c}m\rangle$ is the Berry curvature of band $m$.

Appendix C Evaluation of the anomalous Hall conductivity $\sigma_{xy}^{\mathrm{II}}$

In this Appendix, we calculate explicitly the anomalous Hall conductivity $\sigma_{xy}^{\mathrm{II}}$ which results from the extrinsic (i.e., Fermi surface) contribution, and show that it vanishes when the chemical potential $\mu$ lies sufficiently close to the Weyl points. We have considered the case where an electric field is applied along the $y$ direction as $\bm{E}=E_{y}\bm{e}_{y}$. In this case, the diagonal part of the density matrix induced by the electric field $\langle n_{E}\rangle=\mathrm{diag}[n^{+}_{E\bm{k}},n^{-}_{E\bm{k}}]$ is given by

[TABLE]

where we have used $\varepsilon^{\pm}_{\bm{k}}=\pm\varepsilon_{\bm{k}}=\pm\sqrt{v_{F}^{2}(k_{x}^{2}+k_{y}^{2})+m^{2}}$ and considered the case of $\mu>0$. Then from Eq. (87) we obtain the imaginary part of $[J(\langle n_{E}\rangle)]^{+-}_{\bm{k}}$ as

[TABLE]

where we have used that the contribution from $n_{E\bm{k}}^{+}$ vanishes since $\sin\theta^{\prime}=k_{y}^{\prime}/k_{\perp}^{\prime}$ and $\cos\theta^{\prime}=k_{x}^{\prime}/k_{\perp}^{\prime}$ are odd functions. Similarly, we obtain the real part of $[J(\langle n_{E}\rangle)]^{+-}_{\bm{k}}$ as

[TABLE]

Finally, the extrinsic contribution to the anomalous Hall conductivity $\sigma^{\mathrm{II}}_{xy}$ is calculated from the definition $\sigma^{\mathrm{II}}_{xy}=\mathrm{Tr}[(-e)v_{x}\langle S^{\prime}_{E}\rangle]/E_{y}$. Using Eqs. (84) and (88), in the case of small chemical potential $\mu$ we obtain

[TABLE]

where we have omitted the subscript $-$ in $\varepsilon_{-\bm{k}}$ in the third through last lines. Note that the contributions only from $t=-$ survive since the chemical potential $\mu$ lies at $\varepsilon_{-\bm{k}}^{+}$. Also, we have used the fact that $m_{-}(k_{z})\approx\mp v_{F}^{2}(k_{0}/b)q_{z}$ around the Weyl nodes, where $\bm{q}=\bm{k}-W_{\pm}=(k_{x},k_{y},k_{z}\mp k_{0})$ with $k_{0}=\sqrt{b^{2}-\Delta^{2}}/v_{F}$ and $\bm{q}^{2}\ll 1$.

Appendix D Explicit matrix forms of the magnetic driving term $D_{B}$

In this Appendix, we show the explicit matrix forms of the magnetic driving term (27) obtained from arbitrary diagonal and off-diagonal density matrices. For concreteness, let us consider the case of a two-band model with the two eigenvalues $\varepsilon^{\pm}$ and the density matrix $\langle\rho\rangle$. We define the diagonal and off-diagonal parts of $\langle\rho\rangle$ as $\langle n\rangle$ and $\langle S\rangle$, respectively. The total density matrix of the system is given by $\langle\rho\rangle=\langle n\rangle+\langle S\rangle$. In the eigenstate basis, we have

[TABLE]

where $\mu$ is a chemical potential. We immediately get

[TABLE]

where $\bm{\mathcal{R}}=\sum_{\alpha=x,y,z}\mathcal{R}_{\alpha}\bm{e}_{\alpha}$ is the Berry connection vector.

First, let us calculate explicitly the magnetic driving term originating from the diagonal part of the density matrix $\langle n\rangle$. For concreteness, we consider the case of $\bm{B}=(0,0,B_{z})$. Then the driving term (27) is written as

[TABLE]

After a straightforward calculation, we obtain

[TABLE]

[TABLE]

which results in

[TABLE]

where

[TABLE]

Here, note that $g^{++}_{xy}=g^{++}_{yx}$ and $g^{--}_{xy}=g^{--}_{yx}$. In the case of $n_{\bm{k}}^{\pm}=f_{0}(\varepsilon_{\bm{k}}^{\pm})$ where $f_{0}(\varepsilon_{\bm{k}}^{\pm})$ is the Fermi-Dirac distribution function, we have $\frac{\partial\varepsilon_{\bm{k}}^{\pm}}{\partial k_{y}}\frac{\partial f_{0}(\varepsilon_{\bm{k}}^{\pm})}{\partial k_{x}}-\frac{\partial\varepsilon_{\bm{k}}^{\pm}}{\partial k_{x}}\frac{\partial f_{0}(\varepsilon_{\bm{k}}^{\pm})}{\partial k_{y}}=0$. Namely, the magnetic driving term (161) is purely off-diagonal in this case. Furthermore, in the case of $\varepsilon_{\bm{k}}^{\pm}=\pm\varepsilon_{\bm{k}}$, which applies to Weyl semimetals, we have $g^{+-}=g^{-+}=0$. Then in such a case the magnetic driving term is simplified to be

[TABLE]

Next, let us calculate explicitly the magnetic driving term originating from the off-diagonal part of the density matrix $\langle S\rangle$. For concreteness, we consider the case of $\bm{B}=(0,0,B_{z})$. Then the driving term (27) is written as

[TABLE]

After a straightforward calculation, we obtain

[TABLE]

[TABLE]

where $\Delta=\mathcal{R}_{x}^{++}-\mathcal{R}_{x}^{--}$ and $F=(\mathcal{R}_{x}^{+-}b-\mathcal{R}_{x}^{-+}a)$. In the case of $\varepsilon_{\bm{k}}^{\pm}=\pm\varepsilon_{\bm{k}}$, which applies to Weyl semimetals, we have

[TABLE]

with

[TABLE]

Especially in the case of $a_{\bm{k}}=-ic_{\bm{k}}$ and $b_{\bm{k}}=ic_{\bm{k}}$, i.e., $\langle S\rangle_{\bm{k}}=c_{\bm{k}}\sigma_{y}$, we have

[TABLE]

Finally, we obtain the magnetic driving term which is purely diagonal,

[TABLE]

with $\mathcal{F}_{\bm{k}}=\mathcal{A}_{\bm{k}}-\mathcal{B}_{\bm{k}}$.