# F(R) theories of gravitation

Salvatore Capozziello and Mariafelicia De Laurentis (2015), Scholarpedia, 10(2):31422. | doi:10.4249/scholarpedia.31422 | revision #147843 [link to/cite this article] |

**$f(R)$ gravity** is an extension of **Einstein's General Relativity** derived from relaxing the hypothesis that the Hilbert-Einstein action for the gravitational field is strictly linear in the Ricci curvature scalar $R$, *i.e.* $f(R)= R$. In this sense, **$f(R)$ gravity** represents a class of theories defined as arbitrary functions of $R$. It can be considered as the simplest example of **Extended Theory of Gravity** (Capozziello and De Laurentis, 2011).

## $f(R)$ gravity: Introduction

In **$f(R)$ gravity**, the Hilbert-Einstein action:
\begin{eqnarray}\tag{1}{\cal A}[g]=\int R \sqrt{-g}\,\mathrm{d}^4x\,, \end{eqnarray}
is generalized as:
\begin{eqnarray}\tag{2}{\cal A}[g]=\int f(R) \sqrt{-g} \, \mathrm{d}^4x\,, \end{eqnarray}
where we are using physical units $8\pi G c^{-4}=1$, the determinant of the metric tensor is $g=|g_{\mu\nu}|$ and *f*(*R*) is a function of the Ricci curvature scalar $R$.
This approach is aimed to address problems and shortcomings coming from the Standard Cosmological Model and from the lack of a final theory of Quantum Gravity. Due to these issues, alternative theories of gravity have been considered in order to attempt, at least, a semi-classical scheme where General Relativity and its positive results has to be recovered. One of the most fruitful approaches is that of **Extended Theories of Gravity** (see (Capozziello and De Laurentis, 2011) and references therein) which have become a sort of paradigm in the study of gravitational interaction. They are based on corrections and enlargements of the Einstein theory. The paradigm consists, essentially, in adding higher-order curvature invariants and minimally or non-minimally coupled scalar fields into dynamics related to the effective actions of some fundamental theory (Nojiri and Odintsov, 2011) (Nojiri and Odintsov, 2007) (Capozziello and Francaviglia, 2007) (De Felice and Tsujikawa, 2010) (Sotiriou and Faraoni, 2010). This approach is coherent to the fact that these generalized theories emerge, like Einstein's gravity, from the Gauge Theory, and can be framed into a general bundle structure (Capozziello and De Laurentis, 2010). Other motivations to modify General Relativity come from the Mach Principle which leads to assume a varying gravitational coupling. This principle states that the local inertial frame is determined by some average of the motion of distant astronomical objects (Bondi, 1952). This fact implies that the gravitational coupling can be scale-dependent and related to some scalar field. As a consequence, the concept of inertia and the Equivalence Principle have to be revised.
Besides, any unification scheme as Superstrings, Supergravity or
Grand Unified Theories, takes into account effective actions where
non-minimal couplings to the geometry or higher-order terms of the
curvature invariants are present. Such contributions are due to
one-loop or higher-loop corrections in the strong gravity regimes
near the full (not yet available) Quantum Gravity
(Nojiri and Odintsov, 2007). Specifically, this scheme was adopted in order to
deal with the quantization on curved space-times and the result was
that the interactions among quantum scalar fields and background
geometry or the gravitational self-interactions yield corrective
terms in the Hilbert-Einstein Lagrangian (Birrell and Davies, 1982). Moreover,
it has been realized that such corrective terms are inescapable in
order to obtain the effective action of Quantum Gravity at scales
closed to the Planck one (Vilkovisky, 1992). All these approaches
are not the full Quantum Gravity but are needed as
working schemes towards it.
In summary, higher-order terms in curvature invariants (such as
$R^{2}$, $R_{\mu\nu} R^{\mu\nu}$,
$R^{\mu\nu\alpha\beta}R_{\mu\nu\alpha\beta}$, $R \,\Box R$, or $R
\,\Box^{k}R$) or non-minimally coupled terms between scalar fields
and geometry (such as $\phi^{2}R$) have to be added to the effective
Lagrangian of gravitational field if quantum corrections are
considered. For instance, one can notice that such terms occur in
the effective Lagrangian of strings or in Kaluza-Klein theories,
if the mechanism of dimensional reduction is adopted (Gasperini and Veneziano, 1992).

On the other hand, from a conceptual point of view, there are no *a priori* reason to restrict the gravitational Lagrangian to a
linear function of the Ricci scalar $R$ minimally coupled with
matter (Capozziello and Francaviglia, 2007). Furthermore, the idea that there are
*no exact laws of physics* could be taken into serious account:
in such a case, the effective Lagrangians of physical interactions
are *stochastic functions*. This feature means that the local
gauge invariances (*i.e*. conservation laws) are well approximated
only in the low energy limit and the fundamental physical
constants can vary (Barrow and Ottewil, 1983).

Beside fundamental physics motivations, all these theories
acquired a huge interest in cosmology due to the fact that they
*naturally* exhibit inflationary behaviours able to overcome the
shortcomings of Standard Cosmological Model (based on General Relativity). The
related cosmological models seem realistic and capable of
matching with the Cosmic Microwave Background Radiation observations (Starobinsky, 1980) (Duruisseau and Kerner, 1983) (La and Steinhardt, 1989).
Furthermore, it is possible to show that, via conformal
transformations, the higher-order and non-minimally coupled terms correspond to the Einstein gravity plus one or more than
one minimally coupled scalar fields
(Teyssandier and Tourrenc, 1983) (Maeda, 1989) (Wands, 1994) (Capozziello, de Ritis and Marino, 1998) (Adams, Freese and Guth, 1991).

More precisely, higher-order terms appear as contributions of order two in the derivatives of field equations. For example, a term like $R^{2}$ gives fourth order equations (Ruzmaikina and Ruzmaikin, 1970), $R \ \Box R$ gives sixth order equations (Amendola et al., 1993) (Adams, Freese and Guth, 1991) (Buchdahl, 1951), $R\,\Box^{2}R$ gives eighth order equations (Battaglia-Mayer and Schmidt, 1993) and so on. By a conformal transformation, any 2nd-order derivative term corresponds to a scalar field. The dynamics of any scalar field is given by the corresponding Klein-Gordon equation, which is second order. For example, fourth-order gravity gives Einstein's equations plus one scalar field, sixth-order gravity gives Einstein's equations plus two scalar fields and so on (Schmidt, 1990).

Considering the mathematical point of view, the problem of reducing more general theories to the Einstein standard form has been extensively treated; one can see that, through a Legendre transformation, higher-order theories, under suitable regularity conditions on the Lagrangian, take the form of the Einstein one in which a scalar field (or more than one) is the source of the gravitational field (see for example (Capozziello and Francaviglia, 2007) (Sokolowski, 1989) (Ferraris, Francaviglia and Magnano, 1989) (Magnano and Sokolowski, 1994) (Jakubiec and Kijowski, 1988)).

In any case, the debate on the physical meaning of conformal transformations is far to be solved (see (Capozziello and Faraoni, 2010) (Capozziello and De Laurentis, 2010) (Capozziello, De Laurentis and Faraoni, 2009) (Faraoni, 2008), and references therein for a comprehensive review). Several authors claim for a true physical difference between Jordan frame (higher-order theories and/or variable gravitational coupling) since there are experimental and observational evidences which point out that the Jordan frame could be suitable to better match solutions with data. Others state that the true physical frame is the Einstein one according to the energy theorems (Magnano and Sokolowski, 1994). However, the discussion is open and no definitive statement has been formulated up to now.

The problem can be faced from a more general point of view and the
**Palatini approach** to gravity, (see below) could be useful in this context (Einstein, 1925) (Buchdahl, 1979).
In (Olmo, 2011), this approach is widely discussed for **Extended Theories of Gravity** and several important applications are reported. The fundamental idea of the Palatini formalism is to consider the connection $\Gamma$, entering the
definition of the Ricci tensor, to be independent of the metric
$g$ defined on the space-time ${\cal M}$. The Palatini formulation
for the Hilbert-Einstein theory results to be equivalent
to the purely metric theory: this follows from the fact that the
field equations for the connection $\Gamma$, firstly considered
to be independent of the metric, give the Levi-Civita connection
of the metric $g$. As a consequence, there is no reason to impose
the Palatini variational principle in the
Hilbert-Einstein theory instead of the metric variational
principle.
However, the situation completely changes if we consider theories of gravity
depending on functions of curvature invariants, as $f(R)$, or
non-minimally coupled to some scalar field. In these cases, the
Palatini and the metric variational principle provide different
field equations and the theories thus derived differ (Magnano and Sokolowski, 1994) (Ferraris, Francaviglia and Volovich, 1994). The relevance of the Palatini approach, in
this framework, has been recently proven in relation to
cosmological applications
(Vollick, 2003) (Li and Chu, 2006) (Li, Chan and Chu, 2007) (Capozziello, De Laurentis, Francaviglia and Mercadante, 2009) (Olmo, 2011).

From a physical
point of view, considering the metric $g$ and the connection $\Gamma$
as independent fields means to decouple the metric structure of
space-time and its geodesic structure (being, in general, the
connection $\Gamma$ not the Levi-Civita connection of $g$). The
chronological structure of space-time is governed by $g$ while the
trajectories of particles, moving in the space-time, are governed
by $\Gamma$.
This decoupling enriches the geometric structure of space-time and
generalizes the purely metric formalism. The metric-affine
structure of space-time is naturally translated, by the
Palatini field equations, into a bi-metric structure. Beside the *physical* metric $g$, another metric
$\hat{g}$ is involved. This new metric is related, in the case of **$f(R)$ gravity**, to the connection. As a matter of fact, the
connection $\Gamma$ results to be the Levi-Civita connection of
$\hat{g}$ and thus provides the geodesic structure (Allemandi, Capone, Capozziello and Francaviglia, 2006).

If we consider the case of non-minimally coupled interaction in the gravitational Lagrangian (scalar-tensor theories), the new metric $\hat{g}$ is related to the non-minimal coupling. The new metric $\hat{g}$ can be thus related to different geometric and physical aspects of the gravitational theory. Thanks to the Palatini formalism, the non-minimal coupling and the scalar field, entering the evolution of the gravitational fields, are separated from the metric structure of space-time. The situation mixes when we consider the case of higher-order-scalar-tensor theories. Due to these features, the Palatini approach could greatly contribute to clarify the physical meaning of conformal transformation (Allemandi, Capone, Capozziello and Francaviglia, 2006).

## $f(R)$ gravity in metric formalism

In metric **$f(R)$ gravity**, the field equations are obtained by varying with respect to the metric and not treating the connection independently. The main steps are the same as in the case of the variation of the Hilbert-Einstein action but there are also some important differences.

The variation of the determinant is
\begin{eqnarray}\tag{3}\delta \sqrt{-g}= -\frac{1}{2} \sqrt{-g} g_{\mu\nu} \delta g^{\mu\nu}\,.\end{eqnarray}
The Ricci scalar is defined as
\begin{eqnarray}\tag{4} R = g^{\mu\nu} R_{\mu\nu}.\end{eqnarray}
Therefore the variation with respect to the inverse metric *g ^{μν}* is given by
\[
\begin{align}\tag{5}
\delta R &= R_{\mu\nu} \delta g^{\mu\nu} + g^{\mu\nu} \delta R_{\mu\nu}
= R_{\mu\nu} \delta g^{\mu\nu} + g^{\mu\nu}(\nabla_\rho \delta \Gamma^\rho_{\nu\mu} - \nabla_\nu \delta \Gamma^\rho_{\rho\mu})\,.
\end{align}
\]
Since

*δΓ*is the difference of two connections, it should transform as a tensor. Therefore, it can be written as \[\tag{6}\delta \Gamma^\lambda_{\mu\nu}=\frac{1}{2}g^{\lambda \alpha}\left(\nabla_\mu\delta g_{\alpha\nu}+\nabla_\nu\delta g_{\alpha\mu}-\nabla_\alpha\delta g_{\mu\nu} \right).\] Substituting into the above equation: \[\tag{7}\delta R= R_{\mu\nu} \delta g^{\mu\nu}+g_{\mu\nu}\Box \delta g^{\mu\nu}-\nabla_\mu \nabla_\nu \delta g^{\mu\nu}\,,\] where ∇

^{λ}_{μν}_{μ}is the covariant derivative and $\Box= g^{\mu\nu} \nabla_{\mu}\nabla_{\nu}$ is the D'Alembert operator. We adopt the metric signature $(+,-,-,-)$. Let us consider now a generic analytic function $f(R)$ obeying the variational principle $ \delta \int d^{4}x \sqrt{-g} \, f(R) =0$. We have \begin{eqnarray} && \delta \int d^{4}x \sqrt{-g} \, f(R) = \int d^{4}x \left[ \delta \left( \sqrt{-g} \, f(R) \right) + \sqrt{-g} \, \delta\left(f(R) \right)\right] \nonumber\\ &&\nonumber\\ && = \int d^{4}x \sqrt{-g} \, {\textstyle\left[f'(R)R_{\mu\nu} -\frac{1}{2}g_{\mu\nu}f(R)\right]}\delta g^{\mu\nu} + \int d^{4}x \sqrt{-g} \, f'(R)g^{\mu\nu}\delta R_{\mu\nu} \, ,\nonumber\\ && \tag{8} \end{eqnarray} where the prime denotes differentiation with respect to $R$. These integrals can be evaluated in the local inertial frame. By using \begin{eqnarray}\tag{9} g^{\mu\nu}\delta R_{\mu\nu} =g^{\mu\nu}\partial_{\sigma} \left(\delta G_{\mu\nu}^{\sigma}\right) -g^{\mu\sigma}\partial_{\sigma}\left(\delta G_{\mu\nu}^{\nu}\right) \equiv\partial_{\sigma}W^{\sigma}\,, \end{eqnarray} where \begin{eqnarray}\tag{10} W^{\sigma}\equiv g^{\mu\nu} \delta G_{\mu\nu}^{\sigma} -g^{\mu\sigma}\delta G_{\mu\nu}^{\nu} \,, \end{eqnarray} the second integral in Eq. (8) can be written as \begin{eqnarray}\tag{11} \int d^{4}x \sqrt{-g} \, f'(R)g^{\mu\nu}\delta R_{\mu\nu} =\int d^{4}x \sqrt{-g} \, f'(R) \partial_{\sigma}W^{\sigma} \, . \end{eqnarray} Integration by parts yields \begin{eqnarray}\tag{12} \int d^{4}x \sqrt{-g} \, f'(R)g^{\mu\nu}\delta R_{\mu\nu} &= & \int d^{4}x \frac{\partial}{\partial x^{\sigma} } \left[ \sqrt{-g} \, f'(R)W^{\sigma}\right] - \int d^{4}x \partial_{\sigma}\left[\sqrt{-g}f'(R)\right]W^{\sigma}\, . \end{eqnarray} The first integrand is a total divergence and can be discarded assuming that fields vanish at infinity. We have \begin{eqnarray}\tag{13} \int d^{4}x \sqrt{-g} \, f'(R) g^{\mu\nu}\delta R_{\mu\nu} = -\int d^{4}x \partial_{\sigma} \left[\sqrt{-g} \, f'(R)\right]W^{\sigma} \,. \end{eqnarray}

Let us evaluate $W^{\sigma}$ appearing in Eq. (13). We have \begin{eqnarray}\tag{14} \delta G_{\mu\nu}^{\sigma} = \delta\left[\frac{1}{2}g^{\sigma\alpha}\left(\partial_{\mu}g_{\alpha\nu} +\partial_{\nu}g_{\mu\alpha}-\partial_{\alpha}g_{\mu\nu}\right)\right] = \frac{1}{2}g^{\sigma\alpha}\left[\partial_{\mu}\left(\delta g_{\alpha\nu}\right)+\partial_{\nu}\left(\delta g_{\mu\alpha}\right)-\partial_{\alpha}\left(\delta g_{\mu\nu}\right)\right] \,. \end{eqnarray}

