Self-heating of Strongly Interacting Massive Particles

It was recently pointed out that semi-annihilating dark matter (DM) may experience a novel temperature evolution dubbed as self-heating. Exothermic semi-annihilation converts the DM mass to the kinetic energy. This yields a unique DM temperature evolution, $T_{\chi} \propto 1 / a$, in contrast to $ T_{\chi} \propto 1 / a^{2}$ for free-streaming non-relativistic particles. Self-heating continues as long as self-scattering sufficiently redistributes the energy of DM particles. In this paper, we study the evolution of cosmological perturbations in self-heating DM. We find that sub-GeV self-heating DM leaves a cutoff on the subgalactic scale of the matter power spectrum when the self-scattering cross section is $\sigma_{\rm self} / m_{\chi} \sim {\cal O} (1) \,{\rm cm}^{2} /{\rm g}$. Then we present a particle physics realization of the self-heating DM scenario. The model is based on recently proposed strongly interacting massive particles with pion-like particles in a QCD-like sector. Pion-like particles semi-annihilate into an axion-like particle, which is thermalized with dark radiation. The dark radiation temperature is smaller than the standard model temperature, evading the constraint from the effective number of neutrino degrees of freedom. It is easily realized when the dark sector is populated from the standard model sector through a small coupling.

cutoff on the subgalactic scale of the matter power spectrum when the self-scattering cross section is σ self /m χ ∼ O(1) cm 2 /g.Then we present a particle physics realization of the self-heating DM scenario.The model is based on recently proposed strongly interacting massive particles with pionlike particles in a QCD-like sector.Pion-like particles semi-annihilate into an axion-like particle, which is thermalized with dark radiation.The dark radiation temperature is smaller than the standard model temperature, evading the constraint from the effective number of neutrino degrees of freedom.It is easily realized when the dark sector is populated from the standard model sector through a small coupling.

I. INTRODUCTION
Accumulated observational data begin to test our understanding of how dark matter (DM) is distributed over the Universe and how the DM distribution evolves in time [1].In regard to the structure formation of the Universe, DM is often described as a perfect fluid with zero temperature, which is referred to as cold dark matter (CDM).The CDM paradigm is successful in reproducing the large-scale structure of the Universe observed through cosmic microwave background (CMB) anisotropies [2] and galaxy clustering [3].
Contrary to its success on large scales, CDM predictions appear to be incompatible with the observations on smaller scales [4].Although baryonic physics may play an important role [5][6][7], these small scale issues may be hinting to alternative DM models.One example of the CDM failure is the missing satellite problem, indicating that CDM overpredicts the number of dwarf-size subhalos in a Milky Way-size halo when compared to that of the observed satellite galaxies [8,9].Warm dark matter (WDM) is an interesting possibility in this respect.Gravitational clustering of WDM particles is interrupted on a subgalactic scale because of a sizable thermal velocity of v/c ∼ 10 −3 -10 −4 at the matter-radiation equality.This suppresses dwarf-size halo formation [10,11].A sterile neutrino with a keV mass is a good benchmark model of WDM, where its phenomenology is described by the mass and the mixing angle with an active neutrino in the simplest setup [12].
Another example is the core-cusp problem: Some dwarf-size halos have a cuspy profile as predicted by CDM, while others have a cored profile [13][14][15].Self-interacting dark matter (SIDM) is an interesting solution.The self-scattering cross section of σ self /m χ ∼ 1 cm 2 /g with the DM mass m χ leads to iso-thermalization of DM particles, whose distribution is characterized by a kpc core [16,17].SIDM reproduces the observed diversity of the rotation velocity among similar-size halos by adjusting its distribution sensitively to the baryon distribution [18,19].
A new possibility, which is called self-heating DM, has been proposed recently [20].A characteristic feature of this scenario is that the strength of self-scattering is related to the thermal velocity of DM and, hence, potentially solves several small scale issues simultaneously.The original proposal was based on semi-annihilating DM [21][22][23], χχ → χφ with a light particle φ in the thermal bath.A key observation is that semi-annihilation converts the mass of DM into the kinetic energy, which leads to the novel DM temperature evolution, T χ ∝ 1/a after the freeze-out, instead of T χ ∝ 1/a 2 for free-streaming non-relativistic particles.Self-heating continues as long as self-scattering occurs rapidly.DM self-scattering is an essential ingredient of the scenario because it redistributes the large kinetic energy of the boosted DM particles through the semi-annihilation over the whole DM particles.Stronger self-scattering elongates the duration of self-heating, resulting in a larger thermal motion of DM particles.
In this paper, we investigate the impact of DM self-heating on the matter distribution of the Universe, and propose a viable particle physics realization of self-heating DM.To study the structure formation, we derive the evolution equation of the cosmological perturbations in the self-heating DM scenario.We show that semi-annihilation not only changes the DM temperature evolution but also affects the entropic perturbation of DM.We follow the evolution of the primordial perturbations numerically and show that self-heating DM can leave a subgalactic-scale cutoff in the linear matter power spectrum.We extend the recently proposed strongly interacting massive particle (SIMP) model [24,25] to realize self-heating DM.Pion-like particles in a QCD-like sector semi-annihilate into an axion-like particle (ALP), which is thermalized with dark radiation.We identify a parameter region that is compatible with observational constraints.This paper is organized as follows.In Sec.II, we describe the thermal history of the self-heating DM from its freeze-out to structure formation of the Universe.In Sec.III, we discuss an extension of the SIMP model with a light ALP and dark radiation.We conclude in Sec.IV.In Appendix A, we provide a detailed derivation of the co-evolution equations of the DM number density and temperature.In Appendix B, we present the evolution equations of cosmological perturbations in the self-heating DM scenario.

