The Concentration–Mass relation of massive, dynamically relaxed galaxy clusters: agreement between observations and ΛCDM simulations

ABSTRACT

The relationship linking a galaxy cluster’s total mass with the concentration of its mass profile and its redshift is a fundamental prediction of the Cold Dark Matter (CDM) paradigm of cosmic structure formation. However, confronting those predictions with observations is complicated by the fact that simulated clusters are not representative of observed samples where detailed mass profile constraints are possible. In this work, we calculate the symmetry-peakiness-alignment (SPA) morphology metrics for maps of X-ray emissivity from the three hundred project hydrodynamical simulations of galaxy clusters at four redshifts, and thereby select a sample of morphologically relaxed, simulated clusters, using observational criteria. These clusters have on average earlier formation times than the full sample, confirming that they are both morphologically and dynamically more relaxed than typical. We constrain the concentration–mass–redshift relation of both the relaxed and complete sample of simulated clusters, assuming power-law dependences on mass (κ_m) and 1 + z (κ_ζ), finding κ_m = −0.12 ± 0.07 and κ_ζ = −0.27 ± 0.19 for the relaxed subsample. From an equivalently selected sample of massive, relaxed clusters observed with Chandra, we find κ_m = −0.12 ± 0.08 and κ_ζ = −0.48 ± 0.19, in good agreement with the simulation predictions. The simulated and observed samples also agree well on the average concentration at a pivot mass and redshift providing further validation of the ΛCDM paradigm in the properties of the largest gravitationally collapsed structures observed. This also represents the first clear detection of decreasing concentration with redshift, a longstanding prediction of simulations in data.

galaxies: clusters: general, cosmology: observations, X-rays: galaxies: clusters

1 INTRODUCTION

The ΛCDM paradigm has as one of its core pillars the assumption that the matter budget of the Universe is dominated by dark matter in the form of cold, weakly interacting, massive particles. Under this paradigm we can make detailed predictions for the specific distribution of mass we expect to find within structures that form hierarchically, most immediately the relationship between total mass and the ‘concentration’ parameter describing the distribution of mass within a halo (Navarro, Frenk & White 1997; Navarro et al. 2004; Gao et al. 2008). Verifying these predictions is an important component of testing ΛCDM. However, observational tests require resolved measurements of the total density within dark-matter-dominated gravitational potentials, which can be difficult to obtain. The most common methods for obtaining such measurements are through gravitational lensing of background galaxies (e.g. Okabe et al. 2013; Du et al. 2015; Merten et al. 2015; Shan et al. 2017) or X-ray observations of the intracluster medium, under the assumption of hydrostatic equilibrium (e.g. Vikhlinin et al. 2006; Schmidt & Allen 2007; Amodeo et al. 2016; Mantz, Allen & Morris 2016b). The approaches have complementary strengths and limitations: lensing is a more direct probe of gravitating mass, but is sensitive to all the mass projected along a line of sight, while X-ray measurements can more easily constrain the three-dimensional mass profile but are reliable only for the most dynamically relaxed clusters due to the assumption of equilibrium. In either case, care is needed to ensure that properties of measured cluster samples are compared with predictions for similarly selected samples from simulations as has been discussed previously in works such as Meneghetti et al. (2014).

In the case of X-ray measurements of dynamically relaxed clusters, the ideal case would be to compare with simulated clusters that satisfy an equivalent observational selection for dynamical relaxation based on mock-imaging information, particularly given that mock-observable proxies for relaxation have been found to correlate only weakly with typical dynamical proxies computed directly from simulated particles (Cao, Barnes & Vogelsberger 2021; De Luca et al. 2021). Typically, measurements of the morphology of the X-ray emitting gas are used as a proxy for dynamical relaxation in such a selection; broadly speaking, such measurements aim to distinguish relaxed clusters based on the sharpness of the surface brightness peak at the cluster centre (Vikhlinin et al. 2007; Santos et al. 2008) and/or the overall appearance of symmetry on larger scales (Mohr, Fabricant & Geller 1993; Buote & Tsai 1995; Jeltema et al. 2005; Nurgaliev et al. 2013; Rasia, Meneghetti & Ettori 2013).

In this paper, we aim to replicate a selection based on the SPA algorithm of Mantz et al. (2015) on mock images of simulated clusters. This algorithm in particular was developed to select relaxed clusters for cosmological tests based on the cluster gas-mass fraction (f_gas), and was subsequently shown to produce a sample with smaller intrinsic scatter in f_gas than the previous ‘by-eye’ selection (Mantz et al. 2014, 2022), consistent with identifying more relaxed systems.

In this paper, we apply the SPA methodology to identify observationally relaxed cluster analogs over a range of redshifts from the three hundred project simulations (Cui et al. 2018), a suite of 324 high-resolution hydrodynamical re-simulations of haloes from the dark-matter-only MultiDark simulations. We then determine the concentration–mass–redshift relation of this relaxed subset, and compare it with constraints from 44 hot, dynamically relaxed clusters observed by the Chandra X-ray telescope satisfying the same selection. Section 2 provides a description of both the simulated and observed data sets. Section 3 describes the fitting procedure and gives results for the concentration–mass–redshift relation for both the simulated and observed samples with a brief comparison with literature results. We summarize our results in Section 4.

