Cosmological backreaction within the Szekeres model and emergence of spatial curvature

TL;DR

This paper investigates how non-linear cosmic structure growth within the Szekeres model influences spatial curvature, suggesting that the Universe's low-redshift spatial curvature may deviate from zero due to backreaction effects.

Contribution

It demonstrates that the evolution of spatial curvature depends on the global model, showing that backreaction can lead to a non-zero curvature at low redshifts in the Styrofoam-type Szekeres model.

Findings

In Swiss-Cheese models, the average curvature remains zero, following the background.

In Styrofoam-type models, the average curvature grows to about 0.1 at present.

The results imply low-redshift Universe may have non-zero spatial curvature due to backreaction.

Abstract

This paper discusses the phenomenon of backreaction within the Szekeres model. Cosmological backreaction describes how the mean global evolution of the Universe deviates from the Friedmannian evolution. The analysis is based on models of a single cosmological environment and the global ensemble of the Szekeres models (of the Swiss-Cheese-type and Styrofoam-type). The obtained results show that non-linear growth of cosmic structures is associated with the growth of the spatial curvature $\Omega_{\cal R}$ (in the FLRW limit $\Omega_{\cal R} \to \Omega_k$). If averaged over global scales the result depends on the assumed global model of the Universe. Within the Swiss-Cheese model, which does have a fixed background, the volume average follows the evolution of the background, and the global spatial curvature averages out to zero (the background model is the $\Lambda$CDM model, which is spatially flat). In the Styrofoam-type model, which does not have a fixed background, the mean evolution deviates from the spatially flat $\Lambda$CDM model, and the mean spatial curvature evolves to from $\Omega_{\cal R} =0 $ at the CMB to $\Omega_{\cal R} \sim 0.1$ at $z =0$. If the Styrofoam-type model correctly captures evolutionary features of the real Universe then one should expect that in our Universe, the spatial curvature should build up (local growth of cosmic structures) and its mean global average should deviate from zero (backreaction). As a result, this paper predicts that the low-redshift Universe should not be spatially flat (i.e. $\Omega_k \ne 0$, even if in the early Universe $\Omega_k = 0$) and therefore when analysing low-$z$ cosmological data one should keep $\Omega_k$ as a free parameter and independent from the CMB constraints.

Equations

$$T_{ab}=\rho u_a{u}_b+p{h}_{ab}+{\pi }_{ab},$$

$$h_{ab}=g_{ab}-u_a{u}_b,$$

$$u_{a;b}={\omega }_{ab}+{\sigma }_{ab}+\frac{1}{3}{h}_{ab}\Theta -A_a{u}_b,$$

$$\dot{\Theta }=-\frac{1}{3}{\Theta }^2-\frac{1}{2}(\rho +3p)-2({\sigma }^2-{\omega }^2)+{\mathrm{D}}^a{A}_a+A_a{A}^a+\Lambda ,$$

$${\dot{\sigma }}_{⟨ab⟩}=-\frac{2}{3}\Theta {\sigma }_{ab}-{\sigma }_{c⟨a}{\sigma }^c{}_{b⟩}-{\omega }_{⟨a}{\omega }_{b⟩}+{\mathrm{D}}_{⟨a}{A}_{b⟩}+A_{⟨a}{A}_{b⟩}-E_{ab}+\frac{1}{2}{\pi }_{ab},$$

$${\dot{\omega }}_{⟨a⟩}=-\frac{2}{3}\Theta {\omega }_a-\frac{1}{2}\mathrm{curl}{A}_a+{\sigma }_{ab}{\omega }^b.$$

$$\dot{\rho }=-\Theta (\rho +p)-{\sigma }^{ab}{\pi }_{ab},$$

$$(\rho +p)A_a=-{\mathrm{D}}_{a}p-{\mathrm{D}}^b{\pi }_{ab}-{\pi }_{ab}{A}^b.$$

$${\dot{E}}_{⟨ab⟩}=-\Theta E_{ab}-\frac{1}{2}(\rho +p){\sigma }_{ab}+\mathrm{curl}{H}_{ab}-\frac{1}{2}{\dot{\pi }}_{ab}-\frac{1}{6}\Theta {\pi }_{ab}$$

$$+3{\sigma }_{⟨a}{}^c\left(E_{b⟩c}-\frac{1}{6}{\pi }_{b⟩c}\right)+{\epsilon }_{cd⟨a}\left[2A^c{H}_{b⟩}{}^d-{\omega }^c\left(E_{b⟩}{}^d+\frac{1}{2}{\pi }_{b⟩}{}^d\right)\right],$$

$${\dot{H}}_{⟨ab⟩}=-\Theta H_{ab}-\mathrm{curl}{E}_{ab}+\frac{1}{2}\mathrm{curl}{\pi }_{ab}+3{\sigma }_{⟨a}{}^c{H}_{b⟩c}-{\epsilon }_{cd⟨a}\left(2A^c{E}_{b⟩}{}^d+{\omega }^c{H}_{b⟩}{}^d\right).$$

