On Stability of

Asymptotically Free Mimetic Hořava Gravity

Tobias B. Russ1

__

1 Ludwig-Maxmillians-Universität, Theresienstr. 37, 80333 Munich, Germany

_Contact<EMAIL_ADDRESS>
Abstract

Asymptotically free mimetic gravity has been introduced as a proposal for a
classical limiting curvature theory with the purpose of singularity
resolution. It was found that in a spatially flat universe an initial stage of
exponential expansion with graceful exit is a generic consequence, regardless
of the matter content. In this work I will analyze linear stability of
cosmological perturbations in such a model, considering only the degrees of
freedom of pure mimetic gravity. I show that the addition of Hořava-gravity-
like higher order spatial curvature terms can lift the gradient instability of
scalar perturbations, even when the gradient term has the wrong sign
throughout. Calculating the primordial spectra of tensor and scalar
perturbations in the simplest single component model, I find that the
initially scale invariant spectra turn out to be destroyed later by the
rapidly varying speed of sound at horizon exit.

##  1 Introduction

By the unique status of general relativity, any alternative theory of metric
gravity usually either has to allow for higher derivatives of the metric,
higher dimensions of spacetime or to introduce new fields, separate from the
metric. Another option is to reparametrize the degrees of freedom of the
physical metric itself, e.g. by a disformal transformation [1]. “Mimetic
gravity” stems from the reparametrization of the physical metric $g_{\mu\nu}$
in terms of an auxiliary metric $\tilde{g}_{\mu\nu}$ and the “mimetic field”
$\phi$ as

$g_{\mu\nu}=\tilde{g}_{\mu\nu}\tilde{g}^{\alpha\beta}\phi_{,\alpha}\phi_{,\beta}.$
(1)

This particular disformal transformation is special for two reasons: 1.) It is
singular, explaining how a simple reparametrization can actually lead to new
physics [2], [3], 2.) The physical metric is invariant under Weyl
transformations of the auxiliary metric. This means that the new degree of
freedom introduced by $\phi$ represents what was called a “conformal degree of
freedom of gravity”. Soon after this reparametrization was first introduced in
[4], it was shown in [5] that the mimetic field can be introduced equivalently
as a constrained scalar field, subject to the constraint

$g^{\mu\nu}\phi_{,\mu}\phi_{,\nu}=1.$ (2)

Apart from the dust-like component called “Mimetic Dark Matter” that emerges
as a constant of integration in the modified Einstein equation of mimetic
gravity, the introduction of the mimetic field also enables a wealth of
possible new terms in the gravity action. By breaking shift symmetry in
$\phi$, i.e. introducing a $\phi$ dependent potential, one can produce an
extremely flexible theory where essentially any conceivable background
solution can be realized in a Friedmann universe, cf. [6].

Conversely, if we restrict to shift symmetric theories without higher
derivatives of the metric in the modified Einstein equation, the range of
possibilities for a non-singular universe becomes very narrow. For the most
natural (and arguably only viable) class of modified flat Friedmann universes,
it was shown in [7], [8] that the only thing that can replace the Big Bang
singularity is a smooth transition to a piece of de Sitter spacetime at
limiting curvature. In this work I will show that this class of modifications
also happens to coincide with the class of models that avoid a ghost
instability of scalar metric perturbations.

The mimetic field, by definition (2), provides a global time function whose
gradient is everywhere timelike. In [9] it was shown that this can be used to
covariantly dissect any scalar quantity that is invariant under spatial
diffeomorphisms in the slicing given by $\phi$. In this way it is easy to
write down a Hořava-gravity-like theory with only higher spatial derivatives
but no higher time derivatives or mixed derivatives. In Hořava gravity [10]
such an “asymmetry” between space and time is used to improve the UV behaviour
of the graviton propagator for the purpose of renormalizability. Projectable
Hořava models have been shown to be renormalizable in [11], [12]. Compared to
other covariantized version of Hořava gravity like [13], mimetic Hořava
gravity has the advantage of not having any additional propagating degrees of
freedom in a Minkowski background. Interestingly, also the following
connections between mimetic gravity and Hořava gravity can be drawn: In [14]
it was explored how a dust-like component emerges as a constant of integration
in Hořava-Lifshitz gravity. In [15] an equivalence between the IR limit of
projectable Hořava gravity and a mimetic matter scenario has been found.
Another Hořava-like mimetic model has been presented in [16].

In this paper I will show that higher spatial derivative terms of sixths order
can not only render a power counting renormalizable theory, they also serve to
alleviate the gradient instability of the scalar degree of freedom of mimetic
gravity in an expanding universe.

The paper is organised as follows: In section 2, I merge parts of the
Lagrangians from [8] and [9] to introduce the theory that will be used in the
rest of the paper. In section 3, I re-derive a simplified version of the
background solutions found in [7]. Introducing modified conformal time, these
solutions can be written in closed form. In section 4, I analyse metric
perturbations in a flat Friedmann universe in comoving gauge. I discuss
stability issues and calculate the primordial spectra of tensor and scalar
perturbations for the particular case of a radiation dominated background. The
analysis of the Mukhanov-Sasaki equations with modified dispersion relations
follows similar steps as [17]. In section 5, I summarize my results and give a
brief outlook on possible extensions. In Appendix A, I present the second
order actions for a more general mimetic theory and for the more general case
of perturbations around a non-flat Friedmann universe. In Appendix B, I
perform the linear stability analysis for bouncing solutions driven by higher
order spatial curvature terms in a non-flat universe, as found in [8].
Throughout this paper I use Planck units where $G_{0}=G(\Box\phi=0)=1$,
$\hbar=1$, $c=1$, $k_{B}=1$.

##  2 The theory

Consider the shift-symmetric theory of mimetic gravity defined by

$\mathcal{S}_{g}=-\frac{1}{16\pi}\int\textup{d}^{4}x\sqrt{-g}\left\\{\mathcal{L}_{\textup{nhd}}+\mathcal{L}_{\textup{hd}}+\lambda\left(g^{\mu\nu}\phi_{,\mu}\phi_{,\nu}-1\right)\right\\},$
(3)

where the Lagrangian
$\mathcal{L}=\mathcal{L}_{\textup{nhd}}+\mathcal{L}_{\textup{hd}}$ is divided
into a part $\mathcal{L}_{\textup{nhd}}$ without higher derivatives and a part
$\mathcal{L}_{\textup{hd}}$ with higher spatial derivatives in the
corresponding second order actions.

In [8] we found that the Lagrangian

$\mathcal{L}_{\textup{nhd}}:=f(\Box\phi)R+(f(\Box\phi)-1)\widetilde{R}+2\Lambda(\Box\phi)$
(4)

with

$\widetilde{R}:=2\phi^{,\mu}\phi^{,\nu}G_{\mu\nu}-(\Box\phi)^{2}+\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi,$

leads to a modified Einstein equation that is free of all higher derivatives
of the metric. Using the mimetic constraint, a slight generalization of a
calculation from [9] shows that (apart from considerations involving boundary
terms) it is equivalent to consider the Lagrangian

$\displaystyle\mathcal{L}_{\textup{nhd}}$
$\displaystyle\doteq\tfrac{4}{3}\ell(\Box\phi)-f(\Box\phi)\left(\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi-\tfrac{1}{3}(\Box\phi)^{2}\right)-\widetilde{R},$
(5)

where $\doteq$ denotes equality up to a total covariant derivative and

$\ell(\Box\phi)=\tfrac{3}{2}\Box\phi
F(\Box\phi)-(\Box\phi)^{2}f(\Box\phi)+\tfrac{3}{2}\Lambda(\Box\phi).$ (6)

The function $F(\Box\phi)$ is definded by
$F^{\prime}(\Box\phi)\equiv\mathrm{d}F/\mathrm{d}\Box\phi=f(\Box\phi)$. In the
following, the choice

$\Lambda(\Box\phi)=\tfrac{2}{3}\,\Box\phi\left(\Box\phi
f(\Box\phi)\,-F(\Box\phi)\right)$ (7)

will prove to be particularly simple.

Motivated by the goal of renormalizability along the lines of Hořava gravity
[10], as suggested in [9], I include the sixth order higher spatial derivative
terms

$\mathcal{L}_{\textup{hd}}:=\sigma_{T}^{4}\,\widetilde{C}^{\mu}_{\nu}\widetilde{C}^{\nu}_{\mu}-\frac{\sigma_{S}^{4}}{8}P^{\mu}_{\nu}\nabla_{\mu}\widetilde{R}\nabla^{\nu}\widetilde{R},$
(8)

where the projector $P_{\mu}^{\nu}=\delta_{\mu}^{\nu}-\phi_{,\mu}\phi^{,\nu}$
and the covariant analogues of the spatial Cotton tensor
$\widetilde{C}^{\mu}_{\nu}$ and the spatial Ricci tensor
$\widetilde{R}_{\mu\nu}$ have been introduced in [9] as

$\displaystyle\widetilde{C}_{\nu}^{\mu}$
$\displaystyle:=-\tfrac{1}{\sqrt{-g}}\epsilon^{\mu\rho\kappa\lambda}\nabla_{\lambda}\phi\,\nabla_{\rho}\left(\widetilde{R}_{\nu\kappa}-\frac{1}{4}g_{\nu\kappa}\widetilde{R}\right),$
(9) $\displaystyle\widetilde{R}_{\mu\nu}$
$\displaystyle:=P_{\mu}^{\alpha}P_{\nu}^{\beta}R_{\alpha\beta}+\nabla_{\alpha}\left(\phi^{,\alpha}\nabla_{\mu}\nabla_{\nu}\phi\right).$
(10)

Note that $g^{\mu\nu}\widetilde{R}_{\mu\nu}=-\widetilde{R}$, using the mimetic
constraint.

#### NB.

In a flat Friedmann universe, the Lagrangian studied in [7],

$\displaystyle\widetilde{\mathcal{L}}$
$\displaystyle=f(\Box\phi)R+2\Lambda(\Box\phi)$ (11)
$\displaystyle\doteq\tfrac{4}{3}\ell(\Box\phi)-f(\Box\phi)\left(\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi-\tfrac{1}{3}(\Box\phi)^{2}+\widetilde{R}\right),$

leads to the same background dynamics as (3). However, in general its modified
Einstein equation contains higher mixed derivatives of the metric in the
synchronous frame $t=\phi$. While it was found in [18] that these mixed
derivatives can change the sign of the gradient term of scalar metric
perturbations and help to prevent a gradient instability, it was later
realized in [19], [20] that this comes at the price of introducing an
additional hidden degree of freedom. Even tough this second scalar degree of
freedom does not show up when perturbing around a homogeneous background in
unitary gauge, it was found in [20] that already for perturbations around
Minkowski spacetime with non-homogeneous mimetic field profile it can lead to
instabilities. Since this type of higher mixed derivatives is not present in
$\mathcal{L}_{\textup{nhd}}$, one could hope that no such additional degree of
freedom will appear in this theory. A full Hamiltonian analysis similar to
[19], [20], [21], [22] would of course be a more involved task beyond the
scope of this paper.

##  3 Background dynamics

In the homogeneous, isotropic background given by the flat Friedmann metric

$\mathrm{d}s^{2}=dt^{2}-a^{2}(t)\delta_{ij}\mathrm{d}x^{i}\mathrm{d}x^{j},$
(12)

the only consistent background solution of the mimetic constraint up to shifts
is

$\phi=t,$ (13)

for which we find

$\Box\phi\equiv\kappa=3H\equiv 3\tfrac{\dot{a}}{a}.$ (14)

Either from the equation of motion given in [8] or from analysis of the zeroth
order action (see appendix A), we arrive at the modified Friedmann equation

$\tfrac{2}{3}\left(\kappa\ell^{\prime}(\kappa)-\ell(\kappa)\right)=\frac{c_{\textup{MDM}}}{a^{3}}+8\pi\varepsilon^{m}=:8\pi\varepsilon.$
(15)

The constant of integration $c_{\textup{MDM}}$ describes the contribution of
mimetic matter and $\varepsilon^{m}$ is some general homogeneous, isotropic
matter energy density. Note that this background equation is the same for all
the theories (3), (5) and (11). Using the simplifying choice (7) and the
suggestive notation $G=1/f$ familiar from [7], the modified Friedmann equation
becomes

$H^{2}=\frac{8\pi}{3}G(\kappa)\,\varepsilon.$ (16)

Assuming a monotonically decreasing dependence of $\varepsilon$ on the scale
factor, such a modified Friedmann equation can be understood as an integral
curve of the form $a(H^{2})$ in the phase space spanned by $a$ and $H$, cf.
[8], [23]. The only possible relations of this form which a) are one-to-one b)
have limiting curvature $\kappa<\kappa_{0}$ and c) obey the GR-limit
$G(\kappa)=1+\mathcal{O}\left((\kappa/\kappa_{0})^{2}\right)$, generically
replace the Big Bang singularity by a smooth transition to an initial de
Sitter stage. In (16) such a behaviour can only be realized by “asymptotic
freedom”111Note that in (15) without using (7), such a background behaviour
could be equally well implemented by a choice of $\Lambda$ without asymptotic
freedom. However, as shown in [7], [8], asymptotic freedom becomes unavoidable
for singularity resolution in an anisotropic universe., i.e.
$G(\kappa\to\kappa_{0})\to 0$. For a concrete example with limiting curvature
$\kappa_{0}$, I will take the simple choice

