Self-gravitating Brownian particles in two dimensions: the case of N=2 particles

Self-gravitating Brownian particles in two dimensions: the case of particles

1Introduction

Systems with long-range interactions have recently been the object of considerable interest [1]. One usually considers isolated systems in which the particles evolve according to deterministic Hamiltonian equations. These systems are described by the microcanonical ensemble. Examples of such systems include self-gravitating systems, two-dimensional point vortices, the Hamiltonian mean field (HMF) model etc. However, one may also consider dissipative systems in which the particles, in contact with a thermal bath, evolve according to stochastic Langevin equations. These systems are described by the canonical ensemble. The statistical mechanics of Hamiltonian and Brownian systems with long-range interactions is discussed in [2] at a general level. In this paper, we consider the case of Brownian particles in gravitational interaction. It is known that this system bears deep analogies with simple models of bacterial populations experiencing chemotaxis in biology1 (see, e.g., Ref. [5] for a description of this analogy). In a proper thermodynamic limit , the mean field approximation becomes exact and the dynamics of these systems is described by the Smoluchowski-Poisson system (gravity) [6] or by the Keller-Segel model (chemotaxis) [7]. These equations display rich phenomena such as collapse and evaporation. In particular, in dimensions, the evolution leads to the formation of Dirac peaks if the temperature is below the critical value [8]. These systems have been studied essentially in the mean field approximation, i.e. for a large number of particles. The case of a finite number of particles can be studied numerically by solving the -body stochastic equations. Numerical results will be presented in a companion paper [11]. In the present paper, we consider the extreme case of only Brownian particles in gravitational interaction that can be solved analytically. We analytically obtain the probability density of finding the particles at a distance from each other at time and determine the probability that the particles have coalesced and formed a Dirac peak at time . We also investigate the variance of the distribution and discuss the proper form of the virial theorem for this system. In particular, we show that the virial theorem obtained in [12] is only valid as long as the particles have not formed Dirac peaks.

The paper is organized as follows. In Section 2 we recall the -body coupled stochastic equations describing the evolution of self-gravitating Brownian particles and specifically consider the case . We introduce the center of mass and the reduced particle. We show that the center of mass undergoes a pure Brownian motion and that the reduced particle undergoes a Brownian motion in a central potential . We also recall the “naive” virial theorem obtained in [12] and discuss, with a new light, the distinction between the critical temperatures and . In Section 3, we study the motion of a Brownian particle (reduced particle) in an attractive central potential in . We show that the corresponding Fokker-Planck equation is equivalent to a Schrödinger equation (with imaginary time) with a potential . This equation can be solved analytically in terms of Bessel functions. Then, we can obtain various analytical results such as the probability to find the reduced particle at position at time , the probability that the particle reaches the origin for the first time between and , the probability that the particle has reached the origin at time and the variance of the distribution. We find that the reduced particle has a normal diffusion behavior for small times with a gravity-modified diffusion coefficient and an anomalous diffusion for large times . In particular, the variance increases with time when and tends to zero for when . In Section 4, we consider the case of two self-gravitating Brownian particles in a bounded domain and discuss the differences with the case of an infinite domain. Finally, in Section 5 we show that our study also describes the large time asymptotics of the Smoluchowski-Poisson system (or Keller-Segel model) for . Indeed, in the post-collapse regime, the system is made of a growing central Dirac peak (condensate) surrounded by a dilute halo whose dynamical evolution is eventually described by a Fokker-Planck equation similar to the one studied in the case of particles. We find that the saturation of the mass of the condensate to the total mass is algebraic in an infinite domain and exponential in a bounded domain and we characterize it precisely. In Section 6, we briefly generalize our results to the logarithmic Fokker-Planck equation in dimensions. The Appendices provide complements such as the deterministic limit (Sec. Appendix C), the van Kampen classification (Sec. Appendix F), the correlation functions (Sec. Appendix G) and the general form of the virial theorem for Brownian particles with power law interaction (Sec. Appendix I).

We may note that our study bears some similarities with the stochastic motion (induced by viscosity) of two point vortices studied by Agullo & Verga [13]. However, there also exists crucial differences between the two problems since in our case the interaction is radial leading to the formation of Dirac peaks while in the case of point vortices the interactions is rotational leading to the formation of a spiral structure.

