Cosmology from the Thermal Sunyaev-Zel’dovich Power Spectrum: Primordial non-Gaussianity and Massive Neutrinos

# Cosmology from the Thermal Sunyaev-Zel’dovich Power Spectrum: Primordial non-Gaussianity and Massive Neutrinos

J. Colin Hill and Enrico Pajer Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544
Department of Physics, Princeton University, Princeton, NJ 08544
jch@astro.princeton.edu, epajer@princeton.edu
July 9, 2019
###### Abstract

We carry out a comprehensive analysis of the possible constraints on cosmological and astrophysical parameters achievable with measurements of the thermal Sunyaev-Zel’dovich (tSZ) power spectrum from upcoming full-sky CMB observations, with a particular focus on one-parameter extensions to the CDM standard model involving local primordial non-Gaussianity (described by ) and massive neutrinos (described by ). We include all of the relevant physical effects due to these additional parameters, including the change to the halo mass function and the scale-dependent halo bias induced by local primordial non-Gaussianity. We use the halo model to compute the tSZ power spectrum and provide a pedagogical derivation of the one- and two-halo terms in an appendix. We model the pressure profile of the intracluster medium (ICM) using a parametrized fit that agrees well with existing observations, and include uncertainty in the ICM modeling by including the overall normalization and outer logarithmic slope of the profile as free parameters. We compute forecasts for Planck, PIXIE, and a cosmic variance (CV)-limited experiment, using multifrequency subtraction to remove foregrounds and implementing two masking criteria based on the ROSAT and eROSITA cluster catalogs to reduce the significant CV errors at low multipoles. We find that Planck can detect the tSZ power spectrum with significance, regardless of the masking scenario. However, neither Planck or PIXIE is likely to provide competitive constraints on from the tSZ power spectrum due to CV noise at low- overwhelming the unique signature of the scale-dependent bias. A future CV-limited experiment could provide a detection of , which is the WMAP9 maximum-likelihood result. The outlook for neutrino masses is more optimistic: Planck can reach levels comparable to the current upper bounds eV with conservative assumptions about the ICM; stronger ICM priors could allow Planck to provide evidence for massive neutrinos from the tSZ power spectrum, depending on the true value of the sum of the neutrino masses. We also forecast a % constraint on the outer slope of the ICM pressure profile using the unmasked Planck tSZ power spectrum.

## I Introduction

The thermal Sunyaev-Zel’dovich (tSZ) effect is a spectral distortion of the cosmic microwave background (CMB) that arises due to the inverse Compton scattering of CMB photons off hot electrons that lie between our vantage point and the surface of last scattering Sunyaev & Zeldovich (1970). The vast majority of these hot electrons are located in the intracluster medium (ICM) of galaxy clusters, and thus the tSZ signal is dominated by contributions from these massive objects. The tSZ effect has been used for many years to study individual clusters in pointed observations (e.g., Reese et al. (2011); Plagge et al. (2012); Lancaster et al. (2011); AMI Consortium et al. (2011)) and in recent years has been used as a method with which to find and characterize massive clusters in blind millimeter-wave surveys Marriage et al. (2011); Williamson et al. (2011); Planck Collaboration et al. (2011); Hasselfield et al. (2013); Reichardt et al. (2013). Moreover, recent years have brought the first detections of the angular power spectrum of the tSZ effect through its contribution to the power spectrum in arcminute-resolution maps of the microwave sky made by the Atacama Cosmology Telescope (ACT) and the South Pole Telescope (SPT) Dunkley et al. (2011); Sievers et al. (2013); Reichardt et al. (2012); Story et al. (2012). In addition, three-point statistics of the tSZ signal have been detected within the past year: first, the real-space skewness was detected in ACT data Wilson et al. (2012) using methods first anticipated by Rubiño-Martín & Sunyaev (2003), and second, the Fourier-space bispectrum was very recently detected in SPT data Crawford et al. (2013). The amplitudes of these measurements were shown to be consistent in the SPT analysis, despite observing different regions of sky and using different analysis methods. Note that the tSZ signal is highly non-Gaussian since it is dominated by contributions from massive collapsed objects in the late-time density field; thus, higher-order tSZ statistics contain significant information beyond that found in the power spectrum. Furthermore, the combination of multiple different -point tSZ statistics provides an avenue to extract tighter constraints on cosmological parameters and the astrophysics of the ICM than the use of the power spectrum alone, through the breaking of degeneracies between ICM and cosmological parameters Hill & Sherwin (2013); Bhattacharyaetal2012 (). The recent SPT bispectrum detection used such methods in order to reduce the error bar on the tSZ power spectrum amplitude by a factor of two Crawford et al. (2013).

Thus far, tSZ power spectrum detections have been limited to measurements or constraints on the power at a single multipole (typically ) because ACT and SPT do not have sufficient frequency coverage to fully separate the tSZ signal from other components in the microwave sky using its unique spectral signature. However, this situation will shortly change with the imminent release of full-sky temperature maps from the Planck satellite. Planck has nine frequency channels that span the spectral region around the tSZ null frequency GHz. Thus, it should be possible to separate the tSZ signal from other components in the sky maps to high accuracy, allowing for a measurement of the tSZ power spectrum over a wide range of multipoles, possibly , as we demonstrate in this paper. The proposed Primordial Inflation Explorer (PIXIE) experiment Kogut et al. (2011) will also be able to detect the tSZ power spectrum at high significance, as its wide spectral coverage and high spectral resolution will allow for very accurate extraction of the tSZ signal. However, its angular resolution is much lower than Planck’s, and thus the tSZ power spectrum will be measured over a much smaller range of multipoles (). However, the PIXIE data (after masking using X-ray cluster catalogs — see below) will permit tSZ measurements on large angular scales () that are essentially inaccessible to Planck due to its noise levels; these multipoles are precisely where one would expect the signature of scale-dependent bias induced by primordial non-Gaussianity to arise in the tSZ power spectrum. Assessing the amplitude and detectability of this signature is a primary motivation for this paper.

The tSZ power spectrum has been suggested as a potential cosmological probe by a number of authors over the past few decades (e.g., Cole & Kaiser (1988); Komatsu & Kitayama (1999); Komatsu & Seljak (2002)). Nearly all studies in the last decade have focused on the small-scale tSZ power spectrum () due to its role as a foreground in high angular resolution CMB measurements, and likely because without multi-frequency information, the tSZ signal has only been able to be isolated by looking for its effects on small scales (e.g., using an -space filter to upweight tSZ-dominated small angular scales in CMB maps). Much of this work was driven by the realization that the tSZ power spectrum is a very sensitive probe of the amplitude of matter density fluctuations,  Komatsu & Seljak (2002). The advent of multi-frequency data promises measurements of the large-scale tSZ power spectrum very shortly, and thus we believe it is timely to reassess its value as a cosmological probe, including parameters beyond and including a realistic treatment of the uncertainties due to modeling of the ICM. We build on the work of Komatsu & Kitayama (1999) to compute the full angular power spectrum of the tSZ effect, including both the one- and two-halo terms, and moving beyond the Limber/flat-sky approximations where necessary.

Our primary interest is in assessing constraints from the tSZ power spectrum on currently unknown parameters beyond the CDM standard model: the amplitude of local primordial non-Gaussianity, , and the sum of the neutrino masses, . The values of these parameters are currently unknown, and determining their values is a key goal of modern cosmology.

