Is exponential gravity a viable description for the whole cosmological history?

Here we analyse a particular type of F(R) gravity, the so-called exponential gravity which includes an exponential function of the Ricci scalar in the action. Such a term represents a correction to the usual Hilbert-Einstein action. By using Supernovae Ia, Barionic Acoustic Oscillations, Cosmic Microwave Background and H(z) data, the free parameters of the model are well constrained. The results show that such corrections to General Relativity become important at cosmological scales and at late times, providing an alternative to the dark energy problem. In addition, the fits do not determine any significant difference statistically with respect to the [Formula: see text]CDM model. Finally, such model is extended to include the inflationary epoch in the same gravitational Lagrangian. As shown in the paper, the additional terms can reproduce the inflationary epoch and satisfy the constraints from Planck data.

Introduction

Over the last decade, the study of some modifications of General Relativity have drawn a lot of attention, particularly in the framework of cosmology as an attempt to provide a more natural explanation to the accelerating expansion at early times (inflation) and at late times (dark energy epoch). In this sense, the most simple and natural extension of GR arises as the generalisation of the Hilbert–Einstein action by assuming a non-linear function of the Ricci scalar, which is commonly called f(R) gravity (for a review see [1-8]). Other extensions include curvature invariants such as the Gauss–Bonnet gravity [9-13] or generalisations of the so-called Teleparallel gravity, an equivalent theory to GR constructed as a gauge theory of the translation group leading to a null-curvature theory with non-null torsion (see Ref. [14-17]). Nevertheless, f(R) gravities have been by far the most analysed extension of GR over the last years, also due to the motivation as regards more fundamental theories as string theory [18]. This extensive study has provided a very deep knowledge and comprehension of this type of theories, whose field equations turn out fourth order differential equations instead of second order as in GR. Nevertheless, f(R) gravities can easily be reduced to a type of scalar–tensor theory, i.e. f(R) gravity basically implies the appearance of an extra scalar mode [19-21]. As every theory with extra propagating modes, this may imply the existence of ghosts. Fortunately, this is not the case in f(R) gravities. However, the extra scalar mode may imply violations and deformations of well-known and tested predictions of GR. In order to avoid large corrections at scales where GR is very well tested, f(R) gravities can hide such extra mode through a mechanism known as chameleon mechanism, proposed initially in the framework of scalar–tensor theories [22, 23], but rapidly extended to f(R) gravities [24, 25].

In addition, the versatility of f(R) gravities allows one to reconstruct any cosmological solution with the suitable evolution [26-35]. Then late-time acceleration may arise in a natural way as a consequence of the gravitational theory instead of being the aftermath of any extra unknown field. Moreover, simultaneously f(R) gravities may contribute to the compensation of the large value predicted by quantum field theories for the vacuum energy density, and particularly may play an essential role in the framework of the so-called unimodular gravity theories. In this regard, f(R) gravity scenarios as an alternative to the CDM-cosmology are interesting and attractive, since they are able to describe simultaneously the early-time inflation and the late-time acceleration in the expansion of our universe [19–21, 25, 36–39]. Particularly, some of the most promising inflationary models are constructed within the f(R) gravity scenario, since some of these models can easily reproduce slow-roll inflation by mimicking a cosmological constant at early times and then decaying, leading to a power spectrum for scalar perturbations nearly invariant and a negligible scalar-to-tensor ratio, coinciding with the last data released by the Planck collaboration [40, 41]. This is the case for instance of Starobinsky inflation [42], a quadratic Lagrangian on the Ricci scalar that predicts the correct values for the spectral index and the scalar-to-tensor ratio. Actually, some analysis suggest that any deviation from Starobinsky inflation should be small enough to avoid deviations from its well-established predictions [43]. Keeping this in mind, over the last years some efforts have been focussed on the attempt to unify inflation and dark energy epoch in the framework of f(R) gravities, and particularly within the so-called viable f(R) gravity models [25]. As mentioned above, these viable models pass the well-known local tests, where the scalar mode acquires a large mass through the chameleon mechanism avoiding large corrections with respect to GR. Hence, the local tests or the Solar System tests for viable f(R) theories include correct Newtonian and post-Newtonian limit [24, 25]. In addition, this type of models are capable of reproducing the correct late-time acceleration, in general by simulating an effective cosmological constant that becomes important at late times, while CDM behaviour is recovered at high redshift. Moreover, these models provide good fits when compared with observational data, being almost indistinguishable from CDM [44]. However, viable f(R) gravities contain a type of future cosmological singularity, the so-called sudden singularity, a consequence directly related to the mass of the scalar field that avoids corrections at local scales [45], although such singularity occurs in the future when the right parameters are set and can be avoided by adding some extra terms. Moreover, some extensions of such models are also capable of reproducing inflation at early times, when tends asymptotically to a power Lagrangian, leading to a Starobinsky-like inflation keeping the right predictions [25].

