Nonreciprocal current from electron interactions in noncentrosymmetric crystals: roles of time reversal symmetry and dissipation

TL;DR

This paper demonstrates that electron interactions can induce nonreciprocal current-voltage characteristics in noncentrosymmetric crystals without breaking time-reversal symmetry, highlighting the roles of dissipation and interactions.

Contribution

It reveals that nonreciprocal responses can occur due to electron correlations in inversion-breaking systems without requiring time-reversal symmetry breaking.

Findings

Electron interactions can induce nonreciprocal $I$-$V$ in noncentrosymmetric systems.

Derived a general formula for nonreciprocal response using nonequilibrium Green's functions.

Applied the theory to the Rice-Mele model and discussed potential materials.

Abstract

In noncentrosymmetric crystals with broken inversion symmetry $\mathcal{I}$, the $I-V$ ($I$: current, $V$: voltage) characteristic is generally expected to depend on the direction of $I$, which is known as nonreciprocal response and, for example, found in p-n junction. However, it is a highly nontrivial issue in translationally invariant systems since the time-reversal symmetry ($\mathcal{T}$) plays an essential role, where the two states at crystal momenta $k$ and $-k$ are connected in the band structure. Therefore, it has been considered that the external magnetic field ($B$) or the magnetic order which breaks the $\mathcal{T}$-symmetry is necessary to realize the nonreciprocal $I-V$ characteristics, i.e., magnetochiral anisotropy. Here we theoretically show that the electron correlation in $\mathcal{I}$-broken multi-band systems can induce nonreciprocal $I-V$ characteristics {\it without} $\mathcal{T}$-breaking. An analog of Onsager's relation shows that nonreciprocal current response without $\mathcal{T}$-breaking generally requires two effects: dissipation and interactions. By using nonequilibrium Green's functions, we derive general formula of the nonreciprocal response for two-band systems with onsite interaction. The formula is applied to Rice-Mele model, a representative 1D model with inversion breaking, and some candidate materials are discussed. This finding offers a coherent understanding of the origin of nonreciprocal $I-V$ characteristics, and will pave a way to design it.

Figures

Equations

$$J$$

$${\sigma }_{ij}$$

$$J_i({\omega }_1+{\omega }_2)$$

$${\sigma }_{ijj}(i{\omega }_{n_1},i{\omega }_{n_2})$$

$$G(\omega ,k)$$

$$\Sigma (\omega ,k)$$

$${\left(\begin{array}{cc}{G}_0^R& G_0^K\\ 0& G_0^A\end{array}\right)}^{-1}$$

$$G_E$$

$$\displaystyle\hskip 50.00008pt-(\partial_{k}G_{0}^{-1})G_{0}(\partial_{\omega}G_{0}^{-1})\big{)}\Big{]}G_{0}.$$

$$J_{E^2}=-i\int d\omega dk\text{tr}[v(k)G_{E^2}^{<}(\omega ,k)].$$

$$G_{E^2}$$

$$+\frac{1}{4}{G}_0[({\partial }_{\omega }^2{G}_0^{-1})({\partial }_k^2{G}_0)+({\partial }_k^2{G}_0^{-1})({\partial }_{\omega }^2{G}_0)].$$

$$G_0(\omega ,k)$$

$$G_E(\omega ,k)$$

$$v(k)$$

$$\text{tr}[v(k)G_0(\omega ,k)({\partial }_{\omega }{G}_0^{-1}(\omega ,k))({\partial }_k{G}_E(\omega ,k))]$$

$$=-\text{tr}[v^T(-k)G_0^T(\omega ,-k)({\partial }_{\omega }{G}_0^{-1,T}(\omega ,-k))({\partial }_k{G}_E^T(\omega ,-k))]$$

$$=-\text{tr}[v(-k)G_0(\omega ,-k)({\partial }_{\omega }{G}_0^{-1}(\omega ,-k))({\partial }_k{G}_E(\omega ,-k))],$$

$$\text{tr}[v(\theta )G_{E^2}^{<}(\omega ,\theta )]=-\text{tr}[v(-\theta )G_{E^2}^{<}(\omega ,-\theta )],$$

$$G_{E,11}^{<}$$

$$G_{E,12}^{<}$$

$${\Psi }_{1,k}$$

$$H_{\text{int}}$$

$${\Sigma }_{E,11}^R(k)$$

$$\times [u_{k^{\prime }}{v}_{k^{\prime }}{G}_{E,12}^{<}(k^{\prime })+u_{k^{\prime }}^{\ast }{v}_{k^{\prime }}^{\ast }{G}_{E,21}^{<}(k^{\prime })],$$

$${\Sigma }_{E,11}^R(k)$$

$$J$$

$$\delta J$$

$$=\frac{2aU\tau {{\partial }_k(|u_k{|}^2-|v_k{|}^2)|}_{k=k_F}}{\pi |v_{11,k_F}|E_{g,k_F}}\text{Im}\left[u_{k_F}{v}_{k_F}{v}_{12,k_F}\right]E^2,$$

$$\gamma \equiv \frac{\delta J}{J}\simeq \frac{U}{E_{g,k_F}}\frac{eEa}{W},$$

$$H$$

$$+\Delta \sum _i(-1{)}^i{c}_i^{†}{c}_i.$$

$$H$$

$${ϵ}_{1,k}$$

$$H_{\text{int}}$$

Full text

