Light(ly)-coupled Dark Matter in the keV Range: Freeze-In and Constraints

Dark matter produced from thermal freeze-out is typically restricted to have masses above roughly 1 MeV. However, if the couplings are small, the freeze-in mechanism allows for production of dark matter down to keV masses. We consider dark matter coupled to a dark photon that mixes with the photon and dark matter coupled to photons through an electric or magnetic dipole moment. We discuss contributions to the freeze-in production of such dark matter particles from standard model fermion-antifermion annihilation and plasmon decay. We also derive constraints on such dark matter from the cooling of red giant stars and horizontal branch stars, carefully evaluating the thermal processes as well as the bremsstrahlung process that dominates for masses above the plasma frequency. We find that the parameters needed to obtain the observed relic abundance from freeze-in are excluded below a few tens of keV, depending on the value of the dark gauge coupling constant for the dark photon portal model, and below a few keV, depending on the reheating temperature for dark matter with an electric or magnetic dipole moment. While laboratory probes are unlikely to probe these freeze-in scenarios in general, we show that for dark matter with an electric or magnetic dipole moment and for dark matter masses above the reheating temperature, the couplings needed for freeze-in to produce the observed relic abundance can be probed partially by upcoming direct-detection experiments.


Introduction
Dark matter with mass in the keV to GeV range has been receiving increased attention over the last few years. Numerous mechanisms exist for how such dark matter could have been produced in the early Universe. The mechanism of thermal freeze-out, which is perhaps the best studied mechanism, typically produces dark matter consistent with observations only between about ∼1 MeV to ∼100 TeV, being bounded below by bounds from Big Bang Nucleosynthesis (BBN) [1][2][3] and above by unitarity of the annihilation cross section [4]. However, besides producing dark matter below the GeV scale by traditional thermal freezeout [5], various other related and non-thermal production mechanisms exist, see e.g. .
In this paper, we consider several models for dark matter down to keV masses for which the relic abundance can be produced from freeze-in [9,10]. The couplings between the dark matter and Standard Model (SM) particles needed to obtain the observed relic abundance are typically very small, which naturally allows these models to avoid the BBN bound, since the dark matter particles will not be in chemical equilibrium with SM particles. We show that nevertheless, at least for sufficiently low masses, these models can be probed by constraints from the cooling of various stellar objects.
The models we consider are a dark matter particle coupled to a dark photon that mixes with the photon [29][30][31][32], dark matter with an electric dipole moment (EDM), and dark matter with a magnetic dipole moment (MDM) [33][34][35][36][37][38]. These models can naturally have small interactions with electrically charged SM particles through a small kinetic mixing parameter (for the dark photon portal) or through a higher dimension operator (for the EDM/MDM models). In these models, if the couplings are sufficiently small, the dark matter particles are never in thermal equilibrium with the SM particles, but are produced gradually from the SM thermal bath over time to produce the correct relic abundance. This is called the freeze-in mechanism. Depending on the type of interaction, the production may be dominant at low temperatures (IR freeze-in) or also occur at approved high temperatures (UV freeze-in) [39]. The dark photon model has IR freeze-in, in which case the results do not depend on the reheating temperature, while the models with an electric or magnetic dipole moment have UV freeze-in, where the reheating temperature matters. We calculate the freeze-in parameters for these models, including the contributions from the plasmon decay [40][41][42].
The couplings needed for freeze-in are typically so small that these models cannot be constrained from laboratory experiments. However, constraints from stellar objects, such as red giant stars (RG) and horizontal branch stars (HB), can probe these small couplings [41,[43][44][45][46][47][48][49][50][51]. In stellar objects, SM particles can collide and produce the hypothetical dark sector particles. These dark sector particles can then carry away energy and change the evolution histories of the stellar objects. Since the observed stellar properties are consistent with predictions from the standard stellar models, we can constrain the couplings of the dark sector to the SM particles. Very roughly, since the temperature of the RG and HB stars reach about 10 8 K, dark sector particle masses up to about 10 keV can be probed (more precisely, we will see that both the temperature and the plasma frequency inside the stars set the maximum mass that can be probed).
The remainder of the paper is organized as follows. In Sec. 2, we describe the salient features of the models considered in this paper and some basic constraints on them. In Sec. 3, we describe the freeze-in production in some detail. Sec. 4 discusses the constraints from the RG stars and HB stars. Sec. 5 briefly describes the prospects for probing these models in the laboratory. We present our conclusions in Sec 6. Three appendices provide additional details of our calculations.

Light Dark Matter Models Interacting or Mixed with Photons
In this work we focus on dark matter interacting with photons, either via kinetic mixing through a heavy dark photon or directly due to an electric or magnetic dipole moment.
For the dark photon (A ) portal, we will consider the dark matter candidate to be a fermion (χ) or a complex scalar (φ). The dark photon is the gauge boson of an additional broken U(1) gauge group, and it is kinetically mixed with the photon [29]. The Lagrangian is where is the kinetic mixing parameter, θ W the weak mixing angle, g D = √ 4πα D is the "dark" gauge coupling, and B µν and F µν are the field strength tensors of the hypercharge gauge boson and the dark photon, respectively. If the dark photon is massless, χ or φ is a millicharged particle, for which the stellar cooling constraints have already been discussed in the literature [46,50,52]. Moreover, even if the dark photon is massive but ultralight m m χ , the constraints are similar to the millicharged case [53,54]. We therefore focus here on the "heavy" dark photon case, where m ∼ O(m χ ). We mainly focus on the case with m > 2m χ or m > 2m φ , in which case the dark matter will consist of χ-or φ-particles, and we will only briefly comment on the case with m < 2m χ or m < 2m φ , in which case the dark matter can mostly consist of dark photons.
A model where the dark matter has an electric dipole moment (d χ ) or a magnetic dipole moment (µ χ ) is described by the following term in the Lagrangian

2)
respectively. The d χ and µ χ have mass dimension −1. These effective operators must come from an underlying theory at a larger scale. As an example, the dipole moment can be induced by heavy charged particles (a fermion and a scalar) that couple the dark matter to the SM through a loop [37]. In such a scenario, for charged particles of mass M , the electric dipole moment would be given by where e is the electron charge and g is the coupling between the heavy charged particles and χ. A similar equation holds for µ χ . The mass of the heavy charged particles M could be as light as ∼100 GeV [55] (possibly even slightly smaller [56]) to avoid collider bounds, which can be combined with a limit of g 2 < 4π requiring perturbativity to give an upper bound of d χ 0.5 TeV −1 . Note that this limit is stronger than the LEP limit on dark matter particles with electric or magnetic dipole moment directly, which is about d χ 4 TeV −1 [57] for m χ < 50 GeV. However, in this simple model of an additional scalar and fermion generating the dipole moment, the dark matter mass also receives loop corrections of roughly [37] Since the dark matter mass cannot (trivially) be smaller than its mass correction, we find an upper limit for the electric dipole moment of with a similar equation for the magnetic dipole moment. Of course, one could imagine different UV completions of the dark matter models with a dipole moment that do not have the same strong upper bound on the dipole moment.

