© The Authors 2023

1. Introduction

The hot intracluster medium (ICM) is the largest baryonic component in groups and clusters of galaxies (Gonzalez et al. 2013). It is observable in the X-ray band thanks to its strong Bremsstrahlung continuum emission plus emission lines from highly ionized elements, and it shows temperatures from 1 keV (groups) to more than 10 keV (massive clusters). The density profile of the ICM is usually fitted with the β profile (Cavaliere & Fusco-Femiano 1978), consisting of a flat core and rapidly decreasing outskirts. However, in many clusters, the central electron density is not accurately described by a single β profile, but is observed to be sharply peaked, reaching values significantly larger than the typical n_e ∼ 10⁻³ cm⁻³. In these cases, a double β profile is needed to fit the X-ray surface brightness (as known from ROSAT observations Xue & Wu 2000), effectively defining a central core where the ICM cooling time, t_cool ∝ kT/n_e, is significantly lower than the typical age of the cluster. This implies that a large amount of gas should cool completely on a short timescale (typically less than 1 Gyr) due to radiative losses, leading to a massive cooling flow with a mass deposition rate of the order of several hundreds up to thousands M_⊙ yr⁻¹. Such large values are obtained directly from the brightness profile under the assumption of subsonic flow and constant pressure, as per the so-called “isobaric cooling flow” model (Fabian & Nulsen 1977; Fabian 1994).

The X-ray luminosity in cluster cores is dominated by the hottest ICM component (above a few keV), while the coldest one contributes only a few percent of the total, while being rich in emission lines. At variance with the “isobaric cooling flow” scenario, high-resolution spectroscopy of bright clusters with XMM–Newton have not shown any evidence for the multiphase, line-rich gas predicted by the isobaric cooling model. Instead, it has been observed that the majority of the ICM typically reaches a temperature plateau at about one-third of the virial value (Kaastra et al. 2001; Peterson et al. 2001; Tamura et al. 2001; Donahue & Voit 2004), while the cold gas below this floor is virtually absent. In this framework, the old paradigm of “cooling flow” has been abandoned in favour of the “cool core” scenario (Molendi & Pizzolato 2001). Despite some amount of cooling gas, possibly associated with star-forming episodes in the central galaxy, is still allowed by the observations, the upper limits to the spectroscopic mass-deposition rate values are at least one order of magnitude lower than the rates expected by isobaric cooling flow models (McNamara & O’Connell 1989; Makishima et al. 2001; Edge 2001; Edge & Frayer 2003; McNamara & Nulsen 2007). Recently, spatially resolved spectral analyses of cool cores have confirmed that the central, isobaric mass-deposition rates are significantly lower than the star-formation rates observed in the hosted brightest cluster galaxy (BCG; Molendi et al. 2016). On the other hand, the more convincing cooling flow candidates are limited to only one well-documented case (the Phoenix cluster, see McDonald et al. 2012, 2019; Tozzi et al. 2015; Pinto et al. 2018).

This observational evidence strongly supports the presence of some heating mechanism that prevents the bulk of the ICM from cooling below about one-third of the virial temperature in the cool core. This is now considered a standard condition of the ICM in cool-core clusters, which represent ∼70% of the low-redshift population of X-ray flux-limited sample, according to Hudson et al. (2010). On the other hand, some amount of cold and multiphase gas is now commonly observed in the submm and optical band in star forming regions of the BCG thanks to ALMA (e.g., McNamara et al. 2014; Russell et al. 2017; Temi et al. 2018; Tremblay et al. 2018; Rose et al. 2019; North et al. 2021) and MUSE (e.g., Olivares et al. 2019, 2022; Maccagni et al. 2021), respectively. The presence of multiphase gas suggests that short-lived cooling flows raining all the way down onto the central supermassive black hole (SMBH) in the BCG have time to replenish the cold gas reservoir before being quenched by the ensuing SMBH feedback process. In summary, we are well aware that cool cores are physical systems where the cooling process is counterbalanced by some global heating mechanism that strongly suppresses the mass deposition rate, while still allowing some amount of gas to leak out of the hot phase with a timescale regulated by a complex feeding (cooling) and feedback (heating) cycle (Gaspari & Sądowski, 2017). A comprehensive and systematic comparison of the spectroscopic mass deposition rate to the star formation rate in the BCG and the presence of molecular gas in the cluster core, significantly extending the small sample explored in Molendi et al. (2016), would provide very effective constraints on the baryonic cycle in clusters. In addition, the detection of very diffuse, low-temperature ICM in the center of non-cool-core clusters by the next generation of X-ray bolometers, may provide support for a scenario in which cool-core clusters rapidly switch into the non cool-core phase and vice versa (see Molendi et al. 2023).

Many heating mechanisms have been proposed in the past two decades, among them: thermal conduction (Zakamska & Narayan 2003), viscous dissipation of sound waves (Ruszkowski et al. 2004), supernova feedback (Domainko et al. 2004), turbulence combined with conduction (Dennis & Chandran 2005), cosmic ray–ICM interaction (Guo & Oh 2008; Yang et al. 2019), and feedback from jets and outflows from the central active galactic nucleus (AGN; e.g., McNamara & Nulsen 2007; Gaspari et al. 2012b; Barai et al. 2016; Wittor & Gaspari 2020; McKinley et al. 2022). In particular, the last process is considered the most likely contributor on the basis of the well-documented interactions between radio jets and the surrounding ICM. The large amount of mechanical energy associated with the cavities carved into the ICM by the radio jets may be eventually transformed into thermal energy of the ICM (Blanton et al. 2010, for a review) and may stimulate, at the same time, the cooling of some fraction of the gas (Gaspari et al. 2020, for a review). Therefore, no matter how many mechanisms are contributing, the central AGN is expected to play an important role in regulating cooling, ultimately inducing tight scaling relations between the SMBH and AGN as well as the hot halo properties (Gaspari et al. 2019; Pasini et al. 2021).

In this work we investigate key physical radii that delimit the spherical regions where different phases of the complex baryon cycle are actually taking place, such as the cool-core condensation radius (R_ccc) and the quenched cooling flow radius (R_qcf). In Sect. 2, we introduce and discuss our definition of R_ccc and R_qcf. In Sect. 3, we derive the typical timescales of relevant processes occurring in galaxy clusters. In Sect. 4, we describe the selection of the sample of galaxy clusters observed with Chandra and used in this work. In Sect. 5, we discuss the key parameter represented by the turbulent velocity dispersion of the warm and cold phase of the diffuse baryons. In Sect. 6, we describe data reduction and our analysis strategy. Our results are described in Sect. 7, where we show deprojected timescale profiles in each cluster as a function of radius and we measure the R_ccc and the R_qcf values, followed by an investigation of their distribution across the cluster sample. The physical implication of our findings are discussed in Sect. 8 and our conclusions are summarized in Sect. 9. Throughout the paper, the cosmological model of reference is a ΛCDM with parameters H₀ = 67.8 km s⁻¹ Mpc⁻¹, Ω_Λ = 0.692 and Ω_m = 0.308 (Planck Collaboration XIII 2016). Quoted errors and upper limits correspond to a 1-σ confidence level.

2. Physical definition of cool-core condensation and quenched cooling-flow radius

From the observational point of view, it has been well established that the presence of radio nuclear activity is closely associated with the presence of a cool core (Dunn & Fabian 2006; Sun 2009). The interactions between the relativistic electrons and the thermal electrons of the ICM have been thoroughly studied in spectacular images of few nearby clusters such as Perseus (Fabian et al. 2003a), Hydra A (McNamara et al. 2000), and few other clusters at intermediate redshift (Blanton et al. 2011; Ehlert et al. 2011); in addition, the presence of cavities in the ICM has been explored up to z ∼ 1.2 (Hlavacek-Larrondo et al. 2015). However, while the energy budget associated with cavities is sufficient to switch off the cooling in all the observed cases, the physical mechanism by which the energy of the jet is transferred isotropically to the ICM is still an issue of debate. Possible mechanisms include turbulence (e.g., Gaspari 2015) or weak shocks (e.g., Fabian et al. 2003b; Mathews et al. 2006) driven by the radio-mode activity of the central galaxy. The key issue here is to identify a process that regularly transforms an impulsive and directional energy input of the jet into a gentle heating, smoothly distributed in time and space, to finally shape the regular ICM thermodynamical properties observed in cool-core clusters.

