Cosmological Tests of Gravity
Pedro G. Ferreira


Pedro G. FerreiraCosmological Gravity










Cosmological observations are beginning to reach a level of precision that allows us to test some of the most fundamental assumptions in our working model of the Universe. One such as assumption is that gravity is governed by the General Theory of Relativity. In this review we discuss how one might go about extending General Relativity and how such extensions can be described in a unified way  on large scales. This allows us to describe the phenomenology of modified gravity in the growth and morphology of the large scale structure of the Universe. On smaller scales we explore the physics of gravitational screening and how it might manifest itself in galaxies, clusters and, more generally, in the cosmic web. We then analyse the current constraints from large scale structure and conclude by discussing the future prospects of the field in light of the plethora of surveys currently being planned. 

The key results are:


  * there are a plethora of alternative theories of gravity which are
  restricted by fundamental physics considerations; 

  * there is now a well established formalism for describing
 cosmological perturbations in the linear regime for general
 theories of gravity;

  * gravitational screening can mask modifications to general relativity
 on small scales but may, itself, lead to distinctive signatures in the
 large scale structure of the Universe;

  * current constraints on both linear and non-linear scales may be
 affected by systematic uncertainties which limit our ability to rule
 out alternatives to General Relativity;

  * the next generation of cosmological surveys will dramatically improve
 constraints on General Relativity, by up to two orders of magnitude.





General Relativity, Cosmology, Large Scale Structure, Cosmic Microwave Background, Early Universe

Cosmological Tests of Gravity
    Pedro G. Ferreira^1 
^1Astrophysics, University of Oxford,  Oxford OX1 3RH, UK; email: pedro.ferreira@physics.ox.ac.uk

    Received: 28 September 2018 / Accepted: 5 February 2019
===========================================================================================================================








§ WHY TEST GRAVITY?




 §.§ The success and the uniqueness of General Relativity.

The success of the standard model of fundamental physics is unquestionable and unassailable. There is a complete model for the strong and electro-weak forces (dubbed the standard model – SM) which has passed every single test to which it has been subjected. While some minor anomalies have cropped up, they are not sufficiently significant to overhaul the main picture: a gauge theory of interactions between fermions where electro-weak symmetry is spontaneously broken by a scalar field. The remaining force, gravity, is perfectly described by the general theory of relativity (GR), a theory of a dynamical space time where general covariance  (also called general coordinate invariance or diffeomorphism invariance) plays a fundamental role. Again, observations repeatedly confirm that GR works  in a number of regimes. Observations and experiments on the scale of the Solar System using lunar laser ranging, satellite missions such as Cassini or exquisite measurements of planetary orbits compete with measurements of individual and binary millisecond pulsars as well as triple systems to obtain exquisite constraints of gravity <cit.>.

With the advent of gravitational wave astronomy and black hole imaging, it will be possible to, for the first time, find constraints in the strong gravity regime.


When extrapolated to the largest, cosmological, scales, the SM and GR seem to accurately describe the history and evolution of the universe, as well as the large scale structure of space-time and matter. In that case, the SM and GR have been extrapolated to length scales and energy scales which are well beyond any current experimentally (or astrophysically) acessible regime. In other words,  two theories which pass tests in the laboratory and in the Solar System with flying colours are used to predict what is observed on length scales which are up to fifteen orders of magnitude greater. And the predictions match observations with remarkable consistency: the expansion rate of the universe, the morphology of large scale structure, the abundance of different types of particles, can all be consistently modelled with the current cosmological model. Furthermore, it allows us to predict that there is a dark sector to the Universe, i.e. that in the context of SM and GR, we should expect 95% of the energy density of the Universe to consist of some exotic form of material – dark matter and dark energy. The simplest formulation of this cosmological model has been called Λ- Cold Dark Matter (ΛCDM) where Λ represents the cosmological constant. It is now customary to express cosmological results in terms of constraints on ΛCDM, some of which can reach sub-percent accuracy.

If one delves into the structure of gravity, there are compelling reasons to believe that GR is the unique theory that can explain it. Assuming that gravity is mediated by a spin-2 carrier – a transverse, traceless field – there are a number of arguments that show that GR is the non-linear theory that can correctly explain its dynamics. In a series of lectures and papers <cit.> it was shown that, beginning with the massless, Fierz-Pauli action for a spin-2 field, h_αβ <cit.>, the requirement of general covariance  implies that one needs to consider non-linear corrections and that, when all of them are correctly taken into account, means that the fundamental action is unique and is that of GR.
An alternative, and more classical, viewpoint <cit.> shows that the only second-order, local gravitational field equations derivable from an action containing solely the 4D metric tensor (plus related tensors) are the
Einstein field equations with a cosmological constant. This is backed up by a similar conclusion from the geometrodynamic point of view <cit.>.



 §.§ Hints from the late and early Universe.

While the evidence for the SM and GR as an integral component of the current working model of the Universe is extremely strong, it is important to step back and ponder the way forward. Indeed, with these core theories we can make remarkable predictions on cosmological scales and, in particular, it leads to extremely strong evidence for a dark sector. The evidence is so compelling that the quest to understand this dark sector is the scientific driver of a number of large scale observational programmes. Nevertheless, given that there is no other (non-gravitational) evidence for the dark sector, it is a matter of common sense to question some of the fundamental assumptions that go into the evidence. And the main assumption is that GR is the underlying theory of gravity. An alternative point of view then, is that evidence for the dark sector may signal a break down of GR on cosmological scales, at late time. 

There is tantalizing evidence for deviation from general relativity if one attempts to place constraints on the physics of the very early universe. From observations of the anisotropy of the CMB, it is now possible to place extremely stringent constraints on the statistics of primordial fluctuations. The simplest assumption is that these fluctuations can be characterized in terms of an overall amplitude for both scalar and tensor components along with spectral indices for their spatial morphology. Current constraints on the spectral index place it close to, but not exactly at, scale invariance; they also place an upper-bound on the amplitude of tensor fluctuations (in the form of primordial gravitational waves) relative to scalar fluctuations (in the form of density perturbations) <cit.>. The most favoured theory that can explain these observations is inflation: a period of accelerated expansion in the early universe, amplifying microscopic quantum fluctuations into macroscopic cosmological perturbations.  Constraints from the CMB greatly restrict the range of models which are observationally viable. Remarkably, those that are favoured involve some modification to GR at very early times, either through a non-minimal coupling between the inflaton scalar field and the Ricci tensor or through the inclusion of higher order powers of the Ricci tensor in the fundamental action <cit.>.



 §.§ The (almost) strong field regime.

To some extent it is unsurprising that modifications should arise in a regime of high curvature. Purely from an effective field theory point of view <cit.>, one should expect corrections to GR to emerge in this regime – the Einstein-Hilbert action should have corrections proportional to R^2, R^μνR_μν and so one. Furthermore, in such high curvature regimes, any extra fields (or degrees of freedom) which are frozen at low energies will, most likely, be dynamical and contribute to any non-standard gravitational dynamics. It is then natural to look for modifications away from the weak field regime: while it is practically impossible to access high curvature regimes (apart from the early Universe or in the cores of black holes) one can be slightly less ambitious and look at the regime where the gravitational potential, Φ is of order unity. And, as luck would have it, two brand new windows on gravitational physics have opened up in this regime: gravitational wave detection and black hole imaging.

The advanced Laser Interferometric Gravitational Observatory (aLIGO) has observed a number of coalescing compact objects through their gravitational wave emission <cit.>. A number of merging black hole and neutron star binaries have been used to place stringent constraints on gravitational physics – thus far these constraints have focused on internal consistency tests in the waveform models where the signal is dominated by the weak field regime during the inspiral and after the merger <cit.>. But as the number of events build up, and analytic understanding of the merger improves, it should be possible to dig into the merger where the gravitational potential becomes appreciable. It should also be possible to extract information from the ring-down of the merger where signatures of deviation from GR may become evident in the quasi-normal mode spectrum. From this spectrum it will be possible, for example, to test the no-hair theorem or probe for the presence of extra degrees of freedom <cit.>. Alternatively, by studying mergers which are associated with electromagnetic counterparts, it is possible to place constraints on the propagation of gravitational waves and place constraints on the gravitational sector. Indeed, that has already been done to great effect <cit.> with the recently detected binary neutron star merger GW170817 <cit.>.

A coordinated effort by the consortia of telescopes that make the Event Horizon Telescope to observe Sgr A^* at the center of the Galaxy, will lead to the first image of a black hole event horizon (https://eventhorizontelescope.org) . By analyzing the structure of the accretion flow and how it wraps around the dark inner shadow, it will be possible to probe the structure of space-time close to the Schwarzschild radius and test whether it is consistent with GR, i.e. whether it is has a Kerr geometry. Preliminary measurements have already picked up structure on the scale of the event horizon <cit.> and, while one should expect the intricacies of the physics of accretion flows to complicate the analysis, it will be possible, for the first time, to look for evidence for deviations to GR close to (but not in) the strong gravity regime <cit.>.



 §.§ Modified gravity in cosmology: a roadmap.

Let us now turn our sights onto cosmological scales at late times. A first guess would be that deviations from GR on these scales would be unexpected. While the expanding Universe is, in some sense, a core result of GR or any relativistic form of gravity, and can be thought of as in the strong field regime, any fluctuations away from the homogeneous and isotropic space-time are extremely weak. Naively one would expect that it would be very difficult to construct modifications to GR in the regime where the curvature is small – in fact one might even think it problematic given that Minkowski space is the limit of zero curvature and we believe it to be an excellent approximation in a number of classical and quantum scenarios.

The past few decades have shown that naive expectations about the cosmological regime don't necessarily bear out. For a start, there has been an explosion of theoretical ideas leading to models that modify gravity in the infra-red while leaving it consistent with GR on astrophysical scales <cit.>. These have involved exploring the various loop-holes to Lovelock's theorem – higher dimensions and derivatives, extra fields and non-locality. The various proposals can all be rephrased in terms of extra fields which couple non-minimally to gravity and lead to new fifth forces.  Non-linear, classical effects can lead to screening of gravitational fifth forces in regimes of high gravitational potential and curvature while leaving them unscreened on large scales and in regimes of low curvature <cit.>. Constraints on fifth forces are remarkable <cit.> but are restricted to small scales (compared to cosmology). Hence, looking for modifications to general relativity on large scales is tantamount to looking for the effects of fifth forces in an unscreened regime. This is uncharted territory and one of the few pristine arenas where we might find evidence for new physics.

The search for dark energy has also become a search for modifications to general relativity on large scales. All the main observational programmes currently being planned have, as one part of their science cases, the goal of testing gravity using a range of different observables with a variety of tracers of large scale structure. The purpose of this review is to establish the core ideas underpinning such a programme. I will describe the basic ideas and tools that go into extracting meaningful constraints from cosmological data.

The structure of this review is as follows. In Section <ref> I briefly review current cosmology, describing the key aspects of gravitational physics that play a role in the evolution of the Universe. In Section <ref> I discuss how gravity can be modified on cosmological scales, with a particular emphasis on how this involves extra degrees of freedom and fields. In Section <ref> I work through the basics of gravitational screening and how it might mask modifications to gravity in certain regimes. In Section <ref> I talk through the various, existing constraints on gravity with large scale observables using existing surveys while in Section <ref> I focus on constraints involving the non-linear regime and environmental effects. Finally, in Section <ref>, I map out the future of the field.




§ MODERN COSMOLOGY.




 §.§ The expanding Universe.

 
The evolution of the Universe can be predicted from the field equations of GR, under the assumption that the metric of space-time, g_αβ is homogeneous and isotropic, g_αβ= diag(-1, a^2δ_ij) where a is the scale factor and only depends on time (for simplicity we will assume spatial flatness throughout). The dynamics of a is given by

    H^2≡(ȧ/a)^2=8π G/3ρ̅

where overdot denotes a derivative with regards to time and ρ̅ is the mean density of the Universe. The acceleration of the scale factor is given by

    ä/a=-4π G/3(ρ̅+3P̅)

where P̅ is the total mean pressure of the Universe. Equation <ref> can be rewritten in terms of the fractional energy densities of the various components, Ω_X=8π G ρ̅_X/(3H_0^2) (where X can be: radiation, R, baryons, B, cold dark matter, C, dark energy, DE) such that

    H^2(a)=H^2_0[Ω_ R/a^4+Ω_ B+Ω_ C/a^3+Ω_ DE/a^3(w+1)]

(where for simplicity we have assumed here a constant w). Clearly, a measurement of H as a function of time or scale factor, a, can be used to place constraints on the Ω_X and on the equation of state of dark energy, w=P̅_ DE/ρ̅_ DE.

Distant objects will have their light redshifted: ν_ emitted/ν_ observed=a_0/a≡ 1+z where the ratio of the emitted and observed frequencies depend on the scale factor at emission (a) and at the observer a_0, and where we have defined the redshift z. The relation between distance and redshift can be reconstructed through Equation <ref>. For example, for a Euclidean universe, the angular-diameter distance, D_A, (i.e. the distance reconstructed from knowing the physical size and measuring the angle it subtends in the sky) is given by

    D_A(z)=1/1+z∫_0^zdz/H(z)
 
and the corresponding luminosity distance is given by D_L=(1+z)^2D_A. By measuring redshifts and distances of stars, galaxies, supernovae and matching them to D_A(z) or D_L(z) at z close to 0 we can estimate H(z=0)≡ H_0; by measuring more distant objects, we can infer other cosmological parameters such as the Ω_X's and w.

The expansion of the Universe has an impact on matter and radiation. Throughout its history, the Universe has undergone a series of transitions as the overall temperature decreases. One of the most striking transitions occurs when the temperature is of a few thousand degrees Kelvin, corresponding to a few tenths of an electron-Volt. At that time, when the Universe was about 380,000 years old, at a redshift, z≃ 1100, free electrons and protons combined to form neutral Hydrogen. As a result, the ambient medium which was, until then, optically thick, became transparent; photons which last scattered with residual electrons and protons, remained free to propagate along null geodesics. As a result, we are constantly bombarded by a relic bath of radiation with black-body distribution with temperature, T̅=2.735^o Kelvin and an almost perfectly isotropic distribution. Indeed, it is useful to think of these relic photons as having been emitted  at  time t_* from the surface of sphere – the surface of last scatter – centred on us.




 §.§ Linear cosmological perturbation theory.

It is, of course, with the study of large scale structure that cosmology has really taken off. Going beyond a homogeneous and isotropic universe, one needs to include inhomogeneities. This can be done by perturbing the metric, g_αβ= diag[-(1+2Φ), a^2(1-2Ψ)δ_ij)] (in the Newtonian gauge), the various densities, ρ_X=ρ̅_X(1+δ_X) and pressures, P_X=P̅_X(1+δ_PX) and the temperature of the CMB, T=T̅(1+δ T/T̅); note that Φ, Ψ, δ_X and δ_PX are functions of time, t and position, x⃗ while δ T/T̅ is  function of t, x⃗ and direction n̂.

The evolution of the perturbations can be found through cosmological perturbation theory: taking the full equations for gravity (the Einstein field equations) along with the equations for conservation of mass and momentum (the relativistic Euler equation) and perturbing to linear order. So, for example, in Fourier space, a combination of the 00 and 0i component of the Einstein field equations yields something akin to the relativistic Newton-Poisson equation,

    -k^2Φ=4π G Δρ

