Quantum Transport and Magnetism of Dirac Electrons in Solids

The relativistic Dirac equation covers the fundamentals of electronic phenomena in solids and as such it effectively describes the electronic states of the topological insulators like Bi$_2$Se$_3$ and Bi$_2$Te$_3$. Topological insulators feature gapless surface states and, moreover, magnetic doping and resultant ferromagnetic ordering break time-reversal symmetry to realize quantum anomalous Hall and Chern insulators. Here we focus on the bulk and investigate the mutual coupling of electronic and magnetic properties of Dirac electrons. Without carrier doping, spiral magnetic orders cause a ferroelectric polarization through the spin-orbit coupling. In a doped metallic state, the anisotropic magnetoresistance arises without uniform magnetization. We find that electric current induces uniform magnetization and conversely an oscillating magnetic order induces electric current. Our model provides a coherent and unified description of all those phenomena. The mutual control of electric and magnetic properties demonstrates implementations of antiferromagnetic spintronics. We also discuss the stoichiometric magnetic topological insulator MnBi$_2$Te$_4$.

The relativistic Dirac equation covers the fundamentals of electronic phenomena in solids and as such it effectively describes the electronic states of the topological insulators like Bi2Se3 and Bi2Te3. Topological insulators feature gapless surface states and, moreover, magnetic doping and resultant ferromagnetic ordering break time-reversal symmetry to realize quantum anomalous Hall and Chern insulators. Here we focus on the bulk and investigate the mutual coupling of electronic and magnetic properties of Dirac electrons. Without carrier doping, spiral magnetic orders cause a ferroelectric polarization through the spin-orbit coupling. In a doped metallic state, the anisotropic magnetoresistance arises without uniform magnetization. We find that electric current induces uniform magnetization and conversely an oscillating magnetic order induces electric current. Our model provides a coherent and unified description of all those phenomena. The mutual control of electric and magnetic properties demonstrates implementations of antiferromagnetic spintronics. We also discuss the stoichiometric magnetic topological insulator MnBi2Te4.
Relativistic effects on electrons, exemplified by the spin-orbit coupling, mingle the spin and orbital degrees of freedom and bring the interplay between electric and magnetic properties. Multiferroics is a manifestation in insulators, where a magnetization induces an electric polarization and vice versa [1][2][3][4]. In metals and semiconductors, the spin-orbit coupling enables control of electrons' spin from electric current, and it is essential for spintronics. For example, the Rashba spin-orbit coupling causes the Edelstein effect [5], which produces spin polarization by electric current in inversion-breaking systems. Spintronics conventionally utilizes ferromagnets. Antiferromagnetic spintronics recently has gained more interest owing to various advantages such as fast response, no stray field, and large magnetotransport effects [6][7][8]. However, because the net magnetization vanishes in an antiferromagnet, manipulation and detection remain essential challenges.
A magnetization pattern in general configures a spiral order with strong correlation or with magnetic elements with a fixed magnetic moment. Magnetism breaks timereversal symmetry T even though a spiral magnetic order may have no net magnetic moment. In addition, a spiral order is characterized by a wavevector Q and often breaks inversion symmetry P regardless of the underlying crystalline symmetry.
We study various phenomena related to broken T and P symmetries in magnetic Dirac materials in a unified fashion. The spin-orbit coupling naturally arises from the Dirac equation; as it abides by relativity, the coupling between the electric and magnetic degrees of freedom is contained. There are various materials where the Dirac Hamiltonian becomes the effective model near the chemical potential. Examples are the three-dimensional topological insulators (TIs) Bi 2 Se 3 and Bi 2 Te 3 [9][10][11][12][13][14]. With an insulating bulk, the topologically protected surface states determine the physical properties, which have been extensively studied [15,16]. In the doped case, however, the bulk states dominate the electric and magnetic properties of the system. When a magnetic order is present, the bulk of TIs offer an ideal laboratory to study the Dirac electrons with the exchange coupling to the magnetic moments.
In this work, we consider the electromagnetic response of a gapped Dirac system coupled to local magnetic moments. Our model describes the magnetically doped TIs Bi 2 Se 3 and Bi 2 Te 3 [17][18][19][20][21], where the magnetic dopants couple locally to the Dirac electrons via the exchange coupling. We first show that an inversion-breaking magnetic order can generate a finite electric polarization in the insulating state while the pristine electronic system is centrosymmetric. In a doped metallic state, we reveal that a magnetic order can induce anisotropic resistance. In addition, an electric field produces a uniform magnetization and in reverse an oscillating magnetic order generates direct current. We also discuss the intrinsic magnetic TI MnBi 2 Te 4 [22][23][24][25] and the possibility of inversion-breaking magnetic orders in magnetic TIs.
Model : We consider a three-dimensional isotropic gapped Dirac system. Such an electronic system is realized, for example, in the bulk of TIs. For the TIs Bi 2 Se 3 and Bi 2 Te 3 , the energy bands near the Γ point describe the low-energy behavior, which consists of the spin σ and p orbitals τ from Bi (τ z = +1) and Se/Te (τ z = −1). To linear order in momentum k, the k · p Hamiltonian becomes where the 4 × 4 matrices α = στ x and β = τ z satisfy the anticommutation relations {α a , α b } = {α a , β} = 0 (a = b) and α 2 a = β 2 = I (I: identity matrix) [10]. We set = 1. The pristine system preserves inversion P = τ z and time reversal T = iσ y K with the complex conjugate operator K: PH 0 (k)P −1 = H 0 (−k) and T H 0 (k)T −1 = H 0 (−k). The kinetic term renders the spin and orbital coupling, so that neither is a good quantum number. The sign of the mass can be either positive or negative, which describes the band inversion near the Γ point.
Magnetic dopants such as Mn, Cr, and Fe can substitute the Bi sites of Bi 2 Se 3 and Bi 2 Te 3 [26]. Their  (7). The left panel is a three-dimensional illustration of the anisotropic resistance R , and the center and right panels are the two distinct plane cuts, displaying the anisotropy in the plane perpendicular to the local magnetic moments. (c) Uniform direct current induced by oscillating magnetic orders. For oscillations forming cycloidal and proper screw patterns, the generated direct current is parallel to the wavevector of the magnetic orders, following Eq. (9).
local magnetic moments break time-reversal symmetry and tend to form a magnetic order. In a metallic state, the Ruderman-Kittel-Kasuya-Yosida (RKKY) interaction favors ferromagnetism when the Fermi level is near the Dirac point, but in general a complex magnetic order may occur depending on the Fermi level, anisotropy, and inhomogeneity [27][28][29][30]. An effective spin Hamiltonian reflecting such details of the system determines the magnetic order M (r) = Q M Q e iQ·r , which we take as given in the following analyses.
The exchange coupling yields the local magnetic coupling to the Dirac electrons. We note that the exchange coupling is orbital dependent [31]: Here, we introduceβ = τ z sgn(m) for later convenience. The two coupling constants J and J describe the different strengths of the exchange coupling for the two orbitals (τ z = ±1).
Insulating state: The bulk is insulating when the chemical potential lies inside the mass gap. While the electronic system preserves inversion, the magnetic order may violate it, allowing a finite electric polarization. The calculation of the polarization follows the method by King-Smith and Vanderbilt [32]. We find an inversionbreaking magnetic order produces a finite polarization of the Dirac electrons see Supplemental Material (SM) for details [33]. The result conforms to the analyses of a Ginzburg-Landau model [34] and a microscopic model [35], where certain chiral magnetic orders induce a finite electric polarization. For cycloidal and proper screw orders, only the former induce a finite polarization perpendicular to the wavevector in the magnetization plane according to Eq. (3) [ Fig. 1(a)]. The product JJ implies that the strengths of the exchange coupling should be different for the two orbitals for a finite polarization. The orbitaldependent exchange coupling mixes the conduction and valence bands by the magnetic order to realize a finite polarization.
Effective Hamiltonian in the metallic state: When the system is metallic, we expect various responses to an external electromagnetic field. As charges in the vicinity of the Fermi surface are dominantly responsible to electromagnetic response, it is convenient to derive the effective Hamiltonian for the bands that cross the Fermi energy. We obtain the effective Hamiltonian by following the method by Foldy and Wouthuysen [36], and Tani [37], which we can calculate as a perturbative series in the large mass limit |m| | F | ( F : the Fermi energy measured from a band edge) [33]. In the presence of an external electromagnetic field, the effective Hamiltonian to order m −2 is with Π = p + eA and the momentum operator p = −i∇.
The charge of an electrons is −e. The electric and magnetic fields are E = −∇Φ − ∂A/∂t and B = ∇ × A, respectively, with the scalar potential Φ and the vector potential A. In the effective Hamiltonian,β = ±1 signifies the energy bands:β = +1 corresponds to the conduction band andβ = −1 to the valence band. Although we originally defineβ = τ z sgn(m), it does not precisely label the orbitals after the unitary transformation. The last two terms of the effective Hamiltonian (4) reveal the nontrivial coupling between the Dirac electrons and the magnetic order, which is central to the following results. It manifests the strong spin-orbital coupling embedded in the Dirac Hamiltonian along with the exchange coupling. It also modifies the current density operator at zero frequency. The second term with ∇ × M has a classical analog to the Ampère's circuital law. The third term contains the local magnetic moment M and the spin of the Dirac electrons σ. It implies the possibility of the mutual control of the electric and magnetic degrees of freedom as we will see below.
Current under an electric field : We perform perturbative calculations using functional derivatives to calculate response. We define the action S = T ωn drψ(−iω n + H eff )ψ, where T is the temperature and ω n = (2n + 1)πT is the fermionic Matsubara frequency. Using the partition function Z = DψDψe −S , we obtain the current response in the presence of an external electric field E(ω) = iωA(ω) (Φ = 0) as j a (ω) = ĵ a (ω) = 1 iω We note that it is equivalent to the Kubo formula. We calculate it perturbatively with respect to the exchange couplings J, J , and the inverse mass m −1 , using the unperturbed Green's function G 0 = (ω − H 0 eff − Σ) −1 with H 0 eff = mβ + βk 2 /(2m) and the self-energy Σ. We approximate Σ ≈ −i sgn(ω n )/(2τ ) with a constant τ to describe momentum relaxation in diffusive transport.
The magnetic order alters the current flow. When we focus on a spatially uniform current, the lowest-order corrections by the magnetic order appear as a product of M Q and M * Q . By differentiating Eq. (6) with respect to M Q and M * Q , we obtain the conductivity tensor where the coefficients are given by n( F ) ∝ | F | 3/2 is the carrier density (n > 0 for electrons and n < 0 for holes). We introduce τ ω = τ /(1 − iωτ ) and retain the leading-order contributions in m −1 in the expressions of σ, η, and η . When there is a uniform magnetization M 0 = 0, it yields the anomalous Hall contribution σ AH ab ∝ ε abc M 0,c ; see SM for details [33]. j(ω) depends only on M Q but does not directly depend on Q to this order.
The spatial pattern of the magnetic order modifies the conductivity at second order in M . The first correction with η(ω) reduces the longitudinal conductivity, arising from the exchange coupling −(J + J β )M · σ. The effect is isotropic and it does not require the spin-orbital coupling inherent in the Dirac Hamiltonian. It resembles the magnetoresistance whereas there is no uniform magnetization by assumption. On the other hand, the η term can be traced to the coupling between the magnetic order and current, as we have seen in Eq. (5). It gives rise to anisotropic corrections to the conductivity tensor σ ab (ω) depending on the magnetic order.
The second-order corrections to the conductivity correspond to the anisotropic magnetoresistance and the planar Hall effect [38]. Both cycloidal and screw magnetic orders show the anisotropic resistance R [ Fig. 1(b)]: when the magnetic order lies in the xy plane, the resistance is different in the xy plane and along the z axis. We emphasize, however, that the second-order effect in Eq. (7) appears even without a uniform magnetization. Therefore, when there is no uniform magnetization, namely, σ AH ab = 0, the conductivity tensor is symmetric: σ ab (ω) = σ ba (ω). On the other hand, the anomalous Hall contribution is antisymmetric: σ AH ab (ω) = −σ AH ba (ω) [39]. Magnetization by an electric field : From the coupling between the current and the spin degrees of freedom, we expect that an electric field produces a finite magnetization of Dirac electrons even when the magnetic order has no uniform magnetization. We evaluate the spin expectation value of the Dirac electrons σ in the presence of an external electric field E and the magnetic order M using H eff . The uniform magnetization of the Dirac electrons is given by m Dirac = −gµ B σ /2, where g is the g-factor and µ B is the Bohr magneton. A perturbative calculation finds a finite magnetization under a static and uniform external electric field [33] with Since the magnetization and the electric field transforms differently under inversion, an inversion-breaking magnetic order is necessary to induce magnetization by an electric field [ Fig. 1(a)]. It allows detection of an inversion-breaking magnetic order through the magnetization by applying an electric field. The change of the magnetization under an electric field can be attributed to m Dirac . The effect resembles the Edelstein effect but it appears in the bulk of a TI, where inversion is broken by a magnetic order. The extension to a time-dependence case is straightforward [33]. Current by an oscillating magnetic order : We now investigate whether an external magnetic field induces an electric current. The magnetic field should vary in time as a spatially uniform current cannot exist in the equilibrium. The external magnetic field applied to a metallic system with a magnetic order has the following two effects: it couples to the itinerant electrons to induce cyclotron motion; at the same time, it drives the Rabi oscillation and the Larmor precession of the local magnetic moments.
We first check if a uniform oscillating magnetic field B(ω) induces a uniform current in the presence of a static magnetic order. From the symmetry consideration, the lowest-order contribution should have the form j a (ω) = κ abcd B b (ω)M Q,c M −Q,d with κ abcd linear in Q. However, this mechanism is improbable. The conductivity tensor Eq. (7) is insensitive to inversion breaking, so that the cyclotron motion of the Dirac electrons would not yield a uniform current. We calculate κ ijkl perturbatively and observe that it vanishes to order QJ 2 n( F )/m 2 [33].
We then examine current response by an oscillating magnetic order. If it is finite, an external magnetic field induces an electric current by making the local magnetic moments oscillate. We write the spatial and temporal dependence of the magnetic order as M (r, t) = Qω M Qω e i(Q·r−ωt) . Here, we seek the uniform current response of the form where γ abc is linear in the wavevector Q to capture the inversion breaking by the magnetic order and hence to comply with the symmetry constraint. As a second-order response, the output frequency is the sum of two input frequencies. We can calculate the current response similarly to σ ab (ω) [33]: where we denote M 1 = M Qω1 , M 2 = M −Qω2 , and the coefficients are γ (S) (ω 1 , ω 2 ) and γ (A) (ω 1 , ω 2 ) are symmetric and antisymmetric under the exchange of ω 1 and ω 2 , respectively. γ (A) (ω, −ω) corresponds to zero-frequency response, namely direct current depicted in Fig. 1(c), and γ (S) (ω, ω) to 2ω response. It is worth contrasting the current response in the metallic state (9) with the polarization in the insulating state (3) as they reflect different material properties. First, the current response requires dynamics of the magnetic order whereas the polarization is a thermodynamic quantity defined in the equilibrium. The diffusive nature of the current is manifested in the appearance of the lifetime τ . Second, the current is carried by electric charges near the Fermi energy and it is thus proportional to the carrier density. On the other hand, the polarization only involves the quantities that characterize the system, implying that it requires the information of the entire band structure. Indeed, we cannot obtain Eq. (3) from the effective Hamiltonian (4) but from the original model (1).
As we have discussed, the local magnetic moments oscillate under a time-dependent external magnetic field to induce a uniform electric current. When the oscillation is near resonance, we may expect a larger current response. Since it is a second-order response with respect to the magnetic order, the response should be peaked at the zero frequency and double the resonance frequency. A magnetic order might also be driven by the spin wave spectroscopy technique [40][41][42]. An oscillating magnetic field is induced by periodically aligned wave guides whereby the wavevector of the magnetic field is designed.
Discussions: We have revealed that electromagnetic response of a magnetic TI manifests the entanglement of the spin and orbital degrees of freedom and hence the electric and magnetic properties. Particularly with an inversion-breaking magnetic order, it allows a measurement of electric properties through a magnetic probe and vice versa, and suggest applications in spintronics.
In addition to the magnetically doped TIs, we can consider the stoichiometric magnetic TI MnBi 2 Te 4 . It consists of stacking layers of TI films, bound by the van der Waals interaction [22]. The low-energy effective Hamiltonian is H STI (k) = mτ x +vτ z (ẑ ×σ)·k ⊥ +v z k z τ y , where the stacking direction is set along the z direction and τ z corresponds to the top and bottom TI surface states of a constituent layer [43]. In the SM [33], we confirm that an inversion-breaking magnetic order induces an electric polarization in the insulating state, and derive the effective Hamiltonian for the metallic case to see that the current operator is affected by the magnetic order.
Experimentally, a spiral magnetic order has not yet been reported in magnetic TIs, but yet some experiments reveal noncollinear magnetic orders. Stacking layers of MnBi 2 Te 4 realize a canted antiferromagnetic order [44], and alternating stacks of MnBi 2 Te 4 and Bi 2 Te 3 lead to a variety of heterostructures (MnBi 2 Te 4 ) m (Bi 2 Te 3 ) n [45].
The topological Hall effect is observed in the magnetic/non-magnetic topological insulator heterostructures Cr x (Bi 1−y Sb y ) 2−x Te 3 /(Bi 1−y Sb y ) 2 Te 3 and a theory attributed its origin to a Néel-type skyrmion, consisting of the superposition of the local three spiral orders [46]. The topological Hall effect attributed to skyrmions is also observed in Mn-doped Bi 2 Te 3 topological insulator films [47]. Those observations suggest that various magnetic orders may appear by different stacks and material compositions.
We now estimate the magnitude of the effects that we have discussed using the material parameters of Cr x (Bi 1−y Sb y ) 2−x Te 3 [46]: m = −300 meV, J = −5 meV, J = 1 meV, and the velocity v = 5.0 × 10 5 m/s; see SM for details [33]. The magnetic moment per Cr atom is M ≈ 3µ B . We set F = −100 meV. The RKKY interaction would form a magnetic order in the metallic state with the wavenumber Q = 2k F ≈ 1.5 × 10 9 m −1 . We estimate τ ≈ 5 × 10 −15 s from the longitudinal conductivity 100 Ω −1 cm −1 with |n| ≈ 1.4 × 10 19 cm −3 . Then, the corrections to the conductivity are −2ηM 2 ≈ −8.7 Ω −1 cm −1 and 4η M 2 ≈ −0.4 Ω −1 cm −1 for the isotropic and anisotropic parts, respectively. The magnetization induced by the current density j = 10 8 A/m is m Dirac ∼ 10 −4 A/m. The current densities generated by an oscillating magnetic order at 1 GHz are γ (S) QM 2 ≈ 9.2 A/m 2 for the sum frequency generation and γ (A) JQM 2 ≈ −2.3 × 10 5 A/m 2 . We note that the former grows quadratically with frequency while the latter does linearly. In the insulating state, the electric polarization is ∆P ≈ 1.8 µC/m 2 with the same Q. From those estimates, the electronic response is more likely to be observable that the magnetic one.
In addition to magnetic TIs, we also anticipate similar current response in magnetic Weyl and Dirac semimetals, where an emergent electromagnetic field plays a role as well as the Berry curvature [48][49][50]. A surface, an interface, and a domain wall geometrically break inversion, and thus the existence of a magnetic order can also induce various response. Such structures without inversion support the Dzyaloshinskii-Moriya interaction, which could contribute to a chiral magnetic order to reveal the effects that we have discussed.

