2.1 Reaching Virial Equilibrium

The evolution of galaxy systems at various scales is understood today through the hierarchical formation model (^{Peebles 1980}; ^{Padmanabhan 1993}; ^{Zakhozhay 2018}) in a ΛCDM scenario. According to this model, the smallest systems (galaxy groups and poor clusters) merge -via gravity- to form the largest ones (rich clusters and superclusters), being the observed substructuring a consequence of this process. Subsequently, these galaxy ensembles undergo a non-linear gravitational collapse during which a combination of physical and stochastic mechanisms (^{Lynden-Bell 1967}; ^{Saslaw 1980}; ^{Padmanabhan 1990}) drive the systems into a further state of dynamical relaxation -the virial equilibrium- in which their substructures vanish (^{Araya-Melo et al. 2009}).

Galaxy clustering is an irreversible process where, as galaxies accumulate, different mechanisms tend to increase the number of ways in which, in a statistical sense, internal energy can be distributed within the systems (^{Saslaw 1980}; ^{Saslaw & Hamilton 1984}). Once bound and during collapse, the galaxy system undergoes internal -dissipationless- processes (e.g., violent relaxation, phase mixing, energy equipartition, galaxy mergers, fading of substructures and density or temperature gradients, among others, see, ^{White 1996}; ^{Dehnen 2005}; ^{Binney & Tremaine 2008}) that lead to states of greater dynamical relaxation. All these processes are dominated by gravitational interactions that, during the relaxation time, tend to homogenize the spatial distribution of the member galaxies and to distribute their radial velocities in a Gaussian way -inside the clusters, galaxies are scattered randomly (achieving a quasi-Maxwellian velocity distributions as they tend to virialization, e.g., ^{Saslaw & Hamilton 1984}; ^{Sampaio & Ribeiro 2014}), causing the tendency for macroscopic motions to disappear. This requires the individualization of the motions of the galaxies, so that any substructure (e.g., accreted groups) will be ‘dissolved’ before the cluster virializes. The thermodynamic perturbations and the changes in the distribution of the ICM inside the clusters during their evolution also produce increases in the total entropy of these systems (^{Tozzi & Norman 2001}; ^{Voit 2005}).

The equilibrium state can be understood, in a first approximation, as that in which the gravitational collapse is supported by the effect of inertial -centrifugal or dispersion- forces, achieving relaxed internal configurations (called states of dynamical relaxation), that is, in which there are no unbalanced potentials, such as gravitational-driving forces, within galaxy systems. In this sense, collisionless systems in equilibrium are analogous to self-gravitating fluids because they support gravitational collapse through “pressure gradients” proportional to the velocity dispersion that, at each point, tend to disperse any local increase in particle density (^{Binney & Tremaine 2008}). Furthermore, dynamical equilibrium is characterized by the statistical equality of the cluster mass profiles obtained from different galaxy populations within the cluster (^{Lokas & Mamon 2003}; ^{Lima-Neto et al. 2005}), implying that all these populations are in equilibrium with the cluster potential according to the hydrodynamical equilibrium model by Jeans.

Throughout the evolution of an isolated cluster, its total internal energy U = K + W is conserved so that, as the member galaxies get closer together (spontaneous reduction of inter-particle distance r _ij by mutual attraction), the gravitational potential energy W decreases (becoming more negative) and consequently the total internal kinetic energy K increases. This is reflected observationally as an increase in the velocity dispersion of galaxies, which can stop the collapse. Thus, when the system finds some route to dynamical equilibrium, the increase in K at the expense of the decrease in W does not continue indefinitely, but evolves toward configurations in which the ratio

tends to one (b→1, ^{Saslaw & Hamilton 1984}; ^{Binney & Tremaine 2008}). The asymptotic limit, b = 1, of this tendency is the virial equilibrium:

a state in which dynamical parameters that characterize the global configuration of the system remain, at least temporarily, stationary (relaxed). The virial theorem expresses a statistical equilibrium between the temporal averages of the total internal kinetic and gravitational potential energies, i.e., 〈W〉τ=-2〈K〉τ, but these energies are also the result of an ensemble of interacting galaxies (^{Limber & Mathews 1960}), so it is reasonable to approximate the temporal averages by the average of the ensemble expressed in (2). Virial equilibrium is assumed to be the type of dynamical equilibrium reached by self-gravitating systems.

2.2 Stability

