# Dark matter-electron scattering and freeze-in scenarios in the light of
$Z^{\prime}$ mediation

Basabendu Barman<EMAIL_ADDRESS>Department of Physics, School of
Engineering and Sciences, SRM University-AP, Amaravati 522240, India  Arindam
Das<EMAIL_ADDRESS>Institute for the Advancement of Higher
Education, Hokkaido University, Sapporo 060-0817, Japan  Department of
Physics, Hokkaido University, Sapporo 060-0810, Japan  Sanjoy Mandal
<EMAIL_ADDRESS>Korea Institute for Advanced Study, Seoul 02455, Korea

###### Abstract

We investigate dark matter (DM)-electron scattering in a minimal $U(1)_{X}$
extension of the Standard Model (SM), where the DM can appear as a Majorana
fermion, a complex singlet scalar or a Dirac fermion. To study bounds on the
$U(1)_{X}$ gauge coupling $(g_{X})$ and new gauge boson mass
$(M_{Z^{\prime}})$, from DM-electron scattering, we consider several direct
search experiments like CDMS, DAMIC, SENSEI, PandaX-II, DarkSide-50 and
XENON1T-S2 for different $U(1)_{X}$ charges. In this set-up we consider DM
production via freeze-in both in radiation dominated and modified cosmological
background to project sensitivities on $g_{X}-M_{Z^{\prime}}$ plane satisfying
observed relic abundance. DM-electron scattering could provide comparable, or
even stronger bounds than those obtained from the electron/ muon $(g-2)$, low
energy scattering and intensity frontier experiments within 0.01 GeV $\lesssim
M_{Z^{\prime}}\lesssim$ 0.1 GeV. Constrains from freeze-in could provide
stronger sensitivities for $M_{Z^{\prime}}\gtrsim\mathcal{O}(1)$ GeV, however,
these limits are comparable to those obtained from LHCb, LEP experiments for
$\mathcal{O}(10)$ GeV $\lesssim M_{Z^{\prime}}\lesssim 150$ GeV. In future,
electron-muon scattering (MUonE), proton (FASER, DUNE) and electron/positron
(ILC) beam dump experiments could probe these parameters.

##  I Introduction

Two most prominent pieces of evidence of physics beyond the Standard Model
(SM) are the existence of neutrino mass that oscillation experiments support
ParticleDataGroup:2020ssz  and the dark matter (DM), whose presence is well-
established from several astophysical and cosmological Bertone:2016nfn ;
deSwart:2017heh ; Planck:2018vyg  observations. The origin of neutrino mass
can be explained in several extensions of the SM that employ the so-called
seesaw mechanisms Minkowski:1977sc ; Mohapatra:1979ia ; Gell-Mann:1979vob ;
Magg:1980ut ; Schechter:1980gr ; Wetterich:1981bx ; Lazarides:1980nt ;
Mohapatra:1980yp ; Foot:1988aq ; Ma:1998dn ; Ma:2002pf ; Hambye:2003rt ;
Bajc:2006ia  that can also include the DM itself Ma:2006km , thereby
addressing the two questions in one single framework. In context with DM, the
weakly interacting massive particles (WIMPs) are one of the most prominent
candidates Steigman:1984ac ; Jungman:1995df  (for a review, see, e.g., Refs.
Jungman:1995df ; Bertone:2004pz ; Feng:2010gw ). WIMPs are assumed to be in
thermal equilibrium at temperatures higher than its mass. Corresponding WIMP
relic density is determined by freeze-out, an epoch when the DM annihilation
rate can no longer keep up with the expansion rate of the Universe. However,
strong observational constraints (see, e.g., Refs. Roszkowski:2017nbc ;
Arcadi:2017kky ) on the typical WIMP parameter space motivate searches beyond
this paradigm. As an alternative, feebly interacting massive particles (FIMPs)
have been proposed McDonald:2001vt ; Hall:2009bx ; Bernal:2017kxu . In the
early Universe, FIMPs can be generated from either the decay or annihilation
of states in the visible sector. When the SM bath temperature becomes smaller
than the typical mass scale of the interaction (i.e., the maximum of the DM
and the mediator mass), the generation process becomes Boltzmann suppressed,
giving rise to a constant comoving DM number density. Such a scenario, as
opposed to freeze-out, is referred as the freeze-in Hall:2009bx . The FIMP
paradigm requires very suppressed interaction rates between the dark and
visible sectors, which can be achieved, in its infrared (IR) version, via
couplings several orders of magnitude weaker than that of weak interaction
coupling strength. Such feeble coupling, in general, are beyond the reach of
standard collider and DM scattering experiments.

Direct detection Goodman:1984dc ; Schumann:2019eaa  is a cornerstone of DM
search experiments, serving as a laboratory method to explore particle
interactions of DM in the local DM halo. The objective is to record the rare
instances when DM particles collide with a target material. The sensitivity of
these experiments at low DM masses is fundamentally constrained by the minimum
energy needed to produce an excitation in the material. As it has been argued
in Essig:2011nj ; Lin:2019uvt , direct detection experiments that search for
nuclear recoils caused by DM scattering are ineffective when it comes to DM in
the sub-GeV mass range, as in that case the nuclear recoil energy is far below
the lowest thresholds achieved in typical nuclear recoil experiments. Direct
detection of sub-GeV dark sectors thus requires a different approach, both
theoretically and experimentally. In this case, DM-electron scattering turns
out to be effective as the total energy available in the scattering is
significantly larger to trigger inelastic atomic processes that could lead to
visible signals. Several directions have been explored (see, for example Ref.
Battaglieri:2017aum ) to understand the capability of DM-electron searches in
probing or constraining physics beyond the SM.

From the perspective of particle physics model building, all these new physics
models typically extend the SM particle content with additional fields, and
introduce a discrete symmetry to ensure the stability of the DM over
cosmological timescales. There exist a class of models, where the generation
of neutrino mass and a particle DM can be simultaneously explained under the
same umbrella. The general but minimal $U(1)_{X}$ extension Das:2016zue  of
the SM is such a framework where an SM singlet scalar, charged under
$U(1)_{X}$ gauge group is introduced. The most intriguing aspect of this model
is that incorporating three generations of SM-singlet Right-Handed Neutrinos
(RHNs), as in the seesaw mechanism for generating light neutrino masses after
the breaking of $U(1)_{X}$ and electroweak symmetries, is not merely an option
but rather the simplest solution to eliminate all possible gauge anomalies.
Additionally, there is a new gauge boson $Z^{\prime}$ which interacts
differently with the left and right handed SM charged fermions manifesting a
chiral scenario.

Motivated from these, in this work we have explored the potential of DM-
electron scattering experiments in constraining the $U(1)_{X}$ gauge coupling
$(g_{X})$ and the corresponding mass of the new gauge boson
$(M_{Z^{\prime}})$, that emerges from the spontaneous breaking of $U(1)_{X}$.
We consider three benchmark models, namely, $U(1)_{\rm B-L},\,U(1)_{R}$,
chiral $U(1)_{X}$ and project bounds from several DM-electron scattering
experiments in $g_{X}-M_{Z^{\prime}}$ plane. Since the coupling involved in
this scenario is typically weaker than $SU(2)_{L}$ gauge coupling, we also
discuss DM production in all of the above cases via freeze-in. We take up
three different DM candidates: Majorana, scalar and Dirac, and obtain
corresponding bound in $g_{X}-M_{Z^{\prime}}$ plane to produce the observed DM
abundance via 2-to-2 scattering of the bath particles, mediated by
$Z^{\prime}$. Apart from DM genesis in the standard radiation dominated
background, we also explore the impact of a generalized modified cosmological
background on freeze-in yield. We summarize the status of the
$g_{X}-M_{Z^{\prime}}$ parameter space considering bounds from direct
detection, freeze-in and several other current and future energy and intensity
frontier experiments. Our analysis establishes not only the potential of DM-
electron scattering experiments in constraining parameters of $U(1)_{X}$ gauge
extension of the SM, but also the ability of laboratory and collider tests in
probing modified gravity/cosmological models in future.

The paper is organized as follows. The models are discussed in Sec. II. In
Sec. III, we lay out the formalism of calculating the DM-electron scattering
cross-section and event rates both for atomic and semi-conductor targets. The
freeze-in production of DM is elaborated for the case of radiation domination
and modified cosmology in Sec. IV. In Sec. V we summarize our main results and
show the resulting bounds on $g_{X}-M_{Z^{\prime}}$ plane coming from the
freeze-in and DM-electron scattering studies. Finally we conclude the paper in
Sec. VI.

##  II Models

As advocated in the introduction, we consider a general minimal $U(1)_{X}$
extension of the SM. In $SU(3)_{c}\otimes SU(2)_{L}\otimes U(1)_{Y}\otimes
U(1)_{X}$ framework, the transformation property of the SM lepton doublet
$(\ell_{L}^{i})$ is
$\\{1,2,-\frac{1}{2},x_{\ell}=-\frac{x_{H}}{2}-x_{\Phi}\\}$ and that of the
right handed electron $(e_{R}^{i})$ is $\\{1,1,-1,x_{e}=-x_{H}-x_{\Phi}\\}$.
In this set-up, the SM colored sector follow the transformation property of
quark doublet $(q_{L}^{i})$ as
$\\{3,2,\frac{1}{6},x_{q}=\frac{1}{6}x_{H}+\frac{1}{3}x_{\Phi}\\}$, right
handed up type-quark $(u_{R}^{i})$ as
$\\{3,1,\frac{2}{3},x_{u}=\frac{2}{3}x_{H}+\frac{1}{3}x_{\Phi}\\}$ and right
handed down type quark $(d_{R}^{i})$ as
$\\{3,1,-\frac{1}{3},x_{d}=-\frac{1}{3}x_{H}+\frac{1}{3}x_{\Phi}\\}$,
respectively. The SM-singlet RHNs $(N_{R}^{i})$ are transformed as
$\\{1,1,0,x_{\nu}=-x_{\Phi}\\}$ where $i$ represents three generations of the
fermions. The SM Higgs doublet and SM-singlet BSM scalar transform as
$\\{1,2,-\frac{1}{2},\frac{x_{H}}{2}\\}$ and $\\{1,1,0,2x_{\Phi}\\}$
respectively. These scalars are involved in the breaking of electroweak and
$U(1)_{X}$ symmetries. The Yukawa interactions among the scalars and fermions
in the general $U(1)_{X}$ extension of the SM can be written as

$\displaystyle{\cal L}^{\rm
Y}=-Y_{u}^{\alpha\beta}\overline{q_{L}^{\alpha}}Hu_{R}^{\beta}-Y_{d}^{\alpha\beta}\overline{q_{L}^{\alpha}}\tilde{H}d_{R}^{\beta}-Y_{e}^{\alpha\beta}\overline{\ell_{L}^{\alpha}}\tilde{H}e_{R}^{\beta}-Y_{\nu}^{\alpha\beta}\overline{\ell_{L}^{\alpha}}HN_{R}^{\beta}-\frac{1}{2}Y_{N}^{\alpha}\Phi\overline{(N_{R}^{\alpha})^{c}}N_{R}^{\alpha}+{\rm
h.c.}~{}$ (1)

where $H$ is the SM Higgs doublet and we write $\tilde{H}=i\tau^{2}H^{*}$
considering $\tau^{2}$ as the second Pauli matrix. From the interaction
Lagrangian in Eq. (1) we derive the following relations

$\displaystyle\frac{1}{2}x_{H}$ $\displaystyle=$ $\displaystyle-x_{q}+x_{u}\
=\ x_{q}-x_{d}=\ x_{\ell}-x_{e}=\ -x_{\ell}+x_{\nu}~{};-2x_{\Phi}=2x_{N}.$ (2)

Hence, applying anomaly cancellation conditions from Das:2016zue , and using
Eq. (2) we obtain the $U(1)_{X}$ charges in terms of $x_{H}$ and $x_{\Phi}$ so
that general $U(1)_{X}$ charges of the SM fermions can be expressed as the
linear combination of $U(1)_{Y}$ and B$-$L charges. We find an interesting
feature of this model fixing $x_{\Phi}=1$ without the loss of generality and
varying $x_{H}$. We find, with $x_{H}=-2$, the $U(1)_{X}$ charges of left
handed fermions $\ell_{L}^{i}$ and $q_{L}^{i}$ reduce to zero which further
converts the model into $U(1)_{\rm R}$ scenario. With $x_{\Phi}=1$ and
vanishing $x_{H}$ the $U(1)_{X}$ charge assignment reduces into B$-$L
scenario. Furthermore, taking $x_{H}=-1$, $-0.5$ and $1$, $U(1)_{X}$ charge of
$e_{R}^{i}$, $u_{R}^{i}$ and $d_{R}^{i}$ turns out to be zero. As a result,
the left and right handed fermions interact differently with the BSM neutral
gauge boson $Z^{\prime}$ in this scenario. The renormalizable scalar potential
of this model due to the presence of two scalar fields $H$ and $\Phi$ can be
written as

$\displaystyle V=\
m_{h}^{2}(H^{\dagger}H)+\lambda_{H}(H^{\dagger}H)^{2}+m_{\Phi}^{2}(\Phi^{\dagger}\Phi)+\lambda_{\Phi}(\Phi^{\dagger}\Phi)^{2}+\lambda^{\prime}(H^{\dagger}H)(\Phi^{\dagger}\Phi)\,,$
(3)

and for simplicity we assume that $\lambda^{\prime}$ to be very small. After
the breaking of $U(1)_{X}$ gauge and electroweak symmetries, $H$ and $\Phi$
develop vacuum expectation values (VEVs) as

$\displaystyle\langle H\rangle\ =\ \frac{1}{\sqrt{2}}\begin{pmatrix}v+h\\\
0\end{pmatrix},\,\,\,\langle\Phi\rangle\ =\ \frac{v_{\Phi}+\phi}{\sqrt{2}}\,.$
(4)

The electroweak scale is $v=246$ GeV at the potential minimum with $v_{\Phi}$
being a free parameter. The breaking of general $U(1)_{X}$ symmetry generates
mass of the neutral $Z^{\prime}$ gauge boson as

$M_{Z^{\prime}}=g_{X}\sqrt{4x_{\Phi}^{2}v_{\Phi}^{2}+x_{H}^{2}v^{2}}\,,$ (5)

which further reduces to $M_{Z^{\prime}}\simeq 2g_{X}x_{\Phi}v_{\Phi}$
considering $v_{\Phi}\gg v$ where $g_{X}$ is the general $U(1)_{X}$ gauge
coupling. From Eq. (1) we find that the RHNs interact with $\Phi$ which
generate the Majorana mass $M_{N}\ =\ Y_{N}/\sqrt{2}v_{\Phi}$ for the RHNs
after the breaking of general $U(1)_{X}$ symmetry and the Dirac mass
$m_{D}=Y_{\nu}/\sqrt{2}v$ involving $H$, $\ell_{L}$ and $N_{R}$ is generated
after the electroweak symmetry breaking. These two masses induce the tiny
neutrino mass through the seesaw mechanism followed by the flavor mixing.
Hence we obtain the neutrino mass matrix as