Supplemental Material
In Supplemental Material (SM), we describe the unitary transformation in the large mass limit, the details about the model for magnetically-doped topological insulators (TIs), the calculations of various response in the metallic state and polarization in the insulating state. We also include the model and analysis of the model for the stoichiometric TI MnBi 2 Te 4 . We set = 1 unless otherwise noted.

S1. EFFECTIVE HAMILTONIAN FROM A UNITARY TRANSFORMATION
We derive the effective Hamiltonian for the Dirac system coupled to the magnetic order. To this end, we perform a unitary (Foldy-Wouthuysen-Tani) transformation that diagonalizes the Hamiltonian in the orbital space in the large Dirac mass limit.
We decompose the Hamiltonian as where E and O are the diagonal and off-diagonal terms in the orbital components: H denotes the exchange coupling of the Dirac electrons to the magnetic order. Here, we include the electromagnetic potential (Φ, A) and minimal coupling gives the canonical momentum Π = p + eA. We note that the diagonal part commutes with the matrix β and that the off-diagonal part anticommutes with β: which the following algebra relies on. An unitary transformation e iS1 with a Hermitian operator S 1 converts the Hamiltonian into E 1 and O 1 are the diagonal and off-diagonal terms in the orbital components after the unitary transformation. We determine S 1 to remove the off-diagonal terms at order m 0 , requiring The condition leads to and thus the series expansion of the unitary transformation corresponds to the expansion with respect to the inverse of the Dirac mass m −1 . After performing the unitary transformation e iS1 , we obtain and the diagonal and off-diagonal terms in the orbital components E 1 , O 1 are (S10) E 1 and O 1 again satisfy the same commutation and anticommutation relations We can iterate the same procedure to eliminate off-diagonal components at every order in m. At the j-th repetition, the unitary transformation e iSj leads to the Hamiltonian where the Hermite operator S j is Here we consider H 4 , which is diagonal to order m −3 : To further calculate the expression, we define the matrix Σ as where ε abc is the Levi-Civita symbol. Then, we obtain The last term apparently contains α but can be eliminated; the explicit form depends on the commutation relation between α and H . Therefore, we obtain the formal expression of the effective Hamiltonian as

