A&A 491, 71-87 (2008)
DOI: 10.1051/0004-6361:200809739

Total mass biases in X-ray galaxy clusters

R. Piffaretti - R. Valdarnini

SISSA/ISAS, via Beirut 4, 34014 Trieste, Italy

Received 7 March 2008 / Accepted 31 July 2008

Abstract
Context. The exploitation of clusters of galaxies as cosmological probes relies on accurate measurements of their total gravitating mass. X-ray observations provide a powerful means of probing the total mass distribution in galaxy clusters, but might be affected by observational biases and rely on simplistic assumptions originating from our limited understanding of the intracluster medium physics.
Aims. This paper is aimed at elucidating the reliability of X-ray total mass estimates in clusters of galaxies by properly disentangling various biases of both observational and physical origin.
Methods. We use N-body/SPH simulation of a large sample of $\sim$100 galaxy clusters and investigate total mass biases by comparing the mass reconstructed adopting an observational-like approach with the true mass in the simulations. X-ray surface brightness and temperature profiles extracted from the simulations are fitted with different models and adopting different radial fitting ranges in order to investigate modeling and extrapolation biases. Different theoretical definitions of gas temperature are used to investigate the effect of spectroscopic temperatures and a power ratio analysis of the surface brightness maps allows us to assess the dependence of the mass bias on cluster dynamical state. Moreover, we perform a study on the reliability of hydrostatic and hydrodynamical equilibrium mass estimates using the full three-dimensional information in the simulation.
Results. A model with a low degree of sophistication such as the polytropic $\beta$-model can introduce, in comparison with a more adequate model, an additional mass underestimate of the order of $\sim$$10 \%$ at $r_{{\rm 500}}$ and $\sim$$15 \%$ at $r_{{\rm 200}}$. Underestimates due to extrapolation alone are at most of the order of $\sim$$10 \%$ on average, but can be as large as $\sim$$50 \%$ for individual objects. Masses are on average biased lower for disturbed clusters than for relaxed ones and the scatter of the bias rapidly increases with increasingly disturbed dynamical state. The bias originating from spectroscopic temperatures alone is of the order of $10 \%$ at all radii for the whole numerical sample, but strongly depends on both dynamical state and cluster mass. From the full three dimensional information in the simulations we find that the hydrostatic equilibrium assumption yields masses underestimated by $\sim$10-$15 \%$ and that masses computed by means of the hydrodynamical estimator are unbiased. Finally, we show that there is excellent agreement between our findings, results from similar analyses based on both Eulerian and Lagrangian simulations, and recent observational work based on the comparison between X-ray and gravitational lensing mass estimates.

Key words: galaxies: clusters: general - X-rays: galaxies: clusters - methods: numerical

   
1 Introduction

Galaxy clusters are the largest virialized structure known in the universe. According to the hierarchical clustering model of structure formation, they form by the gravitational collapse of the rare peaks of the primeval density field, on scales of the order of $\sim$10 Mpc. Within this scenario their formation and evolution is a sensitive function of the cosmological matter density parameter $\Omega_{{\rm m}}$ and the mass fluctuation amplitude $\sigma_{\rm 8}$, where $\sigma_{\rm 8}$ is the rms linear fluctuation on scales $8 ~ h^{{\rm -1}}$ Mpc.

Measurements of their evolution rate can be used to asses the growth in mass of such structures, thereby providing a powerful method to constrain the geometry and matter content of the universe (see Voit 2005, and references therein). Moreover, because of the spatial extent of the collapse scale, the cluster baryonic fraction $f_{{\rm b}}$ is expected to be close to the cosmic value $\Omega_{{\rm b}}/ \Omega_{{\rm m}}$. Measurements of $f_{{\rm b}}$ at high redshifts can be used to derive constraints on the equation of state of the dark energy (Allen et al. 2008; Majumdar & Mohr 2004; Haiman et al. 2001).

The importance of cluster of galaxies as cosmological probes will be further strengthened with the forthcoming high redshift surveys. This is of particular relevance in the new era of precision cosmology, in which studies of cluster evolution will provide independent tests with which to compare constraints on cosmological models extracted from observations of the cosmic background radiation (e.g. Spergel et al. 2007) and distance measurements of high redshift supernovae (e.g. Riess et al. 2007; Tonry et al. 2003; Riess et al. 2004).

From the scenario here outlined it follows that in order to use cluster of galaxies as cosmological probes it is crucial to accurately measure their baryonic and total gravitating mass. The methods used to derive cluster masses are mainly based on the velocity dispersion of the optical galaxy populations (e.g. Biviano & Girardi 2003; Rines et al. 2003), observations of the X-ray emitting intracluster medium (ICM) (e.g. Arnaud et al. 2005; Reiprich & Böhringer 2002; Ettori et al. 2002; Finoguenov et al. 2001; Vikhlinin et al. 2006) and on gravitational lensing (e.g. Mahdavi et al. 2008; Smith et al. 2002).

Accurate mass estimates derived from X-ray data are based on the assumptions that both the total potential and the ICM distribution are spherically symmetric and that the ICM is in hydrostatic equilibrium in the cluster potential well. The latter assumption is usually justified by the fact that the estimated ICM sound crossing times are short when compared against cluster ages (Sarazin 1986).

Under these assumptions the ICM is a faithful tracer of the underlying matter distribution and the total mass profile can be determined from the gas radial density and temperature profiles. The density profile is recovered by deprojecting the surface brightness profile, as measured from X-ray flux maps, whereas knowledge of the temperature profile requires the availability of spatially resolved spectroscopy. High quality data taken by the present generation of X-ray telescopes, Chandra and XMM-Newton, allow for nearby clusters accurate measurements of these quantities out to a large fraction of the cluster virial radius (Pratt et al. 2007; Sanderson et al. 2006; Markevitch et al. 1998; Piffaretti et al. 2005; Leccardi & Molendi 2008; Vikhlinin et al. 2005).

The reliability of cluster mass estimated through the X-ray method can be accurately studied using N-body/hydrodynamical simulations. The principal benefit over analytical methods is that the gas evolution can be treated self-consistently. The validity of the numerical approach is supported by X-ray observations, which show the existence of complex thermal structures and of merging activity (see Markevitch & Vikhlinin 2007, for a review).

Since early pioneering studies, hydrodynamical simulations have become a widely used tool to investigate cluster formation and evolution in different cosmological scenarios (cf. Voit 2005). The numerical resolution of the simulations and the modeling of the cluster gas physics has been improved over the years. The latter now incorporates radiative cooling (Davé et al. 2002; Lewis et al. 2000; Yoshikawa et al. 2000; Pearce et al. 2000; Muanwong et al. 2002), metal enrichment of the ICM by supernovae and energy feedback (Kay et al. 2004; Tornatore et al. 2003; Kravtsov et al. 2005; Kay et al. 2003; Borgani et al. 2004; Valdarnini 2003).

The accuracy of cluster X-ray mass estimators has been tested by means of N-body/hydrodynamical simulations in a variety of papers (Kay et al. 2004,2007; Evrard 1990; Nagai et al. 2007a; Evrard et al. 1996; Rasia et al. 2006). Evrard (1990) first pointed out the existence of a significant bias in the binding mass estimates when using the isothermal $\beta$-model. Evrard et al. (1996) confirmed that the source of this discrepancy is related to the isothermal and hydrostatic assumptions. The lack of validity of the hydrostatic assumption is observationally motivated by optical and X-ray maps, which show the existence of substructure with an ongoing merger activity, and is numerically supported by a number of authors (e.g. Kay et al. 2004; Nagai et al. 2007a; Rasia et al. 2004), who found that in the simulations the ICM is not perfectly in hydrostatic equilibrium. This implies the presence of residual gas bulk (laminar) and turbulent (random kinetic) motion, and leads to an underestimate of the masses because of additional non-thermal pressure support which is not accounted for by the hydrostatic equilibrium equation.

With respect to the isothermal $\beta$-model the modeling of the ICM has been significantly improved (e.g. Vikhlinin et al. 2006) through observational progress, which showed the existence of temperature profiles declining with radius (Pratt et al. 2007; Sanderson et al. 2006; Piffaretti et al. 2005; De Grandi & Molendi 2002; Vikhlinin et al. 2005). These features are well reproduced out of the core radii in hydrodynamic simulations which incorporate cooling and feedback (Tornatore et al. 2003; Kay et al. 2003; Borgani et al. 2004; Muanwong et al. 2002; Valdarnini 2003).

In order to properly assess the reliability of cluster X-ray mass estimators it is however necessary to construct mock observations of simulated clusters which must reproduce spectroscopic measurements as expected from X-ray telescopes. This is motivated by the presence of complex thermal structures in the ICM, which bias the (measured) spectral fit temperatures towards lower values than the average emission weighted cluster temperatures defined theoretically (Mazzotta et al. 2004; Mathiesen & Evrard 2001; Valdarnini 2006). The dependence of X-ray mass estimators on spectral biases and other systematics has been investigated through hydrodynamic simulations by a number of authors (Kay et al. 2007; Nagai et al. 2007a; Jeltema et al. 2008; Rasia et al. 2006). Nagai et al. (2007a) argued that mass estimates are biased low ($\sim$5-$20 \%$) even for clusters identified as relaxed.

In order to properly disentangle observational biases that arise from spectroscopic measurements from those due to incomplete relaxation of the gas or from the ones caused by an inaccurate modeling of the radial profiles, it is however necessary to derive X-ray mass estimates from a large sample of simulated clusters. This is the main goal of this paper, in which we apply different X-ray mass estimators to a large set of high-resolution hydrodynamical simulations of galaxy clusters. The physics of the gas includes radiative cooling, star formation, chemical enrichment and energy feedback. The sample comprises $\sim$100 clusters, the size of the sample being a critical quantity in order to extract sub-samples large enough to give meaningful statistics.

We discuss the dependence of the mass bias at different radii upon the adopted analytical models and the chosen radial range used to perform fits of the profiles, the cluster dynamical state as well as the impact on the mass bias which follows from the use of spectral temperatures. As a statistical indicator to quantify the cluster dynamical state through the analysis of X-ray maps we adopt the power ratio method (see Jeltema et al. 2005, and references therein). We also investigate the limit of applicability of the dynamical equilibrium equation when used to recover the cluster true masses in the presence of significant non-thermal gas pressure.

The paper is organized as follows. In Sect. 2 we present the procedure for simulating the cluster sample. In Sect. 3 we describe how we generate and analyze the synthetic X-ray observations that are used in Sect. 4 to recover the total mass distribution. In Sect. 5 we investigate the reliability of the hydrostatic and hydrodynamical mass estimators from the full three-dimensional information provided by the simulations. Finally, we discuss our main results and present our conclusions in Sect. 6.

   
2 N-body/SPH simulations

The considered cosmological model assumes a flat CDM universe, with matter density parameter $\Omega_{\rm m}=0.3$, $\Omega_{\rm\Lambda}=0.7$, $\Omega_{\rm b}=0.0486$ and h = 0.7 is the value of the Hubble constant in units of 100 km ${\rm s}^{-1}~{\rm Mpc}^{-1}$. The power spectrum has been normalized to $\sigma_{\rm 8}=0.9$ on a 8   h^-1 Mpc scale at the present epoch and the primeval spectral index n is set to 1.

The simulation ensemble of galaxy clusters is constructed according to the following procedures. A low-resolution N-body run is first performed starting from an initial reshift $z_{\rm i}=10$in a box of comoving size L, using a P3M code with $N_{\rm p}$ particles. Clusters of galaxies are identified at z=0 using a friends-of-friends (FoF) algorithm. The virial mass and radius are related by $M_{\rm vir}= (4 \pi/3) ~ \Omega_{\rm m} ~ \rho_{\rm c} ~ \Delta_{\rm c} ~ r_{\rm vir}^3$, where $\Delta_{\rm c} \simeq 187 ~ \Omega_{\rm m}^{-0.55}$ for a flat cosmology and $\rho_{\rm c}$ is the critical density. This is defined as $\rho_{\rm c}(z)=3 H(z)^{2} /8 \pi G$, where $H(z)^2=H_{\rm0}^2 E(z)^2$ and $E(z)^2=\Omega_{\rm m} (1+z)^3 + \Omega_{\rm\Lambda}$. In general, the fiducial radius $r_{{\rm\Delta}}$ is defined such that the average of the total density within that radius is $\Delta$ times the critical density, i.e. M(< $r_{{\rm\Delta}})= {4 \pi r_{{\rm\Delta}}^3} ~ \Delta ~ \rho_{\rm c}(z) /3$.

The simulation ensemble is constructed by combining two distinct samples S1 and S2. Table 1 lists the most relevant simulation parameters of the two samples S1 and S2.

