On nonlocally interacting metrics, and a simple proposal for cosmic acceleration

TL;DR

This paper introduces a simple nonlocal modification to general relativity that can explain late-time cosmic acceleration without a cosmological constant, fitting observational data and maintaining physical degrees of freedom.

Contribution

The paper proposes a new nonlocal gravity model inspired by the Deser-Woodard model, which successfully models cosmic acceleration with viable background expansion and healthy degrees of freedom.

Findings

Model fits observational data with m ~ H_0

Cosmic expansion history closely matches standard cosmology

Model avoids physical ghosts despite local formulation issues

Abstract

We propose a simple, nonlocal modification to general relativity (GR) on large scales, which provides a model of late-time cosmic acceleration in the absence of the cosmological constant and with the same number of free parameters as in standard cosmology. The model is motivated by adding to the gravity sector an extra spin-2 field interacting nonlocally with the physical metric coupled to matter. The form of the nonlocal interaction is inspired by the simplest form of the Deser-Woodard (DW) model, $\alpha R\frac{1}{\Box}R$, with one of the Ricci scalars being replaced by a constant $m^{2}$, and gravity is therefore modified in the infrared by adding a simple term of the form $m^2\frac{1}{\Box}R$ to the Einstein-Hilbert term. We study cosmic expansion histories, and demonstrate that the new model can provide background expansions consistent with observations if $m$ is of the order of the Hubble expansion rate today, in contrast to the simple DW model with no viable cosmology. The model is best fit by $w_0\sim-1.075$ and $w_a\sim0.045$. We also compare the cosmology of the model to that of Maggiore and Mancarella (MM), $m^2R\frac{1}{\Box^2}R$, and demonstrate that the viable cosmic histories follow the standard-model evolution more closely compared to the MM model. We further demonstrate that the proposed model possesses the same number of physical degrees of freedom as in GR. Finally, we discuss the appearance of ghosts in the local formulation of the model, and argue that they are unphysical and harmless to the theory, keeping the physical degrees of freedom healthy.

Equations

$$S=\frac{M_{\mathrm{Pl}}^2}{2}\int {\text{d}}^{4}x\sqrt{-g}\left(R+R{h}_1(△)R+R^{\alpha \beta }{h}_2(△)R_{\alpha \beta }+R^{\mu \nu \alpha \beta }{h}_3(△)R_{\mu \nu \alpha \beta }\right)+S_{\text{matter}},$$

$$S=\frac{M_g^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}{R}_g+\frac{M_f^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-f}{R}_f-I_{g,f}+S_{\text{matter}}[g,\Psi ],$$

$$\int {\mathrm{d}}^{D}x\varphi (x){\square }^{-1}\psi (x)$$

$$=\int {\mathrm{d}}^{D}x{\mathrm{d}}^{D}y\psi (y)G^T(y,x)\varphi (x)$$

$$=\int {\mathrm{d}}^{D}y\psi (y)({\square }^{-1}{)}^T\varphi (y)$$

$$\equiv \int {\mathrm{d}}^{D}y\psi (y){\square }^{-1}\varphi (y),$$

$$S=\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}R+\frac{M_f^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-f}{R}_f-\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}\alpha (R_f\frac{1}{\square }R+R\frac{1}{\square }{R}_f)+S_{\text{matter}}[g,\Psi ],$$

$$G_{\mu \nu }+\Delta G_{\mu \nu }$$

$$G_{\mu \nu }^f+\Delta G_{\mu \nu }^f$$

$$\Delta G_{\mu \nu }=$$

$$+{\nabla }_{(\mu }(\frac{1}{\square }{R}_f){\nabla }_{\nu )}(\frac{1}{\square }R)],$$

$$\Delta G_{\mu \nu }^f=$$

$$U$$

$$V$$

$$S=$$

$$+\int {\mathrm{d}}^{4}x\sqrt{-g}{\lambda }_1(R-\square U)+\int {\mathrm{d}}^{4}x\sqrt{-g}{\lambda }_2(R_f-\square V)+S_{\text{matter}}[g,\Psi ],$$

$${\lambda }_1$$

$${\lambda }_2$$

$$S=$$

$$+\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}2\alpha V\square U+S_{\text{matter}}[g,\Psi ].$$

$$\alpha$$

$$f_{\mu \nu }$$

$$V$$

$$\Delta G_{\mu \nu }$$

$$\Delta G_{\mu \nu }^f$$

$$G_{\mu \nu }{\nabla }^{\mu }(\frac{1}{\square }{R}_f)+[1-\frac{1}{2}(\frac{1}{\square }R)]{\nabla }_{\nu }{R}_f-[R_{\rho \nu }{\nabla }^{\rho }(\frac{1}{\square }{R}_f)+{\nabla }_{\nu }{R}_f]+\frac{1}{2}R{\nabla }_{\nu }(\frac{1}{\square }{R}_f)$$

$$=-\frac{1}{2}(\frac{1}{\square }R){\nabla }_{\nu }{R}_f=0.$$

$${\nabla }_{\nu }{R}_f=0.$$

$${\nabla }_f^{\mu }(\sqrt{f^{-1}g}(\frac{1}{\square }R)R_{\mu \nu }^f)+{\nabla }_{\nu }^f{\square }_f(\sqrt{f^{-1}g}(\frac{1}{\square }R))-[R_{\rho \nu }^f{\nabla }^{\rho }(\sqrt{f^{-1}g}(\frac{1}{\square }R))+{\nabla }_{\nu }{\square }_f\sqrt{f^{-1}g}(\frac{1}{\square }R)]$$

$$={\nabla }_f^{\mu }(\sqrt{f^{-1}g}(\frac{1}{\square }R)R_{\mu \nu }^f)-R_{\rho \nu }^f{\nabla }^{\rho }(\sqrt{f^{-1}g}(\frac{1}{\square }R))=\sqrt{f^{-1}g}(\frac{1}{\square }R){\nabla }_f^{\mu }{R}_{\mu \nu }^f=0.$$

$${\nabla }_f^{\mu }{G}_{\mu \nu }^f=0\Rightarrow {\nabla }_f^{\mu }{R}_{\mu \nu }^f-\frac{1}{2}{\nabla }_{\nu }^f{R}_f=0.$$

$${\nabla }_{\nu }^f{R}_f=0.$$

$${\partial }_{\mu }{R}_f=0.$$

$$S=\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}[R+\frac{1}{2}(m^2\frac{1}{\square }R+R\frac{1}{\square }{m}^2)]+S_{\text{matter}}[g,\Psi ],$$

$$S=\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}(R+m^2\frac{1}{\square }R)+S_{\text{matter}}[g,\Psi ].$$

$$U\equiv \frac{1}{\square }R.$$

$$S=\frac{M_{\text{Pl}}^2}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}(R+m^{2}U)+\int {\mathrm{d}}^{4}x\sqrt{-g}\lambda (R-\square U)+S_{\text{matter}}[g,\Psi ].$$

$$(\frac{M_{\text{Pl}}^2}{2}+\lambda )G_{\mu \nu }+\frac{M_{\text{Pl}}^2}{2}{m}^2{g}_{\mu \nu }(1-\frac{U}{2})-{\nabla }_{\mu }{\nabla }_{\nu }\lambda -\frac{1}{2}{g}_{\mu \nu }{\nabla }^{\rho }\lambda {\nabla }_{\rho }U+{\nabla }_{(\mu }\lambda {\nabla }_{\nu )}U=\frac{1}{2}{T}_{\mu \nu },$$

Full text

1 Introduction

The question of why the late-time expansion of the Universe is accelerating is now almost twenty years old, with strong and overwhelming evidence supporting the phenomenon [1] after its initial discovery through the observations of supernovae [2, 3] (see Refs. [4, 5, 6, 7] for recent reviews). The standard model of cosmology, $\Lambda$CDM ($\Lambda$ for the cosmological constant (CC), and CDM for cold dark matter), provides a strikingly successful description of cosmic acceleration in the arguably simplest possible way, i.e., through a single parameter, $\Lambda$. This observationally very well tested model, however, suffers from serious, theoretical issues stemming from a tremendous fine-tuning that is required for compatibility of the observed value of $\Lambda$ with the widely accepted principles of quantum field theory; see, e.g., Ref. [8] for a review. Alongside the other problems with $\Lambda$CDM, different incarnations of the cosmological constant problem have been considered as strong motivations for exploring alternative cosmological models and, therefore, going beyond $\Lambda$CDM [7]. One of the interesting such attempts consists of constructing alternative theories of gravity which would offer mechanisms for cosmic acceleration that are different from a simple cosmological constant in the framework of general relativity (GR); see, e.g., Refs. [9, 10]. Such models must clearly be consistent with the tide of various high-quality cosmological data, at least as well as $\Lambda$CDM, be theoretically well defined and well motivated, (ideally) be simple (i.e. without introducing many free parameters), and offer predictions that are distinguishable from those of $\Lambda$CDM, making the models testable and falsifiable.

One class of interesting models of modified gravity that has attracted significant attention over the past few years is the one with the gravity sector extended by adding an extra rank-2 tensor field similar to the fundamental dynamical field of GR, the metric, describing interacting spin-2 fields. These commonly called bimetric theories have a long history, in connection to massive gravity where gravitons are assumed to possess a nonzero mass, contrary to GR which is the unique theory of massless gravitons [11, 12, 13, 14, 15]. After decades of intensive searches for theoretically consistent theories of nonlinear massive and bimetric gravity, it eventually became possible, after the discovery of ghost-free massive gravity [16, 17, 18, 19, 20, 21, 22, 23, 24, 25] and bigravity [26] (see Refs. [27, 28, 29, 30, 31] for reviews), to explore the cosmological implications of such theories, in particular in connection to cosmic acceleration. It was quickly realized, after the successful construction of bigravity, that it admits Friedman-Lemaître-Robertson-Walker (FLRW) cosmological solutions111See Ref. [32] and references therein for the cosmology of bimetric models with other choices of metrics. which agree with observations at the background level, i.e. they successfully describe the cosmic expansion history even in the absence of an explicit cosmological constant (or vacuum energy) term [33, 34, 35, 36, 37, 38, 39, 40].222See also Ref. [41] for viable background cosmologies of theories with more than two spin-2 fields. It however turned out that the linear cosmological perturbations, investigated extensively in Refs. [42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57], suffer from either ghost or gradient instabilities.333Here, we have referred to theories where matter couples only to one of the two spin-2 fields, the physical metric, and the second metriclike field is considered only as an extra dynamical tensor field interacting with the metric. See Refs. [58, 59, 60, 61, 62, 63, 63, 64, 39, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75] for cases where matter couples to more than one metric directly, or to a composite metric. Although a few potential ways out have been proposed (see, e.g., Refs. [76, 77]), it is still an open question whether any models of only two interacting spin-2 fields with self-accelerating solutions exist that are fully stable linearly at all times, and provide a standard isotropic and homogeneous background evolution. There have been various attempts at constructing interacting-metric theories with kinetic and/or interaction terms other than the ones in the original ghost-free nonlinear theories of massive and bimetric gravity (especially with derivative interactions) [78, 79, 80, 81, 82, 83, 84, 85, 86], but it has proven difficult, if not impossible, to find such new terms that do not revive the so-called Boulware-Deser ghost [87].

Another class of interesting alternative theories of gravity proposed as solutions to the cosmic acceleration problem is the one with nonlocal terms added to the Einstein-Hilbert term of GR action. No matter which definition we choose for general relativity, either the geometrical picture according to Lovelock’s theorem [88] or the quantum-field-theoretical picture in terms of massless spin-2 fields (see Ref. [28]), locality is one of the fundamental assumptions of GR. Clearly, one way of modifying GR is therefore to break this assumption. The appearance of nonlocal terms at low energies (infrared; IR) is a generic feature in effective field theories when massless or light degrees of freedom are integrated out [89, 90, 91, 92, 93, 94, 95, 96, 97], and may also arise more fundamentally in Euclidean quantum gravity [98, 99]. Effective actions with IR nonlocal terms have also been found for theories of massive gravity [100, 101], multimetric gravity [102], and post-Riemannian, affine geometry [103]. IR nonlocalities are usually modelled at the level of the action by adding terms that involve inverse differential operators, such as inverse d’Alembertian $\frac{1}{\Box}$ (or $\Box^{-1}$), which in Fourier space can be considered as a Feynman propagator describing the effects of the integrated-out fields. Such operators modify gravitational interactions at large temporal and spatial scales, and can therefore provide dynamical mechanisms for cosmic acceleration. Another important motivation behind these types of nonlocal modifications in the IR stems from the observation [104, 105] that such operators could provide an appealing solution to the old cosmological constant problem by degravitating a large vacuum energy. Even though no consistent theory of degravitation has been found yet in this context and at the level of the action, degravitation remains an important and inspiring motivation for nonlocal modifications to gravity in the IR. And finally, the fact that many of the no-go theorems for gravity rely on the locality of the action is another motivation to relax this condition. This then opens up a large number of new possibilities for model-building. It is therefore important to try to construct simple forms of nonlocal actions and study their implications in different regimes.

