# Spherical and Plane Symmetric Solutions in Macroscopic Gravity

Anish Agashe<EMAIL_ADDRESS>Department of Physics, The University
of Texas at Dallas,
800 W Campbell Rd, Richardson, TX, 75080, USA Mustapha Ishak
<EMAIL_ADDRESS>Department of Physics, The University of Texas at Dallas,
800 W Campbell Rd, Richardson, TX, 75080, USA

###### Abstract

The theory of macroscopic gravity provides a formalism to average the Einstein
field equations from small scales to largest scales in space-time. It is well
known that averaging is an operation that does not commute with calculating
the Einstein tensor and this leads to a correction term in the field equations
known as back-reaction. In this work, we derive exact solutions to the
macroscopic gravity field equations assuming that the averaged geometry is
plane or spherically symmetric, and the source is taken as vacuum, dust, or
perfect fluid. We then focus on the specific cases of spherical symmetry and
derive solutions that are analogous to the Schwarzschild, Tolman VII, and
Lemaître-Tolman-Bondi solutions. The geodesic equations and curvature
structure are contrasted with the general relativistic counterparts for the
Schwarzschild and Lemaître-Tolman-Bondi solutions.

Keywords: Averaging problem, macroscopic gravity, back-reaction, inhomogeneous
cosmology

##  1 Introduction

An outstanding problem in general relativity (GR) is how the small scale lumpy
universe is averaged to be homogeneous at the largest scales. This problem is
not easily reconcilable due to the difficulties in defining a mathematically
rigorous averaging procedure for non-trivial geometries such as the pseudo-
Riemannian geometry of general relativity. Moreover, the Einstein tensor is a
non-linear function of the metric which would mean that the average of the
Einstein tensor for a given metric is not equal to the Einstein tensor
constructed with the average of that metric, i.e.,
$\left\langle\boldsymbol{E}[\textsl{{g}}]\right\rangle\neq\boldsymbol{E}[\left\langle\textsl{{g}}\right\rangle]$.
These issues are collectively known as the averaging problem in general
relativity [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12].

Therefore, when dealing with large length scales, a set of averaged equations
should be used. In these equations, a correction term should be present to
account for the non-commutativity between averaging and calculating the
Einstein tensor [1, 2, 3, 4, 5, 6]. This term encapsulates the effect of the
microscopic structure of the space-time on the dynamics of the macroscopic
one, and is sometimes referred to as back-reaction. Macroscopic gravity (MG)
is a formalism for covariant space-time averaging [13, 14, 15, 16, 17, 18, 19]
which offers a promising solution to the averaging problem. It is an exact
(non-perturbative) approach which employs averaging bivectors and Lie dragging
of the averaging regions. It is valid for arbitrary classical tensor fields on
a general $n$-dimensional differentiable manifold.

The field equations of MG are the macroscopic analogue of the Einstein field
equations (EFEs). They are characterised by a tensorial correction term which
is a combination of various traces of a quantity called connection
correlation. They can be written in the form of EFEs by taking this correction
term to the right hand side of the equations. Then, the dynamics of a
(assumed) macroscopic space-time geometry can be analysed by determining this
correction term. This can be done by imposing some reasonable assumptions on
the connection correlation [20, 21, 22]. The solutions to the averaged field
equations in MG, with the macroscopic geometry assumed to be FLRW, have been
presented in [20, 21, 22, 23]. Other exact macroscopic geometries like Bianchi
Type-I [23] and static spherically symmetric (Schwarzschild) [24] have also
been explored. Linear perturbations around the FLRW geometry have been
analysed in [22]. In [25, 26, 27], the authors took the microscopic geometry
to be spherically symmetric and wrote it in volume preserving coordinates.
Then, they averaged it using Zalaletdinov’s procedure which becomes trivial in
the volume preserving coordinates, since the averaging operators reduce to
Kronecker delta [15].

In this paper, we derive solutions to the MG field equations where the
macroscopic geometry is assumed to be spherical or plane symmetric (that is,
the less symmetric microscopic geometry averages out to be spherical/plane
symmetric). We start by assuming that the macroscopic geometry admits the
$G_{3}$ group of motion (on $V_{2}$) which contains two special cases – plane
and spherical symmetry111See chapters 15 and 16 in [28] and the references
within for the details. [28], and then, calculate the MG correction term for
such geometries. Our approach to calculate the correction term is similar to
the one taken in [24] to study static spherically symmetric solutions. Here,
we extend the treatment in [24] to include non-static spherically symmetric
and plane symmetric geometries as well. We find that the MG correction term
takes a form of an anisotropic fluid with a qualitative behaviour of an
effective spatial curvature in the field equations. Recently, such geometries
were also analysed in the context of other averaging formalisms [29].

We categorise the solutions based on the source for the space-time – vacuum,
dust, and perfect fluid. Within these three categories, we treat, in detail,
the cases of the static spherically symmetric vacuum solution (Schwarzschild
exterior [30]), the static spherically symmetric perfect fluid solutions
(Schwarzschild interior [30] and Tolman VII [31, 32, 33, 34, 35, 36, 37, 38,
39, 40]), and the non-static spherically symmetric dust solution (Lemaître-
Tolman-Bondi (LTB) [41, 42, 43, 44, 45, 46]). These solutions are the exact
solutions to the MG field equations and do not represent the same geometry as
their counterparts in GR, since they are not the solutions of Einstein field
equations. This is in the same spirit as that of GR – the matter distribution
(at the macroscopic level) should determine the (macroscopic) geometry through
the macroscopic field equations.

The paper is arranged in the following manner: In section 2, we briefly review
the macroscopic gravity formalism. In section 3, we present the spherical and
plane symmetric solutions to the MG field equations. Then, in sections 4, 5
and 6, we analyse specific solutions by assuming the source to be a vacuum,
dust and perfect fluid, respectively. The MG analogue of the Schwarzschild,
Tolman VII and LTB solutions are derived. Finally, we discuss the results in
this paper and make some general remarks in section 7.

The notation and the convention used are as follows: objects associated with
the microscopic geometry are denoted by lowercase letters and those with the
macroscopic geometry by uppercase letters. Greek indices run from $0$ to $3$
and Latin indices from $1$ to $3$. Angular brackets
$\left\langle\cdots\right\rangle$ denote the averaging operation or sometimes
averaged quantities. Covariant differentiation with respect to the macroscopic
connection is denoted by $||$. Indices with round brackets $(~{})$/square
brackets $[~{}]$ are symmetrised/anti-symmetrised; and underlined indices are
not included in (anti-)symmetrisation. The covariant derivative is denoted by
$\nabla_{\mu}$. The sign convention followed is the ‘Landau-Lifshitz Space-
like Convention (LLSC)’ [47]. That is, the signature of the metric is taken to
be Lorentzian $(-,+,+,+)$, the Riemann curvature tensor is defined as,
${r^{\mu}}_{\alpha\nu\beta}~{}=~{}2\partial_{[\nu}{\gamma^{\mu}}_{\underline{\alpha}\beta]}+2{\gamma^{\epsilon}}_{\alpha[\beta}{\gamma^{\mu}}_{\underline{\epsilon}\nu]}$
and $r_{\mu\nu}~{}=~{}{r^{\alpha}}_{\mu\alpha\nu}$ is the Ricci tensor. The
Ricci scalar is defined as
$r~{}=~{}{r^{\mu}}_{\mu}~{}=~{}\textnormal{{g}}^{\mu\nu}r_{\mu\nu}$. Finally,
the units are taken such that $G=1=c$, i.e., $\kappa=8\pi$.

##  2 Macroscopic Gravity Field Equations and their Solutions

Using the concepts of macroscopic electrodynamics [48, 49, 50], a covariant
averaging procedure was developed by Zalaletdinov [13, 14, 15, 16, 17, 18,
19], which can be used in general relativity. It is a generalisation of
averaging on Minkowski space-time and is based on the concept of Lie dragging
of the averaging regions. This procedure is valid for arbitrary classical
tensor fields on any differentiable manifold [13, 14, 15, 16, 17, 18, 19].

Let there be a geometric object (vector, tensor etc.) $p^{\alpha}_{\beta}(x)$
defined on an $n$-dimensional differentiable metric manifold
$(\mathcal{M},\textnormal{{g}}_{\alpha\beta})$. Then, the space-time averaged
value of this object over a compact region $\Sigma\subset\mathcal{M}$ with a
volume $n$-form around a supporting point $x\in\Sigma$, is defined as,

$\left\langle
p^{\alpha}_{\beta}(x)\right\rangle=\dfrac{\int_{\Sigma}\mathcal{A}^{\alpha}_{\mu^{\prime}}(x,x^{\prime})p^{\mu^{\prime}}_{\nu^{\prime}}(x^{\prime})\mathcal{A}^{\nu^{\prime}}_{\beta}(x^{\prime},x)\sqrt{-\textnormal{{g}}^{\prime}}\
\textrm{d}^{n}x^{\prime}}{\int_{\Sigma}\sqrt{-\textnormal{{g}}^{\prime}}\
\textrm{d}^{n}x^{\prime}}$ (2.1)

where $\int_{\Sigma}\sqrt{-\textnormal{{g}}^{\prime}}\
\textrm{d}^{n}x^{\prime}$ is the volume ($V_{\Sigma}$) of the region $\Sigma$
and
$\textnormal{{g}}^{\prime}=\det\left[\textnormal{{g}}_{\alpha\beta}(x^{\prime})\right]$.
The integration is done over all the points $x^{\prime}\in\Sigma$. The
integrand
$\mathcal{A}^{\alpha}_{\mu^{\prime}}(x,x^{\prime})p^{\mu^{\prime}}_{\nu^{\prime}}(x^{\prime})\mathcal{A}^{\nu^{\prime}}_{\beta}(x^{\prime},x)$
is called the _bilocal extension_ of the object $p^{\alpha}_{\beta}(x)$; and
the objects $\mathcal{A}^{\alpha}_{\mu^{\prime}}(x,x^{\prime})$ and
$\mathcal{A}^{\nu^{\prime}}_{\beta}(x^{\prime},x)$ are called the _bilocal
averaging operators_.

Then, bilocal objects $\boldsymbol{\mathcal{F}}$ are defined in terms of these
averaging operators as,

${\mathcal{F}^{\alpha}}_{\beta\rho}=\mathcal{A}^{\alpha}_{\epsilon^{\prime}}\left(\partial_{\rho}\mathcal{A}^{\epsilon^{\prime}}_{\beta}+\nabla_{\sigma^{\prime}}\mathcal{A}^{\epsilon^{\prime}}_{\beta}\mathcal{A}^{\sigma^{\prime}}_{\rho}\right)$
(2.2)

These objects behave as connection coefficients at $x$ and hence, can be
considered as the bilocal extension of the microscopic connection
coefficients. The averages of these objects,
$\left\langle{\mathcal{F}^{\alpha}}_{\beta\rho}\right\rangle$ serve as the
affine connection coefficients of the averaged space-time.

Further, the average of the microscopic Riemann curvature tensor,
$\left\langle{r^{\alpha}}_{\beta\rho\sigma}\right\rangle$ is written as
${R^{\alpha}}_{\beta\rho\sigma}$. A curvature tensor corresponding to the
macroscopic connection can be calculated, and turns out to be related to the
average of the microscopic Riemann curvature tensor through the following
formula,

${M^{\alpha}}_{\beta\rho\sigma}={R^{\alpha}}_{\beta\rho\sigma}+2\left\langle{\mathcal{F}^{\delta}}_{\beta[\rho}{\mathcal{F}^{\alpha}}_{\underline{\delta}\sigma]}\right\rangle-2\left\langle{\mathcal{F}^{\delta}}_{\beta[\rho}\right\rangle\left\langle{\mathcal{F}^{\alpha}}_{\underline{\delta}\sigma]}\right\rangle$
(2.3)

Then, in order to define splitting rules for averaging the differential
Bianchi identities, correlation tensors need to be defined. The connection
correlation is a six rank tensor which takes into account the non-
commutativity of averaging and calculating the Einstein tensor. It is defined
as,