2 DATA

2.1 Observed cluster sample

This work employs the sample of 44 hot, dynamically relaxed clusters observed by Chandra and analysed by Mantz et al. (2022). Each cluster is selected to have an X-ray morphology satisfying the SPA criteria for relaxation, and an ICM temperature ≥5 keV in the isothermal part of the profile (⁠|${\sim}0.5\!-\!1\, r_{2500}$|⁠).¹ The sample spans redshifts 0.018 < z < 1.16. The mass–redshift distribution of the sample can be seen in Fig. 1 (black points). The analysis of Mantz et al. (2022, see also Mantz et al. 2014, 2016a) provides constraints on a parametrized NFW (Navarro et al. 1997) density profile for each cluster,

$$\begin{eqnarray*} \rho (r) \propto \bigg(\frac{r}{r_{\mathrm{s}}}\bigg)^{-1} \bigg(1 + \frac{r}{r_{\mathrm{s}}}\bigg)^{-2}, \end{eqnarray*}$$

(1)

under the assumptions of spherical symmetry and hydrostatic equilibrium. In brief, X-ray spectra extracted from concentric annuli are used to simultaneously fit the NFW mass profile and a non-parametric temperature profile that captures the temperature decrease towards the centre of a cool-core cluster. With small variations due to varying cluster morphologies and data quality, spectral information sufficient to constrain temperatures – and thus the mass profile – is used at radii from ∼50 kpc to >r₂₅₀₀, which for typical concentrations spans the NFW scale radius, r_s. The data thus provide information about the normalization of the mass profile and its shape, which can be expressed as joint constraints (including correlation) on M₂₀₀ and the associated concentration, c = r₂₀₀/r_s, shown in Table 1. The adequacy of the NFW model is verified by comparing the density and temperature profiles with those obtained from non-parametric fits (i.e. not assuming any particular mass model; see e.g. the appendix of Mantz et al. 2016a). For further details, we refer the reader to Mantz et al. (2014), Mantz et al. (2016a), Mantz et al. (2022).

Figure 1.

Mass–redshift for observed sample of 44 relaxed clusters (black points) and relaxed simulated sample at each redshift (blue boxes). The boxplots represent the [0, 0.25, 0.5, 0.75, 1.0] percentiles of relaxed clusters at each redshift slice.

Table 1.

Redshifts, masses, and concentrations of clusters in our observed sample from X-ray data. The final column provides the correlation coefficient between measurement uncertainties in ln M₂₀₀ and ln c, which is accounted for in further analysis.

Cluster	z	M200|$(10^{15}\, \mathrm{M}_\odot)$|	c	ρln
Perseus cluster	0.018	1.08 ± 0.10	4.61 ± 0.33	−0.76
Abell 2029	0.078	1.16 ± 0.04	6.36 ± 0.21	−0.95
Abell 478	0.088	1.41 ± 0.17	4.91 ± 0.48	−0.94
PKS 0745−191	0.103	1.48 ± 0.07	5.20 ± 0.18	−0.93
RX J1524.2−3154	0.103	0.45 ± 0.04	7.81 ± 0.81	−0.91
Abell 2204	0.152	1.31 ± 0.11	7.68 ± 0.60	−0.92
RX J0439.0+0520	0.208	0.50 ± 0.08	8.07 ± 1.50	−0.92
Zwicky 2701	0.214	0.59 ± 0.05	5.41 ± 0.42	−0.92
RX J1504.1−0248	0.215	1.31 ± 0.06	7.28 ± 0.25	−0.86
Zwicky 2089	0.235	0.50 ± 0.07	5.60 ± 0.65	−0.93
RX J2129.6+0005	0.235	0.82 ± 0.14	6.29 ± 1.52	−0.92
RX J1459.4−1811	0.236	0.87 ± 0.13	5.38 ± 1.14	−0.82
Abell 1835	0.252	1.78 ± 0.14	4.55 ± 0.37	−0.94
Abell 3444	0.253	1.18 ± 0.20	4.17 ± 0.73	−0.94
MS 2137.3−2353	0.313	0.53 ± 0.05	7.91 ± 0.82	−0.90
MACS J0242.5−2132	0.314	0.73 ± 0.21	7.82 ± 2.32	−0.83
MACS J1427.6−2521	0.318	0.43 ± 0.09	8.94 ± 2.46	−0.87
MACS J2229.7−2755	0.324	0.56 ± 0.11	6.50 ± 1.28	−0.90
MACS J0947.2+7623	0.345	1.32 ± 0.19	5.88 ± 0.67	−0.93
MACS J1931.8−2634	0.352	1.03 ± 0.10	5.87 ± 0.59	−0.92
MACS J1115.8+0129	0.355	0.95 ± 0.12	6.50 ± 0.97	−0.84
MACS J1532.8+3021	0.363	1.24 ± 0.13	4.79 ± 0.34	−0.95
MACS J0150.3−1005	0.363	0.39 ± 0.07	7.66 ± 1.73	−0.78
RCS J1447+0828	0.376	1.32 ± 0.14	5.00 ± 0.34	−0.85
MACS J0011.7−1523	0.378	0.81 ± 0.17	7.02 ± 1.87	−0.90
MACS J1720.2+3536	0.391	0.85 ± 0.14	6.77 ± 1.51	−0.73
MACS J0429.6−0253	0.399	1.12 ± 0.41	4.58 ± 1.41	−0.91
MACS J0159.8−0849	0.404	1.49 ± 0.31	5.75 ± 1.54	−0.89
MACS J2046.0−3430	0.423	0.45 ± 0.09	7.31 ± 1.60	−0.92
IRAS 09104+4109	0.442	0.80 ± 0.16	6.38 ± 1.11	−0.92
MACS J1359.1−1929	0.447	0.67 ± 0.19	5.73 ± 1.50	−0.91
RX J1347.5−1145	0.451	2.86 ± 0.29	7.49 ± 0.52	−0.87
3C 295	0.460	0.62 ± 0.17	5.85 ± 1.37	−0.94
MACS J1621.3+3810	0.461	0.74 ± 0.13	6.98 ± 1.52	−0.89
MACS J1427.2+4407	0.487	0.66 ± 0.12	6.74 ± 1.37	−0.77
MACS J1423.8+2404	0.539	0.76 ± 0.10	7.05 ± 0.90	−0.89
SPT J2331−5051	0.576	0.66 ± 0.05	5.34 ± 0.29	−0.32
SPT J2344−4242	0.596	1.70 ± 0.10	6.57 ± 0.27	−0.51
SPT J0000−5748	0.702	0.66 ± 0.12	5.16 ± 0.75	−0.86
SPT J2043−5035	0.723	0.87 ± 0.10	3.90 ± 0.25	−0.81
SPT J0615−5746	0.972	2.19 ± 0.26	3.26 ± 0.33	−0.77
CL J1415.2+3612	1.028	0.76 ± 0.19	3.86 ± 0.71	−0.90
3C 186	1.063	0.77 ± 0.25	3.80 ± 1.00	−0.90
SPT J2215−3537	1.160	0.75 ± 0.15	5.40 ± 0.70	−0.83

