Theory for shift current of bosons: Photogalvanic spin current in ferrimagnetic and antiferromagnetic insulators

TL;DR

This paper presents a theoretical framework for understanding how linearly polarized electromagnetic waves generate dc spin currents in noncentrosymmetric magnetic insulators, highlighting a shift current mechanism for magnons.

Contribution

It introduces a non-linear response theory for magnons that explains dc spin current generation without thermal magnons, distinct from spinon-based mechanisms.

Findings

Finite dc spin current conductivity at zero temperature.

Insensitivity of spin current to Gilbert damping.

Feasible electric field intensities for observable spin currents.

Abstract

We theoretically study the optical generation of dc spin current (i.e., a spin-current solar cell) in ordered antiferromagnetic and ferrimagnetic insulators, motivated by a recent study on the laser-driven spinon spin current in noncentrosymmetric quantum spin chains [H. Ishizuka and M. Sato, Phys. Rev. Lett. 122, 197702 (2019)]. Using a non-linear response theory for magnons, we analyze the dc spin current generated by a linearly-polarized electromagnetic wave (typically, terahertz or gigahertz waves). Considering noncentrosymmetric two-sublattice magnets as an example, we find a finite dc spin current conductivity at $T=0$, where no thermally-excited magnons exist; this is in contrast to the case of the spinon spin current, in which the optical transition of the Fermi degenerate spinons plays an essential role. We find that the dc spin-current conductivity is insensitive to the Gilbert damping, i.e., it may be viewed as a shift current carried by bosonic particles (magnons). Our estimate shows that an electric-field intensity of $E\sim10^4-10^6$ V/cm is sufficient for an observable spin current. Our theory indicates that the linearly-polarized electromagnetic wave generally produces a dc spin current in noncentrosymmetric magnetic insulators.

Equations

$$H=\sum _{\mathbf{k}}{\epsilon }_{\alpha }(\mathbf{k}){\alpha }_{\mathbf{k}}^{†}{\alpha }_{\mathbf{k}}+{\epsilon }_{\beta }(\mathbf{k}){\beta }_{\mathbf{k}}{\beta }_{\mathbf{k}}^{†},$$

$$H^{\prime }=$$

$$+\text{h.c.},$$

$$J=$$

$$⟨J⟩(\Omega )=\sum _{\mu ,\nu }\int d\omega {\sigma }_{\mu \nu }(\Omega ;\omega ,\Omega -\omega )h_{\omega }^{\mu }{h}_{\Omega -\omega }^{\nu },$$

$${\sigma }_{\mu \nu }(\Omega ;\omega ,\Omega -\omega )=$$

$$\times \left[\frac{(B_{\mathbf{k}}^{\nu }{)}_{a_2{a}_3}(J_{\mathbf{k}}{)}_{a_3{a}_1}}{\Omega +{\tilde{\epsilon }}_{a_1}(\mathbf{k})-{\tilde{\epsilon }}_{a_3}(\mathbf{k})-i/(2{\tau }_{\mathbf{k}})}-\frac{(J_{\mathbf{k}}{)}_{a_2{a}_3}(B_{\mathbf{k}}^{\nu }{)}_{a_3{a}_1}}{\Omega +{\tilde{\epsilon }}_{a_3}(\mathbf{k})-{\tilde{\epsilon }}_{a_2}(\mathbf{k})-i/(2{\tau }_{\mathbf{k}})}\right],$$

$${\tilde{\epsilon }}_a(\mathbf{k})=$$

$$\mathrm{sgn}(a)=$$

$${\tilde{\rho }}_{\mathbf{k},a}=$$

$$\sigma$$

$$-$$

$$+$$

$$\mathrm{Re}\left[\sigma (0;\omega ,-\omega )\right]=$$

$$-\frac{1}{\pi }\sum _{\mathbf{k}}\mathrm{Re}\left\{\frac{(B_{\mathbf{k}}{)}_{\beta \alpha }(A_{\mathbf{k}}{)}_{\alpha \beta }((B_{\mathbf{k}}{)}_{\beta \beta }+(B_{\mathbf{k}}{)}_{\alpha \alpha })}{{\omega }^2-({\epsilon }_{\alpha }(\mathbf{k})+{\epsilon }_{\beta }(\mathbf{k})+i/2{\tau }_{\mathbf{k}}{)}^2}\right\}.$$

$$H_{\mathrm{tot}}=$$

$$H_0\equiv$$

$$H_{1D}(r_y,r_z)\equiv$$

$$+J(1-\delta ){\mathbf{S}}_A(\mathbf{r}+\hat{x})\cdot {\mathbf{S}}_B(\mathbf{r})$$

$$-(D+D_s){\left[S_A^z(\mathbf{r})\right]}^2-(D-D_s){\left[S_B^z(\mathbf{r})\right]}^2$$

$$-h\left[g_A{S}_A^z(\mathbf{r})+g_B{S}_B^z(\mathbf{r})\right],$$

$$H_Z^{(\omega )}=$$

$$J_x^z=$$

$$+(1-\delta )\left\{S_A^x(\mathbf{r}+\hat{x})S_B^y(\mathbf{r})-S_A^y(\mathbf{r}+\hat{x})S_B^x(\mathbf{r})\right\},$$

$$S_A^z=$$

$$S_A^{+}(\mathbf{r})=$$

$$S_A^{-}=$$

$$S_B^z=$$

$$S_B^{+}(\mathbf{r})=$$

$$S_B^{-}=$$

$$H_0\sim$$

$$+\text{const.}$$

$$h_k^0=$$

$$h_k^x=$$

$$h_k^y=$$

$$h_k^z=$$

$$h_{+}=$$

$$h_{-}=$$

$$a_{\mathbf{k}}=$$

$$b_{-\mathbf{k}}^{†}=$$

$$e^{\mathrm{i}{\Phi }_{\mathbf{k}}}=\frac{h_k^x+\mathrm{i}{h}_k^y}{\sqrt{(h_k^x{)}^2+(h_k^y{)}^2}}$$