with Δρ=∑_Xρ̅_X[δ_X+3(1+w_X)H/aθ_X/k^2] and θ_X=∇⃗·v⃗_X where v⃗_X is the peculiar velocity of component X and k labels the Fourier mode. And the traceless, transverse component (in the absence of anisotropic stress) of the Einstein field equation give

    k^2(Φ-Ψ)=0

The perturbed conservation equations (i.e. the linearized ∇_μ T^μν=0 equation, where ∇ is the covariant derivative and T^μν is the energy momentum tensor) can be used find additional evolution equations for δ_X and θ_X. With a suitable definition of the  sound speed, c^2_X one can replace δ_PX.



 §.§ The observables: power spectra, growth rate and cosmic shear.

When studying the perturbations around a homogeneous and isotropic background, statistics is paramount. First of all, and given that the (ensemble) average (or mean) of all linear perturbations is zero, variance plays a key role. Thus the relevant quantities that need to be considered are the power spectrum of density fluctuations, P(k), and the angular power spectrum of CMB fluctuations, C_ℓ, given by

    ⟨δ^*_M( k')δ_M( k)⟩   =   (2π)^3P(k)δ^(3)( k- k') 
    ⟨ a^*_ℓ' m'a_ℓ m⟩   =   C_ℓδ_ℓℓ'δ_mm'

where δ_M( k) is the Fourier transform of δ_B+δ_C, and a_ℓ m is the spherical harmonic transform of T( n) (or of any two-dimensional observable on a sphere). The various power spectra (and cross-spectra, i.e. the cross  correlation between different fields) are at the heart of the bulk of cosmological analysis.










Of particular importance is the rate of growth of cosmic structure. The linear density contrast can be factorized, δ_M( k,t)=D(t)δ( k,t_0) where t_0 is the physical time today. In the case of a matter dominated Universe (or a Universe with a cosmological constant), the fact that the time dependence of the evolution is scale-independent is strictly satisfied. One can then define two useful quantities, the growth rate,

    f≡dln D/dln a

and the "density weighted growth rate", fσ_8, where σ^2_8=⟨(δ M)^2/M̅^2⟩|_8 h^-1 Mpc, i.e. the mass-variance on 8 h^-1 Mpc scales. These quantities map out the time dependence of gravitational collapse and will play a key role in cosmological tests of gravity. A fair representation of the current observational status of the growth can be found in Figure <ref>: the shaded line shows the expected growth rate from CMB constraints (dominated by density and potential fluctuations at z≃1100) and projected forward assuming ΛCDM, which should be compared with individual and independent constraints (i.e. the various dots with errors bars) from galaxy surveys which are (necessarily) at late redshift. 



The anisotropies in the CMB contain information about the radiation density contrast, δ_γ, the baryon velocity, v⃗_ B and potential fluctuations, Φ at last scattering (known as the Sachs-Wolfe effect) as well as potential fluctuations from last scattering until now (known as the Integrated Sachs-Wolfe effect) through:

    δ T/T(n̂)=1/3δ_γ|_*-v⃗_ B·n̂|_*+Φ|_*-∫_t_*^t_0(Φ̇+Ψ̇)dt/a

where subscript * means evaluated at the time of last scattering, t_*.
 Furthermore, the CMB will be slightly polarized due to the finite thickness of the surface of last scattering and by the reionized electrons at later times. In the left hand panel of Figure <ref> we show the angular power spectra of the anisotropy in CMB intensity and polarization; one can see distinctive features such as a series of oscillations due to the sound waves in the baryon-photon plasma before recombination as well as the damping on small scales (i.e. large ℓ) that arises from the thickness of the surface of last scattering. 

The CMB will also be affected by the intervening structure, from the surface of last scatter until now. As the photons propagate through the perturbed space-time, they will be deflected and the structure of hot and cold spots will be distorted and deformed – the CMB will be lensed. One way the effect will manifest itself is be distorting the angular power spectrum of the CMB; specifically it will very slightly smooth out the peaks and troughs arising from the baryon-photon oscillations generated before recombination. Another way lensing affects the CMB is by modifying the statistics of δ T/T; what was originally a multivariate Gaussian will become non-Gaussian due to lensing and by measuring the 4-point function (or equivalently the trispectrum) and seeing how it deviates from the Gaussian value, it is possible to reconstruct the angular power spectrum of the lensing potential, C^ϕϕ_ℓ.
In the right hand panel of Figure <ref> we show the angular power spectrum of the reconstructed lensing potential, i.e. of C^ϕϕ_ℓ.





Distant galaxies will also be lensed by intervening structure. Measurements of distortions of the shapes of galaxies, i.e. changes in the projected ellipticities, and fluctuations in magnitudes, can be directly related to the distortion tensor

    ( [ 1-κ-γ_1    -γ_2;    -γ_2 1-κ+γ_1 ])= δ_ij+∫_0^χdχ'(Φ+Ψ)_ij[x⃗(χ')]χ'(1-χ'/χ)

where the γ_i are the shears and κ is the convergence, χ=∫_t_0^tdt'/a(t') is the conformal distance and (i,j) are directions perpendicular to the direction of the light ray. One can attempt to reconstruct the gravitational potential probed by the incoming light rays by inverting equation <ref> and, for example, relating it to underlying mass distribution through the Poisson equation. This method has been successfully used for reconstructing the mass distribution of clusters of galaxies as well as for mapping out the large scale distribution of dark matter. Alternatively, one can study the auto and cross correlations of κ and γ_i as well as cross correlations with other measures of the galaxy density fluctuations.



 §.§ Small scales and non-gravitational physics.

While much of what has been discussed is primarily applicable to the largest scales, where linear perturbation theory is applicable, the observables can also be used on smaller scales. Furthermore, we expect the smallest scales to have a larger statistical weight – there are more modes contributing to any particular quantity. But while the statistical weight will improve, the underlying physics becomes more complex as non-linear gravitational collapse plays a role, allied with non gravitational physics. Baryons, through the role of gas and stars or highly energetic processes like supernovae and active galactic nuclei can play a significant role in moulding structure on scales of Megaparsecs and below <cit.>. 

There has been considerable progress in developing both analytic techniques (through different versions of perturbation theory <cit.> and effective field theory <cit.>) and numerical algorithms (with N-body simulations, combined N-body and hydrodynamic simulations and, more recently, machine learning algorithms). In this, more complex regime, one can look at a more varied set of statistics, above and beyond the variances described above. For example, three point and four point statistics (i.e. which are cubic or quartic in the perturbed fields) will have additional information about the physical processes at play.

A different point of view is to look at what seem to be the building blocks of large scale structures, i.e. galaxies and clusters of galaxies or, from the CDM paradigm, halos of dark matter. This more "granular" view of the distribution of matter can be used to constrain a self consistent model – the halo model – which quantifies the structure of the universe in terms of the density profile of individual halos, ρ(r), the number density of halos of a give mass, n(M) and the linear powerspectrum of fluctuations on the largest scales <cit.>.

While there has been tremendous progress in developing efficient quantitative methods for the non-linear regime, and they, arguably, may be of sufficient accuracy for the analysis of current data sets, much work still needs to be done to develop a fully accurate and general approach to this regime.



 §.§ ΛCDM

Putting all the pieces together, it is possible to come up with a consistent picture of what the Universe looks like on large scales. It is accurately described by GR with a cosmological constant (Λ) making up roughly 70% of the energy budget, pressureless cold dark matter (CDM) making up roughly 25% of the energy budget and baryons making up 5%. The evidence for this model has accumulated over two decades. The first compelling evidence came from measurements of the luminosity distance through observations of distant supernovae IA <cit.>. Crucially, a few years later, measurements of the angular power-spectrum of the CMB from ground-based and balloon borne experiments found the first constraints on the geometry of the Universe <cit.>, ruling out curved models and singling out a flat Universe with a cosmological constant. Constraints from the abundance of clusters <cit.> placed constraints on the fractional energy density in dark matter, Ω_M, further cementing ΛCDM as the preferred model for describing the Universe. Subsequent cosmological observations, which we will allude to in Section <ref> have led to precision constraints (in some case at the sub-percent level) on the parameters of ΛCDM.




§ MODIFIED GRAVITY AND COSMOLOGY.




 §.§ What is General Relativity?

Before we delve into what it means to modify gravity, let us establish what we mean by standard gravity. We take it to be Einstein's General Theory of Relativity, where the gravitational field is a metric, g_αβ whose dynamics is described by the Einstein-Hilbert action, and which is minimally coupled to matter. In other words, gravity is encapsulated in the action

    1/16π G∫ d^4x√(-g)(R-2Λ)+∫ d^4x√(-g) L_M(g_αβ,⋯)

where G is Newton's constant, g is the determinant of g_αβ, R is the Ricci scalar, Λ is the cosmological constant and L_M is the Lagrangian density for the matter and radiation fields which are minimally coupled to the metric. It is often convenient to work with the reduced Planck mass, M^2_ Pl=1/8π G. The equations of motions are then given by

    G_μν≡ R_μν-1/2Rg_μν=-Λ g_μν+8π G T_μν.

where we have defined the Einstein tensor, G_μν.

GR is an incredibly special theory, as we hope we made clear in Section <ref>. Before we move venture into the unknown, it is useful to signal potential pitfalls and problems that arise when one tries to consider extensions.  To begin with, GR is endowed with properties that make it a predictive theory. The equations are highly non-linear and, a-priori, there is no reason to assume that it allows a well-posed initial value problem. In other words, given a set of initial conditions, is there a unique and well-defined evolution for the metric, given by equations <ref>? It turns out that there is <cit.> and this means that it is possible to trust numerical solutions to the evolution equations, even in the mildly strong field regime such as, for example, in binary black hole mergers <cit.>. The same cannot be said for the types of extensions we will see below <cit.>. 

Furthermore, as we shall see, extensions to GR almost universally involve considering new degrees of freedom which, quite easily, may be unstable in one form or another. If one considers equations of motions that have higher derivatives, these may suffer from the Ostragradski instability <cit.>; they may also have ghosts or tachyonic instabilities that may render these theories phenomenologically unviable. Any attempt at venturing beyond GR must bear these serious problems in mind.



 §.§ Modified gravity as extra degrees of freedom.

