Observational constraints to a phenomenological f\left(R,\nabla R\right)-model

# Observational constraints to a phenomenological f(R,∇R)-model

R. R. Cuzinatto , C. A. M. de Melo , L. G. Medeiros , P. J. Pompeia cuzinatto@gmail.comcassius.anderson@gmail.comleogmedeiros@gmail.compedropjp@ifi.cta.br
###### Abstract

This paper analyses the cosmological consequences of a modified theory of gravity whose action integral is built from a linear combination of the Ricci scalar and a quadratic term in the covariant derivative of . The resulting Friedmann equations are of the fifth-order in the Hubble function. These equations are solved numerically for a flat space section geometry and pressureless matter. The cosmological parameters of the higher-order model are fit using SN Ia data and X-ray gas mass fraction in galaxy clusters. The best-fit present-day values for the deceleration parameter, jerk and snap are given. The coupling constant  of the model is not univocally determined by the data fit, but partially constrained by it. Density parameter is also determined and shows weak correlation with the other parameters. The model allows for two possible future scenarios: there may be either a premature Big Rip or a Rebouncing event depending on the set of values in the space of parameters. The analysis towards the past performed with the best-fit parameters shows that the model is not able to accommodate a matter-dominated stage required to the formation of structure.

Instituto de Ciência e Tecnologia, Universidade Federal de Alfenas.

Rod. José Aurélio Vilela (BR 267), Km 533, n11999, CEP 37701-970, Poços de Caldas, MG, Brazil.

Rua Bento Teobaldo Ferraz 271 Bloco II, P.O. Box 70532-2, CEP 01156-970, São Paulo, SP, Brazil.

Escola de Ciência e Tecnologia, Universidade Federal do Rio Grande do Norte.

Campus Universitário, s/n - Lagoa Nova, CEP 59078-970, Natal, RN, Brazil.

Inst. de Fomento e Coordenação Industrial, Departamento de Ciência e Tecnologia Aeroespacial.

Praça Mal. Eduardo Gomes 50, CEP 12228-901, São José dos Campos, SP, Brazil.

Instituto Tecnológico de Aeronaútica, Departamento de Ciência e Tecnologia Aeroespacial.

Praça Mal. Eduardo Gomes 50, CEP 12228-900, São José dos Campos, SP, Brazil.

Keywords: Gauge theory, Cosmic acceleration, Higher order gravity, Cosmology.
PACS:98.80.-k, 11.15.-q.

## 1 Introduction

The cosmological constant is the simplest solution to the problem of the present-day acceleration of the universe. It is consistent with all the available observational tests, from the constraints imposed by CMB maps [1, 2] to BAO picks [3, 4] to SN Ia constraints [5, 6, 7, 8, 9] and the estimates of fraction of baryons in galaxy clusters through X-ray detection [10, 11]. Nevertheless, there is an uncomfortable feeling from the fact that  is not satisfactory explained as a vacuum state of some field theory: the energy scale observed for and predicted by quantum field theory badly disagree [12]. This is one of the reasons why alternative explanations for the present-day cosmic acceleration have been put out since 1998, when it was first detected [13, 14].

These proposals can be classified into two broad categories, each one related to a modification of a different side of the Einstein field equations of General Relativity [15].

The modified matter approach introduces a negative pressure (that counteracts gravity) through the energy momentum tensor  on the right-hand side of Einstein equations. Three types of models constitute this approach, namely quintessence models [16, 17, 18, 19, 20], k-essence solutions [21, 22] and perfect fluid models (such as the Chaplygin gas [23, 24]). The quintessence and k-essence models are based on scalar fields whose potential (in the former case) or the kinetic energy (in the later one) drive the accelerated expansion. The perfect fluid models use exotic equations of state (for the pressure as a function of the energy density) to generate negative pressure. The weakness of k-essence and quintessence proposals lies in the fact that they require a scalar field to be dominant today, when the energy scale is low. This requirement demands the mass of such a field  to be extremely small. One would then expect long range interaction of this field with ordinary matter, which is highly constrained by local gravity experiments. The difficulty with the perfect fluid models is the lack of an interpretation to the exotic equations of state based on ordinary physics.

The second approach, called modified gravity, affects the left hand side of Einstein equations [25, 26]. The subtypes of models in this category include braneworld scenarios [27, 28], scalar-tensor theories [29, 30], non-riemannian geometries [31, 32, 33],  gravity [34, 35, 36, 37, 38, 39] and  theories [40, 41, 42]. In the braneworld models there exists an additional dimension through which the gravitational interaction occurs and by which the possibility of a cosmic acceleration is accommodated. Scalar-tensor theories couple the scalar curvature  with a scalar field  in the action , allowing enough degrees of freedom for an accelerated solution. Non-riemannian geometries are another alternative way to explore the geometrization of the gravitational interaction, but allowing new degrees of freedom describing different geometric structures than curvature. The  gravity corresponds to models resulting from non-linear Lagrangian densities dependent on ;  theories build their cosmological models on manifolds whose connection is equipped with torsion instead of curvature.

