# rhapsody-g simulations – I. The cool cores, hot gas and stellar content of massive galaxy clusters

Oliver Hahn, Davide Martizzi, Hao-Yi Wu, August E. Evrard, Romain Teyssier, Risa H. Wechsler

► To cite this version:

Oliver Hahn, Davide Martizzi, Hao-Yi Wu, August E. Evrard, Romain Teyssier, et al., rhapsody-g simulations – I. The cool cores, hot gas and stellar content of massive galaxy clusters. Monthly Notices of the Royal Astronomical Society, 2017, 470 (1), pp.166-186. 10.1093/mnras/stx001. hal-01645183

HAL Id: hal-01645183

https://hal.science/hal-01645183v1

Submitted on 25 Aug 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

# RHAPSODY-G simulations – I. The cool cores, hot gas and stellar content of massive galaxy clusters

Oliver Hahn,$^{1,2*}$ Davide Martizzi,$^3$ Hao-Yi Wu,$^{4,5}$ August E. Evrard,$^5$
Romain Teyssier$^6$ and Risa H. Wechsler$^{7,8}$

Accepted 2017 January 3. Received 2016 December 22; in original form 2015 September 15

## ABSTRACT

We present the RHAPSODY-G suite of cosmological hydrodynamic zoom simulations of 10 massive galaxy clusters at the $M_{vir} \sim 10^{15} M_{\odot}$ scale. These simulations include cooling and subresolution models for star formation and stellar and supermassive black hole feedback. The sample is selected to capture the whole gamut of assembly histories that produce clusters of similar final mass. We present an overview of the successes and shortcomings of such simulations in reproducing both the stellar properties of galaxies as well as properties of the hot plasma in clusters. In our simulations, a long-lived cool-core/non-cool-core dichotomy arises naturally, and the emergence of non-cool cores is related to low angular momentum major mergers. Nevertheless, the cool-core clusters exhibit a low central entropy compared to observations, which cannot be alleviated by thermal active galactic nuclei feedback. For cluster scaling relations, we find that the simulations match well the $M_{500}-Y_{500}$ scaling of Planck Sunyaev–Zeldovich clusters but deviate somewhat from the observed X-ray luminosity and temperature scaling relations in the sense of being slightly too bright and too cool at fixed mass, respectively. Stars are produced at an efficiency consistent with abundance-matching constraints and central galaxies have star formation rates consistent with recent observations. While our simulations thus match various key properties remarkably well, we conclude that the shortcomings strongly suggest an important role for non-thermal processes (through feedback or otherwise) or thermal conduction in shaping the intracluster medium.

**Key words:** methods: numerical – galaxies: clusters: general – dark matter – large-scale structure of Universe – X-rays: galaxies: clusters.

## 1 INTRODUCTION

While simulations of galaxy formation in Milky Way-sized haloes (e.g. Guedes et al. 2011; Christensen et al. 2014) or small cosmic volumes (e.g. Dubois et al. 2014; Vogelsberger et al. 2014; Schaye et al. 2015) are making substantial progress, realizing the population of galaxies that reside in the most massive cosmic haloes, those hosting rich clusters of galaxies, remains a formidable challenge (Kravtsov & Borgani 2012). Compared to a $10^{12} M_{\odot}$ Milky Way

halo, simulating a $10^{15} M_{\odot}$ cluster requires at least an extra order of magnitude in spatial resolution and three orders of magnitude in mass, thus requiring larger and longer computations. In addition, the Gaussian random field nature of the initial conditions biases the progenitors of clusters towards extreme systems at all redshifts, meaning cluster evolution is tied to that of the first stars, the earliest protogalaxies, and the most massive quasars and their supermassive black hole (SMBH) interiors at high redshifts. Finally, the deep nature of the gravitational potential well retains most of the cosmic baryonic content associated with the dark matter, implying that a complex mix of coupled hydrodynamic, magnetohydrodynamic, chemical and radiative processes ought to be solved in order to

<footnote>$^1$Laboratoire Lagrange, Université Côte d'Azur, Observatoire de la Côte d'Azur, CNRS, Blvd de l'Observatoire, CS 34229, F-06304 Nice cedex 4, France</footnote>

<footnote>$^2$Department of Physics, ETH Zurich, CH-8093 Zürich, Switzerland</footnote>

<footnote>$^3$Department of Astronomy, University of California, Berkeley, CA 94720-3411, USA</footnote>

<footnote>$^4$California Institute of Technology, MC 367-17, Pasadena, CA 91125, USA</footnote>

<footnote>$^5$Department of Physics, University of Michigan, Ann Arbor, MI 48109, USA</footnote>

<footnote>$^6$Institute for Computational Science, University of Zurich, CH-8057 Zürich, Switzerland</footnote>

<footnote>$^7$Kavli Institute for Particle Astrophysics and Cosmology & Physics Department, Stanford University, Stanford, CA 94305, USA</footnote>

<footnote>$^8$SLAC National Accelerator Laboratory, Menlo Park, CA 94025, USA</footnote>

<footnote>* E-mail: oliver.hahn@oca.eu</footnote>
gain insights into the cycling and transport of mass, energy/entropy, momentum and metals across ~10 Mpc regions over most of the 13.8 Gyr history of the Universe.

Early hydrodynamic simulations with gas cooling successfully formed multiple galaxies within a cluster (Katz & White 1993; Evrard, Summers & Davis 1994), but the central galaxies were too massive compared to observations. Radio and X-ray observations of cavities near the central galaxy highlighted the need for strong feedback to curtail central cooling (cf. McNamara & Nulsen 2007, 2012). Inclusion of SMBH feedback into semi-analytic methods applied to halo ensembles from N-body simulations led to a reduction in central galaxy stellar masses (Croton et al. 2006; De Lucia & Blaizot 2007). Application within hydrodynamic simulations required estimating accretion rates on to the SMBH of active galactic nuclei (AGNs), a task that in principle requires very high spatial resolution. Schemes such as that developed by Booth & Schaye (2009) were promoted to estimate Bondi–Hoyle accretion rates from simulations with roughly kiloparsec spatial resolution. Gas dynamic simulations using this approach to AGN feedback show improvements to central cluster galaxy morphology (Martizzi, Teyssier & Moore 2012a) as well as improved scaling of hot intracluster medium (ICM) properties with halo mass (e.g. Le Brun et al. 2014; Plabelles et al. 2014, and references therein).

The detailed nature of AGN feedback and its relationship to the phenomenology of the core plasma in clusters remain subjects of active investigation. X-ray observations show that the massive cluster population can be divided into cool-core (CC) and non-cool core (NCC) categories, with the distinction determined by the strength of a surface brightness cusp within the inner ~100 kpc (Allen, Evrard & Mantz 2011). The physics that controls this CC/NCC dichotomy has remained under debate. Major mergers may be capable of driving a CC to NCC transition through shock heating and/or ablation of core gas, but early cosmological gas dynamic simulations with cooling and supernova (SN) feedback found major mergers to be ineffective at heating the core gas (Burns et al. 2008) unless the merger had sufficiently low angular momentum (Poole et al. 2008). The observational finding that morphologically disturbed clusters rarely host CCs offers however strong empirical evidence that mergers play a role in determining core morphology (Pratt et al. 2010).

Still, it remains possible that internal processes, such as AGN jet-driven turbulent heating, may drive the transition from a CC to an NCC state, at least temporarily (Guo & Oh 2009; Parrish, Quataert & Sharma 2009; Guo & Mathews 2010; Parrish, Quataert & Sharma 2010; Ruszkowski & Oh 2010) or may play an important role in pressurizing the cores. Idealized, high-resolution AMR simulations find cold, chaotic accretion on to the SMBH at rates many times the Bondi expectation (Gaspari, Ruszkowski & Oh 2013b), and this has motivated models in which precipitation of turbulent core gas serves to self-regulate AGN accretion and feedback (Gaspari, Brighenti & Ruszkowski 2013a; Li & Bryan 2014; Li et al. 2015). Voit et al. (2015) provide empirical support that such a model, coupled with moderate conductive heating from the external reservoir of ICM plasma, may help explain CC/NCC dichotomy. Ultraviolet (UV) imaging of brightest cluster galaxies in the CLASH (Cluster Lensing And Supernova survey with Hubble) cluster sample (Donahue et al. 2015) reveals knots and extended filamentary structures suggestive of the bipolar streams that emerge in the simulations of Li et al. (2015).

In this paper, we present results from a suite of gas dynamic simulations of 10 high-mass haloes. These simulations extend our previous suite of N-body simulations of massive clusters falling inside a very narrow mass range at $z=0$ (Wu et al. 2013a,b)

to multiphysics adaptive mesh simulations. They include cooling, star formation and SMBHs as well as their respective feedback. These subgrid models have been shown before to reproduce realistic brightest cluster galaxy (BCG) masses (Martizzi et al. 2012a). Here, we investigate how well our simulations reproduce the observed X-ray properties of galaxy clusters, most notably their density, temperature and entropy profiles, and find that, much like in observed systems, our simulations produce a clear CC/NCC dichotomy most clearly seen in the entropy profiles for which we are able to give a physical explanation in terms of major cluster mergers. Reproducing the properties of the hot cluster plasma however, has to go hand-in-hand with achieving also compatible results for the full galaxy population in and around clusters. While the field is still far from predictive simulations across the full multiwavelength range covered by current and future observational studies of clusters, we are able to highlight several successes as well as important shortcomings of such state-of-the-art cluster simulations in a full cosmological context.

The structure of this article is as follows. In Section 2, we introduce the sample of zoom simulations that we study as well as the numerical methods and models that we employ. We then compare the ICM profiles of the simulated clusters with various X-ray observations in Section 3. In Section 4, we investigate the origin of the CC/NCC dichotomy that we find in our sample. We then compare the stellar properties of the galaxies formed in the simulations with observational data in Section 5. Finally, in Section 6, we study the evolution of the simulated clusters along several mass-observable scaling relations important for cosmology. We discuss the impact of specific modelling choices and the influence of numerical effects on our results in Section 7, before we summarize our results and conclude in Section 8.

## 2 DESCRIPTION OF THE SIMULATIONS AND THE NUMERICAL APPROACH

In this section, we describe how we selected a representative sample of massive clusters at fixed mass at $z=0$ from a simulation of a large cosmological volume, combining an average and a tail-biased set of objects selected from a larger sample. We also discuss the details of our numerical and algorithmic approaches, in particular the subgrid models we have employed to account for subresolution physics due to cooling and energy injection by stars and massive black holes.

### 2.1 Initial conditions and general approach

The current RHAPSODY-G simulation suite includes 10 hydrodynamical zoom-in simulations selected from the original RHAPSODY sample of massive galaxy clusters (Wu et al. 2013a). Nine are chosen to have a similar final mass of $M_{vir} \approx 6 \times 10^{14} h^{-1} M_{\odot}$ and the 10th has $M_{vir} \approx 1.3 \times 10^{15} h^{-1} M_{\odot}$. The original RHAPSODY clusters were identified from one of the CARMEN simulations from the LArge Suite of DArk MAter Simulations (LASDAMAS).<fnref n="1" /> The CARMEN simulation on which we base our re-simulations is a cosmological volume of $1 h^{-1} \text{Gpc}$. All simulations are based on a $\Lambda$CDM cosmology with density parameters $\Omega_b = 0.045$ for baryons and $\Omega_m = 0.25$ for total matter, $\Omega_{\Lambda} = 0.75$ for the cosmological constant. The long-wave spectral index is $n_s = 1$, the amplitude normalization is $\sigma_8 = 0.8$ and the Hubble parameter is $h = 0.7$. The original RHAPSODY haloes were

<footnote n="1">¹ http://lss.phy.vanderbilt.edu/lasdamas/</footnote>
identified based on a GADGET-2 (Springel 2005) simulation of the full box using 1120³ particles. The subsequent N-body-only RHAPSODY simulations were also carried out using GADGET-2 but using 'zoom' initial conditions (using MUSIC; Hahn & Abel 2011) for a sphere of 8 h⁻¹ Mpc centred on each selected cluster at z = 0 with an effective resolution of 4096³ (4K) and 8192³ (8K) particles (see Wu et al. 2013a, for details). All initial conditions (both for the original box and all subsequent zooms) were performed using second-order Lagrangian perturbation theory at z = 50.

In this paper, we regenerate the respective initial conditions using MUSIC for 10 clusters (see Section 2.2 below for details on the selection and properties of these haloes) but now including baryon perturbations. We assume here that baryons are fully tracing the dark matter perturbations already at z = 50, which is accurate enough for our purposes here (and common procedure), but strictly speaking not correct in detail (see e.g. Angulo, Hahn & Abel 2013, for a detailed discussion). In addition, we use also Lagrangian perturbation theory for the baryons using the local Lagrangian approximation (see Hahn & Abel 2011, for details on baryon initial conditions (ICs) for grid and particle codes).

## 2.2 The RHAPSODY-G subset of RHAPSODY clusters