${Z^{\alpha}}_{\beta[\gamma}{{}^{\mu}}_{\underline{\nu}\sigma]}=\left\langle{\mathcal{F}^{\alpha}}_{\beta[\gamma}{\mathcal{F}^{\mu}}_{\underline{\nu}\sigma]}\right\rangle-\left\langle{\mathcal{F}^{\alpha}}_{\beta[\gamma}\right\rangle\left\langle{\mathcal{F}^{\mu}}_{\underline{\nu}\sigma]}\right\rangle$
(2.4)

Note that the second term on the right hand side in equation (2.3) is merely a
contraction of the connection correlation and can be written as,
${Q^{\alpha}}_{\beta\rho\sigma}~{}=~{}2\left\langle{\mathcal{F}^{\delta}}_{\beta[\rho}{\mathcal{F}^{\alpha}}_{\underline{\delta}\sigma]}\right\rangle-2\left\langle{\mathcal{F}^{\delta}}_{\beta[\rho}\right\rangle\left\langle{\mathcal{F}^{\alpha}}_{\underline{\delta}\sigma]}\right\rangle=2{Z^{\epsilon}}_{\beta[\rho}{{}^{\alpha}}_{\underline{\epsilon}\sigma]}$.
This tensor is like a curvature tensor in its mathematical construction and
follows the algebraic properties of one. It can be considered as a ‘curvature
deformation’ tensor since it measures the difference between the macroscopic
curvature tensor and the average of the microscopic curvature tensor. Higher
order correlations also exist (restricted up to the dimensionality of the
macroscopic space-time), but can be taken to be zero in a self consistent
manner [13, 14]. This makes the formalism more practical with equations
involving only the connection correlation.

The final form of the averaged Einstein field equations is given by,

${E^{\epsilon}}_{\gamma}=G^{\beta\epsilon}M_{\beta\gamma}-{\frac{1}{2}}\delta^{\epsilon}_{\gamma}G^{\mu\nu}M_{\mu\nu}=\kappa{T^{\epsilon}}_{\gamma}+\left({Z^{\epsilon}}_{\mu\nu\gamma}-{\frac{1}{2}}\delta^{\epsilon}_{\gamma}Q_{\mu\nu}\right)G^{\mu\nu}$
(2.5)

where,
${T^{\epsilon}}_{\gamma}=\left\langle{t^{\epsilon}}_{\gamma}\right\rangle$ is
the average of the energy-momentum tensor for the microscopic matter
distribution, $G^{\mu\nu}$ is the inverse of the macroscopic metric tensor,
$M_{\beta\gamma}$ is the macroscopic Ricci tensor,
${Z^{\epsilon}}_{\mu\nu\gamma}~{}=~{}2{Z^{\epsilon}}_{\mu[\alpha}{{}^{\alpha}}_{\underline{\nu}\gamma]}$
is a Ricci-tensor-like object for the connection correlation and
$Q_{\mu\nu}~{}=~{}{Q^{\epsilon}}_{\mu\epsilon\nu}$.

Therefore, in the theory of macroscopic gravity, averaging out the Einstein
field equations introduces additional terms on right hand side of these
equations, constituting various traces (contractions) of the connection
correlation. We call these additional terms as the MG correction term, denoted
as,

${C^{\epsilon}}_{\gamma}=\left({Z^{\epsilon}}_{\mu\nu\gamma}-\frac{1}{2}\delta^{\epsilon}_{\gamma}Q_{\mu\nu}\right)G^{\mu\nu}$
(2.6)

To derive a solution to the averaged field equations (2.5), we need to
determine the form of this correction term explicitly for a given macroscopic
geometry. This is done by virtue of several algebraic and differential
constraints on the connection correlation [13, 14],

$\displaystyle{Z^{\alpha}}_{\beta\gamma}{{}^{\mu}}_{\nu\sigma}$
$\displaystyle=-{Z^{\alpha}}_{\beta\sigma}{{}^{\mu}}_{\nu\gamma}$ (2.7a)
$\displaystyle{Z^{\alpha}}_{\beta\gamma}{{}^{\mu}}_{\nu\sigma}$
$\displaystyle=-{Z^{\mu}}_{\nu\gamma}{{}^{\alpha}}_{\beta\sigma}$ (2.7b)
$\displaystyle{Z^{\alpha}}_{\beta[\gamma}{{}^{\mu}}_{\nu\sigma]}$
$\displaystyle=0$ (2.7c)
$\displaystyle{Z^{\epsilon}}_{\epsilon\gamma}{{}^{\mu}}_{\nu\sigma}$
$\displaystyle=0$ (2.7d)
$\displaystyle{{{Z^{\alpha}}_{\beta\mu}}^{\gamma}}_{\delta\nu}u^{\nu}$
$\displaystyle=0$ (2.7e)
$\displaystyle{Z^{\alpha}}_{\beta[\gamma}{{}^{\mu}}_{\underline{\nu}\sigma||\lambda]}$
$\displaystyle=0$ (2.7f)
${Z^{\epsilon}}_{\beta[\gamma}{{}^{\mu}}_{\underline{\nu}\sigma}{M^{\alpha}}_{\underline{\epsilon}\lambda\rho]}-{Z^{\alpha}}_{\epsilon[\gamma}{{}^{\mu}}_{\underline{\nu}\sigma}{M^{\epsilon}}_{\underline{\beta}\lambda\rho]}+{Z^{\alpha}}_{\beta[\gamma}{{}^{\epsilon}}_{\underline{\nu}\sigma}{M^{\mu}}_{\underline{\epsilon}\lambda\rho]}-{Z^{\alpha}}_{\beta[\gamma}{{}^{\mu}}_{\underline{\epsilon}\sigma}{M^{\epsilon}}_{\underline{\nu}\lambda\rho]}=0$
(2.7g)

where, $\boldsymbol{u}$ is the unit time-like 4-vector field and
$\boldsymbol{M}$ is the Riemann curvature tensor for the (assumed) macroscopic
space-time geometry.

Equations (2.7a)-(2.7e) are just algebraic constraints and do not depend on
the macroscopic geometry. Solving these leaves us with 121 independent
components in the connection correlation. Equations (2.7e), (2.7f) and (2.7g)
ensure that the higher order correlations are zero. Then, all one has to deal
with is the connection correlation. There are several assumptions on the
structure of the macroscopic geometry and the functional form of the
connection correlation that go into solving the latter two of these [24].
Given these assumptions, the connection correlation can be determined
completely. A systematic way of solving the MG equations has been presented in
[16, 24, 21, 23].

In the following sections, starting with the 121 components, we solve the last
two equations to determine the MG correction term for spherical and plane
symmetric macroscopic geometries.

##  3 Spherical and Plane Symmetric Solutions

For a coordinate system $(x^{0},x^{1},x^{2},x^{3})$, the general form of the
line element for a spherical and plane symmetric space-time is given by [28,
51],

$\textrm{d}s^{2}=G_{\mu\nu}\textrm{d}x^{\mu}\textrm{d}x^{\nu}=-\textrm{e}^{2P}(\textrm{d}x^{0})^{2}+\textrm{e}^{2Q}(\textrm{d}x^{1})^{2}+R^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{3})^{2}\right]$
(3.1)

where, $G_{\mu\nu}$ is the macroscopic metric, $P,Q,R$ are, in general,
arbitrary functions of coordinates $(x^{0},x^{1})$, and just the coordinate
$x^{1}$ for the static case. The time-like unit 4-vector field (the 4-velocity
of time-like particles) admitted by such space-time is given by,
$u^{\nu}=\left[\textrm{e}^{-P},0,0,0\right]$. The function $S(x^{2},k)$ takes
the values $\\{\sin(x^{2}),1\\}$ for $k=\\{1,0\\}$, where $k~{}=~{}KR^{2}$ and
$K$ is the constant Gaussian curvature of the 2-subspace ($x^{0},x^{1}=$
const.). Since the function, $R$, has a fixed value on this 2-subspace, one
can normalise the curvature in this manner [28]. The spherical and plane
symmetric cases are represented by $k=1$ and $k=0$, respectively.

It is conventional to use spherical coordinates $(t,r,\theta,\phi)$ in the
case of spherical symmetry $(k=1)$. Then, the line element takes the form,

$\textrm{d}s^{2}=-\textrm{e}^{2P}\textrm{d}t^{2}+\textrm{e}^{2Q}\textrm{d}r^{2}+R^{2}\left(\textrm{d}\theta^{2}+\sin^{2}\theta\textrm{d}\phi^{2}\right)$
(3.2)

where $P,Q,R$ are arbitrary functions of time, $t$, and the radial coordinate,
$r$, for the non-static case. In the static case, $P,Q$ are functions of only
$r$ and the function $R$ can be chosen to be equal to $r$ [28].

For the case of plane symmetry $(k=0)$, we often use Cartesian coordinates
$(t,x,y,z)$. Then, the line element takes the form,

$\textrm{d}s^{2}=-\textrm{e}^{2P}\textrm{d}t^{2}+\textrm{e}^{2Q}\textrm{d}x^{2}+R^{2}\left(\textrm{d}y^{2}+\textrm{d}z^{2}\right)$
(3.3)

where $P,Q,R$ are arbitrary functions of time, $t$, and the spatial
coordinate, $x$, for the non-static case222Note that, one can choose any one
of the three spatial coordinates to be $x^{1}$. Then, the metric will be
modified accordingly and so would the MG correction term.. In the static case,
$P,Q,R$ are functions of only $x$.

We calculate the MG correction term with the line element for the macroscopic
geometry to be the one in equation (3.1)333We did the same calculations
separately for equations (3.2) and (3.3) for both static and non-static cases.
The results were identical to the respective cases presented here.. We take
the components of connection correlation to be dependent on all the
coordinates, i.e.,
$\boldsymbol{Z}~{}\equiv\boldsymbol{Z}(x^{0},x^{1},x^{2},x^{3})$. We start
with the 121 independent components in $\boldsymbol{Z}$ that remain after
solving constraints (2.7a)-(2.7e). Solving the differential constraint (2.7f)
leaves only 7 independent components. The integrability condition (2.7g) is
trivially satisfied and does not reduce the independent components any
further.

These components are characterised by arbitrary functions of coordinates
$x^{2}$ and $x^{3}$. To make the MG correction term compatible with the
macroscopic geometry, two of these need to be constant. Then, out of the 7
independent components in $\boldsymbol{Z}$, five are arbitrary functions and
two are constants. We label these as $f_{n}(x^{2},x^{3}),\ (n=1,...,5)$ and
$b_{1},b_{2}$. The number of non-zero components in the connection correlation
tensor is 44 and that in the curvature deformation tensor, $\boldsymbol{Q}$,
is 8. The non-Riemannian curvature tensor, $\boldsymbol{R}$, has 44 non-zero
components comprising of three out of the seven variables in $\boldsymbol{Z}$.
We do not list them here as the real important quantity is the MG correction
term, since it enters the field equations. The MG correction term is diagonal
and takes the following form,

$\displaystyle{C^{0}}_{0}$ $\displaystyle=\frac{\mathcal{B}}{R^{2}}$ (3.4a)
$\displaystyle{C^{0}}_{i}$ $\displaystyle=0={C^{i}}_{0}$ (3.4b)
$\displaystyle{C^{1}}_{1}$
$\displaystyle=\frac{\mathcal{B}}{R^{2}}\;;\quad{C^{2}}_{2}=0={C^{3}}_{3}$
(3.4c)

Writing in matrix form, this looks like,

${C^{\epsilon}}_{\gamma}=\frac{\mathcal{B}}{R^{2}}\begin{bmatrix}&1&0&0&0&\\\
&0&1&0&0&\\\ &0&0&0&0&\\\ &0&0&0&0&\end{bmatrix}$ (3.5)

where, $\mathcal{B}=2(b_{1}+b_{2})$ is a constant.

The MG correction term takes the form of an anisotropic fluid. Only the radial
component of the effective pressure is non-zero and is equal in magnitude, but
opposite in sign, to the effective energy density. The MG correction term for
the static spherically symmetric geometry in our case is the same as what was
found in [24], although we have applied fewer restrictions on
$\boldsymbol{Z}$.

###  3.1 The Macroscopic Einstein Field Equations

