Discovery of a New Fundamental Plane Dictating Galaxy Cluster Evolution from Gravitational Lensing

N-body numerical simulations show that the density profile of dark-matter halos in galaxy clusters can be well described by the Navarro–Frenk–White (Navarro et al. 1997; NFW, hereafter) density profile: ${\rho }_{\mathrm{DM}}\propto {(r/{r}_{s})}^{-1}{(1+r/{r}_{s})}^{-2}$, where r is the clustercentric distance and r_s is the characteristic or scale radius. We define the mass inside r_s as M_s. Recent higher-resolution simulations have shown that the internal structure of dark-matter halos reflects the growth history of the halos (Wechsler et al. 2002; Zhao et al. 2003; Ludlow et al. 2013; Correa et al. 2015; More et al. 2015). They show that in the early "fast-rate growth" phase, halos grow rapidly through massive matter accumulation. This growth is often associated with phenomena that erase the previous internal structure of the halos, such as major mergers with other halos. In the subsequent "slow-rate growth" phase, halos gradually grow through moderate matter accretion from their surroundings. There are multiple definitions of the formation time of a halo. One example is the time at which the mass of the main progenitor equals the characteristic mass M_s of the z = 0 halo (Ludlow et al. 2013; Correa et al. 2015), and it approximately represents the end of the fast-growth phase and the transitioning toward the slow-growth phase. Only the outskirts of the halos (r > r_s) gradually grow in the latter phase (Wechsler et al. 2002; Zhao et al. 2003; Ludlow et al. 2013; Correa et al. 2015; More et al. 2015). Thus, dark-matter halos are assembled from the inside out, and the results of the numerical simulations can be interpreted such that the characteristic radius r_s and the mass M_s preserve a memory of the formation time of the halo (Wechsler et al. 2002; Zhao et al. 2003; Ludlow et al. 2013; Correa et al. 2015). In this "inside-out" scenario of halo formation, halos take a range of characteristic densities (${\rho }_{s}\equiv 3{M}_{s}/(4\pi {r}_{s}^{3})$); older halos tend to be more concentrated and have larger characteristic densities, which reflects the higher average density of the universe in the past (Navarro et al. 1997; Wechsler et al. 2002; Zhao et al. 2003; Ludlow et al. 2013; Correa et al. 2015; More et al. 2015). This scenario is in contrast with the classical approach in which halos are continuously modified, even by minor mergers, and constantly changing their profiles so that dark-matter halos lose the memory of their epoch of formation (Gunn & Gott 1972; Press & Schechter 1974).

If the inside-out halo growth scenario is correct, then we would expect that the formation time not only reflects the structural parameters (r_s and M_s) in the form of the characteristic density (ρ_s), but also influences the properties of the X-ray intracluster gas. However, its quantitative dependence is not obvious because the gas is collisional matter in contrast with dark matter. The hot gas is expected to be heated mostly via merger shocks produced when smaller halos fall into the halo (Rasia et al. 2011; Kravtsov & Borgani 2012). However, it is difficult to directly observe the heating process, as the shocks are often located at the outskirts of clusters (Miniati et al. 2000; Ryu et al. 2003), where the gas emission is very faint. We here investigate correlations between the halo parameters (r_s, M_s) and the average X-ray gas temperature T_X because the temperature is supposedly sensitive to the depth of the potential well and the past heating process (Eke et al. 1998). Since the emissivity of the X-ray gas is proportional to the gas density squared, the average measured temperature of a cluster mainly reflects the temperature in the region (r ≲ r_s), where the density is high. While there were previous studies that attempted to investigate correlations among a certain combination of three cluster structural parameters (Schaeffer et al. 1993; Adami et al. 1998; Fujita & Takahara 1999a; Lanzoni et al. 2004; Ota et al. 2006), they adopted parameters such as the galaxy luminosities that are not directly related to the structure of halos.

In the paper, we assume a spatially flat ΛCDM cosmology with Ω_m = 0.27, Ω_Λ = 0.73, and the Hubble constant of H₀ = 70 km s⁻¹ Mpc⁻¹ throughout this paper.

We study the Cluster Lensing And Supernova survey with Hubble (CLASH) observational data set that includes 20 massive clusters, most of which are apparently relaxed, X-ray regular systems⁹ (Postman et al. 2012; Meneghetti et al. 2014). The range of redshifts is 0.187–0.686, and their median redshift is 0.377 (Umetsu et al. 2016). Lensing constraints on the NFW characteristic radius r_s and M_s were obtained from a joint analysis of strong-lensing, weak-lensing shear and magnification data of background galaxies (Umetsu et al. 2016). Their analysis is based on 16-band Hubble Space Telescope observations (Zitrin et al. 2015) and wide-field multi-color imaging taken primarily with Suprime-Cam on the Subaru Telescope (Umetsu et al. 2014). The core-excised X-ray temperatures of the clusters were taken from Postman et al. (2012), in which the temperatures are measured for the region of 50–500 kpc from the cluster centers excluding the cool core at the center of the clusters. We chose the outer radius considering the field of view of the Chandra satellite and the completeness of the data. Since most of the X-ray emissions come from <500 kpc, the increase of the radius does not affect the results. The cluster data are shown in Table 1.

Table 1.  Cluster Data