$$\cosh (2{\Theta }_{\mathbf{k}})=$$

Full text

††preprint: APS/123-QED

I Introduction

Materials subject to an intense incident light shows rich behaviors which are studied in the context of nonlinear response and non-equilibrium phenomena. An example of such is electric shift current in noncentrosymmetric semiconductors and ferroelectrics Belinicher1982 ; Sturman1992 ; Sipe2000 ; Tan2016 ; Morimoto2016 ; Ogawa2017 ; Tokura2018 , where a non-trivial shift of electron position during its optical transition produces a macroscopic electric current. Recent studies revealed that the shift current exhibits strikingly different behaviors from the ordinary photocurrent; the shift current shows unique light-position dependence when it is excited locally Ishizuka2017c ; Nakamura2017 ; Bajpai2019 , and propagates faster than the Fermi velocity of electrons Laman2005 ; Daranciang2012 ; Nakamura2017 . On the other hand, in correlated materials, lower-energy excitations often emerge due to the interaction effect; a typical example is magnetic excitations in Mott insulators. The optical transition of these emergent particles may produce non-trivial phenomena, especially, transport phenomena, related to the nonlinear response of the emergent excitations.

Several recent studies in opto-spintronics and magneto-optics Kirilyuk2010 ; Nemec2018 ; Baltz2018 implies that the intensity and coherence of currently-available electromagnetic waves are sufficient for the control of magnetic excitations or magnetism. Typical results are the following: Magnetization switching by a circularly-polarized laser in ferrimagnets Kimel2005 ; Stanciu2007 ; Khorsand2012 ; Hennecke2019 , laser-driven demagnetization Beaurepaire1996 ; Koopmans2000 ; Stamm2007 , the spin pumping by gigahertz (GHz) or terahertz (THz) waves Kajiwara2010 ; Heinrich2011 , focused-laser driven magnon propagation Satoh2012 ; Hashimoto2017 , intense THz-laser driven magnetic resonance Mukai2016 ; Lu2017 , spin control by THz-laser driven electron transitions Baierl2018 , dichroisms driven by THz vortex beams Cheong2019 , angular momentum transfer between photons and magnons in cavities Haigh2015 ; Osada2016 ; Zhang2016 ; Haigh2016 ; Usami2018 , a ultrafast detection of spin Seebeck effect Kampfrath2018 , a phonon-mediated spin dynamics with THz laser Kampfrath2018b , etc. Moreover, recent theoretical works have proposed several ways of optical control of magnetism: THz-wave driven inverse Faraday effect Takayoshi2014-1 ; Takayoshi2014-2 , Floquet engineering of magnetic states such as chirality ordered states Sato2016 ; Kitamura2017 and a spin liquid state Sato2014 , generation of magnetic defects with laser-driven heat Koshibae2014 ; Fujita2017-1 , applications of topological light waves to magnetism Fujita2017-1 ; Fujita2017-2 ; Fujita2018-1 ; Fujita2018-2 , control of exchange couplings in Mott insulators with high- Mentink2015 and low-frequency Takasan2018 waves, optical control of spin chirality in multiferroic materials Mochizuki2010 , rectification of dc spin currents in magnetic insulators with electromagnetic waves Proskurin2018 ; Ishizuka2019 ; Okuma2019 . These studies are partly supported by recent developments in THz laser science JPhysD2017 ; Cavalleri2017 which realized high-intensity light beams with the photon energy comparable to those of magnetic excitations. Despite these developments, the optical control of the current carried by magnetic excitations is limited to some theoretical proposals.

Among the proposals, a recent theory proposes a mechanism for producing a dc spin current in quantum spin chains without the angular-momentum transfer Ishizuka2019 ; it is distinct from the known mechanisms in which the angular momentum of photons are transferred to the magnet Kajiwara2010 ; Heinrich2011 ; Ohnuma2014 ; Proskurin2018 ; Okuma2019 . The mechanism in Ref. Ishizuka2019, is analogous to that of the shift-current photovoltaic effect Sturman1992 . The close relation between two phenomena are clear from the Jordan-Wigner fermion representation of spin chain; the ground state of the spin chain is a band insulator of Jordan-Wigner fermions, and the photovoltaic response is related to the optical transition of the fermions by the linearly-polarized THz light. However, the relation of this mechanism to the fermion excitations casts doubt on the generality because the low-energy excitations of the ordered magnets are usually magnons, i.e., bosonic excitations.

In this work, we theoretically show that a dc spin current similar to that of the spin chain Ishizuka2019 also appears in ordered antiferromagnetic (AFM) and ferrimagnetic (FRM) insulators by applying a linearly-polarized electromagnetic wave. The symmetry argument in Sec. III shows that the creation of dc spin current with linearly-polarized waves is possible only if both site- and bond-center inversion symmetries are broken. AFM and FRM insulators violate the bond-center inversion symmetry and thereby they naturally satisfy half of the required symmetry condition. The staggered moment is an advantage of considering AFM/FRM insulators for generating a dc spin current. As an example, we consider two-sublattice models with Néel type ground state. Bosonic particles describe the low-energy excitations of these models, i.e., magnons; the ground state is the zero-magnon state. This ground state is very different from that of noncentrosymmetric $S=1/2$ spin chains Ishizuka2019 which are described by a Fermi degenerated state of spinons. Despite the difference, our calculation using a nonlinear response theory finds a finite photovoltaic spin current similar to that of the spinons. We discuss that it is related to the zero-point fluctuation of the quantum magnets. Our theory also indicates that the magnon spin current is shift-current like, i.e., it is insensitive to the magnon lifetime as in the spinon case. This mechanism allows generation of spin current using a linearly-polarized electromagnetic wave and ordinary AFM or FRM insulators.

