Recombination Spectrum

# CosmoSpec: Fast and detailed computation of the cosmological recombination radiation from hydrogen and helium

## Abstract

We present the first fast and detailed computation of the cosmological recombination radiation released during the hydrogen (redshift ) and helium ( and ) recombination epochs, introducing the code CosmoSpec. Our computations include important radiative transfer effects, -shell bound-bound and free-bound emission for all three species, the effects of electron scattering and free-free absorption as well as interspecies () photon feedback. The latter effect modifies the shape and amplitude of the recombination radiation and CosmoSpec improves significantly over previous treatments of it. Utilizing effective multilevel atom and conductance approaches, one calculation takes only seconds on a standard laptop as opposed to days for previous computations. This is an important step towards detailed forecasts and feasibility studies considering the detection of the cosmological recombination lines and what one may hope to learn from the photons emitted per hydrogen atom in the three recombination eras. We briefly illustrate some of the parameter dependencies and discuss remaining uncertainties in particular related to collisional processes and the neutral helium atom model.

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

=1

## 1 Introduction

The energy spectrum of the cosmic microwave background (CMB) allows probing the thermal history of the Universe (Zeldovich & Sunyaev, 1969; Sunyaev & Zeldovich, 1970; Burigana et al., 1991; Hu & Silk, 1993), with several interesting standard and non-standard signals awaiting us (for overview see Chluba & Sunyaev, 2012; Sunyaev & Khatri, 2013; Tashiro, 2014). One of the guaranteed signals is imprinted by the cosmological recombination process (Zeldovich et al., 1968; Peebles, 1968; Dubrovich, 1975) at redshifts . It appears as a spectral distortion of the CMB that today is visible at frequencies (Sunyaev & Chluba, 2009). The cosmological recombination radiation (CRR) constitutes a fundamental signal from the early Universe that can provide an alternative way to measure some of the key cosmological parameters (e.g., Chluba & Sunyaev, 2008b; Sunyaev & Chluba, 2009) and allow testing the physics of the cosmological recombination process. As challenging as the detection of the CRR may be (Sathyanarayana Rao et al., 2015; Desjacques et al., 2015; Balashev et al., 2015), it is a highly desirable target for future CMB studies (Silk & Chluba, 2014), and may be possible with an improved version of PIXIE (Kogut et al., 2011) or PRISM (André et al., 2014).

Nowadays, the detailed recombination history (the free electron fraction as a function of redshift) can be computed with precision around decoupling () using CosmoRec (Chluba & Thomas, 2011) and HyRec (Ali-Haïmoud & Hirata, 2011). This precision is required to ensure that our interpretation of CMB data from Planck gives unbiased results in particular for the spectral index of primordial fluctuation and the baryon density (Rubiño-Martín et al., 2010; Shaw & Chluba, 2011; Planck Collaboration et al., 2015a).

While in terms of standard (atomic) physics the recombination history is expected to be well understood, non-standard process, for instance, due to decaying or annihilating particles (e.g., Peebles et al., 2000; Chen & Kamionkowski, 2004; Padmanabhan & Finkbeiner, 2005; Galli et al., 2009a; Hütsi et al., 2009; Slatyer et al., 2009; Giesen et al., 2012) or variations of fundamental constants (e.g., Kaplinghat et al., 1999; Battye et al., 2001; Rocha et al., 2003; Scóccola et al., 2009; Galli et al., 2009b) can modify the ionization history. These scenarios can be constrained with CMB temperature and polarization data (Planck Collaboration et al., 2015a, b). It is, however, hard to discern between different sources of changes to the recombination history. Similarly, when allowing standard parameters such as the helium mass fraction, , and effective number of neutrino species, , to vary simultaneously, large degeneracies remain (Planck Collaboration et al., 2015a). In all these cases, the CRR is expected to provide extra pieces of the puzzle, which will help to break degeneracies and directly probe the recombination dynamics (Chluba & Sunyaev, 2008b; Sunyaev & Chluba, 2009; Chluba, 2010). In addition, the reaction of the plasma in the pre-recombination era can give hints about early energy release (Lyubarsky & Sunyaev, 1983; Chluba & Sunyaev, 2009b). It is thus important to understand how the CRR depends on different assumptions, and a detailed computation of the signals for the standard scenario is the first step.

Several calculations of the CRR have been carried out in the past using different methods and approximations (e.g., Dubrovich, 1975; Rybicki & dell’Antonio, 1993; Dubrovich & Stolyarov, 1995, 1997; Burgin, 2003; Kholupenko et al., 2005; Wong et al., 2006; Rubiño-Martín et al., 2006; Chluba et al., 2007). In Sect. 3, we dwell deeper into the details, but the biggest problem is that these calculations either included a small number of atomic levels, neglected physical processes and/or took prohibitively long to allow detailed computations for many cosmologies. The situation was very similar back in the early days of CMB cosmology when slow Boltzmann codes (Ma & Bertschinger, 1995; Hu et al., 1995) were used for computations of the CMB power spectra.

In this work, we take the next step for computations of the CRR, developing the software package CosmoSpec1, which allows fast and precise computation of the CRR for a wide range of cosmological models including all important physical processes. The computational method is based on a combination of the effective multilevel atom2 (Ali-Haïmoud & Hirata, 2010) and effective conductance approach (Ali-Haïmoud, 2013), which we briefly review in Sect. 3. These approaches are similar in spirit to the line-of-sight approach for CMB anisotropies (Seljak & Zaldarriaga, 1996) used in modern Boltzmann codes CMBfast, CAMB (Lewis et al., 2000) or CLASS (Lesgourgues, 2011), in the sense that they rely on a factorization of the problem into a slow but cosmology-independent part and a fast, cosmology-dependent part.

With CosmoSpec the calculation of the CRR for one cosmology can be performed in seconds on a standard laptop, a performance that opens the path for detailed parameter explorations and experimental forecasts, applications we leave to future work. We discuss the importance of different physical processes, such as interspecies feedback, collisions and uncertainties in the atomic model for neutral helium (Sect. 4). Overall, CosmoSpec in its current state should be accurate at the level of a few percent for the standard recombination scenario, mainly due to uncertainties of the neutral helium atom model.

