11institutetext: Institut d’Astrophysique de Paris, 98bis Boulevard Arago,
Paris

# Equilibrium solution for cold dynamical systems and self-similarity.

Alard  C

###### Abstract

Aims. Numerical simulations demonstrate a link between dynamically cold
initial solutions and self-similarity. However the nature of this link is not
fully understood. Cold initial conditions alone without further symmetry do
not lead to self-similarity.

Methods. Here we show that when the system approaches equilibrium a new
symmetry appears. The combination of this equilibrium symmetry with the cold
symmetry in the initial conditions leads to full self-similarity. As a
consequence for any initially cold system even if the initial spatial
distribution is not self-similar we will observe an evolution towards self-
similarity near equilibrium.

Results. The case of one dimensional systems or spherically symmetric systems
in 3D are discussed in detail. Systems depending on the energy and other
integrals are also considered. The problem of the degeneracy of the self-
similar solutions at equilibrium is tackled. It is shown that very small
perturbations at the center of the system have the ability to break this
degeneracy and lead to the convergence towards a specific auto-similar
solution.

###### Key Words.:

(Cosmology:) dark matter - Cosmology: theory - Methods: numerical

##  1 Introduction

The cold dark matter (CDM) paradigm in cosmology implies that structures
formed from initial conditions with a very small velocity dispersion. High
resolution numerical simulations of CDM structure formation show that a self
similar regime (see Gunn (1977)) appears in the central region of halos and
extends for about two decades (see Ludlow etal. (2001)). The origin of this
self similar regime and the convergence towards a self similar solution near
equilibrium is well established for specific conjectures (see Bertschinger
(1985) Binney (2004), Lancellotti & Kiessling (2001), Alard (2013b), Halle,
Colombi and Peirani (2019), Schulz etal. (2013)). However a fundamental
problem remains. How is it possible that starting from various non self-
similar initial conditions the system evolves after a number of dynamical
times towards self-similarity ? Furthermore since a large class of self-
similar solutions exists, how does the system makes a choice and evolves
towards a specific similarity index ? Here are the questions we will try to
answer in this work.

##  2 The road to equilibrium for a cold system

As the system evolves from the initial conditions the number of orbits made by
a given point increases. For each new orbits a fold develops and as a
consequence the number of folds in phase space increases quickly. Increasing
the number of folds has the obvious consequence to increase the elongation in
phase space. If we consider the evolution of an initial small area in phase
space as the system evolves this area become more and more elongated. Provided
that the area is small enough the phase space density in the elongated area is
nearly constant. When the elongation become large this small area has enough
extension to cover a whole fold. As a consequence the density in phase space
along the fold is nearly constant. Thus this fold represents an iso-contour of
the density in phase space. Let now imagine that we smooth the solution with
an adaptive filter. The local scale of this filter is of the order of the
inter-fold distance. The smoothing operation using this filter would not
produce a continuous density if the elongated fold with nearly constant
density is not nearly closed. As a consequence with increasing elongation the
fold must become nearly closed. A nearly closed fold is basically a fold near
equilibrium. This reasoning indicates that after a sufficient number of orbits
a cold system should reach a state which is close to equilibrium. A process
that could change the evolution of the system and could stop the increase in
the number of folds and the increase of the local elongation would the
appearance of global oscillations. For instance if an asymmetry is present at
the center a stable oscillation may be produced which would make the
convergence towards equilibrium impossible. An example of such process is
illustrated In Fig.’s (1) and (2). The evolution of a system with two distinct
cores in the initial conditions id presented in Fig. (2). In contrast to Fig.
(1) where a single core is present the system in Fig. (2) does not converge
towards a self-similar solution. The situation may have been different in 3D
where the gravity is stronger and thus where the relaxation process is more
efficient in its ability to erase the details of the initial conditions.

Figure 1: The late evolution of a system with cold initial conditions for 1D
gravity. The density profile of the initial conditions in the $x$ dimension is
a power-law with slope $\frac{-1}{2}$. Figure 2: The late evolution of a
system with cold initial conditions for 1D gravity. The density profile of the
initial conditions in the $x$ dimension is two power-laws with slope
$\frac{-1}{2}$ and different centers. Note that the structure of the self-
similar solution of Fig (1) is completely lost. The dynamic in the central
region is dominated by oscillations between the 2 power-law cores. The
relaxation process and the strength of the gravity in 1D is not strong enough
to overcome this oscillation mode.

