ON ESCAPING A GALAXY CLUSTER IN AN ACCELERATING UNIVERSE

The discovery of the late-time acceleration of the universe is one of the most profound mysteries of physical cosmology. What is at stake with this discovery is the following: either our universe is composed of some exotic "dark energy" whose physics drives the dynamics of acceleration and/or our general relativistic theory of gravity must be extended or modified (Joyce et al. 2015, 2016; Koyama 2016).

Given its profound importance, cosmic acceleration is currently being investigated through a broad constellation of probes that mobilize a range of astrophysical objects and phenomena such as type Ia supernovae, baryon acoustic oscillations, weak gravitational lensing, and galaxy clusters (see Weinberg et al. 2013 for an excellent review of these and other approaches).

Galaxy clusters in particular are vital laboratories that allow us to sensitively probe the physics of large-scale structure formation and thereby constrain cosmological models of our universe (Kravtsov & Borgani 2012). The method most commonly applied is based on the abundance function of clusters, which evolves in shape and amplitude as a function of the cosmological parameters (Vikhlinin et al. 2009; Rozo et al. 2010). The abundance function as a cosmological probe depends not only on robust analytic theory, which is calibrated through simulations, but also on accurate masses as inferred from dynamical, weak-lensing, and X-ray methods.

An alternative way of constraining cosmology that does not depend directly on the abundance function of galaxy clusters was developed in Regoes & Geller (1989) and extended by Regoes (1996). Both of these papers constrain the cosmological matter density parameter through an analysis of the phase (velocity versus distance) space of galaxy clusters. More specifically, as demonstrated by Kaiser (1986), the infall pattern around galaxy clusters forms a trumpet-shaped profile in their phase spaces. This trumpet-shaped profile, also known as a phase-space "caustic," can be inferred from the line-of-sight velocity information. When compared to what is predicted by spherical infall models, an estimate of the matter density parameter of the universe may be inferred. However, the capacity for the amplitude and shape of the caustic to precisely constrain the matter density parameter, at least when considering both linear and nonlinear theory from the spherical infall model, has since been called into question (Diaferio & Geller 1997).

In particular, Diaferio & Geller (1997) demonstrated that the caustic profiles were related to the escape velocity profile of the cluster as mediated by a projection function. This projection function stems from the fact that we observe only the line-of-sight component of a galaxy's velocity and the true velocity vectors can be non-isotropic. Only after projection can the caustics be utilized to infer the mass profiles of galaxy clusters. While this method opened the path for a novel way of estimating the mass profiles of galaxy clusters, the capacity for the caustic to constrain cosmology directly has not been pursued.

In what follows, we argue that a cosmological acceleration term must be included in order to properly model the escape velocity profile of galaxy clusters as inferred from their phase spaces. The reason for this is that the escape velocity profile is often defined by integrating the density profile out to infinity via the integral form of the Poisson equation. However, current analytical expressions for the gravitational potential of galaxy clusters at "infinity" are not well defined. As such, if the potential is not properly normalized it will yield the wrong escape velocity profile.

A recent analysis by Miller et al. (2016) utilized three-dimensional phase spaces from the Millennium simulation (Springel et al. 2005) to demonstrate that the escape velocity edge (e.g., the three-dimensional caustic) can recover the true underlying escape velocity profile of galaxy clusters to high accuracy and precision. The true escape velocity profile is determined from the application of the Poisson equation to the measured cluster density profile. As we discuss below, this technique requires cosmological information pertaining to the acceleration of space.

In their analysis, Miller et al. (2016) also report that the Navarro–Frenk–White (NFW) model overestimates the escape velocity profiles of galaxy clusters by $\sim 10 \%$ (Navarro et al. 1996, 1997). This is because the NFW density profile overestimates the mass beyond the virial radius. Other analytic models of the density profile of dark matter halos, such as the Einasto model (Einasto 1965) and the Gamma model (Dehnen 1993), fare much better. We note that each of these analytical representations of the density profile can be constrained to be identical within ∼r₂₀₀ (the radius at which the average density drops to $200\times$ the critical density of the universe). Where these density profiles differ is in the outskirts, where both the Einasto and Gamma profiles are steeper than the NFW profile.

In this effort, we extend the work of Miller et al. (2016) to include projection effects and also to test the theory and the algorithm on real data. We use projected synthetic data from the light-cone data of Henriques et al. (2012) as well as archival data on 20 galaxy clusters with extensive spectroscopic data and weak-lensing mass profiles. We analyze the cosmology-dependent escape velocity profiles as inferred from their projected phase spaces and assess the viability of our analytic expectations. We note that only when including cosmological effects do we recover values for the velocity anisotropy parameter that are in agreement with ΛCDM simulations and with various other observational studies: Łokas et al. (2006), Benatov et al. (2006), Lemze et al. (2009), Wojtak & Łokas (2010), and Biviano et al. (2013).

