CMB spectral distortions

# The evolution of CMB spectral distortions in the early Universe

J. Chluba and R. A. Sunyaev
Canadian Institute for Theoretical Astrophysics, 60 St. George Street, Toronto, ON M5S 3H8, Canada
Max-Planck Institut für Astrophysik, Karl-Schwarzschild-Str. 1, D-85740 Garching, Germany
Space Research Institute, Russian Academy of Sciences, Profsoyuznaya 84/32, 117997 Moscow, Russia
E-mail: jchluba@cita.utoronto.ca
Accepted 2011 September 8. Received 2011 September 3
###### Abstract

The energy spectrum of the cosmic microwave background (CMB) allows constraining episodes of energy release in the early Universe. In this paper we revisit and refine the computations of the cosmological thermalization problem. For this purpose a new code, called CosmoTherm, was developed that allows solving the coupled photon-electron Boltzmann equation in the expanding, isotropic Universe for small spectral distortion in the CMB. We explicitly compute the shape of the spectral distortions caused by energy release due to (i) annihilating dark matter; (ii) decaying relict particles; (iii) dissipation of acoustic waves; and (iv) quasi-instantaneous heating. We also demonstrate that (v) the continuous interaction of CMB photons with adiabatically cooling non-relativistic electrons and baryons causes a negative -type CMB spectral distortion of in the GHz spectral band. We solve the thermalization problem including improved approximations for the double Compton and Bremsstrahlung emissivities, as well as the latest treatment of the cosmological recombination process. At redshifts the matter starts to cool significantly below the temperature of the CMB so that at very low frequencies free-free absorption alters the shape of primordial distortions significantly. In addition, the cooling electrons down-scatter CMB photons introducing a small late negative -type distortion at high frequencies. We also discuss our results in the light of the recently proposed CMB experiment Pixie, for which CosmoTherm should allow detailed forecasting. Our current computations show that for energy injection because of (ii) and (iv) Pixie should allow to improve existing limits, while the CMB distortions caused by the other processes seem to remain unobservable with the currently proposed sensitivities and spectral bands of Pixie.

###### keywords:
Cosmology: cosmic microwave background – theory – observations

## 1 Introduction

Investigations of the cosmic microwave background (CMB) temperature and polarization anisotropies without doubt have allowed modern cosmology to mature from an order of magnitude branch of physics to a precise scientific discipline, were theoretical predictions today are challenged by a vast amount of observational evidence. Since its discovery in the 60’s (Penzias & Wilson, 1965), the full CMB sky has been mapped by Cobe/Dmr (Smoot et al., 1992) and Wmap (Bennett et al., 2003; Page et al., 2006), while at small angular scales many balloon-borne and ground-based CMB experiments like Maxima (Hanany et al., 2000), Boomerang (Netterfield et al., 2002), Dasi (Halverson et al., 2002), Archeops (Benoît et al., 2003), Cbi (Pearson et al., 2003) and Vsa (Rubiño-Martin et al., 2003) provided important additional leverage, significantly helping to tighten the joint constraints on cosmological parameters. Presently Planck is producing CMB data with unprecedented precision, while both ACT (e.g., see Hajian et al., 2010; Dunkley et al., 2010; Das et al., 2011) and SPT (Lueker et al., 2010; Vanderlinde et al., 2010) are pushing the frontier of CMB power spectra at small angular scales. In the near future SPTpol (McMahon et al., 2009) and (Niemack et al., 2010) will provide additional small scale -mode polarization data, complementing the polarization power spectra obtained with Planck and further increasing the significance of the power spectra. Also Quiet (QUIET Collaboration et al., 2010) already now is observing the CMB in and polarization with the aim to tighten the constraints on the potential level of -modes at multipoles , while and are getting ready to deliver precision CMB polarization data in the near future.

With the new datasets, cosmologists will be able to determine the key cosmological parameters with extremely high precision, also making it possible to distinguish between various models of inflation by measuring the precise value of the spectral index of scalar perturbations, , and constraining its possible running, . Furthermore, with the Planck satellite primordial -mode polarization may be discovered, if the tensor to scalar ratio, , is larger than a few percent (The Planck Collaboration, 2006). This would be the smoking gun of inflation and indicate the existence of gravitational waves, one of the main scientific goals of Planck, Quiet, Spider, Polarbear and other CMB experiments.

However, the anisotropies are not the only piece of information about the early Universe the CMB offers us. In addition, the energy spectrum of the CMB tells the tale about the thermal history of the Universe at very early times, well before any structures had formed and when the baryonic matter and radiation were tightly coupled via Compton scattering. It is well known that any perturbation of the full thermodynamic equilibrium between photons and baryons in the early Universe inevitably leads to spectral distortions in the CMB (Zeldovich & Sunyaev, 1969; Sunyaev & Zeldovich, 1970c, d; Zeldovich et al., 1972; Sunyaev, 1974; Illarionov & Sunyaev, 1974). This for example could be caused by some release of energy by decaying relic particles (Hu & Silk, 1993b) or evaporating primordial black holes (Carr et al., 2010), and depending on their lifetime777For black holes the evaporation timescale is (Hawking, 1974). and the total amount of energy that was transferred to the medium, the residual CMB spectral distortion will be smaller or larger.

The process that erases possible spectral distortions of the CMB in the early Universe and thereby attempts to restores full thermodynamic equilibrium is called thermalization. In the past it has been studied extensively, both analytically and numerically (Illarionov & Sunyaev, 1975a, b; Danese & de Zotti, 1977, 1982; Burigana et al., 1991b, a; Daly, 1991; Barrow & Coles, 1991; Burigana, 1993; Hu & Silk, 1993a, b; Hu et al., 1994; Burigana et al., 1995; Burigana & Salvaterra, 2003; Lamon & Durrer, 2006; Procopio & Burigana, 2009). For very early energy release (), photon production at low frequencies by double Compton emission and the up-scattering of photons to high frequencies by Compton scattering are so efficient that the thermalization process is practically perfect for nearly arbitrary amounts of energy release, such that no distortion should remain today (e.g., see Burigana et al., 1991b). However, below this redshift the CMB spectrum becomes vulnerable, and spectral distortions that could still be observable today can form.

