Unveiling the Dynamics of the Universe

Review Unveiling the Dynamics of the Universe Pedro Avelino 1,2 , Tiago Barreiro 3,4 , C. Sofia Carvalho 5,6 , Antonio da Silva 3,7 , Francisco S.N. Lobo 3,7, *, Prado Martín-Moruno 8 , José Pedro Mimoso 3,7 , Nelson J. Nunes 3,7 , Diego Rubiera-García 3 , Diego Sáez-Gómez 3 , Lara Sousa 1 , Ismael Tereno 5,7 and Arlindo Trindade 1 1 2 3 4 5 6 7 8 * Instituto de Astrofísica e Ciências do Espaço, Universidade do Porto, CAUP, Rua das Estrelas, PT4150-762 Porto, Portugal; [email protected] (P.A.); [email protected] (L.S.); [email protected] (A.T.) Departamento de Física e Astronomia, Faculdade de Ciências, Universidade do Porto, Rua do Campo Alegre 687, PT4169-007 Porto, Portugal Instituto de Astrofísica e Ciências do Espaço, Departamento de Física da Faculdade de Ciências da Universidade de Lisboa, Edifício C8, Campo Grande, P-1749-016 Lisbon, Portugal; [email protected] (T.B.); [email protected] (A.S.); [email protected] (J.P.M.); [email protected] (N.J.N.); [email protected] (D.R.-G.); [email protected] (D.S.-G.) Departamento de Matemática, ECEO, Universidade Lusófona de Humanidades e Tecnologias, Campo Grande, 376, 1749-024 Lisboa, Portugal Instituto de Astrofísica e Ciências do Espaço, Universidade de Lisboa, Tapada da Ajuda, 1349–018 Lisboa, Portugal; [email protected] (C.S.C.); [email protected] (I.T.) Research Center for Astronomy and Applied Mathematics, Academy of Athens, Soranou Efessiou 4, 11–527 Athens, Greece Departamento de Física, Faculdade de Ciências da Universidade de Lisboa, Edifício C8, Campo Grande, P-1749-016 Lisboa, Portugal Departamento de Física Teórica I, Universidad Complutense de Madrid, E-28040 Madrid, Spain; [email protected] Correspondence: [email protected]; Tel.: +351-217-500-484; Fax: +351-217-500-977 Academic Editor: Sergei D. Odintsov Received: 13 May 2016; Accepted: 11 July 2016; Published: 26 July 2016 Abstract: We explore the dynamics and evolution of the Universe at early and late times, focusing on both dark energy and extended gravity models and their astrophysical and cosmological consequences. Modified theories of gravity not only provide an alternative explanation for the recent expansion history of the universe, but they also offer a paradigm fundamentally distinct from the simplest dark energy models of cosmic acceleration. In this review, we perform a detailed theoretical and phenomenological analysis of different modified gravity models and investigate their consistency. We also consider the cosmological implications of well motivated physical models of the early universe with a particular emphasis on inflation and topological defects. Astrophysical and cosmological tests over a wide range of scales, from the solar system to the observable horizon, severely restrict the allowed models of the Universe. Here, we review several observational probes—including gravitational lensing, galaxy clusters, cosmic microwave background temperature and polarization, supernova and baryon acoustic oscillations measurements—and their relevance in constraining our cosmological description of the Universe. Keywords: scalar-tensor gravity; modified theories of gravity; inflationary models; topological defects; observational cosmology; cosmic microwave background radiation; weak lensing; spacetime inhomogeneities Symmetry 2016, 8, 70; doi:10.3390/sym8080070 www.mdpi.com/journal/symmetry Symmetry 2016, 8, 70 2 of 79 1. Introduction During the last few decades Cosmology has evolved from being mainly a theoretical area of physics to become a field supported by high precision observational data. Recent experiments call upon state of the art technology in Astronomy and Astrophysics to provide detailed information on the contents and history of the Universe, which has led to the measurement of parameters that describe our Universe with increasing precision. The standard model of Cosmology is remarkably successful in accounting for the observed features of the Universe. However, a number of fundamental open questions remains at the foundation of the standard model. In particular, we lack a fundamental understanding of the recent acceleration of the Universe [1,2]. What is the so-called “dark energy” that is driving the cosmic acceleration? Is it vacuum energy or a dynamical field? Or is the acceleration due to infra-red modifications of Einstein’s theory of General Relativity (GR)? How is structure formation affected in these alternative scenarios? What are the implications of this acceleration for the future of the Universe? The resolution of these fundamental questions is extremely important for theoretical cosmology. Dark energy models are usually assumed to be responsible for the acceleration of the cosmic expansion in most cosmological studies. However, it is clear that these questions involve not only gravity, but also particle physics. String theory provides a synthesis of these two branches of physics and is widely believed to be moving towards a viable quantum gravity theory. One of the key predictions of string theory is the existence of extra spatial dimensions. In the brane-world scenario, motivated by recent developments in string theory, the observed 3-dimensional universe is embedded in a higher-dimensional spacetime [3]. The new degrees of freedom belong to the gravitational sector, and can be responsible for the late-time cosmic acceleration [4,5]. On the other hand, generalizations of the Einstein-Hilbert Lagrangian, including quadratic Lagrangians which involve second order curvature invariants have also been extensively explored [6–10]. These modified theories of gravity not only provide an alternative explanation for the expansion history of the Universe [11–13], but they also offer a paradigm fundamentally distinct from the simplest dark energy models of cosmic acceleration [14], even from those that perfectly mimic the same expansion history. Nevertheless, it has been realized that a large number of modified gravity theories are amenable to a scalar-tensor formulation by means of appropriate metric re-scalings and field redefinitions. It is, therefore, not surprising that we can think of scalar-tensor gravity theories as a first stepping stone to explore modifications of GR. They have the advantage of apparent simplicity and of a long history of examination. First proposed in its present form by Brans and Dicke for a single scalar field [15], they have been extensively generalised and have maintained the interest of researchers until the present date. For instance, an extensive field of work has been developed in the cosmological dynamics of scalar-tensor theories. This can be elegantly summarised by carrying out a unified qualitative analysis of the dynamical system for a single scalar field. We also have an established understanding of the observational bounds in these models, where we can use the Parametrized Post-Newtonian formalism to constrain the model parameters. Finally, with a conformal transformation, these theories can be recast as matter interacting scalar fields in General Relativity. In this format, they can still play an important role in dark energy modelling, such as in coupled quintessence models. The consideration of a multi-scalar fields scenario, which can be perceived as the possible reflection of a multi-scalar tensor gravity theory allows for cooperative effects between the fields yielding assisted quintessence. Indeed, scalar fields are popular building blocks used to construct models of present-day cosmological acceleration. They are appealing because such fields are ubiquitous in theories of high energy physics beyond the standard model and, in particular, are present in theories which include extra spatial dimensions, such as those derived from string theories. Recently, relative to scalar-tensor theory, much work has been invested in the Galileon models and their generalizations [16]. The latter models allow nonlinear derivative interactions of the scalar field in the Lagrangian and lead to second order field equations, thus removing any ghost-like instabilities. The Lagrangian was first written Symmetry 2016, 8, 70 3 of 79 down by Horndeski in 1974 [17], which contains four arbitrary functions of the scalar field and its kinetic energy. The form of the Lagrangian is significantly simplified by requiring specific self-tuning properties (though it still has four arbitrary functions), however, the screening is too effective, and will screen curvature from other matter sources as well as from the vacuum energy [18]. An alternative approach consists of searching for a de Sitter critical point for any kind of material content [19]. These models might alleviate the cosmological constant problem and can deliver a background dynamics compatible with the latest observational data. A promising alternative to explain the late-time cosmic acceleration is to assume that at large scales Einstein’s theory of GR breaks down, and a more general action describes the gravitational field. Thus, one may generalize the Einstein-Hilbert action by including second order curvature invariants such as R2 , Rµν Rµν , Rµναβ Rµναβ , C µναβ Cµναβ , etc. Some of the physical motivations for these modifications of gravity were inspired on effective models raised in string theory, which indeed may lead to the possibility of a more realistic representation of the gravitational fields near curvature singularities [20]. Moreover, the quantization of fields in curved space-times tell us that the high-energy completion of the Einstein-Hilbert Lagrangian of GR involves higher-order terms on the curvature invariants above [21]. This is in agreement with the results provided from the point of view of treating GR as an effective field theory [22]. Among these extensions of GR the so-called f ( R) gravity has drawn much attention over the last years, since it can reproduce late-time acceleration and in spite of containing higher order derivatives, it is free of the Ostrogradsky instability, as can be shown by its equivalence with scalar-tensor theories (for a review on f ( R) gravity see Refs. [6–9]). Moreover, f ( R) gravities have been also proposed as solutions for the inflationary paradigm [23], where the so-called Starobinsky model is a successful proposal, since it satisfies the latest constraints released by Planck [24]. In addition, the equivalence of f ( R) gravities to some class of scalar-tensor theories has provided an extension of the so-called chameleon mechanism to f ( R) gravity, leading to some viable extensions of GR that pass the solar system constraints [25,26]. Other alternative formulations for these extensions of GR have been considered in the literature, namely, the Palatini formalism, where metric and affine connection are regarded as independent degrees of freedom, which yields an interesting phenomenology for Cosmology [27]; and the metric-affine formalism, where the matter part of the action now depends and is varied with respect to the connection [28]. Recently, a novel approach to modified theories of gravity was proposed that consists of adding to the Einstein-Hilbert Lagrangian an f (R) term constructed a la Palatini [29]. It was shown that the theory can pass the Solar System observational constraints even if the scalar field is very light. This implies the existence of a long-range scalar field, which is able to modify the cosmological and galactic dynamics, but leaves the Solar System unaffected. Note that these modified theories of gravity are focussed on extensions of the curvature-based Einstein-Hilbert action. Nevertheless, one could equally well modify gravity starting from its torsion-based formulation and, in particular, from the Teleparallel Equivalent of General Relativity (TEGR) [30]. The interesting point is that although GR is completely equivalent with TEGR at the level of equations, their modifications (for instance f ( R) and f ( T ) gravities, where T is the torsion) are not equivalent and they correspond to different classes of gravitational modifications. Hence, f ( T ) gravity has novel and interesting cosmological implications, capable of describing inflation, the late-time acceleration, large scale structure, bouncing solutions, non-minimal couplings to matter, etc. [31–33]. Another gravitational modification that has recently attracted much interest is the massive gravity paradigm, where instead of introducing new scalar degrees of freedom, such as in f ( R) gravity, it modifies the graviton itself. Massive gravity is a well-defined theoretical problem on its own and has important cosmological motivations, namely, if gravity is massive it will be weaker at large scales and thus one can obtain the late-time cosmic acceleration. Fierz and Pauli presented the first linear massive gravity. However, it was shown to suffer from the van Dam-Veltman-Zakharov (vDVZ) discontinuity [34,35], namely the massless limit of the results do not yield the massless Symmetry 2016, 8, 70 4 of 79 theory, namely, GR. The incorporation of nonlinear terms cured the problem but introduced the Boulware-Deser (BD) ghost. This fundamental problem puzzled physicists until recently, where a specific nonlinear extension of massive gravity was proposed by de Rham, Gabadadze and Tolley (dRGT), in which the BD ghost is eliminated by a Hamiltonian constraint [16]. This new nonlinear massive gravity has interesting cosmological implications, for instance, it can give rise to inflation, late-time acceleration [16]. However, the basic versions of this theory exhibit instabilities at the perturbative level, and thus suitable extensions are necessary. These could be anisotropic versions, f ( R) extensions, bigravity generalizations, partially-massive constructions. The crucial issue is whether one can construct a massive gravity and cosmology that can be consistent as an alternative to dark energy or other models of modified gravity, and whether this theory is in agreement with high-precision cosmological data, such as the growth-index or the tensor-to-scalar ratio, remains to be explored in detail. Quantum field theory predicts that the universe underwent, in its early stages, a series of symmetry breaking phase transitions, some of which may have led to the formation of topological defects. Different types of defects may be formed depending on the (non-trivial) topology of the vacuum manifold or the type of symmetry being broken. For instance, Domain Walls—which are surfaces that separate domains with different vacuum expectation values—may arise due to the breaking of a discrete symmetry, whenever the vacuum manifold is disconnected. Line-like defects, or Cosmic strings, are formed if the vacuum is not simply connected or, equivalently, if it contains unshrinkable loops. This type of vacuum manifold results, in general, from the breaking of an axial symmetry. Moreover, if the vacuum manifold contains unshrinkable surfaces, the field might develop non-trivial configurations corresponding to point-like defects, known as Monopoles. The spontaneous symmetry breaking of more complex symmetry groups may lead to the formation of textures, delocalized topological defects which are unstable to collapse. The production of topological defects networks as remnants of symmetry breaking phase transitions is thus predicted in several grand unified scenarios and in several models of inflation. Moreover, recent developments in the braneworld realization of string theory suggest that its fundamental objects—p-dimensional D-branes and fundamental strings—may play the cosmological role of topological defects. Topological defects networks, although formed in the early universe, may in most instances survive throughout the cosmological history and leave a variety of imprints on different observational probes. The observational consequences of topological defect networks can be very diverse, depending both on the type of defects formed and on the evolution of the universe after they are generated. Although the possibility of a significant contribution to the dark energy budget has been ruled out both dynamically [36] and observationally [37], light domain walls may leave behind interesting astrophysical and cosmological signatures. For instance, they may be associated to spatial variations of the fundamental couplings of nature (see, e.g., [38]). On the other hand, cosmic strings may contribute significantly to small-scale cosmological perturbations and have consequently been suggested to have significant impact on the formation of ultracompact minihalos [39], globular clusters [40], super-massive black holes [41] and to provide a significant contribution to the reionization history of the Universe [42]. Both cosmic strings and domain walls may be responsible for significant contributions to two of the most significant observational probes: the temperature and polarization anisotropies of the cosmic microwave background and the stochastic gravitational wave background. This fact—alongside the possibility of testing string theory through the study of topological defects—greatly motivates the interest on the astrophysical and cosmological signatures of topological defects. An extremely important aspect of modern cosmology is the synergy between theory and observations. Dark energy models and modified gravity affect the geometry of the universe and cosmological structure formation, impacting the background expansion and leaving an imprint on the statistical properties of the large-scale structure. There are a number of well-established probes of cosmic evolution, such as type Ia supernovae, baryon acoustic oscillations (BAO), weak Symmetry 2016, 8, 70 5 of 79 gravitational lensing, galaxy clustering and galaxy clusters properties [43]. Different methods measure different observables, probing expansion and structure formation in different and often complementary ways and have different systematic effects. In particular, joint analyses with Cosmic Microwave Background (CMB) data are helpful in breaking degeneracies by constraining the standard cosmological parameters. Indeed, CMB has revolutionized the way we perceive the Universe. The information encoded in its temperature and polarization maps provides one of the strongest evidences in favour of the hot Big-Bang theory and has enabled ways to constrain cosmological models with unprecedented accuracy [44]. The CMB also encodes additional information about the growth of cosmological structure and the dynamics of the Universe through the secondary CMB anisotropies. These are originated by physical effects acting on the CMB after decoupling [45], such as the integrated Sachs-Wolfe effect and the Sunyaev-Zel’dovich (SZ) effect, manifest respectively on the largest and arc-minute scales of the CMB. In this review, we will discuss in some detail both a well-established acceleration probe (weak lensing) and a few promising ones related to galaxy cluster properties and the SZ effect. Galaxy clusters are the largest gravitationally bound objects in the Universe and are among the latest bound structures forming in the Universe. For this reason, their number density is highly sensitive to the details of structure formation as well as to cosmological background parameters and cluster abundance is a well established cosmological probe. The Sunyaev-Zel’dovich effect [46,47] is the scattering of CMB photons by electrons in hot reservoirs of ionized gas in the Universe, such as galaxy clusters. In particular, SZ galaxy cluster counts, profiles, scaling relations, angular power spectra and induced spectral distortions are promising probes to confront model predictions with observations. Weak gravitational lensing [48] describes the deflection of light by gravity in the weak regime. Its angular power spectrum is a direct measure of the statistical properties of gravity and matter on cosmological scales. Weak lensing, together with galaxy clustering, is the core method of the forthcoming Euclid mission to map the dark universe. Euclid [49] will provide us with weak lensing measurements of unprecedented precision. To obtain high-precision and high-accuracy constraints on dark energy or modified gravity properties, with both weak lensing and SZ clusters, non-linearities on structure formation must be taken into account. While linear covariant perturbations equations may be evolved with Boltzmann codes, non-linearities require dedicated N-body cosmological simulations [50]. There currently exist a number of simulations for various modified gravity and dark energy models, together with a set of formulas that fit a non-linear power spectrum from a linear one. Hydrodynamic simulations, commonly used in cluster studies, are increasingly needed in weak lensing applications to model various baryonic effects on lensing observables, such as supernova and AGN feedback, star formation or radiative cooling [51]. The ΛCDM framework provides a very good fit to various datasets, but it contains some open issues [52]. As an example, there are inconsistencies between probes, such as the tension between CMB primary signal (Planck) and weak lensing (CFHTLenS) [53], as well as problems with the interpretation of large-scale CMB measurements (the so-called CMB anomalies) [54]. Alternatives to ΛCDM or deviations to General Relativity are usually confronted with data using one of two approaches: model selection or parameterizations. In model selection a specific model is analyzed and its parameters constrained. Such analyses have a narrower scope but may be better physically motivated. Parameterizations are good working tools and are helpful in highlighting a particular feature and in ruling out larger classes of models, however they must be carefully defined in a consistent way. Parameterizations are commonly applied to the dark energy equation of state and to deviations from General Relativity. An example of the latter is the gravitational slip, which provides an unambiguous signature of modified gravity, and can be estimated combining weak lensing measurements of the lensing potential with galaxy clustering measurements of the Newtonian potential. Model parameters in both of the approaches are usually estimated with Monte Carlo techniques, while the viability of different models may be compared using various information criteria. Symmetry 2016, 8, 70 6 of 79 Besides model testing, cosmological data is also useful to test foundational assumptions, such as the (statistical) cosmological principle and the inflationary paradigm. The common understanding is that cosmological structures are the result of primordial density fluctuations that grew under gravitational instability collapse. These primordial density perturbations would have originated during the inflationary phase in the early universe. Most single field slow-roll inflationary models produce nearly Gaussian distributed perturbations, with very weak possible deviations from Gaussianity at a level beyond detection [55–57]. However, non-Gaussianities may arise in inflationary models in which any of the conditions leading to the slow-roll dynamics fail [58], such as the curvaton scenario [59–61], the ekpyrotic inflacionary scenario [62,63], vector field populated inflation [64–66] and multi-field inflation [67–70]. Tests of non-Gaussianity are thus a way to discriminate between inflation models and to test the different proposed mechanisms for the generation of primordial density perturbations. Likewise, the assumption of statistical homogeneity may be tested. Locally, matter is distributed according to a pattern of alternate overdense regions and underdense regions. Since averaging inhomogeneities in the matter density distribution yields a homogeneous description of the Universe, then the apparent homogeneity of the cosmological parameters could also result from the averaging of inhomogeneities in the cosmological parameters, which would reflect the inhomogeneities in the density distribution. The theoretical setup closest to this reasoning is that of backreaction models, where the angular variations in the parameters could also source a repulsive force and potentially emulate cosmic acceleration. Hence the reasoning is to look for these inhomogeneities, not in depth, but rather across the sky and then to use an adequate toy model to compute the magnitude of the acceleration derived from angular variations of the parameters compared to the acceleration driven by a cosmological constant [71]. This work will focus on all of the above-mentioned topics. More specifically, this article is outlined in the following manner: In Section 2, we present scalar-tensor theories, and in Section 3, we consider Horndeski theories and the self-tuning properties. In Section 4, f ( R) modified theories of gravity and extensions are reviewed. In Section 5, an extensive review on topological defects is carried out. The following sections are dedicated to observational cosmology. In particular, in Section 6, cosmological tests with galaxy clusters at CMB frequencies are presented and in Section 7, gravitational lensing will be explored. In Section 8, the angular distribution of cosmological parameters as a measurement of spacetime inhomogeneities will be presented. Finally, in Section 9 we conclude. 2. Scalar-tensor Theories Scalar-tensor theories of gravity stem from the original proposal of Brans-Dicke theory (BD) [15], although theories involving a scalar field with a gravitational role in addition to the tensor metric fields can be traced back to the earlier Kaluza-Klein type theories upon dimensional reduction [72]. One notorious example was Jordan’s proposal for a field theory realization of Dirac’s Large Number Hypothesis [73]. In the course of time, it has been realized that a large number of modified gravity theories are amenable to a scalar-tensor formulation [74] by means of appropriate metric re-scalings and field redefinitions. Interestingly, this possibility of a scalar-tensor description occurs regardless of the field equations being derived through the metric or, alternatively, the Palatini variational prescription [29]. For instance, it is well known that the f ( R) higher-order gravity theories can be cast as a scalar-tensor theory [74–76]. The same happens for the theories based on Lagrangeans with non-linear kinetic terms of a scalar field dubbed k-essence [77]. Symmetry 2016, 8, 70 7 of 79 2.1. General Formalism The general scalar-tensor class of gravity theories encompassing BD’s and a plethora of modified theories that we just alluded to, were first examined by Bergmann [78], Wagoner [79], and Nordvedt [80] (BWN). Their fundamental action can be cast in the following form [81] S= Z p  ω (φ) ,µ − g d x φR − φ,µ φ − 2φλ(φ) + S M , φ 4  (1) where R is the usual Ricci curvature scalar of the spacetime, φ is a scalar field, ω (φ) is a dimensionless function of φ which calibrates the coupling between the scalar field and gravity (in what follows we shall refer to it as the coupling parameter), λ(φ) is another function of the scalar field which can be interpreted both as a potential for φ and as a cosmological parameter and, finally, S M represents the action for the matter fields. It is worth mentioning at this point that the action should include an additional term involving the extrinsic curvature, as in GR [82], to account for a boundary term associated with the variation of Rαβ with respect to gαβ . A discussion of this term in the scalar-tensor theories is provided in [74]. It becomes immediately clear that the scalar field φ plays the role which in Einstein’s GR is confined to the gravitational constant, and has dimensions of a mass squared, φ ∼ 1/G. Since φ is now a dynamical variable, we thus see that the trademark of these theories is the fact that they exhibit a varying gravitational “constant”. In the archetypal Brans-Dicke theory the coupling ω (φ) is constrained to be a constant. The consideration of couplings between a scalar field and the spacetime geometry envisaged by this class of theories allows a generalization to multiple scalar fields so that we may consider multi-scalar-tensor theories. The Lagrangian can then be cast as L g = f (φµ ) R − ωµν (φλ ) αβ g ∇α φµ ∇ β φν + V (φµ ) + 16π Lm , f (φµ ) (2) where R is the usual Ricci curvature scalar of a space-time endowed with the metric gαβ , φµ (µ = 1, . . . , n) represents n scalar fields, f (φµ ) is a function of those scalar fields, ωµν (φ β ) is the metric of the internal space of the φµ fields, V (φµ ) is the scalar potential. The Lagrangian of the matter fields Lm does not carry any explicit dependence on φµ which guarantees that free-falling test particles follow the geodesics of the spacetime (this amounts to keeping valid the conservation of the energy-momentum tensor ∇ β T αβ = 0). We leave aside a detailed consideration of this extended class of theories, which have been addressed in [83,84]. However, in subsection 2.4 we shall analyse a model of quintessence stemming from the cooperative behaviour of several scalar fields coupled to several matter components which as it will be argued can be interpreted as a particular realisation of the models characterised by Equation (2). Taking the variational derivatives of the action (1) with respect to the two dynamical variables gαβ and φ , yields the field equations 1 Rαβ − gαβ R + λ(φ) gαβ 2 φ − =   1 ω (φ) γ ∇α φ∇ β φ − gαβ ∇γ φ∇ φ 2 φ2  Tαβ 1  + ∇ β ∇α φ − gαβ ∇γ φ∇γ φ + 8πGN , φ φ   1 2φ2 λ′ (φ) − 2φλ(φ) 8πGN T − ω ′ (φ)∇γ φ∇γ φ , = 2ω (φ) + 3 2ω (φ) + 3 (3) (4) where T ≡ T γ γ is the trace of the energy-momentum tensor, Tα β , of the matter content of spacetime, and GN is the gravitational constant normalized to its present value (in what follows in this section, we set GN = 1 ). Symmetry 2016, 8, 70 8 of 79 From these equations the role of λ(φ) becomes more apparent. As pointed out in [78], this term enters Equation (3) as a cosmological function, but intervenes in the second Equation (4) as a mass term. A consequence of this latter feature is that the scalar waves propagate with a speed different from the speed of the electromagnetic waves. Another interesting aspect to point out in connection with this term is that the inverse square law for the static gravitational field is modified in the presence of a non-vanishing λ(φ) . This manifests itself in that a Yukawa-type term arises [79,81]. It is important to notice two additional points. First, we assume that the usual relation αβ ∇ β T = 0 establishing the conservation laws satisfied by the matter fields holds. This is ensured by requiring that the matter fields Lagrangian Lm be independent of φ [78] The role of the scalar field is then that of helping to generate the spacetime curvature associated with the metric. Matter may create this field, but the latter cannot act back directly on the matter, which thus responds only to the metric [85]. The second point follows from Equation (4). By allowing the coupling parameter ω (φ) to tend to infinity, it can be seen that the right-hand side of the equation may vanish. This will be the case when ω → ∞ dominates over ω ′ (φ) , λ(φ) and λ′ (φ) [86]. Indeed, in the absence of λ(φ), one requires ω ′ /w3 → 0 [80]. A solution with constant φ is then asymptotically admitted by the Klein-Gordon equation, and thus GR is recovered in the ω → ∞ limit of scalar-tensor theories (see also [83]). However, it is important, on the one hand, to separate the cases where the trace of the energy-momentum tensor is trivially vanishing, i.e., vacuum and radiation. In the latter cases, and again when there is no cosmological potential λ(φ), it has been argued that there are asymptotic behaviors that differ from GR [83,86–88]. On the other hand, in the presence of the potential λ(φ), the asymptotic behaviour is very much dependent on it. In Ref. [86], it was shown that when the potential is of a power-law nature, in particular, when φ λ(φ) ∝ φ2 , it governs the approach to GR, superseding the mechanism due to the divergence of ω (φ). Morevover, as shown below, BD theory can also be an attractor [89]. 2.2. Conformal Picture A convenient device to deal with modified gravity theories is to make use of an appropriate conformal transformation [74,75,82] to bring the theory into the form of GR plus a minimally coupled field. Consider the transformation φ g̃αβ = g , (5) φ0 αβ where φ0 is a constant with the dimensions of mass squared which we introduce to make the scaling factor dimensionless. In order to keep the next results as simple as possible, we shall assume that φ0 = 1 and that φ in the following expressions is dimensionless (i.e., whenever we write φ it corresponds to φ/φ0 where φ0 is some value of reference, for instance φ0 = 1/G). Using the transformed metric g̃αβ and redefining the scalar field as ϕ≡ Z r 2ω (φ) + 3 dφ , 2 φ (6) we find that the action adopts the form S̃ = Z d4 x p Z p   ˜ α ϕ∇ ˜ α ϕ − 2V ( ϕ) + 16π − g̃ R̃ − ∇ d4 x − g̃ L̃m , (7) where V ( ϕ) = U (φ( ϕ)) /φ2 ( ϕ), and also L̃m = φ−2 Lm . Thus, we realize that the action of the STT becomes the action of GR with matter and a massive scalar field. Writing the Einstein field equations in the transformed frame (dubbed the Einstein frame), we obtain R̃αβ −   1 ˜ γ ϕ∇ ˜ γ ϕ − V ( ϕ) g̃αβ , ˜ α ϕ∇ ˜ β ϕ − 1 g̃αβ ∇ g̃αβ R̃ = 8π T̃αβ + ∇ 2 2 (8) Symmetry 2016, 8, 70 9 of 79 r d 2 ˜ ϕ − V ( ϕ) = 8π T̃, dϕ 2ω + 3   r 2 1 φ;α β ˜ ∇ β T̃ α = − . T̃ 2 φ 2ω + 3 (9) (10) It becomes immediately apparent that, when the trace of the energy-momentum tensor of the matter fields in the original frame is vanishing, there is an equivalence between the two pictures related via the conformal transformation. This happens trivially for a vacuum model, and also in the radiation case. In all other cases, the GR-like frame exhibits a coupling between matter and the redefined scalar field. One is then drawn to conclude that the equivalence principle is violated. However, since this variation is induced by the presence of the scalar field which also acts as a source of the Einstein Equation (8) the latter conclusion seems arguable and has been much debated [90–92]. If one considers the scalar field as a component of the total energy-momentum tensor on the right-hand side of the field equations, then the realization that the original matter does not exhibit the geodesic behaviour is a mere consequence of the fact that it couples to the other source component. At this stage we should mention that the conformal transformation mixes the geometric and matter degrees of freedom, which results in many interpretational ambiguities [93]. One may also mention that if one only restricts attention to the Einstein frame, one may also lose sight of the original motivations and modifications of gravity in the geometrical sector [94,95]. However, the Einstein and Jordan frames are physically equivalent, as can be traced back to Dicke’s original paper [96], where the conformal transformation technique was introduced. In fact, both conformal frames are equivalent, in the spirit of Dicke’s paper, provided that in the Einstein frame the units of mass, time and space scale as appropriate powers of the scalar field, and are thus varying. We will discuss the issues of discriminating the geometric and matter sectors in Section 4.1.2. 2.3. Cosmological Dynamics Here we consider the dynamics of cosmological models in general scalar-tensor gravity theories. A number of important, exact solutions of Brans-Dicke models [97,98], and of more general class of theories [99–101] have been derived in the literature. In what follows we resort to a unified qualitative analysis of the major dynamical features of these scalar field cosmological models [86,89,102–105] which has the advantage of relying on a closed, autonomous system of equations where the functional dependence on φ of both the coupling ω (φ) and cosmological potential λ(φ) is left unspecified. The analysis of the relevant features of the dynamics, namely, the set of possible asymptotic behaviors, leads to a classification of the classes of couplings and potentials vis-à-vis to their cosmological impact. Alternative dynamical system analyses that are found in the literature usually envisage specific classes of scalar-tensor theories [106–108] or, when more general [109–112], are not so transparent as the treatment that follows. In order to investigate the cosmological dynamics, it is most convenient to introduce dimensionless expansion normalized variables akin to the square roots of density parameters. This associates to the equilibrium points of the dynamical system (also dubbed singular, critical, or fixed points in the jargon of the theory of dynamical systems) non-static solutions of the cosmological evolution. An alternative qualitative approach based on the bare dimensional variables of the cosmological model, namely, on the scale factor of the universe, Hubble factor and/or on the energy densities of the various source components is also possible, but is of more limited scope in the scrutiny of non-static asymptotic regimes. It is however required whenever the coupling function ω (φ) or the cosmological potential λ(φ) have vanishing minima. Symmetry 2016, 8, 70 10 of 79 In what follows, we shall restrict to the homogeneous and isotropic universes given by the Friedmann-Robertson-Walker (FRW) metric 2 2 2 ds = −dt + a (t)   dr2 2 2 2 2 + r (dθ + sin θ dφ ) , 1 − k r2 (11) where k = 0, ±1 distinguishes the curvature of the spatial hypersurfaces. In the Einstein frame the field equations with a matter content given by a perfect fluid with p = (γ − 1)ρ can be cast as [86] 3 k ã˙ 2 +3 2 ã2 ã 3ã˙ ϕ̇ ϕ̈ + ã ϕ̇2 + V ( ϕ) + ã−3γ m̃( ϕ) 2 dm̃( ϕ) dV ( ϕ) − ã−3γ − , dϕ dϕ 8π Φ∗ = =   (12) (13) √ where ã = Φ/Φ∗ a, the overdots stand for the derivatives with respect to the conformally √ transformed time d̃t̃ = Φ/Φ∗ dt, V ( ϕ) = U (Φ( ϕ))/(Φ( ϕ)/Φ∗ )2 , m̃( ϕ) = µ0 (Φ( ϕ)/Φ∗ )(3γ−4)/2 , where µ0 = 8πρ0 sets the arbitrary initial condition for the energy density of the perfect fluid, and, √ more importantly, ∂ ϕ ln m( ϕ) ∝ (4 − 3γ)α( ϕ) where α = ( 2Ω + 3)−1 is the PPN function that translates the coupling between the scalar field and matter in the Einstein frame as used in Damour and Nordtvedt [83]. These equations are those of GR with a scalar field subject both to a time-independent potential, V ( ϕ), and to a time-dependent potential m̃( ϕ) a−3γ illustrating the fact that in the Einstein frame the masses of particles vary. (Hereafter we shall set 8π/φ∗ = 1). Introducing the new time variable N = ln a and the dimensionless density parameters x2 = ϕ̇2 , 6H 2 y2 = V ( ϕ) , 3H 2 (14) as well as the expansion normalized curvature term K = k/( aH )2 , the Einstein field equations become a fourth order, autonomous dynamical system x ′ = y′ = K′ = ϕ′ = r  i ∂ϕV 3 h y2 + x 2x2 + γ (1 − x2 − y2 + K ) −(3 + K ) x − V 2 r  3 ∂ϕm − (1 − x 2 − y2 + K ) 2 m r   i 3 ∂ϕV 3 h xy + y 2x2 + γ (1 − x2 − y2 + K ) − yK , 2 V 2  h i 3 2 2 2 2x + γ (1 − x − y + K ) − 2K2 , −2K 1 − 2 √ 6x , 3 2   (15) (16) (17) (18) where we have used ρ/3H 2 = Ωm = 1 − x2 − y2 + K. GR models correspond to the case where ∂ ϕ ln m( ϕ) = 0, and Brans-Dicke models are characterised by an exponential coupling m, i.e., by ∂ ϕ ln m( ϕ) = α0 (as well as by V ( ϕ) = 0 in BD original version, but we do not require it here). The crucial point regarding the qualitative study of general models with scalar fields lies in the ϕ′ -equation, since it allows the consideration of arbitrary choices of V ( ϕ) and of m( ϕ) [102]. We compactify the phase space of Equations (15)–(17) by considering ϕ ∈ ℜ ∪ {∞}, and distinguish the fixed points arising at finite values of ϕ, which require x = 0, from those at ϕ = ∞ (which we shall denote ϕ∞ ) that require either x = 0 or ψ = ϕ−1 = 0. The K = 0 and the y = 0 subspaces are invariant manifolds. Symmetry 2016, 8, 70 11 of 79 For K = 0 the system reduces to the Equations (15), (16) and (18) and we find the following fixed points, namely, at finite ϕ: • x = 0, y = 1, ϕ0 where ∂ ϕ V ( ϕ0 ) = 0. This case corresponds to de Sitter solutions, dominated by the scalar field, and arise at non-vanishing maxima or minima of the potential V. Notice that ∂ ϕ m may be different from 0. These solutions are attractors at minima of V and repellors at maxima of V. • x = y = 0, ϕ0 where ∂ ϕ m = 0. The second case corresponds to matter dominated solutions with V = 0 at maxima or minima of m. Note however that the system (15)–(18) is singular when V = 0, translating the fact that the variables in use are not regular. In the original variables, this second class of solutions exists provided that V has also an extremum, or is identically zero. At ϕ∞ , the finite ϕ solutions may exist also at ϕ = ∞ provided that V or m are asymptotically flat. However, other solutions may show up at ϕ = ∞, if both V and m are asymptotically exponential. In fact, q in this case the equilibrium q in ψ = 1/ϕ is trivially satisfied at ψ = 0 [102]. Defining ∂ V ∂ m ϕ ϕ 3 3 The exponential asymptotic behaviour of W = 2 lim ϕ→∞ V and Z = 2 lim ϕ→∞ m . m( ϕ) amounts to having a late-time Brans-Dicke behaviour. One finds the following fixed points (see Figures 1 and 2) [113]: ± ± V = ±1, yV = 0, that lie on the intersection of the invariant lines x2 + y2 = 1 and y = 0. • x∞ These points correspond to the vacuum solutions of Brans-Dicke theory found by [98,99]. BM = −W/3, which lies on the invariant line x2 + y2 = 1, provided |W | < 3. These solutions • x∞ (BM) correspond to those found in [114]. N = −2Z/ (3(2 − γ )), which lies on the invariant line y = 0 and exists if | Z | < 3(2 − γ ) /2. • x∞ These solutions (N) correspond to the matter dominated solutions found by [97,99]. S S = 3γ/2 , ( x S )2 + ( yS )2 = 3x∞ + Z , that lies in the interior of the phase space domain, • x∞ ∞ ∞ Z −W Z −W provided that 9γ/2 + Z ( Z − W ) S 2 < 1. 0 < ( x∞ ) + (yS∞ )2 = ( Z − W )2 These are scaling solutions (S) (see [89,113,115] and references therein) Y   BM S X   V- N V+ Figure 1. The figure represents the fixed points at the ϕ∞ phase-plane for V and m asymptotically exponential. Symmetry 2016, 8, 70 12 of 79 W S 4 V- N 2 BM -1 -2 1 BM Z   2 -2 V+ N -4 S Figure 2. The plot exhibits the regions in parameter space that correspond to the asymptotic Brans-Dicke behaviour. The horizontal axis is the Z-axis and the vertical axis is the W-axis. The regions S and BM are separated by the line 9γ/2 − W 2 + WZ = 0, and the separations between the S and N regions is given by the line 2Z2 − 2WZ + 9γ(2 − γ)/2 = 0. For the K 6= 0 case, as K ′ = 0 for K = 0, we see from the linearization of the system that ∂K ′ = − 2 + 2 Q K =0 , ∂K Q K =0 ≡ i 3 h 2 2x + γ (1 − x2 − y2 ) , 2 (19) which shows that the K = 0 subspace is stable whenever Q < 1 which is precisely the condition for inflationary behaviour. Thus inflation is required to guarantee both the stability of the K = 0 asymptotic solutions and the attraction to GR, whenever the latter applies. As previously referred, when the potential V (φ) yields vanishing minima in the Einstein frame the expansion normalized variables adopted here are not suitable, as they yield a dynamical system that is not defined at those minima. In such a case it is preferable to use variables such as the scale factor, the kinetic energy of the scalar field and its potential, as done in [86], or in some of the references in [106,109]. Summing up the qualitative analysis of the dynamical system (15)–(17) associates exact solutions to a classification of the fixed points, and also allows the consideration of how extended gravity theories dynamically relate to GR [89]. It reveals the interplay of the two main mechanisms that yield general relativity as a cosmological attractor of scalar-tensor gravity theories: (i) the vanishing of the coupling function α, and (ii) the existence of a minimum of the scalar field potential V. The latter mechanism is shown to supersede the former one. Moreover the approach to GR is then characterized by a de Sitter inflationary behavior. This guarantees the stability of the k = 0 GR solutions. We have also shown that at ϕ → ∞ there is a variety of power-law solutions that correspond to Brans-Dicke behavior and that can be attractors of the cosmological dynamics. In what follows we will address a multifield scenario, which can be understood as a particular case of the multi-scalar tensor gravity models of Equation (2). Symmetry 2016, 8, 70 13 of 79 2.4. Assisted Quintessence Assisted inflation is a model of early universe inflation where multiple scalar fields cooperate to fuel the accelerated inflationary expansion [116–118]. Scalar potentials that are not viable to provide inflation as a single field can be seen to work when several such fields are present in parallel. A similar idea can be implemented for our present observed accelerated expansion [119–122]. In the usual quintessence scenario, a single scalar field, uncoupled from dark matter, is responsible for the dark energy component in our universe [123]. In assisted quintessence, we have a scenario where the dark energy sector is composed of multiple scalar fields, with different scalar potentials. Since so little is known about the structure of the dark sector, there is the possibility of a coupling between the scalar fields responsible for inflation and dark matter. Coupled quintessence, where a scalar field interacts with dark matter was introduced in Refs. [124–126]. In [127,128] it was further suggested that the scalar field can couple in a different way to several separate dark matter species. To keep the model as general as possible, we assume the dark matter sector to consist of several matter components with different couplings to the dark energy fields [129]. 2.4.1. General Equations We consider an ensemble of n scalar fields φi cross-coupled to an ensemble of m dark matter components ρα . The cross-couplings are described by the matrix Ciα where indexes i, j identify scalar field indexes and indexes α, β identify the dark matter components. The equation of motion for the fields and the various dark matter components, in a spatially flat Friedmann-Robertson-Walker metric with scale factor a(t), are then written as = κ ∑ Ciα ρα , φ̈i + 3H φ̇i + V,φi (20) α = −κ ∑ Ciα φ̇i ρα . ρ̇α + 3Hρα (21) i Notice that the ensemble of n scalar fields φi cross-coupled to m dark matter components ρα can be perceived as a particular case of the multi-scalar-tensor models of Equation (2). Indeed, if we set the metric of the internal space of the scalar fields to be flat, i.e., ωµν = δµν and allow for a combination of dark matter components in the matter sector, the Ciα coupling coefficients are related to the couplings αi = ∂ ln f /∂φi of the multi-scalar-tensor theory in the Einstein frame, and we obtain the latter equations. The solution for the dark matter component evolution can be given immediately in terms of the values of the fields as ! ρα = ρα 0 exp −3N − κ ∑ Ciα (φi − φi 0 ) . (22) i The rate of change of the Hubble function is κ2 Ḣ = − 2 ! , (23) ! , (24) ∑ ρα + ∑ φ̇i2 α i which is subject to the Friedmann constraint κ2 H = 3 2 where ρφ = ∑i φ̇i2 /2 + V (φ1 , ..., φn ). ∑ ρα + ∑ ρφ α i Symmetry 2016, 8, 70 14 of 79 It is possible to define an effective equation of state parameter for the dark sector, weff = pφ , ρm + ρφ (25) where pφ = ∑i φ̇i2 /2 − V (φ1 , ..., φn ) is the total dark energy pressure and ρm = ∑α ρα the total dark matter energy density. We are interested in looking at the critical points of the evolution to look for viable attractors of late time evolution. For this, we will consider two possible forms for the scalar potential, both leading to scaling solutions: a sum of exponential terms with V1 (φ1 , ..., φn ) = M4 ∑i e−κλi φi and an exponential of a sum of terms with V2 (φ1 , ..., φn ) = M4 e− ∑i κλi φi . 2.4.2. Scalar Field Dominated Solution We start our investigation with the case of a negligible matter contribution (ρα = 0). We can show that the effective equation of state for the attractor is 1 weff = −1 + λeff , 3 (26) where λeff can be interpreted as an effective logarithmic slope for an dynamically equivalent single field exponential potential. For the first potential, V1 , we get that λeff has to obey 1 = λ2eff 1 ∑ λ2 . i (27) i In this case incre asing the number of fields has the effect of decreasing the value of λeff , so that an accelerated expansion can be achieved even in the case of each single field having a slope too large to fuel acceleration. For the second potential, V2 , λeff is given by λ2eff = ∑ λ2i . (28) i Here we see that increasing the number of fields actually increases the value of λeff , making an accelerated expansion more difficult. Both of these results are similar to what is obtained for assisted inflation [116–118]. This was to be expected since the matter sector is being neglected. 2.4.3. Scaling Solution Here we will look for fixed points where all the sectors, kinetic and potential energy densities for the scalar energy and the matter energy density, have non-vanishing contributions. We can again define an effective scalar potential logarithmic slope, λeff . But now we also have to take into account the coupling of the scalar fields to the matter sector. This can also be represented by an effective coupling Ceff . In terms of these parameters, we can obtain the effective equation of state Ceff , λeff − Ceff (29) 2 3 − λeff Ceff + Ceff . (λeff − Ceff )2 (30) weff = and for the dark energy sector Ωφ = Symmetry 2016, 8, 70 15 of 79 For the first potential, the value of λeff is again obtained from Equation (27), whereas the value of Ceff is n C jα , (31) Ceff = λeff ∑ λ j j =1 where α can be taken to be any of the matter species. This is because the sum is constrained to yield the same result for all the matter species, otherwise there is no exact scaling solution. It is instructive to consider these quantities for the simple case of all the potentials being replicas, that is, having the same slopes and couplings. Then λi = λ1 and the couplings C are a diagonal matrix with Cii = C11 . We then get weff = C11 , λ1 − C11 Ωφ = 2 3n − λ1 C11 + C11 , 2 ) (λ1 − C11 (32) where n is the number of scalar fields. The equation of state obtained is the usual result for a single coupled field, but in the scalar field energy we can see a mild dependence in the number of fields. For the second potential, the effective slope λeff obtained in scalar dominance is no longer valid. Instead it is possible to show that we can combine the fields through a rotation so that the kinetic energy of all fields is zero except for the first field. In terms of this rotation matrix Q we get the effective equation of state as λeff = ∑ Q1j λ j , (33) j and the effective coupling to be Ceff = ∑ Q1j Cj1 . (34) j The precise form of Q is hard to get in the general case. But the simple case of a diagonal coupling matrix C can give us an idea of what is to be expected. In this instance, we get Q1i = 1 Cii q ∑ j 1/Cll2 , (35) and so λeff = Ceff ∑ i λi , Cii 1 = Ceff 1 ∑ C2 . i (36) ii Finally, if we consider again the case of n fields all with the same slope and couplings, we get √ √ λeff = n λ1 and Ceff = C11 / n. Here we have a strong dependence in the number of fields in the effective equation of state and the dark energy density, weff = C11 , nλ1 − C11 Ωφ = 2 3n − C11 λ1 n + C11 . (nλ1 − C11 )2 (37) This is clearly in contrast to what we have for the first potential type. 2.4.4. Scalar Potential Independent Solutions Finally we look for solutions where the scalar potential energy is negligible. These solutions, therefore, are independent of the choice we make for the scalar potential. Note that the scalar field still contributes to the overall energy density through its kinetic energy. These solutions can be grouped into three different types: • First, we have the kinetic dominated solutions where all the energy density for the dark sector is due to the kinetic energy of the scalar fields. In this instance, we have weff = Ωφ = 1. Symmetry 2016, 8, 70 16 of 79 • As a second possibility, we have a contribution from the kinetic energy and the matter energy 2 . This will be true for any densities for the dark sector. In this case we have weff = Ωφ = 32 ∑i Ciα matter species α since the sum is constrained to yield the same for all matter species. • Finally, as a third possibility we have a matter dominated scenario. This has been investigated previously in the literature for the case of a single field with two matter components [128,130–132]. If we extend this to more than one field, we will see that there is a consistency condition imposed on the values of the couplings to allow a flat direction in the field space, otherwise this solution will be non existent. This is by its nature a transient solution since the matter contribution will eventually decay away and the scaling solution or the scalar field dominated solution will take over. However, this can be an interesting scenario for model building since it allows the dark matter sector to have a significant contribution before the acceleration of the Universe. If we extend our analysis to the linear perturbations in these models we can show that there is a term proportional to ∑i Ciα Ciβ that can act as a growth or damping term to the density contrast of the matter components. This yields a strong constraint on the types of models that can be compatible with observations [129]. 3. Horndeski Theories and Self-tuning In the previous sections we have considered the simplest scalar-tensor theories of gravity, which are given by the action (1). Those theories can be re-written as a scalar field, which is coupled to the geometry through a term f (φ) R in the action, with a canonical kinetic term and a potential. This is: S= Z p − g d4 x   1 f (φ) R − X + V (φ) + S M , 2 (38) where X = ∇µ φ∇µ φ. Nevertheless, there are more things in heaven and earth for scalar-tensor theories of gravity. K-essence models [133] assume a scalar field minimally coupled to gravity ( f (φ) = 1) with an action that is not necessarily the sum of a kinetic term and a potential but a more general function K (φ, X ). One can go further and consider interaction terms in the Lagrangian containing second order derivatives of the scalar field, as in the kinetic braiding models [134] that have a scalar field Lagrangian given by K (φ, X ) + G (φ, X )φ. These models can be stable because their field equations are only second order. However, theories with Lagrangians containing second order derivatives generically have field equations which contain higher than second order derivatives. Therefore, such theories usually propagate an extra (ghostly) degree of freedom which means that they are affected by the Ostrogradski instability [135]. In fact, in 1974, Horndeski found the most general scalar-tensor action leading to second order equations of motion [17]. Nevertheless, Horndeski work passed mostly unnoticed until Deffayet et al. [136] rediscovered his theory when generalizing the covariantized version [137] of the galileons models [138]. The Horndeski Lagrangian can be written as [17] LH =  4 κ1 (φ, X ) ∇µ ∇α φ R βγ νσ − κ1,X (φ, X ) ∇µ ∇α φ∇ν ∇ β φ∇σ ∇γ φ 3 i µ νσ + κ3 (φ, X ) ∇α φ∇ φ R βγ − 4κ3,X (φ, X ) ∇α φ∇µ φ∇ν ∇ β φ∇σ ∇γ φ h i αβ +δµν F (φ, X ) Rαβ µν − 4F,X (φ, X ) ∇µ ∇α φ∇ν ∇ β φ + 2κ8 (φ, X ) ∇α φ∇µ φ∇ν ∇ β φ   −3 2F,φ (φ, X ) + X κ8 (φ, X ) ∇µ ∇µ φ + κ9 (φ, X ) , αβγ δµνσ (39) where κi (φ, X ) and F (φ, X ) are arbitrary functions, satisfying F,X = κ1,φ − κ3 − 2Xκ3,X . Brans–Dicke theory, k-essence, kinetic braiding and many other models can be seen as particular cases of the most general Horndeski family. Thus, although the Horndeski Lagrangian restricts the number of stable scalar-tensor theories. Actually it has been shown that there are still stable theories beyond the Symmetry 2016, 8, 70 17 of 79 Horndeski Lagrangian [139–141]. Those theories have second order equations of motion only when one particular gauge is fixed. On the other hand, the vacuum energy also gravitates in modified theories of gravity. Thus, the cosmological constant problem [142–144] is still present if, as Weinberg suggested, the extra degree of freedom could screen only a given value of that constant [142]. In this context, the fab four models are based in the observation by Charmousis et al. [18,145] that Weinberg’s no-go theorem can be avoided by relaxing one of the assumptions, that is allowing the field to have a non-trivial temporal dependence once the cosmological constant has been screened. Their dynamical screening is based in requiring Minkowski to be a critical point of the dynamics and, when the critical point is an attractor, it may alleviate the cosmological constant problem. Nevertheless, these models are forced by construction to decelerate its expansion while approaching the Minkowski final state. Thus, a late time accelerating cosmology does not naturally arise in this scenario. When considering the compatibility of dynamical screening the vacuum energy with a late time phase of accelerated expansion, the concept of self-adjustment was extended to non-Minkowskian final states [19]. As we will now show in detail, a scalar field self-tuning to de Sitter can lead to very promising scenarios from a phenomenological point of view [146,147]. In addition, it may alleviate the cosmological constant problem if the field is able to screen the vacuum energy before the material content in a particular model. 3.1. Dynamical Screening In a cosmological context, one can consider a FLRW geometry to express the Horndeski Lagrangian in the minisuperspace. Once the dependence on higher derivatives is integrated by parts the point-like Lagrangian takes the simple form [18] L ( a, ȧ, φ, φ̇) = a3 Zi ( a, φ, φ̇) H i , ∑ where L = V −1 i =0..3 Z d3 x L H , (40) V is the spatial integral of the volume element, H = ȧ/a is the Hubble expansion rate, and an over-dot represents a derivative with respect the cosmic time t. The functions Zi are given by Zi ( a, φ, φ̇) = Xi (φ, φ̇) − k Yi (φ, φ̇) , a2 (41) where Xi and Yi are written in terms of the Horndeski free functions [18]. Taking into account Equation (40), the Hamiltonian density can be written as 1 H ( a, ȧ, φ, φ̇) = 3 a   h i ∂L ∂L φ̇ − L = ∑ (i − 1) Zi + Zi,φ̇ φ̇ H i . ȧ + ∂ ȧ ∂φ̇ i =0..3 (42) Assuming that the matter content, described by ρm ( a), is minimally coupled with the geometry and non-interacting with the field, the modified Friedmann equation is given by H ( a, ȧ, φ, φ̇) = −ρm ( a). (43) Now, one can follow a similar argument as that presented in reference [18] for screening to Minkwoski, the only difference being that we require self-tuning to a more general given on-shell 2 6 = 0 (where the subscript “on” means evaluated on-shell-in-a), as done in solution with H 2 = Hon reference [19]. Thus, assuming that the field is continuous in any phase transition that changes the value of the vacuum energy and noting that the screening of any value of the vacuum energy would lead to the screening of any material content, one has to require Hon to be approached dynamical in order to have acceptable cosmological solutions. That is, we want Hon to be an attractor solution, Symmetry 2016, 8, 70 18 of 79 although this particular adjustment mechanism only ensures that it is a critical point. As it is subtly discussed in reference [18], this requirement leads to the following three conditions: • The field equation evaluated at the critical point has to be trivially satisfied leading the value of the field free to screen. This implies that the minisuperspace Lagrangian density at the critical point takes the form i Zi ( aon , φ, φ̇) Hon = c( aon ) + ∑ i =0..3 1 ζ̇ ( aon , φ) . a3on (44) • The modified Freedman equation evaluated once screening has taken place has to depend on φ̇ to absorb possible discontinuities of the cosmological constant appearing on the r. h. s. of this equation. Thus, taking into account Equations (43) and (44) (and its φ̇-derivative), it leads to ∑ i =1..3 i i Zi,φ̇ ( aon , φ, φ̇) Hon 6= 0. (45) • The full scalar equation of motion must depend on ä to allow for a non-trivial cosmological dynamics before screening. This leads to a condition equivalent to (45) if Hon 6= 0. A particular Lagrangian which satisfies conditions can be written as L ( a, ȧ, φ, φ̇) = a 3 " c( a) + ∑ i =1..3  i Zi ( a, φ, φ̇) H − i Hon  # . (46) This Lagrangian (which is a trivial generalization of that presented in reference [18]) has a critical point at the solution Hon by construction [19]. In order to obtian some information about the Zi ’s, we note that, as it has been explicitly proven in [18] (it must be noted that such a proof is independent of the particular Hon ), two Horndeski theories which self-tune to Hon are related through a total derivative of a function µ ( a, φ). Therefore, one can consider L ( a, ȧ, φ, φ̇) = L ( a, ȧ, φ, φ̇) + µ̇ ( a, φ) , (47) with L ( a, ȧ, φ, φ̇) given by Equation (46). As this equation has to be valid during the whole evolution, we can equate equal powers of H to obtain Z0 = c( a) − ∑ i Zi Hon + i =1..3 φ̇ µ,φ , a3 Z1 = Z1 + 1 µ,a , a2 Z2 = Z2 , Z3 = Z3 , (48) which can be combined to yield [18,19] i Zi ( a, φ, φ̇) Hon = c( a) + φ̇ Hon µ,a ( a, φ) + 3 µ,φ ( a, φ) . 2 a a (49) 2 =Λ 3.2. Requiring the Existence of a de Sitter Critical Point: Hon In the first place, we consider k = 0. In this case the dependence on the scale factor a of the Zi ’s vanishes, as can be easily noted from Equation (41), and the point-like Lagrangian (40) is independent of Yi ’s. For k = 0, therefore, Equation (49) is equal to the Lon ( aon ) /a3on . In the second place, we 2 = Λ, this leads to [19] demand Hon ∑ i =0..3 Zi (φ, φ̇) Λ i/2 = c( a) + √ Λ a2 µ,a + φ̇ µ,φ . a3 (50) Symmetry 2016, 8, 70 19 of 79 As the l. h. s. of this equation is independent of a, the r. h. s. should also be independent of a for any value of φ̇, which allows us to obtain µ ( a, φ). Thus, we have Xi (φ, φ̇) Λi/2 = ∑ i =0..3 √ 1 Λ h(φ) + φ̇ hφ (φ). L a = 3 ( ) on on a3on (51) Therefore, there are three different kind of terms which can appear in the Lagrangian. These are: (i) Xi -terms linear on φ̇; (ii) Xi -terms with non-linear dependence on φ̇, which contribution has to vanish in the on-shell point-like Lagrangian; and, (iii) terms that are not able to self-tune because they contribute through total derivatives or terms multiplied by k in the Lagrangian. Thus, the point-like Lagrangian takes the form [19] L ( a, ȧ, φ, φ̇) = Llinear ( a, ȧ, φ, φ̇) + Lnl ( a, ȧ, φ, φ̇) + d µ( a, φ) + k G ( a, ȧ, φ, φ̇), dt (52) with Llinear and Lnl to be given in the following discussion. 3.2.1. Linear Terms: “The Magnificent Seven” In order to satisfy Equation (51) considering only terms linear on φ̇, it is enough to have √ Xims (φ, φ̇) = 3 Λ Ui (φ) + φ̇ Wi (φ), (53) with the potentials Ui and Wi fulfilling the constraint Wi (φ)Λi/2 = ∑ i =0..3 ∑ Ui,φ (φ)Λi/2 . (54) i =0..3 The Lagrangian of the linear terms is Llinear = a3 ∑ i =0..3 i h √ 3 Λ Ui (φ) + φ̇ Wi (φ) H i , (55) with the potentials satisfying condition (54). As there are a total of eight functions and only one constraint, there are effectively only seven free functions which we coined “the magnificent seven”. The field equation for these models is ε linear   √   H Ḣ = ∑ 3 Λ Ui,φ (φ) − √ Wi (φ) − i Wi (φ) H i , H Λ i =0..3 (56) where we have divided all the equations by a3 for simplicity and we require that the term in the square bracket does not vanish. On the other hand, the Hamiltonian density is Hlinear = ∑ i =0..3 h i √ 3(i − 1) Λ Ui (φ) + i φ̇ Wi (φ) H i . (57) When only this kind of terms are present the modified Friedmann equation is Hlinear = −ρm ( a). These models are able to screen any material content dynamically as long as there is at least one Wi 6= 0 with i 6= 0. The models with only W0 and Ui potentials and those with only Ui ’s would not spoil the screening of other models when combined with them, but they are not able to self-tune by themselves [19]. It is worthy to note that a non-self-tuning Einstein–Hilbert√term is contained in these models and can be explicitly written redefining U2 (φ) as U2 (φ) − 1/(8πG Λ). Symmetry 2016, 8, 70 20 of 79 3.2.2. Non-linear Terms We now consider terms with an arbitrary dependence on φ and φ̇. minisuperspace Lagrangian given by Lnl = a3 These terms have a 3 ∑ Xinl (φ, φ̇) H i . (58) i =0 The functions Xinl are, however, not arbitrary. Taking into account Equation (51), the contribution of these terms has to vanish when Lnl is evaluated at the critical point. Therefore, they have to satisfy the condition 3 ∑ Xinl (φ, φ̇) Λi/2 = 0. (59) i =0 The Hamiltonian density of this models is 3 Hnl = ∑ i =0 h i (i − 1) Xinl (φ, φ̇) + φ̇ Xi,nlφ̇ (φ, φ̇) H i . (60) Taking into account condition (59) and its derivatives, it can be seen that Hnl, os depends on φ̇; thus, these models are always able to self-tune. If only these nonlinear terms are present, the modified Friedmann equation is Hnl = −ρm ( a), and Hlinear + Hnl = −ρm ( a) otherwise. On the other hand, taking into account Lagrangian (58) the field equation can be written as 3 ε nl = ∑ i =0  nl Xi,φ − 3Xi,nlφ̇ H − iXi,nlφ̇  Ḣ nl nl − Xi,φ̇φ φ̇ − Xi,φ̇φ̇ φ̈ H i , H (61) with ε linear + ε nl = 0. 3.3. Cosmology of the Models 3.3.1. Example of a Linear Model The field equation and the Friedmann equation read [146] H′ = 3 ∑i H i √ Λ Ui,φ (φ) − H Wi (φ) ∑i i H i Wi (φ)  , φ′ = √ (1 − Ω) H 2 − 3 ∑i (i − 1) H i Ui (φ) Λ , ∑i i H i+1 Wi (φ) where a prime represents derivative with respect to ln a. As an example of a linear model let us consider the three potentials U2 , U3 and W2 . The constraint equation imposes U2,φ Λ + U3,φ Λ3/2 = W2 Λ, and then U2,φ H′ = −3 H W2 For H ≫ √ Λ, √ ! Λ . 1− H H′ 3 U2,φ , =− H 2 W2 therefore we need: U2,φ /W2 = 1, during a matter domination epoch, and U2,φ /W2 = 4/3, for a radiation domination epoch. For example, the choice of the potentials, U2 = eλφ + 34 e βφ , and W2 = λeλφ + βe βφ , give √ us the desired behaviour, as shown in Figure 3. The de Sitter evolution is attained when H → Λ. Symmetry 2016, 8, 70 21 of 79 ln ( κρ m ) ln ( κρ r ) 30 ln ( κρ ϕ ) 20 10 0 - 10 -5 -4 -3 -2 -1 0 1 - log 10 ( 1+ z ) Figure 3. Energy densities evolution for the models of Section 3.3.1. Unfortunately, the contribution of the field at early times is too large to satisfy current constraints. 3.3.2. Example of a Non-Linear Model We will restrict the analysis to the shift-symmetric cases, which means no dependence on φ, and make use of the convenient redefinition ψ = φ̇. Under these assumptions we obtain the equations of motion [147], H′ = 3(1 + w) Q0 P1 − Q1 P0 , Q1 P2 − Q2 P1 ψ′ = 3(1 + w) Q0 P2 − Q2 P0 , Q2 P1 − Q1 P2 where Q0 , Q1 , Q2 , P0 , P1 , P2 , are non-trivial functions of Xi and H, and the average equation of state parameter of matter fluids is ∑ Ω s (1 + w s ) . 1+w = s ∑s Ωs Let us consider a case involving the three potentials X0 , X1 and X2 such that X2 (ψ) = αψn , X1 (ψ) = −αψn + β , ψm X0 ( ψ ) = − β . ψm We can obtain a model with wψ = w0 + wa (1 − a), such that, w0 = −0.98 and wa = 0.04, which is compatible with current limits and moreover, has a negligible dark energy contribution at early times. The evolution of the energy densities is illustrated in Figure 4. 50 ln (κρ 40 m ) ln (κρ r ) 30 ln (κρ 20 ϕ) 10 0 - 10 - 20 -6 -4 -2 0 -log 10 (1+z) Figure 4. Energy densities evolution for the model with non-vanishing X0 , X1 and X2 . Symmetry 2016, 8, 70 22 of 79 3.4. Summary In this section we have considered a subclass of the Horndeski Lagrangian that leads to a late-time de Sitter evolution of the Universe and that may provide a mechanism to alleviate the cosmological constant problem. We have presented examples of a linear and a non-linear model. In particular, the shift-symmetric non-linear models, the de Sitter critical point is indeed an attractor. We can, thus, understand the current accelerated expansion of our Universe as the result of the dynamical approach of the field to the critical point. The model presented satisfies current observational bounds on the background evolution and appear very promising. 4. f ( R) Modified Theories of Gravity and Extensions Nowadays, a great multidisciplinary effort has firmly established experimentally the standard ΛCDM model of cosmology, which is based on Einstein’s General Theory of Relativity (GR), a cold dark matter source and a tiny cosmological constant. GR itself has been tested on scales ranging from the submilimeter to the solar system [81]. Nonetheless, a number of conceptual issues still refuses to naturally merge into this paradigm. Since the establishment of the theorems by Hawking and Penrose [148–150], space-time singularities both deep inside black holes and in the early universe have been a disturbing issue from a theoretical point of view since they represent regions where Physics looses predictability. On the other hand, introduction of a dark sector in the model allows to fit the data, but implies that up to 96% of the total energy-matter content of the universe is of unknown nature. Though several candidates for the microscopic origin of the dark matter component have been considered, like supersymmetric partners [151], this new physics has not been directly detected at particle accelerators so far. Moreover, the accumulated data in favour of a late-time accelerated expansion of the universe only hints at an exotic form of energy whose equation of state is nearly equal to that of a cosmological constant. In view of these difficulties one may wonder whether this is a signal that the conceptual and mathematical framework represented by GR is reaching its limit of validity and that the increasing complexity of the additional ad hoc hypothesis needed to fit the data could be removed by finding a new framework for gravity where the difficulties (singularities and dark sources) are naturally encompassed. At the more fundamental level, great theoretical progress has been made within the unifying frameworks of both string theory [152] and loop quantum gravity [153] but, given the high energy/curvature scales at which the new effects would take place, they are troubled with large difficulties to make contact with experimental results. In this sense, many authors have focused their attention on phenomenological models where specific issues of the ΛCDM model (singularities, inflation, dark sources, etc) are dealt with. In the high-energy regime, these modified gravity models including higher-order powers of the curvature invariants (see e.g., [8] for a review) are supported by a high-energy completion of the Einstein-Hilbert Lagrangian of GR within the theory of quantized fields in curved space-times [21,154]. 4.1. f ( R) Modified Theories of Gravity 4.1.1. General Formalism Recently, f ( R) modified theories of gravity have been extensively explored in the literature, where the standard Einstein-Hilbert action is replaced by an arbitrary function of the Ricci scalar R [155,156]. These theories can be formulated through different approaches: (i) one could use the metric formalism, which consists in varying the action with respect to the metric gµν only; (ii) in the Palatini formalism [28,157–163], where the metric and the affine connection are treated as independent variables (see Section 4.3 below for more details on this formulation); (iii) and the metric-affine formalism, where the matter part of the action now depends on and is varied with respect to the connection [28]. Symmetry 2016, 8, 70 23 of 79 The action for the f ( R) modified theories of gravity is given by S= 1 2κ Z d4 x p − g f ( R) + S M ( gµν , ψ) , (62) R √ where κ = 8πG. S M ( gµν , ψ) is the matter action, defined as S M = d4 x − g Lm ( gµν , ψ), where Lm is the matter Lagrangian density, in which matter is minimally coupled to the metric gµν and ψ collectively denotes the matter fields. Using the metric approach, by varying the action with respect to gµν yields the following field equation 1 (m) FRµν − f gµν − ∇µ ∇ν F + gµν ∇α ∇α F = κ Tµν , (63) 2 (m) where F ≡ d f /dR. The matter energy-momentum tensor, Tµν , is defined as (m) Tµν = − √ √ 2 δ( − g Lm ) . − g δ( gµν ) (64) Now, considering the contraction of Equation (63), provides the following relationship FR − 2 f + 3 ∇α ∇α F = κ T , (65) which shows that the Ricci scalar is a fully dynamical degree of freedom. f ( R) gravity may be written as a scalar-tensor theory, by introducing a transformation { R, f } → {φ, V } defined as φ ≡ F ( R) , V (φ) ≡ R (φ) F − f ( R (φ)) , (66) so that the action (62) can be written as S= 1 2κ Z [φR − V (φ) + Lm ] p − g d4 x. (67) Note that action is simply a Brans-Dicke type action (1) with the parameter ω = 0. The gravitational field equations of f ( R) gravity, in the scalar-tensor representation, simply reduce to Equations (3) and (4), with a redefinition of the potential V (φ) = 2φλ(φ). 4.1.2. Discriminating between Dark Energy and Modified Gravity Models Due to the fact that subsets of scalar-tensor theories can be obtained from modifications to the gravitational sector through specific mappings, it is fundamental to understand how one may differentiate these modified theories of gravity from dark energy models. This ambiguity requires a practical classification, which is not a trivial issue as one cannot discriminate between the dark energy models and modified gravity, solely from the expansion rate of the Universe. However, as structure formation behaves differently by these two classes of theories, information on the growth of structure, at different scales and redshifts, will break the degeneracy and will serve to discriminate between both models of dark energy and modified gravity [164]. In this context, generic modifications of the dynamics of scalar perturbations, with respect to the ΛCDM background, can be represented by two new degrees of freedom in the Einstein constraint equations, namely, through the functions Q( a, k) and η ( a, k ), where a is the scale factor and k the perturbation scale. In modified theories of gravity, Q( a, k ) results from a mass-screening effect due to local modifications of gravity, and effectively modifies Newton’s constant. In the context of dynamical dark energy models, the function Q( a, k) incorporates the interaction with other fields, due to the perturbations. The function η, absent in ΛCDM, parameterizes the effective stresses due to the modification of gravity or specific dynamical dark energy models. Finally, the scale and time-dependence of both functions, Q and η, can be derived in the considered model and traced on a ( Q, η ) plane. Symmetry 2016, 8, 70 24 of 79 Following [164], a practical manner is to denote the term “modified gravity” if additional contributions to the Poisson equation are presented, which induces Q 6= 1, and if extra effective stresses arise, implying η 6= 1. In this context, “modified gravity” is related to models in which modifications are present in the gravitational sector and in which dark energy clusters or interacts with other fields. Following this practical classification [164], in the context of first order perturbation theory, models with Q = η = 1 are denoted standard dark energy models, such as, a minimally-coupled scalar field with standard kinetic energy [164]. On the other hand, models for Q 6= 1 and η 6= 1 are denoted “modified gravity”, such as scalar-tensor theories, f ( R) gravity, massive gravity and generalized galileons, Horndeski interactions, bi-(multi-) gravity, etc. Thus, in the context of the EUCLID mission [164], the definitions of the functions Q and η are extremely convenient, for instance, EUCLID can distinguish between standard dynamical dark energy and modified gravity by forecasting the errors on Q and η, and several combinations of these functions, such as Q/η. 4.1.3. Late-Time Cosmic Acceleration In this subsection, we show that f ( R) gravity may lead to an effective dark energy, without the need to introduce a negative pressure ideal fluid. Consider the FLRW metric (11). Taking into account the perfect fluid description for matter, we verify that the gravitational field equation, Equation (63), provides the generalised Friedmann equations in the following form [165,166]:  2     ȧ ȧ 1 1 ′ − ṘF ( R) [ f ( R) − RF ( R)] − 3 a 3F ( R) 2 a     ä 1 ȧ ′ 1 ′ 2 ′′ ṘF ( R) + R̈F ( R) + Ṙ F ( R) − [ f ( R) − RF ( R)] + a 2F ( R) a 3 = = κ ρ, (68) 3 κ − (ρ + 3p) . (69) 6 These modified Friedmann field equations may be rewritten in a more familiar form, as  2 ȧ a   ä a = = κ ρtot , 3 κ − (ρtot + 3ptot ) , 6 (70) (71) where ρtot = ρ + ρ(c) and ptot = p + p(c) , and the curvature stress-energy components, ρ(c) and p(c) , are defined as     1 1 ȧ ρ(c) = ṘF ′ ( R) , (72) [ f ( R) − RF ( R)] − 3 κF ( R) 2 a     1 ȧ 1 ṘF ′ ( R) + R̈F ′ ( R) + Ṙ2 F ′′ ( R) − [ f ( R) − RF ( R)] , 2 (73) p(c) = κF ( R) a 2 respectively. The late-time cosmic acceleration is achieved if the condition ρtot + 3ptot < 0 is obeyed, which follows from Equation (71). For simplicity, consider the absence of matter, ρ = p = 0. Now, taking into account the equation of state ωeff = p(c) /ρ(c) , with f ( R) ∝ Rn and a generic power law a(t) = a0 (t/t0 )α [165], the parameters ωeff and α are given by ωeff = − 6n2 − 7n − 1 , 6n2 − 9n + 3 α= −2n2 + 3n − 1 , n−2 (74) respectively, for n 6= 1. Note that a suitable choice of n can lead to the desired value of ωeff < −1/3, achieving the late-time cosmic acceleration. It is interesting to note that observations might slightly favour the fact that the effective dark energy equation-of-state parameter lies in the phantom regime, i.e., ωeff < −1 [37], although the data Symmetry 2016, 8, 70 25 of 79 are still far from being conclusive. In fact, recent fits to supernovae, CMB and weak gravitational lensing data indicate that an evolving equation of state crossing the phantom divide, is mildly favored, and several models have been proposed in the literature. In particular, models considering a redshift dependent equation of state, possibly provide better fits to the most recent and reliable SN Ia supernovae Gold dataset. In a cosmological setting, it has been shown that the transition into the phantom regime, for a single field is probably physically implausible [167], so that a mixture of various interacting non-ideal fluids is necessary. In this context, we refer the reader to [168] for a recent review on specific models in modified gravity. If confirmed in the future, this behavior has important implications for theoretical models of dark energy and modified gravity. For instance, this implies that dark energy is dynamical and excludes the cosmological constant and the models with a constant parameter, such as quintessence and phantom models, as possible candidates for dark energy. All of these models present an extremely fascinating aspect for future experiments focussing on supernovae, cosmic microwave background radiation and weak gravitational lensing and for future theoretical research. 4.2. Curvature-Matter Couplings in f ( R) Gravity Recently, an extension of f ( R) gravity was presented, by considering an explicit curvature-matter coupling. The latter coupling induces a non-vanishing covariant derivative of the energy-momentum tensor, ∇µ T µν 6= 0, which potentially leads to a deviation from geodesic motion, and consequently the appearance of an extra force [169]. Implications, for instance, for stellar equilibrium have been studied in Ref. [170]. The equivalence with scalar-tensor theories with two scalar fields has been considered in Ref. [171], and a viability stability criterion was also analyzed in Ref. [172]. It is interesting to note that nonlinear couplings of matter with gravity were analyzed in the context of the accelerated expansion of the Universe [173–175], and in the study of the cosmological constant problem [176]. The action for f ( R) curvature-matter coupling [169] takes the following form S= Z  1 f ( R) + [1 + λ f 2 ( R)] Lm 2 1  p − g d4 x , (75) where f i ( R) (with i = 1, 2) are arbitrary functions of the curvature scalar R. For notational simplicity we consider κ = 1 throughout this subsection. Varying the action with respect to the metric gµν yields the field equations, given by F1 Rµν − 1 f gµν − ∇µ ∇ν F1 + gµν ∇α ∇α F1 = −2λF2 Lm Rµν 2 1 (m) +2λ(∇µ ∇ν − gµν ∇α ∇α )Lm F2 + (1 + λ f 2 ) Tµν , (76) where we have denoted Fi ( R) = f i′ ( R), and the prime represents the derivative with respect to the scalar curvature. We refer the reader to [177] for the curvature-matter coupling in the Palatini formalism. As mentioned above, an interesting feature of this theory is the non-conservation of the energy-momentum tensor. More specifically, by taking into account the generalized Bianchi identities [169,178,179], one deduces the following conservation equation (m) ∇µ Tµν = i λF2 h (m) gµν Lm − Tµν ∇µ R , 1 + λ f2 (77) where the coupling between the matter and the higher derivative curvature terms describes an exchange of energy and momentum between both. (m) Consider the equation of state for a perfect fluid Tµν = (ρ + p) Uµ Uν + pgµν , where ρ is the energy density and p, the pressure, respectively. The four-velocity, Uµ , satisfies the conditions Symmetry 2016, 8, 70 26 of 79 Uµ U µ = −1 and U µ Uµ;ν = 0. Introducing the projection operator hµν = gµν + Uµ Uν , gives rise to non-geodesic motion governed by the following equation of motion for a fluid element: µ dU µ /ds + Γαβ U α U β = f µ , where the extra force, f µ , is given by   1 λF2 f = (Lm − p) ∇ν R + ∇ν p hµν . ρ + p 1 + λ f2 µ (78) Note that an intriguing feature is that a specific choice for the matter Lagrangian yields different dynamics. For instance, in [180], it was argued that a “natural choice” for the matter Lagrangian density for perfect fluids is Lm = p, based on Refs. [181,182], where p is the pressure. This choice has a particularly interesting application in the analysis of the R−matter coupling for perfect fluids, which implies in the vanishing of the extra force [169]. However, it is important to point out that despite the fact that Lm = p does indeed reproduce the perfect fluid equation of state, it is not unique [183]. Other choices include, for instance, Lm = −ρ [182,184], where ρ is the energy density, or Lm = −na, were n is the particle number density, and a is the physical free energy defined as a = ρ/n − Ts, with T being the temperature and s the entropy per particle (see Ref. [182,183] for details). Hence, it is clear that no immediate conclusion may be extracted regarding the additional force imposed by the non-minimum coupling of curvature to matter, given the different available choices for the Lagrangian density. One may conjecture that there is a deeper principle or symmetry that provides a unique Lagrangian density for a perfect fluid [183]. This has not been given due attention in the literature, as arbitrary gravitational field equations depending on the matter Lagrangian have not always been the object of close analysis. See Ref. [183] for more details. 4.3. The Palatini Approach The discussion outlined in the beginning of this section hints at reconsidering the foundational aspects of gravity as a geometric phenomenon. The key point of such an approach is to make one step back and reconsider the role that basic geometric objects play in gravitational physics. The first such step is to recognize the different role that metric and affine connection have as geometrical structures, and take them as independent objects (Palatini approach [27]), which implies the presence of non-metricity (this is actually the simplest extension towards more general geometric theories of gravitation), i.e., the failure of the independent connection to be metric, Qλµν = ∇λ gµν 6= 0. Note that whether the underlying structure of spacetime is Riemannian or not is a fundamental question of gravity as a geometric phenomenon which has received little attention so far, and its associated phenomenology has been scarcely studied. Indeed, a glance at the modern description of ordered structures with defects in condensed matter physics (like Bravais crystals) teaches us that the presence of point defects induces a continuous description in terms of a geometry with independent metric and affine structures [185–187], which could provide important tools for gravitational physics [188,189]. In the case of GR, one finds that the Palatini approach gives exactly the same result as the standard approach of fixing the Levi-Civita connection a priori and making variations with respect to the metric gµν (metric approach); this result being attached to the particular functional form of the Einstein-Hilbert Lagrangian of GR. Nonetheless, for modified gravity it turns out that it avoids many of the problems of the metric formulation, and generically yields a set of second-order field equations and absence of ghost-like instabilities. To fix ideas let us consider an f (R) extension of GR, defined by the standard action S= 1 2κ 2 Z d4 x p − g f (R) + Sm , (79) where κ 2 = 8πG, g is the determinant of the space-time metric gµν , which is independent of ρ ρ the connection Γλµν . The Ricci tensor, Rµν ≡ Rα µαν = −∂µ Γλλν + ∂λ Γλµν + Γλµν Γρλ − Γλνρ Γµλ , is an object entirely constructed out of the independent connection, while the matter sector, Symmetry 2016, 8, 70 27 of 79 R √ Sm = d4 x − gLm (ψm , gµν ) (where ψm denotes the matter fields) only couples to the metric gµν . We point out that the independence between metric and connection in the Palatini approach also allows for the presence of torsion (the antisymmetric part of the connection). However, for simplicity, in this section, we will assume that torsion vanishes, Γ[µν] = 0, and, in addition, that the Ricci tensor is purely symmetric, R[µν] = 0 (for a derivation of the field equations in full generality under the presence of torsion see [190]). Under these conditions, variation of the action (79) yields δS = 1 2κ 2 Z d4 x p −g h f R Rµν −  i 1 f gµν δgµν + f R gµν δRµν (Γ) + δSm , 2 (80) Since we are working on the Palatini approach, the terms on the variations with respect to δgµν and δRµν (Γ) must vanish independently. This brings two systems of equations of the form f f R Rµν − gµν 2 p ∇Γλ ( − g f R gµν ) = κ 2 Tµν , (81) = 0, (82) where f R ≡ d f /dR and the matter energy-momentum tensor Tµν is defined by Equation (64). Note that in the GR case, f (R) = R, Equation (82) simply expresses the standard metric-connection compatibility condition, while Equation (81) become the Einstein equations. In the f (R) case, tracing with gµν in Equation (81) yields R f R − 2 f = κ 2 T, which is an algebraic relation whose solution, R = R( T ), generalizes the GR relation R = −κ 2 T. The field Equation (81) can be solved by introducing a new rank-two tensor hµν = f R gµν in such a way that the connection Equation (82) √ is rewritten as ∇Γλ ( − hhµν ) = 0, i.e., the independent connection Γλµν becomes the Levi-Civita connection of hµν . This tensor is connected to the physical metric gµν via the conformal transformations given above. With this result, the remaining Equation (81) can be written as [191] κ2 Rµ ( h) = 2 fR ν   f (R) ν ν δµ + Tµ . 2κ 2 (83) where Rµ ν (h) is the Ricci tensor computed with the auxiliary metric hµν . This is a set of second-order Einstein-like equations with all the objects in the right-hand-side depending only on the matter sources (because f (R) ≡ f (R( T )) is a function of the trace of the energy-momentum tensor), which can thus be read as a modified matter source, namely, one can write Gµν = τµν , where τµν accounts for an effective energy-momentum tensor. This property, which results from the special way gravity and matter couple in the Palatini approach, is in sharp contrast with the generic fourth order equations and ghosts that one finds in the (more conventional) metric formulation of these theories (see e.g., [192,193] for a comparison between both approaches). Note that in the absence of matter, Tµ ν = 0, the field equations boil down to those of GR with a cosmological constant term and thus these theories contain no extra propagating degrees of freedom (in particular, the theory is free of ghosts). As a consequence, the dynamics of these theories can be only excited in regions of large curvatures/short distances or for strong fields. Solving the Equation (83) just requires to specify the matter source, Tµ ν , and the object f (R( T )), and to transform the solution hµν back to the physical metric gµν using the conformal transformations above. 4.3.1. Quadratic Cosmology from Gravity-matter Coupling This framework can be extended to include a gravity-matter coupling term g(R) Lm (ψm , gµν ) replacing the matter sector in the action (79), which yields a new set of equations [194] Rµ ν ( h) = κ2 Φ2   f (R) ν ν δ + g (R) T . µ µ 2κ 2 (84) Symmetry 2016, 8, 70 28 of 79 where now Φ ≡ f (R) + 2κ 2 gR Lm , hµν = Φgµν and hµν = Φ−1 gµν . Again, these equations are manifestly second-order due to the same considerations of the previous case. In the following, we consider a FRW cosmology, in the presence of a perfect fluid. Assuming an equation of state P = ωρ (with ω =constant), and inserting these inputs into the field Equation (84) one gets two sets of equations [194] H2 = f ( R) + κ 2 g( R)(1 + 3ω )ρ , 6Φ(1 − ∆)2 Ḣ (1 − ∆) + (85) H2 1 [1 − ∆2 + 6( ρ + p ) ∆ ρ ] = [ f ( R) − κ 2 g( R)(5ρ + 3p)] , 2 12Φ (86) Φ where ∆ ≡ 23 (1 + ω )ρ Φρ and a subindex means a derivative. Both the f (R) and the matter-gravity coupling frameworks can successfully reproduce the background evolution of loop quantum cosmology, whose chief accomplishment is the resolution of the Big Bang singularity applying non-perturbative canonical quantization techniques [195,196]. The singularity is replaced by a non-singular cosmic bounce at an energy of Planck’s density order ρc ≡ c5 /(h̄G2 ) via a quadratic cosmology model of the form   ρ κ2 2 , (87) H = ρ 1+ǫ 6 ρc with ǫ ± 1, where ǫ = +1 corresponds to the solution arising in braneworld scenarios [197] and ǫ = −1 in loop quantum cosmology [198]. In f (R) gravity, quadratic cosmology can be numerically fit with a specific model [199] (which is consistent with the fact that space-time singularities can also be avoided in spherically symmetric charged black hole scenarios of Palatini f (R) gravity [200,201]), while in f ( R) theories with gravity-matter coupling different choices of Φ ≡ Φ(ρ) yields different models able to analytically reproduce this cosmic evolution [194]. This result can further be extended to isotropic and anisotropic cosmologies in Palatini gravity including quadratic Rµν Rµν corrections [202]. 4.3.2. Bouncing Cosmologies in Born-Infeld-like Gravity Models A somewhat different proposal for modified gravity models was introduced by Deser and Gibbons, coined Born-Infeld gravity [203,204], where the action is given by SBI = λ 4 Z 4 d x q gµν + λ −2 R µν − p −g  , (88) where λ is a small parameter, which must be suitably tuned to recover GR at low energies/curvatures. This means that, from an experimental point of view, the action (88) represents a viable theory of gravity, while encompassing new effects at high energies. Born-Infeld gravity has raised great interest due to its many astrophysical and cosmological applications, including the existence of singularity-free cosmologies. To investigate the latter point, let us rewrite, by convenience, the first term in Born-Infeld gravity p √ as − g det M̂, where a hat denotes a matrix and M̂ = ĝ−1 q̂, with qµν ≡ gµν + λ−2 Rµν . Therefore we can rewrite (88) as Z p
SBI = λ4 d4 x − g det M̂ − 1 . (89) Writing in this way, it is clear that Born-Infeld gravity (89) contains two of the invariant polynomials of the matrix M̂ and, therefore, a natural extension of this model would be [205] S = λ4 Z d4 x 4 p −g ∑ β n (en ( M̂)) , n =0 (90) Symmetry 2016, 8, 70 29 of 79 where en are the five invariant polynomials of the matrix M̂ (see [205] for details) and β n is a set of free dimensionless parameters. In this formulation, the original Born-Infeld gravity comes from the contributions e0 ( M̂) = 1 and e4 ( M̂) = det M̂. When testing the non-singular character of solutions in the Palatini approach it is important to establish its robustness, namely, that the singularity resolution is not attached to a particular choice of the functions defining the model, or to some ad hoc procedure to obtain a desired result, like a fine tuning of the parameters. This program has been very successful when applied to spherically symmetric black holes replacing the point-like singularity by a wormhole structure giving continuity to all null and time-like geodesics in a number of gravitational and matter scenarios [206–208]. Similarly, to study bouncing solutions, let us focus on the simple extension of Born-Infeld gravity given by 1 S= 2 κ ǫ Z d4 x p h p i − g f ((det Ω̂)1/2 ) − λ̃ − g + Sm (ψm , gµν ) , (91) where we have redefined the constants to match with GR results at low energies, ǫ ≪ 1, introduced a new matrix Ω̂ ≡ M̂1/2 , for notational simplicity, and added the matter contribution. As usual, in the Palatini formalism, we treat metric and connection as independent entities and obtain the field equations by variation with respect to both of them. After a bit of (lengthy) algebra one gets the field equations for the auxiliary metric tµν as [209] Rµ ν (t) = κ2 ( L G δµ ν + T µ ν ) , 2φ2 (det Ω̂)3/2 (92) with the following definitions: φ ≡ d f /dA and V (φ) = A(φ) f A − f ( A) are functions of the auxiliary scalar (Brans-Dicke-like) field A, the auxiliary metric tµν (associated with a Levi-Civita √ connection, ∇λ [ −ttµν ] = 0) is related to the physical metric gµν as tµν = φ(det Ω̂)1/2 gµα Ω̂αν , and LG ≡ (φ(det Ω̂)1/2 − V (φ) − λ̃)/(ǫκ 2 ) is the Lagrangian density. Note that the equation resulting from variation of Equation (91), namely, det Ω̂ = dV/dφ, establishes an algebraic relation φ = φ(det Ω̂) for every choice of f ( A). Moreover, due to the fact that the equations for the metric gµν can be written as
φ(det Ω̂) + V (φ) + λ̃ µ κ2 ǫ µ −1 µ φ(det Ω̂)[Ω̂ ] ν = δ ν− T ν, (93) 2 2 this implies that Ω̂ = Ω̂( Tµ ν ) and, consequently, φ = φ( Tµ ν ). Therefore, all terms appearing on the right-hand-side of the field Equation (92) are just functions of the matter, like in the case of f (R) theories. Thus, finding a solution of (92) for tµν immediately yields a solution for gµν . Let us assume a perfect fluid energy-momentum tensor with N species of non-interacting fluids, namely, ρ = ∑iN=1 ρi and p = ∑iN=1 pi , with equations of state pi = ωi (ρi )ρi , in a flat FRW spacetime. The Hubble function then reads [209] " # 2 ( ρ + 3p ) + 2( φ (det Ω̂ ) − V − λ̃ ) ǫκ 1 , (94) ǫH 2 =  2 ˙ φ(det Ω̂) + V + λ̃ + ǫκ 2 ρ 6 1+ ∆ 2∆H where ∆ = 2φ2 (det Ω̂)3/2 φ(det Ω̂)+V +λ̃−ǫκ 2 P and ∆˙ 3 =− 2∆H 2∆ N ∑ i =1 " # ∂(det Ω̂) ∂∆ ∂∆ (1 + ωi ) ρ i . + ∂ρi ∂ρi ∂(det Ω̂) (95) To obtain explicit solutions one needs to further specify the family of Born-Infeld gravity extensions. Let us consider the power-law family f (det Ω̂) = (det Ω̂)n where the case n = 1/2 Symmetry 2016, 8, 70 30 of 79 corresponds to the original Born-Infeld gravity. Expanding in series of ǫ ≪ 1, the corresponding Lagrangian yields the result lim S = ǫ ≪1 Z 4 d x p  (1 − λ̃) n −g + 2 R + O ( ǫ ) + Sm , ǫκ 2 κ  (96) which, via a redefinition of constants, is just GR plus a cosmological constant term. Therefore, this kind of theories is in agreement with observational data in the low energy/curvature regime, but provides new physics regarding the early universe. In this sense, at high energies two kinds of solutions are found (see [209] for details): (i) those with ǫρ > 0 and ω > 0 have a vanishing Hubble parameter for some finite energy density ρ B characterized by H 2 ∝ (ρ B − ρ) ∝ ( a − a B )2 , and represent unstable solutions with a minimum value at past infinity, where a(−∞) = a B ; (ii) those with ǫρ < 0 represent genuine non-singular bouncing solutions (as long as ω > −1) with H 2 ∝ (ρ B − ρ) ∝ ( a − a B ) ∝ (t − t B )2 . This is in agreement with the fact that the avoidance of spacetime singularities in black holes within Born-Infeld gravity is achieved when ǫ < 0. As the short discussion above shows, high-energy modifications of GR within the Palatini approach ( f (R), quadratic, Born-Infeld-type, and other generalizations) are able to avoid the big bang singularity via a cosmic bounce followed by a period of de Sitter expansion in a generic and robust manner. At the same time, Palatini theories generically produce an effective cosmological constant at low energies (see Equation (96)). The specific sign and value of this constant depends both on the gravitational action and on the specific choices of parameters. This suggests to consider the combination of these two properties to achieve a model in which a sufficient small cosmological constant compatible with cosmological observations arises at late times as a result of high-energy or low-energy modifications of the theory. This could provide a unified and consistent scenario accommodating the two known de Sitter phases of the cosmic expansion without invoking dark sources of energy. On the other hand, as first suggested by Hawking [210], the existence of electrically charged microscopic black holes originated by large density fluctuations during the early universe (primordial black holes) could have a substantial impact in accounting for the matter density of the hypothesized dark matter particles. The prediction of some Palatini models of the existence of singularity-free, horizonless massive objects below a certain charge/mass scale [211], could offer new avenues for the consideration of black hole remnants [212] as a reasonable source of dark matter. Obtaining precise estimates on the abundance of such objects, as resulting from their production in the early universe, could allow to make detailed cosmological predictions according to this viewpoint. 4.4. Hybrid Metric-Palatini Gravity The hybrid metric-Palatini gravitational theory is a novel approach to modified theories of gravity that consists of adding to the Einstein-Hilbert Lagrangian an f (R) term constructed a la Palatini [29]. Using the respective dynamically equivalent scalar-tensor representation, it was shown that the theory can pass the Solar System observational constraints even if the scalar field is very light. This implies the existence of a long-range scalar field, which is able to modify the cosmological and galactic dynamics, but leaves the Solar System unaffected. The absence of instabilities in perturbations was also verified and the theory provides explicit models which are consistent with local tests and lead to the late-time cosmic acceleration. 4.4.1. General Formalism The hybrid metric-Palatini gravity theory proposed in [29], is given by the following action S= Z d4 x p h i − g R + f (R) + 2κ 2 Lm . (97) Symmetry 2016, 8, 70 31 of 79 can be recast into a scalar-tensor representation [29,213,214] given by the action S= 1 2κ 2 Z d4 x p   3 − g (1 + φ ) R + ∇ µ φ ∇ µ φ − V ( φ ) + Sm , 2φ (98) where Sm is the matter action, κ 2 = 8πG/c3 , and V (φ) is the scalar field potential. Note that the gravitational theory given by Equation (98) is similar to a Brans-Dicke scalar-tensor action with parameter w = −3/2, but differs in the coupling to curvature. The variation of this action with respect to the metric tensor provides the field equations (1 + φ) Gµν = κ 2 Tµν + ∇µ ∇ν φ − ∇α ∇α φgµν − 3 3 1 ∇µ φ∇ν φ + ∇ φ∇λ φgµν − Vgµν , 2φ 4φ λ 2 (99) where Tµν is the matter energy-momentum tensor. The scalar field φ is governed by the second-order evolution equation   1 dV φκ 2 1 µ ∇µ φ∇ φ + φ 2V − (1 + φ) T, = − ∇α ∇ φ + 2φ 3 dφ 3 α (100) which is an effective Klein-Gordon equation [29,213]. 4.4.2. Unifying Local Tests and the Late-Time Cosmic Acceleration In the weak field limit and far from the sources, the scalar field behaves as φ(r ) ≈ φ0 + (2Gφ0 M/3r )e−m ϕ r ; the effective mass is defined as m2ϕ ≡ (2V − Vφ − φ(1 + φ)Vφφ )/3 φ=φ , where 0 φ0 is the amplitude of the background value. The metric perturbations yield (2) h00 (r ) 2Geff M V0 r2 = + , r 1 + φ0 6 (2) hij (r ) =  2γGeff M V0 r2 − r 1 + φ0 6  δij , (101) where the effective Newton constant Geff and the post-Newtonian parameter γ are defined as Geff ≡  G  1 − (φ0 /3) e−m ϕ r , 1 + φ0 γ≡ 1 + (φ0 /3) e−m ϕ r . 1 − (φ0 /3) e−m ϕ r (102) As it is clear from the above expressions, the coupling of the scalar field to the local system depends on φ0 . If φ0 ≪ 1, then Geff ≈ G and γ ≈ 1 regardless of the value of m2ϕ . This contrasts with the result obtained in the metric version of f ( R) theories, and, as long as φ0 is sufficiently small, allows to pass the Solar System tests, even if the scalar field is very light. In the modified cosmological dynamics, consider the spatially flat Friedman-Robertson-Walker (FRW) metric ds2 = −dt2 + a2 (t)dx2 , where a(t) is the scale factor. Thus, the modified Friedmann equations take the form 2 = 2 Ḣ = 3H    1 φ̇ V 2 κ ρ + − 3φ̇ H + , 1+φ 2 4φ   3 φ̇2 1 − φ̈ −κ 2 (ρ + P) + H φ̇ + 1+φ 2 φ (103) (104) respectively. The scalar field Equation (100) becomes φ̈ + 3H φ̇ − φ̇2 φ φκ 2 + [2V − (1 + φ)Vφ ] = − (ρ − 3P). 2φ 3 3 (105) Symmetry 2016, 8, 70 32 of 79 As a first approach, consider a model that arises by demanding that matter and curvature satisfy the same relation as in GR. Taking V (φ) = V0 + V1 φ2 , (106) the trace equation automatically implies R = −κ 2 T + 2V0 [29,213]. As T → 0 with the cosmic expansion, this model naturally evolves into a de Sitter phase, which requires V0 ∼ Λ for consistency with observations. If V1 is positive, the de Sitter regime represents the minimum of the potential. The effective mass for local experiments, m2ϕ = 2(V0 − 2V1 φ)/3, is then positive and small as long as φ < V0 /V1 . For sufficiently large V1 one can make the field amplitude small enough to be in agreement with Solar System tests. It is interesting that the exact de Sitter solution is compatible with dynamics of the scalar field in this model. Relative to the galactic dynamics, a generalized virial theorem, in the hybrid metric-Palatini gravity, was extensively analyzed [215]. More specifically, taking into account the relativistic collisionless Boltzmann equation, it was shown that the supplementary geometric terms in the gravitational field equations provide an effective contribution to the gravitational potential energy. The total virial mass is proportional to the effective mass associated with the new terms generated by the effective scalar field, and the baryonic mass. This shows that the geometric origin in the generalized virial theorem may account for the well-known virial mass discrepancy in clusters of galaxies. In addition to this, astrophysical applications of the model where explored, and it was shown that the model predicts that the effective mass associated to the scalar field, and its effects, extend beyond the virial radius of the clusters of galaxies. In the context of the galaxy cluster velocity dispersion profiles predicted by the hybrid metric-Palatini model, the generalized virial theorem can be an efficient tool in observationally testing the viability of this class of generalized gravity models. Thus, hybrid metric-Palatini gravity provides an effective alternative to the dark matter paradigm of present day cosmology and astrophysics. In a monistic view of physics, one would expect Nature to make somehow a choice between the two distinct possibilities offered by metric and Palatini formalisms. We have shown, however, that a theory consistent with observations and combining elements of these two standards is possible. Hence gravity admits a diffuse formulation where mixed features of both formalisms allow one to successfully address large classes of phenomena. 4.4.3. Future Outlook: More gEneral Hybrid Metric-Palatini Theories The “hybrid” theory space is a priori large. In addition to the metric and its Levi-Civita connection, one has also an additional independent connection as a building block to construct curvature invariants from [216]. Thus one can consider various new terms such as Rµν Rµν , Rµν Rµν , Rµναβ Rµναβ , Rµναβ Rµναβ , RR , etc . (107) Though an exhaustive analysis of such hybrid theories has not been performed, there is some evidence that the so called hybrid class of theories we are exclusively focusing our attention upon here is a unique class of viable higher order hybrid gravity theories. In the more restricted framework of purely metric theories, it is well known that the f (R) class of theories is exceptional by avoiding the otherwise generic Ostrogradski instabilities by allowing for a separation of the additional degrees of freedom in a harmless scalar degree of freedom: as we have already seen, such a separation is possible also for the hybrid theories. Furthermore, it turns out that this feature is a similar exception in the larger space of metric-affine theories, since a generic theory there is inhabited by ghosts, superluminalities or other unphysical degrees of freedom. 4.5. Slow-Roll Inflation in f ( R) Gravity Cosmological inflation has been a great success over the last decades as it solves in an elegant way the problems raised in the Big Bang model: the so-called horizon and flatness problems, Symmetry 2016, 8, 70 33 of 79 basically initial conditions problems. Inflation consists of initial super accelerated phase, during which the above-mentioned problems are solved by requiring a sufficient duration of such a phase (see Ref. [217–219] and references therein). After the end of inflation, the universe reheats and the original hot state of the Big Bang model is recovered. The original idea was proposed by Alan Guth in 1981, which was later improved by Andrei Linde in 1982. Over the years inflation has become a predictable theory, since it provides the seeds, coming from quantum fluctuations, that form the basis of the formation of large scale structures as well as the anisotropies observed in the CMB. Here, we review the essentials of the inflationary paradigm and we focus on the inflationary models within the context of extended theories of gravity, particularly the so-called f ( R) gravity. Furthermore, the usual observables predicted by inflation are obtained and compared with the latest data released by the Plank mission [220]. 4.5.1. Inflationary Parameters Slow-roll inflation has been widely studied in the literature, particularly within the framework of scalar fields, since by the appropriate scalar potential, inflation can be easily reproduced and its exit achieved at the corresponding number of e-folds (see Ref. [221] and references therein). This is probably the simplest model for slow-roll inflation, whose general Lagrangian is given by a single minimally coupled scalar field, Sφ = Z d4 x p   1 − g − ∂µ φ∂µ φ − V (φ) , 2 (108) Then, by assuming a flat FLRW metric, the FLRW equations are obtained 1 3 2 H = φ̇2 + V (φ) , 2 κ2 −  1  2 1 3H + 2 Ḣ = φ̇2 − V (φ) . 2 κ2 (109) Whereas the scalar field equation is given by φ̈ + 3H φ̇ + ∂V (φ) = 0. ∂φ (110) Inflation occurs basically when the scalar field behaves approximately as an effective cosmological constant, and the friction term in (110) dominates H φ̇ ≫ φ̈ while V ≫ φ̇2 . Inflation ends when the kinetic term of the scalar field increases and dominates over the potential, oscillating about the minimum of the potential and reheating the universe. A useful way of describing slow-roll inflation is to define the so-called slow-roll parameters: 1 ǫ= 2 2κ  V ′ (φ) V (φ) 2 , η= 1 V ′′ (φ) , κ 2 V (φ) λ2 = 1 V ′ (φ)V ′′′ (φ) , κ 4 (V (φ))2 (111) where the primes denote derivatives with respect to the scalar field φ. The magnitude of the quantities in (111) have to be small enough during inflation in order to attain the sufficient number of e-folds, or in other words ǫ ≪ 1 and η < 1. When inflation is ending, ǫ ≥ 1 is required. Thus, for the above model the spectral index ns of the curvature perturbations, the tensor-to-scalar ratio r of the density perturbations, and the running of the spectral index αs can be written in terms of the slow-roll parameters as follows ns − 1 = −6ǫ + 2η , r = 16ǫ , αs = dns ∼ 16ǫη − 24ǫ2 − 2ξ 2 . d ln k (112) The values of the above parameters can be compared then with the values provided by the latest Planck results (see [220]). Note that the aim here is to provide the above expressions (112) in terms of a particular underlying model, capable of reproducing slow-roll inflation, in order to compare with Symmetry 2016, 8, 70 34 of 79 the observational data. For the simple model (108), the slow-roll parameters are already expressed in terms of the scalar potential in (111), while the reconstruction of the scalar field Lagrangian (108) given a particular inflationary model H ( N ), is provided by the FLRW equations (109), which yield [222] 2H ′ (φ) , κ 2 H (φ) i 1 h 2 ′ + H ( φ ) H ( φ ) , 3 H ( φ )) ( κ2 ω (φ) = − V (φ) = (113) where ω (φ) is the kinetic term and we have redefined the scalar field as φ = N, where N is the number of e-folds N = log( a/ai ). Hence, by assuming a particular Hubble parameter, the corresponding Lagrangian is easily obtained. Note also that the slow-roll parameters (111) can be expressed in terms of the Hubble parameter by using (113), leading to [223–225] ǫ = η = ξ2 = " 2 #2   −2  ′ H′ (N) H′ (N) H ′′ ( N ) H (N) H(N) − ′ 3+ 6 , + + 4H ( N ) H(N) H(N) H(N) H(N)   −1 "  2 1 H ′′ ( N ) 1 H ′ ( N ) H′ (N) H′ (N) − +3 + 9 3+ 2 H(N) H(N) H(N) 2 H(N) #  ′′ 2 1 H (N) H ′′ ( N ) H ′′′ ( N ) − , + 3 + 2 H′ (N) H′ (N) H′ (N)  ′   −2 " 2 #  1 H (N) H ′′ ( N ) H′ (N) H′ (N) H ( N ) H ′′′ ( N ) + + 6 3+ 3 4 H(N) H(N) H(N) H(N) H ′ ( N )2 (114) (115) 3 H ( N ) H ′′ ( N ) H ′′′ ( N ) H ( N ) ( H ′′ ( N )) H′ (N) H ′′ ( N ) −2 + (116) +4 3 H(N) H(N) ( H ′ ( N )) ( H ′ ( N ))4 #  ′′ 2 2 H ′′′ ( N ) H (N) H ′′ ( N ) H ( N ) ( H ′′ ( N )) H ( N ) H ′′′′ ( N ) , − + 15 ′ −3 + +5 ′ H (N) H′ (N) H (N) ( H ′ ( N ))3 ( H ′ ( N ))2 +9 respectively. By assuming a particular Hubble parameter, the slow-roll conditions are straightforward calculated, while the corresponding Lagrangian for the scalar field is also obtained through the expressions (113). In the next section, we review previous works where instead of considering a scalar field, inflation is driven by modifications of the Hilbert-Einstein action. 4.5.2. Inflation in f ( R) Gravity Over the last years, f ( R) gravity has been widely analyzed, since it can reproduce late-time acceleration as well as inflation with no need of an extra field and/or cosmological constant (for a review see [9]). Particularly, the reconstruction of the corresponding f ( R) action that leads to the appropriate cosmological evolution has drawn much attention (see Ref. [226,227]). Here we review some aspects of f ( R) gravities and specifically its roll within the inflationary paradigm. Let us start with the general action for f ( R) gravity, given in the action (62). The field equations are obtained by varying the above action with respect to the metric tensor, yielding Equation (63). Then, by assuming a spatially flat FLRW metric, the modified FLRW equations are obtained for f ( R) gravity:     f ( R) + 3 H 2 + Ḣ f ′ ( R) − 18 4H 2 Ḣ + H Ḧ f ′′ ( R), 2   ...  f ( R)  − Ḣ + 3H 2 f ′ ( R) + 6 8H 2 Ḣ + 4 Ḣ 2 + 6H Ḧ + H f ′′ ( R) 2
2 +36 4H Ḣ + Ḧ f ′′′ ( R) . 0 = − 0 = (117) (118) Symmetry 2016, 8, 70 35 of 79 As shown in Ref. [228,229], the corresponding slow-roll parameters (117) can be expressed in terms of the f ( R) action as well. Nevertheless, here we are interested in reviewing slow-roll inflation in f ( R) gravity by using the scalar-tensor equivalence of f ( R) gravity and expressing the action (62) in the Einstein frame. Firstly, note that the action (62) in the Jordan frame, can be rewritten by action (67), so that by applying the conformal transformation g̃µν = Ω2 gµν where Ω2 = φ , (119) the action (67) is transformed as follows S̃ = Z 4 d x p   R̃ 1 µ − g̃ − ∂µ φ̃∂ φ̃ − Ṽ (φ̃) , 2 2κ 2 where φ=e κ q 2 φ̃ 3 , Ṽ = e −2κ q 2κ 2 (120) 2 φ̃ 3 V. (121) Hence, the Lagrangian is now written in the form of (108), such that slow-roll inflation occurs when the conditions described above are satisfied. In addition, the slow-roll parameters (111) can be expressed in terms of the original f ( R) action by using the correspondence expressions (66) and (121), leading to 1 (−2 f + f R R)2 , 3 f 4 ( f − R f R )2 R ǫ = η = − (122) 2 f R2 − 4 f RR + R f f R f RR . 3 f RR f R2 ( f − R f R ) Then, given a particular form of the f ( R) action, the corresponding slow-roll parameters can be calculated in terms of the Hubble parameter, which is obtained eventually by solving the FLRW equations. In addition, the corresponding spectral index and scalar-to-tensor ratio are obtained and for a particular number of e-folds, (most of inflationary models require Ne ∼ 55 − 65), both magnitudes can be compared with the constraints from the Planck results [220]. Let us now consider the so-called Starobinsky inflation [230], which has been a very successful inflationary model, taking into account the observational constraints of the Planck satellite. The model has also been analysed from the quantum point of view, where it was shown to provide a correct behaviour and a successful inflationary period [231]. The model corresponds to a particular form of the f ( R) function consisting of the Einstein–Hilbert term and a quadratic correction, given by f ( R) = R + R2 . 6m2 (123) The parameter m2 in the Starobinsky model is a free parameter which is constrained by the observations, and it is in fact related to the mass of the new, dynamical scalar degree of freedom present in the model, usually denoted as the “scalaron”. Although the inflationary dynamics of the model can be analysed either in the Jordan or in the Einstein frame, here we will focus our discussion on the Einstein-frame action, applying the procedure described earlier. Notice that according to the previous discussion, the requirement of the positivity of the second derivative of f , f RR > 0, translates into the requirement that the squared mass m2 remains positive, i.e., m2 > 0. Symmetry 2016, 8, 70 36 of 79 In order to re-express the model in the Einstein frame, we note that for the Starobinsky model the Einstein-frame scalar is given by φ̃ = r   3 −1 φ κ log 1 + . 2 3m2 (124) Given the last relation, and using Equation (121), it is straightforward to derive the Einstein frame potential as Ṽ (φ̃) = 2 √ 3 2 − 2/3κ φ̃ . m 1 − e 4κ 2 (125) Therefore, the Einstein-frame action will be described by the action (120), with the potential defined through Equation (125). During inflation, the background spacetime is described by a FRW spacetime, which after the conformal transformation (119) is expressed as follows ds2 = −dt2 + ã(t)2 dx2 , (126) where ã(t) denotes the scale factor in the Einstein frame. Similarly to Equation (109), the Friedmann equation in this spacetime is described by H̃ 2 ≡  2 ã˙ ã ! φ̃˙ 2 + Ṽ (φ̃) , 2 κ2 = 3 (127) while the equation of motion for the scalar φ takes the form φ̃¨ + 3 H̃ φ̃˙ + Ṽ ′ (φ̃) = 0, (128) where Ṽ ′ (φ̃) ≡ dṼ (φ̃)/dφ̃. The first and second term on the right hand side of the Friedmann equation (127) describe the kinetic and potential energy of the inflaton. Recall that slow–roll inflation occurs in the regime where κφ ≫ 1, i.e., as long as the field φ̃ rolls down the flat part of the potential (125) sufficiently slowly, and the kinetic energy of the inflaton is much smaller than its potential energy, φ̃˙ 2 ≪ Ṽ (φ̃). Then the Friedmann and the scalar field equation, Equations (127) 2 and (128), respectively, become approximately H̃ 2 ≃ κ6 Ṽ (φ̃), and 3 H̃ φ̃˙ ≃ −Ṽ ′ (φ̃). The Hubble function H̃ is approximately constant during inflation, and its time variation is described by the slow–roll parameter ǫ defined as ǫ≡− H̃˙ . H̃ 2 (129) In the slow-roll approximation it is required that ǫ is small, ǫ ≪ 1. We can approximate it as ǫ η ≃ ≃ 1 2κ 2  Ṽ ′ (φ̃) Ṽ (φ̃) 2 4 ≃ 3  q 2 κ φ̃ 3 q 2 κ φ̃ e −1  −2 4 −2 ≃ e 3 2−e 3 4 − 4 1 Ṽ ′′ (φ̃) ≃− e ≃  q  2 2 3 3 κ Ṽ (φ̃) 2 κ φ̃ −1 + e 3 q q 2 κ φ̃ 3 2 κ φ̃ 3 , (130) . The last step in the above relation assumed that κ φ̃ ≫ 1. Inflation ends as long as ǫ ≃ 1, which corresponds to the inflaton’s kinetic energy becoming important. Symmetry 2016, 8, 70 37 of 79 The amount of inflation is measured by the number of e-foldings N which is defined as N≡ Z t end tstart H̃dt , (131) and for agreement with observations one usually requires that N ≃ 55 − 65. Now, under the slow-roll approximation, the above relation takes the following approximate form N ≃ −κ 2 Z φ̃ end Ṽ ( φ̃ ) φ̃start Ṽ ′ (φ̃) dφ ≃ 3 √2/3κ φ̃start . e 4 (132) In the last step, we used the fact that in the slow-roll approximation, κ φ̃ ≫ 1, and we performed the integration with respect to φ̃. Notice that the number of e-foldings is related to the slow-roll parameters as 1 3 1 (133) , η≃− . ǫ≃ 4 N2 N The scalar and tensor fluctuations produced during inflation have a characteristic amplitude and scale dependence, which in the Einstein frame, under the slow-roll approximation, are given by the standard relations for GR minimally coupled to a scalar field [217–219], PR ≃ 1 H̃ 4 N 2 m2 , ≃ 2 3π m2p 4π φ̃˙ 2 (134) 16 H̃ 2 4 m2 ≃ , π m2p π m2p (135) PGW ≃ while the tensor to scalar ratio follows from the above relations as r ≡ PGW /PR ≃ (48/π ) N −2 ≃ 16ǫ. Notice that in the final approximate equalities of the above relations we restricted ourselves to the Starobinsky potential for κ φ̃ ≫ 1, under the slow-roll approximation. The spectral indices corresponding to the scalar and tensor fluctuations, respectively, under the slow-roll regime are given by the relations (112). Note that Planck provides a value for the spectral index ns = 0.968 ± 0.006 and an upper bound for the scalar-to-tensor ratio r < 0.07 coming from a joint analysis of the BICEP2/Keck Array and Planck data [220]. In order to illustrate the power of Starobinsky inflation, let us assume a number of e-folds Ne = 65, which leads to the following values of the inflationary observables: ns = 0.968 , r = 0.00284 . (136) Hence, Starobinsky inflation satisfies quite well the constraints coming from Planck data. Let us consider now a slight modification of Starobinsky inflation to show the uniqueness of the above model, f ( R) = R +  R 6m  2+ λ , (137) where λ is a positive parameter that accounts for the deviation from Starobinsky inflation. In this case, and following the same steps as above, the scalar potential in the Einstein frame yields,  Ṽ (φ̃) = A Be κ q 2 φ̃ 3 −1  12++λλ e −2κ q 2 φ̃ 3 , (138) Symmetry 2016, 8, 70 Here A = λ 2κ 2 (6m)2+λ 38 of 79 and B = (6m)2+λ 2+ λ . This actually does not lead to an effective cosmological constant when κ φ̃ ≫ 1 since at that limit, the scalar potential becomes Ṽ (φ̃) ≃ AB 2+ λ 1+ λ e −κ q 2 3 ( 1+λ λ )φ̃ . (139) Since λ is positive, the potential is an exponential potential at large values of the scalar field which only mimics an effective cosmological constant in the case of λ = 0 when Starobinsky inflation is recovered. Otherwise, for κ φ̃ ≫ 1, the scalar potential is exponentially suppressed, and inflation does not occur. In the case of λ < 0, the scalar field will roll fast, invalidating the slow-roll conditions. Hence, deviations from λ = 0 will not reproduce slow-roll inflation, what may lead to other kinds of inflationary models as chaotic inflation which present deviations from the values of the spectral index and the scalar-to-tensor ratio as pointed by Planck Collaboration [220]). Other extensions of f ( R) gravity can also generate slow-roll inflation. In this sense, the inflationary paradigm has been recently explored within the so-called f ( R, Lφ ) gravities, where Lφ = − 21 ∂µ φ∂µ φ − V (φ) (see [232]). Such kind of theories extends f ( R) gravities by including a nonminimally coupling scalar field. As shown in Ref. [232], this model is equivalent to multifield inflationary scenarios, where isocurvature modes may become important, while adiabatic perturbations can be analytically calculated (see [233]). However, following [232] the adiabatic perturbations can be obtained and the isocurvature modes may be ignored under certain conditions on the action. 5. Topological Defects The Standard Model of Particle Physics describes three of the fundamental interactions between particles — the electromagnetic, weak and strong forces — in a unified quantum field theory. According to this model, the universe underwent, in its early stages, a series of symmetry breaking phase transitions, some of which may have led to the formation of topological defects. 5.1. Spontaneous Symmetry Breaking and Topological Defects To illustrate the process of topological defect formation as a result of spontaneous symmetry breaking, let us consider the simplest model that admits topological defect solutions: the Goldstone model with a single scalar field, φ. In this model, the Lagrangian density is given by L= 2 λ 2 1 φ,µ φ,µ − φ − η2 , 2 4 (140) where λ is a coupling constant. Although this Lagrangian is invariant under transformations of the form φ → −φ — thus having a Z2 symmetry — the Goldstone potential has two degenerate minima at φ = ±η. Hence, once φ acquires a vacuum expectation value (VEV), φ± = h0| φ |0i = ±η , (141) the Z2 symmetry is spontaneously broken, leading to the formation of a kink separating regions where the scalar field has different VEVs (here, |0i is the ground state of the model). In a 1 + 1 dimensional Minkowski spacetime, the equation of motion of the scalar field, ∂2 φ ∂2 φ dV , − 2 =− 2 dφ ∂t ∂x (142) admits the following solution h i φ( x ) = η tanh η (λ/2)1/2 γ( x − x0 − vt) , (143) Symmetry 2016, 8, 70 39 of 79 with γ = (1 − v2 )−1/2 , that corresponds to a kink moving with velocity v in the positive x-direction. Since φ interpolates between φ = −η at x = −∞ and φ = η at x = ∞, there is a region (the domain wall) centred at x = x0 + vt which has non-vanishing energy. This configuration R has a conserved topological charge, Q = dxj0 = φ(+∞) − φ(−∞) = 2η (where jµ = ǫµν φ,ν is the topological current), and, for this reason, it is stable. The constant vacuum state has Q = 0, and, as a result, this configuration is unable to relax into the vacuum while conserving its topological charge. In N + 1 dimensions, the solution in Equation (143) remains valid and represents a ( N − 1)-dimensional planar surface. These surfaces also separate two domains with different VEVs and, for this reason, they are commonly denominated Domain Walls. Domain walls arise due to the breaking of a discrete symmetry, whenever the vacuum manifold is disconnected. However, different types of defects may be formed depending on the (non-trivial) topology of the vacuum manifold (or the type of symmetry being broken). Line-like defects, or Cosmic strings, are formed if the vacuum is not simply connected or, equivalently, if it contains unshrinkable loops. This type of vacuum manifold results, in general, from the breaking of an axial symmetry. Moreover, if the vacuum manifold contains unshrinkable surfaces, the field might develop non-trivial configurations corresponding to point-like defects, known as Monopoles. The spontaneous symmetry breaking of more complex symmetry groups may lead to the formation of textures, delocalized topological defects which are unstable to collapse. 5.2. A VOS Model for Topological Defect Networks of Arbitrary Dimensionality The production of topological defects is predicted in several grand unified scenarios. In most instances, these topological defects may survive throughout cosmological history until late cosmological times, and leave a variety of imprints on different observational probes. A unified phenomenological framework for the description of the cosmological evolution of topological defects of arbitrary dimensionality was developed in [234–236]. This model—which is a generalization of the Velocity-dependent One-Scale (VOS) model for cosmic strings [237]—provides a statistical description of the large-scale dynamics of these networks. This generalized VOS model characterizes quantitatively the macroscopic dynamics of the p-dimensional topological defect network in ( N + 1)-dimensional Friedmann-Robertson-Walker backgrounds by describing the evolution of the characteristic length of the network, L, and of its root-mean-squared velocity (RMS), v̄. The characteristic length of a statistically homogeneous defect network may be defined as ρ̄ = σp N L −p , (144) where ρ̄ is the average defect energy density, and σp is the defect mass per unit p-dimensional area. For infinitely thin and featureless topological defects, the following equations describing the evolution of the characteristic length and velocity of the network (L and v̄, respectively) may be obtained directly from the generalized Nambu-Goto action [236] dv̄ dt dL dt = =  v̄ 1 − v̄ − , L ℓd L 2 c̃ HL + v̄ + v̄ , D ℓd D  2 k (145) (146) −1 1 where D = N − p, and ℓ− d = ( p + 1) H + ℓ f is a damping length scale that includes, not only the effects of cosmological expansion, but also the effect of the frictional forces caused by the scattering of particles off the defects (which is encoded in the friction length scale, ℓ f ). Here, k and c̃ are phenomenological parameters that account, respectively, for the effect of defect curvature and the energy loss associated to defect reconnection on the macroscopic dynamics of the network. The main advantage of this phenomenological model lies precisely on these two parameters. They may either be calibrated using numerical simulations—as was done for cosmic strings [238] and domain Symmetry 2016, 8, 70 40 of 79 walls [239]—to accurately describe the large scale dynamics of standard networks, or be used in a phenomenological description of a large variety of alternative defect scenarios. This simplest version of the VOS model has been extended in order to allow for the description of non-standard defect scenarios, with additional energy-momentum decay channels (e.g., direct energy loss through gravitational radiation [238]) or additional degrees of freedom (e.g., the existence of conserved currents [240]). Other extensions allow for the description of more complex defect scenarios—such as multi-tension cosmic (super-)string networks with junctions [241], semi-local cosmic string networks [242], or hybrid defect networks [243]—by including terms to account for the energy-momentum transfer between different network components. 5.3. Observational Signatures of Topological Defects The observational consequences of topological defect networks can be very diverse, depending both on the type of defects formed and on the evolution of the universe after they are generated. For instance, the average domain wall density grows more rapidly than the background density in the radiation and matter eras and, consequently, domain walls must necessarily be very light in order not to completely dominate the energy density of the universe at late times. Although the possibility of a significant contribution to the dark energy budget has been ruled out both dynamically [36] and observationally [37], it is possible for light domain walls to leave behind interesting astrophysical and cosmological signatures. Other topological defects, such as cosmic strings in a linear scaling regime, have an average density that is roughly proportional to the background density. Hence, standard cosmic string networks do not tend to dominate the energy density of the universe and they naturally generate a roughly scale-invariant spectrum of density perturbations on cosmological scales. Although recent cosmological data, in particular the cosmic microwave background anisotropies, rules out cosmic strings as the main source of large-scale cosmological perturbations (see, however, [244]), it remains plausible that they may play an important role on small scales. In fact, it has been shown that the power spectrum of a string-induced density perturbation component is expected to dominate over a primordial inflationary component on small enough scales [245]. In this context, cosmic strings have been suggested to have significant impact on the formation of ultracompact minihalos [39], globular clusters [40], super-massive black holes [41] and to provide a significant contribution to the reionization history of the Universe [42]. The dynamics of cosmic defects, in particular of domain walls, has also been associated to possible variations of the fundamental couplings of nature (see, e.g., [38]). In this section, we shall focus on the contribution of both cosmic strings and domain walls to two of the most significant observational probes: the temperature and polarization anisotropies of the cosmic microwave background and the stochastic gravitational wave background. 5.3.1. CMB Anisotropies The presence of a topological defect network at late cosmological times would necessarily leave imprints on the Cosmic Microwave Background (CMB). Although current CMB observations seem to be consistent with the inflationary paradigm (in which the fluctuations are seeded in the very early universe) they also allow for a subdominant topological defect contribution [246]. Topological defects source metric perturbations actively throughout the cosmological history. For this reason, their CMB signatures are expected to be fundamentally different from those due to primordial fluctuations. In particular, cosmic defects are expected to generate a significant vector component that would not be present in inflation-seeded scenarios, because vector modes decay rapidly in the absence of a source [247]. The B-mode polarization signal originated by topological defects has then contributions from both tensor and vector modes, and may, for this reason, produce an observationally relevant signal in this channel, despite providing only subdominant contributions to the temperature and E-mode power spectra. The B-mode polarization channel thus offers a relevant observational Symmetry 2016, 8, 70 41 of 79 window for the detection of topological defects. With several CMB experiments planned for the near future—e.g., LITEBIRD, PRISM and CoRE+ (see [248] for a recent overview)—there is the prospect either for the detection of topological defects or for the tightening of current observational constraints in this scenario. Cosmic string and domain wall networks produce temperature anisotropies characterized by significantly different angular power spectra. Domain walls only become cosmologically relevant at late times, when the characteristic length scale of the network is large. For this reason, these networks contribute mostly to the TT angular power spectrum at large scales (or small multipole modes ℓ). Since low-ℓ observational data is highly affected by cosmic variance and has thus large associated uncertainties, current data allows for a significant contribution of standard domain wall networks to the temperature power spectrum. In Ref. [249], the authors have derived a conservative constraint on the domain wall mass per unit area of σ < 3.52 × 10−5 kg m−2 allowed by current observational data (which corresponds to a constraint on the wall-forming symmetry breaking scale of η < 0.92 MeV). Cosmic strings, on the other hand, generate a TT angular power spectrum that has a broad peak at ℓ ∼ 300, and therefore contribute dominantly to the angular power spectrum at higher multipoles (which have small observational uncertainties). For this reason, observational data allow for a significantly larger contribution of domain walls to the temperature power spectrum on large angular scales than it does of cosmic strings. Current CMB observational data constrain the cosmic string tension to be Gµ < 1.3 (2.4) × 10−7 for Nambu-Goto (Abelian-Higgs) strings [37] (which is significantly smaller than the corresponding constraint on domain wall mass per unit area GσL0 = 5.6 × 10−6 ). Theses differences in the shape and amplitude of the temperature power spectra generated by domain wall and cosmic string networks are clearly illustrated in the left panel of Figure 5, where the TT power spectra generated by cosmic string and domain wall networks with the maximum tensions allowed by current observational data are plotted. We chose GσL0 = 5.6 × 10−6 for the domain wall networks, and the weakest constraint on cosmic string tension obtained using Planck data [37], Gµ = 2.4 × 10−7 (which corresponds to the constraint on Abelian-Higgs strings). Figure 5. The total temperature (left panel) and B-mode (right panel) angular power spectra generated by cosmic string and domain networks with the maximum tension allowed by current observational data: GσL0 = 5.6 × 10−6 and Gµ = 2.4 × 10−7 for domain walls and cosmic string, respectively. The spectra were obtained using, respectively, the CMBACT code [250,251], and the extension introduced in [249]. These plots have been adapted from [249]. Although cosmic defect contributions to the temperature power spectrum are necessarily subdominant when compared to the primordial inflationary contribution, they may produce a significant B-mode polarization signal. As for primary temperature anisotropies, cosmic strings and domain walls contribute to the BB angular power spectrum at different scales: domain walls contribute dominantly at low-ℓ, while strings have a significant contribution at small scales. Moreover, the B-mode signal generated by domain wall networks has, in general, smaller amplitude than that of cosmic strings. In the right panel of Figure 5, we plot the total BB power spectrum generated by domain wall networks and cosmic string networks that have the maximum fractional contribution to the TT power spectrum allowed by current observational data. This figure shows that, since the observational constraints on Gµ are more stringent than those on GσL0 , the domain Symmetry 2016, 8, 70 42 of 79 wall contribution may dominate over that of cosmic strings for low multipole modes (particularly if their energy density is close to the maximum allowed by current data). Moreover, the fact that domain walls and cosmic strings have spectra with different shapes and that these defects contribute mostly at different scales should, in principle, allow one to distinguish between their contributions if a signal is detected. 5.3.2. The Stochastic Gravitational Wave Background Cosmic string interactions often result in the formation of loops that detach from the long string network (for instance, in self-intersections or when kinky strings collide). After formation, these loops oscillate relativistically under the effect of their tension, and lose their energy in the form of gravitational waves. These loops, thus, have a finite life span and, as a consequence, the cosmic string network loses energy through this channel. The dynamical effects caused by this energy loss are taken into account in the VOS model through the inclusion of the last (phenomenological) term in Equation (146). The VOS model may be used to estimate the rate of loop production (per comoving volume) as a result of cosmic string interactions dnc c̃ v̄ , (147) = dt α L4 where α is a parameter that characterizes the size of loops, lc = αL(tc ), at the moment of creation, tc . This rate of loop production is strongly dependent on the large-scale properties of the network (through its dependence on v̄ and L) and, therefore, its accurate characterisation requires an accurate description of cosmic string network dynamics. The production of cosmic string loops is expected to occur copiously throughout the evolution of cosmic string networks. At any given time, there are several loops emitting gravitational waves in different directions. The superimposition of these emissions generates a stochastic gravitational wave background (SGWB) that spans a wide range of frequencies [252–255]. The SGWB is often quantified using the spectral density of gravitational waves per logarithmic frequency ( f ) interval, in units of critical density, 1 dρGW , (148) ΩGW = ρc d log f where ρGW is the gravitational radiation energy density. This quantity is strongly dependent on the rate of loop production as a function of time, and, consequently, on the large-scale cosmic string parameters. This spectrum has a characteristic shape: it has a flat portion in the high frequency range, generated by the emissions of loops created during the radiation era, and a pronounced peak at low frequencies that is produced during the matter era. Note however that the amplitude of the spectrum and the characteristics of the peak (its slope, broadness, location and height) vary greatly: they are not only dependent on the cosmic string tension and the large scale properties of the network, but also on the small-scale properties of the loops that generate it. For instance, the size of loops and their emission spectrum dramatically affects the overall shape of the spectrum. The dependence of the SGWB spectrum on these parameters is illustrated in Figure 6, where SGWB spectra generated by cosmic string networks are plotted, for a wide range of values of these macroscopic and microscopic parameters. There is the prospect of probing different parts of the SGWB spectra using current and upcoming gravitational wave detectors in the near future—most notably (Advanced) LIGO [256], (Advanced) VIRGO [257], eLISA [258] and pulsar timing arrays [259,260]. Understanding the shape of the spectrum generated by cosmic strings and its dependence on macroscopic and microscopic parameters is, therefore, crucial. Symmetry 2016, 8, 70 43 of 79 Figure 6. Clockwise, from the top-left Panel: The stochastic gravitational wave background generated by cosmic string networks: 1) For different loop sizes, in the large loop regime (where α is larger than the gravitational back-reaction scale, ΓGµ). Here, cosmic string tension is fixed to Gµ = 10−7 , and we have only considered one harmonic mode of emission (ns ). 2) For different values of α in the small-loop regime (α < ΓGµ). Here, Gµ = 10−7 and ns = 1. 3) For different value of cosmic string tension, Gµ. Here, the loop size parameter is set to α = 10−7 , and ns = 1. 4) For an increasing number of harmonic modes of emission, ns . For larger values of ns , the shape of the spectrum saturates and is identical to that of the spectrum with ns = 104 . Here, α = 10−8 , Gµ = 10−7 , and we have only considered loops with cusps. These plots were adapted from [255]. As in the case of cosmic strings, as domain wall networks evolve and interact, a fraction of their energy is released in the form of gravitational waves. In fact, in [261,262], it has been shown using field theory simulations that domain wall networks produce a spectrum of gravitational waves that peaks at the scale corresponding to the size of the Hubble radius at the time of their decay. Studies of the gravitational wave spectrum generated by domain wall networks are, however, preliminary and limited to the case of metastable networks. 6. Cosmological Tests with Galaxy Clusters at CMB Frequencies The Cosmic Microwave Background (CMB) has revolutionized the way we perceive the Universe. The information encoded in the temperature and polarization of the CMB provides one of the strongest evidences in favour of the hot Big-Bang theory and has enabled ways to constrain cosmological models with unprecedented accuracy [44]. The CMB also provides a unique way to study the formation and evolution of cosmic structure. Physical effects acting on the CMB after decoupling originate secondary CMB anisotropies that encode additional information about the growth of cosmological structure and the dynamics of the Universe [45]. The Sunyaev-Zeldovich (SZ) effect [46,47] is the scatter of CMB photons by plasma in galaxy clusters and filamentary structures. It is among the strongest and most studied sources of secondary anisotropies. Whereas the primordial CMB signal has been widely used to constrain cosmological parameters, the SZ effect is only now giving the first solid steps into adulthood as a probe of Cosmology and structure formation. Largely contributing to this state of affairs are ground-based and satellite observations such as those by SPT [263] and Planck [264] and proposals for future experiments such as in [265]. In this section we review some of the key features of the SZ effect and its application to test cosmological models and the physics of galaxy clusters, in light of Planck cluster observations. We Symmetry 2016, 8, 70 44 of 79 end with an example of the application of cluster data to the study alternative cosmological scenarios involving primordial non-Gaussianities. 6.1. Galaxy Clusters Galaxy clusters are the largest gravitationally bound objects in the Universe. They are formed by gravitational accretion of infalling material at the intersection of filamentary structures. Their typical masses range from 1014 − 1015 h−1 M⊙ and contain hundreds to several thousand galaxies within radii of the order of 1 − 1.5 h−1 Mpc. The intra-cluster medium (ICM) is filled with hot ionized gas, typically at temperatures 1–15 keV, which produces strong X-ray emission and microwave spectral distortions via the SZ effect, see Figure 7. Figure 7. Images of a Galaxy cluster in a hydrodynamic N-Body simulation. The main panel shows a column density map, in logarithmic scale, of a region centered in the cluster. The smaller maps inside the zoomed box show the intensity and polarization of the CMB radiation that emerges from the cluster. The mapped quantities are the cluster thermal SZ effect (top left); kinetic SZ effect (top right); CMB quadrupole induced polarization (bottom left); and tangential velocity polarization (bottom right), [266], due to the interaction of the CMB with the hot gas in the cluster. Because clusters are among the latest bound structures forming in the Universe their number density is highly sensitive to the growth rate and power spectrum of density perturbations as well as to cosmological background parameters like the Universe’s matter density and the Hubble parameter. The central quantity that goes into the computation of the expected abundance of clusters for a given cosmological scenario is the mass function, n( M, z). It gives the comoving number density of virialized halos forming at the redshift z with mass M. As introduced by [267], this can be expressed as: M dn( M, z) f (σ, z) ≡ , (149) ρ̄ d ln σ−1 where f (σ, z) is the so called multiplicity function, ρ̄ is the mean background density and σ( M, z) is the variance of the linear density field at redshift z, smoothed on a scale R, enclosing the mass M, σ2 ( M, z) = D2 (z) 2π 2 Z ∞ 0 k2 P(k )W 2 (k, M )dk. (150) Symmetry 2016, 8, 70 45 of 79 Here, D (z) is the growth factor of linear perturbations, P(k) is the linear density power spectrum extrapolated to redshift zero, and W (k, M) is the Fourier transform of the smoothing kernel in real space. The multiplicity function can be computed analytically or by assessing the number of halos forming in numerical N-body simulations. The following formulae give the multiplicity function for the popular Press–Schecter, f PS , [268] and Jenkins, f J , [267] mass functions: f PS = r   δc2 2 δc exp − 2 , π σ 2σ f J = 0.315 exp(−| ln σ−1 + 0.61|3.8 ) . (151) Here δc is the extrapolated linear overdensity of a spherical perturbation at the time of collapse. The Press–Schecter mass function was the first being derived using analytical arguments, but it overpredicts/underpredicts the abundance of low/high mass halos in simulations. The Jenkins mass function was established using simulations. Several other simulation based mass functions have been proposed to improve the functional form of the multiplicity function, see e.g., [269] for a summary with other multiplicity function expressions. Most of these studies have been focusing on improving systematic uncertainties related to the halo mass definition and quantifying deviations from the universality of f for different cosmologies. This includes extensive mass function studies using primordial non-Gaussian N-body simulations, see e.g., [270], and f ( R) modified gravity simulations, see e.g., [271]. The redshift distribution of observed galaxy cluster abundances can be compared against model predictions using the expression: d2 V d2 N = dz dΩ dzdΩ Z ∞ dn( M, z) 0 dM f survey ( M, z) dM . (152) Here d2 N/dzdΩ is the number of clusters per unit redshift per unit solid angle, d2 V/dzdΩ is the volume element of the underlying model, and f survey ( M, z) is the cluster survey selection function. Setting f survey = 1 and assuming that its possible to observe all cluster masses, Equation (152) gives the full redshift distribution of halos in the Universe, see e.g., [272] for an application to dark energy cluster counts. The f survey function depends on survey parameters such as sky coverage and noise distribution of the observed quantity Sν , which is usually a cluster flux or band luminosity. In a best case scenario, where the sky noise is uniform and the survey is complete above a given flux limit, Slim , R Equation (152) can be computed as an integral over the observed flux, (dn/dM) (dM/dSν ) dSν , with Sν ≥ Sν,lim . Here dM/dSν is the derivative of a scaling relation between M and Sν . Galaxy cluster number counts are therefore a powerful probe of cosmology provided the relation between total mass and the observed survey flux are known. Poor knowledge about the later may imply large uncertainties in the determination of cosmological parameters, as well as leaving a number of unanswered questions regarding the underlying cluster formation mechanism. If gravity is the dominating force driving the formation and evolution of collapsed structures, virialized clusters should have their thermodynamic properties related to the total mass (through the viral theorem). The simplest model describing cluster scaling relations was presented by [273]. This assumes, among other things, that gravitational collapse is scale-free, i.e., clusters of different masses are self-similar replicas of each-other. If non-gravitational effects, such as radiative gas physics and energy feedback, influence the thermodynamic state of clusters one should expect modified scalings. These are usually parameterized as power laws of the form [269], Y = A E(z) βss ( X/X0 )α , (153) where Y, X are cluster properties, A is a normalization parameter giving the amplitude of Y at X = X0 , E(z) = H (z)/H0 , α is the power law index of the (independent) property X, and β ss is a parameter giving the self-similar evolution of the normalization. In general, the amplitude parameter A may itself be a function of redshift. In this case, A evolves in a non self-similar way, usually Symmetry 2016, 8, 70 46 of 79 parameterized by a power-law of redshift, A(z) = A0 (1 + z) β . The parameters α and β ss can be computed analytically within an extended self-similar modelling that also assumes that the baryon gas component in clusters is in hydrostatic equilibrium [269]. For example, under these conditions, the scaling between the SZ integrated flux and mass, Y − M, is given by Equation (153) with the indexes α = 5/3, and β ss = 2/3. Hydrodynamic N-body simulations are the appropriate method to study galaxy cluster scaling relations as they allow us to model the baryon physics throughout the full non-linear evolution of the ICM gas and determine the parameters α, β ss , β and A0 . We refer the reader to references [274–277] for examples of simulations studies showing the impact of various gas physical models, dark energy and primordial non-Gaussianities on several galaxy cluster scalings. Numerical N-body simulations have also been intensively used to study the internal structure of galaxy clusters. These studies demonstrate that galaxy clusters profiles can be used as a sensitive probe of cosmology, cluster physics and structure formation mechanisms. One of the most well known results from simulations, in the context of hierarchical CDM models, is that the radial distribution of dark (collissionless) matter in clusters follows a simple expression known as the Navarro-Frenk-White (NFW) density profile [278] ρNFW (r ) = 4ρs , x (1 + x )2 (154) where x ≡ r/rs , rs is a scale radius, and ρs is the characteristic density at r = rs . Subsequent numerical simulations allowed one to improve this result and to study the impact on cluster profiles of alternative cosmological models such as modified gravity, see [269] for a review. Gas profiles have also been extensively investigated with hydrodynamic N-body simulations. For example, [279] found that the radial pressure profile, P(r ), in SZ observations of resolved clusters follows a universal profile, P0 P (r ) , (155) = γ P500 x (1 + x α )( β−γ)/α where x = r/rs , rs is the scale radius, and P0 , α, β, γ are fitting parameters to derive from observations. P500 is the integrated pressure at r500 (the radius at which the cluster density exceeds 500 times the universe mean background density). Equation (155) is known as the generalized NFW (GNFW) profile of the SZ signal in clusters. 6.2. The SZ Effect: Temperature and Polarization Signal The Sunyaev-Zel’dovich effect is the scattering of CMB photons by electrons in hot reservoirs of ionized gas in the Universe. A particularly effective scattering medium is the ICM where the gas is at temperatures of the order of 108 Kelvin. When the much colder and isotropic CMB radiation propagates through these plasma regions, a fraction of the original CMB photons scatter off the moving electrons and suffer frequency Doppler shifts towards higher energies (inverse Compton scattering) causing a spectral distortion in the CMB spectrum. In the process the total number of photons is conserved. In fact, the SZ effect can be regarded as a composition of two effects: the thermal and kinetic SZ effects. The thermal effect is due to the thermal (random) motion of the electrons in the gas cloud. If the gas has also a bulk peculiar velocity in the CMB radiation frame, there is a kinematic Doppler effect in addition to the thermal effect. The spectral distortions caused by the thermal, ∆Ith , and kinetic, ∆Ik , SZ effects are given by simple formulae [47]: ∆Ith = I0 g( x ) y , ∆Ik = − I0 h( x ) vr τ, c y= kB σT me c2 τ= Z Z Te ne dl , (156) σT ne (l ) dl , (157) Symmetry 2016, 8, 70 47 of 79 where I0 = 2k3B T 3 /h2 c2 , T is the CMB temperature, y is the line-of-sight Compton SZ parameter, τ is the optical depth, vr is the line-of-sight velocity of the center of mass of the gas cloud (positive/negative when the cloud is receding/approaching), Te is the electron density, ne is the electron number density, and g( x ), h( x ) are functions giving the frequency dependence of the effects. Figure 8 shows maps of the thermal and kinetic SZ distortions imprinted on the CMB by cosmic structures in hydrodynamic N-body simulations. The mapped quantities are the y parameter and the temperature fluctuations of the kinetic SZ effect, (∆T/T )k = ∆Ik /( I0 h( x )) = −vr /c τ, see e.g., [280]. Figure 8. Images of the thermal y-comptonization parameter (left panel) and kinetic temperature fluctuations (∆T/T )k = −vr /c τ (right panel) of the SZ effect from hydrodynamic N-body simulations of large scale structure. Maps show the same sky patch, with an angular size of 5 degrees2 . The bright objects in the y map are galaxy clusters. For typical cluster velocities and optical depths, the thermal SZ effect is the dominant effect. Therefore the total SZ signal integrated over the sky angular size of a cluster is well approximate by Y≡ Z 2 y dΩ = d− A Z y dA = k B σT −2 d mc2 A Z V Te ne dV , (158) where d A is the angular diameter distance from the observer to the cluster. This is known as the volume integrated SZ Y-flux or Y-luminosity and is a measure of the total thermal energy density of the cluster. The scattering of CMB photons by hot gas in clusters also generates new (secondary) polarization signatures in the CMB that carry information about the intra-cluster medium. The dominant CMB secondary polarization effect in clusters is due to the existence of an intrinsic quadrupole component from the CMB temperature anisotropy around the scattering media. This effect is known as the CMB quadrupole induced polarization. In addition, the bulk motion of the gas cloud may also induce CMB polarization proportional to the square of the sky transverse velocity of the cloud. This effect is often referred to as induced kinetic polarization. Finally, double scattering events inside the cluster gas cloud also originate polarization, known as double scattering induced polarization. The bottom panels in Figure 7 are detailed simulations of the linear polarization degree and polarization angle of the CMB quadrupole and kinetic induced polarizations in a cluster. These simulations are an essential tool to assess the degree of contamination of the primordial CMB signal by secondary polarization effects and allow to prepare observational strategies to use CMB cluster polarization as an additional Symmetry 2016, 8, 70 48 of 79 probe of cosmology and structure formation, see [266] for a detailed description of cluster polarization simulations and their use as a new observing tool. 6.3. Cosmological Tests Using SZ Cluster Surveys SZ galaxy cluster surveys are recognized as having a great potential to test a wide range of cosmological and astrophysical scenarios. SZ galaxy cluster counts, profiles, scaling relations and SZ angular power spectra are most promising probes to confront model predictions with observations. In this section we will briefly review the current state of cosmological constraints from observations with the Planck satellite. 6.3.1. SZ Cluster Counts The Planck SZ cluster survey is the largest and deepest all-sky cluster survey yet produced at CMB frequencies. Its latest galaxy cluster catalogue (PSZ2, [281]) contains 1653 detections, of which 1203 are confirmed clusters with identified counterparts in other observing bands (X-rays, optical, near infra-red and other SZ experiments). More than a thousand of these clusters possess redshift information. From this catalogue, the authors in [264] selected a sample of 439 clusters with high signal to noise (q = S/N > 6) to infer constraints on cosmological parameters using SZ galaxy cluster number counts. For a given cosmological model, the redshift distribution of clusters abundances detected by Planck above a given noise level, q > qcat can be computed from Equation (152) as: dN (q > qcat ) = dz Z d2 V dΩ dzdΩ Z dn( M, z) χ̂( M, z, l, b, qcat ) dM , dM (159) see [264]. Here M is the mass contained within a cluster radius, R500 , where the cluster density exceeds by 500 times the mean background density of the Universe, χ̂( M, z, l, b, qcat ) = R∞ qcat P [ q | q̄m ( M, z, l, b )] dq is the Planck selection function, (i.e., χ̂ = f survey in Equation (152)) and P[q|q̄m ( M, z, l, b)] is the distribution of q given the mean signal-to-noise value, q̄m ( M, z, l, b), predicted by the model for a cluster of mass M and redshift z located at Galactic coordinates (l, b). This expression reflects the fact that the survey noise level is not uniform on the sky and that clusters are detected above a given noise level, qcat . The distribution P[q|q̄m ] incorporates noise fluctuations and intrinsic scatter in the actual cluster SZ Y signal, (158), computed within R500 , around the mean value, Ȳ ( M, z), predicted from a Y − M scaling relation calibrated using Planck cluster data. This method was used in [264] to impose constraints on cosmological and cluster parameters, such as the variance of the density fluctuations smoothed on a scale of R = 8 h−1 Mpc, σ8 , see Equation (150), the matter density parameter, Ωm0 , the Hubble constant, H0 , and the slope of the Y − M scaling relation, α. Assuming a ΛCDM concordance (base) model their analysis show that Planck SZ cluster constraints on the σ8 − Ωm0 plane are in tension with those obtained from Planck primary CMB anisotropies, unless priors on the hydrostatic mass bias parameter, 1 − b, (which relates the observed baryon mass with the total mass of the cluster) and the effect of gas physics on the growth of structure are changed to values that are in conflict with constraints from other cosmological probes and the predictions from numerical simulations, [264,282]. A way to try to alleviate the tension is to combine cluster and primary CMB data to explore constraints on extensions to the base flat CDM model. This was carried out in [264] by allowing an extra varying parameter (to the base model), namely, the neutrino mass scale, the Thomson optical depth to reionization, the dark energy equation of state, and the Universe curvature parameter. The analysis seems to favor non-minimal neutrino masses (see Figures 11 and 14 in [264]), but this “solution” also implies a decrease on the Hubble constant to values that further deviate from most constraints by astrophysical determinations of H0 . The present Planck analysis also puts in evidence the importance of improving the precision of the cluster mass bias parameter, which in fact may depend on cluster mass and redshift. In the future, large lensing surveys, such as Euclid, WFIRST and LSST, are expected to determine mass biases with Symmetry 2016, 8, 70 49 of 79 one percent level accuracy. With this level of precision on 1 − b, the present cluster combined analysis of Planck data could by itself provide a stringent test of the ΛCDM concordance model, see Figure 13 in [264]. 6.3.2. SZ Power Spectra An important statistic to derive from sky maps of the SZ effect is the angular power spectrum of the thermal and kinetic SZ temperature fluctuations. The theoretical computation of these spectra is straightforward from simulations, such as those in Figure 8, but the computational resources involved to run simulations limit their use if one is required to compute the SZ power spectra for many model realizations. Analytical methods are much faster and allow to explore many model realizations, but they are much less capable of modelling the cluster gas physics and the contribution from the diffuse gas in filamentary structures. Analytically, the angular power spectrum of the thermal SZ signal can be computed as a sum clustering of two terms Cℓtotal = CℓPoisson + Cℓ , where CℓPoisson is the power spectrum resulting from a clustering Poisson distribution of objects (clusters) and Cℓ is a “2-halo” or clustering term giving the power arising from object (cluster two-point) correlations. These terms are given by (see [45] and references therein): CℓPoisson = clustering Cℓ Z zmax = 0 dV dzdΩ Z zmax 0 Z Mmax dn( M, z) dM Mmin dV P(k) × dzdΩ Z Mmax Mmin |yℓ ( M, z)|2 dM dz, dn( M, z) b( M, z)yℓ ( M, z) dM dM (160) 2 dz , (161) where yℓ is the 2D Fourier transform of the projected 3D Compton y-parameter radial profile of individual cluster, P(k, z) is the 3D matter power spectrum at redshift z, and b( M, z) is the linear bias factor that relates the matter power spectrum, P(k, z), to the cluster correlation power spectrum Pcluster (k, M1 , M2 , z) = b( M1 , z)b( M2 , z) D2 (z) P(k, z = 0). The mass and redshift range of the integrations are chosen so as to cover all relevant SZ source contributions up to high enough redshift. Observationally, the SZ all sky anisotropies can be separated from the overall CMB signal by applying component separation techniques that make use of the spectral signature and scale of the SZ effect. The first all-sky maps of the thermal SZ effect were produced from Planck observations [283,284]. These studies allowed to investigate the distribution of the y signal and characterize its higher order statistics, including the power spectrum and bi-spectrum. The observed thermal SZ power spectrum confirmed to be very sensitive to cosmological parameters, especially σ8 and Ωm0 , but these are 0.01 strongly degenerate with the mass bias parameter b (they found that σ8 (Ωm0 /0.28)3/8 = 0.80+ −0.03 if +0.01 3/8 b = 0.2 and σ8 (Ωm0 /0.28) = 0.90−0.03 if b = 0.4). These results are consistent with the constraints on σ8 (Ωm0 /0.28)3/8 from cluster number counts but are again in “soft tension” with the constraints obtained from primary CMB data assuming a ΛCDM concordance base model, see Figure 17 in [284]. 6.3.3. Probing New Physics with the SZ Clusters The present best CMB spectral observations by COBE/FIRAS are unable to determine if the CMB spectrum preserved its planckian shape at all redshifts and if the CMB temperature-redshift scaling is the one predicted by the standard Big-Bang model, T = T0 (1 + z), where T0 is the CMB temperature observed at z = 0. However it is possible to consider alternative models where the CMB has a planckian spectrum at z = 0, but a different temperature-redshift scaling. Observational tests of this scaling require measurements of the CMB temperature at high z, i.e., throughout the Universe’s history. The SZ effect has long been recognized as a most promising method to measure the CMB temperature at the redshift and location of clusters, allowing to test the CMB temperature-redshift relation at high redshifts, see e.g., [285] and references therein. Recent observations by Planck allowed to impose constraints on a simple deviation of the standard scaling, T = T0 (1 + z)1− β , using CMB Symmetry 2016, 8, 70 50 of 79 primary data [44] and SZ effect cluster data, see e.g., [286]. The emerging view is that Planck data impose tight constraints on the T (z) that are consistent with the standard CMB temperature scaling predicted by the standard Big-Bang model (i.e., β = 0). 6.4. Probing Primordial Non-Gaussianities with Galaxy Clusters The formation of the structures observed at large scales in the Universes is one of the most challenging and complex problems still open in modern cosmology. The widespread understanding is that these structures are the result of tiny primordial density fluctuations that might have existed at early times and that eventually grew under gravitational instability collapse. The origin of these tiny primordial density density perturbations is a cause of great debate, yet the inflationary paradigm remains the consensus framework for their generation. The statistical properties of the primordial density perturbations predicted by most standard single-field-slow roll inflationary models are nearly Gaussian in nature, and thus the level of non-Gaussianity is not strong enough to be observationally detected [55–57]. Although recent CMB anisotropies and LSS data seem to support such prediction, more complex models may lead to the generation of larger a non-Gaussianity amplitude, which substantially increases the chances of detection with future high precision observational experiments. For example, under certain conditions, some multi-field inflation models where slow-roll holds, can lead to larger departure from Gaussianity. However, non-Gaussianities may arise in inflationary models in which any of the conditions leading to the slow-roll dynamics fail [58]. The spectrum of models capable of generating larger non-Gaussianities amplitudes include the curvaton scenario [59–61], the ekpyrotic inflacionary scenario [62,63], vector field populated inflation [64–66] and multi-field inflation [67–70]. The current efforts that are being done to detect and constrain primordial non-Gaussianities (hereafter PNG) using a variety of cosmological observables is of the tremendous importance. A positive detection of non-Gaussianity signal would allow to asses the different proposed mechanisms for the generation of primordial density perturbations and in particular, it would make possible to discriminate between different inflation models and rule out some of them. Moreover, our understanding of the key physical processes at the early Universe would benefit greatly from such outcome. There is a wide range of observables to which we can resort to, in order to search for primordial non-Gaussianities. The traditional way to proceed is to measure three-point statistics of the CMB temperature anisotropies maps. Using this observable, there have been in the past reports of a positive detection of non-Gaussianity at a significant level [287]. However, although some contradicting analysis exist (see [288,289] and more recently [290,291]), the interest in testing deviations from Gaussian initial conditions has gained momentum ever since. Nevertheless, the possible existence of primordial non-Gaussianities could also be tested using other cosmological probes such as the characterization of the statistical properties of the structures at large scales. These probes include the the bispectrum and/trispectrum of distributions of galaxies and galaxy clusters (e.g., [292–294]), weak-lensing observations, (e.g., [295,296]) and CMB-LSS [297] as well as CMB-21cm line [298] cross-correlations. Another set of cosmological of probes very discussed in the literature to contain PNG, is the time evolution of the numerical abundances both of massive collapsed objects such as galaxy clusters, and its counterparts, large voids (see e.g [299,300] and [301–304]). 6.4.1. Parametrizing Primordial Non-Gaussianity The cosmological information contained in the power spectrum (or equivalently the two-point correlation function) of the primordial density perturbations is not enough if one wants to look for deviations from Gaussian initial conditions. Therefore we must resort to higher order correlation functions in order to gain further insight on nature of non-Gaussianities. Thus, the bispectrum, BΦ is the lowest order statistics sensible to primordial non-Gaussianities and it is defined as (3) hΦ (k1 ) Φ (k2 ) Φ (k2 )i = (2π )3 δD (k1 + k2 + k3 ) BΦ (k1 , k2 , k3 ) , (162) Symmetry 2016, 8, 70 51 of 79 (3) where Φ is the Bardeen’s potential and δD is the three dimensional Dirac delta. The bispectrum can be be written in the alternative form [305], BΦ (k1 , k2 , k3 ) = f NL (k1 k2 k3 )−2 A (k1 , k2 , k3 ) , (163) where the dimensionless parameter f NL accounts for the amplitude and the shape function A stores the information on the shape of the bispectrum. Due to the large number of inflationary models capable of generating primordial non-Gaussianities, it is common practice to group them according to their resulting bispectrum shape. There are generically four classes of bispectrum shapes frequently studied in the literature: Local, Equilateral, Folded and Orthogonal. However, from these shapes, the most studied in the literature is the Local one, which can be written in terms of Bardeen’s potential and performing a Taylor expansion around an auxiliary Gaussian random field Φ, (see eg. [57])   local Φ (x) = φ (x) + f NL φ2 (x) − hφ2 (x)i . (164) The corresponding bispectrum of the Local shapes is given by,   local Bφloc =2 f NL Pφ (k1 ) Pφ (k2 ) + Pφ (k1 ) Pφ (k3 ) + Pφ (k2 ) Pφ (k3 ) (165) To date, the strongest constraints on primordial non-Gaussianities of the Local type were obtained by the Planck collaboration [291] (see also [290] for previous constraint PNG by Planck local = 2.5 ± 5.7 from temperature alone and f local = 0.8 ± 5 combining collaboration), f NL NL temperature and polarization data. Note however, that although these constraints are consistent with Gaussian initial conditions, they were obtained assuming that the non-linear parameter f NL is scale-independent. Nevertheless, some inflationary models may lead to a scale dependent f NL (see [306] and references therein). 6.4.2. The Non-Gaussian Mass Function There is no unique prescription to incorporate the effects induces by primordial non-Gaussianities in the mass function (see [300,307–312]). The majority of the of the non-Gaussian mass functions are obtained resorting to some expansion of the probability density distribution (PDF) of the primordial density field and either using an extension of the Press-Schecter formalism [268] or the excursion set theory [313,314]. Generically, the modified mass function can be written in the form dn (z, M ) dnNG (z, M ) = G × RNG , dM dM (166) where RNG encodes the correction to the Gaussian mass function due to PNG. In (166), the Gaussian mass function can be simply the one of the PS formalism or any other motivated by N-body numerical simulation with gaussian initial conditions (references of fitting formulas), such as Equation (151). Two of the most discussed and used expression for the RNG correction in the literature are the MVJ [300] and the LV [306]. Both rely on the extension of the PS formalism to non-Gaussian initial conditions. The respective functional forms are, RNG RNG " #   S δ∗ 1 δc2 dS3,M 3,M 3 , + = exp δc 2 × 6 δ∗ d ln σM δc 6σM " 2 1 σM = 1+ S3,M 6 δec δ2 δc4 − 2 2c − 1 4 σM σM ! dS3,M + d ln σM δc2 −1 2 σM (167) !# (168) Symmetry 2016, 8, 70 52 of 79 where δc (z) is the collapse threshold at the redshift z, S3 is the reduced skewness and p δ∗ (z) = δc (z) 1 − S3 δc (z) /3. Figure 9 show the dependence of the non-Gaussian corrective term to the mass function, R NG , the mass and as a function of redshift, for the Equations (167) (blue line) and (168) (red dashed line). In both cases f NL = ±100 was used. It is clear that a higher positive/negative f NL value increases the number of predicted dark matter halos with higher mass. The effects are stronger at high redshift. The MJV prescription seems to predict a higher number of halos than LV for positive values of f NL . The opposite occurs when negative f NL values are considered. z=0 2.0 z=1 RNG 1.5 1.0 0.5 MVJ LV 0.0 1014 1015 M / (M⊙h −1 ) 1014 1015 −1 M / (M⊙h ) Figure 9. The non-Gaussian correction to the mass function due to primordial non-Gaussianities as a function of mass and for two refshifts: on the left panel z = 0 and on the right panel z = 1. Blue line corresponds to the prescription of MVJ (Equation (167)), while red dashed line corresponds to LV formula (Equation (168)). The level of non-Gaussianity of the Local type used was f NL ± 100. 6.4.3. Biased Cosmological Parameter Estimation with Cluster Counts The statistical properties of the primodial density perturbations has a great deal of influence on the redshift distribution of the number of galaxy clusters. Thus, these objects can be used as cosmological probes to search for deviations from Gaussian initial conditions. However, neglecting erroneously the existence of primordial non-Gaussianities, can lead to biases in the estimation of the cosmological parameters using clusters number counts. In [315], the authors investigated how primordial non-Gaussianities influence the estimations of the effective dark energy equation of state parameter, when the information contained in the abundance of galaxy clusters is used to explore different cosmological scenarios. In their work, they computed the effective dark energy equation of state per redshift bin, assuming Gaussian initial conditions, that would allow them to recover the galaxy clusters counts in different non-Gaussian models. From their findings, there resulted a Symmetry 2016, 8, 70 53 of 79 new diagnostic feature for the presence of PNG, in the form of an apparent evolution with time of the effective dark energy equation of state, characterized by the appearance of a discontinuity (more details see [315]). The same authors also estimated the magnitude of the biases that may arise in the determination of a larger set of cosmological parameters with galaxy clusters number counts, again ignoring that PNG may exist (see [316]). Their results indicate that, although the estimation of the present-day dark energy density and respective equation of state parameter do not seem to be sensible to non-Gaussian initial conditions in the primordial matter density field, the biases can be much larger in the case of the amplitude of the primordial density perturbations. Furthermore, the authors also argue that the conflicting constrains on the amplitude of the primordial perturbations obtained by the Planck collaboration using galaxy cluster number counts from the Planck Sunyaev-Zel’dovich Catalog and the primary Cosmic Microwave Background temperature anisotropies, could be alleviated, if a significant non-Gaussianity level at clusters scales exists (see [316] for more details). 6.4.4. The Impact of pRimordial Non-Gaussianities on Galaxy Clusters Scaling Relations In order to use galaxy clusters as cosmological probes, it is essential to have an accurate estimate of their mass. However, this task is not easy, since it is not possible to directly estimate the mass of a galaxy cluster from observational data, and thus one must rely on indirect methods, such as the scaling properties of clusters, to obtain a reliable estimate. The effect of PNG on galaxy clusters scaling relations was investigated for the first time in [277]. The authors resorted to a series of hydrodynamical N-body numerical simulations featuring adiabatic gas physics and allowing for different levels of non-Gaussian initial conditions within the fiducial ΛCDM model. In particular they studied the T − M, S − M, Y − M and YX − M scalings relating the total cluster mass with temperature, entropy and SZ cluster integrated pressure that reflect the thermodynamical state of the intra-cluster medium. Figure 10 shows as an example, the distribution of galaxy clusters for the Y − M, for three different levels of PNG, i.e., f NG ± 500, 0 obtained by [277], while Figure 11 shows the result of fitting Equation (153) for the Y − M scaling to the corresponding distribution of clusters, using the procedure suggested by [276,317]. The main findings of [277], indicate that PNG do have an impact on clusters scaling laws, in particular on the amplitude and redshift evolution of the scalings normalization. In the case of Y − M, the redshift evolution of the its normalization can change as much as 22% when f NL varies from −500 to 500. If f NL is a scale independent parameter and has a small value, as suggested by the current Planck constraints [290,291], then the impact of PNG on clusters scaling laws is negligible when compared with the effect of additional gas physics and other cosmological effects, such as dark energy. On the other hand, if f NL is instead a scale dependent parameter, as suggested by some inflationary models, then PNG could have larger positive/negative amplitude at clusters scales. In this case, the effects of PNG must be taken into account when galaxy cluster data are used to assess the constraining power of future cluster surveys or to estimate the cosmological parameters from observational cluster counts data. Symmetry 2016, 8, 70 54 of 79 Y / (h−2Mpc2) 10−4 10−5 10−6 fNL = −500 fNL = 0 fNL = 500 1014 1015 M / (M⊙h−1) Figure 10. The Y − M cluster scaling relation at redshift zero for values of f NL of −500 (blue squares), 0 (black circles) and 500 (red pentagons). Y−M 1.70 α 1.68 M ≥ 5x1013 Fit : M ≥ 5x1013 1.66 1.64 1.62 log(A) −5.788 −5.790 −5.792 β −0.04 −0.05 −0.06 −0.07 −500 −300 −100 0 100 fNL 300 500 Figure 11. The power-law index α of mass, normalization parameter log10 ( A) (note that here A is the A0 defined in Section 6.1), and power-law index of redshift β for the Y − M scaling, as a function of f NL , for a mass cut of 5 × 1013 M⊙ h−1 , and respective 1σ error bars. The blue solid line and shaded area are respectively, the linear fit to the data points and estimated 95% C.L. confidence interval. Symmetry 2016, 8, 70 55 of 79 7. Testing Gravity with Weak Lensing Gravitational lensing is the name commonly given to the phenomenon of deflection of light by gravitational fields. A lensing system consists of three parts: the source of light, the lens (a inhomogeneous distribution of mass, or a spacetime geometry) and the observer. Gravitational lensing of a point source by intervening matter produces image displacement, or even multiple images of the source. On extended sources, tidal gravitational fields lead to differential deflection, producing shape distortion and size magnification. Gravitational lensing on cosmological scales is a direct probe of Poisson-like constraint equations and probes structure formation independently of mass-to-light properties. It is thus a promising powerful way of testing gravity on cosmological scales and it is a core probe in the Euclid space mission (http://www.euclid-ec.org). 7.1. Gravitational Lensing Gravitational lensing [318–321] involve different types of astrophysical objects that spread in a broad range of mass, mass density and length scale. Various types of systems have been studied, both theoretically and observationally, with the aim of investigating the source, the lens, or the geometrical properties of the spacetime, from the observed images. Sources may range from stars, galaxies, to the last scattering surface, while lenses may be stars, galaxies, galaxy clusters or the large scale structure. The lensing properties are the same in all cases: the effect depends only on the projected two-dimensional mass density distribution of the lens and is independent of the luminosity or composition of the lens; gravitational lenses do not have a well-defined focal point; the effect is achromatic; there is no emission or absorption of photons; the surface brightness is conserved. The conservation of brightness implies that a change in size changes the flux, as if the source was closer to the observer, creating thus a zooming effect without loss of resolution, the so-called natural telescope. The lens equation is a mapping between source ( β) and image (θ ) planes, ~β = ~θ −~α. (169) The deflection field ~α(θ ) is a gradient of the gravitational potential and it is a central quantity of gravitational lensing. It is determined by the lensing potential, i.e., the gravitational potential integrated along the line-of-sight, and it depends thus on the geometry of the system (or the spacetime metric in cosmological-size systems). It is useful to Taylor-expand the lens equation, a procedure valid for small sources in angular size and small deflections (implying small angular size images). To linear order, the lens equation will thus be described by a Jacobian mapping, known as the amplification matrix: ! ∂ αi ∂β i . (170) = δij − Aij (θ ) = ∂ θj ∂ θj When applied to a light ray bundle emitted from an extended source, the symmetrical traceless part of the matrix produces a sheared image, the trace is responsible for an isotropically convergence (or divergence) and the antisymmetrical component introduces a rotation. This decomposition defines the optical components fields in the image plane: the scalar convergence field κ (~θ ) and the shear field γ(θ ). The convergence is given by the Laplacian of the lensing potential, κ= 1 (ψ, +ψ,22 ) , 2 11 (171) and it corresponds thus, through the Poisson equation, to the surface mass density in the effective two-dimensional plane where ψ is defined. The shear is a spin-2 quantity and is given by, γ1 = 1 (ψ, −ψ,22 ) , 2 11 γ2 = ψ,12 . (172) Symmetry 2016, 8, 70 56 of 79 The rotational field is not present, since the deflection field is a gradient and the amplification matrix is thus curl-free. Lensing effects are classified as strong or weak. Strong lensing effects occur at positions of the image plane where the convergence and shear fields take large values (> 1). In these cases multiple images may be formed, as well as strongly deformed multiple images that may merge forming giant arcs or even a full ring, the Einstein ring. The magnification of the images µ(θ ) is given by the inverse of the determinant of the amplification matrix. Images are magnified or demagnified depending on µ being greater or smaller than 1. The sign of the determinant defines the parity of the images, images in regions of det A < 0 are inverted with respect to the source. The points in the image plane where det A = 0 form closed curves named critical curves. The observed giant arcs form near critical curves. Strong lensing occurs at various scales: a quasar (source) - galaxy (lens) system produces multiple images of the quasar with a typical separation of arcsec, while a background star - foreground star system produces multiple images with a typical separation of miliarcsec, which is known as microlensing. If κ . 1 everywhere, multiple images do not form. An example in this intermediate regime are the arclets formed in the images of galaxy clusters, by lensing of background galaxies. The separation of the lensing regimes in an image depends mainly on the position with respect to the source-lens-observer alignment, as illustrated in Figure 12 (see also [322]). OBSERVERS LENS SOURCE (Galaxy Cluster) Strong Lensing regime (non-linear) Intermediate regime (quasilinear) Weak Lensing regime (linear) Figure 12. Illustration of the 3 regimes of lensing. In the linear regime (weak lensing), when κ ≪ 1, the distortion on one source is not noticeable. However, the lensing effect can still be detected statistically when there is a large enough sample of sources suffering a coherent lensing effect. This is the case with images of dense populations of faint galaxies present in deep optical and near-infrared observations. The images of these background galaxies are coherently sheared by a large-scale inhomogeneous matter distribution: the large scale structure of the universe. The shear field of such images is known as cosmic shear and shows angular correlations up to degree scales for sources at z ∼ 1. Cosmic shear correlations are a direct measure of the statistical properties of the lensing potential, which depend on the statistical properties of the cosmological density field and the geometrical properties of spacetime. For this reason cosmic Symmetry 2016, 8, 70 57 of 79 shear is a probe of both geometry and structure formation and it is the lensing effect we deal with in this contribution. 7.2. Cosmic Shear Observationally, a pioneering attempt to measure a cosmological weak lensing signal was made in [323], having produced an upper limit for the signal. A decade later, [324] measured shear correlations with error bars of order 100% and verified the existence of instrumental and atmospheric anisotropical contaminations. A first significant detection, cosmic variance dominated, was made in [325] on a very small field (4 sqarcmin). Finally in 2000, the first definitive detection of the cosmological shear field signal was reported almost simultaneously by four independent teams [326–329]. A large number of surveys from various ground-based observatories was soon undertaken, confirming the detections and using two-point cosmic shear correlations to constrain Ωm and σ8 parameters (see [330] for a compilation). Those surveys culminated on the Canada-France-Hawaii Telescope Legacy Survey (CFHTLS), a dedicated cosmic shear survey covering 170 deg2 that successfully detected a signal on degree scales [331]. In the past 5 years the focus of cosmic shear observations moved from the detection itself to the robust estimate of residual systematics. Currently, the state-of-the-art cosmic shear survey is the Canada-France-Hawaii Telescope Lensing Survey (CFHTLenS) [332], a reanalysis of CFHTLS where new methods have been applied to measure the galaxy distortions [333], verify the source redshift distributions [334] and remove intrinsic alignment correlations [335]. These robust data were used in various cosmological applications, including testing GR on cosmological scales, both alone and in combination with other probes, both for specific gravity theories and parametrized deviations to GR [336–338]. Cosmic shear studies deal with the light from background galaxies that propagates in the perturbed FRW spacetime [48,339,340]. In the Newtonian gauge, considering scalar perturbations and a flat universe, this spacetime is described by the metric ds2 = − (1 + 2Ψ( x, t)) dt2 + a2 (t) (1 − 2Φ( x, t)) δij dxi dx j , (173) where Ψ and Φ are the Newtonian gauge potentials. The lens equation is the solution of the equation of movement for the comoving transverse separation between the perturbed null geodesic and a reference unperturbed one. The deflection at a given point in the image may be seen as the integral of multiple local deflections over cosmological distances. The local deflection of a null geodesic at a given redshift depends on the traveling time of a light ray, and thus on the sum of the metric perturbations, (Ψ + Φ) . The resulting lens equation has a simple form in the Born approximation, i.e., when the line-of-sight integral is performed along the unperturbed geodesic, which is a good approximation for small deflections. In that case the lens equation is β i (~θ, χ) = θi − Z χ 0 dχ′ (χ − χ′ ) χ′ (Ψ + Φ),i (χθ, χ′ ) , χ (174) where χ is the comoving radial coordinate. From here we can read off the local deflections as ~α = ψ,i , defining the lensing potential as ψ(~θ, χ) = Z χ 0 dχ′
(χ − χ′ ) χ′ Ψ(χθ, χ′ ) + Φ(χθ, χ′ ) . χ (175) Note that the dependence on the sum of the two potentials, which is a relativistic feature, is responsible for the well-known factor of 2 discrepancy between General Relativity and Newtonian physics predictions for the deflection of light. The convergence and shear fields may then be computed from the second-order derivatives of the lensing potential, according to Equations (171) and (172), with the following caveat. The convergence and shear fields at an image location θ have Symmetry 2016, 8, 70 58 of 79 contributions from all light rays that may reach that point from a population of galaxies and not necessarily from only a single source at radial coordinate w. The converge then writes 1 κ (~θ ) = 2 Z χH 0   dχ p(χ) ψ,11 (~θ, χ) + ψ,22 (~θ, χ) , (176) where the integration is made over the sources, which are distributed over redshift as p(χ) dχ = p(z) dz. Cosmic shear can only be detected statistically and the most used observables are two-point correlation functions and power spectra. Considering the angular correlation function of the convergence field in the Limber approximation (i.e., neglecting contributions from pairs at large radial separation) and assuming no source clustering (i.e., assuming two independent source distributions), the convergence power spectrum may be written as 1 Pκ (ℓ) = 4 Z  4 l dχ q1 (χ) q2 (χ) PΨ+Φ χ ! |~ℓ| ,χ . χ (177) The lensing power spectrum on an angular multipole scale ℓ has contributions from the power spectra of the Laplacian of the lensing potential at various redshifts z(χ) on a redshift-dependent scale k = l/χ. The functions q include the dependences on the sources distances and distributions. Notice also that the statistical properties of the convergence and shear fields are related and in particular their power spectra are identical: Pκ (ℓ) = Pγ (ℓ). 7.3. Testing Deviations from General Relativity In order to test a cosmological model with cosmic shear data, we need to be able to compute the evolution of the lensing potential. In General Relativity, Einstein equations at perturbation level provide a set of evolution equations for the metric and a set of constraint equations that relate the metric perturbations (the potentials Ψ, Φ) to the perturbations of the cosmological fluid (matter density contrast δ, pressure perturbation δp, fluid velocity v, and anisotropic pressure σ). A combination of the 00 and 0i Einstein equations relates Φ with the comoving matter density perturbation: it is a generalized Poisson equation [341]. In Fourier space, with k2 being the spatial Laplacian, it writes k2 Φ = −4π G a2 ρ∆ , (178) where ∆ is the comoving density perturbation ∆ = δ + 3H (1 + w)v/k2 and w is the equation-of-state parameter w = p/ρ. The ij Einstein equation gives a second constraint, relating the difference Φ − Ψ with the stress σ: it is the anisotropy equation [341], k2 (Φ − Ψ) = 12π G a2 ρ (1 + w)σ. (179) The lensing convergence is determined by the Laplacian of the lensing potential, and thus by, k2 (Φ + Ψ) = −8π G a2 ρ∆. (180) Inserting this in Equation (177), the dependence of the convergence power spectrum on the matter power spectrum becomes explicit. The latter needs then to be worked out from the evolution equations, possibly with the help of Boltzmann codes and N-body simulations. Now, if we consider gravity theories that are modifications of General Relativity, the above equations are modified. Similar to the PPN formalism used in solar systems tests, it is useful to encapsulate the modifications in a set of parameters that can be computed in the various theories and can be easily constrained by experiments [342]. Various simplified versions of this so-called Parameterised Post-Friedmannian formalism have been proposed [343–346], where the Symmetry 2016, 8, 70 59 of 79 modifications are described with two independent parameters. These are the gravitational screening Q = Geff /G, which parametrizes a change in the local gravitational force on cosmological scales, and the gravitational slip η = Ψ/Φ − 1, which parametrizes an anisotropy responsible for a difference between the two potentials. The 2 parameters may be functions of scale and redshift and their values in GR are Q = 1 and η = 0. This simple relations are valid for a large number of theories of gravity. Inserting these parameters in the constraint equations Equations (178) and (179), we can combine them to get  η 2 k2 (Φ + Ψ) = −8π G Q 1 + a ρ∆. (181) 2 Inserting this in Equation (177), the convergence power spectrum may be written as, ij Pγ (ℓ) = 9 2 4 Ω H 4 m 0 Z w h 0 dw gi ( w ) g j ( w ) 2 Σ (k, a) Pδ a2 ( w )   ℓ ,w . w (182) Here it was assumed that ∆ ∼ δ and the Ωm factor comes from the matter density in Equation (180). The parameter Σ( a) = Q( a)(1 + η ( a)/2) is the combination of PPF parameters cosmic shear is sensitive to [346] . The amplitude of the power spectrum depends thus on a degenerate combination of: modified gravity via constraint equations (Σ); distances to the sources (functions g, which depend both on the knowledge of the source redshifts and on the cosmological background); matter density (Ωm h2 ); amplitude of the matter power spectrum (σ8 and structure formation in modified gravity); and external systematic effects. We applied cosmic shear data from the COSMOS survey to constrain the Σ parameter [347]. The COSMOS cosmic shear survey covers an area of 1.64 deg2 where cosmic shear angular correlations were measured on a narrow scale range going up to 20 arcmin [348]. It is however deep, with a limiting magnitude i = 26.7 and a high density of 76 galaxies/arcmin2 , allowing one to cross-correlate populations at different redshifts and perform a tomographic analysis [349]. We assume Σ is scale-independent in the narrow range of scales probed. Following [344], we also assume that the departure from GR originates an effective smooth dark energy Λ with constant w = −1, which implies the background evolution is identical to ΛCDM and the PPF parameters deviate from GR proportionally to ρDE ( a)/ρm ( a). This implies that Q ( a ) = 1 + Q a a3 , η ( a ) = η a a3 . (183) Assuming the PPF parameters remain small, the evolution of Σ also scales as Σ( a) ≈ 1 + Σ a a3 . The last item needed for the evaluation of Equation (182) is the evolution of the matter overdensities, i.e., the growth function. We consider the system of the evolution equations and constraint equations for the case of sub-horizon scales, assuming time derivatives of the perturbations are negligible compared to the spatial derivatives, and ∆ ≈ δ. The evolution equation of the cold dark matter perturbation greatly simplifies, containing only a Hubble drag term and a term on k2 Ψ [350]. In this regime, the growth parameter approximation [351] applies, where d ln(δ/a) = Ωm ( a)γ − 1 . d ln( a) (184) The growth parameter γ has the value γ = 6/11 in GR and   6 Q a + ηa Ω m γ= , 1− 11 2 1 − Ωm (185) in PPF-parameterized modified gravity models [352]. With these assumptions, we attempt to constrain the parameters using only the cosmic shear data in a Bayesian analysis. The cosmic shear data consists of 15 two-point cross-correlation functions, between 6 galaxy redshift bins. All 15 functions have similar shapes but different amplitudes Symmetry 2016, 8, 70 60 of 79 depending on their redshifts. Notice that our parametrization only affects the amplitude of the cosmic shear functions, since we work in the ansatz that Σ and γ are scale-independent. This implies that our constraint will necessarily be a degenerate combination of the 2 modified gravity parameters, the redshift and the standard cosmological parameters affecting the amplitude of the power spectrum. For this reason, we need to use either additional data to constrain the standard cosmology or to introduce strong priors. We introduce a prior in the form of importance sampling [353]. In particular, we take WMAP7 ΛCDM Monte Carlo Markov Chains [354] and evaluate a cube of models in a ( Q a , ηa , f z ) space for each point of the MCMC chain. In this way, the ΛCDM MCMC weights will be changed as a function of the cosmic shear likelihoods in the cubes of models. The ( Q a , ηa ) likelihoods are marginalized over the redshift-calibration nuisance parameter f z that accounts for the source redshift uncertainty. The shear correlation functions are computed using the transfer function of Eisenstein & Hu [355] with the modified linear growth γ and the non-linear halofit prescription [356]. Besides addressing the degeneracies with the standard cosmological parameters and the redshift uncertainties, we also minimise contamination from intrinsic alignments by not using cosmic shear auto-correlation functions in the analysis. The result obtained in this analysis is compatible with GR. We find the following 1-σ constraint for the following combination of the 2 parameters: Σ + 2.3 (γ − 0.55) = 0.99 ± 0.31 . (186) 7.4. Testing Gravity with Future Cosmic Shear Data To take full advantage of future high-precision cosmic shear data, accurate modelling is needed at several levels, and the analyses become much more complex than the one described above. When studying specific modified gravity theories, structure formation needs to be worked out for those models, including non-linear scales, which calls for model-specific perturbative analytical approaches [357] or N-body simulations [50]. Numerical simulations are also required for other purposes. Lensing simulations are needed to estimate the covariance of lensing observables, including cosmic variance. Hydrodynamical simulations are increasingly needed to model various baryonic effects on lensing observables, such as supernova and AGN feedback, star formation or radiative cooling [51]. Working in a model-independent approach requires detailed modelling of the deviations in Einstein equations. Besides applying the PPF parameters, it is also interesting to parametrize the deviations at the observables level. [358] parametrizes deviations of the lensing power spectrum from the non-linear GR matter power spectrum. The latter can be computed from the linear one with a fitting formula [359] or from templates [360]. Scale and time dependencies of the parametrized deviations may also be probed in bins by performing a principal components analysis [361]. Cosmic shear probes deviations from GR through the integrated effect of Ψ and Φ on null geodesics. It has a high signal-to-noise on the scales where the deviations are maximal and is thus more promising than other probes of Φ + Ψ, such as the integrated Sachs-Wolfe effect [362]. However, to be able to break the degeneracy between the two potentials in the modified gravitational effects, cosmic shear needs to be combined with a local probe of matter overdensities, which are sensitive to Ψ alone. One possible combination of cosmic shear and redshift-space distortions is the so-called EG test [363,364], which has the additional advantage of being independent of galaxy bias. A first detection of the EG signal was made in [365]. The combination of weak lensing and galaxy clustering observables is an important asset of the forthcoming Euclid space mission to be launched in 2020, which combines a weak lensing and galaxy clustering wide survey and aims at measuring the rate of cosmic structure growth to a 2% precision [49]. Symmetry 2016, 8, 70 61 of 79 8. Angular Distribution of Cosmological Parameters as a Measurement of Spacetime Inhomogeneities On large scales, the matter distribution is statistically homogeneous. Locally however, matter is distributed according to a pattern of alternate overdensed regions and underdensed regions. Since averaging inhomogeneities in the matter density distribution yields a homogeneous description of the Universe, then the apparent homogeneity of the cosmological parameters could also result from the averaging of inhomogeneities in the cosmological parameters, which would reflect the inhomogeneities in the density distribution. In the context of backreaction models, angular variations in the parameters could also source a repulsive force and potentially emulate cosmic acceleration. Hence the reasoning was to look for these inhomogeneities across the sky and then to use an adequate toy model to compute the magnitude of the acceleration derived from angular variations of the parameters compared to the acceleration driven by a cosmological constant. Recently, a method was proposed to search for inhomogeneities based on a local estimation of cosmological parameters [71]. We demonstrated the method with supernovae (SNe), since they are a sensitive probe of the late–time expansion history of the Universe by means of the dependence of the luminosity distance on both the geometry and the growth of structure. Hence the parameters to be estimated are those that affect the luminosity distance, namely, {Ω M , ΩΛ , H0 }. Supernovae have previously been used to constrain inhomogeneities in the cosmic expansion, namely, by measuring the hemispherical anisotropy of H0 [366], by assessing the dependence of H0 with the position on the cosmic web [367] and by mapping { H0 , q0 } in the sky [368]. Supernovae have also been used to constrain inhomogeneities in the dark energy, namely by constraining fluctuations of a dynamical dark energy [369] and by constraining radial inhomogeneity for a Lemaître–Tolman–Bondi metric with a cosmological constant [370]. We used the type Ia supernova sample compiled by the SNLS–SDSS collaborative effort Joint Light–curve Analysis (JLA) [371] , totalling NSNe = 740 SNe with redshift 0.01 ≤ z ≤ 1.30 and distributed on the sky according to Figure 13. The sample was obtained from http://supernovae.in2p3.fr/sdss_snls_jla/ReadMe.html. Number of SNe per pixel SN in JLA Catalogue 14 15 38 6 7 380 10 6 95 26 73 70 6 380 Figure 13. Left panel: Angular distribution of the JLA type Ia SN sample. The sample consists of type Ia SNe compiled from different surveys as described in Ref. [371] and plotted in Celestial coordinates. Right panel: Angular distribution of the number of type Ia SNe per pixel. The total SN sample is divided into subsamples grouped by pixels in the HEALPix pixelation scheme. The number of SNe in each pixel is printed in the centre of the pixel. The pixel indices {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12} correspond to the SN counts {15, 14, 95, 7, 380, 38, 70, 26, 10, 6, 6, 73}, as indicated by the colour bar in logarithmic scale. In both panels, the dashed line represents the Galactic equator. 8.1. Method of Local Parameter Estimation Our method consists in dividing the SNe over a pixelated map of the sky, with pixels of equal surface area according to the HEALPix pixelation scheme [372]. The SNe subsample that falls in Symmetry 2016, 8, 70 62 of 79 each pixel k is used to estimate the cosmological parameters in that pixel, thus the designation of local estimation. For each cosmological parameter being estimated, there results a map with the same pixelation as that of the SNe subsamples. Each pixel is assumed to be described by a Friedmann–Lemaître–Roberston–Walker metric so that the full sky is an inhomogeneous ensemble of disjoint, locally homogeneous regions. In each pixel, we estimated the cosmological parameters that minimise the chi–square of the fit of the theoretical distance modulus µtheo (z; Ω M , ΩΛ , H0 ) = 5 log[d L (z; Ω M , ΩΛ , H0 )] + 25 , (187) (d L in units of Mpc) computed for the trial values of {Ω M , ΩΛ , H0 }, to the measured distance modulus µdata (z) = m B (z) − ( MB − αx1 + βc) , (188) computed for the estimated light–curve parameters {m B , x1 , c} for each SN. From the local estimation there resulted the maps xik = {Ω M , ΩΛ , H0 }k . From the estimation using the complete SN sample, fid fid designated as the global estimation, we obtained the values xifid = {Ωfid M , ΩΛ , H0 }, which we adopt as fiducial values for the parameters. The fiducial values for the estimated parameters are: fid fid {Ωfid M , ΩΛ , H0 } = {0.256, 0.715, 71.17}. For the derived parameters Ωκ = 1 − Ω M − ΩΛ and fid q0 = Ω M /2 − ΩΛ , the fiducial values are Ωfid κ = 0.029 and q0 = −0.586. For better visualization, we compare the local estimation with the fiducial values by defining the difference maps ∆xik = xik − xifid . Since the SN sample used is a collection from different SN surveys, the SN subsampling in pixels is highly inhomogeneous. As a result, well sampled pixels return robust estimations, whereas poorly sampled pixels return weak constraints. This generates a noise bias. We hypothesised that the noise bias due to the inhomogeneity in the SN subsampling is modelled by the local estimation obtained from randomizing the dependence of redshift with position and then averaging over the different randomizations. For each randomization in z, the SN positions in the sky remain the same as in the original sample. Subtracting the noise bias off the maps previously estimated we obtained unbiased unbias ≡ x − x̄bias + xfid , where x̄bias ≡ xbias ′ maps xik ik i ik ik ij′ k j′ is the average over j randomizations in z. For unbias ≡ xunbias − xfid = x − x̄bias . better visualization, we define the unbiased difference maps as ∆xik ik i ik ik In Figure 14 we show the values of the difference maps at each pixel, before (left panel) and after (right panel) the noise bias subtraction. We also show the corresponding errors. We observe a negative correlation of the number of SNe in each pixel both with the fluctuations about the fiducial values, and with the corresponding errors. These results reflect the poor sampling derived from an inhomogeneous coverage of the sky, both in angular and redshift space, by the SN surveys. A poor sampling into pixels implies a poor sampling in redshift, which limits the constraint of the degeneracies among the parameters and hence compromises the parameters estimation. We also observe that, after the bias subtraction, the fluctuations about the fiducial values are reduced, albeit the errors increase. Consequently, for the correlation of the number of SNe per pixel, the unbiased maps in comparison with the original maps show a smaller correlation with the difference values and a larger correlation with the errors, as expected. These results indicate bias as a measure of the noise bias due to the inhomogeneity of the SN sampling. the validity of x̄ik Comparing the unbiased difference maps with the fiducial values, we measure fluctuations of order 5%–95% for Ω M , 1%–25% for ΩΛ and up to 5% for H0 . Symmetry 2016, 8, 70 No. SNe 400 350 300 250 200 150 100 50 0 ∆Ωi 63 of 79 0.3 0.2 0.1 0.0 0.1 0.2 0.3 Correlation with no. SNe per pixel σ∆Ωi 0.20 0.15 0.10 0.05 ∆ H0 0.00 5 4 3 2 1 0 1 2 2.5 σ∆H 0 2.0 1.5 1.0 0.5 0.0 2 4 6 8 10 12 No. SNe 400 350 300 250 200 150 100 50 0 ∆Ωi Pixel index 0.3 0.2 0.1 0.0 0.1 0.2 0.3 Correlation with no. SNe per pixel: Unbiased σ∆Ωi 0.20 0.15 0.10 0.05 ∆ H0 0.00 5 4 3 2 1 0 1 2 σ∆H 0 2.5 2.0 1.5 1.0 0.5 0.0 2 4 6 8 10 12 Pixel index Figure 14. Correlation of the number of SNe per pixel with the estimated parameters’ fluctuations about the fiducial values, before (left panel) and after (right panel) the noise bias removal. At each pixel we plot: a) in the first panel, the number of SNe; b) in the second panel, the value of the parameters difference {∆Ω M , ∆ΩΛ } respectively as the solid red and solid blue lines, with the dashed black line marking the zero difference; c) in the third panel, the standard deviation {σ∆Ω M , σ∆ΩΛ } respectively as solid red and solid blue lines; d) in the fourth plot, the value of ∆H0 as the solid black line, with the dashed black line marking the zero difference; e) in the fifth panel, the standard deviation σ∆H0 as the solid black line. 8.2. Average over the Local Parameter Estimation When averaging over the homogeneous regions, angular fluctuations in the expansion factor (and consequently in H0 ) induce a backreaction term in the form of an extra positive acceleration [373]. We derived the analytical extra positive acceleration for a toy model of an arbitrary number of disjoint, homogeneous regions and computed the overall deceleration parameter assuming a) no backreaction and b) backreaction for the measured angular distribution of H0 . Symmetry 2016, 8, 70 64 of 79 In the absence of backreaction, the averaging consists in taking the mean weighted by the variance’s inverse wk = 1/Var[ x̄ik ] of parameter xi in each pixel k, Npixel x̄i = x̄ik ∑k = k wk x̄ik Npixel ∑k . (189) wk For the estimated parameters, we find {Ω M , ΩΛ , H 0 } = {0.242, 0.607, 71.06} ± {0.099, 0.115, 0.87}. These values are consistent with, but systematically smaller than, the fiducial values. After the noise bias removal, we find {Ω M , ΩΛ , H 0 }unbias = {0.294, 0.679, 71.37} ± {0.115, 0.120, 0.94}. These values are consistent with the fiducial values. For the derived parameters, we find Ωκ = 0.078 ± 0.131 and q̄0 = −0.451 ± 0.159, and after the noise bias removal we find Ωκ,unbias = 0.002 ± 0.135 and q̄0,unbias = −0.526 ± 0.172. After the noise bias removal, the pixel average values of both Ωκ and q0 become closer to the corresponding values estimated using the complete sample. In the presence of backreaction, the averaging consists in taking the mean weighted by the 3–volume Vk of each pixel k, Npixel x̄i = x̄ik Vk ∑k = Vk x̄ik Npixel ∑k . (190) Vk Since all pixels have the same surface area and hence the angular directions expand the same way today, we can assume that the radial direction also expands the same way today, which implies that the 3–volume is the same today for all pixels. This amounts to identifying a volume as a pixel. Defining vk = Vk / ∑k Vk , then for Npixel disjoint regions, the average of the Hubble parameter is given by Npixel ∑ h H0 iVk = vk H0,k . (191) k Taking the derivative, it follows that   ä = a Vk Npixel ∑ k ä vk + 2 a Npixel Npixel ∑ ∑ l >k k vk vl ( H0,k − H0,l )2 , (192) which decomposes into a linear term in the pixel average of ( ä/a) and a quadratic term in differences of H0 between pairs of pixels. The quadratic (backreaction) term generates an acceleration that is due not to regions speeding up locally, but instead to the slower regions becoming less represented in the average. Note that, in the absence of the quadratic term, the volume average reduces to the pixel average above. Then the volume average of q0 becomes Npixel hq0 iVk = ∑ k vk q0,k −  Npixel Npixel 2 Npixel ∑k vk H0,k 2 ∑ ∑ k l >k vk vl ( H0,k − H0,l )2 . (193) Using the fluctuations in H0 measured in the local parameter estimation, we find q̄0 = −0.452 ± 0.159 and after the noise bias removal we find q̄0,unbias = −0.527 ± 0.173. The volume average of the other parameters is the same as the pixel average. After the noise bias removal, the volume average values become closer to the value for the deceleration estimated from the complete bias as a measure of the sample. The results from both averaging methods strengthen the validity of x̄ik noise bias due to the inhomogeneity of the SN sampling. The quadratic term is of order 10−3 (or lower) times the linear term, which means that the error of the difference is below the standard deviation. Hence for the angular fluctuations in H0 measured Symmetry 2016, 8, 70 65 of 79 with this SN sample, the contribution of the quadratic term in Equation (192) is insignificant, which renders the volume averaging equivalent to the pixel averaging. Hence in the context of this toy model of an inhomogeneous spacetime, backreaction is not a viable dynamical mechanism to emulate cosmic acceleration. A subsequent study using a kinematic parametrisation of the luminosity distance expressed in terms of time derivatives of the scale factor yielded concordant results [374]. 9. Summary and Conclusion In this work, we have explored the dynamics and evolution of the Universe at early and late times, focusing on both dark energy and extended gravity models and their astrophysical and cosmological consequences. More specifically, we presented scalar-tensor theories, focussing on a brief review of the general formalism, the conformal picture, discussed the experimental and observational status of the theory, presented a detailed dynamical system analysis of the cosmological dynamics, and concluded with a quintessence scenario where a multi-scalar-tensor theory is responsible for the present accelerated expansion of the Universe. Furthermore, we considered Horndeski theories and the self-tuning properties, and explored realistic cosmologies that can be constructed, by studying the properties of this intriguing theory. We also reviewed f ( R) modified theories of gravity and extensions, where the foundational questions and astrophysical/cosmological issues of the Palatini approach were extensively explored; we also focussed on recently proposed theories such as on curvature-matter couplings and the hybrid metric-Palatini theory. The inflationary epoch was also presented in f ( R) gravity, where the so-called Starobinsky inflation is known as one of the more reliable candidates for explaining the inflationary paradigm as provided by the recent Planck/Bicep2 constraints. Moreover, we have reviewed the evolution and cosmological consequences of topological defect networks that may be formed in the early universe. In particular, we focused on two of the most significant signatures of these networks: the cosmic microwave background anisotropies and the stochastic gravitational wave background. As mentioned in the Introduction, due to the extremely important aspect of the synergy between theory and observations, the second part of this review was dedicated to observational cosmology. In particular, cosmological tests with galaxy clusters at CMB frequencies were presented. More specifically, galaxy clusters and the thermal and kinetic Sunyaev-Zel’dovich (SZ) effects were reviewed, and the cosmological tests using SZ cluster surveys were presented, focussing on SZ cluster counts, SZ power spectra and the possibility of probing new physics with the SZ clusters; a review on probing primordial non-gausianities with galaxy clusters was also presented, with an emphasis on the parametrization of primordial non-Gaussianities, the non-Gaussian mass function, biased cosmological parameter estimation with cluster counts, the impact of primordial non-Gaussianities on galaxy clusters scaling relations. Furthermore, gravitational lensing was also explored, with a focus on cosmological scales, with the possibility of testing deviations from General Relativity, in particular, with future cosmic shear data. In conclusion, the angular distribution of cosmological parameters as a measurement of spacetime inhomogeneities was presented, using a method of averaging over the local parameter estimation. Thus, cosmological tests were discussed, and their relevance in constraining our cosmological description of the Universe. Acknowledgments: We thank Ippocratis Saltas for helpful comments and suggestions. This work was supported by Fundação para a Ciência e a Tecnologia (FCT) through the research grant UID/FIS/04434/2013. PPA acknowledges financial support from FCT through the Investigador FCT Contract No. IF/00863/2012, funded by FCT/MCTES (Portugal). CSC is funded by FCT, Grant No. SFRH/BPD/65993/2009. FSNL acknowledges financial support from FCT through the Investigador FCT Contract No. IF/00859/2012, funded by FCT/MCTES (Portugal). PMM acknowledges financial support from the Spanish Ministry of Economy and Competitiveness through the postdoctoral training contract FPDI- 2013-16161, and the project FIS2014-52837-P. DRG acknowledges support from FCT, Grant No. SFRH/BPD/102958/2014. DSG acknowledges support from FCT, Grant No. SFRH/BPD/95939/2013. LS is supported by FCT through the grant SFRH/BPD/76324/2011. IT acknowledges support from FCT through the Investigador FCT Contract No. IF/01518/2014 and POPH/FSE (EC) by FEDER funding through the program Programa Operacional de Factores de Competitividade – Symmetry 2016, 8, 70 66 of 79 COMPETE. The authors also acknowledge the COST Action CA15117, supported by COST (European Cooperation in Science and Technology). Author Contributions: All the authors have substantially contributed to the present work. Conflicts of Interest: The authors declare no conflict of interest. References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22. 23. 24. 25. 26. 27. Perlmutter, S.; Aldering, G.; Goldhaber, G.; Knop, R.A.; Nugent, P.; Castro, P.G.; Deustua, S.; Fabbro, S.; Goobar, A.; Groom, D.E.; et al. Measurements of Omega and Lambda from 42 high redshift supernovae. Astrophys. J. 1999, 517, 565–586, doi:10.1086/307221. Riess, A.G.; Filippenko, A.V.; Challis, P.; Clocchiattia, A.; Diercks, A.; Garnavich, P.M.; Craig, G.; Hogan, C.J.; Jha, S.; Kishner, R.P.; et al. Observational evidence from supernovae for an accelerating universe and a cosmological constant. Astron. J. 1998, 116, 1009–1038, doi:10.1086/300499. Maartens, R. Brane world gravity. Living Rev. Rel. 2004, 7, 7. Dvali, G.R.; Gabadadze, G.; Porrati, M. 4-D gravity on a brane in 5-D Minkowski space. Phys. Lett. B 2000, 485, 208. de Rham, C.; Hofmann, S.; Khoury, J.; Tolley, A.J. Cascading Gravity and Degravitation. J. Cosmol. Astropart. Phys. 2008, 0802, 011. Sotiriou, T.P.; Faraoni, V. f(R) Theories of Gravity. Rev. Mod. Phys. 2010, 82, 451. De Felice, A.; Tsujikawa, S. f(R) theories. Living Rev. Rel. 2010, 13, 3. Capozziello, S.; De Laurentis, M. Extended Theories of Gravity. Phys. Rept. 2011, 509, 167. Nojiri, S.; Odintsov, S.D. Unified cosmic history in modified gravity: From F(R) theory to Lorentz non-invariant models. Phys. Rept. 2011, 505, 59. Lobo, F.S.N. The Dark side of gravity: Modified theories of gravity. In Dark Energy-Current Advances and Ideas; Research Signpost: Kerala, India, 2009; pp. 173–204. Capozziello, S. Curvature quintessence. Int. J. Mod. Phys. D 2002, 11, 483. Nojiri, S.; Odintsov, S.D. Modified gravity with negative and positive powers of the curvature: Unification of the inflation and of the cosmic acceleration. Phys. Rev. D 2003, 68, 123512. Carroll, S.M.; Duvvuri, V.; Trodden, M.; Turner, M.S. Is cosmic speed-up due to new gravitational physics? Phys. Rev. D 2004, 70, 043528. Copeland, E.J.; Sami, M.; Tsujikawa, S. Dynamics of dark energy. Int. J. Mod. Phys. D 2006, 15, 1753–1936. Brans, C.; Dicke, R.H. Mach’s principle and a relativistic theory of gravitation. Phys. Rev. 1961, 124, 925–935. De Rham, C. Massive Gravity. Living Rev. Rel. 2014, 17, 7. Horndeski, G.W. Second-order scalar-tensor field equations in a four-dimensional space. Int. J. Theor. Phys. 1974, 10, 363. Charmousis, C.; Copeland, E.J.; Padilla, A.; Saffin, P.M. Self-tuning and the derivation of a class of scalar-tensor theories. Phys. Rev. D 2012, 85, 104040. Martín-Moruno, P.; Nunes, N.J.; Lobo, F.S.N. Horndeski theories self-tuning to a de Sitter vacuum. Phys. Rev. D 2015, 91, 084029. Nojiri, S.; Odintsov, S.D. Where new gravitational physics comes from: M Theory? Phys. Lett. B 2003, 576, 5. Parker, L.; Toms, D.J. Quantum Field Theory in Curved Spacetime: Quantized Fields and Gravity; Cambridge University Press: Oxford, UK, 2009. Cembranos, J.A.R. Dark matter from R2 gravity. Phys. Rev. Lett. 2009, 102, 141301. Bamba, K.; Odintsov, S.D. Inflationary cosmology in modified gravity theories. Symmetry 2015, 7, 220. Ade, P.A.R.; Aghanim, A.N.; Armitage-Caplan, C.; Ashdown, A.M.; Atrio-Barandela, F.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; Bartlett, J.G.; et al. Planck 2013 results. XXII. Constraints on inflation. Astron. Astrophys. 2014, 571, A22. Hu, W.; Sawicki, I. Models of f(R) Cosmic Acceleration that Evade Solar-System Tests. Phys. Rev. D 2007, 76, 064004. Nojiri, S.; Odintsov, S.D. Unifying inflation with LambdaCDM epoch in modified f(R) gravity consistent with Solar System tests. Phys. Lett. B 2007, 657, 238. Olmo, G.J. Palatini Approach to Modified Gravity: f ( R) Theories and Beyond. Int. J. Mod. Phys. D 2011, 20, 413. Symmetry 2016, 8, 70 28. 29. 30. 31. 32. 33. 34. 35. 36. 37. 38. 39. 40. 41. 42. 43. 44. 45. 46. 47. 48. 49. 50. 51. 52. 53. 54. 55. 67 of 79 Sotiriou, T.P.; Liberati, S. Metric-affine f ( R) theories of gravity. Annals Phys. 2007, 322, 935. Harko, T.; Koivisto, T.S.; Lobo, F.S.N.; Olmo, G.J. Metric-Palatini gravity unifying local constraints and late-time cosmic acceleration. Phys. Rev. D 2012, 85, 084016. Linder, E.V. Einstein’s Other Gravity and the Acceleration of the Universe. Phys. Rev. D 2010, 81, 127301. Cai, Y.F.; Capozziello, S.; De Laurentis, M.; Saridakis, E.N. f ( T ) Teleparallel Gravity and Cosmology. 2015, arXiv:1511.07586 [gr-qc]. Harko, T.; Lobo, F.S.N.; Otalora, G.; Saridakis, E.N. Nonminimal torsion-matter coupling extension of f ( T ) gravity. Phys. Rev. D 2014, 89, 124036. Harko, T.; Lobo, F.S.N.;Otalora, G.; Saridakis, E.N. f ( T, T ) gravity and cosmology. J. Cosmol. Astropart. Phys. 2014, 1412, 021. van Dam, H.; Veltman, M.J.G. Massive and massless Yang-Mills and gravitational fields. Nucl. Phys. B 1970, 22, 397. Zakharov, V. I. Linearized gravitation theory and the graviton mass. JETP Lett. 1970, 12, 312. Avelino, P.P.; Martins, C.J.A.P.; Menezes, J.; Menezes, R.; Oliveira, J.C.R.E. Frustrated expectations: Defect networks and dark energy. Phys. Rev. D 2006, 73, 123519. Ade, P.A.R.; Aghanim, N.; Arnaud, M.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; Bartlett, J.G.; Bartolo, N.; et al. [Planck Collaboration]. Planck 2015 results. XIII. Cosmological parameters. 2015, arXiv:1502.01589 [astro-ph.CO]. Avelino, P.P.; Sousa, L. Observational Constraints on Varying-Alpha Domain Walls. Universe 2015, 1, 6–16. Anthonisen, M.; Brandenberger, R.; Scott, P. Constraints on cosmic strings from ultracompact minihalos. Phys. Rev. D 2015, 92, 023521. Barton, A.; Brandenberger, R.H.; Lin, L. Cosmic Strings and the Origin of Globular Clusters. J. Cosmol. Astropart. Phys. 2015, 1506, 022. Bramberger, S.F.; Brandenberger, R.H.; Jreidini, P.; Quintin, J. Cosmic String Loops as the Seeds of Super-Massive Black Holes. J. Cosmol. Astropart. Phys. 2015, 1506, 007. Avelino, P.P.; Liddle, A.R. Cosmological perturbations and the reionization epoch. Mon. Not. Roy. Astron. Soc. 2004, 348, 105. Weinberg, D.H.; Mortonson, M.J.; Eisenstein, D.J.; Hirata, C.; Riess, A.G.; Rozo, E. Observational probes of cosmic acceleration. Phys. Rep. 2013, 530, 87–255. Ade, P.A.R.; Aghanim, N.; Arnaud, M.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Barreiro, R.B.; Bartlett, J.G.; Bartolo, E.; Battaner, R.; et al. Planck Collaboration, Planck 2015 results. XIII. Cosmological parameters. 2015, arXiv:1502.01589. Aghanim, N.; Majumdar, S.; Silk, J. Secondary anisotropies of the CMB. Rep. Prog. Phys. 2008, 71, 066902. Sunyaev, R.A.; Zeldovich, Y.B. The Sunyaev-Zel’dovich effect. Comments Astrophys. Space Phys. 1972 4, 173. Birkinshaw, M. The Sunyaev-Zel’dovich effect. Phys. Rep. 1999, 310, 97. Bartelmann, M.; Schneider, P. Weak Gravitational Lensing. Phys. Rep. 2001, 340, 297–472. Laureijs, R.; Amiaux, J.; Arduini, S.; Brinchmann, J.; Cole, R.; Cropper, M.; Dabin, C.; Duvet, L.; Ealet, A.; Garilli, B.; et al. Euclid Definition Study Report. 2011, arXiv:1110.3193. Winther, H.A.; Schmidt, F.; Barreira, A.; Arnold, C.; Sownak Bose, S.; Claudio Llinares, C.; Marco Baldi, M.; Falck, B.; Hellwing, W.A.; Koyama, K.; et al. Modified gravity N-body code comparison project. Mon. Not. Roy. Astron. Soc. 2015, 454, 4208–4234. Semboloni, E.; Hoekstra, H.; Schaye, J.; van Daalen, M.P.; McCarthy, I.J. Quantifying the effect of baryon physics on weak lensing tomography. Mon. Not. Roy. Astron. Soc. 2011, 417, 2020–2035. Bull, P.; Akrami, Y.; Adamek, J.; Baker, T.; Bellini, E.; Beltran-Jimenez, J.; Bentivegna, E.; Camera, S.; Clesse, S.; Davis, J.H. Beyond Λ CDM: Problems, solutions, and the road ahead. Phys. Dark Universe 2016, 12, 56–99. Joudaki, S.; Blake, C.; Heymans, C.; Choi, A.; Harnois-Deraps, J.; Hildebrandt, H.; Joachimi, B.; Johnson, A.; Mead, A.; Parkinson, D.; et al. CFHTLenS revisited: Assessing concordance with Planck including astrophysical systematics. 2016, arXiv: 1601.05786. Schwarz, D.; Copi, C.; Huterer, D.; Starkman, G. CMB Anomalies after Planck. 2015, arXiv: 1510.07929. Acquaviva, V.; Bartolo, N.; Matarrese, S.; Riotto, A. Gauge-invariant second-order perturbations and non-Gaussianity from inflation. Nucl. Phys. B 2003, 667, 119. Symmetry 2016, 8, 70 56. 57. 58. 59. 60. 61. 62. 63. 64. 65. 66. 67. 68. 69. 70. 71. 72. 73. 74. 75. 76. 77. 78. 79. 80. 81. 82. 83. 84. 68 of 79 Maldacena, J. Non-gaussian features of primordial fluctuations in single field inflationary models. J. High Energy Phys. 2003, 5, 13. Bartolo, N.; Komatsu, E.; Matarrese, S.; Riotto, A. Non-Gaussianity from inflation: Theory and observations. Phys. Rep. 2003, 402, 103. Komatsu, E.; Afshordi, N.; Bartolo, N.; Baumann, D.; Bond, J.R.; Buchbinder, E.I.; Byrnes, C.T.; Chen, X.; Chung, D.J.H.; Cooray, A.; et al. Non-Gaussianity as a Probe of the Physics of the Primordial Universe and the Astrophysics of the Low Redshift Universe. 2003, arXiv:0902.4759 [astro-ph.CO]. Lyth, D.H.; Wands, D. Generating the curvature perturbation without an inflaton. Phys. Lett. B 2002, 524, 5. Bartolo, N.; Matarrese, S.; Riotto, A. Non-Gaussianity in the curvaton scenario. Phys. Rev. D 2004, 69, 043503. Sasaki, M.; Väliviita, J.; Wands, D. Non-Gaussianity of the primordial perturbation in the curvaton model. Phys. Rev. D 2006, 74, 103003. Lehners, J.-L.; Steinhardt, P.J. Non-Gaussian density fluctuations from entropically generated curvature perturbations in ekpyrotic models. Phys. Rev. D 2008, 77, 063533. Lehners, J.-L. Ekpyrotic Nongaussianity: A Review. Adv. Astron. 2010, 903907. Yokoyama, S.; Soda, J. Primordial statistical anisotropy generated at the end of inflation. J. Cosmol. Astropart. Phys. 2008, 8, 005. Karčiauskas, M.; Dimopoulos, K.; Lyth, D.H. Anisotropic non-Gaussianity from vector field perturbations. Phys. Rev. D 2009, 80, 023509. Dimastrogiovanni, E.; Bartolo, N.; Matarrese, S.; Riotto, A. Non-Gaussianity and Statistical Anisotropy from Vector Field Populated Inflationary Models. Adv. Astron. 2010, 752670, doi:10.1155/2010/752670. Polarski, D.; Starobinsky, A.A. Isocurvature perturbations in multiple inflationary models. Phys. Rev. D 1994, 50, 6123. Mukhanov, V.F.; Steinhardt, P.J. Density perturbations in multifield inflationary models. Phys. Lett. B 1998, 422, 52. Rigopoulos, G.I.; Shellard, E.P.S.; van Tent, B.J.W. Large non-Gaussianity in multiple-field inflation. Phys. Rev. D 2006, 73, 083522. Byrnes, C.T.; Choi, K.-Y. Review of Local Non-Gaussianity from Multifield Inflation. Adv. Astron. 2010, 724525. Carvalho C.S.; Marques, K. Angular distribution of cosmological parameters: Measurement of inhomogeneities from type Ia supernovae. 2015, arXiv: 1512.07869 [astro-ph/CO]. Billyard, A.; Coley, A. On the correspondence between exact solutions in Kaluza-Klein theory and in scalar tensor theories. Mod. Phys. Lett. A 1997, 12, 2121. Barrow, J.D.; Tipler, F.J. The Anthropic Cosmological Principle; Cambridge University Press: Cambridge, UK, 1988. Wands, D. Extended gravity theories and the Einstein-Hilbert action. Class. Quant. Grav. 1994, 11, 269. Barrow, J.D.; Cotsakis, S. Inflation and the Conformal Structure of Higher Order Gravity Theories. Phys. Lett. B 1988, 214, 515. Kalara, S.; Kaloper, N.; Olive, K.A. Theories of Inflation and Conformal Transformations. Nucl. Phys. B 1990, 341, 252. Malquarti, M.; Copeland, E.J.; Liddle, A.R.; Trodden, M. A New view of k-essence. Phys. Rev. D 2003, 67, 123503. Bergmann, P.G. Comments on the scalar tensor theory. Int. J. Theor. Phys. 1968, 1, 25. Wagoner, R.V. Scalar tensor theory and gravitational waves. Phys. Rev. D 1970, 1, 3209. Nordtvedt, N. PostNewtonian metric for a general class of scalar tensor gravitational theories and observational consequences. Astrophys. J. 1970, 161, 1059. Will, C.M. The Confrontation between General Relativity and Experiment. Living Rev. Rel. 2014, 17, 4. Wald, R.M., General Relativity; University of Chicago Press: Chicago, IL, USA, 1984; p. 491. Damour, T.; Nordtvedt, K. Tensor-scalar cosmological models and their relaxation toward general relativity. Phys. Rev. D 1993, 48, 3436. Damour, T.; Esposito-Farèse, G. Tensor multiscalar theories of gravitation. Class. Quantum Grav. 1992, 9, 2093. Symmetry 2016, 8, 70 85. 86. 87. 88. 89. 90. 91. 92. 93. 94. 95. 96. 97. 98. 99. 100. 101. 102. 103. 104. 105. 106. 107. 108. 109. 110. 111. 112. 113. 114. 115. 116. 117. 118. 119. 69 of 79 Thorne, K.S.; Will, C.M. Theoretical Frameworks for Testing Relativistic Gravity. I. Foundations. Astrophys. J. 1971, 163, 595. Mimoso, J.P.; Nunes, A.M. General relativity as a cosmological attractor of scalar tensor gravity theories. Phys. Lett. A 1998, 248, 325. Banerjee, N.; Sen, S. Does Brans-Dicke theory always yield general relativity in the infinite omega limit? Phys. Rev. D 1997, 56, 1334. Billyard, A.; Coley, A.; Ibanez, J. On the asymptotic behavior of cosmological models in scalar tensor theories of gravity. Phys. Rev. D 1999, 59, 023507. Mimoso, J.P.; Nunes, A. A Qualitative Analysis of the Attractor Mechanism of General relativity. Astrophys. Space Sci. 2003, 283, 661. Faraoni, V.; Gunzig, E. Einstein frame or Jordan frame? Int. J. Theor. Phys. 1999, 38, 217. Flanagan, E.E. The Conformal frame freedom in theories of gravitation. Class. Quant. Grav. 2004, 21, 3817. Olmo, G.J. Violation of the Equivalence Principle in Modified Theories of Gravity. Phys. Rev. Lett. 2007, 98, 061101. Capozziello, S.; de Ritis, R.; Marino, A.A. Some aspects of the cosmological conformal equivalence between “Jordan frame” and “Einstein frame”. Class. Quant. Grav. 1997, 14, 3243. Capozziello, S.; Lobo, F.S.N.; Mimoso, J.P. Energy conditions in modified gravity. Phys. Lett. B 2014, 730, 280. Capozziello, S.; Lobo, F.S.N.; Mimoso, J.P. Generalized energy conditions in Extended Theories of Gravity. Phys. Rev. D 2015, 91, 124019. Dicke, R.H. Mach’s principle and invariance under transformation of units. Phys. Rev. 1962, 125, 2163. Nariai, H. On the Green’s function in an expanding universe and its role in the problem of Mach’s principle. Progr. Theor. Phys. 1968, 40, 49. O’Hanlon, J.; Tupper, B.O.J. Vacuum-field solutions in the Brans-Dicke theory. Il Nuovo Cimento 1972, 7, 305. Barrow, J.D.; Mimoso, J.P. Perfect fluid scalar-tensor cosmologies. Phys. Rev. D 1994, 50, 3746. Barrow, J.D. Scalar-tensor cosmologies. Phys. Rev. D 1993, 47, 5329. Mimoso, J.P.; Wands, D. Massless fields in scalar—Tensor cosmologies. Phys. Rev. D 1995, 51, 477. Nunes, A.; Mimoso, J.P. On the potentials yielding cosmological scaling solutions. Phys. Lett. B 2000, 488, 423. Charters, T.C.; Nunes, A.; Mimoso, J.P. Stability analysis of cosmological models through Liapunov’s method. Class. Quant. Grav. 2001, 18, 1703. Mimoso, J.P.; Nunes, A. General relativity as an attractor to scalar tensor gravity theories. Astrophys. Space Sci. 1999, 261, 327. Nunes, A.; Mimoso, J.P.; Charters, T.C. Scaling solutions from interacting fluids. Phys. Rev. D 2001 63, 083506. Santos, C.; Gregory, R. Cosmology in Brans-Dicke theory with a scalar potential. Ann. Phys. 1997, 258, 111. Holden, D.J.; Wands, D. Phase-plane analysis of Friedmann-Robertson-Walker cosmologies in Brans-Dicke gravity. Class. Quant. Grav. 1998, 15, 3271. Carloni, S.; Dunsby, P.K.S.; Capozziello, S.; Troisi, A. Cosmological dynamics of Rn gravity. Class. Quant. Grav. 2005, 22, 4839. Kuusk, P.; Jarv, L.; Saal, M. Scalar-tensor cosmologies: General relativity as a fixed point of the Jordan frame scalar field. Int. J. Mod. Phys. A 2009, 24, 1631. Jarv, L.; Kuusk, P.; Saal, M. Potential dominated scalar-tensor cosmologies in the general relativity limit: Phase space view. Phys. Rev. D 2010, 81, 104007. Faraoni, V. Phase space geometry in scalar-tensor cosmology. Ann. Phys. 2005, 317, 366. Tsujikawa, S. Modified gravity models of dark energy. Lect. Notes Phys. 2010, 800, 99. Amendola, L. Scaling solutions in general nonminimal coupling theories. Phys. Rev. D 1999, 60, 043501. Barrow, J.D.; Maeda, K. Extended inflationary universes. Nucl. Phys. 1990, B341, 294. Uzan, J.P. Cosmological scaling solutions of nonminimally coupled scalar fields. Phys. Rev. D 1999, 59, 123510. Liddle, A.R.; Mazumdar, A; Schunck, F.E. Assisted inflation. Phys. Rev. D 1998, 58, 061301. Malik, K.A.; Wands, D. Dynamics of assisted inflation. Phys. Rev. D 1999, 59, 123501. Copeland, E.J.; Mazumdar, A.; Nunes, N.J. Generalized assisted inflation. Phys. Rev. D 1999, 60, 083506. Kim, S.A.; Liddle, A.R.; Tsujikawa, S. Dynamics of assisted quintessence. Phys. Rev. D 2005, 72, 043506. Symmetry 2016, 8, 70 70 of 79 120. Tsujikawa, S. General analytic formulae for attractor solutions of scalar-field dark energy models and their multi-field generalizations. Phys. Rev. D 2006, 73, 103504. 121. Ohashi, J.; Tsujikawa, S. Assisted dark energy. Phys. Rev. D 2009, 80, 103513. 122. Karwan, K. Dynamics of entropy perturbations in assisted dark energy with mixed kinetic terms. J. Cosmol. Astropart. Phys. 2011, 1102, 007. 123. Zlatev, I.; Wang, L.M.; Steinhardt, P.J. Quintessence, cosmic coincidence, and the cosmological constant. Phys. Rev. Lett. 1999, 82, 896. 124. Amendola, L. Perturbations in a coupled scalar field cosmology. Mon. Not. Roy. Astron. Soc. 2000, 312, 521. 125. Holden, D.J.; Wands, D. Selfsimilar cosmological solutions with a nonminimally coupled scalar field. Phys. Rev. D 2000, 61, 043506. 126. Amendola, L. Coupled quintessence. Phys. Rev. D 2000, 62, 043511. 127. Brookfield, A.W.; van de Bruck, C.; Hall, L.M.H. New interactions in the dark sector mediated by dark energy. Phys. Rev. D 2008, 77, 043006. 128. Baldi, M. Multiple Dark Matter as a self-regulating mechanism for dark sector interactions. Ann. Phys. 2012, 524, 602. 129. Amendola, L.; Barreiro, T.; Nunes, N.J. Multifield coupled quintessence. Phys. Rev. D 2014, 90 , 083508. 130. Piloyan, A.; Marra, V.; Baldi, M.; Amendola, L. Supernova constraints on Multi-coupled Dark Energy. J. Cosmol. Astropart. Phys. 2013, 1307, 42. 131. Piloyan, A.; Marra, V.; Baldi, M.; Amendola, L. Linear Perturbation constraints on Multi-coupled Dark Energy. J. Cosmol. Astropart. Phys. 2014, 1402, 45. 132. Baldi, M. Cold dark matter halos in Multi-coupled Dark Energy cosmologies: Structural and statistical properties. Phys. Dark Univ. 2014, 3, 4. 133. Armendariz-Picon, C.; Mukhanov, V.F.; Steinhardt, P.J. Essentials of k-essence. Phys. Rev. D 2001, 63, 103510. 134. Deffayet, C.; Pujolas, O.; Sawicki, I.; Vikman, A. Imperfect Dark Energy from Kinetic Gravity Braiding. J. Cosmol. Astropart. Phys. 2010, 1010, 026. 135. Woodard, R.P. Avoiding dark energy with 1/r modifications of gravity. Lect. Notes Phys. 2007, 720, 403. 136. Deffayet, C.; Gao, X.; Steer, D.A.; Zahariade, G. From k-essence to generalised Galileons. Phys. Rev. D 2011, 84, 064039. 137. Deffayet, C.; Esposito-Farese, G.; Vikman, A. Covariant Galileon. Phys. Rev. D 2009, 79, 084003. 138. Nicolis, A.; Rattazzi, R.; Trincherini, E. The Galileon as a local modification of gravity. Phys. Rev. D 2009, 79, 064036. 139. Zumalacárregui, M.; García-Bellido, J. Transforming gravity: From derivative couplings to matter to second-order scalar-tensor theories beyond the Horndeski Lagrangian. Phys. Rev. D 2014, 89, 064046. 140. Gleyzes, J.; Langlois, D.; Piazza, F.; Vernizzi, F. Healthy theories beyond Horndeski. Phys. Rev. Lett. 2015, 114, 211101. 141. Gleyzes, J.; Langlois, D.; Piazza, F.; Vernizzi, F. Exploring gravitational theories beyond Horndeski. J. Cosmol. Astropart. Phys. 2015, 1502, 18. 142. Weinberg, S. The Cosmological Constant Problem. Rev. Mod. Phys. 1989, 61, 1. 143. Carroll, S.M. The Cosmological constant. Living Rev. Rel. 2001, 4, 1. 144. Kaloper, M.; Padilla, A. Vacuum Energy Sequestering: The Framework and Its Cosmological Consequences. Phys. Rev. D 2014, 90, 084023. 145. Charmousis, C.; Copeland, E.J.; Padilla, A.; Saffin, P.M. General second order scalar-tensor theory, self tuning, and the Fab Four. Phys. Rev. Lett. 2012, 108, 051101. 146. Martín-Moruno, P.; Nunes, N.J.; Lobo, F.S.N. Attracted to de Sitter: Cosmology of the linear Horndeski models. J. Cosmol. Astropart. Phys. 2015, 1505, 33. 147. Martín-Moruno, P.; Nunes, N.J. Attracted to de Sitter II: Cosmology of the shift-symmetric Horndeski models. J. Cosmol. Astropart. Phys. 2015, 1509, 56. 148. Penrose, R. Gravitational collapse and space-time singularities. Phys. Rev. Lett 1965, 14, 57. 149. Penrose, R. Gravitational collapse: The role of general relativity. Riv. Nuovo Cim. Numero Speciale 1969, 1, 252. 150. Hawking, S.W. Singularities in the universe. Phys. Rev. Lett 1966, 17, 444. 151. Jungman, G.; Kamionkowski, M.; Griest, K. Supersymmetric dark matter. Phys. Rept. 1996, 267, 195. 152. Green, M.; Schwarz, J.; Witten, E. Superstring Theory; Cambridge University Press: Cambridge, UK, 1987. Symmetry 2016, 8, 70 71 of 79 153. Thiemann, T. Lectures on loop quantum gravity. Lect. Notes Phys. 2003, 631, 41. 154. Birrell, N.D.; Davies, P.C.W. Quantum Fields in Curved Space; Cambridge University Press: Cambridge, UK, 1982. 155. Buchdahl, H.A. Non-linear Lagrangians and cosmological theory. Mon. Not. Roy. Astron. Soc. 1970, 150, 1. 156. Barrow, J.D.; Ottewill, A.C. The Stability of General Relativistic Cosmological Theory. J. Phys. A 1983, 16, 2757. 157. Ferraris, M.; Francaviglia, M.; Volovich, I. The Universality of vacuum Einstein equations with cosmological constant. Class. Quant. Grav. 1994, 11, 1505. 158. Vollick, D.N. 1/R Curvature corrections as the source of the cosmological acceleration. Phys. Rev. D 2003, 68, 063510. 159. Meng, X.H.; Wang, P. Palatini formation of modified gravity with ln R terms. Phys. Lett. B 2004, 584, 1. 160. Poplawski, N.J. Interacting dark energy in f(R) gravity. Phys. Rev. D 2006, 74, 084032. 161. Li, B.; Chan, K.C.; Chu, M.C. Constraints on f(R) Cosmology in the Palatini Formalism. Phys. Rev. D 2007, 76, 024002. 162. Li, B.; Barrow, J.D.; Mota, D.F. The Cosmology of Ricci-Tensor-Squared gravity in the Palatini variational approach. Phys. Rev. D 2007, 76, 104047. 163. Iglesias, A.; Kaloper, N.; Padilla, A.; Park, M. How (Not) to Palatini. Phys. Rev. D 2007, 76, 104001. 164. Amendola, L.; Appleby, S.; Bacon, D.; Baker, T.; Baldi, M.; Bartolo, N.; Blanchard, A.; Bonvin, C.; Borgani, S.; Branchini, E.; et al. Cosmology and fundamental physics with the Euclid satellite. Living Rev. Rel. 2013, 16, 6. 165. Capozziello, S.; Carloni, S.; Troisi, A. Quintessence without scalar fields. Recent Res. Dev. Astron. Astrophys. 2003, 1, 625. 166. Sotiriou, T.P. Modified Actions for Gravity: Theory and Phenomenology. 2007, arXiv:gr-qc/0710.4438. 167. Vikman, A. Can dark energy evolve to the phantom? Phys. Rev. D 2005, 71, 023515. 168. Cai, Y.F.; Saridakis, E.N.; Setare, M.R.; Xia, J.Q. Quintom Cosmology: Theoretical implications and observations. Phys. Rept. 2010, 493, 1. 169. Bertolami, O.; Boehmer, C.G.; Harko, T.; Lobo, F.S.N. Extra force in f(R) modified theories of gravity. Phys. Rev. D 2007, 75, 104016. 170. Bertolami, O.; Páramos, J. Do f(R) theories matter? Phys. Rev. D 2008, 77, 084018. 171. Bertolami, O.; Páramos, J. On the non-trivial gravitational coupling to matter. Class. Quant. Grav. 2008, 25, 245017. 172. Faraoni, V. A viability criterion for modified gravity with an extra force. Phys. Rev. D 2007, 76, 127501. 173. Nojiri, S.; Odintsov, S.D. Gravity assisted dark energy dominance and cosmic acceleration. Phys. Lett. B 2004, 599, 137. 174. Nojiri, S.; Odintsov, S.D. Dark energy and cosmic speed-up from consistent modified gravity. PoS WC 2004, 2004, 24. 175. Allemandi, G.; Borowiec, A.; Francaviglia, M.; Odintsov, S.D. Dark energy dominance and cosmic acceleration in first order formalism. Phys. Rev. D 2005, 72, 063505. 176. Mukohyama, S.; Randall, L. A Dynamical approach to the cosmological constant. Phys. Rev. Lett. 2004, 92, 211302. 177. Harko, T.; Koivisto, T.S.; Lobo, F.S.N. Palatini formulation of modified gravity with a nonminimal curvature-matter coupling. Mod. Phys. Lett. A 2011, 26, 1467. 178. Koivisto, T. Covariant conservation of energy momentum in modified gravities. Class. Quant. Grav. 2006, 23, 4289. 179. Bertolami, O.; Paramos, J.; Harko, T.; Lobo, F.S.N. Non-minimal curvature-matter couplings in modified gravity. 2008, arXiv:0811.2876 [gr-qc]. 180. Sotiriou, T.P.; Faraoni, V. Modified gravity with R-matter couplings and (non-)geodesic motion. Class. Quant. Grav. 2008, 25, 205002. 181. Schutz, B.F. Perfect Fluids in General Relativity: Velocity Potentials and a Variational Principle. Phys. Rev. D 1970, 2, 2762. 182. Brown, J.D. Action functionals for relativistic perfect fluids. Class. Quant. Grav. 1993, 10, 1579. 183. Bertolami, O.; Lobo, F.S.N.; Páramos, J. Non-minimum coupling of perfect fluids to curvature. Phys. Rev. D 2008, 78, 064036. Symmetry 2016, 8, 70 72 of 79 184. Hawking, S.W.; Ellis, G.F.R. The Large Scale Structure of Spacetime; Cambridge University Press: Cambridge, UK, 1973. 185. Falk, J. Theory of elasticity of coherent inclusions by means of nonmetric geometry. J. Elast. 1981, 11, 359. 186. Kröner, E. The continuized crystal—A bridge between micro and macromechanics? Z. Angew. Math. Mech. 1986, 66, T284–T292. 187. Clayton, J.D. Nonlinear Mechanics of Crystals; Springer: Maryland, MD, USA, 2011. 188. Lobo, F.S.N.; Olmo, G.J.; Rubiera-Garcia, D. Crystal clear lessons on the microstructure of spacetime and modified gravity. Phys. Rev. D 2015, 91, 124001. 189. Olmo, G.J.; Rubiera-Garcia, D. The quantum, the geon, and the crystal. Int. J. Mod. Phys. D 2015, 24, 1542013. 190. Olmo, G.J.; Rubiera-Garcia, D. Importance of torsion and invariant volumes in Palatini theories of gravity. Phys. Rev. D 2013, 88, 084030. 191. Olmo, G.J.; Rubiera-Garcia, D. Palatini f ( R) black holes in nonlinear electrodynamics. Phys. Rev. D 2011, 84, 124059. 192. Exirifard, Q.; Sheikh-Jabbari, M.M. Lovelock gravity at the crossroads of Palatini and metric formulations. Phys. Lett. B 2008, 661, 158. 193. Borunda, M.; Janssen B.; Bastero-Gil, M. Palatini versus metric formulation in higher curvature gravity. J. Cosmol. Astropart. Phys. 2008, 2008, 0811. 194. Olmo, G.J.; Rubiera-Garcia, D. Brane-world and loop cosmology from a gravity-matter coupling perspective. Phys. Lett. B 2015, 740, 73. 195. Bojowald, M. Loop quantum cosmology. Living Rev. Rel. 2008, 11, 4. 196. Ashtekar, A.; Singh, P. Loop Quantum Cosmology: A Status Report. Class. Quant. Grav 2011, 28, 213001. 197. Maartens, R.; Koyama, K. Brane-World Gravity. Living Rev. Rel 2010, 13, 5. 198. Ashtekar, A.; Pawlowski, T.; Singh, P. Quantum nature of the big bang. Phys. Rev. Lett. 2006, 96, 141301. 199. Barragan, C.; Olmo, G.J.; Sanchis-Alepuz, H. Bouncing Cosmologies in Palatini f ( R) Gravity. Phys. Rev. D 2009, 80, 024016. 200. Olmo, G.J.; Rubiera-Garcia, D. Nonsingular black holes in f ( R) theories. Universe 2015, 1, 173. 201. Bambi, C.; Cardenas-Avendano, A.; Olmo, G.J.; Rubiera-Garcia, D. Wormholes and nonsingular space-times in Palatini f ( R) gravity. Phys. Rev. D 2016, 93, 064016. 202. Barragan, C.; Olmo, G.J. Isotropic and Anisotropic Bouncing Cosmologies in Palatini Gravity. Phys. Rev. D 2010, 82, 084015. 203. Deser, S.; Gibbons, G.W. Born-Infeld-Einstein actions? Class. Quant. Grav 1998, 15, L35. 204. Bañados, M.; Ferreira, P.G. Eddington’s theory of gravity and its progeny. Phys. Rev. Lett 2010, 105, 011101. 205. Jiménez, J.B.; Heisenberg, L.; Olmo, G.J. Infrared lessons for ultraviolet gravity: The case of massive gravity and Born-Infeld. J. Cosmol. Astropart. Phys. 2014, 1411, 004. 206. Olmo, G.J.; Rubiera-Garcia, D.; Sanchez-Puente, A. Geodesic completeness in a wormhole spacetime with horizons. Phys. Rev. D. 2015, 92, 044047. 207. Bazeia, D.; Losano, L.; Olmo, G.J.; Rubiera-Garcia, D.; Sanchez-Puente, A. Classical resolution of black hole singularities in arbitrary dimension. Phys. Rev. D 2015, 92, 044018. 208. Olmo, G.J.; Rubiera-Garcia, D.; Sanchez-Puente, A. Classical resolution of black hole singularities via wormholes. Eur. Phys. J. C 2016, 76, 143. 209. Odintsov, O.D.; Olmo, G.J.; Rubiera-Garcia, D. Born-Infeld gravity and its functional extensions. 2014, 90, 044003. 210. Hawking, S. Gravitationally collapsed objects of very low mass. Mon. Not. Roy. Astron. Soc. 1971, 152, 75. 211. Lobo, F.S.N.; Olmo, O.J.; Rubiera-Garcia, D. Semiclassical geons as solitonic black hole remnants. J. Cosmol. Astropart. Phys. 2013, 1307, 011. 212. Chen, P.; Ong, Y.C.; Yeom, D.H. Black Hole Remnants and the Information Loss Paradox. 2014, arXiv:1412.8366 [gr-qc]. 213. Capozziello, S.; Harko, T.; Koivisto, T.S.; Lobo, F.S.N.; Olmo, G.J. Cosmology of hybrid metric-Palatini f(X)-gravity. J. Cosmol. Astropart. Phys. 2013, 1304, 11. 214. Capozziello, S.; Harko, T.; Koivisto, T.S.; Lobo, F.S.N.; Olmo, G.J. Wormholes supported by hybrid metric-Palatini gravity. Phys. Rev. D 2012, 86, 127504. Symmetry 2016, 8, 70 73 of 79 215. Capozziello, S.; Harko, T.; Koivisto, T.S.; Lobo, F.S.N.; Olmo, G.J. The virial theorem and the dark matter problem in hybrid metric-Palatini gravity. J. Cosmol. Astropart. Phys. 2013, 1307, 24. 216. Capozziello, S.; Harko, T.; Koivisto, T.S.; Lobo, F.S.N.; Olmo, G.J. Hybrid metric-Palatini gravity. Universe 2015, 1, 199. 217. Mukhanov, V. Physical Foundations of Cosmology; Cambridge University Press: Cambridge, UK, 2005. 218. Liddle A.R.; Lyth, D.H. Cosmological Inflation and Large-Scale Structure; Cambridge University Press: Cambridge, UK, 2000. 219. Dodelson, S. Modern Cosmology; Academic Press: San Diego, CA, USA, 1999. 220. Ade, P.A.R.; Aghanim, A.N.; Arnaud, M.; Arroja, F.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Ballardini, M.; Banday, A.J.; Barreiro, R.B.; et al. Planck 2015 results. XX. Constraints on inflation. 2015, arXiv:1502.02114 [astro-ph.CO]. 221. Lidsey, J.E.; Liddle, A.R.; Kolb, E.W.; Copeland, E.J.; Barreiro, T.; Abney, M. Reconstructing the inflation potential: An overview. Rev. Mod. Phys. 1997, 69, 373. 222. Elizalde, E.; Nojiri, S.; Odintsov, S.D.; Saez-Gomez, D.; Faraoni, V. Reconstructing the universe history, from inflation to acceleration, with phantom and canonical scalar fields. Phys. Rev. D 2008, 77, 106005. 223. Bamba, K.; Nojiri, S.; Odintsov, S.D. Reconstruction of scalar field theories realizing inflation consistent with the Planck and BICEP2 results. Phys. Lett. B 2014, 737 374. 224. Capozziello, S.; Faraoni, V. Beyond Einstein Gravity; Springer: Dordrecht, The Netherlands, 2010. 225. De la Cruz-Dombriz, A.; Saez-Gomez, D. Black holes, cosmological solutions, future singularities, and their thermodynamical properties in modified gravity theories. Entropy 2012, 14, 1717. 226. Saez-Gomez, D. Modified f ( R) gravity from scalar-tensor theory and inhomogeneous EoS dark energy. Gen. Rel. Grav. 2009, 41, 1527. 227. Nojiri, S.; Odintsov, S.D.; Saez-Gomez, D. Cosmological reconstruction of realistic modified F(R) gravities. Phys. Lett. B 2009, 681, 74. 228. Bamba, K.; Nojiri, S.; Odintsov, S.D.; Sáez-Gómez, D. Inflationary universe from perfect fluid and F ( R) gravity and its comparison with observational data. Phys. Rev. D 2014, 90, 124061. 229. De Felice, A.; Tsujikawa, S.; Elliston, G.; Tavakol, R. Chaotic inflation in modified gravitational theories. J. Cosmol. Astropart. Phys. 2011, 1108, 021. 230. Starobinsky, A.A. A New Type of Isotropic Cosmological Models Without Singularity. Phys. Lett. B 1980, 91, 99. 231. Copeland, E.J.; Rahmede, C.; Saltas, I.D. Asymptotically Safe Starobinsky Inflation. Phys. Rev. D 2015, 91, 103530. 232. De la Cruz-Dombriz, A.; Elizalde, E.; Odintsov, S.D.; Saez-Gomez, D. Spotting deviations from R2 inflation. 2016, arXiv:1603.05537 [gr-qc]. 233. Sasaki, M.; Stewart, E.D. A General analytic formula for the spectral index of the density perturbations produced during inflation. Prog. Theor. Phys. 1996, 95, 71. 234. Avelino, P.P.; Menezes, R.; Sousa, L. p-brane dynamics in N+1-dimensional FRW universes. Phys. Rev. D 2009, 79, 043519. 235. Sousa, L.; Avelino, P.P. p-brane dynamics in (N+1)-dimensional FRW universes: A unified framework. Phys. Rev. D 2011, 83, 103507. 236. Sousa, L.; Avelino, P.P. The cosmological evolution of p-brane networks. Phys. Rev. D 2011, 84, 063502. 237. Martins, C.J.A.P.; Shellard, E.P.S. Quantitative string evolution. Phys. Rev. D 1996, 54, 2535. 238. Martins, C.J.A.P.; Shellard, E.P.S. Extending the velocity dependent one scale string evolution model. Phys. Rev. D 2002, 65, 043514. 239. Leite, A.M.M.; Martins, C.J.A.P.; Shellard, E.P.S. Accurate Calibration of the Velocity-Dependent One-Scale Model for Domain Walls. Phys. Lett. B 2013, 718, 740. 240. Oliveira, M.F.; Avgoustidis, A.; Martins, C.J.A.P. Cosmic string evolution with a conserved charge. Phys. Rev. D 2012, 85, 083515. 241. Avgoustidis, A.; Shellard, E.P.S. Velocity-Dependent Models for Non-Abelian/Entangled String Networks. Phys. Rev. D 2008, 78, 103510. 242. Nunes, A.S.; Avgoustidis, A.; Martins, C.J.A.P.; Urrestilla, J. Analytic Models for the Evolution of Semilocal String Networks. Phys. Rev. D 2011, 84, 063504. 243. Martins, C.J.A.P. Evolution of Hybrid Defect Networks. Phys. Rev. D 2009, 80, 083527. Symmetry 2016, 8, 70 74 of 79 244. Avelino; P.P.; Martins, C.J.A.P. Primordial adiabatic fluctuations from cosmic defects. Phys. Rev. Lett. 2000, 85, 1370. 245. Wu, J.H.P.; Avelino, P.P.; Shellard, E.P.S.; Allen, B. Cosmic strings, loops, and linear growth of matter perturbations. Int. J. Mod. Phys. D 2002, 11, 61. 246. Ade, P.A.R.; Aghanim, N.; Armitage-Caplan, C.; Arnaud, M.; Ashdown, M.; Atrio-Barandela, F.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; et al. [Planck Collaboration], Planck 2013 results. XXV. Searches for cosmic strings and other topological defects. Astron. Astrophys. 2014, 571, A25. 247. Bardeen, J.M. Gauge Invariant Cosmological Perturbations. Phys. Rev. D 1980, 22, 1882. 248. Martin, J.; Ringeval, C.; Vennin, V. How Well Can Future CMB Missions Constrain Cosmic Inflation? J. Cosmol. Astropart. Phys. 2014, 1410, 038. 249. Sousa, L.; Avelino, P.P. Cosmic Microwave Background anisotropies generated by domain wall networks. Phys. Rev. D 2015, 92, 083520. 250. Pogosian, P.; Vachaspati, T. Cosmic microwave background anisotropy from wiggly strings. Phys. Rev. D 1999, 60, 083504. 251. Pogosian, L.; Wasserman, I.; Wyman, M. On vector mode contribution to CMB temperature and polarization from local strings. 2006, astro-ph/0604141. 252. Vilenkin, A. Gravitational radiation from cosmic strings. Phys. Lett. B 1981, 107, 47. 253. Brandenberger, R.H.; Albrecht, A.; Turok, N. Gravitational Radiation from Cosmic Strings and the Microwave Background. Nucl. Phys. B 1986, 277, 605. 254. Sanidas, S.A.; Battye, R.A.; Stappers, B.W. Constraints on cosmic string tension imposed by the limit on the stochastic gravitational wave background from the European Pulsar Timing Array. Phys. Rev. D 2012, 85, 122003. 255. Sousa, L.; Avelino, P.P. Stochastic Gravitational Wave Background generated by Cosmic String Networks: Velocity-Dependent One-Scale model versus Scale-Invariant Evolution. Phys. Rev. D 2013, 88, 023516. 256. Aasi, J.; Abbott, B.P.; Abbott, R.; Abbott, T.; Abernathy, M.R.; Ackley, K.; Adams, C.; Adams, T.; Addesso, P.; Adhikari, R.X.; et al. Advanced LIGO. Class. Quant. Grav. 2015, 32, 074001. 257. Accadia, T.; Agathos, M.; Allocca, A.; Astone, P.; Ballardin, G.; Barone, F.; Barsuglia, M.; Basti, A.; Bauer, T.S.; Bejger, M.; et al. Advanced Virgo Interferometer: A Second Generation Detector for Gravitational Waves Observation. 2015, arXiv:1408.3978 [gr-qc]. 258. Dufaux, J.F. Cosmological Backgrounds of Gravitational Waves and eLISA. ASP Conf. Ser. 2013, 467, 91. 259. Manchester, R.N. The Parkes Pulsar Timing Array. AIP Conf. Proc. 2008, 983, 584. 260. Ferdman, R.D.; van Haasteren, R.; Bassa, C.G.; Burgay, M.; Cognard, I.; Corongiu, A.; D’Amico, N.; Desvignes, G.; Hessels, J.W.T.; Janssen, G.H.; et al. The European Pulsar Timing Array: Current efforts and a LEAP toward the future. Class. Quant. Grav. 2010, 27, 084014. 261. Hiramatsu, T.; Kawasaki, M.; Saikawa, K. On the estimation of gravitational wave spectrum from cosmic domain walls. J. Cosmol. Astropart. Phys. 2014, 1402, 031. 262. Hiramatsu, T.; Kawasaki, M.; Saikawa, K. Gravitational Waves from Collapsing Domain Walls. J. Cosmol. Astropart. Phys. 2010, 1005, 032. 263. Bleem, L.E.; Stalder, B.; de Haan, T.; Aird, K.A.; Allen, S.W.; Applegate, D.E.; Ashby, L.N.; Bautz, M.; Bayliss, M.; Benson, B.A.; et al. Galaxy Clusters discovered via the Sunyaev-Zel’dovic effect in the 2500-square-degree SPT-SZ Survey. Astrophys. J. Suppl. 2015, 216, 27. 264. Ade, P.A.R.; Aghanim, N.; Arnaud, M.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; Bartlett, J.G.; Bartolo, N.; et al. Planck Collaboration, Planck 2015 results. XXIV. Cosmology from Sunyaev-Zeldovich cluster counts. 2015, arXiv:1502.01597 265. André, P.; Baccigalupi, C.; Banday, A.; Barbosa, D.; Barreiro, B.; Bartlett, J.; Bartolo, N.; Battistelli, E.; Battye, R.; Bendo, G.; et al. PRISM (Polarized Radiation Imaging and Spectroscopy Mission): An extended white paper. J. Cosmol. Astropart. Phys. 2014, 2, 006. 266. Ramos, E.P.R.G.; da Silva, A.J.C.; Liu, G.-C. Cosmic Microwave Background Induced Polarization from Single Scattering by Clusters of Galaxies and Filaments. Astrophys. J. 2012, 757, 44. 267. Jenkins, A.; Frenk, C.S.; White, S.D.M.; Colberg, J.M.; Cole, S.; Evrard, A.E.; Couchman, H.M.P.; Yoshida, N. The mass function of dark matter haloes. Mon. Not. Roy. Astron. Soc. 2001, 321, 372. 268. Press, W.H.; Schechter, P. Formation of Galaxies and Clusters of Galaxies by Self-Similar Gravitational Condensation. Astrophys. J. 1974, 187, 425. Symmetry 2016, 8, 70 75 of 79 269. Kravtsov, A.V.; Borgani, S. Formation of Galaxy Clusters. Ann. Rev. Astron. Astrophys. 2012, 50, 353. 270. Wagner, C.; Verde, L. N-body simulations with generic non-Gaussian initial conditions II: Halo bias. J. Cosmol. Astropart. Phys. 2012, 3, 002. 271. Achitouv, I.; Baldi, M.; Puchwein, E.; Weller, J. The Imprint of f ( R) Gravity on Non-Linear Structure Formation. 2015, arXiv:1511.01494 272. Nunes, N.J.; da Silva, A.C.; Aghanim, N. Number counts in homogeneous and inhomogeneous dark energy models. Astron. Astrophys. 2006, 450, 899. 273. Kaiser, N. Evolution and clustering of rich clusters. Mon. Not. Roy. Astron. Soc. 1986, 222, 323. 274. Kay, S.T.; da Silva, A.C.; Aghanim, N.; Blanchard, A.; Liddle, A.R.; Puget, J.-L.; Sadat, R.; Thomas, P.A. The evolution of clusters in the CLEF cosmological simulation: X-ray structural and scaling properties. Mon. Not. Roy. Astron. Soc. 2007, 377, 317. 275. Da Silva, A. C.; Catalano, A.; Montier, L.; Pointecouteau, E.; Lanoux, J.; Giard, M. The impact of dust on the scaling properties of galaxy clusters. Mon. Not. Roy. Astron. Soc. 2009, 396, 849. 276. Aghanim, N.; da Silva, A.C.; Nunes, N.J. Cluster scaling relations from cosmological hydrodynamic simulations in a dark-energy dominated universe. Astron. Astrophys. 2009, 496, 637. 277. Trindade, A.M.M.; da Silva, A. Effect of Priomordial non-Gaussianities on Galaxy Clusters Scaling Relations. 2016, arXiv:1603.09270 278. Navarro, J.F.; Frenk, C.S.; White, S.D.M. A Universal Density Profile from Hierarchical Clustering. Astrophys. J. 1997, 490, 493. 279. Nagai, D.; Kravtsov, A.V.; Vikhlinin, A. Effects of Galaxy Formation on Thermodynamics of the Intracluster Medium. Astrophys. J. 2007, 668, 1. 280. Da Silva, A.C.; Barbosa, D.; Liddle, A.R.; Thomas, P.A. Hydrodynamical simulations of the Sunyaev-Zel’dovich effect: The kinetic effect. Mon. Not. Roy. Astron. Soc. 2001, 326, 155. 281. Ade, P.A.R.; Aghanim, N.; Arnaud, M.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; Barrena, R.; Bartlett, J.G.; et al. Planck Collaboration, Planck 2015 results. XXVII. The Second Planck Catalogue of Sunyaev-Zeldovich Sources. 2015, arXiv:1502.01598. 282. Ade, P.A.R.; Aghanim, N.; Armitage-Caplan, C.; Arnaud, M.; Ashdown, M.; Atrio-Barandela, F.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; et al. Planck Collaboration, Planck 2013 results. XX. Cosmology from Sunyaev-Zeldovich cluster counts. Astron. Astrophys. 2014, 571, A20. 283. Ade, P.A.R.; Aghanim, N.; Armitage-Caplan, C.; Arnaud, M.; Ashdown, M.; Atrio-Barandela, F.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; et al. Planck Collaboration, Planck 2013 results. XXI. Power spectrum and high-order statistics of the Planck all-sky Compton parameter map. Astron. Astrophys. 2014, 571, A21. 284. Aghanim, N.; Arnaud, M.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; Bartlett, J.G.; Bartolo, N.; Battaner, E.; et al. Planck Collaboration, Planck 2015 results. XXII. A map of the thermal Sunyaev-Zeldovich effect. 2015, arXiv:1502.01596. 285. De Martino, I.; Atrio-Barandela, F.; da Silva, A.; Ebeling, H.; Kashlinsky, A.; Kocevski, D.; Martins, C.J.A.P. Measuring the Redshift Dependence of the Cosmic Microwave Background Monopole Temperature with Planck Data. Astrophys. J. 2012, 757, 144. 286. De Martino, I.; Génova-Santos, R.; Atrio-Barandela, F.; Ebeling, H.; Kashlinsky, A.; Kocevski, D.; Martins, C.J.A.P. Constraining the Redshift Evolution of the Cosmic Microwave Background Blackbody Temperature with PLANCK Data. Astrophys. J. 2015, 808, 128. 287. Rudjord, Ø.; Hansen, F.K.; Lan, X.; Liguori, M.; Marinucci, D.; Matarrese, S. An Estimate of the Primordial Non-Gaussianity Parameter f N L Using the Needlet Bispectrum from WMAP. Astrophys. J. 2009, 701, 369. 288. Komatsu, E.; Dunkley, J.; Nolta, M.R.; Bennett, C.L.; Gold, B.; Hinshaw, G.; Jarosik, N.; Larson, D.; Limon, M.; Page, L.; et al. Five-Year Wilkinson Microwave Anisotropy Probe Observations: Cosmological Interpretation. Astrophys. J. S. 2009, 180, 330. 289. Curto, A.; Martínez-González, E.; Barreiro, R.B. Improved Constraints on Primordial Non-Gaussianity for the Wilkinson Microwave Anisotropy Probe 5-Year Data, Astrophys. J. 2009, 706, 399. 290. Ade, P.A.R.; Aghanim, N.; Armitage-Caplan, C.; Arnaud, M.; Ashdown, M.; Atrio-Barandela, F.; Aumont, J.; Baccigalupi, C.; Banday, A.J.; Barreiro, R.B.; et al. Planck Collaboration, Planck 2013 results. XXIV. Constraints on primordial non-Gaussianity. Astron. Astrophys. 2014, 571, A24. Symmetry 2016, 8, 70 76 of 79 291. Ade, P.A.R.; Aghanim, N.; Arnaud, M.; Arroja, F.; Ashdown, M.; Aumont, J.; Baccigalupi, C.; Ballardini, M.; Banday, A.J.; Barreiro, R.B.; et al. Planck Collaboration, Planck 2015 results. XVII. Constraints on primordial non-Gaussianity. 2015, arXiv:1502.01592 [astro-ph.CO]. 292. Sefusatti, E.; Komatsu, E. Bispectrum of galaxies from high-redshift galaxy surveys: Primordial non-Gaussianity and nonlinear galaxy bias, Phys. Rev. D 2007, 76, 083004. 293. Matarrese, S.; Verde, L. The Effect of Primordial Non-Gaussianity on Halo Bias. Astrophys. J. Lett. 2008, 677, L77. 294. Giannantonio, T.; Porciani, C.; Carron, J.; Amara, A.; Pillepich, A. Constraining primordial non-Gaussianity with future galaxy surveys. Mon. Not. Roy. Astron. Soc. 2012, 422, 2854. 295. Schäfer, B.M.; Grassi, A.; Gerstenlauer, M.; Byrnes, C.T. A weak lensing view on primordial non-Gaussianities. Mon. Not. Roy. Astron. Soc. 2012, 421, 797. 296. Hilbert, S.; Marian, L.; Smith, R.E.; Desjacques, V. Measuring primordial non-Gaussianity with weak lensing surveys. Mon. Not. Roy. Astron. Soc. 2012, 426, 2870. 297. Tashiro, H.; Ho, S. Constraining primordial non-Gaussianity with CMB-21 cm cross-correlations? Mon. Not. Roy. Astron. Soc. 2013, 431, 2017. 298. Takeuchi, Y.; Ichiki, K.; Matsubara, T. Application of cross correlations between CMB and large-scale structure to constraints on the primordial non-Gaussianity. Phys. Rev. D 2012, 85, 043518. 299. Robinson, J.; Baker, J.E. Evolution of the cluster abundance in non-Gaussian models. Mon. Not. Roy. Astron. Soc. 2000, 311, 781. 300. Matarrese, S.; Verde, L.; Jimenez, R. The Abundance of High-Redshift Objects as a Probe of Non-Gaussian Initial Conditions. Astrophys. J. 2000, 541, 10. 301. Kamionkowski, M.; Verde, L.; Jimenez, R. The void abundance with non-gaussian primordial perturbations. J. Cosmol. Astropart. Phys. 2009, 1, 010. 302. Lam, T.Y.; Sheth, R.K.; Desjacques, V. The initial shear field in models with primordial local non-Gaussianity and implications for halo and void abundances. Mon. Not. Roy. Astron. Soc. 2009, 399, 1482. 303. D’Amico, G.; Musso, M.; Noreña, J.; Paranjape, A. Excursion sets and non-Gaussian void statistics. Phys. Rev. D 2011, 83, 023521. 304. Sekiguchi, T.; Yokoyama, S. Void bias from primordial non-Gaussianities. 2012, arXiv:1204.2726 [astro-ph.CO]. 305. Fergusson, J.R.; Shellard, E.P.S. Shape of primordial non-Gaussianity and the CMB bispectrum. Phys. Rev. D 2009, 80, 043510. 306. LoVerde, M.; Miller, A.; Shandera, S.; Verde, L. Effects of scale-dependent non-Gaussianity on cosmological structures. J. Cosmol. Astropart. Phys. 2008, 4, 014. 307. Giannantonio, T.; Porciani, C. Structure formation from non-Gaussian initial conditions: Multivariate biasing, statistics, and comparison with N-body simulations. Phys. Rev. D 2010, 81, 063530. 308. Maggiore, M.; Riotto, A. The Halo Mass Function from Excursion Set Theory. III. Non-Gaussian Fluctuations. Astrophys. J. 2010, 717, 526. 309. D’Amico, G.; Musso, M.; Noreña, J.; Paranjape, A. An improved calculation of the non-Gaussian halo mass function. J. Cosmol. Astropart. Phys. 2011, 2, 001. 310. LoVerde, M.; Smith, K.M. The non-Gaussian halo mass function with f N L, g N L and τNL . J. Cosmol. Astropart. Phys. 2011, 8, 003. 311. Achitouv, I.E.; Corasaniti, P.S. Non-Gaussian halo mass function and non-spherical halo collapse: Theory vs. simulations. J. Cosmol. Astropart. Phys. 2012, 2, 2. 312. D’Aloisio, A.; Zhang, J.; Jeong, D.; Shapiro, P.R. Halo statistics in non-Gaussian cosmologies: The collapsed fraction, conditional mass function and halo bias from the path-integral excursion set method. Mon. Not. Roy. Astron. Soc. 2013, 428, 2765. 313. Zentner, A.R. The Excursion Set Theory of Halo Mass Functions, Halo Clustering, and Halo Growth. Intern. J. Mod. Phys. D 2007, 16, 763. 314. Maggiore, M.; Riotto, A. The Halo Mass Function from Excursion Set Theory. I. Gaussian Fluctuations with Non-Markovian Dependence on the Smoothing Scale. Astrophys. J. 2010, 711, 907. 315. Trindade, A.M.M.; Avelino, P.P.; Viana, P.T.P. A new signature of primordial non-Gaussianities from the abundance of galaxy clusters. Mon. Not. Roy. Astron. Soc. 2012, 424, 1442. Symmetry 2016, 8, 70 77 of 79 316. Trindade, A.M.M.; Avelino, P.P.; Viana, P.T.P. Biased cosmological parameter estimation with galaxy cluster counts in the presence of primordial non-Gaussianities. Mon. Not. Roy. Astron. Soc. 2013, 435, 782. 317. Da Silva, A.C.; Kay, S.T.; Liddle, A.R.; Thomas, P.A. Hydrodynamical simulations of the Sunyaev-Zel’dovich effect: Cluster scaling relations and X-ray properties. Mon. Not. Roy. Astron. Soc. 2004, 348, 1401. 318. Schneider, P.; Ehlers, J.; Falco, E.E. Gravitational Lenses; Springer: Berlin, Germany, 1992. 319. Schneider, P.; Kochanek, C.S.; Wambsganss, J. Gravitational Lensing: Strong, Weak and Micro; Springer: Berlin, Germany, 2006. 320. Narayan, R.; Bartelmann, M. Gravitational Lensing. In Formation of Structure in the Universe; Dekel, A., Ostriker, J.P., Eds.; Cambridge University Press: Cambridge, UK, 1999; p. 360. 321. Bartelmann, M. Gravitational Lensing. Class. Quantum Gravity 2010, 27, 233001. 322. Fort, B.; Mellier, Y. Arc(let)s in clusters of galaxies. Astron. Astrophys. Rev. 1994, 5, 239–292 323. Valdes, F.; Jarvis, J.F.; Tyson, J.A. Alignment of faint galaxy images—Cosmological distortion and rotation. Astrophys. J. 1983, 271, 431–441. 324. Mould, J.; Blandford, R.; Villumsen, J.; Brainerd, T.; Smail, I.; Small, T.; Kells, W. A search for weak distortion of distant galaxy images by large-scale structure. Mon. Not. Roy. Astron. Soc. 1994, 271, 31–38. 325. Schneider, P.; van Waerbeke, L.; Jain, B. Kruse, G. A new measure for cosmic shear. Mon. Not. Roy. Astron. Soc. 1998, 296, 873–892. 326. Van Waerbeke, L.; Mellier, Y.; Erben, T.; Cuillandre, J.C.; Bernardeau, F.; Maoli, R.; Bertin, E.; McCracken, H.J.; Le Fèvre, O.; Fort, B.; et al. Detection of correlated galaxy ellipticities from CFHT data: First evidence for gravitational lensing by large-scale structures. Astron. Astrophys. 2000, 358, 30–44. 327. Bacon, D.J.; Refregier, A.R.; Ellis, R.S. Detection of weak gravitational lensing by large-scale structure. Mon. Not. Roy. Astron. Soc. 2000, 318, 625–640. 328. Wittman, D.M.; Tyson, J.A.; Kirkman, D.; Dell’Antonio, I.; Bernstein, G. Detection of weak gravitational lensing distortions of distant galaxies by cosmic dark matter at large scales. Nature 2000, 405, 143–148. 329. Kaiser, N.; Wilson, G.; Luppino, G. Large-Scale Cosmic Shear Measurements. 2000, arXiv:astro-ph/0003338. 330. Hetterscheidt, M.; Simon, P.; Schirmer, M.; Hildebrandt, H.; Schrabback, T.; Erben, T.; Schneider, P. GaBoDS: The Garching-Bonn deep survey. VII. Cosmic shear analysis. Astron. Astrophys. 2007, 468, 859–876. 331. Fu, L.; Semboloni, E.; Hoekstra, H.; Kilbinger, M.; van Waerbeke, L.; Tereno, I.; Mellier, Y.; Heymans, C.; Coupon, J.; Benabed, K.; et al. Very weak lensing in the CFHTLS wide: Cosmology from cosmic shear in the linear regime. Astron. Astrophys. 2008, 479, 9–25. 332. Erben, T.; Hildebrandt, H.; Miller, L.; van Waerbeke, L.; Heymans, C.; Hoekstra, H.; Kitching, T.D.; Mellier, Y.; Benjamin, J.; Blake, C.; et al. CFHTLenS: The Canada-France-Hawaii Telescope Lensing Survey—Imaging data and catalogue products. Mon. Not. Roy. Astron. Soc. 2013, 433, 2545–2563. 333. Miller, L.; Heymans, C.; Kitching, T.D.; van Waerbeke, L.; Erben, T.; Hildebrandt, H.; Hoekstra, H.; Mellier, Y.; Rowe, B.T.P.; Coupon, J.; et al. Bayesian galaxy shape measurement for weak lensing surveys—III. Application to the Canada-France-Hawaii Telescope Lensing Survey. Mon. Not. Roy. Astron. Soc. 2013, 429, 2858–2880. 334. Benjamin, J.; Van Waerbeke, L.; Heymans, C.; Kilbinger, M.; Erben, T.; Hildebrandt, H.; Hoekstra, H.; Kitching, T.D.; Mellier, Y.; Miller, L.; et al. CFHTLenS tomographic weak lensing: quantifying accurate redshift distributions. Mon. Not. Roy. Astron. Soc. 2013, 431, 1547–1564. 335. Heymans, C.; Van Waerbeke, L.; Miller, L.; Erben, T.; Hildebrandt, H.; Hoekstra, H.; Kitching, T.D.; Mellier, Y.; Simon, P.; Bonnett, C.; et al. CFHTLenS: the Canada-France-Hawaii Telescope Lensing Survey. Mon. Not. Roy. Astron. Soc. 2012, 427, 146–166. 336. Simpson, F.; Heymans, C.; Parkinson, D.; Blake, C.; Kilbinger, M.; Benjamin, J.; Erben, T.; Hildebrandt, H.; Hoekstra, H.; Kitching, T.D. et al. CFHTLenS: Testing the laws of gravity with tomographic weak lensing and redshift-space distortions. Mon. Not. Roy. Astron. Soc. 2012, 429, 2249–2263. 337. Dossett, J.N.; Ishak, M.; Parkinson, D.; Davis, T.M. Constraints and tensions in testing general relativity from Planck and CFHTLenS data including intrinsic alignment systematics. Phys. Rev. D 2015, 92, 023003. 338. Harnois-Déraps, J.; Munshi, D.; Valageas, P.; van Waerbeke, L.; Brax, P.; Coles, P.; Rizzo, L. Testing modified gravity with cosmic shear. Mon. Not. Roy. Astron. Soc. 2015, 454, 2722–2735. 339. Munshi, D.; Valageas, P.; Van Waerbeke, L.; Heavens, A. Cosmology with Weak Lensing Surveys. Phys. Rep. 2008, 462, 67–121. Symmetry 2016, 8, 70 78 of 79 340. Kilbinger, M. Cosmology with cosmic shear observations: A review. Rep. Prog. Phys. 2015, 78, 086901. 341. Ma, C.P.; Bertschinger, E. Cosmological Perturbation Theory in the Synchronous and Conformal Newtonian Gauges. Astrophys. J. 1995, 455, 7. 342. Baker, T.; Ferreira, P.G.; Skordis, C.; Zuntz, J. Towards a fully consistent parameterization of modified gravity. J. Cosmol. Astropart. Phys. 2012, 6, 32. 343. Bertschinger, E. On the Growth of Perturbations as a Test of Dark Energy and Gravity. Astrophys. J. 2006, 648, 797–806. 344. Caldwell, R.; Cooray, A.; Melchiorri, A. Constraints on a new post-general relativity cosmological parameter. Phys. Rev. D 2007, 76, 023507. 345. Hu, W.; Sawicki, I. Parametrized post-Friedmann framework for modified gravity. Phys. Rev. D 2007, 76, 104043. 346. Amendola, L.; Kunz, M.; Sapone, D. Measuring the dark side (with weak lensing). J. Cosmol. Astropart. Phys. 2008, 4, 13. 347. Tereno, I.; Semboloni, E.; Schrabback, T. COSMOS weak-lensing constraints on modified gravity. Astron. Astrophys. 2011, 530, 68. 348. Schrabback, T.; Hartlap, J.; Joachimi, B.; Kilbinger, M.; Simon, P. Evidence of the accelerated expansion of the Universe from weak lensing tomography with COSMOS. Astron. Astrophys. 2010, 516, 63. 349. Hu, W. Power Spectrum Tomography with Weak Lensing. Astrophys. J. Lett. 1999, 522, L21–L24. 350. Daniel, S.F.; Linder, E.V.; Smith, T.L.; Caldwell, R.R.; Cooray, A.; Leauthaud, A.; Lombriser, L. Testing general relativity with current cosmological data. Phys. Rev. D 2010, 81, 123508. 351. Linder, E.V. Cosmic growth history and expansion history. Phys. Rev. D 2005, 72, 043529. 352. Linder, E.V.; Cahn, R.N. Parameterized beyond-Einstein growth. Astropart. Phys. 2007, 28, 481–488. 353. Lewis, A.; Bridle, S. Cosmological parameters from CMB and other data: A Monte Carlo approach. Phys. Rev. D 2002, 66, 103511. 354. Komatsu, E.; Smith, K.M.; Dunkley, J.; Bennett, C.L.; Gold, B.; Hinshaw, G.; Jarosik, N.; Larson, D.; Nolta, M.R.; Page, L.; et al. Seven-Year Wilkinson Microwave Anisotropy Probe (WMAP) Observations: Cosmological Interpretation. Astrophys. J. Suppl. 2011, 192, 18. 355. Eisenstein, D.; Hu, W. Baryonic Features in the Matter Transfer Function. Astrophys. J. 1998, 496, 605–614. 356. Smith, R.E.; Peacock, J.A.; Jenkins, A.; White, S.D.M.; Frenk, C.S.; Pearce, F.R.; Thomas, P.A.; Efstathiou, G.; Couchman, H.M.P. Stable clustering, the halo model and non-linear cosmological power spectra. Mon. Not. Roy. Astron. Soc. 2003, 341, 1311–1332. 357. Valageas, P.; Nishimichi, T.; Taruya, A. Matter power spectrum from a Lagrangian-space regularization of perturbation theory. Phys. Rev. D 2013, 87, 083522. 358. Leonard, C.D.; Baker, T.; Ferreira, P.G. Exploring degeneracies in modified gravity with weak lensing. Phys. Rev. D 2015, 91, 083504. 359. Takahashi, R.; Sato, M.; Nishimichi, T.; Taruya, A.; Oguri, M. Revising the Halofit Model for the Nonlinear Matter Power Spectrum. Astrophys. J. 2012, 761, 152–162. 360. Heitmann, K.; Lawrence, E.; Kwan, J.; Habib, S.; Higdon, D. The Coyote Universe Extended: Precision Emulation of the Matter Power Spectrum. Astrophys. J. 2014, 780, 111–129. 361. Zhao, G.B.; Pogosian, L.; Silvestri, A.; Zylberberg, J. Cosmological Tests of General Relativity with Future Tomographic Surveys. Phys. Rev. Lett. 2009, 103, 241301. 362. Schmidt, F. Weak lensing probes of modified gravity. Phys. Rev. D 2008, 78, 043002. 363. Zhang, P.; Liguori, M.; Bean, R.; Dodelson, S. Probing gravity at cosmological scales by measurements which test the relationship between gravitational lensing and matter overdensity. Phys. Rev. Lett. 2007, 99, 141302. 364. Leonard, C.D.; Ferreira, P.G.; Heymans, C. Testing gravity with EG: Mapping theory onto observations. 2015, arXiv:1510.04287. 365. Reyes, R.; Mandelbaum, R.; Seljak, U.; Baldauf, T.; Gunn, J.E.; Lombriser, L.; Smith, R.E. Confirmation of general relativity on large scales from weak lensing and galaxy velocities. Nature 2010, 464, 256–258. 366. Kalus, B.; Schwarz, D. J.; Seikel, M.; Wiegand, A. Constraints on anisotropic cosmic expansion from supernova. Astron. Astrophys. 2013, 553, A56. Symmetry 2016, 8, 70 79 of 79 367. Wotjak, R.; Knebe, A.; Watson, W.A.; Iliev, I.T.; Heß, S.; Rapetti, D.; Yepes, G.; Gottlöber, S. Cosmic Variance of the Local Hubble Flow in Large-Scale Cosmological Simulations. Mon. Not. Roy. Astron. Soc. 2014, 438, 1805. 368. Bengaly, C.A.P.; Bernui, A.; Alcaniz, J.S. Probing cosmological isotropy with Type Ia Supernova. Astrophys. J. 2015, 808, 39. 369. Cooray, A.; Holz, D.E.; Caldwell, R. Measuring Dark Energy Spatial Inhomogeneity with Supernova Data. J. Cosmol. Astropart. Phys. 2010, 11, 15. 370. Valkenburg, W.; Marra, V.; Clarkson, C. Testing the Copernican principle by constraining spatial homogeneity. Mon. Not. Roy. Astron. Soc. 2014, 438, 6. 371. Betoule, M.; Kessler, R.; Guy, J.; Mosher, J.; Hardin, D.; Biswas, R.; Astier, P.; El-Hage, P.; Konig, M.; Kuhlmann, S.; et al. Improved Cosmological Constraints from a Joint Analysis of the SDSS–II and SNLS Supernova Samples. Astron. Astrophys. 2014, 568, 32. 372. Górski, K.M.; Hivon, E.; Banday, A.J.; Wandelt, B.D.; Hansen, F.K.; Reinecke, M.; Bartelmann, M. HEALPix: A Framework for High-Resolution Discretization and Fast Analysis of Data Distributed on the Sphere. Astrophys. J. 2005, 622, 759. 373. Räsänen, S. Accelerated Expansion from Structure Formation. J. Cosmol. Astropart. Phys. 2006, 11, 3. 374. Carvalho, C.S.; Basilakos, S. Angular distribution of cosmological parameters as a probe of inhomogeneities: A kinematic parametrisation. 2016, arXiv:1603.07519 [astro-ph.CO]. c 2016 by the authors; licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC-BY) license (http://creativecommons.org/licenses/by/4.0/).