Freezing in, heating up, and freezing out: predictive nonthermal dark matter and low-mass direct detection

Freeze-in dark matter (DM) mediated by a light (≪ keV) weakly-coupled dark-photon is an important benchmark for the emerging low-mass direct detection program. Since this is one of the only predictive, detectable freeze-in models, we investigate how robustly such testability extends to other scenarios. For concreteness, we perform a detailed study of models in which DM couples to a light scalar mediator and acquires a freeze-in abundance through Higgs-mediator mixing. Unlike dark-photons, whose thermal properties weaken stellar cooling bounds, the scalar coupling to Standard Model (SM) particles is subject to strong astrophysical constraints, which severely limit the fraction of DM that can be produced via freeze-in. While it seems naively possible to compensate for this reduction by increasing the mediator-DM coupling, sufficiently large values eventually thermalize the dark sector with itself and yield efficient DM annihilation to mediators, which depletes the freeze-in population; only a small window of DM candidate masses near the ∼ GeV scale can accommodate the total observed abundance. Since many qualitatively similar issues arise for other light mediators, we find it generically difficult to realize a viable freeze-in scenario in which production arises only from renormalizable interactions with SM particles. We also comment on several model variations that may evade these conclusions.


Introduction
Despite overwhelming evidence demonstrating the existence of dark matter (DM) on galactic and cosmological scales [1], identifying its possible non-gravitational interactions remains one of the greatest challenges in contemporary physics. Historically, much of the DM search program has been guided by the weakly interacting massive particle (WIMP) hypothesis in which ∼ TeV-scale DM carries electroweak charge and acquires its abundance by freezing out of equilibrium from the Standard Model (SM). However, decades of null searches [2] have motivated a new, broader paradigm in which DM is the lightest stable SM singlet in a "hidden sector" with its own dynamics and cosmological history (see [3,4] for a review). In these scenarios the viable DM mass range can span many orders of magnitude, so there is strong motivation for new detection techniques.

JHEP10(2018)136
traditional direct detection experiments, designed for WIMP scattering off heavy nuclei with ∼ keV recoil energy thresholds, many of these techniques exploit observable transitions between small internal energy levels (e.g. atomic ionization) to probe much lighter DM, which typically deposits ∼ meV − few eV per scatter. Furthermore, these targets are particularly sensitive to interactions mediated by light particles φ whose recoil spectra are sharply peaked towards low recoil energies, so the scattering cross section can be greatly enhanced and may be observable even for extremely small couplings. Here g χ/f are the φ couplings to dark/visible particles where f is a SM fermion and, for illustration, we have taken m χ m e and normalized to the characteristic momentum transfer q ∼ αm e for DM scattering off atomic electrons [5]. However, if g f is sufficiently large to thermalize φ with the SM in the early universe, cosmological bounds from BBN and ∆N eff generically require m φ MeV [20], which is outside the range in which eq. (1.1) holds. Thus, in order to exploit this huge enhancement, we demand g f 1 to avoid thermalizing with the photons, which excludes thermal DM production with a sub-MeV mediator, 1 but still allows for non-equilibrium alternatives.
Aside from thermal freeze-out, the only other predictive, UV insensitive production mechanism is "freeze-in" [22,23], in which DM is not initially populated during reheating, but produced later through sub-Hubble interactions with SM particles. Since DM never achieves a thermal number density, its Boltzmann equation is integrable and the final abundance is Ω χ m χ s 0 ρ cr mχ ∞ dT T n 2 SM Hs σv SM→χχ (traditional freeze in), (1.2) where T is the SM temperature, H is the Hubble expansion rate, s is the entropy density, ρ cr is the critical density, n SM is the number density for some SM species, and a 0 subscript represents a present day value. As with freeze-out, the abundance is fixed by the DM-SM interaction and is insensitive 2 to the details of earlier, unknown cosmological epochs (e.g. inflation, reheating). Although the couplings required for freeze-in production are generically very small for light mediators, it may nonetheless be possible to exploit the low thresholds of these experiments to bring the cross section in eq. (1.1) into the observable range. It has been shown that many new direct-detection targets are potentially sensitive to the freeze-in parameter space for DM with interactions mediated by an ultra-light ( keV)

