# Exploring dark sector parameters in light of neutron star temperatures

Guey-Lin Lin<EMAIL_ADDRESS>Institute of Physics, National Yang Ming Chiao
Tung University, Hsinchu 30010, Taiwan  Yen-Hsun Lin
<EMAIL_ADDRESS>Institute of Physics, National Yang Ming Chiao Tung
University, Hsinchu 30010, Taiwan  Institute of Physics, Academia Sinica,
Taipei 11529, Taiwan

###### Abstract

Neutron star (NS) as the dark matter (DM) probe has gained a broad attention
recently, either from heating due to DM annihilation or its stability under
the presence of DM. In this work, we investigate spin-$1/2$ fermionic DM
$\chi$ charged under the $U(1)_{X}$ in the dark sector. The massive gauge
boson $V$ of $U(1)_{X}$ gauge group can be produced in NS via DM annihilation.
The produced gauge boson can decay into Standard Model (SM) particles before
it exits NS, despite its tiny couplings to SM particles. Thus, we perform a
systematic study on $\chi\bar{\chi}\to 2V\to 4{\rm SM}$ as a new heating
mechanism for NS in addition to $\chi\bar{\chi}\to 2{\rm SM}$ and kinetic
heating from DM-baryon scattering. The self-trapping due to $\chi V$
scattering is also considered. We assume the general framework that both
kinetic and mass mixing terms between $V$ and SM gauge bosons are present.
This allows both vector and axial-vector couplings between $V$ and SM fermions
even for $m_{V}\ll m_{Z}$. Notably, the contribution from axial-vector
coupling is not negligible when particles scatter relativistically. We point
out that the above approaches to DM-induced NS heating are not yet adopted in
recent analyses. Detectabilities of the aforementioned effects to the NS
surface temperature by the future telescopes are discussed as well.

##  I Introduction

It has been widely accepted that one-fifth of the total energy of the Universe
consists of dark matter (DM). Though multidisciplinary strategies are employed
to identify its essence, either from direct Aad:2015zva ; Abdallah:2015ter ;
Aalbers:2016jon ; Akerib:2016vxi ; Amole:2017dex ; Akerib:2017kat ;
Aprile:2017iyp ; Aprile:2018dbl ; Aprile:2019xxb ; Aprile:2019jmx  or indirect
detections Aartsen:2014oha ; Choi:2015ara ; Aartsen:2016zhm ; Aguilar:2015ctt
; TheFermi-LAT:2017vmf ; Ambrosi:2017wek , the nature of DM remains a puzzle.
The approach of using neutron star (NS) as the DM probe has been proposed from
the heating effect due to DM Kouvaris:2007ay ; deLavallaz:2010wp ;
Kouvaris:2010vv ; Baryakhtar:2017dbj ; Raj:2017wrv ; Chen:2018ohx ;
Bell:2018pkk ; Acevedo:2019agu ; Joglekar:2019vzy ; Keung:2020teb , the NS
instability caused by DM gravitational collapse Kouvaris:2010jy ; Leung:2011zz
; Kouvaris:2011gb ; McDermott:2011jp ; Guver:2012ba ; Bramante:2013hn ;
Bramante:2013nma ; Kouvaris:2013kra ; Gresham:2018rqo ; Grinstein:2018ptl ;
Garani:2018kkd ; Lin:2020zmm  and gravitation wave emitted from the merger of
binary NS admixed with DM Nelson:2018xtr ; Ellis:2018bkr ; Bauswein:2020kor .
Novel way of constraining long-lived particle through the NSs in the Milky Way
is also investigated recently Leane:2021ihh . In addition, DM self-interaction
naturally arises in various phenomenological models and was proposed to
resolve many issues in the small-scale structure, e.g. core-cusp, missing
satellite, too-big-to-fail and diverse galactic rotation curve, see Ref.
Tulin:2017ara  for a review. Current astrophysical observations constrain DM
self-interaction cross section $\sigma_{\chi\chi}$ in the range Randall:2007ph
; Walker:2011zu ; BoylanKolchin:2011de ; BoylanKolchin:2011dk ; Elbert:2014bma

$0.1\,{\rm cm}^{2}\,{\rm g}^{-1}\leq\sigma_{\chi\chi}/m_{\chi}\leq 10\,{\rm
cm}^{2}\,{\rm g}^{-1}$ (1)

where $m_{\chi}$ is the DM mass.

Without delving into details of model constructions, DM self-interaction can
be understood phenomenologically as an exchange of a vector boson $V$ or a
scalar boson $\phi$ in the dark sector (DS). Here $\phi$ is the dark Higgs
responsible for a spontaneously symmetry breaking of Abelian $U(1)_{X}$ in DS,
and therefore the generation of $V$ boson mass. Assuming DM $\chi$ is a
spin-$1/2$ fermion and charged under $U(1)_{X}$ gauge coupling $g_{d}$, DM
self-interaction is induced by $\mathcal{L}_{{\rm DM-
DM}}=g_{d}\bar{\chi}\gamma_{\mu}\chi V^{\mu}$ and constrained by Eq. (1).
Furthermore, $V$ can mix with SM photon and $Z$ boson through kinetic
Holdom:1985ag ; Galison:1983pa ; Foot:2004pa ; Feldman:2006wd ;
ArkaniHamed:2008qn ; Pospelov:2008jd  and mass mixing terms Babu:1997st ;
Davoudiasl:2012ag ; Davoudiasl:2013aya . These mixing terms appear in the
following Lagrangians:

$\displaystyle\mathcal{L}_{\rm gauge}$ $\displaystyle=$
$\displaystyle-\frac{1}{4}B_{\mu\nu}B^{\mu\nu}+\frac{1}{2}\frac{\varepsilon_{\gamma}}{\cos\theta_{W}}B_{\mu\nu}V^{\mu\nu}-\frac{1}{4}V_{\mu\nu}V^{\mu\nu},$
(2) $\displaystyle\mathcal{L}_{\rm mass}$ $\displaystyle=$
$\displaystyle\frac{1}{2}m_{Z}^{2}Z_{\mu}Z^{\mu}-\varepsilon_{Z}m_{Z}^{2}Z_{\mu}V^{\mu}+\frac{1}{2}m_{V}^{2}V_{\mu}V^{\mu},$
(3)

where $B^{\mu\nu}\equiv\partial_{\mu}B_{\nu}-\partial_{\nu}B_{\mu}$ is the
$U(1)_{Y}$ field strength in SM while $\varepsilon_{\gamma}$ and
$\varepsilon_{Z}$ are the kinetic and $V-Z$ mass mixing parameters
respectively. The electromagnetic (EM) and neutral-current (NC) interactions
between $V$ and SM fermions $f$ resulting from mixing terms in Eqs. (2) and
(3) are given by

$\mathcal{L}_{{\rm DS-SM}}=\left(\varepsilon_{\gamma}eJ^{\rm
EM}_{\mu}+\tilde{\varepsilon}_{Z}\frac{g_{2}}{\cos\theta_{W}}J^{\rm
NC}_{\mu}\right)V^{\mu}$ (4)

where $g_{2}$ is the $SU(2)_{L}$ coupling and $J_{\mu}^{{\rm EM}}$ and
$J_{\mu}^{{\rm NC}}$ are SM electromagnetic and neutral currents,
respectively. The coefficient $\tilde{\varepsilon}_{Z}$ is a linear
combination of two mixing parameters and it reduces to $\varepsilon_{Z}$ for
$m_{V}\ll m_{Z}$. Its general expression is given in Appendix A.

In this paper, we examine the effect of DM heating due to the above
phenomenological setup for a nearby 3 giga-year-old (Gyr-old) and isolated NS.
The associated temperature is around $100\,{\rm K}$ according to the standard
cooling mechanism if there is no other heating source. Therefore, any
temperature deviating from this benchmark value can be potentially due to DM
annihilation in the star. Besides, the dark boson $V$ can be produced via
$\chi\bar{\chi}\to 2V$ for $m_{V}<m_{\chi}$. Thus, $V$ decaying into a fermion
pair inside the star is possible. This implies that the heating from
$\chi\bar{\chi}\to 2V$ cannot be ignored if the decay length of $V$ is smaller
than the star’s radius. Searching the nearby old and cold NS can improve our
understanding on DM. The new dynamics emerging from the above phenomenological
setup will be discussed in the following sections. For completeness, we also
analyze the signal-to-noise ratio (SNR) in the James Webb Space Telescope
(JWST) Gardner:2006ky . Future telescopes such as European Extremely Large
Telescope (E-ELT) and Thirty-Meter Telescope (TMT) Skidmore:2015lga  will
constrain DM properties with unprecedented sensitivities. In the following
sections, we employ the NS mass $M_{0}=1.4M_{\odot}$ and and the radius
$R_{0}=12\,{\rm km}$. We also replace $g_{d}$ with
$\alpha_{\chi}=g_{d}^{2}/4\pi$ and all equations are expressed in terms of
natural units $\hbar=c=k_{B}=1$.

##  II DM capture and NS temperature