##  3 Equilibrium as a symmetry

At equilibrium various quantities like the energy or others integrals of
motion are conserved along individual orbits. Thus the density at equilibrium
is a function of the energy and of other the integrals of motion.

###  3.1 Multi-dimensional systems depending on a single spatial and velocity
variable.

In this section we will consider either one-dimensional gravity or 3
dimensional gravity with spherical symmetry. In such case the density in a
phase space with general coordinates $(x,v)$ can be expressed as a function of
the energy and other integrals like the angular momentum.

####  3.1.1 Systems depending only on the energy.

$f(x,v)\equiv F(E)=F\left(\phi(x)+\frac{v^{2}}{2}\right)$ (1)

For convenience we re-scale the velocity axis, then Eq. 1 reads,

$f(x,v)\equiv F(E)=F\left(\phi(x)+v^{2}\right)$ (2)

It is clear that Eq. 2 is invariant under general rotations in the
$(\sqrt{\phi(x)},v)$ space. Let’s now combine this rotational symmetry in the
$(\sqrt{\phi(x)},v)$ space with the cold symmetry in the initial conditions.
In the cold initial conditions the density is a delta function in the velocity
space. The delta function will be un-changed by a re-scaling of the velocity
axis, the only modification due to the re-scaling of the velocity space will
be a scaling of the amplitude of the function. As a consequence a re-scaling
of the velocity will produce the same final density at equilibrium up to a
scaling factor. Now if we include the rotational symmetry at equilibrium, the
phase-space density must also be invariant by general rotations in the
$\left(\sqrt{\phi(x)},v\right)$ space. Let first make a re-scaling in
velocity,

$v\rightarrow\tilde{v}=\alpha v\\\ $ (3)

We also consider the effect of the velocity re-scaling on the potential
$\phi(x)$,

$\phi(x)\rightarrow\tilde{\phi}(x)\\\ $ (4)

By applying the velocity transformation in Eq. (3) to Eq (2) we obtain,

$f(\tilde{\phi}(x)+\alpha^{2}v^{2})=\beta f(\phi(x)+v^{2})=\beta f(E)$ (5)

Let now consider the special case $\phi(x)=0$, then Eq. (5) reads,

$f(\alpha^{2}v^{2})=\beta f(v^{2})$ (6)

Note that we assumed implicitly that $\phi=0$ implies that $\tilde{\phi}=0$.
In 1 dimension we consider initial conditions with symmetry with respect to
the origin. In radial coordinates in muti-dimension the same symmetry is
implicit. The value $\phi=0$ is observed at the center of the system and is a
consequence of the symmetry of the system. The velocity re-scaling that we
apply conserve this symmetry thus $\tilde{\phi}=0$. Let now apply a rotation
in the $\left(\sqrt{\phi(x)},v\right)$ space. Then $(\phi=0,v)$ is transformed
to $(\sqrt{\phi_{1}(x)},v_{1})$, the norm is conserved in the rotation thus,

$\phi_{1}+v_{1}^{2}=v^{2}.$ (7)

By introducing Eq. (7) in Eq. (6) we obtain,

$f\left(\alpha^{2}\left(\phi_{1}(x)+v_{1}^{2}\right)\right)=\beta
f\left(\phi_{1}(x)+v_{1}^{2}\right)$ (8)

Which leads directly to:

$f\left(\alpha^{2}E_{1}\right)=\beta f\left(E_{1}\right)$ (9)

It is clear that (9) express the full self-similarity of the phase space
density at equilibrium. Note that the transformation of $\phi$ due to the
velocity re-scaling can be obtained using Eq. (8). The left side of the
equation express the phase space density in the transformed coordinates
$(\tilde{\phi},\tilde{v})$ By direct identification we find the transformed
potential $\tilde{\phi}$ (see Eq. 4),

$\tilde{\phi}(x)=\alpha^{2}\phi(x)$ (10)

An equivalent transformation is to consider the effect of the velocity re-
scaling on the variable x. Let’s define $\tilde{x}$ as the transformation of
x. Then,

$\tilde{\phi}(x)=\phi\left(\tilde{x}\right)=\alpha^{2}\phi(x)$ (11)

Let’s now consider the Poisson equation for a n dimensional system,

$x^{1-n}\frac{d}{dx}\left(x^{n-1}\frac{d\phi}{dx}\right)=\rho=\int fd^{n}v$
(12)

