#  Exploring properties of long-lived particles in
inelastic dark matter models at Belle II

Dong Woo Kang1, P. Ko1, Chih-Ting Lu1 1 School of Physics, KIAS, Seoul
130-722, Republic of Korea

###### Abstract

The inelastic dark matter model is one kind of popular models for the light
dark matter (DM) below $O(1)$ GeV. If the mass splitting between DM excited
and ground states is small enough, the co-annihilation becomes the dominant
channel for thermal relic density and the DM excited state can be long-lived
at the collider scale. We study scalar and fermion inelastic dark matter
models for ${\cal O}(1)$ GeV DM at Belle II with $U(1)_{D}$ dark gauge
symmetry broken into its $Z_{2}$ subgroup. We focus on dilepton displaced
vertex signatures from decays of the DM excited state. With the help of
precise displaced vertex detection ability at Belle II, we can explore the DM
spin, mass and mass splitting between DM excited and ground states.
Especially, we show scalar and fermion DM candidates can be discriminated and
the mass and mass splitting of DM sector can be determined within the
percentage of deviation for some benchmark points. Furthermore, the allowed
parameter space to explain the excess of muon $(g-2)_{\mu}$ is also studied
and it can be covered in our displaced vertex analysis during the early stage
of Belle II experiment.

††preprint: KIAS-P21001

##  I Introduction

The nature of dark matter (DM) in our Universe is still a great mysterious
issue. As we know, DM plays a major role in the structure formation
Frenk:2012ph , and its abundance is about 5.5 times larger than the ordinary
matter in the present Universe Aghanim:2018eyx . However, even all robust
evidences to support the existence of DM until now are only connected with
gravitational interactions. Still it is believed that there are non-negligible
couplings between Standard Model (SM) particles and DM in addition to
gravitational interactions. Especially, the Weakly Interacting Massive
Particle (WIMP) DM candidate with the freeze-out mechanism has been
overwhelming in both theoretical and experimental communities during the past
decades as shown in Ref. Plehn:2017fdg  and references therein. This kind of
WIMP DM particles can be searched with direct, indirect detections and also at
colliders experiments. The null signal result from direct detection provides
severe bounds on the cross section of DM and nuclear scattering above a few
GeV Aprile:2018dbl ; Ren:2018gyx ; Aprile:2019dbj . For instance, the spin
independent (spin dependent) DM and nuclear scattering cross section is
restricted to $\sigma_{SI}\lesssim 10^{-46}~{}{\rm cm}^{2}$ Aprile:2018dbl
($\sigma_{SD}\lesssim 2\times 10^{-41}~{}{\rm cm}^{2}$ Aprile:2019dbj ) for
$M_{DM}=100$ GeV. Nevertheless, lighter DM candidates are still less
constrained because the restrictions of fine energy threshold are required
from these direct detection experiments. Therefore, those dark sector models
with MeV to GeV DM candidates become more and more popular for
phenomenological studies and new experimental searches Knapen:2017xzo ;
Lin:2019uvt .

For the GeV scale DM searches, the high-intensity machines such as the BESIII
Ablikim:2009aa , BaBar TheBABAR:2013jta , Belle II Adachi:2018qme ;
Kou:2018nap  and also the fixed target experiments Bjorken:2009mm  can be more
powerful than the high-energy machines such as the Tevatron and LHC
Boveia:2018yeb . On the other hand, DM is not the lonely particle in the dark
sector for most of DM portal models Pospelov:2007mp ; ArkaniHamed:2008qn ;
Pospelov:2008zw ; Baek:2013qwa . In order to connect the SM sector with the
dark sector, a mediator is required. Among these mediator candidates, the dark
photon via the kinetic mixing portal is an attractive type Holdom:1985ag  for
the dark sector with Abelian gauge symmetry. In order to be consistent with
the relic density bound, it’s natural for both DM and dark photon in a similar
energy scale. Therefore, the searches for them at high-intensity machines are
in full swing and relevant constraints can be set up as shown in Ref.
Izaguirre:2015zva ; Duerr:2019dmv . Even though the original purpose for these
high-intensity machines is to study the properties of $J/\psi$ and $B$ mesons,
we can also apply them to explore some models with light DM candidates via the
mono-photon signature Zhang:2019wnz ; Lees:2017lec ; Izaguirre:2015zva ;
Duerr:2019dmv . Furthermore, as pointed out in Ref. Izaguirre:2015zva ;
Mohlabeng:2019vrz ; Duerr:2019dmv , we can not only study the mono-photon
signature, but also the displaced DM signature at B-factories for the
inelastic DM models. Especially, searches for displaced DM signature can cover
some parameter space which the invisible DM signature cannot reach at high-
intensity machines. Hence, we will focus on exploring the properties of long-
lived particles in dark sectors with $U(1)_{D}$ gauge symmetry at Belle II in
this work.

The inelastic (or excited) DM models with extra $U(1)_{D}$ gauge symmetry
TuckerSmith:2001hy ; Baek:2014kna ; Ko:2019wxq ; Baek:2020owl  is one of the
most popular dark sector models with light DM candidates111 In Ref.
TuckerSmith:2001hy , the $U(1)_{D}$ symmetry is explicitly and softly broken
by a dim–2 operator for scalar DM, and there is no dark Higgs boson there. On
the other hand, in the models considered in this paper in that there is dark
Higgs boson which plays important roles in DM phenomenology.. There are at
least two states in the dark sectors and there is an inelastic transition
between them via the new $U(1)_{D}$ gauge boson. If the mass splitting between
these two states is small enough, the co-annihilation channel could be the
dominate one of DM relic density in early Universe Griest:1990kh . It is one
of the unique features of this sort of models. The co-annihilation production
for light DM via thermal freeze out is still consistent with the Cosmic
Microwave Background (CMB) constraint for the amount of parameter space
Baek:2020owl ; Slatyer:2015jla . On the other hand, the constraint from DM and
nuclear inelastic scattering is much weaker than the elastic one in the direct
detection experiments222 The bounds from direct detection experiments are weak
for $\sim O(1)$ GeV DM. Hence we can safely ignore possible contributions to
the DM-nucleon elastic scattering from the SM Higgs boson and dark Higgs boson
in the direct detection experiments.. It makes more allowed parameter space
can be explored in the inelastic DM models. Besides, there are also other rich
phenomenons in this model. For example, it is possible to explain the muon
$g-2$ anomaly Mohlabeng:2019vrz  and XENON1T excess Harigaya:2020ckz ;
Bell:2020bes ; Lee:2020wmh ; Baryakhtar:2020rwy ; Bramante:2020zos ;
An:2020tcg ; Baek:2020owl . Last but not least, the excited DM state can
naturally become long-lived and leave displaced vertex inside detectors after
it has been produced such that these novel signatures can be searched for at
colliders.

The fundamental properties of DM particle include its mass, charge, and spin.
Since we assume DM is electrically neutral or millicharged, only its mass and
spin are left to be determined. Actually, the mass range for WIMP DM can cover
from few MeV to a few hundred TeV and its spin can be $0,1/2,1,3/2$ or even
higher. Most of the time the precise determination of DM spin and mass are
challenging tasks. However, there are already some previous studies about DM
spin and mass determinations. For the DM spin determination, it has been
studied for direct detection experiments in Ref. Catena:2017wzu  and collider
experiments in Ref. Cheng:2010yy ; Melia:2011cu ; Edelhauser:2012xb ;
Christensen:2013sea ; Choi:2018sqc ; Abdallah:2019tpo . For the DM mass
determination, it has been studied for direct detection in Ref. Shan:2009ym ;
Kavanagh:2013wba  and collider experiments in Ref. Cheng:2011ya ;
HarlandLang:2011ih ; HarlandLang:2012gn ; Kobach:2013tga ; Harland-
Lang:2013wxa ; Christensen:2014yya ; Konar:2015hea ; Ko:2016xwd ;
Xiang:2016jni ; Kamon:2017yfx ; Dutta:2017sod ; Kang:2019ukr ; Kim:2019prx ;
Banerjee:2019ktv ; Bae:2020dwf . Especially, using the displaced vertex
information for DM mass determination starts from Ref. Kang:2019ukr . Thanks
to the precise track resolution in the inner detector of ATLAS and CMS
experiments, the DM mass reconstruction can be studied in some BSM models of
the specific cascade processes involving long-lived particles. Furthermore,
the resolutions of displaced vertex, charged lepton and photon momentum can be
even improved at detectors of Belle II experiments Adachi:2018qme . Therefore,
our goal in this paper is to explore the spin and mass properties of inelastic
DM models with dilepton displaced vertex signatures at Belle II. As we will
see in Sec. IV, the scalar and fermion inelastic DM models can be well
discriminated and the DM mass and mass splitting between DM excited and ground
states can be determined within the percentage of deviation for those
benchmark points.

The organization of this paper is as follows. We first review scalar and
fermion inelastic DM models with $U(1)_{D}$ gauge symmetry in Sec. II. Both
analytical representations and numerical results for cross sections of
relevant signal processes are displayed in Sec. III. We also point out how to
distinguish scalar and fermion inelastic DM models via the size and kinematic
distribution of their cross sections in this section. Detailed simulations and
methods to determine DM mass, and mass splitting between DM excited and ground
states are shown in Sec. IV. Finally, we conclude our studies in Sec. V. Some
supplemental formulae for Sec. III can be found in the Appendix A and B.

##  II Inelastic dark matter models

The scalar and fermion inelastic (or excited) DM models with $U(1)_{D}$ gauge
symmetry are reviewed in this section. After the spontaneous symmetry breaking
(SSB) of this $U(1)_{D}$ gauge symmetry, we expect the accidentally residual
$Z_{2}$ symmetry, $\phi_{1}\rightarrow-\phi_{1}$ (scalar) or
$\chi_{1}\rightarrow-\chi_{1}$ (fermion), can be left such that $\phi_{1}$ or
$\chi_{1}$ are stable and become DM candidates in our Universe.

###  II.1 The scalar model Baek:2014kna ; Baek:2020owl

We consider a dark sector with two singlet complex scalars $\Phi$ and
$\phi=(\phi_{2}+i\phi_{1})/\sqrt{2}$ as the dark Higgs and dark matter
sectors, respectively. Both $\Phi$ and $\phi$ are charged under $U(1)_{D}$,
but neutral of the SM gauge symmetry. The $U(1)_{D}$ charges for them are
assigned as $Q_{D}(\Phi)=+2$ and $Q_{D}(\phi)=+1$. Besides, the SM-like Higgs
doublet and other SM particles do not carry $U(1)_{D}$ charges.

The scalar part of the renormalizable and gauge invariant Lagrangian density
is

${\cal
L}_{scalar}=|D_{\mu}H|^{2}+|D_{\mu}\Phi|^{2}+|D_{\mu}\phi|^{2}-V(H,\Phi,\phi),$
(1)

with

$D_{\mu}H=(\partial_{\mu}+i\frac{g}{2}\sigma_{a}W^{a}_{\mu}+i\frac{g^{\prime}}{2}B_{\mu})H,$
(2) $D_{\mu}\Phi=(\partial_{\mu}+ig_{D}Q_{D}(\Phi)X_{\mu})\Phi,$ (3)
$D_{\mu}\phi=(\partial_{\mu}+ig_{D}Q_{D}(\phi)X_{\mu})\phi,$ (4)

where $W^{a}_{\mu}$, $B_{\mu}$, and $X_{\mu}$ are the gauge potentials of the
$SU(2)_{L}$, $U(1)_{Y}$ and $U(1)_{D}$ with gauge couplings $g$, $g^{\prime}$
and $g_{D}$, respectively. The $\sigma_{a}$ is the Pauli matrix and $a$ runs
from 1 to 3. The scalar potential in Eq.(1) is given by

$\displaystyle V(H,\Phi,\phi)=$
$\displaystyle-\mu^{2}_{H}H^{\dagger}H+\lambda_{H}(H^{\dagger}H)^{2}-\mu^{2}_{\Phi}\Phi^{\ast}\Phi+\lambda_{\Phi}(\Phi^{\ast}\Phi)^{2}$
$\displaystyle-\mu^{2}_{\phi}\phi^{\ast}\phi+\lambda_{\phi}(\phi^{\ast}\phi)^{2}+(\mu_{\Phi\phi}\Phi^{\ast}\phi^{2}+H.c.)$
$\displaystyle+\lambda_{H\Phi}(H^{\dagger}H)(\Phi^{\ast}\Phi)+\lambda_{H\phi}(H^{\dagger}H)(\phi^{\ast}\phi)+\lambda_{\Phi\phi}(\Phi^{\ast}\Phi)(\phi^{\ast}\phi),$
(5)

where all parameters are assumed to be real for simplicity.

We then expand $H,\Phi$ fields around the vacuum with the unitary gauge,

$H(x)=\frac{1}{\sqrt{2}}\left(\begin{tabular}[]{c}0\\\
$v+h(x)$\end{tabular}\right)\;\;\;,\;\;\;\Phi(x)=\frac{1}{\sqrt{2}}\left(v_{D}+h_{D}(x)\right),$
(6)

and highlight some important relations in the below :

  * •

The four-points and off-diagonal interactions of the DM sector and new gauge
boson can be obtained from the $|D_{\mu}\phi|^{2}$ term:

$|D_{\mu}\phi|^{2}=(\partial_{\mu}\phi_{2})^{2}+(\partial_{\mu}\phi_{1})^{2}+g^{2}_{D}X_{\mu}X^{\mu}(\phi^{2}_{2}+\phi^{2}_{1})+g_{D}X_{\mu}(\phi_{2}\partial^{\mu}\phi_{1}-\phi_{1}\partial^{\mu}\phi_{2})$
(7)

  * •

The DM sector mass splitting and dark Higgs, DM sector trilinear interactions
are derived from $\mu_{\Phi\phi}\Phi^{\ast}\phi^{2}+H.c.$ terms:

$\mu_{\Phi\phi}\Phi^{\ast}\phi^{2}+H.c.=\mu_{\Phi\phi}(v_{D}+h_{D}(x))(\phi^{2}_{2}-\phi^{2}_{1}).$
(8)

  * •

Finally, the $\phi_{1}$, $\phi_{2}$ masses and their mass splitting can be
represented as

$M_{\phi_{1,2}}=\sqrt{\frac{1}{2}(-\mu^{2}_{\phi}+\lambda_{H\phi}v^{2}+\lambda_{\Phi\phi}v^{2}_{D})\mp\mu_{\Phi\phi}v_{D}},$
(9)

and

$\Delta_{\phi}\equiv
M_{\phi_{2}}-M_{\phi_{1}}=\frac{2\mu_{\Phi\phi}v_{D}}{M_{\phi_{1}}+M_{\phi_{2}}}.$
(10)

Notice the interaction states $(h,h_{D})$ would be rotated to the mass states
$(h_{1},h_{2})$ via the mixing angle $\theta$. The $\sin\theta$ is chosen to
be small enough in our analysis such that it is consistent with the SM-like
Higgs boson data at the LHC333The long-lived dark Higgs phenomenology study in
the inelastic DM model at Belle II can be found in Ref. Duerr:2020muu ..

Finally, since all SM fermions don’t carry $U(1)_{D}$ charges, the only way
for the new $X_{\mu}$ boson and SM fermions to interact is via the kinetic
mixing between $B_{\mu\nu}$ and $X_{\mu\nu}$. The Lagrangian density of this
part can be represented as

${\cal
L}_{X,gauge}=-\frac{1}{4}X_{\mu\nu}X^{\mu\nu}-\frac{\sin\epsilon}{2}B_{\mu\nu}X^{\mu\nu},$
(11)

where $\epsilon$ is the kinetic mixing parameter between these two $U(1)$s. If
we apply the linear order approximation in $\epsilon$, the extra interaction
terms for SM fermions and $Z^{\prime}$ boson can be written as

${\cal L}_{Z^{\prime}f\overline{f}}=-\epsilon
ec_{W}\sum_{f}x_{f}\overline{f}\not{Z}^{\prime}f,$ (12)

where $c_{W}$ is the weak mixing angle and $x_{l}=-1$, $x_{\nu}=0$,
$x_{q}=\frac{2}{3}$ or $-\frac{1}{3}$ depending on the electrical charge of
quark. The dark photon mass can be approximated as

$m_{Z^{\prime}}\simeq g_{D}Q_{D}(\Phi)v_{D}.$ (13)

Notice the correction from the kinetic mixing term is second order in
$\epsilon$ which can be safely neglected here.

###  II.2 The fermion model Ko:2019wxq ; Baek:2020owl

Here we consider a dark sector with the singlet complex scalar $\Phi$ and
Dirac fermion $\chi$ as the dark Higgs and dark matter sectors, respectively.
Both $\Phi$ and $\chi$ are charged under $U(1)_{D}$, but neutral of the SM
gauge symmetry. The $U(1)_{D}$ charges for them are assigned as
$Q_{D}(\Phi)=+2$ and $Q_{D}(\chi)=+1$. Again, the SM-like Higgs doublet and
other SM particles do not carry $U(1)_{D}$ charges.

The scalar part of the renormalizable and gauge invariant Lagrangian density
is

${\cal L}_{scalar}=|D_{\mu}H|^{2}+|D_{\mu}\Phi|^{2}-V(H,\Phi),$ (14)

where $D_{\mu}H$ and $D_{\mu}\Phi$ are the same in Eq.(2) and (3). The scalar
potential in Eq.(14) is given by

$\displaystyle V(H,\Phi)=$
$\displaystyle-\mu^{2}_{H}H^{\dagger}H+\lambda_{H}(H^{\dagger}H)^{2}-\mu^{2}_{\Phi}\Phi^{\ast}\Phi+\lambda_{\Phi}(\Phi^{\ast}\Phi)^{2}$
$\displaystyle+\lambda_{H\Phi}(H^{\dagger}H)(\Phi^{\ast}\Phi),$ (15)

where all parameters are assumed to be real for simplicity.

The Lagrangian density of dark matter sector part is

${\cal L}_{\chi}=\overline{\chi}(i{\hbox to0.0pt{\,/\hss}\partial}+g_{D}{\hbox
to0.0pt{\,/\hss}X}-M_{\chi})\chi-(\frac{f}{2}\overline{\chi^{c}}\chi\Phi^{\ast}+H.c.),$
(16)

where $f$ is assumed to be a real parameter. In order to decompose the Dirac
fermion $\chi$ into a pair of two independent Majorana fermions, $\chi_{1}$
and $\chi_{2}$, we set

$\chi_{1,2}(x)=\frac{1}{\sqrt{2}}(\chi(x)\mp\chi^{c}(x)).$ (17)

After expanding $H,\Phi$ fields around the vacuum with the unitary gauge as
shown in Eq.(6), the Eq.(16) can be written as

$\displaystyle{\cal L}_{\chi}=$
$\displaystyle\frac{1}{2}\overline{\chi_{2}}(i{\hbox
to0.0pt{\,/\hss}\partial}-M_{\chi_{2}})\chi_{2}+\frac{1}{2}\overline{\chi_{1}}(i{\hbox
to0.0pt{\,/\hss}\partial}-M_{\chi_{1}})\chi_{1}$
$\displaystyle-i\frac{g_{D}}{2}(\overline{\chi_{2}}{\hbox
to0.0pt{\,/\hss}X}\chi_{1}-\overline{\chi_{1}}{\hbox
to0.0pt{\,/\hss}X}\chi_{2})-\frac{f}{2}h_{D}(\overline{\chi_{2}}\chi_{2}-\overline{\chi_{1}}\chi_{1}),$
(18)

where $\chi_{1}$, $\chi_{2}$ masses and their mass splitting can be
represented as

$M_{\chi_{1,2}}=M_{\chi}\mp fv_{D},$ (19)

and

$\Delta_{\chi}\equiv(M_{\chi_{2}}-M_{\chi_{1}})=2fv_{D}.$ (20)

Finally, the $Z^{\prime}$ boson mass and its interactions with SM fermions are
the same with Eq.(13) and (12).

###  II.3 The target parameter space and decay width of $\phi_{2}(\chi_{2})$

In this work, we focus on the scenario
$m_{Z^{\prime}}>M_{\phi_{1}}+M_{\phi_{2}}(M_{\chi_{1}}+M_{\chi_{2}})$ such
that the $Z^{\prime}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ decay mode
is kinematically allowed and becomes the dominant one. On the other hand, we
are interested in the co-annihilation dominant channel for DM relic density in
early Universe, so we restrict ourselves to the compressed mass spectrum with
$\Delta_{\phi,\chi}<0.5M_{\phi_{1},\chi_{1}}$ in inelastic DM models
Duerr:2019dmv . Finally, we concentrate on $M_{\phi_{1},\chi_{1}}>100$ MeV
such that the parameter space is free from the Big Bang Nucleosynthesis (BBN)
constraint Boehm:2012gr ; Depta:2019lbe  444Notice that the DM mass can be
extended to lower regions if the light dark Higgs is included as shown in Ref.
Baek:2020owl ..

The $\phi_{1}(\chi_{1})$ is the only DM candidate in our analysis. If the mass
splitting $\Delta_{\phi,\chi}$ is not ignorable compared with the
$M_{\phi_{1}(\chi_{1})}$, $\phi_{2}(\chi_{2})$ is not stable and will decay to
$\phi_{1}(\chi_{1})$ and a SM fermion pair via the off-shell $Z^{\prime}$. The
full analytical formulas for the total width of $\phi_{2}(\chi_{2})$ three-
body decay can be found in Appendix B of Ref. Giudice:2017zke . In the mass
range of our interest, there are three kind of $\phi_{2}(\chi_{2})$ decay
modes 555Since the light $Z^{\prime}$ is isosinglet, it can mix with $\omega$
meson and decay to three pions, $Z^{\prime}\rightarrow\pi^{+}\pi^{-}\pi^{0}$,
if kinematically allowed. However, according to Fig.1 in Ref. Ilten:2018crw ,
the contributions from $Z^{\prime}\rightarrow\pi^{+}\pi^{-}\pi^{0}$ mode are
important only in the adjacent regions of $m_{Z^{\prime}}\approx
m_{\omega}=0.782$ GeV. Hence, we don’t specifically study this regions in our
analysis., $\phi_{2}(\chi_{2})\rightarrow\phi_{1}(\chi_{1})e^{+}e^{-}$,
$\phi_{1}(\chi_{1})\mu^{+}\mu^{-}$ and $\phi_{1}(\chi_{1})\pi^{+}\pi^{-}$. In
our numerical calculations, the total decay width of $\phi_{2}(\chi_{2})$ is
automatically calculated in Madgraph5 aMC@NLO Alwall:2014hca . For the partial
decay width of
$\phi_{2}(\chi_{2})\rightarrow\phi_{1}(\chi_{1})\pi^{+}\pi^{-}$, we rescale
the partial decay width of
$\phi_{2}(\chi_{2})\rightarrow\phi_{1}(\chi_{1})\mu^{+}\mu^{-}$ with the
measured $R(s)$ values in Ref. Zyla:2020zbs .

##  III The relevant cross sections

In this section, we show both analytical representations and numerical results
for cross sections of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ and
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes.
Especially, the method to distinguish spin-0 and spin-1/2 DM candidates in the
inelastic DM models are also discussed here.

###  III.1 The analytical representations

The differential cross section for $e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}$ via
the s-channel $Z^{\prime}$ can be represented as

$\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})}{dt}=\frac{\epsilon^{2}e^{2}g^{2}_{D}}{8\pi
s^{2}}\frac{\left[(M^{2}_{\phi_{1}}-t)(M^{2}_{\phi_{2}}-t)-st\right]}{\left[(s-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}\right]},$
(21)

where $s$, $t$ are Mandelstam variables and $\Gamma_{Z^{\prime}}$ is the total
width of $Z^{\prime}$ boson. In order to study the differential angular
distributions of cross sections, we transfer from $d\sigma/dt$ to
$d\sigma/d\cos\theta$, where $\theta$ is the polar angle of $\phi_{2}$ and it
is defined as the direction of $\phi_{2}$ relative to the positron beam
direction.

The $d\sigma/d\cos\theta$ in the centre-of-mass (CM) frame can be simply
written as

$\left.\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})}{d\cos\theta}\right|_{CM}=\frac{1}{64\pi}\frac{\epsilon^{2}e^{2}g^{2}_{D}E^{2}_{CM}}{[(E^{2}_{CM}-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}]}\xi^{3/2}(1-\cos^{2}\theta),$
(22)

and

$\xi=\sqrt{1-\frac{2(M^{2}_{\phi_{2}}+M^{2}_{\phi_{1}})}{E^{2}_{CM}}+\frac{(M^{2}_{\phi_{2}}-M^{2}_{\phi_{1}})^{2}}{E^{4}_{CM}}},$
(23)

where $E_{CM}$ is the centre-of-mass energy. We can find Eq.(22) is always
proportional to $1-\cos^{2}\theta$. Consequently, it indicates if the
s-channel $Z^{\prime}$ is on-shell produced, it is longitudinally polarized
(helicity $=0$).

However, the formula for $d\sigma/d\cos\theta$ in the Belle II laboratory
(LAB) frame is more tedious, so we show it in Eq.(44) of the Appendix A. As we
will see the numerical results, the differential angular distributions of
cross sections in the LAB frame are just shifted to the initial electron beam
direction for the Belle II machine compared with the ones in the CM frame.

The differential cross section for $e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}$ via
s-channel $Z^{\prime}$ can be represented as

$\frac{d\sigma(e^{+}e^{-}\rightarrow\chi_{1}\chi_{2})}{dt}=\frac{\epsilon^{2}e^{2}g^{2}_{D}}{8\pi
s^{2}}\frac{\left[s^{2}+2(t-M^{2}_{\chi_{1}})(t-M^{2}_{\chi_{2}})+s(2t-(M_{\chi_{2}}-M_{\chi_{1}})^{2})\right]}{\left[(s-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}\right]},$
(24)

and $d\sigma/d\cos\theta$ in the CM frame can be written as

$\displaystyle\left.\frac{d\sigma(e^{+}e^{-}\rightarrow\chi_{1}\chi_{2})}{d\cos\theta}\right|_{CM}=$
$\displaystyle\frac{1}{32\pi}\frac{\epsilon^{2}e^{2}g^{2}_{D}E^{2}_{CM}}{[(E^{2}_{CM}-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}]}\times$
$\displaystyle\left[(1-\frac{(M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}})^{2}}{E^{4}_{CM}}+\frac{4M_{\chi_{1}}M_{\chi_{2}}}{E^{2}_{CM}})\xi+\xi^{3/2}\cos^{2}\theta\right],$
(25)