S2. ISOTROPIC TI MODEL
For the topological insulators Bi 2 Se 3 and Bi 2 Te 3 , the k · p expansion to linear order in momentum k around the Γ point becomes the Dirac Hamiltonian [S1] σ and τ are the Pauli matrices for the spin and orbital degrees of freedom, respectively. τ z = ±1 corresponds to the cation and anion p orbitals. The model satisfy time-reversal, inversion, and three-fold rotational symmetries T = iσ y K, P = τ z , C 3 = exp(iπσ z /3), respectively, where K denotes the complex conjugate operator.
With the rescaling of the momentum k, we can eliminate the coefficients A 1 and A 2 . In addition, we include the electromagnetic potential (Φ, A), corresponding to the electric field E = −∇Φ −Ȧ and the magnetic field B = ∇ × A. The electromagnetic potential replaces the momentum operator p = −i∇ with the gauge-invariant momentum operator Π = p + eA from minimal coupling. Here, the charge of an electron is −e (e > 0). As a result, the Dirac Hamiltonian becomes We do not explicitly write the chemical potential hereafter. β and α are the 4 × 4 matrices satisfying the relations where I is the identity matrix. For the present model, α and β are We suppose doping of magnetic impurities, which couple to the Dirac electrons via the exchange coupling. Since the model consists of the two p orbitals of different origins, the strength of the exchange coupling depends on the orbitals. Thus, the magnetic order M (r, t) affects the Dirac electrons in the form [S2] We include the sign of the mass in the term with J . Since the the matrix β appears with the mass m, its eigenvaluẽ β = ±1 signifies the conduction or valence band in the large mass limit. Therefore, the strength of the exchange coupling is J + J for the conduction band and J − J for the valence band.