Cluster	z	rs	r200	Ms	M200	TX
		(kpc)	(kpc)	(${10}^{14}\,{M}_{\odot }$)	(${10}^{14}\,{M}_{\odot }$)	(keV)
Abell 383	0.187	${304}_{-97}^{+159}$	${1800}_{-189}^{+209}$	${1.4}_{-0.5}^{+1.0}$	${7.9}_{-2.2}^{+3.1}$	6.5 ± 0.24
Abell 209	0.206	${834}_{-192}^{+243}$	${2238}_{-172}^{+161}$	${5.2}_{-1.6}^{+2.2}$	${15.4}_{-3.3}^{+3.6}$	7.3 ± 0.54
Abell 2261	0.224	${682}_{-170}^{+232}$	${2542}_{-188}^{+192}$	${5.8}_{-1.8}^{+2.7}$	${22.9}_{-4.7}^{+5.6}$	7.6 ± 0.30
RX J2129.7+0005	0.234	${294}_{-89}^{+133}$	${1626}_{-154}^{+163}$	${1.1}_{-0.4}^{+0.7}$	${6.1}_{-1.6}^{+2.0}$	5.8 ± 0.40
Abell 611	0.288	${560}_{-172}^{+250}$	${2189}_{-208}^{+204}$	${3.8}_{-1.5}^{+2.3}$	${15.6}_{-4.0}^{+4.8}$	7.9 ± 0.35
MS 2137−2353	0.313	${784}_{-357}^{+557}$	${2064}_{-286}^{+261}$	${4.7}_{-2.6}^{+5.2}$	${13.4}_{-4.9}^{+5.8}$	5.9 ± 0.30
RX J2248.7−4431	0.348	${643}_{-246}^{+422}$	${2267}_{-261}^{+282}$	${4.9}_{-2.3}^{+4.8}$	${18.5}_{-5.7}^{+7.8}$	12.4 ± 0.60
MACS J1115.9+0129	0.352	${738}_{-196}^{+249}$	${2186}_{-174}^{+161}$	${5.1}_{-1.7}^{+2.4}$	${16.6}_{-3.7}^{+4.0}$	8.0 ± 0.40
MACS J1931.8−2635	0.352	${501}_{-221}^{+441}$	${2114}_{-311}^{+355}$	${3.5}_{-1.8}^{+4.6}$	${15.0}_{-5.7}^{+8.9}$	6.7 ± 0.40
RX J1532.9+3021	0.363	${293}_{-114}^{+433}$	${1544}_{-210}^{+191}$	${1.2}_{-0.5}^{+1.6}$	${5.9}_{-2.1}^{+2.5}$	5.5 ± 0.40
MACS J1720.3+3536	0.391	${505}_{-162}^{+248}$	${2055}_{-204}^{+204}$	${3.4}_{-1.4}^{+2.3}$	${14.4}_{-3.9}^{+4.7}$	6.6 ± 0.40
MACS J0416.1−2403	0.396	${642}_{-156}^{+201}$	${1860}_{-154}^{+146}$	${3.4}_{-1.1}^{+1.5}$	${10.7}_{-2.4}^{+2.7}$	7.5 ± 0.80
MACS J0429.6−0253	0.399	${394}_{-143}^{+238}$	${1792}_{-208}^{+225}$	${2.1}_{-0.9}^{+1.8}$	${9.6}_{-3.0}^{+4.1}$	6.0 ± 0.44
MACS J1206.2−0847	0.440	${587}_{-176}^{+248}$	${2181}_{-178}^{+165}$	${4.6}_{-1.7}^{+2.4}$	${18.1}_{-4.1}^{+4.4}$	10.8 ± 0.60
MACS J0329.7−0211	0.450	${254}_{-63}^{+95}$	${1697}_{-127}^{+129}$	${1.4}_{-0.4}^{+0.6}$	${8.6}_{-1.8}^{+2.1}$	8.0 ± 0.50
RX J1347.5−1145	0.451	${840}_{-239}^{+339}$	${2684}_{-230}^{+226}$	${9.8}_{-3.6}^{+5.6}$	${34.2}_{-8.1}^{+9.4}$	15.5 ± 0.60
MACS J1149.5+2223	0.544	${1108}_{-291}^{+404}$	${2334}_{-178}^{+169}$	${10.8}_{-3.7}^{+5.4}$	${25.0}_{-5.3}^{+5.8}$	8.7 ± 0.90
MACS J0717.5+3745	0.548	${1300}_{-271}^{+347}$	${2387}_{-165}^{+154}$	${13.2}_{-3.9}^{+5.3}$	${26.8}_{-5.2}^{+5.6}$	12.5 ± 0.70
MACS J0647.7+7015	0.584	${468}_{-160}^{+254}$	${1884}_{-192}^{+189}$	${3.3}_{-1.3}^{+2.3}$	${13.7}_{-3.8}^{+4.6}$	13.3 ± 1.80
MACS J0744.9+3927	0.686	${574}_{-192}^{+269}$	${1982}_{-185}^{+179}$	${4.9}_{-2.0}^{+3.1}$	${17.9}_{-4.6}^{+5.3}$	8.9 ± 0.80

Download table as:  ASCIITypeset image

