††institutetext: 1 Department of Physical Sciences, Oklahoma State University,
Stillwater, Oklahoma 74078, USA††institutetext: 2 LAPTh, Université Savoie
Mont Blanc, CNRS, B.P. 110, F-74941 Annecy Cedex, France††institutetext: 3
Centre for High Energy Physics, Indian Institute of Science, Bengaluru 560012,
India

# Is the light neutralino thermal dark matter in the MSSM ruled out?

Rahool Kumar Barman1, Genevieve Bélanger2, Biplob Bhattacherjee3, Rohini M.
Godbole3, Rhitaja Sengupta3<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We explore the parameter space of the phenomenological Minimal Supersymmetric
Standard Model (pMSSM) with a light neutralino thermal dark matter
($M_{\widetilde{\chi}_{1}^{0}}\leq m_{h}/2$) for both positive and negative
values of the higgsino mass parameter ($\mu$) that is consistent with current
collider and astrophysical constraints. Our investigation shows that the
recent experimental results from the LHC as well as from direct detection
searches for dark matter by the LUX-ZEPLIN collaboration basically rule out
the $\mu>0$ scenario while only allowing a very narrow region with light
electroweakinos in the $\mu<0$ scenario. These are well within the reach of
the Run-3 of LHC and dedicated efforts to probe this region should be pursued.

##  1 Introduction

The R-parity conserved (RPC) scenario of the minimal supersymmetric extension
of the Standard Model (MSSM) has been among the most favourable choices for
exploring physics beyond the Standard Model (BSM). The RPC-MSSM scenario
alleviates the “naturalness” problem Gildener:1976ih ; PhysRevD.20.2619  in
Standard Model (SM), while also providing a SM-like Higgs boson ($h$) with
mass $m_{h}\sim 125~{}\mathrm{GeV}$ and a stable lightest supersymmetric
particle (LSP), typically the neutralino $\widetilde{\chi}_{1}^{0}$, which can
be a cold dark matter (DM) candidate. The case of the light neutralino
$m_{\tilde{\chi}_{1}^{0}}\leq m_{h}/2$ is of special interest since it is
kinematically feasible for the latter to decay invisibly through
$h\to\tilde{\chi}_{1}^{0}\tilde{\chi}_{1}^{0}$, thus providing an additional
signature for dark matter in the Higgs sector. Several studies have explored
the prospect of a light neutralino DM in the MSSM (in the constrained MSSM,
cMSSM and the phenomenological MSSM, pMSSM) considering the various
experimental constraints of the time PhysRevD.37.719 ; Djouadi:1996mj ;
Belanger:2000tg ; Belanger:2001am ; Belanger:2003wb ; Calibbi:2011ug ;
Dreiner:2012ex ; Ananthanarayan:2013fga ; Calibbi:2013poa ; Belanger:2013pna ;
Han:2014nba ; Belanger:2015vwa ; Hamaguchi:2015rxa ; Barman:2017swy ;
Pozzo:2018anw ; Wang:2020dtb ; KumarBarman:2020ylm ; VanBeekveld:2021tgn .
Collider experiments, like ATLAS and CMS, have made available the latest
results of searches of heavy Higgs bosons ATLAS:2020zms , direct searches of
charginos and neutralinos CMS:2020bfa ; ATLAS:2021moa ; ATLAS:2021yqv ;
CMS:2022sfi , as well as the invisible decay of the SM Higgs boson
ATLAS:2022yvh . The XENON-1T, PICO-60, PandaX-4T, and LUX-ZEPLIN (LZ)
collaborations have also published limits on the DM direct detection (DD)
cross-sections $-$ both spin-dependent (SD) and spin-independent (SI)
XENON:2018voc ; XENON:2019rxp ; PICO:2019vsc ; PandaX-4T:2021bab ;
Aalbers:2022fxq . Among these, the results from the LZ collaboration are the
most stringent ones for the SI DD cross-sections Aalbers:2022fxq . With the
advent of these results, it becomes important that we revisit the MSSM
parameter space of a light neutralino DM which can contribute to the invisible
decay of the SM Higgs boson.

In this paper, we study the current status of the light neutralino DM in MSSM
for both positive and negative values of the higgsino mass parameter, $\mu$.
It should be noted that a positive value of $\mu$ is indicated for a
supersymmetric explanation for the discrepancy between the experimentally
measured value of the $(g-2)_{\mu}$, and the SM prediction Muong-2:2021ojo .
Due to the prevalent uncertainties in the estimation of the hadronic
contributions in the SM prediction, we prefer to have an agnostic attitude
towards the sign of $\mu$. As the Large Hadron Collider (LHC) is gearing up
for Run-3 and will start collecting data soon, a careful study of the overall
status of this scenario is very timely to identify the interesting regions of
the parameter space which can be a focal point of the LHC searches at Run-3.

##  2 Current status of light neutralino dark matter in the MSSM

We consider the pMSSM parameter space with the parameters defined at the
electroweak scale. Our focus is the light neutralino sector with
$m_{\widetilde{\chi}_{1}^{0}}\leq m_{h}/2$ such that it is kinematically
feasible for the SM-like Higgs boson to decay into $\widetilde{\chi}_{1}^{0}$
pair and it can potentially contribute to the invisible decay mode of the
Higgs boson. The input parameters which capture the physics of the Higgs and
electroweakino sectors are: $M_{1}$, the bino mass, $M_{2}$, the wino mass
($M_{1}$ and $M_{2}$ are also collectively referred to as the gaugino masses),
$\mu$, the higgsino mass, $\tan\beta$, the ratio of the Higgs vacuum
expectation value, $M_{A}$, the pseudoscalar mass, $M_{\tilde{Q}_{3l}}$,
$M_{\tilde{t}_{R}}$, $M_{\tilde{b}_{R}}$, the mass of the third generation
squarks, $A_{t}$, trilinear coupling of the stop, and $M_{3}$, the mass of the
gluino. A random scan is performed over the following range of input
parameters:

$\displaystyle 30~{}{\rm GeV}<M_{1}<100~{}{\rm GeV},~{}1~{}{\rm
TeV}<M_{2}<3~{}{\rm TeV},$ $\displaystyle 100~{}{\rm GeV}<|\mu|<~{}2~{}{\rm
TeV},~{}2<\tan{\beta}<50,$ $\displaystyle~{}100~{}{\rm GeV}<M_{A}<5~{}{\rm
TeV},~{}3~{}{\rm TeV}<M_{\tilde{Q}_{3L}}<10~{}{\rm TeV},$ $\displaystyle
3~{}{\rm TeV}<M_{\tilde{t}_{R}}<10~{}{\rm TeV},~{}3~{}{\rm
TeV}<M_{\tilde{b}_{R}}<10~{}{\rm TeV},$ $\displaystyle-10~{}{\rm
TeV}<A_{t}<10~{}{\rm TeV},~{}2~{}{\rm TeV}<M_{3}<5~{}{\rm TeV}$

Since we are interested in light neutralino, it shall dominantly have bino
($\tilde{B}$) component, and therefore, we have scanned $M_{1}$ in the very
low mass region. The coupling of $Z$ and $h$ bosons to a pair of
$\tilde{\chi}_{1}^{0}$ depends also on its higgsino ($\tilde{H}$) and wino
($\tilde{W}$) components, and therefore to circumvent the overabundance of
$\tilde{\chi}_{1}^{0}$ as the DM candidate, we require some $\tilde{H}$ or
$\tilde{W}$ component in it. We are mostly interested in higgsino-like next-to
lightest supersymmetric partner (NLSP) in the present work due to the existing
stronger limits on wino-like NLSPs, and hence $M_{2}$ is varied starting from
1 TeV, and $\mu$ is varied starting from a comparatively lower value of 100
GeV. We cannot push $\mu$ below 100 GeV due to the existing limits on
charginos from the Large Electron Positron (LEP) collider experiments. We have
scanned the parameter space with both positive and negative values of $\mu$.
We have fixed the masses of the first and second generation squarks at 5 TeV,
the masses of all the three generations of sleptons at 2 TeV, and all their
trilinear couplings at zero, in order to keep them decoupled from the particle
spectrum. We perform a dedicated scan where we dynamically tune the $M_{1}$
parameter to keep $m_{\tilde{\chi}_{1}^{0}}$ within a window of 5 GeV of the
$Z$ mass and 3 GeV of the calculated $h$ mass for populating the funnel
regions properly.

We have used FeynHiggs 2.18.1 Heinemeyer:1998yj ; Heinemeyer:1998np ;
Degrassi:2002fi ; Frank:2006yh ; Hahn:2013ria ; Bahl:2016brp ; Bahl:2017aev ;
Bahl:2018qog  to generate the SUSY spectra corresponding to the various sets
of input parameters, and to calculate the Higgs boson mass, and decays in the
Higgs sector. We assume that the lightest CP-even Higgs boson of MSSM is the
SM-like Higgs boson, observed by the ATLAS and CMS collaborations, and the
combined measured mass is quoted as $m_{h}=125.09\pm 0.21({\rm stat})\pm
0.11({\rm syst})$ GeV ATLAS:2015yey . We apply a conservative constraint on
its theoretically calculated mass, $122{\rm~{}GeV}<m_{h}<128{\rm~{}GeV}$.
Since our scan starts from very low ${\rm tan}\beta$ values where satisfying
the observed Higgs boson mass will require high stop masses and $A_{t}$
values, we need to ensure that the latter is not large enough leading to color
and charged breaking minima (CCB) Camargo-Molina:2013sta ; Chowdhury:2013dka ;
Blinov:2013fta . We have implemented the constraint from Ref.
Chowdhury:2013dka  which showed that the electroweak vacuum becomes metastable
and unstable for $|X_{t}|\gtrsim\sqrt{6m_{\tilde{t}_{1}}m_{\tilde{t}_{2}}}$
unless $\mu\sim m_{\tilde{t}_{L},\tilde{t}_{R}}$, where $X_{t}=A_{t}-\mu/{\rm
tan}\beta$.

Next, we apply limits on the partial decay width of the invisible decay of
$Z$-boson from new physics, $\Gamma_{\rm inv}^{\rm new}<2$ MeV ALEPH:2005ab ,
chargino mass, $m_{\chi_{1}^{\pm}}>103$ GeV OPAL:2003wxm , and cross-section
of associated production of neutralinos in final states with jets,
$\sigma(e^{+}e^{-}\rightarrow\tilde{\chi}_{1}^{0}\tilde{\chi}_{2}^{0})\times{\rm
Br}(\tilde{\chi}_{2}^{0}\rightarrow\tilde{\chi}_{1}^{0}+{\rm
jets})+\sigma(e^{+}e^{-}\rightarrow\tilde{\chi}_{1}^{0}\tilde{\chi}_{3}^{0})\times{\rm
Br}(\tilde{\chi}_{3}^{0}\rightarrow\tilde{\chi}_{1}^{0}+{\rm jets})<0.1$ pb
OPAL:2003wxm , as obtained from experiments at the LEP. Subsequently, we add
the flavor physics constraints on various observables, like the branching
fractions of processes $b\rightarrow s\gamma$,
$B_{s}\rightarrow\mu^{+}\mu^{-}$, and $B\rightarrow\tau\nu$ which are required
to satisfy $3.00\times 10^{-4}<{\rm Br}(b\rightarrow s\gamma)<3.64\times
10^{-4}$ HFLAV:2016hnz , $1.66\times 10^{-9}<{\rm
Br}(B_{s}\rightarrow\mu^{+}\mu^{-})<4.34\times 10^{-9}$ CMS:2014xfa , and
$0.78<({\rm Br}(B\rightarrow\tau\nu))_{\rm obs}/({\rm
Br}(B\rightarrow\tau\nu))_{\rm SM}<1.78$ Belle:2010xzn , respectively. We have
used MicrOMEGAS 5.2.13 Belanger:2004yn ; Belanger:2006is ; Belanger:2008sj ;
Belanger:2010gh ; Belanger:2013oya ; Belanger:2020gnr  to calculate both the
LEP and flavor physics observables. Constraints from flavor physics push
$M_{A}$ to be greater than $\sim 800$ GeV irrespective of ${\rm tan}\beta$.