In connection with early energy release two types of CMB distortions, called - and -type distortions, are usually distinguished (e.g., see Illarionov & Sunyaev, 1975a, b). The first type of distortion is characterized by a frequency-dependent chemical potential, , that is very close to constant at high frequencies, and vanishes at very low frequencies. This type of distortion forms at times when the Compton process is still able to achieve full kinetic equilibrium with the electrons (redshifts ), and the CMB spectrum reaches a Bose-Einstein distribution with occupation number , where , with denoting the CMB temperature. However, photon production by double Compton and Bremsstrahlung emission stopped being very fast, so that only at low frequencies full equilibrium can be restored.

At low redshifts (), on the other hand, -type distortions form. These are also known in connection with SZ clusters (Zeldovich & Sunyaev, 1969) and are characterized by a constant temperature decrement at low and an increment at high frequencies. This type of distortion is formed when the Compton process ceases to be extremely efficient, so that full kinetic equilibrium between electrons and photons can no longer be achieved. Photons produced at low frequencies by Bremsstrahlung are therefore not up-scattered very much and remain close to the frequency band they were emitted.

Both types of distortions were tightly constrained by Cobe/Firas (Mather et al., 1994; Fixsen et al., 1996) only leaving room for and at cm and dm wavelength. This imposes strong limits on possible energy release by decaying particles (Hu & Silk, 1993b), which for some particle lifetimes even today are comparable to those derived, for example, from big bang nucleosynthesis (Kusakabe et al., 2006; Kohri & Takahashi, 2010). However, since the 90’s detector technology has advanced very much. Already about ten years ago, a Firas type-II experiment with roughly 100 times better sensitivity in principle would have been possible (Fixsen & Mather, 2002). Additional new technology was used in several flights of the Arcade balloon-borne instrument (Kogut et al., 2004; Fixsen et al., 2011; Kogut et al., 2011), from which improved limits on the CMB temperature at GHz were derived (Seiffert et al., 2011). However, until now no signatures of primordial distortions were found, and the limits on overall CMB spectral distortions derived from Cobe/Firas until today remain most stringent.

However, very recently a new space-based mission, called Pixie (Kogut et al., 2011), was proposed. Pixie will have exquisite spectral capabilities with 400 channels at frequencies ranging from GHz to THz. The claimed sensitivities should allow to tighten the constraints on - and -type distortions by about three orders of magnitude. A detection of and could therefore become possible at the level (see Fig. 12 of Kogut et al., 2011), at least in terms of spectral sensitivity. However, it is not yet totally clear what the foreground limitations to measuring primordial spectral distortions will be.

Under simplifying assumptions analytic approximation in the two extreme cases of and -type distortions can be found, however, for energy release in the redshift range it is clear that a mixture of both types of distortions will form. In that case analytic approximations are more difficult and numerical studies are required. Furthermore, in the past, mainly the special cases of instantaneous energy release at some initial redshift were numerically studied in detail for this regime (e.g. see Burigana et al., 1995; Burigana & Salvaterra, 2003). However, from the physical point of view cases of continuous energy release over extended periods are more natural. In addition, so far only problems with rather large energy injection, close to the allowed upper bounds obtained by Cobe/Firas, were investigated thoroughly. Again from the observational point of view, significantly smaller distortions seem to be favoured, and in addition the sensitivities reached by Pixie demand studies with very small energy release. To allow forecasting the observational possibilities of Pixie with respect to constraints on the thermal history of our Universe, significantly refined computations of the possible spectral distortions with particular physical processes in mind therefore become necessary.

In this work we revisit the cosmological thermalization problem in the light of Pixie. We formulate the thermalization problem allowing to resolve very small spectral distortions at levels well below the current upper limits of Cobe/Firas. We use improved approximations for the double Compton and Bremsstrahlung emissivities, as well as a detailed treatment of the cosmological recombination problem based on CosmoRec (Chluba & Thomas, 2011). We explicitly compute the expected CMB spectral distortions for different energy injection scenarios. In particular, we discuss two scenarios that are even present in the standard cosmological model. One is connected with the interaction of CMB photon with adiabatically cooling electrons and baryons (see Sect. 2.5.1), while the other is caused by the dissipation of acoustic waves in the early Universe (see Sect. 2.5.2). We also demonstrate that the evolution of spectral distortions at very low frequencies GHz is significantly affected by free-free absorption at late times (), when the Compton process no longer is able to keep the electrons and baryons from cooling strongly below the CMB temperature. For large energy release this modification is not very significant, however, for very small energy release, comparable to and at the level , it becomes important.

One of the outcomes of this work is CosmoTherm888CosmoTherm will be made available at www.Chluba.de/CosmoTherm., an open source code that solves the cosmological thermalization problem for different energy injection scenarios in the limit of small distortions. In the future CosmoTherm should be useful for detailed forecasts of possibilities to constrain the thermal history of our Universe by measuring the CMB energy spectrum.

## 2 Formulation of the thermalization problem

A solution of the thermalization problem in the isotropic Universe generally can only be obtained numerically. Under simplifying assumptions several analytic approximations exist, however, these approximations have their limitations. For example, for very large energy injection, non-linear aspects of the problem become important, while for very small energy release, the precise shape of the distortion depends strongly on the small difference in the electron and photon temperature, which is comparable to the distortion itself, and is more difficult to account for in analytic approximations. Below we provide the required equations for the formulation of the problem that allows to include all these aspects consistently.

### 2.1 General aspects of the thermalization problem