All these details cannot be resolved in most cool-core clusters and, therefore, the actual physical processes can hardly be constrained from the macroscopic X-ray quantities such as luminosity and temperature. In recent years, several independent efforts have been devoted to achieve an efficient observational diagnostics to classify clusters according to the presence of a cool core and, at the same time, to understand the dominant physical processes. The many quantities used to define a cool core are all related to the thermodynamical properties of the ICM but they are associated with different physical processes: surface brightness excess or cuspiness (Santos et al. 2008), the temperature gradient (Sanderson et al. 2006; Burns et al. 2008), a steep iron abundance profile (De Grandi et al. 2004; Rasera et al. 2008), steep slopes of central entropy (Pratt et al. 2010), low central cooling time (O’Hara et al. 2006), classical mass deposition rate (Chen et al. 2007), and (clearly) the slope of the electron density profile (Hudson et al. 2010). All these diagnostics are obtained by combining the same three independent X-ray observables: surface brightness, temperature, and line emission.

In a comprehensive overview, Hudson et al. (2010) concluded that the central cooling time is one of the best diagnostics for identifying and characterizing cool cores. Furthermore, Gaspari et al. (2018) showed that the ratio of the cooling time over the turbulence eddy turnover timescale C ≡ t_cool/t_eddy ∼ 1 is a key diagnostic for the condensation extent of the multiphase rain occurring via chaotic cold accretion (hereafter, CCA). In other words, whenever the C-ratio approaches unity there is enough turbulence and quick cooling to directly drive non-linear instabilities. This is confirmed by studies obtained with high-resolution radio/optical telescopes (Gaspari et al. 2018; Olivares et al. 2019, 2022), which show that multiphase filamentary structures can be observed within the region enclosed by t_cool/t_eddy ≈ 1. In this region, the turbulent mixing rate is expected in part to balance the pure cooling flow, in part to drive direct nonlinear thermal instability that then ends up generating a rainfall. As shown by theoretical studies, the ICM can be seen as a hierarchical thermodynamic system that follows a chaotic, top-down multiphase condensation cascade (Gaspari et al. 2017; Voit et al. 2017). Therefore, we assume that an appropriate timescale for feedback is provided by the turbulence timescale, t_eddy. In this work, we leverage these findings to identify the approximately spherical region within which the nonlinear multiphase CCA rain and the triggered feedback response are very effective, via the following condition:

$$\begin{array}{c}\hfill C\equiv t_{\mathrm{cool}}(R_{\mathrm{ccc}})/t_{\mathrm{eddy}}(R_{\mathrm{ccc}})=1.\end{array}$$(1)

This relation effectively defines the cool-core condensation radius, R_ccc. Within R_ccc, we expect direct turbulence instability to be driving localized flickering precipitation and, therefore, the condensation of the light rain that may be responsible for accretion events onto the central SMBH in the BCG, and thus leading to star formation episodes in the BCG.

At the same time, we can identify a region where we assume thermally unstable cooling may ensue from linear perturbations by establishing a threshold in the ratio of the cooling, t_cool, and the free-fall time t_ff (e.g., Field 1965). Voit et al. (2015a) found that the minimum value of t_cool/t_ff fluctuate around values of 10 − 20, concluding that cold clouds start to precipitate out of hot-gas atmospheres when t_cool drops to ten times t_ff. Later, Hogan et al. (2017) showed that the minimum of the t_cool/t_ff ratio in a large sample of observed clusters with constrained nebular emission (tracing the condensed cool gas) is bound between 10 and 40, with few values below 10, which is also supported by hydrodynamical simulations with self-regulated AGN jet feedback (Gaspari et al. 2012b).

Therefore, we argue that an average ratio ≡t_cool/t_ff ∼ 25, despite a large scatter, is a reasonable proxy for tracing the initial growth of linear thermal instability (TI) in heated cooling flows, while a value of 10 traces the lower bound of such a criterion¹. Overall, we define a quenched cooling flow radius (and use it here) when the following condition is met:

$$\begin{array}{c}\hfill t_{\mathrm{cool}}(R_{\mathrm{qcf}})=25\times t_{\mathrm{ff}}(R_{\mathrm{qcf}}).\end{array}$$(2)

In this framework, the quantity R_qcf is expected to be an alternative definition to the “classical” cool-core radius defined on the basis of the cooling time. Indeed, the radius below is often used, whereby the cooling time is shorter than the reference value of 7.7 Gyr², R_classic ≡ R_7.7 Gyr; (Hudson et al. 2010).

Overall, we expect the two newly defined core radii to trace physical transitions from a macro-scale X-ray emitting ICM atmosphere delimited by a “classical” cool core radius to a region where a quenched cooling inflow could potentially develop (< R_qcf), and, eventually, to a region where precipitation and feedback are actively vigorous (< R_ccc). To explore the behaviour of these two spatial scales, we sought to measure R_qcf and R_ccc in an optimally selected sample of massive clusters observed with the Chandra satellite.

3. Timescales in the ICM of massive Galaxy Clusters

In this section, we define the timescale for physical processes relevant to the ICM, which is treated as an optically thin plasma in collisional ionization equilibrium, despite occasional out-of-equilibrium phases that may potentially be reached. However, treating the ICM in steady equilibrium is a fitting approximation and would not bias our results. As previously discussed, we will also assume spherical symmetry, as a requirement that will affect the sample selection in certain ways, as discussed in Sect. 4.

3.1. Cooling time

The ICM X-ray emission is composed of thermal bremsstrahlung (free-free emission) plus line emission from ions of heavy elements. The X-ray luminosity density (energy emitted per unit time at a unit volume) can therefore be written as $L_X\approx n_{\text{e}}^2\mathrm{\Lambda }(T,Z)$, where n_e is the electron density and Λ(T, Z) is the cooling function, which depends on the temperature of both radiative processes and is also related with the abundance, Z, of heavy elements in the ICM (see Fig. 3 of Peterson & Fabian 2006). For the hot ICM (kT > 2 keV), the bremsstrahlung emission dominates and the approximation Λ(T)∝T^1/2 is usually adopted. When temperatures are low, the number of ions (and therefore the number of possible transitions) strongly increases. As a consequence, the enhanced contribution from line emission significantly affects the cooling function. This regime is particularly relevant in cool cores, where low temperatures are always associated with high metallicity, often reaching supersolar values (De Grandi et al. 2004; Liu et al. 2020). Therefore, to describe the cooling efficiency of the ICM accurately at different radii, the full cooling function must be taken into account.

The main energy loss of the ICM is the thermal radiative emission, which is mostly observed in the classic 0.5–10 keV X-ray band for the temperature range we are considering here. Therefore, the cooling time is typically defined as the timescale after which the ICM entirely loses its internal energy via bremsstrahlung radiation, down to the point when the gas eventually recombines and starts loosing energy through other radiative processes. An effective way to estimate the ICM cooling time is obtained by dividing the gas internal energy by the luminosity density of the plasma. We can therefore express the internal energy as (3/2) nkT, obtaining the following for the cooling time:

$$\begin{array}{c}\hfill t_{\mathrm{cool}}\simeq \frac{3}{2}\frac{\mathit{nkT}}{n_{\mathrm{e}}{n}_{\mathrm{i}}\mathrm{\Lambda }(Z,T)},\end{array}$$(3)

where Λ(Z, T) is the cooling function for gas with a specific abundance, Z, and temperature, kT. In this work, we compute the cooling function interpolating the values reported in Sutherland & Dopita (1993). We note that here we adopt a definition of cooling time based on the internal energy rather than the enthalpy (5/2)nkT (Peterson & Fabian 2006). The factor of 5/2 is assumed to account for the inclusion of the extra work-term arising from perfect spherically symmetric isobaric compression. However, the 5/2 value should be considered as an upper limit, since (under realistic conditions) there is no perfect isobaric compression and any contribution from turbulence or AGN heating brings it near the pure 3/2 factor, as shown in Fig. 5 of Gaspari (2015).

For temperatures kT  ≳  2 keV (and metallicity 0.3 Z_⊙), only Bremsstrahlung emission is relevant for such clusters, thus t_cool is well approximated by the expression (see Cavagnolo et al. 2009):

$$\begin{array}{c}\hfill t_{\mathrm{cool}}\simeq {10}^8\mathrm{yr}{\left(\frac{K}{10{\mathrm{keV}\mathrm{cm}}^2}\right)}^{3/2}{\left(\frac{\mathit{kT}}{5\mathrm{keV}}\right)}^{-1},\end{array}$$(4)

where $K\equiv kT/n_{\mathrm{e}}^{2/3}$ is the astrophysical (electron) entropy (Ponman et al. 1999; Tozzi & Norman 2001). We note that using a higher solar metallicity reduces the normalization in Eq. (4) by only 30%. Moreover, hydrodynamical simulations predict that this cooling time ought to have a radial dependence approximated by a power law with a slope of ∼1.3 (Ettori & Brighenti 2008).

3.2. Free-fall time