We may also note that the statistical mechanics of particles in gravitational interaction has been considered by Padmanabhan [14] in (and generalized by Chavanis [12] for the dimensions and ) in the microcanonical and canonical ensembles. However, these authors consider the equilibrium statistical mechanics of self-gravitating particles in a box, and with a small-scale cut-off, while we consider here the dynamical evolution of self-gravitating Brownian particles in a finite or infinite domain without small-scale cut-off. Therefore, we address the time dependent problem and investigate the formation of Dirac peaks.

Finally, the particular character of the dimension in gravity is well-known. We refer for example to [15] for more details and further references.

2The position of the problem

2.1The -body problem

We consider a system of overdamped Brownian particles with mass in gravitational interaction in a space of dimension . Their motion is described by the coupled stochastic equations [12]:

with

for and

for . Here, is the friction coefficient and is a white noise satisfying and where refers to the particles and to the coordinates of space. The diffusion coefficient is given by the Einstein relation

where is the temperature. We assume that the friction is the same for all the particles.

From these stochastic equations, it is possible to derive the Fokker-Planck equation for the -body distribution and then write the BBGKY-hierarchy for the reduced distributions [2]. Let us consider for brevity the single-species system. The proper thermodynamic limit corresponds to in such a way that the normalized temperature is of order unity. In that limit, it can be shown that the mean field approximation becomes exact so that the -body distribution factorizes in a product of one-body distributions: [2]. Furthermore, the one-body distribution, or equivalently the smooth density field , is solution of the Smoluchowski-Poisson system [2]:

The equations generalizing Eqs. (Equation 5)-(Equation 6) for the multi-species case are given in [12].

Up to a change of notations, these equations are isomorphic to a simplified version of the Keller-Segel model of chemotaxis that is valid in the limit of large diffusivity of the chemical and in the absence of degradation [5].

2.2The case : the reduced particle

From now on, we consider only self-gravitating Brownian particles in . In that case, the stochastic equations (Equation 1)-(Equation 3) reduce to

with and . Like for the standard two-body problem, we introduce the center of mass

and the reduced particle

Concerning the motion of the center of mass, we have

where the noise satisfies

Therefore, the center of mass undergoes a pure Brownian motion of the form

with a diffusion coefficient

Concerning the motion of the reduced particle, we have

where the noise satisfies

Therefore, the reduced particle undergoes a Brownian motion in a central potential of the form

with a diffusion coefficient

2.3The naive virial theorem

Let us introduce the total moment of inertia

In [12] (see also Appendix I), an exact closed expression of the virial theorem valid for an arbitrary number of self-gravitating Brownian particles in has been obtained. For , it writes

with the critical temperature

It is instructive to recover this result in a different manner. The positions of the particles and can be expressed in terms of (reduced particle) and (center of mass) as

Substituting these relations in Eq. (Equation 19), we obtain after straightforward algebra

a relation which was of course expected. Now, the Fokker-Planck equation associated with the stochastic motion (Equation 17) of the reduced particle is

Taking the time derivative of

and using simple integrations by parts, we naively2 obtain

This relation exhibits a critical temperature

Introducing the moment of inertia of the reduced particle , we can rewrite Eq. (Equation 26) as

The mean square displacement of the reduced particle satisfies

This is a normal diffusion with a gravity modified diffusion coefficient

The variance increases for and tends to zero in a finite time for . On the other hand, the Fokker-Planck equation associated to the stochastic motion (Equation 13) of the center of mass is simply

and we classically obtain the relation

Finally, summing Eqs. (Equation 26) and (Equation 32) and using Eq. (Equation 23), we recover Eq. (Equation 20). We now clearly see the origin of the two temperatures and that were reported in [12]. In the case , the critical temperature is associated to the dynamics of the reduced particle while the critical temperature enters in the expression of virial theorem for the total moment of inertia (reduced particle and center of mass). This distinction is further discussed in Appendix A in the general case of particles.

2.4The problem

In fact, there is a flaw in the above derivation of the virial theorem because we have naively assumed that the normalization is conserved in time. However, as we shall see, this is not correct. The normalization is not conserved in time because the reduced particle can reach the origin and be “lost” by the system (if it reaches the origin, it remains there for ever). This corresponds to the coalescence of the two particles, resulting in the formation of a Dirac peak, i.e. a new particle of mass . As a result of these “trapping” events

and we must reconsider the problem in more detail.

3Brownian particle in a Newtonian potential in two dimensions

