Spin excitations of magnetoelectric LiNiPO$_4$ in multiple magnetic phases

TL;DR

This study investigates the spin excitations in magnetoelectric LiNiPO$_4$ across various magnetic phases using infrared spectroscopy, revealing diverse magnon modes, electromagnons, and two-magnon continua, and refining magnetic interaction parameters.

Contribution

It provides a detailed analysis of spin excitations in LiNiPO$_4$ under magnetic fields, identifying different active modes and refining exchange and anisotropy parameters with a four-sublattice model.

Findings

Identification of six spin resonances and a two-magnon continuum.

Distinction between magnetic-dipole, electric-dipole, and magnetoelectric excitations.

Refined values of exchange couplings, anisotropies, and Dzyaloshinskii-Moriya interaction.

Abstract

Spin excitations of magnetoelectric LiNiPO$_4$ are studied by infrared absorption spectroscopy in the THz spectral range as a function of magnetic field through various commensurate and incommensurate magnetically ordered phases up to 33\,T. Six spin resonances and a strong two-magnon continuum are observed in zero magnetic field. Our systematic polarization study reveals that some of the excitations are usual magnetic-dipole active magnon modes, while others are either electromagnons, electric-dipole active, or magnetoelectric, both electric- and magnetic-dipole active spin excitations. Field-induced shifts of the modes for all three orientations of the field along the orthorhombic axes allow us to refine the values of the relevant exchange couplings, single-ion anisotropies, and the Dzyaloshinskii-Moriya interaction on the level of a four-sublattice mean-field spin model. This model also reproduces the spectral shape of the two-magnon absorption continuum, found to be electric-dipole active in the experiment.

Tables

Table 1. Table 1: The parameters of the mean-field model used to describe the static magnetic properties and single- and two-magnon excitations in LiNiPO 4 : exchange couplings J i ​ j subscript 𝐽 𝑖 𝑗 J_{ij} and J k subscript 𝐽 𝑘 J_{k} , single-ion anisotropy constants Λ i subscript Λ 𝑖 \Lambda_{i} , Dzyaloshinskii-Moriya coupling D y subscript 𝐷 𝑦 D_{y} , and g 𝑔 g factor g 𝑔 g . Units are in meV except the dimensionless g 𝑔 g .

$J_y$	$J_z$	$J_{xz}$	$J_{xy}$	$J_{yz}$	${\Lambda }_x$	${\Lambda }_y$	$D_y$	$g$	Ref.
0.65	0.16	-0.17	0.16	1.24	0.14	0.74	0.41	2.2	[$\ast$]
0.67	-0.06	-0.11	0.32	1	0.41	1.42	0.32	2.2	[Toft-Petersen et al., 2011]
0.67	-0.05	-0.11	0.3	1.04	0.34	1.82			[Jensen et al., 2009a]
0.59	-0.11	-0.16	0.26	0.94	0.34	1.92			[Li et al., 2009]

Table 2. Table 2: The excitation configurations of ME mode ν 6 subscript 𝜈 6 \nu_{6} . The area of the symbol is approximately proportional to the absorption line area. The color coding is the same as in Fig.  3 .

${\nu }_6$	$E_x^{\omega }$	$E_y^{\omega }$	$E_z^{\omega }$
$H_x^{\omega }$
$H_y^{\omega }$
$H_z^{\omega }$

Equations

$$\mathcal{H}$$

$$\tan \theta \approx \frac{2D_y}{{\Lambda }_x-4(J_{xz}-J_{yz})}.$$

$${\dot{\mathbf{S}}}_i=-\frac{1}{\hslash }{\mathbf{S}}_i\times \frac{\partial \mathcal{H}}{\partial {\mathbf{S}}_i},$$

$$\delta {\dot{\mathbf{S}}}_i=-\frac{1}{\hslash }{\mathbf{S}}_i^0\times \frac{\partial {\mathcal{H}}_{\delta }}{\partial {\mathbf{S}}_i},$$

$$\frac{\partial\,\mathcal{H}_{\delta}}{\partial\,{\bf S}_{i}}=\frac{\partial\,\mathcal{H}}{\partial\,{\bf S}_{i}}\Big{|}_{\{{\bf S}^{0}\}+\{\delta{\bf S}\}}.$$