Noncentrosymmetric crystals exhibit a variety of interesting physical phenomena. These include ferroelectricity Resta , photovoltaic effect (shift current) Grinberg ; Nie ; Shi ; deQuilettes ; Kraut ; Sipe ; Young-Rappe ; Cook17 ; Morimoto-shift16 , and second harmonic generation Boyd ; Bloembergen ; Wu16 . Among them, nonreciprocal dc current response in inversion broken systems has been attracting a keen attention in condensed matter physics. Nonreciprocity (or rectifying effect) is a current response where the $I-V$ characteristic differs when current flows toward left and when it flows toward right (i.e., $I(V)\neq-I(-V)$). The nonreciprocal current response is important both for fundamental physics of inversion broken materials and also for applications such as diode. Conventional example of nonreciprocity is a p-n junction, in which the direction of the current changes the thickness of depletion layer, and hence, the resistivity. Nonlinear current response has been intensely studied in a mesoscopic setup ChristenButtiker96 ; Song99 ; SanchezButtiker04 . In contrast to such artificial heterostructures, nonreciprocity in crystals is a more nontrivial issue. Current responses in crystals are governed by Bloch electrons with good momentum $k$ and their band structure. In the presence of time-reversal symmetry (TRS), the band structure satisfies the relationship $\epsilon_{k\sigma}=\epsilon_{-k{\bar{\sigma}}}$ ($\sigma$ represents the spin and ${\bar{\sigma}}$ the opposite spin to $\sigma$), which indicates that no nonreciprocity appears for noninteracting electrons in the Boltzmann charge transport picture as illustrated in Fig. 1(a). Specifically, the applied electric field causes a shift of the Bloch electrons in the momentum space. The symmetry in the band structure due to the TRS results in symmetric shifts with respect to the direction of the applied electric field $E$, and the conductivity does not depend on the direction of $E$. There are two ways to break TRS: (i) introducing time-reversal breaking term to a microscopic Hamiltonian and (ii) introducing irreversibility at the macroscopic level. The former microscopic TR breaking is achieved with application of an external magnetic field $B$ or introducing magnetic order. Nonreciprocal current response in the presence of magnetic field is known as magnetochiral anisotropy and has been actively studied RikkenNature ; Rikken01 ; RikkenCNT ; RikkenSi ; Pop14 ; Morimoto-chiral16 ; Ideue ; Wakatsuki .

The other way to break TRS is incorporating irreversibility at the macroscopic level, i.e., (ii). Generalizing Onsager’s relation to nonlinear current response, we find that nonreciprocal current may appear due to the effect of dissipation/relaxation even when microscopic Hamiltonian obeys the TRS. Furthermore, there exists a systematic description for the second order nonlinear current responses that is based on the gauge invariant formulation of nonequilibrium Green’s functions under the static $E$ field Onoda06 ; Sugimoto08 . This formulation allows us to show that some electron interaction effects are necessary for nonreciprocal current response in bulk crystals under the TRS, on top of the dissipation effects. Such electron interactions include Coulomb electrons between electrons and electron-phonon interactions. In particular, it turns out that elastic scattering from disorder potential is not able to support nonreciprocal current response. These general symmetry considerations naturally lead us to study nonreciprocal current response with electron interactions in the Boltzmann transport picture that incorporates dissipation effects through relaxation of electron distribution functions.

Indeed, the situation changes in the presence of electron interactions, since electron interactions can modify the effective band structure when the applied electric field changes the electron distributions. In the steady state with nonzero current in noncentrosymmetric crystals, the interaction effect modifies the energy band in an asymmetric way with respect to the direction of $E$, as illustrated in Fig. 1(b), and enables us to circumvent the original constraint of TRS since systems with $E$ and $-E$ are not related with TRS and $\epsilon_{k,\sigma}(E)\neq\epsilon_{-k,\bar{\sigma}}(-E)$ in general. In the Boltzmann transport picture, this asymmetric change of effective band structure leads to nonreciprocal current responses. By using the gauge invariant formulation of nonequilibrium Keldysh Green’s functions, we derive a general formula for the nonreciprocal current in the weak interaction limit. It shows that nonreciprocal current in inversion broken materials is proportional to the strength of electron interaction and inversely proportional to typical band separation and bandwidth. We find that the nonreciprocity in noncentrosymmetric crystals is a quantum mechanical effect that is described by the complex nature of Bloch wave functions and interband matrix element which is unique to inversion broken systems. An estimate of the nonreciprocity shows that doped semiconductors and molecular conductors are good candidate materials for the nonreciprocity from electron correlation. We also discuss possible nonreciprocity in the molecular conductor TTF-CA Tokura88 ; Mitani87 .

In this paper, we focus on the nonlinearity in I-V characteristic of dc transport. Meanwhile, there are other nonlinear current responses in $\mathcal{I}$-broken crystals which have their origins in the complex nature of Bloch wave functions and should be compared with the nonreciprocal response in the present case. One example is a shift current Kraut ; Sipe ; Young-Rappe ; Cook17 ; Morimoto-shift16 , a dc current induced by photoexcitation of electrons beyond the band gap. The shift current is generated from the shift of wave packet centers for valence and conduction bands, and this shift is essentially described by Berry phases of valence and conduction electrons. The present nonreciprocal response and the shift current have similarity in that both rely on the multi-band nature of $\mathcal{I}$-broken systems. Yet, an important difference is that the nonreciprocal response arises from intraband metallic transport that is induced by static electric fields, while the shift current involves optical excitation of interband electron-hole pairs with photon energy larger than the band gap. Other examples are nonlinear Hall effect and low-frequency circular photogalvanic effect (CPGE)  Moore-Orenstein10 ; Sodemann-Fu15 ; Morimoto-semiclassical16 . They are known as geometrical effects described by the Berry curvature dipole of Bloch electrons. They are similar to the present nonreciprocal response in that both are intraband effects. However, the nonlinear Hall effect and the geometrical part of the CPGE are transverse (Hall) responses, in that they are described by off-diagonal components of the nonlinear conductivity tensor ($\sigma_{abb}$ and $\sigma_{aab}$, respectively, with $a\neq b$). In this sense, they are contrasted to the present nonreciprocal current which is a longitudinal current response described by diagonal components $\sigma_{aaa}$ and essentially involves the effect of dissipation.

Results

Time reversal symmetry constrains nonreciprocal current responses in bulk crystals. Based on general symmetry considerations, we show that nonreciprocal current response in crystals generally require two ingredients: (i) dissipation, and (ii) interactions. First, we generalize Onsager’s theorem to nonlinear current responses and show that the effect of dissipation is crucial for nonreciprocal current response. We then show by using gauge invariant formulation of Keldysh Green’s function that nonreciprocal current generally requires some interactions (e.g., electron-electron interactions and electron-phonon interactions). These two conditions suggest that the nonreciprocal current response can be captured by Boltzmann equation picture (that incorporates relaxation of electron distribution function) once we incorporate $E$-linear change of band structure induced by electron interactions.

The nonreciprocal current response is captured by an $E^{2}$ term in the current response. In the Boltzmann transport picture, the current $J$ induced by the applied electric field $E$ is given by

[TABLE]

with the relaxation time $\tau$ and the Fermi velocity $v_{F}$, for a one-dimensional system as depicted in Fig. 1. In noncentrosymmetric systems, the effective band structure with correlation effect can change asymmetrically in an applied electric field, and the Fermi velocity is modified as $v_{F}(E)=v_{F,0}+cE+O(E^{2})$. Therefore, noncentrosymmetric systems can host nonreciprocal current response given by the $E^{2}$ term in $J=(2e^{2}/\hbar)\tau(v_{F,0}E+cE^{2})$. Since the $E$-linear change of the band structure is described by the self energy linear in $E$, we study Green’s function and self energy in the steady state realized with the applied electric field. By using these results, we derive the general formula of nonreciprocal current, and then apply it to Rice-Mele model which is a prototypical model of ferroelectrics.

