Electron spin resonance for the detection of long-range spin nematic order
Shunsuke C. Furuya Tsutomu Momoi
cond-mat.str-el




Condensed Matter Theory Laboratory, RIKEN, Wako, Saitama 351-0198, Japan

Condensed Matter Theory Laboratory, RIKEN, Wako, Saitama 351-0198, Japan
RIKEN Center for Emergent Matter Science (CEMS), Wako, Saitama, 351-0198, Japan


 Spin nematic phase is a quantum magnetic phase characterized by a quadrupolar order parameter.
 Since the quadrupole operators are directly coupled to neither the magnetic field nor the neutron,
 currently, it is an important issue to develop a method for detecting the long-range spin nematic order.
 In this paper we propose that electron spin resonance (ESR) measurements enable us to detect the long-range spin nematic order.
 We show that the frequency of the paramagnetic resonance peak in the ESR spectrum
 is shifted by the ferroquadrupolar order parameter together with other quantities.
 The ferroquadrupolar order parameter is extractable from the angular dependence of the frequency shift.
In contrast, the antiferroquadrupolar order parameter is usually invisible in the frequency shift.
 Instead, the long-range antiferroquadrupolar order yields a characteristic resonance peak in the ESR spectrum,
 which we call a magnon-pair resonance peak.
This resonance corresponds to the excitation of the bound magnon pair at the wave vector
k= 0.
 Reflecting the condensation of bound magnon pairs, the field dependence of the magnon-pair resonance frequency
 shows a singular upturn at the saturation field.
 Moreover, the intensity of the magnon-pair resonance peak shows a characteristic angular dependence
 and it vanishes when the magnetic field is parallel to one of the axes that diagonalize the weak anisotropic interactions.
 We confirm these general properties of the magnon-pair resonance peak in the spin nematic phase
 by studying an S=1 bilinear-biquadratic model on the  square lattice in the linear flavor-wave approximation.
 In addition, we argue applications to the S=1/2 frustrated ferromagnets and
 also the S=1/2 orthogonal dimer spin system SrCu_2(BO_3)_2,
 both of which are candidate materials of spin nematics.
 Our theory for the antiferroquadrupolar ordered phase
 is consistent with many features of the magnon-pair resonance peak experimentally observed in the
 low-magnetization regime of SrCu_2(BO_3)_2.

Electron spin resonance for the detection of long-range spin nematic order
    Tsutomu Momoi
    December 30, 2023
==========================================================================





§ INTRODUCTION


Spin nematic phase is a hidden ordered phase of quantum magnets where
the spin rotation symmetry is spontaneously broken
but, different from the ferromagnetic and the antiferromagnetic phases,
the time reversal symmetry
is kept intact.
The spin nematic phase is characterized by a quadrupolar order parameter made of a symmetric pair
of electron spins <cit.>.
The emergence of the spin nematic phase requires the absence of the spontaneous dipolar orders,
implying interplay and frustration of spin-spin interactions behind the order.
Until today, much
effort has been made to explore the spin nematic phase in S=1/2 frustrated ferromagnets both theoretically <cit.>
and experimentally <cit.>,
in the S=1/2 orthogonal dimer spin system <cit.>,
and also in S≥ 1 spin systems with biquadratic interactions <cit.>.
Recently, the research field of the spin nematic order expands to the field of an iron pnictide superconductor,
FeSe <cit.>.


In the current situation surrounding researches of the spin nematic phase,
one of the most important problems is to develop a method for detecting the spin nematic order
in an experimentally feasible way.
The difficulty in the problem is that the quadrupole operators are directly coupled to neither the magnetic field nor the neutron.
Several theoretical proposals were recently made
by studying theories of the nuclear magnetic resonance
(NMR) <cit.>, inelastic neutron scattering <cit.>,
inelastic light scattering <cit.>,
resonant inelastic x-ray scattering <cit.>, and electron spin resonance (ESR) <cit.>.
As pointed out in Ref. <cit.>,
various quantities related to ESR are naturally coupled to quadrupole operators
and thus ESR is a promising way for detecting the hidden spin nematic order.
In ESR experiments, bound triplon-pair modes were observed in the spin gapped phase of
the orthogonal dimer spin system SrCu_2(BO_3)_2 <cit.> and a bound magnon-pair mode was also
in the fully polarized phase of Sr_2CoGe_2O_7 <cit.>.
Bound magnon pairs can give the instability to the spin nematic ordering when they close the energy
gap <cit.>.
The behavior of these bound pair modes in ESR spectrum is however not yet understood
inside the spin nematic phase from both the theoretical and experimental sides.
In view of the current situations,
an ESR theory for identification of the quadrupolar order in the spin nematic phase is called for.

In this paper, we theoretically study ESR in the spin nematic phase to elucidate how to identify
the quadrupolar order from the ESR spectrum.
We propose two methods of detecting spin nematic orders.
In the first part of the paper, we demonstrate for quite a generic model
that ESR measurements enable us to extract the ferroquadrupolar (FQ) order parameter
(i.e. the spin nematic order parameter
developed at the wave vector k= 0) from the frequency shift of the electron paramagnetic resonance
(EPR) peak in the ESR spectrum,
as one of the authors briefly mentioned in Ref. <cit.>.
Incidentally the antiferroquadrupolar (AFQ) order parameter is not detectable in this method.
In the latter part of the paper,
we propose a complementary ESR measurement suitable for identification of the AFQ order,
taking an example of an S=1 spin model.
In this method, we focus on resonance of bound magnon-pair excitations, which characteristically appear
in quadrupolar order phases signaling the condensation of bound magnon pairs.
The corresponding resonance peak, which we call the magnon-pair resonance peak, appears at a finite frequency
in the ESR spectrum in the AFQ phase.
This resonance peak can be clearly distinguished from others since the peak intensity has a characteristic
angular dependence and the peak frequency also has a distinctive field dependence.
A typical field dependence of the magnon-pair resonance frequency
is shown as the solid curve in Fig. <ref> for the AFQ phase in an S=1 bilinear-biquadratic model.
Having the two methods of ESR for the FQ and the AFQ orders,
we can characterize the FQ and the AFQ phases clearly.
We can apply our theory to S=1/2 candidate materials of spin nematics, such as S=1/2 frustrated ferromagnets and S=1/2
orthogonal dimer spin compound SrCu_2(BO_3)_2.
Our theory for the AFQ ordered phase is consistent with many features of the experimentally observed
magnon-pair resonance peak <cit.> in SrCu_2(BO_3)_2, which suggests the existence of a spin nematic order.



This paper is organized as follows.
In Sec. <ref>, we review briefly an important identity [Eq. (<ref>)] that we rely on in this paper.
Using the identity, we show in Sec. <ref> that the EPR frequency is shifted by the FQ order parameter
depending on the direction of the sample.
We use a little trick in order to extract the FQ order parameter from the frequency shift.
Since the EPR frequency is usually insensitive to the AFQ order parameter, in Sec. <ref>,
we focus on a low-energy bosonic excitation that qualifies as evidence of the presence of the long-range AFQ order.
In Sec. <ref>, we briefly explain the concept of how to detect the boson in ESR experiments in general.
To flesh out the general discussion, we take an example of an S=1 bilinear-biquadratic model on the square lattice (Secs. <ref> and <ref>) and study its ESR with the aid of the linear flavor-wave theory.
The linear flavor-wave theory in the fully polarized phase (Sec. <ref>) and in the AFQ phase
(Sec. <ref>) shows that the boson corresponding to the
bound magnon-pair excitation yields a resonance peak in the ESR spectrum in addition to the EPR one.
The magnon-pair resonance is closely investigated in the fully polarized phase in Sec. <ref>.
Section <ref> is devoted to the investigation of additional resonances in the AFQ phase,
where we find the magnon-pair resonance and another unpaired magnon resonance.
Some related discussions are addressed in Sec. <ref>. Section <ref>
contains
applications to the S=1/2 frustrated ferromagnets and the S=1/2 orthogonal dimer spin system.
Finally, we summarize the paper in Sec. <ref>.









§ FRAMEWORK


Here we describe a generic model that we deal with in the paper and an important identity that we rely on throughout our discussions.
We consider quantum spin systems described by the following Hamiltonian:

    ℋ = ℋ_ SU(2) -H S^z + ℋ',

where ℋ_ SU(2) denotes interactions that respect the SU(2) rotation symmetry of the spin,
the second term is the Zeeman energy,
S^z=∑_i S_i^z is the z component of the total spin S = ∑_i  S_i,
and ℋ' is an anisotropic interaction that breaks weakly the SU(2) symmetry of ℋ_ SU(2).
For simplicity, we take ħ=k_B=a_0=1 hereafter, where a_0 is the lattice constant.


We regard ℋ' as a perturbation to the
Hamiltonian

    ℋ_0 = ℋ_ SU(2) - HS^z.

We emphasize that, in this paper, we consider the cases where the spin nematic phase emerges in the unperturbed system with the Hamiltonian (<ref>) and the perturbation ℋ' affects neither the existence of the
quadrupolar order nor the magnitude of the quadrupolar order parameter.
On the other hand, when we consider ESR of quantum spin systems, however small ℋ' may be,
we must take the anisotropic interaction ℋ' into account.

The ESR absorption spectrum is given by the imaginary part of the retarded Green's function of the total spin.
The ESR spectrum I(ω) in the Faraday configuration is related to 𝒢^R_S^+S^-(ω),
where S^± = S^x ± i S^y are the ladder operators of the total spin,
𝒢^R_O_1O_2(ω) is the retarded Green's function
𝒢^R_O_1O_2(ω)=-i∫_0^∞ dt   e^iω t⟨[O_1(t), O_2(0)]|$⟩.

In fact, when the applied electromagnetic wave is circularly polarized,
the ESR spectrumI(ω)is given by <cit.>
    I(ω) =  ω H_R^2/8[-𝒢^R_S^+S^-(ω)]

within the linear response theory,
whereH_Ris the amplitude of the external oscillating magnetic field.
On the other hand, when the applied electromagnetic wave is completely unpolarized,
the ESR spectrum is given by (appendix <ref>)

    I(ω) = ω H_R^2/8[ -𝒢^R_S^+S^-(ω) - 𝒢^R_S^-S^+(ω) ].

The result (<ref>) is invariant under any rotation around thezaxis, the direction along which the electromagnetic wave propagates.
Since𝒢^R_S^-S^+(ω) = -𝒢^R_S^+S^-(-ω),
the ESR spectrum (<ref>) is fully determined from the retarded Green's function,𝒢^R_S^+S^-(ω).
Here, we note that the frequencyωis positive since it represents the energy of the photon.


In the remainder of the paper, we deal with the ESR spectrum of Eq. (<ref>) because inclusion of𝒢^R_S^-S^+(ω)has no impact on the conclusions of this paper.
The Green's function𝒢^R_S^+S^-(ω)satisfies the identity <cit.>
    𝒢^R_S^+S^-(ω)
        =  2⟨S^z|⟩/ω - H - ⟨[𝒜, S^-]|⟩/(ω-H)^2
            + 1/(ω-H)^2𝒢^R_𝒜𝒜^†(ω),