3.1The Fokker-Planck equation

Let denote the probability density of finding the reduced particle in at time . The evolution of is governed by the Fokker-Planck equation3

The initial distribution is normalized such that . Introducing

the Fokker-Planck equation can be rewritten

In the absence of small and large scale cut-offs, this equation has no steady state since the distribution is not normalizable. We assume that the initial distribution is radially symmetric, so that is radially symmetric for all times. Therefore, we can write the Fokker-Planck equation in the form

As discussed previously, the probability is not conserved because the reduced particle may reach the origin and form a Dirac peak (the two particles coalesce). The probability that the particle has not reached at time is

Taking the time derivative of this quantity and using the Fokker-Planck equation (Equation 37) we obtain

which is non zero since . Therefore, the probability for the particle to form a Dirac peak between and (i.e. to reach for the first time between and ) is

and the probability for the particle to have formed a Dirac peak at time (i.e. to have reached at time ) is

We obviously have . We can now obtain the proper form of the virial theorem associated to the Fokker-Planck equation (Equation 34). Introducing the moment of inertia of the reduced particle

we easily obtain the virial theorem

instead of Eq. (Equation 28). It has to be noted that this relation is not closed since it depends on that must be obtained by solving the Fokker-Planck equation (Equation 34).

3.2The associated Schrödinger equation

Let us consider a general Fokker-Planck equation of the form

For a spherically symmetric distribution in dimensions, it can be rewritten

As is well-known [23], we can transform this Fokker-Planck equation into a Schrödinger equation (with imaginary time) by setting

This yields

with the potential

For a spherically symmetric distribution, the Schrödinger equation (Equation 47) can be rewritten

with

Making the separation of variables

we obtain the eigenvalue equation

Let us note the eigenvalues and the corresponding eigenfunctions. The eigenfunctions are orthogonal with respect to the scalar product

We also normalize them so that . Then, any function can be expanded in the form

If the spectrum is continuous, the sum over must be replaced by an integral over .

3.3The general solution

We consider the Green function which corresponds to the initial condition

The solution on the Fokker-Planck equation (Equation 45) can be expanded on the eigenfunctions in the form

Noting that

and using the initial condition (Equation 55), we finally obtain

3.4The case of a logarithmic potential in

For the Fokker-Planck equation (Equation 37), we have , and . The potential arising in the corresponding Schrödinger equation (Equation 49) is

Therefore, if we assume that initially

the solution of the Fokker-Planck equation (Equation 37) can be written

where is solution of the differential equation

where

Equation (Equation 62) is a Bessel differential equation that can be solved analytically. The solutions that are finite at the origin are of the form

Substituting Eq. (Equation 64) in Eq. (Equation 61), the solution of the Fokker-Planck equation (Equation 37) can be written

where we have made the change of notation for convenience. The normalization constant is determined so as to recover the initial condition (Equation 60) as . Taking in Eq. (Equation 65), we get

Using the closure relation

we obtain

Comparing with Eq. (Equation 60), we find that . Therefore, the solution of the Fokker-Planck equation (Equation 37) with the initial condition (Equation 60) is

Using the identity

valid for , we find that it can finally be written

The distribution is plotted in Figure 1 at different times and for (corresponding to ).

Figure 1: Probability density of finding the particles at a distance {\bf r} from each other, at different times t and for a temperature T=T_*/2 (a=2). From top to bottom: t=0.01, 0.025, 0.05, 0.1, 0.2, 1.
Figure 1: Probability density of finding the particles at a distance from each other, at different times and for a temperature (). From top to bottom: .

Using the identity

we get for :

which tends to Eq. (Equation 60) as expected. On the other hand, for , the probability tends to zero meaning that the particle has been absorbed in after a sufficiently long time so that it is ultimately lost by the system. For , using the identity (Equation 72), we have

For , using the identity

we get

Finally, for , the probability density (Equation 71) becomes

which is the solution of the diffusion equation in . Indeed, in the limit of infinite temperature, the gravity is negligible with respect to diffusion.

3.5The probability to form a Dirac peak

The probability that the particle has not reached at time is given by Eq. (Equation 38). For the distribution (Equation 71), the integral can be performed analytically and we obtain

where

is the incomplete Gamma function. The probability decays because, as time goes on, the particle has more and more chance to reach and form a Dirac. The probability that the particle reaches for the first time between and is given by Eq. (Equation 40). Combining this relation with Eq. (Equation 76), we obtain