On the remaining parameter space points, we apply the limits from the signal
strength measurements of the SM Higgs boson implemented in HiggsSignal 2.6.2
Bechtle:2013xfa ; Stal:2013hwa ; Bechtle:2014ewa , as well as limits from the
collider searches of heavy Higgs bosons in various final-state channels at
collider experiments, like ATLAS and CMS, using the HiggsBounds 5.10.0
Bechtle:2008jh ; Bechtle:2011sb ; Bechtle:2012lvg ; Bechtle:2013wla ;
Bechtle:2015pma  package. The recent search of heavy Higgs bosons decaying to
$\tau$ leptons at ATLAS ATLAS:2020zms  excludes a large part of high ${\rm
tan}\beta$ region for $M_{A}\lesssim 1$ TeV. The parameter space also has to
satisfy the recent limit on the branching of the SM Higgs boson to decay into
invisible particles. We apply the strongest limit which comes from the result
of the recent search for invisible Higgs boson decays, where the latter is
produced in vector-boson fusion (VBF), using an integrated luminosity of 139
fb-1 at ATLAS ATLAS:2022yvh , which restricts
Br($h\rightarrow$invisible)$<0.145$. We refer to all the constraints related
to the Higgs bosons together as the “Higgs constraints” hereafter for
simplicity, which includes the SM-like Higgs mass constraint and constraints
from HiggsSignal 2.6.2, HiggsBounds 5.10.0, and the invisible decay of the SM-
like Higgs boson.

The lightest supersymmetric particle (LSP), $\tilde{\chi}_{1}^{0}$, is a
viable DM candidate in the MSSM, having a thermal freeze-out production in the
early Universe. In the standard cosmology, we require the relic density of the
LSP ($\Omega_{\rm LSP}$) to be equal to the observed DM relic density as
measured by the PLANCK collaboration $\Omega^{\rm obs}_{\rm DM}h^{2}=0.120\pm
0.001$ Aghanim:2018eyx , which assuming a $2\sigma$ interval can vary from
0.118-0.122. Lifting up the requisite that the neutralino LSP forms 100% of
the observed DM relic owing to the possibility of multicomponent DM, we can
modify the relic density constraint to $\Omega_{\rm LSP}\lesssim 0.122$, where
MicrOMEGAS 5.2.13 is used to compute the relic density of
$\tilde{\chi}_{1}^{0}$.

In addition to the relic density constraint, we need to take into
consideration results from the current DD experiments of DM. These experiments
constrain the spin-dependent DM-neutron (SDn) and DM-proton (SDp) as well as
the spin-independent direct detection cross-sections of the lightest
neutralino ($\tilde{\chi}_{1}^{0}$), which is the DM candidate, as a function
of its mass. We use MicrOMEGAS 5.2.13 to compute these cross-sections for the
LSP and then compare them with the 90% confidence level (CL) upper limits
quoted by the PICO-60 (SDp PICO:2019vsc ), XENON-1T (SI XENON:2018voc  and SDn
XENON:2019rxp ), and PandaX-4T (SI PandaX-4T:2021bab ) experiments. The DD
limits are placed assuming that a single DM candidate constitutes the entire
relic. Therefore, if the neutralino DM is underabundant, i.e., $\Omega_{\rm
LSP}<0.122$, then the DD limits are applied on the scaled cross-sections,
where the scaling factor is as follows:

$\xi=\frac{\Omega_{\rm LSP}}{0.122}$ (1)

Being underabundant, the $\chi_{1}^{0}$ component of DM might have evaded the
DD experiments, and therefore, the DD limits on the cross-section of the
lightest neutralino weaken in this scenario.

Moreover, we need to consider the results of direct searches for chargino and
neutralino at the ATLAS and CMS experiments of the LHC. We use the SModelS
2.2.0 Kraml:2013mwa ; Ambrogi:2017neo ; Dutta:2018ioj ; Heisig:2018kfq ;
Ambrogi:2018ujg ; Khosa:2020zar ; Alguero:2020grj ; Alguero:2021dig  package
to implement the electroweakino search constraints on our scanned parameter
space. This version of SModelS includes implementation of results from the
recent search for electroweakinos in the leptonic final states at CMS
CMS:2020bfa  and ATLAS ATLAS:2021moa  and in the hadronic final states at
ATLAS ATLAS:2021yqv , all of which play significant roles in excluding a large
range of $m_{\tilde{\chi}_{2}^{0}}$, especially with the ATLAS analysis
extending the sensitivity to high masses with the hadronic final states.

Figure 1: Scaled SI DM-nucleon cross-section ($\sigma_{SI}\times\xi$) for
$\mu>0$ as a function of the mass of the LSP neutralino DM in the region of
parameter space satisfying LEP, flavor, Higgs, relic density and direct
detection (DD) constraints (grey), along with the regions surviving present
electroweakino searches and invisible branching fraction of Higgs boson with
$m_{\tilde{\chi}_{2}^{0}}$ in the colorbar.

We apply the constraints on our scanned parameter space in two steps $-$ first
“Set A” with constraints from LEP, flavor, Higgs constraints, relic density
and the DD experiments XENON-1T, PICO-60 and PandaX-4T, then “Set B” with
electroweakino constraints from the LHC. Fig. 1 shows the scaled ($\xi$ as
defined in Eq. 1) SI DD cross-sections for $\mu>0$ of the allowed parameter
space after applying all the constraints from “Set A” in grey, and then after
adding constraints from “Set B” in colored points, with the colorbar showing
the mass of $\tilde{\chi}_{2}^{0}$. We observe that the recent electroweakino
constraints have played a crucial role in completely excluding the $Z$-funnel
region of positive $\mu$. However, the colored points in the $h$-funnel region
indicate that the electroweakino searches, as implemented in the SModelS
framework, still allows regions of the parameter space where
$\tilde{\chi}_{2}^{0}$ can be very light ($\sim$ 146-158 GeV) or heavy ($\sim$
855-1330 GeV). Once we overlay the recent 90% CL upper limits on $\sigma_{SI}$
from the LZ experiment Aalbers:2022fxq  and the projected limits from the
XENON-nT experiment on Fig. 1, we find that the current LZ result is strong
enough to exclude the entire allowed parameter region, even in the $h$-funnel,
which survived the present-day electroweakino constraints.

We would like to point out a caveat that one can evade the strong bounds from
LZ if the relic density condition is relaxed such that the coupling can be
small enough to satisfy the strong DD constraint. Assuming that theoretically
we are overestimating the relic density by, say, 20%, we should actually
compare 80% of the value of $\Omega_{\rm LSP}$ from MicrOMEGAS with 0.122,
which leaves a narrow region of parameter space allowed by even the LZ limits.
However these lie just below the current result and might be very easily
accessible to the full exposure of the LZ experiment. Another way to relax the
upper limit on the relic density of the LSP DM can be provided in scenarios of
non-standard cosmology since the observed relic can be satisfied due to the
presence of other physical states which decay later in the Universe, thus,
increasing the entropy density of the Universe. We have observed that there is
a linear relationship between how much we can relax the relic density in the
non-standard cosmological model and the improvement in the DD experimental
limit needed to fully probe that scenario, viz. to completely probe a scenario
in non-standard cosmology where the relic density can be relaxed by a factor
of 5 requires a five-fold improvement from the current LZ limit. We land up
with heavy higgsinos in both these cases, and the electroweakino searches
further push them beyond 850 GeV.

Figure 2: Scaled SI DM-nucleon cross-section ($\sigma_{SI}\times\xi$) (left)
and scaled SD DM-neutron cross-section ($\sigma_{SDn}\times\xi$) (right) for
$\mu<0$ as a function of the mass of the LSP neutralino DM in the region of
parameter space satisfying LEP, flavor, Higgs, relic density and direct
detection (DD) constraints (grey), along with the regions surviving present
electroweakino searches and invisible branching fraction of Higgs boson with
$m_{\tilde{\chi}_{2}^{0}}$ in the colorbar.

Let us now investigate what happens in the negative $\mu$ scenario. Fig. 2
shows a similar plot of the scaled SDn (left) and SI (right) DD cross-sections
for $\mu<0$ of the allowed parameter space, where the colors have the same
meaning as described for Fig. 1. From the left panel of Fig. 2, we observe
that we have parameter region remaining in both the $Z$ and $h$-funnels which
satisfy all the constraints from “Set A” and “Set B”, where the latter
restricts $m_{\tilde{\chi}_{2}^{0}}$ to either very small values ($\sim$
130-142 GeV and $\sim$ 147-157 GeV in the $Z$ and $h$-funnel respectively) or
high values ($\sim$ 856-1385 GeV in the $h$-funnel). The right panel shows the
allowed parameter space in the $m_{\tilde{\chi}_{1}^{0}}$-$\sigma_{\rm SI}$
plane and we observe that the recent LZ result Aalbers:2022fxq  which has
collected data with an exposure of 60 days has excluded a large region in the
$h$-funnel region, leaving behind a narrow marginally allowed region where the
$m_{\tilde{\chi}_{2}^{0}}$ is very light ($\sim$ 149-155 GeV), which we can
expect to be probed in the near future by the LZ experiment with its full 1000
day exposure. On the other hand, the $Z$-funnel is not affected much by the LZ
limits. The projected upper limit on $\sigma_{SDn}$ from the future XENON-nT
experiment (left panel of Fig. 2) shows that it can probe the parameter space
with light $\tilde{\chi}_{2}^{0}$ in both the $Z$ and $h$-funnel regions.

Figure 3: Left: Allowed parameter space for $\mu<0$ after satisfying the LEP,
flavor, Higgs constraints, relic density DM DD constraints from the XENON-1T,
PICO-60 and PandaX-4T experiments (“Set A” in yellow circles), overlayed with
additional constraints from electroweakino searches at the LHC (“Set B” in
light green circles) and the recent constraint on the SI DD cross-section from
the LZ experiment (dark green stars) in the
$m_{\tilde{\chi}_{1}^{0}}$-$m_{\tilde{\chi}_{2}^{0}}$ plane; Right: $R$-values
from SModelS for the allowed parameter space.

The left plot of Fig. 3 shows the parameter space for $\mu<0$ in the
$m_{\tilde{\chi}_{1}^{0}}$-$m_{\tilde{\chi}_{2}^{0}}$ plane, where we can
observe that although the electroweakino searches allow heavy higgsinos in the
$h$-funnel, the very recent LZ DD limit restricts the allowed parameter space
to only light higgsinos, which are particularly interesting to probe at
collider experiments. The right plot of Fig. 3 shows the SModelS $R$-values in
the colorbar, varying from 0 to 1, of the parameter space allowed by all the
current constraints to give the reader an idea of how far these points are
from the present available limits from the LHC experiments $-$ a smaller
$R$-value indicates that the parameter space point lies way outside the
current limit, whereas a $R$-value close to 1 indicates that it lies just on
the border of the limit. The allowed points evade the ATLAS searches for
charginos and neutralinos in the leptonic final states at $\sqrt{s}=8$ TeV
with 20.3 fb-1 of data ATLAS:2014zve ; ATLAS:2014ikz  and at $\sqrt{s}=13$ TeV
with 139 fb-1 of data ATLAS:2021moa , with some of them having very small
$R$-values.

These regions with light charginos and neutralinos, which are evading the
present constraints due to off-shell $Z$ or $h$ bosons, will be very important
for Run-3. To understand their prospect, we perform an analysis of the low
mass higgsino-like electroweakinos in the leptonic final state at
$\sqrt{s}=14$ TeV using the XGBOOST framework in this work. We consider the
process
$pp\rightarrow\tilde{\chi}_{1}^{\pm}\tilde{\chi}_{2}^{0}/\tilde{\chi}_{1}^{\pm}\tilde{\chi}_{3}^{0},~{}\tilde{\chi}_{1}^{\pm}\rightarrow
W^{\pm}\tilde{\chi}_{1}^{0},~{}\tilde{\chi}_{2}^{0}/\tilde{\chi}_{3}^{0}\rightarrow
f\bar{f}\tilde{\chi}_{1}^{0}$ with $m_{\tilde{\chi}_{1}^{\pm}}=142.4$ GeV,
$m_{\tilde{\chi}_{2}^{0}}=149.7$ GeV, $m_{\tilde{\chi}_{3}^{0}}=151.9$ GeV,
and $m_{\tilde{\chi}_{1}^{0}}=61.2$ GeV where $f$ is a SM fermion. We study 11
possible backgrounds for this process $-$ $lll\nu$ ($l\equiv e,\mu,\tau$),
$ZZ$, $t\bar{t}$, $VVV$, $Wh$, $Zh$, ggF and VBF production of $h$ with
$h\rightarrow ZZ^{*}$, $t\bar{t}h$, $t\bar{t}W$, and $t\bar{t}Z$.

We restrict to the leptonic final state which is cleaner for a lighter
benchmark, such as ours. We perform an analysis of the $3l+\rm
E{\\!\\!\\!/}_{T}$ final state where we require exactly three leptons
satisfying $p_{T}>25,25,20$ GeV and $|\eta|<2.4$, and we have put a veto on
$b$-jets with $p_{T}>30$ GeV and $|\eta|<2.5$. In our signal benchmark, since
we do not have any on-shell $Z$-boson, we also veto events where the invariant
mass of a pair of same flavor opposite sign (SFOS) leptons lie within 10 GeV
window of $m_{Z}=91.2$ GeV. After these preselections, we train our signal and
background samples using XGBOOST 111https://xgboost.readthedocs.io/en/stable/
with a set of 21 variables $-$ transverse momenta ($p_{T}$) of the three
leptons, transverse mass ($M_{T}$) and contransverse mass ($M_{CT}$) of each
of the three leptons with the $\rm E{\\!\\!\\!/}_{T}$, minimum and maximum
values of $\Delta R$ between opposite sign lepton pairs along with their
$\Delta\eta$ values, invariant mass of the opposite sign lepton pairs with
minimum and maximum $\Delta R$, missing transverse momentum, number of jets in
the event with the $p_{T}$ of the two leading jets, scalar sum of $p_{T}$ of
all the jets in the event ($H_{T}$), and the invariant mass of the three
leptons. We combine the backgrounds, weighted according to their cross-
sections and use unity weight for the signal 222Following are the
hyperparameters of the XGBOOST model: ‘objective’:‘multi:softprob’,
‘colsample_bytree’:0.3, ‘learning_rate’:0.1, ‘num_class’:12, ‘max_depth’:7,
‘alpha’:5, ‘eval_metric’:‘mlogloss’, ‘num_round’:1000,
‘early_stopping_rounds’:3. After training and validation of the model, we use
it to discriminate the signal benchmark from each of the background classes by
computing the significance of observing the signal over the background events.
At a 14 TeV collider with 137 fb-1 of integrated luminosity, we expect to
observe 310 signal events ($N_{S}$) and a total of 331 background events
($N_{B}$) for a threshold of 0.9 on the XGBOOST output, which upon adding a
20% (50%) systematic uncertainty translates to a significance (using the
formula in Ref. Adhikary:2020cli ) of 3.6 (1.5). We find that the result
sensitively depends on the systematic uncertainty, which might be dominant for
light electroweakinos. These might be either already ruled out by the Run-2
data in analyses which are not yet published and implemented in SModelS, or
can be probed in the Run-3 of LHC, given the systematic uncertainties can be
brought under control.

The future lepton colliders like ILC and CEPC will be crucial for precision
measurements of Higgs boson. The projected upper limit on the invisible
branching of the Higgs boson is 0.4% at ILC Asner:2013psa  and 0.3% at CEPC
An:2018dwb . Although these can probe a significant part of the allowed
parameter space in the $\mu<0$ case, we still have regions with
Br($h\rightarrow$invisible$<0.003$) in both the $Z$ and $h$-funnels. In the
$\mu<0$ case, the partial decay width of the $Z$ boson to
$\tilde{\chi}_{1}^{0}$ ($\Gamma_{\rm inv}^{\rm new}$) is always less than 0.1
MeV for the allowed parameter region that we obtain. Therefore, we do not
expect the Giga-$Z$ option of ILC, which is expected to have a modest
improvement over LEP Carena:2003aj , to be sensitive to this region.

It is interesting to note that if in the future DD experiments, we discover a
light DM in the $Z$-funnel, it will mostly indicate a negative value of $\mu$.
Although we have seen in the present study that the $\mu>0$ scenario is in
severe tension with the experimental results, we have also discussed some
caveats where they can evade these constraints. In the $h$-funnel region, one
major difference between the $\mu>0$ and $\mu<0$ case is that in the former
the experimental limits along with relaxing relic density condition allows the
electroweakinos to be heavy (greater than a TeV), unlike the much lighter
electroweakinos (around 150 GeV) obtained for latter. Thus, an observation of
DM signal in the $h$-funnel from DD experiments would require an additional
observation of a signal in the collider experiments to understand the sign of
the $\mu$ parameter. Further, observation of heavy higgsinos in the collider
experiments might even be a harbinger of non-standard cosmology.

##  3 Conclusion

In summary, this work shows that the current experiments, especially the
recent results from the electroweakino searches at the LHC and the LZ dark
matter DD experiment have severely constrained the positive $\mu$ scenario,
whereas they have squeezed the parameter space to light electroweakinos in the
negative $\mu$ scenario for a light neutralino DM in the MSSM. Such light
higgsinos are also motivated from naturalness arguments. The Run-3 of LHC
shall target this low mass electroweakinos and perform dedicated analyses to
close this narrow gap. The experimental collaborations should look for any
kind of loopholes in the previous analyses of light electroweakino searches
and design analyses to cover them. They should also exhaust the possibility of
the presence of other light degenerate SUSY particles, which might lead to
difficult signatures to observe in the collider experiments. To conclude, at
present we are at a very exciting juncture where the experiments lined up
might exclude the possibility of a light neutralino DM in MSSM altogether, or
we might be very close to start observing the first hints of new physics.

## Acknowledgement

B.B. and R.S. thank Prabhat Solanki and Camellia Bose for useful discussions.
R.K.B. thanks the U.S. Department of Energy for the financial support, under
grant number DE-SC0016013.

## References

  * (1) E. Gildener and S. Weinberg, “Symmetry Breaking and Scalar Bosons,” Phys. Rev. D, vol. 13, p. 3333, 1976.
  * (2) L. Susskind, “Dynamics of spontaneous symmetry breaking in the weinberg-salam theory,” Phys. Rev. D, vol. 20, pp. 2619–2625, Nov 1979.
  * (3) K. Griest and H. E. Haber, “Invisible decays of higgs bosons in supersymmetric models,” Phys. Rev. D, vol. 37, pp. 719–728, Feb 1988.
  * (4) A. Djouadi, P. Janot, J. Kalinowski, and P. M. Zerwas, “SUSY decays of Higgs particles,” Phys. Lett. B, vol. 376, pp. 220–226, 1996.
  * (5) G. Belanger, F. Boudjema, F. Donato, R. Godbole, and S. Rosier-Lees, “SUSY Higgs at the LHC: Effects of light charginos and neutralinos,” Nucl. Phys. B, vol. 581, pp. 3–33, 2000.
  * (6) G. Belanger, F. Boudjema, A. Cottrant, R. M. Godbole, and A. Semenov, “The MSSM invisible Higgs in the light of dark matter and g-2,” Phys. Lett. B, vol. 519, pp. 93–102, 2001.
  * (7) G. Belanger, F. Boudjema, A. Cottrant, A. Pukhov, and S. Rosier-Lees, “Lower limit on the neutralino mass in the general MSSM,” JHEP, vol. 03, p. 012, 2004.
  * (8) L. Calibbi, T. Ota, and Y. Takanishi, “Light Neutralino in the MSSM: a playground for dark matter, flavor physics and collider experiments,” JHEP, vol. 07, p. 013, 2011.
  * (9) H. K. Dreiner, J. S. Kim, and O. Lebedev, “First LHC Constraints on Neutralinos,” Phys. Lett. B, vol. 715, pp. 199–202, 2012.
  * (10) B. Ananthanarayan, J. Lahiri, P. N. Pandita, and M. Patra, “Invisible decays of the lightest Higgs boson in supersymmetric models,” Phys. Rev. D, vol. 87, no. 11, p. 115021, 2013.
  * (11) L. Calibbi, J. M. Lindert, T. Ota, and Y. Takanishi, “Cornering light Neutralino Dark Matter at the LHC,” JHEP, vol. 10, p. 132, 2013.
  * (12) G. Bélanger, G. Drieu La Rochelle, B. Dumont, R. M. Godbole, S. Kraml, and S. Kulkarni, “LHC constraints on light neutralino dark matter in the MSSM,” Phys. Lett. B, vol. 726, pp. 773–780, 2013.
  * (13) T. Han, Z. Liu, and S. Su, “Light Neutralino Dark Matter: Direct/Indirect Detection and Collider Searches,” JHEP, vol. 08, p. 093, 2014.
  * (14) G. Belanger, D. Ghosh, R. Godbole, and S. Kulkarni, “Light stop in the MSSM after LHC Run 1,” JHEP, vol. 09, p. 214, 2015.
  * (15) K. Hamaguchi and K. Ishikawa, “Prospects for Higgs- and Z-resonant Neutralino Dark Matter,” Phys. Rev. D, vol. 93, no. 5, p. 055009, 2016.
  * (16) R. K. Barman, G. Belanger, B. Bhattacherjee, R. Godbole, G. Mendiratta, and D. Sengupta, “Invisible decay of the Higgs boson in the context of a thermal and nonthermal relic in MSSM,” Phys. Rev., vol. D95, no. 9, p. 095018, 2017.
  * (17) G. Pozzo and Y. Zhang, “Constraining resonant dark matter with combined LHC electroweakino searches,” Phys. Lett. B, vol. 789, pp. 582–591, 2019\.
  * (18) K. Wang and J. Zhu, “Funnel annihilations of light dark matter and the invisible decay of the Higgs boson,” Phys. Rev. D, vol. 101, no. 9, p. 095028, 2020.
  * (19) R. Kumar Barman, G. Belanger, and R. M. Godbole, “Status of low mass LSP in SUSY,” Eur. Phys. J. ST, vol. 229, no. 21, pp. 3159–3185, 2020.
  * (20) M. Van Beekveld, W. Beenakker, M. Schutten, and J. De Wit, “Dark matter, fine-tuning and $(g-2)_{\mu}$ in the pMSSM,” SciPost Phys., vol. 11, no. 3, p. 049, 2021.
  * (21) G. Aad et al., “Search for heavy Higgs bosons decaying into two tau leptons with the ATLAS detector using $pp$ collisions at $\sqrt{s}=13$ TeV,” Phys. Rev. Lett., vol. 125, no. 5, p. 051801, 2020.
  * (22) A. M. Sirunyan et al., “Search for supersymmetry in final states with two oppositely charged same-flavor leptons and missing transverse momentum in proton-proton collisions at $\sqrt{s}=$ 13 TeV,” JHEP, vol. 04, p. 123, 2021.
  * (23) G. Aad et al., “Search for chargino–neutralino pair production in final states with three leptons and missing transverse momentum in $\sqrt{s}=13$ TeV pp collisions with the ATLAS detector,” Eur. Phys. J. C, vol. 81, no. 12, p. 1118, 2021.
  * (24) G. Aad et al., “Search for charginos and neutralinos in final states with two boosted hadronically decaying bosons and missing transverse momentum in $pp$ collisions at $\sqrt{s}$ = 13 TeV with the ATLAS detector,” Phys. Rev. D, vol. 104, no. 11, p. 112010, 2021.
  * (25) “Search for electroweak production of charginos and neutralinos at $\sqrt{s}$ =13 TeV in final states containing hadronic decays of WW, WZ, or WH and missing transverse momentum,” 5 2022.
  * (26) G. Aad et al., “Search for invisible Higgs-boson decays in events with vector-boson fusion signatures using 139 $\text{fb}^{-1}$ of proton-proton data recorded by the ATLAS experiment,” 2 2022.
  * (27) E. Aprile et al., “Dark Matter Search Results from a One Ton-Year Exposure of XENON1T,” Phys. Rev. Lett., vol. 121, no. 11, p. 111302, 2018\.
  * (28) E. Aprile et al., “Constraining the spin-dependent WIMP-nucleon cross sections with XENON1T,” Phys. Rev. Lett., vol. 122, no. 14, p. 141301, 2019.
  * (29) C. Amole et al., “Dark Matter Search Results from the Complete Exposure of the PICO-60 C3F8 Bubble Chamber,” Phys. Rev. D, vol. 100, no. 2, p. 022001, 2019.
  * (30) Y. Meng et al., “Dark Matter Search Results from the PandaX-4T Commissioning Run,” Phys. Rev. Lett., vol. 127, no. 26, p. 261802, 2021\.
  * (31) J. Aalbers et al., “First Dark Matter Search Results from the LUX-ZEPLIN (LZ) Experiment,” 7 2022.
  * (32) B. Abi et al., “Measurement of the Positive Muon Anomalous Magnetic Moment to 0.46 ppm,” Phys. Rev. Lett., vol. 126, no. 14, p. 141801, 2021\.
  * (33) S. Heinemeyer, W. Hollik, and G. Weiglein, “FeynHiggs: A Program for the calculation of the masses of the neutral CP even Higgs bosons in the MSSM,” Comput. Phys. Commun., vol. 124, pp. 76–89, 2000.
  * (34) S. Heinemeyer, W. Hollik, and G. Weiglein, “The Masses of the neutral CP - even Higgs bosons in the MSSM: Accurate analysis at the two loop level,” Eur. Phys. J. C, vol. 9, pp. 343–366, 1999.
  * (35) G. Degrassi, S. Heinemeyer, W. Hollik, P. Slavich, and G. Weiglein, “Towards high precision predictions for the MSSM Higgs sector,” Eur. Phys. J. C, vol. 28, pp. 133–143, 2003.
  * (36) M. Frank, T. Hahn, S. Heinemeyer, W. Hollik, H. Rzehak, and G. Weiglein, “The Higgs Boson Masses and Mixings of the Complex MSSM in the Feynman-Diagrammatic Approach,” JHEP, vol. 02, p. 047, 2007.
  * (37) T. Hahn, S. Heinemeyer, W. Hollik, H. Rzehak, and G. Weiglein, “High-Precision Predictions for the Light CP -Even Higgs Boson Mass of the Minimal Supersymmetric Standard Model,” Phys. Rev. Lett., vol. 112, no. 14, p. 141801, 2014.
  * (38) H. Bahl and W. Hollik, “Precise prediction for the light MSSM Higgs boson mass combining effective field theory and fixed-order calculations,” Eur. Phys. J. C, vol. 76, no. 9, p. 499, 2016.
  * (39) H. Bahl, S. Heinemeyer, W. Hollik, and G. Weiglein, “Reconciling EFT and hybrid calculations of the light MSSM Higgs-boson mass,” Eur. Phys. J. C, vol. 78, no. 1, p. 57, 2018.
  * (40) H. Bahl, T. Hahn, S. Heinemeyer, W. Hollik, S. Paßehr, H. Rzehak, and G. Weiglein, “Precision calculations in the MSSM Higgs-boson sector with FeynHiggs 2.14,” Comput. Phys. Commun., vol. 249, p. 107099, 2020.
  * (41) G. Aad et al., “Combined Measurement of the Higgs Boson Mass in $pp$ Collisions at $\sqrt{s}=7$ and 8 TeV with the ATLAS and CMS Experiments,” Phys. Rev. Lett., vol. 114, p. 191803, 2015.
  * (42) J. E. Camargo-Molina, B. O’Leary, W. Porod, and F. Staub, “Stability of the CMSSM against sfermion VEVs,” JHEP, vol. 12, p. 103, 2013.
  * (43) D. Chowdhury, R. M. Godbole, K. A. Mohan, and S. K. Vempati, “Charge and Color Breaking Constraints in MSSM after the Higgs Discovery at LHC,” JHEP, vol. 02, p. 110, 2014.  [Erratum: JHEP 03, 149 (2018)].
  * (44) N. Blinov and D. E. Morrissey, “Vacuum Stability and the MSSM Higgs Mass,” JHEP, vol. 03, p. 106, 2014.
  * (45) S. Schael et al., “Precision electroweak measurements on the $Z$ resonance,” Phys. Rept., vol. 427, pp. 257–454, 2006.
  * (46) G. Abbiendi et al., “Search for chargino and neutralino production at s**(1/2) = 192-GeV to 209 GeV at LEP,” Eur. Phys. J. C, vol. 35, pp. 1–20, 2004.
  * (47) Y. Amhis et al., “Averages of $b$-hadron, $c$-hadron, and $\tau$-lepton properties as of summer 2016,” Eur. Phys. J. C, vol. 77, no. 12, p. 895, 2017.
  * (48) V. Khachatryan et al., “Observation of the rare $B^{0}_{s}\to\mu^{+}\mu^{-}$ decay from the combined analysis of CMS and LHCb data,” Nature, vol. 522, pp. 68–72, 2015.
  * (49) K. Hara et al., “Evidence for $B^{-}->\tau^{-}\bar{\nu}$ with a Semileptonic Tagging Method,” Phys. Rev. D, vol. 82, p. 071101, 2010.
  * (50) G. Belanger, F. Boudjema, A. Pukhov, and A. Semenov, “micrOMEGAs: Version 1.3,” Comput. Phys. Commun., vol. 174, pp. 577–604, 2006.
  * (51) G. Belanger, F. Boudjema, A. Pukhov, and A. Semenov, “MicrOMEGAs 2.0: A Program to calculate the relic density of dark matter in a generic model,” Comput. Phys. Commun., vol. 176, pp. 367–382, 2007.
  * (52) G. Belanger, F. Boudjema, A. Pukhov, and A. Semenov, “Dark matter direct detection rate in a generic model with micrOMEGAs 2.2,” Comput. Phys. Commun., vol. 180, pp. 747–767, 2009.
  * (53) G. Belanger, F. Boudjema, P. Brun, A. Pukhov, S. Rosier-Lees, P. Salati, and A. Semenov, “Indirect search for dark matter with micrOMEGAs2.4,” Comput. Phys. Commun., vol. 182, pp. 842–856, 2011.
  * (54) G. Belanger, F. Boudjema, A. Pukhov, and A. Semenov, “micrOMEGAs$\\_$3: A program for calculating dark matter observables,” Comput. Phys. Commun., vol. 185, pp. 960–985, 2014.
  * (55) G. Belanger, A. Mjallal, and A. Pukhov, “Recasting direct detection limits within micrOMEGAs and implication for non-standard Dark Matter scenarios,” Eur. Phys. J. C, vol. 81, no. 3, p. 239, 2021.
  * (56) P. Bechtle, S. Heinemeyer, O. Stål, T. Stefaniak, and G. Weiglein, “$HiggsSignals$: Confronting arbitrary Higgs sectors with measurements at the Tevatron and the LHC,” Eur. Phys. J. C, vol. 74, no. 2, p. 2711, 2014\.
  * (57) O. Stål and T. Stefaniak, “Constraining extended Higgs sectors with HiggsSignals,” PoS, vol. EPS-HEP2013, p. 314, 2013.
  * (58) P. Bechtle, S. Heinemeyer, O. Stål, T. Stefaniak, and G. Weiglein, “Probing the Standard Model with Higgs signal rates from the Tevatron, the LHC and a future ILC,” JHEP, vol. 11, p. 039, 2014.
  * (59) P. Bechtle, O. Brein, S. Heinemeyer, G. Weiglein, and K. E. Williams, “HiggsBounds: Confronting Arbitrary Higgs Sectors with Exclusion Bounds from LEP and the Tevatron,” Comput. Phys. Commun., vol. 181, pp. 138–167, 2010.
  * (60) P. Bechtle, O. Brein, S. Heinemeyer, G. Weiglein, and K. E. Williams, “HiggsBounds 2.0.0: Confronting Neutral and Charged Higgs Sector Predictions with Exclusion Bounds from LEP and the Tevatron,” Comput. Phys. Commun., vol. 182, pp. 2605–2631, 2011.
  * (61) P. Bechtle, O. Brein, S. Heinemeyer, O. Stal, T. Stefaniak, G. Weiglein, and K. Williams, “Recent Developments in HiggsBounds and a Preview of HiggsSignals,” PoS, vol. CHARGED2012, p. 024, 2012.
  * (62) P. Bechtle, O. Brein, S. Heinemeyer, O. Stål, T. Stefaniak, G. Weiglein, and K. E. Williams, “$\mathsf{HiggsBounds}-4$: Improved Tests of Extended Higgs Sectors against Exclusion Bounds from LEP, the Tevatron and the LHC,” Eur. Phys. J. C, vol. 74, no. 3, p. 2693, 2014.
  * (63) P. Bechtle, S. Heinemeyer, O. Stal, T. Stefaniak, and G. Weiglein, “Applying Exclusion Likelihoods from LHC Searches to Extended Higgs Sectors,” Eur. Phys. J. C, vol. 75, no. 9, p. 421, 2015.
  * (64) N. Aghanim et al., “Planck 2018 results. VI. Cosmological parameters,” 2018.
  * (65) S. Kraml, S. Kulkarni, U. Laa, A. Lessa, W. Magerl, D. Proschofsky-Spindler, and W. Waltenberger, “SModelS: a tool for interpreting simplified-model results from the LHC and its application to supersymmetry,” Eur. Phys. J. C, vol. 74, p. 2868, 2014.
  * (66) F. Ambrogi, S. Kraml, S. Kulkarni, U. Laa, A. Lessa, V. Magerl, J. Sonneveld, M. Traub, and W. Waltenberger, “SModelS v1.1 user manual: Improving simplified model constraints with efficiency maps,” Comput. Phys. Commun., vol. 227, pp. 72–98, 2018.
  * (67) J. Dutta, S. Kraml, A. Lessa, and W. Waltenberger, “SModelS extension with the CMS supersymmetry search results from Run 2,” LHEP, vol. 1, no. 1, pp. 5–12, 2018.
  * (68) J. Heisig, S. Kraml, and A. Lessa, “Constraining new physics with searches for long-lived particles: Implementation into SModelS,” Phys. Lett. B, vol. 788, pp. 87–95, 2019.
  * (69) F. Ambrogi et al., “SModelS v1.2: long-lived particles, combination of signal regions, and other novelties,” Comput. Phys. Commun., vol. 251, p. 106848, 2020.
  * (70) C. K. Khosa, S. Kraml, A. Lessa, P. Neuhuber, and W. Waltenberger, “SModelS Database Update v1.2.3,” LHEP, vol. 2020, p. 158, 2020.
  * (71) G. Alguero, S. Kraml, and W. Waltenberger, “A SModelS interface for pyhf likelihoods,” Comput. Phys. Commun., vol. 264, p. 107909, 2021.
  * (72) G. Alguero, J. Heisig, C. Khosa, S. Kraml, S. Kulkarni, A. Lessa, H. Reyes-González, W. Waltenberger, and A. Wongel, “Constraining new physics with SModelS version 2,” 12 2021.
  * (73) G. Aad et al., “Search for direct production of charginos, neutralinos and sleptons in final states with two leptons and missing transverse momentum in $pp$ collisions at $\sqrt{s}=$ 8 TeV with the ATLAS detector,” JHEP, vol. 05, p. 071, 2014.
  * (74) G. Aad et al., “Search for direct production of charginos and neutralinos in events with three leptons and missing transverse momentum in $\sqrt{s}=$ 8TeV $pp$ collisions with the ATLAS detector,” JHEP, vol. 04, p. 169, 2014.
  * (75) A. Adhikary, N. Chakrabarty, I. Chakraborty, and J. Lahiri, “Probing the $H^{\pm}W^{\mp}Z$ interaction at the high energy upgrade of the LHC,” Eur. Phys. J. C, vol. 81, no. 6, p. 554, 2021.
  * (76) D. M. Asner et al., “ILC Higgs White Paper,” in Community Summer Study 2013: Snowmass on the Mississippi, 10 2013.
  * (77) F. An et al., “Precision Higgs physics at the CEPC,” Chin. Phys. C, vol. 43, no. 4, p. 043002, 2019.
  * (78) M. Carena, A. de Gouvea, A. Freitas, and M. Schmitt, “Invisible Z boson decays at e+ e- colliders,” Phys. Rev. D, vol. 68, p. 113007, 2003.