whereω-His shorthand forω-H + i0and𝒜is an operator defined by

    𝒜 = [ℋ', S^+].

Equation (<ref>) is an exact relation that simply results from the equation of motion ofS^±<cit.>.
Now it is clear that if there was no anisotropy, i.e.ℋ'=0,
the ESR spectrumI(ω)would contain only the single EPR peak,

    I(ω)=I_ EPR(ω) = π H_R^2H/4⟨S^z|δ⟩(ω-H).

The anisotropic interactionℋ', however small it may be, is always present in the materials, giving rise to the finite linewidth
to the EPR peak (<ref>).
It also yields other resonance peaks.




The𝒜operator determines qualitative and quantitative properties of the ESR spectrum.
In particular, Eq. (<ref>) shows that if there is an additional resonance peak well isolated from the EPR one,
it must come from the retarded Green's function𝒢^R_𝒜𝒜^†(ω).
In fact, many interesting resonance peaks are found in ESR experiments and explained on the basis of the
identity (<ref>) <cit.>.

Let us split the ESR spectrum (<ref>) into two parts: the EPR peakI_EPR(ω)and the other peaksI'(ω)well isolated from the former, if they exist,

    I(ω) = I_ EPR(ω) + I'(ω).

In the leading order of the perturbation, the additional peakI'(ω)is given by

    I'(ω)    = ω H_R^2/8(ω-H)^2[ -  G^R_𝒜𝒜^†(ω)],

where the full retarded Green's function𝒢^R_𝒜𝒜^†(ω)in Eq. (<ref>) 
is replaced to  the unperturbed retarded Green's
functionG^R_𝒜𝒜^†(ω)of the unperturbed system (<ref>).
As far asI'(ω)is concerned,
measuring ESR is equivalent to measuring the resonance of the complex operator𝒜given in Eq. (<ref>).
We emphasize that the resonance peak position ofI'(ω)is determined fully in the unperturbed system.
While the intensity ofI'(ω)is proportional to the square of the small coupling constant ofℋ',
the frequency dependence is fully determined from the unperturbed Green's function.




§ FREQUENCY SHIFT BY THE FQ ORDER


In this section, we discuss effects of the FQ order on the EPR absorption peak, showing that the FQ order parameter can be extracted
experimentally from the EPR frequency shift.



 §.§ Generic case

We start with a general discussion.
Let us consider a perturbative expansion of the identity (<ref>) up to the first order of the perturbation,

    𝒢^R_S^+S^-(ω)
        ≃2⟨S^z|⟩/ω- H( 1- ⟨[𝒜, S^-]|_⟩0/2⟨S^z|_⟩01/ω-H),

where⟨O|_⟩0denotes the average taken in the unperturbed system (<ref>)[We note that, when a certain symmetry is broken spontaneously in the unperturbed system and
it is also broken by the anisotropy in the full Hamiltonian ℋ, one must choose a proper direction of the infinitesimal auxiliary
symmetry breaking field in
the unperturbed system so that the ground state is smoothly deformed by imposing weak perturbation.].
The approximation (<ref>) implies that the resonance frequency of the
EPR peak is shifted fromω_r=Hto

    ω_r     ≃ H - ⟨[𝒜, S^-]|_⟩0/2⟨S^z|_⟩0.

The formula (<ref>) is valid
as long as the anisotropic interaction is perturbative,
the paramagnetic resonance peak is not split into plural peaks by
the perturbationℋ'<cit.>, and⟨[𝒜, S^-]|_⟩0is real.
An example of the splitting is found in one-dimensional quantum antiferromagnets with a uniform
Dzyaloshinskii-Moriya (DM)
interaction <cit.>.
While it is nontrivial whether⟨[𝒜, S^-]|_⟩0is real in the presence of the FQ order,
we show in Sec. <ref> that it is indeed real on the basis of a specific model.


In this paper we denote the frequency shift byδω_r=ω_r-H.
The formula (<ref>) is rephrased as

    δω_r = - ⟨[𝒜, S^-]|_⟩0/2⟨S^z|_⟩0.




A key observation is to notice a fact that a quadratic interactionℋ'usually makes
the commutator[𝒜, S^-]quadratic
in the formula (<ref>).
To see it, we take a generic example,

    ℋ' = ∑_p=x,y,z∑_⟨i,j|_⟩nδ^(n)_p S_i^p S_j^p,

whereS_j^p(p=x,y,z) denotes thep-component of spinSoperators onjth site,⟨i,j|_⟩nrepresents a pair of thenth neighbor sitesiandj.
In particular,⟨i,j|_⟩0meansi=j.
The interaction (<ref>) covers quite a wide range of anisotropic spin-spin interactions found in quantum magnets.
In the following discussions, we choose the most dominant anisotropic interaction.
For quantum spin systems withS=1/2, the most likely anisotropic interaction is the anisotropic exchange interaction on the nearest-neighbor
bond (n=1).
For spin systems withS≥1, the most likely one is the single-ion anisotropy
(n=0).
In principle, the anisotropic exchange interaction on the nearest-neighbor bond can be the dominant one even forS≥1.
Because these two types of anisotropy give rise to resemblant frequency
shifts,
it is difficult to judge whether the dominant anisotropic interaction lives in single-spin sites or in bonds connecting two spin sites.
The resemblance was demonstrated in theS=1Heisenberg antiferromagnetic chain <cit.>.


We assume the following inequalities

    |δ^(n)_z|>|δ^(n)_x|≥ |δ^(n)_y|=0,

without loss of generality, for the most dominant anisotropy.
From Eq. (<ref>), the anisotropic interaction (<ref>) leads to the frequency shift

    δω_r
        =- 2δ^(n)_z - δ^(n)_x/2⟨S^z|_⟩0∑_⟨i,j|_⟩n⟨(3S_i^zS_j^z- S_i · S_j)|_⟩0
     
          -
     δ^(n)_x/2⟨S^z|_⟩0∑_⟨i,j|_⟩n⟨S_i^- S_j^-|_⟩0.

This equation contains the operatorS_i^-S_j^-, which creates or annihilates a pair of magnons.
Since bound magnon pairs condense in the spin nematic phase,∑_⟨i,j|_⟩n⟨S_i^-S_j^-|_⟩0is proportional to the condensed amount of bound magnon pairs
at the wave vectork=0, that is, the FQ order parameter.

The average⟨(3S_i^zS_j^z - S_i ·S_j)|_⟩0on the first line of Eq. (<ref>) is
also nonzero simply because of the magnetic field along thezdirection.
In fact, the frequency shift in several one-dimensional quantum magnets, where the long-range spin nematic order is absent,
was explained on the basis of Eq. (<ref>) without the second
line <cit.>.

Let us introduce a little trick in order to get rid of the contribution unrelated to the spin nematic order.
Here, we rotate the material around theyaxis by an angleθ.
The rotation only affects the form of the anisotropic interaction (<ref>) as

    ℋ'
        = ∑_⟨i,j|_⟩n[ (δ^(n)_z cos^2θ + δ^(n)_xsin^2θ) S_i^zS_j^z
     
         
     + (δ^(n)_z sin^2θ + δ^(n)_x cos^2θ) S_i^xS_j^x
     
          - (δ^(n)_z-δ^(n)_x) sinθcosθ (S_i^zS_j^x + S_i^xS_j^z)
     ].

The interactions in the unperturbed systemℋ_0is invariant under the rotation.
The rotated anisotropic interaction leads to the frequency shift,

    δω_r (θ)
        = (δ^(n)_z-δ^(n)_x)(3cos^2θ -1) +δ^(n)_x/2⟨S^z|_⟩0∑_⟨i,j|_⟩n⟨(3S_i^zS_j^z- S_i · S_j)|_⟩0
    + (δ^(n)_z-δ^(n)_x)sin^2θ +δ^(n)_x/2⟨S^z|_⟩0∑_⟨i,j|_⟩n⟨S_i^-S_j^-|_⟩0
     
           
     -(δ^(n)_z-δ^(n)_x)sinθcosθ/2⟨S^z|_⟩0∑_⟨i,j|_⟩n⟨{2(S_i^zS_j^++S_i^+S_j^z) + 3(S_i^zS_j^-+S_i^-S_j^z)}|_⟩0,

where we added the argument toδω_ron the left hand side to clarify that the frequency shift is a function ofθ.
If the rotation does not affect the direction where the FQ order grows,
the angular dependence comes only out of the coefficients of those averages.
The validity of the assumption is to be confirmed in the next subsection for a specific example.

The EPR frequency shift (<ref>) consists of three parts:
(i) the uniaxial part coupled to the average⟨(3S_i^zS_j^z-S_i ·S_j)|_⟩0,
(ii) the FQ order part coupled to the FQ order parameter⟨S_i^-S_j^-|_⟩0,
and (iii) the off-diagonal part coupled to⟨{2(S_i^zS_j^++S_i^+S_j^z) + 3(S_i^zS_j^-+S_i^-S_j^z)}|_⟩0.
Note again that the operator3S_i^zS_j^z-S_i ·S_jis not an FQ order parameter although it is a quadrupole operator.
That operator has a nonzero expectation value when thezaxis is inequivalent to thexyplane by applying the magnetic field
along thezaxis.
Looking at the angular dependence of the frequency shift (<ref>),
we can separate the FQ order parameter from the other terms as follows.


Let us focus onδω_r(0)andδω_r(π/2)because the second line of Eq. (<ref>) vanishes atθ=0 π/2.
Those frequency shifts are rewritten as

    δω_r(0) + δω_r(π/2)    = δ^(n)_z+δ^(n)_x/2⟨S^z|_⟩0(
     ∑_⟨i,j|_⟩n⟨(3S_i^zS_j^z -  S_i · S_j)|_⟩0 + ∑_⟨i,j|_⟩n⟨S_i^-S_j^-|_⟩0
     ), 
    δω_r(0) - δω_r(π/2)    = δ^(n)_z-δ^(n)_x/2⟨S^z|_⟩0(
     3∑_⟨i,j|_⟩n⟨(3S_i^zS_j^z -  S_i · S_j)|_⟩0 - ∑_⟨i,j|_⟩n⟨S_i^-S_j^-|_⟩0
     ).

Sinceδ^(n)_p(p=x,z) are known parameters from ESR measurements at high enough temperatures out of the FQ phase
and the magnetization⟨S^z|_⟩0is also known independently of the ESR experiments.
Thus, combiningδω_r(0)andδω_r(π/2), we can obtain the FQ order parameter∑_⟨i,j|_⟩n⟨S_i^-S_j^-|_⟩0.
In the uniaxially anisotropic case ofδ^(n)_z≠0andδ^(n)_x=δ^(n)_y=0, the procedure is simplified thanks to
the following simple relation,

    1/2δω_r(0) + δω_r(π/2) = δ^(n)_z/2⟨S^z|_⟩0∑_j ⟨S_i^-S_j^-|_⟩0.


 §.§ Specific case


The important remaining tasks in this section are to confirm that the frequency shift (<ref>) is real and that
the averages in Eq. (<ref>) is invariant under the rotation.
To do so, we take a specific example of anS=1bilinear-biquadratic model with the single-ion anisotropy,

    ℋ   = ∑_n=1,2∑_⟨i,j|_⟩n[ J_n  S_i · S_j + K_n ( S_i · S_j)^2]
     
            - H(S^z cosθ + S^x sinθ)
     
            + D ∑_j (S_j^z)^2 + E∑_j {(S_j^x)^2-(S_j^y)^2},

on the square lattice. Rotating the system abouty-axis by angleθ, we can
redefine the Hamiltonian as

    ℋ   =
     ∑_n=1,2∑_⟨i,j|_⟩n[ J_n  S_i · S_j + K_n ( S_i · S_j)^2]
     - HS^z
     
           
     + ∑_j [(Dcos^2θ + Esin^2θ) (S_j^z)^2 
            + (Dsin^2θ + E cos^2θ) (S_j^x)^2 - E (S_j^y)^2
     
            - (D-E)sinθcosθ (S_j^zS_j^x+S_j^xS_j^z)
     ].

The latter representation is easier to handle.
Here we assume, without loss of generality, that bothDandEhave the same sign.
Note that the form of the anisotropic interaction of Eq. (<ref>) is a special case of Eq. (<ref>).
In the language ofδ^(n)_pin Eq. (<ref>), the parametersDandEcorresponds to them as

    δ^(0)_z = D+E,   δ^(0)_x = 2E.

Note that we also imposed the condition (<ref>).
We assume thatJ_n>0andK_n<0forn=1,2and that the single-ion anisotropy can be seen as a perturbation.

It was shown at the mean-field level that the ground state of the unperturbed model
is in the FQ phase when  both|K_1|/J_1and|K_2|/J_2are large enough <cit.>.
The mean-field FQ ground state|ψ_0⟩is represented as a product state,

    |ψ_0⟩ = ∏_j |ϕ_0(θ_H, φ)⟩_j,

of the local state|ϕ_0(θ_H, φ)⟩_j,

    |ϕ_0(θ_H, φ)⟩_j = i (e^iφcosθ_H |1⟩_j - e^-iφsinθ_H |-1⟩_j),

where|m⟩_jis the eigenstate ofS_j^zwith the eigenvaluem.
The anglesθ_Handφ, which are real, are determined so that the ground-state energy is minimized.
The ground state has the FQ order

    ∑_j⟨ψ_0 | (S_j^-)^2 |ψ_0|=⟩ -N e^i2φsin 2θ_H,

which is in general complex.
Here,Nis the number of spins.
In the following, we show thate^i2φis real in the presence of the single-ion anisotropy as long as it is perturbative.


To discuss the ground state energy of the mean-field FQ state (<ref>),
we introduce the quadrupole operators

    Q_j =
     [  Q_j^x^2-y^2; Q_j^3z^2-r^2;       Q_j^xy;       Q_j^yz;       Q_j^zx ]
     =
     [                     (S_j^x)^2 - (S_j^y)^2; [2(S_j^z)^2 -(S_j^x)^2 - (S_j^y)^2 ]/√(3);                 S_j^x S_j^y + S_j^y S_j^x;                 S_j^y S_j^z + S_j^z S_j^y;                 S_j^z S_j^x + S_j^x S_j^z ]

and rewrite the Hamiltonian (<ref>) as

    ℋ   = ∑_⟨i,j|_⟩n[ J_n/2( S_i · S_j +  Q_i · Q_j)
     
            + J_n-K_n/2 ( S_i · S_j -  Q_i · Q_j)
     ]-HS^z
     
            + 2ND/3 + D(3cos^2θ-1)+Esin^2θ/2√(3) Q^3z^2-r^2
            + Dsin^2θ + E(cos^2θ+1)/2 Q^x^2-y^2
            - (D-E)sinθcosθ Q^zx

withQ ≡∑_j Q_j.
Writing the local state as|ϕ(θ_H, φ)⟩_j=e_1|1⟩_j + e_0 |0⟩_j + e_-1|-1⟩_j,
we can represent the ground-state energyE_FQ = ⟨ψ_0|ℋ|ψ_0|$⟩ as <cit.>


    E_ FQ/N   = 2(J_1+J_2)| e ·e̅|^2
     - 2(J_1+J_2-K_1-K_2) |2e_1e̅_-1 - (e_0)^2|^2
     - H (|e_1|^2-|e_-1|^2)
     
            + 2D/3+ D(3cos^2θ-1)+Esin^2θ/6(|e_1|^2+|e_-1|^2-2|e_0|^2)
     
           
     + Dsin^2θ + E(cos^2θ + 1)/2√(2)(e_-1e̅_1+ e_1 e̅_-1)
     -(D-E)sinθcosθ(e_1e̅_0 + e_0 e̅_1 - e_0 e̅_-1 - e_-1e̅_0),


where e̅_a is the complex conjugate of e_a (a=1,0,-1).
Plugging e_1=ie^iφcosθ_H, e_0 = 0, and e_-1=-ie^-iφsinθ_H into Eq. (<ref>),
we obtain

    E_ FQ/N   = 2(J_1+J_2) -2(J_1+J_2-K_1-K_2)sin^2 2θ_H
     
          -H cos 2θ_H
    + 3D(cos^2θ + 1)+Esin^2θ/6
          +Dsin^2θ + E(cos^2θ+1)/√(2)cos2φsin 2θ_H.

Let us determine θ_H and φ in the spirit of the perturbation theory.
First, when D=E=0, the ground-state energy becomes

    E_ FQ^0/N   = 2(K_1+K_2)- H^2/2H_ sat
        + 2(J_1+J_2-K_1-K_2)(cos 2θ_H - H/H_ sat)^2

with the saturation field H_ sat = 4(J_1+J_2-K_1-K_2).
Since K_n<0<J_n for all n=1,2, the angle θ_H is determined to be

    θ_H^0 = 1/2cos^-1(H/H_ sat)

in 0≤θ_H^0≤π/2 when 0≤ H<H_ sat and θ_H^0=0 when H_ sat<H.
We restrict ourselves to the former case where the system is in the FQ phase.
The ground state |ψ_0⟩ has the FQ order,

    ∑_j ⟨ψ_0 |(S_j^-)^2|ψ_0|=⟩ - Ne^i2φ√(1-(H/H_ sat)^2),

where the angle φ determines the direction of quadrupolar directors.
Next, we treat the single-ion anisotropy perturbatively.
Weak anisotropies D and E do not modify the solution θ_H = θ_H^0 at the lowest
order of the perturbation.
The anisotropy, however small it may be, serve as a symmetry breaking field and
determines φ.
From the energy in Eq. (<ref>) the angle φ is chosen, for arbitrary angle θ, as
cos 2φ = 1 when both D and E are negative,
and cos 2φ = -1 when both are positive.
In any case, since sin 2φ=0, we thus find e^i2φ = cos 2φ∈ℝ
and that the FQ order parameter ∑_j ⟨(S_j^-)^2|$⟩ is real, satisfying the relation

    ∑_j ⟨(S_j^-)^2|_⟩0=∑_j ⟨Q^x^2-y^2_j|_⟩0,

for the model (<ref>) at low enough temperatures. The FQ order parameter is thus given by

    ⟨Q^x^2-y^2|_⟩0
         ≡∑_j ⟨Q^x^2-y^2_j|_⟩0 
       = - N cos 2φ√(1-(H/H_ sat)^2).


At the same time, we can also conclude that the obtained solutions ofφandθ_Hare independent of the angleθof the rotation.
That is, the FQ order parameter is invariant under the rotation around theyaxis.
Also, the averages∑_j ⟨(3(S_j^z)^2-2)|_⟩0=Nand∑_j ⟨2(S_j^zS_j^++S_j^+S_j^z)+3(S_j^zS_j^-+S_j^-S_j^z)|_⟩0=0turn out to be independent ofθ.
In particular, the latter average vanishes because it involves the creation and the annihilation of the gapped unpaired magnon.
The assumption made in deriving Eq. (<ref>) in the previous subsection is thus justified.

Finally, the frequency shift (<ref>) of the model (<ref>) at zero temperature becomes

    δω_r
        ≃(D-E)(3cos^2θ-1)+2E/2⟨S^z|_⟩0/N
            - (D-E)sin^2θ + 2E/2⟨S^z|_⟩0⟨Q^x^2-y^2|_⟩0.

We emphasize that the value of the FQ order parameter
is determined independent of perturbative anisotropiesDandE.
Those anisotropies merely fix the angleφ, which was spontaneously determined in the unperturbed system.

Although we do not show that the frequency shift is real for generic cases,
it is reasonable to expect that the result holds true generally as long as the mean field approximation is applicable and
the perturbative anisotropic interaction is given by either the single-ion anisotropy or the anisotropic exchange interaction (<ref>).







 §.§ Experimental applications




For the analysis of experimental data, we comment on the comparison between the frequency shift and the linewidth of the EPR peak.
Should the linewidth be larger than the magnitude of the EPR frequency shift,
the frequency shift would be undetectable.
However, this is not the case as long as the anisotropic interaction is small enough compared to the isotropic ones.
According to the generic perturbation theory (<ref>),
the EPR frequency shift is of the first order of the perturbative anisotropic interaction,
whereas the linewidth is of the second order because it is
proportional toImG^R_𝒜𝒜^†(ω=H)<cit.>.
Therefore, given the long-range FQ order is well developed,
the linewidth is small enough not to mask the frequency shift (<ref>).

We conclude the section, referring to chromium spinel oxides.
Chromium spinelsACr_2O_4(A=Zn, Cd, Hg) are considered as anS=3/2pyrochlore Heisenberg antiferromagnet with biquadratic interactions <cit.>.
Recently it was pointed out that those chromium spinels can have the FQ phase next to the fully polarized phase <cit.>. Indeed, high-field measurements discovered the presence of a classically
unexpected phase just below the fully polarized phase <cit.>.
We emphasize that the result (<ref>) is applicable to those interesting compounds.
In fact, in deriving Eq. (<ref>),
we only specified the anisotropic interactionℋ'and specified neither the spin quantum number nor
the form of the SU(2)-symmetric interaction in the Hamiltonian (<ref>).






§ EMERGENCE OF ANOTHER PEAK BY THE AFQ ORDER: MAGNON-PAIR RESONANCE


The frequency shift (<ref>) is insensitive to the AFQ order
since the AFQ order parameter has alternating sign depending on the position.
In this section, we develop an alternative way for detecting the AFQ order, studying the additional
absorptionI'(ω)in ESR spectrum (<ref>).
Here, the point is that the𝒜operator creates a magnon pair excitation.



 §.§ Concept


Let us explain the concept of our method, taking the following example of the anisotropic interaction:

    ℋ' =  δ' ∑_⟨i,j|_⟩1 (S_i^x S_j^z + S_i^z S_j^x).

Here, we do not specify the precise form ofℋ_SU(2).
We just assume that the unperturbed system (<ref>) has the AFQ ground state.
The𝒜operator for Eq. (<ref>) is given by

    𝒜 = δ' ∑_⟨i,j|_⟩1 (S_i^+S_j^+ -3S_i^z S_j^z +  S_i · S_j).

The first term of Eq. (<ref>) creates or annihilates the magnon pair formed on the nearest-neighbor bond.
This shows that dynamics of the bound magnon pair is directly observable in ESR experiments throungh the
relation (<ref>).


In the spin nematic phase with the AFQ order, the bound magnon pair can acquire an excitation gap, sayΔ, atk=0,
though it becomes gapless at the wave vector of the AFQ ordered ground state.
If so, the ESR spectrum in the spin nematic ordered phase contains a sharp resonance peak whose
frequency corresponds to the excitation gap of the bound magnon pair
atk=0.
As shown in Eq. (<ref>), the ESR spectrum (<ref>) contains a contribution of
the𝒜operator through the imaginary part of the retarded Green's function of𝒜.
Since𝒜involves the creation and annihilation operators of the bound magnon pair,
the ESR spectrum will have a resonance peak atω=Δ, i.e.I'(ω)∝𝒜_MPRδ(ω- Δ).
We call this peak the magnon-pair resonance peak.
In Secs. <ref> and <ref>, we confirm that the ESR spectrum indeed yields the sharp magnon-pair resonance
by taking an example of anS=1spin model on the square lattice.


The interaction (<ref>) results from the rotation
of a spin anisotropy, e.g.δ^(z)∑_⟨i,j ⟩_1S_i^z S_j^z, 
around theyaxis as we did in Eq. (<ref>).
For a general rotation angleθ, the resultant magnon pairing operator of the rotated interaction
(<ref>) shows the angular dependence ofsinθcosθ.
This gives a characteristic angle dependence in the intensity𝒜_MPRof the magnon-pair
resonance peak. [See Eq. (<ref>).]
This angular dependence was a key to characterize qualitatively the quadrupolar liquid state
by using ESR <cit.>.







 §.§ Model


To flesh out the general discussion described in Sec. <ref>,
we take an example of anS=1spin model on the square lattice given by the following unperturbed Hamiltonian

    ℋ_0
     =   ∑_⟨i,j|_⟩1 J_11[ S_i · S_j + ( S_i · S_j)^2]
     
        - ∑_⟨i,j|_⟩2{ J_12[ S_i · S_j + ( S_i · S_j)^2]
     +J_22( S_i · S_j)^2}
        -HS^z

and the perturbative single-ion anisotropy

    ℋ'
     =   ∑_i[(Dcos^2θ + E sin^2 θ) (S_i^z)^2
     
        + (Dsin^2 θ + E cos^2 θ) (S_i^x)^2 - E(S_i^y)^2
     
       - (D-E)sinθcosθ (S_i^z S_i^x + S_i^x S_i^z)
     ].

This form of the anisotropy is the same as in Eq. (<ref>), which
is obtained by rotating the common formD(S_i^z)^2 + E {(S_i^x)^2 - (S_i^y)^2}abouty-axis by angleθ.
We assume that the couplingsJ_11,J_12, andJ_22are all positive.
The bilinear-biquadratic interaction in Eq. (<ref>) has the SU(3) symmetry whenJ_22=0and its low-energy behavior was closely investigated in Ref. <cit.>.
It is easy to confirm that the unperturbed Hamiltonian (<ref>) is written as

    ℋ_0
        = J_11/2∑_⟨i,j|_⟩1  ( Q_i · Q_j +  S_i · S_j)
     
          - ∑_⟨i,j|_⟩2[ J_12/2 ( Q_i · Q_j +  S_i · S_j)
     
          +J_22/2 ( Q_i · Q_j -  S_i · S_j)
     ]
      - HS^z,

by using the quadrupole operatorsQ_jof Eq. (<ref>).



We employ the model (<ref>) for the following reasons.
First, the model (<ref>) does not suffer from a known technical problem of the linear flavor-wave theory, that is,
violation <cit.> of the frequency sum rule  <cit.>
of the dynamical structure factor ∫_0^∞ dωωS^αα(k=0,ω)=0,
which holds true for the unperturbed Hamiltonian (<ref>)at zero magnetic field H=0.
Here,S^αα(k=0,ω)=-Im G^R_S^αS^α(ω)/N.
The linear flavor-wave theory does not always satisfy the exact sum rule but
can be recovered by including three- and four-particle interactions <cit.>.
It is the great advantage of the model (<ref>) that we can omit such a complicated procedure.
This sum rule at zero magnetic field is related to the exact result of the EPR frequencyω=Hin the ESR spectrum (<ref>).
We will briefly comment about an example of the violation of the exact result (<ref>) in the linear flavor-wave theory
in the last part of Sec. <ref>.
Second, the model (<ref>) exhibits the AFQ phase in quite a wide field range
up to the saturation field [H_c^AFQin Eq. (<ref>)].
We will come back to those points in Secs. <ref> and <ref>.
Last but not least, theS=1model is closely related to anS=1/2frustrated ferromagnetic model
on a square lattice <cit.>.
We emphasize that the results about the AFQ phase obtained in the present paper also hold true for theS=1/2model.
Here, the single-ion anisotropy (<ref>) in theS=1model is translated into anisotropic exchange interactions
in theS=1/2model.  See Sec. <ref> for further discussions.


Following the general discussion in Sec. <ref>,
we study the additional peakI'(ω)given by Eq. (<ref>).
Here
we only need to derive the unperturbed retarded Green's functionG^R_𝒜𝒜^†(ω)of the operator𝒜=[ℋ',S^+].
For that purpose,
we use the linear flavor-wave theory <cit.>
for the unperturbed system (<ref>) at zero temperature.




 §.§ Mean-field ground state


