Magnon-phonon coupling effects on the indirect K-edge resonant inelastic X-ray scattering spectrum of a 2D Heisenberg antiferromagnet

TL;DR

This study investigates how magnon-phonon interactions influence the RIXS spectrum of a 2D Heisenberg antiferromagnet, revealing phase-dependent spectral features and the potential for experimental observation.

Contribution

It provides a detailed theoretical analysis of magnon-phonon coupling effects on the RIXS spectrum in both AF and CAF phases, highlighting distinct spectral signatures and sensitivity to anisotropy.

Findings

Maximum RIXS intensity around the K point in AF phase

Broad distribution of intensity around the Y point in CAF phase

Phonon effects can produce distinct peaks below and above the bimagnon peak in CAF phase

Abstract

We compute the effects of magnon-phonon coupling on the indirect K-edge bimagnon resonant inelastic x-ray scattering (RIXS) intensity spectrum of a square lattice Heisenberg antiferromagnet. We analyze the effects of competing nearest and next--nearest magnetic and magnon-phonon coupling interaction in the RIXS spectrum, for both the antiferromagnetic (AF) and the collinear antiferromagnetic (CAF) phases of the model. Utilizing the Dyson-Maleev representation of spin operators, the Bethe-Salpeter ladder approximation scheme for the bimagnon interacting channel, and considering the lowest order magnon-phonon-magnon scattering interaction we highlight distinct features in the X-ray spectrum. Considering damping effects, arising due to the presence of phonons, we find that in the AF phase the RIXS intensity spectrum attains a maximum value primarily localized around the K $\left(\pm\frac{\pi}{2}, \pm \frac{\pi}{2}\right)$ - point. For the CAF phase the intensity is broadly distributed with a significant scattering intensity located around the Y $\left(\pm\frac{\pi}{2}, 0\right)$ - point. Furthermore, in the CAF phase for suitable anisotropy, nearest-, and next-nearest neighbor interaction parameters the phonon effects can manifest itself as a distinct peak both below and above the bimagnon peak. Such a feature is in contrast to the antiferromagnetic spectrum where the effect due to the phonon peak was located consistently beyond the bimagnon peak in the high energy end of the spectrum. Additionally, in the CAF phase we find the RIXS bimagnon-phonon spectrum to be more sensitive to anisotropy compared to its antiferromagnetic counterpart. We conclude that the ultimate effect of magnon-phonon effects in the indirect K-edge RIXS spectrum, in both the antiferromagnetic and the collinear antiferromagnetic phase, is an observable effect.

Figures

Tables

Table 1. Table 1: Definition of the channels v n ​ ( 𝐤 ) subscript 𝑣 n 𝐤 v_{\mathrm{n}}(\mathbf{k}) for AF phase and CAF phase

$n$	AF $v_n\left(𝐤\right)$	CAF $v_n^{\prime }\left(𝐤\right)$
1	$u_{𝐤+𝐪}{u}_{𝐤}\cos k_x$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\cos k_x$
2	$u_{𝐤+𝐪}{u}_{𝐤}\sin k_x$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\sin k_x$
3	$u_{𝐤+𝐪}{u}_{𝐤}\cos k_y$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\cos k_x\cos k_y$
4	$u_{𝐤+𝐪}{u}_{𝐤}\sin k_y$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\sin k_x\cos k_y$
5	$u_{𝐤+𝐪}{v}_{𝐤}$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\cos k_x\sin k_y$
6	$v_{𝐤+𝐪}{u}_{𝐤}$	$u_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\sin k_x\sin k_y$
7	$v_{𝐤+𝐪}{v}_{𝐤}\cos k_x$	$u_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }$
8	$v_{𝐤+𝐪}{v}_{𝐤}\sin k_x$	$v_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }$
9	$v_{𝐤+𝐪}{v}_{𝐤}\cos k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\cos k_x$
10	$v_{𝐤+𝐪}{v}_{𝐤}\sin k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\sin k_x$
11	$u_{𝐤+𝐪}{v}_{𝐤}\cos k_x\cos k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\cos k_x\cos k_y$
12	$u_{𝐤+𝐪}{v}_{𝐤}\sin k_x\cos k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\sin k_x\cos k_y$
13	$u_{𝐤+𝐪}{v}_{𝐤}\cos k_x\sin k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\cos k_x\sin k_y$
14	$u_{𝐤+𝐪}{v}_{𝐤}\sin k_x\sin k_y$	$v_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\sin k_x\sin k_y$
15	$v_{𝐤+𝐪}{u}_{𝐤}\cos k_x\cos k_y$	$u_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\cos k_y$
16	$v_{𝐤+𝐪}{u}_{𝐤}\sin k_x\cos k_y$	$u_{𝐤+𝐪}^{\prime }{v}_{𝐤}^{\prime }\sin k_y$
17	$v_{𝐤+𝐪}{u}_{𝐤}\cos k_x\sin k_y$	$v_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\cos k_y$
18	$v_{𝐤+𝐪}{u}_{𝐤}\sin k_x\sin k_y$	$v_{𝐤+𝐪}^{\prime }{u}_{𝐤}^{\prime }\sin k_y$

Equations

$$H=$$

$$+\frac{1}{2}\sum _{i,{\delta }_2}{J}_2(r_{ii+{\delta }_2}){\mathbb{S}}_i\cdot {\mathbb{S}}_{i+{\delta }_2},$$

$$J(r_{ij})=J({\mathbb{R}}_{ij})+(u_i-u_j)\cdot \nabla J(r){|}_{r={\mathbb{R}}_{ij}},$$

$$H^s=$$

$$+E_1\sum _{\mathbb{k}}{\kappa }_{\mathbb{k}}{ϵ}_{\mathbb{k}}({\alpha }_{\mathbb{k}}^{†}{\alpha }_{\mathbb{k}}+{\beta }_{-\mathbb{k}}^{†}{\beta }_{-\mathbb{k}})+H^4.$$