JHEP10(2018)136
dark-photon [12,24]. However, this particular model has many special features, which do not generalize to other mediators (or even to heavier dark-photons). Most notably, the in-medium thermal suppression of kinetic mixing at finite temperature counterintuitively makes the nearly-massless dark-photon safer from astrophysical and cosmological bounds because all on-shell production vanishes in the massless mediator limit [25,26]. Furthermore, in this limit, the cosmological dark-photon production is also suppressed, so it can be neglected in the Boltzmann system and eq. (1.2) can be straightforwardly applied. However, in models for which mediator production is not parametrically suppressed in the early universe, eq. (1.2) is no longer valid and it is necessary to solve the full Boltzmann system that tracks both χ and φ densities. Similarly, in such models there are nontrivial astrophysical bounds on the mediator-SM coupling g f [25,26]. Thus, it is important to identify other viable light-mediator scenarios that enjoy the enhancement in eq. (1.1) and understand the generic issues they present.
Pioneering earlier work [24] studied DM production through renormalizable portals and its subsequent self-thermalization via dark force interactions. However, this study primarily 3 considered a massless dark-photon mediator, with the above-mentioned special properties. Furthermore, [27] and [28] recently studied the models and constraints for light mediators accessible to new direct detection techniques; however, these analyses leave open the question of whether there exist predictive non-thermal production mechanisms beyond the ultra-light dark-photon mediated scenario.
In this paper we extend this discussion by studying a representative benchmark of relevance for low-mass direct detection: the predictive freeze-in production of dark matter through a light, feebly-coupled scalar mediator φ with Higgs portal mixing. We find that even when g f 1 and all SM-DM interaction rates are slower than Hubble, there is generically a sizable φ population, which can thermalize with χ if g χ is sufficiently large. For fixed, viable choices of g f , there are three qualitatively distinct regimes: • Small coupling: if g χ is sufficiently small, the hidden sector never thermalizes with itself, so hidden-annihilation is always sub-Hubble Γ(χχ → φφ) H, and production approximates traditional freeze-in. In this regime, the χ abundance increases with g χ and eq. (1.2) is approximately valid.
• Large coupling: if g χ is sufficiently large, hidden-annihilation becomes efficient Γ(χχ → φφ) H and the hidden sector thermalizes with itself (but still not with the SM). In this regime DM production initially resembles freeze-in, but eventually the χ-φ system reaches thermal equilibrium (at a lower temperature than the SM) and χ undergoes an additional phase of freeze-out. In this regime, increasing g χ merely converts more of the χ into dark radiation.
• Intermediate coupling: at intermediate values of g χ for which the maximum hiddenannihilation rate is briefly comparable to Hubble Γ(χχ → φφ) ∼ H, production JHEP10(2018)136 interpolates between the above regimes and yields a maximum DM abundance for each choice of mass.
Intriguingly, we find that for viable mediator-SM couplings g f 1, safe from astrophysical and cosmological bounds, it is generically difficult to account for the total freeze-in abundance with only renormalizable couplings to SM particles. Nonetheless, for choices of g f that saturate existing limits, the maximum abundance that can be produced is significantly greater than the naive expectation from eq. (1.2) once mediator-DM interactions are included in the Boltzmann system. This paper is organized as follows: section 2 introduces our benchmark model; section 3 describes the early universe production mechanism; section 4 discusses model variations that evade some of our general conclusions; and section 5 we summarize our findings and suggest future directions of inquiry.

Benchmark model
Our benchmark scenario involves a scalar singlet mediator φ and mixes with the SM through the Higgs portal operators where A and κ are the renormalizable portal couplings. After electroweak symmetry breaking (EWSB) the neutral component of the doublet h mixes with φ and, in the small mixing limit, the system is diagonalized with the shift h → h + sin θ φ, so φ acquires massproportional couplings to SM fermions. Including a φ Yukawa interaction with a fermionic DM candidate χ, the relevant Lagrangian becomes where f is a SM fermion of mass m f and v 246 GeV is the SM Higgs vacuum expectation value. Although this is only one of several possible scalar mediated scenarios, as we will see, it captures much of the essential physics and many of the issues encountered here apply to a much broader class of variations on this simple setup (see section 4 for a discussion and [29,30] for a complementary studies of thermal freeze-out with the same field content). In this model, the cross section for nonrelativistic χe scattering is where q is the three-momentum transfer, µ χe is the reduced mass, and we take the light mediator limit to highlight the parametric enhancement from eq. (1.1) in terms of q ∼ αm e , the characteristic momentum transfer for χ scattering off atomic electrons [5]. Our main results (described in section 3 and displayed in figure 5, right panel) are presented in terms of the convention in eq. (2.3) and are largely independent of the mediator mass so long as m φ αm e , m χ .