$${\dot{\mathbf{S}}}_i=-\frac{1}{\hslash }{\mathbf{S}}_i\times \frac{\partial {\mathcal{H}}_{\delta }}{\partial {\mathbf{S}}_i}+\frac{\alpha }{\hslash S_i}{\mathbf{S}}_i\times {\mathbf{S}}_i\times \frac{\partial {\mathcal{H}}_{\delta }}{\partial {\mathbf{S}}_i},$$

$$\dot{\delta {\mathbf{S}}_i}=$$

$$\gamma \hslash \left[\sum _{i=1}^M\delta {\mathbf{S}}_i(t)\right]=\hat{\chi }(\omega ){\mu }_0{\mathbf{H}}^{\omega }(t).$$

$${\mathcal{N}}_{i,j}\approx \sqrt{{ϵ}_{ii}}\left[1+\frac{{\chi }_{jj}(\omega )}{2}\right].$$

$${\alpha }_{i,j}(\tilde{\omega })=2\pi \tilde{\omega }\sqrt{{ϵ}_{ii}}\mathrm{Im}{\chi }_{jj}(\tilde{\omega }),$$

Full text

1 Introduction

Potential of magnetoelectric (ME) materials in applications relies on the entanglement of magnetic moments and electric polarization Kimura et al. (2003); Fiebig (2005); Spaldin and Fiebig (2005); Eerenstein et al. (2006); Cheong and Mostovoy (2007); Fiebig and Spaldin (2009); Dong et al. (2015); Fiebig et al. (2016). Such an entanglement leads not only to the static ME effect but also to the optical ME effect. One manifestation of the optical ME effect is the non-reciprocal directional dichroism, a difference in the absorption with respect to the reversal of light propagation direction Rikken and Raupach (1997); Rikken et al. (2002); Barron (2004); Arima (2008). The spectrum of non-reciprocal directional dichroism and the linear static ME susceptibility are related via a ME sum rule Szaller et al. (2014). According to this sum rule the contribution of simultaneously magnetic- and electric-dipole active spin excitations to the linear ME susceptibility grows as $\omega^{-2}$ with $\omega\rightarrow 0$. Indeed, strong non-reciprocal directional dichroism has been observed at low frequencies, typically in the GHz-THz range, at spin excitations in several ME materials Kézsmárki et al. (2011); Bordács et al. (2012); Takahashi et al. (2012); Kézsmárki et al. (2014); Kézsmárki et al. (2015); Kuzmenko et al. (2015); Bordács et al. (2015); Okamura et al. (2015); Nii et al. (2017); Yu et al. (2018); Kocsis et al. (2018); Viirok et al. (2019); Okamura et al. (2019). Beside the interest in the non-reciprocal effect, the knowledge of the spin excitation spectrum and selection rules, i.e. whether the excitations are ordinary magnetic-dipole active magnons, electromagnons (electric-dipole active magnons Pimenov et al. (2006)), or ME spin excitations (simultaneously magnetic- and electric-dipole active spin excitations), is crucial in understanding the origin of static ME effect.

It is well established that the static ME effect is present in several olivine-type Li$M$PO4 ($M=\mathrm{Mn}$, Fe, Co, Ni) compounds Mercier and Gareyte (1967); Mercier et al. (1967); Mercier and Bauer (1968); Mercier et al. (1969); Rivera (1994); Kornev et al. (2000); Toft-Petersen et al. (2015); Khrustalyov et al. (2017). LiNiPO4 is particularly interesting due to many magnetic-field-induced phases, some with incommensurate magnetic order, which is unique in the olivine lithium-ortho-phosphate family Toft-Petersen et al. (2011). However, little is known about the spectrum of spin excitations and their selection rules.

THz absorption spectroscopy offers an excellent tool to investigate spin excitation spectra over a broad magnetic field range. As compared to the inelastic neutron scattering (INS), only spin excitations with zero linear momentum are probed, but with a better energy resolution. In addition to excitation frequencies, THz spectroscopy can determined whether the spin excitations are magnetic-dipole active magnons, electromagnons, or ME spin excitations. This information is essential for developing a spin model that would describe the ground and the low lying excited states of the material.