$${\gamma }_1(\mathbb{k})=\frac{\cos k_x+\zeta \cos k_y}{1+\zeta },{\gamma }_2(\mathbb{k})=\cos k_x\cos k_y,$$

$${\kappa }_{\mathbb{k}}=1-\frac{2\eta }{1+\zeta }(1-{\gamma }_2(\mathbb{k})),{\gamma }_{\mathbb{k}}=\frac{{\gamma }_1(\mathbb{k})}{{\kappa }_{\mathbb{k}}},{ϵ}_{\mathbb{k}}=\sqrt{1-{\gamma }_{\mathbb{k}}^2}.$$

$$H^4=$$

$$+\frac{E_1}{SN}\sum _{1,2,3,4}{\delta }_{\mathbb{G}}(1+2-3-4)u_1{u}_2{u}_3{u}_4$$

$$\times (V_{1234}^{(4)}{\alpha }_1^{†}{\beta }_{-4}^{†}{\beta }_{-2}{\alpha }_3+\cdots ),$$

$$A_k$$

$$+A_2\frac{1-{\gamma }_2(k)}{{ϵ}_k},$$

$$A_1$$

$${\Delta }_{A_1}$$

$$A_2$$

$$V_{1234}^{(4)}=$$

$$+{\gamma }_1(4-1)x_1{x}_2-{\gamma }_1(2)x_4-{\gamma }_1(1)x_1{x}_2{x}_4$$

$$-{\gamma }_1(3+4-2)x_3-{\gamma }_1(3+4-1)x_1{x}_2{x}_3]$$

$$+\frac{2\eta }{2(1+\zeta )}[{\gamma }_2(4-2)+{\gamma }_2(4-1)+{\gamma }_2(3-2)$$

$$+{\gamma }_2(3-1)-{\gamma }_2(2)-{\gamma }_2(1)-{\gamma }_2(3+4-2)$$

$$-{\gamma }_2(3+4-1)](x_2{x}_4+{\Phi }_G{x}_1{x}_3).$$

$${\gamma }_1^{\prime }(k)$$

$${\kappa }_k^{\prime }$$

$$A_k^{\prime }$$

$$+A_2^{\prime }\frac{1-{\gamma }_2^{\prime }(k)}{{ϵ}_k^{\prime }},$$

$${\Delta }_{A_1}^{\prime }$$

$$H^{sp}$$

$$+\sum _{i,{\delta }_y}(u_i-u_{i+{\delta }_y})\cdot \nabla J_y(r_{ii+{\delta }_y}){\mathbb{S}}_i^A\cdot {\mathbb{S}}_{i+{\delta }_y}^B$$

$$+\frac{1}{2}\sum _{i,{\delta }_2}[(u_i-u_{i+{\delta }_2})\cdot \nabla J_2(r_{ii+{\delta }_2}){\mathbb{S}}_i^A\cdot {\mathbb{S}}_{i+{\delta }_2}^A$$

$$+(u_{i+{\delta }_x}-u_{i+{\delta }_x+{\delta }_2})\cdot \nabla J_2(r_{i+{\delta }_{x}i+{\delta }_2}){\mathbb{S}}_{i+{\delta }_x}^B\cdot {\mathbb{S}}_{i+{\delta }_x+{\delta }_2}^B]$$

$$u_i=\sum _{q,\lambda }\sqrt{\frac{\hslash }{2Nm{\Omega }_{\lambda }^{ph}(q)}}e(q,\lambda )e^{-iq\cdot {\mathbb{R}}_i}{\phi }_{q\lambda },$$

$$H^{sp}=$$

$$\times [A_{\lambda }({\mathbb{k}}_1,{\mathbb{k}}_2,q){\alpha }_{{\mathbb{k}}_2}^{+}{\alpha }_{{\mathbb{k}}_1}+B_{\lambda }({\mathbb{k}}_1,{\mathbb{k}}_2,q){\beta }_{-{\mathbb{k}}_1}^{+}{\beta }_{-{\mathbb{k}}_2}$$

$$+C_{\lambda }({\mathbb{k}}_1,{\mathbb{k}}_2,q){\alpha }_{{\mathbb{k}}_2}^{+}{\beta }_{-{\mathbb{k}}_1}^{+}+D_{\lambda }({\mathbb{k}}_1,{\mathbb{k}}_2,q){\alpha }_{{\mathbb{k}}_1}{\beta }_{-{\mathbb{k}}_2}].$$

$$g_x(q,\lambda )=\frac{\left|\nabla J_x\right|}{E_1}\sqrt{\frac{\hslash }{2m{\Omega }_{\lambda }^{ph}(q)}}.$$

$${\Lambda }_{\zeta }=g_y(q,\lambda )/g_x(q,\lambda ),{\Lambda }_{\eta }=g_2(q,\lambda )/g_x(q,\lambda ).$$

$${\hat{O}}_q=\sqrt{\frac{2}{N}}\sum _{ij}{e}^{iq\cdot r_i}J(r_{ij}){\mathbb{S}}_i\cdot {\mathbb{S}}_j,$$

$${\hat{O}}_q=\sum _{\mathbb{k}}N(q,\mathbb{k})({\alpha }_{\mathbb{k}+q}^{†}{\beta }_{-\mathbb{k}}^{†}+{\alpha }_{\mathbb{k}}{\beta }_{-\mathbb{k}-q})+\cdots ,$$

$$N(q,\mathbb{k})=$$

$$-{\gamma }_2(q))](u_{\mathbb{k}+q}{v}_{\mathbb{k}}+u_{\mathbb{k}}{v}_{\mathbb{k}+q})+({\gamma }_1(\mathbb{k}+q)+{\gamma }_1(\mathbb{k}))$$