Let then apply the velocity re-scaling transformation to Eq. (12),

$\tilde{x}^{1-n}\frac{d}{d\tilde{x}}\left(\tilde{x}^{n-1}\frac{d\tilde{\phi}}{d{\tilde{x}}}\right)=\tilde{\rho}=\int\tilde{f}d^{n}\tilde{v}$
(13)

Using Eq’s (11), (9) and (12) we obtain,

$\alpha^{2}\left(\left(\frac{dx}{d\tilde{x}}\right)^{2}\frac{d^{2}\phi}{dx^{2}}+{\tilde{x}}^{1-n}\frac{d}{d\tilde{x}}\left({\tilde{x}}^{n-1}\frac{dx}{d\tilde{x}}\right)\frac{d\phi}{dx}\right)=\beta\alpha^{n}x^{1-n}\frac{d}{dx}\left(x^{n-1}\frac{d\phi}{dx}\right)$
(14)

Eq. (14) has two type of solutions for the transformation $\tilde{x}$, (I)
solutions of the type $\frac{dx}{d\tilde{x}}={\rm Constant}$, or (II) general
solutions which depends on $\phi$. As an illustration for $n=1$ the type (II)
general solution reads:

$\frac{dx}{d\tilde{x}}=\beta\alpha^{n-2}\tanh\left[\log\left(\frac{d\phi}{dx}\right)+C_{1}\right]$
(15)

The problem with the solution of type (II) in Eq. (15) is that for a general
potential $\phi(x)$ the identity condition (in the abscence of re-scaling)
$\tilde{x}(x)=x$ for $\beta=1$ and $\alpha=1$ cannot be universally satisfied.
Thus we are left with only the type (I) solution where $\tilde{x}$ is a linear
function of $x$. This type (I) solution reads,

$\tilde{x}={\beta}^{-\frac{1}{2}}\alpha^{-\frac{n}{2}+1}x=\gamma x$ (16)

Eq. (16) implies the similarity of the potential, using Eq. (11) we find,

$\phi\left(\gamma x\right)=\alpha^{2}\phi(x)$ (17)

Similarly the similarity of $f$ may be re-expressed,

$f\left(\gamma x,\alpha v\right)=\beta f\left(x,v\right)$ (18)

###  3.2 System depending on other integrals.

If the phase space density depends on $E$ and the angular momentum $J$ Eq. (2)
now reads,

$f(x,v)\equiv F(E,J)\equiv F\left(\phi(x)+v^{2},|\bf r\times v|\right)$ (19)

Similarly Eq. (5) now reads,

$f\left(\tilde{\phi}(x)+\alpha^{2}v^{2},\tilde{J}\right)=\beta
f(\phi(x)+v^{2},J)=\beta f(E,J)$ (20)

Starting from here we will apply the reasonning we already used in Sec.
(3.1.1). We note the similarity between Eq. (5) and Eq. (20) and proceed
through equations (6) (7) (8) and (9). By this process we arrive at a first
similarity equation for $f$,

$f\left(\alpha^{2}E,\tilde{J}\right)=\beta f(E,J)$ (21)

The Poisson equation is not modified with respect to Sec. (3.1.1) resulting in
the same results for the similarity of the potential and of the $x$ variable
(Eq’s 17 and 18). Using Eq’s (18) and Eq. (3) we find that angular momentum
re-scacle like:

$\tilde{J}=\alpha\gamma J$ (22)

By combining Eq. (21) and Eq. (22) we arrive at a full similarity for $f$,

$f\left(\alpha^{2}E,\alpha\gamma J\right)=\beta f(E,J)$ (23)

##  4 Convergence towards auto-similarity.