In this paper, we focus on the analysis of a type of viable f(R) models that reproduces late-time acceleration by mimicking a cosmological constant but where corrections may have some distinguishable effects. This class of f(R) models are given by a negative exponential of the Ricci scalar in the action, which turns out negligible at large redshifts but becomes important at late times, an effect easily controlled with a free parameter related to current Hubble parameter. Exponential gravity has been previously analysed in Refs. [39, 46–48] as a reliable alternative to other viable f(R) gravities, since GR results are recovered at local scales but reproduce dark energy behaviour at cosmological ones. In addition, the previous analysis has shown the existence of an asymptotically stable de Sitter solution in such exponential Lagrangians, leading to an approximated CDM behaviour at the present time [47, 48]. Moreover, some recent analysis of such type of exponential gravities suggest that observational constraints can be well satisfied from the cosmological point of view, in such a way that f(R) gravity and CDM model turn out nearly indistinguishable, as suggested by previous analysis [49-51]. In addition, exponential gravity can be extended to cover the inflationary stage as well. To do so, an additional exponential is considered in the gravitational action becoming important at large curvature when the inflationary period occurs, and turning out to be negligible as curvature decreases [39, 47]. Hence, in this paper we analyse such a type of f(R) gravities, firstly by fitting the free parameters of the model by using data from Type Ia supernovae, baryon acoustic oscillations (BAO), estimations of the Hubble parameter H(z) and parameters of the cosmic microwave background radiation (CMB) [52-55], and also considering different approaches. Then we analyse how the full gravitational Lagrangian can cover also the inflationary epoch, obtaining the spectral index for scalar perturbations and the tensor-to-scalar ratio.

The paper is organised as follows: Sect. 2 reviews the basics of f(R) gravities, while Sect. 3 is devoted to the introduction of the exponential f(R) gravity model and its dynamical equations. In Sect. 4 the observational data considered in the paper is shown, this includes Union 2.1 observations of Type Ia supernovae, BAO effects, the latest measurements of the Hubble parameter H(z) and CMB parameters. In Sect. 5 we estimate the constraints on the exponential F(R) model from the aforementioned data. In Sect. 6 we investigate the variant of the exponential model with inflation terms in the Lagrangian. Finally Sect. 7 is devoted to the conclusions of the paper.

F(R) gravity

Modified F(R) gravities are described by the following generalisation of the Einstein–Hilbert action [26-35]:where and is the matter action. Einstein General Relativity is very well understood and tested at many scales, so that one should assume the action (2.1) to contain slightly deviations from GR, such that we can rewrite the action in the following way:Here, the function f(R) accounts for the gravitational modifications and should become negligible at scales where GR is well tested. By varying the action (2.1) with respect to the metric tensor , the field equations are obtained:where R and are the Ricci scalar and Ricci tensor, respectively, whereas and is the energy-momentum tensor of matter. Note that the F(R) field equations are fourth order to be compared to the second order of General Relativity. However, the action (2.1) hides an additional scalar mode, such that it can be expressed by the Lagrangian of a type of scalar–tensor theory as follows:where the following relations are found:Hence, in order to avoid large deviations from GR, this additional degree of freedom should be hidden at the appropriate scale, a mechanism commonly known as the chameleon mechanism [22, 23]. In this sense, some F(R) actions which accomplish this requirement have been proposed in the literature [24, 25], particularly some of them with the form of a negative exponential, the type of Lagrangians we are exploring in this manuscript. Nevertheless, let us first analyse the general properties of F(R) gravities, and in particular in the cosmology framework. By assuming a spatially flat Friedman–Lemaître–Robertson–Walker (FLRW) space-time with the metricwhere a(t) is the scale factor of the universe, , and the Ricci scalar is expressed asHere, the Hubble parameter is defined as usual by , where the dot denotes derivatives with respect to the cosmic time. By assuming an energy-momentum tensor as a perfect fluid, where and p are the matter energy density and pressure, the field equations (2.3) turn out to be [26-35]While the divergence of the field equations lead to the energy conservation equation , which in a FLRW metric becomesThe FLRW equations (2.7) can be expressed in terms of other independent variables instead of the cosmic time for convenience. Here we use the number of e-folds, given by with at the present time , to express the above equations in the form of a dynamical system as follows:Here are derivatives with respect to x and we have used the Ricci scalar definition Eq. (2.6) and the continuity equation (2.8). Hence, the analysis of the above system can provide all the information as regards the dynamics produced by a particular action F(R).

Exponential gravity

Let us now introduce the type of exponential F(R) gravity, we are considering in this manuscript [39, 46–51]The model contains just two free parameters and , which may be expressed in a more convenient way as , where is dimensionless [49-51],Note that in principle the model (3.1) can well describe the universe evolution for , including the recombination epoch, the matter-dominated era and the late-time acceleration. This is true as far as , since the exponential becomes negligible and the action (3.1) recovers the usual CDM model at large redshifts, where the curvature becomes much larger than . In Sect. 6 we will also consider corrections to this model such that early-time inflation is also described.