All cosmology-dependent calculations assume a flat ΛCDM model with Hubble constant H₀ = 67.8 km s⁻¹ Mpc⁻¹ and mean matter density Ω_m = 0.307, for consistency with the simulations (see below).

2.2 Simulated cluster sample

The simulated cluster sample is selected from the three hundred project, which is a set of 324 hydrodynamical re-simulations of the most massive clusters (at z = 0) from the (1 Gpc h⁻¹)³ dark-matter-only MDPL2 MultiDark simulation (Klypin et al. 2016).

The SPA selection is based on X-ray emissivity maps generated at four redshifts (z = [0.067, 0.333, 0.592, 0.986]), approximately spanning the redshift range of the observed sample.

For each re-simulation box at each redshift, the map-making procedure identifies the most massive halo in the box (which is usually but not always the same halo at each redshift; Ansarifard et al. 2020). Maps of X-ray emissivity in the rest frame energy band 0.5–2.0 keV, integrated along the line of sight in three orthogonal projections are produced with the code Smac (Dolag et al. 2005). The maps include particles with temperatures above 10⁶ K and densities below the density for star formation. Each particle is weighted by a spline kernel of width equal to the gas particle smoothing length. All maps are 4 × 4 Mpc and 1024 × 1024 pixels, leading to a physical resolution of 3.9 kpc per pixel. The integration length along the line of sight is 10 Mpc. Each of the maps are analysed separately for a total of ∼900 projections at each redshift. Note that these maps contain neither mock-observational artifacts nor emissive sources other than the ICM, though they do contain emission from multiple haloes in each simulation box, if present. Since our intention is to quantify the intrinsic morphology of the target cluster, emissive structures that are unambiguously projections (those appearing in only one projection, since this requires them to be >4 Mpc distant from the cluster centre) were masked in the analysis below.

We computed the SPA metrics independently for each projection of each simulated cluster, as detailed in Appendix A1. Fig. 2 shows that the distributions of the three metrics estimated by the SPA algorithm² – symmetry (s), peakiness (p), and alignment (a) – from the simulations are broadly compatible with the sample of 361 observed clusters analysed by Mantz et al. (2015), although in detail the distributions differ somewhat (Appendix A2). Approximately 8 per cent of projections at each redshift satisfy the SPA criteria for relaxation, with a weak trend in the selected fraction with redshift. This is similar to the relaxed fraction found by Mantz et al. (2015) for clusters selected from Sunyaev–Zel’dovich effect surveys; this is a more fair comparison than with X-ray selected clusters, as peaky clusters (possessing cool, X-ray bright cores) are preferentially detected in X-ray surveys (Rossetti et al. 2017). We consider a simulated cluster to be relaxed if at least two of its projected maps satisfy the SPA criteria. Of all of the clusters, 87.7 per cent are considered relaxed in zero projections, 4.6 per cent are considered relaxed in one projection, 2.4 per cent are considered relaxed in two projections, and 5.3 per cent are considered relaxed in all three projections. Thus, similarly to the individual projections, approximately 8 per cent of the simulated clusters overall are identified as relaxed, with a weak trend in the selected fraction with redshift (Appendix A2).

Figure 2.

