Heating neutron stars with GeV dark matter

An old neutron star (NS) may capture halo dark matter (DM) and get heated up by the deposited kinetic energy, thus behaving like a thermal DM detector with sensitivity to a wide range of DM masses and a variety of DM-quark interactions. Near future infrared telescopes will measure NS temperatures down to a few thousand Kelvin and probe NS heating by DM capture. We focus on GeV-mass Dirac fermion DM (which is beyond the reach of current DM direct detection experiments) in scenarios in which the DM capture rate can saturate the geometric limit. For concreteness, we study (1) a model that invokes dark decays of the neutron to explain the neutron lifetime anomaly, and (2) a framework of DM coupled to quarks through a vector current portal. In the neutron dark decay model a NS can have a substantial DM population, so that the DM capture rate can reach the geometric limit through DM self-interactions even if the DM-neutron scattering cross section is tiny. We find NS heating to have greater sensitivity than multi-pion signatures in large underground detectors for the neutron dark decay model, and sub-GeV gamma-ray signatures for the quark vector portal model.


Introduction
Dark matter (DM) may have a variety of interactions with SM particles and with DM itself, but with strengths that have evaded observation. A neutron star (NS) orbits through large fluxes of halo DM particles which may lose their energy via their interactions with the NS and become gravitationally bound to it. The high density and strong gravity of a NS may be able to compensate the feeble DM interactions and enhance the DM capture rate. The capture of halo DM by a NS had been extensively studied [1][2][3][4][5][6][7]. During the capture process, the strong gravitational potential of the NS accelerates the DM to more than half the speed of light, and DM-neutron scattering releases this kinetic energy to heat up the NS. Consequently, the NS temperature evolution will deviate from the standard cooling profile. A possible observable signal of DM capture by a NS is the detection of unexpectedly hot old neutron stars. The temperature of an old neutron star can be heated by ∼ 100 K to ∼ 2000 K, which is within the near-infrared band of the blackbody spectrum. The thermal emissions from nearby (within 100 pc), faint and isolated NS can be probed by upcoming infrared telescopes such as the James Webb Space Telescope, the Thirty Meter Telescope, and the European Extremely Large Telescope [3].
A DM-neutron cross section of ∼ 2 × 10 −45 cm 2 is large enough to heat up an old neutron star to ∼ 1000 K for DM masses between GeV and PeV. For DM lighter than a GeV, the capture rate is suppressed by Pauli blocking, while for DM heavier than a PeV, multiple scattering is necessary to slow down the halo DM particles. However, the total capture rate must lie below the geometric limit, which corresponds to all the ambient halo DM within the geometric area of the NS being captured.
We study scenarios with three aspects: (1) the DM is of GeV mass, which makes direct detection problematic (2) the DM is a Dirac fermion, so that it matters whether the particle or the antiparticle is the DM, and (3) the DM capture rate can reach the geometric limit. Specifically, we examine NS heating in the neutron dark decay model [8,9] and in a quark vector current portal framework [10,11].
The neutron dark decay model finds its origin in the recent neutron lifetime anomaly which is a ∼ 4σ discrepancy in the neutron lifetimes measured in beam and bottle experiments. If the neutron has the dark decay, n → χ + φ, where χ and φ are dark sector particles, with a partial width of about 7.1 × 10 −30 GeV the discrepancy is alleviated. The scalar φ is almost massless and no heavier than an MeV. The DM particle is very slightly lighter than the neutron and is a Dirac fermion to avoid constraints from neutronantineutron oscillations. Multi-pion signatures in neutron-antineutron oscillation searches by Super-Kamiokande only constrain the model if the DM isχ [12]. The model is interesting in that, as we will see, a NS can be composed of a substantial DM population, so that the DM capture rate can reach the geometric limit through DM self-interactions even if the DM-neutron scattering cross section is small.
As a second example, we consider dark matter that couples to u, d, s quarks through a dimensional-6 vector portal with independent couplings α u,d,s . These couplings can chosen so that the DM capture rate reaches the geometric limit. The NS also gets heated by the annihilation of GeV DM to light mesons (which can be described by chiral perturbation theory [10,11]).
The paper is organized as follows. In section 2, we review the process of DM capture by a NS, and the resultant NS temperature evolution is described in section 3. We study the neutron dark decay model in section 4, and the quark vector current portal model in section 5. We summarize our results in section 6.
2 Dark matter capture in neutron stars DM capture by a NS is primarily governed by DM-nucleon scattering and by DM selfinteractions if a significant DM population is bound by the NS. For weak scale DM, there are stringent upper limits on the DM-nucleon cross section, but constraints on DM selfinteractions are relatively loose. Interestingly, the preferred range for the self-interaction cross section to alleviate the core-cusp problem is 0.1 cm 2 /g σ χχ /m χ 1 cm 2 /g. This corresponds to σ χχ 10 −24 mχ 1 GeV cm 2 , which is much weaker than the upper limit σ χ−nucleon 10 −38 cm 2 from DM direct detection experiments. Therefore, DM selfinteractions may dramatically enhance the capture rate. Other processes, like DM-neutron annihilation, χχ annihilation and neutron decays to DM, also affect DM capture, and are included in our discussion below which is tailored for the neutron dark decay model; the corresponding equations for the quark vector current portal scenario are simpler and obtainable by straightforward modifications.
Because we study scenarios of Dirac fermion DM, the DM particle is either χ orχ. We consider the general case in which the NS is composed of both neutrons and χ, as is the case for the neutron dark decay model we consider. The evolution of the number of DM particles N DM in the neutron star is described by [13] where we distinguish the component N χ produced by neutron decay, n → χ + φ, from the halo DM component N DM because they may have different thermal properties. We assume that the rate of n → χ + φ is large enough to keep the neutrons and χ in thermal equilibrium. Halo DM-neutron elastic scattering contributes to the capture rate, and if DM isχ, halo DM also annihilates with neutrons, which under the assumption of a uniform mass distribution, are respectively given by [1] where the escape velocity of the NS is v esc (R) = 2GM/R 0.63 c,v is the DM dispersion velocity, and ρ DM is the local DM density; the relevant parameter values for the NS and the DM halo are listed in the table below. N n is the total number of neutrons in the NS, and We assume that the neutrons inside the NS behave as a Fermi gas and estimate the Fermi momentum to be p F (3π 2 ρ F /m n ) 1/3 = 437 MeV. DM-neutron scattering only occurs when the momentum exchange δp is larger than p F . We take this Pauli blocking into account by introducing a factor ξ = min(δp/p F , 1) in the above capture rate C s . Note that once the sum of cross sections (ξσ elastic χn for χ DM , or ξσ elastic χn + σ ann χn forχ DM) is larger than critical cross section, σ crit = πR 2 m n /M , and the sum of the capture rate and annihilation rates cannot be larger than the geometric limit, i.e., C c + C ann ≤ C c | geom . This is equivalent to N n (ξσ elastic χn ) ≤ πR 2 if DM is χ, and N n (ξσ elastic χn + σ ann χn ) ≤ πR 2 if DM isχ. For 1 GeV χ DM, the geometric limit C c | geom 8.2 × 10 32 yr −1 corresponds to σ crit 10 −45 cm 2 [5].
The DM capture rate due to scattering on χ from neutron conversion inside the NS or on the trapped DM (whose population is negligible in comparison) is [2] C χχ s = Cχχ s = 3 2 where we have again assumed that the mass density of the NS is uniform. Here, η = 3/2(v N /v), with v N the NS velocity relative to the Galactic center. For these cases, we define the geometric limits, N χ σ χχ→χχ ≤ πR 2 and N χ σχ χ→χχ ≤ πR 2 . The trapped DM with velocity v DM will form its own sphere of radius r DM (t), and the evolution of r DM (t) is derived as follows. The kinetic energy of each DM particle can be expressed in terms of the orbital radius r DM (t) as [2] with the rate of change in kinetic energy given by [1] where a χ is the fractional number of χ in the NS, and 1 − a χ is the fractional number of neutrons in the NS. The first (second) term in brackets corresponds to an energy release δE = 2m r E DM /(m n + m χ ) to the neutron component (χ component) of the NS [2], where m r is the reduced mass of the DM-neutron system. The energy gain, ∆E = 1 2 m χ (v 2 esc − v 2 DM ), results from a drop in the halo DM's potential energy from 1 2 m χ v 2 esc to 1 2 m χ v 2 DM after thermalizing with the trapped DM. Here, Effects of Pauli blocking are included by the factor, ξ = min( √ 2m r v DM /p F , 1). The evolution of r DM (t) is obtained by combining Eqs. (2.4) and (2.5), and and the temperature of the DM sphere T DM is given by 3 2 kT DM (t) = E DM . The last two terms in the second equation in Eq. (2.1) depends on the DM-neutron and DM-antiDM annihilation rates [13] Cχ n a σ ann which depletes the total number of trapped DM.