In Figure 1(a), we show the distribution of clusters in the (log r_s, log M_s, log T_X) space. We see that the data points are closely distributed on a plane that we have determined using a principal component analysis (PCA) to minimize deviations of the data points from the plane (see the Appendix). The arrow P₁ shows the direction on the plane in which the data are most extended, and the arrow P₂ is perpendicular to P₁ on the plane. The plane normal is represented by P₃. The dispersion around the plane or the thickness of the plane is shown in Figure 1(b) and it is only ${0.045}_{-0.007}^{+0.008}$ dex (all uncertainties are quoted at the 1σ confidence level unless otherwise mentioned). The thickness is comparable to that of the well-known fundamental plane for elliptical galaxies in the space of the surface brightness, the effective radius, and the velocity dispersion (∼0.06 dex; e.g., La Barbera et al. 2008; Hyde & Bernardi 2009). In Figure 1(c), we show error bars for individual clusters. In the vertical direction (T_X), we show the temperature errors in Table 1. In the horizontal direction, however, the errors of r_s and M_s are highly correlated, and we display them as a single bar. That is, for each cluster, we draw a bar connecting (r_s^u, M_s^u) and (r_s^l, M_s^l), where the superscripts u and l are the upper and the lower limits shown in Table 1, respectively. Note that we have properly accounted for the correlation for each cluster using the joint posterior probability distribution of the NFW parameters (mass and concentration) when we calculate the plane parameters (see the Appendix). Thus, the actual error is not represented by a single bar in a precise sense. The tight planar distribution in Figure 1 indicates that the structure of the dark-matter halos (r_s and M_s) did make a direct influence on the properties of the intracluster gas (T_X). In the context of the inside-out halo growth scenario, the most natural interpretation of our findings is that the intracluster gas was heated up to around T_X in the fast-rate growth phase when the shape of the potential well (r_s and M_s) was established, and that the gas preserves the memory of the cluster formation as is the case of the dark-matter halo structure.

Zoom In Zoom Out Reset image size

The plane is described by a log r_s + b log M_s + c log T_X =const., with $a={0.76}_{-0.05}^{+0.03}$, $b=-{0.56}_{-0.02}^{+0.02}$, and $c={0.32}_{-0.09}^{+0.10}$. Likelihood contours of the parameters describing the direction of the plane normal, P₃ = (a, b, c), are shown in Figure 2. The estimation of the errors is described in the Appendix. If the intracluster gas at r ≲ r_s simply preserves its pressure equilibrium state at the cluster formation, the gas temperature should reflect the potential depth of the dark-matter halo at the formation. Thus, one may expect that the gas temperature follows the virial theorem in a narrow sense ("virial expectation") at that time, T_X ∝ M_s/r_s, which is one of the main assumptions for the self-similar scaling relations of clusters. The resulting plane, however, is significantly tilted from this virial expectation (Figure 2) and is represented by ${T}_{{\rm{X}}}\,\propto {M}_{s}^{-b/c}{r}_{s}^{-a/c}\propto {M}_{s}^{1.8\pm 0.5}/{r}_{s}^{2.3\pm 0.7}$. Our findings show that the temperature T_X is more sensitive to the depth of the gravitational potential represented by M_s/r_s than the canonical virial expectation because ${T}_{{\rm{X}}}\propto {M}_{s}^{-0.2}{r}_{s}^{-0.3}{({M}_{s}/{r}_{s})}^{2.0}$. In other words, clusters with a deeper potential well tend to have higher temperatures T_X than the virial expectation, or visa versa.

Zoom In Zoom Out Reset image size

Our CLASH sample includes only 20 clusters, and our results may be affected by observational biases. Here, we examine the results of numerical simulations to properly interpret the observations in the context of the CDM cosmology and discuss possible selection bias. First, we analyzed the outputs of MUSIC N-body/hydrodynamical simulations (Meneghetti et al. 2014). These are adiabatic; that is, they do not include any nongravitational effects such as feedback from active galactic nuclei (AGNs) or from supernovae (SNe), and there is no radiative cooling.

The details of the simulations for the MUSIC sample are given in Meneghetti et al. (2014). The MUSIC sample is obtained by resimulating halos selected from the MultiDark cosmological simulation (Prada et al. 2012) in order to achieve a higher resolution. The parallel TREEPM+SPH GADGET code (Springel 2005) is used for the resimulations. The mass resolution for the dark-matter particles is ${m}_{\mathrm{DM}}=9.01\,\times {10}^{8}\ {h}^{-1}\,{M}_{\odot }$ and that for the gas particles is ${m}_{\mathrm{SPH}}=1.9\,\times {10}^{8}\ {h}^{-1}\,{M}_{\odot }$, where the Hubble constant is written as H₀ = 100 h km s⁻¹ Mpc⁻¹ and h = 0.7. The gravitational softening is 6 h⁻¹ kpc for both the gas and dark-matter particles in the high-resolution regions. We select all of the 402 clusters at z = 0.25 with ${M}_{200}\gt 2\times {10}^{14}\ {h}^{-1}\,{M}_{\odot }$ regardless of dynamical state, where M₂₀₀ is the mass enclosed within a sphere of radius (r₂₀₀) within which the mean overdensity equals 200 times the critical density of the universe. We compute the mass-weighted temperature including the core. The mass-weighted formulation is the most appropriate to evaluate the thermal energy of the X-ray gas to be included in the virial theorem. In addition, we kept the core because these simulations are nonradiative and thus do not present cool-core features. The scale radius r_s is obtained by fitting the total density distribution (gas+dark matter) with the NFW profile up to r₂₀₀. The mass M_s is then derived as the mass enclosed by a sphere of radius r_s.