To thermalize spectral distortions after some significant energy release in the early Universe two main ingredients are needed: (i) the number of photons has to be readjusted and (ii) photons have to be redistributed over frequency to restore the shape of a blackbody.

At redshifts electrons and positrons were very abundant and electron-positron, electron-electron and positron-positron Bremsstrahlung were orders of magnitude more efficient in producing soft photons than any other process. Hence thermalization of CMB spectral distortions was extremely rapid. Therefore we shall restrict ourselves to redshifts below . By then the number densities of electrons and positrons had frozen out999More precisely the freeze-out of positrons happens at temperatures or redshift . and only insignificant amounts of positrons were left from the era of electron-positron annihilation. During the subsequent evolution of the (isotropic) Universe the main interactions between the photons and baryonic matter are governed by Compton scattering (CS), double Compton scattering (DC) and normal electron-ion Bremsstrahlung (BR). Furthermore due to the expansion of the Universe the photons suffer from redshifting and the non-relativistic electrons and baryons cool adiabatically.

The photons do not interact directly with the baryonic matter (hydrogen, helium and the nuclei of other light elements), but only indirectly mediated by electrons, which themselves are strongly coupled to the baryons via Coulomb scattering. The timescale, on which the electrons and baryons adjust their energy distributions is much shorter than any other timescale of importance here. Therefore it can be assumed that the electrons and baryons always follow a (relativistic) Maxwell Boltzmann distribution with one temperature . In this case the problem can be formulated with the Boltzmann equations for the evolution of the photon phase space distribution function and a coupled equation describing the time evolution of the electron (and baryon) temperature. In addition the expansion timescale, , where is the Hubble factor, is affected by the energy injection process, however, as we argue in Sect. 2.4 the corresponding correction is negligible. Below we shall write down the photon Boltzmann equation in the expanding isotropic Universe and the evolution equation for the electron temperature with inclusion of heating and cooling.

### 2.2 Evolution of the photons in the expanding Universe

In the early Universe the photons undergo many interactions with the electrons in the Universe. As mentioned earlier, the most important processes are Compton scattering, double Compton scattering, Bremsstrahlung and the adiabatic expansion of the Universe. Among these processes for most times Compton scattering is the fastest. Therefore it is convenient to express all the involved timescales in units of , the Thomson scattering time.

The Boltzmann equation for the evolution of the photon occupation number, , then can be expressed as:

 ∂nν∂τ−HtCν∂nν∂ν=dnνdτ∣∣∣C+dnνdτ∣∣∣DC+dnνdτ∣∣∣BR, (1)

where we introduced the optical depth to electron scattering as the dimensionless time variable. The second term on the left hand side is due to the expansion of the Universe and the right hand side terms correspond to the physical processes quoted above.

For computational purposes it is also convenient to transform to dimensionless frequency. Introducing , where , the expansion term in Eq. (1) can be absorbed so that photons no longer redshift out of the chosen computational domain. Here one point is very important: the temperature can really be any temperature, as long as it remains constant in time. One obvious choice is the observed CMB temperature today, K (Fixsen, 2009), but we will use this freedom to fix our initial conditions conveniently (see Sect. 3.1.1).

For , we have the standard redshift dependent CMB temperature, . However, it is important to mention that with this definition is not necessarily identical to the effective (thermodynamic) temperature of the photon field,

 T∗γ =Tz[G3Gpl3]1/4=Tz[15G3π4]1/4 (2)

which is found by comparing the total energy density of the photon distribution with the one of a pure blackbody. Here is defined by the integral

 G3 =∫x3nxdx=Gpl3+∫x3Δnxdx, (3)

over the photon occupation number, . For a blackbody with thermodynamic temperature one has , and hence , trivially implying that . However, in the presence of spectral distortions after release of energy generally . This for example means that at higher redshifts the effective temperature of the CMB could have been slightly lower than today, and only because of energy injection and subsequent thermalization it has reached its present value. Furthermore, it is important to mention that is known with precision mK (Fixsen & Mather, 2002; Fixsen, 2009), so that deviations of the effective temperature from which are much smaller than this are presently indistinguishable. We will return to this aspect of the problem again in Sect. 3.1.1.

In Eq. (3) we also introduced the distortion, , from a blackbody with temperature . For numerical computations it in general is useful to substitute in this way, as this allows one to cancel the dominant terms and thereby linearize the problem, by only solving for the (small) correction .

#### 2.2.1 Compton scattering

The contribution of Compton scattering to the right hand side of the photon Boltzmann equation (1) can be treated using the Kompaneets equation (Kompaneets, 1956):

 dnxdτ∣∣∣C =θex2∂∂xx4[∂nx∂x+ϕnx(nx+1)] . (4)

Here we introduced the abbreviations and . This equation describes the redistribution of photons over frequency, were the effect of Doppler diffusion and boosting, electron recoil, and stimulated electron scattering are accounted for.

The Kompaneets equation was found as lowest order Fokker-Planck expansion of the collision term for Compton scattering. We neglect higher order relativistic corrections (e.g., see Rephaeli, 1995; Itoh et al., 1998; Challinor & Lasenby, 1998; Sazonov & Sunyaev, 1998; Nozawa et al., 2006). These are expected to affect the results at the level of a few percent at (Chluba, 2005), so that for the purpose of this paper Eq. (4) will be sufficient.

One can easily verify, that the number of photons is conserved under Compton scattering by multiplying Eq. (4) with and integrating over frequency using integration by parts and the fact that the photon distribution vanishes sufficiently fast for and . Furthermore, one can easily check that for photons, which follow a Bose-Einstein distribution, , with constant chemical potential , Eq. (4) vanishes identically if . For in the chosen coordinates the photon distribution in kinetic equilibrium with respect to Compton scattering is with frequency-dependent chemical potential , where we defined .