$G(\kappa)=f(\kappa)^{-1}=1-\left(\frac{\kappa}{\kappa_{0}}\right)^{2},$ (17)

familiar from [7]. Note that in this case $\ell^{\prime\prime}=f^{2}$.
Assuming a single matter component with constant equation of state
$w=p/\varepsilon$, taking a time derivative of (16) we find that during the
inflationary stage

$-\frac{\dot{H}}{H^{2}}=\frac{3(1+w)}{2}G(\kappa)\ll 1,$ (18)

where $\dot{}$ denotes $t$ derivatives. However, the second “slow-roll”
parameter

$\frac{\ddot{H}}{2H\dot{H}}=\frac{3(1+w)}{2}=\mathcal{O}(1)$ (19)

is constant and of order of unity, showing that this background solution does
not fit into a “slow-roll” description.

It is straightforward to obtain the following implicit solution of (16) for
$\kappa(t)$:

$\frac{1+w}{2}\kappa_{0}\,t=\frac{\kappa_{0}}{\kappa}-\mathrm{atanh}\frac{\kappa}{\kappa_{0}}$
(20)

#### Modified conformal time.

Introducing the modified conformal time coordinate $\tilde{\eta}$ by

$\mathrm{d}\tilde{\eta}=\frac{\mathrm{d}t}{a\sqrt{f}},$ (21)

the modified Friedmann equation (16) in modified conformal time $\tilde{\eta}$
looks exactly like the usual Friedmann equation in usual conformal time:

$\left(a_{\tilde{\eta}}\right)^{2}=\frac{8\pi}{3}\varepsilon\,a^{4}$ (22)

A subscript $a_{\tilde{\eta}}=\partial_{\tilde{\eta}}a$ denotes $\tilde{\eta}$
derivatives. Assume that the total energy density is dominated by a component
with equation of state $w$ and parametrized as

$\frac{8\pi}{3}\varepsilon=\left(\frac{c}{H_{l}^{2}}\right)^{\tfrac{1+3w}{2}}\left(\frac{1+3w}{2}\,a\right)^{-3(1+w)},$
(23)

where the prefactors were introduced for later convenience and
$H_{l}=\kappa_{0}/3$ accounts for dimensions such that $c$ is dimensionless.
In this case we find the solution

$a(\tilde{\eta})=\frac{2}{(1+3w)}\frac{\sqrt{c}}{H_{l}}\,\tilde{\eta}^{\frac{2}{1+3w}},$
(24)

where the initial condition $a(\tilde{\eta}=0)=0$ was used. The range of
modified conformal time is then $0<\tilde{\eta}<\infty$.

For the choice (17) the solution for $\kappa(\tilde{\eta})$ is

$\left(\frac{\kappa_{0}}{\kappa}\right)^{2}=\left(\frac{H_{l}}{H}\right)^{2}=1+c\,\tilde{\eta}^{\frac{6(1+w)}{1+3w}}.$
(25)

It describes a smooth transition from exponential expansion at limiting
curvature with $H\sim H_{l}$ to the late time stage dominated by the matter
component with equation of state $w$ where $H\propto t^{-1}$. The end of the
inflationary stage, i.e. the end of accelerated expansion, happens at

$\left(\frac{H_{f}}{H_{l}}\right)^{2}=\frac{1+3w}{3(1+w)},$ (26)

or, in modified conformal time at

$\tilde{\eta}_{f}=\left(\frac{2}{(1+3w)c}\right)^{\frac{1+3w}{6(1+w)}}.$ (27)

##  4 Metric perturbations in comoving gauge

In the same way as the modified Einstein equation of mimetic gravity takes its
simplest form with the choice of time coordinate $t=\phi$, the analysis of
metric perturbations is most readily performed in comoving gauge
$\delta\phi=0$. Metric perturbations in different versions of mimetic gravity
and also in different gauges have been analysed already a number of times,
e.g. in [18, 20, 21, 24, 25, 26, 27, 28, 16]. See appendix A for the
calculation of second order actions around a non-flat Friedmann background in
a more general theory.

Starting with a general ADM metric,

$\textup{d}s^{2}=N^{2}\mathrm{d}t^{2}-\gamma_{ij}\left(\mathrm{d}x^{i}+N^{i}\mathrm{d}t\right)\left(\mathrm{d}x^{j}+N^{j}\mathrm{d}t\right),$
(28)

the mimetic constraint fixes the lapse function as

$N^{2}=\frac{\left(\partial_{0}\phi-N^{i}\partial_{i}\phi\right)^{2}}{\left(1+\gamma^{ij}\partial_{i}\phi\partial_{j}\phi\right)}.$
(29)

In comoving gauge where $\partial_{i}\phi=0$ this implies $N=\dot{\phi}$.
Thus, the spatial slices coincide with slices of constant $\phi$ and the
homogeneous lapse $N$ is determined only by the background. Using the
background solution (13) amounts to setting $N=1$. In this gauge it holds that
the quantities

$\Box\phi=\kappa,\qquad\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi=\kappa^{ij}\kappa_{ij},\qquad\widetilde{R}={{}^{3}\\!R}$
(30)

are given in terms of the extrinsic curvature
$\kappa_{ij}=\frac{1}{2}\left(\dot{\gamma}_{ij}-N_{i,j}-N_{j,i}\right)$ and
intrinsic Ricci curvature ${}^{3}\\!R_{ij}$ of spatial slices and their
traces. Similar, straightforward calculations show that the higher order terms
from (8) take the form

$\displaystyle
P^{\mu}_{\nu}\nabla_{\mu}\widetilde{R}\nabla^{\nu}\widetilde{R}$
$\displaystyle=-\gamma^{kl}\partial_{k}({{}^{3}\\!R})\,\partial_{l}({{}^{3}\\!R}),$
(31) $\displaystyle\widetilde{C}^{\mu}_{\nu}\widetilde{C}^{\nu}_{\mu}$
$\displaystyle={{}^{3}C}^{i}_{j}\,{{}^{3}C}^{j}_{i}.$ (32)

The spatial Cotton tensor is defined as in [9] by

${{}^{3}C}_{j}^{i}=\frac{1}{\sqrt{\gamma}}\epsilon^{ikl}\bar{\nabla}_{k}\left({{}^{3}\\!R}_{jl}-\frac{1}{4}\gamma_{jl}\,{{}^{3}\\!R}\right),$
(33)

where $\bar{\nabla}$ denotes the covariant derivative associated to
$\gamma_{ij}$ and indices on spatial tensors are raised with $\gamma^{ij}$.

Perturbing around the flat Friedmann background (12) in comoving gauge
$\phi=t$, the metric perturbations can be further decomposed as

$\gamma_{ij}=a^{2}(t)\left(e^{-2\psi}{\delta_{ij}}-2E_{,ij}-2F_{(i,j)}+(e^{h})_{ij}\right),\qquad\qquad
N_{i}=\chi_{,i}-aS_{i}$ (34)

where

$F^{i}_{,i}=0,\qquad S^{i}_{,i}=0,\qquad h^{i}_{i}=0,\qquad h^{i}_{j,i}=0.$
(35)

Indices are raised with $\delta^{ij}$ and comma denotes partial derivative.
While the temporal coordinate is completely fixed, there is still a remaining
freedom of choice in the spatial coordinates. In contrast to the synchronous
gauge condition $N_{i}=0$, the choice $E=0$, $F_{i}=0$ fixes the coordinates
uniquely, cf. [29].

###  4.1 Vector perturbations

In mimetic gravity, like in GR, the vector perturbations parametrized by
$F_{i}$ and $S_{i}$ as

$\gamma_{ij}=a^{2}(t)\left({\delta_{ij}}-F_{i,j}-F_{j,i}\right),\qquad\qquad
N_{i}=-aS_{i},$ (36)

are non-dynamical in the absence of sources. Note that the spatial Ricci
tensor ${}^{3}\\!R_{ij}=\mathcal{O}(F^{2})$ is second order in vector
perturbations when perturbing around a spatially flat Friedmann universe.
Hence, higher order spatial curvature terms like in (8) do not contribute to
the linearized equation of motion for vector perturbations. Choosing the gauge
$S_{i}=0$, we can use the equations of motion derived in synchronous
coordinates in [8] and express them in terms of the gauge invariant variable
$V_{i}=S_{i}-a\dot{F}_{i}$. Assuming hydrodynamical matter of the perfect
fluid type, the only non-vanishing vector components of the perturbed energy
momentum tensor $\delta T_{\mu\nu}$ are of the form $\delta
T_{0i}=\left(\varepsilon^{m}+p\right)\delta u_{\perp i}$, [29].

The spatial modified Einstein equations become

$\partial_{t}\left(a^{2}f\left(V_{i,j}+V_{j,i}\right)\right)=0,$ (37)

and have the solution

$V_{i}=\frac{C_{\perp i}}{fa^{2}}.$ (38)

Note that the gauge invariant vector perturbation $V_{i}$ decays both in the
late time limit, like in GR, as well as in the early time limit, as
$f\to\infty$. This is completely in line with the previous results from [8]
that due to asymptotic freedom anisotropies decay during contraction at
limiting curvature. Provided that the constant of integration $C_{\perp i}$
and thus the maximum value of $V_{i}$ is bounded, we see that in this model
vector perturbations never become important. Reinserting this solution into
the $0-i$ modified Einstein equation

$f\Delta V_{i}=16\pi\left(\varepsilon^{m}+p\right)a\,\delta u_{\perp i},$ (39)

with $\Delta=\delta^{kl}\partial_{k}\partial_{l}$, shows that for the physical
velocity $\delta v^{i}=-a^{-1}\delta u_{\perp i}$ it holds that

$\delta v^{i}=-\frac{\Delta C_{\perp i}}{16\pi
a^{4}\left(\varepsilon^{m}+p\right)},$ (40)

where $f$ cancels and which is thus the same equation that one finds in
standard GR. Note that only the matter energy density and pressure are
appearing since mimetic matter does not source vector perturbations.

###  4.2 Tensor perturbations

Tensor perturbations are parametrized by

$\gamma_{ij}=a^{2}\left(\eta_{ij}+h_{ij}+\frac{1}{2}h_{i}^{k}h_{kj}\right),\qquad
N_{i}=0,$ (41)

where indices on $h_{ij}$ are raised with $\delta^{ij}$ and $h^{i}_{i}=0$,
$h^{i}_{j,i}=0$. Expanding (66) to second order in $h_{ij}$ yields the second
order action (see appendix A)

${{}^{(2)}\mathcal{S}}_{h}=\frac{1}{64\pi}\int\\!\mathrm{d}t\,\mathrm{d}^{3}x\,a^{3}\bigg{\\{}f\,(\dot{h}_{ij})^{2}-\frac{1}{a^{2}}(\partial
h_{ij})^{2}-\frac{\sigma_{T}^{4}}{a^{6}}(\partial^{3}h_{ij})^{2}\bigg{\\}},$
(42)

where $\dot{}$ denotes $t$ derivatives, $(\partial
h_{ij})^{2}\equiv\delta^{kl}\partial_{k}h^{ij}\partial_{l}h_{ij}$ and
$(\partial^{3}h_{ij})^{2}\equiv(\partial\Delta h_{ij})^{2}$. Since $f>0$,
tensor perturbations do not exhibit any instabilities provided that
$\sigma_{T}^{4}\geq 0$.

We can read off the propagation speed of gravitational waves,

$c_{T}^{2}=f^{-1}(\kappa)=G(\kappa),$ (43)

and find that it deviates from unity in the early time/high curvature regime.
In particular, note that $c_{T}$ is vanishing at limiting curvature
$\kappa=\kappa_{0}$ in the case of asymptotic freedom. Using (17), right at
the end of inflation $c_{T}^{2}=\tfrac{2}{3(1+w)}$. Already one second later,
in the late time regime where $\kappa\propto t^{-1}$ and $t$ roughly
corresponds to the time since the end of inflation, the propagation speed of
gravitational waves is approximately

$1-c_{T}\approx
1-\sqrt{1-\left(\frac{1}{\kappa_{0}\,t}\right)^{2}}\lesssim\left(\frac{\kappa_{pl}}{\kappa_{0}}\right)^{2}\left(\frac{1\,\textup{sec}}{t}\right)^{2}\times
10^{-86}.$ (44)