All these mechanisms of generating the present-day cosmic dynamics are generally called dark energy models.

Our aim with the present paper is to contribute with the study of dark energy in the branch of modified gravity. We explore the cosmological consequences of a particular model in the context of the gravity testing the model against some of the observational data available. gravity is a theory based on an action integral whose kernel is a function of the Ricci scalar and terms involving derivatives of . Ref. [43] is perhaps an early example of what we call gravity. A more recent example is the (non-local) higher derivative cosmology presented in Refs. [44, 45, 46]. These works intend to eliminate the Big Bang singularity through a bounce produced by a theory based on an action containing the regular Einstein-Hilbert term plus a String-spired term. This last ingredient contains an infinite number of derivatives of , a requirement needed to keep this modified gravity theory ghost-free. References [47, 48] are -gravity prototypes concerned not with the past-incomplete inflationary space-time but with the cosmological constant problem. In this regard, they try to address a similar issue that interests us in the present work: the cosmological acceleration nowadays.

Few years ago, we proposed a gauge theory based on the second order derivative of the gauge potential  following Utiyama’s procedure [49]. This led to a second field strength  (depending on the field strength ) and allowed us to obtain massive gauge fields among other applications, such as the study of gauge invariance of the theory under and groups when . In a subsequent work, gravity was interpreted as a second order gauge theory under the Lorentz group, and all possible Lagrangians with linear and quadratic combinations of and  – or the Riemann tensor, its derivatives and their possible contractions, from the geometrical point of view – were built and classified within that context [50]. In a third paper we selected one of those quadratic Lagrangians for investigating the cosmological consequences of a theory that takes into account - and -terms [51]. That Lagrangian appears in the action  bellow – Eq. (1) – and can be understood as a particular example of an  gravity:

 S=∫d4x√−g(R2χ+β8χ∇μR∇μR−LM) , (1)

where (once we take ), is the coupling constant,  is the covariant derivative, and is the Lagrangian of the matter fields. As usual, .111The authors of Ref.[43] use the same action (1) proposed here but their argument is based on different motivations.

An important feature of is the presence of ghosts (see [46] and references therein). However, here this is not a serious issue because (1) should be considered as an effective action of a fundamental theory of gravity, this last one a ghost-free theory. There are some consistent proposals of extended theories of gravity in the literature. A good example is the String-inspired non-local higher derivative gravitation model presented in Ref. [45]. The ghost-free action of that model is:

 S=∫d4x√−g(M2P2R+12RF(□/M2∗)R) (2)

where  is the Planck mass (, is the gravitational constant), is the mass scale where the higher-derivative terms become important and

 F(□/M2∗)=∞∑n=1fn□n . (3)

In this context, our higher-order model must be treated as a low energy limiting case of (2) where the only relevant contribution of Eq. (3) for the present-day universe is the first term of the series, i.e. , . In this approximation, Eq. (2) is equivalent to the action integral (1) of our -cosmology up to a surface term. The equivalence of both models leads to a relation between the coupling constant and the parameter , namely

 f1=β32πG . (4)

Taking the variation of (1) with respect to the metric tensor (cf. the metric formalism) gives the field equations

 Rμν−12gμνR+βHμν=χTμν , (5)

with

 Hμν =∇μ∇ν[□R]+12∇μR∇νR −Rμν□R−gμν□[□R]−14gμν∇ρR∇ρR , (6)

where and is the Ricci tensor.

A four dimensional spacetime with homogeneous and isotropic space section is described by the Friedmann-Lemaître-Robertson-Walker (FLRW) line element

 ds2=c2dt2−a2(t)[dr21−κr2+r2(dθ2+sin2θdϕ2)] , (7)

where the curvature parameter  assumes the values and is the scale parameter. Substituting the metric tensor from Eq.(7) along with the energy-momentum tensor of a perfect fluid in co-moving coordinates,

 Tμν=(ρ+P)δ0μδ0ν−Pgμν (8)

( is the density and is the pressure of the components fulling the universe), into the gravity field equations Eq.(5) leads to the two higher order Friedmann equations:

 H2+β3F1=8πG3ρ ; (9)
 3H2+2˙H+βF2=−8πGP , (10)

where

 F1 =18HH(4)+108H2...H−18˙H...H+9¨H2+90H3¨H +216H˙H¨H−72˙H3+288(H˙H)2−216H4˙H (11)
 F2 =6H(5)+54HH(4)+138H2...H+126˙H...H +81¨H2+18H3¨H+498H˙H¨H+120˙H3−216H4˙H (12)

which are written in terms of the Hubble function ,

 H≡˙aa , (13)

and where we have set as is largely favored by the observational data. The dot (such as in ) means derivative with respect to the cosmic time and for any arbitrary function of time ,

 f(n)≡dnfdtn ,

i.e., is the fifth-order derivative of the Hubble function with respect to the cosmic time.