If we now insert into Eq. (4) and neglect terms of , we obtain the linear equation

 dΔnxdτ∣∣∣C =Δθex2∂∂xx4ζ+θex2∂∂xx4[∂Δnx∂x+ϕξΔnx], (5)

where for convenience we have introduced

 ζ(x) =−nPl(nPl+1)=−e−x[1−e−x]2≡12∂ξ(x)∂x (6a) ξ(x) =2nPl+1=1+e−x1−e−x (6b)

and , with . One can see from Eq. (5) that the remaining part arising from the CMB background spectrum appears as a source term in the evolution equation, which vanishes identically unless .

Explicitly carrying out part of the derivatives in Eq. (5) and rearranging terms we find

 dΔnxdτ∣∣∣C =De∂2Δnx∂x2+De[4x+ϕξ]∂Δnx∂x +Deϕξ[4x+∂lnξ∂x]Δnx−DeΔϕζ[4x−ξ] (7)

with , , and .

For numerical computations it is convenient to precompute the functions , , , and once the frequency grid is chosen. This accelerates the computation significantly. Furthermore, it is important to use series expansion of the expressions for , in particular when differences are encountered. This is very important to achieve numerical stability and the correct limiting behaviour, and will also be crucial for the formulation of the emission and absorption term in the next section.

#### 2.2.2 Double Compton scattering and Bremsstrahlung

In the cosmological thermalization problem, the double Compton process and Bremsstrahlung provide the source and sink of photons. Their contribution to photon Boltzmann equation, Eq. (1), can be cast into the form (Rybicki & Lightman, 1979; Lightman, 1981; Danese & de Zotti, 1982; Chluba, 2005):

 dnxdτ∣∣∣DC+BR =1x3[1−nx(exe−1)]×K(x,θz,θe) (8)

where the emission coefficient is given by the sum of the contribution due to double Compton scattering and Bremsstrahlung, . Since drops off exponentially for (see below), the main emission and absorption of photons is occurring at low frequencies. At small the photon distribution after a very short time is pushed into equilibrium with , i.e., a blackbody of temperature . If once again we insert and rearrange terms, we find

 dΔnxdτ∣∣∣DC+BR =1−eΔxe1−e−x×K(x,θz,θe)x3 +Δnx(1−exe)×K(x,θz,θe)x3, (9)

with . Note that the first term vanishes if , however, the second is pushed to equilibrium only because of interplay with the first term: if a distortion is present then consequently and hence the first term either leads to net absorption () or emission () until (at particular frequencies) the distortion term is balanced.

#### Double Compton scattering

Due to the large entropy of the Universe at sufficiently high redshifts () DC emission dominates over BR (Danese & de Zotti, 1982). The DC emission coefficient can be given as (Lightman, 1981; Thorne, 1981; Chluba et al., 2007)

 KDC(x,θz,θe) =4α3πθ2z×gdc(x,θz,θe), (10)

where is the fine structure constant and is the effective DC Gaunt factor. In lowest order of the photon and electron energies the DC Gaunt factor factorizes (see Chluba, 2005, for more details). Furthermore, if the photon distribution is not too far from full equilibrium one can approximate using a blackbody ambient radiation field and assuming that . In this case one has (e.g., see Chluba, 2005; Chluba et al., 2007)

 gdc(x,θz,θe) ≈Ipl41+14.16θz×Hdc(x), (11)

where . Here we have included the first order relativistic correction in the photon temperature, however, this term only becomes significant at .

The second factor in Eq. (11) allows to go beyond the soft photon limit, for which was assumed. In lowest order only depends on the ambient photon distribution, but is independent of the electron temperature. It can be computed using (see Chluba, 2005, for more details)

 Hdc(x) ≈1Ipl4∫∞2xx′4nPl(x′)[1+nPl(x′−x)][xx′HG(xx′)]dx′ (12)

where with . The factor was first obtained by Gould (1984) to describe the corrections to the DC emissivity when going beyond the soft photon limit but assuming resting electrons101010Note that is of given by Eq. (27) of Gould (1984). The factor of 2 is to avoid double counting of photons.. In the limit one finds , so that .

Expression (12) was also used in the work of Burigana et al. (1991b). There the approximation was given. However, as mentioned above with the assumptions leading to Eq. (12) the electron temperature is irrelevant, and hence one should set . Furthermore, we reexamined the integral and found that for background photons that follow a blackbody spectrum

 Hpldc(x) ≈e−2x[1+32x+2924x2+1116x3+512x4] (13)

provides a much better approximation to the full numerical result for (cf. Fig 1). This approximation was obtained by replacing and neglecting the induced term in Eq. (12). Furthermore, the resulting expression was rescaled to have the correct limit for . In particular, for Eq. (13) captures the correct scaling . However, since most of the photons are produced at low frequencies we do not expect any big difference because of this improved approximation. Nevertheless, when using the old approximation we found that at early times the spectrum is erroneously brought into full equilibrium at very high frequencies, just by DC emission and absorption.

We note here that if the distortions are not small, then in lowest order the correction to the DC emission can be accounted for by replacing with the solution in the expression for . However, from the observational point of view it seems unlikely that distortions of interest ever exceeded the level , even at . Therefore, the above approximation should be sufficient. Of course this does not include DC emission from very high energy photons that are directly related to the energy injection process. However, in that case the simple approximation used above will anyhow need revision, although the total contribution to the photon production is still expected to be small.

#### Bremsstrahlung

At lower redshifts () Bremsstrahlung starts to become the main source of soft photons. One can define the Bremsstrahlung emission coefficient by (cf. Burigana et al., 1991b; Hu & Silk, 1993a)

 KBR(x,θe) =αλ3e2π√6πθ−7/2ee−xϕϕ3∑iZ2iNigff(Zi,x,θe). (14)

