Free Access
Issue
A&A
Volume 641, September 2020
Article Number A130
Number of page(s) 37
Section Cosmology (including clusters of galaxies)
DOI https://doi.org/10.1051/0004-6361/202038308
Published online 22 September 2020

© ESO 2020

1. Introduction

Weak gravitational lensing (reviewed by Kilbinger 2015) is the name given to small correlated distortions in the observed shapes of galaxies that are caused by coherent gravitational light deflections between the galaxies and the observer. These shape distortions are of a small amplitude compared to the intrinsic shapes of galaxies and therefore one must have large samples of galaxies in order to tease out signal from the background “shape noise”. The correlation of shape distortions can be used to make maps of a galaxy-redshift-weighted-version of the matter distribution, and these can be used to understand the growth and evolution of matter perturbations in the cosmos. Many forthcoming surveys are being designed to measure this effect very precisely, with the ultimate goal of using structure to learn about the accelerated expansion of the Universe.

Weak lensing on its own is a promising tool to infer the distribution of structure. However, in order to extract accurate cosmological-parameter constraints from data sets, it is necessary to have theoretical models that are significantly more accurate than the survey error bars. On very large scales, density fluctuations are small and linear perturbation theory describes the distribution and evolution or the matter accurately enough for all conceivable future data sets (Lesgourgues 2011). However, the weak-lensing signal is predominantly determined by high amplitude, “non-linear”, fluctuations on small scales (Jain & Seljak 1997) and these have not successfully been described by any ab-initio calculation. The standard scheme for modelling the weak-lensing signal is therefore to run N-body simulations over the range of cosmological scenarios under consideration and then to simply measure the fluctuation spectrum as seen in the simulations. One can then either devise fitting functions for the spectra (e.g. Smith et al. 2003; Takahashi et al. 2012; Mead et al. 2015) or build so-called emulation schemes (e.g. Lawrence et al. 2010, 2017; Agarwal et al. 2012; Knabenhans et al. 2019; Angulo et al. 2020) to interpolate between simulation outputs. The lensing signal can then be calculated by integrating over these power spectra along the line-of-sight with weightings to account for the projection and lensing.

Most cosmological simulations only consider the gravitational interaction1, partly because it is simpler and partly because non-gravitational effects are not understood as well. Early analytic work by White (2004) demonstrated that baryons cooling in a halo may contract and core out the inner halo profile, and it is now known that non-gravitational processes can impact the power spectrum of the matter distribution in a way that significantly biases cosmological parameter constraints if not accounted for (e.g. Semboloni et al. 2011; Mohammed et al. 2014; Schneider et al. 2016; Copeland et al. 2018; Huang et al. 2019). The OWLS suite of hydrodynamical simulations of Schaye et al. (2010) were used by van Daalen et al. (2011) to demonstrate that the matter field can be particularly affected by active-galactic nuclei (AGN) feedback and to show that the main effect that reshapes the matter field is the redistribution of gas that is heated as a result of accretion energy that is generated near the central black hole. This heating is extreme enough that gas can even be expelled from the host haloes. Dismayingly, there is a huge range in the possible amplitude of the effects of feedback (Le Brun et al. 2014; McCarthy et al. 2017; Chisari et al. 2018); if this remains unaddressed, then the ability of the next generation of weak-lensing surveys to learn about the accelerated expansion will be severely compromised. The effect of baryonic feedback on matter clustering and on the structure of haloes has also been investigated using other hydrodynamic simulations: Martizzi et al. (2012, 2013) look at the effect of feedback on halo profiles. Velliscig et al. (2015) do the same with the COSMO-OWLS simulations as well as Mummery et al. (2017), using the more recent “BAryons and HAloes of MAssive Systems” (BAHAMAS2) simulations. All authors agree that the statistics of the matter field on scales relevant for forthcoming weak lensing surveys can be affected at the tens of percent level.

Several techniques have been put forward to mitigate the impact that the unknown feedback strength may have on cosmological-parameter constraints from weak-lensing surveys: Eifler et al. (2015) advocate determining a set of “principal components” from libraries of hydrodynamical simulations and then marginalising over these in a data analysis. Mead et al. (2015) use a modified version of the halo model with parameters that determine the effects of feedback on the halo profiles via a change in concentration and an overall bloating of the halo. The relative merits of these two approaches have been explored by Huang et al. (2019). Several authors (Rudd et al. 2008; Semboloni et al. 2011, 2013; Fedeli et al. 2012, 2014; Fedeli 2014; Debackere et al. 2020) have developed halo models that specifically include a gas component that can be used to model the effect of feedback on the matter–matter spectrum. Recently, van Daalen et al. (2020) have shown that the differences in the effect of feedback on the power spectrum can be mainly attributed to different gas fractions in haloes of ∼1014h−1M, and showed that the feedback amplitude is close to “universal” when expressed in this variable. Another promising technique is the “baryonification” method (Schneider & Teyssier 2015; Schneider et al. 2019; Aricò et al. 2020) where haloes in gravity-only simulations are manually deformed in post-processing in such a way that they have profiles appropriate for their gas content. The deformation can be informed either from hydrodynamic simulations or from observational data. A variant of this method is investigated by Dai et al. (2018) who use the actual equations of motion governing gas physics to perform the halo deformation.

A more optimistic position to take is to regard the impact of feedback on observables as a positive (e.g. Harnois-Déraps et al. 2015; Foreman et al. 2016; MacCrann et al. 2017); we may be able to learn about energetic processes within galaxies by analysing the weak-lensing data. An analysis of CFHTLenS data by MacCrann et al. (2015) included the effect of feedback via a single parameter and a functional form extracted from hydrodynamical simulations. Joudaki et al. (2017a) performed a similar analysis and used the model of Mead et al. (2015) to investigate feedback. Both authors found a modest preference in the data for the presence of feedback. Recent Kilo Degree Survey (KiDS) analyses (Hildebrandt et al. 2017, 2020; Joudaki et al. 2017b) use the model of Mead et al. (2015) and have similar results for the feedback amplitude. Recent results from the Deep Lens Survey (Yoon et al. 2019) using the same model show a stronger preference for strong feedback. However, the matter–matter spectrum may not be the ideal tool if one is interested in energetic processes in galaxies since it is not probing the gas directly.

To measure the gas distribution more directly one can use the thermal Sunyaev-Zel’dovich (tSZ) effect: the inverse Compton scattering of the relatively cold cosmic-microwave background (CMB) photons by relatively energetic thermal electrons that predominantly exist in galaxy clusters. This scattering results in an increase in energy of the CMB photons, with a characteristic frequency-dependence that distorts the standard black-body CMB spectrum in the direction of the galaxy cluster on the sky. The amplitude of the effect depends on the product of electron number density and temperature along the line of sight, a quantity with the same units as pressure and that is therefore known as “electron pressure”. As a scattering process, the tSZ amplitude does not decay with distance and therefore traces the large-scale structure out to high redshift. The characteristic frequency dependence of the induced deviation of the CMB spectrum from blackbody can be exploited to make high fidelity maps of the electron pressure distribution in the Universe (e.g. Planck Collaboration XXII 2016). Due to the fact that tSZ emanates from hot gas in dense clusters, hydrodynamic simulations have been an essential tool for modelling and understanding the effect. Early simulations include those by Suto et al. (1998) and Bond et al. (2005). Trac et al. (2011) developed a technique to paint electron pressure onto a gravity-only simulation and therefore to quickly make mock tSZ maps and power spectrum templates. It was later realised that AGN feedback could have a large impact on the tSZ signal (Dolag et al. 2009; Arnaud et al. 2010; Schaye et al. 2010; McCarthy et al. 2010, 2011; Battaglia et al. 2010), and therefore that these violent feedback episodes needed to be carefully included in simulations, which in turn lead to work to investigate the magnitude of this effect (e.g. Martizzi et al. 2012; Battaglia et al. 2012a,b; Puchwein & Springel 2013; Velliscig et al. 2014; Le Brun et al. 2014, 2015; McCarthy et al. 2014; Sijacki et al. 2015; Dolag et al. 2016).

The lensing–tSZ cross correlation has been measured by Van Waerbeke et al. (2014), Hill & Spergel (2014) and Hojjati et al. (2017). Observational work is often conducted using a stacking approach, where the tSZ map is stacked around the locations of suspected tSZ emitters. Makiya et al. (2018) does the same but with galaxies rather than haloes. Hill et al. (2018) and Tanimura et al. (2020) look at tSZ around SDSS galaxies, see a strong signal and both find evidence for a “two-halo term” – evidence of correlation that extends to scales much beyond the halo boundary. Tanimura et al. (2019a) and de Graaff et al. (2019) looked at residual tSZ emission between close pairs of galaxies and found evidence of the filamentary structure that would be expected to exist given the cosmic web. Tanimura et al. (2019b) mask cluster galaxies and find evidence for a residual “intra-cluster” tSZ signal. All evidence that supports tSZ emanating from low-density regions of the Universe as well as from dense clusters. Hojjati et al. (2015) used the COSMO-OWLS simulations to investigate the lensing–tSZ cross correlation measured by Van Waerbeke et al. (2014), similar work was carried out by Battaglia et al. (2015). An advantage of using the cross correlation between lensing and tSZ is that we can ignore many systematic effects that may plague the individual data sets. However, in some cases systematic effects may be correlated between data sets, which may be the case for thermal dust emission from dusty galaxies that will correlate with the lensing signal and may contaminate tSZ maps during the map-making process. Yan et al. (2019) showed that this may affect the cross correlation between CMB lensing and tSZ, more than between lower-redshift galaxy lensing and tSZ, due to the higher redshift contribution to CMB lensing. It may be possible to address issues of contamination from thermal dust using the methods of Addison et al. (2012, 2013). More recently, the cross correlation has been measured and investigated by Baxter et al. (2019) and Omori et al. (2019) in the context of a cross correlation of CMB and galaxy lensing with the Dark Energy Survey, and by Osato et al. (2018, 2020) with the Hyper Suprime-Cam survey.

To understand the tSZ signal Refregier & Teyssier (2002) construct a halo model and compare this model to a simulation in order to investigate the tSZ auto-spectrum, demonstrating good agreement with early simulations. Komatsu & Seljak (2002) construct a similar model and demonstrate how the tSZ spectrum scales with cosmological parameters, particularly demonstrating an extreme sensitivity to σ8, the amplitude of the linear power spectrum. These early works suggested that the tSZ auto-spectrum could contain a wealth of information about cosmological parameters. Holder et al. (2007) examined the relative impacts of feedback (preheating) and changes in cosmology on various tSZ statistics, including the power spectrum. More recent works attempting to utilise tSZ for cosmology are by Shaw et al. (2010), Bolliet et al. (2018) and Hill et al. (2018). In this paper we do not utilise the tSZ auto-spectrum as we are concerned that it may contain poorly understood internal systematics (Planck Collaboration XXIV 2016; although see Horowitz & Seljak 2017). Hill & Spergel (2014) measured and investigated CMB lensing–tSZ correlation while Ma et al. (2015) investigated the galaxy lensing–tSZ detection from Van Waerbeke et al. (2014). In both cases the cross correlation was modelled using the language of the halo model, with NFW profiles being taken for the dark matter and “universal pressure profiles” (UPP; Arnaud et al. 2010) for the electron pressure.

The purpose of this paper is to develop a halo model that may be used for a joint analysis of lensing and tSZ. We demand that the model be able to reproduce the lensing–lensing auto correlation function and the lensing–tSZ cross correlation. We also demand that the physics underpinning the model be the same so that the two signals depend on each other in a physical way. This is different from previous work in which one profile is used for the matter and a completely separate profile is used for the electron pressure. In this way we hope to be able to learn about cosmological parameters and feedback using both data sets simultaneously. Unless otherwise stated, all figures in this paper are made with cosmological parameters inspired by the WMAP9 (Hinshaw et al. 2013) data analysis: Ωm = 0.2793, Ωb = 0.0463, ΩΛ = 1 − Ωm, σ8 = 0.821, ns = 0.972, h = 0.700. These are the same cosmological parameters that are used in the BAHAMAS hydrodynamical simulations.

This paper is set out as follows: in Sect. 2 we show how the halo model can be used to make predictions for the cross spectrum of any field pair, and how we can improve the model by specifically describing a halo in terms of a cold-dark matter (CDM), gas and stellar component. In Sect. 3 we list our ingredient choices for our halo-model calculations. In Sect. 4 we detail the BAHAMAS simulations, against which we compare our model. In Sect. 5 we tune our halo model by fitting its parameters to simulation data, which is the main results of this paper. Our conclusions are in Sect. 6. We also provide a number of appendices: Appendix A examines a common pitfall regarding computation of the two-halo term. Appendix B discusses how we measure power spectra and shot noise from the BAHAMAS simulations where we have particles with different masses and different properties. In Appendix C we present a series of plots that show how changing our halo-model parameters affects our power spectra. Finally, Appendix D details how we compute projected power spectra of tracers of our underlying three-dimensional fields and we show how angular scales of the projected power spectra receive contributions from the underlying three-dimensional scales.

2. Halo model

In this section we review the basic halo-model computation of the cross spectrum between two three-dimensional cosmological fields (Seljak 2000; Peacock & Smith 2000; Scoccimarro et al. 2001; Cooray & Sheth 2002). We start with a general calculation that can be used to calculate the cross spectrum for any cosmological field pair as long as the halo profiles of the fields are known. Some examples of fields for which this calculation could be applied are “matter overdensity”, “gas overdensity”, “electron pressure”, or “galaxy overdensity”. We then specialise to the details of the calculation appropriate for the lensing–lensing and tSZ–lensing power spectra that are the main focus of this paper, where lensing is generated by the matter overdensity field and tSZ by the electron-pressure field3. The modelling presented here draws heavily on work by Refregier & Teyssier (2002), Shaw et al. (2010), Fedeli (2014), Fedeli et al. (2014), Schneider & Teyssier (2015) and Debackere et al. (2020).

2.1. Cross correlations

In the following we omit time dependence from function arguments to make the notation less cluttered. Consider two 3D cosmological fields, u(r) and v(r). The halo model cross spectrum between the two fields at a fixed redshift can be written as the sum of a two- and a one-halo term, respectively

P 2 H , u v ( k ) = P lin ( k ) i = u , v [ 0 b ( M ) W i ( M , k ) n ( M ) d M ] , $$ \begin{aligned}&P_{2\mathrm{H},uv}(k)=P_{\mathrm{lin} }(k) \prod _{i=u,v}\left[\int _0^\infty b(M)W_i(M,k)n(M)\;\mathrm{d}M\right], \end{aligned} $$(1)

P 1 H , u v ( k ) = 0 W u ( M , k ) W v ( M , k ) n ( M ) d M , $$ \begin{aligned}&P_{1\mathrm{H},uv}(k)=\int _0^\infty W_u(M,k)W_v(M,k)n(M)\;\mathrm{d}M, \end{aligned} $$(2)

where Plin(k) is the linear power spectrum of matter fluctuations, M is the halo mass and b(M) is the linear halo bias with respect to the linear matter density (a dimensionless quantity):

δ h ( M , k 0 ) = b ( M ) δ m ( k 0 ) . $$ \begin{aligned} \delta _{\rm h}(M, k\rightarrow 0) = b(M)\delta _{\rm m}(k\rightarrow 0). \end{aligned} $$(3)

n(M) is the halo mass function: the distribution function for the comoving number density of haloes in a mass range (sometimes denoted dn/dM in the literature; units always per mass interval per volume). Equations (1) and (2) contain (spherical) Fourier transforms of the halo profiles of the fields u(r) and v(r) that are being cross correlated:

W u ( M , k ) = 0 4 π r 2 sin ( k r ) kr θ u ( M , r ) d r , $$ \begin{aligned} W_u(M,k)=\int _0^\infty 4\pi r^2\frac{\sin (kr)}{kr}\theta _u(M,r)\;\mathrm{d}r, \end{aligned} $$(4)

where θu(M, r) is the averaged radial profile for the field u(r) in a host halo of mass M. For example, if one were interested in the electron-pressure field then θu would be the pressure profile of free electrons, if one were interested in the matter power spectrum, then θu would be the halo matter overdensity profile. To use these equations to compute the matter–matter power spectrum then one would set θ u = θ v = ρ m ( M , r ) / ρ ¯ $ \theta_{u}=\theta_{v}=\rho_{\mathrm{m}}(M,r)/\bar\rho $ where ρm is the halo matter density profile and ρ ¯ $ \bar\rho $ is the mean matter density. The division by ρ ¯ $ \bar\rho $ ensures that the end result is the overdensity spectrum and we note that W m ( M , k 0 ) = M / ρ ¯ $ W_{\mathrm{m}}(M,k\to0)=M/\bar\rho $. We work using comoving k, so r in the previous equations is also comoving, as is ρ ¯ $ \bar\rho $, which is therefore constant. The units of the halo profile θu are the units of the field u. The units of the halo-Fourier-transform functions Wu (Eq. (4)) are those of the field u(r) multiplied by volume. We interchange between Δ2(k) and P(k) power spectrum definitions:

Δ uv 2 ( k ) = 4 π ( k 2 π ) 3 P uv ( k ) , $$ \begin{aligned} \Delta ^2_{uv}(k)=4\pi \left(\frac{k}{2\pi }\right)^3P_{uv}(k), \end{aligned} $$(5)

where Δ uv 2 (k) $ \Delta^2_{uv}(k) $ has units of the product of fields u and v, while Puv(k) has this and an additional unit of volume.

The dimensionless mass function, g(ν), normalised such that the integral over all ν gives unity, is related to n(M) via

g ( ν ) d ν = M ρ ¯ n ( M ) d M , $$ \begin{aligned} g(\nu )\;\mathrm{d}\nu =\frac{M}{\bar{\rho }}n(M)\;\mathrm{d}M, \end{aligned} $$(6)

and is written in terms of the mass variable ν = δc/σ(M): δc is the critical linear density threshold for halo collapse: δc ≃ 1.686, which has a weak cosmology dependence. σ(M) is the variance in the linear matter field when filtered on a Lagrangian scale R corresponding to mass M = 4πR3/3:

σ 2 ( R ) = 0 Δ lin 2 ( k ) [ 3 ( k R ) 3 ( sin kR k R cos kR ) ] 2 d ln k , $$ \begin{aligned} \sigma ^2(R)=\int _0^{\infty }\Delta _{\rm lin}^2(k)\left[\frac{3}{(kR)^3}(\sin {kR}-kR\cos {kR})\right]^2\;\mathrm{d}\ln {k}, \end{aligned} $$(7)

where the expression in square brackets is the normalised Fourier transform of a real-space top-hat filter.

Note that the adopted halo mass function and bias must satisfy the following properties for matter power spectra to have the correct large-scale limit4:

0 M n ( M ) d M = ρ ¯ , $$ \begin{aligned}&\int _0^\infty Mn(M)\;\mathrm{d}M=\bar{\rho }, \end{aligned} $$(8)

0 b ( M ) M n ( M ) d M = ρ ¯ . $$ \begin{aligned}&\int _0^\infty b(M) M n(M)\;\mathrm{d}M=\bar{\rho }. \end{aligned} $$(9)

In words: these equations enforce that all matter be in a halo of some mass and that, on average, matter is unbiased with respect to itself.

It is worth examining the approximations that lead to the halo-model Eqs. (1) and (2): it has been assumed that haloes trace the underlying linear matter distribution with a linear halo bias, that halo profiles are perfectly spherical with no substructure and that there is no scatter in profile properties at fixed host halo mass. There is also nothing in the modelling to prevent haloes from overlapping. These approximations will break down, and the errors in the eventual power spectrum that they contribute will vary with the fields that are being considered. Unfortunately remedies for these assumptions, where they exist, make the halo model calculation increasingly complex (e.g. Smith et al. 2007; Giocoli et al. 2010; Valageas & Nishimichi 2011; Smith & Markovic 2011; van den Bosch et al. 2013) and are beyond the scope of the work presented here.

2.2. Halo composition and profiles

In this paper we consider each halo of total mass M to be made of separate CDM, gas and stellar components with (possibly time-dependent) mass fraction fi(M) in component i, such that these sum to unity

f c ( M ) + f g ( M ) + f ( M ) = f m ( M ) = 1 . $$ \begin{aligned} f_{\rm c}(M)+f_{\rm g}(M)+f_{*}(M) = f_{\rm m}(M) = 1. \end{aligned} $$(10)

The gas is further separated into that bound to the halo and that ejected from the halo by some feedback process,