As well as the spin-wave theory,
the flavor-wave theory is developed by taking into account quantum fluctuations around the ordered state.
We need to start with deriving the AFQ ground state of the bilinear-biquadratic model (<ref>)
in a site-decoupled semi-classical approximation.
Let us denote its mean-field ground state by|ψ_0⟩.
Under a strong enough magnetic field,|ψ_0⟩is in the fully polarized phase and exactly given by

    |ψ_0⟩ = ∏_j i |1⟩_j.

As the magnetic field is decreased, the ground state of the model (<ref>) enters into a partially polarized phase.
The mean-field ground state of the partially polarized phase has AFQ order, which is given in the form

    |ψ_0⟩= ∏_j |ϕ_0(e^i k_M · r_jθ_H, φ)⟩_j

with

    |ϕ_0(θ_H, φ)⟩_j = i (e^iφcosθ_H |1⟩_j - e^-iφsinθ_H |-1⟩_j)

for0 ≤θ_H ≤π/2and0 ≤φ< π.
Herek_M = (π, π)is the wave vector where the AFQ order grows,r_jspecifies the location of the spinS_j, ande^ik_M ·r_j=±1gives the staggered sign depending on the sublattice.
The factorsinθ_Hrepresents the fraction of the condensed bound magnon pair.
The state|ϕ_0(θ_H,φ)⟩_jis also given by an SU(3) rotation of the polarized state,

    |ϕ_0(θ_H, φ)⟩_j=i exp(i φ S_j^z)exp(i θ_H Q_j^xy)|1⟩_j.

As we did in Sec. <ref>, we determineθ_Hfor the unperturbed system and then studyφdependence perturbatively.


For the unperturbed Hamiltonianℋ_0, the mean-field ground-state energy of the AFQ state
is
given by <cit.>
    E^0_ AFQ/N   =  -2(J_12+J_22)
     + 2(J_11+J_22) cos^2 2θ_H
     
         
     - H cos 2θ_H,

whereNis the number of sites.
This energy is minimized atθ_H=0forH≥H_c^AFQ,
whereH_c^AFQdenotes the saturation field

    H_c^ AFQ = 4(J_11+J_22),

and at

    θ_H =  1/2cos^-1(H/H_c^ AFQ)

forH<H_c^AFQ.


In the AFQ phase, the saturation fieldH_c^AFQis the critical field where bound magnon pairs
start condensing when the field is decreased.
Another mean-field solution is an antiferromagnetically ordered state given by a staggered SU(2) rotation
of the polarized state

    |ϕ_0(e^i k_M · r_jθ_H,φ)⟩_j=i exp(i φ S_j^z) exp( i e^i k_M · r_jθ_H S_j^y)|1⟩_j.

In this solution, the saturation field is

    H_c^ AFM = 2(6J_11+J_22),

where the single magnon closes the energy gap.
WhenHis decreased, ifH_c^AFQ> H_c^AFM, bound magnon pairs condense below the saturation field
before the unpaired magnon does.
The comparison betweenH_c^AFQandH_c^AFMshows that the AFQ phase is realized
for

    4J_11 < J_22

at least near the saturation.
In Sec. <ref>, we see that inclusion of the quantum fluctuation relaxes the condition (<ref>) to

    J_11 < J_22.



The anisotropic perturbationℋ^'[Eq. (<ref>)] can make the ground-state energy depend onφ.
In the mean-field calculation, the first-order perturbation changesE^0_AFQ/Nto

    E_ AFQ/N   = -2(J_12+J_22) +2(J_11+J_22)cos^2 2θ_H
     
       -H cos 2θ_H
    + 3D(cos^2θ + 1)+Esin^2θ/6.

In contrast to the case of the FQ state (<ref>),
the mean-field energy (<ref>) of the AFQ state  does not depend on the angleφdespite
the equivalent form of the anisotropic interaction in these two systems.
Thisφindependence shows that the angleφis undetermined in the AFQ state at the mean-field level.
To determineφ, we need to go beyond the mean-field level.
We leave it undetermined sinceφis insignificant for our calculations below.





 §.§ Linear flavor-wave theory in the fully polarized phase


In this subsection we discuss the flavor-wave theory <cit.> in the fully polarized phase
above the saturation fieldH_c^AFQof the AFQ phase.
Later in Sec. <ref> we discuss the flavor-wave theory in the AFQ phase belowH_c^AFQ.

