Tunable phonon-driven magnon–magnon entanglement at room temperature

We report the existence of entangled steady-states in bipartite quantum magnonic systems at elevated temperatures. We consider dissipative dynamics of two magnon modes in a bipartite antiferromagnet, subjected to interaction with a phonon mode and an external rotating magnetic field. To quantify the bipartite magnon–magnon entanglement, we use entanglement negativity and compute its dependence on temperature and magnetic field. We provide evidence that the coupling between magnon and phonon modes is necessary for the entanglement, and that, for any given phonon frequency and magnon–phonon coupling rate, there are always ranges of the magnetic field amplitudes and frequencies for which magnon–magnon entanglement persists at room temperature.


introduction
Entanglement is a central concept of quantum physics that has developed into a key resource in quantum technology [1].It plays a role in quantum information processing, opens a way towards intercept-resilient quantum communications, and enables increased sensitivity in quantum metrology.If it were possible to overcome the thermal noise and prepare robust entangled quantum states at ambient conditions, the corresponding implications for the future quantum devices would be significant [2,3,4,5,6,7,8,9].
Here, we report evidence for equilibrium phonon-driven magnon-magnon entanglement in bipartite antiferromagnetic materials.We compute entanglement negativity [27] of a pair of magnon modes coupled to a lattice vibration, in the presence of an offresonant external rotating magnetic field.Our aim is to study the stability of bipartite entanglement between magnons to thermal noise.We search for the regime where the magnonic degrees of freedom remain entangled with each other rather than with the environment, and hence can be used as coherent channels in quantum information protocols.We obtain a steady-state entanglement diagram in the temperature and magnetic field plane, which shows the possibility of having magnon-magnon entanglement at room temperature.We analyze the magnon entanglement stability against dissipation caused by magnon-phonon coupling and observe that the two-mode magnon system indeed permits high-temperature entanglement at any value of the magnetic field provided that the magnetic field frequency is adjusted to the proper range.In turn, the existence of crystalline and synthetic antiferromagnetic materials, such as NiO [28], 2D Ising systems like MnPSe 3 [29,30], YIG-based synthetic antiferromagets (SAFs) [31], and perovskite manganites [32,33], particularly SrTcO3 [34], CaTcO3 [35], BiFeO3 [36] that have high Neel temperatures, provide a space for further study of high-temperature entanglement in quantum magnonics.To pursue our analysis, we assume a pair of magnon modes in a bipartite antiferromagnetic (AFM) material subjected to a vibrational mode in an off-resonant rotating magnetic field as shown in Fig. 1.The system Hamiltonian reads Ĥ = Ĥm + Ĥm-m + Ĥz + Ĥph + Ĥm-ph (1) with each term defined as follows.Ĥm and Ĥm-m correspond to single-and two-mode magnon Hamiltonian, respectively, at a given k-vector.To find the explicit form of these terms, we consider an AFM spin Hamiltonian, Ĥs = <i,j> Ŝi I ij Ŝj with nearest neighbor interactions, where Ŝi is the spin vector operator at lattice site i, and I ij is the bilinear interaction tensor matrix between sites i and j, which involves all different types of interactions, e.g., Heisenberg interaction, Dzyaloshinskii-Moriya interaction, anisotropy [23,24,25,26].By using the linearized Holstein-Primakoff transformation Sublattice A :