Pairwise comparison of the distributions of the three SPA metrics found from simulated X-ray maps (blue; contours containing 68.3 and 95.4 per cent of projections) with observed data from Mantz et al. (2015, black points). Red lines mark the threshold values that must all be exceeded for a cluster to be considered relaxed.

A comparison of the mass–redshift distribution for the simulated (blue boxes) and observed (black points) relaxed cluster samples can be seen in Fig. 1. The two samples cover a similar range in mass and redshift at most redshifts. However, the observed sample at high redshifts (z > 0.8) are all at higher masses than the simulated sample. This is likely due to the large available observational volume relative to the simulated sample, combined with the difficulties of observing lower mass objects at the highest redshifts.

A simple way to quantify the dynamical history of a simulated cluster is through the time since formation (difference in lookback time between formation and observation; Δt_form), with ‘formation’ conventionally defined as the point when half of the clusters mass at the time of observation had been accumulated (e.g. Ludlow et al. 2012); clusters with longer times since formation should be more relaxed. We find that the SPA algorithm does indeed preferentially select clusters with longer times since formation. In practice, the distribution of Δt_form (regardless of selection) varies with the redshift of observation, whereas the distribution of differences in redshift between formation and observation, Δz_form, is relatively static for both the SPA-selected and unselected clusters. The latter thus provides a more straightforward way of quantifying the difference between the two samples, as shown in the left-hand panel of Fig. 3. On average, the SPA-selected sample has a higher value of Δz_form = 0.75 ± 0.28 compared with 0.46 ± 0.25 for the unrelaxed sample; a Kolmogorov–Smirnoff test rejects the hypothesis that the two samples are equivalent with a p-value of <10⁻¹⁶. The right-hand panel of Fig. 3 shows the distributions of Δt_form for the SPA-selected clusters in each redshift slice.

Figure 3.

left: Histograms showing distribution of Δz_form for the complete simulated sample (black) and for the relaxed simulated sample only (blue). right: Kernel density estimate (KDE) plot of the distribution of time since formation for the samples of relaxed simulated clusters at each redshift.

For all the simulated clusters, values of r₂₀₀ and M₂₀₀ were determined directly from 3D mass profiles generated about each cluster centre. Concentrations were estimated by fitting NFW models to the 3D profiles between radii of |$0.05\, r_{200}$| and r₂₀₀ (following Meneghetti et al. 2014), minimizing the difference in log-mass between the model and data, with equal weighting. For the relaxed selection, we expect higher than average concentrations (corresponding to earlier formation times, i.e. longer times since the last major merger), as the next section verifies.

3 CONCENTRATION–MASS–REDSHIFT RELATION

Following (Mantz et al. 2016b) (hereafter M16), we consider a power-law model for concentration as a function of both mass and redshift:

$$\begin{eqnarray*} \ln \bigg(\frac{c}{c_0}\bigg) = \kappa _\zeta \ln {\left(\frac{1+z}{1.35} \right)} + \kappa _m \ln {\bigg(\frac{M_{200}}{10^{15} \, \mathrm{M}_\odot } \bigg)}. \end{eqnarray*}$$

(2)

We use the lrgs linear regression code³ to fit this model to both the observed and simulated cluster data, approximating the uncertainty in masses and concentrations determined from the X-ray data as jointly lognormal and accounting for its covariance. We do not include any uncertainty on masses and concentrations from the simulations. The model includes a lognormal intrinsic scatter, σ, which is assumed to be constant with mass and redshift.

Fig. 4 (left-hand panel) shows the observed clusters, along with the best-fitting mass–concentration relation at z = 0.33 (approximately the sample median). Table 2 lists the individual parameter constraints along with the previous results from M16, based on an earlier slightly smaller version of the observed SPA cluster sample. Compared with the M16 sample, the current data incorporate precise measurements of the Perseus Cluster at z = 0.018, additional relaxed clusters at z = 0.376, 0.972, and 1.160, and numerous re-observations of previously included clusters (see Mantz et al. 2022 for details). Owing to this increase in precision at the lowest and highest redshifts, our updated analysis finds significant evidence of evolution (κ_ζ ≠ 0), unlike the M16 analysis which showed an evolution consistent with zero; in contrast the constraints on κ_m and σ are similar to M16. The right-hand panel of Fig. 4 compares the data with the model predictions at low and high redshifts, providing a visual indication of the strength of the observed evolution.

Figure 4.

Left: Mass–concentration data for the observed sample of 44 relaxed clusters, colour coded from blue to red by increasing redshift. Though not shown here, the correlation between mass and concentration measurements is accounted for in our analysis. The dark curve shows the best-fitting power law with mass at approximately the median redshift of the data, with the shaded region indicating the (68.3 per cent) statistical uncertainty in the fit. The lighter curves delineate the 68.3 per cent probability predictive interval, which includes the impact of the intrinsic scatter. Right: the same data are compared with the model predictions and their statistical uncertainty at redshifts of 0.067 and 0.986.

Table 2.

Mass–concentration relation fit to different samples of clusters. Rows 1–4 correspond simulations and observations analysed in this work, as well as the previous generation of observed SPA-selected clusters. Rows 5–7 show fits to other simulated cluster samples. Rows 8–9 show independent measurements in the recent literature from X-ray and gravitational lensing data. The text provides more details about the different sample selections. Uncertainties on the c₀ values derived from the literature should be considered approximate, since we have no way of accounting for the correlations between parameter uncertainties from those fits. For the Dutton & Maccio work, we have estimated κ_ζ by fitting a power law to the normalization parameters of their fits at z = 0, 0.5, and 1.0.

Cluster sample	Selection	κm	κζ	c0	σ
the three hundred simulations	all clusters	−0.06 ± 0.02	−0.31 ± 0.07	4.31 ± 0.13	0.30 ± 0.01
the three hundred simulations	SPA	−0.12 ± 0.07	−0.27 ± 0.19	5.99 ± 0.35	0.24 ± 0.02
Chandra data (this work)	SPA	−0.12 ± 0.08	−0.48 ± 0.19	5.74 ± 0.20	0.18 ± 0.03
Chandra data (Mantz et al. 2016b)	SPA	−0.16 ± 0.07	−0.17 ± 0.26	6.01 ± 0.23	0.16 ± 0.03
MUSIC-2 simulations (Meneghetti et al. 2014)	CLASH	−0.08 ± 0.02	−0.20 ± 0.09	4.10 ± 0.07	0.16
Planck Cosmology simulations (Dutton & Macciò 2014)	Δr and ρrms	−0.09 ± 0.001	−0.58 ± 0.05	3.72 ± 0.009	0.11
‘Q-Continuum/Outer Rim’ simulations (Child et al. 2018)	Δr	−0.08	−0.35	3.51	–
Strong + Weak gravitational lensing data (Merten et al. 2015)	CLASH	−0.32 ± 0.18	0.14 ± 0.52	3.85 ± 0.17	0.07
Chandra data (Amodeo et al. 2016)	by-eye	−0.50 ± 0.20	0.12 ± 0.61	4.63 ± 3.09	0.14 ± 0.09

Fig. 5 shows a comparison between the observed data and the relaxed clusters in the simulations. The average concentration and evolution with mass show good visual agreement. This agreement is also visible in Fig. 6, which shows the mass–concentration data from the simulations, highlighting the relaxed simulated clusters, along with the model fits to the observed data. It is clear again that the fit to the observed data shows qualitative agreement with the simulated sample of relaxed clusters. As expected, the SPA-selected, relaxed clusters (both simulated and observed) have on average higher concentrations at a given mass than the full simulated sample.

Figure 5.

The observed data (grey) restricted to z < 0.6 compared with the relaxed clusters from the simulated data at z = 0.07 (orange), 0.3 (blue), and 0.6 (green).

Figure 6.

Mass–concentration relation of simulated clusters as a function of redshift (top left: z  = 0.067, top right: z  = 0.333, bottom left: z  = 0.592, and bottom right: z  = 0.986). SPA-selected (i.e. relaxed) clusters from the simulations are shown as filled black stars. The black line is the best-fitting model from the analysis of the observed data at each redshift, with the lighter curves showing the predictive interval (accounting for intrinsic scatter). There is good agreement between the relaxed simulated data and the fits to the relaxed observed clusters.

Table 2 shows the parameter constraints from fitting the functional form for equation (2) to the simulated data (both the relaxed subset and all clusters). As expected, the two fits differ significantly in the normalization parameter, with the typical relaxed cluster having a concentration of 5.8, compared with 4.0 for the average cluster. Interestingly, the constraints on the power-law indices with mass and redshift, and the scatter, are consistent between the relaxed subsample and all clusters in the simulated data set.

Fig. 7 compares the joint constraints on the mass and redshift dependence from the fits to the relaxed observed clusters, relaxed simulated clusters, and all simulated clusters. The degeneracy between these parameters in the case of the simulated samples is a natural consequence of the mass distribution as a function of redshift of the simulated haloes, visible in Figs 1 and 6. The agreement between the three results is excellent, with each be consistent with the others within the 68.3 per cent credible region. The normalization parameter from the observed clusters is similarly in good agreement with that from the relaxed simulations, and both are higher than the normalization fitted to all simulated clusters, as expected (Table 2). We note that intrinsic scatter of both the relaxed subset and all simulated clusters is in excess of that inferred from the observed data. The reason for this is not immediately clear, though we have verified that is not solely due to the few clear outliers in Fig. 6. We check whether the radius used to measure the concentration in the simulations is significant by repeating the fit with concentrations measured out to r₅₀₀ rather than r₂₀₀. While we see a slight increase in the median concentrations at each redshift, we find that the slopes as a function of mass and redshift remain consistent within 1σ.

Figure 7.

Credible regions showing the degeneracy between power-law indices for mass (κ_m) and redshift (κ_ζ) dependence of the concentration. blue: shows the fit including all clusters in the simulation, green: shows the fit only including relaxed simulated clusters, and red: shows the fit to the observational data for relaxed clusters.