II. THERMAL HISTORY OF SELF-HEATING DM FROM THE FREEZE-OUT TO STRUCTURE FORMATION
In this section, we discuss the freeze-out of the DM number density, its novel temperature evolution afterwards, and the evolution of cosmological perturbations in the selfheating DM scenario.For our purpose, we consider a scalar DM χ and a light mediator φ.We mainly consider two types of interaction: One is self-scattering χχ ↔ χχ, and the other one is semi-annihilation χχ ↔ χφ, which are the minimal ingredients to realize selfheating of DM.The presence of efficient self-scattering enforces the DM distribution to be where n eq χ = (m 2 χ T χ /2π 2 )K 2 (m χ /T χ ) and K 2 is the 2nd-order modified Bessel function of the second kind.We also assume that a light mediator remains in thermal equilibrium by contacting either to the SM sector or to the dark sector and, hence, f φ = exp[−E φ /T φ ], where T φ = T SM or T φ = T DR , respectively.Its number density is thus given by n eq φ = (m 2 φ T φ /2π 2 )K 2 (m φ /T φ ).In general, annihilation χχ ↔ φφ and also elastic scattering χφ ↔ χφ may exist.

A. Homogeneous and isotropic evolution
The co-evolution equations of n χ and T χ are given, respectively, as where σ semi (σ inv ) is the cross section for the χχ → χφ (χφ → χχ) process.See Appendix A for the derivation.We have added the elastic scattering term for completeness.The expression of momentum exchange rate γ χφ→χφ can be found in Ref. [26] and thus is not repeated here.We remark that self-scattering does not contribute to these equations because it conserves the number and energy of DM particles.The relic abundance of χ coincides with the observed one when where (σv rel ) can (3 × 10 −26 cm 3 /s) is a canonical cross section reproducing the observed relic density of thermal DM [27].The canonical cross section is rescaled by the factor of (T φ /T SM ) fo because T SM, fo (m χ /20)(T SM /T φ ) fo .
We present numerical results in Fig. 1, which show the evolution of DM yield Y χ = n χ /s and the temperature ratio T χ /T φ .Here s = (2π 2 /45) g * s, SM T 3 is the entropy density, and g * s, SM is the effective number of relativistic degrees of freedom in the SM sector.In this figure, we assume that γ χφ→χφ /H 1 and T φ = T SM .Like usual thermal DM with T χ = T φ , the relic abundance is determined around T φ, fo m χ /20 (see Appendix A for a small difference).On the other hand, the DM temperature evolution shows a unique behavior.
Especially, after the freeze-out, the DM temperature scales as T χ ∝ 1/a (see Appendix A for a thermodynamic derivation) despite the fact that DM particles are non-relativistic and that no elastic scattering equilibrates T χ and T φ .Self-heating of DM occurs because a small portion of DM still undergoes semi-annihilation after the freeze-out, and gain the kinetic energy of the order of its mass, which is much larger than the DM temperature.We find that the ratio between the two temperatures approaches1