Onsager’s theorem and its generalization. In this section, we present a general consideration on the nonreciprocal current response in terms of the time reversal symmetry. We generalize Onsager’s relationship to nonlinear current responses, and show that the effect of dissipation is crucial for nonreciprocal current response.

In the linear response, Onsager’s relationship indicates that the conductivity $\sigma_{ij}$ is constrained as

[TABLE]

when the microscopic Hamiltonian preserves time reversal symmetry Onsager-reciprocal . This relationship is derived by considering the time reversal transformation in the Kubo formula for the linear conductivity as explained in Methods. Now we study how Onsager’s theorem can be extended to nonlinear current responses. We consider the second order current response,

[TABLE]

For systems of noninteracting electrons, the nonlinear conductivity $\sigma_{ijj}(i\omega_{n_{1}},i\omega_{n_{2}})$ in the imaginary time formalism satisfies the relationship

[TABLE]

under time reversal symmetry. (For the derivation, see Method section.) Naively, this seems to suggest that the nonlinear conductivity $\sigma_{ijj}(\omega_{1},\omega_{2})$ vanishes in the dc limit ($\omega_{1}\to 0$ and $\omega_{2}\to 0$). However, there is a subtlety in the analytic continuation to real frequencies as follows. We notice that nonlinear conductivity with Matsubara frequencies in the upper half plane is transformed to that with Matsubara frequencies in the lower half plane. Since the real axis is a branch cut in the complex $\omega$ plane, the analytic continuation of $i\omega_{n}\to 0$ for the two quantities, $\sigma_{ijj}(i\omega_{1},i\omega_{2})$ and $\sigma_{ijj}(-i\omega_{2},-i\omega_{1})$, lead to different results in general. This indicates that relationship similar to the Onsager’s relation does not necessarily constrain the dc nonlinear conductivity to vanish.

Interestingly, the extended Onsager’s relation in the above shows that nonreciprocal current response (nonzero $\sigma_{ijj}$) inevitably involves macroscopic irreversibility, i.e., the effect of dissipation, since the branch cut at $\textrm{Im}[\omega]=0$ is associated with macroscopic irreversibility. Specifically, such discontinuity for $\omega\to\pm 0i$ appears in the self energy by incorporating dissipative processes such as impurity scattering. To see this, it is useful to consider the case of linear conductivity. Metallic conductivity $\sigma_{xx}(\omega)$ has a branch cut and the limit of $\omega\to+0i$ gives a dissipative current response which is proportional to the relaxation time $\tau$. In contrast, Hall conductivity $\sigma_{xy}(\omega)$ does not involve such branch cut and corresponds to nondissipative current response (independent of $\tau$). Therefore, the nonreciprocal current response requires dissipation and should be proportional to the relaxation time $\tau$.

In passing, we note that Eq. (4) also indicates that the dissipation is essential for shift current which is a photocurrent caused by an optical resonance at a frequency $\omega$ above the band gap and described by $\sigma_{ijj}(\omega,-\omega)$ Nagaosa-Morimoto . If there is no effect of dissipation, we can naively take analytic continuation of Eq. (4), which leads to $\sigma_{ijj}(\omega,-\omega)=-\sigma_{ijj}(\omega,-\omega)=0$. Thus nonzero shift current requires some irreversibility. This observation is coherent with the fact that shift current essentially relies on optical absorption which is an irreversible process.

Absence of dc nonreciprocal current in noninteracting systems. In this section, we show that dc nonreciprocal current response does not appear when we do not incorporate effects of electron interactions that cause an effective change of the band structure under the applied electric field. We first show that no nonreciprocal current response appears in periodic systems. We then generalize the proof to the systems with static disorder potentials and show that incorporating the effects of elastic scattering does not lead to nonreciprocal current response.

We study systems with an applied electric field by using Keldysh Green’s function and its gradient expansion Rammer86 ; Jauho94 ; Kohler ; Kamenev . In particular, we use its gauge invariant formulation which enables us to treat the effect of $E$ directly Onoda06 . In the presence of a constant external electric field $E$, the Green’s function and the self energy are expanded with respect to $E$ as Onoda06 ; Sugimoto08

[TABLE]

where we set $\hbar=1,e=1$ for simplicity. The unperturbed part of the Green’s function $G_{0}$ is given by

[TABLE]

with the unperturbed Hamiltonian $H$ (without $E$). The linear order correction to the Green’s function $G_{E}$ is given by

[TABLE]

In order to describe the nonequilibrium steady state with applied electric fields, we suppose that the system is coupled to a heat bath. The coupling to the heat bath stabilizes the nonequilibrium electron distribution, and is incorporated through the self energy $\Sigma_{0}$ as $\Sigma_{0}^{R/A}(\omega)=\mp i\Gamma/2$ and $\Sigma_{0}^{K}(\omega)=i\Gamma f(\omega)$, where $\Gamma$ is the coupling strength and $f(\omega)$ is the Fermi distribution function (for details, see Methods) Rammer86 ; Jauho94 .

The second order current response is given by the expectation value,

[TABLE]

In order to show the absence of the second order current response in noninteracting systems, we set $E$-dependent self energy corrections to zero ($\Sigma_{E}=\Sigma_{E^{2}}=0$). Here, we also assumed that the heat bath coupled to the system (and gives $\Sigma_{0}$) is large enough such that it is not modified with applying electric fields. With vanishing $E$-dependent self energies, $G_{E^{2}}$ can be written as Onoda06

[TABLE]

We can show that the expectation value $J_{E^{2}}$ vanishes in the presence of TRS as follows. The TRS defined with $\mathcal{T}=K$ constrains Green’s functions and velocity operator as

[TABLE]

where $T$ denotes transposition with respect to the band index. This transformation law leads to cancellation of the integrand of $J_{E^{2}}$ between $k$ and $-k$. For example, the first term in $G_{E^{2}}$ in Eq. (10) gives the contribution which transforms as

[TABLE]

and cancels out between $k$ and $-k$. (In the last line, we used $\textrm{tr}A=\textrm{tr}A^{T}$.) We can show the cancellation for other terms in $J_{E}$ in a similar way. This indicates that the nonlinear current $\propto E^{2}$ vanishes under the TRS in bulk crystals if we do not incorporate $E$-linear band modification described by $\Sigma_{E}$.