The outline of our paper is as follows: in Section 2 we derive and review the theoretical expectations related to our observable: the escape velocity profile of galaxy clusters. In Section 3 we outline how this observable is inferred from the phase space of galaxy clusters. From thereon, in Section 4, we describe both the synthetic and non-synthetic projected data we utilized to probe our theoretical expectations. In Section 5 we describe projection effects and estimate the value of the velocity anisotropy parameter that we infer by assuming two different theoretical expectations of the escape velocity profile: one that includes the cosmological term and one that does not. In the following section, Section 6, we thoroughly assess likely sources of systematics affecting our analysis. We follow this analysis with a discussion and conclusion in Sections 7 and 8 respectively.

Except for the case of synthetic data in which the cosmological parameters are already defined (Springel et al. 2005), in what follows we assume a flat ΛCDM cosmology with ${{\rm{\Omega }}}_{M}=0.3$, ${{\rm{\Omega }}}_{{\rm{\Lambda }}}=1-{{\rm{\Omega }}}_{M}$, and ${H}_{0}=100\,h\,\mathrm{km}\,{{\rm{s}}}^{-1}\,{\mathrm{Mpc}}^{-1}$ with h = 0.7.

The theory of general relativity and its derivative cosmological models demonstrate that the dynamics of the matter–energy and the universe's expansion dynamics are dialectically entwined. As is often said, matter–energy tells spacetime how to curve and spacetime tells matter–energy how to move. For instance, in the case of galaxy cluster-sized halos, large-scale cosmological dynamics are expressed in the amplitude and shape of the halo mass function (Tinker et al. 2008). Qualitatively, what this tells us is that the dynamics of the galaxies in a gravitationally bound structure such as a galaxy cluster must also necessarily be affected by cosmology.

In this paper we focus on a particular observable, the escape velocity profile of clusters as inferred from their phase space, and test the ways in which analytical formulations of this observable must necessarily introduce a cosmological term in order to accurately describe the escape velocity edges of galaxy clusters. In particular, the velocity profile (v_esc) can be inferred analytically by characterizing the potential (ϕ) that a given test particle must escape from,

$$\begin{eqnarray}&&{v}_{\mathrm{esc}}^{2}=-2\phi .\end{eqnarray} \tag{ 1 }$$

However, as mentioned in Section 1, the cosmological effect on the escape velocity profile has to be included in this potential. We now derive the escape velocity profile that includes cosmology and which has been utilized and tested against both simulations of the ΛCDM universe (Behroozi et al. 2013; Miller et al. 2016) and extensions to general relativity such as Chameleon f(R) gravity (Stark et al. 2016). We then extend those derivations to include projection effects.

2.1. Escape Velocity Profile of a Galaxy Cluster in an Accelerating Universe

Nandra et al. (2012) demonstrated that in the weak-field approximation of general relativity and at sub-horizon scales, a massive particle will still feel a force from the accelerated expansion of space. Following Nandra et al. (2012), then, the effective acceleration experienced by a massive particle with zero angular momentum in the vicinity of a galaxy cluster with gravitational potential Ψ is given by

$$\begin{eqnarray}&&{\boldsymbol{\nabla }}{\rm{\Phi }}={\boldsymbol{\nabla }}{\rm{\Psi }}+{{qH}}^{2}r\hat{r}.\end{eqnarray} \tag{ 2 }$$

The effective potential Φ therefore takes into account both the curvature produced by a density field with potential Ψ and the curvature produced by the acceleration term qH²r. From a Newtonian perspective, then, this last term can be thought of as a repulsive force that opposes the inward pull of the cluster's mass distribution.

In the second term of Equation (2), q is the usual deceleration parameter given by $q\equiv -\tfrac{\ddot{a}a}{{\dot{a}}^{2}}$. We assume a flat universe (${{\rm{\Omega }}}_{k}=0$), and a dark-energy equation-of-state parameter $w=-1$. Also, given that we work in the low-redshift regime (such that ${{\rm{\Omega }}}_{\gamma }(z)\approx 0$) the deceleration parameter can be expressed in terms of the redshift-evolving matter density parameter (${{\rm{\Omega }}}_{M}(z)$) and the redshift-evolving Λ density parameter (${{\rm{\Omega }}}_{{\rm{\Lambda }}}(z)$), $q(z)=\tfrac{1}{2}{{\rm{\Omega }}}_{M}(z)-{{\rm{\Omega }}}_{{\rm{\Lambda }}}(z).$ For our chosen cosmology at the present epoch, we obtain $q(z=0)=-0.55.$ Lastly, the Hubble parameter (H) for this same cosmology is, as usual, $H(z)={H}_{0}E(z)\ ={H}_{0}\sqrt{(1-{{\rm{\Omega }}}_{M})+{{\rm{\Omega }}}_{M}{(1+z)}^{3}}$.