The free-fall time is defined as the timescale of an object falling towards the center of the cluster without any pressure support or any deceleration due to viscosity (Binney & Tremaine 1987). The free-fall time can be written as:

$$\begin{array}{c}\hfill t_{\mathrm{ff}}=\sqrt{\frac{2r}{g(r)}}=\sqrt{\frac{2r{\rho }_g}{\mathrm{d}P/\mathrm{d}r}}=\sqrt{\frac{2r^3}{\mathit{GM}}},\end{array}$$(5)

where M is the total mass within a spherical radius, r. The free-fall time, therefore, depends on the total mass and not on the properties of the ICM, such as density. The total mass profile M(< r) can, in turn, be computed directly from the ICM applying the condition of hydrostatic equilibrium, which requires the knowledge of the density and temperature profiles of the ICM:

$$\begin{array}{c}\hfill M(<r)=-\frac{r^2}{G{\rho }_g}\frac{\mathrm{d}P}{\mathrm{d}r}=-r\frac{kT(r)}{\mu m_{p}G}(\frac{\mathrm{d}log(n_{\mathrm{e}})}{\mathrm{d}log(r)}+\frac{\mathrm{d}log(kT)}{\mathrm{d}log(r)}),\end{array}$$(6)

where kT(r) is the deprojected temperature profile of the ICM, n_e(r) the deprojected electron density profiles, μ is the mean molecular weight of ICM (which is usually assumed to be μ = 0.6). and m_p is the proton mass.

3.3. Turbulence timescale

A number of processes, such as cluster mergers, galaxy motions, and AGN feedback, are capable of producing turbulence in the ICM. The mechanical energy associated with these processes is very high and X-ray observations clearly show that the ICM is significantly affected by it (Iapichino et al. 2010; Gaspari et al. 2014; Liu et al. 2015; Lau et al. 2017; Wittor & Gaspari 2020). Combining X-ray and radio observations, it has been shown that radio jets injects bubbles of relativistic electrons that create cavities within the ICM on scales of 10–100 kpc (McNamara et al. 2005; Diehl et al. 2008; Blanton et al. 2011; Hlavacek-Larrondo et al. 2015; Yang et al. 2019). The expectation is that a relevant fraction of the mechanical energy associated with the buoyant rise of inflated bubbles (or with bulk motions of infalling halos) will transfer into the ICM and eventually produce turbulence. Measurements of turbulence in cluster cores have been reported, in particular, via the power spectra of density fluctuations (Schuecker et al. 2004; Sanders et al. 2010; Gaspari & Churazov 2013; Zhuravleva et al. 2014; Hofmann et al. 2016; Simionescu et al. 2019). The only measurement of turbulence in the ICM at high spectral resolution has been provided by the Hitomi mission (Hitomi & Aharonian 2016), which showed a few 100 km s⁻¹ in turbulent velocities in Perseus cluster. The ultimate observational evidence of the amplitude and distribution of turbulence in the ICM will be provided in the next future by X-ray high-resolution spectra, obtained with the X-ray spectrometers Resolve onboard XRISM³ and X-IFU onboard Athena⁴.

Following Gaspari et al. (2018), we can estimate the turbulence characteristic timescale via the eddy turnover/mixing⁵ time, such as

$$\begin{array}{c}\hfill t_{\mathrm{eddy}}=2\pi \frac{r^{2/3}{L}^{1/3}}{{\sigma }_{v,L}},\end{array}$$(7)

where L is the energy-injection scale and σ_v, L is the typical turbulent velocity measured at the injection scale. The injection scale, L, is related to the AGN feedback influence region and can be approximated via the typical observed size covered by the pair of cavities inflated in the ICM (alternatively, by the extent of the Hα nebula). A phenomenological scaling for L is obtained from the large sample by Shin et al. (2016), with L ∼ 10 (kT/1 keV)² kpc (see also Gaspari et al. 2019). Given the strong dependence on the average temperature, it can reach values as high as ∼200 kpc, as in the case of MS 0735.6+7421 (Vantyghem et al. 2014).

The other relevant parameter, namely, the turbulent velocity dispersion, σ_v, L, is preferentially taken from the literature (Gaspari et al. 2018; Olivares et al. 2019) whenever a direct measurement is available. If no measurements are available, we consider a range of values chosen (as detailed in Sect. 5). Another assumption here is that we assume that the three-dimensional velocity dispersion σ_v, L is obtained multiplying by $\sqrt{3}$ the line-of-sight velocity dispersion of the warm and cold gas measured in (Gaspari et al. 2018; Olivares et al. 2019), which implies isotropic turbulence. This may not be true when violent gas sloshing is present. In fact, during the initial phase, sloshing-induced turbulence may have larger velocity dispersion in the plane of sloshing motion. On the one hand, estimating and removing the effects of gas sloshing is beyond the capability of the current analysis. On the other hand, we argue that its impact is limited. This assumption is based on the fact that sloshing develops on much larger time scales (see ZuHone et al. 2013; ZuHone & Roediger 2016), while the AGN feedback is self-regulated on duty cycles corresponding to ∼1 − 10 Myr (see Gaspari et al. 2012a). Therefore, we expect that the main contribution to the turbulence level is provided by the AGN feedback, while sloshing may boost the turbulence on time scales ranging 0.1 – a few Gyr (for minor and major mergers, respectively; see discussion in Lau et al. 2017). Within this framework, we adopt the assumption of isotropic turbulence bearing in mind that t_eddy can be slightly affected by anisotropic bulk motions of the ICM.

4. Cluster sample selection

The term “cool core” typically refers to the central region of a cluster, approximately in hydrostatic equilibrium, which shows a density profile significantly peaked toward the center, the temperature profile increasing with the radius on a scale of ∼ 50 kpc, an iron abundance profile decreasing with radius on the same scale, and, therefore, a short cooling time. In the absence of a recent major merger, a cool core invariably hosts a BCG at its center.

To achieve a robust measurement of radial timescale profiles down to a few kpc, we decided to perform our study with Chandra data, thanks to its angular resolution unparalleled among current X-ray facilities, which allows us to resolve a scale below 10 kpc virtually at redshifts up to z ∼ 1. In this way, the effective resolution of our analysis will be determined solely by the surface brightness of the ICM emission and the exposure depth of the data. To optimally select a sample of clusters for our study, we considered the Chandra archive to collect all the targets where we can measure the timescale profiles well inside the potential cool core region with an accuracy sufficient to identify the transition radii, R_qcf and R_ccc. As a rule of thumb, we know that the typical size of cool cores is around 40 kpc (Santos et al. 2008), therefore, as the first criterion, we require enough counts within this radius to extract at least two rings with a minimum of 3000 net counts each in the 0.5–7 keV band. This will allow us to measure the temperature with an error lower than ∼20% in the spatial bins within 40 kpc.

On the other hand, we also need to measure temperature and density out to 300–400 kpc, in order to track how timescale profiles behave before approaching the cooling region. Due to the limited Chandra field of view (16′ × 16′ for ACIS-I, 8′ × 8′ for ACIS-S), the requirement to reach ∼ 400 kpc in the FOV translates in a conservative lower limit z > 0.03 and z > 0.1 for ACIS-I and ACIS-S, respectively. We also require z < 0.3 to have at least a width of 3 arcsec for the minimum ring width of ∼15 kpc.

We start from a total sample of 1144 galaxy clusters or groups, consisting of all the publicly available archival data under this category. Our requirement on the redshift range brings this number down to 456 clusters. After applying our criterion on the net number of counts and removing those clusters that strongly depart from a relaxed, a spherical morphology after a visual inspection, we obtained a total of 37 clusters at 0.03 < z < 0.3 satisfying our criteria. With the knowledge that a simple visual inspection does not guarantee a relaxed dynamics, with this step, we excluded all the clusters with disturbed morphologies that would make a spherical deprojection unreliable. We also note that the requirement on the net counts within 40 kpc strongly favors cool-core clusters, since to reach the same angular resolution in our spectral analysis, non-cool-core clusters need to have a significantly deeper exposure to match the same number of photons in the inner 40 kpc as compared to cool-core clusters. Therefore, non cool-core clusters are expected to be underrepresented in our final sample. However, since we are aiming at investigating the correlation between core properties and other observables (and not the distribution across the cluster population), the predominance of cool-core clusters does not strongly affect our conclusions. The list of clusters, with redshift, position, ObsID, and effective total exposure time after data reduction, is shown in Table 1. We note that not all the available exposures for each target are used. In particular, we removed all the off-centered pointings (aiming at the cluster outskirts) and we discarded some pointings with very low exposures (∼5 ks) or with an observation date very far from the bulk of the other observations (which would imply a significantly different effective area). Usually, we keep observations in both ACIS detectors, but in some cases, we discarded the observations in one of the detectors when the corresponding exposure is a minority of the total exposure time; this was done to achieve a simpler spectral analysis at the cost of a negligible loss in signal.