It is easy to generalize the above argument to systems with static disorder potential. We consider a system of the system size $L$ with the periodic boundary condition. We introduce a phase twist at the periodic boundary with the phase $\theta$. In this case, the velocity matrix element $v$ and the nonequilibrium Green’s function $G_{E^{2}}$ become functions of the phase twist $\theta$ instead of the momentum $k$. When the disorder is uniform and the system has translation symmetry on average, physical quantities are obtained by averaging over the phase twist $\theta$. We note that this procedure is very similar to the discussion of Chern number in quantum Hall systems with disorder potential Niu85 . Thus, the nonlinear current response $J_{E^{2}}$ is given by a similar expression to Eq. (9) by replacing $k$ with $\theta$. [The expression for $G_{E^{2}}(\omega,\theta)$ is also obtained by replacing $k$ with $\theta$ in Eq. (10).] Since similar symmetry constraints hold for $G$ and $v$ under the TRS [i.e., $G_{0}(\omega,\theta)=G_{0}(\omega,-\theta)^{T},G_{E}(\omega,\theta)=-G_{E}(\omega,-\theta)^{T},v(\theta)=v(-\theta)^{T}$ ], the integrand of $J_{E^{2}}$ satisfies

[TABLE]

and cancels between $\theta$ and $-\theta$. This proves that elastic scattering from static disorder potential does not induce nonreciprocal current response.

These considerations indicate that $E$-linear change of band structure ($\Sigma_{E}$) is essential for nonreciprocal current response in bulk crystals. The $E$-linear change of band structure requires some kind of electron interactions, such as Coulomb interactions and electron-phonon interactions. Since the current response proportional to $E^{2}$ arises from the $E$-linear change of band structure in the Boltzmann transport picture, it suffices to consider $\Sigma_{E}$ and neglect $\Sigma_{E^{2}}$. Although we can study this nonreciprocal current response by directly looking at $G_{E^{2}}$ with incorporating $\Sigma_{E}$, it is equivalent and more concise to compute $\Sigma_{E}$ and then use the relationship Eq. (1) with the Fermi velocity modified by $E$.

So far, we discussed general conditions to achieve nonreciprocal current response in bulk crystals. In order to proceed to explicit calculations of nonreciprocal current, we need to specify the form of the self energy, i.e., how the self energy $\Sigma$ is expressed in terms of the Green’s function $G$. We consider electron-electron interaction shown in the Feynman diagram Fig. 2(a) and show that it gives rise to nonreciprocal current through $E$-linear band structure change. Incidentally, we also show explicitly that elastic scatterings from isotropic impurity potential [Fig. 2(b)] does not lead to nonreciprocal current, which is consistent with the above general symmetry consideration.

Nonequilibrium steady state under the applied electric field. Now we move on to demonstration of nonreciprocal current responses with electron interactions by performing explicit calculations. We consider the cases of weak interactions and perform Hartree-Fock approximation in the gauge invariant formulation of Keldysh Green’s functions. In order to describe $E$-linear change of the effective band structure, we first study the nonequilibrium steady state under the electric field by looking at $G_{E}^{<}$. Once $G_{E}^{<}$ is obtained, we can compute the $E$-linear change of band structure by studying $\Sigma_{E}^{R}$ that corresponds to the diagram in Fig. 2.

The $E$-linear change of electron occupation has intraband and interband contributions, since the Green’s function for a $\mathcal{I}$-broken system generally has a matrix structure with respect to band index. The intraband contribution is written as

[TABLE]

for the band 1 that we assume crosses with the Fermi energy (for details, see Methods). Here, $k_{F,i}$ are the Fermi momenta for the band 1. This change of the lesser Green’s function linear in $E$ describes the effect of the applied electric field where the electron occupation is shifted in the momentum space as $k\to k+\tau E$ near the Fermi surface (with $\tau=2\pi/\Gamma$). This coincides with the picture of the semiclassical Boltzmann equation as illustrated in Fig. 1(a).

Next, the interband contribution for $G_{E}^{<}$ is given by

[TABLE]

and $G_{E,21}^{<}=-(G_{E,12}^{<})^{*}$, for the bands 1 and 2 (for details of the derivation, see Methods.). Here, we assume that the band 1 is the partially filled valence band and the band 2 is the unoccupied conduction band as illustrated in Fig. 1, and $E_{g,k}$ denotes the band gap at the momentum $k$. This term arises from a quantum mechanical effect that the electric field also modifies the wave function in addition to the shift of the momentum at the Fermi energy. Thus the electron distribution in the steady state effectively has an interband component near the Fermi energy. We note that this interband component of $G^{<}_{E}$ cannot be captured by semiclassical treatment with Boltzmann equation, and is a quantum effect captured by the current approach that uses the gauge invariant formulation of Keldysh Green’s functions. This interband component gives the origin of the nonreciprocity when the electron interaction is incorporated. In contrast, when we consider effects of scattering by short-range impurities within the Born approximation [described by the diagram in Fig. 2(a)], we do not find the $E$-linear change of the effective band structure, as detailed in Methods.

Formula of nonreciprocal current in two band systems. Now we show that nonreciprocal current appears from $E$-linear band structure change once we introduce electron-electron interactions, and derive a general formula for nonreciprocal current in two-band systems. The effect of electron interactions is minimally incorporated by the self energy arising from the Hartree contribution to $\Sigma_{E}^{R}$ as shown in Fig. 2(a).

For simplicity, we consider a two-band model, where the unit cell contains two sites, and the wave functions of valence and conduction bands (labeled by $1$ and $2$, respectively) are represented by

[TABLE]

We then consider two copies of the original system, each labeled by $\uparrow$ and $\downarrow$, and introduce the onsite interaction given by

[TABLE]

with the site index $i$. We treat the effects of the onsite interaction in terms of Hartree-Fock approximation, and study the effective band structure. Since the two copies ($\uparrow$ and $\downarrow$) are decoupled in the noninteracting Hamiltonian, only the Hartree term appears in the present case. (We suppose that the Hartree correction in the equilibrium is already included in the original Hamiltonian.) In the following, we focus on the electronic structure of the $\uparrow$ component, and suppress the label for the two copies for simplicity. By using the momentum space representation of $H_{\textrm{int}}$ (for details, see Methods), the self energy from the Hartree contribution is given by

[TABLE]

with the lattice constant $a$. Now we assume that there are two Fermi momenta at $\pm k_{F}$ with the same Fermi velocity $v_{F}$. By using the Green’s function in the steady state [Eq. (17)], this is expressed as