The first nonlocal mechanism for cosmic acceleration was proposed by Deser and Woodard (DW) [106] (see Ref. [107] for a review),444Strictly speaking, the cosmology of nonlocal models was proposed and studied first in Ref. [98] for models with similar structures, although with no connection to cosmic acceleration, which was not yet discovered at the time. where a simple term of the form $Rf(\frac{1}{\Box}R)$ was added to the standard Einstein-Hilbert term in GR. The function $f$ can have any arbitrary form at the phenomenological level. One very interesting feature of this model is that it does not introduce any new mass scale in the gravity sector, contrary to $\Lambda$CDM where $\Lambda$ is a dimensionful quantity with an observed value far smaller than the other scale in the theory, i.e. the Planck mass, leading to an enormous, unnatural hierarchy that requires extreme fine-tuning. The absence of such a new scale in the DW theory is therefore a highly appealing feature, if the theory would be able to explain the late-time acceleration in a way consistent with observations and without extremely fine-tuned dimensionless parameters. In addition, the model has been proven to not add extra excitations to gravity beyond those of general relativity, i.e. the number of physical degrees of freedom are the same as in GR [108]. Unfortunately, though, the simplest form of the function $f$, i.e. $\alpha\frac{1}{\Box}R$, with $\alpha$ a dimensionless free parameter, does not provide viable cosmic histories, i.e. even at the cosmological background level [109] (see also Ref. [110]). This particular form of $f$ is interesting not only because of its simplicity, but also because the localized formulation of the theory requires only one additional scalar field which is not a ghost; theories with other forms of the function $f$ introduce two scalar fields, one of which is a ghost. Such ghosts have however been argued to be harmless as the fields are only auxiliary and do not add to and do not affect the physical degrees of freedom (which are the same as in GR) by converting them to ghosts [108]. Accepting more complicated forms of $f$, one can show that it is indeed possible to phenomenologically tune it such that any cosmic history can be reconstructed, even an exact $\Lambda$CDM background [109, 111]. The $\Lambda$CDM-equivalent form of the function has however turned out to be highly contrived with several free parameters that need to be fixed observationally, making the model less appealing. The model with the reconstructed $\Lambda$CDM background has been further investigated by studying linear perturbations and structure formation [112, 113, 114, 115], and even though it was originally claimed [113] that the model was strongly ruled out observationally, a counter-claim has recently been made [115] stating that not only is the model consistent with data, it even gives a better fit than $\Lambda$CDM. The origin of the disagreement is not yet known, but there are reasons to believe that it may be related to the framework in which the analysis has been done; the former performs all the calculations in the nonlocal formulation of the theory, while the latter studies the model in its local formulation [116].555In this case, the model contains two auxiliary, scalar fields, and its perturbative analysis resembles that of multi-scalar-tensor theories (see Ref. [117], and references therein). Whether or not the DW model provides an observationally viable model for cosmic evolution, the fact that it requires a highly contrived and ad hoc form of the function $f$ with several free parameters renders it both theoretically and phenomenologically difficult to accept as an interesting alternative to $\Lambda$CDM.

Independently of the DW theory, Maggiore and Mancarella (MM) proposed [118] an alternative model of nonlocal gravity in an attempt to explain cosmic acceleration without a cosmological constant. In the MM model, a term of the form $m^{2}R\frac{1}{\Box^{2}}R$ has been added to the Einstein-Hilbert term. Such a model, contrary to the DW model, requires the introduction of a new, fine-tuned, mass scale into the gravity sector, and in addition, even the simple form of the theory with no complicated function introduced requires two scalar fields in order to localize the theory, one of which is inevitably a ghost. It is however argued, similarly to the DW case, that the presence of the ghost in the theory is not dangerous as the ghostly scalar field is only an auxiliary one with no effects on the physical degrees of freedom of the nonlocal theory, which are the same as in general relativity [119]. This harmlessness of the ghost is guaranteed by fixing the initial conditions of the auxiliary fields such that the localized theory becomes equivalent to the original nonlocal MM theory. One can therefore work with the localized formulation as long as the initial conditions are properly chosen, and therefore the localized theory is used only as a mathematical trick for dealing with computations which otherwise prove difficult in the nonlocal formulation. Although the equivalence of the local and nonlocal formulations has been shown at the cosmological background level [119], it is not fully clear whether it remains so also at the level of the perturbations, as there are already reasons to suspect it to be the case by considering the similar model of DW with the seemingly discrepant results in the perturbative analyses of Refs. [113] and [115]. In addition, one needs to be cautious when trying to quantize the model in the localized form, as if the constraints on the auxiliary fields are not properly taken into account in the quantization procedure, the ghost will render the theory unviable [119, 120]. The cosmological implications of the MM model has been extensively studied and the model has proven to provide cosmic histories, as well as structure formation, consistent with observational data, although being different from that of $\Lambda$CDM [118, 121, 122, 123, 124, 125, 126]. For that, the new mass scale of the model $m$ has to be of a similar order of magnitude as the present value of the Hubble parameter, similarly to the CC term in $\Lambda$CDM. This phenomenologically favored value has however been argued to have not emerged from perturbative quantum loop corrections due to integrating out light fields, and a more complex mechanism must be behind the generation of the nonlocalilty of the MM form [96].

It should be noted that other nonlocal models with structures similar to those of the DW and MM models have also been proposed, where the nonlocal distortion term is built out of the Ricci scalar and exponential functions of $\Box^{-1}$ [127], or tensorial objects, such as the Ricci tensor $R_{\mu\nu}$ or Riemann tensor $R_{\mu\nu\alpha\beta}$, rather than the Ricci scalar $R$. Tensorial nonlocalities are theoretically very interesting, as they could, for example, help alleviate the ultraviolet divergences of GR by modifying the graviton propagator [128], or help implement a consistent degravitation mechanism (see, e.g., Refs. [105, 129] for such attempts in the framework of massive gravity), which is not possible through only the introduction of scalar nonlocalities. Unfortunately, though, tensorial nonlocalities have turned out to generically contain rapidly growing modes that prevent them from providing stable background expansions [130, 131, 132, 133].

Summarizing all the possibilities mentioned above, one can write down the most general nonlocal action, quadratic in the curvature invariants, as [128, 133]

[TABLE]

where $\triangle$ is some differential operator, usually assumed to be the operator $\Box$, $M_{\rm Pl}$ is the reduced Planck mass, and $S_{\text{matter}}$ is the matter action. $h_{1}$, $h_{2}$ and $h_{3}$ are arbitrary functions of the operator $\triangle$ involving negative powers of $\triangle$.666Note that the functions can involve positive powers, in which case the nonlocal modifications to GR are in the ultraviolet. See Refs. [134, 135] for recent progress in the construction of ghost-free, ultraviolet, nonlocal gravity. For $\triangle=\Box$, action (1.1) is the most general parity-invariant quadratic curvature action.777This action contains the MM model with the choice of $h_{1}(\Box)=m^{2}\frac{1}{\Box^{2}}$ and the simplest DW model with $h_{1}(\Box)=\alpha\frac{1}{\Box}$. Although the DW models with sophisticated forms of the function $f(\frac{R}{\Box})$ are not included in action (1.1), we do not consider them here, as we choose to keep the structure of the theory as simple as possible. As discussed above, unfortunately a large (and the most interesting) part of the action has been proven to be problematic, for one reason or another. Although specific forms of the terms are phenomenologically viable and theoretically consistent, such as the MM model of $m^{2}R\Box^{-2}R$, these are only in the scalar sector, and for only restricted forms of the free function $h_{1}$. It is therefore natural, and important, to ask whether there are ways to modify the action and expand its viability to more (and specifically tensorial) models, while keeping the structure of the action and all its interesting features intact.

In this paper, we propose one way to do this, in the framework of bimetric theories of gravity, where the above action is modified by only assuming that one of the two curvature quantities in each pair of $R$, $R_{\mu\nu}$, and $R_{\mu\nu\alpha\beta}$ in the nonlocal terms corresponds to an extra spin-2 field $f_{\mu\nu}$. This is in a sense a minimal modification, as the structure of the action is kept almost the same as the original one, and only the field content is changed through a minimal, single, additional tensor field that accounts, in a unified way, for the scalar and tensorial structures required for the curvature quantities in all the terms. Such models are natural to construct in the framework of interacting spin-2 theories, and they, therefore, reside at the interface of bimetric and nonlocal theories. From the point of view of bimetric theories, given the stringent constraints on the form of possible, consistent, local interactions between two metrics (or, more correctly, a metric and an extra spin-2 field), as we discussed above, here it is natural to ask whether new consistent interactions are possible if we relax the locality condition in constructing theories of gravity. There are various no-go theorems about interacting spin-2 fields, massive or massless, but in all of those theorems locality has been assumed, one way or another. There have already been studies in the literature, see e.g. Refs. [136, 137], where there are hints that such no-go theorems could be evaded by including nonlocal effects. The aim of the present paper is however not to study such no-go theorems in detail and to investigate whether (and how) they can be evaded in nonlocal models. We only take a phenomenological approach here and simply try to study the cosmological implications of a simple, phenomenological model of nonlocally interacting spin-2 fields. We hope though that our work will trigger more theoretical work in the future on such possibilities, where the no-go theorems and the mathematical consistency will be studied rigorously for such theories. Although perhaps the most interesting class of nonlocal models constructed this way are the ones with tensorial terms in the action, as we explained above, we start our investigation of this new direction by only considering the interaction terms that contain only the Ricci scalars, i.e. the terms corresponding to the function $h_{1}$ in the action (1.1), and devote the present paper to only this sector of the full action. We see this as only a first and the simplest step towards the construction of theories of nonlocally interacting metrics with potentially very interesting implications, and leave a detailed investigation of the tensorial terms for future work.

We start the paper with investigating a simple model in which two metrics interact nonlocally through their Ricci scalars only, while each has its own Einstein-Hilbert kinetic term. The structure of the nonlocal interaction is inspired by arguably one of the simplest possible nonlocal terms that one could construct out of the curvature $R$ and the operator $\Box$, i.e., the simplest version of the DW model, $\alpha R\frac{1}{\Box}R$, with all its interesting properties mentioned above. We derive all the field equations for the model, and study the implications of the Bianchi constraints and the conservation of the matter stress-energy tensor. We show that the consistency of the field equations and Bianchi conditions places strong constraints on the properties of the extra spin-2 field in the simplest model considered in this paper, suggesting, strongly, that tensorial interactions must also be added if the theory is to be considered as a bimetric setup with the tensorial properties of an extra spin-2 field involved. That being said, the two-metric theory that we start with leads us to a new, single parameter, single-metric, nonlocal modification to general relativity as another alternative to $\Lambda$CDM, which provides a simple mechanism for cosmic acceleration. The model adds a term of the form $m^{2}\frac{1}{\Box}R$ to GR, and although it is inspired by our first attempt at constructing a nonlocal model of interacting metrics, we can consider it as a purely phenomenological model that could originate from other theoretical frameworks. We show that, contrary to the DW model with also the simple, nonlocal operator $\frac{1}{\Box}$ in its structure, ours provides viable cosmological expansion histories in its simple form. One interesting feature of the model is that, similarly to the DW and MM cases, the nonlocal interaction in the model does not add any new physical degrees of freedom to the noninteracting theory. In addition, even though our model now needs two auxiliary fields to be localized (as opposed to the DW $\alpha R\frac{1}{\Box}R$ model which needs one) and one of the fields is a ghost, they do not affect the physical degrees of freedom by converting them to ghosts. In this respect, our model behaves similarly to the MM nonlocal model of $m^{2}R\frac{1}{\Box^{2}}R$. The model, although resembling the DW model $\alpha R\frac{1}{\Box}R$ with the simplest form of the function $f$, where the quantity $\alpha R$ is replaced by a constant $m^{2}$, has a very different phenomenology, consistent with the observed cosmic evolution while avoiding issues such as sudden future singularities present in the DW models [109]. In this respect, the model is more appealing than the DW models where the $f$ function is constructed in a contrived way, with several free parameters, to describe the cosmic evolution; our model introduces only one free parameter just as in $\Lambda$CDM. In comparison to the MM model, on the other hand, it is arguably a simpler model, as it includes the operator $\frac{1}{\Box}$, rather than $\frac{1}{\Box^{2}}$. Even though the localized formulation of the model requires two auxiliary fields, just as in the MM model, we show that it provides a different cosmic evolution, still different from that of $\Lambda$CDM. It is easier to see the connection of our model to both DW and MM models if we write the nonlocal term in each case in terms of the quantity $X\equiv\frac{1}{\Box}R$. In that case, the “viable" DW model is of the form $Rf(X)$ with $f(X)=a_{1}[\tanh(a_{2}(X+a_{5})+a_{3}(X+a_{5})^{2}+a_{4}(X+a_{5})^{3})-1]$, where $a_{1},...,a_{5}$ are free parameters to be set by observations [109, 111], and the MM model is of the form $m^{2}X^{2}$. The model proposed in this paper is then of the simpler form $m^{2}X$. We show that this simple model is able to provide a cosmic history that is consistent with the observed one while being different from both MM and $\Lambda$CDM models. In the present paper, we study only the background evolution of the Universe, and leave the investigation of perturbations for future work.

This paper is organized as follows. In Sec. 2, we start our investigation of a nonlocal model with two metrics interacting through their Ricci scalars by first presenting the action. We derive the gravitational field equations in Sec. 2.1 in the original, nonlocal formulation. We then localize the theory in Sec. 2.2 by introducing two auxiliary fields, and write down the action, as well as the field equations, in their local forms. Bianchi constraints for both metrics (without assuming any specific forms) are presented in Sec. 2.3, which will be shown to place a strong constraint on the structure of the model, implying that the reference metric (or the second spin-2 field of the model) should have a constant scalar curvature, both spatially and temporally. This leads us to the simple and phenomenologically interesting, single metric, $m^{2}\frac{1}{\Box}R$ model, which we introduce in Sec. 3. After presenting all equations of motion, we begin our investigation of the cosmological implications of the model in Sec. 4, by first deriving in Sec. 4.1 the general evolution equations for background dynamics of the Universe, i.e. the equations equivalent to Friedmann equations in $\Lambda$CDM. We then discuss in Sec. 4.2 the ability of the model in providing cosmic histories and expansion evolutions consistent with observations, in particular how cosmic acceleration can be obtained with no need for a cosmological constant term. We also compare the implications of the model to those of $\Lambda$CDM and Maggiore-Mancarella models, and in Sec. 4.3, we further compare its cosmology to that of the Deser-Woodard $\alpha R\frac{1}{\Box}R$ model. In Sec. 4.4, we discuss the cosmological solutions for the original, two-metric model, focusing on the existence of background solutions for the reference metric, taking into account the implications of the Bianchi constraint. After presenting the cosmology of our model, we discuss in Sec. 5 the appearance of ghosts in the local formulation of the model, and argue that, like in other nonlocal models, such ghosts are of no harm to the single-metric, nonlocal, $m^{2}\frac{1}{\Box}R$ model. In particular, we prove explicitly in Sec. 5.1 that the local formulation contains a ghost, but we show in Sec. 5.2, by analyzing the model in the nonlocal formulation, that the number of physical degrees of freedom is the same as in GR, and their healthiness is not affected by the presence of the nonlocal terms. We also discuss in Sec. 5 the number of degrees of freedom and the issue of ghosts in the nonlocaly-interacting-metric model, and argue that these are more subtle in that case compared to the single-metric, $m^{2}\frac{1}{\Box}R$ model. We argue that this two-metric model possesses 2+2 degrees of freedom linearly and around the cosmological solutions studied in the present work, and it may be the case that the full, nonlinear model contains more degrees of freedom, implying the existence of ghosts in the theory or that it is infinitely strongly coupled around cosmological backgrounds, which may be considered as a major issue for the model. Finally, our conclusions and some discussions are presented in Sec. 6 with suggestions for future work. Appendix A briefly discusses some generalizations of the model, and presents the ghost-free condition in such models.

