Modelling the emergence of cosmic anisotropy from non-linear structures

In order to describe a cosmology with a preferred space-like direction, such as the symmetry axis of an LRS model, it is useful to perform a semi-tetrad decomposition of the field equations with respect to both a preferred time-like vector n and a preferred space-like vector m (satisfying $m_a n^a = 0$, $n_a n^a = -1$, and $m_a m^a = 1$). One can then decompose all vector and tensor quantities with respect to n and m, while computing the Ricci and Bianchi identities for n and m gives rise to a set of equations which is physically equivalent to Einstein's equations. This $1+1+2$ formulation was first introduced by Greenberg [39], developed in [40–46], and worked out in full by Clarkson [47] and Keresztes et al [48]. We will introduce the salient features of this approach here, and explain why they are of use for our approach.

Let us suppose that we have a time-like vector n. One can then define the projection tensor $f_{ab} = g_{ab} + n_a n_b$, as well as derivative operators along and orthogonal to n:

$$\begin{equation} \dot{\phantom{a}} \equiv n^a \nabla_a \qquad \textrm{and} \qquad D_a T_{b c} \equiv f_a^{\phantom{a} d} f_b^{\phantom{b} e} f_c^{\phantom{c} f} \nabla_d T_{ef} \, , \end{equation} \tag{ 1 }$$

with obvious generalizations to tensors of other ranks. These allow us to decompose the covariant derivative of n into kinematic quantities as [49]

$$\begin{equation} \nabla_a n_b = -n_a \dot{n}_b + \frac{1}{3}\Theta f_{ab} + \sigma_{ab} + \eta_{cab}\omega^c, \end{equation} \tag{ 2 }$$

where $\dot{n}_a = n^b \nabla_b n_a$, $\Theta = D_a n^a$, $\sigma_{ab} = D_{\langle a}n_{b\rangle}$ and $\omega_a = \frac{1}{2}\eta_{abc}D^b n^c$. These are, respectively, the acceleration of observers with 4-velocity n, the isotropic expansion of the space-like 3-surfaces orthogonal to n, and the volume-preserving shear and vorticity. In writing these expressions down, we have introduced angled brackets around indices to denote a tensor that is projected orthogonally to n and made symmetric and trace-free, such that e.g.

$$\begin{equation} \sigma_{ab} = D_{\langle a}n_{b\rangle} \equiv \left(f_{(a}^{\ c} f_{b)}^{\ d} - \frac{1}{3}f_{ab}f^{cd}\right)\nabla_c n_d \, . \end{equation} \tag{ 3 }$$

We have also used $\eta_{abc} \equiv n^d \eta_{dabc}$ as the alternating tensor on the 3-surfaces orthogonal to n.

Let us now further suppose that we have a space-like vector m, which is orthogonal to n. The covariant derivative of m can be decomposed as

$$\begin{equation} \nabla_a m_b = - n_a \dot{m}_b + D_a m_b + n_b m^c D_a n_c + n_a n_b m^c \dot{n}_c \, , \end{equation} \tag{ 4 }$$

and a projection tensor onto the two-dimensional 'screen spaces' orthogonal to both n and m can be defined as

$$\begin{equation} M_{ab} = f_{ab} - m_a m_b = g_{ab} + n_a n_b - m_a m_b \, . \end{equation} \tag{ 5 }$$

It can be seen that the projected time and space derivatives of m can be written as

$$\begin{equation} \dot{m}_a = \mathcal{A} \, n_a + \alpha_a \qquad \textrm{and} \qquad D_a m_b = m_a a_b + \frac{1}{2}\phi M_{ab} + \xi\epsilon_{ab} + \zeta_{ab}\, , \end{equation} \tag{ 6 }$$

where $\epsilon_{ab} = \eta_{abc}m^c$ is the alternating tensor on the screen space. The new quantities introduced in equation (6) are $\alpha_a = M_{ac} n^b \nabla_b m^c$ and $\mathcal{A} = -n^a \dot{m}_a = m_a \dot{n}^a$, which are the projections of $\dot{m}_a$ orthogonal and parallel to n, and $a_b = m^c D_c m_b$, $\phi = M^{ab} D_a m_b$, $\zeta_{ab} = \left[M^{\: c}_{(a}M^{\: d}_{b)} - \frac{1}{2}M_{ab}M^{cd}\right]D_{\langle c} m_{d \rangle}$, and $\xi = \frac{1}{2}\epsilon_{ab}M^{ac}M^{bd} D_c m_d$, which are the non-geodesy of m, and the expansion, shear and twist of the screen space along m. We use the names 'expansion' and 'shear' at this point in analogy with the quantites defined for n, despite m being a space-like vector. This is done for linguistic ease.

To go further, we note that all 3-vectors and 3-tensors defined with respect to n can be further projected parallel and orthogonal to m. For 3-vectors $v_a = v_{\langle a \rangle} = f_a^{\ b}v_b$, this can be written as

$$\begin{equation} v_a = V m_a + V_a \, , \end{equation} \tag{ 7 }$$

where $V = v_b m^b$ is the projection of v_a parallel to m_a , and $V_a = M_a^{\ b}v_b$ is the projection of v_a into the screen spaces. Similarly, for projected, symmetric and tracefree 3-tensors $t_{ab} = t_{\langle ab\rangle}$, we can write

$$\begin{equation} t_{ab}= \tau\left(m_a m_b - \frac{1}{2}M_{ab}\right) + 2 \tau_{(a}m_{b)} + \tau_{ab} \, , \end{equation} \tag{ 8 }$$

where $\tau = t_{ab} m^a m^b$, $\tau_a = M_{ab}m_c t^{bc}$ and $\tau_{ab} = [M^{\: c}_{(a}M^{\: d}_{b)} - \frac{1}{2}M_{ab}M^{cd}]t_{cd}$ are the projections of t_ab twice parallel, once parallel and once orthogonal, and twice orthogonal to m, respectively.

Using this approach, we can decompose the acceleration and vorticity 3-vectors associated with n in the following way:

$$\begin{align} \dot{n}_a = \mathcal{A}m_a + \mathcal{A}_a \qquad \textrm{and} \qquad \omega_a = \Omega m_a + \Omega_a \, . \end{align} \tag{ 9 }$$

The only other kinematic quantity is the shear of n, which can be decomposed as

$$\begin{align} \sigma_{ab} = \Sigma \left(m_a m_b - \frac{1}{2}M_{ab}\right) + 2\Sigma_{(a}m_{b)} + \Sigma_{ab} \, . \end{align} \tag{ 10 }$$

The quantities $\mathcal{A}$, Ω and Σ are covariantly defined scalars, while $\mathcal{A}_a$, $\Omega_a$, $\Sigma_a$ and $\Sigma_{ab}$ are covariantly defined 2-vectors and tensors defined in the screen spaces orthogonal to n and m.

We have so far considered the decompositions of the derivatives of n and m. For a complete system we will also require quantities associated with the Ricci and Weyl curvature of the space-time. The former of these can be related to the matter content of the space-time using Einstein's equations, and therefore written in terms of a stress-energy tensor that can be decomposed with respect to n as

$$\begin{equation} T_{ab} = \mu n_a n_b + p f_{ab} + 2 q_{(a}n_{b)} + \pi_{ab}, \end{equation} \tag{ 11 }$$

where $\mu = T_{ab}n^a n^b$, $p = \frac{1}{3}T_{ab}f^{ab}$, $q_a = -T_{bc}n^b f_a^{\ c}$, and $\pi_{ab} = T_{\langle ab \rangle}$ are the energy density, isotropic pressure, momentum density and anisotropic stress, respectively. The last two of these can be further decomposed with respect to m to get

$$\begin{align} q_a = Q m_a + Q_a \qquad \textrm{and} \quad \pi_{ab} = \Pi \left(m_a m_b - \frac{1}{2}M_{ab}\right) + 2\Pi_{(a}m_{b)} + \Pi_{ab} \, , \end{align} \tag{ 12 }$$

where all quantities are defined as in equations (7) and (8), above. Finally, the Weyl tensor can be decomposed with respect to n into an 'electric' and a 'magnetic' part, as

$$\begin{equation} E_{ab} = C_{acbd}n^c n^d \qquad \textrm{and} \qquad H_{ab} = \frac{1}{2}\eta_{efda}C^{ef}_{\quad bc} n^c n^d \, , \end{equation} \tag{ 13 }$$

which can be further decomposed with respect to the space-like vector m as

$$\begin{align*} E_{ab} = \mathcal{E} \left(m_a m_b - \frac{1}{2}M_{ab}\right) + 2\mathcal{E}_{(a}m_{b)} + \mathcal{E}_{ab}, \quad H_{ab} = \mathcal{H} \left(m_a m_b - \frac{1}{2}M_{ab}\right) + 2\mathcal{H}_{(a}m_{b)} + \mathcal{H}_{ab}\, . \end{align*}$$

Again, all quantities in these expressions are defined as in equations (7) and (8). This completes the decomposition of all required quantities with respect to both n and m.

The full and irreducible set of quantities describing the kinematics of n and m, the matter variables and the Weyl curvature is therefore given as follows:

$$\begin{equation} \textrm{2-scalars:} \quad \left \lbrace \Theta, \mathcal{A}, \Sigma, \Omega, \phi, \xi, \mathcal{E}, \mathcal{H}, \mu, p, Q, \Pi \right \rbrace \end{equation} \tag{ 14 }$$

$$\begin{equation} \textrm{2-vectors:} \quad\left \lbrace \mathcal{A}_a, \Sigma_a, \Omega_a, \alpha_a, a_a, \mathcal{E}_a, \mathcal{H}_a, Q_a, \Pi_a \right\rbrace \end{equation} \tag{ 15 }$$

$$\begin{equation} \textrm{2-tensors:} \quad\left \lbrace \Sigma_{ab}, \zeta_{ab}, \mathcal{E}_{ab}, \mathcal{H}_{ab}, \Pi_{ab} \right \rbrace. \end{equation} \tag{ 16 }$$

The reader is referred to [47, 48] for the further details of these decompositions, and to appendix for a full list of all quantities defined above. For the remainder of this paper we will take $\Omega = 0 = \Omega_a$, so that n is hypersurface-orthogonal. This is required for having a foliation of space-like 3-surfaces over which averaging can be performed. It is not a restriction on the properties of matter.

We have so far not specified anything about the time-like vector n, other than it being hypersurface-forming. In the sections that follow we will consider specific choices of this vector, which are orthogonal to foliations with certain special properties. Here, we will simply note that it will be of interest to be able to transform between foliations, and present some of the mathematics that is required to do so.

In order to understand this, let us consider a boost from n to some second time-like vector u. The components of these vectors are related by

$$\begin{equation} u_a = \gamma\left(n_a - v_a\right) \qquad \textrm{and} \qquad n_a = \gamma\left(u_a + w_a\right) \, , \end{equation} \tag{ 17 }$$