# Exact solvability and two-frequency Rabi oscillation in cavity-QED setup
with moving emitter

Mingzhu Weng  Center for Quantum Sciences and School of Physics, Northeast
Normal University, Changchun 130024, China  Zhihai Wang<EMAIL_ADDRESS>Center for Quantum Sciences and School of Physics, Northeast Normal
University, Changchun 130024, China

###### Abstract

In this paper, we investigate the energy spectrum and coherent dynamical
process in a cavity-QED setup with a moving emitter, which is subject to a
harmonic potential. We find that the vibration of the emitter will induce the
effective Kerr and optomechanical interactions. We generalize the Bogliubov
operators approach which dealt with quantum Rabi model, to our cavity-emitter-
vibration system and obtain the energy spectrum exactly. With the assistance
of Bogliubov operators approach, we obtain the energy spectrum of the system
exactly. Furthermore, we show that the dynamics of the system exhibit a two-
frequency Rabi oscillation behavior. We explain such behavior by
optomechanical interaction induced quantum transition between emitter-cavity
dressed states. We hope that the interaction between cavity mode and moving
emitter will provide a versatile platform to explore more exotic effects and
potential applications in cavity-QED scenario.

##  I introduction

Light-matter interaction is a fundamental topic in the modern physics, ranging
from quantum optics and quantum information processing to the condensed matter
physics. The cavity is usually used to adjust the emission of the emitter,
leading to the Purcell effect purcell1946 , which is a central concept in the
field of cavity quantum electrodynamics (QED). In the strong coupling regime,
the cavity-QED setup can be described by the Jaynes-Cummings (JC) model JC1963
, and the single and multiple photon quantum Rabi oscillation have been
studied broadly Agarwal1985 ; Brune1996 ; Garziano2015 .

In the traditional investigations on cavity-QED and waveguide-QED system, the
emitter is usually assumed to be static under the dipole approximation.
However, the vibration degrees of freedom of quantum emitters are recently
received more and more attentions. For example, in waveguide-QED setup, the
waveguide induced interaction between moving emitters has been deeply studied
for both of the cases when the velocity of the emitter is faster and slower
than that of the photons in the waveguide GC2017 ; ES2020 . Also, the motion
of the emitter also leads to the recoil effect QL2013 ; DB2008 ; FD2016 ,
which is predicted by the modulated single-photon scattering line shape. Even
in the cavity-QED setup, the motion of emitter also induces many interesting
phenomena and applications which are absent for the static ones. For example,
the oscillation collapse and revival of atomic transition probability XG1995 ;
LX2000 , the spatial decoherence LZ2005 ; LZ2010 ; LY2001 , the motional
$n$-phonon bundle state YG2021 ; CS2014 , the exotic photon statistics YZ2015
; KM2015 , as well as the dynamical Casimir effect AA2021 ; SS2009 ; OD2019 ;
HW2019 ; WQ2018 ; VM2018 ; SF2015 , just to name a few. On the other hand, in
the recent cavity-QED experiment, the Rydberg atom is usually subject to the
harmonic potential which is generated by the laser or magneto-optical
technology Anderson2011 ; Tikman2016 ; Bounds2018 . Therefore, it naturally
motivates us to investigate the exact energy spectrum and dynamical evolution
of the cavity-QED setup with moving emitter which yields to a harmonic
potential.

In this work, we focus on the quantum effect of the vibration of the two-level
emitter on the energy spectrum and Rabi oscillation of cavity-QED setup. Our
model is similar to the trap ion system which is broadly studied to pursuit
its application in quantum information processing FM2018 ; FM2012 ; LD2018 ;
LD2019 . Instead, we here aim to find the exact energy diagram and study the
coherent dynamics of the system, in order to achieve a basic understand for
the model. With the assistance of a unitary transformation, we find the system
can be effectively described by an emitter-optomechanical cavity Hamiltonian
kippenberg2013 ; liu2014 , with a negligible Kerr term. We find that effective
Hamiltonian possesses a same mathematical structure with quantum Rabi model
DB2011 , and borrow the Bogliubov operators approach QC2012 ; QH2011 ; QH2012
to obtain the exact energy spectrum. It shows that, the optomechanical
interaction will induce a sideband effect, and in each of the sideband, we
observe the Rabi splitting which originates from the emitter-cavity coupling.
We also find that the effective optomechanical interaction leads to the two-
frequency Rabi oscillation and explain it in the dressed state presentation.

##  II Model and Hamiltonian

Figure 1: Schematic diagram of the model: a single mode cavity couples to a
two-level moving emitter which is subject to a harmonic potential.

As schematically shown in Fig. 1, the system we consider is composed by a
single-mode cavity and a movable but spatially confined two-level emitter. The
emitter is characterized by its mass $M$ and the internal energy level spacing
$\Omega$ between the ground state $|g\rangle$ and excited state $|e\rangle$.
We introduce a confinement of the emitter by a harmonic potential of the
oscillator frequency $\omega$. Considering that the spatial motion (vibration)
of the emitter is along the $x$ axis, which is perpendicular to the wall of
the cavity, the Hamiltonian is given by XG1995 ; LX2000

$\displaystyle H$ $\displaystyle=$
$\displaystyle\frac{p^{2}}{2M}+\frac{1}{2}M\omega^{2}x^{2}+\hbar\omega_{a}a^{\dagger}a+\hbar\Omega|e\rangle\langle
e|$ (1) $\displaystyle+\hbar g(a^{\dagger}\sigma_{-}e^{-ikx}+{\rm H.c.}).$

Here, $x$ and $p$ are the emitter’s position and momentum operators,
$\omega_{a}=ck$ is the cavity frequency with $k$ being the photon wave vector,
and $c$ being the velocity of light.
$\sigma_{-}=(\sigma_{+})^{\dagger}=|g\rangle\langle e|$ is the Pauli operator
of the emitter, $a\,({a}^{\dagger})$ is the annihilation (creation) operator
of the cavity field. $g$ is the coupling strength between the emitter and
cavity. In the Hamiltonian Eq. (1), we have applied the rotating wave
approximation by considering $g\ll\\{\omega_{a},\Omega\\}$. It is convenient
to introduce the creation (annihilation) operator $b^{\dagger}\,(b)$, which
satisfies $x=\alpha(b^{\dagger}+b),p=i\hbar(b^{\dagger}-b)/(2\alpha)$
($\alpha=\sqrt{\hbar/2M\omega}$) and the Hamiltonian can be rewritten as

$\displaystyle H$ $\displaystyle=$ $\displaystyle\hbar\omega
b^{\dagger}b+\hbar\omega_{a}a^{\dagger}a+\hbar\Omega|e\rangle\langle e|$ (2)
$\displaystyle+\hbar g[a^{\dagger}\sigma_{-}e^{-ik\alpha(b^{\dagger}+b)}+{\rm
H.c.}]$

by neglecting the constant term.

The operator in the exponential term can be eliminated by performing a unitary
transformation $\tilde{H}=UHU^{\dagger}$ where
$U=e^{ik\alpha(b^{\dagger}+b)a^{\dagger}a}$, and it yields
$\tilde{H}=\tilde{H}_{1}+\tilde{H}_{2}$ with

$\displaystyle\tilde{H}_{1}$ $\displaystyle=$
$\displaystyle\hbar\chi(a^{\dagger}a)^{2}+\hbar\omega_{a}a^{\dagger}a+\hbar\Omega|e\rangle\langle
e|$ (3) $\displaystyle+\hbar g(a^{\dagger}\sigma_{-}+a\sigma_{+})+\hbar\omega
b^{\dagger}b,$ $\displaystyle\tilde{H}_{2}$ $\displaystyle=$ $\displaystyle
i\hbar\eta a^{\dagger}a(b-b^{\dagger}),$ (4)

where

$\chi=k^{2}\alpha^{2}\omega,\,\eta=k\alpha\omega.$ (5)

It is obvious that, the vibrational movement of the emitter induces two
effects. The first one is the Kerr effect as shown by the first term of
$\tilde{H}_{1}$, with the strength $\chi=k^{2}\alpha^{2}\omega=\hbar
k^{2}/(2M)$, which is independent of oscillator frequency $\omega$. Physically
speaking, the movement of the emitter is described by the generation and
absorption of the phonon in second quantization representation, and it is also
accompanied by the generation and absorption of photon in the cavity as shown
by the emitter-photon interaction Hamiltonian $\hbar
g[a^{\dagger}\sigma_{-}e^{-ik\alpha(b^{\dagger}+b)}+{\rm H.c.}]$. Therefore,
it naturally introduces a self-phase modulation to the photon in the cavity.
The other one is the effective coupling between the vibrational degree of
freedom of the emitter and the cavity mode. As given by $\tilde{H}_{2}$, it is
actually an effective optomechanical interaction law1995 ; marquardt2014  with
strength $\eta=k\alpha\omega$, which depends both on the parameters of the
emitter and the harmonic potential. Followed by the typical cavity QED system
with Rydberg atom, we take $\Omega=\omega_{a}=10^{5}$ GHz, $\omega=1\,{\rm
GHz},k=10^{7}\,{\rm m^{-1}},M=10^{-27}\,{\rm kg},g=100\,{\rm MHz}$. Within
these parameters, we will have $\chi=0.05g$, $\eta=g/\sqrt{2}$. Therefore, the
strength of the Kerr effect is much weaker than that of the emitter-cavity
coupling, that is, $\chi\ll g$.

The above cavity QED model can be experimentally realized in the Rydberg atom
platform, in which the parameters can be achieved by
$\Omega=\omega_{a}=10^{5}$ GHz, $k=10^{7}\,{\rm m^{-1}},M=10^{-27}\,{\rm
kg},g=100\,{\rm MHz}$ Anderson2011 ; Tikman2016 ; Bounds2018  . Furhtermore,
the trap of the atom can be realized by the optical tweezers technologies and
the depth of the harmonica trap $\omega$ can be achieve by hundreds of MHz
LTC2012 ; ZY2013 ; LH2017 . Within these parameters, the strength of the Kerr
effect is much weaker than that of the emitter-cavity coupling, that is,
$\chi\ll g$.