General bounds on keV-to-GeV mass dark matter
Our main focus in this paper is on deriving the stellar constraints and freeze-in production of dark matter in the keV to GeV mass range. However, in the remainder of this section we briefly review other bounds on dark matter in or near this mass range.
If the dark matter is in chemical equilibrium with the SM bath in the early Universe, dark matter masses below ∼9.4 MeV (for a Dirac fermion) or ∼6.5 MeV (for a complex scalar) [2,3,58] are in tension with cosmological observables. The reason is that in the cosmological standard model, BBN started at a temperature of about 1 MeV, and the predictions are well confirmed by the measured abundance of light elements; extra relativistic degrees of freedom, N eff , during this evolutionary stage could affect the expansion of the Universe and thus change the temperature at which BBN begins, which would alter the predicted values for the abundance of light elements.
If the particles were never in chemical equilibrium with the SM (which is a necessary condition for freeze-in production), the number density is much smaller than the equilibrium number density, and the contribution to N eff is negligible. We will check this condition below when we compare the parameters needed for freeze-in production to the couplings that would keep the dark matter in chemical equilibrium with the SM bath.
Another constraint on the dark matter mass comes from the existence of small-scale structure. Below dark matter masses of 1 keV, fermionic dark matter cannot account for all of the dark matter in dwarf galaxies [59,60] due to the Pauli principle. Moreover, when the first structures form, the process can be disturbed if (a large component of) the dark matter (either fermionic or bosonic) is too warm. The fast streaming of the particles would then 'wash out' the forming structures. For thermal relics this is the case if the mass is lighter than about 1 keV [61]. In general, the constraint depends on the momentum distribution of the dark sector, which is modified in the non-thermal case by the average momentum |p| compared to the thermal momentum |p| eq [61] m χ/φ ≥ |p| |p| eq × 1 keV . (2.7) While a detailed calculation of the resulting mass bound in our models is beyond the scope of this paper, we do not expect the result to differ significantly from 1 keV (see e.g. [42], which considered a model consisting of dark matter interacting with a very light dark photon and find only an O(1) correction; see also [62]). For the dark-photon portal, large values of α D imply large self-interactions among the dark matter particles mediated by the dark photon. For dark photons with a mass near the dark matter mass, the self-interaction limit on α D from observations of the Bullet cluster, σ SIDM /m χ 2 cm 2 /g [63], become very stringent for small masses. In particular, for m = 3m χ , we need α D < 0.5 for dark matter masses below ∼20 MeV, while we need α D ≤ 10 −6 for dark matter masses below ∼28 keV. We will see below that the stellar constraints disfavor dark matter interacting with a dark photon to constitute a dominant component of dark matter from freeze-in production for dark matter masses below approximately 15 keV for α D = 10 −6 . Below, we will consider values for α D ranging between 10 −6 to 0.5. For dark matter interacting through a dipole moment, the self-interaction limits are not relevant, since the couplings to the mediator-the photon in this case-are very small.
We note that the bounds from N eff , structure formation, and self-interactions may be evaded if the dark matter candidates constitute only a sub-component of the observed dark matter density. Laboratory bounds will be discussed in Sec. 5.

Production via Freeze-In
If the interaction between the dark sector and SM sector is sufficiently small, the dark sector was never in chemical equilibrium with the SM sector throughout the history of the Universe. Excluded from the thermal bath, the dark matter abundance today therefore cannot be set via the typical freeze-out mechanism. Still, as long as there is some small coupling between the SM and dark sector, SM particles in the thermal bath can annihilate to produce dark sector particles. This is called the freeze-in mechanism [9,10], and here we consider the freeze-in production of dark matter interacting with a heavy dark photon mediator, an electric dipole moment, or a magnetic dipole moment, as discussed in Sec. 2.
In general, the number density, n a , of a particle species a produced within the thermal history of the Universe is derived from the Boltzmann equation Here, R(T ) is the number of interactions per unit volume and per unit time in which the particle is produced, where p and E are the momentum and energy, respectively, subscripts i and f correspond to initial and final particles, g is the spin degeneracy, f is the distribution function, + and − correspond to bosons and fermions, respectively, and the final state includes the particle a. R(T ) for 2 → 2 processes is conventionally written as n 2 i σv , and for 1 → 2 processes it is written as n i Γ . The yield Y from freeze-in is found by integrating [64] Here, M Pl 1.22 × 10 22 MeV is the Planck mass, g * (T ) is the effective number of relativistic degrees of freedom at temperature T , g * s (T ) the entropic relativistic degrees of freedom, and g(T ) = 1 + T 3 d ln(g * s ) dT (see e.g. [65]). The first factor of two in Eq. (3.4) accounts for the production of both particles and antiparticles.
For the production of light dark matter coupled to the SM via the dark photon, there are several processes: pair annihilation of the SM particles f withf , plasmon decay, and Z-boson decay, The contribution from Z-boson decay is important for m χ 1 GeV. For m χ 1 GeV, which is our focus in this paper, the Z-boson decay contributes O(10%). For the production of light dark matter coupled to the SM via an electric or magnetic dipole moment, there are pair annihilation of the SM particles f withf and plasmon decay processes that contribute to the production rate We show relevant processes for a dark matter fermion coupled to a dark photon or with an electric or magnetic dipole moment in Fig. 1 (the case of a scalar dark matter particle coupled to a dark photon is similar). We can estimate from dimensional analysis the temperature at which dark matter production via freeze-in is important. Since R(T ) has mass dimension 4, R(T ) ∼ T 4 for dimension-4 operators like the dark photon kinetic mixing term, and R(T ) ∼ T 6 Λ 2 for dimension-5 operators like the electric or magnetic dipole moment interactions. For the dark photon case, this then implies that dY dT ∼ T −2 , and freeze-in production is dominant at low temperature and is not sensitive to the reheating temperature. It is thus said to be infrared-dominated. In contrast, for the dipole moment case, dY dT ∼ T 0 , so processes at all temperatures are relevant, up to the reheating scale, which therefore determines the relic abundance. Figure 1: The dominant processes relevant for the production of dark matter interacting with a dark photon (upper diagrams) and dipole moments (lower diagrams) in the early Universe: pair annihilation (left) and plasmon/Z-boson decay (right).
For the infrared dominated dark photon case, the precise process that dominates the freeze-in production depends on the dark matter mass. If the dark matter is heavier than the electron mass, pair annihilation dominates over plasmon decays (both longitudinal and transverse plasmons). Transverse plasmon decays occur at higher temperatures than the annihilation process, as the plasma frequency has to fulfill 2m χ ω p ≈ T /10. The plasmon production process thus becomes inefficient for T 20m χ , while the pair annihilation process is still very efficient. For dark matter masses above the electron mass, the dominant contributions are expected from those SM particles that are lighter than the dark matter and freeze out only after the freeze-in production has been completed. For example, for m e m χ m µ , the dark matter freeze-in production is dominated by electron-positron pair annihilation. On the other hand, for dark matter masses below the electron mass, the pair annihilation process quickly becomes inefficient, and the decay of (transverse) plasmons yields sizable contributions to the number density. The longitudinal plasmon modes are always suppressed compared to the transverse plasmon modes for infrared dominated production due to lack of available phase space.
For the production of dark matter with an electric or magnetic dipole moment, which is determined by the reheating scale, plasmon decays are always sub-dominant. This is a general conclusion for UV-freeze-in through a high-dimensional operator.
We now consider the two production processes-pair annihilation and plasmon decayin more detail below. For the subdominant Z-boson decay contribution, we provide detailed formulae in Appendix B. Matching the dark matter density to the observed relic abundance today, Ω DM h 2 0.11 [66], will allow us to find the value of the kinetic mixing factor, , or the electric or magnetic dipole moment, d χ or µ χ , that gives the right relic abundance.

