EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 691 (November 2024) A&A, 691 (2024) A301 Full HTML
Open Access
Issue
A&A
Volume 691, November 2024
Article Number A301
Number of page(s) 14
Section Cosmology (including clusters of galaxies)
DOI https://doi.org/10.1051/0004-6361/202449587
Published online 21 November 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model.

Open Access funding provided by Max Planck Society.

1 Introduction

Our understanding of the Universe has been altered significantly with the discovery of its accelerated expansion, first demonstrated by supernovae observations (Riess et al. 1998; Perlmutter et al. 1999). The paradigm of the standard cosmological model widely represents this unknown dark energy with a constant added in Einstein’s gravitational equations, Λ, the so-called cosmological constant. Other ingredients necessary to describe our cosmological model include a nonrelativistic (cold) dark matter (DM) component that represents most of the mass fraction of the large-scale structures (LSSs). These two components constitute the major fraction of the energy content of the present-day Universe and form the backbone of the concordance ΛCDM model.

So far, when considered together with inflation, this parametrization is relatively successful at describing the main observed properties of the Universe (Planck Collaboration VI 2020).

Galaxy clusters represent the most massive virialized structures in the Universe. The number density of clusters of galaxies per unit mass and redshift – the cluster mass function – and their spatial distribution – the two-point correlation function and higher-order statistics – have great potential to constrain the nature of dark energy and test the robustness of the underlying cosmological models (see Clerc & Finoguenov 2023, for a recent review). Cluster counts are currently the most efficient at placing constraints on the matter density parameter (Ωm), the root mean square of the density fluctuations in spheres of 8 Mpc h−1 (σ8), and the sum of the masses of the neutrinos (Σ mν). Additionally, they also break the degeneracies between the cosmological parameters of other probes, such as baryon acoustic oscillations (BAO; Eisenstein et al. 2005; Alam et al. 2017), redshift space distortions (RSD; Percival & White 2009), the weak-lensing (WL) cosmic shear (Troxel et al. 2018; Asgari et al. 2021; Amon et al. 2022; Li et al. 2023), type Ia supernovae and other distance ladders (Riess et al. 2019; Freedman et al. 2019), and the cosmic microwave background (CMB; Bennett et al. 2013; Planck Collaboration VI 2020). Previous and ongoing cluster surveys across a wide wavelength range demonstrate the potential of cluster number counts as a powerful cosmo-logical probe when the biases related to mass calibration are reduced through the utilization of weak gravitational lensing observations (Mantz et al. 2015; Bocquet et al. 2019; Zubeldia & Challinor 2019; Ider Chitham et al. 2020; Garrel et al. 2022; Lesci et al. 2022; Sunayama et al. 2023; Fumagalli et al. 2024; Bocquet et al. 2024a,b, among others). Modern cluster surveys, such as the eROSITA All-Sky cluster survey (see Bulbul et al. 2024; Ghirardini et al. 2024), are revolutionizing the field of cosmology, promising to reach the precision of the CMB experiments on cosmological parameters, and have the potential to measure departures from Einstein’s theory of general relativity (GR) and test theories of dark energy (Wolf & Ferreira 2023). In addition to cluster abundance, eROSITA can probe cosmology through the spatial distribution of galaxy clusters (Seppi et al. 2024) and the X-ray power spectrum.

It also provides an avenue for testing alternative gravity theories, which could, in return, explain the late-time accelerated expansion of the Universe. Indeed, a way of explaining the Universe’s accelerated expansion is to modify the theory of gravitation on large scales, while retaining GR in the high-density Solar System and on galactic scales. Potential departures from GR can be modeled by a set of theories in which the Ricci scalar (R), the quantity describing the infinitesimal volumes in curved space-time, is replaced in the Einstein-Hilbert action by a generic function, f(R). One of the models introducing a physically motivated function is the one proposed by Hu & Sawicki (2007) (HS- f(R) hereafter). This model has direct phenomenological implications. For example, Tsujikawa et al. (2009) and Gannouji et al. (2009) showed that it would significantly increase the structure formation rate, and thus its predictions can be verified. In the context of cluster abundance studies, it increases the expected number density of massive DM halos. Many recent attempts have been made to provide tools for computing summary statistics measuring departures from GR, with the goal of constraining this model observationally. First, numerous codes adapt the standard Boltzmann solvers; for example, CAMB (Lewis & Challinor 2011) and CLASS (Lesgourgues 2011) for modified gravity (MG hereafter) studies with alternative software such as EFTCAMB (Hu et al. 2014), hi_class (Zumalacárregui et al. 2017), MGCAMB (Zucca et al. 2019; Wang et al. 2023), or FRCAMB Xu (2015). These extensions can compute the power spectrum in different MG scenarios.

However, the hierarchical structure formation scenario assumes that galaxy clusters originate from the gravitational collapse of matter overdensities. The cluster abundance can thus be obtained by computing the characteristics of this collapse in the corresponding theory (Kopp et al. 2013; Lombriser et al. 2013), while considering the standard power spectrum. Starting from this fact, several studies have attempted to constrain the f (R) parameters using clusters in the literature. In the first of the studies, Schmidt et al. (2009) used a sample of 49 low-redshift clusters detected by ROSAT and re-observed by Chandra (Vikhlinin et al. 2009) to constrain the HS-f(R) model. This approach was combined with other cosmological probes by Lombriser et al. (2012), showing that the best constraints are obtained when cluster abundance is considered. Cataneo et al. (2015) report a robust upper limit for the value of the background scalar field, log |fR0| < −4.79, by utilizing a sample of 224 ROSAT-detected clusters combined with CMB data. Additionally, Lombriser et al. (2013) and later Cataneo et al. (2016) developed a halo mass function (HMF) formalism to predict accurate cluster abundances, which was followed by the formalisms of von Braun-Bates et al. (2017) and Gupta et al. (2022). These studies explore the implication of f(R) gravity in the halo formation process. In parallel, Hagstotz et al. (2019) developed an alternative method whereby the HMF modeling accounts for the massive neutrinos and MG effects. The modifications of gravity are then encapsulated in the collapse model and the HMF parametrization. In this work, we have used this approach to predict the cluster counts. We emphasize again that all the cosmological quantities related to the matter power spectrum were computed with the standard version of CAMB, without any need for an application of the modified Boltzmann solver.

With 5259 clusters of galaxies selected in its cosmology sample in a wide redshift range, 0.1 < z < 0.8, the western Galactic half of the first eROSITA All-Sky Survey (eRASS1) has the statistical power to reach far beyond constraints available from any cluster survey in the literature. Consequently, this work utilizes the largest intracluster medium (ICM) selected cluster sample to date, in combination with WL follow-up measurements for cluster masses and an accurate selection function to constrain deviations from GR by placing robust limits on the HS parametrization of f(R) gravity.

The paper is organized as follows. In Sect. 2, we describe the datasets used. Section 3 introduces the framework used to describe departures from GR in the HS- f(R) gravity. Section 4 provides a detailed description of the statistical methods. Finally, Sect. 5 presents the constraints on the HS- f(R) gravity and a comparison with the ΛCDM concordance model. Throughout this paper, we use the notation log ≡ log10 and ln ≡ loge.

2 Survey data

The data used in this work are identical to the survey data presented in Bulbul et al. (2024); Kluge et al. (2024); Ghirardini et al. (2024); Grandis et al. (2024). We summarize the multi-wavelength data utilized in this work, including X-ray, optical, and WL observations, in the following section.

2.1 eROSITA All-Sky Survey

The soft X-ray telescope eROSITA (Predehl et al. 2021) on board the Spectrum Roentgen Gamma Mission completed its All-Sky Survey program on June 11, 2020, 184 days after the start of the survey. In this work, we have used the data of the western Galactic half of the eROSITA All-Sky survey (359.9442 deg > l > 179.9442 deg) belonging to the German eROSITA consortium. The raw data were processed through the eSASS software (Brunner et al. 2022). Bulbul et al. (2024); Kluge et al. (2024) compiled two catalogs of clusters detected in this region: the primary galaxy group and cluster catalog and the cosmology sample based on the catalog of the X-ray sources in the soft band (0.2–2.3 keV) provided in Merloni et al. (2024). We utilized the eRASS1 cosmology sample in this work, which was compiled adopting a selection cut of the extent likelihood, ℒext > 6 (see Bulbul et al. 2024, for further details). The optical identification of the clusters selected in the LS DR10-South common footprint produces a final cosmology sample of 5259 clusters of galaxies reaching purity levels of 95%, an ideal sample for the MG studies. Unlike the primary cluster catalog, the redshifts in the cosmology subsample are purely photometric measurements, enabling a consistent assessment of the systematic effects in our analysis. The total survey area covers a 12 791 deg2 region in the western Galactic hemisphere. The count rates were extracted with the 2D-image fitting tool MBProj2D described in Bulbul et al. (2024) and were employed as a mass proxy in the WL mass calibration likelihood (Ghirardini et al. 2024).

2.2 Weak lensing survey data

To perform our mass calibration in a minimally biased way, we utilized the deep and wide-area optical surveys for lensing measurements with the overlapping footprint with eROSITA in the western Galactic hemisphere. The three wide-area surveys used in this work for mass calibration include the Dark Energy Survey (DES), Kilo Degree Survey (KiDS), and the Hyper Suprime-Cam Survey (HSC). We briefly describe the WL survey data here.

The three-year WL data (S19A) covered by the HSC Subaru Strategic Program (Aihara et al. 2018; Li et al. 2022) WL measurements were made around the location of eRASS1 clusters using in the 𝑔, r, i, z, and Y bands in the optical. The total area coverage of HSC is ≈ 500 deg2. The shear profile, 𝑔t (θ), the lensing covariance matrix that serves the measurement uncertainty, and the photometric redshift distribution of the selected source sample were obtained as the HSC WL data products for 96 eRASS1 clusters, with a total signal-to-noise of 40.

We utilized data from the first three years of observations of the Dark Energy Survey (DES Y3). The DES Y3 shape catalog (Gatti et al. 2021) was built from the r, i, z bands using the METACALIBRATION pipeline (Huff & Mandelbaum 2017; Sheldon & Huff 2017). Considering the overlap between the DES Y3 footprint and the eRASS1 footprint, we produced tangential shear data for 2201 eRASS1 galaxy clusters, with a total signal-to-noise of 65 in the tangential shear profile. The analysis and shear profile extraction details are presented in Grandis et al. (2024).

We used the gold sample of WL and photometric redshift measurements from the fourth data release of the Kilo-Degree Survey (Kuijken et al. 2019; Wright et al. 2020; Hildebrandt et al. 2021; Giblin et al. 2021), hereafter referred to as KiDS-1000. We extracted individual reduced tangential shear profiles for a total of 236 eRASS1 galaxy clusters in both the KiDS-North field (101 clusters) and the KiDS-South field (136 clusters), as both overlap with the eRASS1 footprint with a total signal-to-noise of 19.

The WL data in the mass calibration process is detailed and published in Ghirardini et al. (2024); Grandis et al. (2024). We did not reprocess the shear profiles and used the ones used in Ghirardini et al. (2024). This is justified in Sect. 3.3.

3 Modified gravity framework

The evolution of the Universe’s LSSs is sensitive to modification to the standard frameworks of GR and ΛCDM. To investigate potential discrepancies, we need to probe different scales and epochs. Figure 1 shows the redshift and scale sensitivity of the most common cosmological probes. Galaxy cluster number counts efficiently place constraints on the low-redshift or large-scale regime not explored by other probes. This section describes the framework explored in this paper; namely, the HS-f(R) model.

thumbnail Fig. 1

Estimation of the scale and redshift sensitivity of the different cosmological probes, adapted from Preston et al. (2023). The background colors, from blue to green, represent the transition from the linear to the nonlinear regime. We represent cluster abundance with the red rectangle. The cosmology sample of eRASS1 spans the redshift range 0.1 < z < 0.8, as is represented by the solid red rectangle. Future data releases and experiments will increase the upper limit of this interval to the hatched red rectangle. Cluster abundance experiments thus probe a specific regime at LSSs and the nearby Universe, offering a complementary window on structure formation. The constraints from cluster clustering (Seppi et al. 2024) are omitted here. Additionally, the X-ray power spectrum might probe smaller scales.

3.1 HS-f(R) gravity in the context of cluster counts