2 Nonlocally interacting spin-2 fields and $\alpha R_{f}\frac{1}{\Box}R$ term

As described above, our goal is to construct a model of gravity where two metrics (or spin-2 fields) interact nonlocally. This inevitably means that the main ingredients of our model should be two metrics that we call $g_{\mu\nu}$ and $f_{\mu\nu}$. The action of the model should then generically include three main pieces, as in other gravity models: kinetic terms for the metrics, interaction terms between the two, and their couplings to matter. Adhering to the standard recipe for constructing gravity theories, and in regard to our discussion in Sec. 1 on only considering scalar interactions in this paper, we start building our model using three (simplest possible) types of ingredients, i.e. the Ricci scalars $R_{g}$ and $R_{f}$ for $g_{\mu\nu}$ and $f_{\mu\nu}$, respectively, the volume elements $\mathrm{d}^{4}x\sqrt{-g}$ and $\mathrm{d}^{4}x\sqrt{-f}$ with $g$ and $f$ being determinants of the two metrics, and the nonlocal operators $\Box_{g}^{-1}$ and $\Box_{f}^{-1}$, where $\Box_{g}$ and $\Box_{f}$ are the d’Alembertian operators corresponding to $g_{\mu\nu}$ and $f_{\mu\nu}$, respectively. With these elements, we suggest an action of the form

[TABLE]

where the first two terms are the Einstein-Hilbert terms, in charge of giving dynamics to the metrics, and $I_{g,f}$ is the interaction term. Before we introduce the form of $I_{g,f}$, let us discuss our choices for the $g_{\mu\nu}$ and $f_{\mu\nu}$ kinetic terms and the matter action. First of all, following the common recipe in constructing bimetric theories, we assume that $g_{\mu\nu}$ is the physical metric, which couples to matter and is used for measuring distances and time intervals, and that the reference metric $f_{\mu\nu}$ is only an extra spin-2 field which interacts directly only with $g_{\mu\nu}$ and not with matter. That is why the matter action involves only $g_{\mu\nu}$ and matter fields (collectively denoted by $\Psi$). In addition, we have assumed the standard Einstein-Hilbert form for the kinetic terms of $g_{\mu\nu}$ and $f_{\mu\nu}$ with the volume elements $\mathrm{d}^{4}x\sqrt{-g}$ and $\mathrm{d}^{4}x\sqrt{-f}$, respectively. There are various reasons for considering such terms, which follow the arguments in the literature for the consistency and healthiness of kinetic terms [78, 26, 81, 83, 84]. $M_{g}$ and $M_{f}$ are the two (reduced) Planck masses corresponding to $g_{\mu\nu}$ and $f_{\mu\nu}$, respectively, and since we have coupled only $g_{\mu\nu}$ to matter, we assume that $M_{g}$ is equivalent to the standard Planck mass in GR, i.e. $M_{\text{Pl}}$.888We show later in Sec. 2.2 that $M_{f}$ is redundant and is not a free parameter of the model.

Let us now turn to the interaction term $I_{g,f}$. Following our guiding principle of building the simplest possible nonlocal terms for the interaction of the two metrics using the ingredients mentioned above, we construct the $I_{g,f}$ term out of structures of the form $R_{f,g}\frac{1}{\Box_{f,g}}R_{f,g}$. These terms resemble the simplest form of the DW nonlocal theory, and are the simplest, scalar, nonlocal interactions one could construct along the lines of the structure of the action (1.1). We emphasize here again that although the interaction terms involving tensorial quantities $R_{\mu\nu}$ and $R_{\mu\nu\alpha\beta}$ are more interesting and more natural to add for the effects of the reference metric as a tensor field to be fully present, we restrict our investigation in the present paper to only the scalar sector, and leave the analysis of the tensorial nonlocalities for future work. That being said, even the set of simple operators chosen here gives us several possibilities, depending on which combinations of $R_{f,g}$ and $\Box_{f,g}$ we choose. Clearly, one possibility is to include all the terms at this phenomenological level, however since we are interested in the simplest possible model with the least number of free parameters, here we pick only a specific subset of the operators based on some phenomenological reasons.

As discussed above, nonlocal infrared modifications are believed to generically emerge from integrating out light degrees of freedom. For example, if one considers a scalar field with a canonical kinetic term and a nonminimal coupling to gravity, nonlocal structures of the DW and MM forms appear in the effective theory after integrating out the scalar field (see, e.g., Ref. [97]). In order to generate nonlocal terms involving the Ricci scalars of both $g_{\mu\nu}$ and $f_{\mu\nu}$, we can assume that a light scalar field $\varphi$ couples nonminimally to both metrics, giving rise to terms of the above structure when the scalar field is integrated out. In addition, since $g_{\mu\nu}$ is the spin-2 field coupled to matter, it seems more natural to assume that the kinetic term of the scalar field is defined through $g_{\mu\nu}$, i.e. is of the form $g^{\mu\nu}\nabla^{(g)}_{\mu}\varphi\nabla^{(g)}_{\nu}\varphi$. Integrating this field out then leads to the appearance of the inverse of the $\Box_{g}$ operator, rather than $\Box_{f}$. Independently of this particular way of generating nonlocalities, i.e. through integrating out light fields, this choice can also be motivated purely phenomenologically by noticing again that $g_{\mu\nu}$ is the physical metric, suggesting the appearance of differential operators in the action in terms of $g_{\mu\nu}$ rather than $f_{\mu\nu}$.999One may argue that since $g_{\mu\nu}$ is the “physical metric,” it is more natural to trace $(R_{f})_{\mu\nu}$ with $g^{\mu\nu}$ instead of $f^{\mu\nu}$, and therefore construct the action with terms like $(R_{f})_{\mu\nu}g^{\mu\nu}$ instead of $R_{f}$. The argument is based on the fact that the Ricci tensor is the most basic geometrical object as it contains only the connection and not the metric directly. This is certainly a possibility, and perhaps even more natural, but here we have chosen to work with $R$ and $R_{f}$ in the present work purely in order to keep the structure of the gravity sector of the theory as symmetric as possible in terms of the two metrics $g_{\mu\nu}$ and $f_{\mu\nu}$; this clearly need not be the case, and we leave the investigation of such possibilities for future work.’ Similarly, we have assumed that both connections of the metrics are torsion-free, which again is only a simplicity assumption. One should however keep in mind that since here there is a curvature interaction between the two tensors, this is no longer necessarily the natural choice. This then eliminates half of the possibilities for the nonlocal interaction terms. Keeping only the terms that involve both $R_{g}$ and $R_{f}$ (for explicit interactions between the two metrics), we are left with the two terms $R_{g}\frac{1}{\Box_{g}}R_{f}$ and $R_{f}\frac{1}{\Box_{g}}R_{g}$. These terms are equivalent through “integration by parts"101010Note that the variations made at the level of the action treat all Green’s functions as some formal Green’s functions, without yet specifying whether they are retarded or advanced. This choice comes either in in-in computations, or, by choosing the retarded boundary conditions at the end of the calculations (i.e. at the level of the equations of motion) [107]. Then, as explained in detail in e.g. Ref. [138], all $\Box^{-1}$ occurrences inside a nonlocal action should be treated formal, i.e. undetermined linear inverses of $\Box$. When the equations of motion are computed, all the $\Box^{-1}$ should be turned into retarded ones by hand. One implication of treating all the $\Box^{-1}$ as equivalent during the variation is that one can effectively integrate $\Box^{-1}$ by parts, as follows:

$\displaystyle\int\mathrm{d}^{D}x\phi(x)\Box^{-1}\psi(x)$ $\displaystyle\equiv\int\mathrm{d}^{D}x\mathrm{d}^{D}y\phi(x)G(x,y)\psi(y)$

$\displaystyle=\int\mathrm{d}^{D}x\mathrm{d}^{D}y\psi(y)G^{T}(y,x)\phi(x)$

$\displaystyle=\int\mathrm{d}^{D}y\psi(y)(\Box^{-1})^{T}\phi(y)$

$\displaystyle\equiv\int\mathrm{d}^{D}y\psi(y)\Box^{-1}\phi(y),$

(2.2)

since the transposed $(\Box^{-1})^{T}$ is also a right-inverse $\Box(\Box^{-1})^{T}={\rm id}$ [138]. Here $D$ is the number of dimensions, $\phi$ and $\psi$ are arbitrary fields, and $G(x,y)$ is the Green’s function appearing in the equations of motion. if the boundary terms are assumed to vanish,111111Looking at the integration procedure given in the previous footnote, we notice that going from the left-hand side of the first line to the right-hand side we have assumed that the homogeneous solution is vanishing, which corresponds to the minimal choice of the boundary conditions as we will discuss later. If we did not discard this homogeneous solution, we would have the homogeneous solution added to the entire $y$-integral, still within the $x$-integral. and we can therefore construct our action with only one of them and without loss of generality.121212Note, however, that this choice introduces a subtlety in the model when dealing with nonzero initial conditions. Unless explicitly stated, we assume vanishing boundary conditions everywhere in this paper, which is consistent with our choice of initial conditions in the cosmological studies of the single-metric model that we introduce later, as well as with the common choices for the DW and MM models. We comment on nonzero initial conditions later. However, in order to keep the structure of the model symmetric, which significantly simplifies the calculations, we include both terms in our action, which now becomes

[TABLE]

where $\alpha$ is a free, dimensionless parameter, to be constrained observationally. In addition, we have omitted the index $g$ in the operator $\Box_{g}$, as well as in the Ricci scalar $R_{g}$, in order to keep the notation simple. From now on, and throughout the paper, all differential operators and curvature quantities with no metric indices are defined with respect to the physical metric $g_{\mu\nu}$; we use a label $f$ when an operator or a quantity is defined with respect to $f_{\mu\nu}$.

As discussed above, in this paper we study only the minimal action (2.3) for nonlocally interacting metrics, and leave the investigation of the more complete set of possible terms, scalar and tensorial, as well as a detailed and rigorous construction of consistent and theoretically sound models for future work. We believe that such properly constructed models should be different and more sophisticated than the simple and phenomenologically constructed model (2.3) studied here.

2.1 Nonlocal equations of motion

Given our model (2.3), the first step is to derive the modified Einstein field equations, i.e. the equations of motion corresponding to $g_{\mu\nu}$ and $f_{\mu\nu}$, by varying the nonlocal action with respect to the two metrics. We obtain

[TABLE]

where $G_{\mu\nu}$ and $G_{\mu\nu}^{f}$ are Einstein tensors corresponding to the physical and reference metrics $g_{\mu\nu}$ and $f_{\mu\nu}$, respectively. $\Delta G_{\mu\nu}$ and $\Delta G_{\mu\nu}^{f}$ are nonlocal distortion terms, i.e. nonlocal corrections to the Einstein tensors, for both metrics, with the forms131313Here $\nabla$ is a covariant derivative, and $(\mu\nu)$ denotes symmetrization over the indices.

[TABLE]

and $T_{\mu\nu}$ is the stress-energy tensor for matter computed in the usual way through the variation of the matter action with respect to $g_{\mu\nu}$. Note that the $f$-metric equations of motion are expectedly not sourced by matter, as the reference metric $f_{\mu\nu}$ does not couple to matter directly.

2.2 Localization

As mentioned in Sec. 1, a powerful technique for dealing with nonlocal equations is to rewrite them in a localized form, by introducing some auxiliary fields. While this provides the possibility of solving and interpreting the equations using regular local methods, one should be cautious that the local versions of the theory are equivalent to the original nonlocal theory only if some conditions are applied to the fields in such a way that the physical degrees of freedom of the theory are kept intact. Otherwise, the “artificial" local fields can behave as “regular" fields which may affect the implications of the theory, both classically and quantum-mechanically, especially since in most cases some of the extra fields (or their combinations) are of ghost behavior. We discuss this issue in Sec. 5.1, and here only introduce the localized formulation of our two-metric model.

In order to do that, let us introduce the two scalar fields $U$ and $V$,

[TABLE]

The action (2.3) can then be written in the local form

[TABLE]

where we have added the two terms $\lambda_{1}(R-\Box U)$ and $\lambda_{2}(R_{f}-\Box V)$ in order to impose the two conditions (2.8) and (2.9), ensuring that the local and nonlocal actions describe the same equations of motion; $\lambda_{1}$ and $\lambda_{2}$ are the corresponding Lagrange multipliers.

First of all, the variation of action (2.2) with respect to the Lagrange multipliers $\lambda_{1}$ and $\lambda_{2}$ expectedly gives Eqs. (2.8) and (2.9). Let us now vary the action with respect to the fields $U$ and $V$. These give, respectively,

[TABLE]

which fix the two Lagrange multipliers $\lambda_{1}$ and $\lambda_{2}$ in terms of the fields $U$ and $V$. Plugging Eqs. (2.11) and (2.12) back into the action yields