The limiting curvature $\kappa_{0}$ could naturally be taken to lie a few
orders of magnitude below the Planck curvature $\kappa_{pl}=1/t_{pl}$ for the
sake of singularity resolution. Thus, late time experimental constraints like
$1-c_{T}<\mathcal{O}(10^{-15})$, which stems from the multi messenger event
GW170817 ([30], [31]) that happened around $t\sim 10^{9}$ years after the
conjectured inflationary period, are not touched in the slightest by this
deviation of $c_{T}$ from the speed of light in the limiting curvature regime.
They do, however, tightly constrain low curvature / IR modifications due to
additional terms in the Lagrangian like $(\Box\phi)^{2}$ or in particular
$\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi$, see [32], cf.
[33].

#### Mode expansion.

Introducing $y=af^{1/4}$ and substituting the expansion

$h^{i}_{j}(\tilde{\eta},\textbf{x})=\sqrt{32\pi}\int\\!\frac{\textup{d}^{3}k}{(2\pi)^{3/2}}\frac{u_{k}(\tilde{\eta})\,e^{i\textbf{k}\textbf{x}}}{y\sqrt{\mathrm{e}^{m}_{n}\mathrm{e}^{n}_{m}}}\,\mathrm{e}^{i}_{j}(\textbf{k}),$
(45)

where the mode function $u_{k}$ and the contraction
$\mathrm{e}^{m}_{n}\mathrm{e}^{n}_{m}$ of the polarization tensor only depend
on the magnitude of k by isotropy, the second order action in modified
conformal time $\mathrm{d}\tilde{\eta}=\frac{\mathrm{d}t}{a\sqrt{f}}$ becomes

${{}^{(2)}\mathcal{S}}_{h}=\frac{1}{2}\int\\!\textup{d}\tilde{\eta}\,\textup{d}^{3}k\left\\{u_{k,\tilde{\eta}}^{2}-\left(k^{2}+\frac{\sigma_{T}^{4}}{a^{4}}k^{6}-\frac{y_{\tilde{\eta}\tilde{\eta}}}{y}\right)u_{k}^{2}\right\\}.$
(46)

A subscript $u_{k,\tilde{\eta}}\equiv\partial u_{k}/\partial\tilde{\eta}$
denotes derivatives. The modified Mukhanov-Sasaki equation for $u_{k}$ reads

$u_{k,\tilde{\eta}\tilde{\eta}}+\left(\frac{\sigma_{T}^{4}}{a^{4}}k^{6}+k^{2}-\frac{y_{\tilde{\eta}\tilde{\eta}}}{y}\right)u_{k}=0,$
(47)

and the power spectrum of gravitational waves is given by [29]

$\delta^{2}_{h}(k,\tilde{\eta})=\frac{8}{\pi}\frac{k^{3}\left|u_{k}\right|^{2}}{y^{2}}.$
(48)

In the late time limit $f\to 1$, $\tilde{\eta}\to\eta$ is identical to regular
conformal time and (47) becomes just like in GR,

$u_{k,\eta\eta}+\left(k^{2}-\frac{a_{\eta\eta}}{a}\right)u_{k}=0.$ (49)

Initially, at $\tilde{\eta}\to 0$ the $k^{6}$ term is dominating. The
requirement of an initial minimal level of quantum fluctuations determines the
initial conditions at $\tilde{\eta}\to 0$ (up to an unimportant phase factor,
see [29]) as

$u_{k}=\frac{a}{\sigma_{T}k^{3/2}},\qquad
u_{k,\tilde{\eta}}=i\,\frac{\sigma_{T}k^{3/2}}{a}.$ (50)

In the initial region where the $k^{6}$ term is dominating and
$a_{\tilde{\eta}}\,a\ll k^{3}\sigma_{T}^{2}$, the solution of (47) with
initial conditions (50) is well described by the WKB approximation

$u_{k}=\frac{a(\tilde{\eta})}{\sigma_{T}k^{3/2}}\exp\left(-i\int\\!\mathrm{d}\tilde{\eta}\frac{\sigma_{T}^{2}k^{3}}{a(\tilde{\eta})^{2}}\right).$
(51)

Note that the initial spectrum is scale invariant, by virtue of the same
higher order $k^{6}$ term that is needed for a power counting renormalizable
theory. After “horizon” exit, i.e. when $y_{\tilde{\eta}\tilde{\eta}}/y$
becomes dominating, the solution of (47) is

$u_{k}(\tilde{\eta})=y(\tilde{\eta})\left(A^{T}_{k}+B^{T}_{k}\int^{\infty}_{\tilde{\eta}}\\!\frac{\mathrm{d}\tilde{\eta}^{\prime}}{y^{2}(\tilde{\eta}^{\prime})}\right),$
(52)

where the second term is decaying compared to the first term and will be
ignored for the following estimates. The primordial spectrum of tensor
perturbations is hence given by

$\delta^{2}_{h}(k,\tilde{\eta})\approx\frac{8}{\pi}k^{3}\left|A^{T}_{k}\right|^{2}.$
(53)

The initially scale invariant spectrum would only be preserved for modes which
exit the “horizon” before the $k^{2}$ term starts do dominate and in a region
where $y\propto a$, i.e. where the propagation speed $c_{T}$ is almost
constant. This is certainly not satisfied for modes exiting the horizon during
the inflationary stage of an asymptotically free model where $f$ is rapidly
changing. Hence we can expect a primordial spectrum of tensor perturbations
which is far from being scale invariant on large scales.

#### Radiation dominated background.

For a concrete example with the simplest possible background evolution,
consider the case where the total energy density is dominated by a component
with equation of state $w=1/3$. In this case, fixing $c$ by setting
$a(\tilde{\eta}_{f})=1$ at the end of inflation, the solution (24) is
$a(\tilde{\eta})=H_{l}\tilde{\eta}$ and we find

$y=\sqrt[4]{1+(H_{l}\tilde{\eta})^{4}},\qquad\frac{y_{\tilde{\eta}\tilde{\eta}}}{y}=\frac{3H_{l}^{4}\tilde{\eta}^{2}}{(1+(H_{l}\tilde{\eta})^{4})^{2}}.$
(54)

Before the end of inflation at $\tilde{\eta}_{f}=1/H_{l}$, the modified
Mukhanov-Sasaki equation (47) is well approximated by

$u_{k,\tilde{\eta}\tilde{\eta}}+\left(\frac{\tilde{\sigma}_{T}^{4}}{H_{l}^{8}\tilde{\eta}^{4}}k^{6}+k^{2}-3H_{l}^{4}\tilde{\eta}^{2}\right)u_{k}=0,$
(55)

with the dimensionless $\tilde{\sigma}_{T}:=\sigma_{T}H_{l}$. The “horizon”
exit of the mode $k$ happens at $\tilde{\eta}_{k}=k/H_{l}^{2}\sqrt{D}$ where
$D$ solves

$\tilde{\sigma}_{T}^{4}D^{3}+D=3.$ (56)

Note that by the assumption made before, this is only valid for modes
$k^{2}\ll DH_{l}^{2}$ which exit the “horizon” before the end of inflation,
i.e. $\tilde{\eta}_{k}\ll\tilde{\eta}_{f}$.

In the initial region where the $k^{6}$ term is dominating, the solution of
(55) is given by

$u_{k}(\tilde{\eta})=\frac{H_{l}^{2}\,\tilde{\eta}}{k^{3/2}\tilde{\sigma}_{T}}\exp\left(-i\frac{\tilde{\sigma}_{T}^{2}k^{3}}{H_{l}^{4}\tilde{\eta}}\right),$
(57)

where the initial conditions (50) have been taken into account.

It depends only on $\tilde{\sigma}_{T}$ if modes will exit the “horizon”
before or after the gradient term is dominating over the $k^{6}$ term. For
$\tilde{\sigma}_{T}\gg 1$ the gradient term will never become important and
$D=(3/\tilde{\sigma}_{T}^{4})^{1/3}$. By matching the absolute value of
$\left|{u_{k}}\right|$ from (57) and (52) we can estimate

$\left|A^{T}_{k}\right|\approx\frac{H_{l}}{k^{3/2}\tilde{\sigma}_{T}}\frac{a(\tilde{\eta}_{k})}{y(\tilde{\eta}_{k})}\approx\frac{1}{\sqrt{D}\,\tilde{\sigma}_{T}\,k^{1/2}}.$
(58)

On the other hand, if $\tilde{\sigma}_{T}\ll 1$ then there is an intermediate
region where the $k^{2}$ term is dominating and $D=3$. In this region,
starting at

$\tilde{\eta}_{\ast}=\frac{\tilde{\sigma}_{T}}{H_{l}^{2}}\,k<\tilde{\eta}_{k}=\frac{k}{\sqrt{3}H_{l}^{2}},$
(59)

the absolute value $\left|u_{k}\right|$ is approximately constant and we can
estimate

$\left|A^{T}_{k}\right|\approx\frac{H_{l}}{k^{3/2}\tilde{\sigma}_{T}}\frac{a(\tilde{\eta}_{\ast})}{y(\tilde{\eta}_{k})}\approx\frac{1}{k^{1/2}}.$
(60)

Note that this is identical to the result that one would get if
$\tilde{\sigma}_{T}=0$ and the initial conditions were determined from the
$k^{2}$ term.

In summary, the primordial spectrum of large wavelength modes $k^{2}\ll
DH_{l}^{2}$ becomes

$\delta^{2}_{h}(k,\tilde{\eta})\approx\begin{cases}\frac{8}{\pi}\,k^{2}&\text{
if }\tilde{\sigma}_{T}\ll 1\\\
\frac{8}{\pi}\,\frac{k^{2}}{(\sqrt{3}\,\tilde{\sigma}_{T})^{2/3}}&\text{ if
}\tilde{\sigma}_{T}\gg 1\end{cases}$ (61)

Note that here $1/k$ is the physical wavelength at the end of inflation. In
both cases the spectrum is indeed far from being scale invariant, with a large
blue tilt and a spectral index of

$n_{T}=2.$ (62)

For a different background, e.g. $w=0$ the situation becomes more complicated
because $y_{\tilde{\eta}\tilde{\eta}}/y$ initially diverges as
$\sim-1/\tilde{\eta}^{2}$. Hence, after an initial sub-“horizon” $k^{6}$
domination region, there can be an earlier intermediate super-“horizon”
region, followed again by a sub-“horizon” region where $k^{2}$ dominates
before finally exiting the horizon. A similar situation occurs for scalar
perturbations in a dust dominated background, see below.

###  4.3 Scalar perturbations

The scalar metric perturbations in comoving gauge are characterized through
$\zeta$ and $\chi$ by

$\gamma_{ij}=a^{2}(t)e^{2\zeta}\delta_{ij},\qquad\qquad N_{i}=\chi_{,i},$ (63)

where

$\zeta=-\psi-\frac{H}{\dot{\phi}}\delta\phi.$ (64)

is the gauge invariant curvature perturbation in comoving gauge
$\delta\phi=0$.

In contrast to standard GR, even in the absence of any matter fluctuations the
conformal degree of freedom of mimetic gravity can become dynamical. In this
case the action (66) expanded to second order in scalar perturbations (see
appendix A) becomes

$\displaystyle{{}^{(2)}\mathcal{S}}=\frac{1}{8\pi}\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,a^{3}\bigg{\\{}$
$\displaystyle-3\ell^{\prime\prime}\dot{\zeta}^{2}+\tfrac{1}{a^{2}}(\partial\zeta)^{2}-\tfrac{\sigma_{S}^{4}}{a^{6}}\,(\partial^{3}\zeta)^{2}\phantom{\bigg{\\}}}$
(65)
$\displaystyle-\left[\tfrac{1}{3}(\ell^{\prime\prime}-f)\Delta\chi-2\ell^{\prime\prime}a^{2}\dot{\zeta}\right]\tfrac{1}{a^{4}}\Delta\chi\bigg{\\}},$

where $\dot{}$ denotes $t$ derivatives,
$(\partial\zeta)^{2}\equiv\delta^{kl}\partial_{k}\zeta\partial_{l}\zeta$ and
$(\partial^{3}\zeta)^{2}\equiv(\partial\Delta\zeta)^{2}$. In the case
$\ell^{\prime\prime}-f=0$, variation with respect to $\chi$ shows that
$\dot{\zeta}=0$ and thus the conformal degree of freedom is frozen. In the
case $\ell^{\prime\prime}-f\neq 0$, the second order action after integrating
out $\chi$ becomes

${{}^{(2)}\mathcal{S}_{\zeta}}=\frac{1}{8\pi}\int\\!\mathrm{d}t\,\mathrm{d}^{3}x\,a^{3}\bigg{\\{}\frac{3\ell^{\prime\prime}f}{\ell^{\prime\prime}-f}\,\dot{\zeta}^{2}+\frac{1}{a^{2}}\left(\partial\zeta\right)^{2}-\frac{\sigma_{S}^{4}}{a^{6}}\left(\partial^{3}\zeta\right)^{2}\bigg{\\}}.$
(66)

#### No ghost instability.

Since $f>0$, the condition to have no ghost instability reads

$\frac{\ell^{\prime\prime}}{\ell^{\prime\prime}-f}>0\qquad\Rightarrow\quad\ell^{\prime\prime}>f.$
(67)