T /T
< l a t e x i t s h a 1 _ b a s e 6 4 = " U L w 9 a 5 E A < l a t e x i t s h a 1 _ b a s e 6 4 = " U L w 9 a 5 E A  The plot shows the evolution of DM yield (black) and also the evolution of the DM temperature (red).The dashed line corresponds to the equilibrium value.For this plot, we assume no elastic scatterings; thus, the chemical and kinetic decoupling take place at the same time as semi-annihlation decouples from the plasma.One can see that the DM temperature scales as T φ scales even though DM is non-relativistic and kinetically decoupled from the thermal bath.
where γ = (5/4)[1 − m 2 φ /(5m 2 χ )] is the Lorentz boost factor of the final state DM particle.Although we ignored elastic scatterings above, we emphasize that the self-heating is a generic feature of exothermic semi-annihilation even in the presence of elastic scattering.If (γ χφ→χφ /H)| Tχ=T χ, fo 1, the elastic scattering is able to maintain the kinetic equilibrium of DM after the freeze-out, resulting in T χ = T φ .In this case, Eq. ( 1) reproduce a usual discussion found in Refs.[21][22][23], since J (T χ = T φ , T φ ) = 1.Eventually, elastic scattering decouples.After this kinetic decoupling, DM begins to self-heat, and the DM temperature continues to scale as T χ ∝ 1/a.See Fig. 2 for schematic picture of the DM temperature evolution.
We stress that DM self-scattering is necessary for DM self-heating.If self-scattering is absent, semi-annihilation merely produces a small portion of boosted DM particles, which act as a hot component of DM.Thus, the unique temperature evolution continues only The larger the self-scattering cross section is, the longer the self-heating lasts and thus results in a larger thermal motion of DM.The self-scattering decouples at where n = 1/3 for T SM, self > T SM, eq , while n = 1/4 when T SM, self < T SM, eq .Here, T SM, eq 0.8 eV is the SM temperature at the matter-radiation equality.When deriving Eq. ( 5), we the era of self-heating freeze-out & kinetic decoupling of ( semi ) freeze-out of self-interaction ( ) time the era of self-heating freeze-out of self-interaction kinetic decoupling of freeze-out of ( ) ( ) ( semi or ) time 1.The kinetic equilibrium could be maintained after the chemical freeze-out.In this case, the self-heating begins after the decoupling of elastic scattering, and continues until Γ self = H.
and solve Γ self = H for T SM .Remember that σ self /m χ ∼ 1 cm 2 /g forms a kpc core in a subgalactic halo [28]. 2 A resultant large thermal motion of sub-GeV SIDM leaves a subgalactic-scale cutoff in the linear matter power spectrum like keV WDM.We will discuss the structure formation of self-heating DM in the next section.
Before closing this section, we emphasize that the thermal history of φ plays a significant role both in DM searches and in the evolution of the Universe.If φ is massless, it changes the expansion rate of the Universe, as it contributes to the total energy density.Its impact is described by the change in the effective number of neutrino degrees of freedom ∆N eff and is constrained by big bang nucleosynthesis (BBN) and CMB.The latest constraint from the Planck Collaboration is ∆N eff = 3.15 ± 0.24 [2].For massive φ, it may overclose the Universe unless it decays or they annihilate into light particles.If φ decays into visible particles, the late time production of φ through semi-annihilation is severely constrained by the Galactic and extra-Galactic gamma-ray searches as well as the CMB measurement.These constraints require m χ 10 GeV [30,31], which results in a shorter duration of selfheating [see Eq. ( 5)]. 3 The other possibility is that φ decays into dark radiation.In this case, 2 It is claimed that semi-annihilation cooperates with self-interaction to flatten the inner density profile of a subgalactic halo [29].The impact of the self-heating is more significant in smaller halos, and it may alleviate a required strength of self-scattering for solving the core-cusp problem.We do not take this effect into account for simplicity in this paper. 3Thermal relic DM with a sub-GeV mass is still viable if DM particles annihilate into heavier particles [32].
constraints from DM indirect searches cannot be applied, although ∆N eff still constrains the scenario.It will be discussed in Section III.