$$\times (u_{\mathbb{k}+q}{u}_{\mathbb{k}}+v_{\mathbb{k}}{v}_{\mathbb{k}+q})\}.$$

$$I(q,\omega )\propto \sum _n{\left|⟨n|{\hat{O}}_q|0⟩\right|}^2\delta (\omega -{\omega }_{n0}),$$

$$iG(q,\omega )={\int }_0^{\infty }dt{e}^{i\omega t}⟨0|\mathcal{T}{\hat{O}}_q^{†}(t){\hat{O}}_q(0)|0⟩.$$

Full text

PACS number(s)

78.70.Ck, 75.25.-j, 75.10.Jm

I Introduction

Improved X-ray instrumentation resolution coupled with advanced X-ray synchrotron radiation sources have established resonant inelastic X-ray scattering (RIXS) as a significant experimental tool to explore condensed matter systems Ament et al. (2011a); Dean (2015). Experimental  Hill et al. (2008); Braicovich et al. (2009); Schlappa et al. (2009); Braicovich et al. (2010); Guarise et al. (2010); Glawion et al. (2011); Le Tacon et al. (2011); Kim et al. (2012); Bisogni et al. (2012); Dean et al. (2012, 2013); Ulrich et al. (2009); Schlappa et al. (2012); Lee et al. (2013); Monney et al. (2013); Yuan et al. (2017), theoretical van den Brink (2007); Forte et al. (2008a); Nagao and Igarashi (2007); Forte et al. (2008b); Jia et al. (2012); Luo et al. (2014, 2015); Marra et al. (2012, 2013); Klauser et al. (2011); Pakhira et al. (2012), and computational  Johnston et al. (2010); Zhou et al. (2013); Jia et al. (2012); Klauser et al. (2011); Forte et al. (2011); Lee et al. (2013) approaches have investigated RIXS across a wide variety of systems in various dimensions for a range of elementary excitations at different X-ray edges. While the effect of magnon-magnon interactions on the indirect K-edge RIXS spectrum has been investigated both experimentally and theoretically in the antiferromagnetic (AF) and the collinear antiferromagnetic (CAF) phases van den Brink (2007); Nagao and Igarashi (2007); Forte et al. (2008b); Jia et al. (2012); Luo et al. (2014), till date there is no study (theoretical or experimental), which exculsively investigates the role of many body magnon-phonon interaction on the indirect K-edge bimagnon RIXS intensity spectrum. At the K–edge spin angular momentum is conserved in the indirect RIXS process due to the lack of spin–orbit coupling in the 1s electron. Thus the double spin–flip bimagnon excitation is the leading process at the K–edge. For the higher angular momentum shells $(L,M,N,...)$ with finite orbital momentum single spin flip excitations are allowed Ament et al. (2011a).

As a probe RIXS has a high degree of senstivity towards the local environment. Thus it is a natural question to ask - How can phonons which produce lattice modulations (local vibrations) affect the magnetic RIXS spectrum? How can a multimagnon RIXS excitation spectrum, such as that of a bimagnon, be affected by phonons? The interplay of phonons with bimagnons offer the potential to uncover physical phenomena which has been overlooked till date. With next generation beamlines being constructed globally and experimental initiatives likely to probe phonon effects in correlated materials, answers to the above questions are imperative and timely. The theoretical study in this article offers insight on the key experimentally observable signatures which delineate magnon-magnon and magnon-phonon interaction effects in the indirect K-edge RIXS spectrum.

Past investigation on spatial anisotropy and significant frustration within the square lattice Heisenberg magnet has led to the prediction of a two-peak bimagnon structure Luo et al. (2014). The proposed two-peak structure was a consequence of the bimagnon spectrum’s sensitivity to microscopic magnetic interactions. But, in real materials lattice vibrations do matter. Thus, a realistic theoretical model which provides a true depiction of the materials under investigation with a comprehensive account of magnon-phonon coupling is called for.

The coupling between magnetic and lattice degrees of freedom can generate novel physical phenomena. For example, it can effect electronic degrees of freedom Moser et al. (2015); Ament et al. (2011b) induce multiferroic order Toth et al. (2016), create magnon-phonon excitation effects in Raman spectroscopy Johnston et al. (2010), and have an effect on thermal conductivity Chernyshev and Brenig (2015) and optical conductivity Grüninger et al. (2000). Past theoretical studies on magnon-phonon interaction in quantum Heisenberg magnets have alluded to the fact that at low but finite temperatures phonons do play a role on influencing magnetic interaction. This fact has been especially studied within the context of experimental and theoretical studies of spin-phonon interaction on the Raman spectra of Heisenberg antiferromagnets. In these studies it was found that significant broadening effects dominate both the single magnon and two-magnon line shape. Theoretical calculations have attributed the anomalous broadening of the two-magnon spectrum in cuprates to phonon effects Lorenzana and Sawatzky (1995a, b); Sänger (1994); Saenger (1995); Devereaux et al. (1995); Lee and Min (1997).

Unfortunately, Raman spectroscopy is a zero wave vector probe Devereaux and Hackl (2007). Therefore it limits the amount of physical information that can be extracted. But, RIXS is not. The high energy X–ray photons in the RIXS experiments allow for large transferred momenta, with the zero wavevector reproducing the Raman response. Thus, it is appropriate to consider RIXS to explore the full energy-wavevector range to study the effects of magnon-phonon coupling on the bimagnon excitation spectrum. Recently, there has been some experimental  Peng et al. (2015); Yavaş et al. (2010); Moser et al. (2015) and theoretical  Ament et al. (2011c); Devereaux et al. (2016); Yavaş et al. (2010) studies devoted to the study of electron-phonon coupling and its effect on the RIXS spectrum. We note, that our study is different from the existing ones since we are primarily concerned with the role of magnon-phonon interaction on the bimagnon excitation spectrum.