Here we focus on the epoch for , when the content of the universe includes pressureless (non-relativistic) matter and radiation (relativistic particles): , such that the continuity equation (2.8) can be solved and yieldswhere and are the present time values of these components, which can be normalised over the critical density as follows:Here is the Hubble parameter today. Let us first explore the behaviour of model (3.1) during the early universe, when the curvature becomes large as (or for in practical applications). Then the model (3.1) transforms into the CDM model with , so the solutions of the system (2.9) tend asymptotically to CDM at large redshifts, leading toHere the index “” refers to quantities calculated within the CDM model, where and the Hubble parameter today as predicted by the CDM model, while . However, despite the model (3.1) recovering CDM at large redshifts, late-time evolution deviates from CDM, such that the above quantities as measured today, , would differ from CDM unless initial conditions are fixed at , which is not the case of our paper. Note that other viable f(R) models shows a similar behaviour when they are analysed asymptotically [24]. Hence, we have in generalwhere we have denoted by 0 those magnitudes measured today as predicted by our model (3.1). Nevertheless, we can connect the two models through the relation of the physical matter density [24]As will be shown below, this remark is important when performing the fitting analysis for the observable parameters in Sect. 4. Moreover, note that the first FLRW Eq.(2.7) for the CDM model is a constraint equation which evaluated at can be expressed as follows:This expression is very well known in standard cosmology when GR is assumed but breaks down when other gravitational actions beyond GR are considered, like F(R) gravity. In such a case, the first FLRW equation becomes a dynamical equation, since it contains second derivatives of the Hubble parameter. By evaluating the FLRW equation in (2.7) at , the above equation can be expressed asNote that here we have defined , which refers to the cosmological constant term in the action (3.1), while includes the exponential function in (3.1). The smaller is, the closer our model is to CDM at the present time, where Eq. (3.7) is evaluated. Nevertheless, note that our model recovers CDM asymptotically at high redshifts () such that the differences among the relative densities become negligible in both models at high redshifts.

Let us now, for convenience in the calculations, introduce the following dimensionless variables:Hence, the gravitational action (3.1) becomeswhile the system of equations (2.9) takes the form Recall that the variable refers to the number of e-folds. This system can be solved numerically by setting the appropriate initial conditions. As is natural for the model (3.1), we assume initial conditions that match the CDM model at a particular high redshift:which corresponds to the CDM asymptotic solution (3.4) at an initial redshift , or alternatively at . The value of is determined by assuming the following condition:where is a small number in the range , such that our model mimics the CDM solution (3.4) at , and the corresponding solutions practically do not depend on or for all x.

Alternatively, we can also use the following variable [24, 47, 49–51]:Then the equation for H in (2.9) can be rewritten as follows:where the second expression corresponds to the initial condition (3.11). The advantage of the function and its derivative (3.14) lies on their finiteness at (). However, the same does not apply for Eq. (3.10), since is not finite for every redshift. Hence, numerical integration of the system Eq. (3.14) or the corresponding second order differential equation for shows similar difficulties to the system (3.9) and (3.10).

Observational data

Let us now present the data we are using here to fit the free parameters of our model. Besides the late-time evolution data from SNe Ia, BAO and H(z), we are also considering the CMB parameters. Then, to do so, we have to include radiation in our equations, or in other words, assuming Eq. (3.11), or alternatively Eq. (3.14). As our model mimics CDM at high redshifts, we can reduce the number of free parameters by fixing the radiation–matter ratio as provided by Planck [55]:Hence, our model contains four free parameters (3.1):Recall that the Hubble parameter differs from the true Hubble constant , as well as the density parameters . Nevertheless, the Hubble constant can be considered as a nuisance parameter, so that can be marginalized for all fits.

Here, we use the catalogue provided by the Union 2.1 data, which contains 580 points from Type Ia Supernovae (SNe Ia) [52]. BAO data described in Table 1, Refs. [56-74]. We also use 30 estimations of the Hubble parameter H(z) measured from differential ages of galaxies and summarised in Table 2, [75-80]. Finally, the CMB parameters are considered from the Planck mission [55]. In order to proceed with the analysis we use the technique of the minimum , which establishes the best set of the parameters. To do so, we use a two-dimensional grid, such that the number of free parameters (4.2) is reduced to two, either by theoretical considerations or through marginalisation.

Table 1

Values of and A(z) (4.5) with errors and references