Going beyond the description of equilibrium in galaxy systems, we need to discuss if this equilibrium is stable and this topic requires involving the concept of ‘entropy’. As any thermodynamic system, self-gravitating systems progress in the sense of increasing entropy (^{Lifshitz & Pitaevskii 1981}; ^{Tremaine, H´enon & Lynden-Bell 1986}; ^{Pontzen & Governato 2013}) towards the state of dynamical equilibrium described above. Galaxy clustering simulations also confirm this fact (^{Saslaw & Hamilton 1984}; ^{Iqbal et al. 2006}, ²⁰¹¹). Due to the scattering experienced by galaxies in the phase-space of clusters, these become the regions within the LSS where first-order entropy production occurs by increasing the randomness of the motion of galaxies during their gravitational accumulation. Even if galaxy clusters are isolated, some entropy is generated -or produced- due to the presence of internal irreversibilities. The peculiarity here is that the virial equilibrium is not unique, but only a metastable equilibrium state (^{Antonov 1962}; ^{Lynden-Bell & Wood 1968}; ^{Padmanabhan 1990}; ^{Chavanis et al. 2002}). This means that the entropy of self-gravitating systems can grow indefinitely without reaching a global maximum, that is, the virial equilibrium is only a state of local extreme of entropy (^{Padmanabhan 1989}).

The dynamical equilibrium of a galaxy cluster can be disturbed if it actively interacts with its surroundings, for example through mergers with other clusters and/or group accretions or tidal forces. As a result, the cluster takes a route towards a new equilibrium, in a state of higher entropy. Concerning the impact of the interaction, if the accreted groups are very small, the clusters can be kept unperturbed in states close to equilibrium. In more extreme cases, the merger of two massive clusters completely removes the systems from their equilibrium. In dense environments, such as supercluster cores, the accretion of galaxies and groups by the most massive clusters continually disturbs their dynamical states. On the other hand, in less dense environments, such as along filaments or edges near voids, clumpy clusters evolve as quasi-isolated systems, reaching dynamical relaxation possibly faster, without many significant disturbances, but accessing lower entropy levels compared to clusters in “busy” environments. That is, depending on the cosmological environment inside the LSS in which a cluster evolves, its relaxation process may be affected several times -or not- given the amount of matter available in its surroundings.

The most stable states of a galaxy cluster -of mass M and radius R- are favored when the density contrast between its center and its edge is ρ0/ρ(R)<709 and the internal energy is U>-0.335GM2/R (the Antonov instability or gravothermal catastrophe ^{Antonov 1962}; ^{Lynden-Bell & Wood 1968}), and under these conditions the virial equilibrium corresponds to a local maximum of entropy (^{Padmanabhan 1989}, ¹⁹⁹⁰). In fact, the most stable dynamical configurations that self-gravitating systems can access are those in which the particle -or matter- distribution settles on a core-halo structure (^{Binney & Tremaine 2008}; ^{Chavanis et al. 2002}). This theoretical structure, characterized by a collapsed core coexisting with a regular halo, has also been confirmed by simulations (^{Cohn 1980}; ^{Balberg, Shapiro & Inagaki 2002}) and recognized by the distributions of galaxies (^{Sarazin 1988}; ^{Adami et al. 1998}) and the ICM in observed clusters (^{Tozzi & Norman 2001}; ^{Cavagnolo et al. 2009}).

2.3 Estimating the Evolutionary State of Galaxy Systems

Observationally, a cluster close to dynamical equilibrium is distinguished from a non-relaxed one by exhibiting a more regular morphology (^{Sarazin 1988}), both in the optical and X-rays, and a more homogeneous -projected- spatial distribution of member galaxies (without the presence of significant substructures, ^{Caretta et al. 2023}, and references therein), as well as by having a more isotropic galaxy velocity distribution (or Gaussian in the line-of-sight, ^{Girardi & Mezzetti 2001}; ^{Sampaio & Ribeiro 2014}).

Concerning global morphology, the most dynamically relaxed clusters tend to present low ellipticity shapes in the projected distribution of galaxies and X-ray surface brightness maps, with a -possible- single peak at their centers. Several morphological classifications have been proposed following this premise (^{Sarazin 1988}, and references therein). There are also indicators of the internal structure of the galaxy systems, such as the radial profile and the degree of concentration of galaxies (^{Adami et al. 1998}; ^{Tully 2015}; ^{Kashibadze et al. 2020}), as well as a measure of the presence of substructures within them using 1D, 2D and 3D tests (^{Geller & Beers 1982}; ^{Dressler & Shectman 1988}; ^{Caretta et al. 2023}).

In this sense, a significantly substructured system, either from optical observation of galaxy subclumps (e.g., ^{Geller & Beers 1982}; ^{Bravo-Alfaro et al. 2009}; ^{Caretta et al. 2023}) or the detection of multiple peaks in X-ray emission (e.g., ^{Jones & Forman 1984}; ^{Buote & Tsai 1995}; ^{Lagan´a et al 2019}), cannot be considered in dynamical equilibrium. Thus, the presence and significance of substructures reveal how far the galaxy system is from a relaxed and homogeneous global potential. A description of the cluster level of internal substructuring is called its gravitational assembly state (^{Caretta et al. 2023}).

In principle, one can also measure global parameters related to the internal dynamics (^{Carlberg et al. 1996}; ^{Girardi & Mezzetti 2001}) -or dynamical state- of the galaxy system. These parameters are associated, for example, with the mass, radius and velocity dispersion of the system. Such approach supposes that both the observable aspect and structure and the dynamical parameters of the system change in a correlated manner during its evolution, dominated by mechanisms that take it from more irregular and substructured configurations to those with more homogeneous galaxy distributions and more dynamically relaxed (^{Araya-Melo et al. 2009}).