We studied the spin excitation spectra of LiNiPO4 in magnetic field using THz absorption spectroscopy. In the previous INS works two magnon branches were observed below 8 meV Jensen et al. (2009a); Li et al. (2009); Toft-Petersen et al. (2011). Here we broaden the spectral range up to 24 meV, which allows us to observe additional spin excitations and to identify the polarization selection rules for the spin excitations. Using a mean-field model we describe the field dependence of the magnetization and the magnon energies in commensurate phases, from where we refine the values of exchange couplings, single-ion anisotropies, and the Dzyaloshinskii-Moriya interaction. Beside magnons described by the mean-field model few other spin excitations, including two-magnon excitations, are observed.

LiNiPO4 has orthorhombic symmetry with space group $Pnma$. The magnetic Ni2+ ion with spin $S=1$ is inside a distorted O6 octahedron. There are four Ni2+ ions in the structural unit cell forming buckled planes perpendicular to the crystal $x$ axis, as shown in Fig. 1. The nearest-neighbor spins in the $yz$ plane are coupled by strong AF exchange interaction which results in a commensurate AF order below $T_{\textrm{N}}$=20.8 K Kharchenko et al. (2003); Vaknin et al. (2004). The ordered magnetic moments are almost parallel to the crystallographic $z$ axis with slight canting towards the $x$ direction Khrustalyov et al. (2004). On heating above $T_{\textrm{N}}$ the material undergoes a first-order phase transition into a long-range incommensurate magnetic structure. Further heating results in a second-order phase transition into the paramagnetic state at $T_{\mathrm{IC}}=21.7$ K, while short-range magnetic correlations persist up to 40 K Vaknin et al. (2004). Owing to the competing magnetic interactions LiNiPO4 has a very rich $H$–$T$ phase diagram with transitions appearing as multiple steps in the field dependence of the magnetization Khrustalyov et al. (2016); Toft-Petersen et al. (2017). The delicate balance of the nearest-neighbor and the frustrated next-nearest-neighbor exchange interactions puts the material on the verge of commensurate and incommensurate structures, which alternate in increasing magnetic field applied along the $z$ axis as shown in Fig. 2(a) Khrustalyov et al. (2004); Li et al. (2009); Jensen et al. (2009a).

2 Experimental details

LiNiPO4 single crystals were grown by the floating zone method, similarly as described in Ref. [Baker et al., 2011]. Three samples each with a large face normal to one of the principal axes were cut from the same ingot. For optical measurements the slabs with thicknesses from 0.87 mm to 1.09 mm had approximately two degree wedge to suppress interference caused by internal reflections. Samples were mounted on metal discs where the hole depending on the sample size limited the THz beam cross-section to 8 – 16 mm2.

THz measurements up to 17 T were performed in Tallinn with a Martin-Puplett interferometer and a 0.3 K silicon bolometer. High-field spectra from 17 T up to 33 T were measured in Nijmegen High Field Magnet Laboratory using a Bruker IFS 113v spectrometer and a 1.6 K silicon bolometer. The experiments above 17 T were done in the Faraday configuration ($\mathbf{k}\!\parallel\!\mathbf{H}$), while below 17 T both the Faraday and the Voigt ($\mathbf{k}\!\perp\!\mathbf{H}$) configuration experiments were performed. All spectra were measured with an apodized spectral resolution of 0.5  cm*-1*. A linear polarizer was mounted in front of the sample to control the polarization state of the incoming light.

Absorption was determined by using a reference spectrum. The reference spectrum was obtained on the sample in zero magnetic field in the paramagnetic state at $T=30$ K or by measuring a reference hole with the area equal to the sample hole area. In the former case the relative absorption is calculated as $\alpha(H,T)-\alpha(0\text{T},30\text{K})=-d^{-1}\ln\left[I(H,\textit{T})/I(0\text{T},30\text{K})\right]$ where $d$ is the sample thickness and $I$ is the measured intensity. In the latter case the absolute absorption is calculated as $\alpha=-d^{-1}\ln(I/I_{r})$ where $I_{r}$ is the intensity through the reference hole.