The other possible case $\ell^{\prime\prime}<0$ was excluded because for a
smooth low curvature GR limit it must hold that $\ell^{\prime\prime}\to 1$,
$f\to 1$ as $\kappa\to 0$.

The condition (67) constrains the slope of the modified Friedmann equation
(15) by

$\frac{8\pi}{3}\frac{\partial\varepsilon}{\partial{H^{2}}}=\frac{1}{\kappa}\frac{\partial}{\partial{\kappa}}\left(\kappa\ell^{\prime}-\ell\right)=\ell^{\prime\prime}\overset{!}{>}f>0,$
(68)

which can be rewritten as

$0<\frac{\partial H^{2}}{\partial\varepsilon}<\frac{8\pi}{3}G(\kappa).$ (69)

Note that this condition is fully general and does not make use of the
simplifying choice of $\Lambda$ given by (7). Using (7), the condition becomes
$\kappa G^{\prime}<0$ and we see that the running gravitational constant can
only decrease when going to higher curvatures. This condition is always
satisfied for the asymptotically free background solutions from section 3.

#### Instability of spatially flat bouncing solutions.

In a previous mimetic model [34], limiting curvature was realized by a bounce.
In order to achieve a bouncing background solution in a spatially flat
universe, the modified relation $a(H^{2})$ cannot be one-to-one. In fact, the
generic background evolution described in section 3 can only be circumvented
by including multi-valued functions with intricate branch changes in the
Lagrangian, cf. [35]. Moreover, it was discussed in [23] that the bouncing
solution from [34] is unstable under perturbations. From the above it is easy
to see that this is an unavoidable feature of any bouncing solution of (15):
In order to obtain a bounce there must be a region where $H^{2}(\varepsilon)$
decreases until it eventually reaches a zero at some finite $a_{\min}>0$. In
this region the condition (69) is violated and hence there is a ghost
instability.

The stability analysis for bounces in non-flat universes found in [8], driven
by higher order spatial curvature terms, is performed in appendix B.

#### Gradient instability?

As has been noticed in [24], mimetic models without any higher derivatives
typically exhibit a gradient instability. In fact, also for models which
include higher derivatives, a negative square of the speed of sound
$c_{S}^{2}<0$ is something we have to deal with, at least in the low curvature
regime, in any mimetic model with a well behaved GR limit. However, in this
same limit also $\ell^{\prime\prime}\to f$ and thus the scalar degree of
freedom of mimetic gravity stops to propagate.

In a theory with higher mixed derivatives, the sign of the gradient term can
be made negative in the high curvature regime. In fact, for the Lagrangian
(11) the second order action is

$\displaystyle{{}^{(2)}\widetilde{\mathcal{S}}_{\zeta}}=\frac{1}{8\pi}\int\\!\mathrm{d}t\,\mathrm{d}^{3}x\,a^{3}\bigg{\\{}$
$\displaystyle\frac{3\ell^{\prime\prime}f}{\ell^{\prime\prime}-f}\,\dot{\zeta}^{2}+\left[f+\tfrac{1}{a}\partial_{0}\left(a\frac{6f}{\kappa}\right)\right]\tfrac{1}{a^{2}}(\partial\zeta)^{2}+\frac{6f^{\prime}}{\kappa}\tfrac{1}{a^{4}}\left(\partial^{2}\zeta\right)^{2}\bigg{\\}}.$
(70)

A similar second order action has been used in [18] and in [26]. However, the
models considered in these works were not symmetric under shifts of $\phi$.
Using a $\phi$ dependent potential essentially amounts to the introduction of
a time dependent background. In this way it is easy to produce any background
evolution one could wish for, including an inflationary stage, cf. [6]. In
such a highly flexible model, indeed the sign of the gradient term in (70) can
be adjusted more or less independent of the background. However, note that
such an adjustment must restrict to the high curvature regime. Otherwise the
GR limit is violated, as it happens in [26].

Sticking to shift-symmetric mimetic models, the background evolution and the
evolution of perturbations are no longer decoupled, but must be driven by the
same dynamics. Trying to use asymptotic freedom and the background solution
from section 3 in the model (70), we find that the gradient term would
actually become negative during the inflationary stage, however the sign of
the $(\partial^{2}\zeta)^{2}$ term is strictly positive and diverges as
$f^{\prime}\to\infty$. This shows that even though higher mixed derivatives
can make $c_{S}^{2}$ positive at high curvatures, in shift symmetric models
they actually lead to a higher order instability that could be removed only by
a high amount of tuning.

Instead, let us return to our original theory (3) without higher mixed
derivatives but with higher spatial derivatives coming from (8). In the
corresponding second order action (66) the gradient term comes with the
constant prefactor $+1$ and in the case $\sigma_{S}=0$ we would read off the
speed of sound

$c_{S}^{2}=-\frac{\ell^{\prime\prime}-f}{3\ell^{\prime\prime}f}<0.$ (71)

We will continue to use the name $c_{S}^{2}$ to refer to this quantity also in
the case $\sigma_{S}\neq 0$, but it is important to note that in this case due
to the modified dispersion relation the true speed of sound is in general
different from $c_{S}^{2}$. Note that $c_{S}^{2}$ has to be negative
throughout by (67). Using (17), it can be rewritten as

$c_{S}^{2}=-\frac{G(1-G)}{3}.$ (72)

Note that $c_{S}^{2}$ goes to zero in the late time limit as required by the
GR limit $G\to 1$, but it also goes to zero in the early time limit as $G\to
0$. The minimal value $c_{S}^{2}=-1/12$ is reached when $G=1/2$ at
$\kappa=\kappa_{0}/\sqrt{2}$. At the earliest times, the gradient term in (66)
will be dominated by the higher order spatial derivative term
$\tfrac{1}{a^{6}}(\partial^{3}\zeta)^{2}$. In this region there is no
instability, provided that $\sigma_{S}^{4}>0$. Thus, the potential gradient
instability region is sandwiched between the region of domination of the
higher order term and the late time region where the scalar degree of freedom
is “frozen”. These two other regions are without instabilities and hence the
gradient instability gets to act, if at all, only during a limited time. As I
will show below in two concrete examples, provided that $\sigma_{S}\gg
1/H_{l}$, the wrong sign gradient term cannot lead to any dangerous
instability.

#### Mode expansion.

Introducing the time coordinate $\tau$ and the expression $z$,

$\mathrm{d}t=a\sqrt{f}\,\mathrm{d}\tilde{\eta}=a\,\sqrt{\frac{3f\ell^{\prime\prime}}{\ell^{\prime\prime}-f}}\,\mathrm{d}\tau,\qquad
z=\frac{a}{\sqrt{4\pi}}\left(\frac{3f\ell^{\prime\prime}}{\ell^{\prime\prime}-f}\right)^{1/4},$
(73)

the second order action (66) written in terms of the canonically normalized
variable $v=z\zeta$ becomes

${{}^{(2)}\mathcal{S}_{\zeta}}=\frac{1}{2}\int\\!\mathrm{d}\tau\,\mathrm{d}^{3}x\left\\{v_{\tau}^{2}+(\partial
v)^{2}-\frac{\sigma_{S}^{4}}{a^{4}}(\partial^{3}v)^{2}+\frac{z_{\tau\tau}}{z}v^{2}\right\\}.$
(74)

A subscript $v_{\tau}\equiv\partial v/\partial\tau$ denotes derivatives.
Performing the mode expansion of $v$ into Fourier modes $v_{k}$, the modified
Mukhanov-Sasaki equation becomes

$v_{k,\tau\tau}+\left(\frac{\sigma_{S}^{4}}{a^{4}}k^{6}-k^{2}-\frac{z_{\tau\tau}}{z}\right)v_{k}=0.$
(75)

The power spectrum of the curvature perturbation $\zeta$ is given by [29]

$\delta^{2}_{\zeta}(k)=\frac{k^{3}\left|v_{k}\right|^{2}}{2\pi^{2}z^{2}}.$
(76)

In the initial region where the $k^{6}$ term is dominating and
$a_{\tilde{\eta}}\,a\ll k^{3}\sigma_{S}^{2}$, the solution of (75) with
quantum vacuum initial conditions is well described by the WKB approximation

$v_{k}=\frac{a(\tilde{\tau})}{\sigma_{S}k^{3/2}}\exp\left(-i\int\\!\mathrm{d}\tau\frac{\sigma_{S}^{2}k^{3}}{a(\tilde{\eta})^{2}}\right),$
(77)

and it has a scale invariant spectrum. After “horizon” exit, when
$z_{\tau\tau}/z$ is dominating, the solution of (75) is

$v_{k}(\tau)=z(\tau)\left(A^{S}_{k}+B^{S}_{k}\int^{\infty}_{\tau}\\!\frac{\mathrm{d}\tau^{\prime}}{z^{2}(\tau^{\prime})}\right),$
(78)

where the second term is decaying compared to the first term and will be
ignored in the following estimates. At late times the primordial spectrum
after horizon exit is

$\delta_{\zeta}^{2}\approx\frac{k^{3}}{2\pi^{2}}\left|A^{S}_{k}\right|^{2}.$
(79)

The wrong sign gradient term $-k^{2}$ will never get to dominate and cannot
cause instability, provided that at “horizon” exit the condition

$\frac{\sigma_{S}^{4}}{a^{4}}k^{6}\sim\frac{z_{\tau\tau}}{z}\gg k^{2}$ (80)

holds. In other words, the physical wavelength $\lambda_{\textup{phys}}=a/k$
at horizon exit should satisfy

$\lambda_{\textup{phys}}\bigg{|}_{\frac{\sigma_{S}^{4}}{a^{4}}k^{6}\sim\frac{z_{\tau\tau}}{z}}\ll\sigma_{S}.$
(81)

The initial scale invariant spectrum would be preserved only if at horizon
exit $z\propto a$, i.e.

$-\frac{1}{c_{S}^{2}}=\frac{3f\ell^{\prime\prime}}{\ell^{\prime\prime}-f}\approx
const.$ (82)

This condition will in general not be satisfied both at late times where
$\ell^{\prime\prime}-f\to 0$ and at early times where $f\to\infty$,
$\ell^{\prime\prime}\to\infty$. Thus, without a substantial amount of tuning
of the functions $f$ and $\Lambda$ in order for (82) to be satisfied at least
in some intermediate region where the relevant modes exit the “horizon”, one
can already expect a primordial spectrum that will be far from scale
invariant.

Note that all modes $\zeta_{k}$ will at some point exit the “horizon” and
never again re-enter, also after transitioning to the post-inflationary phase.
This is a manifestation of the fact that in the GR limit the scalar degree of
freedom of pure mimetic gravity is dust-like and non-propagating. It is clear
that aiming for a more realistic model, matter perturbations would have to be
taken into account.

#### Radiation dominated background.

In modified conformal time $\tilde{\eta}$ the radiation dominated background
solution from section 3 is given by $a(\tilde{\eta})=H_{l}\tilde{\eta}$, where
$c$ was fixed by setting $a(\tilde{\eta}_{f})=1$, and we can find the explicit
expressions

$z=\sqrt{\tfrac{\sqrt{3}}{4\pi}\left(1+(H_{l}\tilde{\eta})^{4}\right)},\qquad\frac{\mathrm{d}\tilde{\eta}}{\mathrm{d}\tau}=\sqrt{3(1+(H_{l}\tilde{\eta})^{4})},\qquad\frac{z_{\tau\tau}}{z}=18H_{l}^{4}\tilde{\eta}^{2}.$
(83)

The modified Mukhanov-Sasaki equation (75) becomes

$v_{k,\tau\tau}+\left(\frac{\tilde{\sigma}_{S}^{4}}{H_{l}^{8}\tilde{\eta}^{4}}k^{6}-k^{2}-18H_{l}^{4}\tilde{\eta}^{2}\right)v_{k}=0,$
(84)

where the dimensionless $\tilde{\sigma}_{S}=H_{l}\sigma_{S}$ was introduced.
Assuming that $\tilde{\sigma}_{S}\gg\mathcal{O}(1)$, the “horizon” crossing
happens at

$\tilde{\eta}_{k}=\left(\frac{\tilde{\sigma}_{S}^{4}}{18}\right)^{1/6}\frac{k}{H_{l}^{2}}$
(85)

and the condition (80) is satisfied for all modes.

Accelerated expansion ends at $\tilde{\eta}_{f}=1/H_{l}$ and before that we
can approximate $\tau=\tilde{\eta}/\sqrt{3}$ and write the $\tau$ derivatives
in (84) as $\tilde{\eta}$ derivatives. Matching the absolute value of the
initial solution with quantum vacuum initial conditions

$v_{k}(\tilde{\eta})=\frac{H_{l}^{2}\,\tilde{\eta}}{3^{1/4}\tilde{\sigma}_{S}k^{3/2}}\exp\left(-i\frac{\sqrt{3}\tilde{\sigma}_{S}^{2}k^{3}}{H_{l}^{4}\tilde{\eta}}\right),$
(86)