$m_{\nu}=\begin{pmatrix}0&m_{D}\\\ m_{D}^{T}&M_{N}\end{pmatrix}\,,$ (6)

and diagonalizing neutrino mass matrix, the light neutrino mass eigenvalue can
be obtained as $-m_{D}M_{N}^{-1}m_{D}^{T}$. Then the light neutrino flavor
eigenstate $\nu^{i}_{L}$ can be written as $\nu^{i}_{L}\simeq
U_{i\alpha}\nu_{\alpha}+V_{i\alpha}N_{\alpha}$ in terms of light
$(\nu_{\alpha})$ and heavy $(N_{\alpha})$ neutrino mass eigenstates where
$U_{i\alpha}\simeq(U_{\rm PMNS})_{i\alpha}$ where $U_{\rm PMNS}$ is the PMNS
matrix and $V_{i\alpha}\simeq(m_{D}M_{N}^{-1})_{i\alpha}$ being the mixing
between the light and heavy neutrinos through which RHNs interact with the SM
gauge bosons. We consider $m_{D}/M_{N}\ll 1$ and we ignore it when we
calculate $Z^{\prime}$ decay below. The interaction between $Z^{\prime}$ and
the left $(f_{L})$ and right $(f_{R})$ handed fermions under general
$U(1)_{X}$ scenario is given as

$\displaystyle\mathcal{L}^{f}=-g_{X}(\bar{f}_{L}\gamma^{\mu}q_{f_{L}}f_{L}+\bar{f}_{R}\gamma^{\mu}q_{f_{R}}f_{R})Z_{\mu}^{\prime}\,.$
(7)

where $q_{f_{L}}$ and $q_{f_{R}}$ are general $U(1)_{X}$ charges of the SM
left and right handed fermions which manifest the chiral scenario of this
model. We calculate the partial decay widths of $Z^{\prime}$ into different SM
fermions for a single generation as

$\displaystyle\Gamma(Z^{\prime}\to
f\bar{f})=N_{C}\frac{M_{Z^{\prime}}g_{X}^{2}}{24\pi}\left[\left(q_{f_{L}}^{2}+q_{f_{R}}^{2}\right)\left(1-\frac{{m_{f}^{2}}}{{M_{Z^{\prime}}^{2}}}\right)+6\,q_{f_{L}}q_{f_{R}}\frac{m_{f}^{2}}{M_{Z^{\prime}}^{2}}\right]\sqrt{1-\frac{{4m_{f}^{2}}}{{M_{Z^{\prime}}^{2}}}}\,,$
(8)

where $m_{f}$ is the SM fermion mass and $N_{C}=1(3)$ for the SM
leptons(quarks) being the color factor. The partial decay width of
$Z^{\prime}$ into a pair of single generation light neutrinos is given by

$\displaystyle\Gamma(Z^{\prime}\to\nu\nu)=\frac{M_{Z^{\prime}}g_{X}^{2}}{24\pi}q_{f_{L}}^{2}\,,$
(9)

where tiny neutrino mass has been neglected and
$q_{f_{L}}=-\frac{1}{2}x_{H}-x_{\Phi}$ represents general $U(1)_{X}$ charge of
$\ell_{L}$. This vanishes for $x_{H}=-2$ and $x_{\Phi}=1$ and as a result
there will be no tree level decay of $Z^{\prime}$ into a pair of neutrinos in
the $U(1)_{R}$ scenario. In general $U(1)_{X}$ extended SM scenarios,
$Z^{\prime}$ couples with RHNs as following

$\displaystyle\mathcal{L}_{N}=-\frac{1}{2}g_{X}q_{N_{R}}\overline{N}\gamma_{\mu}\gamma_{5}NZ_{\mu}^{\prime}\,.$
(10)

We calculate the corresponding partial decay width of $Z^{\prime}$ into a pair
of single generation of the heavy neutrino as

$\displaystyle\Gamma(Z^{\prime}\to
N_{R}^{\alpha}N_{R}^{\alpha})=\frac{M_{Z^{\prime}}g_{X}^{2}}{24\pi}\,q_{N_{R}}^{2}\,\left(1-4\frac{M_{N}^{2}}{M_{Z^{\prime}}^{2}}\right)^{\frac{3}{2}}\,,$
(11)

with $q_{N_{R}}(=x_{\Phi}=1)$ is the general $U(1)_{X}$ charge of the RHNs as
mentioned before. In the following we discuss few possibilities of having DM
candidates in the generic $U(1)_{X}$ model:

  * (i)

Majorana DM: We consider a scenario, where a Majorana fermions could be
introduced as a potential DM candidate. Among three RHNs in general $U(1)_{X}$
scenario, one generation of the RHNs could have odd $Z_{2}$ parity to ensure
its stability whereas the remaining generations of the RHNs could be
$Z_{2}$-even. We consider $N^{3}_{R}$ as the potential DM candidate while
neutrino mass and flavor mixing will be governed by $N^{1,2}_{R}$. Being
charged under $U(1)_{X}$, DM interacts with the SM sector through $Z^{\prime}$
assuming a tiny mixing in the scalar sector. We write the interaction
Lagrangian of the DM candidate as

$\displaystyle-\mathcal{L}_{\rm DM}^{\rm
M}=\Bigg{\\{}\frac{1}{2}Y_{N}^{3}\Phi\overline{({N_{R}^{3}})^{c}}N^{3}_{R}+h.c.\Bigg{\\}}+\frac{1}{2}g_{X}Q_{\chi}N^{3}\gamma^{\mu}\gamma_{5}N^{3}Z_{\mu}^{\prime}\,,$
(12)

where the first term corresponds to the DM mass after general $U(1)_{X}$
breaking $(M_{N}^{3}=\frac{v_{\Phi}}{\sqrt{2}}Y_{N}^{3}\equiv m_{\rm{DM}})$
and the second term corresponds to its interaction with $Z^{\prime}$.

  * (ii)

Scalar DM: We then take up another scenario, where a complex scalar $\Phi_{1}$
can be introduced as a potential DM candidate with the SM$\otimes U(1)_{X}$
charge assignment $\\{1,1,0,Q_{\chi}\\}$. $\Phi_{1}$ is $Z_{2}$-odd and the
remaining particles are even under $Z_{2}$ parity. $Q_{\chi}$ is the general
$U(1)_{X}$ charge of $\Phi_{1}$ because of which $\Phi_{1}$ interacts with
other fields of the model exchanging $Z^{\prime}$. However, $\Phi_{1}$ can
interact with scalar sector of the model through the following potential

$\displaystyle V\supset\lambda_{\rm
minx_{1}}(H^{\dagger}H)(\Phi_{1}^{\ast}\Phi_{1})+\lambda_{\rm
minx_{2}}(\Phi^{\dagger}\Phi)(\Phi_{1}^{\ast}\Phi_{1})\,.$ (13)

We assume $\lambda_{\rm mix_{1,2}}$ to be very small in the line of scalar
mixing considered in Eq. (3) taken to be small. The interaction between
$\Phi_{1}$ and $Z^{\prime}$ takes the following form

$\displaystyle\mathcal{L}_{\rm DM}^{\rm
s}=Q_{\chi}g_{X}Z^{\prime}_{\mu}\Big{\\{}\Phi_{1}^{\ast}(\partial^{\mu}\Phi_{1})-(\partial^{\mu}\Phi_{1}^{\ast})\Phi_{1}\Big{\\}}+m_{\Phi_{1}}^{2}\Phi_{1}^{\ast}\Phi_{1}\,,$
(14)

where $m_{\Phi_{1}}(=m_{\rm{DM}})$ is the mass of $\Phi_{1}$.

  * (iii)

Dirac DM: Lastly, we consider an alternative scenario where general $U(1)_{X}$
model can be extended with a weakly interacting Dirac fermion field
$\chi_{L,R}$ following the charge assignment $\\{1,1,0,Q_{\chi}\\}$. This
Dirac field could be considered as a potential DM candidate with an assigned
general $U(1)_{X}$ charge $Q_{\chi}$ to ensure its stability. Hence the
interaction Lagrangian can be written as

$\displaystyle\mathcal{L}_{\rm DM}^{\rm
D}=i\overline{\chi}\gamma^{\mu}(\partial_{\mu}+ig_{X}Q_{\chi}Z^{\prime}_{\mu})\chi+m_{\chi}\overline{\chi}\chi\,,$
(15)

with $\chi=\chi_{L}+\chi_{R}$ and $m_{\chi}(=m_{\rm{DM}})$ is the mass of
$\chi$. To ensure the stability of the DM we find that $Q_{\chi}\neq\\{\pm
3x_{\Phi},\pm x_{\Phi}\\}$ $=\\{\pm 3,\pm 1\\}$ assuming $x_{\Phi}=1$ for
$U(1)_{X}$. Except these charges other possibilities could be allowed.

We mention that each DM scenario has been assumed to be present one at a time
and we simply assume that kinetic mixing between $Z^{\prime}$ and SM neutral
gauge boson $Z$ to be very small. With all the interactions and particle
content at our disposal, we now move on to the discussions of quantifying DM-
electron scattering event.

Figure 1: DM-electron elastic scattering relevant for DM direct detection
(left) and SM fermion annihilation into DM applicable for freeze-in (right),
where DM can be a Dirac fermion, a Majorana fermion or a scalar boson.

##  III Dark matter electron scattering

The goal of the direct detection experiments is to observe the energy that is
deposited when a galactic DM particle scatters on the material used for the
detector. The direct detection experiments are in general of two types: DM
scattering on either the atomic nuclei or on the electrons. The main features
of these two types of scattering are determined by kinematics. DM-atomic
nuclei scattering is sensitive to DM mass scale of GeV or higher, because in
order to detect the nuclear recoils, the recoil energy must be above a certain
threshold and sub-GeV DM does not produce collisions with high enough
energies. However, as the DM-electron scattering is inelastic in nature, it
presents a viable framework for the detection of sub-GeV DM particles. This
suggests that larger fraction of the DM particle’s energy is transferred to
the electron in the collision and ionization signals can theoretically be
detected. The kinematics of the DM-electron process where the energy excites
the electron to either a higher energy bound state or an ionized state is far
more complex compare to DM-nuclear elastic scattering as the bound state
electrons does not have definite momentum and could in fact have arbitrarily
high momentum 111The scattering process may take place with any momentum
transfer $q=|{\bf q}|$, however when $q$ deviates too far above the inverse
Bohr radius, $a_{0}^{-1}\approx 3.7$ KeV, scattering rate is suppressed as it
is unlikely for the atomic electron to be found with such a high momentum.. In
this case, the energy transfer to the electron is $\Delta E_{e}=E_{\rm
rec}+E_{b}^{i}$, where $E_{\rm rec}$ is the electron recoil energy and
$-E_{b}^{i}$ is the binding energy of the electron where $i$ referring to an
initial state. Assuming the recoil of the atom is very small, it is possible
to calculate the minimum DM velocity required for electron recoil $E_{\rm
rec}$ Essig:2015cda ; Radick:2020qip ,

$\displaystyle v_{\text{min}}={q\over 2m_{\rm DM}}+{E_{\rm rec}+E_{b}^{i}\over
q}\,.$ (16)

Hence we see that minimum dark matter velocity required to excite an electron
with recoil energy $E_{\rm rec}$ also depends on the electron binding energy.
In order to compute the recoil event rate, one needs the velocity-averaged
differential cross-section, which reads Essig:2011nj ; Kopp:2009et ,

$\frac{d\langle\sigma_{\text{ion}}^{{i}}v\rangle}{d\ln E_{\rm
rec}}=\frac{\bar{\sigma}_{e}}{8\mu^{2}_{e\text{DM}}}\int
dq\,q|f_{\text{ion}}^{{i}}(k^{\prime},q)|^{2}|F_{\text{DM}}(q)|^{2}\eta(v_{\rm
vmin},t)\,,$ (17)

where $\mu_{e\text{DM}}=m_{\rm DM}/(m_{e}+m_{\rm DM})$ is the reduced mass for
DM-electron system, $\eta$ is the mean inverse speed and
$|f_{\text{ion}}^{n\ell}(k^{\prime},q)|^{2}$ is the ionization form factor. We
will discuss about the mean inverse speed and ionization form factor later in
this section. The DM-electron cross-section is conventionally normalized with
a reference cross-section $\bar{\sigma}_{e}$, which is defined as

$\displaystyle\bar{\sigma}_{e}=\frac{\mu_{e\text{DM}}^{2}}{16\pi\,m_{\text{DM}}^{2}\,m_{e}^{2}}\,|\overline{\mathcal{M}_{e\text{DM}}(q)}|^{2}\Big{|}_{q^{2}=\alpha^{2}\,m_{e}^{2}}\,,$
(18)

Here $\bar{\sigma}_{e}$ is the non-relativistic (NR) DM–electron elastic
scattering cross section, but with the 3-momentum transfer $q$ fixed to the
reference value $\alpha\,m_{e}$ [cf. left panel of Fig. 1]. The spin-averaged
squared amplitude for the matrix element encoding DM-electron scattering is
given by $|\overline{\mathcal{M}_{e\text{DM}}(q)}|^{2}$, and are reported in
Appendix. A. In reality we need to include the momentum transfer dependence
and the form factor $F_{\rm DM}(q)$ in Eq. (17) captures this $q$-dependence
of the matrix element:

$\displaystyle|F_{\rm
DM}(q)|^{2}=|\overline{\mathcal{M}_{e\text{DM}}(q)}|^{2}/|\overline{\mathcal{M}_{e\text{DM}}(\alpha
m_{e})}|^{2}.$ (19)

In the NR-limit, we find the following two most interesting limits for the
form factor $F_{\rm DM}(q)$,

$\displaystyle F_{\rm
DM}(q)=\frac{M_{Z^{\prime}}^{2}+\alpha^{2}\,m_{e}^{2}}{q^{2}+M_{Z^{\prime}}^{2}}\simeq\begin{cases}1&M_{Z^{\prime}}\gg\alpha\,m_{e}\\\\[10.0pt]
\alpha^{2}\,m_{e}^{2}/q^{2}&M_{Z^{\prime}}\ll\alpha\,m_{e}\,.\end{cases}$ (20)

Since $\alpha\,m_{e}\simeq 3.7\times 10^{-6}$, for the parameter space of our
interest, we will always have $F_{\rm DM}(q)=1$. Note that the DM-electron
scattering cross-section for fixed $U(1)_{X}$ charges is controlled by the
following set of independent parameters:
$\big{\\{}m_{\text{DM}},\,M_{Z^{\prime}},\,g_{X}\big{\\}}$.