B. Perturbed evolution
The novel evolution of the DM temperature affects the evolution of cosmological perturbations.We derive evolution equations of cosmological perturbations in the self-heating scenario in Appendix B. In the case of self-scattering DM, relevant variables are density contrast δ χ , velocity divergence θ χ , and entropy perturbation π χ .The linearized equations for the cosmological perturbation for self-heating DM are given as where the prime is a derivative with respect to the conformal time and H = a /a.Here we take the conformal Newtonian gauge: We cross-check the result in the synchronous gauge in Appendix B. Semi-annihilation affects δ χ mainly through the relatively large c 2 sχ (4/3)T χ /m χ .One can obtain the matter power spectrum by solving the above equations with equations of the equation of state ω χ and sound speed squared c 2 sχ : with Γ semi = n χ σ semi v rel ∝ 1/a 3 .As the interaction rate for self-scattering becomes smaller than the Hubble expansion rate, the DM temperature behaves as that of free-streaming non-relativistic particles, i.e., T χ ∝ 1/a 2 .
By modifying the publicly available Boltzmann solver CLASS [34], we obtain the present linear matter power spectrum as shown in Fig. 3.The matter power spectrum exhibits In the case of semi-annihilation, this can be realized if φ is heavier than the DM mass [33].In this case, we expect no self-heating of DM, because semi-annihilation is no longer exothermic.
FIG. 3: (Left) Linear matter power spectrum for self-heating DM (red) and CDM (black) at the present Universe.(Right) Linear matter power spectrum normalized with respect to that of CDM.We also show the power spectrum in the thermal WDM models (blue) with m WDM = 4.09 keV and m WDM = 5.3 keV, where the latter is the latest lower bound on the mass of WDM [36].We choose σ self /m χ = 1 cm 2 /g such that self-scattering decouples at T SM, self = 0.5 eV, 0.8 eV, and 1.6 eV for m χ = 0.1 GeV, 1 GeV, and 10 GeV, respectively.We also assume T φ = T SM .
a sharp cutoff around k = O(100) Mpc −1 .This scale turns out to well match the Jeans instability scale at the matter-radiation equality [20,35]: where the temperature ratio r χφ = T χ /T φ should be evaluated at a = min(a self , a eq ).
We compare the resultant matter power spectra in self-heating DM to those in the thermal WDM model (see, e.g., Ref. [37] for details).The Lyman-α forest data constrain the WDM mass as m WDM > 5.3 keV [36].In the case of WDM, the Jeans instability scale appears to be k J,WDM 180 Mpc −1 (m WDM /5.3 keV) 4/3 , while a drop of power in the matter power spectrum takes place around k J,WDM /4, where an additional order one factor can be attributed to the free-streaming of WDM during the radiation-dominated Universe [38].By equating k J k J,WDM /4, we find the following correspondence between WDM and self-heating DM: For the massless mediator φ sharing the temperature with SM particles, the current limit on the thermal WDM mass, m WDM ≥ 5.3 keV, translates into the mass bound of self-heating DM as m χ ≥ 0.1 GeV if the self-interaction decouples after the matter-radiation equality.