Here is the Compton wavelength of the electron, and are the charge, the number density and the BR Gaunt factor for a nucleus of the atomic species i, respectively. Various simple analytical approximations exists (Rybicki & Lightman, 1979), but nowadays more accurate fitting formulae, valid over a wide range of temperatures and frequencies, may be found in Nozawa et al. (1998) and Itoh et al. (2000). In comparison with the expressions summarized in Burigana et al. (1991b) we find differences at the level of for small .

In the early Universe only hydrogen and helium contribute to the BR Gaunt factor, while the other light elements can be neglected. In the non-relativistic case the hydrogen and helium Gaunt factors are approximately equal, i.e., to within a few percent. Therefore assuming that the plasma is still fully ionized the sum in Eq. (14) may be simplified to , where is the baryon number density. However, for percent accuracy one should take the full expressions for and into account, which does not lead to any significant computational burden using the expressions of Itoh et al. (2000).

Furthermore, at redshifts the plasma enters the different epochs of recombination. Therefore, the mixture of the different species (, Hi, Hii, Hei, Heii, Heiii) in the primordial medium has to be followed. We use the most recent computations of the recombination process including previously neglected physical corrections to the recombination dynamics according to Chluba & Thomas (2011).

### 2.3 Evolution of the ordinary matter temperature

As mentioned above, the baryons in the Universe will all follow a Maxwell-Boltzmann distribution with a temperature that is equal to the electron temperature. Caused by the Hubble expansion alone, the matter (electrons plus baryons) temperature would scale as (Zeldovich et al., 1968), implying . However, the electrons Compton scatter many times off background photons and therefore are pushed very close to the Compton equilibrium temperature

 Teqe =I44G3Tz=Tz+[ΔI44G3−ΔG3G3]Tz, (15)

within the (distorted) photon field, where . If the CMB is undistorted, and , however, in the presence of spectral distortions generally . Here we also defined and , which for numerical purposes is better, as the main terms can be cancelled out.

In addition, when spectral distortions are present, the matter cools/heats via DC and BR, and if some significant energy is released, this will in addition heat the medium. Putting all this together the evolution equation for the electron temperature reads (e.g., see Chluba, 2005, for a detailed derivation)

 dρedτ =tC˙Qαhθγ+4~ργαh[ρeqe−ρe]−4~ργαhHDC,BR(ρe)−HtCρe. (16)

Here we introduced ; the heat capacity of the medium111111We neglected relativistic corrections to the heat capacity of the electrons, which at would be of order percent (Chluba, 2005)., ; the energy injection rate, , which for example could be caused by some decaying particles (see Sect. 2.5 for more details); and the energy density of the photon field in units of electron rest mass, , with . Furthermore, we defined . For numerical purposes it is better to rewrite the above equation in , so that this formulation for is very convenient.

The integral arises from the cooling/heating by DC and BR, and reads

 HDC,BR =14G3θz∫[1−nx(exe−1)]×K(x,θz,θe)dx. (17)

For numerical purposes it again is useful to group terms of similar order, when calculating this integral numerically. Usually, does not contribute very significantly to the total energy balance, but for energy conservation over the whole history it is crucial to obtain precise values for it.

#### 2.3.1 Approximate solution for Te

Because the timescale for Compton scattering is extremely short until , at high redshifts one can solve Eq. (16) assuming quasi-stationarity. Even at low redshifts this provides a very good first approximation for the correct temperature. Our full numerical computations confirm this statement. Since the cooling from DC and BR is very small, it is furthermore possible to write , so that

 ρi+1e =ρeq,∗e1+∂ρieHDC,BR(ρie)+Λ (18a) ρeq,∗e =ρeqe+tC˙Q4~ργθγ−HDC,BR(ρie)+∂ρieHDC,BR(ρie)ρie (18b) Λ =HtCαh4~ργ=381+fHe+XeσTXec~ργH (18c)

where is the solution for at time-step . Also note that here we assumed that the spectral distortion from the time-step is used in all the integrals that have to be taken over the photon distribution. This will lead to some small error, which can be controlled by choosing an appropriate step-size. Furthermore, as we will explain below, one can iterate the system at fixed time until the solution converges.

The derivative of with respect to can be computed numerically, however, it is strongly dominated by the term connected with the exponential factor in the definition Eq. (17). For numerical purposes we in general computed all integrals over the photon distribution by splitting the spectrum up into the CMB blackbody part and the distortion. This allows to achieve very high accuracy for the deviations form the equilibrium case, a fact that is very important when tiny distortions are being considered. At low frequencies we again used appropriate series expansions of the frequency dependent functions.

### 2.4 Changes in the expansion history and global energetics

Changes in the thermal history of the Universe generally also imply changes in the expansion rate. For example, considering a massive decaying relic particle, its contribution, , to the total energy density, , will decrease as time goes by. The released energy is then transferred to the medium and because of Compton scattering bulk121212Since the temperature dependent contribution to the total matter energy density is negligibly small, well before recombination virtually all the injected energy will be stored inside the photon distribution. of it ends up in the CMB background, potentially leaving a spectral distortion, and the matter with an increased temperature. Because photons redshift, this implies that the expansion factor changes its redshift dependence, and hence the relation between proper time and redshift is altered.

To include this effect, one has to consider the changes in the energy density of the photons (or more generally all relativistic species) and the relic particles (or more generally the source of the released energy). For simplicity we shall assume that the energy is released by some non-relativistic massive particle, but this can be easily relaxed. We then have

 1a4da4ργdt ≈˙Q (19a) 1a3da3ρXdt =−˙Q, (19b)

where is the energy density of the photons, and is the scale factor.

