THE PROPERTIES OF X-RAY COLD FRONTS IN A STATISTICAL SAMPLE OF SIMULATED GALAXY CLUSTERS (2024)

For the purpose of understanding the impact of additional baryonic physics on the cold front incidence, we have analyzed two additional simulations. Each of the two simulations has identical simulation volumes, initial conditions, and cosmology. The simulated volume is a 128 h⁻¹ comoving Mpc cube, with 256³ root grid zones and 256³ dark matter particles. Each run is allowed to refine to a maximum of five AMR levels, and the refinement criteria are overdensity of 8.0 in either the dark matter or the gas (with respect to the mean on the parent level). This gives the run a peak spatial resolution of 15.6 h⁻¹ kpc. The dark matter particle mass is 9.0 × 10⁹ h⁻¹ M_☉ and the mean baryon mass resolution is 1.4 × 10⁹ h⁻¹ M_☉, a factor of eight improvement over the SFLC run used for our main analysis. For these runs, we use a cosmological model with Ω_m = 0.27, Ω_b = 0.044, Ω_Λ = 0.73, h = 0.7 (in units of 100 km s⁻¹ Mpc⁻¹), σ₈ = 0.9, and using a power spectrum with a spectral index of n_s = 0.97. These values are changes from our previous run, to match Wilkinson Microwave Anisotropy Probe Data Release 5 cosmological parameters. Therefore, there are multiple changes to these simulations from the SFLC, including the mass and spatial resolution, baryonic physics, and cosmological parameters. The most relevant changes to the cosmology are a slight reduction in the value of Ω_m and a slight increase in the value of Ω_Λ. However, in this work, we endeavor to merely compare our two runs where everything is identical except baryonic physics to examine the impact.

The two additional simulations presented here are identical in terms of the initial conditions and simulation properties, except that one has only adiabatic and gravitational physics, and the other includes a prescription for radiative cooling of the gas using non-equilibrium cooling and chemistry for H and He, Cloudy (Ferland etal. 1998) for the metal cooling, and star formation including thermal and metal feedback from supernovae as detailed in Cen & Ostriker (1999). The cooling plus star formation and feedback run includes a spatially uniform but time-varying ultraviolet (UV) radiation background (Haardt & Madau 1996). We turn on the UV background at a redshift z = 7. For the metal line cooling, we interpolate from a grid of data made with Cloudy. For redshifts z > 7, we use a similar grid of heating and cooling data, without the influence of the UV radiation background, assuming collisional ionization only (for more details, see Smith etal. 2008, 2010). Additionally, we have highly time-resolved data outputs for these calculations, allowing us to examine the time history and origin of cold fronts in the simulated clusters in detail. In this paper, we examine only the general incidence of cold fronts in these simulations, and leave a study of the detailed properties of the cold fronts, as well as a full convergence study, to future work.

For our simulation with purely adiabatic physics (CC-Adia), and the identical run with additional baryonic physics (CC-OTH), the data outputs and analysis are identical. We perform the same analyses as done for the SFLC simulation.

The first result of the cold front identification is that in the adiabatic simulation, for the sample of all clusters with M ⩾ 10¹⁴M_☉ out to z = 1, 12% of images have 10 or more cold front pixels, 5% have 20 or more, and roughly 1% have 50 or more. At z = 0 in our adiabatic simulation, we have 116 images with more than 50 cold front pixels. In these images, there are clear, unambiguous cold fronts spanning large areas of the cluster projection (see Figure2). Figure3 shows the result of counting clusters above some minimum pixel count identified as cold fronts. This is equivalent to measuring the comoving total length of cold front features in each image, as each pixel has a constant comoving size in all of these images. In the figure, we plot the trend as a function of redshift for clusters with at least 10 pixels (L_cf>156 h⁻¹ kpc total) identified as a cold front, 20 pixels (L_cf>312 h⁻¹ kpc), and 50 pixels (L_cf>780 h⁻¹ kpc). These numbers are arbitrary, and designed purely to show the general redshift trend. It is clear from this plot that there is only a weak trend in the incidence of cold fronts as a function of redshift out to z= 1.

The choice of a threshold in cold front pixels or comoving length of features is fairly arbitrary, but identification of one pixel with cold front criteria is not particularly interesting. First, it is quite possible to have spurious cold pixels due to line-of-sight superpositions. Additionally, cold fronts are extended features, and so in these images at least several pixels need to be associated in order for a true identification to take place.