Pair Annihilation
For dark matter masses above ∼1 MeV, the freeze-in production is dominated by contributions from annihilation into the dark sector. The rate to produce the dark matter particle f in a thermal bath of temperature T can then be derived from [64,67] where m f is the mass of particle f , K n (x) is the modified Bessel function of the second kind, and s is the usual Mandelstam variable; g i = g j = 2 account for the degrees of freedom of the incoming fermions. Here |M 2 | is averaged over initial state spins and summed over final state spins. We can now consider several cases: • If the dark matter is a fermion interacting with a dark photon, the integrated matrix element of the annihilation of SM fermions f +f with charge q f in units of e and mass m f is given by [64] dΩ with the total width of the dark photon given by [64] 1 These formulae do not include A -Z-mixing, but we show the full formulae that include this mixing in Appendix B. In the second term of Eq. (3.9), the sum is over all fermions that are lighter than the dark photon. If m < 2m e , only the first term contributes. In our calculations, we drop the second term, which is a reasonable approximation as long as α D = g 2 D /4π α 2 . For Eq. (3.8), the major contribution comes from electronpositron annihilations, since the freeze-in is dominated by the lowest temperatures.
• If the dark matter is a boson interacting with a dark photon, the integrated matrix element is and the total width of the dark photon is Again, the second term accounts for the dark photon decays to SM particles, and the first term accounts for the decay to dark matter particles. A factor of 1/3 appears from averaging over the dark photon polarizations. We again drop the second term in our calculations.
• For the case of dark matter with an electric or magnetic dipole moment, the integrated matrix elements are respectively.
With these expressions, Eq. (3.7) can be used to derive the number density averaged interaction rate that enters Eq. (3.4) to yield the relic abundance.

Plasmon decay
Another important production process for dark matter comes from the decay of plasmons in the thermal plasma of the early Universe. In a thermal plasma, the interaction of the photon with charged particles, most dominantly electrons, leads to an effective mass for the photon, which depends on the electron density and temperature of the thermal bath (see Appendix A). At finite temperature, the photon propagator gets renormalized. The additional term acts like a self-energy of the photon, making it effectively massive. These quasi-massive states are called plasmons. The pole of the photon propagator determines the dispersion relations, which are then modified in comparison to the vacuum case. The properties of plasmons differ significantly from photons propagating in vacuum: they move slower than the speed of light and there is a longitudinal mode in addition to the transverse modes. We will refer to the different modes as 'longitudinal' or 'transverse' plasmons. The plasmon production and decay includes all electromagnetic processes where on-shell photons are produced. One thus has to be careful to avoid double counting diagrams that contribute with an on-shell photon. However, here there is no such danger for the annihilation process; the intermediate photon cannot be on-shell, which can be seen by cutting the diagram (see Fig. 1) in the middle. If the photon were on-shell, the inverse process of the right-handside would correspond to a plasmon that decays to an electron-positron pair. This process is kinematically not allowed for plasmons, since the charged SM particles receive corrections to their mass and are too heavy. Thus, in the early Universe, the decay of plasmons and the annihilation of SM particles are two distinct processes that can be treated separately.
The decay of plasmons is an important production mechanisms for neutrinos in stellar objects [40]. Similarly, the plasmons can decay into particles from a dark sector. This process is relevant in the early Universe as well as in stellar objects. In the following, we follow the notation and conventions of [40].
We now want to derive the production rates of dark matter from plasmon decay. In general, the thermally averaged rate is given by Here, K = (ω, k) is the four momentum of the plasmon and P DM , P DM are the four momenta of the outgoing dark matter particles. The sum over the polarizations 'pol' includes one longitudinal (L) and two transverse (T ) modes, thus g L = 1 and g T = 2. The distribution function of the plasmons is a Bose-Einstein distribution f (ω) = 1/(exp(ω/T ) − 1), and Z pol (i.e., Z L and Z T ) is given in Appendix A. To get Eq. (3.15), we follow the calculations in [51] for integrating over the outgoing dark matter phase space. The corresponding f DM (s), with s = ω 2 pol − k 2 , for each model that we consider is The two processes indicated by the expressions in Eq. (3.7) and Eq. (3.14) together then give the total production rate of the freeze-in process that appears in Eq. (3.4), which gives the final dark matter yield. It turns out that for dark matter interacting with a dark photon mediator, the decay of the longitudinal mode contributes only at the percent level compared to the contribution from the transverse mode. However, for the dipole moment models, the contribution from the longitudinal mode can dominate.