In the introduction, we stated Lovelock's theorem and we mentioned that it was a useful guide to how we might deform gravity away from Equation <ref>. Let us now use it to explore a (by no means exhaustive) set of alternatives. For a start, let us consider adding an extra, scalar degree of freedom. One of the most iconic alternatives to GR is Jordan-Brans-Dicke (JBD) gravity <cit.> in which the Einstein-Hilbert term is transformed as M^2_ Pl/2R→φ R and the scalar field is endowed with dynamics with a kinetic energy of the form ω(∂_μφ∂^μφ)/φ. In practice this involves making the Planck mass (or Newton's constant) dynamical, with one free parameter, ω. JBD gravity is a particular example of a scalar-tensor theory for which there is an (almost) complete description if one requires (as in Lovelock's theorem) that the equations of motion are 2^ nd order. Known as the Horndeski (or Covariant Galileon) action <cit.>, it is given by:

    S=∫ d^4x √(-g){∑_i=2^5 L_i[φ,g_μν]+ L_M[g_μν,⋯]},

and the building blocks of the scalar field lagrangian are

    L_2   =    K ,  
     L_3   =     -G_3 φ , 
     L_4   =      G_4R+G_4X{(φ)^2-∇_μ∇_νφ∇^μ∇^νφ}  , 
     L_5   =    G_5G_μν∇^μ∇^νφ
    -1/6G_5X{ (∇φ)^3
    -3∇^μ∇^νφ∇_μ∇_νφφ
           
    +2∇^ν∇_μφ∇^α∇_νφ∇^μ∇_αφ} .

We have that K and G_A are functions of φ and X≡-∇^νφ∇_νφ/2, and the subscripts, X and φ, denote derivatives. The four functions, K and G_A completely characterize this class of theories. Further extensions have been found that preserve the 2^ nd order nature of the equations of motion <cit.>. 

Of particular interest are theories in which additional symmetries are imposed and which will, through naturalness arguments, be protected at the quantum level. Of particular note are Galileon theories <cit.> where the action is invariant under field redefinitions of the form φ→φ+c_μ x^μ +d. Other possibilities are scale invariant theories <cit.>, invariant under joint field/metric transformations of the form φ→λ^-1φ and g_αβ→λ^2 g_αβ.

There are also other possible choices for extra fields. So, for example, one can opt for a four-vector, A_μ, leading to a different class of theories. Of particular note are Einstein-Aether theories <cit.> (and their generalizations –<cit.>) in which one imposes a time-like constraint, A_μ A^μ+1=0; such theories are used as models for spontaneously broken Lorentz symmetry. Other possibilities are extensions of the Maxwell-Proca action, in which the vector field is given an explicit mass <cit.>. 
It is possible to go yet further and include another tensor field, i.e. another metric f_αβ, as the extra degree of freedom. Such theories have garnered a fair amount of attention recently: by constructing a potential, V(g^-1f) one can obtain a theory of a massive graviton which isn't plagued by instabilities, most notably ghosts <cit.>. There are, of course, a myriad of possibilities, adding multiplets of fields or a combination of fields with different spins <cit.>.



 §.§ Alternatives: higher derivatives, higher dimensions and non-locality.

Instead of explicitly adding extra fields, it is possible to do so implicitly in a variety of different ways. One of the simplest ways is to replace R in the Einstein-Hilbert action by a function of other curvature invariants, R→ f(R,R^2, R_αβR^αβ, R_αβμνR^αβμν,⋯) where R_αβ and R_αβμν are the Ricci and Riemann tensors. In doing so, one introduces higher derivative terms in the action which is tantamount to adding new solutions to the equation of motion or, equivalently, adding new degrees of freedom. This can be made clear in the particular case where the f depends solely on the Ricci scalar <cit.>; by defining a new scalar field, φ∼ln df/dR, the theory can be re-expressed as scalar-tensor theory, a subset of the theories discussed above.

One approach which emerges from attempts at unification, such as string theory, is to enlarge the number of dimensions. A notable example is the Kaluza-Klein model <cit.> in which there is one extra dimension and the dynamics is given by the 5-dimensional Einstein-Hilbert action for the metric g_AB where A, B take values from 0 to 4. By redefining the metric as g_AB=e^φ/√(3)g_μν+e^-φ/√(3)A_μ A_ν when A, B take values from 0 to 3, g_μ4=e^-φ/√(3)A_μ and g_44=e^-φ/√(3),  it is possible to trivially equate a higher dimensional theory with a lower dimensional theory with extra fields (albeit in which the fields have an extra coordinate dependence on x^4). By compactifying the extra dimension or by confining the dynamics to a 3-space <cit.>, it is possible to obtain an approximate, but not exactly, 3+1 universe with modified gravitational dynamics (given by extra fields). A particularly popular model – the DGP model – was proposed in <cit.> of a three dimensional "brane" living in a five dimensional space. Through a judicious choice of the brane geometry and tension, vis-a-vis the background space-time it is possible to obtain a "self"-accelerating Universe which mimics one with a cosmological constant (or with dark energy). Furthermore, deviations from GR emerge on large scales (and not on small scales as one expects in "conventional" higher dimensional universes).



The final, more exotic, direction in which we can take this is by considering non-locality <cit.>. It is possible to, for example, include terms in the action of the form R× R/ where the inverse D'Alembertian can now be interpreted as a Greens function (which is by definition non-local). Interesting effects arise in a cosmological context due to the non-locality in time and it can be shown that such theories are equivalent to  scalar tensor theories by replacing φ≡ R and imposing strict constraints on the initial conditions in such a way as to kill all spurious modes arising from the homogenous solutions to this equation. 



 §.§ The speed of gravitational waves and GW170817


Before I delve into the cosmological regime of gravitational theories, it is important to briefly highlight a recent event. The discovery of a binary neutron star (BNS) merger with both gravitational waves and a multitude of electromagnetic signals has heralded the advent of multi-messenger astronomy <cit.>. One of the remarkable bi-products has been the ability to place stringent constraints on the speed of gravitational waves: it is within 10^-15 of the speed of electromagnetic waves. This fact allows us to severely restrict the range of allowed modified gravity theories that have an impact on cosmological scales <cit.> – the extra degrees of freedom will serve as a medium through which the gravitational waves propagate and, depending on their coupling to gravity, may affect their speed. It has been shown that, for example, scalar-tensor theories are strongly constrained by the multi-messenger detection of the BNS, as is a broad swathe of vector-tensor theories.



 §.§ General theories of gravity on cosmological scales.

I briefly covered a variety of ways in which one can extend GR and we have seen that, in general, it can be seen as adding extra degrees of freedom. But we have also seen that the landscape of possible extensions is vast and resides in large functional space of possibilities. This will give pause for thought: given that there is a restricted amount of observational data, the task of comprehensively constraining GR seems difficult at best. The key to progress is to realize that, in this endeavour, one should focus on specific regimes in which certain characteristics of the theory may simplify. Indeed, this is standard practice in high-energy physics, when one focuses on the relevant aspects of a theory at a given energy scale (for example, at the beam luminosity of a particular accelerator). Furthermore, this approach is already standard practice in, for example, the analysis of gravitational data on the scale of the Solar System. There, an efficient formalism for tackling weak field, non-relativistic systems – the Parametrized Post Newtonian (PPN) formalism – has been used to characterize (non)-deviations from GR of various forms <cit.>. It would make sense, then, to apply the same rationale to cosmology.

At its most basic level, cosmology is about the expansion of the Universe, i.e. a detailed description of the time evolution of the scale factor, a(t). As we have seen in GR, this allows us to disentangle the different constituents of the Universe, by constraining the Ω_X in the Friedman equation. Of particular interest is the fact that we may be able to constrain the equation of state of the dark energy component, w_ DE. In other words, observables such as H(z) or D_A(z) can be converted into constraints on w_ DE(z) which in turn can be used to extract the fundamental nature of the dark energy (i.e. if it is a cosmological constant, if it is a scalar field with a particular form of the potential, V(φ), etc).

How does the situation change if one now wishes to constrain gravity itself? Unfortunately, not very much. To see why this is so, we need to realize that any modification to GR can be written in the form

    G_μν=-Λ g_μν+8π G T_μν+U_αβ

where U_αβ can be a functional of g_αβ, the matter fields and any extra degree of freedom (such as scalar, vector and tensor fields). If we now impose homogeneity and isotropy, we have that U_αβ= diag(8π G ρ_U,8π G P_Uδ_ij) where  (ρ_U,P_U) are two free functions of time; note that we have chosen a notation that shows that, at the background level,  modifying GR is essentially equivalent to adding a new energy-momentum tensor, which can also be characterized by an equation of state, w_ U. This means that, at the background level, modifying GR is indistinguishable from adding in a dark energy fluid <cit.> and it is impossible to distinguish between the two with observables such as H(z) or D_A(z).

Given the paucity of information in the background, we need to delve a level deeper in complexity and look at linear perturbations of homogeneous and isotropic background space-time. As we saw, in GR, this leads to a series of  2^ nd order coupled differential equation (as well as constraints) for the perturbed variables δ_X, θ_X, Φ and Ψ. If we now go beyond GR, but remain in the arena of linear perturbations, we can easily guess what the structure of these "extended" equations might be <cit.>. They will be linear and (for now) 2^ nd order but they will include new, extra degrees of freedom. So, for example, if the extra degree of freedom arises from a scalar field, φ, we will now need to include its perturbation, δφ.

Let us now take two of the perturbed Einstein field equations and generalize them. We can modify Equation <ref> to

    -k^2Φ=4π G Δρ+A_0k^2Φ+F_0k^2+α_0k^2δφ+α_1kδφ̇+⋯Γ

and Equation <ref> to

    k^2(Φ-Ψ)=D_0k^2Φ+D_1kΦ̇+K_0k^2Γ+K_1kΓ̇+ϵ_0k^2δφ+ϵ_1k δφ̇+ϵ_2 δφ̈

where kΓ=(Φ̇+ HΨ) and A_0, D_0, D_1, F_0, K_0, K_1, α_0, α_1, ϵ_0, ϵ_1 and ϵ_2 are free functions of time and k; setting these free functions to zero recovers GR but, in principle, one can span the most general class of gravitational theories which involve an extra scalar field. The generality of this approach is appealing yet one is clearly faced with a challenge: is it feasible to, with a finite amount of cosmological data, constrain these multiple functions of time and space?

We have been working at the level of the equations of motion but, in almost full generality, alternative theories to GR will obey an action principle. This means that it would make more sense to extend the action itself. We have already seen that at the full, non-linear level, with the Horndeski action presented in Equation <ref>. But if one limits oneself to the most general theory at the linear level, the exercise simplifies greatly: one simply has to construct the most general action, quadratic in the perturbation variables and which satisfies any desired symmetry principle, the most important of which is linearized diffeomorphism invariance, i.e. invariance under transformations of the type x^α→ x^α+ξ^α. 

There are a few ways in which one can go about constructing the extend linear action. A systematic way is to write the most general quadratic action and then impose symmetries to reduce the number of terms <cit.>. This has been shown to be general and applicable to any type of extra degree of freedom: scalars, vectors and tensors. If one restricts oneself to the case where the extra degree of freedom is a scalar field, a simpler, more efficient approach can be taken with what has been dubbed the Effective Field Theory (EFT) approach <cit.>: starting with a general action which only depends on metric perturbations, one then performs a linearized diffeomorphism transformation, in which ξ^0=δφ and, in such a way one generates the most general linearized scalar-tensor action.

It turns out that, in the case of scalar-tensor theories, there is a straightforward approach of obtaining the most general action for linear perturbations: one takes the Horndeski action of Equation <ref> and Taylor expands to 2^ nd order in the perturbations. If one does that, one finds that, instead of depending on a multitude of free functions of time and scale, as we saw above, the resulting action merely depends on four functions of time, α_A(t) with A=K,B,M,T which can be directly related to the original action <cit.>. These functions are

    M^2_*   ≡   2(G_4-2XG_4X+XG_5φ-φ̇HXG_5X) , 
    
    HM^2_*α_M   ≡   d/dtM^2_* , 
    
    H^2M^2_*α_K   ≡   2X(K_X+2XK_XX-2G_3φ-2XG_3φ X) 
          
    +12φ̇XH(G_3X+XG_3XX-3G_4φ X-2XG_4φ XX) 
          
    +12XH^2(G_4X+8XG_4XX+4X^2G_4XXX)
          
    -12XH^2(G_5X+5XG_5φ X+2X^2G_5φ XX)
          
    +14φ̇H^3(3G_5X+7XG_5XX+2X^2G_5XXX) , 
    
    HM^2_*α_B   ≡   2φ̇(XG_3X-G_4φ-2XG_4φ X) 
          
    +8XH(G_4X+2XG_4XX-G_5φ-XG_5φ X)  
          
    +2φ̇XH^2(3G_5X+2XG_5XX)  , 
     
    M^2_*α_T   ≡   2X[2G_4X-2G_5φ-(φ̈-φ̇H)G_5X]  .

Each of these functions are linked to specific physical properties of the theory: M^2_* and α_M are related to time variations in the background Newton's constant, α_K to the generalized canonical kinetic term of simple DE models, α_B quantifies kinetic mixing between φ and the scalar perturbations of the metric and α_T is associated to modifications to the speed of propagation of tensor modes. 

The fact that general theories of gravity can, at the level of linearized cosmological perturbations, be completely characterized by a handful of free functions is true in general, including, for example, the case where the extra degrees of freedom are 4-vectors or tensors. And it gives hope to the idea that it might be feasible to constrain gravity on large scales. 


 §.§ The quasi-static regime.

It turns out that there is a further simplification if one restricts oneself to sufficiently sub-horizon scales, i.e. scales in which k/(aH)∼ kη≫ 1 where adη=dt allows us to define the conformal time η. Known as the "quasi-static" regime (the cosmological equivalent of the Newtonian limit), this covers scales up to hundreds of Megaparsecs and thus encompasses most of the current and future large scale structure surveys under consideration. In that regime, it is an excellent approximation to effectively freeze out the dynamics of the extra degree of freedom; for example, in the case of a scalar field, this means approximating φ∝ "source" by -∇^2 φ∝ "source". Replacing the extra degree of freedom in the generalized Einstein field equations leads to modified equations <ref> and <ref>:

    -k^2Φ   =   4π G μΔρ
    Φ   =   γΨ

where μ is the modified, effective Newton's constant (and can be rewritten as μ=G_ eff/G_0) and γ is often called the gravitational slip <cit.>. The quasi-static parameters μ and γ are generally functions of time and scale, k; typically the scale dependence either kicks in on very small scales (i.e. deeply in the non-linear regime) or naturally, close to the cosmological horizon when the quasi-static approximation breaks down. This means, that, effectively, when constraining theories of gravity on cosmological scales, one is constraining two free functions of time. An approximate approach is to extend this parametrization to all scales, including super-horizon scales, greatly simplifying any attempts at cosmological parameter estimation.

Note that, it is customary to assume that the quasi-static regime is synonymous with the use of μ and γ on subhorizon scales – we will do so in this review. But we need to point out that this is not entirely general and that in Beyond Horndeski theories <cit.>, there will be additional terms in equations <ref> and <ref> proportional to the matter velocity which also have to be taken into account; this means one would need to extend the number of functions in the quasi-static regime from two to four. 

Over the past few years, quite a few variants of μ and γ have been proposed, it is worth our while to briefly cover them. Some authors choose the following notation 

    -k^2Ψ   =   4π G μΔρ
    Φ   =   ηΨ

If we redefine the μ→μ̃ we have μ̃=μ/γ (where μ is defined in Equation <ref>) and η=γ. μ̃ is of particular use as (we shall see) it sources the growth rate of structure. In the same way, we can define a parameter that sources weak lensing

    Σ=μ(1+1/γ)=μ̃(1+η)

Some authors often define G_M and G_L (where M stands for "Matter" and L stands for "Light") such that

    G_M=μ̃   G_L=Σ.


The quasi-static parameters can be directly connected to more fundamental parameters in underlying theories. For example, in the case of scalar-tensor theories arising from the Horndeski action, these parameters have a specific structure <cit.>

    μ   =   1+h_1k^2/h_2+h_3k^2
    γ   =   h_2+h_3k^2/h_4+h_5k^2

where h_i, i=1,⋯ 5 are functions of time and assembled from the 4 α_X parameters from equations <ref>; thus constraints on the quasi-static parameters lead to constraints on slices of the functional space of free parameters, α_X. And vice versa, by restricting one self to scalar tensor theories, one can establish a restricted set of priors on μ, γ or Σ <cit.>



 §.§ The phenomenology of modified gravity.

The dynamics of linear perturbations in modified versions of GR is sufficiently well understood that there are now codes that accurately model the evolution of perturbations in such theories, calculating cosmological observables. The basic framework is a code that integrates the perturbed Einstein field equations, the coupled set of matter (and other fluid like) conservation equations and the Boltzman equations for radiation and neutrinos. The main, general purpose, Einstein-Boltzman solvers are EFTCAMB <cit.> and MGCAMB <cit.> and ISITGR <cit.>  (derived from CAMB <cit.>) and HiCLASS <cit.> derived from CLASS <cit.> although over a dozen more specific Einstein-Boltzman solvers have been written. Most of these codes have been cross calibrated and agree at the sub-percent level <cit.>.

While we have focused on how we can accurately model general deviations from GR, we now discuss how to constrain these deviations. The most fundamental constraint will come from measurements of the growth rate: from redshift space distortion (RSD), peculiar velocities, the integrated Sachs-Wolfe effect and from tomographic measurements of weak lensing. All of these will constrain the time evolution of δ_ M through f or fσ_8. To see how the growth rate is affected by deviations from GR, it is instructive to focus again on the quasi-static regime <cit.>. There, the evolution equation is now given by

    df/dln a+q(a)f+f^2=3/2Ω_ M(a)μ/γ

where q(a)=1/2[1-3w_ U(a)(1-Ω_ M(a))]. We can clearly see that the growth rate will be affected by both modifications to the background (which is encapsulated here in the equation of state of the effective fluid driving the acceleration, w_ U) and the modified perturbed equations through the quasi-static parameters. As a result there will be a degeneracy between w_ U and μ/γ <cit.> which can only be broken with measurements of cosmological distances. A particular parametrization of the growth is often used, of the form

    f=Ω_M^γ_G

where γ_G is approximately a constant for any specific theory over a narrow range of redshifts <cit.>. For ΛCDM, it can be shown that γ_G=6/11.

Weak lensing will also probe the growth rate – via the modified Newton-Poisson equation it measures the change in the density contrast <cit.>. But it goes further; if we look at equation <ref> we can see that it depends on both Φ and Ψ. If we now focus on convergence we have that κ∝∫ dχ∇^2(Φ+Ψ)∼∫ dχ∇^2(1+1/γ)Φ∼∫ dχ3/2(aH)^2Ω_ Mμ(1+1/γ)δ_ M where we used the Friedman equation to simplify the result. In other words, as well as depending on the growth rate (via δ_M) it probes one of the derived quasi-static parameters described above: Σ. 

Extending the set of cosmological parameters to include extra parameters (such as μ and γ) may lead to degeneracies – we have already seen such a case between w_ U and μ/γ in the case of the growth rate. A particularly problematic parameter is galaxy bias, i.e. the mapping between the galaxy density contrast, δ_ G and the underlying matter density contrast, δ_ M. A simple approximation is that δ_ G≃ b_ Gδ_ M where b_ G is scale and time dependent (although there are now, far more sophisticated models of bias <cit.>). Given that many direct probes of the growth rate are in terms of the galaxy distribution, one often finds that resulting constraints are on f/b_ G. 

A proposal to circumvent the dependence on bias is to consider combinations of statistics that, to some extent, "divide out" the dependence on bias <cit.>. A notable example is

    E_G(R)∝ξ_ GL(R)b_ G/ξ_GG(R)

where ξ_ GG(R) and ξ_GL(R) are the 2-D galaxy-galaxy and galaxy-lensing correlation functions and b_G is an independent measurement of the galaxy bias of the same sample. This statistic was originally formulated in Fourier space but is implemented in real space where a number of practical issues need to be considered <cit.>. Alternatively one can consider the full suite of correlation functions (including the lensing-lensing correlation function) which, combined, can break the degeneracy between b_ G and other parameters.




§ FIFTH FORCES AND GRAVITATIONAL SCREENING.




 §.§ Gravitational regimes and the observational desert.

In the previous section we focused on what one might consider to be a "clean" regime, where dynamics is well described by linear cosmological perturbation theory, and showed how it was relatively straightforward to parameterise general theories of gravity. This regime benefits from the fact that the dynamics is perfectly understood and it is possible to solve the evolution equations with arbitrarily high accuracy. As we saw in Section <ref>, it is where modern cosmology has had a resounding success. But it is also a regime with severe limitations: the larger the scales we have, the less statistical power to constrain parameters. We observe one universe and hence there is a finite amount of information to be extracted from its observation; this cosmic variance is largest for the largest scale modes. It would be desirable to look at smaller scale modes but, as we have seen, non-linear physics (both gravitational and baryonic) play an unavoidable role there.




There is another way of looking at the large scale, cosmological regime in comparison to other regimes <cit.> and this is best illustrated in Figure <ref>. On the horizontal axis one has the typical gravitational potential, Φ (in this plot labelled by ϵ) of a system while on the vertical axis one has the typical curvature  (labelled by ξ), here expressed in terms of the Kretschmann scalar (i.e. a quadratic invariant constructed from the Riemann tensor). Scattered on the plot are a number of different systems (or regimes in which gravity can be constrained), from the Solar System and pulsars to black hole mergers and imaging. Cosmology stands apart in a regime of low curvature.

A striking feature of Figure <ref> is that there seems to be a "desert", a dearth of observation between cosmology and other probes, for values of ξ between 10^-50 and 10^-42. The typical systems found in this desert are galaxies and clusters where the interplay of non-linear gravitational and baryonic physics as well as highly non-linear feedback mechanisms (such as supernovae and Active Galactic Nuclei - AGNs) make them particularly hard to model. Indeed we can see this desert as the transition between two easy regimes: the high curvature part where one is studying isolated systems and the low curvature part where linear theory on a homogeneous background applies.

There is another interesting feature that can be observed. In the high curvature regime, there are strong constraints on gravity - in the case of the Solar System and millisecond pulsars, exceptionally so. One is guaranteed to have modifications to GR in the top right-hand side as one delves into the regime where gravity becomes much stronger. But the current quest for modifications are in the bottom half of the regime where, we argued in the introduction, it is somewhat unnatural for deviations to GR to arise. If it turns out that GR is strictly true in the top half of the plane and violated in the bottom half of the plane, it would seem that some mechanism is shielding regions of high potential from deviations, i.e. from extra fifth gravitational forces and that the transition occurs somewhere in the desert.



 §.§ Gravitational screening.

While the description above is somewhat qualitative, it turns out that there are mechanisms of what is known as gravitational screening which can restore GR in regions of high density, potential or curvature but which allow for the emergence of fifth forces elsewhere <cit.>. A useful, alternative, way of thinking about screening is to realize that a fifth force is always associated to a "charge" and that the strength of the fifth force is proportional to the charge. One can then envisage a mechanism such that  in the screened regime the charge is surpressed and the fifth force is absent while in the unscreened regime, the force has full effect. 

To understand how screening can work, let us revisit one of the simplest scenario, a non-minimally coupled scalar field with a potential <cit.>. It is always possible to transform such a system into a frame in which the scalar field becomes minimally coupled to gravity but non-minimally coupled to the matter sector. In the simplest case, the background scalar field obeys an equation of the form

    φ+dV_ eff/dφ=0

where the effective potential, V_ eff=V-ρ̅φ, i.e. the minimum of the potential will depend on the local matter density. If we solve for dV_ eff/dφ(φ_0)=0, we can then expand φ=φ_0+δφ and we find that, in the quasi-static regime, the perturbations in the scalar field obey an equation of the form

    ∇^2 δφ+m^2(φ_0)δφ∝ " source"

where m^2(φ_0)=d^2V_ eff/dφ^2. If we assume a localized mass, M this means that

    δφ= βMe^-m(φ_0)r/r

where β is a constant of proportionality set by the coupling between the scalar field to matter and we see that m(φ) functions as a cutoff, setting the range of the fifth force. Furthermore, the cut-off depends on the environment so that, if we choose a V∝φ^-n, we have that for large ρ, m(φ_0) is large and the force is short range but if ρ is small, m(φ_0) will be small and the force can be long range. Thus we have that the fifth force can strongly satisfy local constraints yet have an impact on cosmological scales.

We have discussed one particular type of screening, known as Chameleon screening. Interestingly, it has been shown that, if Chameleon screening is to arise in a cosmological context, it has to be on scales less than a Megaparsec, i.e. deep in the non-linear regime; furthermore, the screening scalar field cannot be responsible for accelerated expansion directly through gravitational effects and only as a conventional source of dark energy <cit.>. There are, however, other forms of screening that we should mention. The first one arises in theories in which the scalar mediator of the fifth force has a non-trivial kinetic term, as in, for example bimetric theories of massive gravity or scalar-tensor theories constructed from shift-symmetric term – Galilleons <cit.> – or bimetric theories such as massive gravity. Known as the Vainshtein mechanism <cit.>, the modified kinetic term now transforms the overall strength of the fifth force, the β in equation <ref>, into a function of φ_0 which in turn can depend on the source and environment. Typically one find that near very massive objects, β→ 0 while around lighter or more diffuse objects, β can be non-negligible. Another variant is Symmetron <cit.> where an environmental dependence of the coupling between the scalar field and matter leads to a similar effect on β.



 §.§ The phenomenology of screening.

Screening adds more complexity on to what is already a particularly difficult regime to model. While the payoff for looking at smaller scales is greatly enhanced statistical power (there are more independent samples on smaller scales), one is faced with the non-linear aspects of gravity which are far less well understood than the clean linear regime discussed in the previous section. Furthermore, baryonic physics can now have a significant impact on attempts at making predictions.

There have been some attempts at developing analytic techniques for studying the impact of screening on the non-linear regime. Working in the quasi-static regime, it is possible to extend or modify higher order perturbation theory methods to included the effect of fifth forces <cit.> although the effectiveness of such methods has been somewhat limited. More progress has been made on developing enhanced N-body codes which incorporate the effect of fifth forces <cit.>. There, again using the quasi-static approximation, it is possible to add corrections to the gravitational forces which captures the effect of environment dependent fifth forces. Furthermore, it is possible to couple the modified gravitational evolution with hydrodynamic equations to allow a first foray into the impact of baryonic physics and how it feeds back into gravitational collapse.

Currently, the main codes used for studying the non-linear regime and its interplay with baryons are ECOSMOG <cit.> and ISIS <cit.> which are modified versions of RAMSES <cit.>,  an adaptive-mesh-refinement code, and  MG-GADGET  <cit.> which is a modified version of GADGET <cit.>, an SPH code. There have also been attempts at writing faster N-body solvers with fifth forces, based on particle-mesh (PM) solvers, most notably modified versions of COLA <cit.>. As with the Einstein-Boltzman solvers, there has been a concerted effort to cross-calibrate various codes so that they agree for a range of models on a range of scales <cit.>. 



 §.§ Observational consequences.

If gravitational screening is responsible for the qualitative difference between different regimes, then one needs to take a more nuanced view when undertaking a statistical analysis of large scale structure. In particular, one needs to look at the bigger picture, not only focusing on strictly local values of correlation functions or number counts, but including information about the environment in which these quantities are being estimated. One approach is almost geographic in nature: to construct a gravitational map of the Universe which allows us to identify regimes in which systems may be unscreened.

There have been attempts at constructing a gravitational map of the local Universe (i.e. which extends out to approximately 200h^-1 Mpc) <cit.>. The basic idea in <cit.> is to use a catalogue of galaxies which is as complete as possible in that region – the 2M++ catalogue. One then uses abundance matching to identify the possible dark matter halos associated with those galaxies and to fill in for any smaller halos which will not have been picked up by the survey. Finally, one uses a constrained realization of the long wavelength modes to fill in the rest of the dark matter on large scales. From the combination of all these terms it is possible to construct a map of the gravitational potential, tidal forces and local curvature. It is also possible to model the uncertainty in all the steps leading to an overall local error budget on the map.
This map can then be used to distinguish between the regions that are unscreened and screened.

Once one can identify qualitatively different regions, it is then possible to look for signatures of the unscreened phase, i.e. for evidence of the unsupressed fifth force associated with deviations from GR. The way to go about doing this is to try and identify objects which might behave differently in unscreened environments: some object will either be sufficiently dense, compact or large that they will self screen themselves (i.e. their fifth force charge will be effectively zero) while others will be too diffuse or light to be self screened. The hall mark of the fifth force in unscreened regions will be what seem like effective violations of the equivalence principle.

In this context, there are a number of interesting effects that one might look for <cit.>. For a start, the different components of a galaxy will behave in different ways. For example, in a theory with chameleon (or symmetron screening) the stars in a galaxy will self-screen while the gas will be too diffuse and will not be screened. The net result is that there should be an offset between the optical and radio components of a galaxy if there is screening. By looking for this effect it should be possible to place constraints on the amplitude and range of the a screened fifth force.

Another possibility is by looking at how a galaxy is deformed as it sits in a dark matter halo. The stellar disk might lag behind the halo centre inducing a gradient in the gravitational potential. As a result the galaxy will be warped into a cup shape where the amount of deformation will depend on the strength of the fifth force. The effect will be strongest if the disk is completely perpendicular to the external unscreened force field. 

A more straightforward effect will be that the unscreened component of the galaxy – the gas – will have a faster rotation curve relative to the stellar disk as a result of the additional fifth force having the same direction as the normal gravitational force. Furthermore, in the case of edge-on infall, the stellar and gas disks (and rotation curves) will have greater differences or assymetries with the edges which are closest to the centre of the halo being more compressed than the far side.

Finally, a similar and striking signature may be the relative offset of the galaxy centroid and the super massive black hole residing at its centre <cit.>. If a galaxy falls in an unscreened potential, much of the galaxy will feel a fifth force while the black hole, as a result of the no-hair theorem, will be shielded and won't be subjected to it. 

Another line of attack is to roll out suitably tailored versions of the usual n-point statistics which are used in standard cosmological analysis. In particular, it makes sense to correct estimators so that they downweight screened regions and enhance unscreened regions <cit.>. In the case of chameleon screening, this can easily be achieved by downweighting regions of high density, by "clipping" high density peaks or by introducing some form of nonlinear function of the density field δ̃[δ_M] such that high density values are damped. Examples of such transformations are δ̃=ln(δ+1) and δ̃=[(ρ_*+1)/(ρ_*+ρ_m)]^p where ρ_* and p are free parameters.

An approach that has been gaining traction is to focus on what are, for sure, the most unscreened parts of the Universe: voids. Voids have a scale of about 20-100h^-1 Mpc and are typically surrounded by outflows; they are seen as a new and interesting arena in which to do cosmology, parts of the density field where many of the linear and very mildly linear results still hold. Furthermore, with the vast increase in the size and depth of cosmological surveys, it is now possible to construct meaningful statistics with void catalogues in the same way as one does with galaxy clusters. And given that voids are particularly prone to be affected by unscreened fifth forces, it makes sense to focus on them if one wants to test gravity.

As yet, the exact definition of a void has not yet been completely established and agreed on <cit.>. A number of different factors come into play, from dimensionality to the type and number density of the tracers used, as well as, operationally, the type of void finding algorithm that is used. As a result, void catalogues can be quite different from each other and the resulting statistics can disagree. 

Nevertheless, if one settles on one particular catalogue (or better, considers a range of catalogues) one can proceed to assess a number of statistics. In particular, one looks at the usual statistics applied in cosmological analysis except now in this very particular environment. So, for example, one can assess lensing of background galaxies by the void <cit.>. One can also measure redshift space distortions in the void and compare them to measurements outside the void <cit.>. Another approach is to cross correlate the void profile with the CMB <cit.> in the same direction and in doing so, pull out the integrated Sachs-Wolfe term (the last part of equation <ref>). Finally, one can approach voids in the same way one approaches clusters and quantify their abundance as a function of scale and redshift <cit.>.




§ CONSTRAINTS FROM LARGE SCALE STRUCTURE.




 §.§ The evidence for ΛCDM

We seem to be living in charmed times, at least from the point of view of cosmological observations. With observations of the CMB – with the COBE, WMAP and PLANCK satellite missions and a host of ground-based experiments – complemented by the burgeoning precision probes of LSS, we now have an accurate mathematical model of the origin and evolution of the Universe in which the six parameters that govern it have been measured to exquisite precision <cit.>. Based on GR with a cosmological constant (and described in Section <ref>), ΛCDM fits almost all the data with remarkable consistency. 

The observational apparatus that we have to support ΛCDM is formidable. As well a the CMB experiments, the main observations to contribute to cosmological constraints are: spectroscopic surveys – BOSS <cit.> and WiggleZ <cit.> – measuring galaxy spectra (i.e. redshifts) and angular positions; photometric and weak lensing surveys – KIDS <cit.>, HSC <cit.>  and DES <cit.>  – measuring galaxy images and shapes, angular positions and a crude measurement of redshifts. There are preliminary attempts at line-intensity mapping (for example measuring the 21 cm spin flip transition of neutral hydrogen) but nowhere near at the level required to do cosmology <cit.>.

Just to belabour the point, it is useful to see how well some of the fundamental parameters of ΛCDM are currently constrained <cit.>. The initial conditions of fluctuations are characterized by a spectal scalar index, found to be n_S=0.965±0.004, i.e. almost but not exactly scale invariant; the difference from scale invariance is 8-σ, consistent with (and, indeed, some might say strongly indicative of) a period of primordial inflation. This constraint is supplemented with an upper bound on the tensor to scalar ratio (i.e a measure of the amount of gravitational waves left over from inflation) of r<0.07. Another notable constrain is on the geometry of the universe: the current, tightest constraint is Ω_K=-0.001±0.002. The Hubble parameter is constrained, from cosmological data, to be H_0=67.4±0.5 km s^-1 Mpc. 

A key incentive to consider modifications of gravity is the evidence for accelerated expansion. While ΛCDM is a remarkable fit to the data, it is instructive to consider constraints on the equations of state of the form w_ DE(a)=w_0+w_a(1-a). The most complete analysis using a comprehensive array of data sets <cit.> find that for a constant equation of state (i.e. w_a=0), w_0=-1.028±0.032 if one uses the Planck CMB data combined with distance measures from supernovae and the BAO in the matter power spectrum. If one allows for w_a=0, the same collection of data sets lead to w_0=-0.961±0.077 and w_a=-0.28^+0.31_-0.27 while using the Planck data, 
the BAO and redshift space distortion measurements of growth from galaxy surveys and galaxy weak lensing data lead to w_0=-0.76±0.2 and w_a=-0.72^+0.62_-0.54. There is clearly some lee-way in terms of an evolving equation of state but, as discussed above, this can't necessarily be used to find any evidence for modifications to General Relativity.



 §.§ Minor anomalies.

Given these remarkable constraints and the success of ΛCDM, it is useful to look at the (minor) inconsistencies and anomalies. While I do not yet believe they are at the level that call into question the success of the standard model, they will be useful in interpreting constraints on gravity later on in this section. There have been a number of small internal anomalies flagged in the current CMB data: large scale surpression of power and evidence for mild anisotropy in the low quadropoles <cit.>, mild inconsistency between large angle and small angle angular power spectrum constraints of cosmological parameters <cit.>, a discrepancy in the constraint on σ_8 from primordial and cluster abundances from Sunyaev-Zeldovich measurements <cit.>, etc. 

A persistent anomaly is the discrepancy in the constraint on H_0 from local measurements, such as those using Cepheids in <cit.> as compared to cosmological measurements from, for example, the CMB <cit.>. While both groups claim percent level constraints on the expansion rate, the discrepancy between the two results is much greater, at a level of reasonable significance. Most attempts at coming up with an explanation for this discrepancy from exotic physics do not explore deviations from GR and I will therefore not pursue this anomaly here.

The only anomaly I would highlight here is the inconsistency between the amplitude of lensed CMB fluctuations as inferred directly from the angular power spectrum of temperature fluctuations and the power spectrum of the lensing potential  reconstructed from the maps <cit.> . In particular, if ones free up the amplitude, A_L (such that A_L=1 is the ΛCDM value) of the lensing contribution that smooths the the angular power spectrum of temperature fluctuations (as described in Section <ref>), current constraints from Planck give A_L=1.243±0.095 at the 68% confidence level. Paradoxically, if one attempts to reconstruct the power spectrum of the lensing power spectrum, C^ϕϕ_ℓ, directly from the maps one finds that it is consistent with what one would expect from ΛCDM; i.e. two different methods for estimating the lensed CMB from the same data set give different results.

As mentioned in Section <ref>, a key observable to test gravity is the growth rate of structure, i.e. the speed with which gravitational collapse ensues. The main method for measuring the growth rate is with Redshift Space Distortions (RSD). RSDs can be seen in the redshift space correlation function (or power spectrum): along the line of sight, the proxy for distance – the redshift – will have a correction from peculiar velocities due to the local gravitational potentials <cit.>. This means that correlations along the line of sight will be deformed relative to correlations tangential to the line of sight, i.e. the correlation function (or power spectrum) will be anisotropic. On linear scales – scale of tens of Megaparsecs – the correlation function, or power spectrum is squashed while on smaller, non-linear scales, they are stretched.


A mild (and arguable) anomaly is the evidence for low growth rate coming from RSDs. First highlighted in <cit.> but also thoroughly refuted in <cit.>, the discrepancy can best be assessed in Figure <ref> where we can see the predicted, density weighted growth rate for ΛCDM, marginalized over a host of cosmological parameters and constrained by the Planck 2018 data <cit.>. Two things are of note: first, that current low redshift constraints (from RSDs) are weak compared to predicted CMB constraints and second, that many data points seem to be below the predicted line. I emphasize that the evidence for the low growth rate is not yet significant and, for some, not even indicative. 



A much more significant inconsistency arises in the comparison between galaxy weak lensing constraints and the CMB. If we focus on S_8=σ_8(Ω/0.3)^1/2, we can see in Figure <ref>, that the galaxy weak lensing measurements from the KIDS survey <cit.> leads to a somewhat lower value of S_8 than inferred from the Planck CMB measurements. This discrepancy is present in the DES data <cit.> although at a milder level; it is possible to further soften the difference with a more elaborate likelihood which, one hopes, can capture the statistical and systematic uncertainties more accurately. But, crucially, these two, completely independent, imaging surveys  both find this inconsistency and combining them, or including more data  may only reinforce this discrepancy.


So, in summary, current data consistently points to ΛCDM, which is firmly rooted in GR. There are some mild inconsistencies which are not statistically very significant but might lead one to explore deviations from the standard model. We can now look at what current data has to say about extension of GR. 



 §.§ Parameterizing gravity.

In Section <ref> we discussed how to parametrize deviations from GR, in the linear regime, on large scales. In the case of scalar-tensor theories described by the Horndeski action, this is straightforward to do in terms of the five α_X parameters on both super-horizon and sub-horizon scales. It is also true that most of the statistical power will come from smaller scales where, generally, the quasi-static approximation can be used. This greatly simplifies the analysis as one is restricted to two free function, μ and γ. In most of what follows, we will present constraints in terms of these parameters or combinations of them. But, before we do, a few comments are in order.

Most parameter constraints using μ and γ assume time dependent, but scale free, functions which are valid from sub-horizon to super-horizon scales. This is clearly not correct, in the sense that it does not naturally emerge from the approach described in Section <ref>: one should, at least, expect a scale dependence at or greater than the horizon scale but more generally, the quasi-static approximation breaks down and one needs to, correctly, include the dynamics of the extra degree of freedom (which will, for example, depend on its initial conditions). Having said that, most constraints from LSS are, effectively, on sub-horizon scales (where there is more statistical power) and so one hopes that any inconsistencies in modelling very large scale deviations will be vastly supressed by the data on smaller scales. Effects on the scale of the horizon may be important for large enough surveys (although see <cit.> for how insignificant the constraining power of horizon scale modes is) or in the case of the primordial CMB (although these occur at z≃ 1000 and thus are unaffected by late time deviations from GR). The advantages of focussing on μ, γ outweigh the problems just described. Nevertheless, later on we will look at constraints on scalar-tensor theories where the large scale behaviour is correctly taken into account.

Another, important, point needs to be highlighted. If we restrict ourselves to scale free  μ and γ, one still needs to choose an appropriate time dependence. Strictly speaking, these function will depend on a, ȧ and, in the case of scalar-tensor theories, the background scalar field value, φ_0, all of which depend on t. A common approach has been to assume that deviations from GR are tied to the onset of accelerated expansion which can be encapsulated in the effective value of the background energy density, Ω_ DE or Ω_U <cit.> in the notation introduced in Section <ref>. So, for example, one can take μ=1+α×(Ω_ DE/Ω_ DE0) where Ω_ DE0 is the value of Ω_ DE today, α is a constant and one often refers to μ_0≡μ(z=0)=1+α (and likewise for γ). 
A more general approach is to divide the temporal dependence of μ and γ up into uncorrelated redshift bins – for example, there will be a low redshift bin starting off at z=0 followed by a succession of bins. In fact, this approach is often generalized to include scale dependence, i.e. bins in the fourier mode, k, as well as z.

There are a few drawbacks to any choice of time dependence. First of all, it has been argued that the form μ=1+α(Ω_ DE/Ω_ DE0)  (or of any parameters dependent on Ω_ DE) does not necessarily reflect the evolution of μ and γ one finds in specific models <cit.>. While this may be true in some specific cases,  this time dependence is, generally, a rough approximation of what one might find in theories in which the accelerated expansion is tied to deviations from GR. On the other hand, dividing up the temporal in redshift bins often leads an ill-posed problem (i.e. to many parameters to fit, given the data) and for which results are difficult to determine. Both of theses approximations can be sharpened if, for example, one considers a more systematic expansion of μ and γ in terms of some basis functions, f_n(Ω_ DE) for example, that can capture more elaborate time dependence (for example choosing f_n(x)=x^n or f_n(x)=j_n(x) where the j_n are spherical Bessell functions) or by considering some form of Principal Component Analysis of the time dependence of these functions.  One then has to be careful in choosing priors that adequately reflect the physical behaviour of this functions <cit.>

A second drawback then emerges in that the results one finds are, for current data, strongly dependent on the parametrization one considers (or for example the number of free functions one incorporates into the parametrization). For example, with current data it is possible to constrain a constant μ at the percent level while μ=1+α(Ω_ DE/Ω_ DE0) is constrained, at best, at ten to twenty percent level (an example of this will be seen later on, when we focus on the BOSS data). Having declared the problems, we will stick with the conventional parametrization in the discussion that follows. 

Note that most current analysis combine multiple data sets – e.g. CMB with galaxy weak lensing – so the division between the different results that follow us somewhat artificial. Nevertheless I choose to do so to highlight the impact of the different observables, and their specific effects, on constraints on gravity.





 §.§ Current constraints: the CMB

If I take a selection of constraints on deviations from GR, I naturally begin with the analysis of the Planck data (combined with other data sets). In <cit.> and in <cit.>, a comprehensive analysis of different parametrizations was undertaken with the 2015 and 2018 Planck data releases; the overall take home message was that current cosmological data is broadly consistent with GR. But of particular interest (and representative of the analysis as a whole) are the features seen in Figure <ref> from <cit.>, where constraints on μ̃ (here labelled μ_0) and η are presented.  To begin with we can see that there is a degenerate direction along which constraints are much broader - this direction is set by constraints on the Φ+Ψ from CMB lensing and thus on Σ (if one transforms to the μ, Σ parametrization, one can, to some extent decorrelate the two parameters being constrained). 

Second, we note that the error bars on these parameters are of order 1. While these may seem very weak compared to the precision one obtains on astrophysical scale (for example, the constraints on the analogous parameter to γ from the Solar System is of order 10^-5) they are obtained in a completely different regime, on scales which are fifteen orders of magnitude greater and where the curvature scale is  twenty orders of magnitude smaller (see Figure <ref>).

Finally, and most striking, is that fact that there seems to be a  detection of deviations from GR, away from μ_0=γ_0=1. This push away from GR from the Planck data alone is due to the enhanced amplitude of CMB lensing anomaly (the "A_L anomaly") previously discussed in this section, as can be seen in the left hand panel of Figure <ref>. In the right hand panel we see that,  if one then includes the measurement of the lensing potentials from the maps – C^ϕϕ_ℓ – (which, recall from the discussion above, is an independent measurement to the one used to determine A_L), the constraints are pushed towards GR. Including galaxy weak lensing (in this case from DES) which, as we saw, favours a lower amplitude in the lensing potentials, further pulls the constraints towards GR. Adding the RSD constraints has a small effect in the orthogonal direction.

Is this a significant detection of deviations from GR? Not if one properly accounts for the increase in complexity of the model, even though there is an improvement in the goodness of fit.

GR would have to be excluded at a much higher significance for the detection to be noteworthy. 





 §.§ Current constraints: galaxy weak lensing and growth rate.


We now turn to a different set of constraints, obtained from a combination of galaxy weak lensing from the KIDS-450 with the CMB and RSDs <cit.> . In this case, the parameters being constrained are μ̃ (here dubbed Q) and Σ, binned in both redshift and scale. If we look at the constraints for the lowest redshift bin (labelled by a "2"), on largest scales (assume that we are marginalizing over the other bins) and focus on the case of CMB alone (the red dashed lines in the right hand panel) we see, again, a mild inconsistency with GR; in this case it is milder than in the Planck analysis and comes tied to vastly inflated uncertainties (as is expected from our discussion above, given the number of free parameters – or bins – now included). Weak lensing and RSD data (without CMB) leads to tighter constraints, centered on the GR value. To combine them with CMB data, the authors have introduce a large scale cut to mitigate the effect of non-linear scales. In that case, the constraints including the CMB are completely consistent with GR; again the low lensing potentials from galaxy weak lensing compensate the high A_L constraint from the CMB.
An interesting effect of enlarging the parameter constraints to include deviations from GR is that the "S8" discrepancy can be substantially mitigated. In the left hand panel of Figure <ref> one can see that what were orginally disjunct contours (the dashed lines) can be rendered overlapping (the teal and red shaded contours) by adding in μ and γ. We note that the Dark Energy Survey (DES) has undertaken an analysis of their first year of data <cit.>  to come up with constraints on μ and Σ. Their constraints are of the same order of magnitude as the KIDS-450 constraints although marginally more consistent with the GR values. 





We now briefly look at the SDSS-III
BOSS constraints in <cit.> where the authors concentrate on G_M and G_L which, we know from Section <ref>, can be directly related to μ and γ. Assuming a time dependence for these quantities of the form G_X=1+(G_X^(s)-1)a^s where X=M,L and, in their analysis, s∈{0,1,3}. In Figure <ref> we see how a time independent choice (s=0) leads to very strong constraints while s=1 (and even more so with s=3) leads to weaker constraints as deviations from GR are more suppressed at high redshift. 
 Constraints on the growth index are found to be γ=0.566±0.058 which, again, show strong consistency with GR. 

It is important to note that the analysis of  data to constrain fσ_8 involves a number of steps in which a particular cosmology and theory of gravity is assumed. A concern, then, is that estimates of fσ_8 may be biased by these assumptions. Attempts at estimating the size of this bias <cit.> show that it is insignificant for current surveys but may be an issue when their statistical power greatly improves.
 
 

 A parallel line of attack has been to determine constraints on GR from E_ G. There have been a number of attempts at measuring E_ G using galaxy-galaxy lensing <cit.> and CMB lensing combined with Galaxy velocities <cit.> finding this statistic to be generally consistent with GR, but with a few minor complications (for a thorough analysis <cit.>). First of all, in the case of <cit.>, the large angle contribution of CMB lensing leads to a 2.6σ inconsistency with GR. Second, the choice of background cosmology has a strong impact on the estimate of E_G; for a comfortable consistency with ΛCDM, one needs a slightly lower Ω_M than normally assumed.  The current, state of the art, constraints on E_G are given in Figure <ref>.
 
It is not clear whether E_G is a particularly powerful statistic for testing GR and whether it adds anything more than the other, more standard likelihood analysis currently being undertaken. While, in principle, it should mitigate the contaminating effect of galaxy bias, in practice and due to the way it is implemented, one is still at the mercy of one's best bias model. Furthermore, a proper, full likelihood analysis of all cosmological data will necessarily include cross-correlations between lensing and density measurements and will inevitably work towards breaking any degeneracy between nuisance parameters such as bias and the GR parameters one is attempting to constrain.  
 


 §.§ Model specific constraints: scalar-tensor theories.

Instead of using a general, parametrized approach to constraining GR, it is sometimes useful to focus on specific models. This makes it possible to, to some extent, span different regimes and compare, for example, the constraining power of Solar System probes versus cosmological probes. The workhorse for such constraints is Jordan-Brans-Dicke gravity, the simplest scalar-tensor theory described in Section <ref> with its one free, dimensionless parameter, ω_ BD. Current astrophysical constraints of ω_ BD are remarkable. From the Cassini time delay it is possible to constrain ω_ BD>40,000 <cit.> while from milli-second pulsars the constraints are around ω_ BD>20,000 <cit.>. Current cosmological constraints, using the Planck data are much weaker: ω_ BD>900 and come primarily from the effect of the extra scalar field on the background expansion and thus on the angular diameter distance <cit.>. Measurements of the growth rate are not yet up to the level of precision required to lead to competitive constraints on ω_ BD.

Jordan-Brans-Dicke gravity is not screened on astrophysical scales but more general scalar-tensor theories will be. For example Galileon theories may have Vainshtain screening while conformally coupled theories with potentials may have Chameleon or symmetron screening. If that is the case, then the only way of obtaining tight constraints may be in the cosmological regime (although see Sections <ref>  for other approaches). It would make sense then to look at constraints on a general scalar tensor theory such as the one described by the Horndeski action. Furthermore, as we saw above, the μ, γ parametrization shouldn't be extended to horizon and superhorizon scales (although these scales won't have a strong statistical weight and will contribute little to the final results). The Horndeski action allows a completely self consistent parametrization on all scales in terms of the α parameters.

