Dark Matter freezeout in modified cosmological scenarios

We study the effects of modifying the expansions history of the Universe on Dark Matter freezeout. We derived a modified Boltzmann equation for freeze-out for an arbitrary energy density in the early Universe and provide an analytic approach using some approximations. We then look at the required thermally averaged cross sections needed to obtain the correct relic density for the specific case where the energy density consists of radiation plus one extra component which cools faster. We compare our analytic approximation to a numerical solutions. We find that it gives reasonable results for most of the parameter space explored, being at most a factor of order one away from the measured value. We find that if the new contribution to the energy density is comparable to the radiation density, then a much smaller cross section for Dark Matter annihilation is required. This would lead to weak scale Dark Matter being much more difficult to detect and opens up the possibility that much heavier Dark Matter could undergo freezeout without violating perturbative unitarity.


I. INTRODUCTION
Dark Matter is still one of the biggest unsolved puzzles in modern physics. It continues to escape detection in both direct and indirect detection experiments while collider experiments have yet to be able to identify a statistically significant Dark Matter signal. If one wishes to construct a particle physics model which can explain Dark Matter, then they must ensure that its interactions are quite weak to avoid all these constraints. For recent reviews on the topic of Dark Matter, see [1][2][3][4][5][6] One measurement we do have is the Dark Matter relic density, Ω CDM h 2 = 0.120 ± 0.001 [8]. This measurement has been a key part of model building as a robust model must be able to explain this observation with some production mechanism. A popular production mechanism that has been studied is freeze-out. In this scenario, the Dark Matter starts in thermal equilibrium with the Standard Model at an early time but eventually decouples once the rate of the reactions maintaining chemical equilibrium become comparable to the expansion rate of the universe [1][2][3][4][5][6][7]. This is the usual mechanism used for Weakly Interacting Dark Matter (WIMPs) models [1][2][3][4][5][6][7].
Ensuring that freezeout produces the right amount of Dark Matter can fix some of the couplings in a model. This can lead to a model of Dark Matter being ruled out if the required couplings are too big to simultaneously avoid experimental bounds while producing enough Dark Matter. This has led to models which add more free parameters such as additional interactions between the Dark sector and the Standard Model, including multiple stable Dark Matter species, or considering a different production mechanisms all together.
However, WIMPs are very popular because of the socalled WIMP miracle where one could achieve the correct relic density through freeze-out with weak scale masses * apoulin@physics.carleton.ca and couplings. Many extensions to the Standard Model naturally or can easily accommodate a WIMP like particle. Because of this, it is important to understand how freezeout works and how modifications to the standard cosmological picture could modify the freeze-out results.
In the standard model of cosmology, freezeout occurs during the radiation dominated era. However, because we do not have any observations before Big Bang Nucleosynthesis (BBN), we cannot say for sure what the expansion history of the universe was before BBN. Namely, we cannot exclude the possibility of an additional contribution to the energy density that cools at a faster rate than radiation. This can come about from alternative cosmological models like alternative Dark Energy models [19,20], anisotropic expansion [12,13], brane world cosmology [14][15][16][17][18], some inflaton models [21][22][23], some quintessence models [24-26, 31, 32], or scalar-tensor gravity [33,34]. There has also been work on early matter dominated eras which decay away [27][28][29][30]. All these models have a different expansion rate from the standard Λ CMD model and some previous work has been done to understand the physics in these scenarios [9][10][11]36].
In this work, we develop a general freeze-out equation that can be used for any modified expansion rate and includes the effects of changes in entropy density. In section II, we go through the derivation of the changes to the Boltzmann equations. In section III, we take an analytic approach to obtain a simple way to approximate the relic density. In section IV, we look at numerical results from solving the modified Boltzmann equations for the case of a simple change to the energy density and compare these results to our approximations. Finally, in section V, we summarize our conclusions. Many of the definitions and approximations are the same as those of Kolb and Turner [7], namely CP conservation of the matrix elements, no Fermi degeneracy, no Bose-Einstein condensates, and that the temperature of the Dark Matter will be the same the photon temperature until after freeze-out.