Dark photon + fermion dark matter
In Fig. 2 (left), we show the values of needed to obtain the correct relic abundance from freeze-in for the dark photon portal for a fermionic dark matter candidate. We consider six different scenarios. We see that the value of that yields the correct relic abundance through the freeze-in mechanism does not depend very strongly on the model parameters and mass ratios.  Figure 2: Solid lines show the values of for which the correct dark matter relic abundance is obtained in a model with fermionic (left) and scalar (right) dark matter from freeze-in through a dark photon. We show various choices for the dark-matter-to-dark-photon mass ratio and dark-photon couplings α D . For m < 2m χ/φ and m < 2m e , the dark photon can make up the relic abundance, rather than the fermion χ or scalar φ (see text for details). For parameters above these "freeze-in lines", too much dark matter is produced in the early Universe. Above the dashed lines, the dark matter and SM sector are in chemical equilibrium. Below the chemical equilibrium lines the model is safe from constraints on the number of relativistic degrees of freedom in the early Universe.
We first discuss the case for which the dark photon can decay to dark matter, and choose m = 3m χ . We show lines corresponding to three different values of α D , namely α D = 0.5, α D = α EM 1/137, and α D = 10 −6 . The latter value satisfies the self-interaction bound down to dark matter masses of ∼28 keV, below which production via freeze-in is constrained from stellar cooling, see Sec. 4.
From the matrix element for the annihilation Eq. (3.8) or the plasmon decay Eq. (3.17), we see that far off the resonance, where the decay width is negligible, the production rate is proportional to α D 2 . In contrast, close to the resonance, the coupling α D , which appears also in the decay width, divides out. In the early Universe, a range of temperatures is scanned, so the production always gets large contributions from the resonance at some point. It is for this reason that the scaling of between the different lines is less than a factor of 1/ √ α D .
If the dark matter is heavier than electrons, it is mainly produced through the annihilation process, which rapidly becomes inefficient when the temperature drops below the mass of a dark matter pair. The production from plasmon decays is subdominant, since the plasma frequency, which is the measure for the available phase space of the decay, is much lower than the temperature, ω p ≈ T /10. So the plasmon process contributes less than the annihilation process, since the freeze-in production of dark-photon-mediated dark matter particles is infrared dominated and the plasmon decay process stops at larger temperatures than the annihilation process.
In contrast, for dark matter masses below the electron mass, the dominant contribution comes from plasmon decays. Again, thanks to the infrared domination, the main production happens at late times and stops when the production of dark matter is kinematically forbidden. If this occurs after electrons freeze out, the annihilation process does not contribute anymore. However, plasmon decays can still occur until the plasma frequency ω p ≈ T /10 falls below twice the dark matter mass.
We also present two cases where the dark photon is (partly) lighter than twice the dark matter mass. In this case, the dark photons can make up the entire dark matter abundance for m < 2m e . This scenario is studied in [41], and we find that it does not receive significant corrections from the presence of heavier particles in the dark sector. However, we briefly discuss the composition of the dark relics in these scenarios. The small kink in the orange line (m = 100 keV) in Fig. 2 marks the transition when the dark fermions become heavier than the dark photons. While the dark photons are efficiently produced in the early Universe for all dark fermion masses, they rapidly decay to the dark fermions for dark fermion masses below the kink. On the other hand, for dark fermion masses above the kink, the dark photons are very long-lived and constitute the dark relics. In the same figure, the brown line shows the case where the dark matter is always lighter than the dark photon (m = m χ /10). At dark photon masses below ∼100 keV, only a small number of dark photons is produced directly, such that their initial relic abundance is small. Thus, most of the dark relics that are frozen-in are dark fermions, which may then annihilate into dark photons. However, above roughly m ∼ 100 keV and below m = 2m e , the direct freeze-in production of dark photons becomes sizable compared to dark fermions, and they can constitute the dark relics. For m > 2m e , the dark photon is short-lived and decays rapidly into electron positron pairs. This causes the big kink in the brown line in Fig. 2.
Irrespective of the component that is most dominantly produced by the freeze-in mechanism (dark fermions or dark photons), processes like dark photon annihilation into dark fermions or vice versa could allow for a change in the relative abundance of the two species after freeze-in. This will only occur if the dark gauge coupling is large enough and the two species achieve chemical equilibrium in the dark sector. This has important implications for direct and indirect dark matter searches. This is reminiscent of the "leak-in" scenario discussed in [27,68]. We leave further investigation of this effect to future work.
Finally, we show the values of needed for the ultralight dark photon mediator case [12,13,42,53]. We see that these values are similar to the other cases, especially (at low dark matter masses) the case with m = m χ /10.
When is smaller than the values indicated by the solid curves in Fig. 2, the relic abundance of χ-particles and/or dark photons is less than the observed dark matter relic abundance, so that they form a subdominant dark matter component. Above the freeze-in line, too much dark matter would have been produced through freeze-in, and the Universe would be overclosed. We also show the lines indicating the coupling at which the dark matter would attain chemical equilibrium with the SM bath. This happens if the interaction rate Γ = n eq χ σv for the annihilation process is at any time bigger than the Hubble rate For the interaction rate, the annihilation into electrons is the most relevant process, i.e., the one that gives the most stringent constraint, and the rate is given by Eq. (3.7) divided by the equilibrium number density [67] n eq i = Here, g i is the number of degrees of freedom of the particle species; for a Dirac fermion or a complex scalar, g i = 4. Chemical equilibrium would of course spoil the freeze-in mechanism, and it is thus important to check that the freeze-in line does not get too close to the coupling that allows for chemical equilibration. In Fig. 2 (left), we see that for the dark photon portal dark matter this requirement is fulfilled.

Dark photon + scalar dark matter
The values needed to obtain the correct relic abundance from freeze-in for scalar dark matter coupled to a dark photon are shown with solid lines in Fig. 2 (right). We also show with dashed lines the values above which chemical equilibrium with the SM is reached. We find that the results for scalar dark matter look very similar to the fermion dark matter scenario. Thus, below we will usually consider only fermion dark matter, but we emphasize that our results are approximately applicable to scalar dark matter as well.

Dark matter with an Electric or Magnetic Dipole Moment
The values of the electric (magnetic) dipole moment needed to obtain the correct relic abundance from freeze-in for electric (magnetic) dipole dark matter are shown for different reheating temperatures in the left (right) plot of Fig. 3. Above the lines, too much dark matter is produced in the early Universe, whereas below the lines the relic abundance is lower than the observed amount of dark matter. We find that the dark matter production via freeze-in is dominated for all m χ by the pair-annihilation process rather than from plasmon decay. 2 This is because production via plasmon decays is subdominant compared to production via pair annihilation at high temperatures, while dark matter particles of electron and dipole moments are produced via UV freeze-in. Dark matter with a mass above the reheating temperature cannot be produced efficiently in the early Universe and thus large values of the electric or magnetic dipole moments are needed. However, these values become so large for increasing dark matter masses that the freeze-in line intersects the chemical equilibrium line, so that for even larger dark matter masses, the observed relic abundance cannot be obtained from freeze-in for any value of the dipole moment. electric (magnetic) dipole moment needed to obtain the correct relic abundance from freeze-in for electric (magnetic) dipole dark matter for different reheating temperatures. For parameters above these "freezein lines", too much dark matter is produced in the early Universe. Above the dashed lines, the dark matter and SM sector are in chemical equilibrium, so that this parameter region is not compatible with any form of freeze-in production. Below the chemical equilibrium lines the model is safe from constraints on the number of relativistic degrees of freedom in the early Universe.