When a NS swipes through the space, the DM particles in the halo can scatter
with the baryons and leptons inside the star. Once DMs lost an appreciable
fraction of kinetic energies, they will be gravitationally captured by the NS.
This capture process has been investigated extensively with contributions from
neutrons, protons and leptons as well as relativistic correction included in
Refs. Bell:2020jou ; Bell:2020lmm . In this paper, only neutron contribution
to the capture rate $C_{c}$ is considered. Contributions from other particle
species are ignored due to their small yields. The total DM number $N_{\chi}$
in the star satisfies the differential equation

$\frac{dN_{\chi}}{dt}=C_{c}-C_{a}N_{\chi}^{2},$ (5)

where $C_{a}$ is the DM annihilation rate. Both coefficients $C_{c}$ and
$C_{a}$ are well studied and the expressions can be found in Refs.
Bell:2020jou ; Bell:2020lmm ; Chen:2018ohx  and references therein. We do not
reproduce here. Thus, the exact solution to Eq. (5) is obtained

$N_{\chi}=C_{c}\tau_{{\rm eq}}\tanh(t/\tau_{{\rm eq}})$ (6)

where $\tau_{{\rm eq}}=1/\sqrt{C_{c}C_{a}}$ is the equilibrium timescale. Once
$t>\tau_{{\rm eq}}$, $dN_{\chi}/dt=0$ and $N_{\chi}(t>\tau_{{\rm
eq}})=\sqrt{C_{c}/C_{a}}$ according to Eq. (5). The total annihilation rate at
this stage only depends on the capture rate since
$\Gamma_{a}=C_{a}N_{\chi}^{2}=C_{c}$. Note that $C_{c}$ depends on
$\sigma_{\chi n}$ and $\sigma_{\chi n}\leq\sigma_{\chi n}^{{\rm geom}}\approx
10^{-44}\,{\rm cm}^{2}$ where $\sigma_{\chi n}^{{\rm geom}}$ is the geometric
cross section. In principle, the maximum capture rate is determined by
$C_{c}(\sigma_{\chi n}^{{\rm geom}})=C_{c}^{{\rm geom}}$. Besides, when DM
falls into the NS surface, it is accelerated up to $0.3c-0.5c$. The non-
relativistic (NR) limit for calculating $\sigma_{\chi n}$ is not applicable.
Furthermore one has to consider contributions from axial-vector coupling due
to $V-Z$ mass mixing given by $\mathcal{L}_{\rm mass}$. We have thoroughly
included these effects. A brief discussion on how to compute $\sigma_{\chi n}$
in terms of relativistic kinematics is given in Appendix A.

NS is known to suffer from eternal cooling due to neutrino and photon
emissions. Without extra energy injection, the NS temperature drops until it
releases all its heat. However, if SM particles are produced due to DM
annihilation in the star, these particles can become a heat source and
potentially prevent the star from inevitable cooling. Therefore, the evolution
of NS interior temperature $T_{b}$ is governed by the equation

$\frac{dT_{b}}{dt}=\frac{-\epsilon_{\nu}-\epsilon_{\gamma}+\epsilon_{\chi}}{c_{V}},$
(7)

where $\epsilon_{\nu}\approx 2.1\times 10^{4}\,{\rm erg}\,{\rm cm}^{-3}\,{\rm
s}^{-1}\,(T_{b}/10^{7}\,{\rm K})^{8}$ is the neutrino emissivity,
$\epsilon_{\gamma}\approx 1.8\times 10^{14}\,{\rm erg}\,{\rm cm}^{-3}\,{\rm
s}^{-1}\,(T_{b}/10^{8}\,{\rm K})^{2.2}$ the photon emissivity,
$\epsilon_{\chi}$ the DM emissivity that is responsible for the heating from
DM annihilation and $c_{V}$ the NS heat capacity Kouvaris:2007ay .
Additionally, the surface temperature $T_{s}$ observed by a distant observer
is related to $T_{b}$ by $T_{s}\approx 8.7\times 10^{5}\,{\rm
K}\,(g_{s}/10^{14}\,{\rm cm}\,{\rm s}^{-1})^{1/4}(T_{b}/10^{8}\,{\rm
K})^{0.55}$ where $g_{s}=GM/R^{2}\approx 1.85\times 10^{14}\,{\rm cm\,s}^{-2}$
accounts for the redshift correction from the star’s surface gravity. It is
also pointed out that when $T_{b}<3700\,{\rm K}$, there is no distinction
between $T_{b}$ and $T_{s}$ Chen:2018ohx .

During each annihilation, a pair of DMs release $2m_{\chi}$ of energy in a
form of SM particles or dark bosons depending on which channel are
kinematically allowed. The total energy releasing rate by DM is
$\mathcal{E}_{\chi}=2m_{\chi}\Gamma_{a}\sum_{i}b_{i}$ where $b_{i}$ is the
branching ratio of a specific channel, eg. $e^{\pm},\mu^{\pm},\tau^{\pm}$ or
$q\bar{q}$, and $\sum_{i}b_{i}\leq 1$. Neutrino pair $\nu\bar{\nu}$ is also
part of the annihilation channel in the presence of $V-Z$ mass mixing, but it
cannot contribute to the heating. In addition to the annihilation, DMs also
lose kinetic energies $E_{k}$ to the star through the capturing process. This
has been realized as the kinetic heating Baryakhtar:2017dbj  with the rate
$\mathcal{K}_{\chi}=C_{c}E_{k}=C_{c}m_{\chi}(\gamma-1)$ where
$\gamma=1/\sqrt{1-v^{2}}$ is the Lorentz factor.111Even DM is not captured,
energy deposition still occurs as long as $\chi n$ scattering can happen. On
the other hand, the kinetic heating effect from such uncaptured DM is
relatively small and negligible in our calculation.  Thus, DM emissivity
$\epsilon_{\chi}$ is given by

$\epsilon_{\chi}=\frac{\mathcal{E}_{\chi}+\mathcal{K}_{\chi}}{V},$ (8)

where $V$ is the NS volume.

##  III Decays of dark boson

(a)

(b)

Figure 1: DM heating from dark boson production $\chi\bar{\chi}\to 2V$. Left:
$V$ decays into SM particles before it exits the star. Right: $V$ is self-
trapped due to multiple $\chi V$ scatterings and then decays.

Here we discuss the case of $V$ produced by DM annihilation. $V$ is usually
produced in DM rich environment. If $V$ can subsequently scatter off the
surrounding DMs multiple times, it could lose energy and be self-trapped. It
then decays promptly as shown in Fig. 1b. However, such self-trapping effect
is in general inefficient since $\chi V$ scattering length $\ell_{\chi V}$ is
much larger than the thermal radius $r_{\rm th}$. Hence the scattering rate is
much suppressed and irrelevant to the heating. We leave the detail discussions
in Appendix C. Another trapping is due to the scattering between $V$ and
neutrons. On the other hand the relevant cross section is further suppressed
by the factor $\tilde{\varepsilon}_{Z}^{4}$ and the scattering length is
expected to be much larger than the NS radius. It is safe to omit this effect
in our calculation as well.

However, $V$ can decay into other SM particles before it propagates to the
surface as long as the decay length $\ell_{{\rm dec}}$ is shorter than
$R_{0}$. See Fig. 1a. The decay length is given by $\ell_{{\rm
dec}}=v\gamma\tau_{{\rm dec}}$ with $v\equiv\sqrt{1-m_{V}^{2}/m_{\chi}^{2}}$
the velocity of $V$ and $\tau_{{\rm dec}}\equiv\Gamma_{V}^{-1}$ the lifetime
of $V$ at rest where $\Gamma_{V}$ is the total decay width. Since $V$ is
produced on shell, we do not consider $V$ decaying back to $\chi$ due to
$m_{V}<m_{\chi}$. The probability for $V$ to convert into SM particles after a
propagation distance $r$ is

$F=1-e^{-r/\ell_{{\rm dec}}}.$ (9)

We took $r=R_{0}$ in the calculation. However, if neutrino is the decay
product, it cannot be considered as the heating source and must be subtracted.
By examining the numerical results for $F$, we found that $V$ can decay before
it exits the star in most of our interested parameter space. This implies that
$\chi\bar{\chi}\to 2V$ also plays an important role in NS heating. See
Appendix C for details. Generally speaking, NS contains muons and electrons
that are in degenerate. To enable the decay $V\to e^{+}e^{-}$, $m_{V}$ not
only has to be heavier than $2m_{e}$, the final kinetic energy carried by
$e^{-}$ must also exceed the electron chemical potential to prevent from Pauli
blocking. This condition has been implemented in our study.

Given the information in this section, we summarize that even when
$\chi\bar{\chi}\to 2V$ dominates the annihilation channel for
$m_{V}<m_{\chi}$, the heating effect is still efficient due to $V$ decays.
However, the self-trapping is generally unimportant due to $\ell_{\chi V}\gg
r_{{\rm th}}$ in this paper.

##  IV Implication of DM on NS temperature