Another interesting statistic is the incidence of cold fronts in simulated clusters as a function of mass. For this case, we take clusters with at least 10 identified cold front pixels, which as described earlier is equivalent to the combined length of cold fronts L_cf>156 h⁻¹ kpc. We combine the result for all clusters from 1.0 ⩾ z ⩾ 0.0 in Figure4. What we see here is that the likelihood of finding extended cold fronts in clusters increases with mass. There are several possible interpretations for this result. The first, which is that this is a real physical trend, is that more massive clusters form in more dense local environments, leading to more activity from mergers and accretion, increasing the number of cold fronts generated. A second possibility is that because smaller clusters are smaller in projection, and we have a fixed comoving pixel scale, they are less likely to have extended cold front features purely because of resolution effects. This effect is compounded by the larger cross section for mergers for the larger clusters, generating a higher likelihood for a merger of a given size in that projected volume at any time. It is legitimate to ask how the sample can have only a weak trend in cold front incidence with redshift, yet have a strong mass trend. Given the hard lower mass limit at each redshift, the distribution of clusters as a function of mass in each redshift bin is different. In other words, the clusters selected at each redshift with an identical mass cutoff do not represent the same part of the halo mass function at each redshift. This effect may contribute to the result shown here. We expect that if cold fronts are related to merger activity (as is indicated by simulation studies) that their incidence should trend with redshift, as large-scale structure grows, and the merger rates change as a function of time.

In our above calculation of the fraction of cold front clusters as a function of cluster mass, we include clusters from all redshifts. Does this trend with mass vary as a function of redshift? Figure5 shows the fraction of cold front clusters as a function of mass for four redshift bins in the simulation. The result is that the trend with mass is nearly identical for all redshifts, albeit with slight shifting in the upper mass bins due to small numbers. Effectively, clusters at all redshifts show an identical trend of increasing fraction of clusters with cold fronts with cluster mass.

Using a mass or flux-selected sample of galaxy clusters to study trends as a function of redshift is an inherently flawed method, as has been pointed out by many investigators. The most well-known bias in flux-limited samples is the Malmquist bias, where your survey progressively selects brighter objects at greater distances. For galaxy clusters this has the effect of choosing the higher mass objects at higher redshift. Our sample has a bias working in the opposite direction, which is that we are selecting a mass-limited sample from all redshifts. Effectively we are sampling objects in different parts of the halo mass function at each redshift, not objects that at high redshift would be the precursors to the objects at low redshift in our sample. The result of such a selection is that true evolutionary effects as a function of redshift can easily be masked. One solution is to use a method similar to Hart etal. (2009), selecting objects along an evolutionary "road," such that objects selected at low redshift are, in a statistical sense, the children of the objects at higher z.

To understand our result of the incidence of cold fronts having only weak trending with redshift, we sort our cluster sample in two additional ways. The first method we choose here is to use a Press–Schechter (Press & Schechter 1974) mass function analysis to select a lower mass cutoff at each redshift corresponding to halos of the same rarity in large-scale structure. This threshold choice effectively follows the same cluster population as a function of redshift, allowing us to see evolutionary trends in the incidence of cold fronts.

Our selection is done by calculating the value of σ_M, the rms density fluctuation, as a function of cluster mass and redshift. We then choose a lower mass limit for our sample such that σ_M is a constant with redshift. This choice is motivated by the analytic expressions for the cluster mass function, both derived from linear theory and fit empirically to numerical simulations. The derivations below are taken from the cited references as well as Hallman etal. (2007). The comoving number density of clusters as expressed by Jenkins etal. (2001) is

where σ_M is the rms density fluctuation, computed on mass scale M from the z = 0 linear power spectrum (Eisenstein & Hu 1999), ρ₀ is the mean matter density of the universe, defined as ρ₀ ≡ Ω_mρ_c (with ρ_c being the cosmological critical density, defined as ρ_c ≡ 3H²₀/8πG), and D(z) is the linear growth function, given by this fitting function:

(Carroll etal. 1992), with Ω_m(z) and Ω_Λ(z) defined in the typical way in a flat universe with cosmological constant. For the details of the calculation, see Hallman etal. (2007).

The rms amplitude of the density fluctuations as a function of mass and redshift, smoothed by a spherically symmetric window function with comoving radius R, can be computed from the matter power spectrum using

where is the Fourier transform of the real-space top hat smoothing function:

We choose the value of σ_M for our lower limit from the lower mass limit at z = 1, which is M_limit = 10¹⁴M_☉. Simply put, we use a constant σ_M limit, and take all clusters in the simulation above that lower mass limit at each redshift. The result is shown in Figure6.