## 2 Atomic models

For the computations presented here, we adopt atomic models for hydrogen and helium closely following Rubiño-Martín et al. (2006); Chluba et al. (2007) and Rubiño-Martín et al. (2008). Below we briefly summaries the main aspects of these models, but refer to the above references for details.

### 2.1 Hydrogen and hydrogenic helium

Level energies, electric dipole transitions rates, and photoionization cross sections for hydrogenic ions can be computed analytically using the expressions of Karzas & Latter (1961) and the recursion method of Storey & Hummer (1991) (see also Hey, 2006). Corrections due to the Lamb shift are neglected. We also omit Hi and Heii quadrupole transitions, which coincide with the dipole transition frequencies and have been shown to cause a negligible effect on the recombination dynamics (Grin & Hirata, 2010; Chluba et al., 2012) and thus should not change the Hi and Heii recombination radiation significantly.

We use the hydrogen radiative transfer module of CosmoRec (Chluba & Thomas, 2011), resolving the 2s, 2p, 3s, 3p and 3d states. This includes, Lyman- and resonance scattering as well as two-photon emission and Raman scattering events. Electron scattering is added using a Fokker-Planck approach (Hirata & Forbes, 2009; Chluba & Sunyaev, 2009a) to improve the stability of the transfer calculation, but we neglect modifications caused by the detailed shape of the scattering kernel (Ali-Haïmoud et al., 2010). The results of the transfer calculation are used to obtain the high-frequency distortion caused by hydrogen. The Hi 2s-1s two-photon decay rate is set to (Labzowsky et al., 2005).

For hydrogenic helium, we approximate the 2s-1s two-photon profile using the expression of Nussbaumer & Schmutz (1984):

 ϕ(y)=C[w(1−4γwγ)+αwβ+γ4γ], (1)

where , , , , and . The two-photon profile is normalized as , since two photons are emitted per transition. We set the Heii two-photon decay rate to3 . For Heii, we use an effective multilevel atom approach with Sobolev escape probability formalism. No separate radiative transfer calculation is performed.

### 2.2 Neutral helium

The atomic model for helium significantly more complex than the one of hydrogenic ions and, more importantly, more uncertain. The main complication is due to the fact that the neutral helium atom is a three particle system (nucleus + two electrons), so that no closed analytic solutions exist.

#### Energies and transitions

For the energies and transitions rates among the singlet and triplet levels with principal quantum numbers we use the multiplet tables from Drake & Morton (2007). These tables have a few gaps (e.g., for level ), which we fill with hydrogenic values (see Appendix A for details). To obtain the energies of levels with and we use the quantum defect method (e.g., see Drake, 2006), while for all other levels we use the simple hydrogenic approximation for the level energy, , where , is the electron to helium mass ratio4 and is the Rydberg constant. For levels with , we generally do not include -resolved levels. For transitions between levels with , we use hydrogenic transitions rates, rescaled by the second power5 of the transitions frequency. This approximation is rather crude for transitions involving low- states (i.e., S through D) as illustrated below (Sect. 4.6). No singlet-triplet transitions are allowed among excited states with . For transitions from excited levels to , we use hydrogenic values following Rubiño-Martín et al. (2008), which have the correct average property (Appendix A). We also include a few extra transitions (see Chluba & Sunyaev, 2010; Chluba et al., 2012, for more details), for example, singlet-triplet transitions (Switzer & Hirata, 2008a) and additional Hei quadrupole lines (Cann & Thakkar, 2002) to the ground-state. The radiative transfer problem for Hei is solved using CosmoRec (for details see, Chluba et al., 2012).

We approximate the Hei two-photon profile using (Drake, 1986; Switzer & Hirata, 2008a):

 ϕ(y)=Cw3(w+δ)2[α−βw+γw2], (2)

where , , , , , and . Again, the two-photon profile is normalized as . We use (Drake et al., 1969) for the total two-photon decay rate. The effect of other two-photon and Raman processes (Hirata & Switzer, 2008) on the recombination dynamics and shape of the high-frequency distortion are neglected. We also neglect corrections caused by the detailed photon redistribution kernel of electron and line scattering (Switzer & Hirata, 2008b; Chluba et al., 2012).

#### Photoionization cross section

As pointed out previously (e.g., Rubiño-Martín et al., 2008; Chluba & Sunyaev, 2010; Glover et al., 2014), for the CRR one of the biggest uncertainty in the neutral helium model is due to the lack of data for the photoionization cross sections. For levels with and , we use a combination of the sparse TopBase data (Cunto et al., 1993) and pre-calculated recombination rates from Smits (1996). For all other levels, we use hydrogenic approximations for the photoionization cross section, following Rubiño-Martín et al. (2008). The -resolved photoionization cross section is simply set to using the ionization threshold frequencies for hydrogen and helium, and . For the recombination rate this implies a rescaling by . This approximation is obtained by assuming that the -averaged cross section is equal to the hydrogenic one and using the fact that -resolved cross sections are all equal, up to small corrections due to the fine splitting of -sub-levels (see Appendix A.1).6 It also neglects corrections to the shape of the cross section, in particular due to auto-ionization resonances, which are important for the low- states. This leads to uncertainties in the shape of the free-bound spectrum, as we briefly discuss in Sect. 4.6.

## 3 Computational method

### 3.1 Effective multilevel atom and conductance methods

One of the biggest challenges for detailed and fast computation of the cosmological recombination history and radiation is that in explicit multilevel approaches (e.g., Seager et al., 2000; Chluba et al., 2007) the populations of many levels ( for hydrogen) have to be followed to determine their departures from full equilibrium. For 500-shell calculations of hydrogen this includes some levels, which makes the computations very time-consuming (e.g., Rubiño-Martín et al., 2006; Grin & Hirata, 2010; Chluba et al., 2010). However, the effect of excited states can be integrated out using the effective multilevel atom method for the recombination history (Ali-Haïmoud & Hirata, 2010) and the effective conductance method for the recombination radiation (Ali-Haïmoud, 2013).