f g ( M ) = f bnd ( M ) + f ejc ( M ) . $$ \begin{aligned} f_{\rm g}(M) = f_{\rm bnd}(M)+f_{\rm ejc}(M). \end{aligned} $$(11)

We use the terms “bound” and “ejected” to differentiate gas that is within the halo virial radius from that outside it. In our formalism it is not necessary that the ejected gas component to be gravitationally bound to the host halo. Indeed, it could be that the ejected gas is far outside the virial radius. All that is required is that the gas component was “associated” with the halo, in that it was part of the initial overdensity from which the halo formed. This means that M should not literally be interpreted as the halo mass, because some mass that was initially associated with the halo may be located outside the virial radius. M is therefore a label, and we take this label to coincide with the measured halo mass in gravity-only simulations, thus justifying our use of a mass function and halo bias calibrated on such simulations. When haloes that are identified in hydrodynamic simulations are compared to non-hydrodynamic simulations that are run with the same initial conditions the halo masses differ object-to-object. However, there is almost always an object-to-object correspondence (e.g. van Daalen et al. 2014). Differing halo masses in hydrodynamical simulations arise predominantly because matter has been moved across the (somewhat-arbitrary) halo boundary.

We also separate the stellar mass into that in central and satellite galaxies,

f ( M ) = f cen ( M ) + f sat ( M ) , $$ \begin{aligned} f_{*}(M) = f_{\rm cen}(M)+f_{\rm sat}(M), \end{aligned} $$(12)

where any black hole mass is included in the stellar component.

Each component of the halo is given a spherical density profile, ρi(M, r), that describes how the component is distributed. These are normalised such that

f i ( M ) M = 0 r v 4 π r 2 ρ i ( M , r ) d r . $$ \begin{aligned} f_i(M) M = \int _0^{r_{\rm v}} 4\pi r^2 \rho _i(M,r)\;\mathrm{d}r. \end{aligned} $$(13)

Note that we can write ρm = ρc + ρg + ρ* and a similar equation for Wm.

In this paper we also compute halo-model spectra under the assumption that all matter is CDM, mainly for comparisons with gravity-only simulations. In this case we ignore the gas and stellar components of the halo and set the CDM fraction to unity.

2.3. Large-scale limit for matter fields

It is instructive to consider the k → 0 limit of both the two- and the one-halo terms in Eqs. (1) and (2) for the specific case of matter fields. Let us illustrate this for the example of the cross spectrum of field u(r) with field v(r), in this case θ u ( M , r ) = ρ u ( M , r ) / ρ ¯ $ \theta_u(M,r)=\rho_u(M,r)/\bar\rho $ and Eqs. (4) and (13) give W u ( M , k 0 ) = f u ( M ) M / ρ ¯ $ W_u(M,k\to0)=f_u(M)M/\bar\rho $, therefore:

P 2 H , u v ( k 0 ) P lin ( k 0 ) = i = u , v [ 1 ρ ¯ 0 b ( M ) f i ( M ) M n ( M ) d M ] , $$ \begin{aligned}&\frac{P_{2\mathrm{H},uv}(k\rightarrow 0)}{P_{\rm lin}(k\rightarrow 0)}= \prod _{i=u,v}\left[\frac{1}{\bar{\rho }}\int _0^\infty b(M)f_i(M)Mn(M)\;\mathrm{d}M\right], \end{aligned} $$(14)

P 1 H , u v ( k 0 ) = 1 ρ ¯ 2 0 f u ( M ) f v ( M ) M 2 n ( M ) d M . $$ \begin{aligned}&P_{1\mathrm{H},uv}(k\rightarrow 0)=\frac{1}{\bar{\rho }^2}\int _0^\infty f_u(M) f_v(M) M^2 n(M)\;\mathrm{d}M. \end{aligned} $$(15)

If we adopt the notation

x = 1 ρ ¯ 0 x ( M ) M n ( M ) d M , $$ \begin{aligned} \langle x\rangle =\frac{1}{\bar{\rho }}\int _0^\infty x(M) M n(M)\;\mathrm{d}M, \end{aligned} $$(16)

for an average of property x over halo mass then Eq. (14) becomes

P 2 H , u v ( k 0 ) P lin ( k 0 ) = b f u b f v . $$ \begin{aligned} \frac{P_{2\mathrm{H},uv}(k\rightarrow 0)}{P_{\rm lin}(k\rightarrow 0)}=\langle b f_u\rangle \langle b f_v\rangle . \end{aligned} $$(17)

We see that the amplitude of the two-halo term is governed by the product of the mean-halo-bias-weighted abundances for fields u(r) and v(r). If we consider the auto-spectrum of matter overdensity, then u = v = m, fm(M) = 1 and Eq. (17) tells us that the two-halo term is the linear power spectrum because ⟨b⟩ = 1 (Eq. 9). If fu(M) does not depend on halo mass then the average would reduce to Ωum, where Ωu is the cosmological density parameter for species u. More generally fu(M) will depend on halo mass, and the two-halo term on large scales is then the linear power spectrum weighted by the field abundance multiplied by the halo bias of the haloes in which the field is to be found.

If instead we focus on the one-halo term, Eq. (15) becomes

P 1 H , u v ( k 0 ) = M f u f v ρ ¯ , $$ \begin{aligned} P_{1\mathrm{H},uv}(k\rightarrow 0)=\frac{\langle M f_u f_v\rangle }{\bar{\rho }}, \end{aligned} $$(18)

and we see that the amplitude of the one-halo term is governed by the halo-mass-weighted product of the abundances. For the matter–matter case, u = v = m then fm = 1 and the amplitude of the one-halo term is governed by ⟨M⟩, which is the mean halo mass that a unit of matter is to be found in (not the mean halo mass). For our fiducial cosmology at z = 0 this mass is ∼1013.5h−1M, which gives some indication of the typical halo mass responsible for power in the one-halo term of the matter–matter spectrum.

3. Specific ingredients

The method presented in Sect. 2 has been general, but here we list the specific ingredients used in this paper. We are interested in both the lensing auto correlation function and the cross correlation between tSZ and lensing. This led us to compute the 3D matter–matter and matter–electron pressure power spectra. We are also interested in testing how our model for the different components of a halo compares to those components measured in hydrodynamical simulations, so we shall also compute auto- and cross-spectra of CDM, gas and stars. In this section we take default values for hydrodynamical parameters, but we later free and fit these in Sect. 5.

3.1. Halo demographics

In our calculations we adopt the mass function (g(ν) from Eq. (6)) of Sheth & Tormen (1999)

g ( ν ) d ν = A [ 1 + 1 ( q ν 2 ) p ] e q ν 2 / 2 d ν , $$ \begin{aligned} g(\nu )\;\mathrm{d}\nu =A\left[1+\frac{1}{(q\nu ^2)^p}\right]\mathrm{e}^{-q\nu ^2/2}\;\mathrm{d}\nu , \end{aligned} $$(19)

with p = 0.3, q = 0.707 and A ≃ 0.216. Note that this function has no explicit redshift dependence. We use the halo bias derived from Eq. (19) using the peak-background split formalism (Mo & White 1996; Sheth et al. 2001)

b ( ν ) = 1 1 δ c [ 1 + ν d d ν ln g ( ν ) ] , $$ \begin{aligned} b(\nu )=1-\frac{1}{\delta _{\rm c}}\left[1+\nu \frac{\mathrm{d}}{\mathrm{d}\nu }\ln g(\nu )\right], \end{aligned} $$(20)

which then automatically fulfils the mean-bias condition in Eq. (9). Explicitly for the mass function in Eq. (19) the peak-background split gives

b ( ν ) = 1 + 1 δ c [ q ν 2 1 + 2 p 1 + ( q ν 2 ) p ] · $$ \begin{aligned} b(\nu )=1+\frac{1}{\delta _{\rm c}}\left[q\nu ^2-1+\frac{2p}{1+(q\nu ^2)^p}\right]\cdot \end{aligned} $$(21)

The Sheth & Tormen (1999) relation was fitted to haloes that were identified in N-body simulations with a cosmology-dependent virial overdensity criterion taken from the spherical-collapse model: spherical-overdense regions that are a factor of Δv times denser that the background matter density. We take Δv(z) from the ΛCDM fitting function of Bryan & Norman (1998)

Δ v ( z ) = 1 Ω m ( z ) { 18 π 2 82 [ 1 Ω m ( z ) ] 39 [ 1 Ω m ( z ) ] 2 } . $$ \begin{aligned} \Delta _{\rm v}(z)=\frac{1}{\Omega _{\rm m}(z)}\left\{ 18\pi ^2-82[1-\Omega _{\rm m}(z)]-39[1-\Omega _{\rm m}(z)]^2\right\} . \end{aligned} $$(22)

Δv is known as the virial-collapse density5. The cosmology dependence of the linear-collapse density, δc, calculated according to the spherical-collapse model, was also taken into account in the mass function of Sheth & Tormen (1999). Even though this is a small numerical change it has a larger impact upon the mass function than one might expect due to the fact that it is exponentiated (Courtin et al. 2011; Mead 2017). We use the ΛCDM fitting formula from Nakamura & Suto (1997):

δ c ( z ) = 3 20 ( 12 π ) 2 / 3 { 1 + 0.0123 log 10 Ω m ( z ) } . $$ \begin{aligned} \delta _{\rm c}(z)=\frac{3}{20}(12\pi )^{2/3}\left\{ 1+0.0123\log _{10}\Omega _{\rm m}(z)\right\} . \end{aligned} $$(23)

Note that it is important to use halo structural parameters that are consistent with the way haloes are defined in the mass function. We take haloes defined with the virial definition in Eq. (22) and any conversions between halo mass M and ν use the cosmology-dependent δc from Eq. (23).

3.2. Halo composition

We consider haloes to be made from CDM, gas, and stars. We set the halo CDM fraction to be the constant universal fraction

f c ( M ) = Ω c Ω m , $$ \begin{aligned} f_{\rm c}(M)=\frac{\Omega _{\rm c}}{\Omega _{\rm m}}, \end{aligned} $$(24)

which we consider appropriate because, while feedback may indirectly redistribute CDM within a halo, it will not be able to eject CDM.

Gas is split into bound and ejected gas, with the fraction of bound gas in haloes taken from Schneider & Teyssier (2015):

f bnd ( M ) = Ω b Ω m ( M / M 0 ) β 1 + ( M / M 0 ) β · $$ \begin{aligned} f_{\rm bnd}(M)=\frac{\Omega _{\rm b}}{\Omega _{\rm m}}\frac{(M/M_0)^\beta }{1+(M/M_0)^\beta }\cdot \end{aligned} $$(25)

This function ensures that the bound gas fraction transitions from the universal baryon fraction in high-mass haloes (M ≫ M0) to zero in lower mass systems, with β governing the rate of transition; haloes of mass M0 having lost half of their initial baryon content. By default we take M0 = 1014h−1M and β = 0.6 as in Schneider & Teyssier (2015) as this in in rough accordance with observational results (e.g. Sun et al. 2009; Vikhlinin et al. 2009; Gonzalez et al. 2013). The ejected gas fraction is then fixed such that it accounts for all remaining baryons that were originally associated with the halo assuming that the perturbation was adiabatic

f ejc ( M ) = Ω b Ω m f bnd ( M ) f ( M ) . $$ \begin{aligned} f_{\rm ejc}(M)=\frac{\Omega _{\rm b}}{\Omega _{\rm m}}-f_{\rm bnd}(M)-f_{*}(M). \end{aligned} $$(26)

With these definitions, the requirement that the total halo mass is accounted for (Eq. (10)) is automatically satisfied6.

We follow Fedeli (2014) and set the stellar fraction to be

f ( M ) = A exp [ log 10 2 ( M / M ) 2 σ 2 ] , $$ \begin{aligned} f_{*}(M)=A_*\exp \left[-\frac{\log ^2_{10}(M/M_*)}{2\sigma _*^2}\right], \end{aligned} $$(27)

a form that expresses the fact that star formation efficiency peaks in haloes of mass M* with halo stellar mass fraction A*, while being suppressed for higher and lower halo masses with a logarithmic width of σ*. By default we take A* = 0.03, σ* = 1.2 and M* = 1012.5h−1M, values motivated by Moster et al. (2013) and Kravtsov et al. (2018). We also impose a limit for the high-mass end of f*: f*(M >  M*)≥A*/3, suggested by observational data on the high-halo-mass saturation of the relation between stellar mass and halo mass (e.g. Leauthaud et al. 2012). The stellar halo mass fraction is further split into central and satellite stars, with the stellar content of low-mass haloes assumed to be dominated by a single central galaxy while the stellar content of high-mass haloes is dominated by satellites. For M <  M* we take

f cen ( M ) = f ( M ) , f sat ( M ) = 0 , $$ \begin{aligned} f_{\rm cen}(M) = f_*(M), \quad f_{\rm sat}(M) = 0, \end{aligned} $$(28)

while for M >  M* haloes we have

f cen ( M ) = f ( M ) ( M M ) η , $$ \begin{aligned} f_{\rm cen}(M) = f_*(M)\left(\frac{M}{M_*}\right)^\eta , \end{aligned} $$(29)

f sat ( M ) = f ( M ) [ 1 ( M M ) η ] , $$ \begin{aligned} f_{\rm sat}(M) = f_*(M)\left[1-\left(\frac{M}{M_*}\right)^\eta \right], \end{aligned} $$(30)

with a default of η = −0.3, taken from Moster et al. (2013), such that the stellar mass content of high-mass haloes is dominated by satellites while that of low-mass haloes is dominated by centrals.

In Fig. 1 we show the halo mass fractions of the different components as a function of halo mass: The CDM fraction is constant across halo mass. The bound gas fraction is close to the universal baryon fraction for high mass haloes but depletes below half the universal content for M0 <  1014h−1M. Low-mass haloes have more baryon content in stars than in gas. The stellar fraction peaks at M* = 1012.5h−1M with a stellar mass fraction A* = 0.03. Stellar mass in satellite galaxies dominates higher-mass haloes, while stellar mass in central galaxies dominates the lower-mass haloes.

thumbnail Fig. 1.

Halo mass fractions fu(M) as a function of halo mass for the different components in our model at z = 0: The halo CDM fraction is constant at Ωcm for all halo masses. Star formation efficiency peaks around M* = 1012.5h−1M with a mass fraction A* = 0.03, it saturates at high mass at a value of f* = A*/3. Stars are split into central and satellite, with the halo stellar content being dominated by centrals at low masses and by satellites at high mass. The bound gas fraction is near to the universal baryon fraction for high halo masses, but decreases below half the universal fraction for halo masses lower than M0 = 1014h−1M, the remaining baryons that are neither bound gas or stars are considered to be in unbound gas that is not in the halo. The total mass fraction drops from unity at lower masses due to the ejected gas.

3.3. Halo profiles

In gravity-only simulations it has been known since Navarro et al. (1997; NFW) that gravitational collapse causes dark matter to arrange itself into approximately spherical structures following the “NFW” profile

ρ c ( M , r ) 1 r / r s ( 1 + r / r s ) 2 · $$ \begin{aligned} \rho _{\rm c}(M,r)\propto \frac{1}{r/r_{\rm s}(1+r/r_{\rm s})^2}\cdot \end{aligned} $$(31)

For our halo-model calculation, the profile is truncated at the virial radius, defined such that this encloses an average density of Δv times the mean background density:

M = 4 π 3 r v 3 Δ v ρ ¯ , $$ \begin{aligned} M=\frac{4\pi }{3} r_{\rm v}^3\Delta _{\rm v}\bar{\rho }, \end{aligned} $$(32)

where Δv, the virial-collapse density, is taken from Eq. (22). As we always work in comoving units we take ρ ¯ $ \bar\rho $ to be the (constant) comoving matter density and rv is therefore also comoving. The remaining halo scale radius parameter, rs is then usually specified via the concentration c = rv/rs and we use the concentration-mass relation of Duffy et al. (2008) appropriate for the virial definition of halo radius

c ( M ) = 7.85 ( M 2 × 10 12 h 1 M ) 0.081 ( 1 + z ) 0.71 . $$ \begin{aligned} c(M)=7.85\left(\frac{M}{{2}\,\times \,10^{12}\,h^{-1}\,{M_\odot }}\right)^{-0.081}(1+z)^{-0.71}. \end{aligned} $$(33)

It has been shown (e.g. Rudd et al. 2008; Velliscig et al. 2014; Mummery et al. 2017) that feedback physics changes the dark-matter halo profiles seen in hydrodynamic simulations when compared to profiles in gravity-only simulations. However, the profiles are still well fit by the NFW form, but have slightly altered concentration parameters. In a similar way to Mead et al. (2015), we introduce a change in concentration parameter caused by the ejection of gas in a way that modifies the standard concentration multiplicatively,

c ( M ) c ( M ) [ 1 + ϵ 1 + ( ϵ 2 ϵ 1 ) f bnd ( M ) Ω b / Ω m ] · $$ \begin{aligned} c(M)\rightarrow c(M)\left[1+\epsilon _1 +(\epsilon _2-\epsilon _1)\frac{f_{\rm bnd}(M)}{\Omega _{\rm b}/\Omega _{\rm m}}\right]\cdot \end{aligned} $$(34)

The concentrations of (low-mass) haloes that have lost all of their gas are multiplied by (1 + ϵ1), while those of (high-mass) haloes that have not lost any gas are multiplied by (1 + ϵ2). We take ϵ1 = ϵ2 = 0 as default, which corresponds to no modification. It may be possible to implement an analytical recipe for how gas changes the dark matter profile, such as that presented in White (2004) or Schneider & Teyssier (2015) but we find that this is not necessary for our study.

Gas that is gravitationally bound to a halo is described by the Komatsu & Seljak (2001; KS) density profile, which we take as

ρ bnd ( M , r ) [ ln ( 1 + r / r s ) r / r s ] 1 / ( Γ 1 ) , $$ \begin{aligned} \rho _{\rm bnd}(M,r)\propto \left[\frac{\ln (1+r/r_{\rm s})}{r/r_{\rm s}}\right]^{1/(\Gamma -1)}, \end{aligned} $$(35)

with rs being identical to that in Eq. (31); Γ is the polytropic index for the gas, which we take to be Γ = 1.17 by default. Decreasing Γ makes the gas profile more centrally concentrated. Equation (35) is taken from Martizzi et al. (2013), and is a slightly simplified version of the original KS profile (see Komatsu & Seljak 2002), which we use because it is more convenient to numerically Fourier transform. Rabold & Teyssier (2017) have shown that the KS model provides an extremely good description of haloes that have not undergone violent feedback and Yan et al. (2020) have shown that it also provides a good model for BAHAMAS haloes if Γ is allowed to vary.

We treat the ejected gas in a manner similar to that in Fedeli (2014) and Debackere et al. (2020), which in turn is similar to the way that warm dark matter is treated in Smith & Markovic (2011): we make the approximation that the ejected gas does not contribute to the one-halo term, but only to the two-halo term7. This gives an additive contribution to the two-halo term which is exactly the linear power spectrum multiplied by ⟨bfejc⟩ per gas field, where

b f ejc = 1 ρ ¯ 0 b ( M ) f ejc ( M ) M n ( M ) d M . $$ \begin{aligned} \langle bf_{\rm ejc}\rangle =\frac{1}{\bar{\rho }}\int _0^\infty b(M) f_{\rm ejc}(M) M n(M)\;\mathrm{d}M. \end{aligned} $$(36)

This is equivalent to taking the ejected gas to be a biased tracer of the linear matter density field. It remains to be determined if this is a reasonable assumption, but it does have the virtue of ensuring that we have accounted for all of the mass and thus preserving the physical relations in Eqs. (8) and (9).

We take the stellar density profile for the central galaxies to be a delta function located at the halo centre

ρ cen ( M , r ) δ D ( r ) , $$ \begin{aligned} \rho _{\rm cen}(M,r)\propto \delta _{\rm D}(r), \end{aligned} $$(37)

while we assume that satellite galaxies ρsat(M, r) follow the NFW profile, given in Eq. (31).

To determine the constant of proportionality for the density profile of each component (Eqs. (31), (35) and (37)) we use Eq. (13): assuring that the component makes up the correct mass fraction of a given halo.

To compute any spectra involving the electron pressure we also need to know the temperature of ionised electrons in haloes. We assume all gas to be ionised, and for the bound gas we can use the KS profile to determine the gas temperature:

T g ( M , r ) = T v ( M ) ln ( 1 + r / r s ) r / r s , $$ \begin{aligned} T_{\rm g}(M,r)= T_{\rm v}(M)\frac{\ln (1+r/r_{\rm s})}{r/r_{\rm s}}, \end{aligned} $$(38)