Physical system
and Fourier transformation [11], one obtains Ĥs = Ĥm + Ĥm-m with for a given momentum k-vector.Here, the bosonic annihilation operators âk and b−k represent two counter-propagating magnon modes that describe identifiable magnon modes associated with each sublattice of the AFM.We further assume that the prepared magnon system is subjected to an off-resonant rotating magnetic field with frequency ω 0 < ω a k , ω b −k .An in-plane rotating magnetic field interacts with a spin system through the Zeeman term Ĥz = i µB • S i .By using the same transformations as for Ĥs , we obtain for a given k-vector.B 0 is the strength of the magnetic field, N is the number of spins in the system, S is the total spin per site, and Γ k e −iζ k = 2 N i∈A(B) e −ik•r i with r i being the position vector of site i in the corresponding sublattice.Γ k and ζ k come about by Fourier transforming the Hamiltonian from real space to k-space.It is assumed that the values of Γ k and ζ k for sublattice A are identical to those for sublattice B.
For the magnon-phonon interaction, we consider magnetoelastic coupling, where magnetic degrees of freedom interact with elastic displacements of atoms from their equilibrium positions.Assuming single phonon process [37] with translational and rotational symmetries, the quantized phonon and magnetoelastic spin-phonon coupling terms of the Hamiltonian read [28] where i is summed over all magnetic lattice sites and δ denotes a vector pointing from lattice site i to a neighbouring site.G αβ k ′ δ is the coupling strength between the spins in spatial directions α, β ∈ {x, y, z} and the elastic displacement mode at a given k ′ .The operators ĉk ′ (ĉ † k ′ ) represent annihilation (creation) of an elastic phonon mode with frequency ω c k ′ and crystal momentum k ′ .Taking into account the magnetoelastic coupling only in the transverse direction, the interacting spin-phonon coupling term of the Hamiltonian reduces to By applying the same transformations to the spin degrees of freedom as for the other terms in the Hamiltonian above, the resulting magnetoelastic Hamiltonian takes the form where a term corresponding to a shift of each phonon mode is neglected.This type of magnon-phonon interaction naturally describes the scattering of magnons with emission or absorption of phonons based on momentum conservation [15,38,39,40,41].However, if a given pair of magnon modes âk and b−k are initially prepared in an excited state having large occupation number, the subsequent evolution will be dominated by interactions involving â † k âk and b † −k b−k , thereby forcing the phonon they couple with to have zero momentum, i.e., k ′ = 0.In the light of the dominant term, the effective magnon-phonon interaction Hamiltonian takes the form This implies that the phonon energy in the model system can be written as which physically corresponds to zone-center optical phonon modes.Note that phonon contributions from non-zero k ′ are included implicitly in the dynamics of the system as dissipation and noise terms in the quantum Langevin equations (QLEs) below.It is worth noting that the type of bosonic Hamiltonian expressed in Eq.( 8) is commonly used in optomechanical systems [1,42].By collecting the above terms for a system with magnons of a specific value of k being prepared, the relevant Hamiltonian in the rotating frame, given by applying Û with Ω = µB 0 Γ k ℏ SN 2 .Note that the parameters Ω and ζ k can be tuned by the strength B 0 and azimuthal angle ϕ 0 of the off-resonant rotating magnetic field.Thus, without loss of generality, we assume ϕ 0 + ζ k ≡ 0.
To summarize, the model Hamiltonian takes the form where we drop the subscripts k, −k, 0 for simplicity.We assume that the magnon modes in the AFM have the same frequency ω that are split by an off-resonant external rotating magnetic field with frequency ω 0 .Thus, we have the detuned frequencies ω a = ω +ω 0 and ω b = ω − ω 0 of the two magnon modes.The magnon-magnon coupling with strength g ab is known to give rise to magnon squeezing and entanglement in each magnon eigenstate in AFMs.[23,24,25,26].Physical systems that are relevant for this investigation should have zone center (optical) phonons that couple to magnons with finite crystal momentum (up to the zone boundary).Typically this happens when the frequencies are of similar magnitude.There are many reports in the literature on such materials [28,30,38].Thus, we pay a particular attention to the frequency ratio ω c /ω = 1, although, as Figs. 4 and 5 below show, our results are valid for a wider range of ω c /ω as well.

Dynamics of the system
We consider the system under dissipation caused by thermal fluctuations, uncontrolled coupling to other modes (e.g., phonons with non-zero momentum k ′ ), and Brownian motion.Such dissipative dynamics can be described by nonlinear QLEs [43].
To specify the dynamics of the system, it is handy to introduce dimensionless quadratures for the magnon modes â, b.The dissipative dynamics for magnon modes are described by where κ o is the dissipation rate, Ôin = (ô in + ôin † )/ √ 2 or Ôin = i(ô in † − ôin )/ √ 2 with ô = â, b, is the input noise operator associated with the magnon modes Ô = X, x or Ô = Ŷ , ŷ, respectively.The input noise operators are characterized by the correlation functions with the equilibrium thermal mean magnon occupation numbers In a similar manner we may define to be the dimensionless quadratures associated with the mechanical phonon mode.A Markovian description of quantum Brownian motion for a mechanical mode with large quality factor Q = ω c /γ c ≫ 1 is set by [12,43,44] where the mechanical dissipation rate γ c is associated mainly with phonons of non-zero momentum k ′ , which are not considered in the Hamiltonian Ĥ.The input noise operator ξ is given by the correlation function where is the equilibrium thermal mean phonon occupation number.Therefore, the corresponding nonlinear QLEs, in terms of dimensionless quadratures, take the form

