Variations on the vev flip-flop: instantaneous freeze-out and decaying dark matter

In this work we consider a simple model for dark matter and identify regions of parameter space where the relic abundance is set via kinematic thresholds, which open and close due to thermal effects. We discuss instantaneous freeze-out, where dark matter suddenly freezes-out when the channel connecting dark matter to the thermal bath closes, and decaying dark matter, where dark matter freezes-out while relativistic and later decays when a kinematic threshold temporarily opens. These mechanisms can occur in the vicinity of a one-step or a two-step phase transition. In all cases thermal effects provide this dynamic behaviour, while ensuring that dark matter remains stable until the present day.


Introduction
Understanding the nature of dark matter (DM) is one of the most pressing and elusive problems in physics. For a long time, the best-motivated candidate has been the Weakly Interacting Massive Particle (WIMP). In this picture, DM particles are hypothesised to be heavy (m DM 100 GeV) and to have weak but non-negligible interactions with Standard Model (SM) particles. In the very early universe, these interactions keep the DM and the SM particles in thermal and chemical equilibrium, until eventually the expansion of the universe dilutes the DM density to the extent that DM annihilation into SM particles ceases, and a relic abundance survives to the present day (a process known as freeze-out). However, in spite of a vigorous, multi-pronged search programme [1][2][3][4][5][6], there is still no evidence for WIMPs, and the remaining allowed parameter space has been drastically reduced.
There is currently a large effort studying alternative mechanisms of DM production [7]. In the early universe, the SM particles formed a hot, thermal bath and the effects of this plasma can have a dramatic impact on particle properties and interactions. These finitetemperature effects are known to cause important phenomena in the development of the universe, such as the electroweak phase transition (EWPT). Recently there has been JHEP05(2019)070  interest in understanding the impact of finite temperature effects in mechanisms of DM production [8][9][10][11][12][13][14][15][16][17].
In the present work we examine the mechanism presented in [12,15] where the DM abundance is set not via freeze-out but via a period of DM decay, which occurs when a new scalar field temporarily obtains a vacuum expectation value (vev) just before the EWPT. We find that the key feature required for this period of decay is the opening and closing of kinematic thresholds, due to particle masses obtaining temperature dependent corrections. In the early universe, bosons receive large finite temperature corrections and their masses and vevs should be treated as functions of temperature. In section 2, we propose a simple model which retains the key features of the mechanism in [12]. With this simple model, we identify further regions of parameter space where interesting DM production mechanisms can appear. We first discuss the effective potential at finite temperatures in section 3, the thermal bath in section 4 and Boltzmann equations to track the DM abundance in section 5. With this framework, we describe instantaneous freeze-out in section 6, where DM freezes-out abruptly when a kinematic channel closes. In section 7 we show that a period of DM decay can set the relic abundance, dubbed decaying dark matter, when a new scalar field obtains a vev in a one-step phase transition, while in section 8 we look at a situation similar to that described in [12] where a two-step phase transition (or a vev flip-flop) can lead to a period of decaying dark matter. Finally, in section 9 we briefly consider the experimental constraints on the different scenarios.

The model
To demonstrate the effects as clearly as possible, we consider a simple model consisting of the SM plus a real scalar, S, and two dark sector Dirac fermions, χ and ψ, shown in table 1. In what follows, χ will constitute the dominant component of dark matter. These particles will all be gauge singlets, and will have masses in the GeV to TeV range.
The Lagrangian for these fields is  Table 2. A summary of the mechanisms that operate in various regions of parameter space. In the first column we give the nature of the phase transition of the new scalar S. In the penultimate column we note whether the χ decay and inverse-decay channel closes due to the vev of S becoming sufficiently large, or to the vev of S becoming zero after a 'vev flip-flop'.
In our notation, the SM Higgs potential has µ H 88 GeV and λ H 0.12, with H = (G + , (h + iG 0 )/ √ 2). We choose a basis where the CP-even, neutral component h obtains the vev. We see in eq. (2.2) that the only interaction at dimension-4 between these new particles and the Standard Model is the Higgs portal term (λ p /2) S 2 (H † H). Although in principle there may also be the terms S 3 and S(H † H), we assume these couplings to be negligibly small. In sections 6 and 7 we will consider λ p 1. In this regime the scalar fields decouple and we can consider the evolution of S alone, simplifying the picture. In section 8 we will consider λ p 1. For the new Yukawa couplings we will typically consider y ψ 1, y χψ ∼ 10 −7 and y χ 0. A large y ψ means that ψ will remain in thermal equilibrium throughout the processes affecting the χ abundance. This small y χψ will however make direct and indirect detection of dark matter challenging. Within this model, none of the operators with small couplings will be generated at a significant level via loops.
We will focus on the regions of parameter space where finite temperature effects from the scalar sector have a significant impact on the final abundance of χ. We know from the standard WIMP picture that the observed relic abundance is obtained if DM freezes-out when its mass is 20-30 times greater than the freeze-out temperature. Since we do not have abundances below the equilibrium abundance, this provides a lower limit on m χ . If the DM mass is much greater than the temperature at which the abundance is set, then finite temperature effects will be small compared to zero temperature effects. This tension means that for the mechanisms presented in this paper, we will generally focus on the parameter space where m χ ,m ψ 30 µ S . We choose the signs of the Yukawas so that for positive couplings, S = 0 will give a positive contribution to the mass. In table 2 we summarise the different mechanisms and highlight the relevant parameter space.
Finally we emphasise that many models can accommodate degenerate or nearlydegenerate vector-like fermions, while the hierarchy between the masses of the new fermions and the new scalar may be expected since their mass scales are not directly connected. Similarly, while the strong hierarchy of the new Yukawa couplings y ψ , y χψ and y χ is not explained, its origin could be related to the SM flavour puzzle. Although we focus on this particular region of parameter space to highlight the importance of finite temperature effects, which will prove to dramatically alter the DM relic abundance, we note that these effects may still be numerically important in much wider regions of parameter space.