Eqs. (9) and (10) reduce to the usual Friedmann equations if . From this system of modified Friedmann equations one can obtain the equation expressing the covariant conservation of the energy momentum tensor, , which reads

 ˙ρ+3H(ρ+P)=0 . (14)

In Ref.[51], we considered simultaneously the set formed by Eq. (14) and a combination of Eqs.(9) and (10), namely

 ˙H+β2(F2−F1)=−4πG(ρ+P) . (15)

The system of Eqs.(14)-(15) was solved perturbatively with respect to the coupling constant for pressureless matter. Using the observational data available, we have determined the free parameters of our modified theory of gravity and the scale parameter as a function of the time. The plot of indicates that the decelerated regime peculiar to the universe dominated by matter evolves naturally to an accelerated scenario at recent times.

This makes it promising to perform a more careful analyzes of our  model making a complete maximum likelihood fit to all parameters.

In this paper we will solve the system of Eqs.(14)-(9) exactly, without using any approximation technique. Section 2 presents the script to obtain the complete solution to the higher order cosmological equations. Another improvement we shall do in the present work (when comparing with Ref.[51]) is to use the observational data available to fully test our model. As far as we know, this is the first time this is done for a higher derivative theory of gravity. In Section 3 we construct the likelihood function to be maximized in the process of fitting the model parameters using data from SN Ia observations [8] and from the measurements of the mass fraction of gas in clusters of galaxies [11]. The results of this fit are given in Section 4, where we show the present-day numerical values for the deceleration parameter [52],

 q=−¨aaH2 (16)

the jerk function ,

 j=...aaH3 , (17)

and also of the snap ,

 s=a(4)aH4 . (18)

These last two functions are given in terms of higher-order time derivatives of the scale factor and it is only natural that they show up here. In Section 4.2 we perform an investigation of the dynamical evolution of the model, comparing its behaviour towards the past with the standard solutions for single component universes. Moreover, we show that our -model accommodates a future Big Rip or a Rebouncing scenario depending on the values of the coupling constant . The Rebouncing event is the point where the universe attains its maximum size of expansion and from which it begins to shrink. Our final comments about the features of our model and its validity are given in Section 5.

## 2 Higher order cosmological model

The solution of our higher-order model shall be obtained from the integration of the pair of Eqs.(9) and (14). The definitions

 B≡H40β , (19)

and

 Ωm(t)≡χρ3H20=ρρc , (20)

enables us to express Eq. (9) as

 H2H20+B3H60F1=Ωm(t) (21)

and the conservation of energy-momentum Eq.(14) as

 ˙Ωm+3HΩm=0 , (22)

 Ωm=Ωm0(a0a)3 . (23)

It is standard to take . Moreover, substituting the solution (23) into Eq. (21), leads to a fifth-order differential equation in :

 1H20˙a2a4+BH60F1=Ωm0a3 (24)

with

 F1 =6a(5)˙aa4−6a(4)¨aa4+12a(4)˙a2a3 +3 ...a2a4+18...a¨a˙aa3−78...a˙a3a2 −6¨a3a3−39¨a2˙a2a2+78¨a˙a4a−180˙a6 (25)

In order to solve Eq. (24), we define the non-dimensional variables

 H ≡˙aH0=aHH0 ; (26) Q ≡˙HH0=−aH2H20q ; (27) J ≡˙QH0=aH3H30j ; (28) S ≡˙JH0=aH4H40s , (29)

which are related to the Hubble function, deceleration parameter, jerk and snap. In fact, when these variables are calculated at the present time, , they give precisely the well known cosmological parameters: , ,  and . In the set of variables (26)-(29), one variable is essentially the derivative of the previous one, therefore they can be used to turn the fifth-order differential equation Eq.(24) into five first order differential equations to be solved simultaneously. Following the Method of Characteristics, we choose the scale factor  as the independent variable of the model (instead of the cosmic time ), obtaining the system

 ⎧⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎨⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪⎩dHdz=−1(1+z)2QHdQdz=−1(1+z)2JHdJdz=−1(1+z)2SHdSdz=−1(1+z)2F3(z;H,Q,J,S) , (30)

with

 F3(z;H,Q,J,S) =−2S(1+z)+SQH2 +Q3(1+z)H2+132Q2(1+z)2−13QH2(1+z)3 +30H4(1+z)4+Ωm06B(1+z)H2−16B , (31)

where we have used the relationship between the redshift and the scale factor,

 (1+z)=1a . (32)

The system of Eq.(30) is solved numerically, which is indeed a far more simpler task than dealing directly with Eq. (24). This leads to a function (essentially equal to ) for each set  of cosmological parameters. This function is used to obtain an interpolated function for the luminosity distance  of our model – defined bellow – which is to be fitted to the data from SN Ia observations and from -gas measurements.

## 3 Fitting the model to the SN Ia and f-gas data

### 3.1 SN Ia data

In order to relate the proposed model to the observational data, one has to build the luminosity distance [15],

 (33)