The remaining part of the paper is organized as follows. In Sec. II, we introduce the nonlinear-response theory for two-species magnons, which we will use in the following sections. The main results of this paper are in Secs. III and IV. Section III focuses on the photo-induced spin current in AFM and FRM insulators with a strong one dimensionality, while we study the three-dimensional (3D) magnets in Sec. IV. Effective experimental setups and signatures for investigating the proposed mechanism are discussed in Sec. V. Section VI is devoted to the summary and discussions.

II Nonlinear response theory

We calculate the nonlinear response coefficients for the photo-induced spin current by extending the linear-response theory to the quadratic order in the perturbation. A similar method for fermions is used to calculate the photovoltaic current in semiconductors Kraut1979 ; vBaltz1981 and the spin current of spinons Ishizuka2019 . The derivation of the formula is summarized in Appendix A. We here summarize the outline of the derivation. We also discuss the physical implications.

We consider a two-sublattice AFM/FRM insulator with two species of magnons. The effective Hamiltonian for the magnons is

[TABLE]

where $\alpha_{\bm{k}}$ ($\alpha_{\bm{k}}^{\dagger}$) and $\beta_{\bm{k}}$ ($\beta_{\bm{k}}^{\dagger}$) are the boson annihilators (creators) for the magnons with the momentum $\bm{k}=(k_{x},k_{y},k_{z})$ and $\varepsilon_{a}(\bm{k})$ ($a=\alpha,\beta$) is the energy of the magnons in the $a=\alpha,\beta$ branch with momentum $\bm{k}$. We here consider a general perturbation (spin-electromagnetic-wave coupling)

[TABLE]

and spin-current operator

[TABLE]

Here, $\omega$ is the frequency of ac light, $h^{\mu}_{\omega}$ is the spin-light coupling constant for the $\mu$ direction, and $\psi_{\bm{k}}=(\alpha_{\bm{k}},\beta_{-\bm{k}}^{\dagger})^{T}$.

The nonlinear conductivity is defined by

[TABLE]

where $\langle J\rangle(\Omega)\equiv\int dt\langle J\rangle(t)e^{-{\rm i}\Omega t}$ is the Fourier transform of the expectation value of the spin current $\langle J\rangle(t)$. For the two-sublattice model, the formula for nonlinear spin current conductivity reads

[TABLE]

where

[TABLE]

The relaxation time of magnons, $\tau_{\bm{k}}$, was introduced in Eq. (10), and $\langle\cdots\rangle_{0}$ is the expectation value of $\cdots$ in the equilibrium state of the Hamiltonian in Eq. (1). The conductivity for dc spin current corresponds to the $\Omega=0$ case, $\sigma_{\mu\nu}(0;\omega,-\omega)$. In the rest of this work, we focus on the case $B_{\bm{k}}^{\mu}=B_{\bm{k}}^{\nu}=B_{\bm{k}}$ because we are interested in the response to a linearly polarized light. Hence, we abbreviate the subscripts in the nonlinear conductivity, $\sigma_{\mu\nu}(0;\omega,-\omega)=\sigma(0;\omega,-\omega)$.

We note that the conductivity in Eq. (10) remains non-zero at $T=0$. The substitutions of $\tilde{\rho}_{\bm{k},\alpha}=0$ and $\tilde{\rho}_{\bm{k},\beta}=1$ in Eq. (10) reduce the formula to

[TABLE]

Because of $\tilde{\rho}_{\bm{k},\beta}=1$, the terms involving the off-diagonal component of $B_{\bm{k}}$ remains at $T=0$. In other words, the two-magnon creation/annihilation process plays a crucial role as shown in Fig. 1(d). We focus on the $T=0$ case in the rest of this paper as this process is dominant in the low temperature limit.

From a different viewpoint, Eq. (18) implies the zero-point fluctuation plays a key role in the photovoltaic response of magnons. In our formalism, the zero-point fluctuation is manifested in the Bogoliubov transformation of Holstein-Primakov bosons. This transformation creates $\beta_{\bm{k}}\beta^{\dagger}_{\bm{k}}$ and $\alpha^{\dagger}_{\bm{k}}\beta^{\dagger}_{-\bm{k}}$ terms which contribute to the photovoltaic response in the ground state. $\tilde{\rho}_{\bm{k},\beta}=1$ is another consequence of the Bogoliubov transformation. The importance of the zero-point fluctuation resembles the spinon spin current Ishizuka2019 , in which the Fermi degeneracy of spinons represents the quantum fluctuation of spins. A crucial difference in the current case is the absence of Fermi degeneracy. However, in the case of the AFMs/FRMs, the condensate of Holstein-Primakov bosons plays a similar role to the Fermi degeneracy. The pair-creation process represented by $\alpha_{\bf k}^{\dagger}\beta_{\bf k}^{\dagger}$ generates photovoltaic response of the magnons which is manifested in the denominator of Eq. (18); the sum of eigenenergies, $\varepsilon_{\alpha}(\bm{k})+\varepsilon_{\beta}(\bm{k})$, represents creation/annihilation of a magnon pair. These features implies that the zero-point fluctuation is necessary for the shift current response at $T=0$.

The first term in Eq. (18) vanishes when the ground state has a certain symmetry. For example, collinear magnetic orders with the moments parallel to $S^{z}$ axis are often symmetric with respect to $G={\cal T}M_{x}^{s}$, which is the product of time-reversal operation (${\cal T}$) and the mirror operation for the spin degrees of freedom about $x$ axis ($M_{x}^{s}$). In this case, the real part of $\sigma(0;\omega,-\omega)$ reads

[TABLE]

The conductivities for the models considered in the following sections are calculated using this formula.

III Spatially-anisotropic magnet