We see that 402 simulated MUSIC clusters at the redshift of z = 0.25 form a plane in the (log r_s, log M_s, log T_X) space (Figure 3). In Figure 1(b), we project the simulated clusters on the cross-section of the observed plane, showing that the two sets of data are distributed around the same plane (P₃ = 0), although the band-like distribution of the MUSIC data is slightly tilted (see Figure 2) while passing through the origin. Many of them are found at smaller P₁, because the average radius and mass of the MUSIC clusters are smaller than those of the observed clusters by a factor of a few. The plane angle is consistent with the observed one (90% confidence level) and deviates from the virial expectation (Figure 2). The dispersion around the plane for the simulated clusters is 0.025 dex and is even smaller than the observed one (${0.045}_{-0.007}^{+0.008}$ dex; Table 2). Since unrelaxed clusters tend to have disturbed internal structure, they are expected to increase the dispersion. However, even when we choose the 20% most unrelaxed (UR) clusters in the sample, it is only 0.033. Therefore, the slightly larger observed thickness of the plane is unlikely ascribable to the dynamical state of the systems. The direction of the plane for the UR clusters is also not much different from that of the full sample (Table 3). Although our CLASH clusters are relatively massive and relaxed, these results of the numerical simulations show that the selection bias should not significantly affect the derived plane parameters.

Zoom In Zoom Out Reset image size

Table 2.  Cluster Samples

	Observation	Simulations
		MUSIC	NF0	FB0	FB1	FB0+FB1
Nongravitational Effects	⋯	no	no	yes	yes	yes
Redshift	${0.377}_{-0.190}^{+0.309}$	0.25	0	0	1	0+1
Dispersion around the Plane (dex)	${0.045}_{-0.007}^{+0.008}$	0.025a	0.023	0.031	0.035	0.037

Note.

^aFor the 20% most relaxed (RE) and unrelaxed (UR) clusters in the sample, it is 0.015 and 0.033, respectively. The classification is based on fit residuals to the NFW profile.

Download table as:  ASCIITypeset image

Table 3.  Plane Vectors

Sample	P1	P2	P3
Observation	$({0.55}_{-0.02}^{+0.03},{0.82}_{-0.01}^{+0.01},{0.15}_{-0.06}^{+0.04})$	$(-{0.34}_{-0.07}^{+0.08},{0.07}_{-0.07}^{+0.07},{0.93}_{-0.03}^{+0.03})$	$({0.76}_{-0.05}^{+0.03},-{0.56}_{-0.02}^{+0.02},{0.32}_{-0.09}^{+0.10})$
MUSIC	(0.40,0.81,0.42)	(−0.60, −0.11,0.79)	(0.69, −0.57,0.44)
MUSIC (RE)a	(0.38,0.82,0.42)	(−0.62, −0.11,0.78)	(0.69, −0.56,0.47)
MUSIC (UR)a	(0.39,0.79,0.48)	(−0.58, −0.19,0.79)	(0.71, −0.59,0.39)
NF0	(0.43,0.83,0.36)	(−0.48, −0.13,0.87)	(0.77, −0.55,0.34)
FB0	(0.43,0.84,0.34)	(−0.51, −0.08,0.86)	(0.74, −0.54,0.39)
FB1	(0.37,0.84,0.40)	(−0.55, −0.16,0.82)	(0.75, −0.52,0.40)
FB0+FB1	(0.42,0.82,0.40)	(−0.53, −0.13,0.83)	(0.74, −0.56,0.38)
Observationb	$({0.27}_{-0.01}^{+0.01},{0.90}_{-0.03}^{+0.02},{0.35}_{-0.07}^{+0.06})$	$(-{0.18}_{-0.02}^{+0.03},-{0.31}_{-0.07}^{+0.06},{0.93}_{-0.02}^{+0.02})$	$({0.94}_{-0.00}^{+0.01},-{0.32}_{-0.01}^{+0.01},{0.08}_{-0.02}^{+0.02})$

Notes.

^aThe values for the 20% most relaxed (RE) and unrelaxed (UR) clusters in the MUSIC sample. The classification is based on fit residuals to the NFW profile. ^bFor the plane obtained in the space of (log r₂₀₀, log M₂₀₀, log T_X) (see Section 5.3).

Download table as:  ASCIITypeset image