Figure 2: The energy spectrum diagram for $m=0,1$. The solid lines are the
eigen states of $\tilde{H}_{1}$ and the dashed lines represent the energy-
level transition induced by $\tilde{H}_{2}$.

The Hamiltonian $\tilde{H}_{1}$ is completely solvable due to the conservation
of the excitation number. The eigen values are

$\displaystyle E_{\pm}^{(m,n)}$ $\displaystyle=$
$\displaystyle(m^{2}+m+\frac{1}{2})\hbar\chi+(m+\frac{1}{2})\hbar\omega_{a}+\frac{1}{2}\hbar\Omega+n\hbar\omega$
$\displaystyle\pm\hbar\sqrt{[(m+\frac{1}{2})\chi+\frac{1}{2}\omega_{a}-\frac{1}{2}\Omega]^{2}+(m+1)g^{2}}$

and the corresponding eigen wave function can be obtained as

$\displaystyle|\psi_{+}^{(m,n)}\rangle$ $\displaystyle=$
$\displaystyle\cos\frac{\theta}{2}|m,n,e\rangle+\sin\frac{\theta}{2}|m+1,n,g\rangle,$
(7) $\displaystyle|\psi_{-}^{(m,n)}\rangle$ $\displaystyle=$
$\displaystyle-\sin\frac{\theta}{2}|m,n,e\rangle+\cos\frac{\theta}{2}|m+1,n,g\rangle,$
(8)

where $\tan\theta=2\sqrt{m+1}g/[\Omega-\omega_{a}-(2m+1)\chi]$, and
$|m,n,\sigma\rangle:=|m\rangle_{c}\otimes|n\rangle_{v}\otimes|\sigma\rangle_{a}$
($|\sigma\rangle=|e\rangle,|g\rangle$) represents the state in which the
cavity mode (vibrate mode) is in the bosonic Fock state with $m(n)$
excitations while the emitter is in the state $|\sigma\rangle$. In Fig. 2, we
illustrate the energy diagram for $m=0,1$. Here, the black solid lines are the
eigenstates of $\tilde{H}_{1}$ and the blue dashed lines represent the energy
level transitions between $|\psi_{\pm}^{(m,n)}\rangle$ and
$|\psi_{\pm}^{(m,n\pm 1)}\rangle$, which are induced by $\tilde{H}_{2}$. It
seems that the whole Hamiltonian can only be solved by the perturbation theory
with the presence of $\tilde{H}_{2}$ induced transition. However, thanks to
the excitation number conservation for the internal degree of freedom for the
emitter and the photons in the cavity, that is,
$[a^{\dagger}a+|e\rangle\langle e|,H]=0$, the whole system is still fully
solvable, and the exact energy spectrum can be obtained as what we will
discuss in the follows.

##  III The solution of the Hamiltonian

Now, we derive the exact energy spectrum of the Hamiltonian $\tilde{H}$.
First, we introduce $\tilde{b}=ib$, $\tilde{b^{\dagger}}=-ib^{\dagger}$, then
the Hamiltonian $\tilde{H}$ becomes

$\displaystyle\tilde{H}=$
$\displaystyle\hbar\omega\tilde{b^{\dagger}}\tilde{b}+\hbar\omega_{a}a^{\dagger}a+\hbar\Omega|e\rangle\langle
e|+\hbar g(a^{\dagger}\sigma_{-}+a\sigma_{+})$ (9) $\displaystyle+\hbar
k\alpha\omega a^{\dagger}a(\tilde{b}+\tilde{b^{\dagger}})+\hbar
k^{2}\alpha^{2}\omega(a^{\dagger}a)^{2}.$

In what follows, we will still use the symbol $b$ to represent $\tilde{b}$ for
the sake of simplicity since it does not affect the final result. In the
cavity-emitter basis $\\{|m+1,g\rangle,\,|m,e\rangle\\}$, the Hamiltonian can
be expressed as

$\tilde{H}=\left(\begin{array}[]{cc}H_{11}&\hbar\sqrt{m+1}g\\\
\hbar\sqrt{m+1}g&H_{22}\end{array}\right),$ (10)

where

$\displaystyle H_{11}$ $\displaystyle=$ $\displaystyle\hbar\omega
b^{\dagger}b+(m+1)\hbar\omega_{a}$ $\displaystyle+(m+1)\hbar
k\alpha\omega(b+b^{\dagger})+(m+1)^{2}\hbar k^{2}\alpha^{2}\omega,$
$\displaystyle H_{22}$ $\displaystyle=$ $\displaystyle\hbar\omega
b^{\dagger}b+m\hbar\omega_{a}+\hbar\Omega$ (12) $\displaystyle+m\hbar
k\alpha\omega(b+b^{\dagger})+m^{2}\hbar k^{2}\alpha^{2}\omega.$

The Hamiltonian has the same mathematical structure with that of the quantum
Rabi model (see Eq. (2) in Ref. QC2012 ). It motivates us to apply the
Bogliubov operators approach to solve the eigen spectrum. The basic idea is
that we can introduce two Bogolibov transformations to diagonalize the
Hamiltonian $H_{11}$ and $H_{22}$ respectively, and therefore the wave
function of the whole Hamiltonian $\tilde{H}$ can be obtained two times. Since
they correspond to the same eigenvalue, they should be only different by a
complex constant, and then we can build the transcendental equation for the
eigen energy. Following the process as given in the appendix (the similar
calculation can also be found in Ref. QC2012 ), the transcendental equation is
obtained as $G_{m}(E)=0$, where

$\displaystyle G_{0}(E)$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}\left[\frac{g^{2}}{(-n\omega+\gamma+E/\hbar)(\gamma+E/\hbar)}-1\right]f_{n}(k\alpha)^{n},$

for $m=0$ and

$\displaystyle G_{m}(E)$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}e_{n}[k\alpha(m+1)]^{n}\sum_{n=0}^{\infty}e_{n}(k\alpha
m)^{n}$ (14)
$\displaystyle-\sum_{n=0}^{\infty}f_{n}[k\alpha(m+1)]^{n}\sum_{n=0}^{\infty}f_{n}(k\alpha
m)^{n},$

for $m>0$.

The coefficients $e_{n}$ and $f_{n}$ are defined recursively as

$\displaystyle e_{n}$ $\displaystyle=$
$\displaystyle\frac{-\sqrt{m+1}gf_{n}}{l\omega-\gamma-E/\hbar},$ (15)
$\displaystyle nf_{n}$ $\displaystyle=$ $\displaystyle
K_{n-1}f_{n-1}-f_{n-2},$ (16)

with the initial conditions $f_{0}=1$ , $f_{1}=K_{0}$, and

$K_{n}=\frac{1}{k\alpha\omega}[(n\omega+\beta-E/\hbar)-\frac{(m+1)g^{2}}{n\omega-\gamma-E/\hbar}].$
(17)

Here, $\gamma:=-(m+1)\omega_{a},\beta:=k^{2}\alpha^{2}\omega+(m+1)\omega_{a}$.

Figure 3: (a) $G_{0}$ and (b) $G_{1}$ for $g=100\,{\rm
MHz},\Omega=\omega_{a}=10^{5}\,{\rm GHz},k=2\pi/\lambda=10^{7}\,{\rm m^{-1}}$
and $\omega=1\,{\rm GHz}=10g$.

In Fig. 3 (a) and (b), we plot the $G_{m}(E)$ functions for $m=0$ and $m=1$ by
the blue curves, respectively. Meanwhile, the red curves demonstrate the
divergent behavior at $E/\hbar=n\omega+(m+1)\omega_{a}$, which is implied by
Eq. (17). Therefore, the zero points of the blue curves correspond to the
eigen-energy of the system. As shown in the figure, where we have set the
transition frequency being resonant with the cavity, that is
$\omega_{a}=\Omega$, we can clearly observe the sidebands near $n\hbar\omega$,
which is induced by the vibration of the emitter. Near each sideband, it shows
a Rabi splitting behavior which is given by the emitter-cavity coupling terms
$\hbar g(a^{\dagger}\sigma^{-}+a\sigma^{+})$ in the Hamiltonian.

Recalling that, in Fig. 2, we have plotted the eigenstates of $\tilde{H_{1}}$
by the black solid lines where the energy level spacing between the states
$|\psi_{\pm}^{(m,n)}\rangle$ and $|\psi_{\pm}^{(m,n\pm 1)}\rangle$ is

$\Delta_{m,n}=2\hbar\sqrt{[(m+\frac{1}{2})\chi]^{2}+(m+1)g^{2}}.$ (18)

For the parameter regime considered in Fig. 3, the energy level space achieves
$\Delta_{0,n}\approx 2.11\hbar g$ with $m=0$, which is similar to the space
$\tilde{\Delta}_{0,n}\approx 1.99\hbar g$ in Fig. 3 (a). The similar result
can be also obtained for $m=1$, the result obtained from Eq. (6) is close to
the exact solution given in Fig. 3 (b) as
$|\Delta_{1,n}-\tilde{\Delta}_{1,n}|\approx 0.03\hbar g$. Therefore, the
energy level transitions introduced by $\tilde{H_{2}}$ produce the slight
shift to the energy spectrum of the system.

##  IV The Rabi oscillation

From now on, we will numerically discuss the dynamical evolution of the
system, i.e., to study the Rabi oscillation behavior. Remember that the
effective Hamiltonian $\tilde{H}$ is obtained by a unitary transformation,
correspondingly, we also need to perform the same unitary transformation on
the quantum state. Therefore, preparing the initial pure state as
$|\psi(0)\rangle$, the dynamics of the system is governed by

$|\psi(t)\rangle=U^{\dagger}e^{-i\tilde{H}t}U|\psi(0)\rangle,$ (19)

and the average value for an arbitrary operator $\hat{A}$ reads

$\langle\hat{A}\rangle={\rm Tr}[\hat{A}\rho(t)],$ (20)

where the density matrix $\rho(t)=|\psi(t)\rangle\langle\psi(t)|$.

Figure 4: The Rabi oscillation of the system (a,c) and the corresponding
frequency spectrum (b,d). The initial state is set as
$|\psi(0)\rangle=|m+1,0,g\rangle$ with $m=0$ for (a) and $m=1$ for (c),
respectively. The parameters are set to be same as those in Fig. 3.

As is well known, for the traditional JC model, the emitter and the cavity
field will exchange the excitation, which leads to a perfect Rabi oscillation.
However, for a moving emitter, even when the vibrate mode is in the ground
state, the oscillation behavior is still changed dramatically. In Fig. 4, we
plot the average value $P=\langle\hat{A}\rangle$ with
$\hat{A}=|m,g\rangle\langle m,g|$ with the initial state of the system being
$|\psi(0)\rangle=|m+1,0,g\rangle$. In Fig. 4 (a) and (c), we illustrate the
results for $m=0$ and $m=1$, respectively.

As shown in the figure, for a deep harmonic potential $\omega=10g$, the blue
curves demonstrate perfect Rabi oscillations with periods $T=\pi/g$ and
$T=\pi/(\sqrt{2}g)$ respectively in Fig. 4 (a) and (c). In such a situation,
the emitter is confined tightly by the harmonic potential, and it is similar
to that in standard JC model with static atom. However, for the shallow
harmonic potential, as shown by the red dashed curves in Fig. (a) and (c), the
dynamics diverges from the standard Rabi oscillation, and it shows a two-
frequency oscillation character, which can be obtained by the numerical
Fourier transformation,