Stellar Constraints
For dark matter masses below ∼100 keV, the models under consideration can be constrained from stellar cooling arguments [43]. In some stages of stellar evolution the energy loss into a dark sector is severely constrained from astrophysical observations of globular clusters. We will briefly review these arguments and then discuss the resulting limits.

Critical Stages of Stellar Evolution
A globular cluster is a star cluster with a particularly high density of stars. It is tightly gravitationally bound and has a spherical shape. It is a satellite of a galaxy and most likely it was formed within the star formation process of the parent galaxy. Low mass stars in globular clusters can be used to constrain particle physics properties by the stars' characteristic properties of helium ignition and burning.
To understand the origin of the constraints, it is useful to look at the Hertzsprung-Russeldiagram (HRD) of the globular cluster. Each single star is presented by a point in the plane spanned by the absolute brightness (in magnitudes) and the spectral classes. Note that the latter is correlated with the surface temperature, increasing from the right to the left. Within its lifetime, a star moves through the diagram, starting on the so-called 'main sequence', which is the diagonal from the lower right to the upper left corner. In this stage, it is burning hydrogen to helium in the core. Once the core is transformed into helium, the fusion process moves outwards, building a shell around the core. At this stage, the star becomes a red giant (RG), increasing its magnitude and moving to the red giant branch in the HRD. Depending on the mass of the star, it continues burning helium and eventually heavier elements in the horizontal branch or as a super giant. In its final stage, it becomes either a white dwarf, which is found in the lower left corner of the HRD (faint and hot), or a neutron star or black hole, which are not depicted, since they have no brightness.

Helium Ignition in Red Giants
For low mass stars (0.5M M 2.3M ) helium ignition starts once the core has accumulated to roughly 0.5M . At this stage of stellar evolution the star has reached the tip of the red giant branch. An extra source of cooling would delay helium ignition. The resulting heavier core would imply longer hydrogen burning in the shell and thus a brighter red giant. From the magnitude of the tip of the red giant branch one can thus constrain unknown elementary processes that would enhance the cooling of the star.
Simulations have shown that an extra energy loss of 10 erg g −1 s −1 [43] is consistent with observations. The core density of the red giant is on average 2 × 10 5 g/cm 3 and varies only within a factor of order one. The electrons are degenerate. The temperature is 10 8 K and the electron concentration is Y e = 0.5.

Lifetime on the Horizontal Branch
Once the star is burning helium it moves to the so-called 'horizontal branch' (HB) in the HRD. The stars have a core mass of roughly 0.5M . The stars on the HB differ only in the mass of their hydrogen shell, and hence they have different surface temperatures (or spectral classes) but a similar magnitude; this is why in the HRD they lie on a horizontal line. A globular cluster has hundreds of thousands of stars. This allows one to determine the lifetime of a star on the HB from the ratio of the number of stars on the HB to the number of stars on the red giant branch. It agrees with the prediction from the stellar standard models within 10%. However, an exotic contribution to cooling would result in a faster fuel consumption and has been constrained to be smaller than 10 erg g −1 s −1 [43]. The temperature and electron concentration are the same as for the red giants discussed in Sec. 4.1.1 but the density is slightly smaller, 0.6 × 10 4 g/cm 3 . In this case, the electrons are not degenerate.

Dark Matter Production Mechanisms In Stars
The energy-loss rate or luminosity per unit volume in stellar objects can be written as Here, we use the same notation as in Eq. (3.2), and E out is the energy sum of outgoing dark sector particles. The electron-photon plasma inside stars gives rise to several production channels of light dark matter particles. The core temperature of the stellar objects discussed above reaches up to ∼10 keV. In this energy regime, the dominant production channel for light dark matter particles is through plasmon decays. Additionally, for dark matter masses above the plasma frequency where the plasmon decay is kinematically suppressed, the production via bremsstrahlung processes is relevant. In the following, we will discuss the plasmon decay and the bremsstrahlung processes. The total luminosity will then be given by We note that these processes for the dark matter models with an electric and a magnetic dipole moment have been calculated in detail recently in [51], which appeared during the latter stages of completing our work. We have checked that the production via Compton scattering is subdominant to the production from bremsstrahlung and plasmon decays, so we do not consider it further. We refer the reader to [51] for the detailed calculations also for the production via Compton scattering.

Plasmon Decays
As discussed in Sec , allow decays to massive particles like neutrinos or dark matter. An estimate of the maximum dark matter mass that can be produced from plasmon decays is given by the plasma frequency ω p (see Eq. (A.1)), which reaches roughly 1.6 keV in stars on the horizontal branch and 8.6 keV in red giants before helium ignition (for comparison, it is 0.3 keV in the solar core). For plasmon decays, Eq. (4.1) can be written as [50] dL plasmon The form of this formula is simple to understand: the total energy carried away in the dark sector is given by the energy of the decaying plasmon ω T,L and the rate with which it decays, Γ T,L . The factor of two in the numerator of the first term accounts for the two transverse modes, and the denominator comes from the Bose-Einstein-distribution that is obeyed by the photons in the star of temperature T . The plasmon decay rate to the dark sector is given by Eq. (3.16) The difference between the rate derived here compared to the freeze-in rate Eq. (3.14) is only the factor of the energy ω T /L in Eq. (4.3). To derive the stellar constraints, we are interested in the energy that is taken away from the star, i.e., the energy of the decaying plasmon. For the freeze-in production in the early Universe, the relevant quantity is the number of particles that go into the dark sector.