The most comprehensive analysis of the Horndeski action was undertaken in <cit.>, where the authors used CMB data (a combination of the 2013 Planck power spectra and WMAP polarization), galaxy redshift survey data in the form of the power spectrum (from WiggleZ), the BAO (from BOSS, WiggleZ and 6DFGS) and RSDs (from SDSS, BOSS, WiggleZ, 6DFGS and VIPERS) to place constraints on the αs. In their analysis, the authors used a parametrization of the form α_X=c_XΩ_ DE (where c_X is a constant), along the lines of the one use for μ, γ. On the scales being constrained, α_K was found to be irrelevant and could be fixed to fiducial values without affecting the final results. The typical variance on c_X was found to be O(1) and the parameters were found to be mostly uncorrelated. 

Interestingly, the analysis of <cit.> found that the best fit values of α_B and α_T were away from GR and, for some combination of data sets, at the 95% confidence level. While for CMB alone, the results were perfectly consistent with GR, the inclusion of RSDs pulled the preferred values in such a way that c_T was preferentially less than zero and c_B was preferentially greater than zero. This is somewhat at odds with the analysis of μ, γ in <cit.> where it is the A_L anomaly that is pushing the constraints away from GR and where the RSDs do not play a significant role. The authors do not, however, claim a significant detection of deviations from GR on the grounds that, having introduced more parameters, the evidence (which they compare to GR using the Bayes factor) is not sufficiently significant. In other words, the increase in complexity in choosing Horndeski is not compensated by the improvement in the fit to the data. 