Magnetization up to 32 T was measured in Nijmegen High Field Magnet Laboratory on a Bitter magnet with a vibrating-sample magnetometer (VSM) and additional low field measurements were done using a 14 T PPMS with VSM option (Quantum Design).

3 Mean-field model and magnons

The terms included in the spin Hamiltonian, exchange interactions, single-ion anisotropy terms, and the Zeeman energy, correspond to those also considered in earlier works on LiNiPO4 Li et al. (2009); Jensen et al. (2009a); Toft-Petersen et al. (2011). The model contains four spin variables as classical vectors $\{\mathbf{S}_{1},\mathbf{S}_{2},\mathbf{S}_{3},\mathbf{S}_{4}\}$ in accordance with the four crystallographically non-equivalent positions of the spin $S\,=\,1$ Ni2+ ions in LiNiPO4. The four spins of the magnetic unit cell are connected by five different exchange couplings as presented in Fig. 1. Two of these couplings, $J_{y}$ and $J_{z}$, connect spins at the same crystallographic sites producing, irrespective of the spin state, a constant energy shift in the $\Gamma$ point within the four-sublattice model. Although these terms are omitted in the analysis of single-magnon excitations, they become relevant in the analysis of two-magnon excitations as discussed in Sec. 5.2. The spin Hamiltonian of the magnetic unit cell in the four-sublattice model reads

[TABLE]

Due to the strongly distorted ligand cage of the magnetic ion, the orthorhombic anisotropy of the crystal is taken into account by two single-ion hard-axis anisotropies, $\Lambda_{x},\Lambda_{y}>0$. The parameters in the Zeeman term are the $g$ factor $g$, the Bohr magneton $\mu_{\mathrm{B}}$, and the vacuum permeability $\mu_{0}$. Parameters $J_{xz}$, $J_{xy}$, and $J_{yz}$ are the isotropic Heisenberg exchange couplings as shown in Fig. 1, while $D_{y}$ is the Dzyaloshinskii-Moriya interaction.

According to the neutron diffraction studies Santoro et al. (1966); Vaknin et al. (1999) the ground-state spin configuration of LiNiPO4 in zero magnetic field is a predominantly collinear AF order, where $\mathbf{S}_{1}$ and $\mathbf{S}_{2}$ point in $+z$, while $\mathbf{S}_{3}$ and $\mathbf{S}_{4}$ in $-z$ direction, shown in Fig. 1. Thus, the dominant exchange interaction is the AF $J_{yz}>0$ coupling, while $z$ is an easy axis as $\Lambda_{x}$, $\Lambda_{y}>0$. On top of the collinear order a small alternating canting of spins with net spin along $x$ is superimposed Jensen et al. (2009b). Canting is induced by the Dzyaloshinskii-Moriya coupling $D_{y}$ and breaks the equivalence of $\mathbf{S}_{1}$ and $\mathbf{S}_{2}$ as well as $\mathbf{S}_{3}$ and $\mathbf{S}_{4}$. The canting angle $\theta$ measured from the $z$ axis is approximately

[TABLE]

At each magnetic field, the ground-state spin configuration is obtained by minimizing the energy corresponding to Eq. (1).

The resonance frequencies and amplitudes of modes are calculated using the Landau-Lifshitz equation White (2007)

[TABLE]

where $\dot{{\bf S}}_{i}\equiv{\mathrm{d}}{\bf S}_{i}/{\mathrm{d}}t$.

We solve Eq. (3) for small spin deviations $\{\delta{\bf S}\}\equiv\{\delta{\bf S}_{1},\ldots\delta{\bf S}_{N}\}$ from the equilibrium $\{{\bf S}^{0}\}\equiv\{{\bf S}_{1}^{0},\ldots{\bf S}_{N}^{0}\}$, where $\{{\bf S}\}=\{{\bf S}^{0}\}+\{\delta{\bf S}\}$, with $N$ spins in the magnetic unit cell. It follows from Eq. (3) that $\delta{\bf S}_{i}\perp{\bf S}_{i}$, leaving the spin length constant in the first order of $\delta{\bf S}_{i}$. Inserting $\{{\bf S}\}$ into Landau-Lifshitz equation Eq. (3) and keeping only terms linear in $\delta{\bf S}_{i}$ (terms zero-order in $\delta{\bf S}_{i}$ add up to zero) we get