Furthermore, from X-ray observations one can construct entropy (or temperature or density) profiles that account for the evolutionary history, structure and thermodynamic state of the ICM inside the clusters (^{Tozzi & Norman 2001}; ^{Voit 2005}), which also contributes to the study of the degree of relaxation -or disturbance- of their gravitational potentials (analysis of self-similarity of clusters). Nevertheless, we are not always fortunate enough to detect X-ray emission from clusters nor to have the amount of data necessary to carry out massive studies, so for this we are still limited to optical surveys.

4.1 Observational Data

We use data from [^{Caretta et al.2023}], a sample of 67 galaxy clusters, from Abell/ACO (^{Abell 1958}; Abell et al. 1989) catalogs, with redshifts up to z≈0.15. These clusters were selected because they are among the most well sampled galaxy systems in the nearby Universe, and cover roughly uniformly from poor to rich systems (with ICM-temperatures from 1 to 12 keV), being balanced for including all BM (^{Bautz & Morgan 1970}) types. Three other non-Abell clusters with similar characteristics were included in our sample (AM0227-334, SC1329-313 and MKW03S), which allows us to call (for short) our observational sample Top70 from now on (Table 1). For each cluster, a sample of spectroscopically confirmed member galaxies is available, selected inside the cluster caustics up to a fiducial aperture of 1.3×r200 from the center of the cluster (chosen to be the Central Dominant Galaxy, CDG). The numbers of these presumably virialized members range from 21 to 919 (average 154), while originally at least 100 spectroscopic redshifts were available for each cluster (^{Caretta et al. 2023} for the details of the process for determining membership, caustics and virial radius). One should note that the richness of a cluster is related to both intrinsic and observational conditions -more massive clusters are richer, while nearby clusters tend to be preferentially observed. However, the numbers of members we have are adequate for minimizing observational biases by counting them in bins. Astrometric positions of galaxies have uncertainties of ±0.25'', and radial velocities of ±60 km s^-1 (see, ^{Caretta et al. 2023}, and references therein).

By using different 1D, 2D and 3D methods (e. g., ^{Dressler & Shectman 1988}), applied to the distribution of the members inside the caustics, these authors searched for optical substructures in the cluster sample, supplementing their analysis with X-rays and radio literature data. They found that at least 70% of the clusters in their sample present clear signs of substructuring, with 57% being significantly substructured. The significance of the identified substructures in each cluster was estimated by the fraction of galaxies they contain with respect to the total richness.

The clusters were classified into five assembly state levels according to the presence -or not- and relative importance of substructures: unimodal systems, made up of a regular structure (U); low mass unimodal systems (L); systems with a primary structure and only low significance substructures (P); significantly substructured systems with one main substructure (S); and multimodal conglomerates with more than one main substructure (M). While the L clusters are young poor galaxy systems, representing evolutionary states prior to amalgamation processes but already relatively placidly evolved, the U clusters are old and massive, which have probably grown by mergers and accretions and have already settled close to a relaxation state. P clusters are also old and massive, but still present signs of recent accretions; since these accretions are minor, the relaxation state of the cluster is almost unaffected. Finally, S and M clusters are systems during merging processes, the difference being whether these mergers are minor or major, respectively.

The benefit of using this sample lies in the availability of this discrete classification of gravitational assembly states (see Column 2 of Table 2 below). This will serve to establish correlations between the H _Z estimator applied to the ensemble of galaxies of each cluster and the evolutionary state obtained from direct observational methods.

The observational cluster sample is reported in Table 1: Column 1 shows the name of clusters; Columns 2, 3 and 4 present their right ascension, declination and mean redshift, respectively; Column 5 presents the number of presumably virialized sampled galaxies belonging to the system (inside the caustics and the aperture 1.3×r200); Column 6 presents the radius r ₅₀₀ at which the mean interior overdensity is 500 times the critical density at the corresponding redshift; and Column 7 presents the X-ray temperature of the clusters.

For each observed cluster, the virial mass was estimated using (3) and the virial radius by

where ρvir=18π23H2(z)/8πG is the virialization density assuming a spherical model for nonlinear collapse (e.g., ^{Bryan & Norman 1998}). In all calculations we used α=2.5, assuming a weak anisotropy (e.g., ^{Tully 2015}; ^{Kashibadze et al. 2020}). The line-of-sight velocity dispersions were computed using the Tukey’s biweight robust estimator (^{Beers et al. 1990}). The results of σLOS, Mvir and R _vir are shown, respectively, in Columns 8, 9 and 10 of Table 1.

4.2 Simulation Data