[TABLE]

Before deriving the field equations in the local formulation by varying the localized action with respect to $g_{\mu\nu}$ and $f_{\mu\nu}$, we note that the rescaling

[TABLE]

leaves the action, and therefore the equations of motion, invariant. This means that the quantity $M_{\star}\equiv\frac{M_{f}}{M_{\text{Pl}}}$ is redundant and is not a free parameter. We therefore use this freedom to set $M_{\star}=1$ without loss of generality.

The variation of the action (2.2) with respect to $g_{\mu\nu}$ and $f_{\mu\nu}$ then leads to

[TABLE]

for $\Delta G_{\mu\nu}$ and $\Delta G_{\mu\nu}^{f}$, the nonlocal corrections to Einstein tensors for the two metrics, introduced in Eqs. (2.4) and (2.5). Eqs. (2.17) and (2.18) are identical to Eqs. (2.6) and (2.7) when we localize the latter directly at the level of the equations of motion, as expected. This yields another confirmation of the equivalence of our nonlocal and local formulations of the model, as far as the field equations for $g_{\mu\nu}$ and $f_{\mu\nu}$ are concerned (note that here $M_{\star}=1$).

2.3 Bianchi constraints

In addition to the equations of motion, i.e. the modified Einstein field equations for the metrics, we need to know which extra constraints are imposed on the fields when the Bianchi identities are used for $G_{\mu\nu}$ and $G_{\mu\nu}^{f}$, as well as the conservation of the matter energy-momentum tensor $T_{\mu\nu}$. Imposing $\nabla^{\mu}G_{\mu\nu}=\nabla^{\mu}T_{\mu\nu}=0$ for the $g$-metric field equations (2.4), and using the expression (2.6) for $\Delta G_{\mu\nu}$, we obtain141414Note that here we derive the Bianchi constraints in the nonlocal formulation of the model. Imposing the Bianchi identities and the conservation of $T_{\mu\nu}$ on the equations of motion in the local formulation results in exactly the same Bianchi constraints. The procedure is identical, and we do not repeat it here.

[TABLE]

Assuming $\frac{1}{\Box}R\neq 0$, which we need for the nonlocal modification of gravity in our model, this implies

[TABLE]

By performing similar calculations for the $f$-metric field equations (2.5), imposing $\nabla_{f}^{\mu}G_{\mu\nu}^{f}=0$ as well as using the expression (2.7) for $\Delta G_{\mu\nu}^{f}$, we obtain the constraint

[TABLE]

Note that here a subscript or superscript $f$ indicates that the corresponding quantity or operator is defined in terms of $f_{\mu\nu}$. Now, requiring the prefactor $\sqrt{f^{-1}g}(\frac{1}{\Box}R)$ in Eq. (2.21) to be nonvanishing (otherwise it would yield trivial and uninteresting results) implies $\nabla_{f}^{\mu}R_{\mu\nu}^{f}=0$. On the other hand, we have

[TABLE]

Combining the two conditions, we obtain

[TABLE]

Since $R_{f}$ is a scalar quantity, the covariant derivatives $\nabla_{\mu}$ and $\nabla_{\mu}^{f}$ are independent of the metrics, and both conditions (2.20) and (2.23) imply the Bianchi constraint

[TABLE]

This is a surprising result, as the constraint is very strong; let us understand its implications. The constraint (2.24) tells us that the Ricci scalar of the reference metric must be temporally and spatially constant. This means that the form of the reference metric $f_{\mu\nu}$ is highly constrained, as it has to be a metric with specific dynamics, for example, a Minkowski metric for $R_{f}=0$, or a de Sitter metric for a constant and nonzero $R_{f}$. This would in principle be of no problem if the full set of field equations (2.4) and (2.5) could always be satisfied consistently, for all the interesting choices of the physical metric $g_{\mu\nu}$. The problem however is that the Bianchi constraint (2.24) is a condition independent of the physical system under consideration, and independent of the chosen $g_{\mu\nu}$. This then means that not only should the modified Einstein equations (2.4) for $g_{\mu\nu}$ be solved, but also the extra set of Einstein equations (2.5) for $f_{\mu\nu}$ should additionally be satisfied for the chosen $g_{\mu\nu}$ and with an $f_{\mu\nu}$ that has already been fixed to a metric with a constant curvature. This may overconstrain the system, which means that it may not be possible to always solve the set of equations. An explicit example is cosmology. It may be possible to find a solution for $f_{\mu\nu}$ with a constant curvature for the background evolution, i.e. with $g_{\mu\nu}$ having an FLRW form, but the Bianchi constraint (2.24) may not allow us to perturb such a metric when we perturb $g_{\mu\nu}$, as a perturbed $f_{\mu\nu}$ may not allow a constant curvature. In this paper we will show that a consistent solution for $f_{\mu\nu}$ exists for background cosmology, and we leave the nontrivial question of the existence of perturbative solutions for future work. Although we do not know the answer at this stage, it is possible that the $f$-metric equations cannot be satisfied in all interesting cases. This would then mean that Eqs. (2.5) should be discarded, which in turn would mean that we would not be allowed to vary the action (2.3) with respect to the reference metric $f_{\mu\nu}$. In that case, the kinetic term for $f_{\mu\nu}$ would have no effects on any physical quantities, and can also be dropped. The situation then would be very similar to the ghost-free theory of massive gravity, in which the reference metric is not dynamical, as opposed to bimetric theories with both metrics dynamical, and is fixed for the theory independently of the form of the physical metric $g_{\mu\nu}$. There is however a very important difference between our case and massive gravity, and that is the fact that the reference metric here affects the physical sector only through its Ricci scalar, or in other words, the curvature of $f_{\mu\nu}$ enters the $g_{\mu\nu}$ equations only through the Ricci scalar. In massive gravity, the reference metric, being a rank-2 field, is required in order to give mass to graviton, while in our nonlocal model all of its tensorial properties are lost. Clearly, the situation could change if the structure of the model were extended to include other possible scalar terms. The same could be the case if tensorial nonlocalities were (also) considered, for which either the Bianchi constraint would not be as strong as in our simple model, or the theory would become similar to massive gravity with a fixed reference metric affecting the gravity sector as a full, tensor field, and not only through a scalar contribution. Another possibility in that case would be to add the local, ghost-free interactions of massive or bimetric gravity to our simple nonlocal model, which may also violate the strong Bianchi constraint (2.24). These are all interesting and exciting possibilities, but are beyond the scope of the present paper, and we leave their investigation for future work. For now, we take a closer look at our model (2.3) with only the simple scalar interaction terms, and see whether the structure of the model when the Ricci scalar is assumed to be fixed to a constant would suggest an interesting model. This is the subject of the next section.

3 The $m^{2}\frac{1}{\Box}R$ model

Let us look again at action (2.3) for our model of nonlocally interacting metrics, this time taking into account the condition required for the consistency of the solutions, i.e. the Bianchi constraint (2.24), which, as discussed in the previous section, forces the Ricci scalar of the reference metric, $R_{f}$, to be a constant. Calling the constant combination $2\alpha R_{f}$ simply $m^{2}$,151515Recall that $R_{f}$ has dimension $[M^{2}]$, and $\alpha$ is dimensionless. we can now impose the condition (2.24) at the level of the action. As argued above, as long as we are interested only in the dynamics of the physical metric $g_{\mu\nu}$, we can assume that $f_{\mu\nu}$ is a fixed metric, and therefore, no longer vary the action with respect to $f_{\mu\nu}$.161616The reference metric $f_{\mu\nu}$ is, for example, of a de Sitter form, which is determined purely through the constant curvature $R_{f}$, completely independently of $g_{\mu\nu}$ and matter. We can thus fully ignore the $f$-metric kinetic (Einstein-Hilbert) term in the action. The action then reads

[TABLE]

which is a single-metric model with some nonlocal distortion terms added to GR. Before we discuss the implications of this observation, let us quickly obtain the field equations for action (3.1) and find the relations between the quantities in this model and the ones in our original two-metric model. It is important to note, as we discussed in Sec. 2, that the two terms $m^{2}\frac{1}{\Box}R$ and $R\frac{1}{\Box}m^{2}$ are identical through “integration by parts" if the boundary terms vanish. In that case, the action simplifies to

[TABLE]

Looking at the two actions (3.1) and (3.2), we notice one more subtlety, and that is for the zero value of the parameter $m$. In this case, clearly the action (3.2) reduces to the standard action of GR, as expected, while the second nonlocal term in the action (3.1) can remain nonvanishing. The reason is that an equation like $\Box X=0$ can have nonvanishing solutions for $X$, and therefore, the quantity $\frac{1}{\Box}m^{2}$ is not identically zero for $m=0$. However, this difference between the two actions also stems from the subtleties in choosing the initial and boundary conditions for quantities like $\frac{1}{\Box}m^{2}$ in nonlocal models. If the quantity is set to zero initially and at the boundary, it remains vanishing everywhere and at all times, and the two forms of the action become identical. In this paper, we only consider vanishing initial conditions, and it therefore does not matter for our considerations which form of the action to choose. We therefore use the two forms (3.1) and (3.2) interchangeably. We comment on this again in the next section when we discuss vanishing and nonvanishing values of $m$ in our study of the background cosmology for the model.

Let us now derive the equations of motion for the gravity sector of the model. The procedure is similar to the one for the original two-metric model, and, working in the localized formulation of the model, we first introduce the auxiliary field

[TABLE]

Plugging this into the action, and adding a Lagrange multiplier $\lambda$ in order to impose (3.3), we obtain the localized action

[TABLE]

By varying this local action with respect to $g_{\mu\nu}$, we obtain the modified Einstein equations

[TABLE]

which are identical to Eqs. (2.4) in combination with Eqs. (2.17), after performing the transformations

[TABLE]

The modified Einstein equations (3.5) now take the form

[TABLE]

which can be solved together with the extra equations

[TABLE]

which are obtained through varying the action (3.4) with respect to $\lambda$ and $U$, respectively.

We can see explicitly from the structure of Eq. (3.8) for $V$ (or, equivalently, for the Lagrange multiplier $\lambda$) that this model, although involving only one operator $\Box$ in the nonlocal term, still needs two auxiliary fields for localization, in contrast to the DW model studied in Sec. 4.3. The reason is that Eq. (3.8) does not determine $V$ in terms of the other fields that already exist in the model, or the local operators. It involves the nonlocal operator $\Box^{-1}$ acting on the parameter $m^{2}$, and therefore, after plugging $V$ (or $\lambda$) back into the action (3.4), the nonlocalilty remains. We therefore need an extra auxiliary field to take care of this. The fact that one needs two auxiliary fields for localization, even though the nonlocal operator is of the form $\Box^{-1}$, is interesting also in comparison to the MM model with $\Box^{-2}$. This can be understood by noticing the asymmetric structure of the term $m^{2}\frac{1}{\Box}R$, and the fact that $\Box^{-1}$ acts on both sides when the action is varied.

In addition, we should note that although we obtained the model (3.1) through our original two-metric model, the connection to a bimetric setup is now lost, and from a purely phenomenological point of view the model can simply be taken as a simple nonlocal modification of gravity not necessarily related to a model of interacting metrics, with the nonlocality generated by a completely different mechanism. In this respect, for phenomenologists who are not necessarily interested in the fundamental theory behind the model, the action (3.2) can be considered as an standalone, consistent model of modified gravity, and the starting point for any phenomenological studies. The fact that it provides a viable background cosmology, contrary to the similar DW model $\alpha R\frac{1}{\Box}R$ (as we show in the next sections), while being simpler in structure than the MM model $m^{2}R\frac{1}{\Box^{2}}R$, makes the model appealing. It is however important to note that the $m^{2}\frac{1}{\Box}R$ model introduces a new mass scale, similarly to the MM and differently from the DW models. As discussed in Sec. 1, $m^{2}\frac{1}{\Box}R$ is nothing but a model with the structure $m^{2}X$, where $X\equiv\frac{1}{\Box}R$, in comparison to the MM model with the structure $m^{2}X^{2}$. It is quite interesting that such a simple model provides a well-behaved cosmology, at least at the background level, while introducing only one free parameter, just as in $\Lambda$CDM and the MM model. In the next section, we perform a detailed study of the background cosmology of the model, and leave the exploration of possible theoretical origins of the model, not necessarily based on a theory of interacting metrics, for future work.

4 Cosmology and expansion histories

Now that we have the field equations (3.6) for the metric $g_{\mu\nu}$, as well as the the constraints relating the auxiliary fields $U$ and $V$ to the metric and the parameter $m$, i.e. Eqs. (3.7) and (3.8), we can start the investigation of background cosmology in the $m^{2}\frac{1}{\Box}R$ model.

Let us follow the standard recipe for modelling the background dynamics of the Universe, and assume that the Universe is described by an FLRW metric. Specializing to a spatially flat universe and working in cosmic time $t$, we have

[TABLE]

Here, $a(t)$ is the scale factor, being a function of time only. We furthermore take a perfect-fluid form for the matter source, and therefore assume $T^{\mu}_{~{}\nu}=\mathrm{diag}(-\rho,p,p,p)$, with $\rho$ and $p$ being the matter+radiation energy density and pressure, respectively.

4.1 Background equations

Plugging the FLRW expressions for the metric (4.1) into the $(0,0)$ components of the field equations (3.6), we obtain the (modified) Friedmann equation

[TABLE]

where $H\equiv\frac{\dot{a}}{a}$ is the Hubble rate, and an overdot denotes a derivative with respect to cosmic time. In order to solve this equation, we also need the equations for the auxiliary fields $U$ and $V$, i.e. Eqs. (3.7) and (3.8), which now take the forms

[TABLE]

As is common in similar analyses of modified gravity models, let us introduce a derivative with respect to the number of $e$-folds $N\equiv\text{ln}a$, which we denote by a prime. The cosmic-time derivatives $\dot{U}$, $\ddot{U}$, $\dot{V}$, and $\ddot{V}$ can now be written in terms of the derivatives with respect to $N$,