We see good agreement with previous estimates of the concentration–mass–redshift relation from both simulations (Dutton & Macciò 2014; Meneghetti et al. 2014; Child et al. 2018) and observational studies of X-ray selected clusters (Amodeo et al. 2016) and weak lensing analyses (Merten et al. 2015) (Table 2). Each of these works apply a different criteria for selecting relaxed clusters. Child et al. (2018) selects relaxed clusters between 0 ≤ z ≤ 1 based on a cut on distance between the halo centre and centre of mass (Δr). Dutton & Macciò (2014) also selects clusters using a cut on Δr, but applies an additional cut based on the goodness of the fit of an NFW profile (ρ_rms). Amodeo et al. (2016) selected clusters from Chandra between 0.4 < z < 1.2 after removing clusters with clear evidence of dynamic activity by eye. Lastly Merten et al. (2015) uses the 25 massive clusters between between 0.19 < z < 0.89 in the CLASH data set, 20 of which were chosen by eye for their relaxed appearance in X-ray data, while the sample in Meneghetti et al. (2014) is explicitly chosen from the MUSIC-2 Simulations to match the mass and other X-ray morphology metric distributions of the CLASH sample. However, despite different selection criteria, all of the simulated studies find a dependence on mass with a power-law index of ∼−0.1. The observational studies generally find a slightly stronger mass dependence, but agree with the simulation predictions within ∼2σ. The slope with redshift shows more variation, with Amodeo et al. (2016) and Merten et al. (2015) finding a slope consistent with no or even positive redshift dependence, while Meneghetti et al. (2014) and Child et al. (2018) prefer a negative slope, in agreement with our findings. However, all of these results are consistent with ours within their relatively large uncertainties. We see more variation in the normalization of the various fits, likely being driven by differences in the selection criteria for each study.

4 SUMMARY

We measured the galaxy cluster concentration–mass–redshift relation from Chandra data for a sample of 44 massive, relaxed clusters, as well as for a set of relaxed clusters from the three hundred project simulated data set at z = [0.067, 0.333, 0.592, 0.996]. We applied equivalent selection criteria in both cases, using the SPA X-ray morphology algorithm (Mantz et al. 2015). We find excellent agreement between the observed data and the predictions for similarly relaxed clusters from ΛCDM simulations in the average concentration, as well as its dependence on mass and redshift. These results provide powerful further validation of the ΛCDM paradigm in the properties of the largest gravitationally collapsed structures observed. We find, for the first time, evidence for evolution in the mean concentration from the observed data, with concentrations nominally decreasing by 28 per cent between redshifts 0 and 1. Such evolution is also clearly present in the full simulated sample, although the simulated relaxed subsample is consistent with no evolution at the 2σ level due to its smaller size (and thus larger uncertainties). Constraints from the relaxed simulated clusters are broadly comparable in power to those from the observed data. Particularly as more relaxed clusters are found and observed, larger suites of hydrodynamical simulations (with mock images and appropriate selection) will be needed to provide precise predictions for comparison with observations.

Acknowledgement

We thank Risa Wechsler for useful discussion, and Benjamin T. Floyd, Taweewat Somboonpanyakul, and Michael McDonald for helpful comments on the draft.

The scientific results reported in this article are based on observations made by the Chandra X-ray Observatory, and data obtained from the Chandra Data Archive. Support for this work was provided by the National Aeronautics and Space Administration through Chandra Awards Number GO5-16122X, GO5-16148X, GO6-17112B, GO0-21124B, and AR9-20012X issued by the Chandra X-ray Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of the National Aeronautics Space Administration under contract NAS8-03060. This research has made use of software provided by the Chandra X-ray Center (CXC) in the application packages CIAO. We acknowledge support from the U.S. Department of Energy under contract number DE-AC02-76SF00515 to SLAC National Accelerator Laboratory. The authors acknowledge The Red Española de Supercomputación for granting computing time to complete The300 project in the Marenostrum machine at the Barcelona Supercomputing Center.

DATA AVAILABILITY

Chandra X-ray data are available from the Chandra Data Archive (CDA) at https://cxc.harvard.edu/cda/. Mantz et al. (2022) lists the specific observation IDs employed in that work, which are identical to those employed here.

Footnotes

1

We jointly define a family of characteristic masses and radii of clusters corresponding to ‘overdensities’ Δ in the conventional way, with |$M_\Delta = (4/3) \pi \Delta \rho _\mathrm{c}(z) r_\Delta ^3$|⁠, where ρ_c(z) is the critical density at a cluster’s redshift.

2

For a detailed definition of the SPA metrics see Appendix A1.

3

https://github.com/abmantz/lrgs

4

We have changed notation compared with M15 in order to avoid confusion with f_s, the the fraction of mass in subhaloes (Section A2).

5

The selection of relaxed simulated clusters in this figure is made a posteriori in order to reduce the scatter at large radii due to merging subhaloes, thus visually clarifying the average difference between simulated and observed profiles, but the difference is present irrespective of the selection.

APPENDIX A: SPA METRICS FROM SIMULATED CLUSTERS

A1 Calculations

The raw simulated X-ray maps encode the projected (line-of-sight integrated) X-ray emissivity in the 0.5–2.0 keV band in the cluster rest frame (see Ansarifard et al. 2020 for details). We convert these to maps of mean count rate in the observer-frame 0.6–2.0 keV band (matching the M15 analysis) using the apec model in xspec, assuming a metallicity of 0.3 Solar (according to the Solar abundance table of Asplund et al. 2009) and using the gas temperature at r₂₅₀₀. While the SPA algorithm accounts for the Poisson nature of X-ray data, it can be straightforwardly applied to the mean count rate maps, effectively assuming infinite signal-to-noise. In order to derive observational metrics related as closely as possible to the target simulated cluster, we masked objects in the images that could be unambiguously identified as interlopers, i.e. those appearing in only one projection, since this requires them to be >4 Mpc distant from the cluster centre (Section 2.2). Such masking is only needed for ∼20 of the ∼900 projections at each redshift.