We take the average of the microscopic energy-momentum tensor to be that of an
inhomogeneous relativistic perfect fluid with an isotropic pressure,

${T^{\epsilon}}_{\gamma}=(\rho+p)u^{\epsilon}u_{\gamma}+p\delta^{\epsilon}_{\gamma}$
(3.6)

where, $\rho\equiv\rho(x^{0},x^{1})$ and $p\equiv p(x^{0},x^{1})$ are the
energy density and pressure. Taking the averaged matter distribution to be a
perfect fluid does not reduce the generality of the equations that follow (for
the purposes of this paper), since both vacuum and dust can be considered as
special cases. The conservation of energy-momentum444The differential
constraint (2.7f) of the formalism ensures that the MG correction term has
also vanishing covariant derivative [13, 14]. Hence, the covariant divergence
of the two sides of equation (2.5) vanishes identically – the left side due to
the usual differential Bianchi identity.,
$\nabla_{\epsilon}{T^{\epsilon}}_{\gamma}=0$, gives two equations – one each
for $\gamma~{}=~{}0,1$, while other components are identically satisfied,

$\displaystyle\frac{\dot{\rho}}{\rho+p}+\dot{Q}+2\frac{\dot{R}}{R}$
$\displaystyle=0$ (3.7a) $\displaystyle\frac{p^{\prime}}{\rho+p}+P^{\prime}$
$\displaystyle=0$ (3.7b)

Then, without assuming any equation of state, $p(\rho)$, the macroscopic field
equations for the metric in (3.1) read as

${E^{0}}_{0}=\frac{\textrm{e}^{-2Q}}{R^{2}}\left[2RR^{\prime\prime}+R^{\prime}\left(R^{\prime}-2RQ^{\prime}\right)\right]-\frac{\textrm{e}^{-2P}}{R^{2}}\left[\dot{R}\left(\dot{R}+2R\dot{Q}\right)\right]-\frac{k}{R^{2}}=-8\pi\rho+\frac{\mathcal{B}}{R^{2}}$
(3.8a)
${E^{0}}_{1}=\frac{2\textrm{e}^{-2P}}{R}\left[\dot{R}^{\prime}-\dot{Q}R^{\prime}-P^{\prime}\dot{R}\right]=0$
(3.8b)
${E^{1}}_{0}=-\frac{2\textrm{e}^{-2Q}}{R}\left[\dot{R}^{\prime}-\dot{Q}R^{\prime}-P^{\prime}\dot{R}\right]=0$
(3.8c)
${E^{1}}_{1}=-\frac{\textrm{e}^{-2P}}{R^{2}}\left[2R\ddot{R}+\dot{R}\left(\dot{R}-2R\dot{P}\right)\right]+\frac{\textrm{e}^{-2Q}}{R^{2}}\left[R^{\prime}\left(R^{\prime}+2RP^{\prime}\right)\right]-\frac{k}{R^{2}}=8\pi
p+\frac{\mathcal{B}}{R^{2}}$ (3.8d)
${E^{2}}_{2}={E^{3}}_{3}=\frac{\textrm{e}^{-2Q}}{R}\left[R^{\prime\prime}+R^{\prime}\left(P^{\prime}-Q^{\prime}\right)+R\left(P^{\prime\prime}+{P^{\prime}}^{2}-P^{\prime}Q^{\prime}\right)\right]\\\
-\frac{\textrm{e}^{-2P}}{R}\left[\ddot{R}+\dot{R}\left(\dot{Q}-\dot{P}\right)+R\left(\ddot{Q}+{\dot{Q}}^{2}-\dot{P}\dot{Q}\right)\right]=8\pi
p$ (3.8e)

where, ‘dot’ and ‘prime’ represent differentiation with respect to the
coordinates $x^{0}$ and $x^{1}$, respectively. The plane and spherical
symmetric cases correspond to $k=0$ and $k=1$, respectively. The MG correction
term enters equations (3.8a) and (3.8d) as an additional curvature term. Other
equations are not affected. In the case of static solutions, all the terms
with a dot go to zero.

##  4 Vacuum

In the case of vacuum, the conditions from conservation of energy are
identically satisfied. The MG correction term does not change the field
equations qualitatively since the constant term $\mathcal{B}$ can be absorbed
in the curvature term. Therefore, the usual results for the analysis of vacuum
solutions hold. However, the terms with $k$ would be replaced by
$k+\mathcal{B}$. Therefore, the generalised Birkhoff theorem [28, 52, 53, 54,
55, 56] would still hold and the metric would be necessarily static. This
means that all the terms with ‘dot’ vanish. Moreover, if
$\partial_{\mu}R\partial^{\mu}R>0$, one can choose canonical coordinates,
where $R=x^{1}$ [28, 51]. Then, solving equation (3.8a), we would find,

$\textrm{e}^{2Q}=\left[(k+\mathcal{B})-\frac{2m}{x^{1}}\right]^{-1}$ (4.1)

where, $m$ is a constant that comes from integration and $Q$ is a function of
only $x^{1}$ now. Further, substituting equation (3.8d) in (3.8a), we get,

$P^{\prime}=-Q^{\prime}\ \Rightarrow P=-Q+C_{1}$

where $C_{1}$ is a constant of integration. This gives us,

$\textrm{e}^{2P}=C_{2}\textrm{e}^{-2Q}=C_{2}\left[(k+\mathcal{B})-\frac{2m}{x^{1}}\right]$
(4.2)

where, $C_{2}=\textrm{e}^{2C_{1}}$. For simplicity, we can choose $C_{2}=1$,
as one does in the usual Schwarzschild solution. Then, using equations (4.1)
and (4.2) in equation (3.1), the macroscopic metric takes the following form,

$\textrm{d}s^{2}=-f(x^{1},k,\mathcal{B})(\textrm{d}x^{0})^{2}+\frac{1}{f(x^{1},k,\mathcal{B})}(\textrm{d}x^{1})^{2}+(x^{1})^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{3})^{2}\right]$
(4.3)

where $f=\left[(k+\mathcal{B})-\frac{2m}{x^{1}}\right]$. Therefore, the vacuum
solution is now modified to have the back-reaction term in the line element.
This metric has a coordinate singularity at
$x^{1}=\frac{2m}{(k+\mathcal{B})}$.

###  4.1 The MG-Schwarzschild Exterior Solution

The static spherically symmetric case ($k=1$) gives us the MG-Schwarzschild
(exterior) solution. The line element, in coordinates $(t,r,\theta,\phi)$,
then, looks like,

$\textrm{d}s^{2}=-\left[(1+\mathcal{B})-\frac{2m}{r}\right]\textrm{d}t^{2}+\frac{1}{\left[(1+\mathcal{B})-\frac{2m}{r}\right]}\textrm{d}r^{2}\\\
+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$ (4.4)

It reduces to the usual Schwarzschild exterior solution of GR when we take the
back-reaction term $\mathcal{B}$ to be zero. The line element in equation
(4.4) clearly represents a geometry that is different from that of the
Schwarzschild solution. We know that in the GR-Schwarzschild geometry, the
Ricci tensor and hence the Ricci scalar is zero. But this is not the case in
the MG-Schwarzschild geometry where, ${R^{\mu}}_{\nu}~{}=~{}{\rm
diag}[0,0,-\frac{\mathcal{B}}{r^{2}},-\frac{\mathcal{B}}{r^{2}}]$ and
$R~{}=~{}-\frac{2\mathcal{B}}{r^{2}}$. The Kretschmann scalar is also
modified, and is given by,
$K~{}=~{}\frac{48m^{2}}{r^{6}}~{}+~{}\frac{4\mathcal{B}r(\mathcal{B}r-4m)}{r^{6}}$.
This is exactly the same result as in [24].

As explained above, since the geometry is not Schwarzschild anymore, the
geodesic equations will be modified as well. Using the line element in
equation (4.4), the Lagrangian density will become,

$2\mathcal{L}=-f(r,\mathcal{B})\left(\frac{\textrm{d}t}{\textrm{d}\lambda}\right)^{2}+\frac{1}{f(r,\mathcal{B})}\left(\frac{\textrm{d}r}{\textrm{d}\lambda}\right)^{2}+r^{2}\left(\frac{\textrm{d}\phi}{\textrm{d}\lambda}\right)^{2}$
(4.5)

where, $\lambda$ is an affine parameter and
$f=\left(1+\mathcal{B}-\frac{2m}{r}\right)$. Using the Euler-Lagrange
equations, we get,
$\frac{\textrm{d}t}{\textrm{d}\lambda}~{}=~{}\frac{\gamma}{f(r,\mathcal{B})};\
\frac{\textrm{d}\phi}{\textrm{d}\lambda}~{}=~{}\frac{l}{r^{2}}$, where,
$\gamma$ and $l$ are constants. We know that for geodesics
$2\mathcal{L}=\varepsilon$, where, $\varepsilon=-1,0,+1$ correspond to time-
like, null and space-like geodesics, respectively. Using this, we can obtain
the following orbital equation for the respective type of particles,

$\left(\frac{\textrm{d}u}{\textrm{d}\phi}\right)^{2}=2mu^{3}-(1+\mathcal{B})u^{2}-\frac{2m\varepsilon}{l^{2}}u+\frac{\gamma^{2}+(1+\mathcal{B})\varepsilon}{l^{2}}$
(4.6)

which gives,

$\frac{\textrm{d}^{2}u}{\textrm{d}\phi^{2}}=3mu^{2}-(1+\mathcal{B})u-\frac{m\varepsilon}{l^{2}}$
(4.7)

where, $u\equiv\frac{1}{r}$. The above equation can be written in the form of
the usual equation in GR by rescaling the coordinate $\phi$ and mass $m$ as
$\overline{\phi}=\sqrt{1+\mathcal{B}}\phi$ and
$\overline{m}=\frac{m}{1+\mathcal{B}}$, respectively. The solution to this
equation is, then, given by,

$u=\frac{-m\varepsilon}{(1+\mathcal{B})l^{2}}\left[1+\alpha\cos\left\\{\phi\left(\sqrt{1+\mathcal{B}}+\frac{3m^{2}\varepsilon}{(1+\mathcal{B})^{\frac{3}{2}}l^{2}}\right)\right\\}\right]$
(4.8)

where, $\alpha$ is a small parameter ($|\alpha|<<1$) that quantifies the
deviation from the Newtonian solution.

It is a fruitful exercise to write the metric in equation (4.4) in various
other coordinates that are frequently used to write the GR-Schwarzschild
solution. This will be particularly useful when one wants to study the
behaviour of various geodesics in the MG-Schwarzschild exterior space-time. In
‘isotropic’ coordinates [57], we have,

$\textrm{d}s^{2}=-(1+\mathcal{B})\left(\frac{1-\frac{m}{2(1+\mathcal{B})\overline{r}^{\sqrt{1+\mathcal{B}}}}}{1+\frac{m}{2(1+\mathcal{B})\overline{r}^{\sqrt{1+\mathcal{B}}}}}\right)^{2}\textrm{d}t^{2}+\overline{r}^{2\sqrt{1+\mathcal{B}}-2}\left(1+\frac{m}{2(1+\mathcal{B})\overline{r}^{\sqrt{1+\mathcal{B}}}}\right)^{4}\left(\textrm{d}\overline{x}^{2}+\textrm{d}\overline{y}^{2}+\textrm{d}\overline{z}^{2}\right)$
(4.9)

where,
$r=\overline{r}^{\sqrt{1+\mathcal{B}}}\left(1+\frac{m}{2(1+\mathcal{B})\overline{r}^{\sqrt{1+\mathcal{B}}}}\right)^{2}$
and we have, $\overline{x}~{}=~{}\overline{r}\sin\theta\cos\phi$;
$\overline{y}~{}=~{}\overline{r}\sin\theta\cos\phi$;
$\overline{z}~{}=~{}\overline{r}\cos\theta$. In the Eddington-Finkelstein
coordinates [58, 59, 60], we have,

$\textrm{d}s^{2}=-\left[(1+\mathcal{B})-\frac{2m}{r}\right]\textrm{d}v^{2}+2\textrm{d}v\textrm{d}r+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(4.10)