Having defined the cosmological quantities that compose Equation (2), we can integrate over the physical radius (r) to find the effective potential, and subsequently the escape velocity profile via Equation (1) with effective potential Φ. Integrating Equation (2),

$$\begin{eqnarray}&&{\int }_{r}^{{r}_{\mathrm{eq}}}d{\rm{\Phi }}={\int }_{r}^{{r}_{\mathrm{eq}}}d{\rm{\Psi }}+{{qH}}^{2}{\int }_{r}^{{r}_{\mathrm{eq}}}r^{\prime} {dr}^{\prime} .\end{eqnarray} \tag{ 3 }$$

Note that we are integrating out not to infinity but to a finite radius, r_eq, which is termed the "equivalence radius" in Behroozi et al. (2013). The reason for this is that the escape velocity at infinity is poorly defined. This ambiguity introduces a problem with the normalization of the potential that is used to calculate the escape velocity profile. In particular, as demonstrated in Miller et al. (2016), this offset ends up overestimating the potential by ∼20%. Following Behroozi et al. (2013), we define the equivalence radius to be the point at which the acceleration due to the gravitational potential of the cluster and the acceleration of the expanding universe are equivalent (${\boldsymbol{\nabla }}{\rm{\Phi }}=0$), which yields ${r}_{\mathrm{eq}}={\left(\tfrac{{GM}}{-{{qH}}^{2}}\right)}^{1/3},$ where G is the gravitational coupling constant and M is the mass of the cluster as inferred from a choice of Ψ via the Poisson equation, as detailed in Section 2.2 below. Now, integrating Equation (3) out to r_eq, we have

$$\begin{eqnarray}&&{\rm{\Phi }}(r)={\rm{\Psi }}(r)-{\rm{\Psi }}({r}_{\mathrm{eq}})+\displaystyle \frac{1}{2}{{qH}}^{2}({r}^{2}-{r}_{\mathrm{eq}}^{2})+{\rm{\Phi }}({r}_{\mathrm{eq}}).\end{eqnarray} \tag{ 4 }$$

Setting the boundary condition such that the escape velocity must necessarily be zero at the equivalence radius, $-2{\rm{\Phi }}({r}_{\mathrm{eq}})={v}_{\mathrm{esc}}^{2}({r}_{\mathrm{eq}})=0$, and using Equation (1), we find

$$\begin{eqnarray}&&{v}_{\mathrm{esc}}(r)=\sqrt{-2({\rm{\Psi }}(r)-{\rm{\Psi }}({r}_{\mathrm{eq}}))-{{qH}}^{2}({r}^{2}-{r}_{\mathrm{eq}}^{2})}.\end{eqnarray} \tag{ 5 }$$

This reproduces the result shown in Stark et al. (2016) and Miller et al. (2016). From now on we refer to this escape velocity profile as Einasto qH².

Note that Equation (5) is therefore normalized to yield an escape velocity of zero at the equivalence radius once a Ψ has been chosen. Note also that Equation (5) yields the escape velocity profile in an accelerating universe for any choice of gravitational potential Ψ.

2.2. Gravitational Potential

While it is common to describe the density profile of galaxy clusters with the NFW model of dark matter halos (Navarro et al. 1996, 1997), recent investigations have shown that the NFW potential–density pair overestimates the escape velocity profile by ∼10% within galaxy clusters (Miller et al. 2016). This is because, on average, the NFW profile tends to overestimate the mass in the outskirts of galaxy clusters (see also Diemer & Kravtsov 2015).

Moreover, as is also demonstrated by Miller et al. (2016), in contrast to the NFW model, once the cosmological term (Equation (5)) has been included, both the Gamma (Dehnen 1993) and Einasto (Einasto 1965) gravitational potential profiles can model the radial escape velocity profiles to better than 3% precision. In what follows, we utilize the Einasto profile. However, we emphasize that our analysis is not dependent on the Ψ used. That is, in so far as the gravitational potential profiles are derived from a true density–potential pair, all of our analysis will yield the same results.

Because nearly all published weak-lensing mass profiles utilize the NFW model, we need to map the NFW parameters to the equivalent parameters in the Einasto profile. The Einasto representation of dark matter halos (Einasto 1965) is a three-parameter model ($n,{\rho }_{0},{r}_{0}$) described by the following fitting formula for the density profile:

$$\begin{eqnarray}&&\rho (r)={\rho }_{0}\exp \left[-{\left(\displaystyle \frac{r}{{r}_{0}}\right)}^{1/n}\right].\end{eqnarray} \tag{ 6 }$$

From the density field described by Equation (6) we can derive the gravitational potential using the integral form of the Poisson equation. As demonstrated by Retana-Montenegro et al. (2012), this yields