A restricted analysis of Horndeski gravity taking into account the constraints from GW170817 on the speed of gravitational waves has been undertaken in <cit.> . The authors found that, varying α_B, the values of α_M were discrepant with GR at the 95% confidence level (but fully consistent if α_B=0)  and they argue that this is due to the fact that, with the specific choice of kineticity α_k=0.1, an artificial prior is induced which biases the results.
More recently. the analysis of the Planck 2018 data <cit.>, for the parameter subspace where α_B=-α_M found slight preference for α_M<0 which disappeared once RSD and WL data was included, consistent with their analysis in term of μ, γ.

Of particular interest are the current constraints on Galileon theories. The constraints on the speed of gravitational waves from GW170817 rules out theories with quartic or quintic terms in the action of equation <ref>. Constraints from large scale structure and, in particular, the Integrated Sachs Wolfe effect in equation <ref> can be used to show that the cubic terms lead to a grossly discrepant effect <cit.>. This means that cosmological Galileon gravity is effectively ruled out.

As we have seen, they current situation is that constraints on different parameterizations of deviations from GR on cosmological scales are still weak. Yet, in particular situations there is mild evidence for non GR behaviour. In one case, this evidence is primarily driven by the high A_L anomaly in the angular power spectrum of the CMB and can be compensated by the low amplitude of galaxy weak lensing or RSDs. In another analysis, in this case of Horndeski theories, it is the RSDs that are pulling the constrained values of the αs away from GR. Furthermore, the results depend quite significantly on the time dependence assumed or the cuts to which the data sets have been subjected. Given that the attempts at constraining deviations from GR are so disparate in methodology, as well as few and far between, there is a serious concern that a complete and final set of reliable, calibrated results is still lacking and that can be taken as a basis for more definitive statement about the status of GR on cosmological scales.  Furthermore, the overriding consensus is that there is no evidence for deviations from GR from a Bayesian evidence point of view: the complexity of the extensions outweighs the improvement in the fits.