to the solution (78) after horizon exit, we can estimate the late time
spectrum of large wavelength modes $k^{2}\ll
H_{l}^{2}(18/\tilde{\sigma}_{S}^{4})^{1/3}$ which exit the horizon before the
end of inflation as

$\delta^{2}_{\zeta}(k)\approx\frac{(2/3)^{2/3}}{3\pi}\frac{k^{2}}{\tilde{\sigma}_{S}^{2/3}}.$
(87)

Note that here $1/k$ is the physical wavelength at the end of inflation. It is
far from being scale invariant with a large blue tilt and a spectral index of

$n_{S}-1=2.$ (88)

Combining with (61), the tensor to scalar ratio is given by

$r=\frac{\delta^{2}_{h}}{\delta^{2}_{\zeta}}\approx\begin{cases}24\left(\tfrac{3}{2}\tilde{\sigma}_{S}\right)^{2/3}&\text{
if }\tilde{\sigma}_{T}\ll 1\\\
24\left(\frac{\sqrt{3}}{2}\frac{\tilde{\sigma}_{S}}{\tilde{\sigma}_{T}}\right)^{2/3}&\text{
if }\tilde{\sigma}_{T}\gg 1\end{cases}$ (89)

It can be small only if $\sigma_{T}\gg\sigma_{S}\gg 1/H_{l}$.

#### Dust/MDM dominated background.

If we consider now the background solution dominated by dust or mimetic dark
matter with equation of state $w=0$, the situation gets complicated by an
additional early intermediate super-“horizon” region. Even though in this case
long wavelength modes go through a short gradient instability phase, we will
find that the growth of modes during this stage is completely negligible if
$\sigma_{S}\gg 1/H_{l}$.

Using the background solution in modified conformal time $\tilde{\eta}$ from
section 3 with $w=0$ and fixing $c$ by $a(\tilde{\eta}_{f})=1$ at the end of
inflation, the scale factor is given by
$a(\tilde{\eta})=H_{l}^{2}\tilde{\eta}^{2}/8$ and we calculate

$z=\sqrt{\tfrac{\sqrt{3}}{16\pi}H_{l}\tilde{\eta}\left(1+2a^{3}\right)},\qquad\frac{\mathrm{d}\tilde{\eta}}{\mathrm{d}\tau}=\sqrt{3(1+2a^{3})},\qquad\frac{z_{\tau\tau}}{z}=\frac{3(-1+154a^{3})}{4\tilde{\eta}^{2}}.$
(90)

The modified Mukhanov-Sasaki equation (75) now reads

$v_{k,\tau\tau}+\left(\left(\frac{\tilde{\sigma}_{S}}{H_{l}\,a}\right)^{4}k^{6}+\frac{3}{4\tilde{\eta}^{2}}-k^{2}-\frac{231}{16}\left(H_{l}\,a\right)^{2}\right)v_{k}=0,$
(91)

and we see that compared to (84) there is an additional term coming from
$z_{\tau\tau}/z$.

Assuming $\tilde{\sigma}_{S}\gg 1$, modes with

$k^{2}\gg\mathcal{O}(10^{-1})\frac{H_{l}^{2}}{\tilde{\sigma}_{S}^{4/3}}$ (92)

exit the horizon at

$\tilde{\eta}_{k}\sim\mathcal{O}(1)\frac{\tilde{\sigma}_{S}^{1/3}}{H_{l}^{3/2}}\sqrt{k},$
(93)

and both the gradient term and the $1/\tilde{\eta}^{2}$ term never become
important.

On the other hand, modes with wavelengths larger than (92) already exit the
horizon for a first time at

$\tilde{\eta}_{1,k}\sim\mathcal{O}(1)\frac{\tilde{\sigma}_{S}^{2/3}}{H_{l}^{2}}k.$
(94)

Since $z_{\tau\tau}/z$ is changing sign at
$\tilde{\eta}_{0}=(2^{8}/77)^{1/6}/H_{l}$, modes will at some point shortly
re-enter the horizon due to the gradient term before finally re-exiting again.
Expanding $z_{\tau\tau}/z$ around $\tilde{\eta}_{0}$, we can estimate the
duration $\Delta\tilde{\eta}_{k}$ of the gradient instability region as

$k\Delta\tilde{\eta}_{k}\sim\mathcal{O}(1)\frac{k^{3}}{H_{l}^{3}}\ll\frac{\mathcal{O}(10^{-2})}{\tilde{\sigma}_{S}^{2}}.$
(95)

During this short time span the mode function $v_{k}$ can only grow by a
factor $\sim\exp\left(k\Delta\tilde{\eta}_{k}\right)$ which is completely
negligible for $\tilde{\sigma}_{S}\gg\mathcal{O}(1)$.

Ignoring the effects of the gradient instability region, one would estimate
the primordial spectrum of the longest wavelength modes $k\ll
H_{l}/\tilde{\sigma}_{S}^{2/3}$ which exit the horizon well before the end of
inflation as

$\delta^{2}_{\zeta}\approx\mathcal{O}(1)\frac{k^{3}}{H_{l}}.$ (96)

Again, this is far being scale invariant with a spectral index $n_{S}-1=3$.
Since in this dust/MDM dominated case large wavelength modes go through
several intermediate regions lasting only for a short time, such an analysis
like in the radiation dominated case where the leading order solutions in
different regions were continuously matched at the crossing between regions
should be taken with caution. A full calculation of primordial spectra would
require a numerical study beyond the scope of this paper.

##  5 Conclusions

The initial idea of “Asymptotically Free Mimetic Gravity” was to find a
concrete modified theory of gravity with limiting curvature to address the
singularity problem of GR. It has been successful in achieving this goal in a
variety of settings, including both cosmological as well as black hole
spacetimes. It was found that the concept of “asymptotic freedom”, i.e. the
vanishing of the $\Box\phi$ dependent gravitational constant at limiting
curvature, becomes a crucial ingredient to resolve anisotropic singularities.
Along the way it was realized that Hořava-like higher order spatial curvature
terms can be added in a simple, covariant way to mimetic gravity with the goal
of renormalizability.

Combining both ideas, in this paper I considered “Asymptotically Free Mimetic
Hořava Gravity”. In this theory, an initial stage of exponential expansion
with graceful exit is a necessary feature of any non-singular modified flat
Friedmann universe. It is a natural question to ask whether it could provide a
full-blown inflationary scenario without inflaton. Since the existence of the
inflationary background solution is independent of the matter content, as it
is anyway suppressed by asymptotic freedom, there is no need to assume
vanishing matter density during inflation. Thus, the idea of inflation driven
by asymptotic freedom of gravity could open an interesting possibility to
avoid the necessity of a reheating stage for the sake of populating the
universe with matter after inflation.

In this work I analyzed stability of the inflationary background solutions
under metric perturbations, considering only the degrees of freedom of pure
mimetic gravity. It was found that a ghost instability is naturally avoided by
any model with asymptotic freedom. Although the gradient term of scalar
perturbations is of the wrong sign throughout, its short era of domination is
sandwiched between the domination of higher order spatial curvature terms and
the late time region where the scalar degree of freedom of mimetic gravity
remains frozen forever. Thus, the gradient instability can be circumvented
thanks to higher order spatial curvature terms. Under the condition that
$\sigma_{S}\gg 1/H_{l}$, I showed that the inflationary background solutions
of asymptotically free models are free of any dangerous instability.

After passing stability tests, we have to ask if the primordial spectra
produced by such an inflationary model can agree with CMB observations. In
this paper I showed that for the simplest one-component models the answer to
this question is negative. The initially scale invariant spectra of both
tensor and scalar perturbations can in general not be preserved until the
horizon exit. For the concrete case of a radiation dominated inflationary
background, the primordial spectra of the largest wavelength modes were found
to be far from scale invariant with a large blue tilt and $n_{T}=n_{S}-1=2$.
Moreover, a small tensor-to-scalar ratio $r<\mathcal{O}(10^{-1})$ requires
tuning of the higher order coefficients such that
$\sigma_{T}\gg\mathcal{O}(10^{4})\sigma_{S}$.

However, it is clear that this analysis is incomplete, as any more realistic
model would also have to include matter fluctuations. While the primordial
spectrum of the dust-like scalar degree of freedom of mimetic gravity might be
far from scale invariant, any other spectator field in the inflationary
background solution would still acquire a nearly scale invariant spectrum.
Depending on the details of the conversion process of perturbations between
matter degrees of freedom and the conformal degree of freedom of mimetic
gravity, one could thus speculate that in a “curvaton”-like extension of the
asymptotically free mimetic model, a nearly scale invariant primordial matter
power spectrum would be obtainable. It remains an interesting open question
whether it is possible to construct a viable inflationary scenario from shift-
symmetric mimetic gravity.

## Acknowledgements

I would like to thank my supervisor Slava Mukhanov for his support and
guidance throughout the years of my PhD.

The author’s work is supported by the Deutsche Forschungsgemeinschaft (DFG,
German Research Foundation) under Germany’s Excellence Strategy – EXC-2111 –
390814868.

##  Appendix A Calculation of second order actions

In this appendix I present the calculation of the second order actions used
above. For generality, I will consider the action

$\mathcal{S}_{g}=-\frac{1}{16\pi}\int\textup{d}^{4}x\sqrt{-g}\left\\{L\left(X_{i}\right)+\lambda\left(g^{\mu\nu}\phi_{,\mu}\phi_{,\nu}-1\right)\right\\},$
(97)

where $L$ is a general function of the quantities

$\displaystyle X_{1}$ $\displaystyle=\Box\phi$ $\displaystyle\cong$
$\displaystyle\kappa$ (98) $\displaystyle X_{2}$
$\displaystyle=\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi-\tfrac{1}{3}(\Box\phi)^{2}$
$\displaystyle\cong$ $\displaystyle\tilde{\kappa}^{ij}\tilde{\kappa}_{ij}$
$\displaystyle X_{3}$ $\displaystyle=\widetilde{R}$ $\displaystyle\cong$
$\displaystyle{{}^{3}\\!R}$ $\displaystyle X_{4}$
$\displaystyle=-\widetilde{R}^{\mu\nu}\left(\nabla_{\mu}\nabla_{\nu}\phi-\tfrac{1}{3}\Box\phi\,g_{\mu\nu}\right)$
$\displaystyle\cong$ $\displaystyle{{}^{3}\\!R}^{ij}\tilde{\kappa}_{ij}$
$\displaystyle X_{5}$
$\displaystyle=\widetilde{R}^{\mu\nu}\widetilde{R}_{\mu\nu}-\tfrac{1}{3}\widetilde{R}^{2}$
$\displaystyle\cong$
$\displaystyle{{}^{3}\\!\widetilde{R}^{ij}}\,{{}^{3}\\!\widetilde{R}_{ij}}$
$\displaystyle X_{6}$
$\displaystyle=-P^{\mu}_{\nu}\nabla_{\mu}\widetilde{R}\nabla^{\nu}\widetilde{R}$
$\displaystyle\cong$
$\displaystyle\gamma^{kl}\bar{\nabla}_{k}({{}^{3}\\!R})\bar{\nabla}_{l}({{}^{3}\\!R})$
$\displaystyle X_{7}$
$\displaystyle=-P^{\mu}_{\nu}P^{\alpha}_{\gamma}P^{\beta}_{\delta}\nabla_{\mu}\widetilde{R}^{\gamma\delta}\nabla^{\nu}\widetilde{R}_{\alpha\beta}+\tfrac{1}{3}X_{6}$
$\displaystyle\cong$
$\displaystyle\gamma^{kl}\bar{\nabla}_{k}({{}^{3}\\!\widetilde{R}^{ij}})\,\bar{\nabla}_{l}({{}^{3}\\!\widetilde{R}_{ij}})$

where

$\displaystyle P_{\mu}^{\nu}$
$\displaystyle:=\delta_{\mu}^{\nu}-\phi_{,\mu}\phi^{,\nu}$ (99)
$\displaystyle\widetilde{R}_{\mu\nu}$
$\displaystyle:=P_{\mu}^{\alpha}P_{\nu}^{\beta}R_{\alpha\beta}+\nabla_{\alpha}\left(\phi^{,\alpha}\nabla_{\mu}\nabla_{\nu}\phi\right)$
(100) $\displaystyle\widetilde{R}$
$\displaystyle:=2\phi^{,\mu}\phi^{,\nu}G_{\mu\nu}-(\Box\phi)^{2}+\nabla^{\mu}\nabla^{\nu}\phi\nabla_{\mu}\nabla_{\nu}\phi$
(101)