where $v_a = \gamma^{-1}f_a^{\ b} w_b$ and $\gamma = \left(1-v^2\right)^{-1/2} = \left(1-w^2\right)^{-1/2}$. Here the boost vectors v^a and w^a exist respectively in the surfaces orthogonal to n and u, such that $n_a v^a = u_a w^a = 0$. The projection tensor into the hypersurfaces orthogonal to u (if they exist) can then be defined as $h_{ab} = g_{ab} + u_a u_b$, which allows us to write $w_a = \gamma^{-1} h_a^{\ b} v_b$. The same projection tensor also allows us to define a new vector $\mathbf{\tilde{m}}$ with components

$$\begin{equation} \tilde{m}_a = k \, h_a^{\ b} m_b \, , \qquad \textrm{such that} \qquad m_a = \tilde{k} \, f_a^{\ b} \tilde{m}_b \,, \end{equation} \tag{ 18 }$$

where $k = \left[1+\left(\gamma m_c v^c\right)^2\right]^{-1/2}$ and $\tilde{k} = \left[1+\left(\gamma \tilde{m}_c w^c\right)^2\right]^{-1/2}$. These are the components of the preferred space-like vector m projected orthogonal to u, as illustrated in figure 1.

Zoom In Zoom Out Reset image size

The whole system of $1+1+2$-equations from section 2.1 can now be set up with respect to either n or u, with the preferred space-like vector projected according to equation (18). After we boost from n to u, we wish to express the scalars in the new foliation in terms of those in the old one, plus any terms arising from the boost 3-velocity v that connects them. If $\left\vert \mathbf{v} \right\vert \ll 1$, which is expected if, for example, n corresponds to the canonical foliation of Friedmann–Lemaître–Robertson–Walker (FLRW) spacetime, and u is the matter 4-velocity in that spacetime, then one can expand these expressions order by order in v, leading to vast simplifications. However, in the general case, one need not make this approximation. Performing the full calculations, one finds the following results:

The isotropic expansion $\tilde{\Theta} = h^{ab} \nabla_a u_b$ is related to $\Theta = f^{ab}\nabla_a n_b$ by

$$\begin{equation} \tilde{\Theta} = \gamma\left[\Theta - \kappa - v_a\left(\dot{n}^a + \gamma^2 \dot{v}^a\right) + \gamma^2\left(\frac{1}{3}\kappa v^2 + \beta_{ab} v^a v^b\right)\right], \end{equation} \tag{ 19 }$$

where $\kappa = D_a n^a$ and $\beta_{ab} = D_{\langle a} v_{b\rangle}$. This equation shows how the expansion of the rest spaces defined with respect either n or u are related to each other, and explicitly provides the additional terms that contribute to the expansion after performing a boost.

The scalar shear $\tilde{\Sigma} = \left(h_a^{\ c} h_b^{\ d} - \frac{1}{3}h_{ab}h^{cd}\right)\nabla_c u_d \tilde{m}^a \tilde{m}^b$ is related to Σ by

$$\begin{align} \tilde{\Sigma} &= k^2\gamma \Sigma + \frac{\gamma}{3}\left(\Theta-\kappa\right)\left(k^2\gamma - 1\right) + \frac{\gamma}{3}v_a\left(\dot{n}^a+\gamma^2\dot{v}^a\right) - \frac{\gamma^3}{3}\left(\frac{\kappa}{3}v^2 + \beta_{ab}v^a v^b\right)\nonumber\\ & \quad + k^2\gamma \Bigg[-\beta_{ab}m^am^b -\left(\gamma m_a v^a\right)^2 +\gamma^2 m^a v^b\left(m_c v^c - 1\right)\left(\beta_{ab}+W_{ab}\right) - \gamma^2 m_a v^a\mathcal{A}\nonumber\\ & \quad + \gamma^2 m_b v^b\left(\frac{1}{3}\left(3\Sigma+\Theta-\kappa\right) m_a v^a + \Sigma_a v^a + \dot{v}_a \left(m^a + \gamma^2 m_c v^c v^a\right)\right)\nonumber\\ & \quad - \gamma^4 \left(m_c v^c\right)^2\left(\frac{1}{3}\kappa v^2 + \beta_{ab} v^a v^b\right) \Bigg], \end{align} \tag{ 20 }$$

and the scalar acceleration, $\tilde{\mathcal{A}} = \tilde{m}^a u^b \nabla_a u_b$, and the scalar vorticity $\tilde{\Omega} = \frac{1}{2}\tilde{m}^a \eta_{dabc}u^d \nabla^b u^c$, are

$$\begin{equation} \tilde{\Omega} = -k\gamma^2\left[W_a\left(m_a + \gamma^2 m_bv^b v^a\right) - \frac{1}{2}\gamma^2 m_a v^a \chi_{bc}W^{bc}\right], \end{equation} \tag{ 21 }$$

and

$$\begin{equation} \tilde{\mathcal{A}} = k\gamma^2\left[\mathcal{A} - m_a \dot{v}^a - \gamma^2 m_a v^a v_b \dot{v}^b - \frac{1}{3}m_a v^a\left(3\Sigma + \Theta - \kappa\right) - \Sigma_a v^a + v^a m^b\left(\beta_{ab} + W_{ab}\right)\right] , \end{equation} \tag{ 22 }$$

where we have used that $\Omega = \frac{1}{2} m^a \eta_{dabc}n^d \nabla^b n^c = 0$. One sees immediately that if the vorticity W_a associated with the relative velocity v vanishes, then $\tilde{\Omega} = 0$ as well.

Meanwhile, the twist and expansion of the screen spaces along $\tilde{\mathbf{m}}$, $\tilde{\xi} = \frac{1}{2}\tilde{\epsilon}_{ab} \tilde{M}^{ac}\tilde{M}^{bd}\tilde{D}_c \tilde{m}_d$ and $\tilde{\phi} = \tilde{M}^{ab}\tilde{D}_a \tilde{m}_b$, are related to ξ and φ by

$$\begin{equation} \tilde{\xi} = k^2\gamma\left[\xi + \epsilon_{ab}v^a\left(\alpha^b + \Sigma^b\right)\right]. \end{equation} \tag{ 23 }$$

and

$$\begin{align} \tilde{\phi} &= k\phi + k\mathcal{A} - \tilde{A} - k\gamma m_a v^a \tilde{\Theta} + \left(m^a - \gamma^2 m_b v^b\left(n^a - v^a\right)\right)\nabla_a k\nonumber\\ & \quad + k\gamma^2 \left[-\alpha_a v^a + m_a v^a\left(a_b v^b + \frac{\kappa}{3}\right) - m_a\dot{v}^a + \left(\frac{\phi}{2} M_{ab} + \zeta_{ab}\right)v^a v^b + \left(\beta_{ab} + W_{ab}\right) v^a m^b \right] \nonumber\\ & \quad+ k\gamma^4 m_c v^c\left(-v_a\dot{v}^a + \frac{\kappa}{3}v^2 + \beta_{ab}v^a v^b\right), \end{align} \tag{ 24 }$$

$$\begin{align} \nonumber \textrm{where} \quad \nabla_a k X& = - \gamma^2 k^3 m_c v^c \Bigg[\gamma^2 m_b v^b v^d\left(-n_a \dot{v}_d + \beta_{ad} + W_{ad}\right) + \frac{\gamma^2}{3}m_b v^b \kappa v_a - \alpha_b v^b n_a\nonumber\\ & \quad + v^b\left(\frac{\phi}{2}M_{ab} + \zeta_{ab} + \xi \epsilon_{ab} + m_a a_b\right) - \dot{v}_b m^b n_a + \frac{\kappa}{3}m_a + m^b\left(\beta_{ab} + W_{ab}\right)\Bigg]. \end{align}$$

For the matter variables, the energy density $\tilde{\mu} = T_{ab}u^a u^b$, isotropic pressure $\tilde{p} = \frac{1}{3}T_{ab}h^{ab}$, and scalar heat flux $\tilde{Q} = -T_{bc}u^b h^c_{\ a} \tilde{m}^a$ are given by

$$\begin{align} \tilde{\mu} &= \gamma^2\left(\mu + pv^2 + 2q_a v^a + \pi_{ab}v^av^b\right) , \end{align} \tag{ 25 }$$

$$\begin{align} \tilde{p} &= \left(1+\frac{1}{3}\gamma^2 v^2\right) p + \frac{1}{3}\gamma^2 v^2 \mu + \frac{1}{3}\gamma^2\left(2q_a v^a + \pi_{ab}v^a v^b\right) , \end{align} \tag{ 26 }$$

$$\begin{align} \tilde{Q} &= k\gamma\left[Q + m_c v^c\left(\gamma^2\left(\mu+p+2q_a v^a + \pi_{ab}v^a v^b\right) + \Pi\right) + \Pi_a v^a \right] , \end{align} \tag{ 27 }$$

and the scalar anisotropic stress $\tilde{\Pi} = \left(h_a^{\ c} h_b^{\ d} - \frac{1}{3}h_{ab}h^{cd}\right)T_{cd} \tilde{m}^a \tilde{m}^b$ is given by

$$\begin{align} \tilde{\Pi} &= \gamma^2 \mu\left(\left(k\gamma m_a v^a\right)^2 - \frac{1}{3}v^2\right) + p \left[ k^2 - 1 + 2\left(k\gamma m_a v^a\right)^2 + \frac{1}{3}\gamma^2 v^2\left(3\left(k\gamma m_a v^a\right)^2 - 1\right)\right] \nonumber\\ & \quad+ 2k^2\gamma^2 m_a v^a Q + k^2 \Pi\left(1+2\left(\gamma m_a v^a\right)^2\right) + \frac{\gamma^2}{3}\left(3\left(k\gamma m_c v^c\right)^2-1\right)\left[2q_a v^a + \pi_{ab} v^a v^b\right]. \nonumber\\ \end{align} \tag{ 28 }$$

Finally, one can compute the scalar parts of the Weyl curvature for the new congruence. For the scalar parts of the electric and magnetic Weyl tensors this gives [50]

$$\begin{align} \tilde{\mathcal{E}} &= k^2\gamma^2\left[\left(1+v^2 - \left(m_a v^a\right)^2 - \frac{1}{2}M_{ab}v^a v^b\right)\mathcal{E} + 2\chi^{ab}v_c m_a\mathcal{H}_b + \mathcal{E}_{ab}v^a v^b\right] , \end{align} \tag{ 29 }$$

$$\begin{align} \tilde{\mathcal{H}} &= k^2\gamma^2\left[\left(1+v^2 - \left(m_a v^a\right)^2 - \frac{1}{2}M_{ab}v^a v^b\right)\mathcal{H} - 2\chi^{ab}v_c m_a\mathcal{E}_b + \mathcal{H}_{ab}v^a v^b\right]. \end{align} \tag{ 30 }$$