where,
$v=t+\frac{r}{1+\mathcal{B}}+\frac{2m}{(1+\mathcal{B})^{2}}\ln\left|\frac{(1+\mathcal{B})r}{2m}-1\right|$.
Similarly, we can write,

$\textrm{d}s^{2}=-\left[(1+\mathcal{B})-\frac{2m}{r}\right]\textrm{d}u^{2}-2\textrm{d}u\textrm{d}r+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(4.11)

where,
$u=t-\frac{r}{1+\mathcal{B}}-\frac{2m}{(1+\mathcal{B})^{2}}\ln\left|\frac{(1+\mathcal{B})r}{2m}-1\right|$.
In Kruskal-Szekeres [61, 62] coordinates, we have,

$\textrm{d}s^{2}=-\frac{32m^{3}}{r(1+\mathcal{B})^{4}}\textrm{e}^{\frac{-(1+\mathcal{B})r}{2m}}\textrm{d}U\textrm{d}V+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(4.12)

where, $V=\textrm{e}^{\frac{v(1+\mathcal{B})^{2}}{4m}}$ and
$U=-\textrm{e}^{\frac{-u(1+\mathcal{B})^{2}}{4m}}$. This can be further
written as,

$\textrm{d}s^{2}=-\frac{32m^{3}}{r(1+\mathcal{B})^{4}}\textrm{e}^{\frac{-(1+\mathcal{B})r}{2m}}(\textrm{d}T^{2}-\textrm{d}X^{2})+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(4.13)

where, $T={\frac{1}{2}}(V+U)$ and $X={\frac{1}{2}}(V-U)$.

It should be possible to also write an equivalent form of the MG-Schwarzschild
metric in coordinate systems like that of Lemaître-Novikov [41, 63] and Israel
[64].

##  5 Dust

In the case of dust, the pressure, $p$, in equations (3.7) and (3.8) is zero.
Then, equation (3.7b) means that $P^{\prime}=0$. This implies that $P\equiv
P(x^{0})$. We can divide the solutions into two classes: i) when
$R^{\prime}=0$; $\dot{R}\neq 0$ and, ii) when $R^{\prime}\neq 0$.

###  5.1 Solutions with $R^{\prime}=0,\ \dot{R}\neq 0$

When _$R^{\prime}=0,\ \dot{R}\neq 0$_ , one can choose $R=x^{0}$ [28] and then
integrate equation (3.8d) to get,

$\textrm{e}^{2P}=\left\\{\frac{C_{1}}{x^{0}}-(k+\mathcal{B})\right\\}^{-1}=\frac{x^{0}}{C_{1}-(k+\mathcal{B})x^{0}}$
(5.1)

where, $C_{1}$ is a constant of integration.

Then, assuming $\textrm{e}^{Q}=U(x^{0},x^{1})\textrm{e}^{-P}$, equations
(3.8e) can be written as,

$\ddot{U}+\dot{U}\left(\frac{1}{x^{0}}-3\dot{P}\right)=0$ (5.2)

This can be integrated to find,

$U=F_{1}\int\frac{\textrm{e}^{3P}}{x^{0}}\textrm{d}x^{0}+F_{2}$ (5.3)

where, $F_{1},F_{2}$ are arbitrary functions of the coordinate $x^{1}$ that
arise due to integration. Using equation (5.1), the integration in the above
equation can be solved to get,

$\int\frac{\textrm{e}^{3P}}{x^{0}}\textrm{d}x^{0}=\int\frac{2\textrm{e}^{2P}}{1+(k+\mathcal{B})\textrm{e}^{2P}}\textrm{d}(\textrm{e}^{P})\\\
=\frac{2}{(k+\mathcal{B})}\left\\{\textrm{e}^{P}-\frac{\tan^{-1}\left(\sqrt{k+\mathcal{B}}\
\textrm{e}^{P}\right)}{\sqrt{k+\mathcal{B}}}\right\\}+C_{2}$ (5.4)

where, $C_{2}$ is a constant of integration. Using this, we get,

$\textrm{e}^{Q}=F_{3}-\textrm{e}^{-P}\left\\{\frac{F_{3}\tan^{-1}\left(\sqrt{k+\mathcal{B}}\
\textrm{e}^{P}\right)}{\sqrt{k+\mathcal{B}}}-F_{4}\right\\}$ (5.5)

where, $F_{3}(x^{1})=\frac{2F_{1}}{k+\mathcal{B}}$ and
$F_{4}(x^{1})=C_{2}F_{1}+F_{2}$.

Using equation (5.1) in (3.8a), we get the equation for the energy density,

$\frac{2\textrm{e}^{-2P}}{x^{0}}\left(\dot{P}+\dot{Q}\right)$ (5.6)

If $F_{3}=0$ (i.e., $F_{1}=0$), we recover the vacuum solutions
($\textrm{e}^{P}~{}=~{}\textrm{e}^{-Q}$). Assuming that $F_{3}\neq 0$, we can
rescale the coordinate $x^{1}$ such that $F_{3}=1$. Then the line element in
(3.1) takes the form,

$\textrm{d}s^{2}=-\frac{x^{0}}{C_{1}-(k+\mathcal{B})x^{0}}(\textrm{d}x^{0})^{2}+\left[1-\textrm{e}^{-P}\left\\{\frac{\tan^{-1}\left(\sqrt{k+\mathcal{B}}\
\textrm{e}^{P}\right)}{\sqrt{k+\mathcal{B}}}-F_{4}\right\\}\right]^{2}(\textrm{d}x^{1})^{2}\\\
+(x^{0})^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{3})^{2}\right]$
(5.7)

This reduces to the solution in general relativity when the back-reaction
term, $\mathcal{B}$, is zero. In general relativity, this subclass of
solutions contains generalisations of the Kantowski-Sachs solution (for $k=1$)
[65, 28] and Bianchi Type-I solution (for $k=0$). It would be interesting to
see if the same is true for the MG analogues of these solutions. Given that
the back-reaction term in MG couples with the geometry in a non-trivial
manner, we expect these limits to either change or even not work at all.

###  5.2 Solutions with $R^{\prime}\neq 0$

When $R^{\prime}\neq 0$, we can rescale the coordinate $x^{0}$ such that
$G_{00}=1$. This is equivalent to taking $P=0$ and then defining the ‘dots’ to
be the differentiation with respect to this new rescaled coordinate. Then, we
can integrate equation (3.8b) (or (3.8c)) to get,

$Q=\ln(R^{\prime})+\ln F_{1}$ (5.8)

where, $F_{1}(x^{1})$ is an arbitrary function depending on the boundary
conditions.

By conveniently redefining the function $F_{1}$ as,
$F_{1}(x^{1})~{}\equiv~{}\frac{1}{\sqrt{(k+\mathcal{B})-\varepsilon f^{2}}}$,
we can write,

$\textrm{e}^{2Q}=\frac{{R^{\prime}}^{2}}{(k+\mathcal{B})-\varepsilon
f^{2}(x^{1})}$ (5.9)

The function $f$ is an arbitrary function of the coordinate $x^{1}$ and the
variable, $\varepsilon=(0,\pm 1)$, is to be chosen such that $\textrm{e}^{2Q}$
remains positive. This redefinition makes the other dynamical equations
simpler. These expressions are valid for both static and non-static cases.

In the static case, all the terms with dots go to zero. Then, equation (3.8d)
implies that $\varepsilon f^{2}=0$. When used in equation (3.8a), we get,
$\rho=0$. Therefore, the static solutions in the case of a pressureless fluid
invariably reduce to vacuum solutions which, in the case of spherical symmetry
($k=1$), leads to nothing but the Schwarzschild solution in the Lemaître-
Novikov [28, 42, 63] coordinates.

In the non-static case, using equations (5.8) and (5.9) in (3.8a) gives,

$\frac{{\dot{R}}^{2}+\varepsilon
f^{2}}{R^{2}}+\frac{2\dot{R}\dot{R}^{\prime}+2\varepsilon
ff^{\prime}}{RR^{\prime}}=8\pi\rho$ (5.10)

Then, equation (3.8d) gives the following differential equation for $R$,

$2R\ddot{R}+{\dot{R}}^{2}=-\varepsilon f^{2}$ (5.11)

Performing an integration once, this equation gives,

${\dot{R}}^{2}-\frac{2m(x^{1})}{R}=-\varepsilon f^{2}$ (5.12)

where, $m(x^{1})$ is another arbitrary function coming from integration. This
equation can be completely integrated for both the cases, $\varepsilon=0,\
\varepsilon\neq 0$ (see chapter 15 in [28] for a complete analysis). Using
equation (5.11) and its differentiation, the other two field equations (3.8e)
are identically satisfied.

Using equation (5.12) to substitute values of $\varepsilon f^{2}$ and
$\varepsilon ff^{\prime}$ in equation (5.11), we get555Note that equation
(5.13) can be rearranged to give,  $2m(x^{1})=8\pi\int\rho
R^{2}R^{\prime}\textrm{d}x^{1}$ This equation resembles that for a mass
function and hence justifies the labelling of the function arising from the
integration as ‘$m$’.,

$\frac{2m^{\prime}}{R^{2}R^{\prime}}=8\pi\rho$ (5.13)

Similarly, if we use equation (5.10) to do the same, we get another equation
for the density, in terms of the function $R$ only,

$\frac{\ddot{R}^{\prime}}{R^{\prime}}+\frac{2\ddot{R}}{R}=-4\pi\rho$ (5.14)

Using equation (5.9), the line element in equation (3.1) takes the form,

$\textrm{d}s^{2}=-(\textrm{d}x^{0})^{2}+\frac{{R^{\prime}}^{2}}{(k+\mathcal{B})-\varepsilon
f^{2}}(\textrm{d}x^{1})^{2}+R^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{3})^{2}\right]$
(5.15)

Therefore, the dust solution is modified to have a back-reaction term. The
usual dust solution can be obtained by putting the back-reaction term to be
zero.

####  5.2.1 The MG-LTB Solution

Our analysis until this point is done for both plane and spherical symmetry in
a combined way. Although there has been some interest in the plane symmetric
cosmological models in the past, none of these models has stood the tests of
precision cosmology as well as the spatially homogeneous and isotropic FLRW
model. The only exciting alternative to the FLRW model is the Lemaître-Tolman-
Bondi (LTB) [41, 43, 44] model. The LTB model is a subcase (for
$R^{\prime}\neq 0$) of the spherically symmetric dust solutions to the
Einstein field equations. It is only radially inhomogeneous (due to spherical
symmetry) with the spatial hypersurfaces spherically symmetric about a central
worldline. In this section, we will look at the LTB model in the framework of
macroscopic gravity – the MG-LTB model. To do this, we choose spherical
coordinates $(t,r,\theta,\phi)$ and $k=1$ in equations (5.9) and (5.15). The
line element then looks like,

$\textrm{d}s^{2}=-\textrm{d}t^{2}+\frac{{R^{\prime}}^{2}}{(1+\mathcal{B})-\varepsilon
f^{2}}\textrm{d}r^{2}+R^{2}\left[\textrm{d}\theta^{2}+\sin^{2}\theta\textrm{d}\phi^{2}\right]$
(5.16)

This is a modified form of LTB metric [66, 67] where there is a back-reaction
term in the denominator of the radial metric coefficient.

The dynamical equations governing this model would be given by equations
(5.10) - (5.14). Comparing them to the Friedman equations, we can define local
cosmological parameters within the LTB models [68, 69, 70, 71],

$\displaystyle H(t,r)$ $\displaystyle=\frac{\dot{R}(t,r)}{R(t,r)}$ (5.17a)
$\displaystyle\Omega_{m_{0}}(r)$
$\displaystyle=\frac{2m(r)}{H_{0}^{2}(r)R_{0}^{3}(r)}$ (5.17b)
$\displaystyle\Omega_{c_{0}}(r)$ $\displaystyle=\frac{-\varepsilon
f^{2}}{R_{0}^{2}(r)H_{0}^{2}(r)}$ (5.17c)

where, $H,\ \Omega_{m_{0}},\ \Omega_{c_{0}}$ are the local Hubble parameter,
matter density parameter and curvature density parameter, respectively. The
functions, $H_{0}(r)~{}\equiv~{}H(t_{0},r)$ and
$R_{0}(r)~{}\equiv~{}R(t_{0},r)$, are the boundary values at the present time,
$t_{0}$. Using equation (5.17) in (5.12), we can get an expression for the
expansion history in LTB models,