From the sample of ~100 clusters in the $M_{vir} = 10^{14.8 \pm 0.05} h^{-1} M_{\odot}$ mass bin of the original RHAPSODY sample (Wu et al. 2013a), we selected a subset of 10 clusters that we investigate using multiphysics simulations in this paper. To sample both extreme cases of formation history and more average clusters, we considered the two-dimensional ordering shown in fig. 2 of Wu et al. (2013a), where clusters are ranked first by halo concentration c, and then at similar concentration ranked a second time by the number of subhaloes $N_{sub}$ above $v_{max} > 100 \text{ km s}^{-1}$. The extreme corners of this space are occupied by the systems with IDs 337, 377, 572 and 653. We additionally included system 545 that has similar properties to 337 (high concentration, high substructure fraction). We note that this naturally includes the fossil system 572 that was discussed in more detail already in Wu et al. (2013b) as being a particularly early-forming system and occupying the tail in many halo properties. We complement this subset of 'extreme' clusters with four more clusters taken from the central region of the c-Nsub space. These are the ones with IDs 211, 348, 361 and 448. Additionally, we include cluster 474, which has a mass that is twice as high as the rest of the sample. These clusters, along with various properties discussed in Section 3.1 are listed in Table 1.

## 2.3 Hydrodynamics, N-body and gravity

In order to follow the non-linear multiphysics evolution of our initial conditions, we use here the Eulerian adaptive mesh refinement code RAMSES (Teyssier 2002) and include radiative cooling, as well as subgrid models for star formation and AGN feedback, which we discuss below. RAMSES is based on a MUSCL scheme with an approximate Riemann solver and gravity is solved using the multigrid method directly on the AMR hierarchy. Dark matter as well as stellar and sink particles are evolved using standard N-body techniques.

We employ an overdensity-based (i.e. 'Lagrangian') refinement strategy that splits cells if they reach an overdensity of eight, i.e. refinement of the base grid by n additional levels requires a density of $8^n \bar{\rho}$. The maximum refinement level for the 4K runs (which constitute the main part of our analysis here) is determined by maintaining a smallest cell size of physical 3.8 h⁻¹ kpc. The 4K dark matter N-body particle mass is 8.27 × 10⁸ h⁻¹ M$_{\odot}$, and the

<table><caption>**Table 1.** The haloes selected from the original RHAPSODY sample for the RHAPSODY-G simulations listed along with the $z=0$ masses and the formation redshift $z_{1/2}$. We also list the number of snapshots that enter our stacking analysis and the mean redshift of the stack, as well as whether the cluster is classified as CC in the stack analysis (see Section 3.1 for details). The system RG 572 is a fossil cluster, a significant outlier in almost all its halo properties, but is close to a CC system.</caption><thead><tr><th>Cluster</th><th><i>M</i><sub>vir,0</sub><br>(×10<sup>15</sup> M<sub>☉</sub>)</th><th><i>M</i><sub>500c,0</sub><br>(×10<sup>14</sup> M<sub>☉</sub>)</th><th><i>z</i><sub>1/2</sub></th><th><i>N</i><sub>stack</sub></th><th><i>z̄</i><sub>stack</sub></th><th>CC?</th></tr></thead><tbody><tr><td>RG 211</td><td>1.00</td><td>5.02</td><td>0.27</td><td>7</td><td>0.07</td><td>—</td></tr><tr><td>RG 337</td><td>1.06</td><td>6.59</td><td>0.71</td><td>20</td><td>0.37</td><td>—</td></tr><tr><td>RG 348</td><td>1.15</td><td>6.28</td><td>0.67</td><td>22</td><td>0.30</td><td>—</td></tr><tr><td>RG 361</td><td>1.07</td><td>5.46</td><td>0.65</td><td>19</td><td>0.25</td><td>+</td></tr><tr><td>RG 377</td><td>1.08</td><td>4.89</td><td>0.45</td><td>7</td><td>0.10</td><td>+</td></tr><tr><td>RG 448</td><td>1.03</td><td>5.19</td><td>0.73</td><td>12</td><td>0.18</td><td>+</td></tr><tr><td>RG 474</td><td>2.69</td><td>10.38</td><td>0.47</td><td>4</td><td>0.61</td><td>—</td></tr><tr><td>RG 545</td><td>0.93</td><td>5.12</td><td>0.77</td><td>8</td><td>0.08</td><td>+</td></tr><tr><td>RG 572</td><td>0.96</td><td>5.65</td><td>1.11</td><td>30</td><td>0.37</td><td>(+)</td></tr><tr><td>RG 653</td><td>0.84</td><td>3.75</td><td>0.21</td><td>5</td><td>0.07</td><td>—</td></tr></tbody></table>

initial mass per hydro cell is ~1.82 × 10⁸ h⁻¹ M$_{☉}$, respectively. The 8K run of RG 653, which we consider here for convergence purposes, has twice the mass and force resolution, i.e. 1.9 h⁻¹ kpc, and an eight times smaller N-body particle mass than the 4K runs. The high-resolution Lagrangian patch from which the 8 h⁻¹ Mpc sphere centred on each cluster will form is tagged using a passive scalar colour field that is advected with the gas. Dynamic refinement is restricted to the regions where this colour field is non-zero. We thus focus all computational resources on the forming cluster and its immediate environment.

## 2.4 Modelling cooling, star formation, stellar feedback and chemical evolution

The simulations include optically thin radiative cooling using the cooling rates of Sutherland & Dopita (1993) for H, He and metal line cooling. Different metals are not evolved separately, but the total gas metallicity is advected with the hydrodynamical equations as a passive scalar and is sourced by the SN feedback model. A UV background is included assuming the parametrization of Haardt & Madau (1996) and an instantaneous reionization at $z=10$ taking into account an earlier reionization in the particularly overdense protocluster environment.

The unresolved cold and dense gas that will constitute the interstellar medium (ISM) of galaxies is approximated using a polytropic temperature floor

$$T_{\text{floor}} = T_{*} \left( \frac{n_{\text{H}}}{n_{*}} \right)^{\gamma_{*}-1}, \qquad (1)$$

where $n_{*} = 0.1 \text{ cm}^{-3}$ is the threshold for star formation (see below), $T_{*} = 10^{4}$ K is a characteristic temperature mimicking thermal and turbulent motions in the ISM and $\gamma_{*} = 5/3$ is the effective polytropic exponent. In practice, gas can be heated above the temperature floor, but cannot cool below it.

We assume star formation to occur in a cell when the gas density exceeds $n_{*}$. In this case, a star particle is spawned carrying 20 per cent of the mean baryon mass of the cell. We set the local star formation rate (SFR) per cell as

$$\dot{\rho}_{*} = \epsilon_{*} \frac{\rho_{\text{gas}}}{t_{\text{fl}}}, \qquad (2)$$
with the free-fall time $t_{ff} = (3\pi/32 G\rho)^{1/2}$ of the cell and a star formation efficiency per free-fall time $\epsilon_* = 0.02$. Stellar particles are seeded locally from a Poisson process.

SN feedback is implemented based on the model of Dubois & Teyssier (2008), in which each newly formed star particle releases a fraction $\eta = 0.1$ of its mass and metals with a yield of $y = 0.1$ into its surrounding cells through SNe after 10 Myr (this implies that 1 per cent of the time integrated global SFR is returned as metals to the ISM). In addition, each SN injects an energy of $10^{51}$ erg into the surrounding ISM, which regulates the star formation efficiency at galaxy scale halo masses. The free parameters of the SN feedback have been calibrated to reproduce stellar masses consistent with abundance-matching results (e.g. Behroozi, Wechsler & Conroy 2013c) at halo masses $M_{\text{halo}} \lesssim 10^{12} M_{\odot}$ for haloes resolved with $\gtrsim 1000$ particles (see also fig. 3 of Dubois et al. 2012).

## 2.5 Modelling AGN feedback

The deep potential wells of galaxy clusters require a stronger heating source than SN feedback to prevent a central cooling catastrophe. AGNs provide a natural source of energy in massive galaxies to offset these extreme cooling flows (see e.g. McNamara & Nulsen 2007; Fabian 2012, for reviews). In the simulations discussed in this paper, we include a purely thermal AGN feedback, based on the subgrid models of Springel, Di Matteo & Hernquist (2005) with the additional energy injection thresholding of Booth & Schaye (2009), which is commonly employed in smoothed-particle hydrodynamics (SPH) simulations (e.g. Le Brun et al. 2014; Schaye et al. 2015, for recent examples). We give a brief summary of our AGN feedback model below but refer the reader to Paper 2, Martizzi et al. (2016), for a more detailed exposition.

### 2.5.1 Thermal AGN feedback model

Sink particles, representing SMBHs, are created when contiguous regions of high-density gas that (1) do not already contain a sink particle, (2) are gravitationally bound and (3) have an accretion rate $\dot{M}_{\text{clump}} = \frac{M_4}{t_{ff}} > 30 M_{\odot} \text{yr}^{-1}$, where $M_4$ is the mass contained within a sphere of four computational cells, exceed $10^{-29} \text{ g cm}^{-3}$ (identified on the fly by the clustering method of Bleuler et al. 2015). This guarantees that sink particles are spawned only in the most massive haloes at high redshift. The initial black hole mass is chosen proportional to the clump accretion rate times the Salpeter time.

Sink particles then accrete at a boosted Bondi-Hoyle accretion rate (e.g. Booth & Schaye 2009)

$$ \dot{M}_{\text{BH}} = 4\pi\alpha_{\text{boost}}\tilde{r}_{\text{B}}^2\rho\sqrt{u^2 + c_s^2}, \quad (3) $$

where $u$ is the fluid velocity, $c_s$ the sound speed, $\tilde{r}_B = \min(r_B, 4\Delta x)$ the free-fall limited Bondi radius and

$$ \alpha_{\text{boost}} = \begin{cases} \left(\frac{n_H}{n_*}\right)^2 & \text{if } n_H > n_* = 0.1 \text{ H per cc,} \\ 1 & \text{otherwise} \end{cases} \quad (4) $$

the density-dependent boost factor. We additionally impose the Edgington limit on to $\dot{M}_{\text{BH}}$ assuming that $\epsilon_r = 0.1$ of the accreted rest mass energy is converted into radiation. Numerically, the amount $\dot{M}_{\text{BH}}\Delta t$ up to a maximum of half of the gas mass contained in each cell is removed from cells inside the sink radius. The thermal energy $\Delta E = \epsilon_c \epsilon_r \dot{M}_{\text{BH}} c^2 \Delta t$ would in principle be released into the gas in each time-step, but following Booth & Schaye (2009), we accumulate $E_{\text{AGN}} = \sum \Delta E$ and inject the accumulated $E_{\text{AGN}}$

into the sink radius (i.e. the sphere of radius four cells) only once; $E_{\text{AGN}} > \frac{3}{2} m_{\text{gas}} k_B \Delta T$, where $m_{\text{gas}}$ is the gas mass enclosed by the sink radius. The coupling efficiency $\epsilon_c \approx 0.15$ has been calibrated to reproduce the $M_{\text{BH}}-\sigma$ relation (Teyssier et al. 2011). In our fiducial model, we set $\Delta T = 10^7$ K, but vary this by one order of magnitude up and down when we investigate the effect of trading fewer violent AGN events against more frequent less violent AGN events in Section 7.1. We emphasize however, that the total energy injected by the AGN obviously remains more or less the same. The parameter $\Delta T$ has been shown to allow a tuning of the bulk properties of the ICM (e.g. Le Brun et al. 2014) in SPH simulations.

### 2.5.2 Phenomenological AGN feedback model

Since the AGN energy is thus injected close to the resolution limit into a region containing the densest cells of the simulation, we also consider a more phenomenological model, inspired by the AGN model of Battaglia et al. (2010), in which we distribute $E_{\text{AGN}}$ over a resolved sphere of a radius determined by the black hole mass. Battaglia et al. (2010) used a formula for the injection radius that depends on the halo mass. Since we do not track halo masses on the fly with our algorithm and prefer a local criterion, we decided to express the injection radius used by Battaglia et al. (2010) in terms of the black hole mass by using the relation $M_{\text{BH}} \propto V_c^4$, where $V_c$ is the circular velocity of the halo (Volonteri, Natarajan & Gültekin 2011, noting that for an NFW profile (Navarro et al. 1997) profile the circular velocity is related to halo mass as $V_c \propto M_{200}^{1/3}$). The radius of the injection sphere is then given by

$$ R_{\text{AGN}} = \max \left[ 4\Delta x(1+z), \times 100 \text{ kpc } h^{-1} \times E(z)^{-2/3} \left( \frac{M_{\text{BH}}}{7.6 \times 10^{10} M_{\odot}} \right)^{1/4} \right], \quad (5) $$

where $E(z)$ is the ratio of the Hubble constant at redshift $z$ and $z=0$, i.e. $H_0$:

$$ E(z) = \frac{H(z)}{H_0} = \sqrt{\Omega_{\text{m}}(1+z)^3 + \Omega_{\Lambda}}. \quad (6) $$

$R_{\text{AGN}}$ is larger than four cell sizes at most times. This additionally reduces the thermal energy inserted into each individual cell, reducing immediate loss of energy through efficient cooling at the highest densities and guaranteeing that a large central region is heated by the AGN. Again, we only consider this alternative model, which we term the ‘phenomenological model’, when investigating the impact of changes in the thermal feedback model parameters in Section 7.1.

## 2.6 Halo/galaxy finding and their properties

