Nikolay N. Fimin
Keldysh Institute of Applied Mathematics of Russian Academy of Sciences, Moscow, Russian Federation.
DOI: 10.4236/jhepgc.2026.121021   PDF    HTML   XML   76 Downloads   249 Views   Citations

Abstract

In this paper, the formation criteria of nonstationary pseudoperiodic structures in a cosmological system of gravitating particles are demonstrated. These criteria are based on the properies of solutions of the system of Vlasov-Poisson equations. Furthermore, the properties of solutions of the linearized integral equation for the generalized gravitational potential, leading to the emergence of coherent structures of the quasi-equilibrium type in stationary systems of massive particles, are investigated.

Keywords

Vlasov-Poisson Equations, Dispersion Relation, Liouville-Gelfand Equation

Share and Cite:

1. Introduction

Description of the genesis of cosmological structures of the wall type and pseudo-1D threads that form the so-called “cosmic web” is currently a highly debated topic among research groups around the world. The prevailing concepts in the analysis of complex issues related to evolution distribution of the scales of galaxies and clusters, Currently, the identification of statistical patterns in sets of observations data and attempts at a rational explanation phenomenological approximations of these patterns. At the same time, the study of the mechanism realization of pseudo-ordering on cosmological scales recognized by default unpromising (primarily due to the lack of verified observational material). The generally accepted assumptions here are the results of hydrodynamic and kinetic modeling in the presence of a significant number of model assumptions relative to the initial data, physical assumptions regarding the influence of dark energy, dark matter (scalar fields, presence of exotic particles, etc.).

To a large extent the starting point for calculations (in particular, for the flat metric space-time close to the modern era) are the intermediate and final stages of model Jeans instability, that is, the collapse of the existing structure primary object with a known distribution of elementary substructures in phase space (or, for a hydrodynamic description, with known set of macroparameters and velocity field) [1]. Obviously, in this case, there is considerable arbitrariness in setting additional conditions, which can lead to significantly different results; the only criterion the truth of the modeling, it seems, can only be served by obtaining a physical picture, coinciding with real modern observations (which explains it) special, even excessive, attention to the analysis of distribution statistics). However, based on last requirement, we cannot free ourselves from the influence of a whole set of factors with a fundamentally nonlinear level of influence in the initial-boundary statement tasks that lead to the currently observed set of galaxy distributions. Thus, we cannot fix it based on the facts stationary picture of modern observation eras, specific additional conditions for mathematical modeling problems cosmic web immediately after the end of the inflationary era (and, naturally, in earlier period); this limitation is aggravated introduction into consideration random local perturbations of the parameters of the astrophysical environment with an additional nonlinear interaction of modes (and this interaction should lead to anisotropic self-ordering in the dynamics of the environment on scales significantly exceeding the scale of fluctuations).

The author believes that a special role in the study of the structure and evolution of cosmological structures must acquire approaches, based on the deterministically determined emergence of large-scale structures (coherent in the sense of having translational invariance of substructures in the structure itself, or when compared to similar ones existing nearby) in a cosmological background. That is, “self-assembly” such structures as a result of interaction small-scale disturbances, creating the “order from chaos” (according to the scenario of the formation of “pancakes” by Ya.B. Zeldovich) [2] seems unlikely, and the main ongoing process of structural genesis are the primary formation of a large structure with its further evolution in the form of a transition between states of relative equilibrium (which correspond to the extrema of the entropy of the system with an adiabatic change in its other thermodynamic or topological parameters [3] [4]); system can be multi-connected, exactly like the observed “cosmic web”. In this paper, the author demonstrate the possibility of natural construction of models cosmological structures based on a qualitative analysis of the properties of solutions of the equations of the Vlasov-Poisson system taking into account the modification of the form of the gravitational potential, with the inclusion of an additional member, including cosmological parameters. The possibility of the emergence of multi-connected structure for a flat system gravitating substructures, shown the difference between the kinetic calculation and the hydrodynamic one for the generalized Jeans decay, using methods of bifurcation theory solutions of integral equations the possibility of the appearance of secondary structures between nodes has been demonstrated “cosmic web”.

2. The Vlasov-Poisson System of Equations with a Modified Gravitational Interaction

The Vlasov-Poisson system of equations for describing the dynamics of a multiparticle system (taking into account the Friedmann expansion of the Universe) in $d$ -dimensional ($d=2,3$ ) particles of equal masses $m$ in accordance with the source-like representation of the self-consistent gravitational field proposed in the works [5] [6] is written in the following form:

$\frac{\partial F\left(x,v,t\right)}{\partial t}+{\text{div}}_x\left(vF\right)-\frac{\partial F}{\partial v}{\nabla }_x\Phi \left(F\right)=0,$ (1)

$\Phi \left(F\right)\equiv {\displaystyle {\int }_{{\omega }_{x^{\prime }}}{\displaystyle {\int }_{{\omega }_{v^{\prime }}}{Y}_d\left(x-x^{\prime }\right)F\left(x^{\prime },v^{\prime },t_{*}\right)\text{d}{x}^{\prime }\text{d}{v}^{\prime }}},$