A. Effective Hamiltonian
For the Dirac Hamiltonian H Eq. (S22) with the exchange coupling H Eq. (S25), we perform the unitary transformation to diagonalize the Hamiltonian. We confirm that the exchange coupling commutes with the matrix β: From Eq. (S20), we obtain the effective Hamiltonian to order m −3 as where the operator F is Since the effective Hamiltonian is diagonal in the orbital space,β is regarded as an eigenvalue ±1 hereafter. As the present model is isotropic in the orbital space, we find Σ = σ, which corresponds to the spin of an electron. In the following, we consider the effective Hamiltonian to order m −2 .

S3. RESPONSE FUNCTIONS
We consider the response functions in the metallic state using functional derivatives. We first define the action S using the effective Hamiltonian with k = (k, iω n ). For clarity, we use the simplified notation where T is the temperature and ω n = (2n + 1)πT (n: integer) is the fermionic Matsubara frequency. The partition function Z is given by the path integral Using a functional derivative of the partition function, we can calculate the expectation value of an operatorX. The operatorX should appear in the Hamiltonian in the form where F X (r, t) is regarded as a generalized force that drives the quantity X. Then, the expectation value X satisfies Here, denotes the statistical average in the equilibrium. By further expanding the right-hand side, we can extract the effect of perturbations. The linear response of X to the generalized force F Y (r, t) becomes The result is equivalent to the Kubo formula. We find that the coefficient of the linear response becomes the correlation function This expression hold when the expectation value vanishes in the equilibrium. One can formally extend the expansion to higher orders of perturbations. The operatorX has the formX(r, t) =ψ(r, t)X(r, t)ψ(r, t), where X(r, t) corresponds to the matrix representation of the operator. Then, the calculation of the expectation value becomes the calculation of the connected diagrams of the Green's function. We define the Green's function with the unperturbed Hamiltonian as Here, 0 denotes the statistical average with the unperturbed Hamiltonian H 0 . The unperturbed Hamiltonian for the present system is and we use the empirical self-energy to describe diffusive transport with the lifetime τ . To simplify the notation, we introduce the vertex function Γ X as with q = (q, iΩ m ) and the bosonic Matsubara frequency Ω m = 2mπT (m: integer). Using Wick's theorem, we find that the correlation function becomes where tr stands for the trace of the matrix structure of the Hamiltonian; i.e., the spin matrices from the effective Hamiltonian for the present case.

A. Vertex functions
The analytic continuation requires the Matsubara frequency iΩ m to be replaced with ω + i0 + . In the current model, bosonic Matsubara frequencies correspond to frequencies of external fields. After the analytic continuation iΩ m → ω + i0 + , the vertex functions that we use below are Here, q 0 refers to the temporal component of q as a four-vector, i.e., q 0 = ω with q = (q, ω). Γ A 2 ab (k; q 1 , q 2 ) and Γ AM ab (k; q, Q) are defined from the second-order derivatives of −S: We hereafter omit the variables of Γ M a as it is constant.

B. Current operator
We derive the current operator here. Though we do not use it in the following calculations, it is worth knowing that the magnetic order affects the charge current. The functional derivative of the action gives the current operator J asĴ where J (r, t) becomes The result corresponds to the sum of Γ A , Γ A 2 , Γ AM , and Γ A 2 M . Alternatively, we can obtain J at A = 0 from the relation J = ie[r, H eff ] without inserting the electromagnetic potential.

S4. RESPONSE IN THE METALLIC STATE
In this section, we show the detailed calculations of various response functions that we presented in the main part. The density of states (DOS) repeatedly appears in the following calculations. For the unperturbed Hamiltonian H 0 , the DOS at the chemical potential µ is where we use the step function When we discuss the metallic state, it is convenient to measure the chemical potential from the band edge, which we denote as : Then, the carrier density is The DOS and the carrier density are related by A. Current response under an electric field

