PPrepared for submission to JCAP
Caustic Skeleton & Cosmic Web
Job Feldbrugge a Rien van de Weygaert b Johan Hidding b,c
JoostFeldbrugge d a Perimeter Institute for Theoretical Physics, University of Waterloo,Waterloo, Canada b Kapteyn Astronomical Institute, University of Groningen,Groningen, The Netherlands c Netherlands eScience Center, Amsterdam, The Netherlands d JFA Feldbrugge Studios,Lettelbert, The Netherlands a r X i v : . [ a s t r o - ph . C O ] M a r -mail: [email protected] Abstract.
We present a general formalism for identifying the caustic structure of a dynamically evolvingmass distribution, in an arbitrary dimensional space. The identification of caustics in fluidswith Hamiltonian dynamics, viewed in Lagrangian space, corresponds to the classification ofsingularities in Lagrangian catastrophe theory. On the basis of this formalism we develop atheoretical framework for the dynamics of the formation of the cosmic web, and specificallythose aspects that characterize its unique nature: its complex topological connectivity andmultiscale spinal structure of sheetlike membranes, elongated filaments and compact clusternodes. Given the collisionless nature of the gravitationally dominant dark matter compo-nent in the universe, the presented formalism entails an accurate description of the spatialorganization of matter resulting from the gravitationally driven formation of cosmic structure.The present work represents a significant extension of the work by Arnol’d et al. [11], whoclassified the caustics that develop in one- and two-dimensional systems that evolve accordingto the Zel’dovich approximation. His seminal work established the defining role of emergingsingularities in the formation of nonlinear structures in the universe. At the transition fromthe linear to nonlinear structure evolution, the first complex features emerge at locationswhere different fluid elements cross to establish multistream regions. Involving a complexfolding of the 6-D sheetlike phase-space distribution, it manifests itself in the appearance ofinfinite density caustic features. The classification and characterization of these mass elementfoldings can be encapsulated in caustic conditions on the eigenvalue and eigenvector fields ofthe deformation tensor field.In this study we introduce an alternative and transparent proof for Lagrangian catas-trophe theory. This facilitates the derivation of the caustic conditions for general Lagrangianfluids, with arbitrary dynamics. Most important in the present context is that it allows usto follow and describe the full three-dimensional geometric and topological complexity of thepurely gravitationally evolving nonlinear cosmic matter field. While generic and statisticalresults can be based on the eigenvalue characteristics, one of our key findings is that of thesignificance of the eigenvector field of the deformation field for outlining the entire spatialstructure of the caustic skeleton emerging from a primordial density field.In this paper we explicitly consider the caustic conditions for the three-dimensionalZel’dovich approximation, extending earlier work on those for one- and two-dimensional fluidstowards the full spatial richness of the cosmic web. In an accompanying publication, we applythis towards a full three-dimensional study of caustics in the formation of the cosmic web andevaluate in how far it manages to outline and identify the intricate skeletal features in thecorresponding N -body simulations. ontents A family 134.1.1 The trivial A class 134.1.2 The A caustics 144.1.3 The A caustics 144.1.4 The A ± points 154.1.5 The A caustics 174.1.6 The A ± points 184.1.7 The A caustics 204.2 The D family 214.2.1 The D caustics 214.2.2 The D ± points 244.2.3 The D caustics 244.3 Caustic conditions: physical significance 254.4 Spatial Connectivity: Singularities and Eigenvalue Fields 27 A -family caustics 417.2.2 Evolving D -family caustics 42– i – Discussion & Conclusions 42
A Zel’dovich approximation 50B Non-diagonalizable deformation tensors 51C Shell-crossing conditions: coordinate transformation 52D Lagrangian maps and Lagrangian equivalence 52
D.1 Symplectic manifolds and Lagrangian maps 52D.2 Displacement as Lagrangian map 53D.3 Lagrangian map germs 53D.4 Gradient maps 54D.5 Arnol’d’s classification of Lagrangian catastrophes 54
E Caustic conditions of the normal forms 55
Caustics that emerge in fluid flows are best studied in a Lagrangian space. They are im-portant features, marking the positions where fluid elements cross and multi-stream regionsform. These caustics can be associated to the regions with infinite density, corresponding tolocations where shell-crossing occurs. In the present study, we concentrate specifically on therole of caustics in the formation of the cosmic web. The gravitationally driven formation ofstructure in the universe is dominated by the dark matter component. Given its collisionlessnature, the formalism that we present in this study entails an accurate description of thespatial structure that emanates as a result of its dynamical evolution. The emerging causticseven have a direct physical impact on the baryonic matter, given its accretion into the gravi-tational potential wells delineated by the evolving dark matter distribution. Notwithstandingthis cosmological focus, the caustic conditions and mathematical formalism that we have de-rived for this are of a more generic nature, with a validity that extends to all systems whichallow for a Lagrangian description.The cosmic web is the complex network of interconnected filaments and walls into whichgalaxies and matter have aggregated on Megaparsec scales. It contains structures from afew megaparsecs up to tens and even hundreds of megaparsecs of size. The weblike spatialarrangement is marked by highly elongated filamentary and flattened planar structures, con-necting in dense compact cluster nodes surrounding large near-empty void regions. As borneout by a large sequence of N-body computer experiments of cosmic structure formation (e.g.[60, 65, 69]), these web-like patterns in the overall cosmic matter distribution do represent auniversal but possibly transient phase in the gravitationally driven emergence and evolutionof cosmic structure (see e.g. [4, 23]). In singularity theory, a caustic is the curve of critical values the Lagrangian mapping q → x t ( q ) . – 1 –ccording to the gravitational instability scenario [54], cosmic structure grows from tinyprimordial density and velocity perturbations. Once the gravitational clustering process hasprogressed beyond the initial linear growth phase, we see the emergence of complex patternsand structures in the density field. The resulting web-like patterns, outlined by prominentanisotropic filamentary and planar features surrounding characteristic large underdense voidregions, are therefore a natural manifestation of the gravitational cosmic structure formationprocess.The recognition of the cosmic web as a key aspect of the emergence of structure in theUniverse came with early analytical studies and approximations concerning the emergenceof structure out of a nearly featureless primordial Universe. In this respect the Zel’dovichformalism [72] played a seminal role. The emphasis on anisotropic collapse as agent forforming and shaping structure in the Zel’dovich "pancake” picture [42, 72] was seen as therival view to the purely hierarchical clustering view of structure formation. The successfulsynthesis of both elements in the cosmic web theory of Bond et al. [15] appears to provide asuccesful description of large scale structure formation in Λ CDM cosmology. The cosmic webtheory emphasizes the intimate dynamical relationship between the prominent filamentarypatterns and the compact dense clusters that stand out as the nodes within the cosmicmatter distribution [15, 25, 66]. It also implies that a full understanding of the cosmic web’sdynamical evolution is necessary to clarify how its structural features are connected in theintricate network of the cosmic web. To answer this question we need to turn to a fullphase-space description of the evolving matter distribution and mass flows.The Zel’dovich formalism [72] already underlined the importance of a full phase-spacedescription for understanding cosmic structure formation, however, with the exception of a fewprominent studies [11], the wealth of information in the full 6-D phase-space escaped attention.This changed with the publication of a number of recent publications [1, 30, 53, 58, 62] (foran early study on this observation see [21]) in which it was realized that the morphology ofcomponents in the evolving matter distribution is closely related to its multistream character.This realization is based on the recognition that the emergence of nonlinear structures occursat locations where different streams of the corresponding flow field cross each other.Looking at the appearance of the evolving spatial mass distribution as a 3D phasespace sheet folding itself in 6D phase space, a connection is established between the structureformation process and the morphological classification of the emerging structure. Caustics,which are the subject of this study, mark the regions where the cosmic web begins to form.Based on recent advances and insights, in this study we discuss the role of caustics in theformation of the cosmic web. By tracing the caustics during the formation of the cosmic webwe obtain a skeleton of the current three-dimensional large scale structure.Caustics in fluids with Hamiltonian dynamics, viewed in Lagrangian space, are classi-fied by Lagrangian catastrophe theory [6, 8, 46, 71]. Following this, these results were soonextended to fluids with generic dynamics [17]. For the classification of caustics emergingin the context of a one- and two-dimensional description of cosmic structure formation bythe Zel’dovich approximation, Arnol’d et al. [11] translated this into conditions on the dis-placement field. Following up on this seminal work, Hidding et al. [39] analyzed the overallmorphology and connectivity of caustics that emerge in a displacement field described by theone- and two-dimensional Zel’dovich approximation. The visual illustration of the emergingstructure, for a field of initially Gaussian random density and potential fluctuations, revealedhow the caustics spatially outline the spine of the cosmic web. Feldbrugge et al. [31] elab-orated this into an analytical evaluation of the statistical properties of caustics, assuming a– 2 – igure 1 : Illustration of caustic features in the cosmic web, and its relation with the corre-sponding density field. The N-body simulation is a CDM simulation in an Einstein-de SitterUniverse. On the basis of the initial flow deformation field, the caustics in the matter distri-bution have been identified. The detailed description of these follows in section 4. The redsheets represent the cusps ( A ) singularities which correspond to the walls or membranes ofthe cosmic web. The blue lines and the green points are the swallowtail ( A ) and butterfly( A ) singularities corresponding to the filaments and clusters of the large scale structure.The dark matter distribution in the N-body simulation is represented by a log density colourscheme.random Gaussian initial density field.In the current study we assess the caustics emerging in a one-parameter family of suf-ficiently differentiable maps x t : R (cid:55)→ R , mapping the initial mass distribution to thefinal mass distribution at time t . For practical considerations we consider the evolution ofa collisionless medium of matter in 6-dimensional phase-space. The collisionless Boltzmannequation, known as the Vlasov equation, describes the development of the phase-space density– 3 – ( x, v ) of the medium. In a gravitational field Φ , the phase-space density of mass elementswith velocity v at location x evolves according to ∂f∂t + v k ∂f∂x k − ∂ Φ ∂x k ∂f∂v k = 0 . (1.1)While the medium strictly speaking cannot be considered as a physical fluid, in the sense ofa medium characterized by continuously varying one-valued quantities in Eulerian space, wemight use the term “Lagrangian fluid” or “Vlasov fluid” for the dark matter medium. Forreasons of lucidity, in the remainder of this study we denote a “Vlasov fluid” shortly as “fluid”.Within this context, we give a novel proof of Lagrangian catastrophe theory and thecorresponding caustic conditions for three-dimensional Hamiltonian fluids. These conditionsare expressed in both the eigenvalue and the eigenvector fields of the mass flow deformationtensor. Moreover, our scheme allows us to extend these caustic conditions to fluids with non-Hamiltonian dynamics. Applied to the three-dimensional Zel’dovich approximation, theseconditions on the initial density field lead to a caustic skeleton of the cosmic web. In thisskeleton the walls, filaments and clusters of the large scale structure are directly related tothe A , A , A , D and D caustics of Lagrangian catastrophe theory. See figure 1 for anillustration of the caustic skeleton of the Zel’dovich approximation and a dark matter N -body simulation. A detailed analysis of the caustic skeleton of the Zel’dovich approximationand a comparison with N -body simulations is the subject of a follow-up paper [32].It should be emphasized that the eigenvalue fields of the mass flow deformation tensorhave, for a long time, been successfully used in Lagrangian studies of the cosmic web [24, 49,70]. In these studies, the clusters, filaments and walls are related to the number of eigenvaluesexceeding a threshold. The caustic skeleton here proposed complements their work in that itinclude the information of the eigenvector fields, which so far has been largely neglected.The paper begins section 2 with a concise description of Lagrangian fluid dynamics.The formation of caustics and derivation of the shell-crossing conditions for the occurrenceof multistream regions in a flow field is studied in section 3. These conditions are amongthe main results presented here. In section 4 we apply these shell-crossing conditions to theclassification of catastrophes, described in section 5, to derive the caustic conditions. Section 6discusses the relevance and significance of the caustic structure in the context of the evolvingcosmic mass distribution, and in particular the emergence and morphological structure ofthe cosmic web. Also, it discusses the further application and development of the causticformalism in a cosmological context, outlining the main elements of our project. In section7 we describe the dynamical framework resulting from the considerations above. Finally, insection 8 we summarize the results and discuss possible applications. There exist multiple approaches to fluid dynamics. In the Eulerian approach, the evolutionof the smoothed density and velocity fields is analyzed. The equations of motion of Eulerianfluids are relatively concise and give a reasonably accurate description of the mean flow in afluid element at a given location in the fluid. The Lagrangian view of particle flows is theappropriate tool for following the complex dynamical evolution of fluid elements, includingthe evolution of multi-stream regions and the emergence of caustics, where the caustics arethe critical values of the Lagrangian map. – 4 –n Lagrangian fluid dynamics, we assume every point in space to consist of a mass elementthat is moving with the fluid. Their motion is described by a Lagrangian map x t : L → E ,mapping the initial position q in the Lagrangian manifold L to the position x t ( q ) of the masselement in the Eulerian manifold E at time t . In the context of Lagrangian fluid dynamics, itis most convenient to describe the evolving fluid in terms of the displacement map s t definedby, s t ( q ) = x t ( q ) − q , (2.1)for all q ∈ L . For the Zel’dovich approximation [72] of cosmic structure formation the dis-placement field is given by s t ( q ) = − b + ( t ) ∇ q Ψ( q ) , (2.2)with the growing mode b + and the displacement potential Ψ (appendix A). The displacementpotential is proportional to the linearly extrapolated gravitational potential to the currentepoch φ , i.e. Ψ( q ) = 23Ω H φ ( q ) , (2.3)with H the current Hubble parameter and Ω the current total energy density. In this paperwe always assume the maps x t and s t to be continuous and sufficiently differentiable. Whilein the Lagrangian description a mass element has a constant mass, it may contract, expand,deform and even rotate. This is described in terms of the deformation tensor M , the gradientof the displacement field with respect to the Lagrangian coordinates of a mass element, M = ∂s t ∂q = M , M , M , M , M , M , M , M , M , . (2.4)While mass elements in a Lagrangian fluid are characterized by a few fundamental quantities,which characterize them and remain constant throughout their evolution, most physical prop-erties are basically derived quantities. A good example and illustration of a derived quantityis the density field. The density in a point x (cid:48) ∈ E is defined as the initial mass in the masselement times the ratio of the initial and final volume of the mass element. Formally, this isexpressed as a change of coordinates involving the Jacobian of the map x t , ρ ( x (cid:48) , t ) = (cid:88) q ∈ A t ( x (cid:48) ) ρ i ( q ) (cid:12)(cid:12)(cid:12)(cid:12) ∂x t ( q ) ∂q (cid:12)(cid:12)(cid:12)(cid:12) − = (cid:88) q ∈ A t ( x (cid:48) ) ρ i ( q ) (cid:12)(cid:12)(cid:12)(cid:12) I + ∂s t ( q ) ∂q (cid:12)(cid:12)(cid:12)(cid:12) − . (2.5)This can be written as ρ ( x (cid:48) , t ) = (cid:88) q ∈ A t ( x (cid:48) ) ρ i ( q ) | µ t ( q ) || µ t ( q ) || µ t ( q ) | , (2.6) Note that here we do not explicitly use a distinct notation for vector quantities: q and x t are vectorswhich in conventional cosmology notation are usually written as (cid:126)q and (cid:126)x t . Throughout this paper we use thenotation familiar to the mathematics literature. – 5 –ith A t ( x (cid:48) ) the points q in Lagrangian space L which map to x (cid:48) , i.e., A t ( x (cid:48) ) = { q ∈ L | x t ( q ) = x (cid:48) } , ρ i the initial density field and µ ti the eigenvalues of the deformation tensor M ( q ) , definedby M v i = µ i v i (2.7)with eigenvector v i . The equality in equation (2.6) applies to general deformation tensors ,since the characteristic polynomial of the deformation tensor can be expressed in terms of theeigenvalues χ ( λ ) = det (cid:20) ∂s t ∂q − λI (cid:21) = ( µ t − λ )( µ t − λ )( µ t − λ ) , (2.8)by which det (cid:20) I + ∂s t ∂q (cid:21) = χ ( −
1) = (1 + µ t )(1 + µ t )(1 + µ t ) . (2.9)By substituting derived quantities like density in the, often more familiar, Eulerian fluidequations, we may obtain a closed set of differential equations for the Lagrangian map x t orthe displacement map s t . Note that for practical reasons in this paper we will sometimessuppress the time index of the eigenvalue fields, i.e. µ i = µ ti .Equation (2.6) applies to a fluid with three spatial dimensions. For simplicity, we willrestrict explicit expressions to the -dimensional case . The arguments presented in thispaper straightforwardly generalize to a Lagrangian fluid with an arbitrary number of spatialdimension and it is straightforward to generalize equation (2.6) to d -dimensional fluids in d -dimensional space.The appearance of singularities in equation (2.6) is central to our discussion concerningthe nature of these singularities. They occur when a mass element reaches an infinite density.More formally stated, as we will see in section 3, an infinite density occurs when for at leastone of the i = 1 , . . . , d , µ i = 0 . (2.10)The regions, in which the mapping x t becomes degenerate and the density becomes infiniteare known as foldings , caustics or shocks . They mark important features in the Lagrangianfluid and are the object of study in this paper.While these eigenvalue conditions provide the necessary condition for a mass element topass through a caustic, and reach infinite density, it does not yield the full information neces-sary to infer the geometric structure, spatial connectivity and identity of the caustic. As masselements pass through a multistream region, the spatial properties of the flow will determinethe complexity of the folding of the phase-space sheet in which they are embedded. In thisstudy we demonstrate that the corresponding eigenvectors are instrumental in establishingthe spatial outline and identify of the corresponding caustics. This key realization emanatesfrom the so-called caustic conditions . Note that here we use the general convention to represent the deformation eigenvalue field, with µ i ( q ) the i -th eigenvalue of the deformation tensor, M ( q ) . This differs from the usual convention in cosmology to usethe time-independent representation of the deformation field in the context of the Zel’dovich approximation.Within this formalism, the eigenvalues λ i ( q ) of the deformation field ψ ij = ∂ Ψ( q ) /∂q i ∂q j , are related to theeigenvalues µ i ( q ) via the linear relation µ i ( q, t ) = − b + ( t ) λ i ( q ) , in which b + ( t ) is the growing mode growthfactor. See Appendix A for further details. Formally, it would be appropriate to describe the fluids as ( d + 1) -dimensional fluids, a combination oftheir embedding in a d -dimensional space along with their evolution along time dimension t . – 6 –hroughout our study, we assume that the displacement map s t is continuous and suf-ficiently differentiable. The corresponding eigenvalues are the roots of the characteristicpolynomial of the matrix M = ∂s t /∂q . Since the characteristic equation is a non-linear equa-tion, in principle the eigenvalues could develop singularities and become non-differentiable.However, it can be shown that the eigenvalues can be ordered such that they are continuous.Furthermore the eigenvalues will be assumed to be differentiable whenever the eigenvalues aredistinct. When two eigenvalues coincide, the eigenvalue fields may become non-differentiable. For fluids moving with no dissipation of energy, the Hamiltonian formalism may be applied.Hamiltonian fluids have a potential velocity field v = ∇ φ (2.11)with the velocity potential φ . The mass density ρ and the velocity potential serve as conjugatevariables for the Hamiltonian H , with the equations of motion ∂ρ∂t = + δ H δφ = −∇ · ( ρv ) ,∂φ∂t = − δ H δρ . (2.12)A simple example of a Hamiltonian is H = (cid:90) d x (cid:18) ρ ( ∇ φ ) + e ( ρ ) (cid:19) , (2.13)where e ( ρ ) is the internal energy as a function of density ρ . The first equation of motion inequation (2.12) is equivalent to the continuity equation, while the second equation impliesthe Euler equation ∂v∂t + v · ∇ v = − ρ ∇ p , (2.14)in which p is the pressure of the fluid. For a thorough discussion of fluid mechanics we referto the seminal volumes of [47], and [48]. For detailed and extensive treatments and analysesof Hamiltonian mechanics and Hamiltonian fluids, we refer to the reviews and textbooks by[7], [10], [35], [52], and [59]. The caustics mentioned above result from the folding of the phase space fluid. At the initialtime, t = 0 , the fluid has not yet evolved. The displacement map s is therefore the zero map(eqn. (2.1)), i.e., s ( q ) = 0 (3.1)for all q ∈ L . The map x ( q ) is one-to-one, i.e. each Eulerian coordinate x corresponds to oneLagrangian position q . Throughout the entire volume, the fluid only contains single-streamregions. As the fluid evolves and nonlinearities start to emerge, we see the development of multi-stream regions in the fluid. At the boundary of a multi-stream region, the volume ofa mass element vanishes and its density – following eqn. (2.6) – becomes infinite. At such– 7 – q q s q (cid:48) q (cid:48)(cid:48) C Tx ( q s ) (a) Lagrangian space L x x x t ( q s ) x t ( q (cid:48)(cid:48) ) x t ( q (cid:48) ) x t ( C ) x ( q s ) (b) Eulerian space E Figure 2 : The shell-crossing process of a curve C in a Lagrangian map x t . The left panelshows Lagrangian space, describing the initial positions of the fluid. The right panel showsEulerian space, describing the positions of the fluid at time t . The fluid undergoes shell-crossing in point q s on the curve C (red) at time t . The neighboring points q (cid:48) and q (cid:48)(cid:48) havepassed through the opposing segments of C . The Lagrangian mapping of the curve x t ( C ) (red) develops a non-differentiable point in x t ( q x ) , which is known as a caustic. The arrow T (blue) is the tangent vector of the curve C in point q s .locations in phase space the map x t ( q ) attains a n -to-one character, with n an odd positiveinteger ( n = 3 , , , . . . ). It means that at any one Eulerian location x , streams from n differentLagrangian positions cross.The key question we address here is that of inferring the conditions under which a masselement with Lagrangian coordinate q undergoes shell-crossing. Here we derive the necessaryand sufficient conditions for the process of shell-crossing to occur. These conditions are called shell-crossing conditions . They are the foundation on the basis of which we infer – in section 4– the related conditions on the displacement field for the occurrence of the various classes ofcaustics. These are called the caustic conditions . We infer the caustic conditions for genericas well as Hamiltonian fluid dynamics. A typical configuration resulting from the shell-crossing process – the name by which it isusually indicated – is illustrated in figure 2. It focuses on points q = ( q , q ) that lie ona smooth curve C in Lagrangian space L (fig. 2a). In this context, smooth refers to theassumption that the curve C is C continuous. At time t , the points on the Lagrangian curve C map to the variety x t ( C ) in Eulerian space E (fig. 2b) . The fluid in point q s undergoesshell-crossing at time t . The neighboring points q (cid:48) and q (cid:48)(cid:48) have passed through the opposing In algebraic geometry, a variety is the zero set of a function f , ie. the set of solutions x ∈ E such that f ( x ) = 0 . – 8 – x L q s x ( q s ) q (cid:48)(cid:48) x ( q (cid:48)(cid:48) ) q (cid:48) x ( q (cid:48) ) q x L q s x ( q s ) q (cid:48)(cid:48) x ( q (cid:48)(cid:48) ) q (cid:48) x ( q (cid:48) ) q x L q s x ( q s ) q (cid:48)(cid:48) x ( q (cid:48)(cid:48) ) q (cid:48) x ( q (cid:48) ) Figure 3 : Folding of a one-dimensional fluid in phase space C . The three panels show thetime evolution of the Lagrangian submanifold L (red) of the fluid in phase space. We trackthe evolution of two points ( q (cid:48) , x ( q (cid:48) )) , ( q (cid:48)(cid:48) , x ( q (cid:48)(cid:48) )) forming a multi-stream region and mark thepoint undergoing shell-crossing by ( q s , x ( q s )) . Left panel: the fluid – early in its evolution –consisting of a single-stream region. Middle panel: a fluid during the process of shell-crossing.Right panel: a fluid consisting of a multi-stream region.segments of C . As a result of this, the curve C develops a non-differentiable point in x t ( q x ) ,which is known as a caustic .In a time sequence of three steps, figure 3 illustrates the dynamical process that isunderlying the formation of the caustic at x t ( q s ) . The singularity at x t ( q s ) ∈ x t ( C ) forms asthe result of a folding process in phase space. We may appreciate the emerging structure whenassessing the fate of two neighboring points q (cid:48) , q (cid:48)(cid:48) ∈ C on both sides of q s . While the phasespace sheet x t ( C ) is folded, the points x t ( q (cid:48) ) and x t ( q (cid:48)(cid:48) ) turn around while passing through x t ( q s ) . In figure 3 we observe how the initially single-stream phase space sheet (lefthandpanel) morphs into a configuration marked by shell-crossing as different mass elements q pileup at the same Eulerian position x t ( q s ) (central panel). Subsequently, around x t ( q s ) we noticethe formation of a multi-stream region, with the presence of mass elements q (cid:48) having passedinto a region where mass elements from other Lagrangian locations q are to be found.To infer the shell-crossing conditions, we investigate a curve C in Lagrangian space alongwhich we have points q that will find themselves incorporated in a singularity at Eulerianposition x s ( q s ) . In the case of shell-crossing, points q near the Lagrangian location q s willmap onto the same Eulerian position x ( q s ) . The key realization is that this occurs as points q along a direction T tangential to C are all folded on to a single Eulerian position x s ( q s ) . Thistranslates the question of the shell-crossing condition into one on the identity of a tangentialdirection T ( q ) along which shell-crossing may or will occur. In other words, whether on aparticular curve C – or, more general, a manifold M – there are points q where along oneor more tangential directions T ( q ) to that curve or manifold shell-crossing may or will takeplace.Zooming in on two points q (cid:48) and q (cid:48)(cid:48) in the vicinity of the singularity point q s , we see thatas a result of the folding process the ratio of the distances of the two points in the Lagrangianand Eulerian manifold, must go to zero in the limit that we zoom in on points q (cid:48) and q (cid:48)(cid:48) alongthe Lagrangian curve C at an infinitesimal distance from q s , i.e. ∆ x | ∆ q | = (cid:107) x t ( q (cid:48) ) − x t ( q (cid:48)(cid:48) ) (cid:107)(cid:107) q (cid:48) − q (cid:48)(cid:48) (cid:107) → q (cid:48) , q (cid:48)(cid:48) → q s . (3.2)– 9 –he direct implication of this is, following equation (2.6), that the density in a caustic isinfinite: the volume of the mass element associated to q s vanishes at time t . In essenceit informs us that during shell crossing the points q near Lagrangian location q s , along thetangential direction T to the Lagrangian curve C , map onto the same Eulerian position x ( q s ) .This means that the norm of the directional derivative of x t along the tangential directionvanishes. In other words, along the non-zero tangent vector T along C , (cid:13)(cid:13)(cid:13)(cid:13) ∂x t ∂q T (cid:13)(cid:13)(cid:13)(cid:13) = 0 , (3.3)where ∂x t /∂q is the Jacobian of x t evaluated in q s (see figure 2a). This is equivalent torequiring that ∂x t ∂q T = 0 . (3.4)In terms of the displacement map s t , this condition can be expressed as T + ∂s t ∂q T = 0 , (3.5)with the Jacobian ∂s t ∂q also evaluated in q s . Subsequently consider the eigenvalues µ i andeigenvectors v i of the deformation tensor M = ∂s t ∂q , defined by M v i = µ i v i . (3.6)Under the assumption that the deformation tensor is diagonalizable , we can construct thediagonal matrix M d = diag ( µ , . . . , µ d ) and the eigenvector matrix V = ( v , . . . , v d ) . Foran analysis of the case of non-diagonalizable deformation tensors see appendix B. In threedimensions, with the eigenvalues µ i and eigenvectors v i = ( v i, , v i, , v i, ) , the diagonal matrix M d and eigenvector matrix V are given by M d = µ µ
00 0 µ , V = v , v , v , v , v , v , v , v , v , . (3.7)In terms of V and M d , condition (3.5) reduces to I + M ) VV − T = V ( I + M d ) V − T (3.8)since V is always invertible , using the identity MV = M ( v , . . . , v d ) = ( M v , . . . , M v d ) = ( µ v , . . . µ d v d ) = VM d . (3.9)We thus obtain the condition ( I + M d ) V − T = 0 , (3.10) Unless mentioned otherwise, we will assume all Jacobians to be evaluated in q s . In practice, the assumption of diagonalizability is not really restrictive: non-diagonalizable matrices areunstable, which means that they can be turned diagonalizable by means of a small perturbation in the initialconditions. That is to say, the eigenvectors can always be chosen to be linearly independent. – 10 –hich holds for general diagonalizable deformation tensors. Note that the rows of V − consistof the dual vectors { v ∗ i } of the eigenvectors { v i } , defined by v i · v ∗ j = δ ij for all i and j .Explicitly, this means that V − in three dimensions is given by V − = v ∗ , v ∗ , v ∗ , v ∗ , v ∗ , v ∗ , v ∗ , v ∗ , v ∗ , , (3.11)with v ∗ i = ( v ∗ i, , v ∗ i, , v ∗ i, ) . The product V − T is the vector composed out of the inner productof these dual vectors with the tangent vector T , so that in three dimensions equation (3.10)reduces to (1 + µ ) v ∗ · T (1 + µ ) v ∗ · T (1 + µ ) v ∗ · T = 0 . (3.12)This represents the proof for the shell-crossing condition for one-dimensional subman-ifolds. It states the condition for the tangential direction T along which Lagrangian pointsget folded into an Eulerian singularity point. The obtained condition is a telling expressionfor the central role of both the deformation eigenvalues and eigenvectors in determining theoccurrence of a singularity. Following the proof outlined in the previous subsection 3.1, we arrive at the following twotheorems stipulating the conditions for the formation of singularities by curves C and arbitrarymanifolds M in Lagrangian space L , Theorem: 1 A C continuous curve C ⊂ L forms a singularity under the mapping x t inthe point x t ( q s ) ∈ x t ( C ) ⊂ E , meaning that x t ( C ) is not smooth in x t ( q s ) , if and only if (1 + µ it ( q s )) v ∗ it ( q s ) · T = 0 (3.13) for all i = 1 , , . . . , dim ( L ) , with T a nonzero tangent vector of C in q s . Note that the derived caustic condition is independent of the dynamics of the fluid.In general, both the eigenvalue and eigenvector fields are complex-valued. For Hamiltonianfluids, the relation condition simplifies since the eigenvalue and eigenvector fields are forcedto be real-valued and the eigenvectors can be chosen to coincide with their dual vectors, i.e. v ∗ i = v i .A similar argument holds for higher dimensional submanifolds of L , e.g., sheets andvolumes. These manifolds can be n -dimensional, with n = 1 , . . . , for three-dimensionalfluids. Given an arbitrary manifold M ⊂ L we can consider all curves C ⊂ M passingthrough the point q s ∈ M . The variety x t ( M ) contains a singularity at x t ( q s ) if and only if at– 11 –east one such curve C ⊂ M gets folded under the map x t . Hence for an arbitrary submanifold M , we should consider the one-dimensional shell-crossing condition for the subset of vectors T ∈ M . This proves the general shell-crossing condition: Theorem: 2
A manifold M ⊂ L forms a singularity under the mapping x t in the point x t ( q s ) ∈ x t ( M ) ⊂ E at time t , meaning that x t ( M ) is not smooth in x t ( q s ) , if and only ifthere exists at least one nonzero tangent vector T ∈ T q s M satisfying (1 + µ it ( q s )) v ∗ it ( q s ) · T = 0 (3.14) for all i = 1 , , . . . , dim ( L ) . From this theorem, we immediately observe that the eigenvectors v i are of key importance indetermining the nature of the singularity, in that the shell-crossing condition is not simplythat of µ i = 0 for at least one i . More explicitly, the shell-crossing condition says that µ it ( q s ) = 0 OR v ∗ it ( q s ) · T = 0 for all i , (3.15)indicating that, in addition to one or more eigenvalue constraints µ i = 0 , the shell-crossingcondition consists of complementary constraints. These single out those points q s where theeigenvectors v ∗ j ( q s ) (with j (cid:54) = i ) are orthogonal to a vector T that is restricted to be locatedin the plane tangent to the manifold M in which the singularity emerges. It is this constraintthat is instrumental in defining the area occupied by the corresponding caustic.Note that the shell-crossing conditions are manifestly independent of coordinate choices.While in general the eigenvalue and eigenvector fields generally do depend on the choice ofcoordinates, it can be shown that they are invariant if the corresponding coordinate transfor-mation is orthogonal and global. These transformations include rotations and translations.See appendix C for more details. In section 3, we inferred the general condition for shell-crossing. The condition establishes therelation between the eigenvalue and eigenvector fields of the deformation tensor in Lagrangianspace, and the Lagrangian regions that get incorporated in features of infinite density inEulerian space. Moreover, it allows us to establish the identity of the resulting singularity inEulerian space. T q s M , i.e. all tangent vectors T constrained to be located in the vector space T q s M . In other words,the one-dimensional shell-crossing condition is considered for all vectors T in the vector space of all tangentialvectors to the manifold M in q s ∈ M . – 12 –he stable singularities that emerge can be classified by Lagrangian catastrophe theoryin the A k , D k and E k series (see [9], [34] and [57]). This is described in some detail in section 5 . The A k series is of co-rank , in which co-rank is the number of independent directionsin which the Hessian is degenerate. The A k series corresponds to the caustics for which thedensity diverges due to only one eigenvalue. The D k series is of co-rank and correspondsto the points for which the density diverges due to two eigenvalue fields. The E k series is ofco-rank and corresponds to the points for which three eigenvalue fields. However, for three-dimensional fluids, the points for which all eigenvalues simultaneously satisfy this conditionare degenerate. For this reason we will not discuss them in the context of the present paper.In this section we apply the shell-crossing condition to three-dimensional Lagrangianfluids to obtain the caustics conditions which relate the classification of caustics to the eigen-value and eigenvector field. These conditions have not been derived in earlier work and arenecessary to perform a quantitative study of caustics in large scale structure formation. Insection 5, we summarize the classification of caustics in its traditional form and compare themto the caustic conditions derived here. A family The A family of caustics form when µ i = 0 (4.1)for one i . For diagonalizable deformation tensors, the eigenvector fields { v i } and their dualvector fields { v ∗ i } are linearly independent.For three-dimensional fluids, the A family consists of classes running from the trivial A class, corresponding to the points that never form caustics, the sheetlike A fold, thecurvelike A cusp, the A swallowtail, to the pointlike A butterfly singularity. A class The A class labels the points which never form caustics. According to the shell-crossingcondition, q s will form a singularity at time t if and only if there exists a nonzero tangentvector T ∈ T q s L for which (1 + µ i ( q s )) v ∗ i ( q s ) · T = 0 (4.2)for all i . The point q s will not satisfy this condition if µ i ( q s ) (cid:54) = 0 for all i since the threedual vectors { v ∗ i } of the (generalized) eigenvectors span the tangent space T q s L From the shell-crossing condition we therefore conclude that the three-dimensional variety A , A = { q ∈ L | µ ti ( q ) (cid:54) = 0 for all i and t } , (4.3)consists of the points never forming caustics. In this respect we should note that the displace-ment map at the initial time is the zero map, so that the eigenvalues at the initial time areequal to zero, i.e. µ i ( q ) = 0 for all q ∈ L . Since the eigenvalues are a continuous function oftime, for the cosmologically interesting case of potential flow the requirement for a point q tobelong to A is equivalent to µ ti ( q ) > − . The classification ultimately has its origin in the classification of Coxeter groups – 13 – .1.2 The A caustics Based on the discussion above, we may conclude that for a given i , i = 1 . . . , at time t thepoints A i ( t ) = { q ∈ L | µ ti ( q ) = 0 } (4.4)form a singularity. For three-dimensional fluids, the set A ( t ) forms a two-dimensional sheet,sweeping through space as the fluid evolves. These singularities can be associated to the A fold singularity class.From this, we conclude that the set of points which form a A fold singularity at a time t ∈ [0 , ∞ ) is given by A i = { q ∈ L | µ ti ( q ) = 0 for some t } . (4.5) A caustics Following up on the folding of the fluid to the A i singularity, the A i manifold itself maybe folded into a more complex configuration. The result is a so-called A singularity. Toguide understanding in the emergence of cusps we may refer to the eigenvalue contour mapof figure 4.To infer the identity of the A i caustic, we restrict the criterion for shell-crossing to pointson the A i manifold. In other words, we look for points q s on the surface of the sheetlike variety A i ( t ) that fulfill the criterion for shell-crossing.A point q s ∈ A i ( t ) forms a singularity if there exists a nonzero tangent vector T, T ∈ T q s A i ( t ) , orthogonal to the Span C { v ∗ j | j (cid:54) = i } . Writing the tangent vector T as a linearcombination of the eigenvectors v i , T = α v + α v + α v , (4.6)with α i ∈ C . The caustic conditions tell us that α j = v ∗ j ( q s ) · T = 0 for j (cid:54) = i . (4.7)Given that we know that the i th eigenvalue is real, µ i ∈ R , the eigenvector v i is also real.This means that this condition is satisfied if and only if the tangent vector T is parallel to v i . This is equivalent to the condition that v i is orthogonal to the normal n = ∇ µ ti of themanifold A i ( t ) in the point q s . Explicitly, this means that the inner product of n with v i isequal to 0, µ ti,i ≡ v i · ∇ µ ti = 0 . (4.8)Note that this is the condition that Arnol’d [8] found for the A line for the 2-dimensionalZel’dovich approximation. As we see from the derivation above, the condition is valid in anydimensional space and for general flow configurations.The points q forming a cusp at time t corresponding to eigenvalue field µ i is given by theone-dimensional variety A i ( t ) = { q ∈ L | q ∈ A i ( t ) ∧ µ ti,i ( q ) = 0 } . (4.9)Extrapolating this to the set of all points q that at some time t ∈ [0 , ∞ ) have belonged to orwill be incorporated in a cusp singularity defines a two-dimensional variety A i = { q ∈ L | q ∈ A i ( t ) ∧ µ ti,i ( q ) = 0 for some t } , (4.10)which is the assembly of all A i ( t ) over the time interval t ∈ [0 , ∞ ) .– 14 – q q s A i ( t ) A i v i n (a) Lagrangian space x x x t ( q s ) x t ( A i ( t )) x t ( A i ) (b) Eulerian space Figure 4 : The formation of a cusp ( A ) singularity in a Lagrangian map x t . The left panelshows Lagrangian space, describing the initial positions of the fluid. The right panel showsEulerian space, describing the positions of the fluid at time t . The fluid undergoes shell-crossing along the fold A i ( t ) (red) at time t . The fold gets mapped under the Lagrangianmap to x t ( A ) (red), which is folded into a cusp in the point x t ( q s ) corresponding to q s . Thecusp forms if and only if the normal n of A i ( t ) is orthogonal to the eigenvector field v i in q s .Over time, the cusp traces out the curve A i (blue) which is mapped to x t ( A i ) (blue). A ± points The topology of the sheetlike A i ( t ) variety changes as a function of time. These topologicalchanges occur at critical points of the corresponding eigenvalue field µ ti . It is at these pointswhere in Eulerian space we see the emergence of new features, the disappearance of featuresand/or the merging of features. The critical points are classified as cusp singularities.At minima of the µ i field, a feature gets created. At maxima, a feature gets annihilated.Particularly interesting points are the saddle points. In three-dimensional space, there aretwo classes of saddles in the eigenvalue field µ ti . The index 1 saddles have a Hessian signature ( − − +) , with 1 positive eigenvalue, while the index 2 saddles have a signature ( − + +) .Based on their impact on caustic structure, Arnol’d used a slightly different classificiationscheme, in which the distinguished between A ++3 , A + − and A −− points [8]. The A ++3 pointare identified with the minima , while the A + − points are the saddle points for which the A sheet intersects the two disjoint A sheets. This is illustrated in the upper left panel infigure 6. The additional A −− points correspond to saddle points for which the A sheet doesnot intersect the disjoint A sheets. Because this concerns a non-generic situation, we do nottreat it here. Also note that higher dimensional fluids will have additional A points.In the context of this paper we therefore use a slightly shorter notation for the maxima, Note that in Arnol’d’s notation, related to the Zel’dovich formalism (see appendix A), these are themaxima of the eigenvalue field – 15 – igure 5 : The creation/annihilation of a fold ( A ) sheet in a A +3 point. The upper threepanels show the unfolding of a A +3 singularity in Lagrangian space. The lower three panelsshow the corresponding unfolding in Eulerian space. The two panels on the left show thecusp ( A ) plane on which the cusps form. The middle panels show the appearance of a A +3 singularity in which a fold sheet is formed/removed. The right panels show the resulting fold( A ) sheet. The fold sheet gets folded into a cusp ( A ) curve (red). This configuration isknown as the Zel’dovich pancake (Zel’dovich 1970).
Figure 6 : The merger/splitting of a fold ( A ) sheet in a A − point. The upper three panelsshow the unfolding of a A − singularity in Lagrangian space. The lower three panels show thecorresponding unfolding in Eulerian space. The two panels on the left show two fold ( A )sheets, two cusp ( A ) curves (red) and the cusp ( A ) plane on which the cusps form. Themiddle panels show the merger/splitting of the two fold ( A ) sheets in a A − singularity. Theright panels show the resulting merged fold ( A ) sheet. This configuration is known as the Kissing Lips . – 16 – igure 7 : The merger/splitting of a fold ( A ) sheet in a A −− point. The upper three panelsshow the unfolding of a A −− singularity in Lagrangian space. The lower three panels showthe corresponding unfolding in Eulerian space. The two panels on the left show two fold( A ) sheets, and the cusp ( A ) plane on which the cusps form. The middle panels show themerger/splitting of the two fold ( A ) sheets in a A −− singularity. The right panels show theresulting merged fold ( A ) sheet with the corresponding cusp ( A ) curve.minima and saddles, classifying them as the cusp singularities A +3 an A − , A i +3 = { q ∈ L | q ∈ A i ( t ) ∧ µ ti ( q ) max-/minimum of µ ti at some time t } ,A i − = { q ∈ L | q ∈ A i ( t ) ∧ q saddle point of µ ti at some time t } . (4.11)Note that in this scheme, the saddle points with index and belong to the same singularityclass A i − . For an illustration of the A +3 , A − and A −− singularities, we refer to figures 5, 6and 7. From the caustics conditions we may directly infer that the A i ± points are located onthe A i variety. A caustics In Eulerian space the A i ( t ) variety gets folded in points associated with A swallowtail singu-larities. The identity of the points defining the variety A i ( t ) can be inferred by the applicationof the general shell-crossing condition (eqn. (3.14)) to the A i ( t ) variety (see figure 8). As aconsequence, the A i variety is defined as A i ( t ) = { q ∈ L | q ∈ A i ( t ) ∧ µ ti,ii ( q ) = 0 } , (4.12)with µ ti,ii ( q ) the inner product of the normal n = ∇ µ ti,i with the eigenvector v i , µ ti,ii ≡ v i · ∇ µ ti,i . (4.13)– 17 – q q s A i ( t ) A i A i v i n (a) Lagrangian space x x x t ( q s ) x t ( A i ( t )) x t ( A i ) x t ( A i ) (b) Eulerian space Figure 8 : The formation of a swallowtail ( A ) singularity in a Lagrangian map x t . The leftpanel shows the Lagrangian space describing the initial positions of the fluid. The right panelshows the Eulerian space describing the positions of the fluid at time t . The fluid undergoesshell-crossing along A i ( t ) (red) at time t . The fold gets mapped in Eulerian space, under theLagrangian map, to x t ( A ) (red), which is folded into a cusp in the point x t ( q s ) correspondingto q s . The cusp forms if and only if the normal n of A i ( t ) is orthogonal to the eigenvectorfield v i in q s . Over time, in Lagrangian space the cusp traces out the curve A i (blue) whichin Eulerian space is mapped to x t ( A i ) (blue). Since the cusp ( A i ) curve is tangential to thefold ( A ) curve in q s , the cusp curve x t ( A i ) forms a swallowtail ( A ) singularity. Over time,the swallowtail traces out A i (green), which in Eulerian space is mapped into x t ( A i ) (green).Integrated over time, the points on the varieties A i ( t ) trace out the 1-dimensional variety A i ,i.e. the 1-dimensional line A i is the set of all points A i ( t ) over the time interval t ∈ [0 , ∞ ) , A i = { q ∈ L | q ∈ A i ( t ) ∧ µ ti,ii ( q ) = 0 for some t } . (4.14) A ± points The topology of the variety A i ( t ) changes as a function of time. To this end, we identify thecritical points of the real field µ i,i , µ ti,i ≡ v i · ∇ µ ti . (4.15)Constraining the location of these singularities to the one-dimensional curvelike variety A i ( t ) ,and thus implicitly also to the two-dimensional membrane of the variety A i ( t ) , these A ± pointsmark the locations at which topological changes occur. They represent the sites at which wesee the birth of new singularities in Eulerian space, or the annihilation of and/or merging ofsuch features. These singularities are classified as swallowtail singularities.The birth or death of features on A i ( t ) takes place at maxima and minima of µ ti,i , andis identified with A i +4 singularities. The merging or splitting of features happens at the saddle– 18 – igure 9 : The creation/annihilation of a swallowtail ( A ) singularity in a A +4 point. Theupper three panels show the unfolding of a A +4 singularity in Lagrangian space. The lowerthree panels show the corresponding unfolding in Eulerian space. The two panels on the leftshow a fold ( A ) sheet. The middle panels show a A +4 point on the fold ( A ) sheet. The A +4 point leads to the creation/annihilation of two swallowtail ( A ) singularities. The rightpanels show the resulting cusp ( A ) curves and swallowtail ( A ) singularities. Figure 10 : The merger/splitting of a cusp ( A ) curve in a A − point. The upper threepanels show the unfolding of a A − singularity in Lagrangian space. The lower three panelsshow the corresponding unfolding in Eulerian space. The two panels on the left show a fold( A ) sheet, cusp ( A ) curves and swallowtail ( A ) singularities. The middle panels show themerger/splitting of the cusp ( A ) curves in a A − point. The right panels show the resultingfold ( A ) sheet and cusp ( A ) curves singularities.– 19 – igure 11 : The creation/annihilation of swallowtail singularities in a butterfly ( A ) singu-larity. The upper three panels show the unfolding of a A singularity in Lagrangian space.The lower three panels show the corresponding unfolding in Eulerian space. The two panelson the left show a fold ( A ) sheet, and cusp ( A ) curve. The middle panels show the cre-ation/annihilation of the butterfly ( A ) singularity on the cusp ( A ) curve. The right panelsshow the resulting fold ( A ) sheet, cusp ( A ) curve and swallowtail ( A ) singularities.points of the same field µ ti,i . The latter mark the A i − singularities, A i +4 = { q ∈ L | q ∈ A i ( t ) , µ ti,i ( q ) max-/minimum of µ ti,i | A i ( t ) for some t } ,A i − = { q ∈ L | q ∈ A i ( t ) saddle point of µ ti,i | A i ( t ) for some t } . (4.16)The A ± critical points are constrained to lie on the curvelike variety A i ( t ) . Their identityis therefore determined by the interplay between the geometric properties of two entities. Oneof these is the geometry of the field µ ti,i , the other that of the geometry of the curvelike variety A i ( t ) . For illustrations of the A +4 and A − singularities we refer to figure 9 and 10.From the caustic conditions – as expressed in eqn. (4.12) – we may also immediatelyobserve that the A i ± points belong to the A i variety. In fact, this also represents a conditionon the topology of the field µ ti,i and that of the A i ( t ) variety. A caustics Finally, also the swallowtail curves A i curve get folded in Eulerian space. It leads to theemergence of so-called butterfly singularities, or A singularities. Following the same reason-ing as for the A i and A i varieties, we may infer from the general shell-crossing condition thatthe A i curve gets folded in the points A i . In general this happens when there exists a tangentvector of A parallel to v i , i.e. A i = { q ∈ L | q ∈ A i and v i ∈ T q A i } (4.17)– 20 –n the case three-dimensional case, when the displacement field s t ( q ) is separable into temporaland spatial parts, time evolution can be seen as a progression through a series of surfaces.The folding points can then be found from the relation, A i = { q ∈ L | q ∈ A i ( t ) and µ ti,iii = 0 for some time t } . (4.18)Figure 11 shows an illustration of a A singularity.The butterfly singularity is the highest dimensional singularity that may surface in three-dimensional Lagrangian fluids. It is important to realize that the butterfly singularity onlyexists at one point in space-time. D family The D family of caustics correspond to manifolds for which the condition µ i = 0 , (4.19)holds for two eigenvalue fields simultaneously. From this, we may immediately infer that thesecaustics form at the intersection of two A ( t ) fold sheets, the A i ( t ) and A j ( t ) varieties. Inall, for three-dimensional fluids three classes of D caustics can be identified, the D − elliptic,the D +4 hyperbolic and the D parabolic umbilic caustic. D caustics The D caustics are defined by the points q in Lagrangian space, at which two of the eigen-values have the same value. For instance, the D ij ( t ) caustic, with i (cid:54) = j , is outlined by thepoints q for which at the time t the eigenvalues µ i (t) and µ j ( t ) are equal, µ ti = µ tj . While theeigenvalue µ ti defines the fold sheet A i , and the eigenvalue µ ti the fold sheet A j , the umbilic D ij caustic consist of the set of points q for which D ij ( t ) = { q ∈ L | q ∈ A i ( t ) ∩ A j ( t ) } . (4.20)In three-dimensional space, one would expect that the intersection of the two sheets A i ( t ) and A j ( t ) to consist of one-dimensional curves. This would certainly be true for two sheetsthat would be entirely independent of each other. However, the situation at hand concerns ahighly constrained situation, in which the two eigenvalues µ i and µ j are strongly correlated.Because of the latter, the intersection between the folds A i and A j is considerably morecomplex. Instead of a continuous curve, the intersection consists of isolated, singular points.A telling illustration – and discussion – of this, for the two-dimensional situation, can befound in [39]. D singularities and A varieties To investigate the geometry and structure of the set D ij ( t ) we focus on the particular situationof the set D ( t ) , in which the two first eigenvalues µ and µ have the same value, µ t = µ t = − . Without loss of generality, we transform the coordinate system such that the thirdeigenvector v defines the q axis. This transformation makes the q q -plane the one in whichwe see the folding and collapse of the phase space sheets to the A and A caustics.Assuming that the deformation tensor M is diagonalizable, in this coordinate system ithas the form, M = M M M M
00 0 µ , (4.21)– 21 –n which µ is the third eigenvalue of M . Because the eigenvalues are equal, we get thefollowing 2 conditions for the D caustic. M ( q ) = M ( q ) ,M ( q ) = 0 . (4.22)Hence, the deformation tensor is M = µ µ
00 0 µ . (4.23)As a consequence of the inferred constraints (4.22) for the D singularities is that D ij pointswill always be located on the two corresponding A varieties, A i and A j . We may infer thisfrom the following observation. In the coordinate system introduced above (cf. eq. (4.21)),the eigenvector for the third eigenvalue µ is given by v = (0 , , . The eigenvectors v and v both lie in the q q -plane, and since the matrix upper × matrix is degenerate we havethe freedom to take them to be orthogonal to the gradient of the corresponding eigenvaluefields which will also lay in the q q -plane. This means that v · ∇ µ = µ , = 0 ,v · ∇ µ = µ , = 0 . (4.24)This proves the unfolding D ij → A i and D ij → A j . For the relations between the singularityclasses see section 7.1. For a formal proof see [39]. For the case of a non-diagonalizabledeformation tensor we note that a small perturbation in the initial condition genericallymakes the deformation tensor diagonalizable. The D +4 and D − caustics Shell-crossing for A caustics is a one-dimensional process. A direct implication of this is thatthe related critical points are equivalent up to diffeomorphisms. For the D family this is nolonger true. Shell-crossing for the D -family is two dimensional. As a consequence, the D class consist of hyperbolic ( D +4 ) and elliptic ( D − ) umbilic points, i.e. D ij ( t ) = D + ij ( t ) ∪ D − ij ( t ) . (4.25)In order to infer the corresponding caustic conditions we consider the two constraint quantities Q ( q ) and Q ( q ) (see eq. (4.22)), Q ( q ) = M ( q ) − M ( q )2 ,Q ( q ) = M ( q ) , (4.26)which at the D singularity location vanish, i.e. Q ( q s ) = 0 and Q ( q s ) = 0 . By a Taylorexpansion of Q ( q ) and Q ( q ) in a neighbourhood around the D singularity, we find that forpoints located in the q q -plane, Q ( q ) = a q + b q ,Q ( q ) = c q + d q . (4.27)– 22 – igure 12 : The hyperbolic/elliptic umbilic ( D ± ) singularities. The upper two panels showthe elliptic umbilic ( D +4 ) singularity. The lower panels show the hyperbolic umbilic ( D − )singularity. The two panels on the left are their representations in Lagrangian space andthe two panels on the right their representation in Eulerian space. The black sheets arefold ( A ) sheets corresponding to one eigenvalue field. The green sheets are fold ( A ) sheetscorresponding to a second eigenvalue field. The red lines are cusp ( A ) curves. The point in thecenter depict the hyperbolic/elliptic umbilic ( D ± ) singularities. The hyperbolic umbilic ( D +4 )and elliptic umbilic ( D − ) singularity are also known as the purse and pyramid singularity.In this expansion, we have taken the D singularity to define the origin of the coordinatesystem. The parameters a , b , c and d are the derivatives of Q ( q ) and Q ( q ) at the D location, a = 12 ∂ ( M − M ) ∂q , b = 12 ∂ ( M − M ) ∂q , c = ∂M ∂q , d = ∂M ∂q . (4.28)As proposed by [28], the determinant S M of the corresponding Q Q map, S M = bc − ad = 12 [( M − M ) M − ( M − M ) M ] , (4.29)– 23 –s invariant under rotations in the q q -plane . In the expression above, we have used thenotation M iik = ∂M ii ∂q k , M ikk = ∂M ik ∂q k . (4.30)Using the relations between the matrix elements M , M and M and the eigenvalues µ and µ , we may recast the determinant S M in an explicit expression incorporating theseeigenvalues, S M = 12 [( µ − µ ) , µ , − ( µ − µ ) , µ , ] . (4.31)As [28] pointed out, the transformation can be shown to consist of two branches. Theiridentification surfaces via a rescaling of the determinant via the multiplication by a positive number. We then find that the two branches correspond to two separate singularity classesof the D family, D ± ij ( t ) = { q ∈ L | q ∈ A i ( t ) ∩ A j ( t ) ∧ sign ( S M ) = ± } , (4.32)where the points q ∈ A i ( t ) ∩ A j ( t ) are the points for whom at time t the caustic conditionsare simultaneously valid for two eigenvalues, i.e. µ i = 1 + µ j = 0 . Integrated over time,these D ± ij ( t ) points trace out the curves D ± ij , D ± ij = { q ∈ L | q ∈ A i ( t ) ∩ A j ( t ) ∧ sign ( S M ) = ± , for some time t } . (4.33)For an illustration of the hyperbolic/elliptic umbilic ( D ± ) caustic see figure 12. D ± points The topology of the D ± ij ( t ) variety changes at D ± and D points. The D ± points areanalogous to the A ± points of the A -family. The D ± points occur when i th and j th eigenvaluefield, µ i and µ j , restricted to the points q in the D ± ij variety reaches a minimum or maximum,i.e. D ij +4 = { q ∈ L | q ∈ D + ij ( t ) ∧ µ tk ( q ) max-/min. of µ tk | D + ij ( k = i or k = j ) for some t } D ij − = { q ∈ L | q ∈ D − ij ( t ) ∧ µ tk ( q ) max-/min. of µ tk | D − ij ( k = i or k = j ) for some t } (4.34)Particularly interesting is the fact that the D ± points are always created as a pair. Two D +4 points are created simultaneously, as are D − points. By implication, also the D ± curves(eq. (4.33)) are always created in pairs. This is in contrast to the D points, which go alongwith the creation of a pair consisting of a D +4 and a D − point. D caustics The shell-crossing condition applied to the D ij variety yields the caustic conditions for the D parabolic umbilic singularity. The manifold D ij forms a singularity in the point q s ∈ D ij ( t ) ifand only if the tangent vector T ∈ T q s D ij is normal to v ∗ k , with k (cid:54) = i, j . Hence, the tangentvector T ∈ span C { v i , v k } , i.e. D ij = { q ∈ L | q ∈ D ij and span { v i , v j } ∩ T q D ij (cid:54) = ∅} (4.35) In fact, it can be shown that this determinant is a third-order invariant under rotations [28]. – 24 –or three dimensional fluids in which the deformation tensor is separable in a time factor anda spatial factor, the normal n = ∇ ( µ ti − µ tj ) , is orthogonal to both v i and v j , ( µ i − µ j ) ,i ≡ v i · ∇ ( µ ti − µ tj ) = 0 , ( µ i − µ j ) ,j ≡ v j · ∇ ( µ ti − µ tj ) = 0 . (4.36)The collection of all such points form the variety D ij = { q ∈ L | q ∈ D ij ( t ) ∧ ( µ i − µ j ) ,i = ( µ i − µ j ) ,j = 0 for some time t } . (4.37)The D ij lays on the A i and A j variety. The elliptic and hyperbolic umbilic ( D ± ) pointsmerge in parabolic umbilic ( D ) points, since D ij ( t ) ⊂ D ij ( t ) and S M = 12 { ( µ i − µ j ) ,j µ i,j − ( µ i − µ j ) ,i µ j,i } = 0 . (4.38)The D points are stable singularities in the classification of Lagrangian singularities. Forgeneral dynamics they are unstable and not included in the classification scheme. For a visual appreciation of the process leading to the formation of the various classes ofcaustics identified in the subsections above, it is helpful to consider the phase-space manifoldon which all mass elements are located in 6-D phase-space L × E . This is called the phase-space sheet [see e.g. 1, 62]. The dynamical evolution of a system leads to the folding of thisphase-space sheet. In a sense, we can recognize a hierarchical process in which the phase-spacesheet is wrapped into an increasingly complex pattern. In this process we see the emergenceof a hierarchy of complex spatial folds.The phase-space sheet folding process generates higher order singularities within the A caustic itself. These can only be identified with the help of the complementary eigenvectorconditions. Restricting the manifold M to the points q s located in the A caustic, one mayidentify the subset of points for whom a nonzero vector T exists that (a) is tangent to the A manifold and (b) is orthogonal to the span of dual eigenvectors Span { v ∗ j | j (cid:54) = i } . This subsetfulfils the shell-crossing conditions and maps into a higher order singularity. Proceeding alongthe sequence of caustic conditions leads to the identification of the entire hierarchy of caustics.The classification of A family caustic involves one eigenvalue for which µ i = 0 . It isstraightforward to see that a similar procedure follows for configurations involving more thanone eigenvalue for which µ k = − . For example, if both µ = 0 and µ = 0 , then T will be a vector orthogonal to the dual eigenvector v ∗ . The eigenvalue conditions thereforetrace a line through three-dimensional Lagrangian space. The points q along this line aresingularity points. Along this line we subsequently seek to identify higher-order singularities,by identifying points q s along the line for which a tangent vector T exists fulfilling the shell-crossing conditions.Conversely, note that if µ i (cid:54) = 0 for all i , then there does not exist any T satisfyingthe general shell-crossing condition. – 25 – Figure 13 : Eigenvalue field and singularity points. In the case of the two-dimensionalZel’dovich approximation (see appendix A). The Zel’dovich approximation concerns the spe-cific situation of potential flow, for which the eigenvalues and eigenvectors are real, and showsthe field of the lowest eigenvalue. Top lefthand frame: the contour map illustrates the typicalstructure of the eigenvalue field λ corresponding to a 2-D Gaussian random density field.Indicated are the positions of different A -family singularity points and varieties. The run ofthe A line is particular noteworthy. One may appreciate how the identity of the varioussingularities is determined by the specific geometric character of the eigenvalue field λ ( q ) , asexpressed in its derivatives. Bottom lefthand frame: the panel depicts the run of the eigen-value field along the A curve (in the contour map of top lefthand frame). Note the location ofthe A ± points and A points on the extrema of the curve. The green curve represents the level b + ( t ) − indicating which parts of the A line has formed at three instances depicted in therighthand panels. Righthand panels: the three panels show the evolution, in Eulerian space,of the A line. Note the appearance of the corresponding caustics and the relation betweenthe geometry of the A line in Eulerian space and the A ± and A points corresponding to thethree green lines in the lower lefthand panel. This is a more extensive version of figure 8 inHidding et al. 2014. – 26 – .4 Spatial Connectivity: Singularities and Eigenvalue Fields With the purpose to provide a guide that evokes a visual intuition for the connection be-tween the structure and geometry of the eigenvalue fields and the formation of the varioussingularities, in particular those of the A -family, we include figure 13. It shows a contourmap representing the typical structure of the eigenvalue field µ i . This field corresponds to atwo-dimensional Gaussian random density field. For reasons of convenience, we have assumedhigher eigenvalues to correspond to earlier collapse, and negative ones to no collapse (in otherwords, we have mirrored µ i ). The geometry and topology of the eigenvalue landscape isdecisive for the occurrence of singularities. This may already be inferred from the positionsof different A -family singularity points and varieties, whose positions are indicated on thecontour map.The landscape defined by the eigenvalue contours is varied, characterized by severalpeaks, connected by ridges with lower µ i values. These, in turn, are connected to valleysin which µ i attains negative values that will prevent collapse – along the direction of theeigenvector v i – of the corresponding mass elements at any time. From the density relation(eqn. (2.6)), we know that the region of space that has undergone collapse before the currentepoch (i.e. attained an infinite density) is the superlevel set of the eigenvalue field defined bythe current value µ ti . For each time t , the positive value contours correspond to the A ( t ) fold sheets. Collapse occurs first at the maxima in the field. These mark the birth of newfeatures, and are designated by the label of A +3 points. Evidently, the steepness of the hillaround these maxima, i.e. the gradient ∇ µ i ( q ) , will determine how and which mass elementsaround the hill will follow in outlining the emerging feature around the A +3 points.The run of the A line is particularly noteworthy. The key significance of the A curveis evident from the observation that all A -family singularities are aligned along the ridge.In two-dimensional space, the A curves delineate the points where the eigenvalues µ i aremaximal along the direction of the corresponding local eigenvector. At these points, alongthe eigenvector direction, the gradient of the eigenvalues is zero, i.e. they are the points wherethe eigenvector v i is perpendicular to the local gradient of ∇ µ i of the eigenvalue field. Below,in section 4.1.3, we will see that this follows directly from the shell-crossing conditions thatwere derived in the previous section. Because of this there is a line-up and accumulation ofneighbouring mass elements that simultaneously pass through the singularity. When mappedto Eulerian space, this evokes the formation of an A cusp.To illustrate the connection between A curves and the various singularities even morestrongly, the bottom lefthand panel depicts the run of the eigenvalue field along the A curve.In particular noteworthy is the location of the A ± points and A points on the extrema ofthe curve. A prominent aspect of this is the presence of the A − points at saddle junctions inthe eigenvalue field. These are topologically the most interesting locations, as they evoke themerging of separate fold sheets into a single structure. In other words, they are the pointswhere the topological structure of the field undergoes a transition and where the connectivityof the emerging structural features is established. To establish this even more strongly, thethree righthand panels of figure 13 represent a time sequence of the evolving structure alongthe A line as it is mapped to its appearance in Eulerian space. The evolution follows thelinear Lagrangian Zel’dovich approximation (see [72] and appendix A). We may note theappearance and merging of the corresponding caustics.– 27 – Classification of singularities
The form and morphology in which the various singularities that were inventorized in theprevious section will appear in the reality of a physical system depends on several aspects.The principal influence concern the dynamics of the system, as well as its dimensionality. Thedynamics determines the way the fluid evolves, to a large extent via its dominant influence onthe accompanying flow of the fluid. This affects the morphology of the fluid, and in particularthe occurrence of singularities. Evidently, also the dimensionality of the fluid process will bearstrongly on the occurrence and appearance of singularities. Higher spatial dimensions mayenlarge the number of ways in which a singularity may form. It also influences the ways inwhich singularities can dynamically transform into one another.In this section, we provide an impression of the variety in appearance of singularities.To this end, we will first discuss the generic singularity classification scheme that we follow.It is not the intention of this study to provide an extensive listing of all possible classes offluids. Instead, to make clear in how different physical situations may affect the appearanceof singularities, we restrict our presentation of classification schemes to two different classesof fluids. We also restrict our inventory to fluids in a three- dimensional context. It isthe most representative situation, and at the same time offers a good illustration of otherconfigurations.
To appreciate the role of the dynamics in constraining the evolution and appearance of afluid, and that of the formation and fate of the singularities in the fluid, it is important tounderstand and describe its evolution in terms of six-dimensional phase space.One way of defining phase space C is in terms of the Cartesian product of Lagrangianand Eulerian manifolds L and E , i.e. C = L × E . In this context, the phase space coordinatesof a mass element are ( q, x ) . Every point in phase space ( q, x ) ∈ C represents the initial andfinal position q and x of a mass element at some time t . Evidently, one may also opt for themore conventional definition consisting of space coordinates x and canonical momenta p , inwhich case the phase space coordinate of a mass element are given by ( x, p ) . However, for thedescription of Lagrangian fluid dynamics it is more convenient to follow the first convention.We should note that for this description of phase space Liouville’s theorem does not apply,specifically not for the Euclidean notion of volumes.At the initial time t = 0 , the Lagrangian map is the identity map, i.e. for all q ∈ Lx ( q ) = q . In phase space C , the fluid then occupies the submanifold L = { ( q, x ( q )) ∈C| q ∈ L } . If we equip C with a symplectic structure ω , we can prove this to be a so-calledLagrangian submanifold (for a precise definition of Lagrangian submanifolds see appendix D).Differences in the dynamics of a fluid reveal themselves in particular through majordifferences in the phase space structure and topology of the manifolds delineated by the masselements. To provide an impression of the differences in morphology and classification ofsingularities emerging in fluids of a different nature, specifically that of fluids with a differentdynamical behaviour, we concentrate the discussion on two different classes of Lagrangianfluids:1. Generic Lagrangian fluids .Lagrangian fluids for which the map x t : L → E is a generic continuous and differentiable– 28 –apping from L to E for every time t . The dynamics does not restrict the map x to anyextent. We describe the classification up to local diffeomorphisms, i.e. two singularitiesare considered equivalent if and only if there exist local coordinate transformations,which map them into each other.2. Lagrangian fluids with Hamiltonian dynamics .The evolution of the fluid is governed by a Hamiltonian. This assumption restrictsthe possible evolution of the fluid. Formally, the map x corresponds uniquely to a so-called Lagrangian map. The singularities of Lagrangian maps, known as Lagrangiansingularities, are classified up to Lagrange equivalence.Lagrangian fluids with Hamiltonian dynamics form an important class of fluids: fundamentaltheories of particle physics generally allow for a Hamiltonian description. Nonetheless, in arange of practical circumstances we may encounter fluids that are either more or less con-strained. An example are fluids with effective dynamics. They contain friction terms whichare not described by Hamiltonian systems. Such fluid systems are less restrictive than thosethat are specifically Hamiltonian. On the other hand, there are also Hamiltonian fluids thatare characterized by additional constraints. For the classification of singularities of generic one-family maps x : L → E , with L and E three-dimensional, we follow the classification by Bruce [17]. Bruce showed that the singular-ities that emerge in generic mappings are equivalent to those emerging in the simple linearmaps x t ( q ) = q + t u ( q ) , (5.1)in which u is a vector field on L . In general, the vector field u ( q ) consists of both a longitudinaland a transversal part, u ( q ) = u l ( q ) + u t ( q ) . (5.2)The longitudinal component corresponds to potential motion and has curl zero, ∇ × u l = 0 ,while the transversal component has divergence zero, ∇ · u t = 0 .The classification of singularities in general Lagrangian fluid dynamics is expressed bytheorem 3. We restrict ourselves to listing the classification scheme, in terms of the genericexpressions for the maps x t ( q ) of each of the classified singularities. In appendix E we showthat these normal forms indeed satisfy the corresponding caustic conditions. Note that theclassification was derived using the classification of jet-spaces. It successfully cauterized theproperties of caustics appearing in Lagrangian maps but did not provide a practical way todetect them in realizations. Theorem: 3
A stable singularity occurring in a Lagrangian fluid with generic dynamics is,up to local diffeomorphisms, equivalent to one of the following classes: – 29 – ingularity Singularityclass Map x t ( q ) name A x t ( q ) = q trivial case A x t ( q ) = q + t (cid:0) , , q − q (cid:1) fold A x t ( q ) = q + t (cid:0) , , q q + q − q (cid:1) cusp A x t ( q ) = q + t (cid:0) , , q q + q − q (cid:1) swallowtail A x t ( q ) = q + t (cid:0) , , q q + q q + q − q (cid:1) butterfly D ± x t ( q ) = q + t (cid:0) , q q − q , q ± q + q q − q (cid:1) hyperbolic/elliptic A ± x t ( q ) = q + t (cid:0) , , ( q ± q ) q + q − q (cid:1) A ± x t ( q ) = q + t (cid:0) , , q q ± q q + q − q (cid:1) Note: The normal forms x t ( q ) form the singularity at the origin q = 0 , at t = 1 . Thefirst five singularity classes are the A -family. The subsequent class is the D -family. Thelast two are the normal forms of the A and A points. The A k class has co-rank andco-dimension k − . The D ± singularities have co-rank and are one-dimensional. The evolution of Lagrangian fluids with Hamiltonian dynamics is more constrained than thatof generic Lagrangian fluids. As the fluid develops complex multistream regions, the phasespace submanifold L t = { ( q, x t ( q )) | q ∈ L } for fluids with Hamiltonian dynamics remains aLagrangian submanifold.A key step in evaluating the emerging singularities is that of connecting the displacementmap s t ( q ) to the Lagrangian map. In appendix D.2, we describe in some detail how a givenLagrangian map can be constructed from a Lagrangian submanifold L . A Lagrangian mapcan develop regions in which multiple points in the Lagrangian manifold are mapped to thesame point in the base space.Lagrangian singularities are those points at which the number of pre-images of the La-grangian map undergoes a change. Lagrangian catastrophe theory [5, 13] classifies the stablesingularities. This refers to the stability of singularities with respect to small deformations ofthe Lagrangian manifold of L . This is true up to Lagrangian equivalence, a concept that is ageneralization of equivalence up to coordinate transformation. For a more formal and precisedefinition of Lagrangian equivalence see appendix D.It can be demonstrated [see 13] that every Lagrangian map l : L → C → E is locallyLagrangian equivalent to a so-called gradient map, i.e. the map x t is locally equivalent to x t ( q ) = ∇ q S t , (5.3)for some S t : L → R . By recasting S t in terms of a function Ψ t : L → R , S t = 12 q + Ψ t ( q ) , (5.4)we find that locally the map x can be written in the form x t ( q ) = q + ∇ q Ψ t ( q ) . (5.5)Evidently, this implies that the displacement map is longitudinal, and that the correspondingJacobian ∂s t /∂q is symmetric. – 30 –he classification of singularities of a Lagrangian fluid with Hamiltonian dynamics isexpressed by theorem 4. In appendix E it is shown that these normal forms indeed satisfy thecorresponding caustic conditions. For proofs we refer to Arnol’d [5]. Note that the classifica-tion was derived using the classification of critical points of scalar functions and the theoryof generating functions. It successfully characterized the properties of caustics appearing inLagrangian maps but did not provide a practical way to detect them in realizations. Theorem: 4
A stable Lagrangian singularity of a Lagrangian fluid with Hamiltonian dynam-ics, is locally Lagrange equivalent to one of the following classes:Singularity Singularityclass Map x t ( q ) name A x t ( q ) = q trivial case A x t ( q ) = q + t (cid:0) , , q − q (cid:1) fold A x t ( q ) = q + t (cid:0) q , , q ( q − (cid:1) cusp A x t ( q ) = q + t (cid:0) q , , q q + q − q (cid:1) swallowtail A x t ( q ) = q + t (cid:0) q , q , q q + q q + q − q (cid:1) butterfly D ± x t ( q ) = q + t (cid:0) ± q q − q , hyperbolic/elliptic ± (cid:0) q + q (cid:1) + 2 q q + 2 q − q , q (cid:1) D x t ( q ) = q + t (cid:0) , q − q , q − q (cid:1) parabolic A ± x t ( q ) = q + t (cid:0) q q , ± q q , ( q ± q ) q + q − q (cid:1) A ± x t ( q ) = q + t (cid:0) q , ± q q , q q ± q q + q − q (cid:1) Note: The normal forms x t ( q ) form the singularity at the origin q = 0 , at t = 1 . Thefirst five singularity classes are the A -family. The subsequent two are the D -family. Thelast two are the normal forms of the A and A points. The A k class has co-rank andco-dimension k − . The D k singularities have co-rank and co-dimension k − [5]. Comparing the classification schemes for generic Lagrangian singularities and those forLagrangian fluids with Hamiltonian dynamics, we may note the similarities. Both classifica-tions have an A and a D family. It can be demonstrated that the A singularity classes ofthe scheme for Lagrangian fluids with Hamiltonian dynamics are contained in those corre-sponding to the generic Lagrangian fluid. Concretely, this means that a displacement fieldcorresponding to the Hamiltonian A k class is also an element of the generic A k class.The D families are some what different. The Hamiltonian D class is contained inthe generic D class. However, the Hamiltonian D class has no analogue in the genericclassification scheme. This is a result of the D singularity not being stable under coordinatetransformations.A final remark concerns the singularity classification schemes for higher dimensionalfluids. For these a more elaborate classification scheme applies. This classification scheme isdescribed in appendix D. Singularities generally change their class upon small, but finite, deformations of the displace-ment map s t . The corresponding evolution of a singularity follows the universal unfoldingprocess of singularities. The general behavior is described in the following unfolding diagram,in which the arrows indicate the singularity into which specific singularities can transform.– 31 – A A A A D D For i ≥ , the A i singularities decay into A i − singularities. For i ≥ , the D i singularitiesdecay into either A i − or D i − singularities. In section 7 we will describe how the decay ofsingularities is connected to the evolution of the large-scale structure in the Universe and inoutlining the spine of the cosmic web. The process of formation and evolution of structure in the Universe is driven by the gravita-tional growth of tiny primordial density and velocity perturbations. When it reaches a stageat which the matter distribution starts to develop nonlinearities, we see the the emergenceof complex structural patterns. In the current universe we see this happening at Mega-parsec scales. On these scales, cosmic structure displays a marked intricate weblike pattern.Prominent elongated filamentary features define a pervasive network. Forming the denseboundaries around large tenuous sheetlike membranes, the filaments connect up at massive,compact clusters located at the nodes of the network and surround vast, underdense andnear-empty voids.The gravitational structure formation process is marked by vast migration streams,known as cosmic flows. Inhomogeneities in the gravitational force field lead to the displace-ment of mass out of the lower density areas towards higher density regions. Complex struc-tures arise at the locations where different mass streams meet up. Gravitational collapse setsin as this happens. In terms of six-dimensional phase space, it corresponds to the local foldingof the phase space sheet along which matter – in particular the gravitationally dominant darkmatter component – has distributed itself.
The positions where streams of the dark matter fluid cross are the sites where gravitationalcollapse occurs. The various types of caustics described and classified in our study mark thedifferent configurations in which this process may take place. Their locations trace out a La-grangian skeleton of the emerging cosmic web, marking key structural elements and establish-ing their connectivity (also see the discussion in [39]). In other words, the A , A , A , D , D varieties, in combination with the corresponding A ± , A ± , and D ± points, are the dynamicalelements whose connectivity defines the weaving of the the cosmic web [4, 15, 23, 66, 72].On the basis of this observation, we may obtain the skeleton of the cosmic web by mappingthe caustic varieties defined above to Eulerian space with the Lagrangian map x t . Followingthe identification of the various caustic varieties and caustic points in Lagrangian space, theapplication of the map x t will produce the corresponding weblike structure in Eulerian space.Of central significance in our analysis and description of the cosmic web is the essentialrole of the deformation tensor eigenvector fields in outlining the caustic skeleton and in es-tablishing the spatial connections between the various structural features. So far, Lagrangianstudies of the cosmic web have usually been based on the role of the eigenvalues of the de-formation tensor (for recent work see [24, 49, 70]). Nearly without exception, they ignore theinformation content of the eigenvectors of the deformation tensor. In this work we actuallyemphasize that the eigenvectors are of key importance in tracing the spatial locations of the– 32 – (cid:0) D + -15 -10 -5 0 5 10 15 0.1 1 10 100 ✁ + 1 – 33 – igure 14 : Spatial distribution of singularities in the Lagrangian and Eulerian cosmic web.The figure compares the spine of the cosmic web with the mass distribution in a 2-D N -body simulation. Left panel: initial field of density fluctuations and the skeleton of identifiedsingularities/catastrophes. Right panel: density field of an evolved 2D cosmological N -bodysimulation, in which the Lagrangian skeleton of singularities is mapped by means of theZel’dovich approximation. From Feldbrugge et al. [33].different types of emerging caustic features and, in particular, in establishing their mutualspatial connectivity. This important fact finds its expression in terms of the caustic conditions that we have derived in this study.The study by Hidding et al. [39] illustrated the important role of the deformation fieldeigenvectors in outlining the skeleton of the cosmic web, for the specific situation of A cusplines in the 2-D matter distribution evolving out of a Gaussian initial density field. Thepresent study describes the full generalization for the evolving matter distribution (a) foreach class of emerging caustics in (b) in spaces of arbitrary dimension D . A telling and informative illustration of the intimate relationship between the caustic skeletondefined by the derived caustic conditions and the evolving matter distribution is that offeredby the typical patterns emerging in the two-dimensional situation. Figure 14 provides adirect and quantitative comparison between the caustic skeleton of the cosmic web and thefully nonlinear mass distribution in an N-body simulation. The three panels in the lefthandcolumn show the Lagrangian skeleton for a two-dimensional fluid. The fluid is taken toevolve according to the Zel’dovich approximation [72] (see appendix A), which representsa surprisingly accurate first-order Lagrangian approximation of a gravitationally evolvingmatter distribution [see e.g. 64]. The initial density field of the displayed models is that ofa Gaussian random density field [2, 14], which according to the latest observations and tocurrent theoretical understanding is an accurate description of the observed primordial matterdistribution [27, 45, 55].To enable our understanding of the hierarchical process of structure formation and theresulting multiscale structure of the cosmic web, we assess the caustic structure of the La-grangian matter field at three different resolutions. In figure 14 the field resolution decreasesfrom the top panels to the bottom panels, as the initial density field was smoothed by anincreasingly large Gaussian filter. The contour maps that form the background in these pan-els represent the resulting initial density fields. The red lines trace the A variety, i.e. the A lines, for the largest eigenvalue µ field (also see fig. 13 to appreciate how they are re-lated). Also the A ± points and D ± points are shown, the first as red dots, the latter as blacktriangles.The resulting weblike structure in Eulerian space is depicted in the corresponding right-hand panels. The A lines, A ± points and D ± points are mapped to their Eulerian location bymeans of the Zel’dovich approximation. The red lines, red dots and black triangles representthe Eulerian skeleton corresponding to the Zel’dovich approximation. These are superim-posed on the density field of the corresponding N-body simulations. The comparison betweenthe latter and the Eulerian skeleton reveal that the caustic skeleton – the assembly of A lines, A ± points and D ± points – trace the principal elements and connections of the cosmicweb seen in the N-body simulations remarkably well (see table 1 for the identification of the– 34 – igure 15 : The log density field of a dark matter N -body simulation with Λ CDM cosmol-ogy in a box of h − Mpc with particles and elements of the caustic skeleton of theZel’dovich approximation [51]. Top right panel: the cusp ( A ) sheets (dark blue), the swal-lowtail ( A ) lines (light blue) and the elliptic/hyperbolic umbilic lines (yellow) correspondingto the lowest eigenvalue field of the caustic skeleton. Note that the Zel’dovich approximationconcerns a potential flow, which means that the eigenvalue fields can be ordered. The initialdensity field was smoothed on the scale . h − Mpc. Bottom left panel: the swallowtail ( A )lines (light blue) and the elliptic/hyperbolic umbilic lines (yellow) corresponding to the lowesteigenvalue field of the caustic skeleton. The initial density field is smoothed at . h − Mpc.– 35 –ines and points to the cosmic web). Moreover, by assessing the caustic structure at differentresolutions of the density field, one obtains considerable insight into the multiscale structureand topology of the cosmic web.
One of the unique features facilitated by the caustic conditions that we have derived in theprevious sections is the ability to go beyond the two-dimensional case and construct andexplore the full caustic skeleton of the three-dimensional mass distribution. In the case ofthe skeleton of the cosmic web defined by the three-dimensional mass distribution, the cusp( A ) sheets correspond to the walls or membranes of the large scale structure [15, 23, 51, 66].The swallowtail ( A ) and elliptic/hyperbolic umbilic ( D ± ) lines correspond to the filamentsof the cosmic web and the butterfly ( A ) and parabolic umbilic ( D ) points correspond tothe cluster nodes of the network [3, 15, 23, 51, 66]. The identification of the caustics in thethree dimensional cosmic web is summarized in table 1.To appreciate the impressive level at which the caustic skeleton is outlining the three-dimensional weblike mass distribution, figure 15 provides an instructive illustration. Thefigure depicts elements of the caustic skeleton of the Zel’dovich approximation in a h − Mpc box. The resulting skeleton is superposed on the log density field of a dark matter N -body simulation in a Λ CDM cosmology with particles [51]. We should emphasize thatthe Zel’dovich approximation is linear and that the corresponding skeleton is completely localin the initial conditions. While a full and detailed analysis of these three-dimensional weblikepatterns is the subject of an upcoming accompanying paper [32], the illustrations of figure 15already give a nice impression of the ability of the caustic conditions to outline the spine ofthe cosmic web.The top righthand panel contains the cusp ( A ) sheet (dark blue colour) and the swal-lowtail ( A ) and elliptic/hyperbolic umbilic ( D ± ) lines (light blue colour) corresponding tothe lowest eigenvalue field, superimposed on the density field of the N -body simulation (redshaded log density field values). The pattern concerns the caustics obtained for a displacementfield that is filtered at a length scale of . h − Mpc. Close inspection reveals the close corre-spondence between the cusp sheets of the caustic skeleton and the flattened - two-dimensional- features in the mass distribution of the cosmic web. Notwithstandig this, one may also ob-serve that the two-dimensional skeleton does not capture all the structures present in the N -body simulation. This is predominantly an issue of scale, as the corresponding displace-ment field cannot resolve and trace features whose size is more refined than the . h − Mpcfilter scale.An impression of the more refined structure can be obtained from the bottom left panelof figure 15, which follows the line-like elements of the caustic skeleton at a length scale of . h − Mpc. More specifically, it shows the swallowtail ( A ) and elliptic/hyperbolic umbilic( D ± ) lines of the caustic skeleton. The correspondence of these with the prominent andintricate filamentary pattern in the cosmic mass distribution is even more outstanding thanthat of the A sheets with the membranes in the density field. It is important to realize, andemphasize, that blue curves were generated using only the eigenvalue field corresponding tothe first collapse. This already creates a filament in the network of caustics, without the needto involve the second eigenvalue. In other words, collapse along the second eigenvector is notnecessary to create a filament-like structure (also see [39]). This leads to a radical new insighton structure formation, in that it suggests the different possible late-time morphologies for– 36 –ingularity Singularity Feature in the Feature in theclass name 2D cosmic web 3D cosmic web A fold collapsed region collapsed region A cusp filament wall or membrane A swallowtail cluster or knot filament A butterfly not stable cluster or knot D hyperbolic/elliptic cluster or knot filament D parabolic not stable cluster or knot Table 1 : The identification of the different caustics in the - and -dimensional cosmic webfilaments [40]. We may even relate this to the prominence of the corresponding filamentaryfeatures: as they concern features that have experienced collapse along two directions, theumbilic D ± filaments will have a higher density and contrast than the filigree of more tenuous A ± filaments. An additional observation of considerable interest is that the line-like A and D ± features trace the connectivity of the cosmic web in meticulous detail. Also of decisive interest in their embedding in the cosmic web, is the expected mass distribu-tion in and around the various classes of caustics.Vesilev [67] inferred the density profiles of the various classes of singularities, in casethey emerge as a result of potential motion in a collision-less self-gravitating medium. Foreach of the mass concentrations in and around these singularities, he found scale free power-law profiles. The radially average profiles display the following decrease of density ρ ( r ) as afunction of radius r .Singularity Singularity Profile ρ ( r ) class name A fold ρ ( r ) ∝ r − / A cusp ρ ( r ) ∝ r − / A swallowtail ρ ( r ) ∝ r − / A butterfly ρ ( r ) ∝ r − / D hyperbolic/elliptic ρ ( r ) ∝ r − D parabolic ρ ( r ) ∝ r − log (1 /r ) With respect to these radially averaged profiles, we should realize that the mass distributionin and around the singularities is highly anisotropic. This is true for any dimension in whichwe consider the structure around the singularities.Notwithstanding this, we do observe that the steepest density profiles are those aroundthe point singularities A and D . However, they are mere transient features that will onlyexist for a single moment in time. The point singularities A and D display a less pronouncedbehaviour. However, they move over time. Also, we see that the cusp singularity A possessesa steeper mass distribution that that in and around the sheet singularity A .– 37 – .5 Higher order Lagrangian perturbations Evidently, the details of the dynamical evolution will bear a considerable influence on thedeveloping caustic structure. This not only concerns the dynamics of the system itself, but alsoits description. The examples that we presented in the previous sections showed the causticfeatures developing as the dynamics is predicated on the first-order Lagrangian approximationof the Zel’dovich formalism [72]. The visual comparison with the outcome of the corresponding N -body simulations demonstrated the substantial level of agreement. Nonetheless, given thenature of singularities, the process of caustic formation might be very sensitive to minordeviations of the mass element deformations and hence the modelling of the dynamics. Thismay even strongly affect the predicted population of caustics and their spatial organizationin the skeleton of the cosmic web. Some indications on the level to which the spatial massdistribution is influenced may be obtained from an early series of papers by Buchert andcollaborators [18–22], who were the first to explore the formation of structure in higher-orderLagrangian perturbation schemes and investigate in how far they would effect the occurrencand location of multistream regions. An important finding from their work is that 2nd ordereffects are substantial, while 3rd order ones are minimal. Elaborated and augmented byadditional work [16, 61], 2nd order Lagrangian perturbations – usually designated by thename 2LPT – have been established as key ingredients of any accurate analytical modelingof cosmic structure growth. In a follow-up to the present study, we investigate in detail therepercussions of different analytical prescriptions for the dynamical evolution of the cosmicmass distribution for the full caustic skeleton of the cosmic web.In addition to 2LPT, we will systematically investigate the caustic skeleton in the contextof the adhesion approximation [36–38, 41, 64, 68]. Representing a fully nonlinear extension ofthe Zel’dovich formalism, it includes an analytically tractable gravitational source term for thelater nonlinear stages. It accomplishes this via an artificial viscosity term that emulates theeffects of gravity, resulting in the analytically solvable Burger’s equation. With the effectiveaddition of a gravitational interaction term for the emerging structures, unlike the Zel’dovichapproximation the adhesion model is capable of following the hierarchical buildup of structureand the cosmic web [38, 40, 41]. At early epochs, the resulting matter streams coincidewith the ballistic motion of the Zel’dovich approximation. At the later stages, as the massflows approach multistream regions a solid structure is created at the shell-crossing location.Matter inside these structures is confined to stay inside, while outside collapsed structures theresults from the Zel’dovich approximation and adhesion are identical. The caustics from theZel’dovich approximation are compressed to infinitesimally thin structures, hence unifying theZel’dovich’ idea of collapsed structures in terms of shell crossing with a hierarchical formationmodel. While offering a complete model for the formation and hierarchical evolution of thecosmic web, it does accomplish this by seriously altering the flow pattern involved in thebuildup of cosmic structure. This, in turn, is expected to affect at least to some extent theproperties and evolution of the caustic population and its connectivity. In addition to characterizing the geometric and topological outline of the cosmic web in termsof the caustic skeleton, our study points to another important and related application of theformalism described. The fact that the linear Zel’dovich approximation provides such an ac-curate outline of the skeleton of the cosmic web establishes an important relation between theprimordial density and flow field and the resulting cosmic web. Via the Zel’dovich approxima-tion, we may relate the caustic skeleton directly to the statistical nature and characteristics– 38 –f the primordial density field. In other words, we may directly relate the structure of thecosmic web to the nature of the Gaussian initial density field. This, in turn, establishes adirect link between the geometric and topological properties of the cosmic web and the un-derlying cosmology. Hence a probabilistic analysis of the caustic skeleton may define a pathtowards a solidly defined foundation and procedure for using the structure of the observedcosmic web towards constraining global cosmological parameters and the cosmic structureformation process.The fact that we may invoke Gaussian statistics facilitates the calculation of a widerange of geometric and topological characteristics of the cosmic web, as they are directlyrelated to the primordial Gaussian deformation field, its eigenvalues and eigenvectors. For anexample of such a statistical treatment of -dimensional fluids, we refer to [31]. It describeshow one may not only analytically compute the distribution of maxima, or minima, but alsothe population of singularities and the length of caustic lines. In an accompaying study,we present an extensive numerical analysis of the statistics of - and -dimensional causticskeleton will follow in [32]. This will establish the reference point for the subsequent solidanalytical study of interesting geometric properties of the cosmic web (for the initial stepstowards this program see [33]).This will represent a major extension of statistical descriptions that were solely basedon the eigenvalue fields. The latter would make it possible to study the number density ofclusters and void basins, make predictions on the statistical properties of angular momentum,and even several aspects of the cosmic skeleton (e.g. [29, 56]). As we have argued extensivelyin previous sections, it is only by invoking the information contained in the correspondingeigenvector fields that we may expect to obtain a more complete census of intricate spatialproperties of the cosmic web. The caustic conditions presented in this study reveal the profound relationship between thevarious classes of singularities that may surface in Lagrangian fluids. Besides the aspect of theidentification and classification of singularities, we need to have insight in the transformationand evolution of caustics and caustic networks that accompanies the dynamical evolution ofa fluid. The evolution of the fluid, dictated by the dynamics of the system, generally involvesthe development of ever more distinctive structures and the proliferation of complex structuralpatterns.Tracing the evolution of a fluid starts at an initial time t = 0 . At that time, the dis-placement map s t is the zero map. Amongst others, this implies the fluid does not (yet)contain singularities. Starting from these near uniform initial conditions, the structure in theevolving fluid becomes increasingly pronounced. The phase space sheet that it occupies insix-dimensional space gets increasingly folded. Its projection on Euclidian space follows thisprocess, and it is as a result of the folding process that we see the fluid developing singular-ities. While the dynamical evolution proceeds to more advanced stages, we not only see theappearance of more singularities, but also the transformation of one class of singularities intoanother one. A complementary process that may underlie the changes of local geometry thatof the merging of singularities into a new singularity, itself a manifestation of the hierarchicalbuildup of structural complexity.The eigenvalue landscape in figure 13 offers an instructive tool for facilitating and guidingour understanding and visual intuition for the iterative folding of singularities in phase space– 39 –nd the accompanying caustic transformations. The dynamical evolution of a fluid goes along with a rich palet of local processes. These involvefundamental mutations in the local singularity structure that lead to significant topologicalchanges of the spatial pattern forming in the fluid. In some systems and situations this willbe a key element in the hierarchical buildup of structure.The fundamental notion in these structural mutations in the evolving fluid is that of theruling dynamics of the system evoking changes in the deformation field. Small deformationswill lead to the decay of singularities into different ones belonging to other singularity classes.Conversely, they may get folded according to a rigid order.The sequence of singularity mutations is not random and arbitrary. Due to the strictgeometric conditions and constraints corresponding to the various singularities, expressed inthe caustic conditions discussed extensively in this study, a given singularity is only allowedto transform into a restricted set of other singularities. Conversely, a given singularity mayonly have emanated from a restricted set of other singularities.In most situations a particular singularity can have decayed from only one distinctiveclass of singularities. Some may have descended from two other singularity classes. Likewise,most singularities can decay only into one distinctive other class of singularity. This is truefor all A -family singularities. D -family singularities have a richer diversity of options, withthe D points being able to decay into 3 different ones, while the D points may decay into2 distinct A points. The entire singularity transformation and unfolding sequence may betransparently summarized in the unfolding diagram below. A A i A i A i A i D ij D ij A j A j A j A j The unfolding diagram follows directly from Lagrangian catastrophe theory, although it canalso be derived from the caustic conditions. The unfoldings of an A ik singularities into an A ik − singularities, with k ≥ , follow trivially from the caustic conditions. The same holdsfor the unfolding of the D ij singularities into the D ij singularities. The decay from the D to the A singularities are proven in section 4.2.1. The mutations D ij → A i and D ij → A j follow directly since the shell-crossing of the D ij caustic is analogous to the shell-crossingcondition on the A i and A j caustics. The principal family of singularities – principal in terms of rate of occurrence and spatialdominance – is the A -family. They are induced by singularities in the geometric structure ofone of the eigenvalue fields. In physical terms, they involve one-dimensional collapse on to theemerging singularity. Of a more challenging nature within the evolutionary unfolding of thepatterns emerging in fluid flow is the formation of the D -family of singularities. They occurwhen two fold sheets corresponding to different eigenvalue fields intersect. Amongst others,– 40 –his means that the D singularities connect A singularities corresponding to two eigenvaluefields. A -family caustics The most prominent and abundant singularities are those of the two-dimensional fold sheets A i ( t ) . In Eulerian space, they mark the regions where mass elements are turned inside out asthe density attains infinity. This happens while they represent the locations where separatematter streams are crossing each other. As time proceeds, the fold sheets A i ( t ) sweep over anincreasingly larger Lagrangian region. Ultimately, integrating over time, they mark an entireLagrangian volume, which is labelled as A i . The A i set forms a three-dimensional variety.When we wish to identify where a particular individual fold sheet is born, we turnto the cusp points A i +3 . They are the points on the fold sheets where the correspondingeigenvalue field attains an extremum. Because of this, they mark the sites of birth of thefold singularities. As the A i ( t ) sheets unfold, at the edges their surface gets wrapped in ahigher order singularity, the cusp curves A i ( t ) . In time, these curves move through space andtrace out cusp sheets A i . In the context of the Megaparsec scale matter distribution in theUniverse, the cusp sheets are to be associated with the walls or membranes in the cosmic web[3, 15, 23, 51, 66].A dynamically interesting process occurs at the cusp points A i − , which are the saddlepoints of the corresponding eigenvalue field µ ti that at a given time are encapsulated by thefold sheet A i . At the A i − points, we see the merging or annihilation of fold sheets A i intoa larger structure (cf. figure 13). Mathematically, they mark the key locations where thetopology of the eigenvalue field changes abruptly. Physically, they are associated with themerging of separate structural components, a manifestation of the hierarchical buildup ofstructural complexity [23, 66].Also the cusp curves A i ( t ) can get folded. In Eulerian space, the folding of the cuspcurves manifests itself as A i ( t ) swallowtail points. As time proceeds, these points movethrough space and define the swallowtail curve A i . It is of interest to note that the swallowtailcurve is embedded in the cusp sheet, i.e. A i ⊂ A i . In the context of the cosmic structureformation process, the swallowtail curves outline and trace perhaps the most outstandingfeature of the cosmic web, the pronounced elongated filaments that form the of spine theweblike network [3, 23, 66].Also these features build up in a hierarchical process of small filaments merging intoever larger and more prominent arteries. In the context of the evolving singularity structurethat we study, this process is represented by the A i +4 points and A i − points. They definethe decisive junctions where significant changes in topology occur. For the A i ± points thisconcerns their identity in the gradient of the eigenvalue field, in which the A i +4 are maximaand minima and A i − points are the saddle points. The implication of this is that cusp curvesget created or annihilated at A i +4 points, while they merge or separate at A i − points.The final morphological constituent in this structural hierarchy of singularities is thatof the butterfly points A i . They conclude the A -family of singularities, i.e. the familyof singularities that correspond to the spatial characteristics of the field of one eigenvalue µ i . The swallowtail curves A i get folded at A i butterfly points. In the three-dimensionalstructural pattern that formed in the fluid, these will represent nodes. In the cosmic web,they define the nodal junctions, connecting to the various filamentary extensions that outline– 41 –ts spine [3, 15, 23, 26, 66]. In principle, for a given initial field and dynamical evolution,one might use these identifications to e.g. evaluate how many filaments are connected to thenetwork nodes [4, 56]. D -family caustics The A i ( t ) and A j ( t ) sheets, with i (cid:54) = j , intersect in the elliptic and hyperbolic umbilic points D ± ij ( t ) . In contrast to the A family of singularities, the collapse into D singularities is two-dimensional. It leads to the birth of the socalled umbilic points. Over time, they trace outthe umbilic curve D ± ij . The collapse process may occur in two distinctive ways, indicated bythe labels + and − .The topology of the variety D ± ij ( t ) changes at D ij ± and D points. An interesting char-acteristic of umbilic curves is that they are always created or annihilated in pairs. The D ij ± points correspond to the creation or annihilation of two D ± ij curves of the same signature.By contrast, the D ij points correspond to the creation or annihilation of a pair with one D + ij and one D − ij point. In this study we have developed a general formalism for identifying the caustic structureof a dynamically evolving mass distribution, in an arbitrary dimensional space. Through anew and direct derivation of the caustic conditions for the classification and characterizationof singularities that will form in an evolving matter field, our study enables the practicalimplementation of a toolset for identifying the spatial location and outline of each relevantclass of emerging singularities. By enabling the development of such instruments, and theapplication of these to any cosmological primordial density and velocity field, our studyopens the path towards further insight into the dynamics of the formation and evolution ofthe morphological features populating the cosmic web. In particular significant is that it willenable us to obtain a fundamental understanding of the spatial organization of the cosmicweb, i.e. of the way in which these structural components are arranged and connected.
Caustics are prominent features emerging in advanced stages of dynamically evolving fluids.They mark the positions where fluid elements cross and multi-stream regions form. Theyare associated with regions of infinite density, and often go along with the formation ofshocks. In the context of the gravitationally evolving mass distribution in the universe,caustics emerge in regions in which nonlinear gravitational collapse starts to take place. Assuch, they are a typical manifestation of the structure formation process at the stage whereit transits from the initial linear evolution to that of more advanced nonlinear configurationsinvolving gravitational contraction and collapse. The overall spatial organization of matter atthe corresponding scale is that of the cosmic web, which assembles flattened walls, elongatedfilaments and tendrils and dense, compact cluster nodes in an intricate multiscale weblikenetwork that pervades the Universe.Over the past decades our understanding of the formation and evolution of the cosmicweb has advanced considerably. The availability of large computer simulations have beeninstrumental in this, as they enabled us to follow the cosmic structure formation processin detail (see e.g. [60, 65, 69]). In combination with new theoretical insights [15, 66], this– 42 –as led to the development of a general picture of the emergence of the weblike matter andgalaxy distribution. The full phase-space dynamics of the process and its manifestation in theemerging matter distribution is an instrumental aspect of this that only recently received moreprominent attention. While the study by Zel’dovich [72] already underlined the importance ofa full phase-space description for understanding cosmic structure formation (see also [63, 64]),with the exception of a few prominent studies [11] the wealthy information content of full 6-Dphase-space escaped attention.A series of recent publications initiated a resurgence of interest in the phase-space aspectsof the cosmic structureformation process. They realized that the morphology of components in the evolvingmatter distribution is closely related to its multistream character [1, 30, 53, 58, 62] (foran early study on this observation see [21]). This realization is based on the recognitionthat the emergence of nonlinear structures occurs at locations where different streams of thecorresponding flow field cross each other. Looking at the appearance of the evolving spatialmass distribution as a 3D phase space sheet folding itself in 6D phase space, this establishesa connection between the structure formation process and the morphological classification ofthe emerging structure. Moreover, to further our understanding of the dynamical evolutionand buildup of the cosmic matter distribution, we also need to answer the question in howfar the various emerging structural features connect up in the overall weblike network of thecosmic web.
To be able to answer the questions, we study the emergence of singularities and causticsin a dynamically evolving mass distribution. Our analysis is built on the seminal work byArnol’d, specifically his classification of singularities in Lagrangian catastrophe theory. In athree-dimensional setting we can recognize two series of singularities, the A k and D k series.The 4 classes of A k singularities – A , A , A and A – are the singularities for which thecaustic condition holds for one eigenvalue. The D -family of umbilic singularities – includingthe D +4 , D − and D – are caustics for which the caustic conditions are satisfied by twoeigenvalue simultaneously. In three-dimensional fluids, the case in which all three eigenvaluessimultaneously satisify the caustic conditions, the E -family caustics, is non-degenerate.In order to detect these caustics in practice, we derived the caustic conditions , whichclassify them in terms of both eigenvalue and the eigenvector fields of the deformation tensor.The derivation differs from the classical derivation of catastrophe theory, in terms of gener-ating functions and the classification of its degenerate critical points, in that we work withthe geometry of the system. Moreover, the caustic conditions are not restricted to Hamil-tonian dynamics and apply to all systems which allow for a description with a sufficientlydifferentiable Lagrangian map. On the basis of the derived formalism, we show how the caustics of a Lagrangian fluid form anintricate skeleton of the nonlinear evolution of the fluid. The family of newly derived causticconditions allow a significant extension and elaboration of the work described in Arnold et al.(1982) [11]. Arnol’d et al. classified the caustics that develop in one- and two-dimensionalsystems that evolve according to the Zel’dovich approximation. While [8] did offer a quali-tative description of caustics in the three-dimensional situation, this did not materialize ina practical application to the full three-dimensional cosmological setting. The expressions– 43 –erived in our study, and the specific identification of the important role of the deforma-tion tensor eigenvectors, have enabled us to breach this hiatus. To identify the full spatialdistribution and arrangement of caustics in the evolving three-dimensional cosmic matter dis-tribution, we follow the philosophy exposed in the two-dimensional study by Hidding et al.2014 [31, 39]. By relating the singularity distribution to the spatial properties of the initialGaussian deformation field, [39] managed to identify and show the spatial connectivity ofsingularities and establish how in a hierarchical evolutionary sequence they evolve and mayultimately merge with surrounding structures.When applied to the Zel’dovich approximation for cosmic structure formation, the caus-tic conditions form a skeleton of the caustic web. In the context of the cosmic web, we mayidentify these singularities with different components. This observation by itself leads to someradically new insights into the origin of the structural features in the cosmic web. The A cusp singularities are related to the walls of the skeleton of the comsic web. The A swal-lowtail singularities trace the filamentary ridges and tendrils in the cosmic web. Also the D ± hyperbolic and elliptic umbilic singularities are related to the filamentary spine of the spine,as they define the dense filamentary extensions of the cluster nodes. The butterfly ( A ) andparabolic umbilic ( D ) singularities are both connected with the nodes of the weblike pattern.One immediate observation of considerable interest is that the line-like A and D ± featurestrace the connectivity of the cosmic web in meticulous detail. Perhaps equally or even moreinteresting, and of key importance for our understanding of the dynamical evolution of thecosmic web, is the observation that both filaments and tendrils, as well as nodes, may haveformed due to the folding by the phase-space sheet induced by only one deformation eigen-value: the filamentary A caustics and nodal A caustic belong to the one eigenvalue A familyof caustics. In other words, collapse along the second eigenvector is not necessary to createa filament-like structure, and not even collapse along both second and third eigenvector isneeded for the appearance of nodes (see [39, 40]). This is a new insight as it suggests theexistence of different possible late-time morphologies for filaments and nodes [40].A realization of key importance emanating from our work is that it is not sufficientto limit a structural analysis to the eigenvalues of the deformation tensor field. Usuallyneglected, we argue – and show by a few examples – that it is necessary to include theinformation contained in the (local) deformation tensor eigenvectors, our study has demon-strated and emphasized that for the identification of the full spatial outline of the cosmicweb’s skeleton. In an accompanying numerical study of the caustic skeleton in cosmological N -body simulations, we illustrate how essential it is to invoke the deformation eigenvectorsin the analysis [32]. This study will present a numerical and statistical comparison betweenthe matter distribution in the simulation and the caustic skeleton of the three-dimensionalcosmic web. Amongst the potentially most important applications of the current project is the fact thatthe caustic skeleton inferred from the Zel’dovich approximation adheres closely to the spineof the full nonlinear matter distribution. The direct implication is that we may directlylink the outline of the cosmic web to the initial Gaussian density and velocity field. Onthe basis of the corresponding deformation field, one may then attempt to calculate a rangeof properties analytically. The fact that we may invoke Gaussian statistics facilitates the– 44 –alculation of a wide range of geometric and topological characteristics of the cosmic web,as they are directly related to the primordial Gaussian deformation field, its eigenvaluesand eigenvectors. The first step towards this program were taken by [33]. A few examplesof results of such a statistical treatment for -dimensional fluids are described in [31]. Itdescribes how one may not only analytically compute the distribution of maxima, or minima,but also the population of singularities and the length of caustic lines. This will represent amajor extension of statistical descriptions that were solely based on the eigenvalue fields (seee.g. [29, 56]). Moreover, the ability to infer solid analytical results for a range of parametersquantifying the cosmic web will be a key towards identifying properties of the cosmic webthat are sensitive to the underlying cosmology. This, in turn, would enable the use of theseproperties to infer cosmological parameters, investigate the nature of dark matter and darkenergy, trace the effects of deviations from standard gravity, and other issues of generalcosmological interest.Notwithstanding the observation that the caustic skeleton inferred from the Zel’dovichapproximation appears to closely adhere to the full nonlinear structure seen in N -body sim-ulations, an aspect that still needs to be addressed in detail is the influence of the dynamicalevolution on the the developing caustic structure. This concerns in particular the descrip-tion of the dynamics of the system. Given the nature of singularities, the process of causticformation might be very sensitive to minor deviations of the mass element deformations andhence the modelling of the dynamics. This may even strongly affect the predicted populationof caustics and their spatial organization in the skeleton of the cosmic web. The Zel’dovichformalism [72] is a first-order Lagrangian approximation. A range of studies have shownthat second order Lagrangian descriptions, often named 2LPT, provide a considerably moreaccurate approximation of in particular the mildly nonlinear phases that are critical for un-derstanding the cosmic web [16, 18, 20, 21, 61]. In addition to a follow-up study in which weexplore the caustic structure according to 2LPT and possible systematic differences with thatpredicated by the Zel’dovich approximation, we will also systematically investigate the causticskeleton in the context of the adhesion formalism [36–38, 41]. Representing a fully nonlinearextension of the Zel’dovich formalism through the inclusion of an effective gravitational inter-action term for the emerging structures, it is capable of following the hierarchical buildup ofstructure. While it provides a highly insightful model for the hierarchically evolving cosmicweb, it also affects the flow patterns and hence the multistream structure in the cosmic massdistribution. In how far this will affect the caustic skeleton remains a major question for ourwork.Finally, of immediate practical interest to our project will be identification of the vari-ous classes of singularities that are populating the Local Universe. On the basis of advancedBayesian reconstruction techniques, various groups have been able to infer constrained re-alizations of the implied Gaussian primordial density and velocity field in a given cosmicvolume [43, 44, 49, 50]. From these constrained initial density and deformation fields, wemay subsequently determine the caustic structure in the Local Universe (see e.g. [40]). Theresulting caustic skeleton of the local cosmic web may then be confronted with the structures– clusters, groups and galaxies – that surveys have observed. Ultimately, this will enable usto reconstruct the cosmic history of objects and structures in the local Universe. In summary, the ability to relate the formation and hierarchical evolution of structure in theUniverse to the tale of the emergence and fate of singularities in the cosmic density field– 45 –rovides the basis for a dynamical theory for the development of the cosmic web, includingthat of its substructure. This will be the principal question and subject of the sequel to thework that we have presented here.
Acknowledgements
We thank Sergei Shandarin for having raised our interest in caustics as a key towards thedynamical understanding of the cosmic web. We are very grateful to Bernard Jones for acareful and diligent appraisal of the manuscript, and for the many useful and illuminatingdiscussions and comments. We also thank Adi Nusser, Neil Turok, and Gert Vegter for manyencouraging discussions and the anonymous referee for helpful comments. JF acknowledgesthe Perimeter Institute for facilitating this research through the support by the Government ofCanada through the Department of Innovation, Science and Economic Development Canadaand by the Province of Ontario through the Ministry of Research, Innovation and Science.
References [1] T. Abel, O. Hahn, and R. Kaehler. Tracing the dark matter sheet in phase space.
Mon. Not.R. Astron. Soc. , 427:61–76, November 2012.[2] R. J. Adler.
The Geometry of Random Fields . 1981.[3] M. A. Aragón-Calvo, E. Platen, R. van de Weygaert, and A. S. Szalay. The Spine of theCosmic Web.
Astrophys. J. , 723:364–382, November 2010.[4] M. A. Aragón-Calvo, R. van de Weygaert, and B. J. T. Jones. Multiscale phenomenology ofthe cosmic web.
Mon. Not. R. Astron. Soc. , 408:2163–2187, November 2010.[5] V. I. Arnol’d. Normal forms for functions near degenerate critical points, the Weyl groups of A k , D k , E k and Lagrangian singularities. Functional Anal. Appl , 6:1972, 1972.[6] V. I. Arnol’d. Wave front evolution and equivalent Morse lemma.
Communications in PureApplied Mathematics , 29:557–582, November 1976.[7] V. I. Arnol’d.
Mathematical methods of classical mechanics . 1978.[8] V. I. Arnol’d. Evolution of singularities of potential flows in collisionless media andtransformations of caustics in three-dimensional space.
Trudy Seminar imeni G Petrovskogo ,8:21–57, 1982.[9] V. I. Arnold.
Catastrophe theory . 1984.[10] V. I. Arnol’d and B. A. Khesin. Topological methods in hydrodynamics.
Annual Review ofFluid Mechanics , 24:145–166, 1992.[11] V. I. Arnol’d, S. F. Shandarin, and I. B. Zel’dovich. The large scale structure of the universe. I- General properties One- and two-dimensional models.
Geophysical and Astrophysical FluidDynamics , 20:111–130, 1982.[12] V.I. Arnol’d, S.M. Gusein-Zade, and A.N. Varchenko.
Singularities of Differentiable Maps,Volume 2: Monodromy and Asymptotics of Integrals . Modern Birkhäuser Classics. BirkhauserBoston, 2012.[13] V.I. Arnol’d, A. Varchenko, and S.M. Gusein-Zade.
Singularities of Differentiable Maps:Volume I: The Classification of Critical Points Caustics and Wave Fronts . Monographs inMathematics. Birkhauser Boston, 2012.[14] J. M. Bardeen, J. R. Bond, N. Kaiser, and A. S. Szalay. The statistics of peaks of Gaussianrandom fields.
The Astrophysical Journal , 304:15–61, May 1986. – 46 –
15] J. R. Bond, L. Kofman, and D. Pogosyan. How filaments of galaxies are woven into the cosmicweb.
Nature , 380:603–606, April 1996.[16] F. R. Bouchet, S. Colombi, E. Hivon, and R. Juszkiewicz. Perturbative Lagrangian approach togravitational instability.
Astron. Astrophys. , 296:575, April 1995.[17] J.W. Bruce. A classification of 1-parameter families of map germs r , → r , withapplications to condensation problems. J. Lond. Math. Soc. , II. Ser. 30(2):375–384, 1985.[18] T. Buchert. Lagrangian theory of gravitational instability of Friedman-Lemaitre cosmologiesand the ’Zel’dovich approximation’.
Mon. Not. R. Astron. Soc. , 254:729–737, February 1992.[19] T. Buchert. Lagrangian perturbation theory - A key-model for large-scale structure.
Astron.Astrophys. , 267:L51–L54, January 1993.[20] T. Buchert. Lagrangian Theory of Gravitational Instability of Friedman-Lemaitre Cosmologies- a Generic Third-Order Model for Nonlinear Clustering.
Mon. Not. R. Astron. Soc. , 267:811,April 1994.[21] T. Buchert and J. Ehlers. Lagrangian theory of gravitational instability of Friedman-Lemaitrecosmologies – second-order approach: an improved model for non-linear clustering.
Mon. Not.R. Astron. Soc. , 264, September 1993.[22] T. Buchert, A. L. Melott, and A. G. Weiss. Testing higher-order Lagrangian perturbationtheory against numerical simulations I. Pancake models.
Astron. Astrophys. , 288:349–364,August 1994.[23] M. Cautun, R. van de Weygaert, B. J. T. Jones, and C. S. Frenk. Evolution of the cosmic web.
Mon. Not. R. Astron. Soc. , 441:2923–2973, July 2014.[24] M. S. Chong, A. E. Perry, and B. J. Cantwell. A general classification of three-dimensional flowfields.
Physics of Fluids , 2:765–777, May 1990.[25] J. M. Colberg, K. S. Krughoff, and A. J. Connolly.
Mon. Not. R. Astron. Soc. , 359:272, 2005.[26] J. M. Colberg, K. S. Krughoff, and A. J. Connolly. Intercluster filaments in a Λ CDM Universe.
Mon. Not. R. Astron. Soc. , 359:272–282, May 2005.[27] P. Creminelli, A. Nicolis, L. Senatore, M. Tegmark, and M. Zaldarriaga. Limits onnon-Gaussianities from WMAP data.
Journal of Cosmology and Astroparticle Physics , 5:004,May 2006.[28] T. Delmarcelle.
The Visualization of Second-Order Tensor Fields . PhD thesis, STANFORDUNIVERSITY., 1995.[29] A. G. Doroshkevich. Spatial structure of perturbations and origin of galactic rotation influctuation theory.
Astrophysics , 6:320–330, October 1970.[30] B. L. Falck, M. C. Neyrinck, and A. S. Szalay. ORIGAMI: Delineating Halos UsingPhase-space Folds.
Astrophys. J. , 754:126, August 2012.[31] J. Feldbrugge, J. Hidding, and R. van de Weygaert. Statistics of Caustics in Large-ScaleStructure Formation.
Proceedings of IAU Symposium 308 “The Zeld’ovich Universe: Genesisand Growth of the Cosmic Web” , December 2014.[32] J. Feldbrugge, J. Hidding, and R. van de Weygaert. In preparation.
Mon. Not. R. Astron. Soc. ,2017.[33] Job Feldbrugge. Statistics of caustics in large-scale structure formation. Master’s thesis,Rijksuniversiteit Groningen, the Netherlands, 2014.[34] R. Gilmore.
Catastrophe theory for scientists and engineers . Wiley, 1981.[35] H. Goldstein.
Classical Mechanics . Addison-Wesley, 1980. – 47 –
36] S. N. Gurbatov, A. I. Saichev, and S. F. Shandarin. The large-scale structure of the universe inthe frame of the model equation of non-linear diffusion.
Mon. Not. R. Astron. Soc. ,236:385–402, January 1989.[37] S. N. Gurbatov, A. I. Saichev, and S. F. Shandarin. Large-scale structure of the Universe. TheZeldovich approximation and the adhesion model.
Physics Uspekhi , 55:223–249, March 2012.[38] J.. Hidding.
The Phase-Space Geometry of the Cosmic Web , Ph.D. thesis, University ofGroningen, 2018.[39] J. Hidding, S. F. Shandarin, and R. van de Weygaert. The Zel’dovich approximation: key tounderstanding cosmic web complexity.
Mon. Not. R. Astron. Soc. , 437:3442–3472, February2014.[40] J. Hidding, R. van de Weygaert, and S. Shandarin. The Zeldovich & Adhesion approximationsand applications to the local universe. In R. van de Weygaert, S. Shandarin, E. Saar, andJ. Einasto, editors,
The Zeldovich Universe: Genesis and Growth of the Cosmic Web , volume308 of
IAU Symposium , pages 69–76, October 2016.[41] J. Hidding, R. van de Weygaert, G. Vegter, B. J. T. Jones, and M. Teillaud. The StickyGeometry of the Cosmic Web.
ArXiv e-prints , May 2012.[42] V. Icke.
Astron. Astrophys. , 27:1, 1973.[43] J. Jasche, F. S. Kitaura, C. Li, and T. A. Enßlin. Bayesian non-linear large-scale structureinference of the Sloan Digital Sky Survey Data Release 7.
Mon. Not. R. Astron. Soc. ,409:355–370, November 2010.[44] F.-S. Kitaura. The initial conditions of the Universe from constrained simulations.
Mon. Not.R. Astron. Soc. , 429:L84–L88, February 2013.[45] E. Komatsu, A. Kogut, M. R. Nolta, C. L. Bennett, M. Halpern, G. Hinshaw, N. Jarosik,M. Limon, S. S. Meyer, L. Page, D. N. Spergel, G. S. Tucker, L. Verde, E. Wollack, and E. L.Wright. First-Year Wilkinson Microwave Anisotropy Probe (WMAP) Observations: Tests ofGaussianity.
The Astrophysical Journal Supplement , 148:119–134, September 2003.[46] I. A. Kravtsov and I. I. Orlov. Caustics, catastrophes, and wave fields.
Uspekhi FizicheskikhNauk , 141:591–627, December 1983.[47] L.D. Landau and E.M. Lifshitz.
Volume 1: Mechanics . Pergamom, 1959.[48] L.D. Landau and E.M. Lifshitz.
Volume 6: Fluid Mechanics . Butterworth-Heinemann, 1976.[49] F. Leclercq, J. Jasche, G. Lavaux, B. Wandelt, and W. Percival. The phase-space structure ofnearby dark matter as constrained by the SDSS.
Journal of Cosmology and AstroparticlePhysics , 6:049, June 2017.[50] F. Leclercq, J. Jasche, and B. Wandelt. Bayesian analysis of the dynamic cosmic web in theSDSS galaxy survey.
Journal of Cosmology and Astroparticle Physics , 6:015, June 2015.[51] N. I Libeskind, R. van de Weygaert, M. Cautun, B. Falck, E. Tempel, T. Abel, M. Alpaslan,M. A. Aragoon-Calvo, J. E. Forero-Romero, R. Gonzalez, S. Gottloober, O. Hahn, W. A.Hellwing, Y. Hoffman, B. J. T. Jones, F. Kitaura, A. Knebe, S. Manti, M. Neyrinck, S. E. Nuza,N. Padilla, E. Platen, N. Ramachandra, A. Robotham, E. Saar, S. Shandarin, M. Steinmetz,R. S. Stoica, T. Sousbie, and G. Yepes. Tracing the cosmic web.
ArXiv e-prints , May 2017.[52] P. J. Morrison. Hamiltonian description of the ideal fluid.
Reviews of Modern Physics ,70:467–521, April 1998.[53] M. C. Neyrinck. Origami constraints on the initial-conditions arrangement of dark-mattercaustics and streams.
ArXiv e-prints , February 2012.[54] P. J. E. Peebles.
The large-scale structure of the universe . 1980. – 48 –
55] Planck Collaboration, P. A. R. Ade, N. Aghanim, M. Arnaud, F. Arroja, M. Ashdown,J. Aumont, C. Baccigalupi, M. Ballardini, A. J. Banday, and et al. Planck 2015 results. XVII.Constraints on primordial non-Gaussianity.
Astron. Astrophys. , 594:A17, September 2016.[56] D. Pogosyan, C. Pichon, C. Gay, S. Prunet, J. F. Cardoso, T. Sousbie, and S. Colombi. Thelocal theory of the cosmic skeleton.
Mon. Not. R. Astron. Soc. , 396:635–667, June 2009.[57] T. Poston and I.N. Stewart.
Catastrophe Theory and Its Applications . Pitman, 1978.[58] Nesar S. Ramachandra and Sergei F. Shandarin. Multi-stream portrait of the cosmic web.
Monthly Notices of the Royal Astronomical Society , 452(2):1643–1653, 2015.[59] R. Salmon. Hamiltonian fluid mechanics.
Annual Review of Fluid Mechanics , 20:225–256, 1988.[60] J. Schaye, R. A. Crain, R. G. Bower, M. Furlong, M. Schaller, T. Theuns, C. Dalla Vecchia,C. S. Frenk, I. G. McCarthy, J. C. Helly, A. Jenkins, Y. M. Rosas-Guevara, S. D. M. White,M. Baes, C. M. Booth, P. Camps, J. F. Navarro, Y. Qu, A. Rahmati, T. Sawala, P. A. Thomas,and J. Trayford. The EAGLE project: simulating the evolution and assembly of galaxies andtheir environments.
Mon. Not. R. Astron. Soc. , 446:521–554, January 2015.[61] R. Scoccimarro. Gravitational Clustering from χ Initial Conditions.
Astrophys. J. , 542:1–8,October 2000.[62] S. Shandarin, S. Habib, and K. Heitmann. Cosmic web, multistream flows, and tessellations.
Phys. Rev. D , 85(8):083005, April 2012.[63] S. F. Shandarin and R. A. Sunyaev. The conjecture of the cosmic web. Commentary on:Zel’dovich Ya. B., 1970, A&A, 5, 84.
A&A , 500:19–20, June 2009.[64] S. F. Shandarin and Y. B. Zel’dovich. The large-scale structure of the universe: Turbulence,intermittency, structures in a self-gravitating medium.
Reviews of Modern Physics , 61:185–220,April 1989.[65] V. Springel, S. D. M. White, A. Jenkins, C. S. Frenk, N. Yoshida, L. Gao, J. Navarro,R. Thacker, D. Croton, J. Helly, J. A. Peacock, S. Cole, P. Thomas, H. Couchman, A. Evrard,J. Colberg, and F. Pearce. Simulations of the formation, evolution and clustering of galaxiesand quasars.
Nature , 435:629–636, June 2005.[66] R. van de Weygaert and J. R. Bond. Clusters and the Theory of the Cosmic Web. InM. Plionis, O. López-Cruz, and D. Hughes, editors,
A Pan-Chromatic View of Clusters ofGalaxies and the Large-Scale Structure , volume 740 of
Lecture Notes in Physics, BerlinSpringer Verlag , page 335, 2008.[67] V. A. Vasil’ev. Asymptotic exponential integrals, Newton’s diagram, and the classification ofminimal points.
Plenum Publishing Corporation , 1978.[68] M. Vergassola, B. Dubrulle, U. Frisch, and A. Noullez. Burgers’ equation, Devil’s staircasesand the mass distribution for large-scale structures.
Astron. Astrophys. , 289:325–356,September 1994.[69] M. Vogelsberger, S. Genel, V. Springel, P. Torrey, D. Sijacki, D. Xu, G. Snyder, D. Nelson, andL. Hernquist. Introducing the Illustris Project: simulating the coevolution of dark and visiblematter in the Universe.
Mon. Not. R. Astron. Soc. , 444:1518–1547, October 2014.[70] X. Wang, A. Szalay, M. A. Aragón-Calvo, M. C. Neyrinck, and G. L. Eyink. KinematicMorphology of Large-scale Structure: Evolution from Potential to Rotational Flow.
Astrophys.J. , 793:58, September 2014.[71] I. B. Zeldovich, A. V. Mamaev, and S. F. Shandarin. Laboratory observation of caustics,optical simulation of the motion of particles, and cosmology.
Uspekhi Fizicheskikh Nauk ,139:153–163, January 1983. – 49 –
72] Y. B. Zel’dovich. Gravitational instability: An approximate theory for large densityperturbations.
A&A , 5:84–89, March 1970.
A Zel’dovich approximation
The Zel’dovich approximation is the first order approximation of a Lagrangian pressurelessfluid evolving under self gravity, [72]. The Zel’dovich approximation is the simplest exampleof a Lagrangian fluid with Hamiltonian dynamics and serves as a good illustration of thecaustic conditions. The displacement map of the Zel’dovich approximation factors into aterm depending on time and a term depending on the initial conditions s t ( q ) = − b + ( t ) ∇ q Ψ( q ) , (A.1)with the linearized velocity potential Ψ( q ) and growing mode b + ( t ) . The growing mode can beobtained from linear Eulerian perturbation theory. Up to linear order, the linearized velocitypotential is proportional to the linearly extrapolated gravitational potential at the currentepoch φ ( q ) , i.e. Ψ( q ) = 23Ω H φ ( q ) , (A.2)with current Hubble constant H and current energy density Ω . The linearized velocitypotential Ψ( q ) encodes the initial conditions while the growing mode b + ( t ) encodes the cos-mological evolution of the fluid. For the Zel’dovich approximation it is common to define thedeformation tensor as ψ ij = ∂ Ψ( q ) ∂q i ∂q j (A.3)with eigenvalues λ i ( q ) satisfying µ i ( q, t ) = − b + ( t ) λ i ( q ) . The density in the Zel’dovich ap-proximation can be expressed as ρ ( x (cid:48) , t ) = (cid:88) q ∈ A ( x (cid:48) ,t ) ρ i ( q ) | − b + ( t ) λ ( q ) || − b + ( t ) λ ( q ) || − b + ( t ) λ d ( q ) | , (A.4)with ρ i the initial density field. Caustics occur at q at time t if and only if λ i ( q ) = 1 b + ( t ) (A.5)for at least one i . The eigenvalues λ i are functions determined by the initial gravitationalfield. Equation (A.5) can be pictured as a hyperplane at height /b + ( t ) . Since the Zel’dovichapproximation concerns potential flow, which means that the eigenvalues are real and can beordered such that λ ≥ λ ≥ λ . The intersection of this plane with the graph of the eigen-values undergoes shell-crossing at that time. For the Zel’dovich approximation the causticconditions in terms of the eigenvalues λ i are given by A = { q ∈ L | λ i ( q ) (cid:54) = 1 /b + ( t ) for all t and i } , (A.6) A i ( t ) = { q ∈ L | λ i ( q ) = 1 /b + ( t ) } , (A.7) A i ( t ) = { q ∈ L | q ∈ A i ( t ) and λ i,i ( q ) = 0 } , (A.8) A i ( t ) = { q ∈ L | q ∈ A i ( t ) and λ i,ii ( q ) = 0 } , (A.9) A i ( t ) = { q ∈ L | q ∈ A i ( t ) and λ i,iii ( q ) = 0 } , (A.10) D ± ij ( t ) = { q ∈ L | λ i ( q ) = λ j ( q ) = 1 /b + ( t ) and sign ( S M ) = ± } , (A.11) D ij ( t ) = { q ∈ L | q ∈ D ij ( t ) and ( λ i − λ j ) ,i ( q ) = ( λ i − λ j ) ,j ( q ) = 0 } , (A.12)– 50 –nd the points at which the topology of above sets changes A i +3 = { q ∈ L | q ∈ A i ∧ λ i ( q ) max-/minimum of λ i } , (A.13) A i − = { q ∈ L | q ∈ A i saddle point of λ i } , (A.14) A i +4 = { q ∈ L | q ∈ A i ∧ λ i,ii ( q ) max-/minimum of λ i,ii | A } , (A.15) A i − = { q ∈ L | q ∈ A i saddle point of λ i,ii | A } , (A.16) D ij ± = { q ∈ L | q ∈ D ± ij ∧ λ i ( q ) = λ j ( q ) max-/minimum of λ i | D ± ij = λ j | D ± ij } . (A.17)with the direction derivatives λ i,i = ∇ λ i · v i , λ i,ii = ∇ λ i,i · v i and λ i,iii = ∇ λ i,ii · v i . Notethat the eigenvectors are defined modulo multiplication by a real number and really representlines. B Non-diagonalizable deformation tensors
In sections 3 and 4 we derived the shell-crossing conditions and caustic conditions under theassumption that the deformation tensor ∂s t ∂q = M is diagonalizable. We here extend theseconditions to non-diagonalizable deformation tensors.The eigenvalues µ i are the roots of the characteristic function χ ( λ ) = det ( M − λI ) corresponding to the deformation tensor. The deformation tensor is diagonalizable if and onlyif the algebraic multiplicity – the order of the root – is equal to the geometric multiplicity – thenumber of eigenvectors corresponding to the root – for all eigenvalues. Hence, in order for thedeformation tensor to be non-diagonalizable, two or more eigenvalues need to coincide whilethere are fewer corresponding eigenvectors. In this case we can extend the set of eigenvectorsby adding the necessary generalized eigenvectors to put the deformation tensor in Jordannormal form M = VM J V − (B.1)where V is the generalized modal matrix consisting of the eigenvectors and generalized eigen-vectors and M J is the upper triangular matrix of Jordan normal form containing the eigen-values.In the three-dimensional case, these non-diagonalisable deformation tensors correspondto the hyperbolic/elliptic umbilic ( D ± ) caustics. For simplicity lets restrict to the three-dimensional case where the shell-crossing occurs due to the first and the second eigenvaluefields µ = 1 + µ = 0 . In this case, we extend the set of eigenvectors { v , v } by addingthe generalized eigenvector ¯ v . The Jordan matrix now takes the form M J = µ µ
00 0 µ . (B.2)Substituting equation (B.1) in equation (3.5) we obtain the condition that there needs toexist a non-zero tangent vector of D ± for which (1 + µ ) v ∗ · T + ¯ v ∗ · T = 0 (B.3) (1 + µ )¯ v ∗ · T = 0 (B.4) (1 + µ ) v ∗ · T = 0 , (B.5)– 51 –eplacing equations (3.12). We thus see that the D variety forms a D caustic if and onlyif ¯ v ∗ · T = 0 and v ∗ · T = 0 (for a diagonalizable deformation tensor we obtain the secondcondition). This is equivalent to the condition that T is parallel to the eigenvector v .Finally note that the deformation tensor can only be non-diagonalizable for non-Hamiltoniandynamics for which the parabolic umbilic caustic ( D ) is not stable (see section 5.2). Thiscondition is thus not relevant for Hamiltonian and generic non-Hamiltonian Lagrangian fluidsin three dimensions.This analysis straightforwardly generalizes to the case in which the geometric and alge-braic multiplicity of the eigenvalues differs by more than one for higher dimensional fluids. C Shell-crossing conditions: coordinate transformation
The shell-crossing conditions are manifestly independent of coordinate choices. However,the eigenvalue and eigenvector fields generally do depend on the choice of coordinates. Bythemselves, they do therefore not provide valid descriptions of the dynamics of the fluid.Suppose the displacement field can be written as s = ∇ ψ for some potential ψ . The Hessian H x of ψ , H ij = ∂ ψ∂x i ∂x j , (C.1)transforms non-trivially under the local coordinate transformation x → X ( x ) i.e. H → ˜ H = J T HJ + J T ∇ ( J ) ∇ ψ , (C.2)with J the Jacobian between the coordinate systems X and x , J ij = ∂X i ∂x j . (C.3)From this we immediately infer that the eigenvalue field and eigenvector fields are invariantif the transformation is orthogonal and global, i.e. if J T = J − and ∇ ( J ) = 0 . As may beexpected, these transformations include rotations and translations. D Lagrangian maps and Lagrangian equivalence
We here shortly describe the mathematical background of symplectic manifolds, Lagrangianmanifolds and Lagrangian maps. For a detailed description and derivations we refer to [12, 13].
D.1 Symplectic manifolds and Lagrangian maps A n -dimensional symplectic manifold ( M, ω ) is a smooth n -dimensional manifold M , equippedwith a closed nondegenerate bilinear 2-form ω called the symplectic form. Symplectic mani-folds are always even dimensional for ω to be nondegenerate. In Hamiltonian dynamics thesymplectic form ω can be associated to the Poisson brackets which encodes the dynamics ofthe theory. A Lagrangian manifold L of a n -dimensional symplectic manifold ( M, ω ) is a n -dimensional submanifold of M on which the symplectic form ω vanishes. Let ( B, π ) bea Lagrangian fibration of ( M, ω ) , which is a n -dimensional manifold with a projection map π : M → B for which the fibers π − ( b ) are Lagrangian manifolds for all b ∈ B .An example of a symplectic manifold is phase space consisting of position and canonical– 52 –omenta ( q , . . . , q n , p , . . . p n ) with the symplectic form ω = (cid:80) ni d q i ∧ d p i . An example ofa Lagrangian fibration is { ( q , . . . , q n ) , π } with the projection map π ( q , . . . , q n , p , . . . p n ) =( q , . . . , q n ) .Give a symplectic manifold ( M, ω ) with a Lagrangian fibration ( B, π ) we can for everyLagrangian manifold L define a Lagrangian map ( π ◦ i ) : L → M → B , with i being theinclusion map sending L into M . Two Lagrangian maps ( π ◦ i ) : L → M → B and ( π ◦ i ) : L → M → B are defined to be Lagrangian equivalent if there exist diffeomor-phisms σ, τ and ν such that τ ◦ i = i ◦ σ, ν ◦ π = π ◦ τ and τ ∗ ω = ω , or equivalently thediagram below commutes L ( M , ω ) B L ( M , ω ) B i π i π σ τ ν D.2 Displacement as Lagrangian map
Given a Lagrangian submanifold L we can construct a corresponding Lagrangian map. Firstmap the Lagrangian submanifold L with the inclusion map i : L → C to the correspondingpoints in phase space C , i.e., i : ( q, x ) (cid:55)→ ( q, x ) for all ( q, x ) ∈ L . Subsequently map thesepoints to a base manifold B with the projection map π : C → B . In Lagrangian fluid dy-namics it is convenient to pick the Eulerian manifold E as the base manifold B and definethe projection map as π : ( q, x ) (cid:55)→ x for all ( q, x ) ∈ C . As there will always be an exactcorrespondence between the Lagrangian manifold L and the Lagrangian submanifold L t ⊂ C (there exists a unique point x ∈ E such that ( q, x ) ∈ L t for every q ∈ L ), we can associate theLagrangian map corresponding to L t with the map x t . In summary, the map x t correspondsuniquely to a Lagrangian map for fluids with Hamiltonian dynamics.A Lagrangian map can develop regions in which multiple points in the Lagrangian man-ifold are mapped to the same point in the base space. The points at which the number ofpre-images of the Lagrangian map changes are known as Lagrangian singularities. Lagrangiancatastrophe theory classifies the stable singularities, stable with respect to small deformationsof L , up to Lagrangian equivalence. Lagrangian equivalence is a generalization of equivalenceup to coordinate transformations. For a precise definition of Lagrangian equivalence we referto appendix D. D.3 Lagrangian map germs
In catastrophe theory it is important to consider the Lagrangian map at a point. This isachieved by means of Lagrangian germs. Starting with a point p ∈ M we can consider La-grangian functions F i : U i → B for i = 1 , for small environments U i of p which coincide onthe intersection U ∩ U . The equivalence classes of such Lagrangian functions are Lagrangiangerms. The Lagrange equivalence of Lagrangian maps straightforwardly extends to Lagrangeequivalence of Lagrangian germs. These are the equivalence classes used in the classificationof stable Lagrangian maps, where a Lagrangian germ is stable if and only if every sufficientlysmall fluctuation on the germ is Lagrange equivalent to the germ.– 53 – .4 Gradient maps Every Lagrangian germ is Lagrange equivalent to the germ of a gradient map. That is to say,for every Lagrangian map l = π ◦ i : L → C → E we can for a point ( q, x ) ∈ L locally writethe map as l ( q , . . . , q n , x , . . . , x n ) = (cid:18) ∂S∂q , ∂S∂q , . . . , ∂S∂q n (cid:19) (D.1)for some function S : R n → R . The corresponding map x is given by x ( q , . . . , q n , t ) = (cid:18) ∂S∂q , ∂S∂q , . . . , ∂S∂q n (cid:19) (D.2)for some time t . By writing S = q + Ψ for Ψ : R × R → R we obtain x ( q, t ) = q + ∂ Ψ ∂q , (D.3)with the gradient field s = ∂ Ψ ∂q . (D.4)The Jacobian of the displacement map (cid:20) ∂s∂q (cid:21) ij = ∂ Ψ ∂q i ∂q j (D.5)is symmetric. The set of eigenvectors { v i } can be taken to be orthonormal by which the dualvectors coincide with the eigenvectors, i.e., v ∗ i = v i for all i . A Lagrangian map is locallyequivalent to the Zel’dovich approximation. D.5 Arnol’d’s classification of Lagrangian catastrophes
In section 4, we described the classification of Lagrangian singularities in up to three di-mensions. However the classification extends to higher dimensional singularities. A ( n + 1) -dimensional fluid can contain stable singularities in the A i , D i and E i classes with i ≤ n + 2 ,where the D -class range starts at i = 4 and the E -class is only defined for i = 6 , , . Thesesingularities decompose into lower-dimensional singularities as illustrated in the unfoldingdiagram below. A A A A A A A A A . . .D D D D D A . . .E E E – 54 – : x ( q,
1) = ( q , q , q ) 1 + µ = 1 1 + µ = 1 1 + µ = 1 A : x ( q,
1) = ( q , q , q ) 1 + µ = 1 1 + µ = 1 1 + µ = 2 q µ , = 2 A : x ( q,
1) = ( q , q , q q + q ) 1 + µ = 1 1 + µ = 1 1 + µ = 3 q + q µ , = 6 q µ , = 6 A : x ( q,
1) = ( q , q , q q + q ) 1 + µ = 1 1 + µ = 1 1 + µ = q + 4 q µ , = 12 q µ , = 24 q µ , = 24 A : x ( q,
1) = ( q , q , q q + q q + q ) 1 + µ = 1 1 + µ = 1 1 + µ = q + 2 q q + 5 q µ , = 2 q + 20 q µ , = 60 q µ , = 120 q µ , = 120 A ± : x ( q,
1) = ( q , q , ( q ± q ) q + q ) 1 + µ = 1 1 + µ = 1 1 + µ = q ± q + 3 q µ , = 6 q µ , = 6 A ± : x ( q,
1) = ( q , q , q q ± q q + q ) 1 + µ = 1 1 + µ = 1 1 + µ = q ± q q + 4 q µ , = ± q + 12 q µ , = 24 q µ , = 24 Table 2 : The caustic conditions of the normal forms of the A singularity classes E Caustic conditions of the normal forms
We here verify the caustic conditions for the normal forms in the generic classification ofsingularities given in section 5.2. The normal forms of the the Lagrangian singularities givenin section 5.3 follow analogously.The eigenvalue fields and corresponding derivatives in the direction of the eigenvectorfields are given in table 2. The eigenvalues of the normal form for the trivial ( A ) case equal and thus satisfy the condition µ i (cid:54) = 0 for all i . The third eigenvalue of the normal formof the fold ( A ) singularity equals −1