Temperature evolution
Soon after a NS is formed in a supernova explosion, its core has a temperature of about 10 11 K. It then cools down to 10 8 K through neutrino emission in about 10 5 years. When the core temperature falls below 10 8 K, photon emission dominates the cooling process. Unlike neutrino cooling, whose detailed mechanism is still under debate, photo cooling has less uncertainty, and we focus on this period of a neutron star's life. The interior temperature T int of a NS evolves according to [4] dT int dt where ν,γ,DM are the neutrino, photon and DM emissivities, and c V is the NS heat capacity per unit volume. Treating neutrons and the χ from neutron conversion as ideal Fermi gases, c V is given by [14,15] where the Fermi momenta are The neutrino emissivity is [14,15]  where n 0 = 0.16 fm −3 = 0.16 × 10 39 cm 3 , and n F is the average fermion number density in a NS. 1 Since neutrino emission depends on the eighth power of T int , neutrinos easily escape the NS when it is young. The surface temperature T sur of a NS is related to T int via [16][17][18] where g s = GM/R 2 = 1.85 × 10 14 cm s −2 is the gravitational acceleration at the surface of the NS. Including the effect of gravitational redshift, the observed temperature T obs at infinity is [19] T obs = T sur 1 − 2GM Rc 2 .
The NS luminosity L γ from the outer envelope is given by the Stefan-Boltzmann law: Then the effective photon emissivity is T int 3700 K .
Photon emission dominates the cooling process after 10 5 years, when T obs 10 6 K. Dark matter can inject energy into a NS in several ways. Halo DM-neutron elastic scattering and halo DM-neutron annihilation (if the DM isχ) contribute energy, is the angular average recoil energy transferred from the DM to a neutron in a single collision [19]. Here, For m χ m n we find E R 0.15m χ , which implies that annihilation is more efficient than elastic scattering at heating a NS if the halo DM-neutron annihilation and elastic scattering rates are comparable.
Another source of heat is the annihilation of trapped DM. If the trapped DM isχ, it can annihilate with χ from neutron conversion or with neutrons into SM particles and inject energy, is the efficiency with which energy is absorbed by the NS and depends on the annihilation final states. For instance, f DM = 0 for a purely neutrino final state, and f DM = 1 for a γγ final state. In principle, the contribution fromχ-neutron annihilation also has an efficiency factor, but we approximate this to unity for the final states we consider later; this also applies to the annihilation term in K DM above. The trapped DM also releases its energy via elastic scattering with neutrons and with χ from neutron conversion: From Eqs. (2.5) and (3.6), we see the path of energy conduction. The kinetic energy lost by halo DM to become trapped is transferred to the NS through scattering processes. Summing over the above three contributions, the total DM emissivity is The time evolution of the interior and observed temperatures of a NS without DM heating are shown in Fig. 1. For an old NS of age between 10 8 and 10 9 years, the temperature falls to about 500 K and 150 K, respectively.
In the rest of this section we do not consider the possibility of neutron conversion to χ and DM-neutron annihilation. Neutron star heating by DM capture can compensate the cooling from photon emission once T int falls to ∼ 1000 K. The NS can be heated by two processes: i) kinetic heating by the captured DM, and ii) DM annihilation into SM final states.
In the case of kinetic heating, if the capture rate is at the geometric limit, the observed (surface) temperature increases to 1480 (1660) K after the photon emission and DM kinetic heating processes attain equilibrium, Fig. 2 shows that T obs flattens out at 1480 K after 5 × 10 7 yrs.
DM annihilation consumes the entire DM mass to heat up the NS, and if the annihilation rate is high enough, photon emission and DM heating reach equilibrium earlier.
The observed (surface) temperature increases to 2480 (2780) K, when the photon emission energy-loss rate equals the sum of the DM kinetic and annihilation heating rate: ; see the right panel of Fig. 2. The surface temperature T sur saturates at 2780 K, when the DM annihilation rate equals the DM capture rate, i.e., N 2 DM C a | sat C c . Estimating N DM by multiplying C c = C c | geom with the typical age of an old NS, 5 × 10 8 yr, we find the saturating DM annihilation cross section to be In general, the value of v DM σ ann χχ | sat depends on C c and σ elastic DM−n . For example, consider a smaller capture rate, C c = 10 −4 × C c | geom . Without the heating from DM annihilation, the equilibrium condition, L γ | Tsur=170 K = C c E R , gives a final NS surface temperature T sur = 170 K. Including DM annihilation increases the surface temperature to T sur = 280 K using the criterion, L γ | Tsur=280 K = C c ( E R +m χ ). In this case, v DM σ ann χχ | sat 10 −35 cm 3 /s. In the neutron dark decay model, the trapped DMχ can annihilate with the neutron or χ from neutron conversion to provide additional heating. The observed (surface) temperature can reach 3100 (3440) K, if the photon emission energy-loss rate equals the sum of the DM kinetic and annihilation heating rates:

Neutron dark decay model
The defining feature of the neutron dark decay model is that the neutron decays to dark sector particles χ and φ. In the low energy limit, this can be described as a mixing between the neutron and the Dirac particle χ, which could serve as DM. However, since the DM particle is a Dirac fermion, either χ orχ could be DM, with different interactions with the neutron. Onlyχ can annihilate with the neutron, and only χ is produced from neutron conversion. We separately discuss the phenomenologies of NS heating for these two cases.

Model and NS equation of state
The interaction terms in the model are [8,9] where the heavy scalar Φ = (3, 1) −1/3 (color triplet, weak singlet, hypercharge -1/3) has mass above a TeV, and two Dirac fermionsχ and χ, and a scalar φ, are SM singlets The baryon number assignments for Φ,χ, χ, φ are −2/3, 1, 1, 0, respectively. The annihilation processχχ → φφ produces the observed DM relic abundance if the coupling λ φ 0.04. The first three interaction terms allow the decay n → χφ, which makes the NS unstable [20]. Including the Higgs portal and the g χχ χφ coupling, induces a repulsive χ-neutron interaction, which causes the energy density to increase when converting a neutron into χ, so that the neutron becomes stable inside a NS [9]. Then the interaction g nn nφ is generated from the Higgs portal interaction through the pion with where σ πn = 370 MeV and Higgs mass m h = 125 GeV. Constraints from rapid red giant star cooling [21] require |g n | 10 −14 . The sufficient condition to stabilize the NS is [9] z ≡ m φ |g χ g n | 71 MeV , which puts the NS in the neutron phase, and no χ is produced. Then the NS mass can reach two solar masses with central density of 6n 0 . For very light m φ , the choice, m φ 0.1 eV, g χ 4 × 10 −4 , and µ −0.4 eV, gives z 50 MeV to stabilize the NS, and also provides DM-self scattering cross sections of 0.1 cm 2 /g σ/m χ 1 cm 2 /g which alleviates the tension between N-body simulations of collisionless cold DM and large scale structure observations [9]. However, if m φ > 2.5 eV, g n = −10 −14 , and g χ 4π, z can easily exceed 71 MeV. Therefore, for heavier m φ , the NS is in a mixed phase, and we must solve the equation of state (EoS) equation to obtain the number densities, n χ and n n in the NS. In the mixed phase, the NS can be stabilized by introducing a repulsive DM-self interaction, and achieve a NS mass of about 2M .
We solve the EoS equation as follows. The energy density in a NS in a mixed phase is [9] ε(n n , n χ ) = ε nuc (n n ) + ε χ (n χ ) + n χ n n 2z 2 , where we assume χ is an ideal Fermi gas, and neutrons follow the EoS labeled V 3π + V R in Ref. [22], corresponding to moderately stiff EoSs that incorporate 3-nucleon forces and have been fit to the results of a quantum Monte Carlo. Then, with a (b) = 13.0 (3.21) MeV, α (β) = 0.49 (2.47) [23]. Here, z ≡ m φ /g χ comes from the DM-self interaction, which if mediated by a scalar or vector boson results in an attractive or repulsive force, respectively. A repulsive DM-self interaction can be realized by introducing an additional vector boson into the model; see Ref. [23] for details on the model construction. Here, we simply fix the ratio of z/z = |g χ |/|g n | 2 × 10 5 , although in general, z and z are two independent parameters. The equilibrium condition is  which is used to determine the n and χ compositions of the NS. The total Fermion number density satisfies n F = n n + n χ . The neutron phase is determined by the condition ∂ε/∂n χ | nχ=0 > 0, which requires that no χ be present, because introducing one χ increases the energy density. On the other hand, the condition ∂ε/∂n χ | nχ=n F < 0, transforms the entire NS into a χ star. The mixed phase is defined by ∂ε/∂n χ | 0<nχ<n F = 0. The three phases are shown in the left panel of Fig. 3 in the (z, n F /n 0 ) plane. The shading shows the density ratio a χ ≡ n χ /(n n + n χ ), which is almost independent of z for z > ∼ 0.25 GeV. The minimal composition of χ occurs for n F n 0 , in which case χ contributes about 30% of the total number density. The scenario with DM-self interactions is shown in Fig. 4. The lower panel corresponds to repulsive DM-self interactions which helps to stabilize the neutron star and extends the neutron phase up to z 10 3 GeV. We also solve the Tolman-Oppenheimer-Volkoff equation [24] to check that neutron stars heavier than 2M are obtainable. From the correlation between total pressure P = n 2 F d(ε/n F )/dn F and ε, we find the relations between the NS mass and radius in Fig. 5. From the left and middle panels we see that once z 100 MeV, the NS mass can be larger than 2M for the repulsive case. It is noteworthy that the NS in the repulsive case in the middle panel is in a mixed phase, and can still reach 2M .