$H^{2}(t,r)=H_{0}^{2}(r)\left[\Omega_{m_{0}}(r)\left(\frac{R_{0}}{R}\right)^{3}\right.\\\
+\left.\Omega_{c_{0}}(r)\left(\frac{R_{0}}{R}\right)^{2}\right]$ (5.18)

Integrating this equation, we get an equation for the age of the LTB universe,

$t_{0}-t=\frac{1}{H_{0}}\int_{R}^{R_{0}}\frac{\textrm{d}R}{\left(\Omega_{m_{0}}R_{0}^{2}R^{-1}+\Omega_{c_{0}}R_{0}^{2}\right)^{{\frac{1}{2}}}}$
(5.19)

The next thing to look at, from the point of view of observational cosmology,
would be light propagation. For a radial null geodesic,
$\textrm{d}\theta=0=\textrm{d}\phi$, and from equation (5.16), we get a
constraint equation for the light rays [68, 69, 70, 71],

$\frac{\textrm{d}t}{\textrm{d}\lambda}=-\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{R^{\prime}(t,r)}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}$ (5.20)

where, $\lambda$ is an affine parameter and the negative sign tells that we
are dealing with incoming light rays. Let $t(\lambda)$ and $\tau(\lambda)$ be
two solutions to the equation above. Then, $t+\tau$ would be a solution too.
Then, we can write,

$\displaystyle\frac{\textrm{d}t}{\textrm{d}\lambda}$
$\displaystyle=-\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{R^{\prime}(t,r)}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}$ (5.21a) $\displaystyle\frac{\textrm{d}(t+\tau)}{\textrm{d}\lambda}$
$\displaystyle=-\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{R^{\prime}(t,r)}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}+\frac{\textrm{d}\tau(\lambda)}{\textrm{d}\lambda}$ (5.21b)
$\displaystyle\frac{\textrm{d}(t+\tau)}{\textrm{d}\lambda}$
$\displaystyle=-\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{R^{\prime}(t,r)+\dot{R}^{\prime}(t,r)\tau}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}$ (5.21c)

Using equations (5.21b) and (5.21c), we get,

$\frac{1}{\tau}\frac{\textrm{d}\tau}{\textrm{d}\lambda}=-\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{\dot{R}^{\prime}(t,r)}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}$ (5.22)

Defining the redshift, $z$, as,
$1+z~{}\equiv~{}\frac{\tau(0)}{\tau(\lambda)}$, we can write,

$\frac{\textrm{d}z}{\textrm{d}\lambda}=(1+z)\frac{\textrm{d}r}{\textrm{d}\lambda}\frac{\dot{R}^{\prime}(t,r)}{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}$ (5.23)

This gives us,

$\frac{\textrm{d}r}{\textrm{d}z}=\frac{1}{1+z}\frac{\sqrt{1+\mathcal{B}-\varepsilon
f^{2}}}{\dot{R}^{\prime}(t,r)}$ (5.24)

Further, using equations (5.20) and (5.23), we get,

$\frac{\textrm{d}t}{\textrm{d}z}=-\frac{1}{1+z}\frac{R^{\prime}(t,r)}{\dot{R}^{\prime}(t,r)}$
(5.25)

It is, then, straightforward to calculate the comoving, angular diameter and
luminosity distances, which are given by [68, 69, 70, 71],

$\displaystyle d_{C}(z)$ $\displaystyle=(1+z)R(t(z),r(z))$ (5.26a)
$\displaystyle d_{A}(z)$ $\displaystyle=R(t(z),r(z))$ (5.26b) $\displaystyle
d_{L}(z)$ $\displaystyle=(1+z)^{2}R(t(z),r(z))$ (5.26c)

These equations tell us that the back-reaction term will modify the redshift
dependence of the coordinates, $t$ and $r$, and hence, the distance
calculations.

Looking at equations (5.17) - (5.19), we can see that the equations governing
the MG-LTB solution are the same as the ones for the GR-LTB solution. But the
geometry of MG-LTB model (equation (5.16)) is not the same. This points to an
important feature of models within macroscopic gravity, which is that the
back-reaction (the connection correlation) can be set in such a way that it
only affects either the averaged geometry or the averaged evolution but not
both. For example, we could have chosen the integration function $F_{1}$ in
equation (5.8) to be simply $\frac{1}{\sqrt{k-\varepsilon f^{2}}}$, and then,
equation (5.12) would have had a back-reaction term,
$\frac{\mathcal{B}}{R^{2}}$, in it. The FLRW limit of the LTB model is also
affected in a similar fashion – either we end up with a geometry that is not
FLRW but equation (5.12) reduces to the usual Friedmann equation or we have an
FLRW geometry that evolves differently. In the next section, we will see that
something similar happens in the case of the homogeneous perfect fluid
solutions as well.

##  6 Perfect Fluid

In the case of perfect fluid, the dynamical equations are the same as
equations (3.7) and (3.8). Perfect fluid solutions are more complicated to
solve analytically. In GR, the additional condition due to the isotropy of the
pressure looks like, ${E^{1}}_{1}-{E^{3}}_{3}=0$. However, in MG, this
condition gets modified due to the anisotropy in the effective pressure due to
the back-reaction and we have,
${E^{1}}_{1}-{E^{3}}_{3}=\frac{\mathcal{B}}{R^{2}}$.

###  6.1 Static Solutions

In the static case, all the terms with ‘dots’ in equations (3.7) and (3.8) go
to zero and we can assume $R=x^{1}$ [28, 51]. Then, equation (3.8a) becomes,

${E^{0}}_{0}=\left[x^{1}(k+\mathcal{B}-\textrm{e}^{-2Q})\right]^{\prime}=8\pi\rho(x^{1})^{2}$
(6.1)

We can easily integrate this equation to get,

$\textrm{e}^{2Q}=\left[(k+\mathcal{B})-\frac{2m(x^{1})}{x^{1}}\right]^{-1}$
(6.2)

where, $m(x^{1})$ is a mass function defined as,

$2m(x^{1})=8\pi\int\rho(x^{1})(x^{1})^{2}\textrm{d}x^{1}$ (6.3)

The condition due to isotropy,
${E^{1}}_{1}-{E^{3}}_{3}=\frac{\mathcal{B}}{R^{2}}$, further gives,

$P^{\prime\prime}+(P^{\prime})^{2}-P^{\prime}Q^{\prime}\\\
=\frac{1}{(x^{1})^{2}}\left[1+x^{1}(P^{\prime}+Q^{\prime})-(k+\mathcal{B})\textrm{e}^{2Q}\right]$
(6.4)

Using equations (6.2), (6.3) and (6.4) in (3.8e), we get,

$2x^{1}\left[(k+\mathcal{B})x^{1}-2m\right]P^{\prime}=8\pi p(x^{1})^{3}+2m$
(6.5)

Now, using equation (3.7b) to eliminate $P^{\prime}$ from the equation above,
we get,

$2x^{1}\left[(k+\mathcal{B})x^{1}-2m\right]p^{\prime}\\\ =-(\rho+p)\left[8\pi
p(x^{1})^{3}+2m\right]$ (6.6)

This is the MG analogue of the Tolman-Oppenheimer-Volkoff (TOV) equation [31,
72].

The analysis until here is quite general. However, as in the case of vacuum
solutions, we are interested in spherically symmetric solutions. For that, we
use spherical coordinates $(t,r,\theta,\phi)$ and $k=1$ in equations
(6.1)-(6.6). We will look, in detail, at two static spherically symmetric
perfect fluid solutions – the Schwarzschild interior solution [30] (constant
density) and the Tolman VII solution [31] (variable density).

####  6.1.1 MG-Schwarzschild Interior Solution

The Schwarzschild interior solution is, perhaps, the best known static
spherically symmetric perfect fluid solution in GR. The MG analogue of this
solution has been presented in [24]. It is characterised by a constant energy
density $(\rho~{}=~{}\overline{\rho}~{}=~{}{\rm constant})$ inside a boundary,
say, $r~{}\leq~{}r_{b}$. Then, using equation (6.3), the mass function inside
this boundary becomes,

$2m(r)=8\pi\int^{r}_{0}\overline{\rho}r^{2}\textrm{d}r=\frac{8\pi}{3}\overline{\rho}r^{3}$
(6.7)

and then, using equation (6.2), we get,

$\textrm{e}^{2Q(r)}=\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]^{-1}$
(6.8)

Further, using the constant density term in equation (3.7b), we can integrate
to get,

$p=C_{1}\textrm{e}^{-P}-\overline{\rho}$ (6.9)

Then, using this equation in (6.5), we can integrate to find,

$\textrm{e}^{P}=\frac{3}{2\overline{\rho}}\left(C_{1}-C_{2}\textrm{e}^{-Q}\right)\\\
=\frac{3}{2\overline{\rho}}\left[C_{1}-C_{2}\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right\\}^{{\frac{1}{2}}}\right]$
(6.10)

Putting this equation back in the expression for pressure above (6.9), we get,

$p=\frac{\overline{\rho}}{3}\left(\frac{3\
C_{2}\textrm{e}^{-Q(r)}-C_{1}}{C_{1}-C_{2}\textrm{e}^{-Q(r)}}\right)=\frac{\overline{\rho}}{3}\frac{3\
C_{2}\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]^{{\frac{1}{2}}}-C_{1}}{C_{1}-C_{2}\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]^{{\frac{1}{2}}}}$
(6.11)

The values for the constants $C_{1},C_{2}$ can be determined by matching the
interior solution to the exterior solution. The expression for
$\textrm{e}^{2Q}$ already matches the one for the exterior solution. Requiring
that, for some boundary $r_{b}$, $p(r_{b})=0$, gives us,
$\frac{C_{1}}{C_{2}}~{}=~{}3\textrm{e}^{-Q(r_{b})}$. Using this, the
expression for pressure becomes,

$p=\overline{\rho}\left(\frac{\textrm{e}^{-Q(r)}-\textrm{e}^{-Q(r_{b})}}{3\
\textrm{e}^{-Q(r_{b})}-\textrm{e}^{-Q(r)}}\right)=\overline{\rho}\frac{\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]^{{\frac{1}{2}}}-\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r_{b}^{2}\right]^{{\frac{1}{2}}}}{3\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r_{b}^{2}\right]^{{\frac{1}{2}}}-\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]^{{\frac{1}{2}}}}$
(6.12)

Further, we require that $\textrm{e}^{P(r_{b})}\left|{}_{\rm
int}\right.~{}=~{}\textrm{e}^{P(r_{b})}\left|{}_{\rm ext}\right.$. Using this,
we can completely specify the constants as,

$C_{1}=\overline{\rho}\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r_{b}^{2}\right\\}^{{\frac{1}{2}}}\;;\quad
C_{2}=\frac{\overline{\rho}}{3}$ (6.13)

Then, equation (6.10) becomes,

$\textrm{e}^{P}={\frac{1}{2}}\left[3\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r_{b}^{2}\right\\}^{{\frac{1}{2}}}\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right\\}^{{\frac{1}{2}}}\right]$
(6.14)

Therefore, the MG-Schwarzschild interior solution is characterised by a
constant density $\overline{\rho}$ and a pressure given by equation (6.12).
The metric coefficients are given by equations (6.8) and (6.14). The line
element for the interior solution then looks like,

$\textrm{d}s^{2}=-\frac{1}{4}\left[3\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r_{b}^{2}\right\\}^{{\frac{1}{2}}}-\left\\{(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right\\}^{{\frac{1}{2}}}\right]^{2}\textrm{d}t^{2}\\\
+\frac{1}{\left[(1+\mathcal{B})-\frac{8\pi}{3}\overline{\rho}r^{2}\right]}\textrm{d}r^{2}+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(6.15)

The Schwarzschild interior and exterior solutions derived here match the ones
presented in [24]. However, presenting all the calculations again is justified
since we are working with fewer restrictions on the connection correlation.
This adds to the completeness of the analysis in this paper. In the next
section, we will present an interior solution with non-constant density.

