Origin of longitudinal spin excitations in iron-pnictide parent compounds

We investigate longitudinal spin excitations (LSEs) as a probe of microscopic origin of magnetic ordering in parent pnictides BaFe2As2 and NaFeAs. Currently adopted interpretation of LSEs as bottom of particle-hole continuum points unambiguously toward itinerant-electron magnetism, but is difficult to reconcile with available optical measurements. We study the possibility that the LSEs originate from multi-magnon processes which are not energetically constrained by optical spectroscopy and do not sharply distinguish between local-moment and itinerant scenarios. Two mechanisms, capable of enhancing multi-magnon continuum to the level indicated by neutron scattering experiments, are proposed. The first emphasizes itinerant electrons and is based on electronic transitions between magnetically split bands, while the other relies on purely spin fluctuations close to a magnetic quantum phase transition. Electronic excitations enhance multi-magnon contribution to LSEs for small Fermi surface taking part in the SDW instability, but are insufficient to account for measured intensities. The correct order of LSEs, on the other hand, can be reproduced by the spin fluctuation mechanism for a reasonable set of parameters.


Introduction
Generic proximity of iron-based superconductors (FeSCs) to antiferromagnetism [1] makes identification of the microscopic origin of magnetic ordering an important step toward understanding global phase diagram of these materials. The multi-orbital character of FeSCs allows for three scenarios differing by the role of electronic and spin degrees of freedom. The first, itinerant-electron picture, relies on nesting of the hole-and electron-like Fermi surfaces (FSs) [2,3] which leads to spin density wave (SDW) instability, similarly as in chromium alloys [4]. The second, local-moment scenario, is based on variations of frustrated Heisenberg model [5][6][7][8][9] and emphasizes purely spin fluctuations. The third possibility is that FeSCs should be viewed as a mixture of itinerant electrons and preformed local magnetic moments [10][11][12].
The origin of magnetic ordering in FeSCs is closely related to mechanism of longitudinal spin excitations (LSEs), polarized along the direction of the ordered moments, which have been recently observed in iron-pnictide parent compounds NaFeAs [13] and BaFe 2 As 2 [14] by means of inelastic neutron scattering (INS) experiments. Since, classically, the length of the order parameter does not fluctuate in a local-moment system, the LSEs have been interpreted [13,14] as bottom of the particle-hole (P-H) continuum of an itinerant antiferromagnet (see, e.g., [15]). The latter identification favors a major contria e-mail: maciej.fidrysiak@pwr.edu.pl bution of itinerant electrons to magnetism and calls into question local-moment approaches to these materials.
The P-H excitations are, however, constrained by the details of the electronic structure. Several studies of magnetic dynamics in excitonic SDW systems, based on threeband model with elliptical electron pockets [3], two-band models with perfectly or imperfectly nested FS [16,17], and one-band model [18] agree that the P-H continuum extends above the threshold of twice the SDW gap 2Δ at the ordering wavevector (though it is typically shifted to energies lower than 2Δ at other points in the Brillouin zone). Similarly, more realistic five-orbital model calculations yield substantial suppression of the P-H contribution to the LSEs near the magnetic zone center below the partially opened gap inferred from the density of states [19]. For BaFe 2 As 2 , optical spectroscopy resolves two gap scales 2Δ ≈ 45 meV and 2Δ ≈ 110 meV below the SDW ordering temperature [20,21]. This multi-gap structure is qualitatively consistent with theoretical calculations of optical conductivity, which, however, suggest larger energies ∼100 meV and ∼150 meV [22]. An independent indication of possible opening of the full SDW gaps of 2Δ ∼ 60 meV and 2Δ ∼ 100 meV comes from angle-resolved photoemission spectroscopy (ARPES) [23].
The LSEs in BaFe 2 As 2 near the magnetic zone center appear sharply above the threshold of ≈18 meV, i.e. at energies substantially smaller than the lowest present estimate of 2Δ ∼ 45 meV. While the P-H origin of the LSEs cannot be definitely excluded due to complex electronic structure of iron-pnictides, this observation motivates a search for mechanisms capable of giving rise to in-gap longitudinal magnetic response. Since the experimental threshold for LSEs coincides with the minimal energy required to excite a pair of transverse magnons, one of the candidates is two-magnon (2M) continuum. The latter, contrary to the P-H excitations, is not energetically constrained by the value of the SDW gap, but is sensitive to the structure of the low-energy collective modes instead (e.g., to anisotropy-induced magnon gaps). More importantly, 2M excitations arise both in local-moment and itinerant-electron systems. Sole observation of lowenergy LSEs is hence not sufficient to prove contribution of itinerant electrons to magnetism and further quantitative analysis is required to relate the LSEs to the microscopic mechanism of magnetic ordering.
In this paper we estimate, within both local-moment and band models, fraction of the total low-energy spectral weight consumed by the LSEs originating from multimagnon processes in BaFe 2 As 2 . First, we demonstrate that the 2M continuum, obtained within the anisotropic, local-moment Heisenberg model (extensively used to study spin-waves in iron pnictides), qualitatively reproduces measured longitudinal magnetic response, but calculated INS intensities are too low to account for reported values. Second, we propose two mechanisms capable of enhancing the multi-magnon contribution to LSEs relative to the Heisenberg model prediction and discuss their relevance to iron-pnictides. Within a minimal two-band model, consistent with strictly itinerant-electron magnetism (favored by interpretation of the LSEs as P-H excitations), we find that the 2M contribution to the total low-energy spectral weight becomes substantial for small FS taking part in the SDW instability and scales as ∼k −1 F with Fermi momentum k F . For the typical range of parameters, however, we obtain only twofold enhancement relative to the Heisenberg model result, not enough to match the reported range of intensities. The other proposed mechanism is based on purely spin fluctuations and does not directly involve electronic transitions. The physical motivation to evoke spin fluctuations here is substantial difference between measured ordered and local moments in these materials. To explain the latter, proximity to a frustration-driven magnetic quantum phase transition (QPT) was postulated in early works [8]. More recent approaches usually involve both local moments and itinerant carriers in the spirit of orbital-selective Mottness [24]. Specifically, it has been demonstrated that coupling to itinerant electrons is capable of enhancing quantum fluctuations of the local moments [25]. Here, by employing a phenomenological non-linear σ-model approach with empirical Landau damping term and relating small measured ordered moments to spin fluctuations, we show that the correct order of LSEs can be reproduced by multi-magnon processes for a reasonable set of parameters.

