∎

11institutetext: Institute of Theoretical Physics, Faculty of Mathematics and
Physics, Charles University in Prague, V Holešovičkách 2, 180 00 Prague 8,
Czech Republic
11email<EMAIL_ADDRESS><EMAIL_ADDRESS>
# Modelling Inhomogeneity in Szekeres Spacetime

David Vrba  Otakar Svítek

(Received: date / Accepted: date)

###### Abstract

We study the behaviour of the density contrast in quasi-spherical Szekeres
spacetime and derive its analytical behaviour as a function of $t$ and $r$. We
set up the inhomogeneity using initial data in the form of one extreme value
of the density and the radial profile. We derive conditions for density
extremes that are necessary for avoiding the shell crossing singularity and
show that in the special case of a trivial curvature function, the conditions
are preserved by evolution. We also show that in this special case if the
initial inhomogeneity is small, the time evolution does not influence the
density contrast, however its magnitude homogeneously decreases.

###### Keywords:

inhomogeneous cosmology Szekeres spacetime

###### pacs:

04.20.Jb, 04.20.Dw

##  1 Introduction

Homogeneous cosmological models have successfully explained many important
features of our universe. However, we know that the distribution of matter is
not homogeneous and so these models are only an approximation. In the last
decade modelling of inhomogeneity in cosmology has become a topic of
substantial interest. There are several different approaches to the problem.
Due to the nonlinear nature of Einstein equations one should not rely on
perturbation theory completely therefore exact models with nonuniform
distribution of matter should be considered as well. Among the most studied
exact inhomogeneous solutions of Einstein equations belong Lemaitre-Tolman-
Bondi (LTB) metric 1 , Szekeres solution 2 , Szafron family of solutions 3 ,
Stephani solution 4  or Lemaitre metric 1 , which is a generalization of LTB
for a fluid with nonzero pressure. An overview of inhomogeneous cosmological
solutions can be found in 6 .

LTB metric was extensively studied by Krasinski and Hellaby to model structure
formation 7  \- 9 . Together with Bolejko they used the solution to describe
formation of voids in the universe 10 . LTB metric is a special case of the
Szekeres solution that was discovered by Szekeres 2  and was developed by
Hellaby and Krasinski in a series of papers 11  \- 14 . They gave a
geometrical interpretation of the metric functions in all three different
types of geometry that are quasi-spherical, quasi-pseudospherical and quasi-
planar. The quasi-spherical case is currently the best understood of all three
types. It has found a cosmological application in the study of Bolejko who
focused on structure formation 15 , 16  and constructed models of a void with
an adjourning supercluster. Bolejko also investigated Buchert averaging of the
quasispherical Szekeres metric 17 , constructed a Szekeres-Swiss-cheese model
that he used to estimate the impact of inhomogeneity on the propagation of
light 18  and CMB observations 19 .

Walters and Hellaby showed in a recent paper 20  how to model inhomogeneity in
quasispherical Szekeres metric using initial and final data. They constructed
three models where they specified initial and final radial density profile and
one extreme value of the density on the final time slice in terms of a
deviation function along with the position of the extreme. In this way they
set up all 5 degrees of freedom in the Szekeres metric and worked out an
algorithm to calculate all other metric functions in the model.

For practical purposes it would be desirable to consider the evolution of
inhomogeneities in Szekeres spacetime determined only by initial data. In our
approach we specify the curvature function $f$ and the radial density profile
with one extreme value at the initial time. Thus we set up only 3 degrees of
freedom which is however sufficient to model the density contrast. To make an
appropriate choice for the initial density extreme we investigate the shell
crossing conditions in terms of the density extremes and show that the
resulting constraints are preserved throughout the evolution (for the assumed
value of $f$).

The paper is organized as follows. In the next section we give an overview of
the Szekeres metric and describe its basic properties. In section 3 we derive
conditions that the density extremes have to satisfy in order to avoid a shell
crossing singularity. In section 4 we show that if these conditions are met at
the initial time, then they are fulfilled at later times if the curvature
function vanishes. In section 5 we derive an analytical formula for the
density contrast defined as the difference of the extreme values and for the
case $f=0$ we show that if the initial radial inhomogeneity is small, then the
density contrast is proportional to ${t^{-3}}$ and the radial dependence of
the function does not evolve substantially in time. In section 6 we set up the
model by specifying two more functions and show the time evolution of the
density contrast in a specific example. We conclude in section 7.

##  2 Szekeres spacetime

The Szekeres spacetime is an exact dust solution of Einstein equations without
any symmetries. It was found by Szekeres 2  and generalized by Szafron 3  for
an energy momentum tensor describing a perfect fluid. The LT type of Szekeres
metric can be written as 21

${\rm d}s^{2}=-{\rm
d}t^{2}+\frac{\left(R^{\prime}-R\frac{E^{\prime}}{E}\right)^{2}}{\epsilon+f\left(r\right)}{\rm
d}r^{2}+\frac{R^{2}}{E^{2}}\left({\rm d}p^{2}+{\rm d}q^{2}\right),$ (1)

where

$E\left(r,p,q\right)\equiv\frac{S}{2}\left[\left(\frac{p-P}{S}\right)^{2}+\left(\frac{q-Q}{S}\right)^{2}+\epsilon\right],$
(2)
$E^{\prime}=\frac{S^{\prime}}{2}\left[1-\frac{\left(p-P\right)^{2}}{S^{2}}-\frac{\left(q-Q\right)^{2}}{S^{2}}\right]-\frac{P^{\prime}}{S}\left(p-P\right)-\frac{Q^{\prime}}{S}\left(q-Q\right)$
(3)

and $f$, $P$, $Q$, $S$ are arbitrary functions of $r$. The parameter
$\epsilon$ can have only three values $-1,0,1$ and it determines the geometry
of the two-spaces of constant $t$ and $r$. From Einstein equations it follows
a dynamical equation for the function $R$

${\dot{R}}^{2}=\frac{2M\left(r\right)}{R}+f\left(r\right)+\frac{\Lambda
R^{2}}{3}$ (4)

and an equation for the density evolution

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

where prime denotes a derivative with respect to $r$, $M$ is another arbitrary
function and $\Lambda$ is cosmological constant. So we have two Einstein
equations (4), (5) and by solving the first one we obtain one more arbitrary
function $t_{B}\left(r\right)$, which enters the solution in the form
$t-t_{B}$. From now on we will use a redefined form of the density
$\bar{\rho}\equiv\kappa c^{2}\rho$ and for the sake of simplicity we will drop
the bar.

As mentioned above, the metric does not have any symmetries, in other words it
has no Killing vectors. Nevertheless the three spaces of constant $t$ are
conformally flat 22 .

The parametrization of the metric involves 6 arbitrary functions however the
number of physical degrees of freedom is 5, because we can still rescale the
radial coordinate (the metric is invariant with respect to the transformation
$r^{\prime}=g(r)$).

###  2.1 Interpretation of $f$ and $t_{B}$