Primordial non-Gaussianity is one of the few known probes of the physics of inflation. Models of single-field, minimally-coupled slow-roll inflation predict negligibly small deviations from Gaussianity in the initial curvature perturbations Acquaviva et al. (2003); Maldacena (2003). In particular, a detection of a non-zero bispectrum amplitude in the so-called “squeezed” limit () would falsify essentially all single-field models of inflation Maldacena (2003); Creminelli & Zaldarriaga (2004). This type of non-Gaussianity can be parametrized using the “local” model, in which describes the lowest-order deviation from Gaussianity Salopek & Bond (1990); Gangui et al. (1994); Komatsu & Spergel (2001):

 Φ(→x)=ΦG(→x)+fNL(Φ2G(→x)−⟨Φ2G⟩)+⋯, (1)

where is the primordial potential and is a Gaussian field. Note that , where is the initial adiabatic curvature perturbation. Local-type non-Gaussianity can be generated in multi-field inflationary scenarios, such as the curvaton model Linde & Mukhanov (1997); Lyth & Wands (2002); Lyth et al. (2003), or by non-inflationary models for the generation of perturbations, such as the new ekpyrotic/cyclic scenario Buchbinder et al. (2008); Creminelli & Senatore (2007); Lehners & Steinhardt (2008). It is perhaps most interesting when viewed as a method with which to rule out single-field inflation, however. Current constraints on are consistent with zero Bennett et al. (2012); Giannantonio et al. (2013), but the errors will shrink significantly very soon with the imminent CMB results from Planck. We review the effects of on the large-scale structure of the universe in Section II.

In contrast to primordial non-Gaussianity, massive neutrinos are certain to exist at a level that will be detectable within the next decade or so; neutrino oscillation experiments have precisely measured the differences between the squared masses of the three known species, leading to a lower bound of eV on the total summed mass McKeown & Vogel (2004). The remaining questions surround their detailed properties, especially their absolute mass scale. The contribution of massive neutrinos to the total energy density of the universe today can be expressed as

 Ων≈Mν93.14h2eV≈0.0078Mν0.1eV, (2)

where is the sum of the masses of the three known neutrino species. Although this contribution appears to be small, massive neutrinos can have a significant influence on the small-scale matter power spectrum. Due to their large thermal velocities, neutrinos free-stream out of gravitational potential wells on scales below their free-streaming scale Lesgourgues & Pastor (2006); Abazajian et al. (2011). This suppresses power on scales below the free-streaming scale. Current upper bounds from various cosmological probes assuming a flat CDM+model are in the range eV Sievers et al. (2013); Vikhlinin et al. (2009); Mantz et al. (2010b); de Putter et al. (2012), although a detection near this mass scale was recently claimed in Hou et al. (2012). Should the true total mass turn out to be near eV, its effect on the tSZ power spectrum may be marginally detectable using the Planck data even for fairly conservative assumptions about the ICM physics, as we show in this paper. With stronger ICM priors, Planck could achieve a detection for masses at this scale, using only the primordial CMB temperature power spectrum and the tSZ power spectrum. We review the effects of massive neutrinos on the large-scale structure of the universe in Section II.

In addition to the effects of both known and currently unknown cosmological parameters, we also model the effects of the physics of the ICM on the tSZ power spectrum. This subject has attracted intense scrutiny in recent years after the early measurements of tSZ power from ACT Dunkley et al. (2011) and SPT Lueker et al. (2010) were significantly lower than the values predicted from existing ICM pressure profile models (e.g., Komatsu & Seljak (2001)) in combination with WMAP5 cosmological parameters. Subsequent ICM modeling efforts have ranged from fully analytic approaches (e.g., Shaw et al. (2010)) to cosmological hydrodynamics simulations (e.g., Battaglia et al. (2010, 2012)), with other authors adopting semi-analytic approaches between these extremes (e.g., Trac et al. (2011); Bode et al. (2012)), in which dark matter-only -body simulations are post-processed to include baryonic physics according to various prescriptions. In addition, recent SZ and X-ray observations have continued to further constrain the ICM pressure profile from data, although these results are generally limited to fairly massive, nearby systems (e.g., Arnaud et al. (2010); Planck Collaboration et al. (2013)). We choose to adopt a parametrized form of the ICM pressure profile known as the Generalized NFW (GNFW) profile, with our fiducial parameter values chosen to match the constrained pressure profile fit from hydrodynamical simulations in Battaglia et al. (2012). In order to account for uncertainty in the ICM physics, we free two of the parameters in the pressure profile (the overall normalization and the outer logarithmic slope) and treat them as additional parameters in our model. This approach is discussed in detail in Section III.2. We use this profile to compute the tSZ power spectrum following the halo model approach, for which we provide a complete derivation in Appendix A.

In addition to a model for the tSZ signal, we must compute the tSZ power spectrum covariance matrix in order to forecast parameter constraints. There are two important issues that must be considered in computing the expected errors or signal-to-noise ratio (SNR) for a measurement of the tSZ power spectrum. First, we must assess how well the tSZ signal can actually be separated from the other components in maps of the microwave sky, including the primordial CMB, thermal dust, point sources, and so on. Following THEO () and CHT (), we choose to implement a multi-frequency subtraction technique that takes advantage of the unique spectral signature of the tSZ effect, and also takes advantage of the current state of knowledge about the frequency- and multipole-dependence of the foregrounds. Although there are other approaches to this problem, such as using internal linear combination techniques to construct a Compton- map from the individual frequency maps in a given experiment (e.g., Remazeilles et al. (2013, 2011); Leach et al. (2008)), we find this method to be fairly simple and robust. We describe these calculations in detail in Section IV.

Second, we must account for the extreme cosmic variance induced in the large-angle tSZ power spectrum by massive clusters at low redshifts. The one-halo term from these objects dominates the angular trispectrum of the tSZ signal, even down to very low multipoles Cooray (2001). The trispectrum represents a large non-Gaussian contribution to the covariance matrix of the tSZ power spectrum Komatsu & Seljak (2002), which is especially problematic at low multipoles. However, the trispectrum can be greatly suppressed by masking massive, low-redshift clusters using existing X-ray, optical, or SZ catalogs Komatsu & Kitayama (1999). This procedure can greatly increase the SNR for the tSZ power spectrum at low multipoles. For constraints on it also has the advantage of enhancing the relative importance of the two-halo term compared to the one-halo, thus showing greater sensitivity to the scale-dependent bias at low-. Moreover, even in a Gaussian cosmology, the inclusion of the two-halo term slightly changes the shape of the tSZ power spectrum, which likely helps break degeneracies amongst the several parameters which effectively only change the overall amplitude of the one-halo term; the relative enhancement of the two-halo term due to masking should help further in this regard. We consider two masking scenarios motivated by the flux limits of the cluster catalogs from all-sky surveys performed with the ROSAT X-ray telescope and the upcoming eROSITA555Extended ROentgen Survey with an Imaging Telescope Array, http://www.mpe.mpg.de/erosita/ X-ray telescope. These scenarios are detailed in Section V.0.2; by default all calculations and figures are computed for the unmasked scenario unless they are labelled otherwise.

