Cosmology of modified gravity with a non-local $f(R)$
I. Dimitrijevic\footnote{Faculty of Mathematics
University of Belgrade
Studentski trg 16
Belgrade
Serbia
{\tt ivand@matf.bg.ac.rs
zoranrrr@gmail.com
We consider a modification of GR with a special type of a non-local $f(R)$. The structure of the non-local operators is motivated by the string field theory and $p$-adic string theory. We pay special account to the stability of the de Sitter solution in our model and formulate the conditions on the model parameters to have a stable configuration. Relevance of unstable configurations for the description of the de Sitter phase during inflation is discussed. Special physically interesting values of parameters are studied in details.
1509.04254.tex
\documentclass[a4paper,12pt,titllepage]{article}

%\usepackage{jcappub} % for details on the use of the package, please
                     % see the JCAP-author-manual

\usepackage{verbatim}
\usepackage{amsmath}
\usepackage{amssymb,amsfonts}

\newcommand{\arXiv}[1]{\href{http://www.arXiv.org/abs/#1}{#1}}
\usepackage[colorlinks=true, linkcolor=blue, bookmarks=true]{hyperref}

\newcommand{\be}{\begin{equation}}
\newcommand{\ee}{\end{equation}}
\newcommand{\bea}{\begin{eqnarray}}
\newcommand{\eea}{\end{eqnarray}}
\newcommand{\pd}{\partial}
\DeclareMathOperator{\Hess}{Hess}

\newcommand{\Hc}{\mathcal{H}}
\newcommand{\Gc}{\mathcal{G}}
\renewcommand{\Re}{\mathrm{Re \,}}
\newcommand{\const}{\mathrm{const}}


\newcommand{\dx}{d^4x}
\newcommand{\ii}{\mathrm{i}}
\newcommand{\Real}{\mathbb{R}}
\newcommand{\Natural}{\mathbb{N}}

\newcommand{\pmu}{\partial_{\mu}}
\newcommand{\pnu}{\partial_{\nu}}
\newcommand{\nmu}{\nabla_{\mu}}
\newcommand{\nnu}{\nabla_{\nu}}

\newcommand{\FF}{\mathcal{F}}
\newcommand{\GG}{\mathcal{G}}
\newcommand{\HH}{\mathcal{H}}

\DeclareMathOperator{\ch}{cosh}
\DeclareMathOperator{\sh}{sinh}
\DeclareMathOperator{\tgh}{tanh}
\DeclareMathOperator{\cth}{coth}


\begin{document}
%\begin{titlepage}

\title{Cosmology of modified gravity with a non-local $f(R)$}


\author{I. Dimitrijevic\footnote{Faculty of Mathematics, University of Belgrade, 
Studentski trg 
16,Belgrade, Serbia, {\tt ivand@matf.bg.ac.rs, zoranrrr@gmail.com}},
B. Dragovich\footnote{Institute of Physics, University of Belgrade,
and   Mathematical Institute SANU, Belgrade, Serbia, {\tt dragovich@ipb.ac.rs}}, 
J. Grujic\footnote{Teacher 
Education Faculty, University of Belgrade, Kraljice 
Natalije 43,Belgrade, Serbia, {\tt jelenagg@gmail.com}}, A.S. 
Koshelev\footnote{Departamento  de F\'isica and Centro  de  Matem\'atica  e 
Aplica\c c\~oes,  Universidade  da  Beira  Interior,  6200  Covilh\~a, 
Portugal, {\tt alexey@ubi.pt}}, Z. Rakic$^*$}%\\


\date{~}

\maketitle

\begin{abstract}
We consider a modification of GR with a special type of 
a non-local $f(R)$. The structure of the non-local operators is motivated by 
the  string field theory and $p$-adic string theory. We pay special account to 
the stability of the de Sitter solution in our model and formulate the 
conditions on the model parameters to have a stable configuration. 
Relevance of unstable configurations for the description of the 
de Sitter phase during inflation is discussed. Special physically interesting 
values of parameters are studied in details.
\end{abstract}

%\end{titlepage}



\clearpage

\tableofcontents

\section{Introduction}

In 2015, the gravity community  marks the first century of the General
Relativity (GR), which is viewed as one of the most beautiful
and profound physical theories \cite{wald}. GR is still an acting theory of
gravity usually  presented by Einstein's equation of motion for the
gravitational (metric) field $g_{\mu\nu}$:
\begin{align} R_{\mu\nu} - \frac{1}{2} R g_{\mu\nu} = {8 \pi G} T_{\mu\nu} , \label{eq:1.1}
\end{align}
where $R_{\mu\nu}$ is the
Ricci tensor, $R$ is the Ricci scalar, $T_{\mu\nu}$ is the
energy-momentum tensor of matter, $G$ is the Newtonian constant and the speed 
of light is taken $c = 1$. This
Einstein's equation
can be derived from the Einstein-Hilbert action
\begin{align} S = \frac{1}{16\pi G} \int \sqrt{-g}\, R \, d^4x + \int \sqrt{-g} \mathcal{L}_m \, d^4x ,  \label{eq:1.2}
\end{align}
where $g = \det(g_{\mu\nu})$ and $\mathcal{L}_m$ is the Lagrangian of matter.
The signature of the metric is $(-+++)$ and the dimensionality is manifestly
$4$.

GR has been well tested and confirmed in the Solar system, and
it serves as a theoretical laboratory for gravitational investigations at other
spacetime scales. It has important astrophysical implications predicting
existence of black holes, gravitational lensing and gravitational waves. In
cosmology, GR predicts existence of about $95 \%$
of additional matter, which makes dark side of the universe. Namely, if GR is
the gravity theory for the universe as a whole and if the Universe has the
Friedmann-Lema\^{\i}tre-Robertson-Walker (FLRW) metric (which is homogeneous
and isotropic) at the cosmic scale, then there is about $68\%$ of {\it dark
energy}, $27\%$ of {\it dark  matter}, and only about $5\%$ of {\it visible
matter}  in the universe \cite{planck}.


Despite of remarkable phenomenological achievements and many nice theoretical
properties, GR is not a complete
theory of gravity. It has well known long standing problems both in UV and IR 
regimes. In UV or at short distances GR predicts singularities like the Big 
Bang or black hole ones. In particular, under rather general
conditions, GR contains cosmological solutions which lead to an infinite matter
density at the beginning of the Universe \cite{bordeguthvilenkin}.
When  physical theory contains singularity, it is an evident indication  that
around  it such a theory has to be appropriately  modified. In IR pure GR does 
not provide a decent explanation of the Dark Energy phenomenon since a possible 
cosmological term has an absolutely unnatural and unexplainable tiny value. 
From just a theoretical point of view GR is not a renormalizable theory even if 
being quantized.
The ways to modify GR  usually  come from some more general theories like 
quantum
gravity, string theory as well as astrophysical and
cosmological observations (for a review, see \cite{modreview}).
Unfortunately there is no so far solid fundamental physical principle which 
could 
tell 
us how to find appropriate modification between infinitely many possible 
theoretical constructions.

In the present paper we consider a non-local modification of GR which has 
strong motivation from string field theory (SFT) \cite{ourreview} and $p$-adic 
string theory \cite{padic}. Such a modification was initially introduced in 
\cite{Biswas:2005qr} and is intensively studied recently 
\cite{Biswas:2011ar,Biswas:2012bp,koshelev,Modesto,dragovich-d,anupamnew}. 
After developing the general ideas regarding this modified theory we mainly 
focus on the properties of a de Sitter solution. De Sitter represents one of the 
most important regime in the Universe evolution. At early times de Sitter phase 
with a large expansion rate lead to the rapid inflation. Nowadays slow expansion 
is attributed to the Dark Energy. The principal generalization w.r.t. previous 
papers is a replacement of the non-local term $R\FF(\Box)R$ in the Lagrangian, 
where $\Box$ is the d'Alembertian operator, with a more general one
$$P(R)\FF(\Box)Q(R).$$
The main motivation to go for this more involved structure is the fact that 
just $P=Q=R$ term does not have a stable de Sitter phase w.r.t. linear 
perturbations as will be explained in the next Sections. It is not a problem 
for the de Sitter phase during inflation as it is exactly what is needed to 
exit the inflation. On the other hand, an attempt to describe the present day 
Universe expansion clearly requires a stable configuration and therefore 
demands further model altering.
 It is also
known that a modification with higher order derivatives  improves the situation
with the renormalizability \cite{Stelle,Talaganis:2014ida} during the 
quantization of the model.

In Section 2 we present the model and derive equations of motion. In Section 3 
we develop some general ideas of solution construction and discuss some of 
them. In Section 4 we turn to the de Sitter solution and formulate the 
stability condition for linear perturbations. In Section 5 we explore the 
stability condition in greater details and put down concluding remarks in 
Section 6.

\section{Non-local generalization of GR and equations of motion}

We consider the following  nonlocal gravity action
\begin{equation} \label{lag:1}
S =  \int \dx\sqrt{-g}\left(\frac{M_P^2}2R-\Lambda + \frac{\lambda}2P(R)
\FF(\Box)Q(R) \right) ,
\end{equation}
where $R$ is the scalar curvature, $\Lambda$ is the cosmological constant, $
\mathcal{F}(\Box)= \displaystyle \sum_{n =0}^{\infty} f_{n}\Box^{n}$  is an
analytic function  of the d'Alembert-Beltrami operator $\Box =
\nabla^\mu\nabla_\mu$ where $\nabla_\mu$ is the covariant derivative. The
Planck mass $M_P$ is related to the Newtonian constant $G$ as
$M_P^2=\frac1{8 \pi G}$ and $P$,$Q$ are scalar functions of the scalar
curvature. The dimensionality is fixed to $4$ and our signature is $(-,+,+,+)$.
$\lambda$ is a constant and in principle can be absorbed in the
rescaling of $\FF(\Box)$. However, it is a convenient tool to track the GR
limit which is $\lambda\to 0$.
As it is obvious we are going down the way of generalizing the results obtained 
earlier in the case $P=Q=R$ \cite{Biswas:2012bp}. In the sequel we shall omit 
an explicit further citation of these results referring them rather as 
the \textit{non-local $R^2$} case.

To have physically meaningful expressions and to keep track of the 
non-localities one should introduce the scale of non-locality using a new mass 
parameter $M$. Then the function $\FF$ would be expanded in Taylor series as
$\mathcal{F}(\Box)= \displaystyle \sum_{n =0}^{\infty} \bar 
f_{n}\Box^{n}/M^{2n}$ with all barred constants dimensionless. We will return 
to these notations during the discussion of our results.

Varying the action (\ref{lag:1}) by substituting
\be g_{\mu\nu}\to g_{\mu\nu}+h_{\mu\nu}\label{maindelta}\ee to
the linear order in $h_{\mu\nu}$, removing the total derivatives and
integrating from time to time by parts one gets
\begin{equation} \begin{aligned} \label{EOMPQgravity}
\delta S=&\int
\dx\sqrt{-g}\frac{h^{\mu\nu}}2\bigg[-\tilde G_{\mu\nu}\bigg],
\end{aligned} \end{equation}
where
\begin{equation} \begin{aligned} \label{EOMPQ}
-\tilde G_{\mu\nu}&\equiv-M_P^2G_{\mu\nu}-g_{\mu\nu}\Lambda\\
&+\frac{\lambda}2g_{\mu\nu}P\FF(\Box) Q-\lambda(R_{\mu\nu}-K_{\mu\nu})V\\
&+\frac{\lambda}2\sum_{n=1}^\infty f_n\sum_{l=0}^{n-1}\left(P_\mu^{(l)}
Q_\nu^{(n-l-1)}+P_\nu^{(l)}
Q_\mu^{(n-l-1)}-g_{\mu\nu}(g^{\rho\sigma}P_\rho^{{(l)}}
Q_\sigma^{(n-l-1)}+P^{(l)} Q^{(n-l)})\right).
\end{aligned} \end{equation}
Here $G_{\mu\nu}=R_{\mu\nu}-\frac12g_{\mu\nu}R$ is the Einstein tensor, $
K_{\mu\nu} = \nmu\nnu - g_{\mu\nu} \Box $,
 $V = P_R\FF(\Box)Q + Q_R\FF(\Box)P$ where the subscript $R$ indicates the
derivative w.r.t. $R$ (as many times as it is repeated) and
$$
P^{{(l)}}=\Box^l P,~P_\rho^{{(l)}}=\pd_{\rho}\Box^l P\text{ with the same
for }Q,~P_R,~\dots
$$
Provided a matter source as well the full equations of motion read
\begin{equation} \begin{aligned} \label{EOMPQmatter}
\tilde G_{\mu\nu}=T_{\mu\nu}.
\end{aligned} \end{equation}

Analyzing (\ref{EOMPQ}) we recognize that the first line is the canonical 
EOM for the Einstein's
GR with the cosmological constant, the second line with $\FF(\Box)=1$
represents the extension to the local $f(R)$ type gravities while a
non-constant $\FF(\Box)$ as well as the last line are unique for a higher
derivative (probably non-local) modification of gravity.
The trace equation is of use and we write it separately
\begin{equation} \begin{aligned} \label{tracePQ}
&M_P^2R-4\Lambda + 2\lambda P \FF(\Box)Q -\lambda(R+3\Box)V  \\
-&{\lambda}\sum_{n=1}^{\infty} f_n
\sum_{l=0}^{n-1}\left(g^{\rho\sigma}P_\rho^{{(l)}}
Q_\sigma^{(n-l-1)}+2P^{(l)} Q^{(n-l)}\right)   = -T^\mu_\mu.
\end{aligned} \end{equation}

If either $P$ or $Q$ is a constant then effectively $\FF(\Box)=f_0$ and does
nothing with the non-local theories of interest.
Note
that thanks to the integration by parts there is always the symmetry of an
exchange $P\leftrightarrow Q$.

\section{Solutions construction}

\subsection{Cosmological FRW solutions}

In this paper the primary goal is to attack the cosmological properties of the 
proposed gravity model. This brings us to the cosmologically important metrics 
from which we focus on the Friedmann-Robertson-Walker (FRW) configurations. The 
latter have the following metric
\be
ds^2=-dt^2+a^2(t)\left(\frac{dr^2}{1-Kr^2}
+r^2d\theta^2+r^2\sin^2\theta d\varphi^2\right),\label{frwk}
\ee
where $t$ is the cosmic time, $a(t)$ is the scale factor and $K$ is the spatial 
curvature. Such metrics represent a homogeneous and isotropic universe. The 
resulting tensor $\tilde G_{\mu\nu}$ turns out to be diagonal with two distinct 
components: $(00)$ and $(ii)$ for some spatial index $i=1,2,3$.

The above chosen form of the metric implies that only the matter sources of the 
form $T^\mu_\nu=\mathrm{diag}(-\rho,p,p,p)$ are compatible with the equations 
of motion (\ref{EOMPQmatter}). Here $\rho$ is the energy density and $p$ 
is the pressure density. Such matter sources represent a wide class of 
physically 
important cases and include the most crucial set of perfect fluids which can be 
written as $T_{\mu\nu}=(\rho+p)u_\mu u_\nu+pg_{\mu\nu}$ where $u_\mu$ is the 
4-velocity of the fluid.

It is common that gravity is a theory with constraints thanks to Bianchi 
identities. These identities manifest in the relation $\nabla^\mu 
G_{\mu\nu}\equiv0$. However, one can check that given a generally covariant 
action like (\ref{lag:1}) one ends up with a similar relation
$$\nabla^\mu 
\tilde G_{\mu\nu}\equiv0.$$
This implies like in GR the canonical conservation equation for the matter
$$\nabla^\mu 
T_{\mu\nu}=0.$$

In the sequel we will focus on two special cases: no matter at all and a 
radiation source. The first case represents vacua of our model and technically 
just reads as $T_{\mu\nu}=0$. The second case is more subtle and captures many 
solutions of interest \cite{sols}. Indeed, as already mentioned above, the 
metric 
(\ref{frwk}) leaves us with essentially two equations. For simplicity we choose 
the $(00)$ component and the trace equation. Radiation is characterized by 
$w=p/\rho=1/3$ and is therefore traceless. Given a solution to the trace 
equation with zero RHS we almost solve the whole system of equations since the 
Bianchi identity together with the conservation equation imply that the 
remaining $(00)$ component must be satisfied modulo a radiation-like term at 
most. To be physical, such a radiation must have positive energy density.

\subsection{Non-local $R^2$ case}

It is clear from the complexity of expressions (\ref{EOMPQ}) that 
constructing a general solution is a very ambitious hope. However we remind 
that a significant progress has been achieved in the non-local $R^2$ case by 
considering a simplifying ansatz
\begin{equation}\label{ansatzR}
\Box R=r_1R+r_2  
\end{equation}
with $r_{1,2}$ being 
constants.
It is useful to repeat this procedure here in order to illuminate the steps 
related to our model of interest (\ref{lag:1}).

Consider $P=Q=R$. Then 
equations (\ref{EOMPQmatter}) become
\begin{equation} \begin{aligned} \label{EOMPQR}
T_{\mu\nu}&=-M_P^2G_{\mu\nu}-g_{\mu\nu}\Lambda+\frac{\lambda}2g_{\mu\nu}
R\FF(\Box) 
R-2\lambda(R_{\mu\nu}-K_{\mu\nu})\FF(\Box)R\\
&+\frac{\lambda}2\sum_{n=1}^\infty f_n\sum_{l=0}^{n-1}\left(R_\mu^{(l)}
R_\nu^{(n-l-1)}+R_\nu^{(l)}
R_\mu^{(n-l-1)}-g_{\mu\nu}(g^{\rho\sigma}R_\rho^{{(l)}}
R_\sigma^{(n-l-1)}+R^{(l)} R^{(n-l)})\right).
\end{aligned} \end{equation}
Application of ansatz (\ref{ansatzR}) means
\begin{equation}
\Box^nR=r_1^nR+r_2r_1^{n-1},~n>0\text{ and 
}\FF(\Box)R=\FF_1R+\FF_2,
  \label{ansatzRn}
\end{equation}
where $\FF_1=\FF(r_1)$ and $\FF_2=\frac{r_2}{r_1}(\FF(r_1)-f_0)$.

As explained in the previous Subsection we should proceed with the trace 
equation. For the traceless matter (or for no matter at all) it becomes
\begin{equation}
A_1R-\lambda\FF'(r_1)\left(2r_1R^2+ \partial_\mu R \partial^\mu 
R\right)+A_2=0,
\label{eqEinsteinRonlytraceansatz2}
\end{equation}
whith
\begin{equation*}
\begin{aligned}
A_1&={}M_P^2-\lambda\left(4\FF'(r_1)r_2-2\frac{r_2}{r_1}
(\FF_1-f_0)+6\FF_1r_1\right),\\
A_2&={}-4\Lambda-\lambda\frac{r_2}{r_1}\left(2\FF'(r_1)r_2-2\frac{r_2}{r_1}
(\FF_1-f_0)+6\FF_1r_1\right).
\end{aligned}
\end{equation*}
The above equation is solved provided $A_1=A_2=0$, and \be
\FF'{(r_1)}=0\ .
\label{r1}
\ee
Simple algebra gives
\begin{equation}
 r_2={}-\frac{r_1[M_P^2-6\lambda\FF_1r_1]}{2\lambda[\FF_1-f_0]},
\qquad
\Lambda={}-\frac{r_2M_P^2}{4r_1}. \label{r2lambda}
\end{equation}
After this one is left with the problem of positivity of a possible radiation 
energy density.

\subsection{Ansatz for (\ref{EOMPQ})}

Now we turn to our more generic model. Dramatic simplifications are seen 
provided two conditions are satisfied:
\begin{equation}
  P=Q\text{ and }\Box P=p_1P,
\end{equation}
where $p_{1}$ is a constant.
Application of this ansatz means
\begin{equation}
\Box^nP=p_1^nP,~n>0\text{ and 
}\FF(\Box)P=\FF_{p1}P,
  \label{ansatzRnp}
\end{equation}
where $\FF_{p1}=\FF(p_1)$.
The further explicit substitution in the trace equation (\ref{tracePQ}) yields 
(with $T^{\mu}_{\mu}=0$)
\begin{equation} \begin{aligned} \label{tracePQansatzp}
&M_P^2R-4\Lambda + 2\lambda \FF_{p1}P^2 
-2\lambda(R+3\Box)P_R\FF_{p1}P  -{\lambda}\FF'(p_1)(g^{\rho\sigma}P_\rho
P_\sigma+2p_1P^2)   = 0.
\end{aligned} \end{equation}

Analogously to the non-local $R^2$ case a condition $\FF'(p_1)=0$ 
simplifies the clutter a lot. Then one can see that a specially simple case 
arises if $P=\sqrt{R+R_0}$. For such a choice of function $P$ one gets 
$P_R\sim1/P$ and as the result the latter equation can be solved by adjusting 
the parameters. However, there is still a non-trivial equation to be satisfied
\be
\Box\sqrt{R+R_0}=p_1\sqrt{R+R_0}.\label{sqrt}
\ee
Nevertheless, presence of at least some solutions can be checked numerically. 
Interestingly, such a choice of $P$ has no a known local counterpart as it 
would give just a canonical Einstein-Hilbert term with a cosmological constant.

\subsection{Constant curvature solutions}

This is a special class of solutions which include several important cases.

Substituting $R=\const$ into the trace equation
(\ref{tracePQ}) and using the fact that $P,Q,V$ are now constants one gets
\begin{equation} \begin{aligned} \label{tracePQdS}
&M_P^2R-4\Lambda + 2\lambda P f_0Q -\lambda RV
=0.\end{aligned} \end{equation}
Here $V$ reduces to $(P_RQ+Q_RP)f_0$.
Solving the latter equation is an algebraic rather than differential 
problem.

The very important case is the de Sitter solution as far as it plays a 
crucial role in the description of our Universe. In four dimensions it is 
characterized by
\be
R_{\mu\nu}=\frac R4 g_{\mu\nu}\text{ and } R=\const>0.\label{dS}\ee
Notice that
just a constant $R$ does not mean the space-time
is de Sitter, however. One can check this is a vacuum solution of our model, 
i.e. no matter is needed to support it.

\section{Cosmological expansion}

\subsection{Background}

In this section we look for the cosmological spatially flat de Sitter which can 
be written as
\be
ds^2=-dt^2+e^{2Ht}d\vec x^2
\label{dScosmo}
\ee
with a constant $H$, $t$ the cosmic time and the vector is the 3-dimensional
notion.
This is a particular case of a spatially flat FRW metric
\be
ds^2=-dt^2+a(t)^2d\vec x^2
\label{FRW}
\ee
with $a(t)=a_0\exp(Ht)$. The general definition of $H$ is $H=\dot a/a$ with dot
being the derivative w.r.t. the cosmic time.

Some relevant background quantities (for a general $a$) are
\be
R=12H^2+6\dot H ,~\Gamma_{ij}^0=Hg_{ij},~\Gamma_{j0}^i=H\delta^i_{j},~
 \Box=-\pd_t^2-3H\pd_t+\frac{\delta^{ij}\pd_i\pd_j}{a^2},
\ee
where the indexes $i,j$ range as $1,2,3$.
On the background all quantities are space homogeneous as the metric suggests.

For perturbations we employ the canonical ADM $(1+3)$ decomposition and
introduce the conformal time $\tau$ such that
$$
ad\tau=dt.
$$
Then the general FRW metric
(\ref{FRW}) transforms to
\be
ds^2=a(\tau)^2(-d\tau^2+d\vec x^2).
\label{dScosmoconf}
\ee
For the de Sitter background (\ref{dScosmo})
$$
\tau=-\frac1{a_0H}e^{-Ht}~\Rightarrow~a(\tau)=-\frac1{
H\tau}.
$$
So when $t$ goes from past to future infinity, $\tau$ goes from $-\infty$ to
$0_-$. $t=0$ corresponds to $\tau=-\frac1{a_0H}$.

\subsection{Scalar perturbations}

The metric for the scalar perturbations around a FRW background is defined as
\begin{equation}
ds^2=a(\tau)^2\left[-(1+2\phi)d\tau^2-2\pd_i \beta d\tau
dx^i+((1-2\psi)\delta_{ij}+2\pd_i\pd_j\gamma)dx^idx^j\right].
\label{alex_scalar_pert}
\end{equation}
Out of 4 scalar modes only 2 are gauge invariant.
The convenient gauge invariant variables (Bardeen potentials) are introduced as
\begin{equation}
\Phi=\phi-\frac{1}{a}(a\vartheta)^\prime=\phi-\dot
\chi,\qquad\Psi=\psi+\Hc\vartheta=\psi+H\chi,
\label{GIvars}
\end{equation}
where  $\chi=a\beta+a^2\dot\gamma$, $\vartheta=\beta+\gamma'$, $\Hc(\tau)=a'/a$.
The prime denotes the differentiation with respect to the conformal time
$\tau$ and the dot as before w.r.t. the cosmic time $t$.

The $(1+3)$ structure suggests
to represent the perturbation quantities (which can depend on all 4
coordinates) as
\be
f(\tau, \vec x)=f(\tau,k)Y(k,\vec x),
\label{scalar13}
\ee
where $k=|\vec k|$ comes from the definition of the $Y$-functions as spatial
Fourier modes
\be
\delta^{ij}\pd_i\pd_j Y=-k^2Y.
\ee
Obviously
\be
Y=Y_0e^{\pm i\vec k\vec x}.
\ee
The relevant expressions for the d'Alembertian operator are
\be\Box=-\frac{1}{a^2}\pd_\tau^2-2\frac{a'}{a^3}\pd_\tau-\frac{k^2}{a^2}
=-\partial_t^2-3H\partial_t-\frac{k^2}{a^2}.
\ee
Recall that considering the de Sitter background we do not need the perturbed
d'Alembertian operator to the linear order in perturbation analysis.
However, all the expressions in this subsection are valid for a generic scale
factor $a$.

\subsection{Perturbations}

Perturbations of equations (\ref{EOMPQ}) to the linear order around the dS
vacuum are easy to compute since many terms drop out. What remains after a
careful computation is
\begin{equation} \begin{aligned} \label{deltaEOMPQ}
-&m^2\delta
G^\mu_{\nu}+(R^\mu_{\nu}-K^\mu_{\nu})v(\Box)\delta R=0,
\end{aligned} \end{equation}
where $m^2=M_P^2+\lambda f_0(P_RQ+Q_RP)$ and
$v(\Box)=-\lambda(P_{RR}Q+Q_{RR}P)f_0+2P_RQ_R\FF(\Box)$.
We have used the fact that variation of the $\Box$ acting on a scalar
function is a pure differential operator. Indeed,
varying the box we have
\be(\delta\Box)
f=[-h^{\mu\nu}(\pd_\mu\pd_\nu-\Gamma^\rho_{\mu\nu}\pd_\rho)-g^{\mu\nu}
\gamma^\rho_{\mu\nu}\pd_\rho]f.\ee
Here $\Gamma^\rho_{\mu\nu}$ is the Christoffel symbol.
Also we have used that starting from (\ref{maindelta}) one gets
\be
\delta g^{\mu\nu}=-h^{\mu\nu}
,~\delta \Gamma^\rho_{\mu\nu}=\gamma^\rho_{\mu\nu}=\frac12(\nabla_\mu
h^\rho_\nu+\nabla_\nu h^\rho_\mu-g^{\rho\sigma}\nabla_\sigma h_{\mu\nu}).
\ee
So if $f$ is a constant then $(\delta \Box) f=0$.
The same is true for $K^\mu_\nu$
\be
(\delta K^\mu_\nu)
f=[-h^{\mu\sigma}(\pd_\sigma\pd_\nu-\Gamma^{\rho}_{\sigma\nu}
\pd_\rho)-g^{\mu\sigma}\gamma^\rho_ { \sigma\nu }
\pd_\rho-\delta^\mu_\nu\delta\Box]f,
\ee
which is zero as long as $f$ is a constant.
Hence all the corresponding terms vanish.

Taking the trace of (\ref{deltaEOMPQ}) we get
\begin{equation} \begin{aligned} \label{deltatracePQ}
[m^2+(R+3\Box)v(\Box)]\delta
R=\GG(\Box)\delta R=0.
\end{aligned} \end{equation}
It is a homogeneous equation on $\delta R$. The general method of solving it
is to use the Weierstrass factorization
\be
\GG(\Box)\delta R=\prod_i(\Box-\omega_i^2)e^{\gamma(\Box)}\delta R=0,
\ee
where $\omega_i^2$ are the roots of the algebraic (or perhaps
transcendental) equation $\GG(\omega^2)=0$ and $\gamma(\Box)$ is an entire
function (as a consequence $e^{\gamma(\omega^2)}$ has no roots).
We
assume that there are no multiple roots as they complicate the story slightly.
Then we can solve it for each $\omega_i$ separately
\be \label{eqomega}
(\Box-\omega_i^2)\delta R=0.
\ee
The latter is just a second order equation. Let's denote the corresponding
solution as $\delta R_i$.
 The resulting general solution for
$\delta R$ is a sum
\be
\delta R=\sum_i\delta R_i,
\ee
where each $\delta R_i$ has its arbitrary integration constants. The only thing
we care about is that the total $\delta R$ must be real. Note that in the
sequel the $\delta R=0$ case will also play an essential role.

To proceed we need to specify more concretely the  form of the background
metric and its
perturbations.

\subsection{Solution to (\ref{deltatracePQ})}

As explained after equation (\ref{deltatracePQ}) we essentially must solve
equation (\ref{eqomega}). It is a second order equation which can be written
explicitly as
\be
\left(\pd_\tau^2-\frac2\tau\pd_\tau+{k^2}+\frac{\omega_i^2}{H^2\tau^2}
\right)\delta R=0,
\label{eqdRdS}
\ee
where we taken the dS form of the background.
The solution yields
\be
\delta R_i=(-k\tau)^{ 3/2} \left( C_{1i} J_{\nu_i}(-k\tau) + C_{2i}
Y_{\nu_i}(-k\tau)\right),
\label{dRdS}
\ee
where $J,~Y$ are the Bessel functions
of the first and second kinds,
respectively, with $\nu_i=\sqrt{\frac94-\frac{\omega_i^2}{H^2}}$ and
$C_{1i,2i}$ are the integration constants.
The full answer for $\delta R$ is
\be
\delta R=\sum_i\delta R_i.
\label{dRdSfull}
\ee
A possible case $\delta R=0$ is a viable (and not a trivial) one.

\subsection{Bardeen potentials}

What we have to extract are the Bardeen potentials introduced in
\eqref{GIvars}. To do this we need two equations.

One is given by the
formulation of $\delta R$ in terms of $\Phi$ and $\Psi$ accounting that
the time behavior of $\delta R$ itself is found above. A lengthy computation
gives
\begin{equation}
\delta R= \frac{2}{a^2}\left[k^2
(\Phi-2\Psi)-3\frac{a'}{a}\Phi'-6\frac{a''}{a}\Phi-3\Psi''-9\frac{a'}{a}
\Psi'\right].
\label{deltaRGI}
\end{equation}
This expression is valid for a general $a$ while its dS form is
\begin{equation}
\delta R = -6H^2 \left(4\Phi - \tau(\Phi' + 3\Psi') + \tau^2 \Psi''\right)
+ 2\tau^2 H^2 k^2 \left(\Phi -2\Psi\right).
\label{deltaRGIdS}
\end{equation}


Another equation can be
got as the $i\neq j$ component of the system (\ref{deltaEOMPQ}). To derive the
latter we recite here the not trivially vanishing part (using the notations
introduced in the previous Section)
\begin{equation} \begin{aligned} \label{deltaEOMPQij}
-&m^2\delta G^i_{j}-g^{ik}\pd_k\pd_jv(\Box)\delta
R=0.
\end{aligned} \end{equation}
This upon the explicit sibstitution can be rewritten as
\begin{equation} \begin{aligned} \label{deltaEOMPQijfinal}
-&m^2(\Phi-\Psi)+v(\Box)\delta R=0.
\end{aligned} \end{equation}

\subsection{Solutions for $\Phi$ and $\Psi$}

The equations to solve are (\ref{deltatracePQ}), (\ref{deltaRGIdS}) and
(\ref{deltaEOMPQijfinal}).
Assuming $m^2\neq 0$ and
combining (\ref{deltaEOMPQijfinal}) and (\ref{deltatracePQ}) we yield
\bea
\delta R+(R+3\Box)(\Phi-\Psi)&=&0.\label{eqsimple11}
\eea
The latter equation can be simplified by substituting $\delta
R$ from (\ref{deltaRGIdS}), $\Box$ and
$R=12H^2$ into (\ref{eqsimple11}). This yields a homogeneous equation on
$\Phi+\Psi$
\be
(\tau^2\pd_\tau^2-4\tau\pd_\tau+4+\frac{k^2\tau^2}{3})(\Phi+ \Psi) = 0
\label{phipluspsi}
\ee
and its solution is
\be
\Phi+\Psi=
\eta[c_1(\cos(\eta)+\eta\sin(\eta))+
c_2(-\eta\cos(\eta)+\sin(\eta))]\ ,
\label{phipluspsisol}
\ee
where
$\eta=\frac{k\tau}{\sqrt{3}}
$.

Next, since $\delta R$ is provided independently by (\ref{dRdSfull}) one
readily expresses from (\ref{deltaEOMPQijfinal})
\be
\Phi-\Psi=\frac1{m^2}\sum_i v(\omega_i^2)\delta R_i.\label{phiminuspsisol}
\ee
Now we can find each field algebraically by means of equations
(\ref{phipluspsisol}), (\ref{phiminuspsisol}).

There is a special case $m^2=0$ as we loose the possibility to find out the
Bardeen
potentials separately. The compatibility of system
(\ref{deltatracePQ},\ref{deltaEOMPQijfinal}) requires that either $\delta R=0$
or there is
a root of $v(\Box)$. However, neither of equations coming from
(\ref{deltaEOMPQ}) would help in this case as all of them lack information
about individual Bardeen potentials if $m^2=0$.

Physically this reflects the fact that effectively the Einstein-Hilbert term
vanishes and one reduces the number of propagating degrees of freedom.

%\subsection{Limits}

Note that $\Phi$ and $\Psi$ are bounded if and only if $\Phi+\Psi$ and $\Phi-\Psi$ are bounded.
Let $t\to \infty$, then $\tau \to 0-$.
\be
\lim_{t\to +\infty}(\Phi+\Psi)=\lim_{t\to +\infty} \eta[c_1(\cos(\eta)+\eta\sin(\eta))+
c_2(-\eta\cos(\eta)+\sin(\eta))] = 0,
\ee
\be
\lim_{t\to +\infty}(\Phi-\Psi)=\frac1{m^2}\sum_i v(\omega_i^2)\lim_{t\to +\infty} (-k\tau)^{ 3/2} \left( C_{1i} J_{\nu_i}(-k\tau) + C_{2i}
Y_{\nu_i}(-k\tau)\right).
\ee
For small values of $\tau$ Bessel functions have the following asymptotic
behaviour
\be
J_\nu(z) \sim z^{\Re\nu},~Y_\nu(z) \sim z^{-|\Re\nu|}\text{ for }\Re\nu\neq 0,
Y_\nu(z) \sim \ln z\text{ for }\Re \nu = 0.
\ee
From this we conclude that Bardeen potentials are bounded irrespectively of
the boundary condition provided
\begin{equation}
  |\Re\nu|<\frac32.
  \label{condnu}
\end{equation}


This is in perfect agreement with \cite{Biswas:2012bp} and the comparison should
be done as follows. In that reference a more general class of solutions was 
studied
which asymptote to the de Sitter background at late times while the non-local
part of the Lagrangian is given by
$R\FF(\Box)R$.
In the present paper we stick to the de Sitter solutions from the very
beginning but the non-local part of the Lagrangian is more general. Assuming in
our case $P(R)=Q(R)=R$ we retrieve the results exactly as in Section~4 in
\cite{Biswas:2012bp} where the late time regime is studied.

\section{Understanding condition (\ref{condnu})}

As the main result of the previous Section the question of stability of the de
Sitter vacuum is narrowed to the satisfactory solution to eq. (\ref{condnu}).
$\nu$ in turn depends on the structure of the non-local operator $\GG(\Box)$
such that
\be
\nu=\sqrt{\frac94-\frac{\omega^2}{H^2}}\label{nunu}
\ee
and $\GG(\omega^2)=0$.
Moreover, the requirement that the system does lead to ghosts demands that
there is no more than one such a root $\omega^2$ for the operator
$\GG$.

Tracing back the operator $\GG$ and taking into account that no more than one
root exists one finds
\begin{equation}
\begin{split}
\GG=&m^2+(R+3\Box)v(\Box)\\
=&[M_P^2+\lambda %
f_0(P_RQ+Q_RP)]+(R+3\Box)[-\lambda(P_{RR}Q+Q_{RR}P)f_0+2P_RQ_R\FF(\Box)]\\
=&(\Box-\omega^2)e^{\gamma(\Box)}
\end{split}
\end{equation}
with $\gamma(\Box)$ an entire function. The latter equation can be in principal
inverted to express $\FF(\Box)$ and therefore see how much general it can be.
It is worth to mention that an additional assumption made at the beginning is
that $\FF(\Box)$ itself is an analytic function.

To make things tractable we specialize to monomials
\be
P(R) = R^p, \quad Q(R)=R^q,
\ee
for some nonzero integers $p$ and $q$.
At the background level this results in the following modification of equation
(\ref{tracePQdS})
\begin{equation} \begin{aligned} \label{tracePQdSpq}
&M_P^2R-4\Lambda + \lambda f_0R^{p+q}(2-p-q)
=0.\end{aligned}\end{equation}
This equation can be solved in general w.r.t. $R$ as long as $-3\leq p+q\leq4$.
Also we note the following quantities
\begin{align}
m^2 &=M_P^2+\lambda (p+q)R^{p+q-1} f_0, \\
v(\Box) &=-\lambda R^{p+q-2}((p^2 -p + q^2 -q)f_0+2pq\FF(\Box)).
\end{align}


%\subsection{$p+q\neq0$}

As a useful trick one can analyze 
$\GG(0)$ to see whether the
stability condition can be reached. Indeed, $\GG$ is analytic by 
construction but a compatibility condition must be fulfilled
\be
M_P^2+\lambda R^{p+q-1}(p+q)(2-p-q)f_0=-\omega^2e^{\gamma(0)}.
\label{pqn0}
\ee
It is obvious from (\ref{nunu}) that as long as $\omega^2$ is real it should be
at least positive in order to satisfy (\ref{condnu}). The latter condition
(\ref{pqn0}) clearly shows that $\omega^2$ is indeed real and therefore reduces
to the following necessary inequality
\be
M_P^2+\lambda R^{p+q-1}(p+q)(2-p-q)f_0<0.
\label{pqn01}
\ee
Satisfactory solution to this relation is a necessary stability condition.
From here it is obvious that two special cases, namely $p+q=0$ and $p+q=2$ 
never have a stable de Sitter phase.

In a general situation we have to understand equation 
(\ref{tracePQdSpq}) together with the latter inequality (\ref{pqn01}). 
As a side turn one can simplify {(\ref{pqn01})} using 
the background equation 
(\ref{tracePQdSpq}) to
\be
M_P^2R(p+q-1)>{4\Lambda}(p+q).
\ee
One can see from this representation that for $p+q=1$ (directly related to 
(\ref{sqrt})) one can have a stable de Sitter phase given a negative $\Lambda$. 
The latter condition is not improbable as long as $\lambda f_0<0$.

In attempt to solve the system \eqref{tracePQdSpq} and \eqref{pqn01} 
one can rewrite it as 
\begin{equation}\label{simple}
  1 -s + u = 0, \qquad 1+ u z <0,
\end{equation}
where $s = \frac{4\Lambda}{M_P^2R}, z=p+q, u = \frac{\lambda f_0}{M_P^2} 
R^{z-1}(2-z)$. This latter system looks rather simple but unfortunately does 
not provide immediate new interesting solutions from the physical point of view.


\section{Discussion and conclusion}

The main subject of our investigation was the generalized gravity theory 
represented by (\ref{lag:1}) with the emphasis on the de Sitter solution. This 
theory contains higher derivatives and as such has a potentially dangerous 
property to generate ghosts. However, the ghost-free condition can be clearly 
formulated and is translated to the fact that the operator $\GG(\Box)$ (eq. 
(\ref{deltatracePQ})) has no more than one root. Going further we find that 
simple algebraic conditions \eqref{tracePQdSpq} and \eqref{pqn01} determine
a possibility to have stable configurations. Namely, satisfying the inequality 
\eqref{pqn01} one allows the de Sitter phase to be stable.

Generically one would expect such a non-local model to be a candidate for the 
UV gravity modification. In this case a possible de Sitter phase would be 
attributed to the inflationary stage of the Universe evolution. This way the de 
Sitter phase should not be eternally stable as the inflation should stop. This 
is totally in accord with the study of bouncing Universe models performed in 
\cite{Biswas:2012bp}.

On contrary one can try to use non-local models in an attempt to challenge the 
present day slowly accelerated expansion of the Universe. In this case one has 
to arrange a 
stable evolution for a long period of time to be compatible with the presently 
observed Universe.
Interestingly, for $P=R^p,Q=R^q$ two special cases $p+q=0$ and $p+q=2$ lead to 
no stable de Sitter solutions. The first case is a generalized ``cosmological 
constant'' (since $PQ=1$) but on contrary to just a constant cosmological term 
the stability is lost. The second case is the generalized $R^2$ theory and for 
now is the mainly developed candidate for a renormalizable theory of gravity 
\cite{Talaganis:2014ida}.

One more observed special case is a generalized Hilbert-Einstein term with 
$p=q=1/2$. This leads to a possibility to solve completely the equations of 
motion provided we can find a solution to equation (\ref{sqrt}). However, 
finding a solution turned out to be a complicated problem and we leave it for a 
separate study. In all regimes we can control the order of the non-locality 
scale $M$ using other parameters of the model.

Generically presence or absence of stable de Sitter configurations in model 
(\ref{lag:1}) with monomial $P$ and $Q$ has been reduced to a really simple 
algebraic system (\ref{simple}). In a general situation both stable and 
unstable regimes can be organized based on the parameters of the model.



\section*{Acknowledgments}
This paper was supported in part by  Ministry of 
Education, Science and Technological Development of the Republic of Serbia, 
grant No 174012. AK is supported by the FCT Portugal postdoctoral fellowship 
SFRH/BPD/105212/2014.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\appendix

\begin{thebibliography}{99}

\bibitem{wald}
  R. Wald, \textit{General Relativity} (University of Chicago Press,
1984).

\bibitem{planck}        
  P.~A.~R.~Ade {\it et al.}  [Planck Collaboration],
  %``Planck 2013 results. XVI. Cosmological parameters,''
  Astron.\ Astrophys.\  (2014)
  [\arXiv{1303.5076} [astro-ph.CO]].

\bibitem{bordeguthvilenkin}
  A.~Borde, A.~H.~Guth and A.~Vilenkin,
 % ``Inflationary space-times are incomplete in past directions,''
  Phys.\ Rev.\ Lett.\  {\bf 90} (2003) 151301
  [\arXiv{gr-qc/0110012}].


\bibitem{modreview}  T. Clifton, P. G. Ferreira, A. Padilla and  C. Skordis,  
%Modified gravity and cosmology.
Phys. Rep. \textbf{513}, 1--189 (2012)  [\arXiv{1106.2476v2} [astro-ph.CO]];\\
S. Nojiri and S. D. Odintsov,
%Unified cosmic history in 
%modified gravity: from $F(R)$ theory
%to Lorentz non-invariant models.
Phys. Rep. \textbf{505}, 59--144 (2011)
 [\arXiv{1011.0544v4} [gr-qc]];\\
T. P. Sotirou  and V. Faraoni,
%``$f(R)$ theories of gravity'', Rev. Mod. Phys.
\textbf{82} (2010) 451--497  [\arXiv{0805.1726v4} [gr-qc]].

\bibitem{ourreview}
  I.~Ya.~Aref'eva, D.M.~Belov, A.A.~Giryavets, A.S.~Koshelev, P.B.~Medvedev,
%\textit{Noncommutative Field Theories and (Super)String Field Theories},
\arXiv{hep-th/0111208};\\
 W.~Taylor and B.~Zwiebach,
%  \textit{D-branes, Tachyons, and String Field Theory},
  [\arXiv{hep-th/0311017}].

\bibitem{padic}
V.S.~Vladimirov, I.V.~Volovich, I.V.~Zelenov,
 \textit{$p$-adic
Analysis and Mathe\-matical Physics} (Singapore: World Scientific (1994));\\
B.~Dragovich, A.~Yu.~Khrennikov, S.~V.~Kozyrev and I.~V.~Volovich,
%  \textit{On p-Adic Mathematical Physics},
  Anal.\ Appl.\  {\bf 1} (2009) 1
  [\arXiv{0904.4205} [math-ph]].

\bibitem{Biswas:2005qr}
  T.~Biswas, A.~Mazumdar and W.~Siegel,
%  ``Bouncing universes in string-inspired gravity,''
  JCAP {\bf 0603} (2006) 009
  [\arXiv{hep-th/0508194}].

\bibitem{Biswas:2011ar}
  T.~Biswas, E.~Gerwick, T.~Koivisto and A.~Mazumdar,
%  ``Towards singularity and ghost free theories of gravity,''
  Phys.\ Rev.\ Lett.\  {\bf 108} (2012) 031101
  [\arXiv{1110.5249} [gr-qc]].

  \bibitem{Biswas:2012bp}
  T.~Biswas, A.~S.~Koshelev, A.~Mazumdar and S.~Y.~Vernov,
%  ``Stable bounce and inflation in non-local higher derivative cosmology,''
  JCAP {\bf 1208} (2012) 024
  [\arXiv{1206.6374} [astro-ph.CO]].

\bibitem{koshelev}
%\bibitem{Koshelev:2013ida}
  A.~S.~Koshelev,
%  ``Stable analytic bounce in non-local Einstein-Gauss-Bonnet cosmology,''
  Class.\ Quant.\ Grav.\  {\bf 30} (2013) 155001
  [\arXiv{1302.2140} [astro-ph.CO]];\\
A.~S. Koshelev, S.~Y. Vernov,
%``Cosmological Solutions in Nonlocal Models,"
\arXiv{1406.5887} [gr-qc];\\
%\bibitem{Craps:2014wga} 
  B.~Craps, T.~De Jonckheere and A.~S.~Koshelev,
  %``Cosmological perturbations in non-local higher-derivative gravity,''
  JCAP {\bf 1411}, no. 11, 022 (2014)
  [\arXiv{1407.4982} [hep-th]].

  \bibitem{Modesto}
G. Calcagni,  L. Modesto and P. Nicolini, 
%Super-accelerting bouncing cosmology in assymptotically-free non-local gravity.
 	[\arXiv{1306.5332} [gr-qc]];\\
  L. Modesto, Phys.Rev. D86, 044005 (2012), \arXiv{1107.2403};\\
  L.~Modesto, J.~W.~Moffat and P.~Nicolini,
  %``Black holes in an ultraviolet complete quantum gravity,''
  Phys.\ Lett.\ B {\bf 695}, 397 (2011).

\bibitem{dragovich-d}  B. Dragovich,
%Towards $p$-adic matter in the universe. 
Springer Proceedings in Mathematics and
Statistics \textbf{36}, 13--24 (2013) [\arXiv{1205.4409} [hep-th]].

\bibitem{anupamnew}
  A.~Conroy, A.~Mazumdar, S.~Talaganis and A.~Teimouri,
  %``Non-local gravity in D-dimensions: Propagator, entropy and bouncing 
%Cosmology,''
  \arXiv{1509.01247} [hep-th].
  
 \bibitem{Stelle}
  K. Stelle, Phys.Rev. D{\bf 16}, 953 (1977).

\bibitem{Talaganis:2014ida} 
  S.~Talaganis, T.~Biswas and A.~Mazumdar,
  %``Towards understanding the ultraviolet behavior of quantum loops in 
%infinite-derivative theories of gravity,''
  \arXiv{1412.3467} [hep-th];\\
%  \bibitem{Dona:2015tra} 
  P.~Donà, S.~Giaccari, L.~Modesto, L.~Rachwal and Y.~Zhu,
  %``Scattering amplitudes in super-renormalizable gravity,''
  JHEP {\bf 1508}, 038 (2015)
  [\arXiv{1506.04589} [hep-th]].

\bibitem{sols} 
  I.~Dimitrijevic, B.~Dragovich, J.~Grujic and Z.~Rakic,
  %``Some Cosmological Solutions of a Nonlocal Modified Gravity,''
  Filomat 29:3 (2015), 619-628
  [\arXiv{1508.05583} [hep-th]];\\
  I.~Dimitrijevic, B.~Dragovich, J.~Grujic and Z.~Rakic,
  %``Some Power-Law Cosmological Solutions in Nonlocal Modified Gravity,''
  Springer Proc.\ Math.\ Stat.\  {\bf 111}, 241 (2014).


  \end{thebibliography}

\end{document}
  

1509.04254.extracted_bib
{wald R. Wald, General Relativity (University of Chicago Press, 1984).
1509.04254.extracted_bib
{planck P.~A.~R.~Ade {et al. [Planck Collaboration], %``Planck 2013 results. XVI. Cosmological parameters,'' Astron.\ Astrophys.\ (2014) [\arXiv{1303.5076 [astro-ph.CO]].
1509.04254.extracted_bib
{bordeguthvilenkin A.~Borde, A.~H.~Guth and A.~Vilenkin, % ``Inflationary space-times are incomplete in past directions,'' Phys.\ Rev.\ Lett.\ {90 (2003) 151301 [\arXiv{gr-qc/0110012].
1509.04254.extracted_bib
{modreview T. Clifton, P. G. Ferreira, A. Padilla and C. Skordis, %Modified gravity and cosmology. Phys. Rep. 513, 1--189 (2012) [\arXiv{1106.2476v2 [astro-ph.CO]];\\ S. Nojiri and S. D. Odintsov, %Unified cosmic history in %modified gravity: from $F(R)$ theory %to Lorentz non-invariant models. Phys. Rep. 505, 59--144 (2011) [\arXiv{1011.0544v4 [gr-qc]];\\ T. P. Sotirou and V. Faraoni, %``$f(R)$ theories of gravity'', Rev. Mod. Phys. 82 (2010) 451--497 [\arXiv{0805.1726v4 [gr-qc]].
1509.04254.extracted_bib
{ourreview I.~Ya.~Aref'eva, D.M.~Belov, A.A.~Giryavets, A.S.~Koshelev, P.B.~Medvedev, %Noncommutative Field Theories and (Super)String Field Theories, \arXiv{hep-th/0111208;\\ W.~Taylor and B.~Zwiebach, % D-branes, Tachyons, and String Field Theory, [\arXiv{hep-th/0311017].
1509.04254.extracted_bib
{padic V.S.~Vladimirov, I.V.~Volovich, I.V.~Zelenov, $p$-adic Analysis and Mathe\-matical Physics (Singapore: World Scientific (1994));\\ B.~Dragovich, A.~Yu.~Khrennikov, S.~V.~Kozyrev and I.~V.~Volovich, % On p-Adic Mathematical Physics, Anal.\ Appl.\ {1 (2009) 1 [\arXiv{0904.4205 [math-ph]].
1509.04254.extracted_bib
{Biswas:2005qr T.~Biswas, A.~Mazumdar and W.~Siegel, % ``Bouncing universes in string-inspired gravity,'' JCAP {0603 (2006) 009 [\arXiv{hep-th/0508194].
1509.04254.extracted_bib
{Biswas:2011ar T.~Biswas, E.~Gerwick, T.~Koivisto and A.~Mazumdar, % ``Towards singularity and ghost free theories of gravity,'' Phys.\ Rev.\ Lett.\ {108 (2012) 031101 [\arXiv{1110.5249 [gr-qc]].
1509.04254.extracted_bib
{Biswas:2012bp T.~Biswas, A.~S.~Koshelev, A.~Mazumdar and S.~Y.~Vernov, % ``Stable bounce and inflation in non-local higher derivative cosmology,'' JCAP {1208 (2012) 024 [\arXiv{1206.6374 [astro-ph.CO]].
1509.04254.extracted_bib
{koshelev %
1509.04254.extracted_bib
{Koshelev:2013ida A.~S.~Koshelev, % ``Stable analytic bounce in non-local Einstein-Gauss-Bonnet cosmology,'' Class.\ Quant.\ Grav.\ {30 (2013) 155001 [\arXiv{1302.2140 [astro-ph.CO]];\\ A.~S. Koshelev, S.~Y. Vernov, %``Cosmological Solutions in Nonlocal Models," \arXiv{1406.5887 [gr-qc];\\ %
1509.04254.extracted_bib
{Craps:2014wga B.~Craps, T.~De Jonckheere and A.~S.~Koshelev, %``Cosmological perturbations in non-local higher-derivative gravity,'' JCAP {1411, no. 11, 022 (2014) [\arXiv{1407.4982 [hep-th]].
1509.04254.extracted_bib
{Modesto G. Calcagni, L. Modesto and P. Nicolini, %Super-accelerting bouncing cosmology in assymptotically-free non-local gravity. [\arXiv{1306.5332 [gr-qc]];\\ L. Modesto, Phys.Rev. D86, 044005 (2012), \arXiv{1107.2403;\\ L.~Modesto, J.~W.~Moffat and P.~Nicolini, %``Black holes in an ultraviolet complete quantum gravity,'' Phys.\ Lett.\ B {695, 397 (2011).
1509.04254.extracted_bib
{dragovich-d B. Dragovich, %Towards $p$-adic matter in the universe. Springer Proceedings in Mathematics and Statistics 36, 13--24 (2013) [\arXiv{1205.4409 [hep-th]].
1509.04254.extracted_bib
{anupamnew A.~Conroy, A.~Mazumdar, S.~Talaganis and A.~Teimouri, %``Non-local gravity in D-dimensions: Propagator, entropy and bouncing %Cosmology,'' \arXiv{1509.01247 [hep-th].
1509.04254.extracted_bib
{Stelle K. Stelle, Phys.Rev. D{16, 953 (1977).
1509.04254.extracted_bib
{Talaganis:2014ida S.~Talaganis, T.~Biswas and A.~Mazumdar, %``Towards understanding the ultraviolet behavior of quantum loops in %infinite-derivative theories of gravity,'' \arXiv{1412.3467 [hep-th];\\ %
1509.04254.extracted_bib
{Dona:2015tra P.~Donà, S.~Giaccari, L.~Modesto, L.~Rachwal and Y.~Zhu, %``Scattering amplitudes in super-renormalizable gravity,'' JHEP {1508, 038 (2015) [\arXiv{1506.04589 [hep-th]].
1509.04254.extracted_bib
{sols I.~Dimitrijevic, B.~Dragovich, J.~Grujic and Z.~Rakic, %``Some Cosmological Solutions of a Nonlocal Modified Gravity,'' Filomat 29:3 (2015), 619-628 [\arXiv{1508.05583 [hep-th]];\\ I.~Dimitrijevic, B.~Dragovich, J.~Grujic and Z.~Rakic, %``Some Power-Law Cosmological Solutions in Nonlocal Modified Gravity,'' Springer Proc.\ Math.\ Stat.\ {111, 241 (2014).