The flavor-wave theory is formulated in terms of Schwinger bosons.
Let us denote the creation and annihilation operators of the Schwinger bosons at thejth site byb_j,m^†andb_j,m, respectively.
The flavor indexm=1,0,-1corresponds to the eigenvalue ofS_j^z.
Using these bosons, an operatorO_jat thejth site is written as

    O_j = ∑_m, m'=1,0,-1 b_j,m^†Õ_j^mm' b_j,m'

whereÕ_jdenotes a3 ×3matrix whose element is given byÕ_j^mm'=_j⟨m|O_j|m' ⟩_j.
The Schwinger bosons are subject to the constraint

    ∑_m=1,0,-1 b_j,m^†  b_j,m = 1

for every sitej.


The3×3matrixÕ_jis easily found by writing|1⟩=(1  0   0)^T,|0⟩=(0  1  0)^T, and|-1⟩=(0  0  1)^T.
For example, the matrices for the spin operatorsS_j^α(α=x,y,z)are given by

    S̃_j^x    = 1/√(2)[ 0 1 0; 1 0 1; 0 1 0 ],  S̃_j^y = i/√(2)[  0 -1  0;  1  0 -1;  0  1  0 ], 
    S̃_j^z    =
      [  1  0  0;  0  0  0;  0  0 -1 ].

The matrix representation ofQ_iis easily obtained by combining the matrices ofS_i(<ref>).


The fully polarized phase can be seen as a condensation phase of theb_j,1bosons.
We can replace the operatorsb_j,1andb_j,1^†by ac-number,  the fraction of the condensed boson,

    b_j,1=b_j,1^† = √(1-a_j^† a_j - b_j^† b_j).

Here,  we rewrote boson operators as

    b_j,0 = b_j,    b_j,-1 = a_j,

to simplify the notation.
We thus end up with the Schwinger boson representation of the spin operator,

    S_j^x
        = 1/√(2)(√(1-a_j^† a_j - b_j^† b_j)  b_j + b_j^†√(1-a_j^† a_j - b_j^† b_j)
     + a_j^† b_j + b_j^† a_j
     ),
     
    
     S_j^y
        = -i/√(2)(√(1 - a_j^† a_j - b_j^† b_j)  b_j - b_j^†√(1 - a_j^† a_j - b_j^† b_j)
     + b_j^† a_j - a_j^† b_j),
     
    
     S_j^z
        = 1 - 2a_j^† a_j - b_j^† b_j.

As it is expected, the ladder operator S_j^- = S_j^x - i S_j^y involves the creation operator b_j^†.
Likewise, (S_j^-)^2 = Q_j^x^2-y^2 - i Q_j^xy involves the creation operator a_j^†.
In fact, the Schwinger boson representation of Q_j is
 
    Q_j^x^2-y^2   = √(1-a_j^† a_j^† -b_j^† b_j)  a_j + a_j^†√(1-a_j^† a_j - b_j^† b_j),
      
    
      Q_j^3z^2-r^2   = 1/√(3) (1-3b_j^† b_j),
      
    
      Q_j^xy   = i ( a_j^†√(1-a_j^† a_j - b_j^† b_j) - √(1-a_j^† a_j - b_j^† b_j)  a_j),
      
    
      Q_j^yz   = - i/√(2)( √(1-a_j^† a_j - b_j^† b_j)  b_j - b_j^†√(1-a_j^† a_j - b_j^† b_j)
      + a_j^† b_j - b_j^† a_j),
      
    
      Q_j^zx   = - 1/√(2)[ - (√(1-a_j^† a_j - b_j^† b_j)  b_j + b_j^†√(1-a_j^† a_j - b_j^† b_j))
      +a_j^† b_j + b_j^† a_j
      ].

Up to the linear order of creation and annihilation operators,
onlyQ_j^x^2-y^2≃a_j^†+ a_jandQ_j^xy≃i(a_j^†- a_j)can create the “a” boson that corresponds to the bound magnon pair.


Up to the quadratic order of the creation and annihilation operators,
the unperturbed Hamiltonian (<ref>) turns effectively into

    ℋ_0
        = ∑_ k[ ω_a( k) a_ k^† a_ k + ω_b( k)b_ k^† b_ k],

wherea_kandb_kare the Fourier transforms ofa_jandb_j, respectively, andω_a,b(k)are given by

    ω_a( k)
        = -4 [J_11(1-) - (J_12+J_22)(1-)]
     
          - 8J_22 + 2H,
     
    ω_b( k)
        = -4J_11 (1-) + 4J_12 (1-) + H,

with

    = 1/2(cos k_x + cos k_y),
     
       = cos k_x cos k_y.

In the parameter range

    J_22 > J_11,

the singlea-boson state atk=k_Mhas the lowest eigenenergy when the magnetic field is close to
the saturation field,
while theaboson has the larger excitation energy than thebboson at extremely strong fields.

Two kinds of bosons have the following excitation gaps:

    ω_a( k_M)
        = 2(H-H_c^ AFQ),
    ω_b( k_M)
        =  4(J_22-J_11) + H-H_c^ AFQ.

Thus, the bound magnon pair (aboson) condenses atH=H_c^AFQwhile the unpaired magnon (bboson) remains gapped when
the condition (<ref>) is satisfied.
In contrast, they both have excitation gaps atk=0, which we denote byΔ_aandΔ_b,

    Δ_a    ≡ω_a( 0)= 8J_11 + 2(H-H_c^ AFQ),
     
    Δ_b    ≡ω_b( 0)= H.

Within the framework of the linear flavor-wave theory, the EPR (<ref>) of
the unperturbed system (<ref>)
at temperaturesT≪His understood as the excitation of thebboson from the ground state.





 §.§ Magnon-pair resonance in the fully polarized phase


Here we study the additional ESR spectrumI'(ω), given by Eq. (<ref>), in the fully polarized phase.
Using the linear flavor-wave theory (<ref>), we  evaluate the retarded Green's functionG^R_𝒜𝒜^†(ω).
The𝒜operator determined from the rotated single-ion anisotropy (<ref>) is

    𝒜   =∑_j [ (Dcos^2θ + E sin^2 θ) (Q_j^zx+iQ_j^yz)
     - (Dsin^2θ + Ecos^2θ) Q_j^zx
     +i E Q_j^yz
            - (D-E) sinθcosθ (Q_j^x^2-y^2+iQ_j^xy -√(3)Q_j^3z^2-r^2)
     ].

Up to the linear order of the creation and the annihilation operators,
it is approximated as

    𝒜/√(N)   ≃√(2)(Dcos^2θ + Esin^2θ) b_ k= 0^†
     - 1/√(2)(Dsin^2θ + Ecos^2θ) (b_ k= 0+b_ k= 0^†)
     +E/√(2) (b_ k= 0 - b_ k= 0^†)
     
            - 2(D-E)sinθcosθ  a_ k= 0.

All the terms in the first line of Eq. (<ref>) yield the EPR peak.
The term in the second line, containing theaboson operator, yields the delta-function magnon-pair resonance peakδ(ω-Δ_a)atω=Δ_a.
According to Eq. (<ref>),
the slope of the resonance frequencyω=Δ_aas a function ofHis double of that of the EPR one (<ref>)
because the “a” boson creates the magnon pair and the “b” boson creates the single unpaired magnon.

Equation (<ref>) also indicates that the intensity of the magnon-pair resonance peak shows the angular dependence ofsin^2θcos^2θ:

    I'(ω)
        ≃N(D-E)^2H_R^2/2Δ_a/(Δ_a-H)^2sin^2θcos^2θδ(ω-Δ_a).


 §.§ Linear flavor-wave theory in the AFQ phase


We move on to the discussion of the linear flavor-wave theory in the AFQ phase <cit.>.
Since the angleφof the AFQ directors is not pinned by the anisotropy in the mean-filed approximation,
we consider the AFQ state for the generalφ.
We note that this degeneracy is not lifted even by the first order perturbation of the anisotropy in the linear flavor-wave approximation
as shown in Appendix <ref>. We leave it as an open question to determineφbecause the determination ofφhas little impact on our conclusions in this paper, as shown in Sec <ref>.


In the fully polarized phase, we took into account low-energy excitations from the fully polarized state
by replacing a local base|1⟩with either|0⟩or|-1⟩.
In the AFQ phase, the mean-field ground state [Eq. (<ref>)]
is obtained from the fully polarized state by performing an alternate SU(3) rotation (<ref>),

    |ψ_0⟩ = ∏_j iR̃(e^i k_M · r_jθ_H, φ) |1⟩_j,

where the matrixR̃(e^ik_M ·r_jθ_H, φ)is expressed  as

    R̃(e^i k_M · r_jθ_H, φ) = exp(iφS̃_j^z)
     exp(ie^i k_M · r_jθ_H Q̃_j^xy)

with

    exp(iφS̃_j^z)
        =
     [  e^iφ     0     0;     0     1     0;     0     0 e^-iφ ], 
    exp(ie^i k_M · r_jθ_H Q̃_j^xy)
         =
      [               cosθ_H                    0  e^i k_M · r_jsinθ_H;                    0                    1                    0; -e^i k_M · r_jsinθ_H                    0               cosθ_H ].

In this representation, excitations above the AFQ state are formally described
by local replacements of|1⟩to either|0⟩or|-1⟩, similar to
the case of the fully polarized phase.
As well as the ground state (<ref>), the Schwinger boson representation of an operatorO_jis given by the SU(3) rotation of Eq. (<ref>),

    O_j = ∑_m,m' b_i,m^† [R̃^†(e^i k_M · r_jθ_H) Õ_j R̃(e^i k_M · r_jθ_H)]^mm' b_i,m'.

Forφ=0,
the spin operatorS_jand the quadrupole operatorQ_jin the AFQ phase are related to the  Schwinger boson representation
in the FP phase, given in Eqs. (<ref>)–(<ref>), as follows:

    (
       [  S_j^x; Q_j^zx;        ])    = 1/√(2)(
       [               cosθ_H -e^i k_M · r_jsinθ_H;  e^i k_M · r_jsinθ_H               cosθ_H;                      ]) (
       [ b_j+b_j^† + b_j^† a_j+ a_j^† b_j; b_j+b_j^† - b_j^† a_j- a_j^† b_j;                                  ]),

    (
       [  S_j^y; Q_j^yz;        ])    = 1/√(2)i(
       [               cosθ_H  e^i k_M · r_jsinθ_H; -e^i k_M · r_jsinθ_H               cosθ_H;                      ]) (
       [ b_j - b_j^† + b_j^† a_j - a_j^† b_j; b_j - b_j^† - b_j^† a_j + a_j^† b_j;                                     ]),

    (
       [       S_j^z; Q_j^x^2-y^2;             ])    =
     (
       [                cos 2θ_H   e^i k_M · r_jsin 2θ_H; - e^i k_M · r_jsin 2θ_H                cos 2θ_H;                         ]) (
       [ 1-2a_j^† a_j - b_j^† b_j;              a_j + a_j^†;                          ]),

    (
       [       Q_j^xy; Q_j^3z^2-r^2;              ])    =
    (
       [     -i  a_j + i a_j^†; 1/√(3) (1-3b_j^† b_j);                       ]).

Forφ≠0, they are rotated as

    (
      [ S_j^x (φ); S_j^y (φ);           ])
    =(
      [  cosφ  sinφ; -sinφ  cosφ;       ])
    (
      [ S_j^x (0); S_j^y (0);           ]),
       (
      [        S_j^z (φ); Q_j^3z^2-r^2 (φ);                  ])
    =(
      [        S_j^z (0); Q_j^3z^2-r^2 (0);                  ]),
       (
      [ Q_j^x^2-y^2 (φ);      Q_j^xy (φ);                 ])
    =(
      [  cos 2φ  sin 2φ; -sin 2φ  cos 2φ;         ])
    (
      [ Q_j^x^2-y^2 (0);      Q_j^xy (0);                 ]),
       (
      [ Q_j^zx (φ); Q_j^yz (φ);            ])
    =(
      [  cosφ  sinφ; -sinφ  cosφ;       ])
    (
      [ Q_j^zx (0); Q_j^yz (0);            ]).



Up to the quadratic terms, the unperturbed Hamiltonian (<ref>) is split into two parts,

    ℋ_0
     ≃ℋ_0^a + ℋ_0^b,

where

    ℋ_0^a
        = ∑_ k[ A_ k a_ k^† a_ k + B_ k/2 (a_ k^† a_- k^† +
    a_ ka_- k)],
     
    ℋ_0^b
        = ∑_ k[ C_ k b_ k^† b_ k
     + D_ k/4 (b_ k+ k_M^† b_- k^†
     +b_ k+ k_Mb_- k)].

The parametersA_k,B_k,C_k, andD_kare given by

    A_ k   = 4J_11sin^2 2θ_H - 4J_11 (1-) cos^2 2θ_H 
          + 4J_12(1-) + 4J_22sin^2 2θ_H 
          - 4J_22(1+) cos^2 2θ_H + 2H cos 2θ_H,
     
    
     B_ k   = -4J_11sin^2 2θ_H + 4J_22sin^2 2θ_H,
     
    
     C_ k   = - 4J_11cos^2 2θ_H + 4J_11cos 2θ_H
     
          + 4J_12 (1-) + 4J_22sin^2 2θ_H
     
          +Hcos 2θ_H,
     
    
     D_ k   = -4J_22sin 2θ_H.


To diagonalize the Hamiltonians (<ref>) and (<ref>), we perform the following Bogoliubov transformations,

    [    a_ k; a_- k^† ]   =
     [  coshΘ_ k^a -sinhΘ_ k^a; -sinhΘ_ k^a  coshΘ_ k^a ][    α_ k; α_- k^† ],
     
    [ b_ k+ k_M;   b_- k^† ]   =
     [   coshΘ_ k^b - sinhΘ_ k^b;  -sinhΘ_ k^b   coshΘ_ k^b ][ β_ k+ k_M;   β_- k^† ].

The parametersΘ_k^aandΘ_k^bare determined in order to eliminate the off-diagonal terms:

    Θ_ k^a    = 1/2tanh^-1( B_ k/A_ k),
     
    Θ_ k^b    = 1/2tanh^-1( 2D_ k/C_ k+C_ k+ k_M).