JHEP05(2019)070
3 The effective potential at finite temperatures We first turn to the effective potential of the scalar fields, and describe the finite temperature corrections we include. The effective potential is analogous to the free-energy of a system, and the system (in this case the vacuum scalar field configuration) will move to the minimum of the effective potential. The leading terms in the zero-temperature effective potential are the tree-level potential V tree (h, S), the one-loop Coleman-Weinberg correction V CW (h, S) [18], and a one-loop counterterm V CT (h, S). The leading corrections at finite temperature are the one-loop thermal corrections V T (h, S, T ) [19] and the resummed higher order "daisy" diagrams [20][21][22]. The one-loop effective potential at finite temperature, V eff (h, S, T ), is then The tree level potential V tree is given in eq. (2.2). The T -independent Coleman-Weinberg contribution is [18,21] where the sum is over the eigenvalues of the field-dependent mass matrices of all fields which couple significantly to the scalars and |n i | accounts for their respective numbers of degrees of freedom. In our case n h = n G 0 = n S = 1, n G + = 2, n Z = 3, n W = 6 and n t = −12. We do not include the lighter quarks as they couple only very weakly to the SM Higgs. The coefficient n i is positive for bosons and negative for fermions, due to their different statistics. As a renormalisation scale Λ we use the characteristic scale of the S field, µ S . Using the dimensional regularisation scheme, C i = 5/6 for gauge bosons and C i = 3/2 for scalars and fermions. The field-dependent masses of the CP even neutral scalars, in the basis (h, s), are while the field dependent masses of the remaining bosons are where g and g are the SM u(1) Y and su(2) L coupling constants, respectively. Although the dark sector fermion ψ couples strongly to S, we neglect its contribution since its mass is always much larger than either the temperatures or the scalar field values we consider [23].

JHEP05(2019)070
To ensure that H = µ H / √ λ H , m 2 h = 2µ 2 H and that m S is given by its tree level value at T = 0, we add the counterterm where the factors δ i are (3.12) The one-loop finite temperature correction is [19] V where again we sum over the same eigenvalues as in eq. (3.2). The positive sign in the integrand is for fermions and the negative sign is for bosons. Finally, higher order diagrams containing the bosons give rise to the so-called "daisy" corrections [20][21][22]24]   where the first term should be interpreted as the i-th eigenvalue of the matrix-valued quantity [m 2 (h, S) + Π(T )] 3/2 . Here, m 2 (h, S) is the block-diagonal matrix composed of the individual mass matrices [25]. This sum runs only over the bosonic degrees of freedom. The thermal (Debye) masses [20] in eq. (3.14) are

One-step phase transition
With the effective potential in hand, we can now consider its behaviour in various regions of parameter space. We will first consider the regime where λ p 1, where the new scalar S and the SM Higgs are weakly coupled. In this case the effective potential in h and S decouples, V eff (h, S, T ) = V eff (h, T ) + V eff (S, T ). To consider the evolution of the new scalar, it is sufficient to consider the effective potential as a function of S alone.
In figure 1 (left) we show V eff (S, T ) at several temperatures, for a particular choice of parameters. We will be interested in µ S > 0, so that S obtains a vev at T = 0. We see that at high temperatures, thermal corrections dominate the effective potential and the minimum is at S = 0, so there is no vev. As the universe expands, the temperature reduces until T µ S , at which point the finite T corrections become similar in size to the T = 0 potential. At T c there is a second-order phase transition and S obtains a vev. As the universe cools further, the minima deepen and the vev increases to around 1.7 GeV.
As well as the position of the vev, the effective potential also determines the physical mass of the scalar S in the thermal bath. Finite temperature corrections at one-loop are taken into account by taking the second derivative of the effective potential at a minimum. In figure 1 (right) we show the evolution of the vev along with the physical S mass, as a function of temperature. We see that the mass is large at high temperatures, becomes small through the phase transition (at T c the second derivative of the effective potential goes to zero), and reaches a value µ S at zero temperature.
In sections 6 and 7 we will show that this temperature dependence of the mass and vev can lead to kinematic thresholds opening and closing, allowing dark matter to come into equilibrium or decay in some temperature window.

Two-step phase transition
In the previous section we assumed the new scalar and the SM Higgs to be weakly coupled. If instead they are coupled with λ p ≈ 1, we must consider the effective potential as a function of both S and h. We are in particular interested in the region of parameter space which exhibits a two-step phase transition (also called a vev flip-flop). In figure 2 we show an example of such a transition. If µ S ≈ µ H , it can happen that at high temperatures (top-left) there is one minimum at (h, S) = (0, 0) and neither field has a vev. As the temperature reduces, the finite temperature corrections reduce and minima develop at S = 0 due to the T = 0 potential. As the temperature drops further, further minima appear at h = 0, (top-right). At this point there is a barrier between the S = 0 minima and the h = 0 minima, and there is a period of supercooling while the universe remains in this meta-stable phase. As the temperature drops further, a first order phase transition may take place, when the formation and growth of bubbles of the new phase is energetically favourable [26]. We calculate the temperature of the phase transition (here named the nucleation temperature) using the publicly available code cosmoTransitions [27][28][29][30]. Note that away from the minima the effective potential is gauge dependent, so the nucleation temperature in principle has a residual gauge dependence [25], which we neglect. At the nucleation temperature T n , the universe passes to the phase where h = 0 and S = 0 (bottom-left). As the temperature reduces further, these minima deepen and the universe ends in a phase with h = 246 GeV and S = 0 (bottom-right).
In figure 3 we show the scalar masses and their vevs as a function of temperature. Since at no point do both h and S obtain a vev, there is no mixing between them. As for the one-step phase transition, section 3.1, the scalars have masses similar to T at high temperature. In this example, S obtains its vev in a second order phase transition at T ≈ 230 GeV. As in the one-step phase transition, the mass becomes very small during this transition. After this first transition, S obtains a vev, which grows as the temperature reduces. The mass of S starts to grow while the Higgs continues to get lighter. This continues until the nucleation temperature T n = 116 GeV, when the universe passes to the phase where S = 0 and h = 0 in a strongly first-order phase transition. After this second phase transition, which breaks electroweak symmetry, three components of the Higgs doublet are eaten by the gauge bosons, and the mass of the remaining scalar grows until it reaches 125 GeV at T = 0. Similarly, the h vev grows to 246 GeV at T = 0.