III. EXTENDED SIMP MODEL WITH AN ALP AND DARK RADIATION
In this section, we propose a realization of self-heating DM.We extend the SIMP model [39] with an ALP [33] by introducing dark radiation.The Lagrangian density is given by where H a µν ( H a µν ) is the (dual) field strength of the confining SU(N H ) gauge field and X a µν ( X a µν ) is the (dual) field strength of the quasi-perturbative SU(N X ) gauge field.Gauge charges of the matter contents are summarized in Table I.N and N are N f -flavored (N f ≥ 3).
As the PQ symmetry is spontaneously broken, an ALP arises.The mass of Q and L originates from vacuum expectation value Φ and is assumed to be heavier than the confinement scale of the SU(N H ) gauge group.After integrating out Q and L, we find that the ALP couples to both H a µν H aµν and X a µν X aµν .Once SU(N H ) confines, SU(N f ) L ×SU(N f ) R global symmetry breaks down to the diagonal subgroup SU(N f ) V , and pion-like particles χ emerge.For low energy phenomenology, we obtain the following effective Lagrangian: The expression of L chiral can be found in Ref. [33] and thus is not repeated here.We emphasize that SU(N f ) V is the exact symmetry of the theory; hence, no decay operator of DM pions is allowed.The ALP φ and dark radiation X form a dark plasma with the temperature T φ = T DR .In the following, we take N f = 4, N H = 3, and N X = 2 as a benchmark for the analysis.
Low energy phenomenology is described by three parameters: DM pion mass m χ , DM pion decay constant f χ , and ALP decay constant f .The self-scattering cross section for solving the core-cusp problem, and the semi-annihilation cross section for the observed DM relic density [see Eq. ( 3)] are achieved as [33] Note that we only consider CP-conserving interactions of DM pions and ALPs, because the CP-violaing vacuum angle in the hidden confining sector dynamically vanishes due to an ALP.Although we consider semi-annihilation as a dominant process for the freeze-out of the DM number density, there exist two additional number-changing processes: the 3-to-2 process of χχχ → χχ and annihilation χχ → φφ.Since the annihilation cross section is suppressed by (f χ /f ) 2 compared to that of semi-annihilation, it is sub-domnant.Meanwhile, the 3-to-2 process is sub-dominant for where we have chosen f χ and f such that σ self /m χ 1 cm 2 /g and σ semi v rel = (σv rel ) can (T DR /T SM ) fo .Note that the ALP obtains its mass from explicit breaking of chiral symmetry of the QCD-like sector, and its mass is given as [33] The elastic scattering between the ALP and DM also exists in this model.For the selfheating to occur, the elastic scattering should decouple before the self-scattering decouples.
This can be trivially achieved in our model, because the strength of elastic scattering is suppressed by (f χ /f ) 2 compared to that of semi-annihilation, and by (f χ /f ) 4 compared to that of self-interaction.As a consequence, the momentum exchange rate due to elastic scattering is Thus, elastic scattering is inefficient during and after the freeze-out of DM, and the self-heating begins roughly after the freeze-out of DM.
We assumed that the ALP and dark gauge boson form a thermal bath with T DR during the freeze-out of DM.The XX → XXX rate [40,41] is sufficiently large for α X > O(10 −10 ).
The φX → XX keeps φ in the thermal bath of X until T SM = 1 MeV as [42,43] An ALP decays when it is still semi-relativistic with the thermally averaged decay rate of Non-Abelian dark radiation does not confine until the present Universe for α X (T DR = T DR, fo ) 0.03 (2/N X ), since the confinement scale is given by Λ where α X,0 = α X (T DR = µ 0 ).Neither the ALP nor dark radiation overcloses the Universe.
The model Lagrangian (14) does not contain any interaction that equilibrates the SM and dark sector.Indeed, if two sectors are in thermal equilibrium with each other in the early Universe, ∆N eff tends to exceed unity, which is strongly disfavored by BBN and CMB.The inflaton can decay both into the SM sector and into the dark sector.In this case, (T DM /T SM ) is determined by the branching ratio.Even if the inflaton predominantly decays into the SM sector, the dark sector can be populated through a feeble interaction to the SM sector, while not being completely thermalized with SM particles.In addition to the interactions in Eq. ( 14), we may consider a Higgs portal coupling: The continuous production of the dark sector increases the temperature ratio between the dark sector and the SM plasma until the electroweak phase transition.The temperature ratio between the two sectors at the electroweak phase transition is estimated as , (24) where g * , DR, ew = 83.5 takes into account all the degrees of freedom of particles in Eq. ( 14).
After the freeze-out of DM but before the decay of the ALP, we find g * s, DR = 7.Thus, (T DR /T SM )| ew 0.5 is consistent with ∆N eff < 3.4 [2].