DM-DM scattering cross section
The DM self-scattering cross section arises from the g χχ χφ and λ φχ χφ terms in the Lagrangian. The former is from the t-channel φ exchange diagram, while the later is generated from box diagrams withχ and φ in the loop. Since λ φ 0.04 is much larger than g χ 4×10 −4 , the loop-diagram contribution is comparable with the tree-level one. Since a large fraction of the NS could be composed of χ, DM self-capture is crucial for NS heating.
The DM self-scattering cross section due to the g χχ χφ term has been calculated in Ref. [25]. The velocity-dependent cross section, which is inversely related to the fourth  power of the velocity, was proposed to solve the core-cusp problem. During DM capture by a NS the typical DM velocity reaches v 0.63c, which suppresses this cross section to σ eff χχ→χχ 8.0 × 10 −40 cm 2 . Thus, the DM self-scattering cross section from g χχ χ becomes comparable to that from λ φχ χφ (via box diagrams), as we discuss below.
The DM-self scattering diagrams from the λ φχ χφ term, are shown in Fig. 6. The amplitudes for χ(p 1 )χ(p 2 ) → χ(p 3 )χ(p 4 ) from box-1 and box-2 of Fig. 6 are, respectively, where the relative minus sign arises from Fermi statistics. D µν , D µ,ν and D 0 are loop integration functions defined in LoopT ools [26] as , (4.9) , and µ is the renormalization scale. In order to match the Dirac spinors between box-1 and box-2, we use the Fierz transformation [27] where Dirac spinor w represents either the u or v spinors. Then the crossing operation, p 2 → −p 4 and p 4 → −p 2 , yields the amplitude for the DM self-scattering cross section χχ → χχ from box-3 and box-4. The box diagrams are significantly enhanced by the D µν loop function when the scattering angle in the centre of mass frame approaches θ cm 0 or π. This is due to the nearly massless mediator φ. Fortunately, neither collinear nor head on scattering contribute to the DM captured by DM inside the NS because the net trapped DM number remains the same in both cases. The energy transfer in a DM-DM collision is given by [19] (1 −B)m χ 2B + 2 √B (1 − cos θ cm ) , (4.11) whereB ≡ 1 − 2GM/(c 2 R) for a NS of mass M and radius R. So, the collinear scattering (θ cm 0) cannot slow down the incoming DM enough to be trapped by the NS. On the other hand, head on scattering θ cm π exchanges the momenta of the two initial DM particles such that the incoming DM particle gets trapped and the target particle gets kicked out of the NS. We define an effective DM self-scattering cross section, which is relevant to the DM captured inside the NS: 12) and similarly for χχ → χχ. The cross sections in Fig. 7 are finite. The loop-level contribution from λ φχ χ is comparable with the tree-level contribution from g χχ χφ because λ φ g χ .