The Bogoliubov transformations diagonalize the Hamiltonian (<ref>) to

    ℋ_0
      =∑_ k[ ω_a( k)α_ k^†α_ k
      +ω_b( k) β_ k^†β_ k]

with the following dispersion relations,

    ω_a( k)
        = √(A_ k^2 - B_ k^2),
     
    ω_b( k)
        = C_ k- C_ k+ k_M/2 + √((C_ k+C_ k+ k_M/2)^2 - D_ k^2).

Inheriting the terminology in the fully polarized phase, we call the bosons created byα_k^†andβ_k^†as
the “a” boson and the “b” boson, respectively, also in the AFQ phase.
Atk=k_M, the “a” boson corresponding to the bound magnon pair is gapless,
whereas
the “b” boson corresponding to the unpaired magnon is gapped:

    ω_a( k_M)    = 0, 
    ω_b( k_M)    = J_22-J_11/J_22+J_11H.

The “a” boson is the characteristic Nambu-Goldstone boson that accompanies the AFQ
ordered ground state.
Atk=0, both excitations are
gapped:

    Δ_a    = ω_a( 0) = 8J_11[ 1+  J_22-J_11/J_11{1-(H/H_c^ AFQ)^2}]^1/2,
     
    Δ_b    = ω_b( 0) = H.


The linear flavor-wave theory reproduces the exact EPR frequencyω=H[Eq. (<ref>)]
of the unperturbed system (<ref>) both in the fully polarized phase and in the AFQ phase.
The reproduction of the exact result
is an important criterion of appropriateness of the low-energy effective theory.
The criterion is akin to the sum rule mentioned in Ref. <cit.>.
If we include an SU(2)-symmetric but SU(3)-asymmetric interaction,

    ℋ_ a = -J_21/2∑_⟨i,j|_⟩1(  Q_i · Q_j -  S_i · S_j )

to the unperturbed Hamiltonian (<ref>), the linear flavor-wave theory fails to reproduce the exact EPR frequencyω=HbecauseΔ_bin the AFQ phase is modified to

    Δ_b = 4J_11cos 2θ_H + 4J_22√(1-(J_22-J_21/J_22)^2 sin^2 2θ_H).

This technical problem is an artifact of the linear flavor-wave theory and not essential to our purpose of demonstrating
the general properties of the magnon-pair resonance in the AFQ phase.
Thus, we put aside this probably complicated discussion, restricting ourselves to the model withJ_21=0.





 §.§ Magnon-pair resonance in the AFQ phase


Here we study the additional ESR spectrumI'(ω), given by Eq. (<ref>), in the AFQ phase.
We approximate the operator𝒜up to the quadratic order of the spin operators,
that is, the quadratic order of theβoperators and the linear order of theαoperators.
The𝒜operator (<ref>) is represented as

    𝒜/√(N)≃   3(D-E)cos2θ + D+3E/2√(2) e^-iφ{(cosθ_H coshΘ_ 0^b -sinθ_H sinhΘ_ 0^b)
     β_ 0 + (-cosθ_HsinhΘ_ 0^b + sinθ_H coshΘ_ 0^b)β_ k_M^†}
        -(D-E)cos2θ - (D+3E)/2√(2)e^iφ{ (cosθ_H coshΘ_ 0^b -sinθ_H sinhΘ_ 0^b)
     β_ 0^† + (-cosθ_HsinhΘ_ 0^b + sinθ_H coshΘ_ 0^b)β_ k_M}
        - 2(D-E)sinθcosθ e^-i2φ{cos^2θ_H (α_ 0coshΘ_ 0^a-α_ 0^†sinhΘ_ 0^a)
     -sin^2θ_H(α_ 0^†coshΘ_ 0^a - α_ 0sinhΘ_ 0^a)
     }
        -(D-E)sinθcosθ∑_ k{e^-i2φsin2θ_H (β_ k+ k_M^†coshΘ_ k^b - β_- ksinhΘ_ k^b)(β_ kcoshΘ_ k^b - β_- k+ k_M^†sinhΘ_ k^b)
     
            + 3(β_ k^†coshΘ_ k^b - β_- k+ k_MsinhΘ_ k^b)
     (β_ kcoshΘ_ k^b - β_- k+ k_M^†sinhΘ_ k^b)
     }.

The term containing theαoperators creates
theaboson atk=0,
whereas the linear terms of theβoperators create thebboson
at eitherk=0ork=k_M.
The quadratic term of theβoperators contributes to the two-magnon continuum made of two
scatteringbbosons.
While thebboson excitation atk=0is involved in the EPR peakI_EPR(ω)with the resonance frequencyω= H,
theaboson excitation atk=0gives rise to the magnon-pair resonance peakI_MPR(ω).
Thebboson excitation atk= k_Mand two-unpaired-magnon excitation result in unpaired magnon resonance peakI_k_M(ω)and broad two-magnon continuumI_2-mag(ω), respectively.
In total, the ESR spectrumI(ω)=I_EPR(ω)+I'(ω)contains three sharp peaks
and a broad continuum:

    I(ω) = I_ EPR(ω) + I_ MPR(ω) + I_ k_M(ω) + I_2-mag(ω),

which are given at the leading order of the perturbation by

    I_ EPR(ω)
        ≃π H/4⟨S^z|_⟩0δ(ω-H),
     
    
     I_ MPR(ω)
        ≃𝒜_ MPRδ(ω - Δ_a),
     
    
     I_ k_M(ω)
        ≃𝒜_ k_Mδ(ω-ω_b( k_M)),
     
    
     I_2-mag(ω)
        ≃∑_ kF( k) [sin^22θ_H {δ(ω-ω_b(- k)-ω_b( k))
            +δ(ω-ω_b(- k)+ω_b(- k+ k_M))
     }
            + 18δ(ω-ω_b(- k+ k_M)-ω_b( k))
     ].

The intensities𝒜_MPRand𝒜_k_Mare given by

    𝒜_ MPR/N   =π/8(D-E)^2sin^2θcos^2θΔ_a/(Δ_a-H)^2
    (cos^2θ_H coshΘ_ 0^a +sin^2θ_H sinhΘ_ 0^a)^2,
     
    𝒜_ k_M/N   = π/16ω_b( k_M)/2(ω_b( k_M)-H)^2((D-E)cos2θ-D-3E/2)^2 (sinθ_HcoshΘ_ 0^b-cosθ_HsinhΘ_ 0^b)^2 ,

and the factor F( k) is

    F( k)/N   = 1/8(D-E)^2sin^2θcos^2θ(sinh2Θ_ k^b/2)^2

within the lowest-order perturbation theory.
Note that these are independent of the angle φ of the quadrupolar order.
This φ independence comes as a consequence of the linear flavor-wave approximation and the
first-order perturbation.
Higher-order processes can induce φ dependent corrections to the above results.


The intensities 𝒜_ MPR  and 𝒜_ k_M
can be rephrased as

    𝒜_ MPR/N   = π/8(D-E)^2sin^2 θcos^2 θΔ_a/(Δ_a-H)^2{2H/H_c^ AFQ + [(1+(H/H_c^ AFQ)^2)
     (J_22+J_11/2J_11 - J_22-J_11/2J_11(H/H_c^ AFQ)^2)
     
         
     +4(J_22-J_11)/Δ_a(1-(H/H_c^ AFQ)^2)^2 ]
     },
     
    𝒜_ k_M/N   = π/16((D-E)cos2θ - D-3E/2)^2J_22^2-J_11^2/2J_11^2H(H_c^ AFQ/H - H/H_c^ AFQ).



The presence of the magnon-pair resonance peakI_MPR(ω)is a direct consequence
of the quadrupolar order in the ground state.
We found that the magnon-pair resonance peak appears at the finite
frequencyω=Δ_ain the AFQ phase,
which is continuously connected to the
magnon-pair resonance peak found in the fully polarized phase.
Reflecting the condensation of bound magnon pairs, the field dependence of
the resonance frequencyω=Δ_ashows a singular upturn at the saturation fieldH=H_c^AFQ(Fig. <ref>).
In addition, there is another peakI_k_M(ω)[Eq. (<ref>)] which is absent in the fully polarized phase.
This peakI_k_M(ω)corresponds to creation of the single “b” boson atk=k_M.
Although ESR usually involves excitations atk=0only as Eq. (<ref>) shows,
the AFQ order with the wave vectork_Mmakes the resonance atk=k_Mpossible.


 In general, the magnon-pair resonance peakI_MPR(ω)could be masked by the broad two-magnon continuum.
However there is a better chance to observe this resonance peak near the saturation field.
The lowest energy of the continuum takes the highest value2ω_b (k_M)=2(J_22-J_11)H_c^AFQ/(J_22+J_11)at the saturation field.
For the parameter rangeJ_22 > 1.702 J_11, this lower edge of the continuum is
well above the magnon-pair resonance frequencyΔ_a=8J_11at the saturation.
It is also worth mentioning the field dependence of the intensities.
The intensity of magnon-pair resonance peak inI_MPR(ω)remains finite the
near the saturation field.
In contrast, the intensities of the unpaired magnon peak inI_k_M(ω)and the two-magnon continuumI_2-mag(ω), both of which originate from the unpaired magnon excitations, vanish near
saturation as
  
    𝒜_ k_M   ∝  H_c^ AFQ -  H  
    
       F(k)    ∝  H_c^ AFQ -  H

  forH<H_c^AFQ.
Hence the continuumI_2-mag(ω)disappears
around the saturation fieldH_c^AFQ.
Therefore, our method is more effective under the high magnetic field.


The peak intensity𝒜_MPRof the magnon-pair resonance has a strong field dependence,
showing
a divergence at a certain fieldH_∗belowH_c^AFQ[see Fig. <ref>(a)].
This divergence occurs when the magnon-pair resonance peak merges into the EPR one atH≃H_∗.
The divergence comes from the factorΔ_a/(Δ_a-H)^2in Eq. (<ref>);
as Eq. (<ref>) and Fig. <ref> show,
the excitation gapΔ_aof the bound magnon pair equals toHatH_∗= J_22/[J_22-J_11+J_11(H_c/8J_11)^2].
The intensity𝒜_k_M[Eq. (<ref>)] of the unpaired magnon resonance
also shows the divergence atH=0(Fig. <ref>),
because of merging of the peak into the EPR one atH≃0.


We note that the EPR and the MPR peaks can be mixed under certain interactions.
The EPR peak and the MPR peak are generated by application ofβ_k=0^†andα_k=0^†to a given eigenstate, respectively.
To mix those resonances, an anisotropic interaction including a term such asβ_kα_k'^†orβ_k^†α_k'is necessary.
For example, a uniform DM interaction  withDvector parallel to thexaxis can
generate effectively such an interaction,∑_⟨i,j|_⟩1 D ·S_i ×S_j=∑_⟨i,j|_⟩1 |D| [(b_i-b_i^†)(a_j+a_j^†) -(a_i+a_i^†) (b_j-b_j^†)].
If an anisotropic interaction allows the mixing,
it will be difficut to distinguish the MPR peak from the EPR one when their resonance frequencies are close,
because the EPR peak has a finite linewidth in the presence of anisotropic interactions.
When they are apart, the mixing is not important as long as the anisotropic interactins are perturbative.


Although the intensities suffer from the insignificant divergences,
the resonance frequencies are free from any singular behavior in the linear flavor-wave approximation
except for the singular bent point
at the saturation fieldH_c^AFQdue to physically reasonable characteristic upturn forH>H_c^AFQ(Fig. <ref>).
Using the magnon-pair resonance frequency and the unpaired-magnon resonance frequency atk_M,
we can identify the AFQ order phase in quite a wide field range.



Angular dependence of the intensity of the magnon-pair resonance (<ref>) enables another method of identification free from the technical problems.
The magnon-pair resonance peak shows a characteristic angular dependence

    𝒜_ MPR∝sin^2 θcos^2θ,

which makes the magnon-pair resonance peak vanish when the magnetic field is parallel to thex,y, orzaxes
[as shown in Fig. <ref>(b)].
The angular dependence (<ref>) holds true independent of the model and the theoretical technique,
as we pointed out in Sec. <ref>.
This angular dependence reflects the fact that the operatorsQ_j^yzandQ_j^zxneither create nor annihilate
the bound magnon pair, different fromQ_j^x^2-y^2andQ_j^xy.
Thus, the characteristic angular dependence of𝒜_MPR(<ref>) qualifies as an
evidence of the formation of the bound magnon pair
in the system (<ref>) with the single-ion anisotropy.


We note that the linewidth  of the EPR peak ofS=1/2frustrated ferromagnetic chain compounds is also expected
to show the angular dependence ofsin^2θcos^2θ<cit.>.
That angular dependence of the linewidth comes from the same root as the intensity of the magnon-pair resonance peak.




§ DISCUSSIONS


Here, we discuss some issues related to the magnon-pair resonance peak in the ESR spectrum shown in Sec. <ref>.
We also discuss applications of our theory toS=1/2spin systems.





 §.§ Effects of other anisotropic interactions


In Sec. <ref>, we assumed the single-ion anisotropy (<ref>) as an example of the perturbative anisotropic interaction.
Since ESR depends crucially on the form of the perturbative anisotropic interactionℋ',
it is necessary to confirm that our results obtained in Sec. <ref> is robust against inclusion of other kinds of anisotropic interactions.

The anisotropic exchange interaction (<ref>) on thenth neighbor bond leads to the same result
because the bound magnon pair is not localized at a single site but spread around bonds <cit.>.
If the𝒜operator (<ref>) contains some of operators that have the same symmetry as
the wavefunction of two-magnon bound state, it generates the magnon-pair resonance to the ESR spectrum through the formula (<ref>).
In contrast, the DM interactionℋ'_DM=∑_⟨i,j|_⟩nD_ij ·S_i ×S_jis irrelevant to the magnon-pair resonance (<ref>) and to the unpaired magnon resonance (<ref>) because of the symmetry;
the DM interaction neither create nor annihilate the bound magnon pair on the bond because
it is antisymmetric with respect to the bond-centered inversion, whereas
the wavefunction of the bound magnon pair on the bond is symmetric.




 §.§ Applications to the spin nematic order in S=1/2 spin systems

In Sec. <ref>, we studied theS=1spin model to demonstrate the magnon-pair resonance.
We can apply this result to the spin nematic order inS=1/2spin systems performing an appropriate mapping of low-energy degrees of freedoms.



  §.§.§ S=1/2 frustrated ferromagnets