We emphasize that here it was assumed that the matter temperature re-adjusts, but that the corresponding change in the total energy density of the normal baryonic matter is negligible. Furthermore, Eq. (19) implies that all the injected energy ends up in the photon distribution. For energy injection at very late times, close to the end of recombination, this is certainly not true anymore. In that case, it will be important to add the electron temperature equation to the system, and then explicitly use the Compton heating term for the photons. In this way one has a more accurate description of the heat flow between electrons and photons, and only the spectral distortion is neglected in the problem. However, in terms of the global energetics this should not make any large difference. As we will see below, in the case of adiabatically cooling electrons it turns out that such treatment is necessary in order to define the correct initial condition (see Sect. 3.3).

Equation (19) has the solutions and with

 Δργ(z) =−1a4(z)∫z0a4(z′)˙Qdz′H(z′)[1+z′] ≡−ρ(0)γ(z)∫z0˙Qρ(0)γ(z′)dz′H(z′)[1+z′] (20a) ΔρX(z) =1a3(z)∫z0a3(z′)˙Qdz′H(z′)[1+z′] (20b)

where is the unperturbed photon energy density, and similarly . Here , with , denotes the present day values of the corresponding energy densities. For decaying relic particles with lifetimes shorter than the Hubble time this consequently means . We will use Eq. (20) later to compute the appropriate initial condition for the energy density of the photon field (see Sect. 3.1.1).

The integrals in Eq. (20) themselves depend on , however, assuming that the total change in the energy density is small, one can use , where denotes the unperturbed expansion factor with . Then the modified expansion factor is simply given by with

 ΔH(z)H(0)(z) ≈Δργ(z)+ΔρX(z)2ρ(0)tot(z) =ρ(0)γ(z)2ρ(0)tot(z)∫z0z′−z1+z˙Qρ(0)γ(z′)dz′H(0)(z′)[1+z′] (21)

where is the total energy density of the unperturbed Universe. This equation shows clearly that the change in the total energy density is arising from the fact that energy released by non-relativistic particles is transformed into contributions to the photon energy density: if photons had the same adiabatic index as matter this would not change anything, a conclusion that can be easily reached when thinking about energy conservation. However, the additional redshifting of photons implies a small net change in the expansion rate.

In terms of the thermalization problem, the Hubble factor explicitly only appears in the evolution equation for the matter temperature, Eq. (16), while otherwise it is only present implicitly, via . For it can be neglected, since this will only lead to an extremely small shift in the relation between proper time and redshift that corrects the small spectral distortion, . This conclusion can be reached as well when transforming to redshift as time variable, since then appears in all terms for the photon Boltzmann equation, where is a correction-to-correction. A similar reasoning holds for the matter temperature equation. There one could in principle replace with , however, looking at the quasi-stationary solution for , Eq. (18), it is evident that this will only lead to a minor difference. We checked this statement numerically and found the addition of to be negligible. We conclude that the correction to the Hubble factor enters the problem in second order.

We also comment that above we made several simplifying assumptions about the form of the energy release. For example, if all the energy is released in form of neutrinos, then both the CMB spectrum and the matter temperature remain unaffected, however, the Hubble factor will again change, assuming that the energy density of the source of the energy release itself does not scale like . However, considering this problem in more detail is beyond the scope of this paper and we do not expect major changes in our conclusions for this case.

### 2.5 Energy injection rates for different processes

Having a particular physical mechanism for the energy release in mind, one can specify the energy injection rate associated with this process. Here we assume that most of the energy is going directly into heating of the medium, but only very few photons are produced around the maximum of the CMB energy spectrum, or reach this domain during the evolution of the photon distribution. Also, here we do not treat the possible up-scattering of CMB photons by the decay products, nor do we allow for Bremsstrahlung emission from these particles directly in the regime of the CMB. For a consistent treatment of the thermalization problem it will be very interesting to consider these aspects in more detail, however, this is far beyond the scope of this paper.

With these simplifications it is sufficient to provide in Eq. (16) for each of the processes and then solve the thermalization problem assuming full equilibrium initial conditions. Below we first discuss the effect of adiabatically cooling electrons and the energy injection caused by dissipation of acoustic waves in the early Universe. Furthermore, we give for annihilating and decaying particles, as well as for short bursts of energy injection. The former two processes are present even for the standard cosmological model, while the latter three are subject to strong uncertainties and require non-standard extensions of the physical model.

#### 2.5.1 Adiabatically cooling electrons and baryons

As mentioned above, if the interactions with the CMB are neglected, the temperature of the matter scales as once the baryons and electrons become non-relativistic. At very low redshifts () it is well known that this difference in the adiabatic index of baryonic matter and radiation leads to a large difference in the CMB and electron temperature (Zeldovich et al., 1968) with , which is also very important for the formation of the global 21cm signal from the epoch of reionization (e.g., see Pritchard & Loeb, 2008, and references therein). However, at higher redshifts Compton scattering couples the electrons tightly to the CMB photons, implying that the baryonic matter in the Universe must continuously extract energy from the CMB in order to establish until Compton scattering eventually becomes very inefficient at . Consequently, this should lead to a tiny spectral distortion of the CMB, tiny because the heat capacity of the CMB is extremely large in comparison to the one of matter.

A simple estimate for this case can be given by equating the Compton heating term with the term from the adiabatic Hubble expansion: . This term is already included in Eq. (16), so that no extra has to be given. With this the total extraction of energy from the CMB is roughly given by

 Δργργ∣∣∣cool ≈−αhθγ~ργ,0ln[1+zs1+ze]≈−5.6×10−10ln[1+zs1+ze]. (22)

For the evolution from until this implies a relative energy extraction of , which is only a few times below the claimed sensitivity of Pixie (Kogut et al., 2011).

Here it is interesting to observe that well before the recombination epoch the effective cooling term caused by the electrons has the same redshift dependence as annihilating matter (see below), i.e. . This implies that the type of distortion will be similar in the two cases, however, the distortion caused by the adiabatic cooling of ordinary matter has opposite sign (see Fig. 3), and tend to decrease each other.