IV. CONCLUSION
Self-heating of semi-annihilating DM can suppress subgalactic-scale structure formation when it lasts until the matter-radiation equality.The reduced number of dwarf-size halos can reconcile the possible tension between the CDM paradigm and the observation.
The self-heating of sub-GeV DM is maintained until the matter-radiation equality for σ self /m χ O(0.1-1) cm 2 /g.We have followed the evolution of cosmological perturbations and demonstrated that self-heating sub-GeV DM indeed leaves a cutoff on the subgalactic scale of the linear matter power spectrum.
It is interesting that self-heating DM interrelates a subgalactic cutoff in the linear matter power spectrum and a kpc core of the DM distribution in halos through the thermalization of DM particles.We can take full advantage of astrophysical and cosmological searches of WDM and SIDM to probe self-heating DM.For example, we could tighten the range of the self-interaction strength by analyzing line-of-sight velocity dispersions of dwarf spheroidal galaxies [44] and rotation curves of low-surface brightness galaxies [45] with a larger number of samples.The satellite number counts restrict the cutoff in the linear matter power spectrum [46,47].The matter distribution smoother than the CDM prediction will be tested by the perturbations on strongly lensed systems [48][49][50][51].The top-down structure formation in contrast to the bottom-up one in the CDM paradigm is tested by the further discoveries of high-z galaxies [52][53][54] and by multiple probes of the reionization epoch such as the 21 cm brightness temperature [55][56][57] and its fluctuations due to minihalos [58].
We have proposed an extension of the SIMP model with an ALP and dark radiation as a particle physics realization of sub-GeV self-heating DM.DM pions semi-annihilate into an ALP, which decays into dark radiation.Dark radiation and an ALP forms thermal equilibrium.We have shown that self-heating can be realized for a certain range of model parameters.When the dark sector is populated from the SM sector through a Higgs portal, we can produce the dark sector particles while being consistent with the constraints from BBN and CMB measurements.
We have focused on a velocity-independent self-scattering cross section.On the other hand, the self-scattering cross section diminishing with an increasing velocity may be favored by the constraints from galaxy cluster ellipticities and bullet clusters.One way to realize the velocity-dependent self-scattering cross section is to introduce a light mediator coupling to two DM particles.Extending our discussion to such a case will be intriguing.
T χ , the equation for the number density is not closed by itself.To correctly describe the evolution of the system, one must track the DM temperature evolution as well.
To obtain the evolution equation of the DM temperature, we integrate the Boltzmann equation without any weight.We find the evolution equation of the energy density: with pressure P χ .Since ρ χ = E χ n χ and P χ = p 2 χ / (3E χ ) n χ = T χ n χ , the above equation reads as where ∆E = E φ − E χ Tχ and the function K(T χ , T φ ) is defined as Figure . 1 presents our numerical result of the co-evolution of n χ and T χ .One can see that the freeze-out of DM yield Y χ proceeds in a similar way to a usual discussion found in Refs. [21][22][23].The DM yield is estimated as while determining the freeze-out temperature is ambiguous in our case.In the case of T χ = T φ , the freeze-out temperature is usually determined by ∆(x fo ) = Y χ −Y eq χ = cY eq χ (x fo ) with c being a numerical constant of order unity, where ∆ = (d ln Y eq χ /dx)/(−s σ semi v rel /H).In the self-heating scenario, however, the DM freeze-out is delayed, because the DM temperature shortly increases relative to T φ and enhances the backward semi-annihilation process.This can be seen from Fig. 4, where we present numerically computed ∆ = Y χ − J Y eq χ as a function of m χ /T φ .This delay results only in an O(10)% change of the final DM relic abundance.
the Universe, while the second term is due to the energy injection through semi-annihilation.
If T χ ∝ 1/a 2 as for free-streaming non-relativistic DM particles, then the second term increases with the expansion of the Universe as ∝ (m χ /T χ )(Γ semi /H) ∝ a and thus heats the DM particles.Hence, the DM temperature is determined by the balance between the adiabatic cooling and semi-annihilation heating as T χ ∝ 1/a.