[TABLE]

where we have introduced $v\equiv\alpha h^{2}V$, $u\equiv U$, and $\xi\equiv h^{\prime}/h$, with $h\equiv H/H_{0}$ and $H_{0}$ being the present value of the Hubble rate. Eqs. (4.3) and (4.4) now read

[TABLE]

where we have additionally used $R=6(\xi+2)H^{2}$. The Friedmann equation (4.2) in terms of these new variables takes the form

[TABLE]

where the Universe is assumed to be filled with matter and radiation, with the present density parameters $\Omega_{\text{M}}^{0}$ and $\Omega_{\text{R}}^{0}$, respectively. From this expression we can then read off the effective dark-energylike contribution in our model, with the density parameter

[TABLE]

where $\rho_{\text{NL}}$ is the nonlocal contribution to the total energy density of the Universe $\rho_{\text{tot}}$ at any given time. The evolution of the density parameters $\Omega_{\text{M}}$, $\Omega_{\text{R}}$, and $\Omega_{\text{NL}}=1-\Omega_{\text{M}}-\Omega_{\text{R}}$ can be obtained by solving the continuity equations for matter, radiation, and the nonlocal contribution,

[TABLE]

where $w_{\text{M}}\equiv p_{\text{M}}/\rho_{\text{M}}=0$, $w_{\text{R}}\equiv p_{\text{R}}/\rho_{\text{R}}=1/3$, and $w_{\text{NL}}\equiv p_{\text{NL}}/\rho_{\text{NL}}$ are the equations of state for matter, radiation and the dark-energylike nonlocal contribution, respectively. The density parameters $\Omega_{\text{M}}$, $\Omega_{\text{R}}$, and $\Omega_{\text{NL}}$ then evolve, in terms of $N$, through the equation

[TABLE]

Note that $\xi=-1.5$ and $\xi=-2$ for a universe filled with only matter or radiation, respectively, and we therefore recover the evolution equations for matter- and radiation-domination epochs from the general equation (4.14). The nonlocality equation of state $w_{\text{NL}}$ as a function of $N$ can be obtained using Eq. (4.12) in combination with

[TABLE]

Finally, another important quantity for the background study of any cosmological model is the effective equation of state [9]

[TABLE]

which parametrizes the evolution of the total energy density, and is the key quantity in comparing the model’s predicted background dynamics to observations.

Before we study the implications of these equations for the real Universe, and see whether the model is able to describe the cosmic evolution correctly, let us briefly discuss the curious case of $m=0$. Setting $m$ to zero in the Friedmann equation (4.2), as well as in Eqs. (4.3) and (4.4) for $U$ and $V$, we obtain

[TABLE]

These equations may seem surprising, as $m$ being zero should mean that we are back to standard cosmology, with the standard Friedmann equations, whereas here we seem to have several extra terms remaining in Eq. (4.17). As discussed in Sec. 3, the reason is that we have derived the Friedmann equation using the symmetrized form of the action, i.e. the expression (3.1), in which the second nonlocal term is not necessarily vanishing for a vanishing $m$. However, as we argued before, such a term does vanish if the initial and boundary conditions are set to zero, as in that case the quantity $\frac{1}{\Box}m^{2}$ is always zero. This can be seen here explicitly, as Eq. (4.19) forces $V$ (which is proportional to $\frac{1}{\Box}m^{2}$) to always be zero, if it is initially set to zero. This then in turn forces all the nonstandard terms in Eq. (4.17) to also vanish, and we recover standard cosmology with no cosmological constant term, and therefore with no viable solutions. We can of course use nonzero initial conditions for the auxiliary fields $U$ and $V$, and check whether this $m=0$ case could yield any viable cosmologies, along the lines of the analysis performed in Ref. [126]. The numerical investigation of Eq. (4.17), together with Eqs. (4.18) and (4.19), shows however that the set of equations does not provide an evolution for the Universe consistent with observations. This can be seen qualitatively by looking into the structure of the equations. Note that in Eq. (4.17) all the nonstandard terms are proportional to $V$, $\dot{V}$, or $\ddot{V}$. On the other hand, Eq. (4.19) for $V$ does not have any source, and therefore, whatever initial values we start with for $V$ and its derivatives, the Hubble friction term will reduce $V$ with time, forcing it to vanish eventually. This means that the nonlocal terms in Friedmann equation (4.17) will vanish at late times and, therefore, there will be no dark-energylike behavior in the asymptotic future. One could however think that by appropriately setting the initial values for the auxiliary fields $U$ and $V$, and their derivatives, the nonlocal terms would remain nonvanishing for a sufficiently long period of time over the history of the Universe, with viable behavior. This may seem to work especially because there is a compensating term, $\frac{1}{3}\alpha\dot{V}\dot{U}$, in Eq. (4.17), which might remain nonvanishing, as $U$ is sourced and can compensate for the decaying behavior of $V$. However, a detailed investigation of the dynamical equations shows that $V$ decays very quickly, and the nonlocal terms vanish rapidly, after a short period of time, independently of the initial conditions. We therefore conclude that the model with $m=0$ does not provide a viable cosmological solution, and from now on assume that $m$ is nonzero.

4.2 Cosmic acceleration, and comparison with Maggiore and Mancarella’s $m^{2}R\frac{1}{\Box^{2}}R$ model

Having all the equations needed to study the background solutions for our nonlocal model with a nonzero $m$, here we investigate the implications of the solutions, and explore whether the model can provide a viable cosmological evolution. In particular, we study the possibility of obtaining self-accelerating solutions in the absence of a cosmological constant term.

In order to obtain the cosmic evolution for the model, we need to know how the auxiliary fields $u$ and $v$ evolve with time. For that, we solve Eqs. (4.9) and (4.10), for which we need to know how $\xi$ evolves with time. The evolution of $\xi$ can, on the other hand, be obtained using Eq. (4.11) and its derivative (remember that $\xi=\frac{h^{\prime}}{h}$). This yields the quadratic equation

[TABLE]

for $\xi$. This equation, combined with Eqs. (4.9) and (4.10), gives a closed system of differential equations for $u$, $v$, and $\xi$.

We solve all the equations numerically with the integration initiated at $N=-15$, i.e. well inside the radiation-domination epoch in the early Universe. For the initial conditions we simply set $u_{0}=u^{\prime}_{0}=v_{0}=v^{\prime}_{0}=0$.171717Note that the indices [math] here denote the initial values of the quantities, not their present values. Note that, as discussed before, different initial conditions for the auxiliary fields correspond to different homogeneous solutions of the equations that relate the fields to $R$ and $m^{2}$, defining the operator $\Box^{-1}$ in the nonlocal formulation of the model. The choices we have made here are the simplest ones, commonly made for nonlocal models, and correspond to a vanishing homogeneous solution.181818Let us for example assume that $\phi$ is the homogeneous solution for $U$ in the equation

$U=\frac{1}{\Box}R\iff\Box U=R,$

(4.21)

meaning that $\phi$ satisfies $\Box\phi=0$. It is easy to show that a vanishing homogeneous solution $\phi=0$ corresponds to the minimal choice of the initial conditions that we made, i.e. $U_{0}=U^{\prime}_{0}=0$ in this case. Assuming that the nonlocal effects kicked in at some initial time, and that they were absent before that, terms like $\frac{1}{\Box}R$ are nonvanishing only after that initial time. This means that the homogeneous solution $\phi$ and its derivative are initially zero as they are equal to $U_{0}$ and $U^{\prime}_{0}$, respectively. Since $\Box\phi=0$ does not have a source, if $\phi$ and $\phi^{\prime}$ are initially zero, they will always remain zero. In principle, though, one can relax the constraints on the initial conditions, as done in Ref. [126], and explore the effects of different choices on the dynamics, which corresponds to exploring various local formulations of the model. The evolution we obtain in the present paper corresponds to an attractor solution, and even though we tried a few different initial conditions by hand, we always ended up with the same solution at late times. It however remains to see whether other attractor solutions exist by properly choosing the initial conditions, as has been done in Ref. [126] for the MM model. This requires a proper phase-space analysis of the system, which is beyond the scope of the present paper, and we leave it for future work.

In all the numerical calculations and results presented below, we have set $\Omega_{\text{M}}^{0}\approx 0.31$ and $\Omega_{\text{R}}^{0}\approx 9.2\times 10^{-5}$ for the abundance of matter and radiation at present, consistent with the latest cosmological measurements [139]. One should however note that these values have been obtained for the $\Lambda$CDM model, and one therefore clearly needs to perform a proper statistical analysis of our nonlocal model in order to find the best-fit values for these parameters within the framework of the model. As we will see later, the background evolution of the Universe in this model is close to the $\Lambda$CDM one, and therefore, the chosen values should be close to the values one would obtain through a detailed statistical analysis of the model in comparison to the data. We further set $m^{2}=0.232H_{0}^{2}$, which, as we will see, provides a well-behaved and viable cosmic history.

Our numerical solutions for the auxiliary fields $u$ and $v$ in terms of the number of $e$-folds $N$ are shown in Fig. 1. Fig. 2 then shows the evolution of the normalized Hubble rate $h$ ($\equiv H/H_{0}$) computed through Eq. (4.11). For comparison, we have also shown the evolution of $h$ for the standard $\Lambda$CDM and the MM nonlocal models, where we have assumed the same values for $\Omega_{\text{M}}^{0}$ and $\Omega_{\text{R}}^{0}$ as in our model. We have set $\gamma=0.00891$ for the free parameter of the MM model [118]. The figure shows that $h$ behaves almost identically in the past ($N<0$) for all three models, while behaving significantly differently in the (far) future. In contrast to $\Lambda$CDM, with $h$ becoming almost constant in the future, both MM and our models predict Hubble rates that increase with time, albeit much faster in the MM case. Despite this behavior of $h$ in our model and the fact that it does not become a constant in the future, implying that the Universe does not evolve into a de Sitter phase191919Here, by a de Sitter phase, we mean the Universe being dominated by a nonzero, finite, cosmological constant, with $h$ becoming a finite constant. One could consider cases like ours also a de Sitter phase, with an asymptoticly infinite cosmological constant, and therefore an infinite Hubble rate. as opposed to $\Lambda$CDM, we still obtain an accelerating solution with an effective equation of state $w_{\text{eff}}$ close to $-1$, as we will see later. The reason is the slower increase in $h^{\prime}$ compared to $h$, making $\xi$ and $w_{\text{eff}}$ approach [math] and $-1$, respectively, in the future; see Eq. (4.16). This is similar to what happens in the MM model.

Fig. 3 depicts how the density parameters $\Omega_{\text{M}}$, $\Omega_{\text{R}}$, and $\Omega_{\text{NL}}$ evolve in our nonlocal model. These have been calculated by solving Eq. (4.14) numerically. The figure clearly shows the three epochs in the history of the Universe, i.e. the radiation- and matter-domination eras, as well as the final period in which the nonlocalities dominate. Comparing the same curves with those from $\Lambda$CDM, as well as the MM model, (not shown in the figure) shows an almost exact match between all three models.

Although our model seems to successfully describe the three phases in the expansion history as observations require in terms of the evolution of the density parameters, this is not sufficient for the viability of the model; we further need to study the properties of the energy density at each epoch in terms of the evolution of the effective equation of state $w_{\text{eff}}$, which is shown in Fig. 4. This has been calculated simply through Eq. (4.16) and the evolution of $h$ and $h^{\prime}$. The figure clearly shows that $w_{\text{eff}}$ evolves very closely to its evolution in both $\Lambda$CDM and the MM model. The three epochs of radiation, matter, and dark-energy domination can now be seen with $w_{\text{eff}}$ starting with $w=1/3$ at very early times (radiation domination), decreasing to $w=0$ (matter domination), and then becoming negative at late times (dark-energy domination). Although there are differences in $w_{\text{eff}}$ predicted by the three models at late times, the differences are more pronounced in the future. The asymptotic values of $w_{\text{eff}}$ are the same in all three models ($=-1$), but contrary to $\Lambda$CDM, where $w_{\text{eff}}$ always remains larger than $-1$, both the MM model and ours show phantom behavior in the future, with $w_{\text{eff}}$ decreasing, crossing $-1$, and then increasing again towards $-1$. It is however interesting to notice that our model predicts an effective equation of state that is closer to the $\Lambda$CDM behavior compared to the MM model. In order to understand the source of this observation in an intuitive (but handwaving) way, let us notice that the action of our model includes a nonlocal term of a lower order in $\frac{1}{\Box}R$ compared to the MM one; this makes the effects of the nonlocalities less pronounced. The nonlocal term in our model is of the form $m^{2}X$, with $X\equiv\frac{1}{\Box}R$, while the nonlocality in the MM model is of the $m^{2}X^{2}$ form. $X$ in an FLRW cosmology is a double time-integral of $R$ [118],

[TABLE]

where $t_{*}$ is some initial time at which the nonlocal effects kick in. Since ‘integral’ is a cumulative quantity, $X$ generically increases with time. Considering the cosmological constant in $\Lambda$CDM as a term of zeroth order in $X$, i.e. of the form $m^{2}X^{0}$, it can give a qualitative explanation for why the evolution curves in our model lie between the ones for $\Lambda$CDM and the MM model. One can consider both nonlocal models as cases with some time-dependent and growing cosmological constants, with the growth rate higher in the MM model compared to ours. This suggests that the higher the power in $\frac{1}{\Box}R$ in the action, the larger the deviation from $\Lambda$CDM.