In this section, we apply the above formula to a spin chain with AFM or FRM order, which corresponds to a quasi-one-dimensional (quasi-1D) magnetic compound with a negligible inter-chain interaction. The spins are coupled to the electromagnetic wave through the Zeeman coupling. To make the problem theoretically well-defined, we consider a model which conserves the spin angular momentum $S^{z}$; the model has an easy axis and the applied ac magnetic field is parallel to the ordered moments. The conservation of $S^{z}$ allows us to unambiguously define the spin current operator from the continuity equation. This setup is in contrast to those of usual magnetic resonances and spin pumping Kajiwara2010 ; Heinrich2011 , in which the ac field is perpendicular to the magnetic moment. We use the standard spin-wave approximation to describe magnetic excitations (magnons).

III.1 Model

We consider an ordered noncentrosymmetric spin chain with a two-sublattice unit cell [Fig. 1(a)], whose Hamiltonian is given by

[TABLE]

where $H_{\rm 1D}$ is the spin-chain Hamiltonian with the staggered nearest-neighbor exchange interaction (i.e., dimerization) along the $x$ direction, $H_{0}$ is the bundle of all the chains, and $H_{Z}^{(\omega)}$ is the Zeeman coupling between the spins and the external electromagnetic wave. Here, $\bm{S}_{a}(\bm{r})\equiv(S_{a}^{x}(\bm{r}),S_{a}^{y}(\bm{r}),S_{a}^{z}(\bm{r}))$ ($a=A,B$) is the spin-$S_{a}$ operator on the $a$ sublattice of the unit cell at position $\bm{r}=(r_{x},r_{y},r_{z})$. Symbols $\hat{x}$, $\hat{y}$, and $\hat{z}$ stand for the unit vectors along the $x$, $y$, and $z$ directions, respectively. The parameters in the Hamiltonian $H_{\rm 1D}$ are as follows: $J>0$ is the antiferromagnetic exchange interaction along the spin-chain ($x$) direction, $\delta$ is the dimerization, $D>0$ ($D_{s}$) is the uniform easy-axis (staggered) anisotropy, $g_{A}$ ($g_{B}$) is the $g$ factor for the spins on $A$ ($B$) sublattice, and $h$ is the external static magnetic field along the $S^{z}$ axis. In the spin-light coupling $H_{Z}^{(\omega)}$, $|h_{\omega}|$ and $\arg(h_{\omega})$ are respectively the magnitude and the phase of the ac magnetic field of the linearly-polarized electromagnetic wave. We assume $|D_{s}|<D$, and $|\delta|<1$.

When $S_{A}\neq S_{B}$, the ground state of the model in Eq. (20c) is a FRM-ordered state with magnetization $|S_{A}-S_{B}|$ per a unit cell Marshall1955 ; Lieb1962 . The ground state is a Néel ordered when $S_{A}=S_{B}$. The classical ground state of $H_{0}$ has a collinear order with spins pointing along the $S^{z}$ axis because of the easy axis anisotropy $D$ [Fig. 1(a)]. The anisotropy also produces the spin gap in the excitation spectrum [Fig. 2(a) and 1(b)]. We discuss the effect of the gap and its relation to the frequency dependence of the nonlinear spin conductivity in the next section.

Here we define the spin current for $S^{z}$. Since the model $H_{0}$ conserves the $z$ component of total spin angular momentum, the spin current for $S^{z}$ can be defined from the continuity equations $\partial_{t}S_{A}^{z}=J_{x}^{z}(r_{x}-1,B;r_{x},A)-J_{z}^{z}(r_{x},A;r_{x},B)$ and $\partial_{t}S_{B}^{z}=J_{x}^{z}(r_{x},A;r_{x},B)-J_{z}^{z}(r_{x},B;r_{x}+1,A)$, in which $J_{\beta}^{\alpha}(r,a;r^{\prime},b)$ is the local spin-$S^{\alpha}$ current operator between two neighboring sites $(r,a)$ and $(r^{\prime},b)$ and it flows along the $\beta$ direction. The above continuity equation is obtained from Heisenberg equation of motion for local spins. With these procedures, we find the uniform current operator for $H_{\rm 1D}$ reads

[TABLE]

where $N$ is the total number of unit cells.

III.2 Linear spin-wave approximation

Hereafter, we assume that in the ground state of $H_{\rm 1D}$, the spins on the $A$ sublattice points up while those on $B$ sublattice are down [see Fig. 1 (a)]. The low-energy excitations of $H_{0}$ is calculated by linear spin-wave approximation. Using the Holstein-Primakov bosons, the spin operators are given by

[TABLE]

for the $A$ sublattice and

[TABLE]

for the $B$ sublattice. Up to the linear order in $S_{A}$ and $S_{B}$, $H_{0}$ reads

[TABLE]

where the wave number along the chain ($x$) direction is simply represented by $k$, $a_{\bm{k}}\equiv(1/\sqrt{N})\sum_{\bm{r}}a(\bm{r})e^{{\rm i}\bm{k}\cdot\bm{r}}$, $b_{\bm{k}}\equiv(1/\sqrt{N})\sum_{\bm{r}}b(\bm{r})e^{{\rm i}\bm{k}\cdot(\bm{r}+\hat{x}/2)}$ are the Fourier transformation of Holstein-Primakov bosons. The matrix elements of the magnon Hamiltonian (31) are calculated as

[TABLE]

where

[TABLE]

The quadratic Hamiltonian (31) is diagonalized by the Bogoliubov transformation:

[TABLE]

where $\alpha_{\bm{k}}$ ($\alpha_{\bm{k}}^{\dagger}$) and $\beta_{\bm{k}}$ ($\beta_{\bm{k}}^{\dagger}$) are bosonic annihilation (creation) operators. By choosing

[TABLE]

the Hamiltonian becomes

[TABLE]

where

[TABLE]

Here, we ignored the constant term in $H_{0}$. We note that the dispersions $\varepsilon_{\alpha,\beta}(\bm{k})$ and the phases $(\Theta_{\bm{k}},\Phi_{\bm{k}})$ are all independent of $k_{y}$ and $k_{z}$ because we now consider the 1D model $H_{0}$. Using the same transformation, we find

[TABLE]

and