Phase diagram
Now that we have outlined the two main scenarios of interest, we briefly discuss the parameter space of the effective potential. In figure 4 (left) we show the phase diagram of the tree-level potential for λ S = 3. The parameter space is divided into regions where at zero temperature the global minimum is one with (green) ( h = 0, S = 0), (orange) ( h = 0, S = 0) and (blue) ( h = 0, S = 0). The blue region does not correspond to our universe, since electroweak symmetry is broken, but either the green or orange region may be physical. In these regions we also plot contours of the S mass at T = 0. The orange region shows where a one-step phase transition occurs, discussed in section 3. for the physical (h, 0) and (h, S) minima. In (0, S) electroweak symmetry is not broken, so it is not physical. Right: the region where a two-step phase transition occurs (green) along with regions where S never obtains a vev (red) and where electroweak symmetry is not broken at T = 0 (blue) as a function of m S (T = 0) and λ p . Superimposed on both in dark green is the region where cosmoTransitions gives a two-step phase transition for the one-loop effective potential. simplify calculations we later assume λ p 1, but the mechanisms we discuss would be seen in the whole orange region.
A two-step phase transition, discussed in section 3.2, may occur in the green region. The pixellated region shows points where cosmoTransitions finds that S obtains a vev at some temperature, and then at a later lower temperature, as the SM Higgs obtains its vev, the S vev goes back to zero (also called a vev flip-flop) for the one-loop effective potential. Here, the vev of h gives a contribution to the mass of S and the relation between m S (T = 0) and µ S is We see that when, e.g., m S (T = 0) = 150 GeV, only certain values of λ p are allowed. This is seen to correspond to figure 4 (right), where we again plot the green region, now on the m S (T = 0)-λ p plane. In the lower-right red region, which shows where µ 2 S < 0 at treelevel, the S symmetry is never broken and S never obtains a vev. In the blue region the deepest minima occur when S has a vev and h has no vev, so electroweak symmetry is not broken and this region does not correspond to our universe. The pixellated region again shows where cosmoTransitions calculates that the two-step phase transition successfully completes for the one-loop effective potential. As λ S increases, the two-step phase transition region becomes larger. For a given m S (T = 0), the S vev and the temperature at which S first obtains a vev both increase with λ S . Although not shown in these diagrams, S obtains a vev at higher T for larger portal couplings.

JHEP05(2019)070 4 Finite temperature corrections and the thermal bath
We now turn to a discussion of the thermal bath. In the hot, early Universe, there is enough energy to produce particles with masses T . If there are efficient processes which create and destroy a particle, it will come into equilibrium with a number density n ∼ T 3 . In this work, we will assume that after inflation and above T 100 TeV, there were efficient processes which led to χ, ψ and S to come into equilibrium with the thermal bath. For freeze-in scenarios where χ does not thermalise in the early universe, see [14].
To determine if certain processes are efficient at keeping the new particles in contact with the thermal bath as the temperature cools, we must calculate the rates of these processes and compare them to the rate of Hubble expansion. The masses which enter the Feynman rules are given by the imaginary part of the self energies, which are modified in the presence of thermal corrections [31,32]. Although both bosons and fermions obtain these thermal corrections, their different statistics lead to different boundary conditions in the compactified dimension at finite temperature. The bosonic contributions contain a Matsubara zero-mode while the fermionic contributions do not, so the fermionic contributions are subleading to the bosonic contributions. We therefore only include the thermal corrections to the boson masses. However, from eq. (2.1), we can see that when S has a vev, the effective mass parameter for ψ will receive an extra contribution, so m ψ will still depend on temperature. Since we imagine y χ to be small, we will take m χ to be temperature independent. Since we will mostly be interested in small mass splittings, we introduce the parameter These temperature dependent masses mean that kinematic thresholds can open or close as the temperature reduces. In this work, we focus on scenarios where decay and inverse decay of χ, which is only allowed when m χ > m ψ (T ) + m S (T ), has a dramatic impact on the resulting dark matter relic abundance. In section 6 we consider a scenario where χ remains in equilibrium until the threshold closes, at which point there are suddenly no processes which keep χ in contact with the thermal bath and it immediately freezes-out. We call this process instantaneous freeze-out. In this case, the relic abundance is set by the temperature at which the threshold closes, and m χ 30 T End results in the observed relic abundance. In section 7 we consider a scenario where the kinematic threshold is closed at high temperatures and so χ freezes-out when relativistic. Then, as the finite temperature corrections to m S become smaller, the threshold opens and χ can decay for some time. Finally, the mass of ψ increases due to S obtaining a vev, closing the threshold and stabilising χ, which we call decaying dark matter. The final relic abundance will then be determined by the amount of decay that has occurred when the threshold closes. The abundance of χ will approach but not reach equilibrium and produce the observed relic abundance. Ensuring that ψ and S remain in equilibrium means that χ → ψS efficiently depletes the abundance of χ, with the energy density being passed to the thermal bath. In section 8 we consider a similar scenario, but where S is strongly coupled the SM Higgs and a two-step phase transition can occur. In this case, the mass of ψ reduces when S obtains a vev and the kinematic threshold opens, allowing χ to decay, and closes when the S vev disappears. Again the abundance of χ will approach but not reach equilibrium and can produce the observed relic abundance. We call this the vev flip-flop.
We now consider the constituents of the thermal bath. The new particles we introduce do not change the evolution of the SM particles, which follow the standard cosmology. For the new particles, we will ensure that S and ψ remain in equilibrium when number changing processes of χ are active. Although this is not always necessary it will simplify our calculations. We will discuss the precise processes which keep S and ψ in equilibrium in each section, and here give an overview.
In figure 5 we show three possible ways that S interacts with the SM bath. The processes in figure 5 (a) and (b) can only keep S in thermal equilibrium when electroweak symmetry (EWS) is broken, while process (c) is always active. Process (a) and (b) will mostly be relevant in sections 6 and 7. Process (c) will keep S in thermal equilibrium in section 8, where EWS is not broken during χ decay. The annihilation cross sections and decay rates for the processes in figure 5 are where C f is a colour factor, y f = m f √ 2/v is the Yukawa coupling of the SM fermions to the Higgs, v is the SM Higgs vev, and the factors F γγ W and F γγ f are given in [33,34]. In figure 6 we show processes that may keep ψ in contact with the thermal bath. Since we will be interested in situations where χ is out of equilibrium, we will choose values of y ψ so that ψ can remain in equilibrium through figure 6 (b) until the abundance of χ is fixed. We also want χ to be the dominant dark matter relic, so we will choose y ψ large enough that ψ freezes-out with less than 10% to the relic abundance. To ensure this, we will need y ψ 1 in sections 6 and 7 and |y ψ | 4 in section 8. This large y ψ also has the benefit of providing a large mass contribution to ψ, when S obtains a vev. Process (c) goes as y 4 χψ and will be too small to keep ψ in equilibrium. In what follows, we will see that m ψ (T ) is always smaller than m χ + m S (T ), so the decay process in (a) will always be kinematically forbidden and ψ will only decay via an off-shell S. The annihilation cross section for the processes in figure 6 where we have taken the limit m S → 0. We keep the full m S dependence in our numerical work. Finally, χ may be in equilibrium through the processes shown in figure 7. The decay and inverse decay diagram, figure 7 (a), is proportional to y 2 χψ , while the processes figure 7 (b) and (c) are at least proportional to y 4 χψ . We will choose y χψ small enough that the processes proportional to y 4 χψ do not come into equilibrium at the temperatures of interest. The decay rate of the process in figure 7 (a) is