which amounts to assuming hydrostatic equilibrium. We compute the central temperature Tv(M) as the halo virial temperature for the gas:

3 2 k B T v ( M ) = α G M m p μ p a r v , $$ \begin{aligned} \frac{3}{2}k_{\rm B}T_{\rm v}(M)=\alpha \frac{GMm_{\rm p}\mu _{\rm p}}{a r_{\rm v}}, \end{aligned} $$(39)

where mp is the proton mass, μp is the mean gas particle mass divided by the proton mass, and the factor of a converts the comoving rv to a physical radius. The total halo electron pressure Pe is then given by the ideal gas law

P e ( M , r ) = ρ bnd ( M , r ) m p μ e k B T g ( M , r ) , $$ \begin{aligned} P_{\rm e}(M,r)= \frac{\rho _{\rm bnd}(M,r)}{m_{\rm p}\mu _{\rm e}} k_{\rm B}T_{\rm g}(M,r), \end{aligned} $$(40)

where ρbnd is the halo bound gas density, mp is the proton mass and μe is mean gas particle mass per electron divided by the proton mass8. The first term in this equation is exactly the electron number density. α in Eq. (39) is a parameter that encapsulates deviations from a simple virial relation, for example, unvirialised or unthermalised gas or turbulence, we take the standard α = 1 by default. In our calculations we have neglected relativistic corrections to the effective pressure that contributes to the tSZ spectrum (Itoh et al. 1998; Nozawa et al. 2000), these corrections may be important to consider in the future (Remazeilles et al. 2019).

The ejected gas is treated as if it is distributed with the background linear perturbations, so that it contributes only to the two-halo term. It is given a constant temperature, Tw, which is supposed to be the temperature of the warm-hot intergalactic medium (WHIM; Cen & Ostriker 1999), but we caution the reader against taking this name too literally. We take Tw = 106.5 K as the default temperature of the WHIM (Van Waerbeke et al. 2014).

3.4. Large-scale limit for electron pressure

It is informative to calculate the amplitude scaling of the two- and one-halo terms for the case of the electron pressure, to compare to matter. If we take the k → 0 limit for the halo Fourier transforms (Eq. (4)) of matter we get Wm(M, k → 0) ∝M, but for the electron pressure that originates from bound halo gas this changes to

W p ( M , k 0 ) M 5 / 3 . $$ \begin{aligned} W_{\rm p}(M,k\rightarrow 0) \propto M^{5/3}. \end{aligned} $$(41)

The extra M2/3 for electron pressure comes from the fact that pressure is the product of gas density, which is ∝M, and gas temperature (Eq. (39)), which is ∝M2/3. Using the notation from Sect. 2.3 we can write the large-scale limit of the two-halo term as

P 2 H , mp ( k 0 ) P lin ( k 0 ) = A b f bnd M 2 / 3 + B T w b f ejc , $$ \begin{aligned} \frac{P_{\rm 2H,mp}(k\rightarrow 0)}{P_{\rm lin}(k\rightarrow 0)}= A\langle b f_{\rm bnd} M^{2/3}\rangle + BT_{\rm w}\langle b f_{\rm ejc}\rangle , \end{aligned} $$(42)

where the second term originates from the pressure of non-halo WHIM gas (A and B are constants). For the one-halo term we have

P 1 H , mp ( k 0 ) f bnd M 5 / 3 . $$ \begin{aligned} P_{\rm 1H,mp}(k\rightarrow 0)\propto \langle f_{\rm bnd} M^{5/3}\rangle . \end{aligned} $$(43)

The extra M2/3 in these equations compared to those for the matter–matter power spectrum cause the matter–electron pressure cross spectrum to be dominated by relatively more massive haloes. A corollary of this is that the transition between the two- and one-halo term occurs at larger scales for the matter–electron pressure spectrum compared to the matter–matter spectrum.

3.5. Example power spectra

In Fig. 2 we show the halo-model predictions the spectra of matter–matter, CDM–CDM, gas–gas, stars–stars and matter–electron pressure. We see that at large scales the shape of all curves are identical, a consequence of all components following the linear perturbation distribution (Eq. (14)). For the CDM–CDM spectrum, the ratio at large scales compared to the matter spectrum is exactly (Ωcm)2. For gas–gas and stars–stars it is roughly (Ωgm)2 and (Ω*m)2, but this is not exact because the bias of the haloes in which each field lives also contributes to the two-halo term amplitude (see Eq. (17)). At smaller scales we see that the shapes of the one-halo terms are all quite different, which arises from the different halo profiles of each component. Note that for k ≃ 10 h Mpc−1 the stellar density distribution eventually has more power than the gas, despite the lower abundance: a consequence of the gas being smoothly distributed on small scales while stars are tightly clustered in halo centres. The CDM–CDM spectrum is very similar to the total matter spectrum, which makes sense since CDM dominates the matter budget; however, there are shape differences at smaller scales where the different scale-dependence of contributions of gas–gas and stars–stars are important. The shape of the matter–electron pressure spectrum looks superficially like that of the gas–gas spectrum, but with the two- to one-halo transition taking place at a different scale, due to the different scalings of the two- and one-halo terms. Note that this spectrum cannot be directly compared with the others since it has units of pressure whereas the others are dimensionless.

thumbnail Fig. 2.

Halo-model predictions at z = 0 for the matter–matter (top), CDM–CDM, gas–gas, stars–stars and the matter–electron pressure (bottom) power spectra. We also show the two-halo (long-dashed) and the one-halo (short-dashed) terms for each spectrum apart from matter–matter. At large scales, the shapes of the two-halo terms are all identical, and are that of the linear power spectrum, at small scales the one-halo terms dominates and have significantly different shapes for each spectrum. Note that the transition scale between the two- and one-halo term is at much larger scales for the matter–electron pressure cross spectrum, which is a consequence of this spectrum being dominated by contributions from more massive haloes compared to those that contribute to the matter spectra. The down-turn at small scales in the gas–gas and matter–electron pressure spectra is partly due to the underrepresentation of low-mass haloes and partly due to their relatively smooth halo profiles.

In Fig. 3 we show how the z = 0 matter–matter spectrum and the matter–electron pressure spectrum are built up as a function of halo mass. We show cumulative power spectra computed as we vary the upper limit of halo mass in the integrals in Eqs. (1) and (2) from 1010 to 1016h−1M. In both cases the integration has converged with an upper limit of 1016h−1M. The matter–matter spectrum builds up more gradually as a function of halo mass compared to the matter–electron pressure spectrum, the latter receiving the majority of the power from haloes between 1014 and 1015h−1M and very little power from haloes below 1013h−1M. This, combined with the relatively smooth pressure profiles, contributes to the existence of a peak in the model cross spectrum at k ∼ 7 h Mpc−1, although we caution the reader that this peak may not be physical since we do not model halo substructure. By extension, this implies that the matter–electron pressure cross spectrum is sensitive to higher halo masses than the standard matter–matter spectrum (e.g. Makiya et al. 2018). The reason for this is two fold: first, Eq. (25) tells us that low-mass haloes are deficient in gas and therefore contribute less free-electron number density to the signal. Second, Eqs. (42) and (43) tell us that the overall amplitude of both the two- and the one-halo terms is more sensitive to high-mass haloes. This also ensures that the non-linear scaling of these spectra with σ8 is more extreme, as can be seen in Fig. 4, since increasing the power spectrum amplitude a small amount has a comparatively large impact on the high-mass tail of the halo-mass function. We find that roughly P mm σ 8 3.5 $ P_{\rm mm}\propto \sigma_8^{3.5} $, P mp σ 8 5.8 $ P_{\rm mp}\propto \sigma_8^{5.8} $ and P pp σ 8 8 $ P_{\rm pp}\propto \sigma_8^8 $ with the exact exponent depending on wavenumber. This fact has been noticed in the past (e.g. Komatsu & Seljak 2002; Hill & Pajer 2013; Hill & Spergel 2014) and motivates the use of the tSZ auto-spectrum (Komatsu & Seljak 2002; McCarthy et al. 2014), or the cross correlation of tSZ with lensing (Hill & Spergel 2014; Ma et al. 2015; Hojjati et al. 2017), as a sensitive probe of the power-spectrum amplitude. This statement is equivalent to noting that the high-mass tail of the halo mass function is a sensitive probe of the power spectrum amplitude.

thumbnail Fig. 3.

Cumulative halo model power spectra as a function of the upper-limit of halo mass at z = 0. Left-hand panel: matter–matter power spectrum as the upper limit is raised from 1010 to 1016h−1M (i.e. the lightest coloured curve shows the contribution from only haloes below 1010h−1M), right-hand panel: same but for the matter–electron pressure spectrum. In both cases the converged spectrum from all halo masses is shown in thick black and this has been reached with an upper limit of 1016h−1M. We see that the matter–matter spectrum builds up fairly gradually with halo mass and that all scales receive contributions from a wide range of halo masses. In contrast, the matter–electron pressure spectrum receives very little contribution from haloes less massive than 1013h−1M, with the vast majority of the power coming from haloes between 1014 and 1015h−1M.

thumbnail Fig. 4.

Effect of σ8 variations on the matter–matter power spectrum (left) and on the matter–electron pressure power spectrum (right) at z = 0. We show halo model power spectra for σ8 = 0.7 (blue) to σ8 = 0.9 (red) divided by a central σ8 = 0.8 model. For linear theory the plot would simply be horizontal lines. We see that boosting the amplitude of fluctuations increases the amplitude of the matter–electron power spectrum far more than it does for the matter-matter. At a wavenumber of ≃1 h Mpc−1 the matter-matter power spectrum scales like σ 8 3.5 $ \sigma_8^{3.5} $ while matter–electron pressure scales like σ 8 5.8 $ \sigma_8^{5.8} $. This is because the electron pressure field is dominated by high-mass haloes from the tail of the mass function that are very sensitive to the linear power spectrum amplitude.

4. Hydrodynamic simulations

In this paper we compare auto- and cross-spectra from our halo model to those measured from the BAHAMAS simulations (McCarthy et al. 2017). These are a set of smooth-particle hydrodynamical simulations with separate dark matter and gas particles. The dark matter interacts only gravitationally but the gas also experiences a hydrodynamic pressure force. In addition, processes like star formation, gas heating and cooling, supernovae explosions and AGN formation and evolution are also modelled using “sub-grid” recipes. A full discussion of the sub-grid implementation can be found in Schaye et al. (2010), Le Brun et al. (2014), and McCarthy et al. (2017). Hydrodynamic gas particles have temperature and density as properties, as well as the standard position and velocity, and these can be used to calculate an electron pressure per particle (Appendix B). Star (macro) particles are also created from the gas as the simulation evolves, and these keep track of star formation and death. Black holes are also considered, which can swallow gas particles and grow, although they contribute little to the overall mass.

The BAHAMAS simulation boxes are 400 h−1 Mpc cubes with periodic boundary conditions. Each simulation is of 10243 dark matter and 10243 gas particles with cosmological parameters inspired by the WMAP9 (Hinshaw et al. 2013) data analysis. Initial conditions are created at z = 127 using second-order Lagrangian perturbation theory via the code SGENIC9 (Bird et al. 2020) with appropriate, different transfer functions used for dark matter and baryons. We utilise three different hydrodynamic simulations that differ only in the “strength” of their AGN feedback. This strength is determined by a single parameter, the AGN sub-grid heating temperature, which is the temperature increase given to gas particles that are targeted for feedback. The default AGN model, for which we have 3 realisations, was calibrated to reproduce the amplitude of the observed gas–halo mass relation of groups and clusters, as inferred from X-ray observations (McCarthy et al. 2018). The AGN-LO and AGN-HI models lowered and raised the sub-grid heating temperature so as to approximately bracket the scatter in the observed relation. For comparison purposes there is also a non-hydrodynamic “gravity-only” simulation (DMONLY). This simulation still uses two sets of particles for gas and dark matter, and these different particles start with different initial conditions. The only difference compared to the full hydrodynamic simulations is in the subsequent evolution, where only gravity is considered in DMONLY.

In Fig. 5 we show different fields around a massive halo for the AGN simulation at z = 0. We see that the CDM forms a skeleton around which the other fields are clustered, with the gas density being less tightly clustered than CDM and stellar density being more tightly clustered. Both gas and stars are missing from the lower-mass haloes. The electron pressure follows the gas density, but is restricted to emanate from only the highest gas-density peaks. This is because high-mass haloes also have a higher temperature compared to low-mass haloes, which boosts their electron-pressure signal relative to the density. Low-mass haloes are therefore severely under represented in the electron-pressure distribution.

thumbnail Fig. 5.

Different fields measured from the BAHAMAS AGN simulation at z = 0 around a M ∼ 1015h−1M halo averaged through a square slab of side 20 h−1 Mpc and thickness 5 h−1 Mpc. We show CDM (top left), gas (top right) and star (bottom left) density contrast as well as electron pressure (bottom right). Colour intensity increases with field value logarithmically and we show a dynamic range of 103 for both density contrast and electron pressure. We see that the gas is less clumped than the CDM and that low-mass sub haloes show a large gas depletion compared to the host. Stars are more tightly clustered in the halo centre than the CDM, but are absent in many of the substructures. The total matter density-contrast field would be a sum of the CDM, gas and star fields. The electron pressure broadly follows the gas but can only be seen emanating from the highest gas-density peaks over the dynamic range shown.

In this paper, we are interested in two-point statistics that pertain to weak gravitational lensing and to the tSZ effect. We therefore measure the power spectrum of density fluctuations and the electron-pressure field from the simulations. We work with the full 3D data, rather than with 2D projections, as these are more directly tied to the modelling. Details of how we measure these power spectra from the particles and how we consider the effects of shot noise can be found in Appendix B. We compute the auto- and cross-spectra of all combinations of the fields: total matter overdensity δm, CDM overdensity, δc, gas overdensity, δg, stellar overdensity, δ* and electron pressure Pe. In our work, all overdensities are defined relative to the total matter density (i.e. 1 + δ u = ρ u / ρ ¯ $ 1+\delta_u=\rho_u/\bar\rho $), which ensures

1 + δ m = ( 1 + δ c ) + ( 1 + δ g ) + ( 1 + δ ) . $$ \begin{aligned} 1+\delta _{\rm m}=(1+\delta _{\rm c})+(1+\delta _{\rm g})+(1+\delta _{\rm *}). \end{aligned} $$(44)

The measured power spectra from the BAHAMAS AGN simulation are shown in Fig. 6 together with the halo-model prediction from the model presented in Sect. 3 with the default parameter values. We note that the halo model provides a reasonable model of the data for all spectra shown, but that it is not perfect (the log-log scale can hide some serious defects). At large scales we note that the amplitudes are generally in good agreement, which tells us that we have the overall abundances and halo occupation of each species approximately correct and that we have the mean background electron pressure reasonably well modelled. At smaller scales we note differences that must be due to an incorrect mass function, choice of halo profiles or else due to physical effects that are missing from our simple halo model. For the matter spectra we also note that the standard problem of a power deficit in the transition region between the two- and one-halo terms (e.g. Tinker et al. 2005; Valageas & Nishimichi 2011; Mead et al. 2015) is present for all spectra shown. This defect seems to be less of a problem in the spectra involving the electron pressure. With reference to Fig. 2 we can conjecture that this may be because these spectra are more dominated by the one-halo term and because the transition between the two- and one-halo terms takes place at comparatively larger scales. The existence of a peak in both the matter–electron pressure cross spectrum and electron pressure auto spectrum is hinted at, but by no means confirmed, by the simulation measurements.

thumbnail Fig. 6.

Comparison of the default halo-model prediction with power spectra from the BAHAMAS AGN simulation at z = 0. Points with errors show the measured power spectra with an error-on-the-mean calculated from the finite number of modes that contribute to each k. The upper two panels show spectra for “all matter” (purple), CDM (green), gas (light blue), stars (orange) and electron pressure (dark blue). The electron pressure field has units of 100 eV cm−3. The effect of power aliasing can be seen as an upturn in power at the highest k (∼7 h Mpc−1) shown for each spectrum (see Appendix B). Lines show the halo-model predictions using the default model discussed in Sect. 3: solid lines and points show auto-spectra while broken lines and open points show cross-spectra with the total matter field. We see that the trends in the halo model agree reasonably well with the simulations but that there are disagreements in the details. In the lower two panels we show the “response function” calculated with respect to the matter–matter power in the gravity-only model: P uv hydro ( k ) / P mm gravity ( k ) $ P^{\mathrm{hydro}}_{uv}(k) / P^{\mathrm{gravity}}_{\mathrm{mm}}(k) $. The simulations spectra have been divided by those from the DMONLY simulation and the halo model has been divided by a standard halo-model prediction that assumes all matter to be in NFW haloes. The horizontal-dashed-black lines in the lower left panel show the expected large scale Ωcm, (Ωcm)2, Ωbm, and (Ωbm)2 asymptotes for the CDM and gas cross–matter and auto-spectra. We see that the response functions are generally smoother than the raw power spectra: a consequence of cosmic variance cancellations at large scales and cancellation of aliasing effects at small scales.

The lower panels of Fig. 6 show the power spectrum response functions (e.g. Mead 2017; Cataneo et al. 2019), which we define as the ratio of any power spectrum of any field combination to the matter–matter power spectrum measured in a “gravity-only” model. In the case of the BAHAMAS simulations we calculate this response with respect to the DMONLY simulation. The response has the virtue that it cancels the Gaussian cosmic variance at large scales, which leads to a smooth large-scale response. The measured response functions for the BAHAMAS spectra that involve the electron pressure are less smooth at large scales, which is probably because of the comparative dominance of the one-halo term at large scales and therefore that the large-scale noise is likely to be dominated by Poisson fluctuations in massive halo numbers, rather than the Gaussian variance present in the initial conditions. For the halo model, we calculate this response ratio with respect to a halo model calculation assuming all mass to be in NFW haloes (fc = 1). That the response functions are constant at large scales is indicative that the power on these scales is simply the linear matter spectrum multiplied by some weighting (Eqs. (17) and (42)). For the halo model, the response has the virtue of cancelling some of the standard problems such as the lack of power in the transition region (e.g. Tinker et al. 2005; Mead et al. 2015), as can be seen in Fig. 6, and helps to a lesser extent with the general underestimation of power at smaller scales (Giocoli et al. 2010). Working at the level of the response may also alleviate problems that may arise from not using the most up-to-date ingredients for our halo model. For example, we use the mass function of Sheth & Tormen (1999) whereas there exist more recent mass functions such as Tinker et al. (2008, 2010) or Despali et al. (2016). We prefer Sheth & Tormen (1999) because it was calibrated to a wider range of cosmological parameter space than more modern mass functions.

5. Fitting halo model parameters

The level of agreement between the simulations and our model shown in Fig. 6 in encouraging, but demonstrates that the model is not sufficiently accurate to use to draw robust conclusions from forthcoming survey data. Note that the signal-to-noise is roughly 3:1 in the lensing–tSZ C() measurements of Hojjati et al. (2015) and is expected to be 5:1 for forthcoming KiDS measurements (Tröster et al., in prep.). It is possible that some of this modelling inaccuracy arises from incorrect choices for ingredients, but we think that a substantial amount must also arise from missing features in our basic halo-model calculation. Therefore, following the logic in Mead et al. (2015) we improve our model by fitting parameters that pertain to gas physics to the BAHAMAS simulations at the level of the power spectra. The aim is to find parameters that govern a population of “effective” haloes that describe the power spectra well for the different feedback strengths. We do this at the level of the halo-model “response”, discussed in Sect. 4. Once an accurate model for the response has been developed, we must then multiply this response by an accurate prescription for the matter–matter power spectrum in the gravity-only scenario. In our case we take the accurate prediction to be from HMCODE (Mead et al. 2015, 2016), but we could have used HALOFIT (Smith et al. 2003; Takahashi et al. 2012) or an emulator prediction (e.g. COSMIC EMU; Lawrence et al. 2010, 2017). Note that our model does not generally require us to fit parameters in this way, its utility is more general, fitting just seems to us to be the most obvious way for us to proceed to construct an accurate model.

We present several different models, the utility of each depending on the eventual use case: (1) stars, (2) matter, (3) matter and electron pressure and (4) matter, CDM, gas, and stars. Because the stars are a subdominant component of the total matter budget, we find it necessary to fit parameters that govern the stars separately (1), and then to hold these parameters fixed while we fitted models to the matter (2) and to the matter and electron pressure (3). In case (4), we fitted a model to all the matter fields simultaneously, refitting parameters that govern the star distribution. The free parameters that we have fitted are listed in Table 1 together with descriptions of their physical meaning and their default values. These parameters were chosen by trial-and-error to form minimal set that were able to model the data well, without opening the parameter space too widely.