[TABLE]

III.3 Spin current conductivity

Combining the magnon representation of $(\alpha_{\bm{k}},\beta_{\bm{k}})$ with the formula (19), we compute the nonlinear dc spin-current conductivity for the model $H_{\rm tot}$ under the application of THz laser. We first study the nonlinear conductivity in the clean limit with infinite relaxation time $\tau_{k}\to\infty$. The analytic solution for the conductivity ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]$ obtained from Eq. (19) reads

[TABLE]

when $\omega\in[\omega_{c1},\omega_{c2}]$ and zero otherwise. Here,

[TABLE]

corresponds to the energy for the band bottom of the pair excitation and

[TABLE]

is that for the top of the pair excitation [See Fig. 1(c) and (d)]. The frequency dependence of the conductivity is shown in Fig. 2(a).

Equation (41) is an odd function of $\delta$. This reflects the fact that the inversion-symmetry breaking is necessary for the spin current. $H_{0}$ has two inversion centers when $\delta=0$, $D_{s}=0$, $g_{A}=g_{B}$, and $S_{A}=S_{B}$: one at the center of the bond and the another on the site. The inversion center on the site is broken by the dimerization $\delta$. To see the dependence of $\sigma(0;\omega,-\omega)$ on the model parameters, we explicitly write down the nonlinear conductivity as a function of the parameters, i.e., $\sigma(0;\omega,-\omega)=\sigma(\omega;\delta,g_{A}-g_{B},D_{s},m)$, where $m=\langle S^{z}_{\bm{r}\in A}\rangle-\langle S^{z}_{\bm{r}\in B}\rangle$ is the order parameter of the AFM or FRM insulators. A symmetry argument on the transport coefficient finds $\sigma(\omega;\delta,g_{A}-g_{B},D_{s},m)=-\sigma(\omega;-\delta,g_{A}-g_{B},D_{s},m)$ for the site-center inversion operation. This result is identical to the spinon case in Ref. Ishizuka2019 .

On the other hand, the magnetic order changes the parameter dependence of $\sigma(0;\omega,-\omega)$, which is related to the bond-center inversion operation. The inversion operation about the center of the bonds is broken by the Néel ordering, $D_{s}\neq 0$ or $g_{A}\neq g_{B}$. Therefore, the symmetry operation indicates $\sigma(\omega;\delta,g_{A}-g_{B},D_{s},m)=-\sigma(\omega;\delta,-g_{A}+g_{B},-D_{s},-m)$. In addition, the translation operation about half a unit cell switches $A$ and $B$ sublattices and $m\to-m$; $\sigma(\omega;\delta,g_{A}-g_{B},D_{s},m)=-\sigma(\omega;\delta,g_{A}-g_{B},D_{s},m)$. Hence, the conductivity in the ordered phase is an even function of $g_{A}-g_{B}$ and $D_{s}$. This is a different behavior from the spinon case, in which the conductivity is an odd function of the staggered magnetic field (corresponds to $g_{A}-g_{B}$ in our case).

The conductivity diverges when $\omega$ approaches $\omega_{c1}$. The asymptotic form reads

[TABLE]

where $\delta\omega\equiv\omega-\omega_{c1}$. A similar feature is also found in the spinon case, in which the divergence is related to the singularity of the density of states Ishizuka2019 . On the other hand, the asymptotic form around $\omega=\omega_{c2}$ reads

[TABLE]

The sign of the conductivity is the opposite of that in the lower frequency regime. This is in contrast to the spinon case Ishizuka2019 , in which the sign of the nonlinear conductivity remains the same for all frequencies $\omega\in[\omega_{c1},\omega_{c2}]$.

III.4 Relaxation-time dependence

We next study the damping (relaxation time) dependence of the spin current. In the study of photovoltaic effect, the relaxation-time dependence reflects the microscopic mechanism behind the photovoltaic effect Sturman1992 ; it is called shift current when $\sigma(0;\omega,-\omega)\propto\tau^{0}$ while is an injection current when $\sigma(0;\omega,-\omega)\propto\tau$. In bosonic systems, a slight difference appears in the momentum dependence of the single-particle relaxation time Abrikosov1965 ; it is inversely proportional to the momentum for the Goldstone modes. Therefore, we assume the momentum dependence of damping term as $\tau_{\bm{k}}=1/(\alpha_{0}\varepsilon_{\beta}(\bm{k}))$ so that the momentum dependence is consistent with the field theoretic requirement ($\alpha_{0}$ is the damping factor). Physically, the assumed form of $\tau_{\bm{k}}$ corresponds to the phenomenological Gilbert damping.

We substitute $\tau_{\bm{k}}=1/(\alpha_{0}\varepsilon_{\beta}(\bm{k}))$ in Eq. (19) in order to estimate the relaxation-time dependence of the spin-current conductivity. Figure 2(b) shows the $\alpha_{0}$ dependence of $\sigma(0;\omega,-\omega)$. Our numerical result shows $\sigma(0;\omega,-\omega)$ is insensitive to the damping. A slight difference, however, appears in the high-frequency region, where the smearing due to the damping is more distinct than that in the low-frequency region. This behavior is related to the momentum dependence of $\tau_{k}$, which is inversely proportional to the energy of the magnon. The insensitivity is a signature of a shift-current type photo-induced current Sturman1992 ; this is a similar feature to the spinon case Ishizuka2019 .

IV Three-dimensional magnetic insulators

In this section, we consider a three-dimensional (3D) magnet which consists of coupled spin chains $H_{\rm 1D}$ with a non-negligible inter-chain interaction [See Fig. 1(b)]. We particularly focus on the limit in which $\omega$ is close to the band gap of two-magnon excitations. The procedure of the calculation is the same as the 1D case in the previous section. The static part of the Hamiltonian reads

[TABLE]