[TABLE]

This self energy is an even function with respect to $k$ from TRS (such as $\mathcal{T}=\mathcal{K}$), which is important in obtaining nonreciprocal current response as we will see next.

We now study the nonreciprocal current response by using the self energy $\Sigma_{E}^{R}$. The current induced by an electric field (linearly in $E$) is given by

[TABLE]

from the Boltzmann transport approach. An application of the electric field modifies the band structure as $\epsilon_{1}\to\epsilon_{1}+\frac{E}{2}\Sigma_{E,11}^{R}(k)$, and hence, the Fermi velocity as $v_{11,k_{F}}\to v_{11,k_{F}}+\frac{E}{2}\partial_{k}\Sigma_{E,11}^{R}(k)$. Since the obtained self energy $\Sigma_{E,11}^{R}(k)$ is an even function of $k$, the velocity corrections at $\pm k_{F}$ do not cancel out in evaluating the correction to the current response in Eq. (22). Thus, we obtain the nonlinear current response $\delta J$ (the part of current response proportional to $E^{2}$) as

[TABLE]

which is the general formula for two-band systems in one-dimension. This can be generalized to systems in higher dimensions if we replace the summation over the Fermi points with an integral over the Fermi surface. The above formula indicates that the nonreciprocity ratio $\gamma$ of the nonlinear current to the original current is roughly estimated as

[TABLE]

where $W$ is the band width. Here we used $u_{k},v_{k}\sim 1$ and $v_{11,k_{F}}\sim v_{12,k_{F}}$ for rough order estimates.

The obtained formula indicates that breaking of inversion symmetry is essential for the nonreciprocity. When the system is inversion symmetric, the wave function is expressed with real numbers due to the combination of inversion symmetry $\mathcal{I}$ and TRS ($\mathcal{IT}=\mathcal{K}$). Therefore, we obtain $\textrm{Im}[u_{k_{F}}v_{k_{F}}v_{12,k_{F}}]=0$ in inversion symmetric systems and no reciprocity appears. This clearly shows that the nonreciprocity in the current mechanism essentially relies on the complex nature of wave functions in noncentrosymmetric crystals.

Nonreciprocal current in Rice-Mele model. We study nonreciprocal current in a representative model of ferroelectrics, Rice-Mele model, by taking into account onsite interaction. We show that $E$-linear band structure change is associated with effective modulation of parameters in the Hamiltonian that is induced by the applied electric field $E$.

Rice-Mele model is a representative 1D two-band model with broken inversion symmetry, and is described by a Hamiltonian Rice-Mele ,

[TABLE]

Rice-Mele model is a minimal model for molecular conductors Su ; Nagaosa-Takimoto ; Onoda ; Nagaosa-Morimoto and ferroelectric perovskites Egami . In the momentum representation, the Hamiltonian reads

[TABLE]

where Pauli matrices $\sigma$’s act on two sublattices (A and B) in the unit cell, and $a$ is the lattice constant. For Rice-Mele model, the wave functions in Eq. (18) are given by $u_{k}=-\sin\frac{\theta}{2}$ and $v_{k}=e^{i\phi}\cos\frac{\theta}{2}$, with the parameters $\theta=\cos^{-1}\frac{\Delta}{|\epsilon_{1,k}|}$ and $\phi=\tan^{-1}\delta t$. The energy dispersion for the valence band is given by

[TABLE]

and $\epsilon_{2,k}=-\epsilon_{1,k}$ for the conduction band, as shown in the right panel of Fig. 3 with black line.

We again consider two copies of Rice-Mele model and introduce the onsite interaction given by

[TABLE]

where $\uparrow$ and $\downarrow$ label the two identical copies. By focusing on the electronic structure of the $\uparrow$ component, we suppress the label for the two copies for simplicity.

Now we study the nonequilibrium steady state under the electric field $E$ applied along the 1D chain, by using the nonequilibrium Green’s functions, Eq. (17). The electric field is described by the Hamiltonian, $H_{\textrm{ele}}=-eEa\sum_{i}in_{i}$. The application of the electric field effectively changes the parameters $\delta t$ and $\Delta$, which can be easily obtained within the Hartree approximation since the expectation values are directly computed from the lesser component of the Green’s function. From the Hartree term, the occupation of site A is modified as

[TABLE]

(For details of the derivation, see Methods.) Similarly, the occupation of site B is modified in the opposite way as $\delta n_{B}=-\delta n_{A}$ . Thus the Hartree term effectively changes the staggered potential $\Delta$ as

[TABLE]

Notice that the change of $\Delta$ is opposite in sign depending on the direction of $E$. This situation is schematically illustrated in Fig. 3. Since the parameter changes are asymmetric with respect to the sign of $E$ in the $\delta t-\Delta$ space, the effective band structure $\epsilon_{1,k}(E)$ becomes different for the electric fields $+E$ and $-E$. The nonlinear current in the nonequilibrium steady state is obtained from the conventional Boltzmann equation approach for this modified band structure in the presence of $E$. Namely, the linear conductivity is given by

[TABLE]

with $E$-dependent Fermi velocity $v_{F}(E)=\partial_{k}\epsilon_{1,k}(E)|_{k=k_{F}}$ (where $k_{F}$ is the Fermi momentum). The $E$-linear change of the effective band structure leads to $E$-linear term in $v_{F}(E)$, which results in the nonlinear current response $\propto E^{2}$. Thus the asymmetry in band structure changes leads to the nonreciprocity of the current with respect to the direction of $E$.

The nonreciprocity is quantified by the ratio of the change of electric conductivity $\gamma=[\sigma(E)-\sigma(0)]/\sigma(0)$ in the presence of the applied electric field. We note that $\gamma=\delta J/J$ and the crude approximation is given in Eq. (24). This approximation is also obtained from the $E$-linear change of the parameters in Rice-Mele model in Eq. (30) and Eq. (27). Explicit evaluation of Eq. (23) gives the nonreciprocity ratio of $\gamma=5\times 10^{-7}$ for typical parameters of Rice-Mele model ($\delta t=\Delta=0.3t,U=t,k_{F}=0.1\pi/a$ along with $t=1\textrm{ eV}$ and $a=1\textrm{ \AA}$) and the electric field of $E=10^{5}\textrm{ V/m}$. This order of the nonreciprocity is comparable to those in materials showing magnetochiral anisotropy Rikken01 , as we will discuss further in the discussion section.

Discussions