The right column of (98), denoted by $\cong$, shows the quantities $X_{i}$
evaluated in comoving gauge $\phi=t$ where $N=1$, but the shift $N_{i}$ is
still arbitrary. They are given in terms of the trace and trace-less part
$\kappa=\gamma^{ij}\kappa_{ij}$,
$\tilde{\kappa}_{ij}=\kappa_{ij}-\tfrac{1}{3}\kappa\gamma_{ij}$,
${{}^{3}\\!\widetilde{R}}=\gamma^{ij}\,{{}^{3}\\!R_{ij}}$,
${{}^{3}\\!\widetilde{R}_{ij}}={{}^{3}\\!R_{ij}}-\tfrac{1}{3}{{}^{3}\\!\widetilde{R}}\gamma_{ij}$
of the extrinsic curvature
$\kappa_{ij}=\frac{1}{2N}\left(\dot{\gamma}_{ij}-\bar{\nabla}_{j}N_{i}-\bar{\nabla}_{i}N_{j}\right)$
and of the intrinsic Ricci tensor ${{}^{3}\\!R_{ij}}$ of the spatial slices of
(28), respectively. $\bar{\nabla}_{i}$ denotes the covariant derivative with
respect to the spatial metric $\gamma_{ij}$. Indices on spatial tensors are
raised with $\gamma^{ij}$.

Note that terms depending on the square of the spatial Cotton tensor
$\widetilde{C}^{\mu}_{\nu}$ as introduced in (9) are also covered by the
general Lagrangian ansatz, since

$\widetilde{C}^{\mu}_{\nu}\widetilde{C}^{\nu}_{\mu}\cong{{}^{3}C}^{i}_{j}\,{{}^{3}C}^{j}_{i}\doteq
X_{7}+(\tfrac{1}{3}-\tfrac{3}{8})X_{6}+\tfrac{1}{2}X_{3}X_{5}+3\;{{}^{3\\!}\widetilde{R}}^{i}_{j}{{}^{3}\\!\widetilde{R}}^{j}_{k}{{}^{3}\\!\widetilde{R}}^{k}_{i},$
(102)

where $\doteq$ now denotes equality up to a total covariant spatial derivative
$\bar{\nabla}$. In an isotropic universe the trace-less part of
${}^{3}\\!R_{ij}$ is first order in perturbations, hence the last term in
(102) is always higher order and does not contribute to the second order
action.

Consider now perturbations around the general non-flat Friedmann background

$\mathrm{d}s^{2}=dt^{2}-a^{2}(t)\eta_{ij}\mathrm{d}x^{i}\mathrm{d}x^{j},\qquad\eta_{ij}=\frac{\delta_{ij}}{\left(1+\tfrac{\varkappa}{4}(x^{2}+y^{2}+z^{2})\right)^{2}}$
(103)

in comoving gauge $\phi=t$. The metric perturbations of (28) are then further
decomposed as

$\gamma_{ij}=a^{2}(t)\left(e^{-2\psi}{\eta_{ij}}-2D_{i}D_{j}E-2D_{(i}F_{j)}+(e^{h})_{ij}\right),\qquad
N_{i}=\chi_{,i}-aS_{i}$ (104)

where

$\eta^{ij}D_{j}F_{i}=0,\qquad\eta^{ij}h_{ij}=0,\qquad\eta^{jk}D_{k}h_{ij}=0,\qquad\eta^{ij}D_{j}S_{i}=0$
(105)

and $D_{i}$ denotes the covariant derivative associated with $\eta_{ij}$.

#### Background.

Variation of the zeroth order action
${{}^{(0)}\mathcal{S}}={{}^{(0)}\mathcal{S}_{g}+{{}^{(0)}\mathcal{S}_{m}}}$
with respect to $a^{3}$ yields

$\delta_{a^{3}}{{}^{(0)}\mathcal{S}}=-\frac{1}{16\pi}\int\\!\mathrm{d}^{4}x\,\sqrt{\eta}\left[L-\tfrac{1}{a^{3}}\partial_{0}\left(a^{3}L_{1}\right)-4L_{3}\tfrac{\varkappa}{a^{2}}-16\pi
p\right]\,\delta(a^{3})$ (106)

where subscripts $L_{i}$ denote derivatives of $L$ with respect to $X_{i}$
evaluated on the background. Here it was used that for homogeneous, isotropic
matter

$\delta\mathcal{S}_{m}=\tfrac{1}{2}\int\\!\mathrm{d}^{4}x\sqrt{-g}\,T_{\mu\nu}\delta
g^{\mu\nu}=\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,p\,\delta(a^{3}).$ (107)

The background equation of motion is hence

$L-\frac{1}{a^{3}}\partial_{0}\left(a^{3}L_{1}\right)-4L_{3}\frac{\varkappa}{a^{2}}=16\pi
p.$ (108)

Using the continuity equation
$\frac{1}{a^{3}}\partial_{0}\left(a^{3}\varepsilon\right)=-3\frac{\dot{a}}{a}p$,
the first integral of (108) becomes the modified Friedmann equation

$\frac{1}{2}\left(\kappa
L_{1}-L\right)=\frac{c_{\textup{MDM}}}{u}+8\pi\varepsilon=8\pi\varepsilon,$
(109)

where the constant of integration $c_{\textup{MDM}}$ describes the
contribution of mimetic matter.

#### Tensor perturbations.

Tensor perturbations are parametrized by

$\gamma_{ij}=a^{2}\left(\eta_{ij}+h_{ij}+\frac{1}{2}h_{i}^{k}h_{kj}\right),\qquad
N_{i}=0,$ (110)

where indices on $h_{ij}$ are raised with $\eta^{ij}$. With this
parametrization the inverse spatial metric is

$\gamma^{ij}=\frac{1}{a^{2}}\left(\eta^{ij}-h^{ij}+\frac{1}{2}h^{ik}h_{k}^{j}\right)+\mathcal{O}(h^{3})$
(111)

and it holds that $\sqrt{\gamma}=a^{3}\sqrt{\eta}$ and
$\kappa=3\frac{\dot{a}}{a}$ are still homogeneous up to $\mathcal{O}(h^{3})$.
The extrinsic curvature and intrinsic Ricci scalar up to second order are

$\displaystyle\kappa_{ij}$
$\displaystyle=\frac{\dot{a}}{a}\gamma_{ij}+\frac{1}{2}a^{2}\left(\dot{h}_{ij}+\frac{1}{2}\dot{h}_{ik}h^{k}_{~{}j}+\frac{1}{2}h_{ik}\dot{h}^{k}_{~{}j}\right),$
(112) $\displaystyle{{}^{3}\\!R}$
$\displaystyle\doteq\frac{1}{a^{2}}\left(6\varkappa-\frac{1}{4}D_{k}h^{ij}D^{k}h_{ij}-\frac{1}{2}\varkappa\,h^{ij}h_{ij}\right),$
(113)

where $\doteq$ now denotes equality up to a total covariant spatial derivative
$D$. The spatial Ricci tensor and its trace-less part up to first order are
given by

${{}^{3}\\!R}_{ij}=2\varkappa\eta_{ij}-\frac{1}{2}\Delta h_{ij}+3\varkappa
h_{ij},\qquad{{}^{3}\\!\widetilde{R}}_{ij}=-\frac{1}{2}\Delta h_{ij}+\varkappa
h_{ij},$ (114)

where $\Delta=\eta^{ij}D_{i}D_{j}$. The quantities (98) expanded to second
order are

$\displaystyle X_{1}$ $\displaystyle=3\tfrac{\dot{a}}{a}$ (115) $\displaystyle
X_{2}$ $\displaystyle=\tfrac{1}{4}\dot{h}^{ij}\dot{h}_{ij}$ $\displaystyle
X_{3}$
$\displaystyle\doteq\tfrac{1}{a^{2}}\left(6\varkappa-\tfrac{1}{4}D_{k}h^{ij}D^{k}h_{ij}-\tfrac{1}{2}\varkappa\,h^{ij}h_{ij}\right)$
$\displaystyle X_{4}$
$\displaystyle\doteq\tfrac{1}{8a^{2}}\partial_{0}\left(D_{k}h^{ij}D^{k}h_{ij}+2\varkappa
h^{ij}h_{ij}\right)$ $\displaystyle X_{5}$
$\displaystyle\doteq\tfrac{1}{a^{4}}\left(\tfrac{1}{4}\Delta h^{ij}\Delta
h_{ij}+\varkappa D_{k}h^{ij}D^{k}h_{ij}+\varkappa^{2}h^{ij}h_{ij}\right)$
$\displaystyle X_{6}$ $\displaystyle=0$ $\displaystyle X_{7}$
$\displaystyle\doteq\tfrac{1}{a^{6}}\left(\tfrac{1}{4}D_{k}\Delta
h^{ij}D^{k}\Delta h_{ij}+\varkappa\Delta h^{ij}\Delta
h_{ij}+\varkappa^{2}D_{k}h^{ij}D^{k}h_{ij}\right)$

and we find the second order action

$\displaystyle{{}^{(2)}\mathcal{S}}_{h}=-\frac{1}{64\pi}\int\\!\mathrm{d}^{4}x$
$\displaystyle\sqrt{\eta}\,a^{3}\bigg{\\{}L_{2}\,\dot{h}^{ij}\dot{h}_{ij}-\left(L_{3}+\tfrac{1}{2a}\partial_{0}\left(aL_{4}\right)-\tfrac{2\varkappa}{a^{2}}L_{5}\right)\tfrac{2\varkappa}{a^{2}}\,h^{ij}h_{ij}+$
(116)
$\displaystyle+\left(L_{3}+\tfrac{1}{2a}\partial_{0}\left(aL_{4}\right)-\tfrac{4\varkappa}{a^{2}}\left(L_{5}+\tfrac{\varkappa}{a^{2}}L_{7}\right)\right)\tfrac{1}{a^{2}}h^{ij}\Delta
h_{ij}+\phantom{\bigg{\\}}}$
$\displaystyle+\frac{L_{5}+\tfrac{4\varkappa}{a^{2}}L_{7}}{a^{4}}h^{ij}\Delta^{2}h_{ij}-\frac{L_{7}}{a^{6}}h^{ij}\Delta^{3}h_{ij}\bigg{\\}}.$

Subscripts $L_{i}$ denote derivatives of $L$ with respect to $X_{i}$ evaluated
on the background.

#### Scalar perturbations.

The scalar metric perturbations in comoving gauge are parametrized through
$\zeta$ and $\chi$ by

$\gamma_{ij}=a^{2}(t)e^{2\tilde{\zeta}}\delta_{ij},\qquad\qquad
N_{i}=\chi_{,i},$ (117)

where

$\tilde{\zeta}=\zeta-\ln\left(1+\tfrac{\varkappa}{4}(x^{2}+y^{2}+z^{2})\right).$
(118)

The metric determinant is given by
$\sqrt{\gamma}=a^{3}e^{3\tilde{\zeta}}=a^{3}\sqrt{\eta}\,e^{3\zeta}$ and the
trace of extrinsic curvature is

$\kappa=3\frac{\dot{a}}{a}+3\dot{\zeta}-\bar{\Delta}\chi,$ (119)

where the Laplacian $\bar{\Delta}=\gamma^{ij}\bar{\nabla}_{i}\bar{\nabla}_{j}$
of the covariant derivative $\bar{\nabla}$ with respect to $\gamma_{ij}$ is
given by

$\bar{\Delta}\chi=\tfrac{1}{a^{2}}\Delta\chi-\tfrac{2}{a^{2}}\zeta\Delta\chi+\tfrac{1}{a^{2}}\eta^{ij}\zeta_{,i}\chi_{,j}+\mathcal{O}\left(\chi\zeta^{2}\right).$
(120)

To first order, the traceless part of $\kappa_{ij}$ is given entirely through
$\chi$ as

$\tilde{\kappa}_{ij}=-\bar{\nabla}_{j}\bar{\nabla}_{i}\chi+\tfrac{1}{3}\bar{\Delta}\chi\gamma_{ij}=-D_{j}D_{i}\chi+\tfrac{1}{3}\Delta\chi\eta_{ij}+\mathcal{O}(\zeta\chi).$
(121)

The spatial metric $\gamma_{ij}$ is conformally flat which simplifies the
calculation of

$\displaystyle{{}^{3}\\!R_{ij}}$
$\displaystyle=-\left(\tilde{\zeta}_{,ij}-\tilde{\zeta}_{,i}\tilde{\zeta}_{,j}\right)-\left(\tilde{\zeta}^{,k}_{~{}k}+\tilde{\zeta}^{,k}\tilde{\zeta}_{,k}\right)\gamma_{ij}$
(122)
$\displaystyle=2\varkappa\eta_{ij}-\left(D_{j}D_{i}\zeta+\Delta\zeta\eta_{ij}\right)+\left(\zeta_{,i}\zeta_{,j}-\eta^{kl}\zeta_{,k}\zeta_{,l}\,\eta_{ij}\right)$
$\displaystyle{{}^{3}\\!R}$
$\displaystyle=\frac{e^{-2\zeta}}{a^{2}}\left(6\varkappa-4\Delta\zeta-2\eta^{ij}\zeta_{,i}\zeta_{,j}\right)$
(123)
$\displaystyle=\frac{1}{a^{2}}\left(6\varkappa-4\Delta\zeta-12\varkappa\zeta+12\varkappa\zeta^{2}+8\zeta\Delta\zeta-2\eta^{ij}\zeta_{,i}\zeta_{,j}\right)+\mathcal{O}\left(\zeta^{3}\right)$
$\displaystyle{{}^{3}\\!\widetilde{R}_{ij}}$
$\displaystyle=-\left(D_{j}D_{i}\zeta-\tfrac{1}{3}\Delta\zeta\eta_{ij}\right)+\mathcal{O}(\zeta^{2})$
(124)