Figure 2: DM-electron scattering cross-section for three benchmark points
$m_{\rm DM}=M_{Z^{\prime}}$ (solid black), $m_{\rm DM}=100M_{Z^{\prime}}$
(dotted black) and $m_{\rm DM}=0.01M_{Z^{\prime}}$ (dot dashed black line)
with $g_{X}=10^{-4}$ and $Q_{\chi}=100$. We show existing constraints on the
DM-electron scattering cross-section from Super-CDMS SuperCDMS:2020ymb , DAMIC
DAMIC:2019dcn , SENSEI SENSEI:2020dpa , PandaX-II PandaX-II:2021nsg ,
DarkSide-50 (DS-50) DarkSide:2022knj  and XENON1T-S2 (XE1T-S2) XENON:2019gfn .
Here we do not impose the requirement of right DM relic density.

In Fig. 2 we illustrate the DM-electron scattering cross-section for Majorana
(upper left panel), Scalar (upper right panel) and Dirac DM (bottom panel),
respectively for three benchmark points $m_{\rm DM}=M_{Z^{\prime}}$ (solid
black line), $m_{\rm DM}=100\,M_{Z^{\prime}}$ (dotted black line) and $m_{\rm
DM}=0.01\,M_{Z^{\prime}}$ (dot dashed black line) with $g_{X}=10^{-4}$. We
fixed the $U(1)_{X}$ charges as $x_{H}=-2$ and $x_{\Phi}=1$ and the charge for
the Dirac DM as $Q_{\chi}=100$. In the same plot, we also show the existing
and future constraints on the DM-electron scattering cross-section coming from
experiments such as Super-CDMS SuperCDMS:2020ymb , DAMIC DAMIC:2019dcn ,
SENSEI SENSEI:2020dpa , PandaX-II PandaX-II:2021nsg , DarkSide-50 (DS-50)
DarkSide:2022knj  and XENON1T-S2 (XE1T-S2) XENON:2019gfn .

###  III.1 Calculation of dark matter scattering event rate

We see from Eq. (17) that DM-electron differential scattering rate depends on
three different type of inputs: astrophysical input $\eta(v_{\rm min},t)$,
atomic physics factor $|f^{i}_{\rm ion}(k^{\prime},q)|^{2}$ and the particle
physics input $\bar{\sigma}_{e}|F_{\rm DM}(q)|^{2}$. All the relevant new
physics information are encoded in the particle physics input
$\bar{\sigma}_{e}$ which we already discussed in the above. Here we discuss in
detail about $\eta$ and $f^{i}_{\rm ion}$ which we need to know to compute the
expected event rate in a particular direct detection experiments. The mean
inverse speed $\eta$ is defined as,

$\eta(v_{\rm min},t)\equiv\int_{v_{\rm
min}}^{\infty}\frac{f_{\oplus}\left({\bf v},t\right)}{v}\,d^{3}v,$ (21)

where $f_{\oplus}\left({\bf v},t\right)$ is the DM velocity distribution in
the Erath frame which can be be found simply by applying a Galilean
transformation to velocity distribution in Galactic frame $f_{\infty}({\bf
v})$ 222Assuming the velocity distribution follows the Standard Halo Model, it
can be written as, $f_{\infty}({\bf v})={1\over
N_{\text{esc}}}\left({1\over\pi v_{0}^{2}}\right)^{3/2}e^{-{\bf
v}^{2}/v_{0}^{2}}\text{ for }|{\bf v}|<v_{\rm esc}$ where $N_{\rm esc}$ is a
normalization factor, and $v_{0}\approx 220\text{ km/s}$, whereas the escape
velocity is $v_{\rm esc}\approx 550\text{ km/s}$., so that,

$\displaystyle f_{\oplus}\left({\bf v},t\right)\approx f_{\infty}\left({\bf
v}_{\odot}+{\bf v}+{\bf V}_{\oplus}(t)\right),$ (22)

where ${\bf v}_{\odot}=(11,232,7)\,\,\text{km/s}$ is the velocity of the Sun
in Galactic coordinates and ${\bf V}_{\oplus}(t)$ is the time-dependent
velocity of the Earth in the Solar frame. The atomic physics factor
$|f_{\text{ion}}^{n\ell}(k^{\prime},q)|^{2}$ is the wave function suppression
factor to ionize an electron in the bound state (which is labelled as $i$) to
a final state with momentum $k^{\prime}$ when the momentum transfer is $q$. It
is relatively easy to determine the form factor for the atomic target compare
to the semiconductor target. If the final state of the electron is a plane
wave, so that $k^{\prime}=\sqrt{2m_{e}E_{\rm rec}}$, then for an atomic target
with spherically symmetric full shell with quanta ($n,\ell$), the ionization
factor is given by,

Figure 3: Differential event rate for DM-electron scattering process. In each
panel, the blue dashed, black dot dashed and red solid lines stands for the
target Xenon, Silicon and Germanium, respectively. We fixed the DM mass as
$m_{\rm DM}=10$ MeV and other relevant benchmarks for the parameters as
$g_{X}=10^{-3}$, $m_{\rm DM}=M_{Z^{\prime}}$ such that the cross-section are
$\bar{\sigma}_{e}=7.3\times 10^{-40}\,\text{cm}^{2}$, $2.19\times
10^{-39}\,\text{cm}^{2}$ and $1.46\times 10^{-39}\,\text{cm}^{2}$ for Dirac,
Majorana and Scalar DM, respectively. Note that Dirac DM cross-section is
normalized by the charge $Q_{\chi}^{2}$.

$|f_{\text{ion}}^{n\ell}(k^{\prime},q)|^{2}=\frac{(2\ell+1)k^{\prime
2}}{4\pi^{3}q}\int_{|k^{\prime}-q|}^{|k^{\prime}+q|}dk\,k|\chi_{n\ell}(k)|^{2}\,,$
(23)

where $\chi_{n\ell}$ is the radial part of the bound state wave function in
momentum space which is given as Kopp:2009et ,

$\displaystyle\chi_{n\ell}(p)=2\pi\int\\!dr\,r^{2}R_{n\ell}(r)\int\\!d(\cos\theta)\,P_{\ell}(\cos\theta)\,e^{ipr\cos\theta},$
(24)

where, ${\bf p}$ is a momentum space vector with modulus $p$ and
$P_{\ell}(\cos\theta)$ is a Legendre polynomial. One can approximate the
radial wave functions $R_{n\ell}(r)$ using a linear combination of orbitals
known as Slater type orbitals (STOs) Bunge:1993jsz ,

$\displaystyle R_{n\ell}(r)=\sum_{j}C_{n\ell j}\frac{(2Z_{\ell j})^{n_{\ell
j}+1/2}}{a_{0}^{3/2}\sqrt{(2n_{\ell j})!}}(r/a_{0})^{n_{\ell
j}-1}\exp(-Z_{\ell j}r/a_{0})\,.$ (25)

Here, $a_{0}$ is the Bohr radius, and the parameters $C_{n\ell j}$, $n_{\ell
j}$, and $Z_{\ell j}$ are taken from Bunge:1993jsz . Finally using Eq. (25),
we can evaluate Eq. (24) analytically as follows,

$\displaystyle\chi_{n\ell}(p)=\sum_{j}C_{n\ell j}\,2^{-\ell+n_{\ell
j}}\bigg{(}\frac{2\pi a_{0}}{Z_{\ell
j}}\bigg{)}^{3/2}\bigg{(}\frac{ipa_{0}}{Z_{\ell
j}}\bigg{)}^{\ell}\frac{(1+n_{\ell j}+\ell)!}{\sqrt{(2n_{\ell
j})!}}\,\times{}_{2}F_{1}\bigg{[}\tfrac{1}{2}(2+\ell+n_{\ell
j}),\tfrac{1}{2}(3+\ell+n_{\ell
j}),\tfrac{3}{2}+\ell,-\bigg{(}\frac{pa_{0}}{Z_{\ell
j}}\bigg{)}^{2}\bigg{]}\,,$ (26)

with ${}_{2}F_{1}(a,b,c,x)$ being a hypergeometric function. With these above
information, the differential event rate can now be obtained by summing over
the all possible initial electron states as

$\frac{dR}{d\ln{E_{\rm rec}}}=N_{T}\frac{\rho_{\rm DM}}{m_{\rm
DM}}F(k^{\prime})\sum_{i}\frac{d\langle\sigma_{ion}^{{i}}v\rangle}{d\ln E_{\rm
rec}}\,,$ (27)

where $N_{T}$ is the number of target nuclei and $\rho_{\rm DM}\approx
0.4\text{ GeV}/\text{cm}^{3}$ is the local DM density. $F(k^{\prime})$ is the
Fermi factor which accounts for how the atom itself distorts the scattered
electron’s wave function. In the NR limit this factor is given by

$F(k^{\prime})=\frac{2\pi\nu}{1-e^{-2\pi\nu}}\,\text{ with
}\nu=Z_{\text{eff}}\,(\alpha m_{e}/k^{\prime}),$ (28)

where $\alpha$ is the fine-structure constant and $Z_{\rm eff}$ is the
effective charge that is felt by the scattered electron. Following the Ref.
Essig:2011nj , we set $Z_{\rm eff}=1$ as this is expected to be a good
approximation for outer-shell electrons. In Fig. (3), the blue dotted line
shows the differential event rate, plotted against the electron recoil energy
$E_{\rm rec}$ assuming Xenon atomic target. The upper left, upper right and
bottom panel stands for Majorana, scalar and Dirac DM, respectively. For this
plot we fixed the DM mass as $m_{\rm DM}=10$ MeV and other relevant parameters
as $g_{X}=10^{-3}$, $m_{\rm DM}=M_{Z^{\prime}}$ such that the cross-section
are $\bar{\sigma}_{e}=7.3\times 10^{-40}\,\text{cm}^{2}$, $2.19\times
10^{-39}\,\text{cm}^{2}$ and $1.46\times 10^{-39}\,\text{cm}^{2}$ for Dirac,
Majorana and scalar DM, respectively. The rate was determined using only the
three outermost orbitals (5p, 5s, and 4d), which have binding energies of
approximately 12, 26, and 76 eV, respectively.

Now let us consider the scenario where DM excites electrons in a semiconductor
target above the band gap. In contrast to noble gas targets with binding
energies of $\mathcal{O}(10)$ eV, semiconductor materials allow for electron
ionization energies of $\mathcal{O}(1)$ eV, making them an excellent target
for studying DM-electron scattering 333For example, in Germanium, electron-
hole pairs are ionized to the conduction band by interactions that deposits
energy above the band gap of about 0.67 eV.. Since the electrons in a
semiconductor target are characterized by Bloch wave functions in a periodic
lattice, calculating the ionization form factor is particularly challenging
compare to atomic target. There are many packages to calculate the ionization
form factors such as EXCEED-DM Trickle:2022fwt ; Griffin:2021znd , QEdark
Essig:2015cda  and Quantum-Espresso Giannozzi:2009gev . There are also
analytic and semi-analytic approach to calculate the ionization form factors
Graham:2012su ; Lee:2015qva . Here, we follow Ref. Lee:2015qva , where form
factor was derived using the Roothaan-Hartree-Fock (RHF) wave functions for
the electrons. The delocalized eigenfunctions of the Hamiltonian for an
electron in any periodic potential can be written in terms of localized Bloch
wave functions, $\Psi_{\bf k}({\bf r})$, which may be expressed using Wannier
functions book ,

$\displaystyle\Psi_{\bf k}({\bf r})=\sum_{N}e^{i\,{\bf k}.{\bf
R}_{N}}\phi({\bf r}-{\bf R}_{N}),$ (29)

where $\phi({\bf r})$ is a Wannier function localized at the site ${\bf
R}_{N}$ and ${\bf k}$ are the wave vectors in the first Brillouin zone (BZ).
The wave function of the scattered electron can be roughly represented as a
plane wave for sufficiently high momentum transfers. Hence we can obtain the
scattering cross section by considering the transition of an electron from a
localized initial-state atomic wave function–with a ${\bf k}$-dependent
binding energy–to a final-state wave function with plane-wave solution. The
total differential event rate is then calculated by integrating over all
binding energies and weighting it according to the density of states,

$\frac{dR}{d\ln{E_{\rm rec}}}\approx N_{T}\frac{\rho_{\rm DM}}{m_{\rm
DM}}F(k^{\prime})\int
dE_{b}\,\rho(E_{b})\frac{d\langle\sigma_{\text{ion}}v\rangle}{d\ln E_{\rm
rec}}\,,$ (30)

where $\rho(E_{b})$ is the density of states determined experimentally ssd ;
Graham:2012su , which takes into account that a given set of $E_{b}$ values
corresponds to a different number of $\Psi_{\bf k}$ states. Hence, the density
of states serves as an efficiency factor for scattering at various binding
energies. In Fig. (3), the red-line (black dot-dashed line) shows the
differential event rate, plotted against the electron recoil energy $E_{\rm
rec}$, for DM mass $m_{\rm DM}=10$ MeV assuming the target material is
Germanium (Silicon). In Fig. (3), when we compare the predicted differential
scattering rate for the DM scattering off Germanium/Silicon with that for
Xenon, we see that the rate is larger in Germanium/Silicon for all electron
recoil energies; this is a direct effect of the lower binding energy in the
semiconductor.

##  IV Freeze-in production of dark matter

In this section we obtain constraints on the coupling and mass from the
requirement of producing right DM abundance. Here we are typically interested
in the DM-SM couplings that are feebler than weak-interaction strength,
therefore we prefer freeze-in production of DM over freeze-out. In the
following sections we will consider DM production takes place (i) during
radiation domination (RD) and (ii) in a modified cosmological background prior
to the onset of big bang nucleosynthesis (BBN). In the latter case, as we will
see, the DM-SM coupling $g_{X}$ can be amplified by at least
$\sim\mathcal{O}(10^{2})$ compared to the RD scenario, depending on the choice
of masses and couplings.

###  IV.1 During radiation domination

We first consider the scenario with standard cosmological history, where the
DM genesis takes place in the background of a radiation dominated Universe. We
note, the relevant DM production channels in this model are: (i) on-shell
1-to-2 decay of $Z^{\prime}$ into a pair of DM final state and (ii)
$Z^{\prime}$-mediated $s$-channel 2-to-2 scattering of the bath particles into
a pair of DM. Above the temperature of electroweak (EW) phase transition
$T_{\rm EW}\simeq 160$ GeV, we consider the SM states are absolutely
massless444Here we do not include temperature corrected mass terms for the
particles in the bath., while they all become massive once the EW symmetry is
broken. We therefore consider the total DM yield at present as a sum of the
yield before and after the EW symmetry breaking555This treatment has been
adopted, for example, in Refs. Duch:2017khv ; Barman:2020ifq ;
Bhattacharya:2021rwh .. Now, the Boltzmann equation (BEQ) governing the DM
number density can be written in terms of the DM yield defined as a ratio of
the DM number density to the entropy density in the visible sector, i.e.,
$Y_{\rm DM}=n_{\rm DM}/s$. The BEQ can then be expressed in terms of the
reaction densities as