Finally, we give an estimate of the nonreciprocal response induced by the present mechanism for realistic materials. Typical magnitude of the nonreciprocity is determined by $\gamma=\delta J/J$ in Eq. (24). When the band gap and Coulomb energy are both of the order of 1eV, the ratio $\delta J/J$ reduces to $eEa/W$, which is the ratio between the electric potential in the unit cell and the bandwidth. This allows us to estimate typical nonreciprocity as follows. We consider the current of 1mA that flows in a wire of the area 1mm2, which amounts to a current density of $j=10^{3}\textrm{A/m}^{2}$. For usual metals, conductivity is roughly given by $\sigma\simeq 10^{6}\textrm{A/Vm}$, and hence, the electric field present in the wire is $E=j/\sigma\simeq 10^{-3}\textrm{V/m}$. In this case, the electric potential in the unit cell of $a\simeq 1\AA$ is $eEa\simeq 10^{-13}\textrm{eV}$. Since the bandwidth is typically 1eV, this indicates nonreciprocity ratio is $\delta J/J\simeq 10^{-13}$. This should be compared to the typical order of nonreciprocity for materials showing magnetochiral anisotropy. Bi helix Rikken01 and molecular solids Pop14 show the nonreciprocity measured in resistivity change $\delta\rho$ as $\delta\rho/\rho=\gamma^{\prime}IB$ with $\gamma^{\prime}\simeq 10^{-3}A^{-1}T^{-1}$. For $I=1$mA and $B=1$T, the typical nonreciprocity is $\delta J/J\simeq\delta\rho/\rho\simeq 10^{-6}$. Thus the nonreciprocity induced by electron correlation is very small for good metals. On the other hand, we can expect comparable nonreciprocity for doped semiconductors whose conductivity ranges from $10^{-1}\sim 10^{5}$A/Vm. For example, for the doped Si of $\sigma=10^{-1}$A/Vm and the bandwidth $W\simeq 1$ eV in the presence of the current density $j=10^{3}\textrm{A/m}^{2}$, we obtain the nonreciprocity of $\delta J/J\simeq 10^{-6}$, which becomes comparable with typical materials showing magnetochiral anisotropy.

Another candidate is the molecular conductor TTF-CA which is a strongly correlated insulator. Of course, our theory for weakly correlated metals is not directly applicable. However, it is interesting to estimate the nonreciprocity ratio anyway, since the carriers in TTF-CA (thermally activated or provided by impurity sites) may be treated as electrons having a Fermi surface, and the Hartree approximation sometimes becomes a good approximation at least for the ground states. The typical order of electric field that can be applied is $E\simeq 10^{5}$V/m Tokura88 ; Mitani87 . Since the lattice constant is $a\simeq 1$nm, the electric voltage in the unit cell becomes $eEa\simeq 10^{-4}$ eV, and the band width is given by $W\simeq 0.2$ eV. Thus the nonreciprocity ratio can be $10^{-3}$ which may be comparable with that in magnetochiral anisotropy in Bi helix Rikken01 . We again note that this is a number obtained from a naive application of Eq. (24) to TTF-CA beyond the applicability of our theory, but this suggests that it is an interesting future problem to study TTF-CA as a candidate of strongly correlated materials for nonreciprocity, from both theoretical and experimental points of view.

Our analysis is mostly valid for weakly interacting systems because we adopted Hartree approximation to incorporate the correlation effect. Therefore, the study of nonreciprocal responses of strongly interacting cases remains as an interesting future problem. Meanwhile, our symmetry considerations from generalization of Onsager’s theorem suggests that nonreciprocal current response can generally appear in the presence of dissipation and interactions, regardless of the strength of the interaction. We may also note that Hartree approximation sometimes gives a good description for some ground state properties, even for strong $U$ cases, such as magnetically ordered ground states. Our approach may give a good approximation for nonlinear properties of those states, since the nonreciprocal current response is a nonequilibrium property near the ground state under a moderate electric field.

Methods

Derivation of generalized Onsager’s theorem. In this section, we present general symmetry considerations on the nonreciprocal current response with respect to the time reversal symmetry by extending Onsager’s relationship to nonlinear current.

We consider a system of noninteracting electrons that are described by Green’s function in the Lehmann representation,

[TABLE]

where $|\alpha\rangle$ is a many-body state that satisfies $\hat{H}|\alpha\rangle=E_{\alpha}|\alpha\rangle$ with the many-body Hamiltonian $\hat{H}$, $\beta$ is the inverse temperature, $e^{-\beta\Omega}=\textrm{Tr}[e^{-\beta\hat{H}}]$, and $c_{a}$ and $c_{a}^{\dagger}$ are annihilation and creation operators of an electron with a single particle state $a$. (Here $\alpha,\beta$ are labels for many-body states, whereas $a,b$ are labels for single particle states.) We write the current operator $\hat{v_{i}}$ along the $i$th direction as

[TABLE]

where $v$ is a matrix for a velocity operator in the single particle representation.

In the linear response, Onsager’s relationship indicates that the conductivity $\sigma_{ij}$ is constrained as

[TABLE]

in the presence of time reversal symmetry Onsager-reciprocal . This relationship is derived by considering the time reversal transformation in the Kubo formula for the linear conductivity,

[TABLE]

where $i\omega_{n},i\omega_{m}$ are Matsubara frequencies, and tr is a trace over single particle states (labeled by $a,b$). The time reversal symmetry, $\mathcal{T}=K$, indicates

[TABLE]

These actions of $\mathcal{T}$ in the many-body representation are obtained by using $\mathcal{T}|\alpha\rangle=(|\alpha\rangle)^{*}$ in Eq. (32) and Eq. (33). (We note that this is closely related to symmetry constraint in a single particle Hamiltonian, $H(k)=H^{T}(-k)$, in the momentum representation.) By using these relationships, the Kubo formula can be rewritten as

[TABLE]

and leads to the Onsager’s relationship. Here we rewrote the trace in the reverse order in the second line and used the fact that the transposition in the trace does not change its value.

Next we study how Onsager’s theorem can be extended to nonlinear current responses. We consider the second order current response,

[TABLE]

The nonlinear conductivity $\sigma_{ijj}(i\omega_{n_{1}},i\omega_{n_{2}})$ has a contribution from a triangle diagram which is given by

[TABLE]

since there are no vertex corrections for noninteracting systems. The time reversal symmetry indicates that

[TABLE]