Appendix B: Evolution equations of self-heating DM cosmological perturbations
We derive the evolution equations of cosmological perturbations.Taking the conformal Newtonian gauge [see Eq. ( 9)], we get the linearized Boltzmann equation given by We expand δf as a function of µ = k • q in terms of the Legendre polynomial: Multiplying the Legendre polynomial by the linearized Boltzmann equation and integrating it with respect to µ, we find the Boltzmann hierarchy as In terms of the DM fluid variables [59], one obtains  We have decomposed the pressure perturbation into the isentropic and entropic parts as δP χ = ρχ (c 2 sχ δ χ + π χ ), where the adiabatic sound speed squared is given by The dimensionless constants are defined as scattering of DM.In this limit, we find Eqs.( 6) - (8).We have substituted α 2 = 5, which can be derived as follows.In general, F (τ, k, q) can be expanded by a complete system of functions of q.We can expand F (τ, k, q) = n 1/(2πa 2 m χ T χ ) 2/3 y /2 L +1/2 n (y) e −y F n with y = q 2 /(2πa 2 m χ T χ ) and L α n being the associated Legendre function, as in Refs.[26,60].An advantage of this complete system is that n = 0 gives the dominant contribution in the non-relativistic limit.Since L α 0 = 1, one can find α 2 = 5.One may also calculate other α's in a similar manner, while they do not appear in Eqs. ( 6)-( 8) and thus are not given here.
By taking the synchronous gauge, one can derive the evolution equations given as θ χ = − Hθ χ + k 2 (c 2 sχ δ χ + π χ ) , (B19) where h denotes the trace of h ij .They are equivalent to those in the conformal Newtonian gauge under the gauge transformation [see Eq. ( 27) of Ref. [59]].To check the consistency, one needs to note that c 2 s Hk 2 α = c 2 s H(h + 6η )/2 is negligible when compared to θ χ for non-relativistic DM.
After decoupling of self-scattering [see Eq.( 5)], energetic DM particles through semiannihilation freely stream, while the majority of DM follows the Maxwell-Boltzmann distribution with the temperature T χ ∝ 1/a 2 .After that, Eqs. ( 6)-( 8) still describes the evolution of the majority of DM, although Eqs.(B15)-(B16) are no longer valid.For the decoupling of self-scattering, we take phenomenological approach, i.e., multiplying (1 − e −Γ self /H ) to the equations of ω χ and c 2 sχ [see Eqs. ( 10)- (11)].We may introduce an alternative cutoff, for where a self is the scale factor when Γ self = H.In this case, the decoupling takes place more rapidly compared to the case presented in the main text, and the cutoff in the matter power spectrum appears for larger k.See Fig. 5 for the numerical difference in the two descriptions.
FIG.1: Thermal history of self-heating DM.The plot shows the evolution of DM yield (black) and also the evolution of the DM temperature (red).The dashed line corresponds to the equilibrium value.For this plot, we assume no elastic scatterings; thus, the chemical and kinetic decoupling take place at the same time as semi-annihlation decouples from the plasma.One can see that the DM temperature scales as T φ scales even though DM is non-relativistic and kinetically decoupled from the thermal bath.

FIG. 2 :
FIG.2:The self-heating takes place between the kinetic decoupling and the freeze-out of DM self-interaction.(a) (γ χφ→χφ /H)| Tχ=T χ,fo1.The chemical and kinetic decoupling take place simultaneously.The self-heating of DM begins from the freeze-out of semi-annihilation to the freeze-out of self-interaction.(b) (γ χφ→χφ /H)| Tχ=T χ,fo 1.The kinetic equilibrium could be maintained after the chemical freeze-out.In this case, the self-heating begins after the decoupling of elastic scattering, and continues until Γ self = H.

m = 1
GeV m = 0.1 GeV m = 10 GeV m WDM = 5.3 keV < l a t e x i t s h a 1 _ b a s e 6 4 = " g Z z p l 2 d J 8 8 h y b M u 5 O 6 n W L o s 4 y m A H 7 I E D 4 I A z U A M 3 o A 4 a A I N H 8 A x e w Z v x Z L w Y 7 8 b H p L V k F D P b 4 A + M z x + S Y p c y < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " g Z z p l 2 d J 8 8 h y b M u 5 O 6 n W L o s 4 y m A H 7 I E D 4 I A z U A M 3 o A 4 a A I N H 8 A x e w Z v x Z L w Y 7 8 b H p L V k F D P b 4 A + M z x + S Y p c y < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " g Z z p l 2 d J 8 8 h y b M u 5 O 6 n W L o s 4 y m A H 7 I E D 4 I A z U A M 3 o A 4 a A I N H 8 A x e w Z v x Z L w Y 7 8 b H p L V k F D P b 4 A + M z x + S Y p c y < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " g Z z p l 2 d J 8 8 h y b M u 5 O 6 n W L o s 4 y m A H 7 I E D 4 I A z U A M 3 o A 4 a A I N H 8 A x e w Z v x Z L w Y 7 8 b H p L V k F D P b 4 A + M z x + S Y p c y < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " g Z z p l 2 d J 8 8 h y b M u 5 O 6 n W L o s 4 y m A H 7 I E D 4 I A z U A M 3 o A 4 a A I N H 8 A x e w Z v x Z L w Y 7 8 b H p L V k F D P b 4 A + M z x + S Y p c y < / l a t e x i t > m WDM = 4.09 keV < l a t e x i t s h a 1 _ b a s e 6 4 = " V + N d K S y / q y

FIG. 5 :
FIG.5:(Left) Linear matter matter power spectrum when the decoupling of self-scattering is described by (1 − e −Γ self /H ) (red) and when it is described by e −a/a self (green).(Right) Power spectrum in self-heating DM scenario relative to CDM.

TABLE I :
Gauge charges of matter contents.
The contribution of the dark sector to ∆N eff at the neutrino decoupling is given by