Axion mass in antiferromagnetic insulators

We calculate axion potential in antiferromagnetic insulators from path integral. It is derived from the effective potential for the order parameter of the antiferromagnetic phase in insulators. Static and dynamical axions are defined consistently from the potential. Consequently antiferromagnetic/paramagnetic and topological/normal orders are classified. The dynamical axion is predicted in all phases and its mass turns out to be various values up to eV scale.


I. INTRODUCTION
Axion is a hypothetical elementary particle that is considered to solve the strong CP problem in QCD [1]. It is originally massless due to shift symmetry of the Lagrangian, and it becomes massive after the QCD transition via a nonperturbative effect. Other than the QCD axion, axion particles are known to emerge from superstring theory. Motivated partly by that, axion-like particles (ALPs) are studied in phenomenological, astrophysical, and cosmological contexts. The mass of ALPs is usually considered as a free parameter, leading to more rich phenomenology compared to the QCD axion (see, e.g., Ref. [2] for review).
Besides in the field of particle physics, axion is discussed in condensed matter physics. For example, the magnetoelectric response due to E · B coupling (E and B are electric and magnetic fields, respectively) is already known, and the E · B coupling is the interaction term of axion with the gauge field. Therefore, such magnetoelectric response may be able to understood in terms of "axion." There have been several studies in that context so far. Different from axion in elementary particle physics, there are two types of "axions" in condensed matter physics: static axion and dynamical axion. (Hereafter we use "axion" as an axion in condensed matter physics unless otherwise mentioned.) In terms of Lagrangian, they are written as where α = e 2 /(4π) is the fine-structure constant and we have denoted θ 0 as the static axion (constant) and a as the dynamical axion. The former is pointed out by Refs. [3][4][5][6][7][8] as the E · B coupling. Nonzero θ 0 breaks both time-reversal T and parity P. However, when θ 0 = π those symmetries are intact. This nontrivial state is known as an topological state. Ref. [9], on the other hand, studied the latter to predict that the dynamical axion causes the total reflection of light irradiated to the topological insulators under an antiferromagnetic (AFM) order. They solved equations of motion of photon and the dynamical axion that interact with each other via the E · B interaction term and found an energy gap emerged in the topological insulators. In the similar context, Ref. [10] shows that an electric field is expected to be converted to a magnetic field due to the dynamical axion when the electric field exceeds certain critical value. (See also Refs. [11][12][13][14][15].) Such recent developments in both particle physics and condensed matter physics have motivated particle axion and ALPs search by using axions in insulators [16][17][18][19].
In this paper we calculate the effective potential for the order parameter of the AFM phase from the partition function of the system defined by the path integral. Using the effective potential, axion potential is derived. It is found that the effective potential and axion potential are intuitive and useful in order to distinguish the AFM/paramagnetic (PM) order and the topological/normal phase depending on the model parameters. We clarify the static and dynamical axions in the AFM insulators and how they are related to the AFM order and the topological state. Consequently, the mass of dynamical axion is determined. While it is estimated to be O(meV) in Ref. [9], the axion mass turns out to be larger up to O(eV) or a much suppressed value depending on the property of insulators.
This paper is organized as follows. In Sec. II, the effective potential for order parameter of the AFM phase is derived. Using the result, the AFM/PM phases are determined. The effective potential is interpreted in terms of axion field and axion potential is derived in Sec. III. Sec. IV contains conclusions and discussion.