Linear response
Using the partition function, we can write the linear current response to the external electric field as j a (ω) = 1 iω We note that this equation is equivalent to the Kubo formula. The coefficient corresponds to the electric conductivity: To one-loop order, the conductivity becomes The result coincides with the one obtained from the Drude model.

Effect of a magnetic order without a uniform magnetization
Next, we consider the effect of the magnetic order, which modifies the conductivity. As we assume that the magnetic order does not have uniform magnetization and is characterized by finite wavevectors, the lowest-order correction by the magnetic order appears at second order. Therefore, we can write the correction to the conductivity by the magnetic order as where the coefficient is given by Using the vertex functions Eqs. (S40)-(S44), we obtain the following five contributions: Here we use the notation q = (0, iΩ m ) and Q = (Q, 0). After the momentum integrations, the summation of Matsubara frequencies, and the analytic continuation, we obtain

Effect of a uniform magnetization
When the local magnetic moments have a uniform magnetization M 0 = 0, we expect the anomalous Hall effect j a (ω) = σ AH ab (ω)E b (ω) with the anomalous Hall conductivity σ AH abc (ω) ∝ ε abc M c (0). One may calculate the anomalous Hall conductivity with the Green's function Eq. (S35) and the vertices Eqs. (S40)-(S44); however, one does not find the anomalous Hall conductivity at zero frequency σ AH ab (0) with the same procedure. We need a nonperturbative effect to the model, i.e., a correction to the unperturbed Hamiltonian. Here, we define the unperturbed Hamiltonian as For simplicity but without loss of generality, we assume that the uniform magnetization is oriented along the z axis. m z represents the spin polarization on the Fermi surface, arising from the exchange coupling −(J + J β )M · σ: The unperturbed Green's function uses H mag instead of H 0 : Now we can calculate the anomalous Hall conductivity similarly as Eq. (S57) with G 0 replaced with G 0,mag . For τ −1 , |m z | | F |, we obtain the anomalous Hall conductivity at low frequencies (|ω| Finite magnetization of the local magnetic moments forces spin polarization of the conduction electrons through the Zeeman coupling. We can understand the spin polarization of the conduction electrons as the spin-dependent Fermi energies F ± m z . The spin-polarized conduction electrons with the spin-orbital coupling inherent in the Dirac Hamiltonian lead to the finite anomalous Hall conductivity.

B. Magnetization by current
In the presence of an inversion-breaking magnetic order, the symmetry analysis allows finite uniform magnetization under an external electric field, i.e., electric current. As we have discussed, the lowest-order contributions appear at order M 2 , the induced uniform magnetization of Dirac electrons should have the form The coefficient λ takes the form from which we find the two contributions with q = (0, iΩ m ) and Q = (Q, 0). Using Wick's theorem, we can calculate the two contributions to obtain Here, we should expand the coefficient λ with respect to the wavevector Q and extract odd-order contributions to capture inversion breaking of the magnetic order. In the results above, we retain the terms to linear order in Q. The second term λ (2) is proportional to ω, so that it does not contribute to static uniform magnetization.

C. Current by an external magnetic field
An oscillating external magnetic field may induce electric current if the system breaks inversion. As the model that we consider here does not break inversion without a magnetic order, an inversion-breaking magnetic order is necessary for current response. The uniform current response should have the form In theory, it is convenient to consider Since the uniform magnetic field and the vector potential are related by B(ω) = iq × A(q, ω), we should expandκ with respect to q to findκ abcd (ω, q, Q) = κ abcd (ω, Q) · iε bãb qã.
Using the vertex functions Eqs. (S40)-(S44), we find the five contributions with q = (q, iΩ m ) and Q = (Q, 0). By evaluating the expressions, we can see that κ abcd vanishes at least to order QJ 2 n( )/m 2 ; therefore, we neglect the current directly induced by an oscillating external magnetic field. We cannot exclude the possibility of finite contributions at order n( )J 2 /m 4 or n( )J 2 /(m 3 (|µ| − |m|)) here. We emphasize that careful calculations are necessary, which should satisfy the gauge invariance Eq. (S81).

D. Current response by an oscillating magnetic order
As we have observed that the uniform current directly induced by an external oscillating magnetic field is negligible, we then investigate the current induced by an oscillation of the magnetic order, which has the form where the coefficient is given by γ has the two distinct contributions abc (ω 1 , ω 2 , Q) = kψ (k + Q 1 )Γ AM ab (k + q; −q, Q 1 )ψ(k + q) where we use the notations q = (0, iΩ m ), Q 1 = (Q, iΩ m1 ), and Q 2 = (−Q, iΩ m2 ) with the analytic continuations iΩ m → ω + i0 + , iΩ m1 → ω 1 + i0 + , and iΩ m2 → ω + i0 + . Using Wick's theorem, we can calculate the expressions to obtain γ (1) Like the calculation of the magnetization induced by an external electric field, we should expand the coefficient with respect to the wavevector Q and extract the odd-order contributions in Q to capture inversion breaking by the magnetic order. We can confirm that the uniform current vanishes in the zero-frequency limit (ω 1 , ω 2 → 0) as there must be no uniform current in the equilibrium.

S5. POLARIZATION
This section deals with an insulating case, where the chemical potential µ lies inside the gap (|µ| < |m|). We focus on the the orbital part of the polarization, which reflects the geometric properties of the wave function. Without the magnetic order or the exchange coupling, the electronic system itself preserves inversion and hence there is no polarization. An adiabatic insertion of the exchange coupling may develops finite polarization in the presence of an inversion-breaking magnetic order. We can unambiguously quantify the polarization by measuring from the inversionsymmetric state.
To calculate the electric polarization in an insulator, we follow the method by King-Smith and Vanderbilt [S3]. Suppose that the Hamiltonian H λ describes the electronic system, where the parameter λ (0 ≤ λ ≤ 1) continuously alters the potential. They showed that the change in the polarization per unit volume by an adiabatic change of the parameter λ is Here we write the Bloch wave function as ψ nk = u nk e ik·r with the band index n and the lattice-periodic part u nk , and the charge of an electron is −e (e > 0). The momentum integration is performed in the Brillouin zone and the band index is summed over the all occupied bands. The polarization is defined modulo the lattice period. This expression is concise, but it contains the wave function, which makes an analytic calculation difficult. There is an equivalent expression that uses the Green's function G λ = (ω − H λ ) −1 to calculate the orbital-part of the polarization [S4, S5]: Here, the trace tr includes the summation over all bands. In the following, we consider the massive Dirac Hamiltonian as the unperturbed electronic Hamiltonian and the exchange coupling with the magnetic order as the source of polarization. The Green's function becomes G J = (ω − H − H ) −1 . Note that the Dirac Hamiltonian is defined in the continuum without any lattice structure. The magnetic order with the wavevector Q introduces the periodicity to the system, which could be interpreted as the Brillouin zone. In considering the polarization in an insulating state, we can set the chemical potential µ = 0, so that it lies in the mass gap. When the exchange coupling is much smaller than the Dirac mass (|JM | |m|), we can treat the exchange coupling as a perturbation. Then, the expression of the polarization Eq. (S96) has the integration over the Brillouin zone and the summation of the band index, which can be replaced by the integration of the momentum for −∞ < k a < ∞. Henceforth, the trace means the matrix trace originated from the 4 × 4 Hamiltonian. As a result, the formula for the polarization becomes The polarization is measured from the state with J = 0, where the electronic system remains centrosymmetric and hence the polarization vanishes. The Green's function G J is no longer diagonal in the momentum space in the presence of a magnetic order and the exchange coupling. We now treat the exchange coupling perturbatively for (|JM | |m|), where the system remains insulating at µ = 0. Then, we expand the Green's function G J using the unperturbed one G 0 (k, ω) = [ω − H 0 (k)] −1 . When the magnetic order is written as the equation for the polarization becomes where G 0 = G 0 (k, ω) and M = M (Q).

A. Model calculation
The polarization induced by the exchange coupling is obtained from Eq. (S101) when the exchange coupling is treated perturbatively. The present model contains the two coupling constants for the exchange coupling. To evaluate the formula, we assume that the two coupling constants are proportional during the adiabatic insertion of the exchange coupling. To be more specific, we put J = χJ, where χ remains constant while J evolves. Then, Eq. (S101) leads to the uniform polarization per unit volume (S102) When we use the real-space form M (r), this uniform polarization results in where V is the volume of a unit cell and the spatial integration is performed over a unit cell determined by the magnetic order. We can rewrite the polarization in other forms. Using the relation where we neglect the boundary contribution in the latter, the uniform polarization per unit volume becomes The real-space expressions M (∇ · M ) and M × (∇ × M ) reminds us of the magnetoelectric effect, where magnetic structures that create a magnetic monopole or troidal moment yield finite effect. For example, the former resembles the diagonal magnetoelectric effect, where the polarization is parallel to the magnetization with the coefficient proportional to the charge of the magnetic monopole. We note that we can see the similarity of the magnetic monopole and troidal moment because we now focus on the uniform polarization. Finite polarization requires that the magnetic order realize finite Q Im[M * (Q · M )], which necessarily violates inversion. This factor has the form dr[M (∇ · M )] = − dr[(M · ∇)M ] in the real space and this is compatible with the Ginzburg-Landau theory obtained from the symmetry consideration [S6]. In addition, the strength of the the exchange coupling must be different for the two constituent orbitals of the Dirac electrons, implied by the product JJ . We suppose spiral spin orders with a single Q of the form where e 1 , e 2 , e 3 are orthonormal, and illustrate typical cases in Fig. S1. The condition Eq. (S102) states that there is finite polarization when the plane spanned by e 1 and e 2 (M 1 , M 2 = 0) contains the wavevector Q and that the induced polarization is on the same plane but perpendicular to Q.
The expression for the continuum model is to be contrasted with the atomic model [S7]. In the microscopic cluster model with two magnetic moments S i and S j , the polarization becomes P ∝ e ij × (S i × S j ), where the vector e ij connects the two magnetic moments. Therefore, a noncolinear magnet hosts finite electric polarization. While our expression may be apparently different, we find the same result for uniform polarization for the spin textures listed in Fig. S1. The microscopic model has an empty atomic site between the two magnetic moments and virtual transitions of an electron generate polarization. In our Dirac model, on the other hand, there are multiple bands from different elements that are extended in the bulk. As the conduction and valence bands are related with the nontrivial band topology, a certain magnetic order gives rise to finite polarization.

A. Review
To evaluate the effects that we have obtained with realistic material parameters, we rewrite the results by recovering the Planck constant h and the velocity v; the k · p Hamiltonian reads H 0 (k) = mv 2 β + vα · k. (S109) In the insulating state, the electric polarization generated by an inversion-breaking magnetic order is In the metallic state, we calculated the conductivity, a uniform magnetization, and the current induced by an oscillating magnetic order. The conductivity is with σ 0 (ω) = e 2 |n( F )|τ ω |m| , (S112) The magnetization is The current induced by an oscillating magnetic order is with (S118)

B. Estimates for magnetically-doped TIs
We adopt the values for the magnetically-doped TI Cr x (Bi 1−y Sb y ) 2−x Te 3 in Ref. [S8] presenting experiments of films, which uses the velocity v = 5.0 × 10 5 m/s, the mass m = −0.21m e converted from −300 meV with the electron mass m e = 9.109 kg, the coupling constants for the exchange coupling Jµ B = −5 meV and J µ B = 1 meV. In the following analysis for a metallic state, we set the Fermi energy | F | = 100 meV, which is smaller than half the mass gap |mv 2 | = 300 meV. Then, Eq. (S53) gives the carrier density |n| ≈ 1.4 × 10 19 cm −3 . We estimate the lifetime τ from the longitudinal conductivity σ 0 ≈ 100 Ω −1 cm −1 using the relation σ 0 = e 2 |n|τ /|m| to obtain τ ≈ 5.4 × 10 −15 s, which corresponds to an energy scale /τ ≈ 120 meV. The RKKY interaction presumably plays a dominant role in forming a magnetic order in the metallic state with the wavenumber Q being twice the Fermi wavenumber 2k F . We use Q = 2k F ≈ 1.5 × 10 9 m −1 at | F | = 100 meV. We suppose that each Cr atom has a magnetic moment M = 3µ B from the experimental observation.
We evaluate our results using the material parameters above. In the metallic state, we first calculate the corrections to the conductivity. The isotropic term similar to the magnetoresistance is which reduces the conductivity about a few percent. The reduction depends on the bands labeled byβ, reflecting the strength of the exchange coupling. The magnitude of the η term, which resembles the anisotropic magnetoresistance and the planar Hall effect, is . (S120) It gives rise to the anisotropy in the conductivity.
In the insulating state, the wavenumber of a magnetic order may be different from that in the metallic state, but we use here the same value for the estimate. Then, the electronic polarization induced by a spiral magnetic order becomes ∆P ≈ − eJJ 6π 2 v 3 |m| QM 2 ≈ 1.8 µC/m 2 . (S129) In concluding the estimates, we comment on the surface contributions of magnetic topological materials. In the metallic state, the surface contributions if present should be much smaller than the bulk contributions because of their small volume proportions. In addition, the surface effect is insensitive to the sample thickness, which we could distinguish in experiments. On the other hand, when the bulk is in the insulating state, the surface of a magnetic topological material can be either metallic or insulating depending the Fermi energy. Whether or not the surface contributes to the electric polarization, the electric polarization of the bulk depends on the sample volume or the thickness, whereas the surface contribution does not, which suggests an experimental identification.
We calculate the uniform polarization induced by the exchange coupling using Eq. (S101) to obtain The expressions have the similar structure as the isotropic case in Sec. S5 A, but they reflect the stacking structure of the material.

C. Effective Hamiltonian
By applying the formal expression Eq. (S20), we can obtain the effective Hamiltonian valid for large m: where the matrix Σ is Σ = (−σ x β, −σ y β, σ z ). (S138) We also note α = (ẑ × σ)τ z +ẑτ y . The difference from the previous model Eq. (S26) arise because of the different commutation relation between α and σ, which reflects the planar anisotropy of the present model.

D. Current operator
We can calculate the current response from the oscillation of the magnetic order by following the definition Eq. (S48). As the present model has a complication concerning the anisotropy, it would be useful to write down the general expression first. The current operator J obtained from the Hamiltonian (S20) is The first two terms correspond to the conventional current operator and the contribution from the spin-orbit coupling, respectively. By inserting the explicit form of H , we obtain the current operator for the present model: Now we notice that the current operator couples to the magnetic order M . We can further calculate the Pauli matrices and the differential operators, but we conclude the calculation here; the focus of the section is to see that the electric current couples to the magnetic order.