In this article, we compute the effects of magnon-phonon coupling on the indirect K-edge bimagnon resonant inelastic x-ray scattering (RIXS) intensity spectrum of a square lattice Heisenberg antiferromagnet. The magnetic model theoretically investigated includes both the nearest-neighbor (nn) and next-nearest neighbor (nnn) magnetic and magnon-phonon interaction effects. Since the Heisenberg magnetic model is considered upto the nnn interaction (see Eq. 1), to be consistent in our theoretical formulation we include the magnon-phonon coupling beyond the nn interaction. Using Dyson-Maleev representation of spin operators, the Bethe-Salpeter ladder approximation scheme for the bimagnon interacting channel, and considering the lowest order magnon-phonon-magnon scattering interaction we elucidate the distinct features in the RIXS X-ray spectrum.

Our analysis of the RIXS spectrum, performed for a wide range of model parameters in both the antiferromagnetic and the collinear antiferromagnetic phases of this model, suggest several contrasting behavior in the antiferromagnetic and the collinear antiferromagnetic phase. We consider damping effects in our calculation due to the presence of longitudinal acoustic phonons in our model. We find that in the AF phase the RIXS intensity spectrum attains a maximum value primarily localized around the K $\left(\pm\frac{\pi}{2},\pm\frac{\pi}{2}\right)$ - point. Within the nearest-neighbor model the system is weakly sensitive to the presence of magnon-phonon interactions. For most parameter choices the feature developed is a shoulder in the RIXS spectrum. But, inclusion of the nnn magnetic and magnon-phonon coupling within the isotropic model leads to a splitting of the peak. In contrast, for the CAF phase in the isotropic model the intensity is broadly distributed with a significant scattering intensity located around the Y $\left(\pm\frac{\pi}{2},0\right)$ - point. The rest of the spectral weight appears along the K to M $\left(\pm\pi,0\right)$ path in the Brillouin zone (BZ). Furthermore, in the CAF phase for suitable anisotropy, nearest-, and next-nearest neighbor interaction parameters the phonon effects can manifest itself as a distinct peak both below and above the bimagnon peak. Such a feature is in contrast to the AF spectrum where the effect due to the phonon peak was located consistently beyond the bimagnon peak in the high energy end of the spectrum. Within the anisotropic model the AF RIXS spectrum is merely broadened without any special peak or shoulder development. However, in the CAF phase we find the RIXS bimagnon-phonon spectrum to be more sensitive to anisotropy compared to its antiferromagnetic counterpart. The final RIXS spectra is a result of intricate many body magnon-magnon interactions, influenced by the effect of many-body mangon-phonon interactions.

This article is organized as follows. In Sec. II.1 we introduce the nn and nnn Heisenberg Hamiltonian including the effects of phonons. In Sec. II.2 we write down explicitly the contribution arising from the spin-phonon coupling. In Sec. III we state the RIXS operator, the bimagnon Green function, the corresponding Bethe-Salpeter equation, and the phonon Green function contributions. In Sec. IV we present and discuss our results on the effect of damping and magnon-phonon coupling (Sec. IV.1), phonon contribution to bimagnon RIXS spectrum (Sec. IV.2), frustration and magnon-phonon coupling effects on the AF phase bimagnon phonon spectrum (Sec. IV.3), CAF phase RIXS spectrum (Sec. IV.4), and anisotropy effects in both the AF and CAF phase (Sec. IV.5). Finally, in Sec. V we provide our conclusions.

II Model

II.1 Heisenberg Hamiltonian

The frustrated $J_{1}-J_{2}$ model Yao and Carlson (2008); Uhrig et al. (2009); Goswami et al. (2011) and its spatially anisotropic $J_{x}-J_{y}-J_{2}$ version  Majumdar (2010); Applegate et al. (2010); Schmidt et al. (2010) can support both the $(\pi,\pi)$-AF and the $(\pi,0)$-CAF phase. There are several excellent material realizations of these models in cuprates Dagotto (1994), pnictides Dai (2015); Stewart (2011), and complex vanadium oxide compounds  Tsirlin and Rosner (2009). For our purpose we generalize the model to include the effect of spin-lattice interaction by expanding the spatially dependent exchange interactions $J_{x,y}(\mathbb{r}_{ij})$ and $J_{2}(\mathbb{r}_{ij})$. The presence of lattice vibrations (phonons) cause the ionic distances to vary, which in turn are expected to induce magnon-phonon interactions  Dixon (1980). We write the spin S=1/2 anisotropic $J_{x}-J_{y}-J_{2}$ Heisenberg model on a two dimensional square lattice as

[TABLE]

where $\mathbb{\delta_{2}}=\mathbb{\delta_{x}}+\mathbb{\delta_{y}}$. Assuming first order deviations from the lattice equilibrium positions, the super-exchange coupling $J(\mathbb{r}_{ij})$ can be expanded as Sänger (1994); Lee and Min (1997)

[TABLE]

where $\mathbb{r}_{i}=\mathbb{R}_{i}+\mathbb{u}_{i}$, $\mathbb{R}_{i}$ is the equilibrium position of the site i, and $\mathbb{u}_{i}$ is the displacement operator due to the lattice vibration. Henceforth, we set $J_{x}(\mathbb{R}_{ij})=J_{x}$, $J_{y}(\mathbb{R}_{ij})=J_{y}$ and $J_{2}(\mathbb{R}_{ij})=J_{2}$. The phonon lattice deviations can be quantized resulting in a Hamiltonian which consists of a pure spin contribution Hs and one with a magnon-phonon contribution Hsp. We next analyze the model within the standard Dyson-Maleev spin wave theory approach outlined in Appendix A. Note, the choice of representation Holstein-Primakoff or Dyson-Maleev does not affect the results of our paper. In the AF phase we have

[TABLE]