####  6.1.2 MG-Tolman VII Solution

The obvious solution to look at, after the Schwarzschild solution, is the
Tolman VII solution [31]. The Tolman VII solution is characterised by the
following ansatz on the energy density [39],

$\rho=\rho_{0}\left[1-\beta\left(\frac{r}{r_{b}}\right)^{2}\right]$ (6.16)

where, $\rho_{0}$ is the central density ($\rho_{0}\equiv\rho(r=0)$), $r_{b}$
represents a boundary radius beyond which the solution can be considered to be
the Schwarzschild exterior solution, and $\beta$ is a dimensionless parameter
that takes values between 0 and 1.

Using this ansatz for the energy density, the mass function (equation (6.3))
inside the boundary becomes,

$2m(r)=8\pi\rho_{0}\left(\frac{r^{3}}{3}-\frac{\beta
r^{5}}{5r_{b}^{2}}\right)$ (6.17)

and then, using equation (6.2), we get666In the original paper by Tolman [31],
this form for the metric coefficient was assumed and the expression for the
density was derived using that.,

$\textrm{e}^{2Q}=\left[(1+\mathcal{B})-8\pi\rho_{0}\left(\frac{r^{2}}{3}-\frac{\beta}{5r_{b}^{2}}r^{4}\right)\right]^{-1}$
(6.18)

We define a new variable $x\equiv\frac{r^{2}}{r_{b}^{2}}$. Then, the above
equation can be written as,

$\textrm{e}^{2Q}=\left[(1+\mathcal{B})-A\left(\frac{x}{3}-\frac{\beta}{5}x^{2}\right)\right]^{-1}$
(6.19)

where, $A=8\pi\rho_{0}r_{b}^{2}$, is simply a constant.

Rearranging equation (6.4), we get,

$\textrm{e}^{-2Q}\left[1+rP^{\prime}+rQ^{\prime}(1+rP^{\prime})\right.\\\
\left.-r^{2}(P^{\prime\prime}+{P^{\prime}}^{2})\right]=(1+\mathcal{B})$ (6.20)

Using equation (6.18) and its differentiation with respect to $r$ to find the
value of $Q^{\prime}\textrm{e}^{-2Q}$, the above equation becomes,

$\left[(1+\mathcal{B})-A\left(\frac{x}{3}-\frac{\beta}{5}x^{2}\right)\right]r^{2}(P^{\prime\prime}+{P^{\prime}}^{2})\\\
-\left[(1+\mathcal{B})-A\frac{\beta}{5}x^{2}\right]rP^{\prime}-A\frac{\beta}{5}x^{2}=0$
(6.21)

Then, following [33, 34], we define a function $U(r)$, such that,
$P~{}\equiv~{}\ln U$ and simplify further to get,

$\left(\frac{5(1+\mathcal{B})}{\beta
A}-\frac{5x}{3\beta}+x^{2}\right)r^{2}U^{\prime\prime}\\\
-\left(\frac{5(1+\mathcal{B})}{\beta A}-x^{2}\right)rU^{\prime}-x^{2}U=0$
(6.22)

Writing the differential terms to be with respect to the variable, $x$, we
get,

$\underbrace{\left(\frac{5(1+\mathcal{B})}{\beta
A}-\frac{5x}{3\beta}+x^{2}\right)}_{M}\widetilde{\widetilde{U}}+\underbrace{\left(x-\frac{5}{6\beta}\right)}_{N}\widetilde{U}-\frac{U}{4}=0$
(6.23)

where ‘tilde’ represents differentiation with respect to $x$ and we have used,
$U^{\prime}=\frac{2r}{r_{b}^{2}}\tilde{U}$ and
$U^{\prime\prime}=\frac{2r}{r_{b}^{2}}\tilde{U}+\frac{4r^{2}}{r_{b}^{4}}\tilde{\tilde{U}}$.

Now, we define a variable $w~{}=~{}\ln(\sqrt{M}+N)$. This gives,
$\tilde{w}~{}=~{}\frac{1}{\sqrt{M}}$, and, $\tilde{w}~{}=~{}-NM\sqrt{M}$.
Using this, we get,
$\widetilde{U}~{}=~{}\frac{1}{\sqrt{M}}\frac{\textrm{d}U}{\textrm{d}w}$, and,
$\widetilde{\widetilde{U}}~{}=~{}\frac{1}{M}\frac{\textrm{d}^{2}U}{\textrm{d}w^{2}}~{}-~{}\frac{N}{M\sqrt{M}}\frac{\textrm{d}U}{\textrm{d}w}$.
Substituting these in equation (6.23), we get,

$\frac{\textrm{d}^{2}U}{\textrm{d}w^{2}}+\frac{U}{4}=0$ (6.24)

This is the equation of a simple harmonic oscillator with frequency
${\frac{1}{2}}$ and has the following solution,

$U=C_{1}\cos\left({\frac{1}{2}}w\right)+C_{2}\sin\left({\frac{1}{2}}w\right)$
(6.25)

where, $C_{1},C_{2}$ are constants of integration. Then using the definition
of the function $U$, we can write,

$\textrm{e}^{2P}=\left[C_{1}\cos\left({\frac{1}{2}}w\right)+C_{2}\sin\left({\frac{1}{2}}w\right)\right]^{2}$
(6.26)

where,
$w~{}=~{}\ln\left\\{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{4}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}\right\\}$.

Substituting the expression for the mass function (equation (6.17)) and using
equation (6.25) to find $P^{\prime}=\frac{U^{\prime}}{U}$ in equation (6.5),
we get an expression for the pressure,

$p=\frac{1}{r_{b}}\left(\frac{\beta\rho_{0}}{10\pi}\right)^{{\frac{1}{2}}}\left\\{(1+\mathcal{B})-8\pi\rho_{0}\left(\frac{r^{2}}{3}-\frac{\beta}{5r_{b}^{2}}r^{4}\right)\right\\}^{{\frac{1}{2}}}\left\\{\frac{C_{2}-C_{1}\tan\left({\frac{1}{2}}w\right)}{C_{1}+C_{2}\tan\left({\frac{1}{2}}w\right)}\right\\}-\rho_{0}\left(\frac{1}{3}-\frac{\beta}{5r_{b}^{2}}r^{2}\right)$
(6.27)

The constants $C_{1},C_{2}$ can be determined using the boundary conditions:
$\textrm{e}^{2P(r_{b})}\left|{}_{\rm
tol}\right.~{}=~{}\textrm{e}^{2P(r_{b})}\left|{}_{\rm ext}\right.$ and
$p(r_{b})~{}=~{}0$. The first one gives,

$C_{1}\cos\left({\frac{1}{2}}w_{b}\right)+C_{2}\sin\left({\frac{1}{2}}w_{b}\right)\\\
=\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}$
(6.28a) Using this, the second boundary condition gives,
$C_{2}\cos\left({\frac{1}{2}}w_{b}\right)-C_{1}\sin\left({\frac{1}{2}}w_{b}\right)\\\
=r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)$
(6.28b)

where, $w_{b}\equiv w(r_{b})$. Solving equations (6.28), we get,

$C_{1}=\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}\cos\left({\frac{1}{2}}w_{b}\right)\\\
-r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)\sin\left({\frac{1}{2}}w_{b}\right)$
(6.29a)
$C_{2}=\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}\sin\left({\frac{1}{2}}w_{b}\right)\\\
+r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)\cos\left({\frac{1}{2}}w_{b}\right)$
(6.29b)

where,
$w_{b}=\ln\left\\{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}\right\\}$.

Using these expressions for the constants and simplifying even further,
equations (6.26) and (6.27) become,

$\textrm{e}^{2P}=\left[\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}\cos\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}\right.\\\
\left.+r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)\sin\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}\right]^{2}$
(6.30)

and,

$p=\frac{\beta\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}}{10\pi
r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)}\left\\{(1+\mathcal{B})-8\pi\rho_{0}\left(\frac{r^{2}}{3}-\frac{\beta}{5r_{b}^{2}}r^{4}\right)\right\\}^{{\frac{1}{2}}}\\\
\times\left[\frac{1+r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)\tan\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}}{1-\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}\tan\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}}\right]\\\
-\rho_{0}\left(\frac{1}{3}-\frac{\beta}{5r_{b}^{2}}r^{2}\right)$ (6.31)

The line element for this solution then looks like,

$\textrm{d}s^{2}=\left[\left\\{(1+\mathcal{B})-8\pi\rho_{0}r_{b}^{2}\left(\frac{1}{3}-\frac{\beta}{5}\right)\right\\}^{{\frac{1}{2}}}\cos\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}\right.\\\
\left.+r_{b}\left(\frac{10\pi\rho_{0}}{\beta}\right)^{{\frac{1}{2}}}\left(\frac{1}{3}-\frac{\beta}{5}\right)\sin\left\\{\ln\left(\frac{\frac{r^{2}}{r_{b}^{2}}-\frac{5}{6\beta}+\left(\frac{r^{4}}{r_{b}^{2}}-\frac{5r^{2}}{3\beta
r_{b}^{2}}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}{1-\frac{5}{6\beta}+\left(1-\frac{5}{3\beta}+\frac{5(1+\mathcal{B})}{8\pi\beta\rho_{0}r_{b}^{2}}\right)^{{\frac{1}{2}}}}\right)^{{\frac{1}{2}}}\right\\}\right]^{2}\textrm{d}t^{2}\\\
+\frac{1}{\left[(1+\mathcal{B})-8\pi\rho_{0}\left(\frac{r^{2}}{3}-\frac{\beta}{5r_{b}^{2}}r^{4}\right)\right]}\textrm{d}r^{2}+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}$
(6.32)

The MG-Tolman VII solution is completely specified by equations (6.16),
(6.18), (6.30), (6.31) and (6.32). It reduces to the usual Tolman VII solution
in [33] if we put the back-reaction term to be zero.

Note that the static perfect fluid solutions presented here are not the most
general for their assumed form of density, since we took the constant of
integration in the mass function to be zero. In the more general solutions,
there would be a term $\propto\frac{1}{r}$ in the metric coefficient $G_{11}$
(equations (6.8) and (6.18)) which would lead to a singularity at $r=0$. We
only mention this for the sake of completeness and these cases are beyond the
scope of this paper.

###  6.2 Non-Static Solutions

The non-static plane or spherically symmetric perfect fluid solutions are
usually classified according to the properties of the kinematic variables [28,
73],

$\displaystyle\Theta$
$\displaystyle=\textrm{e}^{-P}\left(\dot{Q}+2\frac{\dot{R}}{R}\right)$ (6.33a)
$\displaystyle\dot{u}_{\alpha}$ $\displaystyle=\left[0,P^{\prime},0,0\right]$
(6.33b) $\displaystyle\sigma^{0}_{0}$
$\displaystyle=0\;;\;-{\frac{1}{2}}\sigma^{1}_{1}=\sigma^{2}_{2}=\sigma^{3}_{3}=\frac{1}{3}\textrm{e}^{-P}\left(\frac{\dot{R}}{R}-\dot{Q}\right)$
(6.33c) And due to the symmetry here, we have,
$\displaystyle\omega_{\alpha\beta}$ $\displaystyle=0$ (6.33d)

where, $\Theta,\dot{u},\sigma$ and $\omega$ are the 4-velocity’s expansion,
acceleration, shear and rotation (or vorticity), respectively [28, 74]. A wide
variety of solutions exist with different assumptions on these variables.

####  6.2.1 Solutions with $\Theta=0$; ${\sigma^{\alpha}}_{\beta}=0$

The vanishing shear and expansion, through equations (6.33a) and (6.33c),
imply that $\dot{Q}~{}=0~{}=~{}\dot{R}$. Then, equation (3.7a) implies that
the energy density does not depend on $x^{0}$ (given that
$\rho~{}+~{}p~{}\neq~{}0$). This will modify the field equations such that no
derivative with respect to $x^{0}$ is present there. This is the same as the
static case. However, one can still have $\dot{P}~{}\neq~{}0$. Then, for the
pressure, we have $\dot{p}~{}\neq~{}0$, if $\dot{P}^{\prime}~{}\neq~{}0$.
Therefore, these solutions are either static or can be generated from static
solutions [28].