where $Y_d\left(r\right)~Y_d^{\left(1\right)}\left(r\right)+Y_d^{\left(2\right)}\cdot r^2+Y_d^{\left(3\right)}$ (this form corresponds to the simultaneous attraction and repulsion between particles); is an integrable (and twice differentiable) function in some bounded domain ${\omega }_d={\omega }_x\subset {ℝ}^d$ , $Y_d^{\left(2\right)}\in {ℝ}^1$ —constant. It can be shown that for an arbitrary finite instant of time $t_{\text{*}}\in {ℝ}_{+}^1$ has place modified Poisson equation:

${{\Delta }_x^{\left(d\right)}\Phi \left(F\right)|}_{t=t_{\ast },\forall t_{\ast }\in {ℝ}_{+}^1}=\frac{2{\pi }^{d/2}}{\Gamma \left(d/2\right)}{\gamma }_d\rho \left(x\right)-C_1,$ (2)

$\rho \left(x\right)~m{\displaystyle \int F\left(x,v,t_{*}\right)\text{d}v}{\gamma }_3\equiv \gamma ,C_1=c^2\Lambda .$

Equation (2) takes the form of the so-called (non-homogeneous) Liouville-Gelfand equation (LGE) [7], if we take the density distribution of particles equal to that in the case in case of integration Maxwell-Boltzmann distribution functions over velocities (here and below we do not take into account the dependence equilibrium distribution function from integrals of motion other than energy): $F=F_0\left(x,v,t_{*}\right)\propto \Omega \left(\epsilon \left[v\left(t_{*}\right);\Phi \left(x\right)\right]\right)$ ; (in this case, the regular branch Φ is selected).

The LGE for the gravitational potential in cosmological systems is essentially local (since it is a differential equation), and with its help it is rather difficult to describe global processes in a self-consistent field for the analysis of the evolution of the entire system of particles (substructures). In the study of non-local effects It is appropriate to consider the integral version of the equation for the potential. However, a natural limitation for the application of this form in analytical and numerical calculations there is an obvious need to take into account the boundaries of the region containing the system of particles under study, on the value of the potential at a given point. Of course, this is due to the fact that for a system of particles interacting according to the modified Newton theorem, it is correct to put on the boundary of the studied region the zero Dirichlet condition ${\Phi \left(x\right)|}_{\partial \omega }=0$ is very difficult (for the general case of asymmetric of a multiparticle system this is unrealistic).

Let us consider the following statement of a physical problem: in the Friedmann flow there are regions of increased density, which are consequences of fluctuation disturbances early in the development of the Universe; could the gravitational interaction between these regions ${\omega }_j$ (each of which contains, perhaps, hierarchy of substructures of various sizes) or in the vicinity each region to form a (quasi)stationary structure that has signs of orderliness? A significant difference compared to charged plasma with Debye screening, allowing pseudo-homogeneous equilibrium distributions of particles is that cosmological substructures in the process of implementation from primary disturbances (local increases) of density macroscale objects have spherically symmetric ($d=3$ ) or radially symmetric ($d=2$ ) forms (in accordance with the Gi-Nidas-Nirenberg theorem [7]); in accordance with the structure of the interparticle potential in the neighborhood spheres with increased density, a radially symmetric layer appears, in in which the attraction to the point which is the center of the said sphere prevails, and semi-infinite spherical layer with repulsion dominance (for a system of two spheres, containing matter with increased density, there is competition between these two centers, leading either to the destruction of one of them, or, in connection with Friedman’s expansion of space, the emergence of repulsion zones between them).

Two approaches can be considered when finding the gravitational potential in the vicinity of a region containing a system of interacting particles:

1) solution of the equation for the potential (non-homogeneous Liouville-Gelfand one) under given Dirichlet conditions on some a priori determined (having considerable arbitrariness in its choice), the boundary of the region $\partial {\omega }_j$ (this boundary must include only one local density concentration); in this case, both the internal and external Dirichlet problems can be considered (which corresponds to attempts to establish the values of the potential—hence, the density of particles—respectively, inside and outside the boundary $\partial {\omega }_j$ ). The influence potentials of neighboring regions can be neglected under certain conditions. The data at the boundaries must be consistent (in the sense of continuity of the solution) with the obtained decision;

2) a similar situation can also be considered for the Neumann problem; however, here it is possible to determine quite accurately the boundary of the region $\partial {\omega }_j$ due to the presence of a maximum of the interaction potential (at this point the derivative of the potential is cancelled out). However, with such a formulation of the problem, the question naturally arises about the validity of the region of description of the particle dynamics only by the zone of attraction of the potential.

The physical aspects of the problem being solved obviously determine the method of setting additional conditions for the nonlinear potential equation. Here, we will limit ourselves to using the Dirichlet conditions (due to the theorem on the equivalence of the gravitational field of a sphere and a point at its center). Thus, we will be interested in whether $\omega$ arise inside (or outside) a fixed region for a given the value of the potential at the boundary of the secondary solutions of the LGE, possessing the property of (quasi)periodicity.

The solution of the Dirichlet problem for the Poisson equation can be represented through the Green’s function of elliptic operator and its meaning potential at the boundaries. In explicit form, this solution has the form:

$\begin{array}{c}\Phi \left(x\right)=-{\displaystyle {\int }_{\omega }\mathcal{G}\left(x,x^{\prime }\right)\left({\gamma }_d\rho \left(x^{\prime }\right)-\frac{c^2\Lambda }{4\pi }\right)\text{d}{x}^{\prime }}\\ -\frac{1}{4\pi }{\displaystyle {\int }_{\partial {\omega }^{\prime }}{\Phi |}_{\partial \omega }n\cdot {\frac{\partial }{\partial x^{\prime }}\mathcal{G}\left(x,x^{\prime }\right)|}_{\partial \omega }\text{d}{S}^{\prime }},\end{array}$ (3)

where $\mathcal{G}\left(x,x^{\prime }\right)$ is the Green’s function of the Dirichlet problem for the inhomogeneous Poisson equation. Choosing a ball $\omega =\left\{\tilde{x};0\le \left|\tilde{x}\right|\le ℛ\right\}$ as the calculation domain, possible in the case of particle distribution close to isotropic assume that the values of the potential on its boundary in accordance with the modified Newton theorem [8] will be equal to ${\Phi |}_{\partial \omega }=-{\gamma }_{d}M/ℛ-c^2\Lambda {ℛ}^2/6$ ($M=Nm$ ). Replacing $\rho \left(x\right)={\left(4\pi \right)}^{-1}{\varkappa }_d\left(T\right){\gamma }_d\exp \left(-\Phi /T\right)$ , and defining the Green’s function for the Dirichlet domain under consideration is:

$\mathcal{G}\left(x,x^{\prime }\right)\equiv 4\pi {\displaystyle \sum }_{\ell =0}^{\infty }{\displaystyle \sum }_{m=-\ell }^{\ell }\frac{Y_{\ell m}^{*}\left({\theta }^{\prime },\phi \right)Y_{\ell m}\left(\theta ,\phi \right)}{2\ell +1}\frac{x_{<}^{\ell }{x}_{>}^{\ell }}{{ℛ}^{2\ell +1}},$

$x_{<}=\min \left(\left|x\right|,\left|x^{\prime }\right|\right),x_{>}=\max \left(\left|x\right|,\left|x^{\prime }\right|\right),$

we obtain a nonlinear integral equation for the potential $\Phi \left(x\right)$ (for $d=3$ ):

$\begin{array}{l}\Phi \left(x\right)=-{\varkappa }_d\left(T\right){\gamma }_d{\displaystyle {\int }_{{\omega }^{\prime }}\left(\frac{1}{\left|x-x^{\prime }\right|}-\mathcal{G}\left(x,x^{\prime }\right)\right)\exp \left(-\Phi \left(x^{\prime }\right)/T\right)\text{d}{x}^{\prime }}\\ -\frac{c^2\Lambda }{6}{x}^2+C_0,C_0=-\frac{{\gamma }_{d}M}{ℛ}-\frac{c^2\Lambda {ℛ}^2}{6}.\end{array}$

After obvious transformations, it can be written as an inhomogeneous equation of the Hammerstein type for the dimensionless potential $U_H$ :

$U_H\left(x\right)=\widehat{\mathfrak{W}}\left(U_H\left(x\right)\right),$ (4)

$\begin{array}{c}\widehat{\mathfrak{W}}\left(U_H\left(x\right)\right)\equiv {\lambda }_H\left(T\right){\displaystyle {\int }_{{\omega }^{\prime }}\underset{\mathcal{K}\left(x,x^{\prime }\right)}{\underbrace{\left(\frac{1}{\left|x-x^{\prime }\right|}-\mathcal{G}\left(x,x^{\prime }\right)\right)}}\exp \left(-U_H\left(x^{\prime }\right)\right)\text{d}{x}^{\prime }}\\ +\alpha \left(\Lambda ,T\right){\left|x\right|}^2;\end{array}$

${\lambda }_H\left(T\right)\equiv \frac{-{\varkappa }_3\left(T\right){\gamma }_3}{T}\exp \left(-\frac{C_0}{T}\right),\alpha \left(\Lambda ,T\right)=-\frac{c^2\Lambda }{6T},U_H=\frac{\Phi -C_0}{T},$

It should be noted that for a spherically symmetric density distribution

${\displaystyle {\int }_{{\omega }^{\prime }}\mathcal{G}\left(x,x^{\prime }\right)\rho \left(x^{\prime }\right)\text{d}{x}^{\prime }}\to C_1\left(=\text{const}\right).$

3. Qualitative Analysis of the Integral Equation for the Potential

Equation (4) contains very complete information about the behavior of the cosmological dynamics of a multiparticle system with gravitational interaction under consideration. First of all, it should be pointed out that the obvious absence of a solution for this equation for $d=3$ in the entire area of research in the form of a constant potential (in contrast to the classical Newtonian case, for which gravity is equivalent attraction at $\Lambda \equiv 0$ ). Thus, in this case, the global uniform distribution of matter, as a condition for the consistency of the modified Poisson equation and the Vlasov equation, can be considered only as some abstract approximation with the aim of further continuation of the solution of a given system of equations with respect to a parameter.