We chose to test the MG model that consists of a modification of the Einstein-Hilbert action in the following form: S=∫​dx4−g(R+f(R)16πG+ℒm),$S = {\mathop \smallint \nolimits^ ^}{\rm{d}}{x^4}\sqrt { - g} \left( {{{R + f(R)} \over {16\pi G}} + {{\cal L}_m}} \right),$(1)

where ℒm represents the Lagrangian of the matter field, R is the Ricci scalar, and f (R) represents the extension to GR, for which f (R) = −2Λ. In the context of f(R) gravity, many functional forms are proposed (see De Felice & Tsujikawa 2010, for a review). Throughout this paper, we focus on the Hu & Sawicki (2007) model, which formalizes f(R) as f(R)=−2ΛRnRn+m2n,$f(R) = - 2{\rm{\Lambda }}{{{R^n}} \over {{R^n} + {m^{2n}}}},$(2)

where m2 represents the curvature scale, and Λ and n (the scaling index) are constants. In the high-curvature regime (R ≫ m2), Eq. (2) becomes f(R)=−2ΛRnRn+m2n=−2Λ11+(m2R)n∼−2Λ(1−(m2R)n).$f(R) = - 2{\rm{\Lambda }}{{{R^n}} \over {{R^n} + {m^{2n}}}} = - 2{\rm{\Lambda }}{1 \over {1 + {{\left( {{{{m^2}} \over R}} \right)}^n}}} \sim - 2{\rm{\Lambda }}\left( {1 - {{\left( {{{{m^2}} \over R}} \right)}^n}} \right).$

To express this function in terms of fR(R) = df /dR, we rewrote it, and obtained f(R)=−2Λ−fR0nR0¯n+1Rn,$f(R) = - 2{\rm{\Lambda }} - {{{f_{{\rm{R}}0}}} \over n}{{{{\overline {{R_0}} }^{n + 1}}} \over {{R^n}}},$

where we have defined the value of the field, fR0=df/dR(R0)=−2nΛm2n/R0¯n+1${f_{{\rm{R}}0}} = {\rm{d}}f/{\rm{d}}R\left( {{R_0}} \right) = - 2n\Lambda {m^{2n}}/{\overline {{R_0}} ^{n + 1}}$, at z=0,R¯0$z = 0,{\bar R_0}$ is the Ricci scalar at z = 0, and overbars denote background quantities. State-of-the-art numerical simulations are required to build a precise HMF suitable for cluster abundance cosmology. In these simulations, it is common to use n = 1 (e.g. Oyaizu 2008; Zhao et al. 2011; Giocoli et al. 2018). For consistency with the HMF calibration method, we thus froze the scaling index to unity. We finally obtained f(R)≈−2Λ−fR0R¯02R.$f(R) \approx - 2{\rm{\Lambda }} - {f_{R0}}{{\bar R_0^2} \over R}.$(3)

In this model, we note that when |fR0| → 0, we recover the usual GR action with a cosmological constant in Eq. (1). Thus, we interpret Λ as the cosmological constant. The only parameter quantifying the deviations from GR is fR0. When |fR0| ≪ 1, the background expansion is not distinguishable from GR. In this paper, our main goal is to measure fR0 or to set an upper limit on its value. In the rest of this section, we briefly describe the behavior of this model.

The variational principle applied to the action in Eq. (1) leads to the modified Einstein’s equation: Gμν−fRRμv−(f2−☐fR)gμv−∇μ∇vfR=8πGTμv.${G_{\mu \nu }} - {f_{\rm{R}}}{R_{\mu v}} - \left( {{f \over 2} - {f_{\rm{R}}}} \right){g_{\mu v}} - {\nabla _\mu }{\nabla _v}{f_{\rm{R}}} = 8\pi G{T_{\mu v}}.$(4)

As usual, Eq. (4) provides the equation of motion of the scalar field, fR: ∇2δfR=a23(δR(fR)−8πGδρm),${\nabla ^2}\delta {f_{\rm{R}}} = {{{a^2}} \over 3}\left( {\delta R\left( {{f_{\rm{R}}}} \right) - 8\pi G\delta {\rho _m}} \right),$(5)

where δ denotes small perturbations of the related quantities. This also provides a Poisson-like equation given by ∇2ψ=16πG3a2ρm−a26δR(fR).${\nabla ^2}\psi = {{16\pi G} \over 3}{a^2}{\rho _m} - {{{a^2}} \over 6}\delta R\left( {{f_{\rm{R}}}} \right).$(6)

For a better understanding of the behavior in this model, we provide the thin-shell condition example first introduced by Khoury & Weltman (2004) and adapted in Li et al. (2012) and Lombriser et al. (2013). Since GR is well constrained in high-density environments like the Solar System, we want to recover it at small scales: this is the so-called “screening mechanism.” As an illustration, we shall consider a top-hat overdensity of radius Rin. Then, the gravitational pull exerted on a particle of mass, m, outside of the overdensity can be expressed as F=GMinmr2(1+13ΔRRin),$F = G{{{M_{{\rm{in}}}}m} \over {{r^2}}}\left( {1 + {1 \over 3}{{{\rm{\Delta }}R} \over {{R_{{\rm{in}}}}}}} \right),$(7)

where Min is the mass inside the top-hat overdensity, r is the distance to the center, and ∆R/Rin is the screening factor defined as ΔRRin=min{ 3| fRin−fRput |2ϕN,1 },${{{\rm{\Delta }}R} \over {{R_{{\rm{in}}}}}} = \min \left\{ {{{3\left| {f_{\rm{R}}^{{\rm{in}}} - f_{\rm{R}}^{{\rm{put}}}} \right|} \over {2{\phi _{\rm{N}}}}},1} \right\},$(8)

where ϕN is the Newtonian gravitational potential, ϕN = GMin/Rin, and fRin $f_{\rm{R}}^{{\rm{in }}}$ and fRout $f_{\rm{R}}^{{\rm{out }}}$ are the values of the scalar filed inside and outside of the top-hat overdensity. Considering that the background expansion is considered in ΛCDM, we can approximate the background Ricci scalar to the one expected in this model. Inside the top-hat overdensity, we consider that the screening mechanism is active; that is, we recover GR. We can then express the scalar field as fRin=fR0(1+4ΩΛ0/Ωm0ρmin/ρm0+4ΩΛ0/Ωm0),$f_{\rm{R}}^{{\rm{in}}} = {f_{{\rm{R}}0}}\left( {{{1 + 4{{\rm{\Omega }}_{{\rm{\Lambda }}0}}/{{\rm{\Omega }}_{{\rm{m}}0}}} \over {\rho _{\rm{m}}^{{\rm{in}}}/{\rho _{{\rm{m}}0}} + 4{{\rm{\Omega }}_{{\rm{\Lambda }}0}}/{{\rm{\Omega }}_{{\rm{m}}0}}}}} \right),$(9)

where ρmin$\rho _{\rm{m}}^{{\rm{in}}}$ is the density inside the considered region. In Eq. (9), we note that the scalar field values become smaller in high-density regions (ρm0≪ρmin)$\left( {{\rho _{{\rm{m}}0}} \ll \rho _{\rm{m}}^{{\rm{in}}}} \right)$. Thus, in Eq. (8), the screening factor becomes ΔRRin∼3| fRout |2ϕN.${{{\rm{\Delta }}R} \over {{R_{{\rm{in}}}}}} \sim {{3\left| {f_{\rm{R}}^{{\rm{out}}}} \right|} \over {2{\phi _{\rm{N}}}}}.$

The condition for the screening to be activated is ∆R/Rin ≪ 1. Thus, we must have | fRout  |≪ϕN$\left| {f_{\rm{R}}^{{\rm{out }}}} \right| \ll {\phi _{\rm{N}}}$. Since fR is a monotonic function of time, a relevant condition would be to have | fR0 |<ϕN.$\left| {{f_{{\rm{R}}0}}} \right| < {\phi _{\rm{N}}}.$

For massive galaxy clusters, the Newtonian potential (in natural units) is on the order of 10−5. Thus, for these gravitationally bound systems, the deviation to GR is screened (i.e., suppressed) when |fR0| < 10−5 (see Sect. 3.3). Overall, the HS-f(R) gravity model is defined with only one additional parameter: fR0. This formalism allows the recovery of GR in the necessary regimes, identified from observational constraints. In the next section, we describe the formalism predicting the DM halo abundance.

3.2 The f (R) halo mass function

The density of clusters per unit mass and redshift was modeled as follows: dnd ln M=ρm,0Md ln σ−1d ln Mf(σ),${{{\rm{d}}n} \over {{\rm{d}}\ln M}} = {{{\rho _{{\rm{m}},0}}} \over M}{{{\rm{d}}\ln {\sigma ^{ - 1}}} \over {{\rm{d}}\ln M}}f(\sigma ),$(10)

where ρm,0 is the matter density at present, f (σ) is the multiplicity function, and σ(R, z) is defined by σ2(R,𝔃)=12π2∫0∞k2P(k,𝔃)(3j1(kR)kR)dk,${\sigma ^2}\left( {R,} \right) = {1 \over {2{\pi ^2}}}\mathop \smallint \limits_0^\infty {k^2}P\left( {k,} \right)\left( {{{3{j_1}\left( {kR} \right)} \over {kR}}} \right){\rm{d}}k,$(11)

where j1(x) = (sin(x) – x cos(x))/x2 is the spherical Bessel function of the first kind of order one, and P(k) is the matter power spectrum. We followed the approach developed in Hagstotz et al. (2019) to describe the specific behavior of the HMF in the HS- f(R) framework when compared to the standard fitting functions and emulators (Tinker et al. 2008; Despali et al. 2016; Bocquet et al. 2020; Euclid Collaboration 2023) obtained from simulations assuming standard GR. The main advantage of this formalism is that the entire framework is derived from the collapse of the structures in f(R) gravity, with standard computations for the power spectrum. Deviations from the GR HMF were modeled as follows: the multiplicity function can be written as in Eq. (10). The goal is to correct the multiplicity function by a factor encapsulating the deviations from GR; therefore, ff(R)(σ)=fGR(σ)⋅η(fR0,σ),${f^{f(R)}}(\sigma ) = {f^{{\rm{GR}}}}(\sigma ) \cdot \eta \left( {{f_{{\rm{R}}0}},\sigma } \right),$(12)

where the correcting factor is given by η(fR0,σ)=fkf(R)(fR0,σ)fkGR(σ).$\eta \left( {{f_{{\rm{R}}0}},\sigma } \right) = {{f_{\rm{k}}^{f(R)}\left( {{f_{{\rm{R}}0}},\sigma } \right)} \over {f_{\rm{k}}^{{\rm{GR}}}(\sigma )}}.$(13)

The multiplicity functions, including the subscript k in the correcting factor, were defined as follows. In the formalism of the multiplicity function, the density field, δ, realizes a random walk (i.e., a Markov process) when smoothed by a simple top-hat filter in k space, with the radius, R, defined in Eq. (11). This means that in the trajectories, δ(R), the next positions depend only on the current position. The multiplicity functions in Eq. (13) are obtained that way and are noted with the subscript k. They are not realistic enough to describe the halo abundance. The f(R) properties are encoded on the collapse of structures and not on the statistics of the density field, and their ratio is well suited to quantify a departure from GR. However, more realistic filters should be considered to derive a proper fitting function, which causes departures from the uncorrelated random walk. Any multiplicity function, including these non-Markovian corrections found in the literature, can be applied for fGR. In this work, we used fGR(σ) from Tinker et al. (2008). We emphasize that the resulting fitting function for HMF we have used in this work was initially calibrated for masses defined as M200m. We converted the HMF to M500c as eRASS1 mass measurements are calibrated for M500,c. Additionally, doing so allows us to stay consistent with the standard ΛCDM cosmology analysis presented in Ghirardini et al. (2024).

In the case of GR, the multiplicity function is given by fkGR(σ)=2aπδcσexp(−a(δc+βσ2)22σ2),$f_{\rm{k}}^{{\rm{GR}}}(\sigma ) = \sqrt {{{2a} \over \pi }} {{{\delta _c}} \over \sigma }\exp \left( { - a{{{{\left( {{\delta _{\rm{c}}} + \beta {\sigma ^2}} \right)}^2}} \over {2{\sigma ^2}}}} \right),$(14)