z	d_z(z)	σ_d	A(z)	σ_A	Refs.	Survey
0.106	0.336	0.015	0.526	0.028	[68]	6dFGS
0.15	0.2232	0.0084	–	–	[73]	SDSS DR7
0.20	0.1905	0.0061	0.488	0.016	[66, 69]	SDSS DR7
0.275	0.1390	0.0037	–	–	[66]	SDSS DR7
0.278	0.1394	0.0049	–	–	[67]	SDSS DR7
0.314	0.1239	0.0033	–	–	[69]	SDSS LRG
0.32	0.1181	0.0026	–	–	[62]	BOSS DR11
0.35	0.1097	0.0036	0.484	0.016	[66, 69]	SDSS DR7
0.35	0.1126	0.0022	–	–	[70]	SDSS DR7
0.35	0.1161	0.0146	–	–	[59]	SDSS DR7
0.44	0.0916	0.0071	0.474	0.034	[69]	WiggleZ
0.57	0.0739	0.0043	0.436	0.017	[60]	SDSS DR9
0.57	0.0726	0.0014	–	–	[62]	SDSS DR11
0.60	0.0726	0.0034	0.442	0.020	[69]	WiggleZ
0.73	0.0592	0.0032	0.424	0.021	[69]	WiggleZ
2.34	0.0320	0.0021	–	–	[65]	BOSS DR11
2.36	0.0329	0.0017	–	–	[64]	BOSS DR11

Table 2

Hubble parameter values H(z) with errors from Refs. [75–80]

z	H(z)	σ_H	Refs.	z	H(z)	σ_H	Refs.
0.070	69	19.6	[78]	0.4783	80.9	9	[80]
0.090	69	12	[75]	0.480	97	62	[76]
0.120	68.6	26.2	[78]	0.593	104	13	[77]
0.170	83	8	[75]	0.6797	92	8	[77]
0.1791	75	4	[77]	0.7812	105	12	[77]
0.1993	75	5	[77]	0.8754	125	17	[77]
0.200	72.9	29.6	[78]	0.880	90	40	[76]
0.270	77	14	[75]	0.900	117	23	[75]
0.280	88.8	36.6	[78]	1.037	154	20	[77]
0.3519	83	14	[77]	1.300	168	17	[75]
0.3802	83	13.5	[80]	1.363	160	33.6	[79]
0.400	95	17	[75]	1.430	177	18	[75]
0.4004	77	10.2	[80]	1.530	140	14	[75]
0.4247	87.1	11.2	[80]	1.750	202	40	[75]
0.445	92.8	12.9	[80]	1.965	186.5	50.4	[79]

Supernovae Ia data

The Union 2.1 compilation provides [52] Sne Ia with their observed (estimated) distance moduli for redshifts in the interval . In order to fit the free parameters of our model, we compare with the theoretical value , where the distance moduli is given byHere is the luminosity distance. The corresponding function is calculated by computing the differences between the SNe Ia observational data and the predictions of a particular model with parameters ,where , is the covariance matrix [52]. The marginalisation over the nuisance parameter is widely described in the literature (see Refs. [81-83]).

BAO data

Baryon acoustic oscillations (BAO) are obtained from galaxy clustering analysis and include measurements of two cosmological parameters [53]where is the sound horizon at the decoupling epoch and is given byThe values (4.5) were estimated for redshifts (and redshift ranges) of galaxies from a peak in the correlation function of the galaxy distribution at the comoving sound horizon scale , which corresponds to the decoupling of the photons . In this paper we use the BAO data from Refs. [56-74] for the parameters (4.5), which provides data points for and 7 data points for A(z), both shown in Table 1. We use the covariance matrices and for correlated data from Refs. [66, 69] described in detail in Ref. [81]. So the function for the values (4.5) yieldswhere and are vector columns with and .

As pointed out above, the Hubble parameter today as predicted by our model (3.1) differs from the one predicted by the CDM model , which is considered here as a free parameter. Both are related by the expression . However, the theoretical values of and A (4.5) do not contain , since the distances (4.3), and are proportional to . In the expression for A(z) we can use the equivalence (3.5) .

All these considerations have to be carefully studied in order to choose the appropriate approach to calculating the sound horizon from different fitting formulae [50, 74, 84, 85]. Here we are considering the following simple fitting formula [81]:with explicit h dependence. For the CDM model, one obtains Mpc as the best fit (see Ref. [81] ). Other approaches give the same predictions as using Eq. (4.7); see Ref. [74].

H(z) data

The Hubble parameter values H at certain redshifts z can be measured with two methods: (1) extraction of H(z) from line-of-sight BAO data [56-65] and (2) making H(z) estimations from differential ages of galaxies [75-80] via the following relation:To avoid additional correlation with the BAO data from Table 1, we use in this paper only values H(z) estimated from differential ages of galaxies, shown in Table 2. The theoretical values naturally depend on . so the function is marginalized over [86]:

CMB data

Unlike the SNe Ia, BAO and H(z) data described above, corresponding to the late-time era , cosmological observations associated with CMB radiation [74, 84, 85] include parameters at the photon-decoupling epoch ( [55]), particularly the comoving sound horizon and the transverse comoving distance,In the present manuscript, we use the CMB parameters in the following form [84, 85]:with the estimations (distance priors) from Ref. [85]Here is the present time baryon fraction. The distance priors (4.11) with their errors and the covariance matrixwere derived in Ref. [85] from the Planck collaboration data [55] with free amplitude of the lensing power spectrum. For the value we use the fitting formula from Refs. [84, 85, 87]; the sound horizon is estimated from Eq. (4.9) as the correction .