JHEP05(2019)070
Note that we do not show χψ ↔ S. Although this channel is kinematically open in the very early universe, when m S (T ) > m χ + m ψ (T ), we always take m χ ,m ψ ≈ (30 − 100) m S (T = 0), so this channel will be closed at the temperatures where the abundance of χ is being set. We also do not show ψ ↔ χS since for the parameter points we consider this decay is always kinematically forbidden. In (d)-(f) we show channels contributing to χS ↔ ψS. The cross section for the process is where we have taken the limit m ψ m χ and m S 0. We have also set v S = 0 since (f) is subdominant to (d) and (e) for the parameters we consider. In our numerical analysis we use the full expressions. We do however always take Γ ψ = Γ S = 0 since the widths of ψ and S are very small. At temperatures near m S but lower than m χ , n S σ(Sχ → Sψ)v ∼ T 3 /m 2 χ ∼ m χ /x 3 where x = m χ /T , hence the rate of this process will generally be suppressed by a factor x 3 compared to Γ(χ → ψS) ∼ m χ . However, the rate can be resonantly enhanced when m S (T ) ≈ 0.
Processes such as χψ → S → SSS and χψ → S → SS (when S has a vev) are also possible, although the rates will be small as the abundance of ψ will be Boltzmann suppressed at the temperatures of interest, and the intermediate S will be a long way offshell. We note that although there is in principle mixing between χ and ψ when S = 0, this mixing will be small and we will ignore its negligible effects. We also neglect diagrams involving λ p . This coupling is taken to be small in section 6 and section 7, while in section 8 these processes lead to a rate significantly smaller than the Hubble rate, due to the x 3 suppression mentioned above.

Dark matter abundance and the Boltzmann equations
We will be interested in the abundance of χ whether or not it is in equilibrium so we will keep track of its abundance using Boltzmann equations. In general the Boltzmann equation for χ isṅ In the parameter space we consider the collision term will have non-negligible contributions from the decay process χ ↔ ψS and the scattering process χS ↔ ψS, The collision term for the decay process χ ↔ ψS is given by

JHEP05(2019)070
where dΠ i is the Lorentz invariant phase space of particle i and f i are their phase space densities. We neglect Pauli-Blocking and Bose-Enhancement, i.e., (1 ± f i ) ≈ 1, and assume Maxwell-Boltzmann distributions for ψ and S. This is a good approximation for ψ, as we will be interested in temperatures much less than the mass of ψ, but since these temperatures are similar to the mass of S, the full S distribution should be used for more detailed calculations. As we argued above, ψ and S will remain in thermal equilibrium, and we assume with zero chemical potential, during the periods of interest. These assumptions are useful as, when combined with the energy-conservation part of the delta-function, it follows that Assuming that the decay process is CP invariant, the collision term then becomes where g χ = 2 is the number of degrees of freedom of χ (χ will satisfy an analogous Boltzmann equation, and will also contribute towards Ω χ ). Note that there is a factor m χ /E = 1/γ. It appears since particles with higher momenta experience increased timedilation in the rest frame of the plasma, and so their lifetime is increased. Since we will be interested in situations where almost all χ decay, this time-dilation in the tail gives an important contribution. To solve this equation in practice, we discretise the integral into bins in momentum space and track the number density in each bin using which gives where For the key results we have checked that changing the number of bins and the highest momentum considered does not significantly change the resulting abundance. We follow an analogous derivation for the collision term of the scattering process χS ↔ ψS, yielding C i χS↔ψS [n i χ ] = −n eq S σv (n i χ − n eq,i χ ) , (5.9)

JHEP05(2019)070
where the thermally averaged cross section for non-identical initial state particles is given by where K 1,2 are modified Bessel functions of the first and second kind, see, for example, [35].
Although we do not need to discretise this collision term in momentum space, we write it in this way to combine it with the collision term for decay. Writing Γ χS→ψS = n eq S σv (5.11) and introducing the yield, Y χ = n χ /s, we can write where x = m χ /T . We will be interested in tracking the relic abundance through changes in the effective number of relativistic degrees of freedom, g * (T ), so we use the conservation of entropy to write We then obtain a differential equation for each momentum bin, Although χ is at times out of equilibrium we check that it never dominates the energy density of the universe, so we can use the usual radiation dominated Hubble expansion of the universe.