$x\,H\,s\,\frac{dY_{\rm DM}}{dx}=\gamma_{\text{22}}+\gamma_{\text{12}}\,,$
(31)

where $x\equiv m_{\text{DM}}/T$ is a dimensionless quantity. The complete
expressions for the reaction densities $\gamma$’s due to 2-to-2 scattering
$(\gamma_{22})$ and 1-to-2 decay $(\gamma_{12})$ can be found in Appendix C.
In the subsequent analysis we shall consider $Z^{\prime}$ to be part of
thermal bath such that its number density follows standard Maxwell-Boltzmann
distribution. The Hubble parameter $H$ in the radiation dominated epoch and
the entropy per comoving volume $s$ are given by

$\displaystyle H\equiv H_{\rm
rad}(T)=\frac{\pi}{3}\,\sqrt{\frac{g_{\star}(T)}{10}}\,\frac{T^{2}}{M_{P}}\,,$
$\displaystyle s(T)=\frac{2\,\pi^{2}}{45}\,g_{\star s}(T)\,T^{3}\,,$ (32)

where $T$ is the temperature of the thermal bath and $g_{\star}(g_{\star s})$
are the relativistic degrees of freedom (DoF) associated with the energy
(entropy) density. Now, it is important to ensure that the DM does not
thermalize with the SM bath in the early Universe, at least till
$T=M_{Z^{\prime}}$. For decay, this can be checked by comparing
$\langle\Gamma_{Z^{\prime}\to\text{DM},\text{DM}}\rangle$ with the Hubble
parameter, where

$\displaystyle\langle\Gamma_{Z^{\prime}\to\text{DM},\text{DM}}\rangle=\Gamma_{Z^{\prime}\to\text{DM},\text{DM}}\times\frac{K_{1}\left(M_{Z^{\prime}}/T\right)}{K_{2}\left(M_{Z^{\prime}}/T\right)}\,,$
(33)

is the thermally averaged $Z^{\prime}$ decay into a pair of DM. Here $K_{1,2}$
are the modified Bessel functions and
$\Gamma_{Z^{\prime}\to\text{DM},\text{DM}}$ is the partial decay width of
$Z^{\prime}$ into a pair of DM final states. On comparison with the Hubble
rate we find

$\displaystyle\frac{\langle\Gamma_{Z^{\prime}\to\text{DM},\text{DM}}\rangle}{H}\Bigg{|}_{T=M_{Z^{\prime}}}\simeq\left(\frac{g_{X}}{10^{-8}}\right)^{2}\,\left(\frac{100\,\text{GeV}}{M_{Z^{\prime}}\,\sqrt{g_{\star}(T)}}\right)\times\begin{cases}7.1\times
10^{-2}&\quad\text{for Dirac DM}\\\\[10.0pt] 3.5\times 10^{-2}&\quad\text{for
Majorana DM}\\\\[10.0pt] 1.7\times 10^{-2}&\quad\text{for Scalar
DM}\,,\end{cases}$ (34)

for $m_{\text{DM}}=M_{Z^{\prime}}/2$. This shows, for the mass window
$M_{Z^{\prime}}$ of our interest, the out of equilibrium condition for DM
production via decay can be ensured for $g_{X}\lesssim 10^{-8}$.

Next, we focus on DM freeze-in through scattering and consider
$M_{Z^{\prime}}<2\,m_{\text{DM}}$ in order to kinematically forbid decay of
$Z^{\prime}$ into DM final states. Here we consider DM production from the SM
bath, via $s$-channel $Z^{\prime}$-mediation [cf. right panel of Fig. 1]. All
relevant production cross-sections are reported in Appendix. B. The condition
for non-thermal production in case of scattering can be obtained by comparing
the 2-to-2 scattering rate $\Gamma_{22}=\gamma_{22}/n_{\rm eq}$ with the
Hubble rate, where $n_{\rm
eq}=T/(2\,\pi^{2})\,m_{\text{DM}}^{2}\,K_{2}\left(m_{\text{DM}}/T\right)$ is
the equilibrium number density of the DM. In order to find an approximate
analytical estimation, we consider $s\gg
m_{\text{DM}}^{2}\,,M_{Z^{\prime}}^{2}$. Using the expressions for scattering
cross-sections in Appendix. B, we obtain

$\displaystyle\frac{\Gamma_{22}}{H}\Bigg{|}_{T=M_{Z^{\prime}}}\simeq\left(\frac{g_{X}}{10^{-4}}\right)^{4}\,\left(\frac{M_{Z^{\prime}}}{100\,\text{GeV}}\right)^{2}\,\frac{1}{\sqrt{g_{\star}(T)}}\times\begin{cases}1.1\times
10^{-5}\,Q_{\chi}^{2}&\quad\text{for Dirac DM}\\\\[10.0pt] 5.9\times
10^{-6}&\quad\text{for Majorana DM}\\\\[10.0pt] 2.9\times
10^{-6}&\quad\text{for Scalar DM}\,,\end{cases}$ (35)

with $m_{\text{DM}}=2\,M_{Z^{\prime}}$ and $x_{H}=-2,\,x_{\Phi}=1$. Note that,
this is rather a conservative estimation, since we consider equilibrium number
density of the DM. Evidently, as $\sigma(s)\propto g_{X}^{4}$, it is possible
to obtain the interaction rate below the Hubble rate for much larger $g_{X}$
compared to the decay case where the interaction rate had a $g_{X}^{2}$
dependence. As a conservative limit, we will always consider $g_{X}\lesssim
10^{-3}$, for freeze-in during RD. To fit the whole observed DM relic density,
it is required that

$Y_{0}\,m_{\text{DM}}=\Omega
h^{2}\,\frac{1}{s_{0}}\,\frac{\rho_{c}}{h^{2}}\simeq 4.3\times
10^{-10}~{}\text{GeV},$ (36)

where $Y_{0}\equiv Y_{\rm DM}(T_{0})$ is the DM yield at the present epoch,
that can be obtained by solving Eq. (31). Here $\rho_{c}\simeq 1.05\times
10^{-5}\,h^{2}$ GeV/cm3 is the critical energy density, $s_{0}\simeq
2.69\times 10^{3}$ cm-3 the present entropy density ParticleDataGroup:2022pth
, and $\Omega h^{2}\simeq 0.12$ the observed abundance of DM relics
Planck:2018vyg . One can find the approximate analytical expression for the
final yield as

$\displaystyle
Y_{0}\simeq\left(\frac{g_{X}}{10^{-6}}\right)^{4}\,\left(\frac{100\,\text{GeV}}{M_{Z^{\prime}}}\right)\,\frac{1}{g_{\star
s}(T)\sqrt{g_{\star}(T)}}\times\begin{cases}1.1\times
10^{-11}\,Q_{\chi}^{2}&\quad\text{for Dirac DM}\\\\[10.0pt] 5.8\times
10^{-12}&\quad\text{for Majorana DM}\\\\[10.0pt] 2.9\times
10^{-12}&\quad\text{for Scalar DM}\,,\end{cases}$ (37)

for $m_{\text{DM}}=2\,M_{Z^{\prime}}$. This shows, a heavier $M_{Z^{\prime}}$
requires larger $g_{X}$ in order to produce the right DM yield such that the
observed DM relic abundance is satisfied. We emphasize that these are only
approximate estimation for the DM yield, we however solve the full BEQ
numerically to obtain the viable parameter space, taking into account the
temperature-dependence of the relativistic degrees of freedom.

###  IV.2 In modified cosmology

Having discussed the freeze-in mechanism during radiation domination, we now
consider a scenario, where we assume the Universe prior to BBN has two
different species: radiation with energy density $\rho_{R}$, and some other
species $\varphi$ with energy density $\rho_{\varphi}$. Our discussion closely
follows Ref. DEramo:2017gpl ; DEramo:2017ecx , where this prescription has
been elaborated666Phenomenological consequences of such faster-than-usual
expansion scenario has been studied, for example, in Refs. Chen:2019etb ;
Mahanta:2019sfo ; Konar:2020vuu ; Arcadi:2021doo ; Barman:2021ifu ;
Mahanta:2022gsi ; Barman:2022njh .. In presence of $\varphi$, total energy
density of the Universe is $\rho=\rho_{R}+\rho_{\varphi}$. In case of a rapid
expansion of the Universe the energy density of $\varphi$ field is expected to
be redshifted faster than the radiation. Accordingly, one can assume
$\rho_{\varphi}\propto a(t)^{-(4+n)}$, where $a$ is the scale factor. Here
$n>0$ implies $\varphi$ energy density dominates over radiation during early
enough times. A general form of $\rho_{\varphi}$ can then be obtained using
the entropy conservation $g_{\star}\left(T\right)^{1/3}aT=\text{constant}$ in
a comoving frame as

$\rho_{\varphi}(T)=\rho_{\varphi}(T_{R})\left(\frac{g_{*s}(T)}{g_{*s}(T_{R})}\right)^{(4+n)/3}\left(\frac{T}{T_{R}}\right)^{(4+n)}\,,$
(38)

where the temperature $T_{R}$ is an unknown variable and can be considered as
the point of equality where $\rho_{\varphi}(T_{R})=\rho_{\text{rad}}(T_{R})$
is achieved. Using this, the total energy density at any temperature $T$ can
be expressed as

$\displaystyle\rho(T)$
$\displaystyle=\rho_{R}(T)+\rho_{\varphi}(T)=\rho_{R}(T)\left[1+\frac{g_{*}(T_{R})}{g_{*}(T)}\left(\frac{g_{*s}(T)}{g_{*s}(T_{R})}\right)^{(4+n)/3}\left(\frac{T}{T_{R}}\right)^{n}\right]\,.$
(39)

The Hubble rate then reads

$\displaystyle H(T)=H_{\rm
rad}\left(T\right)\Biggl{[}1+\Biggl{(}\frac{T}{T_{R}}\Biggr{)}^{n}\Biggr{]}^{1/2}\simeq\begin{cases}H_{\rm
rad}(T)\,,&T\ll T_{R}\\\\[10.0pt] H_{\rm
rad}(T)\,\left(\frac{T}{T_{R}}\right)^{n/2}\,,&T\gg T_{R}\,,\end{cases}$ (40)

where $H_{\rm rad}\left(T\right)$ is the Hubble parameter in the standard
radiation dominated Universe. It is important to note from Eq. (40) that the
expansion rate is larger than what it is supposed to be in the standard
cosmological background for $T>T_{R}$ and $n>0$. The temperature
$T_{R}\gtrsim\left(15.4\right)^{1/n}~{}\text{MeV}$ DEramo:2017ecx ;
DEramo:2017gpl , such that the generation of light nuclei during BBN is not
hampered. For $n=3$, the lower bound on $T_{R}$ turns out to be 2.5 MeV. In
the rest of the analysis we will fix $T_{R}=5$ MeV (which is just above the
lower bound on BBN temperature $T_{\rm BBN}\simeq 4$ MeV) and consider $n\leq
3$. Note that, this bound ceases to exist for $n=0$, as in that case we get
back standard cosmology, without any new extra species driving the expansion
of the Universe. As it has been explained in DEramo:2017gpl , within the
modified cosmological setup, the DM is always under-produced with respect to
the case of a standard history (RD). This can be physically understood from
the fact that as DM production takes place during the epoch of fast expansion,
hence the DM number density easily gets diluted, resulting in under abundance.
Thus, the key point here is: in a faster-than-standard expanding universe,
freeze-in production is suppressed, implying, to produce enough DM to match
observations, larger coupling is required777The same inference also applies to
freeze-in production of DM during an early matter dominated (EMD) era
Hardy:2018bph ; Bernal:2018kcw ; Cosme:2020mck .. This in turn improves the
observational aspect for freeze-in.

Figure 4: DM yield $(Y_{\rm DM})$ as a function of $m_{\text{DM}}/T$, for
standard cosmology (red) and modified scenario (in green, blue and black) for
freeze-in production via $Z^{\prime}$-mediated $\rm{SM~{}SM}\to\rm{DM~{}DM}$
scattering process. We consider $Q_{\chi}=100$ and $T_{R}=5$ MeV, with
$M_{Z^{\prime}}=m_{\text{DM}}/2=100$ GeV. The gray horizontal band corresponds
to observed DM abundance and the vertical dashed black line, $x_{\rm fi}$,
denotes the epoch of freeze-in.

This is shown in Fig. 4, where we illustrate the evolution of DM yield as a
function of $m_{\text{DM}}/T$, for different DM candidates and for a fixed
$T_{R}=5$ MeV. In each case the red curve corresponds to DM production in
standard radiation dominated Universe. As we see, increasing $n$ results in
suppression of DM production. Thus, in order to match the right relic
abundance, one needs to increase $g_{X}$. The behaviour of the DM yield also
clearly shows the IR nature of freeze-in, where freeze-in typically happens at
low temperature, which in this case is given by $x_{\rm fi}\equiv
m_{\text{DM}}/T_{\rm fi}\simeq 4$. It is necessary to mention that although
the evolution of yield seems to be identical in all three cases, that is
because we have fixed the DM mass, which in turn fixes the asymptotic DM yield
$Y_{0}\simeq 2.2\times 10^{-12}\,\left(200\,\text{GeV}/m_{\text{DM}}\right)$,
however the coupling required to produce the observed abundance for different
DM spins is largely different, as we will see in subsequent plots. Due to
faster than usual expansion of the Universe, it is easier to keep the DM out
of equilibrium. Consequently, with larger $n$, it is possible to achieve
larger couplings compared to the standard RD scenario, satisfying freeze-in
requirements.

Although we remain agnostic about the possible cosmological models that can
incorporate such modification to the standard lore, we point towards a few
well-motivated scenarios, where such modification to the Hubble expansion rate
can be realized. For instance, theories with $n=2$ are quintessence fluids
motivated by the accelerated expansion of the Universe Caldwell:1997ii ;
Sahni:1999gb . One possible realization of scalar potential giving rise to
this behaviour is $V(\varphi)\sim\exp\left(-\lambda\,\varphi\right)$
Ratra:1987rm ; Copeland:1997et . For theories with $n>2$ one has to consider
scenarios faster than quintessence. Example of such theories can be found, for
example, in Buchbinder:2007ad , where one assumes the presence of a pre-big
bang “ekpyrotic” phase. Other than these versions of the non-standard
cosmology, there is also the scope for modification of Hubble rate due to
modified gravity theories Dunsby:1997fyr ; Catena:2009tm ; Dent:2009bv ;
Leon:2013qh ; Baules:2019zwk .

##  V Results and discussions