The result clearly changes the sense of the redshift trend from the mass-limited sample. Now we see an increase in the fraction of clusters with cold fronts as a function of time. This selection should allow us to isolate evolutionary effects, as we are sampling the same part of the mass function at each redshift. The fraction generally increases as we move forward in time, which is not particularly surprising, since the result of our sample selection is that the lower mass limit goes up with time. As we have shown, the fraction of clusters with cold fronts is a strong function of mass. What is puzzling is the dip in the fraction at z = 0.2. It is not entirely clear how to explain the appearance of this feature, though it is possible that a varying lower mass limit may create arbitrary features.

The second additional selection mimics an observational selection by a flux limit. Since this simulation has only adiabatic physics, and as such does not predict the X-ray luminosity accurately, we do not use the total X-ray emission calculated from the physical properties of the simulated clusters. Instead, we take a simple approach, using observationally determined X-ray scaling relations to create a mass limit as a function of redshift for our sample.

We take the X-ray scaling relation fit by Mantz etal. (2009), a simple power law

where

and

E(z) is defined in the usual way in a flat ΛCDM universe,

and β^lm₀ and β^lm₁ are the fitted parameters. Mantz etal. (2009) find those parameters to be β^lm₀ = 0.82 ± 0.11, and β^lm₁ = 1.29 ± 0.07 by fitting to the X-ray luminosity function. We use the flux limit from the ROSAT-ESO Flux-Limited X-ray sample (REFLEX; Böhringer etal. 2001). In this case, we are not interested in any particular flux limit, we simply want to show the effect of sample selection on the cold front statistics. The REFLEX flux limit is 2.0 × 10⁻¹² erg s^-1 cm⁻². At each redshift from the simulation, we calculate the X-ray luminosity appropriate to the flux limit, then calculate the cluster mass associated to that luminosity using the scaling relation from Mantz etal. (2009). For the z = 0 clusters, we use z = 0.05 to calculate the value of D_L and the flux. We compare the lower mass limits as a function of redshift in the σ_M- and flux-limited cases in Figure7.

The resulting statistics of cold front incidence are shown in Figure 8. Note that this plot only shows fractions of clusters with cold fronts from z = 0.5 to z = 0, as at this flux limit, we have no clusters in our simulation above the appropriate mass at higher redshifts. In any case, it is clear that a flux limited selection results in a very different picture of the incidence of cold fronts when compared either to the mass-limited or σ_M-limited cases. The flux limit, of course, preferentially selects higher mass objects at higher z. As we have seen in our analysis of the incidence of cold fronts as a function of mass, it is a strong function of mass, so this result is not surprising. In a flux limited sample, we should expect the highest fraction of cold front clusters at high redshifts (since we select the highest mass clusters there), and it should decrease monotonically as we move to lower redshift. In other words, given that the flux limit prefers high mass clusters at high redshift, and we know that there are more cold fronts in high mass clusters at all redshifts from our analysis, we expect the highest fraction of cold front clusters at high z in a flux-limited sample. In this selection, the total fraction of clusters with L_cf>156 h⁻¹ kpc in the full sample is ∼15%, lower than typically found in observed samples.

Here, we study the properties of the identified simulated cold fronts—in particular the jump in temperature and surface brightness across the features. We expect from previous work looking at the shock population (e.g., Ryu etal. 2003; Skillman etal. 2008) that the distribution of number of pixels as a function of ratio of T_sl and 0.3–8 keV X-ray surface brightness should be roughly a power law. As with shocks, in the formation of large-scale structure, we expect more small-ratio jump features than large-ratio ones, meaning that low velocity and temperature contrasts tend to be more plentiful than high contrasts. With our simple cold front identifier, we indeed find that the temperature jump frequency for all the clusters from z = 1 to z = 0 follow a rough power law with a break, as shown in Figure 9. The jumps in surface brightness follow a less simple distribution, flat at low jump ratios and rolling off at large values as shown in Figure 10. A population of pressure continuous features in the simulations should produce a distribution of temperature and density jumps with similar shapes. Though here we show the temperature and surface brightness jumps, under the assumption that the soft X-ray emissivity is only weakly dependent on temperature, the surface brightness and density jumps have very similar shapes (though density is proportional to roughly the square root of X-ray emissivity). It is unclear why the shape of the temperature jump distribution does not qualitatively match the shape of the surface brightness jumps.