This completes the transformation rules for all of our covariantly defined scalars. The 2-vector and tensor quantities that result from the $1+1+2$ decomposition can also be transformed between foliations, but only appear in our averaged description via 'backreaction' terms. We therefore omit presenting their transformation rules here.

In what follows, we will use often use u to refer to the flow lines of an irrotational dust fluid. In the presence of such matter, this vector is uniquely defined and therefore also provides a unique foliation on which to perform averaging. This is the foliation that is most often used in the literature on this subject, and is the one chosen for the standard approach to formulating Buchert's equations [37]. We will then use n to refer to any other foliation, chosen (for example) for observational or geometric reasons. Such freedom is required if we are to allow for the possibility of bulk flows, as the 3-velocity of matter vanishes identically if we choose the foliation to be induced by the flow of dust u. It also allows us the freedom to re-foliate in situations in which the description of space-time might break down if it were specified by the fluid flow, as would happen, for example, in perturbed FLRW cosmologies that contain non-linear structures, or if one wished to consider fluids with non-zero vorticity or caustics in their flow lines.

Averaging in general relativity is a notoriously difficult problem to define, primarily as the theory is defined in terms of tensors, which cannot be straightforwardly compared at different points in a curved space-time [51]. This makes defining averages difficult, as quantities that exist at different points in an averaging domain first have to be transported to the same point for comparison, which is a non-unique and path-dependent process. The problem is further exacerbated by the non-linearity of Einstein's equations, which means that the acts of averaging and evolution do not commute. The evolution of an averaged quantity cannot therefore be guaranteed to be the same as an average taken at some future time.

The first of these problems has been the focus of a large number of different studies, using a variety of different approaches (see e.g. [52–57]). The benefits and drawbacks of each of these approaches are a subject of much debate, but in order to make progress one can consider a more limited case: the averaging of scalars. Such quantities can be readily compared at different points in averaging domains, and it is therefore much more straightforward to compute their averages. Indeed, attempts have been made to construct averaging schemes that completely characterize space-times only in terms of the averages of scalars [58]. In what follows we will proceed by constructing an averaging formalism for the scalars that result from our $1+1+2$ decomposition.

We will be interested in spatial averages of scalar quantities, over the 3-dimensional spatial volumes orthogonal to the irrotational flow lines n. Following Buchert [37], we will define these as follows:

$$\begin{equation} \left\langle S \right \rangle := \frac{1}{V_{\mathcal{D}}}\int_{\mathcal{D}} \mathrm{d}V \: S \: = \: \frac{\int_{\mathcal{D}}\mathrm{d}^3 x\,\sqrt{f} \, S}{\int_{\mathcal{D}}\mathrm{d}^3 x\,\sqrt{f}} \, , \end{equation} \tag{ 31 }$$

where S is the scalar being averaged over the spatial domain $\mathcal{D}$, which has volume $V_{\mathcal{D}} = \int_{\mathcal{D}}\mathrm{d} V$, and where $\mathrm{d} V = \mathrm{d}^3 x\,\sqrt{f}$ is the spatial volume element. In these expressions, we have taken xⁱ to be coordinates on the spatial hypersurface orthogonal to n, and written $f = \mathrm{det}(f_{ab})$ where $f_{ab} = g_{ab}+n_an_b$ is the induced metric. This definition of averaging depends on our choice of foliation and on the selected averaging domain, but is otherwise unambiguous.

The rate of change of this average can be determined by considering that the factor $\sqrt{f}$ represents the volume of an infinitessimal region of space. Assuming the spatial coordinates are fixed by the flow of n, we then have $(\sqrt{f})^{\boldsymbol{\cdot}} = \Theta \sqrt{f}$, where the dot represents the time derivative operator $n^a \nabla_a$ and where $\Theta = \nabla_a n^a$ is the expansion scalar. This allows one to derive Buchert's commutation rule for time derivative and averaging operations [37]:

$$\begin{equation} \left\langle S \right\rangle^{\boldsymbol{\cdot}} - \langle \dot{S}\rangle = \left\langle S \Theta\right\rangle- \left\langle S \right\rangle\left\langle\Theta\right\rangle = \mathrm{Cov}(\Theta,S) \, , \end{equation} \tag{ 32 }$$

where we have defined the covariance as $\mathrm{Cov}(S_1, S_2) : = \left\langle S_1 S_2 \right\rangle - \left\langle S_1\right\rangle \left\langle S_2\right\rangle$. This result relies on the averaging domain $\mathcal{D}$ being propagated by the flow of n. In what follows we will make repeated use of the covariance of scalar quantities, as well as the variance, $\mathrm{Var}(S) : = \mathrm{Cov}(S,S)$.

In equation (32) it is assumed that the lapse function N, between the leaves of the chosen foliation, is independent of spatial position. This is the case if (for example) the time-like 4-vector n is identified with the 4-velocity of a field of irrotational dust u, as in this case $\dot{n}^{a} = \dot{u}^a = 0$. Each leaf in the foliation can then be identified with a single value of the proper time along the world-lines of this dust, and will be the same at all points in each spatial hypersurface. More generally, we will have $\dot{n}_a = D_a \ln N \neq 0$, meaning that the proper time of a family of observers following n is not constant across $\mathcal{D}$. For the purposes of calculating cosmological averages it is desirable to have quantities defined over an averaging domain to have a single value for the time coordinate t. This can be achieved by writing $n_a = - N \partial_a t$ where $N = \dot{t}^{-1}$, such that a derivative with respect to t can be written as $\partial_t S = N \dot{S}$. The commutation rule (32) then becomes

$$\begin{equation} \partial_t\left\langle S\right\rangle - \left\langle \partial_t S\right\rangle = \mathrm{Cov}\left(N\Theta, S\right), \end{equation} \tag{ 33 }$$

which allows us to calculate the evolution of scalar averages as a function of the coordinate time t, once the lapse function N has been determined from the acceleration of n.

One can apply the averaging operation from equation (31) to each of the terms in any scalar equation that is defined locally at each point in space-time. After suitable application of the commutation rule from equation (32) or (33), one can then obtain from that averaged scalar equation a new one that is written entirely in terms averages and their derivatives with respect to time. Equations derived in this way can be thought of as governing the behaviour of averages, rather than the behaviour of locally-defined quantities (as usually occurs in theories of gravity). This is ideal for describing the emergent large-scale properties of a space-time, as is relevant for the study of cosmology. Indeed, the averaged equations that result from this procedure can often be written in a form identical to that of the local gravitational equations of an homogeneous cosmology, with any extra terms that occur being considered to describe the 'back-reaction' effect that results from averaging small-scale structures.

In the standard approach to this problem, three scalar equations are usually considered: the Raychaudhuri equation for the evolution of the expansion, the Hamiltonian constraint equation that relates expansion and energy density, and a contracted Bianchi identity that corresponds to energy conservation. The averages of these are all evaluated in the foliation orthogonal to the flow of dust, and a set of Friedmann-like equations emerge that govern the average expansion [37]. If one wishes to consider the emergence of any anisotropic quantity, however, then this necessarily requires further scalar equations, such as those governing the scalar shear Σ and momentum density Q. This is what we will pursue in the following sections. We note that an interesting step in this direction has already been taken by extending the Buchert formalism to include an evolution equation for the magnitude of the shear, $\sigma^2 = \frac{1}{2} \sigma_{ab}\sigma^{ab}$ [59]. We extend this result further by deriving the full set of averaged scalar equations in the $1+1+2$ formalism. We find a set of equations that can be used to interpret that large-scale cosmological average in terms of an emergent LRS Bianchi model with tilt, and with back-reaction terms that give the effect of small-scale structure on all emergent quantities.

An application of the averaging procedure outlined above requires a choice of foliation, or equivalently a choice of the irrotational time-like 4-vector n. The averages obtained using equation (31) will then correspond to the large-scale properties of the 3-dimensional spaces that constitute the leaves of this foliation. Different foliations will mean that one is considering different 3-dimensional spaces, and hence different averages will be obtained. It is therefore necessary to make sure an appropriate choice of foliation is made, for the situation being considered. This is given further importance by the fact that observers in different frames will infer different cosmological parameters from the Hubble flow around them [21–28].

In general, one might be interested in choosing a foliation that is expected to give results that can be associated with a particular observable [60], that have a particular mathematical or physical meaning associated with them [61], or that are perhaps convenient in some other way [38]. For example, one may wish to construct Hubble diagrams in a frame comoving with the flow of matter [37], or in a frame of 'most uniform Hubble flow' [62–64]. Of course, in order to relate observables corresponding to quantities calculated on different foliations one will need to be able to transform between frames ¹ , as considered in section 2.2.

Let us now consider some specific choices of foliation that one might use:

The comoving foliation exists when the Universe is filled with irrotational dust, and n is chosen to be coincident with the 4-velocity of that dust, u. This choice shares some properties with the 'comoving synchronous' gauge that may be familiar from studies of perturbed FLRW models. It has the distinct benefit of being tied to a physical quantity: the flow of matter. It therefore takes a particularly privileged position in the pantheon of possible choices, but does require the matter flow to be vorticity-free, which is not expected to hold for realistic astrophysical structures. It also cannot account for the existence of any bulk flow, as it corresponds to the choice of frame in which matter is at rest. Nevertheless, this is the standard choice of foliation in much of the literature on mathematical cosmology [49], as well as in many studies of perturbed FLRW.

The gravitational rest frame was defined in [38] to be that in which the magnetic part of the Weyl tensor vanishes:

$$\begin{align*} H_{ab}=0 \,, \end{align*}$$

so that the emergent cosmological model is as Newtonian in character as possible [66]. Cosmologies that obey this condition have been described as 'silent' [67]. Up to second order in perturbations around a Robertson–Walker geometry this choice corresponds to the Poisson gauge [38], but when applied non-perturbatively ends up being a strong restriction on the allowed solutions [68, 69]. It therefore appears to be particularly useful for studies of perturbed FLRW models, as the Poisson gauge is one of the few that is expected to be valid into the regime of non-linear structure formation, where post-Newtonian expansions are required to describe weak gravitational fields [70].

The gravitational wave frame is a more conservative choice defined by the demand that H_ab is divergence-free, such that

$$\begin{align*} D^b H_{ab} = 0 \, . \end{align*}$$

This is a covariant way of stating that the only gravito-magnetic contributions to the curvature measured by an observer comoving with n would be those coming from gravitational waves [71, 72]. It is therefore less restrictive than the condition $H_{ab} = 0$, and typically corresponds to a congruence which is distinct from the flow of dust [73]. It is well-suited to numerical-relativistic cosmological simulations [74], but in a general space-time it is not guaranteed that the frames it picks out are hypersurface forming (though one may be able to take the irrotational part, which would be).