In the local inertial frame, it is
\begin{eqnarray}\tag{15}
\partial_{\alpha}g_{\mu\nu}=\nabla_{\alpha}g_{\mu\nu}= 0\,,
\end{eqnarray}
and
\begin{eqnarray}\tag{16}
\delta G_{\mu\nu}^{\nu}=\frac{1}{2} \, g^{\nu\alpha}\partial_{\mu}\left(\delta
g_{\nu\alpha}\right) \,.
\end{eqnarray}
By combining Eqs. (15) and
(16), one obtains
\begin{eqnarray}
g^{\mu\nu}\delta G_{\mu\nu}^{\sigma}
& = & \frac{1}{2} \, g^{\mu\nu}\left[
-\partial_{\mu}\left(g_{\alpha\nu}\delta
g^{\alpha\sigma}\right)-\partial_{\nu}\left(g_{\mu\alpha}\delta
g^{\sigma\alpha}\right)-g^{\sigma\alpha}\partial_{\alpha}\left(\delta
g_{\mu\nu}\right)\right] = \frac{1}{2} \, \partial^{\sigma}\left(g_{\mu\nu}\delta g^{\mu\nu}
\right)-\partial^{\mu}\left(g_{\alpha\mu}\delta g^{\nu\alpha}\right) \,,
\tag{17} \\
&&\nonumber\\
g^{\mu\sigma}\delta G_{\mu\nu}^{\nu} & = &
-\frac{1}{2} \, \partial^{\sigma}\left(g_{\nu\alpha}\delta g^{\nu\alpha}\right)
\,,\tag{18}
\end{eqnarray}
from which it follows that
\begin{eqnarray}\tag{19}
W^{\sigma}=\partial^{\sigma}\left(g_{\mu\nu}\delta
g^{\mu\nu}\right)-\partial^{\mu}\left(g_{\mu\nu}\delta
g^{\sigma\nu}\right) \,.
\end{eqnarray}
Using this equation, one can write
\begin{eqnarray}
\int d^{4}x \sqrt{-g} \, f'(R)g^{\mu\nu}\delta
R_{\mu\nu}
= \int d^{4}x
\partial_{\sigma}\left[\sqrt{-g} \,
f'(R)\right] \left[\partial^{\mu}\left(g_{\mu\nu}\delta
g^{\sigma\nu}\right)-\partial^{\sigma}\left(g_{\mu\nu}\delta
g^{\mu\nu}\right)\right] \,. \tag{20}
\end{eqnarray}
Integrating by parts and discarding total divergences, one
obtains
\begin{eqnarray}
\int d^{4}x \, \sqrt{-g} \, f'(R)g^{\mu\nu}\delta
R_{\mu\nu} = \int d^{4}x \,
g_{\mu\nu} \partial^{\sigma}\partial_{\sigma}\left[\sqrt{-g}f'(R)\right]\delta
g^{\mu\nu}
- \int d^{4}x \,
g_{\mu\nu}\partial^{\mu}\partial_{\sigma}\left[\sqrt{-g} \, f'(R)\right]
\delta g^{\sigma\nu} \,. \tag{21}
\end{eqnarray}
The variation of the action is then
\begin{eqnarray}\tag{22}
\delta\int d^{4}x \sqrt{-g} \, f(R)
=\int d^{4}x \sqrt{-g}\,\left[
f'(R)R_{\mu\nu}-\frac{1}{2}
f(R) g_{\mu\nu}\right]\delta g^{\mu\nu}+\int d^{4}x \left[
g_{\mu\nu}\partial^{\sigma}\partial_{\sigma}
\left(\sqrt{-g}\, f'(R)\right)-g_{\sigma\nu}\partial^{\mu}
\partial_{\sigma}\left(\sqrt{-g}f'(R)\right)
\right]\delta g^{\mu\nu} \,.
\end{eqnarray}
The vanishing of the variation implies the fourth order
vacuum field
equations
\begin{eqnarray}\tag{23}
f'(R)R_{\mu\nu}-\frac{f(R)}{2} \, g_{\mu\nu} =\nabla_{\mu}
\nabla_{\nu}f'(R)-g_{\mu\nu}\Box
f'(R)\,.
\end{eqnarray}
These equations can be rearranged in the Einstein-like form
\begin{eqnarray}\tag{24}
f'(R) R_{\mu\nu}-\frac{ f'(R)}{2} \, g_{\mu\nu} R
+\frac{f'(R)}{2} \, g_{\mu\nu}R-\frac{f(R)}{2} \, g_{\mu\nu}=
\nabla_{\mu} \nabla_{\nu} f'(R)-
g_{\mu\nu}\Box f'(R)\,,
\end{eqnarray}
and then
\begin{eqnarray}\tag{25}
G_{\mu\nu}=\frac{1}{f'(R)} \left\{
\nabla_{\mu}\nabla_{\nu} f'(R) - g_{\mu\nu}\Box f'(R)
+ g_{\mu\nu} \frac{ \left[ f(R)-f'(R) R \right]}{2}\right\}\,.
\end{eqnarray}
The right hand side of Eq. (25) is then regarded as
an effective stress-energy tensor $T_{\mu\nu}^{(curv)}$, which we call
**curvature energy-momentum tensor**,
sourcing the effective Einstein equations. Although this
interpretation is questionable, in principle, because the field
equations describe a theory different from General Relativity, and one is
forcing upon them the interpretation as effective Einstein
equations, this approach is quite useful in practice.

## The case of $f(R)=R +\alpha R^2$

Quadratic corrections in the Ricci scalar, motivated by the attempts to renormalize General Relativity, constitute a straightforward extension of General Relativity and are particularly relevant in cosmology since they allow the construction of a self-consistent inflationary model (Starobinsky, 1980). Let us derive the field equations for the Lagrangian density \begin{eqnarray} \tag{26} {\cal L}=R+\alpha R^{2}+{\cal L}^{(m)}\,, \end{eqnarray} from the variational principle $\delta\int d^{4}x \sqrt{-g} \, {\cal L}=0 $. We consider first the vacuum case. The variation gives \begin{eqnarray} \tag{27} \int d^{4}x \sqrt{-g} \, G_{\alpha\beta}\delta g^{\alpha\beta} + \alpha \, \delta\int d^{4}x \sqrt{-g} \, R^{2} =0\,, \end{eqnarray} where the variation of $R\sqrt{-g} $ produces the Einstein tensor. The second term on the right hand side of Eq. (27) gives \begin{eqnarray} \tag{28} \delta \int d^{4}x \sqrt{-g} \, R^{2} = -\frac{1}{2}\int d^{4}x \sqrt{-g} \,\, g_{\alpha\beta}\delta g^{\alpha\beta}R^{2} + 2\int d^{4}x \sqrt{-g} \, R\delta R \end{eqnarray} and \begin{eqnarray} \tag{29} \int d^{4}x \sqrt{-g} \, R\delta R = \int d^{4}x \sqrt{-g} \, R\left(\delta g^{\alpha\beta} R_{\alpha\beta} + g^{\alpha\beta}\delta R_{\alpha\beta} \right) \,. \end{eqnarray} By using \begin{eqnarray}\tag{30} g^{\alpha\beta}\delta R_{\alpha\beta} = \nabla_{\alpha}\nabla_{\beta} h^{\alpha\beta} -\Box h \,, \end{eqnarray} where \begin{eqnarray} \tag{31} h^{\alpha\beta} \equiv -\delta g^{\alpha\beta},\ \;\;\;\;\; h \equiv -g_{\alpha\beta}\delta g^{\alpha\beta} \,, \end{eqnarray} one has \begin{eqnarray}\tag{32} \int d^{4}x \sqrt{-g} \, R \, g^{\alpha\beta}\delta R_{\alpha\beta} = \int d^{4}x \sqrt{-g} \, R \left( \nabla_{\alpha}\nabla_{\beta} h^{\alpha\beta} -\Box h\right) \,. \end{eqnarray} Integrating by parts twice, the operators $\nabla_{\alpha}\nabla_{\beta}$ and $\Box$ acting on $ h^{\alpha\beta}$ and $h$, respectively, transfer their action onto $R$ and \begin{eqnarray}\tag{33} \int d^{4}x \sqrt{-g} \, R \, g^{\alpha\beta}\delta R_{\alpha\beta} = \int d^{4}x \sqrt{-g} \, \left(h ^{\alpha\beta} \nabla_{\alpha}\nabla_{\beta} R -h\Box R\right) \,. \end{eqnarray} Using Eq.(31), Eq. (33) becomes \begin{eqnarray}\tag{34} \int d^{4}x \sqrt{-g} \, R \, g^{\alpha\beta}\delta R_{\alpha\beta} =\int d^{4}x \sqrt{-g} \left(-\delta g^{\alpha\beta} \nabla_{\alpha}\nabla_{\beta} R +g_{\alpha\beta}\Box R\delta g^{\alpha\beta}\right) \,. \end{eqnarray} Upon substitution of Eq.(34) into Eq.(29), one obtains \begin{eqnarray} \tag{35} \int d^{4}x \sqrt{-g} \, R\delta R = \int d^4x \, \sqrt{-g} \left(R\delta g^{\alpha\beta}R_{\alpha\beta}-\delta g^{\alpha\beta} \nabla_{\alpha}\nabla_{\beta} R +g_{\alpha\beta}\Box R\delta g^{\alpha\beta}\right)\,, \end{eqnarray} and Eq.(28) takes the form \begin{eqnarray} \delta \int d^{4}x \sqrt{-g} \, R^{2} & = & -\frac{1}{2}\int d^{4}x \sqrt{-g} \, g_{\alpha\beta}\delta g^{\alpha\beta} R^{2}+ 2\int d^{4}x \sqrt{-g} \, \left(R\delta g^{\alpha\beta}R_{\alpha\beta}-\delta g^{\alpha\beta} \nabla_{\alpha}\nabla_{\beta} R +g_{\alpha\beta} \Box R\delta g^{\alpha\beta}\right) \nonumber\\ &&\nonumber\\ & =& \int d^{4}x \sqrt{-g} \, {\textstyle\left(2RR_{\alpha\beta} -\frac{1}{2} \, g_{\alpha\beta}R^{2}\right)}\delta g ^{\alpha\beta} + 2 \int d^{4}x \sqrt{-g} \, \left(g_{\alpha\beta}\Box R- \nabla_{\alpha}\nabla_{\beta} R \right)\delta g^{\alpha\beta}\,. \tag{36} \end{eqnarray} Substituting this equation into Eq.(27) and including the matter part of the Lagrangian ${\cal L}^{(m)}$ which produces the energy-momentum tensor $T_{\mu\nu}^{(m)}$, the field equations \begin{eqnarray}\tag{37} G_{\alpha\beta} + \alpha\left[ 2R\left( R_{\alpha\beta} -\frac{1}{4} g_{\alpha\beta}R\right) + 2\left( g_{\alpha\beta}\Box R-\nabla_{\alpha}\nabla_{\beta} R\right)\right] = T_{\alpha\beta}^{(m)}\,, \end{eqnarray} are obtained.

The trace of Eq. (37) is
\begin{eqnarray}\tag{38}
\Box R -\frac{1}{6\alpha}\left(R+T^{(m)} \right)=0\,.
\end{eqnarray}
One can also define an angular frequency $\omega$
(equivalent to a mass $m$) so that
\begin{eqnarray}\tag{39}
-\frac{1}{6\alpha}=\omega^{2}=m^{2} \,.
\end{eqnarray}
Following this definition,
Eq. (38) becomes
\begin{eqnarray}\tag{40}
\Box R + m^{2}\left(R+T^{(m)} \right)=0 \,.
\end{eqnarray}
Eq. (40) can be seen
as an effective Klein-Gordon equation for the effective scalar
field degree of freedom $R$ (sometimes called *scalaron*).

## Equivalent formalism

We can simplify the analysis of *f*(*R*) theories by introducing an auxiliary field Φ. Assuming $f''(R)\neq 0$ for all *R*, let *V*(Φ) be the Legendre transform of *f*(*R*) so that Φ = *f′*(*R*) and *R* = *V′*(Φ). One obtains the action (O'Hanlon, 1972)

\[\tag{41}{\cal A} = \int d^4x \sqrt{-g} \left[ \Phi R - V(\Phi) + \mathcal{L}_{\text{m}}\right].\]

The Euler–Lagrange equations are
\[\tag{42}V'(\Phi)=R\,,\]
\[\tag{43}\Phi \left( R_{\mu\nu} - \frac{1}{2}g_{\mu\nu} R \right) + \left(g_{\mu\nu}\Box -\nabla_\mu \nabla_\nu \right) \Phi + \frac{1}{2} g_{\mu\nu}V(\Phi) = T_{\mu\nu}\,.\]
Substituting Φ with *f′*(*R*) , we obtain exactly the same equations as before. However, the equations are now second order in the derivatives, instead of fourth order.
By performing a conformal rescaling
\[\tag{44}\tilde{g}_{\mu\nu}=\Phi g_{\mu\nu},\]
we transform to the Einstein frame:
\[\tag{45}R=\Phi^{-1} \left[ \tilde{R} + \frac{3\tilde{\Box} \Phi}{\Phi} +\frac{9}{2}\left(\frac{\tilde{\nabla} \Phi}{\Phi}\right)^2 \right]\,,\]
and then
\[\tag{46}{\cal A} = \int d^4x \sqrt{-\tilde{g}}\left[ -\frac{1}{2}\tilde{R} + \frac{3}{2}\left( \frac{\tilde{\nabla}\Phi}{\Phi} \right)^2 - \frac{V(\Phi)}{\Phi^2} \right]\,,\]
after integrating by parts.

We define \[\tag{47}\tilde{\Phi} = \sqrt{3} \ln{\Phi}\,,\] and substitute \[\tag{48}{\cal A} = \int \mathrm{d}^4x \sqrt{-\tilde{g}}\left[ -\frac{1}{2}\tilde{R} + \frac{1}{2}\left(\tilde{\nabla}\tilde{\Phi}\right)^2 - \tilde{V}(\tilde{\Phi}) \right]\,,\] where \[\tag{49}\tilde{V}(\tilde{\Phi}) = e^{-2/\sqrt{3}\;\tilde{\Phi}}V(e^{\tilde{\Phi}/\sqrt{3}}).\] This is General Relativity minimally coupled to a scalar field.

## $f(R)$ gravity in the Palatini formalism

The fundamental idea of the **Palatini formalism** is to consider the
connection $\Gamma^{\mu}_{\nu\alpha}$,
entering the definition of the Ricci tensor, as a
variable independent
of the space-time metric $ g_{\mu\nu}$. The Palatini
formulation of General Relativity is equivalent
to the metric version of this theory as a consequence of the
fact that the field equations for the connection
$\Gamma^{\alpha}_{\mu\nu}$ give exactly the Levi-Civita connection
of the metric $g_{\mu\nu}$ (Wald, 1984). As a consequence, there
is no particular reason to impose the Palatini variational principle in General Relativity instead of the
metric variational principle.

The situation is different in **Extended Theories of Gravity** depending on functions of
curvature invariants or
non-minimally coupled scalar
fields. In these cases, the
Palatini and the metric variational
principle yield different
field equations (Ferraris, Francaviglia and Volovich, 1994) (Magnano and Sokolowski, 1994).

The **Palatini approach**, in the context of **Extended Theories of Gravity** has been the
subject of much interest in
cosmological applications (Olmo, 2011).
As discussed above, considering the metric $g_{\mu\nu}$ and
the connection $\Gamma^{\alpha}_{\mu\nu}$ as independent fields
amounts to decoupling the metric structure of space-time and its
geodesic structure with the connection
$\Gamma^{\alpha}_{\mu\nu}$ being distinct from the Levi-Civita
connection of $g_{\mu\nu}$. In principle, this decoupling
enriches the geometric structure of space-time and generalizes
the purely metric formalism. By means of the Palatini field
equations, this dual structure of space-time is naturally
translated into a bimetric structure of the theory: instead of a
metric and an independent connection, the **Palatini **
formalism** can **
be seen as containing two independent metrics $g_{\mu\nu}$ and
${\hat g}_{\mu\nu}=f'(R) \, g_{\mu\nu}$. In Palatini **$f(R)$ gravity**,
the new metric
${\hat g}_{\mu\nu}$ determining the geodesics is related to the
connection $\Gamma^{\alpha}_{\mu\nu}$ by the fact that
the latter turns out to be the Levi-Civita connection of
${\hat g}_{\mu\nu}$.
In scalar-tensor gravity, the second metric ${\hat g}_{\mu\nu}$ is
related to the non-minimal coupling
of the Brans-Dicke-like
scalar. In the Palatini formalism the
non-minimal
coupling and
the scalar field are separated from the
metric structure of space-time. Physical consequences of this fact are discussed in (Capozziello, De Laurentis, Francaviglia and Mercadante, 2009) (Amendola, Enqvist and Koivisto, 2011).
However, also other geometrical invariants, besides $R$, can be considered in the Palatini formalism. In (Li, Mota and Shaw, 2008), microscopic and macroscopic behaviors of Palatini modified gravity theories are discussed, in particular a detailed study of $f(R,R^{\mu\nu}R_{\mu\nu})$ models is reported. In (Olmo, Sanchis-Alepuz and Tripathi, 2009) dynamical aspects of Palatini $f(R)$, $f(R^{\mu\nu}R_{\mu\nu})$, and $f(R,R^{\mu\nu}R_{\mu\nu})$ theories are studied. Isotropic and anisotropic bouncing cosmologies in Palatini $f(R)$ and $f(R^{\mu\nu}R_{\mu\nu})$ theories are discussed in (Barragan and Olmo , 2010). A Lagrangian of the type $\mathcal{L}=R+\alpha R_{\mu\nu}R^{\mu\nu} $ is also studied in Palatini formalism in the paper (Buchdahl, 1979).

In the **Palatini formalism**, the Ricci scalar in $f({\cal R})$ is
${\cal R} \equiv {\cal R}( g,\Gamma) \equiv
g^{\alpha\beta} {\cal R}_{\alpha
\beta}(\Gamma )$ and is a generalized Ricci scalar, whereas
$ {\cal R}_{\mu \nu}(\Gamma
)$ is the Ricci tensor of a torsion-free connection
$\Gamma^{\alpha}_{\mu\nu}$ which, *a priori*, has no relations with the space-time metric
$g_{\mu\nu}$ which gives the Ricci scalar $R$. The gravitational sector of the theory is
described by the function $f({\cal R})$, while
$\sqrt{-g}$ denotes the usual scalar density.
The field equations derived with the Palatini variational
principle are
\begin{eqnarray}
&& f^{\prime }({\cal
R}) {\cal R}_{(\mu\nu)}(\Gamma)-\frac{ f({\cal R})}{2}\, g_{\mu
\nu }=T^{(m)}_{\mu\nu}\,, \tag{50} \\
&&\nonumber\\
&& \nabla _{\alpha }^{\Gamma } \left[ \sqrt{-g} \, f' ({\cal
R} )g^{\mu \nu } \right] =0\,, \tag{51}
\end{eqnarray}
where $\nabla^{\Gamma}_{\mu} $ is the covariant derivative of
the non-metric connection $\Gamma^{\alpha}_{\mu\nu}$.

It is important to stress that Eq.(51) is obtained under the assumption that the matter sector described by ${\cal L}^{(m)}$ is functionally independent of the (non-metric) connection $\Gamma^{\alpha}_{\mu\nu}$; however it may contain metric covariant derivatives $\stackrel{g}{\nabla}$ of the matter fields. This means that the matter stress-energy tensor $ T^{(m)}_{\mu\nu} \left[ g,\Psi \right]$ depends on the metric $g_{\mu\nu}$ and on the matter fields collectively denoted by $\Psi$, together with their covariant derivatives with respect to the Levi-Civita connection of $g_{\mu\nu}$. It is easy to see, from Eq.(51), that $\sqrt{-g} \, f' ({\cal R} ) g^{\mu \nu }$ is a symmetric tensor density of weight $1$, which naturally leads to the introduction of a new metric ${\hat g}_{\mu \nu}$ conformally related to $g_{\mu\nu}$ by (Allemandi, Capone, Capozziello and Francaviglia, 2006) (Ferraris, Francaviglia and Volovich, 1994)

\begin{equation}\tag{52} \sqrt{-g} \, f' ({\cal R}) \, g^{\mu \nu} = \sqrt{-{\hat g}} \,\, {\hat g}^{\mu \nu }\,. \end{equation} With this definition $\Gamma^{\alpha}_{\mu\nu} $ is the Levi-Civita connection of the metric ${\hat g}_{\mu\nu}$, with the only restriction that the conformal factor $\sqrt{-g} \, f' ({\cal R}) g^{\mu \nu}$ relating $g_{\mu\nu}$ and ${\hat g}_{\mu\nu}$ be non-degenerate. In the case of the Hilbert-Einstein Lagrangian, it is $f' ({\cal R})=1$ and the statement is trivial.

The conformal transformation
\begin{equation}\tag{53}
g_{\mu\nu} \longrightarrow {\hat g}_{\mu \nu }=f' ({\cal R}) \, g_{\mu
\nu }
\end{equation}
implies that ${\cal R}_{(\mu\nu)}(\Gamma)={\cal R}_{\mu \nu}({\hat g})$.
It is useful to consider the trace of the field
equations (50)
\begin{equation} \tag{54}
f' ({\cal R}) {\cal R} -2f( {\cal R})=g^{\alpha\beta} T^{(m)}_{\alpha\beta}\equiv
T^{(m)} \,,
\end{equation}
which controls the solutions of Eq.(51). We
refer to this scalar equation as the *structural or master equation*
of space-time. *In vacuo* and in the presence of conformally
invariant matter with $T^{(m)}=0$, this scalar equation
admits constant solutions. In these cases, Palatini $f({\cal R})$ gravity reduces to General Relativity with a cosmological constant (Ferraris, Francaviglia and Volovich, 1994).

In the case of interaction with matter fields, the structural equation (53), if explicitly solvable, provides, in principle an expression ${\cal R}=f( T^{(m)})$ and, as a result, both $f({\cal R})$ and $f'({\cal R})$ can be expressed in terms of $T^{(m)}$. This fact allows one to express, at least formally, ${\cal R}$ in terms of $T^{(m)}$, which has deep consequences for the description of physical systems. Matter rules the bi-metric structure of space-time and, consequently, both the geodesic and metric structures which are intrinsically different.

## The Hamiltonian formulation and canonical quantization

**$f(R)$ gravity** can be achieved also in Hamiltonian formulation in some different approaches (Capozziello and Garattini, 2007).
For example, the Ostrogradsky treatment of higher derivative models can be considered or the conformal equivalence between **$f(R)$ gravity** and Einstein’s General Relativity coupled to a scalar field, as we will see below. However, all these formulations are related by canonical transformations as shown in (Deruelle, Sendouda, and Youssef 2009).

In particular, in the framework of the
Arnowitt-Deser-Misner ($\mathcal{ADM}$) formalism (Arnowitt, Deser and Misner, 1963), one
can investigate the possibility to find out cosmological terms
as eigenvalues of $f(R)$-Hamiltonians in a
Sturm-Liouville-like problem (Garattini, 2006). This issue is particularly relevant from several
viewpoints. First of all, the vacuum energy of
gravitational field is not a particular feature of General Relativity where the
cosmological constant has to by added *by hand* into
dynamics. At a classical level, it is well known that **$f(R)$ gravity**, for the Ricci scalar $R$ equal to a constant, exhibit
several de Sitter solutions (Barrow and Ottewill,1983 ) but a definite
discussion, at a fundamental level, considering cosmological terms
as eigenvalues of such theories is necessary. Besides, the
computation of the Casimir energy, the seeking for zero point
energy in different backgrounds (for different $f(R)$, we
expect different zero point energies and, obviously, different
vacuum states), give a track to achieve one-loop energy
regularization and renormalization for this kind of theories
(Elizalde, Odintsov, Romeo, Bytsenko, Zerbini 1994). On the other hand, these issues can be
considered in a multigravity approach to spacetime foam if the
$\mathcal{N}$ spacetimes constituting the foam are supposed to
evolve, in general, with different curvature laws and ground
states (cosmological constants) (Garattini, 2002).

Let us briefly report how to compute the Hamiltonian constraint for General Relativity
considering the standard Hilbert-Einstein theory $f\left(R\right)=R$ and
the Arnowitt-Deser-Misner ($\mathcal{ADM}$) $3+1$ decomposition (Arnowitt, Deser and Misner, 1963). In
terms of these variables, the line element is
\[
ds^{2}=g_{\mu\nu}\left( x\right) dx^{\mu}dx^{\nu}=\left( -N^{2}+N_{i}
N^{i}\right) dt^{2}+2N_{j}dtdx^{j}+g_{ij}dx^{i}dx^{j}.\tag{55}
\]
$N$ is the *lapse function*, while $N_{i}$ the *shift functions*. In this section we adopt the signature $(-,+,+,+)$ and $\kappa =8\pi c^{-4}$ to point out the role of coupling in the various contributions.
In terms of these variables, the gravitational Lagrangian, with the boundary
terms neglected, can be written as
\begin{equation}
\mathcal{L}\left[ N,N_{i},g_{ij}\right] =\sqrt{-g}R=\frac{N{}\,\sqrt{^{3}g}
}{2\kappa}\text{ }\left[ K_{ij}K^{ij}-K^{2}+\,\left( ^{3}R-2\Lambda
_{c}\right) \right],\tag{56}
\end{equation}
where $K_{ij}$ is the second fundamental form, $K=$ $g^{ij}K_{ij}$ is the
trace, $^{3}R$ is the three dimensional scalar curvature and $\sqrt{^{3}g}$ is
the three dimensional determinant of the metric. The conjugate momentum is
simply
\begin{equation}
\pi^{ij}=\frac{\delta\mathcal{L}}{\delta\left( \partial_{t}g_{ij}\right)
}=\left( ^{3}g^{ij}K-K^{ij}\text{ }\right) \frac{\sqrt{^{3}g}}{2\kappa
}.\tag{57}
\end{equation}
By a Legendre transformation, we calculate the Hamiltonian
\begin{equation}
H={\displaystyle\int}
d^{3}x\left[ N\mathcal{H+}N_{i}\mathcal{H}^{i}\right], \tag{58}
\end{equation}
where
\begin{equation}
\mathcal{H}=\left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}-\frac{\sqrt
{^{3}g}}{2\kappa}\left( ^{3}R-2\Lambda_{c}\right) \tag{59}
\end{equation}
and
\begin{equation}
\mathcal{H}^{i}=-2\nabla_{j}\pi^{ji}.\tag{60}
\end{equation}
where $\Lambda_{c}$ is the bare cosmological constant. The equations of motion
lead to two classical constraints
\begin{equation}\tag{61}
\left\{
\begin{array}
[c]{c}
\mathcal{H}=0\\
\mathcal{H}^{i}=0
\end{array}
\right.
\end{equation}
representing invariance under *time* re-parameterization and invariance
under diffeomorphism, respectively. $G_{ijkl}$ is the *supermetric*
defined as
\begin{equation}
G_{ijkl}=\frac{1}{2\sqrt{-g}}(g_{ik}g_{jl}+g_{il}g_{jk}-g_{ij}g_{kl}).\tag{62}
\end{equation}
When $\mathcal{H}$ and $\mathcal{H}^{i}$ are considered as operators acting on
some wave function, we have
\begin{equation}
\mathcal{H}\Psi\left[ g_{ij}\right] =0\tag{63}
\end{equation}
and
\begin{equation}
\mathcal{H}_{i}\Psi\left[ g_{ij}\right] =0.\tag{64}
\end{equation}
Eq. (63) is the Wheeler-DeWitt equation (WDW) (De Witt, 1967). Eqs. (63) and (64) describe the *wave function of the universe* $\Psi\left[
g_{ij}\right] $. The WDW equation represents invariance under *time*
reparameterization in an operatorial form. This standard lore can be applied
to a generic $f(R)$ theory of gravity with the aim to achieve a cosmological
term as an eigenvalue of the WDW equation.

Let us consider now the Lagrangian density describing a generic $f(R)$ model, namely \begin{eqnarray}\tag{65} \mathcal{L}=\sqrt{-g}\left( f\left( R\right) -2\Lambda_{c}\right) ,\qquad\text{with}\;f^{\prime\prime}\neq 0. \end{eqnarray} A cosmological term is added also in this case for the sake of generality. Obviously $f^{\prime\prime}=0$ corresponds to General Relativity. The generalized Hamiltonian density for the $f\left( R\right)$ theory assumes the form \begin{equation}\tag{66} \mathcal{H}=\frac{1}{2\kappa}\left[\frac{\mathcal{P}}{6}\left( {}^{\left( 3\right) }R-2\Lambda_{c}-3K_{ij}K^{ij}+K^{2}\right) +V(\mathcal{P})-\frac {1}{3}g^{ij}\mathcal{P}_{\mid ij}-2p^{ij}K_{ij}\right], \end{equation} where \begin{equation} V(\mathcal{P})=\sqrt{-g}\left[ Rf^{\prime}\left( R\right) -f\left( R\right) \right].\tag{67} \end{equation} Henceforth, the superscript $3$ indicating the spatial part of the metric will be omitted on the metric itself. When $f\left( R\right) =R$, $V(\mathcal{P} )=0$ as it should be. Since \begin{equation}\tag{68} \mathcal{P}^{ij}=-2\sqrt{-g}g^{ij}f^{\prime}\left( R\right) \qquad \Longrightarrow\qquad\mathcal{P=}-6\sqrt{-g}f^{\prime}\left( R\right), \end{equation} we have \begin{equation} \mathcal{H}=\frac{1}{2\kappa}\left[ -\sqrt{-g}f^{\prime}\left( R\right) \left( {}^{\left( 3\right) }R-2\Lambda_{c}-3K_{ij}K^{ij}+K^{2}\right) +V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}-2p^{ij}K_{ij}\right].\tag{69} \end{equation} With the help of Eq. (57), Eq.(83) becomes \begin{equation}\tag{70} \mathcal{H}=f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +\frac{1}{2\kappa}\left[ \sqrt {-g}f^{\prime}\left( R\right) \left( 2K_{ij}K^{ij}\right) +V(\mathcal{P} )+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right)_{\mid ij}-2p^{ij}K_{ij}\right]. \end{equation} However \begin{equation}\tag{71} p^{ij}=\sqrt{-g}K^{ij}, \end{equation} then we obtain \begin{equation}\tag{72} \mathcal{H}=f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +\frac{1}{2\kappa}\left[ 2\sqrt {-g}K_{ij}K^{ij}\left( f^{\prime}\left( R\right) -1\right) +V(\mathcal{P} )+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right)\right)_{\mid ij}\right] \end{equation} and transforming into canonical momenta, one gets \begin{equation} \mathcal{H}=f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +2\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4}^{2}\right] \left( f^{\prime}\left( R\right) -1\right) +\frac{1}{2\kappa}\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right)_{\mid ij}\right].\tag{73} \end{equation} By imposing the Hamiltonian constraint, we obtain \begin{eqnarray}\tag{74} f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right)}R\right]+2\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \left( f^{\prime}\left( R\right) -1\right) +\frac{1}{2\kappa }\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}\right] =-f^{\prime}\left( R\right) \sqrt{-g} \frac{\Lambda_{c}}{\kappa} \end{eqnarray} If we assume that $f^{\prime}\left( R\right) \neq0$ the previous expression becomes \begin{equation}\tag{75} \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g} }{2\kappa}{}^{\left( 3\right) }R\right] +\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4}^{2}\right] \frac{2\left( f^{\prime }\left( R\right) -1\right) }{f^{\prime}\left( R\right) }+\frac{1}{2\kappa f^{\prime}\left( R\right) }\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt {-g}f^{\prime}\left( R\right) \right) _{\mid ij}\right]=-\sqrt{-g} \frac{\Lambda_{c}}{\kappa}{}. \end{equation} Now, we integrate over the hypersurface $\Sigma$ to obtain \begin{eqnarray}\tag{76} \int_{\Sigma}d^{3}x\left\{\left[\left(2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right)}R\right] +\left( 2\kappa\right)\left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \frac{2\left(f^{\prime}\left( R\right)-1\right)} {f^{\prime}\left( R\right)}\right\} +\int_{\Sigma}d^{3}x\frac{1}{2\kappa f^{\prime}\left( R\right) }\left[ V(\mathcal{P})+2g^{ij}\left(\sqrt{-g}f^{\prime}\left( R\right)\right) _{\mid ij}\right] =-\frac{\Lambda_{c}}{\kappa}\int_{\Sigma}d^{3}x\sqrt{-g}. \end{eqnarray} The term \begin{equation}\tag{77} \frac{1}{\kappa}\int_{\Sigma}d^{3}x\frac{1}{f^{\prime}\left( R\right) }g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right)\right)_{\mid ij} \end{equation} appears to be a three-divergence and therefore will not contribute to the computation. The remaining equation simplifies into \begin{equation} \int_{\Sigma}d^{3}x\left\{ \left[ \left( 2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right) }R\right] +\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \frac{2\left( f^{\prime}\left( R\right) -1\right) } {f^{\prime}\left( R\right) }+\frac{V(\mathcal{P})}{2\kappa f^{\prime}\left( R\right) }\right\} =-\frac{\Lambda_{c}}{\kappa}\int_{\Sigma}d^{3}x\sqrt {g}.\tag{78} \end{equation} By a canonical procedure of quantization, one can obtain the vacuum state of a generic $f(R)$ theory.