Our task is to solve the system Eq.(30) for obtaining , then to integrate numerically the non-dimensional luminosity distance,

 dh(z;[M])=H0cdL(z;[M]) ,(c=1) (34)

for each set of cosmological parameters .

Function  appears in the definition of the distance modulus,

 μr=5log10dh+5log10(cH0d0) , (35)

where  is the reference distance of pc. Observational data give of the supernovae and their corresponding redshift .

According to the likelihood method the probability distribution of the data

 LSNIa([M]|{zi,μri})=NSNIae−χ2SNIa([M])2 (36)

( is a normalization constant) is the one that minimizes

 χ2SNIa([M])=∑i(μri−5logdh(zi,[M])−5log10(cH0d0))2σ2μi+σ2ext , (37)

where  is the uncertainty related to the distance modulus of the -th supernova. The data set used in our calculations is composed with the 557 supernovae of the Union 2 compilation [8]. Uncertainty  is given by

 σ2ext=σ2μv+σ2lens , (38)

i.e., it is the combination of the error in the magnitude due to host galaxy peculiar velocities (of order of ),

 σμv=2.172zσvc\ \ \ with\ \ \ σv=300kms−1 , (39)

and the uncertainty from gravitational lensing

 σlens=0.093z (40)

which is progressively important at increasing SN redshifts. Further detail can be found in Refs. [53, 54], for example.

### 3.2 Galaxy cluster X-ray data

Measurements of the fraction of the mass of gas in galaxy clusters in the X-ray band of the spectrum constrain efficiently the ratio and exhibit a lower dependence on the values of other cosmological parameters [10, 11].222 is the non-dimensional baryon density parameter. According to Ref.[11], the mass fraction of gas relates to  as

 fΛCDMgas(z)=KAγb(z)1+s(z)(Ωb0Ωm0)[dΛCDMA(z)dA(z)]1,5 , (41)

where is the angular distance for ,

 dA(z)=dL(z)(1+z)2=1(1+z)1H0z∫0d¯z(1+¯z)H(¯z;[M]) . (42)

The term is the angular distance for the standard CDM model, for which , and Mpc, i.e.

 dΛCDMA(z)=1(1+z)1H0z∫0d¯z[Ωm0(1+¯z)3+ΩΛ]1/2. (43)

The other factors in Eq. (41) are defined bellow:

• is a constant due to the calibration of the instruments observing in X-ray. Following Ref. [11], we take .

• Constant measures the angular correction on the path of an light ray propagating in the cluster due to the difference between CDM model and the model that we want to fit to the data. Ref. [11] argues that for relatively low red-shifts (as those involved here) this factor can be neglected, i.e. taken as .

• accounts for the influence of the non-thermal pressure within the clusters. According to Ref.[11], one may set . Our choice is .

• is the bias factor related to the thermodynamical history of the inter-cluster gas. Once again due to arguments given in Ref.[11], this factor can be modeled as  where and . We set .

• Parameter describes the fraction of baryonic mass in the stars.333Notice that here is not the snap, defined in Eq. (18). Even at risk of confusion, we decided to maintain the notation used in Ref. [11]. Ref. [11] uses and . Here

 h70=H0H700  , (44)

and we will set . Analogously to Eq.(44), we define

 h=H0(100kms−1Mpc−1) , (45)

that will be required bellow.

• The baryon parameter is set to

 Ωb0h2=0.0214±0.002 (46)

following the constraints imposed by the primordial nucleosynthesis on the ratio [55].

The uncertainties in the set of parameters are propagated to build the uncertainty of . It is worth mentioning that the quantities  have Gaussian uncertainties. On the other hand, the parameters  present a rectangular distribution, e.g.  has equal probability of assuming values in the interval . The uncertainty of the parameters obeying the rectangular distribution is estimated as the half of their variation interval divided by . Maybe it is useful to remark that in the case of , we adopted . Moreover, .

Ref. [11] provides data for the fraction of gas and its uncertainty of 42 galaxy clusters with redshift . These data can be used to determine the free parameter of our model by minimizing

 χ2gas([M])=∑i(figas−fΛCDMgas(z,[M]))2σ2fi (47)

Notice that the dependence of on the set of cosmological parameters occurs through the function given by Eq.(41), i.e. ultimately through – Eq. 42.

Function  is related to the distribution

 Lgas([M])=Ngase−χ2gas([M])2 (48)

which should be maximized.

### 3.3 Combining SN Ia and f-gas data sets

The majority of papers that present cosmological models tested through SN Ia data fit – see e.g. Refs. [20, 53] – perform an analitical marginalization with respect to . On the other hand, constraint analyses using -gas data adopt a fixed value for . Our basic reference on -gas treatment, namely Ref.[11], choose Mpc. As we will consider SN Ia and -gas data sets in the same footing, it is convenient to adopt Mpc as a fixed quantity in our calculations (and not the most up-to-date value of Mpc found in [56]).