Instantaneous freeze-out
With this machinery, we can now turn to exploring various mechanisms of dark matter production. In this section, we will focus on the region of parameter space where the new scalar field exhibits a one-step phase transition, section 3.1, and the Yukawa coupling y χψ ∼ 10 −7 is large enough to put χ into equilibrium. We will assume m χ ≈ m ψ ≈ (30-100)×m S (T = 0). This both produces the observed relic abundance of χ and means that at T ≈ µ S both initial state particles in the 2-to-2 processes χχ → ψψ, χχ → SS and χψ → SS are significantly Boltzmann suppressed, dramatically reducing the rate n σv , so these processes can not keep χ in equilibrium. We will further assume that λ p = 10 −5 -10 −3 and y ψ = 2. These couplings ensure that S and ψ remain in equilibrium until the abundance of χ is set. S stays in equilibrium through SS ↔ HH when it is much heavier than the Higgs, and through SS ↔f f , where f is a SM fermion, when it is much lighter than the Higgs and electroweak symmetry is broken at the temperatures of interest. For S lighter than the muon, the Yukawa coupling becomes too small to keep S in equilibrium, putting a lower limit on the S masses we consider. The process S ↔ γγ can also keep S in equilibrium when S has a vev, near T E . The dark sector fermion ψ remains in equilibrium through ψψ ↔ SS.
This mechanism, and the others we present in this paper, depend heavily on the opening and closing of a kinematic threshold, which occurs due to the temperature dependence of the masses of S and ψ. As such, we will be interested in the region of parameter space where χ and ψ have a similar mass, as is natural if they obtain their masses from the same mechanism. In this example we will choose m χ = 77 GeV,m ψ = 74 GeV. We take y χψ = 5 × 10 −8 for the new Yukawa coupling while here and elsewhere we will take y χ ≈ 0. For the new scalar we will first choose m S (T = 0) = 1 GeV and λ S = 1. We can see in figure 1 (right) that with these parameters, S (T = 0) ≈ 1.7 GeV, so m ψ will increase from 74 GeV at high temperatures to 74 GeV + y ψ S (T = 0) ≈ 77.4 GeV at T ≈ 0. At high temperatures, ψ is the lightest dark sector particle, while at T ≈ 0, χ takes this role.
The temperature dependent masses for this set of parameters are shown in figure 8 (left). We see that at high temperatures, m S receives a large correction and the channel χ ↔ ψS is closed. As the temperature reduces, the T dependent corrections become smaller and the mass of S reduces. The channel χ ↔ ψS becomes active at T ≈ 17 GeV. This continues while S goes through its phase transition, where S obtains a vev and its mass starts to increase. The vev contributes to m ψ , so both of these effects act to close the channel. At T E ≈ 3 GeV, m χ becomes smaller than m ψ (T ) + m S (T ), so the process JHEP05(2019)070 χ ↔ ψS becomes kinematically forbidden. At T = 0, m χ < m ψ , so χ is the lightest dark sector particle and cannot decay, even via off-shell processes.
In figure 8 (right) we show the rates of key processes along with the expansion rate of the universe, H, as a function of x = m χ /T . When the rate of a process is much larger than H then it will effectively keep χ or ψ in equilibrium, but when it is much smaller then it will not. We see that between T S to T E , χ ↔ ψS is the dominant process and will act to put χ into equilibrium. The process χS ↔ ψS may also have a rate larger than the expansion rate of the universe, and will contribute. The rate of χS ↔ ψS shows two resonant features, one at x ∼ 20, when m S ∼ 0 and the propagator ψ can be nearly on-shell, and another at x ∼ 40, when m ψ ∼ m χ and diagrams (d) and (e) in figure 7 destructively interfere. Although ψ is kept in chemical equilibrium via ψψ ↔ SS until T E , it decouples soon after at x ∼ 30. Once ψ has departed from chemical equilibrium, the relic abundance of χ cannot be reduced, since any ψ produced will simply decay back to χ. We also show the rate of χS ↔ ψS if finite temperature effects are ignored. We see that in this case the χS ↔ ψS rate is always smaller than the Hubble rate, so χ will freeze-out while relativistic. All other processes involving χ have rates orders of magnitude smaller than the Hubble rate.
In figure 9 (left) we show the evolution of the yield of χ. At small x there are no processes connecting χ to the thermal bath with a rate larger than the Hubble rate, so χ freezes-out. Then, at T S the rate of χ ↔ ψS (and later χS ↔ ψS) becomes larger than the Hubble rate and the χ yield comes back into equilibrium. At T E , however, these channel suddenly becomes inefficient and, since there are no longer any processes keeping χ in equilibrium, it instantaneously freezes-out. We see that in this scenario, the relic abundance of χ is set not by usual, smooth freeze-out, resulting from an interplay between a slowly varying annihilation or decay rate and the Hubble rate, but by a sudden and dramatic change in rates. If these finite temperature effects are ignored, the final abundance is overestimated by ≈ 10 orders of magnitude, since the rate of the most relevant process χS ↔ ψS remains below the Hubble rate throughout.
This means that the resulting yield is not a function of the coupling constant (provided the assumptions above are satisfied), but only of the mass of the dark matter candidate and T E . The final yield is simply given by which is only a function of m χ and T E (which in turn is a function ofm ψ , y ψ and the effective potential). We now consider the parameter space of this mechanism. Rather than considering the full 7-dimensional parameter space, we explore the impact of varying one or two parameters at a time. The final relic abundance is insensitive to y χψ as long as χ ↔ ψS is fast enough to keep χ in equilibrium (providing a lower bound on y χψ ) and as long as χS ↔ ψS does not significantly deplete χ after T E (providing an upper bound on y χψ ). For the benchmark values of the other parameters, this is satisfied for 3 × 10 −8 < y χψ < 1 × 10 −7 . As discussed above, the mechanism occurs while y ψ is large enough to keep ψ in equilibrium while the processes χ ↔ ψS and χS ↔ ψS occur. If ψ goes out of equilibrium before T E , the JHEP05(2019)070 Figure 9. Left: the evolution of the yield of χ as a function of temperature (solid) with the equilibrium yield of χ (dotted) and the yield calculated without accounting for the finite temperature effects (dashed). Also shown is the yield which is required to obtain the observed relic abundance of χ (black, dashed). Right: the required values of m S (T = 0) and λ S for various χ masses, for ∆ = 5% and 2%. In the left light-grey region S is not in equilibrium throughout the decay of χ, in the central grey region ψ decouples before T E , while in the right dark-grey region ψ is not a subdominant relic.
processes will not reduce the energy density in the dark sector and, once m χ < m ψ , ψ will decay to χ. In this variation of the mechanism, the relic abundance is set by the freeze-out of ψ, not via the dramatic change in the rates of χ ↔ ψS and χS ↔ ψS at T E . However, in both scenarios, the Yukawa coupling y χ should be small enough that the 2-to-2 process χχ ↔ SS is inefficient by T E or when ψ freezes-out.
In figure 9 (right) we show the mass of χ required to produce the observed relic abundance, for ∆ = 5% and 2%. We see that m χ is 30 − 100 times larger than m S (T = 0), depending on λ S . Although m S (T = 0) depends on µ S , the S vev also depends inversely on the square root of λ S . This means that as λ S becomes smaller, the S vev gets larger. The temperature of the phase transition, which is close to T E , also increases. These two effects mean that m χ has to increase to compensate. Reducing the mass splitting ∆ has the effect of requiring a larger m χ . This is because a smaller ∆ means that the important 1-to-2 process closes at a higher T , so Y Eq. χ (T = T E ) becomes larger. To yield the observed relic abundance, m χ then has to increase. If the mass splitting is much larger than 5% then the S contribution to m ψ at T = 0 is not large enough to close the mass gap, so m χ > m ψ (T = 0). As mentioned above, if S is lighter than the muon then it is not in equilibrium throughout the process, so the decay does not reduce the energy density of the dark sector. For m µ < m S (T = 0) 1 GeV, either λ p must be reduced below λ p = 10 −3 or there is some tuning between the µ 2 S term and λ p h S 2 /4 to achieve this low mass. If λ p is reduced, care must be taken that S remains in equilibrium throughout. We find that we require λ p 10 −5 in this low m S region of parameter space. On the other hand, as JHEP05(2019)070 m S (T = 0) increases, m ψ also increases (along with m χ ). For the parameters chosen, when ψ is heavy it can freeze-out (at T = T f.o. ψ ) at temperatures higher than T E . In this case, dark sector fermion number conservation implies that the χ abundance is given by With y ψ = 2, if m ψ 1 TeV then ψ will freeze-out with a relic abundance greater than 10% of the observed relic abundance, so the mechanism we are discussing will not be the dominant mechanism in setting the relic abundance. The quartic coupling λ S can become larger than 4 but at some point it becomes non-perturbative and the one-loop analysis breaks down.