The constant mean curvature foliation has a long history [61], and is defined by the condition that the spatial gradient of the expansion scalar (also known as the 'mean curvature') vanishes:

$$\begin{align*} D_a \Theta = 0\, . \end{align*}$$

This foliation has the property of being unique (under certain circumstances), and having a monotonic variation in the expansion scalar between leaves [75]. It therefore provides a plausible candidate for the provision of a universal, cosmological arrow of time. The literature on this particular choice of foliation has been largely resticted to the literature on mathematical cosmology, but it also exists in perturbed FLRW models under the name of 'uniform expansion' gauge.

The list of choices presented above is by no means exhaustive. Other interesting choices of foliation are specified by the 'zero-shear' condition, $\sigma_{ab} = 0$. This is closely related to the gravitational rest frame, in that the Ricci identities provide a constraint that sets H_ab equal to the covariant curl of the shear tensor. Vanishing of the shear therefore corresponds to the vanishing of H_ab , which explains why the Poisson gauge (sometimes referred to as the 'zero shear' gauge) corresponds to the gravitational rest frame. Another is the 'constant density' foliation, in which the spatial gradient of the energy density vanishes, $D_a \mu = 0$. As with the comoving foliation, this has the benefit of being directly tied to the matter content of the space-time. However, it also shares the weakness that it will be highly problematic in the presence of non-linear structures.

In order to extract scalars from anisotropic quantities, we also need to choose a preferred space-like direction m. All 3-vector and tensor quantities can then be decomposed with respect to this direction as prescribed by equations (7) and (8), and the scalar parts of the decomposed quantities can be averaged according to equation (31). Clearly, in order to make sense, this vector must (in some sense) point in the same direction at every point in an averaging domain. Just as in the case of the time-like vector n this presents a choice; there is in general not going to be any single uniquely preferred space-like direction m, and we will need to select a suitable way to define this direction based on the situation at hand. This choice will be important for the outcome of our averaging process.

As in the case of choosing n, we may wish to apply geometric or observational considerations when selecting a preferred direction m. In space-times with symmetry, it may be that a preferred direction is selected by the existence of some Killing vectors. For example, in LRS geometries m can be unambiguously taken to correspond to the spatial axis around which the rotational symmetry exists. This can be identified by searching for the non-degenerate eigenvector of any non-vanishing 3-tensor, such as σ_ab or E_ab [47, 76]. Alternatively, in algebraically special geometries it may be possible to pick out a unique direction from the projections of the canonical null tetrad into our chosen foliation. These will be specified by the properties of the Weyl tensor.

An alternative method for selecting a preferred spatial direction would be to use observational methods at each point within an averaging domain. This method would be better adapted to situations in which anisotropy in a particular observable is being considered, or where the space-time has no explicit symmetries or special algebraic properties (as is the case for the real Universe). For example, one might choose to take the axis of CMB parity asymmetry [77, 78], the direction of greatest asymmetry in the galaxy distribution [14, 15], the quasar or Type Ia supernovae dipole [16, 17], or the direction of greatest variation of a cosmological parameter or coupling constant such as H₀ or the fine structure constant α [9, 79, 80]. As long as these directions line up at different points in space, as one would expect in an anisotropic Universe, then they should provide a well-motivated choice of preferred direction.

As well as considering the choice of foliation and space-like direction, it seems worthwhile to also consider what the model that emerges from averaging the scalars in our $1+1+2$-decomposition will represent. These will be precisely the set of scalars represented in equation (14). It is an extended set compared to the scalars that result from a $1+3$-decomposition, and include (among other things) the scalar parts of 3-vectors such as the momentum density Q, and the scalar part of 3-tensors such as the shear Σ.

The averaged scalars that result from the $1+3$-decomposition, as provided in Buchert's ground-breaking approach [37], are relatively straightforward to interpret. They result in a set of equations that govern the averages of the expansion scalar, the energy density of dust, and the curvature scalar of the 3-spaces:

$$\begin{align*} \{\langle \theta \rangle, \langle \mu \rangle, \langle {}^{(3)}\!R \rangle \} \, , \end{align*}$$

with everything else being collected together into 'back-reaction' terms. Defining an effective scalar factor $a_{\mathcal{D}}$ and an effective curvature parameter $k_{\mathcal{D}}$ of the averaging domain $\mathcal{D}$ by

$$\begin{align*} \langle \theta \rangle = \frac{\dot{a}_{\mathcal{D}}}{a_{\mathcal{D}}} \qquad \textrm{and } \qquad k_{\mathcal{D}} = \frac{\langle {}^{(3)}\!R \rangle}{6 a_{\mathcal{D}}^2} \, , \end{align*}$$

the equations that result bear a striking similarity to the dust-dominated Friedmann equations that govern the behaviour of homogeneous and isotropic cosmologies [37]. This is no accident, as the set of averaged scalars listed above are exactly those required to fully characterise an homogeneous and isotropic Universe. All other quantities will be 3-vectors or 3-tensors, but these must vanish in geometries that are isotropic around every point in space. So, after discarding such quantities in the scalar averaging approach, we are left with exactly the set that is required to describe an FLRW model. One could say that the geometry has been averaged to a Friedmann cosmology.

If we now consider the scalars that result from the $1+1+2$-decomposition, as given in equation (14), then we have a more complicated situation. It can, however, be understood in terms of a well studied class of cosmological models: the LRS models pioneered by Ellis [81, 82]. A space-time is LRS if it admits everywhere a one-dimensional continuous isotropy group, which essentially corresponds to rotations about a local axis of symmetry m. The existence of such a symmetry means that all 3-vectors are proportional to m, and all projected symmetric and trace-free 3-tensors are proportional to $m_a m_b - \frac{1}{2}M_{ab}$, where M_ab is defined in equation (5). In the language of the $1+1+2$-decompostions described in section 2.1, this means that all 2-vectors and 2-tensors must vanish. The entire dynamics of an LRS space-time is therefore described purely by the set of scalars given in equation (14). For this reason, we expect to recover a set of equations for our averaged scalars that can be written in a form similar to the equations that govern LRS cosmologies, with additional back-reaction terms due to the averaging of small-scale structure. We will therefore say that we are averaging to an LRS Bianchi cosmology.

The Bianchi classification contains all spatially homogeneous cosmological models (except Kantowski–Sachs) and divides them into classes labelled I to IX. The classes that admit LRS symmetry are I, II, III, V, VII and IX. Of these, the tilted LRS cosmologies can be found in classes V and VII${}_h$, and the fully isotropic cosmologies can be found in classes I, V, VII and IX. Due to their homogeneity, the only non-zero derivative of any of the scalars S within these models will be a time derivatice. All the governing equations will then necessarily have the form

$$\begin{equation} a_i \, \partial_t \left(N^{p_i} S_i\right) + N^{p_i + 1} \sum_j b_{ij} S_j + N^{p_i + 1}\sum_{j,k} c_{ijk} S_j S_k = 0\,, \end{equation} \tag{ 34 }$$

where S_i labels the allowed scalars from equation (14), b_ij and c_ijk are constants, and $a_i = 0$ or 1 depending on whether the equation is a constraint or an evolution equation. The exponent p_i of the lapse function N depends on the scalar being considered: for the kinematic scalars $\left\lbrace \Theta,\Sigma,\mathcal{A},\phi,\xi\right\rbrace$ we have $p_i = 1$, and for the matter and Weyl curvature scalars, $\left\lbrace \mu, p, Q, \Pi, \mathcal{E}, \mathcal{H}\right\rbrace$ we have $p_i = 2$. These equations are a significant simplification of the full set of equations of motion, and will be the form we expect for the equations that govern the average of the scalars from our $1+1+2$-decompostion. The only difference will be that we will find back-reaction terms on the right-hand side of each equation, instead of zero.

The Ricci identities taken with respect the time-like vector n can be written as

$$\begin{equation*} S_{abc} := 2\nabla_{[a}\nabla_{b]}n_c - R_{abcd}n^d = 0 \, . \end{equation*}$$

Taking the trace of S_abc over its first and third indices, and contracting with n, gives the Raychaudhuri equation ($n^b S^a_{\ \ ba} = 0$):

$$\begin{align} \partial_t \left\langle N \Theta\right\rangle - \left\langle N\mathcal{A}\right\rangle\left(\left\langle N\mathcal{A}\right\rangle + \left\langle N\phi\right\rangle\right) + \frac{1}{3}\left\langle N\Theta\right\rangle^2 + \frac{3}{2}\left\langle N\Sigma\right\rangle^2 +\frac{1}{2}\left\langle N^2\left(\mu + 3 p\right) \right\rangle - \left\langle N^2 \Lambda \right\rangle= \mathcal{B}_1 \,\nonumber\\ \end{align} \tag{ 35 }$$

where $\partial_t S = N \dot{S}$ and $\dot{S} = n^a \nabla_a S$. This is an evolution equation for the expansion, written in terms of coordinate time, and is identical to the corresponding equation from LRS Bianchi cosmologies up to the back-reaction term $\mathcal{B}_1$, and the presence of the lapse function N, which cannot in general be set to unity. We find that the back-reaction term can be written as

$$\begin{align*} \mathcal{B}_1 & = \frac{2}{3}\mathrm{Var}\left( N\Theta \right) - \frac{3}{2}\mathrm{Var} \left(N\Sigma\right) + \mathrm{Var}\left( N\mathcal{A} \right) + \mathrm{Cov}\left(N\mathcal{A},N\phi\right) -2\left\langle N^2\Sigma_a\Sigma^a\right\rangle - \left\langle N^2\Sigma_{ab}\Sigma^{ab}\right\rangle \\[3pt] & \quad + \left\langle N^2 m^a D_a \mathcal{A}\right\rangle + \left\langle N^2 M^{ab} D_a \mathcal{A}_b \right\rangle + \left\langle N^2 \left(\mathcal{A}_a - a_a\right)\mathcal{A}^a\right\rangle + \left\langle \Theta \, \partial_t {N}\right\rangle , \end{align*}$$

which encodes all information about the influence of small-scale inhomogeneities on the acceleration of the expansion of space. This is the only equation that can be derived from S_abc that has a counterpart in the standard approach of averaging to isotropic cosmology. As this is an equation for the acceleration of expansion, a sufficiently large (and positive) $\mathcal{B}_1$ would lead to an accelerating Universe, at least within the spatial domain being considered.

By again taking the trace of S_abc over the first and last indices, but this time contracting the middle index with respect to the space-like vector m we can obtain a momentum constraint equation ( $f_a^{\ b} S^c_{\ \ bc} m^a = 0$):

$$\begin{align} \left\langle N^2 Q \right\rangle + \frac{3}{2}\left\langle N\phi\right\rangle\left\langle N\Sigma\right\rangle = \mathcal{B}_2 \, , \end{align} \tag{ 36 }$$