$$\begin{eqnarray}&&{\rm{\Psi }}(r)=-\displaystyle \frac{{GM}}{r}\left[1-\displaystyle \frac{{\rm{\Gamma }}\left(3n,{\left(\tfrac{r}{{r}_{0}}\right)}^{1/n}\right)}{{\rm{\Gamma }}(3n)}+\displaystyle \frac{r}{{r}_{0}}\displaystyle \frac{{\rm{\Gamma }}\left(2n,{\left(\tfrac{r}{{r}_{0}}\right)}^{1/n}\right)}{{\rm{\Gamma }}(3n)}\right].\end{eqnarray} \tag{ 7 }$$

As in Miller et al. (2016), ${\rho }_{0}$ (or the mass term M) can be thought of as the normalization, r₀ is the scale radius, and n is the index. Following Retana-Montenegro et al. (2012), we use ${\rm{\Gamma }}(3n)=2{\rm{\Gamma }}(3n,{d}_{n})$ where the d_n term is solved for via asymptotic expansion. In particular, we use the d_n term cited therein expanded up to the fifth order. Lastly, ${\rm{\Gamma }}(a,x)$ denotes the upper incomplete gamma function, given as usual by ${\rm{\Gamma }}(a,x)={\int }_{x}^{\infty }{t}^{a-1}{e}^{-t}{dt}.$

As shown in Sereno et al. (2016) and references therein, the mapping between NFW and Einasto profiles is straightforward. In our case, we insist that the two profiles be nearly identical within r₂₀₀. We do this by solving for the Einasto parameters that match analytical NFW density profiles on a cluster-by-cluster basis. As noted by Sereno et al. (2016), recent weak-lensing analyses of stacked clusters cannot distinguish between the Einasto and NFW halo representations of the density profile within this range.

2.3. Comparing Theoretical Escape Velocity Profiles with and without the Cosmological Terms

In Figure 1 we show the resulting profiles for a single cluster-sized halo at z = 0. In particular we plot the escape velocity profiles both for the Einasto qH² profile (Equation (5)) with three different cosmologies (dashed, solid, and dashed–dotted lines) and for the Einasto profile without the cosmological term (dotted line): ${v}_{\mathrm{esc}}(r)=\sqrt{-2{\rm{\Psi }}(r)}.$

Zoom In Zoom Out Reset image size

As set by our aforementioned boundary condition, the Einasto qH² escape velocity profiles all reach an equivalence radius where ${v}_{\mathrm{esc}}(r)\to 0$. We also highlight the significant difference between the Einasto qH² model with ${{\rm{\Omega }}}_{M}=0.3$ and the Einasto model without the cosmological term. Lastly, note that as we increase ${{\rm{\Omega }}}_{M}$ our Einasto qH² profiles converge to the Einasto profile without the cosmological term.

Given that all other data indicate that we do not live in an Einstein–de Sitter universe, we should be able to detect the $\gtrsim 10 \%$ cosmology-dependent effects on the escape velocity profiles of galaxy clusters. We first test our theoretical expectations in N-body simulations in order to thoroughly assess systematic effects (including projection effects, mass scatter, and others, as explained in the subsequent sections) and then utilize archival redshift data and weak-lensing mass estimates of 20 galaxy clusters to test our expectations.

As is clear from our theoretical expectations, our observables are the projected escape velocity profile of galaxy clusters (${v}_{\mathrm{esc}}^{\mathrm{edge}}(r)$) and the observed weak-lensing mass profile. We use the latter in our analytic model of the three-dimensional escape edge and we require the cosmology and the projection function as described in Section 5.2.

We utilize the redshift information of galaxy clusters to generate phase spaces (v_los versus r space) from which we infer the escape velocity profile. More specifically, we do this by transforming the galaxy redshifts at a given angular separation from the cluster center (θ) to line-of-sight velocities (v_los) at the cluster's redshift (z_c) via

$$\begin{eqnarray}&&{v}_{\mathrm{los}}=c\displaystyle \frac{(z-{z}_{c})}{(1+{z}_{c})},\end{eqnarray} \tag{ 8 }$$

where c is the speed of light and z denotes the redshift of individual galaxies. To calculate the physical distance from the cluster's center (r) we calculate the angular diameter distance (d_A) and use the angular separation (θ),

$$\begin{eqnarray}&&r={d}_{A}(z)\theta =\left[\displaystyle \frac{1}{1+z}\displaystyle \frac{c}{{H}_{0}}{\int }_{0}^{z}\displaystyle \frac{{dz}^{\prime} }{E(z^{\prime} )}\right]\theta .\end{eqnarray} \tag{ 9 }$$