Decaying dark matter with a one-step phase transition
Now we continue to consider the region of parameter space where the new scalar field exhibits a one-step phase transition, discussed in section 3.1, but assume a smaller Yukawa coupling, y χψ 10 −7 . With a coupling this small, χ either freezes-out while relativistic or never comes into equilibrium with the thermal bath (as in the freeze-in scenario [14]). Here we will assume that there is some UV physics which connected χ to the thermal bath after reheating, so that χ was in equilibrium at high temperatures. In this scenario χ will have frozen out when relativistic, then χ ↔ ψS, and to a lesser extent χS ↔ ψS, bring χ towards equilibrium after T S . This means that χ will deplete and its abundance will approach the equilibrium abundance. At T E these process become inefficient and the yield of χ stabilises. We will require m χ ≈ m ψ ≈ 30 m S to obtain the observed relic abundance. With these masses and couplings the rates of χ ↔ ψS and χS ↔ ψS are the only processes to ever have a rate greater than the Hubble rate. As discussed in section 6, ψ and S remain in equilibrium throughout the decay process.
To describe the mechanism in this case, we will keep the same scalar parameters chosen above but slightly increase the dark sector fermion masses (we take in this example m χ = 90 GeV andm ψ = 87 GeV) and reduce the Yukawa coupling to y χψ 10 −7 (here we take y χψ = 2.57 × 10 −8 ). The temperature dependence of the masses is shown in figure 10 (left). We see that with these parameter choices, figure 8 (left) is simply shifted to slightly higher masses and the channel χ ↔ ψS is open in the same temperature range, 17 GeV T 3 GeV.
In figure 10 (right) we show the evolution of the yield of χ with x. As expected for y χψ 10 −7 , χ freezes-out while it is still relativistic. However, at T S ≈ 17 GeV the channel χ ↔ ψS becomes kinematically allowed and the χ abundance reduces, approaching the equilibrium curve. The χS ↔ ψS channel also contributes to a lesser extent. However, in this case, the rate is not fast enough to reach equilibrium and Y χ remains larger than Y Eq. χ at all times. This continues until the channel closes at T E ≈ 3 GeV and the yield of χ becomes fixed until the present day, reproducing the observed relic abundance.
In this case the final yield is not set by Y Eq. χ at T E , but by the amount of χ that has decayed in the period from T S to T E . This depends on y χψ and on the temperature JHEP05(2019)070 Figure 10. Left: the masses of χ and ψ as a function of temperature, along with the period where the channel χ ↔ ψS is kinematically allowed. Right: the evolution of the yield of χ as a function of temperature (solid) with the equilibrium yield of χ (dotted) and the yield calculated without accounting for the finite temperature effects (dashed). Also shown is the yield which is required to obtain the observed relic abundance of χ (black, dashed). dependent masses of the particles. As such there is no simple formula for calculating the resulting abundance; the final abundance must be calculated numerically.
We now vary the parameters of the model to see their dependence. We first point out that the relic abundance obtained is exponentially sensitive to the coupling y χψ . In figure 11 (left) we show the evolution of the yield for different values of y χψ . As seen above, for a coupling of y χψ = 2.57 × 10 −8 we obtain Ω χ h 2 = 0.12. With a smaller coupling of y χψ = 2 × 10 −8 we obtain Ω χ h 2 ∼ 340, y χψ = 2.4 × 10 −8 gives Ω χ h 2 = 1.2, while y χψ = 3 × 10 −8 gives Ω χ h 2 = 0.015. We see that a ≈ 10% change in the coupling changes the final abundance by an order of magnitude. As such, this mechanism provides no answer to the coincidence problem (that Ω χ ≈ 5 Ω SM ).
In figure 11 (right) we vary the effective potential, m χ , y χψ and ∆, and plot the curve where the observed relic abundance is reproduced. We see that the observed abundance is reproduced for m χ ≈ (50-200) × m S (T = 0). As in section 6, a smaller quartic coupling λ S results in a larger S vev and phase transition temperature, resulting in a heavier χ. Increasing the mass of S similarly provides a larger vev and a larger phase transition temperature, requiring a heavier χ. A larger y χψ means that χ depletes faster so T E needs to occur at a higher temperature, which is achieved with a heavier S or a smaller λ S . Similarly, a larger ∆ means that T E is naturally lower, so this needs to be compensated with a larger S mass or a smaller λ S . When S is lighter than the muon, S must couple via a first generation SM Yukawa coupling, which is so small that S does not stay in equilibrium throughout the process. The upper bound on m S is again determined by requiring that ψ freezes-out as a subdominant relic, around m ψ ≈ 1 TeV, but we show just a small region of the possible parameter space to demonstrate the impact of varying y χψ and ∆. It may also happen, as in section 6, that ψ freezes-out before T E . This occurs for the dotted line (y χψ = 2.5 × 10 −8 ), since the lower m S (T = 0) reduces T E below the ψ decoupling temperature. For the dashed line (y χψ = 4 × 10 −8 ) the coupling is large enough that χ returns to equilibrium, so the situation is as discussed in section 6. The quartic coupling λ S can also in principle be increased to the perturbativity limit.

