The European Physical Journal C
January 2020 , 80:32 | Cite as
Anisotropic massive Brans–Dicke gravity extension of the standard \(\Lambda \)CDM model
Özgür Akarsu
Nihan Katırcı
Neşe Özdemir
J. Alberto Vázquez
Regular Article - Theoretical Physics
We present an explicit detailed theoretical and observational investigation of an anisotropic massive Brans–Dicke (BD) gravity extension of the standard \(\Lambda \)CDM model, wherein the extension is characterized by two additional degrees of freedom; the BD parameter, \(\omega \), and the present day density parameter corresponding to the shear scalar, \(\Omega _{\sigma ^2,0}\). The BD parameter, determining the deviation from general relativity (GR), by alone characterizes both the dynamics of the effective dark energy (DE) and the redshift dependence of the shear scalar. These two affect each other depending on \(\omega \), namely, the shear scalar contributes to the dynamics of the effective DE, and its anisotropic stress – which does not exist in scalar field models of DE within GR – controls the dynamics of the shear scalar deviating from the usual \(\propto (1+z)^6\) form in GR. We mainly confine the current work to non-negative \(\omega \) values as it is the right sign – theoretically and observationally – for investigating the model as a correction to the \(\Lambda \)CDM. By considering the current cosmological observations, we find that \(\omega \gtrsim 250\), \(\Omega _{\sigma ^2,0}\lesssim 10^{-23}\) and the contribution of the anisotropy of the effective DE to this value is insignificant. We conclude that the simplest anisotropic massive BD gravity extension of the standard \(\Lambda \)CDM model exhibits no significant deviations from it all the way to the Big Bang Nucleosynthesis. We also point out the interesting features of the model in the case of negative \(\omega \) values; for instance, the constraints on \(\Omega _{\sigma ^2,0}\) could be relaxed considerably, the values of \(\omega \sim -1\) (relevant to string theories) predict dramatically different dynamics for the expansion anisotropy.
The base \(\Lambda \)-cold dark matter (\(\Lambda \)CDM) model, relying on the inflationary paradigm [1, 2, 3, 4, 5, 6, 7, 8], is the simplest and most successful cosmological model to describe the dynamics and the large scale structure in agreement with the most of the currently available observational data [9, 10, 11, 12]. It is today credited as the standard cosmological model, yet it is probably not where the story has concluded but the hardest part has just begun. It suffers from severe theoretical issues relating to the cosmological constant \(\Lambda \) being responsible for the late time acceleration of the Universe [13, 14, 15, 16, 17, 18, 19] and, on the observational side, from tensions of various degrees of significance between some existing data sets [20, 21, 22, 23, 24, 25, 26, 27, 28]. Besides, based on the most minimal a priori assumptions, model independent reconstructions of the evolution of the dark energy (DE) equation of state (EoS) parameter w [29, 30, 31] and also model independent diagnoses [32, 33] exhibit a dynamical behaviour of w(z).
Nevertheless, even small deviations from/corrections to the \(\Lambda \)CDM mostly imply/require profound, and sometimes highly non-trivial, modifications to the fundamental theories of physics [34]. Indeed, we still do not have a promising and concrete fundamental theory leading to DE models that are more general than \(\Lambda \) and also can account for the small, but significant, deviations from the \(\Lambda \)CDM model as persistently suggested by the high precision data. Though, depending on the characteristics of the DE models favored by observations, we can decide whether it is more natural to consider DE as an actual physical source [15, 35] or as an effective source originated from a modification [36, 37, 38, 39, 40] to the standard theory of gravity, i.e., general relativity (GR). For instance, eliminating \(\Lambda \) by detecting \(w >-1\) would not be by itself illuminating to the nature of DE, but any detection of \(w < -1\) would be very illuminating to the nature of gravitation. For a perfect barotropic fluid, the adiabatic sound speed \(c^2_a\) is the physical propagation speed of perturbations, and therefore \(w<-1\) (viz., phantom [41] or quintom [42, 43] DE models) typically accompanying by \(c^2_a<0\), implies the instability of perturbations (and/or ghost instabilities, see, e.g., [44]), whereas, as shown in [45, 46], scalar-tensor theories of gravity, such as Brans–Dicke (BD) theory, can lead to an effective source with \(w<-1\) that does not correspond to the change of the sign of \(c^2_a\), and hence perturbations can still be stable. Similarly, any detection of DE with negative energy density would also be hugely informative, such that this for an actual physical source is of course physically ill, whereas it can be an effective source so that negative energy density does not lead to any pathology since this is not the true energy density. On the observational side, the constraints on the EoS parameter of DE persistently indicate that \(w\approx -1\), and so do not exclude \(w < -1\) [12] and it has recently been shown in [21, 32, 47] that DE models with energy densities passing below zero at high redshifts fit the data better and can address the tensions relevant to Lyman-\(\alpha \) forest measurements [20]. Such DE sources are indeed possible, in general, in modified theories of gravity with an effective gravitational coupling strength weaker in the past, e.g., in scalar-tensor theories, if we collect all modifications to the usual Einstein field equations to define an effective DE. Consequently, it can be argued that any detection of a deviation from \(\Lambda \)CDM model implies that the late time acceleration is not driven by \(\Lambda \)1 and of DE yielding \(w< -1\) and/or \(\rho <0\) implies that gravity is not minimally coupled (which eliminate the actual perfect-fluid models of DE and stands as a strong sign in favor of the effective DE models from modified gravity).
Generically, the effective DE models from modified gravity induce non-zero anisotropic stresses,2 which can be realized if we relax the isotropic space assumption of the \(\Lambda \)CDM model along with the modification to GR. Of course, we are not able to observe anisotropic stresses directly, yet they can reveal themselves through their effect on the evolution of the expansion anisotropy as well as on the average expansion rate of the Universe. Then, a natural question is whether the observations allow or suggest an anisotropic space, unless otherwise anisotropic stresses would be irrelevant. Indeed, there are some clues for questioning the spatially maximally symmetric Universe assumption, i.e., Robertson–Walker (RW) background, of the \(\Lambda \)CDM model. This has been mainly motivated by hints of anomalies in the CMB distribution first observed on the full sky by the WMAP experiment [49, 50, 51, 52] and so remained in Planck experiment [53, 54, 55, 56]. So far, the local deviations from the statistically highly isotropic Gaussianity of the CMB in some directions (the so called cold spots) could not have been excluded at high confidence levels [52, 53, 57, 58]. Furthermore, it has been shown that the CMB angular power spectrum has a quadrupole power lower than expected from the best-fit \(\Lambda \)CDM model [59, 60]. Several explanations for this anomaly have been proposed [61, 62, 63, 64, 65, 66, 67, 68, 69, 70] including the anisotropic expansion of the Universe. On the other hand, inflation (canonical) isotropizes the Universe very efficiently [71, 72], leaving a residual anisotropy that is negligible for any practical application in the observable Universe. This could be irrelevant if an anisotropic expansion is developed only well after the matter-radiation decoupling, for example during the domination of DE, say, by means of its anisotropic pressure acting as a late-time source of not insignificant anisotropy [69, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83] (see also, e.g., [84, 85, 86, 87] for constraint studies on the anisotropy of DE). Indeed, the CMB provides very tight constraints on the anisotropy at the time of recombination [88, 89, 90] of the order of the quadrupole temperature, i.e., \((\Delta T/T)_{\ell =2} \sim 10^{-5}\). And, in the simplest anisotropic generalization of the standard \(\Lambda \)CDM (viz., replacing the spatially flat RW metric by Bianchi type I metric), the energy density corresponding to the expansion anisotropy scales as the inverse of the square of the comoving volume, \(\rho _{\sigma ^2}\propto S^{-6}\), which implies an isotropization of the expansion from the recombination up to the present, leading to the typically derived upper bounds on the corresponding density parameter today as \(\Omega _{\sigma ^2,0}\sim 10^{-20}\). However, this is true if the anisotropic expansion is not generated by any anisotropic source, say, by an anisotropic DE (effective or actual), arising after decoupling [69, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83]. Nevertheless, almost all of these strong constraints on the expansion anisotropy in the late Universe, \(\Omega _{\sigma ^2,0}\lesssim 10^{-22}\), assume the non-existence of anisotropic sources (actual or effective) in the late universe, and are derived in fact from the contribution of the expansion anisotropy on the average expansion rate of the Universe [91, 92, 93, 94, 95]. On the other hand, direct observational constraints, e.g., from SNIa data, on the expansion anisotropy of the later Universe are much weaker, namely, on the order of \(\Omega _{\sigma ^2,0}\lesssim 10^{-4}\) for \(z\sim 0\) [96, 97].
It is in fact always possible to imitate any modification to the Einstein field equations as GR with an arbitrary mixture of scalars, vectors and tensors [98, 99]. Therefore, one may think of abandoning any attempt to distinguish between the actual physical DE and the effective DE from modified gravities. However, a modified theory of gravity, e.g., the Brans–Dicke theory, can modify both the average expansion rate of the Universe and the evolution of the expansion anisotropy in a distinctive way depending on its free parameter characterizing the model, e.g., the Brans–Dicke parameter \(\omega \). So, of course, owing to Occam's razor, looking for these two concurrent modifications predicted by a specific modified theory of gravity and testing them against the observational data can provide us a strong reason for favoring or disfavoring the modified theory of gravity under consideration over GR. In this paper, for instance, the effective EoS parameter corresponding to the expansion anisotropy, the present day EoS of the effective DE and the redshift at which it evolves into the phantom region range in specific intervals as \(1< w_{\sigma ^2}\le \frac{5}{3}\), \(-1.14\lesssim w<-1\) and \(0.50 \le z_\mathrm{PDL}< 0.65\) (see [100] for a similar result) for positive values of the BD parameter \(\omega \). Indeed, revealing the origin of the late time acceleration of the Universe considering such distinctive features of the DE models from modified theories of gravity is among the scientific themes underlying some upcoming experiments (e.g., [101]).
Finally, we have also independent motivations for modifying GR by fundamental theoretical physics. Namely, almost all of the attempts to quantize GR introduce deviations from it in the form of extra degrees of freedom, higher powers of the curvature in the action, higher order derivatives in the field equations, or non-local terms. The low-energy limits of the string theories typically yield BD gravity with \(\omega =-1\) and similarly d-brane models yield BD gravity with \(\omega \sim -1\) [102, 103, 104, 105]. Indeed, among all possible alternatives to GR, it is found natural by many to first consider scalar-tensor theories [106, 107, 108, 109, 110], which only add a (usually massive) scalar degree of freedom, the BD-like scalar field, to the two massless, spin-2 polarizations (gravitons) contained in the metric tensor. In relevance with our focus in this paper, the field equations in such models can be regarded as effective Einstein field equations and so the terms originating from the BD-like scalar field and its derivatives, moved to the right hand side, can always be interpreted as an effective fluid. Within this approach, the correspondence between general scalar-tensor theories and this effective fluid has been worked out explicitly first in [111] (see also [112] for the preceding work), and has shown that, in general, the corresponding effective fluid is imperfect. A recent work [113] extended and completed the correspondence showing how a symmetry of BD theory translates into a symmetry of this fluid such that this correspondence is valid for any spacetime geometries. Brans and Dicke's 1961 theory [106, 110], motivated first by the implementation of Mach's principle in gravity, today provides a prototype of scalar-tensor theories [107, 108, 109].
The original BD theory [106, 110] has only one additional constant parameter \(\omega \), determining the deviations from GR, which is recovered in the limit of \(|\omega | \rightarrow \infty \). Theoretically \(\omega \ge -\frac{3}{2}\) is necessary to avoid the Jordan/scalar field (\(\varphi \)) from yielding negative energy density values in the Einstein frame that leads, e.g., the Minkowski vacuum to be unstable. Observations on a wide range of scales further constrain BD theory around GR. The tightest constraints, to date, are imposed by the observations of gravitational phenomena in the Solar System; specifically, \(\omega \gtrsim 40000\) at the 2\(\sigma \) confidence level (C.L.) from observations of radio signals from the Cassini spacecraft as it passes behind the Sun [114]. It can be vastly relaxed for the massive BD theory, i.e., when the Jordan field is accompanied by a potential \(U(\varphi )\), as in this case one should consider the corresponding mass (M) of the Jordan field along with the BD parameter, viz., the \(\{\omega ,M\}\) parameter region [115]. It is shown that \(\omega =O(1)\) is allowed for the Solar System constraints provided that \(M\gtrsim 2\times 10^{-17}\mathrm eV\) [115]. The constraints on M, however, are typically subject to the cosmological observations, since, in particular, it becomes cosmologically significant at substantially lower scales, viz., at the Hubble mass scale \(M_{H_{0}}\simeq 10^{-33}\mathrm eV\), which leaves the constraint from the Cassini spacecraft still valid.
The cosmological constraints on \(\omega \) are complementary to the local constraints as they probe different length and time scales, as well as different epochs of the Universe. The strongest cosmological constraint, which assumes \(U(\varphi )=\mathrm const.\) and considers CMB data from the first Planck release along with the constraints from Big Bang Nucleosynthesis (BBN) light element abundances, gives \(\omega >890\) at the 2\(\sigma \) C.L. [116]. The constraints from Planck 2015 (2013) and a compilation of baryon acoustic oscillations (BAO) data are \(\omega >333\) (\(\omega >208\)) at the 2\(\sigma \) C.L., weakly dependent on the exponent \(n>0\) for a power-law potential \(U(\varphi )\propto \varphi ^n\) mimicking \(\Lambda \) in the late universe [117, 118] (also see [119] for a recent work considering a more complicated potential). The BBN leads to \(\omega \gtrsim 277\) (for a universe containing dust, radiation and conventional vacuum energy), which assumes power-law solutions and, by using general solutions, can be somewhat altered depending on the behaviour of the Jordan field in the early universe [120]. It is forecast that future cosmological experiments will be able to constrain \(\omega \) at a level comparable to the local tests [121, 122, 123].
The typical mass scale of \(M\sim 10^{-33}\mathrm eV\) imposed by cosmological observations is so small that the local constraint \(\omega \gtrsim 40000\) from the Cassini spacecraft must be satisfied at all times in all parts of the Universe, and leads to the conclusion that the massive BD extension of the standard \(\Lambda \)CDM model must be phenomenologically very similar to it throughout most of the history of the Universe, implying that the cosmological features of BD theory would be observationally irrelevant. On the other hand, it is still worth studying in detail and understanding the cosmological constructions within massive BD theory by considering the cosmological constraints on \(\omega \) only, as BD theory is representative of some more general gravity theories that can evade the local constraints (see [37] and references therein for further reading). For instance, straightforwardly, it is possible to consider a further extension of BD theory by allowing the BD parameter to be some functions of the Jordan field, \(\omega =\omega (\varphi )\), in addition to the presence of potential \(U(\varphi )\) and this leads to the possibility of having an \(\omega \) small enough to be significant in the past, while being large enough today to be consistent with the tight constraints imposed by observations of gravitational phenomena in the Solar System. Alternatively, for instance, in the Horndeski theory (the most general scalar–tensor theory having second-order field equations in four dimensions) [124, 125], while on cosmological scales BD theory emerges as an approximation, the derivative self-interactions of the Horndeski scalar can be large enough on small scales (higher curvature than cosmological environments) to screen the scalar resulting at the same time to the recovery GR. For example, it was proposed in [126] that the so called Galileon theories (a close cousin of BD theory) could also explain the late time acceleration of the Universe while evading solar system constraints. A different kind of example is that, in the presence of extra dimensions, it is possible to make BD theory (massive or massless) consistent with gravitational tests (including solar system tests) for \(|\omega |=O(1)\) [127]. See also [128] for a recent overview on cosmological probes of gravity for testing deviations from GR in relevance with such modified gravity theories.
Motivated by the discussions above, we present a detailed theoretical and observational investigation of the anisotropic (LRS Bianchi type I metric) massive Brans–Dicke gravity extension of the standard \(\Lambda \)CDM model, which is characterized by two additional free parameters, namely, the BD parameter \(\omega \) and the corresponding present day density parameter of the expansion anisotropy \(\Omega _{\sigma ^2,0}\). The role of the cosmological constant is taken over by the Jordan field potential of the form \(U(\varphi )\propto \varphi ^2\) and the physical ingredient of the Universe is considered as in the standard \(\Lambda \)CDM. We confine the present work, basically, to the non-negative values of the BD parameter since, here, we mainly consider the extension as a correction to the standard \(\Lambda \)CDM.
2 Anisotropic massive Brans–Dicke gravity extension of the \(\Lambda \)CDM model
We consider the BD action [106, 110] written in the Jordan frame in the following form:
$$\begin{aligned} \begin{aligned} S_\mathrm{JBD}=&\int \mathrm{d}^4 x \sqrt{-g}\bigg [\frac{\varphi ^2}{8}R-\omega \bigg (\frac{1}{2}\nabla _{\mu }\varphi \nabla ^{\mu }\varphi +\frac{1}{2}M^2\varphi ^2\bigg )\bigg ]\\&+S_\mathrm{Matter}, \end{aligned} \end{aligned}$$
where \(\varphi =\varphi (t)\) is the Jordan scalar field (function of cosmic time t only) and \(\omega =\mathrm const.\) is the Brans–Dicke parameter, R is the Ricci scalar, g is the determinant of the metric \(g_{\mu \nu }\), and \(S_\mathrm{Matter}\) is the matter action, which is independent of \(\varphi \) so that the weak equivalence principle is satisfied. It is clear from the way of writing the action that the term \(M^2\) stands as the bare mass-squared of the Jordan field.3 We assume \(M^2=\mathrm{const.}\) so that, as can also be seen from the action, it stands like a cosmological constant as \(2\omega M^2\equiv \Lambda \) and thereby can drive accelerated expansion. Hence, switching to massive BD from GR with a positive cosmological constant provides us with an opportunity to construct \(\Lambda \)CDM-type cosmologies, such that the mass of the Jordan field alone can play the role of positive cosmological constant like in the standard \(\Lambda \)CDM cosmology provided that \(2\omega M^2\equiv \Lambda >0\), and the Jordan field \(\varphi \) varying slowly enough on the top of this can account for small deviations from the standard \(\Lambda \)CDM model in a particular way, which in turn may lead to an improved fit to the observational data w.r.t. the standard \(\Lambda \)CDM model. In line with that, we intend to study the BD extension of the standard \(\Lambda \)CDM model as a correction, and therefore we demand the term \(2\omega M^2\) to be positive definite, which requires \(\omega >0\) as long as we keep \(M^2>0\) to avoid the Jordan field from having an imaginary mass.4 Hence, in this study, unless otherwise is mentioned, we carry out our investigations by assuming \(\omega >0\), which is already stronger than the assumption \(\omega \ge -3/2\) to avoid the Jordan field from yielding negative energy density values in the Einstein frame that leads, for instance the Minkowski vacuum to be unstable.
We consider the simplest anisotropic generalization of the spatially flat and homogeneous spacetime, i.e., locally rotationally symmetric (LRS) Bianchi type I metric, which can be written as follows;
$$\begin{aligned} \mathrm{d}s^2 = -\mathrm{d} t^2+S^2\left[ e^{4\beta }\mathrm{d}x^2+ e^{-2\beta }(\mathrm{d}y^2+\mathrm{d}z^2)\right] , \end{aligned}$$
where \(S=S(t)\) is the mean scale factor. Here, the exponent \(\beta =\beta (t)\) satisfies the relation \(\dot{\beta }^2=\frac{1}{6}\sigma ^2\), where \(\sigma ^2=\sigma _{ij}\sigma ^{ij}\) and \(\sigma ^{ij}\) are shear scalar and tensor, respectively. The scalar curvature for this metric (2) may be written in terms of the mean scale factor and shear scalar as \(R=-6\left( \frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}\right) -\sigma ^2\). Throughout the paper a dot denotes derivative w.r.t. cosmic time t.
We consider all types of matter distribution (namely, the usual cosmological sources such as radiation, baryons, etc.) as isotropic perfect fluids, which are described by the following energy-momentum tensor (EMT)
$$\begin{aligned} {T_{\mu }}^{\nu }={\text {diag}}[-\rho ,p,p,p], \end{aligned}$$
where \(\rho \) and p are the energy density and pressure, respectively.
The field equations for the action (1) within the framework of the metric (2) in the presence of (3) read:
$$\begin{aligned}&3\frac{\dot{S}^2}{S^2}-\frac{1}{2}\sigma ^2-2\omega \frac{\dot{\varphi }^2}{\varphi ^2}+6\frac{\dot{S}}{S} \frac{\dot{\varphi }}{\varphi }-2\omega M^2=\frac{4}{\varphi ^2}\rho , \end{aligned}$$
$$\begin{aligned}&2\frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}+\frac{1}{2}\sigma ^2-\sqrt{\frac{2}{3}}\left[ \dot{\sigma }+\left( 3\frac{\dot{S}}{S}+2\frac{\dot{\varphi }}{\varphi } \right) \sigma \right] \nonumber \\&\quad +(2 \omega +2) \frac{\dot{\varphi }^2}{\varphi ^2} +2\frac{\ddot{\varphi }}{\varphi }+4\frac{\dot{S}}{S}\frac{\dot{\varphi }}{\varphi }-2 \omega M^2=-\frac{4}{\varphi ^2}p, \end{aligned}$$
$$\begin{aligned}&2\frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}+\frac{1}{2}\sigma ^2+\sqrt{\frac{1}{6}}\left[ \dot{\sigma }+\left( 3\frac{\dot{S}}{S}+2\frac{\dot{\varphi }}{\varphi } \right) \sigma \right] \nonumber \\&\quad +(2 \omega +2) \frac{\dot{\varphi }^2}{\varphi ^2}+2\frac{\ddot{\varphi }}{\varphi }+4\frac{\dot{S}}{S} \frac{\dot{\varphi }}{\varphi }-2 \omega M^2=-\frac{4}{\varphi ^2}p, \end{aligned}$$
for which the latter two correspond to the pressure equations along the x-axis and the y- and z-axes, respectively. The Klein–Gordon equation for the Jordan field reads
$$\begin{aligned} \frac{\ddot{\varphi }}{\varphi }+3\frac{\dot{S}}{S}\frac{\dot{\varphi }}{\varphi } -\frac{3}{2\omega }\left( \frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}\right) -\frac{1}{2\omega }\frac{\sigma ^2}{2}+M^2=0. \end{aligned}$$
We consider the standard cosmological fluids, namely, pressureless fluid (\(p_{\mathrm{m}}= 0\)) and radiation/relativistic fluids (\(p_{\mathrm{r}} = \rho _\mathrm{r}/3\)), so that, here, \(\rho =\rho _{\mathrm{m}}+\rho _{\mathrm{r}}\) and \(p=p_{\mathrm{r}}\). If the Jordan field is approximately constant then the first Friedmann Eq. (4) becomes \(3\frac{\dot{S}^2}{S^2}\approx \frac{4}{\varphi ^2}\rho +\frac{\sigma ^2}{2}+2\omega M^2\) (as like in GR), where the effective cosmological gravitational coupling strength is then given by \(G=\frac{\varphi _0^2}{\varphi ^2}\, G_0\), where \(G=\frac{1}{2\pi \varphi ^2}\) and \(G_0=\frac{1}{2\pi \varphi ^2_0}\) with zero subscript denoting today's value (throughout the paper).5
We consider the conventional representation of the true equation for scalar-tensor gravity interacting with matter in the Einsteinian form with the constant \(G_0=G(t_0)\)
$$\begin{aligned} R_{\mu \nu }-\frac{1}{2}R g_{\mu \nu }=8\pi G_0 (T_{\mu \nu ,\mathrm{m}}+T_{\mu \nu ,\mathrm{DE}}), \end{aligned}$$
\(8\pi G_0=4/\varphi ^2_0\) following Refs. [45, 46], where it is discussed that modified gravity theories can be recast in the standard GR form such that (8) where all the new geometrical terms are grouped (on the r.h.s.) to form an effective DE contribution denoted as \(T_{\mu \nu ,\mathrm{DE}}\) (for a detailed discussion, see introduction in Ref. [46]). We note that the shear scalar is kept on the left hand side as it is a part of the Einstein tensor, namely, the metric itself only. Accordingly, we define the effective DE in the following way:
$$\begin{aligned}&\begin{aligned} 3\frac{\dot{S}^2}{S^2}-\frac{1}{2}\sigma ^2=\frac{4}{\varphi ^2_0}\,(\rho +\rho _\mathrm{DE}), \end{aligned} \end{aligned}$$
$$\begin{aligned}&\begin{aligned} 2\frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}+\frac{1}{2}\sigma ^2-\sqrt{\frac{2}{3}}\left( \dot{\sigma }+3\frac{\dot{S}}{S}\sigma \right) = -\frac{4}{\varphi ^2_0}\,(p+p_{\mathrm{DE},x}), \end{aligned} \nonumber \\ \end{aligned}$$
$$\begin{aligned}&\begin{aligned} 2\frac{\ddot{S}}{S}+\frac{\dot{S}^2}{S^2}+\frac{1}{2}\sigma ^2+\sqrt{\frac{1}{6}}\left( \dot{\sigma }+3\frac{\dot{S}}{S}\sigma \right) = -\frac{4}{\varphi ^2_0}\,(p+p_{\mathrm{DE},y}), \end{aligned} \nonumber \\ \end{aligned}$$
where \(\rho _{\mathrm{DE}}\) is the energy density and \(p_{\mathrm{DE},x}\) and \(p_{\mathrm{DE},y}\) are the principal pressures of the effective DE along the x-axis and the y- and z-axes, respectively, and which read
$$\begin{aligned} \rho _{\mathrm{DE}}= & {} \frac{\varphi _0^2}{4}\bigg [\left( \frac{4}{\varphi ^2}-\frac{4}{\varphi _0^2} \right) \rho + 2\omega \frac{\dot{\varphi }^2}{\varphi ^2}-6\frac{\dot{S}}{S}\frac{\dot{\varphi }}{\varphi }+2\omega M^2\bigg ], \nonumber \\ \end{aligned}$$
$$\begin{aligned} p_{\mathrm{DE},x}= & {} \frac{\varphi _0^2}{4}\bigg [\left( \frac{4}{\varphi ^2}-\frac{4}{\varphi _0^2} \right) p-2\sqrt{\frac{2}{3}}\frac{\dot{\varphi }}{\varphi }\sigma \nonumber \\&+(2 \omega +2) \frac{\dot{\varphi }^2}{\varphi ^2} +2\frac{\ddot{\varphi }}{\varphi }+4\frac{\dot{S}}{S}\frac{\dot{\varphi }}{\varphi }-2\omega M^2\bigg ],\end{aligned}$$
$$\begin{aligned} p_{\mathrm{DE},y}= & {} \frac{\varphi _0^2}{4}\bigg [\left( \frac{4}{\varphi ^2}-\frac{4}{\varphi _0^2} \right) p+\sqrt{\frac{2}{3}}\frac{\dot{\varphi }}{\varphi }\sigma \nonumber \\&+(2 \omega +2) \frac{\dot{\varphi }^2}{\varphi ^2} +2\frac{\ddot{\varphi }}{\varphi }+4\frac{\dot{S}}{S}\frac{\dot{\varphi }}{\varphi }-2\omega M^2\bigg ]. \end{aligned}$$
We note that the Jordan field and shear scalar are coupled and they together lead to an anisotropy in the pressure of the effective DE, which may be represented by
$$\begin{aligned} \Delta p_{\mathrm{DE}}=p_{\mathrm{DE},y}-p_{\mathrm{DE},x}=\sqrt{6}\frac{\varphi _0^2}{4}\frac{\dot{\varphi }}{\varphi }\sigma . \end{aligned}$$
This, in turn, leads to a modification in the evolution of the expansion anisotropy provided that \(\varphi \) is dynamical and the space is not exactly isotropic as would be expected from a realistic cosmological model (see, e.g., [132]). We can restate this equation in terms of the cosmological gravitation coupling strength as
$$\begin{aligned} \Delta p_{\mathrm{DE}}=-\sqrt{\frac{3}{2}}\frac{1}{8\pi G_0}\,\frac{\dot{G}}{G}\,\sigma . \end{aligned}$$
This effective anisotropic pressure (in the presence of shear) induced by the Jordan field in the Jordan frame is absent in the Einstein frame (see, e.g., Ref. [132]).
3 Exact solution compatible with standard \(\Lambda \)CDM model
We follow the method, given in [100], for obtaining exact cosmological solutions of a scalar-tensor gravity theory compatible with the \(\Lambda \)CDM model, albeit we consider LRS Bianchi type I spacetime rather than spatially flat Robertson–Walker spacetime.6\(^{,}\)7 To do so, we first re-express the set of differential equations given by Eqs. (4)–(7) using \(\frac{\mathrm{d}z}{\mathrm{d} t}=-H(1+z)\), where \(H=\dot{S}/S\) is the average Hubble parameter and z is the cosmic redshift we define in terms of S as \(z=-1+\frac{1}{S}\). Accordingly, we re-express the energy density equation (4) as
$$\begin{aligned}&3H^2\bigg [1-\frac{2\omega }{3}(1+z)^2\frac{{\varphi '}^2}{\varphi ^2}-2(1+z)\frac{{\varphi '}}{\varphi }\bigg ]\nonumber \\&\quad -\frac{\sigma ^2}{2}-2\omega M^2=\frac{4}{\varphi ^2}\rho , \end{aligned}$$
and obtain the pressure equation as
$$\begin{aligned}&-\frac{\mathrm{d}H^2}{\mathrm{d}z}\left[ 1+z-(1+z)^2\frac{\varphi '}{\varphi }\right] \nonumber \\&\quad +3H^2\bigg \{1+(1+z)^2\bigg [\frac{2}{3}\frac{\varphi ''}{\varphi }+\frac{2}{3}(1+\omega )\frac{{\varphi '}^2}{\varphi ^2}\bigg ]\nonumber \\&\quad -\frac{2}{3}(1+z)\frac{\varphi '}{\varphi }\bigg \}+\frac{\sigma ^2}{2}-2\omega M^2=-\frac{4}{\varphi ^2}p, \end{aligned}$$
using Eqs. (5) and (6), the shear propagation equation as
$$\begin{aligned} \frac{\sigma '}{\sigma }-\frac{3}{1+z}+2\frac{\varphi '}{\varphi }=0, \end{aligned}$$
by subtracting Eq. (5) from Eq. (6), and finally re-express the scalar field equation (7) as
$$\begin{aligned}&-\frac{\mathrm{d}H^2}{\mathrm{d}z}\left[ \omega (1+z)^2\frac{\varphi '}{\varphi }+\frac{3}{2}(1+z)\right] \nonumber \\&\quad +H^2\left[ -2\omega \frac{\varphi ''}{\varphi }(1+z)^2+4\omega (1+z)\frac{\varphi '}{\varphi }+6\right] \nonumber \\&\quad +\frac{\sigma ^2}{2}-2\omega M^2=0, \end{aligned}$$
where \('\) denotes derivative with respect to redshift (\(\mathrm{d}/\mathrm{d}z\)).
We note that Eqs. (18) and (20) have exactly the same mathematical form, that is, a first order linear differential equation in \(H^2\), such as
$$\begin{aligned} A_i(z)\frac{\mathrm{d}H^2}{\mathrm{d}z}+B_i(z)H^2+C_i(z)+D_i=0,\,\, i=1,2, \end{aligned}$$
provided that we set \(p=p_{\mathrm{m}}=0\), viz., the universe is filled with only pressureless matter. We next note that constants \(D_1=D_2=-2\omega M^2\), and that the shear scalar \(\sigma ^2\) is the term that differs in our system of equations from the one given in [100]. Fortunately, it contributes to (18) and (20) in the same way, namely, as \(C_1(z)=C_2(z)=\frac{\sigma ^2}{2}\). In accordance with these points, we assume \(A_1(z)=A_2(z)\) in (18) and (20) and then solve for the rate of change of Jordan field in z as
$$\begin{aligned} \frac{\varphi '}{\varphi }=-\frac{1}{2(1+\omega )}\frac{1}{1+z}, \end{aligned}$$
which in turn renders the coefficients of \(H^2\) identical, i.e, \(B_1(z)=B_2(z)\), as well. Thus, integrating (22), it turns out that the solution of the Jordan field is
$$\begin{aligned} \varphi ^2=\varphi _0^2(1+z)^{-\frac{1}{1+\omega }}, \end{aligned}$$
where \(\varphi _0=\varphi (z=0)\) is the present time value of the Jordan field, and which in turn gives
$$\begin{aligned} G=G_{0} (1+z)^{\frac{1}{1+\omega }}. \end{aligned}$$
The integration of the shear propagation equation given in Eq. (19) gives
$$\begin{aligned} \sigma ^2=\sigma _0^2(1+z)^6\left( \frac{\varphi }{\varphi _0}\right) ^{-4}, \end{aligned}$$
where \(\sigma ^2_0=\sigma ^2(z=0)\) is the present time value of the shear scalar, leading, together with (23), to
$$\begin{aligned} \sigma ^2=\sigma _0^2(1+z)^{6+\frac{2}{1+\omega }}, \end{aligned}$$
whereas it is \(\sigma _\mathrm{GR}^2 \propto (1+z)^6\) in GR (corresponding to the \(| \omega | \rightarrow \infty \) limit of the BD gravity), and thereby the energy density corresponding to the shear scalar can be defined as
$$\begin{aligned} \rho _{\sigma ^2}\equiv \frac{\sigma ^2}{2} \frac{\varphi _0^2}{4}=\frac{\sigma ^2_0}{2} \frac{\varphi _0^2}{4}(1+z)^{6+\frac{2}{1+\omega }}, \end{aligned}$$
by considering the present time value of the Jordan field, i.e., of the cosmological gravitational coupling strength. Because BD theory [110] satisfies the local conservation of EMT, the energy density of the pressureless matter can immediately be written as
$$\begin{aligned} \rho _{\mathrm{m}}=\rho _{\mathrm{m,0}}(1+z)^3, \end{aligned}$$
which is exactly the same as the one for the pressureless matter in the standard \(\Lambda \)CDM model. Finally using \(\varphi ^2\), \(\sigma ^2\) and \(\rho _{\mathrm{m}}\) from Eqs. (23), (26) and (28), respectively, in Eq. (17) we reach the following modified anisotropic Friedmann equation for BD theory in this solution
$$\begin{aligned} \begin{aligned} H^2=&\frac{4}{3}\frac{\gamma }{\varphi _0^2}\left[ \rho _{\mathrm{M}}+\rho _{\mathrm{m},0}(1+z)^{3+\frac{1}{1+\omega }}+\rho _{\sigma ^2,0}(1+z)^{6+\frac{2}{1+\omega }}\right] , \end{aligned} \end{aligned}$$
$$\begin{aligned} \gamma =\frac{6(1+\omega )^2}{(3\omega +4)(2\omega +3)}, \end{aligned}$$
and the energy densities corresponding to the mass of the Jordan field and to the expansion anisotropy today read, respectively,
$$\begin{aligned} \begin{aligned} \rho _{\mathrm{M}}=2M^2\omega \frac{\varphi _0^2}{4}\quad \text {and} \quad \rho _{\sigma ^2,0}=\frac{\sigma _0^2}{2}\frac{\varphi _0^2}{4}. \end{aligned} \end{aligned}$$
One may check for consistency that we recover the H(z) of the standard \(\Lambda \)CDM model when we set \(\omega \rightarrow \infty \) (viz., in the GR limit of BD gravity giving \(\varphi \rightarrow \mathrm{constant}\)) and \(\rho _{\sigma ^2,0}=0\) (viz., isotropic expansion as in the RW spacetime metric). We also note that \(\rho _{\mathrm{M}}\) can be regularized to give a finite positive value consistent with the observations by setting \(M^2\omega =\mathrm{constant}\) so that \(M^2\rightarrow 0\) as \(\omega \rightarrow \infty \) such that the term \(M^2\omega \) arising from their multiplication would always remain unaltered and finite.
3.1 Inclusion of radiation
The solution of the BD extension of standard \(\Lambda \)CDM model, for either isotropic or anisotropic spatially flat spacetimes, can be achieved by using the method given in Ref. [100] provided that \(p=0\), viz., the Universe is filled with only pressureless matter. It covers also the epoch of the Universe when the contributions from the pressureless matter in (29), namely, \(\frac{4}{\varphi ^2}\rho _{\mathrm{m}}\propto (1+z)^{3+\frac{1}{1+\omega }}\), is dominant over the effective cosmological constant \(2\omega M^2\) and one may check that for this epoch, neglecting anisotropy, it reduces to the well known solution with \(S\propto t^{\frac{2+2\omega }{4+3\omega }}\) and \(\varphi ^2\propto t^{\frac{2}{4+3\omega }}\) (see, e.g., [120] and references therein). We could not obtain analytically an exact cosmological solution – other than the trivial solution with \(\varphi =\mathrm const.\) – using the same method when we include relativistic source, \(p_{\mathrm{r}}=\frac{\rho _\mathrm{r}}{3}\). Fortunately the effects of the relativistic source on the background dynamics can be neglected in the late Universe. Namely, because local energy conservation holds in BD gravity in the Jordan frame, we have \(\rho _{\mathrm{r}}\propto (1+z)^4\) and \(\rho _{\mathrm{m}}\propto (1+z)^3\), implying that the effects of the relativistic source on the dynamics of the Universe will be significant for large redshift values, in particular, for the redshift values larger than the matter-radiation equality redshift, \(z_{\mathrm{eq}}=-1+\frac{\rho _{\mathrm{m,0}}}{\rho _{\mathrm{r,0}}}\sim 3380\) [12], as in the standard cosmology [138]. The attractor solution of the BD gravity for radiation domination is well known that it is exactly the GR solution, i.e., \(\varphi = \mathrm{const.}\) [120, 135, 136, 137, 138]. Indeed, in general, BD cosmologies have exact solutions which show that they are dominated by the Jordan field at early times and by the perfect fluid matter sources at late times, which here can be considered as the period all the way down from matter-radiation equality to the earlier times covering the times relevant to Big Bang Nucleosynthesis (BBN). Hence, for \(z>z_{\mathrm{eq}}\), we have \(\varphi = \varphi _1= \mathrm{const.}\), implying \(G=G_1= \mathrm{const.}\) as long as anisotropy is negligible. Note that \(\varphi _1\) (or \(G_1\)) will be different than \(\varphi _0\) (or \(G_0\)) as a consequence of the fact that the Jordan field has been evolving approximately in accordance with (23) for \(z<z_{\mathrm{eq}}\), i.e., between the matter-radiation equality and today. Accordingly, our model will more approach a general relativistic cosmology as the radiation dominates over pressureless matter as we go to higher redshifts, yet the most part of period of the Universe during which \(z<z_{\mathrm{eq}}\) the Jordan field evolves approximately in accordance with (23) and thereby we can estimate that
$$\begin{aligned} \begin{aligned} \varphi _1\sim \varphi _0(1+z_{\mathrm{eq}})^{-\frac{1}{2+2\omega }}\quad&\text {at}\quad z\sim z_{\mathrm{eq}} \\ \quad \quad \quad \quad \text {as well as} \quad \varphi \sim \varphi _1\quad&\text {for}\quad z\gtrsim z_{\mathrm{eq}}, \end{aligned} \end{aligned}$$
which implies
$$\begin{aligned} G_1\sim G_0(1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }} \quad \text {for}\quad z\gtrsim z_{\mathrm{eq}}, \end{aligned}$$
where \(G_1\) is then approximately the cosmological gravitational coupling strength for \(z>z_{\mathrm{eq}}\) (see [120] for further details). Thus, for \(\omega >0\), the cosmological gravitational coupling strength will remain constant (provided that anisotropy is negligible) during the radiation dominated epoch as just like in the standard cosmology, but will be larger than its present time value, e.g., \(G_1/G_0\sim 3380^{\frac{1}{1+\omega }}\). This implies an enhanced expansion rate of the Universe during this epoch, which in turn will affect, for instance, the primordial element abundances. We shall further discuss this point in Sect. 5. We next see from (19) that we approximately have \(\sigma ^2\propto \sigma _0^2(1+z)^6\) like in GR in the presence of isotropic perfect fluid, since we have \(\varphi \sim \mathrm{constant}\) when the radiation is dominant over pressureless matter and expansion anisotropy is still insignificant. As we go to even larger redshift values, say, \(z\gg z_{\mathrm{BBN}}\), expansion anisotropy will dominate over the terms due to the radiation and hence we obtain once again \(\rho _{\sigma ^2}\propto (1+z)^{6+\frac{2}{1+\omega }}\) (which can be obtained simply by setting \(\rho _{\mathrm{M}}=0\) and \(\rho _{\sigma ^2,0}=0\) in our solution).
In the light of the above discussion, although we do not have exact analytical solution when we include the radiation as well, we can find out the contribution from the radiation. First of all, it is obvious that during radiation domination the Jordan field is constant, then the inclusion of radiation to the model will in general slow-down the Jordan field, namely, the model in general will deviate less from GR. Of course, that effect will be insignificant for a viable and realistic cosmological model [namely, cosmological model satisfying, roughly, \(\rho _{\sigma ^2,0}\ll \rho _{\mathrm{r,0}}\ll \rho _{\mathrm{m,0}}\sim \rho _{\mathrm{M}}\) for \(\omega \gtrsim 10\) (see Sect. 3.3.1 for this choice on \(\omega \))] for \(z< z_{\mathrm{eq}}\). Thus, we can write
$$\begin{aligned} H^2\approx & {} \frac{4}{3}\frac{\gamma }{\varphi _0^2}\bigg [\rho _{\mathrm{M}}+\rho _{\mathrm{m},0}(1+z)^{3+\frac{1}{1+\omega }}+\rho _{\mathrm{r},0}(1+z)^{4+\frac{1}{1+\omega }}\nonumber \\&+\rho _{\sigma ^2,0}(1+z)^{6+\frac{2}{1+\omega }}\bigg ]\quad \text {for}\quad z<z_{\mathrm{eq}}, \end{aligned}$$
which can be safely used all the way to the recombination redshift (CMB release), such that it is well known that the Universe should have transited from radiation to matter dominated era when \(z\sim 3380\), and the recombination that leads to photon decoupling should have taken place when \(z\sim 1100\), at which the pressureless matter is still dominant over radiation. Therefore, we note that H(z) given in (34) can safely be used for constraining the model by using CMB radiation data as well.
For the times before the matter-radiation equality, setting \(\varphi \) constant relying on its attractor behavior during the radiation-dominated phase, we can write
$$\begin{aligned} H^2\approx & {} \frac{4}{3}\frac{\gamma }{\varphi _1^2}\left[ \rho _{\mathrm{r},0}(1+z)^{4}+\rho _{\sigma ^2,1}(1+z)^{6}\right] \nonumber \\&\quad \text {for}\quad z_{\mathrm{eq}}<z<z_\mathrm{eq,\sigma ^2,r}. \end{aligned}$$
Note that, for consistency between (34) and (35), here in (35) we have \(\varphi _1\sim \varphi _0(1+z_{\mathrm{eq}})^{-\frac{1}{2+2\omega }}\) in contrast to \(\varphi _{0}\) in (34), in accordance with the above discussions leading to (32), and similarly \(\rho _{\sigma ^2,1}\sim \rho _{\sigma ^2,0}(1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }}\). Using (35), we estimate (at least roughly) the radiation-expansion anisotropy equality redshift as \(z_\mathrm{eq,\sigma ^2,r}\sim -1+ \sqrt{\frac{\rho _{\mathrm{r,0}}}{\rho _{\sigma ^2,1}}}\), which can also be written as follows;
$$\begin{aligned} z_\mathrm{eq,\sigma ^2,r}\sim -1+ \sqrt{\frac{\rho _{\mathrm{r,0}}}{\rho _{\sigma ^2,0}}\left( \frac{\rho _{\mathrm{m,0}}}{\rho _{\mathrm{r,0}}}\right) ^{-\frac{1}{1+\omega }}}. \end{aligned}$$
Avoiding expansion anisotropy from spoiling the physical processes relevant to BBN will then roughly require \(z_\mathrm{eq,\sigma ^2,r}>z_{\mathrm{BBN}}\). Typically, in a viable cosmological model, we expect the matter-radiation equality redshift to be \(z_{\mathrm{eq}}=-1+\frac{\rho _{\mathrm{m,0}}}{\rho _{\mathrm{r,0}}}\sim 3380\), the BBN to take place at \(z_{\mathrm{BBN}}\sim 3\times 10^{8}\) and \(\Omega _{\mathrm{r,0}}\sim 10^{-4}\). Accordingly, taking \(z_\mathrm{eq,\sigma ^2,r}> 3\times 10^{8}\), we obtain \(\frac{\rho _{\sigma ^2, 0}}{\rho _\mathrm{r, 0}}=\frac{\Omega _{\sigma ^2, 0}}{\Omega _\mathrm{r, 0}}\lesssim 10^{-17-\frac{4}{1+\omega }}\), namely, the following stringent constraint on the density parameter of the expansion anisotropy today:
$$\begin{aligned} \Omega _{\sigma ^2,0}\lesssim 10^{-21-\frac{4}{1+\omega }}, \end{aligned}$$
which reduces to \(\Omega _{\sigma ^2,0}\sim 10^{-21}\) in the GR limit (\(\omega \rightarrow \infty \)). It is noteworthy that, in the GR limit, at which anisotropy in the pressure also disappears, our estimation on the expansion anisotropy for today is already of the same order of magnitude with the ones given in [94] from a general test of isotropy using CMB temperature and polarization data from Planck and in [91] from primordial nucleosynthesis, both of them consider GR and only isotropic fluids. This also shows the reliability of our estimation. We note that, for the case \(\omega >0\) we consider in this study, BD theory leads to more stringent constraints on the expansion anisotropy (e.g., \(\Omega _{\sigma ^2,0}\lesssim 10^{-25}\) for \(\omega =0\)), while the case \(\omega <-1\), which we do not consider in the present study, can relax it. For instance, Campanelli et al. [96] have put model-independent upper bounds on the expansion anisotropy, from type Ia SNe observations, in the late Universe, namely, for \(z\lesssim 1.6\), as \(\Omega _{\sigma ^2,0}\lesssim 10^{-4}\) (see also [97] for a recent work). We note that this value is achieved upon choosing \(\omega =-\frac{21}{17}\) in (37), which can be done without worrying for theoretical issues associated with \(\omega <-\frac{3}{2}\), but this will cost large deviations from GR (viz. rapidly varying strength of the gravitational coupling) as well as the standard \(\Lambda \)CDM. We shall further discuss possible modifications on the expansion anisotropy for different values/ranges of \(\omega \) and their comparisons with the standard \(\Lambda \)CDM model in Sect. 3.3.2.
Finally, for the times when the expansion anisotropy dominates over radiation (and hence obviously on pressureless matter), we can write
$$\begin{aligned} H^2\propto (1+z)^{6+\frac{2}{1+\omega }}\quad \text {for}\quad z>z_\mathrm{eq,\sigma ^2,r}, \end{aligned}$$
from our solution (see, e.g., Refs. [132, 139, 140, 141] for anisotropy dominated universes in the context of scalar-tensor theories). This epoch however is out of interest in the present work.
3.2 Effective anisotropic dark energy
The discussions in the previous subsections reveal that the massive BD gravity modification on the gravity sector, manifests itself when radiation becomes subdominant, say, in the relatively late Universe, that is, switching from GR+\(\Lambda \) to the massive BD introduces corrections on the redshift dependencies of both the average expansion rate and expansion anisotropy. Although we consider the presence of radiation when we constrain the model using the cosmological data, while we are discussing the effective DE in what follows, we shall neglect its presence due to the following reasons: (i) Its effect on the effective DE will obviously be negligible in the late universe. (ii) When it is dominant the Jordan field \(\varphi \) freezes, implying the BD gravity mimics GR, and thereby the effective DE will mimic nothing but \(\Lambda \) (not to mention that DE would be subdominant at this epoch). (iii) A practical reason, such that we don't have explicit analytical solution when we include radiation into the model. On the other hand, although we expect anisotropy to be negligible (even with respect to radiation) in the late universe, we keep expansion anisotropy due to the following reasons: (i) To see its effect on the effective DE explicitly. (ii) It leads to an anisotropy in the pressure of the effective DE, which in turn modifies the redshift dependency of the expansion anisotropy w.r.t. GR (\(|\omega |\rightarrow \infty \)). (iii) In relevance with (ii), its modified redshift dependency w.r.t. the one in GR (provided that it can be detected) would be a strong signal in favor of modified gravities (such as BD theory), rather than the presence of a DE (such as \(\Lambda \), scalar fields) ingredient of the Universe assuming that GR is the true theory of gravity.
We implement the solution given in Sect. 3 and accordingly we re-write (9) as follows:
$$\begin{aligned} \begin{aligned} 3H^2=\frac{4}{\varphi ^2_0}\,\left( \rho _{\mathrm{m}}+\rho _{\mathrm{DE}}+\rho _{\sigma ^2}\right) , \end{aligned} \end{aligned}$$
where the energy density of the effective DE, \(\rho _{\mathrm{DE}}\), reads
$$\begin{aligned} \rho _{\mathrm{DE}}=\gamma \rho _{\mathrm{M}}+\rho _{\mathrm{m}}\left[ \gamma \left( \frac{\rho _{\mathrm{m}}}{\rho _{\mathrm{m,0}}}\right) ^{\frac{1}{3(1+\omega )}}-1\right] +\rho _{\sigma ^2}(\gamma -1),\qquad \end{aligned}$$
which gives
$$\begin{aligned} \begin{aligned} \rho _{\mathrm{DE,0}}=&\,\gamma \rho _{\mathrm{M}}+\rho _{\mathrm{m,0}}\left( \gamma -1\right) +\rho _{\sigma ^2,0}(\gamma -1) \end{aligned} \end{aligned}$$
for \(z=0\). Dividing (39) by \(3H_0^2\), we obtain
$$\begin{aligned} \begin{aligned} \frac{H^2}{H_0^2}=\Omega _{\mathrm{DE,0}}f(z)+\Omega _{\mathrm{m,0}}(1+z)^3+\Omega _{\sigma ^2,0}(1+z)^{6+\frac{2}{1+\omega }}, \end{aligned} \end{aligned}$$
where we define \(f(z)=\frac{\rho _{\mathrm{DE}}}{\rho _{\mathrm{DE,0}}}\), reading
$$\begin{aligned} f(z)=\,1+\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}\left( \gamma -1\right) \left[ (1+z)^{6+\frac{2}{1+\omega }}-1\right] \nonumber \\ +\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}\left\{ (1+z)^3\left[ \gamma (1+z)^{\frac{1}{1+\omega }}-1\right] +1-\gamma \right\} , \end{aligned}$$
which gives \(f(z=0)=1\) consistently with \(\Omega _{\mathrm{m,0}}+\Omega _{\mathrm{DE,0}} + \Omega _{\sigma ^2,0}=1\). Note that we define the present time values of the density parameters as \(\Omega _{i,0}=\rho _{i,0}/\rho _\mathrm{c,0}\), where \(\rho _\mathrm{c,0}=3H_0^2\varphi _0^2/4\) is the present time value of the critical energy density, and that we eliminate \(\Omega _\mathrm{M,0}\), that appears in f(z) originally, using the relation
$$\begin{aligned} \begin{aligned} \Omega _{\mathrm{M},0}=\gamma ^{-1}\left[ \Omega _{\mathrm{DE},0}-(\gamma -1)(\Omega _{\mathrm{m,0}}+\Omega _{\sigma ^2,0})\right] , \end{aligned} \end{aligned}$$
which can easily be obtained from (40). Using (44) along with \(\Omega _{\mathrm{m,0}}+\Omega _{\sigma ^2,0}=1-\Omega _{\mathrm{DE,0}}\), we see that the present time density parameter of the effective DE is determined by the Jordan field's mass and the BD parameter \(\omega \) [viz., \(\gamma (\omega )\)] as8
$$\begin{aligned} \Omega _{\mathrm{DE},0}=\omega \frac{2}{3}\left( \frac{M}{H_0}\right) ^2+\frac{\gamma -1}{\gamma }. \end{aligned}$$
We obtain – upon straightforward arrangements in the equations obtained by using (29), (44) in the equations obtained by substituting (23), (30), (31) in (13), (14) – the principal EoS parameters of the effective DE along the x-axis and y- and z-axes in terms of the present day values of the density parameters as follows, respectively,
$$\begin{aligned} w_{\mathrm{DE},x}\equiv&\, \frac{p_{\mathrm{DE},x}}{\rho _{\mathrm{DE}}}=-1+\frac{1}{f(z)}\bigg [\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}\frac{2\omega +2}{2\omega +3}(1+z)^{3+\frac{1}{1+\omega }}\nonumber \\&\quad -\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}(1+z)^{3}-\frac{2}{2\omega +3}\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}(1+z)^{6+\frac{2}{1+\omega }}\nonumber \\&\quad -\frac{2}{3(1+\omega )}\sqrt{\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}} (1+z)^{3+\frac{1}{1+\omega }}\nonumber \\&\quad \times \sqrt{1+\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}(1+z)^{3}+\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}(1+z)^{6+\frac{2}{1+\omega }}} \bigg ], \end{aligned}$$
$$\begin{aligned} w_{\mathrm{DE},y}\equiv&\,\frac{p_{\mathrm{DE},y}}{\rho _{\mathrm{DE}}}=-1\nonumber \\&+\frac{1}{f(z)}\bigg [\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}\frac{2\omega +2}{2\omega +3}(1+z)^{3+\frac{1}{1+\omega }}-\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}(1+z)^{3}\nonumber \\&-\frac{2}{2\omega +3}\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}(1+z)^{6+\frac{2}{1+\omega }}\nonumber \\&+\frac{1}{3(1+\omega )}\sqrt{\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}}(1+z)^{3+\frac{1}{1+\omega }}\nonumber \\&\times \sqrt{1+\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}(1+z)^{3}+\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}(1+z)^{6+\frac{2}{1+\omega }}} \bigg ], \end{aligned}$$
which give
$$\begin{aligned} w_{\mathrm{DE},x,0}&=-1-\frac{\Omega _{\mathrm{m,0}}+2\,\Omega _{\sigma ^2,0}}{(2\omega +3)\Omega _{\mathrm{DE,0}}} -\frac{2\sqrt{\Omega _{\sigma ^2,0}}}{3(1+\omega )\Omega _{\mathrm{DE,0}}},\nonumber \\ \end{aligned}$$
$$\begin{aligned} w_{\mathrm{DE},y,0}&=-1-\frac{\Omega _{\mathrm{m,0}}+2\,\Omega _{\sigma ^2,0}}{(2\omega +3)\Omega _{\mathrm{DE,0}}} +\frac{\sqrt{\Omega _{\sigma ^2,0}}}{3(1+\omega )\Omega _{\mathrm{DE,0}}}\nonumber \\ \end{aligned}$$
for \(z=0\).
The anisotropy of the EoS of the effective DE may be given as the difference between its principal EoS parameters (\(\Delta w_{\mathrm{DE}}=w_{\mathrm{DE},y}-w_{\mathrm{DE},x}\)) as
$$\begin{aligned} \begin{aligned} \Delta w_{\mathrm{DE}}=\,&\frac{1}{(1+\omega )}\frac{(1+z)^{3+\frac{1}{1+\omega }}}{f(z)} \sqrt{\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}} \\&\times \sqrt{1+\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{DE,0}}}(1+z)^{3}+\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}(1+z)^{6+\frac{2}{1+\omega }}}. \end{aligned} \end{aligned}$$
Accordingly, we have
$$\begin{aligned} \Delta w_{\mathrm{DE,0}}=\frac{1}{1+\omega }\,\frac{\sqrt{\Omega _{\sigma ^2,0}}}{\Omega _{\mathrm{DE,0}}}\quad \text {for}\quad z=0, \end{aligned}$$
which is approximately equal to \(\frac{1}{1+\omega }\sqrt{\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{DE,0}}}}\), provided that the effective DE is dominant today (\(z=0\)). Additionally, we can write for dust domination
$$\begin{aligned} \Delta w_{\mathrm{DE}}\sim \frac{1}{(1+\omega )\gamma }\sqrt{\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{m,0}}}}\,(1+z)^{\frac{3}{2}}, \end{aligned}$$
which is positive definite as \(\omega >0\), and for expansion anisotropy domination
$$\begin{aligned} \Delta w_{\mathrm{DE}}\sim -\frac{1}{1+\omega }\,\frac{1}{1-\gamma }, \end{aligned}$$
which depends only on \(\omega \), is negative definite as \(\omega >0\), monotonic in \(\omega >0\) and takes values as \(\Delta w_{\mathrm{DE}}=\{-2,-1.2\}\) for \(\omega =\{0,\infty \}\).
We finally, using \(\dot{\rho }_{\mathrm{DE}}+3H\rho _{\mathrm{DE}}(1+{\bar{w}}_{\mathrm{DE}})=0\), define a mean/volumetric effective EoS parameter for the effective DE as
$$\begin{aligned} \begin{aligned} {\bar{w}}_{\mathrm{DE}}=-1+\frac{1+z}{3}\frac{\rho _{\mathrm{DE}}'}{\rho _{\mathrm{DE}}}=-1+\frac{1+z}{3}\frac{f'(z)}{f(z)}, \end{aligned} \end{aligned}$$
which explicitly reads
$$\begin{aligned} \begin{aligned} {\bar{w}}_{\mathrm{DE}}=&-1+\frac{\Omega _{\mathrm{m},0}(1+z)^3\left[ \gamma \frac{3\omega +4}{3(1+\omega )}(1+z)^{\frac{1}{1+\omega }}-1\right] +\Omega _{{\sigma ^2},0}\frac{2(\gamma -1)(3\omega +4)}{3(1+\omega )}(1+z)^{6+\frac{2}{1+\omega }}}{\Omega _{\mathrm{DE,0}}+\Omega _{\mathrm{m,0}} \left\{ (1+z)^3\left[ \gamma (1+z)^{\frac{1}{1+\omega }}-1\right] +1-\gamma \right\} +\Omega _{\sigma ^2,0}(\gamma -1)\left[ (1+z)^{6+\frac{2}{1+\omega }}-1\right] }, \end{aligned} \end{aligned}$$
$$\begin{aligned} \begin{aligned} {\bar{w}}_{\mathrm{DE,0}}=&-1-\frac{1}{2\omega +3}\,\frac{\Omega _{\mathrm{m},0}}{\Omega _{\mathrm{DE,0}}}-\frac{2(5\omega +6)}{3(2\omega +3)(\omega +1)}\,\frac{\Omega _{\mathrm{\sigma ^2},0}}{\Omega _{\mathrm{DE,0}}} \end{aligned} \end{aligned}$$
for \(z=0\). We note that, given \(\Omega _{\mathrm{DE,0}}>0\) and \(\Omega _{\mathrm{m,0}}\) and \(\Omega _{{\sigma ^2},0}\) are already positive definite, for \(\omega >0\), the contributions from both the pressureless matter and the expansion anisotropy to the present value of the mean EoS of the effective DE, \({\bar{w}}_{\mathrm{DE,0}}\), are negative, i.e., anisotropic BD extension of the standard \(\Lambda \)CDM model predicts phantom like effective DE in the present time Universe.
3.3 Preliminary investigations
3.3.1 Features of effective anisotropic dark energy
In this section, before the observational analysis, we discuss some of the features of the model, in particular those corresponding to the effective DE and the evolution of the expansion anisotropy. To do so, we start with justifying the positivity assumption on \(\omega \) in our study by making use of one of the important parameters about the kinematics of the Universe, the deceleration parameter. We define the mean/volumetric deceleration parameter in terms of the average Hubble parameter as \(q=-1+~\frac{H'}{H}(1+z)\). This reads, as we expect in a viable cosmology \(\Omega _{{\sigma ^2},0}<\Omega _{\mathrm{r,0}}\ll \Omega _{\mathrm{m},0}\),
$$\begin{aligned} \begin{aligned} q\approx -1+\frac{4+3\omega }{2+2\omega }\,\frac{\Omega _{\mathrm{m},0}(1+z)^{3+\frac{1}{1+\omega }}}{1+\gamma \Omega _{\mathrm{m},0}\left[ (1+z)^{3+\frac{1}{1+\omega }}-1\right] } \end{aligned} \end{aligned}$$
for \(z\sim 0\), and
$$\begin{aligned} \begin{aligned} q_0\approx -1+\frac{4+3\omega }{2+2\omega }\Omega _{\mathrm{m},0}, \end{aligned} \end{aligned}$$
for \(z=0\) today. We confine the present study to small deviations from the standard \(\Lambda \)CDM model and hence it is reasonable to expect, in accordance with the most recent Planck results [12], that \( \Omega _{\mathrm{m},0}\sim 0.3\) leading to \(q_0\sim -0.55\) from \(q=-1+\frac{3}{2}\Omega _{\mathrm{m}}\) in the standard \(\Lambda \)CDM, which in turn implies from (58) that we should have \(|\omega | \gtrsim 10\), which leads to \( \omega \gtrsim 10\) since have already assumed \(\omega \ge 0\).9 Unless otherwise, for instance, if \(\omega =0\), then \(q_0\sim -0.4\), which could be decreased to a reasonable value, viz., \(q_0 \sim -0.55\), by decreasing \(\Omega _{\mathrm{m},0}\) considerably, namely \(\Omega _{\mathrm{m},0}\sim 0.22\), which implies very large deviations from the standard \(\Lambda \)CDM. Therefore, our assumption \(\omega \ge 0\) already allows large deviations from the standard \(\Lambda \)CDM, in other words, we do not oversimplify the model but yet basically avoid the cases most likely non-viable that would lead to extended discussions since the model exhibits various complicated dynamics particularly at small negative values of \(\omega \), namely, when \(\omega \sim -1\). Hence, in this section, we shall depict some key/interesting features of the anisotropic BD model and its deviation from the standard \(\Lambda \)CDM model assuming \(\omega \ge 0\).
First of all, one may check easily that, similar to the standard \(\Lambda \)CDM, its anisotropic BD extension asymptotically approaches to de Sitter Universe in the far future,
$$\begin{aligned} H^2\rightarrow \frac{4}{3}\frac{\gamma }{\varphi _0^2}\rho _{\mathrm{M}}=\frac{2\gamma }{3}\omega M^2,\,\,\, \rho _{\sigma ^2}\rightarrow 0\,\,\,\text {and} \,\,\, \rho _{\mathrm{m}}\rightarrow&0 \nonumber \\ \,\,\,\text {for}\,\,\, z\rightarrow&-1, \end{aligned}$$
and accordingly, at this limit, the effective DE mimics exactly the cosmological constant as
$$\begin{aligned} \rho _{\mathrm{DE}}\rightarrow \gamma \rho _{\mathrm{M}},\,\,\, \bar{w}_{\mathrm{DE}}\rightarrow -1\quad \text {and}\quad \Delta w_{\mathrm{DE}}\rightarrow 0. \end{aligned}$$
At low redshifts, viz., at which \(\rho _{\mathrm{m}}\simeq \rho _{\mathrm{m},0}\) and \(\rho _{\sigma ^2}\ll \rho _{\mathrm{m}}\), the effective DE can approximately be described by
$$\begin{aligned} \begin{aligned}&\rho _{\mathrm{DE}}\simeq \gamma \rho _{\mathrm{M}}+(\gamma -1)\rho _\mathrm{m }\quad \text {and}\quad \\&\bar{w}_{\mathrm{DE}}\lesssim -1+\frac{(\gamma -1)\rho _\mathrm{m }}{\gamma \rho _{\mathrm{M}}+(\gamma -1)\rho _\mathrm{m }}\quad \text {for}\quad z\simeq 0. \end{aligned} \end{aligned}$$
Here, we first note that the non-negativity of the effective DE today (\(z=0\)), \(\rho _{\mathrm{DE,0}}\ge 0\), leads to a natural lower bound on the value of the energy density corresponding to the effective cosmological constant (viz., \(\rho _{\mathrm{M}}\)) as \(\rho _{\mathrm{M}}\ge \left( \frac{1}{\gamma }-1\right) \rho _{\mathrm{m},0}\), which allows zero lower bound only at the GR limit, i.e., at \(\omega \rightarrow \infty \) (viz., \(\gamma \rightarrow 1\)). We next note that \(\gamma \) takes values within the range \(\frac{1}{2}\le \gamma < 1\) for \(\omega \ge 0\) and hence the coefficient of \(\rho _{\mathrm{m}}\) is negative definite, \(\gamma -1< 0\), which implies that, because \(\rho _{\mathrm{m}}\) decreases as the Universe expands, \(\rho _{\mathrm{DE}}\) will increase as the Universe expands at \(z\sim 0\), i.e., the effective DE will behave like a phantom fluid (\(\bar{w}_{\mathrm{DE}}\lesssim -1\)) at \(z\sim 0\) [see (56) for \(\bar{w}_{\mathrm{DE}}(z=0)\)]. Namely, we see from (40) or (55) that given \(\rho _{\mathrm{m}}=\rho _{\mathrm{m,0}}(1+z)^3\), the effective DE passes below the phantom-divide-line (PDL, \(w_{\mathrm{DE}}=-1\)) at
$$\begin{aligned} z=z_\mathrm{PDL}\equiv \left( \frac{2\omega +3}{2\omega +2}\right) ^{1+\omega }-1, \end{aligned}$$
and afterwards stays there forever, i.e., \(\bar{w}_{\mathrm{DE}}<-1\) for \(z<z_\mathrm{PDL}\). It is worth noting that the BD gravity predicts almost a specific redshift for the PDL crossing, such that
$$\begin{aligned} \frac{1}{2}\le z_\mathrm{PDL}\le e^{\frac{1}{2}}-1=0.65\quad \text {for}\quad \omega \ge 0, \end{aligned}$$
and \(0.63\le z_\mathrm{PDL}<0.65\) for \(\omega \ge 10\) (see [100] for a similar result). This is indeed an interesting result, such that \(z_\mathrm{PDL}\sim 0.6\) can be taken as the signature of BD gravity. In other words, observations suggesting the presence of a DE source passing below PDL at \(z\sim 0.6\) with high confidence would imply a strong reason for favoring the BD gravity over GR, or vice versa. During this period, \(\Delta w_{\mathrm{DE}}\) can approximately be given by (51) and we note that the effective DE will be almost isotropic, \(\Delta w_{\mathrm{DE}}\approx 0\), for \(\omega >0\) since \(\Omega _{\sigma ^2,0}\ll \Omega _{\mathrm{DE,0}}\) in a realistic cosmology.
At moderate redshifts (e.g., \(z_\mathrm{PDL} \ll z \lesssim z_{\mathrm{eq}}\)), namely, when the contribution from the pressureless matter to the EoS parameter of the effective DE given in (55) becomes significant but those from the expansion anisotropy and the mass of the Jordan field are not, the mean/volumetric EoS parameter of effective DE yields a plateau where
$$\begin{aligned} \begin{aligned} \rho _{\mathrm{DE}}\sim \rho _{\mathrm{m,0}}\gamma \left( \frac{\rho _{\mathrm{m}}}{\rho _{\mathrm{m,0}}}\right) ^{\frac{3\omega +4}{3\omega +3}} \quad \text {as}\quad {\bar{w}}_{\mathrm{DE}}\sim \frac{1}{3\omega +3}\\ \quad \text {for}\quad z_\mathrm{PDL} \ll z \lesssim z_{\mathrm{eq}}. \end{aligned} \end{aligned}$$
This plateau is persistent all the way to the large redshift values at which either the radiation becomes dominant over pressureless matter (\(z\sim z_{\mathrm{eq}}\)), or the contribution from the expansion anisotropy becomes significant. The plateau of \({\bar{w}}_{\mathrm{DE}}\) is located within the range \(0\le {\bar{w}}_{\mathrm{DE}}\le \frac{1}{3}\) depending on the value of \(\omega \ge 0\). For \(\omega =0\), the effective DE behaves like an extra relativistic degree of freedom as \({\bar{w}}_{\mathrm{DE}}=\frac{1}{3}\) during this period and hence requires special attention. On the other hand, as it is discussed above, to obtain a viable model (viz., which deviates from the standard \(\Lambda \)CDM model only slightly), we expect \(\omega \gtrsim 10\) and this places the plateau in a value in the range \(0<{\bar{w}}_{\mathrm{DE}}\lesssim 0.03\). During this period, the evolution of the anisotropy of the effective DE can approximately be given by (52). We note that, although \(\Delta w_{\mathrm{DE}}\) increases as \((1+z)^{\frac{3}{2}}\) during this period, the EoS parameter of the effective DE would not achieve a significant anisotropy during this period even under the weakest observational constraints we give in Sect. 4.
At very/sufficiently large redshifts, viz., when the only dominant contribution is from the expansion anisotropy in (55), the ratio of the energy densities of the effective DE and the expansion anisotropy freezes out;
$$\begin{aligned} \begin{aligned} \rho _{\mathrm{DE}}\sim (\gamma -1)\rho _{\sigma ^2}\quad \text {and}\quad \bar{w}_{\mathrm{DE}}\sim \frac{3\omega +5}{3\omega +3}. \end{aligned} \end{aligned}$$
We note that, except the GR limit \(\omega \rightarrow \infty \) (i.e., \(\gamma \rightarrow 1\)), the mean/volumetric EoS parameter of the effective DE yields a second plateau placing in the range \(1\lesssim \bar{w}_{\mathrm{DE}}\lesssim \frac{5}{3}\) for \(\omega \ge 0\), and the effective DE behaves like the expansion anisotropy with a value of magnitude of ratio in the range \(0<|\rho _{\mathrm{DE}}/\rho _{\sigma ^2}|\lesssim \frac{1}{2}\) but yields negative energy density since the coefficient \(\gamma -1\) is negative for \(\omega \ge 0\). According to this the energy density of the effective DE changes sign and we calculate from (40) that this occurs at
$$\begin{aligned} z=z_\mathrm{DE, pole}\sim \left( \frac{\gamma }{1-\gamma }\,\, \frac{\Omega _{\mathrm{m},0}}{\Omega _{{\sigma ^2},0}}\right) ^ {\frac{\omega +1}{3\omega +4}}-1, \end{aligned}$$
at which its EoS parameter exhibits a pole since at that redshift its energy density becomes zero. For \(z>z_\mathrm{DE, pole}\), \(\Delta w_{\mathrm{DE}}\) asymptotically approaches to (53) with increasing redshift and the EoS of the effective DE will be significantly anisotropic, viz., we have \(1<\bar{w}_{\mathrm{DE}}\le \frac{5}{3}\) and \(-2\le \Delta w_{\mathrm{DE}}<-1.2\) for \(\omega \ge 0\). However, for \(\omega \ge 0\), this interesting behavior is always irrelevant to the observable Universe since it occurs when the Universe is strongly dominated by the expansion anisotropy: Substituting \(z=z_\mathrm{DE, pole}\) into \(\frac{\rho _{\sigma ^2}}{\rho _{\mathrm{m}}}=\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{m,0}}} (1+z)^{3+\frac{2}{1+\omega }}\), which may be read off from (42), we obtain
$$\begin{aligned} \frac{\rho _{\sigma ^2}(z=z_\mathrm{DE,pole})}{\rho _{\mathrm{m}}(z=z_\mathrm{DE,pole})}=\left( \frac{\gamma }{1-\gamma }\right) ^{\frac{3\omega +5}{3\omega +4}}\left( \frac{\Omega _{\mathrm{m,0}}}{\Omega _\mathrm{\sigma ^2,0}}\right) ^{\frac{1}{3\omega +4}}. \end{aligned}$$
We note that, in a realistic cosmological setup, say, \(\omega \gtrsim 10\) (even for \(\omega \ge 0\)) and \(\Omega _{\mathrm{m,0}}\gg \Omega _{\sigma ^2,0}\), this ratio is obviously much larger than unity, implying that \(z_\mathrm{DE,pole}\) occurs always at a redshift at which the expansion anisotropy is already dominant. If we consider the presence of radiation as well, then we have the following: During the radiation domination BD gravity mimics GR (as Jordan field \(\varphi \) freezes out) so that the effective DE will be constant (i.e., mimics cosmological constant) while the expansion anisotropy keeps on growing as \(\rho _{\sigma ^2}\propto (1+z)^6\) as z increases. Eventually, the expansion anisotropy will dominate over radiation and the dynamical effective DE – viz., \(\rho _{\mathrm{DE}}\propto \rho _{\sigma ^2}(\gamma -1)\) as may be seen from (40) by setting \(\rho _{\mathrm{M}}=0=\rho _{\mathrm{m}}\) – will show up again and then its EoS parameter will respectively realize a pole and a plateau that follows that pole as z increases. Thus, in a realistic setup, the pole and the second plateau features of the effective DE occur in the expansion anisotropy dominated universe that takes place before (in terms of time) the radiation dominated universe and hence they are irrelevant to the observable universe.
3.3.2 Modified expansion anisotropy
Almost all of the model dependent observational constraints on the expansion anisotropy in the literature, e.g., [91, 92, 93, 94, 95], consider GR+isotropic sources leading to the steep redshift dependency of the energy density corresponding to the shear scalar as \(\rho _{\sigma ^2} \propto (1+z)^{6}\) and thereby give \(\Omega _{\sigma ^2,0}\lesssim 10^{-21}\). On the other hand, the direct/model independent observational constraints from, e.g., type Ia SNe observations are rather weak, e.g., \(\Omega _{\sigma ^2}\lesssim 10^{-4}\) for \(z\sim 0\) [96, 97]. Such large values may be possible provided that the redshift dependence of \(\rho _{\sigma ^2}\) is modified properly (viz., is made modest), which may be done, in principle, by either introducing an anisotropic source, e.g., anisotropic DE, in GR or considering usual cosmological sources in modified theories of gravity such as BD gravity that leads to an effective anisotropic source.
As we have shown in Sect. 3.1 that switching to the massive BD gravity (1) with \(\omega \ge 0\) results in having more stringent constraints on the expansion anisotropy [see (37) and discussion that follows]. These constraints in fact can be weakened for \(\omega <-1\) (since the expansion anisotropy in this case yields flatter redshift dependency) and even vastly for \(-\frac{5}{3}<\omega <-1\), which however would cost to large deviations from GR, hence from \(\Lambda \)CDM model, such that in this case we will have rapidly varying cosmological gravitational coupling strength [see (24)] or, in the approach realized in (39), to an effective DE significantly deviating from \(\Lambda \) [see (40) along with (27)]. It may be interesting to discuss a bit more on how the expansion anisotropy can be eased down in case \(\omega <-1\) by considering (39) based on the approach that can be described by (8). We see from (65) – or directly from (26) – that the effective EoS parameter corresponding to the expansion anisotropy, viz., shear scalar, is
$$\begin{aligned} w_{\sigma ^2}=1+\frac{2}{3\omega +3}. \end{aligned}$$
We first note that it approaches Zeldovich fluid, i.e., \(w_{\sigma ^2}\rightarrow 1\), at the GR limit \(|\omega |\rightarrow \infty \) and, for \(\omega \ge 0\) (finite), is in the range \(1< w_{\sigma ^2}\le \frac{5}{3}\), which is stiffer than the Zeldovich fluid. Zeldovich fluid yields the most rigid EoS parameter (\(w=1\)) compatible with the requirements of relativity theory [144], which in turn implies that we can not have a minimally interacting source whose energy density redshift dependence is steeper than \((1+z)^{6}\) within GR, whereas we have it in BD gravity from the expansion anisotropy as an effective source provided that \(\omega >-1\) as can be seen from (27). We depict \(w_{\sigma ^2}\) versus \(\omega \) in Fig. 1.
Effective equation of state parameter corresponding to the expansion anisotropy \(w_{\sigma ^2}\) versus Brans–Dicke parameter \(\omega \). \(w_{\sigma ^2}\) exhibits a pole at \(\omega =-1\), and \(w_{\sigma ^2}\rightarrow 1\) as \(|\omega |\rightarrow \infty \)
For \(\omega \le -\frac{4}{3}\), the EoS parameter spanning the range \(-1\le w_{\sigma ^2}<1\), which is familiar from the conventional cosmology: (i) \(\omega =-2\) (\(w_{\sigma ^2}=\frac{1}{3}\)), the expansion anisotropy mimics radiation, viz., plays the role of an extra relativistic species. (ii) \(\omega =-\frac{5}{3}\) (\(w_{\sigma ^2}=0\)), it mimics dust, viz., CDM. (iii) \(\omega =-\frac{3}{2}\) –the scale invariant limit – (\(w_{\sigma ^2}=-\frac{1}{3}\)), it mimics constant negative spatial curvature (hyperboloid space). (iv) \(\omega =-\frac{4}{3}\) (\(w_{\sigma ^2}=-1\)), it mimics \(\Lambda >0\). It might be interesting to remind here that \(-\frac{3}{2}<\omega <-\frac{4}{3}\) corresponds to the range in which the universe exhibits accelerating expansion in the presence of only pressureless matter (see footnote 4). The BD parameter \(\omega \sim -\frac{4}{3}\) appears in d-branes string models [104, 105] and hence it is conceivable that such models would predict \(\sigma ^2\sim \mathrm{const}\). For \(-\frac{4}{3}<\omega <-1\), we have \(w_{\sigma ^2}<-1\), i.e., the expansion anisotropy mimics phantom sources implying that the expansion anisotropy decreases/increases with increasing/decreasing redshift in contrast to the all other cases. The case \(\omega =-1\) (the low energy effective string action limit of BD gravity [104, 105]) is interesting that \(w_{\sigma ^2}\) exhibits a pole, namely, \(\lim _{\omega \, \rightarrow \, -1^+}w_{\sigma ^2}=+\infty \) and \(\lim _{\omega \, \rightarrow \, -1^-}w_{\sigma ^2}=-\infty \), which imply that we must set \(\sigma ^2=0\) in this case (at least in our solution). In relevance with this, in the case of \(\omega >-1\) but \(\omega \approx -1\), \(w_{\sigma ^2}\) will be extremely large, implying that in this case the presence of the expansion anisotropy today or at any moment in the observable Universe will require extreme fine tuning. Finally, for \(\omega =0\) – the lowest value we consider in this study – we have \(\omega _{\sigma ^2}=\frac{5}{3}\) and for \(\omega \gtrsim 10\), which, as we discussed above, is required for small deviations from the \(\Lambda \)CDM model (by keeping isotropic RW metric), we have \(1<w_{\sigma ^2}\lesssim 1.06\). These imply that considering BD gravity with \(w\gtrsim 10\) rather than GR as the law of gravity will additionally lead to small modification in the evolution of the expansion anisotropy (viz., a slightly steeper redshift dependency of the shear scalar compared to the one in GR) in the anisotropic generalization of the \(\Lambda \)CDM model in contrast to DE models that can be described by scalar field/s within GR.
4 Constraints from recent cosmological data
In this section we perform a parameter estimation and provide observational constraints on the free parameters of the model: the pressureless matter density parameter today \(\Omega _{\mathrm{m,0}}\), the baryon density parameter today \(\Omega _\mathrm{b,0}h^2\), and the dimensionless Hubble constant h as well as the two free parameters \(\omega \) and \(\Omega _{\sigma ^2,0}\) that account for the anisotropic BD extension of the \(\Lambda \)CDM model. We consider H(z) given in (34) rather than (39), since the former one contains radiation (viz., \(\Omega _{\mathrm{r,0}}\)) as well and can be safely used all the way to the recombination redshift (covering last scattering surface), so that we could include CMB data in our analysis. For the derived parameters we present the ones relevant to the effective DE. We consider (39), where the effective DE is defined in accordance with the exact explicit solution in redshift given in Sect. 3 and the approach described by (8). In this way, even though this solution ignores the presence of radiation, we could go further and investigate the dynamics of the model under consideration in terms of redshift for \(z<z_{\mathrm{eq}}\) (during which radiation is negligible) analytically in the light of the observational constraints. We keep in mind that H(z) throughout this study is averaged over the volumetric expansion rate of the Universe and hence the method we use constrains the expansion anisotropy through its contribution to the average expansion rate of the Universe (viz., the shear scalar \(\sigma ^2\), which contributes to the H(z) directly via the term \(\rho _{\sigma ^2,0}(1+z)^{6+\frac{2}{1+\omega }}\) and also indirectly via its effect on the \(\rho _{\mathrm{DE}}\) as a result of BD theory).
Constraints on the anisotropic Brans–Dicke extension of the standard \(\Lambda \)CDM model using the combined data sets PLK+BAO+SN+H. For two-tailed distributions, the results are given in \(1\sigma \) and for one-tailed distributions, the results are given in \(2\sigma \). \(\Lambda \)CDM column corresponds to GR limit (Brans–Dicke gravity with \(\omega \rightarrow \infty \)) and the shear scalar \(\sigma ^2=0\) (isotropic). Parameters and ranges of the uniform priors assumed in our analysis and derived parameters are labeled with \(^*\)
Anisotropic BD
\(\Lambda \)CDM
\(\Omega _{\mathrm{m,0}}\)
[0.05, 1.5]
\(\Omega _\mathrm{b,0}h^2\)
0.02241(16)
[0.02, 0.025]
[0.4, 1.0]
\(\log _{10}\omega \)
\(>1.69\)
[\(\omega \rightarrow +\infty \)]
[0, 6]
\(\log _{10}\Omega _{\sigma ^2,0}\)
\(< -8.48\)
\([-12, 0]\)
\(^*{\bar{w}}_{\mathrm{DE,0}}\)
\(>-1.0044\)
\([-1]\)
\(^*\Delta w_{\mathrm{DE,0}}\)
\(<4.23\times 10^{-7}\)
\(^*z_\mathrm{PDL}\)
\(^*\log _{10}z_\mathrm{DE,pole}\)
\(^*z_{\mathrm{eq}}\)
\(3368.62\pm 30.89\)
\(^*2\omega M{^2}\) [\(10^{-66}\)eV\(^2\)]
\(4.48\pm 1.22\)
(\(\Lambda \)) \(4.48 \pm 1.25\)
\(^*M\) [\(10^{-34}\)eV]
\(<1.51\)
\(-2\Delta \ln \mathcal {L}\)
\(-0.62\)
In order to perform the parameter space exploration, we make use of a modified version of the simple and fast Markov Chain Monte Carlo (MCMC) code that computes expansion rates and distances from the Friedmann equation, named SimpleMC [21, 145]. For an extended review of cosmological parameter inference see [146]. The code uses a compressed version of the Planck data (PLK), where the CMB is treated as a "BAO experiment" at redshift \(z=1090\), measuring the angular scale of the sound horizon at that time, a recent analysis of Type Ia supernova (SN) data dubbed Joint Light-curve Analysis (JLA) compressed into a piece-wise linear function fit over 30 bins spaced evenly in \(\log z\), and high-precision Baryon Acoustic Oscillation measurements (BAO), from the comoving angular diameter distance only (viz., growth rate measurement is not considered), the Hubble distance and the volume averaged distance, at different redshifts up to \(z=2.36\). For a more detailed description about the datasets used see [21]. We also include, independent of/separate from other data sets, a collection of currently available data on H(z) obtained from cosmic chronometers (H) (see [143] and references therein). Table 1 displays the parameters used throughout this paper along with the corresponding flat priors; derived parameters labeled with \(^*\).
In the analysis, we assume radiation content by including three neutrino species (\(N_\mathrm{eff}=3.046\)) with minimum allowed mass \(\sum m_{\nu }=0.06\, \mathrm{eV}\) theoretically well determined within the standard model of particle physics and the density parameter of radiation \(\Omega _{\mathrm{r},0}=2.469\times 10^{-5} h^{-2}(1+0.2271 N_\mathrm{eff})\), where dimensionless Hubble constant \(h=H_0/100\, \mathrm{km\,s}^{-1}\,\mathrm{Mpc}^{-1}\) [147]. Throughout the analysis we assume flat priors over our sampling parameters: \(\Omega _{\mathrm{m},0}=[0.05,1.5]\), \(\Omega _{\mathrm{b},0} h^2=[0.02,0.025]\) and \(h=[0.4,1.0]\). Whereas priors for the two parameters that characterize the anisotropic BD extension of the standard \(\Lambda \)CDM model are given by \(\log _{10}\omega =[0,6]\) and \(\log _{10}\Omega _{\sigma ^2,0}= [-12,0]\). The photon energy density today \(\rho _{{\gamma },0}\) is not subject to our analysis since it is well constrained, such that it has a simple \(\rho _{\gamma }=\frac{\pi ^2}{15} T_\mathrm{CMB}^4\) relation with the CMB monopole temperature [148], which is very precisely measured to be \(T_{\mathrm{CMB},0}=2.7255\pm 0.0006\,\mathrm{K}\) [149]. Table 1 summarizes the observational constraints on the free parameters as well as the derived parameters of the model under consideration using the combined datasets PLK+BAO+SN+H; and for comparison those parameters used on the standard \(\Lambda \)CDM model.
(Left panel): 2D marginalized posterior distribution of the parameters \(\log _{10}\Omega _{\sigma ^2}\) and \(\log _{10}\omega \) that determine the extension from the standard \(\Lambda \)CDM model. The 2D constraints are plotted with 1\(\sigma \) and 2\(\sigma \) confidence contours. (Right panel): the 3D posterior distribution in the \(\{\bar{w}_{\mathrm{DE,0}}, ~\Delta w_{\mathrm{DE,0}}, ~\log _{10}\omega \}\) subspace, where the colour code indicates the value of \(\log _{10}\omega \) using the colour bar
(Left panel) 2D marginalised posterior distribution of the parameters corresponding to the mass of the Jordan field M and the Brans–Dicke parameter \(\log _{10}\omega \). (Right panel) 1D constraints corresponding to the redshift of matter-radiation equality (\(z_{\mathrm{eq}}\))
The constraints on g(z) as a result of the data. These show the posterior probability Pr(g|z): the probability of g as normalized in each slice of constant z, with color scale in confidence interval values. The 1\(\sigma \) and 2\(\sigma \) confidence intervals are plotted as black lines and red lines correspond to the best-fit values found over the analysis, while dotted lines to the redshift of matter-radiation equality. (Top) \(\bar{w}_{\mathrm{DE}}(z)\) at low z values (left) and at high redshift values in log-scale (right). (Middle) \(\Delta w_{\mathrm{DE}}(z)\) at low z values (left) and at high redshift values in log-scale (right). (Bottom) \(\Omega _{\mathrm{DE}}(z)\) at low z values (left) and at high redshift values in log-scale (right)
Figure 2 summarizes the constraints on the parameters that characterize the anisotropic BD extension to the standard \(\Lambda \)CDM model. Left panel of this figure displays 2D marginalized posterior distributions on the BD parameter \(\omega \) (deviation from GR) and the density parameter corresponding to the expansion anisotropy today \(\Omega _{\sigma ^2,0}\) (deviation from isotropic expansion). We see from Table 1 that current observations prefer small contributions from the expansion anisotropy \(\log _{10}\Omega _{\sigma ^2,0} < -8.48\) (95% C.L.) with a large BD parameter (small deviation from GR) \(\log _{10}\omega >1.69\) (95% C.L.), signaling that the anisotropic BD extension to the standard \(\Lambda \)CDM model should be studied indeed as a small correction in line with what we have been claiming and basically doing so far. Right panel of Fig. 2 displays 3D posterior distributions in the \(\{\bar{w}_{\mathrm{DE,0}}, \Delta w_{\mathrm{DE,0}}\), \(\log _{10}\omega \}\) subspace: the average/volumetric EoS parameter of the effective DE \(\bar{w}_{\mathrm{DE,0}}\) and the measure of anisotropy of its EoS parameter \(\Delta w_{\mathrm{DE,0}}\) coloured code with \(\log _{10}\omega \). The former two are derived parameters that determine deviations of the effective DE from \(\Lambda \) (\(\bar{w}_{\mathrm{DE}}=-1\) and \(\Delta w_{\mathrm{DE}}=0\)) today and are controlled by \(\omega \). We notice that for large values \(\omega \) (color coded with red) we recover the \(\Lambda \)CDM case at the present time, where \(\Delta w_{\mathrm{DE,0}}=w_{\mathrm{DE},y,0}-w_{\mathrm{DE},x,0} \rightarrow 0\) and \(\bar{w}_{\mathrm{DE,0}} \rightarrow -1\), as expected by our previous analysis. We would like to remind here that in our observational analysis only positive values of \(\omega \) were allowed. The left panel of Fig. 3 shows that current data impose upper bounds on the mass of the Jordan field given by \(M<1.51 \times 10^{-34}\)eV (95% C.L.). The correlation observed on the parameters \(\omega \) and M is an interesting point to note. Even though both parameters have only one tail constraints, the effective cosmological constant \(2\omega M^2\) –a component of the effective DE, see (31) and (40)– is well limited by \(2\omega M{^2}=4.48\pm 1.22\times 10^{-66}\)eV\(^2\) and in agreement with the current value of \(\Lambda =4.48 \pm 1.25\times 10^{-66}\)eV\(^2\), see Table 1. The extra two parameters from the anisotropic BD extension to the \(\Lambda \)CDM model may shift and relax some of the constraints, as it is the case of \(\Omega _{\mathrm{m,0}}=0.3002 \pm 0.0067\), and therefore as a consequence the earlier redshift at matter-radiation equality \(z_{\mathrm{eq}}=3368.62\pm 30.89\), compared to the \(\Lambda \)CDM model with \(z_{\mathrm{eq}} = 3359.45\pm 24.14\) (see right panel of Fig. 3).
We depict, in Fig. 4, the behaviour of the effective DE (top and middle panels) as well as the evolution of its density parameter \(\Omega _{\mathrm{DE}}=\rho _{\mathrm{DE}}/3H^2\) (bottom panel) in redshift, according to the constraints obtained from observational analysis (see Sect. 3.3.1 for the features of the effective DE). Left panel of Fig. 4 displays the behaviors of \(\bar{w}_{\mathrm{DE}}\), \(\Delta w_{\mathrm{DE}}\) and \(\Omega _{\mathrm{DE}}\) for redshifts up to \(z=2\). This panel was drawn by taking random samples from the posterior distribution of the parameter-space and colour coded by its likelihood (bluer colors represent regions of higher probability). Right panel of Fig. 4 displays the extended behaviors of the same functions; from redshift \(z=1\) to \(z=10^6\) for \(\bar{w}_{\mathrm{DE}}\) and from redshift \(z=10\) to \(z=10^6\) for \(\Delta w_{\mathrm{DE}}\) and \(\Omega _{\mathrm{DE}}\). Though, note that the extensions to the redshift values \(z\sim z_{\mathrm{eq}}\) and larger should be seen with a caution that these are obtained by making use of the effective DE derived in Sect. 3.3.1 without considering the presence of radiation. In this panel we depict the probability of a function normalized in each slice of constant z, with colour scale in confidence interval values. The 1\(\sigma \) and 2\(\sigma \) confidence intervals are plotted as black lines and red lines correspond to the best-fit values found over the analysis. At low redshifts, the upper-left panel shows that the crossing of the PDL occurs at \(z_\mathrm{PDL}=0.64839\pm 0.00058\), just as we concluded from the discussion that follows equation (62). For redshift values \(z<z_\mathrm{PDL}\), the effective DE exhibits phantom behaviour \(\bar{w}_{\mathrm{DE}}<-1\) as may be seen from (56) or (61) and for higher redshift values, \(z>z_\mathrm{PDL}\), \(\bar{w}_{\mathrm{DE}}>-1\) as may be seen from (64) and (65). Evolving DE with an EoS parameter being below \(-1\) at present, evolved from \(w>-1\) in the past is named as quintom DE. We obtain the quintom DE as the upper-left panel shows, whereas, the explicit construction of quintom scenario is more difficult than other dynamical DE models, due to a no-go theorem which forbids the EoS parameter of a single perfect fluid or a single scalar field to cross the \(w=-1\) boundary [43]. This property is distinctive that single-scalar-field models with canonical kinetic term are not allowed to satisfy, which also lead to that the Hamiltonian is unbounded from below. Interestingly this property is also achieved throughout some model independent analyses [23, 29, 30]. The middle left panel shows that, as can be deduced from (51) and (52), the effective DE becomes slightly more anisotropic with the increasing redshift, viz., the anisotropy of the EoS parameter \(\Delta w_{\mathrm{DE}}\) increases only about an order of magnitude from its current value \(\Delta w_{\mathrm{DE,0}}\lesssim 10^{-7}\) (see Table 1) to the one at \(z=2\). Looking at the right panel we see that after few redshifts, as pressureless matter becomes dominant (see in the bottom panel that \(\Omega _{\mathrm{DE}}\sim 0.1\) for \(z\sim 1.75\) and \(\Omega _{\mathrm{DE}}\sim 0\) for \(z\gtrsim 10\)), \(\bar{w}_{\mathrm{DE}}\) starts to noticeably climb up and settles in the first plateau of \(\bar{w}_{\mathrm{DE}}\sim 0\) (see (64) and paragraph covering it) lying between \(z\sim 50\) and \(z\sim z_{\mathrm{eq}}\) (viz., throughout the pressureless matter dominated epoch). There is a period in this plateau during which \(\bar{w}_{\mathrm{DE}}>0\) [see (64)], viz., \(\Omega _{\mathrm{DE}}\) increases (i.e., the energy density of the effective DE increases faster than of the pressureless matter) with increasing redshift. However, as can be seen in the bottom right panel, \(\Omega _{\mathrm{DE}}\) during this period can never grow up to considerable values, viz., remains less than a percent at 1\(\sigma \) C.L. and few percents at 2\(\sigma \) C.L.. The anisotropy of the effective DE keeps on increasing with increasing redshift approximately in accordance with (52) during this plateau, but it remains positive definite and insignificant (e.g., \(\Delta w_{\mathrm{DE}}\lesssim 0.07\) at photon decoupling redshift \(z\sim 1100\)) until the effective DE starts to leave this plateau as the expansion anisotropy starts to become dominant. We see that \(\bar{w}_{\mathrm{DE}}\) exhibits a pole at \(\log _{10}z_\mathrm{DE,pole}= 4.80\pm 0.58\) then settles in a new plateau of \(\bar{w}_{\mathrm{DE}}\sim 1\) starting at \(z\sim 10^5\) during which \(\Delta w_{\mathrm{DE}}\) exhibits a similar behavior and settles down at \(\Delta w_{\mathrm{DE}}\sim -1.2\), i.e., the effective DE eventually becomes highly anisotropic. However, note that this last stage of the effective DE starting just after \(z\sim z_{\mathrm{eq}}\) is basically an artifact of omitting radiation in Sect. 3.3.1 to have explicit expression of \(\rho _{\mathrm{DE}}(z)\), strickly speaking, is in fact unlikely to be realized at an observationally relevant past of the Universe (see the discussion starting with equation (65) in Sect. 3.3.1).
Instead, in a realistic setup, the radiation domination should start at \(z=z_{\mathrm{eq}}\) and be maintained all the way to the redshift values at which BBN took place \(z\sim z_{\mathrm{BBN}}\). During which the Jordan field is constant and hence, except an altered cosmological gravitational coupling strength, the Universe evolves exactly the same as in GR (see Sect. 3.1), such that the effective DE mimics \(\Lambda \) (viz., \(\bar{w}_{\mathrm{DE}}= -1\) and \(\Delta w_{\mathrm{DE}}=0 \) implying that it is isotropic and maintains its energy density value at \(z\sim z_{\mathrm{eq}}\) for \(z>z_{\mathrm{eq}}\)) and expansion anisotropy increases as \(\propto (1+z)^6\). However, although the data constrain \(z_\mathrm{DE,pole}\) to be about 1-2 orders of magnitude larger than \(z_{\mathrm{eq}}\), we see in Fig. 4 that, even within the 68% C.L. error region, the abrupt behaviour of the effective DE can start at redshift values less than \(z_{\mathrm{eq}}\), which implies that the expansion anisotropy in this case is too large that it will never allow radiation to be dominant in the Universe, though over pressureless matter. Moreover, even there is a region where abrupt behaviour of the effective DE starts at redshift values larger than \(z_{\mathrm{eq}}\), using the constraint on \(\Omega _{\sigma ^2,0}\) from Table 1 we calculate that the expansion anisotropy dominates over radiation at redshift values smaller than \(z_{\mathrm{BBN}}\), i.e., it will spoil the BBN processes. All these imply that the observational constraint method we performed here is able to put stringent enough constraints on neither \(\Omega _{\sigma ^2,0}\) nor \(\omega \). These two parameters are not independent and hence we should further investigate the upper boundary on \(\Omega _{\sigma ^2,0}\) and lower boundary on \(\omega \), of course, particularly, by using the data providing information about the Universe at \(z\gtrsim z_{\mathrm{eq}}\), such as the peak of the matter-power spectrum relevant to \(z=z_{\mathrm{eq}}\) and BBN relevant to \(z\sim 10^8\).
5 Some further observational consequences and discussions
5.1 Variation of cosmological gravitational coupling strength
The normalized rate of change of the cosmological gravitational coupling strength in our exact solution (neglecting radiation), from (24), reads
$$\begin{aligned} \frac{\dot{G}}{G}=-\frac{H}{1+\omega }, \end{aligned}$$
which is negative definite considering \(\omega \ge 0\) in this study, provided that the Universe is expanding (\(H>0\)). The constraints we found on \(\omega \) and \(H_0\) (see Table 1) predict
$$\begin{aligned} \bigg |\frac{\dot{G}}{G}\bigg | <1.092\times 10^{-12} \, \mathrm{yr}^{-1}\,\, (68\% \,\mathrm{C.L.}) \quad \mathrm{for}\quad z=0, \end{aligned}$$
similar to those given for \(z\sim 0\) from various physical systems in which gravity is not negligible, such as the motion of the bodies of the Solar System, astrophysical and cosmological systems (see [150] for a comprehensive review). Assuming this relation (69) holds all the way to matter-radiation equality, we obtain \(|\dot{G}/G|\lesssim 10^{-8} \, \mathrm{yr}^{-1}\) for \(z\sim z_{\mathrm{eq}}\). On the other hand –given that \(\varphi =\mathrm{const.}\) (and hence \(G=\mathrm{const.}\)) is an attractor solution when radiation is dominant– the redshift dependence of G becomes flatter w.r.t. \(G\propto (1+z)^{\frac{1}{1+\omega }}\) [see (24)] as the radiation becomes more significant w.r.t. pressureless matter and thereby \(|\dot{G}/G|\lesssim 10^{-8} \, \mathrm{yr}^{-1}\) for \(z\sim z_{\mathrm{eq}}\) will never be achieved, but instead G will become almost constant for \(z\sim z_{\mathrm{eq}}\) and remain so for \(z>z_{\mathrm{eq}}\) as long as radiation continues to be dominant. Yet, in accordance with our detailed discussion in Sect. 3.1, (69) is mostly valid all the way to \(z_{\mathrm{eq}}\) and leads to \(G_1\sim G_0(1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }}\) [see (33)]. Accordingly, we find that the constraints from the data predict the upper bound on the relative change in the strength of the cosmological gravitational coupling at \(z\sim z_{\mathrm{eq}}\) w.r.t. its present time value as \(\frac{\Delta G}{G_0}|_{z\sim z_{\mathrm{eq}}}=0.138\,(0.175)\) 68% C.L. (95% C.L.). This, in turn, implies that the strength of the cosmological gravitational coupling during the radiation dominance will be a constant \(G_{1}\) satisfying the following constraints:
$$\begin{aligned} \begin{aligned} G_0\le G_{1}\lesssim 1.138\, G_0\,\,(1.175\, G_0)&\quad 68\% \,\mathrm{C.L.}\, (95\% \,\mathrm{C.L.}) \end{aligned} \end{aligned}$$
for \(z\gtrsim z_{\mathrm{eq}}\). These constraints, of course, are valid also for the epoch of primordial nucleosynthesis that takes place when \(z_{\mathrm{BBN}}\sim 3\times 10^8\), provided that the expansion anisotropy is still insignificant at that redshift. We finally note that these constraints are similar to those obtained from CMB and BBN (see Ref. [150] for a review on CMB and BBN constraints on \(\Delta G/G_0\)).
5.2 Matter-radiation transition
We find using (36) that, along with the constraints on the BD parameter and the radiation density parameter, the upper bound on the expansion anisotropy today, \(\Omega _{\sigma ^2,0}=10^{-8.48}\) (see Table 1), leads to a lower bound on the expansion anisotropy-radiation equality (\(\Omega _{\sigma ^2}=\Omega _{\mathrm{r}}\)) redshift as \(z_\mathrm{eq,\sigma ^2,r}\sim 10^3\), which is obviously not acceptable in a viable cosmological model since it implies that the expansion anisotropy dominates the Universe even at redshift values less than the matter-radiation equality (\(\Omega _{\mathrm{m}}=\Omega _{\mathrm{r}} \)) redshift \(z_\mathrm{eq }\sim 3369\) (see Table 1). Indeed, using the results given in Table 1, we see that the Universe could be dominated by the expansion anisotropy as \(\frac{\Omega _{\sigma ^2,\mathrm eq}}{\Omega _\mathrm{r,eq}+\Omega _\mathrm{m,eq}}=\frac{\Omega _{\sigma ^2,\mathrm eq}}{2\Omega _\mathrm{m,eq}}\lesssim 10^2\) at the matter-radiation equality, and thereby conclude that the upper bound on \(\Omega _{\sigma ^2,0}\) given above should be further reduced to values guaranteeing that the Universe is radiation and matter dominated at the matter-radiation equality. Additionally, transition from radiation to matter domination is one of the most important epochs in the history of the Universe. This transition alters the growth rate of density perturbations: during the radiation era perturbations well inside the horizon are nearly frozen but once matter domination commences, perturbations on all length scales are able to grow by gravitational instability and therefore it sets the maximum of the matter power spectrum in GR as well as BD [138]. Namely, it determines the wavenumber, \(k_{\mathrm{eq}}\), of a mode that enters the horizon, \(H_{\mathrm{eq}}a_{\mathrm{eq}}\), at the matter-radiation transition [138, 151].
In our model, \(k_{\mathrm{eq}}=H_{\mathrm{eq}}a_{\mathrm{eq}}=\frac{H_{\mathrm{eq}}}{(1+z_{\mathrm{eq}})}\) can be estimated analytically by assuming H(z) given in (34) holds all the way to matter-radiation equality. At this point, both the radiation and matter contribute equally to the total energy density. This gives us opportunity to reduce the constraints on \(\Omega _{\sigma ^2,0}\). To do so, using (34), we write
$$\begin{aligned} k_{\mathrm{eq}}=\,&H_0\big [\Omega _\mathrm{M,0}(1+z_{\mathrm{eq}})^{-2}+2\Omega _{\mathrm{m,0}}(1+z_{\mathrm{eq}})^{1+\frac{1}{1+\omega }}\nonumber \\&+\Omega _\mathrm{\sigma ^2,0}(1+z_{\mathrm{eq}})^{4+\frac{2}{1+\omega }}\big ]^{1/2}, \end{aligned}$$
where we use \(a_0=a(z=0)=1\), and \(\Omega _\mathrm{M,0}=1-\Omega _{\mathrm{m,0}}-\Omega _{\mathrm{r,0}}-\Omega _\mathrm{\sigma ^2,0}\). We note that, for \(\omega \gtrsim 50\) and \(z_\mathrm{eq }\sim 3369\) (see Table 1), the contribution to \(k_{\mathrm{eq}}\) from the term with \(\Omega _\mathrm{M,0}\) is negligible, \(k_{\mathrm{eq}}\) is by far more sensitive to \(\Omega _\mathrm{\sigma ^2,0}\) and, for a given set of values of the parameters, \(k_{\mathrm{eq}}\) increases with increasing \(\Omega _\mathrm{\sigma ^2,0}\) and decreases with increasing \(\omega \) (viz., as the BD gravity approaches to GR).
We obtain \(k_{\mathrm{eq}}=0.01034 \pm 0.00012\) for \(\Omega _{\sigma ^2,0}=0\) (isotropic expansion), which is slightly larger than, but yet consistent with the isotropic \(\Lambda \)CDM model value \(k_{\mathrm{eq}}=0.01024 \pm 0.00007\). Accordingly, switching to the massive BD (1) leads only to a slight increase in \(k_{\mathrm{eq}}\), as expected from the constraint \(\omega \gtrsim 50\) (see Table 1). It may be useful to note that these two values are consistent with the recent Planck release [12] value \(k_{\mathrm{eq}}=0.010339\pm 0.000063\) (TT,TE,EE+lowE+lensing+BAO) obtained for the base \(\Lambda \)CDM model. On the other hand, we see that the expansion anisotropy, although negligible today, shifts \(k_{\mathrm{eq}}\) to unrealistically large values, viz., we obtain \(k_{\mathrm{eq}}= 0.15159 \pm 0.00327\) for the upper bound \(\Omega _{\sigma ^2,0}=10^{-8.48}\) (see Table 1). We then work out the values of \(\Omega _{\sigma ^2,0}\) that can shift \(k_{\mathrm{eq}}\) to reasonable values. See Fig. 5 showing explicitly \(k_{\mathrm{eq}}\) with respect to \(\Omega _{\sigma ^2,0}\) with errors to make a simple comparison between the models. We obtain \(k_{\mathrm{eq}}=0.02824 \pm 0.00057\) by setting \(\Omega _{\sigma ^2,0}=10^{-10}\), and \(k_{\mathrm{eq}}=0.01067 \pm 0.00013\) by setting \(\Omega _{\sigma ^2,0}=10^{-12}\), which are still inconsistent with the \(k_{\mathrm{eq}}\) values given for the isotropic \(\Lambda \)CDM model. But then, we obtain \(k_{\mathrm{eq}}=0.01037 \pm 0.00012\) by setting \(\Omega _{\sigma ^2,0}=10^{-13}\), and \(k_{\mathrm{eq}}=0.01034 \pm 0.00012\) by setting \(\Omega _{\sigma ^2,0}=10^{-14}\), where we notice that only the last decimals are different. We observe further that \(k_{\mathrm{eq}}\) does not change for \(\Omega _{\sigma ^2,0}\lesssim 10^{-14}\) anymore and remains consistent with the isotropic \(\Lambda \)CDM model values obtained in this study as well as recent Planck release, which in turn implies that we cannot distinguish \(\Omega _{\sigma ^2,0}\lesssim 10^{-14}\) from \(\Omega _{\sigma ^2,0}=0\) by means of \(k_{\mathrm{eq}}\). We finally calculate using \(\Omega _{\sigma ^2,0}\sim 10^{-14}\) that the Universe is indeed dominated by matter+radiation at matter-radiation equality, viz., we now have \(\frac{\Omega _{\sigma ^2,\mathrm eq}}{2\Omega _\mathrm{m,eq}}\sim 10^{-3}\). Thus, by means of matter-radiation transition, we conclude
$$\begin{aligned} \Omega _{\sigma ^2,0}\lesssim 10^{-14}\quad \text {along with}\quad \omega \gtrsim 50, \end{aligned}$$
where the upper bound on the expansion anisotropy is improved by reducing about six orders of magnitude w.r.t. the ones given in Table 1.
The wavenumber (\(k_{\mathrm{eq}}\)) of a mode that enters the horizon (\(H_{\mathrm{eq}}a_{\mathrm{eq}}\)) at matter-radiation transition (\(z_{\mathrm{eq}}\)), viz., the maximum of matter power spectrum, for some fixed values of the expansion anisotropy today (\(\Omega _{\sigma ^2,0}\)). Anisotropic BD extension of the \(\Lambda \)CDM model, is labeled by green color, the standard base \(\Lambda \)CDM (in this study) is labeled by blue color, while the recent Planck release value for the standard base \(\Lambda \)CDM model is shown with red color
5.3 Big Bang Nucleosynthesis
Big Bang Nucleosynthesis provides a probe of the dynamics of the early Universe, which in turn gives us opportunity to further investigate the constraints on the anisotropic BD extension of the standard \(\Lambda \)CDM model. Such that, in the standard-BBN (SBBN) –assuming the standard model of particle physics is valid and the expansion of the Universe is governed by GR– the processes relevant to BBN take place when the temperature ranges from \(T\sim 1\, \mathrm{MeV}\) to \(T\sim 0.1\, \mathrm{MeV}\) and the age of the Universe from \(t\sim 1\,\mathrm{s}\) to \(t\sim 3\,\mathrm{min}\) corresponding to redshift \(z\sim 3\times 10^8\) at which the Universe is radiation dominated. These, of course, should not be altered significantly in a viable cosmological model. Accordingly, we first looked for the condition of radiation dominance at \(z\sim 3\times 10^8\) (implying \(z_\mathrm{eq,\sigma ^2,r}\gtrsim 3\times 10^8\)) by using the constraints on the relevant parameters from Table 1 and found out that \(\Omega _{\sigma ^2,0}\lesssim 10^{-21}\) for \(\omega \gtrsim 50\) in line with our preliminary investigation in Sect. 3.1 [cf., see Eq. (37)]. On the other hand, \(z_\mathrm{eq,\sigma ^2,r}\) corresponds to the redshift at which \(\rho _{\sigma ^2}/\rho _{\mathrm{r}}=1\), but as may be seen from the investigations in [91, 152] the expansion anisotropy does not lead to a considerable deviation from the SBBN for the \(\rho _{\sigma ^2}/\rho _{\mathrm{r}}\) ratio up to a few percents, viz.,
$$\begin{aligned} \frac{\rho _{\sigma ^2}(z=z_{\mathrm{BBN}})}{\rho _{\mathrm{r}}(z=z_{\mathrm{BBN}})}\lesssim 10^{-2}. \end{aligned}$$
Hence, considering this condition, we can obtain a new constraint on \(\Omega _{\sigma ^2,0}\) stronger than the one given in (37) obtained from the condition \(z_\mathrm{eq,\sigma ^2,r}\gtrsim 3\times 10^8\). To work this out, we first write, from (35),
$$\begin{aligned} \frac{\rho _{\sigma ^2}}{\rho _{\mathrm{r}}}&\sim \frac{\rho _{\sigma ^2,0}(1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }}(1+z)^6}{\rho _{\mathrm{r,0}}(1+z)^4}\nonumber \\&=\frac{\Omega _{\sigma ^2,0}}{\Omega _{\mathrm{m,0}}}(1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }+1}(1+z)^2, \end{aligned}$$
where we also used \(\frac{\Omega _{\mathrm{m,0}}}{\Omega _{\mathrm{r,0}}}=\frac{\rho _{\mathrm{m,0}}}{\rho _{\mathrm{r,0}}}=1+z_{\mathrm{eq}}\). Then, we use in this equation the condition (74) along with the constraints on \(\Omega _{\mathrm{m,0}}\) and \(z_{\mathrm{eq}}\) given in Table 1 by setting \(z=3\times 10^8\) and reach the following constraint
$$\begin{aligned} \Omega _{\sigma ^2,0}\lesssim 10^{-23} \quad \text {along with}\quad \omega \gtrsim 50. \end{aligned}$$
We note that, here, the contribution from the effective anisotropic pressure of the Jordan field during \(z<z_{\mathrm{eq}}\), which is encoded in the term \((1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }}\) in (75), is not significant, such that it leads to only fifteen percent smaller upper bound value for \(\Omega _{\sigma ^2,0}\), viz., \((1+z_{\mathrm{eq}})^{-\frac{1}{1+\omega }}\lesssim 10^{-0.07}=0.85\) for \(\omega \gtrsim 50\). Thus, provided that the condition on the expansion anisotropy today given in (76) is satisfied, the expansion anisotropy will not have considerable effect on BBN, but the altered strength of the cosmological gravitation coupling depending on the BD parameter \(\omega \) will, as we shall investigate in what follows.
Provided that the condition (76) is satisfied, the Universe is radiation dominated during the BBN epoch and hence, as it is discussed in Sect. 3.1, we have \(a\propto t^{\frac{1}{2}}\) as in the SBBN based on GR except that the strength of the cosmological gravitational coupling during this epoch, \(G_1\), can be slightly larger than its present time value, \(G_0\), in our model based on BD gravity, see (71) and (33). Consequently, in accordance with (33), we can write the expansion rate of the Universe during the radiation dominance, hence during the BBN as well, as follows:
$$\begin{aligned} H^2=\frac{G_1 }{G_0}H^2_\mathrm{SBBN}\quad \text {with}\quad H^2_\mathrm{SBBN}=\frac{8\pi G_0 }{3}\frac{\pi ^2}{30}g_{*}T^4, \end{aligned}$$
where T is the temperature and \(g_{*}(T)\) is the effective number of degrees of freedom counting the number of relativistic particle species determining the energy density in radiation as \(\rho _\mathrm{r}=\frac{\pi ^2}{30}g_{*}T^4\). According to this, BD gravity can lead to a larger expansion rate for a given temperature since \(G_1>G_0\) is allowed, see (71). We note that this is analogue of altering the expansion rate of the Universe during BBN by modifying \(g_*\) (for instance, by introducing an extra massless degree of freedom such as sterile neutrino) within GR as
$$\begin{aligned} H^2=\frac{\tilde{g}_* }{g_*}H^2_\mathrm{SBBN} \quad \text {with}\quad H^2_\mathrm{SBBN}=\frac{8\pi G_0 }{3}\frac{\pi ^2}{30}g_{*}T^4. \end{aligned}$$
It is clear from (77) and (78) that we can set the following relation \(\tilde{g}_*=\frac{G_1}{G_0} g_{*}\) implying that the altered G, i.e., \(G_1\), in the BD gravity can equivalently be treated as modified \(g_*\), namely, \(\tilde{g}_*\), in GR. The effect of altered expansion rate of the Universe during BBN due to the modified \(g_*\) within GR is well investigated in the literature [153, 154, 155] and is parametrized in terms of \(S\equiv \frac{H}{H_\mathrm{SBBN}}=\sqrt{\frac{\tilde{g}_*}{g_*}}\), which can be adopted as \(S=\sqrt{\frac{G_1}{G_0}}\) for our work within BD gravity by preserving \(g_*\) as in the SBBN. Such that, deviations from \(S=1\) (SBBN) will modify the neutron abundance and the time available for nuclear production/destruction, changing the BBN-predicted primordial element abundances. In general, it is necessary to access to a BBN code to study the BBN-predicted primordial abundances (viz., mass fractions) as functions of S and \(\eta _{10}\) (the number density of baryons \(n_\mathrm{b}\) to the number density of CMB photons \(n_\mathrm{\gamma }\) defined as \(\eta _{10}=10^{10} n_\mathrm{b}/n_\mathrm{\gamma }\) ). On the other hand, luckily, they have been identified (e.g., in [153, 154, 155]) by extremely simple but quite accurate analytic fits over a limited range in these variables as \(5.5\lesssim \eta _{10}\lesssim 6.5\) and \(0.85\lesssim S \lesssim 1.15\) and, for \(^4\)He (helium) and D (deuterium) these are10 [155]
$$\begin{aligned} Y_{\mathrm{p}}&=0.2381\pm 0.0006+0.0016[\eta _{10}+100 (S-1)], \end{aligned}$$
$$\begin{aligned} y_{\mathrm{DP}}&=45.7(1\pm 0.06)[\eta _{10}-6(S-1)]^{-1.6}, \end{aligned}$$
where, \(Y_{\mathrm{p}}\equiv \frac{n_\mathrm{He}}{n_\mathrm{b}}\) and \(y_{\mathrm{DP}}\equiv 10^5\frac{n_\mathrm{D}}{n_\mathrm{H}}\) are \(^4\)He and D mass fractions, respectively, and here in our study, \(\eta _{10}=273.9\,\Omega _\mathrm{b,0} h^2\) and \( S=\sqrt{\frac{G_1}{G_0}}\simeq (1+z_{\mathrm{eq}})^{\frac{1}{2(1+\omega )}}\), see (33). We see that \(\eta _{10}\approx 6.1\) for both models from the constraints on \(\Omega _\mathrm{b,0} h^2\) given in Table 1 and find that the cosmological data we considered constrain the altered expansion rate of the Universe for \(z>z_{\mathrm{eq}}\), hence during BBN, as \(1<S<1.0665\) at \(68\% \mathrm{\, C.L.}\) and \(1<S<1.0839\) at \(95\% \mathrm{\, C.L.}\). These are well inside the validity intervals of (79) and (80) and hence they can be safely utilized here for obtaining the BBN-predicted \(Y_{\mathrm{p}}\) and \(y_{\mathrm{DP}}\) values, which in turn can be used for a further investigation of the constraint on the BD parameter \(\omega \).
In the standard \(\Lambda \)CDM model [\(S=1\) (GR), \(\Omega _{\sigma ^2,0}=0\) and accommodating SBBN], considering the constraint on \(\Omega _\mathrm{b,0}h^2\) given in Table 1, we obtain \(\eta _{10}=6.1469 \pm 0.0404\) leading, from (79) and (80), to
$$\begin{aligned} Y_{\mathrm{P}}^\mathrm{SBBN}&=0.24793\pm 0.00065\quad (68\%\, \mathrm{C.L.}), \end{aligned}$$
$$\begin{aligned} y_{\mathrm{DP}}^\mathrm{SBBN}&=2.4999 \pm 0.1642\quad (68\%\, \mathrm{C.L.}). \end{aligned}$$
On the other hand, in the BD extension of the standard \(\Lambda \)CDM model, by assuming \(\Omega _{\sigma ^2,0}\lesssim 10^{-23}\) [viz., expansion anisotropy is negligible during BBN, see (76)] and considering the constraints on \(\Omega _\mathrm{b,0}h^2\), \(\omega \) and \(z_{\mathrm{eq}}\) given in Table 1, we obtain \(\eta _{10}=6.1373\pm 0.0448\) and \(1<S<1.0665\) at \(68\% \mathrm{\, C.L.}\) leading, from (79) and (80), to
$$\begin{aligned} Y_{\mathrm{p}}&=0.24898\pm 0.00199\quad (68\%\, \mathrm{C.L.}), \end{aligned}$$
$$\begin{aligned} y_{\mathrm{DP}}&=2.5338\pm 0.1726\quad (68\%\, \mathrm{C.L.}). \end{aligned}$$
The 2D constraints on the BBN predicted \(^4\)He (helium) mass fraction \(Y_{\mathrm{p}}\) (top panel) and D (deuterium) \(y_{\mathrm{DP}}\) mass fraction (bottom panel) in BD gravity are plotted with 1\(\sigma \) (dark green) and 2\(\sigma \) (light green) confidence contours. Grey bands are for the SBBN-predicted mass fractions with 1\(\sigma \). Red bands are independent observational estimates with 1\(\sigma \) from [157, 158]
Observations of helium and hydrogen recombination lines from metal-poor extragalactic H II regions provide an independent method (viz., a direct measurement) for determining the primordial helium abundance and a latest and widely accepted estimate comes from the data compilations of [157] giving \(Y_{\mathrm{p}}=0.2449 \pm 0.0040\) (68 \(\%\) C.L.). Similarly, the most recent estimate of the primordial deuterium abundance comes from the best seven measurements in metal-poor damped Lyman-\(\alpha \) systems studied in [158] giving \(y_{\mathrm{DP}}=2.527\pm 0.030\) (68 \(\%\) C.L.). We note that the BBN-predicted \(Y_{\mathrm{p}}\) and \(y_{\mathrm{DP}}\) for both the \(\Lambda \)CDM model (81) and its anisotropic BD extension (83) obtained by using the cosmological data, led to consistent values with independent observational estimates. We notice, however, slightly larger mean values in the case of BD gravity, as suggested by (79) and (80) when \(S>1\) (assuming \(\eta _{10}\) is fixed). We give a summary of our findings by depicting the 2D marginalized posterior distribution of the BBN-predicted \(Y_{\mathrm{P}}\) via (79)/\(y_{\mathrm{DP}}\) via (80) and \(\omega \) in the top panel/bottom panel of Fig. 6. In which, for a comparison, we depict also the bands of the SBBN-predicted \(Y_{\mathrm{P}}\) via (79)/\(y_{\mathrm{DP}}\) via (80) for the \(\Lambda \)CDM model accommodating SBBN (\(S=1\)) and of their independent observational estimates given in [157, 158]. We note that there is an increasing anti-correlation between \(Y_{\mathrm{P}}\) and \(\omega \) with decreasing \(\omega \) as it approaches its lower bound, whereas for large values of \(\omega \) the \(Y_{\mathrm{P}}-\omega \) contour at 68% C.L. (dark green) approaches the band of the SBNN-predicted \(Y_{\mathrm{P}}\) at 68% C.L. (grey) as it should be. Besides these, more importantly, we see that the \(Y_{\mathrm{P}}-\omega \) contour at 68% C.L. (dark green) stays above the band of the independent observational estimate band (red) for \(\omega \lesssim 250\), namely, its consistency with the independent observational estimate of \(Y_{\mathrm{p}}=0.2449 \pm 0.0040\) from [157] requires \(\omega \gtrsim 250\) at 68% C.L. providing us an improved constraint on the BD parameter \(\omega \). We see that the BBN-predicted \(y_{\mathrm{DP}}\) via (80) is insensitive to \(\omega \), \(y_{\mathrm{DP}}-\omega \) contour at 68% C.L. (dark green) already covers the independent observational estimate band (red) for \(\omega \gtrsim 250\) (the improved constraint from our \(Y_{\mathrm{P}}\) investigation) but is much wider than it (due to the relatively large internal error in (80), viz., \(y_{\mathrm{DP}}\propto 1\pm 0.06\)), and \(y_{\mathrm{DP}}-\omega \) contour for BD gravity (dark green) approaches the SBBN-predicted \(y_{\mathrm{DP}}\) band at 68% C.L. (grey) for large \(\omega \) values as it should be. These show that we are not able to deduce a new constraint on the BD parameter \(\omega \) by comparing the BBN-predicted value of \(y_{\mathrm{DP}}\) via (80) with its independent observational estimate.
Thus, by means of the BBN, we reach the most stringent constraints on \(\Omega _{\sigma ^2,0}\) and \(\omega \), the two free parameters that determine the anisotropic BD extension to the standard \(\Lambda \)CDM model, as follows;
$$\begin{aligned} \Omega _{\sigma ^2,0}\lesssim 10^{-23}\quad \text {and}\quad \omega \gtrsim 250. \end{aligned}$$
Finally, considering these improved constraints, we find that the contribution to this constraint of \(\Omega _{\sigma ^2,0}\) from the anisotropy of the effective pressure of the Jordan field during \(z<z_{\mathrm{eq}}\), which is encoded in the term \((1+z_{\mathrm{eq}})^{\frac{1}{1+\omega }}\) in (75), is quite insignificant, such that, it leads to only a few percent smaller upper bound value for \(\Omega _{\sigma ^2,0}\), viz., \((1+z_{\mathrm{eq}})^{-\frac{1}{1+\omega }}\lesssim 10^{-0.01}=0.977\) for \(\omega \gtrsim 250\).
We have carried out an explicit detailed theoretical and observational investigation of an anisotropic massive Brans–Dicke (BD) gravity extension of the standard \(\Lambda \)CDM model, wherein the extension is characterized by the two additional degrees of freedom: the so called BD parameter \(\omega \) and the present day value of the density parameter corresponding to the shear scalar (a measure of the expansion anisotropy), \(\Omega _{\sigma ^2,0}\). The role of the cosmological constant \(\Lambda \) is taken over by the Jordan field potential of the form \(U(\varphi )=\frac{1}{2}M^2\varphi ^2\) with M being the bare mass of the field. We have considered the LRS Bianchi type I metric, which generalizes the spatially flat RW metric simply by allowing a different scale factor along one of the three orthogonal axes, while preserving the spatial homogeneity and flatness, and the isotropic spatial curvature.11 We have solved the field equations analytically and obtained the average Hubble parameter, H(z), explicitly by extending the method developed in [100]; described the anisotropic effective DE in accordance with the way of defining effective source given in [45, 46], and consistently included the radiation into the model.
The BD parameter \(\omega \), being the measure of the deviation from GR (\(|\omega |\rightarrow \infty \)), by alone characterizes the dynamical behaviour of the effective DE as well as the redshift dependency of the expansion anisotropy. These two affect each other depending on \(\omega \), such that the shear scalar contributes to the dynamics of the effective DE, and its anisotropic stress controls the dynamics of the shear scalar, in particular deviations from its usual form \(\sigma ^2\propto (1+z)^6\) in GR, see H(z) in (34). We planned the current study as an extension in the sense of a correction to the standard \(\Lambda \)CDM model, so we have mainly confined our investigations to small deviations from this model via sufficiently small \(\Omega _{\sigma ^2,0}\) and large \(\omega \) values. We have shown through some preliminary cosmological discussions that \(|\omega |\lesssim 10\) (roughly) cannot be called as a small deviation, yet, for completeness we have extended our investigations to non-negative values of \(\omega \). Indeed, considering the combined data sets PLK+BAO+SN+H, we have obtained \(\omega \gtrsim 50\), \(M<1.51\times 10^{-34}\mathrm eV\) and \(\Omega _{\sigma ^2,0}\lesssim 10^{-8.48}\) (\(\Omega _{\sigma ^2,0}\lesssim 10^{-14}\), when matter-radiation equality is considered).12 Then by means of BBN, we have improved these constraints to \(\omega \gtrsim 250\) (in particular, from the comparison of the helium abundance prediction of the model with the direct measurements) and \(\Omega _{\sigma ^2,0}\lesssim 10^{-23}\). We have also found that the contribution of the anisotropy of the effective DE (viz., \(\Delta w_{\mathrm{DE,0}}<4.23\times 10^{-7}\)) to this constraint on \(\Omega _{\sigma ^2,0}\) is insignificant, namely, led to only a few percent smaller (stronger) upper bound value.13
All these have led us to conclude that, with the observations relevant to the dynamics of the Universe, the simplest anisotropic massive BD gravity extension of the standard \(\Lambda \)CDM model presents no significant deviations from it all the way to the BBN. The strongest cosmological constraints on \(\omega \) present in the literature, e.g., \(\omega \gtrsim 890\) [116], would obviously just strengthen this conclusion. Moreover, we should further consider the strongest local constraint \(\omega >40000\) [114] as well since the constraint \(M<1.51\times 10^{-34}\,\mathrm eV\) obtained here for the mass of the Jordan field leads to no relaxation on it [115]. This in turn implies that, against the local constraints, the cosmological features that arise from replacing GR by massive BD (such as the modified expansion anisotropy due to the anisotropic effective DE) are not significantly sensitive to the current cosmological observations. In other words, in view of the local constraints on \(\omega \), the simplest anisotropic massive BD gravity extension of the standard \(\Lambda \)CDM model cannot be distinguished from its simplest GR based anisotropic extension [95] when we consider the data from cosmological observations.14
It is conceivable that our findings, when we consider only the observations relevant to the background dynamics of the Universe, are representative for the ones that would be obtained in the more general extensions, namely, through a gravity theory more general than the massive BD and therefore can evade the local constraints. For instance, promoting \(\omega \) to be some functions of \(\varphi \), the observations from the Cassini spacecraft place upon strong constraint on \(\omega \) as in the original BD theory. This constraint, however, now only applies to the local value of \(\omega (\varphi )\), i.e., of the present day value of \(\varphi \) in the Solar System. Similarly, it is possible that \(\omega \) has spatial variation, so that, for instance, it could be constrained to be \(\omega >40000\) locally, but could be in line with \(\omega \sim 10^2\) at cosmological scales [159, 160]. See, for instance, [37] and references therein for a further reading on such gravity theories closely related to the BD theory, of which the Horndeski/Galileon [124, 125, 126] theories are the most popular ones. In the anisotropic extensions through such gravity theories, however, one should deal with more complicated field equations, those may not be solved analytically and lead to the explicit expression of H(z), and hence the solutions and findings, presented in this work, are useful as they may be considered as the good approximations to such constructions. Consequently, the theoretical investigations we carried out here are, in general, instructive for the anisotropic extensions of the standard \(\Lambda \)CDM model replacing GR by a modified theory of gravity approximating massive BD at cosmological scales. When the observational constraints are considered as well, our findings against the cosmological constraints (the local constraints) with \(\omega \gtrsim 10^{2}\) (\(\omega \gtrsim 40000\)) are instructive for the extensions through those gravity theories that can (cannot) evade local constraints. Yet, this also implies that the lesson learned from this study based on the massive BD would approximately be valid for such more general constructions, namely, many interesting features of such anisotropic extensions of the \(\Lambda \)CDM model that do not exist within its simple GR based anisotropic extension, would remain insignificant when the model is constrained by the observations.
We note that, in fact, these features become quite significant for \(|\omega |\sim O(1)\), which, interestingly, is in principle the natural order of magnitude for the BD parameter (as, e.g., \(\omega =-1\) appears in the low energy limit of string theories). This by alone, makes it quite interesting, at least theoretically, to further study this region in spite of that, as we have discussed above, it implies large deviations from the standard \(\Lambda \)CDM dynamics, and hence is not expected to be consistent with the real Universe. Such a study, however, is beyond the scope of the current paper, but yet, here, we would like to make a couple of relevant comments. In case of small negative values of \(\omega \), say, \(\omega \sim -1\) but larger than the scale invariant limit \(-\frac{3}{2}\) (below which stability issues appear) the model exhibits several critical points that complicates the investigation in this region. So it is necessary to carry out several separate analyses within this small region. Meanwhile, this region is distinguished that it is relevant to string theories, namely, \(\omega =-1\) corresponds to the low energy effective string action and \(\omega \sim -1\) appears in d-brane constructions [102, 103, 104, 105]. For these reasons in fact, although we have devoted the current work to non-negative \(\omega \) values, we carried out a discussion in Sect. 3.3.2 showing how dramatically the redshift dependence of the expansion anisotropy can alter when \(\omega \sim -1\). For instance, for \(\omega =-\frac{4}{3}\), expansion anisotropy becomes non-dynamical (mimicking \(\Lambda \)) and for \(\omega =-1\), expansion anisotropy should be set to zero. The region \(-\frac{3}{2}<\omega <-\frac{4}{3}\) is also interesting that, in the presence of only dust (without \(\Lambda \) or the mass of Jordan field), it is the region of BD theory leading to accelerated expansion of the Universe, and we notice that the expansion anisotropy also contributes to this acceleration since it also behaves like DE in this region. We do not know, so far, a concrete and successful mechanism making such small values (negative or positive) of \(\omega \) consistent with local experiments and cosmological observations. Yet, involving extra spatial dimensions may be good place for looking for such possibilities as it was shown that it is possible to make BD theory (massive or massless) consistent with the gravitational tests (including solar system tests) for \(|\omega |=O(1)\) in the presence of extra dimensions [127]. Noticing the presence of extra dimensions as part of string theories and the dynamical extra dimensions (internal space) contribute to the evolution of the four dimensional anisotropic spacetime like an anisotropic stress in a similar way that the Jordan field does in BD theory make this option more appealing and interesting [161, 162, 163, 164, 165]. This shows one of the many possibilities of rich behaviors which could be worthwhile to investigate in future studies on the extensions of the standard \(\Lambda \)CDM model wherein the expansion anisotropy exhibits non-trivial behaviors that may have interesting cosmological consequences.
A more rigorous observational investigation of the model may be another direction to follow. We have studied the observational constraints on the parameters of the model by considering the background dynamics of the Universe through the evolution of the comoving volume scale factor. These constraints can be improved when we consider the perturbation sector – as in [116, 117, 118] providing the most stringent cosmological constraints on \(\omega \)– along with, for instance, the full Planck data – as in [93, 94] providing constraints on the anisotropic expansion on the top of the \(\Lambda \)CDM at a level close to the ones from BBN by considering the CMB radiation temperature and polarization data-. Such a further observational investigation of the model requires considerable amount of work beyond the aim of the current paper, as it is necessary to study perturbations on the anisotropic background in the presence of an effective anisotropic source. This would strengthen the constraints on the parameters of the model, but would probably not change our conclusion that the model exhibits no significant deviations from the standard \(\Lambda \)CDM model all the way to the BBN. Yet, we cannot be sure about that unless we undertake this work in future and see the results.
The converse is not true, even if it turns out that DE exactly yields \(w = -1\).
See [48] for a list of well know anisotropic stresses and their effects on the expansion anisotropy.
We consider the bare mass as it is done in [100], where it formally appears as the mass of a minimally coupled canonical scalar field when curvature scalar is dropped. On the other hand, when the effective mass of the Jordan field is considered there are various definitions in the literature (see, e.g., [130] for a recent discussion and further references), though one cannot say that these would not fail from a strict particle physics point of view.
If this is not the case then \(2\omega M^2=\Lambda \) will contribute to the field equations like a negative cosmological constant, which may be compensated by a rapidly changing Jordan field. Within the standard BD gravity (i.e., massless BD gravity) in the presence of pressureless source, there exist cosmological solutions with accelerating expansion if \(-2<\omega <-1\), and Jordan field is real (i.e., the effective cosmological gravitational coupling is positive definite) as well if \(-\frac{3}{2}<\omega <-\frac{4}{3}\) leading to deceleration parameter in the range \(0>q>-1\), and to \(\varphi ^2\propto (1+z)^2\) and \(\varphi ^2\propto (1+z)^3\) in the two boundaries of this range, correspondingly [131]. Hence, obviously, if we consider \(\omega <0\) and \(M^2>0\) leading to \(2\omega M^2<0\) (i.e., effectively negative cosmological constant), then it will be necessary to confine ourselves to the range \(-\frac{3}{2}<\omega <-\frac{4}{3}\), and moreover the decelerating effect of the term \(2\omega M^2<0\) would be compensated by bringing the value of the BD parameter closer to \(-4/3\), which in turn, implies a faster Jordan field [viz., for \(\omega \sim -4/3\) Jordan field changes as fast as the pressureless source, i.e., \(\varphi ^2 \sim (1+z)^3\)] whereas we are looking for slowly changing Jordan field, say, \(\varphi ^2\sim \mathrm{const.}\) since we demand it to do only small modifications on the \(\Lambda \)CDM dynamics.
Note however that for bound systems in the quasi-static regime, e.g., our solar system, the effective Newton's constant is \(G_\mathrm{N}=\frac{2\omega +3}{2\omega +4}\,G_0\) [110], thus observers in a bound system which formed today would measure the same cosmological and local gravitational strength if \(\varphi ^2=\frac{2\omega +4}{2\omega +3}\varphi _0^2\).
We note that in the GR limit (say, \(\varphi =\mathrm{const.}\)) the anisotropic model under consideration here, (4)–(7), reduces mathematically to the standard \(\Lambda \)CDM+stiff matter model in the RW framework [133]. However, the role of the stiff matter with positive energy density in [133] is played by the shear scalar here, so that one can straightforwardly utilize the rich class of solutions given in [133] for the GR limit in our model. This presents a good example that one could find various solutions of the model under consideration here but yet we focus on the solution obtained by extending the method given in [100] to LRS Bianchi type I metric.
We followed the method given in [100] for obtaining an exact cosmological solution of BD gravity theory for the expansion rate H(z), with its \(\Lambda \)CDM model counterpart up to a large redshift, viz., for pressureless matter in a spatially flat, homogeneous and isotropic universe. A more general exact solution of the same setup was given in [134] (much earlier than [100]), was nevertheless given for the scale factor in cosmic time and very complicated for extracting an exact H(z) required for observational analyses.
It may be interesting to note here that for the case \(M=0\) (i.e., the massless Jordan field case), \(\Omega _{\mathrm{DE},0}\) is a function of \(\omega \) only, and if we set \(\Omega _{\mathrm{DE},0}=0.7\) we find that \(\omega =-\frac{5}{4}\) from (45) and \({\bar{w}}_{\mathrm{DE,0}}\approx -1.29\) from (56) assuming the contribution from the expansion anisotropy is insignificant (presence of expansion anisotropy shifts it to even larger negative values). This (\(\omega =-\frac{5}{4}\)), in turn, implies from (26) that the shear scalar (viz., expansion anisotropy) contributes to the modified Einstein field equations interpreted in accordance with (8) like a phantom source with an effective equation of state as \(p=-\frac{5}{3}\rho \).
Measured current values and fitting results of \(q_0\) obtained from different models are given in [142, 143], see also [143] for both the history of the accurate estimations of q and recent q parametrizations without focusing in any concrete cosmological model, using cosmography.
We consider \(Y_{\mathrm{p}}\) and \(y_{\mathrm{DP}}\) equations given in Ref. [155], which were updated with respect to the ones given in [153, 154] in accordance with, e.g, most importantly, the change in the recommended neutron lifetime from \(\tau _\mathrm{n}=885.7\pm 0.8\,\mathrm{s}\) to \(\tau _\mathrm{n}=881.5\pm 1.5\,\mathrm{s}\) by Particle Data Group in 2011 [156]. To adopt these equations from [155] for our work, as we assume the standard particle physics is valid, we set \(\xi =0\) in the original equations given in [155], where \(\xi \) is a parameter to quantify the lepton asymmetry in non-standard particle physics.
In more complicated anisotropic spacetimes, the anisotropic spatial curvature contributes to the shear propagation equation and hence to the redshift dependency of the shear scalar. For instance, the most general spatially flat anisotropic spacetimes are of Bianchi type VII\(_0\), which, in addition to the simple expansion-rate anisotropies, yield anisotropic spatial curvature that mimics traceless anisotropic fluid [48].
We make use of a modified version of the simple and fast MCMC code that computes expansion rates and distances from the Friedmann equation, named SimpleMC [21, 145]. The method we use constrains the expansion anisotropy through its contribution to the average expansion rate of the Universe.
It can be seen from H(z) given in (34) that the sufficiently large negative \(\omega \) values also lead to small deviations from \(\Lambda \)CDM, but in this case we have slightly flatter redshift dependence of the expansion anisotropy, which in turn can relax the constraint on \(\Omega _{\sigma ^2,0}\) at few percents level compared to the one we have obtained along with large positive \(\omega \) values. Nevertheless, for such large negative values of \(\omega \), being less than the scale invariant limit \(-\frac{3}{2}\), one also should cope with stability issues.
For example, the deviation of the effective EoS parameter corresponding to the expansion anisotropy from \(w_{\sigma ^2}=1\) (corresponding to the GR limit) becomes just \(\sim 10^{-5}\) for \(\omega =40000\), whilst it is already very small, viz., \(\sim 10^{-3}\), for the cosmological lower bound \(\omega \sim 250\) (and \(\omega \sim 890\)).
The authors thank to Shahin Sheikh-Jabbari, Mehmet Özkan and Suresh Kumar for valuable discussions. Ö.A. acknowledges the support of the Turkish Academy of Sciences in scheme of the Outstanding Young Scientist Award (TÜBA-GEBİP). Ö.A. is grateful for the hospitality of the Abdus Salam International Center for Theoretical Physics (ICTP) while the part of this research was being carried out. N.K. acknowledges the post-doctoral research support from the İstanbul Technical University (ITU). J.A.V. acknowledges the support provided by the grants FOSEC SEP-CONACYT Investigación Básica A1-S-21925, and UNAM-PAPIIT IA102219.
A.A. Starobinsky, A new type of isotropic cosmological models without singularity. Phys. Lett. B 91, 99 (1980)ADSzbMATHCrossRefGoogle Scholar
A.A. Starobinsky, A new type of isotropic cosmological models without singularity. Adv. Ser. Astrophys. Cosmol. 3, 130 (1987)Google Scholar
A.H. Guth, The inflationary universe: a possible solution to the horizon and flatness problems. Phys. Rev. D 23, 347 (1981)ADSzbMATHCrossRefGoogle Scholar
A.H. Guth, The inflationary universe: a possible solution to the horizon and flatness problems. Adv. Ser. Astrophys. Cosmol. 3, 139 (1987)Google Scholar
A.D. Linde, A new inflationary universe scenario: a possible solution of the horizon, flatness, homogeneity, isotropy and primordial monopole problems. Phys. Lett. B 108, 389 (1982)ADSCrossRefGoogle Scholar
A.D. Linde, A new inflationary universe scenario: a possible solution of the horizon, flatness, homogeneity, isotropy and primordial monopole problems. Adv. Ser. Astrophys. Cosmol. 3, 149 (1987)Google Scholar
A. Albrecht, P.J. Steinhardt, Cosmology for grand unified theories with radiatively induced symmetry breaking. Phys. Rev. Lett. 48, 1220 (1982)ADSCrossRefGoogle Scholar
A. Albrecht, P.J. Steinhardt, Cosmology for grand unified theories with radiatively induced symmetry breaking. Adv. Ser. Astrophys. Cosmol. 3, 158 (1987)Google Scholar
P.A.R. Ade et al., (Planck collaboration), Planck 2015 results. XIII. Cosmological parameters. Astron. Astrophys 594, A13 (2016). arXiv:1502.01589
S. Alam et al., (BOSS Collaboration), The clustering of galaxies in the completed SDSS-III Baryon Oscillation Spectroscopic Survey: cosmological analysis of the DR12 galaxy sample. Mon. Not. Roy. Astron. Soc 470, 2617 (2017). arXiv:1607.03155
T.M.C. Abbott et al., (DES Collaboration), Dark Energy Survey year 1 results: Cosmological constraints from galaxy clustering and weak lensing. Phys. Rev. D 98, 043526 (2018). arXiv:1708.01530
N. Aghanim et al. (Planck Collaboration), Planck 2018 results. VI. Cosmological parameters, [arXiv:1807.06209]
S. Weinberg, The cosmological constant problem. Rev. Mod. Phys. 61, 1 (1989)ADSMathSciNetzbMATHCrossRefGoogle Scholar
V. Sahni, A.A. Starobinsky, The case for a positive cosmological Lambda term. Int. J. Mod. Phys. D 9, 373 (2000). arXiv:astro-ph/9904398 ADSCrossRefGoogle Scholar
P.J.E. Peebles, B. Ratra, The cosmological constant and dark energy. Rev. Mod. Phys. 75, 559 (2003). arXiv:astro-ph/0207347 ADSMathSciNetzbMATHCrossRefGoogle Scholar
T. Padmanabhan, Cosmological constant: the weight of the vacuum. Phys. Rept. 380, 235 (2003). arXiv:hep-th/0212290 ADSMathSciNetzbMATHCrossRefGoogle Scholar
H.E.S. Velten, R.F. vom Marttens, W. Zimdahl, Aspects of the cosmological coincidence problem. Eur. Phys. J. C 74, 3160 (2014). arXiv:1410.2509 ADSCrossRefGoogle Scholar
A. Padilla, Lectures on the cosmological constant problem. arXiv:1502.05296
M. López-Corredoira, Tests and problems of the standard model in Cosmology. Found. Phys 47, 6 711 (2017). arXiv:1701.08720 ADSMathSciNetzbMATHCrossRefGoogle Scholar
T. Delubac et al., (BOSS Collaboration), Baryon acoustic oscillations in the Ly\(\alpha \) forest of BOSS DR11 quasars. Astron. Astrophys. A 574, 59 (2015). arXiv:1404.1801
É. Aubourg et al., Cosmological implications of baryon acoustic oscillation measurements. Phys. Rev. D 92, 123516 (2015). arXiv:1411.1074 ADSCrossRefGoogle Scholar
A.G. Riess et al., A 2.4% determination of the local value of the hubble constant. Astrophys. J. 826, 1 56 (2016). arXiv:1604.01424 ADSCrossRefGoogle Scholar
G.B. Zhao et al., Dynamical dark energy in light of the latest observations. Nat. Astron. 1, 627 (2017). arXiv:1701.08165 ADSCrossRefGoogle Scholar
W.L. Freedman, Cosmology at a crossroads. Nat. Astron. 1, 0121 (2017). arXiv:1706.02739 ADSCrossRefGoogle Scholar
J.S. Bullock, M. Boylan-Kolchin, Small-scale challenges to the \(\Lambda \)CDM paradigm. Ann. Rev. Astron. Astrophys. 55, 343 (2017). arXiv:1707.04256 ADSCrossRefGoogle Scholar
E Di Valentino, Crack in the cosmological paradigm. Nat. Astron. 1, 569 (2017). arXiv:1709.04046 ADSCrossRefGoogle Scholar
A. Bhattacharyya, U. Alam, K.L. Pandey, S. Das, S. Pal, Are \(H_0\) and \(\sigma _8\) tensions generic to present cosmological data? Astrophys. J. 876, 143 (2019). arXiv:1805.04716 ADSCrossRefGoogle Scholar
M. Raveri, W. Hu, Concordance and discordance in cosmology. Phys. Rev. D 99, 043506 (2019). arXiv:1806.04649 ADSMathSciNetCrossRefGoogle Scholar
J.A. Vazquez, M. Bridges, M.P. Hobson, A.N. Lasenby, Reconstruction of the dark energy equation of state. J. Cosmol. Astropart. P 1209, 020 (2012). arXiv:1205.0847 ADSCrossRefGoogle Scholar
S. Hee, J.A. Vazquez, W.J. Handley, M.P. Hobson, A.N. Lasenby, Constraining the dark energy equation of state using Bayes theorem and the Kullback–Leibler divergence. Mon. Not. R. Astron. Soc 466, 369 (2017). arXiv:1607.00270 ADSCrossRefGoogle Scholar
Y. Wang, L. Pogosian, G.B. Zhao, A. Zucca, Evolution of dark energy reconstructed from the latest observations. Astrophys. J. Lett. 869, 1 (2018). arXiv:1807.03772 ADSCrossRefGoogle Scholar
V. Sahni, A. Shafieloo, A.A. Starobinsky, Model independent evidence for dark energy evolution from Baryon acoustic oscillations. Astrophys. J. 793, L40 (2014). arXiv:1406.2209 ADSCrossRefGoogle Scholar
A. Shafieloo, B. L'Huillier, A.A. Starobinsky, Falsifying \(\Lambda \)CDM: model-independent tests of the concordance model with eBOSS DR14Q and Pantheon. Phys. Rev. D 98, 083526 (2018). arXiv:1804.04320 ADSCrossRefGoogle Scholar
R.R. Caldwell, M. Kamionkowski, The physics of cosmic acceleration. Ann. Rev. Nucl. Part. Sci. 59, 397 (2009). arXiv:0903.0866 ADSCrossRefGoogle Scholar
E.J. Copeland, M. Sami, S. Tsujikawa, Dynamics of dark energy. Int. J. Mod. Phys. D 15, 1753 (2006). arXiv:hep-th/0603057 ADSMathSciNetzbMATHCrossRefGoogle Scholar
A. De Felice, S. Tsujikawa, \(f(R)\) theories. Living Rev. Relat. 13, 3 (2010). arXiv:1002.4928 ADSzbMATHCrossRefGoogle Scholar
T. Clifton, P.G. Ferreira, A. Padilla, C. Skordis, Modified gravity and cosmology. Phys. Rept. 513, 1 (2012). arXiv:1106.2476 ADSMathSciNetCrossRefGoogle Scholar
S. Capozziello, M. De Laurentis, Extended theories of gravity. Phys. Rep. 509, 167 (2011). arXiv:1108.6266 ADSMathSciNetCrossRefGoogle Scholar
K. Bamba, S. Capozziello, S. Nojiri, S.D. Odintsov, Dark energy cosmology: the equivalent description via different theoretical models and cosmography tests. Astrophys. Space Sci. 342, 155 (2012). arXiv:1205.3421 ADSzbMATHCrossRefGoogle Scholar
S. Nojiri, S.D. Odintsov, V.K. Oikonomou, Modified gravity theories on a nutshell: inflation, bounce and late-time evolution. Phys. Rep. 692, 1 (2017). arXiv:1705.11098 ADSMathSciNetzbMATHCrossRefGoogle Scholar
R.R. Caldwell, A Phantom menace? Phys. Lett. B 545, 23 (2002). arXiv:astro-ph/9908168 ADSCrossRefGoogle Scholar
Z.K. Guo, Y.S. Piao, X.M. Zhang, Y.Z. Zhang, Cosmological evolution of a quintom model of dark energy. Phys. Lett. B 608, 177 (2005). arXiv:astro-ph/0410654 ADSCrossRefGoogle Scholar
Y.F. Cai, E.N. Saridakis, M.R. Setare, J.Q. Xia, Quintom cosmology: theoretical implications and observations. Phys. Rept. 493, 1 (2010). arXiv:0909.2776 ADSMathSciNetCrossRefGoogle Scholar
S. Dubovsky, T. Gregoire, A. Nicolis, R. Rattazzi, Null energy condition and superluminal propagation. JHEP 0603, 025 (2006). arXiv:hep-th/0512260 ADSMathSciNetzbMATHCrossRefGoogle Scholar
B. Boisseau, G. Esposito-Farese, D. Polarski, A.A. Starobinsky, Reconstruction of a scalar tensor theory of gravity in an accelerating universe. Phys. Rev. Lett. 85, 2236 (2000). arXiv:gr-qc/0001066 ADSCrossRefGoogle Scholar
R. Gannouji, D. Polarski, A. Ranquet, A.A. Starobinsky, Scalar-tensor models of normal and phantom dark energy. J. Cosmol. Astropart. P. 0609, 016 (2006). arXiv:astro-ph/0606287 ADSCrossRefGoogle Scholar
D. Tamayo, J.A. Vazquez, Fourier-series expansion of the dark-energy equation of state. Mon. Not. R. Astron. Soc. 487, 729 (2019). arXiv:1901.08679 ADSCrossRefGoogle Scholar
J.D. Barrow, Cosmological limits on slightly skew stresses. Phys. Rev. D 55, 7451 (1997). arXiv:gr-qc/9701038 ADSCrossRefGoogle Scholar
A. de Oliveira-Costa, M. Tegmark, M. Zaldarriaga, A. Hamilton, The significance of the largest scale CMB fluctuations in WMAP. Phys. Rev. D 69, 063516 (2004). arXiv:astro-ph/0307282 ADSCrossRefGoogle Scholar
H.K. Eriksen, F.K. Hansen, A.J. Banday, K.M. Gorski, P.B. Lilje, Asymmetries in the Cosmic Microwave Background anisotropy field. Astrophys. J. 605, 14 (2004). Erratum: [Astrophys. J. 609, 1198 (2004)]. arXiv:astro-ph/0307507 ADSCrossRefGoogle Scholar
P. Vielva, E. Martinez-Gonzalez, R.B. Barreiro, J.L. Sanz, L. Cayon, Detection of non-Gaussianity in the WMAP 1—year data using spherical wavelets. Astrophys. J. 609, 22 (2004). arXiv:astro-ph/0310273 ADSCrossRefGoogle Scholar
M. Cruz, E. Martinez-Gonzalez, P. Vielva, L. Cayon, Detection of a non-gaussian spot in wmap. Mon. Not. R. Astron. Soc. 356, 29 (2005). arXiv:astro-ph/0405341 ADSCrossRefGoogle Scholar
P.A.R. Ade et al., (Planck Collaboration), Planck 2013 results. XXIII. Isotropy and statistics of the CMB. Astron. Astrophys. A 571, 23 (2014). arXiv:1303.5083
P.A.R. Ade et al., (Planck Collaboration), Planck 2013 Results. XXIV. Constraints on primordial non-Gaussianity. Astron. Astrophys. A 571, 24 (2014). arXiv:1303.5084
P.A.R. Ade et al., (Planck Collaboration), Planck 2013 results. XXV. Searches for cosmic strings and other topological defects. Astron. Astrophys. A 571, 25 (2014). arXiv:1303.5085
P.A.R. Ade et al., (Planck Collaboration), Planck 2013 results. XXVI. Background geometry and topology of the Universe. Astron. Astrophys. A 571, 26 (2014). arXiv:1303.5086
C.L. Bennett, Seven-year Wilkinson Microwave Anisotropy Probe (WMAP) observations: are there cosmic microwave background anomalies? Astrophys. J. Suppl. 192, 17 (2011). arXiv:1001.4758 ADSCrossRefGoogle Scholar
A. Pontzen, H.V. Peiris, The cut-sky cosmic microwave background is not anomalous. Phys. Rev. D 81, 103008 (2010). arXiv:1004.2706 ADSCrossRefGoogle Scholar
G. Efstathiou, A maximum likelihood analysis of the low CMB multipoles from WMAP. Mon. Not. R. Astron. Soc. 348, 885 (2004). arXiv:astro-ph/0310207 ADSCrossRefGoogle Scholar
P.A.R. Ade et al., Planck collaboration Planck 2013 results XVI cosmological parameters. Astron. Astrophys.A 571, 16 (2014). arXiv:1303.5076 CrossRefGoogle Scholar
S. De Deo, R.R. Caldwell, P.J. Steinhardt, Effects of the sound speed of quintessence on the microwave background and large scale structure, Phys. Rev. D 67, 103509 (2003). Erratum: [Phys. Rev. D 69, 129902 (2004)]. arXiv:astro-ph/0301284
J.A. Vazquez, A.N. Lasenby, M. Bridges, M.P. Hobson, A Bayesian study of the primordial power spectrum from a novel closed universe model. Mon. Not. R. Astron. Soc. 422, 1948 (2012). arXiv:1103.4619 ADSCrossRefGoogle Scholar
J.M. Cline, P. Crotty, J. Lesgourgues, Does the small CMB quadrupole moment suggest new physics? J. Cosmol. Astropart. P. 0309, 010 (2003). arXiv:astro-ph/0304558 ADSCrossRefGoogle Scholar
S. Tsujikawa, R. Maartens, R. Brandenberger, Non-commutative inflation and the CMB. Phys. Lett. B 574, 141 (2003). arXiv:astro-ph/0308169 ADSzbMATHCrossRefGoogle Scholar
T. Koivisto, D.F. Mota, Dark energy anisotropic stress and large scale structure formation. Phys. Rev. D 73, 083502 (2006). arXiv:astro-ph/0512135 ADSMathSciNetCrossRefGoogle Scholar
L. Campanelli, P. Cea, L. Tedesco, Ellipsoidal universe can solve the CMB quadrupole problem. Phys. Rev. Lett. 97, 131302 (2006). Erratum: [Phys. Rev. Lett. 97, 209903 (2006)]. arXiv:astro-ph/0606266
A. Gruppuso, A complete statistical analysis for the quadrupole amplitude in an ellipsoidal universe. Phys. Rev. D 76, 083010 (2007). arXiv:0705.2536 ADSMathSciNetCrossRefGoogle Scholar
L. Campanelli, P. Cea, L. Tedesco, Cosmic microwave background quadrupole and ellipsoidal universe. Phys. Rev. D 76, 063007 (2007). arXiv:0706.3802 ADSCrossRefGoogle Scholar
D.C. Rodrigues, Anisotropic cosmological constant and the CMB quadrupole anomaly. Phys. Rev. D 77, 023534 (2008). arXiv:0708.1168 ADSCrossRefGoogle Scholar
L. Campanelli, P. Cea, G.L. Fogli, L. Tedesco, Anisotropic dark energy and ellipsoidal universe. Int. J. Mod. Phys. D 20, 06 (2011). arXiv:1103.2658 zbMATHCrossRefGoogle Scholar
R.M. Wald, Asymptotic behavior of homogeneous cosmological models in the presence of a positive cosmological constant. Phys. Rev. D 28, 2118 (1983)ADSMathSciNetzbMATHCrossRefGoogle Scholar
A.A. Starobinsky, Isotropization of arbitrary cosmological expansion given an effective cosmological constant. JETP Lett. 37, 66 (1983)ADSGoogle Scholar
L.P. Chimento, M.I. Forte, Anisotropic k-essence cosmologies. Phys. Rev. D 73, 063502 (2006). arXiv:astro-ph/0510726 ADSCrossRefGoogle Scholar
R.A. Battye, A. Moss, Anisotropic perturbations due to dark energy. Phys. Rev. D 74, 041301 (2006). arXiv:astro-ph/0602377 ADSCrossRefGoogle Scholar
T. Koivisto, D.F. Mota, Accelerating cosmologies with an anisotropic equation of state. Astrophys. J. 679, 1 (2008). arXiv:0707.0279 ADSCrossRefGoogle Scholar
T. Koivisto, D.F. Mota, Anisotropic dark energy: dynamics of background and perturbations. J. Cosmol. Astropart. P. 0806, 018 (2008). arXiv:0801.3676 ADSCrossRefGoogle Scholar
T. Koivisto, D.F. Mota, Vector field models of inflation and dark energy. J. Cosmol. Astropart. P. 0808, 021 (2008). arXiv:0805.4229 ADSMathSciNetCrossRefGoogle Scholar
Ö. Akarsu, C.B. Kılınç, LRS Bianchi type I models with anisotropic dark energy and constant deceleration parameter. Gen. Rel. Gravit. 42, 119 (2010). arXiv:0807.4867 ADSMathSciNetzbMATHCrossRefGoogle Scholar
A.R. Cooray, D.E. Holz, R. Caldwell, Measuring dark energy spatial inhomogeneity with supernova data. J. Cosmol. Astropart. P. 1011, 015 (2010). arXiv:0812.0376 ADSCrossRefGoogle Scholar
L. Campanelli, A model of universe anisotropization. Phys. Rev. D 80, 3 063006 (2009). arXiv:0907.3703 MathSciNetGoogle Scholar
Ö. Akarsu, T. Dereli, N. Oflaz, Accelerating anisotropic cosmologies in Brans Dicke gravity coupled to a mass-varying vector field. Class. Quant. Gravit. 31, 045020 (2014). arXiv:1311.2573 ADSzbMATHCrossRefGoogle Scholar
T.S. Koivisto, F.R. Urban, Disformal vectors and anisotropies on a warped brane Hulluilla on Halvat Huvit. J. Cosmol. Astropart. P. 1503, 003 (2015). arXiv:1407.3445 ADSMathSciNetCrossRefGoogle Scholar
T.S. Koivisto, F.R. Urban, Doubly-boosted vector cosmologies from disformal metrics. Phys. Script. 90, 095301 (2015). arXiv:1503.01684 ADSCrossRefGoogle Scholar
D.F. Mota, J.R. Kristiansen, T. Koivisto, N.E. Groeneboom, Constraining dark energy anisotropic stress. Mon. Not. R. Astron. Soc. 382, 793 (2007). arXiv:0708.0830 ADSCrossRefGoogle Scholar
S. Appleby, R. Battye, A. Moss, Constraints on the anisotropy of dark energy. Phys. Rev. D 81, 081301 (2010). arXiv:0912.0397 ADSCrossRefGoogle Scholar
S.A. Appleby, E.V. Linder, Probing dark energy anisotropy. Phys. Rev. D 87, 023532 (2013). arXiv:1210.8221 ADSCrossRefGoogle Scholar
L. Amendola, S. Fogli, A. Guarnizo, M. Kunz, A. Vollmer, Model-independent constraints on the cosmological anisotropic stress. Phys. Rev. D 89, 063538 (2014). arXiv:1311.4765 ADSCrossRefGoogle Scholar
E. Martinez-Gonzalez, J.L. Sanz, \(\delta \)T/T and the isotropy of the universe. Astron. Astrophys. 300, 346 (1995)ADSGoogle Scholar
E.F. Bunn, P. Ferreira, J. Silk, How anisotropic is our universe? Phys. Rev. Lett. 77, 2883 (1996). arXiv:astro-ph/9605123 ADSCrossRefGoogle Scholar
A. Kogut, G. Hinshaw, A.J. Banday, Limits to global rotation and shear from the COBE DMR four year sky maps. Phys. Rev. D 55, 1901 (1997). arXiv:astro-ph/9701090 ADSCrossRefGoogle Scholar
J. Barrow, Light elements and the isotropy of the Universe. Mon. Not. R. Astron. Soc. 175, 359 (1976)ADSCrossRefGoogle Scholar
A. Pontzen, Scholarpedia 11, 32340, revision 153016 (2016)Google Scholar
D. Saadeh, S.M. Feeney, A. Pontzen, H.V. Peiris, J.D. McEwen, A framework for testing isotropy with the cosmic microwave background. Mon. Not. R. Astron. Soc. 462, 1802 (2016). arXiv:1604.01024 ADSCrossRefGoogle Scholar
D. Saadeh, S.M. Feeney, A. Pontzen, H.V. Peiris, J.D. McEwen, How isotropic is the Universe? Phys. Rev. Lett. 117, 131302 (2016). arXiv:1605.07178 ADSCrossRefGoogle Scholar
Ö. Akarsu, S. Kumar, S. Sharma, L. Tedesco, Constraints on a Bianchi type I spacetime extension of the standard \(\Lambda \)CDM model. Phys. Rev. D 100, 023532 (2019). arXiv:1905.06949 ADSMathSciNetCrossRefGoogle Scholar
L. Campanelli, P. Cea, G.L. Fogli, A. Marrone, Testing the isotropy of the universe with type Ia supernovae. Phys. Rev. D 83, 103503 (2011). arXiv:1012.5596 ADSCrossRefGoogle Scholar
Y.Y. Wang, F.Y. Wang, Testing the isotropy of the universe with type Ia supernovae in a model-independent way. Mon. Not. R. Astron. Soc. 474, 3516 (2018). arXiv:1711.05974 ADSCrossRefGoogle Scholar
W. Hu, I. Sawicki, A parameterized post-friedmann framework for modified gravity. Phys. Rev. D 76, 104043 (2007). arXiv:0708.1190 ADSCrossRefGoogle Scholar
M. Kunz, L. Amendola, D. Sapone, Dark energy phenomenology (2008). arXiv:0806.1323
B. Boisseau, Exact cosmological solution of a Scalar–tensor gravity theory compatible with the \(\Lambda \)CDM model. Phys. Rev. D 83, 043521 (2011). arXiv:1011.2915 ADSCrossRefGoogle Scholar
L. Amendola et al. (Euclid Theory Working Group), Cosmology and fundamental physics with the Euclid satellite. Living Rev. Relat. 21 (2018)Google Scholar
C.G. Callan, D. Friedan, E.J. Martinez, M.J. Perry, Strings in background fields. Nucl. Phys. B 262, 593 (1985)ADSMathSciNetCrossRefGoogle Scholar
E.S. Fradkin, A.A. Tseytlin, Quantum string theory effective action. Nucl. Phys. B 261, 1 (1985). Erratum: [Nucl. Phys. B 269, 745 (1986)]ADSMathSciNetCrossRefGoogle Scholar
M.J. Duff, R.R. Khuri, J.X. Lu, String solitons. Phys. Rept. 259, 213 (1995). arXiv:hep-th/9412184 ADSMathSciNetCrossRefGoogle Scholar
J.E. Lidsey, D. Wands, E.J. Copeland, Superstring cosmology. Phys. Rept. 337, 343 (2000). arXiv:hep-th/9909061 ADSMathSciNetCrossRefGoogle Scholar
P. Jordan, The present state of Dirac's cosmological hypothesis. Z. Phys. 157, 112 (1959)ADSCrossRefGoogle Scholar
P.G. Bergmann, Comments on the scalar–tensor theory. Int. J. Theor. Phys. 1, 25 (1968)CrossRefGoogle Scholar
K. Nordvedt, Post-Newtonian metric for a general class of scalar–tensor gravitational theories and observational consequences. Astrophys. J. 161, 1059 (1970)ADSMathSciNetCrossRefGoogle Scholar
R.V. Wagoner, Scalar–tensor theory and gravitational waves. Phys. Rev. D 1, 3209 (1970)ADSCrossRefGoogle Scholar
C. Brans, R.H. Dicke, Mach's principle and a relativistic theory of gravitation. Phys. Rev. 124, 925 (1961)ADSMathSciNetzbMATHCrossRefGoogle Scholar
L.O. Pimentel, Energy-momentum tensor in the general scalar–tensor theory. Class. Quant Gravit. 6, 12 (1989)MathSciNetCrossRefGoogle Scholar
M.S. Madsen, Scalar fields in curved spacetimes. Class. Quant. Gravit. 5, 6 (1988)MathSciNetzbMATHCrossRefGoogle Scholar
V. Faraoni, J. Coté, Imperfect fluid description of modified gravities. Phys. Rev. D 98, 084019 (2018). arXiv:1808.02427 ADSMathSciNetCrossRefGoogle Scholar
B. Bertotti, L. Iess, P. Tortora, A test of general relativity using radio links with the Cassini spacecraft. Nature 425, 374 (2003)ADSCrossRefGoogle Scholar
L. Perivolaropoulos, PPN parameter \(\gamma \) and solar system constraints of massive Brans–Dicke theories. Phys. Rev. D 81, 047501 (2010). arXiv:0911.3401 ADSCrossRefGoogle Scholar
A. Avilez, C. Skordis, Cosmological constraints on Brans–Dicke theory. Phys. Rev. Lett. 113, 011101 (2014). arXiv:1303.4330 ADSCrossRefGoogle Scholar
M. Ballardini, F. Finelli, C. Umiltà, D. Paoletti, Cosmological constraints on induced gravity dark energy models. J. Cosmol. Astropart. P 1605, 067 (2016). arXiv:1601.03387 ADSCrossRefGoogle Scholar
C. Umiltà, M. Ballardini, F. Finelli, D. Paoletti, CMB and BAO constraints for an induced gravity dark energy model with a quartic potential. J. Cosmol. Astropart. P 1508, 017 (2015). arXiv:1507.00718 ADSCrossRefGoogle Scholar
M. Rossi, M. Ballardini, M. Braglia, F. Finelli, D. Paoletti, A.A. Starobinsky, C. Umiltà, Cosmological constraints on post-Newtonian parameters in effectively massless scalar–tensor theories of gravity. Phys. Rev. D 100, 103524 (2019). arXiv:1906.10218 ADSCrossRefGoogle Scholar
T. Clifton, J.D. Barrow, R.J. Scherrer, Constraints on the variation of G from primordial nucleosynthesis. Phys. Rev. D 71, 123526 (2005). [arXiv:astro-ph/0504418]ADSCrossRefGoogle Scholar
V. Acquaviva, L. Verde, Observational signatures of Jordan–Brans–Dicke theories of gravity. J. Cosmol. Astropart. P 0712, 001 (2007). arXiv:0709.0082 ADSCrossRefGoogle Scholar
D. Alonso, E. Bellini, P.G. Ferreira, M. Zumalacairregui, Observational future of cosmological scalar–tensor theories. Phys. Rev. D 95, 063502 (2017). arXiv:1610.09290 ADSMathSciNetCrossRefGoogle Scholar
M. Ballardini, D. Sapone, C. Umiltà, F. Finelli, D. Paoletti, Testing extended Jordan–Brans–Dicke theories with future cosmological observations. J. Cosmol. Astropart. P 05, 049 (2019). arXiv:1902.01407 ADSMathSciNetCrossRefGoogle Scholar
G.W. Horndeski, Second-order scalar–tensor field equations in a four-dimensional space. Int. J. Theor. Phys. 10, 363 (1974)MathSciNetCrossRefGoogle Scholar
C. Deffayet, X. Gao, D.A. Steer, G. Zahariade, From k-essence to generalised Galileons. Phys. Rev. D 84, 064039 (2011). arXiv:1103.3260 ADSCrossRefGoogle Scholar
A. Nicolis, R. Rattazzi, E. Trincherini, The Galileon as a local modification of gravity. Phys. Rev. D 79, 064036 (2009). arXiv:0811.2197 ADSMathSciNetCrossRefGoogle Scholar
Ö. Akarsu, A. Chopovsky, V. Shulga, E. Yalçınkaya, A. Zhuk, Weak field limit of higher dimensional massive Brans–Dicke gravity: observational constraints. Phys. Rev. D 101, 024004 (2020). arXiv:1907.04234
M. Ishak, Testing general relativity in cosmology. Living Rev. Relat. 22, 1 (2019). arXiv:1806.10122 ADSCrossRefGoogle Scholar
V. Acquaviva, C. Baccigalupi, S.M. Leach, A.R. Liddle, F. Perrotta, Structure formation constraints on the Jordan–Brans–Dicke theory. Phys. Rev. D 71, 104025 (2005). arXiv:astro-ph/0412052 ADSCrossRefGoogle Scholar
V. Faraoni, Scalar field mass in generalized gravity. Class. Quant. Gravit. 26, 145014 (2009). arXiv:0906.1901 ADSMathSciNetzbMATHCrossRefGoogle Scholar
S. Weinberg, Gravitation and cosmology: principles and applications of the general theory of relativity (Wiley, New York, 1972)Google Scholar
J.P. Mimoso, D. Wands, Anisotropic scalar–tensor cosmologies. Phys. Rev. D 52, 5612 (1995). arXiv:gr-qc/9501039 ADSMathSciNetCrossRefGoogle Scholar
P.H. Chavanis, Cosmology with a stiff matter era. Phys. Rev. D 92, 103004 (2015). arXiv:1412.0743 ADSCrossRefGoogle Scholar
K. Uehara, C.W. Kim, Brans-dicke cosmology with the cosmological constant. Phys. Rev. D 26, 2575 (1982)ADSMathSciNetCrossRefGoogle Scholar
H. Nariai, Gravitational Instability in the Brans–Dicke cosmology. Prog. Theor. Phys. 42, 544 (1969)ADSCrossRefGoogle Scholar
L.E. Gurevich, A.M. Finkelstein, V.A. Ruban, On the problem of the initial state in the isotropic scalar–tensor cosmology of Brans–Dicke. Astrophys. Space Sci. 22, 231 (1973)ADSCrossRefGoogle Scholar
J.D. Barrow, Scalar–tensor cosmologies. Phys. Rev. D 47, 5329 (1993)ADSMathSciNetCrossRefGoogle Scholar
A.R. Liddle, A. Mazumdar, J.D. Barrow, Radiation matter transition in Jordan–Brans–Dicke theory. Phys. Rev. D 58, 027302 (1998). arXiv:astro-ph/9802133 ADSCrossRefGoogle Scholar
J.L. Cervantes-Cota, Bianchi V inflation in the Brans–Dicke theory? Class. Quant. Gravit. 16, 3903 (1999). arXiv:gr-qc/9912047 ADSMathSciNetzbMATHCrossRefGoogle Scholar
J.L. Cervantes-Cota, M. Nahmad, Isotropization of Bianchi type models and a new FRW solution in Brans–Dicke theory. Gen. Relat. Gravit. 33, 767 (2001). arXiv:gr-qc/0005032 ADSzbMATHCrossRefGoogle Scholar
A.Y. Kamenshchik, E.O. Pozdeeva, A.A. Starobinsky, A. Tronconi, G. Venturi, S.Y. Vernov, Induced gravity and minimally and conformally coupled scalar fields in Bianchi-I cosmological models. Phys. Rev. 97, 023536 (2018). arXiv:1710.02681 MathSciNetGoogle Scholar
Ö. Akarsu, T. Dereli, S. Kumar, L. Xu, Probing kinematics and fate of the Universe with linearly time-varying deceleration parameter. Eur. Phys. J. Plus 129, 22 (2014). arXiv:1305.5190 CrossRefGoogle Scholar
A. Gómez-Valent, L. Amendola, \(H_0\) from cosmic chronometers and type Ia supernovae, with Gaussian processes and the novel weighted polynomial regression method. J. Cosmol. Astropart. P 1804, 051 (2018). arXiv:1802.01505 ADSCrossRefGoogle Scholar
YAB Zeldovich, Equation of state at ultra-high densities and its relativistic limitations. Zhur. Eksptl'. i Teoret. Fiz 41 (1961)Google Scholar
https://github.com/slosar/april, version May 2019
L.E. Padilla, L.O. Tellez, L.A. Escamilla, J.A. Vazquez, Cosmological parameter inference with Bayesian statistics. arXiv:1903.11127
E. Komatsu et al., (WMAP Collaboration), Seven-year Wilkinson microwave anisotropy probe (WMAP) observations: cosmological interpretation. Astrophys. J. Suppl. Ser. 192, 18 (2011). arXiv:1001.4538
S. Dodelson, Modern cosmology (Acad Press, New York, 2003)Google Scholar
D.J. Fixsen, The temperature of the cosmic microwave background. Astrophys. J. 707, 916 (2009). arXiv:0911.1955 ADSCrossRefGoogle Scholar
J.P. Uzan, Varying constants, gravitation and cosmology. Living Rev. Relat. 14, 2 (2011). arXiv:1009.5514 ADSzbMATHCrossRefGoogle Scholar
X. Chen, M. Kamionkowski, Cosmic microwave background temperature and polarization anisotropy in Brans–Dicke cosmology. Phys. Rev. D 60, 104036 (1999). arXiv:astro-ph/9905368 ADSCrossRefGoogle Scholar
L. Campanelli, Helium-4 synthesis in an anisotropic universe. Phys. Rev. D 84, 123521 (2011). arXiv:1112.2076 ADSCrossRefGoogle Scholar
J.P. Kneller, G. Steigman, BBN for pedestrians. New J. Phys. 6, 117 (2004). arXiv:astro-ph/0406320 ADSCrossRefGoogle Scholar
G. Steigman, Primordial nucleosynthesis in the precision cosmology Era. Ann. Rev. Nucl. Part. Sci. 57, 463 (2007). arXiv:0712.1100 ADSCrossRefGoogle Scholar
G. Steigman, Neutrinos and big bang nucleosynthesis. Adv. High Energy Phys. 2012, 268321 (2012). arXiv:1208.0032 zbMATHCrossRefGoogle Scholar
K. Nakamura et al., (Particle Data Group), Review of particle physics. J. Phys. G Nucl. Part. Phys 37, 075021 (2010)ADSCrossRefGoogle Scholar
E. Aver, K.A. Olive, E.D. Skillman, The effects of He I \(\lambda 10830\) on helium abundance determinations. J. Cosmol. Astropart. P 1507, 07 (2015). arXiv:1503.08146 Google Scholar
R.J. Cooke, M. Pettini, C.C. Steidel, One percent determination of the primordial deuterium abundance. Astrophys. J. 855, 102 (2018). arXiv:1710.11129 ADSCrossRefGoogle Scholar
J.D. Barrow, C. O'Toole, Spatial variations of fundamental constants. Mon. Not. R. Astron. Soc. 322, 585 (2001). arXiv:astro-ph/9904116 ADSCrossRefGoogle Scholar
T. Clifton, D.F. Mota, J.D. Barrow, Inhomogeneous gravity. Mon. Not. R. Astron. Soc. 358, 601 (2005). arXiv:gr-qc/0406001 ADSCrossRefGoogle Scholar
S.M.M. Rasouli, M. Farhoudi, H.R. Sepangi, Anisotropic cosmological model in modified Brans–Dicke theory. Class. Quant. Gravit. 28, 155004 (2011). arXiv:1105.5086 ADSMathSciNetzbMATHCrossRefGoogle Scholar
S.M.M. Rasouli, M. Farhoudi, P. Vargas Moniz, Modified Brans–Dicke theory in arbitrary dimensions. Class. Quant. Gravit. 31, 115002 (2014). arXiv:1405.0229 ADSMathSciNetzbMATHCrossRefGoogle Scholar
S.M.M. Rasouli, P. Vargas Moniz, Extended anisotropic models in noncompact Kaluza–Klein theory. Class. Quant. Gravit. 36, 075010 (2019). arXiv:1806.03684 ADSMathSciNetCrossRefGoogle Scholar
R. Maartens, V. Sahni, T.D. Saini, Anisotropy dissipation in brane world inflation. Phys. Rev. D 63, 063509 (2001). arXiv:gr-qc/0011105 ADSCrossRefGoogle Scholar
A. Campos, R. Maartens, D. Matravers, C.F. Sopuerta, Brane world cosmological models with anisotropy. Phys. Rev. D 68, 103520 (2003). arXiv:hep-th/0308158 ADSMathSciNetCrossRefGoogle Scholar
Open AccessThis article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
Funded by SCOAP3
1.Department of PhysicsIstanbul Technical UniversityIstanbulTurkey
2.Instituto de Ciencias FísicasUniversidad Nacional Autónoma de México, ApdoMorelosMexico
Akarsu, Ö., Katırcı, N., Özdemir, N. et al. Eur. Phys. J. C (2020) 80: 32. https://doi.org/10.1140/epjc/s10052-019-7580-z
DOI https://doi.org/10.1140/epjc/s10052-019-7580-z
Publisher Name Springer Berlin Heidelberg
Published in cooperation with
EPJC is an open-access journal funded by SCOAP3 and licensed under CC BY 4.0