where $\xi$ can be found in Eq.(23) with $M_{\phi_{1,2}}\rightarrow
M_{\chi_{1,2}}$. For the $M_{\chi_{1,2}}\rightarrow 0$ limit, Eq.(25) is
proportional to $1+\cos^{2}\theta$, which is well known. It shows if the
s-channel $Z^{\prime}$ is on-shell produced, it is transversely polarized
(helicity $=\pm$). However, the mass effects of $M_{\chi_{1,2}}$ spoil parts
of this property. Similarly, we show the formula for $d\sigma/d\cos\theta$ in
the LAB frame in Eq.(45) of the Appendix A. Again, there is a skewing behavior
for the differential angular distributions of cross sections in the LAB frame
compared with the ones in the CM frame.

The initial state radiation (ISR) photon is an useful trigger for signals with
missing energy or soft objects at lepton colliders. Especially, for the
process with on-shell $Z^{\prime}$ production, the ISR photon is used not only
for background rejection, but also for $Z^{\prime}$ invariant mass
reconstruction. To include the ISR photon in
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ process, the
differential cross section can be written as Birkedal:2004xn

$\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma)}{dzd\cos\theta^{\prime}}\simeq{\cal
P}(z,\cos\theta^{\prime})\widehat{\sigma}(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})),$
(26)

and the splitting kernal ${\cal P}(z,\cos\theta^{\prime})$ is

${\cal
P}(z,\cos\theta^{\prime})=\frac{\alpha}{\pi}\frac{1+(1-z)^{2}}{z}\frac{1}{\sin^{2}\theta^{\prime}},$
(27)

where $\alpha$ is the fine structure constant, $z=\frac{E_{\gamma}}{E_{CM}/2}$
is the energy fraction of ISR photon from the initial electron/positron and
$\theta^{\prime}$ is the polar angle of ISR photon and it is defined as the
direction of $\gamma$ to the positron beam direction. The differential form of
$\widehat{\sigma}(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2}))$
can be found in Eq.(21) and (24). In the CM frame, if we assign the four
momentum of the initial electron/positron and ISR photon as
$p_{e^{\pm}}=\left(E_{CM}/2,0,0,\pm E_{CM}/2\right)$ and
$p_{ISR}=\left(E_{\gamma},p_{x,\gamma},p_{y,\gamma},p_{z,\gamma}\right)$, the
four momentum of electron/positron after the radiation is
$p_{e^{\prime\pm}}=\left(E_{CM}-E_{\gamma},-p_{x,\gamma},-p_{y,\gamma},\pm
E_{CM}/2-p_{z,\gamma}\right)$. Let’s consider two limit scenarios to
intuitively catch up the behavior of
$d\widehat{\sigma}(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2}))/d\cos\theta$
distributons. First, if the ISR photon is soft ($z\approx 0$),
$p_{e^{\prime\pm}}\approx p_{e^{\pm}}$ and
$d\widehat{\sigma}(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2}))/d\cos\theta$
distributions are close to Eq.(22) and (25). On the contrary, if the ISR
photon takes almost all energy from the initial electron/positron
($z\rightarrow 1$), the $\phi_{2}(\chi_{2})$ is obvious in the forward or
backward direction. Since
$d\widehat{\sigma}(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2}))/d\cos\theta$
distributions are highly dependent on the $z$ parameter of the ISR photon, the
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes are
not ideal to distinguish scalar and fermion inelastic DM models.

Finally, we focus on the parameter space
$m_{Z^{\prime}}>M_{\phi_{2}(\chi_{2})}+M_{\phi_{1}(\chi_{1})}$ in this work,
so only the three-body decay of $\phi_{2}(\chi_{2})$ via the off-shell
$Z^{\prime}$ is possible. The full analytical representations of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\rightarrow\phi_{1}\phi_{1}(\chi_{1}\chi_{1})f\overline{f}$
are shown in the Appendix B.

###  III.2 The numerical results

Figure 1:  The relations between $M_{\phi_{1},\chi_{1}}$ (GeV) and $\sigma$
(pb) for $e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}$ and
$e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}$ processes (left panel) and
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}\gamma$ and
$e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}\gamma$ processes (right panel).

We first generate both scalar and fermion inelastic DM UFO model files from
FeynRules Alloul:2013bka , and then calculate cross sections of
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})$ and
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})\gamma$ processes via
Madgraph5 aMC@NLO666The $Z^{\prime}$ boson total decay width is automatically
calculated in Madgraph5<EMAIL_ADDRESS>The beam energies are set to be
$E(e^{+})=4.0$ GeV and $E(e^{-})=7.0$ GeV which are consistent with the Belle
II experiment. In order to display the relations of cross sections with
$M_{\phi_{1},\chi_{1}}$ and $m_{Z^{\prime}}$, we study
$m_{Z^{\prime}}=2.5M_{\phi_{1},\chi_{1}}$ and
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ with fixed $\alpha_{D}\equiv
g^{2}_{D}/4\pi=0.1$, $\epsilon=0.01$ and
$\Delta_{\phi,\chi}=0.1M_{\phi_{1},\chi_{1}}$. The relation between
$M_{\phi_{1},\chi_{1}}$ (GeV) and $\sigma$ (pb) for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})$ processes are shown
in the left panel of Fig. 1. Since all of these cross sections are
proportional to $\epsilon^{2}\alpha_{D}$, it is straightforward to rescale
cross sections with different values of $\epsilon$ and $\alpha_{D}$. On the
other hand, the influence from the changes of $\Delta_{\phi,\chi}$ to cross
sections are mild for $\Delta_{\phi,\chi}<0.5M_{\phi_{1},\chi_{1}}$. The
scalar and fermion pair production cross sections can be scaled by
$\beta^{3/2}$ and $\beta^{1/2}$ respectively, where $\beta$ is the velocity of
the final state particle in the CM frame. Because of this extra $\beta$ factor
for the scalar case, cross sections for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}$ are suppressed compared with
$e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}$.

On the other hand, the relation between $M_{\phi_{1},\chi_{1}}$ (GeV) and
$\sigma$ (pb) for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})\gamma$ processes are
shown in the right panel of Fig. 1 with the same parameter settings. Here, the
basic cuts $E(\gamma)>0.1$ GeV and $\eta<2.203$ are applied for the ISR
photon. We find the main contributions of
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})\gamma$ come from
$e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})\gamma$
processes. It explains why scalar and fermion pair cross sections in this
process are very close to each other in the right panel of Fig. 1.

Figure 2:  The $(1/\sigma)(d\sigma/d\cos\theta)$ distributions for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}$ in CM frame (top left panel) and Belle
II LAB frame (bottom left panel) and $e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}$
in CM frame (top right panel) and Belle II LAB frame (bottom right panel).

We then turn to the study of $(1/\sigma)(d\sigma/d\cos\theta)$ distributions
for $e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}$ and
$e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}$ processes with fixed $\alpha_{D}=0.1$,
$\epsilon=0.01$, $\Delta_{\phi,\chi}=0.1M_{\phi_{1},\chi_{1}}$, and
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$. The results in the CM and Belle II
LAB frame are shown in Fig. 2. It is clear that the
$(1/\sigma)(d\sigma/d\cos\theta)$ distributions for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}$ in the CM frame is proportional to
$1-\cos^{2}\theta$ and $e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}$ behaves close
to $1+\cos^{2}\theta$ with some distortions from $M_{\chi_{1,2}}$ effects.
However, since the measurement of time-of-flight for the long-lived particle
is poor at the Belle II machine, we will lose this information and cannot make
the Lorentz transformation for the four-vector of displaced vertex from the
Belle II LAB frame to the CM frame. Hence, the distributions in the CM frame
is only a reference for the comparison with the ones in the Belle II LAB
frame. Hopefully, the real situation for the distributions in the Belle II LAB
frame is just shifted to the electron beam direction for the Belle II machine
compared with the ones in the CM frame. It is still clear to see the different
distributions between scalar and fermion inelastic DM models.

Figure 3:  The $(1/\sigma)(d\sigma/d\cos\theta)$ distributions for
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}\gamma$ in CM frame (top left panel) and
Belle II LAB frame (bottom left panel) and
$e^{+}e^{-}\rightarrow\chi_{2}\chi_{1}\gamma$ in CM frame (top right panel)
and Belle II LAB frame (bottom right panel).

Finally, we fix the same parameter settings as Fig. 2 for the
$(1/\sigma)(d\sigma/d\cos\theta)$ distributions of
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})\gamma$ processes in
the CM and LAB frame in Fig. 3. Again, the distributions in the CM frame is
only for the comparison and it is clear to see the skewing behavior for the
differential angular distributions of cross sections in the LAB frame compared
with the ones in the CM frame. As discussed in the previous discussion,
because these distributions are highly dependent on the $z$ parameter of the
ISR photon, the
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes are
not ideal for the DM spin discrimination. Only when the $Z^{\prime}$ is on-
shell produced and heavy, the ISR photon becomes soft ($z$ is small) such that
differences of $(1/\sigma)(d\sigma/d\cos\theta)$ distributions from scalar and
fermion inelastic DM models show up.

In summary, the size for cross sections of
$e^{+}e^{-}\rightarrow\phi_{2}\phi_{1}(\chi_{2}\chi_{1})$ and their
$(1/\sigma)(d\sigma/d\cos\theta)$ distributions in the Belle II LAB frame can
help us to statistically distinguish fermion inelastic DM model from the
scalar one even for the same parameter settings.

##  IV Search for long-lived particles in inelastic dark matter models at
Belle II

In this section, we study how to search for long-lived particles in inelastic
DM models at Belle II experiment. We first briefly overview the Belle II
experiment and signal signatures, and then show some interesting kinematic
distributions based on four signal benchmark points. We further set up event
selections for dilepton displaced signature and relevant results are
discussed, especially for the discrimination of scalar and fermion inelastic
DM models. Finally, we solve kinematic equations of
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
and $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
processes to determine both $M_{\chi_{1}}$ and $M_{\chi_{2}}$.

###  IV.1 The Belle II experiment and signal signatures

The SuperKEKB accelerator of Belle II experiment is a circular asymmetric
$e^{+}e^{-}$ collider with the nominal collision energy of $\sqrt{s}=10.58$
GeV. The beam parameters are $E(e^{+})=4$ GeV and $E(e^{-})=7$ GeV. The
planned full integrated luminosity for the final dataset is $50~{}{\rm
ab}^{-1}$. In this work, the following Belle II sub-detectors are relevant :
the tracking system including vertex detectors (VXD) and central drift chamber
(CDC), the electromagnetic calorimeter (ECAL), and muon system.

We consider the following two kinds of processes in inelastic DM
models777Because the major contribuitons of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ come from the
on-shell $Z^{\prime}$ production, hence, we only focus on the process
$e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ in the
following analysis unless noted otherwise.:

$\displaystyle
e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\rightarrow\phi_{1}\phi_{1}(\chi_{1}\chi_{1})l^{+}l^{-}\
,$ (28) $\displaystyle e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma\rightarrow\phi_{1}\phi_{1}(\chi_{1}\chi_{1})l^{+}l^{-}\gamma\
.$ (29)

And there are five possible signatures can be searched for at the Belle II
experiment :

  * •

Mono-$\gamma$ : If the excited DM state decays outside the detector or the
decay products are too soft to be detected in Eq.(29).

  * •

Mono-$\gamma$ with prompt lepton pair: The excited DM state promptly decays
and the visible products can be successfully detected and defined in Eq.(29).

  * •

Mono-$\gamma$ with displaced lepton pair : The excited DM state is long-lived
and decays inside the detector leaving the displaced vertex in Eq.(29).

  * •

Only prompt lepton pair : The same as the second one but without ISR photon in
Eq.(28).

  * •

Only displaced lepton pair : The same as the third one but without ISR photon
in Eq.(28).

The Mono-$\gamma$ signature is well-studied as shown in Ref. Zhang:2019wnz ;
Lees:2017lec ; Duerr:2019dmv . Since we focus on searching for long-lived
particles in this paper, only signatures with displaced vertex are studied. On
the other hand, the signature with displaced charged pion pair is also
possible from the long-lived excited DM state decay. However, for simplicity,
here we only study signatures with displaced electron and muon pair.

In order to make our results more realistic, we follow Ref. Adachi:2018qme  to
involve the detector resolution effects with Gaussian smearing at Belle II
experiment. The tracking resolution of electron/muon momentum in the CDC is
given by

$\sigma_{p^{l}_{t}}/p^{l}_{t}=0.0011p^{l}_{t}~{}[\rm{GeV}]\oplus
0.0025/\beta,$ (30)

where $p^{l}_{t}$ is the transverse momentum of electron/muon track and
$\beta$ is its velocity in the natural unit. We conservatively apply
$\sigma_{p_{l}}/p_{l}=0.005$ in our event analysis. On the other hand, the
muon efficiency in the muon system is approximated to $0.98$. The photon
momentum resolution in the ECAL is approximated to
$\sigma_{E_{\gamma}}/E_{\gamma}=2\%$ where $E_{\gamma}$ is the energy of
photon. Finally, we use the resolution of $\sigma_{r_{DV}}=26\mu$m for the
displaced vertex vector of lepton pair in our analysis.