Hence, the function corresponding to the data (4.10)–(4.11) is obtained as follows:which is minimised by marginalizing over the additional parameter , which should be considered as a nuisance parameter, as well as over or . However, for the joint analysis of H(z) and CMB data, the marginalisation over is calculated simultaneously:Let us now present the results for the f(R) model considered here.

Testing exponential F(R) gravity

By considering the SNe Ia, H(z), BAO and CMB data illustrated in the previous section, the above exponential model is well constrained. Here, we calculate these limitations and the best-fitted values of the parameters for the exponential F(R) model (3.1). After marginalizing over (and over for the CMB data in ), the F(R) model (3.1) has three free parameters: , and . Recall that they differ from and , these values are connected and , as shown above in Eq. (3.5). Consequently, the sum as given in (3.7). The sum (3.6) is also not equal 1 in general for the considered F(R) model.

However, firstly let us assume the following condition:This means that the model is assumed to be close to CDM. This assumption relaxes the difficulties to fit the free parameters, since the free parameters of the model can be automatically reduced, leading to two free parameters: and .

The results are depicted in Fig. 1 where the , and regions are shown in the contour plots for the plane (the top-left panel) and for the plane (the top-right panel). The magenta contours correspond to whereas the black lines describe .

Fig. 1

Upper panels: contour plots for the free parameters of the exponential F(R) model (3.1) when assuming , the left panel shows the plane while the plane is depicted in the right panel. Bottom panels: the one-dimensional dependencies of with respect to (left panel) and to (right panel)

At each point in the plane, or in other words, for given values of , (and ), the differential equations (3.9), (3.10) are solved by assuming the CDM model as the initial conditions at high redshift (3.11). Then, once the solution E(x) is determined for each set of the free parameters, the functions: (4.4), (4.6) are obtained. Furthermore, by considering the function , we calculate then the and the by marginalizing over (and over for ), or, in other words, keeping as a nuisance parameter, which avoids further bias on the results, so that we obtain the optimal value (for these fixed and ) and calculate the corresponding values and from Eqs. (3.8) and (3.5); see Ref. [88].

By following this procedure, we can calculate at each point and draw the contour plots in the plane, as shown in the top-right panel of Fig. 1. These calculations were made separately for SNeBAO data (the filled contours) and for SNeCMB data (the black contours). The difference between the panel and the planel is clearly shown for small , while in the limit the model (3.1) tends to the CDM, where and coincide.

In the bottom-left panel of Fig. 1 the one-dimensional dependencies of on are shown for (the black dash-dotted line), for (the solid magenta line) and for (the green dashed line). The latter is not depicted in the upper panels, however, for all three cases, the one-dimensional distributions are calculated from the corresponding two-dimensional matrices describing the contours in the top-left panel. Hence, under the restriction (5.1), the dependence of with respect to the CMB data (4.10)–(4.11) on is very strong (unlike for the late-time SNeBAO data ).

Similarly, the one-dimensional dependency of on (calculated from the two-dimensional matrices ) are depicted at the bottom-right panel of Fig. 1. Under the restriction (5.1) for all functions, the absolute minimum is achieved at the limit , in other words, at the “CDM” limit of the model (3.1).

This conclusion may be illustrated in another way: the curves of the bottom-left panel of Fig. 1 will coincide with the ones for the flat CDM model, since these minima are achieved at large , where coincide with .

While in the bottom-right panel of Fig. 1, one can see the unusual behaviour of the one-dimensional distributions and : these functions have the local minima at . These minima are not shown in the top panels of Fig. 1, because they lie beyond the confidence levels. However, suitable level lines of const show these local minima as “islands” in the or planes. For instance, the corresponding coordinates or optimal values for are , .

Furthermore, let us now consider the model (3.1) for its general case (3.7) beyond the restriction (5.1) . In this case the model has three free parameters: , and (or alternatively and ). So in order to calculate the corresponding , the value of the three parameters has to be given before solving numerically the system (3.9), (3.10) with the initial conditions (3.11), as described above. In this case every function (after marginalisation over and for and ) will depend on , and .

Hence, when we draw the contour plots for and in the plane in the top-left panel of Fig. 2, we minimize these functions over at each point of the plane. In other words, we calculate for and .

Fig. 2

General case where . Upper panels: contour plots for (left panel) and (right panel). Bottom plots: one-dimensional distributions for (solid black lines) and for (solid blue lines and filled contours) are compared with the ones from Fig. 1 for (magenta lines) and for (thin black dash-dotted lines)

At each point of the plane, the minima of the and functions are calculated, also the optimal values of the free parameters , , and are obtained. These values help us to draw the contour plots in the plane in the top-right panel of Fig. 2.