Additionally, we can check the number of identified cold front pixels per simulated cluster image. We have already described this measurement as equivalent to the comoving total length (L_cf) of cold front features in the full image, though we make no effort here to determine whether the features are contiguous. For now, this calculation gives us a sense of the total extent of cold fronts in all the simulated cluster images. Figure 11 shows the result of this calculation, a probability distribution function for the total comoving length of cold fronts in the synthetic cluster maps at three redshifts. Interestingly, the distributions at three redshifts (1.0, 0.5, and 0.0) overlap quite well. The normalization is adjusted by the total number of objects with cold fronts identified, but the relative number of objects with a specific measured length of cold fronts is more or less constant across redshift. The high L_cf end of the distribution grows longer at later times, presumably because the objects progressively get larger, and therefore have larger values of r₅₀₀, thus more area within which to find cold fronts.

Since cold front jumps are observed to be continuous in pressure, we check, in projection, whether these features are consistent with observations. If we take the projected surface brightness jump and use it to estimate (roughly) the projected density jump, and compare that to the corresponding projected temperature jump, we can check to see if pressure continuity is met. The X-ray surface brightness in the 0.3–8.0 keV band is only very weakly dependent on temperature, and is effectively representative of the projected value of n²_e. For one of our simulated cold front clusters, we show the properties of a cold front in Figure 12. We use the square root of the X-ray surface brightness jump as a proxy for the projected density, and multiply it by the value of T_sl to form a projected pressure map. Figure 12 shows the obvious difference between the hot feature (a shock) and the cold feature behind it (a cold front). In the shock case, the jump in temperature corresponds to a pressure jump, and in the cold front case, the temperature jump is effectively continuous in pressure. Note also the opposite direction of the jump in surface brightness and temperature in the cold front. These jumps are typical of the shocks and cold fronts in the simulated clusters.

The cold fronts identified in the images of the simulated clusters have properties consistent with those observed in real galaxy clusters. The jumps in temperature are typically factors of 1.5–3, with a corresponding inverse change in the X-ray surface brightness, such that the features are continuous in pressure. Cold fronts in the simulated images are typically associated with shocks, as is clear in Figure 12. While we have not specifically looked at the observability of these features in the X-ray, it is clear from our data that the surface brightness in the hot shocked region is significantly lower than that in the cluster core. It is possible that many of these shocks would not be detected in X-rays for a typical exposure time, though the cold front would be more easily detectable. In future work, we will explore the statistics of cold fronts given the X-ray observability using more sophisticated synthetic observations including instrument response and X-ray backgrounds.

It is also worth exploring the recent history of the simulated clusters with the most obvious cold front features. Simulations show that cold fronts should be a good diagnostic for recent merger activity (e.g., Mathis etal. 2005; ZuHone & Markevitch 2009), both in the case of sloshing type cold fronts and the more extended merger type cold fronts. A multipanel image of some of the most extended cold fronts in our simulated sample at z = 0 are shown in Figure 13. The left column of this image shows the projected spectroscopic-like temperature (with identified cold fronts shown as white contours), and the right column is the X-ray surface brightness. In these cases it is very clear that the cold front finder has done a sufficiently accurate job of characterizing the cold front regions. We note that these cold fronts are all of the merger type, and not the sloshing type. We do not see the sloshing type cold fronts in these adiabatic clusters, and do not expect them given the lack of radiative cooling and steep central entropy gradients.

For three of these obvious cold front clusters, we have taken snapshots of their recent history to illustrate the process by which these extended cold fronts are formed in the cosmological simulations. Recall that for this large set of simulated clusters, the physics in the simulation is purely adiabatic, meaning that there is no radiative cooling. Therefore any regions where cold and hot gas are in close proximity result either from shocks or adiabatic effects. The histories of these three clusters are shown in Figures14, 15, and 16. In each case, we see that the formation of the cold fronts is a result of a recent merger.

In Figures14 and 16, we show two classic roughly plane-of-the-sky mergers, where in the earlier time steps, the halos of the merger are clearly approaching each other, as evidenced by the shock features between them. At the later epoch, we see the subclusters after core passage, when the cold fronts have developed. We also show the dark matter contours, and there appears to be some mild separation of the dark matter centroids from the gas peaks, though we have not quantified this result in this work. In Figure 15, we see an array of cold fronts throughout the history of the cluster, but in this case the dynamics are more complicated. In the far left panels, we see a situation where there has been an initial core passage of two subclusters. The second panels show the approach of the two objects as they fall back together, with an obvious shock between them, while at the same time a third subcluster approaches from the lower left. In the third panel, the initial merger has had a second core passage, while the extra merging subcluster is now driving a shock as it wraps around the subcluster in the upper left. At the final stage, there are two very obvious cold fronts, but deducing the merging history from the final snapshot could be quite challenging.