In the above equation we have $E_{0}=-\frac{N}{2}J_{x}S^{2}z(1+\zeta)(1-\frac{2\eta}{1+\zeta})$, where N is the total number of sites and the coordination number $z=2$. We also introduce the interaction ratios $\zeta=J_{y}/J_{x}$ and $\eta=J_{2}/J_{x}$, and $E_{1}=J_{x}Sz(1+\zeta)$ in the AF phase. In addition, we define the following functions

[TABLE]

The higher order term $H^{4}$ is given by the expression

[TABLE]

where the $u_{k}$ and $v_{k}$ coefficients arise in the Bogoliubov transformation with $v_{k}=-x_{k}u_{k}$. The momentum labels $k_{1},k_{2},\ldots$ are abbreviated as $1,2,\ldots$. We invoke the conservation of momentum rule upto the reciprocal lattice vector $\mathbb{G}$ with $\delta_{\mathbb{G}}(1+2-3-4)$. The constant and the quadratic terms in $H^{4}$ arising from normal ordering procedure are known as Oguchi corrections. The coefficients are

[TABLE]

Note, the $A_{k}$ equation is different from the previously reported expression Igarashi (1992); Igarashi and Nagao (2005); Majumdar (2010); Majumdar et al. (2012). We find an additional contribution $\Delta_{A_{1}}$ that was previously ignored (please see Appendix A). At the bimagnon RIXS spectrum level the presence or absence of this Oguchi correction term does not change the results quantitatively or qualitatively. However, with the inclusion of anisotropic interactions and magnon-phonon coupling it is important to consider an accurate expression. The AF phase magnon-phonon RIXS spectrum is mildly affected, but, the CAF phase spectrum is unaffected. The quartic interaction vertex in Eq. 6, relevant to our calculation, is given by the expression

[TABLE]

where $\Phi_{\mathbb{G}}=e^{iG_{x}}$. In the CAF phase, $H^{s}$ and $H^{4}$ have the same form as in the AF phase with the new redefined coefficients

[TABLE]

The analytical expressions for $E^{\prime}_{1},A^{\prime}_{1},A^{\prime}_{2}$, and $V^{\prime(4)}_{1234}$ can be obtained by using the replacement $\zeta\leftrightarrow 2\eta$. The same replacement will also generate the coefficients $\gamma^{\prime}_{1}(k),\gamma^{\prime}_{2}(k),\kappa^{\prime}_{k}$,and $\epsilon^{\prime}_{k}$ in $A_{1},A_{2}$, and $V^{(4)}_{1234}$.

II.2 Magnon-phonon Hamiltonian

In this section we focus on the magnon-phonon Hamiltonian contribution $H^{sp}$, generalized to include both the effects of spatial anisotropy and further neighbor interactions. Introducing the Taylor expansion of the exchange coefficients as mentioned earlier, we can write down an expression for the magnon-phonon Hamiltonian in the AF phase as

[TABLE]

where the quantized displacement operator $\mathbb{u}_{i}$ expression is given by

[TABLE]

where $\Omega_{\lambda}^{ph}(\mathbb{q})$ is the dispersion of phonon in branch $\lambda$, $\mathbb{e}(\mathbb{q},\lambda)$ is the phonon’s polarization vector, $\varphi_{\mathbb{q}\lambda}$ is the phonon operator, and m is the reduced ionic mass. To recast Eq. 10 into its spin wave version we use the standard Dyson-Maleev transformation, followed by a Fourier transformation, and a subsequent Bogoliubov transformation to obtain the bosonized magnon-phonon Hamiltonian expression as

[TABLE]

In the above we have introduced the notation

[TABLE]

The spin-phonon coupling coefficients $A_{\lambda}$,$B_{\lambda}$,$C_{\lambda}$, and $D_{\lambda}$ in Eq. 12 are given in Appendix B for both the AF and the CAF phase. Furthermore, in the following discussion we introduce the magnon-phonon coupling ratios for anisotropy and nnn, respectively as

[TABLE]

It is possible to estimate a value for the nn magnon-phonon coupling from experimental data Knoll et al. (1990). For example, in cuprates such as RBa2Cu3O6 (R = Eu, Y), the change in the exchange energy $\nabla J_{x}$ for the Cu-O bond ranges between 2500 - 6000 cm*-1Å-1*. The exchange constant itself varies between 800 - 1000 cm*-1* (96 - 120 meV). With $\hbar\Omega_{ph}$ in a 20 - 40 meV interval and the reduced mass of a Cu-O bond system equal to 2.13 $\times$ 10*-23* g, we find that gx lies between 0.076 - 0.41. For our calculations we have chosen g${}_{x}=0.28$. Reliable experimental data or theoretical estimates on further neighbor magnon-phonnon couplings are either rare to find or difficult to obtain. Thus, for our purposes we make an educated guess of the physically reasonable ratios to simulate the RIXS spectrum. We hope this provides further motivation, both experimentally and theoretically, to investigate the physics of the further neighbor magnon-phonon coupled quantum magnet systems.

III Bimagnon and Magnon-Phonon RIXS

In this section we present the expression for the bimagnon RIXS scattering operator generalized to the cases of spatially dependent exchange interaction. Utilizing the standard definition of the bimagnon RIXS operator, valid within the ultrashort core-hole lifetime (UCL) expansion  van den Brink and van Veenendaal (2005, 2006); Ament et al. (2007), we have

[TABLE]

In its bosonized form the operator reads as

[TABLE]

where we have

[TABLE]

$N^{\prime}(\mathbb{q},\mathbb{k})$ in the CAF phase can be obtained with the replacement $\zeta\leftrightarrow 2\eta$ along with the corresponding coefficients $\gamma^{\prime}_{1}(k),\gamma^{\prime}_{2}(k),u^{\prime}_{k}$, and $v^{\prime}_{k}$. In the following discussion, the energy is in units of $E_{1}$. The scattering intensity is given by