Table 3

Predictions of the exponential F(R) model (3.1) and the CDM for different data sets (, ): and estimates of model parameters

Model	Data	ΩmΛCDM	Ωm0	ΩΛΛCDM	β	Minχ^2/d.o.f
F(R) (3.1)	χΣ32	0.282-0.009+0.010	0.285-0.010+0.012	0.696-0.037+0.025	3.36-2.16+∞	572.07 / 631
F(R) (3.1)	χtot2	0.280-0.001+0.001	0.294-0.007+0.009	0.637-0.062+0.047	2.38-0.80+∞	575.51 / 634
ΛCDM	χΣ32	=Ωm0	0.282-0.009+0.010	0.718-0.010+0.009	∞	572.93 / 633
ΛCDM	χtot2	=Ωm0	0.2772-0.0004+0.0003	0.7228-0.0003+0.0004	∞	583.24 / 636

The same panels and notations of Fig. 1 are used in Fig. 2, but the blue contours corresponds to in Fig. 2 while the blue lines refer to the one-dimensional distributions and in the bottom panels of Fig. 2. In order to compare these results with those obtained under the approximation (5.1), the curves of Fig. 1 are depicted as well, denoted by magenta lines for and by thin black dash-dotted lines for (in the bottom panels).

The black stars in Fig. 2 denote the minimum points of the two-dimensional distributions and . Their coordinates (the optimal values of parameters) are tabulated in Table 3. In the same way, the minimum points for are shown as the blue circles.

The contour plots in Fig. 2 demonstrate that the model (3.1) in absence of the approximation (5.1), and with the free parameters , and the regions of , or confidence level in the plane are essentially enlarged in comparison with Fig. 1. Particularly, for the domains (the blue filled contours) occupies the range , whereas under the restriction (5.1) (the magenta contours in Figs. 1 and 2) the range is . For the joint function (the black contours) these areas are larger in the direction and wider in the direction, especially for the parameter , as shown in the top-right panel of Fig. 2.

These enlarged domains of suitable model parameters include the above-mentioned “islands” or local minima of and functions existed under the restriction (5.1). This effect is hidden in the top panels of Fig. 2 (it is beyond the confidence level), but it is shown in the bottom-right panel, where the local minima at from Fig. 1 (for the magenta and black dash-dotted lines) are naturally included in the general behaviour of (the blue line) and (the solid black line). These one-dimensional distributions determine the optimal values and errors of the parameter ; this information is included in Table 3, where the absolute minima of and the mean of the model parameters are provided.

The optimal values and errors for in Table 3 are deduced from the one-dimensional distributions . They are shown in the bottom-left panel of Fig. 2 as the solid black and blue curves in comparison with the corresponding plots for the case (the dash-dotted and magenta lines; they are taken from Fig. 1). The latter distributions coincide with the predictions of the CDM model.

One can conclude from Table 3 (and Fig. 2) that the absolute minima for the F(R) model (3.1) and are smaller than the ones for CDM model (572.93 and 583.24 respectively). Such a difference lies in the existence of degrees of freedom for the model (3.1).

In Fig. 3 we demonstrate how the numerical solution of the system (3.9), (3.10) behaves (the black solid lines in the phase space diagram) in comparison with the CDM model (the red dashed lines). We show the plots for the Hubble parameter and the Ricci scalar (in the right panels) depending on a in the top panels with the logarithmic scale and the same plots E(z), (with z instead of a) in the usual scale in the bottom panels. The model parameters for both models are taken from Table 3; they are optimal for .

Fig. 3

The top panels in the logarithmic scale illustrate the dependence of , (the blue dots) and on a for the F(R) model (3.1) (the black solid lines) and for the CDM (the red dashed lines). The corresponding plots E(z), and in the usual scale are shown at the bottom panels. For both models the parameters are from Table 3 for

One can see that for the optimal parameters the F(R) and CDM models demonstrate rather close dynamics of the Hubble parameter E with the small future divergence for . The plots for the density (the blue dots) are the same for both models.

For the parameter , the optimal values and errors in Table 3 were found after preliminary calculation of two-dimensional distributions for and . The contour plots for these two-dimensional distributions and the corresponding one-dimensional plot for are depicted in Fig. 4 with the same notations. We see that the CMB data (4.10–4.11) for essentially restrict the bounds of in comparison with the function for SNeBAO data. This is also connected with the factor in the parameter R (4.10) and the narrow restrictions (4.11). The minimum value for (the black star in Fig. 4) is shifted lower from the “CDM line” .

Fig. 4

Contour plot for the plane and the one-dimensional distribution for (solid black lines) and for (solid blue lines and filled contours). The red lines correspond to the CDM model

In the right panel of Fig. 4 we compare the predictions of the F(R) model (3.1) and the CDM model, where . The CDM dependencies coincide with their analogs for the model (3.1) under the restriction (5.1) (after minimising over ). In the CDM case, the errors are essentially smaller. Finally, Fig. 5 shows the evolution of the luminosity distance and the Hubble rate for the exponential gravity model and for CDM for their best fits. As shown, both curves fit the data similarly, such that the two models become indistinguishable.