[TABLE]

where the effective field is

[TABLE]

We solve Eq. (4) by assuming harmonic time dependence $\delta{\bf S}_{i}(t)=\delta{\bf S}_{i}\exp(i\omega t)$. The number of modes is equal to the number of spins in the unit cell.

To calculate the absorption of electromagnetic waves by the magnons we introduce damping. The Landau-Lifshitz-Gilbert equation Gilbert (2004) for the $i$-th spin is

[TABLE]

where $\alpha$ is a positive dimensionless damping parameter and small, $\alpha\ll 1$. Using $\mathbf{A}\times\mathbf{B}\times\mathbf{C}=\mathbf{B}(\mathbf{A}\cdot\mathbf{C})-\mathbf{C}(\mathbf{A}\cdot\mathbf{B})$, and adding a weak harmonic alternating magnetic field, $\mathbf{H}^{\omega}(t)=\mathbf{H}^{\omega}\exp(i\omega t)$, to the effective field yields the following form of the equation of motion up to terms linear in $\delta{\bf S}_{i}$ and $\mathbf{H}^{\omega}$:

[TABLE]

The absorption of electromagnetic waves by the spin system related to magnetic dipole excitations is calculated from the Eq. (3) by inserting $\delta{\bf S}_{i}(t)=\delta{\bf S}_{i}\exp(i\omega t)$ and $\mathbf{H}^{\omega}(t)=\mathbf{H}^{\omega}\exp(i\omega t)$. The frequency-dependent magnetic susceptibility tensor $\hat{\chi}(\omega)$ is obtained by summing up all the magnetic moments in the unit cell, $\mathbf{M}=\gamma\hbar\sum_{i=1}^{N}{\bf S}_{i}$ in Eq. (3), and making a transformation into form

[TABLE]

The absorption coefficient is $\alpha_{i,j}=2\omega c_{0}^{-1}\mathrm{Im}\,\mathcal{N}_{i,j}$, where the complex index of refraction is $\mathcal{N}_{i,j}=\sqrt{\epsilon_{ii}\mu_{jj}}$ assuming small polarization rotation and negligible linear magnetoelectric susceptibilities $\chi^{em}_{ij}$, $\chi^{me}_{ji}$. The magnetic permeability is $\mu_{jj}(\omega)=1+\chi_{jj}(\omega)$ and the background dielectric permittivity is $\epsilon_{ii}$. The polarization of incident radiation is defined as $\{E_{i}^{\omega},H_{j}^{\omega}\}$ where $i$ and $j$ are $x$, $y$, or $z$. If $\chi_{jj}(\omega)\ll 1$,

[TABLE]

Thus, for real $\epsilon_{ii}$ the absorption is

[TABLE]

where units of wavenumber are used, $[\tilde{\omega}_{n}]=$ cm*-1*.

The values of magnetic interactions and anisotropies obtained in this work, see Table 1, reproduce the magnetic field dependence of the magnetization, canting angle $\theta$, and frequencies of four spin excitations in the commensurate magnetic phase of LiNiPO4.

4 Experimental results

The LiNiPO4 samples were characterized by measuring the magnetization along the $x$, $y$, and $z$ directions, shown in Fig. 2. The magnetization increases continuously for $\mathbf{H}\!\parallel\!\mathbf{{x}}$ and $\mathbf{H}\!\parallel\!\mathbf{{y}}$, while for $\mathbf{H}\!\parallel\!\mathbf{{z}}$ there is step at 12 T, 19 T, and 21.5 T. These steps correspond to magnetic field induced changes in the ground state spin structure. The phases I and IV are commensurate, while II, III, and V are incommensurate Toft-Petersen et al. (2011, 2017). The boundary between phases II and III at 16 T, where the periodicity of the incommensurate spin structure changes Toft-Petersen et al. (2011), is hardly visible in the magnetization data Khrustalyov et al. (2004); Toft-Petersen et al. (2017). The size of the magnetic unit cells is the same in phases I and IV Toft-Petersen et al. (2017), i.e., four spins as shown in Fig. 1.

\onecolumngrid@push

