Magnon-mediated Dzyaloshinskii-Moriya torques and pumping of magnons Alexey A. Kovalev, Vladimir A. Zyuzin, and Bo Li Department of Physics and Astronomy and Nebraska Center for Materials and Nanoscience, University of Nebraska, Lincoln, Nebraska 68588, USA (Dated: January 25, 2017) We formulate a microscopic linear response theory of nonequilibrium magnonic torques and magnon pumping applicable to multiple-magnonic-band uniform ferromagnets with Dzyaloshin- skii–Moriya interactions. From the linear response theory, we identify the extrinsic and intrinsic contributionswherethelatterisexpressedviatheBerrycurvatureofmagnonicbands. Weobserve thatinthepresenceofatime-dependentmagnetizationDzyaloshinskii–Moriyainteractionscanact 7 as fictitious electric fields acting on magnons. We study various current responses to this fictitious 1 field and analyze the role of Berry curvature. After identifying the magnon-mediated contribution 0 to the equilibrium Dzyaloshinskii–Moriya interaction, we also establish the Onsager reciprocity be- 2 tweenthemagnon-mediatedtorquesandheatpumping. Weapplyourtheorytothemagnonicheat n pumping and torque responses in honeycomb and kagome lattice ferromagnets. a J PACSnumbers: 85.75.-d,72.20.Pa,75.30.Ds,72.20.My 4 2 I. INTRODUCTION anomalousresponsesinducedbystatisticalforces[29–31]. ] l There is a considerable interest in magnets on lat- l a Itiswellknownthatanelectricfieldcandriveacharge ticeswithnon-trivialgeometryastheyallowobservation h current, whereas in order to understand how to drive a of Berry phase related phenomena such as the thermal s- spincurrentonemightneedtoresorttothefieldofspin- Hall effect of magnons [32–41]. Theoretically, the in- e tronics [1]. Magnetization dynamics generates spin cur- creased magnon damping [40], Dirac magnons [42], and m rents in adjacent normal metal by a phenomenon known the magnon-mediated spin Hall effect [16; 43; 44] have . as spin pumping [2–4]. The discovery of spin pumping been predicted for kagome and honeycomb lattice ferro- t a had a great deal of influence on the development of the magnets. In addition, other manifestations of the Berry m field of spintronics as it led to new insights into the spin phase physics can arise in layered kagome [40] and hon- - Hall [5], spin torque [6; 7], and spin Seebeck effects [8]. eycomb [45] ferromagnets as examined in this work. d The phenomena related to the spin Seebeck effect are In this work, we analyze magnon currents arising in n studied within the field of spincaloritronics [9] in which response to magnetization dynamics (see Fig. 1). In the o c thefocusisoninterplaybetweenthespindegreesoffree- presence of a time-dependent magnetization, DMI can [ dom and heat currents. actasfictitiouselectricfieldsactingonmagnons. Ashas 1 Asheatandspincurrentsarealsocarriedbymagnons, been noted earlier in the introduction, the energy cur- v one naturally arrives at a concept of magnon-mediated rent carried by such magnons contains the ground state 0 spin torques which can lead to thermally induced mo- contribution associated with magnon-mediated equilib- 9 tion of magnetic domain walls [10–12]. Such torques ex- riumDMI.Notethatsuchcorrectionsareimportantonly 9 ist only in noncollinear magnetic structures or when the in systems with non-trivial Berry curvature of magnon 6 Dzyaloshinskii-Moriya interactions (DMI) are present. bands. Here, we concentrate on systems with non-trivial 0 In the latter case, such spin torques have been termed Berrycurvaturebyconsideringvariouscurrentresponses . 1 as DMI torques [13]. Recently, both field-like and in honeycomb and kagome lattice ferromagnets. Our 0 antidamping-likecontributionstoDMItorqueshavebeen linear response calculation of heat currents agrees with 7 studied theoretically [14–18]. It has been noted [13] that thecalculationofmagnon-mediatedthermaltorques[16], 1 : DMItorquescanbeseenasmagnonanalogsofspin-orbit thus confirming the Onsager reciprocity principle (see v torques [19–24]. This suggests that the phenomenology Fig. 1). We also study the feasibility of experimental Xi developed for spin-orbit torques can be readily applied observation of such current responses. to DMI torques [25; 26]. In particular, the intrinsic con- The paper is organized as follows. In section II, we in- r a tribution to DMI torques has been identified [16]. Con- troduce the Hamiltonian describing magnons with mul- tinuing this analogy, one can identify fictitious electric tiple bands and calculate the equilibrium DMI. Next, fields acting on magnons due to time-dependent magne- withinthesamesection,wedescribepumpingofmagnons tization dynamics [11; 27; 28]. One can also identify the in response to magnetization dynamics and thermal magnon-mediated equilibrium contribution to DMI. Due torques within the linear response theory. In the final to such contribution the electron-mediated energy cur- part of section II, we formulate the Onsager relations. rent calculated in response to magnetization dynamics In section III, we apply our theory to honeycomb and fromtheKuboformalismcontainsanunphysicalground- kagome lattice ferromagnets. We conclude our paper in state contribution [26] which needs to be subtracted. section IV. The Appendices A, B, C, and D contain very Similar unphysical contributions have been identified for detailed derivations of our results. 2 of magnons in equilibrium state. The torque operator is introduced as 𝒎 𝒎 z y T =∂ H×m, (3) m x 𝑱 𝑻 𝑻+∆𝑻 where m is a unit vector in the direction of the spin 𝒑 density. WetheninterpretDMIintermsofthemoments of the torque: Figure 1. (Color online) Two effects related by the Onsager ˆ reciprocity principle. Left: Magnetization dynamics pumps 1(cid:28) (cid:29) magnon current Jp and spin current Js = −(cid:126)Jp. This pro- Dαβ = 2 drΨ†(r)(Tαxβ +xβTα)Ψ(r) , (4) cessalsoinvolvesheatcurrentJ carriedbymagnons. Right: eq q A temperature gradient leads to a thermal torque with two whereweassumeafinitesystem. Inordertorepresentan components T and T acting on the uniform magnetization. x y infinite system, we will eliminate the position operator from the final result. The average in Eq. (4) has been calculated in Ref. [16] in a form of a tensor M defined II. THEORY OF MAGNON PUMPING AND β as DMI TORQUES 1 M = Tr[(x ∂ H +∂ Hx )g(E)], (5) In this section, we develop a microscopic linear re- β 2 β m m β sponse theory of magnon pumping and nonequilibrium where g(E) is the Bose distribution function g(E) = magnonic torques applicable to multiple-magnonic-band 1/[exp(βE)−1]. In particular, it has been found that uniform ferromagnets with Dzyaloshinskii–Moriya inter- actions. Tosimplifyformulas,wetakethesystemvolume (cid:26) (cid:27) V = 1 and recover it in the final expressions (19), (20), Mβ =(cid:88) β1 ln(1−e−βεnk)Bm(nβ)(k)−g(Ekn)A(mnβ)(k) , and (28). kn (6) where A. Preliminaries (cid:20) (cid:21) A(n)(k)= (cid:88) Im [η ] 1 [v ] , (7) mβ (cid:101)k nmε −ε (cid:101)kβ mn We consider a noninteracting boson Hamiltonian de- m(cid:54)=n nk mk scribing the magnon fields, which could be, e.g., a result and of the Holstein-Primakoff transformation: ˆ (cid:20) (cid:21) H= drΨ†(r)HΨ(r), (1) Bm(nβ)(k)= (cid:88) Im [η(cid:101)k]nm(ε −2ε )2[v(cid:101)kβ]mn , (8) nk mk m(cid:54)=n where H is a Hermitian matrix of the size N ×N and with the velocity v = ∂ H , the effective field η = Ψ†(r)=[a†1(r),...,a†N(r)]describesN bosonicfieldscor- −∂mHk, and theirkeigenkbaksis representations, v(cid:101)kk = responding to the number of modes within a unit cell T†v T and η = T†η T . Finally, the expression for (or the number of spin-wave bands). The Fourier trans- k k k (cid:101)k k k k the DMI tensor is given by formed Hamiltonian reads D =[M ×m] . (9) H=(cid:88)a†H(k)a , (2) αβ β α k k k It is easy to notice that B(n)(k)=−Ω(n)(k) where now mβ mβ (cid:104)(cid:16) (cid:17) (cid:105) where a†k is the Fourier transformed vector of creation Ω(mnβ)(k)≡i ∂mTk† (∂βTk) nn−(m↔β)isthemixed operators. Hamiltonian in Eq. (2) can be diagonalized spaceBerrycurvatureofthenthband. Thesecondterm bya unitarymatrixT , i.e. E =T†H(k)T and T†T = in Eq. (6) has a clear analogy to the orbital moment [46] k k k k k k 1 where E is the diagonal matrix of band energies, which can be seen after a substitution η →v [25]. N×N k k k and 1 is the N ×N unit matrix. N×N C. Heat and spin pumping by magnetization B. Magnon-mediated Dzyaloshinskii-Moriya dynamics interaction Inthissubsection,wederivethemagnon-mediatedcur- As magnons can exert a torque on magnetization even rent response to slow magnetization dynamics in a sys- in equilibrium, we begin by considering an equilibrium tem with broken inversion symmetry and spin-orbit in- stateofthesystem. SuchequilibriumDMItorquescanbe teractions. TheKubolinearresponseenergycurrentcon- captured by calculating the DMI tensor in the presence tains the ground state energy contribution related to the 3 magnon-mediatedDMIwhichhavebeencalculatedinthe expresstheresultthrougharesponsetensortK contain- aα previoussection. Thus,wewillusetheresultscalculated ing two contributions tK = tI +tII , which are given aα aα aα earlier in order to identify various transport contribu- by tions. ˆ Weareinterestedintheheat,particle,andspincurrent tI = 1 dωg(ω) d ReTr(cid:10)J GRηGA−J GRηGR(cid:11), density responses described by a tensor taα: aα (cid:126)ˆ 2π dω aα aα tII = 1 dωg(ω)ReTr(cid:10)J GRηdGR −J dGRηGR(cid:11), Jaα =−taα·∂tm, (10) aα (cid:126) 2π aα dω aα dω (15) where a is q for the heat current, p for the particle cur- where g(ω) is the Bose distribution function g(ω) = rent, and s for the spin current. Here the spin current 1/[exp((cid:126)ω/k T)−1], GR = (cid:126)((cid:126)ω −H +iΓ)−1, GA = B is related to the magnon particle current density Jp by a (cid:126)((cid:126)ω −H −iΓ)−1, η = −∂ H , J = (vH +Hv)/2, relation J =−(cid:126)J . m q s p andJ =v. Inourcalculations,weadoptaphenomeno- p In the presence of magnetization dynamics, Hamilto- logicaltreatmentandrelatethequasiparticlebroadening nian H acquires a perturbation of the form to the Gilbert damping, i.e. Γ=α(cid:126)ω. ˆ Note that the Kubo response for the energy current H(cid:48) = drΨ†(r)H(cid:48)Ψ(r), (11) density in Eq. (14) contains the bound energy current associated with DMI: where H(cid:48) = ∂ H ·δm(t) and we assume that δm(t) is m JD =Dˆ ·(m×∂ m), (16) small. We are interested in a linear response to the time q t derivative of m(t), thus we write δm(t) = (1/iω)∂ m. t Note that this calculation is similar to the calculation of where tensor Dˆ is given in Eq. (9). This current needs dc current response to electric field with the correspon- to be subtracted from the Kubo current in Eq. (14) in dence A(t)→δm(t) where the perturbation in Eq. (11) order to obtain a transport heat current: leadstoananalogofequilibriumdiamagneticcurrentcor- rection. UsingthelinearresponseKubotheoryweobtain J =JK −JD. (17) qα qα q for the heat and particle current density response: JK =(cid:10)J[0](cid:11) +(cid:10)J[1](cid:11) , (12) To express the response tensor tKaα, we use the Fourier aα aα ne aα eq transformedoperatorsandtheeigenbasisrepresentation for the velocity, (cid:126)v = ∂ E −iA E +iE A , and the or (cid:101)k k k k k k k effective field, −η = ∂ E −iA E +iE A , where (cid:101)k m k m k k m JaKα =ωli→m0(cid:110)−ΠRα(ω)/iω(cid:111)∂tm+(cid:10)Ja[1α](cid:11)eq, (13) Ak =iTk†∂kTk and Am =iTk†∂mTk. We obtain where ΠRα(ω) = Πα(ω + i0) is the retarded corre- tKqα =(cid:88)(cid:110)g(εnk)[−εnkBm(nα)(k)+A(mnα)(k)] lation function related to the follo´wing correlator in kεn g(cid:48)(ε ) (cid:111) (18) Matsubara formalism, Π (iω) = − βdτeiωτ(cid:10)T J[0]h(cid:11) − nk nk (∂ ε )(∂ ε ) , ´ α 0 τ aα 2Γ(n) m nk kα nk with h = − drΨ†(r)∂mHΨ(r) being the nonequilib- k rium field, and j[0](r) = (1/2)Ψ†(r)(vH +Hv)Ψ(r) and ´ q whichaftercombiningwithDMIenergycurrentJD leads J[0] = drj[0](r) being the heat current density and the q q q to the response tensor describing the heat current: macroscopic heat current, respectively. For the particle current we have similar expressions j[0](r)=Ψ†(r)vΨ(r) ´ p and J[p0] = drj[p0](r). Here the velocity operator is 1(cid:88)(cid:88)N 1 given by v = (1/i(cid:126))[r,H]. We also introduce a gra- teqxα =−V 2Γ(n)(∂mεnk)(∂kαεnk)εnkg(cid:48)(εnk), dient correction to the heat ´and particle currents due k n=1 k N to perturbation, i.e., J[a1] = drj[a1](r) where j[q1](r) = tin = 1(cid:88)(cid:88)c (ε )Ω(n) (k), (1/2)Ψ†(r)[(δm(t) · ∂ )(vH + Hv)]Ψ(r) and j[1](r) = qα V 1 nk mkα m s k n=1 Ψ†(r)[δm(t)·∂ v]Ψ(r). Thisanalogofdiamagneticcur- (19) m rent cancels with the term ΠR(0) resulting in the Kubo where ε = [E ] , Γ(n) = αε , g(cid:48)(ε ) = α nk k nn k nk nk contribution of the form (2k T)−1{1−cosh(ε /k T)}−1, c [ε ]=g(ε )ε − B nk B 1 nk nk nk (cid:110) (cid:111) (1/β)ln(1−e−βεnk), V is volume, and we separated the JaKα = lim [ΠRα(0)−ΠRα(ω)]/iω ∂tm. (14) total tensor tqα into the intrinsic and extrinsic contri- ω→0 butions, i.e., t = tex +tin. For the particle current qα qα qα ThecorrelationfunctioninEq.(14)iscalculatedbycon- response only tK tensor needs to be considered, thus pα sidering the simplest bubble diagram for Π and per- we obtain the following expression for the total tensor, α forming the analytic continuation, see e.g. Ref. [16]. We t = tex +tin, divided into the intrinsic and extrinsic pα pα pα 4 contributions: equality H˙(cid:48) = (i/(cid:126))[H,H(cid:48)] = J[0]∂χ and integration by q parts. Following the notations in Ref. [16], we introduce N 1(cid:88)(cid:88) 1 the linear response tensors S and M for the fields h[0] tex =− (∂ ε )(∂ ε )g(cid:48)(ε ), α α pα V 2Γ(n) m nk kα nk nk and h[1] and the total response tensor Lα = Sα +Mα k n=1 k (20) according to equation N tin = 1(cid:88)(cid:88)g(ε )Ω(n) (k). pα V nk mkα h[0]+h[1] =−Lα∂αχ, (26) k n=1 whereM isgivenbyEq.(5)asitfollowsfromEq.(23). The last tensor also describes the spin current response, α For the tensors S we obtain i.e., t =−(cid:126)t . α sα pα S =(cid:88)(cid:110)g(ε )[−ε B(n)(k)+A(n)(k)] α nk nk mβ mβ D. Thermal torques kεn g(cid:48)(ε ) (cid:111) (27) + nk nk (∂ ε )(∂ ε ) . 2Γ(n) m nk kβ nk In this subsection, we derive the magnon-mediated k magnetization torque response to a temperature gradi- We can also separate the total response tensor into the ent in a system with broken inversion symmetry and intrinsic and extrinsic contributions: spin-orbit interactions. The thermal torque is defined N according to equation 1(cid:88)(cid:88) 1 Lex = (∂ ε )(∂ ε )ε g(cid:48)(ε ), α V 2Γ(n) m nk kα nk nk nk T =−β ∂ T, (21) k n=1 k α α N Lin = 1(cid:88)(cid:88)c (ε )Ω(n) (k). whereβα isthethermaltorkancetensorandT describes α V 1 nk mkα torque acting on magnetization and leading to mod- k n=1 (28) ification of the Landau-Lifshitz-Gilbert equation, i.e., For the thermal torkance tensor, we obtain s(1 + αm×)m˙ = m × H + T where H is the ef- eff eff fective magnetic field and s is the spin density. We use β =L ×m/T. (29) α α the Luttinger linear response method [47] in which the temperature gradient is replicated by a perturbation to Hamiltonian H of the form E. Onsager reciprocity relation ˆ 1 H(cid:48) = drΨ†(r)(Hχ+χH)Ψ(r), (22) Wearenowinthepositiontocombinetheresultsfrom 2 previoussubsectionsintooneexpressionthatemphasized where we introduce the temperature gradient as ∂ χ = the Onsager reciprocity relation. In principle, the result i −∂ T/T. The torque response can be found by calculat- ofcalculationofthermaltorquesinthelastsectioncanbe i ing the effective magnon-mediated field: extractedfromtheOnsagerrelationswithoutperforming the calculation. Writing the response tensors in terms of h=h[0]+h[1] =−(cid:10)∂ H(cid:11) −(cid:10)∂ H(cid:48)(cid:11) , (23) the torkance tensors, we obtain m ne m eq   J where for the second term the averaging is done over the pα equilibrium state and for the first term over nonequilib-  Jqα = T rium state induced by the temperature gradient. The magnon-mediated torque acting on the magnetization is given by  σˆ(m) ΠˆT(−m) αα(−m)  −∂αϕ   Πˆ(m) Tλˆ(m) Tβα(−m) ∂αχ , T =m×h. (24) αα(m) Tβα(m) −Λˆ(m) m×∂tm (30) Within the linear response theory, the response h[0] to a where summation over repeated indices is implied, and temperature gradient can be calculated from expression weintroducedtheconductivitytensorσˆ(m),thethermal (cid:110) (cid:111) conductivitytensorλˆ(m),thetensorΠˆ(m)describingthe h[0] = lim [ΠR(Ω)−ΠR(0)]/iΩ ∂ χ, (25) Ω→0 α α α magnonSeebeckandPeltiereffects,andthetensorΛˆ(m) corresponding to LLG equation. The tensor α (m) was α where ΠRα(Ω) = Πα(Ω+i0) is the retarded correlation introduced by analogy with the tensor βα(m) and it is functionrelatedtothefo´llowingcorrelatorinMatsubara giveninEq.(20),i.e.,αα(−m)=tpα×m. Forcomplete- formalism, Π (iΩ) = − βdτeiΩτ(cid:10)T hJ[0](cid:11). Note that nesswealsoaddedaresponsetoananalogofelectricfield α 0 τ qα this correlator differs from the one arising in Eq. (13) in for magnons, −∂ ϕ. Equation (30) immediately follows α the order of operators. In the correlator, we reduce the from Eqs. (19), (20), and (28) given that intrinsic con- perturbation H(cid:48) to the energy current by employing the tributions are odd and extrinsic contributions are even 5 τ second-nearest neighbor DMI is in the z− direction, and d 1 3 the signs of ν are depicted in green in Fig. 2 for the ij directions shown by dashed green arrows. For analytical τ 3 d results, we assume that all DMI are small, i.e. J (cid:29)D[R] 1 andJ (cid:29)D[z]. Inourmodel,weassumethattheorderis τ d 2 in general (mx,my,mz) direction, which can be realized 2 + by application of the magnetic field. Our strategy would a y 1 betofirstunderstandtheroleoftheDMIinthebehavior of magnons for a general direction of the ferromagnetic a - 2 order. After that we will assume that the main order is inthe z− direction, whiletheperturbationsthatdeviate x the order are in the x−y plane. To study the magnons, we perform the Holstein-Primakoff transformation. The Figure2. Schematicsofthegraphenelayerparametersforthe unitcellofthehoneycombferromagnethastwospinsS tight-binding model. Vectors connecting nearest neighbors A aarree uτs1ed=in12(d√e1r3iv,1in),gτth2e=Ha12m(i√l1t3o,n−ia1n),foarndmaτg3no=ns.√13V(−ec1t,o0rs) aannddSb†B(,rh),enbc(re)thcoertrwesopsoentsdionfgbtoosotnheopAeraantdorBs,asu†(brl)a,ttaic(res) √ √ ase1co=nd21-(nea3r,e1s)t,naenidghab2or=D12M(I.3,−1) are used in deriving the areraedisnatrsoudsuucaeld,.SzT=heSH−oals†tae,inan-PdrSim+a=ko√ff2tSra(nSsxfo+rmiSayt)io=n √ A A A A 2S−a†aa (S is the total spin), and the same for B spins. The Fourier image of the Hamiltonian describing under magnetization reversal. The Onsager reciprocity non-interacting magnons written in terms of the Ψ = relation in Eq. (30) is similar to expressions obtained for (a , b )T spinor is k k similar electron-mediated effects in Ref. [25]. Equation (30) can be modified to account for the possibility of (cid:20)3+∆ −γ˜ (cid:21) H =JS k k , (33) magnonaccumulationresultingfromthemagnonmotive −γ˜∗ 3−∆ k k force [28] or temperature gradient [48]. (cid:104) (cid:16) (cid:17) (cid:16)√ (cid:17)(cid:105) where ∆k = 2∆ sin(ky)−2sin k2y cos 32kx , with III. RESULTS FOR HONEYCOMB AND ∆ = m D[z]/J. This type of DMI is a k− dependent z KAGOME FERROMAGNETS mass of magnons. Deriving γ˜ we considered Rashba k DMIinthelowestorderinD[R]/J (cid:28)1parameter. With Inthissection,weapplyourtheorytosinglelayerhon- this assumption eycomb and kagome ferromagnets with DMI. To demon- (cid:32) (cid:33) smtraagtneonexcpulrirceitnltys,hwoewdefiscctriitbioeutsheelheocntreiyccofimelbdssyrsetseumltanin- γ˜k =2ei2k˜√x3 cos k˜2y +e−i√k˜x3, (34) alytically. In principle, our results could also be relevant to layered structures with weakly coupled layers. √ √ where k˜ = k − 3D[R]m , and k˜ = k + 3D[R]m . x x J y y y J x We observe that Rashba DMI plays an effective role of A. Application to honeycomb ferromagnet magnon charge, while order direction (mx,my,0) is an effective vector potential felt by magnons. The eigenvalues of the Hamiltonian are calculated to In this subsection, we study a model of an insulat- be, ing ferromagnet on a honeycomb lattice. This model contains physics discussed above in a transparent and (cid:18) (cid:113) (cid:19) analytical way. We assume a Heisenberg exchange of (cid:15) =JS 3± ∆2 +|γ˜ |2 , (35) k,± k k ferromagnetic sign, in-plane DMI of Rashba type, and second-nearest neighbor DMI. The Hamiltonian is with corresponding eigenfunctions (cid:88) (cid:88) H =−J SiSj + D[R][Si×Sj] (31) vk,+ =[cos(ξ˜k/2)eiχ˜k, −sin(ξ˜k/2)]T, (36) <ij> <ij> (cid:88) +D[z] ν [S ×S ] . (32) and ij i j z <<ij>> vk,− =[sin(ξ˜k/2), cos(ξ˜k/2)e−iχ˜k]T, (37) ThevectorsoftheRashbatypeDMIareshowninFig.2, √ √ where d1 = 21( 3,−1), d2 = 12(− 3,−1), and d3 = where sin(ξ˜k) = |γ˜k|/(cid:112)∆2k+|γ˜k|2, and γ˜k = |γ˜k|eiχ˜k, (0,1), such as D[R] =D[R]d. Note that all vectors, such and the tilde symbol here means that corresponding k as τ and a , are measured in units of lattice spacing a momentaareshiftedbytheRashbaDMI.Unitarymatrix i i 0 which is recovered in the final result. The vector of the that diagonalizes the Hamiltonian is readily constructed 6 and it is given by a0αxeyven/kB a0βxeyven/kB  (cid:16) (cid:17) (cid:16) (cid:17)  0.004 cos ξ˜2k eiχ˜k sin ξ˜2k 0.0006 num naunmal Tk = −sin(cid:16)ξ˜k(cid:17) cos(cid:16)ξ˜k(cid:17)e−iχ˜k . (38) 0.0003 anal 0.002 2 2 Wearenowreadytoderivespinandheatcurrentswhich 0.25 0.5kBT/SJ 0.25 0.5kBT/SJ are driven by magnetization dynamics. We set the dom- inant component of the ferromagnetic order in the z− Figure 3. (Color online) Left: The even component under direction and assume that the magnetization dynamics magnetizationreversalofthetensorα asafunctionoftem- is in the x−y plane. We only focus on the intrinsic con- ij perature. Right: The even component under magnetization tribution to the currents, i.e., due to non-trivial Berry reversal of the torkance tensor β as a function of tempera- ij curvaturesofthemagnonbandstructure. Anexpression ture. In both cases the magnetization is along the z− axis. defining the Berry curvature is For the strength of DMI we use D[z] = D[R] = J/6. Red curvescorrespondtonumericalresultsandbluecurvescorre- Ω =2Im(cid:104)(cid:16)∂ T†(cid:17)(cid:0)∂ T (cid:1)(cid:105)= 1sin(cid:16)ξ˜ (cid:17) (39) spond to analytical results in Eqs. (45) and (46). α,mβ α k mβ k 2 k ×(cid:104)(∂ χ˜ )(cid:16)∂ ξ˜ (cid:17)−(cid:0)∂ χ˜ (cid:1)(cid:16)∂ ξ˜ (cid:17)(cid:105)(cid:20)1 0 (cid:21). α k mβ k mβ k α k 0 −1 SJ (cid:29)T, reads (40) D[R] √3(cid:20) (cid:34)13√3D[z](cid:35) In the following, we focus on the α = x and β = x case, Jpx = J a π sinh z J e−z3 (45) 0 and mention β = y case at the end. Recall that ∆ √ doesnotdependonmβ forβ =(x,y)components,hencke + D[z] 3ζ(3)z3(cid:21)(∂ m) , ∂ ∆ =0. The derivativewithrespect to the direction J 36 t x mβ k of the order m of the remaining functions that depend β whereweintroducedz =T/SJ forbrevity,andsetm = on k˜ is z 1. Similarly, from Eq. (18), the heat current due to the ∂ √ D[R] ∂ √ D[R] Berry curvature at small temperatures, SJ (cid:29)T, reads = 3 ≡ 3 ∂ , (41) ∂mx J ∂k˜y J y D[R]3√3(cid:20) (cid:32)13√3D[z](cid:33) ∂m∂y =−√3DJ[R]∂∂k˜x ≡−√3DJ[R]∂x. (42) Jqx =JS J a0π siDnh[z]√z3I J(cid:21) e−z3 (46) + z4 (∂ m) . This straightforward transformation makes the mixed J 216 t x Berry curvature a regular k− space one, except for the a∂tmβth∆ekK=(cid:48)0=co(cid:0)n0,d4itπio(cid:1)na.nTdhKe B=er(cid:0)ry0,c−ur4vπa(cid:1)tuproeinhtass, eaxntdrecmana IKn bpootihntsc,asweshiale ttehrme r∝emaei−ni3nSTJg iosneduies tdoueK(cid:48)toanΓd 3 3 be approximated as point. We introduced a numerical constant I = ´ (cid:104) (cid:105) 27D[R] ∆ (cid:20)1 0 (cid:21) 0∞dxx2 xexe−x1 −ln(ex−1) = 4π4/45 ≈ 8.65, and Ωx,mx|K(K(cid:48)) ≈− 8 J (cid:0)27∆2+ 3k2(cid:1)3/2 0 −1 . Riemann zeta-function ζ(3)≈1.2. 4 It is straightforward to show that Berry curvature (43) parts of the J and J currents driven by (∂ m) mag- px qx t y netization dynamics vanish. The J and J currents ThecurvatureisthesameforbothK(cid:48)andKpoints. The sy qy √ driven by (∂ m) magnetization dynamics will have the spectrumatthesepointsisfinite,(cid:15) ≈JS(3±3 3|∆|), t y k,± same expressions as in Eqs. (45) and (46). Thus, we cal- but the Berry curvature is of the monopole type. Hence culated even under magnetization reversal components at small temperatures, despite the exponential suppres- αeven = −αeven and βeven = −βeven as it follows from sionofthemagnonnumberattheK(cid:48) andKpoints,there xy yx xy yx Eq. (30). As can be seen from Fig. 3, Eqs. (45) and might be a contribution to the magnon currents due to (46) only qualitatively agree with the numerical results this Berry curvature. At the Γ = (0,0) the spectrum of athighertemperaturesastheBerrycurvaturefromother the lowest band is (cid:15) ≈ 1SJk2, and it will be popu- k,− 4 parts of the Brillouin zone starts to contribute to the re- latedbythemagnonsthemostatlowtemperatures. The sult. Berry curvature is approximated close to this point as D[R] ∆ (cid:20)1 0 (cid:21) Ωx,mx|Γ ≈− J 48ky2kx2 0 −1 . (44) B. Application to kagome ferromagnet AccordingtoEqs.(10)and(20),theparticlecurrentden- Here we apply our theory to the kagome lattice fer- sity due to the Berry curvature at small temperatures, romagnet with the nearest neighbor DMI. The lattice of 7 ν=1 DR 00..2550 α0.a60αxoxdd D=0.3J 0.a002αyexven D=0.3J a2 0 0.3 D=0.2J D=0.2J y x -0.25 D=0.1J 0.01 D=0.1J a1 -0.50 kBT/SJ kBT/SJ 0.5 1. 1.5 0.5 1. 1.5 Figure 5. (Color online) Left: The odd component of the Figure4. (Coloronline)Left: Atwo-dimensionalkagomelat- √ √ tensor α as a function of temperature. The plot is rescaled tice with lattice vectors a = 1( 3,−1) and a = 1( 3,1) ij 1 2 2 2 bymultiplyingitwiththeGilberdampingα. Right: Theeven where atoms are placed in the corners of triangles. Rashba- componentofthetensorα asafunctionoftemperature. In likeDMIvectorsD[R]areshownbybluevectorsperpendicular ij ij both cases the magnetization is along the z− axis. For the tothebonds. Theclockwiseorderingofbondscorresponding strength of the Rashba DMI we use D[R] =D[z] =D. to ν =1 is shown by black arrows. Right: Magnon spectrum of a kagome ferromagnet with DMI D[z] = 0.3J and magne- tizationpointinginthez−direction. Thedistributionofthe αa0βxoxdd/kB a0βyexven/kB BerrycurvatureovertheBrillouinzoneisplottedbythecolor 0.2 D=0.3J D=0.3J coding on top of the spectrum for each subband. 0.04 D=0.2J 0.1 D=0.2J D=0.1J 0.02 D=0.1J thesystemanditsmagnonspectrumareshowninFig.4. kBT/SJ kBT/SJ 0.5 1. 1.5 0.5 1. 1.5 Note that all vectors, such a and a , are measured in 1 2 units of lattice spacing a which is recovered in the final 0 result. We consider a model considered in Ref. [16] with Figure 6. (Color online) Left: The odd component of the a Hamiltonian given by torkance tensor β as a function of temperature. The plot is ij rescaledbymultiplyingitwiththeGilberdampingα. Right: (cid:88) (cid:88) (cid:88) H =−J SiSj −B Siz+ νijDij[Si×Sj], The even component of the torkance tensor βij as a function of temperature. In both cases the magnetization is along the <ij> i <ij> z− axis. For the strength of the Rashba DMI we use D[R] = (47) D[z] =D. whereJ >0correspondstoferromagneticnearestneigh- bor exchange, B is the external magnetic field, and ν ij logicaltreatmentbyrelatingthequasiparticlebroadening describesasignconventionforthenearestneighborDMI, to the Gilbert damping as Γ=α(cid:126)ω. Under a simple cir- i.e., ν = 1 for the clockwise sense of direction and ij cular precession of the magnetization described by angle νij = −1 otherwise (see Fig. 4). Note that vectors θ we have ∂ m=θω[−sin(ωt),cos(ωt),0]T and t D =D[z]zˆ+D[R] have an in-plane Rashba-like compo- ij ij J =θω[αoddcos(ωt)−αevensin(ωt)], nent D[R] directed orthogonally to bonds and outwards px xx yx (48) ij J =θω[αoddsin(ωt)+αevencos(ωt)]. with respect to bond triangles (see Fig. 4). The Rashba- py xx yx like DMI could result from mirror asymmetry with re- Wecannowestimatetheamplitudeofacspincurrentas spect to the kagome planes. At sufficiently low tempera- (cid:113) θ(cid:126)ω (αodd)2+(αeven)2. For a three-dimensional sys- turestheHamiltonianinEq.(47)canbeanalyzedbyap- xx yx plyingtheHolstein-Primakofftransformation. Thecorre- temcontainingweaklyinteractingkagomelayers, wecan sponding magnon spectrum is shown in Fig. 4 where the write αi3jD = αisj/c where c ∝ a0 is the interlayer dis- lower, middle, andupperbandshavetheChernnumbers tance which is comparable to the lattice constant a0. −1, 0, and 1, respectively. For parameters D[z] = 0.1J, D[R] = 0.1J, θ = 0.1◦, We begin by analyzing an effect of magnon pumping ω = 2π×10GHz, kBT = 0.5SJ, and the Gilbert damp- by magnetization dynamics. This effect is characterized ing α = 0.1, we obtain the spin current of amplitude by tensor αα or equivalently by Eq. (10). It is also clear Js ≈10−8J/m2. We suggest to detect such spin currents fromEq.(30)thatthesametensoralsodescribesamag- by the ac inverse spin Hall effect [49]. netizationtorqueinducedbyananalogofelectricfieldfor We also consider an effect of heat pumping by magne- magnons. Weassumeasmall-angleprecessionofmagne- tization dynamics. This effect is characterized by ten- tization around the z− axis. By symmetry considera- sor βα. Here we again assume a small-angle preces- tion, it is sufficient to consider only αeven = −αeven and sion of magnetization around the z− axis. Similar sym- yx xy metry considerations result in relations βeven = −βeven αodd =αodd componentsofthetensorwhereweseparate yx xy texnxsor ααyyinto the parts that are odd and even under and βxoxdd = βyoydd between non-zero components of ten- magnetization reversal, i.e., α = αodd + αeven. The sor β = βodd + βeven separated into the odd and α α α α α α results of our calculations for the two components are even under magnetization reversal parts. The results of shown in Fig. 5. Note that we use a simple phenomeno- our calculations for the two components are shown in 8 Fig. 6. The amplitude of ac heat current is given by ing on the ac inverse spin Hall effect [49]. Additionally, (cid:113) θTω (βodd)2+(βeven)2 whichfortheaboveparameters we obtain an analog of the Hall-like response in systems xx yx with non-trivial Berry curvature of magnon bands. This and T = 50K results in the heat current of amplitude leads to even under magnetization reversal contributions J ≈50kW/m2. q to the response tensors. By the Onsager reciprocity re- After invoking the Onsager relation (30) one can con- lation, this Hall-like response can be related to the anti- firm that estimates obtained in this subsection are com- damping thermal torque [16]. Finally, we identify the parable to estimates for thermal torques obtained in groundstateenergycurrentassociatedwiththemagnon- Ref. [16]. Note also that the phenomenology discussed mediated equilibrium contribution to DMI. This contri- in this paper is similar to Ref. [26], however, the heat bution needs to be subtracted from the Kubo linear re- current is carried by magnons in contrast to electronic sponse result according to our analysis. mechanisms considered before. IV. CONCLUSIONS ACKNOWLEDGMENTS In this work, we explored fictitious electric fields act- ing on magnons in response to time-dependent magne- We gratefully acknowledge useful discussions with tization dynamics in the presence of DMI. We find that K. Belashchenko. This work was supported by the U.S. such fictitious electric fields can drive sizable spin and Department of Energy, Office of Science, Basic Energy energy currents. We suggest a detection scheme rely- Sciences, under Award No. de-sc0014189. Appendix A: Heat current as a response to magnetization dynamics Measurable heat current consists of three parts. Free energy contribution, and non-equilibrium heat current and orbital magnetization heat current carried by magnons. 1. Free energy heat current Magnon mediated Dzyaloshinskii-Moriya interaction contribution to the free energy of the system is (cid:20) (cid:21) ∂m(r) FDMI =D m(r)× (A1) ∂r where D is the Dzyaloshinskii-Moriya tensor we will calculate below. For instance, functionality on x might be DMI due to the boundary or it might be due to spatially dependent magnetization profile. Assuming a time dependence of the magnetization, via a r→r+ωt shift, one can derive the current due to time dependence of DMI part of free energy using continuity equation ∂FDMI +JDMI =0, where ∂t 1 JDMI =− D (∂ m) , (A2) α V αβ t β where V is the volume of the system. The Dzyaloshinskii-Moriya interaction constant is ˆ (cid:28) (cid:29) 1 D = drΨ†(r)(r T +T r )Ψ(r) , (A3) αβ 2 α β β α eq where T =(∂ H ×m) is the torque operator. To calculate the DMI, we introduce β m β (cid:20) dG+ dG−(cid:21) A (η)=iTr v v¯ δ(η−H )−v δ(η−H )v¯ , (A4) αβ αk dη βk k αk k βk dη B (η)=iTr(cid:2)v G+v¯ δ(η−H )−v δ(η−H )v¯ G−(cid:3), (A5) αβ αk βk k αk k βk 9 where v¯ = ∂ H ≡ i[H ,r ] equivalent to the velocity operator definition, with r ≡ i∂ equivalent to the βk mβ k k mβ mβ mβ position operator. It was shown that A − 1dBαβ = 1 Tr(cid:2)r (GA−GR)r −r r (GA−GR)(cid:3)−(α↔β) (A6) αβ 2 dη 4π α mβ α mβ (cid:20) (cid:21) 1 d + Tr (r v¯ −v r ) δ(η−H ) . (A7) 2 α βk αk mβ dη k Also, we derive the Berry curvature parts of A and B . αβ αβ (cid:88)(cid:16) (cid:17) A (η)=−i ∂ T†∂ T δ[η−((cid:15) ) ]−(α↔β), (A8) αβ α k mβ k nn k nn n and a Berry curvature part of the B as αβ (cid:88)(cid:104) (cid:105) B (η)=i ∂ T†(η−H )∂ T δ[η−((cid:15) ) ]−(α↔β). (A9) αβ α k k mβ k nn k nn n Therefore, ˆ ˆ D =(cid:88) ∞ dη˜(cid:20)A (η˜)− 1dBαβ(η˜)(cid:21) η˜dηg(η), (A10) αβ αβ 2 dη˜ k −∞ 0 and it can be shown that ˆ ˆ D =(cid:88) +∞dη˜(cid:20)A (η˜)− 1dBαβ(η˜)(cid:21) η˜dηg(η) (A11) αβ αβ 2 dη˜ n −∞ 0 ˆ (cid:40) ˆ (cid:41) (cid:88) +∞ (cid:16) (cid:17) η˜ = dη˜ −i ∂ T†∂ T δ[η˜−((cid:15) ) ] dηg(η) (A12) n −ˆ∞ α k mβ k nn k nn 0 i (cid:88) +∞ (cid:110)(cid:104) (cid:105) (cid:111) + dη˜ ∂ T†(η˜−H )∂ T g(η˜)δ[η˜−((cid:15) ) ] −(α↔β). (A13) 2 n −∞ α k k mβ k nn k nn 2. Heat current due to magnons We assume that the magnetizaion is varying in time. Next, we assume that due to that there is a time-dependent term in the Hamiltonian. For example, since the DMI depends on the direction of the order, this DMI will be time dependent. The Hamiltonian of the spin waves is then ˆ 1 (cid:104) (cid:105) H = drΨ†(r) Hˆ +Hˆ(cid:48)(t) Ψ(r). (A14) T 2 WedefineHˆ =Hˆ+Hˆ(cid:48)(t). Microscopicexpressionfortheheatcurrentcurrentisderivedviacommutationrelationship T 1 (cid:16) (cid:17) j (r)= Ψ†(r) Hˆ V+VHˆ Ψ(r), (A15) Q 2 T T here V = i[Hˆ ,r] is the full velocity. Velocity has two parts, V = v +v(cid:48), where v = i[Hˆ,r] and v(cid:48) = i[Hˆ(cid:48),r]. T (cid:16) (cid:17) Assumingthatthemagneticorderism(t)=m+δm(t),wewritetheperturbationasHˆ(cid:48)(t)= ∂ Hˆ δm(t). Wewill m useanalogybetweenmagnetizationdynamicsandtheelectromagneticwaves. Thedirectionofthelocalmagnetization canbeseenasavectorpotentialforeffectiveelectromagneticfieldelectricandmagneticfields. Then, ∂m isanalogous ∂t to the electric field, while ∇×m is analogous to the magnetic field. We will then write δm(t)= 1 ∂m(t) ≡ 1∂ m (in ω ∂t ω t Matsubara frequency). The heat current is separated in to two parts 1 (cid:16) (cid:17) j[0](r)= Ψ†(r) Hˆv+vHˆ Ψ(r) (A16) Q 2 10 1 (cid:16) (cid:17) 1 (cid:16) (cid:17) j[1](r)= Ψ†(r) Hˆ(cid:48)v+vHˆ(cid:48) Ψ(r)+ Ψ†(r) Hˆv(cid:48)+v(cid:48)Hˆ Ψ(r) (A17) Q 2 2 1 (cid:104) (cid:16) (cid:17)(cid:105) = Ψ†(r) (δm(t)·∂ ) Hˆv+vHˆ Ψ(r) (A18) 2 m ´ The H(cid:48)(t) will be treated as a perturbation. We will be working with global currents J ≡ 1 drj (r). The heat Q V Q current is conveniently written as (cid:68) (cid:69) (cid:68) (cid:69) (cid:104)J (cid:105)= J[0] + J[1] , (A19) Q Q Q ne eq Where the former one is estimated over non-equilibrium states and is given by Kubo formula, while the later one is due to orbital magnetization of the magnons and is estimated over equilibrium states. a. Non-equilibrium heat current, Kubo formula Kubo formula for an arbitrary operator A(ω), where ω is Matsubara frequency, is ˆ β (cid:104)A(ω)(cid:105) = dτeiωτ(cid:104)T A(0)H(cid:48)(−τ)(cid:105) , (A20) ne τ eq 0 ´ where H(cid:48)(τ)= drΨ†(τ,r)Hˆ(cid:48)Ψ(τ,r) is the perturbing Hamiltonian. ˆ (cid:68) (cid:69) 1/T (cid:68) (cid:69) 1 1 J[0] = dτeiωτ T J[0](0)H(cid:48)(−τ) ≡ S (ω) (∂ m) (A21) Qα ne 0 τ Qα eq V αβ ω t β After all of the transforms, we get S = 1(cid:88)((cid:15) v˜ +v˜ (cid:15) ) (cid:104)T†(∂ H )T (cid:105) g[((cid:15)k)nn]−g[((cid:15)k)mm] , (A22) αβ 2 k αk αk k nm k β k k mn iω+((cid:15)k)nn−(σ3(cid:15)k)mm k where v˜ =T†v T =∂ (cid:15) +A (cid:15) −(cid:15) A , (A23) αk k αk k α k αk k k αk v˜¯ =T†v¯ T =∂ (cid:15) +A¯ (cid:15) −(cid:15) A¯ , (A24) βk k βk k mβ k βk k k βk where A = T†∂ T , and A¯ = T†∂ T ≡ T†∂ T , and where a bar over A¯ symbolizes information that the αk k α k βk k mβ k k β k βk derivative is over β component of the magnetization direction, m . After the transformations we get, β S (ω)= 1(cid:88)g[((cid:15)k)nn]−g[((cid:15)k)mm](v˜ ) [((cid:15) ) +((cid:15) ) ](cid:0)∂ (cid:15) +A¯ (cid:15) +(cid:15) A¯ (cid:1) (A25) αβ 2 iω+((cid:15) ) −((cid:15) ) αk nm k nn k mm mβ k βk k k βk mn k nn k mm kn Expand S (ω) in ω, and get αβ ∂ S (ω)=S[1](0)+S[2](0)+ S[2](ω)| ω, (A26) αβ αβ αβ ∂ω αβ ω=0 where n=m parts of S are αβ S[1](0)=(cid:88)∂g((cid:15))| ((cid:15) ) (∂ (cid:15) ) (∂ (cid:15) ) =−1(cid:88)g[((cid:15) ) ](∂ ∂ (cid:15)2) , (A27) αβ ∂(cid:15) (cid:15)=((cid:15)k)nn k nn α k nn mβ k nn 2 k nn α mβ k nn kn kn where we integrated by parts over k. Term with n(cid:54)=m elements reads S[2](ω)=− 1(cid:88)g[((cid:15)k)nn]−g[((cid:15)k)mm][((cid:15) ) +((cid:15) ) ][((cid:15) ) −((cid:15) ) ]2(A ) (cid:0)A¯ (cid:1) , (A28) αβ 2 iω+((cid:15) ) −((cid:15) ) k nn k mm k nn k mm αk nm βk mn k nn k mm kn