Figure 5: Limits on $g_{X}-M_{Z^{\prime}}$ plane from DM-electron scattering
process in $t$-channel mediated by $Z^{\prime}$ for $x_{H}=-2$ considering
$m_{\rm DM}=10$(100) MeV in left (right) panel for Majorana (upper), Scalar
(middle) and Dirac (lower panel) DM candidates. Shaded regions are ruled out
by existing experimental data whereas prospective bounds could appear from
MUonE, DUNE, FASER2 and ILC-BD experiments. We add constraints obtained from
Freeze-in mechanism considering $m_{\rm DM}=M_{Z^{\prime}}/2$ including
modified cosmology.

Figure 6: Limits on $g_{X}-M_{Z^{\prime}}$ plane from DM-electron scattering
process in $t$-channel mediated by $Z^{\prime}$ for $x_{H}=0$ $(2)$
considering $m_{\rm DM}=100$ MeV in left (right) panel for Majorana (upper),
Scalar (middle) and Dirac (lower panel) DM candidates. Shaded regions are
ruled out by existing experimental data whereas prospective bounds could
appear from MUonE, DUNE, FASER2 and ILC-BD experiments. We add constraints
obtained from Freeze-in mechanism considering $m_{\rm DM}=M_{Z^{\prime}}/2$
including modified cosmology.

We summarize the results in Figs. 5-6, where we show the bounds on
$g_{X}-M_{Z^{\prime}}$ plane, for Majorana, scalar and Dirac DM from electron-
DM scattering and freeze-in scenarios, respectively. In Fig. 5 we show the
constraints for $x_{H}=-2$, the $U(1)_{R}$ case, where in case of electron-DM
scattering we chose $m_{\rm DM}=10$(100) MeV shown in the left (right) panel
for Majorana (upper), Scalar (middle) and Dirac DM (lower panel). We chose
$x_{H}=-2$ because for this charge left handed fermions do not directly
interact with $Z^{\prime}$ whereas right handed fermions do. As a result left
handed neutrinos do not interact with $Z^{\prime}$ forbidding the invisible
decay of $Z^{\prime}$. On the other hand we do the same for the $x_{H}=0$, the
B$-$L scenario where left and right handed fermions interact with $Z^{\prime}$
in the same way. In addition to that we consider another $U(1)_{X}$ charge
where $x_{H}=2$, where all the left and right handed fermions interact
differently with $Z^{\prime}$. The corresponding bounds are shown in Fig. 6
for $x_{H}=0$ (2) in the left (right) panel for the Majorana(upper),
Scalar(middle) and Dirac DM(lower panel) considering $m_{\rm DM}=100$ MeV to
study electron-DM scattering. In these analyses we considered the $U(1)_{X}$
charge of the Dirac DM to be $Q_{\chi}=100$. To estimate the constraints on
$g_{X}-M_{Z^{\prime}}$ plane for different $x_{H}$ from electron-DM scattering
and freeze-in scenarios we apply the perturbativity constraints
$g_{X}(x_{H,\Phi};Q_{\chi})<\sqrt{4\pi}$. We have already mentioned that in
this analysis we chose $x_{\Phi}=1$. Finally we compare our results with the
constraints obtained from different scattering and beam-dump experiments.

For the DM-electron scattering, we consider existing bounds from Super-CDMS
SuperCDMS:2020ymb , DAMIC DAMIC:2019dcn , SENSEI SENSEI:2020dpa , PandaX-II
PandaX-II:2021nsg , DarkSide-50 (DS-50) DarkSide:2022knj  and XENON1T-S2
XENON:2019gfn . For a given DM-type, we choose two benchmark values of DM
mass: 10 MeV and 100 MeV. This is motivated by the fact that scattering
experiments such as CDMS, DAMIC, and SENSEI impose the most stringent bounds
for 10 MeV DM, while PandaX-II, DarkSide-50, and XENON1T-S2 offer greater
sensitivity to DM with a mass of 100 MeV. This is illustrated in Fig. 2.
Therefore, in the $g_{X}-M_{Z^{\prime}}$ plane, these benchmark DM masses
shall provide the strongest constraints. For example, in the left (right)
panel of Fig. 5, we consider $m_{\text{DM}}=10(100)$ MeV, and bounds from
Super-CDMS, DAMIC, SENSEI, PandaX-II, DarkSide-50 (DS-50) and XENON1T-S2 are
shown via green, red, purple, orange, blue and cyan dashed diagonal straight
lines respectively. In the right panel we show the same, but for a DM of mass
100 MeV. We find, irrespective of the DM-spin and mass, XENON1T-S2 provides
most stringent bound for $x_{H}=-2$. In case of Majorana DM, XENON1T-S2
provide strongest limit on the general $U(1)_{X}$ coupling as $2.4\times
10^{-3}\lesssim g_{X}\lesssim 4.3\times 10^{-3}$ with 0.012 GeV $\lesssim
M_{Z^{\prime}}\lesssim$ 0.0203 GeV for $m_{\text{DM}}=10$ MeV and $1.05\times
10^{-3}\lesssim g_{X}\lesssim 1.8\times 10^{-3}$ with 0.02 GeV $\lesssim
M_{Z^{\prime}}\lesssim$ 0.033 GeV for $m_{\text{DM}}=100$ MeV. These limits
are approximately same for the scalar case as well. Bounds below and above
this range of $M_{Z^{\prime}}$ are ruled out by the existing bounds from NA64
NA64:2019auh  and the constraints obtained from the visible decay of
$Z^{\prime}$ in BaBaR BaBar:2014zli ; BaBar:2017tiz  experiment. On the other
hand for Dirac DM because of the larger $U(1)_{X}$ charge, we see that the
limits from different experiments have different end-points. Also, in this
case, corresponding bounds on $g_{X}$ are stronger compared to the other two
cases. We find that anything above $2.7\times 10^{-4}\lesssim g_{X}\lesssim
1.3\times 10^{-3}$ and $1.25\times 10^{-4}\lesssim g_{X}\lesssim 10^{-3}$ for
$m_{\text{DM}}=10$ MeV and 100 MeV within the range 0.01 GeV $\lesssim
M_{Z^{\prime}}\lesssim$ 0.045 GeV and 0.018 GeV $\lesssim
M_{Z^{\prime}}\lesssim$ 0.140 GeV, respectively are ruled out by XENON1T-S2
data. Bounds below and above this range of $M_{Z^{\prime}}$ are ruled out by
E141 Riordan:1987aw  and visible decay of $Z^{\prime}$ in BaBaR BaBar:2014zli
; BaBar:2017tiz , respectively. In case of Dirac-DM scenario we find that the
strongest bound obtained from XENON1T-S2 is slightly stronger than the limits
obtained by the electron $g-2$ limits Asai:2023mzl . We compare our results
for the electron-DM scattering with the bounds estimated form the dark photon
searches in LHCb Foguel:2022unm ; Baruch:2022esd  and CMS experiments
CMS:2023slr  from Asai:2023mzl . These experiments strongly constrain 0.2 GeV
$\lesssim M_{Z^{\prime}}\lesssim$ 70 GeV where limits on the $U(1)_{X}$
coupling vary between $10^{-5}\lesssim g_{X}\lesssim 10^{-4}$. We also find
that constraints from electron $(g-2)$ Asai:2023mzl  can be slightly weaker
than (comparable with) the limits estimated from the DM-electron scattering
process in DS-50 (XENON1T-S2) experiment in case of Dirac-DM scenario for
$m_{\rm DM}=100$ MeV whereas estimated bounds from other experiments are
weaker than the limits from electron $(g-2)$. In case of Majorana and Scalar-
DM scenarios, limits from DM-electron scattering is weaker than the bounds
from electron $(g-2)$ for $m_{\rm DM}=10$ MeV and 100 MeV. Same results can be
obtained for the Dirac-DM case for $m_{\rm DM}=10$ MeV.

If we look at the scenarios for $x_{H}=0$ and $2$ in Fig. 6, we find that
limits obtained from invisible decay of $Z^{\prime}$ involving light neutrinos
in BaBaR BaBar:2014zli ; BaBar:2017tiz  rule out the limits obtained from
electon-DM scattering. In this case we show the bounds estimated considering
neutrino-electron scattering for TEXONO TEXONO:2009knm ; TEXONO:2002pra ;
TEXONO:2006xds , BOREXINO Borexino:2000uvj ; Bellini:2011rx ; Cleveland:1998nv
; Lande:2003ex ; Borexino:2007kvk ; Borexino:2008gab  from Asai:2023mzl  due
to complementarity and notice that these existing bounds rules out the limits
obtained from the electron-DM scattering. This happens because for these
charges light neutrinos directly interact with the $Z^{\prime}$ which are not
possible in case of $x_{H}=-2$ as the $U(1)_{X}$ charge of the SM lepton
doublet vanishes for $x_{\Phi}=1$. In the same line following Asai:2023mzl  we
show the limits on $g_{X}-M_{Z^{\prime}}$ plane for $x_{H}=0$ and 2 studying
the muon neutrino, anti-neutrino scattering with electron from CHARM-II CHARM-
II:1991ydz ; CHARM-II:1989nic ; CHARM-II:1993phx ; CHARM-II:1994dzw
experiment. In this context we also compare the bounds obtained from neutrino-
nucleon scattering from COHERENT COHERENT:2018imc ; COHERENT:2018imc ;
Cadeddu:2017etk ; Cadeddu:2020nbr ; COHERENT:2020iec ; COHERENT:2020ybo
experiment followed by considering the limits obtained from the neutrino
magnetic moment study in GEMMA experiment Beda:2009kx ; Lindner:2018kjo  from
Asai:2023mzl . To show complementarity we present the bounds obtained from the
dark photon searches at LHCb and CMS experiments from Asai:2023mzl  for
$x_{H}=0$ and 2. Limits obtained from DM-electron scattering for these charges
are comparable with the electron $(g-2)$ limits for the Dirac-DM scenario for
$m_{\rm DM}=10$ MeV and 100 MeV, however, it is stronger than the limits
obtained from muon $(g-2)$ Asai:2023mzl . Nature of the bounds for Majorana
and scalar-DM scenarios are same as $x_{H}=-2$. However, it has been already
pointed out that strong constraints from TEXONO and BOREXINO are also stronger
than the $(g-2)$bounds. For $M_{Z^{\prime}}\lesssim 0.04$ GeV beam dump
experiments provide stronger limits.

Finally we show the bounds obtained from the dilepton and dijet searches in
the LEP experiment ALEPH:2006bhb ; ALEPH:1997gvm ;
LEPWorkingGroupforHiggsbosonsearches:2003ing ; ALEPH:2005ab  following
Asai:2023mzl ; KA:2023dyz , with sharp resonance at the $Z-$pole where $g_{X}$
could reach around $0.001$. In the context of comparing scattering scenarios
we mention that there could be a future muon-electron scattering experiment
called MUonE Asai:2023mzl  which could also provide stronger limit for
$x_{H}=-2$ and hence can probe the strongest parameter region obtained from
electron-DM scattering in XENON1T-S2 experiment in case of Dirac-DM scenario.
The limits obtained from XENON1T-S2 for Dirac DM scenario could provide
stonger sensitivity than MUonE experiment for 0.01 GeV $\lesssim
M_{Z^{\prime}}\lesssim 0.025$ GeV which could be tested in future. However, at
the current stage bounds obtained from Majorana and scalar DM scenarios will
provide weaker constraints than the prospective MUonE experiment as shown in
Fig. 5. In case of $x_{H}=0$ and 2 we find from Fig. 6 that future constraints
from MUonE experiment could be weaker than the existing limits from TEXONO,
BOREXINO, COHERENT experiments, respectively.

The freeze-in lines for each DM kind are shown in Fig. 5 for $x_{H}=-2$ and in
Fig. 6 for $x_{H}=0$ (left panel) and $2$ (right panel) via the blue and brown
dashed curves, corresponding to freeze-in production during radiation
domination and during the epoch of modified cosmology, respectively. In all
cases we consider $m_{\text{DM}}=2\,M_{Z^{\prime}}$, such that DM production
via $s$-channel $Z^{\prime}$-mediated scattering dominates. As elaborated in
Sec. IV, a modified cosmological background can improve the DM-SM coupling by
orders of magnitude. This can be prominently seen for $M_{Z^{\prime}}\gtrsim
0.1$ GeV, while for lighter $M_{Z^{\prime}}$ the DM abundance becomes
insensitive to the choice of the mediator mass. Following the prescription in
subsection IV.2, a larger $n$ requires larger $g_{X}$ to produce the right DM
abundance. As an instance, we see, for $n=3$ with $x_{H}=-2$, one needs
$g_{X}\simeq 4\times 10^{-4}$ for $M_{Z^{\prime}}=10$ GeV, while in case of
standard cosmological background, $g_{X}\simeq 2\times 10^{-5}$ is required.
For different $x_{H}$, it is not possible to go beyond $n=3$, as in that case,
$g_{X}\gtrsim 0.1$ is needed to satisfy the observed DM abundance. For such a
large coupling, the out of equilibrium production of DM becomes a question,
and therefore freeze-in remains valid no more. In case of Majorana and scalar
DM, we see, bounds form LHCb Foguel:2022unm ; Baruch:2022esd , BaBaR
BaBar:2014zli ; BaBar:2017tiz , CMS CMS:2023slr  and LEP-II searches are
already probing $n>1$ scenario for different choices of $x_{H}$ when
$M_{Z^{\prime}}\gtrsim 1$ GeV. For Dirac DM, however, $n=3$ remains still
below bounds obtained from LHCb search. This is attributed to the large Dirac
DM charge that in turn demands smaller $g_{X}$ to ensure that the DM does not
overclose the Universe. For lighter $M_{Z^{\prime}}<1$ GeV, electron/positron
beam dump experiments, for example, Orsay Davier:1989wz , NA64 NA64:2019auh ,
KEK Beer:1986qr , E141 Riordan:1987aw , E137 Bjorken:1988as , E774
Bross:1989mp , are providing strong bounds on $g_{X}$ irrespective of the
nature of the DM and for both modified $(n>0)$ and unmodified $(n=0)$
cosmological scenarios. Even after substantial improvement in the DM-SM
coupling because of alternative cosmological scenario prior to BBN, the
freeze-in lines are well beyond the reach of direct search limits in all
cases. We also find that proton beam dump experiments involving Nomad
NOMAD:2001xxt , CHARM CHARM:1985anb , $\nu-$cal Blumlein:2011mv ;
Blumlein:2013cua , NA62 NA62:2023qyn  and FASER FASER:2023tle  could probe the
freeze-in lines for $M_{Z^{\prime}}<1$ GeV. In addition to that, we compare
our results with prospective bounds obtained from proton beam dump in FASER2,
DUNE experiments and electron beam dump scenario at ILC (ILC-BD) experiment
Asai:2022zxw  which could probe the freeze-in scenarios in future.