where a and β were fixed to the values found in Hagstotz et al. (2019) and provided in Table 1, and δc is the critical density for the collapse in GR given in appendix C of Nakamura & Suto (1997), δcGR(𝔃)=3(12π)2/320(1−0.012299 log(1+Ωm−1−1(1+𝔃)3)).$\delta _c^{{\rm{GR}}}() = {{3{{(12\pi )}^{2/3}}} \over {20}}\left( {1 - 0.012299\log \left( {1 + {{\Omega _m^{ - 1} - 1} \over {{{(1 + )}^3}}}} \right)} \right).$(15)

The multiplicity function in the case of HS-f (R) gravity for non-Markovian correction follows fkf(R)(σ)=2aπ1σe−aB¯2/(2σ2)(δcf(R)−3M2∂δcf(R)∂M∂ ln σ∂ ln R),$f_k^{f(R)}(\sigma ) = \sqrt {{{2a} \over \pi }} {1 \over \sigma }{{\rm{e}}^{ - a{{\bar B}^2}/\left( {2{\sigma ^2}} \right)}}\left( {\delta _{\rm{c}}^{f(R)} - {{3M} \over 2}{{\partial \delta _{\rm{c}}^{f(R)}} \over {\partial M}}{{\partial \ln \sigma } \over {\partial \ln R}}} \right),$(16)

where a is a constant. B¯$\bar B$ is the barrier, the threshold beyond which collapse structures are formed. This threshold follows B¯=δcf(R)(fR0,M,𝔃)+βσ2,$\bar B = \delta _{\rm{c}}^{f(R)}\left( {{f_{{\rm{R}}0}},M,} \right) + \beta {\sigma ^2},$

with β being another constant. The key of this formalism is thus to compute the spherical collapse threshold in f(R) gravity, noted as δcf(R)$\delta _{\rm{c}}^{f(R)}$ hereafter. We used the numerical solution provided by Kopp et al. (2013) and given in Eq. (17). For a review, we point the reader to the later reference and the previous ones in this section. The numerical solution for the spherical collapse threshold follows δcf(R)(fR0,M,𝔃)=δcGR(z)×Δ(fR0,M,𝔃),$\delta _{\rm{c}}^{f(R)}\left( {{f_{R0}},M,} \right) = \delta _{\rm{c}}^{{\rm{GR}}}(z) \times \Delta \left( {{f_{R0}},M,} \right),$(17)

where ∆ (fR0, M, z) is given by Δ(fR0,M,𝔃)=1+b2(1+𝔃)−a3(mb−mb2+1)+b3(tanh(mb)−1).$\matrix{ {{\rm{\Delta }}\left( {{f_{R0}},M,} \right) = } \hfill & {1 + {b_2}{{(1 + )}^{ - {a_3}}}\left( {{m_{\rm{b}}} - \sqrt {m_{\rm{b}}^2 + 1} } \right)} \hfill \cr {} \hfill & { + {b_3}\left( {\tanh \left( {{m_{\rm{b}}}} \right) - 1} \right).} \hfill \cr } $(18)

The different parameters of the collapse model are defined as follows, as they were originally fit: { b2=0.0166a3(fR0)=1+exp(−2.08(log| fR0 |+5.56)2)mb(fR0,M,z)=(1+z)a3(logM−m1(1+z)−a4)m1(fR0)=μ1log| fR0 |+μ2a4(fR0)=α4( tanh(0.69(log| fR0 |+6.65))  +1)b3(fR0)=β3(2.41−log| fR0 |). $\left\{ \matrix{ {b_2}\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,{\rm{ = }}\,\,{\rm{0}}{\rm{.0166}} \hfill \cr {a_{(3)}}({f_{R0}})\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = 1 + exp( - 2.08{(log|{f_{(R0)}}| + 5.56)^{(2)}}) \hfill \cr \matrix{ {{m_b}\left( {{f_{R0}},M,} \right)\,\,\,\,\, = } \hfill & {{{(1 + )}^{{a_3}}}\left( {\log M - {m_1}{{(1 + )}^{ - {a_4}}}} \right)} \hfill \cr {} \hfill & {{m_1}\left( {{f_{R0}}} \right) = {\mu _1}\log \left| {{f_{R0}}} \right| + {\mu _2}} \hfill \cr } \hfill \cr \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,{a_4}\left( {{f_{R0}}} \right) = {\alpha _4}\left( {\tanh \left( {0.69\left( {\log \left| {{f_{R0}}} \right| + 6.65} \right)} \right)} \right. \hfill \cr & & & & \,\,\,\left. {\,\,\,\, + 1} \right) \hfill \cr {b_3}\left( {{f_{R0}}} \right)\quad \,\,\,\,\,\,\,\,\,\,\,\,\,\,\, = {\beta _3}\left( {2.41 - \log \left| {{f_{R0}}} \right|} \right) \hfill \cr} \right..$(19)

Three values are then required to describe the collapse threshold in the HS- f(R) parametrization. Those are represented by the parameters µ1, µ2, α4, and β3. In this work, we have used the best-fit values given by Hagstotz et al. (2019), which are summarized in Table 1.

Figure 2 shows the expected distribution of DM halos when we vary log fR0. As was expected, we recover the GR prediction when |fR0| takes small values. At log fR0 < −5, the number of clusters is changed at most by 10%, for objects of a mass of 1014 M⊙. For objects with log fR0 < −6, the expected distribution clusters cannot be distinguished from the GR prediction. We emphasize that for log fR0 = −8, we observe an excess of very massive halos. However, this is related to modeling inaccuracies. Given the small counts of clusters in this mass range, this has no impact on our cosmological inference pipeline.

The HMF derived here was computed for M200m, while the eRASS1 cosmology pipeline considers M500 c. Thus, we followed (Hu & Kravtsov 2003) to convert the mass: we first computed the HMF from Tinker et al. (2008), dn/d ln M500c, for a broad range of halo masses, M500c, and converted these masses assuming the mass concentration relation from Duffy et al. (2008) and the Navarro, Frenk, and White (NFW) density profile (Navarro et al. 1996). The resultant f (R) HMF thus follows dnf(R)d lnM500c=ρm,0M500c| d lnσd lnM500c |fTinker2008(σ(M500c))×fkf(R)(σ(M200m(M500c)))fkGR(σ(M200m(M500c))).$\matrix{ {{{{\rm{d}}{n^{f(R)}}} \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}} = } \hfill & {{{{\rho _{{\rm{m}},0}}} \over {{M_{500{\rm{c}}}}}}\left| {{{{\rm{d}}\ln \sigma } \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}}} \right|{f_{{\rm{Tinker2008}}}}\left( {\sigma \left( {{M_{500{\rm{c}}}}} \right)} \right)} \hfill \cr {} \hfill & { \times {{f_{\rm{k}}^{f(R)}\left( {\sigma \left( {{M_{200{\rm{m}}}}\left( {{M_{500{\rm{c}}}}} \right)} \right)} \right)} \over {f_{\rm{k}}^{{\rm{GR}}}\left( {\sigma \left( {{M_{200{\rm{m}}}}\left( {{M_{500{\rm{c}}}}} \right)} \right)} \right)}}.} \hfill \cr } $(20)

This equation represents the HMF used in this work, although we note that the uncertainties related to the HMF can impact our constraints (Salvati et al. 2020; Artis et al. 2021).

Indeed, due to different numerical simulations, halo finders, various fitting functions used in previous work, and baryonic effects, we expect the departure between the standard GR HMFs in the literature to be as high as 10% (Knebe et al. 2013; Euclid Collaboration 2023). Additionally, Cataneo et al. (2015) derived constraints using a ratio of HMFs computed in GR and f(R) gravity. Their approach was conservative, as an accurate, robust, and computationally efficient modeling of the chameleon screening was lacking at the time. Therefore, it was impossible to take advantage of the full constraining power of the data. The main difference is that the authors used the multiplicity function proposed by Sheth & Tormen (2002). In the ratio of their multiplicity function, they also computed the variance, σ(R), in the respective considered theory. In other words, the term related to the f(R) multiplicity function uses the power spectrum computed in f(R) gravity. Finally, the critical density for the collapse in f(R) gravity is computed if the theory is always unscreened. These modeling choices lead to differences of up to 20% with the predictions considered in this work. As was described earlier, the correction proposed by Hagstotz et al. (2019) is more recent and takes into account more physical effects, while the model of Cataneo et al. (2016) is computationally expensive in the context of a Bayesian analysis. We thus chose to focus on the latest HMF predictions. However, these systematic uncertainties are accounted for following Costanzi et al. (2019); that is, we computed the HMF as dn˜f(R)dlnM500c=dnf(R)d lnM500c(s ln(M/1014M⊙)+q),${{{\rm{d}}{{\mathop n\limits^ }^{f(R)}}} \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}} = {{{\rm{d}}{n^{f(R)}}} \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}}\left( {s\ln \left( {M/{{10}^{14}}{M_ \odot }} \right) + q} \right),$(21)

where q and s represent the errors on the amplitude and slope of the HMF, respectively. Adopting a conservative approach, we used the priors q ~ 𝒩(1, 0.5) and s ~ 𝒩(0, 0.5). Annex A investigates the choice of these priors and their impact on the HMF.

Table 1

Parameters of the HMF used in this analysis.

thumbnail Fig. 2

Relative difference (in percentages) of the HMF at z = 0.5 when we increase fR0. We note that ∆n = dn/d ln M(fR0) – dn/d ln M, where dn/d ln M(fR0) is the HMF described with Eq. (20) and dn/d ln M is obtained with the GR fitting function from Tinker et al. (2008). Increased values of fR0 create more massive halos.

3.3 Weak lensing mass calibration in modified gravity

In this work, we aim to constrain modifications of the theory of gravity with the cluster HMF calibrated with the scaling relation between the X-ray count rate and the WL shear data. Therefore, the predictions for the WL signal of galaxy clusters in the MG theories should be needed for self-consistent analysis, but we show that GR approximation can be valid in this scenario for the eRASS1 sample. As is discussed in the following, we used the prediction from GR for the WL signal of massive clusters by simply utilizing the calibration performed by Grandis et al. (2024). For readability, we provide the basis of the analysis and critical aspects (see Bartelmann & Schneider 2001, for detailed description and review). In the framework of GR, the gravitational potential of a lens cluster deflects the light from distant background sources. The surface brightness, Ss, at position θ of a source with negligible dimensions compared to the lens size is as follows: Sobs(θ)=Ss(Aθ),${S_{{\rm{obs}}}}(\theta ) = {S_{\rm{s}}}(A\theta ),$(22)

where A is the lensing distortion matrix. It is written as A=(1−κ−γ1γ2γ21−κ+γ1).$A = \left( {\matrix{ {1 - \kappa - {\gamma _1}} & {{\gamma _2}} \cr {{\gamma _2}} & {1 - \kappa + {\gamma _1}} \cr } } \right).$(23)

The convergence, κ, relates to isotropic transformations of the image, while the shear components (γ1, γ2) are linked to the image distortion. These two quantities describe the main transformations of the image, and they are computed as follows: the total distribution of mass in clusters is assumed to follow the NFW profile introduced in Navarro et al. (1996), which is well-fit simulated DM halos, and defined as ρ(r)=δ(cΔ)ρc(𝔃)rrΔ/cΔ(1+rrΔ/cΔ)2,$\rho (r) = {{\delta \left( {{c_{\rm{\Delta }}}} \right){\rho _{\rm{c}}}()} \over {{r \over {{r_{\rm{\Delta }}}/{c_{\rm{\Delta }}}}}{{\left( {1 + {r \over {{r_{\rm{\Delta }}}/{c_{\rm{\Delta }}}}}} \right)}^2}}},$(24)

where r∆ is the radius defined such that the mean density of the enclosed volume is ∆ times the mean density of the Universe. ρc = 3H2/8πG is the critical density and c∆ is the concentration parameter. δ is a function that results from the definition of the radius, r∆, and reads: δ(cΔ)=Δ3cΔ3ln(1+cΔ)−cΔ(1+cΔ).$\delta \left( {{c_{\rm{\Delta }}}} \right) = {{\rm{\Delta }} \over 3}{{c_{\rm{\Delta }}^3} \over {\ln \left( {1 + {c_{\rm{\Delta }}}} \right) - {c_{\rm{\Delta }}}\left( {1 + {c_{\rm{\Delta }}}} \right)}}.$(25)