[TABLE]

where $|n\rangle$ represents the excited states in the RIXS intermediate process, and $|0\rangle$ is the ground state. The Fourier transform of the zero temperature time-ordered Green’s function is given by

[TABLE]

Then the scattering intensity can be expressed from the Green’s function as

[TABLE]

The perturbations for our problem are magnon–magnon interaction and magnon–phonon–magnon interaction. The Feynman diagrams are shown in Fig. 2. The Green’s function $G(\mathbb{q},\omega)$ consists of a bimagnon part $G^{\text{bm}}(\mathbb{q},\omega)$ and a magnon-phonon-magnon part $G^{\text{m-ph-m}}(\mathbb{q},\omega)$. The bimagnon part is given by Nagao and Igarashi (2007); Luo et al. (2014)

[TABLE]

where the interacting two magnon correlation function is defined as

[TABLE]

The bimagnon contribution include the effects of magnon-magnon interaction Powalski et al. (2015); Hamer et al. (1992) can be be solved exactly using the Bethe-Salpeter equation. The interaction $u_{1}u_{2}u_{3}u_{4}V^{(4)}_{1234}$ in the spin Hamiltonian can be decomposed into 18 channels (see Appendix C)

[TABLE]

where $\hat{\nu}(\mathbb{k})$ has dimensions of $1\times 18$, $\hat{\Gamma}$ is a matrix of dimension $18\times 18$. The expressions for each of these quantities in the AF and CAF phase are given in Appendix C. After summing the ladder diagrams exactly Luo et al. (2014); Davies et al. (1971); Canali and Girvin (1992); Nagao and Igarashi (2007), the two magnon Green’s function can be expressed as a combination of several matrix products as

[TABLE]

where $\hat{1}$ is a unit matrix of dimensions $18\times 18$ , and we define the non-interacting Green function and the non-interacting polarization propagtor, respectively, as

[TABLE]

and

[TABLE]

The leading order of the magnon-phonon-magnon part (the zeroth order part is already included in bimagnon part), $G^{\text{m-ph-m}}$ takes the form

[TABLE]

In the above Green functions we have introduced the variable $\Delta_{\mathbb{k}}$ as magnon damping due to magnon-phonon-magnon interaction. We use the energy dispersion $\omega_{k}=\kappa_{k}\epsilon_{k}+A_{k}/2S$ in units of $E_{1}$.

IV Bimagnon RIXS Spectra Results

In this section we systematically investigate the effect of damping, magnon-phonon interaction, magnon-magnon interaction, and anisotropy on the K-edge bimagnon phonon indirect RIXS spectrum. We compute the RIXS spectrum along the BZ traversing the path $\Gamma$:(0,0) $\rightarrow$ K$:\left(\frac{\pi}{2},\frac{\pi}{2}\right)$ $\rightarrow$ M:$\left(\pi,0\right)$ $\rightarrow$ $\Gamma$:(0,0).

IV.1 Damping and magnon-phonon coupling effects

We consider longitudinally polarized acoustic phonons, $\mathbb{e}(\mathbf{q},\lambda)\parallel\mathbf{q}$, with a dispersion given by  Sänger (1994)

[TABLE]

where

[TABLE]

To keep the discussion and analysis of the model tractable we first investigate the isotropic version of the model. In spatially isotropic limit of the model, $\zeta=1$, we set $J_{x}=J_{y}=J_{1}$ and $g_{x}=g_{y}=g_{1}$. Thus, the overall energy scale is given by $E_{1}=2J_{1}Sz$. The choice of $(\zeta,\eta)$ parameters are dictated by the magnetization phase diagram of the $J_{x}-J_{y}-J_{2}$ model  Majumdar (2010).

In the absence of a magnetic field magnons in a 2D square lattice are not damped Kopietz (1990); Ty and Halperin (1990). However, with the inclusion of magnon-phonon interaction the lifetime of the magnons are affected Sänger (1994); Saenger (1995). In the regime where the sound velocity is less than the magnon velocity, there exists spontaneously occuring decay processes where a magnon can decompose into a magnon and a phonon. Even at zero temperature, quantum fluctuations arising from magnon-phonon interactions can damp spin wave excitations. Therefore, presence of spin-lattice couplings can have an effect on the RIXS spectrum. Thus, a proper treatment of the magnon-phonon interaction should consider the damping effect.

In Fig. 3 we display our calculations of the energy renormalization and damping effects within our model. In Figs. 3(a) - 3(b) we display the undamped RIXS intensity spectrum at the $M$- and the $K$- point of the BZ in the AF phase. Note, the presence of the van-Hove singularity like sharp peaks in the non-interacting case (dashed lines). These singular structures disappear when we include the two-magnon interaction (red solid lines) Nagao and Igarashi (2007). The scattering of two magnons indeed changes the structure of the response function implying that it is no longer a simple product of the RIXS matrix element and the density of states.

In Fig. 3(c) we study damping on the RIXS spectrum. With increasing damping strength the spectrum is broadened, more so in the low energy regime where the spectrum height decreases with increasing damping strength. The high energy tail of the intensity pattern is not much affected. While damping could potentially arise from various microscopic mechanisms within a 2D square lattice problem, in this article we mainly focus on the effect of magnon-phonon interaction. Also for simplicity, we consider a phenomenological phonon induced constant damping to describe the imaginary part of the self energy. Introducing a constant damping to describe phonon effects is inspired by some previous theories of two magnon Raman spectra in cuprates Weber and Ford (1989); Knoll et al. (1990). Thus, in our calculations we set the damping parameter $\Delta_{\mathbf{k}}=0.1E_{1}$. In Fig. 3(d) we display the trend on the RIXS spectrum that would arise when the nn magnon-phonon coupling $g_{1}$ is increased in strength. Clearly, beyond a threshold magnon-phonon coupling strength a shoulder peak develops. However, the mere inclusion of the magnon-phonon coupling $g_{1}$ is not strong enough to induce any further features. In Figs. 3(e) - 3(f) we display the combined effects of damping and magnon-phonon coupling. At this stage, the nnn interaction is set to zero.