The ideas and calculations developed in Sec. (3) indicates that self-
similarity is a consequence of equilibrium for cold initial conditions.
However self-similar solutions belong to a large class of solutions and the
process leading to the convergence towards a specific solution is not clear.
Here we propose a numerical study of this process for 1D gravity. The
simulations are identical to the 1D simulations produced in Alard (2013a) (see
Sec. 9). A fundamental problem is to asses the number of orbital times
required to reach equilibrium. Since we already know that equilibrium is
equivalent to self-similarity it is also a way to asses the number of orbital
time required to reach self-similarity. A general result obtained from the 1D
simulations is that the convergence towards a stationary regime is obtained
within only one orbital time (see Fig 3). Since the orbital time is a
decreasing function of the distance to the center, the time required for the
convergence towards a stationary regime become infinitely small towards the
center. The energy of the particles on orbits very close to the center is also
very small. This is indeed the case in 1D if the density is a power-law with
exponent less than $-2$. This is also the self-similar solution in 3D observed
in numerical simulations. This solution has an asymptotic power-law slope in
density of $-1$. As a consequence the solution converges first in an area
infinitely close to the center in an infinitely short time and for infinitely
small energies. Once a solution has developed in this area very close to the
center, the degeneracy on the choice of the self-similar solution is broken.
The breaking of the degeneracy itself requires very little energy since it
involves modifying the dynamics of particles with infinitely small energies.
As a consequence the nature of any perturbating process near the center of the
system is very important and will drive the convergence towards a specific
self-similar solution. An example of a process acting in the central region
and leading to a general convergence towards an universal solution is
illustrated in Alard (2013a). In this case if the density in the initial
conditions are less concentrated than a power law of exponent $-\frac{1}{2}$,
a general convergence towards a power-law density with exponent $-\frac{1}{2}$
is observed at equilibrium. The process acting in the central region and
leading to the convergence is due to the density profile of caustics very
close to the center. Another example is the convergence towards an universal
profile for a spherically symmetrical system with constant angular momentum.
This time the convergence is due to the dominating effect of the angluar
momentum in the central region (see Halle, Colombi and Peirani (2019)).

Figure 3: The trajectory in the $(x,v)$ phase space for points within the
central region. The initial density is a power-law profile with exponent
$-0.35$. The presented trajectories start from the initial conditions, the
first orbit is in red color. The value of $X_{0}$ in each graph indicates the
starting position in the initial conditions. In all cases the stationary
regime is reached approximately after one orbital time.

##  5 Conclusion.

This work offers an explanation for the self-similar conjecture observed for
cold initial conditions. Even if the density of these initial conditions is
not self-similar, which would require a power-law density, the system evolves
towards self-similarity. The explanation is that the initial cold-symmetry of
the initial conditions when combined with the symmetry of the equilibrium
state induces a new class of symmetry which is equivalent to self-similarity.
This result is general it was demonstrated in Sec. 3 that self-similarity
occur at equilibrium for cold initial conditions in 1D or spherically
symmetric gravity in 3D. On the other hand if self-similarity is predicted at
equilibrium it does not specify the self-similar solution itself. The self-
similar solutions belong to a large class of solutions and the convergence
towards a specific solution must be induced by some process. An illustration
of this convergence process in 1D gravity is provided in Sec. 4. The
convergence of the solution is operated in a typical orbital time. Since the
orbital time become very small at the center the breaking of the degeneracy
and the convergence towards a specific solution is first operated in the
central region, and then propagates towards the outer regions of the system.
The initial change of the orbit at the center requires only an asymptotically
small change in energy.

## Data Availability

No datasets were generated or analyzed during the current study.

## References

  * Alard (2013a) Alard, C., 2013, MNRAS, 428, 340
  * Alard (2013b) Alard, C., 2013, International workshop Vlasov-Poisson : the numerical approach and its limits, IHP, 2013
  * Binney (2004) Binney J., 2004, MNRAS, 350, 939
  * Bertschinger (1985) Bertschinger, E., 1985, ApJS, 58, 39B
  * Colombi (2017) Colombi, S., Alard, C., 2017, JPlPh, 83, 7002
  * Dehnen & McLaughlin (2005) Dehnen, W., McLaughlin, D., 2005, MNRAS, 363, 1057
  * Fillmore & Goldreich (1984) Fillmore J. A., Goldreich P., 1984, ApJ, 281, 9
  * Gunn (1977) Gunn J. E., 1977, ApJ, 218, 592
  * Halle, Colombi and Peirani (2019) Halle, A., Colombi, S., Peirani, S., 2019, A&A, 621, 8
  * Lancellotti & Kiessling (2001) Lancellotti C., Kiessling M., 2001, ApJ, 549, L93
  * Ludlow etal. (2001) Ludlow A., Navarro J., Springel V., Vogelsberger M., Wang J., White S.,Jenkins A ., Frenk C., 2010, MNRAS, 406, 137
  * Navarro, Frenk & White (1997) Navarro, J., Frenk, C., White, S., 1997, ApJ, 490, 493
  * Schulz etal. (2013) Schulz, A., Dehnen, W., Jungman, G., Tremaine, S.,2013, MNRAS,431,49S