The effective equation of state is the key quantity in constraining the model by background observations, but it is interesting to also investigate the properties of the dark-energy contribution itself through the study of the evolution of its equation of state $w_{\text{DE}}$ and energy density $\rho_{\text{DE}}$. $w_{\text{DE}}$ in our model, i.e. the same quantity as $w_{\text{NL}}$, can be obtained by solving the continuity equation (4.15) using the evolution of $\rho_{NL}$ and its time derivative. Both of these quantities are shown in Fig. 5 for our nonlocal model, as well as for the $\Lambda$CDM and MM models, as functions of the number of $e$-folds $N$. As expected, in both MM model and ours, $\rho_{\text{DE}}$ evolves from zero in the past and during radiation domination to a nonzero value today (for dark-energy domination), in contrast to $\Lambda$CDM with a constant $\rho_{\text{DE}}$. This is because the nonlocal terms behave as integration terms over time, and do not exist in the far past. Another interesting observation is the behavior of $\rho_{\text{DE}}$ in the future. While it remains constant in $\Lambda$CDM, it increases with time in our model, as it does so in the MM model. The rate of increase is however lower in the former, implying again that the future behavior of the Universe is now closer to $\Lambda$CDM. This can be seen also from the evolution of $w_{\text{DE}}$ in Fig. 5, where the nonlocal terms in both nonlocal models have phantomlike behavior over the entire history of the Universe (i.e. with $w_{\text{DE}}<-1$). The deviation from the $\Lambda$CDM value of $w_{\Lambda}=-1$ is however smaller in our nonlocal model compared to MM. In addition, although $w_{\text{DE}}$ changes its behavior around the present time in both models, the change is more dramatic in the MM case.

Using the Chevallier-Polarski-Linder (CPL) parametrization [140, 141],

[TABLE]

where $z$ is redshift, $w_{\text{DE}}$ in our model is best fit by $w_{0}=-1.075$ and $w_{a}=0.045$. This parameterization is however valid only near the present time (in the region $-1<N<0$), and therefore, Eq. (4.23) cannot be used to fit the equation of state at early times or in the future. These values can be compared with the ones for the MM model, i.e. $w_{0}=-1.144$ and $w_{a}=0.084$ [118], which again shows that our nonlocal model gives an expansion history closer to $\Lambda$CDM compared to the MM model.

4.3 Comparison with Deser and Woodard’s $\alpha R\frac{1}{\Box}R$ model

We have seen so far that, similarly to the MM model, our $m^{2}\frac{1}{\Box}R$ nonlocal model provides a viable cosmic expansion history. It is interesting to now take a more detailed look, in terms of the background cosmology, into the other nonlocal model that we mentioned in Sec. 1, i.e. the simplest version of the DW model [106], with $f(\frac{1}{\Box}R)=\frac{1}{\Box}R$, which is very similar to our model in terms of the structure of the nonlocal term. The action for this model is of the form

[TABLE]

where $\alpha$ is a dimensionless, free parameter. Let us first localize the model as usual, by defining the auxiliary field $U$,

[TABLE]

Substituting $U$ into the action and imposing the above condition by introducing a Lagrange multiplier $\lambda$, we obtain the local action

[TABLE]

Varying the action with respect to $\lambda$ gives nothing but the constraint (4.25). By varying the action with respect to the auxiliary field $U$ we obtain

[TABLE]

which can now be used to replace $\lambda$ in the action with a combination of $\alpha$ and $U$. This shows that this model needs only one auxiliary field, $U$, for localization. Finally, we can obtain the modified Einstein equations by varying the action (4.26) with respect to $g_{\mu\nu}$, which gives

[TABLE]

It is easy to also show that the Bianchi identity and conservation of $T_{\mu\nu}$ do not impose any extra constraints on the model and are identically satisfied, as the covariant derivatives of the extra terms in the Einstein equations are identically zero.

Since we are interested in the background cosmology of the model, we assume the FLRW metric

[TABLE]

and obtain the modified Friedmann equation

[TABLE]

with $a$, $h$, $\xi$, $\rho$, $\Omega_{\text{M}}^{0}$, $\Omega_{\text{R}}^{0}$, $N$, and the prime being defined as in previous sections. The constraint (4.25) now reads

[TABLE]

Using this equation, Eq. (4.30) can be written in the simpler form

[TABLE]

We can now use Eq. (4.32) and its time derivative to obtain an expression for $\xi$, which we can then plug into Eq. (4.31) and solve the resulting equation numerically. This gives us the evolution of $\xi$, and therefore the effective equation of state $w_{\text{eff}}=-1-\frac{2}{3}\xi$. As usual, we set the initial conditions $U(N=-15)=U^{\prime}(N=-15)=0$.

Fig. 6 (left panel) shows $w_{\text{eff}}$ as a function of $N$ for the DW model, when we set $\alpha=-0.02285$ (this choice corresponds to the choice made in Ref. [109]).202020Note that the DW model possesses two branches of solutions, one with positive and one with negative values of $\alpha$. The branch corresponding to positive $\alpha$ does not provide a dark-energy-domination phase at late times and in the future, and we therefore do not consider it here. See Ref. [109] for details. For comparison, we have plotted $w_{\text{eff}}$ also for our nonlocal $m^{2}\frac{1}{\Box}R$ model and for $\Lambda$CDM. As it has already been discussed in the literature [109], the simple DW model of $\alpha R\frac{1}{\Box}R$ does not provide a viable background evolution. It does not give a proper matter-domination epoch, and $w_{\text{eff}}$ does not become sufficiently negative at recent times ($w_{\text{eff}}$ obtained for today is $\approx-0.1$). Additionally, the model contains a finite-time singularity in the future, which can be seen from the rapid decrease in $w_{\text{eff}}$, even though it is not necessarily a problem. We have shown in Fig. 6 (right panel) also the evolution of the density parameters for radiation, matter, and the dark-energylike nonlocality, which again shows that the model is not able to produce a proper matter-domination epoch. It is however interesting to note that although the model does not provide a viable cosmic history, the dark-energy-domination epoch occurs at late-times with a choice of the $\alpha$ parameter that is not too small, and therefore not fine-tuned. It has proven difficult to come up with a model of dark energy or modified gravity that is able to provide a viable cosmic acceleration without introducing a new mass scale (normally of the order of $H_{0}$). Therefore, the fact that the DW model provides an evolution, which, although being unviable, is not too far from the observed expansion history, by introducing only one dimensionless parameter with a natural value, encourages searches for viable scale-free models of dark energy.

Finally, we should note that the arguments presented here are valid only for the simple type of DW model, i.e. when the free function $f$ in the full formulation of the model takes the simple form $f(\frac{1}{\Box}R)=\frac{1}{\Box}R$. By tuning the form of the function $f$ it is possible to resolve this problem and reconstruct any arbitrary (and viable) cosmic histories, including $\Lambda$CDM [109, 111], as discussed in Sec. 1. In those cases, however, the model would need a larger number of free parameters, and the form of the function $f$ would become quite contrived, rendering the model less appealing.

4.4 Existence of $f$-metric solutions in the original two-metric framework

Although, as we argued before, our single-metric model (3.1) captures all the phenomenological implications of our original two-metric model, as long as we are interested in only the $g$-metric and matter dynamics, it is interesting to see whether the field equations governing the dynamics of the reference metric $f_{\mu\nu}$ could be satisfied cosmologically if we would still want to keep it dynamical, and vary the action also with respect to the reference metric. A positive (or negative) answer to this question does not affect the validity of our cosmological results, as the effects of the reference metric on physical quantities and their evolution have already been taken into account through the contribution of $f_{\mu\nu}$ to the $g$-metric field equations through $R_{f}$. It is however interesting to see whether the strong constraint imposed on the reference metric, namely $R_{f}$ being a constant, could be satisfied dynamically. We study this question in this section for the background evolution of $f_{\mu\nu}$, and leave the more interesting case of perturbations for future work.

Assuming the reference metric $f_{\mu\nu}$ to be, similarly to $g_{\mu\nu}$, described by an FLRW metric, and again specializing to a spatially flat universe and working in cosmic time $t$, we have

[TABLE]

Here, $X(t)$ and $Y(t)$ are the lapse and scale factor of the reference metric $f_{\mu\nu}$, respectively, both being functions of time only. The choice for the form of the physical metric $g_{\mu\nu}$ is as before, i.e. Eq. (4.1). Because of general covariance, we can still freely choose the cosmic-time coordinate for $g_{\mu\nu}$ and set $g_{00}=-1$. These simple forms for the background metrics significantly simplify the $f$-metric field equations (2.5) together with Eqs. (2.18), for example through the simple forms for quantities like $\sqrt{f^{-1}g}$,

[TABLE]

Let us now consider the implications of the Bianchi constraint (2.24), which for the FLRW form of $f_{\mu\nu}$ implies a constant Ricci curvature,

[TABLE]

with an overdot denoting a derivative with respect to cosmic time. If the Ricci scalar is zero, then $f_{\mu\nu}$ is, e.g., of a Minkowski form,212121Note that $R_{f}=0$ does not imply only a Minkowski form for $f_{\mu\nu}$. For example, an FLRW metric for a universe filled with only radiation also has a zero Ricci curvature. However, since in our model $f_{\mu\nu}$ is not sourced directly by matter, we do not need to consider such dynamical choices for the reference metric in that case, and can simply assume that it is Minkowski. Our cosmological solutions for the physical metric are independent of the actual dynamics of the reference metric, and the only property of the metric that enters our calculations is that $R_{f}$ is zero. We can therefore use a “Minkowski reference metric” and a “reference metric with zero curvature” interchangeably. and if it is constant but nonzero, then the metric is of, e.g, a de Sitter (or anti-de Sitter) form.222222Again, this case can correspond to various types of metrics for $f_{\mu\nu}$. However, what matters for our physical equations of motion is only the Ricci scalar, and not the specific form of the reference metric. Therefore, as long as there are solutions for a de Sitter $f_{\mu\nu}$ consistent with all the equations, we can restrict ourselves to such a form for $f_{\mu\nu}$, and use a “de Sitter reference metric” and a “reference metric with constant and nonzero curvature” interchangeably. We have already shown in the previous sections that $R_{f}=0$ (corresponding to $m^{2}=0$) does not provide a viable background evolution; let us therefore focus on a constant but nonzero $R_{f}$.

Plugging the FLRW expressions for the metrics (4.1) and (4.33) into the $(0,0)$ component of the $f$-metric field equations, we can obtain the Friedmann equation corresponding to the reference metric. We, additionally, explicitly assume a de Sitter form for $f_{\mu\nu}$. For such a maximally symmetric metric with a constant curvature $R_{f}$, we have (in four dimensions)

[TABLE]

for the Riemann, Ricci, and Einstein tensors, respectively. The $(0,0)$ component of the field equations for $f_{\mu\nu}$ then yields

[TABLE]

The Ricci scalar of the $f$-metric, with a de Sitter form, is, on the other hand, given by

[TABLE]

The question we need to answer now is whether Eqs. (4.39) and (4.40) can be simultaneously satisfied, together with the field equations for the physical metric $g_{\mu\nu}$. In order to do this, we plugged the numerical solutions for $u$, $a$, and $h$, found in the previous sections, into the above equations, and solved them for $X$ and $Y$ (keeping $R_{f}$ constant). Our numerical studies show that consistent solutions exist for both $X$ and $Y$, and we therefore conclude that a de Sitter form for $f_{\mu\nu}$ exists for the background dynamics of the model. An important point to note is that, contrary to the $g$-metric equations, which see the $f_{\mu\nu}$ only through the appearance of $R_{f}$ in the interactions, and therefore the quantities $\alpha$ and $R_{f}$ appear in the equations only through the combination $\alpha R_{f}$ (which we denoted collectively by $-\frac{1}{2}m^{2}$), here the two quantities appear separately in the $f$-metric Friedmann equation. The reason is, clearly, that $R_{f}$ appears in both kinetic and interaction terms for $f_{\mu\nu}$, one in combination with $\alpha$ and one independently. This means that, as long as the evolution equations for $f_{\mu\nu}$ are concerned, the model possesses two independent parameters $\alpha$ and $R_{f}$. In our numerical searches of solutions for the background $f_{\mu\nu}$, we fixed $R_{f}$ and $\alpha$ to $\frac{1}{2}m^{2}$ and $1$, respectively. Whether or not a solution exists also at the perturbative level is beyond the scope of the present paper, and we leave its investigation for future work.

5 Auxiliary fields and the problem of ghosts

In this section, we discuss the problem of apparent ghosts in the local formulation of our $m^{2}\frac{1}{\Box}R$ model, which is a generic feature of all nonlocal models, including the DW and MM, as discussed in Sec. 1. We first show that the local formulation contains a ghost, and then argue that the ghost is not harmful to the theory, for similar reasons as in other nonlocal models. We prove, by analyzing the model in its nonlocal formulation, that the physical degrees of freedom of the theory are the same as in GR, and furthermore, that they are not affected by the local (auxiliary) ghost, and remain healthy.

5.1 Ostrogradski ghosts in the local formulation

Let us rewrite the action (3.4) as

[TABLE]

which is formulated in the Jordan frame, as the gravity sector is modified while matter is minimally coupled to gravity. Let us now change the frame to Einstein through the transformations

[TABLE]

by introducing

[TABLE]

Substituting all these into the action, it takes the form

[TABLE]

which represents the action in the Einstein frame, with matter now coupled to both the metric $g_{\mu\nu}$ and the scalar field $\Omega$. In order to write the action in a canonical form for scalar-tensor theories, we introduce the new fields $\phi$ and $\psi$,

[TABLE]

Discarding the boundary terms, the action can now be written in terms of $\phi$ and $\psi$,

[TABLE]

By looking at the kinetic matrix for the scalar fields $\phi$ and $\psi$,

[TABLE]

we notice that its determinant is always negative, meaning that the matrix is negative definite, signalling the presence of an Ostrogradski ghost.