Cold fronts in the simulated clusters result from a variety of dynamical scenarios, ranging from simple binary mergers to multiple simultaneous mergers. The appearance of cold fronts and shocks is however a generic result of all kinds of mergers in these systems.

Recent studies of well-known cold front clusters identify a possible relationship between the presence of X-ray cold fronts and the merger state of the clusters (Owers etal. 2009). There is also some indication that in some fraction of clusters the presence of the cold front may be the only indicator of a recent merger (Ghizzardi etal. 2010). To test these ideas, we have undertaken a study of the X-ray morphology of the clusters in our simulated samples via the measurement of power ratios and centroid shifts of the X-ray images. These quantitative measures are an indicator of whether the X-ray morphologies are obviously disturbed or regular, and are a proxy for deciding whether a cluster is dynamically active or relaxed (Jeltema etal. 2008; Ventimiglia etal. 2008).

In short, the power ratios are the multipole moments of the X-ray surface brightness map from a circular aperture centered on the cluster X-ray centroid (see Buote & Tsai 1995, 1996; Jeltema etal. 2005). As described elsewhere, this method is related to the multipole expansion of the two-dimensional gravitational potential. This set of equations is reproduced from Jeltema etal. (2008). The multipole expansion of the two-dimensional gravitational potential is

and the moments a_m and b_m are

where and Σ is the surface mass density. In the case of X-ray studies, X-ray surface brightness replaces surface mass density in the calculation of the power ratios.

The powers are formed by integrating the magnitude of Ψ_m, the mth term in the multipole expansion of the potential given in Equation (12), over a circle of radius R,

Ignoring factors of 2G, this gives

The higher multipole moments are sensitive to the amount of substructure in the cluster. In most cases, as in this work, we normalize the multipole moments P₂, P₃ and P₄ to P₀ to factor out the overall brightness of the cluster. P₂, the second multipole moment, is sensitive to deviations from circularity, P₃ to deviations from mirror symmetry, and P₄ is similar to P₂, but sensitive to smaller scale asymmetry. Each of the moments has a component of radial decline in the cluster profile which contributes to its value, such that substructure at large radius leads to larger power ratios. Perfectly round clusters (irrespective of the shape of the radial profile) would give all zero power ratios. While each of the power ratios is sensitive to slightly different asymmetries, it is important to note that they are strongly correlated with one another. That means that clusters that are disturbed via merging activity typically show deviations from circularity, as well as from mirror symmetry and small-scale substructure simultaneously. In this work, we use the values of the power ratios inside a radius of r₅₀₀, a radius which is typically observable in X-rays for a large set of galaxy clusters. We note that the value of the power ratio selected is somewhat dependent on the outer radius used, and that merging clusters with disturbances on larger scales than the outer radius may be missed when using r₅₀₀ (Takizawa etal. 2010).

A second type of quantitative measure we employ in the study is the centroid shift (see Mohr etal. 1993; O'Hara etal. 2006; Poole etal. 2006; Maughan 2007). Centroids of the X-ray surface brightness are measured successively in apertures of increasing size. The centroid shift is the variation of the centroid with radius (see above references and Jeltema etal. 2008), and is an indicator of deviations from dynamical equilibrium. Our main interest here is in the change of the mean centroid shifts from the non-cold front clusters to the cold front clusters.

We find that clusters identified as having cold fronts have systematically higher values for the power ratios and centroid shifts than those that do not have cold fronts. This result is statistically significant at >99% when a Kolmogorov–Smirnov test is applied. Table1 shows the median values for the power ratios and centroid shifts for clusters identified as cold front clusters and those without cold fronts. Figures17 and 18 show the distribution of P₃/P₀ and centroid shifts in the cluster sample for both cold front and non-cold front clusters in the sample at z = 0. The other power ratios have very similar distributions.

This result holds for all clusters with even a single pixel identified as a cold front, but it becomes monotonically more significant when we increase the threshold for number of pixels identified as cold fronts. This measure is equivalent to total length of cold front features. Figure 19 shows a histogram of the P₃/P₀ power ratio for the clusters in samples with increasing numbers of cold front pixels. As we raise the threshold for number of pixels identified as cold fronts, the morphology becomes systematically more disturbed, as is expected if cold fronts are related to merging events. This trend is repeated in the other power ratios, as well as the centroid shift measure. In Figure 19, the number of cold front pixels identified in the z = 0 simulated clusters increases as the distribution shifts to the more disturbed (right) side of the distribution. As noted in earlier sections, L_cf can be thought of roughly as a total extent of cold front features, where in these images the pixel scale is 15.6 h⁻¹ kpc.