Furthermore, with Eq. (26) one can estimate the annihilation efficiency that is equivalent to the matter cooling term. It turns out that if some annihilating relic particle is injecting energy with an efficiency then the net distortion introduced by both processes should practically vanish. As we will show below, our computations confirm this estimate (see Fig. 8), demonstrating the precision of our numerical treatment.

#### 2.5.2 Dissipation of acoustic waves

Acoustic waves in the photon-baryon fluid dissipate part of their energy because of diffusion damping (Silk, 1968). As first estimated by Sunyaev & Zeldovich (1970a), and later refined by Daly (1991) and Barrow & Coles (1991), this effect should lead to a spectral distortion in the CMB which depends on the spectral index of scalar perturbation, , and therefore could potentially allow constraining inflationary models (Mukhanov & Chibisov, 1981) using future CMB spectral measurements (Kogut et al., 2011). Following Hu et al. (1994), the heating rate in a radiation dominated Universe can be cast into the form

 ˙Q∣∣ac =F(nS)~ργH(z)[1+z]3(nS−1)/2, (23)

with131313We assumed that all cosmological parameters aside from have their standard values. Furthermore, as in Hu et al. (1994) we simply used the Cobe normalization for our estimate. for . For this implies a total energy dissipation of between and . This is about 10 times larger than the energy extracted in the case of adiabatically cooling matter (Sect. 2.5.1) over the same period. One therefore expects the associated distortion to be about an order of magnitude larger.

Also, like in the case of cooling matter and annihilating particles one has for . The equivalent effective annihilation rate in this case is . This suggests that the bounds on derived from the CMB anisotropies (Slatyer et al., 2009; Galli et al., 2009; Hütsi et al., 2009; Hütsi et al., 2011) already now might be in tension with this heating mechanism, however, the energy release from acoustic waves only goes into heating of the medium, while for annihilating particles extra ionizing and exciting photons are produced, implying that the effect on the ionization history could be much smaller here.

We therefore computed the modification to the recombination history caused by this process using CosmoRec, and found a small change of around the maximum of the Thomson visibility function (Sunyaev & Zeldovich, 1970b) close to , while the electron fraction in the freeze-out tail indeed was affected by several percent. The changes in the freeze-out tail of recombination could have important consequences for the CMB power spectra, however, at late times, close to recombination, the simple formula Eq. (23) should be modified to take into account deviations from tight-coupling and also the assumption of radiation domination is no longer valid, implying that the rate of energy release is overestimated at those epochs. Furthermore, in this case one should more carefully treat the fraction of energy that directly goes into heating of the matter as opposed to terms appearing directly in the photon field. We believe that the overall effect on the recombination dynamics arising from this process is negligible. Nevertheless, refined computations of the energy release rate might become important for precise computations of this problem, but for the purpose of this paper the above simplification will be sufficient.

#### 2.5.3 Annihilating particles

The rate at which annihilating particles inject energy scales like with the number density of the annihilating particles. To parametrize the rate of energy injection we use (compare Padmanabhan & Finkbeiner, 2005; Chluba, 2010)

 dEdt∣∣∣ann =[1−fν]fannNH(1+z)3, (24)

where has dimensions , and is the fraction of the total energy carried away by neutrinos. is related to the abundance of the annihilating particle, its mass, and the thermally averaged annihilation cross section, however, here we treat all these dependencies with one parameter.

Current measurements of the CMB anisotropies already tightly constrain (Slatyer et al., 2009; Galli et al., 2009; Hütsi et al., 2009; Hütsi et al., 2011). This limit is orders of magnitude stronger than the one obtained from Cobe/Firas (e.g., see McDonald et al., 2001; Chluba, 2010), implying that there is only little room for CMB spectral distortions from annihilating matter, unless a large fraction of the energy is actually released in form of neutrinos, i.e. .

With Eq (24) it then follows

 ˙Q∣∣ann =gh(z)mec2dEdt∣∣∣ann, (25)

where we parametrized the fraction of the total energy going into heating using . Similar to Padmanabhan & Finkbeiner (2005) we will adopt . Here, respectively , and denote the free proton, singly and doubly ionized helium ions relative to the total number of hydrogen nuclei, , in the Universe. At high redshifts, well before recombination starts (), one has , such that , while during hydrogen recombination drops rapidly towards . At low redshifts, a significant part of the released energy is going into excitations and ionizations of ions, rather than pure heating, explaining why (Shull & van Steenberg, 1985; Chen & Kamionkowski, 2004).

To estimate the total amount of energy that is injected relative to the energy density of the CMB one can simply assume radiation domination () and compute the integral . In the case of annihilating matter this yields

 Δργργ∣∣∣ann ≈5.2×10−9[[1−fν]fann2×10−23eVs−1][1−Yp0.75][Ωbh20.022]. (26)

Here we have assumed that is independent of redshift. Furthermore, we assumed that before everything thermalizes (Burigana et al., 1991b), while effectively no energy is transferred to the CMB once recombination ends ().

This estimate shows141414Here we neglect corrections because of partial thermalization of distortions at that for one can expect distortions at the level of . Interestingly, part of the energy is released at very high , where one expects a -type distortion (Illarionov & Sunyaev, 1975a, b), however, also at low redshifts some energy is released. This in contrast should lead to a -type distortion (Illarionov & Sunyaev, 1975a, b; Hu & Silk, 1993a), so that in total one expects a mixture of both. We will compute the precise shape of the distortion for some cases below (see Sect. 3.5), confirming these statements.

#### 2.5.4 Decaying relic particles

For decaying particles the energy injection rate is only proportional to the density of the particle. Motivated by the parametrization of Chen & Kamionkowski (2004) we therefore write

 dEdt∣∣∣dec =f∗XΓXNHe−ΓXt. (27)

Again all the details connected with the particle, its mass and abundance (here relative to the number of hydrogen atoms in the Universe) are parametrized by . The exponential factor arises because the relic particles rapidly disappear in the decay process once the age of the Universe is comparable to its lifetime .