JHEP10(2018)136
h f f Figure 1. Leading Feynman diagrams giving rise to non-thermal χ freeze-in production from SM annihilation during the early universe (left) and reannihilation into dark sector mediators φ (middle, right) once the hidden sector thermalizes with itself. Although the SM initiated process must be slower than Hubble expansion to avoid thermalizing the dark sector, the χ − φ interactions can achieve chemical equilibrium at a separate temperature set by the energy density transmitted to the dark sector via the left diagram and those depicted in figure 2.  Figure 2. Representative Feynman diagrams giving contributing to φ freeze in production in the early universe. These processes determine ρ φ which (in addition to ρ χ ) determines the hidden sector's temperature T D . The effective vertex in the right diagram arises from integrating out the top quark for energies in between the weak scale and the QCD confinement scale (see appendix B).
Although our emphasis in this paper is on the freeze-in production for this scenario, it has also been shown that models with light scalar mediators m φ ∼ 10 MeV can yield sufficiently strong DM self interactions σ χχ /m χ ∼ cm 2 g −1 to reduce tension between DMonly N-body simulations and the observed small scale properties of haloes [31]. It may be interesting to explore whether DM with a scalar mediator in the mass range considered here m φ αm e could preserve some of these features (as in [29] for heavier scalars), but this question is beyond the scope of the present work.

Boltzmann equation
Our initial condition assumes a hot, radiation dominated universe in the broken electroweak phase. 4 We further assume that n χ = n φ = 0 after reheating so that sub-Hubble interactions with SM fields are solely responsible for populating the hidden sector. 5 If all SM initiated reactions are slower than Hubble expansion, the Boltzmann equation for χ JHEP10(2018)136 10 -10 Equilibrium Yields, gχ = 10 -3 , sin θ = 10 -10 10 -2 10 -1 10 0 10 1 10 2 10 3 10 -5 Maximum Abundance, eV ≲ mϕ ≪ αme, sin θ = 10 -10 Figure 3. Left: comoving equilibrium yields Y eq χ = n eq χ (T D )/s(T ) plotted against the dimensionless time variable z = m χ /T for a various m χ . Note that the dark temperature is determined primarily by the φ energy density, which is independent of m χ , so each mass point has comparable T D during production; thus, sufficiently large masses are highly Boltzmann suppressed. Right: Maximum DM abundance for a range of m χ with fixed sin θ = 10 −10 , which is still viable, but near the Red Giant exclusion bound (see left panel of figure 5). Below m χ 200 MeV, the blue curve represents the coupling for which the maximum abundance (red curve) can be achieved; larger values merely annihilate away more of the DM as the number density tracks the equilibrium distributions (left) out to larger values of z (see figure 4, left). In the green band between the dashed vertical lines, max Ω χ h 2 > 0.12 so it is possible to achieve the full abundance for two choices of couplings, as shown for a representative mass point in figure 4. production is where H ≡ 1.66 √ g * T 2 /m Pl is the Hubble rate during radiation domination, m Pl = 1.22 × 10 19 GeV is the Planck mass, and g * is the number of relativistic degrees of freedom. Here an eq superscript denotes an equilibrium quantity, and we do not assume equilibrium between visible and hidden sectors (T = T D ); all relevant SM species are in equilibrium with each other, so we omit their superscript. As discussed above, existing bounds on the mixing angle require that the dark sector never thermalize with the SM in the early universe, so n χ n eq χ (T ), but we allow for the possibility that χχ ←→ φφ annihilation (shown in the middle and right diagrams of figure 1) equilibrates the hidden sector so that χ eventually achieves a thermal distribution n χ → n eq Note that if hidden-sector annihilation rate is also sub-Hubble throughout the production process Γ(χχ → φφ) H, eq. (3.1) recovers the usual freeze-in form in eq. (3.1).

Dark temperature evolution
In this scenario, the hidden sector temperature T D is defined in terms of the evolving dark sector energy density, which is also produced through SM interactions. The energy transfer

JHEP10(2018)136
Boltzmann equation for direct χ production is where P χ is the χ pressure density and "dark processes" represents all interactions that exchange φ and χ energy independently of the SM radiation bath. Similarly the contributions from mediator production (see figure 2) are Here we define the thermally averaged quantity where ∆E is the energy transferred to the hidden sector and v mø is the Møller velocity (see appendix A). Note that the definition of the energy transfer differs depending on whether we are considering χ pair production in eq. (3.2) for which ∆E tracks both final state particles, or single φ production in eq. (3.3) where it tracks only the energy transferred to φ. The combined dark sector energy density ρ D = ρ φ + ρ χ defines a dark temperature where ξ D counts the relativistic degrees of freedom in the hidden sector. Note that the dark temperature initially satisfies T D = 0 and depends nontrivially on the visible temperature through eqs. (3.2) and (3.3), which source its growth in eq. (3.5). The dark temperature also determines the equilibrium χ number density in the hidden sector which enters into the r.h.s. of eq. (3.1), where ξ χ is the number of χ polarization states. In our Boltzmann system we have omitted all processes involving lighter SM fermions, which only contribute negligibly to χ and φ production in the early universe; the dark sector is mainly populated while heavy electroweak states are in equilibrium with the SM thermal bath with corrections of order ∼ (m f /m t ) 2 for all lighter SM fermions. For this reason, we can safely neglect longitudinal plasmon mixing effects which can resonantly enhance the φ production at lower temperatures for which m φ is near the plasma frequency [26]; this enhancement does not compensate for the parametric suppression in the coupling relative to continuum production off heavy electroweak states. However, such an effect may be important in model variations in which φ production arises mainly from its coupling to electrons (see section 4 for a discussion). Although freeze-in is insensitive to the precise 10 -3 10 -2 10 -1 10 -9 10 -10 10 -11 10 -12 Comoving Yield, m χ = 1 GeV, sin θ = 10 -10