For each cluster phase space we then identify the edge in radial bins of 0.1 Mpc by finding galaxies with the top 10% of velocities. We follow the interloper removal prescriptions of Gifford et al. (2013, 2016). The latter tested this edge-detection technique on projected data in two different simulations with widely varying values for ${\sigma }_{8}=0.8$ versus ${\sigma }_{8}=0.9$, where ${\sigma }_{8}$ is the normalization of the matter power spectrum on scales of 8 Mpc. They conclude that the edge detection is independent of large variations in the line-of-sight interloper fraction. We also test the robustness of this edge-detection technique as discussed in section 6 below. There are various techniques used to define phase-space edges in the literature (Diaferio & Geller 1997; Lemze et al. 2009; Serra et al. 2011; Geller et al. 2013; Gifford et al. 2013; Rines et al. 2013; Miller et al. 2016). With N-body simulations Miller et al. (2016) have shown that the escape velocity edge (v_esc^edge) can be inferred with $\sim 5 \%$ accuracy. We discuss all of this more thoroughly in Section 6.

The weak-lensing mass profiles are taken from the literature (Sereno 2015). As noted in Section 2.2, these are provided as NFW profiles, and we convert them to Einasto density profiles on an individual basis. When the NFW concentration parameter is not provided, we use the mass concentration from Duffy et al. (2008). We include the weak-lensing mass errors as provided in the literature.

We first test our theoretical expectations with a sample of synthetic clusters from the light-cone data of Henriques et al. (2012). After having assessed relevant systematics (thoroughly described in Section 6 below) with the synthetic data, we conduct our analysis on 20 clusters with weak-lensing and redshift data. We briefly discuss these data sets below.

4.1. Synthetic Data

4.2. Archival Data: Weak-lensing Masses and Galaxy Redshifts

To summarize, with galaxy redshift information for each cluster we create a phase space with its line-of-sight velocities (v_los). From this phase space we infer the escape velocity edge (v_esc^edge), as detailed in Section 3. On the other hand, with mass profile measurements we generate an analytic escape velocity profile through Equation (5) after assuming a gravitational potential with the form of Equation (7) (i.e., the Einasto model of dark matter halos). We expect that if our theoretical expectations can reproduce the edge profile to high precision then the average ratio ${v}_{\mathrm{esc}}^{\mathrm{edge}}/{v}_{\mathrm{esc}}$ should yield unity. What Miller et al. (2016) has shown is that the Einasto model with the additional cosmological term can analytically reproduce velocity escape edges inferred from three-dimensional phase spaces to high precision. Therefore, when comparing our analytic formulation with edges inferred from projected phase spaces, any difference between this ratio and unity should arise from projection effects.

5.1. Projection Effects

As demonstrated by Diaferio & Geller (1997) and Diaferio (1999), the three-dimensional escape velocity profile (v_esc(r)) can be projected by a function of the velocity anisotropy parameter ($g(\beta )$),

$$\begin{eqnarray}&&{v}_{\mathrm{los}}(r)={v}_{\mathrm{esc}}(r)\times {(\sqrt{g(\beta (r))})}^{-1},\end{eqnarray} \tag{ 10 }$$

where $g(\beta (r))$ is given by

$$\begin{eqnarray}&&g(r)=\displaystyle \frac{3-2\beta (r)}{1-\beta (r)}.\end{eqnarray} \tag{ 11 }$$

The anisotropy parameter (β) is given by the ratio of the velocity dispersion in the tangential direction (${\sigma }_{t}$) to the velocity dispersion in the radial direction (${\sigma }_{r}$),

$$\begin{eqnarray}&&\beta (r)=1-\displaystyle \frac{{\sigma }_{t}^{2}}{{\sigma }_{r}^{2}},\end{eqnarray} \tag{ 12 }$$

where ${\sigma }_{t}^{2}=\tfrac{1}{2}({\sigma }_{\theta }^{2}+{\sigma }_{\phi }^{2})$ includes both azimuthal and latitudinal velocity dispersions. The limiting cases are radial infall ($\beta =1$), circular motion ($\beta =-\infty$), and isotropy ($\beta =0$).

Therefore, if we have the three-dimensional velocity information for each of our clusters we can calculate β and directly project our profiles. However, in practice this parameter is difficult to determine and as such we are left to infer what β would be, given an expected theoretical profile for a given cluster, and compare our result with simulations.

5.2. Inferring the Anisotropy Parameter

Given that our theoretical expectation matches the edges in three dimensions to high precision, we should be able to infer the anisotropy parameter by taking the ratio between the escape velocity edge and our theoretical profile via

$$\begin{eqnarray}&&\left\langle \displaystyle \frac{{v}_{\mathrm{esc}}^{\mathrm{edge}}(r)}{{v}_{\mathrm{esc}}(r)}\right\rangle \langle \sqrt{g(r)}\rangle =1.\end{eqnarray} \tag{ 13 }$$

The brackets in Equation (13) signify that the average is calculated over N = 20 clusters. That is, we calculate the ratio for each cluster and then take the average in each radial bin, with a separation of ${\rm{\Delta }}(r/{r}_{200})=0.1.$ Moreover, the averaged ratio is weighted by the error on the ratio of each individual cluster in a given radial bin. The details of our error budget are thoroughly discussed in Section 6 below.