5. Turbulent velocity dispersion estimate

We collected all the measurements available (in the literature) of the turbulent velocity at a particular injection (5–10 kpc) of the warm or cold phase for our selected clusters, converting the line-of-sight measurement into a three-dimensional (3D) value, assuming isotropic turbulence, as discussed in Sect. 3.3. Here, we assume that the turbulent velocity dispersion of the ICM is very close to the value measured for the cold and warm medium. This correlation is shown to have a slope and a normalization close to unity, with a scatter of ∼30% (see Fig. 1 in Gaspari et al. 2018). We retrieved average ensemble measurements of σ_v for 14 clusters from Gaspari et al. (2018) that overlap with our sample. From Olivares et al. (2019), we collected two other new measurements. We note that Olivares et al. (2019) provided other measurements in clusters overlapping with Gaspari et al. (2018); however, they are tied to pencil-beam values, while ensemble values, such as those listed in Table 1 of Gaspari et al. (2018), should be used here, since they better trace the dynamical and thermodynamical properties of the turbulent medium, which we indicate as “macro weather”.

Therefore, 40% of our sample has values of 3D σ_v and associated uncertainties (directly derived from the observed line-of-sight σ_v) that can be directly inserted in Eq. (7) for the corresponding cluster. In Fig. 1 (left panel), we show the 16 clusters in the σ_v − kT plane. We note a large intrinsic scatter with no clear trend with the temperature. Such a large scatter likely dilutes any weak dependence on the (core-excised) temperature, given the current low statistics. The measured values are scattered in the range 200 < σ_v < 400 km s⁻¹, with the average value and rms given by:

$$\begin{array}{c}\hfill ⟨{\sigma }_{\mathrm{v}}⟩=306\pm 78{\mathrm{km}\mathrm{s}}^{-1}.\end{array}$$(8)

Fig. 1. Top panel: Turbulent velocity values measured in the warm or cold phase by Gaspari et al. (2018, blue points) and Olivares et al. (2019, orange points). The reference dotted line σ ∝ T1/3 corresponds to Eq. (9), with the cyan-shaded area showing the 1σ uncertainty. The solid lined corresponds to the constant average value given by Eq. (8), with the orange-shaded area showing its 1σ uncertainty. Bottom panel: the histogram distribution of the clusters with σv measurement (blue) and of the complementary subsample (green).

The constant value ⟨σ_v⟩ and its uncertainty are shown as an orange shaded area in Fig. 1 (left panel). The question is then whether the turbulent velocity in the clusters without direct measurements can be assumed to be within this range. As a first step, we compared the temperature distribution of the 18 clusters with measured σ_v to that of the 25 clusters without σ_v. The two distributions are shown in Fig. 1 (right panel), where it is possible to see that the clusters without σ_v measurement are hotter (hence, more massive) than those with σ_v. Therefore, assuming a constant values for σ_v, despite the large uncertainty, may not fully track the dependence on the temperature-mass scale.

The dependence of σ_v, L on the mass scale has been investigated only in a few works. Gaspari et al. (2018) find ${\sigma }_{\mathrm{v}}\propto M^{0.2}\propto T_{\mathrm{x}}^{0.3}$ in the larger sample of 72 groups and clusters, by leveraging ensemble-beam optical spectra, which are shown to linearly trace hot-gas turbulent velocities within a 0.13 dex scatter. By using a completely different method and band based on X-ray brightness and density fluctuations in 33 groups and clusters, Hofmann et al. (2016) retrieved a slightly negative slope in the hot-gas turbulent Mach numbers (defined as σ_v/c_s, with the sound speed being $c_{\mathrm{s}}\propto T_{\mathrm{x}}^{1/2}$), translating again to ${\sigma }_{\mathrm{v}}\propto T_{\mathrm{x}}^{0.3}$. Overall, if we fit the normalization of a relation with a slope of 0.3 to the available data, we find the following relation

$$\begin{array}{c}\hfill {\sigma }_{\mathrm{v}}=(315\pm 86){\left(\frac{\mathit{kT}}{4.85\mathrm{keV}}\right)}^{0.3}{\mathrm{km}\mathrm{s}}^{-1},\end{array}$$(9)

shown as the light-blue shaded area in Fig. 1. Given the current data, we cannot estimate the 3D σ_v in massive clusters with better approximation than that discussed here. To compute t_eddy(r) in the clusters without σ_v measurement, we preferred to use the σ_sv − kT relation given by Eq. (9); from a theoretical and simulation perspective, we do expect self-regulated AGN feedback (hence, turbulence) scaling at some level with halo mass from poor to rich clusters. Nevertheless, we also discuss the results obtained with both Eqs. (8) and (9) a posteriori to keep the systematics under control, finding that differences due to the different scaling relations are negligible when compared with the statistical uncertainties.

6. Data reduction and analysis

The data reduction was performed using the CIAO software (version 4.12) with CALDB 4.9. We appllied a charge transfer inefficiency correction, time-dependent gain adjustment, grade correction, and pixel randomization. First of all, each observation was reprocessed using the chandra_repro function. The script reads data from the standard data distribution and creates a new bad pixel file, a new level = 2 event file, and a new level = 2 type PHA file for each selected region, with the appropriate response files. We removed high-energy background flares from the event files with the deflare command. In addition, the background is reduced with the VFAINT cleaning whenever possible. The final exposure times are typically lower than the nominal exposure time only by a few percent. The level 2 files obtained in this way are reprojected to match the coordinates of the observation with the longest exposure for each clusters. The merged files are used only for imaging analysis, while spectra are extracted from each single Obsid.

The annular regions that were used to extract the spectra are centered on the emission peak. First, we identified the peak of the surface brightness profile on the total band image smoothed on a 3 arcsec scale and then the center of circle whose photometry maximize the S/N. Usually the two positions differ by few arcsec, given the regular shape of the selected clusters. If the difference is less than 2 arcsec, we used the position of the peak to fit the surface brightness; otherwise, we considered the difference of the two centroids to be significant and then we adopted the maximization of the S/N as a more robust estimate of the cluster center.

We adaptively chose the width of each annulus to guarantee at least 3000 net counts in the 0.5–7 keV band. Thanks to the exquisite angular resolution of Chandra, we do not need to correct any effects caused by the point spread function (PSF) when analyzing the spectra, since the PSF size is much smaller than the bin width. Before creating the spectra, we manually remove all the unresolved sources (mostly foreground and background AGN) visible in the soft, hard and total-band images. For each region, we produce response matrix (RMF) and ancillary response (ARF) files from each Obsid using CIAO. In this way, we keep track of all the differences in the ACIS effective area among different Obsid and, of course, among the different type of ACIS detectors.

Since each galaxy cluster has a different size, we need to define a normalized radius to express our results and compare different clusters. In this work, the radius r₅₀₀, defined as the radius enclosing an average total mass density 500 times larger than the critical density at the cluster redshift, is estimated for each cluster with the relation provided by Vikhlinin et al. (2006):

$$\begin{array}{c}\hfill r_{500}=\frac{0.792}{hE(z)}(\frac{T_{\mathrm{x}}}{5\mathrm{keV}}{)}^{0.527}\mathrm{Mpc},\end{array}$$(10)

where the cosmological evolution factor is $E(z)=\sqrt{\mathrm{\Omega }{(1+z)}^3+{\mathrm{\Omega }}_{\mathrm{\Lambda }}}$ and h ≡ H₀/(100kms⁻¹Mpc⁻¹). Here T_x is the global X-ray temperature, estimated using the emission-weighted average of the temperature measured at r > 50 kpc with a single temperature apec model, which includes thermal bremsstrahlung and line emission. This value, obtained from the core-excised total emission, is considered to be a robust proxy of the more formally-defined virial temperature, since it is not affected by the prominent emission in the core where the temperature may decrease significantly.

We divided the inner ∼ 400 kpc into annuli with about 3000 net counts each in the total (0.5–7 keV) band. The spectra were analyzed with XSPEC v12.8.2 (Arnaud 1996). To model the X-ray emission in each ring, we used a single-temperature apec model where the ratio between the elements refers to the solar elemental abundances as in Asplund et al. (2009). Galactic absorption is modeled with tbabs (Wilms et al. 2000), where the Galactic column density, n_H, at the cluster position is initially set as n_H, tot from Willingale et al. (2013). The temperature, metal abundance, and normalization are left as free parameters. The projected temperature and metallicity profiles are given by the spectral fit. The normalization values of the spectra, K, are linked to the 3D density profile through the relation:

$$\begin{array}{c}\hfill K=\frac{{10}^{-14}}{4\pi D_A^2{(1+z)}^2}{\displaystyle \int }{n}_e{n}_H\mathrm{d}V,\end{array}$$(11)

where D_A is the angular diameter distance, n_e and n_H are the 3D density profiles of electrons and hydrogen atoms, respectively, and the volume integral is performed on the projected annulus and along the line of sight. Assuming spherical symmetry, it is possible to obtain the deprojected electron density averaged in the spherical shell corresponding to the projected radius as:

$$\begin{array}{c}\hfill n_{\mathrm{e}}^2(r)=-\frac{4\times {10}^{14}}{0.82\pi }{D}_A^2{(1+z)}^2{\displaystyle {\int }_r^{\infty }}\frac{\mathrm{d}}{\mathrm{d}s}(\frac{K(s)}{s\mathrm{\Delta }s})\frac{\mathrm{d}s}{\sqrt{s^2-r^2}}.\end{array}$$(12)

where r and s are the physical and projected radii, respectively. The full derivation of Eq. (12) is shown in Appendix A. Here, we assumed that the ratio of hydrogen nuclei to electron density is 0.82, which is appropriate for a fully ionized plasma (Ettori et al. 2002). In order to avoid spurious noise amplification, we applied Eq. (12) to the analytical fit to the function $\sqrt{K(s)/s\mathrm{\Delta }s}$, where Δs is the bin width. The fitting function is assumed (for simplicity) to be the single or double β model that is used for the fit of the deprojected density profile, where a single β model component is expressed as:

$$\begin{array}{c}\hfill n_{\mathrm{e}}(r)=n_0{[1+{\left(\frac{r}{r_c}\right)}^2]}^{-\frac{3}{2}\beta },\end{array}$$(13)

where n₀ is the value of the central density. Cool-core clusters show a pronounced peak in the center, with a plateau typically at r < 2 − 3 × 10⁻²r₅₀₀. On the other hand, non-cool-core clusters have a profile that flattens around r ∼ 0.1 × r₅₀₀, and therefore the values of the central density are significantly lower. The uncertainty in the density profiles, corresponding to a formal 1-σ confidence level, is obtained by computing a large ∼10⁴ set of profiles corresponding to a Monte Carlo sampling of the best-fit parameters given their statistical uncertainty.

The 3D temperature profiles for relaxed galaxy clusters can be modeled with an analytical function obtained as the product of two different regimes, corresponding to the core and the outer region, with opposite slopes (see Vikhlinin et al. 2005, 2006):

$$\begin{array}{c}\hfill T_{3\mathrm{D}}(r)=T_0{t}_{\mathrm{inner}}(r)t_{\mathrm{outer}}(r),\end{array}$$(14)

where the function t_outer(r) describe a gentle decrease at large radii of the form:

$$\begin{array}{c}\hfill t_{\mathrm{outer}}(r)=\frac{{(r/r_t)}^{-a}}{{[1+{(r/r_t)}^b]}^{c/b}},\end{array}$$(15)

with a, b, and c defined positive. The central part of the temperature profile instead, requires a function t_inner(r), which is parameterized as follows:

$$\begin{array}{c}\hfill t_{\mathrm{inner}}(r)=(x+T_{\min }/T_0)/(x+1),x={(r/r_c)}^{a_{\mathrm{inner}}},\end{array}$$(16)

with a_inner defined positive. We note that the temperature profile is described by eight free parameters, which we deem are too many to be meaningfully constrained by our profiles. However, here we are interested mostly in the temperature profiles itself, while we are not directly using the parameters r_t and r_c since they are strongly degenerate with the other parameters.

To derive the deprojected temperature profile T_3D, we fit the projected temperature profile with Eq. (14) and then deprojected the best-fit profile numerically through a standard onion-peeling technique (see Ettori et al. 2002, for details), where the temperature in shells is recovered by correcting the spectral estimate with the emission observed in each ring along the line of sight, assuming a spherically symmetric ICM distribution. We assume that this step does not introduce additional errors in the deprojected temperature profiles, so the uncertainty is equal to that of the projected profile. Clearly, this is correct if we assume that the best-fit analytical profile provides an accurate description, and avoids the typical noise amplification that we would obtain with a straightforward deprojection of the measured projected temperature. We argue that our uncertainty on the temperature is only slightly underestimated since we have selected preferentially spherically symmetric clusters, for which we do expect a smooth projected profile.

The mass profile is then obtained from Eq. (6). The temperature, density, and total mass 3D profiles have been already obtained for several clusters in the literature, and therefore are not presented nor discussed here. We compare our profiles with those in the literature when possible, and find reasonable agreement with some residual uncertainty of the order of 10–20% in a few cases, possibly due to the different calibration, different set of exposures used, and different binning. With respect to a simple review of the literature, our profiles have the advantage to be consistently computed with the same assumptions and calibration. With our deprojected 3D electron density profile, temperature profile, and total mass profile, we are now able to derive the timescale profiles defined in Sect. 3.

7. Results

The profiles of the timescales, t_cool, t_eddy, and t_ff were obtained by plugging the temperature, density, and mass 3D profiles into Eqs. (4), (5), and (7). The timescale profiles are shown in Appendix B for all the clusters in our sample. We note that t_ff profiles can be fairly approximated with a power law, while t_cool profiles show more significant changes in slope. We remark that t_cool is very steep in the presence of a strong cool core (with temperature dropping and density increasing toward the center), while it flattens when the core is less prominent. In some cases, t_cool has a slight increase toward the center. We notice that in some cases this is associated with the presence of a central AGN and we conclude that this effect may be uniquely due to the combination of AGN wings that have not been properly removed and of the size of the spatial bin that tend to smooth the central density peak. Since this does not affect our measurements of R_ccc and R_qcf, we do not further comment on the shape of t_cool and t_ff with regard to the initial few spatial bins.

We also note that in a few systems (A1991, PKS0745, RXCJ1504, CL2089, and A1835) we see TI ratios that drop significantly below 10, usually near values of 7–8. This appear to be in contradiction with previous works, where profiles have never observed below the threshold t_cool/t_ff = 10 (see Fig. 7 of Hogan et al. 2017). While the majority of clusters show a minimum value of the TI-ratio in the 10–40 range, we verify that in all the cases where the profile drops below 10 our results are robust and can be ascribed to our updated analysis and calibration methods. Indeed, we recall that hydrodynamical simulations of self-regulated feeding and feedback in clusters do show the TI ratio falling to 5–10 during major CCA triggering periods (e.g., see Fig. 3 in Gaspari et al. 2012b). We argue that as samples expand and become deeper, we expect to identify several more systems to fill such regions where the TI-ratio falls shortly below 10.

Once we fitted the timescale profiles with Polynomial fitting, we compared the best-fit function for t_cool with that for t_eddy and t_ff. We estimate R_ccc and R_qcf at the radii where t_cool crosses t_eddy and 25 t_ff, respectively (see Sect. 2). The best-fit values of R_ccc and R_qcf and the corresponding errors are computed as follows. We extracted randomly 1000 profiles according to the statistical uncertainties on the best-fit parameters of each profile. For each cluster, we obtained a distribution of intersection radii R_ccc and R_qcf. Then, we fit the distribution of R_ccc and R_qcf with a Gaussian, deriving the best-fit value (the center of the Gaussian) and the 1-σ uncertainty. An example is given in Fig. 2 for the measurement of R_ccc in RXCJ1504. Incidentally, we note that some t_cool profiles significantly change slopes (see, e.g., 2A0335, Hydra, RXCJ1558, and PKS0745), so that, in principle, we have two intersections. However, the one close to the center is typically associated with a region of the profile with larger uncertainties. When looking at the distribution of intersection radii obtained by considering all the intersections between the randomized time-scale profiles, we fit only the most prominent peak, so that we naturally discarded these few potentially ambiguous cases where a secondary peak is present. We visually checked that there were no cases where two peaks of comparable significance are present in the distribution of intersection radii.

Fig. 2. Profiles of teddy, tcool and 25 × tff as a function of the physical radius for RXCJ1504 (top panel). Distribution of intersection radii defined by the condition tcool(Rcc) = teddy(Rcc), obtained by sampling random profiles extracted according to the uncertainty on the best-fit profile (tottom panel). The red dashed line shows the Gaussian fit to the distribution used to derive the central value for Rccc and its 1σ uncertainty.

In some other cases, the best fit profiles do not overlap⁶. Nevertheless, the same procedure allows us to obtain an upper limit at a given confidence level. We set the 1-σ upper limit to R_ccc and R_qcf, measuring the value below which we find the 84% of the of the measured values. The best fit values for R_ccc and R_qcf, together with the core-excised, emission weighted temperature (for r > 40 kpc) and total mass with r₅₀₀, are listed in Table 2.