The SPA metrics are computed by comparing these brightness maps to standard, self-similarly scaled surface brightness levels. The scaling is given by⁴

$$\begin{eqnarray*} B_{\rm S} = \frac{K(z,T,N_\mathrm{H})}{5.7\times 10^8} \frac{E(z)^3}{(1+z)^4} \left(\frac{kT}{\mathrm{keV}}\right) \frac{\mathrm{erg}}{\mathrm{Ms}\, \mathrm{cm}^2\, (0.984\,{\rm arcsec})^2}. \end{eqnarray*}$$

(A1)

Here we depart from M15 by defining K as the conversion factor between bolometric flux in energy units to photon flux (as opposed to energy flux) in the observed energy band. This definition is more appropriate for cooler clusters (kT ≲ 2 keV, given the adopted energy band), for which the energy-to-photon flux ratio varies with temperature; the renormalizing factor of 5.7 × 10⁸ makes the definition approximately equivalent for the hotter clusters studied by M15. The particular brightness levels used in the SPA algorithm are

$$\begin{eqnarray*} S_j = 0.002 \times 10^{0.28j} B_{\rm S}, \end{eqnarray*}$$

(A2)

where j = 0, 1, …, 5.

The motivation for approaching morphology measurement this way is that it allows comparable regions (in terms of enclosed overdensity) of different clusters to be identified based only on observed brightness, given intrinsic self-similarity of the gas density profiles at the relevant radii (outside the core), a roughly self-similar temperature–mass relation (used in deriving equation A1), and the relatively weak cosmological assumptions implicit in the normalized Hubble parameter, E(z). M15 showed that this procedure produces low-scatter scaled surface brightness profiles, and that the isophotes delineated by equation (A2) correspond to approximately comparable radii in units of r₂₅₀₀, without requiring strong assumptions needed to estimate cluster masses (and thus characteristic radii such as r₂₅₀₀) at the outset. However, the regions identified in this way in observed clusters will only be comparable to those identified in simulated clusters if the distribution of gas in the simulated and observed clusters is similar at the radii of interest. In particular, a difference in the gas mass fraction would straightforwardly lead to a difference in the surface brightness associated with a given overdensity, and disproportionately so due to the non-linear dependence of X-ray emissivity on gas density.

We do see evidence for just such a difference between the simulations and the observed data. This can be seen in the left-hand panel of Fig. A1, which compares the measured density profiles for relaxed (SPA-selected) clusters from Mantz et al. (2016a) with the simulated clusters that meet the SPA criteria in at least 2 projections as a function of overdensity.⁵ There is a clear systematic offset at overdensities >1000, with the profiles coming into better agreement at Δ ∼ 500; qualitatively similar behaviour was seen by Li et al. (2020), comparing the same simulations (without the restriction to relaxed clusters) to measurements from McDonald et al. (2017) (see also Campitiello et al. (2022), Sayers et al. (2022)).

Figure A1.

Left: Self-similarly scaled electron density profiles from Chandra observations of relaxed clusters (reproduced from Mantz et al. 2016a, with the Hubble constant adjusted for consistency with the simulations) and simulated clusters satisfying the same morphological selection. Plotting versus the enclosed overdensity, Δ, rather than a simple scaled radius (e.g. r/r₅₀₀) reduces the scatter due to variations in cluster concentration. Blue shading shows the 68 and 95 per cent confidence region from observations, including both statistical uncertainty and intrinsic scatter; red shading shows the 68 per cent region from intrinsic scatter in the simulations. Simulated results are limited to the range of overdensities where most clusters provide data, after excising the central 50 kpc radius. The bottom-left-hand panel shows the ratio to the median observed profile. Right: Comparison of surface brightness between simulated clusters at z = 0.067 (red) and observed clusters (blue; M15) after accounting for self-similar scaling (cf fig. 2 of M15) and adjusting B_S as described in the text. The dashed black levels delineate the five scaled isophotes used in the SPA analysis. Simulated results are limited to the range of overdensities where most clusters provide data, after excising the central 50 kpc radius.

As expected, this translates to an offset in the surface brightness profiles between the simulations and observations in appropriately scaled units, such that directly applying the equations above would mean assessing the morphology of a different region of clusters than the observations we compare to. We therefore adjust B_S by a factor of 2.2 when applied to the simulations, making the scaled surface brightness profiles agree well with the observed data, outside the smallest radii, which are most sensitive to the AGN feedback prescription (right-hand panel of Fig. A1); as expected, this results in the radii identified by the isophotes (see below) resembling those identified in the real data.