In the case of spin-1/2 frustrated ferromagnets, the spin nematic order parameter, defined
on the nearest neighbor bonds,
hask=0wave vector, in which the sign of it alternates inside
the unit cell of the crystal structure <cit.>.
For example, on the square lattice, the two quadrupolar directors
on different bonds along two unit vectorse_1ande_2are orthogonal to each other <cit.>.
Because of this sign change, theS=1/2spin nematic order parameter is not captured
into the frequency shift discussed in Sec. <ref>, even though it hask=0wave vector.

Low-energy degrees of freedom in theS=1/2spin nematic systems are given byS=1spin degrees of freedom
formed on the nearest neighbor bonds <cit.>.
These excitation modes are effectively related to the excitations of theS=1AFQ state through a mapping
between bond degrees of freedom inS=1/2spin systems and on-site spin degrees of freedom inS=1spins <cit.>,
where theS=1spins are assigned on the center points of the nearest-neighbor bonds ofS=1/2spins.
On the square lattice, the gapless excitations withk=k_Min theS=1AFQ state correspond
to the excitations withk=0wave vector andB_1irreducible representation  of the space groupC_4vin
theS=1/2spin nematic states.
Above the saturationH>H_c^AFQ, this mode is the lowest excitation which closes the gap
at the saturation field.
However this mode is inaccessible in ESR measurements, since ESR is directly accessible only to
thek=0wave vector modes with theA_1(trivial) irreducible representation.
Only the gapped excitation modes withk=0can be observed
among the bound magnon pair excitations as same as in theS=1AFQ state discussed in Sec. <ref>.



  §.§.§ S=1/2 orthogonal dimer spin system


The spin nematic phase can also appear in spin-gapped systems when bound magon (triplon) pairs close the energy gap
in an applied field <cit.>.
In theS=1/2Heisenberg model on the Shastry-Sutherland
lattice <cit.>, which is also called an orthogonal dimer spin model,
it was theoretically demonstrated that 
the ground state is an exact dimer state with a finite energy gap <cit.>
and
bound two-triplon excited states <cit.> are
stabilized at zero field by the correlated hopping process <cit.>.
Theoretical calculations <cit.> pointed out that a two-triplon bound state withS^z=2can have a lower energy than two triplon continuum above the gapped ground state and
that the energy-gap closing in an applied magnetic field leads to
the condensation of bound triplon pairs.
Since the lowest energy state of the bound pair in theS^z=2sector has the wave vectork=k_M,
the field-induced condensed phase becomes an AFQ phase [ 
The field-induced antiferroquadrupolar (AFQ) phase in the spin-gapped systems is also described as a low-density condensate of bound magnon pairs, as same as the AFQ phase near the saturation field. The vacuum is a singlet ground state and the bosonic particles are bound triplon pairs in the former case, whereas the vacuum is the fully polarized state and the bosons are bound magnon pairs in the latter. ].

In this system, anisotropic interactions between two orthogonal dimers cause the operatorS_i^+ S_j^+on the inter-dimer
bonds in the𝒜operator. This operator creates a triplon pair on a nearest-neighbor pair of dimers,
which gives rise to a triplon-pair resonance peak in the ESR spectrum.
In the spin gap phase, the resonance frequency behaves as

    ω=Δ_a+2(H_c-H),

whereH_cdenotes the onset-field of the magnetization process andΔ_athe energy gap of the bound triplon
pair atk=0at the critical fieldH=H_c. We note that the bound pair
closes the gap  atk=(π,π)and this excitation is well dispersive, i.e.Δ_a>0<cit.>.
In the magnetic phase aboveH_c, we expect that this peak continuously connects to the triplon pair
resonance peak in the AFQ phase showing a singularity in the field dependence of the frequency at the critical field.

In an ESR study <cit.>, bound triplon-pair resonance peaks were indeed observed in theS=1/2orthogonal dimer spin compound SrCu_2(BO_3)_2.
The lowest-energy resonance peak of the bound triplon pairs shows the field dependence (<ref>), having a strong intensity aroundH=H_c.
Even after the peak frequency changes the slope as a function of a field aroundH=H_c,
the resonance peak remains with strong intensity
inside the magnetic phase between the spin gapped and the 1/8-plateau phases when
the magnetic field is parallel toaaxis.
The implication of this resonance peak inside the magnetic phase has not been properly
understood until now.
Our research elucidates that this ESR result has already suggested the appearance
of a spin nematic order in the ground state
of the field-induced magnetic phase below the 1/8-plateau.
This system deserves further investigations.

To compare the field dependence of the resonance peak with observed results in real compounds, we need
to include mixing between the ground state and excited states. For example,
DM interaction induces a mixing between the singlet ground state and triplon excited states <cit.>.
In the case of the bound two-triplon excited state, anisotropic interactions can induce mixing with the singlet
ground state.
This can be easily seen by considering an anisotropy on the inter-dimer bonds

    2δ (S_i^x S_j^x - S_i^y S_j^y) = δ (S_i^+ S_j^+ + S_i^- S_j^-),

which mixes the bound two-triplon state with the singlet dimer state.
This effect might smear the singularity in the field dependence of the resonance frequency atH=H_c.



 §.§ Magnon-pair resonance in the case of ferroquadrupolar order

Lastly we comment on the additional peaks in the ESR spectrum in the FQ phase.

In the FQ phase, only the EPR peak will be found in the ESR spectrum.
The bound magnon-pair excitation in the FQ phase is gapless atk=0, i.e.Δ_a=0,
but the resonance atω=0is invisible in the ESR spectrum for the factorωin Eq. (<ref>). 
Note that the unpaired magnon resonance peak atk=k_Min the AFQ phase corresponds to the
EPR peak in the FQ phase since the FQ order grows atk = 0.

If the magnetic field is above the saturation fieldH_c^FQ, the bound magnon pair excitation opens a gap,
showing a
characteristic field dependenceΔ_a=2(H-H_c^FQ). This peak is observable in the ESR measurements
as it comes from thek=0modes.
This gives a clear difference from the case of the AFQ order;
the lowest excitation which closes the gap as2(H-H_c^AFQ)above the AFQ phase has
thek=k_Mwave vector [Eq. (<ref>)] and it cannot be observed in ESR measurements.
Thus, the appearance of this peak in the ESR spectrum above the saturation fieldH_c^FQand
the disappearance belowH_c^FQsignal the emergence of the FQ phase belowH_c^FQ.




§ SUMMARY


In this paper we showed that the FQ and the AFQ orders are distinguishable in ESR experiments.
We studied both the frequency shift of EPR resonance and the frequencies of the additional resonance peaks
in the ESR spectrum.

For the generic spin model (<ref>),
the FQ order parameter turned out to shift the resonance frequency of the EPR peak in the ESR spectrum.
The EPR frequency shift shows a characteristic angular dependence [as shown in Eq. (<ref>)]
on a rotation of the material around theyaxis keeping
the magnetic field parallel to thezaxis.
Here we determined theyandzaxes so that the anisotropic spin interactions in these spin components
are, respectively, weakest and strongest.
Thus the FQ order parameter can be extracted from the frequency shift.
For example, as Eqs. (<ref>) and (<ref>) show, the frequency shifts
atθ=0andπ/2enable us to
determine the FQ order parameter⟨S_i^-S_j^-|_⟩0experimentally because only two quantities⟨3S_i^zS_j^z-S_i ·S_j|_⟩0and⟨S_i^-S_j^-|_⟩0are the undetermined variables in these equations.
In particular, when the perturbative anisotropic interaction is uniaxial, the FQ order parameter is simply derived
from the single equation
(<ref>).
The unexplained high field phase of chromium spinels is an interesting target to which this method is applicable.

In the case of the AFQ order, though
the EPR frequency shift is usually insensitive to the order parameter,
fingerprints of the AFQ order appear in the additional resonance peaks other than the EPR peak in the
ESR spectrum.
As far as the resonance peaks well isolated from the EPR one are concerned,
the ESR spectrum, as shown in Eq. (<ref>), is derived from the spectrum of the retarded Green
function of the operator𝒜=[ℋ',S^+]given by the small anisotropic interactionℋ'.
The operator𝒜is usually quadratic containing the magnon pair creation operator
[see Eqs. (<ref>) and (<ref>)].
This is one of the most interesting properties of the ESR spectrum.
Except for the vicinity of the EPR peak atω=H,
measuring the ESR spectrum is effectively equivalent to measuring the spectrum of the  operator𝒜.
We note that, in our pertubative analysis, the anisotropic interaction plays no role of yielding
the spin nematic phase and of making the magnitude of the quadrupolar
order parameter grow. Those are fully determined in the unperturbed system.


The long-range AFQ order yields two additional sharp resonance peaks in the ESR spectrum.
One is attributed to the resonance of the bound magnon-pair excitation,
which we called the magnon-pair resonance.
The magnon-pair resonance is also found in the fully polarized phase adjacent to the AFQ phase,
where the resonance frequency shows the linear field dependence whose slope is double of
that of the EPR frequency, as was experimentally observed in Ref. <cit.>.
(A similar triplon-pair resonance peak was
also experimentally observed in Ref. <cit.>.)
With decreasing the magnetic field, the system enters into the AFQ phase,
where the magnon-pair resonance frequency shows the singular upturn as a function of the magnetic field,
reflecting the condensation of bound magnon pairs (Fig. <ref>).
The other resonance peak is attributed to the excitation of the unpaired magnon at the wave vectork_M=(π,π).
Usually, ESR detects excitations at the wave vectork=0.
In the AFQ phase, since the ground state structure has
the wave vectork_M,
the excitation gap of the magnon atk_Mbecomes
visible in the ESR spectrum as an independent resonance peak.


Our results on the FQ and the AFQ orders are valid as long as (1) the anisotropic interaction is small enough to be seen as a perturbation
to the system
and
(2) the anisotropic interaction is governed mainly by the single-ion anisotropy or the anisotropic exchange interaction.
The weak DM interaction has no impact on the result obtained in this paper because the DM interaction is antisymmetric
with respect to the bond-centered inversion unlike the spin nematic order parameter.
Though we restricted ourselves to the cases of weak anisotropic interactions in this paper,
it is also interesting to investigate cases governed by a large anisotropic interaction such as the case of Ref. <cit.>.
While the formula of the frequency shift [Eq. (<ref>)] is invalid in such cases,
the discussion of the sharp magnon-pair resonance peak isolated from the EPR one will be qualitatively valid even in the case of large anisotropic interactions though we need to derive the full Green's function𝒢^R_𝒜𝒜^†(ω)instead.







§ ACKNOWLEDGMENTS


We thank Akira Furusaki, Masayuki Hagiwara, Hiroyuki Nojiri, Nic Shannon, and Shintaro Takayoshi for helpful discussions.
The present work is supported by JSPS KAKENHI Grant Nos. 16J04731 and 16K05425.



§ POLARIZATION INDEPENDENCE OF THE ESR SPECTRUM

In this Appendix, we derive the ESR spectrum of the unpolarized microwave.
To discuss it, we review the linear response theory of the ESR absorption spectrumI(ω)in a generic spin system described by

    ℋ̃(t) = ℋ - ∫_-π^π dθ∫_-π^π dϕ  A(θ) X(t,θ, ϕ),

whereℋis the Hamiltonian of the spin system of our interest,A(θ)is a total spin operator

    A(θ) = S^xcosθ + S^y sinθ,

andX(t,θ, ϕ)is the oscillating magnetic field,

    X(t,θ, ϕ) = h_R(θ, ϕ) cos(ω t + ϕ).

The amplitudeh_R(θ, ϕ)follows a distribution which can be random or not.
For example, the Hamiltonian under a circularly polarized magnetic field

    ℋ̃(t) = ℋ - H_R/2 (S^+e^iω t + S^- e^-iω t)

is a special case of Eq. (<ref>) with

    h_R(θ, ϕ) = H_R/2[ δ(θ)δ(ϕ) + δ(θ - π2)δ(ϕ - π2)].


The ESR absorption spectrumI(ω)is give by
the energy absorption rate per a period of the oscillating field,

    I(ω) = ω/2π∫_0^2π/ω dt  d/dt[ρ(t)ℋ̃(t)],

whereρ(t)is the density matrix in the Heisenberg picture,

    ρ(t) =  U(t)exp(-ℋ̃(t)/ T)/[exp(-ℋ̃(t)/T)] U^†(t)

with

    U(t) = exp(i∫_0^t dt' ℋ̃(t')).

Ifh_R(θ, ϕ)follows a random distribution, we replace Eq. (<ref>) with

    I(ω) = ω/2π∫_0^2π/ω dt  d/dt[ρ(t)ℋ̃(t)],

whereOdenotes the average of the quantityOwith respect to the random distribution.


Let us first consider the case thath_R(θ, ϕ)is uniquely determined such as Eq. (<ref>).
The energy absorption rate (<ref>) is written as

    I(ω)    =  - ω/2π∫_0^2π/ω dt  ∫_-π^π dθ dϕ [ρ(t) A(θ)] ∂ X(t,θ, ϕ)/∂ t.