In this section, we describe how NS surface temperature $T_{s}$ is affected by
the DM annihilation. If $\epsilon_{\chi}$ is negligible, the standard cooling
mechanism gives $T_{s}\approx 100\,{\rm K}$ for a 3-Gyr-old NS. But when
$\epsilon_{\chi}$ is large enough to counterbalance $\epsilon_{\gamma,\nu}$,
$T_{s}$ could remain in a relatively higher temperature. We present the
numerical results of $T_{s}$ for both $\alpha_{\chi}=1$ and $0.01$ in Figs. 2
and 3 respectively. The adjacent DM density around NS is assumed to be the
same as that of the solar system, $\rho_{\chi}=0.3\,{\rm GeV/cm^{3}}$, since
we aims for the nearby isolated NS. The DM mass scale is shown from $100\,{\rm
MeV}$ to $10^{6}\,{\rm MeV}$. Once $m_{\chi}\lesssim 100\,{\rm MeV}$, all of
the annihilation channels to fermions will be Pauli blocked except neutrinos.
Nonetheless, there is no upper limit for DM mass in NS. But heavier $m_{\chi}$
results in lesser DM number density which makes the NS sensitivity worse. In
addition, Refs. Bramante:2017xlb ; Dasgupta:2019juq  pointed out when
$m_{\chi}\gtrsim\mathcal{O}(10-100)\,{\rm TeV}$, it requires multiple
scatterings to capture the DM and implies that the single-scattering capture
is inefficient. Thus, we restrict our discussions below the TeV DM where NS
has better sensitivity and can be complementary to current DM direct searches.
In the following, we discuss the general trends of the numerical results in
terms of $\alpha_{\chi}=1$, Fig. 2, unless specified otherwise. The
conclusions can be applied to $\alpha_{\chi}=0.01$ directly. A simple
understanding on $\alpha_{\chi}$ is that the dark sector interactions are
proportional to $\alpha_{\chi}^{2}$ and DM-SM interactions are to
$\alpha_{\chi}$. The derivations of such features on the scattering cross
sections for all interactions are given in the appendices.

The values for the parameter $\eta\equiv\varepsilon_{\gamma}/\varepsilon_{Z}$
from top to bottom are $1$ (combined,
$\varepsilon_{Z}=\varepsilon_{\gamma}\neq 0$), 0 (pure $V-Z$ mixing,
$\varepsilon_{\gamma}=0$) and $\infty$ (pure kinetic mixing,
$\varepsilon_{Z}=0$), respectively. From left to right, we have
$m_{V}/m_{\chi}=10$ (heavy mediator), $1$ (equal mass) and $0.1$ (light
mediator). $T_{s}$ is indicated by the color bar placed on the right and the
lowest temperature is $100\,{\rm K}$. Without annihilation, e.g. no anti-DM
exists, solely kinetic heating can raise $T_{s}$ up to $1750\,{\rm K}$. If DM
annihilation is included, $T_{s}$ can maximally reach to $3100\,{\rm K}$.

Various constraints are also plotted, including XENON1T Aprile:2018dbl , XENON
LDM (low mass DM) based on the ionization Aprile:2019xxb  and of Migdal
Aprile:2019jmx  effects, SIDM Randall:2007ph ; Walker:2011zu ;
BoylanKolchin:2011de ; BoylanKolchin:2011dk ; Elbert:2014bma , SN1987A
Sung:2019xie  and beam dump experiments Riordan:1987aw ; Bross:1989mp ;
Abdullah:2018ykz . The parameter curve rendering DM annihilation cross section
at the thermal relic value $6\times 10^{-26}\,{\rm cm}^{3}\,{\rm s}^{-1}$ in
the early Universe is plotted in green on each figure. We have adopted the
method given in Ref. Cirelli:2016rnw  for computing the Sommerfeld enhancement
factor. The DM relative velocity in the early Universe is taken to be $c/3$.
See Appendix B for details. Here we present the thermal relic cross section as
a reference point and refer the readers to Refs. ArkaniHamed:2008qn ;
Cassel:2009wt ; Lin:2011gj  for detailed discussions. In addition, although
the captured DMs can have relatively large Sommerfeld enhancement due to low
velocities,222Assuming DMs are thermalized with the NS core where
$T_{\chi}=T_{b}$. Thus the mean velocity is about $\sqrt{T_{\chi}/m_{\chi}}$.
the enhanced $\sigma v$ only shortens the equilibrium timescale $\tau_{\rm
eq}$. When $t\gg\tau_{\rm eq}$, the total annihilation rate only depends on
the capture rate with $\Gamma_{A}=C_{c}$. NS is generally insensitive to the
Sommerfeld enhancement as long as DMs are in equilibrium.

Figure 2: NS surface temperature $T_{s}$ in the $m_{\chi}-\varepsilon$ plane.
We took the age of NS is 3 Gyrs and the lowest $T_{s}=100\,{\rm K}$ without DM
heating. All figures have $\alpha_{\chi}=1$ and
$\eta=\varepsilon_{\gamma}/\varepsilon_{Z}$. From top to bottom, $\eta=1,0$,
and $\infty$. From left to right, $m_{V}/m_{\chi}=10,1,0.1$. Various
constraints from XENON1T Aprile:2018dbl , XENON LDM Aprile:2019xxb ;
Aprile:2019jmx , SIDM Randall:2007ph ; Walker:2011zu ; BoylanKolchin:2011de ;
BoylanKolchin:2011dk ; Elbert:2014bma , SN1987A Sung:2019xie ; Sung:2021swd ,
beam dump experiments Riordan:1987aw ; Bross:1989mp ; Abdullah:2018ykz  and
the parameter curve rendering thermal relic cross section are shown as well.

Figure 3: The same as Fig. 2 except $\alpha_{\chi}=0.01$.

###  IV.1 Case for $m_{V}\geq m_{\chi}$

When $m_{V}\geq m_{\chi}$, only $\chi\bar{\chi}\to f\bar{f}$ is allowed. A dip
occurs on each plot in Fig. 2 with this mass ordering. The resonant point is
caused by the pole in $\tilde{\varepsilon}_{Z}$ given by Eq. (16) when
$m_{V}=m_{Z}$ with $m_{Z}$ the SM $Z$ boson mass. In fact the value for
$\tilde{\varepsilon}_{Z}$ at this point is
$-i(\varepsilon_{Z}+\varepsilon_{\gamma}\tan\theta_{W})m_{Z}/\Gamma_{Z}$,
which is enhanced by the factor $m_{Z}/\Gamma_{Z}$.

Thus, the DM-neutron scattering cross section $\sigma_{\chi n}$ depends on
$\tilde{\varepsilon}_{Z}$ and is proportional to

$\sigma_{\chi
n}\propto\frac{\alpha_{\chi}\tilde{\varepsilon}_{Z}^{2}}{m_{V}^{4}}\frac{m_{\chi}^{2}m_{n}^{2}}{(m_{\chi}+m_{n})^{2}}\min(\xi,1)$
(10)

in the NR limit. See Eq. (23) for reference.333In the numerical calculation,
we used the general expression for $\sigma_{\chi n}$, Eq. (18), and the
derivation is given in the same appendix. Nonetheless, Eq. (10), or Eq. (23),
is simpler and suitable for our discussions in the main text. The last term
shows the suppression factor due to Pauli blocking where $\xi\sim q/\mu_{F}$
with $q$ the the momentum transfer during the scattering and $\mu_{F}$ the
neutron chemical potential.

In the equilibrium epoch, $t\gg\tau_{\rm eq}$, the total annihilation rate
$\Gamma_{A}=C_{c}\propto\sigma_{\chi n}$.444We found that the equilibrium
condition holds in most of the parameter space in this work. However, in the
calculation we adopted $\Gamma_{A}=C_{a}N_{\chi}^{2}$ with $N_{\chi}$ given by
Eq. (6), instead of simply assuming $\Gamma_{A}=C_{c}$. When
$\tilde{\varepsilon}_{Z}$ is at the resonant point, $\sigma_{\chi n}$ is
enhanced drastically by the factor $m_{Z}^{2}/\Gamma_{Z}^{2}$ so does the DM
heating resulted from DM emissivity $\epsilon_{\chi}$. This accounts for the
dip at $m_{V}=m_{Z}$ in each figure.

On the other hand, DM heating for $m_{\chi}$ in the sub-GeV region is much
stronger. It can be understood that, as $m_{\chi}\ll m_{n}$, $q\propto
m_{\chi}$ while $m_{V}/m_{\chi}$ is held fixed, we have $\sigma_{\chi
n}\propto\tilde{\varepsilon}_{Z}^{2}/m_{\chi}$ according to Eq. (10). Hence a
smaller $m_{\chi}$ leads to a larger $\sigma_{\chi n}$ as well as a more
effective DM heating. However, the effect of DM heating will not grow
indefinitely with $\tilde{\varepsilon}_{Z}$ as $\sigma_{\chi
n}\leq\sigma_{\chi n}^{\rm geom}\approx 10^{-44}\,{\rm cm}^{2}$. The maximum
$T_{s}$ caused by DM heating saturates when $\sigma_{\chi n}=\sigma_{\chi
n}^{\rm geom}$ and is around $3100\,{\rm K}$. This justifies our numerical
results in Fig. 2 that $T_{s}$ does not increase further when
$\tilde{\varepsilon}_{Z}$ is sufficiently large for a given $m_{\chi}$.

