An expanded merger-tree description of cluster evolution
aa r X i v : . [ a s t r o - ph . C O ] O c t Mon. Not. R. Astron. Soc. , 1–11 (2010) Printed 13 November 2018 (MN L A TEX style file v2.2)
An expanded merger-tree description of cluster evolution
Irina Dvorkin ⋆ and Yoel Rephaeli , School of Physics and Astronomy, Tel Aviv University, Tel Aviv, 69978, Israel Center for Astrophysics and Space Sciences, University of California, San Diego, La Jolla, CA 92093-0424
13 November 2018
ABSTRACT
We model the formation and evolution of galaxy clusters in the framework of anextended dark matter halo merger-tree algorithm that includes baryons and incorpo-rates basic physical considerations. Our modified treatment is employed to calculatethe probability density functions of the halo concentration parameter, intraclustergas temperature, and the integrated Comptonization parameter for different clustermasses and observation redshifts. Scaling relations between cluster mass and theseobservables are deduced that are somewhat different than previous results. Modelinguncertainties in the predicted probability density functions are estimated. Our treat-ment and the insight gained from the results presented in this paper can simplifythe comparison of theoretical predictions with results from ongoing and future clustersurveys.
Key words: galaxies: clusters: general - large-scale structure of Universe
Formation of galaxy clusters is of central importance for un-derstanding the evolution of the large scale structure (LSS)of the universe. Statistical properties of clusters - deducedfrom cluster optical, X-ray, and SZ surveys - can be used todetermine the basic cosmological parameters - such as thematter density, the normalization (and spectrum) of primor-dial density fluctuations, and the dark energy equation ofstate - independently of other methods (CMB power spec-trum, galaxy surveys, etc). This can only be achieved if the-oretical tools are developed for a quantitative description ofthe non-linear hierarchical growth of clusters. Hydrodynam-ical numerical simulations of clusters are currently the mostversatile tool for establishing statistical properties of clus-ters, which are commonly specified in terms of mass func-tions and scaling relations between their intrinsic properties.In order to use cluster surveys as cosmological probes itis essential to know how their intrinsic properties - such asformation time, mass and redshift - affect their integratedstatistical properties. Clusters are commonly described asvirialized spherical systems whose masses are dominated bydark matter (DM) and with isothermal intracluster (IC) gasas their main baryonic mass component. Cluster statisticalrelations are usually expressed in terms of the mass and theredshift of observation; for example, the resulting scaling ofthe gas temperature is k B T ∝ M / (cid:2) E ( z )∆ V ( z ) (cid:3) / , (1) ⋆ E-mail: [email protected] where E ( z ) = H ( z ) /H , the Hubble parameter in units of itspresent value and ∆ V ( z ) is the overdensity at virialization.However, this simple description may not be sufficientlyadequate for the description of real clusters. Observationsand numerical simulations indicate that DM density profilesare shallower than isothermal at small radii, and steeperthan isothermal at large radii, and are characterized by ascaling radius r s which marks the transition between thesetwo regions (Navarro, Frenk & White 1995, 1996). Althoughit seems reasonable to assume that the scaling relations arenot much affected by the details of the cluster structure, thestandard description is not sufficiently accurate since themodel parameters - such as r s - depend on the cluster massand redshift. Indeed, numerical simulations show that theformation history of clusters affects the DM scaling radius,such that clusters are systematically denser the earlier theyformed. Cluster mass concentration is usually quantified interms of the parameter c , defined as the ratio of the clustervirial radius to the scaling radius, c = R v /r s . Thus, clus-ters that form earlier have larger concentration parameters(Wechsler et al. 2002); more generally, other physical prop-erties are also influenced by the cluster formation history.The assumption of hydrostatic equilibrium is obviouslyan approximation, since clusters, the largest bound systems,are still forming through mergers of subclumps and accre-tion. Mergers disrupt the state of thermal equilibrium; dur-ing some merger events IC gas temperature and X-ray lumi-nosity are boosted up by factors of up to 3 and 10, respec-tively, as was demonstrated in a series of numerical simu-lations by Ricker and Sarazin (2001). More generally, whenderiving the standard scaling relations one usually ignores c (cid:13) I. Dvorkin and Y. Rephaeli the formation history of the cluster. In the ΛCDM modelthe LSS formed hierarchically, through consecutive merg-ers of smaller structures, and thus clusters with roughly thesame mass and redshift may have very different merger andformation histories. Indeed, both observations and simula-tions reveal scatter in the mass-observable relations at alevel of ∼ − T X − M scatter. How-ever, the amount of scatter may depend on the mass andredshift of the cluster.A meaningful comparison with observational data usu-ally necessitates knowledge of the full probability distribu-tion function (PDF) rather than just the scaling with mass.An example is the case of very large values of the concen-tration parameter and Einstein radius for several clusters(Broadhurst et al. 2008, Zitrin et al. 2009). Based on resultsfrom N-body simulations, the probability of observing clus-ters with the very high measured values was found to bevery low, amounting to a 4 − σ discrepancy with the ΛCDMpredictions. The concentration parameter PDF is a crucialcomponent in this analysis (Sadeh & Rephaeli 2008). It isclear that a better understanding of the origin of this PDFis required.The impact of the formation history on the cluster prop-erties was previously investigated using a series of hydrody-namical simulations (Ricker and Sarazin 2001, Randall andSarazin 2002, Wik et al. 2008). These authors used simu-lations of pairs of merging clusters and analysed the im-pact of recent mergers on IC gas temperature, luminosity,and Comptonization parameter. They found that all thesequantities are boosted for a relatively long time followinga merger, and calculated the effect of this boost on cos-mological parameter estimation. These analyses did not in-clude the full formation history, only the latest merger, andrelied on a small number of simulated clusters with differ-ent masses. Voit and Donahue (1998) showed that the tem-perature evolves less rapidly with mass than in the stan-dard analysis when the recent formation approximation isrelaxed. They assumed gradual mass accretion throughoutthe cluster history, and a one-to-one correspondence betweenthe temperature and the virial energy of the cluster.The statistics of DM halo concentrations and their de-pendence on the halo formation history were investigated using N-body simulations (Bullock et al. 2001, Wechsler etal. 2002, Neto et al. 2007, Gao et al. 2008, Duffy et al. 2008).This approach provides PDFs of the concentration parame-ter which account for different formation histories of differ-ent halos. Results of these works can then be used to inferthe intrinsic scatter in other cluster observables. However, atheoretical approach that is complementary to N-body sim-ulations is needed in order to fully understand the impact ofthe cluster formation history on its properties. The reliabil-ity of the statistical analysis of N-body simulations dependson the simulation volume, which is limited by computationalconstraints. This can be a severe problem if one is interestedin high-mass clusters, which are relatively rare systems. Forexample, the Millennium Simulation (MS, Springel et al.2005), the largest cosmological N-body simulation to date,which follows N = 2160 particles in a periodic box of L = 500 h − Mpc on a side, contains less than 800 relaxedhalos with masses above M = 1 . · h − M ⊙ , and just8 relaxed halos with masses above M = 7 . · h − M ⊙ (Neto et al. 2007, hereafter N07). Moreover, comparisonbetween different simulations is difficult because each uti-lizes different cosmological parameters and halo finding al-gorithms. This difficulty is illustrated by the fact that some-what different mass functions are predicted by different N-body simulations (Jenkins et al. 2001, Sheth & Tormen 2002,Tinker et al. 2008), most likely reflecting the different for-mation histories predicted by these simulations.Predicting the full PDFs of the relevant cluster param-eters in the context of a theoretical model that can be read-ily implemented, would allow a more meaningful statisticalanalysis and the ability to quantify the impact of uncertain-ties in the values of cosmological and cluster parameters onthe main observables deduced from large scale surveys. Inthis paper we develop a model of cluster formation using an-alytically computed DM merger trees, with which we tracethe formation of clusters through major episodal mergersand continuous accretion. We show that our approach pro-vides an improved physical description in comparison withwhat can be obtained from standard scaling relations.This paper is organized as follows. In Section 2 we de-scribe our model of cluster formation. The results of ourgeneralized merger-tree treatment are presented in Section3 and further discussed in Section 4. Throughout the paperwe use the following cosmological parameters: Ω m = 0 . Λ = 0 . H = 73 km/s/Mpc, σ = 0 . Our description of the growth of galaxy clusters is basedon merger trees of DM halos as described numerically bythe modified GALFORM code (Cole et al. 2000; Parkin-son, Cole & Helly 2008), which was successfully employedto construct semi-analytic models of galaxy formation. Thealgorithm implements the excursion set formalism, a key as-pect of which is the conditional mass function: the fractionof mass f ( M | M ) from halos of mass M at redshift z thatconsisted of smaller halos of mass M at an earlier redshift c (cid:13) , 1–11 n expanded merger-tree description of cluster evolution z f ( M | M ) d ln M = r π σ ( δ − δ )( σ − σ ) / exp (cid:18) −
12 ( δ − δ ) ( σ − σ ) (cid:19) d ln σ d ln M d ln M , (2)where σ i = σ ( M i ) is the variance of the linear perturbationfield smoothed on scale M i , and δ i is the critical density forspherical collapse at redshift z i . The original GALFORM al-gorithm is consistent with the Press-Schechter mass function(Press & Schechter, 1974; PS), in the sense that if a grid oftrees is rooted at z = 0, weighted by their mass abundanceaccording to the PS mass function, then the mass functionat higher redshifts again corresponds to PS mass function.The conditional mass function is used to calculate the meannumber of progenitors of mass M at redshift z + dz of ahalo of mass M at z = z : dNdM = 1 M (cid:18) dfdz (cid:19) z = z M M dz (3)The modified GALFORM algorithm was obtained by mak-ing the substitution dNdM → dNdM G ( σ /σ , δ /σ ) , (4)where G is referred to as a perturbing function, to be cal-ibrated by comparison with N-body simulations. Parkin-son et al. (2008) showed that by fitting the outcome ofthe algorithm to the results of the MS they obtained haloabundances which are consistent with the Sheth-Tormenmass function (Sheth & Tormen 2002). Thus, the perturb-ing function G expresses the uncertainty in the choice ofthe correct mass function. In this work we use the fol-lowing parametrization of the perturbing function: G = G ( σ /σ ) γ ( δ /σ ) γ , with the parameters G = 0 . γ = 0 . γ = − .
01 taken from Parkinson et al. (2008).Starting with the specified mass and redshift, the algo-rithm proceeds back in time, checking after each timestepwhether there was a merger; if so, the masses of the mergedhalos are drawn from the distribution (3). Halos with massesbelow some resolution limit M res are not resolved, and areaccounted for as continuously accreted mass. Further detailson the GALFORM algorithm can be found in Parkinson etal. (2008).In the following sections we describe our model of clus-ter formation. Since the DM is the dominant mass com-ponent of clusters we first study the formation history ofcluster-sized DM halos, and then add the IC gas compo-nent, study its properties and related scaling relations. The calculation begins with the construction of a merger treefor a given final halo mass and redshift. For each tree, weuse only the major merger events, that is only those merg-ers for which M > /M < < q , whose value is to be determined.The rationale behind this choice is the assumption that theproperties of a cluster are largely determined by the violentmerger events, during which DM and gas settle in the modi-fied potential well when equilibrium is reestablished. On theother hand, slow accretion of material on the outskirts of the cluster, as well as mergers with low-mass systems (which aretreated in an identical manner in this work) do not cause re-distribution of the cluster components and do not stronglyaffect the physical conditions in the cluster center. Thus, wefollow the major merger events in the cluster history, andrefer to all other processes of mass growth as slow accre-tion. A typical value of the major merger ratio is q = 10;other values will also be considered. The trees were gener-ated up to the redshift z = z fin which depends on the massresolution of the merger tree algorithm, M res , as discussedbelow.The next step is to calculate the density profile of eachhalo in the tree. To this end, the tree is traced down, begin-ning with the smaller masses; for each merger event we useenergy conservation to calculate the density profile of themerged DM halo. For simplicity we assume NFW (Navarro,Frenk & White 1995) profiles for all the halos at all times: ρ DM ( r ) = ρ s x (1 + x ) (5)where x = r/r s = rc/R v is the radial distance expressedin terms of NFW scale radius r s , c is the concentration pa-rameter, and ρ s is the mass density normalization constant.Each halo is completely characterized by its mass, redshiftand concentration parameter, where halo mass is definedas the mass within the virial radius, M = π R v ρ c ∆ V . Allmasses are assumed to be in equilibrium, which is presumedto be attained relatively quickly after each merger event, butwe do not assume the halos are completely virialized withinthe virial radius. In fact, the virial ratio 2 T / | W | approaches1 at the virial radius only for very large c , while for com-monly deduced values of c this ratio is slightly larger than1 at the virial radius (Cole & Lacey 1996; Lokas & Mamon2001).The relation between mass and virial radius depends onredshift, which we identify as the halo formation redshift, z f , defined as follows. If the halo has undergone a majormerger event, we take z f to be the redshift of the last majormerger. This choice is consistent with the main assumptionthat major mergers largely determine the physics of the halo.For halos that did not experience major mergers at all, andwere, according to our interpretation, entirely assembled byminor mergers and continuous accretion, we take z f to bethe redshift at which half of the halo mass had assembled.Difficulty in determining z f by this prescription is encoun-tered only when the branch of the merger tree terminateswhen the halo still has more than half of its mass (recallthat the tree is evolved backwards in time). In GALFORM,the branch is terminated in two cases: either z fin is reached,or the mass of the halo falls below the resolution mass M res .Thus, in order to describe the formation of the smallest ha-los that constitute the tree, z fin should clearly be chosenwell above the expected formation redshift of halos of mass M res . By making this choice we ensure that almost all of thehalos in the tree assemble half of their mass before z fin isreached (that is, later in time), and their formation redshiftcan be traced back by the tree.To start from the earliest halos in the tree and moveforward in time requires specifying their concentration pa-rameters. These are adapted from a fit to a set of N-body c (cid:13) , 1–11 I. Dvorkin and Y. Rephaeli simulations by Bullock et al. (2001) c ( M, z obs ) = 5 (cid:18) M h − M ⊙ (cid:19) − . (1 + z f )(1 + z obs ) , (6)where z f and z obs are the redshifts of formation and obser-vation, respectively. This choice is motivated by the findingof Wechsler et al. (2002) that the concentration parameterscales as c ∼ (1 + z f ) / (1 + z obs ), although note that theirdefinition of the formation redshift is slightly different thanthe one adopted here. Although this choice for the initial c ( M, z ) is somewhat arbitrary, its particular form does notsignificantly influence the results.For each merger event, we calculate the total energyof the system before merging, which depends on the con-centration parameters of the merging halos, and deduce theconcentration parameter of the merged halo from simple en-ergy conservation arguments, motivated by a cluster mergermodel by Sarazin (2002). The total energy of a system oftwo halos before merging is: E tot, = E ( M ) + E ( M ) + U (7)where E ( M i ) are the total energies (potential and kinetic)of each halo and U is the gravitational energy of the twohalos at the point of their largest separation, when they havejust become bound and their relative velocity was negligible: U = − GM M d (8)The distance d , at which the halos became bound, isroughly the mean distance between halos with masses M and M that reside in an overdense region with a scale thatcorresponds to the final mass M f = M + M . We take thisdistance to be d = κd , where d = R + R , and adopt κ = 5as a fiducial value for all halos. This corresponds to a typ-ical distance of several Mpc and a typical relative velocityof several hundreds to a few thousands of km/s, dependingon the masses of the merging clusters, which is in accor-dance with the initial conditions of hydrodynamical clustermerger simulations (Ricker & Sarazin 2001, McCarthy et al.2007, Lee and Komatsu 2010). Very high values of κ produceunrealistically large initial separations and large relative ve-locities, while very low values of κ lead to very small initialdistances, small relative velocities, and consequently, greatlyreduced total energies, which eventually result in very highconcentration parameters of the final halo. In other words, κ was chosen so as to yield realistic values of both the ini-tial separation and the final concentration parameter. Thedependence of the results on κ is discussed below.After the two halos merge, the resultant halo accretesmatter, so when in turn it merges to form a larger halo, ithas more mass than the sum of the masses of its progenitors.We account for the energy of the accreted matter in a veryapproximate way as follows. Given the masses of the twoprogenitor halos, M and M , and the final mass of the halo, M p , the total accreted mass is ∆ M = M p − ( M + M ) = M p − M f . Numerical simulations of galaxy-sized halos (e. g.Wang et al. 2010) seem to indicate that the accreted materialis distributed in the halo outer region. This is quite likelythe case also in cluster-sized halos, so we can estimate the energy due to the accreted mass by writing U acc = − G ( M + M )∆ MR f , (9)where R f is the virial radius of the halo with mass M f = M + M that formed just after the merger. In evaluating U acc we assume that the accreted mass constitutes a rela-tively small fraction of the final halo, and that equation (9)provides a simplified description of the accretion process.We then have for the total energy of the system priorto merging E total = E tot, + U acc (10)This energy is attributed to the resulting merged halo; weassume no mass is lost in the process. By equating E total with the gravitational and potential energy of the resultinghalo, which depend on its concentration parameter, we candeduce the latter. This process is repeated for each halo inthe tree, until arriving at the bottom - the most massivehalo. At the end of this process we obtain the concentrationparameter for the given mass and for one tree, i.e. one re-alization of the halo history. Generating a large number oftrees gives an estimate of the PDF of c ( M ). In our modeling of IC gas we assume that it constitutes asmall fraction of the total cluster mass, and that it does notsignificantly affect the evolution of the cluster. We assumethat the gas has a polytropic equation of state with an adi-abatic index Γ, such that the gas density and pressure arerelated through: P = P ( ρ/ρ ) Γ (11)The solution of the equation of hydrostatic equilibriumfor a polytropic gas inside a potential well of a DM halo withan NFW profile is (Ostriker, Bode, & Babul 2005): ρ ( x ) = ρ (cid:20) − B n (cid:18) − ln(1 + x ) x (cid:19)(cid:21) n , (12)where n = (Γ − − and B is given by: B = 4 πGρ s r s µm p k B T , (13)and µm p is the mean molecular weight. Thus, B depends onthe concentration parameter through ρ s and r s (see equation(5)). The temperature profile is then T ( x ) = T (cid:20) − B n (cid:18) − ln(1 + x ) x (cid:19)(cid:21) (14)As a boundary condition we assume that the gas pres-sure at the virial radius obeys P gas = f g P DM where f g is thegas mass fraction and P DM = ρ DM σ (Ostriker et al. 2005).We obtain σ , the DM (3D) velocity dispersion, by solvingthe Jeans equation for the NFW potential. We expect thisparticular choice for the boundary condition to have a minorinfluence on the results, as discussed in Section 4.For each merger tree we obtain the concentration pa-rameter and the virial radius of the final halo, as describedabove. Taking the observationally deduced value for the adi-abatic index, Γ = 1 .
2, and imposing the above boundary c (cid:13) , 1–11 n expanded merger-tree description of cluster evolution condition to obtain the constant B fully determines the tem-perature profile. By assuming a specific gas mass fraction wealso obtain the full density profile. Repeating this procedurefor a large number of trees provides an estimate of the PDFsof the various physical parameters as a function of the clus-ter mass and the redshift of observation.In order to derive scaling relations we define the meancluster temperature as the emission-weighted value T ew ≡ R Λ( T ) ρ g T dV R Λ( T ) ρ g dV , (15)where integration is over the cluster volume, and the temper-ature dependence of the cooling function is approximatelyΛ( T ) ∝ √ T . Similarly, X-ray luminosity is approximatedby: L X = Z (cid:18) ρ g µm p (cid:19) Λ( T ) dV . (16)Knowledge of the gas temperature allows also calcula-tion of the Comptonization parameter, y , defined as y = Z (cid:18) k B T e m e c (cid:19) n e σ T dl, (17)where the integral is taken along the line of sight to thecluster, σ T is the Thomson scattering cross section (and theelectron temperature T e is assumed to equal the gas tem-perature). The measured SZ intensity (change) is approx-imately proportional to the integrated Y-parameter, givenby the integral of y over the angle the cluster subtends onthe sky Y = Z yd Ω (18)
The model contains several physical parameters, and twonumerical (code-specific) parameters - the number of treerealizations used to estimate the PDFs and the resolutionmass of the merger tree. We shall discuss the latter parame-ters here and defer the discussion of the physical parametersto subsequent sections.A key objective of our model is to determine the PDFof the concentration parameter and IC gas temperature bygenerating a large number of merger tree realizations N .We have found that taking N = 5 · is sufficient to obtainconvergent results - taking larger N does not change thePDF by more than a fraction of a percent. In what follows,we show histograms of 50 binned values obtained from 5 · merger trees.As noted earlier, a resolution mass needs to be selectedfor each tree. This mass is the smallest building block used inthe tree. By sampling different values of the resolution mass,we find that M res has to be at least 3 orders of magnitudebelow M , the final mass for which the tree is built, whiletaking smaller values of M res does not affect the results: forthe mass range we consider, the mean value of c changesby no more than 1% when the resolution mass is loweredfrom M res = 10 − M to M res = 10 − M . The value of M res determines z fin , as discussed above. log (c ) P ( l og ( c )) tree modellognormal fitN07SR080.2 0.4 0.6 0.8 1 1.2 1.40123456 log (c ) P ( l og ( c )) Figure 1.
Concentration parameter probability distributionfunction (PDF) as predicted by the merger-tree model (bars), thelognormal fit to this PDF (solid line), a distribution extractedfrom the MS by N07 (dot-dashed line) and the prediction ofSR08 (dashed line).
Upper panel: M = 4 × h − M ⊙ for themerger tree model and a corresponding range of M = 10 . − . h − M ⊙ from the MS. Lower panel: M = 6 × h − M ⊙ for the merger tree model and a corresponding range of M =10 . − . h − M ⊙ from the MS. The basic outcome of the model is the concentration pa-rameter of the DM halo at a given redshift of observation.Each merger tree results in a slightly different concentra-tion parameter, which depends on the particular structureof the merger tree. Thus, in the limit of a large number oftree realizations the distribution of formation histories pro-vides a PDF of the concentration parameter. The PDF of c for M = 4 × h − M ⊙ at z = 0 is shown in Figure1 (upper panel). A log-normal distribution provides a rea-sonable fit, with < log c > = 0 .
706 and σ log c = 0 . M = 10 . − . h − M ⊙ seen in the MS: < log c > = 0 .
663 and σ log c = 0 . M = 6 × h − M ⊙ is shown in the lower panel, along with a correspondingdistribution for halos in the MS in the mass range of M = 10 . − . h − M ⊙ . For this mass we obtain < log c > = 0 .
758 and σ log c = 0 .
106 from the log-normalfit, compared with < log c > = 0 .
744 and σ log c = 0 . M and M v for a given halo depends on its concentration pa-rameter. c (cid:13) , 1–11 I. Dvorkin and Y. Rephaeli log (c) P ( l og ( c )) G ≠ log (c) P ( l og ( c )) Bullock et al.Neto et al.−0.2 0 0.2 0.4 0.6 0.8 1 1.2012345 log (c) P ( l og ( c )) q=7q=10q=13 −0.2 0 0.2 0.4 0.6 0.8 1 1.2012345 log (c) P ( l og ( c )) κ =4 κ =5 κ =6 Figure 2.
PDF of the concentration parameter from the merger-tree model for M = 10 h − M ⊙ at z = 0 with different modelparameters, as discussed in the text. Shown in the Left upper panel is the dependence on the mass function through the perturbingfunction G (see equation (4)), the dependence on the initial conditions for c ( M, z ) for the earliest halos in the tree in the
Right upperpanel , the dependence on the major merger parameter q = M > /M < in the Left lower panel , and the dependence on κ = d /d (whichparametrizes the initial distance between halos that are about to merge) in the Right lower panel . For comparison, a semi-analytical calculation adaptedfrom Sadeh and Rephaeli (2008; SR08) is also shown. Thislatter treatment was based on an analytical distribution offormation times and a relation between the formation timeand concentration parameter deduced from numerical simu-lations by Wechsler et al. (2002). The merger tree model pre-dicts slightly lower concentration parameters and a slightlybroader distribution function. It is important to note thatboth treatments result in quite similar PDFs that are alsoconsistent with the results of numerical simulations, despiteof the competely different assumptions made in each of theseapproaches.As mentioned earlier, the uncertainty in the correctform of the mass function is quantified by the perturbingfunction G (see equation (4)). Figure 2 (left upper panel)shows how the concentration parameter changes when themerger tree is computed with and without the perturbingfunction G . As expected, the concentration parameter tendsto be larger in the former case, reflecting the earlier forma-tion time of halos in the MS as compared with the extendedPress-Schechter formalism (Wechsler et al. 2002). This re-sult illustrates the rather strong dependence of the PDF of c on the mass function. This dependence has to be accountedfor when comparing results from observations and numericalsimulations.The initial conditions of the tree are the concentrationparameters of the earliest halos. As indicated earlier, theparticular choice of c ( M, z ) for the earliest halos does notappreciably affect the final value of c , as long as this choice isreasonable. For example, Figure 2 (right upper panel) showsthe probability distributions of c with the initial c ( M, z )taken from the fit in equation (6), and a different fit adaptedfrom the results of N07: c = 5 . (cid:18) M h − M ⊙ (cid:19) − . z f z obs (19)It can be seen that there is no significant change in thedistribution function. The influence of these initial conditionon the results is further discussed at the end of section 3.3.The structure of the tree, and hence the calculation of c , depends somewhat on the chosen ratio for major mergers, q = M > /M < . This dependence is shown in Figure 2 (leftlower panel); the choice of q is guided by several physicalconsiderations. On the one hand, it should not be too small,because this would take into account only nearly equal-massmergers. Hydrodynamical simulations (Wik et al. 2008, Mc-Carthy et al. 2007) show that mergers with mass ratios ashigh as q = 10 still lead to strong disruption of equilib-rium in the inner cluster region, and would thus need to betreated as major merger events in our approach. The dy-namical impact of taking higher values of q has not beenexplicitly explored in hydrodynamical simulations. Accord-ingly, we selected this value to be the highest value of q above which the mergers are approximated as continuousmass accretion.The value of κ , which determines the separation atwhich two halos become bound, also influences the resultsquite appreciably. Figure 2 (right lower panel) shows thatdeviations from the fiducial value of κ = 5 can shift the dis-tribution of c due to changes in cluster initial energies. Asdiscussed earlier, the value κ = 5 was chosen so as to pro-duce realistic distances between clusters that are about tomerge, and is consistent with estimates of relative velocitiesof merging clusters (Lee and Komatsu 2010). c (cid:13) , 1–11 n expanded merger-tree description of cluster evolution log (M/M sun h −1 ) l og ( c ) Figure 3. c − M relation at z = 0 (circles), 0 . c indicated by the thin lines. The mass de-pendence is consistent with the results of the MS. The expectation values of the distribution functions fromthe previous sections provide the concentration parameteraveraged over formation histories. It is obviously importantto follow the redshift evolution of c and its distribution withthe final mass. Figure 3 shows the c − M relation for severalredshifts of observation. The results can be well describedby the scaling relation c ∝ M − α , with strong redshift depen-dence of α , ranging from α = 0 .
075 for z = 0 to α = 0 . z = 1, so that the dependence of c on mass is weakerfor higher redshifts. This likely represents the fact that c de-pends on mass through the formation redshift, and the dif-ference in formation redshifts for different masses observedat z = 0 is larger than for different masses observed at ahigher redshift. This flattening of the mass-concentration re-lation at high redshifts is also seen in numerical simulations(Duffy et al. 2008, Gao et al. 2008), although our predictionsfor α are slightly lower at low redshift and slightly higher athigh redshift than those of Duffy et al. Figure 3 also showsthe results of N07 for halos at z = 0 extracted from the MS(thick line) along with the 1 − σ distribution widths (thinlines). These results of the MS are consistent with the pre-dictions of the merger-tree model, although there seems tobe a systematic offset between the respective results fromthese two very different studies.The dependence on the cluster observation redshift,which is often taken to be c ∼ (1 + z obs ) − γ with γ = 1,is also found to be much weaker and mass-dependent, rang-ing from γ = 0 .
38 for M = 10 h − M ⊙ to γ = 0 .
24 for M = 10 h − M ⊙ . This results in slower redshift evolutionthan found by Duffy et al, but is more consistent with thefindings of Gao et al. for massive halos extracted from theMS, especially for masses around M ∼ M ⊙ for whichour result γ = 0 .
31 coincides with the evolution seen by Gaoet al. (note, however, that these authors use the Einasto pro-file to describe DM halos).In general, the scaling relations deduced from numeri-cal simulations are effectively weighted by the mass function,and, since the latter has a sharp cutoff at about the typi-
T [keV] P ( T ) Figure 4.
IC gas temperature PDFs predicted by the merger-treemodel and best fits to a log-normal distribution for 10 h − M ⊙ (circles, solid line), and for 10 h − M ⊙ (asterisks, dashed line).These distributions exhibit a characteristic sharp cutoff at lowtemperatures and a long high-temperature tail, as expected. cal mass of a galaxy cluster, mainly reflect the structure ofsmaller, galaxy-sized halos at low redshifts. Extrapolationsof the results of such simulations cannot faithfully describethe structure of massive halos at high redshifts, as pointedout by Gao et al. Although we use the results of Bullocket al. as the initial conditions for the merger tree - equa-tion (6) - this choice is justified because our final results arenot sensitive to the exact form of these initial conditions.In addition, the initial halos in the merger tree have smallermasses, in the range explored by Bullock et al.Full investigation of the c − M relations and their red-shift evolution neccesitates the use of numerical simulationstargeted at massive, cluster-sized halos. We plan to continueour study in this direction using the hydrodynamical AMRcode Enzo . Since the temperature of IC gas is used as a mass proxyin cluster surveys, its PDF is of great observational im-portance. We have computed this distribution as outlinedabove. Figure 4 shows the PDFs of the emission-weightedtemperature for cluster masses M = 10 h − M ⊙ and M =10 h − M ⊙ at z = 0.As expected, the PDF exhibits a long high-temperaturetail which corresponds to those clusters that were formedatypically early. At low temperatures, on the other hand,there is a sharp cutoff that corresponds to clusters thatformed close to their observation redshift. A log-normal dis-tribution provides a good approximation to the temperaturePDF below M ∼ · h − M ⊙ , as can be seen in Fig-ure 4. The width of the distribution is σ log c = 0 .
07 for10 h − M ⊙ and σ log c = 0 .
08 for 10 h − M ⊙ . Scaling relations of the gas temperature with cluster red-shift, mass, and X-ray luminosity are commonly used in sta-tistical analyses of the cluster population and in the use of c (cid:13) , 1–11 I. Dvorkin and Y. Rephaeli
13 13.2 13.4 13.6 13.8 14 14.2 14.4 14.6 14.8 15−0.200.20.40.60.81 log10(M/h −1 M sun ) l og10 ( T [ k e V ] / ( E ( z ) ∆ V ( z )) / tree modelArnaud et alKotov et al Figure 5.
Mass-temperature scaling relation at z = 0 for themerger-tree model (stars) and from observations (circles and tri-angles). clusters as cosmological probes. Most useful is the T − M relation which can be determined from the probability distri-bution functions. In Figure 5 we show the emission-weightedtemperature versus mass for clusters at z = 0. The tempera-ture was calculated using equation (15). Blue stars representexpectation values of the PDFs, with errorbars indicatingthe distribution variance. The red circles are measurementsof a sample of clusters from Arnaud, Pointecouteau & Pratt(2005), and the black triangles are measurements of anothersample by Kotov & Vikhlinin (2005), where redshift correct-ing factors have been included for both samples. The mergertree results are best-fit with the relation T ∝ M . , which isvery close to the theoretical relation obtained for an isother-mal sphere, T ∝ M / . It can be seen that the results andthe expected scatter are consistent with observations. Notethough the different definitions of mass ( M in Arnaud etal., M in Kotov & Vikhlinin) and temperature (spectraltemperature in both Arnaud et al. and Kotov & Vikhlinin).The observational results suggest that the variance ofthe temperature PDF can be seen to represent the amountof scatter that is expected in observed clusters due to theirdifferent formation history. Note that the error in the mea-sured temperatures is small compared to the scatter, whichis slightly larger than the predicted intrinsic scatter, as ex-pected, since it has additional contributions. For example,not all clusters are fully relaxed and spherical, etc. We findthat the temperature scales as a power-law in mass at allredshifts, T ∝ M ζ ; however ζ varies somewhat with red-shift, from 0 . z = 0 to 0 .
63 for z = 2, which resultsin slower evolution compared to the simple scaling ζ = 2 / minimal possible temperature of a givenmass - the low-temperature endpoint of our PDF - scales as T min ∝ M . − . in our model for all redshifts, in muchbetter agreement with the standard value. Indeed, the stan-dard treatment assumes that the halo is observed immedi-ately after it had formed, which is precisely the situationdescribed by the low-temperature end of the PDF. The ex-pectation value, however, is affected by the width of thePDF, which also depends on mass.The predicted redshift dependence of T is another keyrelation whose knowledge is important as it reflects on clus- log (T [keV]) l og ( L [ e r g / s ] ) Figure 6.
X-ray luminosity vs. emission-weighted temperature(circles with error bars), where the error bars correspond to thedisrtibution width in the logarithm of the luminosity and tem-perature, respectively. A constant gas mass fraction f g = 0 . ter evolution, and its approximate analytic form is neededin comparisons with results of cluster X-ray and SZ surveys.The basic redshift scaling of the temperature is contained inthe relation T ∝ [ E ( z )∆ V ( z )] λ , where λ varies somewhatwith mass, from λ = 0 . M = 10 h − M ⊙ to λ = 0 . M = 10 h − M ⊙ . Thus, T ( z ) is less steep than in thestandard relation (1), where λ = 1 /
3. In addition, the slopeof this scaling relation differs with mass, hinting that thetemperature might not be a separable function in terms ofmass and redshift. The dependence of these results on themodel parameters is discussed in Section 3.8.The luminosity-temperature relation is an importantprobe of the IC gas. In the framework of the presented ap-proach it can be used to test the validity of the simple poly-tropic model. Figure 6 shows the luminosity-temperaturerelation obtained from the merger-tree model, as well as X-ray measurements of a sample of clusters by Pratt et al.(2009). There is reasonable agreement with the data in thehigh-temperature end, with the distribution width approx-imately corresponding to the scatter in the measured val-ues, but the model clearly overpredicts the luminosity oflow-temperature clusters. One reason for this could be non-constant gas mass fraction which, as hinted by observations,is lower in low-mass systems. The dependence of the gasmass fraction on mass and redshift could also be related toadditional physical processes in the IC gas, such as radiativecooling and feedback from supernovae and AGN.
Having determined the IC gas temperature and density pro-files (as outlined above), we can now compute another keyobservable - the integrated Comptonization parameter. Todo so, we also need to specify the gas mass fraction, whichis taken to be f g = 0 . Y ; it exhibits the same general features as the tem-perature distribution, a sharp cutoff at low Y , and a longexponential tail at high values, largely due to clusters that c (cid:13) , 1–11 n expanded merger-tree description of cluster evolution −5 Y ⋅ d A (z) [Mpc /h ] P ( Y ) Figure 7.
PDF of the integrated Comptonization parameter fromthe merger-tree model for M = 10 h − M ⊙ at z = 0 . formed uncharacteristically early. We note that a log-normaldistribution is a poor fit to the outcome of our model.As in the case of IC gas temperature, scaling relationsof the mean values of Y can be computed. The scaling withmass is Y d A ( z ) ∝ M δ where d A ( z ) is the angular diam-eter distance and δ varies with redshift, from δ = 1 .
61 for z = 0 .
01 to δ = 1 .
64 for z = 2. This scaling is close to thestandard result δ = 5 /
3. Similarly, the Comptonization pa-rameter scales with redshift as
Y d A ( z ) ∝ [ E ( z )∆ V ( z )] ε with ε = 0 .
26 for M = 10 h − M ⊙ and ε = 0 .
31 for M = 10 h − M ⊙ . The evolution with redshift is slower thanin the standard description where ε = 1 / The PDFs of cluster observables presented above provide atheoretical basis for comparisons with results of cluster sur-veys. As an example we consider here the predicted temper-ature number counts, which is one of the statistical clusterfunctions that can be used to determine cosmological pa-rameters.The temperature function, that is the cumulative num-ber density of clusters above a certain temperature at a givenredshift (interval) is computed from the following expression n ( T i ) = Z M high M low B ( T i | M, z ) dn ( M, z ) dM dM, (20)where dn ( M, z ) /dM is the mass function, namely the num-ber of halos per unit comoving volume per unit mass. Theselection function B ( T i | M, z ) is usually defined as B = 1 if T ( M, z ) > T i and B = 0 otherwise, where T ( M, z ) is foundaccording to the standard scaling relations (with a sharpcutoff).However, in accord with our treatment here, there is noone-to-one correspondence between temperature and mass,so we need to incorporate the PDF of the temperature inthe calculation of the number counts by using the followingselection function B ( T i | M, z ) = Z ∞ T i P ( T | M, z ) dT . (21)The temperature PDF is described by a log-normal dis-tribution with expectation value and variance taken from −11 −10 −9 −8 −7 −6 −5 −4 T [keV] n [ h M p c − ] sharp cutoffpdf Figure 8.
Temperature number counts: the number density ofclusters above a certain temperature at z = 0 . best fits to the results of our model. The temperature func-tions calculated with our more realistic temperature PDFand that with the standard relation (between temperatureand mass) are shown in Figure 8. In the standard calcula-tion we chose T ( M, z ) to equal the expectation value of therespective PDF. The calculations were performed over themass range 10 − h − M ⊙ using the Sheth-Tormen massfunction.The two calculations coincide for low temperatures, butfor high temperatures our improved treatment yields appre-ciably higher number counts. The reason for this is that ourmore exact treatment takes into account the long tails of thedistribution functions. Thus, low-mass clusters with meantemperatures below T , that do not contribute to n ( T )when the standard scaling is used, can have a significantoverall contribution when the temperature PDF is used. Asdiscussed earlier, this is due to the non-zero probability thatthe formation redshifts of these clusters, and hence also theirtemperatures, were higher than the mean values.As we have mentioned earlier, the log-normal distribu-tion is a mediocre fit to the PDFs of high-mass clusters, anda better understanding of their shapes is required in orderto fully assess their impact on temperature number counts.The above calculation demonstrates the importance of tak-ing temperature PDFs into account in the analysis of clustersurveys. In the previous sections we have shown that our method forthe determination of the PDFs of the various cluster physicalparameters provides a relatively simpler procedure to im-plement than hydrodynamical simulations. The procedureinvolves specifying several free parameters: q - the maximalmajor merger ratio, κ - the parameter that determines theinitial distance between clusters, and the adiabatic index ofthe gas, Γ. We should also add to this list the parametersof the initial c ( M, z ) chosen for the smallest halos in thetree (see equation (6)). These parameters were found not toinfluence the results considerably when chosen reasonably,in accordance with observational results and N-body simu-lations; see the discussion at the end of Section 3.3. c (cid:13) , 1–11 I. Dvorkin and Y. Rephaeli log (c) P ( l og ( c )) Figure 9.
Concentration parameter PDF from the merger-treemodel for M = 10 h − M ⊙ at z = 0 and model parameters inthe range q = 7 −
13 and κ = 4 −
6. The thick curve shows thePDF calculated with the fiducial values.
Since we are mainly interested in the global propertiesof the cluster, such as the emission weighted temperatureand the integrated Comptonization parameter, our resultshave a very weak dependence on a particular choice for theIC gas profile. We have repeated our calculations using the β -profile for the gas: ρ ( r ) = ρ g " (cid:18) rr c (cid:19) − β/ (22)with β = 2 /
3. The gas core radius is given by r c = r s /η where r s is the DM scale radius, and a typical value is η = 2(Ricker and Sarazin, 2001). We have integrated the equationof hydrostatic equilibrium to obtain the gas temperatureprofile, setting the pressure to zero at infinity. With thisprofile, the temperature PDFs change by just a few percentrelative to the polytropic model. We repeated the calculationusing also the β -profile with a different boundary condition,namely setting the gas temperature at the virial radius tothe temperature of the IGM, typically 10 − K. This toohad only a minor impact on the results.In order to estimate the robustness of our model wecompute the errors on the PDFs that result from small de-viations from the fiducial values of the main model param-eters. Figure 9 shows the PDF calculated with parametersin the range q = 7 − , κ = 4 − q = 7 − , κ = 4 − , Γ = 1 . − . < T > = 9 . +0 . − . keV for M = 10 h − M ⊙ at z = 0, while the width ofthe distribution is p < ( T − < T > ) > = 1 . +0 . − . keV. Thislikely is a conservative estimate for the range of the temper-ature uncertainty.The evolution of the variance of the PDF with massand redshift and its uncertainty can also be estimated. Forinstance, if the parameters are varied in the same range q =7 − , κ = 4 − , Γ = 1 . − .
3, and for the same mass of M = 10 h − M ⊙ but for observation redshift of z = 0 . < T > = 12 . +1 . − . keV and p < ( T − < T > ) > = 1 . ± . M = 10 h − M ⊙ at redshift z = 0 are: < T > = 2 . +0 . − . keV and p < ( T − < T > ) > = 0 . +0 . − . keV. The relativeuncertainties in all these cases are similar.Finally, we can estimate the robustness of our resultsfor the evolution of the scaling relations. As an example, wehave checked how the scaling of the temperature with mass( T ∝ M ζ ) and redshift ( T ∝ [ E ( z )∆ V ( z )] λ ) changes whenthe model parameters are varied in the range q = 7 − , κ =4 − , Γ = 1 . − .
3. It turns out that ζ and λ change by nomore than 1% and 6%, respectively, relative to the valuesobtained in Section 3.5. We have presented an expanded merger-tree treatment forthe evolution of galaxy clusters that supplements the sta-tistical description of the dynamical evolution of DM haloswith basic physical considerations that enable us to describealso the properties of IC gas. It should be stressed again thatour approach is statistical by construction and is not meantto provide a prescription for determining the structure of in-dividual halos, but rather to serve as a tool for studying theproperties of a population of clusters. While our treatmentis essentially adiabatic, we have adopted an observationally-deduced value of the polytropic index. By doing so we partlycompensate for the fact that gas cooling is not explicitlytaken into account. Additional justification for the validityof our approach is the fact that we are interested here only instatistical properties of the cluster population, rather thanin detailed spatial profiles of the gas density and tempera-ture in individual (such as cooling-core) clusters.We also assumed that the DM mass profile is not af-fected by the IC gas. Although this approximation is oftenmade in studies of the statistical properties of a popula-tion of clusters (e.g. Bode, Ostriker & Vikhlinin 2009), itis likely to be inaccurate when radiative cooling is impor-tant, or when there is energy exchange between the DM andthe gas components, for example during mergers. Numericalsimulations (Duffy et al. 2010) show that there is a devia-tion of at most 15% in the concentration parameter of groupsand clusters when baryonic physics is included, relative tothe DM only case. The impact of IC gas cooling on the DMdensity profile is often described by adiabatic contractionmodels (e.g. Gnedin et al. 2004). However, the assumptionmade in these models that the baryons initially trace theDM distribution is violated during hierarchical build-up ofhalos. Indeed, Duffy et al. (2010) found that results of thesimulations were not well described by adiabatic contractionmodels beyond 0 . R vir . A natural extension of our modelwould be to incorporate IC gas in the halos that constitutethe merger tree and to follow the joint evolution of bothcomponents.We calculated the PDFs of the cluster concentrationparameter, its IC gas temperature, and integrated Comp-tonization parameter for different masses and redshifts ofobservation. Our deduced PDF of the concentration pa-rameter is well fit with a log-normal distribution, in ac-cord with results from N-body simulations. The tempera- c (cid:13) , 1–11 n expanded merger-tree description of cluster evolution ture PDF for masses below M ∼ · h − M ⊙ can alsobe described with a log-normal distribution. Our deducedmass-observable scaling relations are close to the standardrelations but contain some corrections - notably the evolu-tion of IC gas temperature with redshift is slower than in thesimple model. The results suggest that the gas temperatureis not a separable function of mass and redshift. We show apossible application of our results to the analysis of clustersurveys by calculating IC gas temperature number counts,taking into account the effect of cluster formation history.The probability density functions of the various ob-servables can have important effects on the error estima-tion in the analysis of cluster X-ray and SZ surveys. Asshown by Lima and Hu (2005), large uncertainties in theobservable-mass distributions may substantially degrade theconstraints on cosmological parameters from cluster surveys.The physically-based estimates of the PDFs of the observ-ables considered here provide a tangible basis to begin ad-dressing this aspect.Among the other related applications of the approachpresented here a particularly timely one is the calculation ofthe SZ power spectrum, which will be mapped by the Planck satellite and several ground-based SZ projects. Comparisonsof results from our merger-tree approach and those fromsimulations and semi-analytical treatments (e.g., see Sadeh,Rephaeli & Silk 2007 and references therein) will yield im-portant insight that will help gauging the relative merits anddisadvantages of these very different approaches.
AKNOWLEDGEMENTS
The authors wish to thank the GALFORM team for makingthe code publicly available. This research was supported bya US-Israel Binational Science Foundation grant 2008452.
REFERENCES
Arnaud, M., Pointecouteau, E., Pratt, G. W. 2005, A&A, 441,893Broadhurst, T. J., et al. 2008, ApJ, 695, L9Bullock, J. S., et al. 2001, MNRAS, 321, 559Cole, S., Lacey, C. 1996, MNRAS, 281, 716Cole, S., et al. 2000, MNRAS, 319, 168Cunha, C. E., Evrard, A. E., 2009, Phys. Rev. D, 81, 083509Duffy, A. R., et al. 2008, MNRAS, 390, L64Duffy, A. R., et al. 2010, MNRAS, 405, 2161Gao, L., et al. 2008, MNRAS, 387, 536Gnedin, O. Y., et al. 2004, ApJ, 616, 16Jenkins, A., et al. 2001, MNRAS, 321, 372Kotov, O., Vikhlinin, A. 2005, ApJ, 633, 781Lee, J., Komatsu, E. 2010, ApJ, 718, 60Lima, M., Hu, W., 2005, Phys. Rev. D, 72, 043006 Lokas, E. L., Mamon, G. A., 2001, MNRAS, 321, 155McCarthy, I. G., et al. 2007, MNRAS, 376, 497Navarro, J. F., Frenk, C. A., White, S. D. M., 1995, MNRAS,275, 720Navarro, J. F., Frenk, C. A., White, S. D. M., 1996, ApJ, 462,563Neto, A. F., et al. 2007, MNRAS, 381, 1450Ostriker, J. P., Bode, P., Babul, A., 2005, ApJ, 634, 964Parkinson, H., Cole, S., Helly, J. 2008, MNRAS, 383, 557Press, W. H., Schechter, P. 1974, ApJ, 187, 425Pratt, G. W., et al. 2009, A&A, 498, 361 Randall, S. W., Sarazin, C. L. 2002, ApJ, 577, 579Ricker, P. M., Sarazin, C. L. 2001, ApJ, 561, 621Sadeh, S., Rephaeli, Y. 2008, MNRAS, 388, 1759Sadeh, S., Rephaeli, Y., Silk, J. 2007, MNRAS, 380, 637Sarazin, C. L. 2002, in Merging Processes in Clusters of Galaxies,ed. L. Feretti, I. M. Gioia, G. Giovannini (Dordrecht: Kluwer),1Sheth, R. K., Tormen, G. 2002, MNRAS, 329, 61Springel, V., et al. 2005, Nature, 435, 629Tinker, J., et al. 2008, ApJ, 688, 709Vikhlinin, A., et al. 2006, ApJ, 640, 691Vikhlinin, A., et al. 2009a, ApJ, 692, 1033Vikhlinin, A., et al. 2009b, ApJ, 692, 1060Voit, G. M., Donahue, M., 1998, ApJ, 500, L111Wang, J. et al. 2010, preprint (astro-ph/1008.5114)Wechsler, R. H., et al. 2002, ApJ, 568, 52Wick, D. R., et al. 2008, ApJ, 680, 17Zitrin, A. et al. 2009, ApJ, 707, L102
This paper has been typeset from a TEX/ L A TEX file preparedby the author. c (cid:13)000