The effective multilevel atom7 is an exact generalization of the effective three-level atom model (Peebles, 1968; Zeldovich et al., 1968). In this approach, the case-B recombination coefficient is replaced by effective recombination coefficients to a few low-lying excited states. These coefficients account for stimulated recombinations and the finite probability of photoionization as a captured electron cascades down to the lowest excited states. In the absence of collisional transitions, they are functions of the electron () and photon () temperatures only. These coefficients are supplemented by effective photoionization rates and effective transition rates between the low-lying excited states, . The calculation of the effective rates to high accuracy is computationally demanding, as it requires solving large linear algebra problems. However, it has to be done only once. For any cosmology, the recombination history can then be computed very efficiently by solving an effective few-level system, and interpolating the effective rates whenever they are needed. The effective multilevel approach is the basis for the cosmological recombination codes HyRec (Ali-Haïmoud & Hirata, 2011) and CosmoRec (Chluba & Thomas, 2011) which are now directly used for the analysis of high-precision CMB data from Planck (Planck Collaboration et al., 2015a).

The effective conductance method (Ali-Haïmoud, 2013) allows to factorize the computation of the recombination radiation into a slow, cosmology-independent part, and a fast, cosmology-dependent part. In this approach, the net emissivity (or intensity) in any bound-bound or free-bound transition is rewritten as a sum of temperature-dependent coefficients (effective conductances) times departures from equilibrium of the populations of the lowest-lying excited states (which play the role of tensions). The former are pre-tabulated as a function of transition energy and temperature, while the latter are obtained very efficiently, for any particular cosmology, by solving the effective few-level atom described above.

Here, we extend CosmoRec to allow detailed computation of the recombination history using an effective multilevel approach. Similarly, we update the effective rate coefficients for neutral helium to include 500 shells. While the effective rates and conductances for neutral helium have to be explicitly computed, those for hydrogenic helium can be obtained from those for hydrogen by simple rescaling which we derive in Appendix B.

### 3.2 Electron scattering

Line broadening and shifting caused by scattering with free electrons was found to be important for the Heii recombination radiation (Rubiño-Martín et al., 2008). To include these effects we use the improved approximation for the electron scattering Green’s function of the photon occupation number given by Chluba (2015):

 G(x′,x,y) =exp(−[ln(x/x′)−αy+ln(1+x′y)]2/4yβ)√4πyβx3. (3)

Here, and denotes the scattering -parameter. We also introduced and , with , which were obtained by matching the numerical solution. This approximation improves over the classical solution (, and ) of Zeldovich & Sunyaev (1969) by accounting for the effect of electron recoil at high frequencies and including the blackbody-induced stimulated scattering (Chluba & Sunyaev, 2008a) at low frequencies.

To compute the contribution to frequency bin caused by the scattering of photons emitted at frequency and redshift , we find

 ΔIscj =ΔIi(z)2xi[erf(bu)−erf(bl)]e−ai+yβi ≈ΔIi(z)xie−ai+yβi−b2j√4πyβixu−xlxj. (4)

We integrated across the bin between and with . Here, , and denotes the distortion caused by the process without including electron scattering. The second line in Eq. (3.2) was obtained assuming a narrow bin . The -parameter is computed between the emission redshift and today using the solution for the electron recombination history of CosmoRec. Around the maximum of emission from Heii () we find , so that the line broadening reaches . During neutral helium and hydrogen recombination, electron scattering is much less important and can safely be neglected.

### 3.3 Free-free absorption

At low frequencies (), free-free absorption becomes noticeable (Chluba et al., 2007). So far, this effect was only taken into account for the Hi recombination radiation. It is, however, slightly more important for the helium recombination spectra, as we show below (Sect. 4). To include this process, we modify the emission in each frequency bin and redshift by

 ΔIffi =ΔIi(z)e−τff(xi,z), (5)

where the free-free absorption optical depth is roughly given by, , and is a single redshift-dependent function. This function can be computed using approximations for the free-free Gaunt factors (e.g., Itoh et al., 2000; Draine, 2011). The efficiency of free-free absorption drops exponentially around hydrogen recombination (), but earlier it is important at low frequencies. More details can be found in Chluba (2015).

### 3.4 Feedback due to helium photons

The feedback physics is simple: energetic photons, emitted during the two recombination phases of helium, can ionize and excite the next-lying species at a later time. Specifically, Heii photons feedback onto Hei and Hi atoms, and Hei photons feedback onto Hi atoms. These ionizations and excitations typically take place earlier than the epoch of recombination of the lower-lying species, and are instantaneously compensated by rapid recombinations and decays, generating a pre-recombination radiation from the species that is subjected to feedback. In the case of feedback, the pre-recombination emission from Hei can itself feedback onto hydrogen (Chluba & Sunyaev, 2010; Chluba et al., 2012).

Part of the feedback is already included in previous calculation for the hydrogen recombination radiation (Rubiño-Martín et al., 2008). However, as shown in Chluba & Sunyaev (2010), changes in the recombination dynamics due to Heii feedback cause significant additional modifications of the CRR. In Chluba & Sunyaev (2010), these corrections were only included for 20-shell calculations, which here we extend to 500 shells. Furthermore, changes to the detailed time-dependence of the feedback process because of radiative transfer effects (e.g., broadening by electrons scattering of the Heii Lyman- line) were neglected. These effects can be readily included in our new computation. For this we modify the radiative transfer calculation of neutral helium, including the emission of photons in the Heii Lyman- line and 2s-1s two-photon continuum. These photons then redshift and feedback on neutral helium and hydrogen in their respective pre-recombination eras causing the production of extra photons. Our treatment automatically includes corrections due to line broadening by electron and resonance scattering. The results of these calculations will be presented in Sect. 4.2.

### 3.5 Collisional transitions

Collisional processes can modify the atomic level populations of the highly excited states. This affects the recombination dynamics and the low-frequency emission, originating from transitions involving highly excited states (). For hydrogen, collisional processes were found to have a minor effect for both aspects8 (Chluba et al., 2007, 2010), however, because collisions become more important at higher densities, for helium this may no longer be true, in particular for the pre-recombination emission.