Bremsstrahlung processes
Dark matter particles can be produced by bremsstrahlung processes during the electronproton collisions, see Fig. 4. For dark matter masses above the plasma frequency where the plasmon decay is forbidden, this process dominates. For the calculation of the bremsstrahlung processe, we use the soft radiation approximation (SRA). Then, the energy loss rate Eq. (4.1) for the bremsstrahlung can be separated into the collision part and the radiation part [50]. Although the electron-proton collision is through a t-channel photon, which can have arbitrarily small momentum, the SRA is also valid because of the effective photon mass. Compared to the exact calculations in [51], the SRA gives the equivalent results up to an O(1) factor. For the process e(P 1 ) + p(P 2 ) → e(P 3 ) + e(P 4 ) + χ(P χ ) +χ(Pχ), the rate is found to be dL Brems where the index k corresponds to the radiated off-shell photon with 4-momentum K = (ω, k), s χ = ω 2 − k 2 , and for the solid (dashed) lines. The cooling constraints are derived for stars on the horizontal branch (brown) and red giants (red). In green, we show the parameters for which freeze-in production provides the entire dark matter relic abundance (see also Fig. 2); above the line too much dark matter would have been produced. In blue, we show the parameters for which thermal freeze-out production provides the entire dark matter relic abundance. Above the cyan lines, the dark sector was in chemical equilibrium with the SM bath and is constrained below m χ = 9.4 MeV by N eff . Below ∼1 keV dark matter is constrained from structure formation. Other relevant constraints and some projections from terrestrial searches are shown in Fig. 7. The bounds on scalar dark matter coupling to a dark photon (not shown) are similar.
is the electron-proton elastic scattering amplitude. See Appendix. C for detailed calculations.
Here, we use the approximation that the self energy of the photon is ω 2 p . Note we have used a plus sign in the s-channel propagator to avoid the double counting from the resonance. To get the second line in Eq. (4.4), we again follow the calculations in [51] for integrating over the outgoing dark matter phase space. The factor f DM (s) for each case is shown in Eqs. (3.17)-(3.19).

Results
The total dark luminosity is found by integrating Eq. (4.3) and Eq. (4.4) over the volume of the star and summing the two contributions. If the density and temperature profile is non-trivial, like in the solar case, all quantities depend on the radius, which has to be taken into account for the spatial integration. The constraint on is found by requiring that the dark luminosity does not exceed the limits discussed in Secs. 4.1.
In Fig. 5, we show the stellar constraints for the dark photon portal dark matter for m = 3m χ and for α D = 0.5 (solid lines) and α D = 10 −6 (dashed lines). The brown and red contours show the constraints from the lifetime on the horizontal branch ('HB') and the nondelay of helium ignition in red giants ('RG'), respectively. The plasma frequency in red giant stars is the highest, hence it can probe the largest dark matter masses. When the plasma frequency equals the dark photon mass, the propagator in the cross section is on resonance. The production is enhanced at this parameter point, and the constraint is thus particularly strong. This is seen in the spike-like features in the stellar constraints. While we show the results for a dark matter fermion only, we again note that the bounds on scalar dark matter coupled to a dark photon are very similar.
For the dark-photon-mediated dark matter, we compare in Fig. 5 the stellar constraints to the freeze-in lines (in green), which are also shown in Fig. 2. We find that dark matter that is entirely produced from this mechanism is ruled out below ∼20 keV for α D = 0.5, and below ∼40 keV for α D = 10 −6 . Note that the areas between the respective freeze-in and freeze-out lines (blue) are forbidden in this model, as an overabundance of dark matter would have been produced, overclosing the Universe. Additional decay modes (of the dark photon) beyond the ones assumed in the minimal model setup discussed here, or slight model variations, could open up some of this parameter region (see, e.g., [15][16][17][18][19][20][21][22][23][24][25][26][27]).
Note that we have not derived the constraints from the cooling of white dwarfs and the sun. In the high-mass regime where the bremsstrahlung processes dominate, it is not competitive with the other stellar cooling constraints as the white dwarfs and the sun have a much lower temperature than red giant stars. However, due to the high density the white dwarfs have a high plasma frequency, of ∼23 keV. Thus, a small fraction of the parameter space on the right-hand-side of the red giant tip can in principle be excluded additionally (see e.g. [46], where this was shown for dark photon dark matter).
The stellar constraints for dark matter with a dipole moment are shown in Fig. 6, together with the freeze-in lines. The left (right) plot shows the limits for dark matter with an electric (magnetic) dipole moment, respectively. The limits from the red giants are always stronger than the ones from the horizontal branch stars. No resonant production occurs due to the absence of a mediator in that mass range. We find that for dark matter with a dipole moment, the freeze-in mechanism is constrained for m χ 2 − 5 keV depending on the reheating temperature. We show also the freeze-out parameters from [35], as well as the LEP limit from [57].
The stellar constraints do not have upper boundaries, unlike the supernova 1987A constraints [43,50]. Consider first the supernova 1987A constraint. In this case, the lower boundary of the constrained region is set by the requirement of producing a sufficient number of dark sector particles to carry away more energy than that carried away by neutrinos, which are believed to dominate the energy loss. This would drastically change the cooling of the proto-neutron star, in conflict with observations. The upper boundary of the constrained re-  Figure 6: Stellar cooling constraints derived in this work on dark matter with an electric dipole moment (left) or a magnetic dipole moment (right), from stars on the horizontal branch (brown) and red giants (red). We also show lines for different reheating temperatures along which freeze-in production provides the entire dark matter relic abundance (see also Fig. 3); above the line too much dark matter would have been produced. Above the cyan lines, the dark sector was in chemical equilibrium with the SM bath and is constrained below m χ = 9.4 MeV by N eff . Below ∼ 1 keV dark matter is constrained from structure formation. Above the gray line the models are constrained from LEP data. The blue curve shows the parameters needed to obtain the correct relic abundance from thermal freeze-out.
gion arises from a sufficient number of dark sector particles becoming trapped, thermalizing with the matter inside the proto-neutron star, and failing to carry away sufficient energy. Roughly speaking, if the dark sector particles couple more strongly than neutrinos to the matter inside the proto-neutron star (mostly protons, neutrons, and electrons), the dark sector particles are unable to carry away enough energy, and there is no constraint. This is why there is no supernova 1987A constraint up to arbitrarily high couplings. However, in the case of stellar cooling, the photon dominates the energy loss of the stars. Since it is impossible for the dark sector particles considered in this paper to have stronger couplings to photons than SM particles, the dark sector particles will always carry away more energy than the photon. Moreover, the criteria used for the stellar cooling bounds is that the dark sector particles must carry away less than a fraction of the energy carried away by photons. Thus, only if the dark sector particles interact more strongly than photons, would they fail to carry away sufficient energy. Therefore, there is no upper boundary for the stellar cooling constraints.

