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

# Galaxy clusters in simulations of the local Universe: a matter of constraints

Jenny G. Sorce<fnref n="1" />,<fnref n="2" />,<fnref n="3" />* and Elmo Tempel<fnref n="3" />,<fnref n="4" />

Accepted 2018 February 20. Received 2018 February 12; in original form 2017 November 24

## ABSTRACT

To study the full formation and evolution history of galaxy clusters and their population, high-resolution simulations of the latter are flourishing. However, comparing observed clusters to the simulated ones on a one-to-one basis to refine the models and theories down to the details is non-trivial. The large variety of clusters limits the comparisons between observed and numerical clusters. Simulations resembling the local Universe down to the cluster scales permit pushing the limit. Simulated and observed clusters can be matched on a one-to-one basis for direct comparisons provided that clusters are well reproduced besides being in the proper large-scale environment. Comparing random and local Universe-like simulations obtained with differently grouped observational catalogues of peculiar velocities, this paper shows that the grouping scheme used to remove non-linear motions in the catalogues that constrain the simulations affects the quality of the numerical clusters. With a *less aggressive grouping scheme – galaxies still falling on to clusters are preserved – combined with a bias minimization scheme*, the mass of the dark matter haloes, simulacra for five local clusters – Virgo, Centaurus, Coma, Hydra, and Perseus – is increased by 39 per cent closing the gap with observational mass estimates. Simulacra are found on average in 89 per cent of the simulations, an increase of 5 per cent with respect to the previous grouping scheme. The only exception is Perseus. Since the Perseus-Pisces region is not well covered by the used peculiar velocity catalogue, the latest release lets us foresee a better simulacrum for Perseus in a near future.

**Key words:** methods: numerical – techniques: radial velocities – galaxies: clusters: individual – galaxies: groups: general – large-scale structure of Universe.

## 1 INTRODUCTION

Clusters of galaxies are excellent cosmological probes. Understanding their formation and evolution is thus an entirely logical step in our quest towards understanding the Universe as a whole. However, accessing detailed information about galaxy clusters via observations is far from direct and the extracted information might suffer from observational biases. To complement observational studies, high-resolution simulations of galaxies clusters are now flourishing (e.g. see Wu et al. 2013a,b, 2015; Cui et al. 2016; Elahi et al. 2016; Martizzi et al. 2016; Sembolini et al. 2016a,b; Arthur et al. 2017; Baldi et al. 2017; Barnes et al. 2017; Hahn et al. 2017, for a non-extensive list) and comparisons between observed and simulated clusters emerge (e.g. Jimeno et al. 2017; Natarajan et al. 2017; Suto et al. 2017). Still these comparisons are limited because

of the large variety of cluster types in terms of morphology, mass, evolution stage, etc. (Struble & Rood 1988). Selecting adequately the simulated cluster candidates to be compared with a given observed cluster is not immediately obvious and in a certain sense it is almost impossible to find the exact counterpart (Grossauer et al. 2015). Determining how effectively a numerical cluster represents an observed one is actually subject to uncertainties.

One way to reduce these uncertainties is to use simulations that resemble a portion of the Universe with well-observed clusters. Such simulations host clusters similar to the observed ones in the proper environment and thus make the comparisons between observations and simulations even more legitimate. The most well-observed part of the Universe is undeniably the local Universe. Consequently, efficient simulations of the local Universe down to the cluster scales constitute the optimal choice to perform the detailed comparisons between observations and simulations mentioned earlier. Such simulations unlike typical ones stem from a set of constraints in addition to abiding to a cosmological prior (Bertschinger 1987; Gottlöber,

<footnote n="1">¹University Lyon, University Lyon1, Ens de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon UMR5574, F-69230 Saint-Genis-Laval, France</footnote>

<footnote n="2">²Université de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, F-67000 Strasbourg, France</footnote>

<footnote n="3">³Leibniz-Institut für Astrophysik Postdam (AIP), An der Sternwarte 16, D-14482 Potsdam, Germany</footnote>

<footnote n="4">⁴Tartu Observatory, University of Tartu, Observatoorium 1, 61602 Tõravere, Estonia</footnote>

<footnote>* E-mail: jenny.sorce@univ-lyon1.fr</footnote>
Hoffman & Yepes 2010; Lavaux 2010; Kitaura 2013). These con-
straints can be either peculiar velocities (e.g. Klypin et al. 2003),
equivalently distances (e.g. Lavaux 2016) or redshift surveys (e.g.
Bertschinger & Dekel 1989; Bertschinger et al. 1990; Kitaura et al.
2009; Heß, Kitaura & Gottlöber 2013; Wang et al. 2016). The ini-
tial conditions constrained by the aforementioned measurements
can be produced either forward (e.g. Kitaura & Enßlin 2008; Jasche
& Wandelt 2013; Wang et al. 2013, 2014) or backwards (e.g. Dekel,
Bertschinger & Faber 1990; Ganon & Hoffman 1993; Zaroubi,
Hoffman & Dekel 1999; Lavaux et al. 2008). We use the latter in this
paper and in our previous papers. Our first simulation resembling
the local Universe, obtained successively with the first catalogue
(Tully et al. 2008) of peculiar velocities of the Cosmicflows project
and the second one (Tully et al. 2013), hosted the local large-scale
structure with a remarkable accuracy (Sorce et al. 2014, 2016b).
In addition, large overdensities at the location of prominent local
clusters clearly appeared in these simulations. However, the appli-
cation of a halo finder to these simulations revealed that clusters,
our closest neighbour Virgo cluster excluded (Sorce et al. 2016a),
are not strongly reproduced in the simulation: the largest object in
the large overdensity is not massive enough with respect to expecta-
tions based on observational estimates. In a companion paper (Sorce
& Tempel 2017), we showed that the grouping technique applied
upstream on the catalogue of constraints used to build the initial
conditions might affect the production of massive clusters although
it does not affect overall the local large-scale structure. Namely, the
overdense regions hosting the clusters are more or less pronounced
depending on the grouping scheme used.