To all plots in Fig. 2, DM heating becomes weaker instead of proportional to
$1/m_{\chi}$ for $m_{\chi}\lesssim\mathcal{O}(170)\,{\rm MeV}$. Although DM is
capable of producing $e^{\pm}$ and $\mu^{\pm}$ in this mass range, the
chemical potentials for both particles are
$\mu_{F}^{e}\sim\mathcal{O}(170)\,{\rm MeV}$ and
$\mu_{F}^{\mu}\sim\mathcal{O}(70)\,{\rm MeV}$. All channels are Pauli blocked
and only pions formed by $q\bar{q}$ are allowed until $m_{\chi}<m_{\pi}$.
Nonetheless, in the presence of $V-Z$ mass mixing, neutrinos are also part of
the annihilation products and take a significant branching ratio in the DM
annihilation at such a mass region. But neutrino cannot contribute to the
heating. This explains why $T_{s}$ is much colder when
$m_{\chi}\lesssim\mathcal{O}(170)\,{\rm MeV}$. The DM heating in this region
is mainly due to kinetic heating. As $\sigma_{\chi n}=\sigma_{\chi n}^{\rm
geom}$, the resulted $T_{s}$ is around $1750\,{\rm K}$ from pure kinetic
heating.

Various $\eta$ values in Fig. 2 characterize the contributions from
$\varepsilon_{\gamma,Z}$ to $\tilde{\varepsilon}_{Z}$.555Since neutron is
charge neutral, $Q=0$, the effect of kinetic mixing $Q\varepsilon_{\gamma}$ in
Eq. (15a) has zero contribution to $\sigma_{\chi n}$. Hence $\sigma_{\chi
n}\propto\tilde{\varepsilon}_{Z}^{2}$. Nonetheless, if protons in NS are
considered, then $\varepsilon_{\gamma}$ shall contribute to the DM-proton
cross section $\sigma_{\chi p}$ as a consequence of non-vanishing
$Q\varepsilon_{\gamma}$. Both $\eta=1$ and $0$ are similar because even
$\varepsilon_{\gamma}\neq 0$, its effect to $\tilde{\varepsilon}_{Z}$ is
suppressed by $m_{V}^{2}/m_{Z}^{2}$ as seen from Eq. (16). For $\eta=1$, the
kinetic mixing can contribute comparably to the $V-Z$ mass mixing unless
$m_{V}>m_{Z}$. This can be clearly seen in Fig. 2 that the difference between
$\eta=1$ and $0$ is apparent only in $m_{V}>m_{Z}$ region, which is the region
to the right of the dip. To the left of the dip, the contribution from
$\varepsilon_{\gamma}$ to $\tilde{\varepsilon}_{Z}$ for $\eta=1$ is
negligible.

For $\eta=\infty$, $\varepsilon_{Z}$ vanishes so that the only contribution to
$\tilde{\varepsilon}_{Z}$ comes from $\varepsilon_{\gamma}$. As discussed
earlier, the effect of kinetic mixing term is suppressed by
$m_{V}^{2}/m_{Z}^{2}$ and thus $\sigma_{\chi
n}\propto\tilde{\varepsilon}_{Z}^{2}\propto\varepsilon_{\gamma}^{2}m_{V}^{4}/m_{Z}^{4}$.
The associated DM heating is in general much weaker than the cases with
$\eta=1$ and $0$. However, the advantage of $\eta=\infty$ is that no neutrino
can be produced in the DM annihilation due to the absence of $V-Z$ mass
mixing. The energy released from DM annihilation can be fully deposited into
NS. This accounts for the higher $T_{s}$ than $\eta=1$ and $0$ in terms of the
same $\sigma_{\chi n}$. But the difference is not apparent. Numerical
calculation shows it is around tens to $\mathcal{O}(100)\,{\rm K}$.

###  IV.2 Case for $m_{V}<m_{\chi}$

For light mediator case, the channel $\chi\bar{\chi}\to 2V$ dominates over
$\chi\bar{\chi}\to f\bar{f}$ due to $\alpha_{\chi}\gg\varepsilon_{\gamma,Z}$
in general. As long as $F\sim 1$, $V$ can fully decay into SM particles before
it exits the NS. The resulting heating from $\chi\bar{\chi}\to 2V$ with $V\to
f\bar{f}$ can be appreciable as shown in the rightmost panel of Fig. 2. The
heating region in the $m_{V}<m_{\chi}$ case is much more expanded than the
$m_{V}\gg m_{\chi}$ case since lighter $m_{V}$ induces larger $\sigma_{\chi
n}$ as shown in Eq. (10). The resulting effects from different $\eta$’s are
similar to those in the previous subsection. When DMs mainly annihilate to
$2V$, the thermal relic cross section is controlled by $\alpha_{\chi}$ and
$m_{\chi}$ while it is independent of $\varepsilon_{\gamma,Z}$. Hence the
thermal relic cross section only constrains $m_{\chi}$ when $\alpha_{\chi}$
and $m_{V}$ are fixed. For $\alpha_{\chi}=1$ and $0.01$, the $m_{\chi}$ values
rendering the thermal relic cross section are around $2\times 10^{8}\,{\rm
MeV}$ and $5\times 10^{5}\,{\rm MeV}$, respectively.

##  V Detectability of the future telescope

Figure 4: The exposure time $t_{{\rm exp}}$ for ${\rm SNR}=2$ in JWST. The
region enclosed by the red line indicates $t_{{\rm exp}}\leq 10^{5}\,{\rm s}$.

Since DM annihilation could significantly affect the NS surface temperature
$T_{s}$, we discuss the detectability of $T_{s}$ in the JWST and similar
telescopes in the future. The blackbody spectral flux density with $T_{s}$ at
a given frequency $\nu$ is given by Baryakhtar:2017dbj

$f_{\nu}(\nu,T_{s},d)=\frac{4\pi^{2}\nu^{3}}{e^{2\pi\nu/k_{B}T_{s}}-1}\left(\frac{R_{0}\gamma}{d}\right)^{2}$
(11)

where $d$ is the distance between the NS and the Earth. Taking
$T_{s}=2000\,{\rm K}$ and $d=10\,{\rm pc}$ as an example, we have
$f_{\nu}\approx 0.84\,{\rm nJy}$ at $\nu^{-1}=2\,{\rm\mu m}$. The signal-to-
noise ratio (SNR) for JWST-like telescope is proportional to
$f_{\nu}\sqrt{t_{{\rm exp}}}$ where $t_{{\rm exp}}$ is the exposure time. From
Ref. Gardner:2006ky , JWST covers $0.9\,{\rm\mu m}$ to $2.77\,{\rm\mu m}$
imaging sensitivity in its Near-Infrared Imager (NIRI). A F200W filter
centered at $\nu^{-1}=2\,{\rm\mu m}$ reaches ${\rm SNR}=10$ with
$f_{\nu}=10\,{\rm nJy}$ and $t_{{\rm exp}}=10^{4}\,{\rm s}$.

In Fig. 4, we plot the $t_{{\rm exp}}$ for obtaining ${\rm SNR}=2$ over
$d-T_{s}$ plane. The region enclosed by the red line represents $t_{{\rm
exp}}<10^{5}\,{\rm s}$. There are multiple filters available for NIRI with
$\nu^{-1}$ centered at various different values Gardner:2006ky . We select the
filter with $\nu^{-1}$ most suitably matching the corresponding blackbody
wavelength at $T_{s}$. This explains the zigzag behavior in Fig. 4. In
principle, as $\sigma_{\chi n}\sim\sigma_{\chi n}^{{\rm geom}}$, kinetic
heating can maximally warm the NS up to $1750\,{\rm K}$ without DM
annihilation. For NS that is located within 10 pc, JWST can achieve ${\rm
SNR=2}$ with $t_{{\rm exp}}\leq 10^{5}\,{\rm s}$ for $T_{s}\geq 1750\,{\rm K}$

##  VI Summary and outlook

In this work we have investigated the new dynamics arising from the kinetic
mixing and $V-Z$ mass mixing between the dark gauge boson $V$ of the broken
$U(1)_{X}$ symmetry and neutral gauge bosons in SM. In particular, $V-Z$ mass
mixing induces a resonance at $m_{V}\approx m_{Z}$, which can be seen from the
pole of $\tilde{\varepsilon}_{Z}$ at $m_{V}=m_{Z}$. The axial-vector part of
the coupling between $V$ and SM fermions has been included in our
calculations. As $\chi\bar{\chi}\to 2V$ dominates the annihilation channel for
$m_{V}<m_{\chi}$, $V$ can decay into a pair of SM fermions before it exits NS
and induce NS heating in addition to $\chi\bar{\chi}\to f\bar{f}$. Although
this contribution appears naturally in the dark boson model considered here,
it is usually not included in the model-independent analysis, such as the one
performed in Ref. Chen:2018ohx . We also demonstrated numerically that NS can
provide constraints on sub-GeV DM with feeble coupling to SM particles in
complementary to the current direct search. The detectability with reasonable
$t_{{\rm exp}}$ in JWST telescopes is discussed. Similar conclusion can be
drawn for the future JWST-like telescopes.