Integrating Eq. (Equation 80), we obtain the probability that the particle has formed a Dirac peak at time t:

and we check that as expected. The evolution of for different values of the temperature is shown in Figs. Figure 2 and Figure 3.

Figure 2: Evolution of the function \chi_D(t), giving the probability that the particle has formed a Dirac peak at time t, for different values of the temperature (we have taken T/T_*=1/2, T/T_*=1 and T/T_*=2). The probability increases more rapidly at smaller temperatures.
Figure 2: Evolution of the function , giving the probability that the particle has formed a Dirac peak at time , for different values of the temperature (we have taken , and ). The probability increases more rapidly at smaller temperatures.
Figure 3: Same as Fig.  for larger times.
Figure 3: Same as Fig. for larger times.

For , using the expansion

we obtain

We see that the probability for due to the exponential factor. This tendency is reinforced by the algebraic factor for () while it is reduced for ().

For , using the expansion

we obtain

Therefore, the probability that the particle has not formed a Dirac at time decreases algebraically as . Equation (Equation 85) can be written in the form

where the time

gives an idea of the rapidity at which the Dirac forms as a function of the temperature . The function is represented in Figure 4 and its asymptotic behaviors are given in Appendix Appendix E. We find that for and for . We note that tends to a finite value for while we know that the system does not form a Dirac peak for (indeed for ). Therefore, the physical interpretation of should be considered with care. Another measure of the effect of the temperature on the formation of the Dirac is provided by the quantity

corresponding to evaluated at . This function is represented in Figure 5 and its asymptotic behaviors are given in Appendix Appendix E. We find that for and for .

Figure 4: Evolution of t_{*}(a) as a function of the temperature a=T_*/T.
Figure 4: Evolution of as a function of the temperature .
Figure 5: Evolution of \chi_D(1) as a function of the temperature a=T_*/T.
Figure 5: Evolution of as a function of the temperature .

Finally, the normalized probability density can be written

where is given by Eq. (Equation 71) and by Eq. (Equation 81). We readily check that .

3.6The moment of inertia

The moment of inertia of the reduced particle is defined by

and the variance of the distribution (mean square displacement) is

For the density distribution given by Eq. (Equation 71), the integral can be calculated explicitly yielding

In Appendix B, we check that this relation is consistent with the virial theorem (Equation 43). The evolution of the moment of inertia is represented in Figure 6 for different values of the temperature. For , the moment of inertia increases, for the moment of inertia is constant and for the moment of inertia decreases. This will become clear from the asymptotic behaviors.

For , we get

which can be rewritten

In that case, we have a normal diffusion with a gravity-modified diffusion coefficient

For the variance increases with time while for it decreases. This expression agrees with the naive virial theorem (Equation 29). Indeed, for small times, the probability for the particle to reach is exponentially small so that the probability is conserved and the naive virial theorem holds since there is no Dirac peak.

For , using the expansion (Equation 84), we get

This corresponds to

For , i.e. , the variance increases and goes to for large times. In that case, we have an anomalous diffusion with an exponent . The evolution is always sub-diffusive. The origin of the anomalous diffusion is related to the fact that the particle can be trapped at (and form a Dirac peak). For , the variance decreases and goes to for .

Figure 6: Time evolution of the moment of inertia for different values of the temperature (we have taken T/T_*=1/2, T/T_*=1 and T/T_*=2).
Figure 6: Time evolution of the moment of inertia for different values of the temperature (we have taken , and ).

3.7The most probable position

The most probable value of the distribution is obtained by maximizing , or equivalently , with respect to . This gives

Using the recurrence relation

we obtain

This equation can be rewritten in the parametric form

or equivalently

which gives . Using the asymptotic expansions of , we find that

and

where

Therefore, the radius is decreasing for any temperature and it goes to zero in a finite time depending on the temperature. Some curves are represented in Figure 7.

Figure 7: Time evolution of r_{P} for different values of the temperature (we have taken T/T_*=1/2, T/T_*=1 and T/T_*=2).
Figure 7: Time evolution of for different values of the temperature (we have taken , and ).
Figure 8: Time evolution of r_{*} for different values of the temperature (we have taken T/T_*=1/2, T/T_*=1, T/T_*=2 and T/T_*=4).
Figure 8: Time evolution of for different values of the temperature (we have taken ,