In order to identify haloes, including their stellar, black hole and gaseous content in our simulation, we use a heavily modified version of ROCKSTAR-GALAXIES, which is a special version of the original version of ROCKSTAR (Behroozi, Wechsler & Wu 2013a) adapted to work also for simulations including gas and stars (see e.g. Knebe et al. 2013, where this version was used). In order to interface RAMSES with ROCKSTAR-GALAXIES, we convert all leaf cells of the AMR tree to pseudo-particles of variable mass, which is however typically not varying over more than an order of magnitude due to the Lagrangian refinement strategy.

During the (sub)halo finding, we then calculate all relevant halo and galaxy properties directly inside ROCKSTAR-GALAXIES. Unbinding of gravitationally unbound particles and cells is performed, and we
obtain the masses, radii, centres and bulk velocities of the dark matter, stellar, gaseous and black hole content of each (sub)halo. We use the CONSISTENT-TREES code (Behroozi et al. 2013b) to link the halo/galaxy catalogues across redshifts based on only the dark matter particle IDs.

We additionally compute several galaxy-related quantities, such as the SFR, mean stellar age, surface brightness and magnitude in various photometric bands. For the SFR, we simply sum the mass of all star particles younger than 100 Myr, divided by 100 Myr. For the surface brightness and magnitude calculations, we determine the luminosity of each star particle in a given filter applied to the spectral energy distribution obtained from the simple stellar population models from STARBURST99 (Leitherer et al. 1999), assuming the metallicity and age of the star particle as the mean age and metallicity of the stellar population. In order to derive the surface brightness-limited galaxy masses, we perform a one-dimensional projection and integrate the enclosed stellar mass out to a given surface brightness limit. We note that we do *not* take dust absorption into account, which might affect our comparison with observations.

Finally, we also calculate several properties related to the hot plasma, such as the X-ray luminosity and an X-ray emissivity per cell, in order to allow us to weight quantities by their observability in X-ray observations, based on an APEC emission model (Smith et al. 2001), similar to Biffi et al. (2012). We neglect however, additional sophistications such as actual photon sampling, the effect of the point-spread function (PSF) or spectral fitting in order to determine the X-ray temperature, which are clearly beyond the scope of this paper. Finally, we also compute the Sunyaev–Zeldovich (SZ) flux decrement for each cluster. We detail how we perform X-ray and SZ measurements separately in Sections 6.1 and 6.2, respectively, in the context of the evolution of the clusters along scaling relations.

## 3 COMPARISON OF ICM PROFILES WITH OBSERVATIONS

We perform a detailed comparison of density, temperature, entropy and gas depletion profiles with X-ray observations in this section. To show both the mean and deviations over cosmic time of the simulated profiles, we consider individually stacked profiles over time for each simulated cluster for all snapshots for which the cluster mass falls into a narrow mass range. We discuss this stacking procedure in detail below.

### 3.1 Global cluster properties: masses and cool cores

We first discuss some global properties of our sample of simulated galaxy clusters. To give a visual impression of the diversity of the sample, Fig. 1 shows density-weighted temperature and density maps for all clusters at $z=0$. In Table 1, we list all the simulated clusters including their $z=0$ masses and their half-mass redshift $z_{1/2}$, i.e. the highest redshift at which the most massive progenitor exceeds half the final mass. We quote two $z=0$ mass measurements: the virial mass (calculated using the overdensity criterion of Bryan & Norman 1998) and $M_{500c}$, the mass enclosed inside 500 times the critical density.

For the comparison of the cluster ICM profiles with observational data from the ACCEPT sample (Cavagnolo et al. 2009), in Section 3, we decided to select a very narrow mass bin $4 < M_{500c}(z)/10^{14} M_{\odot} < 6$ in both the ACCEPT sample and the simulated clusters. For the simulated clusters, we take every fifth snapshot out of the ~300 snapshots that we have stored for each cluster for which it falls into this mass bin. The redshift distribution

of both RHAPSODY snapshots as well as ACCEPT clusters selected in this way is in reasonable agreement, as shown in Fig. 2. We see that the redshift sampling is around a very similar median redshift, with the RHAPSODY cluster distribution having a tail to somewhat higher redshifts around $z \sim 0.5$ however. Since this difference is rather modest and the scatter around mean profiles in each stack is also small (see below), we did not apply an additional redshift selection criterion.

We then perform a stack analysis for these snapshots separately for each cluster in order to quantify the mean profiles and the variance around these mean profiles (see Section 3.2 for details). In Table 1, we thus also give the number of snapshots $N_{stack}$ used in the stack analysis of each cluster and the mean redshift $\bar{z}_{stack}$ of the stack. $N_{stack}$ is most of the time, but not always, a larger number if the cluster has a high formation redshift. Finally, we indicate whether the cluster is classified as a CC or NCC cluster based on whether the mean central entropy of the stack entropy profile is below $40 \text{ keV cm}^2$ at $r = 10 \text{ kpc}$. We find that RG 361, RG 377, RG 448 and RG 545 are classified as CC systems according to this classification. Halo RG 572 is close to a CC in principle as well, but we will treat it separately as its core properties are even more extreme. The other systems are NCC clusters. We note that we see no immediate connection between either $z_{1/2}$ or the time the cluster spends in the $4 < M_{500c}(z)/10^{14} M_{\odot} < 6$ mass bin and the CC/NCC distinction. Two CC clusters are from the ‘extreme’ sample and two from the ‘average’. The ratio of four out of nine (excluding the fossil system RG 572 as an extreme outlier) is, within the errorbars of small-number statistics, consistent with the number of CC systems in our ACCEPT subsample (6 out of 28) and with other observational estimates at these mass scales (e.g. Chen et al. 2007; Eckert, Molendi & Paltani 2011) of ~40 per cent.

### 3.2 Comparison with ICM profiles from X-ray data

We next compare the gas temperature, entropy and electron density profiles to observational data. We specifically confront the ICM X-ray profile data from the ACCEPT sample (Cavagnolo et al. 2009). As mentioned above, in order to enable a precise comparison with our RHAPSODY-G clusters, we perform a stringent mass cut of the ACCEPT clusters. We combine X-ray and SZ masses to this end to reduce the influence of hydrostatic mass bias at least to some degree. Specifically, for each cluster in the ACCEPT sample, we find its $M_{500}$ masses from the MCXC catalogue (Piffaretti et al. 2011) and its $M_{SZ}$ mass from the Planck 2015 catalogue (Planck Collaboration XXVII 2016). Cross-matching the three catalogues and selecting all clusters with $4 \times 10^{14} < M/M_{\odot} < 6 \times 10^{14}$ in both $M_{500}$ and $M_{SZ}$, we finally obtain a sample of 29 ACCEPT clusters (see Table B1, all but five from the Abell catalogue). We additionally excluded Abell 2069, which makes the mass cuts but is a strong outlier in both its entropy and electron density profile as it is a merging cluster. As for the RG clusters, we define CC as having a decreasing profile and a central entropy of at most $40 \text{ keV cm}^2$ at $r = 10 \text{ kpc}$. We note that the overall trends of the ACCEPT CC profiles are in good enough agreement with the more recent analysis of Mantz et al. (2016) of CC systems, most notably the slope of the entropy profiles, for our purposes here.

In Fig. 3, we show the comparison between the mass-selected ACCEPT subsample and stacked profiles from the RHAPSODY sample. For each simulated cluster, we show the mean mass-weighted profile as well as the standard deviation around that mean profile, obtaining one ribbon for each RHAPSODY cluster in Fig. 3. For this analysis, we excluded the fossil cluster RG 572 since it shows an extremely
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Mass-weighted temperature and gas density maps at $z=0$ of the simulated clusters used in this paper at the 4K resolution, i.e. $\langle \rho \rangle = \int \rho^2 \, dl / \int \rho \, dl$ and $\langle T \rangle = \int T \rho \, dl / \int \rho \, dl$. Notable extreme cases are as follows: RG 337 is a merging system and RG 572 is a fossil cluster, RG 474 has twice the mass of the other systems; RG 337 and RG 377 are richer in substructure than average and RG 653 is poorer than average. The images are $8\,h^{-1}\,\text{Mpc}$ wide and the projection depth is $4\,h^{-1}\,\text{Mpc}$ centred on each cluster.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Distribution of selected snapshot redshifts from the RHAPSODY-G clusters (blue, solid) that fall into the mass range $4 < M_{500c}(z)/10^{14} M_{\odot} < 6$ as well as the redshift distribution of the subsample of ACCEPT clusters (red, dashed) for which both the Planck SZ mass and the MCXC X-ray mass fall into the same mass range. The small vertical arrows indicate the median redshift of each sample.</figcaption></figure>

fluctuating entropy and core density, caused by a very active and rapidly growing central hypermassive black hole of ~3 × 10¹¹ $M_{\odot}$ at $z=0$. Note that in the full RHAPSODY sample of more than 100 clusters this was still a completely abnormal case. This cluster is so atypical in its core properties already in pure N-body simulations that we only note here that it is also an outlier in its baryonic core properties.

The number of profiles stacked for each cluster is different and given in Table 1. In all cases, we find that the dispersion in each stack around the mean profile is smaller than the difference between the cool- and non-cool-core profiles, indicating a stable bimodality. This can be clearly seen from the small extent of the ribbons in Fig. 3, which indicate a very small amount of scatter around the mean in the stack of each cluster over time. The CCs are thus a much more significant feature than short-term fluctuations in the profiles. The NCC profiles from ACCEPT and RHAPSODY agree well within their respective scatter, for both the electron density profiles and the entropy profiles. The ACCEPT CC clusters however, show a much weaker CC than the RHAPSODY CC clusters: the observed clusters show only a moderate increase in core electron density inside the innermost ~50 kpc, while the RHAPSODY CCs show a strong drop in entropy and increase in electron density already at scales of ≲150 kpc. It is clear that the CC systems are most likely still undergoing overcooling to some degree despite the central AGN that is efficiently fuelled during the CC phase. The CC/NCC dichotomy is thus a long-lived property of our clusters, consistent with the observational constraints of e.g. McDonald et al. (2013). The dichotomy arose naturally in a larger sample of cosmological cluster simulations and can be explained by differences in the assembly history and nature of major mergers of the clusters (see our analysis in Section 4).

Interestingly, outside the core, the temperature profiles show a somewhat discrepant temperature at large radii, specifically an ~30 per cent difference in X-ray temperature at ≳200 kpc, with much better agreement at smaller radii. Upon closer inspection one notices that a similar but slightly weaker offset also exists in the power-law part of the entropy profile, where the simulated profiles are systematically offset to lower entropy. A similar discrepancy can be seen also, e.g. in fig. 7 of Dubois et al. (2011), independently of the much larger range of subgrid models employed in that reference. This may point towards additional physics missing from these simulations, but is arguably to some degree also due to the

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Comparison of ICM profiles from RHAPSODY-G (red and blue ribbons) with observational data from our mass-matched subset of ACCEPT clusters (points with errorbars). Top: electron density profiles; second from top: entropy profiles; third from top: temperature profiles; and bottom: total mass profiles. The CC clusters in the ACCEPT subset are indicated by blue diamonds and the NCC by red circles, while the CC RHAPSODY-G clusters are shown as blue and the NCC as red shaded regions. The simulated haloes were stacked individually for each cluster while occupying the same mass range as the ACCEPT subset. The shaded ribbons indicate the 1σ scatter in each stack, reflecting time variations in the profiles. The errorbars for the ACCEPT subsample show the scatter in the respective profiles in radial bins, where data was available. We excluded the fossil cluster RG 572 from the RHAPSODY sample here (see text).</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Comparison of the cooling time profiles of CC (blue) and NCC (red) clusters with the respective free-fall-time profiles (grey). The profiles have been stacked individually for each cluster as in Fig. 3. All CC clusters have a cooling time below ~1 Gyr in their cores.</figcaption></figure>

15 per cent higher baryon fraction $f_b = \Omega_b/\Omega_m \simeq 0.18$ that we adopted for these simulations compared to a recent Planck value of $f_b = 0.157$ (Planck Collaboration XVI 2014).

In order to investigate better the nature of the CCs, in Fig. 4, we plot cooling and free-fall-time profiles (cf. also fig. 9 of Li et al. 2015). We see that for the NCC systems, the cooling time is at least a factor of 10 above the free-fall time in the core, while for the CC systems, the two time-scales are much closer to each other. Normally, the central AGN should increase the central cooling time dramatically, but it appears to not prevent the formation of CCs in four of our systems. This is despite a dramatic growth that the central black hole undergoes in the CC cases. For one of them, RG 545, we demonstrate in Section 7.1 that the existence of the CC does not depend on the details of the AGN feedback parameters. We are thus led to believe that in our simulations the CC/NCC dichotomy has a cosmological origin, and will investigate this aspect further in Section 4.

### 3.3 Hydrostatic mass and deviations

The cluster is stabilized against gravitational forces by thermal pressure as well as bulk and random flow motions (e.g. Nelson et al. 2014). We next investigate to what degree the thermal pressure alone counteracts gravity in our simulations. In adaptive mesh simulations with a Lagrangian refinement strategy (which is however indispensable if galaxy formation is included), the numerical Reynolds number is relatively low so that no turbulent cascade develops and the energy is quickly dissipated into thermal energy. In contrast to uniform high-resolution simulations with fully developed turbulence (Miniati & Beresnyak 2015), we can focus on thermal effects here inside the outer infall region. We follow Nelson et al. (2014) and define the thermal (or hydrostatic) mass as

$$M_{\text{therm}}(r) = \frac{-r^2}{G \langle \rho \rangle} \frac{\partial \langle P \rangle}{\partial r}, \qquad (7)$$

where the angle brackets denote mass-weighted averages over spherical shells at a given radius $r$. In Fig. 5, we show the ratio $M_{\text{therm}}(r)/M_{\text{grav}}(r)$, where $M_{\text{grav}}(r) = M_{\text{tot}}(<r)$ is the total mass enclosed in radius $r$. All clusters are closest to hydrostatic equilibrium

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Ratio of the hydrostatic and the gravitating mass as a function of radius. Outside of the core and well inside the virial radius, the clusters are very close to hydrostatic equilibrium, with large deviations at radii $\gtrsim 1$ Mpc. There is no obvious difference between the CC systems (blue) and the NCC systems (red). The grey shaded band indicates $\pm 20$ per cent deviation from hydrostatic equilibrium.</figcaption></figure>

at radii of ~300 kpc with a mean deviation of at most 20 per cent inside of ~600 kpc. At large radii, the variance is however substantial, most likely due to incomplete thermalization of the region, which is perturbed by accreting gas. Most importantly for our analysis, we find no significant difference in deviations from hydrostatic equilibrium over radii 100 kpc $\lesssim r \lesssim 1$ Mpc between CC and NCC systems. There is however an indication of systematic differences at $\gtrsim 1$ Mpc, as well as, obviously, at scales of the core. The scale relevant for estimating the hydrostatic mass bias in scaling relations is of course at $R_{500}$. At that scale, we find a median bias $M_{\text{therm}}(r)/M_{\text{grav}}(r)$ of ~1.2 for the CC clusters and of ~1.3 for the NCC clusters indicating the more perturbed nature of the latter. The overestimation of hydrostatic masses, a result contrary to most previous work, may result from our choice of mass-weighted estimators for the pressure and density in equation (7), which is more susceptible to clumping. While a more detailed and robust analysis would be necessary to make a meaningful statement for the bias of observed cluster masses, we here only want to draw the qualitative conclusion that intermediate radii between 200 and 500 kpc are close to hydrostatic equilibrium in our simulations and thus seemingly unaffected by both feedback from the central AGN as well as from accreting gas.

### 3.4 Gas depletion profiles

A crucial observable reflecting the degree to which collisional matter follows the total mass distribution of galaxy clusters is given by the gas depletion profiles. Especially in the cluster cores, cooling and AGN feedback will affect both the gas and the dark matter distribution. In particular, if the central black hole accretes gas quasi-periodically, AGN feedback can transform the cuspy dark matter profiles into cored profiles (Martizzi, Teyssier & Moore 2013, consistent with what we find in the RHAPSODY-G simulations). In Fig. 6, we show the ratio of enclosed gas mass to total mass in units of the universal baryon fraction. We split our sample of clusters into NCC (top) and CC (bottom) systems. We compare our results to the observational relations of Mantz et al. (2014), as indicated by the blue ribbons. We note that the sample of Mantz et al. (2014) specifically only includes relaxed systems, which means that it is clearly biased towards CC systems that have higher baryon fractions in the cores. In addition, we show in the lower part of each panel
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Comparison of cumulative gas (orange ribbons) and baryon fraction (green ribbons) profiles for the same stacked clusters as in Fig. 3, i.e. for the mass range $4 < M_{500c}/10^{14} M_{\odot} < 6$ with the observational results from Mantz et al. (2014). We split again into CC (bottom) and NCC clusters (top), and show stacked results for each simulated cluster individually. The shaded ribbons indicate the $1\sigma$ in the stack of each cluster. In the CC panel, we have indicated the fossil cluster RG 572 with a hatched ribbon.</figcaption></figure>

the contribution of the full baryonic component, i.e. stars and gas, to the total enclosed mass fraction.

The NCC systems follow the observational relation for gas depletion reasonably well. The fossil cluster and the CC systems however show a very large central gas fraction, in many cases even above the universal gas fraction (similar to the results of Burns et al. 2008). This result is clearly in tension with the results of Mantz et al. (2014), and reflects the high electron densities and very low central entropies we have seen for the CC systems above. The conclusion must be that in the case of the CC systems, the AGN feedback is not able to stabilize the core at levels consistent with X-ray observations. The exciting possibility is that other forms of non-thermal feedback (or processes) must be plausibly involved in order to bring these results in line with observations. We note that our results for NCC systems are consistent with published results from SPH simulations (e.g. Battaglia et al. 2013; Planelles et al. 2013; Sembolini et al. 2013) at scales of $0.1 R_{vir}$, but are somewhat higher at larger radii. These observations are clearly aggravated by the fact that the sample of Mantz et al. (2014) is biased towards CCs. From the total baryon fraction profiles, we see that it stays very close to the universal value at all radii $\gtrsim 0.2 R_{200}$, again indicating that

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Assembly history of RG 348 and RG 545. The first is an NCC cluster at $z=0$, the latter is a CC cluster. The top panel shows the cooling time evolution in the core for the two systems, the bottom panel the mass accretion history. Both undergo a major merger between $a \sim 0.5$ and $a \sim 0.6$. While the core of RG 348 is substantially heated during the merger, the effect on RG 545 is not as dramatic. The vertical lines indicate the times for which we show images illustrating the evolution of RG 545 and RG 348 in Figs 8 and 9, respectively.</figcaption></figure>

AGN feedback, while being effective at reducing BCG masses to realistic values (Martizzi et al. 2012a, and see Section 5 below), is not preventing a significant reduction of baryons in the cluster core in our simulations, particularly so for the more relaxed CC systems.

Less visible but as significant is the discrepancy of the gas fraction at $\sim R_{200}/3$, where the profiles are slightly but systematically high with very little scatter. As we see below in Section 7.2, the thermal AGN feedback does not reach these large radii in our simulations, and it is not possible to tune the energy injection parameter $\Delta T$ to deplete the CCs efficiently.

Despite these discrepancies, the use of the gas fraction at $R_{2500}$ as a robust cosmological observable (Allen et al. 2008; Mantz et al. 2014) is strongly supported by our simulations, since $R_{2500}$ appears to be outside the reach of the AGN and shows a virtually unbiased thermal mass (cf. Section 3.3).

## 4 THE ORIGIN OF THE CC/NCC DICHOTOMY

We next inspect the origin of the CC/NCC dichotomy in our simulations. We find that at early times, $z \sim 0.6$ and $z \sim 0.8$, respectively, both RG 348 and RG 545 are CC clusters and have a very similar assembly history, which we show in Fig. 7. Both undergo a major merger with a similar mass ratio at similar times. However, the cooling time of the core (top panel of Fig. 7) rises dramatically more in the case of RG 348 to a value of about 5 Gyr after an initial higher spike. Cluster RG 545 only experiences an increase in cooling time to $\lesssim 2.5$ Gyr. We define the core here as the gas enclosed in an overdensity of $8 \times 10^4$ and evaluate the cooling time at that radius (we find that the free-fall time reduces only slightly over the same time-scale so that roughly $t_{cool}/t_{ff} \propto t_{cool}$). After the merger, both clusters have a more quiescent merger history and both cores start to cool again. The cooling of RG 545 is much more rapid, plausibly since it is quickly experiencing runaway cooling, while RG 348 cools at a slower rate – the time to $z=0$ is still 7 Gyr, shorter than its cooling time – most likely reflecting AGN heating that is able to keep the core hot.
<figure><img src="image_10.png" /><figcaption>**Figure 8.** CC surviving after a major merger (RG 545): we show the entropy (left-hand column) and density (right-hand column). The first row is at $z = 0.82$, the next rows are 322 and 466 Myr later. The top-left panel clearly shows the two CCs of the systems before the merger, but due to enough angular momentum, the cores do not collide (unlike in the case of RG 348, see Fig. 9), and the CC survives. Solid-line and dashed-line circles indicate the cores of the main cluster and the merging cluster, respectively.</figcaption></figure>

After inspecting the time evolution of the systems, the difference in how the mergers proceed is obvious. In the case of RG 545, see Fig. 8, the cores do not collide since the merger has large angular momentum, so that the CC is perturbed but not destroyed and is just sloshing in the centre of the halo. Quite different is the merger that RG 348 experiences: here, the two CCs of the progenitor clusters collide and dissipate within ~100 Myr. We show their time evolution over ~450 Myr during the merger in Fig. 9: the left-hand column shows the minimum entropy along the comoving $1h^{-1}$ Mpc image depth of X-ray-emitting gas with a temperature above 0.1 keV. We also show the density of gas above 0.1 keV, as well as the temperature (in linear scale) and the dark matter density. For RG 545, we only show entropy and density since the merger is less spectacular. From the sequence of RG 348, one sees that the cores pass right through each other (in the moment shown in the second row); after the collision, the CC of the main cluster halo is destroyed, while the one of the less massive progenitor survives weakened. In the last time shown and subsequently, the surviving core mixes quickly with the hot cluster gas due to ram pressure and Kelvin-Helmholtz instabilities. We also note that, as expected, we detect a weak displacement between the BCG as well as the dark matter halo centre and the gas centre after the collision, consistent also with the observational findings of e.g. Semler et al. (2012). Coincidentally, as we submitted this article, Rasia et al. (2015) reported a very similar result from their SPH simulations, which include AGN feedback along with a novel thermal diffusion scheme. While the dichotomy between CC and NCC systems appears to be less pronounced in their simulation, Rasia et al. (2015) are able to reproduce more realistic CCs. They find that the AGN heating and cooling yields a delicate balance in CC systems, which can then be perturbed by cluster mergers at all redshifts leading to a transformation between the NCC and CC states.

The finding of a connection between CC/NCC and mergers in early simulations put forward by Burns et al. (2008) and Planelles & Quilis (2009) thus should be complemented by this new finding. Not only slowly accreting systems can have CCs; the amount of angular momentum in major mergers (i.e. whether the collision is head-on in such a way that the cores do pass through each other) appears to be of high importance as well. This is also consistent with the findings from idealized non-cosmological cluster merger simulations by Poole et al. (2008). Major mergers are thus a necessary but not sufficient condition for the destruction of CCs. The survival of RG 545 as a CC system clearly shows that the cores are not easily disrupted, in some cases, even by a major merger. For all clusters in our sample, it is true that, very similar to the two examples in Fig. 7, the CC clusters never had their cooling time increased beyond at most 20 times the core free-fall time. This is roughly consistent with Meece, O'Shea & Voit (2015) who find that a ratio of at most 10 is needed in order to undergo a significant thermal instability. Since other heat sources are absent, every cluster will slowly cool again after a major merger event (and basically every cluster will experience one or more major mergers throughout its history); the sharp distinction between the CC and NCC systems in terms of their profiles then arises simply due to how much a merger perturbs the core – determined by the merger’s angular momentum – and whether or not rapid cooling can set in subsequently. In the case of RG 348, the CC has been dissipated in no more than a few 100 Myr with a subsequent cooling time that is significantly higher. The disparity in these two time-scales makes this a quick transition, allowing the bimodality apparent in the cluster population to be sustained.

# 5 PROPERTIES OF THE GALAXIES IN AND AROUND THE SIMULATED CLUSTERS

We next discuss the properties of the galaxy population in and around our simulated clusters. Since the high-resolution region is an $8h^{-1}$ Mpc sphere around each cluster at $z=0$, we have a sizeable sample of lower mass haloes ranging from lower mass clusters to group to galaxy scales in the outskirts of our massive clusters. This facilitates a statistical comparison with the more field-dominated observational samples. Nevertheless all galaxies in our simulations do live in a high-density environment, which should be kept in mind when comparing these results to observations.

## 5.1 Comparison with abundance-matching results

Apart from the hot intracluster plasma, galaxy clusters harbour a large population of satellite galaxies. A general shortcoming of galaxy cluster simulations is that in these simulations it is computationally not affordable to resolve the length scales relevant to important aspects of galaxy formation (on say 100 pc and smaller). The hope is that at least the most massive galaxies are resolved well enough. A bona-fide test of the realism of our simulations is thus to investigate whether the central cluster galaxies indeed have stellar masses and other properties that are consistent with observations. Remarkably, despite the fairly poor resolution at galaxy mass halo scales, we find that the star formation efficiency and SN feedback as dictated by the subgrid models is in fact able to reproduce accurate stellar masses across all halo masses that we resolve with $\gtrsim 1000$ particles in our simulations. In Wu et al. (2015), based on the same simulations discussed here, we find stellar mass fractions of $\sim 10$ per cent of the baryon fraction, consistent with observational constraints.

In Fig. 10, we show the comparison of the stellar masses of our galaxies (satellite and central) as a function of their halo mass with
<figure><img src="image_11.png" /><figcaption>**Figure 9.** CC destroyed in a major merger (RG 348). We show the minimum entropy along the line-of-sight (left-hand column) and the density (middle left) of gas with a temperature of at least 0.1 keV, as well as the density-weighted temperature in linear scale (middle right) and the dark matter density (right-hand column). The first row is at $z = 0.61$, the next rows are 171, 330 and 448 Myr later. The top-left panel clearly shows the two CCs of the systems before the merger. During the merger, the cores undergo a strong shock (clearly visible in second temperature panel from top). After the collision, the CCs are dissipated and undergo mixing with higher entropy gas. Solid-line and dashed-line circles indicate the cores of the main cluster and the merging cluster, respectively. The many low-entropy blobs are satellite galaxies that do not carry much mass with them. The depth of the images is comoving $1h^{-1}$ Mpc centred on the more massive cluster.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Comparison between the $M_{*}-M_{500c}$ relation for galaxies from all 10 *RHAPSODY-G* haloes at $z=0$ and abundance-matching constraints. The orange dots represent central galaxies with the mean stellar mass in bins of halo mass shown by bold black circles with errorbars indicating the $1\sigma$ scatter, while the blue dots represent the satellite population, respectively. The solid and dashed black lines show the abundance-matching results from Kravtsov, Vikhlinin & Meshscheryakov (2014), while the dotted red line shows the abundance-matching results from Behroozi et al. (2013c) using the stellar mass functions from Bernardi et al. (2013, remapped from $M_{200c}$ to $M_{500c}$ using the mean relation from the simulation itself). Simulated satellites appear at higher stellar masses as a result of a reduction of their halo mass due to tidal stripping. The drop of stellar masses below the Behroozi relation for masses $\lesssim 10^{12} M_{\odot}$ is a resolution effect (cf. Fig. 11).</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>**Figure 11.** As Fig. 10, but a comparison between the $M_{*}-M_{500c}$ relation for galaxies from a higher resolution simulation of cluster RG 653 (simulated at eight times better mass and twice better spatial resolution than the full sample) at $z=0$ and abundance-matching constraints.</figcaption></figure>

results obtained using the abundance-matching technique (cf. e.g. Conroy, Wechsler & Kravtsov 2006; Moster et al. 2010; Behroozi et al. 2013c). The very massive systems that we consider here occupy however the very tail of the stellar mass functions, and in fact, the older abundance-matching relations (e.g. Moster et al. 2010; Behroozi et al. 2013c) do not capture very well the scaling at the highest masses due to the way stellar mass was counted in the stellar mass functions used. More accurate results have been obtained by Kravtsov et al. (2014) taking explicit care of the intracluster light component; these are fully consistent with the relations of Behroozi et al. (2013c) when updated with the stellar mass functions of Bernardi et al. (2013). We show a comparison of these relations with the stellar masses from our simulations measured using a range of techniques. In particular, we have produced mock V-band images and count the stellar mass inside of isophotes where the surface brightness is above 25 and 24 mag arcsec$^2$, as well as using a simple spherical density threshold of $2.5 \times 10^6 \text{ M}_\odot \text{ kpc}^{-3}$ (panels from left to right in Fig. 10).

We see that, while we largely overproduce stellar mass compared to the original Behroozi et al. (2013c) result, our simulations are consistent with more recent estimates for massive systems. The stellar masses in massive haloes above $M_{500c} \sim 2 \times 10^{12} \text{ M}_\odot$ are still somewhat high, but not dramatically so, and details depend on the exact definition of how stellar mass of central group and cluster galaxies should be measured. Notably, the stellar masses at galactic halo scales seem to undershoot the observational relations. However, this is a resolution effect. In Fig. 11, we show the respective plot for the galaxies in and around cluster RG 653 simulated at twice better spatial and eight times better mass resolution. With higher resolution, the full range of halo masses from $10^{11}$ to $10^{15} \text{ M}_\odot$ is in good agreement, with the slight overproduction of stars in the most massive haloes remaining.

## 5.2 Star formation rates

We next want to investigate the reason for the somewhat high stellar masses in the most massive haloes in our simulation. To this end, we

<figure><img src="image_14.png" /><figcaption>**Figure 12.** SFR as a function of stellar mass for central (orange) and satellite (blue) galaxies from all 10 RHAPSODY-G haloes at $z=0$. Galaxies with no detectable star formation in the previous 500 Myr have been assigned random values between $10^{-3}$ and $10^{-2} \text{ M}_\odot \text{ yr}^{-1}$. For comparison, the main sequence of star-forming galaxies in SDSS from Woo et al. (2013, solid black line) and the division between star-forming and quenched galaxies (black dotted) along with the O II-derived SFRs for the BCG samples of Liu, Mao & Meng (2012) and the UV estimates of McDonald et al. (2015) are shown (note that for the latter, data points with only long lower errorbars are upper limits).</figcaption></figure>

show each of our simulated galaxies (for all clusters) in the stellar-mass versus SFR plane. We estimate our SFRs as the mass of star particles in each simulated galaxy that have an age of less than 100 Myr divided by the 100 Myr. Galaxies that do not have any star particle fulfilling this condition get assigned a random SFR between $10^{-2}$ and $10^{-3} \text{ M}_\odot \text{ yr}^{-1}$. We perform this step in order to allow a visual estimate of the fraction of quenched galaxies in Fig. 12.

From Fig. 12, we see that our galaxies with significant star formation follow a narrow star-forming sequence that is consistent in normalization and slope with the sequence of star-forming galaxies in the Sloan Digital Sky Survey (SDSS) as measured by (Woo et al. 2013, solid black line). In addition, we show the SFRs determined by Liu et al. (2012) for a sample of central cluster galaxies from SDSS, and those of McDonald et al. (2015) for BCGs of clusters detected with the South Pole Telescope (SPT). When comparing with the limit below which Woo et al. (2013) would classify galaxies as ‘quenched’ (black dashed line), we see that many more massive galaxies in this simulation are not quenched – the SDSS quenched fraction exceeds 50 per cent above $M_* \sim 10^{10.5}$ in Woo et al. (2013, private communication). While the SDSS results are averaged over all field galaxies, it may be that some biases exist in these most dense environments around the most massive clusters. However, it seems more plausible that these simulations either miss an additional physical mechanism or lack the resolution to reproduce a realistic quenched fraction at the intermediate masses. While the highest SFRs we observe are roughly consistent with Liu et al. (2012) and perfectly consistent with McDonald et al. (2015, note that the data points with one-sided lower errorbar are upper limits only), it is clear that quenching of intermediate (and maybe even high) mass galaxies is inefficient in our simulations. This is consistent with the somewhat high stellar masses compared to the
abundance-matching results. Finding the mechanism(s) necessary to suppress star formation in massive galaxies is among the most pressing tasks of galaxy formation theory today and it is not surprising that our simulations are not performing significantly better than better-resolved simulations dedicated to study galaxy formation.

# 6 IMPLICATIONS FOR CLUSTER COSMOLOGY: EVOLUTION ALONG SCALING RELATIONS

The use of cluster counts as a cosmological probe exploits population scaling relations that link total system mass to observable properties of the galaxies or hot plasma. In this section, we compare the evolution of the simulated clusters along a range of scaling relations with observational results from X-ray and microwave data and some other published results from simulations. We focus on the SZ effect as well as on X-ray luminosity and temperature. In the absence of a large statistical sample of clusters, tracking a few clusters over time allows us to relate the role of evolutionary processes to their impact on observables.

Since systems in the simulated sample have roughly the same mass at $z=0$, their evolutionary tracks can be used to assess the contribution of their assembly history to scatter in the mass-proxies. In order to study the evolution of the scatter or additional biases relevant to the full population, however, larger samples would be required. We note that, apart from RG572, all clusters grow above $M_{500} \gtrsim 5 \times 10^{13} M_{\odot}$ in the range $1.5 \lesssim z \lesssim 2$, and above $2 \times 10^{14} M_{\odot}$ in the range $0.75 \lesssim z \lesssim 1.25$. While we compare our scalings to observations, we caution that this type of exercise is non-ideal in that (i) we employ true, three-dimensional spherical masses in simulations while observational estimates may be biased with respect to these values; (ii) similarly, the intrinsic properties of the simulations are measured directly in the simulations rather than derived from models applied to mock observations; (iii) the statistics of the simulated sample are not necessarily representative of a broader mass-complete sample; and (iv) the correlated cluster evolution itself might lead to a departure from the scaling relations compared to studies that are performed at fixed redshift using larger samples. The purpose of this exercise is to identify areas of agreement and discrepancy, which in turn should provide insights into future adjustments to the astrophysical feedback model.

## 6.1 Sunyaev-Zeldovich masses

We first compare the Compton-$\gamma$ parameters integrated over our simulated clusters, i.e. calculated as

$$Y_{500} = \int_{|x|<R_{500}} \sigma_T \frac{kT}{m_e c^2} n_e d^3x, \qquad (8)$$

with SZ masses and $Y_{5R500}$ measurements from the Planck high-frequency instrument (HFI) union catalogue 2015 (R2.08), see Planck Collaboration XXXII (2015). We have converted the Planck $Y_{5R500}$ values to $Y_{500}$ using the ratio of 1.814 between the measurement at $5R_{500}$ and at $1R_{500}$ quoted in Melin et al. (2011) and assumed the Planck cosmological parameters when converting $Y_{5R500}$ to units of kpc$^2$. We show the results of this comparison in Fig. 13, along with the unbiased (i.e. setting the hydrostatic mass bias to $b=0$) Planck 2013 mean baseline relation

$$E^{-2/3}(z) \left[ \frac{Y_{500}}{\text{kpc}^2} \right] = 10^{1.81 \pm 0.02} \left[ \frac{M_{500}}{6 \times 10^{14} M_{\odot}} \right]^{1.79 \pm 0.08}, \quad (9)$$

as quoted in Planck Collaboration XX (2014).

<figure><img src="image_15.png" /><figcaption>Figure 13. Comparison of the Planck 2015 SZ cluster sample with the evolutionary tracks of the RG clusters along the M-Y relation. The Planck data is taken from the 2015 union catalogue with rescaled $Y_{500} = Y_{5R500}/1.814$ according to Melin et al. (2011).</figcaption></figure>

Our simulations show good agreement with the Planck 2013 baseline for all systems over the entire range from $\sim 5 \times 10^{13}$ to $10^{15} M_{\odot}$. Despite the large variance between the clusters' assembly histories, we find very little scatter, $\sim 0.2$ dex around the Planck 2013 mean relation (consistent with other simulations, e.g. Sembolini et al. 2013), with only the fossil cluster RG 572 being somewhat of a significant outlier (though still insignificant given the scatter in the Planck M-Y data). In fact, it is remarkable that RG 572 does not stand out. In particular, the scatter is substantially smaller than the scatter around the M-Y scaling relation in the Planck 2015 SZ catalogue. The large observed scatter is mostly due to complications from projection and low angular resolution, which requires joint fitting of $Y_{5R500}$ and cluster size $\theta_{500}$. The large and markedly asymmetric scatter in the Planck M-Y data is an artefact of our simplistic rescaling from $Y_{5R500}$ to $Y_{500}$ with a constant factor and appears in rather stark contrast to the simulation results. Investigating the simulation prediction for $Y_{5R500}$ would require a matched filter approach to reduce interloper effects at larger radii, as employed by the Planck team, and possibly full light-cone projections (e.g. White, Hernquist & Springel 2002), which are beyond the scope of this work.

## 6.2 X-ray observables

In Fig. 14, we show a comparison of the X-ray luminosity of the RG clusters during their evolution as a function of their mass with various scaling relations from the literature as well as a direct comparison with clusters from the MCXC sample (Piffaretti et al. 2011) using the Planck SZ masses and the MCXC luminosity (since otherwise one simply reproduces the scaling relation used to estimate the masses, which is very close to Arnaud et al. 2010). We measure the X-ray luminosity by computing the emissivity per cell using the APEC tabulated emission model and exclude the central core region inside $0.15R_{500}$, i.e.

$$L_X = \int_{0.15 < |x| / R_{500} < 1} kT \epsilon(T, Z) n_e n_H d^3x, \quad (10)$$
<figure><img src="image_16.png" /><figcaption>**Figure 14.** Comparison of the $M_{500}-L_X$ scaling relation with observations and other simulations. The RG clusters are again shown as evolutionary tracks, the luminosity is the integrated 0.1–2.4 keV core-excised luminosity estimated from APEC spectra. We compare this with various other data. The data by Pratt et al. (2009), solid dark blue line, is their fit for the core-excised luminosity in the same band. For the clusters from the MCXC sample, we have used the 0.1–2.4 keV luminosities from MCXC but the respective Planck SZ masses here. The dash–dotted line represents the fit from Arnaud et al. (2010), which is almost identical to the scaling relation used in MCXC, and the dashed line is the best-fitting self-similar solution (i.e. the $\propto M_{500}^{4/3}$ fit to the MCXC/Planck joint data). The dotted line shows the best fit from Biffi et al. (2014), who however use the wider Chandra band. The light blue circles with errorbars show the results from Le Brun et al. (2014) for the 0.1–2.4 keV band in terms of $M_{500}$ including X-ray bias, the squares show their results if instead the true mass $M_{500}$ is used (Le Brun et al., private communication).</figcaption></figure>