Let us consider now the Lagrangian density describing a generic $f(R)$ theory of gravity, namely \begin{equation} \mathcal{L}=\sqrt{-g}\left( f\left( R\right) -2\Lambda_{c}\right) ,\qquad\text{with}\;f^{\prime\prime}\neq0.\tag{79} \end{equation} A cosmological term is added also in this case for the sake of generality. The generalized Hamiltonian density for the $f\left(R\right)$ theory assumes the form \begin{equation}\tag{80} \mathcal{H}=\frac{1}{2\kappa}\left[ \frac{\mathcal{P}}{6}\left( {}^{\left( 3\right) }R-2\Lambda_{c}-3K_{ij}K^{ij}+K^{2}\right) +V(\mathcal{P})-\frac {1}{3}g^{ij}\mathcal{P}_{\mid ij}-2p^{ij}K_{ij}\right], \end{equation} where \begin{equation} V(\mathcal{P})=\sqrt{-g}\left[ Rf^{\prime}\left( R\right) -f\left( R\right) \right].\tag{81} \end{equation} Henceforth, the superscript $3$ indicating the spatial part of the metric will be omitted on the metric itself. When $f\left( R\right) =R$, $V(\mathcal{P} )=0$ as it should be. Since \begin{equation}\tag{82} \mathcal{P}^{ij}=-2\sqrt{-g}g^{ij}f^{\prime}\left( R\right) \qquad \Longrightarrow\qquad\mathcal{P=}-6\sqrt{-g}f^{\prime}\left( R\right), \end{equation} we have \begin{equation} \mathcal{H}=\frac{1}{2\kappa}\left[ -\sqrt{-g}f^{\prime}\left( R\right) \left( {}^{\left( 3\right) }R-2\Lambda_{c}-3K_{ij}K^{ij}+K^{2}\right) +V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}-2p^{ij}K_{ij}\right].\tag{83} \end{equation} With the help of Eq. (57), Eq. (83) becomes \begin{equation}\tag{84} \mathcal{H}=f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +\frac{1}{2\kappa}\left[ \sqrt {g}f^{\prime}\left( R\right) \left( 2K_{ij}K^{ij}\right) +V(\mathcal{P} )+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}-2p^{ij}K_{ij}\right]. \end{equation} However \begin{equation}\tag{85} p^{ij}=\sqrt{-g}K^{ij}, \end{equation} then we obtain \begin{equation}\tag{86} \mathcal{H}=f^{\prime}\left( R\right) \left[ \left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +\frac{1}{2\kappa}\left[ 2\sqrt {g}K_{ij}K^{ij}\left( f^{\prime}\left( R\right) -1\right) +V(\mathcal{P} )+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right)\right) _{\mid ij}\right] \end{equation} and transforming into canonical momenta, one gets \begin{eqnarray} \mathcal{H}=f^{\prime}\left( R\right) \left[\left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}\left( ^{\left( 3\right) }R-2\Lambda_{c}\right) \right] +2\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4}^{2}\right] \left( f^{\prime}\left( R\right) -1\right) +\frac{1}{2\kappa}\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}\right].\tag{87} \end{eqnarray} By imposing the Hamiltonian constraint, we obtain \begin{eqnarray}\tag{88} f^{\prime}\left(R\right)\left[\left(2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right) }R\right]+2\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \left( f^{\prime}\left( R\right)-1\right)+\frac{1}{2\kappa }\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left(R\right) \right) _{\mid ij}\right] =-f^{\prime}\left( R\right)\sqrt{-g} \frac{\Lambda_{c}}{\kappa} \end{eqnarray} If we assume that $f^{\prime}\left( R\right)\neq0$ the previous expression becomes \begin{equation}\tag{89} \left[\left( 2\kappa\right) G_{ijkl}\pi^{ij}\pi^{kl}{}-\frac{\sqrt{-g} }{2\kappa}{}^{\left( 3\right) }R\right] +\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4}^{2}\right] \frac{2\left( f^{\prime }\left( R\right) -1\right) }{f^{\prime}\left( R\right) }+\frac{1}{2\kappa f^{\prime}\left( R\right) }\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt {g}f^{\prime}\left( R\right)\right)_{\mid ij}\right] =-\sqrt{-g} \frac{\Lambda_{c}}{\kappa}{}. \end{equation} Now, we integrate over the hypersurface $\Sigma$ to obtain \begin{eqnarray}\tag{90} \int_{\Sigma}d^{3}x\left\{ \left[ \left( 2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right) }R\right] +\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \frac{2\left( f^{\prime}\left( R\right) -1\right)} {f^{\prime}\left( R\right)}\right\}+\int_{\Sigma}d^{3}x\frac{1}{2\kappa f^{\prime}\left( R\right) }\left[ V(\mathcal{P})+2g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij}\right]=-\frac{\Lambda_{c}}{\kappa}\int_{\Sigma}d^{3}x\sqrt{-g}. \end{eqnarray} The term \begin{equation}\tag{91} \frac{1}{\kappa}\int_{\Sigma}d^{3}x\frac{1}{f^{\prime}\left( R\right) }g^{ij}\left( \sqrt{-g}f^{\prime}\left( R\right) \right) _{\mid ij} \end{equation} appears to be a three-divergence and therefore will not contribute to the computation. The remaining equation simplifies into \begin{equation} \int_{\Sigma}d^{3}x\left\{\left[\left( 2\kappa\right) G_{ijkl}\pi ^{ij}\pi^{kl}{}-\frac{\sqrt{-g}}{2\kappa}{}^{\left( 3\right) }R\right] +\left( 2\kappa\right) \left[ G_{ijkl}\pi^{ij}\pi^{kl}+\frac{\pi}{4} ^{2}\right] \frac{2\left( f^{\prime}\left( R\right) -1\right)} {f^{\prime}\left( R\right) }+\frac{V(\mathcal{P})}{2\kappa f^{\prime}\left( R\right) }\right\} =-\frac{\Lambda_{c}}{\kappa}\int_{\Sigma}d^{3}x\sqrt {-g}.\tag{92} \end{equation} Starting from these considerations, by a canonical procedure of quantization, one obtains the vacuum state of a generic $f(R)$ theory (Capozziello and Garattini, 2007). In this context, the cosmological constant (vacuum state) emerges as a Wheeler-De Witt eigenvalue. The related wave functional can be split by an orthogonal decomposition and then, constructing the transverse-traceless propagator, it is possible to obtain, after a variational minimization, the total one-loop energy density for the TT tensors. Such a quantity explicitly depends on the form of $f(R)$. The one-loop energy regularization and renormalization are achieved by the zeta function regularization method. The resulting renormalized $\Lambda_{0}$ is a running constant which can be set to zero depending on the value of an arbitrary mass scale parameter $\mu$.

## Conformal transformations

Conformal transformations are a useful tool to deal with **Extended Theories of Gravity**.
Let the pair $\{{\cal M}, g_{\mu\nu}\}$ be a space-time, with ${\cal M}$ a smooth
manifold of dimension $n \geq 2$ and $g_{\mu\nu} $ a Lorentzian metric on ${\cal M}$. The point-dependent rescaling of the
metric tensor is
\begin{eqnarray} \tag{93}
g_{\mu\nu} \longrightarrow \tilde{g}_{\mu\nu}=\Omega^2
g_{\mu\nu} \,,
\end{eqnarray}
where the *conformal factor*
$\Omega(x) $ is a nowhere
vanishing, regular
function, is called a *Weyl* or *conformal*
transformation (see (Bronnikov, 2001) (Bronnikov, 2002) (Bronnikov and Shikin, 2002)). Due to this metric rescaling, the
lengths of spacelike and timelike intervals and the norms of
spacelike and timelike vectors are changed, while null vectors
and null intervals of the metric $g_{\mu\nu}$ remain
null in the rescaled metric $\tilde{g}_{\mu\nu}$. The light cones
are left unchanged by the transformation (93)
and the space-times $\{{\cal M}, g_{\mu\nu}\}$ and $\{{\cal M},
\tilde{g}_{\mu\nu}\} $
exhibit the same causal structure (Wald, 1984). A vector that is timelike, spacelike,
or null with respect to the metric $g_{\mu\nu}$ has the same
character with respect to $\tilde{g}_{\mu\nu}$, and *vice-versa*.

The transformation properties of geometrical quantities are (Synge, 1955) (Wald, 1984).