#### Angular-momentum-changing collisions

Angular-momentum-changing collisions are the most important for the recombination history and the CRR. To estimate the effect of -changing collision among energetically degenerate levels (modeling them with hydrogenic wave functions), we use the classical approximation of Vrinceanu et al. (2012), which we extend to general ion and nucleus charge:

 Cnℓ→nℓ′ ≈2.6×10−5√Mme(Z2ion/Z2nuc)Nion√Te/K ×n2[n2(ℓ+ℓ′)−ℓ2<(ℓ+ℓ′+2|Δℓ|)](2ℓ+1)|Δℓ|3cm3s−1, (6)

where is the number density of the projectiles, , is the reduced mass of the ion-target system; the charge of the projectile and that of the nucleus ( for hydrogen and neutral helium, and 2 for singly-ionized helium). We give a simple derivation of Eq. (3.5.1) in Appendix C, valid for . We also argue that the underlying assumption of a straight-line trajectory for the colliding ion is very accurate at the relevant temperatures, even when the target is charged as is the case for singly-ionized helium. For -mixing, proton and -particle collisions are more important than electron collisions (which we neglect).

Equation (3.5.1) automatically obeys detailed balance and provides an estimate for transition with9 , in very good agreement with the full quantum calculation (Vrinceanu et al., 2012). It can be applied for all three species, since -changing collisions are only important for highly-excited states (), which even for neutral helium are close to hydrogenic. It does, however, omit the possibility of collision-induced singlet-triplet state mixing.

For , the quantum-mechanical transition probability scales as for large impact parameters, while the classical transition probability used to obtain Eq. (3.5.1) vanishes beyond a maximum impact parameter (see Fig. 1 in Vrinceanu et al., 2012, as well as our Appendix C). To obtain a finite collision rate with the full quantum description, one therefore needs to cut off the integration at some maximum impact parameter (Vrinceanu & Flannery, 2001). Since the divergence is only logarithmic, the choice of the cutoff is not critical and the result is found to be well approximated by Eq. (3.5.1) (Sadeghpour, 2015, private communication). It is lower by a factor of a few than the previous estimate of Pengelly & Seaton (1964), whose Born approximation calculation required a regularization at both large and small impact parameters.

#### Collisional excitations and de-excitations

For collisional de-excitation rates of ions or atoms by electron impact we use the expression of van Regemorter (1962):

 Ci→j≈20.6Neλ3ijAi→j√Te/KP(hνij/kTe), (7)

where and are the the upper and lower levels, respectively, is the transition wavelength, is the Einstein-A coefficient of the transition and the function depends on whether the target is neutral or positive. At high temperatures, it scales as in both cases. The collisional excitation rate is obtained by detailed balance: .

#### Collisional ionizations and three-body recombinations

For the collisional ionization rate by electron impact, we use the simple expression (see, Mihalas, 1978; Mashonkina, 1996):

 Ci→c≈1.55×1013¯gNeσi→c(νic)√Te/Ke−hνic/kTehνic/kTecms−1, (8)

where is the photoionization cross section at threshold, for hydrogen and neutral helium and for hydrogenic helium. The collisional recombination coefficient follows from detailed balance.

The rate coefficients for collisional ionization and excitation are uncertain. However, they turn out to have a very small effect on the CRR, as we shall discuss in Sect. 4.5.

#### Inclusion in the effective rates and conductances

One can simply generalize the effective multilevel and effective conductance methods to account for collisional transitions. The effective rates and conductances are then not only functions of temperature, but also the collider densities, . In principle this adds more dimensions over which to tabulate the effective rates. However, since the cosmological parameters are known to percent accuracy, so is the evolution of the main collider densities as a function of . It would therefore probably be enough to compute the values and first few derivatives of the effective rates with respect to the densities around their best-fit values . However, given that the uncertainty of the collisional rates is still large compared to the precision of the standard cosmological parameters, in this paper we simply compute the effective rates and conductances at the best-fit collider densities .

## 4 Results

After obtaining all the effective rates and conductances, both the computation of the recombination history and recombination radiation boils down to solving a moderate system of coupled ordinary differential and partial differential equations and numerical integrals. This accelerates the calculation by a large amount, making it feasible to perform one full calculation within seconds on a standard laptop. Here, we illustrate the results for the recombination radiation, highlighting some of the new aspects and comparison with previous computations.

### 4.1 Comparison with previous computations

To validate our new treatment, we first compare our results with previous calculations of the recombination radiation. For hydrogen, detailed multilevel recombination computations for the bound-bound emission, including up to 350 -resolved shells ( individual levels), were presented in Chluba et al. (2010). These calculations also included the effect of Hi Lyman-continuum absorption during the Hei recombination era (Kholupenko et al., 2007; Switzer & Hirata, 2008a), an effect that leads to pre-recombination emission from hydrogen (Rubiño-Martín et al., 2008; Chluba & Sunyaev, 2010). In Ali-Haïmoud (2013), the effective conductance method was used to compute the recombination radiation from hydrogen for 500 -resolved shells, also including the free-bound radiation, which previously was only available for -shell hydrogen calculations (Chluba & Sunyaev, 2006a). Comparing with the results from our new computation for similar settings (in terms of the included physics and atomic model), we find excellent agreement at the level of with these previous computations. While this validates our implementation, here we add corrections due to feedback from Heii photons, which increases the emission caused by hydrogen and neutral helium (see Sect. 4.2).

For singly-ionized helium, we compare with the multilevel computations presented in Rubiño-Martín et al. (2008) for 100 shells. Only results for the bound-bound radiation were previously presented, for which we find excellent agreement both with and without electron scattering included. For the neutral helium recombination emission only the bound-bound emission for up to 30 shells was previously considered (Rubiño-Martín et al., 2008). Extending the multilevel code of Chluba et al. (2010) to include more levels for neutral helium, we confirmed the results of our conductance method. For the bound-bound radiation, at low frequencies () our results deviate from the spectrum presented in Rubiño-Martín et al. (2008). We attribute these differences to variations in the atomic model of neutral helium (in particular the recombination rates to excited states), however, the main high-frequency features agree well with our improved calculations. Here we also add the two free-bound contributions from both helium ions, which have not been obtained before.