To study the evolution in detail, we analyzed another set of simulation data (Rasia et al. 2015). Each of the samples named FB0 and NF0 consists of 29 massive clusters at z = 0, and contain both relaxed and unrelaxed ones. FB0 includes nongravitational effects and NF0 does not. The details of the simulations for samples FB0, FB1, and NB0 are given in Rasia et al. (2015) and Planelles et al. (2017). They are also carried out with the GADGET code (Springel 2005) but including an updated SPH scheme (Beck et al. 2016). The simulations consist in 29 Lagrangian regions around massive clusters with ${M}_{2}00\sim 1\mbox{--}30\times {10}^{14}\,{h}^{-1}\,{M}_{\odot }$ at z = 0. The simulations FB0 and FB1 include phenomena such as heating by AGNs and SNe in addition to radiative cooling, while NF0 is from nonradiative runs. Samples FB0 and NF0 consist of the clusters at z = 0, while sample FB1 refers to the runs at z = 1. The mass resolution for the dark-matter particles is ${m}_{\mathrm{DM}}\,=8.3\times {10}^{8}\ {h}^{-1}\,{M}_{\odot }$ and that for the initial gas particles is ${m}_{\mathrm{SPH}}=1.5\times {10}^{8}\ {h}^{-1}\,{M}_{\odot }$. The gravitational softening is 3.75 h⁻¹ kpc for the both gas and dark-matter particles in the high-resolution regions (Biffi et al. 2017). We derive r_s and M_s for the clusters in the FB0, FB1, and NB0 samples using the same method exploited for the MUSIC simulations. Since this sample, contrary to the previous sample, is built on radiative simulations we do exclude the core. As we did for the observed CLASH sample, the temperatures are obtained in the region between 50 and 500 kpc from the cluster centers. Although we use the mass-weighted temperature in the following discussion, we have confirmed that the results such as the plane angle and thickness are not significantly affected by the choice of the temperature weighting (e.g., spectroscopic-like temperature; Mazzotta et al. 2004) or the choice of the metric radius for temperature measurements (>500 kpc).

We find that each sample forms a plane whose angle is consistent with the observed one (Figure 2). The plane angles for FB0 and NF0 are almost the same, which means that the result is independent of the gas physics. The lack of dependence means that radiative cooling and SNe and AGN feedback counterbalance one another with the effect of not drastically changing the X-ray gas profile on a scale of r_s. Thus, even if our CLASH clusters are affected by some selection bias originating from gas physics (e.g., difference of AGN activities), the bias does not have a significant impact on the plane parameters. The thickness of the plane is, however, increased by the nongravitational effects. In fact, the dispersion around the plane for FB0 is 0.031 dex, which is larger than that for NF0 (0.023 dex), but is still smaller than the observed one even if the observational errors are considered (Table 2). In Figure 2, we also show the plane angle for clusters in the same simulation as FB0 but at z = 1 (sample FB1). Most of the clusters (25/29) are the progenitors of those in FB0. While more clusters should be in the fast-rate growth phase at z = 1, the plane angle is not much different from that at z = 0 (FB0). We find that the clusters in the samples FB0 and FB1 are virtually on the same plane (Figure 4); though, FB1 clusters tend to have smaller physical radii and masses. This indicates that the redshifts of the clusters are unlikely to impact the plane parameters. The dispersion around the common plane is 0.037 dex (FB0+FB1 in Table 2). This indicates that the clusters evolve on this unique plane along the direction of P₁, and that the evolution of cluster halo structure and the thermodynamic history of intracluster gas are strictly regulated by the plane.

Zoom In Zoom Out Reset image size

In general, current numerical simulations are realistically reproducing the observed scaling relations including their slopes (Truong et al. 2018). Any possible small discrepancy of the normalizations does not significantly affect the plane angle and the cluster evolution along the plane.

5.1. The Plane Angle Predicted by a Similarity Solution and Cluster Evolution

In this subsection, we attempt to explain the origin of the peculiar plane angle we found using an analytic solution. Bertschinger (1985) constructed a one-dimensional similarity solution for secondary infall and accretion onto an initially overdense perturbation in an Einstein–de Sitter (Ω₀ = 1) universe. For this solution, an object continues to grow, and the matter density ρ, pressure p at the radius r, and mass m inside r at the cosmological time t can be expressed by nondimensional functions (D, P, and M):

$$\begin{eqnarray}&&\rho (r,t)={\rho }_{{\rm{H}}}D(\lambda ),\\ &&p(r,t)={\rho }_{{\rm{H}}}{({r}_{\mathrm{ta}}/t)}^{2}P(\lambda ),\\ &&m(r,t)=(4\pi /3){\rho }_{{\rm{H}}}{r}_{\mathrm{ta}}^{3}M(\lambda ),\end{eqnarray} \tag{ 1 }$$

where r_ta(t) is the maximum radius that a mass shell reaches (the turnaround radius), ρ_H ∝ t⁻² is the density of the background universe, and λ = r/r_ta is the nondimensional radius. The turnaround radius is represented by r_ta = A_itat^8/9, where A_ita is the coefficient that depends on the overdense perturbation (Bertschinger 1985). The solution has an entropy integral,

$$\begin{eqnarray}&&P{(\lambda )D(\lambda )}^{-\gamma }M{(\lambda )}^{10/3-3\gamma }=\mathrm{const},\end{eqnarray} \tag{ 2 }$$

where γ = 5/3 is the adiabatic index. This relation holds even for a system with a mixture of gas and dark matter (Bertschinger 1985). From Equations (1) and (2), we have $p{\rho }^{-5/3}{m}^{-5/3}\propto {A}_{\mathrm{ita}}^{-3}$, which does not depend on t. The coefficient can be written as ${A}_{\mathrm{ita}}={r}_{\mathrm{ita}}/{t}_{\mathrm{ita}}^{8/9}$, where r_ita and t_ita are the turnaround radius and time of the overdense perturbation, respectively. Note that t_ita is much earlier than "the formation of a cluster" discussed in this study. Assuming that the evolution of the overdense perturbation follows a theory of a spherical collapse, they are represented by