Table 1.

Default hydrodynamical halo model parameters.

To fit we use the Nelder & Mead (1965) simplex algorithm with a large number of initial starting locations so as to avoid local minima. To determine the redshift dependence and relevant parameters, we initially fitted our model to different redshifts separately and we then used this information to parameterise sensible redshift dependencies for some of the parameters. The parameters are fitted to z <  1 BAHAMAS data across redshifts simultaneously with a linear weighting in z and a log weighting in k between 0.015 and 7 h Mpc−1. Our model currently takes ∼0.5 s per z to evaluate all 15 cross spectra we are interested in. The numerical bottle-neck is the computation of the Fourier transforms of the gas density and electron pressure profiles, which must be done numerically.

Before presenting our results we caution the reader against taking the details of our modelling as a serious physical description of the underlying haloes. We are fitting parameters in a halo model to power spectra taken from simulations; we are not fitting the simulations at the level of individual halo profiles. Therefore, the parameters of our hydrodynamical halo model should be thought of as pertaining to “effective” haloes that, through the apparatus of the halo model, provide an accurate model of the spectra we are interested in. They should not be confused with the exact parameters that govern the actual internal structure of actual haloes. Although they may relate to these, it would be necessary to check explicitly that this was so. We remind the reader of the shortcomings of the halo model approach: there may be features in the power spectra of the fields we are interested in that will never be accurately described by a linearly biased population of spherical, virialised haloes. Some of our parameter fitting may account for some of these deficiencies, although working with the response also helps to ameliorate some of these problems. Recall Fig. 3, where we saw that haloes of different masses are more-or-less important to different power spectra. This implies that we may drastically alter the profiles of haloes that contribute negligibly to a particular power spectrum and leave the spectrum almost unchanged. It would obviously be incorrect to over interpret the physical meaning of fitted parameters in this instance. Also remember that computing a power spectrum via the halo model represents a huge compression of the information that is potentially available from the halo profiles individually, which is smeared out across k by integrations. It is perfectly possible for inaccuracies in the modelling of the halo population to add or cancel in various unintuitive ways during these numerical compressions. Huffenberger & Seljak (2003) demonstrated that halo-structural parameters derived directly from haloes extracted from a cosmological simulation were not the same as those derived by fitting the structural parameters in a halo-model power spectrum calculation to the measured power spectrum from the same simulation. Mead et al. (2015) demonstrated that relatively drastic non-physical additions were required in order to improve halo model predictions for the matter–matter power spectrum from the 30% level to the 5% level, and that these were particularly important in the transition region. Tinker et al. (2005) and van den Bosch et al. (2013) both demonstrated that problems exist when trying to model galaxy clustering and galaxy–galaxy lensing using the halo model. To our knowledge, it has never been demonstrated that power spectra from simplistic halo-model calculations agree with those measured from simulations, even if the ingredients for the halo model calculation are taken directly from the simulation of interest.

5.1. Stars

The first model (1) we present is for the star–star power spectrum in Fig. 7. The free parameters are A*, M* and η with their best-fitting values listed in Table 2. We find that we are able to fit a reasonable model to the different simulations with different AGN feedback strengths (temperatures), and that our model follows the data with an accuracy of ≃5% for the response, nicely demarcating the different feedback scenarios. We note a preference for an inverse correlation between both A* and M* and the feedback temperature. This makes physical sense, since increased AGN activity is suppressing star formation. This could result in both an overall suppression in the number of stars in a halo of a given mass and also in pushing the peak star-formation efficiency to lower halo masses. In contrast, the best-fitting η, which governs the split between central and satellite stellar mass, is very similar for each of the different feedback temperatures. Note that there is a scatter in the efficiency of this suppression of star formation by AGN between simulations from different groups. Indeed a correct stellar mass–halo mass relation is one of the targets that hydro simulations strive for, and in practice this can be quite difficult to achieve. We also looked at additionally fitting the parameter σ*, but found that this did not produce a significant improvement in the quality of the fit.

thumbnail Fig. 7.

Halo model (1) for the star–star power spectrum is shown as sold lines while simulation measurements are shown as points with error bars that come from the scatter in power between three different realisations of the AGN BAHAMAS simulation. Note that the power spectra shown here is the measured response multiplied by HMCODE, so cosmic variance is removed at small scales, but we still show the error bar for comparison. Different colours denote the three different AGN feedback temperatures: 107.6 (blue), 107.8 (grey) and 108.0 K (red). Parameters for the model can be found in Table 2. We see that the spectrum is reduced in amplitude as feedback temperature increases, a consequence of AGN feedback suppressing star formation. This is reflected in the modelling as the value of A* and M* decreases as AGN strength increases.

Table 2.

Best-fitting effective halo model parameters for different combinations of power spectra from the BAHAMAS simulations.

5.2. Matter

In Fig. 8 we present model (2), in which we fitted the parameters ϵ1, Γ and M0 to the matter–matter data only. The underlying star model, which contributes to the matter–matter spectrum, is fixed to that shown in Fig. 7 and discussed in the previous subsection. We are able to match the power-spectrum response seen in the simulations at the percent level for each of the feedback models and at all of the redshifts shown. The fitted model has a strong correlation between M0 and AGN strength, which makes physical sense as stronger feedback ejects more halo gas. We also see strong trends in the fitted values of ϵ1 and Γ, which may be due to different amounts of back reaction on the dark matter and different heating of the residual halo gas. The higher AGN heating temperatures favours higher Γ, which corresponds to a less concentrated gas profile, as might be expected from more violent feedback. The upturn in the response function at small scales originates from the stellar contribution to the matter spectrum. In addition, in Fig. 8 we also show how the model fares when pushed beyond the range of AGN heating temperatures to which it was fitted. To do this we linearly interpolate and extrapolate the parameters of our halo model as a function of log TAGN/K. Indeed, this is the reason to use a physically motivated halo model in the first place, as opposed to a simple fitting function or a more blind interpolation between simulation results. This demonstrates some nice physical features of our model that may not be respected by more simplistic fitting functions (such as that in Mead et al. 2015). The suppression in the matter–matter power spectrum in our model originates because gas has been expelled from haloes, with the threshold mass increasing with the AGN heating temperature. Then, to a lesser extent, there are some additional effects from the non-NFW profile of the remaining gas and some back reaction on the concentration of the dark-matter halo component and then some effect from the very centrally concentrated stellar distribution which shows as an upturn at small scales ∼10 h Mpc−1. Because of this, there are some limits to how our model can behave: At one extreme we could raise the AGN heating temperature to a very high value that would result in almost all halo gas being expelled. The maximum suppression that this could have on the power would be to lower the amplitude of the one-halo term by (Ωcm)2 ≃ 0.7. At the other extreme, a low value of the AGN heating temperature would result in almost no gas expulsion and the matter–matter spectrum would be almost unaffected. These limits can be seen to be approached by the extreme temperature values shown in Fig. 8. We also looked at additionally fitting the parameters ϵ2 and β, but found that these did not produce a significant improvement in the quality of the fit.

thumbnail Fig. 8.

Similar to Fig. 7, but for the matter–matter spectrum (2). The power spectra are shown in the top row, while the response functions are shown in the bottom row. We see that the model is nearly perfect (errors at the sub-percent level) across a wide range of scales and for the different feedback scenarios. We see that the amount of power-spectrum suppression increases with the feedback temperature. The critical parameter that governs this in our modelling is the halo mass M0, which governs the halo mass below which more than half of the cosmic baryon density is missing from the halo. In the fitting parameters we see a strong preference for higher values of M0 for stronger feedback temperatures (from blue to red). How our model fares outside the range of AGN heating temperatures to which it was fitted is also shown via the lines that have no corresponding simulation points. To generate these we have interpolated and extrapolated our model parameters as a function of AGN heating temperature, which is shown in 0.1 dex intervals between 107.5 (bluest) and 108.1 K (reddest), indicating that our model is robust to extrapolation.

We also investigated the predictions of our model at fixed AGN heating temperature as we vary the underlying cosmology. As shown in van Daalen et al. (2020), the impact of AGN feedback on the matter–matter spectrum is quite insensitive to the difference between the Planck and WMAP9 cosmologies. We find a similar insensitivity in our model, which is reassuring, but we do find a small difference in the sense that the matter–matter power suppression in the Planck cosmology is predicted to be slightly less than that in WMAP9 at fixed AGN strength, which was shown by van Daalen et al. (2020). This presumably originates from the slightly different baryon fractions between the two different cosmologies, which means that slightly more halo mass is lost in the WMAP9 case.

5.3. Matter and electron pressure

In Fig. 9 we show the main result of this paper, model (3), which is a result of a joint fit to the matter–matter and matter–electron pressure spectra. This model could be used to predict both shear–shear and shear–y tomographic correlations. We do not consider the electron pressure auto-spectrum because while testing we discovered that this auto-spectrum is particularly difficult to fit; we suspect that this is because it receives much of its contributions from only a few of the most massive haloes present in the BAHAMAS simulation, and therefore that the spectra measured from the simulations will have very high Poisson noise contributions. Indeed, we see large variations in this power spectrum between the three different realisations of the BAHAMAS, and this variation dwarfs the differences between the feedback strengths for k <  2 h Mpc−1. The matter–electron pressure spectrum still suffers from this problem somewhat, as can be inferred from the error bars in the second row of Fig. 9, but the wavenumbers for which the variance between the 400 h−1 Mpc boxes is greater than the difference between the feedback models are k <  1 h Mpc−1. We see that we are able to simultaneously match the response functions at the few percent level for the wavenumbers at which the feedback scenarios are clearly demarcated. We see a trend that increasing AGN heating temperature causes an increase in α, the parameter that governs departures from simple virial temperature scaling, suggesting that gas that avoids being ejected is nevertheless heated significantly. However, the general trend is for suppressed small-scale power as AGN strength increases, which in our model arises from the fact that more gas is ejected from the more massive haloes, thus decreasing the pressure overall. In the matter–electron pressure case we note a trend that increasing AGN strength suppresses the power spectrum for scales dominated by the one-halo term; however, particularly at the higher z, we see that the power spectrum is relatively enhanced at larger scales. This could be because the heated gas is pushed out to larger scales by the feedback and this in turn gives an excess large-scale pressure contribution, either directly, or via shock heating the inter-galactic medium. In Fig. 9 we also show how our model for the matter and electron pressure fares when evaluated for AGN heating temperatures outside the range within which the model was constructed. As in Fig. 8 we see that the model behaves sensibly outside of the standard temperature range and thus we suggest that it can be used for cosmological and astrophysical parameter inference. We also looked at additionally fitting the parameters ϵ2 and β, but found that these did not produce a significant improvement in the quality of the fit.

thumbnail Fig. 9.

Same as Figs. 7 and 8, but for model (3) for the matter–matter (first and third row; power and response) and matter–electron pressure spectra (second and fourth row; power and response). The match to the matter–matter spectra is slightly degraded compared to the case when we fitted the matter–matter spectra exclusively, but we still match the response at the few percent level while simultaneously we match the matter–electron pressure spectrum well for k >  1 h Mpc−1. If we compare this to the error bars in the power spectra in the second row we see that our fit is good for scales where the three AGN heating temperature strengths are clearly demarcated in k. At larger scales we suspect the larger errors arise because the power is predominantly coming from few massive haloes. Similar to Fig. 8, we show how the model fares outside the range of AGN strengths to which is was fitted. AGN heating temperature is shown in 0.1 dex intervals between 107.5 K and 108.1 K, while the model was fitted between 107.6 K and 108.0 K, thus demonstrating that our model is robust to extrapolation.

5.4. Matter, CDM, gas and stars

In Fig. 10 we show a model that is the result of a joint fit to all possible combinations of matter, CDM, gas and star spectra. There are 10 distinct power spectra for these 4 fields. In our halo model, the matter field is the sum of the CDM, gas and stars, and so naturally gets upweighted in our fitting when it is included alongside its constituent parts. We see that our model gets reasonable matches to the auto-spectra shown in Fig. 10 and we have also checked that the cross-spectra are well matched (see Fig. 11), with the error always being roughly (but not exactly) the mean of the errors of the auto-spectra of the two contributing fields. We see that the gas spectrum is particularly well matched, which gives us confidence in our choice of the KS profile for the gas density profile. The match to the CDM–CDM response looks less good in Fig. 10, but note that the deviations from this being a constant are very small, of order a few percent. Physically this implies that the CDM is not deformed very significantly when going from a gravity only to a hydrodynamic model. Indeed, the suppression in the matter–matter power spectrum is driven predominantly by the (lack of) contribution from the gas, and then the upturn at smaller scales mainly derives from the contribution of the stars (see Fig. 6, as well as van Daalen et al. 2020). We also looked at additionally fitting the parameters ϵ2, β and σ*, but found that this did not produce a significant improvement in the quality of the fit.

thumbnail Fig. 10.

Same as Figs. 79, but for the best-fitting model (4) to all 10 possible power spectra of matter, CDM, gas and star fields, although we only show the 4 auto-spectra here. The top half of the figure show power spectra while the bottom half show response functions. Columns show different redshifts. Different colours are different feedback temperatures: T = 107.6 (blue), T = 107.8 (grey) and T = 108.0 K (red). Points (with errors) are from BAHAMAS while lines are from our halo model.

5.5. Model accuracy

In Fig. 11 we show the full error matrix for the AGN simulation for each of models (2), (3) and (4) but calculated for every possible power spectra of the 5 fields we have available. The matrix is 5 × 5 but only 15 of the 25 elements are independent due to symmetry. Each cell shows the relative difference for the model cross spectrum compared to the simulation, averaged linearly over z between 0 and 1 and logarithmically over k between ≃0.015 and 7 h Mpc−1. Error matrices for AGN-LO and AGN-HI are not shown, but are similar. We see that when we fitted only the matter–matter spectrum (2) we still have a reasonable match to the CDM, gas and stars spectra and combinations thereof. When we go from (2) to (3) and include the matter–electron pressure spectrum we improve the match to the matter–electron pressure, but at the cost of power spectra that involve the gas. This demonstrates a possible problem in our modelling in how we simultaneously model the gas density and pressure: When we go from case (2) to (4) and try to match matter, CDM, gas and stars power spectra, but ignore electron pressure we match the gas spectrum well, but this success comes at the expense of the pressure, which is ignored in this fit.

thumbnail Fig. 11.

Full error matrix for our fitted models calculated for each field pair for the 107.8 K AGN simulation. The error is averaged linearly over redshift between z = 0 and 1 and over k logarithmically between ≃0.015 and 7 h Mpc−1. This matrix is symmetric because power spectra are symmetric when swapping the field labels. We show the matrix for the three models presented in this work: matter (2, left); matter, electron pressure (3, centre); matter, CDM, gas, stars (4, right). The squares that contain dots are the particular cross-spectra that are fitted in the model. In the left-hand plot, we see that in the matter model we get a reasonable fit to the power spectra all of the constituent fields. In the central plot we get a better match to to the matter–electron pressure spectrum than in the left-hand panel, but this is at the expense of any power spectra that directly involves the gas field, which shows some tension in our model between the gas and pressure modelling. In the right-hand panel, where we fitted matter, CDM, gas and stars, we see that all the matter spectra are reasonable, and better than in the left-hand panel, but achieving this is at the expense of any spectra that involves the pressure field.

6. Discussion

6.1. Summary

We have presented a hydrodynamical halo model that can be used to make analytical calculations of power spectra for combinations of the matter, CDM, gas and star density fields, as well as for cross-spectra of these fields with the electron-pressure field. This model will be used in the future to provide constraints on both cosmology and feedback physics using a combination of lensing and tSZ data. The model has a number of tuneable parameters that pertain to gas physics. We provided four possible sets of values for these parameters that provide good fits to power spectra for: (1) stars, (2) matter, (3) matter and electron pressure and (4) matter, CDM, gas and stars. Parameters of these models are all fitted to power spectra measured from the BAHAMAS simulations for three different AGN feedback strengths, which in this case is parameterised by a single sub-grid parameter: the AGN heating temperature. By interpolating between model parameters as a function of AGN heating temperature we suggest that our model can be used in a cosmological analysis and that it is robust to extrapolation. In particular, model (2) will be useful in lensing-only studies as it is able to reproduce the deficit in matter–matter power seen in hydrodynamic simulations, compared to the DMONLY case, at the percent level. Model (3) provides a slightly worse match to the mater–matter power (accurate to 2%), but captures the temperature variations in the matter–electron pressure spectrum that can be directly integrated to calculate the lensing–tSZ cross correlation (accurate to 15%) as well as the lensing–lensing correlation. This accuracy is sufficient for forthcoming lensing and tSZ data sets (Tröster et al., in prep.).

6.2. Identified problems

6.2.1. Electron pressure errors

Our model is less accurate for matter–electron pressure (15%) compared to matter–matter (2%). We suspect that this is partly because we report the error averaged over log k from 0.015 to 7 h Mpc−1 and linearly over z from 0 to 1. The matter–electron pressure response is noisier at large scales compared to the matter spectra, and this introduces an unavoidable source of error for our necessarily smooth halo model prediction. We may also be seeing the fact that the relatively small box size of BAHAMAS means that the response (or power) has not converged for k <  1 h Mpc−1, and it is possible that our model would provide a better fit to data from a significantly larger simulation.

6.2.2. Response approach

In this paper, we fitted to the power spectrum response, which we created by dividing a measured hydrodynamical power spectrum by a gravity-only spectrum. This response cancels the Gaussian error in power spectra at large scales for the matter spectra, but it may create complicated correlations at smaller scales. In particular, this response may not be appropriate for the electron-pressure spectra, even at large scales, since the large-scale variance does not mimic that of the matter spectra. We fitted parameters of our models to the response using a least-squares approach, with no k or z dependent error. In the future we hope to have access to more and larger simulations that we would need in order investigate this more thoroughly and so that we may derive an appropriate error bar.

6.2.3. Electron pressure auto spectrum

While completing this work we tried to generate a model for the electron pressure auto-spectrum, together with matter–matter and matter–electron pressure; but this proved to be difficult and was abandoned. We discovered that the realisation-to-realisation variance in the electron pressure auto-spectrum measured from BAHAMAS was huge, and dwarfed the difference between feedback strengths for k <  2 h Mpc−1, which left us with very little signal to fit. This in turn suggests that boxes of much greater volume than the BAHAMAS (400 h−1 Mpc)3 must be used to ascertain the feedback dependence of this spectrum on large scales. Therefore, our model should not be used to compute tSZ–tSZ power spectrum at this stage. However, with larger hydrodynamical simulations being run in the future, our approach should still be able to be used to match the power spectrum by refitting parameters. We note that it would also be useful to have more accurate calibration data for the simulations to match, as this is one of the limiting factors for BAHAMAS at the moment.

6.2.4. Joint model

We also tried to generate a joint model for matter, CDM, gas, stars and electron pressure, which has not been presented in this paper because it did not work well. We found that it was particularly difficult to simultaneously model spectra that involved the electron pressure and those that involve the gas density with the same underlying physical model. In the current modelling, these two fields are joined in two different ways: first via M0, which determines the overall gas abundance as a function of halo mass, and second via Γ, the polytropic index of the KS gas profile, which jointly determines the density and pressure profiles. The KS profile is the correct model for undisturbed gas that is in hydrostatic equilibrium in an NFW halo, but this state will clearly be violated in the presence of AGN feedback. It was shown in Yan et al. (2020) that the KS model could provide a reasonable match to gas-density profiles in the BAHAMAS simulations, but it is likely that the polytropic link between density and pressure is not respected in the presence of feedback (Battaglia et al. 2012a,b). It is probable that a better model could be generated by either changing the gas profile, for example by explicitly including heated gas, or by weakening the links between the density and pressure modelling. However, note that some link is necessary if one wants to simultaneously constrain feedback and cosmology with cross-correlation data. It is also possible that some of the inaccuracy stems from our treatment of non-thermal pressure within clusters. We have a single multiplicative parameter, α, which increases or decreases gas temperature independent of halo mass. In future it may be useful to investigate physical models for non-thermal pressure in clusters, such as those presented in Nelson et al. (2014) or Shi et al. (2015).

6.2.5. Non-linear halo bias