###  IV.2 Benchmark points and kinematic distributions

In this study, we use the following four benchmark points (BPs) to display our
analysis :

  * •

(I) $M_{\phi_{1},\chi_{1}}=0.3$ GeV,
$\Delta_{\phi_{1},\chi_{1}}=0.4M_{\phi_{1},\chi_{1}}$,
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ and $\epsilon=2\times 10^{-2}$

  * •

(II) $M_{\phi_{1},\chi_{1}}=3.0$ GeV,
$\Delta_{\phi_{1},\chi_{1}}=0.1M_{\phi_{1},\chi_{1}}$,
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ and $\epsilon=2\times 10^{-3}$

  * •

(III) $M_{\phi_{1},\chi_{1}}=1.0$ GeV,
$\Delta_{\phi_{1},\chi_{1}}=0.4M_{\phi_{1},\chi_{1}}$,
$m_{Z^{\prime}}=2.5M_{\phi_{1},\chi_{1}}$ and $\epsilon=10^{-3}$

  * •

(IV) $M_{\phi_{1},\chi_{1}}=2.0$ GeV,
$\Delta_{\phi_{1},\chi_{1}}=0.2M_{\phi_{1},\chi_{1}}$,
$m_{Z^{\prime}}=2.5M_{\phi_{1},\chi_{1}}$ and $\epsilon=10^{-3}$

with fixed $\alpha_{D}=0.1$. The first BP is inspired from the allowed
parameter space to explain the muon anomalous magnetic moment excess
Mohlabeng:2019vrz . Other three BPs are the ones which cannot be covered from
Mono-$\gamma$ searches from BaBar and Belle II Duerr:2019dmv . In addition,
they can also be searched for in the future long-lived particle experiments,
like CODEX-b Gligorov:2017nwh , SeaQuest Berlin:2018pwi , FASER Feng:2017uoz ,
and MATHUSLA Chou:2016lxi .

Figure 4:  Various kinematic distributions for
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
based on the four BPs in the main text. Top left panel: lepton pair invariant
mass, $M_{ll}$ (GeV), top right panel: projection of the decay vertex distance
on $z$ axis, $r_{z}$ (mm), bottom left panel: lepton pair energy, $E_{ll}$
(GeV), bottom right panel: angular distance between lepton pair, $\Delta
R_{ll}$.

Figure 5:  Various kinematic distributions for $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
based on the four BPs in the main text. Top left panel: lepton pair energy,
$E_{ll}$ (GeV), top right panel: angular distance between lepton pair, $\Delta
R_{ll}$, middle left panel: projection of the decay vertex distance on $z$
axis, $r_{z}$ (mm), middle right panel: the maximum angular distance between
ISR photon and one of the lepton in lepton pair, $\Delta R^{max}_{\gamma l}$,
bottom left panel: photon energy in the CM frame, $E(\gamma)$ (GeV), bottom
right panel : reconstructed $Z^{\prime}$ invariant mass, $m_{Z^{\prime}}$
(GeV).

We show some interesting kinematic distributions based on the above four BPs
for
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
and $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
in Fig. 4 and Fig. 5, respectively. First, the lepton pair invariant mass is
smaller or equal to the mass splitting $\Delta_{\chi}$ (dashed lines),
$2m_{l}\leq M_{ll}\leq\Delta_{\chi}$. Once we have enough signal events, the
threshold value of $M_{ll}$ can help us to roughly determine the mass
splitting $\Delta_{\chi}$ in inelastic DM models. Second, the distance of
displaced vertex of $\chi_{2}$ is not only dependent on $M_{\chi_{2}}$,
$\Delta_{\chi}$, $\epsilon$ and $\alpha_{D}$, but also the boost of $\chi_{2}$
in the LAB frame. Here we show projection of the decay vertex distance on $z$
axis, $r_{z}$. It is clear to see that the $\chi_{2}$ decay length in BP1
(BP2) is the shortest (longest) one on average. The distributions for
projection of the decay vertex distance on the transverse plane, $R_{xy}$, are
similar. On the other hand, $E_{ll}$ and $\Delta R_{ll}$ are energy and
angular distance for the lepton pair. We can find the boost of $\chi_{2}$ is
smaller for the process involving ISR photon by comparing $r_{z}$, $E_{ll}$
and $\Delta R_{ll}$ distributions in Fig. 4 and Fig. 5. The larger boost on
$\chi_{2}$ causes the longer $r_{z}$, larger $E_{ll}$, and smaller $\Delta
R_{ll}$ distributions. Besides, we will see larger $E_{ll}$, and smaller
$r_{z}$ distributions of BP1 suffer from more severe smearing effects from
detector resolutions such that it is more challenge to reconstruct DM mass and
mass splitting between DM excited and ground states in the real experiment.
Third, if the $Z^{\prime}$ is light enough, the ISR photon and $Z^{\prime}$
are almost back-to-back produced such that the lepton from $Z^{\prime}$ decay
can have large angular distance with the ISR photon. Here we show the maximum
angular distance between ISR photon and one of the lepton in lepton pair,
$\Delta R^{max}_{\gamma l}$. Fourth, we boost the $E(\gamma)$ distribution of
ISR photon from the Belle II LAB frame to the CM frame and it is clear to see
its mono-energetic behavior. This feature can largely reduce SM backgrounds
which are continuous in the photon energy spectrum. Finally, the $Z^{\prime}$
invariant mass can be reconstructed from Eq.(39) with measurable observables.
Because of the detector resolution, it is clear to find smearing effects from
the true $Z^{\prime}$ mass (dashed lines), especially for the case of light
$Z^{\prime}$. We will see this phenomenon critically affects abilities to
determine DM mass and the mass splitting between DM excited and ground states
for the process in Eq.(29). Notice those kinematic distributions for the
scalar inelastic DM model are similar to Fig. 4 and Fig. 5, so we do not show
them again.

###  IV.3 Event selections and results

We closely follow Ref. Duerr:2019dmv  to set up event selections for the
displaced signature at Belle II. We only conservatively consider the following
two background-free regions after event selections in our analysis :

  * •

low $R_{xy}$ region ($100\%$ detection efficiency) : $0.2<R_{xy}\leq 0.9$ cm
(electron), $0.2<R_{xy}\leq 17.0$ cm (muon).

  * •

high $R_{xy}$ region ($30\%$ detection efficiency) : $17.0<R_{xy}\leq 60.0$ cm
(electron, muon).

According to the arguments in Ref. Duerr:2019dmv , backgrounds from photon
conversion can be further reduced to a negligible level by requiring the
lepton pair invariant mass $M_{ll}\geq 0.03$ GeV, and their opening angle
larger than $0.1$ rad. Notice backgrounds from the photon conversion to muon
pair in the region $0.9<R_{xy}\leq 17.0$ cm are negligible such that the
analysis for muon in low $R_{xy}$ region can be extended.

Table 1: Event selections for the displaced vertex analysis at Belle II. Objects | Selections
---|---
displaced vertex | (i) $-55\,\mathrm{cm}\leq z\leq 140\,\mathrm{cm}$
(ii) $17^{\circ}\leq\theta^{\text{DV}}_{\text{LAB}}\leq 150^{\circ}$
electrons | (i) both $E(e^{+})$ and $E(e^{-})>0.1\,\mathrm{GeV}$
(ii) opening angle of pair $\theta_{ee}>0.1$ rad
(iii) invariant mass of pair $m_{ee}>0.03\,\mathrm{GeV}$
muons | (i) both $p_{\text{T}}(\mu^{+})$ and $p_{\text{T}}(\mu^{-})>0.05\,\mathrm{GeV}$
(ii) opening angle of pair $\theta_{\mu\mu}>0.1$ rad
(iii) invariant mass of pair $m_{\mu\mu}>0.03\,\mathrm{GeV}$
(iv) veto $0.48\,\mathrm{GeV}\leq m_{\mu\mu}\leq 0.52\,\mathrm{GeV}$
photons | (i) $E^{\gamma}_{\text{LAB}}>0.5\,\mathrm{GeV}$
(ii) $17^{\circ}\leq\theta^{\gamma}_{\text{LAB}}\leq 150^{\circ}$

We summarize all event selections in our analysis in Table 1. For the muon
pair in the final state, we veto the invariant mass region, $0.48\leq
m_{\mu^{+}\mu^{-}}\leq 0.52$ GeV, to reject backgrounds from $K^{0}_{S}$
decay. More discussions on the trigger issue can be found in Ref.
Duerr:2019dmv . Here, event selections in Table 1 are done offline and we
simply assume these events are already triggered and stored.

Table 2: The analysis results of four BPs for process in Eq.(28). The first column is the type of inelastic DM models. The second column is labeled by BP which $e$ and $\mu$ denote the type of lepton pair in the final state. The third column is production cross sections for process in Eq.(28). The fourth and fifth columns are efficiencies of low $R_{xy}$ and high $R_{xy}$ regions defined in the main text, respectively. Final column is the number of signal events with an integrated luminosity of $1ab^{-1}$. Type | BP |  $\sigma$ (fb) | Eff.(low $R_{xy}$) | Eff.(high $R_{xy}$) | $N_{event}$
---|---|---|---|---|---
scalar | BP1$e$ | $948.14$ | $16.98$ | $0\%$ | $1.61\times 10^{5}$
BP2$e$ | $58.39$ | $0.15\%$ | $2.48\%$ | $1.54\times 10^{3}$
BP2$\mu$ | $6.15$ | $0.21\%$ | $3.33\%$ | $217.71$
BP3$e$ | $1.86$ | $10.06\%$ | $0.70\%$ | $200.09$
BP3$\mu$ | $0.61$ | $11.25\%$ | $0.74\%$ | $73.14$
BP4$e$ | $2.23$ | $1.56\%$ | $9.34\%$ | $243.26$
BP4$\mu$ | $0.74$ | $1.72\%$ | $10.78\%$ | $92.50$
fermion | BP1$e$ | $3856.00$ | $14.26\%$ | $0\%$ | $5.50\times 10^{5}$
BP2$e$ | $422.80$ | $0.17\%$ | $2.35\%$ | $1.07\times 10^{4}$
BP2$\mu$ | $44.63$ | $0.22\%$ | $2.97\%$ | $1.42\times 10^{3}$
BP3$e$ | $7.99$ | $10.20\%$ | $0.42\%$ | $848.54$
BP3$\mu$ | $2.69$ | $11.20\%$ | $0.46\%$ | $313.65$
BP4$e$ | $11.71$ | $1.57\%$ | $7.82\%$ | $1.10\times 10^{3}$
BP4$\mu$ | $3.88$ | $1.69\%$ | $8.75\%$ | $405.07$
Table 3: The same as Table 2, but for process in Eq.(29). Type | BP |  $\sigma$ (fb) | Eff.(low $R_{xy}$) | Eff.(high $R_{xy}$) | $N_{event}$
---|---|---|---|---|---
scalar | BP1$e$ | $2472.70$ | $6.70\%$ | $0\%$ | $1.66\times 10^{5}$
BP2$e$ | $159.85$ | $0.16\%$ | $2.27\%$ | $3.88\times 10^{3}$
BP2$\mu$ | $16.85$ | $0.20\%$ | $2.87\%$ | $517.30$
BP3$e$ | $5.13$ | $7.64\%$ | $0.02\%$ | $392.96$
BP3$\mu$ | $1.69$ | $8.83\%$ | $0.03\%$ | $149.73$
BP4$e$ | $7.14$ | $1.86\%$ | $3.29\%$ | $367.71$
BP4$\mu$ | $2.35$ | $2.02\%$ | $2.87\%$ | $114.92$
fermion | BP1$e$ | $2503.60$ | $6.14\%$ | $0\%$ | $1.54\times 10^{5}$
BP2$e$ | $167.10$ | $0.16\%$ | $2.16\%$ | $3.87\times 10^{3}$
BP2$\mu$ | $17.66$ | $0.18\%$ | $2.67\%$ | $503.31$
BP3$e$ | $5.05$ | $7.77\%$ | $0.02\%$ | $393.40$
BP3$\mu$ | $1.70$ | $8.89\%$ | $0.02\%$ | $151.47$
BP4$e$ | $7.14$ | $1.95\%$ | $3.14\%$ | $363.43$
BP4$\mu$ | $2.37$ | $2.05\%$ | $3.44\%$ | $130.11$

The results for four BPs are shown in Table 2 and 3 for processes in Eq.(28)
and (29), respectively. Eff.(low $R_{xy}$) and Eff.(high $R_{xy}$) are
efficiencies of low $R_{xy}$ and high $R_{xy}$ regions after involving the
event selections. Here we use an integrated luminosity of $1~{}{\rm ab}^{-1}$
to calculate number of signal events ($N_{event}$). We can find most events in
BP1 and BP3 (BP2 and BP4) are located in low (high) $R_{xy}$ regions.

Figure 6:  The polar angular distributions in the Belle II LAB frame after
event selections in Table 1 for process in Eq.(28) with electron pair in the
final state and $N_{event}$ in Table 2.

Figure 7:  The same as Fig. 6, but for process in Eq.(29) with electron pair
in the final state and $N_{event}$ in Table 3.