We note that this work only considers $\chi n$ scattering in the capture rate.
This explains why NS is not sensitive to the dark sector when
$\varepsilon_{Z}=0$ ($\eta=\infty$). Neutron interacts with DM only through NC
interaction governed by $\tilde{\varepsilon}_{Z}$. Once $\varepsilon_{Z}=0$,
NC interaction becomes much suppressed since $\varepsilon_{\gamma}$ in
$\tilde{\varepsilon}_{Z}$ is oppressed by $m_{V}^{2}/m_{Z}^{2}$. However, NS
also consists of protons although the fraction of them is rather small. When
protons are included, CC interaction will be involved for the capture of DM
and NS remains sensitive to the dark sector even for $\varepsilon_{Z}=0$. In
general, NS sensitivity will be improved by including proton contributions. We
leave this for future studies.

##  Appendix A DM-neutron interaction

Figure 5: Feynman diagram for DM-neutron scattering. The blob is an effective
vertex that includes both vector and axial-vector contributions from kinetic
mixing and $V-Z$ mass mixing.

When DMs fall into NS, they could scatter with neutrons via exchanging the
dark boson $V$ as shown in Fig. 5. The kinetic mixing and $V-Z$ mass mixing
generate vector and axial-vector interactions between $V$ and SM fermions. The
usual derivation of these interactions proceeds through the diagonalization of
both $\mathcal{L}_{\rm gauge}$ and $\mathcal{L}_{\rm mass}$ in Eqs. (2) and
(3), which gives rise to relations between fields in the gauge basis and those
in mass eigenstate basis. However, since we are only interested in
interactions up to $\mathcal{O}(\varepsilon_{\gamma})$ or
$\mathcal{O}(\varepsilon_{Z})$, we do not need to perform the diagonalization
but rather treating the mixing terms
$\varepsilon_{\gamma}B_{\mu\nu}V^{\mu\nu}/(2\cos\theta_{W})$ and
$\varepsilon_{Z}m_{Z}^{2}Z_{\mu}V^{\mu}$ as perturbations. These two mixing
terms generate the following two-point functions at the tree level

$\displaystyle i\Pi^{\mu\nu}_{V\gamma}$ $\displaystyle=$ $\displaystyle
i\varepsilon_{\gamma}k^{2}g^{\mu\nu},$ $\displaystyle i\Pi^{\mu\nu}_{VZ}$
$\displaystyle=$
$\displaystyle-i(\varepsilon_{\gamma}\tan\theta_{W}k^{2}+\varepsilon_{Z}m_{Z}^{2})g^{\mu\nu},$
(12)

where $k$ is the four-momentum of $V$ entering into kinetic mixing or $V-Z$
mixing vertex. Hence the EM coupling of $V$ to SM fermions results from
multiplying the two-point function $i\Pi^{\mu\nu}_{V\gamma}$, the photon
propagator $iD^{\gamma}_{\alpha\mu}(k)$, and the electromagnetic coupling
$ieA_{\alpha}J^{\alpha}_{\rm EM}$, as shown in Fig. 6. This multiplication
leads to

$\displaystyle ieJ^{\alpha}_{\rm
EM}\frac{-ig_{\alpha\mu}}{k^{2}}i\varepsilon_{\gamma}k^{2}g^{\mu\nu}V_{\nu}=ie\varepsilon_{\gamma}J^{\nu}_{\rm
EM}V_{\nu}.$ (13)

Similarly, NC coupling of $V$ to SM fermions is given by multiplying the two-
point function $i\Pi^{\mu\nu}_{VZ}$, the $Z$ boson propagator
$iD^{Z}_{\alpha\mu}(k)$, and the NC coupling $igZ_{\alpha}J^{\alpha}_{\rm
NC}/\cos\theta_{W}$. This gives rise to

$\displaystyle\frac{ig}{\cos\theta_{W}}J^{\alpha}_{\rm
NC}\frac{-i}{k^{2}-m_{Z}^{2}+im_{Z}\Gamma_{Z}}\left(g_{\alpha\mu}-\frac{k_{\alpha}k_{\mu}}{m_{Z}^{2}}\right)(-i)(\varepsilon_{\gamma}\tan\theta_{W}k^{2}+\varepsilon_{Z}m_{Z}^{2})g^{\mu\nu}V_{\nu}$
(14) $\displaystyle=$ $\displaystyle\frac{-ig}{\cos\theta_{W}}J^{\nu}_{\rm
NC}V_{\nu}\frac{(\varepsilon_{\gamma}\tan\theta_{W}m_{V}^{2}+\varepsilon_{Z}m_{Z}^{2})}{(m_{V}^{2}-m_{Z}^{2}+im_{Z}\Gamma_{Z})}.$

Here we have used the physical conditions $k^{2}=m_{V}^{2}$ and
$k_{\mu}\epsilon^{\mu}_{V}=0$. We have also chosen unitary gauge for the $Z$
boson propagator. Therefore, the interaction vertex between dark boson and
neutron in Fig. 5 has the following Lorentz structure
$ie\bar{\psi}_{n}\gamma^{\mu}(a_{f}+b_{f}\gamma^{5})\psi_{n}$ with

$\displaystyle a_{f}$ $\displaystyle=Q\varepsilon_{\gamma}+\frac{1}{\sin
2\theta_{W}}(I_{3}-2Q\sin^{2}\theta_{W})\tilde{\varepsilon}_{Z},$ (15a)
$\displaystyle b_{f}$ $\displaystyle=-\frac{I_{3}}{\sin
2\theta_{W}}\tilde{\varepsilon}_{Z},$ (15b)

where

$\tilde{\varepsilon}_{Z}=\frac{\varepsilon_{Z}+\varepsilon_{\gamma}\tan\theta_{W}(m_{V}^{2}/m_{Z}^{2})}{(1-m_{V}^{2}/m_{Z}^{2})^{2}+\Gamma^{2}_{Z}/m_{Z}^{2}}\left(1-\frac{m_{V}^{2}}{m_{Z}^{2}}-i\frac{\Gamma_{Z}}{m_{Z}}\right)$
(16)

and $\Gamma_{Z}$ is the $Z$ boson decay width, $Q$ and $I_{3}$ are the
electric charge and the weak isospin respectively. In Tab. 1, we list $Q$ and
$I_{3}$ for various particles. The values for neutron can be obtained by
summing the corresponding quantum numbers of three quarks $udd$ in the low
energy limit.

Figure 6: Feynman diagrams contributing to the coupling of dark boson $V$ to
SM fermions.

Mixing parameters $\varepsilon_{\gamma}$ and $\tilde{\varepsilon}_{Z}$ are
responsible for EM and NC interactions, respectively. EM interaction does not
contribute to $\sigma_{\chi n}$ since $Q=0$ for neutron. On the other hand
$\tilde{\varepsilon}_{Z}$ has a feeble dependence on $\varepsilon_{\gamma}$
with a suppression factor $m_{V}^{2}/m_{Z}^{2}$ when $m_{V}\ll m_{Z}$. This
explains why $\sigma_{\chi n}$ is still nonzero when $\varepsilon_{Z}=0$
($\eta=\infty$).

| $u$ | $d$ | $c$ | $s$ | $t$ | $b$ | $\ell$ | $\nu$
---|---|---|---|---|---|---|---|---
$Q$ | $\frac{2}{3}$ | $-\frac{1}{3}$ | $\frac{2}{3}$ | $-\frac{1}{3}$ | $\frac{2}{3}$ | $-\frac{1}{3}$ | $-1$ | $0$
$I_{3}$ | $\frac{1}{2}$ | $-\frac{1}{2}$ | $\frac{1}{2}$ | $-\frac{1}{2}$ | $\frac{1}{2}$ | $-\frac{1}{2}$ | $-\frac{1}{2}$ | $\frac{1}{2}$
Table 1: Values of $Q$ and $I_{3}$ for quarks, leptons and neutrinos.

The spin-averaged $\chi n$ scattering amplitude is given by

$\displaystyle\overline{|\mathcal{M}_{\chi n}|^{2}}$
$\displaystyle=\frac{8\pi\alpha_{\chi}}{(t-m_{V}^{2})^{2}}\\{-4m_{n}^{2}[(b_{f}^{2}-a_{f}^{2})m_{\chi}^{2}+a_{f}^{2}u+b_{f}^{2}(s+u)]+2(a_{f}^{2}+3b_{f}^{2})m_{\chi}^{4}$
$\displaystyle\quad-4a_{f}^{2}um_{n}^{2}+a_{f}^{2}(t^{2}+2tu+2u^{2})+2(a_{f}^{2}-b_{f}^{2})m_{n}^{4}+b_{f}^{2}(s^{2}+u^{2})\\},$
(17)

where $s$, $t$ and $u$ are the Mandelstam variables. DM scatters with neutrons
with its velocity boosted to $0.3c-0.6c$ by the NS gravity. It must be treated
relativistically. However, neutrons can be treated as at rest since its
chemical potential is $\mathcal{O}(200)\,{\rm MeV}$ in the star.

Therefore, from the method in Ref. Ilisie:2016jta , we are able to write down
the DM-neutron scattering cross section as

$\sigma_{\chi
n}=\frac{1}{16\pi\lambda^{1/2}(s,m_{1}^{2},m_{2}^{2})\lambda^{1/2}(s,m_{3}^{2},m_{4}^{2})}\int_{t_{-}}^{t_{+}}\overline{|\mathcal{M}_{\chi
n}|^{2}}dt$ (18)

where

$\lambda(x,y,z)=x^{2}+y^{2}+z^{2}-2xy-2yz+2xz$ (19)