Decaying dark matter with a two-step phase transition
Finally we consider the region of parameter space where λ p is not negligible, but O(1). As discussed in section 3.2, this can lead to a two-step phase transition where first the new scalar S obtains a vev and later, at the onset of electroweak symmetry breaking (EWSB), the SM Higgs obtains a vev and S goes to zero. This two-step phase transition only occurs in the region of parameter space where µ S ≈ µ H , which means that m S (T = 0) = O(100 GeV). To maintain m χ 30 m S (T = 0), so that the abundance is set by decay, we require m χ , m ψ ≈ 4-5 TeV. Since ψ is now quite heavy, it will naturally freeze-out with a large relic abundance. To counter this and to ensure that ψ is a subdominant relic, we have to take |y ψ | 4 which, while large, is still perturbative. However, corrections from higher orders in the perturbative expansion will not be as small as usually assumed. With this two-step phase transition, we will take y ψ < 0, so that ψ becomes lighter when S obtains a vev. This allows the channel χ ↔ ψS to open when S obtains a vev, and to close when the S vev goes to zero. This situation can lead to the observed relic abundance either by instantaneous freeze-out when y χψ is ∼ 10 −7 , as in section 6, or with a period of DM decay with y χψ 10 −7 , as in section 7. Here we choose to explore the latter, to make a connection with the mechanism described in [12]. In [12] it was shown that a new su(2) L triplet scalar S 3 offers decay channels to a singlet DM candidate when there are two extra dark sector su(2) L triplet fermions, ψ 3 and ψ 3 . In this case, there are the channels χ ↔ ψ + 3 W − and χ ↔ ψ + 3 W − in addition to the process discussed in the present work, χ ↔ ψ 3 S 3 . However, as we demonstrate here, this extra structure is not required and the observed relic abundance can be obtained through the χ ↔ ψS process.
In figure 12 we again show (left) the masses as a function of T . Recall that m ψ (T ) = m ψ + y ψ S . We see that in this case m χ = m ψ = 4 TeV, as long as S = 0, cf., figure 3. However, when S = 0, there is a large negative contribution to m ψ , since y ψ < 0. The mass of S also becomes relatively small between T S and T E , so the channel χ ↔ ψS opens. At T E there is a first-order phase transition and, since S goes to zero, m ψ (T ) returns tõ m ψ and the channel χ ↔ ψS abruptly closes.
As mentioned above, the large Yukawa coupling |y ψ | 4 ensures that ψ remains in equilibrium throughout the depletion of χ. S has a large portal coupling so can easily stay in equilibrium via SS ↔ HH. As discussed in section 4, the rates of all 2-to-2 processes involving χ, except χS ↔ ψS, remain significantly below the Hubble rate at all times. The rate of χS ↔ ψS is suppressed when m χ ≈ m ψ (T ) but enhanced when m S (T ) ≈ 0. This means that the rate of χS ↔ ψS can be larger than the Hubble rate in the period where S obtains a vev. Although it remains subdominant to χ ↔ ψS in this period, we include it in out numerical calculations. Like χ ↔ ψS, the rate of χS ↔ ψS abruptly reduces below the Hubble rate at T E .
In figure 12 (right) we see the yield of χ as a function of x. As in section 7, χ initially freezes-out when relativistic, before beginning to deplete at T S . The depletion continues until it abruptly halts at T E . It stops more abruptly than in section 7 due to  Figure 13. Left: the evolution of the yield for different values of y χ . Right: the y χ required to produced the observed relic abundance in the λ S -m S (T = 0) plane. In the red region S does not obtain a vev, while in the blue region electroweak symmetry is not broken at T = 0. The grey pixels show where ψ freezes-out before T E while at the pink pixels the equilibrium abundance of χ at T E is too large for Y χ to reduce to the observed value.
the first-order phase transition, which quickly changes the mass of ψ at T E . After T E , χ stabilises and the yield remains constant. We see that there is a choice of y ψ such that the observed relic abundance is obtained. Again, once χ has departed from equilibrium at high temperatures, it does not return to equilibrium. We see that properly accounting for these finite temperature effects is essential for calculating the correct abundance.
We now consider the parameter space of this scenario. In figure 13 (left) we show the evolution of the yield for different Yukawa couplings y χψ . We see that the final yield is again exponentially sensitive to this coupling. A small change in the coupling results in an orders of magnitude change in the final yield. When the Yukawa coupling is y χψ 1.2 × 10 −7 we see that the rate is fast enough to bring χ back into equilibrium. We here explicitly see the relationship between instantaneous freeze-out and decaying dark matter. If χ does go back into equilibrium, the final abundance is simply set by the equilibrium abundance at T E , as in section 6. We note that the abundance of χ cannot go below the equilibrium abundance with this mechanism.
In figure 13 (right) we show the value of y χψ required to produce the equilibrium abundance of χ in the m S (T = 0)-λ p plane. We have increased the χ and ψ mass to 4.5 TeV, since otherwise the equilibrium abundance of χ is above the observed abundance at T E for many points, meaning that decay cannot set the abundance (this is still the case for some points, shown in pink). In the red region, µ 2 S is negative so S never obtains a vev and χ does not decay. In the blue region, the deepest minima at zero temperature have S = 0 and h = 0, so electroweak symmetry is not broken (i.e., this region is unphysical). In the region between these curves we show a pixellated region where the universe successfully transitions to the EWSB minimum. We see that larger couplings y χψ JHEP05(2019)070 are required near the red region, where symmetry breaking is weak and S only obtains a vev for a short amount of time. Closer to the blue region the required Yukawa coupling is smaller, as there is a longer amount of time in which χ can decay. The grey pixels show where ψ freezes-out before T E .
We do not show the variation with m χ and m ψ since this mechanism only works in a relatively small window. If these masses become much larger, then ψ is no longer a subdominant relic. On the other hand, since T E ≈ 100 GeV, if the fermions are much lighter, then the equilibrium abundance at T E is too large and χ overcloses the universe.