Potential Reach of Terrestrial Searches
The freeze-in dark matter models discussed in this paper are challenging to detect in the laboratory with direct-detection and accelerator-based experiments. This is not surprising, given the small required couplings. Nevertheless, we illustrate this challenge in Fig. 7 for a "heavy" dark photon mediator with m = 3m χ and α D = 0.5 (unless otherwise indicated), and in Fig. 8 for dark matter interacting with an electric or magnetic dipole moment. We parameterize as usual the reference dark-matter-electron scattering cross section, σ e , and form factor for the dark matter, |F DM (q)| 2 , as [12,53] |M free ( q)| 2 ≡ |M free (αm e )| 2 × |F DM (q)| 2 (5.1) where |M free | 2 is the absolute value squared of the elastic dark-matter-(free)-electron matrix element and q is the magnitude of the three-momentum lost by the dark matter when it scatters off the electron. For each of these models, we can deriveσ e to be where µ χe is the reduced mass between the electron and χ, and v rel is the relative velocity between the incoming dark matter and the incoming electron. A dark photon mediator can be classified as "heavy" and give F DM = 1 once its mass is above the typical momentum transfer, q typ , which varies for different targets. For example, for direct-detection experiments with semiconductor or noble liquid targets, q typ ≡ µ χ,e v rel αm e [53]. So for dark photon masses above a few keV (which is enforced by the stellar constraints), we have F DM = 1. Dark matter interacting with an electric dipole moment has the form factor The form factor for dark matter interacting with a magnetic dipole moment is more complicated,  Figure 7: Solid lines in green, cyan, and red show the values of the dark-matter-electronscattering cross section for which the correct dark matter relic abundance is obtained from freeze-in for Dirac fermion dark matter coupled to a dark photon, for various choices of the dark-photon couplings α D (see Fig. 2). Red and brown-shaded regions show the stellar constraints from red giant and horizontal branch stars (see Fig. 5). Dashed lines show the potential reach of laboratory experiments. The gray shaded regions are excluded from the number of effective relativistic degrees of freedom (see Fig. 5), supernova 1987A, and existing laboratory constraints. The dotted line shows the CMB constraint, which excludes the freezeout line (blue) when the dark matter particle is a Dirac fermion. Projections and constraints for dark matter that is a scalar particle are similar, except with a much weaker CMB bound. If not stated otherwise, the model parameters are m = 3m χ and α D = 0.5. See text for details.
which is a combination of F DM = 1 and F DM = αm e /q. In deriving this form factor, we find an explicit dependence on the relative velocity between the incoming dark matter and the incoming electron in the free 2 → 2 (dark-matter-electron to dark-matter-electron) scattering. A precise calculation of the crystal form factor defined in [53] would need to take this into account. However, here we approximate v rel α and calculate the direct-detection bounds and direct-detection projections using Finally, to convert nuclear recoil cross section sensitivities to σ e , we follow [69]. We show in Fig. 7 and in Fig. 8 the sensitivity (when available in the literature) for a few future planned direct-detection and fixed-target experiments or proposals: a silicon de-   Figure 8: Solid colored lines show the values of the dark-matter-electron-scattering cross section for which the correct dark matter relic abundance is obtained from freeze-in for dark matter interacting with an electric (left) or magnetic (right) dipole moment, for various reheating temperatures (see Fig. 3). Red and brown-shaded regions show the stellar constraints from red giant and horizontal branch stars (see Fig. 6). The dashed line shows the potential reach of a direct-detection experiment using Skipper-CCDs for a 30 kg-year exposure. The gray shaded areas are excluded from the number of effective relativistic degrees of freedom (see Fig. 6), direct detection searches, and supernova 1987A. For the right plot, the region above the dotted line is excluded from the CMB. The freeze-in line for T RH = 10 MeV (orange) stops at the coupling where dark matter would thermalize with the SM sector. See text for details.
tector with a 30-kg-year exposure and single-electron threshold (using, for example, Skipper-CCDs [53,70]), a superfluid helium detector with a 1 kg-year exposure and 10 eV phonon energy threshold [71], an electron-beam fixed-target experiment searching for missing momentum (LDMX, from Fig. 5 in [72]), and an electron-positron collider searching for missing energy (Belle-II) [73] (the latter two do not have sensitivity to dipole moment dark matter in the range of parameters shown in the plot [38]). In Fig. 7, we also show in gray the bound from N eff (also seen in Fig. 5) as well as current laboratory bounds from direct-detection and accelerator-based probes, including XENON10/100/1T, DarkSide-50, DAMIC-SNOLAB, SENSEI, SuperCDMS, E137, LSND, and BaBar [73][74][75][76][77][78][79][80][81][82][83]. In Fig. 8, we show in gray the bound from N eff (also seen in Fig. 6) as well as the direct-detection bounds from [80,83]. At low couplings, the limit reaching to ∼100 MeV is from supernova 1987A (from [50] for the dark photon portal, and from [51] for dark matter coupled to an electric or magnetic dipole moment); we find that the couplings that are probed by the supernova bound lie between the freeze-in and the freeze-out line. The CMB (dotted gray line) sets a strong constraint for Dirac fermion dark matter, but is easily avoided, for example, for scalar dark matter [84,85]. The freeze-out line in Fig. 7 is almost independent of the dark photon mass as long as the dark photon mass is sufficiently far away from 2m χ , so we just present the line for the benchmark case m = 3m χ and α D = 0.5 (although see [86]). We repeat the freeze-in lines from Figs. 2 and 3 as well as the stellar cooling and other bounds from Figs. 5 and 6. As expected, the freeze-in parameters are typically too small to be probed by laboratory searches in the near future. However, interestingly, we see that the freeze-in targets for dark matter interacting with an electric or magnetic dipole moment can be probed for low reheating temperatures with upcoming direct-detection experiments.
Since the dark matter models with an electric or magnetic dipole moment are dimension 5 operators, one can ask how these are UV completed. As discussed in Sec. 2, one simple possibility is to imagine charged scalars and fermions of a common mass M generating the dipole moment operators. In Fig. 8, we show the resulting upper bound on the cross section σ e on this simple UV completion, derived from the upper bound on d χ or µ χ from Eq. (2.6): the black (gray) solid line corresponds to cross sections above which the corrections to the dark matter mass is larger than the tree-level dark matter mass for M at the scale of 100 GeV (1 TeV). Of course, different UV completions may allow for higher cross sections. We do not consider this further in this paper.