$$\begin{eqnarray}&&{r}_{\mathrm{ita}}\propto {m}_{\mathrm{ita}}^{(n+5)/6},\qquad {t}_{\mathrm{ita}}\propto {m}_{\mathrm{ita}}^{(n+3)/4},\end{eqnarray} \tag{ 3 }$$

where m_ita is the mass scale of the perturbation and n is the local slope of the primordial matter power spectrum¹⁰ (Kaiser 1986; Peebles 1993), which is n ∼ −2 at cluster scales (Eisenstein & Hu 1998; Diemer & Kravtsov 2015). Thus, we obtain $p{\rho }^{-5/3}{m}^{-5/3}\propto {m}_{\mathrm{ita}}^{-5/6}$. Assuming that p ∝ ρT_X and $\rho \propto {M}_{s}/{r}_{s}^{3}$ at r ∼ r_s, the relation is ${r}_{s}^{2}{M}_{s}^{-7/3}{T}_{{\rm{X}}}\propto {m}_{\mathrm{ita}}^{-5/6}$. Here, we speculate that the structure of the NFW profile at r ≲ r_s reflects the overdense perturbation that initially collapsed in the similarity solution by Bertschinger (1985). In other words, the fast-rate growth of a dark-matter halo is related to the initial collapse. In fact, Correa et al. (2015) demonstrated that the characteristic density ρ_s of the NFW profile is proportional to the critical density of the background universe ρ_c at the time when the dark-matter halo transits from the fast-rate to the slow-rate growth phase (i.e., the halo formation time). On the other hand, since the initial collapse can be described as a simple spherical collapse of an overdense region, the typical density of the collapsed object is also proportional to ρ_c at the collapse time (Bertschinger 1985). This indicates that both the inner structure of the NFW profile and the overdense perturbation in the similarity solution follow the same evolution and scaling relation. Thus, we assume M_s ∝ m_ita and r_s ∝ r_ita, which leads to

$$\begin{eqnarray}&&{r}_{s}^{2}{M}_{s}^{-3/2}{T}_{{\rm{X}}}=\mathrm{const}.\end{eqnarray} \tag{ 4 }$$

The angle of this plane is shown in Figure 2 (SSol) and it is consistent with the observations. The essential point of the similarity solution is that clusters are not isolated but continuously growing through matter accretion from their outer environments. Therefore, additional contributions, such as the flux of inertia at the cluster surface, should be included in the virial theorem for the complete description of the dynamical state (Bertschinger 1985). In the ΛCDM model (Peebles 1993), the cosmological constant becomes non-negligible at z ≲ ((1−Ω_m)/Ω_m)^1/3 − 1 ∼ 0.39, which is close to the median redshift of our observational sample (Table 2). Thus, about half of the sample may be affected by the cosmological constant. However, although the density profiles of these objects may become steeper in the outskirts (Bertschinger 1985), the effect is not serious because we are interested in the inner region (r ≲ r_s).