A. Single effective component
In this section, we will derive the Boltzmann equation for a general energy density ρ(T ) where T is the photon temperature. We will also use the equation of state p(T ) = w(T )ρ(T ) which relates the energy density to the pressure p(T ). It will be important to not consider w(T ) as simply a constant for this treatment. This energy density can be any well-behaved function of temperature. We start with the Friedman equation: where G is Newton's constant. To take the expansion of the universe into account, it is typical to use the entropy density as a fiducial quantity. It is defined as: We will assume that the there is no entropy injection during freezeout so that the entropy per comoving volume S is conserved. This results in in s ∝ a −3 where a is the scale factor, which givesṡ = −3Hs.
We now introduce the dimensionless quantity x = m/T where m is some appropriate mass scale, typically the mass of the Dark Matter species. It will be important to find an expression forẋ. This can be done by looking aṫ s:ṡ The usual Boltzmann equation for a single species is written as [7]:ṅ where n is the number density, σv is the thermallyaveraged DM annihilation cross-section times velocity, and n eq is the equilibrium number density given by 1 : where g is the number of internal degrees of freedom, m is the particles mass, and K 2 is the modified Bessel function of the second kind. The approximation is valid for large x or small T . 2 We now make the substitutions Y eq = n eq s .
The derivative of Y with respect to time is: where we usedṡ = −3Hs. Combining Eq. 4, 5, and 10, we obtain: Eq. 12 can be used for any function of ρ and w. It is currently written for the case of a single Dark Matter species, but extending it to the case of multiple Dark Matter species is straight forward. See [37,38] for more details on the multi-species case.

B. Multiple energy densities
Equation 12 is valid for any initial ρ(T ). However, in the case where one simply adds an extra contribution to the energy density, we can write the energy density as multiple contributing energy densities which is not only simpler than treating w as a function of temperature, but also makes the equation clearer. In this section, we look at the case where we can write ρ = i ρ i with p i = w i ρ i with w i all constant. Following the same steps used to obtain Eq. 4, we find Combining Eq 5, 10, and 13, we obtain: The radiation energy density in the early universe is given by: where g * (x) is the number of relativistic degrees of freedom. If we set ρ = ρ r in Eq 14, we obtain the usual result The d log g * d log x term is usually ignored because it is small, but we include it in our analysis for completeness.
As a final example, we follow the example from [36]. Let where ρ D , n D > 4, and η = ρ D (T0) ρr(T0) are all constants, T 0 is some constant reference temperature and x 0 = m/T 0 . In our analysis, we will assume that p D = w D ρ D with 3(1 + w D ) = n D . Putting this into a similar form as Eq. 16 which only included radiation, we obtain: The terms in the square bracket are often neglected because in the radiation dominated case, it is approximately 1. However, if ρ D dominates resulting in a large value of η, then it is clear that this term becomes an important contribution.
the same as the mass of the particle as it simplifies some of these equations. However, in the case of multi-species Dark Matter models where not all the species are the same mass, one must distinguish the mass scale in the definition of x from the mass in the definition of neq [38].

III. ANALYTIC SOLUTIONS
In this section, we take an analytic approach to solving equation 12 in terms of ρ and w. We start by defining ∆ = Y − Y eq and combine it with Eq. 11: For small values of x, ∆ and d∆ dx are small as Y tracks Y eq closely. Using these approximations, as well as equations 3 and 6, solving for ∆ gives: The criterion for freezeout is given by ∆(x f ) = cY eq for some c. This c will define the freezeout temperature, the point at which we change from the high temperature limit solution to the low temperature limit solution. We give more details about the choice of c in section IV where we choose c to best fit the numerical results. Using the early time solution, we can obtain an implicit equation for x f : where H is given in Eq 1,ẋ is given in Eq 4, and n eq is given in Eq 6. Without knowing the functional forms of ρ and w, this is as far as we can go in general. Even if Eq. 21 is a complicated function, it can always be solved numerically for a value of x f much quicker than solving 12. Once x f is found, we can look at 11 in the large x case where Y ≈ ∆ and Y eq ≈ 0: Integrating from x = x f to x → ∞ gives: where Y eq = n eq /s, s is given in Eq 2,ẋ is given in Eq 4, and n eq is given in Eq 6. It is convenient to neglect the contribution from (cY eq (x f )) −1 as doing so ensures that Y (x → ∞) is a strictly increasing function of c. Not ignoring this term can lead to numerical issues such as having two values for c orders of magnitude apart which give the correct relic density, or there being no values of c which give the correct relic density. This gives: Using Eq. 24, we can find the final relic density: where is the entropy density today and is the critical energy density today. Here, H 0 denotes the value of the Hubble parameter today and is usually given by H 0 = 100h km/s/Mpc ≈ 70 km/s/Mpc [7].