We want to maximize the probability distributions given in Eqs.(36) and (48) simultaneously, which is done by maximixing

 L([M])=LSNIa([M])Lgas([M]) , (49)

or minimizing

 χ2([M])=χ2SNIa([M])+χ2gas([M]) . (50)

## 4 Results

### 4.1 Estimation of the free parameters of the f(R,∇R)-model

The parameters to be ajusted are of two types. The kinematic parameters are ;444We emphasize that are precisely the values of calculated at the present time . the cosmological parameters are . Maximization of the likelihood function – Eq. (49) – determines all kinematic parameters and one cosmological parameter, namely . Parameter is not fixed by statistical treatment of the data, and the reason for that is evident in the plot of Fig. 1.

The curve of is not Gaussian but rather shows a flat region at high values of . Indeed, values of greater than  are almost equally probable and more favorable than the . This shows it is not possible to determine by maximizing its likelihood function. Instead, we choose four values of to proceed the analysis leading to the determination of the free parameters of the model. These values are , so that we have two values of in the beginning of the plateau () and two of them far from the beginning of the plateau (). These choices are only partially arbitrary: we have done preliminary calculations with other values of  () and the results for the best-fit calculations of are almost the same as those obtained by using .

The standard statistical analysis leads to the best-fit values of the parameters shown in Table 1.

It is clear from Table 1 that our higher order cosmological model presents an acceleration at present time: the deceleration parameter is negative for all the values of the coupling constant . This is a robust result indeed once we have got negative values for with confidence level.

Both SN Ia and -gas data sets give information about astrophysical objects at relatively low redshifts (). That is the reason for the high values of the uncertainties on the parameter  and . Positive values of  mean that the tendency is the increasing of the rate of cosmic acceleration. Because the jerk parameter  does not exhibit the same sign for all values of , it is not possible to argue in favor of an overall behavior for the cosmic acceleration today. This statement is also applicable to the snap obtained with . However, the uncertainties in the value of for  do indicate that it is positive within 68% CL. The negative value of  for  is also guaranteed within the same CL. This difference in the sign of  for and  leads to distinct dynamics for the scale factor towards the future, as we shall discuss in Section 4.2.

The best-fits for  are statistically the same for all the chosen values of . This is linked to the fact that the -gas data produce a sharp determination of this parameter. Notice that the value  of our model is compatible with that of CDM model as reported in Ref.[8] within 95% CL. Inspection of Table 1 indicates that the Gaussian-like distributions centered in each parameter of the set for  overlap the distributions centered in for ; this overlap of distribution curves occurs within 68% CL of each distribution. Moreover, this overlap of distributions is greater than the one observed in the case of the sets for  and , which occurs within 95% CL. Hence, the greater the value of , the greater the overlap of distribution and the less important is our choice of the value of . This emphasizes what was said before about the choice of values for .

The confidence regions (CR) in the and planes for each one of the four values of are shown in Figs. 2-5.

The shape of the confidence regions for  indicates that the deceleration parameter  and the jerk  are correlated parameters; this occurs for all the values of . The four plots point to an accelerated universe at recent times within 95% CL, but the change in the acceleration rate takes positive and negative values with equal probability.

Conversely, the shape of the plots of  in Figs. 25 indicate that  and  are weakly correlated parameters. Notice that the difference between the maximum and the minimum values of  in the range of the plots is only 0.05 with 95% CL: this confirms that the parameter  is sharply determined by the best-fit procedure (something that is consistent with the fact that we have employed -gas data to build the likelihood function).

In order to quantify the correlation between the parameters and , the covariances , and correlation coefficients , were evaluated according to Ref. [57]. The mean values and their respective standard deviations used to calculate the covariances and correlation coefficients were those given in Table 2. The same table shows the values obtained for the covariances and correlation functions for each value of taken in this paper.

The values of the correlation coefficient presented in Table 2 confirm what was stated before: The parameters and are strongly (negatively) correlated while and are weakly correlated for all values of B considered here.

Comparison of the results presented in Table 2 and Table 1 shows that there is no statistical difference between considering the mode and a confidence interval of 68% (Table 1) from taking the mean and standard deviation. This emphasizes the quasi-gaussian behaviour of the probability function.

### 4.2 Behaviour of the f(R,∇R)-model with time

We can investigate the evolution of our higher order model in time by plotting the function as a function of . This is done setting the parameters at their best-fit results for each value of .

There is agreement between our model – regardless the value of – and the standard CDM model in the region . In particular, the curves of for  and  almost coincide: the difference between them is only noticeable at . The transition from a decelerated phase to an accelerated expansion in the higher order model occurs at for all values of .

Towards the future () the redshift takes negative values, . In this region our model present different dynamics depending on the value of the coupling constant. For ,  and , there is an onset of an ever growing acceleration eventually leading to a premature Big Rip, respective to the one taking place in the CDM model. On the other hand, the velocity in the higher order model with  – black curve in Fig. 6 – reaches a maximum value at , the acceleration then changes the sign, leading to a positive deceleration parameter. This deceleration is sufficiently intense to produce a Rebouncing, where the universe reaches a maximum size and then begins to shrink. This happens at .