The spin chains is parallel to the $x$ direction, while the $y$ and $z$ directions are perpendicular to the chains. The ferromagnetic coupling $J_{\perp}>0$ denotes the strength of the inter-chain exchange interaction. We study this model within the linear spin-wave approximation using Holsten-Primakov transformation in Sec. III.2. Focusing on the lower edge of the magnon dispersion, we first expand the matrix elements $h_{\bm{k}}^{a}$ of the magnon Hamiltonian [See Eq. (32)] up to second order in $\bm{k}$:

[TABLE]

We note that the magnon dispersions depend on both intra- and inter-chain wave numbers differently from the 1D case. The dispersion around $\Gamma$ point $\bm{k}=\bm{0}$ is shown in Fig. 3. Using the momentum gradient of the low-energy Hamiltonian with $h_{\bm{k}}^{0,x,y,z}$, we can define the spin current operator; this approximation is essentially equivalent to expanding the lattice spin-current operator in Eq. (40) up to the linear order in $\bm{k}$:

[TABLE]

These equations corresponds to the $k\cdot p$ expansion of the lattice model. Therefore, it should be a good approximation for the lattice model when $\omega$ is close to the gap for two-magnon excitations.

The spin current conductivity is calculated using the formula of Eq. (19). A calculation similar to the 1D model considered in Sec. III gives

[TABLE]

where

[TABLE]

When $\omega$ is close to the lower edge, i.e.,

[TABLE]

$K_{X}$ becomes

[TABLE]

where $\delta\omega=\omega-\omega_{c1}$. Therefore, the asymptotic form of ${\rm Re}\left[\sigma_{ABB}(0;\omega,-\omega)\right]$ is

[TABLE]

Unlike the 1D case, in which the conductivity diverges at the band edge $\omega_{c1}$, the 3D result in Eq. (53) decreases proportionally to $\sqrt{\delta\omega}$ when approaching $\omega_{c1}$. The result is plotted in the inset of Fig. 2(a) with the results for the 1D limit. This difference is a consequence of the difference in the density of states: it diverges in the 1D model while it is proportional to $\sqrt{\delta\omega}$ in the present 3D case.

The approximation we used in this section is accurate when $\omega$ is close to the magnon gap at the $\Gamma$ point in the Brillouin zone. In our model, the band bottom for the two-magnon excitations are at $\Gamma$ point, and the bandwidth of two-magnon excitation along the $x$ and $y$ directions are in the order of $J_{\perp}$ and that for $z$ direction is in the order of $J$. Therefore, our approximation is accurate when $\delta\omega\ll J,J_{\perp}$. This condition is manifested in $J_{\perp}$ in the denominator of Eq. (53), which implies the divergence of ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]$ at $J_{\perp}\to 0$. When $J_{\perp}$ is very small, we expect ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]$ to behave like that of the 1D case. On the other hand, ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]$ looks like Eq. (53) when $J_{\perp}$ is sufficiently large, e.g., when $J_{\perp}\sim J$. Therefore, the 1d result and the result in this section corresponds to the two limits of the 3D magnet.

V Experimental observation

In this section, we discuss experimental methods for detecting signatures of a directional spin current in our mechanism.

V.1 Setup

We here discuss experimental setups for the observation of the spin current generated by linearly-polarized light. The mechanism studied here produces a directional flow of the spin current, which is a distinct feature from the spin pumping Kajiwara2010 ; Heinrich2011 . Therefore, the observation of the directional flow should provide an evidence for our mechanism. We discuss two different mechanisms: First one is an all-optical setup using Kerr rotation or Faraday effect, and the second is a two-terminal setup using inverse spin-Hall effect.

Observation of the spatial distribution of angular momentum in the open circuit setup provides a direct evidence for the optically-generated spin current [See Fig. 4(a)]. In an isolated magnet, the spin current produced by a THz light flows along a direction defined by the magnetic order and the crystal symmetry. Therefore, if the system becomes close enough to a laser-driven non-equilibrium steady state, the angular momentum accumulates at the two ends in an open circuit setup in Fig. 4(a); positive angular momentum on one end and negative on the other end. The angular momentum distribution is anti-symmetric along the direction of the spin current. This distribution is strikingly different from the spin pumping case in which the distribution is symmetric and its difference from the equilibrium state is larger at the focal area of the laser than at the ends.

An all-optical setup using Kerr rotation or Faraday effect would be a useful setup for the observation of such a spatial distribution. Measurement of magnetic moments and its spatial distribution using the optical probe is a commonly used technique for observing the spin current. For instance, this method is used to observe the spin Hall effect Kato2004 . Similarly, observing the magnetization of soft magnet layers attached to the two ends is another possible setup for the experiment [Fig. 4(b)].

The observation of spin current in a two-terminal setup in Fig. 4(c) also enables us to see the directional flow of spin current and to distinguish it from the spin pumping effect. This setup consists of a noncentrosymmetric magnetic insulator which is sandwiched between two metallic leads; the two leads detect spin current via inverse spin Hall effect Saitoh2006 ; Valenzuela2006 ; Kimura2007 . In the photovoltaic mechanism, the spin current in the two leads flows toward the same direction. Therefore, the inverse spin Hall voltage of the two leads has the same sign. In contrast, in the spin pumping, the spin current diffusively flows outward from the magnet; the inverse spin Hall voltage is positive on one side and negative on the other. Therefore, the relative sign of the inverse spin Hall voltage of the two leads can make a distinction between the spin pump and our mechanism.

Finally, we shortly comments on heating effect of applied electromagnetic waves. When we try to detect the photovoltaic spin current with the above setups, spin pumping might also occur due to the heating effect of the applied laser. For such a case, extracting the asymmetric part of the angular-momentum distribution or inverse spin Hall voltage is important to detect an evidence for our mechanism.

V.2 Required intensity of AC field