Table 1: Main parameters of the two cluster samples S1 and S2: comoving box size L of the cosmological N-body runs, number of corresponding particles $N_{\rm p}$, number $N_{\rm cl}$ of the most massive clusters extracted from the sample at z=0, and virial mass $M_{\rm vir}$ of the most (least) massive cluster of sample S1(S2).

A sample is constructed by extracting at z=0 the $N_{\rm cl}$ most massive clusters which are found in a cosmological run of size L and $N_{\rm p}$ particles. For sample S1 the box size and the number of particles are twice and eight times those of sample S2. The value of $N_{\rm cl}$ for sample S1 is chosen such that the last cluster of the sample has its virial mass above that of the first cluster of sample S2. The mass range for the virial masses of the combined samples spans a decade. The random realization of the initial density perturbations are different in the two cosmological simulations.

After completition of the cluster selection, the clusters of the ensemble are resimulated individually using high-resolution hydrodynamic simulations in physical coordinates. The initial conditions of each hydrodynamic simulation are set as follows. The particles of a cluster lying at z=0 within a distance $r_{\rm vir}$ from the cluster center are tagged and located back at an initial redshift $z_{\rm in}=49$ and a cubic region of size $L_{\rm c} \simeq 25$- $50 ~ {\rm Mpc} \propto M_{\rm vir} ^{1/3}$, which contains these particles, is placed at the cluster center. The original Fourier modes of the cosmological simulation are then used to perturb the positions of a uniform lattice of $N_{\rm L}=51^3$ gas particles set inside the cube and high frequency waves are added to the original modes to sample the new Nyquist frequency. The positions of dark matter particles are found similarly. Particles with perturbed positions which lie within a inner sphere of radius $L_{\rm c}/2$ from the cluster centre are kept for the hydrodynamical simulations. For these particles the masses are set in proportion to $\Omega_{\rm b}$ and $(\Omega_{\rm m} - \Omega_{\rm b})$ for gas and dark matter particles, respectively. Tidal fields are modeled by adding to the inner particles an external shell of low-resolution dark matter particles. The shell has an outer radius $L_{\rm c}$ and the mass of a particle is 8 times the sum of the masses of a gas and dark matter particle of the inner region.

The hydrodynamic simulations are run using a multistep TREESPH code in which the gas entropy is explicitly conserved (Goodman & Hernquist 1991). The simulations contain $\simeq$70 000 gas and dark matter particles in the inner region and a similar value of low-resolution dark matter particles in the outer shell. The mass of the gas particles ranges from $m_{\rm g} \simeq 5 \times 10^9~M_{\odot}$ for the most massive cluster of the ensemble, down to $m_{\rm g} \simeq 6 \times 10^8~M_{\odot}$ for the least massive cluster. This mass resolution can be considered adequate for the present purposes, as suggested from analyzing the stability of gas profiles of simulated clusters (Valdarnini 2002). The gravitational softening parameter of gas particles is set to $\varepsilon_{\rm g} =25$ kpc and 15 kpc, for clusters of sample S1 and S2, respectively. For the dark matter particles the softening is rescaled according to $\varepsilon_{i} \propto m_{i}^{1/3}$, where m_i is the mass of the particle i. The softenings are comoving out to z=20, after which they are kept fixed in physical coordinates.

The physical modeling of the gas includes radiative cooling, which depends on the gas temperature and metallicity. Cold gas in high density regions is subject to star formation and gas particles are eligible to form star particles. At each timestep gas particles neighboring a star particle are heated by supernova (SN) explosions of type II and Ia. The gas is also metal enriched through SN explosions. The energy and metal feedback are calculated according to the stellar lifetime and initial mass function. A detailed description of the feedback recipes is given in Valdarnini (2003). The hydrodynamic variables of a simulated cluster are stored at run time at various redshifts in the interval from z=2 down to z=0.

   
3 Simulation and analysis of X-ray observations

At z=0 we extract a total of 153 clusters. For this sample we compute spectroscopic-like global temperatures and select objects with $T_{\rm sl}(<$ $r_{\rm 500}) \geq 2$ keV (see Sect. 3.4 below for the definition of $T_{\rm sl}$). This temperature selection is adopted because when generating spectroscopic temperature profiles we use the approximation derived by Mazzotta et al. (2004), which was developed for continuum-dominated spectra and is therefore not accurate for low temperature systems. Our final sample comprises $\sim$100 temperature-selected clusters spanning a mass range of $8.2 \times 10^{13} ~ h^{-1 }~M_{\odot} \lesssim M(<$ $r_{\rm 200}) \lesssim 1.2 \times 10^{15} ~ h^{-1 }~M_{\odot}$.

For each simulated cluster we generate three independent observations by considering three orthogonal projections, thus obtaining a total of $\sim$300 observed objects. In the following sections we describe how we generate and analyze these mock X-ray observations. Our analysis may be viewed as complementary to those presented in Rasia et al. (2006) and Nagai et al. (2007a). The main difference is the large size of the sample used here, which allows us to achieve results of much higher statistical significance.

   
3.1 X-ray maps

Simulated surface brightness maps are obtained by first choosing a line of sight and then locating the origin at the cluster center in the plane orthogonal to the line of sight; the latter is defined as the location where the gas density reaches its maximum value. Throughout the paper we assume that this position also coincides with the peak of the X-ray emission.

For a given line of sight and redshift z the X-ray surface brightness is defined as:

$$S_{\rm X} (x,y) = \frac{1}{4 \pi (1+z)^4}\int_{-\infty}^{\infty} \varepsilon_{\rm X} {\rm d} l,$$ (1)

where the integral is along the line of sight, x and y are Cartesian coordinates on the chosen plane and $\varepsilon_{\rm X}$ is the X-ray emissivity. The latter quantity can be calculated as $\varepsilon_{\rm X} = \Lambda ~ \rho_{\rm g}^2({\vec x})$, which differs from the true emissivity aside from a constant factor. Here $\rho_{\rm g} ({\vec x})$ is the gas density, $\Lambda = \Lambda(T,Z,{\vec x},E_1,E_2)$ is the cooling function, and T and Z are the gas temperature and metallicity, respectively. The quantities E₁ and E₂ define the energy band [E₁-E₂] used in the X-ray flux measurement. Because of the Lagrangian nature of SPH simulations, the emissivity $\varepsilon_{\rm X}({\vec x})$ is expressed as a summation over particles:

$$\varepsilon_{\rm X} ({\vec x}) = \sum_j \frac{m_j} {\rho_{{\... ...2) \rho_{{\rm g},j}^2 ({\vec x}) W({{\vec x}-{\vec x}_j},h_j);$$ (2)

here the subscript j denotes the values of the quantity at the position of the particle j, W is the SPH smoothing kernel, and h_j the smoothing length of particle j. From this expression for $\varepsilon_{\vec X} ({\bf x})$ it is possible to evaluate the X-ray surface brightness $S_{\rm X}(x,y)$ by performing the integral (1) along the line of sight axis, $\tilde{z}$. For our purposes the integral (1) is replaced by

$$S_X (x,y) = \frac{1}{4 \pi (1+z)^4} \varepsilon_{\rm X}^{2{\rm D}}(x,y),$$ (3)

where $\varepsilon^{2{\rm D}}$ is defined according to Eq. (2) in which the kernel is the 2D Gaussian kernel $W_{2{\rm D}}={\rm exp}[-(\sqrt{ x^2+y^2}/h_j^{\rm G})^2]/ \pi h_{j}^{{\rm G} \;2}$ and $h_j^{\rm G} \simeq 0.82 h_j$.

For each simulated object we generate three surface brightness maps by choosing three orthogonal projections. Motivated by observations we evaluate the maps in the soft energy band [0.5-2] keV. Since particles with temperatures below 10⁵ K do not contribute to the gas emissivity in this energy band, they are not taken into account in the computation.

Cold and dense ICM clumps can affect total mass reconstructions because they produce pronounced irregular features in azimuthally averaged brightness and temperature profiles, and introduce a systematic bias in the spectral temperature determination. We therefore remove small scale clumps from the constructed flux images by performing a masking procedure as done with real X-ray observations. This is accomplished by generating [0.5-2] keV maps on a $1024 \times 1024$ grid of size $2.1 \times r_{200}$. The pixel size of the maps is $\sim$5 and 2 kpc for the most and least massive clusters in our sample, respectively. Pixel sizes are similar to those adopted in observations and are sufficiently small so as to allow a complete identification of all small-scale clumps resolved in our simulations. Clumps are detected and removed following the procedure implemented by Vikhlinin et al. (1998): images are decomposed into wavelets of pixel scale 2ⁿ, with $n=0,1,\ldots 6$, with the significance threshold set to 5 in units of the rms level. Particles that lie in those pixels which are part of a statistically significant structure on a given scale are tagged and identified as part of a clump. These particles are then removed from all the SPH summations. In three dimensions, clumps are identified by those particles that are part of a two dimensional clump in any of the planes. The wavelet algorithm commonly also detects the central peak of the surface brightness distribution. In order to avoid the unwanted masking of the cluster emission in the central region we exclude the inner $0.2 \times r_{200}$ circular aperture from our masking procedure. In addition, we do not impose any limiting flux $f_{\rm X}$ for the detection of small scale clumps. We motivate this choice with the finding that setting $f_{\rm X}$ to $\sim$ $3 \times 10^{-15}~{\rm erg}~/{\rm cm}^{-2}~{\rm s}^{-1}$, as in Nagai et al. (2007a), has negligible effect. Finally, we compute the amount of gas mass removed by the masking procedure and find that even in the most clumpy clusters it is less than a few percent of the total gas mass. If not stated explicitly, all the results presented in the following are derived by including the masking procedure.

   
3.2 Power ratios

Total mass determinations from X-ray observations can be heavily affected by the cluster dynamical state since the method relies on the assumption of hydrostatic equilibrium. The dynamical state of clusters is related to the amount of substructure present in their X-ray surface brightness distribution (e.g. Richstone et al. 1992; Evrard et al. 1993) and various statistical measures have been proposed to quantify cluster substructure (e.g. Buote 2002, and references therein).

In this work we adopt the power ratio method (Buote & Tsai 1995) as a statistical indicator of the cluster dynamical state, since it is widely used to study cluster X-ray morphologies (e.g. Buote & Tsai 1996). This method is expect to be statistically significant when applied to a large cluster sample such as the one studied here. It is in fact unaffected only by mergers along the line of sight, which rarely occur.

According to the power ratio method, the X-ray surface brightness map $S_{\rm X}(\rho,\varphi)$, where $(\rho,\varphi)$ are the conventional polar coordinates, is the source term of the pseudo potential $\Psi(\rho,\varphi)$ which satisfies the 2-D Poisson equation. The pseudo potential is expanded into plane harmonics and the mth coefficients of the expansion are given by:

$$\alpha_m = \int_{R^{'}\leq R_{\rm ap}} {\rm d} ^2 x^{'} S_X(\vec x^{'}) {R^{'}}^m \cos(m\varphi{'}),$$ (4)

$$\beta_m = \int_{R^{'}\leq R_{\rm ap}} {\rm d} ^2 x^{'} S_X(\vec x^{'}) {R^{'}}^m \sin(m\varphi{'}),$$ (5)

where $\vec x{'}=(\rho,\varphi)$ and the integration is over a circular region with aperture radius $R_{\rm ap}$. The mth power ratio is then defined as

$$\Pi_m (R_{\rm ap}) = \log_{10} (P_m/P_0),$$ (6)

where

$$P_m (R_{\rm ap}) {1 \over 2 m^2} (\alpha_m^2 + \beta_m^2) ~~m>0,$$ = (7)

$$[ \alpha_0 \ln(R_{\rm ap}/{\rm kpc}) ]^2.$$ P₀ = (8)

The power ratios $\Pi_{\rm m}(R_{\rm ap})$ are then indicators of the amount of structure present on the scale of the aperture radius $R_{\rm ap}$. The values of $P_{\rm m}$ depend on the choice of the coordinate system. For a fully relaxed configuration $\Pi_{\rm m}\rightarrow - \infty$. As $\Pi _3$ indicates asymmetric distributions, we adopt it here as a measure of the amount of substructure present in a cluster.

For each of the $\sim$300 [0.5-2] keV surface brightness maps (see Sect. 3.1) we set the origin of coordinates at the peak of the X-ray emission. In addition to being fully consistent with the procedure employed in the computation of azimuthally averaged brightness and temperature profiles, this choice of coordinate system enables us to detect $P_{\rm 3}$ values different from zero even in the case of bimodal clusters with nearly equal size components. For the same configuration all the odd moments would vanish if the coordinate system was the frequently adopted centroid frame, which is the coordinate system defined such that $P_{\rm 1}=0$ (Buote & Tsai 1995). We utilize unmasked maps, since the goal of the power ratio analysis is to measure substructure, and for each one we compute $\Pi_3(R_{\rm ap})$ at the same radii where we evaluate the mass biases ( $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$, see Sect. 4). The various integrals involved in the computations are performed according to the SPH prescription. For the simulated clusters of our sample, the values of $\Pi _3$ we found lie in the range between $\sim$-(6-4) for a strongly asymmetric distribution and $\sim$-(12-10) for a cluster with a relaxed configuration.

In order to perform a statistical analysis of mass determination biases as a function of substructure we extract four sub-samples in the following way. For a given overdensity, we first construct the cumulative distribution of the sample values of $\Pi_{\rm 3}$. Then we identify the synthetically observed clusters in four sub-samples ($\Pi _3$ classes, hereafter): first quartile (below $25\%$), second quartile (between $25\%$ and $50 \%$), third quartile (between $50 \%$ and $75\%$), and forth quartile (above $75\%$). In the following these four sub-samples will be referred to as $q_{\rm 1}$, $q_{\rm 2}$, $q_{\rm 3}$, and $q_{\rm 4}$ classes. By construction the most relaxed clusters belong to the $q_{\rm 1}$ class and the most disturbed ones to the $q_{\rm 4}$ class.

   
3.3 Computation and modeling of surface brightness profiles

The azimuthally averaged surface brightness profiles $S_{{\rm X}}(r)$ are computed from the [0.5-2] keV band surface brightness maps as follows. Each map is binned using a grid with circular geometry and coordinates ( $\tilde {\rho}, \tilde{\phi}$); the origin is set at the position of the X-ray emission peak. The grid points are uniformly spaced, linearly in the angular coordinate, and logarithmically in the radial coordinate. The range of the spatial coordinates is between $2 \times 10^{-4}$ and 1.5in units of r₂₀₀. There are $N_{\tilde {\rho}}=140, N_{\tilde{\phi}}=20$ points in the coordinate intervals. The number and spacing of the radial grid points is chosen such that the surface brightness binning in the radial range used in the fits (see Sect. 4 below) is similar to the one adopted in the analysis of real X-ray data. The value of $S_{{\rm X}}(r)$ at a given projected radius r is given by averaging over the azimuthal values of the annulus. These values are also used to define a surface brightness dispersion which is later used to define a weight when fitting the $S_{{\rm X}}(r)$ profile.

In order to explore the influence of the models adopted in the fitting of surface brightness radial profiles on total mass estimates, we adopt two different models. In the simplest and most widely adopted procedure used in the determination of the total gravitational mass of clusters, the surface brightness profile is modeled by a single $\beta$-model (Cavaliere & Fusco-Femiano 1976), i.e.

$$S_X(r) = S_0 \left( 1+ \frac{r^2}{r^2_{\rm c}} \right)^{-3\beta + 0.5}\cdot$$ (9)

The single $\beta$-model is convenient because it is simple to deproject it and obtain the ICM density profile needed to estimate the total mass from the hydrostatic equilibrium.

As a second parametrization we adopt an extended $\beta$-model in which the profiles are modeled as in Vikhlinin et al. (2006). The radial dependence of $\rho_{\rm g}(r)$ is given by the functional form:

$$\rho_{\rm g}^2(r)=\rho_{{\rm g},0}^2\frac{(r/r_{\rm c})^{-\a... ...ac{1}{(1+r^{\gamma}/r_{\rm s}^{\gamma})^{\varepsilon/\gamma}},$$ (10)

where the additional parameters with respect to the standard $\beta$-model allow a much more accurate modeling of the density slope variations. The component present in in Eq. (10) which takes into account a central surface brightness excess is neglected (i.e. we set $\alpha=0$), since, as discussed below, we are excising the very central regions in the fitting procedure. The model parameters are recovered by fitting the projection of Eq. (10) along the line of sight against the surface brightness profiles in the [0.5-2] keV energy band. Following the suggestions of Vikhlinin et al. (2006) we set in Eq. (10) $\gamma=3$.

Notice that for both parametrizations the total mass derived from the hydrostatic equilibrium equation (Eq. (14) below) does not depend on the central gas density value.

As commonly found in simulated clusters, the gas properties in the central regions of our simulated objects are in disagreement with observations. The pronounced steepening of the gas density towards the center in low temperature systems and the sudden temperature drop in the very inner regions are in fact not observed in real clusters. We therefore exclude the inner region within $r_{\rm low}=0.1 \times r_{\rm 200}$ from the surface brightness fits. We emphasize that the precise value for $r_{\rm low}$ is not relevant, since our goal is to quantify the bias in the mass reconstruction at much larger radii.

   
3.4 Computation and modeling of temperature profiles

From the three-dimensional gas temperature distribution $T^{3{\rm D}}(\vec x )$ measured in the simulated clusters we compute both three- and two-dimensional (projected) radial profiles.

A three-dimensional radial temperature profile can be defined as

$$T_{W}^{3{\rm D}}(r)= \frac {\int_{\Delta V} T^{3{\rm D}}(\vec x ) {\it W} {\rm d}^3x } { \int_{\Delta V} {\it W} {\rm d}^3x},$$ (11)

where W is a weight function and the volume integral is over a spherical shell of thickness $\Delta r$ located at distance r from the cluster center. Common choices for $\it W$ are the gas density (mass-weighted temperatures, ${\it W}=\rho_{\rm g}$) and the X-ray emissivity (emission-weighted temperatures, ${\it W}=\varepsilon_{\rm X}$). In this work we consider mass-weighted temperatures $T_{\rm mw}$ and the weight function ${\it W}=\rho^2_{\rm g} T^{-3/4}$, which defines the corresponding spectroscopic-like temperature $T_{\rm sl}$. Mazzotta et al. (2004) found that in the continuum regime ( $T \gtrsim 2$ keV) this choice of the weight function provides accurate approximations of spectroscopic temperatures obtained from X-ray observations (i.e. derived by fitting an integrated spectrum with a single-temperature model). A more sophisticated method would be based on the temperature determinations from either projected spectra (in the case of projected temperature profiles) or from deprojected ones (in the case of deprojected temperature profiles) computed from the simulation outputs. Given the size of our sample and that this procedure is computationally expensive, we adopt spectroscopic-like temperature profiles.

For the cluster under consideration we evaluate the weighted temperature profiles by replacing the integrals in Eq. (11), according to the SPH scheme, by a summation over gas particles. The value of the radial temperature profile $T_{\it W}^{3{\rm D}}(r)$ corresponding to each shell is defined by averaging over a set of $(\theta,\phi)=40~\times~40$grid points uniformly spaced in angular coordinates. The radial coordinate of the shells is binned as in the computation of the surface brightness. For each object we therefore obtain two types of three-dimensional temperature profiles: a mass-weighted profile and a spectroscopic-like one.

The projected temperature profiles $T_{W}^{2{\rm D}}(r)$ are computed in the same way as the three-dimensional profiles, but by considering a cylindrical geometry instead of a spherical one. Furthermore, only spectroscopic-like profiles are computed, since mass-weighted projected profiles are not meaningful in this case. For each cluster we compute three projected spectroscopic-like temperature profiles by consistently considering the same projections used in the computation of surface brightness profiles.

We compare projected to three-dimensional spectroscopic-like profiles and find that in the radial range [0.1-1] $r_{{\rm 200}}$ the difference is generally small. The agreement between the two types of profiles steadily improves from disturbed to relaxed clusters and the difference is negligible for very relaxed clusters.

In the computation of temperatures, in accordance with the considered energy bands, we neglect in the SPH summation particles with temperatures below 10⁵ K. For consistency, we adopt the same temperature cut-off in the computation of mass-weighted temperatures. Moreover, for the spectral weighting we restrict the summation to those gas particles for which kT > 0.5 keV.

It is crucial to remark that spectroscopic-like temperatures are good approximations of the spectroscopic temperatures as would be estimated from X-ray observations, whereas mass-weighted temperatures are a much more accurate indicator of the total binding mass, since they follow the virial relationship.

As for the brightness profiles, we use two different models to fit temperature profiles and exclude the inner aperture within $r_{\rm low}=0.1~\times~r_{\rm 200}$ from the fits. In the simplest case the gas temperature profile is modeled using the commonly adopted polytropic relation

$$T_{\rm g} \propto \rho_{\rm g}^{\gamma - 1},$$ (12)

with $1 \leq \gamma \leq 5/3$.

In the second model we allow greater parametric freedom and follow the modeling proposed by Vikhlinin et al. (2006). We model the temperature profile using

$$T_{\rm g}(r)=T_{{\rm g},0} \frac{(r/r_{\rm t})^{-a}}{[1+(r/r_{\rm t})^b]^{c/b}},$$ (13)

which describes a power-law declining profile with a transition at $r\sim r_{\rm t}$. In the above model the $t_{\rm cool}(r)$ term, which accounts for the temperature decline in the central region because of radiative cooling (Vikhlinin et al. 2006), is absent. This is motivated by our choice of the radial range over which the profiles are fitted.

When modeling temperature profiles, a proper treatment would require us to fit the model projection along the line of sight against the simulated spectral temperature profiles. Moreover, one has to take into account the presence of different temperature components which can bias the single-temperature fit. Given the large size of our sample, we avoid the involved projection and recovery steps and we fit the model in Eq. (13) directly to the 3D spectroscopic-like temperature profiles. These are in fact expected to be in good agreement with deprojected temperature profiles derived from spatially resolved spectroscopy. However, in order to explore the crucial bias introduced by spectroscopic temperatures, we additionally apply the same procedure to 3D mass-weighted temperature profiles, since they provide the azimuthal average of the true gas temperature distribution.

   
4 Total mass determination from X-ray observations

The total gravitating mass of a galaxy cluster is estimated from X-ray observations assuming that the ICM is in hydrostatic equilibrium and that its temperature and density distributions, as well as the total gravitational potential, are spherically symmetric. These assumptions lead to:

$$M_{{\rm }}^{{\rm est}}(<\!\!r)= - \frac{kT_{\rm g}(r) ~ r}{G ... ...n r} + \frac{{\rm d} \ln T_{\rm g}(r)}{{\rm d} \ln r} \right],$$ (14)

where $M_{{\rm }}^{{\rm est}}(<$r) is the estimated total gravitating mass within a cluster-centric distance r, $T_{\rm g}(r)$ and $\rho_{\rm g}(r)$ are the three-dimensional gas temperature and density profiles at the radius r, respectively, and G and $m_{\rm p}$ are the gravitational constant and proton mass. We adopt the value $\mu = 0.58$ for the mean molecular weight of the gas.

   
4.1 The mass bias

In order to quantify the difference between estimated mass $M^{{\rm est}}(<$r) and the actual mass of the simulated object M(<r), we define the mass bias as:

$$b(r)= \frac{M^{\rm est}(<\!\!r)-M^{}(<\!\!r)}{M^{}(<\!\!r)}\cdot$$ (15)

The true mass M(<r) is computed by summing the masses of all the particles (dark matter, gas, and star particles) within the radius r.

In our observational-like analysis the functions $\rho_{\rm g}(r)$ and $T_{\rm g}(r)$ in Eq. (14) are modeled as described in Sects. 3.3 and 3.4. The various options for the choice of fitting functions and temperature profiles allow us to investigate and disentangle biases of different origin. In particular we explore three different cases:

• Polytropic $\beta$-model: for a given cluster and projection we fit the surface brightness profile and the corresponding projected, spectroscopic-like temperature profile using Eqs. (9) and (12), respectively. The derived mass bias is referred to as $b_{{\rm\beta, ~{\it T}_{\rm sl}^{\rm 2D}}}(r)$.

• Extended $\beta$-model with spectroscopic temperatures: for a given cluster and projection we fit the surface brightness profile with Eq. (10). The three-dimensional, spectroscopic-like temperature profile cluster is fitted using Eq. (13). In this case the mass bias is denoted by $b_{{\rm ext \beta, ~{\it T}_{\rm sl}^{\rm 3D}}}(r)$.

• Extended $\beta$-model with mass-weighted temperatures: same as the previous case, but the three-dimensional, mass-weighted temperature profile is adopted instead of the spectroscopic-like profile. The derived mass bias is denoted by $b_{{\rm ext \beta,~{\it T}_{\rm mw}^{\rm 3D}}}(r)$.

For each of the $\sim$300 observations the derived mass bias is evaluated at three different overdensities, corresponding to $\Delta=2500,~500$ and 200. Cluster-centric distances of $r_{{\rm 2500}}$ and $r_{{\rm 500}}$ are usually probed by observations, with $r_{{\rm 2500}}$ being well inside the largest accessible radius and $r_{{\rm 500}}$ typically being the largest distance where the ICM temperature can be reliably estimated through X-ray observations (e.g. Pratt et al. 2007; Piffaretti et al. 2005; De Grandi & Molendi 2002; Leccardi & Molendi 2008; Vikhlinin et al. 2005; Snowden et al. 2008). We additionally choose $r_{{\rm 200}}$ in order to extend our analysis to radii that will be accessible in the near future (e.g. see Reiprich et al. 2008), to quantify the effects of extrapolation to large distances and, most important, to investigate the mass bias in dynamically different regions.

We adopt three different values for the outer boundary of the region considered in the mass reconstruction by using data (i.e. surface brightness and temperature profiles) within $r_{\rm up}=0.5,0.75 ~ {\rm and} ~ 1.05 \times r_{\rm 200}$ (the inner boundary is $r_{\rm low}=0.1 \times r_{\rm 200}$, see Sects. 3.3 and 3.4). The choice of $r_{\rm up}$ is of paramount importance to assess the dependence of the mass bias upon changes in the gas density slope, which can have a significant impact on mass measurement biases. The three different choices allow us to explore a limitation that in practice is imposed by the field of view of the detector and the background level of a given exposure.

To summarize, for a given radial range ([0.1-0.5] $r_{{\rm 200}}$, [0.1-0.75] $r_{{\rm 200}}$, or [0.1-1.05] $r_{{\rm 200}}$, with the latter denoted for simplicity by [0.1-1] $r_{{\rm 200}}$ hereafter), we fit surface brightness and temperature profiles with the analytic functions specific to the adopted model and evaluate Eqs. (14) and (15) at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$. It is important to notice that by averaging over the whole sample we find: $\langle$ $r_{{\rm 2500}}/r_{{\rm 200}}$ $\rangle_{{\rm sample}} = 0.31$ and $\langle$ $r_{{\rm 500}}/r_{{\rm 200}}$ $\rangle_{{\rm sample}}~=~0.67$. As a consequence $M_{{\rm }}^{{\rm est}}(<$r), and therefore $b_{{\rm }}(r)$, is extrapolated using the best fit functions at $r_{{\rm 500}}$ and $r_{{\rm 200}}$ for the fitting range [0.1-0.5] $r_{{\rm 200}}$ and at $r_{{\rm 200}}$ for the fitting range [0.1-0.75] $r_{{\rm 200}}$. No extrapolation is of course necessary when the fitting range [0.1-1] $r_{{\rm 200}}$ is adopted.

The large size of the simulated mock sample ($\sim$300 observations) allows us to use a statistical approach when exploring the dependence of the mass bias on various quantities. For any given mass determination method, overdensity at which the enclosed mass is estimated, radial fitting range, and $\Pi _3$ class q_i (or whole sample) the mass bias values are modeled using distribution fitting. The data set is always modeled with both normal and Weibull distributions in order to describe symmetric and asymmetric data distributions, respectively. We never find, however, strongly skewed distributions. Even when a Weibull distribution yields a better fit to the data than a normal distribution, its mean and standard deviation are only slightly different from those derived from the normal distribution best-fit. Since these small differences are completely irrelevant for our discussion, we quote only results from normal distribution fitting throughout the paper. From the data set of values $b_{i}(r_{{\rm\Delta}})$ ( $i=\beta, T_{\rm sl}^{\rm 2D}$; ${\rm ext} \beta, T_{\rm sl}^{\rm 3D}$; ${\rm ext} \beta, T_{\rm mw}^{\rm 3D}$), distribution fitting yields the mean value $m(b_{i}(r_{{\rm\Delta}}))$ and standard deviation $\sigma(b_{i}(r_{{\rm\Delta}}))$.

Notice that sample S1 is extracted from a volume 8 times larger than the one for S2 (see Sect. 2). This implies that the cluster abundance measured in our sample is not correct. In order to check whether this issue has an important impact on our results we proceed as follows. We create a new sample made of sample S2 and 1/8 of randomly selected objects from sample S1, for which mean values and standard deviations of the mass bias are computed for the all the cases under consideration. The procedure is repeated many times and the derived values are compared with those computed for the original sample. In all the cases we find that the differences are fully negligible and we therefore report results for the original sample in the reminder of the paper.

The four $\Pi _3$ classes q_i (see Sect. 3.2) are defined such that each class contains the same number of objects ($\sim$70). This guarantees a meaningful comparison between mean values and standard deviations derived from different $\Pi _3$ classes.

   
4.2 Results

Total cluster mass determinations from an X-ray observations, in practice might be affected by more that one source of bias. The number of different cases for which we derive total masses (fitting ranges, best fit functions, temperature definitions, etc.) are of course designed in order to investigate the various sources separately, but also to show the effect of their combination. Our results are reported in Table 2 and conveniently shown in three figures (Figs. 1-3). Notice that for clarity the mean values (data points in the figures) and standard deviations (errorbars in the figures) of the various mass bias distributions are given in percent. The reported values show, in most of the cases, the combined effect of different biases and are very useful to estimate the total bias for a given mass reconstruction method and conditions. Nevertheless, some specific cases and the comparison of mass biases derived under different conditions allow us to quantify the effect of single biases separately.

Figure 1: Summary of the results for the mass reconstruction using the polytropic $\beta$-model and projected, spectroscopic-like temperature profiles. The three panels show quantities derived from modeling of the profiles in three different radial ranges. In each panel we show the mean (points) and standard deviation (errorbars) of the mass bias distribution at $r_{{\rm 200}}$ (stars), $r_{{\rm 500}}$ (diamonds), and $r_{{\rm 2500}}$ (triangles). Open/filled symbols indicate quantities derived with/without extrapolation of the mass profile. The results from distribution fitting are shown for the whole sample (all) and when the four $\Pi _3$ classes are used as sub-samples. The results for different overdensities are shifted horizontally to improve clarity. Notice that the mass bias is shown here in percent.

Open with DEXTER

	Figure 2: Same as in Fig. 1, but for the mass reconstruction method adopting the extended $\beta$-model and 3D spectroscopic-like temperature profiles.
Open with DEXTER

   
4.2.1 The modeling bias

In order to disentangle the bias introduced by inaccurate modeling of temperature and density profiles from other biases, we compare the mean mass bias from polytropic $\beta$ and extended $\beta$-models with spectroscopic temperatures only when no extrapolation of the profiles is needed. Notice that the two models are based on projected and three-dimensional spectroscopic-like temperature profiles, respectively. Nevertheless, the difference between the two types of profiles is small (see Sect. 3.4) and does not introduce any significant additional bias, in particular for the relaxed clusters.

We therefore focus on the filled symbols in Figs. 1 and 2 and the corresponding values in Table 2. As ubiquitously found in simulations, the gas density slope in our simulated clusters considerably steepens with radius. In addition, radial temperature profiles also show a non-trivial radial dependence. Thus, the performance of a given model adopted for the surface brightness and temperature fitting strongly depends on the ability to model these slope changes.

From a visual inspection of the fitted profiles we find that, as expected, the extended $\beta$-model is extremely accurate in modeling surface brightness and temperature profiles for any given radial fitting range, especially for the most relaxed clusters. The accuracy of this model is reflected in the fact that, when no extrapolation is involved, the bias mean values and standard deviations do not depend on the radial fitting range (compare the filled symbols in the different panels in Fig. 2 and the values in in Table 2).

The difference between the mean values derived from extended $\beta$- and polytropic $\beta$-models (compare filled symbols in the bottom panels of Figs. 1 and 2 or the corresponding values listed in Table 2) provides a direct measure of the bias due to the inaccurate modeling provided by the polytropic $\beta$-model. The comparison shows that the polytropic $\beta$-model introduces an additional mass underestimate which is on average $\sim$5, 10, and $15 \%$ at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$, respectively.

   
4.2.2 The extrapolation bias

Table 2: Mean values and standard deviations (errors) of the mass bias derived from the observational-like analysis. These are listed in percent for the three adopted models, the three different fitting ranges, and the three overdensities at which they are evaluated. Results are given for the whole sample and for the four $\Pi _3$ classes.

The bias originating from extrapolating of the mass profiles beyond the radial range probed by observation can be extremely large. While this type of bias can be easily kept under control by simply avoiding any mass estimation/comparison beyond a given outer radius, it may be in practice the cause of many disagreements between different mass estimates (i.e. different X-ray analyses, X-ray and lensing, etc.).

The bias introduced by extrapolation of course affects all the explored models. However, in order to evaluate its net effect, we restrict our discussion to the results obtained from the extended $\beta$-model and mass-weighted temperature profiles, since spectroscopic temperatures introduce an additional bias that is discussed below. The polytropic $\beta$-model, in addition to being inaccurate (see Sect. 4.2.1), is also affected by this latter bias. For a given overdensity, the average extrapolation bias can be therefore derived by comparing open with filled symbols in Fig. 3 and the corresponding values in Table 2.

While the extrapolation bias is of course due to the fact that best fit functions for data in a given radial range do not provide a good description of data outside the radial range, it is interesting to notice that it causes a systematic underestimate of the total mass. For the whole sample, the difference between the mean mass bias obtained adopting the the fitting range [0.1- $1] r_{{\rm 200}}$ (i.e. with no extrapolation) and that obtained from the fitting range [0.1- $0.5] r_{{\rm 200}}$ indicates that extrapolation causes an additional $\sim$$10 \%$ average underestimate at $r_{{\rm 200}}$ (see Table 2). For a very large sample the average bias introduced by extrapolation is thus not extremely large. Furthermore, if only the most relaxed clusters (i.e., the sub-sample $q_{\rm 1}$) are taken into account this bias is fully negligible (see Table 2). These considerations are true for average bias values. It is however extremely important to notice that standard deviations of mass biases computed from extrapolated values are at least twice as big as those obtained from non-extrapolated values, independent of the dynamical state of the cluster (see the errorbars in Fig. 3 and errors in Table 2). This clearly shows that, for individual objects, extrapolation can lead to extremely large mass over/underestimates even for the most dynamically relaxed objects.

   
4.2.3 The dynamical state bias

X-ray observations of unrelaxed clusters are expected to yield bias mass estimates because both assumptions of spherical symmetry and hydrostatic equilibrium are not valid.

Here we quantify the mass bias due to the unrelaxed dynamical state of the cluster by considering: estimates derived from the extended $\beta$-model (because of its accuracy, see 4.2.1), 3D mass-weighted temperature profiles (in order to avoid spectroscopic temperature biases), and mass estimates derived without extrapolation (to avoid contamination from extrapolation biases). Our results are shown with filled symbols in Fig. 3 and since for the case under consideration the results do not depend on the adopted fitting range we focus our discussion on the values shown in the lower panel (i.e. the fitting range [0.1- $1] r_{{\rm 200}}$, see Table 2 for the corresponding values).

While at $r_{{\rm 2500}}$ we do not find any appreciable variation of the mean mass bias for the four $\Pi _3$ classes (see Sect. 3.2 for the definition of these classes), at $r_{{\rm 500}}$ and $r_{{\rm 200}}$ we find a clear trend with cluster dynamical state: the mass is on average more underestimated for disturbed clusters than for the relaxed ones. We find, for example, $m(b_{\rm ext \beta, ~{\it T}_{\rm mw}^{\rm 3D}}(r_{\rm 200}))=-3 \%$ and $-23 \%$ for the $q_{\rm 1}$ class (the most relaxed clusters) and $q_{\rm 4}$ class (the most disturbed clusters), respectively.

Furthermore, we find that at all overdensities the individual mass bias values are more scattered around the mean for disturbed objects than for the relaxed ones, as shown by the standard deviations in Fig. 3 and errors in Table 2. This is illustrated more specifically in Fig. 4 for the case of the fitting range [0.1- $1] r_{{\rm 200}}$ and $\Delta =500$.

These results are particularly relevant because they clearly show that mass underestimates in clusters that are identified as relaxed according to the power ratio method are on average small ($\sim$$2-5 \%$), as opposed to the very large bias that affects total mass determinations in disturbed systems. Finally, we remark that, in agreement with other studies (e.g. Jeltema et al. 2008), we do not find any correlation between $\Pi _3$ and true total mass.

	Figure 3: Same as in Fig. 2, but for results derived using mass-weighted temperature profiles.
Open with DEXTER

   
4.2.4 The spectroscopic temperature bias

The difference between emission-weighted, mass-weighted, and spectroscopic ICM temperatures has been addressed by a number of authors (e.g. Rasia et al. 2005; Mazzotta et al. 2004; Gardini et al. 2004; Vikhlinin 2006; Mathiesen & Evrard 2001; Valdarnini 2006). Special emphasis has been placed on global temperatures because of their impact on the scaling relations such as the luminosity-temperature and mass-temperature relations. Recent studies indicate that $T_{\rm mw}(<$ $r_{\rm\Delta}) < T_{\rm sl}(<$ $r_{\rm\Delta}) < T_{\rm ew}(<$ $r_{\rm\Delta})$, where usually $\Delta =500$ (e.g. Ciotti & Pellegrini 2008; Nagai et al. 2007a). This relation is also valid for our simulated sample. By averaging over the whole sample we find the ratios: $T_{\rm mw}(<$ $r_{\rm 500})\!\!\!:\!\!\!T_{\rm sl}(<$ $r_{\rm 500})\!\!\!:\!\!\!T_{\rm ew}(<$ $r_{\rm 500})=1$: $1.13\!\!\!:\!\!\!1.27$. In addition we computed average temperatures from the same sample but without applying the masking procedure described in Sect. 3.1 and find: $T_{\rm mw}(<$ $r_{\rm 500})\!\!:\!\!T_{\rm sl}(<$ $r_{\rm 500})\!\!:\!\!T_{\rm ew}(<$ $r_{\rm 500})=1\!\!:\!\!1.08\!\!:\!\!1.21$. Considering that mass-weighted global temperatures are unaffected by masking, these results show that both emission-weighted and spectroscopic-like global temperatures are on average biased high by $\sim$$5 \%$ if cold ICM clumps are not masked out appropriately.

In the following we focus on the difference between $T_{\rm mw}$ and $T_{\rm sl}$ since emission weighted temperatures are not relevant in our analysis. As noted by Nagai et al. (2007a) global spectroscopic temperatures within large apertures are higher that mass-weighted ones because they are dominated by the inner, hotter cluster region. It is important to notice that accurate total mass determination methods such as the one adopted in this paper rely on spatially resolved temperature determinations and not on global temperatures. Therefore it is crucial to compare $T_{\rm mw}$ and $T_{\rm sl}$ radial profiles and how any difference between the two propagates into total mass estimates.

We compared spectroscopic like and mass-weighted temperature profiles and find, in agreement with Ameglio et al. (2007), that $T_{\rm mw}(r) > T_{\rm sl}(r)$. This inequality, and the fact that the opposite is true when considering global temperatures, can be understood as follows. Global temperatures are computed by averaging within large cluster-centric distances whereas profiles are computed by averaging within radial shells. The density and temperature distributions in the two cases are dramatically different. In the former case the distribution is dominated by the central, hot and dense gas, while in a given radial shell the ICM is obviously distributed in a much narrower range of densities and temperatures. The ICM in our simulations is locally multiphase and made of a dominant component and mostly lower temperature particles which do not substantially contribute to mass-weighted temperature, but bias spectroscopic temperatures low. A simple illustration of the spectroscopic temperature bias for a two-phase ICM can be found in Jia et al. (2008).

In the computation of both $T_{\rm mw}(r)$ and $T_{\rm sl}(r)$ all cold clumps are excluded. This implies that local spectroscopic-like temperatures can underestimate mass-weighted ones also because of cold particles that are not part of clumps. Interestingly, we find that our mass reconstruction method depends very weakly on cold clump removal. Total masses are in fact derived by fitting brightness and temperature profiles with smooth functions of radius, which are in almost all the cases insensitive to the irregularities (spikes or bumps in the brightness profiles and dips or depressions in the temperature profiles) arising from cold clumps when no masking is applied.

In relaxed clusters $T_{\rm mw}$ and $T_{\rm sl}$ agree very well at small and intermediate radii. The agreement becomes gradually worse with increasingly disturbed dynamical state (as measured by the power-ratio method). This reflects the more complex thermal structure of disturbed clusters, where the fraction of cold gas particles stripped from cold sub-structures is relatively high. Furthermore, we find that the relative difference between $T_{\rm mw}(r)$ and $T_{\rm sl}(r)$ increases with radius. This is expected because of the weighting scheme used in the computation of $T_{\rm sl}$ and the presence of cold infalling gas in the cluster outskirts. A similar finding is reported in Rasia et al. (2006). Finally, the relative difference between $T_{\rm mw}$ and $T_{\rm sl}$ is higher in more massive clusters than in low mass systems because in massive objects the spread in temperature is larger (e.g. Ameglio et al. 2007; Valdarnini 2006).

	Figure 4: Mass bias as a function of $\Pi _3$ at $r_{{\rm 500}}$ with estimated masses derived adopting the extended $\beta$-model and 3D mass-weighted temperature profiles. The solid line shows the best fit function $b_{{\rm ext} \beta , T_{\rm mw}^{\rm 3D}}(r_{{\rm 500}})=-0.256$- $0.022 \times \Pi _{\rm 3}(r_{{\rm 500}})$.
Open with DEXTER

The underestimate of the true ICM temperature profile introduced by spectroscopic measurement leads to total mass underestimates because hydrostatic masses (Eq. (14)) depend on both temperature and its logarithmic derivative. In order to disentangle this effect from the other observational biases discussed in the previous subsections, we consider results obtained from the two extended $\beta$-models and compare results obtained from the one adopting spectroscopic-like profiles with those derived using mass-weighted temperature profiles. The net mass bias introduced by spectroscopic temperatures can be readily estimated by computing the difference between the mean bias values listed in Table 2. However, a clearer determination can be performed by defining the spectroscopic bias:

$$b_{\rm spec}(r)= \frac{M^{\rm est}_{{\rm ext} \beta, T_{\rm s... ...(<\!\!r)}{M^{est}_{ext \beta, T_{\rm mw}^{3{\rm D}}}(<\!\!r)},$$ (16)

where $M^{\rm est}_{{\rm ext} \beta, T_{sl}^{3D}}(<$r) and $M^{\rm est}_{{\rm ext} \beta, T_{\rm mw}^{3{\rm D}}}(<$r) are total masses estimated through the extended $\beta$-model with spectroscopic-like and mass-weighted temperatures profiles, respectively. Both are derived adopting the the fitting range [0.1- $1] r_{{\rm 200}}$ in order to avoid any contamination introduced by the extrapolation.

In Table 3 we list the spectroscopic bias mean values and standard deviations computed for the whole sample and for the four $\Pi _3$ classes separately. These values are also shown in Fig. 5 for clarity. The spectroscopic bias is systematic and of the order of $-10 \%$ if the whole sample is considered. However, if the same comparison is performed for the four $\Pi _3$ classes separately one finds that it is larger (in absolute value) in disturbed clusters than in relaxed ones. The mean spectroscopic bias is also greater at larger cluster-centric distances than at small ones. These results are of course directly related to the relative difference between the $T_{\rm mw}$ and $T_{\rm sl}$ profiles discussed above.

Thus, a crucial point of out findings is that also after the removal of all the resolved ICM cold clumps in the simulations, spectroscopic temperature profiles underestimate the true ones, and therefore the total mass, because of a diffuse distribution of cold gas particles.

Table 3: Mean values and standard deviations (errors) of the spectroscopic bias $b_{\rm spec}(r_{\rm\Delta })$ (Eq. (16)). Values are in percent and listed for three different true mass ranges used in their derivation: all masses (the whole sample), lowM-sample (clusters with M(< $r_{\rm 200}) < 1.8 \times 10^{14} ~ h^{-1}~M_{\odot}$), and highM-sample (clusters with M(< $r_{\rm 200}) > 1.8 \times 10^{14} ~ h^{-1}~M_{\odot}$). In each case results are given, at three radii, for the whole sample and for the four $\Pi _3$ classes.

	Figure 5: Mean (points) and standard deviation (errorbars) of the spectroscopic bias $b_{\rm spec}(r)$ (Eq. (16)) at $r_{{\rm 200}}$ (stars), $r_{{\rm 500}}$ (diamonds), and $r_{{\rm 2500}}$ (triangles). The results are shown for the whole sample (all) and when the four $\Pi _3$ classes are used as sub-samples. The results for different overdensities are shifted horizontally to improve clarity.
Open with DEXTER

   
4.2.5 Mass dependence

Our results concerning the influence of the various observational biases allow us to focus on some relevant cases when exploring the dependence of the mass bias on the total, true mass. We can avoid any contamination from biases such as the modeling bias or the extrapolation bias by considering results provided by extended $\beta$-models and without extrapolation of the mass profiles.

As expected, no mass dependence is found for mass biases derived when using the extended $\beta$-model, as shown in Fig. 6 in the case of fitting range [0.1- $1] r_{{\rm 200}}$ and $\Delta =500$. The independence of bias on total mass is found for all fitting ranges and all overdensities. This is, however, only found when mass-weighted temperature profiles are used. Spectroscopic-like temperatures introduce a mass dependence: the average bias is more negative for massive systems, as shown in Fig. 7 for a specific case. This mass dependence arises from the weighting scheme used in the computation of spectroscopic-like temperatures, which is sensitive to the cooler part of the gas distribution, and the scale dependency introduced by cooling on the amount of this cold component versus the cluster mass (see Sect. 4.2.4).

In order to investigate the dependence of the spectroscopic bias on true mass and at the same time on dynamical state, we compute mean values and standard deviations by considering objects in a given (true) mass interval and $\Pi _3$ class. In order to obtain statistically significant values (i.e. a sufficient number of objects in each class), we must adopt a coarse mass binning. After experimenting with different types of mass binning we find that the discussion can be greatly simplified by reporting results for only two sub-samples, which are obtained by ordering the sample by true cluster mass and halving it. In the following these two sub-samples are referred to as the lowM- and highM-sample.

For each of these sub-samples we construct the four $\Pi _3$ classes separately and compute the mean and standard deviation of the spectroscopic bias as done above for the whole sample. The results are reported in Table 3, together with those obtained for the whole sample, and shown in Fig. 8. At $r_{{\rm 2500}}$ we do not find any significant dependence on mass and therefore the average biases reported in Sect. 4.2.4 are valid at any mass. On the other hand, we find a strong mass dependence at $r_{{\rm 500}}$ and $r_{{\rm 200}}$ and that its strength increases with cluster-centric distance. For the lowM-sample the average spectroscopic bias does not significantly vary with overdensity, as opposed to the behavior found for the highM-sample. Because of the findings reported in Sect. 4.2.3, relaxed clusters ($q_{\rm 1}$ and $q_{\rm 2}$ sub-samples) are of particular relevance. The ones in the lowM-sample are very mildly affected by the spectroscopic bias at all cluster-centric distances, whereas for those in the highM-sample the mean spectroscopic bias is $\sim$$-10 \%$ at $r_{{\rm 500}}$ and $\sim$$-15 \%$ at $r_{{\rm 200}}$.

	Figure 6: Mass bias from the extended $\beta$-modeling versus total true mass at $\Delta =500$ for the radial fitting range [0.1- $1] r_{{\rm 200}}$. The estimated total mass is derived using mass-weighted temperature profiles. The solid line shows the best fit function $b_{{\rm ext} \beta, T_{\rm mw}^{3D}}(r_{{\rm 500}})=0.78-0.06 ~ {\rm log} [M(<$ $r_{{\rm 500}})/ h^{-1} ~ M_{\odot}]$.
Open with DEXTER

	Figure 7: Same as Fig. 6, but with the estimated total mass derived using spectroscopic-like temperature profiles. The solid line shows the best fit function $b_{{\rm ext} \beta, T_{\rm sl}^{3{\rm D}}}(r_{{\rm 500}})=2.79-0.21 ~ {\rm log} [M(<$ $r_{{\rm 500}})/ h^{-1} ~ M_{\odot}]$.
Open with DEXTER

	Figure 8: Same as Fig. 5, but for the lowM- (filled symbols) and highM-samples (open symbols), separately.
Open with DEXTER

   
4.2.6 The ideal case

After having examined the various biases affecting the total mass determination, we discuss here the ideal conditions under which the most unbiased mass determination is possible. From our results it is clear that in order to obtain unbiased mass estimates one must select very relaxed clusters (the $q_{\rm 1}$ sub-sample), adopt a model that can accurately fit the ICM density and temperature profiles (the extended $\beta$-model, here), and avoid extrapolation (notice, however, that the extrapolation bias is on average negligible for relaxed clusters. See Sect. 4.2.2). In addition let us assume that the true thermal structure of the ICM (i.e. $T_{\rm mw}(r)$) can be measured. In this case the bias does not depend on the true mass of the system and, at all radii, masses are on average very well reconstructed (see Table 2 and Fig. 3). Notice that there is a residual level of underestimation which might be originating from the still imperfect gas density and temperature modeling or from departure of the gas from hydrostatic equilibrium.

It is of paramount importance, however, to take into account that the true ICM temperature is not accessible through X-ray observations and that we must deal with spectroscopic temperature measurements. As shown in Sect. 4.2.4, the most relaxed clusters are also the less affected by the spectroscopic bias (see Table 3 and Fig. 5). As we have shown in Sect. 4.2.5, the spectroscopic bias depends on cluster mass. It is negligible in low mass clusters, but important for massive ones especially at large radii.

   
5 Total mass determination from the three-dimensional gas distribution

As shown in Sect. 4, mass determinations from X-ray observations suffer from various biases. Of course, in addition to these the assumption of hydrostatic equilibrium (i.e. Eq. (14)) also plays a crucial role in the total mass reconstruction. Moreover, the derivation of Eq.  (14) implicitly assumes spherical symmetry. Whereas the validity of this assumption can be fairly well tested though an X-ray morphological analysis, the robustness of the assumption of hydrostatic equilibrium cannot be directly investigated observationally.

Here, in order to investigate the validity of both assumptions we avoid any bias due to the X-ray reconstruction method and make use of the full three-dimensional gas distribution in the simulation. For consistency with the analysis presented in Sect. 4 we consider only objects with $T_{\rm sl}(<$ $r_{\rm 200}) \geq 2$ keV.

   
5.1 Mass estimators

We introduce the hydrostatic equilibrium mass estimate, $M_{\rm HE}$ by assuming spherical symmetry and rewriting Eq. (14) as

$$M_{\rm HE}(<\!\!r)=-\frac{r^2}{G \rho_{\rm g}} \langle \nabla_r P_{\rm g}(r)\rangle,$$ (17)

where $P_{\rm g}$ is the gas pressure. As a means to avoid any systematic effect generated by estimating gradients from averaged profiles and to achieve maximum accuracy, we evaluate the term $\langle \nabla_r P_{\rm g}(r) \rangle$ directly from the gas particle distribution. To this end, we introduce a spherical shell at the test radius r and a corresponding set of $40 \times 40$ grid points with angular coordinates uniformly spaced in ${\rm cos}\theta, \phi$. At the grid point $\vec x_{\rm gr}$ the pressure gradient is expressed as

$$\vec \nabla P_{\rm g}(\vec x_{\rm gr})=(\gamma-1)\sum m_i u_... ...abla W(\vert\vec x_i -\vec x_{\rm gr}\vert,h_{\vec {\rm g}}),$$ (18)

where u_i is the specific particle internal energy and $\gamma=5/3$. The radial component $\nabla_r$ is then extracted by transforming the gradient vector according to the grid coordinates and finally $\langle \nabla_r P_{\rm g}(r) \rangle$ is obtained by averaging (18) over the set of grid points. We checked the robustness of the estimated average gradients $\langle \nabla_r P_{\rm g}(r) \rangle$ with a two-fold increase in the number of grid points (i.e. $80~\times~80$) and found the results to be unaffected by this choice.

This procedure allows us to define a mass estimate $M_{\rm HE}$ which is influenced by the cluster dynamical state, but in contrast to the mass determinations discussed in Sect. 4 it is not affected by any bias originating from the gas density and temperature profile reconstruction and modeling.

Similarly to Eq. (15) we define the mass bias:

$$b_{\rm HE}(r)= \frac{M_{\rm HE}(<\!\!r)-M^{}(<\!\!r)}{M^{}(<\!\!r)}\cdot$$ (19)

In addition, we also consider the hydrodynamical equilibrium mass estimate (Kay et al. 2004; Rasia et al. 2004):

 

$$M_{\rm DE}(<\!\!r) =M_{\rm HE}(<\!\!r)+M_{\sigma}(<\!\!r)$$

$$=M_{\rm HE}(<\!\!r) - \frac{\sigma_r^2 r}{G}\left( \frac {{\rm d}... ...{\rm d}~ {\rm ln} ~\sigma_r^2}{{\rm d}~ {\rm ln}~ r} + 2 \beta_{\nu}(r)\right),$$ (20)

where $\sigma_{\rm r}$ and $\sigma_{\rm t}$ are the radial and tangential gas velocity dispersion, respectively, and $\beta_{\nu}(r)=1 -{\sigma_{\rm t}^2}/{2\sigma_{\rm r}^2}$ is the gas velocity anisotropy parameter. The estimate of the mass term $M_{\sigma}(<$r) is subject to uncertainties because of the irregularities that are present in the gas velocity dispersion profiles of the simulated clusters. In order to properly estimate the radial derivatives in Eq. (20) we have therefore regularized the measured profiles by applying to them a Savitzky-Golay filter (e.g. Press et al. 1992, Sect. 14.2). This smoothing procedure is effective in removing the small-scale noise and yields mass estimates that are much more regular against the radial dependency. As for $M_{\rm HE}$ we define the mass bias:

$$b_{\rm DE}(r)= \frac{M_{\rm DE}(<\!\!r)-M^{}(<\!\!r)}{M^{}(<\!\!r)}\cdot$$ (21)

Once the mass profiles $M_{\rm HE}(<$r) and $M_{\rm DE}(<$r) for the $\sim$100 objects in our sample are computed, we evaluate the mass biases $b_{\rm HE}(r)$ and $b_{\rm DE}(r)$ at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$. Finally, mean values and standard deviations of the resulting mass bias distributions are computed.

Notice that Eq. (20) is derived from the Boltzmann equation (e.g. Binney & Tremaine 1987, Chap. 4) by assuming steady-state hydrodynamic equilibrium and taking into account the terms from the gas isotropic pressure and anisotropic velocity dispersion. If the latter are neglected (i.e. it is assumed that the velocity dispersion of the gas is much smaller than its temperature), then Eq. (20) yields Eq. (17). In the derivation it is assumed, of course, that the system is spherically symmetric.

The steady-state assumption implies no net radial streaming motions (no mean laminar flow, i.e. $\langle v_r \rangle=0$, which implies $\sigma_r^2= \langle v^2_r\rangle$) and that all time derivatives can be set to zero. Therefore, the applicability of Eq. (20) depends on the degree of spherical symmetry of the system and on the condition $\langle v_r \rangle=0$, while that of Eq. (17) additionally depends on the condition $M_{\sigma}(<$ $r)/M_{\Delta}(<$$r) \ll1$. Hence, the hydrodynamical equilibrium equation naturally provides a more accurate mass estimator than the hydrostatic equilibrium.

The accuracy of spherical symmetry can be assessed by measuring the mean tangential acceleration at the surface of the spherical shell of radius r. To this end, in correspondence with each of the grid points at which we evaluate the pressure gradient in Eq. (18), we also evaluate the gravitation acceleration $\vec a_{g}$ of test particle. We then define a mean torque parameter $\tau _q$ by averaging the acceleration components $a_{\theta }^2+a_{\phi}^2$ and a²_r over the set of grid points:

$$\tau_q^2= \frac{1}{2} \frac{\sum_g a_{\theta}^2+a_{\phi}^2}{\sum_g a^2_r}\cdot$$ (22)

The constraint $\tau_q \ll1$ can be used as a condition to validate the spherical symmetry approximation. Accordingly, the quadrupole moment terms can be neglected in the radial component of the total potential gradient and to the lowest order we have:

$$\vec \nabla_{r} \Phi(\vec r ) \simeq \frac{G M(<\!\!r)}{r^2} ~~~~~~\tau_q \ll1.$$ (23)

It is interesting to notice that the degree of spherical symmetry measured directly from the full 3D information provided by the simulations is highly correlated with that measured from X-ray brightness maps. As shown in Fig. 9 for the specific case of $\Delta=2500$ the torque parameter $\tau _q$ and $\langle \Pi_3 \rangle_{{\rm planes}}$ are correlated. Here $\langle \Pi_3 \rangle_{{\rm planes}}$ is the mean of the three $\Pi _3$ values computed from three orthogonal surface brightness maps (see Sect. 3.2). From a Spearman rank correlation coefficient analysis we find that the correlation is significant at more than the $95 \%$ level at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$. Clusters with small values of the torque parameter $\tau _q$, i.e. whose potential has a high degree of sphericity, are therefore identified as spherical systems when selected by means of the power ratio method.

Figure 9: Mean torque parameter $\tau _q$ (which measures the sphericity of the gravitational potential) as a function of $\langle \Pi_3 \rangle_{{\rm planes}}$ (which measures the sphericity of the gas distribution) at $r_{{\rm 2500}}$. For each cluster the value $\langle \Pi_3 \rangle_{{\rm planes}}$ is the mean of the three $\Pi _3$ values computed from three orthogonal surface brightness maps.

Open with DEXTER

The level of radial streaming motions at the radius r is directly quantified through the radial Mach number $\langle v_r \rangle/c_{\rm s}$, where the gas sound speed $c_{\rm s}$, is computed from $c_{\rm s}^2(r)=5 k T_{\rm g}(r)/(3 \mu m_{\rm p})$.

It follows that the applicability of Eqs. (17) and (20) to measure the cluster gravitating mass within the test radius r is then limited to those systems for which at the radius r the conditions $\tau_q \ll1$ and | $\langle v_r \rangle/c_{\rm s}\vert \ll1$ are satisfied.

   
5.2 Results

The dependence of the mass biases $b_{\rm HE}$ and $b_{\rm DE}$ on the mean torque parameter $\tau _q$ and the radial Mach number $\langle v_r \rangle/c_{\rm s}$ are shown in Figs. 10 and 11, respectively. The figures show that, for a given overdensity, we measure large mass biases (both positive and negative) for clusters with large $\tau _q$ and $\vert\langle v_r \rangle/c_{\rm s}\vert$. As we move from large to small values of the parameters (i.e. to an increasing degree of sphericity of the cluster or to small radial streaming velocities) the scatter of the measured biases substantially decreases. These results demonstrate that, in addition to sphericity, the radial Mach number $\langle v_r \rangle/c_{\rm s}$ is a key parameter in the study of the reliability of the hydrostatic and hydrodynamical mass estimators.

Figure 10: Mass biases $b_{\rm HE}$ (open stars) and $b_{\rm DE}$ (filled triangles) as a function of mean torque parameter $\tau _q$ at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$. In each panel we show the mean values derived from the whole sample for $b_{\rm HE}$ (dashed line) and $b_{\rm DE}$ (dot-dashed line), and the mean values (plotted only up to 0.1) derived from the sub-sample with $\tau _q < 0.15$ and $\vert\langle v_r \rangle /c_{\rm s}\vert < 0.1$ for $b_{\rm HE}$ (short-dashed line) and $b_{\rm DE}$ (solid line). In each panel are also reported these mean values along with the standard deviations. Notice that in order to focus on the bulk of the data we have omitted the results for a few clusters with very large bias and large $\tau _q$.

Open with DEXTER

	Figure 11: Mass biases as a function of radial Mach number $\langle v_r \rangle/c_{\rm s}$. Symbols, lines, and values are as explained in Fig. 10. Also here few clusters with very large bias and large $\vert\langle v_r \rangle/c_{\rm s}\vert$ are omitted from the plot to improve clarity.
Open with DEXTER

It is important to notice that the parameters $\tau _q$ and $\langle v_r \rangle/c_{\rm s}$ are not uncorrelated. We find in fact that objects with large $\tau _q$ also tend to have large $\vert\langle v_r \rangle/c_{\rm s}\vert$, as shown in Fig. 12 for $\Delta=2500$. The correlation between these parameters is found to be significant at more than the $95 \%$ level at all the three considered overdensities.

	Figure 12: Mean torque parameter $\tau _q$ as a function of Mach number $\vert\langle v_r \rangle/c_{\rm s}\vert$ at $r_{{\rm 2500}}$.
Open with DEXTER

In addition to the considerable decrease of the scatter with decreasing values of $\tau _q$ and $\vert\langle v_r \rangle/c_{\rm s}\vert$, a visual inspection of Figs. 10 and 11 indicates that on average the mass estimators $M_{\rm HE}$ is biased low and that masses computed through the more accurate estimator $M_{\rm DE}$ are on average well recovered. This is quantitatively shown by the mean values derived from the whole sample, which are listed in Table 4.

It is crucial that the average biases derived from the hydrostatic mass estimator $M_{\rm HE}$ are in agreement with those obtained from the X-ray mass reconstruction procedure presented in Sect. 4. A comparison between the values in Table 4 and those listed in Table 2 for the extended $\beta$-model with mass-weighted temperature profiles (since these follow the virial relationship more accurately that spectroscopic ones) shows that in the two cases the true mass is similarly underestimated. Furthermore, we find that, as for the mass bias derived using mass-weighted temperatures in Sect. 4, both $b_{\rm HE}$ and $b_{\rm DE}$ do not depend on cluster true mass.

Most important, the mean values of the $b_{\rm HE}$ and $b_{\rm DE}$ distributions show that if in addition to the thermal pressure also the ICM non-thermal pressure component (i.e. the gas anisotropic velocity dispersion) is taken into account, the total mass of the system is on average better reconstructed (see also Kay et al. 2004; Rasia et al. 2004).

Table 4: Mean values and standard deviations (errors) of the mass biases $b_{\rm HE}$ and $b_{\rm DE}$ derived from the whole sample and for the sub-sample with $\tau _q < 0.15$ and $\vert\langle v_r \rangle /c_{\rm s}\vert < 0.1$. Values are given in percent at three different overdensities.

In oder to investigate the accuracy of the mass estimators $M_{\rm HE}$ and $M_{\rm DE}$ more precisely we proceed as follows. At a given overdensity, we create sub-samples of clusters with $\tau _q$ and | $\langle v_r \rangle$ $/c_{\rm s}\vert$ smaller that some given threshold values and in each case we compute mean values and standard deviations of $b_{\rm HE}$ and $b_{\rm DE}$. It is important to find a compromise between small threshold values for $\tau _q$ and $\vert\langle v_r \rangle/c_{\rm s}\vert$ (in order to guarantee as much as possible the applicability of the mass estimators) and a sub-sample with enough objects so as to allow a meaningful distribution fitting. We experimented with various cuts and find that the optimal choice is to select clusters with $\tau_q \langle~0.15$ and $\vert\langle v_r \rangle /c_{\rm s}\vert < 0.1$, which yields sub-samples of at least $\sim$40 objects at each overdensity. Results derived from this sub-sample are listed in Table 4. We remark that the results are not sensitive to exact choice of the threshold values. As expected the standard deviations are much smaller than those derived for the whole sample. The standard deviations for the sample with $\tau_q > 0.15$ and $\vert\langle v_r \rangle /c_{\rm s}\vert \rangle 0.1$ are $\sim$2-3 times larger that those obtained for the sample with $\tau_q \langle~0.15$ and | $\langle v_r \rangle$ $/c_{\rm s}\vert < 0.1$.

The analysis of simulated objects for which the applicability criteria of Eqs. (17) and (20) are most valid allows us, therefore, to draw more robust conclusions on the reliability of the mass estimators $M_{\rm HE}$ and $M_{\rm DE}$. Hydrostatic masses are biased low by $\sim$$5 \%$ at $r_{{\rm 2500}}$ and by $\sim$$10 \%$ at $r_{{\rm 500}}$ and $r_{{\rm 200}}$, whereas the hydrodynamical equilibrium mass estimator is unbiased at any overdensity.

   
6 Discussion and conclusions

This work is aimed at elucidating the reliability of X-ray total mass estimates in clusters of galaxies using N-body/SPH simulation of a large sample of clusters. The physical modeling of the gas includes radiative cooling, star formation, supernova heating, and metal enrichment. The large number of simulated clusters enables us to derive very robust conclusions through a statistical analysis of the sample. The total mass is recovered adopting an observational-like approach and is compared with the true mass in the simulations. Surface brightness and temperature profiles that we generate from the simulations are used to estimate the cluster mass at different overdensities ( $r_{{\rm 2500}}, r_{{\rm 500}}$, and $r_{{\rm 200}}$) by means of the hydrostatic equilibrium equation. We explore various models and conditions under which the mass is reconstructed in order to disentangle different mass biases. In addition, a power ratio analysis of the surface brightness maps allows us to assess the dependence of the mass bias on cluster dynamical state. Moreover, we study the reliability of hydrostatic and hydrodynamical equilibrium mass estimates using the full three-dimensional gas distribution in the simulation.

In the following we list and discuss our main findings.

1.

Our analysis shows that it is very important to use analytical models with a large amount of parametric freedom when modeling the shape of the ICM temperature and surface brightness radial profiles. A model with a low degree of sophistication such as the polytropic $\beta$-model can introduce a very large modeling bias. Compared to the more sophisticated extended $\beta$-model, which is found to be extremely accurate in following the slope changes of the gas profiles, it additionally leads to average mass underestimates of the order of $\sim$5, 10, and $15 \%$ at $r_{{\rm 2500}}$, $r_{{\rm 500}}$, and $r_{{\rm 200}}$, respectively.

2.

The bias originating from extrapolating of the mass profiles beyond the radial range probed by observation can be extremely large. We find that the underestimate from extrapolation alone is of the order of $\sim$$10 \%$ at $r_{{\rm 200}}$, and lower at smaller cluster-centric distances, when considering an average over the whole sample. However, for individual objects, the extrapolation bias can be as large as $\sim$$50 \%$.

3.

The unrelaxed dynamical state of a cluster can also lead to mass underestimates. The total mass is on average biased lower for disturbed clusters than for the relaxed ones. Furthermore, we find that the bias values are much more scattered around the mean for disturbed objects than for the relaxed ones. If mass-weighted temperature profiles are adopted, the mean mass bias is at most $\sim$$-5 \%$ at all cluster-centric distances for relaxed clusters, but can be as large as $\sim$$-20 \%$ for the most disturbed ones.

Our results are in excellent agreement with the findings of Jeltema et al. (2008), as shown by the comparison of Fig. 4 in this work and Fig. 9 (top panel) in their paper, where the correlation between the same quantities are shown. This agreement is of particular importance, considering that the analysis by Jeltema et al. (2008) is based on simulations performed with the adaptive mesh refinement code Enzo. Our findings are also in good agreement with the results of Kay et al. (2007).

4.

Mass estimates derived using spectroscopic temperatures are lower that those derived from mass-weighted temperatures (i.e. the true gas temperatures). We find that the spectroscopic temperature bias, i.e. the bias originating from spectroscopic temperatures alone, is of the order of $-10 \%$ for the whole numerical sample. A similar value is found by Kay et al. (2007).

Mass underestimates derived adopting spectroscopic-like temperature profiles (7 and $17 \%$ on average at $\Delta=2500$ and 500 for the whole sample) are in good agreement with the values found by Nagai et al. (2007a) (12 and $16 \%$ are the corresponding values), whose results are based on mock X-ray observations derived from simulations performed with an Eulerian code. We also find very good agreement by performing the same comparison for relaxed and unrelaxed clusters separately. We notice that total mass biases derived from the simulations by Nagai et al. (2007a) and adopting mass-weighted profiles are on average $-7 \%$ at $r_{{\rm 500}}$ for relaxed clusters (Lau et al. 2007). The corresponding value that we find from our Lagrangian simulations is $-5 \%$. Considering the different nature of the numerical codes and the different implementation of the various physical processes, the agreement is extremely relevant.

Our results are in tension with the findings of Rasia et al. (2006), who find, on average, much stronger biases. We notice, however, that a fair comparison is not possible since the work by Rasia et al. (2006) focused only on 5 clusters.

Our analysis also shows that the spectroscopic bias depends on dynamical state. We find that the spectroscopic bias is rather small (-2,-5, and $-7 \%$ at $\Delta= 2500, 500$ and 200) for the most relaxed clusters and of the order of $-12 \%$ for the most disturbed ones.

5.

While the mass bias derived from mass-weighted temperature profiles does not depend on true cluster mass, the one derived from spectroscopic-like temperature exhibits a strong mass dependence. For clusters with M(< $r_{\rm 200}) < 1.8 \times 10^{14} ~ h^{-1}~M_{\odot}$ we find that, at all radii, the spectroscopic bias contributes on average less than $7 \%$ of the total bias, and less than $4 \%$ if one considers only the most relaxed clusters. On the other hand, for clusters with M(< $r_{\rm 200}) > 1.8 \times 10^{14} ~ h^{-1}~M_{\odot}$ we find that the spectroscopic bias is larger: e.g. $-7 \%$ and $-12 \%$ for the relaxed clusters at $r_{{\rm 500}}$ and $r_{{\rm 200}}$, respectively.

6.

Even in the ideal case, i.e. when the mass of relaxed clusters is estimated without extrapolation and adopting a very accurate model for the ICM density and temperature profiles, total masses are affected by the spectroscopic temperatures. In this case the mean total mass bias is $\sim$(-4,-3), (-12,-9), and $(-15,-3) \%$ at $r_{{\rm 2500}}, r_{{\rm 500}}$, and $r_{{\rm 200}}$ for the (most, least) massive clusters in our simulated sample.

7.

Even if we assume that the true (mass-weighted) ICM temperature can be reconstructed, masses are biased low (by -2, -9, and $-12 \%$ at $\Delta=2500$, 500, and 200, respectively, when taking the average over the whole sample).

8.

The latter finding prompted us to investigate the possible violation of the hydrostatic equilibrium by using the full, three-dimensional information provided by the simulations. In agreement with the results from our observational-like analysis, we find that the hydrostatic equilibrium assumption yields masses underestimated by $\sim$10-$15 \%$. This implies that the origin of the bias is not the X-ray reconstruction method. The same level of mass bias has been found by Rasia et al. (2004) and Burns et al. (2008).

9.

In order to elucidate the origin of this bias we compute masses using both hydrostatic and dynamical equilibrium mass estimators. The dynamical equilibrium takes into account both thermal and non-thermal pressure of the ICM. We find that the masses estimated through dynamical equilibrium are on average well estimated. The mean bias profiles of the ``feedback'' clusters simulated by Kay et al. (2004) are in very good agreement with our average values at all overdensities.

10.

We explore the conditions of applicability of both estimators, i.e. spherical symmetry (torque parameter $\tau_q \ll1$) and small radial streaming velocities (radial Mach number $\vert\langle v_r \rangle/c_{\rm s}\vert \ll1$). We find that the biases $b_{\rm HE}$ and $b_{\rm DE}$ do not depend on cluster true mass and that their scatter decreases with decreasing $\tau _q$ and | $\langle v_r \rangle/c_{\rm s}\vert$.

11.

For clusters with small values of $\tau _q < 0.15$ and | $\langle v_r \rangle/c_{\rm s}\vert < 0.1$ we find that the average mass underestimate found from the hydrostatic equilibrium estimator ($5 \%$ at $r_{{\rm 2500}}$ and $10 \%$ at $r_{{\rm 500}}$ and $r_{{\rm 200}}$) is extremely well corrected by adopting the dynamical equilibrium estimator.

Our results of course depend on the physics included in the simulation. Our modeling of the gas physical processes is incomplete and does not include, for example, energy feedback from active galactic nuclei (AGN). The inclusion of more realistic physics is particularly relevant in the cluster cores, since it is supposed to provide additional heating to the gas which could help to solve the cooling flow problem as well as the overcooling of baryons actually present in current simulations (e.g. Voit 2005).

The effect on the thermal status of the ICM of energy feedback from the AGN are however not expected to modify in a significant way the ICM properties outside the cluster cores. The validity of this assumption is justified from the radial behavior of the measured temperature profiles, which are in good agreement at $r \gtrsim 0.1 \times r_{\rm 200}$ with the present simulations (Valdarnini 2006) and the ones of other authors (e.g. Nagai et al. 2007b).

Importantly, we find that, in addition to cold gas clumps, also the diffuse cold gas component which is left after the removal of all resolved clumps substantially biases spectroscopic temperatures low. The good agreement of our temperature profiles with those found by other authors using different codes suggests that the amount of the cool gas component present in our simulations is not affected by numerical issues.

There is another issue connected to the gas physical modeling in the simulations that could be relevant for our analysis. In the runs performed here the artificial viscosity is treated according to the standard SPH formulation (Monaghan 2005), which is a numerical scheme comparatively viscous (Morris & Monaghan 1997) and inadequate to follow the development of fluid turbulence. It has been shown that the level of kinetic energy in random gas motion could be as high as $\sim$$30 \%$ of the gas thermal energy (Vazza et al. 2006; Dolag et al. 2005), when the numerical viscosity scheme is generalized as in Morris & Monaghan (1997) to properly treat fluid turbulence.

It is worth noticing that very high levels of random motion are inconsistent with the upper limits recently derived from the comparison of X-ray and weak lensing mass estimates (Zhang et al. 2008; Mahdavi et al. 2008). Direct measurements of the various sources of non-thermal pressure (gas motions, cosmic rays, magnetic fields, etc.) are extremely difficult, and the comparison between X-ray and lensing masses presently provides the most significant constraints on the level of non-thermal pressure in the ICM. Zhang et al. (2008) obtain a ratio of $1.09 \pm 0.08$ between the weak lensing and X-ray mass estimates extracted at $r_{{\rm 500}}$ from a sample of 19 clusters. At the same overdensity, (Mahdavi et al. 2008) find that the ratio between X-ray and lensing is $0.78 \pm 0.09$ ( $0.85~\pm~0.10$ after correction for excess structure along the line of sight) for a sample of 18 clusters. The fairly good agreement between these recent observational results and the values found in the present analysis and other similar studies therefore implies a low level of turbulence present in the ICM, thereby suggesting that the gas physics outside the cluster cores is well approximated by the simulations presented here.

While we have shown that at present both theoretical models and observations seem to indicate a deviation from the hydrostatic equilibrium of the order of $\sim$$10 \%$, further investigations are needed in order to draw robust conclusions on this important issue. Future measurements of X-ray and lensing masses for large samples and measurements of the ICM velocity structure will provide valuable information on the validity of the hydrostatic equilibrium and therefore on the reliability of X-ray clusters as cosmological probes.

Moreover, these measurements will also be of particular relevance for constraining gas physics in galaxy clusters. The amount of physical viscosity is in particular a key parameter which quenches the level of turbulence present in the ICM. Physical viscosity has been incorporated in hydrodynamic simulations of galaxy clusters by Sijacki & Springel (2006) and it was shown that even a modest amount of physical viscosity has significant consequences on ICM properties. Therefore X-ray mass estimates, extracted from a statistically meaningful sample of hydrodynamical SPH simulations which include physical viscosity, could profitably be used to indirectly measure the level of ICM viscosity when contrasted against X-ray and weak lensing mass measurements.

These constraints will also likely have a significant impact on those scenarios in which the cooling flow problem is solved by providing additional heating to the gas through energy dissipation.

Acknowledgements

We acknowledge the anonymous referee for detailed and valuable comments that helped in improving the presentation of the manuscript. We would like to thank Daisuke Nagai and Alexey Vikhlinin for useful comments. R.P. thanks Stefano Borgani for stimulating discussions.

References

• Allen, S. W., Rapetti, D. A., Schmidt, R. W., et al. 2008, MNRAS, 383, 879 [NASA ADS] [CrossRef]
• Ameglio, S., Borgani, S., Pierpaoli, E., & Dolag, K. 2007, MNRAS, 382, 397 [NASA ADS] [CrossRef] (In the text)
• Arnaud, M., Pointecouteau, E., & Pratt, G. W., 2005, A&A, 441, 893 [NASA ADS] [CrossRef] [EDP Sciences]
• Binney, J., & Tremaine, S. 1987 (Princeton, NJ: Princeton University Press) (In the text)
• Biviano, A., & Girardi, M. 2003, ApJ, 585, 205 [NASA ADS] [CrossRef]
• Borgani, S., Murante, G., Springel, V., et al. 2004, MNRAS, 348 [CrossRef],1078
• Buote, D. A. 2002, in Merging Processes in Cluster of Galaxies, ed. L., Feretti, I. M., Gioia, & G., Giovannini (Dordrecht Kluwer) [arXiv:astro-ph/0106057] (In the text)
• Buote, D. A., & Tsai, J.C., 1995, ApJ, 452, 522 [NASA ADS] [CrossRef] (In the text)
• Buote, D. A., & Tsai, J. C. 1996, ApJ, 458, 27 [NASA ADS] [CrossRef] (In the text)
• Burns, J. O., Hallman, E. J., Gantner, B., Motl, P. M., & Norman, M. L. 2008, ApJ, 675, 1125 [NASA ADS] [CrossRef] (In the text)
• Cavaliere A., & Fusco-Femiano R. 1976, A&A, 49, 137 [NASA ADS] (In the text)
• Ciotti, L., & Pellegrini, S. 2008, MNRAS, 387, 902 [NASA ADS] [CrossRef]
• Davé, R., Katz, N., & Weinberg, D. H. 2002, ApJ, 579, 23 [NASA ADS] [CrossRef]
• De Grandi, S., & Molendi, S. 2002, ApJ, 567, 163 [NASA ADS] [CrossRef]
• Dolag, K., Vazza, F., Brunetti, G., & Tormen, G. 2005, MNRAS, 364, 753 [NASA ADS] [CrossRef]
• Ettori, S., De Grandi, S., & Molendi, S. 2002, A&A, 391, 841 [NASA ADS] [CrossRef] [EDP Sciences]
• Evrard, A. E. 1990, ApJ, 363, 349 [NASA ADS] [CrossRef]
• Evrard, A. E., Mohr, J. J., Fabricant, D. G., & Geller, M. J. 1993, ApJ, 419, L9 [NASA ADS] [CrossRef]
• Evrard, A. E., Metzler, C. A., & Navarro, J. F. 1996, ApJ, 469, 494 [NASA ADS] [CrossRef]
• Finoguenov, A., Reiprich, T. H., & Böhringer, H. 2001, A&A, 368, 749 [NASA ADS] [CrossRef] [EDP Sciences]
• Gardini, A., Rasia, E., Mazzotta, P., et al. 2004, MNRAS, 351, 505 [NASA ADS] [CrossRef]
• Goodman, J., & Hernquist, L. 1991, ApJ, 378, 637 [NASA ADS] [CrossRef] (In the text)
• Haiman, Z., Mohr, J. J., & Holder, G. P. 2001, ApJ, 553, 545 [NASA ADS] [CrossRef]
• Jeltema, T. E., Canizares, C. R., Bautz, M. W., & Buote, D. A. 2005, ApJ, 624, 606 [NASA ADS] [CrossRef] (In the text)
• Jeltema, T. E., Hallman, E. J., Burns, J. O., & Motl, P. M. 2008, ApJ, 681, 167 [NASA ADS] [CrossRef]
• Jia, S. M., Boehringer, H., Pointecouteau, E., Chen, Y., & Zhang, Y. Y. 2008, [arXiv:0806.1575] (In the text)
• Kay, S. T., Thomas, P. A., & Theuns, T. 2003, MNRAS, 343, 608 [NASA ADS] [CrossRef]
• Kay, S. T., Thomas, P. A., Jenkins, A., & Pearce, F. R. 2004, MNRAS, 355, 1091 [NASA ADS] [CrossRef]
• Kay, S. T., da Silva, A. C., Aghanim, N., et al. 2007, MNRAS, 377, 317 [NASA ADS] [CrossRef]
• Kravtsov, A. V., Nagai, D., & Vikhlinin, A. A. 2005, ApJ, 625, 588 [NASA ADS] [CrossRef]
• Lau, E., Kravtsov, A. V., & Nagai, D. 2007, in prep. (In the text)
• Leccardi, A., & Molendi, S. 2008, A&A, 486, 359 [NASA ADS] [CrossRef] [EDP Sciences]
• Lewis, G. F., Babul, A., Katz, N., et al. 2000, ApJ, 536, 623 [NASA ADS] [CrossRef]
• Mahdavi, A., Hoekstra, H., Babul, A., & Henry, J. P. 2008, MNRAS, 128
• Markevitch, M., & Vikhlinin, A. 2007, Phys. Rep., 443, 1 [NASA ADS] [CrossRef] (In the text)
• Markevitch, M., Forman, W. R., Sarazin, C. L., & Vikhlinin, A. 1998, ApJ, 503, 77 [NASA ADS] [CrossRef]
• Mathiesen, B., & Evrard, A. E., 2001, ApJ, 546, 100 [NASA ADS] [CrossRef]
• Mazzotta, P., Rasia, E., Moscardini, L., & Tormen, G. 2004, MNRAS, 354, 10 [NASA ADS] [CrossRef]
• Majumdar, S., & Mohr, J. J. 2004, ApJ, 613, 41 [NASA ADS] [CrossRef]
• Monaghan, J. J. 2005, Rep. Prog. Phys., 68, 1703 [NASA ADS] [CrossRef] (In the text)
• Morris, J. P., & Monaghan, J. J. 1997, J. Comp. Phys., 136, 41 [NASA ADS] [CrossRef] (In the text)
• Muanwong, O., Thomas, P. A., Kay, S. T., & Pearce, F. R. 2002, MNRAS, 336, 527 [NASA ADS] [CrossRef]
• Nagai, D., Vikhlinin, A., & Kravtsov, A. V. 2007, ApJ, 655, 98 [NASA ADS] [CrossRef]
• Nagai, D., Kravtsov, A. V., & Vikhlinin, A. 2007, ApJ, 668, 1 [NASA ADS] [CrossRef] (In the text)
• Pearce, F. R., Thomas, P. A., Couchman, H. M. P., & Edge, A. C. 2000, MNRAS, 317, 1029 [NASA ADS] [CrossRef]
• Piffaretti, R., Jetzer, Ph., Kaastra, J. S., & Tamura, T. 2005, A&A, 433, 101 [NASA ADS] [CrossRef] [EDP Sciences]
• Pratt, G. W., Böhringer, H., Croston, J. H., et al. 2007, A&A, 461, 71 [NASA ADS] [CrossRef] [EDP Sciences]
• Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992, Numerical Recipe (Cambridge: Cambridge Univ. Press) (In the text)
• Rasia, E., Tormen, G., & Moscardini, L. 2004, MNRAS, 351, 237 [NASA ADS] [CrossRef]
• Rasia, E., Mazzotta, P., Borgani, S., et al. 2005, ApJ, 618, L1 [NASA ADS] [CrossRef]
• Rasia, E., Ettori, S., Moscardini, L., et al. 2006, MNRAS, 369, 2013 [NASA ADS] [CrossRef]
• Reiprich, T. H., & Böhringer, H. 2002, ApJ, 567, 716 [NASA ADS] [CrossRef]
• Reiprich, T. H., Hudson, D. S., Zhang, Y. Y., et al. 2008, [arXiv:0806.2920] (In the text)
• Richstone, D., Loeb, A., & Turner, E.L. 1992, ApJ, 393, 477 [NASA ADS] [CrossRef]
• Riess, A. G., Strolger, L. G., Tonry, J., et al. 2004, ApJ, 607, 665 [NASA ADS] [CrossRef]
• Riess, A. G., Strolger, L. G., Casertano, S., et al. 2007, ApJ, 659, 98 [NASA ADS] [CrossRef]
• Rines, K., Geller, M. J., Kurtz, M. J., & Diaferio, A. 2003, AJ, 126, 2152 [NASA ADS] [CrossRef]
• Sanderson, A. J. R., Ponman, T. J., & O'Sullivan, E. 2006, MNRAS, 372, 1496 [NASA ADS] [CrossRef]
• Sarazin, C. L. 1986, Rev. Mod. Phys., 58, 1 [NASA ADS] [CrossRef] (In the text)
• Sijacki, D., & Springel, V. 2006, MNRAS, 371, 1025 [NASA ADS] [CrossRef] (In the text)
• Smith, G. P., Smail, I., Kneib, J. P., et al. 2002, MNRAS, 330, 1 [NASA ADS] [CrossRef]
• Snowden, S. L., Mushotzky, R. F., Kuntz, K. D., & Davis, D. S. 2008, A&A, 478, 615 [NASA ADS] [CrossRef] [EDP Sciences]
• Spergel, D. N., Bean, R., Doré, O., et al. 2007, ApJS, 170, 377 [NASA ADS] [CrossRef] (In the text)
• Tonry, J. L., Schmidt, B. P., Barris, B., et al. 2003, ApJ, 594, 1 [NASA ADS] [CrossRef]
• Tornatore, L., Borgani, S., Springel, V., et al. 2003, MNRAS, 342, 1025 [NASA ADS] [CrossRef]
• Valdarnini, R. 2002, ApJ, 567, 741 [NASA ADS] [CrossRef] (In the text)
• Valdarnini, R. 2003, MNRAS, 339, 1117 [NASA ADS] [CrossRef]
• Valdarnini, R. 2006, New Astron., 12, 71 [NASA ADS] [CrossRef]
• Vazza, F., Tormen, G., Cassano, R., Brunetti, G., & Dolag, K. 2006, MNRAS, 369, L14 [NASA ADS]
• Vikhlinin, A. 2006, ApJ, 640, 710 [NASA ADS] [CrossRef]
• Vikhlinin, A., McNamara, B. R., Forman, W., et al. 1998, ApJ, 502, 558 [NASA ADS] [CrossRef] (In the text)
• Vikhlinin, A., Markevitch, M., Murray, S. S., et al. 2005, ApJ, 628, 655 [NASA ADS] [CrossRef]
• Vikhlinin, A., Kravtsov, A., Forman, W., et al. 2006, ApJ, 640, 691 [NASA ADS] [CrossRef]
• Voit, G. M. 2005, RMP, 77, 207 [CrossRef] (In the text)
• Yoshikawa, K., Jing, Y. P. & Suto, Y. 2000, ApJ, 535, 593 [NASA ADS] [CrossRef]
• Zhang, Y. -Y., Finoguenov, A., Boehringer, H., et al. 2008, A&A, 482, 451 [NASA ADS] [CrossRef] [EDP Sciences]