We further compare polar angular distributions in the Belle II LAB frame after
event selections in Table 1 for processes in Eq.(28) and (29) with electron
pair in the final state in Fig. 6 and 7 with $N_{event}$ in Table 2 and 3,
separately. In Fig. 6, since the selected regions can be approximately
background free, we can directly count the events bin by bin to distinguish
the signals are either from scalar or fermion inelastic DM models for the
above four BPs even only ${\cal L}=1~{}{\rm ab}^{-1}$ is used. However, in
Fig. 7, it is hard to distinguish the signals are either from scalar or
fermion inelastic DM models except for the BP2. We have discussed this
behavior in Sec. III that we only have the chance to separate these two models
if the ISR photon is soft enough. Therefore, in addition to the usual search
process in Eq.(29) at B-factories, we propose to consider the process in
Eq.(28) with the displaced vertex trigger suggested in Ref. Duerr:2019dmv
which is more sensitive to distinguish the spin of DM in inelastic DM models.

Figure 8:  The future bounds from
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ and
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes for
event selections in Table 1 with the integrated luminosity of $50~{}{\rm
ab}^{-1}$. Here parameters $\alpha_{D}=0.1$,
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ and
$\Delta_{\phi,\chi}=0.1M_{\phi_{1},\chi_{1}}$ are fixed and $90\%$ C.L.
contours which correspond to an upper limit of 2.3 events with the assumption
of background-free are applied. The model-independent LEP bound Hook:2010tw
and BaBar mono-$\gamma$ bound Lees:2017lec  are also shown.

Finally, we use the most aggressive integrated luminosity with $50~{}{\rm
ab}^{-1}$ for event selections in Table 1 at Belle II for scalar and fermion
inelastic DM models and predict the future bounds from
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ and
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes in
Fig. 8. Notice the off-shell $Z^{\prime}$ production in the second process is
also considered. Here we fix the parameters, $\alpha_{D}=0.1$,
$m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ and
$\Delta_{\phi,\chi}=0.1M_{\phi_{1},\chi_{1}}$, and apply $90\%$ C.L. contours
which correspond to an upper limit of 2.3 events with the assumption of
background-free. Notice constraints from model-independent LEP bound
Hook:2010tw  and BaBar mono-$\gamma$ bound Lees:2017lec  are added for the
comparison. We closely follow Ref. Duerr:2019dmv  for the recasting of BaBar
mono-$\gamma$ constraints in inelastic DM models. For
$M_{\phi_{1},\chi_{1}}<0.4$ GeV, even dilepton displaced vertex is located in
our target regions, the lepton pair in the final state is too soft and
collinear to pass the event selections in Table 1. On the contrary, for
$M_{\phi_{1},\chi_{1}}\gtrsim 4.9$ GeV, the lepton pair in the final state
becomes more energetic and well-separated, but dilepton displaced vertex is
too short to reach our target regions and production cross sections are also
highly suppressed. We find our results in the second process for fermion
inelastic DM model are consistent with the ones in Ref. Duerr:2019dmv  and we
further include detector resolution effects in our analysis. Besides, we
believe our results for the first process in both scalar and fermion inelastic
DM models are first shown in the literature. Most importantly, the future
bounds from the first process can be even stronger than the usual ones in the
second process.

Figure 9:  The same as Fig. 8 but for $m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$
and $\Delta_{\phi,\chi}=0.1M_{\phi_{1},\chi_{1}}$. The green shaded region
bounded by the green dashed lines is the $2\sigma$ allowed region for the
$(g-2)_{\mu}$ excess and the lighter gray region excluded by the $(g-2)_{\mu}$
at $5\sigma$ C.L.

Moreover, inspired by BP1 which can explain the muon $(g-2)_{\mu}$
Mohlabeng:2019vrz , we perform the same analysis in Fig. 9 for
$\alpha_{D}=0.1$, $m_{Z^{\prime}}=3M_{\phi_{1},\chi_{1}}$ and
$\Delta_{\phi,\chi}=0.4M_{\phi_{1},\chi_{1}}$ fixed. We show the 2$\sigma$
allowed and the 5$\sigma$ excluded regions for the muon $(g-2)_{\mu}$ as well
as the model-independent LEP bound Hook:2010tw  and the BaBar mono-$\gamma$
bound Lees:2017lec . Once the mass of DM is increasing, the mass splitting
between DM ground and excited states is also enhanced, hence, the final state
lepton pair becomes more energetic. Therefore, one can find the BaBar
mono-$\gamma$ bound is slightly weakened when the displaced channel is open.
Since the recasting in Ref Mohlabeng:2019vrz  does not include the requirement
on the angle $\theta_{LAB}^{DV}$ in Table 1, the results are weaker than the
ones in Ref. Duerr:2019dmv . That is the reason why this parameter space is
still allowed in Ref Mohlabeng:2019vrz . In order to make a more precise
recasting, we follow Ref. Duerr:2019dmv  for the BaBar mono-$\gamma$ bound and
it is shown in dotted gray line in Fig. 9 which can already cover the
$(g-2)_{\mu}$ allowed window in the parameter space. Even the recasting BaBar
mono-$\gamma$ bound can close this region, we find our result can give much
stronger bound, $\epsilon\sim\mathcal{O}(10^{-5}-10^{-4})$, especially in the
parameter space where the mono-$\gamma$ searches are suppressed. Therefore we
can explicitly cover this area using displaced lepton search at the Belle II
in a very first stage of the run. We found the future bounds from the process
without ISR photon can be stronger than the one with ISR photon in Fig. 8. On
the other hand, in the lower $M_{\phi_{1}}$ region, the bound for process with
ISR photon get stronger because its cross sections are larger than the ones
without ISR photon as shown in Fig. 1.

###  IV.4 Determination of the DM mass and mass splitting of DM sector

In general, we cannot uniquely determine the DM mass at colliders because of
lacking enough constraints for invisible particles in the final state. We can
take a mono-$\gamma$ search for DM pair production at B-factories as an
example. If the DM is denoted by $\chi$, the process for mono-$\gamma$ search
at B-factories is $e^{+}e^{-}\rightarrow\chi\chi\gamma$ (or
$e^{+}e^{-}\rightarrow\chi\overline{\chi}\gamma$). There are eight unknown
values from four-momentum of two DMs in the final state. Unfortunately, only
five constraints in this process : four from the four-momentum conservation
and one from the same mass for the DM pair. We still need three extra
conditions to uniquely determine the DM mass for each event.

Now we turn our attention to the case of inelastic DM models. For the process
in Eq.(28) or (29), if the $\phi_{2}(\chi_{2})$ is long-lived and leave the
displaced vertex at the Belle II detectors, we will have two extra
constraints. Again, there are still eight unknown values from four-momentum of
two $\phi_{1}^{\prime}s(\chi_{1}^{\prime}s)$ in the final state. However,
because of the charge neutrality of the $\phi_{2}(\chi_{2})$, a three-momentum
vector of $\phi_{2}(\chi_{2})$ is proportional to the direction of displaced
vertex (DV) Kang:2019ukr ; Bae:2020dwf

$\overrightarrow{p}_{\phi_{2}(\chi_{2})}=|\overrightarrow{p}_{\phi_{2}(\chi_{2})}|\widehat{r}_{DV},$
(31)

where $\overrightarrow{p}_{\phi_{2}(\chi_{2})}$ is the three-momentum vector
of $\phi_{2}(\chi_{2})$ and $\widehat{r}_{DV}$ is the unit vector of displaced
vertex from $\phi_{2}(\chi_{2})$. Therefore, we have two more constraints, and
there are seven constraints for this kind of processes in total. We still need
one more condition to uniquely determine the DM mass. Let’s first study the
kinematic equation for processes in Eq.(28) and (29) and solve them event-by-
event from our Monte Carlo samples and then involve the detector resolution
effects888Here we use the fermion inelastic DM model as an example. The same
method can be applied to the scalar inelastic DM model, so we will not repeat
it thereafter..

For
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$,
we first set the direction of DV and following four momentum of $e^{+}$,
$e^{-}$, $\chi_{2}$ and $\chi_{1}$ as

$\widehat{r}_{DV}=(\sin\theta\cos\phi,\sin\theta\sin\phi,\cos\theta),$ (32)

and

$\displaystyle p_{e^{+}}=(E_{+},0,0,E_{+}),$ $\displaystyle
p_{e^{-}}=(E_{-},0,0,-E_{-}),$ $\displaystyle
p_{\chi_{2}}=(E_{\chi_{2}},|\overrightarrow{p_{\chi_{2}}}|\widehat{r}_{DV}),$
$\displaystyle p_{\chi_{1}}=(E_{\chi_{1}},\overrightarrow{p_{\chi_{1}}}),$
(33)

with $E_{-}\neq E_{+}$ in the Belle II LAB frame, then according to energy and
momentum conservation, we can write down

$\displaystyle E_{\chi_{1}}=E_{+}+E_{-}-E_{\chi_{2}},$
$\displaystyle\overrightarrow{p_{\chi_{1}}}=(-|\overrightarrow{p_{\chi_{2}}}|\sin\theta\cos\phi,-|\overrightarrow{p_{\chi_{2}}}|\sin\theta\sin\phi,E_{+}-E_{-}-|\overrightarrow{p_{\chi_{2}}}|\cos\theta).$
(34)

After applying the above two equations to
$E^{2}_{\chi_{i}}=|\overrightarrow{p_{\chi_{i}}}|^{2}+M^{2}_{\chi_{i}}$ for
$i=1,2$, the $E_{\chi_{2}}$ can be written as

$\displaystyle E_{\chi_{2}}=$
$\displaystyle\frac{1}{2\left[\sin^{2}\theta(E^{2}_{+}+E^{2}_{-})+2(1+\cos^{2}\theta)E_{+}E_{-}\right]}[(E_{+}+E_{-})(4E_{+}E_{-}+M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}})$
$\displaystyle\pm|(E_{-}-E_{+})\cos\theta|((M^{2}_{\chi_{1}}-4E_{+}E_{-})^{2}-2[2\sin^{2}\theta(E^{2}_{+}+E^{2}_{-})+4\cos^{2}\theta
E_{+}E_{-}$
$\displaystyle+M^{2}_{\chi_{1}}]M^{2}_{\chi_{2}}+M^{4}_{\chi_{2}})^{1/2}].$
(35)

Figure 10:  The solutions for Eq.(37) with $E_{\chi_{2}}$ in Eq.(35) of
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
process on the $(M_{\chi_{1}},\Delta_{\chi})$ plane for four BPs. Here we
display three arbitrary Monte Carlo events for each BP with red, green and
blue lines. On the other hand, the black line shows the kinematic endpoint
measurement of $M^{max}_{l^{+}l^{-}}$.

We then consider the energy and momentum conservation for
$\chi_{2}\rightarrow\chi_{1}f\overline{f}$,

$\displaystyle E_{\chi_{2}}=E_{\chi^{\prime}_{1}}+E_{f}+E_{\overline{f}}\equiv
E_{\chi^{\prime}_{1}}+E_{V^{\prime}},$
$\displaystyle|\overrightarrow{p_{\chi_{2}}}|\widehat{r}_{DV}=\overrightarrow{p_{\chi^{\prime}_{1}}}+\overrightarrow{p_{f}}+\overrightarrow{p_{\overline{f}}}\equiv\overrightarrow{p_{\chi^{\prime}_{1}}}+\overrightarrow{p_{V^{\prime}}}.$
(36)

Finally, we can receive the following equation for inputs of $E_{V^{\prime}}$,
$\overrightarrow{p_{V^{\prime}}}$ and $\widehat{r}_{DV}$,

$M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}}-2E_{\chi_{2}}E_{V^{\prime}}+E^{2}_{V^{\prime}}-|\overrightarrow{p_{V^{\prime}}}|^{2}+2\sqrt{E_{\chi_{2}}^{2}-M^{2}_{\chi_{2}}}(\widehat{r}_{DV}\cdot\overrightarrow{p_{V^{\prime}}})=0,$
(37)

where $E_{\chi_{2}}$ is shown in Eq.(35). In Fig. 10, we display three
arbitrary Monte Carlo events on the $(M_{\chi_{1}},\Delta_{\chi})$ plane for
four BPs according to Eq.(37). On the other hand, the black line shows the
kinematic endpoint measurement $M^{max}_{l^{+}l^{-}}$. We can find that the
lines from all events and $M^{max}_{l^{+}l^{-}}$ cross to the same point which
is the true $(M_{\chi_{1}},\Delta_{\chi})$ in our four BPs. This is a simple
application of the Kinematic Focus Point Method proposed in Ref. Kim:2019prx .

Figure 11:  The same as Fig. 10, but the solutions for Eq.(37) with
$E_{\chi_{2}}$ in Eq.(41) of $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
process.

Similarly, for $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$,
we use the same settings in Eq.(32),(33) and set the four momentum of ISR
photon in the Belle II LAB frame as

$p_{\gamma}=(E_{\gamma},\overrightarrow{p_{\gamma}}).$ (38)

According to energy and momentum conservation, the four momentum of on-shell
$Z^{\prime}$ can be written as

$p_{Z^{\prime}}=(E_{Z^{\prime}},\overrightarrow{p_{Z^{\prime}}})=(E_{+}+E_{-}-E_{\gamma},-p_{x,\gamma},-p_{y,\gamma},E_{+}-E_{-}-p_{z,\gamma}),$
(39)

where $\overrightarrow{p_{\gamma}}=(p_{x,\gamma},p_{y,\gamma},p_{z,\gamma})$
and $p^{2}_{Z^{\prime}}=m^{2}_{Z^{\prime}}$. Furthermore, we represent the
four momentum of $\chi_{1}$ as