\begin{eqnarray} \tag{94}
\tilde{g}^{\mu\nu}=\Omega^{-2}\, g^{\mu\nu} \,,
\;\;\;\;\;\;\;\;\;\;
\tilde{g}=\Omega^{2n} \, g \,,
\end{eqnarray}
for the inverse metric and the metric determinant,
\begin{eqnarray} \tag{95}
\tilde{\Gamma}^{\alpha}_{\beta \gamma}=
\Gamma^{\alpha}_{\beta\gamma}+\Omega^{-1}\left(
\delta^{\alpha}_{\beta} \nabla_{\gamma}\Omega +
\delta^{\alpha}_{\gamma} \nabla_{\beta} \Omega
-g_{\beta\gamma}\nabla^{\alpha} \Omega \right) \,,
\end{eqnarray}
for the Christoffel symbols,
\begin{eqnarray}
\tilde{R}_{\alpha\beta\gamma}^{\delta} &=&
{R_{\alpha \beta \gamma}}^{\delta}+2 \,
\delta^{\delta}_{[\alpha} \nabla_{\beta]}\nabla_{\gamma} ( \ln
\Omega )
-2g^{\delta \sigma} g_{\gamma [ \alpha}\nabla_{\beta
]}\nabla_{\sigma}
( \ln \Omega )+ 2 \nabla_{[ \alpha} ( \ln \Omega ) \, \delta^{\delta}_{\beta
]}
\nabla_{\gamma}( \ln \Omega )+\nonumber\\&&
-2\nabla_{[ \alpha}( \ln \Omega ) \, g_{\beta ]
\gamma} \, g^{\delta \sigma} \, \nabla_{\sigma}( \ln \Omega )
- 2g_{\gamma [ \alpha} \delta^{\delta}_{\beta ]} \,
g^{\sigma \rho }
\nabla_{\sigma} ( \ln \Omega ) \, \nabla_{\rho} ( \ln \Omega ) \,,
\tag{96}
\end{eqnarray}
for the Riemann tensor,
\begin{eqnarray}
{\tilde R}_{\alpha\beta } & =& R_{\alpha\beta }
-(n-2) \nabla_{\alpha}\nabla_{\beta }
( \ln \Omega )
-g_{\alpha \beta } g^{\rho\sigma } \, \nabla_{\sigma}
\nabla_{\rho} ( \ln \Omega )+\left( n-2 \right) \nabla_{\alpha} ( \ln \Omega )
\nabla_{\beta}( \ln \Omega )
-\left( n-2 \right) g_{\alpha\beta }\, g^{\rho\sigma} \,
\nabla_{\rho}( \ln \Omega ) \nabla_{\sigma}( \ln \Omega ) \,,
\tag{97}
\end{eqnarray}
for the Ricci tensor, and
\begin{eqnarray}
\tilde{R} & \equiv & \tilde{g}^{\alpha\beta}
\tilde{R}_{\alpha\beta }=
\Omega^{-2} \left[ R-2 \left( n-1 \right) \Box \left( \ln \Omega
\right) - \left( n-1 \right) \left( n-2 \right) \,
\frac{g^{\alpha\beta} \nabla_{\alpha} \Omega \, \nabla_{\beta}
\Omega}{\Omega^2}
\right] \,, \tag{98}
\end{eqnarray}
for the Ricci scalar. In the case of $n=4$ space-time dimensions,
the transformation
property of the Ricci scalar can be written as
\begin{eqnarray}
\tilde{R} &=& \Omega^{-2} \left[ R-\frac{6 \Box \Omega}{\Omega}
\right] = \Omega^{-2} \left[ R-\frac{12 \Box (
\sqrt{\Omega})}{\sqrt{\Omega}}
+ \frac{3g^{\alpha\beta} \nabla_{\alpha} \Omega \nabla_{\beta}
\Omega}
{\Omega^2}\right]\,. \tag{99}
\end{eqnarray}
The Weyl tensor
${C_{\alpha\beta\gamma}}^{\delta}$, with the last index
controvariant, is conformally
invariant
\begin{eqnarray} \tag{100}
{\tilde C}_{\alpha\beta\gamma}^{\delta}=
{C_{\alpha\beta\gamma}}^{\delta}\,,
\end{eqnarray}
but the same tensor with indices raised or lowered with respect to
${C_{\alpha\beta\gamma}}^{\delta}$ is not. This property
explains the name *conformal tensor* used for
${C_{\alpha\beta\gamma}}^{\delta}$ (Lorentz, 1937).
If the original metric $g_{\alpha\beta} $ is Ricci-flat *i.e.*, $R_{\alpha\beta}=0$), the
conformally transformed metric $\tilde{g}_{\alpha\beta}$ is not given by
Eq.(97). In the conformally transformed
world, the conformal
factor $\Omega$ plays the role of an
effective
matter and
this fact has consequences for the physical interpretation of
the theory. A vacuum metric in the Jordan frame is not such in
the Einstein frame, and the interpretation
of what is matter and
what is gravity becomes frame-dependent (Sotiriou, Liberati, and Faraoni, 2008). However, if the
Weyl tensor
vanishes in one frame, it also vanishes
in the conformally related frame. Conformally flat metrics are
mapped into conformally flat metrics.
Since, in general, tensorial quantities are not
invariant under conformal transformations, neither are the tensorial
equations describing geometry and physics. An equation involving
a tensor field $\psi$ is said to be *conformally invariant*
if there exists a number $w$ (the *conformal weight* of
$\psi$) such that, if $\psi$ is a solution of a tensor equation
with the metric $g_{\mu\nu}$ and the associated geometrical
quantities, $\tilde{\psi} \equiv \Omega^w \psi$ is a solution of
the corresponding equation with the metric $\tilde{g}_{\mu\nu}$
and the associated geometry.
In addition to geometric quantities, one needs to consider
the behavior of standard
forms of matter under conformal transformations. Most forms of matter and fields are not conformally
invariant: invariance under
conformal transformations is a very
special property.
In general, the covariant conservation equation for a
(symmetric) stress-energy tensor $T_{\alpha\beta}^{(m)} $
representing
ordinary matter,
\begin{eqnarray} \tag{101}
\nabla^{\beta} \, T_{\alpha\beta}^{(m)} =0\,,
\end{eqnarray}
is not conformally invariant (Wald, 1984). The conformally
transformed $\tilde{T}_{\alpha\beta}^{(m)} $ satisfies the
equation
\begin{eqnarray} \tag{102}
\tilde{\nabla}^{\beta} \, \tilde{T}_{\alpha\beta}^{(m)} = -
\tilde{T}^{(m)}
\,\tilde{\nabla}_{\alpha} \left( \ln \Omega \right)\,.
\end{eqnarray}
Clearly, the conservation equation (101) is conformally
invariant only for a matter
component that has vanishing trace
$T^{(m)}$ of the energy-momentum tensor. This feature is
associated with light-like behavior; examples are the
electromagnetic field and a radiative fluid with equation of
state $p^{(m)}=\rho^{(m)}/3$.
Unless $T^{(m)}=0$, Eq.(102) describes an exchange of
energy and momentum between matter and the scalar field $\Omega$,
reflecting
the fact that matter and the geometric factor $\Omega$ are
directly coupled in the Einstein frame
description.
Since the geodesic equation ruling the motion of free particles
in General Relativity can be derived from the conservation
equation
(101) (*geodesic hypothesis*), it
follows that timelike geodesics of the original metric
$g_{\alpha\beta}$
are not geodesics of the rescaled metric $\tilde{g}_{\alpha\beta}$ and *vice-versa*. Particles in free fall in the world
$\left\{{\cal M}, g_{\alpha\beta} \right\}$ are subject to a force proportional to the
gradient $\tilde{\nabla}^{\alpha} \Omega$ in the rescaled world
$\left\{ {\cal M}, \tilde{g}_{\alpha\beta} \right\}$ (this is often
identified as a *fifth force* acting on all massive particles
and, therefore, it can be said that no massive test particles exist in the
Einstein frame). The stress-energy tensor definition in terms
of the matter action ${\cal A}^{(m)}=\int d^4x\, \sqrt{-g} \, {\cal L}^{(m)} $ is
\begin{eqnarray}
\tag{103}
\tilde{T}_{\alpha\beta}^{(m)}=\frac{-2}{\sqrt{ -\tilde{g} } } \,
\frac{ \delta \left( \sqrt{-\tilde{g}} \, \, {\cal L}^{(m)} \right) }{\delta
\tilde{g}^{\alpha\beta} } \,,
\end{eqnarray}
together with the rescaling (93) of the metric, yields
\begin{equation}\tag{104}
\tilde{T}_{\alpha\beta}^{(m)}= \Omega^{-2}\, T_{\alpha\beta}^{(m)}\,,\quad {\tilde{T}_{\alpha}^{\beta}}^{(m)} =\Omega^{-4}\,{T_{\alpha}^{\beta}}^{(m)}\,,\quad
{\tilde{T}^{\alpha\beta}}^{(m)}=\Omega^{-6}\, {T^{\alpha\beta}}^{(m)}\,,\quad \tilde{T}^{(m)}=\Omega^{-4}\, T^{(m)}\,.
\end{equation}
The last relation makes it clear that the trace vanishes
in the Einstein frame if and only if it
vanishes in the Jordan frame.
Conformal invariance
corresponds to the absence of a characteristic length (or mass)
scale in the physics. In general, the effective potential $V(\phi)$, coming from conformal transformations,
contains dimensional parameters (such as a mass $m$, that is a further *characteristic gravitational length*).
This means that the further degrees of freedom coming from **Extended Theories of Gravity** give rise to features that could play a fundamental role in the dynamics of astrophysical structures addressing, in principle, the dark matter problem.

## $f(R)$ cosmology

Cosmology is one of the straightforward applications that motivates the introduction of **$f(R)$ gravity**. As we will discuss below, **$f(R)$ gravity** could constitute a geometrical explanation for the dark side of the universe (Capozziello, 2002). The main ingredient of such an interpretation lies on the fact that the further degrees of freedom of **$f(R)$ gravity** can be assembled into an effective curvature stress-energy tensor (Capozziello, De Laurentis and Lambiase, 2012) which could give rise to dark energy effects.

### Cosmological Equations

Let us start from the action (2) where we add also matter. We have \begin{equation}\tag{105} {\cal A}=\int d^4x \sqrt{-g} \left[f(R)+{\cal L}_{(m)} \right]\,{.} \end{equation}

The field are

\begin{equation}\tag{106} f'(R)R_{\alpha\beta}-\frac{1}{2}f(R)g_{\alpha\beta}= f'(R)^{;\alpha\beta}(g_{\alpha\mu}g_{\beta\nu}-g_{\alpha\beta}g_{\mu\nu})+ \tilde{T}^{(m)}_{\alpha\beta}\,, \end{equation} which can be recast in the more expressive form \begin{equation}\tag{107} G_{\alpha\beta}=R_{\alpha\beta}-\frac{1}{2}g_{\alpha\beta}R=T^{(curv)}_{\alpha\beta}+T^{(m)}_{\alpha\beta}\,, \end{equation} where \begin{equation} \tag{108} T^{(curv)}_{\alpha\beta}=\frac{1}{f'(R)}\left\{\frac{1}{2}g_{\alpha\beta}\left[f(R)-Rf'(R)\right]+ f'(R)^{;\alpha\beta}(g_{\alpha\mu}g_{\beta\nu}-g_{\alpha\beta}g_{\mu\nu}) \right\} \end{equation} and \begin{equation}\tag{109} T^{(m)}_{\alpha\beta}=\frac{1}{f'(R)}\tilde{T}^{(m)}_{\alpha\beta}\,, \end{equation} is the stress-energy tensor of matter where we have taken into account the nontrivial coupling to geometry.

Reducing the action to a point-like, Friedmann-Robertson-Walker one, we can write for the geometric part \begin{equation}\tag{110} {\cal A}_{(curv)}=\int dt {\cal L}(a, \dot{a}, R, \dot{R})\,{,} \end{equation} where dot means derivative with respect to the cosmic time. The scale factor $a$ and the Ricci scalar $R$ can be considered as canonical variables. This position could seem arbitrary since $R$ depends on $a, \dot{a}, \ddot{a}$, but it is generally used in canonical quantization. The definition of $R$ in terms of $a, \dot{a}, \ddot{a}$ introduces a constraint which eliminates second and higher order derivatives in action (110), and gives a system of second order differential equations in $\{a, R\}$. Action (110) can be written as \begin{equation}\tag{111} {\cal A}_{(curv)}=2\pi^2\int dt \left\{ a^3f(R)-\lambda\left [ R+6\left ( \frac{\ddot{a}}{a}+\frac{\dot{a}^2}{a^2}+\frac{k}{a^2}\right)\right]\right\}\,{,} \end{equation} where the Lagrange multiplier $\lambda$ is derived by varying with respect to $R$. It is \begin{equation}\tag{112} \lambda=a^3f'(R)\,{.} \end{equation} The point-like Lagrangian is then \begin{equation} \tag{113} {\cal L}={\cal L}_{(curv)}+{\cal L}_{(m)}=a^3\left[f(R)-R f'(R)\right]+6a\dot{a}^2f'(R)+6a^2\dot{a}\dot{R}f''(R)-6ka f'(R)+a^3p_{(m)}\,, \end{equation} where we have taken into account also the fluid matter contribution which is a pressure term.

The Euler-Lagrange equations are \begin{equation} \tag{114} 2\left(\frac{\ddot{a}}{a}\right)+\left(\frac{\dot{a}}{a}\right)^2+ \frac{k}{a^2}=-p_{(tot)}\,, \end{equation} and \begin{equation} \tag{115} f''(R)\left[R+6\left(\frac{\ddot{a}}{a}+\frac{\dot{a}^2}{a}^2+\frac{k}{a^2}\right)\right]=0\,. \end{equation} The dynamical system is completed by the energy condition \begin{equation} \tag{116} \left(\frac{\dot{a}}{a}\right)^2+\frac{k}{a^2}=\frac{1}{3}\rho_{(tot)}\,. \end{equation}

### Curvature Quintessence

Combining Eq.(114) and Eq.(116), we obtain the Friedmann equation \begin{equation} \tag{117} \left(\frac{\ddot{a}}{a}\right)=-\frac{1}{6}\left[\rho_{(tot)}+3p_{(tot)} \right]\,, \end{equation} where it is clear that the accelerated or decelerated behaviour depends on the r.h.s.. It is \begin{equation} \tag{118} p_{(tot)}=p_{(curv)}+p_{(m)}\;\;\;\;\;\rho_{(tot)}=\rho_{(curv)}+\rho_{(m)}\,, \end{equation} where we have put in evidence the curvature and matter contributions.