This mass profile creates a gravitational potential noted as ϕ. We note the lens cluster, l, and background sources (typically galaxies). In the thin lens approximation (the typical size of the lens is negligible compared to the distance), we can write the lensing potential of a cluster as ψ(θ)=2c2DlsD1Ds∫​ϕ(D1θ,𝔃)d𝔃,$\psi (\theta ) = {2 \over {{c^2}}}{{{D_{{\rm{ls}}}}} \over {{D_1}{D_{\rm{s}}}}}{\mathop \smallint \nolimits^ ^}\phi \left( {{D_1}\theta ,} \right){\rm{d}},$(26)

where θ = (θ1,θ2) is the position on the sky. Dls, Dl, and Ds are, respectively, the distance between the lens cluster and the source galaxy, between the observer and the lens, and between the observer and the galaxy. The convergence and shear component can be expressed as a function of the derivatives of the lensing potential through { κ(θ)=12(∂2ψ∂θ12+∂2ψ∂θ12)γ1(θ)=12(∂2ψ∂θ12−∂2ψ∂θ12).γ2(θ)=∂2ψ∂θ1∂θ2$\{ \matrix{ {\matrix{ {\kappa (\theta )} & { = {1 \over 2}\left( {{{{\partial ^2}\psi } \over {\partial \theta _1^2}} + {{{\partial ^2}\psi } \over {\partial \theta _1^2}}} \right)} \cr {{\gamma _1}(\theta )} & { = {1 \over 2}\left( {{{{\partial ^2}\psi } \over {\partial \theta _1^2}} - {{{\partial ^2}\psi } \over {\partial \theta _1^2}}} \right).} \cr } } \hfill \cr {{\gamma _2}(\theta ) = {{{\partial ^2}\psi } \over {\partial {\theta _1}\partial {\theta _2}}}} \hfill \cr } $(27)

These expressions, combined with the definition of the lensing potential (26), lead to the following: if we define the surface mass density as ∑ (R)=2∫0∞ρ(R,𝔃)d𝔃,$\mathop \sum \nolimits^ \left( R \right) = 2\mathop \smallint \limits_0^\infty \rho \left( {R,} \right){\rm{d}},$(28)

where the matter density follows the NFW profile (Eq. (24)), the convergence profile is expressed as κ(R)=Σ(R)Σc,$\kappa (R) = {{{\rm{\Sigma }}(R)} \over {{{\rm{\Sigma }}_{\rm{c}}}}},$(29)

and the radial dependence of the tangential shear follows γt(R)=Σ(<R)−Σ(R)Σc,${\gamma _{\rm{t}}}(R) = {{{\rm{\Sigma }}( < R) - {\rm{\Sigma }}(R)} \over {{{\rm{\Sigma }}_{\rm{c}}}}},$(30)

where Σ(< R) is he mean surface mass density in he radius, R, and Σc = c2/4πG · Ds/DlDls. These quantities were derived in Wright & Brainerd (2000). The WL mass calibration was performed through the “reduced tangential shear,” gt(R)=γ(R)1−κ(R).${g_{\rm{t}}}(R) = {{\gamma (R)} \over {1 - \kappa (R)}}.$(31)

The reduced tangential shear profiles are related to the underlying mass through P(g^t∣MWL,z^)$P\left( {{{\hat g}_t}\mid {M_{{\rm{WL}}}},\hat z} \right)$ (see Sect. 4.2). Since it is the consequence of the gravitational potential, it carries the necessary information on the underlying mass. To summarize, we assumed a standard lens equation resulting from GR to constrain our masses.

In this work, we consider one specific way of testing potential deviations: we probe the Hu-Sawicki parameterization of the f(R) gravity.

The theory of gravity impacts cluster WL in two ways:

  1. The paths of light through the cluster’s gravitational potential might be altered by the modified theory, leading to the same lensing potential (Bekenstein & Sanders 1994) causing a different distortion in the shape of background galaxies.

  2. The effective Newtonian potential of galaxy clusters in a modified theory might deviate from the one predicted in GR.

The first point is of no concern in the modification we explore. In the case of f(R) gravity, the cosmic expansion remains the same (see Sect. 3.1), but the geodesics are potentially affected. However, Schmidt (2010) showed that they remain virtually unchanged, up to a negligible factor on the order of ≲ 10−4. We can thus use the generic lensing equations in HS- f(R) gravity.

Halo matter profiles in f(R) gravity have been extensively studied by Mitchell et al. (2018); Ruan et al. (2023), focusing on the mass concentration relation and the 3D density profiles. Neither of these two observables exactly matches the 2D projected matter density maps needed to float the WL mass calibration in a cluster cosmological context (Grandis et al. 2021a). Nonetheless, the general findings of Mitchell et al. (2018) about the maximal halo mass at which halos become fully screened can still be applied. That work finds that halos with a mass larger than 10p2h−1M⊙${10^{{p_2}}}{h^{ - 1}}{{\rm{M}}_ \odot }$ are screened, and thus behave like they would in GR. The parameter is given by p2=1.50(log| fRO |−log(1+𝔃))+21.6,${p_2} = 1.50\left( {\log \left| {{f_{{\rm{RO}}}}} \right| - \log (1 + )} \right) + 21.6{\rm{,}}$(32)

where we rounded before the numerical errors reported in Mitchell et al. (2018). We evaluated this expression at z = 0.1, the lower redshift limit of our cluster sample. Given that we are working with an X-ray-selected sample, this is also the red-shift at which we reach the lowest limiting mass. For the values log10 |fR0| = 5.25, 5.5, 6, we find p2 = 13.69, 13.31, 12.56. The first value of p2 corresponds to a mass of 1013.85 M⊙. Considering that 95.4% of the objects in the cosmology sample have an estimated mass greater than this value, we can argue that the cluster one halo regime that we analyze is also unaltered in the f(R) gravity parameter space we explore.

4 Statistical inference

This section provides a brief description of the cluster cosmology inference pipeline. A detailed review can be found in the ΛCDM eRASS1 cosmology results (Ghirardini et al. 2024). The main element is that the cluster abundance statistic follows a Poisson statistic (i.e., the cluster abundance as a function of their observable is a Poisson process), coupled with WL shear measurements, calibrating the scaling relation between mass and X-ray observables. A mixture model for eliminating the contamination due to active galactic nuclei (AGNs) and background fluctuations misclassified as clusters in the eRASS1 cosmology subsample was employed to perform a precision cosmology experiment.

We express λ as the intensity of the Poisson process (the number density of objects per unit of observables) and x as the vector of the observables. This intensity depends on the model, and we note as Θ the set of cosmological parameters. By definition, the expected number of objects whose observable properties belong to the subsample, Ω, of the observable parameter space follows N{x∈Ω}(Θ)=∫Ωλ(x|Θ)dx.${N_{\left\{ {x \in \Omega } \right\}}}\left( \Theta \right) = \mathop \smallint \limits_\Omega \lambda \left( {x|\Theta } \right){\rm{d}}x.$(33)

In our model, the global observable vector is defined as x={ C^R,z^,λ^,ℋ,g+ },$x = \left\{ {{{\hat C}_{\rm{R}}},\hat z,\hat \lambda ,{\cal H},{g_ + }} \right\},$

where C^R${\hat C_{\rm{R}}}$ are the observed count rate, z^$\hat z$ are the observed photometric redshifts, λ^$\hat \lambda $ is the observed richness, ℋ is the sky positions, which need to be considered due to the uneven exposure of the survey, and 𝑔+ is the shear profile obtained from the different WL surveys described in Sect. 2.2. The intensity of the process thus describes the number density of the detected clusters in the observables parameter space. The observed value is the Poisson realization of the expected value, N{x∊Ω}(Θ). In the framework we adopt in this paper, we need to describe the intensity of the different species we are considering in order to predict the abundance of the objects. Then, the general shape of the Poisson log-likelihood will follow In ℒ(Θ)= ∑iIn(λ(xi|Θ))−∫xλ(x|Θ)dx,${\rm{In}}\,{\cal L}\left( \Theta \right) = \,\mathop \sum \limits_i {\rm{In}}\left( {\lambda \left( {{x_i}|\Theta } \right)} \right) - \mathop \smallint \limits_x \lambda \left( {x|\Theta } \right){\rm{d}}x,$(34)

where the index, i, runs over the observed objects. Using the Poisson likelihood, we neglect the cosmic variance (Hu & Kravtsov 2003) and the super sample variance (Lacasa & Grain 2019) as these effects are subdominant given the area covered by eRASS1 and the limited number of clusters included in this analysis (see Ghirardini et al. 2024, for discussion). Our catalog contains three classes of objects: galaxy clusters (C), which are of interest, and contaminants, such as AGNs, and background fluctuations misclassified as clusters (NC). Our model simultaneously accounts for the cluster counts and the contaminant fractions through the Poisson mixture model. The total density is the sum of the three-component model: λtot(x∣Θ)=λC(x∣Θ)+λAGN(x)+λNC(x).${\lambda _{{\rm{tot}}}}(x\mid {\rm{\Theta }}) = {\lambda _{\rm{C}}}(x\mid {\rm{\Theta }}) + {\lambda _{{\rm{AGN}}}}(x) + {\lambda _{NC}}(x).$(35)

In terms of the total number of objects in each class, we obtain Ntot(Θ)=NC(Θ)+fAGNNtot(Θ)+fNCNtot(Θ),${N_{{\rm{tot}}}}({\rm{\Theta }}) = {N_{\rm{C}}}({\rm{\Theta }}) + {f_{{\rm{AGN}}}}{N_{{\rm{tot}}}}({\rm{\Theta }}) + {f_{{\rm{NC}}}}{N_{{\rm{tot}}}}({\rm{\Theta }}),$(36)

where fAGN and fNC are the respective fractions of contaminants. Consequently, we can express the number of AGNs, false detections, and the total number of objects in the catalog as a function of the cosmology-dependent number of the cluster, NC : { Ntot(Θ)=(1/(1−fAGN−fNC))NC(Θ)NAGN(Θ)=(fAGN/(1−fAGN−fC))NC(Θ)NNC(Θ)=(fNC/(1−fAGN−fNC))NC(Θ) .$\left\{ {\matrix{ {{N_{{\rm{tot}}}}({\rm{\Theta }}) = \left( {1/\left( {1 - {f_{{\rm{AGN}}}} - {f_{{\rm{NC}}}}} \right)} \right){N_{\rm{C}}}({\rm{\Theta }})} \hfill \cr {{N_{{\rm{AGN}}}}({\rm{\Theta }}) = \left( {{f_{{\rm{AGN}}}}/\left( {1 - {f_{{\rm{AGN}}}} - {f_{\rm{C}}}} \right)} \right){N_{\rm{C}}}({\rm{\Theta }})} \hfill \cr {{N_{{\rm{NC}}}}({\rm{\Theta }}) = \left( {{f_{{\rm{NC}}}}/\left( {1 - {f_{{\rm{AGN}}}} - {f_{{\rm{NC}}}}} \right)} \right){N_{\rm{C}}}({\rm{\Theta }})} \hfill \cr } .} \right.$(37)

Starting from the formalism above, we can write the total number of objects as Ntot(Θ)=11−fAGN−fNC∫xλC(x∣Θ)dx.${N_{{\rm{tot}}}}({\rm{\Theta }}) = {1 \over {1 - {f_{{\rm{AGN}}}} - {f_{{\rm{NC}}}}}}\int_x {{\lambda _{\rm{C}}}(x\mid {\rm{\Theta }}){\rm{d}}x.} $(38)

The number density of contaminants follows λAGN(x|θ) = NAGN(Θ)𝒫AGN(x), and λNC(x|θ) = NNC(Θ)𝒫NC(x). In this equation, 𝒫 describes the probability distribution function of the respective object, depending on the observables. These terms are described in Sect. 4.3. Finally, the likelihood (Eq. (34)) becomes lnℒ(Θ)=∑iln( λC(xi∣Θ)                       +NAGN(Θ)PAGN(xi)+NNC(Θ)PNC(xi) )                        −11−fAGN−fNC∫xλC(x∣Θ)dx.$\eqalign{ & \ln {\cal L}(\Theta ) = \sum\limits_i {\ln } \left( {{\lambda _{\rm{C}}}\left( {{x_i}\mid \Theta } \right)} \right. \cr & \left. {\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, + {N_{{\rm{AGN}}}}(\Theta ){{\cal P}_{{\rm{AGN}}}}\left( {{x_i}} \right) + {N_{{\rm{NC}}}}(\Theta ){{\cal P}_{{\rm{NC}}}}\left( {{x_i}} \right)} \right) \cr & \,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\, - {1 \over {1 - {f_{{\rm{AGN}}}} - {f_{{\rm{NC}}}}}}\int_x {{\lambda _{\rm{C}}}} (x\mid \Theta ){\rm{d}}x. \cr} $(39)

The method described in the appendix of Bocquet et al. (2015) separates the mass calibration likelihood from the number count likelihood. For the clusters that belong to the footprints of the WL surveys, the number density becomes λC(x\{ ɡ+ },ɡ+∣Θ)=λC(x\{ ɡ+ }∣Θ)(ɡ+∣x\{ ɡ+ },Θ).${\lambda _{\rm{C}}}\left( {{x_{\backslash \left\{ {{g_ + }} \right\}}},{g_ + }\mid \Theta } \right) = {\lambda _{\rm{C}}}\left( {{x_{\backslash \left\{ {{g_ + }} \right\}}}\mid \Theta } \right){\cal P}\left( {{g_ + }\mid {x_{\left\{ {\left\{ {{g_ + }} \right\}} \right.}},\Theta } \right).$(40)

This formalism is applied in Eq. (39) and allows for the calibration of X-ray and optical scaling relations (see Sect. 4.2 and Ghirardini et al. (2024) for details of the mass calibration).

In the following sections, we briefly describe the components of Eq. (39). The cluster abundance term, λC, is described in Sect. 4.1, the WL likelihood is given in Sect. 4.2, and the mixture model accounting for contamination is provided in Sect. 4.3.

4.1 Cluster abundance likelihood

The term λC in Eq. (39) represents the expected number density of clusters. It is obtained through λC=dNC dC^R d d^ dℋ^=∫z∫M∫CRP(I∣CR,z,ℋ^)P(C^R∣CR)                                              P(CR∣M,z)P(z^∣z)dn dM dV dz dM dz dCR.$\eqalign{ & {\lambda _{\rm{C}}} = {{{\rm{d}}{N_{\rm{C}}}} \over {{\rm{d}}{{\hat C}_{\rm{R}}}{\rm{d}}\hat {\cal H}}} = \int_ {\int_M {\int_{{C_R}} {{\cal P}\left( {I\mid {C_{\rm{R}}},,\hat {\cal H}} \right){\cal P}\left( {{{\hat C}_{\rm{R}}}\mid {C_{\rm{R}}}} \right)} } } \cr & & & \,\,\,\,\,\,\,{\cal P}\left( {{C_{\rm{R}}}\mid M,} \right){\cal P}(\mid ){{{\rm{d}}n} \over {{\rm{d}}M}}{{{\rm{d}}V} \over {{\rm{d}}}}{\rm{d}}M\,{\rm{d}}{\rm{d}}{C_{\rm{R}}}. \cr} $(41)

dn/dM is the HMF, dV/dz is the comoving volume, and (I∣CR,ɀ,ℋ^)${\cal P}\left( {I\mid {C_{\rm{R}}},,\hat {\cal H}} \right)$ is the probability of detecting a cluster; that is, the selection function, the fraction of detected objects at a given observed sky position, redshift, and count rates, shown in Clerc et al. (2024). We note that the fraction of detected objects is estimated for a ΛCDM cosmology. We assume that this quantity is not significantly modified by the models considered in this work. The other terms are explained below. Consistently with (Ghirardini et al. 2024), the count rates are assumed to be linked to the underlying mass of the objects through (Ghirardini et al. 2024) 〈 lnCR¯CR,p|M,ɀ 〉=lnAX+bX(M,ɀ)·lnMMp+eX(ɀ),$\left\langle {\left. {\ln {{\overline {{C_{\rm{R}}}} } \over {{C_{{\rm{R}},{\rm{p}}}}}}} \right|{\mkern 1mu} M,} \right\rangle = \ln {A_{\rm{X}}} + {b_{\rm{X}}}(M,\left. {\ln {{\overline {{C_{\rm{R}}}} } \over {{C_{{\rm{R}},{\rm{p}}}}}}} \right|{\mkern 1mu} M,) \cdot \ln {M \over {{M_{\rm{p}}}}} + {e_{\rm{X}}}(\left. {\ln {{\overline {{C_{\rm{R}}}} } \over {{C_{{\rm{R}},{\rm{p}}}}}}} \right|{\mkern 1mu} M,),$(42)

where CR,p = 0.1 cts/s and Mp = 2 × 1014M⊙ are the pivot count rates and mass. All the masses are expressed in solar masses, M⊙. The other terms are bX(M,ɀ)=(BX+CX⋅lnMMp+FX⋅ln1+ɀ1+ɀp),${b_{\rm{X}}}(M,) = \left( {{B_{\rm{X}}} + {C_{\rm{X}}} \cdot \ln {M \over {{M_{\rm{p}}}}} + {F_{\rm{X}}} \cdot \ln {{1 + } \over {1 + {_{\rm{p}}}}}} \right),$(43)

where zp = 0.35 is the pivot redshift, and eX(ɀ)=DX⋅lndL(ɀ)dL(ɀp)+EX⋅lnE(ɀ)E(ɀp)+GX⋅ln1+ɀ1+ɀp,${e_{\rm{X}}}() = {D_{\rm{X}}} \cdot \ln {{{d_{\rm{L}}}()} \over {{d_{\rm{L}}}\left( {{_{\rm{p}}}} \right)}} + {E_{\rm{X}}} \cdot \ln {{E()} \over {E\left( {{_{\rm{p}}}} \right)}} + {G_{\rm{X}}} \cdot \ln {{1 + } \over {1 + {_{\rm{p}}}}},$(44)

where dL is the luminosity distance and E(z) = H(z)/H0 is the evolution function. The parameters {AX, BX, CX, DX, EX, FX} were fit jointly with the cosmological parameters to provide consistency and avoid astrophysical biases. We assume that the true count rates are the following, CR~Log−N(CR(M,ɀ)¯,σX),${C_{\rm{R}}}\~\log - {\cal N}\left( {\overline {{C_{\rm{R}}}(M,)} ,{\sigma _{\rm{X}}}} \right),$

where Log-𝒩 characterizes the log-normal distribution, CR(M,ɀ)¯$\overline {{C_{\rm{R}}}(M,)} $ is the output of the scaling relation in Eq. (42), and σX is the dispersion on the selection function, fit together with the cosmological parameters. We also include the measurement errors on count rates supplied by the MultiBand Projector in 2D (MBProj2D) tool (Sanders et al. 2018), and the observed count rates are C^R~Log−N(CR,σ^),${\hat C_{\rm{R}}}\~\log - {\cal N}\left( {{C_{\rm{R}}},\hat \sigma } \right),$

where σ^$\hat \sigma $ is fixed at the value given by the X-ray processing pipeline, MBProj2D (see Bulbul et al. 2024, for details). The observed redshifts were obtained through eROMaPPer, the optical identification tool (see Kluge et al. 2024), as in the following: P(ɀ^∣ɀ)=(1−cɀ)⋅N(bɀ⋅ɀ,σɀ(1+ɀ)) +cɀ⋅N(bɀ⋅ɀ+cshift ,ɀ,σɀ(1+ɀ)),$\eqalign{ & P(\hat \mid ) = \left( {1 - {c_}} \right) \cdot {\cal N}\left( {{b_} \cdot ,{\sigma _}(1 + )} \right) \cr & & + {c_} \cdot {\cal N}\left( {{b_} \cdot + {c_{{\rm{shift }},}},{\sigma _}(1 + )} \right), \cr} $(45)

where cz, bz, σz, and cshift,z are the fit parameters of the distribution. The richness, a proxy for the number of galaxies belonging to a given cluster, was computed starting from Grandis et al. (2021b) as 〈lnλ∣M,ɀ〉=lnAλ+bλ(ɀ)ln(MMp)+Cλln(1+ɀ1+ɀp),$\langle \ln \lambda \mid M,\rangle = \ln {A_\lambda } + {b_\lambda }()\ln \left( {{M \over {{M_{\rm{p}}}}}} \right) + {C_\lambda }\ln \left( {{{1 + } \over {1 + {_{\rm{p}}}}}} \right),$(46)

where all the masses are expressed in solar mass (M⊙), and the redshift dependence of the mass slope follows bλ(ɀ)=Bλ+Dλln(1+ɀ1+ɀp).${b_\lambda }() = {B_\lambda } + {D_\lambda }\ln \left( {{{1 + } \over {1 + {_{\rm{p}}}}}} \right).$(47)

The optical richness is not one of the main observables used in the mass calibration. It is only used to eliminate the contamination in the sample through the mixture method. We kept the richness limit low (λ^>3)$(\hat \lambda > 3)$ to ensure that no additional optical selection effects were introduced in the analysis.

Table 2

Priors on the cosmological parameters used in this analysis when they are involved in the corresponding model.

4.2 Weak-lensing likelihood

An essential step between the cluster number counts and X-ray observables obtained from the eROSITA survey and the HMF is calibrating the cluster mass scaling relations. We utilized the WL shear measurements in calibrating cluster scaling relations as it is currently the most reliable method with minimal bias available (Euclid Collaboration 2024). We used the surveys presented in Sect. 2 – DES, KiDS, and HSC – to achieve this goal. Following Grandis et al. (2024), we assumed a relation between WL inferred masses and true mass in the following form: 〈 lnMWLMp M,ɀ 〉=b(ɀ)+bMln(MMp).$\left\langle {\ln {{{M_{{\rm{WL}}}}} \over {{M_{\rm{p}}}}}{\rm{ }}M,} \right\rangle = b() + {b_M}\ln \left( {{M \over {{M_{\rm{p}}}}}} \right).$(48)

The DES, KiDS, and HSC surveys were used simultaneously to calibrate the count rate to mass (Eq. (42)) and richness to mass (Eq. (46)) relations. The scatter around this scaling is logσMWL2=s(ɀ)+sMlog(MMp),$\log \sigma _{{M_{{\rm{WL}}}}}^2 = s() + {s_M}\log \left( {{M \over {{M_{\rm{p}}}}}} \right),$

where the mass component is consistent with zero, and we fixed sM = 0 in the rest of the analysis, as is prescribed by Grandis et al. (2024).

4.3 Contamination fraction

Although the eRASS1 cosmology subsample is a high-fidelity catalog with a purity level reaching 95%, the remaining 5% of contaminants should be eliminated from our analysis. We applied a mixture model fitting for the fraction of AGN and random sources. We refer the reader to Ghirardini et al. (2024) for a detailed method review.

5 Results

This section provides our results for the HS- f (R) gravity based on cluster counts from the eRASS1 cosmology catalog with WL mass calibration. We first show that our mass scale is consistent with the base ΛCDM analysis from Ghirardini et al. (2024). The priors on the cosmological parameters are shown in Table 2.

5.1 Comparison of the mass scale with GR-ΛCDM

One of the key ingredients of cluster abundance cosmology is the relation between the observables used to detect the objects and the underlying cluster mass. In this work, the main observables are the X-ray count rates of the clusters observed in the first eROSITA All-Sky Survey. The modeling of the scaling relations between observables and the cluster mass is described in Sect. 4. The parameters of this count rate to mass relation are fit jointly with the cosmological parameters. This means that when considering extensions to GR that do not affect cluster masses due to screening mechanisms, we should find consistent scaling relations, as the fit cosmology is similar to the one found in ΛCDM.

The redshift dependence of the mass slope, CX, was fixed to zero, the luminosity distance dependence, DX, to -2, and normalized Hubble parameter dependence, EX, to 2. Figure 3 shows the fit parameters, AX, BX, FX, and GX, when compared to the standard ΛCDM analysis presented in Ghirardini et al. (2024) and the WL mass calibration part of the likelihood presented in Grandis et al. (2024). The latter represents the information from the WL mass calibration of DES clusters.

Overall, there is good agreement between our different models whenever we consider deviation from GR. The mass calibration parameters are indeed not affecting our results. This means that the mass calibration is robust and that changes in the cosmological model primarily affect the underlying distribution of clusters represented by the HMF.

thumbnail Fig. 3

Comparison of the posteriors on the scaling relation obtained in the f (R) framework (the mass of the neutrinos is fixed) and the standard ΛCDM analysis presented in Ghirardini et al. (2024) (gray). The different models are fully consistent with one another. We overplot the contours obtained from WL only from Grandis et al. (2024) (green).
Table 3

Constraints on the HS- f (R) gravity model in different scenarios.

thumbnail Fig. 4

Comparison of the posteriors obtained on the Ω – σ8 plane. We represent the ΛCDM analysis introduced in Ghirardini et al. (2024) in gray, while the posteriors of HS- f (R) case are shown in red.

5.2 Constraints on f(R) gravity

We applied two different models to constrain the HS- f (R) gravity parameterization using the Hagstotz et al. (2019) formalism. In the first one, we assumed massless neutrinos, while in the second one, we simultaneously constrained the sum of the neutrino masses and the parameter log |fR0| (Eq. (3)). We summarize our results in Table 3. In both cases, eRASS1 cluster number counts provide upper limits on log |fR0|. In the first case, we find the following constraints, summarized in Figure 4: Ωm=0.28±0.02σ8   =0.89±0.04,S8   =0.85±0.04$\eqalign{ & {\Omega _{\rm{m}}} = 0.28 \pm 0.02 \cr & {\sigma _8}\,\,\, = 0.89 \pm 0.04, \cr & {S_8}\,\,\, = 0.85 \pm 0.04 \cr} $(49)

at a 68% confidence level. These constraints are in good agreement with the results of the standard ΛCDM analysis presented in Ghirardini et al. (2024). For the first time, we have obtained constraints on the f (R) parameters with cluster abundance only, obtaining log| fR0 |<−4.31 at 95% confidence (∑ mv=0 eV).$\log \left| {{f_{{\rm{R}}0}}} \right| < - 4.31{\rm{ at }}95\% {\rm{ confidence }}\left( {\sum {{m_v}} = 0\,{\rm{eV}}} \right){\rm{. }}$(50)

Figure 5 presents the measurements in the Ωm – log |fR0| plane. Our constraints show a consistency with GR, as has been suggested by previous surveys using cluster abundance to probe HS- f (R) gravity (see Sect. 5.3). Our results suggest that large unscreened departures from GR can be ruled out by cluster counts alone in the near future, providing an independent test of gravity on scales ~10 Mpc. Additionally, Hagstotz et al. (2019) show that massive neutrinos can counteract the effect of MG by reducing the abundance of massive structures. For instance, massive neutrinos increase the collapse barrier presented in Eq. (17). For Σmν = 0.3 eV, the collapse barrier is not distinguishable from the GR barrier if log |fR0| = –5. It is thus necessary to investigate the case in which the sum of the mass of the neutrinos is fit jointly with log |fR0|. We used the prescription provided in Costanzi et al. (2013) and derived the linear quantities from the cold DM power spectrum instead of the total matter power spectrum. When we freed the masses of the neutrinos, we obtained log| fR0 |<−4.08 at 95% confidence (∑ mv free ).$\log \left| {{f_{{\rm{R}}0}}} \right| < - 4.08{\rm{ at }}95\% {\rm{ confidence }}\left( {\sum {{m_v}} {\rm{ free }}} \right).$(51)

We additionally find the following constraints on the cosmological parameters, including the neutrino masses: Ωm=0.31±0.03σ8=0.87±0.05S8=0.84±0.04∑ mv<0.49eV at 95% confidence .$\matrix{ {{\Omega _{\rm{m}}}} \hfill & { = 0.31 \pm 0.03} \hfill \cr{{\sigma _8}} \hfill & { = 0.87 \pm 0.05} \hfill \cr{{S_8}} \hfill & { = 0.84 \pm 0.04} \hfill \cr{\sum {{m_v}} } \hfill & { < 0.49{\rm{eV at }}95\% {\rm{ confidence }}} \hfill \cr } .$(52)

These constraints are again fully consistent with those obtained in the base ΛCDM (Ghirardini et al. 2024) (see Table 3). We investigated the consistency between our S8 measurements and the one given in the standard ΛCDM analysis in Ghirardini et al. (2024), who finds S8=0.86±0.01.${S_8} = 0.86 \pm 0.01.$(53)

Figure 6 shows the different models, compared to the results obtained by Planck Collaboration VI (2020). In all cases, our S8 value is in good agreement with the CMB. Additionally, the three models are fully compatible. However, the upper limits on the neutrino mass are higher due to the increased free parameters in the fits. It is worth noting that we solely used cluster abundances with no priors on the cosmological parameters when constraining fR0 . The ability to find upper limits with cluster counts alone can primarily be explained by the increased statistics provided by the eRASS1 cluster sample. This allows us to provide a self- consistent physical framework. Moreover, the constraints benefit from the shear measurement from high-quality DES, KIDS, and HSC mass calibration employed in this work (Grandis et al. 2024). It is therefore clear that in both cases, whether or not the neutrino masses are allowed to be free, the constraints from the eRASS1 cluster count measurements remain consistent with GR with a cosmological constant, for which fR0 = 0, as is shown in Eq. (3).

thumbnail Fig. 5

Constraints obtained on log fR0 with eRASS1 clusters. No external cosmological datasets are combined with our data, which is shown in this figure. Constraints with a fixed neutrino mass are shown in red, while constraints for which the neutrino mass is free are shown in magenta.
thumbnail Fig. 6

Comparison of the S 8 parameter obtained using the different models. We show the 1σ uncertainties. The case introducing massless neutrinos is in red, while the one considering massive neutrinos is in pink. Results from the CMB are shown in green. The gray area represents the 68% and 95% confidence regions of the ΛCDM analysis from Ghirardini et al. (2024).

5.3 Comparison with other constraints

Hereafter, all the constraints are shown at the 95% confidence level. In this work, we provide constraints on f (R) gravity with cluster abundance only for the first time. This allows us to test gravity with a self-consistent physical framework at a given physical scale and redshift range. Indeed, recent forecasts (Mitchell et al. 2021; Vogt et al. 2024) predict that future cluster abundance studies alone will be able to distinguish GR from alternative scenarios. Similarly, an alternative way of probing cosmology with the peaks of the matter density field is to consider WL peaks. Those also show promising results (Davies et al. 2024). This section compares our results with the ones obtained by other surveys and studies, from the smallest scales to the cosmological constraints. Indeed, due to the unknown nature of the potential modification of GR, it is important to test the theory at all physical scales.

Due to the nature of chameleon screening, the most stringent constraints are obtained from dwarf galaxies in cosmic voids. Desmond & Ferreira (2020) found |fR0| < –7.85 using galaxy morphology, which practically rules out the Hu-Sawicki model of f (R) gravity at astrophysical scales. Despite these results, testing the theory at larger scales is still relevant.

In the case of galaxy clusters, previous studies obtained consistent results with those presented in this work, but combining cluster abundances with the CMB and other LSS probes (Schmidt et al. 2009; Lombriser et al. 2012). When combined with CMB data, the most recent study by Cataneo et al. (2015) obtained log |fR0| < –4.79. We reach a similar precision with clusters only, primarily due to the statistical power of the eRASS1 cosmology sample, the largest ICM-selected sample to date.

Other probes at the cosmological scale typically obtain weaker constraints than the ones probing the small scales. They include the WL 2pt statistic. Tröster et al. (2021), using the KIDS 3 × 2 pt measurements of cosmic shear from KiDS-1000, did not find significant constraints on f (R) gravity using their data. Conversely, Liu et al. (2016) obtained log(fR0) < –5.15 combining measurements of the WL peaks with CMB priors. Hu et al. (2016) used supernovae data from JLA, WiggleZ, and CFHTLenS, in combination with CMB and BAO data, to find log(fR0) < –3.2. More recently, Kou et al. (2024) found an upper limit of log |fR0| < –4.61, cross-correlating CMB lensing with galaxy clustering.

Consequently, our results are consistent with the latest cosmological constraints. Astronomical probes on astrophysical and cosmological scales have found no evidence of a deviation from GR. Recent forecasts established for the Euclid survey (Casas et al. 2023) show that next-generation surveys will have the ability to distinguish GR from f (R) gravity with extreme precision, potentially ruling out the model definitively.

6 Discussion and conclusions

In this work, we present results on potential deviations from standard GR and ΛCDM using X-ray detected eRASS1 clusters, in combination with LS DR10-South for optical confirmation and redshift measurements, and DES, KiDS, and HSC for WL mass calibration. We have derived constraints on the HS model of f (R) gravity. We have considered both the case of massive and massless neutrinos.

All the models we have studied are statistically consistent with GR. Given the constraining power of eRASS1, we have chosen to present constraints from the cluster abundance only. This choice allows us to measure our parameters with a consistent physical framework for the first time. Indeed, this work shows that cluster counts alone have a great potential to constrain MG at the scale at which they are sensitive. Deeper surveys with eROSITA will significantly increase the constraints that we obtain here and might rival the results obtained at smaller scales. Overall, combined with the follow-up WL observations, the first eROSITA’s All-Sky survey has the statistical power to constrain the cosmological parameters with percent-level precision and to test the GR at large scales.

Lastly, we have investigated the consistency between our S8 measurements and the one given in the standard ΛCDM analysis and found good agreement. Modifying gravitation at the scale covered by cluster abundance does not significantly affect the S8 measurement. Our error bars are increased, but all posteriors are fully compatible. This emphasizes that the S 8 measurement provided in the ΛCDM analysis is robust against modifications and in good agreement with Planck Collaboration VI (2020).

Finally, the precision of our analysis mostly relies on the available HMF models. We have used the framework proposed by Hagstotz et al. (2019), while others like Cataneo et al. (2016) propose alternative approaches. Additional work will be necessary to increase the precision of our constraints by improving the reliability of the HMF, not only for f (R) gravity but also for any additional model of interest.

eROSITA, thus far, has collected data for more than four all-sky surveys. The final eRASS will detect about 105 galaxy clusters and groups (projected from the eFEDS results in Liu et al. 2022; Bulbul et al. 2022). The deeper eROSITA data and extensive cluster catalogs, in combination with state-of-the-art WL mass calibration, will allow us to tighten our constraints on the f (R) gravity models.

Acknowledgement

The full acknowledgments are available in Appendix B.

Appendix A Impact of the errors on the halo mass function

In this work, we use correcting factors for the HMF errors. Indeed, several methods yield different results for predicting the HMF in the HS- f (R) scenario. We use the prediction of Hagstotz et al. (2019), but alternative approaches exist like the one introduced by Cataneo et al. (2016). These different approaches yield up to 20% differences in the prediction of the HMF.

We use the analytic form introduced in Costanzi et al. (2019), where dn˜f(R)dlnM500c=dnf(R)dlnM500c(sln(M/1014M⊙)+q),${{{\rm{d}}{{\tilde n}^{f(R)}}} \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}} = {{{\rm{d}}{n^{f(R)}}} \over {{\rm{d}}\ln {M_{500{\rm{c}}}}}}\left( {s\ln \left( {M/{{10}^{14}}{M_ \odot }} \right) + q} \right),$

We compare the posterior distribution of our parameters s and q and their impact on the HMF. This is shown in figure A.1. We fix the cosmology to the results of Planck Collaboration VI (2020). Then, we compute the HMF for GR (using the fit provided by Tinker et al. (2008)) and for different values of fR0. Assuming that the GR HMF is correct, we then use the posterior distribution of the parameters s and q to quantify their impact on the HMF. We show that including the slope and amplitude correction of the HMF corrects most of the expected difference between the different models of f (R), i.e., 20% difference. In the future, a more comprehensive modeling of the HMF uncertainties will be necessary to obtain high-precision measurements on the MG parameters.

thumbnail Fig. A.1

Relative difference of the HMF in different MG scenarios. ∆n = dn/dln M(fR0) – dn/d ln M, where dn/d ln M(fR0) is the HMF described with Eq. 20 and dn/d ln M is obtained with the GR fitting function from Tinker et al. (2008). We use the posterior distribution of the correcting parameters s and q and compare their impact to one of the MG parameters. The dark gray area corresponds to the 68% confidence level area, while the light gray one corresponds to the 95% confidence area.

Appendix B Acknowledgements

The authors thank Prof. Catherine Heymans for her help and availability for the implementation of the blinding strategy and her helpful comments on the draft.

This work is based on data from eROSITA, the soft X-ray instrument aboard SRG, a joint Russian-German science mission supported by the Russian Space Agency (Roskosmos), in the interests of the Russian Academy of Sciences represented by its Space Research Institute (IKI), and the Deutsches Zentrum für Luft und Raumfahrt (DLR). The SRG spacecraft was built by the Lavochkin Association (NPOL) and its subcontractors and is operated by NPOL with support from the Max Planck Institute for Extraterrestrial Physics (MPE).

The development and construction of the eROSITA X-ray instrument was led by MPE, with contributions from the Dr. Karl Remeis Observatory Bamberg & ECAP (FAU Erlangen- Nuernberg), the University of Hamburg Observatory, the Leibniz Institute for Astrophysics Potsdam (AIP), and the Institute for Astronomy and Astrophysics of the University of Tübingen, with the support of DLR and the Max Planck Society. The Arge- lander Institute for Astronomy of the University of Bonn and the Ludwig Maximilians Universität Munich also participated in the science preparation for eROSITA.

The eROSITA data shown here were processed using the eSASS software system developed by the German eROSITA consortium (Brunner et al. 2022).

V. Ghirardini, E. Bulbul, A. Liu, C. Garrel, S. Zelmer, and X. Zhang acknowledge financial support from the European Research Council (ERC) Consolidator Grant under the European Union’s Horizon 2020 research and innovation program (grant agreement CoG DarkQuest No 101002585). CNES financially supported N. Clerc. M. Cataneo acknowledges the financial support provided by the Alexander von Humboldt Foundation through the Humboldt Research Fellowship program, as well as support from the Max Planck Society and the Alexander von Humboldt Foundation in the framework of the Max Planck-Humboldt Research Award endowed by the Federal Ministry of Education and Research. T. Schrabback and F. Kleinebreil acknowledge support from the German Federal Ministry for Economic Affairs and Energy (BMWi) provided through DLR under projects 50OR2002, 50OR2106, and 50OR2302, as well as the support provided by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under grant 415537506. The Innsbruck group also acknowledges support provided by the Austrian Research Promotion Agency (FFG) and the Federal Ministry of the Republic of Austria for Climate Action, Environment, Energy, Mobility, Innovation and Technology (BMK) via the Austrian Space Applications Programme with grant numbers 899537 and 900565.

The Legacy Surveys consist of three individual and complementary projects: the Dark Energy Camera Legacy Survey (DECaLS; Proposal ID #2014B-0404; PIs: David Schlegel and Arjun Dey), the Beijing-Arizona Sky Survey (BASS; NOAO Prop. ID #2015A-0801; PIs: Zhou Xu and Xiaohui Fan), and the Mayall z-band Legacy Survey (MzLS; Prop. ID #2016A-0453; PI: Arjun Dey). DECaLS, BASS, and MzLS together include data obtained, respectively, at the Blanco telescope, Cerro Tololo Inter-American Observatory, NSF’s NOIRLab; the Bok telescope, Steward Observatory, University of Arizona; and the Mayall telescope, Kitt Peak National Observatory, NOIRLab. Pipeline processing and analyses of the data were supported by NOIRLab and the Lawrence Berkeley National Laboratory (LBNL). The Legacy Surveys project is honored to be permitted to conduct astronomical research on Iolkam Du’ag (Kitt Peak), a mountain with particular significance to the Tohono O’odham Nation.

Funding for the DES Projects has been provided by the U.S. Department of Energy, the U.S. National Science Foundation, the Ministry of Science and Education of Spain, the Science and Technology FacilitiesCouncil of the United Kingdom, the Higher Education Funding Council for England, the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign, the Kavli Institute of Cosmological Physics at the University of Chicago, the Center for Cosmology and Astro-Particle Physics at the Ohio State University, the Mitchell Institute for Fundamental Physics and Astronomy at Texas A&M University, Financiadora de Estudos e Proje- tos, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvi- mento Científico e Tecnológico and the Ministério da Ciência, Tecnologia e Inovação, the Deutsche Forschungsgemeinschaft, and the Collaborating Institutions in the Dark Energy Survey.

The Collaborating Institutions are Argonne National Laboratory, the University of California at Santa Cruz, the University of Cambridge, Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas-Madrid, the University of Chicago, University College London, the DES-Brazil Consortium, the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH) Zürich, Fermi National Accelerator Laboratory, the University of Illinois at Urbana-Champaign, the Institut de Ciències de l’Espai (IEEC/CSIC), the Institut de Física d’Altes Energies, Lawrence Berkeley National Laboratory, the Ludwig-Maximilians Universität München and the associated Excellence Cluster Universe, the University of Michigan, the National Optical Astronomy Observatory, the University of Nottingham, The Ohio State University, the OzDES Membership Consortium, the University of Pennsylvania, the University of Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the University of Sussex, and Texas A&M University.

Based on observations made with ESO Telescopes at the La Silla Paranal Observatory under program IDs 177.A-3016, 177.A-3017, 177.A-3018, and 179.A-2004, and on data products produced by the KiDS consortium. The KiDS production team acknowledges support from Deutsche Forschungsgemeinschaft, ERC, NOVA, and NWO-M grants; Target; the University of Padova; and the University Federico II (Naples).

This paper makes use of software developed for the Large Synoptic Survey Telescope. We thank the LSST Project for making their code available as free software at http://dm.lsst.org

This work made use of the following Python software packages: SciPy1 (Virtanen et al. 2020), Matplotlib2 (Hunter 2007), Astropy3 (Astropy Collaboration 2022), NumPy4 (Harris et al. 2020), CAMB (Lewis & Challinor 2011), pyCCL5 (Chisari et al. 2019), GPy6 (GPy 2012), climin7 (Bayer et al. 2015), ultranest8 (Buchner 2021)

References

  1. Aihara, H., Armstrong, R., Bickerton, S., et al. 2018, PASJ, 70, S8 [NASA ADS] [Google Scholar]
  2. Alam, S., Ata, M., Bailey, S., et al. 2017, MNRAS, 470, 2617 [Google Scholar]
  3. Amon, A., Gruen, D., Troxel, M. A., et al. 2022, Phys. Rev. D, 105, 023514 [NASA ADS] [CrossRef] [Google Scholar]
  4. Artis, E., Melin, J.-B., Bartlett, J. G., & Murray, C. 2021, A&A, 649, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  5. Asgari, M., Lin, C.-A., Joachimi, B., et al. 2021, A&A, 645, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Astropy Collaboration (Price-Whelan, A. M., et al.) 2022, ApJ, 935, 167 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bartelmann, M., & Schneider, P. 2001, Phys. Rep., 340, 291 [Google Scholar]
  8. Bayer, J., Osendorfer, C., Diot-Girard, S., Rueckstiess, T., & Urban, S. 2015, TUM, Tech. Rep. [Google Scholar]
  9. Bekenstein, J. D., & Sanders, R. H. 1994, ApJ, 429, 480 [NASA ADS] [CrossRef] [Google Scholar]
  10. Bennett, C. L., Larson, D., Weiland, J. L., et al. 2013, ApJS, 208, 20 [Google Scholar]
  11. Bocquet, S., Saro, A., Mohr, J. J., et al. 2015, ApJ, 799, 214 [NASA ADS] [CrossRef] [Google Scholar]
  12. Bocquet, S., Dietrich, J. P., Schrabback, T., et al. 2019, ApJ, 878, 55 [Google Scholar]
  13. Bocquet, S., Heitmann, K., Habib, S., et al. 2020, ApJ, 901, 5 [Google Scholar]
  14. Bocquet, S., Grandis, S., Bleem, L. E., et al. 2024a, Phys. Rev. D, submitted [arXiv:2310.12213v2] [Google Scholar]
  15. Bocquet, S., Grandis, S., Bleem, L. E., et al. 2024b, Phys. Rev. D, submitted [arXiv:2401.02075] [Google Scholar]
  16. Brunner, H., Liu, T., Lamer, G., et al. 2022, A&A, 661, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Buchner, J. 2021, The Journal of Open Source Software, 6, 3001 [NASA ADS] [CrossRef] [Google Scholar]
  18. Bulbul, E., Liu, A., Pasini, T., et al. 2022, A&A, 661, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Bulbul, E., Liu, A., Kluge, M., et al. 2024, A&A, 685, A106 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Casas, S., Cardone, V. F., Sapone, D., et al. 2023, A&A, submitted [arXiv:2306.11053] [Google Scholar]
  21. Cataneo, M., Rapetti, D., Schmidt, F., et al. 2015, Phys. Rev. D, 92, 044009 [NASA ADS] [CrossRef] [Google Scholar]
  22. Cataneo, M., Rapetti, D., Lombriser, L., & Li, B. 2016, J. Cosmology Astropart. Phys., 2016, 024 [CrossRef] [Google Scholar]
  23. Chisari, N. E., Alonso, D., Krause, E., et al. 2019, CCL: Core Cosmology Library, Astrophysics Source Code Library [record ascl:1901.003] [Google Scholar]
  24. Clerc, N., & Finoguenov, A. 2023, in Handbook of X-ray and Gamma-ray Astrophysics, eds. C. Bambi, & A. Santangelo, 123 [Google Scholar]
  25. Clerc, N., Comparat, J., Seppi, R., et al. 2024, A&A, 687, A238 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Costanzi, M., Villaescusa-Navarro, F., Viel, M., et al. 2013, J. Cosmology Astropart. Phys., 2013, 012 [CrossRef] [Google Scholar]
  27. Costanzi, M., Rozo, E., Simet, M., et al. 2019, MNRAS, 488, 4779 [NASA ADS] [CrossRef] [Google Scholar]
  28. Davies, C. T., Harnois-Déraps, J., Li, B., et al. 2024, MNRAS, 533, 3546 [NASA ADS] [CrossRef] [Google Scholar]
  29. De Felice, A., & Tsujikawa, S. 2010, Liv. Rev. Relativ., 13, 3 [NASA ADS] [CrossRef] [Google Scholar]
  30. Desmond, H., & Ferreira, P. G. 2020, Phys. Rev. D, 102, 104060 [NASA ADS] [CrossRef] [Google Scholar]
  31. Despali, G., Giocoli, C., Angulo, R. E., et al. 2016, MNRAS, 456, 2486 [NASA ADS] [CrossRef] [Google Scholar]
  32. Duffy, A. R., Schaye, J., Kay, S. T., & Dalla Vecchia, C. 2008, MNRAS, 390, L64 [Google Scholar]
  33. Eisenstein, D. J., Zehavi, I., Hogg, D. W., et al. 2005, ApJ, 633, 560 [Google Scholar]
  34. Euclid Collaboration (Castro, T., et al.) 2023, A&A, 671, A100 [CrossRef] [EDP Sciences] [Google Scholar]
  35. Euclid Collaboration (Giocoli, C., et al.) 2024, A&A, 681, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Freedman, W. L., Madore, B. F., Hatt, D., et al. 2019, ApJ, 882, 34 [Google Scholar]
  37. Fumagalli, A., Costanzi, M., Saro, A., Castro, T., & Borgani, S. 2024, A&A, 682, A148 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  38. Gannouji, R., Moraes, B., & Polarski, D. 2009, J. Cosmology Astropart. Phys., 2009, 034 [CrossRef] [Google Scholar]
  39. Garrel, C., Pierre, M., Valageas, P., et al. 2022, A&A, 663, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  40. Gatti, M., Sheldon, E., Amon, A., et al. 2021, MNRAS, 504, 4312 [NASA ADS] [CrossRef] [Google Scholar]
  41. Ghirardini, V., Bulbul, E., Artis, E., et al. 2024, A&A, 689, A298 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Giblin, B., Heymans, C., Asgari, M., et al. 2021, A&A, 645, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Giocoli, C., Baldi, M., & Moscardini, L. 2018, MNRAS, 481, 2813 [Google Scholar]
  44. GPy 2012, GPy: A Gaussian process framework in python, http://github.com/SheffieldML/GPy [Google Scholar]
  45. Grandis, S., Bocquet, S., Mohr, J. J., Klein, M., & Dolag, K. 2021a, MNRAS, 507, 5671 [NASA ADS] [CrossRef] [Google Scholar]
  46. Grandis, S., Mohr, J. J., Costanzi, M., et al. 2021b, MNRAS, 504, 1253 [NASA ADS] [CrossRef] [Google Scholar]
  47. Grandis, S., Ghirardini, V., Bocquet, S., et al. 2024, A&A, 687, A178 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Gupta, S., Hellwing, W. A., Bilicki, M., & García-Farieta, J. E. 2022, Phys. Rev. D, 105, 043538 [NASA ADS] [CrossRef] [Google Scholar]
  49. Hagstotz, S., Costanzi, M., Baldi, M., & Weller, J. 2019, MNRAS, 486, 3927 [Google Scholar]
  50. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  51. Hildebrandt, H., van den Busch, J. L., Wright, A. H., et al. 2021, A&A, 647, A124 [EDP Sciences] [Google Scholar]
  52. Hu, W., & Kravtsov, A. V. 2003, ApJ, 584, 702 [Google Scholar]
  53. Hu, W., & Sawicki, I. 2007, Phys. Rev. D, 76, 064004 [CrossRef] [Google Scholar]
  54. Hu, B., Raveri, M., Frusciante, N., & Silvestri, A. 2014, Phys. Rev. D, 89, 103530 [Google Scholar]
  55. Hu, B., Raveri, M., Rizzato, M., & Silvestri, A. 2016, MNRAS, 459, 3880 [Google Scholar]
  56. Huff, E., & Mandelbaum, R. 2017, arXiv e-prints [arXiv:1702.02600] [Google Scholar]
  57. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [NASA ADS] [CrossRef] [Google Scholar]
  58. Ider Chitham, J., Comparat, J., Finoguenov, A., et al. 2020, MNRAS, 499, 4768 [NASA ADS] [CrossRef] [Google Scholar]
  59. Khoury, J., & Weltman, A. 2004, Phys. Rev. D, 69, 044026 [CrossRef] [Google Scholar]
  60. Kluge, M., Comparat, J., Liu, A., et al. 2024, A&A, 688, A210 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  61. Knebe, A., Pearce, F. R., Lux, H., et al. 2013, MNRAS, 435, 1618 [NASA ADS] [CrossRef] [Google Scholar]
  62. Kopp, M., Appleby, S. A., Achitouv, I., & Weller, J. 2013, Phys. Rev. D, 88, 084015 [NASA ADS] [CrossRef] [Google Scholar]
  63. Kou, R., Murray, C., & Bartlett, J. G. 2024, A&A, 686, A193 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Kuijken, K., Heymans, C., Dvornik, A., et al. 2019, A&A, 625, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  65. Lacasa, F., & Grain, J. 2019, A&A, 624, A61 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  66. Lesci, G. F., Nanni, L., Marulli, F., et al. 2022, A&A, 665, A100 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Lesgourgues, J. 2011, arXiv e-prints [arXiv:1104.2932] [Google Scholar]
  68. Lewis, A., & Challinor, A. 2011, CAMB: Code for Anisotropies in the Microwave Background, Astrophysics Source Code Library [record ascl:1102.026] [Google Scholar]
  69. Li, B., Zhao, G.-B., & Koyama, K. 2012, MNRAS, 421, 3481 [NASA ADS] [CrossRef] [Google Scholar]
  70. Li, X., Miyatake, H., Luo, W., et al. 2022, PASJ, 74, 421 [NASA ADS] [CrossRef] [Google Scholar]
  71. Li, X., Zhang, T., Sugiyama, S., et al. 2023, Phys. Rev. D, 108, 123518 [CrossRef] [Google Scholar]
  72. Liu, X., Li, B., Zhao, G.-B., et al. 2016, Phys. Rev. Lett., 117, 051101 [NASA ADS] [CrossRef] [Google Scholar]
  73. Liu, A., Bulbul, E., Ghirardini, V., et al. 2022, A&A, 661, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Lombriser, L., Slosar, A., Seljak, U., & Hu, W. 2012, Phys. Rev. D, 85, 124038 [NASA ADS] [CrossRef] [Google Scholar]
  75. Lombriser, L., Li, B., Koyama, K., & Zhao, G.-B. 2013, Phys. Rev. D, 87, 123511 [NASA ADS] [CrossRef] [Google Scholar]
  76. Mantz, A. B., von der Linden, A., Allen, S. W., et al. 2015, MNRAS, 446, 2205 [Google Scholar]
  77. Merloni, A., Lamer, G., Liu, T., et al. 2024, A&A, 682, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  78. Mitchell, M. A., He, J.-h., Arnold, C., & Li, B. 2018, MNRAS, 477, 1133 [NASA ADS] [CrossRef] [Google Scholar]
  79. Mitchell, M. A., Arnold, C., & Li, B. 2021, MNRAS, 508, 4157 [NASA ADS] [CrossRef] [Google Scholar]
  80. Nakamura, T. T., & Suto, Y. 1997, Progr. Theor. Phys., 97, 49 [NASA ADS] [CrossRef] [Google Scholar]
  81. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1996, ApJ, 462, 563 [Google Scholar]
  82. Oyaizu, H. 2008, Phys. Rev. D, 78, 123523 [NASA ADS] [CrossRef] [Google Scholar]
  83. Percival, W. J., & White, M. 2009, MNRAS, 393, 297 [NASA ADS] [CrossRef] [Google Scholar]
  84. Perlmutter, S., Aldering, G., Goldhaber, G., et al. 1999, ApJ, 517, 565 [Google Scholar]
  85. Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  86. Predehl, P., Andritschke, R., Arefiev, V., et al. 2021, A&A, 647, A1 [EDP Sciences] [Google Scholar]
  87. Preston, C., Amon, A., & Efstathiou, G. 2023, MNRAS, 525, 5554 [NASA ADS] [CrossRef] [Google Scholar]
  88. Riess, A. G., Filippenko, A. V., Challis, P., et al. 1998, AJ, 116, 1009 [Google Scholar]
  89. Riess, A. G., Casertano, S., Yuan, W., Macri, L. M., & Scolnic, D. 2019, ApJ, 876, 85 [Google Scholar]
  90. Ruan, C.-Z., Cuesta-Lazaro, C., Eggemeier, A., et al. 2023, MNRAS, 527, 2490 [NASA ADS] [CrossRef] [Google Scholar]
  91. Salvati, L., Douspis, M., & Aghanim, N. 2020, A&A, 643, A20 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  92. Sanders, J. S., Fabian, A. C., Russell, H. R., & Walker, S. A. 2018, MNRAS, 474, 1065 [NASA ADS] [CrossRef] [Google Scholar]
  93. Schmidt, F. 2010, Phys. Rev. D, 81, 103002 [NASA ADS] [CrossRef] [Google Scholar]
  94. Schmidt, F., Vikhlinin, A., & Hu, W. 2009, Phys. Rev. D, 80, 083505 [NASA ADS] [CrossRef] [Google Scholar]
  95. Seppi, R., Comparat, J., Ghirardini, V., et al. 2024, A&A, 686, A196 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  96. Sheldon, E. S., & Huff, E. M. 2017, ApJ, 841, 24 [NASA ADS] [CrossRef] [Google Scholar]
  97. Sheth, R. K., & Tormen, G. 2002, MNRAS, 329, 61 [Google Scholar]
  98. Sunayama, T., Miyatake, H., Sugiyama, S., et al. 2023, arXiv e-prints [arXiv:2309.13025] [Google Scholar]
  99. Tinker, J., Kravtsov, A. V., Klypin, A., et al. 2008, ApJ, 688, 709 [Google Scholar]
  100. Tröster, T., Asgari, M., Blake, C., et al. 2021, A&A, 649, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  101. Troxel, M. A., MacCrann, N., Zuntz, J., et al. 2018, Phys. Rev. D, 98, 043528 [Google Scholar]
  102. Tsujikawa, S., Gannouji, R., Moraes, B., & Polarski, D. 2009, Phys. Rev. D, 80, 084044 [NASA ADS] [CrossRef] [Google Scholar]
  103. Vikhlinin, A., Burenin, R. A., Ebeling, H., et al. 2009, ApJ, 692, 1033 [Google Scholar]
  104. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261 [CrossRef] [Google Scholar]
  105. Vogt, S. M. L., Bocquet, S., Davies, C. T., Mohr, J. J., & Schmidt, F. 2024, Phys. Rev. D, 109, 123503 [NASA ADS] [CrossRef] [Google Scholar]
  106. von Braun-Bates, F., Winther, H. A., Alonso, D., & Devriendt, J. 2017, J. Cosmology Astropart. Phys., 2017, 012 [CrossRef] [Google Scholar]
  107. Wang, Z., Mirpoorian, S. H., Pogosian, L., Silvestri, A., & Zhao, G.-B. 2023, J. Cosmology Astropart. Phys., 2023, 038 [CrossRef] [Google Scholar]
  108. Wolf, W. J., & Ferreira, P. G. 2023, Phys. Rev. D, 108, 103519 [NASA ADS] [CrossRef] [Google Scholar]
  109. Wright, C. O., & Brainerd, T. G. 2000, ApJ, 534, 34 [Google Scholar]
  110. Wright, A. H., Hildebrandt, H., van den Busch, J. L., & Heymans, C. 2020, A&A, 637, A100 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  111. Xu, L. 2015, arXiv e-prints [arXiv: 1506.03232] [Google Scholar]
  112. Zhao, G.-B., Li, B., & Koyama, K. 2011, Phys. Rev. D, 83, 044007 [NASA ADS] [CrossRef] [Google Scholar]
  113. Zubeldia, Í., & Challinor, A. 2019, MNRAS, 489, 401 [NASA ADS] [CrossRef] [Google Scholar]
  114. Zucca, A., Pogosian, L., Silvestri, A., & Zhao, G. B. 2019, J. Cosmology Astropart. Phys., 2019, 001 [CrossRef] [Google Scholar]
  115. Zumalacárregui, M., Bellini, E., Sawicki, I., Lesgourgues, J., & Ferreira, P. G. 2017, J. Cosmology Astropart. Phys., 2017, 019 [CrossRef] [Google Scholar]

All Tables

Table 1

Parameters of the HMF used in this analysis.

In the text
Table 2

Priors on the cosmological parameters used in this analysis when they are involved in the corresponding model.

In the text
Table 3

Constraints on the HS- f (R) gravity model in different scenarios.

In the text

All Figures

thumbnail Fig. 1

Estimation of the scale and redshift sensitivity of the different cosmological probes, adapted from Preston et al. (2023). The background colors, from blue to green, represent the transition from the linear to the nonlinear regime. We represent cluster abundance with the red rectangle. The cosmology sample of eRASS1 spans the redshift range 0.1 < z < 0.8, as is represented by the solid red rectangle. Future data releases and experiments will increase the upper limit of this interval to the hatched red rectangle. Cluster abundance experiments thus probe a specific regime at LSSs and the nearby Universe, offering a complementary window on structure formation. The constraints from cluster clustering (Seppi et al. 2024) are omitted here. Additionally, the X-ray power spectrum might probe smaller scales.
In the text
thumbnail Fig. 2

Relative difference (in percentages) of the HMF at z = 0.5 when we increase fR0. We note that ∆n = dn/d ln M(fR0) – dn/d ln M, where dn/d ln M(fR0) is the HMF described with Eq. (20) and dn/d ln M is obtained with the GR fitting function from Tinker et al. (2008). Increased values of fR0 create more massive halos.
In the text
thumbnail Fig. 3

Comparison of the posteriors on the scaling relation obtained in the f (R) framework (the mass of the neutrinos is fixed) and the standard ΛCDM analysis presented in Ghirardini et al. (2024) (gray). The different models are fully consistent with one another. We overplot the contours obtained from WL only from Grandis et al. (2024) (green).
In the text
thumbnail Fig. 4

Comparison of the posteriors obtained on the Ω – σ8 plane. We represent the ΛCDM analysis introduced in Ghirardini et al. (2024) in gray, while the posteriors of HS- f (R) case are shown in red.
In the text
thumbnail Fig. 5

Constraints obtained on log fR0 with eRASS1 clusters. No external cosmological datasets are combined with our data, which is shown in this figure. Constraints with a fixed neutrino mass are shown in red, while constraints for which the neutrino mass is free are shown in magenta.
In the text
thumbnail Fig. 6

Comparison of the S 8 parameter obtained using the different models. We show the 1σ uncertainties. The case introducing massless neutrinos is in red, while the one considering massive neutrinos is in pink. Results from the CMB are shown in green. The gray area represents the 68% and 95% confidence regions of the ΛCDM analysis from Ghirardini et al. (2024).
In the text
thumbnail Fig. A.1

Relative difference of the HMF in different MG scenarios. ∆n = dn/dln M(fR0) – dn/d ln M, where dn/d ln M(fR0) is the HMF described with Eq. 20 and dn/d ln M is obtained with the GR fitting function from Tinker et al. (2008). We use the posterior distribution of the correcting parameters s and q and compare their impact to one of the MG parameters. The dark gray area corresponds to the 68% confidence level area, while the light gray one corresponds to the 95% confidence area.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.