$p_{\chi_{1}}=(E_{Z^{\prime}}-E_{\chi_{2}},\overrightarrow{p_{Z^{\prime}}}-|\overrightarrow{p_{\chi_{2}}}|\widehat{r}_{DV}),$
(40)

and write down the $E_{\chi_{2}}$ as

$\displaystyle E_{\chi_{2}}=$
$\displaystyle\frac{1}{2\left[E^{2}_{Z^{\prime}}-(\widehat{r}_{DV}\cdot\overrightarrow{p_{Z^{\prime}}})^{2}\right]}[E_{Z^{\prime}}(m^{2}_{Z^{\prime}}+M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}})$
$\displaystyle\pm|\widehat{r}_{DV}\cdot\overrightarrow{p_{Z^{\prime}}}|\sqrt{(m^{2}_{Z^{\prime}}+M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}})^{2}-4M^{2}_{\chi_{2}}\left[E^{2}_{Z^{\prime}}-(\widehat{r}_{DV}\cdot\overrightarrow{p_{Z^{\prime}}})^{2}\right]}].$
(41)

Finally, the kinematic equation for $\chi_{2}\rightarrow\chi_{1}f\overline{f}$
is the same as Eq.(37) with $E_{\chi_{2}}$ in Eq.(41). Again, we show three
arbitrary Monte Carlo events on the $(M_{\chi_{1}},\Delta_{\chi})$ plane for
four BPs according to Eq.(37) in Fig. 11. All events and
$M^{max}_{l^{+}l^{-}}$ cross to the true $(M_{\chi_{1}},\Delta_{\chi})$ in our
four BPs.

Figure 12:  The solutions for Eq.(37) with $E_{\chi_{2}}$ in Eq.(35) of
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
process on the $(M_{\chi_{1}},M_{\chi_{2}})$ plane for four BPs after
involving detector resolution effects and event selections. Here 100 signal
events are taken into account for each BP and the bin size for both $x$ and
$y$ axes are set to be $10$ MeV.

Figure 13:  The same as Fig. 12, but the solutions for Eq.(37) with $E_{\chi_{2}}$ in Eq.(41) of $e^{+}e^{-}\rightarrow Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$ process.  BP | $N_{phys}$ | $(M_{\chi_{2}},\,M_{\chi_{1}})^{\textrm{true}}$ | rms
---|---|---|---
$(M_{\chi_{2}},\,M_{\chi_{1}})^{\textrm{peak}}$
BP1 | $4473$ | $(0.42,\,0.30)$ | $(0.168,\,0.175)$
$(0.43,\,0.32)$
BP2 | $4915$ | $(3.30,\,3.00)$ | $(0.175,\,0.190)$
$(3.30,\,3.00)$
BP3 | $4856$ | $(1.40,\,1.00)$ | $(0.172,\,0.192)$
$(1.40,\,1.00)$
BP4 | $4918$ | $(2.40,\,2.00)$ | $(0.155,\,0.170)$
$(2.40,\,2.00)$
Table 4: The peak measured values and root mean square (rms) for $(M_{\chi_{2}},M_{\chi_{1}})$ of four benchmark points in $e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$ process. $N_{phys}$ is the number of physical solutions from 100 signal events. All numbers are in GeV unit.  BP | $N_{phys}$ | $(M_{\chi_{2}},\,M_{\chi_{1}})^{\textrm{true}}$ | rms
---|---|---|---
$(M_{\chi_{2}},\,M_{\chi_{1}})^{\textrm{peak}}$
BP1 | $901$ | $(0.42,\,0.30)$ | $(0.114,\,0.138)$
$(0.47,\,0.35)$
BP2 | $4914$ | $(3.30,\,3.00)$ | $(0.121,\,0.128)$
$(3.30,\,3.00)$
BP3 | $377$ | $(1.40,\,1.00)$ | $(0.216,\,0.402)$
$(1.41,\,1.01)$
BP4 | $2824$ | $(2.40,\,2.00)$ | $(0.126,\,0.173)$
$(2.40,\,2.00)$
Table 5: The same as Table. 4, but in $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
process.

In order to make the Kinematic Focus Point Method fit to reality, we involve
the detector resolution effects from Sec. IV.1 and event selections in Table.
1. Here we conservatively assume there are 100 signal events for each BP can
be recorded. Since we can solve the above kinematic equations for any two
signal events, we will get $C^{100}_{2}=4950$ solutions from 100 signal
events. After removing the unphysical solutions, we show the distributions of
these solutions on the $(M_{\chi_{1}},M_{\chi_{2}})$ plane for four BPs of
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
and $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
processes in Figs. 12 and 13, respectively. Here the bin size for both $x$ and
$y$ axes are set to be $10$ MeV. Moreover, the most probable mass values and
their statistical errors for four BPs in our simulation are summarized in
Table. 4 and 5. The statistical errors are estimated as the root mean square
(rms) value with respective to the most probable values
$rms=\sqrt{\sum^{N_{phys}}_{i=1}(M_{i}-M^{peak})^{2}/N_{phys}}$, where
$N_{phys}$ is the number of physical solutions from 100 signal events and
$M^{peak}$ is the most probable value.

In Fig. 12, except for BP1 on top left panel, $M_{\chi_{1}}$ and
$M_{\chi_{2}}$ can be determined within the percentage of deviation for other
three BPs. As shown in Fig. 4, both the decay length of $\chi_{2}$ is relative
short and $E_{ll}$ is relative large for BP1 compared with other three BPs
such that events in BP1 cannot be avoided to severely suffer from detector
resolution effects. This explains why the physical solutions are reduced and
the reconstructed $M_{\chi_{1}}$ and $M_{\chi_{2}}$ are shifted to larger
values with higher smearing for BP1. In Fig. 13, $M_{\chi_{1}}$ and
$M_{\chi_{2}}$ can still be determined within the percentage of deviation for
BP2 and BP4. However, apart from detector resolution effects from charged
leptons and displace vertex, the ability to precisely pin down $M_{\chi_{1}}$
and $M_{\chi_{2}}$ for the process in Eq.(29) also relies on how well the on-
shell $Z^{\prime}$ can be reconstructed. As shown in Fig. 5, the peak of
$m_{Z^{\prime}}$ is rather spread out (or even shifted) for BP1 and BP3 than
for BP2 and BP4. According to Eq.(39), the reconstructed $Z^{\prime}$ four-
momentum is highly dependent on the detector resolution effects of ISR photon.
Once the $E(\gamma)$ increases, its smearing is also enhanced such that the
reconstruction of the on-shell $Z^{\prime}$ becomes poor. In the end, the
physical solutions of BP1 and BP3 are largely reduced and shifting behaviors
of $M_{\chi_{1}}$ and $M_{\chi_{2}}$ become severe in Fig. 13 compared with
the ones in Fig. 12. This is another reason to encourage our experimental
colleagues to search for not only the usual process in Eq.(29) at B-factories,
but also the process in Eq.(28) with the displaced vertex trigger which is
more sensitive to determine the mass and mass splitting of DM sector for some
parameter space in the inelastic DM models.

##  V Conclusion

The dark sectors with light DM candidate are still all the rage as the null
signal result from DM direct detetion has been reported until now. On the
other hand, the quietness of beyond the SM searches at the LHC also hints new
physics signatures may hide in the elusive corner. New kind of signatures such
as long-lived particles at colliders become more and more popular and may
guide a royal road for new physics evidences. Among these dark sector models,
the inelastic DM model is an appealing example which can both allow light DM
candidate and predict long-lived particle which can be searched for at
colliders. Therefore, we focus on the displaced vertex signatures in inelastic
DM models at Belle II in this work.

In the inelastic DM models, if the mass splitting between DM excited and
ground states is small enough, the co-annihilation becomes the dominant
channel for thermal relic density and the DM excited state can be long-lived
at the collider scale. We first review scalar and fermion inelastic DM models
with $U(1)_{D}$ gauge symmetry in Sec. II and point out the dark Higgs sector
caused the mass splitting between DM excited and ground states and provided a
mechanism to the $Z^{\prime}$ mass. Besides, the off-diagonal interaction
between $Z^{\prime}$ and DM sector is derived from the covariant derivative
term of DM sector.

The analytical representations and numerical results for cross sections of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ and
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$ processes are
studied in Sec. III. In the first process, there is an extra $\beta$ factor
suppression for boson pair cross sections compared with fermion pair ones. In
the second process, the dominant channel comes from $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\gamma$.
Therefore, cross sections in this process for scalar and fermion inelatic DM
models are very close to each other. In addition, the polar angle
distributions of $\phi_{2}(\chi_{2})$ at the Belle II LAB frame as shown in
Fig. 2 are helpful to distinguish these two models.

The novel dilepton displaced vertex signatures at Belle II are the main
targets of this study. We include detector resolution effects in Sec. IV.1 and
event selections in Table. 1 for signal processes in Eq.(28) and (29). In Fig.
6, we display four BPs of scalar and fermion inelastic DM models and these two
models can be discriminated even with an integrated luminosity less than
$1~{}\rm{ab}^{-1}$. Our analysis results are summarized in Figs. 8 and 9,
which indicate that the future bounds from dilepton displaced vertex searches
in inelastic DM models at Belle II are stronger than previous constraints for
two benchmark mass windows $0.4\lesssim M_{\phi_{1},\chi_{1}}\lesssim 4.9$ GeV
and $0.1\lesssim M_{\phi_{1},\chi_{1}}\lesssim 1.1$ GeV, respectively. For the
latter case, our result can cover the allowed parameter space which can
explain the muon $(g-2)_{\mu}$. Therefore the early stage of the Belle II
experiment can explicitly close off this area in the parameter space. We
further apply the Kinematic Focus Point method for
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}$
and $e^{+}e^{-}\rightarrow
Z^{\prime}\gamma\rightarrow\chi_{1}\chi_{2}\gamma\rightarrow\chi_{1}\chi_{1}l^{+}l^{-}\gamma$
processes to determine both $M_{\chi_{1}}$ and $M_{\chi_{2}}$ masses.
Especially, we find the mass and mass splitting of DM sector can be determined
within the percentage of deviation as shown in Table. 4 and 5.

## Acknowledgment

We thank Myeonghun Park and Youngjoon Kwon for very helpful discussions and
useful information. The work is supported in part by KIAS Individual Grants,
Grant No. PG076201 (DWK), PG021403 (PK) and PG075301 (CTL) at Korea Institute
for Advanced Study, and by National Research Foundation of Korea (NRF) Grant
No. NRF-2019R1A2C3005009 (PK), funded by the Korean government (MSIT).

##  Appendix A The analytical representations for cross sections of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ in the LAB frame

In this appendix, we show the analytical representations for cross sections of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})$ in the LAB frame
which have been mentioned in Sec. III. We first set the four momentum of the
initial $e^{+}$, $e^{-}$ and $\phi_{2}(\chi_{2})$, $\phi_{1}(\chi_{1})$ in the
LAB frame as

$\displaystyle p_{e^{+}}=(E_{+},0,0,E_{+}),$ $\displaystyle
p_{e^{-}}=(E_{-},0,0,-E_{-}),$ $\displaystyle
p_{\phi_{2}(\chi_{2})}=(E_{\phi_{2}(\chi_{2})},p_{x,\phi_{2}(\chi_{2})},p_{y,\phi_{2}(\chi_{2})},p_{z,\phi_{2}(\chi_{2})}),$
$\displaystyle
p_{\phi_{1}(\chi_{1})}=(E_{\phi_{1}(\chi_{1})},p_{x,\phi_{1}(\chi_{1})},p_{y,\phi_{1}(\chi_{1})},p_{z,\phi_{1}(\chi_{1})}).$
(42)

According to the four momentum conservation for this scattering process,
$p_{\phi_{1}(\chi_{1})}$ can be written as

$p_{\phi_{1}(\chi_{1})}=(E_{+}+E_{-}-E_{\phi_{2}(\chi_{2})},-p_{x,\phi_{2}(\chi_{2})},-p_{y,\phi_{2}(\chi_{2})},E_{+}-E_{-}-p_{z,\phi_{2}(\chi_{2})}).$
(43)

The $d\sigma/d\cos\theta$ for $e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}$ via
s-channel $Z^{\prime}$ in the LAB frame can be written as

$\displaystyle\left.\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})}{d\cos\theta}\right|_{LAB}=$
$\displaystyle\frac{\epsilon^{2}e^{2}g^{2}_{D}}{256\pi
E_{+}E^{2}_{-}[(4E_{+}E_{-}-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}]}\times$
$\displaystyle[8E_{+}E_{-}\left(E_{\phi_{2}}(E_{+}+E_{-})+(E_{-}-E_{+})|\overrightarrow{p_{\phi_{2}}}|\cos\theta\right)$
$\displaystyle-4E_{+}E_{-}(M^{2}_{\phi_{1}}+3M^{2}_{\phi_{2}})+(M^{2}_{\phi_{2}}-M^{2}_{\phi_{1}})^{2}$
$\displaystyle+\left(4E_{+}(E_{\phi_{2}}-E_{-}-|\overrightarrow{p_{\phi_{2}}}|\cos\theta)+M^{2}_{\phi_{1}}-M^{2}_{\phi_{2}}\right)\times$
$\displaystyle\left(4E_{-}(E_{\phi_{2}}-E_{+}+|\overrightarrow{p_{\phi_{2}}}|\cos\theta)+M^{2}_{\phi_{1}}-M^{2}_{\phi_{2}}\right)]\times$
$\displaystyle\left(-\frac{dE_{\phi_{2}}}{d\cos\theta}+\cos\theta\frac{d|\overrightarrow{p_{\phi_{2}}}|}{d\cos\theta}-\sin\theta|\overrightarrow{p_{\phi_{2}}}|\right),$
(44)