We next estimate the required ac electromagnetic field for generating an observable spin current. We here assume a spin current of $J_{s}=10^{-16}$ J/cm2 is observable. This estimate is based on a Boltzmann theory calculation for spin Seebeck effect in a ferromagnet Hirobe2017 ; Ishizuka2019 . The details of the estimate is briefly explained in Appendix B. We use the following parameters as a typical value for 1D insulating magnets: $J=100k_{B}$ J, $\delta=0.1$, $S_{A}=S_{B}=1$, $g_{A}-g_{B}=0.1\mu_{B}$ J/T, $h_{+}=10k_{B}$ J, and a the light with a frequency which is $\hbar\delta\omega=6\pi\hbar\times 10^{11}$ Hz above the band gap. Here, $\hbar$ is the Planck constant. With these parameters, the conductivity for the 1D AFM/FRM chain is ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]\sim 10^{-14}$ J/(cm2T2). Therefore, the required magnitude of oscillating magnetic field to produce a spin current of $J_{s}=10^{-16}$ J/cm2 is $B\sim\sqrt{\frac{J_{s}}{|{\rm Re}\left[\sigma(0;\omega,-\omega)\right]|}}\sim 0.1$ T. This corresponds to the electric field $E=cB\sim 10^{4}-10^{5}$ V/cm under the assumption of $c=10^{8}$ m/s which is a typical value of speed of light in insulators. Similar estimate for the 3D magnet with $J=100k_{B}$ J, $J_{\perp}=10k_{B}$ J, $\delta=0.1$, $S_{A}=S_{B}=1$, $g_{A}-g_{B}=0.1\mu_{B}$ J/T, $h_{+}=10k_{B}$ J, and $\omega=2\pi\times 10^{12}$ Hz gives ${\rm Re}\left[\sigma(0;\omega,-\omega)\right]\sim 10^{-11}$ J/(cm2T2) and $E=cB\sim 10^{5}-10^{6}$ V/cm. Our estimate predicts that the photovoltaic spin current is experimentally observable by using a moderate-intensity THz light.

V.3 Candidate material

We believe the photovoltaic spin current should be seen generically in noncentrosymmetric magnets. In a recent work Ishizuka2019 , the authors find three kinds of spin-light couplings induce the spin current in a spin chain, and this work presents photovoltaic spin current in ordered magnets. These results imply the generation of photovoltaic spin current is a universal phenomenon in noncentrosymmetric magnetic insulators. One such material is ferrimagnetic diamond chains Okamoto2003 ; Shores2005 ; Vilminot2006 . These materials often have a distortion associated with trimerization, which breaks the inversion symmetry Okamoto2003 . Also, a large density of states for the magnon excitations is expected in this material because it is a quasi-1D magnet. Thus the ferrimagnetic phase of the diamond chain is a promising candidate for studying the spin current.

VI Summary and Discussion

To summarize, we studied the spin current generation through the shift current mechanism in ferrimagnetic/antiferromagnetic insulators. Our theory uses a nonlinear response theory, which is a natural generalization of the linear response theory. Based on this method, we find that the illumination of a linearly-polarized light produces the magnon current in noncentrosymmetric magnets with antiferromagnetic/ferrimagnetic order. The photovoltaic spin current appears even at the zero temperature where no magnon excitation exists in the equilibrium; the current is related to the two-magnon excitation process and not to the optical transition of existing (thermally-excited) magnons. We stress that the photo-induced spin current in our mechanism is carried by electrically-neutral particles. The relaxation-time dependence of the spin current indicates that our photovoltaic effect is a “shift current”, i.e., the nonlinear conductivity is insensitive to the damping. Our theory clearly shows that the shift current mechanism, which is well known in electron (fermion) systems, is also relevant to systems with bosonic excitations whose the ground state is the vacuum of bosons (zero boson state).

Our result implies the zero-point quantum fluctuation is a key for the shift-current type photocurrent. In the spinon spin current Ishizuka2019 , the optical transition of a fermionic excitation plays a crucial role for the photocurrent. In contrast to these cases, the ground state of the ordered magnets is the zero-magnon state. Therefore, no optical transition of the existing magnons. Despite the crucial difference, we find a finite photovoltaic spin current at the zero temperature. The magnon photocurrent we found is ascribed to the optical transition of the “condensed” Holstein-Primakov bosons. In the antiferromagnets/ferrimagnets, the ground state is a condensate of Holstein-Primakov bosons, which is technically represented by the Bogoliubov transformation. The optical transition of the condensed Holstein-Primakov bosons allows generation of the shift-current type photocurrent even at the zero temperature. On the other hand, we find that the nonlinear conductivity is zero at $T=0$ for the ferromagnetic version of the model considered here. From this viewpoint, the two-magnon creation is similar to the particle-hole pair creation in semiconductors; the optical transition of fermions from the valence band to the conduction band is equivalent to the pair creation. As the condensation of the Holstein-Primakov bosons is a manifestation of zero-point fluctuation, the zero-point fluctuation is the essence for the shift-current type photovoltaic effects in the magnetic insulators.

Our results implies that the dc spin current generation using linearly polarized light is generally possible in the magnets without inversion symmetry.

Appendix A Derivation of Kraut-von Baltz formula for Bosons

Here, we shortly explain the derivation of the nonlinear conductivity in two-band boson systems. We used the formula in Eq. (75) for the analytic calculations and Eq. (66) for numerical results with a finite Gilbert damping.

We calculate the nonlinear response coefficients using a formalism similar to the linear response theory. We assume a system with a time-dependent perturbation $H^{\prime}=-\sum_{\mu}\hat{B}_{\mu}F_{\mu}(t)$, where $\hat{B}_{\mu}$ is an operator and $F_{\mu}(t)$ is a time-dependent field; the Hamiltonian reads $H=H_{0}+H^{\prime}$. The expectation value of an observable $\hat{A}$ reads $\langle\hat{A}\rangle(t)={\rm Tr}\left[\hat{\rho}(t)\hat{A}\right]/Z,$ where $\rho(t)$ is the density matrix at time $t$ and $Z\equiv{\rm Tr}\rho(t)$. By expanding $\rho(t)$ up to the second order in $F_{\mu}(t)$, the Fourier transform of $\langle A\rangle(t)$, $\langle A\rangle(\Omega)$, reads