A. Solving the Boltzmann equation numerically
In this section, we will investigate numerical solutions to equation 12. To do so, we will assume that the Dark Matter is composed of scalars and that the thermally averaged cross section is approximately constant with σv = σ 0 .
We will assume that the Dark Matter was in thermal equilibrium with the Standard Model well before it froze out. 3 The parameters that we are interested in varying are the mass of the Dark Matter (m), the value of η, and the thermally averaged cross section σ 0 , as well as investigating the cases of n D = 6, 8. The n D = 6 scenario corresponds to some quintessence models with a kination phase while the n D = 8 corresponds to brane world cosmology or some late inflaton decay models [36]. We are not interested in changing x 0 because any such change can be absorbed into a change in η. For this numerical analysis, we took x 0 = 25. The required thermally averaged cross section to obtain the observed relic density ΩCDMh 2 = 0.120 as a function of mass. The black, red, orange, green, blue, and purple (from the top to the bottom) curves represent η = 0, 1, 10, 10 2 , 10 3 , 10 4 , respectively. The solid colored curves represent nD = 6 while the dashed curves represent nD = 8. These coincide for the black (top) curve.
In figure 1, we show the required thermally averaged cross section to obtain the observed relic density Ω CDM h 2 = 0.120 for various values of η. First, we should comment on the shape of the η = 0 curve, namely why it starts by increasing, followed by a sharp decrease, followed by another increase. It should be noted that the typical value of x for freeze-out is in the range x = 10 to x = 30, so the important temperatures to consider are on the order of m/10. Since Ω CDM ∼ mn where m is the Dark matter mass and n is the Dark Matter number density, we see that the general trend for increasing the mass should be to decrease the number density n by increasing the thermally averaged cross section. In the low mass and high mass regimes, this is what happens since g * (x) does not change significantly for the important values of x. However, around m = 50 to m = 200 GeV, the Dark Matter is freezing out around the same time as the QCD phase transition, resulting in a rapidly changing g * (x). The increased numbers of degrees of freedom as we increase the mass results in requiring a decrease in the thermally averaged cross section to obtain the correct relic density. The derivative term 1 3 d log g * d log x is small at all temperatures and does not significantly influence the above trends. Although it is usually ignored, we include it in our analysis.
As expected from equation 12, we see from figure 1 that increasing the value of η results in requiring a smaller thermally averaged cross section. Even modest values of η can have large effects on this value. We also see that the n D = 8 case always results in needing a smaller thermally averaged cross section.
In figure 2, we show the required thermally averaged cross section as a function of η to obtain the correct relic density Ω CDM h 2 = 0.120 for various masses. One thing to notice is that the masses do not follow any nice pattern in terms of where they start off for low η. This is understood by imagining picking points from the η = 0 curve in figure 1 at various masses and projecting them onto the σ 0 axis. Because of the irregular shape of the curve, the starting points for the masses will also be irregular.
One thing we do see is that at very small η, there is no effect as expected, but at around η = 0.1 for n D = 6 and η = 0.01 for n D = 8, the effects of the modified expansion start to be seen. Between η = 0.1 and η = 10 for both cases, we see that the required thermally averaged cross section starts to decreases according to a power law. From equation 18, we can deduce that this power law will have the form Figure 1 and 2 give us an idea of how the cross section needs to change to get the correct relic density. In figure  3, we show a full map which includes both effects. The various colors represent the value of log 10 σ 0 . Here, we see that the cross section gets much smaller as η increases and does not differ much as the mass increases. At first glance, it may be difficult to see any obvious differences between the n D = 6 case on the top and the n D = 8 case on the bottom. In figure 4, we plot the ratio σ (6) 0 /σ (8) 0 where the superscript indicates the value of n D . This shows us that the differences are more important at low mass and only differs by a factor of order 1. At this point it should be noted that if the cross section required to obtain the correct relic density becomes too small, then it will bring the thermal equilibrium assumption into question. Although we do not worry about that in our case, if one was to consider much larger values of η, one would need to keep this in mind. Furthermore, our analysis will only hold if there is no entropy injection at the end of the new cosmological era defined by the dominance of the new energy density. If instead of a rapidly cooling energy density with n D > 4 we choose, for example, a new matter dominated era with n D = 3, the new energy density would need to decay away which would create an entropy injection diluting the DM. Further details about these scenarios can be found in [27][28][29][30].

