3 Model

A spatially isolated, spherical N-body system is studied. The N particles of the system are accelerated by gravity, and forces induced by boundary conditions and perturbations, respectively. In all, a particle can be accelerated by five force types: 1) confinement 2) energy dissipation due to velocity dependent friction 3) external forcing 4) long-range attraction (self-gravity) and 5) short-range repulsion. The five force types are discussed in detail in Sects. 3.2-3.5.

3.1 Units

A steep potential well confines the particles to a sphere of radius $R\approx1$ (see Sect. 3.2). The total mass M and the gravitational constant G equal one as well, G=M=1. The dimensionless energy is the energy measured in units of GM²/R, that is in our model $\approx1$. Thus energies and other quantities with energy dimension (like the temperature) can be considered as dimensionless.

   
3.2 Confinement

In order to keep the model simple and to enable a comparison with established theoretical results of canonical isothermal spheres we apply a confinement that prevents gravitationally unbound particles from escaping and keeps thus the particle number constant.

The confinement is realized through a steep potential well

$$\Phi_{\rm conf} \propto R^{16},$$ (2)

where R is the distance from the center of our spatially isolated system. The walls of the potential well are steep, so that the particles are only significantly accelerated by the confinement when they approach R=1. However, the potential walls must not be too steep, to avoid unphysical accelerations near the wall, due to a discrete time-step.

Unlike non-spherical reflecting enclosures, the applied potential well conserves angular momentum.

   
3.3 Energy dissipation

The effect of different dissipation schemes is studied: local dissipation, global dissipation and a scheme similar to dynamical friction. For convenience we call the latter hereafter "dynamical friction''.

   
3.3.1 Local dissipation

The local dissipation scheme represents inelastic scattering of interstellar gas constituents. That is, friction forces are added that depend on the relative velocities and positions of the neighboring particles.

A particle is considered as a neighbor if its distance is $r<\lambda\,\epsilon$, where $\epsilon$ is the softening length and $\lambda$ is a free parameter. The friction force has the form:

$$F_i=\left\{\begin{array}{r@{\;\;:\;\;}l} \Lambda\frac{m^2 r_i... ...ec{V}\cdot\vec{r}<0\\ 0 & {\rm otherwise}, \end{array}\right.$$ (3)

where V is the relative velocity of two neighboring particles, $\zeta$, $\eta\ge 0$ are free parameters and i=x,y,z. The term $\vec{V}\cdot\vec{r}/r=V\cos\vartheta<0$ ensures a convergent flow and that dissipation affects only the linear momentum. As a consequence, angular momentum is locally conserved.

For $\epsilon\ll r < \lambda\,\epsilon$ the friction force is

$$F_i\propto \frac{(V \cos\vartheta)^{\zeta}}{r^{\eta}}\, e_i,$$ (4)

where e_i is the ith component of the unity vector. Thus the friction force increases with the relative velocity and decreases with the distance, provided that $\zeta,\;\eta>0$. The condition $r<\lambda\,\epsilon$ ensures the local nature of the energy dissipation.

3.3.2 Global dissipation

The global dissipation depends not on the relative velocity V, but on the absolute particle velocity v. The global friction force is

$$F_i=-\alpha\, v_i,$$ (5)

where $\alpha$ is a free parameter. The same friction force was already used in the shearing box experiments of Toomre & Kalnajs (1991) and Huber & Pfenniger (2001a, 2001b).

3.3.3 Dynamical friction

Setting the velocity-dispersion $\sigma=1$, which is the dispersion of a virialized self-gravitating system with M=R=1, Chandrasekhar's dynamical friction formula can be parameterized to a good approximation by,

$$F_i=\Gamma\frac{v_i}{(4+v^3)},$$ (6)

where $\Gamma$ is a free parameter (Chandrasekhar 1943; Binney & Tremaine 1994). In contrast to global dissipation this scheme has the feature that F_idoes not increase asymptotically with v. There is a maximum at $v\approx 1.3$. Thus high speed particles in collapsed regions no longer dissipate energy and gravitational runaway is impeded.

   
3.4 Forcing scheme

If the dissipated energy is replenished by a forcing scheme the system can be subjected to a continuous energy-flow. Here the energy injection is due to time-dependent boundary conditions, that is, due to a perturbation potential. This perturbation should on the one hand be non-periodic and quasi-stochastic, on the other hand it should provide an approximately regular large-scale energy injection in time-average. A simple perturbation potential, meeting these conditions, has the linear form,

$$\Phi_{\rm pert}(x,y,z,t)=\gamma [B_x(t) x + B_y(t) y + B_z(t) z],$$ (7)

where $x,\ y$ and z are Cartesian coordinates, $\gamma$ is a free parameter determining the strength of the perturbation, and

 

$$\sum^3_{i=1} A_{x,i}\sin(\omega_{x,i} t + \phi_{x,i})$$ B_x(t) =

$$\sum^3_{i=1} A_{y,i}\sin(\omega_{y,i} t + \phi_{y,i})$$ B_y(t) =

$$\sum^3_{i=1} A_{z,i}\sin(\omega_{z,i} t + \phi_{z,i}),$$ B_z(t) = (8)

where t is the time. The amplitudes A_j,i and the phases $\phi_{j,i}$ are arbitrary fixed constants and remain unchanged for all simulations (j=x,y,z). The frequencies are given through $\omega_{j,i}=\delta\Omega_{j,i}$, where $\Omega_{j,i}$ are arbitrary, but not rationally dependent constants (to avoid resonances) as well. A_j,i, $\phi_{j,i}$ and $\Omega_{j,i}$ lay down the form of the potential. Then, the amplitude and the frequency of the perturbation are controlled by two parameters, namely $\gamma$ and $\delta$, respectively (Huber 2001).

The perturbations are a linear combination of stationary waves that do not inject momentum.

If our system represents a molecular cloud then a perturbation similar to those described above can be due to star clusters, clouds or other massive objects passing irregularly in the vicinity. Indeed such stochastic encounters must be quite frequent in galactic disks and we assume that the average time between two encounters is,

$$\tau_{\rm pert}=\frac{1}{f_{\rm pert}} \mathrel{\mathchoice {... ...$\scriptscriptstyle ...$$ (9)

where $f_{\rm pert}$ is the mean frequency of the encounters and $\tau_{\rm dyn}$ is the dynamical time. Thus, these encounters can provide a continuous low frequency energy injection on large scales.

Such a low frequency forcing scheme induces at first ordered particle motions, which are transformed in the course of time to random thermal motion, due to gravitational particle interaction.

   
3.5 Self-gravity, repulsion and force computation

Interaction forces include self-gravity and repulsive forces whose strength can be adjusted by a parameter. Repulsive forces may be a result of sub-resolution processes such as star-formation or quantum degeneracy. The particle interaction potential reads,

$$\Phi_p (r)=-\frac{Gm}{\sqrt{r^2+\epsilon^2}}\left (1-\xi\frac{\epsilon^2}{(r^2+\epsilon^2)}\right),$$ (10)

where m is the particle mass and r is here the distance from the particle center. On larger scales $\epsilon)$ the potential is a Plummer potential $\Phi_{\rm Pl}$ (Plummer 1911; Binney & Tremaine 1994). On small scales $\epsilon)$ the deviation from a Plummer potential and the strength of the repulsive force is determined by the parameter $\xi$. The interaction potential becomes for instance repulsive in the range $\vert r\vert<\epsilon \sqrt (3\xi -1)$ if $\xi > 1/3$ (see Fig. 1). If $\xi =0$, $\Phi_{\rm p}=\Phi_{\rm Pl}$ holds.

Figure 1: The potential as a function of r in units of Gm/R, where m is the particle mass. A test particle at $\vert r\vert<\epsilon \sqrt (3\xi -1)$ feels a repulsive force if $\xi > 1/3$. The size of the softening parameter $\epsilon$ is indicated, as well.

The interaction forces are computed on the Gravitor Beowulf Cluster at the Geneva Observatory with a parallel tree code. This code is based on the Fortran Barnes & Hut (1986, 1989) tree algorithm, and has been efficiently parallelized for a Beowulf cluster. It is available on request.

The time-integration is the leap-frog algorithm with uniform time-step, which ensures the conservation of the symplectic structure of the conservative dynamics. The time-step is,

$$\Delta t \le 0.1\frac{\epsilon}{\sigma_{\rm v}},$$ (11)

where $\sigma_{\rm v}$ is the velocity-dispersion of the initial state.

The accuracy of the force computation is given through the tolerance parameter, which for the studies presented in this paper is $\theta\le 0.58$.

3.6 Code testing

In order to test the code, the evolution of the angular momentum is checked. We find that the angular momentum Jof a system with $10\,000$ particles and local dissipation scheme, that is initially $J=3.4\times 10^{-5}$ in units such that energy is dimensionless, remains small, $J<4\times 10^{-5}$. This is illustrated by means of two examples. Two simulations with equal dissipation strength are carried out. One with short range repulsion, $\xi =2/3$, and the other without, $\xi =0$. After 10crossing times the angular momentum is $J=3.3\times 10^{-5}$ and $J=3.8\times 10^{-5}$, respectively. The deviation from the initial value is larger for the simulation without short-range repulsion. This is because the dissipation leads in this case to a stronger mass concentration, i.e., to shorter particle distances. Indeed, the potential energy after 10 crossing times is $U\approx-12$and $U\approx-17$ for the simulation with and without short-range repulsion, respectively.

3.7 Parameters

The different dissipation schemes, the forcing scheme and the interaction potential with the short-distance repulsive force, have several free parameters. In order to do a reasonably sized parameter study the particle number has to be limited to a maximum of $N=160\,000$. However, the absolute particle number is not important, but the effect of changing it is. Thus, N also is varied.

The model parameters and their ranges are indicated in Table 1.

  

Table 1: Model paramaters, characterizing energy injection (perturbation potential), dissipation, interaction potential (self-gravity and repulsion) and mass resolution.

$$\begin{tabular}{c\vert r@{$,\ldots,\,$}l\vert l\vert l} \hlin... ...0$ & Particle number & Mass resolution\\ \hline \end{tabular}$$