Semiclassical dynamics and nonlinear charge current

Abstract

Electron conductivity is an important material property that can provide a wealth of information about the underlying system. Especially, the response of the conductivity with respect to electromagnetic fields corresponds to various nonlinear charge currents, which have distinct symmetry requirements and hence can be used as efficient probes of different systems. To help the band-structure engineering of such nonlinear currents, a universal treatment of electron dynamics up to second order expressed in the basis of the unperturbed states are highly useful. In this work, we review the general semiclassical framework of the nonlinear charge currents.

Graphical abstract

Keywords

nonlinear charge current / nonlinear electron dynamics / nonlinear anomalous Hall effect / linear magnetoresistance / negative longitudinal magnetoresistance / Berry curvature / quantum metric / positional shift

1 Introduction

Electron transport has long been a central topic in solid state physics, not only due to its fascinating complexity but also because it can provide rich information about the underlying electronic structure and dynamics. In the most simple Drude theory of metals for the diffusive transport, the conductivity is determined by the density, effective mass of carriers, and the phenomenological relaxation time [1]. Such simple picture can be generalized using the dynamics of Bloch electrons, to take into account the real energy spectrum and the Fermi surface geometry, as well as the microscopic scattering process [1].

In recent years, this classic paradigm of electron transport receives profound modifications, due to a deeper understanding of the geometric structure in the momentum space induced by the wave functions instead of the energy spectrum. The most well-known examples are various Hall effects [2, 3]. Central to the Hall effect is the concept of the Berry curvature. Its integration yields the Chern number, responsible for a topological classification of matter [4]. Moreover, it directly modifies the electron motion and hence is indispensable in describing the static and optical transport phenomenon [3, 5-7]. Berry curvature is also the analog of the magnetic field in the momentum space. As such, the source of Berry curvature is the momentum-space analog of magnetic monopoles, which can be realized by Weyl/Dirac points [8, 9]. In three dimensions, such isolated monopoles in momentum space can response differently to electromagnetic fields and cause a characteristic strong negative longitudinal magnetoresistance [8, 9].

On the other hand, the charge current can be expanded with respect to the electromagnetic fields as follows:

\begin{array}{}\text{(1)}& J_i\text{}=\text{}{\sigma }_{ij}{E}_j+{\sigma }_{ijk}{E}_j{E}_k+{\sigma }_{(ij\mathrm{,}k)}{E}_j{B}_k+{\sigma }_{(ij\mathrm{,}kk)}{E}_j{B}_k{B}_k+\cdots \mathrm{.}\end{array}

Berry curvature offers a first order correction to the electron dynamics. As such, it can give an accurate description of {\sigma }_{ij} whose antisymmetric part is the anomalous Hall effect. In addition to this linear current, there are various nonlinear currents with coefficients satisfying different symmetry requirements and hence can be used to probe distinct systems. To further demonstrate the symmetry requirement, we can expand the last three coefficients in terms of transport relaxation time τ in the semiclassical diffusive transport regime [10]

\begin{array}{}\text{(2)}& {\sigma }_{ijk}\text{}=\text{}{\sigma }_{ijk\mathrm{,0}}+{\sigma }_{ijk\mathrm{,1}}\tau +{\sigma }_{ijk\mathrm{,2}}{\tau }^2\mathrm{,}\end{array}

\begin{array}{}\text{(3)}& {\sigma }_{(ij\mathrm{,}k)}\text{}=\text{}{\sigma }_{(ij\mathrm{,}k)\mathrm{,0}}+{\sigma }_{(ij\mathrm{,}k)\mathrm{,1}}\tau +{\sigma }_{(ij\mathrm{,}k)\mathrm{,2}}{\tau }^2\mathrm{,}\end{array}

\begin{array}{}\text{(4)}& {\sigma }_{(ij\mathrm{,}kk)}\text{}=\text{}{\sigma }_{(ij\mathrm{,}kk)\mathrm{,0}}\text{}+\text{}{\sigma }_{(ij\mathrm{,}kk)\mathrm{,1}}\tau \text{}+\text{}{\sigma }_{(ij\mathrm{,}kk)\mathrm{,2}}{\tau }^2\text{}+\text{}{\sigma }_{(ij\mathrm{,}kk)\mathrm{,3}}{\tau }^3\mathrm{.}\end{array}

The symmetry of different coefficients is summarized in Table 1. It is thus tempting to find how different coefficients depend on the geometric properties of Bloch electrons in the momentum space. This can provide experimental guide and help the band structure engineering of larger nonlinear charge currents. For this purpose, the electron dynamics should be extended to second order.

In this work, we will focus on the semiclassical treatment of the electron dynamics and review the effort of generalizing such theory up to second order and deriving the nonlinear currents in the Boltzmann transport framework. Our paper is organized as follows. In Section 2, we discuss the generalization of the semiclassical equations of motion up to second order through an appropriately constructed wave packet up to first order. We then discuss the validity of the semiclassical theory though the magnetoelectric coefficient and magnetic susceptibility. With the correct semiclassical theory up to second order, in Section 3, we extend the Boltzmann transport theory and derive its solution in steady-state up to third order. We then use this universal treatment to derive various nonlinear currents, including the electric-field-correction to the nonlinear Hall effect ({\sigma }_{jii}), linear magnetoresistance ({\sigma }_{(ij\mathrm{,}k)}), and quadratic magnetoresistance ({\sigma }_{(ij\mathrm{,}kk)}).

2 Semiclassical theory of Bloch electrons in electromagnetic fields

To derive various nonlinear charge currents, one will need the correct electron dynamics beyond linear order under electromagnetic fields. For this purpose, certain perturbation technique should be employed. However, the periodic nature of crystals causes difficulties for such process. We note that the perturbation from electromagnetic fields reads eE\cdot r+\frac{1}{2}e[A(r)\cdot \widehat{v}+\widehat{v}\cdot A(r)]. It contains the position operator, which is unbounded in extended systems. Therefore, direct application of the time-independent perturbation is inappropriate.

To circumvent such difficulty, for magnetoconductivty one can first obtain the eigenstate of Bloch electrons under magnetic field, which are Landau levels for effective continuum models or Hofstadter spectrum for tight-binding models. Then the static magnetoconductivity is obtained using the standard Kubo formula under an appropriate limit. The dependence of the conductivity on magnetic field is encoded in the quantization of original Bloch states under magnetic field. Although this method is standard and has wide applicability, it cannot give universal understanding of the physics behind the magnetoconductivity, as the quantization into Landau levels or the Hofstadter spectrum varies for different systems.

Alternatively, one can construct a phase space using an appropriate position variable R and momentum variable P from the Bloch states. The electron dynamics is then contained in the equations of motion in the phase space. If R and P are canonical variables, their dynamics is solely determined by the energy. Therefore, the whole procedure amounts to deriving an effective energy in terms of R and P under external fields. If R and P are non-canonical, their Poisson brackets are also needed, which determine the structure of the equations of motion. The magnetoconductivity is then derived in combination with the Boltzmann equation for the distribution function. Although this method treats external fields as perturbations, and is only derived up to the second order currently, it has the benefit that it yields analytical results independent of the model detail. In fact, all terms in the equations of motion are expressed using the unperturbed Bloch functions, and hence can directly access both the spectrum and the geometric structure in the unperturbed Hilbert space.

In this section, we follow the second method by treating Bloch electrons as wave packets, whose center of mass position and momentum naturally span the phase space. We will first sketch the construction of the wave packet and show how its internal structure modifies the electron dynamics in phase space at first order of external fields. A detailed discussion of both aspects can be found in Ref. [3]. We then extend such theory up to second order and derive the second order correction to equations of motion. The validity of the semiclassical theory can be confirmed through the magnetoelectric effect and the quantization of effective energies into Landau levels under magnetic field. Since all elements in the semiclassical dynamics are expressed using unperturbed Bloch states, they can be evaluated using the output from the first-principles codes.

2.1 Properties of the wave packet

The dynamics of Bloch electrons under uniform electromagnetic fields is governed by the following crystal Hamiltonian in non-relativistic limit

\begin{array}{}\text{(5)}& {\widehat{H}}_{\text{f}}={\widehat{H}}_0[\widehat{p}+eA(r);r]+\frac{g{\mu }_B}{\hslash }B\cdot \widehat{s}+eE\cdot r\mathrm{,}\end{array}

where {\widehat{H}}_0(\widehat{p};r) is the unperturbed Hamiltonian for periodic crystals with \widehat{p} being the momentum operator, g is the gyromagnetic ratio, {\mu }_B is the Bohr magneton, and \widehat{s} is the spin operator. The magnetic field enters through both the minimal coupling and the Zeeman coupling. Generally speaking, {\widehat{H}}_{\text{f}} does not respect the translational symmetry, unless the magnetic flux through a unit cell is a rational number times the flux quantum such that the magnetic translational symmetry is recovered [11-14].

Although {\widehat{H}}_{\text{f}} exactly determines the dynamics, to directly diagonalize it is a difficult task. Instead, we can simplify the problem by assuming that the solution to the Schrödinger equation (\text{i}\hslash {\partial }_t-{\widehat{H}}_{\text{f}})\psi =0 takes the ansartz of a wave packet, which is the essence of the semiclassical theory. To construct the wave packet, we first evaluate the electromagnetic potential in the full Hamiltonian at the center of mass position r_c of the wave packet, and hence recover the translational symmetry of the original unperturbed crystal. The resulting local Hamiltonian reads

\begin{array}{}\text{(6)}& {\widehat{H}}_c={\widehat{H}}_0[\widehat{p}+eA(r_c);r]+\frac{g{\mu }_B}{\hslash }B\cdot \widehat{s}+eE\cdot r_c\mathrm{,}\end{array}

Its eigenenergy forms the local Bloch bands {\epsilon }_n[p+(e/\hslash )A(r_c)]+eE\cdot r_c with n being the band index and \hslash p being crystal momentum. Its eigenfunction is the Bloch function {\text{e}}^{\text{i}p\cdot r}|u_n[p+eA(r_c)/\hslash ]\rangle.

In the following, we focus on a single band with index 0 and assume that it is well separated from all the other bands. When electromagnetic fields are weak, the dynamics of a Bloch electron starting from some state in this band will still be confined in the same band. Therefore, the wave packet is constructed as the superposition of Bloch states from the band 0 [15]

\begin{array}{}\text{(7)}& |W\rangle ={\displaystyle \int }\text{d}p{C}_0(p){\text{e}}^{\text{i}p\cdot r}|u_0[p+eA(r_c)/\hslash ]\rangle \mathrm{.}\end{array}

|W\rangle should be normalized, indicating {\displaystyle \int }\text{d}p|C_0{|}^2=1.

There are two constraints on the coefficient C₀. First, the wave packet is assumed to be sharply localized in the momentum space. This requires that the magnitude of C₀ satisfies |C_0{|}^2\approx \delta (p-p_c). Here \hslash p_c is the center of mass momentum of the wave packet, i.e.,

\begin{array}{}\text{(8)}& p_c=\langle W|p|W\rangle \mathrm{.}\end{array}

Secondly, the construction of the wave packet requires a local Hamiltonian around the center of mass position r_c. This r_c has to be determined in a self-consistent manner, which exerts a constraint on the phase of C₀ [15]:

\begin{array}{}\text{(9)}& r_c=\langle W|r|W\rangle ={\frac{\partial \gamma }{\partial p}|}_{p=p_c}+{\mathcal{A}}_0(k_c)\mathrm{,}\end{array}

where \gamma =-arg(C_0), {\mathcal{A}}_0(p)=\langle u_0(p)|\text{i}{\partial }_p|u_0(p)\rangle is the intraband Berry connection, and \hslash k_c=\hslash p_c+eA(r_c) is the gauge-independent physical momentum. These two constraints complete the construction of the wave packet. Although the wave packet is for Bloch states from a single band, it can be generalized to account for the general multiband case [16, 17].

2.2 Semiclassical dynamics up to first order

With the knowledge of the wave packet, we now derive its evolution. The coefficient C₀ in the wave packet can be determined through the variational principle, i.e., the least action principle. Since C₀ is a complex function, with independent magnitude and phase, determined by k_c (or equivalently, p_c) and r_c respectively, applying the variational principle to the Lagrangian will yield the dynamics of r_c and k_c, i.e., the phase space equations of motion.

By evaluating the Lagrangian under the wave packet up to first order, we obtain [15]

\begin{array}{}\text{(10)}& \begin{array}{ll}L\hfill & =\langle W|\text{i}\hslash {\partial }_t-{\widehat{H}}_0-{\widehat{H}}_1|W\rangle \hfill \\ \hfill & =-(r_c-{\mathcal{A}}_0)\cdot \hslash {\dot{k}}_c-\frac{1}{2}eB\times r_c\cdot {\dot{r}}_c-{\tilde{\epsilon }}_0\mathrm{,}\hfill \end{array}\end{array}

where {\widehat{H}}_1=\frac{e}{4}B\cdot [(r-r_c)\times \widehat{v}-\widehat{v}\times (r-r_c)]+eE\cdot (r-r_c) is the first order correction to the local Hamiltonian {\widehat{H}}_c and \tilde{\epsilon }={\epsilon }_0-B\cdot m+eE\cdot r_c with m=m_{orb}-(g{\mu }_B/\hslash )\langle u_0|\widehat{s}|u_0\rangle, is the modified Band energy after taking account of both the magnetic moment and electric dipole of the wave packet, which will be discussed later. Here we use the symmetric gauge for the vector potential.

The Berry connection {\mathcal{A}}_0 is an essential ingredient in the Lagrangian. In the absence of {\mathcal{A}}_0, the above Lagrangian reduces to the conventional canonical form for electrons under electromagnetic fields:

\begin{array}{}\text{(11)}& L=\left(\hslash k_c-\frac{1}{2}eB\times r_c\right)\cdot {\dot{r}}_c-{\tilde{\epsilon }}_0\mathrm{.}\end{array}

The appearance of {\mathcal{A}}_0 indicates that r_c and k_c are not canonical variables. To make this statement clearer, we make the following transformation [18]:

\begin{array}{}\text{(12)}& r_c\text{}=\text{}q+{\mathcal{A}}_0+\frac{1}{2}\frac{e}{\hslash }(B\text{}\times \text{}{\mathcal{A}}_0\cdot {\partial }_p){\mathcal{A}}_0\text{}+\text{}\frac{1}{2}\frac{e}{\hslash }{\text{Ω}}_0\text{}\times \text{}(B\text{}\times \text{}{\mathcal{A}}_0)\mathrm{,}\end{array}

\begin{array}{}\text{(13)}& k_c\text{}=\text{}p+\frac{1}{2}\frac{e}{\hslash }B\times q+\frac{e}{\hslash }B\times (r_c-q)\mathrm{,}\end{array}

where {\text{Ω}}_0={\nabla }_p\times {\mathcal{A}}_0 is the Berry curvature, which can be viewed as the geometrical Berry phase per unit area [19]. Here the argument of {\mathcal{A}}_0 and {\text{Ω}}_0 is p+\frac{e}{2\hslash }B\times q. Then the Lagrangian in Eq. (10) recovers the canonical form for p and q

\begin{array}{}\text{(14)}& L=\hslash p\cdot \dot{q}-{\tilde{\epsilon }}_0\mathrm{.}\end{array}

Eqs. (12) and (13) thus describe the connection between physical variables and canonical variables.

One consequence of the noncanonicality is that the phase space measure for the volume element \text{d}{r}_c\text{d}{k}_c has to change based on the Jacobian \partial (r_c\mathrm{,}{k}_c)/\partial (q\mathrm{,}p). The resulting phase space density of states reads [20]

\begin{array}{}\text{(15)}& \mathcal{D}=1+\frac{e}{\hslash }B\cdot {\text{Ω}}_0(k_c)\mathrm{.}\end{array}

Since the physical position and momentum take the crystal volume and the Brillouin zone as their range, it is natural to evaluate the statistical average of any operator in the physical phase space spanned by (r_c\mathrm{,}{k}_c). Therefore, the modified density of states \mathcal{D} should always be used if the magnetic field is involved.

The noncanonicality of r_c and k_c will also lead to a nontrivial symplectic form, and hence affects the structure of the equations of motion. From Eq. (10), the Euler-Lagrangian equations of motion yield the following phase space dynamics [15]

\begin{array}{}\text{(16)}& {\dot{r}}_c=\frac{\partial {\tilde{\epsilon }}_0}{\hslash \partial k_c}-{\dot{k}}_c\times {\text{Ω}}_0(k_c)\mathrm{,}\end{array}

\begin{array}{}\text{(17)}& \hslash {\dot{k}}_c=-eE-e{\dot{r}}_c\times B\mathrm{.}\end{array}

The second term in the velocity equation is the anomalous velocity arising from the noncanonicality. The above equations of motion have wide applicability in the static spin, heat, and charge transport phenomena [3] and optical phenomena [5-7, 21-24].

2.3 Electromagnetic dipoles

One direct application of the semiclassical theory is to derive the electromagnetic multipoles in crystals. As the wave packet is sharply localized in the momentum space, it has a finite width in the real space. The charge and current distributions across such spread need not to be uniform. To represent such internal anisotropy of the wave packet, it is natural to use electromagnetic multipoles. For illustration purpose, we will consider the electromagnetic dipoles.

The electric dipole of the wave packet reads \langle W|-er|W\rangle =-e{r}_c, with r_c given in Eq. (9). Under the periodic gauge [3, 25], the Bloch wave function is periodic in the momentum space, so is the coefficient C₀. Therefore, by integrating r_c in the Brillouin zone, the first term in Eq. (9) from the phase of C₀ can at most yield an integer multiples of 2\pi and hence can be ignored. The remaining term yields the electric polarization in crystals

\begin{array}{}\text{(18)}& P=-e{\displaystyle \int }\frac{\text{d}{k}_c}{{(2\pi )}^3}{\mathcal{A}}_0(k_c)\mathrm{.}\end{array}

Originally, the electric polarization is derived by first considering the adiabatic charge pumping current and then identifying it as the displacement current from the temporal variation of the polarization [26, 27]. Here our method yields the same electric polarization under the periodic gauge for Bloch functions.

The magnetic dipole can be evaluated as follows [15, 28]:

\begin{array}{}\text{(19)}& \begin{array}{ll}{m}_{\text{orb}}\hfill & =-e\langle W|\frac{1}{2}[(r-r_c)\times \widehat{v}-\widehat{v}\times (r-r_c)]|W\rangle \hfill \\ \hfill & =-\frac{e}{2}{\displaystyle \sum _{n\ne 0}}{\mathcal{A}}_{0n}\times v_{n0}\mathrm{,}\hfill \end{array}\end{array}

where {\mathcal{A}}_{0n}=\langle u_0|\text{i}{\partial }_{k_c}|u_n\rangle is the interband Berry connection and v_{n0}=\langle u_n|\widehat{v}|u_0\rangle is the interband velocity element. In this work, we use the convention that for intraband matrix elements, we only keep one band index. The magnetic moment m_orb represents how a wave packet couples to an external magnetic field, i.e., shifting the band energy through an effective Zeeman coupling -B\cdot m_{\text{orb}}. As such, m_orb contributes to the orbital magnetization as shown later. m_orb can usually be accessed through the optical activity of crystals. For example, in ferromagnetic materials, the cross-gap part of m_orb (with 0 standing for occupied bands and n for the unoccupied bands in Eq. (19)) can be measured through the f-sum rule for the circular dichroism [29, 30]. Moreover, in noncentrosymmetric metals, the natural optical activity in the semiclassical regime is determined by the integration of the momentum space dipole of the magnetic moment (defined as v_i{(m_{\text{orb}})}_j with v being the band velocity) on the Fermi surface [6, 7, 22].

Interestingly, if one starts from the relativistic Dirac Hamiltonian and construct a coherent wave packet for the upper two bands, the orbital magnetic moment of this wave packet naturally recovers the spin magnetic moment [18]. Mathematically, this is equivalent to the Foldy–Wouthuysen transformation that reduces the Dirac Hamiltonian to the non-relativistic Schrödinger Hamiltonian [31, 32]. Here the semiclassical theory offers a heuristic and clear picture for the emergence of spin.

To derive the orbital magnetization in crystals, one should evaluate the free energy under magnetic field. On one hand, the magnetic moment of the wave packet couples to the magnetic field and modifies the band energy. On the other hand, the phase space density of states is changed by the magnetic field. Taking both corrections into consideration, one obtains the following free energy

\begin{array}{}\text{(20)}& F={\displaystyle \int }\mathcal{D}\frac{\text{d}k}{{(2\pi )}^3}(-k_{B}T)\ln \left[1+\exp \left(\frac{{\tilde{\epsilon }}_0-\mu }{k_{B}T}\right)\right]\mathrm{.}\end{array}

Therefore, the orbital magnetization can be derived by taking the derivative of the free energy with respect to the magnetic field [33]

\begin{array}{}\text{(21)}& M=-\frac{\partial F}{\partial B}={\displaystyle \int }\frac{\text{d}k}{{(2\pi )}^3}\left(f_{0}m-\frac{e}{\hslash }{g}_e{\text{Ω}}_0\right)\mathrm{,}\end{array}

where g_e=-(k_{B}T)\ln [1+\exp (({\epsilon }_0-\mu )/(k_{B}T))] is the grand potential density, and f₀ is the equilibrium Fermi distribution function. Here and hereafter we will ignore the subscript c in the integration over k_c. This is the contribution of the orbital magnetization from band 0. The total contribution can be obtained by summing over the band index. The orbital magnetization thus contains two contributions. The first one is due to the relative motion inside the wave packet, or the self-rotation of the wave packet. The second one is due to the rotation of the wave packet as a whole, or the revolution of the wave packet. This same magnetization can also be obtained through the Wannier function approach [34, 35], the exact Hofstadter spectrum [36, 37], or the linear response theory [38, 39].

2.4 Semiclassical dynamics up to second order

Although the semiclassical dynamics in Eqs. (16) and (17) is very powerful, it cannot fully account for the second order response functions such as the magnetoresistance, magnetoelectric coefficient, and so on. For this purpose, we need to extend the above semiclassical dynamics up to second order.

The construction of the wave packet has to be modified accordingly. The external electromagnetic fields generally affect the eigenfunction. As a result, the wave packet should be the superposition of the true eigenstates instead of the unperturbed ones. We can expand the true eigenstate in the basis of the unperturbed ones, and obtain the following wave packet [40]

\begin{array}{}\text{(22)}& |W\rangle ={\displaystyle \int }\text{d}p{\text{e}}^{\text{i}p\cdot r}\left[C_0(p)|u_0\rangle +{\displaystyle \sum _{n\ne 0}}{C}_n(p)|u_n\rangle \right]\mathrm{.}\end{array}

Here for simplicity we drop the argument p+(e/\hslash )A(r_c) of |u_0\rangle and |u_n\rangle. If |W\rangle satisfies the Schrödinger equation, we have

\begin{array}{}\text{(23)}& \langle u_n|{\text{e}}^{-\text{i}p\cdot r}(\text{i}\hslash {\partial }_t-{\widehat{H}}_{\text{f}})|W\rangle =\mathrm{0,}\forall n\ne 0.\end{array}

Therefore, we use this as an additional constraint on the wave packet in Eq. (22), which determines the connection between C_n and C₀. Specifically, solving this constraint at the linear order of electromagnetic fields yields [40]

\begin{array}{}\text{(24)}& C_n=\frac{G_{n0}}{{\epsilon }_0-{\epsilon }_n}{C}_0-\frac{\text{i}}{2}\frac{e}{\hslash }[B\times (\text{i}{\partial }_p+{\mathcal{A}}_0-r_c)C_0]\cdot {\mathcal{A}}_{n0}\mathrm{,}\end{array}

where G_{n0}=-B\cdot M_{n0}+eE\cdot {\mathcal{A}}_{n0} with M_{n0}=\frac{e}{2}{\displaystyle {\sum }_{m\ne 0}}(v_{nm}+v_0{\delta }_{mn})\times {\mathcal{A}}_{m0}-g{\mu }_B{s}_{n0}/\hslash being the interband element of the total (spin plus orbital) magnetic moment. We comment that different from the discussion here, in Refs. [40] and [41], the spin magnetic moment has not been explicitly added.

Equation (24) indicates that the correction to the wave packet from the magnetic field has two contributions with distinct structures and hence different origins. It is useful to first analyze the nature of these two perturbations as they prelude the form of the second order correction to the phase space dynamics. The first term has the conventional form of the perturbation correction with the energy gap as the denominator. It is called the vertical mixing as it mixes Bloch states from different bands but at the same p point in the Brillouin zone. It also contains both adiabatic and nonadiabatic perturbation, as the wave packet |W\rangle is time-dependent through the argument p+(e/\hslash )A(r_c) of the Bloch state.

In comparison, the second term in C_n modifies the wave packet in the following way:

\begin{array}{}\text{(25)}& \begin{array}{l}{\displaystyle \int }\text{d}p{\text{e}}^{\text{i}p\cdot r}\frac{-\text{i}}{2}\frac{e}{\hslash }[B\times (\text{i}{\partial }_p+{\mathcal{A}}_0-r_c)C_0]\cdot {\displaystyle \sum _{n\ne 0}}{\mathcal{A}}_{n0}|u_n\rangle \\ ={\displaystyle \int }\text{d}p{\text{e}}^{\text{i}p\cdot r}\frac{e}{2\hslash }\left(B\times (r-r_c)\cdot \widehat{D}|u_0\rangle +B\cdot {\text{Ω}}_0|u_0\rangle \right)\mathrm{.}\end{array}\end{array}

where \widehat{D}={\partial }_p+\text{i}{\mathcal{A}}_0 is the covariant derivative that ensures the gauge invariance. Notice that the momentum argument of |u_0\rangle is p+(e/\hslash )A(r_c). The first term in Eq. (25) has the meaning of shifting the momentum to p+(e/\hslash )A(r). This is a unique property of the perturbation from the magnetic field, and it suggests that the correction to the wave function also obeys the Peierls substitution. The same property is also derived using the Moyal product and phase space formulation of quantum mechanics [42]. The second term in Eq. (25) makes the total correction normal to the original state {\text{e}}^{\text{i}p\cdot r}|u_0\rangle to eliminate redundancy in the perturbative correction. Since the second term in C_n effectively mixes Bloch states in the same band but at neighbouring p points, it is referred to as the horizontal mixing. Both corrections can be easily visualized in Fig. 1.

Using the modified wave packet, we obtain the following Lagrangian

\begin{array}{}\text{(26)}& L=-(r_c-{\mathcal{A}}_0^{\text{t}})\cdot \hslash {\dot{k}}_c-\frac{1}{2}eB\times r_c\cdot {\dot{r}}_c-{\tilde{\epsilon }}_0\mathrm{,}\end{array}

where {\mathcal{A}}_0^{\text{t}} is the Berry connection evaluated in the true eigen state instead of the unperturbed one. It can be put in the form of the original Berry connection plus a positional shift correction: {\mathcal{A}}_0^{\text{t}}={\mathcal{A}}_0+{\mathcal{A}}_0^{\text{'}}. Before discussing the expression of the positional shift {\mathcal{A}}_0^{\text{'}} and modified energy \tilde{\epsilon }, we comment that the Lagrangian up to second order has the exact same form with the previous one in Eq. (10). Consequently, the connection between physical and canonical variables also has the same form except that the Berry connection {\mathcal{A}}_0^{\text{t}} in the exact eigenstate should be used

\begin{array}{}\text{(27)}& r_c=q+{\mathcal{A}}_0^{\text{t}}+\frac{e}{2\hslash }(B\text{}\times \text{}{\mathcal{A}}_0^{\text{t}}\cdot {\partial }_p){\mathcal{A}}_0^{\text{t}}+\frac{e}{2\hslash }{\text{Ω}}_0^{\text{t}}\text{}\times \text{}(B\text{}\times \text{}{\mathcal{A}}_0^{\text{t}})\mathrm{,}\end{array}

\begin{array}{}\text{(28)}& k_c=p+\frac{e}{2\hslash }B\times q+\frac{e}{\hslash }B\times (r_c-q)\mathrm{,}\end{array}

where {\text{Ω}}_0^{\text{t}}={\nabla }_p\times {\mathcal{A}}_0^{\text{t}} is the Berry curvature evaluated using the exact eigenstate. As a result, the phase space density of states reads

\begin{array}{}\text{(29)}& \mathcal{D}=1+\frac{e}{\hslash }B\cdot {\text{Ω}}_0^{\text{t}}\mathrm{.}\end{array}

More importantly, the equations of motion up to second order also keep the same form [40]

\begin{array}{}\text{(30)}& {\dot{r}}_c=\frac{\partial {\tilde{\epsilon }}_0}{\hslash \partial k_c}-{\dot{k}}_c\times {\text{Ω}}_0^{\text{t}}(k_c)\mathrm{,}\end{array}

\begin{array}{}\text{(31)}& \hslash {\dot{k}}_c=-eE-e{\dot{r}}_c\times B\mathrm{.}\end{array}

The positional shift is gauge-independent. It represents the additional shift of the wave packet center after applying external electromagnetic fields