The halo model provides a good description of the power spectra for two-halo dominated regimes and the one-halo dominated regimes, but it is frustrating that it fails in the transition region between these two regimes. In retrospect, the reason for this is obvious: The assumption made is that haloes are linearly tracing an underlying linear matter field, and this assumption clearly breaks down at scales that correspond to this transition: k ∼ 0.5 h Mpc−1 at z = 0 for matter–matter. Linear bias is not valid at this wavenumber (Desjacques et al. 2018) and the two-halo term predicted by the model at these and smaller scales will be, quite simply, wrong. We suggest that a fruitful line of future inquiry would be to improve the treatment of halo bias within the standard halo-model framework (e.g. Smith et al. 2007). We have somewhat sidestepped this issue in this paper by fitting the BAHAMAS response functions with halo-model responses, rather than the power spectra directly. This is certainly of benefit because we cancel the Gaussian noise present at large-scales in the simulations. This makes the halo model response an almost perfect tool for understanding matter–matter power response for dark energy models measured in simulations without hydrodynamics (Mead 2017; Cataneo et al. 2019) when the underlying linear spectrum is fixed. This would also be true in hydrodynamic simulations if we could be assured that the troubles in the transition region manifest themselves in the same way, and with the same scale dependence, for other power spectra as they do for matter–matter. Unfortunately, from Fig. 2 we can see that this is not exactly the case, since the transition region occurs at slightly different scales for the different power spectra (particularly the scale for matter–electron pressure is quite different from matter–matter): there simply cannot be some universal k-dependent prescription that one could apply to solve all problems in all cases. It is therefore possible that the response is not the optimal quantity to consider when one considers the different power spectra measured in hydrodynamic simulations.

6.2.6. Expelled gas modelling

One further shortcoming in our modelling is in our treatment of the gas expelled from haloes. In this paper we do not consider any sort of one-halo contribution arising from expelled gas, and instead account for the effects of this gas only in the two-halo term. In reality, we know that gas expelled from a halo will still be correlated with that halo, and will be found in an extended shell outside the virial radius. Since this shell is larger than the halo virial radius, it may provide a significant contribution to the power at scales that correspond to the transition between the two- and one-halo terms in the classic halo model, and we have so far ignored this. In the early stages of this work we attempted to explicitly include a gas component outside the halo virial radius, but we found that this introduces some problems with the standard halo-model formalism: In the calculation, it is necessary to truncate the density profile of each NFW halo at the virial radius, otherwise they each have an infinite mass. If the gas halo is allowed to extend beyond the virial radius then we introduce an explicit anti-correlation between gas and CDM because we have gas in regions where we explicitly have no CDM. This anti-correlation allows the CDM–gas cross spectrum to be negative in some regions of parameter space. While negative cross-spectra are not problematic in general, in cosmology it would be unrealistic due to the general structure of the cosmic web. Indeed, such a negative cross spectrum is not seen in simulations. It may be possible to improve this by using a different CDM halo profile with a finite mass as the radius tends to infinity, but this introduces difficulty in deciding which mass to associate with which halo when identifying haloes in simulations. The fact that we treat the ejected gas in the two-halo term as a simple multiplicative bias also generates another problem: the two-halo term for ejected gas is not suppressed at small scales relative to the shape of linear theory. For other fields, it is always suppressed by a weight that depends on the halo profile and halo bias (Eq. (1)). The details of this two-halo term suppression cannot be taken seriously because they occur at scales where the halo bias is no longer linear, however, this suppression does have the virtue of ensuring that the two-halo term is always much smaller than the one-halo term at small scales. Because the gas profiles themselves lack small-scale structure, the two-halo term for gas (and electron pressure) can be larger than the one-halo term at very small scales (k ∼ 100 h Mpc−1). Clearly this is unphysical and would need to be addressed in further work.

6.2.7. Important halo masses

A further potential problem for this approach is the differing halo masses that are important for the tSZ effect compared to those important for galaxy weak lensing. As shown in Fig. 3, the matter–electron pressure spectrum derives most contribution from higher-mass haloes than matter–matter. It is therefore these haloes that we learn most about from the lensing-tSZ cross correlation. However, it is lower-mass haloes that mostly affect the suppression in the power spectrum that is important for the lensing-lensing spectrum. A fruitful avenue for future research may be to learn more directly about the gas distribution from the haloes that are most relevant, those of 1013 to 1014h−1M. It is difficult to do this using the tSZ-lensing cross correlation directly, and instead one may have to use some sort of clipping (e.g. Simpson et al. 2011, 2013; Hill & Pajer 2013; Giblin et al. 2018) technique to mask the influence of the highest mass haloes. We note that there is already a great deal known about the impact of baryons on the matter distribution from tSZ and X-ray studies of galaxies, groups, and clusters from targeted observational studies. One could imagine using the model presented in this paper to look at both the lensing and tSZ signal (or lensing–tSZ) around haloes of specific masses. This would add information about the signal as it arises from different mass bins, somehow akin to the information in Fig. 3. However, this would necessitate having reliable measurements of individual halo masses and would require modelling effects such as halo mis-centring (e.g. Yan et al. 2020). It would also need to be ascertained whether the halo model can be trusted for calculations that are binned in halo mass and how problems with these calculations compare to problems in the standard when integrated over all halo masses.

6.3. Conclusion

We intend to use this model in the near future to analyse the cross correlation between the tSZ map from Planck and weak lensing from the KiDS data, as well as the lensing auto correlation data (Tröster et al., in prep). This work has the potential to simultaneously learn about cosmology and the magnitude of AGN feedback in the Universe and therefore to inform galaxy formation modelling and future hydrodynamical simulations as well as to shrink error bars on cosmological parameters. As alluded to in this paper, the tSZ data itself also has the power to provide tight constraints on σ8 and initial tests show that this could break degeneracies that exist in a traditional lensing analysis.

So far, we have only presented the best-fitting model parameters for each scenario and we have not investigated errorson, or degeneracies between, our model parameters, but we hope to be able to do this in future. We note that this would require us to create meaningful error bars on our response measurements.

It is clear that the fitted parameters presented here are very BAHAMAS specific, which may worry some readers who would prefer a model that could encompass a wider range of hydrodynamical simulations. However, we think that our model is general enough that it could be refitted to other simulations in the future. We also note van Daalen et al. (2020) has demonstrated a strong correlation between the matter–matter response function and the gas fraction in haloes of ∼1014h−1M and that the BAHAMAS simulations reproduce observations of gas fractions in haloes of this mass to within the observational error bar. Given that the gas content in haloes of this mass also is important for tSZ suggests that tuning our model to BAHAMAS was a reasonable first choice. In the future, it may be possible to rewrite our results in terms of more physical quantities, such as the halo gas fraction (e.g. Chisari et al. 2019; van Daalen et al. 2020), rather than the unphysical AGN heating temperature. It would also be desirable to confront our model with different hydrodynamical simulations, run with different underlying assumptions and methods, for example those run with the Adaptive Refinement Tree code (Kravtsov et al. 2002; Nagai et al. 2007) or those of ILLUSTRIS (Vogelsberger et al. 2014).

In the future, it may be possible to use this halo model to make predictions for higher-order statistics and these may lead to increased constraints on the model parameters. We also envisage the model being used to make predictions for other “hydrodynamical” observable quantities, such as n-point statistics that involve the X-ray field, which is a different function of gas temperature and density compared to the tSZ and provides complimentary information about the state of gas in the universe. The utility of the halo model that we present does not require its parameters to be fitted to simulation data, this was simply one option, and the option that we decided to pursue in this paper. In the future the model can be used to provide quick estimates for the effect of new feedback scenarios that have yet to be simulated, or to explore non-standard cosmological scenarios that have yet to be hydrodynamical simulated. It would even be possible to investigate the non-linear coupling of feedback with non-standard cosmological scenarios.


1

Often simulations that only consider the gravitational interaction are called “dark-matter only”. We prefer “gravity only” because these simulations do contain baryons, and the initial conditions are set up to account for baryons, but they exclusively consider the gravitational interaction between particles when they are subsequently run.

3

We use matter overdensity but straight electron pressure, not the pressure overdensity, because these are the exact quantities probed by lensing and tSZ respectively.

4

Achieving these limits is difficult numerically because of the large amount of mass contained in low-mass haloes according to most popular mass functions. Special care must be taken with the two-halo integral in the case of power spectra that involve the matter field. See Appendix A.

5

The non-standard factor of Ωm(z) in the denominator of Eq. (22) arises because we work with respect to the matter density, rather than critical density.

6

In the unusual situation that the halo mass fraction in bound gas plus that in stars would be greater than the universal baryon fraction we subtract the excess mass from the stars. This only happens for very high halo masses or unusual sets of parameters.

7

Numerically this is achieved by setting the ejected gas profile to zero in the one-halo term, and to a delta function in the two-halo term.

8

If the gas were comprised only of ionised hydrogen and helium, with a hydrogen mass fraction fH, we would have μp = 4/(3 + 5fH) and μe = 2/(1 + fH). For the BAHAMAS with fH = 0.752, μe ≃ 0.59 and μe ≃ 1.14. However, we adopt the values μp = 0.61 and μe = 1.17 to be compatible with gas metallicity in the simulations and the way the electron pressure field was measured in post processesing.

Acknowledgments

AJM acknowledges support from a CITA National Fellowship and (together with LVW) from NSERC. AJM and TT have received funding from the Horizon 2020 research and innovation programme of the European Union under Marie Skłodowska-Curie grant agreements Nos. 702971 (AJM) and 797794 (TT). TT and CH acknowledge support from the European Research Council under grant number 647112, while IGM acknowledges the same funding under grant number 769130. CH also acknowledges 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. We thank an anonymous referee for useful comments that improved the quality of this paper.