In addition, we use data from the IllustrisTNG simulation (^{Springel et al. 2018}; ^{Nelson et al. 2019}), a set of cosmological simulations that assumes initial conditions consistent with the [^{Collaboration 2016}] results and takes into account magnetohydrodynamical effects. In particular, we use the TNG300-1 data cube that contemplates the largest volume allowing the study of the distribution of galaxies and massive objects such as clusters. We take halos (equivalent to galaxy clusters) from the TNG300-1 simulation at redshift z = 0,² which has a resolution of 2500³ dark matter particles and 2500³ baryonic matter particles. We sampled 248 halos with masses greater than 1014h-1M⊙ with 370,119 member subhalos, of which 366,940 are classified as galaxies. Information related to the peculiar velocity, mass and position of the center of each halo and subhalo was extracted.

Each halo in the simulation contains, on average, 1,300 subhalos within its virial radius, greatly outnumbering the tracers (galaxies) in our observational sample of clusters. This large difference is due to the number of low-mass systems (dwarf galaxies subhalos) entering the count of TNG300-1; these galaxies are hard to detect in real clusters due to their low luminosity. Thus, to avoid statistical differences between the parameters estimated from observational and simulated data, we limit the selection of member subhalos by taking only those with mass greater than 2.0×1010h-1M⊙. For each halo, all member subhalos up to a distance of 1.3×r200 from the center (the particle with the least gravitational potential energy) were taken. Also, σLOS, Mvir and R _vir were calculated reproducing the observational procedure.

Figure 1 shows the projected distribution of member subhalos for two halos in the sample, one high-entropy (left) and one low-entropy (right). Note that, as expected, the distribution of subhalos is more random and homogeneous in a high entropy halo, while more substructured and elongated in a low entropy halo. This also happens in the observed clusters, reinforcing our hypothesis that evolutionary changes in galaxy systems progress in the direction in which the system dissolves substructures, becoming observationally more regular and homogeneous, with higher entropy values.

Fig. 1 Two examples of sampled TNG300 halos. Left: subhalo distribution for a high entropy halo (TNG-halo-34). Right: subhalo distribution for a low entropy halo (TNG-halo-87). Each dot represents a member subhalo: the small dots are subhalos with masses less than 2.0×1010M⊙ while the big dots are the subhalos taken for our analysis. 

4.3 Results on the Assembly State of Observed Clusters

To appreciate the correlation between the H _Z entropy (shown in Column 3 of Table 2 below) and the gravitational assembly level -and therefore with the evolutionary state- of galaxy systems, we present in Figure 2 the distributions of the H _Z values with respect to the assembly state classes. These are displayed in the form of boxplots, which allow graphically describing the locality, dispersion and asymmetry of data classes (groups) of a quantitative variable through its quartiles (^{Heumann & Shalabh 2016}). The classes were ordered from more to less relaxed systems. Although the box overlap does not allow one to unambiguously discriminate the class to which an arbitrary individual cluster should belong, the statistical correlation is clear: the sequence U-P-S-M-L follows directly the decreasing median values of entropy. The only marginal difference case is between U and P systems, which is understandable if we consider that the primary systems (with only low significance substructures) are dynamically very similar to the unimodal ones -both tend to be relatively massive and evolved systems. This discussion will be resumed later.

Fig. 2. Boxplots of H _Z -entropy for the five assembly classifications of clusters. The lower and upper extremes of each box are the 25th and 75th percentiles, respectively, while the central red line marks the median. Whiskers extend to the most extreme non-outlier data, while outliers are represented by red ‘+’ symbols. The color figure can be viewed online. 

5.1 Shannon Entropy of Galaxy Distributions in Phase-Spaces

Another way to evaluate the entropy estimates that H _Z can provide is by using a method that does not take into account equilibrium assumptions, but allows us to characterize the internal states of a galaxy system from the raw distribution of its member galaxies. In information theory, for example, entropy is a measure of the uncertainty of a random variable (or source of information, e.g., ^{Shannon 1948}; ^{Cover & Thomas 2006}). If X is a discrete random variable with possible observable values x∈X, which occur with probability p(x)=PrX=x, the Shannon (or information) entropy of X is defined as

where the sum is performed over all possible values of the variable, and the base λ of the logarithm is chosen according to the entropy ‘units’ to use (nats, λ=e; bits, λ=2; bans, λ=10). One of the criteria used by ^{Shannon (1948)} to define H _S ensures that it increases as the possible values of X begin to appear with equal frequency, taking higher values when they become equiprobable, i.e., when there are no special configurations in the data distribution that provide more information, increasing the uncertainty.

Now, the raw coordinates of the galaxies in a cluster, i.e. the observable set of triples (RA,Dec,z) of right ascension, declination and redshift, are distributed within a solid angle, which can be approximated by a cylinder with a circular base in the plane of sky and depth along the line-of-sight (see, Figure 3, left panel). Inside the cylinder, the position of each galaxy can be expressed in the form x=(r,θ,z), where r is its projected distance from the cluster center, θ its -azimuthal- angle with respect to the local north direction in the projected sky distribution (see, Figure 3, right panel), and the redshift z is a measure of its radial velocity. In fact, the cylinder of x-coordinates (two spatial and one velocity) may be considered a -projected- phase-space for the galaxy ensemble, and the distribution of the variables (r,θ,z) inside it depends on the dynamical state of the galaxy cluster.