##  VI Conclusions

We consider a general $U(1)_{X}$ extension of the SM where three generations
of RHNs are introduced to generate the neutrino mass through seesaw mechanism
after spontaneous breaking of $U(1)_{X}$. On cancellation of gauge and mixed
gauge-gravity anomalies we find, the $U(1)_{X}$ charge assignments of the left
and right handed SM fermions are different. As a result, they interact
differently with beyond the SM neutral gauge boson $Z^{\prime}$ of the model.
We consider three different charges $x_{H}=-2,$ 0 and $2$, while fixing
$x_{\Phi}=1$. In the first case, $U(1)_{X}$ charges of the left handed SM
fermions vanish manifesting the $U(1)_{R}$ scenario, in the second case, left
and right handed SM fermions interact with $Z^{\prime}$ in the same way
manifesting B$-$L scenario; the third case is a chiral one where left and
right handed SM fermions interact differently with $Z^{\prime}$. To explore
the phenomenological implications of this property, we consider electron-DM
scattering processes mediated by $Z^{\prime}$ for Majorana, scalar and Dirac
type DM, to estimate constraints on $g_{X}-M_{Z^{\prime}}$ plane and compare
with existing bounds from several electron-DM scattering experiments,
considering two benchmark choices of DM masses [cf. Fig. 2]. In order to check
the consistency between theoretical estimation and experimental predictions,
we have also estimated events rates under different experimental set-up for
electron-DM scattering [cf. Fig. 3]. Finally, we find, in case of $x_{H}=-2$
electron-DM scattering provides strong bounds for different DM models. Among
them, strongest bound is obtained from the XENON1T-S2 experiment. However, in
case of $x_{H}=0$ and 2, we find that neutrino-electron and neutrino-nucleon
scattering experiments provide strong constrains on $g_{X}-M_{Z^{\prime}}$
plane which are stronger than the limits estimated from electron-DM
scattering. Therefore depending on the electron-DM scattering experiment,
strongest limits could be provided in $U(1)_{R}$ scenario while the parameter
space is compared to the scattering and beam-dump experiments approximately
for the range 0.01 GeV $\lesssim M_{Z^{\prime}}\lesssim 1$ GeV for Dirac-DM.
However, this window of $M_{Z^{\prime}}$ is narrow in case of Majorana and
Scalar type DM models.

In the same framework, we also study freeze-in production via
$\rm{SM~{}SM}\to\rm{DM~{}DM}$ scattering (mediated by $Z^{\prime}$) in a
radiation dominated and in a modified cosmological background, reproducing
observed DM relic abundance forbidding the decay of $Z^{\prime}$ into DM
candidates [cf. Fig. 4]. In the modified cosmological scenario, we consider
the Universe expands faster compared to the standard radiation domination
prior to the onset of big bang nucleosynthesis (BBN). As before, we estimate
bounds on $g_{X}-M_{Z^{\prime}}$ plane for $x_{H}=-2$, 0 and $2$ taking
different DM candidates into account. We see, existing bounds from electron
and proton beam dump experiments rules out freeze-in predictions for
$M_{Z^{\prime}}\lesssim 0.1$ GeV, whereas depending on $x_{H}$ and nature of
the DM candidates these bounds become weaker for $M_{Z^{\prime}}\gtrsim 10$
GeV. In case of Majorana and scalar DM candidates, bound from freeze-in
production in modified cosmological background is ruled out by existing LHCb
constraints, while the bounds obtained from freeze-in production during
radiation domination stay beyond the current LHCb reach. However, in case of a
Dirac DM, the estimated bounds from freeze-in are always strong for 0.8 GeV
$\lesssim M_{Z^{\prime}}\lesssim 150$ GeV and region for
$M_{Z^{\prime}}\lesssim 0.8$ GeV is ruled out by existing limits form beam
dump experiments. We thus infer, present bounds from beam dump, LHCb, CMS,
BaBaR and LEP can potentially put bound on pre-BBN modification to standard
cosmology, thereby constraining different models of modified gravity and
cosmology. However, freeze-in couplings required to produce the observed DM
abundance remains well beyond the reach of DM-electron scattering experiments
[cf. Fig. 5 and 6]. Proton beam dump experiments FASER2, DUNE and electron
beam dump experiment in ILC could probe parameter regions on
$g_{X}-M_{Z^{\prime}}$ plane for different general $U(1)_{X}$ charges,
providing a complementarity in future.

Acknowledgments.– The work of S.M. is supported by KIAS Individual Grants
(PG086002) at Korea Institute for Advanced Study.

##  Appendix A Kinematics and scattering amplitudes

We consider the DM$-e$ elastic scattering process DM $e^{-}\to e^{-}$ DM in
the $t-$channel mediated by $Z^{\prime}$. We define the four momenta in the
laboratory frame as

$\displaystyle
p_{1}=\left(E_{\text{DM}},\,\vec{p}_{1}\right),~{}~{}p_{2}=\left(m_{e},\,\vec{0}\right),~{}~{}p_{3}=\left(E_{\text{DM}}^{\prime},\,\vec{p}_{3}\right)p_{4}=\left(E_{4},\,\vec{p}_{4}\right)$
(41)

From 4-momenta conservation we further obtain

$\displaystyle E_{\text{DM}}^{\prime}=E_{\text{DM}}+m_{e}-E_{4}\text{ with
}E_{4}=m_{e}+E_{r},$ (42) $\displaystyle p_{1}\cdot p_{2}=p_{3}\cdot
p_{4}=E_{\text{DM}}\,m_{e},~{}~{}p_{3}\cdot p_{2}=p_{4}\cdot
p_{1}=m_{e}\,\left(E_{\text{DM}}+m_{e}-E_{4}\right)$ (43) $\displaystyle
p_{1}\cdot p_{3}=m_{\text{DM}}^{2}-m_{e}^{2}+m_{e}\,E_{4},~{}~{}p_{2}\cdot
p_{4}=m_{e}\,E_{4},~{}~{}q^{2}=\left(p_{4}-p_{2}\right)^{2}=-2m_{e}\,E_{\rm
rec}\,.$ (44)

The differential cross-section in the laboratory frame reads

$\displaystyle\frac{d\sigma}{dE_{\rm
rec}}=\frac{1}{32\,|\vec{p}_{1}|^{2}\,m_{e}}\,\left|\overline{\mathcal{M}}\right|^{2}_{\text{DM}e\to
e\text{DM}}\,,$ (45)

where

$\displaystyle\left|\overline{\mathcal{M}}\right|^{2}_{\chi e\to
e\chi}=\frac{g_{X}^{4}\,Q_{\chi}^{2}\,m_{e}}{\left(2\,E_{\rm
rec}\,m_{e}+M_{Z^{\prime}}^{2}\right)^{2}+M_{Z^{\prime}}^{2}\,\Gamma_{Z^{\prime}}^{2}}\Big{[}\left(5x_{H}^{2}+12x_{H}x_{\Phi}+8x_{\Phi}^{2}\right)\Big{\\{}E_{\text{DM}}^{2}m_{e}+m_{e}(E_{\text{DM}}-E_{\rm
rec})^{2}-$ $\displaystyle m_{\chi}^{2}(E_{\rm
rec}+m_{e})\Big{\\}}-4\,m_{e}\,(x_{H}+x_{\Phi})\,(x_{H}+2x_{\Phi})\left(E_{\rm
rec}m_{e}+m_{\chi}^{2}\right)+8m_{e}m_{\chi}^{2}(x_{H}+x_{\Phi})(x_{H}+2x_{\Phi})\Big{]}\,,$
(46)

is for Dirac DM,

$\displaystyle\left|\overline{\mathcal{M}}\right|^{2}_{Ne\to
eN}=\frac{g_{X}^{4}\,x_{\Phi}^{2}\,m_{e}}{\left(2\,E_{\rm
rec}\,m_{e}+M_{Z^{\prime}}^{2}\right)^{2}+M_{Z^{\prime}}^{2}\,\Gamma_{Z^{\prime}}^{2}}\Big{[}\left(5x_{H}^{2}+12x_{H}x_{\Phi}+8x_{\Phi}^{2}\right)\times\Big{\\{}E_{\text{DM}}^{2}m_{e}+m_{e}(E_{\text{DM}}-E_{\rm
rec})^{2}+$ $\displaystyle M_{N}^{2}(E_{\rm
rec}+m_{e})\Big{\\}}-4m_{e}(x_{H}+x_{\Phi})(x_{H}+2x_{\Phi})\left(E_{\rm
rec}m_{e}+M_{N}^{2}\right)-8m_{e}M_{N}^{2}(x_{H}+x_{\Phi})(x_{H}+2x_{\Phi})\Big{]}\,,$
(47)

is for Majorana DM and

$\displaystyle\left|\overline{\mathcal{M}}\right|^{2}_{\Phi e\to
e\Phi}=\frac{g_{X}^{4}\,x_{\Phi}^{2}\,m_{e}}{\left(2\,E_{\rm
rec}\,m_{e}+M_{Z^{\prime}}^{2}\right)^{2}+M_{Z^{\prime}}^{2}\,\Gamma_{Z^{\prime}}^{2}}\Big{[}x_{H}^{2}\Big{\\{}20E_{\text{DM}}\,m_{e}\,(E_{\text{DM}}-E_{\rm
rec})-E_{\rm rec}m_{e}^{2}-2m_{\Phi}^{2}(5E_{\rm rec}+m_{e})\Big{\\}}+$
$\displaystyle
24x_{H}x_{\Phi}\Big{\\{}2E_{\text{DM}}m_{e}(E_{\text{DM}}-E_{\rm rec})-E_{\rm
rec}m_{\Phi}^{2}\Big{\\}}+16x_{\Phi}^{2}\Big{\\{}2E_{\text{DM}}m_{e}(E_{\text{DM}}-E_{\rm
rec})-E_{\rm rec}m_{\Phi}^{2}\Big{\\}}\Big{]}\,,$ (48)

in case of scalar DM.

##  Appendix B Scattering cross-section

Here we report expressions for the scattering cross-sections in the center of
mass frame of a pair of massless SM fermions $(f)$ going into a pair of DM
particles as $f\bar{f}\to$ DM DM, mediated by $s$-channel $Z^{\prime}$ where
$f=\ell$ (charged lepton), $\nu$ (neutrino), $u$ (up quark) and $d$ (down
quark), respectively. For charged lepton initial states we have

$\displaystyle\sigma(s)_{\ell^{+}\ell^{-}\to\text{DM}\text{DM}}\simeq\frac{g_{X}^{4}}{96\pi\,\left[(s-M_{Z^{\prime}}^{2})^{2}+\Gamma_{Z^{\prime}}^{2}\,M_{Z^{\prime}}^{2}\right]}\,\sqrt{1-\frac{4\,m_{\text{DM}}^{2}}{s}}\begin{cases}Q_{\chi}^{2}\,\left(s+2\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}+12\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)&\text{for
Dirac}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}+12\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{2}&\text{for
Majorana}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}+12\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{4}&\text{for
Scalar}\,.\end{cases}$ (49)

For light neutrino initial states

$\displaystyle\sigma(s)_{\nu\nu\to\text{DM}\text{DM}}\simeq\frac{g_{X}^{4}}{48\pi\,\left[(s-M_{Z^{\prime}}^{2})^{2}+\Gamma_{Z^{\prime}}^{2}\,M_{Z^{\prime}}^{2}\right]}\,\sqrt{1-\frac{4\,m_{\text{DM}}^{2}}{s}}\begin{cases}Q_{\chi}^{2}\,\left(s+2\,m_{\text{DM}}^{2}\right)\,\left(x_{H}+2\,x_{\Phi}\right)^{2}&\text{for
Dirac}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(x_{H}+2\,x_{\Phi}\right)^{2}}{4}&\text{for
Majorana}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(x_{H}+2\,x_{\Phi}\right)^{2}}{4}&\text{for
Scalar}\,.\end{cases}$ (50)

For up-quark initial states,

$\displaystyle\sigma(s)_{uu\to\text{DM}\text{DM}}\simeq\frac{g_{X}^{4}}{2592\pi\,\left[(s-M_{Z^{\prime}}^{2})^{2}+\Gamma_{Z^{\prime}}^{2}\,M_{Z^{\prime}}^{2}\right]}\,\sqrt{1-\frac{4\,m_{\text{DM}}^{2}}{s}}\begin{cases}Q_{\chi}^{2}\,\left(s+2\,m_{\text{DM}}^{2}\right)\,\left(17\,x_{H}^{2}+20\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)&\text{for
Dirac}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(17\,x_{H}^{2}+20\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{2}&\text{for
Majorana}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(17\,x_{H}^{2}+20\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{4}&\text{for
Scalar}\,.\end{cases}$ (51)

For down-quark initial states,

$\displaystyle\sigma(s)_{dd\to\text{DM}\text{DM}}\simeq\frac{g_{X}^{4}}{2592\pi\,\left[(s-M_{Z^{\prime}}^{2})^{2}+\Gamma_{Z^{\prime}}^{2}\,M_{Z^{\prime}}^{2}\right]}\,\sqrt{1-\frac{4\,m_{\text{DM}}^{2}}{s}}\begin{cases}Q_{\chi}^{2}\,\left(s+2\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}-4\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)&\text{for
Dirac}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}-4\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{2}&\text{for
Majorana}\\\\[10.0pt]
x_{\Phi}^{2}\,\frac{\left(s-4\,m_{\text{DM}}^{2}\right)\,\left(5\,x_{H}^{2}-4\,x_{H}\,x_{\Phi}+8\,x_{\Phi}^{2}\right)}{4}&\text{for
Scalar}\,.\end{cases}$ (52)

##  Appendix C Reaction densities

The reaction density for to $1\to 2$ decay process is given by

$\displaystyle\gamma_{12}$
$\displaystyle=\int\,\prod_{i=1}^{3}\,(2\,\pi)^{4}\,\delta^{(4)}\,\left(p_{a}-p_{1}-p_{2}\right)\,f_{a}^{\text{eq}}\,\left|\mathcal{M}\right|^{2}_{a\to
1,2}=\frac{g_{a}}{2\,\pi^{2}}\,m_{a}^{2}\,\Gamma_{a\to
1,2}\,T\,K_{1}\left(\frac{m_{a}}{T}\right)\,.$ (53)

The reaction density corresponding to 2-to-2 processes reads