Fig. 5

Luminosity distance (left panel) and Hubble parameter (right panel) for the best fit when considering the general case described by (3.7) for the model (3.1). The best fit for CDM is also depicted

Exponential model and inflation

The exponential F(R) model (3.1) considered in the previous sections, describes all observational manifestations of the late-time acceleration. However, such a model can also explain the early-time inflation when introducing some suitable modifications in the form of F(R) as follows [47]:These additional terms can generate the expected cosmological constant during the inflationary era, when R is near or larger than . The natural number helps to avoid the effects of inflation during the matter era when and the last term in Eq. (6.1) is necessary for a successful exit from inflation. As pointed out in Ref. [47], the model (6.1) has the following properties: a de Sitter phase naturally arises during inflation in the high-curvature regime, the inflationary terms become negligible in the small curvature era , inflation ends successfully and avoids anti-gravity effects and instabilities during the matter era. To satisfy all these properties, the following requirements are obtained over the free parameters:where can vary in the rangeIn addition, inflation occurs when , where an unstable de Sitter point arises in the phase space, if the condition [47]for the functionis satisfied. The condition (6.4) may be deduced, if we search a de Sitter solution const, const of the system (2.9) before the matter era.

If , we can neglect the terms in G(R) and the condition (6.4) for de Sitter solutions takes the formHere the constant in Eq. (6.1) is expressed as . In Ref. [47] the de Sitter solution with parameters , , , was considered. The corresponding values here are , with satisfy the condition (6.5). However, below we analyse a wider set of inflationary solutions which obey the observational limitations for the slow-roll parameters.

Let us now focus on the realisation of slow-roll inflation within the model (6.1) and its predictions. In order to do so, the scalar–tensor counterpart of F(R) gravities is more convenient than its original action, such that the F(R) gravities can be expressed in terms of a scalar field as shown in (2.4) and (2.5). By varying the action (2.4) with respect to the scalar field , it yieldsHere recall that the relations of Eq. (2.5) for the scalar field and its potential hold:In order to analyse slow-roll inflation for the model (6.1), the action (2.4) can be transformed into the Einstein frame by the following conformal transformation:The action (2.4) turns out to beHere, we have redefined the scalar field to keep the kinetic term in a canonical form:The FLRW equations for the action (6.8) become simpler than working directly on the F(R) action:whereas the scalar field equation is given byDuring slow-roll inflation the scalar field mimics an effective cosmological constant, what basically means that and . Both conditions can also be expressed through the so-called slow-roll parametersThese quantities remain very small when inflation occurs such that and , while at the end of inflation . In addition, the slow-roll parameters (6.12) are related to the amplitude and scale dependence of the perturbations originated during inflation, such that the spectral index of the perturbations and the tensor-to-scalar ratio, are given bySince both values are very well constrained by the last data from Planck and Bicep2 collaborations [40, 41], which giveThen we can test whether the model (6.1) is capable to satisfy such constraints. Firstly, let us analyse the action (6.1), as during inflation, , the action (6.1) can be approximated as follows:Then, by Eq. (6.9), the relation among the scalar field and the curvature is obtained:Here recall that . During inflation , and consequently , so Eq. (6.16) can be approximated byThe scalar potential yieldsInflation usually requires a number of e-foldings –65, which is defined asApplying the slow-roll approximation to the above relation (6.19), the number of e-foldings yieldswhere . For the model analysed here, this expression cannot be solved analytically even by taking some approximations. Hence, in order to illustrate the powerful of the model, let us consider an example for the parameter of the model that satisfy the conditions (6.2):Note that the value corresponds to Starobisnky inflation. Then, by integrating (6.20), the following results are obtained:As shown above in (6.14), these values lie within the allowed ranges provided by Planck, such that the model (6.1) can reproduce well inflation and then recover late-time acceleration, leading to a unified description of the universe evolution.

Let us now analyse the system of equations (2.9) during the inflationary epoch, which are reduced toHere , , , .

In the top-left panel of Fig. 6, the F(R) function (6.1) is depicted for the valuesand compared to the exponential gravity (3.1) and the CDM model. As shown, the two F(R) models are not distinguishable at low curvature regimes but they are when the curvature becomes very large, as during inflation. Moreover, the model (6.1) practically coincides with the CDM model in the range (or ), while differs for because of the term. Finally for , both F(R) models (3.1) and (6.1) behave similarly and deviate from CDM.

Fig. 6

Plots of F(R) (top-left), (top-right), and (the bottom-right panel) for the models: (3.1), (6.1) and the CDM model with the parameters from Table 3. In the bottom-left panel the level lines of are shown for the model (6.1)

