Abstract
The cosmic web structure is studied with the concepts and methods of fractal geometry, employing the adhesion model of cosmological dynamics as a basic reference. The structures of matter clusters and cosmic voids in cosmological N-body simulations or the Sloan Digital Sky Survey are elucidated by means of multifractal geometry. A nonlacunar multifractal geometry can encompass three fundamental descriptions of the cosmic structure, namely, the web structure, hierarchical clustering, and halo distributions. Furthermore, it explains our present knowledge of cosmic voids. In this way, a unified theory of the large-scale structure of the universe seems to emerge. The multifractal spectrum that we obtain significantly differs from the one of the adhesion model and conforms better to the laws of gravity. The formation of the cosmic web is best modeled as a type of turbulent dynamics, generalizing the known methods of Burgers turbulence.
1. Introduction
The evolution of the universe at large is ruled by gravity. Although the Einstein equations of the theory of general relativity are nonlinear and difficult to solve, the early evolution of the universe can be described by an exact solution, namely, the FLRW solution, plus its linear perturbations. These perturbations grow, and they grow faster on smaller scales, becoming nonlinear. Then, large overdensities arise that decouple from the global expansion. This signals the end of a dust-like description of the matter dynamics and the need for a finer description of it, in which dissipative processes, in particular, play a prominent role. The dissipation is essentially a transfer of kinetic energy from one scale to another smaller scale, due to nonlinear mode coupling, as is characteristic of fluid turbulence.
Indeed, the successful adhesion model of large-scale structure formation is actually a model of turbulence in irrotational pressureless flow [1, 2]. This model generates the well-known cosmic web structure. The geometry and formation of this web structure are the subject of this study. The cosmic web is a foam-like structure, formed by a web of sheets surrounding voids of multiple sizes. In fact, the range of sizes is so large that we can speak of a self-similar structure. This motivates its study by means of fractal geometry.
In fact, fractal models of the universe predate the discovery of the cosmic web structure and arose from the idea of a hierarchy of galaxy clusters that continues indefinitely towards the largest scales, an idea championed by Mandelbrot [3]. In spite of the work of many cosmologists along several decades, e.g., Peebles [4, 5], the debate about the scale of transition to homogeneity is not fully settled [4–10]. However, the magnitude of the final value of the scale of homogeneity shall not diminish the importance of the fractal structure on lower scales, namely, in the highly nonlinear clustering regime [11]. At any rate, the notion of a hierarchy of disconnected matter clusters has turned out to be naive. Mandelbrot [3] already considered the early signs of cosmic voids and of filamentary structure in the universe and, with this motivation, he proposed the study of fractal texture. However, only the study of the cosmological dynamics of structure formation can unveil the precise nature of the intricate cosmic web geometry.
The geometry of the cosmic web has been observed in galaxy surveys [12–14] and cosmological -body simulations [15–17]. In addition, the fractal geometry of the large-scale structure has been also studied in the distribution of galaxies [18–20] and in cosmological -body simulations [21–23]. Recent -body simulations have very good resolution and they clearly show both the morphological and self-similar aspects of the cosmic web. However, while there seems to be a consensus about the reality and importance of this type of structure, a comprehensive geometrical model is still missing. In particular, although our present knowledge supports a multifractal model, its spectrum of dimensions is only partially known, as will be discussed in this paper.
We must also consider halo models [24], as models of large-scale structure that are constructed from statistical rather than geometrical principles. In the basic halo model, the matter distribution is separated into a distribution within “halos”, with given density profiles, and a distribution of halo centers in space. Halo models have gained popularity as models suitable for analyzing the results of cosmological -body simulations and can even be employed for designing simulations [25]. Halo models are also employed for the study of the small-scale problems of the Lambda cold dark matter (LCDM) cosmology [26]. Remarkably, the basic halo model, with some modifications, can be combined with fractal models [27–31]. At any rate, halo models have questionable aspects and a careful analysis of the mass distribution within halos shows that it is too influenced by discreteness effects intrinsic to -body simulations [31, 32].
We review in this paper a relevant portion of the efforts to understand the cosmic structure, with a bias towards methods of fractal geometry, as regards the description of the structure, and towards nonlinear methods of the theory of turbulence, as regards the formation of the structure. So this review does not have as broad a scope as, for example, Sahni and Coles’ review of nonlinear gravitational clustering [33], in which one can find information about other topics, such as the Press-Schechter approach, the Voronoi foam model, and the BBGKY statistical hierarchy. General aspects of the cosmic web geometry have been reviewed by van de Weygaert [34]. Here, we shall initially focus on the adhesion approximation, as a convenient model of the cosmic web geometry and of the formation of structures according to methods of the theory of fluid turbulence, namely, the theory of Burgers turbulence in pressureless fluids [35, 36]. However, the adhesion model is, as we shall see, somewhat simplistic for the geometry as well as for the dynamics. For the former, we need more sophisticated fractal models, and, for the latter, we have to consider, first, a stochastic version of the adhesion model and, eventually, the full nonlinear gravitational dynamics.
2. The Cosmic Web Geometry
In this section, we study the geometry of the universe on middle to large scales, in the present epoch, and from a descriptive point of view. As already mentioned, there are three main paradigms in this respect: the cosmic web that arises from the Zeldovich approximation and the adhesion model [1, 2], the fractal model [3], and the halo model [24]. These paradigms have different origins and motivations and have often been in conflict. Notably, a fractal model with no transition to homogeneity, as Mandelbrot proposed [3], is in conflict with the standard FLRW cosmology and, hence, with the other two paradigms. Indeed, Peebles’ cosmology textbook [5] places the description of “Fractal Universe” in the chapter of Alternative Cosmologies. However, a fractal nonlinear regime in the FLRW cosmology does not generate any conflict with the standard cosmology theory [11] and actually constitutes an important aspect of the cosmic web geometry. Furthermore, we propose that the three paradigms can be unified and that they represent, from a mature point of view, just different aspects of the geometry of the universe on middle to large scales. Here we present the three paradigms, separately but attempting to highlight the many connections between them and suggesting how to achieve a unified picture.
2.1. Self-Similarity of the Cosmic Web
Although the cosmic web structure has an obviously self-similar appearance, this aspect of it was not initially realized and it was instead assumed that there is a “cellular structure” with a limited range of cell sizes [37]. However, the web structure is generated by Burgers turbulence, and in the study of turbulence and, in particular, of Burgers turbulence, it is natural to assume the self-similarity of velocity correlation functions [35, 36]. Therefore, it is natural to look for scaling in the cosmic web structure. The study of gravitational clustering with methods of the theory of turbulence is left for Section 3, and we now introduce some general notions, useful to understand the geometrical features of the cosmic web and to connect with the fractal and halo models.
The FLRW solution is unstable because certain perturbations grow. This growth means, for dust (pressureless) matter following geodesics, that close geodesics either diverge or converge. In the latter case, the geodesics eventually join and matter concentrates on caustic surfaces. In fact, from a general relativity standpoint, caustics are just a geometrical phenomenon that always appears in the attempt to construct a synchronous reference frame from a three-dimensional Cauchy surface and a family of geodesics orthogonal to this surface [38]. But what are, in general, nonreal singularities, owing to the choice of coordinates, become real density singularities in the irrotational flow of dust matter.
One can consider the study and classification of the possible types of caustics as a purely geometrical problem, connected with a broad class of problems in the theory of singularities or catastrophes [39]. To be specific, the caustic surfaces are singularities of potential flows and are called Lagrangian singularities [39, ch. 9]. A caustic is a critical point of the Lagrangian map, which maps initial to final positions. After a caustic forms, the inverse Lagrangian map becomes multivalued; that is to say, the flow has several streams (initially, just three). Although this is a general process, a detailed study of the formation of singularities in cosmology has only been carried out in the Zeldovich approximation, which is a Newtonian and quasi-linear approximation of the dynamics [1, 33]. Recently, Hidding, Shandarin, and van de Weygaert [40] have described the geometry of all generic singularities formed in the Zeldovich approximation, displaying some useful graphics that show the patterns of folding, shell crossing, and multistreaming. Since the classification of Lagrangian singularities is universal and so is the resulting geometry, those patterns are not restricted to the Zeldovich approximation.
Unfortunately, the real situation in cosmology is more complicated: the initial condition is not compatible with a smooth flow, so the theory of singularities of smooth maps can only be employed as an approximation, suitable for smoothened initial conditions. This approach yields some results [33]. However, a nonsmooth initial condition gives rise to a random distribution of caustics of all sizes, with an extremely complex distribution of multistreaming flows, whose geometry is mostly unexplored. Besides, the Zeldovich approximation fails after shell crossing. The real multistreaming flow generated by Newtonian gravity has been studied in simplified situations, for example, in one dimension (which can describe the dimension transverse to caustics). Early -body simulations showed that the thickness of an isolated multistreaming zone grows slowly, as if gravity makes particles stick together [1]. The analytical treatment of Gurevich and Zybin [41], for smooth initial conditions, proves that multistreaming gives rise to power-law mass concentrations. One-dimensional simulations with cosmological initial conditions (uniform density and random velocities) [42–44] show that particles tend to concentrate in narrow multistreaming zones pervading the full spatial domain. Furthermore, both the space and phase-space distributions tend to have self-similar properties.
At any rate, before any deep study of multistreaming was undertaken, the notion of “gravitational sticking” of matter suggested a simple modification of the Zeldovich approximation that suppresses multistreaming. It is natural to replace the collisionless particles by small volume elements and hence to assume that, where they concentrate and produce an infinite density, it is necessary to take into account lower-scale processes. The simplest way to achieve this is embodied in the adhesion model, which supplements the Zeldovich approximation with a small viscosity, giving rise to the (three-dimensional) Burgers equation [1]. One must not identify this viscosity with the ordinary gas viscosity but with an effect of coarse graining of the collisionless Newtonian dynamics [45]. Even a vanishing viscosity is sufficient to prevent multistreaming. Actually, the vanishing viscosity limit is, in dimensionless variables, the high Reynolds-number limit, which gives rise to Burgers turbulence. Although we postpone the study of cosmic turbulence to Section 3, we summarize here some pertinent results of the adhesion model.
A deep study of the geometry of the mass distribution generated by the adhesion model has been carried out by Vergassola et al. [46]. Their motivation was to compare the prediction of the adhesion model with the result of the Press-Schechter approach to the mass function of collapsed objects (this approach is explained in [33]). But Vergassola et al. actually show how a self-similar mass distribution arises (including the analytical proof for a particular case). The adhesion model (with vanishing viscosity) has an exact geometrical solution, in terms of the convex hull of the Lagrangian potential [46]. In one dimension and with cosmological initial conditions, the adhesion process transforms the Lagrangian map to a random devil’s staircase, with an infinite number of steps, whose lengths correspond to the amounts of stuck mass (Figure 1, left-hand side). To be precise, if the initial condition consists of a uniform density and a Gaussian velocity field with a power-law energy spectrum, then, at any time , singularities arise in the velocity field and mass condensates on them, giving rise to the steps of the devil’s staircase. Such condensations initially contain very small mass but grow as mass sticking proceeds, preserving a self-similar distribution in space and time.
Let us explain in more detail the properties of the self-similar mass distribution that arises from a scale-invariant initial velocity field, such that namely, a fractional Brownian field with Hurst exponent [46]. These fields are Gaussian random fields that generalize the Brownian field, with and independent increments, in such a way that, for , same-sign increments are likely (persistence), whereas, for , opposite-sign increments are likely (antipersistence) [3, ch. IX]. If , then the field tends to be smooth, whereas if , then the field tends to be discontinuous. A fractional Brownian velocity field corresponds to Gaussian density fluctuations with power-law Fourier spectrum of exponent , in dimensions. There is a scale , called “coalescence length” by Vergassola et al. [46], such that it separates small scales, where mass condensation has given rise to an inhomogeneous distribution, from large scales, where the initial conditions hold and the distribution is still practically homogeneous ( is a redefined time variable, namely, the growth factor of perturbations [1, 46]). In one dimension, we have a devil’s staircase with steps that follow a power-law, so that the cumulative mass function isHowever, the number of steps longer than , that is to say, of large mass condensations, has an exponential decay (like in the Press-Schechter approach). These results can be generalized to higher dimensions, in which the web structure is manifest (Figure 1, right-hand side). However, in higher dimensions, the full geometrical construction is very complex and hard to visualize [47].
Even though the results of the adhesion model are very appealing, we must notice its shortcomings. On the one hand, it is based on the Zeldovich approximation, which is so simple that the mass function is directly given by the initial energy spectrum. The Zeldovich approximation is not arbitrary and is actually the first order of the Lagrangian perturbation theory of collisionless Newtonian dynamics, as thoroughly studied by Buchert (e.g., [48]); but some nice properties of the Zeldovich approximation are not preserved in higher orders. At any rate, the Zeldovich approximation is exact in one dimension, yet the one-dimensional adhesion model is not, because the instantaneous adhesion of matter is not realistic and differs from what is obtained in calculations or -body simulations, namely, power-law mass concentrations instead of zero-size mass condensations [41–44, 49]. The calculations and simulations also show that the final mass distribution is not simply related to the initial conditions.
Nevertheless, it can be asserted that the locations of power-law mass concentrations in cosmological -body simulations, in one or more dimensions, form a self-similar distribution that looks like the cosmic web of the adhesion model. Indeed, all these mass distributions are actually particular types of multifractal distributions [27–30, 43]. We now study fractal models of large-scale structure from a general standpoint.
2.2. Fractal Geometry of Cosmic Structure
The concept of cosmic web structure has its origin in the adhesion model, which is based on the study of cosmological dynamics. In contrast, fractal models of large-scale structure have mostly observational origins, combined with ancient theoretical motivations (the solution of the Olbers paradox, the idea of nested universes, etc.). Actually, the ideas about the structure of the universe have swung between the principle of homogeneity and the principle of hierarchical structure [5]. The cosmological principle is the modern formulation of the former, but it admits a conditional formulation that supports the latter [3, §22]. There is no real antithesis and the synthesis is provided by the understanding of how a fractal nonlinear regime develops below some homogeneity scale, consequent to the instability of the FLRW model. This synthesis is already present in the adhesion model, in which the homogeneity scale is the above-mentioned coalescence length.
The notion of hierarchical clustering, namely, the idea of clusters of galaxies that are also clustered in superclusters and so onwards, inspired the early fractal models of the universe. This idea can be simply realized in some three-dimensional and random generalization of the Cantor set, with the adequate fractal dimension [3]. In mathematical terms, the model places the focus on the geometry of fractal sets, characterized by the Hausdorff dimension. A random fractal set model for the distribution of galaxies has been supported by measures of the reduced two-point correlation function of galaxies, which is long known to be well approximated as the power-lawwith exponent and correlation length Mpc, as standard but perhaps questionable values [3–10]. In fact, the name correlation length for is essentially a misnomer [9], and is to be identified with the present-time value of in the adhesion model, that is to say, with the current scale of transition to homogeneity. Although the value of has been very debated, (3) implies, for any value of , that there are very large fluctuations of the density when and indeed that the distribution is fractal with Hausdorff dimension [3–11].
It would be convenient that fractals were described only by its Hausdorff dimension, but it was soon evident that one sole number cannot fully characterize the richness of complex fractal structures and, in particular, of the observed large-scale structure of the universe. Indeed, fractal sets with the same Hausdorff dimension can look very different and should be distinguishable from a purely mathematical standpoint. This is obvious in two or higher dimensions, where arbitrary subsets possess topological invariants such as connectivity that are certainly relevant. Furthermore, it is also true in one dimension, in spite of the fact that any (compact) fractal subset of is totally disconnected and has trivial topology.
In an effort to complement the concept of fractal dimension, Mandelbrot was inspired by the appearance of galaxy clusters and superclusters and introduced the notion of lacunarity, without providing a mathematical definition of it [3, §34]. Lacunarity is loosely defined as the property of having large gaps or voids (lacuna). This property gained importance in cosmology with the discovery of large voids in galaxy redshift surveys, which came to be considered a characteristic feature of the large-scale structure of the universe [13, 14]. The precise nature of cosmic voids is indeed important to characterize the geometry of the cosmic web and will be studied below.
In general, for sets of a given Hausdorff dimension, there should be a number of parameters that specify their appearance, what Mandelbrot calls texture [3, ch. X]. Lacunarity is the first one and can be applied to sets in any dimension, including one dimension. Others are only applicable in higher dimension; for example, parameters to measure the existence and extension of subsets of a fractal set that are topologically equivalent to line segments, an issue related to percolation [3, §34 and 35].
In fact, the methods of percolation theory have been employed in the study of the cosmic web [1]. Those methods are related to the mathematical notion of path connection of sets, which is just one topological property among many others. Topological properties are certainly useful, but topology does not discriminate enough in fractal geometry, as shown by the fact that all (compact) fractal subsets of are topologically equivalent. According to Falconer [50], fractal geometry can be defined in analogy with topology, replacing bicontinuous transformations by a subset of them, namely, bilipschitz transformations. The Hausdorff dimension is invariant under these transformations, although it is not invariant under arbitrary bicontinuous transformations. Therefore, the set of parameters characterizing a fractal set includes, on the one hand, topological invariants and, on the other hand, properly fractal parameters, beginning with the Hausdorff dimension. Unfortunately, there has been little progress in the definition of these parameters, beyond Mandelbrot’s heuristic definition of texture parameters.
Besides the need for a further development of the geometry of fractal sets, we must realize that the concept of fractal set is not sufficient to deal with the complexity of cosmic geometry. This is easily perceived by reconsidering the mass distribution generated by the adhesion model, in particular, in one dimension. In the geometry of a devil’s staircase, we must distinguish the set of lengths of the steps from the locations of these steps, which are their ordinates in the graph on the left-hand side of Figure 1. These locations constitute a dense set; that is to say, they leave no interval empty [46]. The set is, nevertheless, a countable set and therefore has zero Hausdorff dimension. But more important than the locations of mass condensations is the magnitude of these masses (the lengths of the steps). To wit, we are not just dealing with the set of locations but with the full mass distribution. Mandelbrot [3] does not emphasize the distinction between these two concepts and actually gives no definition of mass distribution; but other authors do; for example, Falconer includes an introduction to measure theory [50], which is useful to precisely define the Hausdorff dimension as well as to define the concept of mass distribution. A mass distribution is a finite measure on a bounded set of , a definition that requires some mathematical background [50], although the intuitive notion is nonetheless useful.
The geometry of a generic mass distribution on is easy to picture, because we can just consider the geometry of the devil’s staircase and generalize it. It is not difficult to see that the cumulative mass that corresponds to a devil’s staircase interval is given by switching its axes; for example, the total mass in the interval is given by the length of the initial interval that transforms to under the Lagrangian map. The switching of axes can be performed mentally in Figure 1 (left). The resulting cumulative mass function is monotonic but not continuous. A general cumulative mass function needs only to be monotonic, and it may or may not be continuous: its discontinuities represent mass condensations (of zero size). Given the cumulative mass function, we could derive the mass density by differentiation, in principle. This operation should present no problems, because a monotonic function is almost everywhere differentiable. However, the inverse devil’s staircase has a dense set of discontinuities, where it is certainly nondifferentiable. This does not constitute a contradiction, because the set of positions of discontinuities is countable and therefore has zero (Lebesgue) measure, but it shows the complexity of the situation. Moreover, in the complementary set (of full measure), where the inverse devil’s staircase is differentiable, the derivative vanishes [46].
In summary, we have found that our first model of cosmological mass distribution is such that the mass density is either zero or infinity. We might attribute this singular behavior to the presence of a dense set of zero-size mass condensations, due to the adhesion approximation. Therefore, we could expect that a generic mass distribution should be nonsingular. Of course, the meaning of “generic” is indefinite. From a mathematical standpoint, a generic mass distribution is one selected at random according to a natural probability distribution (a probability distribution on the space of mass distributions!). This apparently abstract problem is of interest in probability theory and has been studied, with the result that the standard methods of randomly generating mass distributions indeed produce strictly singular distributions [51]. This type of distributions has no positive finite density anywhere, like our discontinuous example, but the result also includes continuous distributions. In fact, the mass distributions obtained in cosmological calculations or simulations [41, 43] seem to be continuous but strictly singular nonetheless; namely, they seem to contain dense sets of power-law mass concentrations.
Let us consider the example of random mass function displayed in Figure 2, which is generated with one of the methods referred to in [51]. Actually, the real mass distribution, namely, the mass density, is obtained as the derivative of the mass function and is not a proper function, because it takes only the values zero or infinity. The graphs in Figure 2 are the result of coarse graining with a small length, but it is evident that the values of the density are either very large or very small. Furthermore, the large or small values are clustered, respectively. They form mass clusters or voids, the former occupying vanishing total length and the latter occupying the full total length, in the limit of vanishing coarse-graining length. Unlike the inverse devil’s staircase, the present mass distribution function is continuous, so its graph contains no vertical segments, in spite of being very steep and actually nondifferentiable at many points. Therefore, there are no zero-size mass condensations but just mass concentrations. These mass concentrations are clustered (in our example, as a consequence of the generating process). In view of its characteristics, this random mass function illustrates the type of one-dimensional mass distributions in cosmology, which can be considered generic. Some further properties of generic strictly singular mass distributions are worth noticing.
To any mass distribution is associated a particular set, namely, its supporting set, which is the smallest closed set that contains all the mass [50]. This set is the full interval in the examples of Figure 1 or Figure 2; but it can be smaller in strictly singular mass distributions and it can actually be a fractal set, namely, a set with Hausdorff dimension strictly smaller than one. In fact, if Mandelbrot [3] seldom distinguishes between fractal sets and fractal mass distributions, it is probably because, in most of his examples, the support is a fractal set and the mass is uniformly distributed on it. In this simple situation, the distinction is superfluous. However, we have to consider singular mass distributions in general for the cosmic geometry and take into account the fact that the observed clustering may not be a property of the supporting set but a property of the mass distribution on it. This aspect of clustering is directly related to the nature of voids: in a fractal set, the set of voids is simply the complementary set, which is totally empty, whereas the voids in a mass distribution with full support have tiny mass concentrations inside (see Figures 1 and 2). The latter option describes better our present knowledge of cosmic voids [28, 52, 53].
In conclusion, the appropriate geometrical setting for the study of the large-scale structure of the universe is the geometry of mass distributions, as a generalization of the geometry of sets. Naturally, a generic mass distribution is more complex in three-dimensional space than in one-dimensional space. In particular, one has to consider properties such as connectivity, which are important for the web structure (Figure 1). Since the supporting set of a geometrical model of cosmic structure is probably trivial, namely, the full volume, one cannot formulate the morphological properties in terms of the topology of that set. So it is necessary to properly define the topology of a mass distribution. This can be done in terms of the topology of sets of particular mass concentrations but is beyond the scope of this work. The classification of mass concentrations is certainly useful, as a previous step and by itself. The standard way of classifying mass concentrations constitutes the multifractal analysis of measures [50, ch. 17].
2.3. Multifractal Analysis of Mass Distributions
The original idea of hierarchical clustering of galaxies considers galaxies as equivalent entities and constructs a sort of fractal set as the prolongation of the hierarchy towards higher scales (up to the scale of homogeneity) and towards lower scales (the prolongation to lower and indeed infinitesimal scales is necessary to actually define a mathematical fractal set). The chief observable is the fractal dimension, which has been generally obtained through (3), namely, from the value of estimated through the correlation function of galaxy positions [3–6]. By using the galaxy positions, the density in (3) is the galaxy number density instead of the real mass density. Both densities should give similar results, if the mass were distributed uniformly over the galaxy positions. However, Pietronero [18] noticed that galaxies are not equivalent to one another and that, in fact, their masses span a broad range. Therefore, he argued that the fractal dimension based on galaxy positions would not be enough and one should extend the concept of fractal to that of multifractal, with a spectrum of different exponents [18].
This argument was a good motivation for considering mass distributions other than Mandelbrot’s typical example, namely, a self-similar fractal set with mass uniformly distributed on it. This type of distribution, characterized by only one dimension, is usually called unifractal or monofractal, in contrast with the general case of multifractal distributions, characterized by many dimensions. However, one can have a sample of equivalent points from an arbitrary mass distribution, as is obvious by considering this distribution as a probability distribution. Actually, several researchers were soon making multifractal analyses of the distribution of galaxy positions, without considering their masses, which were not available [19, 20]. This approach is not mathematically wrong, but it does not reveal the properties of the real mass distribution. We show in Section 2.4 the result of a proper multifractal analysis of the real mass distribution.
Let us introduce a practical method of multifractal analysis that employs a coarse-grained mass distribution and is called “coarse multifractal analysis” [50]. A cube that covers the mass distribution is divided into a lattice of cells (boxes) of edge-length . Fractional statistical moments are defined aswhere the index runs over the set of nonempty cells, is the mass in the cell , is the total mass, and . The restriction to nonempty cells is superfluous for but crucial for . The power-law behavior of for is given by the exponent such that where is the box-size. The function is a global measure of scaling.
On the other hand, mass concentrates on each box with a different “strength” , such that the mass in the box is (if , then the smaller is, the larger is and the greater the strength is). In the limit of vanishing , the exponent becomes a local fractal dimension. Points with larger than the ambient space dimension are mass depletions and not mass concentrations. It turns out that the spectrum of local dimensions is related to the function by [50]. Besides, every set of points with common local dimension forms a fractal set with Hausdorff dimension given by the Legendre transform .
Of course, a coarse multifractal analysis only yields an approximated , which however has the same mathematical properties as the exact spectrum. The exact spectrum can be obtained by computing for a sequence of decreasing and establishing its convergence to a limit function.
For the application of coarse multifractal analysis to cosmology, it is convenient to assume that the total cube, of length , covers a homogeneity volume; that is to say, is the scale of transition to homogeneity. The integral moments , , are connected with the -point correlation functions, which are generally employed in cosmology [4] (we can understand the nonintegral moments with as an interpolation). A straightforward calculation, starting from (4), shows that , for , can be expressed in terms of the coarse-grained mass density aswhere is the average corresponding to coarse-graining length . In particular, is connected with the density two-point correlation function, introduced in Section 2.2. In fact, we can use (3) to integrate in a cell and calculate Assuming that , we obtain Therefore, , and, comparing with (5), . In a monofractal (or unifractal), this is the only dimension needed and actually , which implies and . In this case, we can write (5) asThis relation, restricted to integral moments, is known in cosmology as a hierarchical relation of moments [4, 7]. Actually, if all are equal for some in (4), then , for all , and therefore . Let us notice that, when , (6) implies that ; that is to say, we have a “monofractal” of dimension 3. In cosmology, this happens in the homogeneous regime, namely, for , but not for .
The multifractal spectrum provides a classification of mass concentrations, and is easy to compute through (4) and (5), followed by a Legendre transform. Its computation is equally simple for mass distributions in one or several dimensions. Actually, we know what type of multifractal spectrum to expect for self-similar multifractals: it spans an interval , is concave (from below), and fulfills [50]. Furthermore, the equality is reached at one point, which represents the set of singular mass concentrations that contains the bulk of the mass and is called the “mass concentrate.” The maximum value of is the dimension of the support of the mass distribution. We deduce that is tangent to the diagonal segment from to and to the horizontal half line that starts at and extends to the right. The values of at these two points of tangency are 1 and 0, respectively. Beyond the maximum of , we have . Notice that the restriction to and, in particular, to the calculation of only integral moments, misses information on the structure corresponding to a considerable range of mass concentrations (possibly, mass depletions). An example of this type of multifractal spectrum, with , is shown in Figure 3 (a full explanation of this figure is given below).
The maximum of , as the dimension of the support of the mass distribution, provides important information about voids [53]. If equals the dimension of the ambient space, then either voids are not totally empty or there is a sequence of totally empty voids with sizes that decrease so rapidly that the sequence does not constitute a fractal hierarchy of voids. This type of sequence is characterized by the Zipf law or Pareto distribution [3, 53, 54]. It seems that in cosmology and, furthermore, that there are actually no totally empty voids, as discussed below. In mathematical terms, one says that the mass distribution has full support, and we can also say that it is nonlacunar. We have already seen examples of nonlacunar distributions in Section 2.2, which are further explained below. Let us point out a general property of nonlacunar distributions. When equals the dimension of the ambient space, the corresponding is larger, because (it could be equal, for a nonfractal set of points, but we discard this possibility). The points with such values of are most abundant [ is maximum] and also are mass depletions that belong to nonempty voids. The full set of points in these voids is highly complex, so it is not easy to define individual voids with simple shapes [53] and may be better to speak of “the structure of voids” [52].
As simple examples of multifractal spectrum calculations, we consider the calculations for the one-dimensional mass distributions studied in Section 2.2. The case of inverse devil’s staircases is special, because the local dimension of zero-size mass condensations is . Furthermore, the set of locations of mass condensations is countable and hence has zero Hausdorff dimension; that is to say, (in accord with the general condition ). The complementary set, namely, the void region, has Hausdorff dimension and local dimension [deduced from the step length scaling (2)]. Since the set of mass concentrations fulfills , for , it is in fact the “mass concentrate.” The closure of this set is the full interval, which is the support of the distribution and is dominated by voids that are not totally empty, a characteristic property of nonlacunar distributions. This type of multifractal spectrum, with just one point on the diagonal segment from to and another on the horizontal half line from to the right, has been called bifractal [55]. However, the multifractal spectrum that is calculated through the function involves the Legendre transform and, therefore, is the convex hull of those two points, namely, the full segment from to . The extra points find their meaning in a coarse-grained approach to the bifractal.
Let us comment further on the notion of bifractal distribution. We first recall that a bifractal distribution of galaxies was proposed by Balian and Schaeffer [20], with the argument that there is one fractal dimension for clusters of galaxies and another for void regions. In fact, the bifractal proposed by Balian and Schaeffer [20] is more general than the one just studied, because they assume that the clusters, that is to say, the mass concentrations, do not have but a positive value. Such a bifractal can be considered as just a simplification of a standard multifractal, because any multifractal spectrum has a bifractal approximation. Indeed, any multifractal spectrum contains two distinguished points, as seen above, namely, the point where , which corresponds to the mass concentrate, and the point where is maximum. The former is more important regarding the mass and the second is more important regarding the supporting set of the mass distribution (the set of “positions”), whose dimension is . However, this supporting set is riddled by voids in nonlacunar multifractals.
The second example of multifractal spectrum calculation corresponds to the one-dimensional random mass distribution generated for Figure 2. This mass distribution is just a mathematical construction that is not connected with physics but is a typical example of continuous strictly singular mass distribution. The multifractal spectrum is calculable analytically, but we have just employed coarse multifractal analysis of the particular realization in Figure 2. This is a good exercise to assess the convergence of the method for a mass distribution with arbitrary resolution. In Figure 3 are displayed four coarse spectra, for , and , which show perfect convergence in the full range of , except at the ends, especially, at the large- end, corresponding to the emptiest voids. The limit is a standard spectrum of a self-similar multifractal, unlike the bifractal spectrum of the inverse devil’s staircase. Since the mass distribution of Figure 2 is nonlacunar, reaches its highest possible value, . The long range of values with suggests a rich structure of voids.
Finally, let us generalize the multifractal spectrum of the inverse devil’s staircase generated by the one-dimensional adhesion model to the three-dimensional case, for a later use. In three dimensions, in addition to point-like masses, there are filaments and sheets. Therefore, the (concave envelope of the) graph of is the union of the diagonal segment from to and the segment from to , with some for voids, as mass depletions (the conjecture that, along a given line, there is essentially a one-dimensional situation [46] suggests that ). There are three distinguished points of the multifractal spectrum in the diagonal segment, namely, , , and , corresponding to, respectively, nodes, filaments and sheets, where mass concentrates. Each of these entities corresponds to a Dirac-delta singularity in the mass density: nodes are simple three-dimensional Dirac-delta distributions, while filaments or sheets are Dirac-delta distributions on one- or two-dimensional topological manifolds. There is an infinite but countable number of singularities of each type. Moreover, the set of locations of each type of singularities is dense, which is equivalent to saying that every ball, however small, contains singularities of each type. The mathematical description is perhaps clumsy but the intuitive geometrical notion can be grasped just by observing Figure 1 (it actually plots a two-dimensional example, but the three-dimensional case is analogous).
2.4. Multifractal Analysis of the Large-Scale Structure
Of course, the resolutions of the available data of large-scale structure are not as good as in the above algorithmic examples, but we can employ as well the method of coarse multifractal analysis with shrinking coarse-graining length. As regards quality, we have to distinguish two types of data: data from galaxy surveys and data from cosmological -body simulations.
The first multifractal analyses of large-scale structure employed old galaxy catalogues. These catalogues contained galaxy positions but no galaxy masses. Nevertheless, Pietronero and collaborators [6, 8] managed to assign masses heuristically and carry out a proper multifractal analysis, with (4) and (5). The results were not very reliable, because the range of scales available was insufficient. This problem also affected contemporary multifractal analyses that only employed galaxy positions. Recent catalogues contain more faint galaxies and, therefore, longer scale ranges. Furthermore, it is now possible to make reasonable estimates of the masses of galaxies. In this new situation, we have employed the rich Sloan Digital Sky Survey and carried out a multifractal analysis of the galaxy distribution in the data release 7, taking into account the (stellar) masses of galaxies [56]. The results are reasonable and are discussed below. At any rate, galaxy survey data are still poor for a really thorough multifractal analysis.
Fortunately, the data from cosmological -body simulations are much better. On the one hand, the state-of-the-art simulations handle billions of particles, so they afford excellent mass resolution. Several recent simulations provide a relatively good scaling range, namely, more than two decades (in length factor), whereas it is hard to get even one decade in galaxy redshift surveys. Of course, the scale range is still small if compared to the range in our example of random mass distribution in one dimension, but it is sufficient, as will be shown momentarily. On the other hand, the bodies simulate the full dark matter dynamics, whereas galaxy surveys are restricted to stellar matter, which only gives the distribution of baryonic matter (at best). Of course, the distribution of baryonic matter can be studied in its own right. A comparison between multifractal spectra of large-scale structure in dark matter -body simulations and in galaxy surveys is made in [56], using recent data (the Bolshoi simulation, Figure 4, and the Sloan Digital Sky Survey). We show in Figure 5 the shape of the respective multifractal spectra.
Both the spectra of Figure 5 are computed with coarse multifractal analysis, so each plot displays several coarse spectra. The plot in the left-hand side corresponds to the Bolshoi simulation, which has very good mass resolution: it contains particles in a volume of . The coarse spectra in the left-hand side of Figure 5 correspond to coarse-graining length Mpc/ and seven subsequently halved scales. Other cosmological -body simulations yield a similar result and, in particular, one always finds convergence of several coarse spectra to a limit function. Since the coarse spectra converge to a common spectrum in different simulations, we can consider the found multifractal spectrum of dark matter reliable. Notice however that the convergence is less compelling for larger , corresponding to voids and, specifically, to emptier voids.
In the right-hand side of Figure 5 is the plot of four coarse spectra of a volume-limited sample (VLS1) from the 7th data release of the Sloan Digital Sky Survey, with 1765 galaxies in a volume of , for coarse-graining volumes (the coarse-graining length can be estimated as ). Four coarse spectra are the maximum number available for any sample that we have analyzed [56]. These four coarse spectra converge in the zone , of relatively strong mass concentrations, but there are fewer converging spectra for larger and the shape of the spectrum of mass depletions () is questionable. While it seems that the multifractal spectrum reaches its maximum value possible, , so that the distribution can be nonlacunar, the values for larger are quite uncertain [56]. In particular, it is unreliable that , a quite large value, to be compared with the considerably smaller and more reliable value of derived from -body simulations.
At any rate, we can assert a good concordance between the multifractal geometry of the cosmic structure in cosmological -body simulations and galaxy surveys, to the extent that the available data allow us to test it, that is to say, in the important part of the multifractal spectrum up to its maximum (the part such that ). The common features found in [56] and visible in Figure 5 are as follows: (i) a minimum singularity strength ; (ii) a “supercluster set” of dimension where the mass concentrates; and (iii) a nonlacunar structure (without totally empty voids). It is to be remarked that corresponds to the edge of diverging gravitational potential.
As regards the visual morphological features, a web structure can be observed in both the -body and galaxy distributions. In particular, the slice of the Bolshoi simulation displayed in Figure 4 definitely shows this structure and, furthermore, it looks like the two-dimensional distribution from the adhesion model displayed in the right-hand side of Figure 1. However, when we compare the real multifractal spectrum of the large-scale structure with the multifractal spectrum of the three-dimensional mass distribution formed according to the adhesion model, as seen at the end of Section 2.3, we can appreciate considerable differences. Therefore, in spite of the fact that the real cosmic web looks like the web structure generated by the adhesion model, the differences in their respective multifractal spectra reveal differences in the distributions. Indeed, the cosmic web multifractal spectrum, well approximated by the left-hand side graph in Figure 5, contains no point-like or filamentary Dirac-delta distributions and it does not seem to contain sheet distributions, although the graph passes close to the point . The differences can actually be perceived by carefully comparing the right-hand side of Figure 1 with Figure 4. Notice that the features in the adhesion model cosmic web appear to be sharper, because they correspond to zero-size singularities, namely, nodes having zero area and filaments zero transversal length. In Figure 4, point-like or filamentary mass concentrations are not so sharp, because they correspond to power-law concentrations.
In fact, the notion of gravitational sticking in the adhesion model is a rough approximation that can only be valid for sheets, because the accretion of matter to a point () or a sharp filament () involves the dissipation of an infinite amount of energy (in Newtonian gravity). It can be proved that any mass concentration with generates a diverging gravitational potential. Therefore, it is natural that the real cosmic web multifractal spectrum starts at (Figure 5). Moreover, Figure 5 shows that ; that is to say, the set of locations of singularities with has dimension zero, which suggests that the gravitational potential is finite everywhere.
The above considerations are based on the assumption of continuous matter, and even the adhesion model initially assumes a fluid model for the dynamics, even though it gives rise to point-like mass condensations (at nodes). However, -body dynamics is intrinsically discrete and this can modify the type of mass distribution on small scales. In fact, cosmological -body simulations have popularized halo models, in which the matter distribution on small scales is constructed on a different basis.
2.5. Halo Models
The idea of galactic dark matter halos has its roots in the study of the dynamics in the outskirts of galaxies, which actually motivated the introduction of the concept of dark matter. These dark halos are to be related to the visible shapes of galaxies but should be, in general, more spherical and extend to large radii, partially filling the intergalactic space and even including several galaxies. Therefore, the large-scale structure of dark matter should be the combination of quasi-spherical distributions in halos with the distribution of halo centers. This picture is analogous to a statistical model of the galaxy distribution introduced by Neyman and Scott in 1952 [57], in which the full distribution of galaxies is defined in terms of clusters of points around separate centers combined with a distribution of cluster centers. Neyman and Scott considered a “quasi-uniform” distribution of cluster centers and referred to the Poisson distribution [57]. Of course, modern halo models consider the clustering of halo centers [24].
The halo model is apparently very different from the models that we have seen above, but it is actually connected with the basic fractal model and its cluster hierarchy. Indeed, we can consider Neyman and Scott’s original model as consisting of a hierarchy with just two levels, namely, point-like galaxies and their clusters, with little correlation between cluster centers. A specific formulation of such model is Fry’s hierarchical Poisson model, which he employed to study the statistics of voids [58]. This two-level hierarchy can be generalized by adding more levels (as already discussed [53]). An infinite cluster hierarchy satisfying Mandelbrot’s conditional cosmological principle constitutes the basic fractal model. Of course, halo models are normally defined so that they have patterns of clustering of halos that are independent of the distributions inside halos, departing from self-similarity.
In fact, the distribution of halo centers is often assumed to be weakly clustered and is sometimes treated perturbatively. However, it tends to be recognized that halos must preferentially lie in sheets or filaments and actually be ellipsoidal [25, 31]. Neyman and Scott required the distribution of points around the center to be spherical, without further specification. Modern halo models study in detail the halo density profiles, in addition to deviations from spherical symmetry. For theoretical reasons and also relying on the results of cosmological -body simulation, density profiles are assumed to have power-law singularities at their centers. There are a number of analytic forms of singular profiles, with parameters that can be fixed by theory or fits [24]. A detailed account of this subject is beyond the scope of this work. In our context, the power-law singularity and its exponent, which normally has values in the range , are especially important.
We can compare the density singularities in halos with the singular mass concentrations already studied in multifractal distributions. A coarse-grained multifractal distribution indeed consists of a collection of discrete mass concentrations and can be formulated as a fractal distribution of halos [27, 28]. Models based on this idea constitute a synthesis of fractal and halo models. Furthermore, the idea naturally agrees with the current tendency to place halos along filaments or sheets and hence assume strong correlations between halo centers. The range of the density power-law exponent is equivalent to the range of the mass concentration dimension, defined by . In view of Figure 5, this range is below the value for the mass concentration set [fulfilling ], which is . This difference can be attributed to a bias towards quasi-spherical rather than ellipsoidal halos, because the former are more concentrated. In the web structure generated by the adhesion model, nodes are more concentrated than filaments or sheets.
It is to be remarked that the theories of halo density profiles are consistent with a broad range of power-law exponents and even with no singularity at all (Gurevich-Zybin’s theory [41] is an exception). In fact, the quoted density-singularity range is mainly obtained from the analysis of -body simulations. For this and other reasons, it is necessary to judge the reliability of -body simulations in the range of scales where halos appear [32]. Notice that the reliability of -body simulations on small scales has been questioned for a long time [59–61]. The problem is that two-body scattering spoils the collisionless dynamics, altering the clustering properties. Although the mass distribution inside an -body halo is smooth (except at the center), this can be a consequence of discreteness effects: an -body halo experiences a transition from a smooth distribution to a very anisotropic and nonsmooth web structure over a scale that should vanish in the continuum limit [31].
To conclude this section, let us try to present an overview of the prospects of the halo model. The original idea of the halo model, namely, the description of the cosmic structure in terms of smooth halos with centers distributed in a simple manner and, preferably, almost uniformly, seems somewhat naive, inasmuch as it is an ad hoc combination of two simple types of distributions, which can only be valid on either very small or very large scales. Moreover, smooth halos seem to be the result of -body simulations in an unreliable range of scales. In spite of these problems, the notion of galactic dark matter halos, that is to say, of small-scale dark matter distributions that control the dynamics of baryonic matter and the formation of galaxies, is surely productive. However, the cold dark matter model has problems of its own, because the collisionless dynamics tends, anyway, to generate too much structure on small scales, and it is expected that the baryonic physics will change the dynamics on those scales [62]. The halo model provides a useful framework to study this question [26]. The ongoing efforts to unify the halo and multifractal models [27–31] will hopefully help to give rise to a better model of the large-scale structure of the universe.
3. Formation of the Cosmic Web
The halo or fractal models are of descriptive nature, although they can be employed as frameworks for theories of structure formation. In contrast, the Zeldovich approximation and the adhesion model constitute a theory of structure formation. Unfortunately, this theory is quite simplistic, insofar as the nonlinear gravitational dynamics is very simplified, to the extent of being apparently absent. At any rate, in the correct equations, the equation for the gravitational field is linear and the nonlinearity is of hydrodynamic type (see the equations in Section 3.2). This nonlinearity gives rise to turbulence and, in the adhesion model in particular, to Burgers turbulence, which is highly nontrivial.
So we place our focus on the properties of Burgers turbulence. However, we begin by surveying other models of the formation of the cosmic web that are not directly based on approximations of the cosmological equations of motion (Section 3.1). Next, we formulate these equations in Section 3.2 and proceed to the study of turbulence in Burgers dynamics, first by itself (Section 3.3) and second in the setting of stochastic dynamics (Section 3.4). Finally, we consider the full nonlinear gravitational dynamics (Section 3.5).
3.1. Models of Formation of Cosmic Voids
Let us recall two models of the formation of the cosmic web, namely, the Voronoi foam model [33, 34] and Sheth and van de Weygaert’s model of the hierarchical evolution of voids [63]. Both these models focus on the formation and evolution of cosmic voids, but the Voronoi foam model is of geometric nature whereas the hierarchical void model is of statistical nature. They are both heuristic and are not derived from the dynamical equations of structure formation, but the Voronoi foam model can be naturally connected with the adhesion model.
The Voronoi foam model focuses on the underdense zones of the primordial distribution of small density perturbations. The minima of the density field are the peaks of the primordial gravitational potential and, therefore, the seeds of voids. Indeed, voids form as matter flows away from those points. The basic model assumes that the centers of void expansion are random; namely, they form a Poisson point field, and it also assumes that every center generates a Voronoi cell. The formation of these cells can be explained by a simple kinematic model that prescribes that matter expands at the same velocity from adjacent centers and, therefore, a “pancake” forms at the middle; that is to say, a condensation forms in a region of the plane perpendicular to the joining segment through its midpoint [34, Fig. 40]. The Voronoi foam results from the simultaneous condensation on all these regions.
Let us describe in more detail the formation of the cosmic web according to this kinematic model [33, 34]. After the matter condenses on a wall between two Voronoi cells, it continues to move within it and away from the two nuclei, until it encounters, in a given direction, the matter moving along the other two walls that belong to the cell of a third nuclei. The intersection of the three walls is a Voronoi edge, on which the matter from the three walls condenses and continues its motion away from the three nuclei. Of course, the motion along the edge eventually leads to an encounter with the matter moving along three other edges that belong to the cell of a fourth nuclei, and it condenses on the corresponding node. All this process indeed follows the rules of the adhesion model, but the Zeldovich approximation of the dynamics is replaced by a simpler kinematic model, which produces a particular cosmic web, namely, a Voronoi foam. It is easy to simulate matter distributions that form in this way, and those distributions are suitable for a correlation analysis [34]. This analysis reveals that the two-point correlation function is a power-law, with an exponent that grows with time, as matter condenses on lower dimensional elements of the Voronoi foam. Naturally, is close to two in some epoch (or close to the standard value ).
The only parameter of the Voronoi foam model with random Voronoi nuclei is the number density of these nuclei. In fact, this construction is well known in stochastic geometry and is known as the Poisson-Voronoi tessellation. The number density of nuclei determines the average size of Voronoi cells, that is to say, the average size of cosmic voids. The exact form of the distribution function of Voronoi cell sizes is not known, but it is well approximated by a gamma distribution that is peaked about its average [64]. Therefore, the Voronoi foam model reproduces early ideas about cosmic voids, namely, the existence of a “cellular structure” with a limited range of cell sizes [37]. This breaks the self-similarity of the cosmic web found with the adhesion model by Vergassola et al. [46] (see Section 2.1). The self-similarity of cosmic voids, in particular, is demonstrated by the analysis of galaxy samples or cosmological -body simulations [53]. Notice that gamma distributions of sizes are typical of voids of various shapes in a Poisson point field [53]. In a Voronoi foam with random nuclei, the matter points would be the nodes of condensation, namely, the Voronoi vertices, and these do not constitute a Poisson point field but are somehow clustered. At any rate, a distribution of sizes of cosmic voids that is peaked about its average is not realistic.
Moreover, a realistic foam model must account for the clustering of minima of the primordial density field and for the different rates of void expansion, giving rise to other types of tessellations. In fact, such models seem to converge towards the adhesion model [33]. The adhesion model actually produces a Voronoi-like tessellation, characterized by the properties of the initial velocity field (or the initial density field) [47]. In cosmology, the primordial density field is not smooth and, hence, has no isolated minima. Although a smoothened or coarse-grained field can be suitable for studying the expansion of cosmic voids, the distribution of voids sizes is very different from the one predicted by the basic Voronoi foam model [53].
Let us now turn to Sheth and van de Weygaert’s hierarchical model of evolution of voids [63]. It is based on an ingenious reversal of the Press-Schechter and excursion set approaches to hierarchical gravitational clustering. In fact, the idea of the Press-Schechter approach is better suited for the evolution of voids than for the evolution of mass concentrations (a comment on the original approach can be found at the end of Section 3.5). The reason is the “bubble theorem”, which shows that aspherical underdense regions tend to become spherical in time, unlike overdense regions [33, 34].
Sheth and van de Weygaert’s hierarchical model of voids is more elaborate than the Voronoi foam model but it also predicts a characteristic void size. To be precise, it predicts that the size distribution of voids, at any given time, is peaked about a characteristic void size, although the evolution of the distribution is self-similar. Actually, the lower limit to the range of void sizes is due to a small-scale cut-off, which is put to prevent the presence of voids inside overdense regions that are supposed to collapse (the “void-in-cloud process”) [63]. The resulting distribution of void sizes looks like the gamma distribution of the Voronoi foam model, and Sheth and van de Weygaert [63] indeed argue that the hierarchical model, which is more elaborate, somehow justifies the Voronoi foam model, which is more heuristic.
If the small-scale cut-off is ignored, by admitting that voids can form inside overdense regions, then the distribution of void sizes in Sheth and van de Weygaert’s model has no characteristic void size and, for small sizes, it is proportional to the power of the void volume [63]. Power-law distributions of void sizes characterize lacunar fractal distributions, which follow the Zipf law for voids, but the exponent is not universal and depends on the fractal dimension [53]. At any rate, the results of analyses of cosmic voids and the cosmic web itself, presented in [53] and Section 2.4, show that the matter distribution is surely multifractal and nonlacunar. This implies that cosmic voids are not totally empty but have structure inside, so that the self-similar structure does not manifest itself in the distributions of void sizes but in the combination of matter concentrations and voids (as discussed in [53]). In fact, the cosmic web generated by the adhesion model with scale-invariant initial conditions is a good example of this type of self-similarity (Section 2.1).
Therefore, the conclusion of our brief study of models of formation of cosmic voids is that just the adhesion model is more adequate than these other models for describing the formation of the real cosmic web structure. The adhesion model derives from the Zeldovich approximation to the cosmological dynamics, which we now review.
3.2. Cosmological Dynamics
We now recall the cosmological equations of motion in the Newtonian limit, applicable on scales small compared to the Hubble length and away from strong gravitational fields [4, 5, 33]. The Newtonian equation of motion of a test particle in an expanding background is best expressed in terms of the comoving coordinate and the peculiar velocity , where is the Hubble constant (for simplicity, we can take in the present time). So the Newton equation of motion, , with total gravitational field , can be written aswhere the net gravitational field, , is defined by subtracting the background field . The equation for this background field is that is to say, the dynamical FLRW equation for pressureless matter (with a cosmological constant). This equation gives and hence . To have a closed system of equations, we must add to (9) the equation that gives in terms of the density fluctuations,and the continuity equation,These three equations form a closed nonlinear system of equations.
The equations can be linearized when and where the density fluctuations are small. Within this approximation, one actually needs only (9), which is nonlinear because of the convective term in but becomes linear when the density fluctuations and hence and are small. The linearized motion iswhere is the growth rate of linear density fluctuations. Redefining time as , we have a simple motion along straight lines, with constant velocities given by the initial net gravitational field.
3.3. Burgers Turbulence
The Zeldovich approximation prolongs the linear motion (13) into the nonlinear regime, where is not small. To prevent multistreaming, the adhesion model supplements the linear motion with a viscosity term, as explained in Section 2.1, resulting in the equationwhere is the peculiar velocity in -time. To this equation, it must be added the condition of no vorticity (potential flow), , implied by . It is to be remarked that the viscosity term is in no way fundamental and the same job is done by any dissipation term, in particular, by a higher-order hyperviscosity term [65]. Eq. (14) is the (three-dimensional) Burgers equation for pressureless fluids (tildes are suppressed henceforth). Any prevents multistreaming and the limit is the high Reynolds-number limit, which gives rise to Burgers turbulence. Whereas incompressible turbulence is associated with the development of vorticity, Burgers turbulence is associated with the development of shock fronts, namely, discontinuities of the velocity. These discontinuities arise at caustics and give rise to matter accumulation by inelastic collision of fluid elements. The viscosity measures the thickness of shock fronts, which become true singularities in the limit .
A fundamental property of the Burgers equation is that it is integrable; namely, it becomes a linear equation by applying to it the Hopf-Cole transformation [1, 35, 36]. This property shows that the Burgers equation is a very special case of the Navier-Stokes equation. In spite of it, Burgers turbulence is a useful model of turbulence in the irrotational motion of very compressible fluids. At any rate, the existence of an explicit solution of the Burgers equation makes the development of Burgers turbulence totally dependent on the properties of the initial velocity distribution, in an explicit form. The simplest solutions to study describe the formation of isolated shocks in an initially smooth velocity field. Fully turbulent solutions are provided by a self-similar initial velocity field, such that , as advanced in Section 2.1. Indeed, (14) is scale invariant in the limit ; namely, it is invariant under simultaneous space and time scalings, , , for arbitrary , if the velocity is scaled as . Therefore (14) has self-similar solutions, such thatA self-similar velocity field that is fractional Brownian on scales larger than the “coalescence length” contains a distribution of shocks that is self-similar on scales below . Such velocity field constitutes a state of decaying turbulence, since kinetic energy is continuously dissipated in shocks [1, 35, 36]. In the full space , as is adequate for the cosmological setting, the transfer of energy from large scales, where the initial conditions hold, to the nonlinear small scales, proceeds indefinitely.
The simplicity of this type of self-similar turbulence allows one to analytically prove some properties and make reasonable conjectures about others [35, 36, 46]. A case that is especially suitable for a full analytical treatment is the one-dimensional dynamics with , namely, with Brownian initial velocity. As regards structure formation, the dynamics consists in the merging of smaller structures to form larger ones, in such a way that the size of the largest structures is determined by the “coalescence length” and, therefore, they have masses that grow as a power of [1, 46]. The merging process agrees with the “bottom-up” picture for the growth of cosmological structure [5]. This picture is appropriate for the standard cold dark matter cosmology.
From the point of view of the theory of turbulence, centered on the properties of the velocity field, the formed structure consists of a dense distribution of shock fronts with very variable magnitude. Therefore, the kinetic energy dissipation takes place in an extremely nonuniform manner, giving rise to intermittency [1]. A nice introduction to the phenomenon of intermittency as a further development of Kolmogorov’s theory of turbulence can be found in Frisch’s book [66, ch. 8]. Kolmogorov’s theory assumes that the rate of specific dissipation of energy is constant. Intermittency is analyzed by modeling how kinetic energy is distributed among structures, either vortices in incompressible turbulence or shocks in irrotational compressible turbulence. Mandelbrot proposed that intermittency is the consequence of dissipation in a fractal set rather than uniformly [3, §10]. Mandelbrot’s idea is embodied in the -model of the energy cascade [66]. However, this model is monofractal and a bifractal model is more adequate, especially for Burgers turbulence. Although this bifractality refers to the velocity field [66], it corresponds, in fact, to the bifractality of the mass distribution in the one-dimensional adhesion model already seen in Section 2.3. Let us briefly study the geometry of intermittency in Burgers turbulence.
Kolmogorov’s theory of turbulence builds on the Richardson cascade model and is based on universality assumptions, valid in the limit of infinite Reynolds number () and away from flow boundaries. As explained by Frisch [66], it is possible to deduce, by assuming homogeneity, isotropy, and finiteness of energy dissipation in the limit of infinite Reynolds number, the following scaling law for the moments of longitudinal velocity increments:where is the specific dissipation rate, and . To introduce the effect of intermittency, the Kolmogorov scaling law can be generalized towhere and does not depend on (and does not have to be related to ). Eq. (18) expresses that the fractional statistical moments of have a power-law behavior given by the exponent . Therefore, the equation is analogous to the combination of (4) and (5) for multifractal behavior. The functions and are indeed connected in Burgers turbulence (see below).
If in (18), as occurs in (16), then the probability is Gaussian. The initial velocity field that we have chosen for Burgers turbulence is Gaussian, namely, fractional Brownian, with (Section 2.1). Therefore, if (which would fit (16) for ). In general, intermittency manifests itself as concavity of [66]. In the Burgers turbulence generated by the given initial conditions, a type of especially strong intermittency develops for any and . This intermittency is such that , in addition to being concave, has a maximum value equal to one. This can be proved by generalizing the known argument for isolated shocks [35] to a dense distribution of shocks, as is now explained.
Let us replace in (18) the ensemble average with a spatial average (notice that fractional statistical moments of in (4) are defined in terms of a spatial average). The average in (18), understood as a spatial average, can be split into regular and singular points (with shocks). At regular points, the inverse Lagrangian map is well defined and, from (13), Self-similarity in time allows one to set without loss of generality. Hence, we deduce the form of the spatial velocity increment over at a regular point ,Since the properties of the Lagrangian map are known, one can deduce properties of the velocity increment. For example, at a regular point , the term is subleading for , because the derivatives of the inverse Lagrangian map vanish.
To proceed, we restrict ourselves to one dimension, in which the analysis is simpler and we have, for , where the first term is due to the regular points and the second term is due to the set of -discontinuities (shocks) such that . Clearly, putting back, corresponds to the self-similar expansion of voids, whereas and corresponds to shocks. If , the first term dominates as , and vice versa. Therefore, , for , and , for . Notice that the number of terms in the sum increases as , but the series is convergent for , due to elementary properties of devil’s staircases: if the spatial average is calculated over the interval that corresponds to an initial interval , then and when (actually, the series converges for , where is the Besicovitch-Taylor exponent of the sequence of gaps, which is in this case equal to , because of the gap law [3, p. 359]).
The form of , with its maximum value of one, implies that the probability is very non-Gaussian for , as is well studied [35, 36] and we now show. In one dimension, the ratio is a sort of “coarse-grained derivative”, related by (20) to the coarse-grained density (the length of the initial interval that is mapped to the length is the coarse-grained density normalized to the initial uniform density). Hence, we derive a simple relation between the density and the derivative of the velocity, namely, . We can also obtain the probability of velocity gradients from , in the limit . However, this limit is singular and, indeed, the distribution of velocity gradients can be called strictly singular, like the distribution of densities, because, at regular points, where , then , and at shock positions, where , then . Naturally, it is the shocks that make very non-Gaussian and forbid that grows beyond . The intermittency is so strong that the form of is very different from Kolmogorov’s Eq. (16) and the average dissipation rate seems to play no role.
As the mass distribution can be expressed in terms of the velocity field, the intermittency of the latter is equivalent to the multifractality of the former. Indeed, the bifractal nature of , in one dimension, is equivalent to the dual form of for or . Let us notice, however, that is universal, whereas the bifractal studied in Section 2.3 depends on the value of (at the point for voids), and so does . Actually, the full form of also depends on .
For the application to cosmology, we must generalize these results to three dimensions. Although this generalization is conceptually straightforward, it brings in some complications: in addition to point-like mass condensations, there appear another two types of mass condensations, namely, filaments and sheets, the latter actually corresponding to primary shocks (furthermore, there is mass flow inside the higher-dimensional singularities, and this flow must be determined [47]). The form of is still determined by the primary shocks and is, therefore, unaltered [35, 36]. Nevertheless, the mass distribution itself is no longer bifractal, as is obvious. The form of has already been seen in Section 2.3; namely, it is the union of the diagonal segment from to and the segment from to , with .
In summary, the adhesion model, with initial velocity and density fields having Gaussian power-law fluctuations, namely, being fractional Brownian fields, leads to a self-similar and strongly intermittent Burgers turbulence, such that various quantities of interest can be calculated analytically. This dynamics gives rise to a self-similar cosmic web that is quite realistic. Therefore, it constitutes an appealing model of structure formation in cosmology [1, 46]. However, it has obvious shortcomings: (i) in cosmology, the initial but “processed” power spectrum is not a power-law [5, 33]; (ii) the Burgers equation is integrable and therefore the adhesion model may be too simplified model of the gravitational dynamics, which actually is chaotic; and (iii) the condensation of matter in shock fronts is considered as an inelastic collision but the dissipated kinetic energy is assumed to disappear without further effect. In general, chaos and dissipation are connected: chaotic dynamics erases memory of initial conditions, giving rise to an irreversible process in which entropy grows. To deal with these problems, we need a more general framework. A first step is taken by connecting the Burgers equation with a well-studied stochastic equation, namely, the Kardar–Parisi–Zhang (KPZ) equation, as we study next.
3.4. The Stochastic Burgers Equation
After studying freely decaying Burgers turbulence in the adhesion model, we proceed to the study of forced Burgers turbulence [35, 36]. It is understood that the forcing in cosmology is not external but is just the feedback of dissipated energy. Indeed, we intend to achieve a more complete understanding of the dynamics that is based on a general relation between dissipation and fluctuating forces, such as the relation embodied in the classical fluctuation-dissipation theorem [67]. So we rely on the theory of stochastic dynamics, with origin in the Langevin equation for Brownian motion. The Langevin equation is derived by simply adding a fluctuating force to Newton’s second law. This force is due to the interaction with the environment and must contain a slowly varying dissipative part and a rapidly varying part with null average (a “noise”).
The Langevin equation has a limited scope. A qualitative development in stochastic dynamics consisted in considering partial differential equations, with fluctuating forces that depend on space and time. A simple example is the Edwards-Wilkinson equation, which results from adding noise to the diffusion equation [68]. Assuming that the scalar field of the Edwards-Wilkinson equation is the “height” of some substrate, the equation can be employed in the study of stochastic growth. The Edwards-Wilkinson equation is linear but a more suitable equation for stochastic growth, namely, the Kardar–Parisi–Zhang (KPZ) equation, was obtained by adding to it a nonlinear term. For several reasons, this equation has become a paradigm in the theory of stochastic nonlinear partial differential equations. In fact, by taking the scalar field in the KPZ equation as a velocity potential, the equation becomes equivalent to the stochastic Burgers equation, in which a noise term, that is to say, a random force, is added to the right-hand side of (14), so that it contains dissipation plus noise.
The presence of noise in the stochastic Burgers equation alters the evolution of the velocity and density fields. For example, mass condensations can gain energy and hence fragment, so adhesion is no longer irreversible. It may even not happen at all. The different kind of evolution induced by noise is best appreciated if the spatial region is bounded, with size (typically, the use of periodic boundary conditions makes space toroidal). Then, the total mass and kinetic energy are bounded and the free decay of Burgers turbulence eventually leads to a simple state. For example, in one dimension, all the kinetic energy is dissipated and all the mass condensates on a single point (for generic initial conditions) [35, 36]. As noise injects energy, the dissipative dynamics cannot relax to the state of minimum energy but to a fluctuating state, which may or may not have mass condensations, according to the characteristics of the noise.
In the Langevin dynamics, the particle motion relaxes until it is in thermal equilibrium with the medium and has the Maxwell velocity distribution [67]. The strength of the noise is given by the temperature , and the velocity fluctuations fulfillwhere is the mass of the particle and is the relaxation time. This equation says that the variance of is independent of the absolute time , and this holds for all , provided that . In contrast, the value of explicitly depends on and on the initial condition . The proportionality of the variance of to is characteristic of Brownian motion.
The thermal noise in the Langevin equation is Gaussian and macroscopically uncorrelated [67]; that is to say, it has white frequency spectrum. However, there is no reason to assume thermal equilibrium, in general. Indeed, various “colors” are considered in stochastic dynamics. In the theory of stochastic partial differential equations the noise is Gaussian and has a general power spectrum . In particular, a noise with power-law spectrum of spatial correlations but still uncorrelated in time is favored, namely, with (see below).
The Edwards-Wilkinson equation for a scalar field is equivalent to the stochastic Burgers equation for with the convective term suppressed. The solution of the Edwards-Wilkinson equation is more difficult than the solution of the Langevin equation, but it is a linear equation and hence soluble, by decoupling of the Fourier spatial modes. The dynamics also consists in a relaxation to a stationary state [68]. However, the relaxation takes place at a rate that depends on the “roughness” of at . Indeed, Fourier spatial modes have relaxation times inversely proportional to , so relaxation is faster on the smaller scales.
The KPZ dynamics also begins to relax on the small scales. On these scales, the nonlinear term, which corresponds to the convective term in the Burgers equation, becomes important. As said above, the combination of nonlinearity and dissipation tends to build structure whereas the noise tends to destroy it. In fact, the KPZ dynamics is ruled by the competition between nonlinearity and noise. To proceed, let us introduce the dynamic scaling hypothesis [68], as a generalization of (22).
Given a general stochastic partial differential equation for a scalar field , it is reasonable to assume thatwhere and are critical exponents and is a scaling function. This function expresses the crossover between the power laws of purely spatial and purely temporal correlations: approaches a constant for and behaves as for . Eq. (23) generalizes (22) to both space and time and, hence, is analogous to the space-time similarity already considered in (15) [the redefinition makes the connection with (15) or (22) more apparent]. (Let us notice that the averages considered in this section are taken with respect to realizations of the random noise, unlike in the preceding section, in which they are taken with respect to realizations of the random initial conditions.) Eq. (15) implies that the initial conditions are modified below a scale that grows as a power of time [with exponent ] and, likewise, (23) implies that relaxation to the noise-induced asymptotic distribution takes place below a scale that grows as a power of , with exponent .
In the case of the Edwards-Wilkinson equation in dimensions for a field , dynamic scaling is borne out by the explicit solution. In the solution, the Fourier mode is a linear combination of the noise Fourier modes, so a Gaussian noise produces a Gaussian [68]. This field fulfills (23), so and , for fixed , are actually fractional Brownian fields with Hurst exponents and , respectively. With noise spectrum , the explicit solution yields and , as can also be deduced from the scaling properties of the equation. Therefore, and , at fixed , become less rough as grows. For example, in one dimension and with white noise, , so describes, for fixed , a Brownian curve, whereas is uncorrelated (white noise like). Moreover, the dynamics consists in a relaxation to thermal equilibrium, in which fluid elements acquire the Maxwell velocity distribution, in accord with the fluctuation-dissipation theorem.
The KPZ equation is nonlinear and its solution is a non-Gaussian field. Therefore, (23) is just an assumption, supported by theoretical arguments and computer simulations. The one-dimensional case, with white noise, is special, insofar as the effect of the nonlinear term is limited, the fluctuation-dissipation theorem holds, and , like in the Edwards-Wilkinson equation. For general dimension and noise uncorrelated in time, the Burgers equation has Galilean invariance (invariance under addition of a constant vector to plus the consequent change of ), which enforces the relation [68]. In one dimension and with white noise, this relation implies that , which shows that the KPZ and Edwards-Wilkinson equations differ in temporal correlations; that is to say, the convective term of the stochastic Burgers equation does have an effect. Nevertheless, the dynamics also consists in a relaxation to the Maxwell velocity distribution. In the state of thermal equilibrium, there is no structure.
At any rate, for cosmology, we must consider the three-dimensional KPZ equation. The study of the KPZ equation with “colored” noise has been carried out by Medina et al. [69], employing renormalized perturbation theory, namely, perturbation of the linear equation by the nonlinear term and subsequent renormalization. Within this approach, the types of noise that give rise to dynamical scaling are found as fixed points of the dynamical renormalization group. For small , the noise power spectrum assumes the universal form that is to say, the noise consists of a white noise plus a power-law long-range noise [69]. Unfortunately, in three dimensions, renormalized perturbation theory is of little consequence: the only stable renormalization group fixed point is the trivial one (vanishing nonlinear convective term), and it is stable only if [69]. This means that a noise with sufficient power on large scales inevitably leads to a strong-coupling stationary state and questions the application of perturbation theory to the stochastic Burgers equation with spatially correlated noise in three dimensions. In spite of it, the stochastic adhesion model has been long studied in cosmology and mostly with perturbative treatments, which are still being applied [70]. Certainly, nonperturbative methods seem more sensible [71].
In fact, renormalized perturbation theory, in any dimension , is valid only for [69]. For larger , the long-range noise correlations induce analogous long-range correlations in and and give rise to intermittency of (as studied for in [72]). In particular, the field becomes strongly non-Gaussian. This is a purely nonlinear effect and may seem counterintuitive, because, in the Edwards-Wilkinson equation, and , at fixed , become less rough as grows. The progressive smoothness of turns it into a real hydrodynamical field (in contrast with the white noise like velocity in with , for example). Although the nonlinearity, namely, the convective term of the Burgers equation, does not spoil this smoothing, generally speaking, it gives rise to shocks and, hence, intermittency. The combination of smoothing and formation of shocks can be seen, for example, in simulations in [72, Fig. 4]. Therefore, for , we have the regime in which structures can form in spite of the noise, that is to say, the properly turbulent regime. As we focus on the field in this regime, we replace (23) with (18) [we have ].
The turbulent regime of forced Burgers turbulence has been well studied, with a combination of computer simulations and theoretical arguments related to the Kolmogorov theory. A simple argument of Chekhlov and Yakhot [73], for , shows that corresponds to an “almost constant” (logarithmic) energy flux in Fourier space (a balanced Richardson energy cascade). Furthermore, their numerical simulations show that develops shocks (with the corresponding mass condensations) which give rise to power-law correlations with exponents and . The latter corresponds to the Kolmogorov scaling, in accord with (16). Let us also notice that the noise exponent is such that the noise strength has the dimensions of the dissipation rate . But Chekhlov and Yakhot’s argument is incomplete, so that larger values of are suitable [72]. A breakthrough in the study of intermittency and Kolmogorov scaling in forced Burgers turbulence (in one dimension) has been the use by Polyakov [74] of nonperturbative methods borrowed from quantum field theory. Polyakov [74] assumes that the Burgers noise correlation function (in ordinary space) is twice differentiable; that is to say, . Boldyrev [75] has extended Polyakov’s methods to the range .
For cosmology, we must consider . Kolmogorov’s universality and its consequences for cosmology are considered as a guiding principle in [71]. The conclusion is that the interesting range of the exponent for the cosmic structure, in , is . The lower limit is such that the noise strength has the dimensions of and the Burgers noise correlation function is proportional to . Moreover, diverges for if . The reason for the upper limit is that the -dependent part of the noise correlation function depends explicitly on the system size for . The range corresponds to large-scale forcing that is such that the dissipation rate only depends on what happens on large scales and such that the -dependent part of the noise correlation function, proportional to , is universal.
Furthermore, there is intermittency for , in particular, in the range . To measure it, we can use (18), where now the exponent depends on the noise exponent (and the average is with respect to the random noise). The intermittency increases with increasing and progressively affects exponents with lower values of [76]. The analysis of the range obtains strong intermittency, namely, for , but . In particular, while in the range [76]. At , the exponent reaches its maximum , and it remains fixed for . The minimum value , for , corresponds to the Kolmogorov scaling (16).
Especially interesting for cosmology is that the density two-point correlation function can be obtained in terms of correlation functions of the velocity field; in fact, it can be obtained in terms of just the two-point correlation function [71]. Furthermore, the density two-point correlation function is a power-law, taking the form in (3), with an exponent determined by , namely, . Given the range of , we deduce that which is a reasonable range (although the preferred value of is somewhat larger).
To summarize, the stochastic Burgers equation with spatially correlated noise constitutes an interesting model of structure formation, in which the cosmic web is not determined by ad hoc initial conditions but is the result of the interplay of the inelastic gravitational condensation of matter with the consequent dynamical fluctuations, so that energy can be conserved on the average. The relaxation to the asymptotically stable state takes place over a scale that grows with time, namely, , where is determined by the dynamics instead of the initial conditions. The asymptotically stable state has simple scaling properties and, furthermore, it achieves Kolmogorov’s universality. However, this universality does not lead to uniqueness, in the sense of total independence of the large scales: indeed, there is an allowed range for the noise exponent . At any rate, the model predicts a stable fractal cosmic web, such that the exponent of the density two-point correlation function is within a reasonable range. Furthermore, the formation of structures that are stable in the highly nonlinear regime connects, in spirit, with famous Peebles’ stable clustering hypothesis [4], to be studied in the next section.
The rationale for a fluctuating force is not restricted to the cosmological equations in the Zeldovich approximation. One can employ the original equations, namely, (9), (11), and (12). This system of equations can be supplemented with viscous and random forcing terms [77]. At any rate, it behooves us to take a more general standpoint and study the full dynamics of gravitational clustering.
3.5. Nonlinearity, Chaos, and Turbulence in Gravitational Clustering
The first comprehensive attempts to treat analytically the theory of large-scale structure formation are due to Peebles and collaborators [4]. They employed general principles of statistical mechanics to formulate the problem and added a scaling hypothesis to obtain definite solutions. This interesting work did not lead to a geometrical picture of structure formation until the connection with Mandelbrot’s ideas [3] and the development of the Zeldovich approximation and the adhesion model [1]. Let us recall the basics of Peebles’ approach (a complete introduction to it, with some novelties, is provided by Shani and Coles [33]).
Surely, the most general approach to the nonequilibrium statistical mechanics of a system of particles is based on the Liouville equation and the derived hierarchy of equations for phase-space probability functions involving increasingly more particles, known as the BBGKY hierarchy (BBGKY stands for Bogoliubov-Born-Green-Kirkwood-Yvon). Every successive equation of this hierarchy gives the evolution of an -particle correlation function but involves the -particle correlation function. Therefore, to have a manageable set of equations, it is necessary to close the hierarchy at some level, by assuming some relation between the -particle correlation function and lower order correlation functions. Closure approaches are also much employed in the statistical theory of turbulence, and indeed the Kolmogorov theory can be understood as a closure approach with an additional scaling hypothesis [66]. Peebles and collaborators [4] follow a similar path, considering density correlation functions in addition to velocity correlation functions, which are the only ones present in incompressible turbulence.
In fact, the scaling hypothesis of Peebles et al. focuses on the reduced two-point density correlation function , which fulfills a simple equation that just expresses the conservation of particle pairs:where is the mean relative peculiar velocity of particle pairs. Of course, the presence of makes this equation not closed. To proceed, it is necessary to assume closure relations not only for (26) but also for the full set of equations. Nevertheless, it is possible to obtain information about only from the conservation of particle pairs, using scaling arguments [4]. However, these arguments are restricted to the Einstein-de Sitter cosmological model, with null curvature and therefore without characteristic scale, so that and are power laws, namely, and . Furthermore, the initial spectrum of density fluctuations must be a power-law.
Given the absence of characteristic scales, Peebles and collaborators propose similarity solutions of the equations, namely, solutions that are functions of the sole variable , with an exponent that is to be determined. In particular, they seek a similarity solution of (26). As this equation is not closed, one has to assume a relation between the expectation value of the relative velocity and and, furthermore, that this relation adopts a simple scaling form in the linear or very nonlinear regimes, namely, in the cases or . In the linear case, one knows that , which is in accord with the Zeldovich approximation, , with for the Einstein-de Sitter cosmology ( is the exponent of the Fourier power spectrum of density fluctuations). Therefore, one deduces that .
Of course, the problematic case is the very nonlinear regime. To deal with it, Peebles and collaborators proposed the stable clustering hypothesis, namely, that the average relative velocity of particle pairs vanishes in physical (not comoving) coordinates [4]. This hypothesis may seem natural for tightly bound gravitational systems. It plays a role that is somewhat analogous to Kolmogorov’s hypothesis of constant dissipation rate in turbulence. Indeed, the constancy of provides a link between large and small scales, and so does the stable clustering hypothesis. In particular, substituting in (26) and, furthermore, taking into account the fact that , it is easy to solve the equation and obtain , where . In terms of variables or , where is the physical distance,Naturally, as regards the -dependence, (27) is just a particular case of (3), in which is expressed in terms of , which can be derived from the initial power spectrum exponent (in the Einstein-de Sitter cosmology), giving . As regards the dependence on , that is to say, on time, (27) has new content, which can be generalized.
First, let us remark that (27) can be deduced without invoking the equation for conservation of particle pairs or closure relations, by means of a heuristic argument also due to Peebles [5, §22]. This argument is based on the estimation of the density of a bound system, which can be assumed to give the average conditional density, namely, the average density at distance from an occupied point. This statistic is expressed as where is the background cosmic density. For a bound system, , and the density of a system of size is the conditional density . The average conditional density is especially useful to define the extremely nonlinear limit, in which and while its product stays finite. This limit actually allows one to dispense with the transition to homogeneity and study directly a fractal mass distribution [3, 18].
To estimate the density of a bound system, one can follow the evolution of an overdensity of comoving size , which initially grows as . Over a time , the overdensity is of the order of magnitude of the background density and the linear theory is no longer valid [notice that this time is the “coalescence time” that can be deduced in the adhesion model by inverting the “coalescence length” expression ]. In terms of the physical size , the time for nonlinearity fulfills the relation , which, solving for , gives . At this time, the background density is , and our overdensity, as approximately twice the background density, follows the same rule. In fact, a factor of twice the background density seems hardly sufficient for the formation of a bound system; for example, the spherical collapse model favors a factor that is about 180 [33]. The exact value of this factor is irrelevant and the important point is that the self-similar evolution of an overdense blob of size leads to a bound system with density . Therefore, interpreting this density as a conditional density, one obtains in accord with (27).
We can draw more general conclusions by only appealing to the stable clustering hypothesis and without assuming any scaling in the cosmological evolution. Indeed, if particle pairs keep, on the average, their relative physical distances, then the average conditional density is a function that does not depend on time. Therefore, , where is an arbitrary function. This form of is also obtained as the general solution of (26) under stable clustering conditions, namely, and , which allow us to write (26) as One can check that solves the equation by direct substitution.
Of course, the stable clustering hypothesis can be questioned: Peebles says that “merging might be expected to dissipate the clustering hierarchy on small scales” [5, §22], and indeed the effect of merging is often presented as an argument against the hypothesis. However, this argument disregards that stable clustering occurs only on the average, that is to say, that merging can be accompanied by splitting, keeping the average physical distance of particle pairs constant. Actually, the conservation of energy in the local inertial frame demands it so. We have used this argument for introducing the stochastic adhesion model in Section 3.4, a model in which the merging and subsequent fragmentation of mass condensations take place in a stable manner. The stable clustering hypothesis in cosmology is actually related to the possibility of studying, in the local inertial frame and without a cosmological model, the gravity laws and, specifically, the laws of bound systems, a possibility of which Newton certainly took advantage. (The problem of the influence of the global cosmological expansion on local dynamics is still being discussed [78] and is connected with another polemic issue, namely, back-reaction [79].)
Accepting the stable clustering hypothesis, we can further argue that is singular as and that it must have, in particular, a power-law form . But this form of is independent of the hypothesis. Notice that the average conditional density has dimensions of mass density, so that a finite limit of in the limit would define a universal mass density. However, there is no universal mass density in the theory of general relativity. (Adding quantum mechanics, we have the universal mass density kg/m3 ( is the Planck constant). In this regard, a large but reasonable universal mass density would constitute a naturalness problem analogous to but less serious than the problem of the cosmological constant (or vacuum energy density) [80]. Of course, reasonable universal mass densities can be obtained by adding other universal constants, e.g., the nucleon mass, the electron charge) On the other hand, a universal mass density could also be given by the initial conditions, but the only density in them is , which is time dependent and too small: it is actually the limit of the average conditional density for . Therefore, it is natural to assume that is singular as .
A power-law singularity is not the only possibility but is the most natural possibility. In fact, the power-law form of is required by the notion of hierarchical clustering (Section 2.2), which is related to the full form of Peebles’ closure approach. If we dispense with the full form and just keep the stable clustering with a power-law form of , there is no way to obtain the value of the exponent , because there is no reason to relate it to the initial conditions, which do not need to be self-similar. The evolution of the FLRW cosmological model does not need to be self-similar either. Therefore, if could be derived from the initial conditions, this derivation would present an awfully difficult problem. A possible option that is conceptually different is that should be determined by a nonlinear eigenvalue problem [81]. One example of the solution of a problem of this type is provided by Polyakov’s nonperturbative theory of forced Burgers turbulence [74]. Another example, which is even more relevant, is provided by Gurevich and Zybin’s theory [41], described below.
Let us remark that the presence of length scales, in the FLRW cosmological model or in the initial spectrum of density fluctuations, does not rule out a sort of space-time similarity in the strongly nonlinear regime and let us remark that such a situation has already been shown to arise in the stochastic adhesion model, in Section 3.4. Naturally, the zero-size mass condensations of the adhesion model are not realistic, so it is necessary to analyze what type of mass concentrations can form under the effect of Newton’s gravity.
The formation of power-law density singularities in collisionless gravitational dynamics has been proved by Gurevich and Zybin’s analytic theory [41]. This study is especially interesting, since it brings to the fore the role of nonlinear dynamics and chaos in the formation of gravitational structure. They take an isolated overdensity of generic type, namely, a generic local maximum of the density field (assumed to be differentiable) and study its evolution. This evolution leads to a density singularity in a dynamical time and, after that, it consists in the subsequent development of a multistreaming flow with an increasing number of streams that oscillate about the origin, within a region defined by the first caustic. In this region, the multiplicity of streams in opposite directions gives rise to strong mixing, which leads to a steady state distribution with an average density. This average density has a power-law singularity at the origin, with an exponent that can be calculated (under reasonable assumptions). This exponent is independent of the particular form of the initial overdensity.
In the case of spherical symmetry, Gurevich and Zybin find It corresponds to a diverging gravitational potential, but the divergence is weak and does not induce large velocities. Let us recall from Section 2.4 that the local mass concentration exponent is found in cosmological -body simulations and the SDSS stellar mass distribution and notice that it corresponds to in the case of spherical symmetry. Arguably, the strongest mass concentrations have spherical symmetry. The density is also the singular isothermal density profile, where thermal equilibrium is represented by an -independent average velocity and the Maxwell distribution of velocities. An -independent circular velocity of stars is observed in the outskirts of spiral galaxies (this observation is, of course, one of the motivations for the existence of dark matter); so the isothermal density profile may have a role on these small scales [5, §3]. However, the collisionless gravitational dynamics does not lead to thermal equilibrium. In particular, the distributions of velocities in the steady states found by Gurevich and Zybin are very anisotropic [41].
For the sake of generality, Gurevich and Zybin treat the case without full spherical symmetry, that is to say, with low ellipticity, so that transverse velocities are small. Remarkably, this case leads to a singularity that is not exactly a power-law but can be well approximated by one, namely, with 1.7–1.9. Furthermore, this form is not very sensitive to the degree of ellipticity. The opposite case, of very large ellipticity, can be approximated by one-dimensional or two-dimensional solutions. Gurevich and Zybin only consider one-dimensional flows and derive . In this case, it is possible to study in some detail the nonlinear caustic oscillations and the role of transfer of energy successively to higher harmonics [41]. The picture is analogous to the Richardson energy cascade in the theory of incompressible turbulence [66]. As regards the one-dimensional power-law exponent 4/7, it has to be noticed that Aurell et al. [49] analyze a soluble example of one-dimensional collapse, with different initial conditions, and obtain a less singular value.
In general, that is to say, in the collapse in dimensions higher than one and without symmetry, it is commonly accepted that power-law singularities must appear, but there is no agreement about the exponent. This question is related, of course, to the problem of halo density profiles in the halo model, examined in Section 2.5. We have seen that a range of exponents has been considered, relying on various arguments and on the results of simulations. Let us notice that halo models usually assume a bias towards quasi-spherical mass concentrations and, therefore, towards more singular values of the mass concentration exponent . Here we have to consider all singularities on the same footing.
If Gurevich and Zybin’s theory were to be completed by adding the case of two-dimensional collapse, then one could think of obtaining a universal exponent by averaging the exponents over suitable positions of singularities at nodes, filaments, and sheets of the cosmic web. Gurevich and Zybin indeed consider the problem of how to integrate their results for isolated singularities in a hierarchical structure. They divide this structure into an intermediate nonlinear range, set to 50–100 Mpc, where the cosmic web (“cellular”) structure appears, and a lower-scale range, where steady nondissipative gravitational singularities appear [41]. However, this is an artificial restriction of the cosmic web structure, amounting to a separation into a regular “cellular” structure and a sort of halos, which are smooth, except for their central singularities. In fact, a distribution of isolated and smooth halos is only an approximation of the real cosmic web [27]. The halos are part of the cosmic web structure rather than separate entities, and they may not be smooth [31].
In particular, Gurevich and Zybin’s treatment of the collapse of isolated overdensities ignores that maxima of the initial density field are not isolated: in the “bottom-up” model of growth of cosmological structure, isolated density maxima of the initial density field only exist provided that this field is coarse grained, which introduces an artificial scale [82]. The resulting halos, with central singularities but smooth profiles, belong in the model of fractal distribution of halos, which is just a coarse-grained approximation of a multifractal model [27, 31].
In fractal geometry, the density is not a suitable variable, because it only takes the values zero or infinity, as explained in Section 2.2. Besides, in a multifractal model, with a range of local exponents , the average of local exponents over spatial positions takes a special form: for global magnitudes such as the moments , each is dominated by a specific set of singularities, with local dimension . The bulk of the mass belongs to the mass concentrate, with exponent , but the singularities of this set are not particularly strong in the real cosmic web (, see Figure 5). According to the argument, given above, that the density of bound systems is to be identified with the average conditional density and, hence, is connected with , we may think that the important exponent is . Certainly, the focus is usually placed on the exponent in (3), which is the one most accessible to measures and the one that Peebles’ closure approach tries to determine. However, the correlation dimension does not directly determine the corresponding local dimension , unless the multifractal spectrum is known [in general, ].
These difficulties do not arise if we have only one fractal dimension. In fact, Peebles’ closure approach only employs the three-point correlation function, in addition to the two-point correlation function, and assumes the hierarchical form for it [4, 33], which is valid for a monofractal but not for a general multifractal. However, the present data support full multifractality (Figure 5) and so does Gurevich and Zybin’s analytic theory, unless it is restricted to quasi-spherical mass concentrations. Moreover, even the intuitive image of the cosmic web, consisting of sheets, filaments, and knots, shows that it could hardly involve only one dimension. But the three-dimensional morphology of the cosmic web is by no means the cause of multifractality: one-dimensional cosmological -body simulations already show multifractality [43]. It seems that we are obliged to accept that we have to deal with a range of fractal dimensions, namely, with a full multifractal spectrum that we cannot quite model yet.
Multifractality of the mass distribution formed by gravitational clustering is surely a universal phenomenon, linked to the intermittency of highly nonlinear processes. Nonlinear gravitational clustering is indeed a phenomenon very related to turbulence, which began to be studied somewhat earlier but is in a similar state of development. The specific nature of gravity leads to some helpful constraints, chiefly, the stable clustering hypothesis, but this is not sufficient for a solution and specific nonlinear methods are needed [33]. The guiding idea for studying nonlinear phenomena like turbulence or gravitational clustering is surely the scaling symmetry, which is naturally motivated by the observed multifractality. The renormalization group prominently features among the analytic methods based on scaling, as has been advocated before [83]. The efficacy of this method is attested by its results for the stochastic adhesion model, specifically, the results obtained employing nonperturbative formulations of the renormalization group [84]. Stochastic gravitational clustering demands a random noise with an arbitrary exponent, as proposed by Antonov [77], who only applies perturbation theory. One could try to determine the arbitrary exponent from universality constraints, like in [71]. The stochastic approach deserves to be explored further.
After considering sophisticated methods of nonlinear dynamics, let us briefly comment on a simple method of obtaining the mass function of objects formed by gravitational collapse, namely, the Press-Schechter method [33]. This method is based on the spherical collapse model, which is not realistic in three dimensions. In fact, Vergassola et al. [46] find that the Press-Schechter mass function does not agree with the mass function derived from the adhesion model in more than one dimension. They define this mass function referring it to mass condensations in balls of fixed size, given that there are extended objects of arbitrary size, such as filaments, in more than one dimension. With a similar definition of the mass function, the form of this function has been studied in the distributions generated by cosmological -body simulations [28, 29]. The result is a mass function of Press-Schechter type but with a fixed power-law exponent that is independent of (the initial power spectrum exponent). This is contrary to the Press-Schechter approach.
One can see that the Press-Schechter approach is, in a way, analogous to the Peebles approach, in the sense that the latter also finds, under questionable assumptions, a simple relation, in this case, a relation between and another important exponent, namely, the -exponent . However, it seems that simple approaches to the nonlinear gravitational dynamics that attempt to bypass its true complexity have limited scope.
4. Discussion
We have reviewed the main ideas and theories that have led to the current understanding of the geometry of the cosmic structure. The appealing denomination of it as a “cosmic web” is most appropriate in regard to its three-dimensional appearance, but some important properties, related to its fractal geometry, are independent of the web morphology and already appear in one-dimensional cosmological models, in which the geometry is much simpler.
We have also considered halo models, very briefly. Halo models are based on the statistical properties of discrete point distributions rather than on the geometrical analysis of continuous mass distributions. However, we have shown that the correlations between points cannot be too simple and must consist of a hierarchy that allows us to naturally connect them with fractal models. Of course, this connection can actually take place only in the limit of an infinite number density of points, in which we obtain a continuous mass distribution.
The abstract study of the geometry of mass distributions has been essentially a mathematical subject. Since it is well developed, we propose to take advantage of this body of knowledge and we argue that the cosmic web structure must be studied as a strictly singular but continuous mass distribution. The geometry of these distributions is not simple even in one dimension, because they have nonisolated singularities. However, in the one-dimensional case, every mass distribution can be described in terms of the mass distribution function, which is just a monotonic function that gives the mass distribution by differentiation (it must be differentiable almost everywhere). If the function is continuous, so is the mass distribution. It is not easy to generalize this construction to higher dimensions, where various geometrical and topological issues arise. At any rate, a study of the type of singularities in a continuous mass distribution is always possible. This is the goal of multifractal analysis, which is based on the analysis of the behavior of fractional statistical moments when . In cosmology, only the integral moments are normally considered, but they do not provide sufficient information.
The simplest strictly singular and continuous mass distribution consists of a uniform mass distribution on a fractal set, namely, on a self-similar set of the type of the Cantor set. This mass distribution has just one kind of singularities. Therefore, it is a monofractal, described by just one dimension, the Hausdorff dimension of the fractal set. This type of fractal has a sequence of individual empty voids that is characterized by a particular form of the Zipf law.
Next in complexity, we consider a bifractal distribution, specifically, the type of distribution generated by the one-dimensional adhesion model. It is actually a noncontinuous distribution with some interesting properties; namely, it does not contain any properly fractal set, with nontrivial Hausdorff dimension, but it is obviously self-similar. If we disregard its void structure, it just consists of a collection of point-like masses, whose magnitudes follow a power-law, which is one cause of its self-similarity. Another cause is the spatial distribution of the masses, which gives rise to a self-similar void structure. In fact, this structure is crucial for the fractal nature, or rather bifractal, of this mass distribution, because the set of point-like masses, with , has nothing fractal to it, and the fractality is given by the other point of the bifractal spectrum, at .
The transition from the preceding bifractal to a full multifractal is best perceived by focusing on the mass distribution function: let us substitute the points of discontinuity, that is to say, the points with vertical segments in the graph, which correspond to , by points of nondifferentiability, such that the derivative grows without limit. We so obtain a continuous monotone function that has weaker singularities, with . This may seem an artificial procedure that is difficult to be realized in practice. On the contrary, any generic mass distribution is of this type. Moreover, these distributions normally have a nontrivial multifractal spectrum. If the graph of the function, besides not having vertical segments, does not have horizontal segments, then there are no individual empty voids and the multifractal is nonlacunar. In this case, the voids look like the voids of the adhesion model distribution, but they have a more complex structure, because the multifractal spectrum contains a full interval with .
Naturally, the large-scale mass distribution is three-dimensional and, therefore, more complex, with notable morphological features. These features, in particular the shape of voids, inspired the early models of the cosmic web, e.g., the Voronoi foam model. Nowadays, important information on the structures of matter clusters and voids can be obtained from the multifractal analysis of -body simulations of the dynamics of cold dark matter alone or with baryonic matter, in combination with the multifractal analysis of galaxy catalogues. The most obvious facts are, first, the self-similarity of the structures and, second, that the multifractal spectra are incompatible with monofractal distributions. We also notice the concordance of the multifractal spectra of cold dark matter and baryonic matter, to the extent allowed by the quality of the data. This concordance supports the hypothesis of a universal multifractal spectrum for the structure generated by gravitational clustering of both dark and baryonic matter.
Of course, this hypothesis challenges us to explain how such universality arises. Some features of the multifractal spectrum can be easily related to natural properties of the mass distribution. For example, is related to the nonlacunar structure or full support of the mass distribution. It can be explained by the adhesion model, which is a reliable approximation for the rough structure of voids (although not good enough for predicting the detailed multifractal spectrum for ). On the opposite side of the multifractal spectrum, we have that . This value disagrees with the adhesion model prediction (), and, indeed, the adhesion model is a poor approximation for the formation of strong mass concentrations. However, the basic condition of boundedness of the gravitational potential imposes . In fact, a monofractal of dimension equal to one was proposed by Mandelbrot for a similar reason, namely, the moderate value of cosmic velocities [3, §9]. But certainly the large-scale mass distribution is not monofractal, and the argument, properly considered, only sets a lower bound to the local fractal dimension.
To go beyond these two properties of the mass distribution, we need to understand better the dynamics of structure formation in cosmology. It is useful to begin with a further analysis of the adhesion model, that is to say, with a study of Burgers turbulence, because the hydrodynamic type of nonlinearity of the dynamical equations proves that the methods of the theory of turbulence are appropriate. We consider the Richardson energy cascade, the Kolmogorov hypothesis, and the intermittency of the velocity field. The Burgers turbulence with scale-invariant initial conditions lends itself to a fairly complete analysis and one can derive the form of the function that encodes the type of intermittency.
However, a treatment of the dynamics of structure formation that only involves the velocity field (because it is based on the Zeldovich approximation) is bound to be insufficient. A crucial problem of the adhesion model is that decaying Burgers turbulence does not conserve energy. This can be amended by studying forced Burgers turbulence, in the context of stochastic dynamics. The stochastic Burgers (or KPZ) equation gives rise to a scale-invariant asymptotic state that is independent of the initial conditions. The techniques based on dynamical scaling are powerful, but perturbation theory cannot deal with the three-dimensional problem and one has to apply nonperturbative methods. Simple universality arguments lead to a constraint on and hence to a reasonable range for the two-point correlation exponent . But it is doubtful that these arguments can be pushed much further while keeping within the scope of Burgers turbulence.
Therefore, we must keep in mind the methods of turbulence but deal with the full gravitational dynamics. A lot of good work has been done on nonlinear gravitational dynamics and we cannot even mention most of it, so we have limited ourselves to two relevant approaches: the classic Peebles closure approach and Gurevich-Zybin’s theory of nondissipative gravitational singularities. The full form of Peebles’ closure approach involves some questionable hypotheses, but we have focused precisely on the stable clustering hypothesis for the highly nonlinear regime. This most natural hypothesis gives a prominent role to the conditional mass density and, hence, to fractal geometry. For the calculation of power-law exponents, namely, fractal dimensions, the stable clustering hypothesis falls short and additional assumptions are necessary.
Definite power-law exponents of gravitational singularities are provided by the Gurevich-Zybin theory, which consists in an appropriate treatment of multistreaming in collisionless gravitational collapse. However, this treatment can only obtain exponents for singularities that form in isolation. We find no simple way to relate the exponents obtained by this theory with the dimensions in the multifractal spectrum of the real cosmic web. At any rate, the derivation of power-law singularities from the full dynamics of gravitational collapse is certainly an advance over the zero-size singularities of the adhesion model.
In conclusion, the problem of explaining the geometry of the cosmic web on a dynamical basis has only been solved partially. This problem seems to have a similar status to the problem of the geometry of the flow in fluid turbulence. Hopefully, the current developments in the study of nonlinear dynamics will bring further progress.
Conflicts of Interest
The author declares that there are no conflicts of interest regarding the publication of this paper.