$$\displaystyle{\rm FLRW~{}universe~{}~{}}\left\{\begin{array}[]{llllll}E_{ab}\equiv 0\\ H_{ab}\equiv 0\\ \omega_{ab}\equiv 0\\ \sigma_{ab}\equiv 0\\ \pi_{ab}\equiv 0\\ D_{a}p\equiv 0\end{array}\right.$$

$$\dot{\Theta }=-\frac{1}{3}{\Theta }^2-\frac{1}{2}(\rho +3p)+\Lambda ,$$

$$\dot{\rho }=-\Theta (\rho +p).$$

$$3\frac{\ddot{a}}{a}=-4\pi G(\rho +3p)+\Lambda ,$$

$$3\frac{{\dot{a}}^2}{a^2}=8\pi G\rho -3\frac{k}{a^2}+\Lambda ,$$

$$3\frac{{\ddot{a}}_{\mathcal{D}}}{a_{\mathcal{D}}}=-4\pi G⟨\rho {⟩}_{\mathcal{D}}+\Lambda +{\mathcal{Q}}_{\mathcal{D}},$$

$$3\frac{{\dot{a}}_{\mathcal{D}}^2}{a_{\mathcal{D}}^2}=8\pi G⟨\rho {⟩}_{\mathcal{D}}-\frac{1}{2}⟨\mathcal{R}{⟩}_{\mathcal{D}}+\Lambda -\frac{1}{2}{\mathcal{Q}}_{\mathcal{D}},$$

$$\displaystyle\big{(}{\cal Q}_{\cal D}a^{6}_{\cal D}\dot{\big{)}}+a^{4}_{\cal D}\big{(}\langle{\mathcal{R}}\rangle_{\cal D}a^{2}_{\cal D}\dot{\big{)}}=0,$$

$$⟨A{⟩}_{\mathcal{D}}=\frac{{\int }_{\mathcal{D}}{d}^{3}x\sqrt{|h|}A}{{\int }_{\mathcal{D}}{d}^{3}x\sqrt{|h|}},$$

$$a_{\mathcal{D}}={\left(\frac{V_{\mathcal{D}}}{V_{\mathcal{D},i}}\right)}^{1/3},$$

$${\mathcal{Q}}_{\mathcal{D}}=\frac{2}{3}\left(⟨{\Theta }^2{⟩}_{\mathcal{D}}-⟨\Theta {⟩}_{\mathcal{D}}^2\right)-2⟨{\sigma }^2{⟩}_{\mathcal{D}}.$$

$$H_{\mathcal{D}}^2=\frac{8\pi G}{3}⟨\rho {⟩}_{\mathcal{D}}-\frac{1}{6}⟨\mathcal{R}{⟩}_{\mathcal{D}}+\frac{1}{3}\Lambda -\frac{1}{6}{\mathcal{Q}}_{\mathcal{D}},$$

$${\Omega }_m^{\mathcal{D}}$$

$${\Omega }_{\Lambda }^{\mathcal{D}}$$

$${\Omega }_{\mathcal{R}}^{\mathcal{D}}$$

$${\Omega }_{\mathcal{Q}}^{\mathcal{D}}$$

$$\mathrm{d}{s}^2=\mathrm{d}{t}^2-\frac{1}{\epsilon -K}{\left[R^{\prime }+\frac{R}{S}\left(S^{\prime }\cos \theta +N\sin \theta \right)\right]}^2\mathrm{d}{r}^2-R^2\left(\mathrm{d}{\theta }^2+{\sin }^2\theta \mathrm{d}{\varphi }^2\right)$$

$$+{\left(\frac{R}{S}\right)}^2\left\{{\left[S^{\prime }\sin \theta +N(1-\cos \theta )\right]}^2+{\left[({\partial }_{\varphi }N)(1-\cos \theta )\right]}^2\right\}\mathrm{d}{r}^2$$

$$+2{\left(\frac{R}{S}\right)}^2\left[S{S}^{\prime }\sin \theta +SN(1-\cos \theta )\right]\mathrm{d}r\mathrm{d}\theta -2{\left(\frac{R}{S}\right)}^{2}S({\partial }_{\varphi }N)\sin \theta (1-\cos \theta )\mathrm{d}r\mathrm{d}\varphi ,$$

$${\sigma }_{ab}=\Sigma {\mathrm{e}}_{ab},E_{ab}=\mathcal{W}{\mathrm{e}}_{ab},$$

$$\rho =\frac{2\left(M^{\prime }-3M{E}^{\prime }/E\right)}{R^2\left(R^{\prime }-R{E}^{\prime }/E\right)},$$

$$\Theta =\frac{{\dot{R}}^{\prime }+2\dot{R}{R}^{\prime }/R-3\dot{R}{E}^{\prime }/E}{R^{\prime }-R{E}^{\prime }/E},$$

$$\Sigma =-\frac{1}{3}\frac{{\dot{R}}^{\prime }-\dot{R}{R}^{\prime }/R}{R^{\prime }-R{E}^{\prime }/E},$$

$$\mathcal{W}=\frac{M}{3R^3}\frac{3R^{\prime }-R{M}^{\prime }/M}{R^{\prime }-R{E}^{\prime }/E},$$

$$\frac{1}{6}\mathcal{R}=\frac{1}{3}\rho +6{\Sigma }^2-\frac{1}{9}{\Theta }^2+\frac{1}{3}\Lambda ,$$

$$\mathcal{R}=2\frac{K}{R^2}\left(1+\frac{R{K}^{\prime }/K-2R{E}^{\prime }/E}{R^{\prime }-R{E}^{\prime }/E}\right).$$

$${\dot{R}}^2=-K+\frac{2M(r)}{R}+\frac{1}{3}\Lambda R^2.$$

$$r=R_i,$$

Full text

1 Introduction

Backreaction is a process that describes feedback of structure formation on the mean global evolution of the Universe. As long as the evolution of cosmic structures is well within the linear regime, the Universe should successfully be approximated by perturbations around the FLRW model and its global evolution should follow the Friedmann solution. This part of our Universe’s evolution seems to be well understood. What is less understood and still debatable among cosmologists is the current epoch of the Universe’s evolution with its non-linear growth of cosmic structures. Some cosmologists argue that in the non-linear regime the average evolution of the Universe may deviate from the Friedmannian evolution.

The Friedmannian evolution is at the core of the standard cosmological model. Observational data are mostly interpreted within this framework. Also this framework is embedded within $N$-body simulations. The $N$-body simulations are employed to trace the evolution of our Universe in the non-linear regime. Mostly this is done by employing the Newtonian physics on small scales, and assuming that the global evolution follows the Friedmann solution and is unaffected by local interactions. In principle it is possible to construct a Newtonian cosmology based on a system that expands homothetically and obey the Friedmannian evolution [1]. Also the linear perturbations of such a system are known to be consistent with linear perturbations of the FLRW model [2]. However, in the non-linear regime the situation may be different [3, 4], for example one of the recent studies shows that the relativistic effects lead to Yukawa-type interactions between matter particles [5]. Studies based on the weak-filed limit (i.e. applying post-Newtonian corrections) suggest that the mean evolution is well approximated by the Friedmannian evolution [6, 7, 8]. On the other hand, studies that try to implement backreaction-type effects into the Newtonian Cosmology and $N$-body codes find large effects [9, 10] — such an approach has sparked recently a debate, see Ref. [11] and the rebuttal response in Ref. [12].

The presence of backreaction is a mathematical consequence of the non-linear structure of equations that govern the evolution of the Universe [13]. Cosmological backreaction has been debated for the last 30 years [14, 15]. The nature of this debate focuses on the magnitude of backreaction. On one hand, predictions derived from models based on the FLRW framework and $N$-body simulations are consistent with observational data; this is inductive evidence suggesting that the backreaction could be negligibly small. On the other hand, presence of tensions between various data and experiments could be related to backreaction effects. Finally, some people argue that the presence of dark energy, which dominates the energy budget of the standard cosmological model, is just an artefact and nothing else as a manifestation of strong backreaction effects.

This paper aims to be another voice in the debate on the backreaction and presents the study of the backreaction within the Szekeres model, which is an exact inhomogeneous cosmological solution of the Einstein equations. The structure of the paper is as follows: Sec. 2 continues with the introduction to the phenomenon of backreaction; Sec. 3 introduces the Szekeres model; Secs. 4 and 5 present the results of the analysis of the backreaction phenomenon within the Szekeres model; Sec. 6 concludes the analysis.

2 Evolution of matter in the Universe and backreaction

The energy momentum tensor of a viscous fluid with no energy transfer can be written as

[TABLE]

where $\rho$ is energy density, $p$ is pressure, $\pi_{ab}$ is the anisotropic stress tensor, and $h_{ab}$ is the spatial part of the metric in 3+1 split

[TABLE]

and $u_{a}$ is the matter velocity flow, whose gradient can be decomposed as

[TABLE]

where $\omega_{ab}=u_{[a;|\sigma|}h^{\sigma}{}_{b]}$ is rotation, $\sigma_{ab}=u_{(a;|\sigma|}h^{\sigma}{}_{b)}-\frac{1}{3}h_{ab}\Theta$ is shear, $\Theta=u^{a}{}_{;a}$ is expansion, and $A^{a}=u^{a}{}_{;b}u^{b}$ is acceleration.

Evolution of expansion, shear, and rotation are given by [16, 17]

[TABLE]

The equations for density, pressure, and anisotropic stress are [16, 17]

[TABLE]

Finally, the evolution of the electric ($E_{ab}$) and magnetic ($H_{ab}$) parts of the Weyl curvature are given by [16, 17]

[TABLE]

The homogeneous and isotropic FLRW models form a special subset of all possible solutions of the above equations. The FLRW solution is characterised by vanishing Weyl curvature, vanishing rotation and shear, zero anisotropic stress and pressure gradients [18]

[TABLE]

In such a case all the above given evolution equations reduce only to 2 equations that fully describe the evolution of a spatially homogeneous and isotropic system

[TABLE]

After some algebra, it can be shown that the above equations are equivalent to the Friedmann equations

[TABLE]

where the relation between the scale factor $a(t)$ and the expansion rate is $\Theta=3\dot{a}/a$, and the spatial curvature is ${\cal R}=6k/a^{2}$.

If the universe is homogeneous, then the average evolution is exactly the same as the evolution of an individual worldline, and is given by the Friedmann equations. If the universe is inhomogeneous, then the average over all individual worldlines, may deviate from the solution of a uniform universe, and therefore deviate from the Friedmannian evolution. This is what backreaction describes. When studying backreaction, one focuses on all neglected (in the FLRW case) terms [cf. eq. (2.3)–(2.9)] and investigates if it is possible that all these terms can affect the mean global evolution of the inhomogeneous system. In other words, if global (mean) evolution of the volume of the Universe (i.e. $\Theta$) and/or matter (i.e. $\rho$) is the same as prescribed by the Friedmann solution [i.e. eqs. (2.11) and (2.12) or equivalently by eqs. (2.13) and (2.14)] or if the contribution from the shear, rotation, and Weyl curvature can affect the expansion rate $\Theta$, and subsequently change its global (mean) evolution compared to the FLRW case.

Because the evolution equations (2.3)–(2.9) are complicated, we still lack a satisfactory description of backreaction for a real Universe. If we limit the analysis only to irrotational and pressureless fluids (the Szekeres model discussed below can only describe irrotational dust) then the scalar parts of the above equations can be averaged and reduced to [19]

[TABLE]

where the dot $\dot{}\,$ denotes partial time derivative $\dot{}\equiv\partial_{t}$, $\langle{\mathcal{R}}\rangle_{\cal D}$ is an average of the spatial Ricci scalar $\mathcal{R}$, $\langle{\ }\rangle_{\cal D}$ is the volume average over the hypersurface of constant time

[TABLE]

the scale factor $a_{\cal D}$ is defined as

[TABLE]

where VD is the volume of the domain ${\cal D}$, and $V_{{\cal D},i}$ is its initial value. Finally the function ${\cal Q}_{\cal D}$ is

[TABLE]

Equations (2.16) can be used to define the Hubble parameter $H_{\cal D}$

[TABLE]

which then can be used to introduce the cosmic quartet [20]

[TABLE]

It follows from eq. (2.17) that iff ${\cal Q}_{\cal D}=0$ at all times then the average spatial curvature reduces to the Friedmannian curvature, i.e. $\langle{\mathcal{R}}\rangle_{\cal D}\sim a_{\cal D}^{-2}$. This implies that iff ${\cal Q}_{\cal D}$ vanishes at all times then the Buchert equations reduce to the Friedmann equations, meaning that the non-linear effects associated with structure formation cannot affect the average evolution of the Universe. The question and issue of the ongoing debate on backreaction is related to the amplitude of backreaction. Some studies suggest large backreaction [21, 22, 23, 24, 25, 26, 27, 28, 29, 30], while some suggest otherwise [31, 32, 33]. The debate on whether the backreaction in realistic models of the Universe is strictly zero has already been settled [14]. However, it is still debated whether a small deviation from zero also means almost Friedmannian evolution, or not. This question is addressed in this paper using the Szekeres solution of the Einstein equations. The investigation is based on a realistic model of the local cosmological environment [34] and the ensemble of Szekeres models that is based on this type of structures.

3 Szekeres model

The Szekeres model [35] is one of the most general inhomogeneous, exact, cosmological solutions of the Einstein equations [18, 36]. The metric of the Szekeres model [35, 37] in the spherical coordinates is [38]

[TABLE]

where ${}^{\prime}\equiv\partial/\partial r$, $N(r,\phi)\equiv P^{\prime}\cos\phi+Q^{\prime}\sin\phi$, $\varepsilon=\pm 1,0$ and $K=K(r)\leq\varepsilon$, $S=S(r)$, $P=P(r)$, and $Q=Q(r)$ are arbitrary functions of $r$.

Systems that can be described using the Szekeres solution have no vorticity $\omega_{ab}=0$, no viscosity $\pi_{ab}=0$, no pressure $p=0$, and no gravitational radiation $H_{ab}=0$. The Szekeres models are of Petrov type D, the shear and electric Weyl tensor can be written as

[TABLE]

where ${\rm e}_{ab}=h_{ab}-3z_{a}z_{b}$ where $z^{a}$ is a space-like unit vector aligned with the Weyl principal tetrad. As a result the fluid equations (2.3)–(2.9) reduce only to 4 scalars [39, 36, 40]

[TABLE]

where $E^{\prime}/E=-(S^{\prime}\cos\theta+N\sin\theta)/S$. The spatial curvature follows from the “Hamiltonian” constraint

[TABLE]

and is given by

[TABLE]

Thus, the whole evolution of the system, is reduced only to a single equation for the function $R$

[TABLE]

To define the Szekeres model 5 arbitrary functions of radial coordinate $r$ need to be specified. In this paper, the radial coordinate has been chosen as

[TABLE]

where $R_{i}$ is the value of $R$ at the last scattering instant, and the arbitrary functions were chosen to be

[TABLE]

Thus, the model is prescribed by following parameters: $r_{0}$ (size of perturbation), $\Delta r$ (transition zone), $m_{0}$ (central density contrast of the initial perturbation), $\eta$ (dipole’s size), and $\rho_{i}$ (background density at the initial instant). The initial instant is set to be the last scattering instant, so $\rho_{i}=(1+z_{CMB})^{3}\Omega_{m}3H_{0}^{2}/(8\pi G)$, where $z_{CMB}$ is the CMB’s redshift $z_{CMB}=1090$.

4 Local non-linear evolution and emergence of local spatial curvature within a single cosmological environment

The model of a local cosmological environment is specified by fixing the free parameters to

[TABLE]

The evolution of the system is calculated by solving eq. (3.8) from the last scattering instant till the present-day instant. The present-day density contrast presented in the right panel of Fig. 1 is defined as

[TABLE]

where $\rho_{0}$ is the present-day background density, $\rho_{0}=3H_{0}^{2}\Omega_{m}/(8\pi G)$. The evolution of the density contrast at 4 different locations is also presented in Fig. 1. The evolution of the density contrast evaluated within this Szekeres model is compared to the evolution of the density contrast $\delta_{lin}$ which follows from the linear approximation [41]

[TABLE]

where the initial conditions for $\delta_{i}$ and $\dot{\delta}_{i}$ have been chosen to be the same as in the Szekeres model at the last scattering instant.

Two upper left panels in Fig. 1 show places where the models is quite homogeneous, and where the evolution is well described by the linear approximation. The other two panels show places where the inhomogeneity is large and the growth is non-linear. A more detailed map on local inhomogeneity and averages is presented in Fig. 2. Figure 2 presents the volume averages evaluated within spherical domains. Each point in Fig. 2 presents the value of the averaging within a domain centred at that point and of radius 5 Mpc (cf. the matter horizon in Ref. [42]). The volume averages were evaluated using the code SzReD111The code SzReD is not yet publicly available under a free licence. It is accessible on collaboration basis and anyone interested in using the code is advised to contact the author.. The code uses the Healpix grid to evaluate integrals around any given point in space and time. The average expansion rate and shear at the present instant normalised by $H_{0}$ are presented in the upper panels. The kinematic backreaction $\Omega_{\cal Q}^{\cal D}$ and spatial curvature $\Omega_{\mathcal{R}}^{\cal D}$ are presented in the lower panels of Fig. 2.

By comparing Figs. 1 and 2 it is apparent that the non-linear evolution is associated with places where the present-day spatial curvature $\Omega_{\cal R}$ is not negligible. It should be pointed out, that at the initial instant the spatial curvature is negligibly small, which means that the non-linear growth is related to the emergence of the spatial curvature. The emergence of the spatial curvature is better depicted in Fig. 3. Figure 3 shows the evolution of backreaction and emergence of the spatial curvature at 4 regions presented in Fig. 1. The evolution of the cosmic quartet (2.20) at these 4 different locations (presented in Fig. 3) shows that regions which evolve into the non-linear regime have also a substantial build-up of the spatial curvature, and the contribution of the spatial curvature to the evolution of the system is comparable with the contribution from matter.

There is nothing new or unexpected about this result. The initial perturbations of the spatial curvature, as seen from the Hamiltonian constraint (3.6), are of the same order as density perturbations. These perturbations are then enhanced in the course of evolution. This also happens in the FLRW regime, where the Friedmannian evolution also allows for the change of the spatial curvature

[TABLE]

Thus initial perturbations of curvature are enhanced also by the Friedmannian evolution and its evolution is comparable (but not exactly the same) to the one presented in Fig. 3. This shows that a naive expectation that within a spatially flat (globally) universe, the spatial curvature is everywhere the same, is simply not accurate. A more realistic expectation is that within the inhomogeneous regions (i.e. locally) $\Omega_{\cal R}$ (or equivalently $\Omega_{k}$) is of a comparable amplitude as $\Omega_{m}$.

Another significant result, presented in Fig. 3 is that in all cases the contribution from the kinematic backreaction is much smaller than the contribution from the spatial curvature $\Omega_{\cal R}$. This result is consistent with findings reported in Refs. [21, 26, 43, 44, 29], which show that even small perturbations of the kinematic backreaction $\Omega_{\cal Q}$ can lead to non-Friedmannian evolution, or as in the case of this Section to a highly non-linear evolution. This shows that the kinematic backreaction ${\cal Q}$ is merely a part of the backreaction phenomenon and non-linear evolution. Thus, care should be exercised when debating the relevance of backreaction using only arguments based on the amplitude of the kinematic backreaction.

The results reported in this Section were obtained based on a model of a single cosmological environment: a pair of a cosmic void and overdensity, with the size of inhomogeneity below 100 Mpc. Beyond that scale the system is homogeneous. The reason for this is that the size of inhomogeneous structures within the Szekeres model increases with radius $r$, so in order to model realistic structures one has to limit the number of structures to 2 or 3 [45] (but see [46, 47] for an alternative approach). So beyond 100 Mpc in order to eliminate extremely elongated inhomogeneities the model is set to be homogeneous. As a consequence on the scale of 100 Mpc (and beyond) this model is not suitable to study backreaction.

5 Global evolution and emergence of global spatial curvature within the ensemble of Szekeres models

The model considered in Sec. 4 can only be applied to study inhomogeneities on small scales, i.e. much smaller than the scale of homogeneity (i.e. $\sim 100$ Mpc). In order to study backreaction on a much larger, global scale one needs to implement a different approach to the Szekeres model, such as for example the ensemble approach. In this approach we consider an ensemble of Szekeres worldlines, where each worldline is a separate Szekeres model specified by functions (3.9)–(3.13) with the free parameters set to

[TABLE]

where ${\cal U}_{[0-1]}$ is a random number between 0 and 1 (uniform distribution).

For a large number of realisations, the above prescription provides a wide range of inhomogeneous structures from void-like structures (as in Sec. 4) to systems with central overdensity and an adjacent void (as for example in Ref. [45]).

Two configurations are investigated:

Swiss-Cheese-type configuration

The Swiss-Cheese configuration has a lattice of a homogeneous FLRW regions and inhomogeneities smoothly approach the FLRW background (as shown in Fig. 1). Through the homogeneous lattice inhomogeneities are joined to other inhomogeneities. The position of a worldline (with respect to the inhomogeneity defined by (5.1)) is randomly selected by randomly generating the pseudo-Cartesian coordinates $x,y,z$

[TABLE]

which are then used to evaluate the Szekeres coordinates

[TABLE]

With a sufficiently large number of worldlines, this approach reproduces a Swiss-Cheese-type configuration, i.e. each type of inhomogeneity with its asymptotically-approaching FLRW region is well mapped. 2. 2.

Styrofoam-type configuration

The Styrofoam-type configuration consists of densely packed closed-cell structures that do not exhibit any fixed FLRW lattice. The inhomogeneity is still defined in the same way as above, i.e. by (5.1), but the asymptotically-approaching FLRW region is excluded from the Monte Carlo simulation by selecting only worldlines that are close to the central inhomogeneity, with coordinates selected as

[TABLE]

where

[TABLE]

The difference between this type of configuration and the Swiss-Cheese-type configuration is that the inhomogeneous regions are only mapped in their central parts and there is no asymptotically-approaching FLRW region. On one hand this is an advantage – no FLRW lattice. On the other hand, this means that inhomogeneities are just stuck together without any attempt to smoothly join them together, which in the Swiss-Cheese configuration is obtained via the FLRW lattice.

The ensemble consists of $10^{7}$ different and independent worldlines – since the model is silent ($H_{ab}=0$, $\nabla_{a}p=0$, $\pi_{ab}=0$, and $q_{a}=0$) there is no commutation between the worldlines and therefore each worldline evolves independently. The volume average of a function $A$ within the ensemble of the Szekeres models is

[TABLE]

where $A_{n}$ is the value of a function for a specific worldline, and $v_{n}$ is the volume around each worldline

[TABLE]

where ${\rm\det}\,g$ is a determinant of the metric (3.1). Finally, the size of the domain is

[TABLE]

For $10^{7}$ worldlines, the comoving size of the domain of the ensemble defined by (5.1) is approximately $1050$ Mpc.

As in Sec. 4, the initial instant is set to be the last scattering instant. The state of the ensemble at the initial instant presented in Fig. 4, shows that initially the model is quite homogeneous, with average density $\langle{\rho}\rangle$ and domain expansion rate $H_{\cal D}$ quickly approaching the $\Lambda$CDM values, and the spatial curvature and kinematic backreaction approaching $\,0\,$ when averaged over a sufficiently large domain. As seen from Fig. 5, after 13.8 Gyr of evolution the statistical homogeneity is still present at scales beyond 100 Mpc — the averaging over random domains of size ${\cal D}>100{\rm~{}Mpc}$ produces the same results (i.e. the cosmic variance is negligible small at these scales). However, the mean values depend on the global model, and the results differ between the Swiss-Cheese-type and Styrofoam-type models.

The mean evolution within the Swiss-Cheese-type model follows the background $\Lambda$CDM model: average density $\langle{\rho}\rangle=\rho_{0}$, domain expansion parameter $H_{\cal D}=H_{0}$, the model is practically spatially flat with $\Omega_{\cal R}=2.4\times 10^{-4}$, and the kinematic backreaction is negligibly small $\Omega_{\cal Q}=-1.6\times 10^{-5}$. However, within the Styrofoam-type model the average density falls below the $\Lambda$CDM, i.e. $\langle{\rho}\rangle=0.88\,\rho_{0}$ and the domain expansion parameter is faster than in $\Lambda$CDM background $H_{\cal D}=1.03\,H_{0}$. In addition, the mean spatial curvature is $\Omega_{\cal R}=0.11$, and the kinematic backreaction is $\Omega_{\cal Q}=-1.1\times 10^{-2}$, which shows significant deviation from the Friedmannian evolution of the background model.

The explanation for the difference between the global behaviour of the Swiss-Cheese-type model and the Styrofoam-type model is presented in Fig. 6. Figure 6 shows the volume fraction of various types of structures within the Styrofoam-type and Swiss-Cheese-type models. The underdense fraction is defined as

[TABLE]

where $V_{underdense}$ is volume occupied by regions with $\rho<0.9\,\rho_{\Lambda CDM}$. The overdense fraction is defined as

[TABLE]

where $V_{overdense}$ is volume occupied by regions with $\rho>1.1\,\rho_{\Lambda CDM}$. Finally, the lattice fraction is defined as

[TABLE]

where $V_{lattice}$ is volume occupied by regions with $0.99\,\rho_{\Lambda CDM}<\rho<1.01\,\rho_{\Lambda CDM}$.

As seen the Swiss-Cheese model is dominated by the asymptotically-homogeneous regions that only deviate by less than $1\%$ from the $\Lambda$CDM model. Such regions occupy more than $50\%$ of the total volume of the Swiss-Cheese-type model, while the contribution from the underdense ($\rho<0.9\,\rho_{\Lambda CDM}$) and overdense ($\rho>1.1\,\rho_{\Lambda CDM}$) is merely at the percent level, and the rest of the volume is occupied by almost $\Lambda$CDM-like regions. Consequently, the averages and the mean global evolution follows closely the $\Lambda$CDM model. Contrary to the Swiss-Cheese model, within the ensemble of the Szekeres models of the Styrofoam type, the Monte Carlo simulation only probes the central parts of inhomogeneity and avoids the asymptotic FLRW regions. As a result, the volume is dominated by underdense regions with $f_{u}=0.75$ at the present-day instant, and the overdense regions occupy $20\%$ of the total volume. The reason why underdense regions dominate the volume is simply because voids expand faster than the overdense regions, and so they quickly start to occupy larger factions of the volume. Interestingly, by comparing the time scales of Fig. 6 to Fig. 1 one can see that the instant when voids start to dominate the total volume is when they enter the non-linear regime.

The evolution of the cosmic quartet (2.20) is presented in Fig. 7. Not surprisingly, since the volume of the Swiss-Cheese-type model is dominated by the $\Lambda$CDM-like regions, the evolution of the cosmic quartet follows the $\Lambda$CDM behaviour. As for the Styrofoam-type model, which volume is dominated by voids, the evolution of the cosmic quartet, when compared with Fig. 3, resembles the void-like behaviour. It is worth pointing out that within the Styrofoam-type model there is no fixed background. The mean density evolves slightly differently than in the $\Lambda$CDM model, which was used to specify the model at the initial instant. In particular, the evolution of the mean spatial curvature is unlike in the $\Lambda$CDM which is spatially flat.

The emergence of the mean global spatial curvature (within the Styrofoam-type model) has an intuitive explanation. Within this model underdense regions dominate from the start, but since density contrast of overdense regions increase at a slightly higher rate – this is because, as seen from eq. (2.3) both density and shear negatively contribute to $\dot{\Theta}$ – thus overdense regions more quickly pass the $10\%$-threshold used in the definition of $f_{o}$ and $f_{u}$, and thus in the left panel of Fig. 6 it looks like overdense regions slightly dominate in the first 500 My of evolution. At that instant the spatial curvature of the Styrofoam-type model is $\Omega_{\cal R}\approx 5\times 10^{-3}$. After that instant $f_{u}>f_{o}$ and soon afterwards the amplitude of $\Omega_{\cal R}$ becomes non-negligible. This phenomenon can also be understood in terms of the analysis presented in Ref. [48], which focused on the dynamical system of eqs. (2.15)–(2.17). The analysis of the dynamical system (2.15)–(2.17) showed that even a tiny perturbation in ${\cal Q}_{\cal D}$ can drive the system into another basin of attraction that is dominated by the averaged curvature [48]. This property of a global gravitational instability of the Friedmannian model has been identified as the reason for large curvature deviations. The Styrofoam-type model based on the ensemble of the Szekeres models provides an explicit realisation of this scenario.

6 Conclusions

The models considered in this paper (local environment and global ensemble) addressed the issue of the amplitude of backreaction, and its impact on the evolution of a cosmological system. The analysis was based on the exact, cosmological, and non-symmetrical solutions of the Einstein equations, i.e. the Szekeres model, and it provided examples of relativistic models with non-vanishing backreaction (cf. [14]). The obtained results show that the non-linear evolution and backreaction are closely associated with the spatial curvature. The growth of inhomogeneities cannot be separated from the growth of the spatial curvature, which variation across the Universe is comparable with the variation in the matter field, and in addition the global mean spatial curvature does not necessarily average out to zero. This implies that the emergence of the Cosmic Web in the real Universe should also be associated with the emergence of the mean spatial curvature.

There is nothing “exotic” about the obtained results, and the second-order effects do not “magically” appear with an amplitude a few orders of magnitude larger than the first order effects. These findings can easily be understood in a simple and logical manner. Initially, at the last scattering instant, the Universe is fairly homogeneous with only small perturbations present. As long as the growth of structures is linear and perturbations are small, the differences in the expansion rates are negligible and both types of regions: underdense and overdense expand at a similar, background rate. Once the growth of cosmic structure is non-linear (e.g. matter shear $\sigma^{2}$ is no longer negligible) the expansion rate of overdense regions efficiently slows down (cf. eq. (2.3)). When that happens the volume of the Universe becomes dominated by voids (cf. Fig. 6). Once the volume of the Universe starts to be dominated by voids then the total volume (eg. within the cosmic horizon) increases faster than in the FLRW model, and since matter is conserved, the mean density is lower than in the $\Lambda$CDM model. This explains why, as seen in Fig. 5, the mean density is below the $\Lambda$CDM matter density, and the expansion rate is slightly higher than $H_{0}$. As for the spatial curvature, the density perturbations are coupled with expansion rate perturbations and curvature perturbations (cf. eq. 3.6)). Thus, initially at the last scattering instant, not only density but also tiny curvature perturbations are present. In the course of evolution, these perturbations also grow (in the FLRW regime, $\Omega_{\cal R}\to\Omega_{k}=-k/\dot{a}^{2}$), which leads to large variations of spatial curvature across cosmic structures (cf. Fig. 2). In addition, since the volume of the Universe in the non-linear regime is dominated by voids, the mean global spatial curvature evolves from $\Omega_{\cal R}=0$ to $\Omega_{\cal R}\approx 0.1$ (cf. Fig. 7).