Using Equation (13), for a given average ratio, we can find what average $g(\beta (r))$ is needed to make that ratio unity. Then by inverting Equation (11) we can find the anisotropy parameter β via

$$\begin{eqnarray}&&\beta (r)=\displaystyle \frac{3-g(r)}{2-g(r)}.\end{eqnarray} \tag{ 14 }$$

Lastly, the theoretical expectation (v_esc) can be either the escape velocity profile with the cosmological term (i.e., Einasto qH², Equation (5)) or without it (i.e., just the Einasto potential, Equation (7): ${v}_{\mathrm{esc}}(r)=\sqrt{-2{\rm{\Psi }}(r)}$). The differences between the inferred β parameters for these two analytic profiles are detailed in Sections 5.2.1 and 5.2.2, respectively, for the synthetic sample. For the archival data we compare only our inferred β assuming the Einasto qH² theory (see Section 5.2.3).

5.2.1. Synthetic Data and Theory with the Cosmological Term

Assuming v_esc(r) in Equation (13) is the Einasto qH² model (Equation (5)), we show the weighted average ratio in Figure 2 (black stars) for a single set of 20 clusters in our synthetic sample.

Zoom In Zoom Out Reset image size

With Equations (13) and (14) we then calculate the ${\chi }^{2}$ for $0.3\lt r/{r}_{200}\lt 1$ in order to infer the most likely average β for this set of 20 clusters. Note that we focus on the region $0.3\lt r/{r}_{200}\lt 1$ because simulation results have shown that the anisotropy parameter is on average constant across different redshifts within this radial range (Serra et al. 2011; Lemze et al. 2012; Munari et al. 2013). After inferring β from this method we project the average profile. The results for this sample are also shown in Figure 2 (black dots).

Moreover, with the ${\chi }^{2}$ we calculated, and assuming a Gaussian likelihood ${ \mathcal L }\propto \exp [-{\chi }^{2}/2]$, we can generate a likelihood plot of our inferred average β for each of our 10 sets of 20 synthetic clusters. The result is shown by the gray band in Figure 3. The band represents the 1σ error on the distribution of likelihoods for all 10 sets of averaged cluster ratios and their inferred β's. In the same figure we compare our inferred β with two other results from synthetic data.

Zoom In Zoom Out Reset image size

The red vertical dashed line is the average anisotropy parameter (also for $0.3\lt {r}_{200}\lt 1$) as directly measured in three dimensions using Equation (12) calculated by Iannuzzi & Dolag (2012). The $2\sigma$ bootstrap error on the mean is shown by red dotted lines. The sample used by Iannuzzi & Dolag (2012) is composed of the 1000 clusters in the Millennium simulation at z = 0 with virial masses greater than $2\times {10}^{14}\,{M}_{\odot }$. This is a much larger sample size than we use, and as such we consider this β to be a robust estimate of the true anisotropy parameter.

The black vertical dashed lines in Figure 3 come from a direct calculation of the anisotropy parameter (using Equation (12)) for the superset of 200 clusters with masses ${M}_{200}\gt 4\times {10}^{14}\,{M}_{\odot }$ that we use in this work. Our superset has a slightly larger average β, but our measurement and that of Iannuzzi & Dolag (2012) are within $\sim 2\sigma$ of each other. We hypothesize that this small difference could be attributed to the fact that our simulated data include the orphan galaxies in the catalog of Guo et al. (2011), whereas Iannuzzi & Dolag (2012) exclude these. Note, for instance, that the β calculated in Lemze et al. (2012) utilized particles and is also larger than that of Iannuzzi & Dolag (2012).

The gray line (and band) is the likelihood, which represents β using only the projected phase-space profiles of 20 clusters (i.e., no three-dimensional information). As expected, this likelihood is much larger than the error bounds on β from the three-dimensional information for larger samples. The reason it is larger is because the ${\chi }^{2}$ analysis includes representative errors on both the weak-lensing masses (25%) and the escape edges (5%). However, this likelihood fully captures the true underlying three-dimensional radially averaged velocity anisotropy.

5.2.2. Synthetic Data and Theory without the Cosmological Term

In both of our synthetic data samples shown in Figures 2 and 3 we took the ratio between escape velocity edge and the Einasto qH² analytic profile. Now, assuming we do not actually need the cosmology-dependent term to accurately reproduce the escape edge, we test whether or not we can recover the true anisotropy parameter. More specifically, the v_esc(r) in the ratio now utilizes Equation (1) with the Einasto potential of (7): ${v}_{\mathrm{esc}}(r)=\sqrt{-2{\rm{\Psi }}(r)}$ (rather than Equation (5)). The result is shown in Figure 3 (purple band).