Results
The salient feature of this model is that neutrons can convert into χ inside the NS, which makes the NS composed of n and χ in most of the interesting parameter space. Then the DM self-scattering cross sections from the box diagrams in Fig. 6, that are significantly larger than the critical cross section σ crit , enhance the capture rate above the geometric limit. Consequently, the NS can be heated up to 1500 K. If furtherχ − n andχ − χ annihilations are allowed, the NS temperature might reach 3100 K depending on the final state particles from annihilation. We are interested in the parameter regions which can explain the neutron lifetime anomaly. The masses m χ , m φ , and mχ in this model need to satisfy the relations [8] 937.992 MeV < m χ + m φ < 939.565 MeV ,  We choose three benchmark points of Ref. [12], within the region. We fix λ φ = 0.04 to give the correct DM relic density [8], and g χ = 4 × 10 −4 to alleviate the core-cusp problem [9]. For the neutron dark decay model, the DM can be eitherχ or χ, so we separately discuss these cases below.

χ is DM
In this subsection, we consider the case in which χ is DM, so there are no DM-neutron and DM-antiDM annihilation processes involved. Figure 8 shows the temperature increase in neutron stars older than 10 9 years in the parameter region of Eq. (4.13). The panels from left to right correspond to attractive DM self-interaction, no DM self-interaction, and repulsive DM self-interaction scenarios. For each panel, the higher temperature region corresponds to a mixed phase of NS, and the lower temperature region corresponds to the neutron phase. A dramatic temperature change occurs at the boundary of these two phases. For attractive DM self-interactions and no DM self-interactions, the boundary occurs for m φ 0.2 eV, which corresponds to z 100 MeV. For repulsive DM self-interactions, the phase transition gradually occurs for 10 eV m φ 100 eV, which corresponds to 5 GeV z 50 GeV.
In the neutron phase, DM capture relies primarily on DM-neutron scattering. We can see that the NS temperature is always below 200 K. Because the DM-neutron cross section is too small to saturate the geometric limit, the kinematic recoil energy of halo DM cannot heat up the NS. In the mixed phase, there are substantial χ from neutron conversion inside the NS, and so, the DM self-capture kicks in and dramatically enhances the halo DM capture rate to the geometric limit. This results in an observed NS temperature of 1580 K, when the equilibrium condition L γ | Tsur=1660 K = C c | geom ( E R ) is satisfied.