These findings were obtained based on the Styrofoam-type model. In the Swiss-Cheese model, at the local scales the picture is similar (cf. Figs. 1–3) — i.e. non-linear evolution leads to large differences between expansion rates of underdense and overdense regions, as well as, large variations of the spatial curvature — however on global scales the results are different, which can be linked to a fact that in the Swiss-Cheese model the volume is dominated by the $\Lambda$CDM-like regions. Since the volume of the real Universe is in fact dominated by cosmic voids [49], it is reasonable to conclude that the Styrofoam-type model is more realistic than the Swiss-Cheese model. As a results, one can expect that the global mean spatial curvature of our Universe should also be dominated by voids, and subsequently it should deviate from zero in the low-redshift Universe, leading to its evolution from spatial flatness ($\Omega_{\cal R}=\Omega_{k}=0$) in the early Universe to a negative spatial curvature ($\Omega_{\cal R}\sim\Omega_{k}>0$) at the present day epoch, in a similar manner as presented in Fig. 7.

It is interesting to point out that, in fact, the analysis of low-redshift data, such as a supernova data alone (without combining it with the CMB) implies large spatial curvature, i.e. $\Omega_{k}\approx 0.2$ and only after inclusion of the CMB reduces to $\Omega_{k}=0.005\pm 0.009$ [50]. Also, supernova data alone point towards slightly lower values of $\Omega_{m}$ compared to the CMB constraints, which could be understood, not just in terms of various systematics, but also partly in terms of findings presented in Fig. 7. In addition, there is a known tension between the values of $H_{0}$ derived from the CMB [51] and the local measurements [52] which again, apart from various systematics [53], could also be partly explained in terms of findings presented in Fig. 5.