where $n_e(x)$ and $n_H(x)$ are the electron and hydrogen number density, $T(x)$ is the electron temperature and $\epsilon(T, Z)$ is the (temperature- and metallicity-dependent) emissivity taken from the APEC tables and integrated over the same energy range 0.1–2.4 keV as the MCXC clusters. The exclusion of the core is necessary as it introduces a strong variability of the X-ray luminosity due to AGN events where the central density fluctuates and unrealistically high luminosities can be obtained. The fit given by Pratt et al. (2009, solid blue line), is using the same core exclusion radii as for our simulation data. The other samples that we compare against however, do include the core, so that the comparisons have to be taken with a grain of salt. For comparison, we also show the results of Le Brun et al. (2014), who use the same AGN model in SPH simulations. We show their results in terms of the true mass $M_{500}$ (as for our data), and for a more ‘observational’ $M_{500}$ that includes X-ray bias (as published in their paper). We note that we show their data for the same 0.1–2.4 keV band and in terms of the true mass $M_{500}$ that they have kindly made available to us. While their data appears thus low in this plot compared to the observations, this discrepancy can be alleviated by taking into account hydrostatic mass bias as done in their actual paper (not shown here). For completeness, we also show the best-fitting relation from Arnaud et al. (2010), which is basically identical to the scaling relation used for the $M_{500}$ measurements in MCXC, as well as the best fit from

Biffi et al. (2014). The latter authors use the Chandra bolometric band (0.3–12 keV), which explains the higher X-ray luminosity compared to the MCXC data. In comparison, our clusters evolve along a scaling relation that is somewhat shallower than both Arnaud et al. (2010) and Le Brun et al. (2014), but consistent with Pratt et al. (2009) and (Biffi et al. 2014, who however used the larger spectral window) but slightly steeper than a simple self-similar scaling $\propto M^{4/3}$ (shown as a green dashed line in the figure and whose normalization has been obtained by fitting to the MCXC data). In fact, the best-fitting slope in our case is closer to 1.2, but we do not want to make a more quantitative statement at this stage since the comparison of evolving small (autocorrelated) samples should be taken with some caution for making rigorous predictions. It suffices to observe here that the simulation X-ray luminosities are consistently higher (by about a factor of 2 around $10^{14} M_{\odot}$ and at most $\sim 20$ per cent at $10^{15} M_{\odot}$) than the both the relation given by Pratt et al. (2009) and the MCXC luminosities, with a scaling relation that is slightly steeper than the self-similar expectation of 4/3. There is a weak indication that the normalization of the evolutionary track of clusters is persistent, i.e. RG 211 tends to be lower than RG 361, which is lower than RG 348 for a range of masses. This is plausibly indicative of a connection between X-ray luminosity and assembly history at fixed mass. We will investigate the origin of scatter in future work. Once again, it is remarkable that RG 572 does not stand out dramatically, its core-excised X-ray luminosity is only slightly lower than the rest of the sample. Since the use of Planck $M_{500}$ shifts the MCXC clusters away from the Arnaud et al. (2010) relation, we expect that hydrostatic mass bias will play an important role when comparing the RG clusters to observations. From our analysis in Section 3.3, it is however clear that the 20–30 per cent bias we estimated there cannot account for the offset between the RG scaling relation and observations.

To conclude our comparison with X-ray scaling relations, in Fig. 15, we show the mass–temperature relation for a subset of the MCXC clusters for which we could take the X-ray temperatures from the ACCEPT cluster catalogue along with the data from Vikhlinin et al. (2009), which agrees perfectly with lower scatter with the MCXC sample we show. Due to the known biases in estimating temperatures from simulated clusters that reflect temperatures measured from X-ray spectroscopy (see e.g. Biffi et al. 2014, for in-depth comparisons), we show a range of differently weighted temperature estimates in the different panels of the figure. First, we compare with the spectroscopic-like temperature as introduced by Mazzotta et al. (2004) as follows:

$$T_{\text{sl}} = \frac{\int n^2 T^{\alpha-1/2} dV}{\int n^2 T^{\alpha-3/2} dV}, \quad \text{with } \alpha = 0.75, \qquad (11)$$

as a fit to spectroscopic temperature estimates with Chandra or XMM–Newton. We additionally do not include cells with a temperature below 0.5 keV in this average, which is necessary in the case of cooling to avoid the inclusion of non-X-ray-emitting gas. Furthermore, we again exclude the core inside $0.15 R_{500}$ in order to avoid even larger fluctuations in the estimated cluster temperature due to the influence of the AGN model on this region. This is also common practice in observations to reduce scatter in the scaling relations (e.g. Mantz et al. 2010). We found that a temperature estimated by computing the mean energy of all photons emitted using an APEC model (see discussion above) above 0.1 keV leads to basically identical results within the scope of this first analysis, so that we do not consider these ‘mock’ X-ray observations in this paper further.
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Comparison of the ACCEPT/MCXC clusters $M500-T_X$ relation with the evolutionary tracks of the RG clusters. From left to right, we show comparisons with different estimates of the gas temperature in the simulations: (left) the spectroscopic-like temperature, with core excluded, (center) volume-weighted temperature and (right) the mass-weighted temperature. The grey dots are a combination of X-ray masses from the MCXC catalogue with temperature estimates taken from the ACCEPT catalogue, and the black long-dashed line is the best-fitting power law to this data, while the short-dashed line indicates a power-law fit to RG 474, which covers the largest range in mass evolution. The open squares represent the measurements of Vikhlinin et al. (2009).</figcaption></figure>

The other two temperature estimates we consider are a simple mass and a volume-weighted average of all cells with no core excluded and no cut in minimum gas temperature, so that they reflect the true mean temperatures, not including possible observational biases.

We note that the spectroscopic-like (as an APEC emissivity weighted temperature not shown here) temperature has a much larger variance compared to the volume- and mass-weighted averages. The variance substantially increases if the core is also included due to the $n^2$-dependence of the emissivity. In all cases, and thus plausibly not affected by potential biases compared to observational results, the $M-T$ scaling that we observe is slightly steeper than the one exhibited by the MCXC/ACCEPT data. In particular, we find that the cluster temperature tends to be somewhat low for masses above $10^{14} M_\odot$. The discrepancy, while systematic, is not dramatic. It is consistent with our comparison of the temperature profiles with the ACCEPT data, in which we found a similar temperature discrepancy at large radii.

To summarize, in no case did we find a very strong dependence of any cluster scaling relation on the assembly history when the core is excluded. A weak dependence might be possible for the $M-L_X$ relation in our data. This result is consistent with earlier analysis in this direction, by e.g. Rasia et al. (2011), who show how explicit minor and major mergers evolve along the overall population scaling relation, and Jeltema et al. (2008), who however found a dependence on the dynamical state once the scaling relations are not compared against the true mass but against, e.g. the hydrostatic mass. This is very plausibly the case for our simulations as well, since the hydrostatic mass bias can be large close to $R_{500c}$ (see our analysis in Section 3.3). It is thus not entirely clear if and how the ACCEPT comparison may be biased by hydrostatic mass errors that could alleviate also the discrepancies we observed between the profiles in Section 3.2.

# 7 DISCUSSION OF THE NUMERICAL MODELLING APPROACHES

We have established already in the previous sections that our results are numerically converged by comparing the 10 runs at 4K

resolution to the higher resolution 8K run of RG 653. Next, we investigate the dependence of our results on particular choices of the AGN feedback model parameters and discuss discrepancies with published results based on Lagrangian methods.

## 7.1 AGN model dependence: X-ray properties

In this section, we investigate the robustness of our results to changes in the thermal AGN feedback model. In particular, we want to see whether the cold cores we find are stable to such changes. We furthermore will investigate the degree of change in the scaling relations that can result from such astrophysical changes beyond the resolution of the simulation in the subgrid model.

While it is clear that the very short cooling times in the centres of massive clusters necessitate a mechanism to prevent efficient cooling in their cores in order to reproduce both the observed gas fractions and realistic galaxy masses, what this mechanism is, or how it should be modelled on the scales accessible to three-dimensional cosmological simulations, is less clear. Virtually all current state-of-the-art cosmological simulations that tap into this mass regime invoke some form of AGN feedback to resolve this problem. As shown in Fig. 6, the AGN model we adopted fails however to reduce the gas fraction at radii $\gtrsim 0.5 R_{vir}$, where our simulations predict somewhat high baryon fractions compared to e.g. Mantz et al. (2014). If AGNs however have such a dramatic effect on the entire ICM, most of cluster cosmology would have to rely on a very tight relation between cluster mass and AGN feedback energy. The analysis of Le Brun et al. (2014) e.g. shows that the entire baryon gas profile can be varied by tuning the energy accumulation threshold $\Delta T$ of the feedback model of Booth & Schaye (2009). While this is possible, it should be expected that AGN feedback could be a significant source of scatter, relating processes at the pc scale to the Mpc scale of the ICM. In this section, we thus repeat this analysis and investigate how sensitive our results are to variations of the energy accumulation parameter $\Delta T$.

As a benchmark, we show in Appendix A also the profiles obtained for clusters RG 348 and RG 545 in non-radiative runs (see Fig. A1). A small difference between the CC and NCC system can already be seen in these parameter-free simulations and without
<figure><img src="image_18.png" /><figcaption>**Figure 17.** Impact of the AGN model on ICM profiles: as the top three panels of Fig. 3, but only for the CC cluster RG 545 using two different implementations and a range of parameters of the thermal AGN injection threshold. The phenomenological model has the most significant impact by raising the entropy also at small radii, changing the electron density to a power law in the core, and flattening the temperature profile. The lowest injection threshold $\Delta T = 10^6$ K produces frequent outward travelling blast waves that show up in the entropy and temperature profile but only have a modest effect on the density. In no case radii outside $\sim 100$ kpc are affected, however. For comparison, we also show the profiles of the ACCEPT CC subset discussed above (blue diamonds).</figcaption></figure>

Mantz et al. (2014). We thus also want to investigate the impact of the AGN model parameters on the gas depletion profiles. In Fig. 18, we show essentially the same plot as in Fig. 6, varying again the injection threshold parameter $\Delta T$. While we saw that the lowest injection threshold $\Delta T = 10^6$ K was able to impact the entropy and density profiles, it also reduces the central gas fraction, but only very slightly so, and only at radii $\lesssim R_{vir}/10$. Similarly, the larger injection regions of the phenomenological model reduce the central gas fraction below the universal value. Once again, in no case could we affect regions outside the core leaving the gas fraction inconsistent with observations at $\sim 0.2 R_{vir}$. The energy available from the thermal blast wave, even after accumulating energy to heat to higher temperatures is completely insufficient to reduce the gas fraction outside the core to the $f_{gas} \sim 0.5 \Omega_b / \Omega_m$ observed by Mantz et al. (2014). Plausibly, thermal conduction or kinetic feedback might play a role here to distribute energy better towards the outskirts.

<figure><img src="image_19.png" /><figcaption>**Figure 18.** Impact of AGN model on gas fraction: as Fig. 6, but only for the CC cluster RG 545 using two different implementations and a range of parameters of the thermal blast-wave AGN model. The phenomenological model has the most significant impact to keep the gas fraction below the universal baryon fraction also at small radii, but in all cases the gas fraction at $\gtrsim 0.1 R_{vir}$ is high compared to the constraints from Mantz et al. (2014). In contrast to the NCC case, the thermal blast wave here is insufficient to substantially reduce the gas fraction in the inner parts of the halo.</figcaption></figure>

## 7.3 A note on fundamental differences between AMR and SPH

Various authors using SPH simulations have argued that the parameter $\Delta T$ has a strong impact on the physical properties of both the galaxies and the intracluster gas. As we have discussed above, this is in contrast to our own findings. Le Brun et al. (2014) find changes in the gas density, Pike et al. (2014) in the gas pressure out to $R_{500}$. Furthermore, Le Brun et al. (2014) found that $\Delta T = 10^8$ K allowed them to make cluster properties compatible with a range of observables. We discussed in more detail the impact of the two flavours of energy injection as well as the energy accumulation threshold on the various observables and various physical properties in Section 7.1 above. As we demonstrated there, in our Eulerian AMR simulations, the particular choice of $\Delta T$ and even the size of the injection region does not affect larger radii. It is thus plausible that thermal feedback has a different impact in Lagrangian and Eulerian simulations. Unfortunately, to the best of our knowledge, no direct comparison exists beyond the indirect estimate of Chaudhuri, Majumdar & Nath (2013) who find that the amount of energy per unit mass that the AGN has to provide is larger for SPH than for AMR. While the lack of entropy cores in purely adiabatic SPH simulations has been known for a long time (see Frenk et al. 1999; Power, Read & Hobbs 2014, for the original result and a recent explanation of its origin), more modern formulations of SPH are able to alleviate this discrepancy with Eulerian methods (see e.g. Rasia et al. 2014; Sembolini et al. 2016, for recent method comparisons). At the same time, the effect of feedback models in controlled experiments is not well documented. We hope that a comparison of the numerical discrepancies, as suggested by our results, will be undertaken by the community in the near future. The results we have obtained in this paper might however suggest that in SPH simulations, purely thermal feedback could affect much larger scales than in AMR simulations. It is at least plausible that in Lagrangian schemes a purely thermal feedback always ultimately leads to either an associated kinetic component or an outward convection of heated particles. Such a resolution-scale convection would of course be strongly suppressed in Eulerian methods.
## 7.4 A note on metal enrichment