Earlier studies have investigated the consequences of primordial non-Gaussianity for the tSZ power spectrum Sadeh et al. (2007); Roncarelli et al. (2010), though we are not aware of any calculations including the two-halo term (and hence the scale-dependent bias) or detailed parameter constraint forecasts. We are also not aware of any previous work investigating constraints on massive neutrinos from the tSZ power spectrum, although previous authors have computed their signature Shimon et al. (2011). Other studies have investigated detailed constraints on the primary CDM parameters from the combination of CMB and tSZ power spectrum measurements Taburet et al. (2010). Many authors have investigated constraints on and from cluster counts, though the results depend somewhat on the cluster selection technique and mass estimation method. Considering SZ cluster count studies only, Shimon et al. (2011) and Shimon et al. (2012) investigated constraints on from a Planck-derived catalog of SZ clusters (in combination with CMB temperature power spectrum data). The earlier paper found a uncertainty of eV while the later paper found eV; the authors state that the use of highly degenerate nuisance parameters degraded the results in the former study. In either case, the result is highly sensitive to uncertainties in the halo mass function, as the clusters included are deep in the exponential tail of the mass function. We expect that our results using the tSZ power spectrum should be less sensitive to uncertainties in the tail of the mass function, as the power spectrum is dominated at most angular scales by somewhat less massive objects (Komatsu & Seljak (2002). Finally, a very recent independent study Mak & Pierpaoli (2013) found eV for Planck SZ cluster counts (with CMB temperature power spectrum information added), although they estimated that this bound could be improved to eV with the inclusion of stronger priors on the ICM physics.

Our primary findings are as follows:

• The tSZ power spectrum can be detected with a total SNR using the imminent Planck data up to , regardless of masking;

• The tSZ power spectrum can be detected with a total SNR between and 22 using the future PIXIE data up to , with the result being sensitive to the level of masking applied to remove massive, nearby clusters;

• Adding the tSZ power spectrum information to the forecasted constraints from the Planck CMB temperature power spectrum and existing data is unlikely to significantly improve constraints on the primary cosmological parameters, but may give interesting constraints on the extensions we consider:

• If the true value of is near the WMAP9 ML value of , a future CV-limited experiment combined with eROSITA-masking could provide a detection, completely independent of the primordial CMB temperature bispectrum; alternatively, PIXIE could give evidence for such a value of with this level of masking;

• If the true value of is near 0.1 eV, the Planck tSZ power spectrum with eROSITA masking can provide upper limits competitive with the current upper bounds on ; with stronger external constraints on the ICM physics, Planck with eROSITA masking could provide evidence for massive neutrinos from the tSZ power spectrum, depending on the true neutrino mass;

• Regardless of the cosmological constraints, Planck will allow for a very tight constraint on the logarithmic slope of the ICM pressure profile in the outskirts of galaxy clusters, and may also provide some information on the overall normalization of the pressure profile (which sets the zero point of the relation).

The remainder of this paper is organized as follows. In Section II, we describe our models for the halo mass function and halo bias, as well as the effects of primordial non-Gaussianity and massive neutrinos on large-scale structure. In Section III, we describe our halo model-based calculation of the tSZ power spectrum, including the relevant ICM physics. We also demonstrate the different effects of each parameter in our model on the tSZ power spectrum. In Section IV, we consider the extraction of the tSZ power spectrum from the other components in microwave sky maps via multifrequency subtraction techniques. Having determined the experimental noise levels, in Section V we detail our calculation of the covariance matrix of the tSZ power spectrum, and discuss the role of masking massive nearby clusters in reducing the low- cosmic variance. In Section VI, we use our tSZ results to forecast constraints on cosmological and astrophysical parameters from a variety of experimental set-ups and masking choices. We also compute the expected SNR of the tSZ power spectrum detection for each possible scenario. We discuss our results and conclude in Section VII. Finally, in Appendix B, we provide a brief comparison between our forecasts and the Planck tSZ power spectrum results that were publicly released while this manuscript was under review Planck Collaboration et al. (2013).

The WMAP9+eCMB+BAO+ maximum-likelihood parameters Hinshaw et al. (2012) define our fiducial model (see Section III.3 for details). All masses are quoted in units of , where and is the Hubble parameter today. All distances and wavenumbers are in comoving units of . All tSZ observables are computed at GHz, since ACT and SPT have observed the tSZ signal at (or very near) this frequency, where the tSZ effect leads to a temperature decrement in the CMB along the line-of-sight (LOS) to a galaxy cluster.

## Ii Modeling Large-Scale Structure

In order to compute statistics of the tSZ signal, we need to model the comoving number density of halos as a function of mass and redshift (the halo mass function) and the bias of halos with respect to the underlying matter density field as a function of mass and redshift. Moreover, in order to extract constraints on and from the tSZ power spectrum, we must include the effects of these parameters on large-scale structure. We describe our approach to these computations in the following.

### ii.1 Halo Mass Function

We define the mass of a dark matter halo by the spherical overdensity (SO) criterion: () is the mass enclosed within a sphere of radius () such that the enclosed density is times the critical (mean matter) density at redshift . To be clear, subscripts refer to masses referenced to the critical density at redshift , with the Hubble parameter at redshift , whereas subscripts refer to masses referenced to the mean matter density at redshift , (this quantity is constant in comoving units).

We will generally work in terms of a particular SO mass, the virial mass, which we denote as . The virial mass is the mass enclosed within a radius  Bryan & Norman (1998):

 rvir=(3M4πΔcr(z)ρcr(z))1/3, (3)

where and . For many calculations, we need to convert between and various other SO masses (e.g., or ). We use the NFW density profile Navarro et al. (1997) and the concentration-mass relation from Duffy et al. (2008) in order to do these conversions, which require solving the following non-linear equation for (or ):

 ∫rδ,c04πr′2ρNFW(r′,M,cvir)dr′=43πr3δ,cρcr(z)δ (4)

where is the concentration parameter ( is the NFW scale radius) and we replace the critical density with the mean matter density in this equation in order to obtain instead of . After solving Eq. (4) to find , we calculate via .

The halo mass function, describes the comoving number density of halos per unit mass as a function of redshift. We employ the approach developed from early work by Press and Schechter Press & Schechter (1974) and subsequently refined by many other authors (e.g., Sheth & Tormen (1999); Sheth et al. (2001); Jenkins et al. (2001); Tinker et al. (2008)):

 dn(M,z)dM = ¯ρmMdln(σ−1(M,z))dMf(σ(M,z)) (5) = −¯ρm2M2R(M)3σ2(M,z)dσ2(M,z)dR(M)f(σ(M,z)),

where is the variance of the linear matter density field smoothed with a (real space) top-hat filter on a scale at redshift :

 σ2(M,z)=12π2∫k3Plin(k,z)W2(k,R(M))dlnk, (6)

where is the linear theory matter power spectrum at wavenumber and redshift . Note that the window function is a top-hat filter in real space, which in Fourier space is given by

 W(k,R)=3x2(sinxx−cosx), (7)

where . In Eq. (5), the function is known as the halo multiplicity function. It has been measured to increasingly high precision from large -body simulations over the past decade Jenkins et al. (2001); Warren et al. (2006); Tinker et al. (2008); Bhattacharya et al. (2011). However, many of these calibrated mass functions are specified in terms of the friends-of-friends (FOF) mass rather than the SO mass, hindering their use in analytic calculations such as ours. For this reason, we use the parametrization and calibration from Tinker et al. (2008), where computations are performed in terms of the SO mass with respect to the mean matter density, , for a variety of overdensities. The halo multiplicity function in this model is parametrized by

 (8)

where are (redshift- and overdensity-dependent) parameters fit from simulations. We use the values of these parameters appropriate for the halo mass function from (Tinker et al., 2008) with the redshift-dependent parameters given in their Eqs. (5)–(8); we will hereafter refer to this as the Tinker mass function. Note that the authors of that study caution against extrapolating their parameters beyond the highest redshift measured in their simulations () and recommend setting the parameters equal to their values at higher redshifts; we adopt this recommendation in our calculations. Also, note that our tSZ power spectrum calculations in Section III are phrased in terms of the virial mass , and thus we compute the Jacobian using the procedure described in Eq. (4) in order to convert the Tinker mass function to a virial mass function .

We compute the smoothed matter density field in Eq. (6) by first obtaining the linear theory matter power spectrum from CAMB at and subsequently rescaling by , where is the linear growth factor. We normalize by requiring that deep in the matter-dominated era (e.g., at ). The resulting is then used to compute the mass function in Eq. (5).

Note that we assume the mass function to be known to high enough precision that the parameters describing it can be fixed; in other words, we do not consider to be free parameters in our model. These parameters are certainly better constrained at present than those describing the ICM pressure profile (see Section III.2), and thus this assumption seems reasonable for now. However, precision cosmological constraints based on the mass function should in principle consider variations in the mass function parameters in order to obtain robust results, as has been done in some recent X-ray cluster cosmology analyses Mantz et al. (2010a). However, we leave the implications of these uncertainties for tSZ statistics as a topic for future work.

#### ii.1.1 Effect of Primordial non-Gaussianity

The influence of primordial non-Gaussianity on the halo mass function has been studied by many authors over the past two decades using a variety of approaches (e.g., Lucchin & Matarrese (1988); Colafrancesco et al. (1989); Chiu et al. (1998); Robinson et al. (2000); Koyama et al. (1999); Verde et al. (2001); Matarrese et al. (2000); Lo Verde et al. (2008); D’Amico et al. (2011)). The physical consequences of the model specified in Eq. (1) are fairly simple to understand for the halo mass function, especially in the exponential tail of the mass function where massive clusters are found. Intuitively, the number of clusters provides information about the tail of the probability distribution function of the primordial density field, since these are the rarest objects in the universe, which have only collapsed recently. For positive skewness in the primordial density field (), one obtains an increased number of massive clusters at late times relative to the case, because more regions of the smoothed density field have , the collapse threshold ( in the spherical collapse model). Conversely, for negative skewness in the primordial density field (), one obtains fewer massive clusters at late times relative to the case, because fewer regions of the smoothed density field are above the collapse threshold. As illustrated in recent analytic calculations and simulation measurements LoVerde & Smith (2011); Pillepich et al. (2010); Grossi et al. (2009); Wagner et al. (2010); Grossi et al. (2007), these changes can be quite significant for the number of extremely massive halos () in the late-time universe; for example, the abundance of such halos for can be times larger than in a Gaussian cosmology. These results have been used as a basis for recent studies constraining by looking for extremely massive outliers in the cluster distribution (e.g., Hoyle et al. (2011); Cayón et al. (2011); Mortonson et al. (2011); Enqvist et al. (2011); Harrison & Hotchkiss (2012); Hoyle et al. (2012)).

We model the effect of on the halo mass function by multiplying the Tinker mass function by a non-Gaussian correction factor:

 (dndM)NG=dndMRNG(M,z,fNL), (9)

where is given by Eq. (5). We use the model for given by Eq. (35) in LoVerde & Smith (2011) (the “log-Edgeworth” mass function). In this approach, the density field is approximated via an Edgeworth expansion, which captures small deviations from Gaussianity. The Press-Schechter approach is then applied to the Edgeworth-expanded density field to obtain an expression for the halo mass function in terms of cumulants of the non-Gaussian density field. The results of LoVerde & Smith (2011) include numerical fitting functions for these cumulants obtained from -body simulations. We use both the expression for and the cumulant fitting functions from LoVerde & Smith (2011) to compute the non-Gaussian correction to the mass function. This prescription was shown to accurately reproduce the non-Gaussian halo mass function correction factor measured directly from -body simulations in LoVerde & Smith (2011), and in particular improves upon the similar prescription derived in Lo Verde et al. (2008) (the “Edgeworth” mass function).

Note that we apply the non-Gaussian correction factor to the Tinker mass function in Eq. (9), which is an SO mass function, as mentioned above. The prescription for computing makes no assumption about whether is an FOF or SO mass, so there is no logical flaw in this procedure. However, the comparisons to -body results in LoVerde & Smith (2011) were performed using FOF halos. Thus, without having tested the results of Eq. (9) on SO mass functions from simulations, our calculation assumes that the change in the mass function due to non-Gaussianity is quasi-universal, even if the underlying Gaussian mass function itself is not. This assumption was tested in Wagner et al. (2010) for the non-Gaussian correction factor from Lo Verde et al. (2008) (see Fig. 9 in Wagner et al. (2010)) and found to be valid; thus, we choose to adopt it here. We will refer to the non-Gaussian mass function computed via Eq. (9) using the prescription from LoVerde & Smith (2011) as the LVS mass function.

#### ii.1.2 Effect of Massive Neutrinos

It has long been known that massive neutrinos suppress the amplitude of the matter power spectrum on scales below their free-streaming scale,  Lesgourgues & Pastor (2006):

 kfs≈0.082H(z)H0(1+z)2(Mν0.1eV)h/Mpc. (10)

Neutrinos do not cluster on scales much smaller than this scale (i.e., ), as they are able to free-stream out of small-scale gravitational potential wells. This effect leads to a characteristic decrease in the small-scale matter power spectrum of order in linear perturbation theory Lesgourgues & Pastor (2006); Abazajian et al. (2011). Nonlinear corrections increase this suppression to for modes with wavenumbers  Abazajian et al. (2011).

The neutrino-induced suppression of the small-scale matter power spectrum leads one to expect that the number of massive halos in the low-redshift universe should also be decreased. Several papers in recent years have attempted to precisely model this change in the halo mass function using both -body simulations and analytic theory Brandbyge et al. (2010); Marulli et al. (2011); Ichiki & Takada (2012). In Brandbyge et al. (2010), -body simulations are used to show that massive neutrinos do indeed suppress the halo mass function, especially for the largest, latest-forming halos (i.e., galaxy clusters). Moreover, the suppression is found to arise primarily from the suppression of the initial transfer function in the linear regime, and not due to neutrino clustering effects in the -body simulations. This finding suggests that an analytic approach similar to the Press-Schecter theory should work for massive neutrino cosmologies as well, and the authors subsequently show that a modified Sheth-Tormen formalism Sheth & Tormen (1999) gives a good fit to their simulation results. Similar -body simulations are examined in Marulli et al. (2011), who find generally similar results to those in Brandbyge et al. (2010), but also point out that the effect of on the mass function cannot be adequately represented by simply rescaling to a lower value in an analytic calculation without massive neutrinos. Finally, Ichiki & Takada (2012) study the effect of massive neutrinos on the mass function using analytic calculations with the spherical collapse model. Their results suggest that an accurate approximation is to simply input the -suppressed linear theory (coldbaryonic-only) matter power spectrum computed at to a CDM-calibrated mass function fit (note that a similar procedure was used in some recent X-ray cluster-based constraints on  Vikhlinin et al. (2009)). The net result of this suppression can be quite significant at the high-mass end of the mass function; for example, eV leads to a factor of decrease in the abundance of halos at as compared to a massless-neutrino cosmology Ichiki & Takada (2012). We follow the procedure used in Ichiki & Takada (2012) in our work, although we input the suppressed linear theory matter power spectrum to the Tinker mass function rather than that of Bhattacharya et al. (2011), as was done in Ichiki & Takada (2012). We will refer to the -suppressed mass function computed with this prescription as the IT mass function.

### ii.2 Halo Bias

Dark matter halos are known to cluster more strongly than the underlying matter density field; they are thus biased tracers. This bias can depend on scale, mass, and redshift (e.g., BBKS (); Mo & White (1996); Smith et al. (2007)). We define the halo bias by

 b(k,M,z)=√Phh(k,M,z)P(k,z), (11)

where is the power spectrum of the halo density field and is the power spectrum of the matter density field. Knowledge of the halo bias is necessary to model and extract cosmological information from the clustering of galaxies and galaxy clusters. For our purposes, it will be needed to compute the two-halo term in the tSZ power spectrum, which requires knowledge of .

In a Gaussian cosmology, the halo bias depends on mass and redshift but is independent of scale for , i.e. on large scales (e.g., Tinker et al. (2010)). We compute this linear Gaussian bias, , using the fitting function in Eq. (6) of Tinker et al. (2010) with the parameters appropriate for SO masses (see Table 2 in Tinker et al. (2010)). This fit was determined from the results of many large-volume -body simulations with a variety of cosmological parameters and found to be quite accurate. We will refer to this prescription as the Tinker bias model.

Although the bias becomes scale-dependent on small scales even in a Gaussian cosmology, it becomes scale-dependent on large scales in the presence of local primordial non-Gaussianity, as first shown in Dalal et al. (2008). The scale-dependence arises due to the coupling of long- and short-wavelength density fluctuations induced by local . We model this effect as a correction to the Gaussian bias described in the preceding paragraph:

 b(k,M,z)=bG(M,z)+ΔbNG(k,M,z), (12)

where the non-Gaussian correction is given by Dalal et al. (2008)

 ΔbNG(k,M,z)=2δc(bG(M,z)−1)fNLα(k,z). (13)

Here, (the spherical collapse threshold) and

 α(k,z)=2k2T(k)D(z)c23ΩmH20 (14)

relates the linear density field to the primordial potential via . Note that is the linear matter transfer function, which we compute using CAMB. Since the original derivation in Dalal et al. (2008), the results in Eqs. (13) and (14) have subsequently been confirmed by other authors Slosar et al. (2008); Matarrese & Verde (2008); Giannantonio & Porciani (2010) and tested extensively on -body simulations (e.g., Desjacques et al. (2009); Dalal et al. (2008); Pillepich et al. (2010); Smith et al. (2012)). The overall effect is a steep increase in the large-scale bias of massive halos, which is even larger for highly biased tracers like galaxy clusters. We will refer to this effect simply as the scale-dependent halo bias.

The influence of massive neutrinos on the halo bias has been studied far less thoroughly than that of primordial non-Gaussianity. Recent -body simulations analyzed in Marulli et al. (2011) indicate that massive neutrinos lead to a nearly scale-independent increase in the large-scale halo bias. This effect arises because of the mass function suppression discussed in Section II.1.2: halos of a given mass are rarer in an cosmology than in a massless neutrino cosmology (for fixed ), and thus they are more highly biased relative to the matter density field. However, the amplitude of this change is far smaller than that induced by , especially on very large scales. For example, the results of Marulli et al. (2011) indicate an overall increase of % in the mean bias of massive halos at for eV as compared to . Our implementation of the scale-dependent bias due to local yields a factor of increase in the large-scale () bias of objects in the same mass range at for . Clearly, the effect of is much larger than that of massive neutrinos, simply because it is so strongly scale-dependent, while only leads to a small scale-independent change (at least on large scales; the small-scale behavior may be more complicated). Moreover, the change in bias due to is larger at higher redshifts (), whereas most of the tSZ signal originates at lower redshifts. Lastly, due to the smallness of the two-halo term in the tSZ power spectrum compared to the one-halo term (see Section III), small variations in the Gaussian bias cause essentially no change in the total signal. For all of these reasons, we choose to neglect the effect of massive neutrinos on the halo bias in our calculations.

## Iii Thermal SZ Power Spectrum

The tSZ effect results in a frequency-dependent shift in the CMB temperature observed in the direction of a galaxy group or cluster. The temperature shift at angular position with respect to the center of a cluster of mass at redshift is given by Sunyaev & Zeldovich (1970)

 ΔT(→θ,M,z)TCMB = gνy(→θ,M,z) = gνσTmec2∫LOSPe(√l2+d2A|→θ|2,M,z)dl,

where is the tSZ spectral function with , is the Compton- parameter, is the Thomson scattering cross-section, is the electron mass, and is the ICM electron pressure at location with respect to the cluster center. We have neglected relativistic corrections in Eq. (III) (e.g., Nozawa et al. (2006)), as these effects are relevant only for the most massive clusters in the universe (). Such clusters contribute non-negligibly to the tSZ power spectrum at low-, and thus our results in unmasked calculations may be slightly inaccurate; however, the optimal forecasts for cosmological constraints arise from calculations in which such nearby, massive clusters are masked (see Section VI), and thus these corrections will not be relevant. Therefore, we do not include them in our calculations.

Note that we only consider spherically symmetric pressure profiles in this work, i.e. in Eq. (III). The integral in Eq. (III) is computed along the LOS such that , where is the angular diameter distance to redshift and is the angular distance between and the cluster center in the plane of the sky (note that this formalism assumes the flat-sky approximation is valid; we provide exact full-sky results for the tSZ power spectrum in Appendix A). In the flat-sky limit, a spherically symmetric pressure profile implies that the temperature decrement (or Compton-) profile is azimuthally symmetric in the plane of the sky, i.e., . Finally, note that the electron pressure is related to the thermal gas pressure via , where is the primordial hydrogen mass fraction. We calculate all tSZ power spectra in this paper at GHz, where the tSZ effect is observed as a decrement in the CMB temperature (). We make this choice simply because recent tSZ measurements have been performed at this frequency using ACT and SPT (e.g., Wilson et al. (2012); Crawford et al. (2013); Story et al. (2012); Sievers et al. (2013)), and thus the temperature values in this regime are perhaps more familiar and intuitive. All of our calculations can be phrased in a frequency-independent manner in terms of the Compton- parameter, and we will often use “y” as a label for tSZ quantities, although they are calculated numerically at GHz.

In the remainder of this section, we outline the halo model-based calculations used to compute the tSZ power spectrum, discuss our model for the gas physics of the ICM, and explain the physical effects of each cosmological and astrophysical parameter on the tSZ power spectrum.

### iii.1 Halo Model Formalism

We compute the tSZ power spectrum using the halo model approach (see Cooray & Sheth (2002) for a review). We provide complete derivations of all the relevant expressions in Appendix A, first obtaining completely general full-sky results and then specializing to the flat-sky/Limber-approximated case. Here, we simply quote the necessary results and refer the interested reader to Appendix A for the derivations. Note that we will work in terms of the Compton- parameter; the results can easily be multiplied by the necessary factors to obtain results at any frequency.

The tSZ power spectrum, , is given by the sum of the one-halo and two-halo terms:

 Cyℓ=Cy,1hℓ+Cy,2hℓ. (16)

The exact expression for the one-halo term is given by Eq. (76):

 Cy,1hℓ=∫dzχ(z)d2VdzdΩ∫dMdndM∣∣ ∣∣∫kdkJℓ+1/2(kχ(z))~y3D(k;M,z)∫cdz′H(z′)(1+z′)√χ(z′)Jℓ+1/2(kχ(z′))∣∣ ∣∣2, (17)

where is the comoving distance to redshift , is the comoving volume element per steradian, is the halo mass function discussed in Section II.1, is given in Eq. (68), and is a Bessel function of the first kind. In the flat-sky limit, the one-halo term simplifies to the following widely-used expression (given in e.g. Eq. (1) of Komatsu & Seljak (2002)), which we derive in Eq. (81):

 Cy,1hℓ≫1≈∫dzd2VdzdΩ∫dMdn(M,z)dM|~yℓ(M,z)|2, (18)

where

 ~yℓ(M,z)≈4πrsℓ2s∫dxx2sin((ℓ+1/2)x/ℓs)(ℓ+1/2)x/ℓsy3D(x;M,z). (19)

Here, is a characteristic scale radius (not the NFW scale radius) of the profile given by and is the multipole moment associated with the scale radius. For the pressure profile from Battaglia et al. (2012) used in our calculations, the natural scale radius is . In our calculations, we choose to implement the flat-sky result for the one-halo term at all — see Appendix A for a justification of this decision and an assessment of the associated error at low- (the only regime where this correction would be relevant).

The exact expression for the two-halo term is given by Eq. (82):

 (20)

where is the linear theory matter power spectrum at (which we choose to set equal to 30), is the halo bias discussed in Section II.2, and refers to the expression for given in Eq. (19) evaluated with . This notation is simply a mathematical convenience; no flat-sky or Limber approximation was used in deriving Eq. (82), and no appears in . In the Limber approximation Limber (1954), the two-halo term simplifies to the result given in Komatsu & Kitayama (1999), which we derive in Eq. (84):

 Cy,2hℓ≫1≈∫dzd2VdzdΩ[∫dMdn(M,z)dMb(k,M,z)~yℓ(M,z)]2Plin(ℓ+1/2χ(z);z). (21)

We investigate the validity of the Limber approximation in detail in Appendix A. We find that it is necessary to compute the exact expression in Eq. (20) in order to obtain sufficiently accurate results at low-, where the signature of the scale-dependent bias induced by is present (looking for this signature is our primary motivation for computing the two-halo term to begin with). In particular, we compute the exact expression in Eq. (20) for , while we use the Limber-approximated result in Eq. (21) at higher multipoles.

The fiducial integration limits in our calculations are for all redshift integrals, for all mass integrals, and for all wavenumber integrals. We check that extending the wavenumber upper limit further into the nonlinear regime does not affect our results. Note that the upper limit in the mass integral becomes redshift-dependent in the masked calculations that we discuss below, in which the most massive clusters at low redshifts are removed from the computation.

We use the halo mass functions discussed in Section II.1 (Tinker, LVS, and IT) and the bias models discussed in Section II.2 (Tinker and scale-dependent bias) in Eqs. (18), (20), and (21). The only remaining ingredient needed to complete the tSZ power spectrum calculation is a prescription for the ICM electron pressure profile as a function of mass and redshift. Note that this approach to the tSZ power spectrum calculation separates the cosmology-dependent component (the mass function and bias) from the ICM-dependent component (the pressure profile). This separation arises from the fact that the small-scale baryonic physics that determines the structure of the ICM pressure profile effectively decouples from the large-scale physics described by the background cosmology and linear perturbation theory. Thus, it is a standard procedure to constrain the ICM pressure profile from cosmological hydrodynamics simulations (e.g., Battaglia et al. (2012); Bode et al. (2012)) or actual observations of galaxy clusters (e.g., Arnaud et al. (2010); Planck Collaboration et al. (2013), which are obtained for a fixed cosmology in either case (at present, it is prohibitively computationally expensive to run many large hydrodynamical simulations with varying cosmological parameters). Of course, it is also possible to model the ICM analytically and obtain a pressure profile (e.g., Komatsu & Seljak (2001); Shaw et al. (2010). Regardless of its origin (observations/simulations/theory), the derived ICM pressure profile can then be applied to different background cosmologies by using the halo mass function and bias model appropriate for that cosmology in the tSZ power spectrum calculations. We follow this approach.

Note that because the tSZ signal is heavily dominated by contributions from collapsed objects, the halo model approximation gives very accurate results when compared to direct LOS integrations of numerical simulation boxes (see Figs. 7 and 8 in Battaglia et al. (2012) for direct comparisons). In particular, the halo model agrees very well with the simulation results for , which is predominantly the regime we are interested in for this paper (on smaller angular scales effects due to asphericity and substructure become important, which are not captured in the halo model approach). These results imply that contributions from the intergalactic medium, filaments, and other diffuse structures are unlikely to be large enough to significantly impact the calculations and forecasts in the remainder of the paper. Contamination from the Galaxy is a separate issue, which we assume can be minimized to a sufficient level through sky cuts and foreground subtraction (see Section IV).

### iii.2 Modeling the ICM

We adopt the parametrized ICM pressure profile fit from Battaglia et al. (2012) as our fiducial model. This profile is derived from cosmological hydrodynamics simulations described in Battaglia et al. (2010). These simulations include (sub-grid) prescriptions for radiative cooling, star formation, supernova feedback, and feedback from active galactic nuclei (AGN). Taken together, these feedback processes typically decrease the gas fraction in low-mass groups and clusters, as the injection of energy into the ICM blows gas out of the cluster potential. In addition, the smoothed particle hydrodynamics used in these simulations naturally captures the effects of non-thermal pressure support due to bulk motions and turbulence, which must be modeled in order to accurately characterize the cluster pressure profile in the outskirts.

The ICM thermal pressure profile in this model is parametrized by a dimensionless GNFW form, which has been found to be a useful parametrization by many observational and numerical studies (e.g., Nagai et al. (2007); Arnaud et al. (2010); Plagge et al. (2012); Planck Collaboration et al. (2013)):

 Pth(x)P200,c=P0(x/xc)γ[1+(x/xc)α]β,x≡r/r200,c, (22)

where is the thermal pressure profile, is the dimensionless distance from the cluster center, is a core scale length, is a dimensionless amplitude, , , and describe the logarithmic slope of the profile at intermediate (), large (), and small () radii, respectively, and is the self-similar amplitude for pressure at given by Kaiser (1986); Voit (2005):

 P200,c=200GM200,cρcr(z)Ωb2Ωmr200,c. (23)

In Battaglia et al. (2012) this parametrization is fit to the stacked pressure profiles of clusters extracted from the simulations described above. Note that due to degeneracies the parameters and are not varied in the fit; they are fixed to and , which agree with many other studies (e.g., Nagai et al. (2007); Arnaud et al. (2010); Plagge et al. (2012); Planck Collaboration et al. (2013). In addition to constraining the amplitude of the remaining parameters, Battaglia et al. (2012) also fit power-law mass and redshift dependences, with the following results:

 P0(M200,c,z) = 18.1(M200,c1014M⊙)0.154(1+z)−0.758 (24) xc(M200,c,z) = 0.497(M200,c1014M⊙)−0.00865(1+z)0.731 (25) β(M200,c,z) = 4.35(M200,c1014M⊙)0.0393(1+z)0.415. (26)

Note that the denominator of the mass-dependent factor has units of rather than as used elsewhere in this paper. The mass and redshift dependence of these parameters captures deviations from simple self-similar cluster pressure profiles. These deviations arise from non-gravitational energy injections due to AGN and supernova feedback, star formation in the ICM, and non-thermal processes such as turbulence and bulk motions Battaglia et al. (2012, 2012). Eqs. (22)–(26) completely specify the ICM electron pressure profile as a function of mass and redshift, and provide the remaining ingredient needed for the halo model calculations of the tSZ power spectrum described in Section III.1, in addition to the halo mass function and halo bias. We will refer to this model of the ICM pressure profile as the Battaglia model.

Although it is derived solely from numerical simulations, we note that the Battaglia pressure profile is in good agreement with a number of observations of cluster pressure profiles, including those based on the REXCESS X-ray sample of massive, clusters Arnaud et al. (2010), independent studies of low-mass groups at with Chandra Sun et al. (2011), and early Planck measurements of the stacked pressure profile of clusters Planck Collaboration et al. (2013).

We allow for a realistic degree of uncertainty in the ICM pressure profile by freeing the amplitude of the parameters that describe the overall normalization () and the outer logarithmic slope (). To be clear, we do not free the mass and redshift dependences for these parameters given in Eqs. (24) and (26), only the overall amplitudes in those expressions. The outer slope is known to be highly degenerate with the scale radius  (e.g., Battaglia et al. (2012); Plagge et al. (2012)), and thus it is only feasible to free one of these parameters. The other slope parameters in Eq. (22) are fixed to their Battaglia values, which match the standard values in the literature. We parametrize the freedom in and by introducing new parameters and defined by:

 P0(M200,c,z) = CP0×18.1(M200,c1014M⊙)0.154(1+z)−0.758 (27) β(M200,c,z) = Cβ×4.35(M200,c1014M⊙)0.0393(1+z)0.415. (28)

These parameters thus describe multiplicative overall changes to the amplitudes of the and parameters. The fiducial Battaglia profile corresponds to . We discuss our priors for these parameters in Section VI.

### iii.3 Parameter Dependences

Including both cosmological and astrophysical parameters, our model is specified by the following quantities:

 {Ωbh2,Ωch2,ΩΛ,σ8,ns,CP0,Cβ,(fNL,Mν)}, (29)

which take the following values in our (WMAP9+BAO+ Hinshaw et al. (2012)) fiducial model:

 {0.02240,0.1146,0.7181,0.817,0.9646,1.0,1.0,(0.0,0.0)}. (30)

As a reminder, the CDM parameters are (in order of their appearance in Eq. (29)) the physical baryon density, the physical cold dark matter density, the vacuum energy density, the rms matter density fluctuation on comoving scales of at , and the scalar spectral index. The ICM physics parameters and are defined in Eqs. (27) and (28), respectively, is defined by Eq. (1), and is the sum of the neutrino masses in units of eV. We have placed and in parentheses in Eq. (29) in order to make it clear that we only consider scenarios in which these parameters are varied separately: for all cosmologies that we consider with , we set , and for all cosmologies that we consider with , we set . In other words, we only investigate one-parameter extensions of the CDM concordance model.

For the primary CDM cosmological parameters, we use the parametrization adopted by the WMAP team (e.g., Hinshaw et al. (2012)), as the primordial CMB data best constrain this set. The only exception to this convention is our use of , which stands in place of the primordial amplitude of scalar perturbations, . We use both because it is conventional in the tSZ power spectrum literature and because it is a direct measure of the low-redshift amplitude of matter density perturbations, which is physically related more closely to the tSZ signal than . However, this choice leads to slightly counterintuitive results when considering cosmologies with , because in order to keep fixed for such scenarios we must increase (to compensate for the suppression induced by in the matter power spectrum).

For the fiducial model specified by the values in Eq. (30), we find that the tSZ power spectrum amplitude at is at GHz. This corresponds to at GHz (the relevant ACT frequency) and at GHz (the relevant SPT frequency). The most recent measurements from ACT and SPT find corresponding constraints at these frequencies of  Sievers et al. (2013) and  Crawford et al. (2013) (using their more conservative error estimate). Note that the SPT constraint includes information from the tSZ bispectrum, which reduces the error by a factor of 2. Although it appears that our fiducial model predicts a level of tSZ power too high to be consistent with these observations, the results are highly dependent on the true value of , due to the steep dependence of the tSZ power spectrum on this parameter. For example, recomputing our model predictions for gives at GHz and at GHz, which are consistent at with the corresponding ACT and SPT constraints. Given that is within the error bar for WMAP9 Hinshaw et al. (2012), it is difficult to assess the extent to which our fiducial model may be discrepant with the ACT and SPT results. The difference can easily be explained by small changes in and is also sensitive to variations in the ICM physics, which we have kept fixed in these calculations. We conclude that our model is not in significant tension with current tSZ measurements (or other cosmological parameter constraints), and is thus a reasonable fiducial case around which to consider variations.

Figs. 1 and 2 show the tSZ power spectra for our fiducial model and several variations around it, including the individual contributions of the one- and two-halo terms. In the fiducial case, the two-halo term is essentially negligible for , as found by earlier studies Komatsu & Kitayama (1999), and it only overtakes the one-halo term at very low- (). However, for , the influence of the two-halo term is greatly enhanced due to the scale-dependent bias described in Section II.2, which leads to a characteristic upturn in the tSZ power spectrum at low-. In addition, induces an overall amplitude change in both the one- and two-halo terms due to its effect on the halo mass function described in Section II.1.1. While this amplitude change is degenerate with the effects of other parameters on the tSZ power spectrum (e.g., ), the low- upturn caused by the scale-dependent bias is a unique signature of primordial non-Gaussianity, which motivates our assessment of forecasts on using this observable later in the paper.

Fig. 2 shows the results of similar calculations for . In this case, the effect is simply an overall amplitude shift in the one- and two-halo terms, and hence the total tSZ power spectrum. The amplitude shift is caused by the change in the halo mass function described in Section II.1.2. Note that the sign of the amplitude change is somewhat counterintuitive, but arises due to our choice of as a fundamental parameter instead of , as mentioned above. In order to keep fixed while increasing , we must increase , which leads to an increase in the tSZ power spectrum amplitude. Although this effect is degenerate with that of and other parameters, the change in the tSZ power spectrum amplitude is rather large even for small neutrino masses (% for eV, which is larger than the amplitude change caused by ). This sensitivity suggests that the tSZ power spectrum may be a useful observable for constraints on the neutrino mass sum.

We demonstrate the physical effects of each parameter in our model on the tSZ power spectrum in Figs. 411, including the effects on both the one- and two-halo terms individually. Note that the limits on the vertical axis in each plot differ, so care must be taken in assessing the amplitude of the change caused by each parameter. Except for and , the figures show % variations in each of the parameters, which facilitates easier comparisons between their relative influences on the tSZ power spectrum. On large angular scales (), the most important parameters (neglecting and ) are , , and . On very large angular scales (), the effect of is highly significant, but its relative importance is difficult to assess, since the true value of may be unmeasurably small. Note, however, that is important over the entire range we consider, even if its true value is as small as eV. Comparison of Figs. 4 and 8 indicates that the amplitude change induced by eV (for fixed ) is actually slightly larger than that caused by a % change in around its fiducial value.

We now provide physical interpretations of the effects shown in Figs. 411:

• (4): The change to the halo mass function discussed in Section II.1.1 leads to an overall increase (decrease) in the amplitude of the one-halo term for (). This increase or decrease is essentially -independent, is also seen at in the two-halo term, and is % for . More significantly, the influence of the scale-dependent halo bias induced by is clearly seen in the dramatic increase of the two-halo term at low-. This increase is significant enough to be seen in the total power spectrum despite the typical smallness of the two-halo term relative to the one-halo term for a Gaussian cosmology.

• (4): The presence of massive neutrinos leads to a decrease in the number of galaxy clusters at late times, as discussed in Section II.1.2. This decrease would lead one to expect a corresponding decrease in the amplitude of the tSZ signal, but Fig. 4 shows an increase. This increase is a result of our choice of parameters — we hold constant while increasing , which means that we must simultaneously increase , the initial amplitude of scalar fluctuations. This increase in (for fixed ) leads to the increase in the tSZ power spectrum amplitude seen in Fig. 4. The effect appears to be essentially -independent, although it tapers off slightly at very high-.

• (6): Increasing (decreasing) the amount of baryons in the universe leads to a corresponding increase (decrease) in the amount of gas in galaxy clusters, and thus a straightforward overall amplitude shift in the tSZ power spectrum (which goes like ).

• (6): In principle, one would expect that changing should change the tSZ power spectrum, but it turns out to have very little effect, as pointed out in Komatsu & Seljak (2002), who argue that the effect of increasing (decreasing) on the halo mass function is cancelled in the tSZ power spectrum by the associated decrease (increase) in the comoving volume to a given redshift. We suspect that the small increase (decrease) seen in Fig. (6) when decreasing (increasing) is due to the fact that we hold constant when varying . Thus, is also held constant, and thus is decreased (increased) when is increased (decreased). This decrease (increase) in the baryon fraction leads to a corresponding decrease (increase) in the tSZ power spectrum amplitude, as discussed in the previous item. The slight -dependence of the variations may be due to the associated change in required to keep constant, which leads to a change in the angular diameter distance to each cluster, and hence a change in the angular scale associated with a given physical scale. Increasing (decreasing) requires increasing (decreasing) in order to leave unchanged, which decreases (increases) the distance to each cluster, shifting a given physical scale in the spectrum to lower (higher) multipoles. However, it is hard to completely disentangle all of the effects described here, and in any case the overall influence of is quite small.

• (8): An increase (decrease) in has several effects which all tend to decrease (increase) the amplitude of the tSZ power spectrum. First, is decreased (increased), which leads to fewer (more) halos, although this effect is compensated by the change in the comoving volume as described above. Second, for fixed , this decrease (increase) in leads to fewer (more) baryons in clusters, and thus less (more) tSZ power. Third, more (less) vacuum energy leads to more (less) suppression of late-time structure formation due to the decaying of gravitational potentials, and thus less (more) tSZ power. All of these effects combine coherently to produce the fairly large changes caused by seen in Fig. 8. The slight -dependence may be due to the associated change in required to keep and constant, similar (though in the opposite direction) to that discussed in the case above. Regardless, this effect is clearly subdominant to the amplitude shift caused by , which is only slightly smaller on large angular scales than that caused by (for a % change in either parameter).

• (8): Increasing (decreasing) leads to a significant overall increase (decrease) in the amplitude of the tSZ power spectrum, as has been known for many years (e.g. Komatsu & Seljak (2002)). The effect is essentially -independent and appears in both the one- and two-halo terms.

• (10): An increase (decrease) in leads to more (less) power in the primordial spectrum at wavenumbers above (below) the pivot, which we set at the WMAP value (no ). Since the halo mass function on cluster scales probes much smaller scales than the pivot (i.e., much higher wavenumbers ), an increase (decrease) in should lead to more (fewer) halos at late times. However, since we require to remain constant while increasing (decreasing) , we must decrease (increase) in order to compensate for the change in power on small scales. This is similar to the situation for described above. Thus, an increase (decrease) in actually leads to a small decrease (increase) in the tSZ power spectrum on most scales, at least for the one-halo term. The cross-over in the two-halo term is likely related to the pivot scale after it is weighted by the kernel in Eq. (21), but this is somewhat non-trivial to estimate. Regardless, the overall effect of on the tSZ power spectrum is quite small.

• (10): Since sets the overall normalization of the ICM pressure profile (or, equivalently, the zero-point of the relation), the tSZ power spectrum simply goes like .

• (11): Since sets the logarithmic slope of the ICM pressure profile at large radii (see Eq. (22)), it significantly influences the total integrated thermal energy of each cluster, and thus the large-angular-scale behavior of the tSZ power spectrum. An increase (decrease) in leads to a decrease (increase) in the pressure profile at large radii, and therefore a corresponding decrease (increase) in the tSZ power spectrum on angular scales corresponding to the cluster outskirts and beyond. On smaller angular scales, the effect should eventually vanish, since the pressure profile on small scales is determined by the other slope parameters in the pressure profile. This trend is indeed seen at high- in Fig. 11. Note that a % change in leads to a much larger change in the tSZ power spectrum at nearly all angular scales than a % change in , suggesting that simply determining the zero-point of the relation may not provide sufficient knowledge of the ICM physics to break the long-standing ICM-cosmology degeneracy in tSZ power spectrum measurements. It appears that constraints on the shape of the pressure profile itself will be necessary.

## Iv Experimental Considerations

In this section we estimate the noise in the measurement of the tSZ power spectrum. The first ingredient is instrumental noise. We describe it for the Planck experiment and for an experiment with the same specifications as the proposed PIXIE satellite Kogut et al. (2011). The second ingredient is foregrounds777To be precise we will consider both foregrounds, e.g. from our galaxy, and backgrounds, e.g. the CMB. On the other hand, in order to avoid repeating the cumbersome expression “foregrounds and backgrounds” we will collectively refer to all these contributions as foregrounds, sacrificing some semantic precision for the sake of an easier read.. We try to give a rather complete account of all these signals and study how they can be handled using multifrequency subtraction. Our final results are in Fig. 12. Because of the several frequency channels, Planck and to a much larger extent PIXIE can remove all foregrounds and have a sensitivity to the tSZ power spectrum mostly determined by instrumental noise.

### iv.1 Multifrequency Subtraction

We discuss and implement multifrequency subtraction888We are thankful to K. Smith for pointing us in this direction. along the lines of THEO (); CHT (). The main idea is to find a particular combination of frequency channels that minimize the variance of some desired signal, in our case the tSZ power spectrum. We hence start from

 ^aSZℓm=∑νiwiaℓm(νi)gνi, (31)

where refers to our estimator for the tSZ signal at 150 GHz (the conversion to a different frequency is straightforward), are the different frequency channels relevant for a given experiment, are the weights for each channel, are spherical harmonic coefficients of the total measured temperature anisotropies at each frequency and finally is the tSZ spectral function defined in Section III, allowing us to convert from Compton- to . We can decompose the total signal according to with enumerating all other contributions. We will assume that , i.e. different contributions are uncorrelated with each other. Dropping for the moment the and indices, the variance of is then found to be

 ⟨^aSZ^aSZ⟩=CSZ(∑νiwi)2+∑νiνjwiwj∑fCf(νi,νj)gνigνj, (32)

where is the tSZ power spectrum at 150 GHz as given in Eq. (16) and (again the index is implicit) is the cross-correlation at different frequencies of the of each foreground component (we will enumerate and describe these contributions shortly). To simplify the notation, in the following we will use

 C(νi,νj)=Cij≡∑fCf(νi,νj)gνigνj. (33)

We now want to minimize with the constraint that the weights describe a unit response to a tSZ signal, i.e., . This can be done using a Langrange multiplier and solving the system

 ∂i[⟨(^aSZ)2⟩+λ(∑iwi−1)]=∂λ[⟨(^aSZ)2⟩+λ(∑iwi−1)]=0. (34)

Because of the constraint , the term in is independent of (alternatively one can keep this term and see that it drops out at the end of the computation). Then the solution of the first equation can be written as

 wi=−λ(C−1)ijej=0, (35)

where is just a vector with all ones and is the inverse of in Eq. (33). This solution can then be plugged back into the constraint to give

 wi=(C−1)ijejek(C−1)klel, (36)

which is our final solution for the minimum-variance weights. From Eq. (32) we see that the total noise in each