The grouping applied to the constraints is an absolute require-
ment. Indeed, the technique we use to build the initial conditions
is linear but constraints include galaxies and their velocities in
various environments including dense environments like clusters.
Therefore, non-linear motions find their way into the catalogue of
constraints and are passed along to the reconstruction technique
that cannot handle them in an appropriate manner. This same phe-
nomenon is visible in redshift surveys in the form of fingers of
god (Jackson 1972) and Kaiser's effect (Kaiser 1987). These effects
need to be suppressed since they affect reconstructions based on
redshift surveys (Kitaura et al. 2012). A different modelling can be
applied to the small scales with respect to the large scales: Spherical
collapse versus Lagrangian perturbation (Heß et al. 2013). However,
a certain balance is necessary: while it is necessary to group to sup-
press non-linear motions, galaxies in the field and a fortiori galaxies
infalling on to clusters are essential to retrieve the proper density
field and obtain an optimal reconstruction of the local Universe
(Sorce, Hoffman & Gottlöber 2017). In this paper, we show that as
expected what is true for the reconstruction of the local Universe
(Sorce & Tempel 2017) is also valid for its simulations.

This paper starts with a brief description of the catalogue of
constraints, the grouping algorithm applied to it, and the different
steps to build the constrained initial conditions. Then, the resulting
simulations of the local Universe are analysed and compared to
those obtained with the earlier released version of the grouped
catalogue as well as to random simulations. Finally, a conclusion
closes the paper.

## 2 BUILDING CONSTRAINED INITIAL CONDITIONS

The different steps to produce constrained initial conditions used in
the project have been widely described and summarized in previous

papers (e.g. Sorce et al. 2016b). In the following, they are briefly
reminded.

2.1 The catalogue

The second catalogue of radial peculiar velocities or more pre-
cisely of direct distance measurements of the Cosmicflows project
constitutes our set of constraints. Published in Tully et al. (2013),
it contains more than 8000 accurate galaxy distances mostly
(~88 per cent) obtained with the Tully–Fisher (Tully & Fisher 1977)
and the Fundamental Plane (Colless et al. 2001) methods. Cepheids
(Freedman et al. 2001), Tip of the Red Giant Branch (Lee, Freed-
man & Madore 1993), Surface Brightness Fluctuation (Tonry et al.
2001), supernovae of type Ia (Jha, Riess & Kirshner 2007), and other
miscellaneous methods constitute the remaining ~12 per cent. It ex-
tends up to about 250 $h^{-1}$ Mpc and about 50 per cent of the data are
within 70 $h^{-1}$ Mpc and 90 per cent within 160 $h^{-1}$ Mpc.

2.2 The grouping scheme

The grouping scheme is widely described in Tempel et al. (2016) and
our application to the catalogue of constraints is detailed in Sorce
& Tempel (2017). A brief description is given here as a reminder.

Tempel et al. (2016) introduced a new grouping method (here-
after, Tempel grouping scheme). This method is based on a widely
used Friends-of-Friends (FoF) percolation method, where different
linking lengths in radial (along the line of sight) and in transversal (in
the plane of the sky) directions are used but the conventional FoF
groups are refined using multimodality analysis. More precisely,
Tempel et al. (2016) use a model-based clustering analysis to check
the multimodality of groups found by the FoF algorithm and they
separate nearby/merging systems. In Sorce & Tempel (2017), we
tested different linking lengths and settled for the default one (0.25
$h^{-1}$ Mpc at redshift zero) so as to group sufficiently to remove non-
linear motions without large residuals and not too much so as to
preserve the infall on to the clusters.

The grouping scheme thus provides the groups to which the dif-
ferent galaxies that populate the second catalogue of Cosmicflows
belong to as well as their total velocity. This information is com-
bined with the galaxy distance estimates given by the second cat-
alogue of Cosmicflows to access galaxy radial peculiar velocities
(the constraints).

Furthermore, the constrained simulations obtained with this
grouping scheme are to be compared with the first generation of
constrained simulations obtained with the second Cosmicflows cat-
alogue of radial peculiar velocities and the grouping version (here-
after, Tully grouping scheme) released via the Extragalactic Dis-
tance Database<fnref n="1" /> (Tully et al. 2009). We remind that this earlier
scheme is based on literature groups and thus is not a systematic
scheme: within 30 Mpc, groups are those identified by Tully (1987),
further away groups are those given in the literature like Abell's cat-
alogue (Abell, Corwin & Olowin 1989).

## 2.3 Bias minimization, Reconstruction, Reverse Zel'dovich Approximation, Constrained Realizations, and Rescaling

Five more steps are required to complete the construction of the
constrained initial conditions:

<footnote n="1">¹ http://edd.ifa.hawaii.edu/</footnote>
(i) Minimization of the biases (Sorce 2015) inherent to any observational radial peculiar velocity catalogue. This minimization permits removing the spurious infall on to the local Volume and gives a proper Virgo cluster in the simulations and larger masses for the other nearby clusters.

(ii) Reconstruction of the cosmic displacement field with the Wiener-Filter (WF) technique (linear minimum variance estimator, in an abridged form WF; Zaroubi et al. 1995, 1999) applied to the peculiar velocity constraints.

(iii) Relocation of the constraints to the positions of their progenitors using the Reverse Zel'dovich Approximation and the reconstructed cosmic displacement field (Doumler et al. 2013) and replacing noisy radial peculiar velocities by their WF 3D reconstructions (Sorce et al. 2014) to ensure that structures are at the proper position at redshift zero.