[TABLE]

with the nonlinear conductivity

[TABLE]

Here, $E_{n}$ is the eigenenergy of the many-body eigenstate $n$, $\tau_{mn}$ is the relaxation time, and $O_{nm}$ ($O=A,B_{\mu},B_{\nu}$) is the matrix element of $\hat{O}$ in the eigenstate basis of $H_{0}$.

We here consider a periodic free-boson system in which all matrices $A$, $B^{\mu}$, and $B^{\nu}$ have the following form:

[TABLE]

where $\alpha_{\bm{k}}$ ($\alpha_{\bm{k}}^{\dagger}$) and $\beta_{\bm{k}}$ ($\beta_{\bm{k}}^{\dagger}$) are the annihilation (creation) operators of the boson eigenstates with momentum $\bm{k}$, and $O_{\bm{k}}=A_{\bm{k}},B^{\mu}_{\bm{k}},B^{\nu}_{\bm{k}}$. The theory for spinwave excitations of many antiferromagnetic models with a Néel-type order reduces to the above form by using Holstein-Primakov and Bogoliubov transformations.

For the two-band system, we can express Eq. (55) using single-particle eigenstates. We note that $A$, $B_{\mu}$, and $B_{\nu}$ for the two-band system above do not conserve the particle number. However, all operators are quadratic in the annihilation/creation operators and consists of only for terms: $\alpha_{\bm{k}}^{\dagger}\alpha_{\bm{k}}$, $\beta_{-\bm{k}}\beta_{-\bm{k}}^{\dagger}$, $\beta_{-\bm{k}}\alpha_{\bm{k}}$, and $\alpha_{\bm{k}}^{\dagger}\beta_{-\bm{k}}^{\dagger}$. Therefore, only few terms out of the possible Wick decomposition remain nonzero, similar to that of the systems with conserved particle number. Using these features, we find

[TABLE]

Here,

[TABLE]

and we assumed the relaxation time only depends on $\bm{k}$. It is worth noting that the conductivity remains finite at $T=0$ despite there are no excitations. Technically, this is a consequence of $\tilde{\rho}_{\bm{k},\beta}$, which is 1 at $T=0$. Physically, this is because the pair creation/annihilation processes contribute to the spin current even at $T=0$.

We here focus on the $T=0$ limit. In this limit, $\tilde{\rho}_{\bm{k},\alpha}=0$ and $\tilde{\rho}_{\bm{k},\beta}=1$. Using these results, we obtain

[TABLE]

As we discussed in the main text, certain symmetries restricts the first term to be zero; this is the case for the models we consider in the main text. Assuming the first term vanishes, we find

[TABLE]

We used this formula for the calculation of nonlinear conductivity in the main text.

Appendix B Boltzmann theory for spin Seebeck effect

The magnitude of spin current $J_{s}=10^{-16}$ J/cm2 is the estimate for the spinon spin current produced by the spin Seebeck effect in a recent experiment Hirobe2017 . We here summarize the method and result discussed in a supplemental material of a recent work Ishizuka2019 .

The spin current is estimated from the Seebeck effect of magnons whose dispersion is given by

[TABLE]

near the $\Gamma$ point of $\bm{k}=\bm{0}$. This magnon dispersion corresponds to that of a ferromagnetic heisenberg model with exchange interaction $J$, uniaxial anisotropy $D$, and the magnetic field $h$ parallel to the anisotropy. The current is calculated using the semiclassical Boltzmann theory, in which the current reads

[TABLE]

Here, $f_{\bm{k}}(\bm{r})$ is the density of magnons with momentum $\bm{k}$ at position $\bm{r}$ and $v_{z}\equiv\partial_{k_{z}}\varepsilon(\bm{k})$ is the group velocity of magnons. $f_{\bm{k}}(\bm{r})$ is calculated from the Boltzmann equation with temperature gradient

[TABLE]

where $f^{(0)}_{\bm{k}}(\bm{r})$ is the density at the equilibrium. Here, the relaxation-time approximation is used to simplify the calculation of collision integral on the right hand side. The spin current induced by the spin Seebeck effect is estimated by substituting the solution of $f_{\bm{k}}(\bm{r})$ in Eq. (78) into the current formula in Eq. (77)

In the Boltzmann theory, the spin current by the spin Seebeck effect reads

[TABLE]

where $\Lambda=(6\pi^{2})^{1/3}/a$ is the cutoff for magnon dispersion and

[TABLE]

Using a set of typical parameters $S=1$, $J_{H}=100k_{B}$ J, $D=0$ J, $h=\mu_{B}$ J, $a=4\times 10^{-10}$ m, $\alpha=10^{-2}$, $T=100$ K, $\Delta T=3\times 10^{4}$ K/m, we find $J_{s}\sim 10^{-12}$ J/cm2 for the ferromagnet. We assume this value as the typical spin current density in the insulating ferromagnets.

A recent experiment on quasi-one-dimensional magnets observed a spin current which is 10*-4* of what is typically observed in a ferromagnetic phase Hirobe2017 . Therefore, we assume $J_{s}\sim 10^{-16}$ J/cm2 as the experimental resolution for the spin current.

Acknowledgements.

We thank Ryosuke Matsunaga and Youtarou Takahashi for fruitful discussions. We also thank Wataru Murata for providing Fig. 4. H.I. was supported by JSPS KAKENHI Grant Numbers JP18H04222, JP19K14649, and JP18H03676, and CREST JST Grant Numbers JPMJCR16F1. M.S. was supported by JSPS KAKENHI (JP17K05513), and Grant-in-Aid for Scientific Research on Innovative Area “Nano Spin Conversion Science” (Grant No.17H05174) and “Physical Properties of Quantum Liquid Crystals” (Grant No. 19H05825).