Within the linear response, the trace [ρ(t)A(θ)] is approximated as

    [ρ(t)A(θ)]
        ≃⟨A(θ)|+⟩i∫_0^∞ dt' ∫_-π^π dθ' dϕ'  
     X(t-t',θ', ϕ')⟨[A(t',θ), A(0,θ')]|.⟩

Here, the average ⟨·|$⟩ is taken with respect to the HamiltonianℋandA(t, θ)is defined as

    A(t', θ) = e^it'ℋ A(θ) e^-it'ℋ.

Taking these relations into account, we find that the ESR energy absorption rate (<ref>) is given by

    I(ω)
        ≃ -iω/2π∫_0^2π/ω dt  ∫_0^∞ dt'
     ∫_-π^π dθ dϕ dθ' dϕ'  ∂ X(t,θ, ϕ)/∂ tX(t-t',θ',ϕ') ⟨[A(t',θ), A(0,θ')]|⟩
       =  iω/8∫_0^∞ dt' ∫_-π^π dθ dϕ dθ' dϕ'   h_R(θ, ϕ)h_R(θ', ϕ')
     sin(ω t' + ϕ - ϕ')
     {⟨[S^+(t'), S^-(0)]|e⟩^-i(θ - θ')
                       
    + ⟨[S^-(t'), S^+(0)]|e⟩^i(θ-θ')
     + ⟨[S^+(t'), S^+(0)]|e⟩^-i(θ + θ') + ⟨[S^-(t'), S^-(0)]|e⟩^i(θ + θ')}.

Whenh_R(θ, ϕ)is given by Eq. (<ref>), the ESR absorption rate (<ref>) becomes

    I(ω)
        = ω H_R^2/8[ -𝒢^R_S^+S^-(ω)],

which reproduces Eq. (<ref>).

We next consider the case thath_R(θ, ϕ)follows a random distribution.
Applying the linear response theory (<ref>) to the random averaged energy absorption rate
(<ref>), we find
 
    I(ω)
         =  - ω/2π∫_0^2π/ω dt  ∫_-π^π dθ dϕ [ρ(t) A(θ)] ∂ X(t,θ, ϕ)/∂ t
       ≃iω/8∫_0^∞ dt' ∫_-π^π dθ dϕ dθ' dϕ'  h_R(θ, ϕ)h_R(θ', ϕ')sin(ω t' + ϕ - ϕ')
     {⟨[S^+(t'), S^-(0)]|e⟩^-i(θ - θ')
                       
    + ⟨[S^-(t'), S^+(0)]|e⟩^i(θ-θ')
     + ⟨[S^+(t'), S^+(0)]|e⟩^-i(θ + θ') + ⟨[S^-(t'), S^-(0)]|e⟩^i(θ + θ')}.

If the applied electromagnetic wave is “white”, that is, if
the random averageh_R(θ, ϕ)h_R(θ', ϕ')satisfies

    h_R(θ, ϕ)h_R(θ', ϕ') = H_R^2/(2π)^2δ(θ-θ')δ(ϕ-ϕ'),

the ESR energy absorption rate (<ref>) is simplified as

    I(ω)
        = ω H_R^2/8[ -𝒢^R_S^+S^-(ω) - 𝒢^R_S^-S^+(ω)].


§ PERTURBATIONS OF THE GROUND STATE ENERGY IN THE AFQ PHASE


This appendix is devoted to estimation of the ground state energy in the AFQ phase of theS=1model on the square lattice
described by the Hamiltonian

    ℋ = ℋ_0 + ℋ',

whereℋ_0is the unperturbed Hamiltonian (<ref>).
We take the single-ion anisotropy (<ref>)
as the perturbationℋ'.



As we discussed in Sec. <ref>,  the angleφthat specifies the direction of the AFQ order growing is not determined at the mean-field level.
Here, we estimate theφdependence of the ground-state energy using the linear flavor-wave theory explained in Sec. <ref>.
The perturbative single-ion anisotropy shifts the ground state energy from its unperturbed value by an amount

    δ E_0 = ⟨ GS|ℋ'| GS|,⟩

up to the first order of the perturbation.
Here,| GS⟩is the ground state in the AFQ phase of the unperturbed Hamiltonian, or the vacuum annihilated
byα_ kandβ_ kfor all k:

    α_ k| GS⟩ = β_ k| GS⟩ = 0.




The energy shift (<ref>) in the Schwinger boson language is given by

    δ E_0
        = const. - ( D - 3(D-E)/2sin^2θ) ∑_ k⟨ GS|b_ k^† b_ k| GS|⟩
          +1/2{ Dsin^2θ + E(1+cos^2θ)}cos 2φsin 2θ_H ∑_ k⟨ GS|(2a_ k^† a_ k+ k_M + b_ k^† b_ k+  k_M)| GS|.⟩

Using the Bogoliubov transformations (<ref>) and (<ref>) and
the property (<ref>) in the ground state,
we can further reduceδ E_0to

    δ E_0
        = const. - ( D - 3(D-E)/2sin^2θ) ∑_ ksinh^2Θ_ k^b.

The shift (<ref>) is independent of the angleφand hence leavesφundetermined.


73
fxundefined [1]
 ifx#1
fnum [1]
  #1firstoftwo
 secondoftwo
 
fx [1]
  #1firstoftwo
 secondoftwo
 





noop [0]secondoftwo

ref[1]@startlink#1@href
href[1]#1@endlink
anitize@url [0] `
12 `$12
  `&12 `#12 `1̂2 `_12 `%12
startlink[1]
endlink[0]

rl [1]href #1












@bib@innerbibempty

 [Blume and Hsieh(1969)]blume_hsieh
  authorauthorM. Blume and authorY. Y. Hsieh,  http://dx.doi.org/10.1063/1.1657616journaljournalJ. Appl. Phys. volume40, pages1249 (year1969)NoStop
 [Chen and Levy(1971)]chen_leby
  authorauthorH. H. Chen and authorP. M. Levy,  10.1103/PhysRevLett.27.1383journaljournalPhys. Rev. Lett. volume27, pages1383 (year1971)NoStop
 [Andreev and Grishchuk()]andreev_grishchuk
  authorauthorA. Andreev and authorI. Grishchuk, @noop journalZh. Eksp. Teor. Fiz. 87, 467 (1984) [Sov. Phys. JETP  60, 267 (1984)] NoStop
 [Chubukov(1991)]chubukov_1dnematic
  journalauthorauthorA. V. Chubukov, 
  10.1103/PhysRevB.44.4693journaljournalPhys.
  Rev. B volume44, pages4693
  (year1991)NoStop
 [Shannon et al.(2006)Shannon, Momoi, and Sindzingre]shannon_nematic_sq
  authorauthorN. Shannon, authorT. Momoi,   and authorP. Sindzingre, 10.1103/PhysRevLett.96.027213journaljournalPhys. Rev. Lett. volume96, pages027213 (year2006)NoStop
 [Momoi et al.(2006)Momoi,
  Sindzingre, and Shannon]momoi2006
  authorauthorT. Momoi, authorP. Sindzingre,   and authorN. Shannon,  10.1103/PhysRevLett.97.257204journaljournalPhys. Rev. Lett. volume97, pages257204 (year2006)NoStop
 [Ueda and Totsuka(2007)]ueda2007
  authorauthorH. T. Ueda and authorK. Totsuka, 10.1103/PhysRevB.76.214428journaljournalPhys. Rev. B volume76, pages214428 (year2007)NoStop
 [Shindou and Momoi(2009)]shindou2009
  authorauthorR. Shindou and authorT. Momoi,  10.1103/PhysRevB.80.064410journaljournalPhys. Rev. B volume80,pages064410 (year2009)NoStop
 [Shindou et al.(2011)Shindou, Yunoki, and Momoi]shindou2011
  authorauthorR. Shindou, authorS. Yunoki,   and authorT. Momoi,  10.1103/PhysRevB.84.134414journaljournalPhys. Rev. B volume84, pages134414 (year2011)NoStop
 [Heidrich-Meisner et al.(2006)Heidrich-Meisner, Honecker, and Vekua]heidrichmeisner2006
  authorauthorF. Heidrich-Meisner, authorA. Honecker,  and authorT. Vekua,  10.1103/PhysRevB.74.020403journaljournalPhys. Rev. B volume74,pages020403 (year2006)NoStop
 [Vekua et al.(2007)Vekua,
  Honecker, Mikeska, and Heidrich-Meisner]vekua_1dnematic
  authorauthorT. Vekua, authorA. Honecker,
  authorH.-J. Mikeska,  andauthorF. Heidrich-Meisner, 10.1103/PhysRevB.76.174420journaljournalPhys. Rev. B volume76, pages174420 (year2007)NoStop
 [Kecke et al.(2007)Kecke,
  Momoi, and Furusaki]kecke_1dnematic
  authorauthorL. Kecke, authorT. Momoi,  andauthorA. Furusaki,  10.1103/PhysRevB.76.060407journaljournalPhys. Rev. B volume76, pages060407 (year2007)NoStop
 [Hikihara et al.(2008)Hikihara, Kecke, Momoi, and Furusaki]hikihara_1dnematic
  authorauthorT. Hikihara, authorL. Kecke,
  authorT. Momoi,  and authorA. Furusaki, 
  10.1103/PhysRevB.78.144404journaljournalPhys.
  Rev. B volume78, pages144404
  (year2008)NoStop
 [Sudan et al.(2009)Sudan,
  Lüscher, and Läuchli]sudan_1dnematic
  authorauthorJ. Sudan, authorA. Lüscher,   and authorA. M. Läuchli, 10.1103/PhysRevB.80.140402journaljournalPhys. Rev. B volume80, pages140402 (year2009)NoStop
 [Ueda and Totsuka(2009)]ueda2009
  authorauthorH. T. Ueda and authorK. Totsuka, 10.1103/PhysRevB.80.014417journaljournalPhys. Rev. B volume80, pages014417 (year2009)NoStop
 [Zhitomirsky and Tsunetsugu(2010)]zhitomirsky_j1j2
  authorauthorM. E. Zhitomirsky and authorH. Tsunetsugu, http://stacks.iop.org/0295-5075/92/i=3/a=37001journaljournalEuro. Phys. Lett. volume92, pages37001 (year2010)NoStop
 [Sato et al.(2013)Sato,
  Hikihara, and Momoi]sato_q1dnematic
  authorauthorM. Sato, authorT. Hikihara,   and authorT. Momoi,  10.1103/PhysRevLett.110.077206journaljournalPhys. Rev. Lett. volume110, pages077206 (year2013)NoStop
 [Starykh and Balents(2014)]starykh_q1dnematic
  authorauthorO. A. Starykh and authorL. Balents,  10.1103/PhysRevB.89.104407journaljournalPhys. Rev. B volume89, pages104407 (year2014)NoStop
 [Momoi et al.(2012)Momoi,
  Sindzingre, and Kubo]momoi2012
  authorauthorT. Momoi, authorP. Sindzingre,   and authorK. Kubo,  10.1103/PhysRevLett.108.057206journaljournalPhys. Rev. Lett. volume108, pages057206 (year2012)NoStop
 [Ueda and Momoi(2013)]UedaMomoi2013
  authorauthorH. T. Ueda and authorT. Momoi, 10.1103/PhysRevB.87.144417journaljournalPhys. Rev. B volume87, pages144417 (year2013)NoStop
 [Janson et al.(2016)Janson,
  Furukawa, Momoi, Sindzingre,
  Richter, and Held]janson2016
  authorauthorO. Janson, authorS. Furukawa,
  authorT. Momoi, authorP. Sindzingre, authorJ. Richter,  and authorK. Held, 
  10.1103/PhysRevLett.117.037206journaljournalPhys. Rev. Lett. volume117, pages037206 (year2016)NoStop
 [Nawa et al.(2014)Nawa,
  Okamoto, Matsuo, Kindo,
  Kitahara, Yoshida, Ikeda,
  Hara, Sakurai, Okubo,
  Ohta, and Hiroi]nawa_nmr_exp
  authorauthorK. Nawa, authorY. Okamoto,
  authorA. Matsuo, authorK. Kindo, authorY. Kitahara, authorS. Yoshida, authorS. Ikeda, authorS. Hara, authorT. Sakurai, authorS. Okubo,
  authorH. Ohta,  and authorZ. Hiroi, 
  10.7566/JPSJ.83.103702journaljournalJ. Phys.
  Soc. Jpn. volume83, pages103702
  (year2014)NoStop
 [Büttgen et al.(2014)Büttgen, Nawa, Fujita, Hagiwara, Kuhns, Prokofiev, Reyes, Svistov, Yoshimura, andTakigawa]buttgen_nmr_exp
  authorauthorN. Büttgen, authorK. Nawa,
  authorT. Fujita, authorM. Hagiwara, authorP. Kuhns, authorA. Prokofiev, authorA. P.Reyes, authorL. E. Svistov, authorK. Yoshimura,
   and authorM. Takigawa, 10.1103/PhysRevB.90.134401journaljournalPhys. Rev. B volume90, pages134401 (year2014)NoStop
 [Nawa et al.(2017)Nawa,
  Takigawa, Krämer, Horvati ćć, Berthier,
  Yoshida, and Yoshimura]nawa2017
  authorauthorK. Nawa, authorM. Takigawa,
  authorS. Krämer, authorM. Horvati ćć, authorC. Berthier,
  authorM. Yoshida,  andauthorK. Yoshimura,  10.1103/PhysRevB.96.134423journaljournalPhys. Rev. B volume96, pages134423 (year2017)NoStop
 [Grafe et al.(2017)Grafe,
  Nishimoto, Iakovleva, Vavilova, Alfonsov, Sturza, Wurmehl, Nojiri, Rosner, Richteret al.]grafe_nmr_exp
  authorauthorH.-J. Grafe, authorS. Nishimoto,
  authorM. Iakovleva, authorE. Vavilova, authorA. Alfonsov, authorM.-I. Sturza, authorS. Wurmehl, authorH. Nojiri, authorH. Rosner, authorJ. Richter,  et al., 
  10.1038/s41598-017-06525-0journaljournalSci.
  Rep. volume7, pages6720 (year2017)NoStop
 [Orlova et al.(2017)Orlova,
  Green, Law, Gorbunov,
  Chanda, Krämer, Horvati ćć, Kremer,
  Wosnitza, and Rikken]orlova
  authorauthorA. Orlova, authorE. L. Green,
  authorJ. M. Law, authorD. I. Gorbunov, authorG. Chanda, authorS. Krämer, authorM. Horvati ćć, authorR. K. Kremer, authorJ. Wosnitza,  and authorG. L. J. A. Rikken, 
  10.1103/PhysRevLett.118.247201journaljournalPhys. Rev. Lett. volume118, pages247201 (year2017)NoStop
 [Yoshida et al.(2017)Yoshida, Nawa, Ishikawa, Takigawa, Jeong, Krämer, Horvati ćć, Berthier,
  Matsui, Goto, Kimura,
  Sasaki, Yamaura, Yoshida,
  Okamoto, and Hiroi]yoshida2016
  authorauthorM. Yoshida, authorK. Nawa,
  authorH. Ishikawa, authorM. Takigawa, authorM. Jeong, authorS. Krämer, authorM. Horvati ćć, authorC. Berthier, authorK. Matsui, authorT. Goto, authorS. Kimura, authorT. Sasaki,
  authorJ. Yamaura, authorH. Yoshida, authorY. Okamoto,  and authorZ. Hiroi, 
  10.1103/PhysRevB.96.180413journaljournalPhys.
  Rev. B volume96, pages180413
  (year2017)NoStop
 [Momoi and Totsuka(2000a)]momoi2000b
  authorauthorT. Momoi and authorK. Totsuka,  10.1103/PhysRevB.62.15067journaljournalPhys. Rev. B volume62, pages15067 (year2000a)NoStop
 [Nojiri et al.(2003)Nojiri,
  Kageyama, Ueda, and Motokawa]nojiri2003
  authorauthorH. Nojiri, authorH. Kageyama,
  authorY. Ueda,  and authorM. Motokawa, 
  10.1143/JPSJ.72.3243journaljournalJ. Phys.  Soc. Jpn. volume72, pages3243
  (year2003)NoStop
 [Papanicolaou(1988)]papanicolaou2
  authorauthorN. Papanicolaou, 
  http://dx.doi.org/10.1016/0550-3213(88)90073-9journaljournalNucl. Phys. B volume305,pages367  (year1988)NoStop
 [Tanaka et al.(2001)Tanaka,
  Tanaka, and Idogaki]tanaka2001
  authorauthorK. Tanaka, authorA. Tanaka,   and authorT. Idogaki, http://stacks.iop.org/0305-4470/34/i=42/a=304journaljournalJ. Phys. A: Math. Gen. volume34, pages8767 (year2001)NoStop
 [Harada and Kawashima(2002)]harada
  authorauthorK. Harada and authorN. Kawashima,  10.1103/PhysRevB.65.052403journaljournalPhys. Rev. B volume65, pages052403 (year2002)NoStop
 [Tsunetsugu and Arikawa(2006)]tsunetsugu_nematic_tri
  authorauthorH. Tsunetsugu and authorM. Arikawa,  10.1143/JPSJ.75.083701journaljournalJ. Phys. Soc. Jpn. volume75, pages083701 (year2006)NoStop
 [Läuchli et al.(2006)Läuchli, Mila, and Penc]lauchli_blbq_tri
  authorauthorA. Läuchli, authorF. Mila,   and authorK. Penc,  10.1103/PhysRevLett.97.087205journaljournalPhys. Rev. Lett. volume97, pages087205 (year2006)NoStop
 [Takata et al.()Takata,
  Momoi, and Oshikawa]takata_pyrochlore
  authorauthorE. Takata, authorT. Momoi,   and authorM. Oshikawa,@noop  http://arxiv.org/abs/arXiv:1510.02373arXiv:1510.02373NoStop
 [Wang et al.(2015)Wang,
  Kivelson, and Lee]wang_fese_nat
  authorauthorF. Wang, authorS. A. Kivelson,
   and authorD.-H. Lee,  10.1038/nphys3456journaljournalNat. Phys. volume11, pages959
  (year2015)NoStop
 [Yu and Si(2015)]yu_fese_afq
  authorauthorR. Yu and authorQ. Si,  10.1103/PhysRevLett.115.116401journaljournalPhys. Rev. Lett. volume115, pages116401 (year2015)NoStop
 [Wang et al.(2016)Wang,
  Hu, and Nevidomskyy]wang_fese_fq
  authorauthorZ. Wang, authorW.-J. Hu,  andauthorA. H. Nevidomskyy, 10.1103/PhysRevLett.116.247203journaljournalPhys. Rev. Lett. volume116, pages247203 (year2016)NoStop
 [Sato et al.(2009)Sato,
  Momoi, and Furusaki]sato_nmr_1dnematic_1
  authorauthorM. Sato, authorT. Momoi,  andauthorA. Furusaki,  10.1103/PhysRevB.79.060406journaljournalPhys. Rev. B volume79, pages060406 (year2009)NoStop
 [Sato et al.(2011)Sato,
  Hikihara, and Momoi]sato_nmr_1dnematic_2
  authorauthorM. Sato, authorT. Hikihara,   and authorT. Momoi,  10.1103/PhysRevB.83.064405journaljournalPhys. Rev. B volume83, pages064405 (year2011)NoStop
 [Smerald and Shannon(2016)]smerald_nmr_nematic
  authorauthorA. Smerald and authorN. Shannon,  10.1103/PhysRevB.93.184419journaljournalPhys. Rev. B volume93, pages184419 (year2016)NoStop
 [Shindou et al.(2013)Shindou, Yunoki, and Momoi]shindou2013
  authorauthorR. Shindou, authorS. Yunoki,   and authorT. Momoi,  10.1103/PhysRevB.87.054429journaljournalPhys. Rev. B volume87, pages054429 (year2013)NoStop
 [Smerald et al.(2015)Smerald, Ueda, and Shannon]smerald_nematic
  authorauthorA. Smerald, authorH. T. Ueda,   and authorN. Shannon,  10.1103/PhysRevB.91.174402journaljournalPhys. Rev. B volume91, pages174402 (year2015)NoStop
 [Onishi(2015)]onishi_1dnematic
  authorauthorH. Onishi,  10.7566/JPSJ.84.083702journaljournalJ. Phys. Soc. Jpn. volume84, pages083702 (year2015)NoStop
 [Michaud et al.(2011)Michaud, Vernay, and Mila]mila2011
  authorauthorF. Michaud, authorF. Vernay,   and authorF. Mila,  10.1103/PhysRevB.84.184424journaljournalPhys. Rev. B volume84, pages184424 (year2011)NoStop
 [Savary and Senthil()]savary_rixs
  authorauthorL. Savary and authorT. Senthil, @noop  http://arxiv.org/abs/arXiv:1506.04752arXiv:1506.04752NoStop
 [Furuya(2017)]furuya_esr_1dnematic
  authorauthorS. C. Furuya,  10.1103/PhysRevB.95.014416journaljournalPhys. Rev. B volume95, pages014416 (year2017)NoStop
 [Akaki et al.(2017)Akaki,
  Yoshizawa, Okutani, Kida,
  Romhányi, Penc, and Hagiwara]akaki_esr_quad
  authorauthorM. Akaki, authorD. Yoshizawa,
  authorA. Okutani, authorT. Kida, authorJ. Romhányi, authorK. Penc,  and authorM. Hagiwara, 
  10.1103/PhysRevB.96.214406journaljournalPhys.
  Rev. B volume96, pages214406
  (year2017)NoStop
 [Kubo and Tomita(1954)]kubo_tomita
  authorauthorR. Kubo and authorK. Tomita, 10.1143/JPSJ.9.888journaljournalJ. Phys. Soc. Jpn. volume9, pages888 (year1954)NoStop
 [Oshikawa and Affleck(2002)]oshikawa_esr
  authorauthorM. Oshikawa and authorI. Affleck,  10.1103/PhysRevB.65.134410journaljournalPhys. Rev. B volume65, pages134410 (year2002)NoStop
 [Ozerov et al.(2015)Ozerov,
  Maksymenko, Wosnitza, Honecker, Landee, Turnbull, Furuya, Giamarchi, and Zvyagin]ozerov_esr_dimpy
  authorauthorM. Ozerov, authorM. Maksymenko,
  authorJ. Wosnitza, authorA. Honecker, authorC. P. Landee, authorM. M. Turnbull, authorS. C. Furuya, authorT. Giamarchi,  and authorS. A. Zvyagin, 
  10.1103/PhysRevB.92.241113journaljournalPhys.
  Rev. B volume92, pages241113
  (year2015)NoStop
 [Furuya and Oshikawa(2012)]furuya_esr_boundary
  authorauthorS. C. Furuya and authorM. Oshikawa,  10.1103/PhysRevLett.109.247603journaljournalPhys. Rev. Lett. volume109, pages247603 (year2012)NoStop
 [Note1()]Note1
  noteWe note that, when a certain symmetry is broken
  spontaneously in the unperturbed system and it is also broken by the
  anisotropy in the full Hamiltonian ℋ, one must choose a
  proper direction of the infinitesimal auxiliary symmetry breaking field in
  the unperturbed system so that the ground state is smoothly deformed by
  imposing weak perturbation.Stop
 [Maeda and Oshikawa(2005)]maeda_shift_jpsj
  authorauthorY. Maeda and authorM. Oshikawa,  10.1143/JPSJ.74.283journaljournalJ. Phys. Soc. Jpn. volume74, pages283 (year2005)NoStop
 [Povarov et al.(2011)Povarov, Smirnov, Starykh, Petrov, and Shapiro]povarov_1d_dm
  authorauthorK. Y. Povarov, authorA. I. Smirnov, authorO. A. Starykh, authorS. V. Petrov,
   and authorA. Y. Shapiro, 10.1103/PhysRevLett.107.037204journaljournalPhys. Rev. Lett. volume107,pages037204 (year2011)NoStop
 [Furuya et al.(2011)Furuya,
  Suzuki, Takayoshi, Maeda,  and Oshikawa]furuya_esr_haldane
  authorauthorS. C. Furuya, authorT. Suzuki,
  authorS. Takayoshi, authorY. Maeda,  and authorM. Oshikawa, 
  10.1103/PhysRevB.84.180410journaljournalPhys.
  Rev. B volume84, pages180410
  (year2011)NoStop
 [Maeda et al.(2005)Maeda,
  Sakai, and Oshikawa]maeda_esr_shift
  authorauthorY. Maeda, authorK. Sakai,  andauthorM. Oshikawa,  10.1103/PhysRevLett.95.037602journaljournalPhys. Rev. Lett. volume95, pages037602 (year2005)NoStop
 [Furuya et al.(2012)Furuya,
  Bouillot, Kollath, Oshikawa,  and Giamarchi]furuya_esr_bpcb
  authorauthorS. C. Furuya, authorP. Bouillot,
  authorC. Kollath, authorM. Oshikawa,  and authorT. Giamarchi, 
  10.1103/PhysRevLett.108.037204journaljournalPhys. Rev. Lett. volume108, pages037204 (year2012)NoStop
 [Mori and Kawasaki(1962)]mk
  authorauthorH. Mori and authorK. Kawasaki,  10.1143/PTP.27.529journaljournalProg. Theor. Phys. volume27, pages529 (year1962)NoStop
 [Penc et al.(2004)Penc,
  Shannon, and Shiba]penc_pyrochlore
  authorauthorK. Penc, authorN. Shannon,   and authorH. Shiba,  10.1103/PhysRevLett.93.197203journaljournalPhys. Rev. Lett. volume93, pages197203 (year2004)NoStop
 [Kimura et al.(2011)Kimura,
  Hagiwara, Takeuchi, Yamaguchi, Ueda, Ueda, and Kindo]Hg
  authorauthorS. Kimura, authorM. Hagiwara,
  authorT. Takeuchi, authorH. Yamaguchi, authorH. Ueda, authorY. Ueda,  and authorK. Kindo,  10.1103/PhysRevB.83.214401journaljournalPhys. Rev. B volume83, pages214401 (year2011)NoStop
 [Miyata et al.()Miyata,
  Ueda, and Takeyama]Cd
  authorauthorA. Miyata, authorH. Ueda,  andauthorS. Takeyama, @noop
  http://arxiv.org/abs/arXiv:1302.3664arXiv:1302.3664NoStop
 [Miyata et al.(2012)Miyata,
  Ueda, Ueda, Motome,
  Shannon, Penc, and Takeyama]Zn
  authorauthorA. Miyata, authorH. Ueda,
  authorY. Ueda, authorY. Motome, authorN. Shannon, authorK. Penc,  and authorS. Takeyama,  10.1143/JPSJ.81.114701journaljournalJ. Phys. Soc. Jpn. volume81, pages114701 (year2012)NoStop
 [Smerald and Shannon(2013)]smerald_nematic_tri
  authorauthorA. Smerald and authorN. Shannon,  10.1103/PhysRevB.88.184430journaljournalPhys. Rev. B volume88, pages184430 (year2013)NoStop
 [Hohenberg and Brinkman(1974)]Hohenberg1974
  authorauthorP. C. Hohenberg and authorW. F. Brinkman,  10.1103/PhysRevB.10.128journaljournalPhys. Rev. B volume10,pages128 (year1974)NoStop
 [Papanicolaou(1984)]papanicolaou1
  authorauthorN. Papanicolaou, 
  http://dx.doi.org/10.1016/0550-3213(84)90268-2journaljournalNucl. Phys. B volume240,pages281  (year1984)NoStop
 [Shastry and Sutherland(1981)]ShastryS
  authorauthorB. S. Shastry and authorB. Sutherland, 
  https://doi.org/10.1016/0378-4363(81)90838-XjournaljournalPhysica B+C volume108, pages1069  (year1981)NoStop
 [Miyahara and Ueda(1999)]MiyaharaU
  authorauthorS. Miyahara and authorK. Ueda,  10.1103/PhysRevLett.82.3701journaljournalPhys. Rev. Lett. volume82,pages3701 (year1999)NoStop
 [Knetter et al.(2000)Knetter, Bühler, Müller-Hartmann, andUhrig]knetter2000
  authorauthorC. Knetter, authorA. Bühler,
  authorE. Müller-Hartmann,   and authorG. S. Uhrig,  10.1103/PhysRevLett.85.3958journaljournalPhys. Rev. Lett. volume85, pages3958 (year2000)NoStop
 [Totsuka et al.(2001)Totsuka, Miyahara, and Ueda]totsuka2001
  authorauthorK. Totsuka, authorS. Miyahara,   and authorK. Ueda,  10.1103/PhysRevLett.86.520journaljournalPhys. Rev. Lett. volume86, pages520 (year2001)NoStop
 [Momoi and Totsuka(2000b)]momoi2000a
  authorauthorT. Momoi and authorK. Totsuka,  10.1103/PhysRevB.61.3231journaljournalPhys. Rev. B volume61,pages3231 (year2000b)NoStop
 [Note2()]Note2
  noteThe field-induced antiferroquadrupolar (AFQ) phase in the
  spin-gapped systems is also described as a low-density condensate of bound
  magnon pairs, as same as the AFQ phase near the saturation field. The vacuum
  is a singlet ground state and the bosonic particles are bound triplon pairs
  in the former case, whereas the vacuum is the fully polarized state and the
  bosons are bound magnon pairs in the latter.Stop
 [Romhányi et al.(2011)Romhányi, Totsuka, and Penc]romhanyi2011
  authorauthorJ. Romhányi, authorK. Totsuka,
   and authorK. Penc,  10.1103/PhysRevB.83.024413journaljournalPhys. Rev. B volume83, pages024413 (year2011)NoStop