where the back-reaction term is given by

$$\begin{align*} \mathcal{B}_2 &= \frac{2}{3}\left\langle N^2 m^a D_a \Theta \right\rangle - \left\langle N^2 m^a D_a \Sigma \right\rangle - \frac{3}{2}\mathrm{Cov}\left(N\phi,N\Sigma\right)- \left\langle N^2 M^{ab} D_a \Sigma_b\right\rangle\nonumber\\ & \quad + 2\left\langle N^2 a_b \Sigma^b\right\rangle + \left\langle N^2\Sigma_{ab}\zeta^{ab}\right\rangle. \end{align*}$$

This term describes the direct contribution from inhomogeneity to the large-scale momentum density in the direction of m.

There are two remaining equations that can be derived from S_abc . The first of these is an evolution equation for the scalar part of the shear ($f_{\langle a}^{\ \ d} f_{b\rangle}^{\ \ e} n^c S_{dec} m^a m^b = 0$):

$$\begin{align} & \partial_t\left\langle N\Sigma\right\rangle + \frac{2}{3}\left\langle N\Theta\right\rangle\left\langle N\Sigma\right\rangle + \frac{1}{2}\left\langle N\Sigma\right\rangle^2+\left\langle N^2\mathcal{E}\right\rangle - \frac{1}{2}\left\langle N^2\Pi\right\rangle\nonumber\\ & \quad - \frac{1}{3}\left(2\left\langle N\mathcal{A}\right\rangle - \left\langle N\phi\right\rangle \right)\left\langle N\mathcal{A}\right\rangle = \mathcal{B}_3 \, , \end{align} \tag{ 37 }$$

where the back-reaction term is given by

$$\begin{align*} \mathcal{B}_3 & = \frac{1}{3}\mathrm{Cov}\left(N\Theta,N\Sigma\right) + \frac{2}{3}\mathrm{Var}\left( N\mathcal{A} \right) - \frac{1}{3}\mathrm{Cov}\left(N\phi, N\mathcal{A}\right) + \frac{2}{3}\left\langle N^2 m^a D_a \mathcal{A}\right\rangle - \frac{1}{2}\mathrm{Var}\left(N\Sigma \right) \\ & \quad \nonumber - \frac{1}{3}\left\langle N^2 M^{ab} D_a \mathcal{A}_b \right\rangle - \frac{1}{3}\left\langle N^2\Sigma_a\Sigma^a\right\rangle + \frac{1}{3}\left\langle N^2 \mathcal{A}_a\mathcal{A}^a\right\rangle + \frac{1}{3}\left\langle N^2\Sigma_{ab}\Sigma^{ab}\right\rangle + 2\left\langle N^2\alpha_a\Sigma^a\right\rangle \\ & \quad \nonumber - \frac{2}{3}\left\langle N^2 a_a\mathcal{A}^a\right\rangle + \left\langle \Sigma\, \partial_t {N} \right\rangle \, . \end{align*}$$

This drives the generation of anisotropy in the expansion of space. The last equation to be derived from S_abc is a constraint that is usually used to set either the acceleration of n (i.e. $\mathcal{A}$) or the twist of m (i.e. ξ) to zero. In this case, however, we find ($f_a^{\ d} f_b^{\ e} n^c S_{dec}\eta^{abf} m_f = 0$):

$$\begin{align} \left\langle N\mathcal{A}\right\rangle\left\langle N\xi\right\rangle = \mathcal{B}_{4}\,, \end{align} \tag{ 38 }$$

where the back-reaction term is given by $\mathcal{B}_{4} = -\frac{1}{2}\left\langle N^2 \epsilon_{ab}D^a \mathcal{A}^b\right\rangle - \mathrm{Cov}\left(N\mathcal{A},N\xi\right)$, In this case we see that any non-zero effect from the small-scale inhomogeneity will lead to a more complicated result, with $\left\langle N \mathcal {A}\right\rangle \neq 0 \neq \left\langle N\xi\right\rangle$.

We can also consider the Ricci identities defined with respect to the preferred space-like direction m:

$$\begin{align*} R_{abc} := 2\nabla_{[a}\nabla_{b]}m_c - R_{abcd}m^d =0 \, . \end{align*}$$

This time, the only equation that has a counterpart in isotropic cosmology is obtained by contracting over the last two indices with the screen space projection tensor, and contracting the first index with m. This gives the Hamiltonian constraint equation ($m^a M^{bc} R_{abc} = 0$):

$$\begin{align} & \frac{2}{9}\left\langle N\Theta\right\rangle^2-\frac{1}{2}\left\langle N\phi\right\rangle^2+\frac{1}{3}\left\langle N\Theta\right\rangle\left\langle N\Sigma\right\rangle -\left\langle N\Sigma\right\rangle^2 - \frac{2}{3}\left(\left\langle N^2\mu\right\rangle + \left\langle N^2 \Lambda \right\rangle\right)\nonumber\\ & \quad - \frac{1}{2}\left\langle N^2 \Pi\right\rangle-\left\langle N^2\mathcal{E}\right\rangle+2\left\langle N\xi\right\rangle^2 = \mathcal{B}_5 . \end{align} \tag{ 39 }$$

The back-reaction term in this case is given by the following expression:

$$\begin{align*} \mathcal{B}_5 &= \left\langle N^2 m^a D_a \phi \right\rangle + \frac{1}{2}\mathrm{Var}\left( N\phi \right) - \frac{2}{9}\mathrm{Var}\left( N\Theta \right) + \mathrm{Var}\left( N\Sigma \right)- \frac{1}{3}\mathrm{Cov}\left(N\Theta,N\Sigma\right) - 2\,\mathrm{Var}\left(N\xi \right) \\ \nonumber & \quad + \left\langle N^2\zeta_{ab}\zeta^{ab}\right\rangle + \left\langle N^2 \Sigma_a\Sigma^a\right\rangle + \left\langle N^2 a_b a^b\right\rangle - \left\langle N^2 M^{ab} D_a a_b\right\rangle. \end{align*}$$

This term acts in the same way as an additional effective energy density, in what would reduce to the Friedmann equation in isotropic cosmology. Comparing (39) to the trace of the Gauss embedding equation, $^{(3)}R = 2\mu - \frac{2}{3}\Theta^2 + 2\sigma^2 + 2\Lambda$, and using $\sigma^2 = \frac{3}{4}\Sigma^2 + \Sigma_a \Sigma^a + \frac{1}{2}\Sigma_{ab}\Sigma^{ab}$, one finds that the average (lapse-weighted) Ricci curvature of the space-like hypersurfaces is given by

$$\begin{align*} \left\langle N^2 \: ^{(3)} R \right\rangle & = -\frac{3}{2}\left\langle N\phi\right\rangle^2 - \frac{3}{2}\left\langle N\Sigma\right\rangle^2 + \left\langle N\Theta\right\rangle\left\langle N\Sigma\right\rangle - 3\left\langle N^2\mathcal{E}\right\rangle - \frac{3}{2}\left\langle N^2\Pi\right\rangle + 6\left\langle N\xi\right\rangle^2 \nonumber \\ & \quad - 3\left\langle N^2 m^a D_a \phi\right\rangle - \frac{3}{2}\,\mathrm{Var}\left( N\phi \right) - \frac{3}{2}\,\mathrm{Var}\left( N\Sigma \right) + \mathrm{Cov}\left(N\Theta,N\Sigma\right) + 6\,\mathrm{Var}\left( N\xi \right)\nonumber \\ & \quad - \left\langle N^2\Sigma_a\Sigma^a\right\rangle - 3\left\langle N^2 a_b a^b\right\rangle + 3\left\langle N^2 M^{ab}D_a a_b\right\rangle - 3\left\langle N^2\zeta_{ab}\zeta^{ab}\right\rangle + \left\langle N^2 \Sigma_{ab}\Sigma^{ab}\right\rangle. \end{align*}$$

The second equation to be gained from contracting R_abc with the screen space projector is to be found using the time-like vector n, and is an evolution for the expansion of m in the preferred space-like direction ($n^a M^{bc} R_{abc} = 0$):

$$\begin{align} \partial_t \left\langle N\phi\right\rangle - \left\langle N^2 Q\right\rangle -\left(\frac{2}{3}\left\langle N\Theta\right\rangle - \left\langle N \Sigma\right\rangle\right)\left(\left\langle N\mathcal{A}\right\rangle - \frac{1}{2}\left\langle N\phi\right\rangle\right) = \mathcal{B}_{6} \, , \end{align} \tag{ 40 }$$

where we find

$$\begin{align*} \mathcal{B}_{6} &= \frac{2}{3}\mathrm{Cov}\left(N\Theta,N\phi\right) + \frac{1}{2}\mathrm{Cov}\left(N\Sigma,N\phi\right) + \frac{2}{3}\mathrm{Cov}\left(N\Theta,N\mathcal{A}\right)-\mathrm{Cov}\left(N\Sigma,N\mathcal{A}\right) \\ \nonumber & \quad + \left\langle N^2 M^{ab} D_a \Sigma_b\right\rangle - \left\langle N^2 \zeta_{ab}\Sigma^{ab}\right\rangle + \left\langle N^2 \mathcal{A}^a\left(\alpha_a - \Sigma_a\right)\right\rangle + \left\langle N^2 a^a \left(\Sigma_a - \mathcal{A}_a\right)\right\rangle + \left\langle \phi \, \partial_t {N}\right\rangle. \end{align*}$$

This equation has no counterpart in isotropic cosmological modelling.

The first of the three further equations that can be derived from R_abc are a constraint on the magnetic part of the Weyl tensor ($\epsilon^{ab} n^c R_{abc} = 0$):

$$\begin{align} \left\langle N^2 \mathcal{H}\right\rangle - 3\left\langle N\xi\right\rangle\left\langle N\Sigma\right\rangle = \mathcal{B}_{7} \, , \end{align} \tag{ 41 }$$

where we find $\mathcal{B}_{7} = 3\,\mathrm{Cov}\left(N\xi,N\Sigma\right) + \left\langle N^2 \epsilon_{ab}\left(D^a\Sigma^b - \zeta^{ac}\Sigma^b_{\:c}\right)\right\rangle$. For our averaged cosmology to be 'silent' we therefore require $\mathcal{B}_{7}$ to vanish along with either $\left\langle N\xi\right\rangle$ or $\left\langle N\Sigma\right\rangle$, or to conspire to cancel $3\left\langle N\xi\right\rangle\left\langle N\Sigma\right\rangle$. The second equation is an evolution equation for the twist of m ($n^a \epsilon^{bc}R_{abc} = 0$):