IV.2 Phonon contribution to RIXS spectrum

The origin of magnon-phonon interaction is dynamical in nature. Thus summing up an infinite set of diagrams, as in the bimagnon RIXS intensity case, can pose a serious calculation challenge. To proceed with the calculation we note that higher order terms generated by the magnon-phonon interaction decrease rapidly. Thus, akin to the celebrated Migdal theorem used within the context of electron-phonon scattering in superconductivity we consider only the leading order magnon-phonon diagram, see Fig. 2(b).

In Fig. 4 we showcase the effects of $g_{1}$ on the magnon-phonon-magnon intensity as the nnn frustration parameter is varied from 0.1 to 0.3. It appears that for a given nnn interaction strength, wthin the isotropic model, the magnon-phonon interaction causes a rearrangement of the spectral strength. While at the $M-$ point the reshuffling of the spectral weight is quite prominent, that at the $K$- point is minimally affected. However, for both the locations in the BZ the spectral intensity scales in proportion to the magnon phonon coupling.

In Figs. 5(a) and  5(b) we investigated the magnon-phonon-magnon RIXS intensity spectrum as a function of $\eta$. It is observed that with increasing interaction the spectra has a downshift. Note, the miniscule unphysical negative contributions in the intensity is an artifact of the exclusion of the zeroth order contribution from the magnon-phonon-magnon vertex function Lee and Min (1997). The final total RIXS intensity which includes the contribution from the phonon induced ladder interactions restores this term and naturally yields an overall physical positive RIXS intensity spectrum (see Fig. 6 and Fig. 7). The observed spectral downshift is a characteristic feature that is also noticed in neutron  Coldea et al. (2001); Hayden et al. (2004) and Raman spectrum  Chen et al. (2011) with frustrated interactions. Our present calculations confirm that such an effect can also occur even within the magnon-phonon-magnon channel. Interestingly, while the trend towards downshift itself is robust, the peak-dip-peak structure appears to be dependent on the BZ location. At the characteristic $(\pi,0)$ (M point in BZ) the intensity pattern is complimentary to that observed at the K-point. In Figs. 5(c) and  5(d) the effect of $\Omega^{ph}$ on the magnon-phonon-magnon intensity at M point and K point is displayed to give a sense of how the RIXS spectrum may be affected. We observe within our choice of parameters $\Omega^{ph}$ does not affect the intensity spectrum.

From Figs. 5(e) - 5(f) it is clearly evident that with increasing nnn interaction $\eta$ the phonon RIXS spectrum undergoes a shift in the spectral weight. While most of the weight is localized around the K point, with enhanced nnn interaction some of the weight disperses towards the BZ corner. The 2D plot suggests that the phonon effects are primarly localized along the $(\pi/2,\pi/2)\rightarrow(\pi,0)$ line in BZ. However, with the inclusion of the nnn phonon interaction the spectral leaking is subdued, with the spectrum becoming more localized at the $K$ point, see Figs. 5(g) - 5(h).

IV.3 Effect of frustration and magnon-phonon coupling

In Fig. 6 we show the combined effect of including nnn interaction and phonon coupling in the calculation. From Fig. 6(a) - 6(b) it is clearly evident that nnn interaction introduces a downshift of the spectral weight. In Figs. 6(c) - 6(d) we show the dependence on the RIXS spectrum as the magnon–phonon coupling is introduced. Comparing the line plots for the $\eta=0.3$ case with that in Figs. 6(c) - 6(d) we notice that including magnon-phonon coupling causes a peak development, in addition, to introducing further broadening. The broadening effect is more prominent at the M- point. We notice from Figs. 6(e) - 6(h) that the effect of next-nearest neighbor magnon-phonon coupling is minimal at the M point, however, has a significant effect at the K- point. Further pronounced effects of the nnn magnon-phonon coupling contribution is clearly visible in Fig. 6(h). A 2D RIXS intensity pattern tracking the evolution of the total interacting RIXS peak development is shown in Fig. 7. We observe that the most prominent signal is at the K- point with satellite peaks due to magnon-phonon coupling developing as the strength of the nnn contribution is enhanced. The effect of the nn magnon-phonon coupling appears to be minimal within this model. This fact is evident by observing the relatively unchanged RIXS spectrum features when the magnon-phonon coupling is changed.

IV.4 RIXS intensity in CAF phase

The 2D quantum Heisenberg AF can support a collinear AF phase for relevant magnetic interaction parameters  Majumdar (2010). While the AF phase provides information on the parent magnetic compunds of cuprates, the CAF phase magnets are relevant in understanding how pnictide superonductivity can arise. To provide a comprehensive understanding of the model under study, in Fig. 8, we display our findings for the CAF phase. Compared to Fig. 7 for the AF phase where the RIXS intensity peaks at the K$:\left(\frac{\pi}{2},\frac{\pi}{2}\right)$ point, in the CAF phase the RIXS spectral intensity weight attains its maximum value near $\left(\frac{\pi}{2},0\right)$. Futher comparison of the magnon-phonon-magnon RIXS plots in the two phases, Fig. 8(b) and Figs. 5(e) - 5(h), suggest that the magnon-phonon interaction creates additional satellite structures in the CAF phase for suitable interaction parameters. The spectral intensity of the magnon-phonon crorrelation is shown in Fig. 8(b) where we notice that the intensity is spread over a much wider region in the BZ compared to the AF phase. The additional ripples in the CAF phase spectrum arises from the magnon-phonon-magnon spectrum. Systematic studies of the damping effect are displayed in Fig. 9. Similar to the AF phase, with increasing damping the RIXS bimagnon curves are broadened out. For the magnon-phonon-magnon plots, strong damping tends to suppress the satellite structures in the spectrum more at the K point, compared to the one at the Y point.