The sign of the function $f$ effects the solution of the Einstein equation
(4). It is a dynamical equation for $R$ and it looks very similar to the
Friedman equation, except that here the functions $f$ and $M$ depend on $r$.
We will assume only the case when $\Lambda=0$. There are three different types
of evolution 23 . If $f<0$ the evolution is elliptic which means that the
universe is first expanding and at some point the expansion stops and the
universe is collapsing to a final singularity. The solution is given in a
parametric form

$R=\frac{M}{\left(-f\right)}\left(1-\cos\eta\right),\quad\eta-\sin\eta=\frac{\left(-f\right)^{\frac{3}{2}}\left(t-t_{B}\right)}{M}.$
(6)

For $f>0$ the evolution is hyperbolic in the sense that the sign of the
expansion does not change, the universe either expands or collapses depending
on initial conditions and the parametric solution to the equation (4) is

$R=\frac{M}{\left(-f\right)}\left(\cosh\eta-1\right),\quad\sinh\eta-\eta=\frac{f^{\frac{3}{2}}\left(t-t_{B}\right)}{M}.$
(7)

For $f=0$ the evolution is parabolic and is given by equation

$R=\left(\frac{9}{2}M\right)^{\frac{1}{3}}\left(t-t_{B}\right)^{\frac{2}{3}}.$
(8)

Since $f$ is a function of $r$, the universe may have different evolution in
different regions. The parabolic evolution can be on the boundary between two
regions, one having elliptic evolution and the other one hyperbolic. So the
sign of $f$ determines the sign of scalar curvature of the three-spaces of
constant $t$ and when $f=0=f^{\prime}$, the three-spaces are flat which
corresponds to a presence of pure decaying modes as was shown in sussman .

The function $t_{B}$ is called the bang time function, because $t=t_{B}$ is
the moment when big bang happened. So unlike in homogeneous models, the
initial moment of evolution is position dependent.

The interpretation of the other metric functions ($R$, $M$, $P$, $Q$, $S$ )
depends on $\epsilon$ and from now on we will only consider the case
$\epsilon=+1$ which is often called the quasi-spherical case.

###  2.2 Coordinate transformation and interpretation of $R$

We can make a coordinate transformation 11

$\frac{p-P}{S}=\cot\frac{\theta}{2}\cos\phi,\quad\frac{q-Q}{S}=\cot\frac{\theta}{2}\sin\phi$
(9)

and rewrite the induced two-metric of the surfaces of constant $t$ and $r$
into coordinates $\theta$ and $\phi$. After applying the transformation we get

${\rm d}s^{2}=\frac{4}{S^{2}}\sin^{4}\frac{\theta}{2}\left({\rm d}p^{2}+{\rm
d}q^{2}\right)={\rm d}\theta^{2}+\sin^{2}\theta{\rm d}\phi^{2},$ (10)

which is a metric on a unit sphere. The transformation (9) is nothing but a
stereographic projection and the function $E$ describes how the $(p,q)$ plane
is mapped onto a unit sphere. Every sphere is multiplied by $R^{2}$ and the
function is sometimes called the areal radius, because it is actually the
radius of a sphere on a comoving coordinate $r$.

###  2.3 Geometrical meaning of $E$

In order to understand more about the geometrical properties of the metric (1)
it is important to investigate the function $E$. We can see that this function
appears in the metric and in the equation for the density (5) in the form
$\frac{E^{\prime}}{E}$. Particularly, we can investigate when the function is
equal to zero and what are its extreme values. Using the transformation (9) in
(2) and (3) we can write $\frac{E^{\prime}}{E}$ as

$\frac{E^{\prime}}{E}=-\frac{S^{\prime}\cos\theta+\sin\theta\left(P^{\prime}\cos\phi+Q^{\prime}\sin\phi\right)}{S}.$
(11)

The equation $\frac{E^{\prime}}{E}=0$ now becomes

$S^{\prime}\cos\theta+P^{\prime}\sin\theta\cos\phi+Q^{\prime}\sin\theta\sin\phi=0$
(12)

and after realizing that $\cos\theta=z$, $\sin\theta\cos\phi=y$ and
$\sin\theta\sin\phi=x$, the equation (12) becomes

$S^{\prime}z+P^{\prime}x+Q^{\prime}y=0,$ (13)

which is an equation of a plane that goes through the origin of the spherical
coordinate system and intersects the $r=t=const$ sphere in a great circle. Now
calculating the derivatives of (11) with respect to $\theta$ and $\phi$ and
putting it equal to zero we can find that there are two extreme values and
they are located at opposite sites on the sphere, in other words if one
extreme is at the coordinates $\left(\theta_{1},\phi_{1}\right)$, the other
one is at $\left(\pi-\theta_{1},\phi_{1}+\pi\right)$ and they are located
symmetrically with respect to the plane $\frac{E^{\prime}}{E}=0$, because the
coordinates of the extremes are exactly the components of the unit normal to
the plane $\frac{E^{\prime}}{E}=0$. The values of the extremes are

$\displaystyle\left(\frac{E^{\prime}}{E}\right)_{max}$ $\displaystyle=$
$\displaystyle\frac{\sqrt{\left(S^{\prime}\right)^{2}+\left(P^{\prime}\right)^{2}+\left(Q^{\prime}\right)^{2}}}{S},$
$\displaystyle\left(\frac{E^{\prime}}{E}\right)_{min}$ $\displaystyle=$
$\displaystyle-\frac{\sqrt{\left(S^{\prime}\right)^{2}+\left(P^{\prime}\right)^{2}+\left(Q^{\prime}\right)^{2}}}{S},$
(14)

where max and min refers to maximum and minimum respectively. So the extremes
have opposite values and the function $\frac{E^{\prime}}{E}$ behaves on the
sphere like a dipole 11 .

##  3 Shell crossing conditions in terms of the density extremes

We are now going to study the formula for the density (5). We can see that
under certain circumstances the density can change sign or possibly diverge if
the denominator becomes zero. Those points where this happens are called the
shell crossing singularity. The conditions that the metric functions $M$, $f$,
$t_{B}$ and $R$ have to satisfy in order to avoid shell crossing singularity
can be found in 11 . In this section we will derive the conditions that the
extreme values of the density $\rho_{max}$ and $\rho_{min}$ have to satisfy in
order to avoid the shell crossing singularity.

We can rewrite (5) as

$\rho\left(t,r,\theta,\varphi\right)=2\frac{M^{\prime}-3M\frac{E^{\prime}}{E}}{R^{2}\left(R^{\prime}-R\frac{E^{\prime}}{E}\right)}=\frac{R^{\prime}\rho_{LT}-R\frac{E^{\prime}}{E}\rho_{AV}}{R^{\prime}-R\frac{E^{\prime}}{E}},$
(15)

where we define

$\rho_{LT}\left(t,r\right)\equiv\frac{2M^{\prime}}{R^{2}R^{\prime}}$ (16)

and

$\rho_{AV}\left(t,r\right)\equiv\frac{6M}{R^{3}}.$ (17)

Here $\rho_{LT}$ is the density that we get from (5) if we set
$\frac{E^{\prime}}{E}=0$. It is also a radial density in the sense that on a
given $r$ it is the value of the density around the great circle that lies in
the plane that defines the dipole on the sphere. The index $LT$ refers to
Lemaitre-Tolman-Bondi metric, because its density is given with exactly the
same formula as (16). If we make the choice