For particles with very short lifetime a large amount of energy could be injected, without violating any of the CMB constraints, since at early times the thermalization process is very efficient (e.g. see Sunyaev & Zeldovich, 1970c; Hu & Silk, 1993b). This implies that in particular early energy release (at ) could still lead to interesting features in the CMB spectrum. Again it will be important to see how the characteristics of the distortion change from -type to -type as the lifetime of the particle increases.

Under the above circumstances the constraints derived from Cobe/Firas to date are still comparable to those from other probes (Kusakabe et al., 2006; Kohri & Takahashi, 2010; Kogut et al., 2011), at least for particles with intermediate lifetimes close to . However, these constraints could be significantly improved with future spectral measurement of the CMB (Kogut et al., 2011). In addition, constraints obtained with the CMB temperature and polarization anisotropies strongly limit the possible amount of energy injection by decaying particles with long lifetimes (Chen & Kamionkowski, 2004; Zhang et al., 2007).

To estimate the total amount of energy injected by the particle we again take the integral . This then leads to

 Δργργ∣∣∣dec ≈10−5[f∗X8×105eV][1−Yp0.75][Ωbh20.022][1+zX5×104]−1J, (28)

where for simplicity we once again assumed radiation domination (implying ). Furthermore, we introduced the redshift corresponding to the lifetime of the particle, , with . In addition, we defined the integral , which for is very close to unity. Ignoring corrections because of the thermalization process151515We will include these in our estimates of Sect. 3.6.1. (e.g., see Hu & Silk, 1993b; Chluba, 2005), Eq. (28) implies that for particles with lifetimes of the limits from Cobe/Firas are for our parametrization. Note that can be found by replacing in Eq. (25) with the expression Eq. (27).

#### 2.5.5 Quasi-instantaneous energy release

As next case we discuss the possibility of short bursts of energy release. We model this case using a narrow Gaussian in cosmological time around the heating redshift with some width :

 dEdt∣∣∣δ =f∗δ[1+zh]4e−[t−th]2/2σh2√2πσh2. (29)

With this parametrization, assuming sufficiently small , one has

 Δργργ∣∣∣δ ≈10−5[f∗δ2.6×10−6eV]. (30)

Computationally it is demanding to follow very short bursts of energy injection, however, our code allows to treat such cases by setting the step-size appropriately. Again can be found by replacing in Eq. (25) with the expression Eq. (29).

With the above formula we can also study the effect of going from quasi-instantaneous to more extended energy release by changing the value of . For short bursts we chose , but we also ran cases with to demonstrate the transition to more extended energy release (see Fig. 19).

## 3 Results for different thermal histories

In this section we discuss some numerical results obtained by solving the coupled system of equations described above. We first give a few details about the new thermalization code, CosmoTherm, which we developed for this purpose. It should be possible to skip this section, if one is not interested in computational details. In Sect. 3.3-3.7 we then discuss several physical scenarios and the corresponding potential CMB spectral distortion.

### 3.1 Numerical aspects

To solve the cosmological thermalization problem we tried two approaches. In the first we simplified the computation with respect to the dependence on the electron temperature. At each time-step the solution for the photon distribution is obtained using the partial differential equation (PDE) solver developed in connection with the cosmological recombination problem (Chluba & Thomas, 2011), which allows us to setup a non-uniform grid spacing using a second order semi-implicit scheme in both time and spatial coordinates. However, we assume that the differential equation for the electron temperature can be replaced with the quasi-stationary approximation, Eq. (18), and solved once the spectrum at is obtained. Since both the evolution of the electron temperature and the spectral distortions is usually rather slow one can iterate the system of equations for fixed in a predictor-corrector fashion, until convergence is reached, and then proceed to the next time.

In the second approach, we modified the PDE solver of Chluba & Thomas (2011) to include the additional integro-differential equation for the electron temperature. This results in a banded matrix for the Jacobian of the system, which has one additional off-diagonal row and column because of the integrals over the photon field needed for the computation of the electrons temperature. Such system can be solved in operations, and hence only leads to a small additional computational burden in comparison to the normal PDE problem. We typically chose grid-points logarithmically spaced in frequency over the range to . However, for precise conservation of energy we implemented a linear grid at frequencies above . Furthermore, we distributed more points in the range to , where most of the total energy in the photon field is stored. Even for about points we achieved very good conservation of energy and photon number. We tested the convergence of the results increasing up to times, as well as widening/narrowing the frequency range.

We then included the coupled PDE/ODE system into the ODE stepper of CosmoRec (Chluba et al., 2010). This stepper is based on an implicit Gear’s method up to six order and allows very stable solution of the problem in time, with no need to iterate extensively. Furthermore, this approach in principle enables us to include the full non-linearity of the problem in , however, for the formulation given above this did not make any difference.

We found that the two approaches described above give basically identical results however the second is significantly faster and more reliable for short episodes of energy injection. In the current version of CosmoTherm this approach was adopted. For typical settings one execution takes a few minutes on a single core, however, for convergence tests one execution with and small step-size () took several hour. Also for cases with quasi-instantaneous energy release a small step-size was necessary during the heating phase, such that the runtime was notably longer.

We also would like to mention that in the current implementation of CosmoTherm we assume that the recombination history is not affected by the energy injection process and therefore can be given by the standard computation carried out with CosmoRec 161616CosmoRec is available at www.Chluba.de/CosmoRec.. This assumption is incorrect at low redshifts for cases with late energy release, as small amounts of energy can have a significant impact on the ionization history. We plan to consider this aspect of the problem in some future work.

#### 3.1.1 Initial condition

We performed the computations assuming full equilibrium at the starting redshift, . Usually we ran our code starting at . However, in order to end up with a blackbody that has an energy density corresponding to the measured CMB temperature, , today, we had to modify the temperature of the photons and baryons at