is the Källén function,

$t_{\pm}=\frac{1}{2}\sum_{i=1}^{4}m_{i}^{2}-\frac{s}{2}-\frac{1}{2s}(m_{1}^{2}-m_{2}^{2})(m_{3}^{2}-m_{4}^{2})\pm\frac{\lambda^{1/2}(s,m_{1}^{2},m_{2}^{2})\lambda^{1/2}(s,m_{3}^{2},m_{4}^{2})}{2s},$
(20)

and

$s=m_{1}^{2}+m_{2}^{2}+2E_{1}m_{2},$ (21)

where $m_{1}=m_{3}=m_{\chi}$ and $m_{2}=m_{4}=m_{n}$ for the $\chi n$
scattering. In Eq. (21), the energy $E_{1}=\gamma m_{1}$ is the total energy
carried by particle 1, which is DM.

###  A.1 Pauli blocking in the $\chi n$ scattering

Note that if the momentum transfer $\sqrt{-t}$ in Eq. (18) is smaller than the
Fermi momentum, the suppression by Pauli blocking takes effect. We include
this in the numerical calculation by incorporating the method in Ref.
Bell:2020jou . Our result agrees with Ref. Bell:2020jou  in the three
benchmark scenarios that $\overline{|\mathcal{M}_{\chi n}|^{2}}$ are constant,
$t$-dependent and $t^{2}$-dependent.

###  A.2 Axial-vector contribution in the NR limit

If $\chi$ can be treated non-relativistically as well, we have
$s=m_{\chi}^{2}+m_{n}^{2}+2m_{\chi}m_{n}$,
$u=m_{\chi}^{2}+m_{n}^{2}-2m_{\chi}m_{n}$ and $t=0$. Therefore the amplitude
and the cross section become,

$\overline{|\mathcal{M}_{\chi n}^{{\rm NR}}|^{2}}=\frac{64\pi
a_{f}^{2}\alpha_{\chi}m_{\chi}^{2}m_{n}^{2}}{m_{V}^{4}}$ (22)

and

$\sigma_{\chi n}^{{\rm
NR}}=\frac{4a_{f}^{2}\alpha_{\chi}}{m_{V}^{4}}\frac{m_{\chi}^{2}m_{n}^{2}}{(m_{\chi}+m_{n})^{2}}$
(23)

which are independent of $b_{f}$ where it determines the strength of axial-
vector coupling.

##  Appendix B DM annihilation

(a) To SM particles

(b) To dark bosons

Figure 7: Various channels for DM annihilation

We can divide the DM annihilation into two categories, which are $m_{V}\geq
m_{\chi}$ and $m_{V}<m_{\chi}$, respectively. For the prior case, DM can only
annihilate into SM particles as shown in Fig. 7a. For the later one, as long
as $g_{d}\gg e\varepsilon_{\gamma,Z}$, the dominant annihilation products are
two dark boson $V$ as shown in Fig. 7b. The amplitude for $\chi\bar{\chi}\to
f\bar{f}$ is given by

$\displaystyle\overline{|\mathcal{M}_{\chi\bar{\chi}\to f\bar{f}}|^{2}}$
$\displaystyle=\frac{8\pi\alpha_{\chi}}{(s-m_{V}^{2})^{2}+m_{V}^{2}\Gamma_{V}^{2}}\\{a_{f}^{2}[-2(m_{f}^{2}+m_{\chi}^{2})(-m_{f}^{2}-m_{\chi}^{2}+2u)+s^{2}+2su+2u^{2}]$
$\displaystyle\quad+b_{f}^{2}[-4m_{\chi}^{2}(m_{f}^{2}+t+u)-2m_{f}^{4}+6m_{\chi}^{4}+t^{2}+u^{2}]\\}$
(24)

where $\Gamma_{V}$ is the $V$ decay width. Assuming DM is at rest in the star,
the amplitude can be simplified into,

$\overline{|\mathcal{M}_{\chi\bar{\chi}\to
f\bar{f}}|^{2}}=\frac{128\pi\alpha_{\chi}}{(4m_{\chi}^{2}-m_{V}^{2})^{2}+m_{V}^{2}\Gamma_{V}^{2}}m_{\chi}^{4}\left[a_{f}^{2}\left(1+\frac{1}{2}\frac{m_{f}^{2}}{m_{\chi}^{2}}\right)+b_{f}^{2}\left(1-\frac{m_{f}^{2}}{m_{\chi}^{2}}\right)\right].$
(25)

The partial decay widths of $V$ are given by

$\Gamma_{V}=\frac{m_{V}}{12\pi}\sqrt{1-4\frac{m_{f}^{2}}{m_{V}^{2}}}\left[a_{f}^{2}\left(1+2\frac{m_{f}^{2}}{m_{V}^{2}}\right)+b_{f}^{2}\left(1-4\frac{m_{f}^{2}}{m_{V}^{2}}\right)\right]$
(26)

for $V\to f\bar{f}$ and

$\Gamma_{V}=\frac{\alpha_{\chi}}{3}m_{V}\sqrt{1-4\frac{m_{\chi}^{2}}{m_{V}^{2}}}\left(1+2\frac{m_{\chi}^{2}}{m_{V}^{2}}\right)$
(27)

for $V\to\chi\bar{\chi}$. Note that we have omitted the Heaviside theta
function $\theta(m_{V}-2m_{\chi,f})$ in the above expressions but it is always
implemented when we perform the calculation to ensure the energy conservation.
Besides, when $m_{\chi}>m_{V}$, the channel $\chi\bar{\chi}\to 2V$ is allowed
and the amplitude is

$\displaystyle\overline{|\mathcal{M}_{\chi\bar{\chi}\to 2V}|^{2}}$
$\displaystyle=-\frac{32\pi^{2}\alpha_{\chi}{2}}{(t-m_{\chi}^{2})^{2}(u-m_{\chi}^{2})^{2}}\\{m_{V}^{4}[6m_{\chi}^{2}(t+u)-6m_{\chi}^{4}+t^{2}-8tu+u^{2}]$
$\displaystyle\quad+4m_{V}^{2}[m_{\chi}^{4}(t+u)-4m_{\chi}^{2}tu+tu(t+u)]-m_{\chi}^{4}(3t^{2}+14tu+3u^{2})$
$\displaystyle\quad+m_{\chi}^{2}(t^{3}+7t^{2}u+7tu^{2}+u^{3})+6m_{\chi}^{2}-tu(t^{2}+u^{2})\\}.$
(28)

In the NR limit,

$\overline{|\mathcal{M}_{\chi\bar{\chi}\to
2V}|^{2}}=256\pi^{2}\alpha_{\chi}^{2}\frac{m_{\chi}^{2}(m_{\chi}^{2}-m_{V}^{2})}{(m_{V}^{2}-2m_{\chi}^{2})^{2}}.$
(29)

Thus, the general expression for annihilation cross section is obtained by
using the Fermi golden rule,

$\sigma v=\frac{\overline{|\mathcal{M}|^{2}}}{32\pi
m_{\chi}^{2}}\sqrt{1-\frac{m_{f}^{2}}{m_{\chi}^{2}}}\theta(m_{\chi}-\mu_{F}^{f}),$
(30)

where $m_{f}$ is the final state particle mass and $\mu_{F}^{f}$ the chemical
potential of fermion $f$ in the star. There is no chemical potential for dark
boson $V$. Therefore, we arrive at

$(\sigma v)^{f\bar{f}}=4\alpha_{\chi}\kappa
m_{\chi}^{2}\sqrt{1-\frac{m_{f}^{2}}{m_{\chi}^{2}}}\theta(m_{\chi}-\mu_{F}^{f}),$
(31)

where

$\kappa=\frac{1}{(4m_{\chi}^{2}-m_{V}^{2})^{2}+m_{V}^{2}\Gamma_{V}^{2}}\left[a_{f}^{2}\left(1+\frac{1}{2}\frac{m_{f}^{2}}{m_{\chi}^{2}}\right)+b_{f}^{2}\left(1-\frac{m_{f}^{2}}{m_{\chi}^{2}}\right)\right]$
(32)

for $\chi\bar{\chi}\to f\bar{f}$ and

$(\sigma
v)^{2V}=8\pi\alpha_{\chi}^{2}\sqrt{1-\frac{m_{V}^{2}}{m_{\chi}^{2}}}\frac{(m_{\chi}^{2}-m_{V}^{2})}{(m_{V}^{2}-2m_{\chi}^{2})^{2}}$
(33)

for $\chi\bar{\chi}\to 2V$. The total annihilation cross section is the sum of
both

$\sigma v=(\sigma v)^{f\bar{f}}+(\sigma v)^{2V}.$ (34)

We note that the second term contributes when $m_{V}<m_{\chi}$.

##  Appendix C Dark boson in the star

Dark bosons can be produced from DM annihilation once $m_{V}<m_{\chi}$. This
channel is thought to have feeble effect on the heating since $V$ interacts
with the NS medium weakly and escapes without any trace. However, we found
that, depending on the strength of $\varepsilon_{\gamma,Z}$, $V$ can decay
into SM particles before it reaches the surface of the star. In the case that
the decay length $\ell_{{\rm dec}}$ is much smaller than the star’s radius,
the total energy released from the annihilation can be fully deposited to the
star. See Fig. 1a. We also examine the case that $V$ is produced in the DM
rich region in the star’s center. $V$ could undergo multiple scattering with
the surrounding DMs and self-trapped until it decays. See Fig. 1b. This is
another way to extract energy from $V$. We discuss both effects in the
following.