Conclusions
In this work, we discussed the freeze-in production of dark matter in the keV-to-GeV mass range as well as the constraints from stellar cooling. We considered two distinct scenarios: fermionic and bosonic dark matter that is coupled to the SM through kinetic mixing between the photon and a dark photon, as well as fermionic dark matter interacting with SM photons through an electric or magnetic dipole moment.
When the dark matter interactions with SM particles are small, the dark sector is not in thermal equilibrium with the SM in the early Universe, and dark matter production can occur through freeze-in from fermion-antifermion annihilation and from the decay of plasmons. The latter dominates for sub-MeV masses and pushes the couplings needed to obtain the observed relic abundance from freeze-in to very small values. For the dark photon portal models, the production in the early Universe is infrared dominated. In contrast, the dark matter models with a dipole moment are described by dimension-five operators, so that the freeze-in production occurs at all temperatures. Consequently, the freeze-in parameters depend on the reheating temperature. We also checked that the dark matter does not thermalize with the SM thermal bath, such that the bounds from BBN and N eff , which usually constrain dark matter with masses in the sub-MeV-range, are avoided.
In addition to deriving the freeze-in production, we also calculated the stellar cooling constraints, and find that the strongest limits are from the non-delay of helium ignition at the red giant tip. This bound excludes freeze-in production of dark photon portal dark matter with masses below, for example, 20 keV for a dark gauge coupling of α D = 0.5. For dark matter with an electric or magnetic dipole moment with T RH > 1 GeV, dark matter masses below 4 keV are disfavored from stellar cooling constraints.
Finally, we discussed the potential to probe these models in laboratory experiments. In the case of dark matter coupled via a dark photon, some part of the parameter space that can be probed by future experiments is already ruled out by the red giant constraint. Towards larger dark matter masses, the freeze-in lines are too low to be probed in the foreseeable future, but present potential targets for future, very ambitious, experiments. However, for dark matter with an electric or magnetic dipole moment, and for dark matter masses above the reheating temperature, the freeze-in production in the early Universe is suppressed; relatively large couplings are required to then obtain the correct relic abundance, so that these scenarios can be partially probed with upcoming direct-detection experiments.

A Properties Of Photons In A Thermal Plasma
In a plasma, electrons can move freely and thus affect the propagation of electromagnetic waves. They become a combination of coherent vibrations of not only the electromagnetic field, but also the electron density. Quantization leads to a spin-1 field with one longitudinal and two transverse polarization modes. We review the material needed to derive our results in this paper, basing our discussion on [40].
The effectively massive photon modes are caused by a modified dispersion relation for the photon in a plasma. For a photon in vacuum, the relation between its frequency ω and wave vector k is simply given by ω 2 = k 2 . For plasmons, this relation is subject to modifications depending on the electron density n e and temperature T . The modified dispersion relations give rise to a non-zero phase-space ω 2 − k 2 allowing for decays to massive particles. Note that in principle free protons and nuclei could also contribute to the plasma effect. However, they are much heavier than the electrons and thus more inert, so their contribution turns out to be negligible.
A characteristic quantity of a plasma is its plasma frequency (n e (E) +n e (E)) .
It is in general a function of the temperature T , as the electron (positron) density follows the Fermi distribution n e/ē = [e (E∓µ)/T + 1] −1 with the chemical potential µ. For the explicit computation it is helpful to replace v = p E . Defining allows the definition of the quantity v = ω 1 /ω p , which intuitively is the typical electron velocity. With these ingredients, the general dispersion relations valid at all temperatures and densities up to first order in the electromagnetic fine structure constant α are given by [40] The transverse mode satisfies ω T > k for all values of k. In contrast, the dispersion relation for the longitudinal mode can cross the light cone if k becomes larger than ω L . This prevents the longitudinal plasmon from propagating and constrains the longitudinal wave vector to a maximal value The renormalization of the propagator determines the propagation of plasmons. However, when interactions are considered it is useful to change from the mass to the interaction basis. The coupling to the electromagnetic current then gets renormalized.
such that the dressed polarization vectors are [51]  We now want to discuss specific limits that are helpful for our numerical implementation of the calculations. In general, as k → 0 the dispersion relations ω t/l approach the plasma frequency. For large wave numbers k T and small electron density, the situation of the vacuum is restored, ω T → k and the longitudinal mode disappears.
In the relativistic limit, T m e or µ m e , Eqs. (A.3) and (A.4) simplify as v = 1 and k max → ∞. The plasma frequency reduces to ω 2 p,rel. = 4α 3π µ 2 + π 2 T 2 3 . (A.10) In the degenerate limit, T µ − m e , the plasma frequency can be expressed in terms of the Fermi momentum p F In the dispersion relations, v can be replaced by the Fermi velocity v F = p F E F with the Fermi energy E F = p 2 F + m 2 e . In the classical limit, the electrons are non-relativistic and non-degenerate, T m e −µ. The plasma frequency is given by Most contributions to the freeze-in for the dark photon portal comes from late stages of the thermal history of the early Universe. For dark matter masses below the electron mass, the classical limit is important as production occurs partly when the electrons are nonrelativistic. At temperatures of tens of keV, the lepton asymmetry becomes important, such that the sum of the electron and positron number densities is given by In the last term, η B ≈ 6×10 −10 is the baryon to photon ratio and n γ = 2ζ 3 T 3 /π 2 is the photon number density with the Riemann zeta function value ζ 3 ≈ 1.2. Since the baryon number density seems to coincide with the number density of electrons, the last term accounts for the asymmetry.

B Inclusion of A − Z Mixing In Freeze-In Calculations
The mixing term between the dark photon and the hypercharge gauge boson in the Lagrangian reads L ⊃ − 2 cos θ W F µν B µν . (B.1) After the electroweak symmetry breaking, this term can be written with gauge boson mass eigenstates, where Z µν is the field strength of the Z boson. The second term is negligible at low energies, but can be relevant for energies larger than the GeV scale. In this work, we mainly focus on the sub-GeV scale, so the contribution from Z-mixing is less than O(10%). However, we include the contribution from Z-mixing in our calculations, and briefly summarize the relevant formula in this Appendix. For the Z-mixing contribution, we ignore plasma effects because the effects do not open a new production channel, and the correction is not significant. Also, we do not include Z-mixing for the stellar bounds, as the temperature of the stellar objects are very small compared to the Z-boson mass.

B.1 Z-Boson Decay
The last term in Eq. (3.5) describes the contribution from the Z-boson decay to a dark matter pair, which dominates for 10 GeV m DM < m Z /2. The term for the case of fermionic dark matter χ can be written as

B.2 Annihilation through the Z boson
In Eq. (3.8), we only show the amplitude for production through the photon. Here, we show the full amplitude with the Z-boson: where θ W is the weak mixing angle, g Z = e cos θ W sin θ W , C f V = T 3 f − 2q f sin 2 θ, C f A = T 3 f , and Γ Z 2.5GeV is the decay width of the Z-boson.