Sami-Matias Niemi


Turku 2011

From the Department of Physics and Astronomy

University of Turku

Turku, Finland

Supervised by

Dr. Pekka Heinämäki and Dr. Pasi Nurmi

Department of Physics and Astronomy

Tuorla Observatory

University of Turku

Turku, Finland

Reviewed by

Prof. Heikki Salo

Division of Astronomy

Department of Physics

University of Oulu

Oulu, Finland


Dr. Antti Tamm

Tartu Observatory

Tõravere, Estonia


Prof. Volker Müller

Leibniz-Institut für Astrophysik

Potsdam (AIP)

Potsdam, Germany

One could not be a successful scientist without realizing that,
in contrast to the popular conception supported by newspapers
and mothers of scientists, a goodly number of scientists are
not only narrow-minded and dull, but also just stupid.

Dr. James Dewey Watson
Nobel Laureate


Countless warm and well deserved thanks are due to all the people who contributed, and more importantly, forced me to push through the hard times when quitting felt like the only right thing to do. These thanks are intended especially for the benefit of those who are actually not interested in the rest of the thesis, but are driven to read this section and this section only to find out if their names are listed. To their convenience, I have deliberately minimised the number of names spelled out.

First words of appreciation must go to my supervisors Pekka Heinämäki and Pasi Nurmi who exposed me to the mysteries of cosmological simulations, formation of large-scale structure, and cosmology. I may not have been a good student - I was hardly ever present during my years of studies - thus, warm thanks for sticking with me. I also want to express my gratitude to Mauri Valtonen for his ideas in galaxy group studies. I am in debt to Henry Ferguson and Rachel Somerville for guiding me through the last paper of the thesis. Your guidance truly convinced me that galaxy formation and evolution is worth spending countless long, but also fruitful, hours when others sleep.

For financial support, that made this work possible, I thank the following contributors: the Väisälä foundation, Nordic Optical Telescope Science Association, Space Telescope Science Institute, and G.J. Wulff foundation. Special thanks go to my mother Liisa Misukka for always providing extra funds when times weren’t so great.

I would like to thank the warm staff of Nordic Optical Telescope for making my stay in la Palma so enjoyable. I did not only enjoy eating pata negra and drinking wine (both of which I did more than I care to admit), but I also learned invaluable lessons about observational astronomy. For a hardcore theoretician this might sound odd, but astronomy is actually observation driven science and we should all make sure that we have at least seen the night sky. I would also like to thank all the individuals who made me feel welcome to STScI and taught me more about observations. Special thanks go to Stefano Casertano, Danny Lennon, Harry Ferguson, and Rachel Somerville for countless career advices and extremely useful science discussions. I would also like to thank Michael Wolfe for reading the early draft of this manuscript and standing for my Finglish. Last but most importantly I would like to thank Carolin Villforth for the countless science arguments and all the wonderful moments I have been able to share with you.

I thank the official reviewers of this thesis, Drs. Heikki Salo and Antti Tamm, for timely reading and useful criticism. I am greatly honoured that Prof. Volker Müller has agreed to act as my opponent in the public disputation of this dissertation. Finally, to those who wish to know the tips and tricks of accomplishing a PhD, I confidently say that all the people mentioned here are part of it and to them I owe my deepest gratitude.

Sami-Matias Niemi

1] ]1 ]

List of publications

I Are the nearby groups of galaxies gravitationally bound objects?
Niemi S.-M., Nurmi P., Heinämäki P. and Valtonen M.
MNRAS 382, 1864 (2007)
II The origin of redshift asymmetries: how CDM explains anomalous redshift
Niemi S.-M. and Valtonen M.
A&A 494, 857 (2009)
III Formation, evolution and properties of isolated field elliptical galaxies
Niemi S.-M., Heinämäki P., Nurmi P. and Saar E.
MNRAS 405, 477 (2010)
IV Physical properties of Herschel selected galaxies in a semi-analytical galaxy formation model
Niemi S.-M., Somerville R.S., Ferguson H.C., Huang K.-H., Lotz J. and Koekemoer A.M.
submitted to MNRAS

Chapter 1 Introduction

How did the Universe begin? How will it evolve? How did all the structures such as galaxies, groups, and clusters we observe in the night sky begin to form? How did they grow and how will they evolve? These are some of the most profound questions the mankind have and shall seek to answer.

Observations of the cosmic microwave background (CMB) radiation from the time when the Universe was only 111Cold Dark Matter model has been assumed, see Section 2.4.1. years old, and significantly smaller than today, have shown that the early Universe was extremely smooth containing only small temperature anisotropies (Komatsu et al., 2009; Hinshaw et al., 2009). Even so, we can observe a vast amount of different type of structure in visible light and in other frequencies (e.g. Blanton & Moustakas, 2009, and references therein) when looking at the night sky. The smooth early Universe must therefore have gone through radical changes when evolving from the smooth primordial gas and dark matter density field to inhomogeneous structures such as dark matter haloes, galaxies, groups and clusters observed, directly or indirectly, today.

In physical cosmology the global evolution of the temperature, pressure, and density fields can be studied using Friedmann’s world models. These models can describe the evolution of an homogeneous and isotropic universe. However, as we do not live in an empty universe, a proof for that is obvious - you are reading this write up - we must also describe how the small anisotropies, seen in the CMB, evolve. During the early times, in the so-called linear regime, the evolution of the density field can be followed using linear perturbation theory and Newtonian gravity. However, after the density perturbations grow enough they start to collapse and their evolution turns nonlinear. Their evolution can still be followed using analytical approximations, however, their validity is limited. Fortunately, more accurate methods have also been developed.

Due to the inherently nonlinear nature of gravity cosmological -body simulations have become an invaluable tool when the growth of structure is being studied and modelled closer to the present epoch. Large simulations with high dynamical range (e.g. Springel et al., 2005; Boylan-Kolchin et al., 2009) have made it possible to model the formation and growth of cosmic structure with unprecedented accuracy. Moreover, galaxies, the basic building blocks of the Universe, can also be modelled to good accuracy in cosmological context and studied from their initial formation down to the present time (e.g. Naab et al., 2007). For example, semi-analytical models of galaxy formation (e.g. Somerville & Primack, 1999; Somerville et al., 2001, 2008) allow us to populate dark matter haloes with galaxies that are formed from baryonic matter when lacking hydrodynamical simulations (e.g. Springel et al., 2005; Croton et al., 2006; Lucia & Blaizot, 2007). Reassuringly, both semi-analytical models as well as hydrodynamical simulations, which both model for example an inflow of gas, how gas can cool and heat up again, how stars are formed within galaxies, and how stellar populations evolve, are in reasonable agreement with observational data at lower redshifts (e.g. Kitzbichler & White, 2007). Instead, at high redshifts the small number of candidate galaxies (e.g. Yan et al., 2010; Labbé et al., 2010; Bouwens et al., 2011) still complicates more detailed comparisons to simulations (e.g. Dayal et al., 2011; Razoumov & Sommer-Larsen, 2010). Even so, simulations can be used to make predictions for different observables and aid when interpreting observational results. The growth of cosmic structures, cosmological -body simulations, and formation of galaxies are briefly reviewed in Chapter 2.

Despite all the simulations and successes in recent years and decades, there are still many unanswered questions in the field of galaxy formation and evolution. One of the longest standing issues in galaxy evolution is the significance of the formation place and thus initial conditions to a galaxy’s evolution in respect to environment, often formulated simply as “nature versus nurture” like in human development and psychology. We are therefore left to ponder if the galaxies we see today are simply the product of the primordial conditions in which they formed, or whether experiences in the past change the path of their evolution. Unfortunately, our understanding of galaxy evolution in different environments is still limited, albeit the morphology-density relation (e.g. Oemler, 1974; Dressler, 1980) has shown that the density of the galaxy’s local environment can affect its properties. For example, on average, luminous, non-starforming elliptical and lenticular galaxies have been found to populate denser regions than star forming spiral galaxies. Consequently, the environment should play a role in galaxy evolution, however, despite the efforts, the exact role of the galaxy’s local environment remains open.

A group of galaxies is the most common galaxy association in the Universe (e.g. Holmberg, 1950; Humason et al., 1956; Turner & Gott, 1976; Huchra & Geller, 1982; Ramella et al., 1995; Zabludoff & Mulchaey, 1998). As such, more than half of all galaxies are found in groups and small clusters. They are therefore important cosmological indicators of the distribution of matter in the Universe. Moreover, groups and clusters can also provide important clues for galaxy formation and evolution physics as the properties of galaxies can be studied as a function of the local environment. As a result, the environmental dependency in galaxy evolution can be understood to some extent in terms of the group environment (e.g. Moore et al. 1996, 1998; Mihos 2004; Fujita 2004, but see also Kauffmann et al. 2004; Blanton 2006). Note, however, that the debate over the exact role of group environment is far from over. Groups of galaxies, their properties and the group environment in the context of galaxy evolution is discussed in Chapter 3.

One fundamental question concerning groups of galaxies, and closely related to galaxy evolution in groups, is whether the identified systems are gravitationally bound or not. This is a valid concern because from the observational point of view, groups and their member galaxies are not well defined. Early studies based on simulations gave hints that not all observed systems of galaxies are dynamically relaxed, but might be in the process of formation (e.g. Diaferio et al., 1994; Frederic, 1995a, b). Despite this, many observational studies, even today, treat all identified groups like they were gravitationally bound structures. The work presented in Papers I and II tackles this issue and is summarised in Chapters 3.5 and 3.6, respectively. Results of this work show that a significant fraction of systems of galaxies are gravitationally unbound when the most often used grouping algorithm, namely Friends-of-Friends, is being applied to simulated data. This result has several important implications for the studies of galaxy groups and for the evolution of galaxies in groups. For example, Paper II shows that groups with a large excess of positive redshifts are more often gravitationally unbound than groups that do not show any significant excess. Fortunately, this prediction can be used, together with the methods developed for and presented in Paper I, to assess whether observed groups are likely to be gravitationally bound or not.

Elliptical galaxies are most often found in dense environments like the cores of groups and clusters (e.g. Dressler, 1980). Yet, observations have shown that there is a significant population of isolated elliptical galaxies that are found in under-dense regions with no bright nearby companions (e.g. Aars et al., 2001; Reda et al., 2004; Smith et al., 2004; Denicoló et al., 2005; Collobert et al., 2006). Whether these galaxies originally formed in under-dense regions or if the local environment has impacted their evolution, in form of a collapsed group, is a profound question with long reaching implications. Key observations of galaxy evolution and the significance of the environment for galaxy evolution are briefly reviewed in Chapter 4.

The results of a theoretical case study of isolated field elliptical galaxies, Paper III, are also presented and summarised in Chapter 4. These results show that three different yet typical formation mechanisms can be identified, and that isolated field elliptical galaxies reside in relatively light dark matter haloes excluding the possibility that all of them are collapsed groups as suggested earlier. Additionally, also another case study concerning luminous infrared galaxies, Paper IV, is discussed. The results of this study imply a strong correlation, such that more infrared-luminous galaxies are more likely to be merger-driven. However, the results also imply that a significant fraction (more than half) of all high redshift infrared-luminous galaxies detected by Herschel Space Observatory are able to attain their high star formation rates without enhancement by a merger. These and other results discussed in Chapter 4 imply that both “nature” and “nurture” play a role in galaxy evolution.

Chapter 2 Formation of Structure

“…the biggest blunder of my life.”

Albert Einstein

2.1 Historical perspective

The study of the structure formation of the Universe dates back to 1610 and Galileo Galilei who realised that the Galaxy can be resolved into stars when observed through a telescope. Galilei’s observations can also be considered as a starting point for early observational cosmology. However, it was René Descartes and Thomas Wright who were likely the first ones to speculate and publish their cosmological views. Already around 1760s Immanuel Kant and Johann Lambert developed the first hierarchical model of the Universe, although since these early models it took almost two hundred years before cosmology developed into a physical science and before the formation of structure in the Universe could be truly appreciated and studied in a physical context. Below I briefly mention a few key moments from the history that have lead to the structure formation theory, as we know it today. However, many important events are not mentioned, thus I refer the interested reader to more comprehensive reviews, see e.g. Ratra & Vogeley (2008).

It was Albert Einstein and his General Theory of Relativity (GR) in 1915 (Einstein, 1915) that gave birth to the physical cosmology as we know it today. GR is a framework that explains one of the four fundamental forces of the Universe, namely gravity, and enabled cosmologists to predict the behaviour of a model universe. As a result, it became possible, for the first time in the history of mankind, to formulate self-consistent models that describe the Universe and large-scale structure. As a consequence, in 1917 Einstein derived the first fully self-consistent model of the Universe (Einstein, 1917). However, soon after he realised that without modifications his field equations predicted that a static Universe was not stable. At the time the Universe was assumed to be static, thus, Einstein introduced a cosmological constant , that enabled a static universe, to solve the issue.

In 1924 Aleksander Friedmann and in 1927 Georges Lemâitre derived solutions for expanding universes, paving the way for evolving universe models. Later, Friedmann’s models became the standard models describing the dynamics of the Universe. In 1927 Lemâitre first proposed what has come to be known as the Big Bang theory of the origin of the Universe, albeit the name was introduced by Fred Hoyle. The framework for the Big Bang model relies on the Einstein’s GR, the Cosmological Principle, Friedmann’s equations, and it is also the standard theory for the origin of the Universe.

In 1935 Howard Percy Robertson and Arthur Geoffrey Walker derived independently the space-time metric for all isotropic, homogeneous, uniformly expanding models of the Universe. However, the Friedmann world models are isotropic and homogeneous, thus, all observable structures such as galaxies and groups are absent. Consequently, the next step towards developing more realistic models of the Universe was to include small density perturbations and to study their development under gravity, namely the formation of structure. Fortunately, already in 1902, well before GR, sir James Hopwood Jeans had shown that the stability of a perturbation depends on the competition between gravity and pressure (Jeans, 1902): gas pressure prevents gravitational collapse on small spatial scales and gives rise to acoustic oscillations. Jeans showed that density perturbations can grow only if they are heavier than a characteristic mass111Now referred to as the Jeans’ mass. scale, while below this scale dissipative fluid effects remove energy from the acoustic waves, which dampens them. The application of the Jeans criterion and the growth of spherically symmetric perturbations in an expanding universe were worked out by Lemâitre and Richard Tolman in the 1930s. Albeit it was not before 1946 when Evgenii Lifshitz worked out relativistic perturbation theory and started applying it to the linear growth of cosmic structure. Finally, a general scheme for structure formation was first outlined by Lev Davidovich Landau and Lifshitz in the 1950s, and developed further by Phillip James Edwin Peebles during the 1970s.

However, to truly appreciate the study of the formation of structure of the Universe, observational constrains on the initial density perturbations were required. The serendipitous discovery of Cosmic Microwave Background (CMB) radiation by Arno Penzias and Robert Wilson in (Penzias & Wilson, 1965), predicted already in by Ralph Alpher, Robert Herman, and George Gamow (Gamow, 1948; Alpher & Herman, 1948), paved the way for understanding the initial conditions of large-scale structure formation. Measurements of the CMB describe the initial conditions after recombination of the very early Universe and led way to the standard model of cosmology. Furthermore, large galaxy surveys of recent years have helped to set constrains for the models of structure formation at low redshift.

2.2 Observational background of cosmology

According to the Big Bang model, the background radiation from the sky measured today comes from the so-called last scattering surface. As the name implies the surface of last scattering is a spherical surface where the Cosmic Microwave Background photons were scattered for the last time222This however is not true for all CMB photons: some have scattered from free electrons that have become available due to reionization. before arriving at our microwave detectors. This decoupling of photons from matter happened years after the Big Bang during the epoch of recombination when the rate of Thomson scattering became slower than the expansion of the Universe (e.g. Weymann, 1966; Peebles, 1968). At that moment, photon interactions with matter became insignificant, leading to the CMB radiation. This moment also defines the “optical” horizon; the largest volume from which we can receive information via photons.

On the largest scales the most robust evidence for the isotropy of the Universe comes from the CMB measurements, while galaxy surveys compliment the CMB information by probing later epochs and smaller angular scales. Therefore, in the next two Sections I will briefly review what the CMB and galaxy survey observations can tell us about the formation of structure.

2.2.1 Cosmic microwave background radiation

The epoch when the ionisation state of the intergalactic gas changed from being a fully ionised plasma to a neutral gas is known as the epoch of recombination. This is the redshift when the detailed anisotropy structure of the early Universe was imprinted onto the Cosmic Microwave Background (CMB). It can therefore provide information on the initial density perturbations and describe the Universe on the epoch well before galaxies, groups and the formation of large-scale structure.

Ever since the discovery of the CMB radiation, observations have played a key role in shaping and constraining the standard model of cosmology (for a review, see e.g. Hu & Dodelson, 2002; Bartelmann, 2010). The CMB observations probe the earliest observable Universe and hence the initial conditions of the Universe, such as its homogeneity, isotropy, and flatness. Based on the CMB observations it has been established, first by Cosmic Background Explorer (COBE) and later by Wilkinson Microwave Anisotropy Probe (WMAP), that the electromagnetic spectrum of the CMB is extremely close to a thermal blackbody with a temperature K. Moreover, COBE and WMAP has established that the largest temperature anisotropies in the CMB are of the order of (Smoot et al., 1992; Hinshaw et al., 2009), as shown in Figure 2.1. The absence of K fluctuations alone show that the matter in the Universe must be dominantly something that does not interact electromagnetically (Peebles, 1982), i.e., “dark”.

Figure 2.1: The foreground-reduced Internal Linear Combination map based on the five year WMAP data. Image from Hinshaw et al. (2009).

The angular temperature fluctuations associated with the primordial density perturbations are assumed to originate in a rather narrow range of redshifts. If this holds then the pattern of the angular temperature fluctuations in the CMB map (Fig. 2.1) gives us a direct snapshot of the distribution of radiation and energy at the moment of recombination. The angular scale corresponds to the Hubble radius at recombination, which can be taken as a dividing line between the small-scale perturbations that have been substantially modified by gravity and the large-scale inhomogeneities that have not changed much. The fluctuations on large angular scales arise from inhomogeneities with wavelengths exceeding the Hubble radius at recombination. As a result, they provide pristine information about the primordial inhomogeneities. On the other hand, sub-horizon perturbations are formed by primordial sound waves.

The anisotropy of the CMB can be divided into two types: primary anisotropy, due to effects which occur at the last scattering surface and earlier; and secondary anisotropy, due to effects such as interactions of the CMB photons with hot gas or gravitational potentials (for example, the Sunyaev-Zel’dovich and integrated Sachs-Wolfe effects), between the last scattering surface and the observer. The structure of the CMB primary anisotropies is mainly determined by two effects: acoustic oscillations and, on small angular scales, photon diffusion (also known as Silk damping (Silk, 1968)). The photon diffusion damping arises from the fact that the photon-baryon-electron (PBE) fluid is not tightly coupled and the photons can diffuse through the fluid, while the acoustic oscillations result from the constructive and deconstructive interference. Overdensities in the dark matter compresses the fluid due to their gravity until the rising pressure in the coupled PBE fluid is able to counteract gravity. The cosmological importance of this is that the PBE fluid underwent acoustic oscillations, while the dark matter, being decoupled, did not.

The density fluctuations in the early Universe are assumed to be critical for structure formation, because they can provide the seeds from which the structures within the Universe can grow and eventually collapse to form the first stars and galaxies (e.g. Bromm et al., 1999; Barkana & Loeb, 2001; Abel et al., 2002; Bromm et al., 2009, and references therein). The density perturbations of the early Universe are thought to have a very specific character when inflation is assumed: they form a Gaussian random field (Bardeen et al., 1986), which is nearly scale-invariant according to the spectral index measured by WMAP (Hinshaw et al., 2009; Komatsu et al., 2009; Dunkley et al., 2009; Jarosik et al., 2011). I will return to this in Section 2.5 where a more detailed discussion of the initial density perturbations is presented.

Power Spectrum of the CMB

To maximise the information a CMB map (such as Fig. 2.1) can provide, the CMB information is most often presented in the form of a power spectrum in terms of the angular scale or multipole moment as shown in Figure 2.2. The power spectrum, which is a spherical harmonic transform333Fourier transformation is not possible on a sphere, thus, spherical harmonics which are analogous are used instead. of the CMB map, and polarisation of the CMB radiation provide a wealth of information (see e.g. Tristram & Ganga, 2007, for data analysis methods) for both constraining cosmological parameters (Table 2.2) and for understanding the formation of the large-scale structure in the Universe. The CMB anisotropy is a powerful cosmological probe because the parameters which determine the spectrum can all be directly related to the basic cosmological parameters such as the energy densities , the dark energy equation of state , and the Hubble parameter .

Figure 2.2: The temperature power spectrum for the seven-year WMAP data. The solid line shows the predicted spectrum for the best-fit flat CDM model (Section 2.4.1). The error bars on the data points represent measurement errors while the shaded region indicates the uncertainty in the model spectrum arising from cosmic variance. Image from Jarosik et al. (2011).

The general shape of the power spectrum (Fig. 2.2) - a plateau at large angular scales (small ) and acoustic peaks at small angular scales (large ) - confirms that the spectrum is predominantly nearly scale-invariant and adiabatic in agreement with the basic predictions of the Big Bang and inflationary paradigm. The dominant acoustic peaks in the CMB power spectra are caused by the collapse of dark matter over-densities and the oscillation of the photon-baryon fluid into and out of these over-densities (Lineweaver, 2003). The underlying physical notion is that the pressure of photons can erase anisotropies, whereas the gravitational attraction of baryons makes them to collapse and to form dense haloes. As a result, these effects can create acoustic oscillations, which give the CMB its characteristic peak structure (see Fig. 2.2). The first acoustic peak is associated with perturbations on the scale of the sound horizon at the last scattering surface Mpc, while following peaks are on the scales less than the sound horizon. Combining the information about the heights and locations of the peaks, many cosmological parameters can be determined with good accuracy independent of other observations such as galaxy surveys.

2.2.2 Large galaxy surveys

The large galaxy surveys of today such as the 2dF Galaxy Redshift Survey (2dFGRS; Colless, 1999; Colless et al., 2001; Percival et al., 2001) and the Sloan Digital Sky Survey (SDSS; York et al., 2000; Stoughton et al., 2002; Abazajian et al., 2003) have provided large, statistically significant samples of different types of galaxies. The ability to explore many dimensions of galaxy properties and scaling relations simultaneously and homogeneously has been greatly beneficial. Moreover, spectroscopic observations and multi-wavelength imaging allows galaxies to be sorted in classes and sub-populations by, e.g., morphology, environment and luminosity, while large sky coverage allows galaxies to be grouped in groups and clusters enabling studies of galaxy evolution as a function of environment (e.g. Blanton & Berlind, 2007; Mateus et al., 2007, see also Chapter 4). However, large galaxy surveys not only provide good statistics for galaxy properties but they can also be used to study the large-scale structure and cosmology.

Despite the fact that the CMB radiation is very smooth, the visible Universe that is dominated by the light from galaxies looks highly inhomogeneous and consists of structures from the scales of isolated galaxies and voids, through groups (Chapter 3) and clusters to superclusters and to large filaments between them (see Fig. 2.3). One aim of large galaxy redshift surveys is therefore to map the three dimensional distribution of galaxies, in order to understand the properties of this distribution and what it implies about the contents and evolution of the Universe. Because the topology of the distribution of galaxies is closely related to the initial conditions of the Universe and to the assumption that the initial perturbations were Gaussian fluctuations with random phases on large scales this mapping can also provide information concerning the conditions in the early Universe independent from the CMB measurements. Unfortunately, the spatial distribution of galaxies, groups and clusters, depends not only on the matter distribution in the Universe, but also on how they form in the matter density field. It is therefore important to understand galaxy formation (Section 2.8) and evolution (Chapter 4) in detail when studying the clustering of galaxies.

A representation of the large-scale distribution of galaxies on the sky in the 2dFGRS is shown in Fig. 2.3. From this figure alone it is obvious that galaxies form larger structures such as clusters and filaments and that the visible light is unevenly distributed in the Universe even on relatively large scales. Even though the distribution of galaxies becomes smoother and smoother when larger and larger scales are considered, non-random structure is still present in forms of superclusters and filaments between them. This non-random structure is sometimes called the cosmic web, as the long filaments of dark and baryonic matter seem to form a “threaded" structure.

The filaments seen in large galaxy surveys are the largest known structures in the Universe and can be up to Mpc444Here refers to the dimensionless Hubble parameter defined such that . Also, see equation 2.3.2 for a definition of the Hubble parameter. long (Einasto et al., 1980; Batuski & Burns, 1985; White et al., 1987; Bahcall, 1988). Filaments are important structures for galaxy formation as it is assumed that they are the channels that carry baryons to the nodes of the filaments where clusters of galaxies are formed. Filaments can also help cool gas to avoid shock-heating (e.g. Dekel & Birnboim, 2006; Dekel et al., 2009a; Kereš et al., 2009), while over-densities in filaments can form gravitationally bound dark matter haloes. There is however still a debate in how gas can enter and coalesce into dark matter haloes and how it cools down to form stars and eventually galaxies (see e.g. Kereš et al., 2005; Kaufmann et al., 2006, and references therein). The exact role of filaments in the formation and evolution of galaxies is therefore currently unclear.

While filaments are the largest known structures, superclusters (e.g. Araya-Melo et al., 2009) are the largest non-percolating galaxy systems (Oort, 1983; Bahcall, 1988; Einasto et al., 2007, 2008). Unlike super clusters, the scales of the largest voids are in general to times the scale of a regular relaxed cluster, i.e., up to Mpc, although the size measurements vary greatly (Zeldovich et al., 1982; Rood, 1988; Vogeley et al., 1994; Lindner et al., 1995; El-Ad & Piran, 1997; Hoyle & Vogeley, 2004; Ceccarelli et al., 2006; von Benda-Beckmann & Müller, 2008; Tinker & Conroy, 2009). Because of all the structure in the cosmic web, the galaxy distribution seems to have a sponge-like topology, with both high- and low-density regions forming an interconnected network, where voids are separated from high-density regions by flattened structures called “walls”. This raises an obvious question: how can we quantify the clustering of different types of objects and what does this tell us about the formation of the large-scale structure?

Figure 2.3: The projected distribution of galaxies in the nearby Universe as a function of redshift and Right Ascension. Earth is at the centre, and each blue point represents a galaxy. Courtesy of the 2dFGRS website.
Two-point correlation function and the power spectrum

Among the simplest methods to measure clustering properties of galaxies (or of other objects like quasars, groups, clusters, etc.) is with the spatial two-point correlation function . It describes the excess probability above Poisson of finding an object at distance from another object selected at random over that expected in a uniform, random distribution (see Peebles, 1980, for a complete discussion). We can now write the probability to find galaxies in infinitesimal small volumes and as follows


where is the mean galaxy density. In practise, it is however often more convenient to derive the two-point correlation function using, for example, the Landy & Szalay (1993) estimator.

Large datasets provided by galaxy surveys have been used to study, for example, the spatial correlation functions (e.g. Connolly et al., 2002; Scranton et al., 2002; Zehavi et al., 2004; Masjedi et al., 2006), clustering of matter, galaxies (e.g. Coil et al., 2007) and groups (e.g. Coil et al., 2006a) in the Universe. Consequently, allowing to set constrains for the structure formation and cosmological parameters (Lahav & Suto, 2004; Tegmark et al., 2004b). In galaxy surveys, redshifts of galaxies are usually used as distances, thus the correlation function is said to work in redshift-space. However, because of peculiar velocities555The peculiar velocity is the velocity that remains after subtracting off the contribution due to the Hubble expansion., an isotropic distribution in real-space will appear anisotropic in redshift-space and vice versa. The redshift-space correlation function therefore differs from the real-space correlation function. This effect is known as the redshift-space distortion (for observational studies, see e.g. Hamilton, 1998; Tegmark et al., 2004a). It is important to note that redshift-space distortions due to peculiar velocities along the line of sight will introduce systematic effects to the estimate of . For example, at small separations, random motions within a virialized overdensity cause an elongation along the line of sight (dubbed as “fingers of God"). On the other hand, on large scales, coherent infall of galaxies into forming structures causes an apparent contraction of structure along the line of sight (dubbed as the “Kaiser effect").

On scales smaller than Mpc the real-space correlation function is well approximated by a power law


where the slope and Mpc is the correlation length. This shows that galaxies are strongly clustered on scales Mpc, and the amplitude of clustering becomes weak on scales much larger than Mpc. Note, however, that the exact values of and are found to depend on the properties of the galaxies. Particularly, brighter and redder galaxies are more strongly clustered than fainter and bluer ones (Norberg et al., 2001; Zehavi et al., 2005; Coil et al., 2006b; Wang et al., 2008). Additionally, early-type galaxies have been found to be much more clustered on small scales, leading to a morphology-density relation (Peacock, 2002, see also Chapter 4).

The galaxy correlation function is a measure of the degree of clustering in either the spatial or the angular distribution of galaxies. The spatial two-point correlation (or autocorrelation) function and the power spectrum forms a Fourier-transform pair, i.e.


Here is the volume within which is defined. Assuming isotropy and that the two-point correlation function is spherically symmetric leads to


Note that the function allows only wave-numbers to contribute to the amplitude of the fluctuations on the scale .

Figure 2.4 shows a matter power spectrum at the present time. According to the figure on large scales (small wave-numbers Mpc) the current matter power spectrum still has its primordial shape. This shape corresponds to a power law dependence on scale; , where is the so-called spectral index, assumed to be close to unity (for theoretical background, see Section 2.5.2 and for the reference value, see Table 2.2). The horizon scale at the epoch when the matter and radiation densities are equal (the matter-radiation equality, , see Table 2.2 for a value) is imprinted upon the power spectrum as the scale at which the spectrum turns over. Hence, the peak position in the spectrum corresponds to the Jeans length (Eq. 2.5.1) at matter-radiation equality (I will return to this in Section 2.5). The position of this turnover corresponds to a physical scale determined by the matter and radiation densities . Moreover, the shape of the observed power spectrum depends on the amount and the nature of the matter in the Universe, providing constrains for cosmology. For example, if all of the dark matter were hot then the matter power spectrum would fall off sharply to zero to the right of the peak.

Figure 2.4: The matter power spectrum as a function of wave-number at the present time. Note the turnover at Mpc. The plot combines data from different scales: CMB, large galaxy surveys, weak lensing and Ly forest, in order of decreasing co-moving wavelength. In addition, there is a single data point for galaxy clusters. Figure from Tegmark et al. (2004a).

Large redshift surveys can be used not only to study the power spectrum of galaxy clustering but also the presence of the acoustic oscillations of baryons (e.g. Tegmark et al., 2004a; Cole et al., 2005; Eisenstein et al., 2005; Percival et al., 2007). The physics of these oscillations are analogous to those of the CMB acoustic oscillations. The amplitude of baryonic acoustic oscillations (BAOs) is however suppressed in comparison to CMB because not all the matter in the Universe is composed of baryons. In practise, many studies of BAOs have taken advantage of the clustering of luminous red galaxies (LRGs; e.g. Padmanabhan et al., 2007; Sánchez et al., 2009). The clustering of LRGs also allows the values of cosmological parameters to be derived independently from the CMB measurements. Furthermore, the observed power spectrum allows to constrain on both the amplitude and the scale dependence of the galaxy bias (e.g. Padmanabhan et al., 2007). This ultimately links the galaxy power spectrum to the matter power spectrum. Despite this connection, large galaxy surveys can provide constrains for structure formation models independent of the CMB measurements.

A viable structure formation model has to therefore simultaneously explain both the smoothness of the CMB and the clear evidence from the galaxy surveys that the assumptions of isotropy and homogeneity do not hold on smaller scales, but galaxies and other inhomogeneities such as groups do form. Moreover, large galaxy surveys together with the CMB observations can be used to set strict constrains for cosmological parameters (see Section 2.4.1 and Table 2.2). Importantly, these observations can provide information about structures of different sizes and times from the recombination to the present epoch. A successful structure formation model must therefore be able to explain all the current observations. However, before we start looking into structure formation models in more detail, some tools to study the dynamical evolution of an expanding universe that is isotropic and homogeneous must be introduced. In the following Sections I shall give a minimalistic overview of the dynamics of an expanding universe on top of which the structure formation theory can be build upon. For more detailed treatment, I refer the interested readers to great textbooks of e.g. Peebles (1980); Dodelson (2003); Longair (2008); Mo et al. (2010).

2.3 Dynamics of an expanding universe

Being able to model the dynamical evolution of an expanding universe is a basic requirement for any structure formation model. Because space-time can be curved and is not static, we must rely on General Relativity when deriving the equations that govern the evolution of the background universe.

2.3.1 Einstein equation

Einstein’s General Relativity (GR) enabled self-consistent models of the Universe to be constructed as it relates matter and energy to the geometrical properties of the Universe. In GR, the gravity field is described by Einstein’s field equation:


where is the Einstein tensor, is Newton’s gravitational constant and is the energy-momentum tensor. Note that matter is incorporated in Einstein’s equation through the energy-momentum tensor. On large scales, matter can be approximated as a perfect fluid characterised by an energy density , pressure and four velocity . Now the energy-momentum tensor may be written as


where the equation of state depends on the properties of matter. Often in cosmologically interesting cases (or more general as in Eq. 2.3.3). Equation 2.3.1 shows that in GR the strength of the gravitational field depends not only on the energy density , but also on the pressure .

It is also possible to write the Einstein equation in a form that explicitly shows the cosmological constant , now


where is the Ricci tensor, is the scalar curvature, is the cosmological constant, and is the unit tensor666Note that later in this Chapter will refer to the density perturbation field. defined by the metric such that . As the cosmological constant can be interpreted as the contribution of vacuum energy to the Einstein equation it can also be included in the energy-momentum tensor.

According to the Cosmological Principle (Milne, 1933) the Universe is homogeneous and isotropic, at least on large enough scales, thus space-time can be described by the Robertson-Walker metric:


where the spatial positions are described by spherical coordinates , is the speed of light and is the scale factor. The scale factor and the changing scale of the Universe causes the cosmological redshift , which can be defined as


Here is the present value of the scale factor and is the value of the scale factor at the time when the light was emitted. As a result, the redshift can be used to parameterise the history of the Universe; a given corresponds to a time when our Universe was times smaller than now. The importance of the metric 2.3.1 is that it allows to define the invariant interval between events at any epoch or location in an expanding universe, and thus determines the metric, Riemann, and Ricci curvature tensors.

The Einstein equation describes the geometry of space which is curved by matter and energy. It is also the basic equation of GR that the dynamical variables characterising the gravitational field must follow. The cosmological evolution of relativistic matter can therefore be derived from the Einstein equation (2.3.1) when the metric of the space-time and the energy-momentum tensor are fixed.

2.3.2 Friedmann equations

Alexander Friedmann was the first to derive a pair of equations that can describe the expansion rate of a homogeneous and isotropic universe. The first Friedmann equation can be derived directly from Einstein’s field equation (2.3.1), which can be simplified using the Robertson-Walker metric (Eq. 2.3.1). Now, Friedmann’s first equation can be written as


where is Newton’s gravitational constant, is the spatial curvature ( or ), is the sum over all energy densities (e.g. baryons, photons, neutrinos, dark matter and dark energy), and is the Hubble parameter that describes the rate of the expansion:


Here denotes the time derivative of the scale factor. Note that the dynamical content of the metric is encoded in the function . Friedmann’s second equation can be derived from the trace of Einstein’s field equation and can be written as:


Now the is the second time derivative of the scale factor and is the pressure.

The importance of the two Friedmann equations (2.3.2 and 2.3.2) is that they determine the two unknown functions; the scale factor and the energy density . Consequently, the Friedmann equations can describe the dynamics of a homogeneous, isotropic and expanding universe, because the scale factor completely describes the time evolution of such a universe. The Friedmann equations therefore determine the expansion rate of the Universe, based on the density of material within it and the curvature of space. Note that, a Robertson-Walker metric (2.3.1) whose scale factor satisfies Friedmann’s equations is called the Friedmann-Lemait̂re-Robertson-Walker metric and the cosmological standard model upholds that the Universe at large is described by such a metric.

Construction of a cosmological model requires solving the Friedmann equations, resulting to the expansion rate as a function of time, and hence the size of the Universe. However, solving the Friedmann equations alone is not enough for a cosmological model, it is also essential to know how the energy density changes as a function of time.

2.3.3 Evolution of energy density: the fluid equation

The evolution of energy density in an expanding universe can be described with the so-called fluid equation. The fluid equation, which holds only for adiabatic processes, can be obtained from the Friedmann’s equations when equation 2.3.2 is rewritten using equation 2.3.2 resulting in:


where denotes the energy density and is the pressure, as defined earlier. The fluid equation expresses the conservation of mass-energy, thus it is also known as the energy conservation or the continuity equation. The second term in the fluid equation corresponds to the loss in energy because the pressure of the material has done work as the volume of the Universe increased. However, because the energy is conserved the energy lost from the fluid via the work done goes into gravitational potential energy. To generalise, for an adiabatically expanding volume the entropy per unit co-moving volume is conserved, and the expansion of the Universe causes an increase or decrease of its internal energy depending on whether the pressure is smaller or larger than zero.

For a given equation of state, , the fluid equation gives the density and pressure as a function of the scale factor . The equation of state, which describes the energy content of the Universe, is often parametrized with in the following way:


Here the subscript denotes the species of the material. Note, however, that this so-called perfect fluid hypothesis holds only for material whose pressure is directly related to its density. Finally, if is time-independent, then substituting Eq. 2.3.3 into 2.3.3 gives the time evolution of the mean energy density of the Universe as follows:


The fluid equation together with the equation of state allow the derivation of the mean energy density, temperature and pressure of the Universe at any redshift from their values at the present time. At early times the Universe is assumed to be radiation dominated thus we can approximate that the energy content of the Universe is dominated by an ultra-relativistic radiation fluid for which . As a result, the mean energy density evolves proportional to . After the nucleosynthesis, but before the recombination, at 777The exact time depends on the matter density of the Universe and the Hubble constant. The redshift given is for the reference values, see Table 2.2., the densities of non-relativistic matter and relativistic radiation are equal. However, after this point the Universe turns into a matter dominated. Now, a non-relativistic gas can be approximated with a fluid of zero pressure888sometimes referred to as dust and we take . As a consequence, the mean energy density evolves . Finally, at recent epochs the energy density seems to have become dominated by vacuum energy. In order to keep a constant energy density as the Universe expands, the pressure must be negative, therefore, for vacuum energy we take . Note that can also be taken for the cosmic inflation.

Table 2.1 summarises the evolution of energy density, pressure, and temperature as a function of the scale parameter . Note, however, that these scaling relations only hold if the equation of state remains constant with respect to time, while in reality such a simplification may not hold on all times. It should also be kept in mind that although the contribution of baryons and photons to the present day energy budget is small, they make an important contribution to shaping the matter power spectrum. Moreover, in realistic cases the Universe is not made out of a single material component as presented above. Fortunately, each material, baryons, photons, dark matter, neutrinos, dark energy, etc., obey their own fluid equation containing the appropriate expression for its pressure if the fluids are non-interacting. One can therefore take a linear combination of the terms and substitute that into the Friedmann equations in case of more realistic models.

Dominant component Energy density Pressure Temperature
Vacuum energy
Table 2.1: Thermodynamics of a homogenous and isotropic universe.

2.4 Ingredients of structure formation

The ultimate goal of a structure formation theory is to describe how the phase transition progressed from almost perfectly homogeneous initial fields to all the structure we observe. To get closer to achieving this goal - to model the formation of structure in an evolving background universe - we must next concentrate on small density perturbations. A realistic structure formation model must be able to describe the evolution of the density field in the Universe with time when the field contains small fluctuations. The usual approach is to model the fluctuations as a perturbation to a smooth background which we assume is homogeneous and isotropic.

The key idea of any structure formation model is that if there are small perturbations, i.e., fluctuations in the energy density of the early Universe, then gravitational instability can amplify them leading to virialized structures such as the galaxies, groups, and clusters we observe today. To model the formation of structure in a realistic, self-consistent, and physical way several ingredients are required (Coles & Lucchin, 2002):

  1. a background cosmology (Section 2.3),

  2. an initial fluctuation spectrum (Section 2.5.2),

  3. a choice of fluctuation mode and a statistical distribution of fluctuations,

  4. a Transfer function (Section 2.5.3),

  5. a recipe for the nonlinear evolution (Sections 2.6 and 2.7), and

  6. a prescription to relate mass fluctuations to observable light (Section 2.8).

As can be seen from the comprehensive list above, detailed modelling of structure formation involves several different ingredients. All of which interact in a complicated manner. To complicate the matter even further, most of the above items involve one or more assumptions (see e.g. Coles & Lucchin, 2002). It should, however, be kept in mind that most of the assumptions are physically motivated, albeit this does not exclude the possibility that they are inaccurate or even incorrect.

Figure 2.5: A logic flow chart for the formation of structure. The exponential expansion refers to the expansion of the Universe during inflation. The radiation, matter and dark energy domination refers to the dominant energy content during a given time period, while the linear and nonlinear regimes refer to the evolution of the density perturbation field. The gastrophysics indicate the regime where hydrodynamics and also forces other than gravity play a significant role. Red hexagons show the initial and the observed power spectrum. Galaxy formation is described in more details in Figure 2.8.

Figure 2.5 shows a logical flow chart describing the formation of structure. The chart shows how structures start to form from small initial density field fluctuations after the Big Bang and proceed through various steps to galaxy formation and to the matter power spectrum observed, for example, in large galaxy surveys (Section 2.2.2). At present, the physical mechanism that can best describe the initial density field is inflation (Guth, 1981; Mukhanov & Chibisov, 1981; Linde, 1982; Narlikar & Padmanabhan, 1991). The initial conditions of the Universe are thought to arise from the scale-invariant quantum-mechanical zero-point fluctuations of the scalar field that drove the inflation in the very early Universe (Guth & Pi, 1982; Hawking, 1982; Baumann, 2007). Inflation predicts, for example, that the initial fluctuations are adiabatic (i.e. perturbations are in thermal equilibrium) and behave as a Gaussian random field with a nearly scale invariant spectrum. In the following Sections I will therefore concentrate on adiabatic fluctuations with a Gaussian random phase and leave out isocurvature perturbations and non-Gaussian random fields (this is the structure formation ingredient III).

The first ingredient of the structure formation describes the global evolution of the background universe (Section 2.3), and is usually done using Friedmann equations (Section 2.3.2). The rest of the ingredients are related to the formation and evolution of density perturbations under gravity in an expanding background universe. These will be described in the following sections in more detail: the spectrum of the initial fluctuations and the Transfer function are discussed in Section 2.5, while the nonlinear evolution using both analytical methods and cosmological body simulations is explored in Sections 2.6 and 2.7, respectively. Finally, to form realistic galaxies that can be compared to the galaxies observed in large galaxy surveys at least hydrodynamics, baryonic physics, and star formation (e.g. McKee & Ostriker, 2007) must be considered. One realisation for modelling gas, star formation, and feedback processes is the semi-analytical models of galaxy formation, which are discussed in Section 2.8.4. However, before exploring the theory of structure formation, lets introduce the most successful structure formation model so far in more detail.

2.4.1 The Cold Dark Matter model

The Cold Dark Matter (CDM) -model is nowadays accepted by the majority of astronomers as a standard model of Big Bang cosmology and cosmological structure formation. The success of CDM is widely recognised and is due to its simplicity, yet at the same time, it has the capability to simultaneously explain several profound observations of the Universe. CDM can explain the structure and existence of the cosmic microwave background, the large-scale structure of galaxy groups and clusters, weak and strong gravitational lensing, and the accelerated expansion of the Universe inferred from type Ia supernovae observations (e.g. Narlikar & Padmanabhan, 2001). The statistical analysis of observations (Section 2.2) strongly support the flat CDM cosmological model with the total energy density equal to the critical density. Table 2.2 summarises the energy budget and the values of the basic CDM parameters based on WMAP results (Hinshaw et al., 2009; Komatsu et al., 2009; Dunkley et al., 2009; Jarosik et al., 2011).

The term of the standard cosmology stands for the Einstein’s cosmological constant, often assumed to be the vacuum energy of space, and dubbed as dark energy due to its unknown origin (Narlikar & Padmanabhan, 2001; Frieman et al., 2008). The energy budget of the model is dominated by this unknown dark energy; the -term constitutes almost per cent of the total energy composition. The remaining per cent of the energy budget is matter, however, about per cent of the matter is assumed to be in form of non-baryonic dark matter. Because the dark matter particles are assumed to be non-relativistic, the dark matter is said to be “cold”. The term cold dark matter was introduced by Peebles and Richard Bond in to cover the wide range of particles that were (and have been) suggested for the origin of this unknown gravitating material. Note that the coldness of the dark matter particles is actually required by the large-scale structure: hot dark matter, i.e., relativistic particles, does not predict enough structure on small scales.

The shape of the CDM power spectrum is such that structures form from smaller to larger structures, i.e., “bottom-up”, with galaxies forming first followed by the formation of groups and clusters. More general, in all CDM-models, independent of the -term, the initial density fluctuations have larger amplitudes on smaller scales, thus the CDM-models are hierarchical; larger structures form by clustering of smaller objects via gravitational instability (e.g. Davis et al., 1985; Frenk et al., 1988; White & Frenk, 1991; Bullock et al., 2001). In CDM-models the collapse of matter happens when a local perturbation starts to turn around while the Universe is expanding. The process of collapsing continues until the internal velocity of system’s components are large enough to hold the system against more collapse. As a result, a dark matter halo is formed. Note that unlike baryonic matter, the behaviour of dark matter does not depend on the scale of the system since dark matter only interacts gravitationally. Thus, in the CDM-model dark matter haloes with different sizes and masses are scaled versions of each other.

Quantity Symbol Value Error
Total Density
Dark Energy Density
Matter Energy Density
Dark Matter Density
Baryonic Matter Density
Hubble Parameter
Power Spectrum Normalisation
Scalar Spectral Index
Redshift of Matter-radiation Equality
Redshift of Decoupling
Age of Decoupling yr
Sound Horizon at Decoupling Mpc
Redshift of Reionization
Reionization Optical Depth
Age of the Universe Gyr
Table 2.2: Values of the basic CDM parameters based on the WMAP 7-year results (Jarosik et al., 2011).

Even though the CDM-model is widely accepted, there are still many unknowns. For example, several particle candidates exist for dark matter (see e.g. Baltz, 2004; Muñoz, 2004; Bertone et al., 2005), however, none has been observed thus far. One of the leading candidates for the dark matter particle is the lightest stable supersymmetric particle called neutralino, which is weakly interacting and massive, but several other candidates, for example Axions, exist. Hence, the type of the dark matter particles is yet to be confirmed. The nature of dark energy is even more mysterious, though the leading candidate is the vacuum energy of space (Peebles & Ratra, 1988; Frieman et al., 2008). As both dark matter and energy reveal themselves only via gravity all attempts to detect them directly have been unsuccessful thus far. A significant amount of work remains therefore to be done before the formation of structure and cosmology can be considered as fully understood.

2.5 Evolution of initial perturbations: structure formation

The present consensus in cosmology, as seen in the previous Sections, is that the observed structures developed from small initial perturbations of the physical fields (density, velocity, gravitational potential, etc.) resulting from the instability of the Friedmann models for small perturbations. Hence, to model structure formation we must model how the density perturbation field evolves. A full treatment would require GR and would proceed by perturbing the background metric and the energy-momentum tensor (see e.g. Weinberg, 2008). However, a fully relativistic treatment is beyond the scope of this introduction to the formation of structure, thus I will describe a Newtonian method, which gives an excellent approximation. Even so, it should be kept in mind that the following analysis holds only for perturbations on scales much smaller than the Hubble radius (i.e. on sub-horizon scale), because the Newtonian description assumes instantaneous gravity (i.e. the speed of gravity has been assumed to be infinite).

2.5.1 Growth of small perturbations in an expanding universe

The hierarchy of cosmic structures is assumed to have grown from primordial density seed fluctuations, which can be described with the density contrast . The density contrast, as a function of the co-moving coordinates , motivates the study of the density perturbation field. It can be defined as


where is the mean density. A critical feature of the field is that it inhabits a universe that is isotropic and homogeneous in its large-scale properties.

In order to describe the structure formation in an expanding universe we must follow the evolution of the initial perturbation field as a function of time, while gravitation magnifies the perturbations in both the baryonic and dark matter distribution. After recombination the amplitude of the density fluctuations is . Thus, the onset of structure formation happens well within the linear regime where . Consequently, a linear perturbation theory can be used as long as the density field does not turn nonlinear (for a comprehensive presentation, see e.g. Coles & Lucchin, 2002).

On large scales matter can be described with a perfect fluid approximation. Thus, at any given time matter can be characterised by the energy density distribution , the entropy per unit mass , and the vector field of three-velocities . These quantities satisfy the hydrodynamical equations that allow the study of the behaviour of small perturbations in a homogeneous, isotropic background. The equations of motion for a non-relativistic fluid are the continuity, Euler, and Poisson equation. The continuity equation states that the change in the mass inside an element of the fluid equals to the mass convected into the element, i.e., it defines the conservation of mass. On the other hand, the Euler equation states that the acceleration of a small fluid element is due to the difference in pressure acting on opposite sides of the element, while the Poisson’s equation describes the relation between the potential fluctuations and the density perturbations causing them.

In their basic form (see e.g. Coles & Lucchin, 2002; Longair, 2008) the continuity, Euler and Poisson equation hold for a smooth background. However, in this Section we derive the evolution of the density perturbations and hence we must consider small perturbations to the background. All the quantities involved must be then written as a sum of the smooth background and the perturbed quantity, e.g., in case of pressure: . Here, corresponds to the smooth background pressure field, while is a small perturbation to this smooth component. The perturbed quantities can then be substituted to the basic continuity, Euler and Poisson equation. After ignoring terms higher than the first order in the perturbation and subtracting the zero order equations one obtains the linear perturbation forms of these equations. Thus, it is possible to describe the evolution of the density, velocity, potential and pressure fields in an expanding universe with Newtonian gravity by using the perturbed continuity, Euler, and Poisson equations, together with the conservation of entropy (for a full derivation, see e.g. Peebles, 1980; Longair, 2008; Mo et al., 2010).

The linear perturbation theory, briefly described above, holds that, during the matter-dominated era, the density field of sub-horizon perturbations can be described with the growth equation as follows


Here is the Hubble parameter (defined in Eq. 2.3.2), is the speed of sound, denotes the differentials with respect to co-moving coordinates, and is a partial time derivative. The above second-order growth equation has been written in a general999Setting gives the often shown form of the growth equation. but single-fluid form as a function of cosmic time and co-moving coordinate . This gives a linear approximation for the growth of density perturbations in an expanding universe. Note that the second term on the left-hand side is the so-called Hubble drag term, which tends to suppress perturbation growth due to the expansion of the Universe. On the other hand, the first term on the right-hand side is the gravitational term, which causes perturbations to grow via gravitational instability, while the last term on the right-hand side is a pressure term and is due to the spatial variations in density.

The solution to the growth equation (2.5.1) apply to the evolution of a single Fourier mode of the density field. However, in the linear regime, the equations governing the evolution of the perturbations are all linear in perturbation quantities. It is then useful to expand the perturbation fields in chosen mode functions. If the curvature of the Universe can be neglected, the mode function can be chosen to be plane waves. Now the perturbation fields can be represented by their Fourier transforms. We therefore seek a wave solution for of the form


We can now write a wave equation for after taking a Fourier transform of equation 2.5.1. Because each mode is assumed to evolve independently, we can write


The equation 2.5.1, sometimes called the Jeans equation, describes the evolution of each of the individual modes , corresponding to . Note that because (in Eq. 2.5.1) is in co-moving units, the wave-vectors are also and that the derivatives (denoted by dots) are time derivatives, because does not explicitly depend on a spatial position.

If the dark matter is cold and collisionless we can neglect the pressure term in Eq. 2.5.1. This allows us to write a general solution in a form of two linearly independent power laws, i.e.,


where and are constants to be determined by initial conditions. The growth (or Jeans) equation therefore has two solutions: a growing and decaying mode. The latter is hardly interesting for structure formation, thus, hereafter we concentrate on the growing mode. The growing mode is described by the growth factor defined such that the density contrast at the time is related to the density contrast today by


It is useful to note that the solution to equation 2.5.1 can either grow or decrease depending on the sign of the


The density perturbations can grow only if the second term in the above equation dominates, while the transition takes place at the wave-number for which the two terms are equal, at


We can now write 2.5.1 in terms of the physical wavelength using the simple relation , and by doing so obtain the Jeans length:


which defines a scale length on which structures can grow. In general, on scales smaller than the Jeans length, i.e., or , the solution to the Jeans equation corresponds to a sinusoidal sound wave; so pressure can counter gravity. Due to the damping caused by the Hubble drag term, there is no growth of structure for sub-Jeans scales, but the solution is oscillatory. Instead, on scales longer than the Jeans length, but smaller than the horizon, the pressure can no longer support the gravity and the solution can grow. As the Jeans length is time dependent in an expanding universe, for example, before the recombination Mpc, while after the Jeans length is only kpc, a given mode may switch between periods of growth and stasis governed by the evolution of . But what does all this mean for the growth of perturbations?

The growth rate of a density perturbation depends on epoch or more precisely on what component dominates the global expansion, whether a perturbation -mode is super- or sub-horizon, and the Jeans length. As already noted in the case of the evolution of the energy density (Section 2.3.3), the early Universe was assumed to be radiation dominated until the time of matter-radiation equality . During the epoch of radiation domination the Universe can be taken as flat and for the growth equation can be solved by (for a detailed derivation, see e.g. Coles & Lucchin, 2002). Thus, for all perturbations, the growing mode (here we ignore the decaying mode) outside the horizon (on super-horizon scales) grows as . Instead, on the sub-horizon101010Note that when the universe was radiation dominated and thus the Jeans length is always close to the size of the horizon. scales, the cold and collisionless dark matter, which has no pressure of its own and is not coupled to photons, grows at most logarithmically . After the matter-radiation equality, matter begins to dominate the dynamics. On the super-horizon scale all perturbations (dark matter, baryons, and photons) grow as . Dark matter, being pressureless, grows with the same rate (e.g. Coles & Lucchin, 2002) also on sub-horizon scales. However, baryons are still coupled to the radiation until the time of decoupling . As a result, the sub-horizon perturbations in the baryons cannot grow but instead oscillate. Finally, at the time of decoupling , the rate of collisional ionization does not dominate any longer and the baryons can decouple from the photons. At this point the baryonic perturbations can start to grow as on scales . On the smaller scales they instead continue to oscillate. Finally, at the latest times when the -term is assumed to be dominant, the growing mode solution is .

Note, however, that in general the presentation in this Section applies only to adiabatic perturbations in a non-relativistic fluid with a single component. Fortunately though, the Newtonian perturbation theory is valid even with the presence of relativistic energy components, such as radiation and dark energy, as long as they can be considered smooth and their perturbations can be ignored. In this case they contribute only to the background solution.

2.5.2 Statistical description of the initial fluctuations

In the previous Section an equation (2.5.1) describing the evolution of the density perturbation field (Eq. 2.5.1) was presented. However, as was noted, it is often convenient to consider the density perturbation field by the superposition of many modes. The natural tool for achieving this is via Fourier analysis in case the co-moving geometry is flat or can be approximated as such (Eq. 2.5.1). In such case the power spectrum can be defined as


where the angle brackets indicate an average. As a result, the power spectrum describes how much the density field varies on different scales. Note that in an isotropic universe, the density perturbation spectrum cannot contain any preferred direction. Thus, we must have an isotropic power spectrum and we can write it simply as a function of wave-number rather than a vector. Even with such simplification, the power spectrum provides a complete statistical characterisation of a particular kind of stochastic process: a Gaussian random field (Bardeen et al., 1986).

Thus, the power spectrum can characterise the statistical properties of the cosmological perturbations. This is highly useful in order to be able to relate theory to observations. For example, results of large galaxy surveys (e.g., correlation functions, Section 2.2.2) suggest that the spectrum of the initial fluctuations must have been very broad with no preferred scales. Because power spectrum is related to a two-point correlation function by a Fourier transform (Eqs. 2.2.2 and 2.2.2), it is then natural to assume that the power spectrum of the initial fluctuations generated in the early phases of the Big Bang is of a power-law form. Thus,


where is the amplitude, is a wave-number in physical units Mpc and is a free parameter. Hence, the power spectrum also describes the normalisation of the spectrum of density perturbations on large physical scales.

Because there is yet no theory for the origin of the cosmological perturbations, the amplitude of the power spectrum has to be fixed by observations. The amplitude is therefore set using either the so-called COBE normalisation111111In the COBE normalisation the amplitude of the large scale temperature anisotropies in the CMB are used to constrain the amplitude. (e.g. Bunn & White, 1997) or by the variance of the density fluctuations within spheres of Mpc radius, (for the reference model value, see Table 2.2). The parameter can be measured, for example, using the cosmic-shear autocorrelation function, the abundance and evolution of the galaxy-cluster population, the statistics of Lyman- forest lines (Seljak et al., 2006), or by counting the number of hot X-ray emitting clusters in the local universe (see e.g. Bartelmann, 2010, and references therein).

According to equation 2.2.2 the power-law form of the power spectrum

corresponds to a two-point correlation function of form


The mass within a fluctuation is , thus, the spectrum of a mass fluctuation is


Finally, the root-mean-square (rms) density fluctuation at mass scale can be written as


As can be seen from the equations above, the spectral index has a significant role in the structure formation. If , the power spectrum is called tilted: a tilted spectrum is called “red” if and “blue” if (for a review, see e.g. Abbott & Wise, 1984; Lucchin & Matarrese, 1985). A red spectrum shows that there is more structure at large scales, while a blue spectrum describes that there is more structure at small scales. The special case is found when the spectral index equals unity.

The Harrison-Zel’dovich power spectrum

Simple inflationary theories predict that right after inflation the matter power spectrum would have a simple power-law form. Consequently, the primordial, or Harrison-Zel’dovich (Harrison, 1970; Zeldovich, 1972), power spectrum can be written as with the spectral index equals unity. The simple power-law form, , now results the spectrum of density perturbations to have a following form


while the two-point correlation function takes the form:


The importance of the Harrison-Zel’dovich spectrum is the property that it is scale-invariant: the density contrast had the same amplitude on all scales when the perturbations came through their particle horizons during the radiation dominated era. Interestingly, the scale-invariant spectrum corresponds to a metric that is a fractal, leading to a fractal nature of the Universe (e.g. Jones et al., 1988; Balian & Schaeffer, 1989).

The current large-scale observations (see Section 2.2.2 and the citations therein) are reasonably well fit by an scale-invariant primordial spectrum of perturbations. Theoretically the spectral index would be precisely unity if inflation lasts forever. As this is obviously not the case, the spectral index must however deviate slightly from the unity. It can be shown (for a detailed discussion, see Liddle & Lyth, 2000) that must be slightly smaller than unity, in agreement with the CDM -model value ( in Table 2.2). It can also be shown that the spectral index of the temperature fluctuations (Fig. 2.1) as a function of angular scale depends only upon the spectral index of the initial power spectrum. Thus, for the Harrison-Zel’dovich power spectrum, the amplitude is independent of the angular scale.

2.5.3 The Transfer function

If we wish to model how the form of the power spectrum evolves as a function of time, the statistical description of the initial fluctuations described by the initial power spectrum must be evolved. During the evolution the form can be modified by several physical phenomena. For example, radiation and relativistic particles can cause kinematic suppression of growth of the initial perturbations. Moreover, the imperfect coupling of photons and baryons may also cause dissipation of perturbations. On the other hand, gravity will amplify the perturbations and eventually leads to collapsed and bound structures. Thus, real power spectra result from modification of any primordial power by a variety of processes: growth under self-gravity, the effects of pressure, and dissipative processes. In general, however, modes of short wavelengths have their amplitude reduced relative to those of long wavelengths.

A possible way to quantify how the shape of the initial power spectrum is modified by different physical processes as a function of time is to use a simple function of a wave-number, namely the Transfer function . For statistically homogeneous initial Gaussian fluctuations, the shape of the original power spectrum is changed by physical processes and the processed power spectrum is related to its primordial form via the Transfer function as follows


Here is the solution of the linearised density perturbations equation (2.5.1), i.e., the growth factor. Hence, once the Transfer function is known, one can calculate the post-recombination power spectrum from the initial conditions.

The form of the Transfer function is a function of the amount and type of the dark matter particles. As Section 2.4.1 described, the currently favoured dark matter particles are non-relativistic. Damping processes can also effect during the linear evolution. As noted already, the cold dark matter does not suffer from strong dissipation, but on scales less than the horizon size at matter-radiation equality there is a kinematic suppression of growth on small scales. Additional complication to the form of the Transfer function arises from having a mixture of matter (both collisionless dark matter and baryonic plasma) and relativistic particles (collisional photons and collisionless neutrinos). One more complication for the shape of the Transfer function arises from the fact that sub-horizon perturbations grow differently during the radiation and matter dominated eras (Section 2.5.1). Due to the complicated form of the Transfer function it therefore must, in general, be calculated using an approximation formula, e.g., by Bond & Efstathiou (1984); Bardeen et al. (1986) or more precisely numerically using publicly available programs such as CMBfast121212 (Seljak & Zaldarriaga, 1996).

2.5.4 Evolution of the initial power spectrum to the present time

Section 2.2 showed that the current matter power spectrum is far from its initial form, even though at scales Mpc a Harrison-Zel’dovich power law is a good approximation. Thus, if one assumes that the initial power spectrum has a Harrison-Zel’dovich form after the inflation, it must have evolved significantly. As described above, the Transfer function can describe how the shape of the initial power spectrum evolves through the epochs of horizon crossing and radiation-matter equality. For the largest scales Mpc the perturbations are still small even today, and one can use the Transfer function. However, for smaller scales such as galaxies, groups and clusters, the inhomogeneities have become so large at later times that the physics of structure growth has become nonlinear.

The caveats in mind, we can still provide an approximation for the form of the present-day power spectrum. If the mass fluctuations inside the co-moving horizon radius at matter-radiation equality are independent of time, and assuming that the dark matter is cold, the present day power spectrum can be approximated with a functional form

where the co-moving wave-number Mpc. By and large this form is in agreement with Figure 2.4, however, the real power spectrum has a smooth maximum, which is caused by the different rates of growth before and after the matter-radiation equality. Note that the co-moving wave-number describes the location of the peak and is set by the co-moving horizon radius at matter-radiation equality. The steep decline for structures smaller than the horizon radius reflects the suppression of structure growth during radiation domination.

The evolution of the matter power spectrum from its initial to the current form (see Fig 2.4) can be summarised as follows. Before the matter-radiation equality, the matter distribution followed mostly that of the radiation. Because radiation has significant pressure perturbations were forced to oscillate on sub-horizon scales. Instead, the largest perturbations were too large for radiation pressure to be able to hold back the collapse, thus, perturbations on scales larger than the horizon scale were able to grow (with the rate given in Section 2.5). This caused the matter power spectrum to increase in height and started to induce a bump at small scales (larger ). As the times went on, the horizon scale increased and larger scales were able to oscillate, causing the bump to shift to larger scales (smaller ). Instead, after matter-radiation equality, , the dark matter was able to grow with a rate given in Section 2.5. Because the dark matter is assumed to be cold and collisionless it does not have significant pressure, and hence there were practically no more acoustic oscillations. During the matter dominated epoch the lack of pressure in dark matter allows it to continue to collapse, causing the whole of the matter power spectrum to increase. As a result, the turn-over point of the matter power spectrum is frozen into the power spectrum, which corresponds to the co-moving horizon radius at matter-radiation equality. Finally, closer to the current time, the small scale perturbations (high ) have turned nonlinear. Consequently, their growth is fast causing the matter power spectrum to rise on smaller scales. Quantitatively, such an evolution has been noted to lead to the current form of the matter power spectrum.

Unfortunately, the power spectrum can only describe the statistical properties of the density contrast, not the evolution of the individual density fluctuations. As a result, the evolution of individual perturbations as well as the nonlinear evolution on small scales must be studied using some other techniques. In such a case one must resort to, for example, analytical techniques or cosmological -body simulations.

2.6 Nonlinear evolution: I. Analytical methods

As the density contrast (Eq. 2.5.1) approaches unity, the evolution of the density fluctuations becomes nonlinear. In the course of nonlinear evolution, overdensities contract, causing matter to flow from larger to smaller scales causing the power spectrum to deform. Even though the linear perturbation theory (Section 2.5.1) fails for , the onset of nonlinear evolution can still be described analytically.

2.6.1 Spherical top-hat model

The simplest analytical model for the nonlinear evolution of a discrete perturbation is called the spherical top-hat model (for a textbook review, see e.g. Padmanabhan, 1993). In this approximation the perturbation evolves according to Birkhoff’s theorem; in a spherically symmetric situation, matter external to the sphere will not influence its evolution. An evolving density perturbation will therefore eventually stop expanding, turns around, and collapses. Because in this model the perturbation has no internal pressure, it collapses to infinite density. Interpreting this literally leads to a conclusion that all spherical perturbations would result to black holes. In reality, this however has obviously not happened, and thus, one should be aware of its limitations.

In the spherical top-hat model a perturbation reaches a maximum size at the time of a turnaround . According to the model, the perturbation will then collapse at . In a realistic case during the collapse the gravitational potential energy must be converted into kinetic energy of the particles involved in the collapse because a collisionless system cannot dissipate energy. This can be achieved, for example, via the process of violent relaxation (Lynden-Bell, 1967). The collapsed object will therefore eventually relax to a structure supported by random motions and satisfy the virial theorem, in which the internal kinetic energy of the system is equal to half of its gravitational potential energy (see Section 3.5.2 and equation 3.0).

If one assumes that the relaxed object virialises at , then a mean overdensity within a virial radius can be derived using the virial theorem. The mean overdensity within the virial radius at the time of virialisation can now be written as


where is the background density at time . In the case of a non-zero cosmological constant and for a flat universe, an approximation for the mean overdensity can be written


where (Bryan & Norman, 1998). This simple approximation can be used to derive an average density of a virialized object formed through gravitational collapse in an expanding universe. In the simplest approximation , this leads to . Note also that the virial theorem and the spherical collapse model can also be used to estimate the redshift at which the object became virialized.

One of the shortcomings of the spherical top-hat model is however the assumption that the perturbations were exactly spherically symmetric. Hence, a more general approximation should be considered.

2.6.2 Zel’dovich approximation

Given the fact that fluctuations of early times were small, it is reasonable to assume that at later epochs only the growing mode has a significant amplitude. Now, if one assumes that the density field grows self-similarly with time, the onset of nonlinear evolution can be described by the so-called Zel’dovich approximation (Zel’Dovich, 1970). The Zel’dovich approximation is a form of the linear perturbation theory and is applicable to a pressureless fluid. The basic assumptions of the Zel’dovich approximation are as follows: the scales of interest are much smaller than the size of horizon; the universe is dominated by the matter component; and the curvature of the universe is zero. The Zel’dovich approximation does therefore not assume spherical symmetry, like the top-hat model.

The Zel’dovich approximation is a Lagrangian description for the growth of perturbations and it specifies the growth of structure by giving the displacement and the peculiar velocity of each mass element in terms of the initial position (Zel’Dovich, 1970; Shandarin & Zeldovich, 1989). Furthermore, the Zel’dovich approximation is a kinematic approximation by nature: particle trajectories are straight lines. The first nonlinear structures to form in this approximation will be two-dimensional sheets, also called the Zel’dovich pancakes. However, the approximation is not valid after the formation of the pancakes when shell crossing will start to occur. Thus, other techniques to follow the evolution of density perturbations further into the nonlinear regime are required.

2.6.3 The Press-Schechter formalism

The Press-Schechter (PS; Press & Schechter, 1974) theory, which was derived heuristically using the linear growth theory and the spherical top-hat model (Section 2.6.1) provides an analytical description for the evolution of gravitational structure in a hierarchical universe. In the PS formalism, an early universe is assumed to be well-described by an isotropic random Gaussian field of small density perturbations. Moreover, the phases of fluctuations are assumed to be random so that the field is entirely defined by its power spectrum (Bower, 1991). The basic idea of the PS theory is to imagine smoothing the cosmological density field at any epoch on a given scale so that the mass scale of virialized objects of interest satisfies . However, as noted earlier, the growth of the density perturbations can only be followed with simple analytical techniques until they become nonlinear. The PS formalism circumvents this difficulty by assuming that the region collapses rapidly and independently of its surroundings once it has turned nonlinear. As a result, the collapsed region can be described as a single large body to the rest of the universe. This simplification allows the linear equations to be applied, however, one must still take into account the nonlinear single body objects when modelling the formation of large-scale structure.

The PS formalism allows the modelling of the growth of cosmic structure in a highly simplified universe (Cole, 1991). Perhaps more importantly, it also allows to estimate the mass function of the collapsed objects. The PS formalism leads to a halo mass function, which has the form of a power law multiplied by an exponential. While the PS formalism gives a reasonable approximation to the numerical data, it has been shown to underestimate the number of massive systems, while over-predicting the number of “typical” mass objects (e.g. Governato et al., 1999; Sheth & Tormen, 1999; Jenkins et al., 2001). Thus, a more realistic and detailed modelling of the nonlinear growth of structure in an expanding universe is required for more vigorous comparisons. This can be achieved, for example, by using numerical simulations.

2.7 Nonlinear evolution: II. Cosmological simulations

2.7.1 Background

Cosmological -body simulations provide a robust method to study large-scale structure and the formation and growth of cosmic structure of a universe on the nonlinear regime. This is possible because the equations of motions are integrated numerically (e.g. Springel et al., 2001, and the references therein). The basic idea of cosmological simulations was founded in the 1960s and it owes its existence to early few body simulations. The first cosmological simulations were run with a modest number of particles and only the formation of a few galaxies were followed. Instead, today large cosmological simulations use billions of particles, however, many simulations still use pure dark matter and no baryons. To overcome the issue that in dark matter only simulations galaxies must be placed by hand using, e.g., semi-analytical methods (Section 2.8.4) hydrodynamics can also be modelled. Hydrodynamical simulations (Katz & Gunn, 1991; Navarro & Benz, 1991; Katz, 1992; Cen, 1992; Rosswog, 2009) that contain gas particles are gaining popularity with increasing computing power, but their volumes are still modest compared to dark matter only simulations (Frenk et al., 1999; Teyssier, 2002; Sommer-Larsen et al., 2003; Springel & Hernquist, 2003; Davé et al., 2010; Razoumov & Sommer-Larsen, 2010; Tantalo et al., 2010). Moreover, even in hydrodynamical simulations some key (sub-grid) physics, such as star formation, is not modelled directly from the first principles but using similar prescriptions as in semi-analytical models.

Observations of the cosmic microwave background, discussed in Section 2.2.1, show that the perturbations of the gravitational potential are caused by non-relativistic material (e.g. Hu & Dodelson, 2002). Therefore, cosmological -body simulations can usually operate on Newtonian limit without the framework of General Relativity. Even so, the expansion of the Universe must be taken into account. This is often done using a co-moving coordinate system, which moves as a function of time as the Universe ages and expands. Consequently, cosmological -body simulation codes essentially follow the evolution of the density field in an expanding background by following the motions of particles caused via gravity by integrating the equations of motions numerically.

The basic assumptions of a modern cosmological -body simulation code can be summarised as follows:

  • mass content of the Universe is build up mainly from dark matter;

  • gravity is the only notable force on large scales;

  • each dark matter particle in the simulation volume represents several particles of the real Universe and they are collisionless;

  • a simulation starts from initial conditions with all modes well within the linear regime;

  • periodic boundary conditions are adopted, resulting in that no particle can be lost during the simulation.

These notions provide the basic assumptions that most simulations obey. It is also noteworthy that the assumption that the dark matter particles are collisionless also means that the evolution of the Universe is driven by the mean gravitational potential rather than two-body interactions. In what follows, I briefly describe a general idea of generating initial conditions and how to follow the motions of particles.

2.7.2 Setting up the initial conditions

For most galaxy formation and large-scale structure problems, setting up the initial conditions of a cosmological -body simulation can be split into three parts:

  1. generating a power spectrum;

  2. generating a Gaussian random density field using the power spectrum;

  3. imposing density perturbation field on the particle distribution.

The first step, generation of a power spectrum, defines the dark matter model. The power spectrum can be generated, for example, by taking a primordial power spectrum (Section 2.5.2) and then multiplying it with the Transfer function (Section 2.5.3) of a chosen cosmology. The second step sets up a “smooth” distribution of particles (for a technical description, see e.g. Martel, 2005) by generating a single realisation of the density field in -space. The third step is to impose density perturbations with the desired characteristics, i.e., the assignment of displacements and velocities to particles. A suitable particle distribution, i.e., a linear fluctuation distribution can be generated using, for example, the Zel’dovich approximation (Section 2.6.2, but see also Efstathiou et al., 1985). Note that when following this technique the matter density and velocity fluctuations are initialised at the starting redshift chosen usually such that all modes in the simulation volume are still within the linear regime.

After the initial conditions have been set and the starting redshift has been chosen, the simulation can be evolved towards the current epoch by using a cosmological -body code. In general, the code allows the time evolution of the simulated particles to be followed by integrating the equations of motions.

2.7.3 Equations of motions

Several different techniques to follow the gravitational evolution of the density field in cosmological -body simulations exist. In this Section the basics are briefly introduced, while more detailed descriptions of different techniques can be found from the literature (see e.g. Efstathiou et al., 1985; Barnes & Hut, 1986; Couchman, 1991; Cen, 1992; Xu, 1995; Kravtsov et al., 1997; Teyssier, 2002; Aarseth, 2003; Springel & Hernquist, 2003; Bagla & Padmanabhan, 1997, and references therein).

In usual cases, dark matter (and stars if applicable) can be modelled as a self-gravitating collisionless fluid in cosmological simulations. Since the number of dark matter particles is large, two-body scattering events are assumed to be seldom. As a result, it is convenient to describe the system in terms of the single particle distribution function in phase space. Now, if we make a reasonable assumption that there are no collisions between particles, the evolution of the distribution function of the fluid follows, in the co-moving coordinates , the collisionless Boltzmann equation:


where the self-consistent potential is the solution of Poisson’s equation


Here is the mass density of the single-particle phase space. Unfortunately though, the coupled equation pair consisting of the collisionless Boltzmann and Poisson equation is difficult to solve directly. Thus, simulations often follow the so-called -body approach, where the smooth phase fluid is represented by particles which are integrated along the characteristic curves of the collisionless Boltzmann equation. Consequently, the problem is conveniently reduced to a task of following Newton’s equations of motion for a large number of particles under their own self-gravity (Springel et al., 2001).

The dynamics of particles can be described by the Hamiltonian:


where and are the co-moving coordinate vectors, the corresponding canonical momenta are given by , and the is the interaction potential. Note that the time dependency of the Hamiltonian is caused by the time dependency in the scale parameter . Before the equations of motions of simulated particles can be derived the interaction potential has to be solved. When periodic boundary conditions are assumed the interaction potential can be solved from equation:


where is the side length of the simulation volume, , and is the particle density distribution function. Finally, after the interaction potential has been solved the Hamilton’s equations of motions


can be derived. To follow the time evolution of the simulated particles the derived equations of motions must be integrated, after making a small variation to time, via e.g. the “leapfrog integration scheme” (e.g. Dolag et al., 2008). Note that in most cases, the particle motion integrals are time-integrals and require integrating the scale parameter .

2.7.4 Identifying dark matter haloes

Far in the nonlinear regime, towards the end of a simulation run, bound structures start to form (for an example, see Fig. 2.6). After their formation they grow in mass either by accretion or by merging with other bound structures. These bound structures can be identified from simulations by using the so-called halo finders that search for collections of dark matter particles that are gravitationally bound. The bound structures of particles are called dark matter haloes due to their relatively spherical nature.

The most popular algorithm to identify virialized haloes is likely the so-called Friends-of-Friends halo finder (Davis et al., 1985). This simple algorithm links all particles with distances less than a linking length to a single halo. In general, the linking length is set to correspond to the mean virialisation overdensity (Eq. 2.0) derived using the spherical top-hat model (Section 2.6.1). The linear linking length of the Friends-of-Friends halo finder is a free parameter, often taken to be a fraction, e.g. , of the mean particle separation. For other halo finding algorithms, see for example Eisenstein & Hut (1998); Neyrinck et al. (2005); Kim & Park (2006); Knollmann & Knebe (2009) and references therein.

Figure 2.6: Example of structure formation as a function of time in the Millennium-II simulation. Courtesy of Michael Boylan-Kolchin and the Millennium-II.
Figure 2.7: A zoom sequence from to Mpc into the most massive halo in the Millennium-II simulation volume at redshift zero. Courtesy of Michael Boylan-Kolchin and the Millennium-II.

Figure 2.7 shows an extraction from the Millennium-II simulation (Boylan-Kolchin et al., 2009). The Figure shows a zoom sequence from to Mpc into the most massive halo in the simulation at redshift zero visualising the dynamical range of the simulation. Figure 2.6 shows the structure formation as a function of time. This set of images shows the growth of the most massive halo over the cosmic time. The left column is Mpc, the centre column is Mpc, and the right is Mpc in co-moving units. From top to bottom, the regions plotted are at redshifts and .

2.7.5 Resolution effects

The simplified scheme of collisionless particles used in cosmological -body simulations leads to finite mass- and force-resolution. In an ideal case the number of simulation particles should be as large as possible to enable a detailed study of formation and growth of cosmic structure in all scales ranging from the smallest dwarf galaxies to the largest clusters and filaments. However, in reality this choice is in general limited by the available computing resources. Moreover, the number of particles must also be balanced with the choice of a simulation volume size to compete against the cosmic variance.

The mass resolution of a simulation can be derived when the simulation volume and the number of dark matter particles have been chosen. The mass resolution, in units of solar mass , describes the mass of a single dark matter particle and can be derived from


Here is the side length of the simulation volume in Mpc, is the number of particles in the volume , and is a constant. Note that the mass resolution sets a hard limit: no object below the mass resolution can form in a simulation. In reality, however, the smallest structures to form and which are identifiable must be made out of several tens of particles. This renders the effective mass resolution at least an order of magnitude worse than indicated by equation 2.0. Furthermore, a finite mass resolution also limits the ability to study the internal structures of dark matter haloes with masses close to the resolution limit.

On the other hand, a finite force resolution arises from the fact that the gravitational force between two particles diverges as their distance approaches to zero. In reality, however, the gravitational force between two extended objects is finite. The force resolution, which is more subtle effect than mass resolution, is more complicated to quantify because it depends on the simulation code used. In a mesh-based simulation code, the force is automatically softened on the scale of a chosen mesh. Instead, in particle-particle algorithms, the force softening is often applied artificially by modifying Newton’s gravity law by writing it as:


where is the mass of a single particle, is the distance between the two particles, and is the gravitational softening length. Non-zero softening length now guarantees that the force does not diverge even when . However, by doing so the simultaneously sets a limit for the highest density contrast that can be resolved.

2.8 Galaxy formation and evolution

2.8.1 Background

So far I have briefly shown how the formation and growth of large scale dark matter structure can be modelled from the CMB down to the current epoch. However, observations such as large galaxy surveys (Section 2.2.2) cannot yet directly probe dark matter, thus one should also try to model the luminous component or baryonic matter, i.e. the galaxy formation and evolution. A detailed discussion of the theory of galaxy formation, is however beyond the scope of this introduction. Instead, in what follows I try briefly to summarise the main concepts of galaxy formation, illustrated in Figure 2.8, to provide background for the following chapters. For more detailed presentations, I refer the interested reader to the great text books of Longair (2008) and Mo et al. (2010).

The previous sections showed how to model dark matter haloes. If we now assume that galaxies form and reside in dark matter haloes, it becomes obvious that the properties of the galaxy population are related to the cosmological density field and to the dark matter halo population. One can therefore try to link the properties of dark matter haloes to the properties of observed galaxies by using statistical arguments. As a result, it can be shown that the properties of the galaxy population depend on the properties of the dark matter halo and subhalo populations (see e.g. Mo et al., 2010, for detailed discussion). This allows, for example, the galaxy luminosity function to be compared to the dark matter halo mass function. The correspondence of light and mass is important, for example, when the mass power spectrum is being derived from observations (Section 2.2.2), because a correlation between the observed light and the underlying mass must be assumed.

2.8.2 Linking halo mass to galaxy luminosity

To overcome the difficulty of linking dark matter haloes to luminous galaxies, Vale & Ostriker (2004) (see also Oguri, 2006) proposed that a galaxy’s luminosity can be related to its host dark matter halo’s virial mass via a simple relation as follows:


where , , , and are free parameters and is scaled such that


However, Cooray & Milosavljević (2005) showed that the relation between the mass of a dark matter halo and the luminosity of the galaxy it hosts is not straightforward due to the complicated baryonic physics involved.

The baryonic content of a dark matter halo becomes dynamically important in the nonlinear regime when dark matter haloes are forming. Consequently, to model realistic galaxies, such as the Galaxy we live in, baryons must be modelled, albeit they do not contribute to the structure formation as much as dark matter. Hydrodynamical effects such as heating and cooling processes of gas, shocks, star formation, and feedback processes must all be taken into account when the formation of baryonic structure is being considered. One solution is to use hydrodynamical simulations to model gas directly, but the lack of fundamental theories for physical processes involved in the formation and evolution of galaxies, such as star formation, render them less than optimal. As a result, the hydro-simulations also require simple prescriptions for the so-called sub-grid physics, which cannot be modelled directly. What can be modelled then?

2.8.3 Modelling of galaxy formation and evolution

Galaxy formation (for a comprehensive view, see e.g. Longair, 2008; Mo et al., 2010) is expected to proceed via a two-stage process originally outlined by White & Rees (1978), but see also Hoyle (1953); Binney (1977); Rees & Ostriker (1977); Silk (1977) for early development. In this paradigm, the gravitational instability acting on collisionless dark matter results in the formation of self-gravitating dark matter haloes (as already noted earlier). Because baryons, initially well mixed with the dark matter, are assumed to “feel” the dark matter via gravity, they also participate in this collapse after the dark matter haloes have started to form. However, unlike the dark matter, the gas is not collisionless, but can dissipate. As a result, in a very simplified picture, the gas can be assumed to be heated by shocks to the virial temperature of the dark matter halo during this infall. After which, the hot gas can cool radiatively, on a time scale set by atomic physics.

During the collapse and cooling, the gas is assumed to condense to the cores of collapsed dark matter haloes. However, it is assumed that this process is not the same in all haloes. In smaller structures such as galaxy host haloes, the dominant physical process is cooling, which allows baryons to be more centrally concentrated than dark matter. In contrast, in larger structures, baryons experience a deeper gravitational potential and can therefore gain potential energy as they fall to the centre of a halo. This process heats baryonic matter and increases its temperature via shocks (e.g. Birnboim & Dekel, 2003). As a result, baryonic matter experiences pressure forces which do not let them to be as concentrated as its host dark matter. In large-scale structures the process of cooling will therefore be highly important in the dense cores of dark matter where the gas can cool. Finally, the cold gas can fragment into stars, and a galaxy is born.

To simplify, in hierarchical models, such as the CDM (Section 2.4.1), the galaxy formation involves at minimum the following three stages:

  1. the hierarchical formation of dark matter haloes,

  2. the accretion of gas into the haloes, and

  3. the cooling and fragmentation of the hot gas into stars.

Figure 2.8 shows a logic flow chart for galaxy formation131313The flow chart is by no means a complete description of all gas physics that may play a role in galaxy formation, but rather tries to capture the main aspects.. The paths leading to the formation of various galaxies (green ellipses) are drawn from the initial conditions set by the cosmological framework, discussed in earlier Sections and outlined in Figure 2.5. Note that the flow chart does not include any feedback effects, which have been found to be significant and will be discussed later.

Figure 2.8: A logic flow chart for galaxy formation. The dominantly gravitational processes have been coloured with blue, while mainly hydrodynamical processes are coloured with yellow. The red hexagons indicate a binary choice, while dark green ellipses note the end products of the simplified processes described in the flow chart. The figure is an adaptation of Fig. 1.1 from Mo et al. (2010).

2.8.4 Semi-analytical models of galaxy formation

Due to the complicated physics related to galaxy formation and evolution (as seen in Fig. 2.8), simple rules that can be easily varied to study the importance of different physical processes are highly useful. Semi-analytical models (SAMs) of galaxy formation try to fill this void by encoding simplistic rules for the formation and evolution of galaxies within a cosmological framework. A SAM is a collection of physical recipes that describe an inflow of gas, how gas can cool and heat up again, how stars are formed within galaxies, how stellar populations evolve and how black holes grow using simplified physics (Cole, 1991; White & Frenk, 1991; Lacey & Silk, 1991; Kang et al., 2005; Baugh, 2006; Lucia & Blaizot, 2007). SAMs can also easily include different feedback effects: stellar winds, active galactic nuclei (AGN) or supernovae (SNe) feedbacks (e.g. Croton et al., 2006; Somerville et al., 2008; Ricciardelli & Franceschini, 2010), for example. Hence, SAMs try to describe all the gas physics that goes into galaxy formation and evolution, but is not modelled in the dark matter only simulation.

Due to their nature SAMs can be used to explore ideas of galaxy formation and evolution and to understand which physical processes are the most important in the life of a galaxy by changing the recipes describing the physics. SAMs can also be applied to the so-called sub-grid physics that operates below the resolution of a simulation. As it is not yet possible to simulate all star formation processes (McKee & Ostriker, 2007; Krumholz, 2011) in a cosmological context, sub-grid physics must be modelled with simplified physics even when hydrodynamics is involved.

The backbone of a SAM is the evolution of dark matter haloes. Often, this evolution is parameterised with dark matter halo merger trees (Fig. 2.9) that allow the hierarchical nature of gravitational instabilities to be explicitly taken into account (Baugh et al., 1998). Dark matter merger trees describe how the dark matter haloes form via mergers of smaller haloes. They provide the backdrop for the introduction of the baryonic component which reacts gravitationally to the growing network of dark matter potential wells. Even though modern studies derive merger trees directly from simulations, this is by means not necessary as they can be derived also by using Monte Carlo techniques. In such case the extended Press-Schechter (EPS) theory (Bond et al., 1991; Lacey & Cole, 1993) is often used. The Monte Carlo techniques provide a fast method of generating merger trees, however, they have been found to be less than reliable in some cases (e.g. Cole et al., 2008).

Figure 2.9: Example of a merger tree as a function of the lookback time. Symbols are colour-coded as a function of the colour and their area scales with the stellar mass. Only progenitors more massive than M are shown.

After a cosmological model has been chosen and the merger trees have been generated the baryonic processes must be taken into account. A SAM typically consists of the following steps: 1) follow the three baryonic components: hot and cold gas, and stars and adopt a recipe for disk formation; 2) specify a recipe for the conversion rate between the three components, including star formation and feedback effects; 3) keep track of the metallicity of each component; 4) convert the star formation history and metallicity of the stellar populations into luminosities; and 5) adopt a recipe for galaxy-galaxy mergers. In what follows, I briefly describe these steps in turn.

Heating and cooling of gas

Modern galaxy formation theories assume that the gas density profile follows that of dark matter. When gas falls into the potential well of a dark matter halo it is assumed to be shock-heated (White & Rees, 1978) to the virial temperature of the halo, given by


Here is the mean molecular mass of the gas, is the mass of a hydrogen atom, is the Boltzmann’s constant, is the mass of the halo and is the virial radius, within which the mean density is 200 times the critical density.

Before stars can form from the cold K molecular clouds the shock-heated hot gas must cool radiatively (e.g. Helly et al., 2003). The cooling time can be defined, for example, as the ratio of the thermal energy density and the cooling rate per unit volume. In this case, the cooling time can be written as


where is the mean particle mass, refers to the Boltzmann constant, is the hot gas density, and is the cooling function. Note that the cooling rate is a function of metallicity of the gas and the virial temperature of the halo , and thus, the cooling is in general more effective in higher density regions. Additionally, in highly simplified scenarios, the more metal-enriched gas tends to cool faster.

Different cooling mechanisms (inverse Compton scatter, molecular and atomic cooling, and bremsstrahlung) can however be dominant at different times and temperatures complicating the gas cooling modelling significantly. For example, in massive haloes, where K, gas is fully collisionally ionised and cools mainly through bremsstrahlung emission from free electrons. In the temperature range K K excitation and de-excitation mechanisms dominate, while in haloes with K gas is mainly neutral and the cooling processes are therefore suppressed. However, in a simplified scenario the gas is assumed to be able to cool if cooling time is shorter than some characteristic timescale, which is model dependent.

As the cold condensed gas accumulates in the central regions of the dark matter haloes it can be identified as the ISM of the protogalaxy. What is assumed to follow after the cold gas has settled down is a disk formation.

Disk formation

The disk formation may be modelled after, for example, the amount of cold gas and the properties of the host halo. The underlying physical notion is that when structures grow and collapse in the early universe, they exert tidal torques on each other. This provides each collapsing dark matter halo with some angular momentum. As a result, SAMs typically assume that the cold gas will form a disk with the same specific angular momentum as the dark matter halo, while the size and rotation velocity are determined by the spin parameter:


Here is the angular momentum of the halo, refers to the total energy of the halo, and is the mass of the dark matter halo. When adopting this simple prescription, what follows is that the mass and spin of the disk is tightly coupled to the mass and spin of the dark matter halo.

The dimensionless angular momentum is a measure of the degree of rotational support of the galaxy. Note, however, that the typical values of of collapsed dark matter haloes have been found to be significantly smaller than that of the largely flattened centrifugally supported disk galaxies we observe today with . Hence, a considerable amount of dissipation must have occurred to produce the observed disks.

Star formation

In a typical SAM star formation (for a general review, see Kennicutt, 1998a) is assumed to take place in the disks of galaxies, while the actual onset of star formation is assumed to occur once the surface density of cold gas exceeds a critical density (Kennicutt, 1998b, 1989) in a molecular cloud (Krumholz, 2011). Ideally, the star formation law should be derived from the first principles as a function of the physical conditions, such as density, temperature, metallicity, radiation and magnetic fields of the ISM, however, the detailed physics involved in the fragmentation of the cold gas, collapse and onset of a protostar, and the physical conditions of the ISM are not yet well understood. I will, however, return to the importance of the ISM in Section 4.5.

Due to the complicated physics involved, SAMs often derive the star formation rate (SFR) of a galaxy using a simple relation


where is the density of the cold gas, is the characteristic timescale, and is a measure for the efficiency of star formation. Note, however, that several different forms of the above equation have been developed. These recipes for the star formation efficiency range from simple models that assume that is a constant to models that are proportional to the dynamic time of the galaxy and take into account, for example, the circular velocity and/or the radius of the disc. In any event, the above star formation law is a variant of the empirical Schmidt (1959) law, for which it has been assumed that the SFR is controlled by the self-gravity of the gas.

Closely tied to the star formation in galaxies is the number of stars of a given mass that forms, that is the the initial mass function (IMF; e.g. Kroupa, 2001; Chabrier, 2003). SAMs typically assume that the IMF of stellar populations is universal when modelling star formation. Note, however, that theoretical arguments (e.g. Davé, 2008) and indirect observational evidence suggest that the stellar IMF may evolve with, e.g., time (van Dokkum, 2008, and references therein) or environment, casting a shadow over the assumption of universality.

Feedback effects

The early SAMs, using recipes similar to the above ones, were however, not able to reproduce the observed form of the galaxy luminosity function (LF; Eq. 3.0). They often over-predicted the number of both faint and bright galaxies, especially in the infrared (e.g. Benson et al., 2003; Croton et al., 2006; Benson & Devereux, 2010). To alleviate the discrepancy i.e. to limit the number of faint and bright galaxies a feedback mechanism was introduced. To regulate the star formation in both light and massive dark matter haloes the feedback was divided to two separate mechanisms that operate on different mass regimes.

Modern SAMs model active galactic nuclei (AGN) feedback, which can suppress the cooling flow in high mass systems (e.g. Silk & Rees, 1998; Croton & Farrar, 2008; del P Lagos et al., 2008). The AGN feedback (Sijacki et al., 2007) provides additional energy that can suppress the cooling of hot gas generating a sharp cut-off to the high-luminosity end of the LF. In many models the strength of the AGN feedback depends directly on the mass accretion of the black hole , thus, the modified cooling rate can be written, for example, as


Here, refers to the black hole accretion efficiency. The additional energy from the AGN can prevent gas from cooling and is more important in later times when galaxies have more massive black holes. The AGN feedback can therefore help to regulate the SFR at later epochs and to prevent the overproduction of very massive galaxies. Another effect of AGN feedback is seen on the ages of high stellar mass systems, which are significantly older for AGN feedback models (Khalatyan et al., 2008).

While the AGN feedback affects mainly massive galaxies, the supernovae (SN) feedback is important for lighter galaxies and for the metal enrichment of the inter-stellar and possibly even the inter-galactic medium (Vecchia & Schaye, 2008). The SN feedback helps to self-regulate the process of star formation throughout the galaxies’ history. In the absence of the SN feedback star-formation rates (SFRs) are extremely high in early times and fairly low in more recent epochs. However, with the SN feedback, SFRs are initially lower so more gas is available for later periods of star formation (see also Section 4.2.3).

The SN feedback can function because in the models it is assumed that a SN can blow gas out of a star forming disc. Moreover, it is also assumed that the rate of mass ejection is proportional to the total mass of stars formed. As a result, the re-heating of gas due to the feedback can then be modelled as follows