Full Calculation
Naive 'Freeze-In' value of the reheat temperature, we assume it to be above m t for all SM electroweak states to be relativistic when dark sector production begins; if the reheat temperature is lower, DM production will proceed through lighter fermions and the resulting abundance will be significantly suppressed, but the qualitative features of our discussion are unaffected. We note that a full thermalization also requires that φ-number changing reactions (e.g. 2φ → 3φ) reach equilibrium in addition to demanding that the χχ ←→ φφ annihilation rate exceeds Hubble expansion. To see that this is generically true, note that prior to thermalization both χ and φ number densities are determined by SM initiated processes in eqs. (3.2) and (3.3). However, direct χ production in eq. (3.2) always involves an additional penalty of g 2 χ 1, so n χ ∼ g 2 χ n φ just prior to thermalization. Thus, when the χχ → φφ annihilation rate thermalizesn χ σ(χχ → φφ) > H -the 2φ → 3φ rate is also in equilibriumn φ σ(2φ → 3φ) > H -because the greater φ density easily compensates for the O(10 − 100) phase space penalty for the 3-body final state as long as g χ 0.1 and the φ quartic coupling λ φ g χ . The condition on g χ is trivially satisfied in our parameter space of interest in section 3.3 and the condition on λ φ can easily be realized without loss of essential generality.

Numerical results
In this section we present our main findings as shown in figures 3, 4, and 5. Before commenting on specifics of these plots, a few general comments are in order. The χ production mechanism described in this paper can best be understood by considering two extreme limits:

JHEP10(2018)136
• Freeze-in regime (g χ 1) If the DM-mediator coupling is very small, the dark sector never thermalizes with itself, so χ production is very similar to the more familiar freeze-in mechanism. In this regime χφ interactions are negligible, so the χ abundance arises from SM initiated processes. Since all such reactions require virtual φ exchange, the rates for these reactions scale as g 2 χ sin 2 θ, so increasing g χ increases the χ abundance at late times.
• Dark thermalization regime (g χ ∼ 1): if the DM-mediator coupling is large, both χ and φ abundances initially grow as a function of time since n χ = n φ = 0 as our initial conidtion. However, once the χχ ←→ φφ annihilation rate exceeds Hubble expansion, the dark sector thermalizes with itself even though neither χ nor φ ever thermalizes with the SM due to the feeble Higgs-φ mixing angle. This thermalized dark sector has its own (lower) temperature T D < T , so around T D ∼ m χ , the χ population freezes out via χχ → φφ annihilation. In this regime, increasing g χ decreases the χ abundance at late times.
In between these limiting regimes, there is a critical value of g χ which maximizes Ω χ for a given choice of m χ . For such critical values of g χ the χ − φ thermalize with each other and freeze out quasi-relativistically T D ∼ m χ just before Boltzmann suppression depletes the χ abundance (see figure 3) If this theoretical maximum abundance exceeds the observed DM abundance for a given value of m χ , Ω max χ h 2 > 0.12, then there are two choices of g χ that can yield Ω χ h 2 = 0.12; one smaller value in the freeze-in regime and the other in the dark thermalization regime. However, as we will see below, for most choices of m χ the maximum achievable abundance falls far short of Ω DM so χ can only be a subdominant fraction of the total DM density.
For the relevant parameter space considered in this paper m φ αm e , the hidden sector temperature is primarily set by ρ φ whose growth rate is proportional to sin 2 θ in eq. (3.3), whereas the collision terms for ρ χ in eq. (3.2) all depend on g 2 χ sin 2 θ and are further suppressed. Furthermore, ρ φ continues to accumulate through gg → gφ processes long after direct χ production processes (e.g. tt → χχ) have frozen out. Thus, the cosmological evolution of T D is insensitive to the DM mass, so sufficiently large values for which m χ T D encounter hidden-sector Boltzmann suppression throughout their production, φ-thermalization, and re-annihilation phases.
The left panel of 3 presents the equilibrium χ yield Y eq χ ≡ n eq χ (T D )/s(T ) normalized to the total entropy, where n eq χ (T D ) is computed using eq. (3.6) and depends on the visible temperature implicitly through eq. (3.5). Here we see that, as m χ increases for a fixed choice of g χ , the heights of the equilibrium curves begin to fall, which limits the total abundance that can be generated for a given mass point; increasing g χ increases the abundance in the freeze-in regime, but once the dark sector equilibrates with itself, the yield tracks these equilibrium distributions and eventually diminishes as more of the DM annihilates away into φ radiation. Thus, for each choice of mass there is a maximum abundance as represented by the red curve in figure 3 (right panel), which corresponds to the value of g χ from the blue curve. Note that inside the green band, the blue curve splits into