\onecolumngrid@pop

The zero-field THz absorption spectra measured at 3.5 K are shown in Fig. 3. Three absorption lines are identified as magnetic-dipole active magnons: $\nu_{1}=16.1$ cm*-1*, $\nu_{2}=36.2$ cm*-1*, and $\nu_{3}=48.4$ cm*-1*. The excitation $\nu_{5}=56.4$ cm*-1* is an $E_{y}^{\omega}$-active electromagnon. The excitations $\nu_{4}=54.8$ cm*-1* and $\nu_{6}=66.4$ cm*-1* are ME spin excitations; $\nu_{4}$ is $\{E^{\omega}_{x},H^{\omega}_{z}\}$-active, while $\nu_{6}$ is present in five different combinations of oscillating electric and magnetic fields with the strongest intensity in $E^{\omega}_{z}$ polarization, see Table 2. There is an $E^{\omega}_{x}$-active broad absorption band $\nu_{7}$.

All seven modes $\nu_{1},\ldots,\nu_{7}$ are absent above $T_{\textrm{N}}$. Since no sign of structural changes has been found in the neutron diffraction Jensen et al. (2009b) and in the spectra of Raman-active phonons at $T_{\textrm{N}}$ Fomin et al. (2002), the lattice vibrations can be excluded and all new modes are assigned to spin excitations of LiNiPO4.

The magnetic field dependence of resonance frequencies and absorption line areas is presented in Fig. 4 as obtained from the fits of the absorption peaks with the Gaussian line shapes. When the magnetic field is applied in $\mathbf{H}\!\parallel\!\mathbf{{x}}$ or $\mathbf{H}\!\parallel\!\mathbf{{y}}$ directions, Fig. 4(a) or (b), we found a continuous evolution of modes up to the highest field of 33 T. On contrary, for $\mathbf{H}\!\parallel\!\mathbf{{z}}$ we observed discontinuities in the spin excitation frequencies, approximately at 12 T, 19 T and 21.5 T. These fields correspond to the field values where the steps are seen in the magnetization in Fig. 2. The boundary between II and III at 16 T is not visible in the THz spectra. Apparently the spin excitation spectra are rather insensitive to the change of the magnetic unit cell size within the incommensurate phase.

The mean-field model (Sec. 3) predicts four magnon modes for a four sub-lattice system and they are assigned to $\nu_{1}$, $\nu_{2}$, $\nu_{3}$ and $\nu_{5}$. The magnetic field dependence and the selection rules of the magnetic-dipole active magnons $\nu_{1}$, $\nu_{2}$ and $\nu_{3}$ are reproduced well by the mean-field model, Fig. 4. However, only the energy of the magnon $\nu_{5}$ is reproduced by the model and not the intensity as this excitation is found to be an electromagnon in the experiment.

The resonances $\nu_{4}$, $\nu_{6}$, and the band $\nu_{7}$ cannot be described within the four-sublattice mean-field model. The weak $\nu_{4}$ mode is a ME spin excitation, $\{E^{\omega}_{x},H^{\omega}_{z}\}$ active, which might be related to a spin-stretching excitation allowed for $S>1/2$ Penc et al. (2012). The $\nu_{6}$ mode is a ME two-magnon excitation and $\nu_{7}$ an $E^{\omega}_{x}$-active two-magnon excitation band, as will be discussed below. The excitations $\xi_{8},\xi_{9},\xi_{10}$, and $\theta_{12},\theta_{13}$ are only present in the incommensurate phases II, III, and in phase V with more than four spins per magnetic unit-cell and thus cannot be explained by the present four-sublattice model. The field dependence and the selection rules of $\eta_{11}$, the only mode found experimentally in the four-sublattice commensurate phase IV, are described by the model, Fig. 4(c) and 4(f).

There are two resonances in the vicinity of the $\nu_{1}$ mode as indicated by blue symbols in Fig. 4(b) and (e). These two modes have a well-defined selection rule, $\mathit{H}^{\omega}_{x}$. Because they are at low frequency but not described by the mean field model we assign them to impurity modes.