One can conclude that the difference between the models (3.1) and (6.1) is negligible during the radiation/matter era, so their observational predictions coincide. In order to illustrate such behaviour for both F(R) models as well as for the CDM model, the top-right panel in Fig. 6 depict the evolution of the scale factor , where is dimensionless time. Here the differences between the F(R) models are not visible, because it takes place for , but the late time and future behaviour of the CDM model differs from that for the F(R) models.

The early-time evolution and the inflationary epoch ( and ) of the model (6.1 in comparison with the CDM model are shown in the bottom-right panel of Fig. 6. Here the dashed and dash-dotted lines correspond to the Hubble parameter and dots describe the Ricci scalar in their normalized form (3.8) as functions of the scale factor.

The inflationary solution is unstable [47]: after e-foldings the de Sitter solution decays and the evolution transforms into the CDM behaviour.

Conclusions

Exponential gravity may be considered as an alternative to the so-called concordance model. The main gravitational action studied along this manuscript and described in (3.1) represents an slightly correction to the usual Hilbert–Einstein action with a cosmological constant. Such correction is modelled by a free parameter, which has been called , such that it controls the scale at which corrections to GR become important. As shown in some previous works [39, 47], such a type of F(R) models can reproduce a late-time acceleration epoch. However, the aim of this work was to show in an accurate way whether such type of gravitational actions fulfil the necessary cosmological constraints leading to a reasonable bound on the crucial parameter . Note that an equivalent analysis was performed in [50] and more recently in [51]. In our paper we have updated the constraints obtained in previous works by assuming new released data from the last years. In addition, we have also analysed the exponential model at late times by following two approaches: the first one by assuming the condition (5.1) and the second one by assuming a more general approach. The former provides a more restricted case as we are forcing the model to mimic CDM at the present time, while the latter keeps the model free, except for the initial conditions which are the same for the two cases. As expected, the more restrictive on the theoretical conditions are, the better the constraints on the free parameter turn out, as shown in the bottom panels of Fig. 2. On the other hand, whereas the parameter is not upper bounded (recall that GR is recovered for ), the fits realised in the paper, where SNe Ia, BAO, CMB and H(z) data were used, provides a sufficiently small lower bound on , which may have consequences at the perturbation level, an aspect that should be studied in the future. In addition, the values for the matter density do not differ too much among the one given by the CDM model and the one from the exponential gravity model, either when some approximations are assumed or when the general case is considered. Moreover, the is a bit smaller for the exponential gravity case than in CDM, such that the best fit does not correspond to CDM, although the difference is not statistically significant. From a qualitative point of view, the tiny differences among exponential gravity and CDM can be shown by looking at the form of the action in (3.1). By a sufficient small exponent, i.e. a parameter large enough, the Lagrangian mimics quite well the Hilbert–Einstein action with a cosmological constant, such that there is not significant difference on the cosmological evolution among both models, as shown in Fig. 4.

In addition, the exponential gravity action (3.1) can be extended in such a way that the new action (6.1) can reproduce inflation as well. Note that inflationary models within the framework of F(R) theories have been widely studied in the literature, as shown by one of the most successful inflationary models, the so-called Starobinsky inflation [42]. In this sense, we have proposed here a model where the exponential term responsible for the late-time acceleration is suppressed at the large curvature regime and consequently its induced cosmological constant while two additional terms may become important: a different effective cosmological constant (much larger than ) and a power term . Recalling that Starobinsky inflation is described by a term, the inclusion of just generalised the latter and ensures a successful exit from inflation [47]. However, as shown in some previous papers, such exponent has to be in order to avoid instabilities [47]. Here, we have extended the previous analysis by using the usual scalar–tensor counterpart and obtaining the explicit form of the potential for the scalaron. Then we computed explicitly the spectral index of curvature fluctuations during inflation. An example fulfilling all the requirements provides an spectral index that leads to a nearly invariant power spectrum and an almost null amplitude for the tensor modes, predictions in agreement with the last data released by the Planck collaboration. Note that the constant establishes the energy scale at which the last terms in the action (6.1) become important, such that then the action also provides a quasi-de Sitter inflationary expansion, similar to Starobinsky model, where the term guarantees a successful exit from inflation. As shown, the values for the free parameters which provide the correct values for the spectral index and the scalar-to-tensor ratio, also avoid further corrections when inflation ends. Such model is then able to reproduce inflation and successfully exit.

Hence, we can conclude that the full gravitational Lagrangian (6.1) is capable of reproducing inflation and late-time acceleration in such a way that no other fields are required. Recently, one more extension of this type of exponential gravity with log-corrected term was proposed to explain the unified universe history (see Ref. [89]). As shown here, the Lagrangian satisfies the observational constraints, with a no statistical significant difference with respect to the CDM model, what means that they cannot be distinguished. Simultaneously, the F(R) model (6.1) provides the right predictions for the inflationary epoch. The next steps should be focussed on the analysis of cosmological perturbations and the possible effects of such Lagrangian at the astrophysical level.