As with the gray band described above, the purple band represents the likelihood on the inferred β with the $1\sigma$ error representing the standard deviation in the 10 sets of 20 synthetic clusters. We see that if we remove the cosmological term, the average anisotropy parameter is much larger and therefore we cannot recover the simulation results. More specifically, assuming the red lines are correct in Figure 3, i.e., that they truly describe the anisotropy parameter, we can rule out the non-cosmology-dependent theory (purple band) at the $6.3\sigma$ level.

In summary, these tests on the simulations provide two important results: first that our algorithm can recover the average β given only projected data and second that the cosmological model plays a significant role.

5.2.3. Archival Data and Theory with the Cosmological Term

After getting a sense of what we should expect for a given sample size of 20 clusters, and after making concrete the abstract necessity of including a cosmology-dependent term in our escape velocity profile by studying the underlying velocity anisotropy distributions with 200 synthetic clusters, we perform the same analysis on an archival data set of 20 galaxy clusters.

The ratio for our sample of 20 clusters, similar to Figure 2, is shown in Figure 4. As with the synthetic data, we calculate the ${\chi }^{2}$, and after assuming a Gaussian we can infer the likelihood for the average anisotropy parameter for our set of 20 clusters. The result is shown in Figure 5 where we compare the likelihood of β inferred for the archival set of 20 clusters (black line) to the distribution of $\beta$ inferred from the synthetic set of clusters (gray band, as in Figure 3).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We find that for our sample of just 20 clusters we recover the peak of the likelihood as inferred from simulations to high accuracy. In particular we find that the value of the most likely anisotropy parameter for our archival data on 20 galaxy clusters is $\beta ={0.248}_{-0.360}^{+0.164}$ at the 68% confidence level. Note that we calculate this interval by assuming that the distribution to the right of the peak is Gaussian and then find where $Q\gt -2\mathrm{ln}[{ \mathcal L }(\beta )/{{ \mathcal L }}_{\max }]$ (where Q = 1 for our single-parameter model). After this error is calculated then we find the error to the left of the peak by integrating the likelihood from the rightmost error up to the value of β that yields 0.68 times the total area of the likelihood.

In Figure 5 the data likelihood (black line) is larger than the simulation likelihood for β (gray band) and with a longer tail to negative β. This is because in the simulations we applied a representative error on the cluster masses fixed at 25%, whereas the data have errors that vary from 15% to 44%. The larger errors on the weak-lensing masses allow for more negative anisotropy. In fact, as we discuss below, the weak-lensing mass errors are the dominant component of our error budget.

The agreement between the anisotropy parameter inferred with the archival data of 20 clusters and the ΛCDM simulation results (assuming the Einasto qH² theory) is clear. To understand the robustness of this result, we carefully consider possible systematics.

As mentioned in the previous section, we carefully consider the uncertainties that make up our averaged ratio and we weight our data with those errors. In particular, we weight the ratio for a given cluster at a given radius by its error at that radius. More specifically, we propagate the error on the ratio by considering the error on the numerator (i.e., ${v}_{\mathrm{esc}}^{\mathrm{edge}}(r)$) and the error on the denominator (i.e., the theoretical escape velocity profile, v_esc(r)). We consider various uncertainties plaguing these two components of the average ratio below.

Throughout, we use a ${\chi }^{2}$ statistic and assume that the errors on the ratios are Gaussian. The fact that our likelihood is centered near the truth for the simulated data indicates that the assumption of Gaussianity is not a bad choice for our analysis.

6.1. Escape Velocity Edge

6.2. Weak-lensing Masses

6.3. Mass–Concentration Relation

Another component that is implicit in our calculation of v_esc(r) is the utilization of a mass–concentration relation to obtain the NFW density profiles. Recalling from previous sections, we utilize this NFW density to infer the Einasto parameters as described above. In particular, we utilize the mass–concentration relation of Duffy et al. (2008) for both the synthetic and archival data samples. Most importantly, this is the relation also used in the meta-catalog we utilize (Sereno 2015). The relation is given by

$$\begin{eqnarray}&&{c}_{200}{\left.({M}_{200},z)={A}_{200}{\left(\displaystyle \frac{{M}_{200}}{{M}_{\mathrm{piv}}}\right)}^{{B}_{200}}(1+z\right)}^{{C}_{200}},\end{eqnarray} \tag{ 15 }$$

where ${A}_{200}=5.71$, ${B}_{200}=-0.084$, ${C}_{200}=-0.47$, and ${M}_{\mathrm{piv}}=2\times {10}^{12}\,{h}^{-1}\,{M}_{\odot }$. We recalculated our inferred β from our samples by employing the $1\sigma$ error variations on the relation's parameters (${A}_{200},{B}_{200},{C}_{200}$) and found that the inferred β varied by only $\sim 1 \%$. Perhaps this is expected given the relative flatness of the mass–concentration relation at the high-mass end of the spectrum, which is where most of our clusters lie.