The distributions of R_ccc and R_qcf in our sample are shown in Fig. 3. To obtain the distributions of R_qcf and R_ccc, first we resampled the R_ccc and R_qcf values of each cluster 1000 times by randomly varying the profiles of t_cool, t_eddy, and t_ff according to their uncertainty. Eventually, we averaged the number of values falling in each bin and obtained the final histogram distributions of R_qcf and R_ccc that properly take into account the uncertainties in the timescale profiles.

Fig. 3. Histogram of the novel cool-core radii Rccc and Rqcf, compared with the classical cool-core radius (within which tcool < 7.7 Gyr), as fractions of r500 (top) and in kpc unit (bottom).

We find that the distribution of R_ccc is peaked at ∼0.01 r₅₀₀, and is entirely included within ∼0.07 r₅₀₀. We also notice an apparent bimodality in the distribution of R_qcf, which is peaked at ∼0.04 r₅₀₀, and shows a second peak at 0.07 r₅₀₀. In principle the bimodality in the distribution of R_qcf is expected to be significant, since it is visible after the randomization of the R_qcf values which would have erased the dip at ∼0.05 r₅₀₀ – as this was only due do discreteness effects. However, as we show in the lower panel of Fig. 3, the bimodality disappear when the distribution is plotted as a function of the physical radius. Therefore, we do not investigate further this feature. In addition, we recall that, since our sample selection is biased towards cool-core clusters, some features in the distribution of R_ccc and R_qcf in our sample should not be ascribed to general properties of the cluster population. Finally, we note that the distribution of R_ccc is shifted at lower values by a factor of ∼3 with respect to the distribution of R_qcf.

We note that, historically, the classic cooling radius is often defined as the radius, where t_cool = t_age, and t_age is a typical “age” of the object, roughly estimated as t(z_o)−t(z ∼ 2)∼1/H₀, t is the age of the universe at different redshifts z and z ∼ 2 is assumed to mark the epoch of galaxy cluster formation. Another arbitrary definition of cooling radius is obtained by directly comparing the cooling time to some absolute reference timescale, for instance, 1 Gyr, or a specific value such as 7.7 Gyr (Hudson et al. 2010). These criteria allow one to derive a simple and immediate order-of-magnitude estimate of the actual cooling radius, at variance with our definition based on a local, physical criterion, independent from the cosmic epoch. Therefore, we will refer to the “classical” cool-core definition adopting the fixed threshold t_cool < 7.7 Gyr as the radius R_7.7, and compare it to the physically motivated values of R_ccc and R_qcf. The “classical” cool-core R_7.7 distribution is shifted to larger values with respect to R_ccc and R_qcf, distributed in the broad range of 0.1 − 0.2 r₅₀₀, as shown in Fig. 3. Despite the wide distribution, R_7.7 is clearly disconnected from the physically motivated values for R_ccc and R_qcf. We argue that these values describe more accurately the condensation and the intermittent cooling flow-regions, as we further discuss in the next section (Sect. 8).

In Fig. 4, we investigate in more details the relation between R_qcf and R_ccc. We can directly verify that R_qcf > R_ccc for all the clusters, as expected, with a few cases for which R_qcf ∼ R_ccc within the errors. With the bootstrap fitting method, and using the relation:

$$\begin{array}{c}\hfill \log R_{\mathrm{qcf}}=A\times \log R_{\mathrm{ccc}}+\gamma ,\end{array}$$(17)

Fig. 4. Quenched cooling flow radius Rqcf plotted against the cool-core condensation radius Rccc, including the upper limits in Rqcf and Rccc (green arrows). The solid line is the best-fit power law relation (with 1-σ uncertainty shown by the shaded area), while the dashed line shows Rqcf ∝ Rccc with arbitrary normalization. The color of each point is color coded according to the temperature of each cluster. Error bars correspond to 1-σ confidence level.

we obtain $A=0.{46}_{-0.03}^{+0.02}$ and $\gamma =-0.{41}_{-0.05}^{+0.04}$. We included also the upper limits in this analysis, shown with arrows in Fig. 4. The correlation is statistically very significant, but also shows a considerable intrinsic scatter in addition to the statistical noise. The interesting result is that the relation is significantly different from a linear relation, R_qcf ∝ R_ccc. From this, we conclude that larger quenched cool cores have a larger probability to host a larger cool-core condensation region, implying more vigorous precipitation events. In other words, if we assume that a large cool core (R_qcf ∼ 0.1 r₅₀₀) is the signature of a cooling and feedback activity ongoing for a long time, we may expect that the precipitation condition (and therefore the condensation region) is extended to the entire cool-core region (R_ccc ∼ R_qcf), while in smaller cool cores, the condensation region is progressively smaller (we discuss this further in Sect. 8).

In Fig. 5, we plot the normalized (left panel) and physical (right panel) quenched cool-core radius, R_qcf, and cool-core condenstation radius, R_ccc, versus the total mass. The plots show a substantial scatter with no clear dependence on the mass over the range spanned by our sample, except for a small hint of lower values at larger masses (relatively to r₅₀₀, but not as absolute values). The large intrinsic scatter may reflect the different history of each cluster, with no clear dependence on the total mass (the integrated accretion history) but, rather, on the occurrence of recent major mergers, which have the effect of erasing the cool core and restarting the cooling process and the formation of a new cool core. The smaller R_ccc and R_qcf values at M₅₀₀ ∼ 10¹⁵ M_⊙ are affected by small number statistics and probably by a selection bias: smaller cool cores survive the criterion on the net counts within 40 kpc because they are substantially brighter. On the other hand, smaller clusters must have larger cool cores to pass the selection criteria. Extending our sample to the lower mass, galaxy group will be important to understand potential mass trends, which are nevertheless expected to be weak (e.g., Gaspari et al. 2018).

Fig. 5. Quenched cool core radius Rqcf (red) and cool-core condensation radius Rccc (blue) plotted versus the total mass M500. Error bars correspond to 1-σ confidence level. In the top (bottom) panel the radii are in units of r500 (physical kpc). No dominant trend is observed.

The dependence of R_ccc and R_qcf on redshift is shown in Fig. 6. In both cases, there is an increasing trend with redshift. We do not attempt to interpret this relation in terms of evolution of the cluster population due to the limited size of our sample and our selection technique.

Fig. 6. Quenched cooling flow radius Rqcf (red) and the cool-core condensation radius Rccc (blue) plotted versus redshift. The black dots show the size of the innermost bins, to visualize of how well we are able to resolve the inner regions compared with the Rqcf values.

8. Discussion

Cool cores are one of the most important features of galaxy clusters. As introduced in Sect. 1, they represent the heart of the ICM atmosphere that fill the galaxy cluster potential well. They roughly demark the reservoir of gas out of which the central SMBH can feed recursively over the cosmic evolution (see the multiscale diagram in the review by Gaspari et al. 2020). However, they are often arbitrarily defined, by comparing the cooling time with ad-hoc timescales such as 1 or 7.7 Gyr (see Hudson et al. 2010 for a detailed discussion). Previous studies also focused on visual definitions, such as the peak in the temperature profile (e.g., Vikhlinin et al. 2006), or classical pure or unheated cooling flow rates (e.g., McDonald et al. 2018). Inspired by the results of Gaspari et al. (2018), we investigated here a more physically-driven way to establish the “cooling-flow” region, namely, where the ICM plasma actually condenses, subsequently forming multiphase filaments and clouds that rain onto the SMBHs via CCA (e.g., Tremblay et al. 2018; Rose et al. 2019; Juráňová et al. 2019; Olivares et al. 2019, 2022; Temi et al. 2022). As a result of our investigation, we propose that the ratio of crucial ICM feeding and feedback timescales provides a dimensionless and unbiased way to assess such a cooling region, which can potentially differ from the above definition, namely, that of the classical ad hoc “cool core”. Here, we discuss the implications of our main findings, namely: (1) definition and meaning of the condensation cool-core radius R_ccc; (2) definition and meaning of the quenched cooling flow radius R_qcf; (3) R_ccc − R_qcf relation and constraints on the AGN feeding-feedback duty cycle.

8.1. Definition and meaning of th condensation cool-core radius