This difference between Big Rip and Rebouncing scenarios deserves a more careful statistical analysis. We shall consider this in a broader space of parameters, and not using only the best-fit values for the parameters of our model, as we have done so far. A careful analysis shows that the deceleration parameter  influences weakly on the possible future regimes of our model. The parameters that decisively decide in favor of a Big Rip scenario or a Rebouncing event are the jerk  and the snap . Fig. 7 displays the confidence regions plots for .

The straight red lines passing through the confidence regions of the plots Fig. 7 separate the sector of the space of parameters where the Rebouncing occurs (lower part) from that region of values corresponding to a Big Rip (upper part). A numerical analysis shows that 39% of the confidence region allows for the Rebouncing scenario in the higher order model with and 42% of the confidence region for the model with  leads to Rebouncing. A similar analysis shows that the model with  produces a premature Big Rip for all the values of the pair in the 95% confidence region; converselly, the model with  engenders Rebouncing within the 95% CR.

At last, it is interesting to explore the features of the model for . Extrapolating the results presented in Fig. 6 we verified that our model does not reproduce the decelerated phase typical of a matter dominated Friedmann universe. Besides, for the higher order model exhibits a behaviour compatible with555The result (51) is still valid even when radiation is added to the model.

 13<ω<1.253, (51)

where the equation of state is , i.e. the equation for a perfect fluid in the standard Friedmann cosmology. On the other hand, a linear perturbative analysis shows that the perturbations in gravitation and matter fields do not grow for the equation of state with as given in (51). Therefore, the model proposed here can not describe the structure formation in the early universe.

## 5 Final Remarks

In this paper we have studied the phenomenological features of a particular modified gravity model built with an action integral dependent on the Ricci scalar and also on the contraction  of its covariant derivatives. The Friedmann equations resulting from the field equation and the FLRW line element involve sixth-order time-derivatives of the scale factor . We solved those differential equations numerically – but without any approximation – assuming a null cuvature parameter () and pressureless matter. Our higher order model depends on five parameters, namely . The present day value of the Hubble constant is taken as given – and not constrained by data fit – the reason for that being consistency with the -gas data analysis. Four of the five parameters  are determined by the best-fit to the data for the SN Ia and fraction of gas in galaxy clusters  available from Refs.[8] and [11], respectively.

The shape of the likelihood function  in terms of shows that it is not possible to determine this parameter by maximization. We were forced to select some values of and proceed to the best-fit of the other parameters in the set . It was clear that the  is statistically independent of the initial conditions, regardless the values of and the other best-fit parameters. In addition the sign of both  and  depend on , so that the behaviour of the model towards the future is not unique.

The confidence regions for the pair of parameters , and were built. They showed that  and  are strongly correlated parameters, whereas  and  are weakly correlated quantities. This was confirmed by the results in Table 2. The plot of the confidence region for the pair  is relevant to understand the future evolution of our higher-order model, which can exhibit either a premature Big Rip or a Rebouncing depending on the value of and on the values of  and .

If additional data for high objects were used, as Gamma Ray Bursts, for instance, they might reduce the uncertainties in our estimations of the jerk and the snap, and help us to decide if the Rebouncing scenario is favored with respect to the Big Rip, or conversely.

A positive feature of the particular -model presented here is its consistency with the standard CDM dynamics at recent times. Maybe this is anticipated because the model has a great number of free parameters to be fit to the data. The high number of kinematic parameters (i.e., those involving acceleration and its first and second derivatives) is a consequence of the model being of the fifth-order in derivatives of , a fact that could always make room for agreement with the observations. However, the freedom for bias is not unlimited. Even with all these free parameters there is a constraint to the dynamics of the model towards the past.

Besides this consistency with the standard CDM dynamics (at recent times) the model does not reproduce a large scale structure formation phase. This shows that just fitting a model to actual data is not enough to falsify or prove the model. General essential characteristic must be present too, as the instabilities allowing the large scale structure formation. Therefore, the extra term presented in the action (1) should be ruled out as viable alternative for dark energy.

#### Acknowledgements

This paper is dedicated to Prof. Mario Novello on the ocasion of his 70th birthday. RRC thanks FAPEMIG-Brazil (grant CEX–APQ–04440-10) for financial support. CAMM is grateful to FAPEMIG-Brazil for partial support. LGM acknowledges FAPERN-Brazil for financial support.

## References