where $\Delta=\eta^{ij}D_{i}D_{j}$. Combining these results, the quantities
(98) expanded to second order read

$\displaystyle X_{1}$
$\displaystyle=3\tfrac{\dot{a}}{a}+3\dot{\zeta}-\tfrac{1}{a^{2}}\Delta\chi+\tfrac{2}{a^{2}}\zeta\Delta\chi-\tfrac{1}{a^{2}}\left(D\zeta\right)^{2}$
(125) $\displaystyle X_{2}$
$\displaystyle\doteq\tfrac{1}{a^{4}}\left(\tfrac{2}{3}(\Delta\chi)^{2}-2\varkappa\left(D\chi\right)^{2}\right)$
$\displaystyle X_{3}$
$\displaystyle=\tfrac{1}{a^{2}}\left(6\varkappa-4\Delta\zeta-12\varkappa\zeta+12\varkappa\zeta^{2}+8\zeta\Delta\zeta-2\left(D\zeta\right)^{2}\right)$
$\displaystyle X_{4}$
$\displaystyle\doteq\tfrac{1}{a^{4}}\Delta\chi\left(\tfrac{2}{3}\Delta\zeta+2\varkappa\zeta\right)$
$\displaystyle X_{5}$
$\displaystyle\doteq\tfrac{1}{a^{4}}\left(\tfrac{2}{3}(\Delta\zeta)^{2}-2\varkappa\left(D\zeta\right)^{2}\right)$
$\displaystyle X_{6}$
$\displaystyle\doteq\tfrac{16}{a^{6}}\left(\left(D\Delta\zeta\right)^{2}-6\varkappa\left(\Delta\zeta\right)^{2}+9\varkappa^{2}\left(D\zeta\right)^{2}\right)$
$\displaystyle X_{7}$
$\displaystyle\doteq\tfrac{1}{a^{6}}\left(\tfrac{2}{3}\left(D\Delta\zeta\right)^{2}-6\varkappa\left(\Delta\zeta\right)^{2}+12\varkappa^{2}(D\zeta)^{2}\right)$

with the notation
$\left(D\zeta\right)^{2}\equiv\eta^{ij}\zeta_{,i}\zeta_{,j}$.

Neglecting matter fluctuations of
$\mathcal{S}_{m}=\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,a^{3}\,e^{3\zeta}p$ and
using the background equation of motion (108), the total second order action
${{}^{(2)}\mathcal{S}}={{}^{(2)}\mathcal{S}_{g}+{{}^{(2)}\mathcal{S}_{m}}}$
becomes

$\displaystyle{{}^{(2)}\mathcal{S}}$
$\displaystyle=-\frac{1}{32\pi}\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,a^{3}\bigg{\\{}-\left[L_{3}-\tfrac{3}{a}\partial_{0}\left(aL_{13}\right)-\tfrac{12\varkappa}{a^{2}}L_{33}\right]\,\tfrac{12\varkappa}{a^{2}}\zeta^{2}+$
(126)
$\displaystyle+L_{11}9\dot{\zeta}^{2}-\left[L_{3}-\tfrac{3}{a}\partial_{0}\left(aL_{13}\right)-\tfrac{\varkappa}{a^{2}}(24L_{33}+L_{5})+\tfrac{6\varkappa^{2}}{a^{4}}(12L_{6}+L_{7})\right]\tfrac{4}{a^{2}}\,\zeta\Delta\zeta+\phantom{\bigg{\\}}}$
$\displaystyle+\left[16L_{33}+\tfrac{4}{3}L_{5}-\tfrac{12\varkappa}{a^{2}}\left(16L_{6}+L_{7}\right)\right]\tfrac{1}{a^{4}}\,\zeta\Delta^{2}\zeta-\left[32L_{6}+\tfrac{4}{3}L_{7}\right]\tfrac{1}{a^{6}}\,\zeta\Delta^{3}\zeta\phantom{\bigg{\\}}}$
$\displaystyle+\left[(8L_{13}+\tfrac{4}{3}L_{4})\left(3\varkappa\zeta+\Delta\zeta\right)-6L_{11}a^{2}\dot{\zeta}+(L_{11}+\tfrac{4}{3}L_{2})\Delta\chi+4\varkappa
L_{2}\chi\right]\tfrac{1}{a^{4}}\Delta\chi\bigg{\\}}.$

Subscripts $L_{i}$ denote derivatives of $L$ with respect to $X_{i}$ evaluated
on the background.

Variation with respect to $\Delta\chi$ yields

$(L_{11}+\tfrac{4}{3}L_{2})\Delta\chi+4\varkappa
L_{2}\chi=3L_{11}a^{2}\dot{\zeta}-(4L_{13}+\tfrac{2}{3}L_{4})\left(3\varkappa\zeta+\Delta\zeta\right).$
(127)

In the spatially flat case the second order action for $\zeta$ after
integrating out $\chi$, assuming $L_{11}+\tfrac{4}{3}L_{2}\neq 0$, becomes

$\displaystyle{{}^{(2)}\mathcal{S}_{\zeta}}=-$
$\displaystyle\frac{1}{32\pi}\int\\!\mathrm{d}^{4}x\,a^{3}\bigg{\\{}\frac{12L_{11}L_{2}}{L_{11}+\tfrac{4}{3}L_{2}}\,\dot{\zeta}^{2}-\left[L_{3}-\tfrac{1}{a}\partial_{0}\left(a\frac{4L_{2}L_{13}-\tfrac{1}{2}L_{11}L_{4}}{L_{11}+\tfrac{4}{3}L_{2}}\right)\right]\tfrac{4}{a^{2}}\,\zeta\Delta\zeta\phantom{\bigg{\\}}}$
$\displaystyle+\left[16L_{33}+\tfrac{4}{3}L_{5}-\frac{(4L_{13}+\tfrac{2}{3}L_{4})^{2}}{L_{11}+\tfrac{4}{3}L_{2}}\right]\tfrac{1}{a^{4}}\,\zeta\Delta^{2}\zeta-\left[32L_{6}+\tfrac{4}{3}L_{7}\right]\tfrac{1}{a^{6}}\,\zeta\Delta^{3}\zeta\bigg{\\}}.$
(128)

#### NB.

Note that in the case $L_{33}=L_{5}=L_{6}=L_{7}=0$ it is tempting to make the
choice, $L_{4}=-6L_{13}$ such that

$\displaystyle{{}^{(2)}\mathcal{S}_{\zeta}}=-\frac{1}{32\pi}\int\\!\mathrm{d}^{4}x\,a^{3}\bigg{\\{}$
$\displaystyle\frac{12L_{11}L_{2}}{L_{11}+\tfrac{4}{3}L_{2}}\,\dot{\zeta}^{2}-\left[L_{3}-3\tfrac{1}{a}\partial_{0}\left(aL_{13}\right)\right]\tfrac{4}{a^{2}}\,\zeta\Delta\zeta\bigg{\\}}.$
(129)

However, note that in this case a gradient instability can only be prevented
at the expense of introducing a gradient instability for tensor perturbations.

##  Appendix B Stability analysis of higher order spatial curvature bounces

As seen in [8], the goal of singularity resolution in spatially non-flat
universes requires the introduction of a scalar spatial curvature dependent
potential. By isotropy of the background, there is a degeneracy of higher
order spatial curvature terms that will lead to the same background dynamics.
In the following I will consider the generalization

$\mathcal{L}=\mathcal{L}_{\textup{nhd}}+\mathcal{L}_{\textup{hd}}+V(\widetilde{R})+\alpha^{4}\,\widetilde{R}\left(\widetilde{R}^{\mu\nu}\widetilde{R}_{\mu\nu}-\tfrac{1}{3}\widetilde{R}^{2}\right),$
(130)

where now $\sigma_{T}$, $\sigma_{S}$ and $\alpha$ can depend on
$\widetilde{R}$. Note that in the spatially flat case this Lagrangian is
equivalent to (3).

#### Background.

The modified Friedmann equation with non-vanishing spatial curvature becomes

$\tfrac{2}{3}\left(\kappa\ell^{\prime}-\ell\right)=\tfrac{1}{2}\left(V\left({{}^{3}\\!R}\right)-{{}^{3}\\!R}\right)+8\pi\varepsilon.$
(131)

It is independent of the higher order terms $\sigma_{T}$, $\sigma_{S}$ and
$\alpha$. This modified Friedmann equation has been discussed in [8]. A bounce
is made possible if the term $V\left({{}^{3}\\!R}\right)$ is negative and
dominating over $\varepsilon$ at small scale factor. Assuming a cubic
potential,

$V\left({{}^{3}\\!R}\right)=-6\,\delta^{4}\left(\frac{\varkappa}{6}{{}^{3}\\!R}\right)^{3}=-\frac{6\,\delta^{4}}{a^{6}},$
(132)

this will be satisfied, provided that the matter equation of state satisfies
$w<1$. The bounce then happens at the minimal value of the scale factor

$a_{\min}=\left(\frac{1+3w}{2}\right)^{\tfrac{1+w}{1-w}}\left(\frac{\sqrt{c}}{H_{l}}\right)^{-\tfrac{1+3w}{3(1-w)}}\delta^{\tfrac{4}{3(1-w)}},$
(133)

where (23) was assumed. In the case $\varkappa=1$ a re-collapse happens at

$a_{\max}=\left(\frac{1+3w}{2}\right)^{-\tfrac{3(1+w)}{1+3w}}\frac{\sqrt{c}}{H_{l}}.$
(134)

#### Second order actions.

The second order action for tensor perturbations is

$\displaystyle{{}^{(2)}\mathcal{S}}_{h}=\frac{1}{64\pi}\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,a^{3}\bigg{\\{}$
$\displaystyle
f\,\dot{h}^{ab}\dot{h}_{ab}-\left(1-V^{\prime}+\tfrac{6\varkappa^{2}}{a^{4}}(\sigma_{T}^{4}+2\alpha^{4})\right)\tfrac{2\varkappa}{a^{2}}\,h^{ab}h_{ab}$
(135)
$\displaystyle-\left(1-V^{\prime}+\tfrac{8\varkappa^{2}}{a^{4}}(2\sigma_{T}^{4}+3\alpha^{4})\right)\tfrac{1}{a^{2}}(Dh_{ab})^{2}$
$\displaystyle-(7\sigma_{T}^{4}+6\alpha^{4})\tfrac{\varkappa}{a^{6}}(D^{2}h_{ab})^{2}-\frac{\sigma_{T}^{4}}{a^{6}}(D^{3}h_{ab})^{2}\bigg{\\}},$

where $(Dh_{ab})^{2}\equiv\eta^{cd}D_{c}h^{ab}D_{d}h_{ab}$,
$(D^{2}h_{ab})^{2}\equiv\Delta h^{ab}\Delta h_{ab}$,
$(D^{3}h_{ab})^{2}\equiv\eta^{cd}D_{c}\Delta h^{ab}D_{d}\Delta h_{ab}$.

Neglecting matter perturbations, the second order action for scalar
perturbations, after integrating out $\chi$ is

$\displaystyle{{}^{(2)}\mathcal{S}_{\zeta}}=\frac{1}{8\pi}\int\\!\mathrm{d}^{4}x\sqrt{\eta}\,a^{3}\bigg{\\{}$
$\displaystyle\dot{\zeta}\,\frac{3\ell^{\prime\prime}f\left(\Delta+3\varkappa\right)}{(\ell^{\prime\prime}-f)\Delta-3f\varkappa}\,\dot{\zeta}-\left[\tfrac{1}{3}(1-V^{\prime})+4V^{\prime\prime}\tfrac{\varkappa}{a^{2}}\right]\tfrac{\varkappa}{a^{2}}\,9\zeta^{2}$
(136)
$\displaystyle+\left[1-V^{\prime}+\tfrac{3\varkappa}{a^{2}}(8V^{\prime\prime}+\tfrac{\varkappa}{a^{2}}(2\alpha^{4}-3\sigma_{S}^{4})\right]\tfrac{1}{a^{2}}\,(D\zeta)^{2}+\phantom{\bigg{\\}}}$
$\displaystyle-\left[4V^{\prime\prime}+\tfrac{2\varkappa}{a^{2}}(\alpha^{4}-3\sigma_{S}^{4})\right]\tfrac{1}{a^{4}}\,(D^{2}\zeta)^{2}-\tfrac{\sigma_{S}^{4}}{a^{6}}\,(D^{3}\zeta)^{2}\bigg{\\}}.$

#### Ghost instability?

The condition for no ghost instability of the mode characterized by the
eigenvalue $-k^{2}$ of the curved Laplacian $\Delta$ can be written as