### 4.2 Contributions from the three recombination eras

In Fig. 1, we illustrate the individual contributions from bound-bound and free-bound transitions of the three atomic species. We only show the emission among levels with , and do not display the high-frequency distortion from the Lyman lines and and 2s-1s continuum emission, though we compute them consistently as as part of the radiative transfer problem. The effect of collisions and Heii feedback were neglected to produce Fig. 1 and will be discussed below.

The hydrogen recombination radiation shows several feature (Balmer, Paschen and Brackett lines) related to the bound-bound transitions (), with the emission merging to a quasi-continuum at low frequencies. Free-free absorption is noticeable at . Although on average at a times lower level, the distortions from the two helium eras also show significant structure. For the Heii spectrum, the effect of electron scattering is important, but nevertheless distinct features remain visible. The total Hei spectrum, shows absorption features at and . The first feature is related to the 10 830 line () and the second to absorption in the 3D-2P triplet-triplet transitions (see Rubiño-Martín et al., 2008, for more details). For the helium contributions, free-free absorption becomes significant at , exponentially cutting off all recombination emission at .

In Figure 2, we show the total (bound-bound and free-bound) recombination radiation from hydrogen and helium for our detailed 500-shell calculations. We also included the high-frequency distortion, which for the helium spectra is strongly absorbed and reprocessed, reappearing as pre-recombination hydrogen Lyman- emission (Rubiño-Martín et al., 2008; Chluba et al., 2012). Overall, the effect of helium is visible in several bands, introducing features, shifts in the line positions and enhanced emission. In particular, when adding all feedback corrections (see below for more discussion), the effect of helium becomes significant, in some parts well in excess of the expected average contribution. This will help significantly when attempting to use the CRR to determine the pre-stellar (un-reprocessed) helium abundance.

### 4.3 Effect of Heii feedback

As explained in Sect. 3.4, feedback caused by high-frequency photons emitted during the recombination era changes the CRR of hydrogen and neutral helium by ionizing neutral atoms in their pre-recombination eras (Chluba & Sunyaev, 2010). In Fig. 3, we illustrate the effect on the hydrogen bound-bound spectrum. Since the number of helium atoms is lower than the number of hydrogen atoms, the effect is not as large, however, additional features appear in the recombination spectrum, e.g., around and . The dominant additional feedback is due to the reprocessing of the Heii Lyman- line by Hei continuum absorption. This causes pre-recombination emission of Hei which then feeds back in the Hi Lyman-continuum around (Chluba & Sunyaev, 2010). For comparison, we also show the case without any feedback at all (dashed/blue) and including primary Hei photon feedback (dotted/violet). The latter case already goes beyond the calculations of (Rubiño-Martín et al., 2008), which neglected the feedback from the Hei two-photon continuum.

The effect on the Hei recombination spectrum is much more dramatic (see Fig. 4). This is because in contrast to the feedback for hydrogen there is roughly one feedback photon per helium atom. The feedback occurs mainly around and is dominated by the Heii Lyman- (Chluba & Sunyaev, 2010). Overall, Heii feedback enhances the total number of photons emitted by helium itself from to . We included all free-bound and bound-bound emission of helium atoms in this estimate. Without Heii feedback, hydrogen emits roughly , which with feedback increases to . In total, are emitted in the standard recombination process, with contributed by helium alone.

### 4.4 Derivatives with respect to Yp, Neff and Ωbh2

As illustrated in Chluba & Sunyaev (2008b), the CRR directly depends on parameters like the CMB monopole temperature, , the baryon density, , and the expansion rate today, . Similarly, the helium mass fraction, , and other parameters affecting the expansion rate, such as the cold dark matter density, , and effective number of neutrino species, , modify the CRR. Although detailed parameter forecasts and consideration of optimal experimental design will be left for future work, here we illustrate how the CRR depends on and . Collisions are neglected in the calculation.

In Fig. 5 we show the logarithmic derivatives of the CRR with respect to (keeping a flat Universe by varying ), and . For comparison, we also show the CRR itself. The logarithmic derivative with respect to has a similar amplitude and shape as the CRR. This is expected since to leading order the total emission is simply proportional to the number of baryons in the Universe, i.e., . The differences relative to the CRR are due to changes in the recombination dynamics and at low frequencies the efficiency of free-free absorption.

The logarithmic derivatives with respect to and show a much richer structure, with both positive and negative features. The CCR is much less sensitive to changes of , since it only enters through the expansion rate. The helium abundance changes both the relative contribution of the two helium spectra to the hydrogen spectrum but also affects the recombination dynamics through feedback and modifications of the recombination rate. Naively, neglecting these dynamical aspects one would simply expect the total CRR to be the sum of the hydrogen and helium templates, which gives fewer features in the derivative (Fig. 6). Since the derivatives with respect to , and have very different shapes, it is clear that the CRR will help break degeneracies between these parameters. To what extent and with which strategy this may be possible will be studied in detail in a future paper.

### 4.5 Effect of collisions

The main effect of -changing collisions is to push the angular-momentum substates within a given shell into statistical equilibrium. This enhances the very low-frequency emission () while slightly lowering the emission at higher frequencies. This is because radiative recombinations are more efficient to the low- states (), which in the absence of collisions are slightly overpopulated in comparison to the high- states. Low- states depopulate via large transitions, emitting at high frequencies, whereas transitions among high- states, due to dipole selection rules, are restricted to smaller jumps, emitting at lower frequencies. The -changing collisions tend to level the low-/high- population imbalance and hence enhance the low-frequency emission (see also, Chluba et al., 2007).

Collisional excitations and de-excitations push the relative populations of excited levels closer to Boltzmann equilibrium at the gas temperature, which is very close to the CMB temperature at the relevant redshifts. This therefore lowers the net radiative transition rate, which is proportional to the departure of the level populations from Boltzmann equilibrium with one another.

At low , the populations are closer to Boltzmann equilibrium relative to the lowest () excited states, while at high they are closer to Saha equilibrium with the continuum. Collisional ionizations and three-body recombinations drive the highly-excited states even closer to Saha equilibrium with the continuum, steepening the level of departures from Saha equilibrium from low (largest departure) to high (smallest departure), leaving the low- states mostly unaffected. This has two effects. First, it reduces the net free-bound emission, which is proportional to the population departure from Saha equilibrium. Second, it enhances the bound-bound emission, which depends on the rate of change of population departures from equilibrium as a function of .

All these effects are most significant at the low frequencies where the highly excited states radiate. They are also slightly more important for helium than for hydrogen, due to higher densities at the relevant times. In Fig. 7, we illustrate these effects using 100-shell calculations. We find -changing collisions and collisional ionization to be subdominant, even for the Heii recombination radiation. Increasing these collision rates by two orders of magnitude makes their effect visible at low frequencies (see figure). The mixing of sub-levels produces a correction to the low-frequency recombination radiation, which for the Hei spectrum becomes slightly more noticeable (in the total CRR of Hei as opposed to for Hi and Heii). However, in all cases the emission at remains practically unaffected.

In Fig. 8, we show the correction to the total CRR caused by collisions. We vary the number of shells in the calculation, to illustrate the changes at low and high frequencies. At low frequencies (), the effect reaches percent-level, while at higher frequencies it remains , nearly independent of . This shows that collisional processes need to be taken into account for precise computations of the recombination radiation. For baseline calculations we neglect collisions, which implies an uncertainty at percent level for the total CRR.

### 4.6 Uncertainties due to the neutral helium model

As mentioned above, the atomic model for neutral helium is not as accurate as those of hydrogen and singly-ionized helium. This implies that the neutral helium contributions to the CRR is also uncertain. Here, we illustrate some of the crucial aspects and estimate the level of the uncertainty.

At this point, the level energies are probably the most certain. For , we compared the hydrogenic approximation and quantum defect energies with the data from Drake & Morton (2007). The biggest deviations are found for the S and P states, which at depart from the hydrogenic values by . The energy values are matched very well by the quantum defect calculation, which is used to extrapolate to higher levels. For , we find levels with to depart from the hydrogenic values by .

A larger uncertainty is found for the dipole transition rates among low- level involving S, P and D states. In Fig. 9 we compare the transitions rates from Drake & Morton (2007) for the singlet P-S and S-P series with the hydrogenic approximation, . This approximation works quite well for the P-1S sequence, converging to the hydrogenic approximation at the level for . However, for the other cases, the hydrogenic approximation is significantly off, overestimating the transition rates for the P-S sequence roughly twice and underestimating the S-P sequence by a factor of around . For the triplet levels, we find a similar mismatch. For singlet transitions involving D states, the hydrogenic approximation becomes accurate at the level of . Transition series involving only higher states approach hydrogenic values at the percent level. Overall, our analysis implies that low- transitions starting with and ending at are not well approximated using hydrogenic values and an improved neutral helium model including many shells is ultimately needed.

As mentioned above, the photonionization cross sections are probably the largest source of uncertainty in the neutral helium model. This is because the hydrogenic approximation neglects corrections to the shape of the cross section from the non-hydrogenic character of the wave-functions and due to auto-ionization resonances, which are important for the low- states (S, P and D). This is illustrated in Fig. 10, where we compare the TopBase cross section with the hydrogenic approximation for some examples. We find that even for , the cross sections for the S and P states show clear non-hydrogenic character already close to the threshold frequency. Although we do not expect auto-ionization resonances to significantly contribute to the total photoionization and recombination rates (which could have an effect on the recombination dynamics) they affect the precise shape of the Hei free-bound continuum. Similarly, for the S and P states we find the slopes of the cross section to be more shallow than the hydrogenic approximation. For the D states, the match is found to be much better already, but auto-ionization resonances are still present even there. For the triplet P states, TopBase cross sections only exist for and . Thus, the shape of the free-bound continuum for is likely affected at a significant level. For , we also expect the free-bound continuua to show dependence on the shortcomings of the hydrogenic approximation. A improved quantum mechanical model of neutral helium will therefore eventually be needed for highly accurate predictions of the CRR.

## 5 Conclusions and future work

We carried out the first comprehensive computation of the CRR from hydrogen and helium including all relevant processes for 500-shell atoms. Our computations improve over previous calculations in terms of the total number of included shells and the modeling of feedback processes. For the first time, we include the free-bound emission from helium and model the effect of collisional processes at high redshift. We furthermore overcome the performance-bottleneck of previous treatments using effective multilevel atom and conductance methods. With CosmoSpec one calculation of the CRR now only takes seconds on a standard laptop. This opens a way to more detailed forecasts and parameter explorations with the CRR, applications that we will address in a future paper.

The precision of our calculations is currently limited by the atomic model for neutral helium (see Sect. 4.6 for more detailed discussion) and uncertainties due to collisions. For neutral helium, uncertainties related to the photoionization cross section in particular are expected to exceed the level of 10%-30% for the Hei free-bound continuum. This should enter the total CRR at the level of a few percent. Similarly, transition rates, multiplet ratios and radiative/collisional singlet-triplet transitions among excited levels should be considered in more detail. For the total CRR, this implies an uncertainty at the level of a few percent. We also show that collisional mixing of -substates is noticeable at low frequencies (). At higher frequencies, corrections due to collisions remain at the sub-percent level (Fig. 8).

Our calculations clearly highlight the importance of helium feedback for the detailed shape of the CRR (Fig. 3 and 4), which introduces additional features into the derivative of the CRR with respect to the helium abundance (Fig. 6). This process is thus required for a detailed modeling of the dependence of the CRR on the helium abundance, which may be used to extract this parameter from future spectral distortion measurements. The CRR may also help break degeneracies between , and due to the differing sensitivity of the signal on these parameters. In combination with future CMB anisotropy and large-scale structure measurements (Abazajian et al., 2015), this could deliver additional constraints on neutrino physics.

Our calculations still use a few approximations that could affect the results at the level of . At this point, we did not separately include the effect of Heii Lyman-. This will change the dynamics of Heii recombination by Lyman- feedback (Chluba & Sunyaev, 2007; Switzer & Hirata, 2008a; Ali-Haïmoud et al., 2010) and also modify the interspecies feedback slightly. Similarly, for the Heii recombination history we neglected detailed radiative transfer effects, stimulated Heii 2s-1s two-photon emission and Raman processes, all effects that were studied in detail for the hydrogen recombination history (Chluba & Sunyaev, 2006b; Kholupenko & Ivanchik, 2006; Hirata & Switzer, 2008; Chluba & Sunyaev, 2008c; Hirata, 2008). Corrections due to charge transfer between Heii and Hei are expected to be small, however, in more detailed computations this processes should be considered. Also, one should treat the direct ionizations from the ground state and continuum escape more carefully.

The CRR probes the Universe at earlier stages than the CMB anisotropies, providing a new window for complementary constraints on non-standard physics. For instance, variations of fundamental constants, such as the fine-structure constant or electron-to-proton mass ratio, would shift the position and relative amplitudes of the recombination lines. We highlight that the effect of variations of and can be treated self-consistently using the effective rate and conductance method (see Appendix B). We find the effective recombination and photoionization rates to scale as and , respectively. These effects have been incorporated by HyRec for a while10, but they differ from previously adopted scalings (e.g., Kaplinghat et al., 1999; Scóccola et al., 2009). The differences may change the constraints derived from Planck data (Planck Collaboration et al., 2015b).

Another interesting effect is the extra ionization and excitation from annihilating or decaying dark matter particles, which could significantly enhance the CRR. It may therefore be possible to derive interesting upper limits on dark matter properties with future spectral distortion measurements, even in the absence of a detection. Tests of spatial variations in the composition of the Universe, for example, due to non-standard BBN models or cosmic bubble collisions (Dai et al., 2013; Silk & Chluba, 2014; Chary, 2015) may be another line of research. We look forward to exploring these directions in the future.

## Acknowledgments

JC is supported by the Royal Society as a Royal Society University Research Fellow at the University of Cambridge, UK. YAH is supported by NSF Grant No. 0244990, NASA NNX15AB18G, the John Templeton Foundation, and the Simons Foundation. Use of the GPC supercomputer at the SciNet HPC Consortium is acknowledged. SciNet is funded by: the Canada Foundation for Innovation under the auspices of Compute Canada; the Government of Ontario; Ontario Research Fund - Research Excellence; and the University of Toronto.

## Appendix A Improved J-resolved transitions rates from hydrogenic values

In the helium model of Rubiño-Martín et al. (2008), transitions among triplet states with which are not contained in the multiplet tables of Drake & Morton (2007) are filled with hydrogenic values using the approximation

 AHenLJ→n′L′J′ ≈(2J′+1)(2L′+1)(2S′+1)AHnL→n′L′. (9)

This approximation assumes that the partial rates for individual sub-levels within a multiplet are the same and that the total transition rate averaged over all and is given by the hydrogenic (non--resolved) values. This can be readily confirmed by computing the -averaged transition rate11

 AHenL→n′L′ =∑JJ′(2J+1)(2L+1)(2S+1)AHenLJ→n′L′J′≈AHnL→n′L′ (10)

after inserting Eq. (9). From the atomic physics point of view, this approximation is incorrect. It is well-known (e.g., Goldberg, 1935; Condon & Shortley, 1963) that the transitions within a given multiplet differ in strength and that their relative ratio deviates by much more than expressed with Eq. (9). Using the relation , between the line strength and transition rate together with the sum rule (e.g., see Edmonds, 1960)

 SnLJ→n′L′J′=(2J+1)(2J′+1)∣∣∣{LJSJ′L′1}∣∣∣2SnL→n′L′, (11)

where denotes the Wigner-6-symbol, one can show that the transition rates within a multiplet should scale like

 AHenLJ→n′L′J′≈(2L+1)(2J′+1)∣∣∣{LJSJ′L′1}∣∣∣2AHnL→n′L′. (12)

For , we therefore have the three cases

 rJ−1 =(L+J+2)(L+J+1)(L+J−1)(L+J−2)4J(2J+1)L(2L−1) (13a) rJ =(2J+1)(J−L−1)(L−J−2)(J+1)(L+J+1)(L+J−2)rJ−1 (13b) rJ+1 =J(J−L)(L−J−3)(2J+1)(L+J−1)(L+J+2)rJ (13c)

for the line ratios . Within these multiplets, the strongest line is the one for , which for roughly equals the hydrogenic value, . For large , the intensity of the second line () drops like , while for the third line we have . Clearly, this is a very different behavior than obtained with Eq. (9). For , we similarly find the line ratios

 rJ−1 =(J+1)(J−L−3)(L−J)(2J+1)(L+J)(L+J+3)rJ (14a) rJ =(2J+1)(L−J−1)(J−L−2)J(L+J+1)(L+J+4)rJ+1 (14b) rJ+1 =(L+J+4)(L+J+3)(L+J+1)(L+J)4(2L+3)(L+1)(2J+1)(J+1), (14c)

where now the strongest line in the multiplet is the one for .

When considering transitions from levels with (non--resolved) to levels with (-resolved), we average Eq. (12) over the initial level , giving us the approximation

 AHenL→n′L′J′≈(2J′+1)(2L′+1)(2S′+1)AHnL→n′L′, (15)

which also follows from Eq. (9). Notice that to obtain the hydrogenic approximation we also scale by the average transition frequency ratio to the third power. This usually is a small correction.

We find that the changes to the recombination history and spectrum caused by the above improvements are small (). The main reason is that only a few transitions among the high- states at are missing from the multiplet tables of Drake & Morton (2007). In addition, for transitions between levels with (non--resolved) to levels with (-resolved) the approximation of Rubiño-Martín et al. (2008) already represented the correct value, as shown in Eq. (15). This explains why the overall correction remains small.

### a.1 J-resolved photoionization cross sections

The total photoionization cross section from level is the sum of the photoionization cross sections to continuum states (label “c”, with an implicit integration over final electron momentum) with angular momentum quantum numbers :

 σnLJ=∑L′J′σnLJ→cL′J′, (16)