$$\begin{align} \partial_t\left\langle N\xi\right\rangle + \frac{1}{2}\left(\frac{2}{3}\left\langle N\Theta\right\rangle - \left\langle N\Sigma\right\rangle\right)\left\langle N\xi\right\rangle - \frac{1}{2}\left\langle N^2\mathcal{H}\right\rangle = \mathcal{B}_{8} \,, \end{align} \tag{ 42 }$$

with back-reaction term

$$\begin{align*} \mathcal{B}_{8} &= \frac{2}{3}\mathrm{Cov}\left( N\Theta, N\xi\right) + \frac{1}{2}\mathrm{Cov}\left(N\Sigma,N\xi\right) + \frac{1}{2}\left\langle N^2 \epsilon_{ab}\left(a^a + \mathcal{A}^a\right)\left(\alpha^b+\Sigma^b\right)\right\rangle \\ & \quad + \frac{1}{2}\left\langle N^2 \epsilon_{ab}D^a\alpha^b\right\rangle + \frac{1}{2}\left\langle N^2 \epsilon_{ac}\zeta_b^{\ c}\Sigma^{ab}\right\rangle + \left\langle \xi\,\partial_t {N}\right\rangle\,. \end{align*}$$

The final equation involving R_abc is a contraint on the product of the expansion of m in the preferred space-like direction and its twist ($m^a \epsilon^{bc} R_{abc} = 0$):

$$\begin{align} \left\langle N\phi\right\rangle\left\langle N\xi\right\rangle = \mathcal{B}_{9} \, . \end{align} \tag{ 43 }$$

The back-reaction term in this equation is

$$\begin{align*} \mathcal{B}_{9} = -\left\langle N^2 m^a D_a \xi \right\rangle - \mathrm{Cov}\left(N\phi,N\xi\right) + \frac{1}{2}\left\langle N^2 \epsilon_{ab}\left(D^a a^b + \Sigma^a a^b\right)\right\rangle\,. \end{align*}$$

As in equation (38), the back-reaction term in this last equation prevents one obtaining the usual result that either $\left\langle N\phi\right\rangle$ or $\left\langle N\xi\right\rangle$ vanishes, with $\left\langle N\phi\right\rangle \neq 0 \neq \left\langle N\xi\right\rangle$ being required if $\mathcal{B}_{9}\neq 0$.

Our final set of equations are derived from the Bianchi identities:

$$\begin{align} W_{abcde} := \nabla_{[a}R_{bc]de} = 0 \, . \end{align} \tag{ 44 }$$

The only equation that could be derived from these identities in the case of isotropic cosmologies would be the energy conservation equation, which in our case takes the form ($n^a W_{abc}^{\ \ \ bc} = 0$):

$$\begin{align} \partial_t \left\langle N^2 \mu\right\rangle + \left\langle N\Theta\right\rangle\left(\left\langle N^2\mu\right\rangle+\left\langle N^2 p \right\rangle \right) + \frac{3}{2}\left\langle N\Sigma\right\rangle\left\langle N^2 \Pi\right\rangle +\left(\left\langle N\phi\right\rangle + 2\left\langle N\mathcal{A}\right\rangle\right)\left\langle N^2 Q \right\rangle = \mathcal{B}_{10} \, , \end{align} \tag{ 45 }$$

where the back-reaction term is given by

$$\begin{align*} \mathcal{B}_{10} &= -\mathrm{Cov}\left(N^2 p,N\Theta\right) - \left\langle N^3 M^{ab} D_a Q_b \right\rangle - \left\langle N^3 m^a D_a Q \right\rangle - \mathrm{Cov}\left(N\left(\phi+2\mathcal{A}\right),N^2Q\right)\\ \nonumber & \quad - \frac{3}{2}\mathrm{Cov}\left(N\Sigma,N^2\Pi\right) - 2\left\langle N^3 \mathcal{A}_a Q^a\right\rangle - 2\left\langle N^3\Sigma_a\Pi^a\right\rangle - \left\langle N^3\Sigma_{ab}\Pi^{ab}\right\rangle \\ \nonumber & \quad + \left\langle N^3 a_a Q^a\right\rangle + 2\left\langle N \mu\,\partial_t {N}\right\rangle\,. \end{align*}$$

The back-reaction term here provides an additional source that drives the rate of change of the average energy density, $\left\langle N^2 \mu\right\rangle$. All other equations correspond to anisotropic cosmologies only, and vanish in the isotropic limit. It can be seen that if one chooses the foliation such that n is the 4-velocity of pressureless dust, then $\mathcal{B}_{10} = 0$. For a generic foliation, however, this need not be the case, as an observer that is not comoving with the dust will typically measure non-zero pressure, momentum density and anisotropic stress.

The first of the remaining equations gives conservation of momentum density ($m^a W_{abc}^{\ \ \ bc} = 0$):

$$\begin{align} & \partial_t\left\langle N^2 Q\right\rangle + \left(\frac{4}{3}\left\langle N\Theta\right\rangle + \left\langle N\Sigma\right\rangle\right) \left\langle N^2 Q\right\rangle - \left\langle N\mathcal{A}\right\rangle\left\langle N^2\left(\mu + p + \Pi\right)\right\rangle\nonumber\\ & \quad + \frac{3}{2}\left\langle N^2 \Pi\right\rangle\left\langle N\phi\right\rangle = \mathcal{B}_{11}. \end{align} \tag{ 46 }$$

This equation contains the back-reaction term

$$\begin{align*} \mathcal{B}_{11} &= -\frac{1}{3}\mathrm{Cov}\left(N\Theta, N^2 Q\right) - \left\langle N^3 m^a D_a \left(p + \Pi\right) \right\rangle - \left\langle N^3 M^{ab} D_a \Pi_b\right\rangle - \mathrm{Cov}\left(\frac{3}{2}N\phi+N\mathcal{A},N^2\Pi\right) \\ \nonumber & \quad - \mathrm{Cov}\left(N\Sigma, N^2 Q\right) - \mathrm{Cov}\left(N^2\left(\mu+p\right), N \mathcal{A}\right) - \left\langle \mathcal N^3 {A}_a\Pi^a\right\rangle + \left\langle N^3 \zeta_{ab}\Pi^{ab}\right\rangle \\ \nonumber & \quad + \left\langle N^3 \left(\alpha_a-\Sigma_a\right)Q^a\right\rangle + 2\left\langle N^3 a_a \Pi^a\right\rangle + 2\left\langle N Q \,\partial_t {N}\right\rangle\,. \end{align*}$$

This term drives the evolution of the momentum density, and if non-zero can therefore source a bulk flow in the averaged cosmology.

There are four remaining equations to be derived from W_abcde . Two of these are evolution equations for the scalar parts of the electric and magnetic Weyl tensor, and two are constraints for these quantities. The evolution equation for the electric part of the Weyl tensor is given by ( $f\,^{de} n^c W_{cdabe} m^a m^b = 0$):

$$\begin{align} & \partial_t\left\langle N^2 \mathcal{E}\right\rangle + \frac{1}{2}\partial_t\left\langle N^2 \Pi\right\rangle + \left(\left\langle N\Theta\right\rangle - \frac{3}{2}\left\langle N\Sigma\right\rangle \right)\left\langle N^2\mathcal{E}\right\rangle + \frac{1}{2}\left(\frac{1}{3}\left\langle N\Theta\right\rangle + \frac{1}{2}\left\langle N\Sigma\right\rangle \right)\left\langle N^2\Pi\right\rangle\nonumber\\ & \!\!\quad- \frac{1}{3}\left(\frac{1}{2}\left\langle N\phi\right\rangle-2\left\langle N\mathcal{A}\right\rangle\right)\left\langle N^2 Q \right\rangle + \frac{1}{2}\left(\left\langle N^2 \mu\right\rangle+\left\langle N^2 p \right\rangle \right) \left\langle N\Sigma\right\rangle - 3\left\langle N\xi \right\rangle\left\langle N^2 \mathcal{H}\right\rangle = \mathcal{B}_{12} , \end{align} \tag{ 47 }$$

with back-reaction

$$\begin{align*} \mathcal{B}_{12} &= -\frac{1}{3}\left\langle N^3 m^a D_a Q \right\rangle + \frac{1}{3}\mathrm{Cov}\left(N\Theta,N^2\Pi\right) + \left\langle N^3\epsilon_{ab}D^a \mathcal{H}^b\right\rangle + \frac{1}{6}\left\langle N^3 M^{ab} D_a Q_b\right\rangle \\ \nonumber & \quad + \frac{3}{2}\mathrm{Cov}\left(N\Sigma,N^2\mathcal{E}\right) + \left\langle N^3\epsilon_{ab}\mathcal{H}^{bc}\zeta^a_{\ c}\right\rangle -\frac{1}{2}\mathrm{Cov}\left(N\Sigma,N^2\left(\mu+p+\frac{1}{2}\Pi\right)\right) \\ \nonumber & \quad+ \frac{1}{3}\mathrm{Cov}\left(N^2 Q,\frac{1}{2}N\phi-2N\mathcal{A}\right) - \frac{1}{6}\left\langle N^3 \Sigma_a\Pi^a\right\rangle + \frac{1}{3}\left\langle N^3 \left(\mathcal{A}_a + a_a\right)Q^a\right\rangle \\ \nonumber & \quad + 2\left\langle N^3 \epsilon_{ab}\mathcal{A}^a\mathcal{H}^b\right\rangle - \left\langle N^3 \Sigma_{ab}\left(\mathcal{E}^{ab}+\frac{1}{2}\Pi^{ab}\right)\right\rangle + 3\,\mathrm{Cov}\left(N\xi, N^2\mathcal{H}\right) + \left\langle N^3\alpha_a\Pi^a\right\rangle \\ \nonumber & \quad + \left\langle N^3\left(2\alpha_a+\Sigma_a\right)\mathcal{E}^a\right\rangle + 2\left\langle N \mathcal{E}\,\partial_t {N}\right\rangle\,. \end{align*}$$

The evolution equation for the magnetic part of the Weyl tensor is ($\eta_a^{\ cd} n\,^{e} n\,^{f} W_{ecdfb}m^a m^b = 0$):

$$\begin{align} \partial_t\left\langle N^2\mathcal{H}\right\rangle + \left(\left\langle N\Theta\right\rangle - \frac{3}{2}\left\langle N\Sigma\right\rangle\right)\left\langle N^2\mathcal{H}\right\rangle + 3\left\langle N\xi\right\rangle\left(\left\langle N^2\mathcal{E}\right\rangle-\frac{1}{2}\left\langle N^2\Pi\right\rangle\right) = \mathcal{B}_{13} \, , \end{align} \tag{ 48 }$$

which has back-reaction term