(iv) Production of the density fields constrained by the modified observational peculiar velocities combined with a random realization to restore statistically the missing structures using the Constrained Realization technique (CR; Hoffman & Ribak 1991, 1992; van de Weygaert & Bertschinger 1996).

(v) Rescaling of the density fields to build constrained initial conditions and increasing the resolution by adding small-scale features (e.g. GINNUNGAGAP code<fnref n="2" />).

To enrich our comparisons, a set of random (typical) initial conditions is prepared. All the initial conditions are built within the Planck cosmology framework ($\Omega_m = 0.307$, $\Omega_\Lambda = 0.693$, $H_0 = 67.77$, $\sigma_8 = 0.829$; Planck Collaboration XVI et al. 2014) in 500 $h^{-1}$ Mpc boxes with 512³ particles (particle mass: $8 \times 10^{10} h^{-1}$ M$_\odot$). The simulations are run with the N-body code GADGET (Springel 2005).

# 3 LOCAL UNIVERSE-LIKE SIMULATIONS

In this section, nine simulations of each type (Tully grouping scheme, Tempel grouping scheme, and random) are used for further comparisons. A total of 27 simulations are thus run and every group of three simulations (one simulation per type in a given group) is based on the exact same random realization (cf. CR). For each type, nine simulations permit studying the stability of constrained simulations and the same random seeds for each group allow us to study the effect of the grouping technique. First the local large-scale structure is the object of focus before directing the efforts towards studying the clusters (dark matter haloes). For simplification, throughout the rest of the paper, 'Tully and Tempel constrained simulations' are adopted as a shortened notation.

## 3.1 The local large-scale structure

Fig. 1 shows the local large-scale structure obtained for two constrained simulations. The top panel presents two supergalactic slices of the local Universe obtained with Tully grouped version of the catalogue while the bottom panel gives the local structures resulting from Tempel grouping scheme. The solid contours stand for the overdensities while the dotted ones represent the underdensities. The green colour stands for the mean field. A few structures are identified with blue names. Overall, the local large-scale structure is well reproduced in both cases. It is very similar and the differences appear only in the details, at the cluster scale level. For instance, Coma appears more distinctly with Tempel grouping scheme than

with Tully's. In addition, if Virgo is well delimited in both cases, Centaurus appears more clearly for Tempel grouping scheme than for Tully's.

Sorce et al. (2016b) showed that the cosmic variance is reduced by a factor of 2–3 in the inner part of the box for the first generation of constrained simulations based on the second catalogue of Cosmicflows and Tully grouping scheme. It is interesting to quantify it when changing the grouping scheme for Tempel’s. Fig. 2 shows the average of the variances (filled circles) and their standard deviation (error bars) obtained when comparing pairs of random (R) and pairs of constrained ($C_{\text{Tully}}$, $C_{\text{Tempel}}$) simulations.

The variance is defined as the scatter around the 1:1 relation obtained when comparing every cell from a simulation to its exact counterpart in the other simulation of the pair. Once all variances are derived, their mean and standard deviation are derived and plotted as a filled circle with error bars. Since most of the constraints are within the inner part of the box, it is reasonable to compare not only the entire boxes but also their inner parts. Consequently, the process is repeated cutting the boxes to compare smaller and smaller regions.

The first panel of Fig. 2 shows the variance between pairs of simulations of the same nature. Clearly and as expected, the random simulations differ on average by ~32 per cent more from each other (black) than the constrained ones (blue). An identical reduction of the cosmic variance by a factor of 1.5 on average and 2.5–3 within the inner part of the box is observed for both grouping schemes. There is a limit to the method used to derive the cosmic variance. This is visible when comparing only the inner, and thus smaller, parts of the boxes: the mean variance decreases by 20 per cent for the pairs of random simulations. This is entirely due to the higher probability of finding small empty regions than large empty regions. Note however that since it is not improbable to find dense regions (even if the probability is low), the standard deviation is on average about 10 times larger when considering the inner part of the random box than when considering the entire random box. This inconvenience happens only for the pairs of random simulations. It is indeed well known that there are structures in the inner part of the box for the constrained simulations by construction (the local Universe has structures there).

The second and third panels of Fig. 2 show the variance between pairs of random and constrained simulations as well as between pairs of constrained simulations obtained with different grouping schemes. The third panel averages only the variances obtained for pairs of simulations sharing the same random realization, while the second panel averages the variances of the other pairs (not sharing the same random realization). Four points are worth noticing:

(i) The average variance obtained for pairs of constrained simulations based on different grouping schemes is smaller by about 20 per cent than that obtained for pairs of random and constrained simulations. In addition when comparing large volumes, the shape of the curve drawn by the variances is identical to that obtained when comparing constrained simulations obtained with the same grouping scheme. This confirms that overall the grouping scheme does not affect the simulation of the local large-scale structure.

(ii) However, when reaching the inner part of the box, the mean variance between the simulations increases by up to 20–25 per cent with respect to its minimum rather than continuing its decrease. While the result is completely expected when comparing random and constrained simulations – higher probability of finding a small empty region in the random simulation to be compared to the known

<footnote n="2">² https://github.com/ginnungagapgroup/ginnungagap</footnote>
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Supergalactic slices of the local large-scale structure obtained in constrained simulations. A 5 $h^{-1}$ Mpc smoothing scale has been applied to the fields. Different grouping schemes are used to remove non-linear motions from the constraint catalogue of galaxy radial peculiar velocities. Top: Tully grouping scheme, released with the catalogue and used for the first generation of constrained simulations based on the second catalogue of Cosmicflows. Bottom: Tempel grouping scheme, tested in this paper. The contours stand for the density. Solid lines show overdensities while dotted lines represent underdensities. The green colour is the mean field. A few structures are named in blue. Overall, the local large-scale structure is properly reproduced in both simulations. The differences appear only at the cluster scale. For instance, Coma and Centaurus are more clearly defined in the simulation obtained with Tempel grouping scheme.</figcaption></figure>

structures in the very nearby Universe – in the comparison between constrained simulations based on different grouping schemes, the finding is entirely due to the small differences noted at the cluster scale in Fig. 1.

(iii) In the third panel, if the smallest by 40–45 per cent average variance is still that obtained when comparing constrained simulations, it is worth noticing that the mean variance increases by a factor of up to 3–4 with the decrease in size of the compared regions. While this is entirely due to the fact that the weakly constrained part of the box – hence the random realization – dominates to a large extent when comparing the totality of the random and constrained boxes, when comparing constrained simulations it emphasizes that the large-scale structure is quite unaffected by the grouping scheme down to volumes of ~$(100 \ h^{-1} \ Mpc)^3$. The latter affects the simulations only at the cluster scale. The shape of the curve drawn by the variances when comparing only the inner part of the box is indeed similar when comparing only constrained simulations sharing the

same random realization and when comparing constrained simulations whatever random realization they have been constructed of.

(iv) Within the inner part of the box, Tempel grouping scheme results in constrained simulations that differ by 20 per cent more from the random simulations than those obtained with Tully grouping scheme. This is in agreement with the results found in Sorce & Tempel (2017): the densities are more pronounced with Tempel grouping scheme than with Tully's; hence, the constrained simulations differ more from the random ones in the former case than in the latter: a majority of underdensities are compared with higher overdensities.

Before focusing on the clusters (dark matter haloes) in a detailed way in the two different types of constrained simulations, it is worth comparing the power spectra and mass functions of the simulations. Fig. 3 shows the $1\sigma$ confidence interval of the ratios of the power spectra and mass functions of the entire box and for a 160 $h^{-1}$ Mpc
<figure><img src="image_3.png" /><figcaption>**Figure 2.** Average variance (filled circle) and its standard deviation (error bar) between density fields of simulations as a function of the size of the compared sub-box. From left to right: comparisons of pairs of random (R, black) and constrained ($C_{\text{Tully}}$ dark blue, $C_{\text{Tempel}}$ light blue) simulations, comparisons between random and constrained simulations (green and yellow) as well as between constrained simulations obtained with different grouping schemes (red) that do not share the same random realization (middle) and that share the same random realization (right). The goal of the constrained simulations is fulfilled: the cosmic variance is reduced with respect to that of random simulations.</figcaption></figure>

radius sphere of constrained and random simulations. The Amiga's halo finder is used to find the dark matter haloes in all the simulations (Knollmann & Knebe 2009).

