{'fraction_non_alphanumeric': 0.052311561089581786, 'fraction_numerical': 0.027070646411379993, 'mean_word_length': 4.033010625737899, 'pattern_counts': {'":': 0, '<': 35, '<?xml version=': 0, '>': 54, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 30, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Exact and Perturbed Friedmann-Lemaitre Cosmologies

In this thesis we first apply the 1+3 covariant description of general relativity to analyze n-fluid Friedmann-Lemaitre (FL) cosmologies; that is, homogeneous and isotropic cosmologies whose matter-energy content consists of n non-interacting fluids. We are motivated to study FL models of this type as observations suggest the physical universe is closely described by a FL model with a matter content consisting of radiation, dust and a cosmological constant. Secondly, we use the 1+3 covariant description to analyse scalar, vector and tensor perturbations of FL cosmologies containing a perfect fluid and a cosmological constant. In particular, we provide a thorough discussion of the behaviour of perturbations in the physically interesting cases of a dust or radiation background.


Cosmological Models
The goal of cosmology is to describe the large-scale structure, origins and evolution of the universe as a whole. Within the framework of general relativity, a cosmological model, denoted (M, g, u) consists of a differential manifold M, a Lorentzian metric g (which describes the geometry of the space-time) and a unit timelike vector field u, g ij u i u j = −1, (1.1) that defines the world-lines of a family of fundamental observers. In cosmology one represents the matter content as a smooth distribution (a continuum) whose particles represent galaxies in the present epoch.
The time-like vector field u is the 4-velocity of the particles in this continuum. One member of the congruence of time-like curves defined by u, namely the world-line of our galaxy, plays a particularly important role, since it is the position from which we make observations. The evolution of the universe is assumed to be governed by the Einstein field equations (EFEs), where R ab is the Ricci tensor, R is the Ricci scalar (R = g ab R ab ) and T ab is the energy-momentum tensor of the distribution of matter and energy in the universe. Here and throughout we will use geometerized units so that c = 1 and 8πG = 1, where c is the speed of light in vacuum and G is the gravitational constant.
There are two fundamental assumptions that are used to construct models of the universe: I. The universe is isotropic about our galaxy. More precisely, at any event on the world-line of our galaxy, all directions orthogonal to the world-line are physically and geometrically equivalent, when the universe is viewed on a sufficiently large scale.
A variety of cosmological observations (e.g. number counts of distant galaxies and peculiar velocity surveys) provide support for assumption I in the present epoch. The strongest support is provided by the temperature of the cosmic microwave background radiation (CMBR) which is isotropic to within 10 −5 (i.e. ∆T T ∼ 10 −5 ).
II. The Copernican Principle. This principle states that our galaxy is not in a preferred position in the universe. Unlike assumption I, this assumption is philosophical, rather than observational.
Assumptions I and II together imply that the universe satisfies the so-called Cosmological Principle, namely i) the universe is isotropic about all fundamental world-lines, and ii) the universe is spatially homogeneous, i.e. all physical and geometric quantities do not vary in space.
It is a well-known consequence of assumptions I and II (we give a concise proof in Chapter 2; see Theorem 2.1.1) that the metric g ab can be written in the form of the Robertson-Walker (RW) metric, 1 whereg αβ is a 3-metric of constant curvature. The preferred timelike vector field u is given by (1.4) The Robertson-Walker metric provides the basis for the Friedmann-Lemaître (FL) cosmological models that play a central role in modern cosmology. A Friedmann-Lemaître cosmology 2 is a cosmological model (M, g, u) that satisfies the Einstein field equations (1.2) with a suitable energy-momentum tensor T ab , with metric g and fundamental congruence u having the Robertson-Walker form (1.3) and (1.4).
In this thesis we study FL cosmologies and their perturbations within a general mathematical framework. In section 1.2 we provide the necessary geometrical background for this study. 1 Here we use Greek indices to denote spatial indices. They can take on the values 1, 2 or 3. 2 The first explicit FL cosmological models were first discovered by Friedmann in 1929 (for some historical details, see Ellis (1989)).

Projection and Differentiation
Projection parallel to u is achieved by simply contracting a vector or tensor with u a . Projection into the 3-space orthogonal to u is defined by the projection tensor which satisfies the identities Corresponding to these two types of projection, two derivatives are defined.
The overdot˙denotes differentiation along the fundamental congruence, defined for scalars aṡ and for tensors asṪ where ∇ a denotes covariant differentiation.
Following Maartens (1997), we use angle brackets to denote the orthogonal projections of vectors, the orthogonally projected symmetric trace-free part of tensors, and their time derivatives, i.e.
Further, we will use the projection symbol ⊥ to identify orthogonal projection over all tensor indices, i.e.
(T ab cd ) ⊥ = h a p h b q h r c h s d T pq rs . (1.10) The spatial covariant derivative∇ a is defined by projecting all free indices of the covariant derivative (compare Ellis (1973), ). For example, (1.11) The covariant derivative of a scalar can then be decomposed into a derivative along the fundamental congruence u and a spatial covariant derivative via ∇ a f =∇ a f − u aḟ . (1.12) The totally skew pseudotensor is defined by where g = det(g ab ). (1.14) The projected permutation tensor is then defined by (compare Maartens et al. (1995)), ǫ abc = η abcd u d . (1.15) and is used to define the spatial curl of a vector by curl X a = ǫ abc∇ b X c .

Kinematic Quantities
We now focus on the decomposition of the covariant derivative of the fundamental 4-velocity. Following Ellis (1973, p8-11), we decompose the covariant derivative of a unit timelike vector field u as ∇ b u a = σ ab + ω ab + Hh ab −u a u b , (1.18) with σ ab symmetric and trace-free, ω ab antisymmetric, and σ ab u a = 0, ω ab u a = 0,u a u a = 0.
(1. 19) It follows from (1.18) and (1.19) that Here, H is the Hubble scalar 5 ,u a is the acceleration vector, σ ab is the rate of shear tensor and ω ab is the vorticity tensor. In practice we will usually replace ω ab by the vorticity vector, defined by ω a = 1 2 η abcd u b ω cd (1.24) which satisfies u a ω a = 0. We can write the vorticity tensor in terms of the vorticity vector using the relation ω ab = η abcd ω c u d . (1.25) The magnitude of the rate of shear tensor and the vorticity tensor are defined by σ 2 = 1 2 σ ab σ ab , ω 2 = ω a ω a .

Proof:
Given ω a =u a = 0, equation (1.23) implies u [a;b] = u [a,b] = 0. It then follows from Proposition 1.2.1 that f ,[a t ,b] = 0. This in turn implies that f = g(t) and hence that one can redefine t to set f = 1.
In this case, the scalar t that labels the hypersurfaces is clock time along the normal congruence.

Decomposition of the Energy-Momentum Tensor
The energy-momentum tensor T ab is a symmetric tensor that specifies the total matter-energy content of a cosmological model. We use the standard decomposition of the energy-momentum tensor with respect to a timelike vector field u (see Ellis (1973), p7): T ab = µu a u b + q a u b + q b u a + ph ab + π ab , (1. 26) where q a u a = 0, π ab u b = 0 π a a = 0, π ab = π ba .
(1. 27) In this equation, µ is the total energy density measured by an observer moving with fluid 4-velocity u a , q a is the energy flux relative to u a (representing heat conduction and diffusion), p is the isotropic pressure and π ab is the trace-free anisotropic pressure (due to processes such as viscosity).
A perfect fluid with 4-velocity u is defined by q a = 0, π ab = 0, (1.28) in which case the energy-momentum tensor (1.26) simplifies to T ab = µu a u b + p(g ab + u a u b ). (1.29)

The Weyl Curvature Tensor
The Weyl conformal curvature tensor (or simply the Weyl tensor) is defined in terms of the Riemann tensor and the Ricci tensor according to C abcd = R abcd − 1 2 (g ac R bd + g bd R ac − g bc R ad − g ad R bc ) + 1 6 R(g ac g bd − g ad g bc ).
(1.30) (see, for example, Ellis 1973, p7). We note that the Weyl curvature tensor has all the symmetry properties of the Riemann tensor and is in addition trace-free (C a bad = 0).
It is helpful to view (1.30) as a decomposition of the Riemann tensor R abcd , with its 20 independent components, into a tracefree part C abcd and the trace R ab : each of which have 10 independent components. The Ricci tensor is determined at each point by the energy-momentum tensor through the field equations (1.2). On the other hand, as we shall see later, the Weyl tensor governs the propagation of changes in the gravitational field, i.e. gravitational waves.
Within the 1 + 3 formalism the Weyl tensor can be decomposed further into an "electric part" E ab and a "magnetic part" H ab by projecting C abcd and its dual orthogonal to u: (1. 32) We note that E ab and H ab are symmetric, trace-free and orthogonal to u: (1. 33) In addition, this decomposition describes the Weyl tensor completely since 7 E ab = 0 = H ab ⇐⇒ C abcd = 0. (1.34) Note that through (1.30) E ab , H ab and R ab provide a full description of the space-time curvature as described by the Riemann tensor.

Spatial Curvature
If u is irrotational, we have seen that there exists a family of hypersurfaces (i.e. 3-dimensional manifolds) orthogonal to u. In this situation the projection tensor h ab can be viewed as the 3-metric induced on these hypersurfaces and the spatial covariant derivative∇ a can be viewed as the covariant derivative determined by this 3-metric.
Analogous to the spacetime Ricci identities, the 3-space Ricci identities are given bŷ where (3) R abcd is the spatial Riemann tensor and X a is any vector field in the hypersurfaces, i.e. X a u a = 0. It can be shown by applying (1.35) to an arbitrary vector field X that is orthogonal to u, expanding the left hand side using (1.11) and (1.5) and introducing the spacetime Riemann tensor via that (3) R abcd is related to R abcd via (3) R abcd = (R abcd ) ⊥ − Θ ac Θ bd + Θ ad Θ bc , (1.37) where Θ ab is the expansion tensor of the congruence u, defined by Θ ab = σ ab + Hh ab .
The spatial Ricci tensor and spatial Ricci scalar are then defined in the usual manner: Further, the trace-free spatial Ricci tensor is defined as (3) S ab = (3) R ab − 1
(1. 40) We can obtain an expression for the spatial Ricci scalar and trace-free spatial Ricci tensor in terms of the spacetime quantities from (B.5), incorporating the Weyl tensor through (1.30) and (1.32): S ab = E ab + 1 2 R ab − Hσ ab + σ c a σ b c . (1.42) Then using the field equations (1.2) in the form (A.2)-(A.5) to eliminate R and R ab , we obtain S ab = E ab + 1 2 π ab − Hσ ab + σ c a σ b c .
(1. 44) We conclude this section with a brief digression on spaces of constant curvature. A differentiable manifold M of dimension n with metric g is said to be a space of constant curvature if and only if R abcd = C(g ac g bd − g ad g bc ).
(1. 45) It follows by contraction that R = n(n − 1)C, and R ab − 1 n Rg ab = 0, (1. 46) i.e. the tracefree Ricci tensor S ab ≡ R ab − 1 n Rg ab is zero. In addition it follows from the twice contracted Bianchi identities that if n > 2 then R is constant (see proposition 1.2.3). If n = 3, the Riemann tensor is completely determined by the Ricci tensor according to 8

Evolution and Constraint Equations
In the 1 + 3 covariant formalism, the information contained in the EFEs is displayed in an indirect manner by using the Ricci identities applied to the fundamental 4-velocity u, The latter can be written in divergence form by contracting once and then twice and using (1.30) to introduce the Weyl tensor: 9 8 See, for example, Ellis (1973) eq. (80). 9 Also see Wainwright and Ellis (1997), p27. and The Riemann tensor in (1.50) is likewise expressed in terms of the Weyl tensor and Ricci tensor using (1.30). In all these identities, the Ricci tensor is expressed in terms of the energy-momentum tensor (1.26) using the EFEs in the form Finally, one projects all three identities parallel to and orthogonal to u. The result is a set of evolution equations giving the derivatives of H, σ ab , ω a , E ab + 1 2 π ab , H ab , µ and q a along the congruence u, and a set of constraint equations that contain only the spatial covariant derivative∇ a . We refer the reader to Appendix A for the complete list of equations that arise from this process.

Overview of the Thesis
The main goal of this thesis is twofold. Firstly, we analyze the dynamics of so-called n-fluid FL cosmologies, that is, homogeneous and isotropic cosmologies whose matter-energy content consists of n non-interacting fluids. Secondly, we give a unified coordinate-independent discussion of perturbations of FL cosmologies containing a perfect fluid and a cosmological constant. Some parts of this thesis are new, while other parts provide a unified presentation of known results for these cosmologies. Our treatment throughout is completely coordinate-independent since it is based on the so-called 1 + 3 covariant description of general relativity and in this respect it differs from much of the current literature.
In Chapter 1, we introduce the 1 + 3 covariant description of general relativity and in Appendix A we give the full system of evolution equations and constraint equations. In this sense the thesis is selfcontained. However, a reader who has no prior knowledge of the 1+3 covariant description will encounter an initially steep learning curve. We assume the reader has a basic knowledge of general relativity and some familiarity with dynamical systems theory (see, for example, Wainwright and Ellis (1997), Chapter 4).
Chapter 2 plays a supporting role. In this chapter we describe the geometrical, i.e. coordinate independent properties of the Robertson-Walker metric, and give two characterizations of the Friedmann-Lemaître models that highlight their very special nature.
In Chapter 3, we give a complete qualitative description of the dynamics of n-fluid FL cosmologies.
A central role is played by a new formulation of the evolution equations for these models as a dynamical system on a compact n-dimensional state space. The most important example is the 3-fluid model in which the three fluids are radiation, pressure-free matter and a cosmological constant, since it provides a model of the physical universe. Our discussion of this model is one of the highlights of this chapter.
In Chapter 4, we give a unified derivation of the governing equations for scalar, vector and tensor perturbations of FL cosmologies containing a perfect fluid and a cosmological constant, using the so-called geometrical approach pioneered by . We then analyze the behaviour of the perturbations in various regimes using both exact and approximate solutions. In particular we determine how the perturbations affect the so-called observational quantities that determine the extent to which a perturbed model deviates from an exact FL model.
The Appendices contain a wealth of additional material relating to perturbations of FL universes, not all of which is needed for the purposes of this thesis.

COSMOLOGIES
In this chapter we first derive the geometrical (i.e. coordinate independent) properties of the Robertson-Walker metric within the framework of the 1 + 3 description of general relativity. We then give two characterizations of the Friedmann-Lemaître cosmologies, firstly using the kinematic quantities and secondly using the Weyl tensor.

Properties of the RW Metric
In this section we show that the high symmetry of the RW line-element leads to strong restrictions on the kinematic quantities of the fundamental congruence, on the Weyl tensor and on the Ricci tensor. We then derive the RW line-element assuming that the Cosmological Principle holds.
Since the projection tensor h ab has components h αβ = g αβ , h a0 = 0, it follows that The decomposition (1.18) immediately implies σ ab = 0, ω ab = 0,u a = 0 and Thus, H is purely a function of t, which implies∇ a H = 0.
We now consider the restrictions on the spatial curvature associated with the Robertson-Walker metric.
We note that H ab = 0 is an immediate consequence of (2.10), on applying proposition 2.
Proof: We note that R ab = 0 follows immediately from (2.11). In order to show R ab u a h b c = 0, we make use of the contracted Ricci identities in the form 2 where t is a scalar field. We then introduce co-moving coordinates x α , i.e. coordinates x α that are constant along the world lines of the fundamental congruence. Then the metric and 4-velocity assume the form Next, using the metric (2.17), it can be quickly shown that  Proof: The Cosmological Principle implies that any geometrically defined vector orthogonal to u must be zero, since otherwise it would define a preferred direction orthogonal to u. In particular, we have ω a =u a =∇ a H = 0. Likewise, any geometrically defined rank two symmetric tracefree tensor X ab that is orthogonal to u (X ab u b = 0) must be zero, since otherwise its eigenvectors would define a preferred direction orthogonal to u. In particular we have σ ab = (3) S ab = 0. The desired result thus follows from proposition 2.1.5.

Characterizations of the FL Cosmologies
In this section we present two characterizations of FL cosmologies in terms of the kinematic quantities and Weyl curvature.
We can now present the main result of this section, which we shall refer to as the Weyl characterization of the Friedmann-Lemaître models.

COSMOLOGIES
The main goal of this chapter is a description of the dynamics of n-fluid FL cosmologies. In section 3.1, using the results from chapter 1, we derive the basic evolution equations for n-fluid FL cosmologies and describe their behaviour using dynamical systems techniques. In section 3.2, we describe the degrees of freedom of the n-fluid FL cosmologies in two ways, firstly, using intrinsic parameters and secondly, using

Behaviour and Classification of n-Fluid FL Cosmologies
In this section we derive the evolution equations for n-fluid FL cosmologies without introducing local coordinates and a line-element, by specializing the general evolution and constraint equations (appendix A) using the restrictions derived in chapter 2. We formulate the evolution equations as ODEs for a set of dimensionless variables in a compact state, which enables us to use dynamical systems methods to describe the evolution of n-fluid FL cosmologies qualitatively.

Basic Evolution Equations
We know from theorems 2.2.1 and 2.2.2 and proposition 2.1.4 that FL models have the following properties: i) the kinematic quantities of the fundamental congruence u satisfy ii) the Weyl curvature is zero: and iii) the stress-energy tensor satisfies π ab = 0, q a = 0, and hence are purely functions of clock time t along the fundamental congruence. With the above restrictions, the general evolution and constraint equations reduce to two evolution equations for µ and H, namely (A.8) and (A.14), which readμ where the overdot denotes differentiation with respect to t. The spatial geometry of the model is determined by the Gauss equation in the form (1.43) and (1.44). It follows that the tracefree Ricci tensor is zero, and that the spatial Ricci scalar is determined by H and µ according to In order to obtain the standard evolution equations, we introduce the length scale function ℓ(t), which is defined by H, up to a constant multiple, according to This definition is consistent with the length scale defined by the RW metric (1.3).
To proceed we need the following result: Since ℓ is only defined by (3.9) up to a multiplicative constant, we may choose it to satisfy 1 1 6 (3) where K is the spatial curvature parameter. The value for K indicates whether (3) R is positive, negative or zero. Models where K = −1, 0 and +1 are then referred to, respectively, as open FL models, flat FL models and closed FL models.
It is important to note that if K = ±1, then ℓ(t) is completely determined via eq (3.11), by the spatial geometry, represented by (3) R. On the other hand, if K = 0, equation (3.11) does not restrict ℓ(t), which means that if K = 0 then ℓ(t) can be rescaled with a multiplicative constant.

n-Fluid FL Cosmologies
In this chapter, we consider FL cosmological models whose matter content consists of n fluids with the same four-velocity u. The total energy density and pressure are equal to a sum of the energy density and pressure of the individual components: The fluids are assumed to be non-interacting, which means that each fluid component individually satisfies the conservation equation (3.12), i.e.μ We further assume that each fluid obeys a linear equation of state with equation of state parameters γ i (i = 1, . . . , n), i.e.
The most important physical values are γ i = 1 for pressure-free matter and γ i = 4 3 for radiation. The value γ i = 0 can be used to represent a cosmological constant and γ i = 2 is sometimes considered to correspond to a so-called stiff fluid. For definiteness we order the equation of state parameters according to γ 1 > γ 2 > · · · > γ n , (3.19) and assume that γ 1 ≥ 1, i.e. at least one matter component is matter with a non-negative pressure. The conservation equations (3.16) then assume the forṁ As well, the Friedmann equation (3.14) and Raychaudhuri equation (3.13) assume the forṁ In summary, we will use the term n-fluid FL cosmology to refer to an FL cosmological model whose matter content consists of n non-interacting fluids with the same 4-velocity u and equations of state given by (3.17)- (3.19). Following Jantzen and Uggla (1992), we will use abbreviations such as RDC-universe and RCΛ-universe to identify n-fluid FL cosmologies, where R denotes radiation, D denotes dust, C denotes curvature and Λ denotes a cosmological constant.
We have shown that the state of an n-fluid FL cosmology as a function of time t is determined by the functions (ℓ, µ 1 , . . . , µ n ) which satisfy the system of ODEs (3.20) and (3.21). The state of an n-fluid cosmology can equivalently be described by the variables (H, µ 1 , . . . , µ n ), which satisfy the evolution The first of these is a restatement of (3.20) using the definition (3.9) of H. The second follows from (3.6) on using (3.15) and (3.17). Further, upon applying (3.9), the Friedmann equation (3.21) gives an additional constraint that must be satisfied, namely In the case of dust (γ = 1), this solution is the Einstein-de Sitter universe (Einstein and de Sitter (1932)).

Hubble-Normalized Scalars
In this section we introduce so-called Hubble-normalized scalars, which form a set of dimensionless variables that, together with the Hubble scalar H, describe the state of n-fluid FL cosmologies. These dimensionless scalars also play an important role in relating the models to observations (see section 3.2.2).
Firstly, the density parameters Ω i , which describe the dynamical effect of the matter densities, are defined by The total density parameter Ω, which describes the dynamical effect of the total matter content, is defined by It follows from (3.15) and (3.31) that Lastly, the dynamical effect of the spatial curvature is described by Ω k , defined by It follows from (3.11) that In terms of the parameters Ω and Ω k , the Friedmann equation (3.21) assumes the simple form which may be further expanded in terms of Ω i as We note that (3.36) implies that Ω determines the spatial curvature as follows: There is an additional dimensionless variable that plays a fundamental role, namely the deceleration parameter q, defined by Clearly, if q < 0 then the expansion is accelerating (l > 0) and if q > 0 then the expansion is decelerating (l < 0). In an n-fluid FL model, q can be expressed in terms of the density parameters by using the Raychaudhuri equation (3.22) with (3.9) and (3.31). We write The quantities Ω i , Ω, Ω k and q are dimensionless, whereas H has dimensions of (length) −1 . Together, these variables are fundamental quantities in the FL models and more general models (see, for example, Wainwright and Ellis (1997), section 5.2 and chapters 12 and 13).

A Compact State Space for n-fluid FL Cosmologies
The one drawback of using H-normalized variables is that they cannot give a unified description of models that have both an expanding (H > 0) and a contracting (H < 0) epoch. In particular, when H = 0 the density parameters (3.31) and (3.34) are undefined (i.e. Ω i → +∞ as H → 0). In order to circumvent this difficulty, we introduce a set of variables which are bounded throughout the evolution of the model.
For an n-fluid FL cosmology, with total energy density we define the dimensionless matter variables, which describe the relative significance of the different matter components: i.e. the χ i are bounded.
Following Wainwright (1996, p124-5), we now introduce a curvature parameterΩ so as to compactify the state space:Ω and that sign(Ω) = sign(H). gives the flat FL models. We now derive evolution equations forΩ and the χ i . It is convenient to introduce a quantityq that plays the role of the deceleration parameter q, but unlike q, is bounded when H = 0. Motivated by (3.42), we defineq Together, (3.19) and (3.50) imply thatq is restricted by The form of (3.55) suggests defining a new time variableτ by It is essential to note that this change of variable is well-defined, since cosΩ > 0 on account of (3.48). In terms ofτ , (3.55) reads We differentiate (3.54) with respect toτ and use (3.58). After simplifying using (3.54), we obtain dΩ dτ = −q cos 2Ω cosΩ. The state space is not compact due to the strict inequality in (3.48) and so we seek to compactify this space, since then we can apply theorems from the theory of dynamical systems.
In order to compactify the space, we simply observe that the right hand sides of the ODEs (3.57) and (3.60) are well-defined, in fact analytic, at the pointsΩ = ± π 2 , which form part of the boundary of the state space. We can thus extend the state space to include the points withΩ = ± π 2 , thereby compactifying it. Since lim as follows from (3.54), we interpret the setsΩ = ± π 2 as representing FL cosmologies with µ = 0. It follows from (3.14) that the setΩ = ± π 2 represents the Milne universe, given by (3.28).

State Space Representation of n-Fluid FL Cosmologies
State Space: (3.62) Evolution Equations:

Comments:
i) The ODEs (3.63) and (3.64) describe the essential dynamics of n-fluid FL cosmologies in terms of the dimensionless variables (Ω, χ 1 , . . . , χ n ) and the dimensionless timeτ . Since the state space, as defined by (3.62) is compact, the solutions of the ODEs are defined for allτ ∈ R, and as a result the ODEs define a dynamical system on the state space (see, for example, Wainwright and Ellis (1997), p87, Corollary 4.1). The evolution of an n-fluid FL cosmology is thus described by an orbit of the dynamical system (a solution curve of the ODEs). As a consequence, we are guaranteed the existence of a past attractor A − of the ODEs (i.e. an invariant set to which all orbits, except possibly a set of measure zero, are past asymptotic to asτ → −∞) and a future attractor A + (i.e. as a past attractor, except asτ → +∞). The past attractor describes the asymptotic regime at early times and the future attractor describes the asymptotic regime at late times.
ii) This description of the evolution of n-fluid FL cosmologies as a dynamical system on a compact state space is new. It generalizes the formulation given by Wainwright (1996) for the case n = 2.
iii) The physical state of an n-fluid FL cosmology is determined by the variables (H, µ 1 , . . . , µ n ), which are related to the dimensionless variables via the total density µ, according to  This freedom implies that each orbit determines a 1-parameter family of cosmological models, whose physical variables are related by a multiplicative constant. Since an n-dimensional dynamical system has an (n − 1)-parameter family of orbits, it follows that the class of n-fluid FL cosmologies is labelled by n parameters. We will specify these parameters explicitly in section 3.2.
iv) Clock time t is related toτ via the ODE (3.56) and can thus be expressed as an integral involving Ω(τ ) and µ(τ ), the latter determined from (3.67). Clock time will take on values in an interval where t i is finite or −∞ and t f is finite or +∞, depending on the model.

Invariant Sets and Equilibrium Points in the Compact State Space
Our goal in the rest of this section is to describe the qualitative behaviour of the orbits in the compact state space, thereby understanding the dynamics of n-fluid FL cosmologies. Since the state space is compact, it follows that every orbit has a past and future attractor (see, for example, Wainwright and Ellis (1997), p99-100).
In general the attractors of a dynamical system are complicated invariant sets. But in the present situation we will show that they are determined by the equilibrium points of the dynamical system (section 3.1.7), which we now present.
2) Milne set. This set is the n − 1 dimensional subset given byΩ = ± π 2 (which corresponds to Ω = 0 and Ω k = 1 on account of (3.69) and (3.36)). The orbits in this set describe the Milne vacuum solution.
3) m-fluid sets (1 ≤ m < n). These sets are determined by setting n − m of the χ i variables to be zero, with the m remaining χ i variables being non-zero. The orbits in these sets describe m-fluid FL cosmologies.
3) Einstein static points (E) are given bỹ and so constitute an (n − 2)-dimensional set of equilibrium points. These equilibrium points correspond to the n-fluid interpretation of the Einstein static model (3.29) and only exist when γ n < 2/3.

Classification of n-Fluid FL Cosmologies
In this section, we will prove that the n-fluid FL cosmologies that have an epoch of expansion (i.e. H > 0 in some time interval) form three qualitatively different generic 4 subclasses, defined as follows: 1) Ever-expanding models satisfy H > 0 for all t.
The main result of this section now follows immediately. To the best of the author's knowledge, the classification theorem is new for n > 2, although it is a well known result for n = 2, when the matter content is a perfect fluid (γ 1 > 2 3 ) and a cosmological constant (γ 2 = 0) (see Rindler (1977), p234-8).

Asymptotic Behaviour of Solutions
In this section we determine the asymptotic behaviour of the three generic classes of n-fluid FL cosmologies identified in the previous section. We will first describe the behaviour of the length scale ℓ(t) as a function of clock time t in each of the three generic classes of FL cosmologies, and then use these results to determine the past attractor and future attractor. Important in this analysis is the set of orbits asymptotic to an Einstein static solution. In the following proposition, we demonstrate that orbits of this type are non-generic. Proof: We can write the evolution equations (3.57) and (3.60) as a dynamical system of the form The physical state space is the hyperplane given by (3.62), and hence is n dimensional. Consider the Jacobian matrix J ij = ∂F j ∂X i . Upon evaluating J ij on the Einstein static set, given by (3.74), we obtain where i, j = 1, . . . , n and χ i,E denotes the value of χ i at the Einstein equilibrium point. The Jacobian matrix, when restricted to the physical state space, has two non-zero eigenvalues and n − 2 zero eigenvalues. 7 It follows 8 that each Einstein static equilibrium point has a 1-dimensional stable manifold and a 1-dimensional unstable manifold. Thus the set of orbits future asymptotic to the Einstein static set has dimension exactly one larger than the dimension of the Einstein static manifold, i.e.
In the next proposition, we establish the behaviour of the length scale ℓ as a function of timeτ , for ever-expanding models.
We now determine the asymptotic behaviour of clock time t (i.e. asτ → ±∞). Let where the dependence of t onτ is determined via equation (3.56). Then t i and t f are finite or infinite depending on the type of universe and during the evolution t i < t < t f . Clock time is related to the length scale via the Friedmann equation. In particular, for any n-fluid FL cosmology, we can integrate (3.20) to where λ i is an integration constant whose inverse can be interpreted as the characteristic length scale of fluid i. On substituting (3.88) into the Friedmann equation (3.21), we obtain iii) For a bouncing n-fluid FL model, t i and t f are infinite and ℓ(t) satisfies ℓ(t) → +∞ as t → ±∞.

Proof: Case i)
In this case, we will make use of the following elementary results from calculus: Recall that ever-expanding models satisfy H(τ ) > 0 ∀τ . Thus by Proposition 3.1.5, solutions which are not asymptotic to the Einstein static solution satisfy ℓ → 0 asτ → −∞ and ℓ → ∞ asτ → ∞. Let For an ever-expanding model, this function is defined and continuous (for 0 < a < ∞) since H > 0 implies that the denominator is strictly positive (see (3.89)). Hence, we can integrate (3.89) to obtain Similarly, we can integrate (3.89) to obtain Further, the integrand satisfies For the following two cases we make use of the fact thatq andl are related bÿ Case ii) For a recollapsing model, Proposition 3.1.6 impliesq(τ ) > 0 for allτ ∈ R. Then (3.95) implies that ℓ(t) is concave down on its domain and so we can conclude there exists t i and t f such that ℓ(t i ) = 0 and ℓ(t f ) = 0.
We must now show that ℓ does not go to infinity in finite time. Choose ℓ 0 so that H > 0 for ℓ ≥ ℓ 0 .
The past and future attractors of the dynamical system on the compact state space can now be determined using theorem 3.1.2.

97)
and Proof: We begin by obtaining algebraic expressions for χ i andΩ. From (3.88) we have Using (3.47), we write the Friedmann equation (3.25) in the form Then using (3.45) and (3.99) we have To determine the limits of (3.100) and (3.102) as ℓ → 0, we write (3.103) and, using (3.44), (3.99) and (3.100), It follows that lim and (3.108) It now follows from Theorem 3.1.2 and equations (3.72) and (3.73) that orbits describing models in the three generic classes approach the equilibrium points given in table 3.1, asτ → −∞ and asτ → +∞.

Comment:
Since the orbits describing models in the three generic classes approach the equilibrium points in table 3.1 asτ → ±∞, it follows that the density parameters Ω i and Ω k approach their values at the equilibrium points, which we give in table 3.2.

Parameterization of FL Models
In this section we introduce the n essential parameters that label the family of n-fluid FL cosmologies, namely a conformal parameter λ and n − 1 dimensionless mass parameters m 1 , . . . , m n−1 , that we shall refer to as intrinsic parameters. In section 3.2.2, we introduce an equivalent set of parameters, namely the observational parameters Ω i,0 and H 0 and show how they determine the intrinsic parameters. Finally, we show that the evolution equations (3.57) and (3.60) admit n − 1 independent conserved quantities, which characterize the orbits in the n-dimensional state space. In particular, we will show that the values of the conserved quantities are determined in terms of the intrinsic parameters.

The Intrinsic Parameters
We now integrate the conservation equation and demonstrate that the constants of integration of this equation can be used to parameterize n-fluid FL cosmologies. As we have seen previously in (3.88), on The constants of integration λ i have dimension (length) −1 and so each can be thought of as defining a characteristic length scale for each fluid via λ −1 i . When (3.109) is substituted into (3.21), there results a first order ODE for ℓ(t): dℓ dt In order to define a dimensionless length scale and time variable for n-fluid FL cosmologies, one needs to choose a preferred λ i that will play the role of an overall conformal factor. Typically we will choose λ to be λ n , which is a natural choice if γ n = 0. In this case, it follows from (3.109) that The dimensionless length L and time T are then defined as L = λℓ, and T = λt.  10 For i = n, this equation states that mn = 1.
so that In dimensionless form, the Friedmann equation (3.110) then reads It then follows from (3.112) that ℓ(t) can be written as The line element then takes the form where L is given by (3.117) and dΣ 2 is a metric of constant curvature. 11 From (3.9), (3.112) and (3.115) we obtain Further, using (3.31), (3.114) and (3.120) the density parameters Ω i are written as In summary, equations (3.114) and (3.118) show that the n-fluid FL cosmologies form an n-parameter family, labelled by a conformal factor λ and n − 1 mass parameters m 1 , . . . , m n−1 . We will refer to these parameters as the intrinsic parameters. then so isgij = λ −2 gij for any constant λ, with the matter terms scaling appropriately, e.g.μ = λ 2 µ.

The Observational Parameters
In this section we introduce the observational parameters associated with n-fluid FL cosmologies, based on the current epoch t 0 , which represents the age of the universe. We denote the values of H, q, Ω i and Ω k at t 0 by where H 0 is the Hubble constant. Then, the set of constants are collectively referred to as the observational parameters (see, for example, Wainwright and Ellis (1997), p55-8).
In addition, we require the parameter ℓ 0 , which denotes the length scale at the present time, i.e ℓ 0 = ℓ(t 0 ). This quantity is not an observational parameter, but if K = 0, ℓ 0 can be determined in terms of observational parameters using (3.35) and (3.36), i.e. via .

(3.123)
In the case of zero curvature (K = 0), ℓ 0 can be chosen arbitrarily in order to fix the multiplicative factor that determines ℓ (also see section 3.1.1).
We now present the observational formulation of the n-fluid FL cosmologies. In this formulation, the dimensionless length scale and dimensionless time variable are given by Using (3.124) and (3.125), the Friedmann equation (3.14) assumes the form We use uniqueness of the solution to (3.126) under some appropriate initial conditions (see section C. 2) to conclude that in an expanding epoch the length scale is given by It follows from (3.9), (3.124) and (3.126) that the Hubble scalar assumes the form 12 and further that the density parameters (3.31) and (3.45) assume the form and The parameters Ω k,0 , q 0 , t 0 in (3.122) can be related to the parameters {Ω i,0 , H 0 }, as follows: First, Ω k,0 is determined in terms of Ω i,0 via the Friedmann equation (3.37): Third, knowing H 0 , we can obtain t 0 by integrating (3.126) as and Ω k,0 has been eliminated using (3.131). As a consequence, the only independent observational parameters are the n + 1 quantities defined by It is useful to note that the normalized scale factor can be expressed in terms of the redshift z, defined by (see Wainwright and Ellis (1997), eq. (2.41)). One can thus use (3.128) and (3.129) to express H and Ω i in terms of z.
We now present the link between the observational parameters and the intrinsic parameters. On sub- (3.137)

Conserved Quantities
We now show that the evolution equations (3.57) and (3.60) admit n − 1 conserved quantities.

From (3.45) and (3.114) we have that
We now eliminate λ 2 /µ and L by taking products of powers of three distinct matter variables (we require n ≥ 3). Indeed, it follows from (3.138) that In the case of K = 0, we also obtain a conserved quantity describing the curvature of a model. We use (3.9) and (3.47) to rewrite the Friedmann equation (3.14) as Then from (3.138) and (3.140), we have Since there can be only n − 1 independent conserved quantities, we must fix some of the values of i, j and k in (3.139) and (3.141) to obtain a minimal independent set. We choose i = 1, k = n in (3.139) (with 2 ≤ j ≤ n − 1) and i = 1 and j = n in (3.141) and so define 13 (3.142) 13 Observe that if we choose j = 1 or j = n in (3.142) we simply obtain M1 = 1 and Mn = 1. and Upon recalling that m n = 1, the conserved quantities are then determined by the intrinsic parameters m 1 , . . . , m n−1 according to and (3.145)

Comment:
We have derived these conserved quantities algebraically, and it follows immediately from It involves a lengthy calculation to verify (3.146) using the ODEs (3.63) and (3.64).
If K = 0 (i.e.Ω = π 4 ) then K = 0. As described in section 3.2.2, we can use the freedom to rescale ℓ 0 to set The conserved quantities (3.142) and (3.143) can also be expressed in terms of the density parameters, as follows: The conserved quantities are thus related to the observational parameters by and upon evaluating (3.149) and (3.150) at t = t 0 . We note that conserved quantities of this type have been previously discovered, for example, by Ehlers and Rindler (1989) for RDCΛ FL universes, and more generally by Lake (2006) for n-fluid FL universes, using a different method.

2-Fluid FL Cosmologies -Qualitative Analysis
We now use the results of sections 3.1 and 3.2 to give a qualitative analysis of the dynamics of 2-fluid FL cosmologies. The results of this section, in particular figures 3.1-3.3, are contained in Wainwright (1996), although the use of the conserved quantity K is new.

General Features
The state space is two-dimensional, being described by the variables (Ω, χ 1 , χ 2 ) subject to the restrictions (3.46). We use (3.46) to eliminate χ 1 via The state space is then the rectangle in the (Ω, χ 2 )-plane given by In accordance with theorem 3.1.1, there are two qualitatively different classes of 2-fluid FL cosmologies, depending on whether γ 2 > 2 3 or γ 2 < 2 3 (recall that γ 1 and γ 2 are assumed to satisfy 0 ≤ γ 2 < γ 1 ≤ 2). We will discuss these two classes in sections 3. An important difference between these two cases is that if γ 2 < 2 3 the Einstein static universe, represented by an equilibrium point E, is included in the family of models. The point E is given by 15  15 The restriction γ2 < 2 3 ensures that χ2,E < 1, as required. γ 2 < 2 3 thenq is negative on a subset of the state space. Figure 3.2 shows that the orbits of ever-expanding models eventually enter this subset, and as a result the models are accelerating at late times.
The portraits of the orbits are drawn by first sketching the various invariant sets and equilibrium points, as given in section 3.1.5, in particular the boundary of the state space (3.154) and the flat FL invariant set Ω = ± π 4 . In the case γ 2 < 2 3 , the orbits that are past and future asymptotic to the Einstein static point E should be drawn. At this stage, the orbits that have been drawn represent exceptional models and form a "skeleton" for the state space. The orbits of typical models can now be drawn by noting the past and future attractors as given in theorem 3.1.3. In drawing the portraits it is helpful to note that the evolution equations (3.63) and (3.64) are invariant under the interchange i.e. changing the sign ofΩ reverses the direction of time. In particular, the invariant set − π 2 ≤Ω ≤ − π 4 represents cosmologies that evolve in the same way as those with π 4 ≤Ω ≤ π 2 , with the difference that the direction of time is reversed. These models are contracting throughout their evolution (Ω < 0 implies H < 0) and hence are not potential models of the real universe.
In practice the state space can be sketched numerically 16 using the conserved quantity K that is defined by (3.143) for the family of 2-fluid cosmologies. On substituting (3.153) into (3.155) and choosing n = 2, we obtain (3.157) The orbits are given by The value of K at E is denoted K E and is given by It follows that the equation describes the orbits that are past and future asymptotic to E, i.e. it describes the stable and unstable manifolds of E.
The value of K determines whether a model is ever-expanding or not, as shown in the following table. 16 For example, using the plots[implicitplot] command in Maple.
Restriction on K Qualitative Behaviour Further, the sign of K determines whether the model is open, flat or closed. We note that all generic classes of FL models are covered by the given values of K.
In the case γ 2 > 2 3 the critical value of K is K crit = 0, indicating that the flat FL models play an exceptional role as regards to future evolution. In the case γ 2 < 2 3 , K crit = K E and the exceptional models are those that are past or future asymptotic to the Einstein static universe. This distinction is clearly indicated by comparing figure 3.1 (RDC-universes) with figure 3.2 (DCΛ-universes).
We note that the conserved quantity K can also be expressed in terms of the density parameters. It follows from (3.150) with n = 2 that with

FL Cosmologies with Radiation and Dust (RDC-universes)
We now consider 2-fluid FL cosmologies containing radiation (γ 1 = 4 3 ) and dust (γ 2 = 1). The portrait of the orbits, shown in figure 3.1 is drawn as described in the previous subsection.
For the RDC-universes, the expression (3.157) for the conserved quantity K becomes Equivalently, the expression (3.161) gives K in terms of the density parameters as The value of K is infinite or zero on the skeleton of the state space: Note that K is indeterminate at the past attractor A − = F + r (χ d = 0,Ω = π 4 ) and at the future attractor . This indeterminacy arises since infinitely many orbits meet at A − and A + .
Discussion: Typical models differ qualitatively depending on the value of K, which determines where the orbit lies in the state space. Firstly, for ever-expanding models if the orbit passes close to M + r , the dust component will never be significant dynamically (i.e. Ω d ≪ 1) and this will correspond to K being large and negative. On the other hand if the orbit passes close to F + d then Ω d will attain a maximum value close to 1, corresponding to K being negative and close to zero. Secondly, recollapsing models differ quantitatively by how close the orbit comes to F ± d , i.e. by the maximum value of χ d , which in turn will depend on K. One can determine the dependence on K explicitly, as follows: i) For ever-expanding models, we obtain ii) For recollapsing models an expression for χ r | min and χ d | max can be obtained by settingΩ = 0 in (3.163) and solving the resulting second degree polynomial. We obtain

FL Cosmologies with Dust and Cosmological Constant
We now consider 2-fluid FL cosmologies containing radiation (γ 1 = 4 3 ) and dust (γ 2 = 1). The portrait of the orbits, shown in figure 3.2 is drawn as described in section 3.3.1.
For the DCΛ-universes, the expression (3.157) for the conserved quantity K becomes (3.167) Equivalently, the expression (3.161) gives K in terms of the density parameters as The critical value of χ Λ , i.e. the value of χ Λ at the Einstein static equilibrium point, is given by (3.155) as It then follows from (3.159) that the critical value of K is The value of K is infinite or zero on the skeleton of the state space: Note that K is indeterminate at the past attractor and at the future Discussion: Firstly, ever-expanding models differ qualitatively depending on the magnitude and the sign of the curvature scalar Ω k . If the orbit stays close to flat FL (Ω = π 4 ), the curvature will never be significant dynamically (i.e. Ω k ≪ 1) and this will correspond to K being close to zero. Secondly, recollapsing models (and bouncing models) differ qualitatively by how close the orbit comes to the Einstein static model, i.e. by the maximal (or minimal) value of χ Λ . Since χ Λ < χ Λ,E in recollapsing models (see (3.155)), the cosmological constant will never be the dominant fluid. Similarly, since χ Λ > χ Λ,E in bouncing models, there will be no dust dominant epoch, i.e. χ d ≈ 1. One can determine the dependence on K explicitly, as follows: i) For ever-expanding models, we obtain ii) For recollapsing and bouncing models, an expression for χ Λ | max and χ Λ | min can be obtained by settingΩ = 0 in (3.167) and solving the resulting third degree polynomial.

3-Fluid FL Cosmologies -Qualitative Analysis
We now use the results of sections 3.1-3.3 to give a qualitative analysis of the dynamics of 3-fluid FL cosmologies in which the matter content is radiation, dust and a cosmological constant (Λ > 0) (labelling the χ i as χ r , χ d and χ Λ ). The analysis given in this section is new.

General Features
The state space is three-dimensional, being described by the variables (Ω, χ r , χ d , χ Λ ) subject to the restriction (3.46). We use (3.46) to eliminate χ r via The physical state space is then a solid triangular prism in R 3 given by We can construct a "skeleton" for the state space using various invariant subsets, as shown in figure   3.4. The choice of χ i = 0 for any one fluid gives a 2-fluid FL cosmology, identified as follows: Further, the orbits in the planesΩ = ± π 4 represent the flat 3-fluid universes. The Einstein static universe is represented by a line segment of equilibrium points denoted E r E d , which is the intersection of the planesΩ = 0 andq = 0. The stable manifold S + of the set of Einstein static points is a surface through the line segment E r E d terminating at deS − when χ Λ = 1 and at the line This surface can be described by an equation of the form We know from theorem 3.1.1 that there are three generic classes of RDCΛ FL universes. Their orbits are contained in invariant subsets bounded by the surfaces S ± , as follows: Ever-expanding universesΩ > |f (χ d , χ Λ )| (the subset to the right of S + and S − ), Recollapsing universes (the subset between S + and S − and below the planeq = 0), Bouncing universes (the subset between S + and S − and above the planeq = 0).
The planesΩ = 0 andq = 0, although not invariant sets, play an important role in describing the dynamics. The planeΩ = 0 divides the state space into a region of expansion (Ω > 0) and a region of contraction (Ω < 0). The variableq is given by (3.65) which here specializes tõ On using (3.172) to eliminate χ r , we find thatq = 0 on the plane Equation (3.177) then divides the state space into a region of acceleration (q < 0) and a region of deceleration (q > 0) as follows: It follows from the general discussion in section 3.2.3 that there are two conserved quantities for the RDCΛ FL universes, defined by The orbits of typical models are described by the equations i.e. the orbits are the intersections of these two families of surfaces in the state space. As a special case the flat universes satisfy K = 0 (i.e.Ω = π 4 ) and then the equation describes the orbits in the invariant setΩ = π 4 and can be used to numerically sketch them, as in figure 3.5.

Discussion:
Knowing the past attractor A − = F + r and the future attractor A + = F − r ∪ deS + one can visualize the orbits of typical universes. For example, the orbits of ever-expanding models join F + r to deS + . The models will differ qualitatively depending on whether they come close to F + d (Ω d ≈ 1 in a neighbourhood of F + d ) and whether they come close to the Milne set (Ω = π 2 ) or to the surfaces S ± , in which case there will be an epoch in which spatial curvature is dominant (i.e. |Ω k | ≈ 1). One can determine the points in state space at which |Ω k | and χ d attain a maximum value, as follows. If K = 0, it follows from (3.36) and (3.69) thatΩ and |Ω k | are simultaneously extremal. Since in an ever-expanding modelΩ is extremal whenq = 0, we conclude that |Ω k | is extremal on the plane given by (3.177). The matter variable χ d will attain a maximum when χ ′ d = 0, which impliesq = 1 2 , as follows from (3.64) and (3.65). After substituting (3.172) to eliminate the χ r term in (3.176), we conclude that χ d will attain a maximum on the plane χ d + 4χ Λ = 1. (3.184)

A Model of the Physical Universe
In this section we discuss the RDCΛ FL universes as viable models for the large scale dynamical behaviour of the physical universe. Current observations of distant galaxies, type I a supernovae and the cosmic microwave background radiation support the following conclusions: i) The universe is expanding in the present epoch (H 0 > 0).
ii) The universe is accelerating in the present epoch (q 0 < 0).
iii) The early universe was radiation-dominated (hot big-bang).
iv) Pressure-free matter (baryonic matter and dark matter) is dynamically significant in the present epoch.
v) The spatial geometry is close to flatness in the present epoch.

Explicit FL Cosmologies
In We are primarily interested in expanding models, which begin at a big-bang singularity, i.e. the length scale satisfies ℓ(t s ) = 0 at some time t s , and dℓ dt > 0 in some interval. Since the DE (3.196) is autonomous, it is invariant under the change t → t + C. We will use this freedom to set t s = 0, thereby fixing the constant of integration that arises in solving (3.196). We will thus solve (3.196) subject to the requirement dℓ dt > 0 and the initial condition ℓ(0) = 0.
The solution ℓ(η) can in principle be expressed in terms of clock time t using the equation (D.1). As with (3.196), when solving (3.198) for ℓ(η), we will impose the requirement dℓ dη > 0 and the initial condition In the following sections we consider five classes of models in which either (3.196) or (3.198) can be transformed into (3.195). The solution of (3.195) subject to the requirement dL dT > 0 and the initial condition L(0) = 0 is given in Appendix C. In each case, we give the solution for ℓ(t) or ℓ(η) and note that the corresponding expressions for the energy densities and pressure are given by equation

Single Fluid FL Cosmologies
In this case we consider solutions to the Friedmann equation with n = 1, 2 3 < γ ≤ 2 and arbitrary curvature. We choose the conformal factor to be λ = λ 1 . The Friedmann DE (3.198) assumes the form We perform the change of variable  The result is The parameter λ is not essential, since by rescaling ℓ it can be assigned any value. Equation (3.209) thus represents the unique flat FL universe with a single fluid. 17 In the case of radiation (γ = 4 3 , β = 1) one can integrate (D.1) to obtain, In the case of dust (γ = 1, β = 2) one can integrate (D.1) in closed form, but cannot invert the expression in terms of elementary functions. In this case t and η are related by

FL Vacuum Cosmologies
In this case we solve the Friedmann equation where the only fluid component is a cosmological constant, i.e. n = 1, γ 1 = 0. We choose the conformal factor to be λ = λ 1 , which leads to (3.220)

Flat FL with a Single Fluid and Λ
In this case we solve the Friedmann equation with n = 2, γ 1 = γ arbitrary, γ 2 = 0 and zero curvature (i.e. K = 0). We choose the conformal factor to be λ = λ 2 , which leads to (3.213). As described in section and so . (3.225) Using (3.114), the energy density is given by These solutions were first given by Harrison (1967) in equation (23). They represent ever-expanding models and correspond to the orbit F + → deS + along the K = 0 invariant set, as in figure 3.3.

FL Cosmologies with Radiation and Dust (RDC-Universes)
In this case we solve the Friedmann equation with n = 2, γ 1 = 4 3 and γ 2 = 1. We choose λ = λ 2 as the conformal parameter. Then, by (3.113), the single mass parameter is given by (3.227) The Friedmann equation (3.198) assumes the form We apply the change of variable which gives dL dT and so where S K is given by (3.205). Using (3.114), the energy density is given by In this case, one can integrate (D.1) so as to obtain an expression for the clock time t = t(η), but cannot invert this expression to obtain ℓ = ℓ(t): (3.234) The solution (3.232) contains the one-fluid dust and radiation solutions given in (3.206) as special cases.
One can see this by writing (3.232) in terms of parameters λ r = λ 1 and λ d = λ 2 , and letting λ r → ∞ or λ d → ∞.
These solutions were first given 18 by Chernin (1966). They describe all orbits in the state space for RDC-universes with H > 0 in some epoch (see figure 3.1) and hence describe two qualitatively different classes of models: ever-expanding models if K = −1 (the orbits F + 1 → M + ) or K = 0 (the orbit F + 1 → F + 2 ) and recollapsing models if K = +1 (the orbits F + 1 → F − 1 ). The dimensionless parameter m determines key physical properties of the models. First, it determines the epoch of matter-radiation

FL Cosmologies with Radiation and Λ (RCΛ-Universes)
In this case we solve the Friedmann equation with n = 2, γ 1 = 4 3 , γ 2 = 0 and arbitrary curvature. We choose the conformal factor to be λ = λ 2 so that (3.213) holds. Then, by (3.113), the single mass parameter is given by (3.238) The Friedmann equation (3.196) then takes the form Upon making a change of variable according to we obtain dL dT (3.241) 18 Harrison reports that in the case K = +1 this solution was given by Lemaître (1927), but we have not verified this reference.
and so Using (3.114), the energy density is given by  This result follows upon rewriting (3.243) using standard hyperbolic identities as We note that at the critical value of m, denoted m * and given by m * = 1 4 , the solution (3.243) simplifies to This solution is future asymptotic to the Einstein static solution and corresponds to the orbit F + → E in figure 3.2.

Explicit Solutions in Terms of Observational Parameters
In sections 3.5.1-3.5.5, we have presented the explicit solutions using intrinsic parameters, namely, a conformal factor λ, the curvature indicator K and, for the 2-fluid solutions with non-zero curvature, a dimensionless mass parameter m, which describes the relative significance of the two fluids.
In order to link with observations it is desirable to write the solutions in terms of the observational parameters. In order to do so we need to use the formula that relates the intrinsic parameters λ i to the observational parameters, derived in section 3.2.2. In the case of K = 0, we can use equations (3.123) and (3.137) to write (3.250) One can think of (3.250) as relating λ −1 i , the length scale determined by the i th fluid, to ℓ 0 , the length scale determined by the spatial curvature. In particular, we will need the following special cases of (3.250): Note that the Friedmann equation (3.36), evaluated at t 0 gives where Ω 0 is the total density parameter evaluated at t 0 . In the case of a cosmological constant (i.e. γ = 0), we can obtain an expression for λ Λ directly from (3.137), as follows: In the case of K = 0 and n > 1, we use the freedom to scale ℓ in order to set λ 1 equal to λ n . Then, as a result of (3.137), the length scale ℓ 0 is determined explicitly in terms of Ω i,0 and H 0 according to The intrinsic parameters λ i are then obtained via (3.137).
There is no work to be done in writing the densities in observational form, since they are given in general by equation (3.125), which we repeat for convenience: where a = ℓ/ℓ 0 .
In order to illustrate this process we give the length scale function for the RDC-universes in terms of the observational parameters. We first rewrite ℓ(η) in (3.232) in terms of λ r and λ d : Then (3.251) and (3.252) immediately give where Ω 0 = Ω r,0 + Ω d,0 . Then Ω r,0 = 0 gives the DC-universe and Ω d,0 = 0 gives the RC-universe.

Peacock (1999) refers to the terms
, and Ω r,0 2|1 − Ω 0 | as the "dimensionless masses" (p79), when discussing RDC-universes. Note that these quantities are closely related to those given in (3.251) and (3.252), upon applying (3.253). These quantities also appear in the DC-and RC-universes as given by Coles and Lucchin (1995) (see p39, eq (2.4.2) for dust in terms of η and p41, eq (2.5.2) for radiation in terms of t).

Discussion
We summarize the principal features of the dynamics of n-fluid FL universes. i) n-fluid FL universes have at most one epoch of expansion (H > 0) and at most one epoch of acceleration (q < 0, or equivalently, w < − 1 3 ). These results depend crucially on the assumption that the fluids are non-interacting. We refer to Tolman (1934), p402 and p429-431, and to Clifton and Barrow (2007) for a discussion of models in which this assumption is not made.
iii) The n-fluid FL cosmologies admit n − 1 conserved quantities formed from the density parameters Ω 1 , . . . , Ω n . In particular, they admit one curvature quantity and n − 2 matter quantities.
We now comment on an important difference between n-fluid FL cosmologies with n > 2 and 2fluid cosmologies. For a 2-fluid FL model, there are two essential parameters, a conformal parameter λ and a mass parameter m 1 , or equivalently, two observational parameters Ω 1,0 and Ω 2,0 . There is also one conserved quantity K that keeps track of the spatial curvature. In an RDCΛ-model of the physical universe, the maximum value of χ dust is close to unity, which ensures that there is an extended matter-dominated epoch.
The dynamical systems analysis of the n-fluid FL universes leads to an interesting conclusion regarding the behaviour of the spatial curvature, as follows. For any ǫ > 0, there is an open set of n-fluid FL universes with γ n < 2 3 whose spatial curvature scalar Ω k satisfies |Ω k | < ǫ throughout the evolution. In other words, if γ n < 2 3 , in particular, if fluid n is a cosmological constant, the total matter-energy content "controls" the spatial curvature throughout the evolution. This result is a consequence of the fact that generic ever-expanding models satisfy Ω k → 0 as ℓ → 0 and since γ n < 2 3 , Ω k → 0 as ℓ → +∞. It follows that |Ω k | attains a maximum value ǫ which depends on the orbit in question. By restricting to orbits that lie within a sufficiently small neighbourhood of the flat FL invariant set Ω k = 0, the value ǫ can be made arbitrarily small. The above result has been found using different methods by Lake (2005), for DCΛ FL universes. He indicates that the result will also hold for RDCΛ universes (see also Lake (2006)). Lake comments on the significance of this result in connection with the so-called flatness problem (Lake (2005)).
One of the main results of this chapter is the analysis of the RDCΛ FL universes and the comparisons with observations. Discussion of these models dates back to the 1930's. Our analysis is novel in that it gives new insight into the dynamics of these models, by representing the evolution as orbits in a compact state space. We now give a brief survey of previous work on these models. The work that is closest to ours is that of Ehlers and Rinder (1989), who give a dynamical systems analysis of RDCΛ-universes, using Hubble-normalized variables. They give the explicit expressions for the density parameters as functions of the scale factor [equation (19)] and also show the existence of two conserved quantities [equation (29)]. Their state space representation is inevitably incomplete since the Hubble-normalized variables are undefined at the instant of maximum expansion in recollapsing and bouncing models.

PERTURBATIONS OF FL MODELS
The high degree of isotropy of the cosmic microwave background (CMB), when combined with the Copernican Principle, provides strong support for the belief that the large scale structure of the observable universe is very well described by the Friedmann-Lemaître universes, which are exactly isotropic and spatially homogeneous. The real universe is, of course, not exactly isotropic and spatially homogeneous, since there is complex structure associated with the observed distribution of galaxies, the overall expansion may not be exactly isotropic, and there may be primordial gravitational waves. But the current belief of most cosmologists is that the deviations from an exact FL cosmological model are sufficiently small that they can be described by considering linear perturbations of the FL models.

Historical Development
In this section we briefly discuss the history of the theory of linear perturbations of FL cosmologies. There are two main approaches, which we shall refer to as the metric approach and the geometrical approach.

Metric Approach
The metric approach to perturbations of the FL models was introduced by Lifshitz (1946) and subsequently discussed in greater detail by Lifshitz and Khalatnikov (1963). In the metric approach, one distinguishes a background spacetime, which is a FL model, and a perturbed spacetime, which represents the physical universe. The metric for the perturbed spacetime is written in the form where g (0)ab is the metric of the FL model and δg ab is a small perturbation, i.e. |δg ab | ≪ |g (0)ab |. The fluid 4-velocity and energy density in the perturbed spacetime are likewise written in the form The density contrast is defined by

The problem of gauge
In defining a perturbation, one is effectively establishing a one-to-one correspondence between points in the physical spacetime and points in the background spacetime. However, this correspondence is not unique: one may make an infinitesimal change in the correspondence keeping the background coordinates The most important ambiguity concerns the density contrast δ, as defined by (4.4). Bardeen defines two gauge-invariant quantities in terms of δ, denoted by ǫ m and ǫ g . The first equals the density contrast δ when one uses the so-called co-moving gauge, while the second equals δ in the so-called Newtonian gauge (see Bardeen 1980, p22).
ii) Second, as emphasized by Hawking (1966), the metric tensor is not a physically significant quantity since one cannot measure it directly.

Geometrical Approach
The geometrical approach has its origins in a paper by Hawking (1966). Instead of using the components of a perturbed metric as basic variables, he proposed using the evolution equations for the kinematic quantities, the Weyl curvature and the matter density to study how perturbations evolve. Although this paper initiated a significant new approach to studying linear perturbations of FL, the analysis of density perturbations was flawed and the paper had little impact. Subsequently, motivated by Hawking's paper, Lyth and Mukherjee (1988) used the evolution equations for µ and H to study density perturbations in an FL cosmology, although they did not do so in a fully gauge-invariant way. The decisive step was taken by , when they proposed using the spatial gradient of the matter density (∇ a µ) as the basic variable to describe density perturbations. Unlike the density contrast δ, the spatial gradient is a In this chapter we give an introduction to the geometrical approach to linear perturbations of FL models.

The Linearized Einstein Field Equations
The geometrical approach is based on the conviction that one should use variables that have a direct physical or geometric meaning, and are both coordinate-independent and gauge-invariant. One thus uses tensorial quantities that are identically zero in an FL universe, and hence gauge-invariant, on account of the Stewart-Walker lemma.
We shall see that a complete set of tensorial quantities is the following: i) The density and Hubble gradients:∇ a µ,∇ a H ii) The fluid kinematic quantities: σ ab , ω a iii) The Weyl curvature: E ab , H ab These variables have two additional advantages: a) they directly describe the deviation of the physical universe from an idealized FL universe, and b) they are, in principle, observable, when normalized to be dimensionless using a power of H. We refer in particular to the fundamental paper of Sachs and Wolfe (1967), who analyzed the observations of distant galaxies from this point of view, and Maartens et al. (1995b), who used the cosmic microwave background observations to bound these quantities.
In the geometrical approach one thus replaces the EFEs for the metric tensor components by the equivalent system of evolution equations and constraints, satisfied by the tensorial quantities i)-iii), as given in Appendix A.
For simplicity we will assume the matter content consists of a perfect fluid with barotropic equation of state p = p(µ) and a cosmological constant Λ. We will use the standard notation, where c s is the speed of sound.

The Linearization Process
In order to derive the governing equations for linear perturbations of a FL cosmology, we must distinguish between quantities which are non-zero in the background (referred to as zero-order quantities) and those which are zero in the background (referred to as first-order quantities We note that the accelerationu a is not an independent quantity, since it can be expressed in terms of the density gradient X a via (A.17), on using (4.5) and (4.8): (4.10) The evolution equations for the zero-order quantities are simply the evolution equations for a FL cosmology with a single fluid and cosmological constant, given by (3.5) and (3.6), with µ replaced by µ + Λ and p by p − Λ to make the cosmological constant explicit. We also include the evolution equation for w, which follows from (4.12) and (4.5). terms of the zero-order quantities. We use (3.11) to write (3.8) as follows:

Zero-Order Evolution EquationsḢ
(4.14) The evolution equations for the first order quantities are obtained from the general system of evolution equations and constraints in Appendix A, by a process of linearization. There are three aspects to the linearization procedure: Drop all products of the first order quantities (4.7) and their spatial derivatives in the general nonlinear evolution equations.
L 2 ) The resulting equations form a coupled system of linear PDEs for the first order quantities, whose coefficients depend on the zero order quantities. For example, there are terms such as HX a , HE ab and µσ ab . The linearization process involves replacing the zeroth order quantity by its value in the background FL model in each of these products.
The resulting equations also contain the projected spatial covariant derivative operator∇ a acting on first order quantities. The linearization process involves replacing∇ a , which is defined on the physical spacetime, by the corresponding operator on the background FL spacetime.

The Linearized Evolution and Constraint Equations
We now derive the governing equations for linearized perturbations by applying the linearization procedure L 1 to the general system of evolution equations and constraint equations in Appendix A.
One aspect of the linearization process, namely dealing with the terms∇ aub , curl(u a ) and∇ a (∇ bu b ) involving the spatial derivative ofu a , requires particular attention. First, it follows from (4.10) and L 1 that 3∇ where F (µ) is a function of µ whose specific form is unimportant. Then the term F (µ)XaX b is dropped.
Second, it follows in a similar manner from (4.10), (4.8) and L 1 that Finally, it again follows from (4.10) and repeated application of L 1 that  to eliminate terms involving spatial derivatives ofu a .

First-Order Constraint Equations∇
There are two additional constraints involving the curls of the spatial gradients X a and Z a . Indeed, using (4.10), (4.16) and (4.17) leads to curl(X a ) = 6(1 + w)µHω a .

Comments:
i) The evolution equations (4.20)-(4.25) and constraints (4.26)-(4.30) form a closed system of linear first order PDEs for the first order quantities (4.7). The coefficients µ, w, c 2 s and H are zero-order quantities and are evaluated in the background FL model. Note that if K = 0 the term K/ℓ 2 in (4.21) is given in terms of zero order quantities by the constraint (4.14).
ii) It can be shown that the system of evolution equations and constraints is consistent: that is, the time derivatives of each constraint equation (4.26)-(4.30) are identically satisfied as a consequence of the other equations. As a result, once a set of initial data is chosen to satisfy the constraint equations, these equations will hold at all later times (see Bruni et al. (1992), p42). Most authors, however, do not make use of the full system of evolution and constraint equations.

Harmonic Decomposition
Before proceeding, we must discuss the fundamental role of spatial harmonics in the geometrical theory.
The harmonic decomposition is analogous to a Fourier decomposition of a quantity that is defined along the background 3-spaces. Hence, harmonics are defined so as to be constant along flow lines (i.e., independent of proper time) and orthogonal to the fluid 4-velocity u. We refer the reader to Appendix F for a complete discussion of spatial harmonics in the geometrical approach.
The key idea behind the harmonic decomposition is that a scalar, vector or tracefree symmetric rank 2 tensor orthogonal to u a can be written as a linear combination of scalar, vector and tensor harmonics with purely time-dependent coefficients. In the literature it is tacitly assumed that the harmonics are complete in some appropriate space of tensor fields.
The expansions are as follows: ab . Using this decomposition, it is natural to consider the kinematic quantities and Weyl tensor components of the real universe to be the fundamental description of the linear perturbations (since they are each gauge-invariant, each of these components can be expanded using the harmonic basis functions).
The decoupling of the perturbation types is extremely important in the linearized theory since it allows us to consider scalar, vector and tensor perturbations independently. For each perturbation type, we will proceed as follows: First, we first define a basic variable. Second, we use the linearized evolution equations (4.11)-(4.30) to derive a second order differential equation which describes the evolution of the basic variable, which we shall call the governing DE. Finally, we use the linearized evolution and constraint equations to express all kinematic quantities and Weyl tensor components in terms of the basic variable, hence providing a complete solution to the linearized equations for each type of perturbation.

Scalar Perturbations
We first focus on scalar perturbations, which are defined as solutions of the linearized equations (4.11)- or equivalently, using (4.8), We also rescale Z a according to Z a = ℓZ a , i.e.
Z a = 3ℓ∇ a H.

The Governing DE
In this section we derive the governing DE for scalar perturbations, following Ellis et al. (1990) (also see Wainwright and Ellis (1997, p290-3)). Since we are considering scalar harmonics, we can write Expressing these equations in terms of the rescaled variables D a and Z a leads tȯ These equations can be combined to give a second order evolution equation for D a . On differentiating (4.40) and using (4.11)-(4.13), (4.14) and (4.41), we obtain We now expand each term in terms of scalar harmonics. The derivatives of D a along the fundamental congruence are obtained by differentiating (4.38) and notingQ (0,k) a = 0 (see (F.10)), as follows: (4. 45) In order to remove the Laplacian term from the evolution equation (4.42), we use the identity (F.13), which can be written as On applying (4.45) and (4.46) to (4.42), we obtain

Kinematic Quantities and Weyl Curvature
We now write the kinematic quantities and Weyl curvature in terms of the density gradient coefficients D (k) and the scalar harmonics Q (0) . Since we are considering scalar harmonics, Z a and σ ab can be expanded as We can use (4.40) to express Z (k) in terms of D (k) andḊ (k) , and since ω a = 0 it follows from (4.26) and (F.16) that Expanding E ab in terms of scalar harmonics via ab , (4.51) and using (F.16) leads to (4.52) The coefficients E (k) can then be expressed in terms of D (k) via the constraint (4.29). We thus obtain the following expressions for the harmonic coefficients of the non-zero Hubble-normalized first-order quantities: 4

Vector Perturbations
In this section we focus on vector perturbations, defined as solutions of the linearized equations (4.20)-

The Governing DE
In this section, we derive the governing DE for vector perturbations. On differentiating (4.56) along the fundamental congruence and using (4.23), we obtaiṅ We expand (4.57) in terms of the associated vector harmonics P (4.58) On substituting (4.58) into (4.57), we obtaiṅ This equation is the governing DE for vector perturbations, in terms of clock time t.

Kinematic Quantities and Weyl Curvature
We now write the kinematic quantities and Weyl curvature in terms of the co-moving vorticity coefficients W (k) and vector harmonics Q (1) . By the definition of vector perturbations, all quantities can be expanded in terms of vector harmonics and associated vector harmonics, as follows: ab .
We can now use (4.26), (4.29) and (4.30) to determine the coefficients σ (k) , E (k) and H (k) in terms of W (k) . We thus obtain the following expressions for the harmonic coefficients of the first-order quantities:

Tensor Perturbations
In this section we focus on tensor perturbations, defined as solutions of the linearized equations (4.11)-(4.30) that can be expanded in terms of the tensor harmonics Q (2) ab . Tensor perturbations are interpreted as gravitational waves, described by a coupling of the electric and magnetic Weyl curvature.
As regards to the choice of basic variable, we shall see that for tensor perturbations the only non-zero first order quantities are σ ab , E ab and H ab . On the grounds of mathematical simplicity we choose the comoving shear tensor as the basic variable, denoted X ab and defined by X ab = ℓσ ab . (4.69)

The Governing DE
In this section, we derive the governing DE for tensor perturbations. 5 It is a consequence of applying tensor harmonics that rank two tensors have zero spatial divergence, i.e. We now use (4.24) to calculate (ℓE ab )˙⊥, noting thatl = Hℓ. We use the constraint (4.28) to eliminate H ab and then apply (E.31) to the "curl curl X ab " term. The result is We can now combine these equations to give a second order DE for X ab (differentiate (4.72) using (4.11) and (4.73), and then use (4.72) to eliminate E ab from the resulting equation). After expressing the coefficients of X ab in terms of the density parameters Ω and Ω Λ (see 3.31), we obtain 6 We now expand X ab in terms of tensor harmonics according to Substituting (4.75) into (4.74) and using (F.38) yields This is the governing equation for tensor perturbations, in terms of clock time t.

Kinematic Quantities and Weyl Curvature
We now write the kinematic quantities and Weyl curvature in terms of the dimensionless shear coefficients X (k) , tensor harmonics Q (2) ab and associated tensor harmonics P (2) ab . Since (4.70) confirms that σ ab , E ab and H ab can be expanded in terms of tensor harmonics, we can write 6 It should be noted that the evolution equation ( where ∆σ ab =σ ⊥ab −∇ 2 σ ab , so as to be consistent with the remainder of that paper. This equation is then equivalent to (4.74) for Λ = 0 on recalling that X ab = ℓσ ab and Θ = 3H.
where 7 P ab = ℓcurl(Q ab ). (4.78) Then using (4.69), (4.22) and (4.28) we can write σ ab , E ab and H ab in terms of X (k) . We thus obtain the following expressions for the harmonic coefficients of the non-zero first-order quantities:

Dynamics of the Linear Perturbations: General Features
In this section we discuss the dynamics of the linear perturbations. We first discuss the importance of conformal time in analyzing perturbations of FL and rewrite the governing DEs for each perturbation type in terms of this time variable. Finally, we examine the dependence of the first-order quantities on the solutions to the governing DEs.

Conformal Time and the Particle Horizon
It is a consequence of propositions 4.4.1, 4.5.1 and 4.6.1 that the spatial gradients, kinematic quantities and Weyl tensor components can be determined directly from the set of solutions of the governing DEs.
As a consequence, we turn our attention to the behaviour of solutions to these equations. 7 See definition in section F.3.
We make the choice of conformal time η as the independent quantity, defined in terms of the length whereg αβ is a 3-metric of constant curvature.
Conformal time is directly linked to the notion of the particle horizon, defined as the boundary of the part of the universe that is visible to us. In terms of conformal time, the distance to the particle horizon (see figure 4.1) is given by Since the wavelength λ for a perturbation of mode k at some fixed time η is given by Hence, if kη ≪ 1 the wavelength of the perturbation mode is large compared to the distance to the particle horizon (since λ ≫ d H ). In this situation it is customary to say that the perturbation is outside the horizon.

Governing DEs Using Conformal Time
We now present the governing DEs for scalar, vector and tensor perturbations in terms of conformal time and specialized to a γ-law equation of state, i.e. p = (γ − 1)µ. It follows from (4.5) that w = c 2 s = γ − 1. (4.87) The three governing DEs, given by (4.47), (4.59) and (4.76) in terms of clock time t, assume the following form in terms of conformal time. The change of variable can be performed simply by using these identities: which follow from (4.82).

Governing DEs for Perturbations using Conformal Time
Scalar In the case of vector perturbations, the governing DE (4.90) can be integrated directly, on applying 8 We then obtain W (k) = C(λℓ) 3γ−4 . In the case of scalar and tensor perturbations, the second-order governing DEs can be transformed into normal form by making a change of dependent variable. In both cases the resulting DE 9 has the following form, in terms of a constant frequency ω and a potential U (η): as in Bardeen (1980).

Scalar
Change of variable: Governing DE: with

General Solutions for Scalar and Tensor Perturbations with Λ = 0
In the case Λ = 0, the potential U (η) in the governing DEs is given by for both scalar and tensor perturbations (see equations (4.99) and (4.103)). The background solution is for K = (0, −1, +1). On noting that ΩH 2 ℓ 2 = 1 3 µℓ 2 , it follows from (3.207) that The value of ω, as given by where a = (c s k, k) (4.110) for scalar and tensor perturbations, respectively.

Functional Dependence of First-Order Quantities
We are interested in describing to what extent a perturbed FL cosmology deviates from an exact FL cosmology.
It is known that cosmological observations place bounds on the first order quantities, specifically on dimensionless quantities formed by dividing the first order quantities by an appropriate power of the The first analysis of this type was due to Sachs and Wolfe (1967) who showed that observations of distant galaxies could, in principle, place bounds on these quantities. More recently, Maartens et al (1995a,b) have shown that the high isotropy of the CMBR ( ∆T T ≈ 10 −5 ) places strong bounds on these quantities. 10 It is natural to use these Hubble-normalized quantities as a measure of how close the physical model is to an idealized FL model, and we will thus calculate their time dependence for the three perturbation types. A growing mode is an indication of an instability in the background FL model.
We now briefly mention the connection between the geometrical approach and Bardeen's gaugeinvariant metric approach to cosmological perturbations (Bardeen (1980)). In Bardeen's approach, the basic variables are gauge-invariant linear combinations of the metric perturbations. However, it turns out that Bardeen's metric potentials Φ H , Ψ and H (2) T have the same time dependence as certain expressions 10 The analysis is complicated and subject to assumptions being imposed on the time derivatives of the multipoles of the CMBR.
formed from the basic variables in the geometrical approach. It can be shown that 11

Perturbations of FL Models with Pressure-Free Matter
For the case where the background FL model contains dust and possibly a cosmological constant, the analysis of scalar perturbations simplifies significantly since one can avoid performing the harmonic decomposition. The evolution equation (4.42) for the fractional density gradient D a simplifies tö Namely, we use the fact that the fluid 4-velocity is a geodesic in the case of dust. Thus, we can drop the spatial projections from (4.120) and so obtain Although it follows that D a does not depend on the perturbation mode k, a harmonic decomposition is required in order to obtain expressions for the shear and electric Weyl tensor (see, for example, (4.54) and (4.55)).

The General Solution in Integral Form for Scalar Perturbations
A second simplification is that the evolution equation (4.122) has a first integral, which we derive as follows, using the spatial gradients D a , Z a and R a .
In the case of dust, the evolution equations, given in (4.40), (4.41) and (B.13) simplify significantly, and can be written asḊ Note that these DEs are not independent, since the 3 spatial gradients are linearly dependent. This result can be seen upon taking the spatial gradient of (1.43) and linearizing so as to obtain Observe that the quantityC a , defined by 12  On noting that, in the case of dust, (4.11) can be rewritten aṡ   If H → 0 or H → ∞ as η → η * (where η * = ±∞ is allowed), we can rewrite D (+) as On applying L'Hôpital's rule, in the limit as η → η * , we obtain Then (4.130) gives an expression for the asymptotic behaviour of the growing mode in terms of zero-order quantities: Although the general solution (4.133) cannot always be obtained in closed form, one can analyze the asymptotic behaviour of these solutions. We are primarily interested in three asymptotic epochs, namely flat FL (η ≈ 0), Milne (η ≫ 1, Ω k ≈ 1) and de Sitter (η ≈ η f , Ω Λ ≈ 1). The general behaviour of perturbations within these epochs shall be discussed in detail in section 4.10.

Scalar Perturbations with Λ = 0
We now use the general solution (4.135) to obtain the solution for a dust background with Λ = 0. The length scale (4.105), when specialized to the case of dust, is given by λℓ = (S K (η/2)) 2 , (4.139) and the Hubble parameter by H = ℓ ′ /ℓ 2 and (4.139), where C K (T ) is defined analogous to S K (T ) to be for K = (0, −1, +1).
We now consider the case of K = 0 separately. In this case, the growing mode can be determined from (4.135), (4.139) and (4.140), on integrating This integral is elementary and can be expressed in terms of H and ℓ after a short calculation: It follows from (4.107) with β = 2 that which allows us to rewrite the Friedmann equation (3.36) as We can now eliminate the H 2 ℓ 2 term in (4.143), obtaining The solution of (4.122) for K = 0 is then where H and ℓ are given by (4.139) and (4.140).
In the case of K = 0, the integral (4.142) evaluates to where the coefficient 1 10 is chosen so as to be asymptotic to the growing mode solution for K = 0 when η → 0. The solution for a flat FL background is then We substitute (4.152) into (4.151) and use T as the independent variable. After applying a standard identity, we obtain a .
The growing mode defined in this way will satisfy lim t→0 D (+) = 0.  where and n 2 = k 2 + 3K. Here, the coefficients c + and c − are constants that depend on k.
Tensor perturbations describe the propagation of gravitational waves, which are represented by the electric and magnetic Weyl tensor, E ab and H ab . We specialize the solution (4.160) to K = −1 and use (3.206) and H = ℓ ′ /ℓ 2 to determine ℓ and H. On applying (4.80) and (4.81), we can determine the asymptotic behaviour of E ab and H ab in three asymptotic regimes: The first regime (η ≪ 1 and kη ≪ 1) describes long wavelength perturbations that occur outside the particle horizon close to the singularity; the second regime (η ≪ 1 and kη ≫ 1) describes perturbations within the particle horizon close to the singularity; lastly, the long term behaviour of the perturbation, when spatial curvature is significant, is described in the third regime (η ≫ 1).
η ≪ 1 and kη ≪ 1 η ≪ 1 and kη ≫ 1 Here c ± , A and φ are arbitrary constants that depend on k. In the first regime we can see that the tensor perturbations do not represent gravitational waves since they have a power law dependence on η.
We depict the behaviour of E (k) /H 2 in figure 4.4, for three given choices of n.

Perturbations of Radiation-Filled FL Models
We  The coefficients c + and c − are constants that depend on k.
The solution for X (k) (η) is obtained by choosing β = 1 in the general solution (4.109), (4.111) and (4.113). Since the index has the value µ = 3 2 , these solutions are elementary functions, given by (H.14)-(H. 19). Upon comparing the functional dependence in each case, we obtain the following unified form of the solution: and n 2 = k 2 + 3K. The coefficients c + and c − are constants that depend on k.
As with perturbations of a dust background, we specialize the solution (4.163) to K = −1 and use (3.206) and H = ℓ ′ /ℓ 2 to determine ℓ and H. On applying (4.80) and (4.80), we can determine the asymptotic behaviour of E ab and H ab in three asymptotic regimes.
η ≪ 1 and kη ≪ 1 η ≪ 1 and kη ≫ 1 Here c ± , A and φ are arbitrary constants that depend on k. We depict the behaviour of E (k) /H 2 in

Perturbations in Asymptotic Epochs of FL Models
In this section we give approximate solutions of the governing ODEs for scalar and tensor perturbations in three asymptotic epochs for FL cosmologies. We also give the behaviour of the spatial gradients, kinematic quantities and Weyl tensor components for these approximate solutions. Further, to facilitate comparisons, we give the Bardeen (1980) gauge-invariant metric potentials for the approximate solutions.
The three asymptotic epochs are as follows: Flat FL. The flat FL background is described by This solution approximates the dynamics of any FL cosmology with an initial singularity with cosmological constant and arbitrary curvature at early times, since these models satisfy

The Flat FL Asymptotic Epoch
The solutions of the governing equations in the flat FL background are obtained by substituting (4.171) into (4.109) and (4.93) and are given by In the expressions for D (k) (η) and X (k) (η) we can identify two asymptotic regimes, namely kη ≪ 1 and kη ≫ 1. The restriction kη ≪ 1 means that the perturbation mode corresponding to k is outside the particle horizon, while kη ≫ 1 means the perturbation mode is inside the particle horizon (see sec- i.e. the spatial curvature is dynamically negligible in the regime kη ≪ 1.
The approximate solutions and the resulting time dependence of the first order quantities are given in

The Milne Asymptotic Epoch
The governing DEs in the Milne background are obtained by substituting (4.173) into (4.97), (4.90) and (4.101) and are given byD where holds for scalar perturbations and holds for tensor perturbations. For scalar perturbations there are two cases, namely ω 2 > 0 and ω 2 ≤ 0.
The solutions are if γ = 1, (4.187) where A is an arbitrary constant amplitude.
In order to calculate the approximate expression for E ab we need an approximate expression for Ω subject to η ≫ 1. A simple calculation using (3.206) and (3.207) yields the approximate asymptotic behaviour Ω = exp − 2 β η ≪ 1.

Perturbations of Milne.
Quantity Scalar Perturbations Vector Perturbations Tensor Perturbations
In order to calculate the approximate expression for E ab we need an approximate expression for Ω

Discussion
We now summarize the principal features of linear perturbations of FL cosmologies.
i) For each perturbation type the complete system of linearized evolution and constraint equations is reduced to a single ODE for the harmonic modes, together with expressions for the Hubble-  We contract (A.1) with h ab and u a and use (1.26) in order to obtain the following relations: (see Ellis (1973), p24). Using (1.30) we also obtain two important relations relating the Riemann curvature tensor and Weyl tensor components: The twice-contracted Bianchi identities (1.53) lead to the following evolution equations:  where T = −3H 2 + µ + Λ +∇ au a +u au a − 2σ 2 + 2ω 2 .
(A. 23) In order to proceed, we require an identity which permits the interchange of the˙derivative and the spatial derivative∇ a when applied to a scalar field f . A straightforward calculation using (1.18) yields for any scalar f (also see equation (92) in Ellis et al. (1990)).
We now derive evolution equations for X a and Z a , following Ellis and van Elst (1998, p56). Substi-

SPATIAL CURVATURE
In this appendix we introduce the 3-curvature when the fundamental congruence has non-zero vorticity.
This analysis leads to a generalization of the Gauss equation (1.37), the 3-Ricci tensor (1.39) and 3-Ricci scalar.

B.1 The Generalized Gauss Equation
We note that the projected covariant derivative∇ a , defined in (1.11), is a non-commutative operator on the space-time manifold when vorticity is present, i.e. ω a = 0. After a short calculation using (1.5), (1.11) and (1.18), we obtain∇ We now derive a similar identity for tensors with one or more index.
Upon performing a similar computation as in (B.1) and using (1.50), we obtain the expression Ellis et al. (1990), eq. (79)). We note that this definition reduces to the 3-Riemann tensor given in (1.37) when u a is irrotational. A similar identity to (B.4) can be obtained for tensors with two or more indices, generalized according tô On contracting (B.5), we can obtain an expression for the spatial Ricci scalar and trace-free spatial Ricci tensor in terms of the spacetime quantities from (B.5), incorporating the Weyl tensor through (1.30) and (1.32): Using the field equations (1.2) in the form (A.2)-(A.5) to eliminate R and R ab , we obtain S ab = E ab + 1 2 π ab − Hσ ab − Hω ab − 2 3 (σ 2 − ω 2 )h ab +σ c a σ cb + ω c a ω cb + ω c a σ cb + σ c a ω cb . (B.10)

B.2 The Spatial Gradient of the 3-Ricci Scalar
The spatial gradient of the 3-Ricci scalar is defined by R.
(B.11) An evolution equation for R a is obtained in the same manner as for D a and Z a , defined in (4.35) and   iii) 2-fluid FL plus arbitrary curvature. The equation of state parameters must satisfy where ∆γ is constant.
iv) 3-fluid FL. The equation of state parameters must satisfy where ∆γ is constant.

C.1 Solutions with a Stiff Fluid
There are three known cases where we can obtain a solution to the Friedmann equation where one fluid is a so-called stiff fluid, i.e. γ s = 2. In this section we briefly present each of these solutions.

SD-Universes.
We choose the conformal parameter to be λ = λ d , which, on using ( SRC-Universes. We choose the conformal parameter to be λ = λ r , which, on using (3.113), leads to a single mass parameter defined by (C.14) In this case the Friedmann equation reads We apply the change of variable Subject to the IC L(0) = 0, the solution is SDΛ-Universes. We choose the conformal parameter to be λ = λ Λ , which, on using (3.113), leads to two mass parameters defined by In this case the Friedmann equation reads We apply the change of variable L = λ 2 s λℓ 3 , T = 3λt, (C.21) Subject to the IC L(0) = 0, the solution is

C.2 Existence-Uniqueness for the Friedmann Equation
In this section we prove existence-uniqueness of solutions to the Friedmann DE (3.115), given by for any choice of parameters m i > 0, i = 1, . . . , n − 1 and subject to the conditions We can rewrite (C.24) as Since γ 1 > γ 2 > · · · > γ n and γ 1 > 2 3 , this DE is of the form Since dL dT > 0, we can rewrite (C.27) as The general solution to (C.29) is then, for all m i > 0, where G(L) is the antiderivative of the RHS of (C.29). The requirement that T = 0 when L = 0 fixes C, and hence gives a unique solution.
APPENDIX D
Proof. It follows from the definition of an ever-expanding model (see section 3.1.6) that Using the chain rule we expand (D.1) so as to obtain It then follows from (D.2) and (D.3) that η is a strictly increasing function of ℓ for t > 0.
Let η 0 and η * be the value of conformal time at ℓ 0 and ℓ * , respectively. We now integrate (D.3) over (η 0 , η * ) so as to obtain Since η is strictly increasing, we note that η * is bounded above if and only if F (ℓ * ) is bounded as ℓ * → +∞. Similarly, we note that η * is bounded below if and only if F (ℓ * ) is bounded as ℓ * → 0 + . We now prove the two results separately: (a) Assume that condition (a) is not met. It then follows from (3.196) that ∃ a constant C > 0 such that dℓ dt < C ∀ t < t 0 . The integral (D.4) then takes the form which is unbounded.
If condition (a) is met it then follows from (3.196) that ∃ constants C > 0 and ǫ > 0 such that dℓ dt > Cℓ −ǫ ∀ t < t 0 . The integral (D.4) then takes the form which is bounded.
(b) Assume that condition (b) is not met. It then follows from (3.196) that ∃ a constant C > 0 such that dℓ dt < C ∀ t > t 0 . The integral (D.4) then takes the form which is unbounded.
If condition (b) is met it then follows from (3.196) that ∃ constants C > 0 and ǫ > 0 such that dℓ dt > Cℓ ǫ ∀ t < t 0 . The integral (D.4) then takes the form which is bounded.
Thus, we conclude that the constraints hold, completing the proof.
We note that this theorem also holds in the case of closed FL as long as dℓ dt > C, where C is a constant satisfying C > 0. This constraint prevents models which recollapse or whose length scale tend to a constant as t → +∞.

COVARIANT DIFFERENTIAL IDENTITIES
In this appendix, we present a list of useful differential identities involving the derivative operators used in this thesis. This list represents an extension of the list given in Ellis et al. (1990) Appendix A by incorporating identities involving the Laplacian and curl operators.
Zero-order 3-curvature quantities for FL models: Before proceeding, we must first obtain zero order expressions for the background curvature. It quickly follows from (3.11), (1.40) and (1.47) that, to zero order, the background 3-curvature quantities satisfy the following relationships: We now present a list of linearized identities, derived using (1.7), (1.11), (1.16)-(1.17), (B.1) and (B.6). In all of the following relations, we shall assume that∇ a f , X a and T ab and all their derivatives are first order quantities. Further, we assume that T ab is symmetric and trace-free.

First-order time derivatives:
The following set of linearized identities are used when interchanging the derivative along the fundamental congruence (˙) with∇, curl(·) or∇ 2 : :

COVARIANT HARMONICS
In this appendix, we present an overview of general identities and results related to the covariant harmonics. We extend the results of Bruni et al. (1992) Appendix B by introducing the notion of associated vector and tensor harmonics.
As in Hawking (1966), we define the spatial harmonics Q Hereafter, as we are working in the geometrical formalism, the use of the term harmonics will refer to the covariant harmonics defined by (F.1) and (F.2). For the problem at hand, we only require scalar harmonics, vector harmonics and two-index tensor harmonics (which, for purposes of brevity will be referred to as tensor harmonics).
We briefly discuss notation for harmonics and quantities derived from harmonics. The 0 superscript of the harmonic quantity Q (0,k) denotes that it is a scalar harmonic. Similarly, we use 1 and 2 to denote vector and tensor harmonics, respectively. It is common in the literature to drop the wavenumber k when writing a harmonic quantity, so writing Q (0,k) as Q (0) instead. The wavenumber of the harmonic being used is instead implicitly defined by the equation where the harmonic appears. Further, here and elsewhere, we shall use Q to denote the whole set of harmonics and the standard quantities constructed from Q: we have scalar harmonics Q (0) , Q ab and tensor harmonics Q ab . We now give the definitions and properties of each of these three classes of harmonic types. it can be quickly shown that Q a and Q ab satisfy the curl-free constraint, curl(Q a ) = 0, curl(Q ab ) = 0, (F.8) congruence-orthogonality, u a Q a = 0, u a Q ab = 0, u b Q ab = 0, (F.9) and time invariance,Q a = 0,Q ab = 0. (F.10) We are restricted in our choice of k 2 to the eigenvalues determined by Lifshitz (1946), which in turn are determined by the spatial curvature according to k 2 = n 2 − K n > 0 if K = −1, 0, n = 3, 4, . . . if K = +1. (F.11) We note that it is common to include n = 2 as a possibly harmonic mode in the case of K = +1, but it can be shown via the constraint equation (4.29) and identity (F.16) that this mode is unphysical. In fact, this mode is gauge-dependent and associated with the freedom of choice of the velocity field u a .
Some important identities involving the vector quantity Q a can be derived using the differential identities in Appendix E and (F.3)-(F.5): ℓ∇ a Q a = kQ, (F.12) Similarly, some important identities involving the tensor quantity Q ab can be derived using the differential identities in appendix E and (F.3)-(F.6):

F.2 Vector Harmonics
Vector harmonics are defined as solutions of∇ where Q a is divergence-free and time-independent, ∇ a Q a = 0,Q a = 0. The tensor-valued vector harmonics are related using (F.25), (F.26) and (E.24), as follows: ℓ curl(P ab ) = 1 2 (k 2 + 2K)Q ab , (F. 30) ℓ curl(Q ab ) = 1 2 P ab . (F.31) As with scalar harmonics, we are restricted in our choice of k 2 to the eigenvalues determined by Lifshitz (1946), according to Again, it is common to include n = 2 as a possibly harmonic mode in the case of K = +1, but it can be shown via the constraint equation (A.21) that this mode is unphysical.
Using the differential identities in Appendix E and (F.21)-(F.23) we can derive the following identity for the vector harmonics Q a and associated vector harmonics P a (in this case, we may replace Q a with P a to obtain a second identity): Some identities involving the tensor-valued vector harmonics Q ab and associated tensor-valued vector harmonics P ab can be derived using the differential identities in appendix E and (F.21)-(F.26). In each of the following identities, the pair (Q ab , Q a ) may be replaced with (P ab , P a ):

F.3 Tensor Harmonics
Tensor harmonics are defined as solutions of where Q ab is divergence-free, trace-free and time independent: ∇ b Q ab = 0, Q a a = 0,Q ab = 0. (F.39) As with vector harmonics, every tensor harmonic Q ab defines an associated tensor harmonic P ab via P ab = ℓcurl(Q ab ).

(F.40)
It can be shown that P ab satisfies the definition of a tensor harmonic (F.38)-(F.39) using (E.24), (E. 29) and (E.13). It then follows immediately from (F.38) that the quantities Q ab and P ab satisfy congruence orthogonality: u a Q ab = 0, u b Q ab = 0, u a P ab = 0, u b P ab = 0. (F.41) Upon taking the curl of (F.40) and applying (F.38) and (E.31), we obtain Q ab in terms of P ab : Q ab = (k 2 + 3K) −1 ℓcurl(P ab ). (F.42) As with scalar and vector harmonics, we are restricted in our choice of k 2 to the eigenvalues determined by Lifshitz (1946): We can obtain an important identity for Q ab and P ab using (B.6) and (E.1), as follows (the identity for P ab is identical, except with Q ab replaced with P ab ): ℓ 2∇ [a∇b] Q cd = 1 2 K (h ac Q bd − h bc Q ad ) + (h ad Q bc − h bd Q ac ) .

THE 3-COTTON-YORK TENSOR
In this appendix we introduce the 3-Cotton-York tensor, which is useful in the analysis of tensor perturbations. Specifically, we are interested in using the 3-Cotton-York tensor in place of the co-moving shear X ab in developing a governing DE for tensor perturbations. It is of interest in this regard due to its relation with Bardeen's metric potential H T . We will further show that the coefficients of the shear and Weyl tensor components can be rewritten in terms of C (k) . To the best of the author's knowledge, this analysis of tensor perturbations in terms of the 3-Cotton-York tensor is new.
In terms of the trace-free spatial Ricci tensor (1.40), the 3-Cotton-York tensor is defined as (see, for example, van Elst and Uggla (1997, p2680-1)) (3) C ab = curl( (3) S ab ), (G.1) and hence is symmetric and tracefree. In the linearized context, we can use (B.10) and (4.28) to obtain an alternative first-order expression for (3) C ab in terms of the Weyl tensor components, namely C ab = −HH ab + curl(E ab ).

(G.2)
We note that an important property of this quantity is that it has zero spatial divergence in the linear theory.
Indeed, upon taking the divergence of (G.

G.1 Evolution Equation for the 3-Cotton-York Tensor
Instead of choosing the co-moving shear to be the basic quantity for tensor perturbations (as in section 4.6), we can instead choose the basic quantity to be the dimensionless 3-Cotton-York tensor, defined by Upon recalling from (G.3) that the 3-Cotton-York tensor has zero divergence, it quickly follows that C ab also has zero divergence, i.e.∇ b C ab = 0. (G.6) Hence, C ab can be expanded purely in terms of associated tensor harmonics, according to ab , (G.7) The derivatives of C ab along the fundamental congruence are obtained by differentiating (G.7), givinġ This is then the second-order DE governing the evolution of the 3-Cotton-York tensor. In terms of conformal time, this DE reads C ′′ (k) + 2HℓC ′ (k) + (k 2 + 2K)C (k) = 0, (G. 12) which is identical to Bardeen's DE for the metric potential H T , given in Bardeen (1980) equation (4.14).

G.2 Tensor Perturbations
We now show the relation between (G.11) and the governing DE for tensor perturbations (4.76). We note that it follows by (G.3) that the 3-Cotton-York tensor is purely tensorial, i.e. the scalar and vector contributions to (3) C (k) are exactly zero. Expanding (3) C (k) in terms of associated tensor harmonics yields

H.2 Special Solutions of the Governing DEs
For certain values of γ, µ and ν, the Bessel functions and associated Legendre P functions reduce to closed-form trigonmetric expressions. Further, in each of these cases we can derive a second linearly independent solutionQ of the associated Legendre DEs (H.3) and (H.5) that is defined so as to reduce to J µ (bz) in the limit z → 0. This solution is a preferred choice, since it emphasizes the relation between solutions for K = ±1 and solutions with K = 0 in the early time limit. In this section, we present the first few instances of this behaviour and its importance in modelling perturbations with a given fluid background.
For brevity, we shall make use of the notation C(z) = cos(z), S(z) = sin(z), (H.7) and use the tilde ∽ to denote dependence up to a multiplicative constant.

H.3 Resonant Solutions of the Legendre DE
For certain special cases, the associated Legendre functions P and Q are no longer linearly independent (see Abramowitz and Stegun p333, eq. (8.1.8)). One such example is the case of scalar perturbations with a dust background (γ = 1), where a second solution must be determined through alternative means. We

H.4 Asymptotic Behaviour
We now present the asymptotic behaviour of the Bessel functions and associated Legendre functions, both for early times and late times. as t → ∞.