where the spin quantum number is implicit and held fixed throughout. Each angular-momentum-resolved transition is proportional to the line strength (Condon & Shortley, 1963):

 σnLJ→cL′J′∝SnLJ,cL′J′2J+1. (17)

Summing over final states, and using the sum rule (Edmonds, 1960)

 ∑J′SnLJ,cL′J′=2J+12L+1SnL,cL′, (18)

we arrive at

 σnLJ∝12L+1∑L′SnL,cL′, (19)

which is independent of .

## Appendix B Rescaling of the effective rates for hydrogenic ions

Noting that the dipole transition rates, (between levels and , where is the upper level) and photoionization cross sections12, , scale as

 Ai→j ∝σi→cν2dν∝αν3i→j/ca20Z−2=α5meZ4, (21)

we find that the radiative rate coefficients scale as (compare, Ali-Haïmoud & Hirata, 2011)

 Ri→j =Ai→j[1+nνij(Tγ)]=α5meZ4F(kTγhνij) (22a) βi→c =8πc2∫∞νicν2σi→c(ν)nν(Tγ)dν=α5meZ4B(kTγhνic) (22b) αc→i =8πc2feqi(Te)∫∞νicν2σi→c(ν)[1+nν(Tγ)]e−hν/kTedν =α2Zm2eA(kTγhνic,kTehνic) (22c) feqi(Te) =gi2gch3ehνic/kTe(2πmekTe)3/2. (22d)

Here, is the fine structure constant; the charge of the hydrogenic ion; its reduced mass; the Bohr radius; the bound-bound / bound-free transition frequency, respectively; the blackbody photon occupation number; and the statistical weights of the level and the continuum particle.

The inverse dipole transitions rate, , is simply given by . Because the probabilities of an electron reaching a specific level from another level within the excited levels just depends on ratios and sums of and , the effective rate coefficients , and (Eq. 17, 18 and 19 of Ali-Haïmoud & Hirata, 2011) for any hydrogenic ion can be obtained from those of hydrogen. Defining the factors , and , we find

 Ai(Tγ,Te) =fAAHi(fTTγ,fTTe) (23a) Bi(Tγ) =fBBHi(fTTγ) (23b) Ri→j(Tγ) =fBRHi→j(fTTγ), (23c)

where and are the reduced masses of hydrogen and the hydrogenic ion with charge in units of the electron rest mass. We also included the explicit dependence on , where is the valued used for the effective rate computation.

The effective conductances of hydrogenic helium can be similarly obtained with

 Gen′n(Tγ) =fBGe,Hn′n(fTTγ) (24a) Gin′n(Tγ) =fBGi,Hn′n(fTTγ) (24b) dGifb(Tγ)dν =fBdGi,Hfbdν(fTTγ) (24c)

using those of hydrogen (superscript ’H’) with the notation of Ali-Haïmoud (2013) for the conductances. In addition, we need to shift the energies of each frequency bin in the conductance table. These simple scalings are no longer exactly valid when collisions are included (Sect. 3.5). They also neglect smaller corrections due to the Lamb-shift, which depend on as well.

## Appendix C Classical derivation of ℓ-changing collisional transition rates

In this Appendix, we give a simple derivation of the classical result of Vrinceanu et al. (2012) in the Born approximation. We recover the leading term of their classical result up to corrections of order .

Let us consider an electron moving in a Keplerian orbit around a nucleus with charge , and subject to a perturbing force per unit mass acting on the electron. We have for singly-ionized helium and 1 for hydrogen and neutral helium, for which we approximate the helium nucleus and inner electron as a point charge, and the orbiting electron is the outermost, excited one.

The Gaussian perturbation equation (see e.g. Roy, 1982) for the osculating orbital elements are

 dcidt=∂ci∂\varv⋅f. (25)

We emphasize that these equations are exact. We define , where is the elementary charge and is the reduced mass of the electron-nucleus system. The energy per unit mass is , and evolves under the influence of the perturbing force according to

 dϵdt=\varv⋅f, (26)

and the angular momentum per unit mass evolves according to

 dLdt=r×f. (27)

We now make two assumptions:
We assume that the perturbing force is nearly constant accross the orbit of the electron. For the problem of interest ( generated by a passing ion), this is equivalent to the electric dipole approximation, where we assume that the dimensions of the orbit are much smaller than the impact parameter, .
We assume that the force varies on a timescale long compared to the orbital period. This is the orbital adiabatic approximation, and it requires that , where is the transition frequency from to the neighboring energy states.

With these assumptions, one can average the rate of change of orbital elements over one orbit, and obtain, first, that the energy is conserved, and second, that the secular rate of change of angular momentum is

 ⟨dLdt⟩=⟨r⟩×f, (28)

where is the time-average of the position vector over one orbit. Note that this reduces to the standard formula for a torque on an electric dipole moment when .

The time-averaged position vector for a Keplerian orbit can be shown to be , where is the semi-major axis and is the eccentricity vector, pointing from the orbital focus to percienter, and having magnitude , the orbital eccentricity. The secular rate of change of angular momentum is therefore

 ⟨dLdt⟩=32af×e. (29)

One can derive a similar expression for the rate of change of the eccentricity vector itself,

 ⟨dedt⟩=32κf×L, (30)

though we will not need it in the Born approximation.

We now specify to the perturbing force created by a passing ion of charge , moving on a trajectory :

 f=Zionq2eme^RR2. (31)

The specific angular momentum of the ion-nucleus system is conserved in the collision (up to corrections of order ). Defining as the polar angle of the passing ion, which starts at , we have , where is the impact parameter and the velocity of the ion at infinity. We may therefore rewrite the perturbing force as

 Missing or unrecognized delimiter for \left (32)

where form a fixed orthonormal basis in the plane of the ion’s orbit, with pointing in the incoming direction of the ion.

We now define the dimensionless angular momentum

 l≡L√κa, (33)

whose magnitude is . In terms of quantum numbers for a hydrogenic atom, . The secular equation can be rewritten as

 dldΦ=α[cosΦ