6.4. Cosmological Parameters

By utilizing archival data on just 20 galaxy clusters and by picking a fiducial cosmology within the range of what is expected from cosmological probes, we are able to recover the average anisotropy parameter in agreement with ΛCDM simulations. In this sense, then, we are already implicitly constraining cosmology by picking a reasonable choice of values for h, ${{\rm{\Omega }}}_{m}$, and w. What remains to be seen, however—and this much we defer to future work—is how precisely we can constrain cosmology with Equation (5) once β is independently inferred for each cluster, given the scatter on weak-lensing masses (our dominant systematic). If β can be inferred for each cluster through an independent technique (e.g., via the Jeans equation), we can leverage this to constrain cosmological parameters in the near future.

We note that our resulting average velocity anisotropies are well in agreement with other analyses (Benatov et al. 2006; Łokas et al. 2006; Lemze et al. 2009; Wojtak & Łokas 2010; Biviano et al. 2013). For example, with a sample of only six nearby relaxed Abell clusters, Łokas et al. (2006) constrains the anisotropy parameter to $-1.1\lt \beta \lt 0.5$ at the 95% confidence level. Our results with 20 clusters basically reproduce this constraint on $\beta .$

Furthermore, we note that our treatment would not affect observables that are either first or second derivatives of the potential, such as analyses of Jeans mass or inferences of X-ray masses. The reason for this is that the dominant term in Equation (5) is the normalization constant ${\rm{\Psi }}({r}_{\mathrm{eq}})$. The techniques that are affected by our theoretical expectations are those that deal directly with the escape velocity profile as such (or in other words, with the cluster's potential as inferred from dynamics). Therefore our analysis matters when the normalization of the potential matters. This is not the case for the two aforementioned methods nor for weak-lensing masses. In short, cosmology matters for masses inferred from escape velocity, as shown in a three-dimensional ΛCDM simulation (Miller et al. 2016), two-dimensional simulations (Figure 3), and real data (Figure 5, black line).

One technique that deals directly with the escape velocity profile as such is the caustic technique. Hoekstra et al. (2015) has recently demonstrated that mass estimates by the caustic technique (see Rines et al. 2013) are on average underestimated when compared to the weak-lensing masses by $\sim 22 \% .$ If we were to drop the cosmological terms in our theoretical profile, the overall escape velocity profile would increase (as in Figure 1). As such, in order to match this Einasto profile without the cosmological term, the Einasto qH² would always infer a higher mass. Interestingly, we find that M₂₀₀(Einasto ${{qH}}^{2})/{M}_{200}$(Einasto) = 1.22, exactly reproducing the result of Hoekstra et al. (2015). We note this as a possible explanation for this discrepancy and as a call to return to reflect upon the cosmological dependence of the escape velocity edges as already argued long ago by Kaiser (1986) and Regoes & Geller (1989).

With archival data on just 20 galaxy clusters with extensive redshift information and weak-lensing estimates of mass profiles we demonstrate the need to include a cosmology-dependent term in the analytic model of the escape velocity profile of galaxy clusters. We conduct our analysis also utilizing 10 sets of 20 synthetic galaxy clusters to study underlying systematics and projection effects. We find that our analytic formulation provides remarkable agreement with both sets of data.

More specifically, we leverage the complications involved in projecting the line-of-sight velocities related to the anisotropy parameter (β) and utilize this information to quantify the necessity of including a cosmological term in our analytic theory. We find that if we do not include a cosmological term in our analytic theory we infer velocity anisotropies that are inconsistent with both numerical results and observational data.

Throughout our analysis we picked a fiducial cosmology to probe our theoretical expectations and showed that there is a degeneracy between the velocity anisotropy parameter and cosmology. However, by independently inferring the anisotropy parameter, and combining this with the cosmology-dependent Einasto qH² theoretical profiles, one can in principle constrain cosmology.

As such, we have briefly advanced the capacity for the escape velocity profile of galaxy clusters to become a novel probe of cosmology in the near future given its sensitivity to the physics of cosmic acceleration. What is required to realize this is the next generation of weak-lensing data for galaxy clusters as well as deep spectroscopic follow-up (e.g., the Dark Energy Survey (Diehl et al. 2014) and the Dark Energy Spectroscopic Instrument—DESI). We defer a systematic of study of this kind to future work.

This material is based upon work supported by the National Science Foundation under Grant No. 1311820 and 1256260 and the Department of Energy grant DE-SC0013520. The Millennium Simulation databases used in this paper and the web application providing online access to them were constructed as part of the activities of the German Astrophysical Virtual Observatory (GAVO). This research has made use of the VizieR catalog access tool, CDS, Strasbourg, France. Lastly, we would like to thank the anonymous referee whose comments have improved both the quality of presentation and level of rigor of our analysis.