Results and discussion
We now show that it is possible to maintain magnon-magnon entanglement in an equilibrium setting.We focus on the steady-state regime, where d dt ⟨ Ô⟩ = 0 for Ô = X, Ŷ , x, ŷ, p, q.Any operator can then be written as a steady-state expectation value plus an additional quantum fluctuation term, i.e., Ô = ⟨ Ô⟩ + δ Ô(t) [45].Imposing the steady-state conditions on Eq. ( 18) and retaining fluctuations up to linear order, we obtain a set of linearized QLEs where u(t) = (δ X, δ Ŷ , δx, δ ŷ, δ q, δ p) T and n(t) = ( Xin , Ŷ in , xin , ŷin , 0, ξ) T .Information about steady-state expectation values and coupling constants is encoded into the drift matrix: Here, ∆ ac = ω a +g ac ⟨q⟩ and ∆ bc = ω b +g bc ⟨q⟩ are the effective magnon frequency detunings induced by the magnon-phonon interaction.The effective magnon-phonon couplings are given by M X = g ac ⟨ X⟩, M Y = g ac ⟨ Ŷ ⟩, and M x = g bc ⟨x⟩, M y = g bc ⟨ŷ⟩.Assuming low dissipation of the magnon modes compared to their frequencies in typical AFM structures [46,47], i.e., κ a , κ b ≪ ∆ ac , ∆ bc , we obtain the steady-state solutions as Once the frequencies and exchange parameters are set as parameters of the model Hamiltonian in Eq. ( 11), the explicit form of the drift matrix A mainly depends on the steady-state value of ⟨q⟩, which is given as a root of the polynomial of degree five specified by the last algebraic expression in Eq. (21).As a frequency shift in magnon modes, only the real-valued roots of ⟨q⟩ are physically meaningful.For a given set of parameter values, the polynomial has only one real-valued root.
The linear form of QLEs in Eq. (19) ensures that the dynamics of quantum fluctuations is Gaussian and can be completely characterized by the corresponding covariance correlation matrix (CCM) [48].By using the fact that different components of the input noise vector n(t) are uncorrelated, the steady-state CCM denoted by V can be obtained through the Lyapunov equation [12,44] where ] is the diffusion matrix.Solving Eq. ( 22) for V , one can evaluate the bipartite magnonmagnon entanglement by computing the logarithmic negativity of two-mode Gaussian states defined as [12,27,48,49] where ν− is the minimum symplectic eigenvalue of a matrix, V 4 , being the 4 × 4 reduced CCM obtained by projecting V onto the two magnon modes.This eigenvalue can be computed as ν− = min eig |iΥ Ṽ4 |, for the symplectic matrix Υ = ⊕ 2 j=1 iσ y (with σ y being the Pauli-y matrix) and the partial transpose Ṽ4 = P 1|2 V 4 P 1|2 with P 1|2 = diag(1, −1, 1, 1).
The magnon-magnon logarithmic negativity in the steady-state limit as a function of temperature T and external magnetic field strength Ω is shown in Fig. 2 for ω c /ω = 1.
Figure 2: (Color online) Phonon-driven magnon-magnon entanglement E N as a function of temperature and the frequency Ω.Note that Ω is associated with the strength of the off-resonant rotating magnetic field, B 0 (see discussion below Eq.( 10) ).The color bar displays a qualitative illustration of the entanglement, where the white area indicates the separable zone while the colored zones shows non-zero magnon-magnon entanglement E N in the Ω − T plane.Double peaks are associated with the degeneracy of two polarized magnons in AFM systems, which is broken in the presence of non-zero magnetic field frequency ω 0 .The gap between the peaks is proportional to the frequency difference of the two magnons.The ratio ω 0 /ω with respect to the magnon frequency at the degeneracy point, i.e., ω, allows one to control the magnon-magnon entanglement region in the Ω − T parameter space.The lower the ratio ω 0 /ω, the weaker the external magnetic field required to achieve high entanglement at a given temperature (see also Fig. 3).Here, the remaining parameters of the Hamiltonian in Eq. ( 11) satisfy ω c /ω = 1 and g ac /g ab = g bc /g ab = 0.001.
Note that in the following numerical analysis, γ c is taken to be 10 2 Hz, and κ a,b ≈ 1 GHz.As can be seen, our analysis reveals magnon-magnon entanglement at room temperature under the phonon mode effect.The high-temperature result for magnon-magnon entanglement is reproduced across a wide range of Hamiltonian parameters through complementary plots in Figs. 3, 4, and 5.For any given rotation frequency ω 0 of the magnetic field, there always exists a certain range of the field amplitude Ω, where magnon-magnon entanglement can be realized at high temperature.Moreover, the entangled region including the high-temperature entanglement zone in the Ω − T parameter space can be efficiently placed around practically feasible parameter values by controlling the frequency ratio ω 0 /ω.In other words, magnon-magnon entanglement can be realized for any positive value of Ω at a given temperature T , provided the frequency ratio ω 0 /ω is properly adjusted.The results depicted in Fig. 2 were qualitatively replicated across a range of magnon-phonon coupling ratio g ac /g ab = g bc /g ab ∈ [10 −5 , 10 −3 ].This shows that, regardless of the strength of the magnon-phonon interaction in the reasonable range relative to the coupling strength between the two magnon modes, a similar outcome as presented in Fig. 2 can be achieved.Thus, the strength of magnon-phonon interaction, as treated here, does not have a detrimental effect on the phase transition of the magnon-magnon entanglement.It is instead so that, in the steady-state limit, the phonon degrees of freedom have a key role in establishing magnon-magnon entanglement, in the sense that in the absence of phonons, the two magnon modes are completely disentangled.We demonstrate this by noting that the Hamiltonian in Eq.( 11) reduces to when the phonon mode is absent.Through the similar method as shown above, we obtain the corresponding drift matrix and numerically confirm that E ab N = 0, for any choice of physical parameters of the spin Hamiltonian.We employ the drift matrix from Eq. ( 25) to solve Eq. ( 22) for the CCM, i.e., V .Subsequently, we determine the minimum eigenvalue ν − for evaluating for evaluating the logarithmic negativity.The procedure is repeated for a wide range of distinct parameter values within the drift matrix described in Eq. ( 25).In all cases, ν− consistently exceeds 1  2 at any temperature, indicating a complete absence of entanglement.This indicates that the nontrivial entanglement diagram obtained in Fig. 2 is dissipationand phonon-driven.To elaborate further on this point, in Fig. 3 we show how some cuts of the entangled region in Fig. 2 at three different temperatures move along the parameter domain of Ω upon varying the control parameter ω 0 /ω.The plots indicate that for a given value of Ω, there is a narrow interval for ω 0 /ω where magnon-magnon entanglement is non-zero and that this interval depends on temperature.The appearance of two lines in these plots at higher values of ω 0 /ω and temperatures corresponds to the two peaks in Fig. 2.These peaks result from the detuning frequencies of the two magnon modes, where ω a = ω + ω 0 and ω b = ω − ω 0 .As ω 0 approaches zero, causing the two magnon modes to become degenerate, the two peaks converge and coincide at ω 0 = 0, as indicated in Fig. 3.
Figures 2 and 3 are plotted for fixed value of the frequency ratio ω c /ω = 1 and fixed ratios g ac /g ab and g bc /g ab , and shows segments of entanglement in a parameter space of magnetic field strength and frequency ratio ω 0 /ω.In a search for a real material that can serve as an optimal host for robustly entangled magnon modes, it is instructive to also understand the dependence of entanglement on these material parameters.This is shown in Fig. 4 and 5.One can see that, in a wide range of couplings and phonon frequencies, it is possible to have non-zero entanglement at rather large (and even ambient) temperatures if frequency and amplitude of the rotating magnetic field are tuned to proper values.Fig. 4 demonstrates that for any given magnon-phonon coupling g ac = g bc , there always exists a narrow window of phonon frequencies, where the same magnon-magnon entanglement as shown in the Fig. 2 occurs.Figure 4: (Color online) Segments of magnon-magnon entanglement E N in Fig. 2 are shown by the effect of two factors, the magnon-phonon coupling ratios, g ac /g ab (= g bc /g ab ), and the phonon-magnon frequency ratio ω c /ω.The blue area marks a range in the g ac /g ab − ω c /ω parameter domain, which gives the same magnon-magnon entanglement pattern as in Fig. 2, all the other uncolored regions depict zero-entanglement (separable).Although the plots correspond to the non-zero magnon-magnon entanglement at Ω = 0.1932 THz, ω 0 /ω = 0.07, as well as temperatures T = 100 K (left) and T = 300 K (right), the same range of phonon-magnon frequency ratio and magnon-phonon coupling ratio is obtained for other points of the entanglement diagram in Fig. 2. The higher the temperature the narrower the range of relative parameters.
We further examine the effect of asymmetric magnon-phonon couplings on the magnonmagnon entanglement region.For given ω c /ω, Fig. 5 illustrates a wide range of different coupling strengths g ac ̸ = g bc , which allows the same entanglement phase diagram as in Fig. 2. Figure 5 shows that, even up to as high temperatures as T = 300 K, for different values of coupling rates varying in a wide range, there is a value of phonon frequency ω c /ω that guarantees non-zero magnon-magnon entanglement.Besides, it implies that a broad family of antiferromagnetic materials are likely to have a similar magnon-magnon entanglement diagram, as shown in Fig. 2.There are numerous suitable classes of AFM compounds that can maintain the steadystate entanglement of magnon modes within the AFM structure, even up to room temperature.Among them, NiO is a concrete example of a crystalline AFM [50,51,52,53] with high Néel temperature of T NiO N = 523K and large magnon lifetime.Notably, it has been shown that the zone center phonon frequency matches the magnon frequency at finite k frequencies ∼ 11.3 THz and ∼ 17.3 THz [28].The perovskite manganites also form a promising class of materials in this regard, with many compounds of varying chemical composition that have been synthesized (see, e.g., Refs.[31,32]).For these systems, upon chemical modulation, both the magnetic and lattice properties can be tuned.This has particular relevance for finding compositions that allow for finer adjustments of the magnon and phonon frequencies relevant for the results of the present investigation.
We end our discussion with a remark that the dissipative magnon entanglement studied here is of the mixed-state type and is achieved in the asymptotic steady-state limit.This is different from the entanglement observed in the pure energy eigenstates of ideal closed magnon systems studied previously in Ref. [23].