References

  1. Addison, G. E., Dunkley, J., & Spergel, D. N. 2012, MNRAS, 427, 1741 [NASA ADS] [CrossRef] [Google Scholar]
  2. Addison, G. E., Dunkley, J., & Bond, J. R. 2013, MNRAS, 436, 1896 [Google Scholar]
  3. Agarwal, S., Abdalla, F. B., Feldman, H. A., Lahav, O., & Thomas, S. A. 2012, MNRAS, 424, 1409 [CrossRef] [Google Scholar]
  4. Angulo, R. E., Zennaro, M., Contreras, S., et al. 2020, ArXiv e-prints [arXiv:2004.06245] [Google Scholar]
  5. Aricò, G., Angulo, R. E., Hernández-Monteagudo, C., et al. 2020, MNRAS, 495, 4800 [CrossRef] [Google Scholar]
  6. Arnaud, M., Pratt, G. W., Piffaretti, R., et al. 2010, A&A, 517, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  7. Bartelmann, M., & Schneider, P. 2001, Phys. Rep., 340, 291 [NASA ADS] [CrossRef] [Google Scholar]
  8. Battaglia, N., Bond, J. R., Pfrommer, C., Sievers, J. L., & Sijacki, D. 2010, ApJ, 725, 91 [NASA ADS] [CrossRef] [Google Scholar]
  9. Battaglia, N., Bond, J. R., Pfrommer, C., & Sievers, J. L. 2012a, ApJ, 758, 74 [NASA ADS] [CrossRef] [Google Scholar]
  10. Battaglia, N., Bond, J. R., Pfrommer, C., & Sievers, J. L. 2012b, ApJ, 758, 75 [NASA ADS] [CrossRef] [Google Scholar]
  11. Battaglia, N., Hill, J. C., & Murray, N. 2015, ApJ, 812, 154 [CrossRef] [Google Scholar]
  12. Baxter, E. J., Omori, Y., Chang, C., et al. 2019, Phys. Rev. D, 99, 023508 [NASA ADS] [CrossRef] [Google Scholar]
  13. Bird, S., Feng, Y., Pedersen, C., & Font-Ribera, A. 2020, JCAP, 2020, 002 [CrossRef] [Google Scholar]
  14. Bolliet, B., Comis, B., Komatsu, E., & Macías-Pérez, J. F. 2018, MNRAS, 477, 4957 [NASA ADS] [CrossRef] [Google Scholar]
  15. Bond, J. R., Contaldi, C. R., Pen, U.-L., et al. 2005, ApJ, 626, 12 [NASA ADS] [CrossRef] [Google Scholar]
  16. Bryan, G. L., & Norman, M. L. 1998, ApJ, 495, 80 [NASA ADS] [CrossRef] [Google Scholar]
  17. Cacciato, M., Lahav, O., van den Bosch, F. C., Hoekstra, H., & Dekel, A. 2012, MNRAS, 426, 566 [NASA ADS] [CrossRef] [Google Scholar]
  18. Cataneo, M., Lombriser, L., Heymans, C., et al. 2019, MNRAS, 488, 2121 [NASA ADS] [CrossRef] [Google Scholar]
  19. Cen, R., & Ostriker, J. P. 1999, ApJ, 514, 1 [NASA ADS] [CrossRef] [Google Scholar]
  20. Chisari, N. E., Richardson, M. L. A., Devriendt, J., et al. 2018, MNRAS, 480, 3962 [NASA ADS] [CrossRef] [Google Scholar]
  21. Chisari, N. E., Mead, A. J., Joudaki, S., et al. 2019, Open J. Astrophys., 2, 4 [CrossRef] [Google Scholar]
  22. Cooray, A., & Sheth, R. 2002, Phys. Rep., 372, 1 [NASA ADS] [CrossRef] [Google Scholar]
  23. Copeland, D., Taylor, A., & Hall, A. 2018, MNRAS, 480, 2247 [CrossRef] [Google Scholar]
  24. Courtin, J., Rasera, Y., Alimi, J.-M., et al. 2011, MNRAS, 410, 1911 [NASA ADS] [Google Scholar]
  25. Dai, B., Feng, Y., & Seljak, U. 2018, JCAP, 11, 009 [CrossRef] [Google Scholar]
  26. Debackere, S. N. B., Schaye, J., & Hoekstra, H. 2020, MNRAS, 492, 2285 [CrossRef] [Google Scholar]
  27. de Graaff, A., Cai, Y.-C., Heymans, C., & Peacock, J. A. 2019, A&A, 624, A48 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  28. Desjacques, V., Jeong, D., & Schmidt, F. 2018, Phys. Rep., 733, 1 [NASA ADS] [CrossRef] [Google Scholar]
  29. Despali, G., Giocoli, C., Angulo, R. E., et al. 2016, MNRAS, 456, 2486 [NASA ADS] [CrossRef] [Google Scholar]
  30. Dolag, K., Borgani, S., Murante, G., & Springel, V. 2009, MNRAS, 399, 497 [Google Scholar]
  31. Dolag, K., Komatsu, E., & Sunyaev, R. 2016, MNRAS, 463, 1797 [CrossRef] [Google Scholar]
  32. Duffy, A. R., Schaye, J., Kay, S. T., & Dalla Vecchia, C. 2008, MNRAS, 390, L64 [NASA ADS] [CrossRef] [Google Scholar]
  33. Eifler, T., Krause, E., Dodelson, S., et al. 2015, MNRAS, 454, 2451 [NASA ADS] [CrossRef] [Google Scholar]
  34. Fedeli, C. 2014, JCAP, 4, 28 [Google Scholar]
  35. Fedeli, C., Dolag, K., & Moscardini, L. 2012, MNRAS, 419, 1588 [NASA ADS] [CrossRef] [Google Scholar]
  36. Fedeli, C., Semboloni, E., Velliscig, M., et al. 2014, JCAP, 8, 28 [Google Scholar]
  37. Foreman, S., Becker, M. R., & Wechsler, R. H. 2016, MNRAS, 463, 3326 [NASA ADS] [CrossRef] [Google Scholar]
  38. Giblin, B., Heymans, C., Harnois-Déraps, J., et al. 2018, MNRAS, 480, 5529 [NASA ADS] [CrossRef] [Google Scholar]
  39. Giocoli, C., Bartelmann, M., Sheth, R. K., & Cacciato, M. 2010, MNRAS, 408, 300 [NASA ADS] [CrossRef] [Google Scholar]
  40. Gonzalez, A. H., Sivanandam, S., Zabludoff, A. I., & Zaritsky, D. 2013, ApJ, 778, 14 [NASA ADS] [CrossRef] [Google Scholar]
  41. Harnois-Déraps, J., van Waerbeke, L., Viola, M., & Heymans, C. 2015, MNRAS, 450, 1212 [NASA ADS] [CrossRef] [Google Scholar]
  42. Hildebrandt, H., Viola, M., Heymans, C., et al. 2017, MNRAS, 465, 1454 [Google Scholar]
  43. Hildebrandt, H., Köhlinger, F., van den Busch, J. L., et al. 2020, A&A, 633, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  44. Hill, J. C., & Pajer, E. 2013, Phys. Rev. D, 88, 063526 [NASA ADS] [CrossRef] [Google Scholar]
  45. Hill, J. C., & Spergel, D. N. 2014, JCAP, 2, 030 [NASA ADS] [CrossRef] [Google Scholar]
  46. Hill, J. C., Baxter, E. J., Lidz, A., Greco, J. P., & Jain, B. 2018, Phys. Rev. D, 97, 083501 [CrossRef] [Google Scholar]
  47. Hinshaw, G., Larson, D., Komatsu, E., et al. 2013, ApJS, 208, 19 [NASA ADS] [CrossRef] [Google Scholar]
  48. Hojjati, A., McCarthy, I. G., Harnois-Deraps, J., et al. 2015, JCAP, 10, 047 [CrossRef] [Google Scholar]
  49. Hojjati, A., Tröster, T., Harnois-Déraps, J., et al. 2017, MNRAS, 471, 1565 [CrossRef] [Google Scholar]
  50. Holder, G. P., McCarthy, I. G., & Babul, A. 2007, MNRAS, 382, 1697 [NASA ADS] [Google Scholar]
  51. Horowitz, B., & Seljak, U. 2017, MNRAS, 469, 394 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Huang, H.-J., Eifler, T., Mandelbaum, R., & Dodelson, S. 2019, MNRAS, 488, 1652 [NASA ADS] [CrossRef] [Google Scholar]
  53. Huffenberger, K. M., & Seljak, U. 2003, MNRAS, 340, 1199 [NASA ADS] [CrossRef] [Google Scholar]
  54. Itoh, N., Kohyama, Y., & Nozawa, S. 1998, ApJ, 502, 7 [NASA ADS] [CrossRef] [Google Scholar]
  55. Jain, B., & Seljak, U. 1997, ApJ, 484, 560 [NASA ADS] [CrossRef] [Google Scholar]
  56. Joudaki, S., Blake, C., Heymans, C., et al. 2017a, MNRAS, 465, 2033 [NASA ADS] [CrossRef] [Google Scholar]
  57. Joudaki, S., Mead, A., Blake, C., et al. 2017b, MNRAS, 471, 1259 [NASA ADS] [CrossRef] [Google Scholar]
  58. Kaiser, N. 1992, ApJ, 388, 272 [NASA ADS] [CrossRef] [Google Scholar]
  59. Kilbinger, M. 2015, Rep. Prog. Phys., 78, 086901 [Google Scholar]
  60. Kilbinger, M., Heymans, C., Asgari, M., et al. 2017, MNRAS, 472, 2126 [NASA ADS] [CrossRef] [Google Scholar]
  61. Knabenhans, M., Stadel, J., Marelli, S., et al. 2019, MNRAS, 484, 5509 [NASA ADS] [CrossRef] [Google Scholar]
  62. Komatsu, E., & Seljak, U. 2001, MNRAS, 327, 1353 [Google Scholar]
  63. Komatsu, E., & Seljak, U. 2002, MNRAS, 336, 1256 [NASA ADS] [CrossRef] [Google Scholar]
  64. Kravtsov, A. V., Klypin, A., & Hoffman, Y. 2002, ApJ, 571, 563 [NASA ADS] [CrossRef] [Google Scholar]
  65. Kravtsov, A. V., Vikhlinin, A. A., & Meshcheryakov, A. V. 2018, Astron. Lett., 44, 8 [NASA ADS] [CrossRef] [Google Scholar]
  66. Lawrence, E., Heitmann, K., White, M., et al. 2010, ApJ, 713, 1322 [NASA ADS] [CrossRef] [Google Scholar]
  67. Lawrence, E., Heitmann, K., Kwan, J., et al. 2017, ApJ, 847, 50 [CrossRef] [Google Scholar]
  68. Leauthaud, A., Tinker, J., Bundy, K., et al. 2012, ApJ, 744, 159 [NASA ADS] [CrossRef] [Google Scholar]
  69. Le Brun, A. M. C., McCarthy, I. G., Schaye, J., & Ponman, T. J. 2014, MNRAS, 441, 1270 [NASA ADS] [CrossRef] [Google Scholar]
  70. Le Brun, A. M. C., McCarthy, I. G., & Melin, J.-B. 2015, MNRAS, 451, 3868 [NASA ADS] [CrossRef] [Google Scholar]
  71. Lemos, P., Challinor, A., & Efstathiou, G. 2017, JCAP, 5, 014 [CrossRef] [Google Scholar]
  72. Lesgourgues, J. 2011, ArXiv e-prints [arXiv:1104.2934] [Google Scholar]
  73. Lo Verde, M., Miller, A., Shandera, S., & Verde, L. 2008, JCAP, 4, 14 [Google Scholar]
  74. Ma, Y.-Z., Van Waerbeke, L., Hinshaw, G., et al. 2015, JCAP, 9, 046 [NASA ADS] [CrossRef] [Google Scholar]
  75. MacCrann, N., Zuntz, J., Bridle, S., Jain, B., & Becker, M. R. 2015, MNRAS, 451, 2877 [NASA ADS] [CrossRef] [Google Scholar]
  76. MacCrann, N., Aleksić, J., Amara, A., et al. 2017, MNRAS, 465, 2567 [CrossRef] [Google Scholar]
  77. Makiya, R., Ando, S., & Komatsu, E. 2018, MNRAS, 480, 3928 [NASA ADS] [CrossRef] [Google Scholar]
  78. Martizzi, D., Teyssier, R., Moore, B., & Wentz, T. 2012, MNRAS, 422, 3081 [NASA ADS] [CrossRef] [Google Scholar]
  79. Martizzi, D., Teyssier, R., & Moore, B. 2013, MNRAS, 432, 1947 [NASA ADS] [CrossRef] [Google Scholar]
  80. McCarthy, I. G., Schaye, J., Ponman, T. J., et al. 2010, MNRAS, 406, 822 [NASA ADS] [Google Scholar]
  81. McCarthy, I. G., Schaye, J., Bower, R. G., et al. 2011, MNRAS, 412, 1965 [NASA ADS] [CrossRef] [Google Scholar]
  82. McCarthy, I. G., Le Brun, A. M. C., Schaye, J., & Holder, G. P. 2014, MNRAS, 440, 3645 [NASA ADS] [CrossRef] [Google Scholar]
  83. McCarthy, I. G., Schaye, J., Bird, S., & Le Brun, A. M. C. 2017, MNRAS, 465, 2936 [NASA ADS] [CrossRef] [Google Scholar]
  84. McCarthy, I. G., Bird, S., Schaye, J., et al. 2018, MNRAS, 476, 2999 [NASA ADS] [CrossRef] [Google Scholar]
  85. Mead, A. J. 2017, MNRAS, 464, 1282 [CrossRef] [Google Scholar]
  86. Mead, A. J., Peacock, J. A., Heymans, C., Joudaki, S., & Heavens, A. F. 2015, MNRAS, 454, 1958 [NASA ADS] [CrossRef] [Google Scholar]
  87. Mead, A. J., Heymans, C., Lombriser, L., et al. 2016, MNRAS, 459, 1468 [Google Scholar]
  88. Mo, H. J., & White, S. D. M. 1996, MNRAS, 282, 347 [NASA ADS] [CrossRef] [Google Scholar]
  89. Mohammed, I., Martizzi, D., Teyssier, R., & Amara, A. 2014, ArXiv e-prints [arXiv:1410.6826] [Google Scholar]
  90. Moster, B. P., Naab, T., & White, S. D. M. 2013, MNRAS, 428, 3121 [Google Scholar]
  91. Mummery, B. O., McCarthy, I. G., Bird, S., & Schaye, J. 2017, MNRAS, 471, 227 [CrossRef] [Google Scholar]
  92. Nagai, D., Kravtsov, A. V., & Vikhlinin, A. 2007, ApJ, 668, 1 [NASA ADS] [CrossRef] [Google Scholar]
  93. Nakamura, T. T., & Suto, Y. 1997, Prog. Theor. Phys., 97, 49 [NASA ADS] [CrossRef] [Google Scholar]
  94. Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [NASA ADS] [CrossRef] [Google Scholar]
  95. Nelder, J. A., & Mead, R. 1965, Comput. J., 7, 308 [Google Scholar]
  96. Nelson, K., Lau, E. T., & Nagai, D. 2014, ApJ, 792, 25 [NASA ADS] [CrossRef] [Google Scholar]
  97. Nozawa, S., Itoh, N., Kawana, Y., & Kohyama, Y. 2000, ApJ, 536, 31 [NASA ADS] [CrossRef] [Google Scholar]
  98. Omori, Y., Giannantonio, T., Porredon, A., et al. 2019, Phys. Rev. D, 100, 043501 [CrossRef] [Google Scholar]
  99. Osato, K., Flender, S., Nagai, D., Shirasaki, M., & Yoshida, N. 2018, MNRAS, 475, 532 [CrossRef] [Google Scholar]
  100. Osato, K., Shirasaki, M., Miyatake, H., et al. 2020, MNRAS, 492, 4780 [CrossRef] [Google Scholar]
  101. Peacock, J. A., & Smith, R. E. 2000, MNRAS, 318, 1144 [NASA ADS] [CrossRef] [Google Scholar]
  102. Planck Collaboration XXII. 2016, A&A, 594, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  103. Planck Collaboration XXIV. 2016, A&A, 594, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  104. Puchwein, E., & Springel, V. 2013, MNRAS, 428, 2966 [NASA ADS] [CrossRef] [Google Scholar]
  105. Rabold, M., & Teyssier, R. 2017, MNRAS, 467, 3188 [NASA ADS] [CrossRef] [Google Scholar]
  106. Refregier, A., & Teyssier, R. 2002, Phys. Rev. D, 66, 043002 [NASA ADS] [CrossRef] [Google Scholar]
  107. Remazeilles, M., Bolliet, B., Rotti, A., & Chluba, J. 2019, MNRAS, 483, 3459 [NASA ADS] [CrossRef] [Google Scholar]
  108. Roncarelli, M., Moscardini, L., Tozzi, P., et al. 2006, MNRAS, 368, 74 [NASA ADS] [CrossRef] [Google Scholar]
  109. Roncarelli, M., Moscardini, L., Borgani, S., & Dolag, K. 2007, MNRAS, 378, 1259 [NASA ADS] [CrossRef] [Google Scholar]
  110. Rudd, D. H., Zentner, A. R., & Kravtsov, A. V. 2008, ApJ, 672, 19 [NASA ADS] [CrossRef] [Google Scholar]
  111. Schaye, J., Dalla Vecchia, C., Booth, C. M., et al. 2010, MNRAS, 402, 1536 [Google Scholar]
  112. Schmidt, F. 2016, Phys. Rev. D, 93, 063512 [CrossRef] [Google Scholar]
  113. Schneider, A., & Teyssier, R. 2015, JCAP, 12, 049 [NASA ADS] [CrossRef] [Google Scholar]
  114. Schneider, A., Teyssier, R., Potter, D., et al. 2016, JCAP, 2016, 047 [NASA ADS] [CrossRef] [Google Scholar]
  115. Schneider, A., Teyssier, R., Stadel, J., et al. 2019, JCAP, 2019, 020 [CrossRef] [Google Scholar]
  116. Scoccimarro, R., Sheth, R. K., Hui, L., & Jain, B. 2001, ApJ, 546, 20 [NASA ADS] [CrossRef] [Google Scholar]
  117. Seljak, U. 2000, MNRAS, 318, 203 [NASA ADS] [CrossRef] [Google Scholar]
  118. Semboloni, E., Hoekstra, H., Schaye, J., van Daalen, M. P., & McCarthy, I. G. 2011, MNRAS, 417, 2020 [NASA ADS] [CrossRef] [Google Scholar]
  119. Semboloni, E., Hoekstra, H., & Schaye, J. 2013, MNRAS, 434, 148 [NASA ADS] [CrossRef] [Google Scholar]
  120. Shaw, L. D., Nagai, D., Bhattacharya, S., & Lau, E. T. 2010, ApJ, 725, 1452 [NASA ADS] [CrossRef] [Google Scholar]
  121. Sheth, R. K., & Tormen, G. 1999, MNRAS, 308, 119 [NASA ADS] [CrossRef] [Google Scholar]
  122. Sheth, R. K., Mo, H. J., & Tormen, G. 2001, MNRAS, 323, 1 [NASA ADS] [CrossRef] [Google Scholar]
  123. Shi, X., Komatsu, E., Nelson, K., & Nagai, D. 2015, MNRAS, 448, 1020 [NASA ADS] [CrossRef] [Google Scholar]
  124. Sijacki, D., Vogelsberger, M., Genel, S., et al. 2015, MNRAS, 452, 575 [NASA ADS] [CrossRef] [Google Scholar]
  125. Simpson, F., James, J. B., Heavens, A. F., & Heymans, C. 2011, Phys. Rev. Let., 107, 271301 [NASA ADS] [CrossRef] [Google Scholar]
  126. Simpson, F., Heavens, A. F., & Heymans, C. 2013, Phys. Rev. D, 88, 083510 [NASA ADS] [CrossRef] [Google Scholar]
  127. Smith, R. E., & Markovic, K. 2011, Phys. Rev. D, 84, 063507 [CrossRef] [Google Scholar]
  128. Smith, R. E., Peacock, J. A., Jenkins, A., et al. 2003, MNRAS, 341, 1311 [Google Scholar]
  129. Smith, R. E., Scoccimarro, R., & Sheth, R. K. 2007, Phys. Rev. D, 75, 063512 [NASA ADS] [CrossRef] [Google Scholar]
  130. Sun, M., Voit, G. M., Donahue, M., et al. 2009, ApJ, 693, 1142 [NASA ADS] [CrossRef] [Google Scholar]
  131. Suto, Y., Sasaki, S., & Makino, N. 1998, ApJ, 509, 544 [NASA ADS] [CrossRef] [Google Scholar]
  132. Takahashi, R., Sato, M., Nishimichi, T., Taruya, A., & Oguri, M. 2012, ApJ, 761, 152 [NASA ADS] [CrossRef] [Google Scholar]
  133. Tanimura, H., Hinshaw, G., McCarthy, I. G., et al. 2019a, MNRAS, 483, 223 [NASA ADS] [CrossRef] [Google Scholar]
  134. Tanimura, H., Aghanim, N., Douspis, M., Beelen, A., & Bonjean, V. 2019b, A&A, 625, A67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  135. Tanimura, H., Hinshaw, G., McCarthy, I. G., et al. 2020, MNRAS, 491, 2318 [Google Scholar]
  136. Tinker, J. L., Weinberg, D. H., Zheng, Z., & Zehavi, I. 2005, ApJ, 631, 41 [NASA ADS] [CrossRef] [Google Scholar]
  137. Tinker, J., Kravtsov, A. V., Klypin, A., et al. 2008, ApJ, 688, 709 [NASA ADS] [CrossRef] [Google Scholar]
  138. Tinker, J. L., Robertson, B. E., Kravtsov, A. V., et al. 2010, ApJ, 724, 878 [NASA ADS] [CrossRef] [Google Scholar]
  139. Trac, H., Bode, P., & Ostriker, J. P. 2011, ApJ, 727, 94 [NASA ADS] [CrossRef] [Google Scholar]
  140. Valageas, P., & Nishimichi, T. 2011, A&A, 527, A87 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  141. van Daalen, M. P., Schaye, J., Booth, C. M., & Dalla Vecchia, C. 2011, MNRAS, 415, 3649 [NASA ADS] [CrossRef] [Google Scholar]
  142. van Daalen, M. P., Schaye, J., McCarthy, I. G., Booth, C. M., & Vecchia, C. D. 2014, MNRAS, 440, 2997 [CrossRef] [Google Scholar]
  143. van Daalen, M. P., McCarthy, I. G., & Schaye, J. 2020, MNRAS, 491, 2424 [CrossRef] [Google Scholar]
  144. van den Bosch, F. C., More, S., Cacciato, M., Mo, H., & Yang, X. 2013, MNRAS, 430, 725 [NASA ADS] [CrossRef] [Google Scholar]
  145. Van Waerbeke, L., Hinshaw, G., & Murray, N. 2014, Phys. Rev. D, 89, 023508 [NASA ADS] [CrossRef] [Google Scholar]
  146. Velliscig, M., van Daalen, M. P., Schaye, J., et al. 2014, MNRAS, 442, 2641 [NASA ADS] [CrossRef] [Google Scholar]
  147. Velliscig, M., Cacciato, M., Schaye, J., et al. 2015, MNRAS, 453, 721 [NASA ADS] [CrossRef] [Google Scholar]
  148. Vikhlinin, A., Burenin, R. A., Ebeling, H., et al. 2009, ApJ, 692, 1033 [NASA ADS] [CrossRef] [Google Scholar]
  149. Vogelsberger, M., Genel, S., Springel, V., et al. 2014, Nature, 509, 177 [NASA ADS] [CrossRef] [Google Scholar]
  150. White, M. 2004, Astropart. Phys., 22, 211 [NASA ADS] [CrossRef] [Google Scholar]
  151. Yan, Z., Hojjati, A., Tröster, T., Hinshaw, G., & van Waerbeke, L. 2019, ApJ, 884, 139 [CrossRef] [Google Scholar]
  152. Yan, Z., Raza, N., Van Waerbeke, L., et al. 2020, MNRAS, 493, 1120 [CrossRef] [Google Scholar]
  153. Yoo, J., Tinker, J. L., Weinberg, D. H., et al. 2006, ApJ, 652, 26 [NASA ADS] [CrossRef] [Google Scholar]
  154. Yoon, M., Jee, M. J., Tyson, J. A., et al. 2019, ApJ, 870, 111 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Evaluation of the two-halo term

There is some confusion in the literature about how exactly to numerically evaluate the integral for the halo model two-halo term:

P 2 H , u v ( k ) = P lin ( k ) i = u , v [ 0 b ( M ) W i ( M , k ) n ( M ) d M ] . $$ \begin{aligned} P_{2\mathrm{H},uv}(k)=P_{\mathrm{lin} }(k) \prod _{i=u,v}\left[\int _0^\infty b(M)W_i(M,k)n(M)\;\mathrm{d}M\right]. \end{aligned} $$(A.1)

We shall define the integral appearing in Eq. (A.1) to be

I u ( M 1 , M 2 , k ) = M 1 M 2 b ( M ) W u ( M , k ) n ( M ) d M . $$ \begin{aligned} I_u(M_1,M_2,k)=\int _{M_1}^{M_2} b(M)W_u(M,k)n(M)\;\mathrm{d}M. \end{aligned} $$(A.2)

In this appendix we specialise to the two-halo term for the matter–matter spectrum, u = v = m, as this provides a good illustration of the source of confusion. Equation (A.2) should be evaluated over all halo masses, M1 → 0 to M2 → ∞. Confusion arises because most numerical implementations evaluate this integral only over a finite range of halo mass, which seems sensible given that usually only a finite range of masses are thought to be relevant. However, note that for the matter spectrum we must have the following large-scale limit respected on physical grounds:

P 2 H , mm ( k 0 ) = P lin ( k 0 ) . $$ \begin{aligned} P_{2\mathrm{H},\mathrm{mm}}(k\rightarrow 0)=P_{\mathrm{lin} }(k\rightarrow 0). \end{aligned} $$(A.3)

So the integral in Eq. (A.2) must tend to unity when u = m in the k → 0 limit. In this case, W m ( M , k 0 ) M / ρ ¯ $ W_{\mathrm{m}}(M,k\to0)\to M/\bar\rho $ and Eq. (A.2) becomes

I m ( M 1 , M 2 , k 0 ) = 1 ρ ¯ M 1 M 2 b ( M ) M n ( M ) d M . $$ \begin{aligned} I_{\rm m}(M_1,M_2,k\rightarrow 0)=\frac{1}{\bar{\rho }}\int _{M_1}^{M_2} b(M)Mn(M)\;\mathrm{d}M. \end{aligned} $$(A.4)

If we pick sensible-seeming limits of M1 = 1010h−1M and M2 = 1016h−1M, then we find, for a standard cosmological model and a standard mass function at z = 0, that Im(M1 = 1010h−1M, M2 = 1016h−1M, k → 0)≃0.67. The problem is that, if taken literally, most popular mass functions have a large amount of the total matter contained in very low mass haloes, and so the convergence of Eq. (A.2) is slow as M1 is lowered. The upper limit of M2 is not usually a problem as long as it is reasonably high, say 1016h−1M at z = 0, which is effectively infinite here.

To address this issue, one thing that should absolutely not be done is to “renormalise” the mass function by multiplying it by a constant to enforce that Im(M1, ∞, k → 0) = 1 for a particular choice of M1. Physically, this amounts to removing low-mass haloes and putting all their mass into higher-mass haloes. In an example where Im(M1, ∞, k → 0) = 0.67 this would amount to increasing by ∼50% the numbers of haloes between the integration limits! Using a mass function “renormalised” in this way will drastically change the scale dependence in the two- and one-halo terms because it increases the number of haloes above M1.

There are two choices for how to proceed: the first (e.g. Yoo et al. 2006; Cacciato et al. 2012) is to assume that W m ( M < M 1 , k ) = M / ρ ¯ $ W_{\mathrm{m}}(M < M_1,k)=M/\bar\rho $, which is equivalent to taking low-mass haloes to have delta-function profiles. One can then numerically evaluate Eq. (A.2) above the limit of M1 and define a new function

A ( M 1 ) = 1 I m ( M 1 , , k 0 ) , $$ \begin{aligned} A(M_1)=1-I_{\rm m}(M_1,\infty ,k\rightarrow 0), \end{aligned} $$(A.5)

which is the missing part of the integral from below M1 (this only needs to be computed once). One can then add A(M1) every time Eq. (A.2) is evaluated:

I m ( M 1 , , k ) M 1 b ( M ) W m ( M , k ) n ( M ) d M + A ( M 1 ) . $$ \begin{aligned} \begin{split} I_{\rm m}(M_1,\infty ,k)\rightarrow&\int _{M_1}^{\infty } b(M)W_{\rm m}(M,k)n(M)\;\mathrm{d}M \\&+A(M_1). \end{split} \end{aligned} $$(A.6)

Note that this correction is additive, not multiplicative. This approximation is sufficient as long as the physical shapes of haloes at and below M1 do not influence the calculation; however, it is only the appropriate correction for Eq. (A.2) for matter. The second choice is to alter the mass function so that the missing signal from haloes below M1 is considered to be in haloes of mass exactly M1 (Schmidt 2016). This amounts to making the substitution

n ( M ) n ( M ) = n ( M ) + A ( M 1 ) δ D ( M M 1 ) b ( M 1 ) M 1 / ρ ¯ , $$ \begin{aligned} n(M)\rightarrow n^{\prime }(M)=n(M)+\frac{A(M_1)\delta _{\rm D}(M-M_1)}{b(M_1)M_1/\bar{\rho }}, \end{aligned} $$(A.7)

in Eq. (A.2). This reduces to a different additive correction to Eq. (A.2) of

I u ( M 1 , , k ) M 1 b ( M ) W u ( M , k ) n ( M ) d M + A ( M 1 ) W u ( M 1 , k ) M 1 / ρ ¯ · $$ \begin{aligned} \begin{split} I_u(M_1,\infty ,k)\rightarrow&\int _{M_1}^{\infty } b(M)W_u(M,k)n(M)\;\mathrm{d}M \\&+A(M_1)\frac{W_u(M_1,k)}{M_1/\bar{\rho }}\cdot \end{split} \end{aligned} $$(A.8)

This approximation is sufficient provided the physical shapes of haloes below M1 do not influence the calculation, but it maintains scale-dependence in the two-halo term compared to Eq. (A.6). We prefer the second option because it applies evaluations of Eq. (A.2) for any field u(r), not just for matter (note the “general field” subscript u in Eq. (A.8), compared to the “matter” subscript m in Eq. (A.6)). We apply this mass-function alteration only to the two-halo term in all of our calculations.

Note that the problem of slow convergence with respect to the lower integration limit does not usually affect the one-halo term for matter spectra because low-mass haloes contribute very little to the one-halo term, whose amplitude and shape depends on ⟨M⟩ (notation from Sect. 2.3) which is ∼1013.5h−1M at z = 0 for a standard cosmological model, and therefore is usually included in a sensible integration range for Eq. (A.2).

In cases where the power spectrum is unaffected by low-mass haloes then the discussion in this appendix is not relevant. For example, for the galaxy–galaxy number density power you only need to compute the integrals down to the halo mass below which the haloes contain no galaxies (obviously). In the case of the electron-pressure field, so little pressure comes from low mass haloes that this correction can be ignored. Considering the discussion in this appendix will, however, be relevant for any spectrum that involves the matter field because the low-mass haloes contribute significantly to the overall mass.