4.4.2χ is DM
In this subsection, we assumeχ is the DM candidate. Therefore, additional DM-neutron and DM-antiDM annihilation processes enhance the NS heating.
The DM-antiDM annihilation is through the χχ → φφ process. Whether or not χχ → φφ enhances the NS temperature, depends on whether or not the decay products of φ can be absorbed by the NS. If φ mixes with SM Higgs according to Ref. [9], φ → γγ is the dominant channel, so that NS heating can be further enhanced. For scenarios in which φ decays into neutrinos or dark sector particles, DM-antiDM annihilation does not contribute to the heating process. In the upper and lower rows of Fig. 9, we separately show the two scenarios in which the final state particles are absorbed or not absorbed by the NS.
In Fig. 9, for each panel, there are higher and lower temperature regions respectively corresponding to the mixed and neutron phases. The upper row of Fig. 9, which shows the neutron phase, has an additional DM-neutron annihilation process (compared to the χ DM case) to heat up the NS. However, its contribution is insignificant and the observed temperature is below 200 K. In the mixed phase, again, the substantial component of χ in the NS and large DM-antiDM scattering help the capture rate to reach the geometric limit, but the additional DM-χ annihilation cannot heat up the NS, because the annihilation final states cannot been absorbed. The result is that kinetic heating raises the NS temperature to 1580 K.
In the mixed phase, the DM-antiDM annihilation process enhances the NS observed temperature up to 3100 K corresponding to a surface temperature of 3440 K; see the lower panel of Fig. 9. This occurs when the equilibrium condition L γ | Tsur=3440 K = C c | geom ( E R + 2m χ ) is satisfied. But in the neutron phase, the temperature is lower than 200 K because there is no χ component from neutron conversion to annihilate with DMχ.