The exchange parameters obtained by fitting the mean-field model to THz spectra are presented in Table 1. Our model also reproduces the magnetization for commensurate phases I and IV, Fig. 2. The canting angle of spins given by the parameters of the current work, Table 1 and Eq. (2), is $\pm\theta=8.1^{\circ}$ in zero field, in good agreement with the value determined by elastic neutron scattering, $(7.8\pm 2.6)^{\circ}$, as reported in Ref. [Toft-Petersen et al., 2011; Jensen et al., 2009b].

5 Discussion

5.1 One-magnon excitations

The four sub-lattice mean-field model describes four magnons $\nu_{1}$, $\nu_{2}$, $\nu_{3}$, and $\nu_{5}$, among which $\nu_{1}$ and $\nu_{2}$ can be identified as $\Gamma$-point magnon modes observed in the INS spectra Jensen et al. (2009a), whereas the $\nu_{5}$ resonance has also been detected by the Raman spectroscopy Fomin et al. (2002).

The zero-field frequencies of $\nu_{1}$ and $\nu_{2}$ are related to the single-ion anisotropies $\Lambda_{x}$ and $\Lambda_{y}$, respectively. Furthermore, the selection rules for the $\nu_{1}$ and $\nu_{2}$ suggest that they are anisotropy-gapped magnons, since in both cases the magnetic dipole moment oscillates perpendicular to the corresponding anisotropy axis, along $y$ for the $\nu_{1}$ mode and along $x$ for $\nu_{2}$ in zero field. Moreover, the mean-field model reproduces the rotation of the magnetic dipole moment of $\nu_{1}$ ($\nu_{2}$) towards the $z$ axis in increasing magnetic field $\mathbf{H}\!\parallel\!\mathbf{{y}}$ ($\mathbf{H}\!\parallel\!\mathbf{{x}}$). The reappearance of $\nu_{1}$ in phase IV, marked as $\eta_{11}$, is also predicted by the model.

The frequencies of $\nu_{3}$ and $\nu_{5}$ depend strongly on the weak $J_{xy}$ and $J_{xz}$ exchange interactions connecting the two AF systems, $\{\mathbf{S}_{1},\mathbf{S}_{4}\}$ and $\{\mathbf{S}_{2},\mathbf{S}_{3}\}$. While the FM $J_{xz}$ only shifts the average frequency of $\nu_{3}$ and $\nu_{5}$, the AF $J_{xy}$ affects the difference frequency. The zero-field selection rules of these excitations – magnetic dipole moment along $z$ for $\nu_{3}$ and the absence of magnetic-dipole activity of $\nu_{5}$ – are reproduced by the model.

Our model does not describe the incommensurate phases II, III and V. However, it reproduces the frequency of the lowest $\eta_{11}$ mode in the commensurate phase IV, Fig. 4(c).

5.2 Two-magnon excitations

Two-magnon excitations appear in the absorption spectra when one absorbed photon creates two magnons with the total $\mathbf{k}$ vector equal to zero Richards (1967). The two-magnon absorption is the strongest where the density of magnon states is the highest, usually at the Brillouin zone boundary. Since the product of the two spin operators has the same time-reversal parity as the electric dipole moment, the two-magnon excitation by the electric field is allowed; this mechanism usually dominates over the magnetic-field induced two-magnon excitation Richards (1967).

The broad absorption band between 60 and 115 cm*-1*, shown in Fig. 5(a), appears below $T_{\textrm{N}}$ and is $\mathit{E}^{\omega}_{x}$ active. A similar excitation band observed by Raman scattering was attributed to spin excitations Fomin et al. (2002). In another olivine-type crystal LiMnPO4, a broad band in the Raman spectrum was assigned to a two-magnon excitation with the line shape reproduced using the magnon density of states (DOS) Filho et al. (2015).

We calculated the magnon DOS numerically on a finite-size sample of $4\times 4\times 4$ unit cells with 256 spins using the model represented by Eq. (1) but extended by the $J_{y}$ and $J_{z}$ couplings shown in Fig. 1. The two-magnon DOS was obtained by doubling the energy scale of the single-magnon DOS and is shown in Fig. 5(b). Since the observed broad absorption band emerges in the energy range of the high magnon DOS we assign this absorption band to a two-magnon excitation.