Let us consider the linearized (near the point $U_H\left(0\right)=U_0$ ) version of Equation (4) in the homogeneous and inhomogeneous cases:

$U^{†}={{\widehat{\mathfrak{W}}}^{\prime }}_0\left(U^{†}\right),U^{‡}={{\widehat{\mathfrak{W}}}^{\prime }}_0\left(U^{‡}\right)+\alpha {\left|x\right|}^2.$(5)

The operator $\widehat{I}-{{\widehat{\mathfrak{W}}}^{\prime }}_0$ , where ${{\widehat{\mathfrak{W}}}^{\prime }}_0\equiv {\widehat{\mathfrak{W}}}^{\prime }\left[U_0\right]$ is the Frechet derivative of the integral operator Hammerstein from the right side of the first equation belongs class of Noetherian operators of zero index with a weak singularity and an alternating kernel; since the Green’s function is symmetric with respect to its coordinates, for matter distributions that deviate weakly from spherically symmetric, ${{\widehat{\mathfrak{W}}}^{\prime }}_0$ can also be considered self-adjoint. For small amplitudes of deviation $\delta U^{†}$ from the selected solution, we get the opportunity to apply the well-known mathematical apparatus of Fredholm operator analysis to the homogeneous linear equation $\delta U^{†}-\lambda {\displaystyle {\int }_{{\omega }^{\prime }}\mathcal{K}\left(x-x^{\prime }\right)\text{exp}\left(-U_0^{†}\right)\delta U^{†}\text{d}{x}^{\prime }}=0$ . For simplicity of calculations, without loss of generality, we can take $U_0^{†}=0$ , and look for periodic solutions of the last equation (in the form of an expansion in eigenfunctions ${\mathfrak{b}}_j=c_{\omega }\exp \left(iqx\right)$ of the kernel $\mathcal{K}$ in the domain $\omega \subseteq {ℝ}^3$ ) in the form $\delta U^{†}={\displaystyle {\sum }_j{\mathfrak{a}}_j\left(c_{\omega }\right)\text{exp}\left(iqx\right)}$ ($q_{\ell =1,2,3}=2\pi /d$ , the cubic case is obviously generalized to $d_{\ell }\ne d_k$ ). Substituting this expression into the first equation of the set (5) gives

$1=\lambda {\displaystyle {\int }_{{\omega }^{\prime }}\mathcal{K}}\left(x-x^{\prime }\right)\text{exp}\left(-iq\left(x-x^{\prime }\right)\right)\text{d}\left|x-x^{\prime }\right|.$ (6)

We introduce the critical value of the parameter $\lambda ={\lambda }_c$ (corresponding to the case $d\to \infty$ , $q_{\ell }\equiv 0$ ), with its help we can write criterion for the existence of periodic solutions for the linearized integral homogeneous Poisson equation:

$\lambda ={\left({\displaystyle {\int }_{\omega }\mathcal{K}\left(r\right)\frac{\sin \left(qr\right)}{qr}{r}^2\text{d}r\text{d}\theta \text{d}\varphi }\right)}^{-1}\ge {\lambda }_c\equiv {\left({\displaystyle {\int }_{\omega }\mathcal{K}\left(r\right)r^2\sin \left(\theta \right)\text{d}r\text{d}\theta \text{d}\varphi }\right)}^{-1}$ (7)

Obviously, this criterion is only suitable for the case $\Lambda \equiv 0$ , and only for the 2–partial approximation of interaction. This implies, that the root $q$ of equation (6) is unique (in this case the potential distribution will be purely periodic); if there are several (incommensurate) roots of $q_s$ , the potential distribution will belong to the class of almost-periodic functions. Taking into account collective functional interactions in a system of $N$ particles leads to a (homogeneous) equation of the form

$\begin{array}{c}\delta U^{†}\left(x\right)={\displaystyle \sum _{k=1,\cdots ,N}{\displaystyle {\int }_{{\omega }_1}\cdots }{\displaystyle {\int }_{{\omega }_k}{\mathcal{K}}_k}}\left(x,x_1,\cdots ,x_k\right)\\ \times \left(\text{exp}\left(-\delta U^{†}\left(x_1\right)\right)\cdots \text{exp}\left(-\delta U^{†}\left(x_k\right)\right)\right)\text{d}{x}_1\cdots \text{d}{x}_k.\end{array}$

After linearizing this equation we obtain

$\delta U^{†}\left(x\right)={\displaystyle \sum _k{\lambda }_k}{\displaystyle {\int }_{{\omega }_1}\cdots }{\displaystyle {\int }_{{\omega }_k}{\mathcal{K}}_k}\left(x,x_1,\cdots ,x_k\right){\displaystyle \sum _{\ell }^k\delta U^{†}}\left(x_{\ell }\right){\displaystyle \prod _{s=1}^k\text{d}{x}_s}.$

The corresponding to (7) criterion for the emergence of three-dimensional periodic solutions can be represented in the following form:

$\begin{array}{l}{\displaystyle {\int }_{\omega }{\displaystyle \sum _k{\displaystyle \sum _s^k{\lambda }_k}}}\frac{\sin \left(qr\right)}{qr}\left({\displaystyle {\int }_{{\omega }_1}\cdots }{\displaystyle {\int }_{{\omega }_{s-1}}{\displaystyle {\int }_{{\omega }_{s+1}}\cdots }}{\displaystyle {\int }_{{\omega }_k}{\mathcal{K}}_s}\left(x,x_1,\cdots ,x_k\right)\frac{{{\displaystyle \prod }}_{n=1}^k\text{d}{x}_n}{\text{d}{x}_s}\right)\\ \times r^2\sin \left(\theta \right)\text{d}r\text{d}\theta \text{d}\varphi =1.\end{array}$ (8)

Therefore, taking into account cluster interactions in a multiparticle system is consistent with taking into account two-particle interactions: the emergence of a periodic distribution structure potential (and density of matter in space) in space (in a linear homogeneous approximation) occurs abruptly upon reaching a certain (quasi-equilibrium kinetic [5]) temperature in the system under study. The periods of the structures are determined from conditions (6) and (8).

For a non-homogeneous linear equation (the second equation (5)) one can, using the Hilbert-Schmidt theorem, obtain an explicit form (unique) solutions in the form of a resolvent series (uniformly and absolutely convergent Fourier series in the eigenfunctions of the kernel $\mathcal{K}$ ), for $\lambda \ne {\lambda }_{\ell }$ ($\ell =1,2,\cdots$ ):

$\delta U^{‡}\left(x\right)=\lambda {\displaystyle {\int }_{\omega }{\displaystyle \sum _{\ell =1}^{\infty }\frac{{\mathfrak{b}}_{\ell }\left(x\right){\mathfrak{b}}_{\ell }\left(x^{\prime }\right)}{{\lambda }_{\ell }-\lambda }\alpha {\left|x^{\prime }\right|}^2\text{d}{x}^{\prime }}}+\alpha {\left|x\right|}^2,$

where ${\lambda }_{\ell }$ are the characteristic numbers of the homogeneous equation (corresponding to the eigenfunctions ${\mathfrak{b}}_{\ell }\left(x\right)$ ). Thus, the periodicity of the structure, when taking into account the repulsion (i.e., when $\Lambda \ne 0$ ) degenerates—in the simplest case—into a composition of a sinusoidal function and an increasing branches of the parabola, which leads to a smoothing of the periodicity and the actual dominance of repulsion at a certain distance (in the interaction channel between two previously formed heterogeneities in a pseudo-uniform distribution of matter); since this process can be considered as a pairwise interaction, the mechanism of formation becomes clear voids in the initially uniformly distributed material continuum: the anisotropy of the directions of the particle momenta in the channels between the inhomogeneities I and II forms non-cubic lattice ($d_{1,2,3}=d$ ), and quasi-one-dimensional layers ($d_3\ll d_{1,2}$ ), the velocity profile of which in the mentioned channel is slowed down for initially fast particles (and thus synchronized as a uniform distribution in velocity space) due to repulsion in the far zone the second component of the macrosystem (heterogeneity II at the other end of the channel). It can be assumed that a significant part of the mass heterogeneity turns into a flat formation corresponding to the first minimum potential, and this situation in the channel of pairwise interaction is mirrored: repulsion, understood as an external third-party force, from the side II this leads to quasi-stationarity of the wall, formed by the inhomogeneity I, and vice versa; Friedmann’s expansion of space leads to the fact that the initially attracted massive inhomogeneities appear at distances where repulsion becomes more significant due to the influence of the cosmological term.

4. The Problem of Initial Perturbation for the Linearized Integral Version of the Vlasov-Poisson Equation

Above, we considered the equation for the potential from the Vlasov-Poisson system of equations. The second equation in this system yields even more interesting results. Let’s transform the equation under study from an integro-differential equation to an integral equation. To do this, we represent its solution in the following form: $F\left(x,v,t\right)=F_0\left(x,v\right)+f\left(x,v,t\right)$ , $\left|F_0\right|\gg \left|f\right|$ . Substituting this expression into the Vlasov equation, we obtain:

$\frac{\partial f}{\partial t}+v\frac{\partial f}{\partial x}=m^{-1}\frac{\partial F_0\left(x,v\right)}{\partial v}\cdot \frac{\partial }{\partial x}{\displaystyle \iint {\mathfrak{Y}}_3^{\left(\vartheta \right)}\left(\left|x-x^{\prime }\right|\right)f\left(x^{\prime },v^{\prime },t\right)\text{d}{v}^{\prime }\text{d}{x}^{\prime }}.$ (9)

Using the inverse of the total derivative operator (or the theory of semigroups of linear operators), we arrive at the form of the Vlasov equation in the representation of shifts along trajectories:

$\begin{array}{l}f\left(x,v,t\right)=f\left(x-v\left(t-t_0\right),v,t_0\right)\\ +{\displaystyle {\int }_{t_0}^t\text{d}{t}^{\prime }}\tilde{B}{\nabla }_v\exp \left(-\frac{m{v}^2}{2T}-\frac{1}{T}{\Phi }_{-}\left(t-t^{\prime }\right)\right)\cdot {\widehat{\mathcal{K}}}_d^{(-)}\left[\vartheta \right]\varrho ,\end{array}$