Let us remind ourselves of the fatal consequences of ghosts (see, e.g., Refs. [142, 143]). The resulting unboundedness of the Hamiltonian of the system from below can make the classical theory fully unstable. Although such instabilities could be acceptable if the unstable modes do not grow too rapidly at both the background and perturbative levels, keeping the theory consistent with observations, the ghosts are definitely fatal at the quantum level. They instantaneously reach states with arbitrarily large negative energies, and therefore decay into matter particles very quickly, filling the Universe with an unacceptably large amount of particles. The theory, therefore, does not possess a stable and well-defined vacuum, and should be rejected. We should note however that there are ways to avoid such a disastrous scenario by modifying the decay rate of the ghost fields through the violation of Lorentz invariance above some energy scale where new physics appears, making the decay time larger than the age of the Universe [144, 145, 146, 147]. Since we do not violate Lorentz invariance in our model, the appearance of the ghost in the local formulation may seem fatal, rendering the model excluded. As we discussed in Sec. 1, and detail in the next subsection, the ghost in our model is only an auxiliary field and not a physical degree of freedom. It does not affect the healthiness of the theory as long as we keep in mind that the localized theory must be equivalent to the original nonlocal one by imposing appropriate initial conditions on the auxiliary fields. In that case, they do not affect the physical degrees of freedom and the theory remains viable.

5.2 Nonlocal formulation and the number of physical degrees of freedom

Let us now analyze the $m^{2}\frac{1}{\Box}R$ model in its original formulation and without localization. The Einstein field equations for the model are identical to those corresponding to the physical metric $g_{\mu\nu}$ presented in Sec. 2.1 for the two-metric model, when the quantity $-2\alpha R_{f}$ is replaced by $m^{2}$. The field equations then read

[TABLE]

where $G_{\mu\nu}$ is the Einstein tensor, and $\Delta G_{\mu\nu}$ is the nonlocal distortion term, with the form

[TABLE]

Our goal here is to count the number of physical degrees of freedom in the model, and prove that they are all healthy. Similar discussions for the Deser-Woodard and Maggiore-Mancarella models have been presented in Refs. [108] and [119], respectively.

Let us follow the procedures of Refs. [108] and [119], and choose the synchronous gauge to write $g_{\mu\nu}$ as

[TABLE]

In general relativity, the dynamical equations of motion for $h_{ij}$, i.e. those that are second-order in time derivatives, are the $(i,j)$ components of the Einstein field equations, and $G_{ij}$ reads

[TABLE]

The $(\mu,0)$ components contain at most first-order time derivatives ${\cal O}(\partial_{t})$, and hence are constraints on the initial data of the metric, i.e. $h_{ij}(t_{0}),\dot{h}_{ij}(t_{0})$, the number of which determines the number of propagating degrees of freedom of the theory.232323For a more rigorous explanation, see, e.g., footnote 12 of Ref. [119]. Considering $\Delta G_{\mu 0}$ we now argue that it vanishes when evaluated at the hyper-surface of the initial conditions.

In $d+1$-dimension spacetimes, the action of the d’Alambertian operator is specified by the Green’s function

[TABLE]

where $x$ and $y$ are two spacetime points. Using this the solution to an inhomogeneous equation $\Box_{x}f(x)=F(x)$ for a function $f(x)$ sourced by another function $F(x)$ can be found by

[TABLE]

Therefore, all we need to compute is the Green’s function $G(x,y)$ from Eq. (5.14) and then solve the integral (5.15) in order to compute $f(x)$. Eq. (5.14) has two solutions, one causal (retarded) and the other one acausal (advanced). In order to maintain causality, we need to make sure that only the retarded Green’s function is used, by assuming242424It has however been argued, e.g. in Ref. [120], that by interpreting the nonlocal action as a quantum effective one, causality is automatically satisfied.

[TABLE]

This requirement, however, is not sufficient to fully specify the Green’s function, because the defining equation is a second-order differential equation and therefore one needs to specify the initial conditions for it. Assuming that our nonlocal model is valid only below some energy scale as an effective field theory, and therefore after some initial time $t_{0}$, we therefore further assume

[TABLE]

which mean that the nonlocality effects begin at $t_{0}$ and are absent before that. The immediate implication of these conditions is that, taking into account Eq. (5.15), the quantities $\Box^{-1}F$ and their first-order time derivatives vanish at the initial time $t_{0}$. Now Eqs. (5.11) for $\Delta G_{\mu 0}$ tell us that the nonlocal distortion terms in the $(\mu,0)$ Einstein equations vanish initially except, potentially, for the term

[TABLE]

Here, the piece that can potentially contain second-order time derivatives is

[TABLE]

and the rest of the term includes only spatial derivatives or one time derivative. Working again in the synchronous gauge as given in Eq. (5.12), we then have $g_{00}g^{00}=1$, and therefore, (5.20) vanishes, implying that (5.19) does not contain any second-order time derivatives, and $\Delta G_{\mu 0}$ then vanishes at the initial time $t=t_{0}$. We can therefore conclude that our nonlocal model possesses the same number of propagating degrees of freedom as in GR (similarly to the DW and MM nonlocal models).

Our arguments above have been based on the single-metric formulation of our nonlocal model, where the dynamics of the reference metric are fully ignored, for the reasons discussed in the previous sections. It is though interesting to also count the number of physical degrees of freedom for the $f_{\mu\nu}$ sector of the theory and see whether they would be affected by nonlocalities.

It turns out, however, that the situation is more subtle in this case. Let us perform a similar study of the $(\mu,0)$ components of the Einstein field equations for the reference metric. Repeating the same procedure as for $g_{\mu\nu}$, it is easy to show that the potentially dangerous terms are now

[TABLE]

which we need to show to vanish initially. Similarly to the $g_{\mu\nu}$ case, the piece that can potentially include second-order time derivatives is

[TABLE]

In the two-metric theory, we only have one set of diffeomorphism invariance and one set of transformations which act on the metric like a diffeomorphism but do not transform the coordinates. The situation is similar to bimetric theories, where one does not have two independent sets of diffeomorphism invariance. This means that the only set of diffeomorphism invariance can be used to take one of the metrics (say, $g_{\mu\nu}$) to a synchronous form, but we do not have another set of diffeomorphisms to do the same for the other metric ($f_{\mu\nu}$ here). This in turn means that we cannot use the same arguments as for $g_{\mu\nu}$ for counting the number of degrees of freedom for $f_{\mu\nu}$ in general and for arbitrary solutions. We can however restrict ourselves to the specific cosmological solutions and the specific form of the reference metric that we have considered in this work. We have seen that the consistency conditions for the solutions of the two-metric model force $f_{\mu\nu}$ to be a metric with constant and nonzero Ricci scalar, nonperturbatively. Even though it is possible for such a metric to, in general, not satisfy $f_{00}f^{00}=1$, we have shown that there is at least one solution for $f_{\mu\nu}$ that possesses this property, and that is the de Sitter reference metric studied here, for which $f_{00}f^{00}=1$. This therefore shows that (5.22) vanishes in our case, and, consequently, $\Delta G_{\mu 0}^{f}$ also vanishes at $t=t_{0}$, similarly to the $g_{\mu\nu}$ case. This means that the number of propagating degrees of freedom for the $f$ metric around the cosmological solution considered here is the same as in a GR-like theory for the reference metric. In addition, given that $f_{\mu\nu}$ is a fixed, unphysical262626Note that here by $f_{\mu\nu}$ being “fixed” we simply mean that we can choose it to be of any arbitrary form, with a constant and nonzero $R_{f}$, for any physical system that we are interested in, as long as we can satisfy its own equations of motion. Clearly, the Bianchi constraints, forcing a constant Ricci scalar for the reference metric, do not fix its form to a specific metric (see Footnote 22), and we always have the freedom to choose any form for it, either universally and independently of the particular physical system under investigation, or differently in different cases, as long as it has a constant and fixed $R_{f}$ in all cases. Here, therefore, by “fixed” we simply mean that we “choose” the form of $f_{\mu\nu}$ before working with the model, for example for cosmology, and then fix it to that form for the entire analysis. Also note that by calling $f_{\mu\nu}$ “unphysical” we simply mean that it is decoupled from matter and is unobservable, for the reasons explained in the text. Although, strictly speaking, “unphysical” is not a correct word, we adhere to it in this paper as it is commonly used in the field of multi-metric gravity for describing reference metrics. metric which does not couple to matter, we conclude that our original nonlocally-interacting-metric model possesses the same number of “physical" degrees of freedom as in GR (corresponding to the physical metric $g_{\mu\nu}$) around the cosmological backgrounds studied here. Counting the number of total propagating degrees of freedom (including those of the $f$ metric) around more general solutions requires a more detailed and careful analysis, and is beyond the scope of the present paper.

Let us however point out that even though there seem to be only $2+2$ degrees of freedom in the linear spectrum of the two-metric model for the cosmological solutions studied here, it is likely that more degrees of freedom are present in the nonlinear theory due to the usual no-go theorems for constructing ghost-free theories of interacting massless gravitons (see the discussion in the introduction for possibilities of evading the no-go theorems when locality is violated). This will then imply that around the cosmological backgrounds, these modes are infinitely strongly coupled since their kinetic term does not show up at the linear level where the number of dynamical degrees of freedom is reduced. In general, this renders the use of linear perturbation theory invalid. If further detailed studies of the the model shows that in fact the full theory either has ghosts or is infinitely strongly coupled around cosmological backgrounds, then this will serve as a major drawback to the two-metic, nonlocal model studied in this paper. These are very interesting and important questions to address and need to be investigated in a dedicated and separate work, which we leave for future. Note however that the phenomenological, simple, single-metric model of $m^{2}\frac{1}{\Box}R$, which all our cosmological studies in the present paper have been based on, is immune from these problems. This model can have a completely different origin, fully unrelated to interacting spin-2 fields.

Up until now, we have only shown that the nonlocal terms do not add additional propagating degrees of freedom to the theory for our physical modes corresponding to the physical metric $g_{\mu\nu}$ in general, and in our specific cosmological setup for the reference metric $f_{\mu\nu}$, if certain initial conditions are imposed. However, there is still the possibility that nonlocalities turn the graviton modes into ghosts, and therefore render the theory unstable. It is therefore important to study the propagation of graviton modes and ensure that they always stay healthy. Although a full and rigorous (perturbative) stability analysis of the solutions is necessary for our model, similar to any other modifications to GR, it goes beyond the scope of the present paper, and we leave it for future work. We, however, again follow Refs. [108] and [119], and study the question of whether there could be any ghost instabilities from the change in the sign of the kinetic terms for gravitons, coming from nonlocalities.

Let us again analyze the $m^{2}\frac{1}{\Box}R$ model first. Graviton modes can in principle become ghosts if the $(i,j)$ components of the Einstein equations are allowed to contain second-order time derivatives of the metric with an overall negative sign. In order to explicitly check whether this can happen, we first note that the Ricci scalar $R$ contains at most first-order temporal derivatives on the spatial components of the metric.272727In principle, $R$ contains a term $-\Box\log h$, with $h$ being the determinant of $h_{ij}$ in Eq. (5.12), but since $h$ is already constrained by the $(0,0)$ component of Einstein equations, this term does not involve any unconstrained second-order time derivative terms (see Ref. [119] for more discussions). In addition, in an expression that contains $\frac{1}{\Box}R$, $R$ is integrated over twice temporally, leaving no time derivatives in the expression. Hence, in order for terms including $\frac{1}{\Box}R$ to contain second-order time derivatives, they need to involve two time derivatives acting on $\frac{1}{\Box}R$. Let us now look at the $(i,j)$ components of the Einstein equations,

[TABLE]

Since in Eqs. (5.2) there are no higher than first-order time derivatives acting on $\frac{1}{\Box}R$, the only second-order time derivatives acting on the spatial components of $g_{\mu\nu}$ come from the Einstein tensor $G_{ij}$. We therefore need to check only the sign of the term multiplied by $G_{ij}$ in order to see whether it can be altered by nonlocal effects. Specifically, we need

[TABLE]

in order to ensure the absence of ghostly gravitons. Let us now remind ourselves that the mass parameter $m$ should be of the order of $H_{0}$ for the model to provide a viable cosmic evolution. In addition, similarly to the other nonlocal models of DW and MM, we do not expect our model to exhibit a vDVZ discontinuity [148, 149], meaning that no screening mechanism is needed for the model to reduce to GR on solar system scales if the mass scale of the theory is as low as the Hubble rate today. The term $\frac{1}{\Box}m^{2}$ in the condition (5.24) is therefore completely negligible on small scales. This all means that any potential danger for gravitons to turn into ghosts can only be on cosmological scales. Now, in order to see whether the condition (5.24) is satisfied cosmologically, we use its equivalent expression in terms of the auxiliary field $v$ introduced in Sec. 4.1, namely282828Here, we only consider the background dynamics, and ignore the effects of perturbations, as long as they remain small and stable. A detailed investigation of such effects is beyond the scope of this paper, and we leave it for future work.

[TABLE]

The evolution of $1-2\frac{v}{h^{2}}$ is presented in Fig. 7 (left panel). We clearly see that it is always positive, and, therefore, the condition (5.25) is always satisfied.292929Note that the behavior of $v$ at late times, depicted in Fig. 1, is independent of the initial conditions for the auxiliary fields. It is because the solution is an attractor, and even by setting the initial conditions differently, $v$ moves very rapidly to its negative values. It is however possible for some choices of initial conditions to violate the condition (5.25) at early times, and one should therefore be careful with such choices.

For completeness, let us end this section by looking also at the ghost instabilities for the $f$-metric gravitons, in the original two-metric model and in our specific cosmological setup. We know that $f_{\mu\nu}$ has fixed dynamics, as it is forced to be of, e.g., a de Sitter form. In fact, $f_{\mu\nu}$ only introduces a scale in the model, $R_{f}$, which is manifest in the single-metric formulation of the model, i.e. $m^{2}\frac{1}{\Box}R$. The fact that $f_{\mu\nu}$ can always be fixed to a metric like de Sitter (at least at the background level) is sufficient to guarantee that $f$-gravitons are healthy. In addition, the consequence of $f$ having a fixed $R_{f}$, namely that our model can be fully formulated in terms of only one metric, implies, on its own, that the $f$-gravitons do not affect the properties of the physical $g$-gravitons and matter. For the completeness of our discussions, however, we repeat an analysis similar to that of $g_{\mu\nu}$ for the $f_{\mu\nu}$ Einstein equations in order to show explicitly that the conditions for the healthiness of $f$-gravitons are satisfied.