$$\begin{align*} \mathcal{B}_{13} &= -\left\langle N^3\epsilon_{ab}D^a\left(\mathcal{E}^b-\frac{1}{2}\Pi^b\right)\right\rangle + \frac{3}{2}\mathrm{Cov}\left(N\Sigma,N^2\mathcal{H}\right) + 2\left\langle N^3\epsilon_{ab}\mathcal{E}^a \mathcal{A}^b\right\rangle\nonumber \\ & \quad + \left\langle N^3\left(2\alpha_a + \Sigma_a\right)\mathcal{H}^a\right\rangle + \frac{1}{2}\left\langle N^3\epsilon_{ab}\Sigma^a Q^b\right\rangle - \left\langle N^3\Sigma_{ab}\mathcal{H}^{ab}\right\rangle + \frac{1}{2}\left\langle N^3\epsilon_{ab}\mathcal{E}^{ac}\zeta^b_{\ c}\right\rangle \\ \nonumber & \quad - 3\,\mathrm{Cov}\left(N\xi,N^2\mathcal{E}-\frac{1}{2}N^2\Pi\right) + \left\langle N \mathcal{H} \, \partial_t {N}\right\rangle\,. \end{align*}$$

The remaining two equations are a constraint equation for the electric Weyl tensor ($f_a^{\ b} n^c n^d W_{db \ \ ec}^{\ \ e} m^a = 0$):

$$\begin{align} \frac{3}{2}\left(\left\langle N^2\mathcal{E}\right\rangle + \frac{1}{2}\left\langle N^2\Pi\right\rangle\right)\left\langle N\phi\right\rangle + \left(\frac{1}{3}\left\langle N\Theta\right\rangle-\frac{1}{2}\left\langle N\Sigma\right\rangle\right)\left\langle N^2 Q\right\rangle = \mathcal{B}_{14}\,, \end{align} \tag{ 49 }$$

which we find to have a back-reaction term given by

$$\begin{align*} \mathcal{B}_{14} &= -\left\langle N^3 m^a D_a \mathcal{E}\right\rangle + \frac{1}{3}\left\langle N^3 m^a D_a \mu \right\rangle - \frac{1}{2}\left\langle N^3 m^a D_a \Pi\right\rangle - \left\langle N^3 M^{ab} D_a \left(\mathcal{E}_b+\frac{1}{2}\Pi_b\right)\right\rangle \\ \nonumber & \quad - \frac{3}{2}\mathrm{Cov}\left(N\phi,N^2\mathcal{E}+\frac{1}{2}N^2\Pi\right)+\frac{1}{2}\mathrm{Cov}\left(N^2 Q,N\Sigma-\frac{2}{3}N\Theta\right)+ \frac{1}{2}\left\langle N^3\Sigma_a Q^a\right\rangle \\ \nonumber & \quad + \left\langle N^3 \epsilon_{ab}\Sigma^{ac}\mathcal{H}_c^{\ b}\right\rangle + \left\langle N^3\left(\mathcal{E}_{ab}+\frac{1}{2}\Pi_{ab}\right)\zeta^{ab}\right\rangle + \left\langle N^3\left(2\mathcal{E}_a+\Pi_a\right)a^a\right\rangle, \end{align*}$$

and a constraint for the magnetic part of the Weyl tensor ($\eta^{bc}_{\ \ a} n^d W_{db \ \ ec}^{\ \ e} m^a = 0$):

$$\begin{align} \frac{3}{2}\left\langle N\phi\right\rangle\left\langle N^2\mathcal{H}\right\rangle + \left\langle N^2 Q \right\rangle \left\langle N\xi\right\rangle = \mathcal{B}_{15}. \end{align} \tag{ 50 }$$

The final back-reaction term for this equation is given by

$$\begin{align*} \mathcal{B}_{15} &= \left\langle N^3 m^a D_a \mathcal{H} \right\rangle - \frac{3}{2}\mathrm{Cov}\left(N\phi,N^2\mathcal{H}\right) - \mathrm{Cov}\left(N^2Q,N\xi\right) - \frac{1}{2}\left\langle N^3\epsilon_{ab}D^a Q^b\right\rangle \\ \nonumber & \quad - \left\langle N^3 M^{ab} D_a \mathcal{H}_b\right\rangle + 2\left\langle N^3 a_b \mathcal{H}^b\right\rangle + \left\langle N^3 \zeta_{ab}\mathcal{H}^{ab}\right\rangle - \left\langle N^3\epsilon_{ab}\Sigma^a_{\ c}\left(\mathcal{E}^{bc}+\frac{1}{2}\Pi^{bc}\right)\right\rangle. \end{align*}$$

Equations (35)–(50) provide a complete set that can be used to describe the large-scale properties of an anisotropic cosmological model after averaging.

These equations are a considerable complication on those that govern exact LRS cosmologies, due to the presense of the back-reaction terms, $\mathcal{B}_i$. If all of these terms are sufficiently small when calculated for a domain $\mathcal{D}$, then the expansion of $\mathcal{D}$ should be expected to be well approximated by an exact LRS Bianchi model. If this is not the case, then we will be in a situation where back-reaction of small-scale structures on the large-scale properties of the cosmology is no longer negligible. The reader may note that the equations above are all written in terms of the coordinate time t, which (by definition) takes the same value at every spatial averaging domain $\mathcal{D}$, as required. In the next section we will consider a family of exact cosmological models, to help develop our understanding of these terms by explicitly calculating them for an example geometry.

Let us first consider the foliation composed of leaves that coincide with the surfaces of transitivity of the Killing vectors from equation (52). These are spatially homogeneous 3-dimensional spaces with a time-like normal n that must satisfy $n_a X_3^a = 0$. This immediately implies $v = -C/X$, and results in a set of spaces of constant η. The set of non-vanishing scalars in this foliation is larger than in the foliation orthogonal to the matter flow, as it contains a bulk flow that gives rise to non-zero pressure, heat flux and anisotropic stress. The full set of scalars is $\left\lbrace \Theta, \Sigma, \mathcal{E}, \phi, \mu, p, Q,\Pi \right\rbrace$, which can be immediately obtained from equations (58)–(62). For each of these scalars (S) we also know that their projected covariant derivative ($D_a S$) must vanish, due to homogeneity. This immediately implies that this foliation is both constant-density and constant-mean-curvature.

The line-element in equation (51) can be recast into coordinates adapted to this foliation by making the transformation [83]

$$\begin{equation} t = Cx + \int dT \, \beta(T) \qquad \textrm{and} \qquad r = - x + \int dT \frac{\alpha(T)}{\beta(T)} \,, \end{equation} \tag{ 64 }$$

where $\alpha = {C}/({X^2-C^2})$ and $\beta^2 = 1 - {C^2}/({X^2-C^2})$. This results in

$$\begin{equation} ds^2 = -dT^2 + A(T)dx^2 + B(T)e^{2x}\left(dy^2+dz^2\right), \end{equation} \tag{ 65 }$$

where $X^2 = C^2 + A$ and $Y^2 = B e^{4x} \exp{\left\{2\int dT \, {\alpha}/{\beta}\right\}}$. All quantities that were functions of t + Cr are now purely functions of T, and moreover we have that $\mathbf{n} = \partial_T$ and $\mathbf{X}_3 = -\partial_x + y\partial_y + z\partial_z$. This makes clear the orthogonality of n with the homogeneous spaces spanned by the Killing vectors from equation (52). The lapse function is equal to unity in these coordinates.

The behaviour of the expansion rates of space, and the expansion-normalised variables $\left\lbrace \Sigma^{(n)}, \phi^{(n)}, \mathcal{E}, \mu^{(n)} \right\rbrace$, are shown in figures 2 and 3. Figure 2 shows the expansion rates parallel and orthogonal to the direction picked out by the bulk flow, and displays significant anisotropy at early times (small η). As the bulk flow decays at late times (large η), the difference in expansion rates in different directions also decays, demonstrating that space isotropises in this foliation. Figure 3 shows the other kinematic scalars as a function of η, normalised by the appropriate power of the expansion rate Θ to make the resultant quantity dimensionless. One sees shear Σ and and electric Weyl curvature $\mathcal{E}$ both increase at early times, consistent with growing anisotropy in this limit. In both plots η is made dimensionless by normalising it relative to $\eta_\mathrm{eq}$, its value when $Y = -D/3k^2$ (such that the two terms on the right-hand side of equation (54) are equal). In figure 2 the expansion rates are made dimensionless by normalising them with respect to the mean expansion rate $\bar{H}^{(n)} = \frac{1}{3} \Theta^{(n)}$ at the equality time $\eta_\mathrm{eq}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In this foliation it can be seen that all of back-reaction scalars $\mathcal{B}_i$ from section 3.3 must be equal to zero. This happens because all variances and covariances must vanish (as $\int_{\mathcal{D}} d^3 x\,\sqrt{f}\, S(T) = S(T) V_{\mathcal{D}}$), and because all terms containing projected covariant derivatives (such as $m^a D_a S$) must also vanish. The cosmology one observes in this foliation can therefore be entirely described by the locally defined scalar quanitites, without any need for averaging. This is as exactly as expected for a space-time that admits a homogeneous foliation. Of course, in most geometries such a simplifying set of symmetries will not exist, and in more realistic cosmologies one would need to perform explicit averages in order to gain a set of quantities that could be used to describe the large-scale properties of space. Let us now consider an inhomogeneous foliation of the Farnsworth solutions, to show how this would work in this simple example space-time.

Let us now foliate the exact same space-time into hypersurfaces orthogonal to the matter flow u, with induced metric $h_{ab} = g_{ab} + u_au_b$. This corresponds to v = 0 in equations (58)–(62), and to surfaces of constant coordinate time t. The lapse is again equal to unity, so covariant and coordinate time derivatives are equivalent. The only non-zero quantities in this case are the scalars $\left\lbrace \Theta, \Sigma, \mathcal{E}, \mu, \phi \right\rbrace$, which in terms of the parameters from the solution in equation (55) are

$$\begin{align} \Theta &= \frac{2\left[3W\sinh{\eta}\left(\cosh{\eta}-1\right)+2C\left(2\cosh{\eta}+1\right)\right]}{W\left(\cosh{\eta}-1\right)\left[W\left(\cosh{\eta}-1\right)^2-2C\sinh{\eta}\right]} \end{align} \tag{ 66 }$$

$$\begin{align} \Sigma &= \frac{8C\left(\sinh^2{\eta}+\cosh{\eta}-1\right)}{3W\left(\cosh{\eta}-1\right)^2\left[W\left(\cosh{\eta}-1\right)^2-2C\sinh{\eta}\right]} \end{align} \tag{ 67 }$$

$$\begin{align} \phi &= \frac{4}{W\left(\cosh{\eta}-1\right)}, \end{align} \tag{ 68 }$$

and the dust density

$$\begin{equation} \mu = \rho = \frac{24}{W\left(\cosh{\eta}-1\right)\left[W\left(\cosh{\eta}-1\right)^2-2C\sinh{\eta}\right]}. \end{equation} \tag{ 69 }$$