where
$|\overrightarrow{p_{\phi_{2}}}|=\sqrt{E^{2}_{\phi_{2}}-M^{2}_{\phi_{2}}}$ and
$E_{\phi_{2}}$ is the same as Eq.(35) with $\chi_{1,2}\rightarrow\phi_{1,2}$.

On the other hand, the $d\sigma/d\cos\theta$ for
$e^{+}e^{-}\rightarrow\chi_{1}\chi_{2}$ via s-channel $Z^{\prime}$ in the LAB
frame can be written as

$\displaystyle\left.\frac{d\sigma(e^{+}e^{-}\rightarrow\chi_{1}\chi_{2})}{d\cos\theta}\right|_{LAB}=$
$\displaystyle\frac{\epsilon^{2}e^{2}g^{2}_{D}}{32\pi
E_{+}E^{2}_{-}[(4E_{+}E_{-}-m^{2}_{Z^{\prime}})^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}]}\times$
$\displaystyle[2(E^{2}_{+}+E^{2}_{-})E^{2}_{\chi_{2}}(1+\cos^{2}\theta)-(E_{-}-E_{+})(M^{2}_{\chi_{2}}-M^{2}_{\chi_{1}})|\overrightarrow{p_{\chi_{2}}}|\cos\theta$
$\displaystyle+(E_{+}+E_{-})\left(4(E_{-}-E_{+})|\overrightarrow{p_{\chi_{2}}}|\cos\theta+M^{2}_{\chi_{1}}-M^{2}_{\chi_{2}}\right)E_{\chi_{2}}$
$\displaystyle-2(E^{2}_{+}+E^{2}_{-})M^{2}_{\chi_{2}}\cos^{2}\theta+4E_{+}E_{-}M_{\chi_{1}}M_{\chi_{2}}]\times$
$\displaystyle\left(-\frac{dE_{\phi_{2}}}{d\cos\theta}+\cos\theta\frac{d|\overrightarrow{p_{\phi_{2}}}|}{d\cos\theta}-\sin\theta|\overrightarrow{p_{\phi_{2}}}|\right),$
(45)

where
$|\overrightarrow{p_{\chi_{2}}}|=\sqrt{E^{2}_{\chi_{2}}-M^{2}_{\chi_{2}}}$ and
$E_{\chi_{2}}$ can be found in Eq. (35).

##  Appendix B The full analytical representations of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\rightarrow\phi_{1}\phi_{1}(\chi_{1}\chi_{1})f\overline{f}$

In this appendix, we display the full analytical representations of
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}(\chi_{1}\chi_{2})\rightarrow\phi_{1}\phi_{1}(\chi_{1}\chi_{1})f\overline{f}$
as mentioned in Sec. III. The differential cross section for
$e^{+}(p_{1})e^{-}(p_{2})\rightarrow\phi_{1}(p_{3})\phi_{1}(p_{4})f(p_{5})\overline{f}(p_{6})$
can be represented as

$d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{1}f\overline{f})=\frac{(2\pi)^{4}|{\cal
M}_{2\rightarrow 4}|^{2}}{2s}\times
d\Phi_{4}(p_{1}+p_{2};p_{3},p_{4},p_{5},p_{6}),$ (46)

where ${\cal M}_{2\rightarrow 4}$ is the scattering amplitude for
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{1}f\overline{f}$, $s$ is the square of
centre-of-mass energy and $\Phi_{4}(p_{1}+p_{2};p_{3},p_{4},p_{5},p_{6})$ is
the four-body phase space. Since $\phi_{2}$ is on-shell produced in this
process, we apply the narrow width approximation for on-shell $\phi_{2}$ to
separate $e^{+}e^{-}\rightarrow\phi_{1}\phi_{1}f\overline{f}$ to
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}$ with
$\phi_{2}\rightarrow\phi_{1}f\overline{f}$. The differential phase space
$d\Phi_{4}(p_{1}+p_{2};p_{3},p_{4},p_{5},p_{6})$ can be expanded as

$d\Phi_{4}(p_{1}+p_{2};p_{3},p_{4},p_{5},p_{6})=d\Phi_{2}(p_{1}+p_{2};p_{3},l)\times
d\Phi_{3}(l;p_{4},p_{5},p_{6})(2\pi)^{3}dl^{2},$ (47)

where $l$ is the four momentum of on-shell $\phi_{2}$. Therefore, we can
disassemble the differential cross section in the following form,

$d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{1}f\overline{f})=d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})dB(\phi_{2}\rightarrow\phi_{1}f\overline{f}),$
(48)

where

$d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})=\frac{(2\pi)^{4}|{\cal
M}_{2\rightarrow 2}|^{2}}{2s}\times d\Phi_{2}(p_{1}+p_{2};p_{3},l),$ (49)

and

$\displaystyle dB(\phi_{2}\rightarrow\phi_{1}f\overline{f})$
$\displaystyle=\frac{2\pi}{\Gamma_{\phi_{2}}}\delta(l^{2}-M^{2}_{\phi_{2}})\frac{1}{2M_{\phi_{2}}}|{\cal
M}_{1\rightarrow 3}|^{2}d\Phi_{3}(l;p_{4},p_{5},p_{6})(2\pi)^{3}dl^{2}$
$\displaystyle=\frac{d\Gamma_{\phi_{2}\rightarrow\phi_{1}f\overline{f}}}{\Gamma_{\phi_{2}}}\delta(l^{2}-M^{2}_{\phi_{2}})(2\pi)^{4}dl^{2},$
(50)

where $\Gamma_{\phi_{2}}$ is the total width of $\phi_{2}$, ${\cal
M}_{2\rightarrow 2}$ is the scattering amplitude for
$e^{+}e^{-}\rightarrow\phi_{1}\phi_{2}$ and ${\cal M}_{1\rightarrow 3}$ is the
decay amplitude for $\phi_{2}\rightarrow\phi_{1}f\overline{f}$. After applying
the narrow width approximation for on-shell $\phi_{2}$, $|{\cal
M}_{2\rightarrow 4}|^{2}$ can be transferred to

$|{\cal M}_{2\rightarrow 4}|^{2}=|{\cal M}_{2\rightarrow
2}|^{2}\times\frac{\pi\delta(l^{2}-M^{2}_{\phi_{2}})}{M_{\phi_{2}}\Gamma_{\phi_{2}}}|{\cal
M}_{1\rightarrow 3}|^{2}.$ (51)

Furthermore, we integrate out some phase space variables and leave the
following kinematic variables :

  * •

$\cos\theta$ : $\theta$ is the direction of $\phi_{2}$ relative to the
positron beam direction.

  * •

$x_{-}$ : $x_{-}=2E_{f}/M_{\phi_{2}}$ in the $\phi_{2}$ rest frame.

  * •

$x_{+}$ : $x_{+}=2E_{\overline{f}}/M_{\phi_{2}}$ in the $\phi_{2}$ rest frame.

for the representation of differential cross section,

$\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{1}f\overline{f})}{d\cos\theta
dx_{-}dx_{+}}=\frac{d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})}{d\cos\theta}\frac{dB(\phi_{2}\rightarrow\phi_{1}f\overline{f})}{dx_{-}dx_{+}}.$
(52)

Notice the opening angle between the fermion pair can be written as

$\cos\theta_{f\overline{f}}=1-2(x_{-}+x_{+}-1+M^{2}_{\phi_{1}}/M^{2}_{\phi_{2}})/(x_{-}x_{+}).$
(53)

The $d\sigma(e^{+}e^{-}\rightarrow\phi_{1}\phi_{2})/d\cos\theta$ can be found
in Eq.(22) in the CM frame and Eq.(44) in the LAB frame and
$d\Gamma(\phi_{2}\rightarrow\phi_{1}f\overline{f})/dx_{-}/dx_{+}$ can be
written as

$\frac{d\Gamma(\phi_{2}\rightarrow\phi_{1}f\overline{f})}{dx_{-}dx_{+}}=\frac{g^{2}_{D}\epsilon^{2}\alpha
M^{3}_{\phi_{2}}}{8\pi^{2}}\frac{(1-x_{+})(1-x_{-})M^{2}_{\phi_{2}}-(M_{\phi_{2}}-\Delta_{\phi})^{2}}{\left[M^{2}_{\phi_{2}}(x_{+}+x_{-}-1)-m^{2}_{Z^{\prime}}\right]^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}}.$
(54)

One can find $x_{+}$ and $x_{-}$ are symmetric in this formula.

The associated formulae for
$e^{+}(p_{1})e^{-}(p_{2})\rightarrow\chi_{1}(p_{3})\chi_{1}(p_{4})f(p_{5})\overline{f}(p_{6})$
can be obtained from Eq.(46) to (52) by replacing $\phi_{1,2}$ to
$\chi_{1,2}$. The $d\sigma(e^{+}e^{-}\rightarrow\chi_{1}\chi_{2})/d\cos\theta$
can be found in Eq.(25) in the CM frame and Eq.(45) in the LAB frame and
$d\Gamma(\chi_{2}\rightarrow\chi_{1}f\overline{f})/dx_{-}/dx_{+}$ can be
written as

$\displaystyle\frac{d\Gamma(\chi_{2}\rightarrow\chi_{1}f\overline{f})}{dx_{-}dx_{+}}=$
$\displaystyle\frac{g^{2}_{D}\epsilon^{2}\alpha
M^{3}_{\chi_{2}}}{16\pi^{2}}\frac{1}{\left[M^{2}_{\chi_{2}}(x_{+}+x_{-}-1)-m^{2}_{Z^{\prime}}\right]^{2}+m^{2}_{Z^{\prime}}\Gamma^{2}_{Z^{\prime}}}$
$\displaystyle\times\\{-M^{4}_{\chi_{2}}\left[x^{2}_{+}+x^{2}_{-}+2(x_{+}+x_{-})\right]+4M^{3}_{\chi_{2}}\Delta_{\chi}(1+x_{+}+x_{-})$
$\displaystyle-M^{2}_{\chi_{2}}\left[(6+x_{+}+x_{-})\Delta^{2}_{\chi}+2(6+x_{+}+x_{-})m^{2}_{f}\right]$
$\displaystyle+2M_{\chi_{2}}(\Delta^{3}_{\chi}+4\Delta_{\chi}m^{2}_{f})-m^{2}_{f}(\Delta^{2}_{\chi}+4m^{2}_{f})\\}.$
(55)