Conclusion
We have examined bipartite magnon-magnon entanglement in a general setting of an AFM material subjected to a dissipative elastic displacement (phonon) mode, thermal noise, and an external rotating magnetic field.We have provided evidence that it is possible to observe phonon-driven entanglement between the magnon modes of the two sublattices even at room temperature.The high-temperature magnon-magnon entanglement for a given material and magnon-phonon coupling, can be ensured by tuning the external magnetic field frequency and amplitude in a physically reasonable range of scales.The presence of a non-zero magnon-phonon coupling is necessary for maintaining entanglement between magnons in antiferromagnets in a dissipative steady-state regime, while the strength of the coupling does not significantly affect the existence of this entanglement.It would be interesting to examine whether AFMs can be used as a resource for different quantum information applications, for instance quantum state transfer [54], and to study the potential role of magnon-magnon entanglement in such processes.

Figure 1 :
Figure 1: (Color online) Schematic illustration of bipartite collective spin waves in an antiferromagnetic material subjected to a vibrational phonon mode in an off-resonant rotating magnetic field.The rotating magnetic field is illustrated in the lower left corner while the coupling between magnons in sublattice A and B to the phonons are illustrated in the lower right part of the figure.The upper figure depicts schematically the spins interacting with the lattice.

Figure 3 :
Figure 3: (Color online) Magnon-magnon entanglement E N as a function of the frequency ratio ω 0 /ω and external magnetic field strength Ω ∝ B 0 at three temperatures T = 100 K (left), T = 200 K (middle) and T = 300 K (right).The dark blue color specifies nonzero entanglement (entangled), all the other uncolored regions depict zero-entanglement (separable).We keep g ac /g ab = g bc /g ab = 0.001 and ω c /ω = 1, which are the same as the parameter values used for Fig. 2.

Figure 5 :
Figure 5: (Color online) The range of asymmetric magnon-phonon couplings, which allows for the same magnon-magnon entanglement pattern in Fig. 2, for a given phonon frequency.The contour plots correspond to the non-zero magnon-magnon entanglement segment at Ω = 0.1932 THz, ω 0 /ω = 0.07 and temperature T = 300 K (T = 100 K for inset) for ω c /ω = 0.75 (yellow), ω c /ω = 1 (pink) and ω c /ω = 1.25 (blue).The uncolored regions depict zero-entanglement (separable).Similar patterns occur for any other point in the entangling diagram of Fig. 2 and any other value of ω c /ω.The higher the temperature the narrower the range of asymmetric magnon-phonon couplings for a given ω c /ω.