The two-magnon process in iron-pnictides
The 2M contribution to LSEs can be discussed in the most straightforward manner within the Heisenberg model given by the Hamiltonian previously used to study spin waves (SWs) in BaFe 2 As 2 [14,26,27]. The summations are performed over Fe sites and the in-plane exchange constants read SJ 1a = 59.2 ± 2.0 meV, SJ 1b = −9.2 ± 1.2 meV, SJ 2 = 13.6 ± 1.0 meV [26]. In reference [14] spin anisotropies SD x = 0.196 meV, SD y = −0.311 meV, and the out-of-plane exchange SJ c = 0.333 meV have been chosen so that the calculated SW energies correspond to 90% of the maximum magnetic response. Since the SWs acquire significant width in BaFe 2 As 2 and the threshold for the 2M process is determined by the minimal energy required to excite a pair of transverse magnons, we identify SW masses with the energy for which the transverse response becomes non-zero at the magnetic Brillouin zone center. This procedure gives somewhat lower estimates SD x = 0.117 meV, SD y = −0.166 meV, and SJ c = 0.231 meV which are, however, in good agreement with SJ c = 0.22 meV reported in reference [27].
To the leading order of the spin-wave theory (SWT) [28], imaginary parts of the transverse (χ yy 0 , χ zz 0 ) and longitudinal (χ xx 0 ) spin susceptibilities (for ω > 0) are given by: and Here q SDW = (π, 0, π) is the ordering wavevector corresponding to the columnar antiferromagnetic phase of BaFe 2 As 2 , and N counts the lattice sites. In Figure  1 the longitudinal component Imχ xx 0 (ω, q SDW ) is plotted as a function of energy. Two main qualitative features of the INS data [14] are reproduced, i.e. threshold of ≈18 meV for LSEs and suppression of the scattering intensity around ≈28 meV (seen as a cusp in Fig. 1) which appears at slightly higher Energy ω (meV)  (1) at the magnetic Brillouin zone center. The observed [14] threshold of ∼18 meV for the LSEs is reproduced. energy than experimental ≈25 meV though. Quantitative analysis is not straightforward due to lack of absolute intensity measurements. The LSEs can be, nonetheless, estimated to account for ≈ 5.7% of the total low-energy response integrated up to ≈25 meV at the magnetic zone center. We define I αα ≡ 25 meV 0 dωImχ αα 0 (ω, q SDW ) and obtain I xx /(I xx + I yy + I zz ) ≈ 0.05-0.09%, depending on the value of effective spin S = 0.28-0.54 measured [29,30] for BaFe 2 As 2 . To make comparison with INS, one should take into account broadening of the SWs [26], not reproduced by the SWT susceptibilities (2) and (3). This would give larger ratio I xx /(I xx +I yy +I zz ), but is unlikely to explain ≈5.7% reported for BaFe 2 As 2 .