####  6.2.2 Solutions with $\Theta\neq 0$; ${\sigma^{\alpha}}_{\beta}=0$

The most interesting and well studied case is that of expanding, shear free
solutions. We will look at this case in some detail since it has the FLRW
solution as a subcase. Vanishing shear, through equation (6.33c), implies
$\dot{Q}=\frac{\dot{R}}{R}$. This can be integrated to give,

$R=X(x^{1})\textrm{e}^{Q}$ (6.34)

where, $X$ is an arbitrary function arising from integration. Using this in
equation (3.8b) and then integrating, we get,

$\textrm{e}^{P}=\dot{Q}\textrm{e}^{-Y(x^{0})}$ (6.35)

where, $Y$ is another arbitrary function coming from integration. Then, using
this and equation (6.34), the expansion scalar becomes,
$\Theta=3\textrm{e}^{Y}$.

Using equations (6.34) and (6.35) in equations (3.8a) and (3.8d), the energy
density and pressure can be written as,

$3\textrm{e}^{2Y}-\frac{\textrm{e}^{-2Q}}{X^{2}}\left(2XX^{\prime\prime}+{X^{\prime}}^{2}+2X^{2}Q^{\prime\prime}\right.\\\
\left.+4XX^{\prime}Q^{\prime}+X^{2}{Q^{\prime}}^{2}+k+\mathcal{B}\right)=8\pi\rho$
(6.36)

and,

$\frac{\textrm{e}^{-3Q}}{\dot{Q}X^{2}}\partial_{0}\left[\textrm{e}^{Q}\left\\{\left(X^{\prime}+XQ^{\prime}\right)^{2}-(k+\mathcal{B})\right.\right.\\\
\left.\left.-\textrm{e}^{2(Q+Y)}\right\\}\right]=8\pi p$ (6.37)

The isotropy condition, ${E^{1}}_{1}-{E^{3}}_{3}=\frac{\mathcal{B}}{R^{2}}$,
gives,

$\partial_{0}\left[\textrm{e}^{Q}\left\\{{X^{\prime}}^{2}-X\left(X^{\prime\prime}+X^{\prime}Q^{\prime}\right)\right.\right.\\\
\left.\left.+X^{2}\left({Q^{\prime}}^{2}-Q^{\prime\prime}\right)-(k+\mathcal{B})\right\\}\right]=0$
(6.38)

Further, the metric becomes,

$\textrm{d}s^{2}=-\dot{Q}^{2}\textrm{e}^{-2Y}(\textrm{d}x^{0})^{2}+\textrm{e}^{2Q}(\textrm{d}x^{1})^{2}\\\
+\textrm{e}^{2Q}X^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{2})^{2}\right]$
(6.39)

Equation (6.38) can be integrated and then further solved to produce several
classes of solutions [28], especially in the case of spherical symmetry. The
one that is of interest to us is the spherically symmetric solution $((k=1);\
x^{\mu}=(t,r,\theta,\phi))$ with a homogeneous distribution of matter density
and pressure $(\rho~{}\equiv~{}\rho(t);\ p~{}\equiv~{}p(t))$. The homogeneity
of the pressure, through equation (3.7b), implies $P^{\prime}~{}=~{}0$. As
explained in section 5, this allows one to put $P=0$ and write the metric
coefficient $Q$ as in equation (5.9). Further equation (6.35), then, implies
$\dot{Q}^{\prime}~{}=~{}0$. Using this, we can write,

$Q~{}=~{}Q_{1}(t)+Q_{2}(r)$ (6.40)

Using this in equation (6.34), we get,

$R^{\prime}=Z^{\prime}a(t)$ (6.41)

where, $Z(r)=X\textrm{e}^{Q_{2}(r)}$ and $a(t)=e^{Q_{1}(t)}$.

Then, using equations (5.9) and (6.41), the metric in equation (6.39) (for the
spherically symmetric case) becomes,

$\textrm{d}s^{2}=-\textrm{d}t^{2}+\frac{{Z^{\prime}}^{2}a^{2}}{F^{2}(r)}\textrm{d}r^{2}+Z^{2}a^{2}\left(\textrm{d}\theta^{2}+\sin^{2}\theta\textrm{d}\phi^{2}\right)$
(6.42)

where, the function $F$ is an arbitrary function that comes from integration.
Now, since $Z,Z^{\prime}$ are functions of only $r$, we can perform a
coordinate transformation $r=Z(r)$ such that the above metric becomes,

$\textrm{d}s^{2}=-\textrm{d}t^{2}+a^{2}\left(\frac{\textrm{d}r^{2}}{F^{2}(r)}+r^{2}\textrm{d}\theta^{2}+r^{2}\sin^{2}\theta\textrm{d}\phi^{2}\right)$
(6.43)

This is nothing but the spatially homogeneous and isotropic FLRW metric.

The function $F$ can be (conveniently) redefined such that
$F^{2}(r)~{}=~{}(1~{}+~{}\mathcal{B})~{}-~{}Kr^{2}$. Given this redefinition,
the MG field equations reduce to the usual Friedmann equations. On the other
hand, if $F^{2}$ is taken to be simply $1~{}-~{}Kr^{2}$, we recover the usual
FLRW geometry, but now, its evolution is affected by the back-reaction term.
This highlights the obvious, albeit important, point that we discussed briefly
in the case of LTB solution, which is, in MG, the FLRW geometry evolves
differently than in GR and the geometry that evolves like a FLRW geometry is
different from FLRW.

####  6.2.3 Solutions with ${\sigma^{\alpha}}_{\beta}\neq 0$;
$\dot{u}^{\alpha}=0$

The vanishing acceleration, through equation (6.33b), implies that
$P^{\prime}~{}=~{}0$. Then, equation (3.7b) would mean that
$p^{\prime}~{}=~{}0$. Further, equation (3.8b) or (3.8c) gives,

$\dot{R}^{\prime}~{}=~{}\dot{Q}R^{\prime}$ (6.44)

Since the integration of this equation depends on what kind of a function $R$
is, the solutions with shear and vanishing acceleration are divided into the
following three categories:

i) When $R=\textrm{constant}=C(\textrm{say})$, equation (3.8a) and (3.8d)
gives us,

$-8\pi\rho=\frac{k+\mathcal{B}}{C^{2}}=8\pi p$ (6.45)

Using this, equation (3.8e) can be solved to find $Q$. Since
$P~{}\equiv~{}P(x^{0})$, the time coordinate can be rescaled such that $P=0$
(similar to what we had in the case of dust). The line element, in this case,
looks like,

$\textrm{d}s^{2}=-(\textrm{d}x^{0})^{2}+S^{2}\left\\{\frac{x^{0}}{C},-(k+\mathcal{B})\right\\}(\textrm{d}x^{1})^{2}\\\
+C^{2}\left[(\textrm{d}x^{2})^{2}+S^{2}(x^{2},k)(\textrm{d}x^{3})^{2}\right]$
(6.46)

ii) When $R~{}\neq~{}\textrm{constant}$; $R^{\prime}~{}=~{}0$, one can choose
$R~{}=~{}x^{0}$. The field equations then give,

$\displaystyle\frac{\textrm{e}^{-2P}}{(x^{0})^{2}}\left\\{1+2x^{0}\dot{Q}+(k+\mathcal{B})\textrm{e}^{2P}\right\\}$
$\displaystyle=8\pi\rho$ (6.47a)
$\displaystyle\frac{\textrm{e}^{-2P}}{(x^{0})^{2}}\left\\{1-2x^{0}\dot{P}-(k+\mathcal{B})\textrm{e}^{2P}\right\\}$
$\displaystyle=8\pi p$ (6.47b)
$\displaystyle\frac{\textrm{e}^{-2P}}{(x^{0})^{2}}\left\\{1-(x^{0})^{2}\left(\ddot{Q}+\dot{Q}^{2}-\dot{P}\dot{Q}\right)\right\\}$
$\displaystyle=4\pi(\rho+3p)$ (6.47c)

These equations can be solved to find one of the two functions $P,Q$ by
prescribing the other one. We know that due to zero acceleration,
$P^{\prime}~{}=~{}0$. However, it should be noted that, even though there is
no $Q^{\prime}$ in the field equations, it is not zero in general.

iii) When $R^{\prime}~{}\neq~{}0$, one can choose $P~{}=~{}0$. Then, this case
becomes similar to the corresponding case with dust. Equation (6.44) can be
integrated to find,

$\textrm{e}^{2Q}=\frac{{R^{\prime}}^{2}}{(k+\mathcal{B})-\varepsilon
f^{2}}\qquad\quad\varepsilon=0,\pm 1$ (6.48)

The energy density and pressure can be determined from the field equations,
which reduce to,

$\displaystyle\frac{{\dot{R}}^{2}+\varepsilon
f^{2}}{R^{2}}+\frac{2\dot{R}\dot{R}^{\prime}+2\varepsilon
ff^{\prime}}{RR^{\prime}}$ $\displaystyle=8\pi\rho$ (6.49a)
$\displaystyle\frac{\dot{R}^{2}+\varepsilon f^{2}}{R^{2}}+\frac{2\ddot{R}}{R}$
$\displaystyle=-8\pi p$ (6.49b)
$\displaystyle\frac{\ddot{R}^{\prime}}{RR^{\prime}}+\frac{2\ddot{R}}{R}$
$\displaystyle=-4\pi(\rho+3p)$ (6.49c)

The first equation above is the same as the one in the case of dust and the
third equation comes from differentiating the second one with respect to
$x^{1}$. The metric takes the same form as the one on equation (5.15).

##  7 Conclusion

In this paper, we have analysed the effects of averaging on plane and
spherically symmetric macroscopic geometries within the framework of
macroscopic gravity. We found that the back-reaction takes a form of an
anisotropic fluid and enters the field equations in the form of an additional
spatial curvature. This was consistent with the findings in [25, 26, 24] for
spherically symmetric geometries. Here, we have extended the analysis further
to include plane symmetric geometries as well. We categorised the solutions
based on the type of source and then analysed various subcases within.

Taking the source to be vacuum led to the MG-Schwarzschild solution.
Similarly, for dust, the non-static solutions led to the MG-LTB solution and
for a perfect fluid, we derived Schwarzschild interior and Tolman VII (static
case) and some non-static solutions. Our approach in this work was different
from how back-reaction has been treated in the literature in that we
considered the influence of back-reaction on the geometry instead of the
dynamics. Since in MG, the back-reaction modifies the field equations, it will
influence the geometry and not just the dynamical evolution.

## Acknowledgements