$R\left(t_{i}\right)\equiv r$ (18)

and use the fact that the function $M$ depends only on $r$, we can express it
on the initial time slice when $\frac{E^{\prime}}{E}=0$ as

$M=\frac{1}{2}\int_{0}^{r}\rho_{LT0}r^{\prime 2}{\rm d}r^{\prime},$ (19)

where the index $0$ refers to the value of the function at the initial time
$t_{i}$. In that case the definitions $(\ref{35})$ and $(\ref{36})$ can be
further rewritten as

$\rho_{LT}\left(t,r\right)=\frac{1}{R^{2}R^{\prime}}\rho_{LT0}r^{2},\quad\quad\rho_{AV}\left(t,r\right)=\frac{3}{R^{3}}\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}.$ (20)

The initial condition (18) and formulas (19) - (20) will however not be used
in the following calculations (with the exception of equation (44) giving
critical point position) and the main results of this section are independent
of it, thus giving us generally applicable criteria.

It is good to take a look at the derivative of the density with respect to
$\frac{E^{\prime}}{E}$

$\frac{\partial\rho}{\partial\frac{E^{\prime}}{E}}=RR^{\prime}\frac{\rho_{LT}-\rho_{AV}}{\left(R^{\prime}-R\frac{E^{\prime}}{E}\right)^{2}}.$
(21)

From the derivative we can see, that if
$R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)>0$ the derivative is positive and
the density is growing as $\frac{E^{\prime}}{E}$ increases. On the other hand
if $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)<0$ the density will decrease as
$\frac{E^{\prime}}{E}$ increases. Either way we can see that the density
behaves on each sphere also like a dipole in the sense that it has two extreme
values, maximum and minimum and they are located at the same position as
extreme values of the function $\frac{E^{\prime}}{E}$. In the case when
$R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)<0$ we can write for the density
maximum

$\rho_{max}=\frac{R^{\prime}\rho_{LT}-R\left(\frac{E^{\prime}}{E}\right)_{min}\rho_{AV}}{R^{\prime}-R\left(\frac{E^{\prime}}{E}\right)_{min}}$
(22)

and for the minimum

$\rho_{min}=\frac{R^{\prime}\rho_{LT}-R\left(\frac{E^{\prime}}{E}\right)_{max}\rho_{AV}}{R^{\prime}-R\left(\frac{E^{\prime}}{E}\right)_{max}}.$
(23)

In the case when $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)>0$ the role of
$\left(\frac{E^{\prime}}{E}\right)_{max}$ interchanges with
$\left(\frac{E^{\prime}}{E}\right)_{min}$ in the last two equations.

We split the investigation of the density into three parts, in the first part
we investigate when both the numerator and the denominator are positive, in
the second case we check when they are both negative and in the last case we
take a look at the special case when they are both zero. It will be convenient
to split the first part into three more subcases depending on the value of
$R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)$. For the investigation we will
also need a formula for the function $\frac{E^{\prime}}{E}$ given in terms of
$\rho_{LT}$ and $\rho_{AV}$ which we can derive from (15)

$\frac{E^{\prime}}{E}=\frac{R^{\prime}}{R}\frac{\rho-\rho_{LT}}{\rho-\rho_{AV}}.$
(24)

First case: numerator and denominator of the density are both positive.
Subcase A: $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)<0$

From the positivity of the denominator in (15) we get

$\frac{R^{\prime}}{R}>\left(\frac{E^{\prime}}{E}\right)_{max}=\frac{R^{\prime}}{R}\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}},$
(25)

because as follows from (21) if $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)<0$
we can see that $\frac{E^{\prime}}{E}$ has maximum where the density has
minimum. From the dipole property

$\left(\frac{E^{\prime}}{E}\right)_{max}=-\left(\frac{E^{\prime}}{E}\right)_{min}$
(26)

we know that $\left(\frac{E^{\prime}}{E}\right)_{max}$ is non-negative and
(25) implies that $R^{\prime}>0$ and therefore we have to consider
$\rho_{LT}<\rho_{AV}$. We can now set conditions for $\rho_{min}$ so that the
right hand side in (25) is nonnegative. Inequality (25) simplifies to

$1>\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}}$ (27)

It is reasonable to require $\rho_{min}\leq\rho_{LT}$, because the minimum
should be the smallest value on each sphere ($\rho_{min}=\rho_{LT}$
corresponds to the case when the density is homogeneously distributed on the
whole sphere). Since we consider $\rho_{LT}<\rho_{AV}$ it follows that
$\rho_{min}<\rho_{AV}$ and the right hand side of (25) is nonnegative. It
means that besides the obvious condition $\rho_{min}\leq\rho_{LT}$ we don’t
have any other restriction for the minimum.

Due to the dipole property we can try to derive a condition for the maximum

$\frac{R^{\prime}}{R}>\left(\frac{E^{\prime}}{E}\right)_{max}=-\frac{R^{\prime}}{R}\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}}.$
(28)

Since $\rho_{max}$ should be greater than or equal to $\rho_{LT}$ it implies
that $\rho_{max}<\rho_{AV}$ in order for the right hand side in (28) to be
non-negative. This condition is reasonable since $\rho_{AV}$ is not true value
of the density as follows from its definition (17). Taking this into
consideration we can solve the inequality (28) and the solution is

$\rho_{max}<\frac{1}{2}\left(\rho_{LT}+\rho_{AV}\right)$ (29)

and this condition is clearly more restrictive than $\rho_{max}<\rho_{AV}$.

For the numerator of (15) we need

$\frac{R^{\prime}}{R}\rho_{LT}>\frac{E^{\prime}}{E}\rho_{AV}$ (30)

which can be rewritten as

$\frac{R^{\prime}}{R}\frac{\rho_{LT}}{\rho_{AV}}>\left(\frac{E^{\prime}}{E}\right)_{max}=\frac{R^{\prime}}{R}\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}}.$
(31)

It again follows that $R^{\prime}>0$ and $\rho_{LT}<\rho_{AV}$. The right-hand
side of (31) has to be non-negative and we know that $\rho_{min}\leq\rho_{LT}$
and since $\rho_{LT}<\rho_{AV}$ we also have $\rho_{min}<\rho_{AV}$ so the
inequality (31) is always true and we don’t get any new condition for
$\rho_{min}$. We now use the property (26) and look for a condition for
$\rho_{max}$

$\frac{\rho_{LT}}{\rho_{AV}}>-\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}}.$
(32)

We know that $\rho_{max}\geq\rho_{LT}$ so to make the right hand side non-
negative we need $\rho_{max}<\rho_{AV}$. The solution to the inequality (32)
is

$\rho_{max}<2\frac{\rho_{AV}\rho_{LT}}{\rho_{LT}+\rho_{AV}}.$ (33)

The conditions for the denominator and the numerator have to be satisfied
simultaneously, so for the maximum we have (29) and (33) and (as can be seen
from figure 1) the condition (33) is stronger in the case of
$\rho_{LT}<\rho_{AV}$ and has to be fulfilled in order for the density to be
positive. As far as the minimum is concerned, apart from the obvious condition
$\rho_{min}\leq\rho_{LT}$ it is not constrained at all.

First case: numerator and denominator of the density are both positive.
Subcase B: $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)>0$

The subcase B will be investigated similarly. The only difference is that from
(21) we know that if $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)<0$ the
function $\frac{E^{\prime}}{E}$ has maximum where the density has maximum and
this modifies (22) and (23) accordingly. From the denominator of (15) we have

$\frac{R^{\prime}}{R}>\left(\frac{E^{\prime}}{E}\right)_{max}=\frac{R^{\prime}}{R}\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}},$
(34)

which again implies $R^{\prime}>0$ and $\rho_{LT}>\rho_{AV}$ because
$\left(\frac{E^{\prime}}{E}\right)_{max}\geq 0$ and it simplifies to

$1>\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}}.$ (35)

The analysis is basically the same, we again need the right-hand side to be
non-negative. We know that $\rho_{max}\geq\rho_{LT}$ and therefore
$\rho_{max}>\rho_{AV}$ as well, so the inequality (35) will be always true as
long as $\rho_{LT}>\rho_{AV}$ and we don’t get any further condition for the
maximum. For the minimum we have

$\frac{R^{\prime}}{R}>\left(\frac{E^{\prime}}{E}\right)_{max}=-\frac{R^{\prime}}{R}\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}}.$
(36)

In order for the right-hand side to be non-negative we need
$\rho_{min}>\rho_{AV}$ which is a valid condition for the density minimum,
because as mentioned above $\rho_{AV}$ is not true value of the density. The
inequality (36) has the solution

$\rho_{min}>\frac{1}{2}\left(\rho_{AV}+\rho_{LT}\right).$ (37)

From the numerator of (15) we have the condition

$\frac{R^{\prime}}{R}\frac{\rho_{LT}}{\rho_{AV}}>\left(\frac{E^{\prime}}{E}\right)_{max}=\frac{R^{\prime}}{R}\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}}.$
(38)

It again implies $R^{\prime}>0$ and $\rho_{LT}>\rho_{AV}$ and it does not give
us any constraint for the maximum, because we know that
$\rho_{max}\geq\rho_{LT}$ and therefore $\rho_{max}>\rho_{AV}$. So (38) holds
as long as $R^{\prime}>0$ and $\rho_{LT}>\rho_{AV}$. For the minimum we can
write

$\frac{R^{\prime}}{R}\frac{\rho_{LT}}{\rho_{AV}}>\left(\frac{E^{\prime}}{E}\right)_{max}=-\frac{R^{\prime}}{R}\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}}.$
(39)

The right-hand side will be positive if $\rho_{min}>\rho_{AV}$ and the
solution to (39) is

$\rho_{min}>2\frac{\rho_{LT}\rho_{AV}}{\rho_{LT}+\rho_{AV}}.$ (40)

All together in the case of $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)>0$ we
have two conditions (37) and (40) for the minimum and as can be seen from
figure 1 the first one is stronger. For the maximum, except for the obvious
$\rho_{max}\geq\rho_{LT}$, we don’t have any constraint so the density will be
positive as long as (37) holds and $\rho_{max}\geq\rho_{LT}$.

Figure 1: In this figure we can see 4 different surfaces that constrain the
density extremes. In the right part of the figure they are from the bottom to
the top: $\rho_{AV}$, $2\frac{\rho_{AV}\rho_{LT}}{\rho_{AV}+\rho_{LT}}$,
$\frac{1}{2}\left(\rho_{AV}+\rho_{LT}\right)$, $\rho_{LT}$.

First case: numerator and denominator of the density are both positive.
Subcase C: $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)=0$

Before we consider this special case we can first rewrite (15) as

$\rho=\frac{R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)+\rho_{AV}\left(R^{\prime}-R\frac{E^{\prime}}{E}\right)}{R^{\prime}-R\frac{E^{\prime}}{E}}$
(41)

Now using $R^{\prime}\left(\rho_{LT}-\rho_{AV}\right)=0$ we can see that the
formula for the density simplifies

$\rho=\rho_{AV}\frac{R^{\prime}-R\frac{E^{\prime}}{E}}{R^{\prime}-R\frac{E^{\prime}}{E}}=\rho_{AV},$
(42)

in other words the value of the density is independent of $\theta$ and $\phi$
and the density is homogeneously distributed on the whole sphere as follows
also from (21) because the derivative of the density is now zero. Using the
same argument as in subcase A and B we have $R^{\prime}>0$ and therefore
$\rho_{LT}=\rho_{AV}$. For the extreme values of the density we obviously need
in this case

$\rho_{max}=\rho_{min}=\rho_{AV}=\rho_{LT}.$ (43)

The point where $\rho_{LT}=\rho_{AV}$ is especially interesting because the
shell-crossing conditions for $\rho_{max}$ and $\rho_{min}$ change here and
the location of $\rho_{max}$ and $\rho_{min}$ interchanges by passing through
the plane $\frac{E^{\prime}}{E}=0$. The position of this critical point
$r_{c}$ generally evolves and we can get its value as the solution to the
equation

$3\frac{R^{\prime}}{R}\int_{0}^{r_{c}}\rho_{LT0}r^{2}{\rm
d}r=\rho_{LT0}r_{c}^{2},$ (44)

where we used the initial condition (18) and formulas (20) in order to express
$\rho_{LT}$ and $\rho_{AV}$. So the time evolution of the critical point will
depend on the time evolution of the function $\frac{R^{\prime}}{R}$.

Second case: numerator and denominator of the density are both negative.

In order for the denominator of the density (15) to be negative, we need

$\frac{R^{\prime}}{R}<\left(\frac{E^{\prime}}{E}\right)_{min}.$ (45)

This condition will not be fulfilled unless $R^{\prime}<0$ because
$\left(\frac{E^{\prime}}{E}\right)_{min}\leq 0$. The investigation is than
similar as in the first case and leads to the same conditions for $\rho_{min}$
and $\rho_{max}$ as can be easily verified.

On the other hand if $R^{\prime}>0$ and we allow the denominator to be
negative on a given $r$, it will be also positive on that $r$ for some
specific $\theta$ and $\phi$. This behaviour may or may not be correct,
depending on the behaviour of the numerator and specifically depending on
whether or not the numerator changes sign at the same $\theta$, $\phi$. For
the numerator we have inequality

$R^{\prime}\rho_{LT}-R\frac{E^{\prime}}{E}\rho_{AV}<0.$ (46)

This can be rewritten as

$\frac{R^{\prime}}{R}\frac{\rho_{LT}}{\rho_{AV}}<\left(\frac{E^{\prime}}{E}\right)_{min}.$
(47)

The same argument as in the case of denominator tells us that this will be
satisfied only if $R^{\prime}<0$. So if the numerator is negative on a given
$r$ but $R^{\prime}>0$, there is a region on this $r$ where it is positive
too. For the positiveness of the denominator in a region $(\theta,\phi)$ we
have

$\frac{R^{\prime}}{R}>\frac{E^{\prime}}{E}$ (48)

and for the numerator the condition is

$\frac{R^{\prime}}{R}\frac{\rho_{AV}}{\rho_{LT}}>\frac{E^{\prime}}{E}.$ (49)

Clearly the conditions are not the same unless $\rho_{LT}=\rho_{AV}$, so
except for this special case, there will always be a region on a given $r$
where the numerator and the denominator will have different signs. The places
where the sign changes are in both cases circles on the sphere, that are
parallel to the great circle that defines the plane of the dipole. The order
in which the circles go are seen in figure 2, where the parallel circles are
mapped as parallel lines.

So if in the second case $R^{\prime}>0$ there are no conditions for the
density extremes that would prevent shell crossing.

(a)  (b)

Figure 2: The behaviour of the density for $\rho_{AV}>\rho_{LT}$ in panel (a)
and $\rho_{AV}<\rho_{LT}$ in panel (b) when the conditions for no shell-
crossing are not met. The circle represents a sphere on a given $r$, the solid
line represents the plane $\frac{E^{\prime}}{E}=0$. The dashed line represents
a plane
$\frac{E^{\prime}}{E}\frac{\rho_{LT}}{\rho_{AV}}=\frac{R^{\prime}}{R}$, that
is parallel to the plane $\frac{E^{\prime}}{E}=0$, the density becomes zero
here. The dotted line is the plane $\frac{E^{\prime}}{E}=\frac{R^{\prime}}{R}$
and it is the plane where the density diverges. The density is negative
between the dashed and dotted line.

Third case: numerator and denominator of the density are both zero.

The denominator is zero when

$\frac{R^{\prime}}{R}=\frac{E^{\prime}}{E},$ (50)

but this equation will be satisfied on a given $r$ for all $\theta$ and $\phi$
only if both $\frac{R^{\prime}}{R}$ and $\frac{E^{\prime}}{E}$ are zero. This
consequently means that
$E^{\prime}=P^{\prime}=Q^{\prime}=S^{\prime}=R^{\prime}=0$. The same arguments
are valid for the numerator of the density. As was discussed in 11  and in the
context of LTB metric in 23  this point will not be a singularity if also
$M^{\prime}=0$. Since $M^{\prime}$ and $E^{\prime}$ are independent of time it
implies that $R^{\prime}$ cannot evolve and we need stationary configuration
in order to avoid shell crossing.

Conditions for $\rho_{min}$ | Conditions for $\rho_{max}$
---|---
$\rho_{AV}>\rho_{LT}$
$\rho_{min}\leq\rho_{LT}<\rho_{AV}$ | $\rho_{max}\in\langle\left.\rho_{LT};2\frac{\rho_{AV}\rho_{LT}}{\rho_{LT}+\rho_{AV}}\right)$
$\rho_{AV}<\rho_{LT}$
$\rho_{min}\in\left(\frac{\rho_{AV}+\rho_{LT}}{2};\rho_{LT}\right.\rangle$ | $\rho_{max}\geq\rho_{LT}>\rho_{AV}$
$\rho_{AV}=\rho_{LT}$
$\rho_{min}=\rho_{LT}$ | $\rho_{max}=\rho_{LT}$
Table 1: The list of shell crossing conditions for density extremes in terms
of $\rho_{LT}$ and $\rho_{AV}$.

The table 1 summarizes the conditions for the extreme values of the density
that we obtained so that shell crossing would be avoided. From the dipole
property (26) it follows that there is a constraint between the density
extremes

$\frac{\rho_{min}-\rho_{LT}}{\rho_{min}-\rho_{AV}}=-\frac{\rho_{max}-\rho_{LT}}{\rho_{max}-\rho_{AV}},$
(51)

which means that knowing one extreme value allows us calculate the other one
with a formula

$\rho_{max}=\frac{\rho_{min}\left(\rho_{AV}+\rho_{LT}\right)-2\rho_{LT}\rho_{AV}}{2\rho_{min}-\rho_{LT}-\rho_{AV}}$
(52)

that follows from (51). And for the other extreme we just interchange min for
max in the last equation. Consequently it is sufficient to ensure that just
one of the extremes satisfies the conditions in table 1 and the constraint
(51) ensures that the other density extreme has a value that does not break
the conditions in table 1.

##  4 Time dependence of the shell crossing conditions, case $f=0$

In the previous section we derived conditions that $\rho_{max}$ and
$\rho_{min}$ have to satisfy in order to avoid shell crossing singularity. To
find out if a shell crossing occurs at a given time, we need to calculate
$\rho_{max}$ or $\rho_{min}$ at that time and then check if the conditions are
met. It would be useful to have conditions in terms of the initial
$\rho_{LT0}$ and $\rho_{min0}$ or $\rho_{max0}$ that would ensure no shell
crossing at any time during the evolution. In this section we show that in the
special case when $f=0$, if we avoid shell crossing at the initial time, it is
guaranteed that no shell crossing occurs during the time evolution.

First we show that if on the initial time slice $\rho_{AV0}>\rho_{LT0}$, then
this condition holds at any later time. We will use the initial condition (18)
and start with the inequality that we want to prove $\rho_{AV}>\rho_{LT}$

$\frac{3}{R^{3}}\int_{0}^{r}\rho_{LT0}r^{\prime 2}{\rm
d}r^{\prime}>\frac{1}{R^{2}R^{\prime}}\rho_{LT0}r^{2},$ (53)

this can be rewritten as

$\frac{R^{\prime}}{R}r\,\frac{3\int_{0}^{r}\rho_{LT0}r^{\prime 2}{\rm
d}r^{\prime}}{\rho_{LT0}r^{3}}>1$ (54)

and the inequality sign depends on the sign of $R^{\prime}$ and here we
assumed $R^{\prime}>0$. We can see that the time dependence in the last
inequality is hidden in the function $\frac{R^{\prime}}{R}r$. The choice $f=0$
has parabolic evolution and the solution to the equation (4) is given by (8).
If we substitute in (8) for $M$ from (19) we get

$R\left(t,r\right)=\left(\frac{9}{4}\right)^{\frac{1}{3}}\left(\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}\right)^{\frac{1}{3}}\left(t-t_{B}\right)^{\frac{2}{3}}.$
(55)

We can fix the bang time function $t_{B}$ using (18) to be

$t_{B}=\left[t_{i}-\frac{2}{3}\left(\frac{r^{3}}{\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}}\right)^{\frac{1}{2}}\right],$ (56)

where $t_{i}$ is the initial moment of the evolution. By calculating the
radial derivative of (55) we can express $\frac{R^{\prime}}{R}r$ as

$\frac{R^{\prime}}{R}r=\frac{\rho_{LT0}r^{3}}{3\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}}-\frac{2}{3}\frac{t_{B}^{\prime}}{t-t_{B}}r.$ (57)

After calculating radial derivative of the bang time function (56) and
substituting it in the last formula we obtain

$\frac{R^{\prime}}{R}r=\frac{\rho_{LT0}r^{3}}{3\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}}+\frac{2}{3}\frac{\frac{3\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm
d}r^{\prime}}{\rho_{LT0}r^{3}}-1}{\left(\frac{3\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm
d}r^{\prime}}{\rho_{LT0}r^{3}}\right)^{\frac{3}{2}}}\,\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}.$
(58)

Now we define

$A\equiv\frac{\rho_{AV0}}{\rho_{LT0}}=\frac{3\int_{0}^{r}\rho_{LT0}r^{\prime
2}{\rm d}r^{\prime}}{\rho_{LT0}r^{3}}$ (59)

and using this definition we can rewrite (58) as

$\frac{R^{\prime}}{R}r=\frac{1}{A}+\frac{2}{3}\frac{A-1}{A^{\frac{3}{2}}}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}.$
(60)

From this formula it can be shown that $\frac{R^{\prime}}{R}r$ is a monotonic
function of $t$. For $t=t_{i}$, its value is 1, which follows from our initial
condition $R\left(t_{i}\right)=r$ and for $t\to\infty$,
$\frac{R^{\prime}}{R}r\to\frac{1}{A}$ so it is increasing if
$\rho_{AV0}>\rho_{LT0}$ and it is decreasing if $\rho_{AV0}<\rho_{LT0}$. In
any case the function is positive which implies $R^{\prime}>0$ so our
assumption of the inequality sign in (54) is justified. Using the definition
(59) and substituting for $\frac{R^{\prime}}{R}r$ into (54) we get

$\frac{2}{3}\frac{A-1}{\sqrt{A}}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}>0.$
(61)

Since $t>t_{B}$ the last inequality will be satisfied as long as $A>1$, which
is equivalent to $\rho_{AV0}>\rho_{LT0}$ as can be seen from the definition of
$A$. So we can see exactly what we wanted to prove, if $\rho_{AV0}>\rho_{LT0}$
than $\rho_{AV}>\rho_{LT}$ at any time. Similarly we could prove that if
$\rho_{AV0}<\rho_{LT0}$ than $\rho_{AV}<\rho_{LT}$ at any time. Next we will
assume that $\rho_{AV0}>\rho_{LT0}$ and we will express the shell crossing
condition $\rho_{min}<\rho_{LT}$ in terms of initial data. For the function
$\frac{E^{\prime}}{E}$ we have formula (24) and since the function does not
depend on time, we can express it on the initial time slice

$\frac{E^{\prime}}{E}_{max}=\frac{1}{r}\frac{\rho_{min0}-\rho_{LT0}}{\rho_{min0}-\rho_{AV0}}.$
(62)

We plug the last formula into (23) and after some calculations we get

$\rho_{min}=\frac{r^{3}}{R^{3}}\rho_{LT0}\frac{\frac{\rho_{min0}}{\rho_{LT0}}\left(1-\frac{\rho_{AV0}}{\rho_{LT0}}\right)}{\frac{R^{\prime}}{R}r\left(\frac{\rho_{min0}}{\rho_{LT0}}-\frac{\rho_{AV0}}{\rho_{LT0}}\right)-\frac{\rho_{min0}}{\rho_{LT0}}+1}$
(63)

and we used
$\frac{\rho_{AV}}{\rho_{LT}}=\frac{R^{\prime}}{R}r\frac{\rho_{AV0}}{\rho_{LT0}}$
which follows from its definitions (16) and (17). We will now define

$C\equiv\frac{\rho_{min0}}{\rho_{LT0}}$ (64)

and using this definition with (59) we can rewrite (63) as

$\rho_{min}=\frac{r^{3}}{R^{3}}\rho_{LT0}\frac{C\left(1-A\right)}{\frac{R^{\prime}}{R}r\left(C-A\right)+\left(1-C\right)}.$
(65)

The inequality $\rho_{min}<\rho_{LT}$ can now be rewritten as

$\frac{R^{\prime}}{R}r\frac{\left(1-A\right)C}{\frac{R^{\prime}}{R}r\left(C-A\right)+\left(1-C\right)}<1.$
(66)

We need to multiply the last inequality by the denominator, but in order to do
that we need to find out its sign. From (59) and (64) it follows

$C-A=-\frac{1}{\rho_{LT0}}\left(\rho_{AV0}-\rho_{min0}\right)<0,$ (67)

next using the minimal value of $\frac{R^{\prime}}{R}r$ to be $\frac{1}{A}$ we
can write for the denominator of (66)

$\left(1-C\right)-\frac{R^{\prime}}{R}r\left(A-C\right)\leq
1-C-\frac{1}{A}\left(A-C\right)=C\left(\frac{1}{A}-1\right)<0,$ (68)

because we assume $\rho_{AV0}>\rho_{LT0}$, i.e $A>1$. So the denominator is
negative and when we multiply by it in (66), we get after some calculations

$\frac{R^{\prime}}{R}rA>1,$ (69)

which is inequality (54) written in terms of $A$ and we already proved that it
holds as long as $A>1$. So we can see that if $\rho_{AV0}>\rho_{LT0}$ than
$\rho_{min}<\rho_{LT}$ is fulfilled and the shell crossing will be avoided at
any time. Similarly it could be shown that in the case
$\rho_{AV0}<\rho_{LT0}$, i.e. $A<1$, if $\rho_{max0}>\rho_{LT0}$ than
$\rho_{max}>\rho_{LT}$ at any time.

We can also see that if $\rho_{AV0}=\rho_{LT0}$ at some point $r_{c}$, it
means that we have $A=1$ and from (60) it follows $\frac{R^{\prime}}{R}r=1$.
So the position of the critical point $r_{c}$ as defined in (44) is time
independent.

##  5 Density contrast, case $f=0$

We will choose $\rho_{LT0}$ in such a way that $\rho_{LT0}<\rho_{AV0}$. We
already expressed $\rho_{min}$ in terms of $\rho_{LT}$ and $\rho_{min0}$ in
(65). Similarly we can express $\rho_{max}$ in terms of $\rho_{LT0}$ and
$\rho_{max0}$ as

$\rho_{max}=\frac{r^{3}}{R^{3}}\rho_{LT0}\frac{D\left(1-A\right)}{\frac{R^{\prime}}{R}r\left(D-A\right)-\left(D-1\right)},$
(70)

where we define

$D\equiv\frac{\rho_{max0}}{\rho_{LT0}}.$ (71)

We are now interested in the time evolution of the difference of the extreme
values of the density, so we define

$\displaystyle\Delta\rho\equiv\rho_{max}-\rho_{min}=\frac{r^{3}}{R^{3}}\rho_{LT0}\left(1-A\right)\times$
$\displaystyle\left[\frac{D}{\frac{R^{\prime}}{R}r\left(D-A\right)-\left(D-1\right)}-\frac{C}{\frac{R^{\prime}}{R}r\left(C-A\right)-\left(C-1\right)}\right]\
$ (72)

We can now substitute for $R$ from (55) and after some more calculations it
can be rewritten in the form

$\Delta\rho\left(t,r\right)=\frac{4}{3}\frac{1}{\left(t-t_{B}\right)^{2}}\frac{\Delta\rho_{0}}{\rho_{LT0}}h\left(t,r\right)\,,$
(73)

where we have defined $\Delta\rho_{0}\equiv\rho_{max0}-\rho_{min0}$,

$h\left(t,r\right)\equiv\frac{1-A}{A}\
\frac{1-\frac{R^{\prime}}{R}rA}{\left(\frac{R^{\prime}}{R}r\right)^{2}\left(D-A\right)\left(C-A\right)+\left(D-1\right)\left(C-1\right)},$
(74)

and we used

$2\left(CD+A\right)-\left(C+D\right)\left(A+1\right)=0,$ (75)

which can be derived from the constraint (52). We will now investigate the
behaviour of the function $h$.

We set $A=1+\epsilon$ and take a look at how the function $h$ behaves in the
case when $\epsilon<<1$, which corresponds to the situation when the radial
derivative of $\rho_{LT0}$ is small. First we approximate
$\frac{R^{\prime}}{R}r$ that is given by (60),

$\displaystyle\frac{R^{\prime}}{R}r$ $\displaystyle=$
$\displaystyle\frac{1}{1+\epsilon}\left(1+\frac{2}{3}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}\frac{\epsilon}{\sqrt{1+\epsilon}}\right)\approx$
$\displaystyle\approx$ $\displaystyle
1+\left(\frac{2}{3}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}-1\right)\epsilon+\left(2-\frac{1}{3}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}\right)\epsilon^{2}+o\left(\epsilon^{3}\right).$

So in the numerator on the right hand side of (74) we get

$\left(1-A\right)(1-\frac{R^{\prime}}{R}rA)=-\left(1-\frac{R^{\prime}}{R}r\right)\epsilon+\frac{R^{\prime}}{R}r\epsilon^{2}\approx\frac{2}{3}\frac{1}{t-t_{B}}\sqrt{\frac{3}{\rho_{LT0}}}\epsilon^{2}+o\left(\epsilon^{3}\right).$
(77)

In order to simplify the denominator in (74) we set $D=1+\delta$ and $C=1-\xi$
and we assume that $\Delta\rho_{0}<<\rho_{LT0}$ in which case $\delta<<1$ and
$\xi<<1$. In this case the approximation of the denominator reads

$\displaystyle
A\left[\left(\frac{R^{\prime}}{R}r\right)^{2}\left(C-A\right)\left(D-A\right)+\left(C-1\right)\left(D-1\right)\right]\approx$
$\displaystyle\quad\approx\left(C-A\right)\left(D-A\right)+\left(C-1\right)\left(D-1\right)$
(78)

and for the function $h$ we get

$h\approx\frac{2}{\sqrt{3}}\frac{1}{t-t_{B}}\frac{1}{\sqrt{\rho_{LT0}}}\frac{\left(A-1\right)^{2}}{\left(C-A\right)\left(D-A\right)+\left(C-1\right)\left(D-1\right)}.$
(79)

From (75) it follows

$2CD=AC+AD-2A+C+D.$ (80)

Using (71) we can rewrite the denominator in (79)

$\left(C-A\right)\left(D-A\right)+\left(C-1\right)\left(D-1\right)=A^{2}-2A+1=\left(A-1\right)^{2}.$
(81)

So for $h$ we can write

$h\approx\frac{2}{\sqrt{3}}\frac{1}{t-t_{B}}\frac{1}{\sqrt{\rho_{LT0}}}$ (82)

and for $\Delta\rho$ we have

$\Delta\rho\approx\frac{8\sqrt{3}}{9}\frac{1}{\left(t-t_{B}\right)^{3}}\frac{\Delta\rho_{0}}{\rho_{LT0}^{\frac{3}{2}}}.$
(83)

From the last equation we can see that if the radial derivative of
$\rho_{LT0}^{-\frac{3}{2}}$ is small then it does not effect much the shape of
the initial difference of the density extremes. Also if at late time $t$ the
bang time function is small compared to $t$ we can see that $\Delta\rho$ is
proportional to $\frac{1}{t^{3}}$. We can try to evaluate (83) at the initial
time $t=t_{i}$

$\Delta\rho\left(t=t_{i}\right)\approx\frac{8\sqrt{3}}{9}\frac{27}{8}\frac{\left(\int\rho_{LT0}r^{2}{\rm
d}r\right)^{\frac{3}{2}}}{r^{\frac{9}{2}}}\frac{\Delta\rho_{0}}{\rho_{LT0}^{\frac{3}{2}}}=A^{\frac{3}{2}}\Delta\rho_{0},$
(84)

so we can see that if A is close to 1 the approximation gives us what we
expect.

##  6 Model specification

The Szekeres spacetime has 5 degrees of freedom, so to fully specify the model
we need to set up 5 functions. But since we are only interested in the density
contrast, it is sufficient to specify just three functions, which is the
curvature function $f$ that we assume to be zero, the initial radial profile
$\rho_{LT0}$ and one extreme value of the density. We will set up two models
$a$ and $b$. For the initial radial profile we choose

$\displaystyle a:$
$\displaystyle\quad\rho_{LT0}=\rho_{b0}\left(1+\frac{1}{10}e^{-\frac{r^{2}}{500}}\right),$
(85) $\displaystyle b:$
$\displaystyle\quad\rho_{LT0}=\rho_{b0}\left(1+\frac{1}{5}e^{-\frac{r^{2}}{500}}\right).$

Both of them are peaked at the origin and $\rho_{b0}$ is the background
density at the initial time. The $\rho_{LT0}$ is chosen in such a way that
besides the origin it is less then $\rho_{AV0}$, therefore for the extreme
value we specify the density minimum since the only requirement for it is that
it has to be equal to $\rho_{LT0}$ at the origin and less everywhere else. We
choose $\rho_{min0}$ as

$\displaystyle a:$
$\displaystyle\quad\rho_{min0}=\rho_{b0}\left(1+\frac{1}{10}e^{-\frac{r^{2}}{400}}\right),$
(86) $\displaystyle b:$
$\displaystyle\quad\rho_{min0}=\rho_{b0}\left(1+\frac{1}{5}e^{-\frac{r^{2}}{400}}\right).$

The models differ by the size of the radial inhomogeneity that is in the model
$b$ twice as big as in the model $a$ and we want to demonstrate, that the
approximation formula will work better for the model $a$, because the function
$A$ is closer to 1 than in the model $b$. The functions $\rho_{LT0}$,
$\rho_{min0}$, $\rho_{max0}$ and $\rho_{AV0}$ for the model $a$ are shown in
figure 3a. The function $A$ for both models is shown in figure 3b. The initial
time is chosen as $t_{i}=5\cdot 10^{5}y$ that approximately corresponds to the
time of last scattering. The final time $t_{f}=13.7\cdot 10^{9}y$ which is
approximately present time. The initial density contrast
$\frac{\Delta\rho_{0}}{\rho_{b0}}$ and the final density contrast
$\frac{\Delta\rho}{\rho_{b}}$ for both models are shown in figures 4,
$\rho_{b}$ denotes the background density at the final time. The impact of the
time evolution on the initial shape is minimal as we expected since $A$ is
close to $1$, on the other hand the magnitude of the density contrast drops
significantly because of the factor $t^{-3}$ in (83). We can also see in the
figure 4 that the formula (83) approximates the density contrast better for
the model $a$, that corresponds to the lower peak, because we chose the radial
inhomogeneity lower than in model $b$.

(a)  (b)

Figure 3: Panel (a): The functions $\rho_{min0}$, $\rho_{max0}$, $\rho_{LT0}$
and $\rho_{AV0}$ for the model $a$. The dotted line represents the chosen
initial density minimum $\rho_{min0}$ according to (86). The dashed line is
the initial density maximum $\rho_{max0}$ computed according to (52). The
solid line between the dotted and dashed lines is the chosen initial radial
density profile $\rho_{LT0}$ as specified in (85). The top solid line is the
computed $\rho_{AV0}$. All values on the vertical axis are divided by
$\rho_{b}$. Panel (b): The function $A$ as defined in (59) for both models.
The lower curve corresponds to the model $a$, the upper curve corresponds to
the model $b$.

(a)  (b)

Figure 4: Panel (a): The behaviour of the initial density contrast
$\frac{\Delta\rho_{0}}{\rho_{b0}}$ for both models. Panel (b): the behaviour
of the final density contrast $\frac{\Delta\rho}{\rho_{b}}$ for both models.
The solid curve represents the exact formula (73). The dotted curve is the
density contrast as calculated according to the approximation formula (83).
The lower and upper curves correspond to the model $a$ and $b$ respectively.
The lower curve is apparently approximated better, because the radial
inhomogeneity was chosen smaller.

##  7 Conclusion

We studied model of inhomogeneity in quasispherical Szekeres model. We set up
only 3 of 5 degrees of freedom, which is sufficient for studying the evolution
of the density contrast. The lack of specification of the last two degrees of
freedom means that we do not have the detailed information about the density
distribution, particularly we do not know the position of the extreme values
on the spheres. We specify the initial radial density profile $\rho_{LT0}$
which is the value of the density around the great circle that lies in the
plane that defines the dipole. Next, we choose one extreme value of the
density at the initial time, either $\rho_{min0}$ or $\rho_{max0}$, and the
last function that we choose is the curvature function $f$. In order to choose
an appropriate value for the density extremes, we investigated the shell
crossing conditions in terms of density $\rho_{max}$, $\rho_{min}$,
$\rho_{LT}$ and $\rho_{AV}$. We derived conditions that $\rho_{max}$ and
$\rho_{min}$ have to satisfy in order to avoid shell crossing and we showed
that in the special case $f=0$, if the conditions are fulfilled on the initial
time slice, then they will hold at any time.

Next, we derived an analytical formula for the density contrast $\Delta\rho$
as a function of $t$ and $r$. In the special case $f=0$, we derived an
approximation formula that is valid if the initial inhomogeneity is small and
we showed that in this approximation the density contrast is proportional to
the initial density contrast and depends on time as $t^{-3}$, so there is a
decrease in magnitude during time evolution, however the shape of the function
is preserved. It shows that the dynamics is very simple and close to
homogeneous one for small inhomogeneity confirming the expected behaviour. In
this sense one may argue that a small inhomogeneity can be successfully
treated in perturbation theory and the influence of nonlinearity is
negligible.

The next research will be focused on the situation when the curvature function
is chosen more generally and is not equal to zero.

###### Acknowledgements.

D. V. was supported by grant GAUK 398911 and project SVV-267301. O. S. was
supported by grant GAČR 14-37086G.

## References

  * (1) Lemaitre, G.: L’Univers en expansion. Ann. Soc. Sci. Bruxelles A 53, 51-85 (1933)
  * (2) Szekeres, P.: Quasispherical gravitational collapse. Phys. Rev. D 12, 2941-8 (1975)
  * (3) Szafron, D.A.: Inhomogeneous cosmologies: New exact solutions and their evolution. J. Math. Phys. 18, 1673 (1975)
  * (4) Stephani, H.: Über Lösungen der Einsteinschen Feldgleichungen, die sich in einen fünfdimensionalen flachen Raum einbetten lassen. Commun. Math. Phys. 4, 137 (1967)
  * (5) Krasinski, A.: Inhomogeneous Cosmological Models. CUP, Cambridge (1997)
  * (6) Krasinski, A., Hellaby, Ch.: Formation of a galaxy with a central black hole in the Lemaître-Tolman model. Phys. Rev. D 69, 043502 (2004)
  * (7) Krasinski, A., Hellaby, Ch.: More examples of structure formation in the Lemaître-Tolman model. Phys. Rev. D 69, 023502 (2004)
  * (8) Hellaby, Ch., Krasinski, A.: Alternative methods of describing structure formation in the Lemaitre-Tolman model. Phys. Rev. D 73, 023518 (2006)
  * (9) Bolejko, K., Krasinski, A., Hellaby, Ch.: Formation of voids in the Universe within the Lemaitre-Tolman model. Mon. Not. Roy. Astron. Soc. 362, 213-228 (2005)
  * (10) Hellaby, Ch., Krasinski, A.: You cannot get through Szekeres wormholes: Regularity, topology, and causality in quasispherical Szekeres models. Phys.Rev. D 66, 084011 (2002)
  * (11) Hellaby, Ch., Krasinski, A.: Physical and geometrical interpretation of the Szekeres models. Phys.Rev. D 77, 023529 (2008)
  * (12) Krasinski, A.: Geometry and topology of the quasiplane Szekeres model. Phys. Rev. D 78, 064038 (2008); erratum Phys. Rev. D 85, 069903(E) (2012)
  * (13) Krasinski, A., Bolejko, K.: Geometry of the quasihyperbolic Szekeres models. Phys. Rev. D 86, 104036 (2012)
  * (14) Bolejko, K.: Structure formation in the quasispherical Szekeres model. Phys. Rev. D 73, 123508 (2006)
  * (15) Bolejko, K.: Evolution of cosmic structures in different environments in the quasispherical Szekeres model. Phys. Rev. D 75, 043508 (2007)
  * (16) Bolejko, K.: Volume averaging in the quasispherical Szekeres model. Gen. Rel. Grav. 41, 1585-1593 (2009)
  * (17) Bolejko, K., Celerier, M.N.: Szekeres Swiss-cheese model and supernova observations. Phys. Rev. D 82, 103510 (2010)
  * (18) Bolejko, K.: The Szekeres Swiss Cheese model and the CMB observations. Gen. Rel. Grav. 41, 1737-1755 (2009)
  * (19) Walters, A., Hellaby, Ch.: Constructing Realistic Szekeres Models from Initial and Final Data. JCAP 1212, 001 (2012)
  * (20) Hellaby, Ch.: The nonsimultaneous nature of the Schwarzschild R=0 singularity. J. Math. Phys. 37, 2892-905 (1996)
  * (21) Berger, B.K., Eardley, D.M., Olson, D.W.: Note on the spacetimes of Szekeres. Phys. Rev. D, 16, 3086-9 (1977)
  * (22) Hellaby, Ch.: A Kruskal-Like Model with Finite Density. Class. Quantum Grav., 4, 635-50 (1987)
  * (23) Sussman, R., A., Bolejko, K.: A novel approach to the dynamics of Szekeres dust models. Class. Quantum Grav., 29, 065018 (2012)