Fig. 3 Sample of galaxies for A2199 cluster in ^{Caretta et al. (2023)}. Color bars represent the redshift (z) distribution. Left panel: the cylinder along the line-of-sight direction. Right panel: projected distribution in the sky plane with RA and Dec coordinates transported to the origin (the FRG, see the text). The color figure can be viewed online. 

By considering the cylinder of sampled galaxies as a source of information for the x-distribution, the probability of finding a galaxy in the neighborhood of position x (i.e., the probability that the position random variable X takes a value very close to x inside the cylinder) can be approximated as p(x)≃f-rθz(x)Δx, where f-rθz must be the observed -or empirical- joint probability density function (PDF) that represents the actual distribution of the galaxy ensemble in the variable x. Replacing this in (14), one can compute the Shannon entropy of the galaxy distribution for a cluster in the form

where λ=e has been chosen, and the sum is performed over all the n discrete partitions Δx made to the domain X=0,Rvir×0,360×zmin,zmax of the observed joint PDF. Here, z _min and z _max are the minimum and maximum values of redshift in the sample of galaxies of each cluster.

Shannon entropy is a well-established measure for quantifying disorder or uncertainty in a system (^{Cover & Thomas 2006}). In this context, H _S provides a measure of the degree of randomness in the distribution of galaxies in the phase-space of a cluster (or any other galaxy system). By calculating this entropy, we are evaluating the information contained in the galaxy ensemble and how galaxies are distributed in different regions of phase-space. Close to equilibrium, the memory (information) of the initial conditions of cluster formation is lost (^{Binney & Tremaine 2008}; ^{Araya-Melo et al. 2009}). The general character of H _S comes from the fact that it is not restricted only to thermodynamic variables, but to any type of data X that contains information about the state of a system, so it can be used to study characteristics of non-equilibrium systems as they evolve. In addition, it does not consider the microstates of a system as equiprobable, keeping a certain relationship in mathematical form and meaning with the Gibbs entropy (^{Jaynes 1957}). .

5.2 Shannon Entropy Estimations

To calculate the Shannon entropy, we first established the observed joint PDF f-rθz(x) of each cluster in two different ways. First, by counting galaxies independently in bins of width Δr, Δθ and Δz in the radial-r, azimuthal-θ and redshift-z directions, respectively, we constructed 1D-histograms that describe the observed distribution of that variable in the galaxy ensemble. Then, by a smoothing technique, the observed PDFsf-r(r), f-θ(θ) and f-z(z) were obtained for the variables from their respective normalized histograms (see, red solid lines on graphs in Figure 3). For this, we used a kernel density estimator that allows a non-parametric fit of PDFs of random variables, adapting directly to the data (e.g., ^{Heumann & Shalabh 2016}). This method is particularly useful when the actual distribution of a data set is unknown, as is the case for the galaxy distributions inside the (phase-space) cylinders. A standard Gaussian smoothing kernel was used with bins of widths Δr=0.15 Mpc, Δθ=12∘, cΔz=200 km s^-1 and supports in the intervals [R _vir], [0,360] and zmin,zmax for the f-r, f-θ and f-z distributions, respectively. Finally, we took f-rθz(r,θ,z)=f-r(r)f-θ(θ)f-z(z) assuming statistical independence of the variables r, θ and z in the galaxy distributions.

On the other hand, as is evident in the right panel of Figure 3, galaxies with the highest and lowest redshifts prefer the center of the projected distribution of the cluster. This correlation between the variables r and z is physically justified since the galaxies acquire a greater speed during their transit through the central regions of the clusters (i.e., where the gravitational potential well is the deepest). Projection effects can also occur, especially for halo galaxies that are located in front or behind the core and close to the line-of-sight. No significant correlation was detected between the variables r or z with θ instead. Thus, to construct the observed joint PDF in the second way we only took into account the correlation between r and z, such that f-'rθz(r,θ,z)=f-rz(r,z)f-θ(θ), where the f-rz functions were obtained by counting galaxies in bi-dimensional bins of “area” ΔrΔz and smoothed by a multivariate Gaussian surface kernel, using the same bin width values as above.

Two statistical tests, Kolmogorov-Smirnov and Rank-Sum, showed no significant differences between both methods, confirming the null hypothesis that the PDFs f-rθz and f-'rθz represent the same distribution with a confidence level of 95% in both tests. This is probably because the fraction of galaxies that present a strong rz-correlation is very small, implying a low statistical weight. Then, for simplicity we choose the first option (i.e., that with the three variables assumed to be statistically independent, see Figure 4) to compute the probability values