Naively, this seems to suggest that the nonlinear conductivity $\sigma_{ijj}(\omega_{1},\omega_{2})$ vanishes in the dc limit ($\omega_{1}\to 0$ and $\omega_{2}\to 0$). However, we notice that nonlinear conductivity with Matsubara frequencies in the upper half plane is transformed to that with Matsubara frequencies in the lower half plane. Since the real axis is a branch cut in the complex $\omega$ plane, the analytic continuation of $i\omega_{n}\to 0$ for the two quantities, $\sigma_{ijj}(i\omega_{1},i\omega_{2})$ and $\sigma_{ijj}(-i\omega_{2},-i\omega_{1})$, lead to different results in general. This indicates that relationship similar to the Onsager’s relation does not necessarily constrain the dc nonlinear conductivity to vanish. Instead, this extended Onsager’s relation indicates that nonreciprocal current necessarily involves irreversibility such as dissipation and relaxation.

In a similar manner, we can also derive an extended Onsager’s relation for shift current. Shift current is dc current induced by optical absorption above the band gap and photoexcitation of electron-hole pairs that have finite polarization Sipe ; Morimoto-shift16 . It is described by a nonlinear current response, $J_{i}(\omega_{1}+\omega_{2})=\sigma_{ijj}^{\textrm{shift}}(\omega_{1},\omega_{2})E(\omega_{1})E(\omega_{2})$ with $\omega_{2}\approx-\omega_{1}$. The nonlinear conductivity $\sigma^{\textrm{shift}}$ has two contributions as $\sigma^{\textrm{shift}}=\sigma^{\textrm{tr}}(\omega_{1},\omega_{2})+\sigma^{\textrm{bubble}}(\omega_{1},\omega_{2})$, where the latter piece is a correlation function of paramagnetic current $\hat{v}_{i}$ and diamagnetic current $\hat{v}_{\textrm{dia,ij}}\equiv\sum_{ab}(v_{\textrm{dia},ij})_{ab}c_{a}^{\dagger}c_{b}$ Kim17 . Time reversal symmetry leads to the same relation,

[TABLE]

since $\sigma^{\textrm{bubble}}$ also obeys the same transformation law under the TRS with $\sigma^{\textrm{tr}}$ as follows. In the momentum representation, matrix elements for diamagnetic current are given by

[TABLE]

Accordingly, TRS constrains diamagnetic current operator as

[TABLE]

due to an extra $k$ derivative. The nonlinear conductivity for shift current is written as

[TABLE]

Under TRS, this transforms as

[TABLE]

where we used the symmetry between $i\omega_{n_{1}}$ and $i\omega_{n_{2}}$ to fit the transformation law with that for $\sigma^{\textrm{tr}}$. Therefore, nonzero shift current also requires irreversibility that introduces a branch cut at the real axis in the $\omega$ space and makes two limits $\omega\to\pm i0$ different. In this case, the irreversibility comes from optical transition and creation of electron hole pairs across the band gap.

Keldysh Green’s function. In this section, we summarize basic notations of Keldysh Green’s functions that we need for our discussion Jauho94 ; Kamenev ; Aoki-RMP14 . In the Keldysh Green’s function formalism, we consider the Keldysh component of the Green’s function in addition to the retarded and advanced Green’s function. Keldysh component describes the electron occupation in the nonequilibrium state, while the retarded and advanced components describe the spectrum of the system. The Dyson equation for the Green’s function is given by

[TABLE]

with the Hamiltonian $H$.

In the thermal equilibrium, the Keldysh Green’s function is obtained by solving the Dyson equation. We suppose that the system is weakly coupled to a heat bath with broad spectrum, which determines the electron distribution of the system. The coupling to the heat bath (such as electron reservoirs) is described by the self energy given by Jauho94

[TABLE]

where $\Gamma$ is the strength of the coupling to the bath, and $f(x)=1/[1+\exp(x/k_{B}T)]$ is the Fermi distribution function with the temperature $T$.

The observables in the nonequilibrium steady state is obtained from Keldysh Green’s function. We define the lesser component of the Green’s function as

[TABLE]

By using $G^{<}$, we can write the expectation value of a general fermion bilinear as

[TABLE]

The lesser Green’s function is concisely obtained from the equation

[TABLE]

where the lesser component of the self energy encodes the information of the electron distribution and is given by

[TABLE]

Keldysh Green’s function under the applied electric field. In this section, we study the nonequilibrium electron distribution realized under the applied electric field. We compute the $E$-linear part of the lesser Green’s function $G_{E}^{<}$ in Eq. (6) in the gauge invariant formulation. In doing so, we use the diagram in Fig. 2(b) to specify the form of self energy $\Sigma^{<}_{E}$ in Eq. (6). (We note that the electron interaction in Fig. 2(a) does not change electron distribution and does not contribute to $\Sigma^{<}_{E}$. Furthermore, it turns out in the end that the contribution of impurity scattering to $\Sigma^{<}_{E}$ is actually negligible under TRS.) Specifically, we consider the delta function type impurity [$V(r)=u\delta(r-r_{0})$ with density $n$]. In the second-order Born approximation, the self energy is given by

[TABLE]

which corresponds to the diagram in Fig. 2(b). In the right hand side, $G_{E}$ denotes the bare Green’s function that does not include the effect of impurity scattering. We note that the impurity scattering also modifies the self energy $\Sigma_{0}$ in the zeroth order in $E$, but this correction only changes the coupling $\Gamma$ in Eq. (48) and can be absorbed by redefining $\Gamma$ accordingly.

The current response in the nonequilibrium steady state under the electric field is captured by the lesser Green’s function $G_{E}^{<}$ which gives a contribution linear in $E$. We consider a multiband system and suppose that the Bloch wave functions are given by $\Psi_{i,k}$ which satisfy $H\Psi_{i,k}=\epsilon_{i,k}\Psi_{i,k}$ with the energy dispersion $\epsilon_{i,k}$ (where $i$ is the band index). First we start with the intraband component of $G_{E,ii}^{<}$ for the band $i$ (where we omit the band index $i$ in the following, for simplicity). By assuming that $G_{0}$ and $v$ have a single component, equation (8) gives

[TABLE]

where we used $\partial_{\omega}f(\omega)=-\delta(\omega-\epsilon_{F})$. This expression is simplified by using the relationship

[TABLE]

that holds for a positive integer $n$, and the Fermi momenta $k_{F,i}$, where we only keep the leading order in terms of $1/\Gamma$. By using Eq. (54) with Eq. (53), the impurity scattering gives rise to $\Sigma_{E}^{<}(\omega)$ given by

[TABLE]