Let us therefore look at the $(i,j)$ components of Eq. (2.7). We note that the only second-order time derivatives acting on $f_{\mu\nu}$ come from either the second-order derivatives of the determinant of $f_{\mu\nu}$, or the Einstein and Riemann tensors. The determinant is constrained through the $(0,0)$ component of the Einstein equations, and its second-order derivatives do not affect the propagation of the $f$-gravitons. Remembering again that $R_{f}$ is a constant, the potentially dangerous term is therefore (setting $M_{f}=M_{\text{Pl}}$)

[TABLE]

The condition for the kinetic sign of the $f$-gravitons to stay positive is therefore

[TABLE]

with $u$ being defined in Sec. 4.1. The cosmological evolution of $1-2\alpha u\sqrt{f^{-1}g}$ given in Fig. 1 (right panel), corresponding to the solution found for a de Sitter $f_{\mu\nu}$ discussed in Sec. 4.4, demonstrates that the quantity is always positive,303030Note again that this has been shown only for the specific, de Sitter solution of $f_{\mu\nu}$ considered in this paper, and does not necessarily apply to all possible solutions. It therefore demonstrates the absence of ghosts only around a certain class of backgrounds. and the condition (5.27) is always satisfied.313131Similarly to the comments in Footnote 29, here also one should be careful with choosing the initial conditions for the auxiliary fields in order to guarantee the condition (5.27). The fact that $f_{\mu\nu}$ has a fixed and healthy dynamical form, e.g. de Sitter, requires an extra care to ensure that the initial conditions are compatible with the healthiness of $f$-gravitons.

6 Conclusions and outlook

We began the paper with the question of whether two or more spin-2 fields could interact nonlocally, resulting in new models of modified gravity in the infrared, avoiding the cosmological problems present in locally interacting theories, and providing an explanation for the late-time acceleration of the Universe. Aiming at constructing models which would serve as a bridge between two classes of modifications to general relativity, namely, multimetric and nonlocal theories, we assumed nonlocal interactions between the metric of spacetime and an extra spin-2 field with no direct interaction with matter. We focused only on cases with scalar curvatures in the action, and further restricted our studies to a specific, minimal subset of possible forms for the interaction between the two metrics, which is, arguably, the simplest possible choice, inspired by the nonlocal model proposed by Deser and Woodard. We derived the field equations for both metrics, as well as the constraints imposed on the model by applying Bianchi identities and the conservation of matter energy-momentum tensor, and found that the latter would place a very strong constraint on the form of the reference metric, forcing it to have a (spatially and temporally) constant Ricci scalar, independent of the form of the physical metric. This might imply that such a simple nonlocality, with only the Ricci scalar present in the interaction, would not provide a consistent bimetric theory with a dynamical reference metric with respect to which the action is varied. We however showed in this paper that such a consistent dynamical metric does exist at the cosmological background level, and it is still an open question whether such a metric can exist also at the level of perturbations. In case the answer to this question turns out to be negative, one could then fix the reference metric to a specific form, with no kinetic term in the action and without varying the action with respect to it, similarly to the nondynamical reference metric in theories of massive gravity. However, the fact that the Bianchi constraint forces the reference metric to have a constant Ricci scalar, which is the only quantity appearing in the nonlocal interaction terms, makes the introduction of a tensorial field of little motivation, as the scalar quantity does not represent the tensorial structure of such a field. We expect, however, that the entire story will change if one includes tensorial interactions in the theory, i.e. terms involving the Ricci and Riemann tensors, $R_{\mu\nu}$ and $R_{\mu\nu\alpha\beta}$. In addition, one could add the ghost-free local interaction terms to the theory. In either case the strong Bianchi constraint is expected to be violated, and therefore the theory could possess nontrivial and consistent implications, with potentially very rich phenomenologies. It is also interesting to study possible ways of obtaining nonlocally interacting multimetric theories from underlying local theories, for example by properly integrating out light fields interacting with spin-2 fields. This would provide models with nontrivial and interesting phenomenologies. These are all exciting routes to explore, and we leave them for future work.

Inspired by our first and simple attempt at building models of nonlocally interacting metrics, we have proposed in this paper a new, simple, single-parameter model of modified gravity on cosmological scales that successfully provides a mechanism for the late-time accelerated expansion of the Universe, without an explicit cosmological constant. The model has a simple structure, similar to the simplest version of the Deser and Woodard’s nonlocal model. We have derived the equations of motion, and the model has then been localized by introducing two auxiliary scalar fields, providing a framework in which the equations of motion can be handled more easily, and potentially solved for cosmology and any other systems of interest. We have then studied the cosmology of the model at the background level, investigating whether viable expansion histories could be achieved. Our detailed studies have shown that one can indeed obtain a cosmic evolution consistent with the observed one. In particular, we have demonstrated that the model provides all three epochs of radiation, matter and dark energy domination, in agreement with observations. The effective equation of state has also been calculated, and shown to be consistent with its observed evolution. The cosmic history in our model closely resembles that of $\Lambda$CDM over the entire evolution (even more so than the nonlocal model of Maggiore and Mancarella), and the effective equation of state approaches $-1$, although the Hubble rate increases in the future in contrast to the future de Sitter phase of $\Lambda$CDM with the Hubble rate becoming constant. The effective equation of state today is slightly smaller than $-1$, and first decreases for about one $e$-fold in the future and then continuously increases towards $-1$. Even though the present dark-energy equation of state is very close to the $\Lambda$CDM value, with the CPL parameters $w_{0}\sim-1.075$ and $w_{a}\sim 0.045$, the difference between the two is sufficiently large (about $7.5\%$ for $w_{0}$) that the model can observationally be distinguished from $\Lambda$CDM using the forthcoming cosmological surveys. We have also compared the cosmic background evolution to that of the simplest form of the model of Deser and Woodard with no viable cosmic history, and have shown explicitly how the two are different. In order to know whether the model proposed in this paper is consistent with all the observational data, one needs to study the formation and evolution of cosmic structure for the model. This requires a full perturbative analysis, and we leave it for future work. In addition, we have restricted ourselves in this paper to zero initial conditions for the auxiliary fields in the local formulation of the model, and it is interesting to know whether new types of viable solutions exist if the initial conditions are set differently. This would require a full phase-space analysis of the model, and remains to be investigated.

We have finally discussed the question of apparent ghosts in the local formulation of our model, where one of the two auxiliary scalar fields is always a ghost. We have demonstrated that the number of physical degrees of freedom in our single-metric, nonlocal model is the same as in general relativity, when proper initial conditions necessary for the completeness of the model are imposed. This shows that the scalar fields appearing in the local formulation are auxiliary and unphysical, and therefore the ghosts are harmless to the theory. Additionally, we have shown that the physical degrees of freedom are not affected by the ghosts, and always stay healthy. The issue of ghosts has also been discussed for some generalizations of the single-metric model, where we have shown that there could be cases, with a particular ghost-free condition satisfied, that possess no ghosts in their local formulations.

We have also discussed the issue of ghosts and the number of degrees of freedom in the original, nonlocaly-interacting-metric model, and argued that this two-metric model contains 2+2 degrees of freedom linearly and around the cosmological solutions that we have studied in this paper. We have emphasized that the full, nonlinear model may possess more degrees of freedom if the known no-go theorems in constructing ghost-free theories of interacting, massless, spin-2 fields hold true for these nonlocal theories as well, which will then imply that the theory has ghosts or is infinitely strongly coupled around cosmological backgrounds, invalidating the application of linear cosmological perturbation theory. This may be considered as a major issue for the model, implying that this nonlocally-interacting-metric model as a full theory might be either inconsistent or not applicable to cosmology in its present form. We consider this as another motivation for extending or modifying the model. From the theoretical point of view, these are pressing issues that have not been addressed in the present paper, and the model has therefore not yet been shown to be consistent. We leave an investigation of these questions to future work, where a fully theoretical study of the model and its consistency will be performed. We however emphasize again that the single-metric, nonlocal model whose cosmological implications have been studied in this paper is immune from these issues.

Even though we have been led to our phenomenological, nonlocal, single-metric model through an investigation of nonlocally interacting metrics, it can be used as a standalone, nonlocal model for all phenomenological studies, with no reference to any two-metric scenarios. It could be considered as an effective model originating from a completely different, more fundamental, underlying theory. We leave the investigation of such possibilities for future work. It is interesting that even though we tried to construct our model through two metrics interacting nonlocally, the consistency constraints on the solutions of the simplest, scalar-based models were so stringent that the final model became effectively single-metric with a simple structure. Whether the model could be obtained in a completely different way, and whether other models with similarly simple structures exist, are currently two open questions that need to be investigated. It is still an open question whether a more sophisticated model of modified gravity can be constructed by letting the physical metric interact with another spin-2 field nonlocally, in an attempt to obtain a viable cosmic history with self-acceleration at late times and in the future. In particular, such models have the very important and interesting potential of providing viable, scale-free models of modified gravity with no fine-tuning of parameters. All theories of multimetric gravity proposed so far have been constructed through local interactions, and most of them suffer from problems of ghost and gradient instabilities. We believe that relaxing the assumption of locality could open up new directions for extending the class of such theories, with interesting cosmological implications. It also remains to be seen, through rigorous theoretical investigations, whether (and if so, why) nonlocal theories can evade the well-known no-go theorems on the construction of consistent theories of interacting (massive or massless) spin-2 fields.

Acknowledgments

We thank Ana Achúcarro, Nima Khosravi, Tomi S. Koivisto, Alexey S. Koshelev, Henrik Nersisyan, Javier Rubio, Angnis Schmidt-May, Adam R. Solomon, and Nico Wintergerst for helpful discussions. We are also grateful to the anonymous referee for thoroughly and critically reading the manuscript and for invaluable comments. V.V. is supported by a de Sitter fellowship of the Netherlands Organization for Scientific Research (NWO). Y.A. and A.S. acknowledge support from the NWO and the Dutch Ministry of Education, Culture and Science (OCW), and also from the D-ITP consortium, a program of the NWO that is funded by the OCW. L.A. acknowledges support from DFG through the project TRR33 “The Dark Universe.”

Appendix A Ostrogradski ghosts in generalizations to $m^{2}G(\frac{1}{\Box}R)+RF(\frac{1}{\Box}m^{2})$

nonlocalities

In this appendix we briefly study a generalization of our $m^{2}\frac{1}{\Box}R$ model to cases that resemble the DW model with a general form of the function $f(\frac{1}{\Box}R)$. Let us therefore consider interaction terms of the form

[TABLE]

where $F$ and $G$ can be any arbitrary functions. As usual, defining the auxiliary fields

[TABLE]

we can localize the action for the generalized model, which then takes the form

[TABLE]

For analytical functions $G(U)$ and $F(V)$, such as exponential or polynomial, with well-defined Taylor expansions, and using integration by parts, we can show $m^{2}G(U)=RG(V)$ and $RF(V)=m^{2}F(U)$ at the level of the action. This means that we can rewrite the terms involving the functions $G$ and $F$ in the action (A) in terms of only $U$ or $V$. Without loss of generality, and by redefining the functions $G$ and $F$, we can then rewrite the model in such a way that $G=F$.

We can now vary the action with respect to the fields appearing in action (A). The variation with respect to the Lagrange multipliers $\lambda_{1}$ and $\lambda_{2}$ gives the constraints (A.2) and (A.3), while the variation with respect to the auxiliary fields $U$ and $V$ gives

[TABLE]

respectively, where a prime here denotes a derivative with respect to the argument of the function. In order to solve these equations for $\lambda_{1}$ and $\lambda_{2}$, we introduce two new auxiliary fields $\tilde{U}$ and $\tilde{V}$,

[TABLE]

In principle, we need to introduce two new Lagrange constraints $\tilde{\lambda}_{1}(m^{2}F^{\prime}(U)-\Box\tilde{U})+\tilde{\lambda}_{2}(RF^{\prime}(V)-\Box\tilde{V})$. However, by doing that and performing the variation of the action with respect to $\tilde{U}$ and $\tilde{V}$ we get $\tilde{\lambda}_{1}=\tilde{\lambda}_{2}=0$ for the Lagrange multipliers $\tilde{\lambda}_{1}$ and $\tilde{\lambda}_{2}$. This implies that $\lambda_{1}$ and $\lambda_{2}$ in the action (A) are themselves the two extra auxiliary fields that we need in order to localize the action. For the simplicity of notation, we keep $\tilde{U}$ and $\tilde{V}$ for these extra fields, instead of $\lambda_{1}$ and $\lambda_{2}$. We therefore need four auxiliary fields for localizing our generalized, nonlocal model, instead of two for the simple choice of $F(U)=U$ studied in the paper. This is consistent with what happens in the DW model with general forms for $f(\frac{1}{\Box}R)$, where one needs two auxiliary fields for localization, in contrast to the simple case of $f(\frac{1}{\Box}R)=\frac{1}{\Box}R$.

The action (A) now reads

[TABLE]

Let us now repeat the procedure of Sec. 5.1, and study the apparent ghosts in the localized formulation of our generalized model. We therefore use the transformations

[TABLE]

by introducing

[TABLE]

to write the action in the Einstein frame, which now takes the form

[TABLE]

We now introduce the new field $\phi$,

[TABLE]

and rename the other auxiliary fields as

[TABLE]

Discarding the boundary terms, the action becomes

[TABLE]

In order for the kinetic matrix of the auxiliary fields $\{\phi,\psi,\chi,\zeta\}$ to be positive definite, we require

[TABLE]

which then requires

[TABLE]

This means that, depending on the form of $F$ and the dynamics of $\phi$ and $\psi$, it is possible for the generalized model to be free of ghosts in the local formulation, in contrast to the simple model with $F(\psi)=\psi$ which always has a ghost.