The electric Weyl curvature is given by equation (63). The normalised expansion rates as a function of η are shown in figure 4, and the other expansion-normalised scalars are shown in figure 5. These plots, displaying quantities calculated in the foliation generated by the matter flow u, may be directly compared with figures 2 and 3.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In order to average these quantities we need an averaging domain $\mathcal{D}$ on each of the surfaces orthogonal to u. If we take this domain to have coordinate extents $\left\lbrace \Delta y, \Delta z, \Delta r\right\rbrace$, then the spatial volume is given by

$$\begin{align} V_{\mathcal{D}}(t) = \int_{\mathcal{D}} \mathrm{d}y\mathrm{d}z\mathrm{d}r \sqrt{h(\eta)} = \frac{\Delta y \Delta z W^3 k^2}{16 C} \,e^{2t/C} \, \int_{\eta(r_\mathrm{min},t)}^{\eta(r_\mathrm{max}, t)} d\eta \, J(\eta), \end{align} \tag{ 70 }$$

where $J(\eta) = \left(\cosh{\eta}-1\right)^2 \left[W\left(\cosh{\eta}-1\right)^2 - 2C\sinh{\eta}\right]\exp{\left\{-{W}\left(\sinh{\eta} - \eta\right)/C\right\}}$. We can then compute all the desired scalar averages as ²

$$\begin{equation} \left\langle S(t)\right\rangle = \frac{\int_{\eta(r_\mathrm{min},t)}^{\eta(r_\mathrm{max}, t)} d\eta \, J(\eta) S(\eta)}{\int_{\eta(r_\mathrm{min},t)}^{\eta(r_\mathrm{max}, t)} d\eta \, J(\eta)}\, , \end{equation} \tag{ 71 }$$

where $\Delta r = r_\mathrm{max} - r_\mathrm{min}$. In choosing the integration domain (i.e. $r_\textrm{min/max}$), it is important to avoid the big bang singularity that is located at $r_\textrm{BB}(t)$ such that $t + C r_\textrm{BB}(t) = 0$. For any C > 0 this is satisfied by $r_\mathrm{min} \gt0$. The value of $r_\mathrm{max}$ needs only to be sufficiently large that $\eta_\mathrm{min}$ and $\eta_\mathrm{max}$ are relatively far apart, but once this is ensured to be true the backreaction effects are only weakly dependent on its value. If we choose $r_\mathrm{min} = 1$, $r_\mathrm{max} = 20$ and C = 2 as an example, then one obtains the integration domain shown in figure 6. At early times the distance from $r_\mathrm{min}$ to $r_\mathrm{max}$ corresponds to a large separation in η, while at late times the integration domain shrinks (as the space isotropises).

Zoom In Zoom Out Reset image size

The evolution with t of our averaged scalars, normalised by the averaged expansion $\left\langle\Theta\right\rangle$, is shown in figure 7. The inhomogeneity in these spaces has been averaged out to give a set of purely time-dependent scalar functions. These averages provide a cosmological description of the large-scale properties of the spaces orthogonal to the matter flow u, including anisotropy. The magnitude of the back-reaction scalars that produce them are displayed in figure 8, where we focus on $\left\lbrace\mathcal{B}_1,\mathcal{B}_3,\mathcal{B}_6,\mathcal{B}_{12}\right\rbrace$ as they appear in evolution rather than constraint equations. These quantities are shown as ratios of the time derivative of the relevant average, as per the equations in section 4, in order to show the fractional contribution they make to their evolution. In each case, the averaging domain used is the one displayed in figure 6, and we have again used $\left\{C,W\right\} = \{2,125\}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

It can be seen from figure 8 that the back-reaction effect becomes large at early times, and dissipates to zero at late times. This is as expected in the Farnsworth solution, as anisotropy (and hence inhomogeneity in the matter-rest-space foliation) is entirely due to the bulk flow $v = -C/X$, which becomes small as the scale factor X grows and the space isotropises. This example demonstrates that the large-scale averaged properties of a space depend on the chosen foliation, and that one requires in addition a choice of averaging domain (and in general a spatial direction m). The kinematic and matter quantities that result can display different behaviours when calculated on different foliations, and so can the magnitude of back-reaction scalars. This clearly demonstrates the need for these choices to be made carefully, and for the interpretation of scalar averaging results in cosmology to be understood in terms of these choices.

• Preferred (hypersurface-forming) time-like vector: n.
• Preferred space-like vector: m.
• Projection tensor into space-like hypersurfaces orthogonal to n: $f_{ab} = g_{ab} + n_a n_b$.
• Projected covariant derivative on space-like hypersurfaces: $D_a T_{b c} = f_a^{\phantom{a} d} f_b^{\phantom{b} e} f_c^{\phantom{c} f} \nabla_d T_{ef}$.
• Proper time derivative along n: $\dot{\phantom{a}} = n^a \nabla_a$.
• Projection tensor into screen spaces orthogonal to both n and m: $M_{ab} = f_{ab} - m_a m_b$.
• Alternating tensor on space-like hypersurfaces: $\eta_{abc} = n^d \eta_{dabc}$.
• Alternating tensor on screen space: $\epsilon_{ab} = \eta_{abc}m^c$.

• Expansion of space-like hypersurfaces: $\Theta = D_a n^a$.
• Volume-preserving shear: $\sigma_{ab} = D_{\langle a}n_{b\rangle} = \left(f_{(a}^{\ c} \,f_{b)}^{\ d} - \frac{1}{3}f_{ab}f\,^{cd}\right)\nabla_c n_d$.
• Scalar part of shear: $\Sigma = \sigma_{ab} m^a m^b$.
• Vector part of shear: $\Sigma_a = M_{ab} m_c \sigma^{bc}$.
• Tensor part of shear: $\Sigma_{ab} = \left(M^{\ c}_{(a}M^{\ d}_{b)} - \frac{1}{2}M_{ab}M^{cd}\right)\sigma_{cd}$.
• Vorticity: $\omega_a = \frac{1}{2}\eta_{abc} D^b n^c$ (choosen to be zero throughout, for n).
• Scalar part of vorticity: $\Omega = \omega_a m^a = 0$.
• Vector part of vorticity: $\Omega_a = M_{ab} \omega^b = 0$.
• Acceleration: $\dot{n}_a = n^b \nabla_b n_a$.
• Scalar part of acceleration: $\mathcal{A} = \dot{n}_a m^a$.
• Vector part of acceleration: $\mathcal{A}_a = M_{ab}\dot{n}^b$.

• Expansion of screen space along m: $\phi = M^{ab} D_a m_b$.
• Area-preserving shear of screen space along m: $\zeta_{ab} = M^{\ c}_{(a}M^{\ d}_{b)} D_c m_d - \frac{1}{2}\phi M_{ab}$.
• Twist of screen space: $\xi = \frac{1}{2}\epsilon_{ab} M^{ac} M^{bd} D_c m_d$.
• Non-geodesy of m: $a_b = m^c D_c m_b$.
• Scalar part of $\dot{m}_a$: $\mathcal{A} = -n_a\dot{m}^a$, since $n_a m^a = 0$.
• Vector part of $\dot{m}_a$: $\alpha_a = M_a^{\ b} \dot{m}_b$.

• New time-like vector (not necessarily hypersurface-forming): u.
• Projection tensor orthogonal to u: $h_{ab} = g_{ab} + u_a u_b$.
• Relation between time-like vectors: $n_a = \gamma\left(u_a + w_a\right) \ \Leftrightarrow \ u_a = \gamma\left(n_a - v_a\right)$, where $\gamma = \left(1 - v_a v^a\right)^{-1/2} = \left(1 - w_a w^a \right)^{-1/2}$.
• Boost vectors relating n and u: $w_a = \gamma^{-1} h_a^{\ b} v_b$, and $v_a = \gamma^{-1} f_a^{\ b} w_b$.
• Relation between space-like vectors, projected into surfaces orthogonal to the time-like vectors: $\tilde{m}_a = k h_a^{\ b} m_b$, and $m_a = \tilde{k}f_a^{\ b} \tilde{m}_b$, where $k = \left[1+\left(\gamma m_c v^c\right)^2\right]^{-1/2}$ and $\tilde{k} = \left[1+\left(\gamma \tilde{m}_c w^c\right)^2\right]^{-1/2}$.
• Alternating tensor on two-dimensional surfaces orthogonal to n and v: $\chi_{ab} = \eta_{abc} v^c$.
• Expansion of boost v (i.e. velocity divergence): $\kappa = D_a v^a$.
• Shear of v: $\beta_{ab} = D_{\langle a}v_{b \rangle}$.
• Vorticity of v: $W_a = \frac{1}{2}\eta_{abc} D^{a}v^{b}$.

• Electric part of Weyl tensor: $E_{ab} = C_{acbd} n^c n^d$.
• Scalar part of E_ab : $\mathcal{E} = E_{ab} m^a m^b$.
• Vector part of E_ab : $\mathcal{E}_a = M_{ab} m_c E^{bc}$.
• Tensor part of E_ab : $\mathcal{E}_{ab} = \left(M^{\ c}_{(a}M^{\ d}_{b)} - \frac{1}{2}M_{ab}M^{cd}\right)E_{cd}$.
• Magnetic part of Weyl tensor: $H_{ab} = \frac{1}{2}\eta_{efda}C^{ef}_{\quad bc} n^c n^d$.
• Scalar part of H_ab : $\mathcal{H} = H_{ab} m^a m^b$.
• Vector part of H_ab : $\mathcal{H}_a = M_{ab} m_c H^{bc}$.
• Tensor part of H_ab : $\mathcal{H}_{ab} = \left(M^{\ c}_{(a}M^{\ d}_{b)} - \frac{1}{2}M_{ab}M^{cd}\right)H_{cd}$.

• Energy density: $\mu = T_{ab} n^a n^b$.
• Isotropic pressure: $p = \frac{1}{3} T_{ab} f^{ab}$.
• Heat flux/momentum density: $q_a = - T_{bc} n^b f_a^{\ c}$.
• Scalar part of heat flux: $Q = q_a m^a$.
• Vector part of heat flux: $q_a = M_{ab} q^b$.
• Anisotropic stress: $\pi_{ab} = T_{\langle ab \rangle} = \left(f_{(a}^{\ c} f_{b)}^{\ d} - \frac{1}{3}f_{ab}f^{cd}\right)T_{cd}$.
• Scalar part of anisotropic stress: $\Pi = \pi_{ab} m^a m^b$.
• Vector part of anisotropic stress: $\Pi_a = M_{ab} m_c \pi^{bc}$.
• Tensor part of anisotropic stress: $\Pi_{ab} = \left(M^{\ c}_{(a}M^{\ d}_{b)} - \frac{1}{2}M_{ab}M^{cd}\right)\pi_{cd}$.