Nevertheless, these various analysis' have shown the way forward in how to constrain gravity on cosmological scales. A clear approach to the analysis is emerging and, more importantly, a sense of the systematic problems and pitfalls that must be dealt with if one is to obtain significant constraints of GR with current and future data.




§ ALTERNATIVE CONSTRAINTS - CLUSTERS, GALAXIES AND VOIDS.


As we have seen, there has been some progress in applying the techniques for constraining GR on large scales to existing data, with mixed results. The situation is less developed on smaller scales where, on the one hand, non-linear gravitational collapse and baryonic physics can play a significant role and, on the other hand, the distinctive signatures of screening can emerge. Existing constraints are fewer, more scattered and, unlike on large scales, very model specific. In this section we will undertake a brief survey of the more durable results.




 §.§ Cluster masses.

To begin with, and as we saw discussed in Sections <ref> and <ref>, it makes sense to look for evidence for gravitational slip by comparing the way relativistic particles (i.e. photons) and non-relativistic particles (e.g. baryons or dark matter) respond to a gravitational potential and mass. One such approach is to compare different ways of measuring cluster masses: X-ray, Sunyaev-Zeldovich and weak lensing observations. Each one of these observations measures different things, from temperature and pressure profiles to gravitational potentials, and it is necessary to construct a model (using, for example the fluid equations and assuming hydrostatic equilibrium) to bring them all together. Such an approach has been undertaken in <cit.> and used to place constraints on the different free parameters in an individual cluster, including one that controls deviations from gravity. The calculation is done in the non-linear regime and has been done for one particular model, f(R); the authors found df/dR ≲ 5× 10^-5. This method has also been applied to constraints on the Vainshtein mechanism where no evidence for screening was found <cit.>.

Instead of looking at individual clusters, one can look at the statistics of ensembles of clusters and, for example, the cluster mass function, n(M). Arguably, this statistic will not be overwhelmingly contaminated by baryonic physics given the size of the objects and the scales probed. Attempts at constraining this statistic necessarily involve N-body simulations (although progress has been made developing semi-analytic methods) and are, again, model specific. Again, most effort has gone into testing f(R) (and some other alternatives) using  estimates of cluster abundances from a range of cluster catalogues <cit.>. The constraints are at the level of df/dR ≲ 1-10× 10^-5. 


 §.§ Screened fifth forces in galaxies.

Some progress has been made in looking for screening within galaxies. The latest, most detailed search for chameleon screening  has, arguably, lead to evidence of a screened fifth force  of the form F=-Δ G e^-mr/r^2 with a strength of Δ G/G≃ 0.01 on scales of 1/m∼ 1 Mpc. This evidence comes from two very different types of analysis. In the first case, a systematic offset between the screened stellar component and the unscreened HI component of a selection of galaxies is consistent with the presence of a weak fifth force <cit.>. A number of systematic effects have been included and are taken into account in the error budget although two caveats most be considered. If one considers a (overly) conservative error on the centre of the HI gas (effectively doubling the uncertainty that is presented in the data  set), the detection goes away and one obtains a strong upper bound,  Δ G/G<0.01-0.1 depending on the screening scale. But, and more importantly, little is know about how baryonic effects may lead to a signal that mimics the fifth force. While the signature that is constrained is very specific, a detailed analysis of ultra-high resolution simulations of galaxies in a large enough volume is required to completely confirm or rule out the detection.

The second type of the analysis looks at another signature described in Section <ref>: warped stellar disks as the galaxies fall into their (displaced) dark matter halos <cit.>. There the authors used r-band images of mostly late time galaxies in the NASA Sloan Atlas to compile a catalogue of approximately 10^3 warps with a specific focus on the "cup" shaped signature of a screened fifth force. Again, they found a strong indication of a fifth force with Δ G/G ≃ 0.01 on scales of around 2 Mpc, surprisingly consistent with the measurement of star-gas offsets. Remarkably, this new analysis is with an altogether different data set and involving a very different signature. Nevertheless it is potentially plagued by its own, as yet, unquantified systematics: the effects of baryonic physics is poorly understood but more importantly, there is no clear idea yet of how warps arise in galaxy. While there is a consensus that warps are induced by environmental effects, the fact that warps have been found in isolated galaxies calls this environmental origin into question. A better understanding of how warps form would greatly help in firming up this detection of a screened fifth force.

A different, but related, constraint on screened forces has been found by attempting an Eotvos like test on galactic/cosmological scales. In essence, the idea is to compare the motion of a "small" stellar system and a supermassive black hole in the gravitational potential of another, larger system. A black hole is chosen as one of the test bodies because in some theories <cit.>, it will satisfy a no-hair theorem and therefore will be unscreened (the authors are inspired by, and very specifically address, the case of Galileon theories which are endowed with Vainshtein screening).  As a result, and for a reasonable choice of the screening parameters, one should expect an offset of O( kpc) between a black hole and its host galaxy.

In <cit.> and <cit.>, the authors looked at the black hole centred on M87 as it falls under the influence of a sub-clump of the Virgo cluster around M84 and M86. They found that the fifth force field strength of the Galileon must be less than 10^3 (km/s)^2/kpc. The authors had to make a number of assumptions about the mass to light ratio, concentration, halo profile, etc, but ultimately found that their constraint is stronger than current constraints on Galileons from the Solar System. Again, much is to be done to better understand the systematics and it is clear that more robust constraints will be obtained with more objects but the route is promising for constraining (or detecting) Vainshtein screening.


 §.§ Cosmic voids.

The arena, par excellence, where constraints on the screening mechanism should be most effective, is in cosmic voids. As explained in section <ref>, one should expect forces to be unscreened as they are regions of low density, curvature or potential. There has been a reasonable amount of effort in thinking about what signatures to look for, primarily based on analytic calculations and numerical simulations. The use of voids as a tool of more standard cosmological analysis is also growing with first attempts at measuring the void-correlation function already in hand. 

Surprisingly, there are no de facto constraints on deviations from GR with current data. The only attempt, thus far, is embedded in efforts to constrain cosmological parameters from the void correlations function. The authors of <cit.> have found a constraint on the linear f/b=0.417±0.089 at median redshift, z̅=0.57, which is entirely consistent with Λ CDM. Thus, they claim, there is no evidence for deviations from GR within voids.




§ FUTURE PROSPECTS. 




 §.§ Theoretical status.

Testing gravity with cosmology has come a long way over the last decade. And, in some aspects, one might consider it a mature field. For a start, we have a much clearer understanding of what it means to construct alternatives to General Relativity. Primarily it involved adding extra degrees of freedom but do so one has to be careful that instabilities, ghosts or other pathologies aren't unwittingly introduced. And even if they are, they may turn out to be benign. So, while in the past, there were a set prejudices and a disparate set of rules of what one could or could not do when tampering with GR, there is now a much more coherent and sensible approach. One of the by products is that we have a clearer understanding of the landscape of gravitational theories and how they link together. For example, scalar-tensor theories take pride of place and are consumately well understood even though that understanding is not necessarily complete. And their structure can be mimicked in vector-tensor theories leading to a burgeoning exploration of this different part of the landscape. Tensor-tensor theories – massive gravity or bigravity theories – will reduce, in a particular sector to a type of scalar tensor theories but also have properties that are uniquely their own. And so on as one moves into higher dimensions and non-locality.

From the point of view of cosmology, and in particular on the largest scales, there is now a completely well understood theory of a linearly perturbed expanding universe. So much so that it can be (and has been) implemented in precise Einstein-Boltzman solvers that allow us to predict (linear) cosmological observables with exquisite precision. To a lesser extent, we also have a good idea of how smaller scales may be subjected to gravitational screening mechanisms that can mask fifth forces; we have a reasonable complete view of the different possible mechanisms and have numerical implementation of some of them. 



 §.§ Observational status.

When it comes to obtaining constraints from current data, the situation is more uncertain. The first steps have been taken and there are now a number of pipelines that can take current observations of the large scale structure of the universe – the CMB, galaxy and weak lensing surveys – and find a set of constraints on gravitational parameters, i.e. parameters that encapsulate deviations from GR on linear scales. The results are intriguing in that anomalies in the CMB power spectrum due to what seems like excess gravitational lensing push the constraints away from GR. Paradoxically, the low amplitude of lensing in galaxy gravitational lensing surveys can compensate when combined with the CMB data. Measurements of growth seem marginally lower than what is expected from GR yet do not yet have enough statistical power to play significant role.

An overall concern is the disparate set of methods and parametrizations being used by the (few) different groups which makes it difficult to compare and, more importantly, cross calibrate the different results. This differs from cosmological constraints of parameters within ΛCDM; there, a large number of groups have found constraints on, for example, the fractional densities of the different energy densities, the expansion rate, the spectral index of fluctuations, etc and there is a very clear and accurate consensus of what they are. Until this level of rigor is applied to the gravitational parameters, it is difficult to interpret what the current constraints actually mean or whether they are robust.

Nevertheless, where there are constraints, they are still weak, roughly of O(1). This means that, even in the case where there is a 3σ discrepancy with GR, the increase in complexity by including extra gravitational parameters still outweighs the improvement in the likelihood. This is typically characterized in terms of Bayes factors, comparing the Bayesian evidence between the two different scenarios. The fact that, from this point of view, deviations from GR are disfavoured relative to GR parallels other attempts at extending ΛCDM and merely confirms that the base model is still the best candidate for explaining the Universe.

On smaller scales, deeply into the non-linear regime where more complex, baryonic physics can play a role, tests of gravity are still in their infancy <cit.>. While it is an extraordinarily difficult regime to work in and prone to countless systematics, it is also, potentially immensely rewarding. Data on these scales is accruing at a phenomenal rate and, if systematics are properly modelled and under control, they should allow for constraints (or detections) with a high statistical power. The first steps have been taken with a few detections and constraints of an altogether new type of effect – gravitational screening – which have emerged. There is clearly a vast and uncharted territory to be explored.



 §.§ Planned surveys and prospects.