$\displaystyle\gamma_{22}=\int\prod_{i=1}^{4}d\Pi_{i}\left(2\pi\right)^{4}\delta^{(4)}\biggl{(}p_{a}+p_{b}-p_{1}-p_{2}\biggr{)}f_{a}{{}^{\text{eq}}}f_{b}{{}^{\text{eq}}}\left|\mathcal{M}_{a,b\to
1,2}\right|^{2}$
$\displaystyle=\frac{T}{32\pi^{4}}\,g_{a}g_{b}\,\int_{s_{\text{min}}}^{\infty}ds\,\frac{\biggl{[}\bigl{(}s-m_{a}^{2}-m_{b}^{2}\bigr{)}^{2}-4m_{a}^{2}m_{b}^{2}\biggr{]}}{\sqrt{s}}\sigma\left(s\right)_{a,b\to
1,2}\,K_{1}\left(\frac{\sqrt{s}}{T}\right)\,,$ (54)

with $a,b(1,2)$ as the incoming (outgoing) states and $g_{a,b}$ are
corresponding degrees of freedom. Here
$f_{i}{{}^{\text{eq}}}\approx\exp^{-E_{i}/T}$ is the Maxwell-Boltzmann
distribution. The Lorentz invariant 2-body phase space is denoted by:
$d\Pi_{i}=\frac{d^{3}p_{i}}{\left(2\pi\right)^{3}2E_{i}}$. The amplitude
squared (summed over final and averaged over initial states) is denoted by
$\left|\mathcal{M}_{a,b\to 1,2}\right|^{2}$ for a particular 2-to-2 scattering
process. The lower limit of the integration over $s$ is
$s_{\text{min}}=\text{max}\left[\left(m_{a}+m_{b}\right)^{2},\left(m_{1}+m_{2}\right)^{2}\right]$.