${\Phi }_{-}\left(t-t^{\prime }\right)\equiv {\Phi }_{\min }\left(x-v\left(t-t^{\prime }\right)\right),$

${\widehat{\mathcal{K}}}_d^{(-)}\left[\vartheta \right]\varrho \equiv \frac{\partial }{\partial x}{\displaystyle \int {\mathfrak{Y}}_d^{\left(\vartheta \right)}\left(\left|x-x^{\prime }-v\left(t-t^{\prime }\right)\right|\right)\varrho \left(x^{\prime },t^{\prime }\right)\text{d}{x}^{\prime }}.$(10)

Here $\varrho \left(x^{\prime },t^{\prime }\right)\equiv {\displaystyle \int f\left(x^{\prime },v^{\prime },t^{\prime }\right)\text{d}{v}^{\prime }}$ , $t_0$ is some initial instant, selected on the time axis, for example, by the occurrence of some perturbation in a system described using retarded potentials, $\tilde{B}=B/m$ . To simplify the calculations, we assume the temperature to be constant in the part of the system that we are interested in.

We will seek the solution $f\left(x,v,t\right)$ in the form of Fourier-Laplace integrals of spatial harmonic waves with decreasing amplitude (or waves with gaps in the coefficients), that is, in the form ${\displaystyle {\int }_{{ℝ}_w^3}{f}_w\left(v,t\right)\text{exp}\left(iwx\right)\text{d}w}$ , where $w={\left(w_{\ell }\right)}_{\ell =1,2,3}^{\text{T}}$ , $w_{\ell }=w_{\ell }^{\left(Re\right)}+i{w}_{\ell }^{\left(Im\right)}$ (respectively, for the density $\varrho \left(x,t\right)$ we have $\varrho \left(x,t\right)={\displaystyle {\int }_{{ℝ}_w^3}{\varrho }_w\left(t\right)\text{exp}\left(iwx\right)\text{d}w}$ . Thus, the kernel on the right-hand side of the definition of ${\widehat{\mathcal{K}}}_3\varrho$ can, together with the exponentially decreasing part of the density, turn into an integrable (in a limited region of space ${ℝ}^3$ due to the aforementioned “problem of infinite mass”) function, obviously, under the fulfillment of the appropriate conditions on the imaginary exponents of the exponential function $w_{\ell }^{\left(Im\right)}$ : the elements $f_w$ and ${\varrho }_w$ belong to the set of coefficients of the generalized Fourier integral, and the conditions for the existence of these integrals are the a priori absence of exponential growth of the distribution function $f\left(x,v,t\right)$ .

Before moving on to the Fourier transforms, we must integrate expression (10) over the velocities (to form the equation for the density). On the left-hand side is the Fourier coefficient ${\varrho }_w\left(t\right)$ , the first term on the right-hand side of the resulting integral equation has the form:

${}^{\left[1\right]}\mathcal{D}{f}_w\left(t\right)={\displaystyle \int f_w\left(v,t_0\right)\text{exp}\left(-iwv\left(t-t_0\right)\right)\text{d}v},$

the second term on the right-hand side of (10) can be transformed to the form ${\displaystyle {\int }_{t_0}^t{}^{\left[2\right]}{\mathcal{D}}_w\left(t-t^{\prime }\right)\text{d}{t}^{\prime }}$ , where:

$\begin{array}{l}{}^{\left[2\right]}{\mathcal{D}}_w\left(t-t^{\prime }\right)\\ ={\displaystyle {\int }_{{ℝ}_v^3}\text{d}v}{\left[\gamma \tilde{B}{\nabla }_v\exp \left(-\frac{m{v}^2}{2T}-\frac{{\Phi }_{-}\left(t-t^{\prime }\right)}{T}\right)\right]}_1\exp \left(iw\left(x-v\left(t-t^{\prime }\right)\right)\right)\\ \times \frac{\partial }{\partial x}{\left[{\displaystyle \int \text{d}{x}^{\prime }}{\left|x-v\left(t-t^{\prime }\right)-x^{\prime }\right|}^{-1}\text{exp}\left(-iw\left(x-v\left(t-t^{\prime }\right)-x^{\prime }\right)\right)\varrho \left(x^{\prime },t^{\prime }\right)\right]}_2.\end{array}$

It should be noted that in the vicinity of point $\Phi ={\overline{\Phi }}_{\min }$ on the phase plane, the fixed branch of the minimum gravitational potential function is close to constant or has a locally parabolic structure (which corresponds to the inertial motion of the particle or motion near a state of stable or unstable equilibrium, corresponding to a generalized Lagrange point).

We rename $\left|x-v\left(t-t^{\prime }\right)-x^{\prime }\right|=\xi$ and integrate in brackets ${\left[\mathrm{...}\right]}_2$ :

$4\pi {\displaystyle {\int }_0^{\infty }\text{exp}\left(w^{Im}\xi \right)\frac{\text{sin}\left(w^{Re}\xi \right)}{w^{Re}\xi }\xi \text{d}\xi }\equiv \mathscr{H}\left(w^{Re},w^{Im}\right).$

Therefore, we finally have:

${}^{\left[2\right]}{\mathcal{D}}_w\left(t-t^{\prime }\right)=\mathscr{H}\left(w^{Re},w^{Im}\right){\displaystyle {\int }_{{ℝ}_v^3}\text{d}v}\exp \left(iw\left(x-v\left(t-t^{\prime }\right)\right)\right)i\left(w\cdot {\left[\mathrm{...}\right]}_1\right)$

Thus, to obtain the Fourier transforms of the density ${\varrho }_w$ , the following equation is valid, the Volterra integral equation of the second kind:

${\varrho }_w\left(t\right)={}^{\left[1\right]}\mathcal{D}{f}_w\left(t\right)+{\displaystyle {\int }_{t_0}^t{}^{\left[2\right]}{\mathcal{D}}_w\left(t-t^{\prime }\right){\varrho }_w\left(t^{\prime }\right)\text{d}{t}^{\prime }}.$ (11)

The solution to this type of equation can be obtained in terms of the one-sided Laplace transform (with respect to the time variable):

$\widehat{ℒ}{f}_w\left(\omega \right)={\displaystyle {\int }_0^{\infty }{f}_w\left(t\right)\text{exp}\left(-\omega t\right)\text{d}t},$

$f_w\left(t\right)=Heav\left(t\right){\left(2\pi i\right)}^{-1}{\displaystyle {\int }_{{\zeta }_0-i\infty }^{{\zeta }_0+i\infty }\widehat{ℒ}{f}_w\left(\omega \right)\text{exp}\left(\omega t\right)\text{d}\omega },$

${\varrho }_w\left(t\right)={\left(2\pi i\right)}^{-1}{\displaystyle {\int }_{{\zeta }_0-i\infty }^{{\zeta }_0+i\infty }\text{exp}\left(\omega t\right)}\frac{\widehat{ℒ}\left({}^{\left[1\right]}\mathcal{D}{f}_w\right)\left(\omega \right)}{1-\widehat{ℒ}\left({}^{\left[2\right]}{\mathcal{D}}_w\right)\left(\omega \right)}\text{d}\omega ,$(12)

$\widehat{ℒ}\left({}^{\left[k\right]}{\mathcal{D}}_w\right)\left(\omega \right)={\displaystyle {\int }_0^{\infty }\text{exp}{\left(-\omega t\right)}^{\left[k\right]}{\mathcal{D}}_w\left(t\right)\text{d}t}$

(the integral is taken in the right half-plane of the complex plane $\omega =-\zeta -i\infty$ along a line parallel to the $\Re \left(\omega \right)=0$ axis, $Heav\left(t\right)$ is a Heaviside function). The poles of the integrand are determined by the equation $\widehat{ℒ}{(}^{\left[2\right]}{\mathcal{D}}_w)=1$ , which is reduced to the form (for $\Re \left(\omega \right)>0$ ):

$\gamma \mathscr{H}\left(w_{\ell }^{Re},w_{\ell }^{Im}\right){\displaystyle {\int }_{{ℝ}^3}\frac{iw{\nabla }_v{F}_0}{\omega +iwv}}=-1.$

This dispersion equation in the general case, of course, can only be solved numerically; However, in limiting cases, the dependence $\omega =\omega \left(w\right)$ is fairly easy to determine, using asymptotic estimates. In particular, for small values of the modulus of the “wave vector” it can be expanded in a Taylor series in $\left|w\right|$ (up to the first non-vanishing term) the factor in the integrand ($\omega /\left(\omega +ivw\right)\approx w/\omega -iv{w}^2/{\omega }^2$ ), and choose a coordinate system with the abscissa axis aligned with the vector $w$ (to obtain a scalar dependence), then:

$\begin{array}{l}\gamma \mathscr{H}\left(w\right)\frac{w^2}{{\omega }^2}{\displaystyle \int [B{T}^{-1}\left(-m{v}_1+{{\Phi }^{\prime }}_{\min }\cdot \left(t-t^{\prime }\right)\right)\begin{array}{c}\\ \end{array}}\\ \times {\exp \left(-\frac{m{v}_1^2}{2T}-\frac{{\Phi }_{-}\left(t-t^{\prime }\right)}{T}\right)]}_3\text{d}{v}_1=1.\end{array}$

Thus, we obtain an approximate dispersion relation (integrated over time), which can be briefly written as

${\omega }^2=\gamma w^2\mathscr{H}\left(w\right){\displaystyle \int {\left[\mathrm{...}\right]}_3\text{d}{v}_1}.$ (13)

The quantity $\mathscr{H}\left(w^{Re},w^{Im}\right)>0$ , so the nature of the system’s dynamics depends on the sign of the expression in brackets ${\left[\mathrm{...}\right]}_3$ : if ${\left[\mathrm{...}\right]}_3<0$ , then $\omega \in ℂ$ (a purely imaginary number); otherwise, $\omega$ is real. This means that a perturbation in our system under study, qualitatively following the law $~\text{exp}\left(\omega t+iwv\right)$ , will have an oscillatory character in time if the system, located in the vicinity of a conditional equilibrium state, is subjected at the initial instant of time to a perturbation with a long wavelength and criterion (13) is satisfied for imaginary $\omega$ ; for real frequencies, the primary disturbance is damped or growing (unstable) depending on the sign of $\omega \lessgtr 0$ . Furthermore, the integrand can change sign, transitioning from oscillatory behavior of the system to oscillatory (this is influenced by both an increase in the observation time interval and a transition to a different branch of the multivalued function ${\Phi }_{\min }$ ).

Let us turn to a particular solution of the Vlasov-Poisson equation (for the 3-dimensional case) of the form

$f^{+}\left(x,v,t\right)~c\left(k\right)f_w^{†}\left(v\right)\text{exp}\left(i\omega t-iwx\right).$

As a simplifying assumption, let’s assume that in the expression

$F_0\left(x,v\right)=B\exp \left(-\frac{m{v}^2}{2T}-\frac{1}{T}{\Phi }_{\min }\left(x\right)\right)$

the value ${\Phi }_{\min }\approx \phi =\text{const}$ (which is due to the system’s state being close to the COP/conditional equilibrium, in which by definition ${{\Phi }^{\prime }}_{\min }=0$ ); we have not yet used the delay approach here—this will happen a little later, so we can locally accept this assumption. The substitution of the ansatz into the equation yields

$f_w^{†}\left(v\right)=\frac{\gamma }{m}\mathscr{H}\left(w\right)\frac{\left(w{\nabla }_v\mathfrak{M}\left(v\right)\right)}{wv-\omega }{\displaystyle \int f_w^{†}\left(v\right)\text{d}{v}_1\text{d}{v}_2\text{d}{v}_3},$

$\mathfrak{M}\left(v\right)=B\exp \left(-\phi /T\right)\text{exp}\left(-m{v}^2/2\right).$

The condition for the existence of a nontrivial solution in the class of functions of the form $f^{+}$ is obtained after integrating over the velocities of both parts (dispersion relation $\omega =\omega \left(w\right)$ ):

$\gamma m^{-1}\mathscr{H}\left(w\right)B\exp \left(-\phi /T\right){\displaystyle \int \frac{\left(wgra{d}_v\right)}{wv-\omega }\text{d}v}=1.$

The magnitudes of the components of the wave vector $w$ are complex, and the integral is taken in the sense of the principal value symmetrically about the $v_1$ -axis, if the $w$ -axis is directed along the radial coordinate axis; The sum of the two terms when going around the pole will correspond to the retarded and advanced potentials. From this we have:

$f_w^{†}\left(v\right)=c\left(w\right)\gamma m^{-1}\mathscr{H}\left(w\right)B\exp \left(-\phi /T\right)\left(w{\nabla }_v/\left(\omega -wv\right)\right).$

Thus, we have obtained the $\omega \left(w\right)$ dependence for a system of massive particles with self-gravity. It provides information about the nature of the system’s dynamics (oscillatory in time and/or space with decreasing or constant amplitude, exponentially decreasing and unstable). Of extreme importance for the astrophysical background of dispersion solution analysis is the fact that in this case there is no need to consider a set of random fluctuations in the system (which, in turn, somehow coherently self-order over colossal distances in one or two directions to form the observed cosmological structures). The mechanism for the emergence of a pseudo-oscillatory (in time and in selected spatial directions) realization of the matter distribution is inherent in the nature of the system of Vlasov-Poisson equations.

5. Conclusion

The paper proposes a kinetic model the emergence of periodicity of gravitational strata—voids separated by two-dimensional surfaces, due to the presence of Poisson structures in the Hubble one-dimensional flow, associated with quasi-oscillatory equation that is a consequence of Poisson’s own equation. This approach is a development of models the emergence of flat structures on cosmological scales, and has certain advantages over the mentioned models: in accordance with developed approach two-dimensional structures do not arise depending on random disturbances density of the environment, and are causally determined by physical reality in the form of the presence between particles of gravitational interaction, which is what causes do not remain uniform during the evolution of the Hubble flow distributed throughout space, and form a cellular macrostructure (as a consequence of the intersection of strata).

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Bisnovatij-Kogan, G.S. (2011) Relativistic Astrophysics and Physical Cosmology. Krasand.
[2]	Zeldovich, Ya.B. and Novikov, I.D. (1975) Structure and Evolution of the Universe. Nauka.
[3]	Chavanis, P.H. (2003) Gravitational Instability of Isothermal and Polytropic Spheres. Astronomy & Astrophysics, 401, 15-42.[CrossRef]
[4]	Tremaine, S., Hénon, M. and Lynden-Bell, D. (1986) H-Functions and Mixing in Violent Relaxation. Monthly Notices of the Royal Astronomical Society, 219, 285-297.[CrossRef]
[5]	Vlasov, A.A. (1967) Statistical Distribution Functions. Nauka.
[6]	Vlasov, A.A. (1978) Nonlocal Statistical Mechanics. Nauka.
[7]	Dupaigne, L. (2011) Stable Solutions of Elliptic Partial Differential Equations. Springer.
[8]	Gurzadyan, V.G. (1985) The Cosmological Constant in the McCrea-Milne Cosmological Scheme. Observatory, 105, 42.