Although current constraints are weak and inconsistent, the hope is that with the new generation of surveys, the situation will change dramatically. Quite possibly, over the next decade, we will transition from O(1) constraints on deviations from GR to O(10^-1) or even O(10^-2). A number of planned experiments dominate the landscape. The Euclid satellite (http://sci.esa.int/euclid/), primarily funded by the European Space Agency, will survey 15,000 square degrees of the sky, building up a catalogue with spectra of upto 10^7 galaxies. In parallel, a ground-based spectroscopic survey, DESI (http://desi.lbl.gov), based in Arizona will collect spectra of galaxies of a similar number and area. The Large Synoptic Survey Telescope, LSST,  (http://www.lsst.org) based in the Atacama Desert will undertake a deep, imaging, survey (r>27)  with an area of up to 20,000 square degrees. The Square Kilometre Array, SKA, (https://www.skatelescope.org) will produce a radio survey of HI emission and galaxies in the frequency range 50-1760 MHz, mapping structure out to z≃ 3-4. The Simons observatory (https://simonsobservatory.org/index.php), also in the Atacama Desert will map the CMB to exceptionally high resolution and sensitivity, while the Stage IV CMB survey (https://cmb-s4.org) will, hopefully, be a coordinated set of multiple ground based observatories which will cover upto 40% of the sky.

One can forecast the precision with which one will be able to constrain the gravitational parameters. For example, with Euclid or LSST, we expect to reduce the uncertainty in μ and γ from O(1) today down to O(10^-2) . If we restrict ourselves to scalar-tensor theories described by the Horndeski action, again, and consider combinations of LSST, SKA and Stage IV experiments, we find the uncertainty in the α_X is reduced to O(10^-1) to O(10^-2) <cit.>. An interesting case is that of JBD gravity where we found that current constraints place ω_ BD>10^3; future constraints, using combinations of LSST, SKA and Stage IV will push ω_ BD> few × 10^4, comparable with constraints on the Solar System and astrophysical scales. Once we reach this level, we will truly have entered a new age in gravitational physics, with constraints on cosmological scales playing a significant role in understanding gravity beyond GR. In fact, if the values of μ and γ (or equivalently the current values of the α_X) which current data seem to find discrepant with GR hold up, then there will be uncontrovertible evidence of new physics.

Future surveys will also allow us to greatly improve searches for screening, or fifth forces more generally, in galaxies and other astrophysical structures. For example, with the SKA it will be possible to obtain constrains on Δ G/G≃ 10^-9 using HI-optical offsets of galaxies centres <cit.>.



 §.§ Synergies with other gravitational probes.

With the parallel development of observational cosmology, gravitational wave physics and black hole imaging, interesting (and inevitable) synergies have developed. For a start, and as more and better detections of compact object mergers are recorded, it will become possible to undertake cosmological constraints with such events. A first attempt with GW170817 at constraining the Hubble constant has shown the tremendous effectiveness of this new window on the expansion rate of the Universe <cit.>. The constraint on the speed of gravitational waves has already had a devastating impact on the fauna of modified gravity theories, leading to a severe cull of acceptable models <cit.>. Furthermore, the absence (or presence) of friction in the propagation of gravitational waves can have an impact on what deviations from GR are allowed <cit.>. And the promise of detecting a stochastic background of gravitational waves is leading to renewed scrutiny of novel effects that might arise in deviations from GR.

A particular example of synergy between cosmological and gravitational wave physics has recently been invoked in <cit.> where it was shown that many (primarily scalar-tensor) theories in which the extra fields had a cosmological impact, yet in which the speed of gravitational waves was equivalent to that of light, would obey no-hair theorems to some degree. This means that, for example, any smoking gun for deviations of gravity would have to be on large, cosmological, scales or in the, still dynamical, regime post-merger of binary collision events. I.e. these theories would not manifest themselves during the inspiral or in black hole imaging.

This emerging and exciting synergy is the hallmark of what I believe is a new and exciting phase in gravitational physics, a new "golden age" for general relativity. Until now the focus, and quest, for new physics has been in the high energy regime, of particle accelerators and direct (and indirect) dark matter searches. With the multiple new windows on the gravitational universe, in which cosmology takes pride of place, one would hope that new forces and phenomena are on the verge of discovery. In the very least, we will end up with a cast iron theory for gravity, tested over an enviable range of scales and regimes.






§ DISCLOSURE STATEMENT

The author is not aware of any affiliations, memberships, funding, or financial holdings that
might be perceived as affecting the objectivity of this review. 




§ ACKNOWLEDGMENTS

Much of my understanding of this field was developed with a wide number of collaborators, to all of whom I am grateful. I would like to thank specific input from Alex Barreira, Harry Desmond, Shahab Joudaki, Kazuya Koyama, Martin Kunz, Macarena Lagos, Danielle Leonard, Antony Lewis and Alessandra Silvestri. I am also extremely grateful to Eva-Maria Mueller, Shahab Joudaki and Sukhdeep Singh for producing "bespoke" images for this review. This work was supported by ERC Grant No: 693024, STFC and the Beecroft Trust.






00
natexlab#1#1

ar-style2.bst


[Adelberger, Heckel & Nelson (2003)]Adelberger2003
Adelberger E, Heckel BR, Nelson AE 2003
Ann. Rev. Nucl. Part. Sci 53:77

[Alam et al. (2016)]Alam2016
BOSS collaboration, Alam A, Ata M, Bailey S, Beutler F et al 2016 
M.N.R.A.S. 470:2617

[Allwright & Lehner (2018)]Allwright2018
Allwright G., Lehner L 2018
Submitted ArXiv:1808.07897

[Alonso et al. (2017)]Alonso2017
Alonso D, Bellini E, Ferreira PG, Zumalacarregui M 2017
Phys. Rev. D95, 063502

[Amendola  et al. (2008)]Amendola2008
Amendola L, Kunz M, Sapone D 2008
J.C.A.P. 0804:013

[Amendola et al. (2018)]Amendola2018
Amendola L, Kunz M, Saltas I 2018
J.C.A.P. 1808:030

[Amin et al. (2009)]Amin2009
Amin MA, Wagoner RV, Blandford RD, 2009 
M.N.R.A.S. 390:131

[Amon et al. (2017)]Amon2017
Amon A, Blake C, Heymans C, Leonard CD et al 2017 
M.N.R.A.S. 479:3422


[Amon et al. (2018)]Amon2018
Amon A, Blake C, Heymans C, Leonard CD, Asgari, M et al.  2018
M.N.R.A.S 479:3422

[Arkani-Hamed  et al. (1998)]Arkani1998
Arkani-Hamed N, Dimopoulos S, Dvali G 1998
Phys. Lett. B 429:263

[Asvathaman et al. (2017)]Asvathaman2017
Asvathaman A, Heyl J, Hui L 2017
M.N.R.A.S. 465:3261

[Avilez & Skordis (2014)]Avilez2014
Avilez A, Skordis C 2014
Phys. Rev. Lett. 113:011101

[SDSS (2002)]Bahcall2002
SDSS collaboration, Bahcall N, Dong F, Bode P, Kim R at al  2002
Astroph. Journ. 585:182

[Baker et al. (2012)]Baker2012
Baker T, Ferreira PG, Skordis C 2012
Phys. Rev. D87:024015

[Baker et al. (2014)]Baker2014
Baker T, Ferreira PG, Skordis C 2014
Phys. Rev. D89:024026
[Baker & Bull (2015)]Baker2015
Baker T, Bull P 2015
Astroph. Journ. 811:116

[Baker et al. (2015)]Baker2015a
Baker T, Psaltis D, Skordis C
Astroph. Journ.  802:63



[Baker et al. (2017)]Baker2017
Baker T, Bellini E, Ferreira PG, Lagos M, Noller J, Sawicki I 2017
Phys. Rev. Lett. 119:251301

[Baker et al. (2018)]Baker2018
Baker T, Clampitt J, Jain B, Trodden M 2018
Phys. Rev. D98:023511

[Barreira et al. (2015)]Barreira2015
Barreira A, Cautun M, Li B, Baugh C, Pascolo, S 2015
J.C.A.P. 1508:028

[Barreira et al. (2016)]Barreira2016
Barreira A, Sanchez A, Schmidt F 2016
Phys. Rev. D94:084022

[Battye & Pearson (2011)]Battye2011
Battye RA, Pearson JA 2011 
Phys. Rev. D 84:125032

[Belgacem et al. (2018)]Belgacem2018
Belgacem E, Dirian Y, Foffa S, Maggiore M 2018
Phys. Rev. D97:104066

[Bellini & Sawicki (2014)]Bellini2014
Bellini E, Sawicki I 2014
J.C.A.P. 1407:050

[Bellini et al. (2016)]Bellini2016
Bellini E, Cuesta A, Jimenez R,  Verde L 2016
J.C.A.P. 1602:053

[Bellini et al. (2017)]Bellini2017
Bellini E, Barreira A, Frusciante N, Hu B, Peirone S et al. 2017
Phys. Rev. D97:023520

[Bernardeau et al. (2002)]Bernardeau2002
Bernardeau F, Colombi S, Gaztanaga E, Scoccimarro R 2002
Phys. Rept. 367:1

[Berti et al. (2015)]Berti2015
Berti E, Barrausse E, Cardoso V, Gualtieri L, Pani P et al. 2015
Class. Quant. Grav. 32:243001

[Bertotti et al. (2003)]Bertotti2003
Bertotti B, Iess L, Tortora P 2003
Nature 425:374

[Bertschinger & Zukin (2008)]Bertschinger2008
Bertschinger E, Zukin P 2008 
Phys. Rev. D78:024015

[WiggleZ (2011)]Blake2011
WiggleZ collaboration, Blake C, Kazin E, Beutler F, Davis T, Parkinson D et al. 2011
M.N.R.A.S. 418:1707

[Blake et al. (2016)]Blake2016
Blake C, Joudaki S, Heymans C, Choi A et al 2016 
M.N.R.A.S. 456:2806

[Blas et al. (2011)]Blas2011
Blas D, Lesgourgues J, Tram T 2011
J.C.A.P. 1107:034

[Blas et al. (2011)]Blas2011a
Blas D, Shaposhnikov M, Zenhausern D 2011
Phys. Rev. D84:044001

[Bloomfield et al. (2012)]Bloomfield2012
Bloomfield J, Flanagan E, Park M, Watson S 2012
J.C.A.P. 1308:010

[Bose et al. (2017)]Bose2017
Bose B, Koyama K, Hellwing W, Zhao GB, Winther H  2017
Phys. Rev. D96:023519

[Brans & Dicke (1961)]JBD
Brans C, Dicke RH 1961
Phys. Rev. 124:925

[Broderick et al.(2013)]Broderick2013
Broderick A, Johannsen T, Loeb A, Psaltis D 2013
Astroph. J. 784:7


[Burgess (2003)]Burgess2003
Burgess C 2003
Living Rev. Rel. 7:5

[Burrage & Sakstein (2017)]Burrage2017
Burrage C, Sakstein J 2017
Living. Rev. Rel. 21:1

[Cabre et al. (2012)]Cabre2012
Cabre A, Vikram V, Zhao GB, Jain B, Koyama K 2012
J.C.A.P. 1207:034

[Cai et al. (2014)]Cai2014
Cai YC, Neyrink M, Szapudi I, Cole S, Frenk C  2014
Astroph. Journ. 786:110

[Carrasco et al. (2012)]Carrasco2012
Carrasco JJ, Hertzberg M, Senatore L 2012
J.H.E.P. 1209:082

[Cataneo et al. (2015)]Cataneo2015
Cataneo M, Rapetti D, Schmidt F, Mantz A, Allen S et al. 2015
Phys. Rev. D92:044009

[Cataneo et al (2018)]Cataneo2018
Cattaneo M, Lombriser L, Heymans C et al, 2018 
Submitted to M.N.R.A.S
ArXiv:1812.05594

[Chisari et al. (2018)]Chisari2018
Chisari NE, Richardson M, Devriendt J, Dubois Y, Schneider A et al. 2018
M.N.R.A.S. 480:3962

[Choquet-Bruhat & Geroch (1969)]Choquet1969
Choquet-Bruhat Y, Geroch R  1969
Comm. Matt. Phys 14:329

[Clampitt & Cai (2012)]Clampitt2012
Clampitt J, Cai Y-C 2012
M.N.R.A.S. 431:749C

[Clifton et al. (2012)]Clifton2012
Clifton T, Ferreira PG, Padilla A, Skordis C 2012
Phys. Rep. 513:1

[Cooray & Sheth (2002)]Cooray2002
Cooray A, Sheth R 2002
Phys. Rept. 372:1

[Creminelli et al. (2009)]Creminelli2009
Creminelli P, D'Amico G, Norena J, Vernizzi F 2009
J.C.A.P. 0902:018

[Creminelli & Vernizzi (2017)]Creminelli2017
Creminelli P, Vernizzi F 2017
Phys. Rev. Lett. 119:251301

[Crittenden et al. (2012)]Crittenden2012
Crittenden R, Zhao GB, Pogosian L, Samushia L, Zhang X  2012
J.C.A.P. 1202:048

[Cusin et al. (2018)]Cusin2018
Cusin G, Lewandowksi M, Vernizzi F 2018
J.C.A.P. 1804:061

[BOOMERanG (2000)]deBernardis1999
BOOMERanG Collaboration, de Bernardis P, Ade P, Bock, JJ, Boscaleri A et al. 1999
Nature 404:955

[De Felice et al. (2011)]DeFelice2011
De Felice A, Kobayashi T, Tsujikawa S 2011
Phys. Lett. B706:123

[Deffayet et al. (2009)]Deffayet2009
Deffayet C, Esposito-Farese G, Vikman A 2009
Phys. Rev. D79:084003

[DeLaTorre et al. (2016)]delatorre2016
de la Torre S, Jullo E, Giocoli C, Pezzotta A et al 2016
Astron. Astrophys. 608:44

[ De Rahm et al. (2011)]deRahm2011
de Rahm C, Gabadadze G, Tolley A 2011
Phys. rev. Lett 106:231101

[DES (2018a)]Abbott2018
DES Colaboration, Abbott TC, Abdalla FB, Avila A, Banerji M et al.  2018a
Phys. Rev. D98:043526

[DES (2018b)]Abbott2018DES
DES Colaboration, Abbott TC, Abdalla FB, Avila A, Banerji M et al.  2018b
Submitted to Phys. Rev. D ArXiv:1810:02499

[Deser (1970)]Deser1970
Deser S 1970
Gen. Rel. Grav. 1:9

[Deser & Woodard (2007)]Deser2007
Deser S, Woodard R 2007
Phys. Rev. Lett. 99:111301

[Desjacques et al. (2016)]Desjacques2016
Desjacques V, Jeong D, Schmidt F 2016
Phys. Rept. 733:1

[Desmond et al. (2018a)]Desmond2018
Desmond H, Ferreira PG,Lavaux G, Jasche J  2018a
M.N.R.A.S. 474:3152

[Desmond et al. (2018b)]Desmond2018a
Desmond H, Ferreira PG, Lavaux G, Jasche J 2018b
Submitted ArXiv: 1802.07206

[Desmond et al. (2018c)]Desmond2018b
Desmond H, Ferreira PG, Lavaux G, Jasche J 2018c
Accepted to Phys. Rev. ArXiv: 1807.01486

[Desmond et al. (2018)]Desmond2018d
Desmond H, Ferreira PG, Lavaux G, Jasche J 2018d
Submitted ArXiv: 1807.11742

[Doeleman et al. (2008)]Doeleman2008
Doeleman S, Weintroub J, Rogers A, Plambeck, Freund R et al. 2008
Nature 455:78

[Donoghue (1994)]Donoghue1994
Donoghue J 1994
Phys. Rev. D50:3874

[Dossett et al. (2011)]Dossett2011
Dossett J, Ishak M, Moldenhauer J 2011 
Phys. Rev. D84 123001

[Dvali et al. (2000)]Dvali2000
Dvali G, Gabadadze G, Porrati M 2000
Phys. Lett. B485:208

[Espejo et al. (2018)]Espejo2018
Espejo J, Peirone S, Raveri M, Koyama K, Pogosian et al. 2018
Submitted ArXiv:1809.01121

[Ezquiaga & Zumalacarregui (2017)]Ezquiaga2017
Ezquiaga JM, Zumalacarregui M 2017
Phys. Rev. Lett. 119: 251304




[Fierz & Pauli (1939)]FP1939
Fierz M, Pauli W 1939
Proc. Roy. Soc.  A173:211

[Ferreira et al. (2016)]Ferreira2016
Ferreira PG, Hill C,  Ross G 2016
Phys. Lett. B763:174

[Ferreira & Skordis (2010)]Ferreira2010
Ferreira PG, Skordis C 2010
Phys. Rev. D81:104020

[Feynman et al. (1995)]Feynman1995
Feynman RP, Mornigo FB, Wagner WG 1995
Feynman lectures on gravitation Addison Wesley

[Freire et al. (2012)]Freire2012
Freire PC, Wex N, Esposito-Farese G, Verbiest JPW, Bailes M et al. 2012
M.N.R.A.S 423:3328

[Garcia-Bellido et al. (2011)]GarciaBellido2011
Garcia-Bellido J, Rubio J, Shaposhnikov M, Zenhausern D 2011
Phys. Rev. D84:123504

[Gleyzes et al. (2013)]Gleyzes2013
Gleyzes J, Langlois D, Piazza F, Vernizzi F 2013
J.C.A.P. 1308:025

[Gleyzes et al. (2015a)]Gleyzes2015a
Gleyzes J, Langlois D, Vernizzi F 2015a
Int. J. Mod. Phys D23:1443010

[Gleyzes et al. (2015b)]Gleyzes2015b
Gleyzes J, Langlois D, Piazza F, Vernizzi F 2015b
Phys. Rev. Lett. 114:211101

[Gubitosi et al. (2012)]Gubitosi2012
Gubitosi G, Piazza F, Vernizzi F 2012
J.C.A.P. 1302:032

[Hamaus et al. (2015)]Hamaus2015
Hamaus N, Sutter PM, Lavaux G, Wandelt B 2015
J.C.A.P. 1511:036

[Hamaus et al. (2016)]Hamaus2016
Hamaus N et al. 2016
Phys. Rev. Lett. 117:091302

[MAXIMA (2000)]Hanany2000
MAXIMA Collaboration, Hanany S,  Ade P, Balbi A, Bock J 2000
Astroph. Journ. 545:L5

[Hassan & Rosen (2012)]Hassan2012
Hassan F, Rosen R 2012
J.H.E.P. 1202:126

[Heisenberg (2014)]Heisenberg2014
Heisenberg L 2014
J.C.A.P. 1405:015


[KIDS (2016)]Hildebrandt2016
KIDS Collaboration, Hildebrandt H, Viola M, Heymans C, Joudaki S et al. 2016
M.N.R.A.S. 465:1454

[Hinterbichler & Khoury  (2010)]Hinterbichler2010
Hinterbichler K, Khoury J 2010
Phys. Rev. Lett. 104:231301

[Hinterbichler (2012)]Hinterbichler2012
Hinterbichler K 2012
Rev. Mod. Phys 84:671

[Hojiati et al. (2011)]Hojiati2011
Hojiati A, Pogosian L, Zhao GB 2011
J.C.A.P. 1108:005


[Hojman et al. (1976)]Hojman1976
Hojman SA, Kuchar K, Teitelboim C 1976
Annals. Phys. 96:88

[Horndeski (1974)]Horndeski1974
Horndeski G W 1974
Int. J. Theor. Phys. 10:363

[Hu et al. (2013)]Hu2013
Hu B, Raveri M, Frusciante N, Silvestri A 2013
Phys. Rev. D89:103530

[Hu & Sawicki (2007)]Hu2007
Hu W, Sawicki I 2007 
Phys. Rev. D76:104043


[Hui & Nicolis (2012)]Hui2012
Hui L, Nicolis A 2012
Phys. Rev. Lett. 109:051304

[Jain & Vanderplas (2011)]Jain2011
Jain B, Vanderplas J 2011
J.C.A.P. 10:032

[Jacobson & Mattingley (2000)]Jacobson2000
Jacobson T, Mattingley D 2000
Phys. Rev. D64:024028

[Johannsen & Psaltis (2010)]Johannsen2010
Johannsen T, Psaltis D (2010)
Astroph. Journ. 718:446

[KIDS (2018)]Joudaki2018
KIDS collaboration, Joudaki S, Blake C, Johnson A, Amon A et al.    2018
M.N.R.A.S. 474:4894

[Joyce et al. (2014)]Joyce2014
Joyce A, Jain B, Khoury J, Trodden M 2014
Phys. Rept. 568:1

[Kaiser & Sfakianakis (2014)]Kaiser2014
Kaiser D, Sfakianakis E 2014
Phys. Rev. Lett 112:011302

[Kaiser (1987)]Kaiser1987
Kaiser N 1987
M.N.R.A.S. 227:1

[Khoury & Weltman (2004)]Khoury2004
Khoury J, Weltman A 2004
Phys. Rev. D69:044026

[Koboyashi et al. (2011)]Koboyashi2011
Koboyashi T, Yamaguchi M, Yokoyama J 2011 
Prog. Theor. Phys. 126:511

[Koyama et al. (2009)]Koyama2009
Koyama K, Taruya A, Hiramitsu T 2009
Phys. Rev. D79:123512

[Kreisch & Komatsu (2017)]Kreisch2017
Kreisch CD, Komatsu E 2017
Submitted ArXiv:1712.02710

[Kunz (2012)]Kunz2012
Kunz M 2012
Compt. Rend. Phys. 13:539

[Lagos at al (2016)]Lagos1
Lagos M, Baker T, Ferreira PG, Noller J  2016
J.C.A.P. 1608:007

[Lagos & Ferreira (2016)]Lagos2
Lagos M, Ferreira PG 2016
J.C.A.P. 1701:047

[Lagos at al (2016)]Lagos3
Lagos M, Bellini E, Noller J, Ferreira PG, Baker T  2016
J.C.A.P. 1803:021


[Lehner & Pretorius (2014)]Lehner2014
Lehner L, Pretorius F 2014
Ann. Rev. Astron. Astroph. 552:661

[Leonard et al. (2015)]Leonard2015
Leonard D, Ferreira PG, Heymans C
J.C.A.P. 1512:051

[LIGO (2016a)]Abbott2016a
LIGO and Virgo Scientific Collaborations, Abbott BP, Abbott R, Abbott TD, Abernathy MR et al.  2016a
Phys. Rev. Lett. 116:061102

[LIGO (2016b)]Abbott2016b
LIGO and Virgo Scientific Collaborations, Abbott BP, Abbott R, Abbott TD, Abernathy MR  et al. 2016b
Phys. Rev. Lett. 116:221101

[LIGO (2017a)]Abbott2017
LIGO and Virgo Scientific Collaborations, Abbott BP, Abbott R, Abbott TD, Acernese F al  2017a
Phys. Rev. Lett. 119:161101

[LIGO  (2017b)]Ligo2017
LIGO and Virgo Scientific Collaborations, Abbott BP, Abbott R, Abbott TD, Acernese F al 2017b
Nature 551:85

[Linder et al. (2016)]Linder2016
Linder E, Sengor G, Watson S 2016
J.C.A.P. 1605:053

[Llinares et al. (2014)]Llinares2014
Llinares C, Mota D, Winther H 2014
A.& A. 562:A78

[Leonard et al. (2018)]Leonard2018
Leonard C, Ferreira PG, Heymans C 2018
J.C.A.P. 12:051

[Lewis et al. (1999)]Lewis1999
Lewis A, Challinor A, Lasenby A 1999
Astroph. Jour. 538:473

[Li et al. (2012)]Li2012
Li B, Teyssier R, Koyama K 2012
J.C.A.P. 1:51

[Linder (2005)]Linder2005
Linder E 2005
Phys. Rev. D72:043529

[Lombriser et al. (2012)]Lombriser2012
Lombriser L, Slosar A, Seljak U, Hu W 2012
Phys. Rev. D85:124038

[Lovelock (1972)]Lovelock1972
Lovelock D 1972
J. Math. Phys 13: 874

[Macaulay et al. (2013)]Macaulay2013
Macaulay E, Wehus I, Eriksen HK 2013
Phys. Rev. Lett. 111:161301

[Maggiore & Mancarella (2014)]Maggiore2014
Maggiore M, Mancarella M 2014
Phys. Rev. D90:023005

[Mancarella (2017)]Mancarella2017
Mancarella M 2017
PhD. Thesis 

[HSC (2018)]Mandelbaum2018
HSC collaboration, Mandelbaum R, Miyatake H, Hamana T, Oguri M et al. (HSC collaboration) 2018
PASJ 70:S25

[Melchior et al. (2014)]Melchior2014
Melchior P, Sutter P, Sheldon ES, Krause E, Wandelt, BD et al. 2014
M.N.R.A.S. 440:2922

[Miller et al. (2002)]Miller2002
Miller A, Beach J, Bradley S, Caldwell R, Chapman H et al. 2002
Astroph. Journ. Suppl. 140:115

[Mueller et al. (2018)]Mueller2018
Mueller EM, Percival W, Linder E, Alam S, Zhao GB et al. 2018
M.N.R.A.S. 475:2122

[Nadathur et al. (2015)]Nadathur2015
Nadathur S, Hotchkiss S, Diego JM, Iliev IT, Gottlober S et al. 2015
M.N.R.A.S. 449:3997

[Nicolis et al. (2009)]Nicolis2009
Nicolis A, Rattazzi R, Trincherini E 2009
Phys. Rev. D79:064036

[Peirone et al. (2018)]Peirone2018
Peirone S, Koyama K, Pogosian L, Raveri M, Silvestri A  2018
Phys. Rev. D97:043519

[SCP (1998)]Perlmutter1998
Supernova Cosmology Project (SCP) Collaboration, Perlmutter S, Aldering G, Goldhaber G, Knop, RA et al. 1998
Astroph. Jour. 517:565

[Planck (2014)]Ade2014
Planck Collaboration, Ade P, Aghanim N, Armitage-Caplan C, Arnaud M et al. 2014
A.& A. 571:A20

[Planck (2016)]Ade2016
Planck Collaboration, Ade P, Aghanim N, Arnaud M, Ashdown M et al. 2016
A.&A. 594:A14

[Planck (2018a)]Aghanim2018
Planck Collaboration, Aghanim N, Akrami Y, Ashdown M, Aumont J, Baccigalupi C  et al.  2018a
Submitted  ArXiv:1807.06209

[Planck (2018b)]Akrami2018a
Planck Collaboration,  Akrami Y, Arroja F, Ashdown M, Aumont J, Baccigalupi C. al  2018b
Submitted 1807.06211

[Pogosian & Silvestri (2016)]Pogosian2016
Pogosian L, Silvestri A 2016
Phs. Rev. D94:104014

[Puchwein et al. (2013)]Puchwein2013
Puchwein E, Baldi M, Springel V 2013
M.N.R.A.S. 436:348

[Pullen et al. (2015)]Pullen2015
Pullen AR, Alam S, Ho S 2015 
M.N.R.A.S. 449:4326

[Pullen et al. (2016)]Pullen2016
Pullen A,  Alam S, He S, Ho S 2016
M.N.R.A.S. 460:4098

[Randall & Sundrum (1999)]Randall1999
Randall L, Sundrum R. 1999
Phys. Rev. Lett 83:3370

[Riess et al. (1998)]Riess1998
Riess A, Filippenko A, Challis P, Clocchiatti A et al. 1998
Astron. J 116:1009

[Riess et al. (2016)]Riess2016
Riess A, Macri LM, Hoffmann S, Scolnic D, Casertano S et al. 2016
Astroph. Journ. 826:56

[Renk et al. (2017)]Renk2017
Renk J, Zumalacarregui M, Montanari M, Barreira A 2017
J.C.A.P. 1710:020

[Reyes, R. et al. (2010)]Reyes2010
Reyes R, Mandelbaum R, Seljak U, Baldauf T, Gunn JE et al. 2010
Nature 464:256

[Sakstein et al. (2016)]Sakstein2016
Sakstein J, Wilcox H, Bacon D, Koyama K, Nichol R 2016
J.C.A.P. 1607:019

[Sakstein et al. (2017)]Sakstein2017
Sakstein J, Jain B, Heyl J, Hui L 2017
Astroph. Journ. 844:L14

[Samushia et al. (2014)]Samushia2014
Samushia L, Reid B, White M, Percival WJ, Cuesta AJ et al. 2014
M.N.R.A.S. 439:3504

[Schmidt et al. (2011)]Schmidt2011
Schmidt F, Vihlinin A,  Hu W 2011
Phys. Rev. D83:063503

[Schwarz et al. (2015)]Schwarz2015
Schwarz D, Copi C,  Huterer D, Starkman G 2015
Class. Quant. Grav. 33:184001

[Simpson & Peacock (2010)]Simpson2010
Simpson F, Peacock J 2010
M.N.R.A.S. 402:2009

[Singh et al. (2018)]Singh2018
Singh S, Alam S, Mandelbaum R, Seljak U, Rodriguez-Torres S et al. 2018
Submitted ArXiv:1803.08915

[Spergel et al. (2015)]Spergel2015
Spergel D, Flauger R, Hlozek R 2015
Phys. Rev. D91:023518

[Springel (2005)]Springel2005
Springel V 2005
M.N.R.A.S. 364:1105

[Sotiriou & Faraoni (2010)]Sotiriou2010
Sotiriou T, Faraoni V 2010
Rev. Mod. Phys. 82:451

[Sutter et al. (2012)]Sutter2012
Sutter P, Lavaux G, Wandelt B, Weinberg D 2012
Astroph. Journ. 761:44

[Switzer et al. (2013)]Switzer2013
Switzer E, Masui KW, Bandura K, Calin L-M, Cheng T-C et al. 2013
M.N.R.A.S. 434:L46

[Tattersall et al. (2018)]Tattersall2018
Tattersall O, Ferreira P, Lagos M 2018
Phys. Rev. D97:084005

[Terukina et al. (2014)]Terukina2014
Terukina A, Lombriser L, Yamamoto K, Bacon D, Koyama K et al. 2014
J.C.A.P. 4:013

[Teyssier (2002)]Teyssier2002
Teyssier R 2002
A.& A. 385:337

[Thorne & Will (1971)]Thorne1971
Thorne K, Will C 1971
Astroph. Jour. 163:595

[Vainshtein (1972)]Vainshtein1972
Vainshtein AI 1972
Phys. Lett. 39B:393

[Valogiannis & Bean (2017)]Valogiannis2017
Valogiannis G,  Bean R 2017
Phys. Rev. D95:103515

[Valogiannis & Bean (2018)]Valogiannis2018
Valogiannis G,  Bean R 2018
Phys. Rev. D97:023535

[Wang et al. (2012)]Wang2012
Wang J, Hui L, Khoury J
Phys. Rev. Lett. 109:241301

[Weinberg (1964)]Weinberg1964
Weinberg S 1964
Phys. Rev. 135:B1049

[White (2016)]White2016
White M 2016
J.C.A.P. 1611:057

[Will (2014)]Will2014
Will C
Living. Rev. Rel. 4:4

[Winther et al. (2015)]Winther2015
Winther H, Schmidt F, Barreira A, Arnold C, Bose S et al. 2015
M.N.R.A.S. 454:4208

[Winther et al. (2017)]Winther2017
Winther H, Koyama K, Manera M, Wright B, Zhao GB et al. 2017
J.C.A.P. 1708:006

[Woodard (2007)]Woodard2007
Woodard R 2007
Lect.Notes Phys. 720: 403

[Zhang et al. (2007)]Zhang2007
Zhang P, Liguori M, Bean R,  Dodelson S 2007
Phys. Rev. Lett. 99:141302

[Zlosnik et al. (2007)]Zlosnik2007
Zlosnik T, Ferreira PG, Starkman G 2007
Phys. Rev. D75:044017

[Zumalacarregui & Garcia-Bellido (2014)]Zuma2014
Zumalacarregui M, Garcia-Bellido J 2014
Phys. Rev. D89:064046

[Zumalacarregui et al. (2017)]Zuma2017
Zumalacarregui M, Bellini E, Sawicki I, Lesgourgues J, Ferreira PG 2017
J.C.A.P. 1708:019