This simple model may also explain the vector P₁. Since we assumed r_s ∝ r_ita and M_s ∝ m_ita, we obtain ${r}_{s}\propto {M}_{s}^{1/2}$ from Equation (3). The direction of the line ${r}_{s}\propto {M}_{s}^{1/2}$ on the log r_s–(log M_s plane is almost the same as that of P₁ projected on the log r_s–log M_s plane, especially for P₁ for the simulation samples (Table 3). Figure 5(a) shows that the MUSIC data points are actually distributed along the line ${r}_{s}\propto {M}_{s}^{1/2}$ (thick solid arrow) on the log r_s–log M_s plane. This means that the cluster's elongated distributions along P₁ (Figures 3 and 4) reflect the evolution of the typical overdense perturbation along ${r}_{\mathrm{ita}}\propto {m}_{\mathrm{ita}}^{1/2}$, which can also be interpreted as the evolution of clusters during the fast-rate growth phase. The direction of P₁ for the observational sample (first line of Table 3) is slightly different from those for the simulations. This may be because the observational sample is biased toward high temperature clusters.

Zoom In Zoom Out Reset image size

Cluster formation time is associated with the characteristic density (e.g., Fujita & Takahara 1999b; Ludlow et al. 2013; Correa et al. 2015). The dashed lines in Figure 5(a) are isochrones or ρ_s = const. The MUSIC data points are widely distributed along each isochrone, which reflects a variety of cluster masses for a given formation time. In other words, it reflects a variety of peak densities of initial density fluctuations of the universe. Individual clusters evolve approximately in the direction of P₁ (Figures 3 and 4), or in the direction to which ${M}_{s}{r}_{s}^{1/2}$ increases (Figure 5(a)), that is, along the line ${r}_{s}\propto {M}_{s}^{1/2}$ (${r}_{s}\propto {M}_{s}^{1/1.65}$ in a detailed study; Zhao et al. 2009). Figure 5(b) shows that ${M}_{s}{r}_{s}^{1/2}$ and T_X are correlated, which reflects that T_X evolves according to the structure evolution of dark halos, although the projection of the band-like distribution of clusters in Figure 3(a) onto the ${M}_{s}{r}_{s}^{1/2}$–T_X plane disperses the relation. In Figure 5(a), clusters become denser (having larger ρ_s) or older in the direction of the thin dotted arrow. In Figure 5(c), the correlation between ρ_s and T_X is not clear, because a possible correlation is obscured by the projection of clusters with various masses along the dashed lines in Figure 5(a) (see also Figure 3(a)). However, each cluster moves mainly along the direction of P₁, which shows that the cluster temperature increases as the density ρ_s decreases, although major mergers (A, B, and E in Figure 5(c)) derail the cluster significantly from the evolution along P₁. We emphasize that the MUSIC data distribution on the plane does not follow a single line but has a finite spread, and M_s has a distribution for a given formation time. In general, the correlation between two parameters of clusters is not necessarily represented by a line but is often represented by a broad band (Figure 5). This is because the mass scale of the initial density fluctuation of the universe has a distribution for a given peak height (e.g., Barkana & Loeb 2001). Since the peak height also has its own distribution, these result in a two-dimensional band-like distribution of clusters. In general, clusters move along the band on the plane in the direction of P₁ (see also the next subsection).

5.2. Stability of the Plane Against Mergers

The similarity of clusters is generally good because most of them are well represented by the NFW profile. The thinness of the plane reflects the excellent similarity because cluster structure is well described by the similarity solution discussed in Section 5.1. However, clusters occasionally experience mergers that might break the similarity.

In Figure 6(a), we follow the evolution of a cluster that undergoes three major mergers in the sample FB0+FB1. This cluster experiences mergers around the times indicated by A, B, and E. While M_s tends to correlate with r_s, the temperature T_X tends to inversely correlate with r_s. This can be explained as follows. When a large substructure is merging with a cluster, it does not dissolve immediately when it touches the viral radius of the cluster. The 3D cluster gas density profiles will, therefore, include the gas of the substructure as it moves within the cluster external atmosphere toward the inner regions. Very large substructures can even reach the center and cross it. This is the main origin of changing r_s and M_s: when the substructure is in the outskirts, the profile is flatter and thus r_s and M_s are larger (less concentrated cluster), while when the substructure is closer to the center, the object appears more concentrated, i.e., with smaller r_s and M_s. On the other hand, T_X increases for a moment, after the shock has time to propagate throughout a large region of the cluster.

Zoom In Zoom Out Reset image size

Thanks to the behavior of the three parameters r_s, M_s, and T_X, the cluster does not substantially deviate from the plane (P₃ = 0) even during a merger (Figure 6(b)), which contributes to the thinness of the plane. Figure 6 also shows that after the end of the last major merger indicated by the letter F, the three parameters do not change much and the cluster remains in almost the same position on the plane. This indicates that the current cluster structure is determined at the last major merger.

5.3. Cluster Distribution in the Space of (log r₂₀₀, log M₂₀₀, log T_X)

Contrary to the inside-out scenario, dark-matter halos possess only a common global density that is related to the virial overdensity in the classical picture of cluster formation (Gunn & Gott 1972; Press & Schechter 1974; Lacey & Cole 1993). A representative density is ρ₂₀₀, which is 200 times the critical density of the universe at the cluster redshift. Individual halos can thus be characterized by M₂₀₀, the mass enclosed within a sphere of radius r₂₀₀ within which the mean overdensity equals ρ₂₀₀. The mass M₂₀₀ is often regarded as the total mass of a cluster. The values of r₂₀₀ and M₂₀₀ for our observed clusters are shown in Table 1.

As an alternative parameter combination, we check the cluster distribution in the space of (log r₂₀₀, log M₂₀₀, log T_X). In Figure 7(a), the data points are distributed on a plane almost vertical to the log r₂₀₀–log M₂₀₀ plane, which reflects the obvious relation of ${M}_{200}=4\pi {\rho }_{200}{r}_{200}^{3}/3\propto {r}_{200}^{3}$ regardless of T_X. The redshift dependence of ρ₂₀₀ is small because the majority of the clusters in our sample are distributed in a relatively narrow range at low-intermediate redshifts, although the small redshift dependence slightly slants the plane. The dispersion in the direction of the plane normal (P₃) is small (Figure 7(b); ${0.016}_{-0.001}^{+0.001}$ dex). A weak correlation (P₁), in which clusters with a larger M₂₀₀ tend to have a larger T_X, is seen on the plane, but the dispersion (P₂) is relatively large (Figure 7(a); ${0.099}_{-0.011}^{+0.011}$ dex). Although this correlation P₁ has been studied as a mass–temperature relation (e.g., Sun et al. 2009; Connor et al. 2014), the large dispersion suggests that it is not a primary relation, and that the combination of parameters (r_s, M_s, T_X) rather than (r₂₀₀, M₂₀₀, T_X) is appropriate for studying the connection between the X-ray gas and the dark halo structure and the interplay between gas heating and gravitational collapse in detail. The well-known large scatter of the concentration (c₂₀₀ ≡ r₂₀₀/r_s)–mass (M₂₀₀) relation (e.g., Bullock et al. 2001) may be because r_s has less to do with r₂₀₀ and M₂₀₀.

Zoom In Zoom Out Reset image size

In this study, we showed that observational data of high-mass galaxy clusters form a plane in the three-dimensional logarithmic space of their characteristic radius r_s, mass M_s, and X-ray temperature T_X with a very small orthogonal scatter. Since the evolution history of a cluster is encoded in r_s and M_s, the tight correlation suggests that the gas temperature was determined at a specific cluster formation time. We also found that the plane is tilted with respect to T_X ∝ M_s/r_s, which is the plane expected in the case of simplified virial equilibrium. This strange plane angle can be explained by a similarity solution, which indicates that clusters are not isolated but continuously growing through matter accretion from their outer environments. In other words, the effects of the growth must be considered when the internal structure of clusters is discussed. We have shown that numerical simulations reproduce the observed plane and its angle, regardless of the gas physics implemented in the code. The simulations show that clusters evolve along the plane and they do not deviate much from the plane even during major mergers, which contributes to the overall thinness of the plane.

Further work is needed to understand the relation between the similarity solution and the formation of the NFW halo structure (Salvador-Solé et al. 1998; Williams et al. 2004). Large-scale high-resolution numerical simulations from high to low redshifts would be useful to study these topics. It would also be interesting to extend the observational sample to objects other than massive clusters in the nearby universe. Those objects include clusters at higher redshifts and less massive systems, such as galaxy groups and elliptical galaxies that present a break from the cluster self-similarity behavior (Ponman et al. 1999). It would also be interesting to study relationships between the cluster plane and galactic-scale ones (e.g., Goulding et al. 2016). Galaxy clusters have been used to constrain cosmological parameters, such as the amount of matter and dark energy, and to investigate the growth of large-scale structure. In particular, the evolution of the abundance of rare massive clusters above a given mass threshold is highly sensitive to the expansion history and the growth rate of mass density fluctuations (e.g., Vikhlinin et al. 2009). Thanks to the thinness of the plane we discovered, precise determinations of the fundamental plane can be used to calibrate cluster mass–observable relations, a key ingredient of the cluster cosmology (Fujita et al. 2018).

This work was supported by MEXT KAKENHI No. 15K05080 (Y.F.). K.U. acknowledges support from the Ministry of Science and Technology of Taiwan through the grants MoST 103-2112-M-001-030-MY3 and MoST 106-2628-M-001-003-MY3. E.R. acknowledges support from the ExaNeSt and EuroExa projects, funded by the European Union's Horizon 2020 research and innovation programme under grant agreements No. 671553 and No. 754337, respectively.

The plane that represents the planarly distributed data points in three-dimensional space can be obtained so that the deviations of the points from the plane are minimized. Here, we define three vectors (principal components) in the space of (x, y, z) ≡ (log r_s, log M_s, log T_X). The first component P₁ is defined as the direction to which the points have the largest variance. The second component P₂ is orthogonal to P₁ and is the direction to which the points have the largest variance under the orthogonal condition between P₁ and P₂. The third component P₃ is orthogonal both to P₁ and P₂, which means that the points have the least variance to that direction. That is, P₃ is the normal to the plane that represents the planarly distributed points.

We find the plane through a PCA. Assuming that each data point is given by (x_i, y_i, z_i), a covariant matrix can be defined as

$$\begin{eqnarray*}A=\displaystyle \sum _{i}\left(\begin{array}{ccc}{({x}_{i}-\bar{x})}^{2} & ({x}_{i}-\bar{x})({y}_{i}-\bar{y}) & ({x}_{i}-\bar{x})({z}_{i}-\bar{z})\\ ({x}_{i}-\bar{x})({y}_{i}-\bar{y}) & {({y}_{i}-\bar{y})}^{2} & ({y}_{i}-\bar{y})({z}_{i}-\bar{z})\\ ({x}_{i}-\bar{x})({z}_{i}-\bar{z}) & ({y}_{i}-\bar{y})({z}_{i}-\bar{z}) & {({z}_{i}-\bar{z})}^{2}\end{array}\right),\end{eqnarray*}$$

where $\bar{x}$, $\bar{y}$, and $\bar{z}$ are the average of x_i, y_i, and z_i, respectively. If the eigenvalues for the matrix A are designated as λ₁, λ₂, and λ₃ (${\lambda }_{1}\gt {\lambda }_{2}\gt {\lambda }_{3}$), the corresponding eigenvectors are the principal components P₁, P₂, and P₃, respectively. The dispersions of data in the directions of P₁, P₂, and P₃ are represented by $\sqrt{{\lambda }_{1}}$, $\sqrt{{\lambda }_{2}}$, and $\sqrt{{\lambda }_{3}}$, respectively.

For the observational sample, we estimate the uncertainties of the plane parameters (P₁, P₂, P₃, λ₁, λ₂, λ₃) by accurately propagating the errors in the observed three cluster parameters (r_s, M_s, T_X) using Monte-Carlo simulations. Among the three parameters (r_s, M_s, T_X), the errors in the NFW halo parameters (r_s, M_s) are highly correlated with each other. Thus, when we estimate the errors, we directly use the joint posterior probability distribution P(M₂₀₀, c₂₀₀) of the NFW parameters for each cluster, accounting for the full uncertainty in both mass and concentration. The posterior probability distributions for the individual clusters were obtained by Umetsu et al. (2016) using Markov-chain Monte-Carlo sampling. We perform Monte-Carlo realizations of the temperature T_X assuming random Gaussian errors. Finally, we generate a total of 10⁶ realizations of the data and derive errors of the plane parameters.