From the curvature-stress-energy tensor, we can define a curvature pressure \begin{equation} \tag{119} p_{(curv)}=\frac{1}{f'(R)}\left\{2\left(\frac{\dot{a}}{a}\right)\dot{R}f''(R)+\ddot{R}f''(R)+\dot{R}^2f'''(R) -\frac{1}{2}\left[f(R)-Rf'(R)\right] \right\}\,, \end{equation} and a curvature density \begin{equation} \tag{120} \rho_{(curv)}=\frac{1}{f'(R)}\left\{\frac{1}{2}\left[f(R)-Rf'(R)\right] -3\left(\frac{\dot{a}}{a}\right)\dot{R}f''(R) \right\}\,. \end{equation}

From Eq. (117), the accelerated behaviour is achieved if
\begin{equation}
\tag{121} \rho_{(tot)}+ 3p_{(tot)}< 0\,,
\end{equation}
which means
\begin{equation}
\tag{122} \rho_{(curv)}> \frac{1}{3}\rho_{(tot)}\,,
\end{equation}
assuming that all matter components have non-negative pressure.
In other words, conditions to obtain acceleration depends on the
relation
\begin{equation}
\tag{123}
\rho_{(curv)}+3p_{(curv)}=\frac{3}{f'(R)}\left\{\dot{R}^2f'''(R)+\left(\frac{\dot{a}}{a}\right)\dot{R}f''(R)
+\ddot{R}f''(R)-\frac{1}{3}\left[f(R)-Rf'(R)\right]\right\}\,,
\end{equation}
which has to be compared with matter contribution. However, it
has to be
\begin{equation}
\tag{124}
\frac{p_{(curv)}}{\rho_{(curv)}}=w_{(curv)}\,,
\qquad -1\leq w_{(curv)}<0\,.
\end{equation}
The form of $f(R)$ is the main ingredient to obtain **curvature**
quintessence**.**

### Examples of exact solutions

As simple choice in order to fit the above prescriptions, we ask for solutions of the form \begin{equation} \tag{125} f(R)=f_0 R^n\,,\qquad a(t)=a_0\left(\frac{t}{t_0}\right)^{\beta}\,. \end{equation} The interesting cases are for $n\neq 1$ and $\beta\geq 1$ (accelerated behaviour). Inserting Eqs. (125) into the above dynamical system, we obtain the exact solutions \begin{equation} \tag{126} \beta=2\,;\qquad n=-1,\frac{3}{2}\,;\qquad k=0\,. \end{equation} In both cases, the deceleration parameter is \begin{equation} \tag{127} q_0=-\frac{1}{2}\,, \end{equation} in agreement with the cosmic accelerated behavior.

The case $n=3/2$ deserves a further discussion. Considering conformal transformation from Jordan frame to Einstein frame, it is possible to give explicit form for the scalar field potential. It is \begin{equation}\tag{128} \tilde{g}_{\alpha\beta}\equiv f'(R)g_{\alpha\beta}\,{,}\qquad \varphi=\sqrt{\frac{3}{2}}\ln f'(R)\,{.} \end{equation} The conformal equivalence of the Lagrangians gives \begin{equation}\tag{129} {\cal L}=\sqrt{-g}\,f_0R^{3/2}\longleftrightarrow \tilde{\cal L}=\sqrt{-\tilde{g}}\left[-\frac{\tilde{R}}{2}+ \frac{1}{2}\nabla_{\mu}\varphi\nabla^{\mu}\varphi-V_0\exp\left( \sqrt{\frac{2}{3}}\varphi\right)\right]\,{,} \end{equation} in our physical units. This is the so-called Liouville field theory and it is one of the few cases where a fourth-order Lagrangian can be expressed, in the Einstein frame, in terms of elementary functions under a conformal transformation. It is possible to obtain the general cosmological solution

\begin{equation}\tag{130} a(t)=a_0\sqrt{c_4t^4+c_3t^3+c_2t^2+c_1t+c_0}\,. \end{equation}

The constants $c_i$ are combinations of the initial conditions.
Their values determine the cosmological evolution. For
example, $c_4\neq 0$ gives a power law inflation while, if the
regime is dominated by the linear term in $c_1$, we get a
radiation-dominated stage.
The above one is just a simple model showing how the dark energy issue can be recovered into the framework of **$f(R)$ gravity**. However more realistic models can be worked out as reported in literature (Capozziello and Francaviglia,2008) (Nojiri and Odintsov,2007)(Starobinsky,1979). In particular, in (Nojiri and
Odintsov, 2003), it has been proposed the possibility of unifying the
early-time inflation with late-time acceleration including quintessence.
A detailed review in this topic is given in (Nojiri
and Odintsov, 2011).
It is also worth noticing that, as for other dark energy
models, **$f(R)$ gravity** may show up the finite-time future singular
behaviour of one of the known four types as classified in
(Nojiri, Odintsov and Tsujikawa, 2005).
Several works discussed the finite-time future singularities in F(R)
gravity (for review, see (Nojiri and Odintsov,
2011)).
Furthermore (Bamba, Nojiri and Odintsov, 2008)
have demonstrated
that the removal of future singularities in **$f(R)$ gravity** is easily
achieved by the addition of $R^2$ term. Simultaneously with such removal,
this term straightforwardly leads to the unification of inflation with
dark energy.
In other words, a universe model which consistently contains inflation
and dark energy, as predicted by **$f(R)$ gravity**, does not present future
singularities.

### The phase-space view of $f(R)$ gravity

The phase space portrait of cosmological models, deduced from **$f(R)$ gravity**, can be discussed by analytical and numerical methods (see (Capozziello, Occhionero and Amendola, 1993) (Carloni, Dunsby, Capozziello and Troisi,2005 ) (Carloni and Dunsby, 2009) (Abdelwahab, Goswami and Dunsby,2012 ). According to this approach, peculiar structures, such as attractors and singular points, emerge.
Trajectories of interest are those undergoing an inflationary/dark energy expansion and then reaching Friedmannian asymptotic stable phases. These features can be moreover discussed through effective potentials deduced starting from the form of **$f(R)$ gravity** (see in particular (Capozziello, Occhionero and Amendola, 1993)). Different kinds of potential regions of phase space can be recognized distinguishing various physical behaviors. There are, in general, allowed regions where trajectories reach the Friedman phases after undergoing an inflationary period, disconnected regions where trajectories, although physical, never reach the Friedman stages and forbidden regions where there are no physical solutions. A survey of the global phase space can be given via Poincaré projections where suitable variables are adopted. This picture allows to represent singularities and attractors at finite and at infinite.
The phase-space view can be considered one of the most powerful approaches to deal with dynamical systems in cosmology in order to extract physical information from solutions (see *e.g.* (Wainwright and Ellis, 1997)).

## Cosmological perturbations

Cosmological perturbations can be considered also in the framework of **$f(R)$ gravity**. The evolution of the scalar perturbations since the radiation epoch until the present time can be compared with the evolution of perturbations in the concordance model, *i.e.*, the $\Lambda$CDM model. This method allows to compare the effects due to $f(R)$ modifications in the matter power spectrum as measured today. Specifically, while cosmography is an approach to test the background universe, the matter power spectrum allows to test the clumpy universe. The matter power spectrum provides information on clustering processes coming from the early radiation up to the current epoch. Furthermore, it is a tool to discriminate between General Relativity and modified gravity at the background level. Here, we will briefly sketch the scalar perturbations for metric **$f(R)$ gravity** giving the main quantities that have to be considered in the perturbation analysis.

### Scalar perturbations

The scalar part of cosmological perturbations can be written in the so called Newtonian gauge assuming the conformal invariance, that is
(Mukhanov, Feldeman and Branderberger, 1992) (Dodelson, 2003)
\begin{eqnarray}\tag{131}
ds^2 = a^2\left[\left(1+2\Phi\right)d\eta^2 - \left(1-2\Psi\right)\delta_{ij}dx^idx^j\right],
\end{eqnarray}
where $\Phi\left(\eta,x^i\right)$ and $\Psi\left(\eta,x^i\right)$ are the gauge invariant Bardeen potentials (Bardeen, 1980) and $\eta$ the conformal time. At first order, the energy-momentum tensor consists of the quantities
(Mukhanov, Feldeman and Branderberger, 1992) (Dodelson, 2003)
\begin{eqnarray}\tag{132}
T^0_0 &= (\rho+\delta\rho),
\qquad
T^0_i &= -(\rho+p)v_i,
\qquad
T^i_j &=-(p+\delta p)\delta^i_j,
\end{eqnarray}
where $\delta\rho$ is the energy density perturbation, $\delta p$ is the pressure perturbation and $v$ is the peculiar velocity. In the simplest approach, we discard the effects of anisotropic stress terms. In the metric **$f(R)$ gravity**, the $(i,j)$ components relate the potentials $\Phi$ and $\Psi$ with the perturbation of $f'(R)$ (Pogosian and Silvestri, 2007)
\begin{align}\tag{133}
\left(\Phi-\Psi\right)f'(R) = \delta f'(R) = f''(R)\delta R.
\end{align}
It is important to stress that the equality between the two potentials, valid in General Relativity, is no longer true in **$f(R)$ gravity**; the derivative $ f'(R)$ can be interpreted as a new degree of freedom in the background evolution. This means that $\delta f'(R)$ is a new degree of freedom for perturbations. Let us define the new variables $\Psi^+$ and $\Xi$,
\begin{eqnarray}\tag{134}
\Psi^+ &\equiv \frac{\Phi+\Psi}{2}, \qquad \Xi &\equiv \frac{\delta f'(R)}{ f'(R)} = \Phi-\Psi.
\end{eqnarray}
$\Psi^+$ corresponds to $\Phi^+$ in Ref. (Pogosian and Silvestri, 2007), while $\Xi$ corresponds to the variable $\chi$ of the same reference divided by $ f'(R)$. We choose the variable $\Xi$ instead of $\chi$, because it represents the difference between the two potentials and allows us to control the integration when the theory deviates from General Relativity ($ f'(R)\gg 1$). Using these variables, we can write the $(0,0)$ and $(i,0)$ components of the perturbed Einstein equations as (de la Cruz-Dombriz, Dobado and Maroto, 2008)
\begin{eqnarray}
\tag{135}
&&\left[1+\frac{1}{2}\frac{\left(f'(R)\right)_N}{ f'(R)}\right]\left(\Psi^+\right)_N
+\left[1+\frac{\left( f'(R)\right)_N}{ f'(R)} + \frac{k^2}{3\mathcal{H}^2}\right]\Psi^+
+\frac{1}{4}\frac{\left( f'(R)\right)_N}{ f'(R)}\left(\Xi\right)_N
-\frac{1}{2}\left[1 - \mathcal{H}_N +2\frac{\left( f'(R)\right)_N}{ f'(R)}\right]\Xi=-\frac{a^2\kappa^2}{6\mathcal{H}^2 f'(R)}\delta\rho,
\end{eqnarray}
\begin{eqnarray}\tag{136}
\left(\Psi^+\right)_N
+ \left[1+\frac{1}{2}\frac{\left( f'(R)\right)_N}{ f'(R)}\right]\Psi^+
-\frac{3}{4}\frac{\left( f'(R)\right)_N}{ f'(R)}\Xi
&= - \frac{a^2\kappa^2\left(\rho+p\right)}{2 f'(R)}\frac{v}{\mathcal{H}}.
\end{eqnarray}
The $N$-subscript indicates the derivative with respect to $N\equiv\log(a)$. ${\mathcal{H}}$ is the Hubble parameter written in the conformal time $\eta$. The above equations can be combined to obtain the evolution equations for $\Psi^+ $ and $\Xi$, that is
\begin{align}
\tag{137}
\left(\Psi^+ \right)_N &= -\Psi^+ - \frac{1}{4}\frac{\left( f'(R)\right)_N}{ f'(R)}\left(2\Psi^+ - 3\Xi\right) - \frac{a^2\kappa^2(\rho+p)}{2 f'(R)}\frac{v}{\mathcal H},
\\
\tag{138}
\left(\Xi\right) _N &= \Xi + \frac{a^2\kappa^2(\rho+p)}{ f'(R)}\frac{v}{\mathcal H} + \frac{1}{2}\frac{\left( f'(R)\right)_N}{ f'(R)}\left(2\Psi^+ - 3\Xi\right) - \frac{2}{3\mathcal H^2}\frac{f_R}{\left( f'(R)\right)_N}\left[2k^2\Psi^+ + \frac{a^2\kappa^2}{ f'(R)}\left(\delta\rho-3\mathcal{H}(\rho+p)v\right)\right]
+\frac{2}{\mathcal{H}}\frac{ f'(R)}{\left( f'(R)\right)_N}\left(\mathcal H - \mathcal H_N\right)\Xi.
\end{align}
Eqs. (137) and (138) are equivalent to those obtained in Ref. (Pogosian and Silvestri, 2007). The evolution equation for $\Psi^+$ (137) can be used to eliminate the dependence of $(\Xi)_N$ on $(\Psi^+)_N$ in Eq.~(138).
This is, finally, the evolution equation for the metric perturbations.

We need now the evolution equations for perturbed matter quantities that are $\delta\rho$ and $v$. For $L$ perfect fluids, each with an energy-momentum tensor as (132), we can define the total, $\delta$, and individual, $\delta^{(i)}$, relative density perturbation, respectively, as \begin{eqnarray}\tag{139} \delta = \sum_{i=1}^L\frac{\rho^{(i)}}{\rho}\delta^{(i)}, \qquad \delta^{(i)} \equiv \frac{\delta\rho^{(i)}}{\rho^{(i)}}, \end{eqnarray} while the total and individual peculiar velocities are related by \begin{align}\tag{140} v &= \sum_{i=1}^L \frac{\rho^{(i)}+p^{(i)}}{\rho+p}v^{(i)}. \end{align} By perturbing the conservation equations of the energy momentum tensor, we find that, for adiabatic and non-interacting fluids, each pair $\delta^{(i)}$ and $v^{(i)}$ satisfies the equations (Mukhanov, Feldeman and Branderberger, 1992) \begin{align} \tag{141} &\left(\delta^{(i)}\right)_N + 3\left(c_s^{(i)2}-w^{(i)}\right) - (1+w^{(i)})k^2\frac{v^{(i)}}{\mathcal H}= 3(1+w^{(i)})\left(\Psi^+ + \frac{1}{2}\Xi\right)_N, \\ \tag{142} &\left(v^{(i)}\right)_N + \left(1-3c_s^{(i)2} \right)v^{(i)} + \frac{c_s^{(i)2}}{1+w^{(i)}}\frac{\delta^{(i)}}{\mathcal H} = - \frac{1}{\mathcal H}\left(\Psi^+ -\frac{1}{2} \Xi\right). \end{align} Here, $w^{(i)}\equiv p^{(i)}/\rho^{(i)}$ and $c_s^{(i)2}\equiv p^{(i)}_N/\rho^{(i)}_N$ are, respectively, the state parameter and the squared speed of sound of the $(i)$-fluid.

Eqs. (137), (138), (141), and (142) form a closed set of equations that allows the evolution of perturbations since the early radiation epoch until the present time. Once the present day value of the perturbation variables is computed, they can be related with observable quantities.

However also other approaches are possible for perturbation theory in **$f(R)$ gravity**. For example, the evolution of density perturbations can be dealt with the Covariant and Gauge Invariant theory. In (Carloni, Dunsby and Troisi, 2008), $f(R)\sim R^n$ power-law models are are dealt with general perturbations equations obtaining exact solutions for scales much bigger than the Hubble radius. These solutions have a some interesting features. In particular, all values of $n$ show growing modes for the density contrast, even if the universe undergoes an accelerated expansion. Such a behavior does not occur in General Relativity, where as soon as dark energy dominates, the density contrast experiences a decay. This peculiarity is a particular feature of fourth order gravity models.

### Matter Power Spectrum and Large Scale Structure constraints on $f(R)$ gravity

The evolution of scalar perturbations in**$f(R)$ gravity** can be tracked since the radiation epoch until the present time, and compared with the evolution of perturbations in the $\Lambda$CDM model. In particular, the modifications in the matter power spectrum due to **$f(R)$ gravity** can be compared, in principle, with today measured cosmological parameters.

While cosmography is one of the several approaches to test the background universe, the matter power spectrum is a tool to test the clumpy universe. Such a spectrum provides information about the process of clustering of matter in the universe at all scales, *i.e.*, from those corresponding to the early radiation dominated epoch up to those which have recently reentered the horizon. In addition, it is a formidable tool to break the possible degeneracy between General Relativity and modified gravity at the background level. In order to obtain the matter power spectrum, we need the evolution of the scalar perturbations as deduced above. In other words, the perturbative analysis carried above allows to see how the $f(R)$ corrections affect the matter power spectrum $\mathcal{P}_{\delta_{m}}$.
Notice however that the correct definition of $\mathcal{P}_{\delta_{m}}$ uses the energy density perturbation $\delta_m$ in the comoving gauge (for a discussion on this point see *e.g.* Refs. (Wands and Slozar,2009) (Bruni, Crittenden, Koyama, Maartens, Pitrou and Wands,2011), (Clifton, Ferreira,Padilla and Skordis2012)), while the variables used here correspond to the Newton gauge. In order to account for this gauge difference, the matter power spectrum can be calculated as
\begin{equation}\tag{143}
\mathcal{P}_{\delta_{m}} \equiv \langle|\delta_{m}^{(com)}|^2\rangle = \langle|\delta_{m} - 3\mathcal{H}v_{m}|^2\rangle.
\end{equation}

As an example, let as consider a cosmological model filled with dust and an effective energy density given by **$f(R)$ gravity**. In principle, the $f(R)$ model can account for the radiation content of the universe and drives its late time acceleration.
The initial conditions can be chosen as
\begin{align}
\tag{144}
\Psi^+(N_{ini})=1,
~~~~~~
\left(\Psi^+\right)_N(N_{ini})=0,
~~~~~~
\Xi(N_{ini})=0,
~~~~~~
\left(\Xi\right)_N(N_{ini})=0,
\end{align}
to match the initial conditions in the $\Lambda$CDM model with a radiation component and a single inflationary field, *i.e.*, the field $\Psi^+$ is initially constant and the modes are well outside of the horizon. During the radiation dominated epoch, the $f(R)$ corrections, in general, increase the value of the metric perturbations by several orders of magnitude, in contrast with what happens in General Relativity. Furthermore, from the above analysis, the initial value of the density perturbation $\delta$ are related to $\Psi^+(N_{ini})$ as
\begin{equation}
\tag{145}
\delta_m(N_{ini}) =
- \frac{6\mathcal{H}^2f_R}{a^2\kappa^2\rho}\left[1+\frac{\left(f_R\right)_N}{f_R} + \frac{k^2}{3\mathcal{H}^2}\right]\Psi^+ (N_{ini})\,.
\end{equation}
Initially one can have $(f_R)_N/f_R\ll1$, furthermore, given that all the relevant modes are outside of the horizon, we have that $k^2/(3\mathcal{H}^2)\ll1$. Therefore, the relation between $\delta$ and $\Psi^+$ is ruled by the factor $6\mathcal{H}^2f_R/a^2\kappa^2\rho$. In General Relativity, the Friedmann equation gives $3\mathcal{H}^2=a^2\kappa^2\rho$ ($f_R=1)$ and the matter and metric perturbations have the same order of magnitude. However, in **$f(R)$ gravity**, if the matter energy density is initially sub-dominant with respect to the $f(R)$ effective energy density, the factor on the right hand side of Eq.~(145) can become large, and so does $\delta$. Considering the General Relativity limit on Eq.~(145) we obtain $\delta_m=-2\Psi^+$. This relation refers to the total matter energy density perturbation $\delta$. For a General Relativity model, beside matter, we have also radiation, therefore we need two initial conditions for $\delta_m$ and $\delta_r$. Those are obtained by imposing the adiabaticity conditions, see for example (Dodelson, 2003), ( Amendola and Tsujikawa, 2010)
\begin{equation}\tag{146}
\frac{\delta_m}{1+w_m} = \frac{\delta_r}{1+w_r} = \frac{\delta}{1+w},
\end{equation}
where $w_m=0$, $w_r=1/3$, and initially $w\approx1/3$. Therefore, in a General Relativity model with matter, radiation and initial adiabatic perturbations, we have the initial relation between the matter and the metric perturbations
\begin{align}\tag{147}
\delta_m = \frac{3}{4}\delta_r=\frac{3}{4}\delta = -\frac{3}{2}\Psi^+.
\end{align}
The presence of radiation and the adiabaticity of the initial perturbations thus corrects the General Relativity limit of Eq.~(145) by a factor of $3/4$.

In **$f(R)$ gravity**, the initial value of $\delta_m$ can be orders of magnitude different than the one in General Relativity due to the effective stress-energy tensor. As a result, the matter power spectrum can be very different from the one of the $\Lambda$CDM model. Besides the overall increase/decrease in the amplitude of the spectrum, the $f(R)$ effects change the shape of the spectrum, in particular on the high $k$ regime. This is consistent with the scale-dependent effects of $f(R)$ on the matter perturbation (Pogosian and Silvestri, 2007) (de la Cruz-Dombriz, Dobado and Maroto, 2008).

In summary, one can study the evolution of the first order metric and matter perturbations since the radiation dominated epoch until the present time and obtain the theoretical power spectrum of the matter perturbation. The numerical methods to achieve results can be based, for example, on the so called "quasi static approximation" (de la Cruz-Dombriz, Dobado and Maroto, 2008). In general, the effects of the $f(R)$ corrections on the matter perturbation are stronger on the modes with higher wave-numbers, as is expected from the corrections to the closed evolution equation of $\delta_m$ (de la Cruz-Dombriz, Dobado and Maroto, 2008). More precisely, this can be seen directly on the matter power spectra obtained: for low $k$ end of the spectrum, the shape of the spectra obtained is consistent with the one of the $\Lambda$CDM model (up to a multiplicative factor); for modes with $k> k_{eq}\approx0.01$ Mpc$^{-1}$ the shape of the spectra obtained starts to diverge from the one in the $\Lambda$CDM model. Notice however that the matter power spectrum can be obtained assuming that the initial amplitude profile of the initial conditions is that of $\Phi^+$ for single field inflation, that is
\begin{equation}\tag{148}
\mathcal P _{\Psi^+} = \langle |\Psi^+|^2\rangle = \frac{8\pi^2}{9k^3}A_s\left(\frac{k}{k_0}\right)^{n_s-1},
\end{equation}
which could be modified due to $f(R)$ effects. Here, $A_s$ is the curvature $\mathcal R$ power spectrum at some pivot scale $k_0$ ; $n_s$ is the scalar spectrum index that can be obtained by observations (*e.g.* the PLANCK mission (Ade et al.,2014). The relation $\mathcal{R}=(3/2)\Psi^+$ can be used during the radiation epoch and for modes outside of the horizon.
As a final remark, Large Scale Structure and perturbation analysis vs observations can put, in principle, strong constraints on the form and the evolution of any **$f(R)$ gravity** model.

## Dark Matter and Dark Energy as Curvature Effects

The above considerations can be worked out in view of addressing the problem of the dark side of the universe.
Dark matter issues come from dynamical mass estimates of self-gravitating systems.
In several astrophysical observations, there is more matter dynamically inferred than that can be accounted for from luminous components. This mass discrepancy is usually attributed to additional (missing or dark) matter, assuming the validity of Newton law of gravity at astrophysical scales.
Oort was the first that posed the *missing matter* problem (Oort, 1932) (Oort, 1960). By observing the Doppler red-shift values of stars moving near the plane of our Galaxy, he asserted that he could calculate how fast the stars were moving. He found that there had to be enough matter inside the Galaxy such that the central gravitational force was strong enough to keep the stars from escaping, much as Sun's gravitational pull keeps a planet in its orbit. But when the calculation was made, it turned out that there was not enough mass in the Galaxy. The discrepancy was not small: the Galaxy had to be at least twice as massive as the sum of the mass of all its visible components combined.
In addition, in the 1960's the radial profile of the tangential velocity of stars in their orbits around the Galactic Center, as a function of their distance from that center, was measured. It was found that typically, once we get away from the Galactic Center, all the stars travel with the same velocity independent of their distance out from the Galactic Center.

There were problems, too, at a larger scale. In 1933, Zwicky announced that when he measured the individual velocities of a large group of galaxies, known as the Coma Cluster, he found that all of the galaxies that he considered were moving so rapidly relative to one another that the cluster should have come apart long ago. The visible mass of the galaxies making up the cluster was far too little to produce enough gravitational force to hold the cluster together. So not only our own Galaxy was lacking mass, but also the whole Coma Cluster of galaxies was suffering the same problem at a different scale (Zwicky, 1933). Initially, the problem was only approached by leaving Newton's law inviolated and postulating the existence of some invisible dark entities to make up the missing mass. At the beginning, it has never came to mind anyone to go back and examine the basic assumption that only gravity was at work in these cases. It was easier to patch up the theory introducing invisible entities. Many names have been coined to define this invisible entity, a bit as in the days of ether.

There are the Massive Compact Halo Objects (MACHOs), objects like black holes, and neutron stars that purportedly populate the outer reaches of galaxies like the Milky Way. Then there are the Weakly Interacting Massive Particles (WIMPs), which possess mass, yet do not interact with ordinary matter (baryons such as protons and neutrons) because they are composed by something unknown. Dark (missing) matter (DM) even comes in two flavors, hot (HDM) and cold (CDM). The CDM is supposedly to be in dead stars, planets, brown dwarfs ("failed stars") etc., while HDM is postulated to be fast moving in particles floating throughout the universe. It should be constituted by neutrinos, tachyons etc. But where is all of this missing matter? The truth is that after many years of looking for it, there is still no definitive proof that WIMPs exist, or that MACHOs will ever make up more than five percent of the total reserve of missing dark stuff.

Besides, by adding a further ingredient, the cosmological constant $\Lambda$, such a model (now $\Lambda$CDM) has become the new cosmological paradigm usually called the *concordance model*.
In fact, high quality data coming from the measurements of
cluster properties as the mass, the correlation function and the
evolution with redshift of their abundance, the Hubble diagram of Type Ia
Supernovae, the optical surveys of
large scale structure, the anisotropies
of the cosmic microwave background, the cosmic
shear measured from weak lensing surveys and
the Lyman-$\alpha$ forest absorption are
evidences toward a spatially flat universe with a subcritical
matter content and undergoing a phase of accelerated expansion.
Interpreting all this information in a self-consistent model is
the main task of modern cosmology and $\Lambda$CDM model provides
a good fit to the most part of the data giving a
reliable picture of the today observed universe.

Nevertheless, it is affected by serious theoretical shortcomings
that have motivated the search for alternative candidates
generically referred to as *dark energy* or *quintessence*. Such models range from scalar fields rolling down
self interaction potentials to phantom fields, from
phenomenological unified models of dark energy and dark matter to
alternative gravity theories.

Essentially, dark energy (or any alternative component) has to act as a negative pressure fluid which gives rise to an overall acceleration of the Hubble fluid. Despite of the clear mechanisms generating the observed cosmological dynamics, the nature and the fundamental properties of dark energy remain essentially unknown notwithstanding the great theoretical efforts made up to now.

The situation for dark matter is similar: its clustering and distribution properties are fairly well known at every scale but its nature is unknown, up to now, at fundamental level.

On the other hand, the need of unknown components as dark energy
and dark matter could be considered nothing else but as a signal
of the breakdown of Einstein General Relativity at astrophysical
(galactic and extragalactic) and cosmological scales.
In this context, **Extended Theories of Gravity** could be, in
principle, an interesting alternative to explain cosmic
acceleration and large scale structure
without any missing components. In their simplest version, the Ricci
curvature scalar $R$, linear in the Hilbert-Einstein action, could
be replaced by a generic function $f(R)$ whose true form could be
"reconstructed" by the data. In fact, there is no a priori reason
to consider the gravitational Lagrangian linear in the Ricci
scalar while observations and experiments could contribute to
define and constrain the "true" theory of gravity (Capozziello and Faraoni, 2010) (Capozziello and De Laurentis, 2010).

Coming to the weak-field limit, which essentially means considering Solar System scales, any alternative relativistic theory of gravity is expected to reproduce General Relativity which, in any case, is firmly tested only in this limit and at these scales (Will,1993). Even this limit is a matter of debate since several relativistic theories do not reproduce exactly the Einstein results in their Newtonian limit but, in some sense, generalize them.

In general, all these efforts can be included in the so called Modified Newtonian Dynamics (MOND), first proposed by Milgrom (Milgrom, 1983) (Milgrom and Sanders, 1983), in the attempt to explain missing matter without dark matter but assuming a change into dynamics at scales larger than Solar System's ones.
In general, any relativistic theory of gravitation yields
corrections to the weak-field gravitational potentials which, at the post-Newtonian level and
in the Parametrized Post-Newtonian formalism,
could constitute a test of these theories (Will,1993).
This point deserves a deep discussion. Beside the fundamental physics motivations coming from Quantum Gravity and unification theories (Capozziello and De Laurentis, 2011), **Extended Theories of Gravity** pose the problem that there are further gravitational degrees of freedom (related to higher order terms, non-minimal couplings and scalar fields in the field equations) and gravitational interaction is *not* invariant at any scale. This means that, besides the Schwarzschild radius, other characteristic gravitational scales come out from dynamics. Such scales, in the weak field approximation, should be responsible of characteristic lengths of astrophysical structures that should result *confined* in this way (Capozziello and De Laurentis,2012).

## The Weak Field Limit of $f(R)$ gravity

We can deal with the Newtonian and the post-Newtonian limit of **$f(R)$ gravity**
adopting the spherical symmetry. The solution of field equations
can be obtained considering the general spherically symmetric metric:
\begin{eqnarray}\tag{149}
ds^2\,=\,g_{\sigma\tau}dx^\sigma
dx^\tau=g_{00}(x^0,r)d{x^0}^2-g_{rr}(x^0,r)dr^2-r^2d\Omega
\end{eqnarray}
where $x^0\,=\,ct$ and $d\Omega$ is the solid angle.
In order to develop the Newtonian limit, let us consider the
perturbed metric with respect to a Minkowskian background
$g_{\mu\nu}\,=\,\eta_{\mu\nu}+h_{\mu\nu}$. The metric entries can
be developed as:
\begin{eqnarray}\tag{150}
\left\{\begin{array}{ll} g_{tt}(t,
r)\simeq1+g^{(2)}_{tt}(t,r)+g^{(4)}_{tt}(t,r)\,,
\\\\
g_{rr}(t,r)\simeq-1+g^{(2)}_{rr}(t,r)\,,\\\\
g_{\theta\theta}(t,r)=-r^2\,,\\\\
g_{\phi\phi}(t,r)=-r^2\sin^2\theta\,,
\end{array}\right.\,
\end{eqnarray}
where we are assuming $c\,=\,1$,$x^0=ct\rightarrow t$ and expansion is at order $c^{-2}$ and $c^{-4}$ .
Since we want to obtain the most general result, we does not
provide any specific form for $f(R)$ Lagrangian. We assume,
however, analytic Taylor expandable $f(R)$ functions with respect
to a certain value $R\,=\,R_0$:

\begin{eqnarray}\tag{151} f(R)=\sum_{n}\frac{f^n(R_0)}{n!}(R-R_0)^n\simeq f_0+f_1R+f_2R^2+f_3R^3+...\,, \end{eqnarray} where we have assumed $R_0=0$. In order to obtain the weak field approximation, one has to insert expansions (150) and (151) into field Eqs. (23) and expand the system up to the orders ${\mathcal O}(0)$, ${\mathcal O}(2)$ e ${\mathcal O}(4)$ (that is, as stated above at order $c^{-2}$ and $c^{-4}$). This approach provides general results and specific (analytic) theories are selected by the coefficients $f_i$ in Eq.(151). It is worth noticing that, at the order ${\mathcal O}(0)$, the field equations give the condition $f_0 =0$ and then the solutions at further orders do not depend on this parameter. If we consider the ${\mathcal O}(2)$-order approximation, the field equations in vacuum, results to be \begin{eqnarray}\tag{152} \left\{\begin{array}{ll} f_1rR^{(2)}-2f_1g^{(2)}_{tt,r}+8f_2R^{(2)}_{,r}-f_1rg^{(2)}_{tt,rr}+4f_2rR^{(2)}=0\,, \\\\ f_1rR^{(2)}-2f_1g^{(2)}_{rr,r}+8f_2R^{(2)}_{,r}-f_1rg^{(2)}_{tt,rr}=0\,, \\\\ 2f_1g^{(2)}_{rr}-r\left[f_1rR^{(2)}-f_1g^{(2)}_{tt,r}-f_1g^{(2)}_{rr,r}+4f_2R^{(2)}_{,r}+4f_2rR^{(2)}_{,rr}\right]=0\,, \\\\ f_1rR^{(2)}+6f_2\left[2R^{(2)}_{,r}+rR^{(2)}_{,rr}\right]=0\,, \\\\ 2g^{(2)}_{rr}+r\left[2g^{(2)}_{tt,r}-rR^{(2)}+2g^{(2)}_{rr,r}+rg^{(2)}_{tt,rr}\right]=0\,. \end{array} \right.\end{eqnarray} It is evident that the trace equation (the fourth in the system (152)), provides a differential equation with respect to the Ricci scalar which allows to solve exactly the system (152) at ${\mathcal O}(2)$-order. Finally, one gets the general solution:

\begin{eqnarray}\tag{153}
\left\{\begin{array}{ll}
g^{(2)}_{tt}=\delta_0-\frac{Y}{f_1r}-\frac{\delta_1(t)e^{-r\sqrt{-\xi}}}{3\xi
r}+\frac{\delta_2(t)e^{r\sqrt{-\xi}}}{6({-\xi)}^{3/2}r}
\\\\
g^{(2)}_{rr}=-\frac{Y}{f_1r}+\frac{\delta_1(t)[r\sqrt{-\xi}+1]e^{-r\sqrt{-\xi}}}{3\xi
r}-\frac{\delta_2(t)[\xi r+\sqrt{-\xi}]e^{r\sqrt{-\xi}}}{6\xi^2r}
\\\\
R^{(2)}=\frac{\delta_1(t)e^{-r\sqrt{-\xi}}}{r}-\frac{\delta_2(t)\sqrt{-\xi}e^{r\sqrt{-\xi}}}{2\xi
r}\end{array} \right.
\end{eqnarray}
where $\xi\doteq\displaystyle\frac{f_1}{6f_2}$,
and $Y$ is an arbitrary integration
constant.
When we consider the limit $f(R)\rightarrow R$, in the case
of a point-like source of mass $M$, we recover the standard
Schwarzschild solution. Let us notice that the
integration constant $\delta_0$ is dimensionless, while the two
arbitrary functions of time $\delta_1(t)$ and $\delta_2(t)$ have
respectively the dimensions of $length^{-1}$ and $length^{-2}$.
The functions of time $\delta_i(t)$ ($i=1,2$) are completely arbitrary since the
differential equation system (152) contains only spatial
derivatives and can be fixed to constant values. Besides, the
integration constant $\delta_0$ can be set to zero since it
represents an unessential additive quantity for the potential.
It is possible to write the general solution of the problem considering
the previous expression (149). In order to match at infinity
the Minkowskian prescription for the metric, one can discard the
Yukawa growing mode and then we have:
\begin{eqnarray}\tag{154}
\left\{\begin{array}{ll}ds^2\,=\,\biggl[1-\frac{2GM}{f_1r}-\frac{\delta_1(t)e^{-r\sqrt{-\xi}}}{3\xi
r}\biggr]dt^2- \biggl[1+\frac{2GM}{f_1r}-\frac{\delta_1(t)(r\sqrt{-\xi}+1)e^{-r\sqrt{-\xi}}}{3\xi
r}\biggr]dr^2-r^2d\Omega\\\\R\,=\,\frac{\delta_1(t)e^{-r\sqrt{-\xi}}}{r}\end{array}\right.
\end{eqnarray}
At this point, one can provide the
gravitational potential. The first of (153) gives the
second order solution in term of the metric expansion (see the
definition (150)). This term coincides with the
gravitational potential at the Newton order. In particular, since
$g_{tt}\,=\,1+2\Phi_{grav}\,=\,1+g_{tt}^{(2)}$, the
gravitational potential of $f(R)$-gravity, analytic in the Ricci
scalar $R$, is
\begin{eqnarray}\tag{155}
\Phi_{grav}\,=\,-\left[\frac{GM}{f_1r}+\frac{\delta_1(t)e^{-r\sqrt{-\xi}}}{6\xi
r}\right]\,.
\end{eqnarray}
This general result means that the standard Newton potential is
achieved only in the particular case $f(R)=R$ while it is not so
for any analytic $f(R)$ models. Eq.(155) deserves
some comments. The parameters $f_{1,2}$ and the function
$\delta_1$ represent the deviations with respect the standard
Newton potential. To test these theories of gravity inside the
Solar System, we need to compare such quantities with respect to
the current experiments, or, in other words, Solar System
constraints should be evaded fixing such parameters. On the other
hand, these parameters could acquire non-trivial values (*e.g.*
$f_1\neq 1,\,\delta_1(t)\neq 0,\,\xi\neq 1$) at scales different
from the Solar System ones.
We note that the $\xi$ parameter can be related to an effective mass being
$
m^2= (3\xi)^{-1}
$
and can be interpreted also as an effective length $L$.
Eq. (155) can be recast as
\begin{equation}
\tag{156} \Phi(r) = -\frac{G M}{
(1+\delta) r}\left(1+\delta e^{-\frac{r}{L}}\right)\,,
\end{equation}
where the first term is the Newtonian-like part of the potential
for a point-like mass $\displaystyle{\frac{M}{1+\delta}}$ and the second
term is a modification of gravity including a scale
length $L$ associated to the above coefficients of the Taylor expansion.
If $\delta=0$ the Newtonian potential and the standard gravitational coupling are recovered.
Comparing Eqs.(155) and (156), we obtain that
$\displaystyle{1+\delta=f_1}$, and $\delta$ is related to $\delta_1(t)$
through
\begin{equation}\tag{157}
\delta_1=-\frac{6GM}{L^2}\left(\frac{\delta}{1+\delta}\right)
\end{equation}
where $\displaystyle{\frac{6GM}{L^2}}$ and $\delta_1$ can be assumed quasi-constant.
Under these assumptions, the scale length $L$ could naturally reproduce several phenomena that range from Solar System to cosmological scales. Understanding
on which scales the modifications to General Relativity are working or what is the weight of corrections
to the Newtonian gravitational potential is a crucial point that could confirm or rule out these alternative approaches (Capozziello and De Laurentis, 2012).

## Solar system and Equivalence Principle constraints

Solar System tests for relativistic theories of gravity include gravitational redshift, deflection of light by the Sun, planetary orbit precession at perihelion. General Relativity is consistent with these experimental tests. Considering them gives a useful framework where predictions of different relativistic theories can be parameterized in a systematic way (Capozziello and De Laurentis, 2011). In particular, the Post Newtonian limit framework (Will,1993) has become a basic tool to connect gravitational theories with Solar System experiments.
Furthermore, Post Newtonian formalism can be applied to all metric theories of gravitation where self-gravitating bodies satisfy the Einstein Equivalence Principle. In general, it is possible to find constraints on $f(R)$ models (and any **Extended Theories of Gravity**) from Solar System experiments combined with experiments on the violation of Equivalence Principle (Capozziello and Tsujikawa, 2008).

In an environment of high density such as Earth or Sun, the Ricci scalar $R$ is larger than the cosmological background. If the outside of a spherically symmetric body is the vacuum, the metric can be described by a Schwarzschild-like exterior solution with $R=0$. In the presence of non-relativistic matter with a matter-energy density $\rho_m$, this gives rise to a contribution to the Ricci scalar $R$ of the order $\rho_m$. The task is now to test if the above result involving the Yukawa correction to the potential, is compatible or not with Solar System tests and Equivalence Principle. We will deal with the further degrees of freedom, coming from $f(R)$, that can be dealt under the standard of an additional scalar field coming into dynamics. This approach is allowed by conformal transformations.

When the mass of an equivalent scalar field degree of freedom is heavy in a region with high density, a spherically symmetric body can be represented as a thin-shell, so that an effective coupling due to the fifth force (the Yukawa correction) is suppressed through a Chameleon Mechanism ( Khoury and Weltman, 2004). In this mechanism, the mass of the effective scalar field depends on the density environment. If the matter density is sufficiently high, the field acquires a heavy mass about the potential minimum. On the other hand, the field has lighter mass in a low-density cosmological environment relevant for dark energy. In this case, it can propagate freely. As long as a spherically symmetric body has thin-shell around its surface, the effective coupling between the field and matter becomes much smaller than the coupling (De Felice and Tsujikawa, 2010).

Starting from the action for $f(R)$, introducing the new metric $\tilde{g}_{\mu \nu}$ and a scalar field $\phi$, conformal transformations allow to write

\begin{eqnarray}\tag{158} \tilde{g}_{\mu \nu}=\psi g_{\mu \nu}\,, \quad \phi=\sqrt{\frac{3}{2}}\,{\rm ln}\,\psi\,, \end{eqnarray}

where, clearly, the conformal transformation is specified by $\psi=\partial f/\partial R$. We can rewrite the $f(R)$ action in the Einstein frame as (Maeda, 1989) (Capozziello and Faraoni, 2010) (Capozziello and De Laurentis, 2010), in physical units \begin{eqnarray} \tag{159} {\cal S} &=&\int{\rm d}^{4}x\sqrt{-\tilde{g}}\left[ -\frac{1}{2}\tilde{R}+\frac{1}{2} g^{\mu\nu} \phi_{;\mu}\phi_{;\nu}-V(\phi) \right] +{\cal L}^{(m)} (\tilde{g}_{\mu \nu} e^{2\beta \phi}, \Psi_m)\,, \end{eqnarray}

where

\begin{eqnarray}\tag{160} \beta=-\frac{1}{\sqrt{6}}\,,\quad V(\phi(\psi))=\frac{R(\psi)\psi-f}{2\psi^2}\,. \end{eqnarray}

The field $\phi$ is directly coupled to the non-relativistic standard matter by a coupling $\beta$ that we will specify below. In a spherically symmetric space-time, the variation of the action (159) with respect to the scalar field $\phi$ gives

\begin{eqnarray} \tag{161} \frac{d^2 \phi}{d {\tilde{r}}^2}+ \frac{2}{{\tilde{r}}} \frac{ d\phi}{d {\tilde{r}}}= \frac{d V_{eff}}{d\phi}\,, \end{eqnarray}

where ${\tilde{r}}$ is the distance from the center of symmetry and

\begin{eqnarray} \tag{162} V_{eff}(\phi)=V(\phi)+ e^{\beta \phi}\rho^*\,. \end{eqnarray}

Here $\rho^*$ is a conserved quantity in the Einstein frame. It is related to the matter-energy density $\rho$, in the Jordan frame, via the relation $\rho^*=e^{3\beta \phi}\rho$.

Assuming that a spherically symmetric body has a constant density $\rho^*=\rho_A^*$ inside the body (${\tilde r}<{\tilde {r}}_c$) and that the energy density outside the body (${\tilde r}>{\tilde {r}}_c$) is $\rho^*=\rho_B^*$, the mass $M_c$ of the body and the gravitational potential $\Phi_c$ at the radius ${\tilde r}_c$ are given by $M_c=(4\pi/3){\tilde r}_c^3 \rho_A^*$ and $\Phi_c=M_c/8\pi {\tilde{r}}_c$, respectively. The effective potential $V_{eff} (\phi)$ has two minima at the field values $\phi_A$ and $\phi_B$ satisfying $V_{eff}' (\phi_A)=0$ and $V_{eff}' (\phi_B)=0$, respectively. The former corresponds to the region with a high density that gives rise to a heavy mass squared $m_A^2 \equiv V_{eff}''(\phi_A)$, whereas the latter to the lower density region with a lighter mass squared $m_B^2 \equiv V_{eff}''(\phi_B)$. In the high-density regime (with a heavy field mass), the spherically symmetric body acquires a thin-shell under the Chameleon Mechanism. When the thin-shell develops inside the body, the following thin-shell parameter is much smaller than the order of unity (see (Capozziello and Tsujikawa, 2008) and references therein):

\begin{eqnarray} \tag{163} \frac{\Delta {\tilde{r}}_c}{{\tilde{r}}_c} =\frac{\phi_B-\phi_A}{6\beta \Phi_c}\,. \end{eqnarray}

Solving Eq.(161) with appropriate boundary conditions, the field profile outside the body (${\tilde{r}}>{\tilde{r}}_c$) is given by

\begin{eqnarray} \tag{164} \phi({\tilde{r}}) \simeq -\frac{\beta_{eff}}{4\pi} \frac{M_c e^{-m_B({\tilde{r}}-{\tilde{r}}_c)}}{{\tilde{r}}}+\phi_B\,, \end{eqnarray}

where the magnitude of the effective coupling, $\beta_{eff}= (3\beta) (\Delta {\tilde{r}}_c/{\tilde{r}}_c)$, is much smaller than unity when the thin-shell is formed. The bound on the thin-shell parameter from experimental tests of the Post Newtonian parameter in Solar Systems can be derived. The spherically symmetric metric in Einstein frame is given by

\begin{eqnarray} \tag{165} d{\tilde s}^2=\psi d s^2={\tilde{g}}_{\mu\nu} d{\tilde{x}}^\mu d{\tilde{x}}^\nu=\left[1-2{\cal {\tilde{A}}}({\tilde{r}})\right]dt^2-\left[1+2{\cal {\tilde{B}}}({\tilde r})\right]d{\tilde{r}}^2-{\tilde{r}}^2d\Omega^2\,, \end{eqnarray} where ${\tilde{A}}({\tilde{r}})$ and ${\tilde B}({\tilde r})$ are functions of ${\tilde r}$ and, as above, $d\Omega^2=d\theta^2+(\sin^2 \theta)d\phi^2$. In the weak field limit, it is ${\tilde{A}}({\tilde{r}})\ll 1$ and ${\tilde B}({\tilde{r}})\ll 1$. The metric outside the spherically symmetric body with mass $M_c$ is given by ${\tilde A}({\tilde {r}})\simeq{\tilde B}({\tilde{r}})\simeq \frac{GM_c}{{\tilde {r}}}$. Let us transform back the metric (165) into the Jordan frame. It is \begin{eqnarray}\tag{166} {\rm d}s^2=[1-2{\cal{A}}(r)]{\rm d}t^2 -[1+2{\cal{B}}(r)]{\rm d}r^2-r^2 {\rm d} \Omega^2\,, \end{eqnarray} Comparing the two metrics under the condition $|\phi| \ll 1$, we find that the following relations hold (Faulkner et al, 2007)

\begin{eqnarray}\tag{167} r = \psi^{-1/2}\tilde{r}\,, \quad {\cal{A}}(r) \simeq \tilde{{\cal{A}}}(\tilde{r})+ \frac{\phi(\tilde{r})}{\sqrt{6}}\,, \quad {\cal{B}}(r) \simeq \tilde{{\cal{B}}}(\tilde{r})+ \frac{\tilde{r}}{\sqrt{6}} \frac{{\rm d}\phi(\tilde{r})}{{\rm d} \tilde{r}}\,. \end{eqnarray}

Under the weak-field approximation in the Einstein frame, the metric components $\tilde{{\cal{A}}}(\tilde{r})$ and $\tilde{{\cal{B}}}(\tilde{r})$ outside a spherically symmetric body with mass $M_c$ are given by $\tilde{{\cal{A}}}(\tilde{r}) \simeq \Phi_c=M_c/(8\pi \tilde{r})$. Provided that $|\phi| \ll 1$, one has $\psi \simeq 1$ and hence $\tilde{r} \simeq r$. In the following, we omit a tilde for the quantities in the Einstein frame because the condition $|\phi| \ll 1$ always holds for the analytic $f(R)$ models we are studying. Using the thin-shell solution (164), we find that the metrics ${\cal{A}}(r)$ and ${\cal{B}}(r)$ in the Jordan frame are

\begin{eqnarray}\tag{168} {\cal{A}}(r)=\frac{G_{\rm eff}M_c}{r}\,,\quad {\cal{B}}(r)=\gamma \frac{G_{\rm eff}M_c}{r}\,, \end{eqnarray}

where the effective gravitational coupling $G_{\rm eff}$ and the post-Newtonian parameter $\gamma$ are given by

\begin{eqnarray} \tag{169} G_{\rm eff} &\simeq& G \left[ 1-\frac{\sqrt{6}}{3}\beta_{\rm eff} e^{-m_B (r-r_c)} \right]\,,\\ \\ \tag{170} \gamma &\simeq& \frac{1+(\sqrt{6}\beta_{\rm eff}/3)(1+m_{B} r) e^{-m_B (r-r_c)}}{1-(\sqrt{6}\beta_{\rm eff}/3) e^{-m_B (r-r_c)}}\,. \end{eqnarray}

If the condition $m_Br \ll 1$ is satisfied in an environment where local gravity experiments are carried out, we obtain

\begin{eqnarray}\tag{171} G_{\rm eff} &\simeq& G \left( 1+\frac{\Delta r_c}{r_c} \right)\,, \\ \\ \tag{172} \gamma &\simeq& \frac{1-\Delta r_c/r_c} {1+\Delta r_c/r_c}\,. \end{eqnarray}

As long as the body has a thin-shell ($|\Delta r_c/r_c| \ll 1$), one can satisfy the experimental bounds $|\gamma-1| \ll 1$ (Will,1993). In the thick-shell regime, we just need to change $\beta_{\rm eff}$ in Eqs.(169)-(170) to $\beta=-1/\sqrt{6}$. Under the condition $m_Br \ll 1$, we obtain $G_{\rm eff}=4G/3$ and $\gamma=1/2$, which contradicts with the experimental bound mentioned above. If the mass $m_B$ is heavy so that $m_Br \gg 1$ for $r \gtrsim r_c$, the mass $m_A$ inside the body is typically much heavier to satisfy the relation $m_A r_c \gg 1$. In such cases the body has a thin-shell, which means that the post-Newtonian parameter is given by Eq.(172).

In Ref. (Faulkner et al, 2007), it was shown that under the Chameleon Mechanism the post-Newtonian parameter, $\gamma={\cal B}(r)/{\cal A}(r)$, is approximately given by (172) provided that the condition $m_B r \ll 1$ holds at Solar System scales. The present tightest constraint on $\gamma$ is $|\gamma-1| <2.3 \times 10^{-5}$ (Will,1993), which translates into

\begin{eqnarray} \tag{173} \frac{\Delta r_c}{r_c}<1.15 \times 10^{-5}\,. \end{eqnarray} Let us now place experimental bounds from a possible violation of Equivalence Principle. The thin-shell condition around the Earth under the Chameleon Mechanism has to be considered ( Khoury and Weltman, 2004). The Earth has a radius $r_{\oplus}=6 \times 10^{3}$ km with a mean density $\rho_{\oplus} \simeq 5.5$ g/cm$^3$. The atmosphere exists in the region $r_{\oplus}<r<r_{atm}$ with a homogeneous density $\rho_{atm} \simeq 10^{-3}$ g/cm$^3$. The region outside the atmosphere ($r>r_{atm}$) has a homogenous density $\rho_G \simeq 10^{-24}$ g/cm$^3$. The atmosphere exists in the region $r_{\oplus}<r<r_{atm}$ with a homogeneous density $\rho_{atm} \simeq 10^{-3}$ g/cm$^3$. The region outside the atmosphere ($r>r_{atm}$) has a homogenous density $\rho_G \simeq 10^{-24}$ g/cm$^3$. Defining the gravitational potentials as $\Phi_{\oplus}=\rho_{\oplus}r_{\oplus}^2/6$ and $\Phi_{\rm atm}=\rho_{atm}r_{atm}^2/6$, we have that $\Phi_{\oplus} \simeq 5.5 \times 10^3 \Phi_{\rm atm}$ because $\rho_{\oplus} \simeq 5.5 \times 10^3 \rho_{atm}$ and $r_{\oplus} \simeq r_{\rm atm}$. Recalling the relation $\Delta r_{atm}/r_{atm}= (\phi_G-\phi_{atm})/(6\beta \Phi_{atm})$, where $\phi_G$ and $\phi_{atm}$ correspond to the field values at the local minima of the effective potential (162) in the regions $r>r_{atm}$ and $r_{\oplus}<r<r_{atm}$ respectively, we find $\Delta r_{\oplus}/r_{\oplus} \equiv -(\phi_G-\phi_{atm}) /\sqrt{6} \Phi_{\oplus} \simeq 2.0 \times 10^{-4} (\Delta r_{atm}/r_{atm})$. When the atmosphere has a thin-shell, then the thickness of the shell ($\Delta r_{atm}$) is smaller than that of the atmosphere: $r_s=10$-10$^2$km. Taking the value $r_s=10^2$ km and $r_{atm} =6.5 \times 10^3$ km, we obtain $\Delta r_{\rm atm}/r_{atm}<1.6 \times 10^{-2}$. Hence the condition for the atmosphere to have a thin-shell is estimated as

\begin{eqnarray} \tag{174} \frac{\Delta r_{\oplus}}{r_{\oplus}} \lesssim 10^{-6}\,. \end{eqnarray}

Solar System tests of Equivalence Principle consider, essentially, free-fall acceleration of the Moon and the Earth toward the Sun (see (Capozziello and De Laurentis, 2011) for a discussion on the Equivalence Principle). The constraint on the difference of the two accelerations is given by

\begin{eqnarray} \tag{175} \eta \equiv 2\frac{|a_{\rm Moon}-a_{\oplus}|} {a_{\rm Moon}+a_{\oplus}}<10^{-13}\,. \end{eqnarray}

The Sun and the Moon have the thin-shells like the Earth (Capozziello and Tsujikawa, 2008), where the field profiles outside the spheres are given as in Eq.(164) with the replacement of corresponding quantities. We note that the acceleration induced by a fifth force with the field profile $\phi(r)$ and the effective coupling $\beta_{eff}$ is $a^{fifth}= |\beta_{eff}\phi(r)|$. Then the accelerations $a_{\oplus}$ and $a_{Moon}$ are (De Felice and Tsujikawa, 2010) (Capozziello and Tsujikawa, 2008)

\begin{eqnarray}\tag{176} a_{\oplus} &\simeq& \frac{GM_{\odot}}{r^2} \left[ 1+3 \left(\frac{\Delta r_{\oplus}} {r_{\oplus}}\right)^2 \frac{\Phi_{\oplus}}{\Phi_{\odot}} \right]\,,\\ \\\tag{177} a_{\rm Moon} &\simeq& \frac{GM_{\odot}}{r^2} \left[ 1+3 \left(\frac{\Delta r_{\oplus}} {r_{\oplus}}\right)^2 \frac{\Phi_{\oplus}^2} {\Phi_{\odot}\Phi_{\rm Moon}}\right]\,, \end{eqnarray}

where $\Phi_{\odot} \simeq 2.1 \times 10^{-6}$, $\Phi_{\oplus} \simeq 7.0 \times 10^{-10}$ and $\Phi_{Moon} \simeq 3.1 \times 10^{-11}$ are the gravitational potentials of Sun, Earth and Moon, respectively. Hence the condition (175) translates into

\begin{eqnarray} \tag{178} \frac{\Delta r_{\oplus}}{r_{\oplus}}<2 \times 10^{-6}\,, \end{eqnarray}

which gives the same order of the upper bound as in the thin-shell condition (174) for the atmosphere. The constraint coming from the violation of strong equivalence principle (Will,1993) provides a bound $\Delta r_{\oplus}/r_{\oplus}<10^{-4}$ (Capozziello and Tsujikawa, 2008), which is weaker than (178).

In summary, any correction that one can take into account for the Newtonian potential has to satisfy the above constraints on Equivalence Principle and Solar System. This statement has to hold also for **$f(R)$ gravity**.

## Massive gravitational modes

**$f(R)$ gravity ** can be linearized producing a third
massive mode of gravitational radiation. Taking the trace of (25) one gets

\begin{equation} 3\square f'(R)+Rf'(R)-2f(R)=0,\tag{179}\end{equation}

and, with the identifications (Capozziello, Corda and De Laurentis, 2008)

\begin{equation} \begin{array}{ccccc} \Phi\rightarrow f'(R) & & \textrm{and } & & \frac{dV}{d\Phi}\rightarrow\frac{2f(R)-Rf'(R)}{3}\end{array}\tag{180}\end{equation}

a Klein-Gordon equation for the effective $\Phi$ scalar field is obtained

\begin{equation} \square\Phi=\frac{dV}{d\Phi}.\tag{181}\end{equation} We assume $\Phi_{0}$ to be a minimum for $V$

\begin{equation} V\simeq\frac{1}{2}\alpha\delta\Phi^{2}\Rightarrow\frac{dV}{d\Phi}\simeq m^{2}\delta\Phi,\tag{182}\end{equation}

and the constant $m$ has the mass dimension.

Assuming

\begin{equation} \begin{array}{c} g_{\mu\nu}=\eta_{\mu\nu}+h_{\mu\nu}\,, \qquad \Phi=\Phi_{0}+\delta\Phi.\end{array}\tag{183}\end{equation}

to first order in $h_{\mu\nu}$ and $\delta\Phi$. Considering $\widetilde{R}_{\mu\nu\rho\sigma}$, $\widetilde{R}_{\mu\nu}$ and $\widetilde{R}$ as the linearized quantities that correspond to $R_{\mu\nu\rho\sigma}$, $R_{\mu\nu}$ and $R$, the linearized field equations are

\begin{equation}\tag{184} \begin{array}{c} \widetilde{R}_{\mu\nu}-\frac{\widetilde{R}}{2}\eta_{\mu\nu}=(\partial_{\mu}\partial_{\nu}h_{f}-\eta_{\mu\nu}\square h_{f})\\ \\{}\square h_{f}=m^{2}h_{f},\end{array}\end{equation}

where

\begin{equation} h_{f}\equiv\frac{\delta\Phi}{\Phi_{0}}.\tag{185}\end{equation}

$\widetilde{R}_{\mu\nu\rho\sigma}$ and Eqs. (184) are invariants for gauge transformations

\begin{equation}\tag{186} \begin{array}{c} h_{\mu\nu}\rightarrow h'_{\mu\nu}=h_{\mu\nu}-\partial_{(\mu}\epsilon_{\nu)}\,, \quad\delta\Phi\rightarrow\delta\Phi'=\delta\Phi;\end{array}\end{equation}

then

\begin{equation} \bar{h}_{\mu\nu}\equiv h_{\mu\nu}-\frac{h}{2}\eta_{\mu\nu}+\eta_{\mu\nu}h_{f}\tag{187}\end{equation}

can be defined. Considering the transformation for the parameter $\epsilon^{\mu}$

\begin{equation} \square\epsilon_{\nu}=\partial^{\mu}\bar{h}_{\mu\nu},\tag{188}\end{equation} a gauge can be chosen, that is

\begin{equation} \partial^{\mu}\bar{h}_{\mu\nu}=0.\tag{189}\end{equation}

In this way, field equations read like

\begin{equation} \square\bar{h}_{\mu\nu}=0\tag{190}\end{equation}

\begin{equation} \square h_{f}=m^{2}h_{f}\tag{191}\end{equation}

Solutions of Eqs. (190) and (191) are plane waves

\begin{equation} \bar{h}_{\mu\nu}=A_{\mu\nu}(\overrightarrow{p})\exp(ip^{\alpha}x_{\alpha})+c.c.\tag{192}\end{equation}

\begin{equation} h_{f}=a(\overrightarrow{p})\exp(iq^{\alpha}x_{\alpha})+c.c.\tag{193}\end{equation}

where

\begin{equation} \begin{array}{ccc} k^{\alpha}\equiv(\omega,\overrightarrow{p}) & & \omega=p\equiv|\overrightarrow{p}|\\ \\q^{\alpha}\equiv(\omega_{m},\overrightarrow{p}) & & \omega_{m}=\sqrt{m^{2}+p^{2}}.\end{array}\tag{194}\end{equation}

In Eqs. (190) and (192) the equation and
the solution for the standard waves of General Relativity (Landau and Lifsits, 1999) (Misner, Thorne and Wheeler, 1973)
have been obtained, while Eqs. (191) and (193)
are respectively the equation and the solution for the massive mode related to **$f(R)$ gravity **
(see also (Bogdanos, Capozziello, De Laurentis and Nesseris,2010)).

The fact that the dispersion law for the modes of the massive field $h_{f}$ is not linear has to be emphatized. The velocity of any massles mode (arising from General Relativity) $\bar{h}_{\mu\nu}$ is the light speed $c$. On the other hand, the dispersion law (the second of Eq. (194)) for the modes of $h_{f}$ is that of a massive wave-packet. The group-velocity is

\begin{equation} \overrightarrow{v_{G}}=\frac{\overrightarrow{p}}{\omega},\tag{195}\end{equation}

that is the velocity of a massive particle with mass $m$ and momentum $\overrightarrow{p}$.

From Eqs. (194) and Eq. (195), we obtain

\begin{equation} v_{G}=\frac{\sqrt{\omega^{2}-m^{2}}}{\omega}.\tag{196}\end{equation}

The mass is

\begin{equation} m=\sqrt{(1-v_{G}^{2})}\omega.\tag{197}\end{equation}

From Eq. (187) it is

\begin{equation} h_{\mu\nu}=\bar{h}_{\mu\nu}-\frac{\bar{h}}{2}\eta_{\mu\nu}+\eta_{\mu\nu}h_{f}.\tag{198}\end{equation}

Considering the massles case one obtains

\begin{equation} \begin{array}{c} \square\epsilon^{\mu}=0\\ \\\partial_{\mu}\epsilon^{\mu}=-\frac{\bar{h}}{2}+h_{f},\end{array}\tag{199}\end{equation}

which gives the total transversality of the field. In the massive case

\begin{equation} \square\epsilon^{\mu}=m^{2}h_{f}.\tag{200}\end{equation}

In the same way, it is possible to show that it does not exist any linear relation between the tensorial field $\bar{h}_{\mu\nu}$ and the massive field $h_{f}$. After some gauge considerations it is possible to show that the gravitational waves can be written as

\begin{equation} h_{\mu\nu}(t,z)=A^{+}(t-z)e_{\mu\nu}^{(+)}+A^{\times}(t-z)e_{\mu\nu}^{(\times)}+h_{f}(t-v_{G}z)\eta_{\mu\nu}.\tag{201}\end{equation}

The term $A^{+}(t-z)e_{\mu\nu}^{(+)}+A^{\times}(t-z)e_{\mu\nu}^{(\times)}$
describes the two standard polarizations of gravitational waves which
arise from General Relativity, while the term $h_{f}(t-v_{G}z)\eta_{\mu\nu}$
is the massive field arising from **$f(R$) gravity**.
In other words, the derivative function $f'(R$) generates
a third massive mode for gravitational waves which is not
present in standard General Relativity.

## Discussion

**$f(R)$ gravity** is a class of effective theories representing a new approach to the gravitational interaction. The paradigm is that Einstein's General Relativity has to be extended in order to address several shortcomings emerging at ultraviolet and infrared scales. These are essentially due to the lack of a final, self-consistent theory of quantum gravity. From the astrophysical and cosmological viewpoints, the goal is to encompass phenomena like dark energy and dark matter under a geometric standard related to the possibility that gravitational interaction depends on the scales.

This *geometric* view, in principle, does not need the introduction of further particle ingredients and preserves all the well-posed results of General Relativity, being based on the same fundamental principles (Equivalence Principle, diffeomorphism invariance, gauge invariance, etc.).

The main criticism to this approach is that, until now, no $f(R)$ model, or any **Extended Theory of Gravity**, succeeds in addressing the whole phenomenology ranging from quantum to cosmological scales. Besides, the $f(R)$ description of dark side of the universe is substantially equivalent to that related to the hypothesis of dark material constituents. The need of one or more than one *experimentum crucis*, capable of retaining or ruling out one of the two concurring pictures, is pressing to solve the debate.

## References

- Capozziello, S. and De Laurentis, M. (2011). Extended Theories of Gravity.
*Physics Reports*509: 167. arXiv:1108.6266 [gr-qc]. doi:10.1016/j.physrep.2011.09.003. Bibcode: 2011PhR...509..167C.

- Nojiri, S. and Odintsov, S.D. (2011). Unified cosmic history in modified gravity: from F(R) theory to Lorentz non-invariant models.
*Phys.Rept.*505: 59. arXiv:1011.0544 [gr-qc]. doi:10.1016/j.physrep.2011.04.001.

- Nojiri, S. and Odintsov, S.D. (2007). Introduction to modified gravity and gravitational alternative for dark energy.
*Int. J. Geom. Meth. Mod. Phys.*4: 115. arXiv:hep-th/0601213. doi:10.1142/S0219887807001928.

- Capozziello, S. and Francaviglia, M. (2008). Extended Theories of Gravity and their Cosmological and Astrophysical Applications.
*Gen.Rel.Grav.*40: 357. arXiv:0706.1146 [astro-ph]. doi:10.1007/s10714-007-0551-y.

- De Felice, A. and Tsujikawa, S. (2010). f(R) Theories.
*Living Reviews in Relativity*13: 3. arXiv:arXiv:1002.4928 [gr-qc]. doi:10.12942/lrr-2010-3. Bibcode: 2010LRR....13....3D.

- Sotiriou, T. P. and Faraoni, V. (2010). f(R) Theories of Gravity.
*Reviews of Modern Physics*82: 451. arXiv:0805.1726 [gr-qc]. doi:10.1103/RevModPhys.82.451. Bibcode: 2010RvMP...82..451S.

- Capozziello, S. and De Laurentis, M. (2010). A Review about the Invariance Induced Gravity:Gravity and Spin from Local Conformal-Affine Symmetry.
*Foundations of Physics*40: 867. arXiv:0910.2881 [hep-th]. doi:10.1007/s10701-009-9367-9.

- Bondi, H. (1956).
*Cosmology*. Cambridge Univ. Press, Cambridge, UK. ISBN 9780521141185.

- Birrell, N.D. and Davies, P.C.W. (1982).
*Quantum Fields in Curved Space*. Cambridge Univ. Press, Cambridge, UK. ISBN 9780521278584.

- Vilkovisky, G.A. (1992). Effective action in quantum gravity.
*Class. Quant. Grav.*9: 895. doi:10.1088/0264-9381/9/4/008.

- Gasperini, M. and Veneziano, G. (1992). Boosting Away Singularities from Conformal String Backgrounds.
*Phys. Lett. B*227: 256. arXiv:hep-th/9209052. doi:10.1016/0370-2693(92)90802-B.

- Barrow, J.D and Ottewill, A.C. (1983). The stability of general relativistic cosmological theory.
*J. Phys. A*16: 2757. doi:doi:10.1088/0305-4470/16/12/022.

- Starobinsky, A.A. (1980). A new type of isotropic cosmological models without singularity.
*Phys. Lett. B*91: 99. doi:10.1016/0370-2693(80)90670-X. Bibcode: 1980PhLB...91...99S.

- Duruisseau, J.P.; Kerner, P. and Eysseric (1983). Non-Einsteinian lagrangians assuring non-singular cosmological solutions.
*Gen. Rel. Grav.*15: 797.

- La, D. and Steinhardt, P.J. (1989). Extended inflationary cosmology.
*Phys. Rev. Lett.*62: 1066. doi:http://dx.doi.org/10.1103/PhysRevLett.62.376.

- Teyssandier, P. and Tourrenc, P. (1983). The Cauchy problem for the $R+R^2$ theories of gravity without torsion.
*J. Math. Phys.*24: 2793. doi:http://dx.doi.org/10.1063/1.525659.

- Maeda, K. (1989). Towards the Einstein-Hilbert action via conformal transformation.
*Phys. Rev. D*39: 3192. doi:http://dx.doi.org/10.1103/PhysRevD.39.3159.

- Wands, D. (1994). Extended gravity theories and the Einstein-Hilbert action.
*Class. Quantum Grav.*11: 269. arXiv:gr-qc/9307034. doi:10.1088/0264-9381/11/1/025.

- Capozziello, S.; de Ritis, A.A. and Marino (1998). Recovering the Effective Cosmological Constant in Extended Gravity Theories.
*Gen. Rel. Grav.*30: 1247. arXiv:gr-qc/9804053. doi:10.1023/A:1026651129626.

- Adams, F.; Freese, K. and Guth, A. (1991). Constraints on the scalar-field potential in inflationary models.
*Phys. Rev. D*43: 965. doi:http://dx.doi.org/10.1103/PhysRevD.43.965.

- Ruzmaikina, T.V. and Ruzmaikin, A.A. (1970). Quadratic Corrections to the Lagrangian Density of the Gravitational Field and the Singularity.
*Sov. Phys. JETP*30: 372. Bibcode: 1969JETP...30..372R.

- Amendola, L.; Battaglia-Mayer, A.; Capozziello, S.; Gottlober, S.; Muller, V. and Schmidt, H.J. (1993). Generalized sixth order gravity and inflation.
*Class. Quantum Grav.*10: L43. doi:10.1088/0264-9381/10/5/001.

- Buchdahl, H.A. (1951). Uber die Variationsableitung von fundamental-invarianten beliebig hoher Ordnung.
*Acta Math*85: 63.

- Battaglia-Mayer, A. and Schmidt, H.J. (1993). The de Sitter space-time as attractor solution in eighth order gravity.
*Class. Quantum Grav.*10: 2441. doi:10.1088/0264-9381/10/11/026. Bibcode: 1993CQGra..10.2441M.

- Schmidt, H.J. (1990). Variational derivatives of arbitrarily high order and multiinflation cosmological models.
*Class. Quantum Grav.*7: 1023. doi:10.1088/0264-9381/7/6/011. Bibcode: 1990CQGra...7.1023S.

- Sokolowski, L.M. (1989). Uniqueness of the Metric Line Element in Dimensionally Reduced Theories.
*Class. Quantum Grav.*6: 2045. doi:10.1088/0264-9381/6/1/006.

- Ferraris, M.; Francaviglia, M. and Magnano, G. (1988). Do non–linear metric theories of gravitation really exist?.
*Class. Quantum Grav.*5: L95. doi:10.1088/0264-9381/5/6/002.

- Magnano, G. and Sokolowski, L.M. (1994). On physical equivalence between nonlinear gravity theories and a general relativistic selfgravitating scalar field.
*Phys.Rev. D*50: 5039. arXiv:gr-qc/9312008. doi:10.1103/PhysRevD.50.5039.

- Jakubiec, A. and Kijowski, J. (1988). On Theories of Gravitation With Nonlinear Lagrangians.
*Phys.Rev. D*37: 1406. doi:10.1103/PhysRevD.37.1406.

- Capozziello, S. and Faraoni, V. (2010).
*Beyond Einstein gravity: A Survey of gravitational theories for cosmology and astrophysics*. Fundamental Theories of Physics.**170**. Springer. ISBN 978-94-007-0164-9.

- Capozziello, S. and De Laurentis, M. (2010).
*Invariance Principles and Extended Gravity: Theories and Probes*. Nova Science Publishers, Inc.. ISBN 978-1-61668-500-3.

- Capozziello, S.; De Laurentis, M. and Faraoni, V. (2009). A bird's eye view of f(R)-gravity.
*The Open Astronomy Journal*2: 1874. arXiv:arXiv:0909.4672v2 [gr-qc]. doi:10.2174/ 1874381101003010049.

- Faraoni, V. (2008).
*f(R) gravity: Successes and challenges*. arXiv:0810.2602v1 [gr-qc].

- Einstein, A. (1925).
*Sitzung-ber. Preuss. Akad.Wiss.*414.

- Buchdahl, H.A. (1979). Quadratic Lagrangians and Palatini's device.
*J. Phys. A*12 (8): 1229. doi:doi:10.1088/0305-4470/12/8/017. Bibcode: 1979JPhA...12.1229B.

- Olmo, G.J. (2011). Palatini Approach to Modified Gravity: f(R) Theories and Beyond.
*Int. J. Mod. Phys. D*20: 413. arXiv:arXiv:1101.3864v1 [gr-qc]. doi:10.1142/S0218271811018925.

- Ferraris, M.; Francaviglia, I. and Volovich (1994). The Universality of Einstein Equations.
*Class. Quantum Grav.*11: 1505. arXiv:gr-qc/9303007v2. doi:10.1088/0264-9381/11/6/015.

- Vollick, D.N. (2003). 1/R Curvature Corrections as the Source of the Cosmological Acceleration.
*Phys.Rev.D*68: 063510. arXiv:astro-ph/0306630. doi:10.1103/PhysRevD.68.063510.

- Li, B. and Chu, M.C. (2006). CMB and Matter Power Spectra of Early f(R) Cosmology in Palatini Formalism.
*Phys.Rev.D*74: 104010. arXiv:astro-ph/0610486v1. doi:10.1103/PhysRevD.74.104010.

- Li, B.; Chan, K.C. and Chu, M.C. (2007). Constraints on f(R) Cosmology in the Palatini Formalism.
*Phys.Rev.D*76: 024002. arXiv:astro-ph/0610794v2. doi:10.1103/PhysRevD.76.024002.

- Capozziello, S.; De Laurentis, M.; Francaviglia, M. and Mercadante, S. (2009). From Dark Energy and Dark Matter to Dark Metric.
*Foundations of Physics*39: 1161. arXiv:0805.3642v4. doi:10.1007/s10701-009-9332-7.

- Allemandi, G.; Capone, M.; Capozziello, S. and Francaviglia, M. (2006). Conformal aspects of Palatini approach in Extended Theories of Gravity.
*Gen. Rel. Grav.*38: 33. arXiv:hep-th/0409198v1. doi:10.1007/s10714-005-0208-7.

- O'Hanlon, J. (1972).
*Phys. Rev. Let.*29: 137.

- Wald, R.M. (1984).
*General Relativity*. Chicago University Press. ISBN 978-0226870335.

- Amendola, L.; Enqvist, K. and Koivisto, T. (2011). Unifying Einstein and Palatini gravities.
*Phys. Rev. D*83: 044016. arXiv:1010.4776 [gr-qc]. doi:10.1103/PhysRevD.83.044016.

- Li, B.; Mota, D.F. and Shaw, D.J. (2008). Microscopic and Macroscopic Behaviors of Palatini Modifed Gravity Theories.
*Phys. Rev. D*78: 064018. arXiv:0805.3428 [gr-qc]. doi:10.1103/PhysRevD.78.064018.

- Olmo, G.J.; Sanchis-Alepuz, H. and Tripathi, S. (2009). Dynamical Aspects of Generalized Palatini Theories of Gravity.
*Phys. Rev. D*80: 024013. arXiv:0907.2787 [gr-qc]. doi:10.1103/PhysRevD.80.024013.

- Barragan, C. and Olmo, G. (2010). Isotropic and Anisotropic. Bouncing Cosmologies in Palatini Gravity.
*Phys. Rev. D*82: 084015. arXiv:1005.4136 [gr-qc].

- Capozziello, S. and Garattini, R. (2007). The cosmological constant as an eigenvalue of f(R)-gravity Hamiltonian constraint.
*Class. Quantum Grav.*24: 1627. arXiv:gr-qc/0702075. doi:10.1088/0264-9381/24/6/016.

- Deruelle, N.; Sendouda, Y. and Youssef, A. (2009). Various Hamiltonian formulations of f(ℛ) gravity and their canonical relationships.
*Phys. Rev. D*80: 084032. arXiv:0906.4983 [gr-qc]. doi:10.1103/PhysRevD.80.084032.

- Arnowitt, R.; Deser, S. and Misner, C. W. (1962).
*Gravitation: An Introduction to Current Research*. Edit by L. Witten, John Wiley \& Sons, Inc., New York.

- Garattini, R. (2006). The Cosmological Constant as an Eigenvalue of a Sturm-Liouville Problem and its Renormalization.
*J. Phys. A*39: 6393. arXiv:gr-qc/0510061. doi:10.1088/0305-4470/39/21/S33.

- Elizalde, E.; Odintsov, S.; Romeo; Bytsenko and Zerbini (1994).
*Zeta regularization techniques with applications*. Edit by L. Witten, John Wiley \& Sons, Inc., New York. ISBN 978-981-02-1441-8.

- De Witt, B.S. (1967). The Cosmological Constant as an Eigenvalue of a Sturm-Liouville Problem and its Renormalization.
*Phys. Rev.*160: 1113. arXiv:gr-qc/0510061.

- Bronnikov, K.A. (2001). Spherically symmetric false vacuum: No-go theorems and global structure.
*Phys. Rev. D*64: 064013. arXiv:gr-qc/0104092. doi:10.1103/PhysRevD.64.064013.

- Bronnikov, K.A. (2002). Scalar-tensor gravity and conformal continuations.
*J. Math. Phys.*43: 6096. arXiv:gr-qc/0204001. doi:10.1063/1.1519667.

- Bronnikov, K.A. and Shikin, G.N. (2002). Spherically symmetric scalar vacuum: no-go theorems, black holes and solitons.
*Grav.Cosmol.*8: 107. arXiv:gr-qc/0109027.

- Synge, J.L. (1955).
*Relativity: The General Theory*. Amsterdam: North Holland.

- Lorentz, M.A. (1937).
*Collected Papers*. 5.**363**. Nijhoff, The Hague.

- Sotiriou, T.; Liberati, S. and Faraoni, V. (2008). Theory of gravitation theories: a no-progress report.
*Int.J.Mod.Phys.D*17: 393. arXiv:0707.2748 [gr-qc]. doi:10.1142/S0218271808012097.

- Capozziello, S. (2002). Curvature and quintessence.
*Int.J.Mod.Phys.D*11: 483. arXiv:gr-qc/0201033v1. doi:10.1142/S0218271802002025.

- Capozziello, S.; De Laurentis, M. and Lambiase, G. (2012). Cosmic relic abundance and f(R) gravity.
*Phys. Lett. B*715 (1): 1. arXiv:1201.2071 [gr-qc]. doi:10.1016/j.physletb.2012.07.007.

- Starobinsky, A.A. (1979). Spectrum of relict gravitational radiation and the early state of the universe.
*JETP Lett.*30: 682. doi:10.1088/0305-4470/16/12/022. Bibcode: 1979ZhPmR..30..719S.

- Nojiri, S. and Odintsov, S.D. (2003). Modified gravity with negative and positive powers of the curvature: unification of the inflation and of the cosmic acceleration.
*Phys.Rev.D*68: 123512. arXiv:hep-th/0307288v4. doi:10.1103/PhysRevD.68.123512.

- Nojiri, S.; Odintsov, S.D. and Tsujikawa, S. (2005). Properties of singularities in (phantom) dark energy universe.
*Phys.Rev.D*71: 063004. arXiv:hep-th/0501025. doi:10.1103/PhysRevD.71.063004.

- Bamba, K. and Odintsov, S.D. (2008). Inflationary cosmology and the late-time accelerated expansion of the universe in non-minimal Yang-Mills-F(R) gravity and non-minimal vector-F(R) gravity.
*Phys.Rev. D*77: 123532. arXiv:arXiv:0803.3384 [hep-th]. doi:10.1103/PhysRevD.77.123532.

- Capozziello, S. and Amendola, L. (1993). The Phase space view of inflation. 2: Fourth order models.
*Int.J.Mod.Phys. D*1: 615. arXiv:arXiv:0803.3384 [hep-th]. doi:10.1142/S0218271892000318.

- Carloni, S.; Dunsby, P.K.S.; Capozziello, S. and Troisi, A. (2005). Cosmological dynamics of R^n gravity.
*Class.Quant.Grav.*22: 4839. doi:10.1088/0264-9381/22/22/011.

- Carloni, S.; Troisi, A. and Dunsby, P.K.S. (2009). Some remarks on the dynamical systems approach to fourth order gravity.
*Gen.Rel.Grav.*41: 1757. doi:10.1088/0264-9381/22/22/011.

- Abdelwahab, M.; Goswami, R. and Dunsby, P.K.S. (2009). Cosmological dynamics of fourth order gravity: A compact view.
*Phys.Rev. D*85: 083511. doi:10.1103/PhysRevD.85.083511.

- Wainwright, J. and Ellis, G.F.R. (1997).
*Dynamical System in Cosmology*. Eds.Cambridge University Press,Cambridge .

- Brandenberger, R.H. (1992). Theory of cosmological perturbations. Part 1. Classical perturbations. Part 2. Quantum theory of perturbations. Part 3. Extensions.
*Phys.\ Rept.*215: 203.

- Dodelson, S. (2003).
*Gauge Invariant Cosmological Perturbations*. Academic Press .

- Bardeen, J.M. (1932). Gauge Invariant Cosmological Perturbations.
*Phys. Rev. D*22: 1882. doi:10.1103/PhysRevD.22.1882.

- Pogosian, L. and Silvestri, A. (2008). The pattern of growth in viable f(R) cosmologies.
*Phys.\ Rev.\ D*77: 023503.

- de la Cruz-Dombriz, A.; Dobado, A. and Maroto, A.L. (2008). On the evolution of density perturbations in f(R) theories of gravity.
*Phys.\ Rev.\ D*77: 123515.

- Carloni, S.; Dunsby, P.K.S. and Troisi, A. (2008). The evolution of density perturbations in f(R) gravity.
*Phys. Rev. D*77: 024024. arXiv:arXiv:0707.0106 [gr-qc]. doi:10.1103/PhysRevD.77.024024.

- Wands, D. and Slozar, A. (2009). Scale-dependent bias from primordial non-Gaussianity in general relativity.
*Phys. Rev. D*79: 123507. arXiv:0902.1084 [astro-ph.CO]. doi:10.1103/PhysRevD.79.123507.

- Bruni, D.; Crittenden, R.; Koyama, K.; Maartens, R.; Pitrou, C. and Wands, D. (2012). Disentangling non-Gaussianity, bias and GR effects in the galaxy distribution.
*Phys. Rev. D*85: 041301. arXiv:1106.3999 [astro-ph.CO]. doi:10.1103/PhysRevD.85.041301.

- Clifton, T.; Ferreira, P.G.; Padilla, A. and Skordis, C. (2012). Modified Gravity and Cosmology.
*Phys.Rept.*513: 1-189. arXiv:1106.2476 [astro-ph.CO]. doi:10.1016/j.physrep.2012.01.001.

- Amendola, L. and Tsujikawa, S. (2010).
*Dark Energy: Theory and Observations*. Cambridge University Press, Cambrige .

- Ade, P.~A.~R. and et al. (2014). Planck 2013 results. XVI. Cosmological parameters.
*Astron.\ Astrophys.*571: A16. arXiv:1303.5076 [astro-ph.CO]. doi:10.1016/j.physrep.2012.01.001.

- Oort, J.H. (1932). The Force Exerted by the Stellar System in the Direction Perpendicular to the Galactic Plane and Some Related Problems.
*Bull. Astr. Neth.*6: 249. Bibcode: 1932BAN.....6..249O.

- Oort, J.H. (1960). Note on the Determination of Kz and on the Mass Density Near the Sun.
*Bull. Astr. Neth.*15: 45. Bibcode: 1960BAN....15...45O.

- Zwicky, F. (1933). Die Rotverschieb ung von extragalaktischen Nebeln.
*Helv. Phys. Acta*6: 110.

- Will, C.M. (1993).
*Theory and Experiment in Gravitational Physics*. 2nd ed. Cambridge University Press, Cambridge, UK.

- Milgrom, M. (1983). A modification of the Newtonian dynamics as a possible alternative to the hidden mass hypothesis.
*The Astrophysical Journal*270: 389. doi:10.1086/161130. Bibcode: 1983ApJ...270..365M.

- Milgrom, M. and Sanders, R.H. (1983).
*The Astrophysical Journal Letters*315: 493.

- Capozziello, S. and De Laurentis, M. (2012). The dark matter problem from f(R) gravity viewpoint.
*Annalen der Physik*524: 545. doi:10.1002/andp.201200109.

- Capozziello, S. and Tsujikawa, S. (2008). Solar system and equivalence principle constraints on f(R) gravity by chameleon approach.
*Phys. Rev. D*77 (1): 1. arXiv:0712.2268 [gr-qc]. doi:10.1103/PhysRevD.77.107501.

- Khoury, J. and Weltman, A. (2004). Chameleon Cosmology.
*Phys. Rev. D*69 (1): 044026. arXiv:astro-ph/0309411. doi:10.1103/PhysRevD.69.044026.

- Maeda, K.I. (1989). Towards the Einstein-Hilbert Action via Conformal Transformation.
*Phys. Rev. D*39: 3159. doi:10.1103/PhysRevD.39.3159.

- Faulkner, T.; Tegmark, M.; Bunn, E.F. and Mao, Y. (2007). Constraining f(R) gravity as a scalar tensor theory.
*Phys. Rev. D*76: 063505. arXiv:astro-ph/0612569. doi:10.1103/PhysRevD.76.063505.

- Capozziello, S.; Corda, C. and De Laurentis, M. (2008). Massive gravitational waves from f(R) theories of gravity: Potential detection with LISA.
*Physics Letters B*669 (5): 255–259. arXiv:0812.2272. doi:10.1016/j.physletb.2008.10.001. Bibcode: 2008PhLB..669..255C.

- Landau, L. and Lifsits, E. (1999).
*Teoria dei campi*. Editori riuniti edition III . ISBN 978-8835909002.

- Misner, C.W.; Thorne, K.S. and Wheeler, J.A. (1973).
*Gravitation*. W.H.Freeman and Company . ISBN 978-0716703440.

- Bogdanos, C.; Capozziello, S.; De Laurentis, S. and Nesseris (2010). Massive, massless and ghost modes of gravitational waves from higher-order gravity.
*Astroparticle Physics*34: 236. arXiv:0911.3094 [gr-qc]. doi:10.1016/j.astropartphys.2010.08.001.