###  C.1 Decay length

Figure 8: Fraction $F$ of dark boson decay into SM particles that contributes
to the heating effect.

The dark boson decay length with time dilation effect is given by

$\ell_{{\rm dec}}=v\gamma\tau_{{\rm dec}},$ (35)

where $v=\sqrt{1-m_{V}^{2}/m_{\chi}^{2}}$ is the $V$ velocity and $\tau_{{\rm
dec}}=\Gamma_{V}^{-1}$ the $V$ lifetime at rest. Let us assume that $V$ is
produced in the center of the star and its propagation distance is $R_{0}$.
Fig. 8 presents $F$ defined in Eq. (9), i.e., the fraction of $V$ converting
into SM particles after traveling a distance $r=R_{0}$, as functions of
$\varepsilon_{Z,\gamma}$ and $m_{\chi}$ for $m_{V}=0.1m_{\chi}$. We have
subtracted neutrino contributions from $F$ since they cannot generate heat.
Since the branching ratio of $V$ decays to neutrinos is nonzero in the case of
$V-Z$ mass mixing, $F$ is generally smaller than $1$ for $\eta=1$. For
$\eta=\infty$, no neutrinos can be produced, thus $F$ can reach unity.

In these figures, the chemical potential for electron $\mu_{F}^{e}$ is about
$\mathcal{O}(170)\,{\rm MeV}$. For a dark boson at rest with
$m_{V}\leq\mu_{F}^{e}$, $V\to e^{+}e^{-}$ can be Pauli blocked even for
$m_{V}\geq 2m_{e}$. On the other hand, if $V$ is highly boosted as a result of
heavy DM annihilation, $V\to e^{+}e^{-}$ is not Pauli blocked as long as
$m_{\chi}\geq\mu_{F}^{e}$. Therefore, to enable $V\to f\bar{f}$ decays, two
conditions are required. The first is $m_{\chi}\geq\mu_{F}^{f}$ and the second
is $m_{V}\geq 2m_{f}$.

###  C.2 Dark boson-DM interaction length

Figure 9: $\chi V$ scattering via $s$ and $t$ channels.

Feynman diagrams contributing to $\chi V$ scattering are shown in Fig. 9 and
the amplitude is given by

$\displaystyle\overline{|\mathcal{M}_{\chi V}|^{2}}$
$\displaystyle=\frac{64\pi^{2}}{3}\frac{\alpha_{\chi}^{2}}{(s-m_{\chi}^{2})^{2}(t-m_{\chi}^{2})^{2}}\\{m_{V}^{4}[6m_{\chi}^{2}(s+t)-6m_{\chi}^{4}+s^{2}-8st+t^{2}]$
$\displaystyle\quad-
m_{\chi}^{4}(3s^{2}+14st+3t^{2})+m_{\chi}^{2}(s^{3}+7s^{2}t+7st^{2}+t^{3})$
$\displaystyle\quad+4m_{V}^{2}[m_{\chi}^{4}(s+t)-4stm_{\chi}^{2}+st(s+t)]+6m_{\chi}^{8}-st(s^{2}+t^{2})\\}.$
(36)

To compute the scattering cross section $\sigma_{\chi V}$, it is fair to
assume DM at rest. However, $V$ is produced with relativistic velocity since
$m_{\chi}>m_{V}$. We follow the procedure given in Eqs. (18)-(21) and set
$m_{1}=m_{3}=m_{V}$ and $m_{2}=m_{4}=m_{\chi}$. Thus,

$\sigma_{\chi
V}=\frac{1}{16\pi\lambda(s,m_{V}^{2},m_{\chi}^{2})}\int_{t_{-}}^{t_{+}}\overline{|\mathcal{M}_{\chi
V}|^{2}}dt.$ (37)

Note that $\chi V$ scattering is not subject to Pauli blocking since DMs do
not become degenerate in the presence of annihilation.

Figure 10: The ratio $\ell_{\chi V}/r_{{\rm th}}$ for $\eta=1$ and $\infty$.
We take $\alpha_{\chi}=1$ and $T_{\chi}=1000\,{\rm K}$ in the calculation.

The $\chi V$ scattering length $\ell_{\chi V}$ is given by

$\ell_{\chi V}=(n_{\chi}\sigma_{\chi V})^{-1},$ (38)

with $n_{\chi}\equiv N_{\chi}/V_{\chi}$ the average DM number density. The
volume characterizing DMs in NS is $V_{\chi}=4\pi r_{{\rm th}}^{3}/3$ where

$r_{{\rm th}}\approx 2.4\times 10^{3}\,{\rm
cm}\,\left(\frac{T_{\chi}}{10^{5}\,{\rm K}}\frac{10\,{\rm
MeV}}{m_{\chi}}\right)^{1/2}$ (39)

is the thermal radius. If $\ell_{\chi V}\ll r_{{\rm th}}$, $V$ can scatter
with surrounding DMs multiple times and gradually lose its kinetic energy.
However, our numerical result shows that $\ell_{\chi V}\gg r_{\rm th}$ in all
of our interested parameter space. In Fig. 10, we take $\alpha_{\chi}=1$ and
$T_{\chi}=1000\,{\rm K}$. The choice $\alpha_{\chi}<1$ makes $\ell_{\chi V}$
even longer due to a weaker $\chi V$ interaction. For $\eta=1$, the region for
$\ell_{\chi V}/r_{{\rm th}}<1$ happens when $m_{\chi}\lesssim 300\,{\rm MeV}.$
However, even $V$ can be self-trapped, it hardly decays into particles other
than neutrinos because the allowed channels, eg. $e^{\pm}$ and $\mu^{\pm}$,
are Pauli blocked. For $\eta=\infty$, only a very small parameter space leads
to $\ell_{\chi V}/r_{{\rm th}}<1$. Therefore, we conclude that the self-
trapping of $V$ is insignificant, hence only $V$ decays contribute to the
energy injection.