With this preprocessing completed, the SPA calculations proceed as described by M15, with the exception that the simulated maps do not require estimation or propagation of measurement uncertainties. Briefly, peakiness (p) is defined as the average brightness in the circular region bounded by the radius where the cluster’s surface brightness profile exceeds S₅ (equation A2), weighted by a factor of 1 + z (Equation 5 of M15). Symmetry (s) and alignment (a) depend on the centres of elliptical shapes fitted to the distribution of pixels falling in the ranges S₀–S₁, S₁–S₂, etc.: s measures how well the centres of these ellipses agree with a global centre, and a measures how well they agree with one another (Equations 6–7 of M15). Clusters with particularly faint centres may never reach the surface brightness thresholds corresponding to the brighter of these isophotes (this can occur in both the simulated and real data), and only those isophotes for which an ellipse can be fit contribute to s and a.

A2 Comparison with observed clusters

SPA values were calculated for each of the cluster projections, which we consider independently unless otherwise noted. The distributions of the SPA metrics for the simulated clusters can be seen in Fig. A2. The top row shows the distributions of parameters at each redshift separately. The bottom row compares a combination, approximately weighted by the redshift distribution of the M15 clusters, with the observed distributions from M15 (z ∈ [0.08, 1.1]). The simulated distributions of s and a are somewhat broader than those measured from data, and offset to smaller values. We see hints of evolution in the mode of the s distributions, but otherwise the s and a distributions are quite consistent in shape across redshifts.

Figure A2.

Plot of s, p, and a distributions for simulated clusters. top: KDE plots for distributions split by redshift. bottom: KDE plot for weighted combination of redshifts. Weights were chosen to approximate the same redshift distribution as the data used in the M15 analysis.

With p, we see a mild but consistent trend towards higher values in the mode and a larger spread with increasing redshift. The low p tail at z = 1.0 is primarily driven by clusters with kT < 5 keV (Fig. A3). Nevertheless, the upper tail of the peakiness distribution, exceeding the SPA selection threshold, is in good agreement with the data; it is also remarkably consistent for the simulated clusters across redshift, except for the lowest-redshift slice.

Figure A3.

The fraction of the projections for our simulated cluster sample identified as relaxed (green), peaky (orange), and undisturbed (blue) as a function of redshift (left) and mass-weighted temperature within r₂₅₀₀ (right). The y-errors represent the 68 per cent confidence interval. right: the x–errors represent the bin sizes in temperature. Only z = 0.3, 0.6 are included in the temperature plot to control for the correlation between redshift and temperature.

We apply the criterion for relaxation defined by M15, and also define ‘peaky’ and ‘undisturbed’ selections analogous to that work. The fraction of the sample identified as relaxed (s > 0.87, p > −0.82, and a > 1.0), peaky (p > −0.82), and undisturbed (s > 0.87 and a > 1.0) as a function of redshift can be seen in Fig. A3 (cf. Fig. A2). We see a slight trend towards increasing relaxed fraction with increasing redshift. The fraction of relaxed clusters for z > 0.3 is roughly consistent with the 11 per cent of observed clusters classified as relaxed by M15.

Comparing with fig. 18 of M15, we see that the relaxed, peaky, and undisturbed fractions are in good agreement with their SZ-selected cluster samples (respectively: 8.5, 14, and ∼20 per cent), tending to be lower on average than the measured fractions for their X-ray selected samples. In total, 11 per cent of the observed and 8 per cent of the simulated sample is identified as relaxed. This corresponds to 44 observed and 96 simulated clusters (15 at z = 0.07, 23 at z = 0.3, 28 at z = 0.6, and 30 at z = 1.0) included in the relaxed sample.

We also show a comparison between the SPA criteria for relaxation and selection criteria for simulated clusters based on the virial ratio (η), the centre-of-mass offset (Δr) and the fraction of mass in subhaloes (f_s) in Table A1. These metrics are commonly used in studies of simulated clusters as simply measurable proxies for relaxation (Cui et al. 2018; Cao et al. 2021; De Luca et al. 2021). As can be seen, the criteria used in previous works bracket the selectivity of the SPA cuts, with the selection criteria from Cui et al. (2018) identifying only 4.5 per cent of clusters as relaxed with limited overlap, with the SPA-selected sample. On the other end, the De Luca et al. (2021) cuts select ∼40 per cent of clusters as relaxed and include almost 80 per cent of the SPA selected clusters. While all selected samples show large scatter in Δz_form, the cuts based on the SPA relaxed sample show a higher value of Δz_form than the relaxed samples defined using previous cuts.

Table A1.

Criteria for cluster relaxation taken from previous works. (1) gives the cuts used to select the sample of simulated clusters, (2) is the total fraction of the sample selected as relaxed across all redshifts, (3) gives the mean and standard deviation of Δz_form for the relaxed sample, and (4–5) show the completeness and purity of the sample relative to the SPA selection.

Selection Criteria	Selected fraction (per cent)	Δzform	Completeness (per cent)	Purity (per cent)
0.85 < η2500 < 1.15, Δr2500 < 0.04, fs, 2500 < 0.1 (Cui et al. 2018)	4.5	0.63 ± 0.29	9.2	17.7
Δr500 < 0.1, fs, 500 < 0.1 (De Luca et al. 2021)	43.3	0.62 ± 0.30	78.5	15.7
Δr500 < 0.07, fs, 500 < 0.1 (Cao et al. 2021)	40.2	0.63 ± 0.30	77.2	16.7
s > 0.87, a > 1.0, p > −0.82 (M15)	8.7	0.75 ± 0.28	–	–