II. POTENTIAL OF ANTIFERROMAGNETIC ORDER PARAMETER
The topological state in insulators is consider to be realized by a (quasi) gapless state via the spin-orbit coupling. For modeling such a state, a minimum setup is to consider 4-by-4 matrix, i.e., two energy states (or orbits or sublattices) with up and down spins. Such a model is often described by the gamma matrices (even in lower or higher dimension space [3]). In the momentum space, it is given by where k = (k x , k y , k z ) is the wavenumber and Γ a are the gamma matrices satisfying The representation of the gamma matrices depends on the basis. In Refs. [9,20], for example, the basis of the antibonding and bonding are chosen, i.e., (|P 1 + z , σ , where σ j and s j (j = x, y, z) are the Pauli matrices. The coefficients d a , on the other hand, are given by theoretical model and/or the first-principles computation. Among models, the so-called Dirac model is a low energy effective model to describe the topological insulators expanding a Hamiltonian around a Dirac point. In the Dirac model, d a are given as where k 2 = |k| 2 and the gamma matrices are given in Eq. (5). 0 , A, B, m, and m 5 have a mass dimension of one. ( 0 is not important in our later discussion, which is explained in Appendix C.) The Dirac model can be derived from a certain model, such as an extended Fu-Kane-Mele-Hubbard model on a diamond lattice at halffilling, discussed in Refs. [19,21], or the first-principles calculation of the layered, stoichiometric crystals Sb 2 Te 3 , Sb 2 Se 3 , Bi 2 Te 3 and Bi 2 Se 3 [20]. In general, the coefficient A corresponds to the spin-orbit coupling. 1 If it is zero, then θ is exactly zero and no topological phase appears, which will be easily understood from Eq. (28). The coefficient m and B describe the energy gap. They originate in, for instance, the spin-orbit coupling [20], or hopping strength anisotropy due to the lattice distortion in the diamond lattice [21]. m (and B) effectively determines whether the insulator is topological or normal. Namely, m/B < 0 (m/B > 0) corresponds to the topological (normal) state for B = 0. When B = 0, negative (positive) m gives the topological (normal) state. Finally m 5 is an order parameter of the AFM phase and plays an important role for describing axion in insulators, which will be discussed in detail soon. While in the numerical study we adopt the Dirac model, the following analytical calculation is applied for generic model described by Eq. (2). Besides, in Appendix A we give another numerical results in the effective model for 3D topological insulators given in Eq. (A1) [9].
1 To be precise, the spin-orbit interaction itself is local and it is constant in momentum space. The momentum dependence comes from nonlocal interaction coupled to spin-orbit coupling.
The m 5 Γ 5 term appears due to the Hubbard-Stratonovich transformation or the mean field approximation of the Hubbard interaction term, which causes the AFM order written in a sublattice basis ( |A, σ , |B, σ ), where n Aσ = ψ † Aσ ψ Aσ and n Bσ = ψ † Bσ ψ Bσ are number density of electron with spin σ =↑, ↓ at sublattice A and B, respectively. U is a parameter with mass dimension one, V is the volume of the material and N is the number of sublattice A (B) in the material. This interaction term involves four electrons and it is difficult to analyze. However, it can be rewritten by the Hubbard-Stratonovich transformation with introducing a scalar field φ in the path integral. (See Appendix B for details.) As a result, the interaction term becomes where n I− = n I↑ − n I↓ (I =A, B), and Here the momentum integral is executed in the first Brillouin zone unless otherwise mentioned, which is applied to other momentum integral in the later discussion. The term proportional to n A− − n B− gives m 5 Γ 5 term in the sublattice basis. (See Appendix D for the explicit expression of Γ a .) Here we have renamed m 5 as When m 5 = 0 the Hamiltonian is T and P conserving while a finite value of m 5 breaks both. In addition, it should be noticed that the mass term for φ emerges, which must be included in the effective Lagrangian for φ.
For later convenience we separate the Hamiltonian as where Here d a 0 is defined by From the Hamiltonian, the partition function is defined by functional integral of four-component electron field ψ and φ, where Here H is the Hamiltonian in the coordinate space. Now we integrate out ψ in Eq. (15) to derive the effective Lagrangian for φ. It is straightforward to get Here "Tr" includes that trace over matrices and spacetime integral. To further rewrite the above expression, let us define Green's function as where H 0 is H 0 (k) written in the coordinate space. Then, separating H as H = H 0 + δH and using we get The effective potential for φ is derived as polynomial expansion from Eq. (21). Using the Fourier transformation of the propagator, where k · x = k 0 x 0 − k · x and |d 0 | 2 ≡ 4 a=1 d a 0 d a 0 , it is obtained as (see Appendix C for detail) A similar calculation is done in Refs. [15,18,22]. However, the expansion is truncated at n = 2. On top of that, the mass term, the second term in Eq. (24) is missing. 2 In fact, both the negative term derived from the last term in Eq. (21) and the mass term are crucial to understand the AFM/FM phase and the topological/normal state, which will be shown below. Fig. 1 shows V φ (φ) as function of φ. In the calculation we introduce dimensionless wavenumberk ≡ ka where a is the lattice size (typicallyÅ) and take A/a = 1 eV and B/a 2 = 0.5 eV 3 and the integral is executed in [−1, 1] region. This is because the Dirac model is low energy effective model and the (normalized) momentum cannot take as large as π. (Instead if the effective model for 3D topological insulators given in Ref. [9] is used, then the momentum integral is done in [−π, π]. See Appendix A for the results.) It is seen that the curvature at φ = 0 becomes negative when U gets larger. This situation can be understood by expanding V φ around φ = 0, Here we have used Eq. (9). It is seen that the coefficient of φ 2 term becomes negative when U > U crit where U crit satisfies Namely, for U < U crit the minimum of the potential located at φ = 0, corresponding to the PM phase. When U > U crit , on the other hand, the phase transition occurs to get nonzero φ 0 , which satisfies and the AFM order appears. In the next section, we show how the situation is described in terms of the axion field.

III. AXION POTENTIAL
So far we have derived the effective potential for φ. It can be described by axion field because φ is related to the axion field θ. To be concrete, θ is given by d a [9], 2 There are other typos in Ref. [22]. We thank A. Sekine for confirming this. 3 A larger value of B is obtained in first-principles calculation [20].
This may indicate it is difficult to have the AFM order in the materials such as Fe-doped Bi 2 Te 3 (see also later discussion). In our analysis, we choose O(1) eV for A and B. where |d| 2 = 5 a=1 d a d a and ijkl is Levi-Civita symbol with i, j, k, l being 1, 2, 3, 5. It is checked that θ = 0 or ±π when d 5 = φ = 0 and the value depends on d 4 . (We note that θ = +π and −π are identical and hereafter we focus on θ = π unless otherwise stated. See also Appendix D for another derivation of θ using chiral anomaly. However, it is not equivalent to the above expression, except for a specific case.) For instance, in the Dirac model given in Eq. (6), m/B > 0 (< 0) gives θ = 0 (π) when d 5 = 0 [23]. We have confirmed that the results are consistent with Ref. [23] (see Figs. 4 and 5 in Appendix A). Using Eq. (28), we can derive φ as function of θ. Let us write Then, the potential as function of θ can be obtained by Fig. 2 shows V θ as function of θ/(2π) in the Dirac model. It is seen the axion potential nontrivially depends on the sign of m. When U is smaller than U crit , the potential minimum is located at θ = 0 (π) for m > 0 (m < 0). V θ is multivalued function of θ since φ is multivalued function of θ for m > 0. For m < 0, on the other hand, θ is monovalent function, i.e., θ → ±π for φ → 0±. (See Fig. 4 in Appendix A.) In terms of potential minimum of V φ , both θ = 0 and π correspond to φ (or m 5 ) = 0, i.e., the PM phase. In terms of θ, on the other hand, they are different topologically. Namely, θ = 0 (π) is classified as the normal insulator (topological insulator). When U gets larger than U crit , nonzero φ 0 becomes a minimum of the potential and insulator is in the AFM order. In terms of axion, both θ = 0 and π becomes unstable and a new stable point appears. It should be noticed that the effective potential given in Eq. (24) is crucial to derive this result, which consists of the negative potential term derived from integrating out the electron field and the positive mass term that comes consequently from the Hubbard-Stratonovich transformation of the Hubbard interaction term (7). Now the definition of static and dynamical axions are clear. The static axion θ 0 is defined as a stationary point of V θ (θ), and the dynamical axion a as a fluctuation around θ 0 . Denoting φ 0 as corresponding value of θ 0 , we expand φ and the potential V φ , where J is called stiffness and obtained as (see Appendix C for derivation) Using Eqs. (31) and (32), the action for the dynamical axion a is obtained, Therefore, the axion mass is obtained by for U < U crit and for U > U crit . Eqs. (33) and (36) are consistent with the results in Ref. [9], except for a factor of 4 larger compared to ones in the literature. In the literature, however, how the value of φ 0 is determined is not clearly stated. Rather it is treated as another parameter that is independent of the model parameter in d a . By using the effective potential, on the contrary, φ 0 is determined uniquely as a stationary point of the potential when the model parameters are fixed.
To demonstrate the result, we plot the axion mass on phase plane in Fig. 3. There are phases classified by the AFM/PM and the topological/normal orders, and they are determined by U and m for the given parameters A and B of the Dirac model. The AFM (PM) phase corresponds to φ 0 =nonzero (0), and the topological (normal) phase corresponds to m < 0 (m > 0). In the PM phase, θ 0 = 0 (π) in the normal (topological) phase meanwhile θ 0 takes a value of 0 < θ 0 < π in the AFM phase depending on φ 0 . It is found that the axion mass is typically O(1) eV in the four phases. (The same result is obtained in the effective model for 3D topological insulators. See Fig. 8 in Appendix A.) This result seems to be different from the value argued in the literature [9] where the O(meV) axion in the topological AFM insulators is discussed. Let us elaborate on this issue closely.
First of all, it is clear that the existence of dynamical axion does not depend on whether insulator is topological or not. This point is mentioned in Ref. [10]. Fig. 3 shows that the axion mass can be much suppressed near the critical point of the AFM phase in both the topological and normal insulators. In Ref. [9], on the other hand, m 5 (= φ 0 ) = 1 meV is taken, which leads to meV scale axion mass. At a glance this looks possible; a small value of U can give a suppressed m 5 since m 5 ∝ U M − z in mean field approximation (M − z in Ref. [9] corresponds to m z in Appendix B 3). However, it is no longer the AFM order if U < U crit . Therefore, suppressed U cannot be naively considered in this simple setting. This relates to an intrinsic issue that Bi 2 Se 3 does not have an AFM phase. A possible way out suggested in Ref. [9] is doping Fe to create the AFM order. In this way, the AFM order could be realized. However, doping Fe may not be successful due to the first-principles calculation [24]. A current candidate is Mn-doped topological insulator, Mn 2 Bi 2 Te 5 suggested in Ref. [23]. The first-principles calculation done by Ref. [25] shows that Mn 2 Bi 2 Te 5 has rich magnetic topological quantum states, including the AFM, FM and other magnetic states. The energy gap of the AFM phase indicates m 5 ∼ O(10 -100 meV). The total energy of each magnetic states are also found to be the same order, which might indicate that controlling the realization of the magnetic states would be challenging. Therefore, the next step is to confirm these magnetic states experimentally. 4 In a nutshell, the axion mass is roughly O(eV) but it can be suppressed near the critical boundary between the AFM and PM orders. Such dynamical axion with various mass range might be utilized in the future particle axion or ALPs search. On the other hand, in order to precisely predict the axion mass in such specific circumstances, the first-principles calculation to determine the 4 Another interesting candidate is MnBi 2 Te 4 -family materials. Although they are not supposed to have dynamical axion, they are found to have controllable AFM and ferromagnetic orders. Therefore it is worth studying their properties for the future application. See, for example, Refs. [26][27][28][29][30][31].
parameters is required. Besides, we need experimental evidence of having the AFM order in Fe-or Mn-doped Bismuth Selenide or Bismuth Telluride. It is worth noticing that the dynamical axion should exist in the PM phase. 5 On the other hand, there would be other degree of freedom in the PM phase. For example, φ is introduced in the circumstance where the AFM order is assumed (see, e.g., Eq. (B18)). If this assumption is relaxed, another scalar field should appear. Therefore, it would be nontrivial to find out dynamical axion mode under such contaminated circumstances. If the possible obstacles are overcome, however, the dynamical axion in the PM phase might be another tool for the future particle axion search.

IV. CONCLUSIONS AND DISCUSSION
We have investigated axion in AFM insulators. In the analysis, the effective potential for the order parameter of the AFM phase is derived from the path integral. From the effective potential, the order parameter of the AFM is determined as the stationary point of the potential, which distinguishes the AFM or PM order. Axion potential is obtained from the effective potential. It is shown that the axion field value θ = 0 or π is the minimum of the axion potential in the PM order, meanwhile in the AFM order different values of θ are obtained as the minimum of the potential. The field value of the stationary point of axion is the static axion, and the dynamical axion is defined as a fluctuation around stationary value. The dynamical axion is predicted in four types of phases distinguished by the AFM/PM and the topological/normal states. The mass of axion turns out to be typically O(eV). On the other hand, it can be much suppressed near the phase boundary between the AFM and PM orders. In either case, it is crucial to determine the value of the order parameter in a possible candidate, e.g., Mn 2 Bi 2 Te 5 , to predict the axion mass.
Another possible direction to investigate the static/dynamical axion is to utilize the magnetoelectric effect. The magnetoelectric couplings in generic insulators are obtained by the first-principles calculation by Ref. [32]. For example, θ = 1.3 × 10 −3 and 1.07π are obtained for Cr 2 O 3 and Bi 2 Se 3 . In addition, if Bi is magnetized, the value of θ is shifted to θ = π ± 0.55. On the experimental side, quantized Faraday and Kerr rotation due to Bi 2 Se 3 film is observed precisely [33]. Recently, anisotropic topological magnetoelectric effect is observed in axion insulators [34]. Those are examples of "axion electrodynamics" as discussed in Refs. [9,10]. Such methods using optical responses would make it possible to directly probe the order parameter of the magnetic states and the properties of the static/dynamical axion. The formalism given here would be a help for finding and modeling realistic materials for the future particle axion or ALPs search. Searching a wide range of materials including ones that have no topological states or have ferromagnetic/paramagnetic phases might open a new possibility for particle axion search.

ACKNOWLEDGMENTS
We are grateful to Makoto Naka for valuable discussions in the early stage of this project and careful read-ing of the manuscript. We also thank Akihiko Sekine for useful discussions. This work was supported by JSPS KAKENHI Grant Numbers JP17K14278, JP17H02875, JP18H05542, and JP20H01894.

Appendix A: Effective model for 3D topological insulators
We give numerical results in the effective model for 3D topological insulators studied in Ref. [9]. The model is parametrized as, where the gamma matrices are the same as Eq.  In this section, we discuss a possible origin of δH term given in Eq. (13). For pedagogical purpose, we provide three types of derivations; derivations in the continuum space (Appendix B 1), in a discrete space (Appendix B 2) and by the mean field approximation (Appendix B 3).

In continuum space
To begin with, we consider one of sublattices since it is enough to get the basic idea of emerging a scalar degree of freedom that corresponds to m 5 (= φ). The interaction term to cause the AFM order is given by following Hamiltonian: 6 where arrow stands for the spin, V is the volume of the material, N is the number of site, and It is convenient to define the following quantities: and take since the total number of electrons are unchanged. (This corresponds to Eq. (B46) in the mean field approximation.) Therefore, only n − is a dynamical degree of freedom, and H int becomes 7 6 This term corresponds to the continuum limit of Eq. (B23) by making a replacement, 1 For later discussion we define the Fourier expansion of ψ σ and n σ as Consequently, H int in the momentum space is given by In terms of Euclidean time τ = it, action of the interaction term is written by Now we execute the Hubbard-Stratonovich transformation. It is also used in Ref. [35], which studies topological superconductors and superfluids. (See Ref. [36] for the renormalization group approach.) It is effectively the same as the inverse procedure of integrating out the heavy scalar particles to give the effective Lagrangian. (See Appendix E for details.) Introducingφ(k) and using the following identity, the Euclidean action gives Here we have omitted irrelevant constant term in the last step and we have used The second term on the rhs of Eq. (B13) is given in the momentum space as Now we apply the above result to two sublattices, A and B. The interaction term to start with is where N is the number of the sublattice A (B). As in the previous discussion, we define n I± (x) (I=A, B), Since we are interested in the AFM order, n A− and n B− are not independent and related by 8 Therefore, we obtain and Here the summation over spin indices are implicit. This term corresponds to d 5 Γ 5 (= φΓ 5 ) term in Eq. (2) after changing the basis from (|P 1 + z , σ , |P 2 − z , σ ) to (|A, σ , |B, σ ) where σ = ↑, ↓. In this basis, Γ 5 = diag{1, −1, −1, 1} and Eq. (B19) becomes where c k = (c Ak↑ , c Ak↓ , c Bk↑ , c Bk↓ ) and 9

In discrete space
As in the previous subsection, we first consider one of the sublattices and extend the result to the other later. In the discrete space, we use annihilation operator c iσ where i is a label of site instead of the wavefunction ψ σ (x). The interaction Hamiltonian is where While defining n i± as only n i− is dynamical since n i+ =const. as in Eq. (B4). Consequently, the interaction Hamiltonian is given by It is useful to give Fourier transformation of c iσ and n iσ , which leads to the Hamiltonian in the momentum space, As in the continuum case, the Hubbard-Stratonovich transformation is done by using Eq.(B12). Then the Euclidean action of the interaction term is given by Here we have defined the Fourier expansion of φ i as and the second term on the rhs of Eq. (B32) is written in the momentum space as Extension to the other sublattice is trivial. As in the previous subsection, let us call the two sublattice as A and B. Then the interaction Hamiltonian is Repeating the discussion so far, the Euclidean action is obtained as The interaction part in momentum space is Therefore we have obtained the same result as Eq.(B20) butφ 2 (0) is replaced byφ 0 . Finally Eq. (B37) becomes where c i = (c Ai↑ , c Ai↓ , c Bi↑ , c Bi↓ ). When we make a continuum limit for space, i.e., We show how to derive the effective potential (24) and stiffness J (33) from Tr(G 0 δH) n term in Eq. (21). It corresponds to calculation of one-loop diagram with external scalar fields those three-momenta are zero. It is similar to the Coleman-Weinberg potential [37], but the the result turns out to be non-logarithmic function. In the following calculation, the continuum limit in the coordinate and momentum space is taken. We define the Fourier expansion of G 0 and φ as First we show how to derive the effective potential. Note that the trace vanishes when n is odd. Thus, let us begin with n = 2: Here trace part gives where index a is contracted. To give the effective potential we can take k = 0. As a result, it is given in a simple form, In the similar manner, higher polynomials are calculated as where n ≥ 2. To get the φ n terms, we can take k i = 0. Then the trace part is where n = 2m (m = 1, 2, 3, · · · ). Using this formula, the last term in Eq. (21) gives Here we assume that q 0 integral is from −∞ to +∞, and we use the Wick rotation of the the integral path as Finally the effective potential for φ (= m 5 ) is obtained as The stiffness J is obtained as an coefficient of the time-derivative term of the dynamical field. It is determined after the potential minimum is found. Therefore, instead of using Eqs . (12) and (13), we take where now d 5 = φ 0 (and the others are the same) and δφ = φ − φ 0 . The propagator is given similarly to Eq. (23) as where |d| 2 ≡ 4 a=1 d a d a + φ 2 0 . To derive the time-derivative kinetic term, the quadratic term of δH should be expanded in terms of k 0 (but k = 0 can be taken). The trace part is given by After expanding with respect to k 0 , we get the kinetic term where δφ is Fourier transformation of δφ and J = d 3 q (2π) 3 |d 0 | 2 4(|d 0 | 2 + φ 2 0 ) 5/2 .

(C16)
Appendix D: θ term as chiral anomaly In a special case of the Dirac model, θ can be derived as the path integral measure due to the chiral rotation of the electron field ψ [19,21]. In this section we apply the technique to the present model and clarify an issue claimed in Ref. [21] that θ derived in this technique can not be applied in some circumstances.
In the ( |A, σ , |B, σ ) (σ =↑, ↓) basis, the Gamma matrix is given as To be concrete, Γ a corresponds toα µ rotated by U 3 in Ref. [19], except for an extra sign inα 5 . Or if one makes another transformation, the Dirac α matrices in the Dirac representation are obtained, Here the Dirac α matrices corresponds to (α, β) = (Γ 1 ,Γ 2 ,Γ 3 ,Γ 4 ). Therefore, where we have changed variable k as Ak = k , and introduced the Dirac gamma matrix γ µ = (β, βα) and γ 5 = iγ 0 γ 1 γ 2 γ 3 . Consequently, this Hamiltonian leads to the following action, whereψ = ψ † γ 0 and we have added vector potential A µ of U (1) em . Then chiral rotation of ψ gives the chiral anomaly term in the QED [38,39], whereF µν = 1 2 µνρσ F ρσ ( 0123 = +1) and Θ = π 2 [1 − sgn(m)]sgn(φ) + tan −1 φ m . (D11) Here we have set the domain of tan −1 as [−π, π]. This is consistent with Ref. [21], except that the sign of m and we have used single degree of freedom for electron field ψ. In the model considered in Refs. [19,21], there are three Dirac points and they expand the Hamiltonian at each Dirac points. That is why they have three degree of freedom for electron field. (The sign of M f e iκ f γ5 in Eq. (30) of Ref. [21] is supposed to be a plus. Comparing to the result in Ref. [19] (v1 on arXiv), their result is a factor of two smaller. In addition, it is claimed that first term is 1/2, instead of π/2, in the topological phase. 10 ) The authors of Ref. [21] claim that the expression of θ cannot be applied when U m f /δt 1.
(U m f /δt corresponds to φ/m in our notation.) This is misleading because the derivation by chiral rotation in the path integral is exact and there is no approximation. However, it should be noted that in this calculation the cutoff of the momentum integral is taken to be infinity. To put it more correctly, therefore, this result is reliable when U m f /δt is smaller than the cutoff momentum. Although the Dirac model is the low energy effective model, if the cutoff momentum is taken to be infinity, then we get from Eq. (28) which is the same result with Eq. (D11).