## References

  * (1) G. Aad et al. [ATLAS Collaboration], Eur. Phys. J. C 75, 299 (2015) [Erratum ibid 75, 408 (2015)] [arXiv:1502.01518 [hep-ex]].
  * (2) J. Abdallah et al., Phys. Dark Univ. 9-10, 8 (2015) [arXiv:1506.03116 [hep-ph]].
  * (3) J. Aalbers et al. [DARWIN Collaboration], JCAP 1611, 017 (2016) [arXiv:1606.07001 [astro-ph.IM]].
  * (4) D. S. Akerib et al. [LUX Collaboration], Phys. Rev. Lett. 118, 021303 (2017) [arXiv:1608.07648 [astro-ph.CO]].
  * (5) C. Amole et al. [PICO Collaboration], Phys. Rev. Lett. 118, 251301 (2017) [arXiv:1702.07666 [astro-ph.CO]].
  * (6) D. S. Akerib et al. [LUX Collaboration], Phys. Rev. Lett. 118, 251302 (2017) [arXiv:1705.03380 [astro-ph.CO]].
  * (7) E. Aprile et al. [XENON Collaboration], Phys. Rev. Lett. 119, 181301 (2017) [arXiv:1705.06655 [astro-ph.CO]].
  * (8) E. Aprile et al. [XENON Collaboration], Phys. Rev. Lett. 121, 111302 (2018) [arXiv:1805.12562 [astro-ph.CO]].
  * (9) E. Aprile et al. [XENON], Phys. Rev. Lett. 123, 241803 (2019) [arXiv:1907.12771 [hep-ex]].
  * (10) E. Aprile et al. [XENON], Phys. Rev. Lett. 123, 251801 (2019) [arXiv:1907.11485 [hep-ex]].
  * (11) M. G. Aartsen et al. [IceCube PINGU Collaboration], arXiv:1401.2046 [physics.ins-det].
  * (12) K. Choi et al. [Super-Kamiokande Collaboration], Phys. Rev. Lett. 114, 141301 (2015) [arXiv:1503.04858 [hep-ex]].
  * (13) M. G. Aartsen et al. [IceCube Collaboration], Eur. Phys. J. C 77, 146 (2017) [arXiv:1612.05949 [astro-ph.HE]].
  * (14) M. Aguilar et al. [AMS Collaboration], 211101 (2015).
  * (15) M. Ackermann et al. [Fermi-LAT Collaboration], Astrophys. J. 840, 43 (2017) [arXiv:1704.03910 [astro-ph.HE]].
  * (16) G. Ambrosi et al. [DAMPE Collaboration], Nature 552, 63 (2017) [arXiv:1711.10981 [astro-ph.HE]].
  * (17) C. Kouvaris, Phys. Rev. D 77, 023006 (2008) [arXiv:0708.2362 [astro-ph]].
  * (18) A. de Lavallaz and M. Fairbairn, Phys. Rev. D 81, 123521 (2010) [arXiv:1004.0629 [astro-ph.GA]].
  * (19) C. Kouvaris and P. Tinyakov, Phys. Rev. D 82, 063531 (2010) [arXiv:1004.0586 [astro-ph.GA]].
  * (20) M. Baryakhtar, J. Bramante, S. W. Li, T. Linden and N. Raj, Phys. Rev. Lett. 119, 131801 (2017) [arXiv:1704.01577 [hep-ph]].
  * (21) N. Raj, P. Tanedo and H. B. Yu, Phys. Rev. D 97, 043006 (2018) [arXiv:1707.09442 [hep-ph]].
  * (22) C. S. Chen and Y. H. Lin, JHEP 1808, 069 (2018) [arXiv:1804.03409 [hep-ph]].
  * (23) N. F. Bell, G. Busoni and S. Robles, JCAP 1809, 018 (2018) [arXiv:1807.02840 [hep-ph]].
  * (24) J. F. Acevedo, J. Bramante, R. K. Leane and N. Raj, arXiv:1911.06334 [hep-ph].
  * (25) A. Joglekar, N. Raj, P. Tanedo and H. B. Yu, arXiv:1911.13293 [hep-ph].
  * (26) W. Y. Keung, D. Marfatia and P. Y. Tseng, JHEP 07, 181 (2020) [arXiv:2001.09140 [hep-ph]].
  * (27) C. Kouvaris and P. Tinyakov, Phys. Rev. D 83, 083512 (2011) [arXiv:1012.2039 [astro-ph.HE]].
  * (28) S. C. Leung, M. C. Chu and L. M. Lin, Phys. Rev. D 84, 107301 (2011) [arXiv:1111.1787 [astro-ph.CO]].
  * (29) C. Kouvaris, Phys. Rev. Lett. 108, 191301 (2012) [arXiv:1111.4364 [astro-ph.CO]].
  * (30) S. D. McDermott, H. B. Yu and K. M. Zurek, Phys. Rev. D 85, 023519 (2012) [arXiv:1103.5472 [hep-ph]].
  * (31) T. Güver, A. E. Erkoca, M. Hall Reno and I. Sarcevic, JCAP 1405, 013 (2014) [arXiv:1201.2400 [hep-ph]].
  * (32) J. Bramante, K. Fukushima and J. Kumar, Phys. Rev. D 87, 055012 (2013) [arXiv:1301.0036 [hep-ph]].
  * (33) J. Bramante, K. Fukushima, J. Kumar and E. Stopnitzky, Phys. Rev. D 89, 015010 (2014) [arXiv:1310.3509 [hep-ph]].
  * (34) C. Kouvaris and P. Tinyakov, Phys. Rev. D 90, 043512 (2014) [arXiv:1312.3764 [astro-ph.SR]].
  * (35) M. I. Gresham and K. M. Zurek, Phys. Rev. D 99, 083008 (2019) [arXiv:1809.08254 [astro-ph.CO]].
  * (36) B. Grinstein, C. Kouvaris and N. G. Nielsen, Phys. Rev. Lett. 123, 091601 (2019) [arXiv:1811.06546 [hep-ph]].
  * (37) R. Garani, Y. Genolini and T. Hambye, JCAP 1905, 035 (2019) [arXiv:1812.08773 [hep-ph]].
  * (38) G. L. Lin and Y. H. Lin, JCAP 08, 022 (2020) [arXiv:2004.05312 [hep-ph]].
  * (39) A. Nelson, S. Reddy and D. Zhou, JCAP 07 (2019), 012 [arXiv:1803.03266 [hep-ph]].
  * (40) J. Ellis, G. Hütsi, K. Kannike, L. Marzola, M. Raidal and V. Vaskonen, Phys. Rev. D 97 (2018), 123007 [arXiv:1804.01418 [astro-ph.CO]].
  * (41) A. Bauswein, G. Guo, J. H. Lien, Y. H. Lin and M. R. Wu, [arXiv:2012.11908 [astro-ph.HE]].
  * (42) R. K. Leane, T. Linden, P. Mukhopadhyay and N. Toro, [arXiv:2101.12213 [astro-ph.HE]].
  * (43) S. Tulin and H. B. Yu, Phys. Rept. 730, 1 (2018) [arXiv:1705.02358 [hep-ph]].
  * (44) S. W. Randall, M. Markevitch, D. Clowe, A. H. Gonzalez and M. Bradac, Astrophys. J. 679, 1173 (2008) [arXiv:0704.0261 [astro-ph]].
  * (45) M. G. Walker and J. Penarrubia, Astrophys. J. 742, 20 (2011) [arXiv:1108.2404 [astro-ph.CO]].
  * (46) M. Boylan-Kolchin, J. S. Bullock and M. Kaplinghat, Mon. Not. Roy. Astron. Soc. 415, L40 (2011) [arXiv:1103.0007 [astro-ph.CO]].
  * (47) M. Boylan-Kolchin, J. S. Bullock and M. Kaplinghat, Mon. Not. Roy. Astron. Soc. 422, 1203 (2012) [arXiv:1111.2048 [astro-ph.CO]].
  * (48) O. D. Elbert, J. S. Bullock, S. Garrison-Kimmel, M. Rocha, J. Oñorbe and A. H. Peter, Mon. Not. Roy. Astron. Soc. 453, 29 (2015) [arXiv:1412.1477 [astro-ph.GA]].
  * (49) B. Holdom, Phys. Lett. 166B, 196 (1986)
  * (50) P. Galison and A. Manohar, Phys. Lett. 136B, 279 (1984)
  * (51) R. Foot, Int. J. Mod. Phys. D 13, 2161 (2004) [astro-ph/0407623].
  * (52) D. Feldman, B. Kors and P. Nath, Phys. Rev. D 75, 023503 (2007) [hep-ph/0610133].
  * (53) N. Arkani-Hamed, D. P. Finkbeiner, T. R. Slatyer and N. Weiner, Phys. Rev. D 79, 015014 (2009) [arXiv:0810.0713 [hep-ph]].
  * (54) M. Pospelov and A. Ritz, Phys. Lett. B 671, 391 (2009) [arXiv:0810.1502 [hep-ph]].
  * (55) K. S. Babu, C. F. Kolda and J. March-Russell, Phys. Rev. D 57, 6788 (1998) [hep-ph/9710441].
  * (56) H. Davoudiasl, H. S. Lee and W. J. Marciano, Phys. Rev. D 85, 115019 (2012) [arXiv:1203.2947 [hep-ph]].
  * (57) H. Davoudiasl, H. S. Lee, I. Lewis and W. J. Marciano, Phys. Rev. D 88, no. 1, 015022 (2013) [arXiv:1304.4935 [hep-ph]].
  * (58) J. P. Gardner, J. C. Mather, M. Clampin, R. Doyon, M. A. Greenhouse, H. B. Hammel, J. B. Hutchings, P. Jakobsen, S. J. Lilly and K. S. Long, et al. Space Sci. Rev. 123, 485 (2006) [arXiv:astro-ph/0606175 [astro-ph]] and JWST pocket guide
  * (59) W. Skidmore et al. [TMT International Science Development Teams & TMT Science Advisory Committee], Res. Astron. Astrophys. 15, 1945-2140 (2015) [arXiv:1505.01195 [astro-ph.IM]].
  * (60) N. F. Bell, G. Busoni, S. Robles and M. Virgato, JCAP 09, 028 (2020) [arXiv:2004.14888 [hep-ph]].
  * (61) N. F. Bell, G. Busoni, S. Robles and M. Virgato, [arXiv:2010.13257 [hep-ph]].
  * (62) J. Bramante, A. Delgado and A. Martin, Phys. Rev. D 96 (2017), 063002 [arXiv:1703.04043 [hep-ph]].
  * (63) B. Dasgupta, A. Gupta and A. Ray, JCAP 08 (2019), 018 [arXiv:1906.04204 [hep-ph]].
  * (64) E. M. Riordan, M. W. Krasny, K. Lang, P. De Barbaro, A. Bodek, S. Dasu, N. Varelas, X. Wang, R. G. Arnold and D. Benton, et al. Phys. Rev. Lett. 59, 755 (1987)
  * (65) A. Bross, M. Crisler, S. H. Pordes, J. Volk, S. Errede and J. Wrbanek, Phys. Rev. Lett. 67, 2942 (1991)
  * (66) M. Abdullah, J. B. Dent, B. Dutta, G. L. Kane, S. Liao and L. E. Strigari, Phys. Rev. D 98, 015005 (2018) [arXiv:1803.01224 [hep-ph]].
  * (67) A. Sung, H. Tu and M. R. Wu, Phys. Rev. D 99, 121305 (2019) [arXiv:1903.07923 [hep-ph]].
  * (68) A. Sung, G. Guo and M. R. Wu, [arXiv:2102.04601 [hep-ph]].
  * (69) M. Cirelli, P. Panci, K. Petraki, F. Sala and M. Taoso, JCAP 05 (2017), 036 [arXiv:1612.07295 [hep-ph]].
  * (70) S. Cassel, J. Phys. G 37 (2010), 105009 [arXiv:0903.5307 [hep-ph]].
  * (71) T. Lin, H. B. Yu and K. M. Zurek, Phys. Rev. D 85 (2012), 063503 [arXiv:1111.0293 [hep-ph]].
  * (72) V. Ilisie, Concepts in Quantum Field Theory, Springer (2016)