The models like the one presented in this paper, i.e. Styrofoam-type model, are not perfect realisations of our Universe. There are a number of limitations, for example the lack of rotation $w_{ab}$ excludes presence of virialised regions; the lack of magnetic Weyl tensor $H_{ab}$ excludes multiple eigenvalues of shear; and lack of heat flow $q_{a}$ excludes energy transfers from one cosmic cell to another. This all means that one should exercise caution when drawing conclusion as to the properties of our Universe. Still, the obtained results should encourage further studies of the cosmological backreaction and development of numerical cosmology towards realistic models of the Universe [54, 55, 56]. Also, these results should encourage cosmologists to think outside the box, especially when dealing with low-redshift data, and allow for example for $\Omega_{k}\neq 0$, even if at the CMB the Universe was spatially flat222Separation of the low-$z$ from high-$z$ data has already been implemented by some cosmologists, for example in Ref. [57] the supernova analysis was conducted using different sets of cosmological parameters for low-$z$ and high-$z$ data; in Ref. [58] the CMB data was analysed independently from the low-$z$ cosmology that enters via the distance to the last scattering surface.. Allowing for $\Omega_{k}$ to be a free parameter at low-$z$ and independent from high-$z$ constraints, could in principle reduce some of the tensions and inconsistencies in the observed data [59, 60], but the actual analysis remains to be done. In addition, cosmologists should also aim at directly measuring the curvature of the low-redshift Universe to check if it indeed deviates from the CMB constraints [61].

Acknowledgments

I would like to thank Thomas Buchert, Jan Ostrowski, Boudewijn Roukema, David Wiltshire, and anonymous referee for their comments and suggestions. This work was supported by the Australian Research Council through the Future Fellowship FT140101270.