We find (see Paper 2, Martizzi et al. 2016, for more details) that both our ICM as well as the galaxy population have a metallicity that is slightly lower than what is observed in the gas (e.g. De Grandi et al. 2004; Matsushita 2011; Werner et al. 2013) and stellar metallicity (e.g. Gallazzi et al. 2005). A metallicity distribution more consistent with observations has been reproduced in various SPH simulations (e.g. Sijacki et al. 2007; Planelles et al. 2014; Rasia et al. 2015). Also, in the ILLUSTRIS simulation (not SPH) realistic metallicities were produced (Vogelsberger et al. 2014), albeit at the price of baryon fractions that are completely inconsistent with observational constraints (Haider et al. 2016, their fig. 1). On the other hand, our results are consistent with the low metallicities found by other AMR simulations, e.g. Dubois et al. (2014). This aspect is another possible systematic discrepancy between the methods and should be investigated in more detail. We note that standard SPH does not include any mixing of metals (or any tracers, although explicit metal diffusion can be added; see e.g. Shen, Wadsley & Stinson 2010), while Eulerian schemes are forced to diffuse metals at the resolution scale. Planelles et al. (2014) included a smoothed metallicity estimate when calculating cooling times and state that their unsmoothed metallicity field is very noisy. We investigate the aspect of metal enrichment in our simulations further in Paper 2, Martizzi et al. (2016).

# 8 SUMMARY AND CONCLUSIONS

We present simulations of nine clusters of mass $M_{vir} \sim 10^{15} M_{\odot}$ and one of twice that mass, including cooling, star formation, as well as SN and AGN feedback with a physical resolution of $3.8 h^{-1}$ kpc. The simulations include the environment of the clusters inside of spheres of $8 h^{-1}$ Mpc around each cluster at $z=0$. In this paper, we compare in detail the ICM density, temperature, entropy and gas depletion profiles with X-ray data by performing a time-ensemble analysis for each cluster over a narrow mass bin and comparing with observed clusters in the same mass range. Next, we investigate the evolution of our simulated clusters over cosmic time with a range of cosmological observables that serve as mass proxies. In this paper, we focus on the mass versus Compton-y, X-ray luminosity and X-ray temperature scaling relations that are of particular importance for cluster cosmology. We also establish the numerical convergence of our results with resolution and their robustness against changes in the AGN feedback parameters.

We summarize our findings as follows.

(i) We find a persistent CC/NCC dichotomy in our clusters. The CCs are insensitive to changes in the thermal AGN feedback model parameters.

(ii) We link the disruption of CCs to low angular momentum major mergers. Major mergers with enough angular momentum leave the CCs intact. Core disruption occurs on time-scales of at most a few 100 Myr, with much increased core-cooling times after the disruption, leading to a quick transition and thus a stable bimodality.

(iii) Our simulations agree with the Planck $M_{500}-Y_{500}$ scaling relation with very little scatter. We do not identify a strong dependence of the scatter on the accretion history or the AGN model parameters.

(iv) The RG clusters are more X-ray luminous than a comparison sample from the MCXC catalogue. The clusters evolve along scaling relations in the $M-L_X$ plane that are consistent with self-similar scaling. There is a slight indication that the scatter in this relation correlates with details of the assembly history.

(v) The NCC clusters reproduce density, entropy and mass profiles of an ACCEPT comparison sample well, and are roughly consistent with the observed gas depletion profiles. The CC systems have excess central gas and a too low central entropy compared to the ACCEPT clusters. In addition, there is a general indication that at large radii, the simulated clusters have a slightly too low entropy and temperature, and a slightly too high density compared to observations.

(vi) The galaxies forming in our simulations have realistic masses and are consistent with abundance-matching results across three decades in halo mass. At higher masses, the simulated galaxies are slightly more massive than observed at a given halo mass, although we caution that observational issues complicate a detailed comparison. The SFRs at the high-mass end are consistent with recent observational constraints for BCGs.

The discrepancies we observed and listed above are plausibly related to shortcomings of the simplistic central thermal blast-wave model:

(i) In our AMR simulations, we find that thermal AGN feedback does not affect the ICM at significantly large radii. We see no effect on gas at scales of ~$R_{vir}/2$, nor is the AGN able to mildly stabilize the CC systems. Once a CC forms, it cools below observed core entropies and leads to a core baryon fraction inconsistent with observational constraints, which then explains the discrepancies we see in comparison to observations.

(ii) The above finding is discrepant with SPH results in the literature. In our simulations, details of the AGN energy injection are irrelevant for global ICM properties. In particular, the X-ray and SZ scaling relations are unaffected by details of the thermal AGN model. This is quite in contrast to the findings of, e.g. Le Brun et al. (2014) and points to a possible discrepancy between SPH/Lagrangian and Eulerian methods and how feedback couples to gas in such simulations.

The inability of the thermal AGN model to shape larger scales in our simulations plausibly points to other forms of energy injection (e.g. through kinetic feedback; see Choi et al. 2015; Li et al. 2015) or additional processes in shaping the ICM. Several published results suggest that thermal conduction might play a central role both in stabilizing CCs and at larger scales (e.g. Guo, Oh & Ruszkowski 2008; Parrish et al. 2009; Ruszkowski & Oh 2010; Arth et al. 2014). We will investigate these aspects in future research.

# ACKNOWLEDGEMENTS

We thank Yohan Dubois, James Bartlett, Pawel Biernacki, Elena Rasia and Amandine Le Brun for discussions and comments. We are grateful to Peter Behroozi for kindly making his ROCKSTAR-GALAXIES code available to us for further modifications, and for providing data shown in Figs 9 and 10. Finally, we thank the anonymous referee for a careful reading of our draft as well as insightful and valuable suggestions.

OH acknowledges support from the Swiss National Science Foundation (SNSF) through the Ambizione fellowship. H-YW acknowledges the support by the US Department of Energy under contract number DE-FG02-95ER40899. DM acknowledges support from the Swiss National Science Foundation. RHW received support from the US Department of Energy contract to SLAC no. DE-AC02-76SF0051. This work was supported by a grant from the Swiss National Supercomputing Centre (CSCS) under project ID s416.
REFERENCES

Allen S. W., Rapetti D. A., Schmidt R. W., Ebeling H., Morris R. G., Fabian A. C., 2008, MNRAS, 383, 879

Allen S. W., Evrard A. E., Mantz A. B., 2011, ARA&A, 49, 409

Angulo R. E., Hahn O., Abel T., 2013, MNRAS, 434, 1756

Arnaud M., Pratt G. W., Piffaretti R., Böhringer H., Croston J. H., Pointecouteau E., 2010, A&A, 517, A92

Arth A., Dolag K., Beck A. M., Petkova M., Lesch H., 2014, MNRAS, preprint (arXiv:1412.6533)

Battaglia N., Bond J. R., Pfrommer C., Sievers J. L., Sijacki D., 2010, ApJ, 725, 91

Battaglia N., Bond J. R., Pfrommer C., Sievers J. L., 2013, ApJ, 777, 123

Behroozi P. S., Wechsler R. H., Wu H.-Y., 2013a, ApJ, 762, 109

Behroozi P. S., Wechsler R. H., Wu H.-Y., Busha M. T., Klypin A. A., Primack J. R., 2013b, ApJ, 763, 18

Behroozi P. S., Wechsler R. H., Conroy C., 2013c, ApJ, 770, 57

Bernardi M., Meert A., Sheth R. K., Vikram V., Huertas-Company M., Mei S., Shankar F., 2013, MNRAS, 436, 697

Biffi V., Dolag K., Böhringer H., Lemson G., 2012, MNRAS, 420, 3545

Biffi V., Sembolini F., De Petris M., Valdarnini R., Yepes G., Gottlöber S., 2014, MNRAS, 439, 588

Bleuler A., Teyssier R., Carassou S., Martizzi D., 2015, Comput. Astrophys. Cosmol., 2, 5

Booth C. M., Schaye J., 2009, MNRAS, 398, 53

Bryan G. L., Norman M. L., 1998, ApJ, 495, 80

Burns J. O., Hallman E. J., Gantner B., Motl P. M., Norman M. L., 2008, ApJ, 675, 1125

Cavagnolo K. W., Donahue M., Voit G. M., Sun M., 2009, ApJS, 182, 12

Chaudhuri A., Majumdar S., Nath B. B., 2013, ApJ, 776, 84

Chen Y., Reiprich T. H., Böhringer H., Ikebe Y., Zhang Y.-Y., 2007, A&A, 466, 805

Choi E., Ostriker J. P., Naab T., Oser L., Moster B. P., 2015, MNRAS, 449, 4105

Christensen C. R., Brooks A. M., Fisher D. B., Governato F., McCleary J., Quinn T. R., Shen S., Wadsley J., 2014, MNRAS, 440, L51

Conroy C., Wechsler R. H., Kravtsov A. V., 2006, ApJ, 647, 201

Croton D. J. et al., 2006, MNRAS, 365, 11

De Grandi S., Ettori S., Longhetti M., Molendi S., 2004, A&A, 419, 7

De Lucia G., Blaizot J., 2007, MNRAS, 375, 2

Donahue M. et al., 2015, ApJ, 805, 177

Dubois Y., Teyssier R., 2008, A&A, 477, 79

Dubois Y., Devriendt J., Teyssier R., Slyz A., 2011, MNRAS, 417, 1853

Dubois Y., Devriendt J., Slyz A., Teyssier R., 2012, MNRAS, 420, 2662

Dubois Y. et al., 2014, MNRAS, 444, 1453

Eckert D., Molendi S., Paltani S., 2011, A&A, 526, A79

Evrard A. E., Summers F. J., Davis M., 1994, ApJ, 422, 11

Fabian A. C., 2012, ARA&A, 50, 455

Frenk C. S. et al., 1999, ApJ, 525, 554

Gallazzi A., Charlot S., Brinchmann J., White S. D. M., Tremonti C. A., 2005, MNRAS, 362, 41

Gaspari M., Brighenti F., Ruszkowski M., 2013a, Astron. Nachr., 334, 394

Gaspari M., Ruszkowski M., Oh S. P., 2013b, MNRAS, 432, 3401

Guedes J., Callegari S., Madau P., Mayer L., 2011, ApJ, 742, 76

Guo F., Mathews W. G., 2010, ApJ, 717, 937

Guo F., Oh S. P., 2009, MNRAS, 400, 1992

Guo F., Oh S. P., Ruszkowski M., 2008, ApJ, 688, 859

Haardt F., Madau P., 1996, ApJ, 461, 20

Hahn O., Abel T., 2011, MNRAS, 415, 2101

Haider M., Steinhauser D., Vogelsberger M., Genel S., Springel V., Torrey P., Hernquist L., 2016, MNRAS, 457, 3024

Jeltema T. E., Hallman E. J., Burns J. O., Motl P. M., 2008, ApJ, 681, 167

Katz N., White S. D. M., 1993, ApJ, 412, 455

Knebe A. et al., 2013, MNRAS, 428, 2039

Kravtsov A. V., Borgani S., 2012, ARA&A, 50, 353

Kravtsov A., Vikhlinin A., Meshscheryakov A., 2014, Astrophys. J., preprint (arXiv:1401.7329)

Le Brun A. M. C., McCarthy I. G., Schaye J., Ponman T. J., 2014, MNRAS, 441, 1270

Leitherer C. et al., 1999, ApJS, 123, 3

Li Y., Bryan G. L., 2014, ApJ, 789, 54

Li Y., Bryan G. L., Ruszkowski M., Voit G. M., O'Shea B. W., Donahue M., 2015, ApJ, 811, 73

Liu F. S., Mao S., Meng X. M., 2012, MNRAS, 423, 422

McDonald M. et al., 2013, ApJ, 774, 23

McDonald M. et al., 2015, ApJ, 811, 111

McNamara B. R., Nulsen P. E. J., 2007, ARA&A, 45, 117

McNamara B. R., Nulsen P. E. J., 2012, New J. Phys., 14, 055023

Mantz A., Allen S. W., Ebeling H., Rapetti D., Drlica-Wagner A., 2010, MNRAS, 406, 1773

Mantz A. B., Allen S. W., Morris R. G., Rapetti D. A., Applegate D. E., Kelly P. L., von der Linden A., Schmidt R. W., 2014, MNRAS, 440, 2077

Mantz A. B., Allen S. W., Morris R. G., Schmidt R. W., 2016, MNRAS, 456, 4020

Martizzi D., Teyssier R., Moore B., 2012a, MNRAS, 420, 2859

Martizzi D., Teyssier R., Moore B., Wentz T., 2012b, MNRAS, 422, 3081

Martizzi D., Teyssier R., Moore B., 2013, MNRAS, 432, 1947

Martizzi D., Hahn O., Wu H.-Y., Evrard A. E., Teyssier R., Wechsler R. H., 2016, MNRAS, 459, 4408

Matsushita K., 2011, A&A, 527, A134

Mazzotta P., Rasia E., Moscardini L., Tormen G., 2004, MNRAS, 354, 10

Meece G. R., O'Shea B. W., Voit G. M., 2015, ApJ, 808, 43

Melin J.-B., Bartlett J. G., Delabrouille J., Arnaud M., Piffaretti R., Pratt G. W., 2011, A&A, 525, A139

Miniati F., Beresnyak A., 2015, Nature, 523, 59

Moster B. P., Somerville R. S., Maulbetsch C., van den Bosch F. C., Macciò A. V., Naab T., Oser L., 2010, ApJ, 710, 903

Navarro J. F., Frenk C. S., White S. D. M., 1997, ApJ, 490, 493