Quark vector current portal dark matter
We consider Dirac DM with mass around a GeV that couples to quarks through a vector current interaction. It is difficult for current DM direct detection experiments to probe this scenario because the recoil energy is much lower than the typical detector threshold. However, the leading DM annihilation final state is π + π − , which produces MeV photons that can be observed by near future instruments that will fill in the "MeV-gap" in the cosmic photon spectrum [10,11]. Through the quark vector current, we also expect substantial DM-neutron scattering that will enable a NS to capture halo DM, which in turn will heat the NS.

DM-nucleon scattering cross section
Consider a Dirac fermion DM particle χ that couples to quarks through a vector-vector current, where α q are the coupling strengths and Λ is a cutoff scale. To describe DM capture by a NS, the DM-neutron scattering cross section should be calculated in the relativistic limit, since the DM particles are accelerated close to the speed of light. The DM-neutron and DM-proton cross sections are given by [19,28] dσ χn,p (s, t) d cos θ cm = c χn,p Λ 4 where θ cm is the scattering angle in the center mass frame andμ ≡ m χ /m n m χ /m p . Here, c χp,n = (α u B p,n u + α d B p,n d ) 2 , with the integrated nuclear form-factors, B p u = B n d = 2 and B n u = B p d = 1. The nucleon form factor is |F n (E R )| 2 = exp[−E R /(0.114 GeV)] [28], where E R is the recoil energy in the initial n or p rest frame. For DM capture by a NS, in the initial nucleon rest frame, the energy of DM due to gravitational acceleration is , where we have neglected the thermal motion of the DM. The expressions for the other kinematic variables are where | − → p 0 | = √ s 2 λ 1/2 (1, m 2 χ /s, m 2 n /s) and λ(x, y, z) ≡ x 2 + y 2 + z 2 − 2xy − 2xz − 2yz. An example of the DM-neutron scattering cross section for DM capture by a NS is provided in Fig. 10. GeV, the DM-neutron cross section is larger than the critical cross section. Therefore, we expect the corresponding DM capture rate to reach the geometric limit and an old NS temperature can be heated up to 1500 K. The sensitivity provided by NS heating is significantly greater than that from future observations of MeV cosmic photons, which are sensitive to α q ∼ O(1) [19,28].

