A halo model approach for the 21-cm power spectrum at cosmic dawn
AA halo model approach for the 21-cm power spectrum at cosmic dawn
Aurel Schneider ∗ and Sambit Giri † Center for Cosmology and Theoretical Astrophysics,Institute for Computational Science, University of Zurich, Switzerland.
Jordan Mirocha ‡ McGill University, Department of Physics & McGill Space Institute, Montr´eal, Canada. (Dated: November 26, 2020)Prior to the epoch of reionisation, the 21-cm signal of the cosmic dawn is dominated by theLyman- α coupling and gas temperature fluctuations caused by the first sources of radiation. Whileearly efforts to model this epoch relied on analytical techniques, the community quickly transitionedto more expensive semi-numerical models. Here, we re-assess the viability of simpler approaches thatallow for rapid explorations of the vast astrophysical parameter space. We propose a new analyticalmethod to calculate the 21-cm power spectrum based on the framework of the halo model. Both theLyman- α coupling and temperature fluctuations are described by overlapping radiation flux profilesthat include spectral red-shifting and source attenuation due to look-back (light-cone) effects. The21-cm halo model is compared to the semi-numerical code exhibiting generally goodagreement, i.e., the power spectra differ by less than a factor of three over a large range of k -modesand redshifts. We show that the remaining differences between the two methods are comparableto the expected variations from modelling uncertainties associated with the abundance, bias, andaccretion rates of haloes. While these current uncertainties must be reduced in the future, our worksuggests that inference at acceptable accuracy will become feasible with very efficient halo modelsof the cosmic dawn. I. INTRODUCTION
Current and upcoming radio interferometer telescopessuch as LOFAR [1, 2], MWA [3, 4], HERA [5], or SKA[6] are expected to detect for the first time the stronglyredshifted 21-cm clustering signal during and prior to theepoch of reionisation. These measurements will open upa new window on the universe [7, 8], providing insightsinto both astrophysics [9–11] and cosmology [12, 13].Next to the very first stars and galaxies [14, 15], the21-cm signal has the potential to find evidence for exoticsources of radiation [16, 17], new signs from the dark mat-ter sector [18–22] or other deviations from the standardΛCDM cosmological model [23–25].Predicting the 21-cm signal of the cosmic dawn, how-ever, remains a challenging task. On the one hand, thecomplicated physics of radiation-hydrodynamics com-bined with the enormous ranges of relevant scales makebrute-force simulations extremely difficult. On the otherhand, the poorly known characteristics of early star-forming sources as well as the complicated interplay be-tween gas cooling and feedback add important uncer-tainties that need to be either understood or properlyparametrised.First attempts to predict the 21-cm clustering be-fore reionisation were based on analytical techniques, us-ing a combination of cosmological perturbation theoryand excursion-set modelling prescriptions for the sources ∗ [email protected] † [email protected] ‡ [email protected] [14, 26, 27] (for a more recent attempt see also Ref. [28]).While these calculations were able to predict many im-portant features of the 21-cm power spectrum – such asfor example the characteristic double-peaked shape dueto the Lyman- α coupling and heating epochs – it remainsunclear how well they agree with more detailed calcula-tions [30].A further important step towards more realistic pre-dictions of the 21-cm signal at cosmic dawn was takenwith the development of semi-numerical methods such as [31, 32] or simfast21 [33] (see also Refs. [34–36] for other semi-numerical methods). These modelsnumerically evolve the matter perturbations, Lyman- α coupling, and temperature fluctuations on a grid, wherethe source distributions are either obtained via a halo cat-alogue from N-body simulations or via an excursion-setrecipe. Semi-numerical techniques consist of a major im-provement with respect to analytical approaches, mainlybecause they follow the evolution of the spin temperaturein configuration space and are therefore able to producemaps of the 21-cm signal.In principle, more accurate predictions can be obtainedby post-processing numerical N -body simulations usingradiative-transfer (RT) calculations. This is the strategyfollowed e.g. by Grizzly [37, 38], C -Ray [39, 40], and CRASH [41]. While the former code is based on approx-imate but faster one-dimensional RT calculations, thelatter two follow a full three dimensional ray-tracing ap-proach. First results for the epoch of pre-reionisation canbe found in Refs. [38, 42] as well as Ref. [43] and [44].Full radiation-hydrodynamical simulations, includingself-consistent formation of sources, are becoming an in-creasingly common tool to study the process of reionisa-tion [e.g. Refs. 45–49]. However, for the pre-reionisation a r X i v : . [ a s t r o - ph . C O ] N ov epoch of cosmic dawn, such simulations remain rare[50]. This is because they require narrowly binned multi-frequency RT calculations in a simulation box that re-solves the mini-haloes hosting the first sources and thataccounts for the large distances travelled by Lyman- α and X-ray radiation.In this paper we propose an analytical approach topredict the 21-cm global signal and power spectrum atcosmic dawn. While being potentially less accurate thanthe different numerical techniques mentioned above, sucha model has the advantage of providing fast predictionsbased on a well defined framework. Different parametri-sations of the source modelling or effects from vary-ing cosmology can be readily implemented and tested.An analytical method is also particularly well suitedfor parameter inference, which is an important aspectof 21-cm cosmology due to the large uncertainties re-lated to early-universe galaxy formation that have to beparametrised and marginalised over. Comparing multi-parameter models with observations requires a very largenumber of calculations that can only be performed withfast prediction routines.Our method is inspired by earlier work of Holzbauerand Furlanetto [51] who used the halo model to predictthe clustering from the Lyman- α coupling between thegas and the first star-light at cosmic dawn. We extendthis approach adding a description for the temperaturefluctuations as well as an improved modelling of the haloaccretion and the star formation rates. With this athand, we obtain a complete prediction of the 21-cm powerspectrum at cosmic dawn.The paper is structured as follows. In Sec. II we discussthe source modelling including halo mass function, bias,gas accretion, and star-formation efficiency. Sec. III sum-marises the derivation of the differential brightness tem-perature. In Sec. IV we present the halo model of fluxprofiles and we show results for the 21-cm power spec-trum assuming different star-formation efficiencies. Wethen go on and compare our model to other analyticaland semi-numerical methods in Sec. V before concludingin Sec. VI. II. SOURCE MODELLING
Accurately quantifying the abundance, distribution,and emission of sources is a crucial step for any methodpredicting the 21-cm signal. The source distribution, forexample, can be readily quantified using prescriptions forthe halo mass function and the halo bias. Regardingthe source emission, analytical and semi-numerical mod-els usually rely on estimates for the halo accretion ratescombined with a parametrisation of the star-formationefficiency. While the former describes the accretion ofprimordial gas onto the parent halo, the latter genericallyparametrises all star formation and feedback processes.In the present section, we will first discuss the halomass function and biasing model, before we describe and compare different estimates for halo mass accretion rates.At the end we present our parametrisation for the starformation efficiency parameter.
A. Halo mass function and bias
Both the halo mass function and the halo bias arecentral ingredients of the halo model which our methodis based on. For the halo mass function, we assumethe standard form motivated by the extended Press-Schechter (PS) approach [52, 53], i.e., dnd ln M = 12 ¯ ρ m M f ( ν ) d ln νd ln M (1)where the peak-height of perturbations is defined as ν = δ c ( z ) /σ with δ c ( z ) = 1 . /D ( z ), D ( z ) being thecosmological growth factor. The first crossing distribu-tion us given by f ( ν ) = A (cid:114) qνπ (cid:2) qν ) − p (cid:3) exp( − qν/ , (2)which comes with two free model parameters q and p [see Ref. 54]. The factor A is fixed via the relation A =(1 + 2 − p Γ(1 / − p ) / √ π ) − .Using the peak-background split model introduced inRef. [55], the halo bias can be described within the sameextended PS formalism, leading to the relation [see e.gRef. 56] b ( M ) = 1 + qν − δ c ( z ) + 2 pδ c ( z )[1 + ( qν ) p ] . (3)In the following, we use the values q = 0 .
85 and p = 0 . N body simulations from Refs. [20, 57] than theoriginal Press-Schechter [52] or Sheth-Tormen [54] massfunctions (which are characterised by q = 1, p = 0 and q = 0 . p = 0 .
3, respectively).
B. Halo mass accretion
There exists different ways to analytically quantify thegrowth of haloes. Here we discuss three models that areused in the literature, and we compare them to resultsfrom N -body simulations.
1. Exponential (EXP) halo growth
The simplest way to obtain a halo accretion history isto assume that haloes grow exponentially. This assump-tion has been shown to provide a good match to simu-lations, especially at high redshifts [see e.g. Refs 59–61].We consider the following equation M ac ( M, z ) = M exp [ α ( z − z )] (4) M a c [ M / h ] EXP accretionAM accretionEPS accretionSimulations 5 10 15 20 25 30z10 d M a c / d t [( M / h ) / y r ] M = 10 M /h M = 10 M /h M = 10 M /h M = 10 M /h 5 10 15 20 25 30z10 S F R D [ M / M p c h / y r ] EXP accretionAM accretionEPS accretion
FIG. 1. Halo growth, mass accretion rate, and star-formation rate density (SFRD) for the three different models discussed inSec. II B.
Left:
Halo mass as a function of redshift for four different final masses M (green, brown, red, and blue lines). Theorange lines correspond to the mean halo growth from the B20 [58] simulations. Centre:
Mass accretion rates for the samemodels and halo masses.
Right:
Resulting star-formation rate densities (SFRD) for the three different mass-accretion ratemodels (see Eq. 13) where M is the halo mass at the final redshift z . Thefactor α can be determined using simulations and hasbeen shown to be only weakly mass dependent at highredshift. Following Ref. [60] we therefore assume α =0 .
79 independent of the halo mass M .In the left-hand and centre panels of Fig. 1 we showthe redshift evolution of the halo mass ( M ac ) and theaccretion rate ( ˙ M ac ). The dash-dotted lines correspondto the exponential (EXP) growth model of Eq. (4). Themodel is in very good agreement with the numerical N -body simulations from Behroozi et al. [58, B20] (orangelines). Note, however, that the simulations do not coverall redshifts and mass ranges of interest to the presentstudy.
2. Extended Press-Schechter (EPS) method
Halo accretion rates can also be obtained by meansof the extended Press-Schechter (EPS) formalism. Forexample, Neistein and van den Bosch [62] calculated anensemble of EPS merger trees and determined the aver-age growth of their main branches. Based on this, theyproposed the equation dM ac dz = − (cid:114) π M ac (cid:112) S ( M ac /Q ) − S ( M ac ) dδ c ( z ) dz , (5)where S ≡ σ at redshift zero [see also Ref. 63]. This dif-ferential equation can be solved assuming M ac ( M, z ) ≡ M , where z designates the final redshift of interest.The value Q has to be selected empirically, and we use Q = 1 . Q = 2 − . Q -parameter according tothese simulations.The advantage of the EPS model compared to themuch simpler exponential model lies in the fact thatEq. (5) is, in principle, sensitive to changes in cosmology.Whether the true cosmology dependence of the mass ac-cretion rate can be accurately reproduced by the EPSmodel, remains, however, to be tested.
3. Abundance matching (AM) method
The final method we are investigating here is inspiredby the abundance matching (AM) technique [64] and wasfirst applied in Ref. [65] as a measure of high-redshiftgalaxy growth. The method aims to connect haloes be-tween different redshift bins ( z n ) by matching (cid:90) ∞ M n dM n ( M, z n ) = (cid:90) ∞ M n − dM n ( M, z n − ) , (6)where n ( M, z ) is the cumulative halo mass function de-fined by n ( M, z ) = (cid:90) ∞ M dnd ln M (cid:48) (cid:18) dM (cid:48) M (cid:48) (cid:19) . (7)Connecting all masses M n at different redshifts z n allowsus to estimate the accretion rate of haloes. Note thatthe AM method connects haloes in a strictly hierarchicalway, which means that the i -largest halo at the final red-shift z will be assumed to stay the i largest halo at allhigher redshifts. As a consequence, the model implicitlyassumes haloes to exclusively grow via smooth accretion,since halo mergers cannot be accounted for.In Fig. 1 the halo growth and accretion rates of theAM model are shown as dashed lines. Compared tothe other methods, the AM model predicts significantlyslower halo growth over all redshifts and mass ranges.The AM results do not match the B20 simulations verywell. Note that our findings are in qualitative agreementwith Ref. [66], where the AM model was compared to afitting function of Ref. [61]. C. Star-formation efficiency and stellar-to-halomass ratio
The halo growth rate discussed above provides an esti-mate of the amount of total matter accreted onto a halo.However, we still need to parametrise the star-formationefficiency, i.e., how much of that accreted matter willbe transformed into stars that emit radiation. FollowingRefs. [67–69], we define f ∗ ( M ) ≡ ˙ M ∗ / ˙ M ac , (8)where ˙ M ∗ and ˙ M ac are the stellar and halo accretionrates. Note that Eq. (8) is different from the stellar-to-halo mass ratio ˜ f ∗ ( M ) ≡ M ∗ /M ac (9)often used in the literature [see e.g. Ref. 11, 70]. Thestar-formation efficiency and the stellar-to-halo mass ra-tio can be connected by the integral˜ f ∗ ( M ) = 1 M ac (cid:90) f ∗ ˙ M ac dt. (10)Note that for the simple case of a redshift-independentstar-formation efficiency, Eq. (10) leads to ˜ f ∗ = f ∗ . Amore general study of the relation between star-formationefficiency and stellar-to-halo mass ratio, including ex-plicit redshift dependences motivated by feedback pro-cesses, can be found in Refs. [65, 71].In this paper, we use the functional form f ∗ ( M ) = 2(Ω b / Ω m ) f ∗ , ( M/M p ) γ + ( M/M p ) γ × S ( M ) (11)as a parametrisation for the star-formation efficiency.Eq. (11) consists of a double-power law, multiplied witha small-scale function S ( M ) = [1 + ( M t /M ) γ ] γ , (12)that may provide either a suppression or a boost at thetruncation mass scale M t . A suppression at small scalescould naturally occur at scales where atomic cooling pro-cesses become inefficient. A small-scale boost, on the other hand, could emerge due to the presence of popula-tion III stars in mini-haloes.Note that Eq. (11) does not depend on redshift, and wecan therefore set f ∗ = ˜ f ∗ in this paper. We will nonethe-less formally distinguish between f ∗ and ˜ f ∗ in order toavoid confusion and to acknowledge the fact that a morerealistic description of the star-formation efficiency maywell include an explicit redshift dependence. D. Star-formation rate density and collapsefraction
With the halo mass function, accretion rate, and starformation efficiency at hand, it is possible to calculatethe mean star-formation-rate density (SFRD), which isan important ingredient for calculating the global 21-cmsignal (see Sec. III). We define the SFRD as the inte-gral over the halo mass function weighted by the star-formation efficiency and the accretion rate, i.e.,˙ ρ ∗ ( z ) = (cid:90) dM dndM f ∗ ( M ) ˙ M ac ( M, z ) . (13)The definition above differs from the relation ˙˜ ρ ∗ ( z ) =¯ ρ m df coll /dt , which is based on the global collapse fraction( f coll ) and is often used in the literature [e.g. Refs. 14,72]. The collapse fraction is given by the integral f coll ( z ) = 1¯ ρ m (cid:90) dM ˜ f ∗ ( M ) dndM M, (14)providing the total ratio of stars to matter in the uni-verse. Note that Eq (14) includes the stellar-to-halo massratio inside of the integral and is therefore slightly dif-ferent from the standard definition in the literature [seee.g. Refs. 73, 74].The star-formation rate densities based on the threedifferent mass accretion models introduced above areshown as black lines in the right-hand panel of Fig. 1.While they have similar general trends, the exponentialgrowth (EXP) model is about a factor of two larger thanthe method based on abundance matching (AM). Theextended Press-Schechter (EPS) model lies in between,being closer to the AM method at very high and closerto the EXP model at lower redshifts.The differences between the star-formation rate densi-ties shown in Fig. 1 directly affect the modelling of theLyman- α coupling, gas heating, and ionisation. This isone of the main reasons why comparisons between dif-ferent methods to calculate the 21-cm signal are hard tointerpret when they are using different implementationsfor the star-formation rate (as we will see in Sec. V). III. 21-CM BRIGHTNESS TEMPERATURE
The 21-cm differential brightness temperature ( T ) isa function of the background radiation ( T γ ), the spintemperature of the gas ( T s ), the neutral hydrogen frac-tion ( x HI ) and the gas density field ( δ b ) which all dependon redshift z and position x . Following e.g. Ref. [75], thebrightness temperature can be written as T ( x , z ) = 27 x HI (1 + δ b ) × (cid:18) Ω b h . (cid:19) (cid:18) . m h (1 + z )10 (cid:19) (cid:18) − T γ T s (cid:19) (15)in milli-Kelvin [mK]. Assuming a standard ΛCDM modelwithout exotic radio sources, the background tempera-ture is dominated by the cosmic microwave background(CMB) radiation. The rightmost expression in Eq. (15)can be written as (cid:18) − T γ T s (cid:19) (cid:39) x tot x tot (cid:18) − T γ T k (cid:19) , (16)where x tot ≡ x α + x c ( x α and x c denoting the radiativeand collisional coupling coefficients). Throughout thispaper, we set the collisional coupling to zero, as it is onlyimportant at very high redshifts beyond z ∼ T k ) is obtained via the differen-tial equation32 dT k ( x , z ) dz = T k ( x , z ) ρ ( x , z ) dρ ( x , z ) dz − Γ h ( x , z ) k B (1 + z ) H , (17)where ρ is the matter density and Γ h the heating sourceterm. The latter is given by the sumΓ h ( x , z ) = 4 π (cid:88) i f i f X,h × (cid:90) ∞ ν th dν (cid:0) ν − ν thi (cid:1) h P σ i ( ν ) J X,ν ( x , z ) (18)with i = { H , He } , i.e. the hydrogen and helium com-ponents with fractions f i = n i /n b , threshold energies ν th i h P = { } eV, and cross sections σ i ( ν ). For thefraction of X-ray energy deposited as heat, we assume f X,h = ¯ x . e , where ¯ x e is the free electron fraction [seeFig. 4 in Ref. 76] that is calculated according to Eqs. 12and 13 in Ref. [74].The radiation coupling coefficient, induced by theWouthuysen-Field effect [77, 78], can be written as x α ( x , z ) = 1 . × (1 + z ) S α J α ( x , z ) , (19)where we set S α = 1 for simplicity [following e.g. Refs.79, 80].The gas temperature and Lyman- α coupling depend onthe flux terms J X and J α . While the former is dominatedby X-ray radiation between a few hundred to a few thou-sands eV, the latter stems from the narrow spectral rangebetween the Lyman- α and the Lyman-limit frequencies.We parametrise the spectral energy distributions of theLyman- α and X-ray flux as simple power laws I s ( ν ) = A s ν − α s , (20) with s = { α, X } . The normalisation A s is defined so thatintegrating I s ( ν ) over the corresponding energy range be-comes unity [see Ref. 81]. The number emissivity of UVphotons between the Lyman- α and Lyman-limit range isgiven by ε α ( ν ) = N α m p I α ( ν ) , (21)where N α is the number of photons per baryon emittedin the range between ν α and ν LL . The energy emissivityof X-ray photons is ε X ( ν ) = f X c X I X ( ν ) νh P (22)where f X is a free parameter of order unity and c X isa normalisation factor constrained by observations. Itis set to c X = 3 . × erg yr s − M − (cid:12) based on thefindings of Ref. [82]. A. Global Signal
The global differential brightness temperature is di-rectly obtained by averaging Eqs. (15) and (16). Wethereby set δ b to zero and assume spatially averaged val-ues for the gas temperature and coupling coefficient. TheLyman- α coupling is obtained via Eq. (19), where themean Lyman- α flux is given by¯ J α ( z ) = (1 + z ) π n m (cid:88) n =2 f n (cid:90) z ( n )max z dz (cid:48) cε α ( ν (cid:48) ) H ( z (cid:48) ) ˙ ρ ∗ ( z (cid:48) ) , (23)with ν (cid:48) = ν (1+ z (cid:48) ) / (1+ z ). The recycling fractions f n aretaken from Ref. [72] with the sum truncated at n m = 23.The integration limit is given by z ( n )max = (1 + z )[1 − ( n + 1) − ] / (1 − n − ) − , (24)designating the maximum redshift from which photonscan Doppler-shift into the Lyman resonances.The global X-ray number flux per frequency is givenby the relation [70]¯ J X,ν ( z ) = (1 + z ) π (cid:90) ∞ z dz (cid:48) cε X ( ν (cid:48) ) H ( z (cid:48) ) e − τ ν ˙ ρ ∗ ( z (cid:48) ) . (25)Here we have introduced the optical depth parameter τ defined as τ ν ( z, z (cid:48) ) = (cid:90) z (cid:48) z dz (cid:48)(cid:48) dldz (cid:48)(cid:48) (cid:88) i n i σ i ( ν (cid:48)(cid:48) ) , (26)where i = { HI , HeI } and where ν (cid:48)(cid:48) is the frequency red-shifted from the source at z (cid:48) to z (cid:48)(cid:48) . Plugging Eq. (25)into Eq. (18) leads to the mean heating rate ¯Γ h ( z ). Theglobal temperature evolution is finally obtained by solv-ing 32 d ¯ T k dz = 3 ¯ T k (1 + z ) − ¯Γ h ( z ) k B (1 + z ) H , (27)which corresponds to Eq. (17) at order zero in densityperturbations.The solution of Eq. (27) together with ¯ x α from Eq. (23)allows us to obtain the global differential brightness tem-perature ( ¯ T ). Examples of the ¯ T signal are shown inthe bottom-left panel of Fig. 4. B. Power spectrum
The 21-cm brightness temperature defined in Eq. (15)is a function of both redshift and position. At linearorder, the 21-cm perturbations are given by δ ( x , z ) = β b δ b + β α δ α + β h δ h + β p δ p − δ dv , (28)This expansion is identical to the one proposed inRef. [73], except that, for reasons that will become ev-ident later on (see Sec. IV C), we furthermore separatethe temperature fluctuations into a heating ( δ h ) and aprimordial ( δ p ) term: δ T = f T δ h + (1 − f T ) δ p , (29)where f T ≡ ( ¯ T k − ¯ T p ) / ¯ T k . In this context primor-dial means prior to any heating from sources, i.e., theregime where the gas cools adiabatically and the tem-perature fluctuations are seeded by the matter perturba-tions. Note furthermore that Eq. (29) is a direct conse-quence of the assumption T k = T h + T p .The other terms of Eq. (28) designate the baryonperturbations ( δ b ), the Lyman- α coupling perturbations( δ α ), and the perturbations due to the line-of-sight ve-locity gradient ( δ dv ) caused by redshift-space distortioneffects. At linear order and in Fourier space, δ dv is simplygiven by δ dv = µδ m , where µ is the cosine of the anglebetween the wave vector k and the line of sight [83].The pre-factors of the individual perturbations inEq. (28) are given by β b (cid:39) , (30) β α = ¯ x α ¯ x tot (1 + ¯ x tot ) , (31) β h (cid:39) f T ¯ T γ ( ¯ T k − ¯ T γ ) , (32) β p (cid:39) (1 − f T ) ¯ T γ ( ¯ T k − ¯ T γ ) , (33)and only depend on redshift z but not on the positionvector x .Based on Eqs. (28) and (29), it is straight-forward tocalculate the power spectrum and to sort all componentswith respect to their power of µ . Taking the average overthe angle then leads to [84] P = P αα + P hh + P pp + P bb + 2 (cid:16) P αh + P αp + P αb + P hp + P hb + P pb (cid:17) + 23 (cid:16) P αm + P hm + P pm + P bm (cid:17) + 15 P mm . (34) In the next section, we propose a method to calculate allindividual auto and cross power spectra of Eq. (34) inorder to obtain a fast, analytical estimate for the powerspectrum of the 21-cm brightness fluctuations. Note fur-thermore, that when talking about the 21-cm power spec-trum, we will mean either P or the expression ¯ T ∆ depending on the context. The dimnesionless quantity∆ is defined as ∆ ≡ k P / (2 π ). IV. HALO MODEL
In the framework of the original halo model haloes areconsidered as the building blocks of the universe [see e.g.Ref. 56]. The halo abundance, distribution, and inter-nal profiles are used to calculate the matter power spec-trum at both linear and nonlinear scales. In the contextof 21-cm clustering, however, the focus is not so muchon haloes per se, but rather on the sources inhabitinghaloes, which emit radiation, thereby affecting the gascells around them. Instead of halo profiles, we thereforeconsider radiation profiles around sources that extendinto the intergalactic space far beyond the halo limits.This means that there will typically be several overlap-ping radiation profiles from different sources affecting anysingle gas volume.Although the picture of overlapping radiation profilesis rather different from the original halo model, it turnsout that the 21-cm power spectrum can be described ina very similar way. In this section, we first introducethe formalism, before going into the details of the sourceprofiles and the description of temperature fluctuations.At the end we show the resulting 21-cm power spectrumassuming three benchmark models with different astro-physical parameters.
A. Power spectrum description
In the context of the 21-cm halo model, the powerspectra of different components can be calculated in thefollowing way: P XY ( k, z ) = β X β Y (¯ ρ b f coll ) (cid:90) dM dndM ˜ f ∗ M | u X || u Y | ,P XY ( k, z ) = β X (¯ ρ b f coll ) (cid:90) dM dndM ˜ f ∗ M | u X | b X (35) × β Y (¯ ρ b f coll ) (cid:90) dM dndM f ∗ M | u Y | b Y × P lin ,P XY ( k, z ) = P XY ( k, z ) + P XY ( k, z ) , where P lin ( k ) the linear matter power spectrum, u X,Y ( k, M, z ) the Fourier transformed flux profile, and b X,Y ( M, z ) the halo bias. The subscripts X and Y referto the mass ( m ), baryon ( b ), Lyman- α ( α ), and heat-ing ( h ) components. Note furthermore that Eq. (35) de-scribes both auto and cross spectra depending on whether FIG. 2. Normalised Lyman- α coupling (left), X-ray energy deposition (centre), and heating (right) profiles. In the top-row wevary the redshift at fixed halo mass and in the bottom-row we vary the halo mass at fixed redshift. X = Y or X (cid:54) = Y . This means that the halo model pro-vides all components of Eq. (34) except the ones thatinclude the primordial gas temperature ( p ). We will de-rive these in Sec. IV C.In the halo model framework, the one-halo term ( P XY )describes the clustering within one single source profile,while the two-halo term ( P XY ) accounts for the signalinduced by different sources. As a consequence P XY and P XY dominate at small and large scales, respectively,with a transition region corresponding to the typical sizeof the source profiles. Not surprisingly, only the two-haloterm carries information about the spatial distribution ofsources via the components b X,Y and P lin . The one-haloterm, on the other hand, carries information about theshot-noise of sources, and its shape is only controlled bythe Fourier transformed radiation profiles. B. Flux profiles
The key components of the halo model are the radi-ation profiles around sources. Following Holzbauer andFurlanetto [51], the profile of the Ly- α radiation can be written as ρ α ( r | M, z ) = 14 πr n m (cid:88) n =2 f n ε α ( ν (cid:48) ) f ∗ ˙ M ac ( z (cid:48) | M, z ) (36)with ν (cid:48) = ν (1 + z (cid:48) ) / (1 + z ). The profile is proportional tothe halo mass accretion ˙ M ac and exhibits the character-istic r decrease with radius r . It furthermore dependson the look-back redshift z (cid:48) = z (cid:48) ( r ) that corresponds tothe redshift when a photon at radius r has been emittedat the source. At the emission redshift z (cid:48) , the source isin an earlier stage of evolution (compared to the redshiftof the signal z ), which means that its accretion rate issmaller as well. The look-back redshift is obtained byinverting the co-moving distance r ( z (cid:48) | z ) = (cid:90) z (cid:48) z cH ( z (cid:48)(cid:48) ) dz (cid:48)(cid:48) (37)which has to be done numerically.The shape of the Lyman- α flux profile is plotted in theleft-hand panels of Fig. 2, where the top and bottom pan-els specifically highlight the dependencies on redshift ( z )and halo mass ( M ). Here we have used the EPS mass ac-cretion model, but the plots look very similar if another FIG. 3. Normalised Fourier transforms of the Lyman- α (solid) and the heating profiles (dashed). In the left- and right-handpanels we vary redshift and halo mass while keeping the other one constant. model is applied instead. All Lyman- α profiles are char-acterised by a 1 /r decrease close to the source, whichbecomes gradually steeper towards the outer parts. Thesteepening is a direct result of the radiation originatingfrom the source at a higher redshift ( z (cid:48) ) when the ac-cretion rate onto the source was smaller. At around 200Mpc/h, the profiles exhibit a steep drop that is due tothe photons having redshifted out of the Lyman- α series.At very high redshift beyond z=25, the drop-off shifts to-wards much smaller radii below 100 Mpc/h, because thesource is so young that the radiation did not have timeto expand further.A closer inspection of the Lyman- α flux profiles inFig. 2 reveals small discontinuities in the form of step-like features between r = 0 . −
100 Mpc/h. These smallsteps are a consequence of the sum in Eq. 36 and havebeen predicted in earlier work [see Ref. 51].The heating of the gas depends on the flux of X-rayradiation and can be defined similarly to the Lyman- α profile, i.e., ρ xray ( r | M, z ) = 1 r (cid:88) i f i f X,h × (cid:90) ∞ ν th dν (cid:0) ν − ν th i (cid:1) h P σ i ( ν ) ε X ( ν (cid:48) ) e − τ ν f ∗ ˙ M ac ( z (cid:48) | M, z ) , (38)where i = { H , He } . Compared to Eq. (36) there is an ad-ditional attenuation term due the optical depth definedin Eq. (26).The shape of the energy deposition profile ( ρ xray ) isshown in the middle panels of Fig. 2, the top and bottompanels again highlighting its redshift and mass depen-dencies. As a general trend, the profile extends furtherin radius (to around 500 Mpc/h) and is more graduallysuppressed than the Lyman- α coupling profile. This isdue to the fact that hard X-ray radiation can travel largedistances until it is deposited as heat, and X-ray photons do not redshift out of a well defined spectral range. Sincethe large distances traveled by X-ray photons also leadto increased look-back redshifts (the difference between z and z (cid:48) ), the energy deposition profile shows a strongerattenuation towards large radii.The energy deposition from X-ray emission leads to astrong increase of the temperature fluctuations aroundsources. The corresponding profile can be obtained viathe differential equation32 dρ h ( r | M, z ) dz = 3 ρ h ( r | M, z )(1 + z ) − ρ xray ( r | M, z ) k B (1 + z ) H . (39)The first term on the right-hand-side of this equation de-scribes the cooling due to the expansion of space, thesecond term corresponds to the energy deposition fromX-ray radiation. The solution of Eq. (39) is called theheating profile ( ρ h ) in agreement with the notation ofSec. III B, where we separated the temperature fluctua-tions into a primordial and a heating term (see Eq. 29).The heating profile is illustrated on the right-hand-side of Fig. 2, where the top and bottom panels againshow the dependency with redshift and halo mass. Theshape of the heating profile ( ρ h ) is similar to the energydeposition profile ( ρ xray ), which is not a surprise, since X-ray radiation is assumed to be the only source of heatingin Eq. (39).The source profiles of Eqs. (36) and (39) are key ingre-dients of the halo model. However, in order to substitutethem into Eq. (35), we need to first calculate their nor-malised Fourier transforms. They are given by u i ( k | M, z ) = (cid:82) drr ρ i ( r | M, z ) j ( kr ) (cid:82) drr ρ i ( r | M, z ) (40)where j ( x ) = sin( x ) /x is the spherical Bessel function oforder 0. The subscript i either stands for α or h , denotingthe Lyman- α emission and heating profiles.In Fig. 3 we show the Fourier transformed profiles withtheir redshift and halo mass dependencies. All profiles goto unity at low k -modes, a key characteristic that is guar-anteed by the normalisation of Eq. (40). Towards highervalues of k , the profiles become strongly suppressed. Ingeneral, the Lyman- α profiles are more suppressed thanthe heating profiles, which is a consequence of the factthat soft X-ray photons do not travel far before being ab-sorbed by the gas. This early absorption leads to an ex-cess of small-scale clustering compared to the more freelyemitted Lyman- α flux. One exception to this behaviourbecomes visible at z = 30 (see red lines in the left-handpanel), where both profiles have a very similar shape.This is not surprising because at very high redshifts, theshape of the profiles is driven by the emission of a veryyoung source with photons that did not have time totravel far.Fig. 3 shows a clear redshift dependence of the nor-malised Fourier profiles with more small scale clusteringtowards higher redshifts. The mass dependence, on theother hand, is very weak. This means that, in principle,the halo model of Eq. (35) could be simplified consider-ably by assuming the profiles u i not to depend on halomass. Note, however, that for the sake of completenesswe keep the full mass dependence in our model. C. Temperature fluctuations
In Sec. III B, we have separated the temperature fluc-tuations into a heating ( δ h ) and a primordial ( δ p ) term.While the heating term is sourced by the X-ray flux emis-sion and can therefore be readily described by the halomodel, the primordial fluctuations are driven by the mat-ter fluctuations and can be directly solved via Eq. (17)at the linear level. Since we have separated out the con-tribution from the sources, we can set Γ = 0. After lin-earising T = ¯ T (1 + δ p ) as well as ρ = ¯ ρ (1 + δ ) we obtainthe solution δ p = (1 + δ ) / − . (41)for the primordial (prior to sources) temperature fluctu-ations. Note that in Fourier space, the matter perturba-tions are readily obtained by setting δ = √ P mm . Withthis at hand, we can now derive the reminding auto andcross power spectra of Eq. (34) which include the primor-dial heating term and are thus not covered by Eq. (35).They are given by P pp ( k, z ) = β p ( z ) δ p ( k, z ) , (42) P Xp ( k, z ) = β p ( z ) (cid:112) P XX ( k, z ) δ p ( k, z ) , (43)where X again stands for { α, h, b, m } . D. Results
Based on the formalism derived above, it is now possi-ble to calculate the 21-cm power spectrum for a given setof model parameters. In this section, we will first present our choices of parametrisation before showing the powerspectrum for a selected set of redshift and k -modes.For the Lyman- α sources, we assume N α = 10000evenly distributed ( α α = 0) over the energy range be-tween the Lyman- α and Lyman-limit frequencies (seeEq. 21). The X-ray energy emission is defined by f X = 1plus a power-law spectral energy distribution with α X =1 . E = 0 . − q = 0 .
85. The halo bias is modelled based thepeak-background split approach (Eq. 3) with the sameparameters than what is used for the halo mass function.In order to highlight the sensitivity of the results tothe source parametrisation, we focus on three benchmarkmodels (A, B, and C) that are characterised by differentstar formation efficiencies ( f ∗ ). All models have the samelarge-scale behaviour given by the double power law ofEq. (11) with f ∗ , = 0 . γ = 0 . γ = − .
61, and M p = 2 × M (cid:12) /h [see Ref. 69]. At small mass scales,however, Model A assumes a strong additional suppres-sion (with M t = 5 × M (cid:12) /h, γ = 1, γ = −
4) mim-icking the effects of inefficient cooling processes. ModelC, on the other hand, is characterised by a boost of f ∗ to-wards very small masses (with M t = 10 M (cid:12) /h, γ = 1, γ = 1). Such a behaviour can be motivated by thepresence of Population-III stars in mini-haloes. ModelB finally shows neither additional suppression nor boostof f ∗ but a continuation of the power-law decrease downto the smallest masses (i.e. γ = 0). All models aretruncated at M min = 5 × M (cid:12) /h.The star-formation efficiencies of the three benchmarkmodels are plotted in the top-left panel of Fig. 4. Thetypical shape from the double-power law prescription isvisible at large masses above ∼ M (cid:12) /h. At smallerhalo masses the models diverge showing the characteristicsuppression, power-law continuation, and boost of thebenchmark models A, B, and C described above.The effect of the different star-formation efficiencieson the global differential brightness signal is illustratedin the top-right panel of Fig. 4. While model A leadsto a narrow absorption signal at a rather low redshiftof z ∼
14, model B and C show wider troughs shiftedtowards z ∼
16 and z ∼
18, respectively[85].The 21-cm power spectrum as a function of k -valuesare plotted in the three middle panels of Fig. 4 with in-creasing redshift from left to right. Many of the lineshave a wave like feature with a local flattening or mini-mum at k ∼ . − k -modes. At k = 0.01 h/Mpc (left) and 0.1 h/Mpc (mid-dle) the power spectra show the characterised double-peak feature, which indicate the two characteristic epochswhere the Lyman- α coupling and the X-ray heating dom-inate. At k = 1 h/Mpc the two peaks have merged into0 FIG. 4.
Leftmost panels:
Star-formation efficiency (top) and resulting global signal (bottom) for the three benchmark modelswith suppressed, unchanged and boosted small-scale behaviour (blue, cyan, and magenta).
Remaining panels:
Power spectraas a function of k -modes (top) and redshift (bottom) for the same models. one single broader peak in qualitative agreement withother work from the literature [see e.g. Refs. 30, 43].The results presented in Fig. 4 are based on the EPSmass accretion model plus the halo mass function andbias parameters described above. However, it is impor-tant to notice that these modelling choices introduce sig-nificant uncertainties regarding the 21-cm power spec-trum. In the Appendices A and B, we quantify the ef-fects due to the choice of the mass accretion model, thehalo mass function and halo bias. The bottom line ofthis analysis is that modelling choices, such as switchingfrom Sheth-Tormen to Press-Schechter halo prescription,or using the AM instead of the EPS accretion rate mod-elling, has an effect on the 21-cm power spectrum thatcan in some cases be as large as a factor of ∼
10. Thedifference can be significantly reduced if the flux param-eters N α are re-adapted so that the global signals areforced to match, but a remaining difference of the powerspectrum of about a factor of ∼ V. COMPARISON WITH OTHERAPPROACHES
In this section we compare the 21 cm halo modelwith the analytical model from Refs. [14, 27] (abbrevi-ated to BLPF model) and with the semi-numerical code [31, 32]. While the former is more of a gen-eral consistency check, the latter consist of a true testfor our model. Note, however, that the comparison with is not straight-forward, mainly because thereare subtle differences in the parametrisation of sourceproperties which cannot be fully accounted for withoutchanging the code itself, something we postpone to futurework.
A. Comparing to the analytical approach of BLPF
A first analytical calculation of the 21-cm power spec-trum at cosmic dawn has been performed by Barkana andLoeb [14], focusing on the perturbations induced by theLyman- α coupling. Their model has been extended to in-clude temperature fluctuations by Pritchard and Furlan-etto [27] which is why we abbreviate it as BLPF. In con-trast to the halo model approach, which is centred on theradiation sources as building blocks, the BLPF approachfocuses on the gas, calculating the light-cone effects fromthe surrounding sources. A good summary of the modelcan be found in Ref. [70].For this comparison we have implemented the BLPFmodel, mostly following the descriptions in Ref. [27]. Thesource parametrisation, however, which includes the halomass function, the bias, and the spectral energy distri-bution, has been adapted to the description presented inthis paper. We furthermore omit any shot-noise correc-tions, since there is no fully worked-out model that in-cludes both shot-noise of the Lyman- α and temperaturefluctuations [27]. Finally, the star-formation rate den-sity is based on the EPS accretion rate modelling (see1 FIG. 5. Comparison between the 21-cm halo model developed in this paper (solid lines) and the BLPF method from Refs. [14, 27](coloured bands) assuming astrophysical parameters from the benchmark model B (see Fig 4).
Top panels: global signal (left)and the total 21-cm power spectrum with respect to k -modes (centre) and redshift (right). Middle panels:
Auto power spectraof the Lyman- α coupling (left), temperature (centre), and mass (right) components. Bottom panels:
Corresponding cross powerspectra of the individual components.
Eq. 13) and not on the time derivative of the collapsefraction. These changes with respect to the original workof Refs. [14, 27] allow us to carry out a fair comparison,where any resulting discrepancies can be fully attributedto differences between the methods and are not the resultof different source descriptions.The comparison between the BLPF and the halo modelis performed using the source parametrisation of bench-mark model B. This means we assume the radiation fluxparameters N α = 10000 and f X = 1, a non-truncateddouble-power law for the star formation efficiency, as wellas power-law spectra with indices α α = 0 and α X = 1 . k -modes for a selection of fourdifferent redshifts. The redshift values are chosen to lie inthe heating dominated regime ( z = 13), in the transitionregime of maximum absorption ( z = 16), in the Lyman- α regime ( z = 19), and at the epoch of the very first stars( z = 22). The solid lines correspond to the 21-cm halomodel, while the coloured bands show the results from2the BLPF model. In general, the BLPF model predictsless power, especially during the heating (black) and, toa lesser extend, during the Lyman- α epochs (dark green).The differences between the models are typically of theorder of a few, but they can grow to about an order ofmagnitude for specific redshifts and k -ranges. At verysmall scales above k ∼ k ∼ .
01, 0 .
06, 0 .
2, and0 .
95 h/Mpc. As before, we observe a good qualitativeagreement between the models, both of them showing acharacteristic double-peak feature for low k -modes merg-ing into one single peak at k ∼ k -values in the halo model compared to the BLPFmodel. In general, we conclude that the two models dif-fer by no more than a factor of a few in their redshiftevolution, with some exceptions where the difference cangrow to about an order of magnitude at most.In the middle and bottom rows of Fig. 5 we show theauto and cross power spectrum for the individual compo-nents α , T , and m . We observe as a general rule that thehalo model and BLPF power spectra are well convergedat the largest scales (lowest k values) before they startto slowly diverge towards higher k -values. Beyond k ∼ P αα and P T T auto spectra as well asthe P αT cross spectrum. The P mm auto spectrum, onthe other hand, shows very little differences except at k > k and z -values the differences can grow to about an order ofmagnitude. We want to emphasise, however, that at thispoint we do not know which of the two models is moreaccurate. Although former findings by Ref. [30] suggestthat the BLPF model lacks power with respect to semi-numerical calculations as well, we want to remind thatthese findings were based on a different source parametri-sation and are therefore not directly comparable. B. Comparing to
The code [32, 86] is based on a semi-numerical approach to predict 21-cm maps by solving thespin temperature evolution and the reionization processon a three-dimensional grid. The matter field is evolvedaccording to a first-order Zel’Dovich displacement [87]. The sources are not resolved individually, but their dis-tribution is calculated using an excursion-set method di-rectly applied to the matter field.A one-to-one comparison between our model and is not straight-forward because of subtle dif-ferences in the parametrisation and implementation thatmay significantly affect the results. For example, we donot exactly know which halo bias and halo mass functionagrees best with the excursion set implementation usedin . Furthermore, there are small differences inthe parametrisations of the spectral energy distributionsand the star-formation efficiency, that might be of rele-vance.In the code, the stellar-to-halo mass ra-tio is defined as a power law, followed by an exponen-tial cutoff towards small scales (with default parame-ters f ∗ , = 0 .
05 for the amplitude, α ∗ = 0 . M turn = 5 × M (cid:12) for the massscale of the exponential downturn, see Ref. [11]). Weattempt to reproduce the same functional form withoutchanging the parametrisation described in Eq. (11). Thebest agreement is found with the parameters f ∗ , = 0 . M p = 10 M (cid:12) , γ = − . γ = − . M t = 2 × M (cid:12) , γ = 1 .
4, and γ = −
4. A comparison of the twofunctions is shown in the top-left panel of Fig. 6.Regarding the halo mass function, we use the Sheth-Tormen model with q = 0 .
707 and p = 0 .
3. Althoughthis mass function is not such a good fit to high-redshift N -body simulations (see Ref. [20, 57] and discussions inSec. II A), it is used as the reference for the collapse frac-tion calculated in [see e.g. Eq. 14 in 32]. Forthe halo bias, on the other hand, we apply Eq. (3) withthe Press-Schechter (PS) parameters q = 1 and p = 0.This is because in , sources are populated withrespect to the PS conditional mass function. Finally, wedo not know which prescription for the halo-accretionrate is supposed to match best with the algorithm of . For this reason we decide to show the EPSas well as the EXP approach, since both are in goodagreement with high-redshift N -body simulations (seeSec. II B).In order to allow for a meaningful comparison, we se-lect the flux parameters N α and f X so that we obtain agood match to the global signal. The spectralenergy range and power-law index of the X-ray radia-tion, on the other hand, is kept fix at the default valuesof (which are E min = 0 . E max = 2 keV,and α X = 1).In the top-right panel of Fig. 6, we plot the globaldifferential brightness temperature. The results fromthe 21-cm halo model with EPS and EXP accretion areshown as solid and dash-dotted lines, while the globalsignal from is plotted as broad blue band. Al-though the agreement between the lines and the bandis very good (which is not so surprising, since we haveselected the flux parameters N α and f X to obtain thebest fit to the global signal), the curves arenot identical. The differences are of order 10 percent or3 FIG. 6. Comparison between the halo model and 21cmFAST.
Top-left:
Star formation efficiency of Eq. 11 fitted to the one from21cmFAST (blue line and band, respectively).
Top-right:
Global differential brightness temperature obtained by assuming theEXP and EPS accretion rates (dash-dotted and solid lines) and matched to the 21cmFast results (band) by refitting the totalLyman- α and X-ray photon numbers. Centre and Bottom:
Corresponding power spectra as a function of k -modes and redshift.Note that redshifts below z ∼
10 should not be trusted since effects related to reionisation become important. less, and could either stem from unaccounted deviationsin the source parametrisation or from the different waysthe global signal is calculated.The central panel of Fig. 6 shows the power spectrumas a function of k -modes for three different redshifts, rep-resenting the regime dominated by heating ( z ∼ . z ∼ .
9, blue), and theregime dominated by the Lyman- α coupling ( z ∼ .
0, cyan). The results from are again shown asbroad coloured bands. Note that we plot two runs withbox-size L = 300 and 450 Mpc, respectively (where thenumber of low- and high-resolution cells are kept constantat 300 and 1200 ). The resulting power spectra from the21-cm halo model are plotted as solid and dashed lines.The bottom panel of Fig. 6 illustrates the power spec-trum, this time as a function of redshift. At large scales(low k -values) the two characteristic peaks due to the4Lyman- α coupling and the heating epochs are clearly vis-ible at z ∼
11 and 14. Arounf k ∼ . andthe 21-cm halo model is best at very large and very smallscales. In between, at the transition scale (see k = 0 . .
44 h/Mpc), there is some visible discrepancies be-tween the two models, with predicting morepower by a factor of a few between z ∼
11 and 13. Fur-thermore, there are increasing differences towards verylow and very high redshifts. Note, however, that at red-shifts below z ∼
10 the halo model results cannot betrusted because effects from reionisation are ignored. Atvery high redshifts, on the other hand, both predictionsfrom and the halo model become increasinglyuncertain due to the poor knowledge of source numbersand distributions.As a summary, let us emphasise that we find goodqualitative agreement between and the 21-cmhalo model. At the quantitative level, the differencesdo not exceed a factor of ∼
3, except at very low andvery high redshifts, where the models cannot be fullytrusted. Note that this level of disagreement is of thesame order than the expected systematic effects due tomodelling choices regarding the halo accretion, the haloabundance, and the halo bias (see Appendix A and B formore details). We therefore conclude that, at the currentstage, it is impossible to know if the observed differencesbetween the 21-cm halo model and are a resultof different modelling choices or if they hint towards morefundamental issues regarding the halo model approach.
VI. CONCLUSIONS
Fast and sufficiently accurate predictions for the 21-cmclustering signal are important in order to develop a bet-ter understanding of the cosmic dawn, the epoch of thehigh-redshift universe right before the phase transitionfrom neutral to ionised hydrogen. Many important as-pects of the prediction pipeline, especially related to thesource modelling, remain unknown, and fast models mayhelp to explore the vast parameter space of possibilities.In this paper we present a new analytical method basedon the framework of the halo model, where Lyman- α cou-pling and temperature fluctuations are described with thehelp of overlying flux profiles that properly include red-shifting and source attenuation due to the expansion ofthe universe and the finite speed of light. The model pro-vides a natural framework to predict all auto and crosspower spectra of the Lyman- α coupling, the temperaturefluctuations, and the matter perturbations. The temper-ature fluctuations are separated in a primordial compo-nent sourced by the matter perturbations and a heatingcomponent induced by the first sources. The effects formthe process of reionisation could, in principle, be addedto the framework, but this is left as future work.A distinctive advantage compared to other analytical methods is that the halo model approach naturally in-cludes shot-noise effects induced by the low number ofsources in the early universe. This is important for highredshifts and small scales, especially when investigatingmodels with very luminous sources.Next to presenting the framework of a 21-cm halomodel, we investigate the effects of the stellar accretion,the source abundance, and the halo bias. We show thatdifferent choices regarding these model ingredients canstrongly affect both the global signal absorption troughand the 21-cm power spectrum. If the flux parametersare re-fitted to correct for the shift in the global signal,the effect on the power spectrum becomes smaller, butthe difference can still be as large as a factor of a few.This means that analytical and semi-numerical methodsrequire a precise understanding of structure formation atvery early times in order to avoid significant systematicerrors.In order to check the general validity of the 21-cm halomodel, we compare it to the earlier analytical model byBarkana, Loeb, Pritchard, and Furlanetto (BLPF) intro-duced in Refs. [14, 27]. In general, we find good quali-tative agreement between the two approaches regardingboth the k -mode and redshift evolution. At the quanti-tative level, the halo model predicts a somewhat largerclustering signal compared to the BLPF model. However,the difference is no more than a factor of a few, exceptfor specific redshifts and k -values where it can grow toabout an order of magnitude.Furthermore, we compare the 21-cm halo model to thesemi-numerical code [32, 86]. Although an ex-act comparison between the two methods is currently un-feasible due to small differences in the source parametri-sations, we nevertheless find a very encouraging agree-ment. At redshifts before the onset of reionisation, thedifference between the halo model and staysbelow a factor of ∼ k -modes and redshifts). Whether the re-maining differences between the models are due to thespecifics of the source modelling, or whether they are asign for a more fundamental failure of the halo modelremains to be investigated in the future.In general, a detailed apple-to-apple comparison in-cluding analytical, semi-numerical, and full simulation-based calculations of the 21-cm signal would be extremelybeneficial for the community. Only a combined use ofaccurate but very expensive simulations together withmuch faster approximate methods will lead to a morecomplete understanding of the complex and rich signalfrom the epoch of the cosmic dawn. ACKNOWLEDGMENTS
We thank Romain Teyssier and Kai Hofmann for veryhelpful suggestions related to the model. This work issupported by the Swiss National Science Foundation viathe grant
PCEFP2 181157. FIG. 7. Mass accretion modelling and how it affects the 21-cm signal assuming the source parameters of fiducial model B.
Top-left:
Halo growth for the three accretion models EPS, EXP, and AM (see Sec. II B).
Top-centre:
Star-formation ratedensity for the same models.
Top-right:
Resulting global 21-cm signal, where the flux parameters of the EXP and AM models( N α , f X ) are either kept the same (solid lines) or where they are modified so that the minimum of the absorption signal lie atexactly the same redshift the the one from the EPS model (dashed lines). Centre and Bottom:
Resulting power spectra as afunction of k -modes and redshift. Appendix A: Effects from the mass accretionmodelling
In Sec II B we have investigated different modellingchoices for the halo growth rate based on an extendedPress-Schechter (EPS) prescription, an exponential ac- cretion rate (EXP), and an abundance-matching tech-nique (AM). We have compared these cases to simula-tions of B20 [58] before selecting the EPS model as ourstandard method for the paper.In this appendix we have a closer look at the EPS,EXP, and AM halo accretion rates, focusing on how they6affect the 21-cm global signal and power spectrum. Wethereby assume the fiducial model B for the source pa-rameters. More details about model B can be found inSec. IV D.The top panels of Fig. 7 show the halo growth (left),the star-formation rate density (centre), and the global21-cm signal (right) for the EPS, EXP, and AM models.Note that the absorption trough of the global signal isshifted towards higher redshifts when going from the AMto the EXP prescription.In order to discriminate between effects originatingfrom the global signal and the ones affecting the powerspectrum, we allow the Lyman- α and X-ray flux param-eters N α and f X to be renormalised so that the globalsignal absorption troughs of the AM and EXP modelsare aligned with the ones from the EPS model. The re-sulting signal from such a renormalised flux analysis aregiven by the dashed lines in Fig. 7.The power spectrum as a function of k -modes andredshift obtained from the different accretion ratemodels are plotted in the middle and bottom panels ofFig. 7. There are large differences of up to an orderof magnitude between the models, clearly exceedingthe differences observed in the global signal alone.The differences are smaller but still significant whenlooking at the models with matching global signalsdue to renormalised flux parameters (see dashed anddark-green solid lines). The power spectrum from thesemodels can still differ by about a factor of 3 at most. Wetherefore conclude that the mass accretion rate has tobe modelled at high accuracy in order to obtain reliablepredictions of the 21-cm power spectrum. Appendix B: Effects from halo bias and massfunction prescriptions
In Appendix A we have shown that modelling choicesregarding the mass accretion rate may significant affectthe 21-cm signal. Here we focus on another central mod- elling component: the halo mass function and corre-sponding bias prescription introduced in Eqs. (1-3).In order to test the sensitivity of the 21-cm global sig-nal and power spectrum, we vary the free model param-eters q and p of the first crossing distribution f ( ν ) ofEq. (2). We thereby investigate three main cases, theoriginal Press and Schechter [52] model ( q = 1, p = 0),the Sheth and Tormen [54] model ( q = 0 . p = 0 . q = 0 . p = 0 .
3) that isin best agreement with results from high-redshift simu-lations (see discussion in Sec. II A).The corresponding halo bias and mass functions of thethree models are shown in Fig. 8 (top-left and top-centrepanels). Note that the difference between the models be-come quite large in both cases, especially towards largehalo masses. As a result, the global 21-cm signal is shiftedsubstantially, the maximum of the absorption troughmoving from z ∼
16 to z ∼
18 for the case of the STinstead of the PS prescription. In order to counteractthis effect at the level of the global signal, we show twomore models where the flux parameters N α and f X aremodified in order to bring the absorption trough backin line with the intermediate model (dashed lines). Note,however, that there are remaining differences in the shapeof the absorption signal that cannot be absorbed by sucha simple recalibration.In the remaining panels of Fig. 8, we show the 21-cmpower spectrum as a function of k -modes and redshift.There are significant differences, especially between theST model (brown line) and the other two cases (solid redand purple lines). The recalibrated models, on the otherhand, are much closer together, but differences of abouta factor of two remain.We conclude that next to the halo accretion rate, thehalo bias and mass functions need to be known to goodaccuracy in order to avoid large modelling errors with thehalo model approach. Note, that similar uncertaintiesalso exist for semi-numerical methods or even simulationresults that may only be trusted if they have realistichalo abundance and distributions. [1] M. V. Haarlem et al. , Electronics , 54 (2006).[2] F. Mertens et al. , Mon. Not. Roy. Astron. Soc. , 1662(2020), arXiv:2002.07196 [astro-ph.CO].[3] S. Tingay et al. , Publications of the Astronomical Societyof Australia , 1 (2013).[4] C. M. Trott et al. , Mon. Not. Roy. Astron. Soc. , 4711(2020), arXiv:2002.02575 [astro-ph.CO].[5] D. R. DeBoer et al. , Publ. Astron. Soc. Pac. , 045001(2017), arXiv:1606.07473 [astro-ph.IM].[6] L. Koopmans et al. , PoS AASKA14 , 001 (2015),arXiv:1505.07568 [astro-ph.CO].[7] P. Madau, A. Meiksin, and M. J. Rees, Astrophys. J. , 429 (1997), arXiv:astro-ph/9608010.[8] P. Shaver, R. Windhorst, P. Madau, and A. de Bruyn, Astron. Astrophys. , 380 (1999), arXiv:astro-ph/9901320.[9] A. Ewall-Wice, T.-C. Chang, J. Lazio, O. Dore, M. Seif-fert, and R. Monsalve, Astrophys. J. , 63 (2018),arXiv:1803.01815 [astro-ph.CO].[10] J. Mirocha et al. , (2019), arXiv:1903.06218 [astro-ph.CO].[11] J. Park, A. Mesinger, B. Greig, and N. Gillet, Mon.Not. Roy. Astron. Soc. , 933 (2019), arXiv:1809.08995[astro-ph.GA].[12] M. McQuinn, O. Zahn, M. Zaldarriaga, L. Hernquist,and S. R. Furlanetto, Astrophys. J. , 815 (2006),arXiv:astro-ph/0512263.[13] A. Liu et al. , (2019), arXiv:1903.06240 [astro-ph.CO]. FIG. 8. Effects from the mass function and bias modelling on the 21-cm signal, assuming source parameters from the fiducialmodel B.
Top-left and centre:
Halo bias and mass functions for different values of q and p . The purple and brown linescorrespond to the Press-Schechter ( q = 1, p = 0) and Sheth-Tormen ( q = 0 . p = 0 .
3) mass functions, while the red lineshows an case in-between ( q = 0 . p = 0 .
3) which is in better agreement with high redshift simulations.
Top-right:
Global21-cm signal resulting from these mass function and bias prescriptions, where the flux parameters ( N α , f X ) are either kept thesame (solid lines) or where they are changed so that the absorption troughs lie at the same redshift (dashed lines). Centre andBottom:
Resulting power spectra as a function of k -modes and redshsifts.[14] R. Barkana and A. Loeb, Astrophys. J. , 1 (2005),arXiv:astro-ph/0410129 [astro-ph].[15] B. Greig et al. , (2020), arXiv:2006.03203 [astro-ph.CO].[16] H. Liu and T. R. Slatyer, Phys. Rev. D , 023501 (2018),arXiv:1803.09739 [astro-ph.CO]. [17] O. Mena, S. Palomares-Ruiz, P. Villanueva-Domingo,and S. J. Witte, Phys. Rev. D , 043540 (2019),arXiv:1906.07735 [astro-ph.CO].[18] M. Sitwell, A. Mesinger, Y.-Z. Ma, and K. Sigurd-son, Mon. Not. Roy. Astron. Soc. , 2664 (2014), arXiv:1310.0029 [astro-ph.CO].[19] L. Lopez-Honorez, O. Mena, A. Molin´e, S. Palomares-Ruiz, and A. C. Vincent, JCAP , 004 (2016),arXiv:1603.06795 [astro-ph.CO].[20] A. Schneider, Phys. Rev. D , 063021 (2018),arXiv:1805.00021 [astro-ph.CO].[21] A. Lidz and L. Hui, Phys. Rev. D , 023011 (2018),arXiv:1805.01253 [astro-ph.CO].[22] O. Nebrin, R. Ghara, and G. Mellema, JCAP , 051(2019), arXiv:1812.09760 [astro-ph.CO].[23] J. C. Hill and E. J. Baxter, JCAP , 037 (2018),arXiv:1803.07555 [astro-ph.CO].[24] W. Yang, S. Pan, S. Vagnozzi, E. Di Valentino, D. F.Mota, and S. Capozziello, JCAP , 044 (2019),arXiv:1907.05344 [astro-ph.CO].[25] L. Lopez-Honorez, O. Mena, S. Palomares-Ruiz,P. Villanueva-Domingo, and S. J. Witte, JCAP , 026(2020), arXiv:2004.00013 [astro-ph.CO].[26] S. Furlanetto, M. Zaldarriaga, and L. Hernquist, Astro-phys. J. , 1 (2004), arXiv:astro-ph/0403697.[27] J. R. Pritchard and S. R. Furlanetto, Mon. Not. Roy.Astron. Soc. , 1680 (2007), arXiv:astro-ph/0607234.[28] J. Raste and S. Sethi, Astrophys. J. , 55 (2018),arXiv:1711.03827 [astro-ph.CO].[29] J. Raste and S. Sethi, Astrophys. J. , 56 (2019),arXiv:1810.12746 [astro-ph.CO].[30] M. G. Santos, A. Amblard, J. Pritchard, H. Trac,R. Cen, and A. Cooray, Astrophys. J. , 1 (2008),arXiv:0708.2424 [astro-ph].[31] A. Mesinger and S. Furlanetto, Astrophys. J. , 663(2007), arXiv:0704.0946 [astro-ph].[32] A. Mesinger, S. Furlanetto, and R. Cen, Monthly Noticesof the Royal Astronomical Society , 955 (2010).[33] M. Santos, L. Ferramacho, M. Silva, A. Amblard, A. C.C. . Ist, and U. Irvine, Monthly Notices of the RoyalAstronomical Society , 2421 (2009).[34] A. Fialkov, R. Barkana, A. Pinhas, and E. Visbal, Mon.Not. Roy. Astron. Soc. , 36 (2014), arXiv:1306.2354[astro-ph.CO].[35] S. Mutch, P. M. Geil, G. Poole, P. Angel, A. Duffy,A. Mesinger, and J. Wyithe, Monthly Notices of theRoyal Astronomical Society , 250 (2016).[36] A. Hutter, Mon. Not. Roy. Astron. Soc. , 1549 (2018),arXiv:1803.00088 [astro-ph.CO].[37] R. M. Thomas et al. , Mon. Not. Roy. Astron. Soc. ,32 (2009), arXiv:0809.1326 [astro-ph].[38] R. Ghara, T. R. Choudhury, and K. K. Datta, Mon.Not. Roy. Astron. Soc. , 1806 (2015), arXiv:1406.4157[astro-ph.CO].[39] G. Mellema, I. Iliev, M. Alvarez, and P. Shapiro, NewAstron. , 374 (2006), arXiv:astro-ph/0508416.[40] M. Friedrich, G. Mellema, I. Iliev, and P. Shapiro,Monthly Notices of the Royal Astronomical Society ,2232 (2012).[41] L. Graziani, A. Maselli, and B. Ciardi, Monthly Noticesof the Royal Astronomical Society , 722 (2013).[42] R. Ghara, G. Mellema, S. K. Giri, T. R. Choudhury,K. K. Datta, and S. Majumdar, Mon. Not. Roy. Astron.Soc. , 1741 (2018), arXiv:1710.09397 [astro-ph.CO].[43] H. E. Ross, K. L. Dixon, R. Ghara, I. T. Iliev, andG. Mellema, Mon. Not. Roy. Astron. Soc. , 1101(2019), arXiv:1808.03287 [astro-ph.CO].[44] M. B. Eide, B. Ciardi, L. Graziani, P. Busch, Y. Feng,and T. Di Matteo, Mon. Not. Roy. Astron. Soc. , 6083 (2020), arXiv:2009.06631 [astro-ph.CO].[45] P. Ocvirk et al. , Mon. Not. Roy. Astron. Soc. , 1462(2016), arXiv:1511.00011 [astro-ph.GA].[46] A. Doussot, H. Trac, and R. Cen, Astrophys. J. , 18(2019), arXiv:1712.04464 [astro-ph.CO].[47] P. Ocvirk et al. , Mon. Not. Roy. Astron. Soc. , 4087(2020), arXiv:1811.11192 [astro-ph.GA].[48] J. Rosdahl, H. Katz, J. Blaizot, T. Kimm, L. Michel-Dansac, T. Garel, M. Haehnelt, P. Ocvirk, andR. Teyssier, Monthly Notices of the Royal AstronomicalSociety , 994 (2018).[49] H. Katz et al. , Mon. Not. Roy. Astron. Soc. , 2200(2020), arXiv:1905.11414 [astro-ph.GA].[50] B. Semelin, E. Eames, F. Bolgar, and M. Cail-lat, Mon. Not. Roy. Astron. Soc. , 4508 (2017),arXiv:1707.02073 [astro-ph.CO].[51] L. N. Holzbauer and S. Furlanetto, Monthly Notices ofthe Royal Astronomical Society , 718 (2011).[52] W. H. Press and P. Schechter, Astrophys. J. , 425(1974).[53] J. R. Bond, S. Cole, G. Efstathiou, and N. Kaiser, As-trophys. J. , 440 (1991).[54] R. K. Sheth and G. Tormen, Mon. Not. Roy. Astron. Soc. , 61 (2002), arXiv:astro-ph/0105113 [astro-ph].[55] R. K. Sheth and G. Tormen, Mon. Not. Roy. Astron. Soc. , 119 (1999), arXiv:astro-ph/9901122.[56] A. Cooray and R. K. Sheth, Phys. Rept. , 1 (2002),arXiv:astro-ph/0206508.[57] I. Iliev, G. Mellema, P. Shapiro, U.-L. Pen, Y. Mao,J. Koda, and K. Ahn, Monthly Notices of the RoyalAstronomical Society , 2222 (2011).[58] P. Behroozi et al. , (2020), 10.1093/mnras/staa3164,arXiv:2007.04988 [astro-ph.GA].[59] J. McBride, O. Fakhouri, and C. Ma, Monthly Noticesof the Royal Astronomical Society , 1858 (2009).[60] A. Dekel, A. Zolotov, D. Tweed, M. Cacciato, D. Cev-erino, and J. Primack, Mon. Not. Roy. Astron. Soc. ,999 (2013), arXiv:1303.3009 [astro-ph.CO].[61] H. Trac, R. Cen, and P. Mansfield, Astrophys. J. ,54 (2015), arXiv:1507.02685 [astro-ph.CO].[62] E. Neistein and F. C. van den Bosch, Mon. Not. Roy.Astron. Soc. , 933 (2006), arXiv:astro-ph/0605045.[63] C. A. Correa, J. S. B. Wyithe, J. Schaye, and A. R.Duffy, Mon. Not. Roy. Astron. Soc. , 1514 (2015),arXiv:1409.5228 [astro-ph.GA].[64] A. Vale and J. Ostriker, Mon. Not. Roy. Astron. Soc. , 189 (2004), arXiv:astro-ph/0402500.[65] S. R. Furlanetto, J. Mirocha, R. H. Mebane, andG. Sun, Mon. Not. Roy. Astron. Soc. , 1576 (2017),arXiv:1611.01169 [astro-ph.GA].[66] R. H. Mebane, J. Mirocha, and S. R. Furlan-etto, Mon. Not. Roy. Astron. Soc. , 4544 (2018),arXiv:1710.02528 [astro-ph.GA].[67] C. Mason, M. Trenti, and T. Treu, Astrophys. J. ,21 (2015), arXiv:1508.01204 [astro-ph.GA].[68] G. Sun and S. R. Furlanetto, Monthly Notices of theRoyal Astronomical Society , 417 (2016).[69] J. Mirocha, S. R. Furlanetto, and G. Sun, Mon. Not.Roy. Astron. Soc. , 1365 (2017), arXiv:1607.00386.[70] J. R. Pritchard and A. Loeb, Rept. Prog. Phys. ,086901 (2012), arXiv:1109.6012 [astro-ph.CO].[71] J. Mirocha, (2020), 10.1093/mnras/staa3150,arXiv:2008.04322 [astro-ph.GA].[72] J. R. Pritchard and S. R. Furlanetto, Mon. Not. Roy. Astron. Soc. , 1057 (2006), arXiv:astro-ph/0508381[astro-ph].[73] R. Barkana and A. Loeb, Astrophys. J. , 474 (2004),arXiv:astro-ph/0310338.[74] J. Mirocha, Mon. Not. Roy. Astron. Soc. , 1211(2014), arXiv:1406.4120 [astro-ph.GA].[75] S. Furlanetto, Mon. Not. Roy. Astron. Soc. , 867(2006), arXiv:astro-ph/0604040.[76] S. Furlanetto and S. J. Stoever, Monthly Notices of theRoyal Astronomical Society , 1869 (2010).[77] S. A. Wouthuysen, AJ , 31 (1952).[78] G. B. Field, Proceedings of the IRE , 240 (1958).[79] X.-L. Chen and J. Miralda-Escude, Astrophys. J. , 1(2004), arXiv:astro-ph/0303395.[80] C. M. Hirata, Mon. Not. Roy. Astron. Soc. , 259(2006), arXiv:astro-ph/0507102 [astro-ph]. [81] J. Mirocha, G. J. A. Harker, and J. O. Burns, Astrophys.J. , 11 (2015), arXiv:1509.07868 [astro-ph.CO].[82] M. Gilfanov, H. Grimm, and R. Sunyaev, Mon. Not. Roy.Astron. Soc. , L57 (2004), arXiv:astro-ph/0301331.[83] S. Bharadwaj and S. S. Ali, Mon. Not. Roy. Astron. Soc. , 142 (2004), arXiv:astro-ph/0401206.[84] R. Barkana and A. Loeb, Astrophys. J. Lett. , L65(2005), arXiv:astro-ph/0409572.[85] Note that although our parametrisation is similar tothat of Mirocha et al. [69], the global signal is shiftedto slightly higher redshift due to differing assumptionsabout the star-formation efficiency normalisation and theminimum halo mass.[86] S. G. Murray, B. Greig, A. Mesinger, J. B. Mu˜noz,Y. Qin, J. Park, and C. A. Watkinson, (2020),10.21105/joss.02582, arXiv:2010.15121 [astro-ph.IM].[87] Y. B. Zel’Dovich, Astron. Astrophys.500