Overall, the power spectra of the constrained simulations are below those of the random simulations on large scales as already noticed by Sorce et al. (2016b). Tests conducted on mock catalogues are in favour of the data as the most likely culprit, either as an intrinsic property or/and because of their modelling via, for instance, their grouping (rather than the succession of well-established mathematical procedures including Wiener filtering, Reverse Zel'dovich Approximation, and Constrained Realizations). This paper focuses on studying the impact of the data grouping modelling. Indeed, before any possibility of concluding that this behaviour is an intrinsic property of our local environment, any data modelling must be investigated. The power spectra of the simulations obtained with Tempel grouping scheme have on average slightly higher (10 per cent) values than those of the simulations obtained with Tully grouping scheme on large scales: the light blue zone tends to be above the dark blue zone. Consequently, the grouping scheme is partly responsible for the observation made by Sorce et al. (2016b). Namely, the power spectra of the constrained simulations obtained with Tempel grouping scheme have smaller values than those of the random simulations but to a lesser extent than those obtained with Tully grouping scheme; alternatively the red zone is above 1.0 meaning that the power spectra of Tempel constrained simulations have higher values than those of Tully constrained simulations on the large scales. This improves the probability of the local power spectrum given the Planck power spectrum.

Regarding the mass functions in the entire box or in the sphere, constrained simulations tend to have less massive haloes than random simulations as already observed by Sorce et al. (2016b). The same discussion as above is here also valid. Tempel constrained simulations have on average 1.5 more massive (above $2 \times 10^{14} h^{-1} M_{\odot}$) haloes than Tully constrained simulations within the $160 h^{-1}$ Mpc radius sphere as shown by the light blue zone that is on average above the dark blue one or by the red zone that is clearly above 1.0 on average at the high-mass end. These observations reinforce our expectations that Tempel scheme produces constrained simulations with more massive haloes than Tully scheme.

Although the grouping scheme alleviates the tension between power spectra and mass functions of constrained and random simulations, further investigations (an investigation of the data uncertainty modelling is currently underway) are necessary to conclude as to the reason for the observed residual. It could be either an intrinsic property of our local environment or another data modelling that needs improvement or both. Still, the next section proves that the constrained simulations are completely valid at least within $30 h^{-1}$ Mpc where Virgo, Hydra, and Centaurus clusters are perfectly reproduced as well as for zoom-in simulations of these clusters.

## 3.2 Local clusters of galaxies

In this section, the dark matter haloes counterparts of local observed clusters are looked for and studied in the constrained simulations obtained with the two different grouping schemes. Five local clusters of different masses and at various distances from us are selected for further studies: Virgo, Coma, Perseus, Centaurus, and Hydra. Their unique counterpart in each one of the constrained simulations is searched for in the list of dark matter haloes extracted from the simulations with the halo finder. Note that only simulacra with masses higher than $10^{14} h^{-1} M_{\odot}$ are considered. In addition, distances between simulacra and observed clusters cannot exceed 30 per cent of the distance estimate of the clusters. Regardless, the most important point is that if simulacra are slightly shifted in positions with respect to the observed cluster, their shifts are consistent so that their locations do not differ significantly from each other as shown hereafter.

Fig. 4 gives the percentage of simulations in which a simulacrum of the observed clusters is found (top) as well as the average mass of the simulacra (bottom). Overall, Tempel scheme does not increase significantly the percentage of success in getting a simulacrum: 87±14 against 89±19 per cent. This is in agreement with the fact that the local large-scale structure is well simulated in both cases, namely there are overdensity regions at the location of clusters. However, Tempel scheme increases the average mass of all the simulacra especially those of Virgo and Centaurus. While the mass of Virgo candidates is doubled (a factor of 2.2) that of Centaurus candidates is more than fivefold (a factor of 5.1) to reach an
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Top: $1\sigma$ confidence interval of the ratio of the power spectra of constrained and random simulations (blue) and of constrained simulations (red). Middle and bottom: the same as the top panel but for the mass functions of the entire box and of a 160 $h^{-1}$ Mpc radius sphere.</figcaption></figure>

excellent agreement with recent observational estimates (within $2\sigma$ and $1\sigma$, respectively). Indeed, Tully (2015) published recently the virial masses of these local clusters in $M_{\odot}$ with distances consistent with $H_0 = 75 \text{ km s}^{-1} \text{ Mpc}^{-1}$. The only uncertainties related to the virial masses that are provided are those of the bi-weight projected virial radii. A propagation of uncertainty using the sole

bi-weight project virial radii is far from optimal. Since Tully (2015) also supplies us with the luminosity masses that follow a 1:1 relation for clusters more massive than $10^{14} h^{-1} M_{\odot}$, the difference between the luminosity and the virial masses gives a rough estimate of the virial mass uncertainty. Table 1 summarizes these masses to be compared with $M_{200}$ (i.e. the mass enclosed in a sphere with a mean density of 200 times the critical density of the Universe) also included in the table. This mass derived by the halo finder is known to be proportional to the virial mass (given by the halo finder) via a factor of $0.80 \pm 0.03$ (e.g. Sorce et al. 2016a). Assuming the virial masses given by both observational estimates and the halo finder to be roughly similar, we overplot them for comparisons on Fig. 4 with blue thick dashed lines as well as the $1\sigma$ uncertainty of the conversion with thinner lines. The orange dot-dashed lines represent a dynamical mass estimate of the Virgo cluster (Lee, Kim & Rey 2015) that can be assumed to be roughly similar to a virial mass estimate for an unrelaxed cluster. Although this value is higher than the general values found in the literature for the mass of the Virgo cluster (e.g. Karachentsev et al. 2014, to give another reference), it is interesting to mention this value since it is based on the reconstruction of the dynamics of galaxies in filaments around the Virgo cluster. This value is thus obtained via both observation and numerical reconstructions like our value. The red dotted lines stand for a very recent estimate of Virgo's mass via the first turn around radius by Shaya et al. (2017) with $M_{200} = 4.9 \pm 0.7 M_{\odot}$ (Sorce et al. 2016a). The average mass of the first Virgo clusters is within $2\sigma$ of this estimate validating our previous study of the Virgo cluster with the constrained simulations that stated the good quality of the simulacra. Centaurus simulacra are now in extraordinary agreement with observations: their average mass is within $1\sigma$ of the estimated mass. Hydra simulacra have now masses in excellent agreement with observational estimates as well: in Table 1, the means are quasi-identical ($2.6$ and $2.5 \times 10^{14} h^{-1} M_{\odot}$) and the standard deviation of the simulated halo masses is almost equal to the uncertainty of the observational mass estimate ($0.8$ and $0.7 \times 10^{14} h^{-1} M_{\odot}$). While the mass of the Coma cluster is only increased by 24 per cent for Tempel scheme with respect to Tully's, it is worth noticing that it is now present in 100 per cent of the simulations. The only cluster that has less efficient simulacra is Perseus. However Perseus-Pisces region is poorly constrained with the second catalogue of Cosmicflows. The newly released third catalogue (Tully, Courtois & Sorce 2016) that contains more constraints in that region allows us to foresee good simulacra for Perseus in a near future.

It is important to note that it is the combination of the grouping scheme and the *bias minimization* that allows us to get such results. Without the bias minimization scheme, Virgo's success rate drops to 60 per cent with masses barely above $10^{14} h^{-1} M_{\odot}$ as for Centaurus although its success rate is 100 per cent, the most massive simulacrum is about $4 \times 10^{14} h^{-1} M_{\odot}$ and the majority of the simulacra have masses barely above $10^{14} h^{-1} M_{\odot}$.

The top of Fig. 5 shows the relative change between the properties of the simulacra obtained in the two sets of constrained simulations. We define the relative change as the difference between the parameter value of the dark matter halo in Tempel constrained simulation and that in Tully's divided by the value in Tully's. Interestingly, Virgo already showed to be very stable (Sorce et al. 2016a) is quasi unchanged in Tempel constrained simulations with respect to Tully's. The only exception is the z supergalactic coordinate. This is expected as the z direction is the less constrained because of the zone of avoidance. The other clusters present simulacra those relative changes are the most important for the velocity components. Perseus presents also relative changes of 5–10 times
<table><caption>**Table 1.** Clusters from Tully (2015) with $H_0 = 75$ km s$^{-1}$ Mpc$^{-1}$: (1) cluster name, (2) supergalactic longitude, (3) supergalactic latitude, (4) distance, (5) virial mass, (6) virial mass converted. Simulated clusters from this paper: (7) $M_{200}$ for Tully grouping, (8) $M_{200}$ for Tempel grouping.</caption><thead><tr><th rowspan="2">Cluster</th><th colspan="4">Observation</th><th colspan="3">Simulation</th></tr><tr><th>(2)<br>sgl<br>(°)</th><th>(3)<br>sgb<br>(°)</th><th>(4)<br>d<br>(Mpc)</th><th>(5)<br>M<br>(10<sup>14</sup> M<sub>☉</sub>)</th><th>(6)<br>M<sub>200</sub><br>(10<sup>14</sup> h<sup>-1</sup> M<sub>☉</sub>)</th><th>(7)<br>M<sub>200</sub><br>(10<sup>14</sup> h<sup>-1</sup> M<sub>☉</sub>)</th><th>(8)<br>M<sub>200</sub><br>(10<sup>14</sup> h<sup>-1</sup> M<sub>☉</sub>)</th></tr></thead><tbody><tr><td>Virgo</td><td>103.0008</td><td>-2.3248</td><td>14.9</td><td>7.01 ± 1.7</td><td>4.21 ± 1.0</td><td>3.05 ± 0.5</td><td>6.6 ± 0.5</td></tr><tr><td>Centaurus</td><td>156.2336</td><td>-11.5868</td><td>38.7</td><td>10.8 ± 3.9</td><td>6.48 ± 2.3</td><td>1.51 ± 0.4</td><td>7.58 ± 1.3</td></tr><tr><td>Hydra</td><td>139.4478</td><td>-37.6063</td><td>41.0</td><td>4.39 ± 1.3</td><td>2.63 ± 0.8</td><td>1.55 ± 0.3</td><td>2.50 ± 0.7</td></tr><tr><td>Perseus</td><td>347.7159</td><td>-14.0594</td><td>52.8</td><td>16.3 ± 4.2</td><td>9.78 ± 2.5</td><td>2.32 ± 1.6</td><td>2.46 ± 1.8</td></tr><tr><td>Coma</td><td>89.6226</td><td>8.1461</td><td>73.3</td><td>15.9 ± 1.2</td><td>9.54 ± 0.7</td><td>1.78 ± 0.5</td><td>2.2 ± 0.9</td></tr></tbody></table>

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Top: percentage of simulations with a simulacrum of the observed cluster whose name is given in the top right corner of the panel for both types of grouping schemes. Bottom: average mass (histogram), standard deviation (error bar) of the different simulacra for both types of grouping schemes. Clearly the average masses of the different simulacra are higher when using Tempel grouping scheme. Blue thick dashed lines show the virial mass estimates from Tully (2015) converted to $M_{200}$. The blue thinned dashed lines show the $1\sigma$ uncertainty in the conversion from $M_{vir}$ to $M_{200}$. The red thick (thin) dotted lines show the latest observational mass estimate of the Virgo cluster from the first turn around radius converted to $M_{200}$ ($\pm 1\sigma$ uncertainty in the conversion) while the orange dot-dashed lines stand for a dynamical mass estimate – equivalent to a virial mass for unrelaxed clusters – obtained studying galaxies in filaments falling into the Virgo cluster.</figcaption></figure>

the y and z supergalactic coordinates confirming that it is the less constrained clusters in terms of positions. Still, all the maximum relative changes are held below 10 times the parameters.

The bottom of Fig. 5 gives the mean variation (filled circle) and standard deviation (error bar) of the properties of the dark matter haloes in Tempel (red) and Tully (blue) constrained simulations taken separately. The variation is defined as the standard deviation of a parameter divided by the parameter value of a given simulacrum. Virgo appears very stable for both grouping schemes and even more stable (by a factor of 5) in terms of the z supergalactic coordinate when using Tempel grouping scheme rather than Tully's. The y and z velocity components of Centaurus counterparts are more stable when using Tempel grouping scheme. For the y component, the standard deviation of the variation decreases from 18 per cent to less than 3 per cent. The other clusters have more mitigated variation and standard deviation values.

Overall, the efficiency of changing the grouping scheme for a less aggressive scheme (see Sorce & Tempel 2017, for a detailed discussion) results in galaxy clusters more:

(i) present: on average, clusters are simulated with a success rate increased by 5 per cent with respect to Tully grouping scheme;

(ii) stable: on average, the stability of the parameters is increased by a factor of 3 with respect to the parameters of haloes obtained with Tully grouping scheme;

(iii) massive: on average, haloes are 39 per cent more massive than those obtained with Tully grouping scheme.

These assertions are reinforced when considering Centaurus and Virgo: Centaurus simulated mass is now within $1\sigma$ of the observation mass estimate. Virgo is perfected: its z-component is in particular more constrained than before with a standard deviation decreased by a factor of 5.

## 4 CONCLUSION

Galaxy clusters are excellent cosmological probes whose formation and evolution still need to be understood in detail. Numerical simulations of clusters constitute a formidable complementary
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Top: relative change and standard deviation (filled circle and error bar) between the parameters of the cluster simulacra obtained in the simulations produced with the two different grouping schemes. Bottom: variation (filled circle) and standard deviation (error bar) of the parameters of the cluster simulacra found in the simulations obtained with Tully (blue) and Tempel (red) grouping schemes.</figcaption></figure>

approach to their observations. However, the diversity of galaxy clusters complicates comparisons with their numerical counterparts on a one-to-one basis down to the simulated and observed galaxy populations.

Such detailed comparisons are feasible in the context of simulations that resemble the local Universe provided that the latter reproduce the local large-scale structure down to the cluster scales. In our first generation of constrained simulations made with the second catalogue of galaxy peculiar velocity of the Cosmicflows project combined with a bias minimization scheme, large overdensities were present at the location of local clusters but massive enough dark matter simulacra of the latter were found only for the Virgo cluster.

In a previous study, we showed that one of the essential steps in the process of building the constrained initial conditions affects the overdensity values of the reconstructed field. This step consists in grouping the catalogue of constraints (galaxies and their peculiar velocities) to remove non-linear motions that would affect the linear reconstruction. However, this first study demonstrated that the grouping must be made with parsimony to preserve the infall on clusters and thus to increase the local densities. This study goes further as it probes the impact of the grouping scheme (Tully's based on the literature and Tempel's based on an advanced FoF algorithm) on the final product: the simulations that resemble the local Universe.

Overall, the same large-scale structure is simulated with both grouping schemes. However, a slight increase (10 per cent) of the power spectrum on large scales is observed with the less aggressive grouping scheme. This is ideal as it improves the probability of the power spectrum of the local Universe given Planck power spectrum. Tempel grouping scheme also increases the mass function at the high end, more precisely the most massive haloes are heavier than with Tully grouping scheme. These new observations imply the real need to inquire further on the impact of the data modelling on the resulting simulations. It is essential to determine whether the

residual is intrinsic to the data or due to some data modelling or both. This study has shown that the grouping is partly responsible; in an ongoing study, we are investigating the impact of the data uncertainty modelling.

A thorough study of five of the local clusters (Virgo, Centaurus, Hydra, Coma, and Perseus) reveals that their simulacra are better representative. The Virgo simulacrum is still very stable (present in 100 per cent of the simulations) and its mass is increased by 50 per cent with respect to the first generation of constrained simulations we produced. In both cases, the masses are within $2\sigma$ of the observational mass estimate, the previous one on the low side, the new one on the upper side. The z supergalactic coordinate of the simulacra presents a standard deviation divided by a factor of 5, implying that although the z direction is that of the zone of avoidance, it is possible to constrain further the z position of the Virgo dark matter haloes with a moderate grouping scheme. The most incredible advantage of Tempel grouping scheme is visible for the Centaurus cluster. In the new set of constrained simulations, Centaurus simulacra are five times more massive than before and are within $1\sigma$ of the recent observational estimates. Coma is also improved in the sense that a simulacrum is now present in 100 per cent of the simulation with a mass increased by nearly 25 per cent. Hydra's and Perseus simulacra are also slightly more massive than in the first generation of simulations. The formers, present in 80 per cent of the simulations, have the quasi same mean (within 5 per cent) and standard deviation (within 12 per cent) as the observational estimate.

All in all, using Tempel grouping scheme improves considerably the simulacra of Centaurus and perfects those of Virgo *provided that it is combined with the bias minimization scheme*. Indeed, without the latter, all the advantages of using Tempel rather than Tully grouping scheme disappear. The combination of Tempel grouping scheme and the bias minimization scheme ameliorates the simulacra of Coma, Hydra, and Perseus although there is still room
for improvements. First, the third catalogue of peculiar velocities of the Cosmicflows project will offer us more data especially in the direction of Perseus-Pisces. Secondly, a better modelling of the uncertainties in the bias minimization scheme (so far a 5 per cent uncertainty is applied to all the distances obtained after minimization and is propagated to the velocities) is under study. Thirdly, a new grouping algorithm based on point processes with interactions is investigated.

Now, the constrained simulations of the local Universe produced via the method described in this paper that still uses *only peculiar velocity data sets as constraints* (in the sense that no additional density constraints are added at the positions of the clusters, the velocity constraints contain both the position and mass information of the clusters) not only resemble the local large-scale structure and have Virgo dark matter simulacra but also stable Centaurus haloes with masses within $1\sigma$ of observational estimates as well as better representatives of Coma, Hydra, and Perseus.

A large number of zoom-in dark matter simulations of these haloes will permit making statistical studies of these local clusters regarding their formation, their substructures, etc. In addition, further zoom-in hydrodynamical simulations of these haloes are planned to study the galaxy populations of these various local clusters to be compared with their observational counterparts. Links between properties of galaxy populations in local clusters of different types (various masses, formation histories, substructures, etc.) will be highlighted to further refine our understanding of galaxy formation and evolution in clusters.

## ACKNOWLEDGEMENTS

The authors would like to thank the referee Francisco-Shu Kitaura and Stefan Gottlöber for providing comments that helped clarified the paper as well as Sergey Pilipenko for very useful discussions. JS acknowledges support from the Astronomy ESFRI and Research Infrastructure Cluster ASTERICS project, funded by the European Commission under the Horizon 2020 Programme (GA 653477) as well as from the ‘l’Oréal-UNESCO Pour les femmes et la Science’ and the ‘Centre National d’études spatiales (CNES)’ postdoctoral fellowship programmes. ET was supported by ETAg grants IUT40-2, IUT26-2 and by EU through the ERDF CoE grant TK133. The authors gratefully acknowledge the Gauss Centre for Supercomputing e.V. (http://www.gauss-centre.eu) for providing computing time on the GCS Supercomputers SuperMUC at LRZ Munich and Jureca at JSC Juelich.

## REFERENCES

* Abell G. O., Corwin H. G., Jr, Olowin R. P., 1989, ApJS, 70, 1

* Arthur J. et al., 2017, MNRAS, 464, 2027

* Baldi A. S., De Petris M., Sembolini F., Yepes G., Lamagna L., Rasia E., 2017, MNRAS, 465, 2584

* Barnes D. J., Kay S. T., Henson M. A., McCarthy I. G., Schaye J., Jenkins A., 2017, MNRAS, 465, 213

* Bertschinger E., 1987, ApJ, 323, L103

* Bertschinger E., Dekel A., 1989, ApJ, 336, L5

* Bertschinger E., Dekel A., Faber S. M., Dressler A., Burstein D., 1990, ApJ, 364, 370

* Colless M., Saglia R. P., Burstein D., Davies R. L., McMahan R. K., Wegner G., 2001, MNRAS, 321, 277

* Cui W. et al., 2016, MNRAS, 458, 4052

* Dekel A., Bertschinger E., Faber S. M., 1990, ApJ, 364, 349

* Doumler T., Hoffman Y., Courtois H., Gottlöber S., 2013, MNRAS, 430, 888

* Elahi P. J. et al., 2016, MNRAS, 458, 1096

* Freedman W. L. et al., 2001, ApJ, 553, 47

* Ganon G., Hoffman Y., 1993, ApJ, 415, L5

* Gottlöber S., Hoffman Y., Yepes G., 2010, preprint (arXiv:1005.2687)

* Grossauer J. et al., 2015, ApJ, 807, 88

* Hahn O., Martizzi D., Wu H.-Y., Evrard A. E., Teyssier R., Wechsler R. H., 2017, MNRAS, 470, 166

* Heß S., Kitaura F.-S., Gottlöber S., 2013, MNRAS, 435, 2065

* Hoffman Y., Ribak E., 1991, ApJ, 380, L5

* Hoffman Y., Ribak E., 1992, ApJ, 384, 448

* Jackson J. C., 1972, MNRAS, 156, 1P

* Jasche J., Wandelt B. D., 2013, MNRAS, 432, 894

* Jha S., Riess A. G., Kirshner R. P., 2007, ApJ, 659, 122

* Jimeno P., Broadhurst T., Lazkoz R., Angulo R., Diego J.-M., Umetsu K., Chu M.-c., 2017, MNRAS, 466, 2658

* Kaiser N., 1987, MNRAS, 227, 1

* Karachentsev I. D., Tully R. B., Wu P.-F., Shaya E. J., Dolphin A. E., 2014, ApJ, 782, 4

* Kitaura F. S., Enßlin T. A., 2008, MNRAS, 389, 497

* Kitaura F. S., Jasche J., Li C., Enßlin T. A., Metcalf R. B., Wandelt B. D., Lemson G., White S. D. M., 2009, MNRAS, 400, 183

* Kitaura F.-S., 2013, MNRAS, 429, L84

* Kitaura F.-S., Erdodu P., Nuza S. E., Khalatyan A., Angulo R. E., Hoffman Y., Gottlöber S., 2012, MNRAS, 427, L35

* Klypin A., Hoffman Y., Kravtsov A. V., Gottlöber S., 2003, ApJ, 596, 19

* Knollmann S. R., Knebe A., 2009, ApJS, 182, 608

* Lavaux G., 2010, MNRAS, 406, 1007

* Lavaux G., 2016, MNRAS, 457, 172

* Lavaux G., Mohayaee R., Colombi S., Tully R. B., Bernardeau F., Silk J., 2008, MNRAS, 383, 1292

* Lee M. G., Freedman W. L., Madore B. F., 1993, ApJ, 417, 553

* Lee J., Kim S., Rey S.-C., 2015, ApJ, 807, 122

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

* Natarajan P. et al., 2017, MNRAS, 468, 1962

* Planck Collaboration XVI et al., 2014, A&A, 571, A16

* Sembolini F. et al., 2016a, MNRAS, 459, 2973

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

* Shaya E. J., Tully R. B., Hoffman Y., Pomarède D., 2017, ApJ, 850, 207

* Sorce J. G., 2015, MNRAS, 450, 2644

* Sorce J. G., Tempel E., 2017, MNRAS, 469, 2859

* Sorce J. G., Courtois H. M., Gottlöber S., Hoffman Y., Tully R. B., 2014, MNRAS, 437, 3586

* Sorce J. G., Gottlöber S., Hoffman Y., Yepes G., 2016a, MNRAS, 460, 2015

* Sorce J. G. et al., 2016b, MNRAS, 455, 2078

* Sorce J. G., Hoffman Y., Gottlöber S., 2017, MNRAS, 468, 1812

* Springel V., 2005, MNRAS, 364, 1105

* Struble M. F., Rood H. J., 1988, S&T, 75, 16

* Suto D., Peirani S., Dubois Y., Kitayama T., Nishimichi T., Sasaki S., Suto Y., 2017, PASJ, 69, 14

* Tempel E., Kipper R., Tamm A., Gramann M., Einasto M., Sepp T., Tuvikene T., 2016, A&A, 588, A14

* Tonry J. L., Dressler A., Blakeslee J. P., Ajhar E. A., Fletcher A. B., Luppino G. A., Metzger M. R., Moore C. B., 2001, ApJ, 546, 681

* Tully R. B., 1987, ApJ, 321, 280

* Tully R. B., 2015, AJ, 149, 171

* Tully R. B., Fisher J. R., 1977, A&A, 54, 661

* Tully R. B., Shaya E. J., Karachentsev I. D., Courtois H. M., Kocevski D. D., Rizzi L., Peel A., 2008, ApJ, 676, 184

* Tully R. B., Rizzi L., Shaya E. J., Courtois H. M., Makarov D. I., Jacobs B. A., 2009, AJ, 138, 323

* Tully R. B. et al., 2013, AJ, 146, 86

* Tully R. B., Courtois H. M., Sorce J. G., 2016, AJ, 152, 50

* van de Weygaert R., Bertschinger E., 1996, MNRAS, 281, 84

* Wang H., Mo H. J., Yang X., van den Bosch F. C., 2013, ApJ, 772, 63
Wang H., Mo H. J., Yang X., Jing Y. P., Lin W. P., 2014, ApJ, 794, 94

Wang H. et al., 2016, ApJ, 831, 164

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

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

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

Zaroubi S., Hoffman Y., Fisher K. B., Lahav O., 1995, ApJ, 449, 446

Zaroubi S., Hoffman Y., Dekel A., 1999, ApJ, 520, 413

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