varying one of the variables in its respective domain (e.g., r∈0,Rvir, θ∈0,360 and z∈zmin,zmax), while keeping the other two constant. Thus, p (x) was obtained in about 1,000 positions x=(r,θ,z) inside the data cylinder of each real (and simulated) cluster and these values were used to compute the respective H _S -entropies by (15) for the galaxy ensembles. The results for Top70 clusters are shown in Column 4 of Table 2.

Fig. 4 Probability density functions for the radial-r, azimuthal-θ and redshift-z distributions of member galaxies in the A2199 cluster. The red solid lines represent the observed PDFs obtained by smoothing (normal kernel) from observational data histograms, while dashed blue lines represent the relaxed PDFs of equilibrium distributions. The error bars represent the difference between the height of the smoothed curve and the height of the corresponding bin in each 1D-histogram. The color figure can be viewed online. 

5.3 Results on the Shannon Entropy for Observed and Simulated Clusters

The relation between the H _Z -entropy estimator, proposed in the present work, and the Shannon entropy H _S of the real galaxy distributions can be seen in Figure 5, where the dashed line represents the best fit curve, with a coefficient of determination (R ², ^{Heumann & Shalabh 2016}) of 0.886. This figure shows a high degree of correlation (almost linear) between the H _Z and H _S entropies that, when measured by the Pearson and Spearman coefficients (e.g., ^{Gibbons & Chakraborti 2003}; ^{Heumann & Shalabh 2016}), gives the values 0.932 and 0.922, respectively. Below (see, §7), we offer what we believe to be a possible explanation for such an explicit correlation. In addition, the figure also shows the association of these H _Z and H _S entropy values with the assembly state of the clusters (represented by the U-P-S-M-L classes). Although the correlation is not evident, it is noticeable that all U and P clusters are placed in the locus of points with higher values of both H _Z and H _S , while all L are located in the region with lower values of these entropies.

Fig. 5 Scatter plot of H _Z vs. H _S made with the entropy values estimated for the Top70 clusters. The symbols and color scale of the points represent the U-P-S-M-L assembly level classification of clusters performed by ^{Caretta et al. (2023)}. The dashed line represents the best -power law- fit with R2=0.886. The color figure can be viewed online. 

A procedure similar to that performed with Top70 clusters was applied to the TNG300 sample to compare the H _Z and H _S entropies estimated for the simulated halos. Figure 6 shows a significant correlation, with Pearson and Spearman coefficients of 0.709 and 0.678, respectively, between the dynamical and Shannon entropies, which are very similar to those of real clusters. We did not carry out a classification of simulated clusters in their different assembly levels (A) that allows us to analyze the distribution of H _Z in each class, as in the case of real clusters. We hope to do this in future work.

Fig. 6 Scatter plot of H _S vs. H _Z entropies made with the data of the complete TNG300 sample. The dashed line represents the best power law fit with R2=0.505. 

6.1 The Relaxation Probability of Galaxy Systems

Several studies reveal that clusters close to dynamical equilibrium have distributions of galaxies whose radial density and LOS velocity profiles tend to ones that can be represented by specific mathematical functions (^{Saslaw & Hamilton 1984}; ^{Sarazin 1988}; ^{Adami et al. 1998}; ^{Sampaio & Ribeiro 2014}). The above is also true for the ICM entropy profiles in X-ray observations (e. g., ^{Voit 2005}) or for the distribution of DM subhalos in cosmological simulations (^{Davis et al. 1985}; ^{Springel et al. 2005}). Thus, we consider that it is also possible to characterize the evolutionary state of a galaxy system by measuring how far the currently observed PDFs are from their expected equilibrium functional shapes. For this, we need to use (or choose) a consistent reference equilibrium model that describes the distribution when the galaxy ensemble is relaxed, as well as a metric that allows us to estimate the distance of the observed distribution from that of the equilibrium model.

We limit ourselves to a simple reference equilibrium model, considering clusters represented by a spherical distribution of galaxies, with a homogeneous core-halo spatial configuration and an isotropic velocity distribution without net angular momentum.

In principle, the radial distribution of galaxies can be well described [?, see]]ad1998 by different single mass profiles proposed for clusters in the literature, both core- (^{King 1962}; ^{Einasto 1965}) and cuspy- (^{Hernquist 1990}; ^{Navarro, Frenk & White 1996}) dominated ones. Nonetheless, there has been a tendency, especially for fitting DM halos, to use more complex functions, with a larger number of free parameters (^{Dehnen 1993}; ^{Fielder et al. 2020}; ^{Diemer 2023}) -although they are essentially double and/or truncated power laws- in order to better accommodate the inner and outer slopes. On the other hand, since core-dominated profiles usually perform slightly better for real clusters (^{Sarazin 1988}; ^{Adami et al. 1998}; ^{García-Manzanárez 2022}) -as commented before, the core-halo structure is the one expected for self-gravitating systems at equilibrium- and considering the simplest possible analytical form, we choose the King density profile at this point. Such a profile can be approximated analytically by the form