\begin{array}{}\text{(32)}& r_c=\langle W|r|W\rangle ={\frac{\partial \gamma }{\partial p}|}_{p=p_c}+{\mathcal{A}}_0(k_c)+{\mathcal{A}}_0^{\text{'}}(k_c)\mathrm{.}\end{array}

The positional shift at first order reads

\begin{array}{}\text{(33)}& {\mathcal{A}}_0^{\text{'}}=2{\displaystyle \sum _{n\ne 0}}\text{R}\text{e}\frac{{\mathcal{A}}_{0n}{G}_{n0}}{{\epsilon }_0-{\epsilon }_n}-\frac{1}{2}\frac{eB}{\hslash }{ϵ}_{ijk}{\gamma }_{ji\mathcal{l}}{\widehat{e}}_{\mathcal{l}}\mathrm{,}\end{array}

where {ϵ}_{ijk} is the Levi-Civita symbol and {\gamma }_{ji\mathcal{l}}=\frac{1}{2}({\partial }_{\mathcal{l}}{g}_{ij}+{\partial }_i{g}_{j\mathcal{l}}-{\partial }_j{g}_{i\mathcal{l}}) is the Christoffel symbol. Here and hereafter, we use the Einstein summation convention for repeated indices and the partial derivative is with respect to the crystal momentum unless otherwise specified. g_{ij} is a Fubini-Study metric tensor, usually called quantum metric [43-46]. It measures the distance between two neighbouring Bloch states in the momentum space. In fact, it is a part of a more general concept, the quantum geometrical tensor, which, for band 0, is defined as

\begin{array}{}\text{(34)}& {\mathcal{G}}_{ij}=\langle {\partial }_i{u}_0|{\partial }_j{u}_0\rangle -{({\mathcal{A}}_i)}_0{({\mathcal{A}}_j)}_0\mathrm{.}\end{array}

The real and imaginary part of {\mathcal{G}}_{ij} yield the quantum metric and Berry curvature, respectively [43, 44].

In the positional shift, the first term is due to the vertical mixing and the second term is due to the horizontal mixing. It is interesting that the horizontal mixing yields a purely geometric correction to the Berry connection. This is not coincidental, as the horizontal mixing addresses neighbouring Bloch states in the momentum space, whose difference is proportional to the Berry connection and whose distance is determined by the quantum metric. We comment that this positional shift is also envisioned in Ref. [47], using a similar technique.

The positional shift reflects the change in the electric dipole moment of the wave packet. When external fields are applied, the electric polarization has to pick up the following change. First, the density of states has to be changed to \mathcal{D}. Second, the Berry connection {\mathcal{A}}_0\to {\mathcal{A}}_0^{\text{t}}+\frac{e}{2\hslash }(B\times {\mathcal{A}}_0^{\text{t}}\cdot {\partial }_p){\mathcal{A}}_0^{\text{t}}+\frac{e}{2\hslash }{\text{Ω}}_0^{\text{t}}\times (B\times {\mathcal{A}}_0^{\text{t}}) according to Eq. (27). As a result, the first order change in electric polarization reads

\begin{array}{}\text{(35)}& \delta P=-e{\displaystyle \int }\frac{\text{d}k}{{(2\pi )}^3}\left[\frac{e}{2\hslash }({\text{Ω}}_0\cdot {\mathcal{A}}_0)B+{\mathcal{A}}_0^{\text{'}}\right]\mathrm{.}\end{array}

The first term is the Abelian Chern–Simons 3-form, which yields the topological part of the orbital magnetoelectric coefficient. The second term yields the electric polarizability and the cross-gap part of the orbital magnetoelectric polarizability, consistent with the calculation using the linear response theory [48]. Such consistency confirms the validity of the first order correction to the Berry phase.

The modified band energy in Eq. (26) should contain the correction up to second order, i.e., it can be put in the following form: {\tilde{\epsilon }}_0={\epsilon }_0-B\cdot m+eE\cdot r_c+{\epsilon }_{(2)}. The second order correction describes the change in the electric dipole and magnetic moment of the wave packet in response to external fields. It reads [41]

\begin{array}{}\text{(36)}& \begin{array}{ll}{\epsilon }_{(2)}=\hfill & \frac{1}{4}\frac{e}{\hslash }(B\cdot {\text{Ω}}_0)(B\cdot m)-18\frac{e^2}{{\hslash }^2}{ϵ}_{sik}{ϵ}_{tj\mathcal{l}}{B}_s{B}_t{g}_{ij}{\alpha }_{k\mathcal{l}}\hfill \\ \hfill & +\nabla \cdot P_{\text{E}}+{\displaystyle \sum _{n\ne 0}}\frac{G_{0n}{G}_{n0}}{{\epsilon }_0-{\epsilon }_n}\hfill \\ \hfill & -eB\cdot ({\mathcal{A}}_0^{\text{'}}\times v_0)+\frac{e^2}{8m}(B^2{g}_{ii}-B_i{g}_{ij}{B}_j)\mathrm{,}\hfill \end{array}\end{array}

where {\alpha }_{k\mathcal{l}}={\partial }_k{\partial }_{\mathcal{l}}{\epsilon }_0/{\hslash }^2 is the inverse effective mass tensor, m is the free electron mass, G_{0n}=(G_{n0}{)}^{\star }, and

\begin{array}{}\text{(37)}& P_{\text{E}}=\frac{1}{2}e{\displaystyle \sum _{n\ne 0}}\text{R}\text{e}[(B\times {\mathcal{A}}_{0n})G_{n0}]\mathrm{.}\end{array}

In the energy correction, the first two terms are purely due to the horizontal mixing and hence are geometrical corrections involving Berry curvature and quantum metric. The third term is due to the cross effect of the horizontal and vertical mixing. P_E is called energy polarization as it accounts for the energy dipole in the momentum space. The fourth term is purely due to the horizontal mixing and has the conventional form of the perturbation correction to the energy. The fifth term shifts the momentum argument of {\epsilon }_0. The last term is the expectation of the second order correction to the local Hamiltonian.

We comment that if one starts from an effective tight-binding Hamiltonian, the last term of Eq. (36) will change. Instead of the free electron mass, one has to use the Hessian operator {\widehat{\text{Γ}}}_{ij}={\partial }^2{\widehat{H}}_0/\partial {\widehat{p}}_i\partial {\widehat{p}}_j. Then the last term in Eq. (36) should be replaced by the following two terms

\begin{array}{}\text{(38)}& \begin{array}{ll}{\epsilon }_{\text{H}}=\hfill & -\frac{e^2}{16}{(B\times \partial )}_i{(B\times \partial )}_j\langle u_0|{\widehat{\text{Γ}}}_{ij}|u_0\rangle \hfill \\ \hfill & +\frac{e^2}{8}{\displaystyle \sum _{(m\mathrm{,}n)\ne 0}}{(B\text{}\times \text{}{\mathcal{A}}_{0m})}_i{({\widehat{\text{Γ}}}_{ij})}_{mn}{(B\text{}\times \text{}{\mathcal{A}}_{n0})}_j\mathrm{.}\hfill \end{array}\end{array}

Here {({\widehat{\text{Γ}}}_{ij})}_{mn}=\langle u_m|{\widehat{\text{Γ}}}_{ij}|u_n\rangle.

2.5 Landau level quantization and magnetic susceptibility

If one starts from a continuum tight-binding Hamiltonian, electronic states will be quantized into Landau levels under magnetic field. Therefore, as long as the semiclassical theory can successfully reproduce the Landau level quantization up to second order, it is enough to account for various response functions at second order in magnetic field.

There are various quantization rules for Landau levels, such as the Onsager’s rule [49], the Maslov canonical operator method [50], the Gutzwiller trace formula [51], and so on. For semiclassical theory, the most relevant one is the Onsager’s rule [49], which reads

\begin{array}{}\text{(39)}& S({\epsilon }_n)=2\pi \left(n+\frac{1}{2}\right)\frac{eB}{\hslash }\mathrm{,}\end{array}

where S({\epsilon }_n) is the area in the momentum space enclosed by an equal-energy contour. Later, this quantization rule is generalized up to first order to take account of corrections from Berry phase and magnetic moment [3, 28, 52-56].

To further generalize the Onsager’s rule, one observes that the area S({\epsilon }_n) is proportional to the total electron density below {\epsilon }_n. Since the area S as a function of energy is a smooth function, such electron density should be understood as the smooth semiclassical electron density {\rho }_{\text{semi}}. Then it has been proved that the Onsager’s rule can be replaced by the following density quantization rule [57]

\begin{array}{}\text{(40)}& {\rho }_{\text{semi}}({\epsilon }_n)=\left(n+\frac{1}{2}\right)\frac{eB}{h}\mathrm{.}\end{array}

This suggests that the smooth semiclassical electron density always intersects with the true stepwise electron density at the half-filling points, as shown in Fig. 2.

In general {\rho }_{\text{semi}} can be expanded near B=0 in power series of B. Up to second order we have

\begin{array}{}\text{(41)}& \left(n+\frac{1}{2}\right)eBh=\frac{S({\epsilon }_n)}{4{\pi }^2}+B{\frac{\partial M}{\partial \mu }|}_{\mu ={\epsilon }_n}+\frac{1}{2}{B}^2{\frac{\partial \chi }{\partial \mu }|}_{\mu ={\epsilon }_n}\mathrm{,}\end{array}

where M and χ are the magnetization and magnetic susceptibility respectively, evaluated at zero magnetic field and zero temperature. It is interesting to note that Eq. (41) indicates the nonlinearity in the Landau level fan diagram when plotted as B against 1/n due to the appearance of the susceptibility. Such nonlinearity is consistent with experiments [58-64]. Based on Eq. (41), as long as the semiclassical theory can yield the correct susceptibility, it can give correct second order responses in magnetic field.

To evaluate the magnetic susceptibility, we need to calculate the free energy up to second order. In the semiclassical framework, the free energy can be expressed as follows:

\begin{array}{}\text{(42)}& F={\displaystyle \int }\frac{\text{d}{k}_c}{8{\pi }^3}\mathcal{D}(g_e(\tilde{\epsilon })+g_{\text{L}})\mathrm{.}\end{array}

Here the second term g_{\text{L}} is the Peierls–Landau magnetic free energy: g_{\text{L}}=-(e^2{f}_0^{\prime }/48{\hslash }^2)B_s{B}_t{ϵ}_{sik}{ϵ}_{tj\mathcal{l}}{\alpha }_{ij}{\alpha }_{k\mathcal{l}}, where f_0^{\prime } is the energy derivative of the Fermi distribution function f₀. For isotropic bands, the effective mass tensor α is diagonal, and g_{\text{L}} will reduce to its familiar form [1]. This term originates from higher order corrections to the replacement of the quantum mechanical commutator with the classical Poisson bracket [42].

After collecting terms at second order in B, one can find the following free energy F={\displaystyle \int }{g}^{″}[\text{d}k/(8{\pi }^3)] with

\begin{array}{}\text{(43)}& \begin{array}{ll}{g}^{″}\hfill & =g_{\text{L}}+\frac{f_0^{\prime }}{2}{(B\cdot m)}^2-f_0^{\prime }{v}_0\cdot P_{\text{E}}\hfill \\ \hfill & +f_0{\displaystyle \sum _{n\ne 0}}\frac{G_{0n}{G}_{n0}}{{\epsilon }_0-{\epsilon }_n}+\frac{e^2{f}_0}{8m}(B^2{g}_{ii}-B_i{g}_{ij}{B}_j)\hfill \\ \hfill & -\frac{3e{f}_0}{4\hslash }(B\cdot {\text{Ω}}_0)(B\cdot m)-\frac{e^2{f}_0}{8{\hslash }^2}{ϵ}_{sik}{ϵ}_{tj\mathcal{l}}{B}_s{B}_t{g}_{ij}{\alpha }_{k\mathcal{l}}\mathrm{.}\hfill \end{array}\end{array}

The magnetic susceptibility can be readily obtained through its definition

\begin{array}{}\text{(44)}& {\chi }_{ij}=-{\frac{{\partial }^{2}F}{\partial B_i\partial B_j}|}_{\mu \mathrm{,}T\mathrm{,}B=0}\mathrm{.}\end{array}

The first two terms in Eq. (43) are Peierls–Landau diamagnetic and Pauli paramagnetic contributions. The physical meaning of the other terms in Eq. (43) can be illustrated by re-expressing it in the Wannier function basis and then taking the atomic insulator limit. In particular, the two terms in the second line reduce to the familiar form of the Van-Vleck paramagnetic free energy and Langevin diamagnetic free energy in atomic physics, respectively. The last term is purely geometrical as it solely comes from the horizontal mixing. The energy polarization and geometrical magnetic free energy approach zero in the atomic limit, indicating that they have no analog in atomic physics and are novel terms in crystals.

The susceptibility obtained from Eq. (43) in the semiclassical framework is consistent with the result derived using other methods. For example, the susceptibility has been calculated in the honeycomb lattice using Hofstadter spectrum or linear response theory [65, 66]. Eq. (43) yields the exactly same result as shown in Fig. 3. In Ref. [67], it is also proved that Eq. (43) is consistent with the susceptibility from linear response theory [68]. This confirms the validity of the semiclassical dynamics under magnetic field.

2.6 Evaluating positional shift and energy correction in tight-binding models

The positional shift and energy correction can be implemented in first-principles codes. For this purpose, besides the energy spectrum, one needs three additional matrix elements: the velocity matrix element v_{mn}=\langle m|\widehat{v}|n\rangle, the Hessian matrix element {({\text{Γ}}_{ij})}_{mn}=\langle u_m|{\widehat{\text{Γ}}}_{ij}|u_n\rangle, and the spin Pauli matrix element {\sigma }_{mn}=\langle m|\widehat{\sigma }|n\rangle, all of which can be evaluated in first-principles codes in principle [69-72]. In the following, we sketch the process of re-expressing the positional shift and energy correction in terms of the energy spectrum and these three matrix elements.

We first examine the positional shift. The interband part of the Berry connection is related to the interband velocity element:

\begin{array}{}\text{(45)}& v_{mn}=\frac{\text{i}}{\hslash }({\epsilon }_m-{\epsilon }_n){\mathcal{A}}_{mn}\mathrm{,}\forall m\ne n\mathrm{.}\end{array}

Using this identity, we manipulate the interband mixing element G_{n0}

\begin{array}{}\text{(46)}& G_{n0}=\frac{\text{i}eB\hslash }{2}\cdot {\displaystyle \sum _{m\ne 0}}\frac{(v_{nm}+v_0{\delta }_{mn})\times v_{m0}}{{\epsilon }_m-{\epsilon }_0}+\frac{g{\mu }_B}{2}B\cdot {\sigma }_{n0}\mathrm{.}\end{array}

Then the first term in the positional shift can be put in the desired form

\begin{array}{}\text{(47)}& 2\hslash \text{I}m{\displaystyle \sum _{n\ne 0}}\frac{v_{0n}{G}_{n0}}{{({\epsilon }_0-{\epsilon }_0)}^2}\mathrm{.}\end{array}

For the remaining term we note that

\begin{array}{}\text{(48)}& \begin{array}{lll}{ϵ}_{ijk}{\gamma }_{ji\mathcal{l}}\hfill & =\hfill & {ϵ}_{ijk}{\partial }_i{g}_{j\mathcal{l}}\hfill \\ \hfill & =\hfill & -2{\hslash }^3{ϵ}_{ijk}{\displaystyle \sum _{n\ne 0}}\mathrm{Re}\frac{[(v_i{)}_0-{(v_i)}_n](v_j{)}_{0n}{(v_{\mathcal{l}})}_{n0}}{{({\epsilon }_0-{\epsilon }_n)}^3}\hfill \\ \hfill & \hfill & +{\hslash }^2{ϵ}_{ijk}{\displaystyle \sum _{n\ne 0}}\mathrm{Re}\frac{{\partial }_i{(v_j)}_{0n}{(v_{\mathcal{l}})}_{n0}+(j\leftrightarrow \mathcal{l})}{{({\epsilon }_0-{\epsilon }_n)}^2}\mathrm{.}\hfill \end{array}\end{array}

The only unknown quantity has the form {\partial }_i{(v_j)}_{mn}, which can be manipulated as follows:

\begin{array}{}\text{(49)}& \begin{array}{ll}{\partial }_i{(v_j)}_{mn}\text{}=\hfill & {\displaystyle \sum _{m^{\prime }\ne m}}\frac{\hslash {(v_i)}_{m{m}^{\prime }}{(v_j)}_{m^{\prime }n}}{{\epsilon }_m-{\epsilon }_{m^{\prime }}}\text{}-\text{}{\displaystyle \sum _{m^{\prime }\ne n}}\frac{\hslash {(v_j)}_{m{m}^{\prime }}{(v_i)}_{m^{\prime }n}}{{\epsilon }_{m^{\prime }}-{\epsilon }_n}\hfill \\ \hfill & +\hslash {({\text{Γ}}_{ij})}_{mn}+\text{i}[({\mathcal{A}}_i{)}_m-{({\mathcal{A}}_i)}_n](v_j{)}_{mn}\mathrm{.}\hfill \end{array}\end{array}

The first three terms are in the desired form. The last term contains the intraband Berry connection which is gauge-dependent. Therefore, it cannot be effectively evaluated in the first-principles codes. However, it can be shown that it does not contribute to Eq. (48) and hence need not to be evaluated. Equations (46) to (49) are enough to transform the positional shift in the desired form.

We now examine the energy correction. Using the following identities

\begin{array}{}\text{(50)}& {\text{Ω}}_0=-{\hslash }^2{\displaystyle \sum _{n\ne 0}}\mathrm{Im}\frac{v_{0n}\times v_{n0}}{{({\epsilon }_0-{\epsilon }_n)}^2}\mathrm{,}\end{array}

\begin{array}{}\text{(51)}& m=-\frac{e\hslash }{2}{\displaystyle \sum _{n\ne 0}}\mathrm{Im}\frac{v_{0n}\times v_{n0}}{{\epsilon }_0-{\epsilon }_n}-\frac{g{\mu }_B{\sigma }_0}{2}\mathrm{,}\end{array}

\begin{array}{}\text{(52)}& g_{ij}={\hslash }^2{\displaystyle \sum _{n\ne 0}}\mathrm{Re}\frac{{(v_i)}_{0n}{(v_j)}_{n0}}{{({\epsilon }_0-{\epsilon }_n)}^2}\mathrm{,}\end{array}

\begin{array}{}\text{(53)}& {\alpha }_{ij}=({\text{Γ}}_{ij}{)}_0+2{\displaystyle \sum _{n\ne 0}}\text{R}\text{e}\frac{{(v_i)}_{0n}{(v_j)}_{n0}}{{\epsilon }_0-{\epsilon }_n}\mathrm{,}\end{array}

together with Eqs. (45) to (49), it is straightforward to put all terms except the energy polarization contribution in Eq. (36) in the desired form. We comment that for the tight-binding Hamiltonian, one may have to substitute the last term in Eq. (36) with Eq. (38) which has realistic Hessian matrix. In this case, we have

\begin{array}{}\text{(54)}& \begin{array}{ll}{\epsilon }_{\text{H}}=\hfill & \frac{e^2{\hslash }^2}{8}{\displaystyle \sum _{(m\mathrm{,}n)\ne 0}}\mathrm{Re}\frac{{(B\times v_{0m})}_i{({\widehat{\text{Γ}}}_{ij})}_{mn}{(B\times v_{n0})}_j}{({\epsilon }_0-{\epsilon }_m)({\epsilon }_n-{\epsilon }_0)}\hfill \\ \hfill & +\frac{e^2}{8}{ϵ}_{sik}{ϵ}_{tj\mathcal{l}}{B}_s{B}_t{\displaystyle \sum _{n\ne 0}}\mathrm{Im}[({\text{Γ}}_{ij}{)}_{n0}{\partial }_{\mathcal{l}}{({\mathcal{A}}_k)}_{0n}]\hfill \\ \hfill & +\frac{e^2}{8}{ϵ}_{sik}{ϵ}_{tj\mathcal{l}}{B}_s{B}_t{\displaystyle \sum _{n\ne 0}}\mathrm{Im}[({\mathcal{A}}_k{)}_{0n}{\partial }_{\mathcal{l}}{({\text{Γ}}_{ij})}_{n0}]\mathrm{.}\hfill \end{array}\end{array}

To evaluate the last two terms, one notes that

\begin{array}{}\text{(55)}& {\partial }_{\mathcal{l}}{({\mathcal{A}}_k)}_{0n}=-\frac{\text{i}{\hslash }^2[(v_{\mathcal{l}}{)}_0-{(v_{\mathcal{l}})}_n](v_k{)}_{0n}}{{({\epsilon }_0-{\epsilon }_n)}^2}+\frac{\text{i}\hslash {\partial }_{\mathcal{l}}{(v_k)}_{0n}}{{\epsilon }_0-{\epsilon }_n}\mathrm{,}\end{array}

\begin{array}{}\text{(56)}& \begin{array}{ll}{\partial }_{\mathcal{l}}{({\text{Γ}}_{ij})}_{n0}\hfill & ={\displaystyle \sum _{m\ne n}}\frac{\hslash {(v_{\mathcal{l}})}_{nm}{({\text{Γ}}_{ij})}_{m0}}{{\epsilon }_n-{\epsilon }_m}-{\displaystyle \sum _{m\ne 0}}\frac{\hslash {({\text{Γ}}_{ij})}_{nm}{(v_{\mathcal{l}})}_{m0}}{{\epsilon }_m-{\epsilon }_0}\hfill \\ \hfill & +{({\partial }_{\mathcal{l}}{\text{Γ}}_{ij})}_{n0}+\text{i}[({\mathcal{A}}_{\mathcal{l}}{)}_n-{({\mathcal{A}}_{\mathcal{l}})}_0]({\text{Γ}}_{ij}{)}_{n0}\mathrm{.}\hfill \end{array}\end{array}

Using Eq. (55) and Eq. (49), one can put the second term in {\epsilon }_{\text{H}} in the desired form plus extra terms from intraband Berry connection. After inserting Eq. (56) into {\epsilon }_{\text{H}}, the first term in the second line of Eq. (56) does not contribute due to the anti-symmetrization of the spatial indices. The last term in Eq. (56) will cancel the extra gauge-dependent terms in Eq. (55). This completes the manipulation of {\epsilon }_{\text{H}}.

The manipulation of energy polarization contribution in Eq. (36) can be done as follows:

\begin{array}{}\text{(57)}& \nabla \cdot P_{\text{E}}=\frac{1}{2}e{B}_i{ϵ}_{\mathcal{l}ij}\text{R}\text{e}{\displaystyle \sum _{n\ne 0}}[{\partial }_{\mathcal{l}}{({\mathcal{A}}_j)}_{0n}{G}_{n0}+{({\mathcal{A}}_j)}_{0n}{\partial }_{\mathcal{l}}{G}_{n0}]\mathrm{,}\end{array}

where

\begin{array}{}\text{(58)}& \begin{array}{ll}{\partial }_{\mathcal{l}}{G}_{n0}\hfill & =-\frac{eB}{2}{\displaystyle \sum _{m\ne 0}}[{\partial }_{\mathcal{l}}(v_{nm}+v_0{\delta }_{nm})]\times {\mathcal{A}}_{m0}\hfill \\ \hfill & -\frac{eB}{2}{\displaystyle \sum _{m\ne 0}}(v_{nm}+v_0{\delta }_{nm})\times ({\partial }_{\mathcal{l}}{\mathcal{A}}_{m0})\hfill \\ \hfill & +\frac{g{\mu }_B\hslash B}{2}\cdot \left[{\displaystyle \sum _{m\ne n}}\frac{{(v_{\mathcal{l}})}_{nm}{\sigma }_{m0}}{{\epsilon }_n\text{}-\text{}{\epsilon }_m}\text{}-\text{}{\displaystyle \sum _{m\ne 0}}\frac{{(\sigma )}_{nm}{(v_{\mathcal{l}})}_{m0}}{{\epsilon }_m\text{}-\text{}{\epsilon }_0}\right]\hfill \\ \hfill & +\frac{g{\mu }_B\hslash B\cdot {\sigma }_{n0}}{2}\text{i}[({\mathcal{A}}_{\mathcal{l}}{)}_n-{({\mathcal{A}}_{\mathcal{l}})}_0]\mathrm{.}\hfill \end{array}\end{array}

Again, terms contain intraband Berry connection can be ignored as they will cancel each other in the end.

We comment that if one starts from the Schrödinger Hamiltonian, the Hessian matrix can be directly evaluated without using wave functions in first-principles codes. For Schrodinger Hamiltonian, {\widehat{\text{Γ}}}_{ij}=1/(m{\delta }_{ij}). Therefore, one has {({\text{Γ}}_{ij})}_{mn}=1/(m{\delta }_{ij}{\delta }_{mn}). Then for the semiclassical dynamics, one only needs to evaluate the velocity and spin matrix elements using output from first-principles codes.

As a concrete example, we consider the following tight-binding Hamiltonian describing a 3D TI in the Bi₂Se₃ family on a simple cubic lattice [4, 73, 74]

\begin{array}{}\text{(59)}& {\widehat{H}}_{\text{tb}}=2\lambda {\sigma }_z(s_x\sin k_y-s_y\sin k_x)+2{\lambda }_z{\sigma }_y\sin k_z+{\sigma }_x{M}_k\mathrm{,}\end{array}

where σ and s are Pauli matrices in orbital and spin space respectively, and M_k=ϵ-2t(\cos k_x+\cos k_y+\cos k_z). At ϵ=6t, the above Hamiltonian is at the phase transition point between trivial phase and topological insulator phase, and supports a Dirac point at the Γ point. The Dirac point can split into Weyl point if the following additional term is added to the Hamiltonian

\begin{array}{}\text{(60)}& {\widehat{H}}_{\text{break}}=b_0{\sigma }_y{s}_z+b_z{s}_z\mathrm{.}\end{array}

These two terms break inversion and time reversal symmetry, respectively. The resulting Weyl semimetal phase is displayed in Fig. 4.

Here we choose the semiconductor phase and calculate the first and second order correction to the band energy, as shown in Fig. 5. It can be seen that both first and second order correction vary drastically near the band crossing points and the small band gap. In fact, both corrections diverge at the Γ point due to the band crossing. Across the small global band gap as shown in Fig. 4, both first and second order corrections change sign. It is also interesting to note that near the Γ point, the lowest two bands have similar first order corrections, but are opposite to those of the remaining two bands. This property is consistent with the calculation of magnetic moment for the low-energy Dirac model. In comparison, the second order correction have different properties. The lowest and the second highest band have similar second order corrections, opposite to those of the remaining two bands. This can be most easily seen from Eq. (36). The Berry curvature, effective mass tensor and energy gap reverse sign for the two bands forming a Dirac cone, leading to the sign change in their second order energy correction.

3 Nonlinear charge current

With the semiclassical theory up to second order at hand, we can derive the nonlinear charge current up to third order, i.e. conductivity up to second order in fields. We focus on the semiclassical regime with weak electromagnetic fields. There is one important issue that we want to discuss before the detailed derivation of nonlinear currents. In the previous section, the semiclassical theory is generally derived for a single Bloch band that is well separated from all the other bands. In reality, this condition is hardly met as band-crossing points are generally present. In this case, the previous semiclassical theory can still be used to derive the conductivity as long as the Fermi surface is not close to those band-crossing points. The reason is as follows. As a perturbation theory, the semiclassical theory may fail when the band gap is small, such as near the band crossing points. However, the conductivity is a Fermi surface property. As a result, as long as the band-crossing points are deep inside the Fermi sea, the perturbation theory still works near the Fermi surface.

In this section, we will first establish the general theory of the charge current beyond the linear order. Then we will discuss several important examples, including the response of the anomalous Hall conductivity to electromagnetic fields, linear magnetoresistance in time-reversal-broken materials, chiral anomaly in Weyl semimetals and two different mechanisms for the negative longitudinal magnetoresistance: the intrinsic contribution and the current jetting.

3.1 General theory of nonlinear charge currents

The current in the semiclassical theory reads

\begin{array}{}\text{(61)}& J=-e{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\mathcal{D}\dot{r}f\mathrm{,}\end{array}

where f is the electron distribution function, which is the Fermi function in equilibrium. Here and hereafter, we drop the subscript c in the center of mass position r_c and momentum k_c for simplicity. From the semiclassical dynamics in Eq. (30) and (31), one can show that

\begin{array}{}\text{(62)}& \mathcal{D}\dot{r}=\tilde{v}+eE\times {\text{Ω}}^{\text{t}}+\frac{e}{\hslash }\left(\tilde{v}\cdot {\text{Ω}}^{\text{t}}\right)B\mathrm{,}\end{array}

where \tilde{v}=\partial \tilde{\epsilon }/\hslash \partial k is the modified band velocity. Here and hereafter, we ignore the band index 0 in relevant intraband quantities for simplicity and the subscript 0 only has the meaning of zeroth order.

The remaining factor in the current is the distribution function f, which is typically solved from the Boltzmann equation. If the sample is homogeneous and reaches a steady state under external fields, under the relaxation time approximation, the Boltzmann equation reads

\begin{array}{}\text{(63)}& \dot{k}\cdot \frac{\partial f}{\partial k}={\frac{\text{d}f}{\text{d}t}|}_{\text{collision}}=-\frac{f-f_0}{\tau }\mathrm{.}\end{array}

Here the argument of equilibrium distribution f₀ is the modified band energy \tilde{\epsilon }. Solving the Boltzmann equation perturbatively, one has [75]

\begin{array}{}\text{(64)}& f={\displaystyle \sum _{m=0}^{\infty }}{(-\tau \dot{k}\cdot {\partial }_k)}^m{f}_0(\tilde{\epsilon })\mathrm{.}\end{array}

Equations (61), (62) and (64) are enough to derive the current except the anomalous Hall current up to third order of the external fields, i.e., the conductivity up to second order. Other than the anomalous Hall current, the equilibrium distribution f₀ alone does not contribute to the current. Therefore, at least the first order correction to f₀ is required, which automatically contains an external field due to \dot{k}. As a result, one only needs to evaluate the factor \mathcal{D}\dot{r} up to second order, which is given in the semiclassical dynamics. For the anomalous Hall conductivity, the semiclassical theory can only yields its first order correction, as only the first order correction to Berry curvature is obtained in the semiclassical theory up to second order.

To explicitly write out the solution in Eq. (64) up to third order in fields, we note that

\begin{array}{}\text{(65)}& -\tau \dot{k}=\frac{\tau }{\hslash }\frac{eE+e\tilde{v}\times B+\frac{e^2}{\hslash }(E\cdot B){\text{Ω}}^{\text{t}}}{1+\frac{e}{\hslash }B\cdot {\text{Ω}}^{\text{t}}}\mathrm{.}\end{array}

Therefore, we write the distribution function as follows:

\begin{array}{}\text{(66)}& f=f_0(\tilde{\epsilon })+f_1+f_2+f_3\mathrm{,}\end{array}

where the subscript 1, 2, and 3 stand for different orders in τ, and

\begin{array}{}\text{(67)}& f_1=\frac{\tau }{\hslash }\frac{eE+\frac{e^2}{\hslash }(E\cdot B){\text{Ω}}^{\text{t}}}{1+\frac{e}{\hslash }B\cdot {\text{Ω}}^{\text{t}}}\cdot \partial f_0(\tilde{\epsilon })\mathrm{,}\end{array}

\begin{array}{}\text{(68)}& \begin{array}{ll}{f}_2=\hfill & \frac{e^2{\tau }^2[E_i+{(v\times B)}_i]E_j}{{\hslash }^2}{\partial }_i{\partial }_j{f}_0(\epsilon -B\cdot m)\hfill \\ \hfill & +\frac{2e^2{\tau }^2[E_j+{(v\times B)}_j]{\partial }_i{\partial }_j{f}_0(\epsilon )}{{\hslash }^2}\left[\frac{E_i+{(v\times B)}_i}{1+\frac{e}{\hslash }B\cdot \text{Ω}}+\frac{e}{\hslash }(E\cdot B){\text{Ω}}_i\right]\hfill \\ \hfill & +\frac{e^2{\tau }^2[E_j+{(v\times B)}_j]{\partial }_i{f}_0(\epsilon )}{{\hslash }^2}{\partial }_j\left[\frac{E_i+{(v\times B)}_i}{1+\frac{e}{\hslash }B\cdot \text{Ω}}+\frac{e}{\hslash }(E\cdot B){\text{Ω}}_i\right]\mathrm{,}\hfill \end{array}\end{array}

\begin{array}{}\text{(69)}& f_3=\frac{e^3{\tau }^3{E}_i}{{\hslash }^3}{[(E+v\times B)\cdot \partial ]}^2{\partial }_i{f}_0(\epsilon )\mathrm{.}\end{array}

There is another contribution to the nonlinear conductivity. The system can reach a thermodynamic equilibrium together with a uniform magnetic field. During such process, if the sample remains charge neutral, the carrier density is fixed. As the magnetic field modifies the band energy, the chemical potential has to vary accordingly.

Without loss of generality, we assume B is along the z direction, i.e. B=B_z. At first order in B, the change in chemical potential can be fixed using the following expansion of the electron density n

\begin{array}{}\text{(70)}& \text{d}n={\frac{\partial n}{\partial B}|}_{\mu \mathrm{,}T\mathrm{,}B=0}B+{\frac{\partial n}{\partial \mu }|}_{T\mathrm{,}B=0}{\mu }_{(1)}=0.\end{array}

Note that {\frac{\partial n}{\partial \mu }|}_{T\mathrm{,}B=0}=g(\mu ), which is simply the density of states at the chemical potential for the unperturbed band. Moreover, {\frac{\partial n}{\partial B}|}_{\mu \mathrm{,}T} can be connected to the magnetization through the Maxwell equation:

\begin{array}{}\text{(71)}& {\frac{\partial n}{\partial B}|}_{\mu \mathrm{,}T}={\frac{\partial M_z}{\partial \mu }|}_{T\mathrm{,}B}\mathrm{.}\end{array}

We comment that for insulators at zero temperature, -\frac{\partial M_z}{\partial \mu }\to {\sigma }_{xy}, and Eq. (71) reduces to the familiar Streda formula. Combining Eqs. (70) and (71), we obtain

\begin{array}{}\text{(72)}& {\mu }_{(1)}=-\frac{{(\partial M_z/\partial \mu )}_{T\mathrm{,}B=0}}{g(\mu )}\mathrm{.}\end{array}

At second order in B, the expansion of n reads

\begin{array}{}\text{(73)}& \begin{array}{ll}\text{d}n=\hfill & {\frac{\partial n}{\partial \mu }|}_{T\mathrm{,}B=0}{\mu }_{(2)}+\frac{1}{2}{\frac{{\partial }^{2}n}{\partial B^2}|}_{T\mathrm{,}\mu \mathrm{,}B=0}{B}^2\hfill \\ \hfill & +\text{}\frac{1}{2}{\frac{{\partial }^{2}n}{\partial {\mu }^2}|}_{T\mathrm{,}B=0}{[{\mu }_{(1)}]}^2\text{}+\text{}{\frac{{\partial }^{2}n}{\partial B\partial \mu }|}_{T\mathrm{,}B=0}B{\mu }_{(1)}\mathrm{.}\hfill \end{array}\end{array}

{\mu }_{(2)} can be solved by setting the above equation equal to zero. The result can be simplified using Eq. (71) and the definition of the isothermal magnetic susceptibility {\chi }_{zz}=\partial M_z/\partial B, yielding

\begin{array}{}\text{(74)}& {\mu }_{(2)}=-\frac{1}{2}\frac{B^2\frac{\partial {\chi }_{zz}}{\partial \mu }+\frac{\partial g}{\partial \mu }{[{\mu }_{(1)}]}^2+2\frac{{\partial }^2{M}_z}{\partial {\mu }^2}B{\mu }_{(1)}}{g(\mu )}\mathrm{.}\end{array}

These change in chemical potential will also contribute to the nonlinear conductivity:

\begin{array}{}\text{(75)}& \begin{array}{ll}{J}_{\mu }=\hfill & -e{\mu }_{(1)}\frac{\partial }{\partial \mu }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\mathcal{D}\dot{r}f\hfill \\ \hfill & +\frac{e^2\tau }{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}v(E\cdot v)\left\{{\mu }_{(2)}{f}_0^{″}-\text{}\frac{1}{2}{[{\mu }_{(1)}]}^2{f}_0^{‴}\right\}\mathrm{,}\hfill \end{array}\end{array}

where f_0^{″} and f_0^{‴} are second and third order energy derivatives of f₀.

Equation (61) with Eqs. (62) and (67)–(69), combined with Eqs. (72), (74), and (75) yields a complete description of the nonlinear conductivity up to second order in the framework of the semiclassical theory. This procedure gives a complete account for the drifting part in the Boltzmann equation. However, it ignores the dependence of the collision integral on electromagnetic fields. Nevertheless, the semiclassical theory can yield a qualitatively valid result and offers a fresh perspective of how the band properties beyond the spectrum affect the conductivity.

Finally, we comment that to get the corresponding change in the resistivity, one can invert the conductivity tensor. If change in conductivity \delta {\sigma }_{xx} is obtained, the corresponding change in resistivity is

\begin{array}{}\text{(76)}& \frac{\delta {\rho }_{xx}}{{\rho }_{xx}}=-\frac{\delta {\sigma }_{xx}}{{\sigma }_{xx}}\mathrm{,}\end{array}

provided that the Hall conductivity {\sigma }_{xy} is much smaller than the longitudinal conductivity {\sigma }_{xx}.

3.2 Nonlinear anomalous Hall conductivity

The anomalous Hall effect refers to a Hall-type current in ferromagnets solely driven by an electric field. It is a topic under extensive studies [2]. In the past, there are three mechanisms identified: the intrinsic contribution, the skew-scattering contribution, and the side-jump contribution. Here we focus on the intrinsic contribution, which is due to the anomalous velocity in Eq. (16) [3]. From the perspective of the semiclassical theory, the Berry curvature acts as the magnetic field in the momentum space and bends the electron trajectory likewise, leading to a Hall-type current. The anomalous Hall effect has also been studied in noncolinear antiferromagnets, in which the net spin magnetization vanishes but the orbital magnetization is still present [76-85].

From symmetry consideration, the anomalous Hall effect requires broken time-reversal symmetry. In other words, it requires either a net spin magnetization as in ferromagnets or orbital magnetization as in noncolinear antiferromagnets. It is forbidden if time reversal symmetry is present. However, the nonlinear Hall effect may still be present, due to the manipulation of the magnetization through electromagnetic fields.

In systems with both time reversal and inversion symmetry, a net magnetization can be induced through the magnetic susceptibility by a magnetic field. This magnetization can then lead to a Hall-type current. The analytical expression for this correction can be obtained by plugging the anomalous velocity (second term in Eq. (62)) into the current in Eq. (61). The result reads [40]

\begin{array}{}\text{(77)}& J=\frac{e^2}{\hslash }E\times {\displaystyle \int }[\hslash v\times {\mathcal{A}}^{\prime }(B)+\text{Ω}(B\cdot m)]f_0^{\prime }\frac{\text{d}k}{8{\pi }^3}\mathrm{,}\end{array}

where {\mathcal{A}}^{\prime }(B) stands for the part of {\mathcal{A}}^{\prime } solely dependent on B and m is the spin plus orbital magnetic moment. It can be explicitly checked that compared with the susceptibility in Eq. (43), the second term in Eq. (77) is part of the geometric contribution to the magnetic susceptibility. Moreover, the first term is also part of the Van-Vleck and energy polarization contributions to the susceptibility. This current yields the conductivity {\sigma }_{(ij\mathrm{,}k)\mathrm{,0}} in Eq. (3).

This current has exactly the same dependence on electromagnetic fields with the ordinary Hall current. However, they have different origins. The ordinary Hall current is due to the Lorentz force and has the form of {\omega }_c\tau with {\omega }_c being the cyclotron frequency. In contrast, the current in Eq. (77) requires a nontrivial structure in the momentum space and does not involve the relaxation time. In fact, the ratio of the resistivity {\rho }_{xy}^{\prime } from Eq. (77) and the ordinary Hall resistivity {\rho }_{xy}^{\text{ord}} can be put in the following form

\begin{array}{}\text{(78)}& \frac{{\rho }_{xy}^{\prime }}{{\rho }_{xy}^{\text{ord}}}={\left({\rho }_{xx}\frac{e^2}{4h}\right)}^{2}S(\mu )\mathrm{,}\end{array}

where the first factor is universal and S(\mu ) is a model-dependent factor but independent of the scattering process. The universal scaling factor means that {\rho }_{xy}^{\prime } will dominate the ordinary Hall effect when {\rho }_{xx} is large, i.e., in dirty metals/semiconductors with relatively small relaxation time. In a typical Hall-bar measurement set-up, both {\rho }_{xy}^{\prime } and {\rho }_{xy}^{\text{ord}} will contribute to the Hall current. Therefore, the total Hall resistivity should be dependent on the relaxation time. To differentiate one contribution from the other, one should change the universal scaling factor through temperature, film thickness, or doping, and measure the scaling behaviour based on Eq. (78).

In systems that simultaneously break the time reversal and inversion symmetry but preserve the combined symmetry, the anomalous Hall current also vanishes. This symmetry is exactly the magnetoelectric symmetry which allows a magnetization to be induced by an electric field through the magnetoelectric coefficient [85]. This magnetization can then lead to a Hall-type current. The analytical expression can be obtained in a similar way with the magnetic-field-correction to anomalous Hall effect. The result reads [40]

\begin{array}{}\text{(79)}& J=e^{2}E\times {\displaystyle \int }[v\times {\mathcal{A}}^{\prime }(E)]f_0^{\text{'}}\frac{\text{d}k}{8{\pi }^3}\mathrm{,}\end{array}

where {\mathcal{A}}^{\prime }(E) stands for the part of {\mathcal{A}}^{\prime } solely dependent on E. This current yields the conductivity {\sigma }_{ijk\mathrm{,0}} in Eq. (2). The connection between this Hall-type current and the magnetoelectric effect can be mostly seen for a two-band model, as shown in Fig. 6.

In noncentrosymmetric but time-reversal-invariant materials, a net magnetization can be induced in the nonequilibrium steady state through the Edelstein effect. Such magnetization can further lead to a Hall-type current. The analytical expression can be obtained by plugging Eq. (62) and (67) into Eq. (61) and keeping terms at second order in E, and reads [86]

\begin{array}{}\text{(80)}& J=\frac{e^3}{\hslash }\tau E\times {\displaystyle \int }\text{Ω}(E\cdot v)f_0^{\text{'}}\frac{\text{d}k}{8{\pi }^3}\mathrm{.}\end{array}

This current yields the conductivity {\sigma }_{ijk\mathrm{,1}} in Eq. (2). Interestingly, the requirement of the mirror symmetry in Fig. 6(b) also works for this current. But one should keep in mind that for the current in Eq. (79), as the time reversal symmetry is broken, the inherent spin texture in the sample is also subjected to the mirror operation. The observation of this nonlinear Hall current has been reported in several recent experiments [87-89]. We comment that in Eq. (80), v\text{Ω} together constitutes a pseudotensor which represents the first order moment of the Berry curvature in the momentum space. Hence it is referred to as the Berry curvature dipole. This Berry curvature dipole, together with the magnetic moment dipole also play essential roles in the study of the natural optical activity [6, 7, 22, 85, 90, 91].

Finally, using the quantum kinetic equations for the evolution of the density matrix under both electric field and disorders, it has been found that in noncentrosymmetric materials, besides the Berry curvature dipole contribution, there are also nonlinear side jump and skew scattering mechanisms responsible for the nonlinear anomalous Hall effect [92]. It can compete with the nonlinear anomalous Hall effect due to Berry curvature dipole.

3.3 Linear magnetoresistance

Linear magnetoresistance, referring to the linear dependence of the resistance on magnetic field, was first observed in simple metals [93-95] and under continuous studies afterwards. There is a surge of interest in recent years, as it appears in various novel materials, such as multilayer graphene [96], topological insulators [97-100], and Weyl/Dirac semimetals [101-106]. It has been shown that the linear magnetoresistance can arise in the ultra-quantum regime for Dirac systems when only the lowest Landau level is partially filled [107-109], as well as in inhomogeneous samples with mobility fluctuation [110-113]. In the semiclassical regime, linear magnetoresistance is subject to a symmetry requirement. Generally speaking Onsager relation requires that {\sigma }_{xx}(B)={\sigma }_{xx}(-B), which forbids odd power dependence of the conductivity on magnetic field. However, if the time reversal symmetry is broken and a magnetization is present, Onsager relation becomes {\sigma }_{xx}(B\mathrm{,}M)={\sigma }_{xx}(-B\mathrm{,}-M), which can be satisfied in principle for terms with odd powers of B.

To derive the linear magnetoresistance, or equivalently the linear magnetoconductivity in the semiclassical regime, one collects terms linear in B in the current in Eq. (61). To simplify the result, we first consider the transverse magnetoresistance and assume that the electric field is along \widehat{x} direction and B is along \widehat{z} direction. This configuration will yield the following conductivity (here we use the definition in Eq. (3)) [114]

\begin{array}{}\text{(81)}& \begin{array}{ll}{\sigma }_{(xx\mathrm{,}z)\mathrm{,1}}=\hfill & \frac{e^2}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\left[2v_x\frac{\partial m_z}{\partial k_x}+\frac{e}{\hslash }{(v_x)}^2{\text{Ω}}_z\right]f_0^{\prime }\hfill \\ \hfill & +\frac{e^2}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}{(v_x)}^2(m_z+\delta {\mu }_{(1)}/B_z)f_0^{″}\mathrm{,}\hfill \end{array}\end{array}

\begin{array}{}\text{(82)}& {\sigma }_{(xx\mathrm{,}z)\mathrm{,2}}=-e^3{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(v_x{v}_y{\alpha }_{xx}-v_x^2{\alpha }_{xy})f_0^{\prime }=0.\end{array}

For longitudinal magnetoresistance, we assume that the electric field is parallel to the magnetic field, i.e., it is along \widehat{z} direction. Then the magnetoconductivity reads

\begin{array}{}\text{(83)}& \begin{array}{l}\text{}\\ \text{}\text{}\\ \begin{array}{ll}{\sigma }_{(zz\mathrm{,}z)\mathrm{,1}}\text{}=\hfill & \frac{e^2}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\left[2v_z\frac{\partial m_z}{\partial k_z}\text{}+\text{}\frac{e}{\hslash }{v}_z(v_z{\text{Ω}}_z\text{}-\text{}2v\cdot \text{Ω})\right]f_0^{\prime }\hfill \\ \hfill & +\frac{e^2}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}{(v_z)}^2(m_z+\delta {\mu }_{(1)}/B_z)f_0^{″}\mathrm{,}\hfill \end{array}\end{array}\end{array}

\begin{array}{}\text{(84)}& {\sigma }_{(zz\mathrm{,}z)\mathrm{,2}}\text{}=\text{}-e^3{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(v_z{v}_y{\alpha }_{xz}-v_z{v}_x{\alpha }_{yz})f_0^{\prime }=0.\end{array}

Both {\sigma }_{(xx\mathrm{,}z)\mathrm{,1}} and {\sigma }_{(zz\mathrm{,}z)\mathrm{,1}} are clearly violations to the Kohler’s rule [115]. In other words, the magnetic field affects the conductivity not through the Lorentz force, but through the two extra corrections in the first order semiclassical equations of motion. On one hand, it couples to the k-dependent magnetic moment to modifies the band energy and hence band velocity. On the other hand, it couples to the Berry curvature and changes the phase of space density of states, or equivalently, the carrier density. In the longitudinal magnetoconductivity, there is an additional term with the factor v\cdot \text{Ω}. This actually measures the flux strength of the momentum space magnetic field \text{Ω}. This unique contribution to the longitudinal magnetoresistance is due to the coupling between the real and momentum space lorentz force, i.e., the magnetic field along \widehat{z} can first bend the velocity along \widehat{x} (or \widehat{y}) to the direction \widehat{y} (or \widehat{x}), which is further bent to the \widehat{z} direction due to the momentum space magnetic field {\text{Ω}}_x (or {\text{Ω}}_y). Finally, we comment that both {\sigma }_{(xx\mathrm{,}z)\mathrm{,1}} and {\sigma }_{(zz\mathrm{,}z)\mathrm{,1}} are linear in magnetic moment and Berry curvature, which are odd in time reversal operation. As a result, these two conductivities only appear in materials that break time reversal symmetry.

In the above, we show that {\sigma }_{(xx\mathrm{,}z)\mathrm{,2}} and {\sigma }_{(zz\mathrm{,}z)\mathrm{,2}} vanish identically in general. We comment that if a detailed consideration of scattering process is performed, there is a contribution to {\sigma }_{(xx\mathrm{,}z)} that is of second order in τ, due to the skew scattering process [114]. Therefore, it requires a nontrivial antisymmetric part of the scattering probability W_{k{k}^{\prime }} and hence the breaking of time reversal symmetry.

When the chemical potential approaches the hot spot of the Berry curvature and the magnetic moment, the latter can have a large contribution to the linear magnetoresistance, as shown in Fig. 7. As the chemical potential approaches the band edge, the Berry curvature and magnetic moment increases and so does the linear magnetoresistance. However, the semiclassical theory is only valid when (\mu -\text{Δ})\tau /\hslash >1. Otherwise, the scattering from different impurities will be correlated and a fully quantum mechanical treatment should be used.

3.4 Negative longitudinal magnetoresistance induced by chiral anomaly

It is well known that the chemical potential difference in real space can lead to a current. Its counterpart in momentum space, i.e., the current due to the chemical potential difference in momentum space, has completely different origins and attracts great attention recently as it manifests in topological semimetals.

To understand this effect, we note that from the semiclassical theory, the Berry curvature is the analog of the magnetic field in the momentum space. The flux of the Berry curvature, defined as (1/2\pi ){\displaystyle \oint }\text{d}S\cdot \text{Ω}, can be nonzero through a closed surface, indicating the appearance of a monopole charge in the momentum space. This is most easily realized by a pair of Weyl points

\begin{array}{}\text{(85)}& \widehat{H}={\lambda }_{i}vk\cdot \sigma \mathrm{.}\end{array}

where {\lambda }_i=\pm 1 stands for a positive or negative charge, or equivalently, the chirality of the Weyl node. For such pair of Weyl nodes, there is a current driven by the magnetic field. This is the chiral magnetic current, derived using the last term in Eq. (62)

\begin{array}{}\text{(86)}& J_{\text{CME}}=-\frac{e^2}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(v\cdot \text{Ω})B{f}_0=-\frac{e^{2}B}{4{\pi }^2{\hslash }^2}{\displaystyle \sum _i}{\lambda }_i{\mu }_i\mathrm{,}\end{array}

where {\mu }_i is the chemical potential for each Weyl node. This current has to vanish in equilibrium [73], even when the two Weyl nodes are shifted in energy and the relative chemical potential is different.

However, away from equilibrium, the chemical potential for different valleys can respond to external fields in different ways, leading to an inhomogeneous chemical potential distribution in the momentum space and hence a net current through J_{\text{CME}} [7]. To see this, we first write down the Boltzmann equation:

\begin{array}{}\text{(87)}& \dot{k}\cdot \partial f_0=-{\displaystyle \int }\text{d}{k}^{\prime }\mathcal{D}(k^{\prime })W_{k{k}^{\prime }}[f(k)-f(k^{\prime })]\mathrm{,}\end{array}

where W_{k{k}^{\prime }} is the scattering probability. Using the semiclassical equations of motion in Eqs. (16) and (17), one finds that

\begin{array}{}\text{(88)}& \mathcal{D}\dot{k}=-\frac{e}{\hslash }E-\frac{e}{\hslash }\tilde{v}\times B-\frac{e^2}{{\hslash }^2}(E\cdot B)\text{Ω}\mathrm{.}\end{array}

Therefore, the Boltzmann equation becomes

\begin{array}{}\text{(89)}& \begin{array}{l}\frac{e}{\mathcal{D}(k)\hslash }\left[E+\frac{e}{\hslash }(E\cdot B)\text{Ω}\right]\cdot \partial f_0\\ ={\displaystyle \int }\text{d}{k}^{\prime }\mathcal{D}(k^{\prime })W_{k{k}^{\prime }}[f(k)-f(k^{\prime })]\mathrm{.}\end{array}\end{array}

In multivalley systems such as a pair of Weyl nodes, a special solution to the above equation may emerge. To see this, we note that usually the transport current is determined by the fastest scattering rate which is the intravalley scattering. For the intravalley scattering with f_0(k) and f_0(k^{\prime }) from the same type of Weyl node, one can multiply Eq. (89) by \mathcal{D}(k) and integrate it around the Weyl node. The result is

\begin{array}{}\text{(90)}& \frac{e}{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}E\cdot \left[v{f}_0^{\prime }-\partial (B\cdot m{f}_0^{\prime })+\frac{e}{\hslash }B(\text{Ω}\cdot v)f_0^{\prime }\right]=0.\end{array}

However, for a Weyl node, although the integration of the first two terms vanish, the integration of the third term is proportional to the charge of the Weyl node and hence does not vanish. It indicates that charge can be pumped from the positive Weyl node to the negative Weyl node. This contradiction shows that the intervalley scattering is indispensable in solving the above Boltzmann equation.

To explicitly derive this solution, one can take f_0(k) to be from the positive Weyl node and f_0(k^{\prime }) from the negative Weyl node. Then under the relaxation time approximation, Eq. (90) for the positive Weyl node becomes

\begin{array}{}\text{(91)}& \frac{e^2}{{\hslash }^2}{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(E\cdot B)(\text{Ω}\cdot v)f_0^{\prime }=\frac{\text{Δ}{N}_{+}-\text{Δ}{N}_{-}}{{\tau }_{\text{int}}}\mathrm{,}\end{array}

where \text{Δ}{N}_{+} and \text{Δ}{N}_{-} is the density of particles near the Fermi surface, and {\tau }_{\text{int}} is the intervalley scattering time. Eq. (91) means that the pumping of particles from one Weyl node to another through drifting has to be balanced by the diffusive collision through the intervalley scattering process. It leads to an inhomogeneous chemical potential distribution in the momentum space

\begin{array}{}\text{(92)}& \delta {\mu }_{+}-\delta {\mu }_{-}=-\frac{e^2}{4{\pi }^2{\hslash }^2}\frac{{\tau }_{\text{int}}}{g(\mu )}(E\cdot B)\mathrm{,}\end{array}

where g(\mu )={\mu }^2/(2{\pi }^2{\hslash }^3{v}^3) is the density of states at the Fermi surface and {\tau }_{\text{int}} is the intervalley scattering time. This chemical potential difference for Weyl nodes with opposite chiralities is the chiral anomaly effect [8, 9, 116-118].

This imbalance in the chemical potential can lead to a current through Eq. (86)

\begin{array}{}\text{(93)}& J_{\text{CME}}=\frac{e^4{v}^3{\tau }_{\text{int}}}{2{\pi }^3\hslash {\mu }^2}(E\cdot B)B\mathrm{.}\end{array}

This current is only present when E is parallel to B. The corresponding conductivity is quadratic in magnetic field and increases with it, indicating a negative longitudinal magnetoresistance. It corresponds to {\sigma }_{(ij\mathrm{,}kk)\mathrm{,1}} in Table. 1 and does not need to break time reversal or inversion symmetry.

Besides the semiclassical theory, the negative magnetoresistance has also been confirmed using chiral kinetic theory [119-123] and exact quantized Landau levels [124]. In the semiclassical regime, the negative magnetoresistance for one pair of Weyl nodes is also extended in other situations. It has been shown that samples with Dirac points also have negative magnetoresistance, provided that the chiral charge is approximately conserved and the chiral relaxation time is large [125]. Moreover, when the Dirac point opens a gap due to the suppression of the lattice symmetry, the chiral relaxation time is still present, so is the negative magnetoresistance [126]. Experimentally, the negative longitudinal magnetoresistance in Weyl and Dirac semimetals has been reported and interpreted as the chiral magnetic effect [102, 106, 127-136], but similar signals could also arise from two different mechanisms discussed in the following two sections.

Finally, we comment that the Weyl nodes is gapless and near those band-crossing points, the semiclassical theory may fail. However, we can set the Fermi surface to be well above the Weyl point. The flux through the Fermi surface does not change with the Fermi energy. As a result, the discussion using the semiclassical theory still works.

3.5 Intrinsic quadratic longitudinal magnetoresistance

The magnetoconductivity due to chiral anomaly and chiral magnetic effect is one contribution to {\sigma }_{(ij\mathrm{,}kk)} in Eq. (4). It is not a complete description of {\sigma }_{(ij\mathrm{,}kk)}. To fully account for the drifting contribution to {\sigma }_{(ij\mathrm{,}kk)}, one has to use the semiclassical theory up to second order, as presented in Section 3.1. Specially, one can plug Eqs. (67) and (62) into Eq. (61) and collect terms up to second order in magnetic field. The current responsible for {\sigma }_{(ij\mathrm{,}kk)\mathrm{,1}} is complicated but it can be put in several groups:

\begin{array}{}\text{(94)}& J=J_1+J_2+J_3+J_4+J_5\mathrm{.}\end{array}

The first contribution has the form

\begin{array}{}\text{(95)}& J_1=\frac{e^4}{{\hslash }^2}\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}[(v\text{}\times \text{}B)\cdot (E\times \text{Ω})](v\text{}\times \text{}B)\text{}\times \text{}\text{Ω}{f}_0^{\prime }\mathrm{.}\end{array}

This term is purely due to the anomalous velocity as both E×Ω and (v×B)×Ω are parts of the anomalous velocity. This current is particularly interesting. If the electromagnetic fields are along the same direction, which can be assume to be the \widehat{z} direction without loss of generality, this current becomes

\begin{array}{}\text{(96)}& J_{\mathrm{1,}z}=-\frac{e^4}{{\hslash }^2}E{B}^2\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}{(v_x{\text{Ω}}_x+v_y{\text{Ω}}_y)}^2{f}_0^{\prime }\mathrm{.}\end{array}

This will always lead to a negative longitudinal magnetoresistance.

The second contribution has the form

\begin{array}{}\text{(97)}& \begin{array}{ll}{J}_2=\hfill & \frac{e^2\tau }{2}{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}{\widehat{e}}_i{\alpha }_{ij}{E}_j{(B\cdot m)}^2{f}_0^{″}\hfill \\ \hfill & -\frac{e^2\tau }{{\hslash }^2}{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\frac{\partial (B\cdot m)}{\partial k}(E\cdot {\partial }_k)[(B\cdot m)f_0^{\prime }]\mathrm{.}\hfill \end{array}\end{array}

This is solely due to the Zeeman coupling between magnetic field and the magnetic moment, containing both the spin and orbital part.

The third contribution has the form

\begin{array}{}\text{(98)}& \begin{array}{ll}{J}_3=\hfill & \frac{e^3}{{\hslash }^2}\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(E\cdot B)v(\text{Ω}\cdot {\partial }_k)[(B\cdot m)f_0^{\prime }]\hfill \\ \hfill & +\frac{e^3}{{\hslash }^2}\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(E\cdot B)(\text{Ω}\cdot v)\frac{\partial (B\cdot m)}{\partial k}{f}_0^{\prime }\hfill \\ \hfill & +\frac{e^3}{{\hslash }^2}\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}[(v\times B)\times \text{Ω}](E\cdot {\partial }_k)[(B\cdot m)f_0^{\prime }]\hfill \\ \hfill & +\frac{e^3}{{\hslash }^2}\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}\left\{\left[\frac{\partial (B\cdot m)}{\partial k}\times B\right]\times \text{Ω}\right\}[(E\cdot v)f_0^{\prime }]\mathrm{.}\hfill \end{array}\end{array}

This current is due to the coupling between the anomalous velocity and the Zeeman energy correction.

The fourth contribution has the form

\begin{array}{}\text{(99)}& \begin{array}{ll}{J}_4=\hfill & e^2\tau {\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(2{\widehat{e}}_i{\alpha }_{ij}{E}_j{f}_0^{\prime }+v(E\cdot v)f_0^{″}){\epsilon }_{(2)}\hfill \\ \hfill & +\frac{e^3\tau }{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}[(v\times B)\cdot (E\times {\text{Ω}}^{\prime })]v_0{f}_0^{\prime }\hfill \\ \hfill & -\frac{e^3\tau }{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}(v\cdot {\text{Ω}}^{\text{'}})B(E\cdot v)f_0^{\prime }\mathrm{.}\hfill \end{array}\end{array}

The first line is due to the second order correction to the band energy and the remaining terms are due to the first order correction to the Berry curvature, i.e., {\text{Ω}}^{\prime }.

From the expressions of J₁ to J₄, we find that J₁ is purely due to the orbital motion of Bloch electrons, while J₂ to J₄ can contain both spin and orbital contribution. The orbital part has been calculated in Ref. [137] and the spin part has been calculated in Ref. [138]. Here J₁ to J₄ are complete expressions that contain the spin part, the orbital part, and the spin-orbital part contribution.

The last contribution is due to the correction to the chemical potential:

\begin{array}{}\text{(100)}& J_5=\frac{e^2\tau }{\hslash }{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}v(E\cdot v)\left\{{\mu }_{(2)}{f}_0^{″}-\frac{1}{2}{[{\mu }_{(1)}]}^2{f}_0^{‴}\right\}\mathrm{.}\end{array}

In nonmagnetic metals, {\sigma }_{(ij\mathrm{,}kk)\mathrm{,3}} is also nonzero. The corresponding current reads

\begin{array}{}\text{(101)}& J^{\prime }=-\frac{e^4}{{\hslash }^2}{\tau }^3{\displaystyle \int }\frac{\text{d}k}{8{\pi }^3}v[(v\times B\cdot {\partial }_k{)}^2(E\cdot v)]f_0^{\prime }\mathrm{.}\end{array}

This current requires a peculiar Fermi surface such that the integration does not vanish [75].

The current J and J’ complete the description of the quadratic magnetoresistance. In J’ the magnetic field always couples to the relaxation time. Therefore, the corresponding magnetoresistance follows Kohler’s rule [115]. In comparison, J violates Kohler’s rule and its corresponding magnetoresistance does not depend on the relaxation time. Moreover, J does not need a highly anisotropic Fermi surface. As long as the Berry curvature, magnetic moment, and magnetic susceptibility does not vanish at the same time, J will be nonzero. This intrinsic contribution to magnetoresistance is also confirmed by a fully quantum mechanical treatment [139].

The current J can yield a negative longitudinal magnetoresistance. In fact, for the Dirac Hamiltonian in Eq. (85), the ratio of the quadratic longitudinal magnetoresistance \delta \rho from J to that from chiral anomaly reads

\begin{array}{}\text{(102)}& \frac{\delta \rho }{\delta {\rho }_{\text{CA}}}=+\frac{8}{45}\frac{\tau }{{\tau }_{\text{int}}}\mathrm{.}\end{array}

The global positive sign means that \delta \rho is negative, like \delta {\rho }_{\text{CA}}. Usually \tau \ll {\tau }_{\text{int}}. Therefore, in Weyl semimetals, the negative magnetoresistance should be dominated by the contribution from the chiral anomaly in the semiclassical regime. However, the negative magnetoresistance from J can persist in metals without any Dirac/Weyl point, as shown in Fig. 8. Panel (a) suggests that the longitudinal magnetoresistance is negative while the transverse one is positive. This causes the sign change of the effective magnetoresistance as the angle between E and B continuously changes from 0 to \pi /2, as shown in Panel (b).

It is interesting to further check the behaviour of the magnetoconductivity in a lattice model. We use the model Hamiltonian in Eqs. (59) and (60) and choose the parameters such that it is in the semiconductor phase as shown in Fig. 4(d). We then calculate the quadratic magnetoconductivity for the transverse and longitudinal configuration, using the expression J₁ to J₅. Both spin and orbital part of the magnetic moment is considered. The result is presented in Fig. 9. It can be found that both magnetoconductivity show strong dependence on the chemical potential. Consistent with the result in the low-energy model, the magnetoresistance will be strongly enhanced when chemical potential is near the band edge, as both the anomalous velocity and energy correction are enhanced. Moreover, the longitudinal magnetoconductivity is generally positive, indicating a negative longitudinal magnetoresistance, while the transverse magnetoconductivity is negative. Near the band edge, the longitudinal magnetoconductivity is larger than the transverse one, indicating a sizeable transition angle θ from negative to positive magnetoresistance, as the relative angle between E and B changes from 0 to \pi /2. This is similar to the negative magnetoresistance induced by the chiral anomaly.

3.6 Current jetting

Beyond the Boltzmann transport theory, the magnetoresistance is also greatly affected when the current distribution in the sample is inhomogeneous. For example, for a matchstick sample with two leads for current input and output, as shown in Fig. 10, the current in the sample can be highly inhomogeneous provided that the transverse magnetoresistance is sufficiently large [115, 140]. This is the current jetting problem.

The current jetting problem can be solved through the charge continuity equation: \partial \rho /\partial t+\nabla \cdot J=0. In the steady state \partial \rho /\partial t=0. Therefore, in the sample, the following equation for the electrostatic potential is satisfied

\begin{array}{}\text{(103)}& {\sigma }_{ij}{\partial }_i{\partial }_{j}V=0.\end{array}

Since the Hall conductivity is antisymmetric, it does not affect the above equation. However, it can affect the boundary condition. If the Hall conductivity can be ignored, and the longitudinal conductivity is determined by the Drude formula: {\sigma }_{xx}={\sigma }_{yy}={\sigma }_0/(1+{\omega }_c^2{\tau }^2) and {\sigma }_{zz}={\sigma }_0, the above equation can be put in the form of the Poisson equation

\begin{array}{}\text{(104)}& \frac{{\partial }^{2}V}{\partial x^2}+\frac{{\partial }^{2}V}{\partial y^2}+R^2\frac{{\partial }^{2}V}{\partial z^2}=\mathrm{0,}\end{array}

where R=\sqrt{{\sigma }_{zz}/{\sigma }_{xx}}. The boundary condition should also be determined by the charge continuity equation.

A large R will greatly stretch the equipotential contour along \widehat{z} direction. It can cause a great portion of current lies in a sharp cone centered around the \widehat{z} axis [115]. This is similar to the behaviour of the longitudinal magnetoresistance due to chiral anomaly. In fact, if the injection of the current only localized near a small region in the surface and does not spread to the whole surface, the current injection may cause a strong negative magnetoresistance [106, 141, 142].

4 Outlook

We have shown how to derive various nonlinear currents in the framework of the semiclassical theory and Boltzmann transport theory. Band properties beyond the spectrum play essential roles in those nonlinear currents. However, the semiclassical theory ignores the electron interaction and the modification to electron-impurity scattering from electromagnetic fields. The semiclassical theory has been extended in Fermi liquid using Keldysh formalism up to first order [143, 144]. Moreover, the Boltzmann transport theory is also extended to treat the drifting and collision in equal-footing up to first order [123, 145]. It is tempting to extend such theories up to second order to give a complete treatment of various nonlinear currents.

For the three nonlinear anomalous Hall effects, the one due to Berry curvature dipole has attracted several experimental efforts, while no experiments are reported to measure the other two corrections. The correction due to magnetic field should be generally present in any system and competes with the ordinary Hall effect. Distinguishing the magneto nonlinear Hall effect from the ordinary Hall effect is similar to extracting the intrinsic anomalous Hall from total anomalous Hall signal [146-151].

Moreover, in samples that have combined time reversal and inversion symmetry (but break each one separately), the linear anomalous Hall effect as well as the nonlinear Hall effect due to Berry curvature dipole vanishes identically and only the intrinsic nonlinear Hall current contributes. This type of materials is ideal to test the intrinsic nonlinear Hall effect. If the time reversal and inversion symmetry, as well as their combined symmetry are all broken, the anomalous Hall and its two electric nonlinear version should all be present. There are still methods to distinguish them. One can first perform a scaling of the relaxation time to distinguish the intrinsic contribution from the extrinsic contribution. This is similar to extracting intrinsic, skew-scattering and side-jump contribution from the total anomalous Hall signal [146-151]. The nonlinear part can be further extracted by reversing the direction of the electric field and finding the difference in the Hall signal.

For the negative longitudinal magnetoresistance, the difficulty is to distinguish the contribution from chiral anomaly from the intrinsic and current jetting contributions. The current jetting effect may be reduced by attaching the lead across the whole end of the sample so that the current is uniformly injected. In comparison, the main difference between the chiral anomaly and the intrinsic contribution is a hierarchy of relaxation time. One needs to control the ratio of the intervalley or chiral relaxation time to the transport relaxation time to distinguish these two contributions.

There is another issue in evaluating nonlinear currents in first-principles codes. Depending on the type of the unperturbed Hamiltonian, the second order energy has different forms for Schrodinger Hamiltonian or effective tight-binding Hamiltonian, mainly due to different Hessian matrices. Then there are two methods to evaluate the second order energies and related effects: one can start from the Schrödinger Hamiltonian and use the wave function from first-principles codes as approximation to evaluate the analytical expression derived for the Schrödinger Hamiltonian; on the other hand, one can start from the effective tight-binding Hamiltonian and uses its response to electromagnetic fields to approach the real one. The second one requires to evaluate the effective Hessian matrices while that of the first one is simply 1/(m{\delta }_{ij}). It is tempting to study which method is more accurate.

References

[1] N. W. Ashcroft and M. Mermin, Solid State Physics, Harcourt, Orlando, New York, 1976

[2] N. Nagaosa, J. Sinova, S. Onoda, A. H. MacDonald, and N. P. Ong, Anomalous Hall effect, Rev. Mod. Phys. 82(2), 1539 (2010)

[3] D. Xiao, M. C. Chang, and Q. Niu, Berry phase effects on electronic properties, Rev. Mod. Phys. 82(3), 1959 (2010)

[4] X. L. Qi and S. C. Zhang, Topological insulators and superconductors, Rev. Mod. Phys. 83(4), 1057 (2011)

[5] T. Morimoto and N. Nagaosa, Topological nature of nonlinear optical effects in solids, Sci. Adv. 2(5), e1501524 (2016)

[6] S. Zhong, J. E. Moore, and I. Souza, Gyrotropic magnetic effect and the magnetic moment on the Fermi surface, Phys. Rev. Lett. 116(7), 077201 (2016)

[7] J. Ma and D. A. Pesin, Chiral magnetic effect and natural optical activity in metals with or without Weyl points, Phys. Rev. B 92(23), 235205 (2015)

[8] D. T. Son and B. Z. Spivak, Chiral anomaly and classical negative magnetoresistance of Weyl metals, Phys. Rev. B 88(10), 104412 (2013)

[9] N. Armitage, E. Mele, and A. Vishwanath, Macroscopic polarization in crystalline dielectrics: The geometric phase approach, Rev. Mod. Phys. 66, 899 (2018)

[10] The highest order of τ is the same as the highest order of the field, which is the property of the asymptotic solution to the Boltzmann equation, as implied in Eq. (64).

[11] J. Zak, Magnetic translation group, Phys. Rev. 134(6A), A1602 (1964)

[12] J. Zak, Magnetic translation group (II): Irreducible representations, Phys. Rev. 134(6A), A1607 (1964)

[13] D. R. Hofstadter, Energy levels and wave functions of Bloch electrons in rational and irrational magnetic fields, Phys. Rev. B 14(6), 2239 (1976)

[14] M. C. Chang and Q. Niu, Berry phase, hyperorbits, and the Hofstadter spectrum, Phys. Rev. Lett. 75(7), 1348 (1995)

[15] G. Sundaram and Q. Niu, Wave-packet dynamics in slowly perturbed crystals: Gradient corrections and Berry-phase effects, Phys. Rev. B 59(23), 14915 (1999)

[16] D. Culcer, Y. Yao, and Q. Niu, Coherent wave-packet evolution in coupled bands, Phys. Rev. B 72(8), 085110 (2005)

[17] R. Shindou and K. I. Imura, Noncommutative geometry and non-Abelian Berry phase in the wave-packet dynamics of Bloch electrons, Nucl. Phys. B 720(3), 399 (2005)

[18] M.-C. Chang and Q. Niu, Berry curvature, orbital moment, and effective quantum theory of electrons in electromagnetic fields, J. Phys.: Condens. Matter 20, 193202 (2008)

[19] M. V. Berry, Quantal phase factors accompanying adiabatic changes, Proc. R. Soc. Lond. A Math. Phys. Sci. 392, 45 (1984)

[20] D. Xiao, J. Shi, and Q. Niu, Berry phase correction to electron density of states in solids, Phys. Rev. Lett. 95(13), 137204 (2005)

[21] J. E. Moore and J. Orenstein, Confinement-Induced Berry phase and helicity-dependent photocurrents, Phys. Rev. Lett. 105(2), 026805 (2010)

[22] S. Zhong, J. Orenstein, and J. E. Moore, Optical gyrotropy from axion electrodynamics in momentum space, Phys. Rev. Lett. 115(11), 117403 (2015)

[23] T. Morimoto, S. Zhong, J. Orenstein, and J. E. Moore, Semiclassical theory of nonlinear magneto-optical responses with applications to topological Dirac/Weyl semimetals, Phys. Rev. B 94(24), 245121 (2016)

[24] E. Deyo, L. E. Golub, E. L. Ivchenko, and B. Spivak, Semiclassical theory of the photogalvanic effect in non-centrosymmetric systems, arXiv: 0904.1917 (2019)

[25] R. Resta, Manifestations of Berry’s phase in molecules and condensed matter, J. Phys.: Condens. Matter 12(9), R107 (2000)

[26] R. D. King-Smith and D. Vanderbilt, Theory of polarization of crystalline solids, Phys. Rev. B 47, 1651(R) (1993)

[27] R. Resta, Macroscopic polarization in crystalline dielectrics: The geometric phase approach, Rev. Mod. Phys. 66(3), 899 (1994)

[28] M. C. Chang, and Q. Niu, Berry phase, hyperorbits, and the Hofstadter spectrum: Semiclassical dynamics in magnetic Bloch bands, Phys. Rev. B 53(11), 7010 (1996)

[29] I. Souza and D. Vanderbilt, Dichroic f-sum rule and the orbital magnetization of crystals, Phys. Rev. B 77(5), 054438 (2008)

[30] W. Yao, D. Xiao, and Q. Niu, Valley-dependent optoelectronics from inversion symmetry breaking, Phys. Rev. B 77(23), 235406 (2008)

[31] L. L. Foldy and S. A. Wouthuysen, On the Dirac theory of spin 1/2 particles and its non-relativistic limit, Phys. Rev. 78(1), 29 (1950)

[32] E. Blount, Extension of the Foldy–Wouthuysen transformation, Phys. Rev. 128(5), 2454 (1962)

[33] D. Xiao, Y. Yao, Z. Fang, and Q. Niu, Berry-phase effect in anomalous thermoelectric transport, Phys. Rev. Lett. 97(2), 026603 (2006)

[34] T. Thonhauser, D. Ceresoli, D. Vanderbilt, and R. Resta, Orbital magnetization in periodic insulators, Phys. Rev. Lett. 95(13), 137205 (2005)

[35] D. Ceresoli, T. Thonhauser, D. Vanderbilt, and R. Resta, Orbital magnetization in crystalline solids: Multi-band insulators, Chern insulators, and metals, Phys. Rev. B 74(2), 024408 (2006)

[36] O. Gat and J. E. Avron, Magnetic fingerprints of fractal spectra and the duality of Hofstadter models, New J. Phys. 5, 44 (2003)

[37] O. Gat, and J. E. Avron, Semiclassical analysis and the magnetization of the Hofstadter model, Phys. Rev. Lett. 91(18), 186801 (2003)

[38] J. Shi, G. Vignale, D. Xiao, and Q. Niu, Quantum theory of orbital magnetization and its generalization to interacting systems, Phys. Rev. Lett. 99(19), 197202 (2007)

[39] T. Qin, Q. Niu, and J. Shi, Energy magnetization and the thermal Hall effect, Phys. Rev. Lett. 107(23), 236601 (2011)

[40] Y. Gao, S. A. Yang, and Q. Niu, Field Induced Positional Shift of Bloch electrons and its dynamical implications, Phys. Rev. Lett. 112(16), 166601 (2014)

[41] Y. Gao, S. A. Yang, and Q. Niu, Geometrical effects in orbital magnetic susceptibility, Phys. Rev. B 91(21), 214405 (2015)

[42] E. Blount, Bloch electrons in a magnetic field, Phys. Rev. 126(5), 1636 (1962)

[43] J. P. Provost, and G. Vallee, Riemannian structure on manifolds of quantum states, Commun. Math. Phys. 76(3), 289 (1980)

[44] T. Neupert, C. Chamon, and C. Mudry, Measuring the quantum geometry of Bloch bands with current noise, Phys. Rev. B 87(24), 245103 (2013)

[45] J. Anandan and Y. Aharonov, Geometry of quantum evolution, Phys. Rev. Lett. 65(14), 1697 (1990)

[46] R. Resta, The insulating state of matter: A geometrical theory, Eur. Phys. J. B 79(2), 121 (2011)

[47] M. V. Berry, Quantum phase corrections from adiabatic iteration, Proc. R. Soc. Lond. A 414(1846), 31 (1987)

[48] A. M. Essin, A. M. Turner, J. E. Moore, and D. Vanderbilt, Orbital magnetoelectric coupling in band insulators, Phys. Rev. B 81(20), 205104 (2010)

[49] L. Onsager, Interpretation of the de Haas–van Alphen effect, Philos. Mag. 43(344), 1006 (1952)

[50] K. Reijnders, T. Tudorovskiy, and M. Katsnelson, Semiclassical theory of potential scattering for massless Dirac fermions, Ann. Phys. 333, 155 (2013)

[51] M. C. Gutzwiller, Periodic orbits and classical quantization conditions, J. Math. Phys. 12(3), 343 (1971)

[52] M. Wilkinson, An example of phase holonomy in WKB theory, J. Phys. A 17(18), 3459 (1984)

[53] R. Rammal, Landau level spectrum of Bloch electrons in a honeycomb lattice, J. Phys. France 46(8), 1345 (1985)

[54] G. P. Mikitik and Y. Sharlai, Manifestation of Berry’s phase in metal physics, Phys. Rev. Lett. 82(10), 2147 (1999)

[55] P. Carmier and D. Ullmo, Berry phase in graphene: Semiclassical perspective, Phys. Rev. B 77(24), 245413 (2008)

[56] J. N. Fuchs, F. Piéchon, M. O. Goerbig, and G. Montambaux, Topological Berry phase and semiclassical quantization of cyclotron orbits for two dimensional electrons in coupled band models., Eur. Phys. J. B 77(3), 351 (2010)

[57] Y. Gao and Q. Niu, Zero-field magnetic response functions in Landau levels, Proc. Natl. Acad. Sci. USA 114(28), 7295 (2017)

[58] A. A. Taskin, K. Segawa, and Y. Ando, Oscillatory angular dependence of the magnetoresistance in a topological insulator Bi_1-xSb_x, Phys. Rev. B 82(12), 121302 (2010)

[59] J. G. Analytis, R. D. McDonald, S. C. Riggs, J. H. Chu, G. S. Boebinger, and I. R. Fisher, Two-dimensional surface state in the quantum limit of a topological insulator, Nat. Phys. 6(12), 960 (2010)

[60] Z. Ren, A. A. Taskin, S. Sasaki, K. Segawa, and Y. Ando, Large bulk resistivity and surface quantum oscillations in the topological insulator Bi₂Te₂Se, Phys. Rev. B 82(24), 241306 (2010)

[61] B. Sacépé, J. B. Oostinga, J. Li, A. Ubaldini, N. J. Couto, E. Giannini, and A. F. Morpurgo, Gate-tuned normal and superconducting transport at the surface of a topological insulator, Nat. Commun. 2, 575 (2011)

[62] C. Brüne, C. X. Liu, E. G. Novik, E. M. Hankiewicz, H. Buhmann, Y. L. Chen, X. L. Qi, Z. X. Shen, S. C. Zhang, and L. W. Molenkamp, Quantum Hall effect from the topological surface states of strained bulk HgTe, Phys. Rev. Lett. 106(12), 126803 (2011)

[63] F. Xiu, L. He, Y. Wang, L. Cheng, L. T. Chang, M. Lang, G. Huang, X. Kou, Y. Zhou, X. Jiang, Z. Chen, J. Zou, A. Shailos, and K. L. Wang, Manipulating surface states in topological insulator nanoribbons, Nat. Nanotechnol. 6(4), 216 (2011)

[64] F. X. Xiang, X. L. Wang, M. Veldhorst, S. X. Dou, and M. S. Fuhrer, Observation of topological transition of Fermi surface from a spindle torus to a torus in bulk Rashba spin-split BiTeCl, Phys. Rev. B 92(3), 035123 (2015)

[65] G. Gómez-Santos, and T. Stauber, Measurable lattice effects on the charge and magnetic response in graphene, Phys. Rev. Lett. 106(4), 045504 (2011)

[66] A. Raoux, M. Morigi, J. N. Fuchs, F. Piéchon, and G. Montambaux, From dia- to paramagnetic orbital susceptibility of massless fermions, Phys. Rev. Lett. 112(2), 026402 (2014)

[67] M. Ogata and H. Fukuyama, Orbital magnetism of Bloch electrons (I): General formula, J. Phys. Soc. Jpn. 84(12), 124708 (2015)

[68] In Ref. [67], it is found that except the energy polarization contribution, all the other terms are consistent. For the energy polarization contribution, Ref. [41] contains a typo. When inserting the second order energy in Eq. (42), the energy polarization in Ref. [41]has an additional 1/4 factor by mistake. After removing such factor as given in Eq. (43), the energy polarization has the same expression with Eq. (2.31) in Ref. [67].

[69] N. Marzari and D. Vanderbilt, Maximally localized generalized Wannier functions for composite energy bands, Phys. Rev. B 56(20), 12847 (1997)

[70] I. Souza, N. Marzari, and D. Vanderbilt, Maximally localized Wannier functions for entangled energy bands, Phys. Rev. B 65(3), 035109 (2001)

[71] J. R. Yates, X. Wang, D. Vanderbilt, and I. Souza, Spectral and Fermi surface properties from Wannier interpolation, Phys. Rev. B 75(19), 195121 (2007)

[72] N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza, and D. Vanderbilt, Maximally localized Wannier functions: Theory and applications, Rev. Mod. Phys. 84(4), 1419 (2012)

[73] M. M. Vazifeh and M. Franz, Electromagnetic response of Weyl semimetals, Phys. Rev. Lett. 111(2), 027201 (2013)

[74] L. Fu, and E. Berg, Odd-parity topological superconductors: Theory and application to Cu_xBi₂Se₃, Phys. Rev. Lett. 105(9), 097001 (2010)

[75] H. K. Pal and D. L. Maslov, Necessary and sufficient condition for longitudinal magnetoresistance, Phys. Rev. B 81(21), 214438 (2010)

[76] K. Ohgushi, S. Murakami, and N. Nagaosa, Spin anisotropy and quantum Hall effect in the kagomé lattice: Chiral spin state based on a ferromagnet, Phys. Rev. B 62(10), R6065 (2000)

[77] R. Shindou and N. Nagaosa, Orbital ferromagnetism and anomalous Hall effect in antiferromagnets on the distorted fcc Lattice, Phys. Rev. Lett. 87(11), 116801 (2001)

[78] M. Taillefumier, B. Canals, C. Lacroix, V. K. Dugaev, and P. Bruno, Anomalous Hall effect due to spin chirality in the Kagomé lattice, Phys. Rev. B 74(8), 085105 (2006)

[79] A. Kalitsov, B. Canals, and C. Lacroix, Anomalous Hall effect due to magnetic chirality in the pyrochlore lattice, J. Phys. Conf. Ser. 145, 012020 (2009)

[80] H. Takatsu, S. Yonezawa, S. Fujimoto, and Y. Maeno, Unconventional anomalous Hall effect in the metallic triangular-lattice magnet PdCrO₂, Phys. Rev. Lett. 105(13), 137201 (2010)

[81] M. Udagawa and R. Moessner, Anomalous Hall effect from frustration-tuned scalar chirality distribution in Pr₂Ir₂O₇, Phys. Rev. Lett. 111(3), 036602 (2013)

[82] H. Chen, Q. Niu, and A. MacDonald, Anomalous Hall effect arising from noncollinear antiferromagnetism, Phys. Rev. Lett. 112(1), 017205 (2014)

[83] M. T. Suzuki, T. Koretsune, M. Ochi, and R. Arita, Cluster multipole theory for anomalous Hall effect in antiferromagnets, Phys. Rev. B 95(9), 094406 (2017)

[84] G. Y. Guo and T. C. Wang, Large anomalous Nernst and spin Nernst effects in the noncollinear antiferromagnets Mn₃X (X=Sn, Ge, Ga), Phys. Rev. B 96(22), 224415 (2017)

[85] L. Landau, E. Lifshitz, and L. Pitaevskii, Electrodynamics of Continuous Media, New York: Pergamon Press, 1984

[86] I. Sodemann and L. Fu, Quantum nonlinear Hall effect induced by berry curvature dipole in time-reversal invariant materials, Phys. Rev. Lett. 115(21), 216806 (2015)

[87] S. Y. Xu, Q. Ma, H. Shen, V. Fatemi, S. Wu, T. R. Chang, G. Chang, A. M. M. Valdivia, C. K. Chan, Q. D. Gibson, J. Zhou, Z. Liu, K. Watanabe, T. Taniguchi, H. Lin, R. J. Cava, L. Fu, N. Gedik, and P. Jarillo-Herrero, Electrically switchable Berry curvature dipole in the monolayer topological insulator WTe₂, Nat. Phys. 14(9), 900 (2018)

[88] K. Kang, T. Li, E. Sohn, J. Shan, and K. F. Mak, Observation of the nonlinear anomalous Hall effect in 2D WTe₂, arXiv: 1809.08744 (2018)

[89] Q. Ma, S.Y. Xu, H. Shen, D. Macneill, V. Fatemi, T.-R. Chang, A. M. M. Valdivia, S. Wu, Z. Du, C.-H. Hsu, S. Fang, Q. D. Gibson, K. Watanabe, T. Taniguchi, R. J. Cava, E. Kaxiras, H.-Z. Lu, H. Lin, L. Fu, N. Gedik, and P. Jarillo-Herrero, Observation of the nonlinear Hall effect under time reversal symmetric conditions, Nature 565, 337 (2019)

[90] A. Malashevich and I. Souza, Band theory of spatial dispersion in magnetoelectrics, Phys. Rev. B 82(24), 245118 (2010)

[91] J. Orenstein and J. E. Moore, Berry phase mechanism for optical gyro ropy in stripe-ordered cuprates, Phys. Rev. B 87(16), 165110 (2013)

[92] S. Nandy and I. Sodemann, Symmetry and quantum kinetics of the non-linear Hall effect, arXiv: 1901.04467 (2019)

[93] J. R. Reitz, and A. W. Overhauser, Magnetoresistance of potassium, Phys. Rev. 171(3), 749 (1968)

[94] P. A. Penz and R. Bowers, Strain-dependent magnetoresistance of potassium, Phys. Rev. 172(3), 991 (1968)

[95] B. K. Jones, Strain-dependent magnetoresistance of sodium and potassium, Phys. Rev. 179(3), 637 (1969)

[96] A. L. Friedman, J. L. Tedesco, P. M. Campbell, J. C. Culbertson, E. Aifer, F. K. Perkins, R. L. Myers-Ward, J. K. Hite, C. R. Jr Eddy, G. G. Jernigan, and D. K. Gaskill, Quantum linear magnetoresistance in multilayer epitaxial graphene, Nano Lett. 10(10), 3962 (2010)

[97] D. X. Qu, Y. S. Hor, J. Xiong, R. J. Cava, and N. P. Ong, Quantum oscillations and Hall anomaly of surface states in the topological insulator Bi₂Te₃, Science 329(5993), 821 (2010)

[98] H. Tang, D. Liang, R. L. J. Qiu, and X. P. A. Gao, Two-dimensional transport-induced linear magneto-resistance in topological insulator Bi₂Se₃ nanoribbons, ACS Nano 5(9), 7510 (2011)

[99] X. Wang, Y. Du, S. Dou, and C. Zhang, Room temperature giant and linear magnetoresistance in topological insulator Bi₂Te₃ nanosheets, Phys. Rev. Lett. 108(26), 266806 (2012)

[100] J. Tian, C. Chang, H. Cao, K. He, X. Ma, Q. Xue, and Y. P. Chen, Quantum and classical magnetoresistance in ambipolar topological insulator transistors with gate-tunable bulk and surface conduction, Sci. Rep. 4(1), 4859 (2015)

[101] L. He, X. Hong, J. Dong, J. Pan, Z. Zhang, J. Zhang, and S. Li, Quantum transport evidence for the three-dimensional Dirac semimetal phase in Cd₃As₂, Phys. Rev. Lett. 113(24), 246402 (2014)

[102] T. Liang, Q. Gibson, M. N. Ali, M. Liu, R. J. Cava, and N. P. Ong, Ultrahigh mobility and giant magnetoresistance in the Dirac semimetal Cd₃As₂, Nat. Mater. 14(3), 280 (2015)

[103] J. Feng, Y. Pang, D. Wu, Z. Wang, H. Weng, J. Li, X. Dai, Z. Fang, Y. Shi, and L. Lu, Large linear magnetoresistance in Dirac semimetal Cd₃As₂ with Fermi surfaces close to the Dirac points, Phys. Rev. B 92(8), 081306 (2015)

[104] M. Novak, S. Sasaki, K. Segawa, and Y. Ando, Large linear magnetoresistance in the Dirac semimetal TlBiSSe, Phys. Rev. B 91(4), 041203 (2015)

[105] A. Narayanan, M. Watson, S. Blake, N. Bruyant, L. Drigo, Y. Chen, D. Prabhakaran, B. Yan, C. Felser, T. Kong, P. Canfield, and A. Coldea, Linear magnetoresistance caused by mobility fluctuations in n-doped Cd₃As₂, Phys. Rev. Lett. 114(11), 117201 (2015)

[106] X. Huang, L. Zhao, Y. Long, P. Wang, D. Chen, Z. Yang, H. Liang, M. Xue, H. Weng, Z. Fang, X. Dai, and G. Chen, Observation of the chiral-anomaly-induced negative magnetoresistance in 3D Weyl semimetal TaAs, Phys. Rev. X 5(3), 031023 (2015)

[107] A. A. Abrikosov, Quantum magnetoresistance, Phys. Rev. B 58(5), 2788 (1998)

[108] A. A. Abrikosov, Quantum linear magnetoresistance, Europhys. Lett. 49(6), 789 (2000)

[109] C. M. Wang and X. L. Lei, Linear magnetoresistance on the topological surface, Phys. Rev. B 86(3), 035442 (2012)

[110] C. Herring, Effect of random inhomogeneities on electrical and galvanomagnetic measurements, J. Appl. Phys. 31(11), 1939 (1960)

[111] M. M. Parish and P. B. Littlewood, Non-saturating magnetoresistance in heavily disordered semiconductors, Nature 426(6963), 162 (2003)

[112] N. A. Porter and C. H. Marrows, Linear magnetoresistance in n-type silicon due to doping density fluctuations, Sci. Rep. 2(1), 565 (2012)

[113] N. Kozlova, N. Mori, O. Makarovsky, L. Eaves, Q. Zhuang, A. Krier, and A. Patanè, Linear magnetoresistance due to multiple-electron scattering by low-mobility islands in an inhomogeneous conductor, Nat. Commun. 3(1), 1097 (2012)

[114] H. Chen, Y. Gao, D. Xiao, A. H. MacDonald, and Q. Niu, Semiclassical theory of linear magnetoresistance in crystalline conductors with broken time-reversal symmetry, arXiv: 1511.02557 (2015)

[115] A. B. Pippard, Magnetoresistance in Metals, Cambridge University Press, Cambridge, England, New York, 1989

[116] H. B. Nielsen and M. Ninomiya, The Adler–Bell–Jackiw anomaly and Weyl fermions in a crystal, Phys. Lett. B 130(6), 389 (1983)

[117] V. Aji, Adler–Bell–Jackiw anomaly in Weyl semimetals: Application to pyrochlore iridates, Phys. Rev. B 85(24), 241101 (2012)

[118] A. A. Burkov, Chiral anomaly and transport in Weyl metals, J. Phys.: Condens. Matter 27(11), 113201 (2015)

[119] M. A. Stephanov, and Y. Yin, Chiral kinetic theory, Phys. Rev. Lett. 109(16), 162001 (2012)

[120] D. T. Son and N. Yamamoto, Kinetic theory with Berry curvature from quantum field theories, Phys. Rev. D 87(8), 085016 (2013)

[121] B. Z. Spivak and A. V. Andreev, Magnetotransport phenomena related to the chiral anomaly in Weyl semimetals, Phys. Rev. B 93(8), 085107 (2016)

[122] Y. Hidaka, S. Pu, and D. L. Yang, Relativistic chiral kinetic theory from quantum field theories, Phys. Rev. D 95(9), 091901 (2017)

[123] A. Sekine, D. Culcer, and A. H. MacDonald, Quantum kinetic theory of the chiral anomaly, Phys. Rev. B 96(23), 235134 (2017)

[124] A. A. Burkov, Chiral anomaly and diffusive magnetotransport in Weyl metals, Phys. Rev. Lett. 113(24), 247203 (2014)

[125] A. A. Burkov, Negative longitudinal magnetoresistance in Dirac and Weyl metals, Phys. Rev. B 91(24), 245157 (2015)

[126] A. Andreev and B. Spivak, Longitudinal negative magnetoresistance and magnetotransport phenomena in conventional and topological conductors, Phys. Rev. Lett. 120(2), 026601 (2018)

[127] C. Shekhar, A. K. Nayak, Y. Sun, M. Schmidt, M. Nicklas, I. Leermakers, U. Zeitler, Y. Skourski, J. Wosnitza, Z. Liu, Y. Chen, W. Schnelle, H. Borrmann, Y. Grin, C. Felser, and B. Yan, Extremely large magnetoresistance and ultrahigh mobility in the topological Weyl semimetal candidate NbP, Nat. Phys. 11(8), 645 (2015)

[128] X. Yang, Y. Li, Z. Wang, Y. Zheng, and Z. A. Xu, Chiral anomaly induced negative magnetoresistance in topological Weyl semimetal NbAs, arXiv: 1506.03190 (2015)

[129] Z. Wang, Y. Zheng, Z. Shen, Y. Lu, H. Fang, F. Sheng, Y. Zhou, X. Yang, Y. Li, C. Feng, and Z. A. Xu, Helicity-protected ultrahigh mobility Weyl fermions in NbP, Phys. Rev. B 93(12), 121112 (2016)

[130] C. L. Zhang, S. Y. Xu, I. Belopolski, Z. Yuan, Z. Lin, B. Tong, G. Bian, N. Alidoust, C. C. Lee, S. M. Huang, T. R. Chang, G. Chang, C. H. Hsu, H. T. Jeng, M. Neupane, D. S. Sanchez, H. Zheng, J. Wang, H. Lin, C. Zhang, H. Z. Lu, S. Q. Shen, T. Neupert, M. Z. Hasan, and S. Jia, Signatures of the Adler–Bell–Jackiw chiral anomaly in a Weyl fermion semimetal, Nat. Commun. 7(1), 10735 (2016)

[131] H. J. Kim, K. S. Kim, J. F. Wang, M. Sasaki, N. Satoh, A. Ohnishi, M. Kitaura, M. Yang, and L. Li, Dirac versus Weyl fermions in topological insulators: Adler-Bell-Jackiw anomaly in transport phenomena, Phys. Rev. Lett. 111(24), 246603 (2013)

[132] J. Xiong, S. K. Kushwaha, T. Liang, J. W. Krizan, M. Hirschberger, W. Wang, R. Cava, and N. Ong, Evidence for the chiral anomaly in the Dirac semimetal Na₃ Bi, Science 350(6259), 413 (2015)

[133] J. Feng, Y. Pang, D. Wu, Z. Wang, H. Weng, J. Li, X. Dai, Z. Fang, Y. Shi, and L. Lu, Large linear magnetoresistance in Dirac semimetal Cd₃As₂ with Fermi surfaces close to the Dirac points, Phys. Rev. B 92(8), 081306 (2015)

[134] C. Z. Li, L. X. Wang, H. Liu, J. Wang, Z. M. Liao, and D. P. Yu, Giant negative magnetoresistance induced by the chiral anomaly in individual Cd₃As₂ nanowires, Nat. Commun. 6(1), 10137 (2015)

[135] H. Li, H. He, H. Z. Lu, H. Zhang, H. Liu, R. Ma, Z. Fan, S. Q. Shen, and J. Wang, Negative magnetoresistance in Dirac semimetal Cd₃As₂, Nat. Commun. 7(1), 10301 (2016)

[136] C. Zhang, E. Zhang, W. Wang, Y. Liu, Z. G. Chen, S. Lu, S. Liang, J. Cao, X. Yuan, L. Tang, Q. Li, C. Zhou, T. Gu, Y. Wu, J. Zou, and F. Xiu, Room-temperature chiral charge pumping in Dirac semimetals, Nat. Commun. 8, 13741 (2017)

[137] Y. Gao, S. A. Yang, and Q. Niu, Intrinsic relative magnetoconductivity of nonmagnetic metals, Phys. Rev. B 95(16), 165135 (2017)

[138] X. Dai, Z. Du, and H. Z. Lu, Negative magnetoresistance without chiral anomaly in topological insulators, Phys. Rev. Lett. 119(16), 166601 (2017)

[139] H. W. Wang, B. Fu, and S. Q. Shen, Intrinsic magnetoresistance in three-dimensional Dirac materials with low carrier density, Phys. Rev. B 98(8), 081202 (2018)

[140] K. Yoshida, Transport of spatially inhomogeneous current in a compensated metal under magnetic fields (III): A case of bismuth in longitudinal and transverse magnetic fields, J. Appl. Phys. 51(8), 4226 (1980)

[141] R. D. Reis, M. O. Ajeesh, N. Kumar, F. Arnold, C. Shekhar, M. Naumann, M. Schmidt, M. Nicklas, and E. Hassinger, On the search for the chiral anomaly in Weyl semimetals: The negative longitudinal magnetoresistance., New J. Phys. 18(8), 085006 (2016)

[142] Z. Yuan, H. Lu, Y. Liu, J. Wang, and S. Jia, Large magnetoresistance in compensated semimetals TaAs₂ and NbAs₂, Phys. Rev. B 93(18), 184405 (2016)

[143] R. Shindou and L. Balents, Artificial electric field in Fermi liquids, Phys. Rev. Lett. 97(21), 216601 (2006)

[144] R. Shindou and L. Balents, Gradient expansion approach to multiple-band Fermi liquids, Phys. Rev. B 77(3), 035110 (2008)

[145] D. Culcer, A. Sekine, and A. H. MacDonald, Interband coherence response to electric fields in crystals: Berry-phase contributions and disorder effects, Phys. Rev. B 96(3), 035106 (2017)

[146] Y. Tian, L. Ye, and X. Jin, Proper scaling of the anomalous Hall effect, Phys. Rev. Lett. 103(8), 087206 (2009)

[147] S. H. Chun, Y. S. Kim, H. K. Choi, I. T. Jeong, W. O. Lee, K. S. Suh, Y. S. Oh, K. H. Kim, Z. G. Khim, J. C. Woo, and Y. D. Park, Interplay between carrier and impurity concentrations in annealed Ga_1-xMn_xAs: Intrinsic anomalous Hall effect, Phys. Rev. Lett. 98(2), 026601 (2007)

[148] W. L. Lee, S. Watauchi, V. L. Miller, R. J. Cava, and N. P. Ong, Dissipationless Anomalous Hall Current in the Ferromagnetic Spinel CuCr₂Se_4-xBr_x, Science 303(5664), 1647 (2004)

[149] R. Mathieu, A. Asamitsu, H. Yamada, K. S. Takahashi, M. Kawasaki, Z. Fang, N. Nagaosa, and Y. Tokura, Scaling of the anomalous Hall effect in Sr_1-xCa_xRuO₃, Phys. Rev. Lett. 93(1), 016602 (2004)

[150] B. C. Sales, R. Jin, D. Mandrus, and P. Khalifah, Anomalous Hall effect in three ferromagnetic compounds: EuFe₄Sb₁₂, Yb₁₄MnSb₁₁, and Eu₈Ga₁₆Ge₃₀, Phys. Rev. B 73(22), 224435 (2006)

[151] C. Zeng, Y. Yao, Q. Niu, and H. H. Weitering, Linear magnetization dependence of the intrinsic anomalous Hall effect, Phys. Rev. Lett. 96(3), 037204 (2006)

Rights & permissions

Higher Education Press and Springer-Verlag GmbH Germany, part of Springer Nature 2019