We find that a more realistic definition for the roughly spherical region where condensation and precipitation are strongly effective, is provided by the condensation cool-core radius R_ccc. Here, we find typical values r ∼ 0.01 < R_ccc/r₅₀₀ < 0.05, that is, about 5–10× smaller than the above classical definitions of the cool core radius (see upper panel of Fig. 3). This has been achieved by leveraging the C-ratio crossing around unity, which is a simple and clear physical threshold that characterizes the triggering of the turbulent non-linear thermal instability and related top-down multiphase condensation (e.g., Gaspari et al. 2017; Voit 2018; Olivares et al. 2022). In our sample of massive galaxy clusters, C ∼ 1 yields a distribution of R_ccc peaked at ∼20 kpc (see lower panel of Fig. 3), which is comparable to the values found in the upper envelope of the Gaspari et al. (2018) sample (including poor clusters groups as well). It is important to note that such a novel cool-core condensation radius is not an hypothetical region over which a massive pure cooling flow would ensue, but the actual physical region set by the observed balance of feeding and feedback processes. The scattered but rather flat R_ccc − M₅₀₀ relation (see Fig. 5) suggests that more local properties (such as the BCG mass) may drive the condensation core evolution. We defer this issue to a subsequent study.

8.2. Definition and meaning of the quenched cooling flow radius

We further explored a complementary criterion t_cool/t_ff < 25 to trace the region over which a quenched cooling flow may potentially develop from linear TI (Gaspari et al. 2012b; Sharma et al. 2012; Voit et al. 2015b). Despite the phenomenological non-unitary threshold that is still hard to understand in a comprehensive theoretical framework, the related quenched cooling flow radius, R_qcf, is another valuable indication of the potentially condensing region out of the heated macro cooling flow. Although it is affected by larger fluctuations in the distribution function, ranging from 0.02 to 0.13R₅₀₀ (or from 20 to 130 kpc; see Fig. 3), R_qcf tends to be on average larger than R_ccc by factor up to ∼3. We thus find that R_ccc is the inner part of the condensation region affected by direct turbulent precipitation and CCA, while R_qcf represents the more extended “weather” over which we expect the secular development of a quenched cooling flow. We expect each cluster to be in different weather stages of the self-regulated feeding and feedback cycle, oscillating between extended and more localized condensation rainfalls (< 0.05 r₅₀₀). On the other hand, R_qcf only traces the global quenched cooling flow region, up to ≥0.1, r₅₀₀, which does not necessarily imply actual condensation. Indeed, most nebular emission and warm gas is typically contained within smaller radii, which are comparable to R_ccc (e.g., Gaspari et al. 2018). Therefore, in comparison to the quenched cooling flow radius, R_qcf, R_ccc represents a more reliable and stable indicator for the effective cool core in terms of condensed matter, which can be leveraged in simulations and semi-analytical models, as well as in the interpretation of observations.

Incidentally, we also explored an alternative definition of R_qcf, corresponding to the condition t_cool/t_ff = 10. In this case, we find that less than ten clusters would have a well-defined crossing value, while the majority would instead have an undefined value, due to their profiles slowly approaching the condition t_cool/t_ff = 10 in the flat part of the profile. We note, however, that the condition t_cool/t_ff = 10 must be considered as a lower bound to the R_qcf criterion, as we already mentioned. We are aware that the actual threshold is more likely a range of values, rather than a single-value threshold, and that as a consequence, the actual R_qcf distribution may be more scattered than that presented in this work. Nevertheless, we argue that an average ratio t_cool/t_ff ∼ 25 is a reasonable proxy to trace the initial growth of linear thermal instability (TI) in heated cooling flows (see Hogan et al. 2017). Finally, considering a range of values for the threshold, apart from introducing a large scatter, would not change the qualitative aspects of our results.

8.3. R_ccc-R_qcf relation and constraints on the AGN duty cycle

The R_qcf − R_ccc relation can give us further insight into the AGN feeding-feedback duty cycle. Notably, there is no linear relation between the two condensation radii: they tend to become comparable at large 100 kpc values but non-linearly diverge as they approach smaller sizes (Fig. 4). The divergence of such radii toward smaller scales can be interpreted through the micro and meso precipitation having a more flickering duty cycle than the macro (ensemble) weather. Such a trend is qualitatively consistent with the CCA variability found in high-resolution 3D hydrodynamical simulations (Gaspari et al. 2017), showing that the characteristic power spectral frequency of the rains has a negative slope (flicker noise). Notably, such feature has been probed in other recent cluster surveys and observations (e.g., McDonald et al. 2021; Somboonpanyakul et al. 2021). Extrapolating below the cluster regime, the sublinear slope implies that smaller halos are expected to have a larger R_qcf − R_ccc separation; indeed, recent multiwavelength constraints suggest that the raining region is substantially more compact in isolated galaxies (Temi et al. 2022) compared with the global cooling flow zone. Furthermore, we suggest that such differences in “cool-core” size could allow to reconcile better the differences between the cooling rates Ṁ_cool retrieved via imaging out of the classical large cool core with those constrained via spectroscopy, with the latter usually limited within smaller regions similar to R_ccc (e.g., Molendi et al. 2016).

To recap, the above combined evidences suggest a global core radii picture based on R_ccc ≲ R_qcf ≪ R_classic. As shown in Fig. 3, each cluster has a “classical” cool-core region R_classic, which envelopes a sphere with radius up to ∼0.2 − 0.3 r₅₀₀ (or several 100 kpc). However, such R_classic is purely based on an arbitrary cosmological threshold t_cool < 7.7 Gyr (Hudson et al. 2010). Within such classical cool core where most of the X-ray radiation is emitted (but the hot gas is not able to actually condense into the lower gas phases) we find the long-term quenched cool core R_qcf, where feeding and feedback processes balance out secularly, akin to an extended macro weather. Inside such cool core, we find the effective condensation rain and flickering CCA traced via R_ccc and directly tied to the nebular emission and warm or cold gas detections.

9. Conclusions

We computed the profiles of the cooling time t_cool(r), free-fall time t_ff(r), and turbulence cascade timescale t_eddy(r) in 37 massive (2 × 10¹⁴ < M₅₀₀ < 2 × 10¹⁵ M_⊙) galaxy clusters observed with Chandra with high S/N at 0.03 < z < 0.29. The profiles have been obtained from observables entirely derived from X-ray datasets, such as temperature, electron density, and metal abundance profiles. We derive for each cluster a condensation core radius R_ccc defined as the radius where t_cool = t_eddy, and a (quenched) cooling flow radius R_qcf defined as t_cool = 25 × t_ff. We accurately evaluate the statistical uncertainties on these timescales and radii, and explore the distribution of their values across our cluster sample. Our main results are as follows.

• The distribution of the condensation core R_ccc peaks at ∼0.01 − 0.02 r₅₀₀ and entirely included within ∼0.07 r₅₀₀. The distribution of R_qcf is broader and shifted on average to ∼ 3× larger values with respect to the R_ccc distribution.