Appendix B: Power spectrum measurements from hydrodynamical simulations

In this paper we measure the power spectrum of overdensity fluctuations and electron pressure using particle data from hydrodynamic simulations. We do this using a Fast-Fourier Transform (FFT) algorithm applied to Cartesian fields, defined on a regular mesh, that are created from particle data. Generally, we consider a mesh-defined field u(r) where ui is the contribution to u(r) per particle located at ri

u ( r ) = i = 1 N u i δ D ( r r i ) , $$ \begin{aligned} u(\boldsymbol{r})=\sum _{i=1}^N u_i \delta _{\rm D}(\boldsymbol{r}-\boldsymbol{r}_i), \end{aligned} $$(B.1)

with N being the total number of particles contributing to u(r).

To create density-contrast fields from particles in hydrodynamical simulations it is necessary to assign a different density per particle to the mesh because hydrodynamical particles may have different masses mi. We therefore write

ρ ( r ) ρ ¯ = i = 1 N m i m ¯ δ D ( r r i ) , $$ \begin{aligned} \frac{\rho (\boldsymbol{r})}{\bar{\rho }}={\sum _{i=1}^N \frac{m_i}{\bar{m}}} \delta _{\rm D}(\boldsymbol{r}-\boldsymbol{r}_i), \end{aligned} $$(B.2)

where m ¯ $ \bar{m} $ is the expected mean mass in a cell. In our definitions of overdensity, this division to create a density contrast is always by the mean mass of all matter, not only of the specific species being considered. The situation is more complicated when creating the electron pressure fields from the particles. Hydrodynamic gas particles (labelled with i) come tagged with an internal energy (kBTi) and a local density (ρi). To convert the particle gas temperature to a contribution to the total electron pressure per particle, Pe, i, a quantity that can be summed to create an electron-pressure field on mesh via

P e ( r ) = i = 1 N P e , i δ D ( r r i ) , $$ \begin{aligned} P_{\rm e}(\boldsymbol{r})=\sum _{i=1}^N P_{\mathrm{e},i} \delta _{\rm D}(\boldsymbol{r}-\boldsymbol{r}_i), \end{aligned} $$(B.3)

we use the ideal gas law:

P e , i V = N e , i k B T i $$ \begin{aligned} P_{\mathrm{e},i}V=N_{\mathrm{e},i}k_{\rm B}T_i \end{aligned} $$(B.4)

where Ne, i is the total number of free electrons in the gas particle in the mesh-cell volume V:

N e , i = m i m p 1 μ e , $$ \begin{aligned} N_{\mathrm{e},i}=\frac{m_i}{m_{\rm p}}\frac{1}{\mu _{\rm e}}, \end{aligned} $$(B.5)

where mi is the gas particle mass, mp is the proton mass and μe is the number of free electrons per proton mass. In creating the electron-pressure field from the particles we also calculate the particle hydrogen number density

n H , i = f H ρ i m p , $$ \begin{aligned} n_{\mathrm{H},i}=\frac{f_{\rm H} \rho _i}{m_{\rm p}}, \end{aligned} $$(B.6)

where fH ≃ 0.752 is the hydrogen gas mass fraction. We exclude particles with nH >  0.1 cm−3 as these particles are presumably cold, neutral and should eventually form stars. We checked that imposing this exclusion does not greatly affect our results.

Pe, i from Eq. (B.4) is directly related to the quantity Υi used in Roncarelli et al. (2006, 2007) and McCarthy et al. (2018):

Υ i = σ T k B T i m e c 2 m i μ e m p , $$ \begin{aligned} \Upsilon _i=\sigma _{\rm T}\frac{k_{\rm B}T_i}{m_{\rm e}c^2}\frac{m_i}{\mu _{\rm e} m_{\rm p}}, \end{aligned} $$(B.7)

where me is the electron mass and σT is the Thompson scattering cross section. Υi differs from Pe, i only by a factor of volume and some constants and has dimensions of area rather than pressure.

In N-body simulations it is well known that the raw measured power spectrum receive an unphysical additive contribution from the finite number of particles – so-called discreteness or “shot” noise. This contribution arises due to the automatic correlation of particles with themselves at zero separation. In the correlation function this contribution is confined to zero separation, but in P(k) the contribution is spread evenly over all wavenumbers by the Fourier transform (the Fourier transform of a delta function is a constant). In a standard N-body simulation with equal-mass particles the shot noise contribution to the matter–matter P(k) spectrum is:

S = L 3 N , $$ \begin{aligned} S=\frac{L^3}{N}, \end{aligned} $$(B.8)

where N is the total number of particles in the simulation and L is the box size. As expected, the shot noise contribution is larger when there are fewer particles. It is common practice to subtract this contribution from power spectra measured from N-body simulations since in the limit of an infinite number of particles it would vanish.

From hydrodynamical simulations we are interested in calculating auto- and cross-spectra between between different fields generated from the particles; for example, density, pressure or temperature fields. In the case that the fields to be correlated are generated from distinct sets of particles then the shot noise contribution is zero (e.g. gas–star density). However, in other cases the fields might be generated using (a subset of) same particles (e.g. matter–star density) or from different properties of the same particles (e.g. gas density–electron pressure, which both come from the same hydrodynamical gas particles). There is also the confounding issue that each particle does not provide the same contribution to the eventual field (e.g. individual stellar and gas particles have different masses, gas particles contribute different pressures).

The additive shot-noise contribution for the cross spectrum, Puv(k), between two fields u(r) and v(r), comprised of particles, where a subset N $ \tilde{N} $ of those particles contribute to both fields u(r) and v(r), can be shown to be

S uv = L 3 M 2 i = 1 N u i v i , $$ \begin{aligned} S_{uv}=\frac{L^3}{M^2}\sum _{i=1}^{\tilde{N}} u_i v_i, \end{aligned} $$(B.9)

where L is the simulation box size and M is the total number of mesh cells. Note that the shot-noise contribution is not really a function of M, as this will always cancel with factors of M that come from ui and vi. Equation (B.9) can be derived by taking the Fourier transform of Eq. (B.1) and then considering the part of the power spectrum that would still arise even if the particle positions ri are completely uncorrelated, so that the only contribution is from the self-correlation of each particle with itself.

Let us demonstrate that Eq. (B.9) reduces to more standard formulae: if we are interested in shot noise in the matter power spectrum we start from the particle masses from which we create a density field. The density-contrast field generated from the particles is given by Eq. (B.2). We can write m ¯ $ \bar{m} $ as the total mass in the simulation divided by the total number of mesh cells:

m ¯ = 1 M i = 1 N m i . $$ \begin{aligned} \bar{m}=\frac{1}{M}\sum _{i=1}^{N} m_i. \end{aligned} $$(B.10)

Using Eq. (B.9) and considering shot-noise in the matter–matter power spectrum m we get:

S mm = L 3 i = 1 N m i 2 [ i = 1 N m i ] 2 · $$ \begin{aligned} S_{\rm mm}=L^3\frac{\sum _{i=1}^N m_i^2}{\left[\sum _{i=1}^N m_i \right]^2}\cdot \end{aligned} $$(B.11)

If all the particle have equal masses then Eq. (B.11) reduces to Eq. (B.8) as required.

In Fig. B.1 we show the contribution of shot noise to the raw auto-spectra measured for the different fields from the AGN BAHAMAS simulations. Subtracting shot noise is most important for the gas spectrum, which makes sense because gas is the most diffuse component and therefore the one for which discreteness effects are most obvious. At z = 2 at k = 7 h Mpc−1 shot noise can makes up 8% of the raw measured gas–gas spectrum. Shot noise is also important at a similar level in some of the cross-spectra that we measure. For all of the results presented in this paper shot noise has been subtracted from the power spectra.

thumbnail Fig. B.1.

Shot-noise contributions to auto-spectra of the various fields measured from the AGN BAHAMAS simulation. We show the ratio of the shot-noise to the auto-spectra as a function of scale for matter (purple), CDM (green), gas (blue), stars (orange) and electron pressure (yellow). The dashed line shows the cutoff for a one percent contribution. We see that shot noise is most important for the gas power spectrum at higher redshifts, where it can contribute as much as 8% of the raw measured power for the smallest scales shown. Generally shot noise becomes less important as the simulation evolves and structure develops.

Appendix C: Variations in halo-model parameters

In this appendix we present figures that demonstrate the effect of varying our fitted parameters on the halo profiles and on the eventual power spectra that we are interested in. We show the effect of varying A* (Fig. C.6), M* (Fig. C.7), η (Fig. C.9), ϵ1 (Fig. C.1), Γ (Fig. C.5), M0 (Fig. C.3), α (Fig. C.10) and Twhim (Fig. C.11) on the power spectra, response, halo density profiles and halo mass fractions. These plots are extremely useful if one wants to gain physical insight on the effect of changing the hydrodynamical parameters on the eventual Fourier Space statistics, which can otherwise be quite difficult to understand.

We see that increasing A* increases the amplitude of the stellar content of haloes independently of their mass, and this has an almost scale-independent effect on any power spectra that includes the star field. There is a small back reaction on the gas because stars can only be added at the expense of gas, but this is very small because stars are a small fraction of the available halo baryons, even in the most extreme cases. Changing M* changes the halo mass at which star formation peaks, which has a more scale-dependent effect on the eventual power spectra as stellar mass is shifted into haloes with different structural properties. Altering η changes the distribution of stellar mass between the extended satellite galaxy part and the central galaxy. This has no effect at the largest scales, but effects the one-halo term. Shifting more stellar mass into central galaxies makes this term more and more shot-noise like.

ϵ1 governs the concentration modification of lower-mass galaxies that have ejected most of their gas. This has no effect on the power spectra at large scales, but changes them at small scales where there is sensitivity to the halo concentration. Γ governs the gas density and pressure profiles and has large effects on the halo profiles. Decreasing Γ makes the gas and pressure profiles more centrally concentrated and therefore boosts the associated power spectra at small scales. Altering M0 changes the mass corresponding to haloes that have lost half of their gas (with lower-mass haloes having lost more than half). This has a comparatively large effect on spectra, with a strong scale dependent effect on the matter–gas (no change at large scales) and a more scale-independent effect on spectra involving the pressure.

Increasing α increases the temperature of gas that is bound to the host halo. This boosts the pressure signal in a perfectly scale-independent way. There is no effect whatsoever on any matter spectra since these are completely insensitive to gas temperature. Conversely, increasing Tw boosts the temperature of the ejected gas, which is correlated only on large scales and thus only affects the large-scale, two-halo portions of power spectra that involve the pressure.

thumbnail Fig. C.1.

Effect of varying ϵ1, which governs deviations in halo concentration for low-mass haloes that have ejected their gas, on halo profiles (top right two panels), halo mass fractions (bottom right) and power spectra and responses (left column; top left two panels for matter spectra, bottom left two panels for pressure). Increasing ϵ1 increases halo concentration and this generally boosts power. From the top right we see this has the largest impact for lower mass haloes because these have ejected most of their gas. From the left-hand panels, this has comparatively more of an impact on power spectra involving gas than on the other spectra. We can also see similar effects on spectra that involve the electron pressure.

thumbnail Fig. C.2.

Same as Fig. C.1 but for ϵ2, which governs deviations in halo concentration for high-mass haloes that have retained their gas. Increasing ϵ2 increases halo concentration and this generally boosts power. From the top right we see this has the largest impact for higher mass haloes because these have retained most of their gas. From the left-hand panels, this has comparatively more of an impact on power spectra involving gas than on the other spectra. We can also see some effects on spectra that involve the electron pressure.

thumbnail Fig. C.3.

Same as Fig. C.1 but for M0, which governs the halo mass below which haloes have lost most of their gas. As can be seen in the bottom-right panel, changing M0 changes the split of gas between the bound component, which dominates high-mass haloes and the unbound component, which is the case for low-mass haloes that have jettisoned most of their gas. The effects of this on the bound halo profiles can be seen in the top-right panel. Decreasing M0 adds to the amplitude of the gas density profiles and this boosts power in power spectra involving gas. We see a similar, but more limited effect on spectra involving the electron pressure. This limitation arises because low-mass haloes are lower temperature, so contribute much less to spectra involving the electron pressure than they do to those involving the gas density. There is also an overall amplitude affect on the pressure spectra because the unbound gas is taken to have a fixed Tw = 106 K that contributes power on large scales.

thumbnail Fig. C.4.

Same as Fig. C.1 but for β, the power-law index governing the decline in halo gas fraction for low-mass haloes. As can be seen from the bottom-right panel, increasing β makes the decline in gas fraction sharper, and this has a comparatively large effect on low-mass haloes as can be seen in the other panels on the right-hand side of the figure. High mass haloes are affected in the opposite sense because the gas fraction is forced to be 0.5 at 1014h−1M. The effects on the power spectra can be seen in the left-hand panels, where the matter–gas spectrum can be seen to be most affected.

thumbnail Fig. C.5.

Same as Fig. C.1 but for Γ, the polytropic index for the gas. Increasing Γ makes the gas profile less concentrated, as can be seen in the halo-profile plots (top right). If Γ is lowered then gas is more concentrated and the power spectrum of gas is boosted at small scales. Effects on power spectra that do not involve gas are much smaller. The same trends are seen in the electron-pressure profile and spectra involving the electron pressure, which are also sensitive to this parameter as the gas pressure is partly determined by the gas density.

thumbnail Fig. C.6.

Same as Fig. C.1, but for A*, which governs the peak-star-formation efficiency in haloes. Lower-right panel: A* effects the overall abundance of stellar mass, but not the ratio in central and satellite stellar mass. There is some back reaction on the gas abundance that occurs in order to keep the baryon fraction universal, so that halo gas is slightly depleted as a result of this boost, so there is a small effect on gas profiles. The effect on the halo profiles can be seen in the top right panel. Increasing A* boosts star-formation efficiency, which has a fairly uniform effect on the amplitude of any power spectra involving the star density as can be seen in the top two panels in the left-hand column. There is an effect on the matter–matter spectrum at small scales as stars start to dominate the power. The effect on electron pressure is minimal, but not non-existent due to the gas abundance being affected.

thumbnail Fig. C.7.

Same as Fig. C.1 but for M*, the halo mass at which star formation efficiency peaks. Changing this parameter shifts the halo-stellar-mass fraction to lower and higher masses, as can be seen in the bottom-right panel, there is a small back-reaction effect on the gas abundance. Changing M* indirectly affects the halo profiles at different masses, as can be seen in the top right panel, but note that the delta-function component of the stellar halo profile is not shown in this plot. The power spectra involving stars are affects in a more scale-dependent manner as can be seen in the top-left two panels.

thumbnail Fig. C.8.

Same as Fig. C.1 but for σ*, which governs the width of the stellar-mass distribution around the peak (which is at M* = 1012.5h−1M here). As can be seen in the bottom-left panel, decreasing σ* shrinks the stellar-mass distribution, resulting in fewer stars in both higher and lower mass haloes than M*. The top-right panel shows that this has the biggest impact on 1014h−1M haloes, not higher-mass haloes due to the saturation we impose at high masses. In the left-hand panels we see that the effects of this are constrained to be on power spectra involving the stars. There is only a tiny effect on gas spectra which is due to baryonic mass conservation.

thumbnail Fig. C.9.

Same as Fig. C.1 but for η, which determines the split of stellar halo mass between central and satellite galaxies for M >  M*, where M* = 1012.5h−1M in this figure. In the lower-right panel we see that increasing η from the fiducial −0.3 means that more mass is put into central galaxies, whereas decreasing it places more mass in satellites. As a result of this, the amplitude of the stellar halo is changed, as can be seen in the top-right panel, with decreasing η putting more mass into the NFW part of the profile and thus boosting the halo amplitude. This has a scale-dependent effect on the matter–stars power spectrum, as can be seen in the top-left two panels.

thumbnail Fig. C.10.

Same as Fig. C.1 but for α. α governs deviations of the halo temperature from the virial relation. We see that this affects the pressure profiles (middle right) by uniformly scaling the amplitude, because the pressure is directly proportional to temperature; increasing the temperature increases the pressure. We can see the effect of this on the pressure power spectra in the lower left two panels where the effect is not uniform due to the complicated dependence of the power spectra on the underlying halo populations in the two- and one-halo terms, which have different scalings with temperature.

thumbnail Fig. C.11.

Same as Fig. C.1 but for Tw. Varying Tw only affects the temperature of the unbound gas that contributes to the two-halo term only. We see that as we increase the gas temperature we boost the two-halo term. Note that there is a floor to the effect that lowering Tw can have, which is because the one-halo term contributes at large scales in the pressure spectra too, and lowering the temperature of this gas beyond a certain point means that only the two-halo term becomes dwarfed by the one-halo term at large scales.

Appendix D: Projected fields

In this paper we are primarily interested in the cross correlation between tSZ and various lensing quantities. We use the metric convention of Bartelmann & Schneider (2001):

d s 2 = d t 2 a 2 [ d r 2 + f k 2 ( r ) ( d θ 2 + sin 2 θ d ϕ 2 ) ] , $$ \begin{aligned} \mathrm{d}s^2=\mathrm{d}t^2-a^2[\mathrm{d}r^2+f_k^2(r)(\mathrm{d}\theta ^2+\sin ^2\theta \mathrm{d}\phi ^2)], \end{aligned} $$(D.1)

where a is the scale factor, r is the comoving distance, and fk is the comoving, angular-diameter distance.

The weak-lensing signal is produced by the bending of light by the large-scale structure of the Universe. With some approximations, gravitational lensing is conveniently summarised by the convergence field κ, which can be written as

κ ( θ ) = 3 2 Ω m ( H 0 c ) 2 0 r H f k ( r ) a ( r ) q ( r ) δ m ( r , θ ) d r . $$ \begin{aligned} \kappa (\boldsymbol{\theta })=\frac{3}{2}\Omega _{\rm m}\left(\frac{H_0}{c}\right)^2\int _0^{r_{\rm H}}\frac{f_{k}(r)}{a(r)}q(r)\delta _{\rm m}(r,\boldsymbol{\theta })\;\mathrm{d}r. \end{aligned} $$(D.2)

Here q(r) is the lensing-efficiency kernel, which weights redshifts along the line-of-sight according to their contribution to the total lensing

q ( r ) = z ( r ) p ( z ) f k [ r ( z ) r ] f k [ r ( z ) ] d z , $$ \begin{aligned} q(r)=\int _{z(r)}^\infty p(z^{\prime })\frac{f_k[r^{\prime }(z^{\prime })-r]}{f_k[r^{\prime }(z^{\prime })]}\;\mathrm{d}z^{\prime }, \end{aligned} $$(D.3)

where p(z) is the normalised redshift distribution of lensed galaxies. If one is interested in CMB lensing then p(z) = δD(z − z*) where z* is the redshift of the last-scattering surface. In this case:

q ( r , r ) = f k ( r r ) f k ( r ) , $$ \begin{aligned} q(r,r_*)=\frac{f_k(r_*-r)}{f_k(r_*)}, \end{aligned} $$(D.4)

where r* is the comoving distance to z*.

The tSZ signal in the CMB is produced when high energy, free electrons increase the energy of passing, cool CMB photons via inverse Compton scattering. This results in a boost in the effective temperature of the CMB at the location of the electrons that is in proportion to the electron “pressure”, usually expressed in terms of the dimensionless Compton “y” parameter:

y ( θ ) = σ T m e c 2 0 r H P e ( r , θ ) a 3 ( r ) a ( r ) d r , $$ \begin{aligned} y(\boldsymbol{\theta })=\frac{\sigma _{\rm T}}{m_{\rm e}c^2}\int _0^{r_{\rm H}} \frac{P_{\rm e}(r,\boldsymbol{\theta })}{a^3(r)}\;{a(r)\mathrm{d}r}, \end{aligned} $$(D.5)

where σT is the Thompson scattering cross section, me is the electron mass, c is the speed of light. The integral is taken over comoving distance between us and the particle horizon, rH. The electron pressure can be written as Pe = nekBTe where ne is the comoving electron number density and Te is the physical electron temperature. The factors of a convert our comoving r and electron pressure to a physical quantities.

Given a general 2D field U(θ) that is the projected version of a 3D field u(r, θ) via some projection kernel XU(r):

U ( θ ) = 0 r H X U ( r ) u ( r , θ ) d r $$ \begin{aligned} U(\boldsymbol{\theta })=\int _0^{r_{\rm H}} X_U(r) u(r,\boldsymbol{\theta })\;\mathrm{d}r \end{aligned} $$(D.6)