Two-band model
Since the Heisenberg model (1) significantly underestimates the LSEs, we consider the possibility that electronic transitions between magnetically split bands soften the magnetic excitations and enhance the contribution of the 2M continuum to the total low-energy spectral weight. We employ a minimal two-dimensional (d = 2) model of multiband magnetism involving two species of fermions with perfectly nested electron-and hole-like FSs. While, given the complex electronic structure of iron-pnictides, this approach is certainly an oversimplification, it shares several relevant features with more realistic models and hence can provide insight into the qualitative trends. Specifically, including ellipticity of electron pockets and increasing the number of bands does not change the threshold behavior of the P-H contribution to LSE continuum at the magetic zone center, even if the FS is not fully gapped upon magnetic transition [3,17]. This threshold is expected to be important as it sets the scale for virtual electronic transitions which generate interactions between low-energy collective modes. Moreover, the structure of the low-energy spin-waves in the SDW phase (whose multiple excitations give rise to multi-magnon continuum) is robust to variations of the microscopic details as a consequence of the SDW coherence [18].
The model is given by the Hamiltonian Here "ex" stands for "excitonic", species of fermions are indexed by (l) with l = 1, 2, (1,2) k = 2t(cos k x + cos k y ) ± e 0 are singleparticle dispersions controlled by the hopping integral t and the energy shift e 0 , and H int is the interaction term. At half-filling, hole-and electron-like FSs are perfectly nested with the wavevector q 0 = (π, π), which allows to define a single function e k ≡ 2t(cos k x + cos k y ) − e 0 so that In Appendix A we show that the relative one-and two-magnon intensities are not sensitive to the precise form of the single-particle dispersions Out of five independent local coupling channels contributing to H int , we retain a combination controlled by two parameters U 1 and U 2 with U 2 U 1 . This choice ensures that SDW is favored over different possible ground states [31,32] and leads to technical simplifications.
Discussion of the 2M contribution to LSEs requires going beyond random phase approximation (RPA) which underestimates fluctuations and does not capture the 2M continuum. Here we employ a bosonization scheme based on decoupling of H int by two sets of Hubbard-Stratonovich fields φ (1) and φ (2) (one for each coupling constant U 1 and U 2 ), followed by integrating out fermions. Technical details are shifted to Appendix B (see also Refs. [33,34]). This approach reproduces RPA spin susceptibilities in the saddle-point approximation and allows to systematically include corrections responsible for the 2M processes by means of loop expansion. In particular, the zero tempera- P,α sector, which does not take part in the SDW instability.
The low-energy 2M contribution to longitudinal dynamical spin susceptibility is given by the one-loop diagram shown in Figure 2 which, physically, describes decays of the longitudinal excitation mode into a pair of transverse magnons. For 0 < ω 2Δ and T → 0, we get where R ex 2M = (2dxΔ · U 1 · J ex eff ) −2 is 2M intensity renormalization factor.
The dependence of R ex 1M , R ex 2M , and J ex eff on t/U 1 for fixed e 0 /t = 3.0 is summarized in Figure 3a. Both renormalization factors undergo reduction for increasing t/U 1 , but R ex 2M /R ex 1M rapidly saturates and has a finite limit as t/U 1 → ∞. For t/U 1 ≈ 0.3, where the asymptotic behavior sets in, there is a small hump in J ex eff (thin dashed line) which signals crossover from strong-to weak-coupling regime. Figure 3b shows R ex 2M /R ex 1M in the weak coupling limit (t/U 1 → ∞) as a function of 4 − e 0 /t (the expression 4 − e 0 /t is a monotonic function of the FS size and, for small Fermi

Effects of quantum spin fluctuations
We now turn to a qualitatively different mechanism supporting sizable multi-magnon contribution to low-energy LSEs, based on enhanced quantum spin fluctuations. The latter, contrary to the itinerant-electron scenario, discussed in previous section, can be realized in strictly localized spin systems. The physical motivation for turning attention to spin fluctuations is apparent discrepancy between large ordered moments m † > 2 μ B obtained from first-principle calculations [36,37] and small m † < 1 μ B measured for these compounds [30]. Early attempts to resolve this issue were based on an assumption of proximity to a frustration-induced magnetic quantum phase transition (QPT) in local-moment J 1 -J 2 antiferromagnets [8]. While the relevant QPT in this scheme is of the firstorder, in two-and quasi-two dimensional systems it is presumably very close to second-order [38,39]. Substantial softening of the magnetic fluctuations in iron-pnictides is independently supported by combined first-principle and Ginzburg-Landau studies [40]. Finally, more recent studies, consistent with proximity to an orbitally selective Mott transition [24], usually adopt the picture of coupled itinerant electrons and preformed local moments. In particular, it has been demonstrated [25] that coupling to the band electrons can boost quantum fluctuations of the localized spins in iron-pnictides.
We note that the Heisenberg model (1) is not suitable for the present discussion as (for parameters relevant to BaFe 2 As 2 ) weak quantum fluctuations are signaled by small SWT correction to the ordered moment Δm † ≈ 0.286 μ B . The two-band model of Section 3 does not support sizable spin fluctuations as well. Here we employ an effective description, based on the non-linear σ-model (NLσM) [41,42] given by the action The latter does not require detailed knowledge of the microscopic mechanism of the spin fluctuation enhancement and (in contrast to frustrated J 1 -J 2 models [39,43]) can be readily tuned by the coupling constant g and SW velocities c i to reproduce both the measured ordered moment m † and the low-energy SW dispersion. Since the NLσM (9) with single order parameter n (subjected to the non-linear constraint n 2 = 1) does not include degeneracy between (π, 0) and (0, π) points in the BZ in the high-temperature tetragonal phase, we restrict to low temperatures, where this symmetry is broken [44]. The SW velocities c i and magnon gaps ω α can be obtained directly from the low-energy SW dispersion of BaFe 2 As 2 E 2 α (k) ≈ ω 2 α + c 2 i k 2 i (in this section wavevectors are measured relative to the magnetic Brillouin zone (BZ) center). We get c x ≈ 173 meV, c y ≈ 112 meV, c z ≈ 8.9 meV, c 0 ≡ (c x c y c z ) 1/3 ≈ 56 meV. The empirical damping parameter Γ (k) turns out to be anisotropic, but energy-independent. Here, for simplicity, we set it to a constant Γ (k) ≡ 59.8 meV [26] corresponding to the vicinity of the magnetic zone center (we have checked that the results are not sensitive to the details of the damping at large |k|). While extension of the energy-independent Γ down to the lowest frequencies is difficult to reconcile with expected coherence-driven sharpening of the SWs close to k = 0 [3,18], it might be viewed as a way of modeling the partial transfer of the magnetic spectral weight to the high-energy P-H continuum, inevitable in the presence of coupling to itinerant electrons. The spin anisotropy parameters ω α are determined from the requirement [14] that Imχ 11 0 (ω, 0) and Imχ 22 0 (ω, 0) attain 90% of their maximal values at ω = 18.9 meV and ω = 11.6 meV, respectively.
To treat the NLσM quantitatively, we extend the number of n-field components from physical three to N keeping gN = const. and perform formal large-N expansion which has proven effective near magnetic QPTs (see, e.g., [45]). The NLσM undergoes a QPT from the ordered phase with non-zero expectation value of n to disordered phase as the coupling constant g exceeds the critical value g c given (in the large-N limit) by: where G β (iω n , k) = (ω 2 n + 2Γ |ω n | + c 2 i k 2 i + ω 2 β ) −1 are magnon propagators. The right hand side of equation (10) implicitly depends on the cutoff Λ regularizing the integral at large wavevectors. An estimate for Λ can be obtained from the requirement of the magnetic BZ volume conservation, formally stated as 4/3πΛ 3 ≡ (2π) 3 /2 [46]. This procedure gives Λc 0 ≈ 173 meV, which is consistent with local dynamical susceptibility measurements [47] and is close to other estimates ∼200 meV [41]. The coupling constant g can be obtained from the large-N scaling of the order parameter | n | = 1 − g/g c by identification of | n | with the ordered moment m † normalized to its bare value m † 0 [46]. While m † = 0.87μ B is known from experiments for BaFe 2 As 2 [30], m † 0 is not easily determined due to inherently multi-orbital character of iron-pnictides. The Hund's rule for fully localized Fe 2+ ions yields m † 0 = 4 μ B [47], DFT calculation gives m † 0 ∼ 2.6 μ B [36], frustrated J 1 -J 2 models suggest m † 0 between 2 and 4 μ B [39,41], and anisotropic Heisenberg models (cf. Eq. (1)) are consistent with m † 0 ∼ 1 μ B . We hence compute the longitudinal response in entire regime m † 0 = 0.87-4.0 μ B . The large-N transverse and longitudinal susceptibilities χ αα 0 (α = 1, . . . , N − 1) and χ NN 0 read χ αα 0 (iω n , k) = gc 0 · G α (iω n , k), where is the magnon polarization operator. In Figure 4a the imaginary parts of the analytically continued susceptibilities (11)-(12) at the magnetic BZ center for m † 0 = 3.0 μ B and physical N = 3 are plotted as a function of energy. The longitudinal component Imχ 33 0 (ω, 0) is pronounced, but does not exhibit a sharp threshold behavior as in the case of Heisenberg model (Fig. 1). The reason is that the form of empirical damping yields linear ω dependence of Imχ 11 0 (ω, 0) and Imχ 22 0 (ω, 0) in the low-energy sector, instead of more appropriate threshold behavior. The multi-magnon continuum is sensitive to the structure of transverse excitations, hence it extends to ω → 0 as well as an artifact of the approximation. Since fine tuning of the transverse magnon propagators would simultaneously affect the LSEs, it is unlikely to change substantially the relative intensities I 33 /(I 11 + I 22 + I 33 ). The ratio I 33 /(I 11 + I 22 + I 33 ) for N = 3 is shown in Figure 4b as a function of m † 0 . The regime of m † 0 ≈ 0.87 μ B , where the magnitude of quantum fluctuations is small, has been labeled "classical spin" and is realized by the anisotropic Heisenberg model (1). For a fully localmoment system with m † 0 = 4.0 μ B one gets a large ratio I 33 /(I 11 + I 22 + I 33 ) ≈ 13.6%. Finally, the experimental value I 33 /(I 11 + I 22 + I 33 ) ∼ 5 % is reproduced for m † 0 ≈ 3 μ B which slightly exceeds m † 0 = 2.6 μ B , suggested by first-principle calculations, but is smaller than m † 0 = 4.0 μ B relevant to fully localized spins. This range of bare moment would be qualitatively consistent with the picture of local moments coupled to fluctuating itinerant electrons discussed in reference [25]. As I 33 /(I 11 + I 22 + I 33 ) turns out to be sensitive to the precise value of m † 0 (Fig. 4b), the multi-magnon contribution to the low-energy LSEs deserves further attention within more microscopic models involving both localized and itinerant electrons.

Summary
We have studied the possibility that the LSEs, observed in BaFe 2 As 2 , originate from multi-magnon processes. The latter, contrary to the P-H continuum of an itinerant antiferromagnet (AF), are not energetically constrained by the value of the SDW gap and can be naturally reconciled with optical spectroscopy. Two mechanisms capable of enhancing the multi-magnon contribution to the total low-energy spectral weight have been proposed. The first is based on renormalization of the 2M continuum by virtual electronic transitions between magnetically split bands of an itinerant AF. We have shown that, although enhancement is possible for small FS taking part in the SDW instability, it is not sufficient to explain the data for the range of parameters relevant to BaFe 2 As 2 . The other discussed mechanism relies on spin fluctuations in the vicinity of a magnetic QPT that we have studied within the effective NLσM with Landau damping term. We have calculated fraction of the total low-energy spectral weight consumed by LSEs for BaFe 2 As 2 as a function of bare moment which varies between different classes of microscopic models. It has been demonstrated that the correct order of relative LSE intensities can be reproduced within the spin fluctuation driven scenario, provided that the bare moment is sufficiently large.

Appendix A: Two-band model: robustness to the details of the band structure
In this Appendix, we show that the ratio of renormalization factors R ex 2M /R ex 1M , estimated in Section 3, is insen-sitive to the details of the band structure. This check is essential, because both R ex 1M and R ex 2M strongly depend on the precise form of (l) k , while the minimal model (4) cannot match simultaneously all relevant energy scales of BaFe 2 As 2 for employed dispersion. Specifically, size of the FS, SDW gap Δ, and the magnetic ordering temperature T SDW ≈ 140 K of BaFe 2 As 2 can be reproduced at the mean field level by choosing e 0 /t = 3.3, U 1 = 961 meV, and t = 1000 meV. The scale of J ex eff ≈ 780 meV, however, with no more adjustable quantities, exceeds significantly the experimental in-plane exchange ∼100 meV.
We consider modified dispersions with one more parameter 0 < λ < 1 and discuss the evolution of relevant energy scales and intensity renormalization factors as a function of λ. This choice of˜ (l) k preserves the bandwidth, shape of the FS, and perfect nesting, but the Fermi velocity v F strongly depends on λ. The characteristic energy scales for fixed Δ = 21.3 meV, e 0 /t = 3.3, and t = 1000 meV are shown in Figure A.1a as a function of λ. For λ < 0.9 the exchange coupling is determined by Fermi velocity as J ex eff ≈ v F /2, which is expected of an itinerant AF [15]. The magnetic ordering temperature in this regime is insensitive to λ. The SDW coupling constant U 1 decreases as λ is increased which signals enhancement of the SDW order characterized by Δ/U 1 . A qualitative change occurs when J ex eff approaches 2Δ, where J ex eff no longer follows v F and T SDW becomes λdependent. Figure A.1b shows the dependence of R ex 1M and R ex 2M on λ, both of which undergo sizable enhancement as λ increases. The ratio R ex 2M /R ex 1M is, however, remarkably robust to λ except for a narrow range around λ ≈ 1 (coinciding with the crossover at J ex eff ≈ 2Δ seen in Fig. A.1a). For J ex eff ∼ 100 meV, R ex 2M /R ex 1M is barely altered compared to the analysis with the simple dispersion (l) k .

Appendix B: Two-band model: calculation of the spin susceptibilities
In this Appendix we present derivation of the low-energy form of the transverse and longitudinal spin susceptibilities given by equations (6) and (8).
The action of the two-band model (4) is given by S = dτ¯ k,lc (B.1) Above we have employed a compact notation P ≡ (iω n , p) (ω n denote Matsubara frequencies) and introduced background currents J (l) , coupled linearly to the generalized spin operatorsŜ α lm = S ·c (l) σ α c (m) as i,α (one for each term proportional to U l in H int ), where i runs over the lattice sites, and α = 0 . . . 3, are now used to decouple the interaction Hamiltonian with the use of identity By integrating out fermions and performing a shift φ (l) T N ), we arrive at the effective action (1,2) K,K = φ (1,2) K−K ,α σ α τ 1,0 , andδ K ≡ N T δ K . The Pauli matrices τ α and σ α act on the band (isospin) and spin indices, respectively, and τ 0 ll = σ 0 ll = δ ll . We have also defined Q 0 ≡ (0, q 0 ) and τ 0± ≡ 1 2 (τ 0 ± τ 3 ) to simplify the expressions.
To obtain the lowest order contribution to the longitudinal spin susceptibility for 0 < ω 2Δ, one needs to consider the one-loop diagram shown in Figure 2. The threepoint vertex Γ