## References

  * (1) C. S. Frenk and S. D. M. White, Annalen Phys. 524, 507-534 (2012) doi:10.1002/andp.201200212 [arXiv:1210.0544 [astro-ph.CO]].
  * (2) N. Aghanim et al. [Planck], Astron. Astrophys. 641, A6 (2020) doi:10.1051/0004-6361/201833910 [arXiv:1807.06209 [astro-ph.CO]].
  * (3) M. Bauer and T. Plehn, Lect. Notes Phys. 959, pp. (2019) doi:10.1007/978-3-030-16234-4 [arXiv:1705.01987 [hep-ph]].
  * (4) E. Aprile et al. [XENON], Phys. Rev. Lett. 121, no.11, 111302 (2018) doi:10.1103/PhysRevLett.121.111302 [arXiv:1805.12562 [astro-ph.CO]].
  * (5) X. Ren et al. [PandaX-II], Phys. Rev. Lett. 121, no.2, 021304 (2018) doi:10.1103/PhysRevLett.121.021304 [arXiv:1802.06912 [hep-ph]].
  * (6) E. Aprile et al. [XENON], Phys. Rev. Lett. 122, no.14, 141301 (2019) doi:10.1103/PhysRevLett.122.141301 [arXiv:1902.03234 [astro-ph.CO]].
  * (7) S. Knapen, T. Lin and K. M. Zurek, Phys. Rev. D 96, no.11, 115021 (2017) doi:10.1103/PhysRevD.96.115021 [arXiv:1709.07882 [hep-ph]].
  * (8) T. Lin, PoS 333, 009 (2019) doi:10.22323/1.333.0009 [arXiv:1904.07915 [hep-ph]].
  * (9) M. Ablikim et al. [BESIII], Nucl. Instrum. Meth. A 614, 345-399 (2010) doi:10.1016/j.nima.2009.12.050 [arXiv:0911.4960 [physics.ins-det]].
  * (10) B. Aubert et al. [BaBar], Nucl. Instrum. Meth. A 729, 615-701 (2013) doi:10.1016/j.nima.2013.05.107 [arXiv:1305.3560 [physics.ins-det]].
  * (11) I. Adachi et al. [Belle-II], Nucl. Instrum. Meth. A 907, 46-59 (2018) doi:10.1016/j.nima.2018.03.068
  * (12) E. Kou et al. [Belle-II], PTEP 2019, no.12, 123C01 (2019) [erratum: PTEP 2020, no.2, 029201 (2020)] doi:10.1093/ptep/ptz106 [arXiv:1808.10567 [hep-ex]].
  * (13) J. D. Bjorken, R. Essig, P. Schuster and N. Toro, Phys. Rev. D 80, 075018 (2009) doi:10.1103/PhysRevD.80.075018 [arXiv:0906.0580 [hep-ph]].
  * (14) A. Boveia and C. Doglioni, Ann. Rev. Nucl. Part. Sci. 68, 429-459 (2018) doi:10.1146/annurev-nucl-101917-021008 [arXiv:1810.12238 [hep-ex]].
  * (15) M. Pospelov, A. Ritz and M. B. Voloshin, Phys. Lett. B 662, 53-61 (2008) doi:10.1016/j.physletb.2008.02.052 [arXiv:0711.4866 [hep-ph]].
  * (16) N. Arkani-Hamed, D. P. Finkbeiner, T. R. Slatyer and N. Weiner, Phys. Rev. D 79, 015014 (2009) doi:10.1103/PhysRevD.79.015014 [arXiv:0810.0713 [hep-ph]].
  * (17) M. Pospelov, Phys. Rev. D 80, 095002 (2009) doi:10.1103/PhysRevD.80.095002 [arXiv:0811.1030 [hep-ph]].
  * (18) S. Baek, P. Ko and W. I. Park, JHEP 07, 013 (2013) doi:10.1007/JHEP07(2013)013 [arXiv:1303.4280 [hep-ph]].
  * (19) B. Holdom, Phys. Lett. B 166, 196-198 (1986) doi:10.1016/0370-2693(86)91377-8
  * (20) E. Izaguirre, G. Krnjaic and B. Shuve, Phys. Rev. D 93, no.6, 063523 (2016) doi:10.1103/PhysRevD.93.063523 [arXiv:1508.03050 [hep-ph]].
  * (21) M. Duerr, T. Ferber, C. Hearty, F. Kahlhoefer, K. Schmidt-Hoberg and P. Tunney, JHEP 02, 039 (2020) doi:10.1007/JHEP02(2020)039 [arXiv:1911.03176 [hep-ph]].
  * (22) Y. Zhang, W. T. Zhang, M. Song, X. A. Pan, Z. M. Niu and G. Li, Phys. Rev. D 100, no.11, 115016 (2019) doi:10.1103/PhysRevD.100.115016 [arXiv:1907.07046 [hep-ph]].
  * (23) J. P. Lees et al. [BaBar], Phys. Rev. Lett. 119, no.13, 131804 (2017) doi:10.1103/PhysRevLett.119.131804 [arXiv:1702.03327 [hep-ex]].
  * (24) G. Mohlabeng, Phys. Rev. D 99, no.11, 115001 (2019) doi:10.1103/PhysRevD.99.115001 [arXiv:1902.05075 [hep-ph]].
  * (25) D. Tucker-Smith and N. Weiner, Phys. Rev. D 64, 043502 (2001) doi:10.1103/PhysRevD.64.043502 [arXiv:hep-ph/0101138 [hep-ph]].
  * (26) S. Baek, P. Ko and W. I. Park, Phys. Lett. B 747, 255-259 (2015) doi:10.1016/j.physletb.2015.06.002 [arXiv:1407.6588 [hep-ph]].
  * (27) P. Ko, T. Matsui and Y. L. Tang, JHEP 10, 082 (2020) doi:10.1007/JHEP10(2020)082 [arXiv:1910.04311 [hep-ph]].
  * (28) S. Baek, J. Kim and P. Ko, Phys. Lett. B 810, 135848 (2020) doi:10.1016/j.physletb.2020.135848 [arXiv:2006.16876 [hep-ph]].
  * (29) K. Griest and D. Seckel, Phys. Rev. D 43, 3191-3203 (1991) doi:10.1103/PhysRevD.43.3191
  * (30) T. R. Slatyer, Phys. Rev. D 93, no.2, 023527 (2016) doi:10.1103/PhysRevD.93.023527 [arXiv:1506.03811 [hep-ph]].
  * (31) K. Harigaya, Y. Nakai and M. Suzuki, Phys. Lett. B 809, 135729 (2020) doi:10.1016/j.physletb.2020.135729 [arXiv:2006.11938 [hep-ph]].
  * (32) N. F. Bell, J. B. Dent, B. Dutta, S. Ghosh, J. Kumar and J. L. Newstead, Phys. Rev. Lett. 125, no.16, 161803 (2020) doi:10.1103/PhysRevLett.125.161803 [arXiv:2006.12461 [hep-ph]].
  * (33) H. M. Lee, [arXiv:2006.13183 [hep-ph]].
  * (34) M. Baryakhtar, A. Berlin, H. Liu and N. Weiner, [arXiv:2006.13918 [hep-ph]].
  * (35) J. Bramante and N. Song, Phys. Rev. Lett. 125, no.16, 161805 (2020) doi:10.1103/PhysRevLett.125.161805 [arXiv:2006.14089 [hep-ph]].
  * (36) H. An and D. Yang, [arXiv:2006.15672 [hep-ph]].
  * (37) R. Catena, J. Conrad, C. Döring, A. D. Ferella and M. B. Krauss, Phys. Rev. D 97, no.2, 023007 (2018) doi:10.1103/PhysRevD.97.023007 [arXiv:1706.09471 [hep-ph]].
  * (38) H. C. Cheng, Z. Han, I. W. Kim and L. T. Wang, JHEP 11, 122 (2010) doi:10.1007/JHEP11(2010)122 [arXiv:1008.0405 [hep-ph]].
  * (39) T. Melia, JHEP 01, 143 (2012) doi:10.1007/JHEP01(2012)143 [arXiv:1110.6185 [hep-ph]].
  * (40) L. Edelhauser, K. T. Matchev and M. Park, JHEP 11, 006 (2012) doi:10.1007/JHEP11(2012)006 [arXiv:1205.2054 [hep-ph]].
  * (41) N. D. Christensen and D. Salmon, Phys. Rev. D 90, no.1, 014025 (2014) doi:10.1103/PhysRevD.90.014025 [arXiv:1311.6465 [hep-ph]].
  * (42) S. Y. Choi, Phys. Rev. D 98, no.11, 115037 (2018) doi:10.1103/PhysRevD.98.115037 [arXiv:1811.10377 [hep-ph]].
  * (43) W. Abdallah, A. Hammad, S. Khalil and S. Moretti, Phys. Rev. D 100, no.9, 095006 (2019) doi:10.1103/PhysRevD.100.095006 [arXiv:1907.08358 [hep-ph]].
  * (44) C. L. Shan, New J. Phys. 11, 105013 (2009) doi:10.1088/1367-2630/11/10/105013 [arXiv:0903.4320 [hep-ph]].
  * (45) B. J. Kavanagh and A. M. Green, Phys. Rev. Lett. 111, no.3, 031302 (2013) doi:10.1103/PhysRevLett.111.031302 [arXiv:1303.6868 [astro-ph.CO]].
  * (46) H. C. Cheng and J. Gu, JHEP 10, 094 (2011) doi:10.1007/JHEP10(2011)094 [arXiv:1109.3471 [hep-ph]].
  * (47) L. A. Harland-Lang, C. H. Kom, K. Sakurai and W. J. Stirling, Eur. Phys. J. C 72, 1969 (2012) doi:10.1140/epjc/s10052-012-1969-2 [arXiv:1110.4320 [hep-ph]].
  * (48) L. A. Harland-Lang, C. H. Kom, K. Sakurai and W. J. Stirling, Phys. Rev. Lett. 108, 181805 (2012) doi:10.1103/PhysRevLett.108.181805 [arXiv:1202.0047 [hep-ph]].
  * (49) A. C. Kobach, Phys. Rev. D 88, 116001 (2013) doi:10.1103/PhysRevD.88.116001 [arXiv:1308.5671 [hep-ph]].
  * (50) L. A. Harland-Lang, C. H. Kom, K. Sakurai and M. Tonini, JHEP 06, 175 (2014) doi:10.1007/JHEP06(2014)175 [arXiv:1312.5720 [hep-ph]].
  * (51) N. D. Christensen, T. Han, Z. Qian, J. Sayre, J. Song and Stefanus, Phys. Rev. D 90, 114029 (2014) doi:10.1103/PhysRevD.90.114029 [arXiv:1404.6258 [hep-ph]].
  * (52) P. Konar and A. K. Swain, Phys. Rev. D 93, no.1, 015021 (2016) doi:10.1103/PhysRevD.93.015021 [arXiv:1509.00298 [hep-ph]].
  * (53) P. Ko and H. Yokoya, JHEP 08 (2016), 109 doi:10.1007/JHEP08(2016)109 [arXiv:1603.04737 [hep-ph]].
  * (54) Q. F. Xiang, X. J. Bi, Q. S. Yan, P. F. Yin and Z. H. Yu, Phys. Rev. D 95, no.7, 075037 (2017) doi:10.1103/PhysRevD.95.075037 [arXiv:1610.03372 [hep-ph]].
  * (55) T. Kamon, P. Ko and J. Li, Eur. Phys. J. C 77, no.9, 652 (2017) doi:10.1140/epjc/s10052-017-5240-8 [arXiv:1705.02149 [hep-ph]].
  * (56) B. Dutta, T. Kamon, P. Ko and J. Li, Eur. Phys. J. C 78, no.7, 595 (2018) doi:10.1140/epjc/s10052-018-6071-y [arXiv:1712.05123 [hep-ph]].
  * (57) Z. Flowers, Q. Meier, C. Rogan, D. W. Kang and S. C. Park, JHEP 03, 132 (2020) doi:10.1007/JHEP03(2020)132 [arXiv:1903.05825 [hep-ph]].
  * (58) D. Kim, K. T. Matchev and P. Shyamsundar, JHEP 10, 154 (2019) doi:10.1007/JHEP10(2019)154 [arXiv:1906.02821 [hep-ph]].
  * (59) S. Banerjee, B. Bhattacherjee, A. Goudelis, B. Herrmann, D. Sengupta and R. Sengupta, [arXiv:1912.06669 [hep-ph]].
  * (60) K. J. Bae, M. Park and M. Zhang, Phys. Rev. D 101, no.11, 115036 (2020) doi:10.1103/PhysRevD.101.115036 [arXiv:2001.02142 [hep-ph]].
  * (61) M. Duerr, T. Ferber, C. Garcia-Cely, C. Hearty and K. Schmidt-Hoberg, [arXiv:2012.08595 [hep-ph]].
  * (62) C. Boehm, M. J. Dolan and C. McCabe, JCAP 12, 027 (2012) doi:10.1088/1475-7516/2012/12/027 [arXiv:1207.0497 [astro-ph.CO]].
  * (63) P. F. Depta, M. Hufnagel, K. Schmidt-Hoberg and S. Wild, JCAP 04, 029 (2019) doi:10.1088/1475-7516/2019/04/029 [arXiv:1901.06944 [hep-ph]].
  * (64) G. F. Giudice, D. Kim, J. C. Park and S. Shin, Phys. Lett. B 780, 543-552 (2018) doi:10.1016/j.physletb.2018.03.043 [arXiv:1712.07126 [hep-ph]].
  * (65) P. Ilten, Y. Soreq, M. Williams and W. Xue, JHEP 06, 004 (2018) doi:10.1007/JHEP06(2018)004 [arXiv:1801.04847 [hep-ph]].
  * (66) J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli and M. Zaro, JHEP 07, 079 (2014) doi:10.1007/JHEP07(2014)079 [arXiv:1405.0301 [hep-ph]].
  * (67) P. A. Zyla et al. [Particle Data Group], PTEP 2020, no.8, 083C01 (2020) doi:10.1093/ptep/ptaa104
  * (68) A. Birkedal, K. Matchev and M. Perelstein, Phys. Rev. D 70, 077701 (2004) doi:10.1103/PhysRevD.70.077701 [arXiv:hep-ph/0403004 [hep-ph]].
  * (69) A. Alloul, N. D. Christensen, C. Degrande, C. Duhr and B. Fuks, Comput. Phys. Commun. 185, 2250-2300 (2014) doi:10.1016/j.cpc.2014.04.012 [arXiv:1310.1921 [hep-ph]].
  * (70) V. V. Gligorov, S. Knapen, M. Papucci and D. J. Robinson, Phys. Rev. D 97, no.1, 015023 (2018) doi:10.1103/PhysRevD.97.015023 [arXiv:1708.09395 [hep-ph]].
  * (71) A. Berlin, S. Gori, P. Schuster and N. Toro, Phys. Rev. D 98, no.3, 035011 (2018) doi:10.1103/PhysRevD.98.035011 [arXiv:1804.00661 [hep-ph]].
  * (72) J. L. Feng, I. Galon, F. Kling and S. Trojanowski, Phys. Rev. D 97, no.3, 035001 (2018) doi:10.1103/PhysRevD.97.035001 [arXiv:1708.09389 [hep-ph]].
  * (73) J. P. Chou, D. Curtin and H. J. Lubatti, Phys. Lett. B 767, 29-36 (2017) doi:10.1016/j.physletb.2017.01.043 [arXiv:1606.06298 [hep-ph]].
  * (74) A. Hook, E. Izaguirre and J. G. Wacker, Adv. High Energy Phys. 2011, 859762 (2011) doi:10.1155/2011/859762 [arXiv:1006.0973 [hep-ph]].