AA acknowledges that a code to find exact-FLRW solutions in MG was written by
Tharake Wijenayake and Mustapha Ishak [23], using the computer algebra
software Maple (https://www.maplesoft.com) and the openly available package
GRTensor (https://github.com/grtensor/grtensor). All the results in this paper
were obtained by modification and expansion of this code.

## References

  * [1] G.F.R. Ellis.  Relativistic Cosmology: Its Nature, Aims and Problems.  Fundam. Theor. Phys., 9:215–288, 1984.
  * [2] G.F.R. Ellis and W. Stoeger.  The ’fitting problem’ in cosmology.  Class. Quant. Grav., 4:1697–1729, 1987.
  * [3] G.F.R. Ellis.  83 years of general relativity and cosmology: Progress and problems.  Class. Quant. Grav., 16:A37–A75, 1999.
  * [4] M. F. Shirokov and I. Z. Fisher.  Isotropic Space with Discrete Gravitational-Field Sources. On the Theory of a Nonhomogeneous Isotropic Universe.  Astron. Zh., 39:899, January 1962.
  * [5] M. F. Shirokov and I. Z. Fisher.  Isotropic Space with Discrete Gravitational-Field Sources. On the Theory of a Nonhomogeneous Isotropic Universe.  Sov. Astro., 6:699, April 1963.
  * [6] M. F. Shirokov and I. Z. Fisher.  Isotropic Space with Discrete Gravitational-Field Sources. On the Theory of a Nonhomogeneous Isotropic Universe.  General Relativity and Gravitation, 30(9):1411–1427, September 1998\.
  * [7] Reza K. Tavakol and Roustam Zalaletdinov.  On the domain of applicability of general relativity.  Found. Phys., 28:307–331, 1998.
  * [8] Timothy Clifton.  Back-Reaction in Relativistic Cosmology.  Int. J. Mod. Phys. D, 22:1330004, 2013.
  * [9] Thomas Buchert and Syksy Räsänen.  Backreaction in Late-Time Cosmology.  Annual Review of Nuclear and Particle Science, 62(1):57–79, November 2012.
  * [10] George F.R. Ellis and Thomas Buchert.  The Universe seen at different scales.  Phys. Lett. A, 347:38–46, 2005.
  * [11] George F.R. Ellis.  Inhomogeneity effects in Cosmology.  Class. Quant. Grav., 28:164001, 2011.
  * [12] R.J. van den Hoogen.  Averaging Spacetime: Where do we go from here?  In 12th Marcel Grossmann Meeting on General Relativity, pages 578–588, 3 2010.
  * [13] Roustam M. Zalaletdinov.  Averaging out the Einstein equations.  General Relativity and Gravitation, 24(10):1015–1031, October 1992\.
  * [14] Roustam Zalaletdinov.  Towards a theory of macroscopic gravity.  Gen. Rel. Grav., 25:673–695, 1993.
  * [15] Marc Mars and Roustam M. Zalaletdinov.  Space–time averages in macroscopic gravity and volume-preserving coordinates.  Journal of Mathematical Physics, 38(9):4741–4757, Sep 1997.
  * [16] Roustam M. Zalaletdinov.  Averaged Lagrangians and MacCallum–Taub’s Limit in Macroscopic Gravity.  General Relativity and Gravitation, 28(8):953–979, August 1996\.
  * [17] ROUSTAM ZALALETDINOV.  The averaging problem in cosmology and macroscopic gravity.  International Journal of Modern Physics A, 23(08):1173–1181, Mar 2008.
  * [18] R M Zalaletdinov.  Space-time Averages of Classical Physical Fields.  Ann. Eur. Acad. Sci., (gr-qc/0411004):344–371. 28 p, Oct 2004.
  * [19] Roustam M. Zalaletdinov.  Averaging problem in general relativity, macroscopic gravity and using Einstein’s equations in cosmology.  Bull. Astron. Soc. India, 25:401–416, 1997.
  * [20] A. A. Coley, N. Pelavas, and R. M. Zalaletdinov.  Cosmological solutions in macroscopic gravity.  Physical Review Letters, 95(15), Oct 2005.
  * [21] R. J. van den Hoogen.  A complete cosmological solution to the averaged einstein field equations as found in macroscopic gravity.  Journal of Mathematical Physics, 50(8):082503, Aug 2009.
  * [22] Timothy Clifton, Alan Coley, and Robert van den Hoogen.  Observational cosmology in macroscopic gravity.  Journal of Cosmology and Astroparticle Physics, 2012(10):044–044, Oct 2012.
  * [23] Tharake Wijenayake and Mustapha Ishak.  Expansion and Growth of Structure Observables in a Macroscopic Gravity Averaged Universe.  Phys. Rev. D, 91(6):063534, 2015.
  * [24] R.J. Van Den Hoogen.  Spherically Symmetric Solutions in Macroscopic Gravity.  Gen. Rel. Grav., 40:2213–2227, 2008.
  * [25] A. A. Coley and N. Pelavas.  Averaging spherically symmetric spacetimes in general relativity.  Phys. Rev. D, 74(8):087301, October 2006.
  * [26] A. A. Coley and N. Pelavas.  Averaging in spherically symmetric cosmology.  Phys. Rev. D, 75(4):043506, February 2007.
  * [27] S. Ph. Tegai and I. V. Drobov.  Averaging schwarzschild spacetime.  Phys. Rev. D, 96:024041, Jul 2017.
  * [28] Hans Stephani, D. Kramer, Malcolm A.H. MacCallum, Cornelius Hoenselaers, and Eduard Herlt.  Exact solutions of Einstein’s field equations.  Cambridge Monographs on Mathematical Physics. Cambridge Univ. Press, Cambridge, 2003.
  * [29] Timothy Clifton and Roberto A. Sussman.  Cosmological Backreaction in Spherical and Plane Symmetric Dust-Filled Space-Times.  Class. Quant. Grav., 36(20):205004, 2019.
  * [30] Karl Schwarzschild.  Über das Gravitationsfeld eines Massenpunktes nach der Einsteinschen Theorie.  Sitzungsberichte der Königlich Preußischen Akademie der Wissenschaften (Berlin, pages 189–196, January 1916.  English translation republished as: Schwarzschild, K. “Golden Oldie”: On the Gravitational Field of a Mass Point According to Einstein’s Theory. General Relativity and Gravitation 35, 951–959 (2003), arXiv:physics/9905030v1, Correction by: Senovilla, J.M.M. The Schwarzschild solution: corrections to the editorial note. Gen Relativ Gravit 39, 685–693 (2007).
  * [31] Richard C. Tolman.  Static solutions of Einstein’s field equations for spheres of fluid.  Phys. Rev., 55:364–373, 1939.
  * [32] Max Wyman.  Radially Symmetric Distributions of Matter.  Physical Review, 75(12):1930–1936, June 1949.
  * [33] A. L. Mehra.  Radially symmetric distribution of matter.  Journal of the Australian Mathematical Society, 6(2):153–156, 1966\.
  * [34] M C Durgapal and G L Gehlot.  Spheres with varying density in general relativity.  Journal of Physics A: General Physics, 4(6):749–755, nov 1971.
  * [35] M. C. Durgapal and P. S. Rawat.  Non-rigid massive spheres in general relativity.  Mon. Not. R. Astro. Soc., 192:659–662, September 1980.
  * [36] M. S. R. Delgaty and Kayll Lake.  Physical acceptability of isolated, static, spherically symmetric, perfect fluid solutions of Einstein’s equations.  Computer Physics Communications, 115(2-3):395–415, December 1998\.
  * [37] Nicholas Neary, Mustapha Ishak, and Kayll Lake.  Tolman type VII solution, trapped null orbits, and w-modes.  Phys. Rev. D, 64(8):084001, October 2001.
  * [38] Nicholas Neary and Kayll Lake.  r-modes in the Tolman VII solution.  arXiv e-prints, pages gr–qc/0106056, June 2001.
  * [39] Ambrish M. Raghoonundun and David W. Hobill.  Possible physical realizations of the Tolman VII solution.  Phys. Rev. D, 92(12):124005, December 2015.
  * [40] P. C. Vaidya and G. K. Patwardhan.  Relativistic distributions of matter of radial symmetry.  J. Univ. Bombay, 12N3:23, 1943.
  * [41] Georges Lemaître.  L’Univers en expansion.  Annales de la Société Scientifique de Bruxelles, 53:51, January 1933.
  * [42] George Abbé Lemaître and M. A. H. MacCallum.  The Expanding Universe.  General Relativity and Gravitation, 29(5):641–680, May 1997.
  * [43] Richard C. Tolman.  Effect of inhomogeneity on cosmological models.  Proceedings of the National Academy of Sciences, 20(3):169–176, 1934.
  * [44] H. Bondi.  Spherically Symmetrical Models in General Relativity.  Monthly Notices of the Royal Astronomical Society, 107(5-6):410–425, 12 1947.
  * [45] B. Datt.  Über eine Klasse von Lösungen der Gravitationsgleichungen der Relativität.  Zeitschrift fur Physik, 108(5-6):314–321, May 1938.
  * [46] B. Datt.  On a Class of Solutions of the Gravitation Equations of Relativity.  General Relativity and Gravitation, 31:1619, October 1999.
  * [47] Charles W. Misner, K.S. Thorne, and J.A. Wheeler.  Gravitation.  W. H. Freeman, San Francisco, 1973.
  * [48] G. Russakoff.  A Derivation of the Macroscopic Maxwell Equations.  American Journal of Physics, 38(10):1188–1195, October 1970.
  * [49] H.A. Lorentz.  The Theory of Electrons and Its Applications to the Phenomena of Light and Radiant Heat.  Dover books on physics. Dover Publications, 2003.
  * [50] John David Jackson.  Classical Electrodynamics.  Wiley, 1998.
  * [51] Hubert Goenner and John Stachel.  Einstein Tensor and 3-Parameter Groups of Isometries with 2-Dimensional Orbits.  Journal of Mathematical Physics, 11(12):3358–3370, December 1970\.
  * [52] George David Birkhoff and Rudolph Ernest Langer.  Relativity and modern physics.  Harvard University Press, 1923.
  * [53] M. Cahen and R. Debever.  Sur le théorème de Birkhoff.  C. R. Acad. Sci. (Paris), 260:815, 1965.
  * [54] A. Barnes.  On Birkhoff’s theorem in general relativity.  Communications in Mathematical Physics, 33(1):75–82, March 1973\.
  * [55] Hubert Goenner.  Einstein tensor and generalizations of Birkhoff’s theorem.  Communications in Mathematical Physics, 16(1):34–47, March 1970\.
  * [56] C. Bona.  A new proof of the generalized Birkhoff theorem.  Journal of Mathematical Physics, 29(6):1440–1441, June 1988.
  * [57] Arthur Stanley Eddington.  The Mathematical Theory of Relativity.  Cambridge University Press, 1930.
  * [58] A. S. Eddington.  A Comparison of Whitehead’s and Einstein’s Formulæ.  Nature, 113(2832):192, February 1924.
  * [59] David Finkelstein.  Past-Future Asymmetry of the Gravitational Field of a Point Particle.  Physical Review, 110(4):965–967, May 1958.
  * [60] Roger Penrose.  Gravitational Collapse and Space-Time Singularities.  Phys. Rev. Lett., 14(3):57–59, January 1965.
  * [61] M. D. Kruskal.  Maximal Extension of Schwarzschild Metric.  Physical Review, 119(5):1743–1745, September 1960.
  * [62] G. Szekeres.  On the singularities of a Riemannian manifold.  Publ. Math. Debrecen, 7:285–301, 1960.
  * [63] I. D. Novikov.  R- and T-regions in a spacetime with a spherically symmetric space (in Russian).  Comm. State Sternberg Astron. Inst., 132:3, 1964.
  * [64] Werner Israel.  New interpretation of the extended schwarzschild manifold.  Phys. Rev., 143:1016–1021, Mar 1966.
  * [65] R. Kantowski and R. K. Sachs.  Some spatially homogeneous anisotropic relativistic cosmological models.  Journal of Mathematical Physics, 7(3):443–446, 1966.
  * [66] Andrzej Krasinski.  Inhomogeneous cosmological models.  Cambridge Univ. Press, Cambridge, UK, 3 2011.
  * [67] A. Molina and J.M.M. Senovilla, editors.  Inhomogeneous cosmological models. Proceedings, Spanish Relativity Meeting, Mao, ERE 94, Spain, September 12-14, 1994, Singapore, 1995\. World Scientific.
  * [68] Kari Enqvist.  Lemaitre–tolman–bondi model and accelerating expansion.  General Relativity and Gravitation, 40(2-3):451–466, Dec 2007\.
  * [69] Kari Enqvist and Teppo Mattsson.  The effect of inhomogeneous expansion on the supernova observations.  Journal of Cosmology and Astroparticle Physics, 2007(02):019–019, Feb 2007.
  * [70] Juan Garcia-Bellido and Troels Haugbølle.  Confronting Lemaitre Tolman Bondi models with observational cosmology.  Journal of Cosmology and Astroparticle Physics, 2008(4):003, April 2008.
  * [71] Håvard Alnes and Morad Amarzguioui.  Cmb anisotropies seen by an off-center observer in a spherically symmetric inhomogeneous universe.  Phys. Rev. D, 74:103520, Nov 2006.
  * [72] J. R. Oppenheimer and G. M. Volkoff.  On massive neutron cores.  Phys. Rev., 55:374–381, Feb 1939.
  * [73] Eduard Herlt.  Spherically Symmetric Nonstatic Perfect Fluid Solutions with Shear.  General Relativity and Gravitation, 28(8):919–934, August 1996\.
  * [74] Eric Poisson.  A Relativist’s Toolkit: The Mathematics of Black-Hole Mechanics.  Cambridge University Press, 12 2009.