B. Comparison with analytic approach
When deriving the analytic approach in section III, the value of c was ambiguous. The standard lore suggests that if σv ≈ σ 0 x −n for x 3, then one should choose c such that c(c + 2) = n + 1 [7]. The value of n is related to the velocity dependence of the cross section and, in We found that if we make fewer simplifications while using this approach, such as using Eq. 21 to find x f and Eq. 24 to find Y (x → ∞) without neglecting terms, we need much different values of c from the typical c(c+2) = n + 1 to obtain the same relic density. See figure 5 for the values of c needed to obtain the correct relic density using the cross section values obtained in figure 3. We see that the value of c changes drastically when we change η but does not change much as the mass changes. Using this information, we can fit for c to obtain the simple Using this approximation, we can calculate the relic density and see how it compares to the actual value. In figure  6, we plot This approach works well for most of the parameter space explored but can be off by a factor of order one in some places. In particular at low η and low mass, the analytic approximation tends to be too small while at high eta and low mass, the result n D = 8 tends to be too large. This does allow one to have the correct order of magnitude for the relic density which, depending on the purpose, can be good enough. In the range considered, the biggest deviations in the n D = 6 case were Ωh 2 = 0.0715 on the low end and Ωh 2 = 0.161 on the high end. For the n D = 8 case, the biggest deviations were Ωh 2 = 0.0812 on the low end and Ωh 2 = 0.255 on the high end.
One may be concerned about the large range of c obtained and the validity of the approximations made in section III. Recall that in the small x regime, we assumed that ∆ = Y − Y eq was small compared to Y eq and for large x, we assume that Y eq is much smaller than ∆. Inevitably, when we go from one to the other, the assumption will break down and both values will be comparable. Where to choose this point becomes a question of what best fits the data.

V. SUMMARY AND CONCLUSIONS
In this work, we derived a general Boltzmann equation which can be used for any modified expansion history of the universe. This was done by characterizing the energy density as a general function and letting w from the equation of state vary with temperature. We then looked at the specific example where we could write the energy density as a sum of components where all the values of w from the equations of state for each component were constant. Finally, we looked at the case of a single extra energy density component parameterized by how fast it cools, n D , and how abundant it was compared to radiation, η. Using these results, we found an analytic approach to approximating the relic density.
We then solved the modified Boltzmann equation for this last case to compare the effects to the standard picture. We found that for n D = 6 and n D = 8, the required cross section to obtain the measured relic does not changes significantly for η 10 −1 and follows a power law for η 10 with some transition region between them. Importantly, a larger value of η required a smaller cross section. This is quite simple to understand if we were to imagine what the expansion history of the universe would be if we ignored radiation. In this hypothetical situation, the scale factor would go as a ∝ t 2/n D , which gives a Hubble expansion rate of H = 2t/n D , meaning that at the same temperature, the universe would be expanding slower than the radiation dominated case where H = t/2. Because of this, it takes more time for the rate of expansion to be comparable to the rate of the reaction keeping the Dark Matter in equilibrium, so freezeout would occur later. A smaller cross section for the equilibrium reaction, which normally results in the Dark Matter freezing out earlier, is required to balance this effect.
We also showed that in the n D = 8 case, the cross section needed is comparable to the n D = 6 case. The only significant difference being in the low mass, high η case where they start to differ by a factor of order 1.
This result has some important consequences for model building and Dark Matter detection. If one of these modified cosmological scenarios occurred in our Universe and Dark Matter went through freezeout, then it may well be the case that at weak scale masses, the cross sections are orders of magnitude smaller than current bounds. However, this opens the possibility that heavier Dark Matter particles underwent freezeout. The bound of the mass of Dark Matter undergoing freezeout is about 100 TeV from perturbative unitarity [6], but this bound could be relaxed in these scenarios. This again poses a problem for Dark Matter detection as they are not designed to look for Dark Matter at these masses. All in all, these results suggests that Dark Matter detection may be a bigger challenge than previously expected if there is a significant changes to the expansion rate of the early Universe.

ACKNOWLEDGMENTS
The author thanks Stephen Godfrey, Catarina Cosme, and Mara Dutra for guidance and helpful conversation. The author also thanks Keith Dienes and Brooks Thomas for introducing him to the subject of Dark Matter and freezeout. This work was supported by the Natural Sciences and Engineering Research Council of Canada.