Nelson K., Lau E. T., Nagai D., Rudd D. H., Yu L., 2014, ApJ, 782, 107

Parrish I. J., Quataert E., Sharma P., 2009, ApJ, 703, 96

Parrish I. J., Quataert E., Sharma P., 2010, ApJ, 712, L194

Piffaretti R., Arnaud M., Pratt G. W., Pointecouteau E., Melin J.-B., 2011, A&A, 534, A109

Pike S. R., Kay S. T., Newton R. D. A., Thomas P. A., Jenkins A., 2014, MNRAS, 445, 1774

Planck Collaboration XVI, 2014, A&A, 571, A16

Planck Collaboration XX, 2014, A&A, 571, A20

Planck Collaboration XXVII, 2016, A&A, 594, A27

Planck Collaboration XXXII, 2015, A&A, 581, A14

Panelles S., Quilis V., 2009, MNRAS, 399, 410

Panelles S., Borgani S., Dolag K., Ettori S., Fabjan D., Murante G., Tornatore L., 2013, MNRAS, 431, 1487

Panelles S., Borgani S., Fabjan D., Killedar M., Murante G., Granato G. L., Ragone-Figueroa C., Dolag K., 2014, MNRAS, 438, 195

Poole G. B., Babul A., McCarthy I. G., Sanderson A. J. R., Fardal M. A., 2008, MNRAS, 391, 1163

Power C., Read J. I., Hobbs A., 2014, MNRAS, 440, 3243

Pratt G. W., Croston J. H., Arnaud M., Böhringer H., 2009, A&A, 498, 361

Pratt G. W. et al., 2010, A&A, 511, A85

Rasia E., Mazzotta P., Evrard A., Markevitch M., Dolag K., Meneghetti M., 2011, ApJ, 729, 45

Rasia E. et al., 2014, ApJ, 791, 96

Rasia E. et al., 2015, ApJ, 813, L17

Ruszkowski M., Oh S. P., 2010, ApJ, 713, 1332

Schaye J. et al., 2015, MNRAS, 446, 521

Sembolini F., Yepes G., De Petris M., Gottlöber S., Lamagna L., Comis B., 2013, MNRAS, 429, 323

Sembolini F. et al., 2016, MNRAS, 457, 4063

Semler D. R. et al., 2012, ApJ, 761, 183

Shen S., Wadsley J., Stinson G., 2010, MNRAS, 407, 1581

Sijacki D., Springel V., Di Matteo T., Hernquist L., 2007, MNRAS, 380, 877
Smith R. K., Brickhouse N. S., Liedahl D. A., Raymond J. C., 2001, ApJ, 556, L91

Springel V., 2005, MNRAS, 364, 1105

Springel V., Di Matteo T., Hernquist L., 2005, MNRAS, 361, 776

Sutherland R. S., Dopita M. A., 1993, ApJS, 88, 253

Teyssier R., 2002, A&A, 385, 337

Teyssier R., Moore B., Martizzi D., Dubois Y., Mayer L., 2011, MNRAS, 414, 195

Vikhlinin A. et al., 2009, ApJ, 692, 1033

Vogelsberger M. et al., 2014, Nature, 509, 177

Voit G. M., Donahue M., Bryan G. L., McDonald M., 2015, Nature, 519, 203

Volonteri M., Natarajan P., Gültekin K., 2011, ApJ, 737, 50

Werner N., Urban O., Simionescu A., Allen S. W., 2013, Nature, 502, 656

White M., Hernquist L., Springel V., 2002, ApJ, 579, 16

Woo J. et al., 2013, MNRAS, 428, 3306

Wu H.-Y., Hahn O., Wechsler R. H., Mao Y.-Y., Behroozi P. S., 2013a, ApJ, 763, 70

Wu H.-Y., Hahn O., Wechsler R. H., Behroozi P. S., Mao Y.-Y., 2013b, ApJ, 767, 23

Wu H.-Y., Evrard A. E., Hahn O., Martizzi D., Teyssier R., Wechsler R. H., 2015, MNRAS, 452, 1982

## APPENDIX A: COMPARISON WITH NON-RADIATIVE RUNS

In Fig. A1, we show the $z=0$ profiles of electron density, entropy, temperature and gas depletion for the CC cluster RG 545 (blue lines) and the NCC cluster RG 348 (red lines) when radiative cooling, star and black hole formation and all feedback processes are disabled (solid lines). For comparison, we also show the $z=0$ results using the fiducial full physics runs discussed in the main part of the paper as dashed lines. We note that even in the non-radiative runs, the CC cluster has a higher density and lower entropy inside ~200 kpc. However, obviously no CC can emerge in the non-radiative runs, but it is clear that the merger that is occurring for RG348 leaves also a distinct signature even when no cooling is present.

<figure><img src="image_20.png" /><figcaption>**Figure A1.** Comparison between the full physics runs ('CSF', dashed lines) and non-radiative runs (solid lines) for the CC halo RG545 and the NCC halo RG348. We show the same profiles as in Fig. 3 as well as the gas depletion profiles (as in Fig. 6), but at $z=0$ only and thus not averaged over different snapshots.</figcaption></figure>
APPENDIX B: ACCEPT CLUSTER
COMPARISON SUBSAMPLE

In Table B1, we list the key properties of the subset of the ACCEPT clusters that we use for comparison in the main text of this article.

<table><caption>**Table B1.** Catalogue of ACCEPT clusters that have been selected for comparison with our simulated clusters. The clusters occupy the mass range between 4 and $6 \times 10^{14} M_{\odot}$ in $M_{500}$ based on both X-ray mass estimates as given in the MCXC catalogue and SZ mass estimates as given in the Planck 2015 SZ union catalogue. In this table, we list additional properties of these clusters as compiled from the ACCEPT, MCXC and Planck 2015 union catalogues. Compton-y $Y_{5R500}$ have been converted from the Planck SZ catalogue quantities given in arcmin$^2$ to kpc$^2$ assuming $h = 0.7$, $\Omega_m = 0.3$ and $\Omega_\Lambda = 0.7$ in the angular diameter distance relation.</caption><thead><tr><th>Cluster name</th><th>z</th><th>M<sup>MCXC</sup><sub>500</sub><br>(×10<sup>14</sup> M<sub>☉</sub>)</th><th>R<sup>MCXC</sup><sub>500</sub><br>(Mpc)</th><th>M<sup>Planck</sup><sub>SZ</sub><br>(×10<sup>14</sup> M<sub>☉</sub>)</th><th>T<sub>cl</sub><br>(keV)</th><th>K<sub>0</sub><br>(keV cm<sup>2</sup>)</th><th>L<sup>MCXC</sup><sub>X,500</sub><br>(10<sup>44</sup> erg s<sup>-1</sup>)</th><th>E<sup>-2/3</sup>Y<sub>500</sub><br>(kpc<sup>2</sup>)</th><th>Cool core?</th></tr></thead><tbody><tr><td>ABELL 85</td><td>0.0558</td><td>5.32</td><td>1.21</td><td>4.92</td><td>6.9</td><td>12.5</td><td>5.10</td><td>1.70</td><td>+</td></tr><tr><td>ABELL 141</td><td>0.23</td><td>4.72</td><td>1.10</td><td>5.67</td><td>5.31</td><td>205.03</td><td>5.16</td><td>19.97</td><td>—</td></tr><tr><td>ABELL 267</td><td>0.23</td><td>4.93</td><td>1.11</td><td>5.04</td><td>6.79</td><td>168.56</td><td>5.53</td><td>13.99</td><td>—</td></tr><tr><td>ABELL 399</td><td>0.0716</td><td>4.25</td><td>1.12</td><td>5.24</td><td>5.8</td><td>153.2</td><td>3.59</td><td>4.08</td><td>—</td></tr><tr><td>ABELL 586</td><td>0.171</td><td>5.20</td><td>1.16</td><td>5.17</td><td>8.7</td><td>94.75</td><td>5.62</td><td>6.23</td><td>—</td></tr><tr><td>ABELL 611</td><td>0.288</td><td>4.60</td><td>1.06</td><td>5.50</td><td>6.69</td><td>124.93</td><td>5.33</td><td>64.94</td><td>—</td></tr><tr><td>ABELL 907</td><td>0.1527</td><td>5.03</td><td>1.14</td><td>5.41</td><td>5.04</td><td>23.38</td><td>5.30</td><td>5.85</td><td>—</td></tr><tr><td>ABELL 963</td><td>0.2056</td><td>4.73</td><td>1.11</td><td>5.83</td><td>6.6</td><td>55.77</td><td>5.03</td><td>13.02</td><td>—</td></tr><tr><td>ABELL 1413</td><td>0.1426</td><td>5.55</td><td>1.19</td><td>5.95</td><td>8.9</td><td>64.03</td><td>6.04</td><td>6.71</td><td>—</td></tr><tr><td>ABELL 1650</td><td>0.0843</td><td>4.12</td><td>1.10</td><td>4.45</td><td>5.89</td><td>37.96</td><td>3.47</td><td>1.75</td><td>—</td></tr><tr><td>ABELL 1651</td><td>0.084</td><td>4.39</td><td>1.13</td><td>5.07</td><td>7.0</td><td>89.46</td><td>3.85</td><td>2.64</td><td>—</td></tr><tr><td>ABELL 1664</td><td>0.1276</td><td>4.06</td><td>1.08</td><td>4.28</td><td>3.5</td><td>14.4</td><td>3.57</td><td>3.53</td><td>+</td></tr><tr><td>ABELL 1795</td><td>0.0625</td><td>5.53</td><td>1.22</td><td>4.47</td><td>7.8</td><td>18.99</td><td>5.48</td><td>0.94</td><td>+</td></tr><tr><td>ABELL 1995</td><td>0.3186</td><td>5.87</td><td>1.14</td><td>4.92</td><td>8.6</td><td>374.35</td><td>8.28</td><td>21.18</td><td>—</td></tr><tr><td>ABELL 2034</td><td>0.113</td><td>4.07</td><td>1.09</td><td>5.85</td><td>7.15</td><td>232.64</td><td>3.51</td><td>3.17</td><td>—</td></tr><tr><td>ABELL 2069<sup>a</sup></td><td>0.116</td><td>4.57</td><td>1.13</td><td>5.31</td><td>7.9</td><td>453.25</td><td>4.26</td><td>8.29</td><td>—</td></tr><tr><td>ABELL 2104</td><td>0.1554</td><td>4.42</td><td>1.10</td><td>5.74</td><td>9.31</td><td>160.61</td><td>4.23</td><td>15.39</td><td>—</td></tr><tr><td>ABELL 2111</td><td>0.23</td><td>4.66</td><td>1.09</td><td>5.73</td><td>8.02</td><td>107.36</td><td>5.05</td><td>16.17</td><td>—</td></tr><tr><td>ABELL 2244</td><td>0.0967</td><td>4.49</td><td>1.13</td><td>4.38</td><td>5.57</td><td>57.58</td><td>4.05</td><td>3.15</td><td>—</td></tr><tr><td>ABELL 2294</td><td>0.178</td><td>4.23</td><td>1.08</td><td>5.98</td><td>7.1</td><td>156.31</td><td>4.05</td><td>6.87</td><td>—</td></tr><tr><td>ABELL 2409</td><td>0.1479</td><td>4.63</td><td>1.12</td><td>5.06</td><td>5.5</td><td>73.81</td><td>4.53</td><td>6.75</td><td>—</td></tr><tr><td>ABELL 3364</td><td>0.1483</td><td>4.50</td><td>1.11</td><td>4.89</td><td>6.59</td><td>268.55</td><td>4.32</td><td>4.99</td><td>—</td></tr><tr><td>ABELL 3571</td><td>0.0391</td><td>4.51</td><td>1.15</td><td>4.63</td><td>7.6</td><td>79.31</td><td>3.82</td><td>0.63</td><td>—</td></tr><tr><td>ABELL 3827</td><td>0.0984</td><td>4.59</td><td>1.14</td><td>5.77</td><td>8.05</td><td>164.58</td><td>4.20</td><td>3.79</td><td>—</td></tr><tr><td>CL J1226.9+3332</td><td>0.89</td><td>4.39</td><td>0.83</td><td>5.70</td><td>10.4</td><td>166.03</td><td>11.25</td><td>181.25</td><td>—</td></tr><tr><td>MS 0735.6+7421</td><td>0.216</td><td>4.60</td><td>1.09</td><td>5.02</td><td>5.45</td><td>15.96</td><td>4.87</td><td>10.28</td><td>+</td></tr><tr><td>MS 0906.5+1110</td><td>0.163</td><td>4.74</td><td>1.12</td><td>5.42</td><td>8.1</td><td>104.23</td><td>4.86</td><td>7.21</td><td>—</td></tr><tr><td>RXCJ0331.1-2100</td><td>0.188</td><td>4.05</td><td>1.06</td><td>4.34</td><td>4.61</td><td>11.4</td><td>3.82</td><td>5.36</td><td>+</td></tr><tr><td>ZWCL 1358+6245</td><td>0.328</td><td>4.61</td><td>1.05</td><td>4.81</td><td>7.2</td><td>20.67</td><td>5.62</td><td>29.16</td><td>+</td></tr></tbody></table>

Note. $^a$Excluded.

This paper has been typeset from a TeX/LA TeX file prepared by the author.