suggesting a finite central density ρ0=ρ(0) and the existence of a core of radius r _c . The parameters r _c , ρ0 and γ can be determined by the best fit of the model (17) to the spatial distribution of galaxies. For three-dimensional distributions it has been found that γ=3/2 and ρ0=9σLOS2/4πGrc2, while for observed two-dimensional (or projected) distributions γ=1 and ρ0=9σLOS2/2πGrc c (see, ^{Rood et al. 1972}; ^{Schneider 2015}), where G is the gravitational constant. Thus, here we consider the projection (on the RA-Dec plane) of a cluster in equilibrium to have a King-type radial distribution of galaxies described by the (1D) radial-PDF of the form

in which the count of galaxies must be performed in bins of length Δr instead of rings of area ΔA=2πrΔr, that is, counting the number of galaxies located between r and r+Δr for each distance from the cluster center. Integrating (18) from 0 to R _vir, the normalization condition to find the ‘relaxed’ core radius (r’ _c ) in virial equilibrium is

For the azimuthal distribution of galaxies on the RA-Dec projected plane of a regular cluster we expect a continuous uniform azimuthal-PDF of the form

since, under the assumed equilibrium model, the probability of finding galaxies in any direction of the plane must be the same if there are no deformities (flattening or elongation) in the cluster morphology and no substructures when counting galaxies in slices of width Δθ at fixed radius.

Finally, for the (3D) galaxy velocities inside clusters in equilibrium we expect a quasi-Maxwellian distribution n (isotropic, e.g., ^{Sarazin 1988}; ^{Sampaio & Ribeiro 2014}), so that the LOS components of velocities have a normal distribution described by the redshift-PDF

where c is the speed of light and cz- is the mean LOS velocity of the cluster.

The construction of the reference PDFs associated to the data is done in three steps. First, we determine numerically the relaxed core radius r’ _c (see, Column 5 of Table 2 for observational sample), taking into account the normalization condition (19), to be used in freq. For R _vir, also used in freq, and for cz- and σLOS, the last ones for fzeq, we take the previously calculated parameters (see, Columns 10, 4 and 8, respectively, of Table 1 for observational sample). The fθeq distribution is trivial and does not require observational parameters of clusters.

Next, for each cluster we construct its corresponding relaxed mock cluster, which has as many particles as observed galaxies but distributed according to (or following) the corresponding freq, fθeq and fzeq equilibrium PDFs.

The third step is done by fitting the tuned equilibrium models to the particle distributions of the mock clusters in the radial, azimuthal, and redshift components, using bin widths and smoothing levels equal to those used in the observed PDFs. With this, we are able to construct the relaxed PDFs, f~r, f~θ and f~z, i.e., the ones expected in the dynamical relaxation state. Figure 3 shows an example of the observed and relaxed PDFs for the A2199 cluster.

Now, the relaxation probability of a galaxy system can be defined as a distance between the observed PDFs and the relaxed ones, i.e., between the current dynamical state of the galaxy ensemble (characterized by f-r, f-z and f-θ) and its most probable equilibrium state (characterized by f~r, f~θ and f~z), in the probability space. For this, we use the ^{Hellinger distance (Hellinger 1909)}, a metric used to quantify the similarity between two distributions in the same probability space so that, if f1 and f2 represent two PDFs for the same variable, the Hellinger distance between them is defined as

where the integration must be carried out over the domain of the functions, and the property 0≤H(f1,f2)≤1 allows us to define the relaxation probability, Prelax≡1-H, of a galaxy ensemble. Thus, the closer the f-i and f~i functions are to each other, the more relaxed a galaxy system can be considered. For a system close to equilibrium H⟶0 and Prelax⟶1.

Under the same assumption of statistical independence between the variables r, θ and z used for the observed joint PDFs, the relaxed joint PDFs were then defined as f~rθz(r,θ,z)=f~r(r)f~θ(θ)f~z(z). For calculations of H, we take f1=f-rθz and f2=f~rθz in (22) for each galaxy cluster. The calculated values for the relaxation probability, Prelax, for the observed clusters are shown in Column 6 of Table 2.

In addition, the top panel of Figure 6 shows the distribution of Prelax values with respect to the assembly state classes. Like what happened with the H _Z -entropy, there is a clear correlation between the relaxation probability and the level of gravitational assembly of a system: the sequence U-P-S-M-L considered here goes from the most likely relaxed to the least relaxed systems with respect to the chosen equilibrium model. The triple correlation between the dynamical entropy, the relaxation probability and the assembly level of clusters can be seen in the bottom panel of Figure 7, a scatter plot of Prelaxvs.HZ built with the mean and median values (see Table 3) obtained for H _Z and Prelax in each assembly classification.