• We find no significant correlation between total mass and the two radii R_ccc and R_qcf, with a hint of larger values (relatively to r₅₀₀, but not as absolute values) toward low masses (below 3 × 10¹⁴ M_⊙, which shall be further explored with larger samples including galaxy groups. Instead, both and R_ccc appear to moderately increase with redshift.

• Supported by theoretical models and high-resolution hydrodynamical simulations of the multiphase condensation rain and CCA (e.g., Gaspari et al. 2017, 2018), we find that R_qcf is a tracer of the extended quenched cooling flow or macro weather. The smaller cool-core condensation radius R_ccc traces instead the actual inner rain, which also drives CCA feeding episodes onto the SMBH. Both proposed core radii are an order of magnitude smaller than the classical and ad-hoc cool-core definition R_7.7 Gyr, where the halo is emitting strong X-ray radiation, but not condensing or inflowing.

• We find that the correlation between R_qcf and R_ccc is non-linear, diverging at lower values up to a factor of 3 over an order of magnitude in r₅₀₀. This relation between the two scales allows us to infer some features of the duty cycle and related appearance or disappearance of the local cooling flow. Specifically, the slope of the R_qcf − R_ccc relation is measured to be $0.{46}_{-0.03}^{+0.02}$, and suggests that the micro CCA rain is flickering on and off in the macro weather, in agreement with hydrodynamical simulations (Gaspari et al. 2017).

As shown above, to describe the fate of the multiscale cooling gas and its interplay with the cluster atmosphere, we need a combination of high spatial resolution and good spectroscopic capabilities. For instance, among the future generation of X-ray instruments, those onboard Athena will play an important role to help in this direction. It is expected to achieve a PSF below 10 arcsec over the entire field of view and an effective area at 1 keV a factor > 5 larger than the current detectors, along with the Wide Field Imager (WFI) providing sensitive wide field imaging and spectroscopy and the X-ray Integral Field Unit (X-IFU) delivering spatially resolved high-resolution X-ray spectroscopy – altogether, these instruments will allow for an unprecedented level of accuracy in mapping the state of the gas within 0.1 r₅₀₀ and cool cores in systems at T > 2 keV up to a redshift of ≲0.4. Further, high angular-resolution X-ray facilities such as Lynx⁷ and AXIS⁸ are proposing to investigate with subarcsecond resolution over a FoV of 400–500 arcmin² the X-ray sky, improving this capability of a factor ∼100 with respect to Chandra ACIS-I. All such core missions will enable crucial steps forward in the characterization and understanding of the weather in galaxy cluster cores by linking the micro, meso, and macro scales as well as above key condensation processes.

Acknowledgments

We thank the Referee for their constructive criticism and helpful suggestions. This work was supported by the Bureau of International Cooperation, Chinese Academy of Sciences under the grant GJHZ1864. We acknowledge financial contribution from the agreement ASI-INAF n.2017-14-H.0. P.T. acknowledges financial support through grant PRIN-MIUR 2017WSCC32: “Zooming into dark matter and proto-galaxies with massive lensing clusters”. M.G. acknowledges partial support by NASA HST GO-15890.020/023-A; this work is part of the broader BlackHoleWeather program. S.E. acknowledges financial contribution from the contracts ASI-INAF Athena 2019-27-HH.0, INAF mainstream project 1.05.01.86.10, and funding from the European Union Horizon 2020 Programme under the AHEAD2020 project (grant agreement n. 871158). The Chandra raw data used in this paper are available to download at the Chandra Data Archive website (https://cxc.cfa.harvard.edu/cda/). The reduced data are also available upon request from the corresponding author.

References

Appendix A: Derivation of the deprojected electron density (Equation 12)

When fitting the spectra of the projected annuli, we obtain the projected, emission-weighted temperature and metallicity, as well as the projected normalization of the emission measure, corresponding to the normalization of the spectrum. The normalization parameter K of of the mekal spectrum in Xspec is defined as:

$$\begin{array}{c}\hfill K=\frac{{10}^{-14}}{4\pi D_A^2{(1+z)}^2}{\displaystyle \int }{n}_e(r)n_H(r)\mathrm{d}V\sim \frac{{10}^{-14}0.82}{4\pi D_A^2{(1+z)}^2}{\displaystyle \int }{n}_e{(r)}^2\mathrm{d}V,\end{array}$$(A.1)

where 0.82 is the cosmic ratio, n_H/n_e, and n_e(r) is the radial density profile of the electron density, which is the quantity we want to derive from the observed (projected) values. In the assumption of spherical symmetry, the differential volume element dV corresponding to a projected annulus of finite width, Δs, depends on the projected radius, s, and on the coordinate z along the line of sight, simply as dV = 2πsΔsdz. Since $z=\sqrt{(r^2-s^2)}$, we obtain:

$$\begin{array}{c}\hfill \mathrm{d}V=2\pi s\mathrm{\Delta }s\frac{r\mathrm{d}r}{\sqrt{r^2-s^2}}.\end{array}$$(A.2)

If we define the function f(r) as

$$\begin{array}{c}\hfill f(r)\equiv \frac{0.82{10}^{-14}}{4D_A^2{(1+z)}^2}{n}_e^2(r),\end{array}$$(A.3)

then we have rewritten K as a function of the projected radius, s, in a compact expression:

$$\begin{array}{c}\hfill \frac{K(s)}{s\mathrm{\Delta }s}=2{\displaystyle {\int }_s^{\infty }}\frac{rf(r)\mathrm{d}r}{\sqrt{r^2-s^2}}.\end{array}$$(A.4)

If we apply the Abel transform, we can write:

$$\begin{array}{c}\hfill f(r)=-\frac{1}{\pi }{\displaystyle {\int }_r^{\infty }}\frac{d}{\mathit{ds}}(\frac{K(s)}{s\mathrm{\Delta }s})\frac{\mathrm{d}s}{\sqrt{s^2-r^2}}.\end{array}$$(A.5)

This allows us to write the square of the electron density with an expression based on observables. Combining Equations A.3 and A.5, we obtain

$$\begin{array}{c}\hfill n_e^2(r)=-\frac{4\times {10}^{14}}{0.82\pi }{D}_A^2{(1+z)}^2{\displaystyle {\int }_r^{\infty }}\frac{d}{\mathit{ds}}(\frac{K(s)}{s\mathrm{\Delta }s})\frac{\mathrm{d}s}{\sqrt{s^2-r^2}},\end{array}$$(A.6)

which is the same as Equation 12 given in Section 6.

Appendix B: Timescale profiles

In Fig. B.1, we show the profiles of the cooling time, t_cool, turbulence eddy time, t_eddy, and free-fall time, t_ff. Each row refers to a cluster following the same order of Table 2. In the left panels, we show the points corresponding to the spatial bins of the spectral fits, with the points obtained from the best-fit deprojected temperature and density profiles for t_cool and t_ff. In addition, we plot the polynomial best fit to the time scale profiles as a continuous line, showing the 1-σ uncertainty with a shaded area. In the right column, we show the ratio of the best-fit function of t_cool, t_eddy, and t_ff. In particular, t_cool/t_ff defines R_qcf when t_cool/t_ff = 25, while t_cool/t_eddy defines R_ccc when t_cool/t_eddy = 1. In Fig. B.2, we can see that the results of using Equation 8 and Equation 9 to find R_ccc are basically the same, and the error between Equation 9 and the real data are also small, so we use Equation 9 to calculate R_qcf for those galaxy clusters without measurement data.

Fig. B.1. Profiles of the cooling time tcool, the turbulence eddy turnover time teddy, and the free fall time tff (left). Ratio of tcool/tff that defines Rqcf, when tcool/tff = 25, and tcool/teddy, which defines Rccc when tcool/teddy = 1 (right).

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.1. continued.

Fig. B.2. Comparison chart of the results(16 galaxy clusters) of calculating Rccc with Equation 8 and measured data (left). Comparison chart of the results(16 galaxy clusters)of calculating Rccc with Equation 9 and measured data (right). Comparison chart of the results of calculating Rccc with Equations 8 and 9 (bottom).

All Tables

All Figures

Fig. 1. Top panel: Turbulent velocity values measured in the warm or cold phase by Gaspari et al. (2018, blue points) and Olivares et al. (2019, orange points). The reference dotted line σ ∝ T1/3 corresponds to Eq. (9), with the cyan-shaded area showing the 1σ uncertainty. The solid lined corresponds to the constant average value given by Eq. (8), with the orange-shaded area showing its 1σ uncertainty. Bottom panel: the histogram distribution of the clusters with σv measurement (blue) and of the complementary subsample (green).

In the text

Fig. 2. Profiles of teddy, tcool and 25 × tff as a function of the physical radius for RXCJ1504 (top panel). Distribution of intersection radii defined by the condition tcool(Rcc) = teddy(Rcc), obtained by sampling random profiles extracted according to the uncertainty on the best-fit profile (tottom panel). The red dashed line shows the Gaussian fit to the distribution used to derive the central value for Rccc and its 1σ uncertainty.

In the text

	Fig. 3. Histogram of the novel cool-core radii Rccc and Rqcf, compared with the classical cool-core radius (within which tcool < 7.7 Gyr), as fractions of r500 (top) and in kpc unit (bottom).
In the text

Fig. 4. Quenched cooling flow radius Rqcf plotted against the cool-core condensation radius Rccc, including the upper limits in Rqcf and Rccc (green arrows). The solid line is the best-fit power law relation (with 1-σ uncertainty shown by the shaded area), while the dashed line shows Rqcf ∝ Rccc with arbitrary normalization. The color of each point is color coded according to the temperature of each cluster. Error bars correspond to 1-σ confidence level.

In the text

	Fig. 5. Quenched cool core radius Rqcf (red) and cool-core condensation radius Rccc (blue) plotted versus the total mass M500. Error bars correspond to 1-σ confidence level. In the top (bottom) panel the radii are in units of r500 (physical kpc). No dominant trend is observed.
In the text

	Fig. 6. Quenched cooling flow radius Rqcf (red) and the cool-core condensation radius Rccc (blue) plotted versus redshift. The black dots show the size of the innermost bins, to visualize of how well we are able to resolve the inner regions compared with the Rqcf values.
In the text

	Fig. B.1. Profiles of the cooling time tcool, the turbulence eddy turnover time teddy, and the free fall time tff (left). Ratio of tcool/tff that defines Rqcf, when tcool/tff = 25, and tcool/teddy, which defines Rccc when tcool/teddy = 1 (right).
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.1. continued.
In the text

	Fig. B.2. Comparison chart of the results(16 galaxy clusters) of calculating Rccc with Equation 8 and measured data (left). Comparison chart of the results(16 galaxy clusters)of calculating Rccc with Equation 9 and measured data (right). Comparison chart of the results of calculating Rccc with Equations 8 and 9 (bottom).
In the text