## References

  * (1) Particle Data Group Collaboration, P. A. Zyla et al., “Review of Particle Physics,” PTEP 2020 no. 8, (2020) 083C01.
  * (2) G. Bertone and D. Hooper, “History of dark matter,” Rev. Mod. Phys. 90 no. 4, (2018) 045002, arXiv:1605.04909 [astro-ph.CO].
  * (3) J. de Swart, G. Bertone, and J. van Dongen, “How Dark Matter Came to Matter,” Nature Astron. 1 (2017) 0059, arXiv:1703.00013 [astro-ph.CO].
  * (4) Planck Collaboration, N. Aghanim et al., “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641 (2020) A6, arXiv:1807.06209 [astro-ph.CO]. [Erratum: Astron.Astrophys. 652, C4 (2021)].
  * (5) P. Minkowski, “$\mu\to e\gamma$ at a Rate of One Out of $10^{9}$ Muon Decays?,” Phys. Lett. B 67 (1977) 421–428.
  * (6) R. N. Mohapatra and G. Senjanovic, “Neutrino Mass and Spontaneous Parity Nonconservation,” Phys. Rev. Lett. 44 (1980) 912.
  * (7) M. Gell-Mann, P. Ramond, and R. Slansky, “Complex Spinors and Unified Theories,” Conf. Proc. C 790927 (1979) 315–321, arXiv:1306.4669 [hep-th].
  * (8) M. Magg and C. Wetterich, “Neutrino Mass Problem and Gauge Hierarchy,” Phys. Lett. B 94 (1980) 61–64.
  * (9) J. Schechter and J. W. F. Valle, “Neutrino Masses in SU(2) x U(1) Theories,” Phys. Rev. D 22 (1980) 2227.
  * (10) C. Wetterich, “Neutrino Masses and the Scale of B-L Violation,” Nucl. Phys. B 187 (1981) 343–375.
  * (11) G. Lazarides, Q. Shafi, and C. Wetterich, “Proton Lifetime and Fermion Masses in an SO(10) Model,” Nucl. Phys. B 181 (1981) 287–300.
  * (12) R. N. Mohapatra and G. Senjanovic, “Neutrino Masses and Mixings in Gauge Models with Spontaneous Parity Violation,” Phys. Rev. D 23 (1981) 165.
  * (13) R. Foot, H. Lew, X. G. He, and G. C. Joshi, “Seesaw Neutrino Masses Induced by a Triplet of Leptons,” Z. Phys. C 44 (1989) 441.
  * (14) E. Ma, “Pathways to naturally small neutrino masses,” Phys. Rev. Lett. 81 (1998) 1171–1174, arXiv:hep-ph/9805219.
  * (15) E. Ma and D. P. Roy, “Heavy triplet leptons and new gauge boson,” Nucl. Phys. B 644 (2002) 290–302, arXiv:hep-ph/0206150.
  * (16) T. Hambye, Y. Lin, A. Notari, M. Papucci, and A. Strumia, “Constraints on neutrino masses from leptogenesis models,” Nucl. Phys. B 695 (2004) 169–191, arXiv:hep-ph/0312203.
  * (17) B. Bajc and G. Senjanovic, “Seesaw at LHC,” JHEP 08 (2007) 014, arXiv:hep-ph/0612029.
  * (18) E. Ma, “Verifiable radiative seesaw mechanism of neutrino mass and dark matter,” Phys. Rev. D 73 (2006) 077301, arXiv:hep-ph/0601225.
  * (19) G. Steigman and M. S. Turner, “Cosmological Constraints on the Properties of Weakly Interacting Massive Particles,” Nucl. Phys. B 253 (1985) 375–386.
  * (20) G. Jungman, M. Kamionkowski, and K. Griest, “Supersymmetric dark matter,” Phys. Rept. 267 (1996) 195–373, arXiv:hep-ph/9506380 [hep-ph].
  * (21) G. Bertone, D. Hooper, and J. Silk, “Particle dark matter: Evidence, candidates and constraints,” Phys. Rept. 405 (2005) 279–390, arXiv:hep-ph/0404175.
  * (22) J. L. Feng, “Dark Matter Candidates from Particle Physics and Methods of Detection,” Ann. Rev. Astron. Astrophys. 48 (2010) 495–545, arXiv:1003.0904 [astro-ph.CO].
  * (23) L. Roszkowski, E. M. Sessolo, and S. Trojanowski, “WIMP dark matter candidates and searches—current status and future prospects,” Rept. Prog. Phys. 81 no. 6, (2018) 066201, arXiv:1707.06277 [hep-ph].
  * (24) G. Arcadi, M. Dutra, P. Ghosh, M. Lindner, Y. Mambrini, M. Pierre, S. Profumo, and F. S. Queiroz, “The waning of the WIMP? A review of models, searches, and constraints,” Eur. Phys. J. C 78 no. 3, (2018) 203, arXiv:1703.07364 [hep-ph].
  * (25) J. McDonald, “Thermally generated gauge singlet scalars as selfinteracting dark matter,” Phys. Rev. Lett. 88 (2002) 091304, arXiv:hep-ph/0106249.
  * (26) L. J. Hall, K. Jedamzik, J. March-Russell, and S. M. West, “Freeze-In Production of FIMP Dark Matter,” JHEP 03 (2010) 080, arXiv:0911.1120 [hep-ph].
  * (27) N. Bernal, M. Heikinheimo, T. Tenkanen, K. Tuominen, and V. Vaskonen, “The Dawn of FIMP Dark Matter: A Review of Models and Constraints,” Int. J. Mod. Phys. A 32 no. 27, (2017) 1730023, arXiv:1706.07442 [hep-ph].
  * (28) M. W. Goodman and E. Witten, “Detectability of Certain Dark Matter Candidates,” Phys. Rev. D 31 (1985) 3059.
  * (29) M. Schumann, “Direct Detection of WIMP Dark Matter: Concepts and Status,” J. Phys. G 46 no. 10, (2019) 103003, arXiv:1903.03026 [astro-ph.CO].
  * (30) R. Essig, J. Mardon, and T. Volansky, “Direct Detection of Sub-GeV Dark Matter,” Phys. Rev. D 85 (2012) 076007, arXiv:1108.5383 [hep-ph].
  * (31) T. Lin, “Dark matter models and direct detection,” PoS 333 (2019) 009, arXiv:1904.07915 [hep-ph].
  * (32) M. Battaglieri et al., “US Cosmic Visions: New Ideas in Dark Matter 2017: Community Report,” in U.S. Cosmic Visions: New Ideas in Dark Matter.  7, 2017.  arXiv:1707.04591 [hep-ph].
  * (33) A. Das, S. Oda, N. Okada, and D.-s. Takahashi, “Classically conformal U(1)’ extended standard model, electroweak vacuum stability, and LHC Run-2 bounds,” Phys. Rev. D 93 no. 11, (2016) 115038, arXiv:1605.01157 [hep-ph].
  * (34) R. Essig, M. Fernandez-Serra, J. Mardon, A. Soto, T. Volansky, and T.-T. Yu, “Direct Detection of sub-GeV Dark Matter with Semiconductor Targets,” JHEP 05 (2016) 046, arXiv:1509.01598 [hep-ph].
  * (35) A. Radick, A.-M. Taki, and T.-T. Yu, “Dependence of Dark Matter - Electron Scattering on the Galactic Dark Matter Velocity Distribution,” JCAP 02 (2021) 004, arXiv:2011.02493 [hep-ph].
  * (36) J. Kopp, V. Niro, T. Schwetz, and J. Zupan, “DAMA/LIBRA and leptonically interacting Dark Matter,” Phys. Rev. D 80 (2009) 083502, arXiv:0907.3159 [hep-ph].
  * (37) SuperCDMS Collaboration, D. W. Amaral et al., “Constraints on low-mass, relic dark matter candidates from a surface-operated SuperCDMS single-charge sensitive detector,” Phys. Rev. D 102 no. 9, (2020) 091101, arXiv:2005.14067 [hep-ex].
  * (38) DAMIC Collaboration, A. Aguilar-Arevalo et al., “Constraints on Light Dark Matter Particles Interacting with Electrons from DAMIC at SNOLAB,” Phys. Rev. Lett. 123 no. 18, (2019) 181802, arXiv:1907.12628 [astro-ph.CO].
  * (39) SENSEI Collaboration, L. Barak et al., “SENSEI: Direct-Detection Results on sub-GeV Dark Matter from a New Skipper-CCD,” Phys. Rev. Lett. 125 no. 17, (2020) 171802, arXiv:2004.11378 [astro-ph.CO].
  * (40) PandaX-II Collaboration, C. Cheng et al., “Search for Light Dark Matter-Electron Scatterings in the PandaX-II Experiment,” Phys. Rev. Lett. 126 no. 21, (2021) 211803, arXiv:2101.07479 [hep-ex].
  * (41) DarkSide Collaboration, P. Agnes et al., “Search for Dark Matter Particle Interactions with Electron Final States with DarkSide-50,” Phys. Rev. Lett. 130 no. 10, (2023) 101002, arXiv:2207.11968 [hep-ex].
  * (42) XENON Collaboration, E. Aprile et al., “Light Dark Matter Search with Ionization Signals in XENON1T,” Phys. Rev. Lett. 123 no. 25, (2019) 251801, arXiv:1907.11485 [hep-ex].
  * (43) C. F. Bunge, J. A. Barrientos, and A. V. Bunge, “Roothaan-Hartree-Fock Ground-State Atomic Wave Functions: Slater-Type Orbital Expansions and Expectation Values for Z = 2-54,” Atom. Data Nucl. Data Tabl. 53 (1993) 113–162.
  * (44) T. Trickle, “Extended calculation of electronic excitations for direct detection of dark matter,” Phys. Rev. D 107 no. 3, (2023) 035035, arXiv:2210.14917 [hep-ph].
  * (45) S. M. Griffin, K. Inzani, T. Trickle, Z. Zhang, and K. M. Zurek, “Extended calculation of dark matter-electron scattering in crystal targets,” Phys. Rev. D 104 no. 9, (2021) 095015, arXiv:2105.05253 [hep-ph].
  * (46) P. Giannozzi et al., “QUANTUM ESPRESSO: a modular and open-source software project for quantumsimulations of materials,” J. Phys. Condens. Matter 21 no. 39, (2009) 395502.
  * (47) P. W. Graham, D. E. Kaplan, S. Rajendran, and M. T. Walters, “Semiconductor Probes of Light Dark Matter,” Phys. Dark Univ. 1 (2012) 32–49, arXiv:1203.2531 [hep-ph].
  * (48) S. K. Lee, M. Lisanti, S. Mishra-Sharma, and B. R. Safdi, “Modulation Effects in Dark Matter-Electron Scattering Experiments,” Phys. Rev. D 92 no. 8, (2015) 083517, arXiv:1508.07361 [hep-ph].
  * (49) A. Bohm, A. Mostafazadeh, H. Koizumi, Q. Niu, and J. Zwanziger, The Geometric Phase in Quantum Systems: Foundations, Mathematical Concepts, and Applications in Molecular and Condensed Matter Physics.  04, 2003.
  * (50) D. J. Chadi and M. L. Cohen, “Tight-binding calculations of the valence bands of diamond and zincblende crystals,” physica status solidi (b) 68 no. 1, (1975) 405–419. https://onlinelibrary.wiley.com/doi/abs/10.1002/pssb.2220680140.
  * (51) M. Duch, B. Grzadkowski, and D. Huang, “Strongly self-interacting vector dark matter via freeze-in,” JHEP 01 (2018) 020, arXiv:1710.00320 [hep-ph].
  * (52) B. Barman, S. Bhattacharya, and B. Grzadkowski, “Feebly coupled vector boson dark matter in effective theory,” JHEP 12 (2020) 162, arXiv:2009.07438 [hep-ph].
  * (53) S. Bhattacharya, S. Chakraborti, and D. Pradhan, “Electroweak symmetry breaking and WIMP-FIMP dark matter,” JHEP 07 (2022) 091, arXiv:2110.06985 [hep-ph].
  * (54) Particle Data Group Collaboration, R. L. Workman et al., “Review of Particle Physics,” PTEP 2022 (2022) 083C01.
  * (55) F. D’Eramo, N. Fernandez, and S. Profumo, “When the Universe Expands Too Fast: Relentless Dark Matter,” JCAP 05 (2017) 012, arXiv:1703.04793 [hep-ph].
  * (56) F. D’Eramo, N. Fernandez, and S. Profumo, “Dark Matter Freeze-in Production in Fast-Expanding Universes,” JCAP 02 (2018) 046, arXiv:1712.07453 [hep-ph].
  * (57) S.-L. Chen, A. Dutta Banik, and Z.-K. Liu, “Leptogenesis in fast expanding Universe,” JCAP 03 (2020) 009, arXiv:1912.07185 [hep-ph].
  * (58) D. Mahanta and D. Borah, “TeV Scale Leptogenesis with Dark Matter in Non-standard Cosmology,” JCAP 04 no. 04, (2020) 032, arXiv:1912.09726 [hep-ph].
  * (59) P. Konar, A. Mukherjee, A. K. Saha, and S. Show, “A dark clue to seesaw and leptogenesis in a pseudo-Dirac singlet doublet scenario with (non)standard cosmology,” JHEP 03 (2021) 044, arXiv:2007.15608 [hep-ph].
  * (60) G. Arcadi, J. P. Neto, F. S. Queiroz, and C. Siqueira, “Roads for right-handed neutrino dark matter: Fast expansion, standard freeze-out, and early matter domination,” Phys. Rev. D 105 no. 3, (2022) 035016, arXiv:2108.11398 [hep-ph].
  * (61) B. Barman, P. Ghosh, F. S. Queiroz, and A. K. Saha, “Scalar multiplet dark matter in a fast expanding Universe: Resurrection of the desert region,” Phys. Rev. D 104 no. 1, (2021) 015040, arXiv:2101.10175 [hep-ph].
  * (62) D. Mahanta and D. Borah, “WIMPy leptogenesis in non-standard cosmologies,” JCAP 03 (2023) 049, arXiv:2208.11295 [hep-ph].
  * (63) B. Barman and A. Ghoshal, “Probing pre-BBN era with scale invariant FIMP,” JCAP 10 (2022) 082, arXiv:2203.13269 [hep-ph].
  * (64) E. Hardy, “Higgs portal dark matter in non-standard cosmological histories,” JHEP 06 (2018) 043, arXiv:1804.06783 [hep-ph].
  * (65) N. Bernal, C. Cosme, T. Tenkanen, and V. Vaskonen, “Scalar singlet dark matter in non-standard cosmologies,” Eur. Phys. J. C 79 no. 1, (2019) 30, arXiv:1806.11122 [hep-ph].
  * (66) C. Cosme, M. Dutra, T. Ma, Y. Wu, and L. Yang, “Neutrino portal to FIMP dark matter with an early matter era,” JHEP 03 (2021) 026, arXiv:2003.01723 [hep-ph].
  * (67) R. R. Caldwell, R. Dave, and P. J. Steinhardt, “Cosmological imprint of an energy component with general equation of state,” Phys. Rev. Lett. 80 (1998) 1582–1585, arXiv:astro-ph/9708069.
  * (68) V. Sahni and A. A. Starobinsky, “The Case for a positive cosmological Lambda term,” Int. J. Mod. Phys. D 9 (2000) 373–444, arXiv:astro-ph/9904398.
  * (69) B. Ratra and P. J. E. Peebles, “Cosmological Consequences of a Rolling Homogeneous Scalar Field,” Phys. Rev. D 37 (1988) 3406.
  * (70) E. J. Copeland, A. R. Liddle, and D. Wands, “Exponential potentials and cosmological scaling solutions,” Phys. Rev. D 57 (1998) 4686–4690, arXiv:gr-qc/9711068.
  * (71) E. I. Buchbinder, J. Khoury, and B. A. Ovrut, “New Ekpyrotic cosmology,” Phys. Rev. D 76 (2007) 123503, arXiv:hep-th/0702154.
  * (72) P. K. S. Dunsby, B. A. C. C. Bassett, and G. F. R. Ellis, “Covariant analysis of gravitational waves in a cosmological context,” Class. Quant. Grav. 14 (1997) 1215–1222, arXiv:gr-qc/9811092.
  * (73) R. Catena, N. Fornengo, M. Pato, L. Pieri, and A. Masiero, “Thermal Relics in Modified Cosmologies: Bounds on Evolution Histories of the Early Universe and Cosmological Boosts for PAMELA,” Phys. Rev. D 81 (2010) 123522, arXiv:0912.4421 [astro-ph.CO].
  * (74) J. B. Dent, S. Dutta, and R. J. Scherrer, “Thermal Relic Abundances of Particles with Velocity-Dependent Interactions,” Phys. Lett. B 687 (2010) 275–279, arXiv:0909.4128 [astro-ph.CO].
  * (75) G. Leon, J. Saavedra, and E. N. Saridakis, “Cosmological behavior in extended nonlinear massive gravity,” Class. Quant. Grav. 30 (2013) 135001, arXiv:1301.7419 [astro-ph.CO].
  * (76) V. Baules, N. Okada, and S. Okada, “Braneworld cosmological effect on freeze-in dark matter density and lifetime frontier,” Eur. Phys. J. C 82 no. 7, (2022) 643, arXiv:1911.05344 [hep-ph].
  * (77) NA64 Collaboration, D. Banerjee et al., “Improved limits on a hypothetical X(16.7) boson and a dark photon decaying into $e^{+}e^{-}$ pairs,” Phys. Rev. D 101 no. 7, (2020) 071101, arXiv:1912.11389 [hep-ex].
  * (78) BaBar Collaboration, J. P. Lees et al., “Search for a Dark Photon in $e^{+}e^{-}$ Collisions at BaBar,” Phys. Rev. Lett. 113 no. 20, (2014) 201801, arXiv:1406.2980 [hep-ex].
  * (79) BaBar Collaboration, J. P. Lees et al., “Search for Invisible Decays of a Dark Photon Produced in ${e}^{+}{e}^{-}$ Collisions at BaBar,” Phys. Rev. Lett. 119 no. 13, (2017) 131804, arXiv:1702.03327 [hep-ex].
  * (80) E. M. Riordan et al., “A Search for Short Lived Axions in an Electron Beam Dump Experiment,” Phys. Rev. Lett. 59 (1987) 755.
  * (81) K. Asai, A. Das, J. Li, T. Nomura, and O. Seto, “Probing for chiral Z’ gauge boson through scattering measurement experiments,” Phys. Rev. D 109 no. 7, (2024) 075026, arXiv:2307.09737 [hep-ph].
  * (82) A. L. Foguel, G. M. Salla, and R. Z. Funchal, “(In)Visible signatures of the minimal dark abelian gauge sector,” JHEP 12 (2022) 063, arXiv:2209.03383 [hep-ph].
  * (83) C. Baruch, P. Ilten, Y. Soreq, and M. Williams, “Axial vectors in DarkCast,” JHEP 11 (2022) 124, arXiv:2206.08563 [hep-ph].
  * (84) CMS Collaboration, “Search for prompt production of a GeV scale resonance decaying to a pair of muons in proton-proton collisions at $\sqrts=13~\mathrmTeV,”.
  * (85) TEXONO Collaboration, M. Deniz et al., “Measurement of Nu(e)-bar -Electron Scattering Cross-Section with a CsI(Tl) Scintillating Crystal Array at the Kuo-Sheng Nuclear Power Reactor,” Phys. Rev. D 81 (2010) 072001, arXiv:0911.1597 [hep-ex].
  * (86) TEXONO Collaboration, H. B. Li et al., “Limit on the electron neutrino magnetic moment from the Kuo-Sheng reactor neutrino experiment,” Phys. Rev. Lett. 90 (2003) 131802, arXiv:hep-ex/0212003.
  * (87) TEXONO Collaboration, H. T. Wong et al., “A Search of Neutrino Magnetic Moments with a High-Purity Germanium Detector at the Kuo-Sheng Nuclear Power Station,” Phys. Rev. D 75 (2007) 012001, arXiv:hep-ex/0605006.
  * (88) Borexino Collaboration, G. Alimonti et al., “Science and technology of BOREXINO: A Real time detector for low-energy solar neutrinos,” Astropart. Phys. 16 (2002) 205–234, arXiv:hep-ex/0012030.
  * (89) G. Bellini et al., “Precision measurement of the 7Be solar neutrino interaction rate in Borexino,” Phys. Rev. Lett. 107 (2011) 141302, arXiv:1104.1816 [hep-ex].
  * (90) B. T. Cleveland, T. Daily, R. Davis, Jr., J. R. Distel, K. Lande, C. K. Lee, P. S. Wildenhain, and J. Ullman, “Measurement of the solar electron neutrino flux with the Homestake chlorine detector,” Astrophys. J. 496 (1998) 505–526.
  * (91) K. Lande and P. Wildenhain, “The Homestake chlorine solar neutrino experiment: Past, present and future,” Nucl. Phys. B Proc. Suppl. 118 (2003) 49–54.
  * (92) Borexino Collaboration, C. Arpesella et al., “First real time detection of Be-7 solar neutrinos by Borexino,” Phys. Lett. B 658 (2008) 101–108, arXiv:0708.2251 [astro-ph].
  * (93) Borexino Collaboration, G. Alimonti et al., “The Borexino detector at the Laboratori Nazionali del Gran Sasso,” Nucl. Instrum. Meth. A 600 (2009) 568–593, arXiv:0806.2400 [physics.ins-det].
  * (94) CHARM-II Collaboration, D. Geiregat et al., “An Improved determination of the electroweak mixing angle from muon-neutrino electron scattering,” Phys. Lett. B 259 (1991) 499–507.
  * (95) CHARM-II Collaboration, K. De Winter et al., “A Detector for the Study of Neutrino - Electron Scattering,” Nucl. Instrum. Meth. A 278 (1989) 670.
  * (96) CHARM-II Collaboration, P. Vilain et al., “Measurement of differential cross-sections for muon-neutrino electron scattering,” Phys. Lett. B 302 (1993) 351–355.
  * (97) CHARM-II Collaboration, P. Vilain et al., “Precision measurement of electroweak parameters from the scattering of muon-neutrinos on electrons,” Phys. Lett. B 335 (1994) 246–252.
  * (98) COHERENT Collaboration, D. Akimov et al., “COHERENT Collaboration data release from the first observation of coherent elastic neutrino-nucleus scattering,” arXiv:1804.09459 [nucl-ex].
  * (99) M. Cadeddu, C. Giunti, Y. F. Li, and Y. Y. Zhang, “Average CsI neutron density distribution from COHERENT data,” Phys. Rev. Lett. 120 no. 7, (2018) 072501, arXiv:1710.02730 [hep-ph].
  * (100) M. Cadeddu, N. Cargioli, F. Dordei, C. Giunti, Y. F. Li, E. Picciau, and Y. Y. Zhang, “Constraints on light vector mediators through coherent elastic neutrino nucleus scattering data from COHERENT,” JHEP 01 (2021) 116, arXiv:2008.05022 [hep-ph].
  * (101) COHERENT Collaboration, D. Akimov et al., “First Measurement of Coherent Elastic Neutrino-Nucleus Scattering on Argon,” Phys. Rev. Lett. 126 no. 1, (2021) 012002, arXiv:2003.10630 [nucl-ex].
  * (102) COHERENT Collaboration, D. Akimov et al., “COHERENT Collaboration data release from the first detection of coherent elastic neutrino-nucleus scattering on argon,” arXiv:2006.12659 [nucl-ex].
  * (103) A. G. Beda, E. V. Demidova, A. S. Starostin, V. B. Brudanin, V. G. Egorov, D. V. Medvedev, M. V. Shirchenko, and T. Vylov, “GEMMA experiment: Three years of the search for the neutrino magnetic moment,” Phys. Part. Nucl. Lett. 7 (2010) 406–409, arXiv:0906.1926 [hep-ex].
  * (104) M. Lindner, F. S. Queiroz, W. Rodejohann, and X.-J. Xu, “Neutrino-electron scattering: general constraints on $Z^{\prime}$ and dark photon models,” JHEP 05 (2018) 098, arXiv:1803.00060 [hep-ph].
  * (105) ALEPH, DELPHI, L3, OPAL, LEP Electroweak Working Group Collaboration, J. Alcaraz et al., “A Combination of preliminary electroweak measurements and constraints on the standard model,” arXiv:hep-ex/0612034.
  * (106) ALEPH Collaboration, R. Barate et al., “Study of the muon pair production at center-of-mass energies from 20-GeV to 136-GeV with the ALEPH detector,” Phys. Lett. B 399 (1997) 329–341.
  * (107) LEP Working Group for Higgs boson searches, ALEPH, DELPHI, L3, OPAL Collaboration, R. Barate et al., “Search for the standard model Higgs boson at LEP,” Phys. Lett. B 565 (2003) 61–75, arXiv:hep-ex/0306033.
  * (108) ALEPH, DELPHI, L3, OPAL, SLD, LEP Electroweak Working Group, SLD Electroweak Group, SLD Heavy Flavour Group Collaboration, S. Schael et al., “Precision electroweak measurements on the $Z$ resonance,” Phys. Rept. 427 (2006) 257–454, arXiv:hep-ex/0509008.
  * (109) S. K. A., A. Das, G. Lambiase, T. Nomura, and Y. Orikasa, “Probing chiral and flavored $Z^{\prime}$ from cosmic bursts through neutrino interactions,” arXiv:2308.14483 [hep-ph].
  * (110) M. Davier and H. Nguyen Ngoc, “An Unambiguous Search for a Light Higgs Boson,” Phys. Lett. B 229 (1989) 150–155.
  * (111) G. A. Beer et al., “Emission of Muonium Into Vacuum From a Silica Powder Layer,” Phys. Rev. Lett. 57 (1986) 671–674.
  * (112) J. D. Bjorken, S. Ecklund, W. R. Nelson, A. Abashian, C. Church, B. Lu, L. W. Mo, T. A. Nunamaker, and P. Rassmann, “Search for Neutral Metastable Penetrating Particles Produced in the SLAC Beam Dump,” Phys. Rev. D 38 (1988) 3375.
  * (113) A. Bross, M. Crisler, S. H. Pordes, J. Volk, S. Errede, and J. Wrbanek, “A Search for Shortlived Particles Produced in an Electron Beam Dump,” Phys. Rev. Lett. 67 (1991) 2942–2945.
  * (114) NOMAD Collaboration, P. Astier et al., “Final NOMAD results on muon-neutrino —$>$ tau-neutrino and electron-neutrino —$>$ tau-neutrino oscillations including a new search for tau-neutrino appearance using hadronic tau decays,” Nucl. Phys. B 611 (2001) 3–39, arXiv:hep-ex/0106102.
  * (115) CHARM Collaboration, F. Bergsma et al., “Search for Axion Like Particle Production in 400-GeV Proton - Copper Interactions,” Phys. Lett. B 157 (1985) 458–462.
  * (116) J. Blumlein and J. Brunner, “New Exclusion Limits for Dark Gauge Forces from Beam-Dump Data,” Phys. Lett. B 701 (2011) 155–159, arXiv:1104.2747 [hep-ex].
  * (117) J. Blümlein and J. Brunner, “New Exclusion Limits on Dark Gauge Forces from Proton Bremsstrahlung in Beam-Dump Data,” Phys. Lett. B 731 (2014) 320–326, arXiv:1311.3870 [hep-ph].
  * (118) NA62 Collaboration, E. Cortina Gil et al., “Search for dark photon decays to $\mu$+$\mu$- at NA62,” JHEP 09 (2023) 035, arXiv:2303.08666 [hep-ex].
  * (119) FASER Collaboration, H. Abreu et al., “Search for dark photons with the FASER detector at the LHC,” Phys. Lett. B 848 (2024) 138378, arXiv:2308.05587 [hep-ex].
  * (120) K. Asai, A. Das, J. Li, T. Nomura, and O. Seto, “Chiral Z’ in FASER, FASER2, DUNE, and ILC beam dump experiments,” Phys. Rev. D 106 no. 9, (2022) 095033, arXiv:2206.12676 [hep-ph].