Here, the numerator in the right hand side vanishes since the TRS leads to $\sum_{k_{F,i}}\frac{v_{k}}{|v_{k}|}=0$, and hence, $\Sigma_{E}^{<}(\omega)=0$ follows. Thus we obtain

[TABLE]

This change of the lesser Green’s function linear in $E$ describes the effect of the applied electric field where the electron occupation is shifted in the momentum space as $k\to k+\tau E$ near the Fermi surface (with $\tau=2\pi/\Gamma$). This corresponds to the picture from the semiclassical Boltzmann equation as illustrated in Fig. 1(a).

Next we consider the interband component, $G_{E,12}^{<}$, by focusing on the valence and conduction bands which are labeled by 1 and 2, respectively. Equation (8) gives

[TABLE]

We assume that the Fermi energy is located within the band 1 and does not cross the band 2. In this case, the second term in the right hand side reduces to

[TABLE]

with $E_{g,k}=\epsilon_{2,k}-\epsilon_{1,k}$, where we only kept the leading term with respect to $1/E_{g,k}$. (Here we used Eq. (55) for $G_{0,11}^{R}G_{0,11}^{A}$ and discarded the second term.) Since the right hand side is inversely proportional to the band gap $E_{g,k}$, the self energy $\Sigma_{E,12}^{<}$ obtained from Eq. (53) is proportional to $\Gamma/E_{g,k_{F}}$, which is negligible in the left hand side of the above equation given that $G_{0,22}^{A}\propto 1/E_{g,k_{F}}$. Therefore the lesser part of the Green’s function is given by

[TABLE]

We note that $G_{E,21}^{<}$ is obtained from the relationship

[TABLE]

as a consequence of the hermiticity of expectation values in Eq. (50).

Effective band dispersion with impurity scattering. In this section, we study the effective band dispersion in the presence of $E$ and impurity scattering by looking at $\Sigma_{E}^{R}$. We show that impurity scattering is insufficient for nonreciprocal current response because the change of the band dispersion turns out to be the same for positive and negative electric fields.

From Eq. (8), the retarded part of the equation for $G_{E}$ reads

[TABLE]

with $v_{k}=\partial_{k}H$, where we used $\partial_{\omega}\Sigma^{R}=0$. For simplicity, we consider a two-band system, where the Green’s function is given by

[TABLE]

where $i,j=1,2$ are labels for valence and conduction bands, respectively. For the diagonal components, we obtain

[TABLE]

since the second term in Eq. (62) vanishes trivially. The diagonal part of the self energy is momentum independent and vanishes as

[TABLE]

Off-diagonal part is determined from

[TABLE]

By integrating over the momentum, we obtain

[TABLE]

which leads to nonzero $\Sigma_{E,21}^{R}$ in general. Therefore, the effective Hamiltonian is given by

[TABLE]

and the effective band structure of the valence band in the presence of $E$ is obtained by diagonalizing $H$ as

[TABLE]

This is an even function with respect to $E$; the effective band structure depends on the strength of electric field $|E|$, but is independent of the direction of the applied field. Therefore, no reciprocal current appears when we use Boltzmann equation approach based on this modified band structure.

We note that this conclusion is not changed even when we treat the impurity scattering by self-consistent Born approximation. In the self-consistent Born approximation, $G_{E}$ in Eq. (53) is taken as a full Green’s function including the effect impurity scattering. In this case, the self energy $\Sigma_{E}$ is obtained by repeating the above calculation and taking convergence. In the every step of the repetition, the energy dispersion is modified as Eq. (68) and still gives a symmetric dispersion in $k$. After repeating this many times, the dispersion remains symmetric in $k$. Therefore, self-consistent treatment of impurity scattering still gives no nonreciprocal current response.

Electron-electron interaction in two-band model. In this section, we derive Eq. (20) for the self energy that arises from the electron-electron interaction in the case of a two-band model. We also derive Eq. (29) for the expectation value of density operators. These expressions are obtained by using the momentum space representation of the interaction Hamiltonian.

We consider the onsite interaction that is given by

[TABLE]

with the site index $n$. Expressing the Hartree contribution to the self energy requires momentum representations of the density operators $n_{A,i}$ and $n_{B,i}$, where we omit the indices for two copies ($\uparrow$ and $\downarrow$) since the expressions are identical for two copies. For the wave functions in Eq. (18), the creation operators of Bloch states are written as

[TABLE]

where $N$ is the system size. By using inverse Fourier transformation, the creation operators in the site basis are expressed with Bloch states as

[TABLE]

where $k$ runs momenta in the first Brillouin zone (e.g., $k=2\pi j/Na$ for $j=0,\ldots,N-1$ with lattice constant $a$). Now the density operators are given by

[TABLE]

The retarded part of the self energy is given by Hanai16

[TABLE]

by using the momentum space representation for the interaction $H_{\textrm{int}}$ which is given by

[TABLE]

with the band index $m_{i}$. The first term in Eq. (76) is the Hartree term, and the second is the Fock term. The momentum representations of the interaction that are relevant for the Hartree term contributing to $\Sigma_{E,11}^{R}$ are given by

[TABLE]

By using Eq. (76), the self energy $\Sigma_{E,11}^{R}$ is written as

[TABLE]

The first two terms in the integral vanishes due to TRS. Specifically, $G_{E,ii}(k)$ is an odd function of $k$ due to TRS as in Eq. (16), and $|u_{k}|^{2}$ and $|v_{k}|^{2}$ are even functions of $k$, which indicates that the first two terms vanish after integrating over $k^{\prime}$. Thus we end up with

[TABLE]

where we replaced the sum $\sum_{k^{\prime}}$ with the integral $Na\int\frac{dk}{2\pi}$.

Next, we derive the changes of the density $\delta n_{A}$ and $\delta n_{B}$ caused by the electric field $E$. By using Eq. (75) and Eq. (50), the change of the density at A site is given by

[TABLE]

Since the first and second terms vanish due to TRS, we obtain

[TABLE]

Similarly, the change of the density at A site is given by

[TABLE]

which is opposite in sign compared to $\delta n_{A}$.

Acknowledgements. This work was supported by the Gordon and Betty Moore Foundation’s EPiQS Initiative Theory Center Grant (TM), and by Grants-in-Aid for Scientific Research from the Ministry of Education, Science, Sports and Culture No. 24224009 and 25400317, CREST, Japan Science and Technology (grant no. JPMJCR16F1), and ImPACT Program of Council for Science, Technology and Innovation (Cabinet office, Government of Japan, 888176) (NN).