Cosmology
(gχ = 0)  Higgs-Portal Mediator, sin θ = 10 -10 , eV ≲ mϕ ≪ α me Figure 5. Left: available parameter space for a light, Higgs-mixed scalar particle plotted alongside constraints from 5th forces [10,[35][36][37][38][39], stellar cooling [26], and late time φ → γγ decays [40] labeled "Cosmology". The vertical line represents the approximate boundary where the χe scattering rate undergoes a parametric shift for scattering off atomic electrons, which is taken to be the reference value for presenting different techniques in [4]. For other materials with different characteristic momenta, this shift takes place for different values of m φ . Note that the cosmological bounds from [40] assume only a single φ without any additional contributions from the presence of an additional dark sector; including these contributions would only serve to strengthen these bounds, which further motivates values of m φ < αm e . Right: constraints and projections for the χe cross section for sub-GeV DM. Some representative projections for the SENSEI 2e − , Silicon 1e − , and superconductor targets are taken from [4], and the Dirac metals projection is from [15]; all projections assume Ω χ Ω DM . For points outside the green vertical band, the blue curve depicts (max Ω χ /Ω DM )×σ e where sub-Hubble production as described in section 3 generates the maximum fractional abundance (max Ω χ < Ω DM ) for each m χ . Inside the vertical green band, the blue curves represent σ e evaluated at g χ for which the full DM abundance can be accommodated; note there are two such values as depicted in figure 4 and on the right panel of figure 3. Although the Higgs portal mixing sin θ = 10 −10 may seem ad-hoc, this is chosen as a representative viable value near the stellar cooling exclusion region depicted on the left panel of this figure.
two segments corresponding to the distinct values that yield the observed DM abundance Ω χ h 2 0.12, and not the maximum abundance (as it does outside this band).
This counterintuitive behavior can be understood from the plots in figure 4 which shows the production history for a representative mass point m χ = GeV and sin θ = 10 −10 as the coupling g χ is varied. On the left panel, the red curves are solutions to the Boltzmann equation in eq. (3.1) for two choices of coupling that yield the observed DM abundance; the blue curves at the bottom of the plot also show the "naive" freeze-in result for the same choices of g χ ; as in eq. (1.2), the blue curves are computed by merely integrating the SM collision terms in eq. (3.1) without including the additional χχ ←→ φφ dynamics. The right panel of figure 4 presents the DM abundance as a function of g χ as computed with the full Boltzmann system (red curve) and with the naive freeze-in integration (blue curve).
On the right panel of figure 5 we present our main result: the predictive parameter space for χ production in the early universe plotted in terms of the fiducial χe scattering where A(q 0 ) is the matrix element for a free-electron target evaluated at reference momentum q 0 = αm e following the convention in [4], which also applies to the other projections shown. For the blue curve in this plot, this cross section is given by eq. (2.3) in the m φ αm e limit and evaluated at g χ values, which either yield the full DM abundance or the largest possible subdominant fraction where this is not possible as shown in figure 3 (right); in latter case the cross section is rescaled by the fractional abundance. The vertical green band in plot highlights the mass window for which the total abundance can be achieved and, as in figure 3 (right), the blue curve splits into two segments in the shaded green band to reflect the two values of g χ that can accommodate the total abundance in; the top segment corresponds to χ that yield freeze-out in the hidden sector, whereas the bottom segment presents the freeze-in like value. Note that for most of the parameter points shown here max Ω χ 10 −2 Ω DM (see figure 3, right), so DM self interaction bounds do not constrain the viable parameter space; for a discussion of this issue see [27]. Relatedly, there is no perturbative unitary violation in χχ scattering for the g χ along this curve since the values that yield the maximum possible abundance for viable choices of sin θ are all well below unity (see figure 4, right).
For completeness, we also show the viable parameter space for a Higgs portal mediator on the left panel of figure 5 to motivate our benchmarks. The cosmological bound presented here is taken from [40], which studied a minimal extension of the SM with only the additional scalar particle φ, which mixes with the SM Higgs. In our scenario, the dark sector also contains χ which can increase the φ production through additional χχ → φφ annihilation not considered in [40]. Thus, the cosmology region shown in figure 5 (left) is conservative; larger values of g χ can make this bound stronger, but a full treatment of this bound is beyond the scope of this paper. This plot shows that, even in the least constrained scenario where the bounds on sin θ are unaffected by additional φ production, thereby maximally allowed choices for m φ < αm e must still satisfy sin θ 10 −9 , thereby severely limiting the DM abundance that can be achieved via freeze-in like dynamics.
Finally, we also note that for the couplings considered here, χχ → χχ elastic scattering may also constrained by observations of galaxy halo shapes, cores, and cluster collisions (see [41] for a review). However, for most of parameter space considered in this paper, the maximum achievable χ abundance (see right panel of figure 3) only constitutes a subdominant fraction of the total DM density. Thus, the usual σ χχ /m χ 0.1cm 2 /g bound does not apply to this subdominant DM fraction. However, inside the green band in figures 3 (right) and 5 (right), there is a window around m χ ∼ GeV where the total abundance can be achieved for g χ ∼ 10 −3 . To estimate the constraint in this window, we invoke the momentum transfer cross section which removes the forward scattering contributions from the overall result

JHEP10(2018)136
where v is the DM relative velocity and in the last equality we use the result from mediator in the m χ v/m φ 1 regime [41]. For reference values m χ = GeV, m φ = αm e , v = 10 −3 and g χ = 10 −3 we find that σ T /m χ = 1.6 × 10 −3 cm 2 /g, which is safe from self interaction bounds. For smaller masses where χ can only be a small subcomponent of the DM, this bound does not apply.

Model variations
Thus far, we have restricted our attention to the scenario in which both χ and φ are produced through h-φ mixing during the early universe. We have found that the bounds on sin θ 5 × 10 −10 (see left panel of figure 5), severely limit the total energy density transferred to the hidden sector, thereby making it impossible to achieve the total DM abundance for m χ few 100 MeV for all values of g χ . Here we consider some possible solutions to this problem that still feature sub-Hubble DM production and light scalar mediators.

Different SM current
Since much of the discussion above depends on the Higgs-portal coupling relations, it may be possible to evade some of the above complications with a different pattern of couplings. For the scenario studied in this paper, the main impediment to generating the observed DM abundance over the eV m φ αm e parameter space of interest is the limit sin θ 5 × 10 −10 , which is driven primarily by resonantly enhanced φ production through its coupling to electrons in Red Giants (RG) [26] -see figure 5 (left). For a minimal Higgs-mixing scenario, the same mixing angle θ rescales the φ coupling to electroweak states and thereby also suppresses the total energy energy/number density transferred to the hidden sector via electroweak production in the early universe. However, the early universe production is dominated by the heaviest SM states (t, h, W/Z) at T ∼ v, so if the φ coupling to these particles can be enhanced relative to the electron coupling in a minimal Higgs-mixing scenario it is possible to increase the overall abundance at late times.

Two Higgs doublets
One simple extension, which can accommodate such a variation involves mixing φ with the the neutral CP-even states of a two-Higgs doublet model (2HDM) to break the usual SM relations between quark and lepton couplings. For instance, in a Type II 2HDM, φ could mix with the predominantly up-type doublet, thereby increasing the relative ratio of t/e couplings and enhancing dark sector particle production in electroweak processes. However, it is not clear whether a concrete realization of such a scenario can simultaneously accommodate the observed DM abundance while satisfying the existing limits on 2HDMs, which typically force the 125 GeV CP even scalar close to alignment limit, which approximately reproduces the usual SM coupling patters [42,43], so a careful study is necessary to see if this is possible.

Vectorlike fermions
Alternatively, the scalar mediator can couple to heavy vectorlike quarks, which are integrated out to induce effective interactions with SM particles (e.g. φ G µν G µν , φ F µν F µν , or φ qq if vectorlike states mix with SM quarks). If the vectorlike quarks are sufficiently decoupled, all dark sector production proceeds through effective SM interactions and therefore this class of models makes a firm prediction for the DM-SM scattering rate. Furthermore, in this variation there is no tree-level electron coupling, so Red Giant bounds can be significantly weakened [26] and more DM can be produced via freeze-in, but this scenario is harder to test with new detection techniques, most of which probe DM-electron scattering in various materials.
If, instead, the mediator couples to heavy vectorlike leptons, the IR theory will contain the φ F µν F µν operator (and also φ¯ if the new heavier states mix with right handed leptons). However, this variation does not ameliorate the Red Giant bounds, which now constrain all early-universe production. Unlike in the Higgs-mixing scenario, where production off heavier electroweak states in eqs. (3.2) and ( 3.3) is parametrically enhanced relative to the electron scattering cross section, here the production and detection depend on the same couplings, so dark sector production is greatly suppressed. Although the dominant φ production in this leptophilic scenario occur at lower temperatures T ∼ m φ where there is a resonant enhancement from φ-longitudinal-plasmon mixing, this effect has been found to be modest (order-few) relative to non-resonant continuum production in astrophysical contexts [26], so it is unlikely to compensate for the large reduction in the overall rate. Nonetheless, a full calculation of resonant φ production at these lower temperatures is worth further investigation. 6

Asymmetric freeze-in
For the the scenario considered in earlier sections, the abundance is ultimately limited by the additional χχ → φφ annihilation phase that becomes efficient at depleting the χ population for sufficiently large g χ . However, if the DM carries a conserved global quantum number and has an asymmetric population at late times, then increasing g χ will not lead to further depletion once all antiparticles are annihilated away. In order to maintain the qualitative features of the scenario considered here (light mediator and sub-Hubble production), the freeze-in interaction that produces the DM must also satisfy the Sakharov conditions [45]. Although models of asymmetric freeze-in have been proposed in the literature [46], the asymmetry they produce is suppressed by subtle cancellations [47] and it is not clear whether these limitations can be overcome, but this problem deserves further study (see [48] for a discussion).

Additional production source
One simple way to enhance DM production is to introduce a new source of χ or φ production beyond just the Higgs portal coupling to SM particles in eq. (2.2). For instance φ could JHEP10(2018)136 also mix with a new, heavy singlet scalar S which thermalizes with the SM in the early universe. Since the S-φ mixing angle is a new free parameter, it is possible to increase the χ and φ abundances via S initiated scattering and decay processes, which still realize freeze-in production if this mixing is suitably small. However, this modification is no longer predictive and manifestly depends on presently unknown UV physics, thereby detracting somewhat from the appeal of the freeze-in mechanism.
In addition to changing the SM-φ flavor structure, the SM extensions discussed in section 4.1 can provide additional production sources if the new heavier particles are populated on shell after reheating. New interactions from these states can enhance φ and χ production in the early universe, but the magnitude of this enhancement now depends on the details of the extended sector(s), which adds UV sensitivity to the production mechanism.

Concluding remarks
The recent burst of creativity in devising new techniques sub-GeV direct detection may soon make it possible to probe feebly coupled "freeze-in" dark matter coupled to a very light keV mediator. However, the only predictive, UV complete example in which this has been demonstrated (DM coupled to an ultra-light dark-photon) has subtly special properties that may not generalize to other mediators. To understand how robustly such scenarios can be probed, we have studied the freeze-in production of dark matter through a light scalar mediator with Higgs portal mixing.
In our analysis we find that the recently updated stellar cooling bounds on light scalars have severe, generic implications for any such freeze-in scenario. Most notably, for any sufficiently light ( keV) scalar mediator, compensating for the suppression from this bound requires a large, order-one coupling to dark matter, which quickly thermalizes the dark sector with itself. Even though all interactions between dark and visible sectors remain slower than Hubble expansion, the dark sector's separate thermal bath enables efficient dark matter annihilation into mediators, thereby depleting its abundance once the mediator-darkmatter coupling increases beyond the threshold for hidden sector thermalization. Thus, for each choice of dark matter mass and Higgs mixing-angle, there is a maximum cosmological abundance that can be accommodated if the dark sector is only produced through its interactions with the Standard Model. We find that for Higgs portal mixing-angles near the limit imposed by astrophysical bounds, it is impossible to generate the observed dark matter abundance except for a narrow window in the ∼ 100 MeV -few GeV range. Furthermore, even inside this range, the cross sections for direct detection off electron targets are several orders of magnitude below future sensitivity projections for proposed experimental techniques.
We have also considered some possible model variations, which may alter these conclusions and allow for successful freeze-in with a light ( keV) scalar mediator. In order to produce more DM during the early universe, while satisfying astrophysical bounds (which mostly constrain the mediator-electron coupling), it may be possible to generalize the minimal Higgs-portal scenario sector to mix the mediator with a larger electroweak sector (possibly a 2HDM) and thereby enhance the top-mediator coupling. This modification JHEP10(2018)136 increases early universe freeze-in production while heavy electroweak states are still in equilibrium, but maintains an appreciable coupling to electrons for direct detection at late times. Similarly, it may be possible to enhance DM production by preferentially coupling the light mediator to a heavy vectorlike fourth generation of quarks, which are integrated out to generate higher dimension operators with SM quarks and gluons. Such an extension may alleviate astrophysical bounds on the mediator coupling and thereby enable viable freeze in production. Other possible variations include either additional (beyond Standard Model) sources of mediator or DM production and/or additional CP violation to realize a light-mediator variation on asymmetric freeze-in, but these possibilities are beyond the scope of the present work and deserve further study.

A Thermal averaging
In this appendix we apply the methods used in [49] to compute thermal averages for annihilation cross sections and energy transfer rates in the early universe. H ≡ȧ/a = 1.66 √ g * T 2 /m P l is the Hubble expansion rate during radiation domination and the collision term for the number density is where |M| 2 is the spin averaged squared matrix element and we have omitted bose(pauli) enhancement(blocking) factors. Using detailed balance f eq χ (p 1 )f eq χ (p 2 ) = f eq A (p 3 )f eq B (p 4 ) and following the derivation in [49], the Botlzmann equation becomes where the thermally averaged annihilation cross section is defined to be

JHEP10(2018)136
where K 1,2 are Bessel functions of the first and second kinds, respectively. In the traditional "Freeze In" regime, the initial condition is n χ (0) = 0 and the χ production rate in SM annihilation processes is slower than Hubble expansion, so the n 2 χ term can be neglected and the Boltzmann equation becomes integrable and we can define comoving yields Y i ≡ n i /s and change the time variable to z ≡ m χ /T to obtain which can be integrated to obtain the asymptotic abundance at late times where s 0 2969 cm −3 is the present day CMB entropy, ρ cr = 8.1h 2 × 10 −47 GeV 4 is the critical density.
Energy density collision terms. Here we generalize the argument in [49] to calculate the thermally averaged energy transferred to particle 1 (which we will later identify with χ or φ) in a 1 + 2 ←→ 3 + 4 process governed by the Boltzmann equation where P 1 and ρ 1 are respectively the pressure and energy densities for species 1. Ignoring quantum statistical factors and assuming the other particles are in equilibrium with the radiation bath, the collision term is j p j |M| 2 E 1 f eq 3 (p 3 )f eq 4 (p 4 ) − f 1 (p 1 )f eq 2 (p 2 ) , (A.9) where we make no assumption about particle 1 being in thermal equilibrium with the others. Using detailed balance f eq 3 (p 3 )f eq 4 (p 4 ) = f eq 1 (p 1 )f eq 2 (p 2 ) and the definition of the cross section where F = (p 1 · p 2 ) 2 − m 2 1 m 2 2 is the usual flux factor, so we get Since the production/absorption rate is always slower than Hubble expansion for our purposes, the f 1 term can be neglected, so we get C ρ = d 3 p 1 2E 1 (2π) 3 d 3 p 2 2E 2 (2π) 3 4F σ 12→34 E 1 f eq 1 (p 1 )f eq 2 (p 2 ) = P 12→34 n eq 2 n eq 1 , (A.12)

JHEP10(2018)136
where the thermally averaged power-transfer rate is the Møller velocity is v m ø ≡ F/E 1 E 2 , and the denominator factor can be evaluated analytically (A.14) To evaluate the numerator of eq. (A.13), we follow the procedure in [49] we perform all trivial angular integrations and perform the coordinate transformation (E 1 , E 2 cos θ) → (E + , E − , s) where θ is the CM angle between the incoming three-vectors, E ± = E + ± E − , and s = (p 1 + p 2 ) 2 to get where s * ≡ (m 1 + m 2 ) 2 and we have used ε = 1 − s * /s E 2 + − s, so the power transfer becomes 16) which is valid even in the m i → 0 limit provided that the Bessel functions are properly expanded to cancel the mass dependence in the prefactor. Using detailed balance, we can relate the resulting collision term to the reverse process P 12→34 n eq 1 n eq 2 = P 34→12 n eq 3 n eq 4 , (A. 17) which describes the SM+SM → SM φ processes in eq. (3.3) where we identify φ with particle 1 of this argument. For processes in which 2χ particles are produced via SM annihilation, a similar argument yields the expression which tracks the energy transfer to both hidden sector particles produced in this process. whose form is appropriate for the collision terms in eq. (3.2).

B Cross sections
DM hidden sector annihilation χχ → φφ. In the hidden sector, the χχ → φφ annihilation cross section is

JHEP10(2018)136
where again we have summed (but not averaged) over all colors. Note that the r.h.s. is always positive in the physical phase space where t < 0 and that there is a t-channel, forward scattering singularity as t → 0, which is regulated by the QCD Debye mass m 2 g (T ) = 8πα s T 2 at finite temperature.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.