$f(\omega_{0})=\frac{1}{\sqrt{2\pi}}\int dtP(t)e^{-i\omega_{0}t}$ (21)

and the spectrum strength corresponds to (a) and (c) are given in (b) and (d),
respectively. Here, we clearly observe the spectrum splitting, which is
represented by the red dashed lines.

The splitting can be understood from the energy-level diagram in Fig. 2, which
shows the energy-level transition between different sidebands. Taking the
subspaces with $m=0$ and $n=0,1$ as an example, the states
$|\psi_{\pm}^{(0,0)}\rangle$ will couple to states
$|\psi_{\pm}^{(0,1)}\rangle$ simultaneously, that is, it forms four transition
channels as shown by the dashed lines. However, in the parameter regime we
consider, the coupling between $|\psi_{+}^{(0,0)}\rangle$ and
$|\psi_{-}^{(0,1)}\rangle$ will play the most important role due to their
smallest energy spacing. As a result, a simplified energy-level diagram can be
given by Fig. 5 (a). Neglecting the effective Kerr interaction, whose strength
$\chi$ is much smaller than that of the emitter-cavity coupling $g$, the
$|\psi_{+}^{(0,0)}\rangle\leftrightarrow|\psi_{-}^{(0,1)}\rangle$ transition
intensity $\mu$ is

$\displaystyle\hbar\mu\approx\langle\psi_{+}^{(0,0)}|\hbar\eta
a^{\dagger}a(b-b^{\dagger})|\psi_{-}^{(0,1)}\rangle=\frac{1}{2}\hbar\eta.$
(22)

with

$\displaystyle|\psi_{+}^{(0,0)}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}(|0,0,e\rangle+|1,0,g\rangle),$
$\displaystyle|\psi_{-}^{(0,1)}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}(-|0,1,e\rangle+|1,1,g\rangle),$ (23)

where we have considered $\omega_{a}=\Omega$. As a result, it forms another
two dressed states $|\Psi_{+}\rangle$ and $|\Psi_{-}\rangle$, which are the
superposition of $|\psi_{+}^{(0,0)}\rangle$ and $|\psi_{-}^{(0,1)}\rangle$ as
shown in Fig. 5 (b). Therefore, the Rabi oscillation process can be
approximately considered as the oscillation between the states
$|\psi_{-}^{(0,0)}\rangle$ and $|\Psi_{\pm}\rangle$, the corresponding
transition frequency $\omega_{\pm}$ can be obtained as what follows.

Figure 5: (a) The original simplified energy-level diagram. (b) The
interpretation for the two-frequency Rabi oscillation behavior.

Neglecting the Kerr interaction and considering the situation with
$\omega_{a}=\Omega$, we will have

$\displaystyle E_{+}^{(0,0)}$ $\displaystyle=$
$\displaystyle\hbar(\omega_{a}+g),$ (24a) $\displaystyle E_{-}^{(0,1)}$
$\displaystyle=$ $\displaystyle\hbar(\omega_{a}+\omega-g),$ (24b)
$\displaystyle E_{-}^{(0,0)}$ $\displaystyle=$
$\displaystyle\hbar(\omega_{a}-g).$ (24c)

Therefore, the eigen energies of $|\Psi_{\pm}\rangle$ are obtained as

$\displaystyle E_{\pm}$ $\displaystyle=$
$\displaystyle\hbar\omega_{a}+\frac{1}{2}\hbar\omega\pm\hbar\sqrt{(g-\frac{1}{2}\omega)^{2}+\mu^{2}}.$
(25)

As a result, the energy level transition frequencies of the system shown in
Fig. 5(b) are

$\displaystyle\hbar\omega_{\pm}$ $\displaystyle=$ $\displaystyle
E_{\pm}-E_{-}^{(0,0)}$ (26) $\displaystyle=$ $\displaystyle\hbar
g+\frac{1}{2}\hbar\omega\pm\hbar\sqrt{(g-\frac{1}{2}\omega)^{2}+\mu^{2}}.$

For the considered parameters $\omega=2g$ in Fig. 4 (a) and (b), the coupling
strength achieves $\mu\approx 0.15g$, and $\omega_{\pm}\approx(2\pm 0.15)g$,
which coincides with the two peaks in Fig. 4 (b) (see the red dashed curve).
The similar results can also be obtained for $m=1$, and the two-peak structure
for the spectrum strength which is given by the red dashed curve in Fig. 4 (d)
can be predicted and the transition frequencies are approximately obtained as
$\omega_{\pm}^{\prime}\approx(2.9\pm 0.23)g$.

##  V Conclusion

In this paper, we investigate the energy spectrum and the Rabi oscillation
behavior in a light-matter interaction model with a moving emitter. We
introduce a harmonic potential to confine the vibration degree of the emitter,
and show that the vibration of the emitter will induce an effective Kerr
interaction and opto-mechanical coupling. With the assistance of Bogliubov
operators approach, we obtain the exact energy spectrum of the system.
Furthermore, with a shallow potential, we find that the Rabi oscillation will
exhibit a two-frequency character, which is dramatically different from that
of a static emitter.

In the previous studies, it was shown that the mechanical squeezing can be
realized in the optomechanical system wollman2015 . Therefore, we hope our
study about the vibration induced optomechanical interaction can be applied in
the squeezed state preparation and is furthermore beneficial for quantum
precision measurement and sensing.

Note- During the preparation of this work, we find a similar investigation
about the optomechanical strong coupling between a single cavity photon and a
single atom chang2021 .

###### Acknowledgements.

We thank Profs. X. X. Yi, Y. Li and L.-P. Yang for helpful discussion. This
work is supported by the funding from Ministry of Science and Technology of
China (No. 2021YFE0193500) and the Natural Science Foundation of China (Nos.
11875011 and 12047566).

##  Appendix A Exact solution

In this appendix, we give the detailed derivation of the $G$-function, whose
zero points yield the eigen energy of the system. The same approach to deal
with quantum Rabi model can be found in Ref. QC2012 .

Based on Eq. (10) in the main text, we introduce the Bogoliubov operators,

$b^{\dagger}=B^{\dagger}-k\alpha(m+1),\,b=B-k\alpha(m+1),$ (27)

to generate the new bosonic operators $B$ and $B^{\dagger}$. Thus, we can
remove the linear term of the diagonal elements of the Hamiltonian matrix and
simplify it to

$\displaystyle\tilde{H}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cc}\hbar\omega
B^{\dagger}B-\hbar\gamma&\hbar\sqrt{m+1}g\\\ \hbar\sqrt{m+1}g&\hbar\omega
B^{\dagger}B-\hbar k\alpha\omega(B^{\dagger}+B)+\hbar\beta\end{array}\right).$
(30)

where

$\gamma=-(m+1)\omega_{a},\beta=k^{2}\alpha^{2}\omega+m\omega_{a}+\Omega.$

The wave function can then be assumed as

$\displaystyle|\Phi\rangle_{B}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}\sum_{n=0}^{\infty}\sqrt{n!}e_{n}|n\rangle_{B}\\\
\sum_{n=0}^{\infty}\sqrt{n!}f_{n}|n\rangle_{B}\end{array}\right),$ (34)

where $e_{n}$ and $f_{n}$ are the expansion coefficients. $|n\rangle_{B}$ is
an extended coherent state. It has the following properties

$\displaystyle|n\rangle_{B}$ $\displaystyle=$
$\displaystyle\frac{(B^{\dagger})^{n}}{\sqrt{n!}}|0\rangle_{B}=\frac{(b^{\dagger}+k\alpha(m+1))^{n}}{\sqrt{n!}}|0\rangle_{B},$
(35) $\displaystyle|0\rangle_{B}$ $\displaystyle=$ $\displaystyle
e^{-\frac{1}{2}k^{2}\alpha^{2}(m+1)^{2}-k\alpha(m+1)b^{\dagger}}|0\rangle_{b}.$
(36)

Here the vacuum state represented by the Bogoliubov operator is defined as the
eigenstate of the annihilation operator $b$.

By the Schrödinger equation, we will have

$\displaystyle\sum\limits_{n=0}^{\infty}\hbar(n\omega-\gamma)\sqrt{n!}e_{n}|n\rangle_{B}+\hbar\sqrt{m+1}g\sum\limits_{n=0}^{\infty}\sqrt{n!}f_{n}|n\rangle_{B}$
$\displaystyle=E\sum\limits_{n=0}^{\infty}\sqrt{n!}e_{n}|n\rangle_{B},$
$\displaystyle\hbar\sqrt{m+1}g\sum\limits_{n=0}^{\infty}\sqrt{n!}e_{n}|n\rangle_{B}+\sum\limits_{n=0}^{\infty}\hbar(n\omega+\beta)\sqrt{n!}f_{n}|n\rangle_{B}$
$\displaystyle-\hbar
k\alpha\omega\sum\limits_{n=0}^{\infty}(\sqrt{n}f_{n}\sqrt{n!}|n-1\rangle_{B}+\sqrt{n+1}f_{n}\sqrt{n!}|n+1\rangle_{B})$
$\displaystyle=E\sum\limits_{n=0}^{\infty}\sqrt{n!}f_{n}|n\rangle_{B}.$

Left-multiplying both sides of the above equation by ${}_{B}\langle l|$ gives

$(l\omega-\gamma-E/\hbar)e_{l}=-\sqrt{m+1}gf_{l},$ (38)
$(l\omega+\beta-E/\hbar)f_{l}-k\alpha\omega(l+1)f_{l+1}-k\alpha\omega=-\sqrt{m+1}ge_{l}.$
(39)

The coefficients $e_{l}$ and $f_{l}$ have the following relationship

$\displaystyle e_{l}$ $\displaystyle=$
$\displaystyle\frac{-\sqrt{m+1}gf_{l}}{l\omega-\gamma-E/\hbar},$ (40)
$\displaystyle lf_{l}$ $\displaystyle=$ $\displaystyle
K_{l-1}f_{l-1}-f_{l-2},$ (41)

where

$K(l)=\frac{1}{k\alpha\omega}[(l\omega+\beta-E/\hbar)-\frac{(m+1)g^{2}}{l\omega-\gamma-E/\hbar}].$
(42)

with $f_{0}=1$ and $f_{1}=K_{0}$.

Similarly, we can define another Bogoliubov operator A
($b^{\dagger}=A^{\dagger}-k\alpha m$). The transformed Hamiltonian then reads

$\displaystyle\tilde{H}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cc}\hbar\omega A^{\dagger}A+\hbar
k\alpha\omega(A^{\dagger}+A)+\hbar\beta^{\prime}&\hbar\sqrt{m+1}g\\\
\hbar\sqrt{m+1}g&\hbar\omega
A^{\dagger}A-\hbar\gamma^{\prime}\end{array}\right),$ (45)

where

$\gamma^{\prime}=-m\omega_{a}-\Omega,\
\beta^{\prime}=k^{2}\alpha^{2}\omega+(m+1)\omega_{a}.$

The wave function can also be written as

$\displaystyle|\Phi\rangle_{A}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{c}(-1)^{n}\sum_{n=0}^{\infty}\sqrt{n!}f_{n}^{\prime}|n\rangle_{A}\\\
(-1)^{n}\sum_{n=0}^{\infty}\sqrt{n!}e_{n}^{\prime}|n\rangle_{A}\end{array}\right),$
(49)

and they obey the properties

$\displaystyle|n\rangle_{A}$ $\displaystyle=$
$\displaystyle\frac{(A^{\dagger})^{n}}{\sqrt{n!}}|0\rangle_{A}=\frac{(b^{\dagger}+k\alpha
m)^{n}}{\sqrt{n!}}|0\rangle_{A},$ (50) $\displaystyle|0\rangle_{A}$
$\displaystyle=$ $\displaystyle e^{-\frac{1}{2}k^{2}\alpha^{2}m^{2}-k\alpha
mb^{\dagger}}|0\rangle_{b}.$ (51)