$\frac{\ell^{\prime\prime}f\left(k^{2}-3\varkappa\right)}{(\ell^{\prime\prime}-f)k^{2}+3f\varkappa}>0.$
(137)

In the case $\ell^{\prime\prime}\neq f$, the condition for short wavelength
modes $k^{2}\gg\varkappa$ is the same as (67) in the spatially flat case. On
the other hand, in the limit $\ell^{\prime\prime}-f\to 0$, $f\to 1$ (which
applies also in the region around the bounce where $H$ vanishes) the condition
becomes

$\varkappa k^{2}>3.$ (138)

In an open universe $\varkappa=-1$ the spectrum of the curved Laplacian is
continuous and bounded by $k^{2}\geq 1$. Hence all modes suffer from a ghost
instability. In a closed universe $\varkappa=1$, however, the eigenvalues of
$\Delta$ are discrete and given by [36], [37]

$k^{2}=n^{2}-\varkappa,\qquad n\in\mathbb{N}_{\geq 1}.$ (139)

The discreteness of the spectrum is due to the periodic boundary conditions
that eigenfunctions have to satisfy. Note that only the modes $n=1,2$ would
suffer from a ghost instability. However, these two longest wavelength modes
can be shown to correspond to pure gauge terms [38], [39]. In conclusion, in
the closed case $\varkappa=1$ there is no ghost instability.

#### Gradient instability?

Let us now consider $\varkappa=1$ and take for simplicity the cubic potential
(132) and constant $\sigma_{S}^{4}\geq 0$, $\sigma_{T}^{4}\geq 0$,
$\alpha^{4}\geq 0$. In this case tensor perturbations do not exhibit any
instability and the sixth order term in (136) does not lead to an instability
of scalar perturbations. The condition following from the right sign of the
forth order term reads

$\alpha^{4}\geq 3\sigma_{S}^{4}+2\delta^{4}.$ (140)

From the gradient term we can read off that there is a gradient instability
wherever

$1+\tfrac{3}{a^{4}}(2\alpha^{4}-3\sigma_{S}^{4}-7\delta^{4})>0.$ (141)

Using (140), we see that if $\delta^{4}\leq 2\sigma_{S}^{4}$ this condition is
satisfied at all times and there is a gradient instability throughout.
Conversely, if $\delta^{4}>2\sigma_{S}^{4}$ it is possible to avoid the
gradient instability in the region around the bounce, provided that

$a_{\min}^{4}<3\left(7\delta^{4}+3\sigma_{S}^{4}-2\alpha^{4}\right).$ (142)

In order to avoid the gradient instability at all times, it would be necessary
to have

$a_{\max}^{4}<3\left(7\delta^{4}+3\sigma_{S}^{4}-2\alpha^{4}\right)\overset{!}{<}9\left(\delta^{4}-2\sigma_{S}^{4}\right),$
(143)

where for the last inequality (140) was used. Note, however, that in this case

$\left(\frac{a_{\max}}{a_{\min}}\right)\sim\left(\frac{a_{\max}}{\delta}\right)^{\tfrac{4}{3(1+w)}}\overset{!}{<}\mathcal{O}(10),$
(144)

and hence a gradient instability region is unavoidable for any universe that
undergoes any significant amount of expansion.

## References

  * Bekenstein [1993] J. D. Bekenstein,  _The Relation between physical and gravitational geometry,_Phys. Rev. D 48, 3641 (1993), arXiv:gr-qc/9211017 .
  * Deruelle and Rua [2014] N. Deruelle and J. Rua,  _Disformal Transformations, Veiled General Relativity and Mimetic Gravity,_JCAP 09, 002 (2014), arXiv:1407.0825 [gr-qc] .
  * Domènech  _et al._ [2015] G. Domènech, S. Mukohyama, R. Namba, A. Naruko, R. Saitou, and Y. Watanabe,  _Derivative-dependent metric transformation and physical degrees of freedom,_Phys. Rev. D 92, 084027 (2015), arXiv:1507.05390 [hep-th] .
  * Chamseddine and Mukhanov [2013] A. H. Chamseddine and V. Mukhanov,  _Mimetic Dark Matter,_JHEP 11, 135 (2013), arXiv:1308.5410 [astro-ph.CO] .
  * Golovnev [2014] A. Golovnev,  _On the recently proposed Mimetic Dark Matter,_Phys. Lett. B 728, 39 (2014), arXiv:1310.2790 [gr-qc] .
  * Chamseddine  _et al._ [2014] A. H. Chamseddine, V. Mukhanov, and A. Vikman,  _Cosmology with Mimetic Matter,_JCAP 06, 017 (2014), arXiv:1403.3961 [astro-ph.CO] .
  * Chamseddine  _et al._ [2019a] A. H. Chamseddine, V. Mukhanov, and T. B. Russ,  _Asymptotically Free Mimetic Gravity,_Eur. Phys. J. C 79, 558 (2019a), arXiv:1905.01343 [hep-th] .
  * Chamseddine  _et al._ [2020] A. H. Chamseddine, V. Mukhanov, and T. B. Russ,  _Non-Flat Universes and Black Holes in Asymptotically Free Mimetic Gravity,_Fortsch. Phys. 68, 1900103 (2020), arXiv:1912.03162 [hep-th] .
  * Chamseddine  _et al._ [2019b] A. H. Chamseddine, V. Mukhanov, and T. B. Russ,  _Mimetic Hořava gravity,_Phys. Lett. B 798, 134939 (2019b), arXiv:1908.01717 [hep-th] .
  * Horava [2009] P. Horava,  _Quantum Gravity at a Lifshitz Point,_Phys. Rev. D 79, 084008 (2009), arXiv:0901.3775 [hep-th] .
  * Barvinsky  _et al._ [2016] A. O. Barvinsky, D. Blas, M. Herrero-Valea, S. M. Sibiryakov, and C. F. Steinwachs,  _Renormalization of Hořava gravity,_Phys. Rev. D 93, 064022 (2016), arXiv:1512.02250 [hep-th] .
  * Barvinsky  _et al._ [2019] A. O. Barvinsky, M. Herrero-Valea, and S. M. Sibiryakov,  _Towards the renormalization group flow of Horava gravity in $(3+1)$ dimensions, _Phys. Rev. D 100, 026012 (2019), arXiv:1905.03798 [hep-th] .
  * Germani  _et al._ [2009] C. Germani, A. Kehagias, and K. Sfetsos,  _Relativistic Quantum Gravity at a Lifshitz Point,_JHEP 09, 060 (2009), arXiv:0906.1201 [hep-th] .
  * Mukohyama [2009] S. Mukohyama,  _Dark matter as integration constant in Horava-Lifshitz gravity,_Phys. Rev. D 80, 064005 (2009), arXiv:0905.3563 [hep-th] .
  * Ramazanov  _et al._ [2016] S. Ramazanov, F. Arroja, M. Celoria, S. Matarrese, and L. Pilo,  _Living with ghosts in Hořava-Lifshitz gravity,_JHEP 06, 020 (2016), arXiv:1601.05405 [hep-th] .
  * Cognola  _et al._ [2016] G. Cognola, R. Myrzakulov, L. Sebastiani, S. Vagnozzi, and S. Zerbini,  _Covariant Hořava-like and mimetic Horndeski gravity: cosmological solutions and perturbations,_Class. Quant. Grav. 33, 225014 (2016), arXiv:1601.00102 [gr-qc] .
  * Bianco  _et al._ [2018] S. Bianco, V. N. Friedhoff, and E. Wilson-Ewing,  _Modified dispersion relations, inflation and scale invariance,_Phys. Rev. D 97, 046006 (2018), arXiv:1609.06891 [gr-qc] .
  * Zheng  _et al._ [2017] Y. Zheng, L. Shen, Y. Mou, and M. Li,  _On (in)stabilities of perturbations in mimetic models with higher derivatives,_JCAP 08, 040 (2017), arXiv:1704.06834 [gr-qc] .
  * Zheng [2021] Y. Zheng,  _Hamiltonian analysis of Mimetic gravity with higher derivatives,_JHEP 01, 085 (2021), arXiv:1810.03826 [gr-qc] .
  * Ganz  _et al._ [2019a] A. Ganz, N. Bartolo, and S. Matarrese,  _Towards a viable effective field theory of mimetic gravity,_JCAP 12, 037 (2019a), arXiv:1907.10301 [gr-qc] .
  * Takahashi and Kobayashi [2017] K. Takahashi and T. Kobayashi,  _Extended mimetic gravity: Hamiltonian analysis and gradient instabilities,_JCAP 11, 038 (2017), arXiv:1708.02951 [gr-qc] .
  * Ganz  _et al._ [2019b] A. Ganz, P. Karmakar, S. Matarrese, and D. Sorokin,  _Hamiltonian analysis of mimetic scalar gravity revisited,_Phys. Rev. D 99, 064009 (2019b), arXiv:1812.02667 [gr-qc] .
  * de Cesare [2019] M. de Cesare,  _Reconstruction of Mimetic Gravity in a Non-SingularBouncing Universe from Quantum Gravity,_Universe 5, 107 (2019), arXiv:1904.02622 [gr-qc] .
  * Ijjas  _et al._ [2016] A. Ijjas, J. Ripley, and P. J. Steinhardt,  _NEC violation in mimetic cosmology revisited,_Phys. Lett. B 760, 132 (2016), arXiv:1604.08586 [gr-qc] .
  * Firouzjahi  _et al._ [2017] H. Firouzjahi, M. A. Gorji, and S. A. Hosseini Mansoori,  _Instabilities in Mimetic Matter Perturbations,_JCAP 07, 031 (2017), arXiv:1703.02923 [hep-th] .
  * Hosseini Mansoori  _et al._ [2021] S. A. Hosseini Mansoori, A. Talebian, and H. Firouzjahi,  _Mimetic inflation,_JHEP 01, 183 (2021), arXiv:2010.13495 [gr-qc] .
  * Arroja  _et al._ [2016] F. Arroja, N. Bartolo, P. Karmakar, and S. Matarrese,  _Cosmological perturbations in mimetic Horndeski gravity,_JCAP 04, 042 (2016), arXiv:1512.09374 [gr-qc] .
  * Gorji  _et al._ [2018] M. A. Gorji, S. A. Hosseini Mansoori, and H. Firouzjahi,  _Higher Derivative Mimetic Gravity,_JCAP 01, 020 (2018), arXiv:1709.09988 [astro-ph.CO] .
  * [29] V. Mukhanov, Physical Foundations of Cosmology, Cambridge University Press, 2005.
  * Abbott  _et al._ [2017a] B. P. Abbott _et al._ (LIGO Scientific, Virgo),  _GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral,_Phys. Rev. Lett. 119, 161101 (2017a), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ [2017b] B. P. Abbott _et al._ (LIGO Scientific, Virgo, Fermi-GBM, INTEGRAL),  _Gravitational Waves and Gamma-rays from a Binary Neutron Star Merger: GW170817 and GRB 170817A,_Astrophys. J. Lett. 848, L13 (2017b), arXiv:1710.05834 [astro-ph.HE] .
  * Casalino  _et al._ [2019] A. Casalino, M. Rinaldi, L. Sebastiani, and S. Vagnozzi,  _Alive and well: mimetic gravity and a higher-order extension in light of GW170817,_Class. Quant. Grav. 36, 017001 (2019), arXiv:1811.06830 [gr-qc] .
  * Emir Gümrükçüoğlu  _et al._ [2018] A. Emir Gümrükçüoğlu, M. Saravani, and T. P. Sotiriou,  _Hořava gravity after GW170817,_Phys. Rev. D 97, 024032 (2018), arXiv:1711.08845 [gr-qc] .
  * Chamseddine and Mukhanov [2017] A. H. Chamseddine and V. Mukhanov,  _Resolving Cosmological Singularities,_JCAP 03, 009 (2017), arXiv:1612.05860 [gr-qc] .
  * Brahma  _et al._ [2018] S. Brahma, A. Golovnev, and D.-H. Yeom,  _On singularity-resolution in mimetic gravity,_Phys. Lett. B 782, 280 (2018), arXiv:1803.03955 [gr-qc] .
  * Akama and Kobayashi [2019] S. Akama and T. Kobayashi,  _General theory of cosmological perturbations in open and closed universes from the Horndeski action,_Phys. Rev. D 99, 043522 (2019), arXiv:1810.01863 [gr-qc] .
  * Abbott and Schaefer [1986] L. F. Abbott and R. K. Schaefer,  _A General, Gauge Invariant Analysis of the Cosmic Microwave Anisotropy,_Astrophys. J. 308, 546 (1986).
  * Lifshitz and Khalatnikov [1963] E. M. Lifshitz and I. M. Khalatnikov,  _Investigations in relativistic cosmology,_Adv. Phys. 12, 185 (1963).
  * Bardeen [1980] J. M. Bardeen,  _Gauge Invariant Cosmological Perturbations,_Phys. Rev. D 22, 1882 (1980).