IV.5 Anisotropy, frustration, and magnon-phonon coupling

In Fig. 10, we compare and contrast the effects of spatial anisotropy between the AF and the CAF phase to reveal subtle differences between the RIXS response. A prominent unique two-shoulder peak is seen to develop in the CAF phase, which for the AF phase has typically been absent for the physical parameter region that was investigated. While in the AF phase the effects of anisotropy overpower the magnon-phonon effects, in the CAF phase the phonon peak structures are evident. These differences are crucial in distinguishing between the two different types of magnetic ordering. We also point out that at the Y-point the CAF phase spectrum, the magnon-phonon coupling introduces peaks both at the low and high energy spectrum. This feature is a consequence of the multi-satellite spectral feature of the magnon-phonon spectrum as seen in Fig. 9.

V Conclusion

Although the theory of magnon-magnon interactions on RIXS spectrum is well established there has been no thorough investigation on the role of phonons on the bimagnon magnetic spectrum of RIXS in the insulating 2D quantum Heisenberg magnet. In particular, microscopic investigations of the interaction between spin and lattice degrees of freedom and its effects on the bimagnon excitation is completely missing. In this paper, we provide a comprehensive theory of magnon-phonon interaction effects in both the AF and CAF phase at the K-edge. We compute the RIXS spectrum including the effects of damping, magnon-phonon coupling both at the nn and nnn level in the AF and the CAF phase. A detailed comparison between the two spectra reveal stark differences in the location of the maximum intensity points, the nature of the effect of phonon on the two-peak structure, and the distribution of the indirect RIXS spectral weight. From the perspective of experiments, the local sensitivity of magnon-phonon correlation on the bimagnon spin dynamics should be clearly evident in the nature of the varying spectral shapes. The unique differences manifest when anisotropy and frustration is taken into account. The final outcome of the magnon-phonon effects in the indirect K-edge RIXS spectrum, in both the AF and the CAF phase, is an experimentally observable feature. Based upon our analysis we infer that the spectrum is a delicate balance between competing nearest and next-nearest neighbor magnon and magnon-phonon coupling strengths. Finally, we hope our work will inspire future experimental investigation on the hitherto unexplored consequences of phonon effects on bimagnon RIXS dynamics at the K-edge.

Acknowledgements.

We thank Mark Dean for useful discussions. T.D. acknowledges invitation, hospitality, and kind support from Sun Yat-Sen University. T. D. acknowledges funding support from Augusta University Scholarly Activity Award and from Sun Yat-Sen University Grant No. OEMT–2017–KF–06. Z. X. and D. X. Y. are support by NKRDPC-2017YFA0206203, NSFC-11574404, NSFC-11275279, NSFG-2015A030313176, Special Program for Applied Research on Super Computation of the NSFC-Guangdong Joint Fund (the second phase), Leading Talent Program of Guangdong Special Projects.

Appendix A Dyson Maleev, Fourier and Bogoliubov transformation

The standard two sublattice Dyson-Maleev representation used to bosonize the spin operator is given by

[TABLE]

for the sublattice A, and

[TABLE]

for the sublattice B. In the above $a(a^{\dagger})$ represents the bosonic annihilation (creation) operator on sublattice A and $b(b^{\dagger})$ represents the bosonic annihilation (creation) operator on sublattice B. Introducing the Fourier transform of bosonic operator we have

[TABLE]

where N is the total number of sites. Using the Bogoliubov transformations

[TABLE]

we diagonalize the quadratic terms in $H^{s}$. The expressions for $u_{\mathbb{k}}$ and $v_{\mathbb{k}}$ are given by

[TABLE]

Further useful relations used in the calculation include $u_{k}^{2}+v_{k}^{2}=1/\epsilon_{k}$ ad $u_{k}v_{k}=-\gamma_{k}/(2\epsilon_{k})$. The additional $\Delta_{A_{1}}$ Oguchi correction term $A_{k}$ ($A^{{}^{\prime}}_{k}$) in the AF phase is given by

[TABLE]

and in the CAF phase by

[TABLE]

Appendix B Magnon-phonon vertex

The spin-phonon vertex $A_{\lambda},\,B_{\lambda},C_{\lambda},\,D_{\lambda}$ are short for $A_{\lambda}(k_{1},k_{2},q),\,B_{\lambda}(k_{1},k_{2},q),C_{\lambda}(k_{1},k_{2},q),\,D_{\lambda}(k_{1},k_{2},q)$. Their analytical expressions are given by

[TABLE]

For the AF phase we define

[TABLE]

where $\Delta_{x}(k),\Delta_{y}(k),\Delta_{2}(k)$ are short for $\Delta_{x}(q,\lambda,k),\Delta_{y}(q,\lambda,k)$ and $\Delta_{2}(q,\lambda,k)$ respectively,

[TABLE]

and

[TABLE]

For the CAF phase, we define

[TABLE]

where $\Delta_{x}(k),\Delta_{y}(k),\Delta_{2}(k)$ are the same as in AF phase. But,

[TABLE]

Appendix C Vertices and $\Gamma$ matrix

The vertex $u_{1}u_{2}u_{3}u_{4}\,V^{(4)}_{1234}$ can be transformed into a separable form (see Eq. 23) with 18 channels and a q-dependent $\Gamma(q)$ matrix. The definition of these channels are given in Table 1 Luo et al. (2015).

The non-zero matrix elements of $\Gamma(q)$ in AF phase are given by

[TABLE]

where

[TABLE]

The non-zero matrix elements of $\Gamma^{\prime}(q)$ in CAF phase are given by

[TABLE]

where

[TABLE]