There is another dominantly electric-dipole active spin resonance $\nu_{6}$ at 66.5 cm*-1* not reproduced by our mean-field model. Since the frequency of the electric-dipole active $\nu_{6}$ mode is at the maximum of the two-magnon DOS and in a magnetic field, $\mathbf{H}\parallel\mathbf{z}$, splits into a lower and an upper resonances with effective $g$ factors $g_{-}=4.24\pm 0.07$ and $g_{+}=4.00\pm 0.04$, i.e two times larger than that expected for one spin-flip excitation, we interpret the $\nu_{6}$ resonance as a two-magnon excitation. The singular behavior in the DOS coinciding with the $\nu_{6}$ resonance peak corresponds to the flat magnon dispersion along the R–T line in the Brillouin zone. This mode is weakly magnetic-dipole active as well and is therefore a ME resonance.

The magnetic field dependence of the two-magnon excitation $\nu_{6}$ was modeled by calculating the field dependence of the magnon DOS. The result is shown in Fig. 4. The splitting of the resonance in magnetic field is observed only for $\mathbf{H}\!\parallel\!\mathbf{{z}}$ and is reproduced by the model calculation. In the calculation $J_{y}$ was set to 0.65 meV to reproduce the instability of phase I at 12 T, while $J_{z}$ =0.16 meV was used to reproduce the zero-field frequency of $\nu_{6}$. The magnitude of $J_{y}$ and $J_{z}$ is similar to the ones in INS studies Toft-Petersen et al. (2011); Jensen et al. (2009a); Li et al. (2009) but $J_{z}$ has the opposite sign. In high-symmetry cases the electric-dipole selection rules of two-magnon excitations can be reproduced by group theoretical analysis Loudon (1968), but the low magnetic symmetry of LiNiPO4 hinders such an analysis. Nevertheless, it is expected that the $\Delta S=0$ two-magnon continuum $\nu_{7}$ has different selection rules than the $\Delta S=2$ two-magnon excitation $\nu_{6}$ due to their different symmetry.

6 Conclusions

We measured the magnetic field dependence of THz absorption spectra in various magnetically ordered phases of LiNiPO4. We have revealed a variety of spin resonance modes: three magnons, an electromagnon, an electric-dipole active two-magnon excitation band, and a magnetoelectric two-magnon excitation. The abrupt changes in the magnon absorption spectra coincide with the magnetic phase boundaries in LiNiPO4. The magnetic dipole selection rules for magnon absorption and the magnetic field dependence of magnon frequencies in the commensurate magnetic phases are described with a mean-field spin model. With this model the additional information obtained from the magnetic field dependence of mode frequencies allowed us to refine the values of exchange couplings, single-ion anisotropies, and Dzyaloshinskii-Moriya interaction. The significant differences found in magnetic interaction parameters compared to former studies are the opposite sign of $J_{z}$ exchange coupling, the smaller values of the $J_{xy}$ exchange coupling and the $\Lambda_{x}$ and $\Lambda_{y}$ anisotropies. The mean-field model did not explain the observed magneto-electric excitation $\nu_{4}$ and the spin excitations in the incommensurate phases. In the future, more about the magneto-electric nature of LiNiPO4 spin excitations can be learned from non-reciprocal directional dichroism studies as in LiCoPO4 Kocsis et al. (2018).

7 Acknowledgments

The authors are indebted to László Mihály and Karlo Penc for valuable discussions. This project was supported by institutional research funding IUT23-3 of the Estonian Ministry of Education and Research, by the European Regional Development Fund project TK134, by the bilateral program of the Estonian and Hungarian Academies of Sciences under the Contract No. NKM-47/2018, by the Hungarian NKFIH Grant No. ANN 122879, by the BME Nanotechnology and Materials Science FIKP grant of EMMI (BME FIKP-NAT), and by the Deutsche Forschungsgemeinschaft (DFG) via the Transregional Research Collaboration TRR 80: From Electronic Correlations to Functionality (Augsburg-Munich-Stuttgart). D. Sz. acknowledges the FWF Austrian Science Fund I 2816-N27 and V. K. was supported by the RIKEN Incentive Research Project. High magnetic field experiments in HFML were funded by EuroMagNET under the EU Contract No. 228043.