Fig. 7 Top panel: Boxplots of Prelax values for the five assembly states of clusters from ^{Caretta et al. (2023)}. Bottom panel: Scatter plot of Prelaxvs.HZ made with the mean (black markers) and median (blue markers) H _Z and Prelax values presented in Table 3. The symbols of the points represents the U-P-S-M-L assembly level classes. The black dashed and blue solid lines represents the linear fits to the mean and median data, with R2 of 0.964 and 0.950, respectively. The color figure can be viewed online. 

Finally, we applied a Kolmogorov-Smirnov (KS-) test to evaluate, with a confidence level set at 90% (or significance level of 0.1), the null hypothesis that the H _Z or Prelax values for the galaxy systems classified into two different assembly classes (U, P, S, M and L) come from the same continuous distribution. Table 4 shows the p-values, in the range [0,1], that resulted from the KS-test. The closer a p-value is to 1, the more similar are the distributions of H _Z and Prelax values in two different assembly (A) classes.

6.2 Relation Between H _Z and the Cluster Concentration Indices

Apart from the discrete characterization of the assembly state used in the previous sections, one can also probe some continuous parameters, estimated from optical or X-ray data (^{Carlberg et al. 1996}; ^{Girardi & Mezzetti 2001}; ^{Zhang et al. 2011}; ^{Caretta et al. 2023}, and references therein), that are expected to correlate with the evolutionary state of the galaxy system. Here we present, for instance, the concentration indices of a cluster, both from optical galaxy distributions and X-ray from ICM: more relaxed clusters tend to be denser at their centers and have higher values for these indices.

For the optical data, we use the Maximum Likelihood Estimation method (MLE) to determine the optimal concentration indices for the King and NFW radial profiles fitted to the observed projected distributions of galaxies in the clusters of the Top70 sample. The King profile employed for the fitting is the one described by equation (17), while the NFW profile takes the form

where r _s is the scale radius of the galaxy cluster. The concentration indices of each cluster are related to its respective virial and characteristic radii (r _c or r _s obtained by MLE) in the form cK=Rvir/rc and cNFW=Rvir/rs (^{Binney & Tremaine 2008}), and the obtained values are compiled in Columns 7 and 8 of Table 2. We also probed the ratio c500=Rvir/r500 (using the available X-ray data in Table 1); the obtained values are shown in Column 9 of the Table 2.

Figure 8 shows the relationship between the H _Z -entropy and the estimated -cluster- concentration indices for the Top70 sample. The Pearson correlation coefficients between the H _Z -entropy values and those corresponding to the indices c _K, c _NFW and c ₅₀₀ are 0.448, 0.512, 0.429, respectively. Despite the not strong statistical significance and the remarkable dispersion, an evident growing trend can be seen in the concentration indices with the growth of the H _Z -entropy in the systems. This tendency is also subtly manifested in the assemblage classes (represented by the marker symbols in Figure 8): the U and P clusters tend to have higher concentration indices, while the L clusters tend to present lower values of these. S and M clusters appear with a more extended range of concentration indices. It is interesting to note that if we leave only U and P clusters in the above relation, for example for c ₅₀₀, the correlation increases from 0.429 to 0.603 -these are the typical clusters that appear in studies based on X-rays emission of ICM.

Fig. 8 Scatter plots of H _Z -entropy vs. concentration indices (c _K, c _NFW, and c ₅₀₀) for the Top70 cluster sample. The symbols used for the markers have the same meaning as in Figure 5. The dashed line represents the best linear fit to the data, with R2 of 0.203, 0.264 and 0.215, respectively from top to bottom. The color figure can be viewed online. 

We have also probed other continuous dynamical parameters, obtained from X-rays emission of ICM, associated to the evolutionary state of galaxy systems: the concentration c and the luminosity concentration c _L (resepctively from ^{Parekh et al. 2015}; ^{Zhang et al. 2017}), the Gini and M ₂₀ morphological coefficients (from ^{Parekh et al. 2015}; ^{Lovisari et al. 2017}), the substructure level S _C (from ^{Andrade-Santos et al. 2012}), and the C _SB and C _SB4 concentration parameters (from ^{Andrade-Santos et al. 2017}). Although the statistics are usually poor, and the aspects captured from X-ray data are not necessarily similar to the ones captured by optical data (which include H _Z ), the results are all similar to the ones showed for c500, for example increasing correlation when only U and P clusters are considered.

It is remarkable that the only parameter that correlates very well with H _Z is r’ _c (with a Pearson coefficient of 0.972). The theory states that the natural tendency of a gravitational system towards its evolution is establishing a structure core-halo. This is exactly what the correlation of H _Z and r’ _c may be suggesting: the core radius grows with the entropy of the cluster.

We consider that H _Z is an efficient parameter for estimating the relaxation/dynamical state of galaxy clusters, maybe better than the others we have compared in this incomplete exercise, because it can capture more important aspects of such an evolutionary snapshot. However, it is important to note that no parameter can, alone, account for all the evolutionary aspects, such as interaction with the cosmological environment, collapse and mass growth, galaxy formation, AGN activity, feedback processes, among others, making paramount to consider different observations and analyses for constructing the most complete picture possible.