Following the previous steps, the relationship between the two coefficients
can be obtained as

$\displaystyle
e_{l}^{\prime}=\frac{-\sqrt{m+1}gf_{l}^{\prime}}{l\omega-\gamma^{\prime}-E/\hbar},$
(52)

The corresponding recursive relationship is

$\displaystyle lf_{l}^{\prime}$ $\displaystyle=$ $\displaystyle
K_{l-1}^{\prime}f_{l-1}^{\prime}-f_{l-2}^{\prime},$ (53) $\displaystyle
K^{\prime}(l)$ $\displaystyle=$
$\displaystyle\frac{1}{k\alpha\omega}[(l\omega+\beta^{\prime}-E/\hbar)-\frac{(m+1)g^{2}}{l\omega-\gamma^{\prime}-E/\hbar}],$

with $f_{0}^{\prime}=1$ and $f_{1}^{\prime}=K_{0}^{\prime}$.

Since both the wave functions (34) and (49) are the true eigenfunction for a
nondegenerate eigenvalue $E$, they should be proportional to each other, that
is, $|\Phi\rangle_{B}=r|\Phi\rangle_{A}$, where $r$ is a complex constant.
Projecting both sides of this identity onto the original vacuum state
${}_{b}\langle 0|$ , we have

$\displaystyle\sum_{n=0}^{\infty}\sqrt{n!}e_{n}{}_{b}\langle 0|n\rangle_{B}$
$\displaystyle=$ $\displaystyle
r(-1)^{n}\sum_{n=0}^{\infty}\sqrt{n!}f_{n}^{\prime}{}_{b}\langle
0|n\rangle_{A},$ $\displaystyle\sum_{n=0}^{\infty}\sqrt{n!}f_{n}{}_{b}\langle
0|n\rangle_{B}$ $\displaystyle=$ $\displaystyle
r(-1)^{n}\sum_{n=0}^{\infty}\sqrt{n!}e_{n}^{\prime}{}_{b}\langle
0|n\rangle_{A},$

and from (36) and (51), we obtain

$\displaystyle\sqrt{n!}_{b}\langle 0|n\rangle_{B}$ $\displaystyle=$
$\displaystyle(k\alpha(m+1))^{n}e^{-\frac{1}{2}k^{2}\alpha^{2}(m+1)^{2}},$
$\displaystyle(-1)^{n}\sqrt{n!}_{b}\langle 0|n\rangle_{A}$ $\displaystyle=$
$\displaystyle(-k\alpha m)^{n}e^{-\frac{1}{2}k^{2}\alpha^{2}m^{2}}.$ (56)

Then we have to consider the situations with $m=0$ and $m\neq 0$,
respectively.

When $m=0$, eliminating the ratio constant $r$ gives

$\sum_{n=0}^{\infty}e_{n}(k\alpha)^{n}\sum_{n=0}^{\infty}e_{n}^{\prime}0^{n}=\sum_{n=0}^{\infty}f_{n}(k\alpha)^{n}\sum_{n=0}^{\infty}f_{n}^{\prime}0^{n},$
(57)

which yields

$\sum_{n=0}^{\infty}e_{n}(k\alpha)^{n}e_{0}^{\prime}=\sum_{n=0}^{\infty}f_{n}(k\alpha)^{n}f_{0}^{\prime}.$
(58)

With (40) and (52), we get

$\sum_{n=0}^{\infty}\frac{-gf_{n}}{n\omega-\gamma-E/\hbar}(k\alpha)^{n}\frac{gf_{0}^{\prime}}{\gamma^{\prime}+E/\hbar}-\sum_{n=0}^{\infty}f_{n}(k\alpha)^{n}f_{0}^{\prime}=0.$
(59)

Setting $\Omega=\omega_{a}$, we obtain the transcendental equation for the
eigen energy $E$ as

$\displaystyle G_{0}(E)$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}\frac{g^{2}f_{n}}{(-n\omega+\gamma+E/\hbar)(\gamma+E/\hbar)}(k\alpha)^{n}$
(60) $\displaystyle-\sum_{n=0}^{\infty}f_{n}(k\alpha)^{n}=0.$

For $m\neq 0$, we will similarly reach

$\displaystyle G_{m}(E)$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}e_{n}[k\alpha(m+1)]^{n}\sum_{n=0}^{\infty}e_{n}(k\alpha
m)^{n}$
$\displaystyle-\sum_{n=0}^{\infty}f_{n}[k\alpha(m+1)]^{n}\sum_{n=0}^{\infty}f_{n}(k\alpha
m)^{n}=0.$

which are Eq. (LABEL:Gfunction0) and Eq. (14) in the main text for $m=0$ and
$m\neq 0$, respectvely.

## References

  * (1) E. M. Purcell, Phys. Rev. 69, 681 (1946).
  * (2) E. T. Jaynes and F. W. Cummings, Proc. IEEE 51, 89 (1963).
  * (3) G. S. Agarwal, J. Opt. Soc. Am. B 2, 480 (1985)
  * (4) M. Brune, F. S.-Kaler, A. Maali, J. Dreyer, E. Hagley, J. M. Raimond, and S. Haroche, Phys. Rev. Lett. 76, 1800 (1996).
  * (5) L. Garziano, R. Stassi, V. Macrí, A. F. Kockum, S. Savasta and F. Nori, Phys. Rev. A 92, 063830 (2015).
  * (6) G. Calajó and P. Rabl, Phys. Rev. A 95, 043824 (2017).
  * (7) E. S.-Burillo, A. G.-Tudela, and C. G.-Ballestero, Phys. Rev. A 102, 013726 (2020).
  * (8) Q. Li, D. Z. Xu, C. Y. Cai, and C. P. Sun, Sci. Rep. 3, 3144 (2013).
  * (9) D. Braun and J. Martin, Phys. Rev. A 77, 032102 (2008).
  * (10) F. Damanet, D. Braun, and J. Martin, Phys. Rev. A 93, 022124 (2016).
  * (11) X. G. Wang and C. P. Sun, J. Mod. Optics 42, 515 (1995).
  * (12) L. X. Cen and S. J. Wang, J. Phys. A: Math. Gen. 33, 3697 (2000).
  * (13) L. Zheng, C. Li, Y. Li, and C. P. Sun, Phys. Rev. A 71, 062101 (2005).
  * (14) L. Zheng, C. P. Yang, and F. Nori, Phys. Rev. A 82, 062106 (2010).
  * (15) L. You, Phys. Rev. A 64, 012302 (2001).
  * (16) Y. G. Deng, T. Shi, and S. Yi, Photon. Res. 9, 1289 (2021).
  * (17) C. S. Muñoz, E. d. Valle, A. G. Tudela, K. Müller, S. Lichtmannecker, M. Kaniber, C. Tejedor, J. J. Finley, and F. P. Laussy, Nat. Photon. 8, 550 (2014).
  * (18) Y. Zhang, J. Zhang, S. X. Wu, and C. S. Yu, Ann. Phys. 361, 563 (2015).
  * (19) K. M. Birnbaum, A. Boca, R. Miller, A. D. Boozer, T. E. Northup, and H. J. Kimble, Nature 436, 87 (2005).
  * (20) A. Agustí, L. G. Álvarez, E. Solano, and C. Sabín, Phys. Rev. A 103, 062201 (2021).
  * (21) S. Scheel and S. Y. Buhmann, Phys. Rev. A 80, 042902 (2009).
  * (22) O. D. Stefano, A. Settineri, V. Macrí, A. Ridolfo, R. Stassi, A. F. Kockum, S. Savasta, and F. Nori, Phys. Rev. Lett. 112, 030402 (2019).
  * (23) H. Wang, M. P. Blencowe, C. M. Wilson, and A. J. Rimberg, Phys. Rev. A 99, 053833 (2019).
  * (24) W. Qin, V. Macrí, A. Miranowicz, S. Savasta, and F. Nori, Phys. Rev. A 100, 062501 (2019).
  * (25) V. Macrí, A. Ridolfo, O. D. Stefano, A. F. Kockum, F. Nori, and S. Savasta, Phys. Rev. X 8, 011031 (2018).
  * (26) S. Felicetti, C. Sabín, I. Fuentes, L. Lamata, G. Romero, and E. Solano, Phys. Rev. B 92, 064501 (2015).
  * (27) S. E. Anderson, K. C. Younge, and G. Raithel, Phys. Rev. Lett. 107, 263001 (2011).
  * (28) Y. Tikman, I. Yavuz, M. F. Ciappina, A. Chacón, Z. Altun, and M. Lewenstein, Phys. Rev. A 93, 023410 (2016).
  * (29) A. D. Bounds, N. C. Jackson, R. K. Hanley, R. Faoro, E. M. Bridge, P. Huillery, and M. P. A. Jones, Phys. Rev. Lett. 120, 183401 (2018).
  * (30) F. Zhou, L. Yan, S. Gong, Z. Ma, J. He, T. Xiong, L. Chen, W. Yang, M. Feng and V. Vedral, Sci. Adv., 2, e1600578 (2018).
  * (31) L. Chen, W. Wan, Y. Xie, F. Zhou and M. Feng, Chin. Phys. Lett., 29, 033701 (2012).
  * (32) C. J. Trout, M. Li, M. Gutiérrez, Y. Wu, S.-T. Wang, L. Duan and K. R. Brown, New J. Phys., 20, 043038 (2018).
  * (33) K. A. Landsman, Y. Wu, P. H. Leung, D. Zhu, N. M. Linke, K. R. Brown, L. Duan and C. Monroe, Phys. Rev. A, 100, 022332 (2019).
  * (34) T. Ramos, V. Sudhir, K. Stannigel, P. Zoller, and T. J. Kippenberg, Phys. Rev. Lett. 110, 193602 (2013).
  * (35) H. Wang, X. Gu, Y.-x. Liu, A. Miranowicz, and F. Nori, Phys. Rev. A 90, 023817 (2014).
  * (36) D. Braak, Phys. Rev. Lett. 107, 100401 (2011).
  * (37) Q. H. Chen, C. Wang, S. He, T. Liu, and K. Wang, Phys. Rev. A 86, 023822 (2012).
  * (38) Q. H. Chen, T. Liu, Y. Y. Zhang, and K. L. Wang, Europhys. Lett. 96, 14003 (2011).
  * (39) Q. H. Chen, L. Li, T. Liu, and K. L. Wang, Chin. Phys. Lett. 29, 014208 (2012).
  * (40) T. Li, Z.-X. Gong, Z.-Q. Yin, H. T. Quan, X. Yin, P. Zhang, L.-M. Duan and X. Zhang, Phys. Rev. Lett., 109, 163001 (2012).
  * (41) Z.-Q. Yin, T. Li, X. Zhang and L.M. Duan, Phys. Rev. A, 88, 033614 (2013) .
  * (42) H.-K. Li, E. Urban, C. Noel, A. Chuang, Y. Xia, A. Ransford, B. Hemmerling, Y. Wang, T. Li, H. Häffner and X. Zhang, Phys. Rev. Lett., 118, 053001 (2017).
  * (43) C. K. Law, Phys. Rev. A 51, 2537 (1995).
  * (44) M. Aspelmeyer, T. J. Kippenberg, and F. Marquardt, Rev. Mod. Phys. 86, 1391 (2014).
  * (45) E. E. Wollman, C. U. Lei, A. J. Weinstein, J. Suh, A. Kronwald, F. Marquardt, A. A. Clerk, and K. C. Schwab, Science 349, 952 (2015).
  * (46) J. A.-Luengo and D. E. Chang, arXiv: 2108.03526 (2021).