• [1] D. N. Spergel et al., First year Wilkinson Microwave Anisotropy Probe (WMAP) observations: Determination of cosmological parameters, Astrophys. J. Suppl. 148 (2003), 175; E. Komatsu et al., Seven-year Wilkinson Microwave Anisotropy Probe (WMAP) observations: Cosmological interpretation, Astrophys. J. Suppl. 192:18 (2011).
• [2] Planck Collaboration: P. A. R. Ade et al., Planck 2013 results. I. Overview of products and scientific results, arXiv: 1303.5062 (2013); Planck Collaboration: P. A. R. Ade et al., Planck 2013 results. XVI. Cosmological parameters, arXiv: 1303.5076 (2013).
• [3] D. J. Eisenstein et al., Detection of the Barion Acoustic Peak in the large-scale correlation function of SDSS luminous red galaxies, Astrophys. J. 633 (2005), 560; W. J. Percival et al., Baryon Acoustic Oscillations in the Sloan Digital Sky Survey Data Release 7 Galaxy Sample, Mon. Not. Roy. Astron. Soc. 401 (2010), 2148.
• [4] S. Cole et al., The 2dF Galaxy Redshift Survey: Power-spectrum analysis of the final dataset and cosmological implications, Mon. Not. Roy. Astron. Soc. 362 (2005), 505; F. Beutler, The 6dF Galaxy Survey: Baryon Acoustic Oscilations and the local Hubble constant, Mon. Not. Roy. Astron. Soc. 416 (2011), 3017.
• [5] P. Astier et al., The Supernova legacy Survey: Measurement of ,  and  from the first year data set, Astron. Astrophys. 447 (2006), 31.
• [6] A. G. Riess el al., Type Ia Supernova discoveries at from the Hubble Space Telescope: Evidence for past deceleration and constraints on dark energy evolution, Astrophys. J. 607 (2004), 665; A. G. Riess el al., New Hubble Space Telescope discoveries of Type Ia Supernovae at : Narrowing constraints on the early behavior of dark energy, Astrophys. J. 659 (2007), 98;
• [7] W. M. Wood-Vasey et al., Observational constraints on the nature of the dark energy: first cosmological results from the ESSENCE Supernova Survey, Astrophys. J. 666 (2007), 694.
• [8] R. Amanullah et al., Spectra and Hubble Space Telescope light curves of six type Ia supervonae at  and the Union2 Compilation, Astrophys. J. 716 (2010), 712.
• [9] N. Suzuki et al., The Hubble Space Telescope Cluster Supernova Survey: V. Improving the Dark Energy Constraints Above and Building an Early-Type-Hosted Supernova Sample, Astrophys. J. 746 (2012), 85.
• [10] S. Schindler, -different ways to determine the matter density of the universe, Space Science Reviews 100 (2002), 299.
• [11] S. W. Allen et al., Improved constraints on dark energy from Chandra X-ray observations of the largest relaxed galaxy clusters, Mon. Not. Roy. Astron. Soc. 383 (2008), 879.
• [12] S. Weinberg, The cosmological constant problem, Rev. Mod. Phys. 61 (1989) 1.
• [13] A. G. Riess et al., Observational evidence from supernovae for an accelerating universe and a cosmological constant, Astron. J. 116 (1998), 1009.
• [14] S. Perlmutter et al., Measurements of  and  from 42 high-redshift supernovae, Astrophys. J. 517 (1999), 565.
• [15] L. Amendola e S. Tsujikawa, Dark Energy: Theory and Observations, Cambridge University Press, New York (2010).
• [16] B. Ratra and P. J. E. Peebles, Cosmological consequences of a rolling homogeneous scalar field, Phys. Rev D 37 (1988), 3406.
• [17] R. R. Caldwell, R. Dave, and P. J. Steinhardt, Cosmological imprint of an energy component with general equation-of-state, Phys. Rev Lett. 80 (1998), 1582.
• [18] S. M. Carrol, Quintessence and the rest of the world, Phys. Rev Lett. 81 (1998), 3067.
• [19] A. Hebecker and C. Wetterich, Natural quintessence?, Phys. Lett. B 497 (2001), 281.
• [20] L. Amendola, G. C. Campos, and R. Rosenfeld, Consequences of dark matter-dark energy interaction on cosmological parameters derived from SN Ia data, Phys. Rev. D 75 (2007), 083506.
• [21] T. Chiba, T. Okabe, and M. Yamaguchi, Kinetically driven quintessence, Phys. Rev. D 62 (2000), 023511.
• [22] C. Armendariz-Picon, V. F. Mukhanov, and P. J. Steinhardt, Essentials of k-essence, Phys Rev. D 63 (2001), 103510.
• [23] A. Y. Kamenshchik, U. Moschella, and V. Pasquier, An alternative to quintessence, Phys. Lett. B 511 (2001), 265.
• [24] M. C. Bento, O. Bertolami, and A. A. Sen, Generalized Chaplygin gas, accelerated expansion and dark energy-matter unification, Phys. Rev. D 66 (2002), 043507.
• [25] S. Capozziello and M. De Laurentis, Extended theories of gravity, Phys. Rep. 509 (2011), 167.
• [26] S. Nojiri and S. D. Odintsov, Unified cosmic history in modified gravity: from theory to Lorentz non-invariant models, Phys. Rep. 505 (2011), 59.
• [27] G. R. Dvali, G. Gabadadze, and M. Porrati, 4D gravity on a brane in 5D minkowski space, Phys. Lett. B 485 (2000), 208.
• [28] V. Sahni and Shtanov, Braneworld models of dark energy, JCAP 0311 (2003), 014.
• [29] N. Bartolo and M. Pietroni, Scalar tensor gravity and quintessence, Phys. Rev. D 61 (2000), 023518.
• [30] F. Perrota, C. Baccigalupi, and S. Matarrese, Extended quintessence, Phys. Rev. D 61 (2000), 023507.
• [31] Z. Chang, X. Li, Modified Friedmann model in Randers-Finsler space of approximate Berwald type as a possible alternative to dark energy hypothesis, Phys. Lett. B 676 (2009) 173.
• [32] K. S. Adhav, LRS Bianchi Type-I Universe with Anisotropic Dark Energy in Lyra Geometry, Int. J. Astr. Astrophys., 1 No. 4, (2011) 204. doi: 10.4236/ijaa.2011.14026.
• [33] R. Casana, C. A. M. de Melo, B. M. Pimentel, Massless DKP field in a Lyra manifold, Class. Quantum Gravity 24 (2007) 723.
• [34] S. Capozziello, Curvature quintessence, Int. J. Mod. Phys. D 11 (2002), 483.
• [35] A. De Felice and S. Tsujikawa,  theories, Living Rev. Rel. 13 (2010), 3.
• [36] T. P. Sotiriou and V. Faraoni,  theories of gravity, Rev. Mod. Phys. 82 (2010), 451.
• [37] J. Santos et al., Latest supernovae constraints on cosmologies, Phys. Lett. B 669 (2008), 14.
• [38] N. Pires, J. Santos and J. S. Alcaniz, Cosmographic constraints on a class of Palatini gravity, Phys. Rev. D 82 (2010), 067302.
• [39] S. Nojiri and S. D. Odintsov, Introduction to modified gravity and gravitational alternative for dark energy, Int. J. Geom. Meth. Mod. Phys. 4 (2007), 115.
• [40] G. Bengochea and R. Ferraro, Dark torsion as the cosmic speed-up, Phys. Rev. D 79 (2009), 124019.
• [41] E. V. Linder, Einstein’s other gravity and the acceleration of the universe, Phys. Rev. D 81 (2010), 127301.
• [42] K. Bamba et al., Equation of state for dark energy in gravity, JCAP 1101 (2011), 021.
• [43] S. Gottlöber, H-J Schmidt, A. A. Starobinsky, Sixth-order gravity and conformal transformations, Class. Quantum Gravity 7 (1990), 893.
• [44] T. Biswas, A. Mazumdar and W. Siegel, Bouncing universes in string-inspired gravity, JCAP 0603 (2006), 009.
• [45] T. Biswas, T. Koivisto and A. Mazumdar, Towards a resolution of the cosmological singularity in non-local higher derivative theories of gravity, JCAP 1011 (2010), 008.
• [46] T. Biswas et al., Stable bounce and inflation in non-local higher derivative cosmology, JCAP 1208 (2012), 024.
• [47] N. Arkani-Hamed et al., Non-local modification of gravity and the cosmological constant problem, hep-th/0209227.
• [48] A. O. Barvinsky, Nonlocal action for long-distance modifications of gravity theory, Phys. Lett. B 572 (2003), 109.
• [49] R. R. Cuzinatto, C. A. M. de Melo, and P. J. Pompeia, Second order gauge theory, Ann. Phys. 322 (2007), 1211.
• [50] R. R. Cuzinatto, C. A. M. de Melo, and P. J. Pompeia, Gauge formulation for higher order gravity, Eur. Phys. J. C 53 (2008), 99.
• [51] R. R. Cuzinatto, C. A. M. de Melo, and P. J. Pompeia, Cosmic acceleration from second order gauge gravity, Astrophys. Space Sci. 332 (2011), 201.
• [52] M. Visser, Jerk, snap, and the cosmological equation of state, Class.Quant.Grav. 21 (2004) 2603.
• [53] L. G. Medeiros, Realistic Cyclic Magnetic Universe, Int. J. Mod. Phys. D 21 (2012), 1250073.
• [54] D. E. Holz and E. V. Linder, Safety in numbers: Gravitational Lensing Degradation of the Luminosity Distance-Redshift Relation, Astrophys. J. 631 (2005) 678.
• [55] D. Kirkman et al., The cosmological baryon density from the deuterium to hydrogen ratio towards QSO absorption systems: D/H towards Q1243+3047, Astrophys. J. Suppl. 149 (2003) 1.
• [56] A. G. Riess el al., A 3% solution: determination of the Hubble constant with the Hubble Space Telescope and Wide Field Camera 3, Astrophys. J. 730 (2011), 119.
• [57] J. Beringer et al., Review of Particle Physics, Phys. Rev. D 86 (2012), 010001.

You are adding the first comment!
How to quickly get a good reply:
• Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
• Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
• Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters