Theoretical dark matter halo kinematics and triaxial shape

In a recent paper, Salvador-Sol\'e et al. (2012) have derived the typical inner structure of dark matter haloes from that of peaks in the initial random Gaussian density field, determined by the power-spectrum of density perturbations characterising the hierarchical cosmology under consideration. In the present paper, we extend this formalism to the typical kinematics and triaxial shape of haloes. Specifically, we establish the link between such halo properties and the power-spectrum of density perturbations through the typical shape of peaks. The trends of the predicted typical halo shape, pseudo phase-space density and anisotropy profiles are in good agreement with the results of numerical simulations. Our model sheds light on the origin of the power-law-like pseudo phase-space density profile for virialised haloes.


INTRODUCTION
Virialised haloes in N -body simulations of cold dark matter (CDM) cosmologies show a wide variety of ellipsoidal shapes. On the contrary, their structural and kinematic properties are remarkably similar from one object to another. They are little sensitive not only to the mass, redshift, environment and even specific cosmology, but also to their individual shape. Only their scaling shows a mild dependence on some of these properties. As shown by gravo-hydrodynanic simulations, baryons introduce a larger scatter in the properties of haloes at their central region. However, in the present paper, we will concentrate on pure dark matter haloes and we will not deal with such secondary baryonic effects.
The typical spherically averaged halo density profile, ρ (r), is well-fitted, down to about one hundredth the virial radius, by the so-called NFW profile (Navarro et al. 1997) as well as by the Einasto (1965) profile, which gives slightly better fits down to smaller radii (Navarro et al. 2004;Merritt et al. 2005Merritt et al. , 2006Stadel et al. 2009;Navarro et al 2010). The velocity dispersion profile, σ(r), is reasonably well-fitted by the solution of the Jeans equation for spherically symmetric isotropic systems with null value at infinity (Cole & Lacey 1996;Merritt et al. 2006). More remarkably, Taylor & Navarro (2001) showed that the pseudo phase-space density profile is very nearly a pure power-law, with index ν ≈ −1.875 (see also Ascasibar et al. 2004;Rasia et al. 2004;Bett et al. 2007;Faltenbacher et al. 2007;Vass et al. 2009;Navarro et al 2010) and a similar relation arises from the radial velocity dispersion component, σr(r) (Dehnen & McLaughlin 2005). In equation (1), σ(r) is the velocity dispersion over the whole spherical shell with r, so the variance coincides with the spherical average of the local value at points with r, σ 2 (r) = σ 2 loc (r). Finally, Hansen & Moore (2006) found that the velocity anisotropy profile β(r) behaves linearly with the logarithmic derivative of the density, with a and b respectively equal to about −0.2 and 0.8 (Hansen & Stadel 2006; see also Ludlow et al. 2011 for an alternative expression), although with a substantial scatter this time.
The origin of all these trends is certainly related with the way dark matter clusters. In hierarchical cosmologies, haloes grow through continuous mergers with notably different dynamic effects according to the relative mass of the captured and capturing objects. For this reason, it is usually distinguished between major mergers, with a dramatic effect each, and minor mergers, contributing together with the capture of diffuse matter (if any) to the so-called accretion, responsible of just a smooth secular evolution of the system. Some authors have attempted to explain the typical halo density profile as the result of repeated major (or intermediate) mergers (Syer & White 1998;Salvador-Solé et al. 1998;Subramanian et al. 2000;Dekel et al. 2003). Others have concentrated instead in the effects of pure accretion (PA) (Avila-Reese et al. 1998;Nusser & Sheth 1999;Del Popolo et al. 2000;Manrique et al. 2003;Ascasibar et al. 2004;Salvador-Solé et al. 2007). Both extreme scenarios have also been investigated regarding the possible origin of the pseudo phase-space density and velocity anisotropy profiles (Hansen & Moore 2006). The PA scenario has received much support from the results by Wang & White (2009) showing that all typical halo trends are already set in the first generation haloes formed by monolithic collapse (i.e. no major merger; only accretion of diffuse matter) in warm dark matter cosmologies.
Regarding the shape, CDM haloes are found to be triaxial ellipsoids, with a trend towards prolate rather than oblate shapes (e.g. Frenk et al 1988;Dubinski & Carlberg 1991;Warren et al. 1992;Cole & Lacey 1996;Springel et al. 2004;Allgood et al. 2006;Hayashi et al. 2007;Macciò et al. 2007;Stadel et al. 2009;Vera-Ciro et al. 2011). Inside each individual object, the typical minor to major axial ratio takes a roughly uniform value of about 0.6, with a slight trend to an outward-decreasing triaxiality (Frenk et al 1988;Bullock 2002;Jing & Suto 2002;Springel et al. 2004;Kasun & Evrard 2005;Bailin & Steinmetz 2005;Allgood et al. 2006;Hayashi et al. 2007;Bett et al. 2007;Stadel et al. 2009;Vera-Ciro et al. 2011). The main axis is preferentially aligned, at all scales, along with the filament feeding the halo (e.g. Lemson & Kauffmann 1999;Basilakos et al. 2006;Patiri et al. 2006;Macciò et al. 2007;Ragone-Figueroa et al. 2010;Vera-Ciro et al. 2011). This indicates that the memory of the preferred direction of major mergers and accretion is not erased during virialisation (Vera-Ciro et al. 2011) or, equivalently, that the shape of virialised haloes depends on that of their seeds. Moreover, as haloes are not supported by rotation but by the local anisotropic velocity tensor, the fact that their triaxial shape is related to the shape of protohaloes automatically implies that their kinematics must also be related to it.
The seeds of haloes are believed to be peaks (secondary maxima) in the primordial random Gaussian density field filtered at the scale of the halo. The isodensity contours in the immediate vicinity of peaks are triaxial (Doroshkevich 1970) and rather prolate with a trend to become more spherical for very high peaks (Bardeen et al. 1986;hereafter BBKS). As well-known, the monolithic collapse of non-spherical systems is highly non-radial (Zeldovich 1970), giving rise to filaments and triaxial virialised objects. Thus, it is natural to believe that the shape of peaks is somehow translated into that of haloes, in agreement with the above mentioned alignments. A few authors (Lee et al. 2005;Rossi et al. 2011) have tried to make the link between the shape of haloes and that of peaks through the modelling of ellipsoidal collapse. Unfortunately, these models do not account for the highly non-linear effects of shell-crossing during virialisation, which play a crucial role in setting the final properties of virialised haloes. On the other hand, there is in the literature no attempt to relate the kinematics of haloes with the shape of their seeds.
In a recent paper, Salvador-Solé et al. (2012;hereafter SVMS) have shown that the kinematics of virialised haloes in (bottom-up) hierarchical cosmologies with dissipationless collisionless dark matter depends on their triaxial shape, contrarily to their spherically averaged density profile which does not. This allowed SVMS to infer, under the assumption of PA, the typical spherically averaged density profile for haloes from that of peaks in the primordial density field, determined by the power-spectrum of density perturbations. Furthermore, SVMS showed that the density profile for haloes having undergone major mergers is indistinguishable from that of haloes grown by PA, so the model actually holds for all haloes regardless of their individual aggregation history.
In the present paper, we extend the SVMS model to the kinematics and triaxial shape of virialised objects. Under the PA assumption and neglecting any possible rotation tidally induced by surrounding matter, we derive the halo shape, velocity anisotropy profile and pseudo phase-space density profile from the triaxial shape of peaks, taking into account the full virialisation process. We first assume the simple case of PA and then study the foreseeable effects of major mergers. This allows us to establish the link between those typical halo properties and the power-spectrum of density perturbations. The theoretical predictions obtained when this formalism is applied to CDM haloes are in good agreement with the results of numerical simulations.
The paper is organised as follows. In Section 2, we derive some general relations valid for triaxial systems, regardless of whether they are in equilibrium or not. Assuming PA, these relations are used, in Section 3, to make the link between the triaxial shape of a virialised object formed by PA and that of its seed. The typical velocity anisotropy and velocity dispersion profiles for virialised objects are derived in Section 4 from their triaxial shape. In Section 5, we apply the model to CDM haloes. The origin of the power-law-like form of the pseudo phase-space density profile and the effects of major mergers are discussed in Section 6. Our results are summarised in Section 7. A package with the numerical codes used in the present paper is publicly available from www.am.ub.es/∼cosmo/haloes&peaks.tgz.

MEAN SQUARED AND CROSSED FLUCTUATION PROFILES
There are in the literature several ways to characterise the shape of a triaxial system with semiaxes a b c: the axial ratios, the ellipticity and prolateness (e.g. BBKS) and the triaxiality parameter (e.g. Franx et al. 1991), among others. In the present paper, we use the primary and secondary eccentricities, respectively defined as As the eccentricities vary, in general, over the system (we assume from now on that the centre of symmetry remains the same at all scales), we will deal with the eccentricity profiles ep(r) and es(r), being r the radius of a sphere associated with each ellipsoid (see details below). Numerical studies of virialised haloes usually work with spherically averaged profiles. Of course, triaxial systems show non-null departures δρ(r, θ, ϕ) and δΦ(r, θ, ϕ) of the local density ρ(r, θ, ϕ) and gravitational potential Φ(r, θ, ϕ) from their respective spherical averages ρ (r) and Φ (r). Hence, the (scaled) mean squared or crossed density and potential fluctuation profiles, (δρ/ ρ ) 2 (r), (δΦ/ Φ ) 2 (r) and (δρ/ ρ )(δΦ/ Φ ) (r), also give a measure of the triaxial shape of those systems. Below, we relate these profiles to the eccentricity profiles defined above.
Take the Cartesian axis z aligned along with the major axis of the ellipsoidal isodensity contour ρ iso (r), labelled by the radius The local density at the point (r, θ, ϕ) then takes the form where eps stands for one of the two eccentricities and esp for the other one 1 . The spherically averaged density at r is then equal to where we have defined the function G(r) = 1 + 1 1 − e 2 s (r) + 1 1 − e 2 p (r) = a 2 (r) 1 a 2 (r) Dividing equation (5) by ρ (r) and replacing ρ iso (r) by its value given in equation (6), we obtain 1 + δρ ρ (r, θ, φ) = 3 G(r) sin 2 θ cos 2 φ + sin 2 θ sin 2 φ 1 − e 2 sp (r) + cos 2 θ 1 − e 2 ps (r) .
The mean squared density fluctuation over the sphere with radius r, can then be readily obtained from δρ/ ρ given by equation (8). The result is The relations between the eccentricity profiles and the mean squared potential and crossed density-potential fluctuation profiles follow from the relation (10) and the relations between such fluctuation profiles and the mean squared density fluctuation profile. To derive them we need first to develop the member on the left of the Poisson equation, Taking into account that, by the Gauss theorem, Φ satisfies the usual Poisson integral relation for spherically symmetric systems we arrive after some algebra at the exact relation where ξ(r) = (d ln M/d ln r) −1 and ζ(r) = −(d ln | Φ |/d ln r) −1 . Equation (13) shows that it is possible to infer δΦ/ Φ from δρ/ ρ for the appropriate boundary conditions. This requires, however, the full characterisation of δρ(r, θ, ϕ), which is, in general, a complicate function of the spatial coordinates. For simplicity, we will concentrate in ellipsoidal systems (i) with non-rotating symmetry axes and (ii) with rms density fluctuation profile of the power-law form, (δρ/ ρ ) 2 1/2 = Qr κ , or, equivalently, satisfying ∂[δρ(r, θ, ϕ)/ ρ (r)]/∂r = κr −1 δρ(r, θ, ϕ)/ ρ (r). 2 Condition (i) is essentially satisfied by simulated haloes (see Jing & Suto 2002;Springel et al. 2004;Kasun & Evrard 2005;Bailin & Steinmetz 2005;Hayashi et al. 2007;Vera-Ciro et al. 2011), while condition (ii) can be seen as a reasonable approximation, at least within some finite radial range 3 . In these conditions, multiplying equation (13) by the scaled density fluctuation, performing the spherical average and taking into account the divergence theorem, we are led to being K(r) = 1 + κξ(r)[2 + (1 − κ)ζ(r)]. Equation (14) is an ordinary differential equation that can be solved, for trivial consistency boundary condition at r = 0, for the mean crossed density-potential fluctuation. Similarly, multiplying equation (13) by the scaled potential fluctuation, the same derivation as above leads to Equation (15) cannot yet be solved for the mean squared potential fluctuation because of the presence of the unknown mean squared radial derivative of the scaled potential (the last term on the right). But this drawback can be bypassed by multiplying equation (13) being and I(r) = I(r) − κ/2. Once again, equation (17) can be solved, for trivial consistency boundary conditions at r = 0, for the mean squared potential fluctuation. Note that the mean ∂(δρ/ ρ )/∂r δΦ/ Φ profile entering the Jeans equation for anisotropic spherically averaged triaxial systems (eq. [42] below) can also be inferred from the radial derivative of the mean crossed density-potential fluctuation and the condition ∂[δρ(r, θ, ϕ)/ ρ (r)]/∂r = κr −1 δρ(r, θ, ϕ)/ ρ (r) above. One particular case of interest is that of homologous triaxial systems, i.e. with κ = 0. Equations (14), (15) and (17) That is, all squared or crossed fluctuation profiles are uniform and equal to each other. Another particular case of interest is that found in asymptotic regimes, i.e. with all the profiles behaving as power-laws. Equations (15), (16) and (17) then imply that both the rms density and potential fluctuation profiles have the same form, proportional to r κ . This does not mean that the isodensity and isopotential contours are equally non-spherical at any given radius: the positive proportionality factors Q and P of the rms density and potential fluctuation profiles, respectively, are different in general. Specifically, for an object with power-law spherically averaged density profile ρ (r) ∝ r −α , the ratio Q/P depends on κ and α in the following way, as inferred from equations (14) and (17). Figure 1 shows the ratio Q/P that results for α and κ respectively in the ranges 2 < α < 3 and −2 < κ < 2. The function is well-behaved everywhere except in the line κ = 4α − 10 where it diverges by alternating positive and nonphysical negative values on opposite sides of the line, in a series of segments separated by a few well-behaved points. Far enough from the strip, Q/P tends to take values greater than one for any value of κ. Notice that, as ρ(r, θ, ϕ) is positive, the negative values of δρ necessarily satisfy the condition |δρ| < ρ . Moreover, in rather prolate ellipsoidal systems the solid angle where δρ is negative is less than the solid angle where it is positive and, as the mean δρ is null, we also necessarily have that the positive values of δρ are less than ρ . Consequently, we have |δρ/ ρ | < 1, in general. Then, provided the mass distribution is smooth enough, which is guaranteed in triaxial virialised objects, equation (13) also implies |δΦ/ Φ | < 1. The angular cross-correlation of these two dimensionless quantities is even smaller, so any physical quantity reporting to the triaxial system can be written as a series expansion of the dimensionless fluctuation profiles measuring the deviation from spherical symmetry. We will take advantage of this fact in Sections 3 and 4.
All the previous relations hold regardless of whether the triaxial system is in equilibrium or not. Thus, they apply both to the virialised object and its seed. From now on, all quantities referring to the protoobject are denoted by subindex p.

ECCENTRICITY PROFILES
We want to relate the eccentricities of a virialised triaxial object with those of its seed. As there are two eccentricities, we need two equations. These equations must obviously report on properties involving the shape of these systems.
One of such properties is the volume of ellipsoidal isodensity contours. The ratio between the volumes of the ellipsoids encompassing a given mass in the final and initial systems must be equal, to leading order in the deviation from spherical symmetry, to the ratio between the volumes of the corresponding spheres encompassing identical mass, To write the member on the right of equation (21) we have taken into account that the mass of an ellipsoidal isodensity contour with semiaxes, a, b and c, labelled by R coincides, owing to the relations (5) and (6) arising from the definition (4) of the labelling radius, with the mass of the sphere with radius R, On the other hand, the function rp(r) in equation (21) is the solution, for the boundary condition rp(0) = 0, of the differential equation that follows from differentiation of equation (23) holding for spheres with identical mass in the initial and final systems. For simplicity in the notation, we omit from now on the explicit dependence of rp on r.
We stress that equation (21) is only valid to leading order in the deviations from spherical symmetry. This does not mean that the volume of each ellipsoid is approximated by that of the corresponding sphere. (By doing this, we would lose the information on the shape of the system.) What we approximate is the whole ratio between the volumes of both ellipsoids by the ratio between the volumes of the corresponding spheres. Were the axial ratios of the ellipsoid conserved over the evolution of the system, the relation (21) would be exact. Actually, the axial ratios vary during virialisation. But this variation is of higher order in the deviation from spherical symmetry, so the ratio between the volumes of the two ellipsoids is kept equal to that between the two spheres to leading order.
Taking into account the relation (3), equation (21) leads to A second property involving the triaxial shape of the system is the energy of the sphere with radius R encompassing a fixed mass M (see SVMS), This can be written in the form E(R) = E (R) + δE (R). The so-called "spherical" total energy, measures the total energy in the sphere with M , were the mass distribution inside it spherically symmetric and endowed with the velocity variance profile s 2 (r) that would result taking into account the energy lost by shell-crossing during virialisation but not the potential energy exchanged with the rest of the system owing to its non-spherical mass distribution. In the protoobject, s 2 p (rp) can be taken equal to the real velocity variance profile σ 2 p (r), but in the virialised object s 2 (r) differs from σ 2 (r). Nonetheless, the profile s 2 (r) can be obtained, by differentiation of equation (27), from the mass and spherical energy profiles, M (r) and E (r), respectively, derived in SVMS. The result was and where we have defined the "spherical" potential energy, This leads to the following explicit expression showing that s 2 (r) can be calculated from the spherically averaged density profile ρ (r). As this latter profile is independent of the shape of the system, so must the function E (r) (see eq. [27]). Consequently, all the effects of triaxiality in E(r) are included in the residual Specifically, the part of it harbouring the difference σ 2 (r) − s 2 (r) corrects the spherical total energy in the sphere from the gravitational energy exchanged (it may be positive or negative) at a distance with the rest of the system, and the other part with the crossed density-potential fluctuation profile corrects it for the different potential energy due to the actually non-spherically symmetric mass distribution of the system. As virialisation is driven by the energy lost through shell-crossing (not by the energy exchanged at a distance between shells, which is a consequence of the former effect and much less important), the ratio between the final and initial total energies in the sphere is, to leading order in the deviation from spherical symmetry, equal to the ratio between the corresponding spherical total energies also accounting for that energy loss We stress, once again, that equation (34) is valid to leading order in the deviation from spherical symmetry. This does not mean that the total energy in the sphere is approximated by its spherical counterpart. (By doing this we would lose all the information on the shape of the system.) It is the whole ratio between the total energies in the final and initial spheres which is approximated by the ratio between their spherical counterparts. Thus, taking into account the relation E(R) = E (R) + δE (R), equation (34) implies to leading order in the deviation from spherical symmetry. Substituting the expressions for δE and δEp in the virialised object and its seed (eq. [33]) into the relation (35) and differentiating it, we are led to where D(R) ≡ E (R)/Ep(R), the so-called "spherical" energy dissipation factor, is given by (eqs. [29] and [30]) Then, taking into account equations (10) and (14), both for the virialised object and its seed, equation (36) establishing the link between the respective mean crossed density-potential fluctuation profiles leads to where and S(r) = 5 2 being and γ(r) the logarithmic derivative of the mean crossed density-potential fluctuation profile.
Equations (25) and (38) determine the eccentricities of the virialised object grown by PA from those of its seed, as wanted. From equations (38) and (39), we see that this relation involves the profile σ 2 (r), which depends on the shape of the final system and, hence, on that of its seed. But such a dependence cannot be determined without previously determining the anisotropy profile. This will be done in Section 4. Only at small r, the term S(r) in equation (38) becomes negligible 4 . This means that, except for the factor U (r), equation (38) takes there the form of an identity relation, just as equation (25). Interestingly, the set of algebraic equations (38) and (25) is solvable only for a very narrow range of U (r) values around unity. This is shown in Figure 2, where we plot the solution space for varying values of U (r) inferred by means of the algebraic solving procedure provided in the Appendix. Consequently, the eccentricities of the virialised object near the centre are necessarily close to those of the protoobject at the same location.

KINEMATIC PROFILES
As mentioned, the kinematics of a non-rotating triaxial virialised object must depend on its shape and conversely. This relation is however hard to establish from first principles, even under the assumption of PA as made here. It can nonetheless be guessed from comparison with the simplest idealised case of strict spherical symmetry. See also Section 5 for a more fundamental justification.
In that idealised case, the final velocity dispersion would be purely radial because the system would expand, collapse and virialise radially (shell-crossing would take place between spherical shells). The total energy within spheres is, to leading order in the dimensionless fluctuation profiles measuring the deviation from spherical symmetry, the same in spherically symmetric systems as in triaxial ones (see eq. [26]). Thus, non-radial infall due to non-spherical shells (or to any possible instability; see below) should cause the transfer from radial to tangential kinetic energy together with a deviation of the gravitational potential from its spherical average without altering, to leading order, the total energy of the system. Assuming one-to-one correspondence between the two effects, we are led to the conclusion that, to leading order, the fractional velocity variance transferred from the radial to the tangential direction (equal to half the fractional 1-D tangential velocity variance generated) should be equal to half the typical fractional deviation of the potential from its spherical average, Equation (40) might seem to imply that it is impossible to have tangential velocities in spherically symmetric virialised objects (with null δΦ), which would be manifestly wrong. What it actually implies is this very conclusion but for objects formed by PA 5 . Indeed, as discussed in SVMS, such objects are always ellipsoidal (they have non-null δΦ). Even if the seeds are spherically symmetric, the gravitational effect of the surrounding matter (non-null shear tensor) or radial orbit instability, in the artificial case of isolated seeds, automatically leads to non-radial collapse and the formation of triaxial virialised objects (see Sec. 5 for more details). Given the relation (40), the anisotropy profile, is a function of the rms potential fluctuation profile for the virialised object. As shown in Section 3, such a fluctuation profile is related to its counterpart in the seed. However, such a relation involves the velocity dispersion profile itself (eq. [38]). Thus, σ 2 (r) must be determined before the rms potential fluctuation profile.
To do this we will make use of the (exact) generalised Jeans equation for steady non-spherical virialised objects (see SVMS), which, taking into account equation (41), can be written as . (42) The member on the right of equation (42) is of second order in the dimensionless fluctuation profiles. Like β(r), these profiles are functions of their counterparts in the seed and of σ(r). Consequently, equation (42) is a differential equation for σ(r), which can be solved for the usual boundary condition of null dispersion at infinity 6 . Moreover, the member on the right of equation (42) is of second order in the deviations from spherical symmetry, while all remaining terms are of lower order; the anisotropy profile itself is correct only to first order. Therefore, it can be neglected in equation (42), which notably simplifies the solution of this equation. Once σ(r) has been determined, the eccentricities of the virialised object can be calculated from equations (25) and (38) and β(r) from equation (41).

APPLICATION TO CDM HALOES
In SVMS, the spherically averaged density profile for CDM haloes was shown to arise from the spherically averaged density profile of halo seeds. As this latter profile is determined by the power-spectrum of density perturbations, it was possible to infer the typical halo density profile for CDM haloes directly from the power-spectrum of density perturbations. The extension of the model presented in the present paper relates the shape and kinematics of CDM haloes to the shape of their seeds, which is also determined, of course, by the power-spectrum of density perturbations. Thus, it should also be possible to infer the typical shape and kinematics of haloes directly from the power-spectrum of density perturbations.
To do this we need to know the eccentricity profiles for protohaloes, while what is only known is the typical eccentricities at peaks (i.e. at rp = 0) in the initial density field filtered by a Gaussian window, calculated by BBKS. Nonetheless, the former can be accurately inferred from the latter, following the procedure developed in SVMS for the spherically averaged density profile.
As pointed out in SVMS, the fact that haloes grown by PA develop from the inside out implies that the density profile for each halo ancestor in the continuous series leading to the final object exactly matches that for the halo within one triaxial isodensity contour. All these halo ancestors evolve from peaks in the initial density field filtered at the ancestor scale. Thus, the typical spherically averaged density contrast profile for the protohalo, δp (rp), convolved by a Gaussian window with any radius R f must be equal to the typical density contrast of a peak, δ pk (R f ), in the density field filtered at scale R f , The typical trajectory δ pk (R f ) followed by peaks evolving by PA into the series of ancestors leading to a typical halo with M at t is the solution of the differential equation (Manrique & Salvador-Solé 1996), for the appropriate boundary condition (through the δ pk (t) and R f (M ) relations given by eqs.
[40] and [41] of SVMS). In equation (44), the inverse of xe(R f , δ pk ) is the mean inverse curvature x (i.e. minus the Laplacian over the second order spectral moment) of peaks with δ pk at R f . The differential equation (44) can be solved, which then allows one to inverse the Fredholm integral equation (43) and obtain the wanted density contrast profile for the typical protohalo (see SVMS for details). This procedure was followed in SVMS to infer the unconvolved spherically averaged density profile for the seed of a typical halo.
In the present case, we must take into account that the squared semiaxes of each peak over the peak trajectory leading to the typical halo with M at t coincide (except for an arbitrary scaling) with minus the second order spatial derivatives λ 2 j of the filtered density contrast at the peak, in Cartesian coordinates xj with x1 aligned along the major semiaxis. They are therefore related to the density contrast profile, δp(r), for the protohalo through Writing the density contrast profile in the halo in terms of its spherical average (this relation is similar to that between the respective total densities, given by eqs. [5] and [6]), performing the second order spatial derivatives of the Gaussian window inside the integral on the right of equation (45) and taking x = 0, changing to spherical coordinates, integrating over θ and φ and averaging over the peak curvature, we arrive at where x (R f , δ pk ) is the mean curvature of peaks with δ pk and R f and Aj(R f , δ pk ) (j = 1÷3) are its associated dimensionless semiaxes λj scaled to the square root of the Laplacian, so that they satisfy A 2 1 + A 2 2 + A 2 3 = 1. These semiaxes depend on δ pk and R f through x in a well-known form calculated by BBKS. Thus, we can calculate them as well as the curvature x over the peak trajectory δ pk (R f ) leading to the halo with M at t and invert equation (46) by means of the same procedure as used for equation (43). In addition, the integral on the left of equation (46) can also be calculated from the known spherically averaged protohalo density contrast profile. Thus, taking the ratios between the solutions for different j obtained from inversion of equation (46), we can infer the eccentricity profiles for the seeds of typical haloes grown by PA.
Once the eccentricity profiles, (ep)p(rp) and (es)p(rp), are known, the typical halo shape and kinematics can be inferred following the steps described in Sections 3 and 4: 1) infer (through eqs.
[25] and [38]) the eccentricity profiles ep(r) and es(r) of the halo as a function of σ(r); 2) determine from them the mean squared density fluctuation profile for the halo (eq. [10]) as a function of σ(r); 3) infer the corresponding mean squared potential fluctuation profile (eqs [14] and [17]) as a function of σ(r); 4) solve the generalise Jeans equation (42) for the velocity dispersion profile σ(r); and 5) determine (from eqs. [25] and [38]) the typical eccentricity profiles, ep(r) and es(r), now using the explicit values of σ(r), and the anisotropy profile β(r) (eq. [40] after deriving the mean squared potential fluctuation profile).
Unfortunately, the accurate theoretical eccentricities so inferred could not be compared with the results of numerical simulations because there is so far no such typical profiles drawn from simulations. We just know some main trends of the typical shape of haloes. Thus, it is not worth at this stage to carry out such a complex derivation. In fact, the purpose of the present paper is not to infer accurate eccentricity and kinematic profiles for CDM haloes directly from the power-spectrum of the concordance model but rather to verify the validity of the model and try to understand the origin of the universal trends shown by such halo properties in numerical simulations. And these two objectives are much better obtained by means of the following simpler and comprehensive approximate procedure.
Peaks are triaxial with a rather prolate shape (BBKS). As mentioned, their eccentricities depend on δ pk and R f through their curvature x in such a way that, the larger x, the more spherical is the peak. The quantity x increases progressively with increasing R f over typical peak trajectories, meaning that peaks become increasingly spherical over those tracks, very slowly first and much more rapidly at the end. Thus, the same trend must be found in protohaloes for increasing rp as well as in typical haloes for increasing r. (It is true that the relation between the eccentricities in the seed and the halo also depends on σ(r), but this dependence can be neglected at small radii, while at large radii the decrease with increasing rp of protohalo eccentricities is so marked that the dependence on σ(r) cannot reverse it.) Therefore, according to the present model, haloes should be rather prolate and approximately homologous at small and intermediate radii and tend to become more spherical near the halo edge. This behaviour is fully consistent with the results of numerical simulations (e.g. Jing & Suto 2002;Bailin & Steinmetz 2005;Allgood et al. 2006;Stadel et al. 2009).
Moreover, we can have an idea on the overall values of halo eccentricities because, as mentioned, they are expected to be rather uniform except near the edge. The typical ratio between the major and minor axes (the major and intermediate axes) in peaks with low and moderate average curvature is ∼ 1.7 (∼ 1.3) (BBKS), implying typical values of the eccentricities ep (es) of ∼ 0.81 (∼ 0.64). Such typical values, independent of the filtering scale (they only show a moderate dependence on the curvature) were derived by BBKS for peaks in the old standard CDM model (i.e. in the Einstein-de Sitter cosmology), so they are valid for any CDM cosmology as all cosmologies converge for high enough redshifts as those corresponding to primordial peaks to the Einstein-de Sitter model. As U (r) is close to one with a small bias towards two at small r, the expected values of  ep (es) in haloes at small radii should be close to, perhaps a little larger than, those of peaks, say, ∼ 0.9 (∼ 0.8). This is also consistent with the values found in numerical simulations (Frenk et al 1988;Bullock 2002;Jing & Suto 2002;Springel et al. 2004;Kasun & Evrard 2005;Bailin & Steinmetz 2005;Allgood et al. 2006;Hayashi et al. 2007;Bett et al. 2007;Stadel et al. 2009).
At large r, the dependence on σ 2 (r) of the relation between the eccentricities in the halo and the seed becomes important. Yet, we can still determine the expected outer asymptotic behaviour of the eccentricities in the halo, taking only into account that, as mentioned, the object becomes more spherical at very large r. In that asymptotic regime, the power indexes of the mean squared density, squared potential and crossed density-potential fluctuation profiles coincide (see Sec. 3) and the more or less marked departure of isodensity and isopotential contours from spheres, given by the ratio Q/P , depends on the index κ and the asymptotic logarithmic slope −α of the density profile (eq. [20]). For values of α approaching 3 as near the outer asymptotic regime of virialised haloes, the only possible negative values of κ leading to outward-decreasing eccentricities with isopotential contours moderately more spherical than the isodensity contours (Q/P not much greater than one) as found in the literature Kasun & Evrard 2005;Hayashi et al. 2007) are in the range −0.2 κ < 0 (see the regions in dark blue, or quite strong grey if black and white, in Fig. 3). Thus, the outer asymptotic logarithmic slope, κ, is also severely constrained in this approximate treatment.
And what about the typical halo kinematics? The typical mean squared potential fluctuation profile for CDM haloes in log-log should be approximately given by the natural spline between a uniform value of ∼ 0.11, as implied by the above mentioned typical eccentricities, below about one hundredth the virial radius and a straight line with logarithmic slope κ beyond about one tenth the virial radius. We will consider two values of κ, −0.175 and −0.1, in order to better sample the allowed range −0.2 κ < 0 and to see how robust the theoretical typical halo kinematics are against variations in the individual shapes of haloes. (Different triaxial shapes must give rise to different mean squared potential fluctuation profiles.) The two squared potential fluctuation profiles so built are shown in Figure 4.
Once the form of the squared potential fluctuation profile has been fixed, we can calculate the ratio σ 2 t (r)/σ 2 (r) (eq. [40]) and, from it, the anisotropy profile (eq. [41]). Figure 5 shows the velocity anisotropy profiles inferred from the two different squared potential fluctuation profiles. As can be seen, they are both in good agreement with the results of numerical simulations. Specifically, for κ = −0.1, we find an anisotropy profile that closely follows the typical 'universal' law proposed by Hansen & Stadel (2006) (the red line in the right panel). While, for κ = −0.175, we find it closer to the anisotropy profile found by Navarro et al (2010) in a realisation of a Milky Way mass halo (the red line in the left panel). The only noticeable disagreement in this latter case is that the theoretical profiles show no cutoff near the halo edge as found by Navarro et al But this is a very unusual feature (it does not follow the 'universal' trend) and likely reflects an incomplete relaxation of this particular simulated halo at those large radii.
Using the previous theoretical velocity anisotropy profiles, we can solve the generalised Jeans equation (42) to find the corresponding velocity dispersion profiles and, using the typical spherically averaged halo density profile derived in SVMS,   the corresponding pseudo phase-space density profiles, ρ (r)/σ 3 (r). In Figure 6, we show the results obtained from the two mean squared potential fluctuation profiles above. In both cases, the theoretical pseudo phase-space density profile is very nearly a power-law with a logarithmic slope of −1.875, in full agreement with the results of numerical simulations. The two solutions essentially overlap; there is just a slight undulation at intermediate radii in the solution obtained from κ = −0.175, which indicates that the spline carried on between the inner and outer asymptotes of the squared potential fluctuation profile is somewhat deficient in this case. Nonetheless, the fact that so distinct mean squared potential fluctuation profiles as adopted here give so similar solutions demonstrates that the pseudo phase-space density profile is very robust against variations in halo shape (see Sec. 6 for the possible origin of this behaviour). Note instead the clear dependence of the zero-point of the pseudo phase-space density profile on the mass of the halo predicted by the model. This cannot be compared with simulations as this particular dependence has not yet been analysed. It would be worth trying to confirm this prediction of the model. According to our results, the velocity anisotropy profile should show an important dependence on the shape of the system, which is consistent with the rather large scatter shown by this profile in numerical simulations. The velocity dispersion profile follows from the Jeans equation (42) for such an anisotropy profile and the density profile ρ (r). As the latter profile is independent of the shape of the object, it cannot balance in the Jeans equation the dependence of β(r) on that property. Consequently, σ(r) must also depend on the shape of the object. Yet, the pseudo phase-space density profile turns out to be very robust against variations in that property; there are just very slight undulations depending on it. Consequently, σ(r) must also be similarly robust. This is once again in agreement with the small scatter shown by this profile in simulated haloes.
To end up we want to remark that we have made so far no reference to the CDM cosmology used; all the preceding (approximate) results hold, regardless of the cosmology provided only haloes grow hierarchically from the bottom-up (see SVMS). This is once again in agreement with the results of numerical simulations. Of course, the finer details of the predicted profiles will depend on the specific power-spectrum used, but this will not modify the main properties discussed here.

DISCUSSION
As just shown, the present model recovers the right form of the anisotropy and pseudo phase-space density profiles of simulated haloes. The origin of the former profile is clear: it arises from the tight relation (40) between the tangential velocity dispersion and the mean potential fluctuation profile and the typical form (essentially the well-constrained inner and outer asymptotic behaviours) of this latter profile. But the origin of the power-law-like pseudo phase-space density profile is less obvious. On the other hand, we have assumed so far that haloes grow by PA while they grow both through PA and major mergers. In this Section, we address these two important points.

Origin of the Pseudo Phase-Space Density Profile
The form of the pseudo phase-space density profile of simulated haloes, recovered by the present model, coincides with that predicted by Bertschinger's (1985) model for collisionless spherically symmetric, self-similar, systems evolving by PA. This implies that the ultimate reason for such a coincidence cannot be neither the shape (spherically symmetric vs. triaxial) of the system nor the form (power-law vs. non-power-law) of the density profile. The only more fundamental aspect in the formation of virialised objects grown by PA directly or indirectly included in both models is the way the coarse-grained phase space density increases via phase mixing (entropy generation) when the system loses energy via shell-crossing and contracts.
Bertschinger's model is dynamical and can properly follow the effects of shell-crossing, while the SVMS model is steady and focuses on the equilibrium state resulting from that process. However, the inside-out growth condition the SVMS model relies on directly follows from the lack of apocentre-crossing during virialisation, a condition that is also fulfilled by Bertschinger's model. Thus, despite the distinct approach, the effects of shell-crossing are similarly accounted for in the two models.
Certainly, the power-law form of the pseudo phase-space density profile predicted by Bertschinger's model is tightly related to the self-similarity assumption, while the present model makes no such an assumption. However, PA is always close to self-similar. It is driven by gravitation, a scale-free force, and the initial conditions are also very approximately scale-free: the initial mass distribution is essentially uniform (with essentially the critical density) and the velocity field corresponds to the scale-free unperturbed Hubble-flow form. Thus, the solution of the present model at small z should not be very far from self-similar, except for the effects of the increasing departure from the Einstein-de Sitter universe owing to the specific density profile of the seed. In this sense, the predicted profiles should not be far from power-laws, particularly the pseudo phase-space density profile that appears to be so little sensitive to the details of the seed. Another fundamental assumption of Bertschinger's model is the radial infall of the system, while the SVMS model takes into account that the collapse and virialisation is actually non-radial. However, if the tangential velocities develop at the expense of the initial radial velocities without altering the total velocity dispersion, as assumed in equation (40), the velocity dispersion that emerges as a consequence of shell-crossing should be the same in spherically symmetric as in triaxial systems.
The situation is therefore as follows. The power-law form, with index equal to -1.875, of the pseudo phase-space density profile found by Bertschinger is the direct consequence of the shell-crossing produced in self-similar radial PA. Thus, the conditions that: i) PA is always quite close to self-similar and ii) the velocity dispersion generated by shell-crossing is not modified, according to the condition (40), by the appearance of tangential velocities in non-radial PA, should guarantee that the evolution of the coarse-grained phase-space density found in Bertschinger's model is kept essentially unaltered in the general (non-striclty self-similar and non-radial) case of PA. This would explain why the present model recovers the pseudo phase-space density profileà la Bertschinger in agreement with the results of numerical simulations. Note that such a pseudo phase-space density profile does not necessarily implies that the density profile must also be roughly proportional to r −2.25 as found by Bertschinger. This would be the case if PA were strictly self-similar, so to end up with a power-law density profile, and, more importantly, radial. As stated by Bertschinger (1985) when comenting on his density profile ρ ∝ r −2.25 : "Self-similar relaxation is not complete [as it proceeds along one direction only] and thus heuristically should not lead to a density profile as flat as ρ ∝ r −2 [as predicted by Lynden-Bell (1967) in the case of complete 3D relaxation]. It is uncertain to what extent the addition of angular momentum (i.e. non-radial orbits) will change this result,...". In other words, if PA is non-radial, than the power index of the density profile is expected to change, in agreement with the results of SVMS, despite the fact that the pseudo-phase space density profile, driven by the evolution of the phase-space density during virialisation in the absence of apocentre-crossing, should be kept roughly unaltered as a consequence of condition (40).
As mentioned, a key point in the present model, arising from the particular way shell-crossing proceeds in accreting virialised haloes, is their inside-out growth. Bertschinger's model is instead based on the self-similar evolution of such systems and does not care about the implications this has on the growth of the steady object. But, if our explanation of the power-lawlike pseudo phase-space density profile is correct, the two models should essentially coincide (except for the symmetry and the strict self-similarity condition), meaning that Bertschinger's solution should also approximately satisfy the inside-out growth condition. This is indeed the case. As stated by Bertschinger, as a consequence of shell-crossing, particle orbits rapidly become approximately periodic, that is, the mass inside them stabilises and the inner system does not change anymore. In other words, the system develops a steady core and grows from the inside-out as new shells virialise. Of course, such a behaviour can only be approximately satisfied in Bertschinger's model. As a consequence of self-similarity, steadiness can only be strictly achieved in the limit of infinite time or vanishing radii; particle orbits can never become exactly periodic, they can only tend to become so (see Betschinger's Fig. 9 at the base of his comment regarding the inside-out growth). To obtain a fully steady solution, the self-similarity assumption must be replaced by the inside-out growth condition as in the present model.

Major Mergers
As has long been known, there is in radial PA a one-to-one correspondence between the density profiles of virialised objects and their seeds (Bertschinger 1985;Del Popolo et al. 2000). The reason for this is that in radial PA, there is no phase-mixing of particles turning around at the same time (Bertschinger 1985)) or, equivalently, there is no apocentre-crossing of particle orbits during virialisation. This result was extended in SVMS to non-radial PA.
An important consequence of this result was that the spherically averaged density profile ρ (r) for a virialised halo does not allow one to tell whether or not it has suffered major mergers 7 . The reason for this is that, given a virialised halo with any arbitrary aggregation history, we can always think about one peak with appropriate density profile (or equivalently, the spherical energy distribution Ep(M )), leading by PA to a virialised object with identical spherically averaged density profile. Such a putative accreting seed of the halo really exists; it is the peak tracing the current evolution by smooth accretion of the halo according to the peak formalism. Thus, it is a normal peak contributing to the peak number density just as the halo associated with it is a normal halo contributing to the halo mass function. Therefore, the halo is indistinguishable from one evolved by PA, regardless of its real aggregation history.
But what about the shape and kinematics of haloes? Do they allow one to tell between objects formed by PA and having suffered major mergers? The answer is always the same: major mergers go unnoticed in these properties as well. The proof is similar to that carried out in SVMS for the density profile. Using the one-to-one relations inferred in Section 3 thanks to the results of SMVS, given a halo with any arbitrary aggregation history, one can always think about one peak not only with the appropriate density profile, independent of the ellipsoidal shape of the peak, leading by PA to the spherically averaged profile of the halo (SVMS), but also with a shape also leading to its shape and kinematics. The shape of this putative accreting seed, coincides, by continuity, with that of the peak tracing the halo in the peak formalism (see SVMS) and is thus a normal peak contributing to the peak number density. As a consequence, the halo is indistinguishable from the one that would grow by PA from that seed not only regarding the spherically averaged density profile, but also the shape and kinematics. This result can be summarised by saying that major mergers go unnoticed in virialised objects regarding all those properties.
As discussed in SVMS, the previous conclusion reflects the fact that virialisation is a real relaxation process yielding the memory loss of the halo past history. This obviously affects not only the inner structure but also the shape and kinematics of virialised objects. Certainly, the one-to-one correspondence between virialised haloes and their seeds in PA implies that there is no full memory loss in PA due to virialisation, in agreement with the results of simulations (e.g. Vera-Ciro et al. 2011). However, the fact that it is impossible to reconstruct the seed of any simulated virialised halo by running the simulation backwards is evidence of the existence of a time arrow indicating the direction of virialisation. What makes both aspects consistent is the fact that, as just mentioned, the properties of virialised haloes do not allow one to tell whether they have suffered major mergers or they have evolved by PA. Therefore, we cannot unambiguously determine their initial conditions. In other words, we can reconstruct their putative accreting seeds, but these may or may not be their real (possibly multinode) seeds.
An important corollary of this result is that the present model for the shape and kinematics of virialised objects formed by PA also holds for objects having suffered major mergers. One must simply infer the properties of the virialised object from those of its putative accreting seed, that is the peak tracing its current evolution by smooth accretion.

SUMMARY AND CONCLUSIONS
The model developed in SVMS for the inner structure (spherically averaged density profile) of virialised haloes has been extended to deal with their triaxial shape (eccentricity profiles) and kinematics (velocity dispersion and anisotropy profiles). To do this we have considered the simple scenario of haloes evolving by pure accretion. Under this assumption we have derived the shape and kinematics of the final objects from the shape of their seeds. The reason why dealing with the case of pure accretion is sufficient is that all virialised halos, even those having undergone major mergers, can be seen to arise in this way from the peak which, according to the peak formalism, traces their current evolution. The fact that in pure accretion the final objects keeps the memory of the initial conditions is a consequence of virialisation being achieved, in this case, through shell-crossing with no crossing of particle apocentres. This is not contradictory with the idea that virialisation is a real relaxation process. On the contrary, the fact that major mergers do not leave any particular imprint in the shape and kinematics of virialised objects implies that one cannot be sure about the real aggregation history or, equivalently, about the real seed of any given virialised object.
Applied to virialised haloes in hierarchical cosmologies, the model allows one to make the link between the typical (mean) shape and kinematics of these objects and the power-spectrum of random Gaussian density perturbations. Following a simpler comprehensive approximate approach, we have shown that the overall shape and kinematics predicted by the model are consistent with the results of simulations of CDM cosmologies. In particular, the theoretical anisotropy profile has the same universal form as empirically found, with a substantial scatter due to the variety of peak (halo) shapes. Likewise, the theoretical pseudo phase-space density profile has very approximately a power-law form with the logarithmic slope of −1.875 shown by simulated haloes. This profile is much more robust, that is, it does not essentially depend on the peak (halo) shape.
The present model suggests that the origin of the power-law-like form of the pseudo phase-space density profile of virialised haloes formed by pure accretion is due to the particular evolution of the coarse-grained phase-space density during virialisation. This takes place through the phase mixing produced by shell-crossing in one direction with no crossing of particle apocentres. The fact that PA is always very nearly self-similar and that the velocity dispersion is conserved where a tangential velocity component develops in non-radial PA would explain why the solution found by Bertschinger (1985) for self-similar, radial PA, with the same kind of phase mixing, is very nearly preserved. which can be readily solved for the quantities re(r)/r and r/a(r). Then, the eccentricities ep(r) and es(r) can be obtained by means of equations (A5) and (A6).