we can write the angular power spectrum between two such projected fields, U and V, in terms of the 3D power spectrum of u and v using the Limber approximation (Kaiser 1992) as

C UV ( ) = 0 r H X U ( r ) X V ( r ) f k 2 ( r ) P uv ( k ( r ) , z ( r ) ) d r , $$ \begin{aligned} C_{UV}(\ell )=\int _0^{r_{\rm H}}\frac{X_U(r)X_V(r)}{f^2_k(r)} P_{uv}\left(k(r),z(r)\right)\;\mathrm{d}r, \end{aligned} $$(D.7)

where k(r) = ( + 1/2)/fk(r) and z(r) corresponds to the redshift at comoving distance r. We have increased the low- accuracy of the approximation to 𝒪(−2) by including the lowest-order correction  →  + 1/2 (Lo Verde et al. 2008). For the lensing projection kernel we have

X κ ( r ) = 3 2 Ω m ( H 0 c ) 2 f k ( r ) q ( r ) a ( r ) , $$ \begin{aligned} X_\kappa (r)=\frac{3}{2}\Omega _{\rm m}\left(\frac{H_0}{c}\right)^2\frac{f_k(r)q(r)}{a(r)}, \end{aligned} $$(D.8)

and we must project the 3D matter–matter power spectrum, Pmm(k). For y projection we have

X y ( r ) = σ T m e c 2 1 a 2 ( r ) , $$ \begin{aligned} X_y(r)=\frac{\sigma _{\rm T}}{m_{\rm e}c^2}\frac{1}{a^2(r)}, \end{aligned} $$(D.9)

and we must project the 3D matter–electron pressure power spectrum, PmP(k), to obtain the lensing-tSZ cross correlation.

There has been some recent discussion on the validity of the Limber approximation for the interpretation of cosmological observations. We note that comparisons between this approximation and a full calculation show that the Limber approximation is valid given the scales we are interested in, and that we will be interested in in the near future. We point the interested reader towards Hill & Pajer (2013), Kilbinger et al. (2017) and Lemos et al. (2017). Note that the accuracy of the Limber approximation is generally better for Cκy() compared to for Cκκ() due to the broader projection kernel for y compared to that of κ.

The contributions to Cκκ() and Cκy() as a function of k, z and r are shown in Fig. D.1 for different . Note that, for a given , Cκy() receives its contributions from lower z and r compared to Cκκ(), and also from a much broader range of k, z and r. To generate this figure we used a p(z) distribution (Eq. (D.3)) taken from KiDS galaxies for redshifts z = 0.1 → 0.9 (Joudaki et al. 2017b).

thumbnail Fig. D.1.

Normalised contribution to each of C() for κκ (top) and κy (bottom) as a function of scale k, redshift z and (comoving) distance r. The lensing is taken from KiDS galaxies for redshifts z = 0.1 → 0.9. Each of these curves is the normalised integrand in Eq. (D.7) re-expressed as a function of either k, z or r (so dC()/dk, dC()/dz, dC()/dr). In practice, we carry out our integration in r between the observer and z(r) = 9 and we checked that our integrations are robust to this upper limit. The power spectra being integrated are set to zero for k <  10−4h Mpc−1 and k >  102h Mpc−1 and again our results are robust to this. Note the generally broader and lower–z contributions for κy compared to those for κκ.

All Tables

Table 1.

Default hydrodynamical halo model parameters.

Table 2.

Best-fitting effective halo model parameters for different combinations of power spectra from the BAHAMAS simulations.

All Figures

thumbnail Fig. 1.

Halo mass fractions fu(M) as a function of halo mass for the different components in our model at z = 0: The halo CDM fraction is constant at Ωcm for all halo masses. Star formation efficiency peaks around M* = 1012.5h−1M with a mass fraction A* = 0.03, it saturates at high mass at a value of f* = A*/3. Stars are split into central and satellite, with the halo stellar content being dominated by centrals at low masses and by satellites at high mass. The bound gas fraction is near to the universal baryon fraction for high halo masses, but decreases below half the universal fraction for halo masses lower than M0 = 1014h−1M, the remaining baryons that are neither bound gas or stars are considered to be in unbound gas that is not in the halo. The total mass fraction drops from unity at lower masses due to the ejected gas.

In the text
thumbnail Fig. 2.

Halo-model predictions at z = 0 for the matter–matter (top), CDM–CDM, gas–gas, stars–stars and the matter–electron pressure (bottom) power spectra. We also show the two-halo (long-dashed) and the one-halo (short-dashed) terms for each spectrum apart from matter–matter. At large scales, the shapes of the two-halo terms are all identical, and are that of the linear power spectrum, at small scales the one-halo terms dominates and have significantly different shapes for each spectrum. Note that the transition scale between the two- and one-halo term is at much larger scales for the matter–electron pressure cross spectrum, which is a consequence of this spectrum being dominated by contributions from more massive haloes compared to those that contribute to the matter spectra. The down-turn at small scales in the gas–gas and matter–electron pressure spectra is partly due to the underrepresentation of low-mass haloes and partly due to their relatively smooth halo profiles.

In the text
thumbnail Fig. 3.

Cumulative halo model power spectra as a function of the upper-limit of halo mass at z = 0. Left-hand panel: matter–matter power spectrum as the upper limit is raised from 1010 to 1016h−1M (i.e. the lightest coloured curve shows the contribution from only haloes below 1010h−1M), right-hand panel: same but for the matter–electron pressure spectrum. In both cases the converged spectrum from all halo masses is shown in thick black and this has been reached with an upper limit of 1016h−1M. We see that the matter–matter spectrum builds up fairly gradually with halo mass and that all scales receive contributions from a wide range of halo masses. In contrast, the matter–electron pressure spectrum receives very little contribution from haloes less massive than 1013h−1M, with the vast majority of the power coming from haloes between 1014 and 1015h−1M.

In the text
thumbnail Fig. 4.

Effect of σ8 variations on the matter–matter power spectrum (left) and on the matter–electron pressure power spectrum (right) at z = 0. We show halo model power spectra for σ8 = 0.7 (blue) to σ8 = 0.9 (red) divided by a central σ8 = 0.8 model. For linear theory the plot would simply be horizontal lines. We see that boosting the amplitude of fluctuations increases the amplitude of the matter–electron power spectrum far more than it does for the matter-matter. At a wavenumber of ≃1 h Mpc−1 the matter-matter power spectrum scales like σ 8 3.5 $ \sigma_8^{3.5} $ while matter–electron pressure scales like σ 8 5.8 $ \sigma_8^{5.8} $. This is because the electron pressure field is dominated by high-mass haloes from the tail of the mass function that are very sensitive to the linear power spectrum amplitude.

In the text
thumbnail Fig. 5.

Different fields measured from the BAHAMAS AGN simulation at z = 0 around a M ∼ 1015h−1M halo averaged through a square slab of side 20 h−1 Mpc and thickness 5 h−1 Mpc. We show CDM (top left), gas (top right) and star (bottom left) density contrast as well as electron pressure (bottom right). Colour intensity increases with field value logarithmically and we show a dynamic range of 103 for both density contrast and electron pressure. We see that the gas is less clumped than the CDM and that low-mass sub haloes show a large gas depletion compared to the host. Stars are more tightly clustered in the halo centre than the CDM, but are absent in many of the substructures. The total matter density-contrast field would be a sum of the CDM, gas and star fields. The electron pressure broadly follows the gas but can only be seen emanating from the highest gas-density peaks over the dynamic range shown.

In the text
thumbnail Fig. 6.

Comparison of the default halo-model prediction with power spectra from the BAHAMAS AGN simulation at z = 0. Points with errors show the measured power spectra with an error-on-the-mean calculated from the finite number of modes that contribute to each k. The upper two panels show spectra for “all matter” (purple), CDM (green), gas (light blue), stars (orange) and electron pressure (dark blue). The electron pressure field has units of 100 eV cm−3. The effect of power aliasing can be seen as an upturn in power at the highest k (∼7 h Mpc−1) shown for each spectrum (see Appendix B). Lines show the halo-model predictions using the default model discussed in Sect. 3: solid lines and points show auto-spectra while broken lines and open points show cross-spectra with the total matter field. We see that the trends in the halo model agree reasonably well with the simulations but that there are disagreements in the details. In the lower two panels we show the “response function” calculated with respect to the matter–matter power in the gravity-only model: P uv hydro ( k ) / P mm gravity ( k ) $ P^{\mathrm{hydro}}_{uv}(k) / P^{\mathrm{gravity}}_{\mathrm{mm}}(k) $. The simulations spectra have been divided by those from the DMONLY simulation and the halo model has been divided by a standard halo-model prediction that assumes all matter to be in NFW haloes. The horizontal-dashed-black lines in the lower left panel show the expected large scale Ωcm, (Ωcm)2, Ωbm, and (Ωbm)2 asymptotes for the CDM and gas cross–matter and auto-spectra. We see that the response functions are generally smoother than the raw power spectra: a consequence of cosmic variance cancellations at large scales and cancellation of aliasing effects at small scales.

In the text
thumbnail Fig. 7.

Halo model (1) for the star–star power spectrum is shown as sold lines while simulation measurements are shown as points with error bars that come from the scatter in power between three different realisations of the AGN BAHAMAS simulation. Note that the power spectra shown here is the measured response multiplied by HMCODE, so cosmic variance is removed at small scales, but we still show the error bar for comparison. Different colours denote the three different AGN feedback temperatures: 107.6 (blue), 107.8 (grey) and 108.0 K (red). Parameters for the model can be found in Table 2. We see that the spectrum is reduced in amplitude as feedback temperature increases, a consequence of AGN feedback suppressing star formation. This is reflected in the modelling as the value of A* and M* decreases as AGN strength increases.

In the text
thumbnail Fig. 8.

Similar to Fig. 7, but for the matter–matter spectrum (2). The power spectra are shown in the top row, while the response functions are shown in the bottom row. We see that the model is nearly perfect (errors at the sub-percent level) across a wide range of scales and for the different feedback scenarios. We see that the amount of power-spectrum suppression increases with the feedback temperature. The critical parameter that governs this in our modelling is the halo mass M0, which governs the halo mass below which more than half of the cosmic baryon density is missing from the halo. In the fitting parameters we see a strong preference for higher values of M0 for stronger feedback temperatures (from blue to red). How our model fares outside the range of AGN heating temperatures to which it was fitted is also shown via the lines that have no corresponding simulation points. To generate these we have interpolated and extrapolated our model parameters as a function of AGN heating temperature, which is shown in 0.1 dex intervals between 107.5 (bluest) and 108.1 K (reddest), indicating that our model is robust to extrapolation.

In the text
thumbnail Fig. 9.

Same as Figs. 7 and 8, but for model (3) for the matter–matter (first and third row; power and response) and matter–electron pressure spectra (second and fourth row; power and response). The match to the matter–matter spectra is slightly degraded compared to the case when we fitted the matter–matter spectra exclusively, but we still match the response at the few percent level while simultaneously we match the matter–electron pressure spectrum well for k >  1 h Mpc−1. If we compare this to the error bars in the power spectra in the second row we see that our fit is good for scales where the three AGN heating temperature strengths are clearly demarcated in k. At larger scales we suspect the larger errors arise because the power is predominantly coming from few massive haloes. Similar to Fig. 8, we show how the model fares outside the range of AGN strengths to which is was fitted. AGN heating temperature is shown in 0.1 dex intervals between 107.5 K and 108.1 K, while the model was fitted between 107.6 K and 108.0 K, thus demonstrating that our model is robust to extrapolation.

In the text
thumbnail Fig. 10.

Same as Figs. 79, but for the best-fitting model (4) to all 10 possible power spectra of matter, CDM, gas and star fields, although we only show the 4 auto-spectra here. The top half of the figure show power spectra while the bottom half show response functions. Columns show different redshifts. Different colours are different feedback temperatures: T = 107.6 (blue), T = 107.8 (grey) and T = 108.0 K (red). Points (with errors) are from BAHAMAS while lines are from our halo model.

In the text
thumbnail Fig. 11.

Full error matrix for our fitted models calculated for each field pair for the 107.8 K AGN simulation. The error is averaged linearly over redshift between z = 0 and 1 and over k logarithmically between ≃0.015 and 7 h Mpc−1. This matrix is symmetric because power spectra are symmetric when swapping the field labels. We show the matrix for the three models presented in this work: matter (2, left); matter, electron pressure (3, centre); matter, CDM, gas, stars (4, right). The squares that contain dots are the particular cross-spectra that are fitted in the model. In the left-hand plot, we see that in the matter model we get a reasonable fit to the power spectra all of the constituent fields. In the central plot we get a better match to to the matter–electron pressure spectrum than in the left-hand panel, but this is at the expense of any power spectra that directly involves the gas field, which shows some tension in our model between the gas and pressure modelling. In the right-hand panel, where we fitted matter, CDM, gas and stars, we see that all the matter spectra are reasonable, and better than in the left-hand panel, but achieving this is at the expense of any spectra that involves the pressure field.

In the text
thumbnail Fig. B.1.

Shot-noise contributions to auto-spectra of the various fields measured from the AGN BAHAMAS simulation. We show the ratio of the shot-noise to the auto-spectra as a function of scale for matter (purple), CDM (green), gas (blue), stars (orange) and electron pressure (yellow). The dashed line shows the cutoff for a one percent contribution. We see that shot noise is most important for the gas power spectrum at higher redshifts, where it can contribute as much as 8% of the raw measured power for the smallest scales shown. Generally shot noise becomes less important as the simulation evolves and structure develops.

In the text
thumbnail Fig. C.1.

Effect of varying ϵ1, which governs deviations in halo concentration for low-mass haloes that have ejected their gas, on halo profiles (top right two panels), halo mass fractions (bottom right) and power spectra and responses (left column; top left two panels for matter spectra, bottom left two panels for pressure). Increasing ϵ1 increases halo concentration and this generally boosts power. From the top right we see this has the largest impact for lower mass haloes because these have ejected most of their gas. From the left-hand panels, this has comparatively more of an impact on power spectra involving gas than on the other spectra. We can also see similar effects on spectra that involve the electron pressure.

In the text
thumbnail Fig. C.2.

Same as Fig. C.1 but for ϵ2, which governs deviations in halo concentration for high-mass haloes that have retained their gas. Increasing ϵ2 increases halo concentration and this generally boosts power. From the top right we see this has the largest impact for higher mass haloes because these have retained most of their gas. From the left-hand panels, this has comparatively more of an impact on power spectra involving gas than on the other spectra. We can also see some effects on spectra that involve the electron pressure.

In the text
thumbnail Fig. C.3.

Same as Fig. C.1 but for M0, which governs the halo mass below which haloes have lost most of their gas. As can be seen in the bottom-right panel, changing M0 changes the split of gas between the bound component, which dominates high-mass haloes and the unbound component, which is the case for low-mass haloes that have jettisoned most of their gas. The effects of this on the bound halo profiles can be seen in the top-right panel. Decreasing M0 adds to the amplitude of the gas density profiles and this boosts power in power spectra involving gas. We see a similar, but more limited effect on spectra involving the electron pressure. This limitation arises because low-mass haloes are lower temperature, so contribute much less to spectra involving the electron pressure than they do to those involving the gas density. There is also an overall amplitude affect on the pressure spectra because the unbound gas is taken to have a fixed Tw = 106 K that contributes power on large scales.

In the text
thumbnail Fig. C.4.

Same as Fig. C.1 but for β, the power-law index governing the decline in halo gas fraction for low-mass haloes. As can be seen from the bottom-right panel, increasing β makes the decline in gas fraction sharper, and this has a comparatively large effect on low-mass haloes as can be seen in the other panels on the right-hand side of the figure. High mass haloes are affected in the opposite sense because the gas fraction is forced to be 0.5 at 1014h−1M. The effects on the power spectra can be seen in the left-hand panels, where the matter–gas spectrum can be seen to be most affected.

In the text
thumbnail Fig. C.5.

Same as Fig. C.1 but for Γ, the polytropic index for the gas. Increasing Γ makes the gas profile less concentrated, as can be seen in the halo-profile plots (top right). If Γ is lowered then gas is more concentrated and the power spectrum of gas is boosted at small scales. Effects on power spectra that do not involve gas are much smaller. The same trends are seen in the electron-pressure profile and spectra involving the electron pressure, which are also sensitive to this parameter as the gas pressure is partly determined by the gas density.

In the text
thumbnail Fig. C.6.

Same as Fig. C.1, but for A*, which governs the peak-star-formation efficiency in haloes. Lower-right panel: A* effects the overall abundance of stellar mass, but not the ratio in central and satellite stellar mass. There is some back reaction on the gas abundance that occurs in order to keep the baryon fraction universal, so that halo gas is slightly depleted as a result of this boost, so there is a small effect on gas profiles. The effect on the halo profiles can be seen in the top right panel. Increasing A* boosts star-formation efficiency, which has a fairly uniform effect on the amplitude of any power spectra involving the star density as can be seen in the top two panels in the left-hand column. There is an effect on the matter–matter spectrum at small scales as stars start to dominate the power. The effect on electron pressure is minimal, but not non-existent due to the gas abundance being affected.

In the text
thumbnail Fig. C.7.

Same as Fig. C.1 but for M*, the halo mass at which star formation efficiency peaks. Changing this parameter shifts the halo-stellar-mass fraction to lower and higher masses, as can be seen in the bottom-right panel, there is a small back-reaction effect on the gas abundance. Changing M* indirectly affects the halo profiles at different masses, as can be seen in the top right panel, but note that the delta-function component of the stellar halo profile is not shown in this plot. The power spectra involving stars are affects in a more scale-dependent manner as can be seen in the top-left two panels.

In the text
thumbnail Fig. C.8.

Same as Fig. C.1 but for σ*, which governs the width of the stellar-mass distribution around the peak (which is at M* = 1012.5h−1M here). As can be seen in the bottom-left panel, decreasing σ* shrinks the stellar-mass distribution, resulting in fewer stars in both higher and lower mass haloes than M*. The top-right panel shows that this has the biggest impact on 1014h−1M haloes, not higher-mass haloes due to the saturation we impose at high masses. In the left-hand panels we see that the effects of this are constrained to be on power spectra involving the stars. There is only a tiny effect on gas spectra which is due to baryonic mass conservation.

In the text
thumbnail Fig. C.9.

Same as Fig. C.1 but for η, which determines the split of stellar halo mass between central and satellite galaxies for M >  M*, where M* = 1012.5h−1M in this figure. In the lower-right panel we see that increasing η from the fiducial −0.3 means that more mass is put into central galaxies, whereas decreasing it places more mass in satellites. As a result of this, the amplitude of the stellar halo is changed, as can be seen in the top-right panel, with decreasing η putting more mass into the NFW part of the profile and thus boosting the halo amplitude. This has a scale-dependent effect on the matter–stars power spectrum, as can be seen in the top-left two panels.

In the text
thumbnail Fig. C.10.

Same as Fig. C.1 but for α. α governs deviations of the halo temperature from the virial relation. We see that this affects the pressure profiles (middle right) by uniformly scaling the amplitude, because the pressure is directly proportional to temperature; increasing the temperature increases the pressure. We can see the effect of this on the pressure power spectra in the lower left two panels where the effect is not uniform due to the complicated dependence of the power spectra on the underlying halo populations in the two- and one-halo terms, which have different scalings with temperature.

In the text
thumbnail Fig. C.11.

Same as Fig. C.1 but for Tw. Varying Tw only affects the temperature of the unbound gas that contributes to the two-halo term only. We see that as we increase the gas temperature we boost the two-halo term. Note that there is a floor to the effect that lowering Tw can have, which is because the one-halo term contributes at large scales in the pressure spectra too, and lowering the temperature of this gas beyond a certain point means that only the two-halo term becomes dwarfed by the one-halo term at large scales.

In the text
thumbnail Fig. D.1.

Normalised contribution to each of C() for κκ (top) and κy (bottom) as a function of scale k, redshift z and (comoving) distance r. The lensing is taken from KiDS galaxies for redshifts z = 0.1 → 0.9. Each of these curves is the normalised integrand in Eq. (D.7) re-expressed as a function of either k, z or r (so dC()/dk, dC()/dz, dC()/dr). In practice, we carry out our integration in r between the observer and z(r) = 9 and we checked that our integrations are robust to this upper limit. The power spectra being integrated are set to zero for k <  10−4h Mpc−1 and k >  102h Mpc−1 and again our results are robust to this. Note the generally broader and lower–z contributions for κy compared to those for κκ.

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.