Chiral Lagrangian and DM annihilation
We now calculate NS heating due to DM-antiDM annihilation. At the GeV scale, DMquark vector current interactions can be described by Chiral perturbation theory, such that the DM annihilate into pseudoscalar or vector mesons. We focus on √ s 1.15 GeV, so that we only need to include the χχ → K + K − , K L K S , ρπ, ωπ channels.
The Feynman rules for GeV DM couplings to low-energy QCD pseudoscalar meson and vector meson can be found in appendix B of Ref. [10]. Then the vector meson propagator < 0|T (ρ µν , ρ αβ )|0 > is [29] g µα g νβ (m 2 , (5.4) and the polarization of ρ µν is [k µ ν (k) − k ν µ (k)] /m ρ . The polarization sum between ρ µν and ρ µ ν is given by k µ k ν g νµ + k ν k µ g µν − k µ k µ g νν + k ν k ν g µµ /m 2 ρ . Using the ρ propagator in Eq. (5.4) and the χχρ, K + K − ρ vertices from Appendix B of Ref. [10], the amplitude squared for and the values for f V , h p , and F can be found in Ref. [10] In terms of the Mandelstam variables, where θ is the angle between − → p and − → k , and | − In the threshold limit, s → 4m 2 χ ⇒ u − t = 0 and u + t = −2(m 2 χ − m 2 K ). Then the amplitude squared can be simplified to The total and partialχχ annihilation cross sections are shown in Fig. 10 including interference effects. For √ s > 1. 15 GeV, other channels are kinematically viable, like a glueball with neutral pions. Because the calculation of glueball emission is beyond the scope of this work, we only consider the DM annihilation cross section for √ s < ∼ 1.15 GeV. Moreover, as long as the DM annihilation rate is large enough to maintain equilibrium between the DM capture rate and depletion rates, including the new channels do not further increase the temperature of the NS. Without including the DM annihilation channels above √ s =1.15 GeV, we still obtain a conservative estimate of NS heating for DM masses above 0.575 GeV.

Results
In Fig. 11, we shown the observed temperature of the NS due to the vector-vector current couplings to quarks in Eq. (5.1). For α u or α d larger than O(10 −4 ), DM capture DM heats up the NS to more than 1480 K, which is shown by the black curve. However, for α u = −2α d the DM-neutron scattering cross section vanishes, and the NS does not get heated. This feature is indicated by the dashed line in Fig. 11.
In Fig. 12, we vary α s and α u + α d , and fix α u = α d , Λ = 100 GeV and m χ = 1 GeV. Clearly, T obs is insensitive to the parameter α s , which modifies vσ ann , but not σ elastic χn ; α s only affects O(100) K temperatures. We may understand the features of Fig. 12 as follows. First, focus on the region of T obs above 1000 K, where σ elastic χn is close to σ crit 2×10 −45 cm 2 and the DM capture rate C c reaches the geometric limit. This corresponds to α u α d 4 × 10 −5 which gives a DM annihilation cross section, vσ ann O(10 −33 ) cm 3 /s, which is six orders of magnitude larger than vσ ann | sat O(10 −39 ) cm 3 /s. α s only alters vσ ann within a similar magnitude, but cannot suppress it down to vσ ann | sat . Thus, for T obs around 1000 K, T obs is insensitive to α s . However, the situation is different when the final T obs is of O(100) K, which corresponds to much smaller values of C c and σ elastic χn . Take C c = 10 −4 × C c | geom as an example. This corresponds to σ elastic χn = 2 × 10 −49 cm 2 and α u α d 4 × 10 −7 , which gives vσ ann O(10 −37 ), which is much smaller than the saturating annihilation cross section, Figure 11. T obs (in K) in the vector portal DM framework by varying α u and α d . We fix α s = 0, Λ = 100 GeV and m χ = 1 GeV. vσ ann | sat O(10 −35 ). This means that increasing vσ ann by varying α s enhances the final NS temperature T obs . This behavior at O(100) K is evident from the dark blue region in Fig. 12. For α u + α d = 2 × 10 −7 , increasing |α s | from 10 −8 to 10 −5 raises T obs , which plateaus for |α s | > 10 −5 . The little spike at α s 3 × 10 −7 is due to destructive interference between the DM annihilation channels.

Summary
We have investigated NS heating by the capture of GeV-mass DM. We discussed the generic scenario that the NS could be in a mixed phase composed of both neutrons and a substantial population of DM from neutron conversion. In this case, the geometric limit of the DM capture rate can be saturated through DM self-interactions without DM-neutron interactions.
A NS can be in a mixed phase in the neutron dark decay model (that explains the neutron lifetime anomaly), because neutrons are able to convert to DM. We demonstrated that a NS in mixed phase can be stable and its mass can be as heavy as 2M by solving the equation of state and Tolman-Oppenheimer-Volkoff equations.
To illustrate the effect of DM capture on NS heating, we chose the above mentioned neutron dark decay model and the quark vector current portal framework. For the neutron dark decay model, since the DM self-scattering cross section is crucial to estimate the DM capture rate, we calculated the tree-level and one-loop box diagram contributions. In the mixed phase of a NS, DM self-scattering can enhance the DM capture rate up to the geometric limit without DM-neutron interactions. We find that for m φ 100 eV, the sensitivity of near future infrared instruments is greater than afforded by multi-pion signatures at Super-Kamiokande, Hyper-Kamionkande, and DUNE.
For quark vector portal DM, since the NS is in the neutron phase, halo DM is captured only via DM-neutron interactions. We find that the capture rate is close to geometric limit for α u,d O(10 −4 ), in which case the NS is heated to ∼ 1500 K. This is four order of magnitude more sensitive than the detection of MeV cosmic gamma-rays which is sensitive to α u,d O(1). We also find that NS heating is not sensitive to α s , unless future telescopes can observe NS temperatures of around 100 K.