Experimental constraints
We now briefly consider experimental tests of the model in different regions of its parameter space. The usual probes of dark matter, direct and indirect detection of the galactic population and direct dark matter production at colliders, may be effective for the parameter space for instantaneous freeze-out, considered in section 6. However, since in this section the relic abundance is not determined by y χψ , there is no expectation of a coupling of a certain size (as there is in typical freeze-out scenarios). In sections 7 and 8, where we consider decaying dark matter, these searches are hindered by the small coupling y χψ , as is common for freeze-in models. However, direct detection is beginning to probe even these small couplings [36].
Aside from directly searching for χ, the model may first be probed via its other particle content. In all cases the energy density in χ must be passed to the SM thermal bath. This means that S and ψ must be in equilibrium while the χ ↔ ψS channel is open, which in turn means that the S 2 |H| 2 portal coupling must be ≈ 10 −3 (although for masses less than a GeV it could be as small as 10 −5 ). If S is lighter than m h /2, it may be detected by precision measurements of the Higgs. The current limit excludes λ p 10 −2 [37,38] while future colliders such as the ILC and CEPC will probe λ p 10 −3 [39][40][41]. Pairs of heavy S particles may be produced via an off-shell Higgs at proton colliders. A future 100 TeV machine with 30 ab −1 of integrated luminosity will be able to exclude the two-step phase transition region of section 8 at 2σ but not discover an S in this region at 5σ [42][43][44]. It will not be able to probe the scalar sector in sections 6 and 7.
In sections 6 and 7, any remaining ψ and S can decay, so these will not make any of the dark matter relic abundance. However, in section 8 there is a subdominant population of both ψ and S (S has no vev at T = 0 so the Z 2 symmetry may be intact). Assuming that the Z 2 symmetry is not softly broken by effects from particles beyond this simple model, the subdominant populations of S and ψ may be detected by direct or indirect detection experiments. However, the rates seen in these experiments depend on their abundances, which are unconnected to the χ abundance.
In the scenario in section 8 there is a first-order phase transition, as h obtains a vev and S loses its vev. If this phase transition is strong enough there is the interesting possibility of detecting a stochastic gravitational wave background resulting from the first order phase transition in upcoming space based interferometers such as LISA [45].
In this work we have explored the impact of finite temperature corrections on the relic abundance of a simple model of dark matter, which consists of two dark sector fermions (one being a dark matter candidate, χ, and a further dark sector fermion ψ) and a new scalar, S. We calculate the leading finite temperature effects, using the effective potential of the new scalar and the SM Higgs, and track the dark matter abundance using Boltzmann equations. We emphasise that in different regions of parameter space the new scalar can either simply obtain a vev (which we focus on in sections 6 and 7), or pass through a twostep phase transition during EWSB (section 8). These finite temperature effects provide corrections to the scalar and fermion masses. We assume that χ was in equilibrium at high temperatures and focus on regions of parameter space where decays and inverse decays χ ↔ ψS are kinematically allowed for a period of time. We also ensure that S and ψ remain in equilibrium with the SM thermal bath during the time of interest. We have explored the parameter space of the model and identified different regions where finite temperature effects lead to non-standard mechanisms of dark matter production.
In section 6 we discussed instantaneous freeze-out and considered the parameter space where χ can quickly come into equilibrium through the processes χ ↔ ψS and χS ↔ ψS. When these channels abruptly become inefficient, χ instantaneously freezes-out, resulting in a tight relationship between the temperature of decoupling and the resulting relic abundance. We find viable dark matter candidates in the range m χ ≈ (3 GeV-300 GeV), although χ and ψ should have similar mass (at the 5% level) to avoid significant fine tuning. However, the final abundance depends exponentially on the temperature at which the χ ↔ ψS and χS ↔ ψS channels become inefficient.
In section 7 we discuss the case where χ is only weakly coupled to the thermal bath and consider decaying dark matter via a one-step phase transition. In this case, rather than quickly coming into equilibrium, the yield of χ slowly approaches the equilibrium yield. The final abundance of χ depends on all of the model parameters in a complicated manner. We again find viable dark matter candidates in the range m χ ≈ (3 GeV-300 GeV). In this scenario, the final abundance is exponentially sensitive to theψχS coupling.
In section 8 we again focus on the case where χ is weakly coupled and describe decaying dark matter via a two-step phase transition (although we point out that instantaneous freeze-out can also occur during a two-step phase transition). We show that by changing the sign of a Yukawa coupling, we can again achieve the observed dark matter abundance with a period of decay. In this case S must have a mass similar to the SM Higgs, resulting in a smaller region of viable parameter space.
Finally we survey the experimental prospects for detecting particles in these scenarios. In all cases the energy density in χ must be passed to the SM thermal bath. This means that S must be in equilibrium while the dominant χ ↔ ψS channel is open. The S 2 |H| 2 portal coupling cannot be too small, so a light S may be detected by precision measurements of the Higgs, while a heavy S may be detected at a future 100 TeV proton collider. Direct and indirect detection of χ is usually hindered by its small couplings. In section 8 the subdominant population of ψ may be detected by these means. In this case there are also JHEP05(2019)070 interesting possibilities in detecting a stochastic gravitational wave background resulting from a first order phase transition.
Although we study a simple model of dark matter, these mechanisms for setting the dark matter relic abundance, some which we discuss here for the first time, may be relevant in a wide range of Beyond the Standard Model theories.