text:
# Slow evolution of elliptical galaxies induced by dynamical friction I. Capture of a system of satellites

## 1 Introduction

If we consider a globular cluster of mass $`m_{gc}10^5M_{}`$ orbiting inside a spherical galaxy, given the mass ratio $`10^5`$ with respect to the mass of a typical star, we find that it should suffer dynamical friction on a time-scale $`T_{fr}10^9`$yr, through scattering of the stars. \[The term dynamical friction refers to the slow relaxation process associated with discrete two-body encounters investigated in the pioneering studies of Chandrasekhar (cha43 (1943)). In these classical analyses, estimate is given of the cumulative effect of discrete encounters on a test particle passing through a homogeneous system of field particles. This effect, associated with the granularity of the system in phase space, is ignored in pure Vlasov descriptions of stellar systems. The term dynamical friction is often extended to describe loosely the physical origin of orbital decays of satellites of finite mass also in real inhomogeneous systems and simulations, which fall obviously outside the classical description of Chandrasekhar.\] If the cluster starts from a circular orbit located at $`r_i`$, in a time $`T_{fr}`$ it will reach a lower orbit close to $`r_f`$, while a fraction of the energy initially associated with such a cluster $`E_i=m_{gc}\mathrm{\Phi }(r_i)+(1/2)m_{gc}r_i(d\mathrm{\Phi }/dr)_{r_i}`$ is released into the distribution function of the hosting galaxy. If we refer to the realistic possibility of the presence of several thousands of such globular clusters we may thus argue that on a fraction of the Hubble time the underlying stellar system may well absorb a non-negligible fraction of its binding energy and thus should slowly evolve (see also Lotz et al. lot01 (2001)).

Under the suggestive scenario just outlined, the issue we would like to address in this paper is the following conceptually challenging problem. In a two-component system, where the more massive component (the “galaxy”) is basically collisionless and the minority component (the globular cluster system or a large number of mini-satellites caught in minor mergers) is slowly dragged in toward the galaxy center by processes of dynamical friction with the stars of the galaxy, how does the distribution function of the galaxy evolve during the process? A resolution of this challenging problem requires a clarification of the basic mechanism of dynamical friction under inhomogeneous conditions and a discussion of the collective effects involved in the reaction of the underlying distribution function. Note that while the infall of a single cluster (or satellite) is a problem with a preferred direction (given by the angular momentum of the satellite), the capture of a large number of randomly distributed clusters (or mini-satellites) can preserve spherical symmetry, if present to begin with.

It might be argued that the slow evolution in the quasi-isotropic many-cluster case described above could be modeled in terms of an adiabatic evolution of a distribution function within spherical symmetry. For the different problem of the adiabatic growth of a central massive black hole, such a procedure has indeed been developed (e.g., see Cipollina & Bertin cip94 (1994) and references therein). Here it is not clear how a semi-analytic procedure to calculate such slow evolution, taking advantage of the conservation of the relevant actions, can be formulated. In addition, given the relation between dynamical friction and resonant interactions (e.g., see Weinberg wei89 (1989)), in spite of its slowness the process is likely to be inherently non-adiabatic. Therefore, it seems that the only viable way to an answer is that of approaching the problem by means of suitable N-body simulations.

In the long run, results from this project could also shed light on another theoretically interesting scenario, which so far has found major applications only in the context of the internal dynamics of globular clusters. This is a picture, pioneered by Bonnor (bon56 (1956)) in his analysis of gas spheres, by which stellar systems may be subject to a process of gravitational catastrophe (Lynden-Bell & Wood lyn68 (1968)), when their central concentration exceeds a certain threshold value. The catastrophe consists in an instability process, where evolution leads to cores that become more and more concentrated. The source of this instability is thermodynamical, in the sense that it results from a counter-intuitive property of self-gravitating systems, which can be described in terms of negative specific heat. Recently, this general scenario has been demonstrated to hold for an astrophysically interesting family of partially relaxed stellar systems (Bertin & Trenti 2003).

The study of relaxation in N-body simulations and by means of N-body simulations is a vast field of research (e.g., see Hohl hoh73 (1973) and Huang et al. hua93 (1993)). The study of the sinking of a single satellite inside a galaxy, as a result of dynamical friction, was undertaken systematically about two decades ago, initially through simplified simulations, within a restricted three-body problem and with the center of mass of the primary galaxy being kept fixed, and then by means of more realistic, self-consistent simulations (Bontekoe & van Albada bon87 (1987); Zaritsky & White zar88 (1988); Bontekoe bon88 (1988)). The results showed that the general scaling predicted by the classical formulae derived in the “homogeneous limit” (or “local limit”) (Chandrasekhar cha43 (1943)) was basically applicable, although the sinking time scale turned out to be significantly longer (by a factor of $`23`$ for the regimes studied by the numerical experiments) if the self-consistent response of the galaxy was properly incorporated. A convincing interpretation of the sinking process can be made, at least in the limit where the mass of the sinking satellite is small, by investigating the linear response in the galaxy distribution function through stellar dynamical techniques (Fridman & Polyachenko fri84 (1984); Palmer & Papaloizou pal85 (1985); Weinberg wei86 (1986), wei89 (1989); Bertin et al. ber94 (1994), hereafter BPRV94; Palmer pal94 (1994)). Such response consists of an “in-phase” contribution, which has little relevance to the mechanism responsible for dynamical friction, and an “out-of-phase” contribution, associated with the “resonant” interaction, which produces a wake able to set a significant torque on the satellite. These techniques demonstrate that the collective wake in the galaxy, largely dominated by the dipolar ($`l=1`$) response, can trail significantly away from the orbiting object (when the satellite is outside the galaxy, the response is concentrated around a location opposite to the object with respect to the global center of mass), so that the resulting torque on the satellite that is responsible for dynamical friction can be significantly reduced with respect to other more naive calculations. Discrepancies between earlier simulations and the predictions of the stellar dynamic semi-analytic theory were ascribed to non-linearity effects (in the sense that the satellite-to-galaxy mass ratio $`M_s/M_G0.1`$ considered in the simulations was judged to be excessive; see Hernquist & Weinberg her89 (1989)), but the general picture appears to have been clarified. It is noted that the differences in behavior with respect to the naive expectations of the traditional formulae for dynamical friction become substantial when the size of the satellite is finite (Weinberg (wei89 (1989)) refers to the length ratio $`R_s/R_G`$; the multipoles associated with the perturbation with $`l`$ above a certain threshold are bound to be unimportant), while the differences should become negligible for very compact satellites (in which case all multipoles contribute, with a divergence that is removed by the Coulomb logarithm).

Several aspects of the general problem of dynamical friction have been revisited recently. The issue of the contrast between “local” and “non-local” effects has been re-discussed (Prugniel & Combes pru92 (1992); Maoz mao93 (1993); Cora, Muzzio & Vergne cor97 (1997)), with the confirmation that the description of the friction process in terms of the “homogeneous limit” (Chandrasekhar cha43 (1943)) should be adequate in the limit of low-mass, relatively compact satellites. A very recent investigation (Cora, Vergne & Muzzio cor01 (2001)) also suggests that the process depends very little on the regularity or the stochasticity of the stellar orbits in the galaxy, in contrast with some earlier conjectures (Pfenniger pfe86 (1986)). Other interesting problems are the problem of circularization (of the satellite orbit), the relation of the friction effect to the past history of the dynamical event under consideration, and the issue of a direct calculation of the relevant dissipative force from the response of the stellar system to the perturbing satellite (e.g., see Séguin & Dupraz seg94 (1994, 1996); see also Colpi col98 (1998); Colpi & Pallavicini col98a (1998); Colpi, Mayer & Governato col99 (1999); Nelson & Tremaine nel99 (1999)).

Note that so far investigations of processes of dynamical friction have mostly focused on the issue of the fate of the sinking objects (and mostly for the case of one object; but see Tremaine, Ostriker & Spitzer tre75 (1975)), while the problem of the reaction of the distribution function of the underlying stellar system (especially in the quasi-spherical limit of a system of many heavy objects outlined above) remains practically unexplored.

Our paper is organized as follows. In Sect. 2 we summarize the theoretical framework that allows us to understand the mechanisms of dynamical friction within the context of a single sinking satellite. In Sect. 3 we describe the model adopted. For the present exploratory investigation, the choice of the model (both for the galaxy and for the satellite or the fragments) is mostly dictated by convenience (for the study of dynamical mechanisms and for an easier comparison with previous analyses) rather than by realistic requirements. In Sect. 4 we describe the diagnostic tools that will be used in the simulations and some preliminary tests. In Sect. 5 we proceed to illustrate the results of the simulations for the case of one sinking satellite and for the quasi-spherical case for which the sinking of many smaller mini-satellites is considered. In Sect. 6 we briefly draw our conclusions and set out the plans for future investigations.

## 2 Theoretical framework for the problem of a single sinking satellite

Here we briefly recall the main issues involved in the problem of a single sinking satellite of “diameter” $`2R_s`$, initially located on a circular orbit at distance $`r_sR_G`$ from the center of mass of the entire system. Here $`R_G`$ denotes the “radius” of the galaxy. Let $`\mathrm{\Omega }_s>0`$ be the angular velocity of the satellite along its circular orbit. The case of a very compact satellite is likely to be within the reach of the classical formulae for dynamical friction (Chandrasekhar cha43 (1943)). In the case where the satellite is sufficiently extended, for the description of its density distribution and the associated gravitational potential, only the first harmonics in a multipole expansion are important, up to $`l_{max}\pi r_s/(2R_s)`$ (see Weinberg wei89 (1989)).

A “friction” torque acting on the satellite is expected to result from the action of suitable resonances with star orbits. To understand this, we may refer to the action of the satellite as a perturbation over a basic state that is non-rotating, spherically symmetric, and characterized by isotropic distribution function $`f_0(E)`$, where $`E=v^2/2+\mathrm{\Phi }_0(r)`$ is the star energy, per unit mass, in the absence of the perturbation, and we may then follow the analysis developed earlier (BPRV94).

Analytical results can be derived in the limit when the mass of the satellite is vanishingly small (with respect to the mass of the galaxy, so that the resulting perturbation in the galaxy distribution function can be treated by a linear theory; see subsection 2.2). In practice, for real systems or for numerical simulations, some of the evolution effects are going to be associated with the finite mass of the satellite. Once we move to consider the possibility of a quasi-equilibrium configuration for a composite system (galaxy plus satellite or a population of mini-satellites), we then have to worry about genuine instabilities that may occur, independent of the granularity present in the system. All of this will contribute to evolution although, in principle, it is not related to the original gedanken experiment at the basis of the classical definition of dynamical friction.

Another subtle aspect of the problem is related to geometry. Broadly speaking, we may think that the difference between the case of straight orbits (considered in the classical description of dynamical friction) and that of orbits in the inhomogeneous potential of the galaxy is taken into account in the following discussion of the resonant response. Yet, the cumulative effects on a single star, because of multiple encounters with the satellite (for example, for stars on quasi-circular orbits), grow in time differently from the classical picture, even before dealing with the collective behavior of the entire set of stars that characterizes the resonant response.

### 2.1 Single-star Resonances

The frequency $`\omega `$ associated with the $`m`$-th component of a multipolar term of order $`l`$ in the expansion of the potential $`\mathrm{\Phi }_s`$ of the satellite is readily related to the angular frequency of the satellite $`\omega =m\mathrm{\Omega }_s`$. Therefore, if we refer to Eq. (80) of BPRV94, we find that the appropriate resonance condition is $`m\mathrm{\Omega }_sp\mathrm{\Omega }_r(E,J)+s\mathrm{\Omega }_\theta (E,J)=0`$; here $`s`$, $`m`$, $`p`$ are integers ($`ls,m+l`$; $`\mathrm{}<p<\mathrm{}`$) and $`\mathrm{\Omega }_r(E,J)>0`$, $`\mathrm{\Omega }_\theta (E,J)>0`$ are the two frequencies (dependent on the star specific energy $`E`$ and angular momentum $`J`$) that characterize the star orbit and depend on the properties of the basic state. As discussed in BPRV94, only the components with $`s=l,l2,..l`$ contribute to the density perturbation.

For an extended satellite, the $`p=0`$ component, which corresponds to an average along the star radial orbit excursion, is expected to dominate so that the resonance condition should reduce to $`\mathrm{\Omega }_\theta (E,J)=m\mathrm{\Omega }_s/s`$. Note that in the course of time, because of dynamical friction, the angular frequency of the satellite is going to change. Clearly, the resonance conditions that we have just stated are applicable provided that $`\dot{\mathrm{\Omega }}_s/\mathrm{\Omega }_s\mathrm{\Omega }_\theta (E,J)`$, which we can rewrite as $`\dot{\mathrm{\Omega }}_s\mathrm{\Omega }_s^2/l`$.

### 2.2 “Shielding” and global resonances

If we consider the linearized Vlasov-Poisson equation for our problem, we see that Eq. (54) of BPRV94 is modified into $`f_1=(\mathrm{\Phi }_s+\mathrm{\Phi }_1)`$, where $``$ is a linear operator that basically describes one integration along the unperturbed characteristics. The quantity $`f_1`$ can be called the linear response of the distribution function. This response has a “driven” part, associated with $`\mathrm{\Phi }_s`$, and a “self-gravitating” part, related to the potential $`\mathrm{\Phi }_1`$, that must be consistent with the response itself, following the Poisson equation $`^2\mathrm{\Phi }_1=4\pi Gf_1d^3v`$.

In plasma physics, the latter contribution is often called the “shielding” associated with the collective behavior of the plasma. Such shielding is often ignored, although it should be emphasized that this step is improper, because the omitted term is a linear contribution. A relatively simple analysis can then be made, by neglecting the self-gravity of the response, showing the connection between single-star resonances and dynamical friction. Note that this is the only ground where a comparison with Chandrasekhar’s (cha43 (1943)) work can be made.

The general case, which includes the possibility of resonance with discrete global modes of the system, is currently beyond the reach of analytical investigations (see also Vesperini & Weinberg 2000).

## 3 Model

As a model of the primary stellar system (the “galaxy”), following Bontekoe & van Albada (bon87 (1987)) we consider an isotropic polytropic model of index $`n=3`$. The distribution function is thus of the form $`f_G(r,v)=A(E_0E)^{n3/2}`$ for $`E<E_0`$, with $`f_G(r,v)=0`$ for $`E>E_0`$, where $`E=\mathrm{\Phi }_G(r)+v^2/2`$ is the single-star energy, $`E_0=\mathrm{\Phi }_G(R)=GM/R`$ is the gravitational potential at the boundary $`r=R,`$ and $`M`$ is the total mass of the system. Polytropes can be described in terms of the Lane-Emden function $`\theta (\xi ),`$ which satisfies the differential equation:

$$\frac{1}{\xi ^2}\frac{d}{d\xi }(\xi ^2\frac{d\theta }{d\xi })=\theta ^n,$$
(1)
where $`n`$ is the polytropic index, with the boundary conditions $`\theta =1`$, $`d\theta /d\xi =0`$ at $`\xi =0`$. From the function $`\theta (\xi )`$ one can find the potential $`\mathrm{\Phi }_G(r)`$ and the associated density $`\rho _G(r)`$.

The secondary object (the “satellite”) is described by a rigid Plummer density distribution (corresponding to a polytrope of index $`n=5`$), characterized by $`\theta =(1+0.5\xi ^2)^{1/2}`$. The potential and the mass distribution of the satellite are then defined by two scales, $`R_s`$ and mass $`M_s`$, so that $`\mathrm{\Phi }_s(r)=GM_s(R_s^2+r^2)^{1/2}`$, with cumulative mass $`M_s(r)=M_sr^3(R_s^2+r^2)^{3/2}`$. This distribution extends to infinity, but 90% of the mass is contained within $`3.7R_s.`$

As a result of the interaction with the satellite, the distribution function of the galaxy is a time-evolving $`f`$ that will be associated with a mean field $`\mathrm{\Phi }`$, to be derived from the Poisson equation $`\mathrm{\Delta }\mathrm{\Phi }(r,t)=4\pi Gfd^3v=4\pi G\rho (r,t)`$. We solve the Poisson equation on a spherical polar mesh of $`24\times 12\times 24`$ cells in the $`r,\theta `$, and $`\phi `$ directions. The radial mesh is logarithmic. The angular mesh is fixed, with 12 equal divisions in $`\mathrm{cos}\theta `$ and 24 in $`\phi .`$ On this mesh the density, the potential, and the force field are expanded in associated Legendre functions $`P_l^m(\mathrm{cos}\theta )`$ up to $`P_6^6.`$ At each time-step the center of the spherical mesh used to calculate the potential and the force field in the galaxy is repositioned so as to coincide with the center of mass of the galaxy. The equations of the motion (for the particles representative of the galaxy and, separately, for the satellite) are integrated in Cartesian coordinates with the so called leapfrog scheme. This is basically the N-body code described by van Albada & van Gorkom (van77 (1977)) and then improved and used by van Albada (van82 (1982)) and by Bontekoe & van Albada (bon87 (1987)). Therefore, we do not provide here a full description of the method employed, since it is available in the papers quoted.

### 3.1 Sinking of a single satellite

A first set of numerical experiments addresses the sinking of a single satellite, initially located on a circular orbit at $`\theta =\pi /2`$.

The galaxy is considered as a collisionless stellar system, so that its representative particles $`(i=1,\mathrm{}..,N)`$ interact with each other through the mean gravitational potential $`\mathrm{\Phi }`$ and directly with the satellite:

$$\frac{d^2r_i}{dt^2}=\mathrm{\Phi }(r_i,t)\frac{GM_s(r_ir_s)}{[R_s^2+(r_ir_s)^2]^{3/2}}.$$
(2)
Note that $`\mathrm{\Phi }`$ is updated separately by the Poisson-solver scheme, following the original method described by van Albada & van Gorkom (van77 (1977)) (this step effectively closes the system of equations), which guarantees a pure time-reversible evolution.

In turn, the satellite is taken to interact with all the representative particles of the galaxy directly via two-body forces:

$$\frac{d^2r_s}{dt^2}=\underset{i=1}{\overset{i=N}{}}\frac{Gm(r_ir_s)}{\left[R_s^2+(r_ir_s)^2\right]^{3/2}},$$
(3)
In the equations above, $`r_i,m`$ are the position vectors and the mass of the representative particles of the galaxy, while $`r_s`$ is the position vector of the satellite. All position vectors are measured with respect to an inertial frame of reference with its origin at the center of mass of the combined (galaxy + satellite) system.

### 3.2 Sinking of many fragments

Another set of experiments studies an initial configuration where the satellite is fragmented into several ($`N_f`$) pieces (“mini-satellites” or “fragments”). For simplicity, we take the fragments to be characterized by the same Plummer density distribution with mass $`M_f`$ (such that $`M_s=N_fM_f`$) and the same lengthscale $`R_f=R_s`$ (a smaller lengthscale $`R_f`$ would be more realistic for the astronomical system we have in mind, but would be less easy to simulate).

The mutual interactions between mini-satellites can be kept “on” by considering

$$\frac{d^2r_i}{dt^2}=\mathrm{\Phi }(r_i,t)\underset{j=1}{\overset{j=N_f}{}}\frac{GM_f(r_ir_{fj})}{[R_f^2+(r_ir_{fj})^2]^{3/2}},$$
(4)
for the representative particles of the galaxy (labeled by $`i=1,\mathrm{}.,N`$), and

$$\frac{d^2r_{fj}}{dt^2}=a_j^{(1)}+a_j^{(2)}$$
(5)
for the mini-satellites (labeled by $`j=1,\mathrm{},N_f`$), where

$$a_j^{(1)}=\underset{i=1}{\overset{i=N}{}}\frac{Gm(r_ir_{fj})}{\left[R_f^2+(r_ir_{fj})^2\right]^{3/2}}$$
(6)
is the acceleration produced by the galaxy and $`a_j^{(2)}`$ denotes the gravitational acceleration exerted on the $`j`$-th mini-satellite by the other fragments. The gravitational interaction between two extended bodies is complicated to describe. For the purposes of the present paper, where we are not interested in the phenomena associated with the details of such interaction, we simplify the simulations by referring to the following prescription:

$$a_j^{(2)}=\underset{k=1,kj}{\overset{k=N_f}{}}\frac{GM_f(r_{fk}r_{fj})}{\left[R_f^2+(r_{fk}r_{fj})^2\right]^{3/2}}.$$
(7)
All position vectors are measured with respect to an inertial frame with its origin at the center of mass of the entire system.

The mutual interaction between mini-satellites can be ignored and turned “off” simply by taking $`a_j^{(2)}=0`$. This procedure may be the more realistic procedure if we wish to simulate a system of globular clusters inside a galaxy. In fact, for such system, the few-body-problem effects associated with such interaction are expected to play a secondary role.

As to the initial distribution of the fragments we have considered three cases, to be described below. The first is the most natural generalization of the study of the sinking of a single satellite. A posteriori, we have found that in such a case the effects related to the lack of equilibrium in the initial conditions are too strong and confuse the picture of evolution. (Similar problems, related to the lack of equilibrium, were present also in the case of a single satellite, but those are not emphasized in this paper because that case has been well studied in the past and we prefer to focus here on the new quasi-spherical problem in the presence of many fragments). Therefore we have devised quieter starts so as to better disentangle the effects associated with the granularity of the system from those associated with the lack of equilibrium.

Given the relatively small number of fragments involved, we anticipate that different realizations of the physical cases described below may be associated with different evolutions.

#### 3.2.1 I: Fragments quasi-isotropically distributed on circular orbits on a thin shell

We first consider the case where all the fragments are distributed initially with random positions on a thin shell, located at radius $`r_{shell}(0)`$. The fragments are given initial velocities appropriate for circular orbits, with random orientations, based on the gravitational field generated by the galaxy alone.

The purpose of simulations of this type is to derive information relevant for a more realistic situation where the fragments have a full three-dimensional space distribution, so that the interaction between the galaxy and the minority population is differential with radius. Such ideally desired simulations are not feasible with current codes. To move one step further in the desired direction we have then considered the case of two shells of fragments initially located at different radii.

Unfortunately, due to the finite mass of the set of fragments, these simulations are characterized by violent epicyclic oscillations that shake the system on the fast dynamical time scale and tend to persist long after the initially expected transient.

#### 3.2.2 II: Fragments quasi-isotropically distributed on circular orbits on a thick shell with improved initial conditions

In order to set up a smoother quasi-equilibrium initial condition, we have decided to proceed as follows. We refer to a spherically symmetric density distribution for a shell of the form $`\rho _{shell}=\rho _0\mathrm{exp}[(rr_{shell}(0))^2/a^2]`$ for $`|rr_{shell}(0)|2a`$ (and vanishing for $`|rr_{shell}(0)|>2a`$). The normalization constant $`\rho _0`$ is taken in such a way that the total mass associated with the shell is $`M_s`$. In view of the choice of model for the fragments, we have decided to take $`2a=R_s`$.

Such density distribution generates a potential $`\mathrm{\Phi }_{shell}(r)`$, which can be calculated numerically. We now consider a modified distribution function for the galaxy $`f_G(r,v)`$ taken to be of the polytropic form as in Sect. 3, but with $`E=\mathrm{\Phi }_G(r)+\mathrm{\Phi }_{shell}(r)+v^2/2`$. We then have to integrate a Lane-Emden equation similar to Eq. (1) but modified by the presence of the external potential $`\mathrm{\Phi }_{shell}(r)`$. This yields a potential $`\mathrm{\Phi }_G(r)`$ different from that considered in Sect. 3, to be inserted into the definition of $`f_G(r,v)`$.

Finally, we consider the clumpy realization of the shell density distribution by reducing it to a set of $`N_f`$ fragments of mass $`M_f`$. The $`j`$th fragment is initially distributed in $`r`$ at a location $`r_{fj}`$ (with random angular coordinates) so that the set of fragments reproduces the assumed density distribution of the shell; its initial velocity is that appropriate for the circular velocity of a test particle in the combined potential $`\mathrm{\Phi }_G(r)+\mathrm{\Phi }_{shell}(r)`$.

This procedure ensures a quieter start while leaving the basic picture of the single shell case unchanged. It can be generalized to the case of two shells of fragments initially located at different radii.

#### 3.2.3 III: Fragments constructed by clumping part of the galaxy distribution function

The cases considered above are aimed at describing the interaction between the galaxy and a minority component with phase-space properties distinct from those of the galaxy. As we have seen while addressing the need for a quieter start, such distinction is likely to be accompanied by effects that are due to the lack of equilibrium in the initial conditions, difficult to separate from the secular effects that are traditionally called dynamical friction. In this assessment, we can go even one step further and we may argue that some of these initial conditions (with the different phase-space properties considered) may actually be characterized by dynamical instabilities. If this were the case, it would be very hard or even impossible to distill the true effects of dynamical friction out of the simulations. In order to address this point, we will run parallel simulations of case II with a smooth, Vlasov realization of the shell density, which should allow us to single out the effects of possible Vlasov instabilities of such system (see Subsect. 3.2.4).

In view of this discussion, we have decided to address a third set of simulations that is less relevant for the astrophysical case that has motivated this paper but is definitely more interesting for the study of the problem of dynamical friction.

Here we consider the initial equilibrium state to be basically that described by the galaxy model alone, without satellites or shells or fragments, as defined in the first paragraph of Sect. 3. The mini-satellites here are then introduced by clumping together part of the distribution function $`f_G`$ in such a way that it corresponds to a diffuse “shell” in space. Formally, we separate out a piece of the distribution function $`f_{shell}`$ (for simplicity, we keep the notation “shell” even if we are considering a situation different from that of case I and case II). Then the system made of $`f_Gf_{shell}`$ is simulated by the Vlasov code, while the contribution $`f_{shell}`$ is reduced to clumps that are brought back to the form of the fragments and treated by means of direct interaction with the galaxy representative particles (as described in subsection 3.2). Note that in this case the “clumps” or “fragments” are distributed in phase space exactly as the particles of the galaxy (thus minimizing the lack of equilibrium in the initial conditions and the risk of introducing undesired sources of instability); in particular, the clumps are initially characterized by an isotropic distribution in velocity space.

We should emphasize that this third class of simulations, while addressing the settling of heavy masses dragged inwards by dynamical friction, is conceptually distinct from the study of the mass stratification process in self-gravitating systems, which may be performed by means of direct or Fokker-Planck codes (e.g., see Spitzer spi87 (1987)).

#### 3.2.4 Purely collisionless reference models

For the three cases just outlined, the simulations carried out following the description given in subsection 3.2 will be accompanied by reference simulations of similar but purely collisionless models. These can be seen as cases for which the number of fragments formally diverges, while the total mass contained in the fragments remains constant. In the simulations the case of very large number of fragments cannot be treated in terms of direct interaction; instead, for these cases, the fragments are treated as simulation particles of the collisionless component. It is clear that for case III, in such reference collisionless models, the identification of $`f_{shell}`$ is meaningless. In practice, to monitor the properties of our integration scheme, we have labeled the representative particles associated with $`f_{shell}`$ by a different “color”.

The purpose of these parallel reference runs is to check, as much as possible, that we are not confusing evolution effects associated with the initial conditions with the effects that are genuinely associated with the granularity of the system.

#### 3.2.5 Initial set up

In order to derive the initial distribution of the particles in the simulation we start from the mass distribution $`M(r)`$ which is known on a discrete set of points labelled by the index $`j`$, with $`M_j=M(r^j)`$.

First we determine the distance of a particle from the galaxy center. If all the particles of the galaxy have the same mass $`m=M/N`$, then for the $`i^{th}`$ particle we define $`xx=m(i1/2),`$ and then find the cell $`j`$ such that $`M_j>xx>M_{j1}.`$ The distance from the center $`r_i`$ of the $`i^{th}`$ particle is then found by an interpolation between $`r^{j1},r^j`$, and $`r^{j+1}`$. Next we choose two uniformly distributed random numbers: $`q_\phi [0,1]`$ for the $`\phi `$ coordinate of a particle, such that $`\phi =2\pi q_\phi `$, and $`q_\theta [1,1]`$ for the $`\theta `$ coordinate, such that $`\mathrm{cos}\theta =q_\theta .`$ The Cartesian coordinates are then found from the spherical coordinates by standard formulae. The choice of the initial velocities is performed in a similar way, using a rejection method (e.g., see Press et. al. (pre92 (1992))) to obtain the desired distribution function. For the potential $`\mathrm{\Phi }_G(r)`$ we proceed as for the mass distribution; the value of the potential at the radius of the particle location is found by interpolation between the values at $`r^{j1},r^j`$, and $`r^{j+1}.`$

In Case III, of a “subtracted” shell, the particles of the galaxy are divided in two “species”: light particles, with mass $`m=M/N=(MM_s)/N_1`$, and heavy particles with mass $`M_f=M_s/N_f`$ respectively, where $`M_s`$ is the mass of the shell, $`N_f`$ the number of heavy particles, and $`N_1`$ the number of light particles. Initially, the heavy particles are all located inside a thin shell. Inside this shell there are no light particles. The $`i^{th}`$ light particle, for $`i[1,i1=NM(r_{shell})/M],`$ where $`r_{shell}`$ is the position of the “shell”, is located as in the non-subtracted case at the distance $`r_i`$ obtained by interpolation between the $`r^{j1},r^j`$, and $`r^{j+1},`$ where $`j`$ is such that the condition $`M_j>m(i1/2)>M_{j1}`$ is satisfied. Then, for the $`k^{th}`$ heavy particle, with $`k[1,N_f]`$, we adopt $`xx=m(i11/2)+M_f(k1/2)`$ and then find the cell $`j`$ where $`M_j>xx>M_{j1}.`$ The distance from the center $`r_k`$ of the $`k^{th}`$ heavy particle is then found by the interpolation rule mentioned above. After that the remaining light particles are distributed according the previous rule. Finally, in the simulations, the light particles interact with one another through their mean field, while the heavy particles interact directly with the light particles and among themselves, as described at the beginning of Sect. 3.2.

#### 3.2.6 Other possible initial conditions

Another possible initial condition is to turn on slowly the gravity associated with the massive objects or massive shells. In a collisionless simulation, it would lead to the “empirical construction” of equilibrium systems in which dynamical friction is absent and would thus be, for our purposes, essentially equivalent to the procedure described in Case II. Such construction of equilibrium models can also be approached analytically (see Cipollina & Bertin cip94 (1994) and references therein).

The case where the satellite or the fragments are treated as discrete objects would provide information on the friction force, different from the information provided in this paper, but such approach would have no special merits, in terms of physical justification and physical interpretation, with respect to the one adopted here.

### 3.3 The choice of the code

The choice of code that we have made (basically, the code used by Bontekoe & van Albada (bon87 (1987))) appears to be appropriate for the study of the effects of dynamical friction, which are dominated by low-$`l`$ wakes (see Sect. 2). Still, one may ask the reason why we have resorted to that code instead of codes that are commonly utilized for a variety of studies in galactic dynamics or in the cosmological context. In particular, one popular class of codes, the tree-codes (e.g., see Barnes & Hut (bar86 (1986))), has found modern realizations (e.g., GADGET; see Springel, Yoshida, & White spr01 (2001)) that are particularly appealing since they are flexible, widely used, and have gone through many tests. Tree-codes treat near particles in terms of direct interactions via a softened potential, because the use of the true Newtonian potential would introduce a generally undesired collisionality among the simulation particles. Tree codes are often used to study the evolution of collisionless systems; for this purpose, the softening length has to be chosen in an optimal way, depending on the physical characteristics of the system that is being investigated.

The use of a tree code in our project would then pose severe interpretation problems, because the non-physical effects associated with the introduction of the softening parameter would show up precisely in relation to the issues that we are trying to investigate.

## 4 Diagnostics and test cases

The results of our simulations can be studied in terms of the energy and angular momentum budget as a function of time and by means of figures that illustrate the spatial structure of the density and of the pressure distribution within the galaxy in the course of time. Another set of figures will provide plots of the azimuthally averaged density profile for the galaxy and of the location $`r_s`$ of the satellite as a function of time; for the case of many mini-satellites, we will show a plot of the radius $`r_x`$ of the sphere that contains different fractions ($`x`$ percent) of the total mass in the form of mini-satellites. We first proceed by specifying units and definitions for the relevant quantities.

### 4.1 Units

Unless specified otherwise, the units adopted for mass, length, and time are the mass of the galaxy ($`M`$), the radius of the galaxy ($`R`$), and the natural crossing time ($`t_{cr}=GM^{5/2}(2K_{gal})^{3/2}`$ evaluated at the beginning of the simulation; here $`K_{gal}`$ is the total kinetic energy associated with the stars in the galaxy). For the galaxy model considered in this paper, the revolution period relative to a circular orbit at the periphery (at $`R`$) is $`11t_{cr}`$. To return to physical units, we may refer to the case where $`M=210^{11}M_{}`$, $`R=10`$ kpc, so that $`t_{cr}=0.1810^8`$ yrs and the revolution period at the periphery is $`210^8`$ yrs.

### 4.2 Conservation laws

In the course of the simulation the total energy and the total angular momentum should remain constant at their initial values $`E_{tot}`$ and $`J_{tot}`$.

In the case of a single satellite we write the energy conservation as

$$E_{tot}=E_{gal}+K_{sat}+E_{int}=$$
(8)
$$\underset{i=1}{\overset{i=N}{}}\frac{m}{2}\left[v_i^2+\mathrm{\Phi }(r_i)\right]+\frac{M_s}{2}v_s^2\left[\underset{i=1}{\overset{i=N}{}}\frac{GM_sm}{(R_s^2+(r_ir_s)^2)^{1/2}}\right],$$

where $`E_{gal}=K_{gal}+W_{gal}`$ is the total (kinetic plus gravitational) energy of the galaxy in the absence of the satellite. From the point of view of the satellite, it is natural to introduce the energy $`E_{sat}=K_{sat}+E_{int}`$, so that we may split the total energy into $`E_{tot}=E_{gal}+E_{sat}`$. The angular momentum conservation is given simply by $`J_{tot}=J_{gal}+J_{sat}=_{i=1}^{i=N}mr_i\times v_i+M_sr_s\times v_s`$.

In the case of several mini-satellites, following the model outlined in Sect. 3, the energy takes the form

$$E_{tot}=E_{gal}+(K_{fra}+W_{fra})+E_{int}=\underset{i=1}{\overset{i=N}{}}\frac{m}{2}\left[v_i^2+\mathrm{\Phi }(r_i)\right]+$$
(9)
$$\underset{j=1}{\overset{j=N_f}{}}\frac{M_f}{2}\left[v_{fj}^2\underset{kj}{\overset{k=N_f}{}}\frac{GM_f}{(R_f^2+(r_{fj}r_{fk})^2)^{1/2}}\right]$$

$$\underset{i=1}{\overset{i=N}{}}\underset{j=1}{\overset{j=N_f}{}}\frac{GmM_f}{(R_f^2+(r_ir_{fj})^2)^{1/2}},$$

while the angular momentum can be written as $`J_{tot}=J_{gal}+J_{fra}=_{i=1}^{i=N}mr_i\times v_i+_{j=1}^{j=N_f}M_fr_{fj}\times v_{fj}`$. Note that the energy associated with the system of mini-satellites now includes a term $`W_{fra}`$ that describes the gravitational energy internal to the system of fragments. By analogy with the case of the single satellite, we may wish to refer to the energy $`E_{fra}=(K_{fra}+W_{fra})+E_{int}`$ as the total energy associated with the fragments, so that $`E_{tot}=E_{gal}+E_{fra}`$.

### 4.3 Density and pressure distributions

In order to describe the galaxy density perturbation, we refer to the density response defined as $`\rho _1(t,r)=\rho (t,rr_0(t))\rho _0(r)`$; here $`\rho _0(r)`$ is the adopted density distribution of the galaxy in the absence of satellites, $`r_0(t)`$ is the position of the center of mass of the galaxy in an inertial frame of reference. To diagnose the structure of the response it will be useful to consider the characteristics of the first multipoles.

In the simulations the density on the spherical mesh is expanded in associated Legendre functions $`P_l^m(\mathrm{cos}\theta )`$ up to $`l=6`$

$$\rho (r,\theta ,\phi )=\underset{l=0}{\overset{l=6}{}}[A_{l,0}(r)Y_l^0(\mathrm{cos}\theta )+\underset{m=1}{\overset{m=l}{}}\sqrt{2}A_{l,m}(r)\mathrm{Re}(Y_l^m(\mathrm{cos}\theta ,\phi ))$$
(10)
$$+\underset{m=1}{\overset{m=l}{}}\sqrt{2}B_{l,m}(r)\mathrm{Im}(Y_l^m(\mathrm{cos}\theta ,\phi ))],$$

where the coefficients of the expansion into real spherical harmonics are given by suitable integrations over the relevant solid angle

$$A_{l,m}(r)=\sqrt{2}\rho (r,\theta ,\phi )\mathrm{Re}(Y_l^m(\mathrm{cos}\theta ,\phi ))𝑑\mathrm{\Omega },$$
(11)
$$A_{l,0}(r)=\rho (r,\theta ,\phi )Y_l^0(\mathrm{cos}\theta )𝑑\mathrm{\Omega },$$
(12)
$$B_{l,m}(r)=\sqrt{2}\rho (r,\theta ,\phi )\mathrm{Im}(Y_l^m(\mathrm{cos}\theta ,\phi ))𝑑\mathrm{\Omega }.$$
(13)
In particular, the dipole ($`l=1`$) component of the density response $`\rho _1`$ is given by

$$\rho _{l=1}(r,\theta ,\phi )=A_{1,0}\sqrt{3/4\pi }\mathrm{cos}\theta A_{1,1}\sqrt{3/4\pi }\mathrm{sin}\theta \mathrm{cos}\phi B_{1,1}\sqrt{3/4\pi }\mathrm{sin}\theta \mathrm{sin}\phi .$$
(14)
The galaxy pressure tensor is defined in terms of the star distribution function in the standard way. In the case of a single satellite sinking into a galaxy along the equatorial plane we have $`v_r0`$ and $`v_\theta 0.`$

### 4.4 Tests

We have performed many preliminary experiments aimed at testing the numerical code.

The numerical implementation of the polytropic equilibrium configuration of Sect. 3 with $`N=510^3÷510^5`$ particles has been tested by generating an unperturbed galaxy that has remained quiescent for about 50 crossing times. The radii of the spheres containing $`0.5\%,1\%,1.5\%,\mathrm{}.,5\%`$ of the total mass have been checked to remain constant in time to within $`1\%`$.

Then, we have tested that a satellite with $`R_s=0.1R`$ and very small mass, $`M_s=10^8M`$, behaves as a light “test particle” without disturbing the main galaxy. This satellite, placed initially at the periphery of the galaxy $`r_s=R`$ on an orbit with $`v_{s\phi }=0.536R/t_{cr}`$ and $`v_{sr}=0.95410^3R/t_{cr}`$, remains very close to its initial circular path, never exceeding an eccentricity of $`e=0.05`$. At the end of the run, at time $`t50t_{cr}`$ the orbital energy of this light satellite is conserved to within $`1\%`$.

#### 4.4.1 Energy and momentum conservation

In addition to the specific tests mentioned above, we have checked that globally energy and angular momentum are conserved, at $`t50t_{cr}`$, to better than $`1\%`$.

#### 4.4.2 Comparison with earlier work by Bontekoe and van Albada (1987)

The main results of the simulations of the sinking of a single satellite by Bontekoe & van Albada (1987) have been checked in detail, with the use of simulations with $`N`$ up to $`510^5`$ (see Sect. 5.1 below). This adds much confidence, in relation both to the overall structure of the code used in this paper and to the modeling of dynamical friction used in these investigations (since results have been checked to be largely independent of $`N`$ over a range significantly larger than that available earlier).

## 5 Results of the simulations

We have performed a number of runs, mostly based on $`N=510^5`$. Runs labeled by $`A`$ refer to experiments with a single sinking satellite, initially located on a circular orbit at $`r_s=R`$; various combinations of satellite mass $`M_s/M`$ and satellite size $`R_s/R`$ have been considered. For these runs (some are listed in Table 1), we have recorded the time $`t_{fall}/t_{cr}`$ at which the satellite reaches the central regions; for some experiments (in particular, for those labeled $`A_8`$ and $`A_9`$) the satellite is unable to reach the center and we have thus studied the radius of minimum radial approach $`r_{min}`$.

Runs labeled by $`B`$ are listed in Table 2 and refer to the case of the slow evolution of the system in the presence of many fragments, initially located in a single shell, following the procedure of case I described in Sect. 3.2.1. Runs of type $`BS`$ are based on the smoother initialization of case II outlined in Sect. 3.2.2. Runs of type $`BT`$ refer to the procedure of case III of Sect. 3.2.3. Here the fourth column records the initial radius of the shell and the fifth column lists the time at the end of the experiment; all runs are made with the mutual interaction, among fragments, on, but we have also performed a few runs for which the interaction is turned off. Runs of type $`C`$ involve two shells of fragments.

### 5.1 The case of one sinking satellite

We have studied the energy transfer between the galaxy and the sinking satellite and the overall energy balance, using the total energy conservation as a diagnostic of the code accuracy. In particular we have considered the quantity $`\mathrm{\Delta }E_{sat}=K_{sat}+E_{int}K_{sat}(t=0)E_{int}(t=0)`$. For example, for run $`A_1`$ (characterized by $`M_s/M=0.1`$, $`R_s/R=0.1`$, $`t_{fall}/t_{cr}=48.125`$) the energy lost by the satellite is initially very small (on the order of $`0.02`$, in the units following the conventions described in Sect. 4.1, at time $`t38.5t_{cr}`$), and then rapidly rises (to $`0.075`$ at time $`t49.5t_{cr}`$); during the run the total energy is conserved to better than $`0.002`$. Similarly, we have considered the angular momentum balance. Here the relatively large amount of angular momentum already associated with the galaxy right after the beginning of the simulation ($`J_{gal}0.0054`$ at $`t5.5t_{cr}`$) is not yet due to the scattering of star orbits by the sinking satellite, but is mostly related to the overall orbital motion of the galaxy in the inertial frame of reference of the center of mass; eventually, the satellite is dragged in, down to the center, and loses completely its angular momentum, which is acquired by the galaxy.

During the process, following the steps indicated by Bontekoe & van Albada (1987), we can “measure” the Coulomb logarithm. The measurement of $`\mathrm{ln}\lambda `$ consists of the following steps: (i) Determine the velocity $`\widehat{v}_s`$ of the satellite with respect to the local streaming motion of the stars (the streaming motion is determined as a mean velocity of the stars in the “vicinity” of the satellite); (ii) Measure the energy $`E_{sat}`$ of the satellite at many instants, then determine $`dE_{sat}/dt;`$ (iii) Determine the density $`\rho _G`$ of the galaxy at the position of the satellite by linear interpolation of densities between the 8 nearest mesh corners; (iv) Calculate $`F(\widehat{v}_s)`$, the fraction of particles with velocity smaller than $`\widehat{v}_s`$. Finally the Coulomb logarithm is estimated from the relation $`4\pi G^2\mathrm{ln}\lambda =v_s(dE_{sat}/dt)/(M_s^2\rho _GF(\widehat{v}_s))`$.

The main results that we have obtained by examining the experiments for the case of a single sinking satellite are the following:

* Qualitatively, the general conclusions derived by Bontekoe & van Albada (1987) have been checked to hold. However, we have found that, while the classical formula of dynamical friction can be used to fit the orbital decay of the satellite over a significant time interval, the determined value of the Coulomb logarithm changes significantly in the course of the simulation (see Fig. 1); during the initial transient and the final part of the orbital evolution, when the satellite reaches the central regions, the change in the nominal value of the Coulomb logarithm can be dramatic, as already noticed by Bontekoe & van Albada (1987).
* The process of dynamical friction has been checked to operate via the generation of wakes, of which the dipolar $`l=1`$ component appears to play a dominant role (see Figs. 2 and 3). This generally agrees with the scenario outlined by Weinberg (1989); see also Sect. 2.
* As a result of the scattering of stars due to the sinking satellite, the galaxy acquires not only some angular momentum (the induced systematic motions are approximately those of solid body rotation, see Fig. 4), but also a significantly anisotropic pressure tensor. Figure 5, relative to run $`A_1`$, suggests that the final configuration reached is reasonably well characterized by $`p_{rr}p_{\theta \theta }`$, while the toroidal pressure $`p_{\phi \phi }`$ is significantly larger even in the central regions. As a result, except for the outer parts, the evolved distribution function of the galaxy may turn out to be accessible to models with $`f=f(E,J_z)`$.
* In some runs, for which the diameter of the satellite is reduced, the satellite turns out to fall inwards faster, but to be unable to reach the center of the galaxy (see Fig. 6). We have checked that this effect does not depend on the number N by running a simulation of type $`A_9`$ with N up to $`210^6`$ (in such a way that the number of simulation particles in the region inside the radius of minimum approach increases from about 600 to about 2500; correspondingly, the observed radius of minimum approach changes by less than 1 %). We should recall that simulations with a code of the type we are using, with a finite expansion in spherical harmonics, cannot deal realistically with spatially small satellites (see comment at the end of Sect. 2.1). Still, a “saturation” of dynamical friction of the kind noted in our simulations is likely to be real, because the modeling problem associated with the spherical harmonic expansion is less severe for the late phases of the simulation, when the satellite is close to the galaxy center. The ultimate fate of a single sinking object might be related to other interesting issues, such as those met in the discussion of the formation of massive black hole binaries (e.g., see Quinlan & Hernquist qui97 (1997), Milosavljevich & Merritt mil01 (2001)), which go well beyond the scope of the present paper (that addresses the slow evolution of a galaxy on the large scale as a result of the capture of a large number of small objects); what we have noted here is just meant to say that on this specific feature we confirm the results of previous investigations (in particular, those of Bontekoe & van Albada 1987) carried out with a much smaller number of particles.
* Finally we note that the overall density distribution of the galaxy, at the end of the simulation, turns out to be less concentrated with respect to the initial distribution (see Fig. 7). In passing, we note that a figure such as Fig. 7 naively suggests that the total mass is not conserved, but we have checked that in this respect this is only a false impression left by our way of representing the density profile; in particular, one should recall that the weight of the volume element rapidly increases with radius.

### 5.2 The case of many fragments

#### 5.2.1 Case I

For various runs of type $`B`$ and $`C`$ we have checked that the total energy conservation generally holds, over a period of $`100t_{cr}`$, to better than one part over 500. Yet, it is difficult to estimate the impact of dynamical friction on the galaxy evolution, because secular effects are masked by effects related to the initial lack of equilibrium. From this class of numerical experiments we have obtained the following results:

* We have observed a process of shell diffusion, which is a collective effect qualitatively different from what might have been imagined by simply superposing the effects of orbital decay of the individual fragments. The phenomenon persists, presumably via interaction with the wakes induced by the fragments in the galaxy, even when the direct fragment-fragment interaction is turned off. Sharp coherent oscillations are noted at early times; they are recognized as epicyclic oscillations generated by our choice of initial velocities for the fragments. We have also performed tests on the process of shell diffusion and oscillations on the basis of 400 fragments initially located close to $`r=0.6R`$.
* We have noted systematic changes of the galaxy density distribution. The general trend emerging from the simulations is that significant peaking of the density distribution may occur, provided the initial radius of one shell of fragments is well inside the galaxy. However, much of the change is apparently associated with the initial lack of equilibrium.
* The process of density profile modification is not accompanied by the development of significant pressure anisotropy.

In conclusion, simulations of Case I turn out to be of little direct use for the main objectives of our project. However, they are extremely useful and instructive, because they demonstrate that, in order to properly single out slow relaxation processes, one has to be very careful about the choice of initial conditions and, in particular, because they provide concrete examples of some misleading effects that can occur as a result of the use of improper models.

#### 5.2.2 Case II

The use of the smoother initialization described as case II in Sect. 3.2.2 leads to runs where the secular effects associated with the granularity of the fragment population are more convincingly identified. The orbital decay of the fragments is illustrated in Fig. 8, where we also show for comparison a parallel run where the shell is treated as collisionless. The interaction which generates dynamical friction of the galaxy on the fragments is responsible for the slow, but systematic, development of a tangential bias in the pressure tensor, as illustrated in Fig. 9 and for an overall change in the density profile of the galaxy (see Fig. 10).

The check against the collisionless run is worth a special digression. In order to diagnose the collisionality present in our system better, we have performed the following test. We have studied the correlation between single particle energy (for the particles simulating the galaxy) at two different times. If the simulation is truly collisionless, the single particle energy $`E_p=mv^2/2+m\mathrm{\Phi }(r)`$, where $`\mathrm{\Phi }(r)`$ is the mean field potential, should be strictly conserved. Indeed, by means of a plot similar to the lower right frame of Fig. 8, we have checked that for the run illustrated by the lower left frame of Fig. 8 more than $`90\%`$ of the particles have an energy correlation such that $`|[E_p(t=118t_{cr})E_p(t=14t_{cr})]/E_p(t=14t_{cr})|<0.06`$. In turn, for run $`BS_2`$, otherwise similar, but characterized by the presence of 100 fragments, a similar level of correlation for single particle energy (at 6 %) is preserved within the same time interval only by 50 % of the particles; the lower right frame of Fig. 8 thus illustrates the “damage” associated with dynamical friction.

#### 5.2.3 Case III

The even smoother procedure of case III leads to slower orbital decays, as illustrated in the lower frames of Fig. 11. For this type of simulations, the initial distribution of fragment orbits replicates the isotropic star distribution of the underlying galaxy. Because of dynamical friction, the fragment orbits are subject to a slow process of circularization, as illustrated in the upper frames of Fig. 11.

For each of the fragments that are involved in the simulation described as Case II, we can repeat the analysis that has led to the measurement of the Coulomb logarithm from the study of the decay of the orbit of a single satellite. For example, for run $`BS_2`$ we have 100 independent such determinations. In Fig. 12 we record, for run $`BS_2`$ a measurement of an average value of the Coulomb logarithm. We note that the fact that the value of the Coulomb logarithm thus determined remains of the same order of magnitude as that obtained in Fig. 1 proves that the classical expression for the effects of dynamical friction captures the correct scaling with respect to the mass of the object subject to friction.

## 6 Discussion and conclusions

In this paper we have set up a numerical laboratory for the study of dynamical friction in inhomogeneous quasi-spherical galaxies and made several experiments that have elucidated the physical mechanisms that participate in the process.

The general character of the simulations presented here is distinctly different from that adopted in previous studies in this general research area, because we have addressed a physical scenario in which the slow collapse induced by dynamical friction approximately preserves the spherical symmetry present initially. This configuration is inspired by the astrophysical problem of studying the evolution of a spherical galaxy in the presence of a substantial globular cluster system or as a result of the capture of a large set of mini-satellites dragged in along random directions.

From the theoretical point of view, the framework adopted here has an important advantage with respect to the traditional case of the study of the decay of the orbit of a single satellite, because it allows us to run smoother simulations that are basically free from other effects unrelated to dynamical friction, such as those associated with lack of equilibrium in the initial configuration.

In the past, most of the attention has been focused on the effect of the host galaxy on the object subject to dynamical friction. Here, we have been able not only to test the adequacy of the classical formulae of dynamical friction (derived for ideal homogeneous models) in the inhomogeneous galaxy environment, but also to measure the evolution of the density and pressure tensor profiles in the galaxy induced by the presence of friction processes on a minority population of heavier particles. We have made quantitative tests that convince us that we are detecting genuine effects associated with secular evolution.

### 6.1 Plans for future investigations

Starting from the encouraging quantitative results obtained in this paper, there are now two major issues that we would like to address in investigations planned for the near future, which appear to be approachable by the use of the numerical laboratory developed here:

* In this paper we have limited our study to the presence of one or two shells of heavy particles, so as to probe the reaction of the galaxy distribution function in a controlled way and to test whether the interaction involved can be described as “local”. The first immediate objective that we have in mind is to devise a set of more realistic simulations, able to represent the distribution of the globular cluster system in galaxies better.
* So far we have limited our discussion to the study of the orbital decay in an $`n=3`$ polytropic, isotropic galaxy model. Now we would like to test systematically the process of slow collapse induced by dynamical friction as a function of the concentration properties of the equilibrium model, for the case of more realistic galaxy models.

###### Acknowledgements.

Part of this work was supported by the Italian MIUR. We thank Tjeerd van Albada for providing us with a copy of his code and for a number of useful conversations.


format: text

text:
# An XML-based document suite

## Introduction

We have designed and implemented the XDOC document suite as a workbench for the flexible processing of electronically available documents in German. We have decided to exploit XML \[xml-standard\] and its accompanying formalisms (e.g. XSLT \[XSL\]) and tools (e.g. xt \[clarkSite\] ) as a unifying framework. All modules in the XDOC system expect XML documents as input and deliver their results in XML format.

XML – and ist precursor SGML – offers a formalism to annotate pieces of (natural language) texts. To be more precise: If a text is (as a simple first approximation) seen as a sequence of characters (alphabetic and whitespace characters) then XML allows to associate arbitrary markup with arbitrary subsequences of contiguous characters. Many linguistic units of interest are represented by strings of contiguous characters (e.g. words, phrases, clauses etc.). To use XML to encode information about such a substring of a text interpreted as a meaningful linguistic unit and to associate this information directly with the occurrence of the unit in the text is a straightforward idea. The basic idea is further backed by XMLs demand that XML elements have to be properly nested. This is fully concordant with standard linguistic practice: complex structures are made up from simpler structures covering substrings of the full string in a nested way.

The end users of our applications are domain experts (e.g. medical doctors, engineers, …). They are interested in getting their problems solved but they are typically neither interested nor trained in computational linguistics. Therefore the barrier to overcome before they can use a computational linguistics or text technology system should be as low as possible.

This experience has consequences for the design of the document suite. The work in the XDOC project is guided by the following design principles that have been abstracted from a number of experiments and applications with ”realistic” documents (i.a. emails, abstracts of scientific papers, technical documentation, …):

* The tools shall be usable for ‘realistic’ documents.
One aspect of ‘realistic’ documents is that they typically contain domain-specific tokens that are not directly covered by classical lexical categories (like noun, verb, …). Those tokens are nevertheless often essential for the user of the document (e.g. an enzyme descriptor like EC 4.1.1.17 for a biochemist).
* The tools shall be as robust as possible.
In general it can not be expected that lexicon information is available for all tokens in a document. This is not only the case for most tokens from ‘nonlexical’ types – like telephone numbers, enzyme names, material codes, … –, even for lexical types there will always be ‘lexical gaps’. This may either be caused by neologisms or simply by starting to process documents from a new application domain with a new sublanguage. In the latter case lexical items will typically be missing in the lexicon (‘lexical gap’) and phrasal structures may not or not adequately be covered by the grammar.
* The tools shall be usable independently but shall allow for flexible combination and interoperability.
* The tools shall not only be usable by developers but as well by domain experts without linguistic training.

Here again XML and XSLT play a major role: XSL stylesheets can be exploited to allow different presentations of internal data and results for different target groups; for end users the internals are in many cases not helpful, whereas developers will need them for debugging.

The tools in the XDOC document suite can be grouped according to their function:

* preprocessing
* structure detection
* POS tagging
* syntactic parsing
* semantic analysis
* tools for the specific application: e.g. information extraction

In all tools the results of processing is encoded with XML tags delimiting the respective piece of text. The information conveyed by the tag name is enriched with XML attributes and their resp. values.

## Preprocessing

Tools for preprocessing are used to convert documents from a number of formats into the XML format amenable for further processing. As a subtask this includes treatment of special characters (e.g. for umlauts, apostrophes, …).

## Structure detection

We accept raw ASCII texts without any markup as input. In such cases structure detection tries to uncover linguistic units (e.g. sentences, titles, …) as candidates for further analysis. A major subtask is to identify the role of interpunction characters.

If we have the structures in a text explicitly available this may be exploited by subsequent linguistic processing. An example: For a unit classified as title or subtitle you will accept an NP whereas within a paragraph you will expect full sentences.

In realistic texts even the detection of possible sentence boundaries needs some care. A period character may not only be used as a full stop but may as well be part of an abbreviation (e.g. ‘z.B.’ – engl.: ‘e.g.’ – or ‘Dr.’), be contained in a number (3.14), be used in an email address or in domain specific tokens. The resources employed are special lexica (e.g. for abbreviations) and finite automata for the reliable detection of token from specialized non-lexical categories (e.g. enzyme names, material codes, …).

These resources are used here primarily to identify those full stop characters that function as sentence delimiters (tagged as IP). In addition, the information about the function of strings that include a period is tagged in the result (e.g. ABBR).

###### Example 1

results of structure detection

```
Anwesend<IP>:</IP>
<ABBR>Univ.-Prof.</ABBR>
<ABBR>Dr.</ABBR><ABBR>med.</ABBR>Dieter Krause<IP>,</IP>
Direktor des Institutes fuer Rechtsmedizin

```

## POS tagging

To try to assign part-of-speech information to a token is not only a preparatory step for parsing. The information gained about a document by POS tagging and evaluating its results is valuable in its own right. The ratio of token not classifiable by the POS tagger to token classified may e.g. serve as an indication of the degree of lexical coverage.

In principle a number of approaches is usable for POS tagging (e.g. \[brill:92\]). We decided to avoid approaches based on (supervised) learning from tagged corpora, since the cost for creating the necessary training data are likely to be prohibitive for our users (especially in specialized sublanguages).

The approach chosen was to try to make best use of available resources for German and to enhance them with additional functionality. The tool chosen is not only used in POS tagging but serves as a general morpho-syntactic component for German: MORPHIX.

The resources employed in XDOC’s POS tagger are:

\- the lexicon and the inflectional analysis from the morphosyntactic component MORPHIX

\- a number of heuristics (e.g. for the classification of token not covered in the lexicon)

For German the morphology component MORPHIX \[finkler.neumann:88\] has been developed in a number of projects and is available in different realisations. This component has the advantage that the closed class lexical items of German as well as all irregular verbs are covered. The coverage of open class lexical items is dependent on the amount of lexical coding. The paradigms for e.g. verb conjugation and noun declination are fully covered but to be able to analyze and generate word forms their roots need to be included in the MORPHIX lexicon.

We exploit MORPHIX - in addition to its role in syntactic parsing \- for POS tagging as well. If a token in a German text can be morphologically analysed with MORPHIX the resulting word class categorisation is used as POS information. Note that this classification need not be unique. Since the tokens are analysed in isolation multiple analyses are often the case. Some examples: the token ‘der’ may either be a determiner (with a number of different combinations for the features case, number and gender) or a relative pronoun, the token ‘liebe’ may be either a verb or an adjective (again with different feature combinations not relevant for POS tagging).

In addition since we do not expect extensive lexicon coding at the beginning of an XDOC application some tokens will not get a MORPHIX analysis. We then employ two techniques: We first try to make use of heuristics that are based on aspects of the tokens that can easily be detected with simple string analysis (e.g. upper-/lowercase, endings, …) and/or exploitation of the token position relative to sentence boundaries (detected in the structure detection module). If a heuristic yields a classification the resulting POS class is added together with the name of the employed heuristic (marked as feature SRC, cf. example 3). If no heuristics are applicable we classify the token as member of the class unknown (tagged with XXX).

To keep the POS tagger fast and simple the disambiguation between multiple POS classes for a token and the derivation of a possible POS class from context for an unknown token are postponed to syntactic processing. This is in line with our general principle to accept results with overgeneration when a module is applied in isolation (here: POS tagging) and to rely on filtering ambiguous results in a later stage of processing (here: exploiting the syntactic context).

###### Example 2

domain-specific tagging

```

<PRODUCT Method="Sandguss" Material="CC333G">
    <N>Gussstueck</N>
    <NORM>
         <N>EN</N>
         <NR>1982</NR>
    </NORM>
    <IP>-</IP>
    <MAT-ID>CC333G</MAT-ID>
    <IP>-</IP>
    <METHODE>GS</METHODE>
    <IP>-</IP>
    <MODELLNR>XXXX</MODELLNR>
</PRODUCT>

```

The example above is the result of tagging a domain-specific identifier. The token is annotated as a PRODUCT with description of the used method and material. It is a typical token in the domain of casting technology.

## Syntactic parsing

For syntactic parsing we apply a chart parser based on context free grammar rules augmented with feature structures.

Again robustness is achieved by allowing as input elements:

* multiple POS classes,
* unknown classes of open world tokens and
* tokens with POS class, but without or only partial feature information.

###### Example 3

unknown token classified as noun with heuristics

```
<NP TYPE="COMPLEX" RULE="NPC3" GEN="FEM"
        NUM="PL" CAS="_">
  <NP TYPE="FULL" RULE="NP1" CAS="_"
            NUM="PL" GEN="FEM">
       <N SRC="UNG">Blutanhaftungen</N>
  </NP>
  <PP CAS="DAT">
    <PRP CAS="DAT">an</PRP>
    <NP TYPE="FULL" RULE="NP2" CAS="DAT"
            NUM="SG" GEN="FEM">
      <DETD>der</DETD>
      <N SRC="UC1">Gekroesewurzel</N>
    </NP>
  </PP>
</NP>

```

The latter case results from some heuristics in POS tagging that allow to assume e.g. the class noun for a token but do not suffice to detect its full paradigm from the token (note that there are ca two dozen different morphosyntactic paradigms for noun declination in German).

For a given input the parser attempts to find all complete analyses that cover the input. If no such complete analysis is achievable it is attempted to combine maximal partial results into structures covering the whole input.

A successful analysis may be based on an assumption about the word class of an initially unclassified token (tagged XXX). This is indicated in the parsing result (feature AS) and can be exploited for learning such classifications from contextual constraints. In a similar way the successful combination from known feature values from closed class items (e.g. determiners, prepositions) with underspecified features in agreement constraints allows the determination of paradigm information from successfully processed occurrences. In example 4 features of the unknown word ”Mundhoehle” could be derived from the features of the determiner within the PP.

###### Example 4

unknown token classified as adjective and features derived through contextual constraints

```
<NP TYPE="COMPLEX" RULE="NPC3" GEN="MAS" NUM="SG"
    CAS="NOM">
  <NP TYPE="FULL" RULE="NP3" CAS="NOM" NUM="SG"
        GEN="MAS">
    <DETI>kein</DETI>
    <XXX AS="ADJ">ungehoeriger</XXX>
    <N>Inhalt</N>
  </NP>
  <PP CAS="DAT">
    <PRP CAS="DAT">in</PRP>
    <NP TYPE="FULL" RULE="NP2" CAS="DAT" NUM="SG"
        GEN="FEM">
      <DETD>der</DETD>
      <N SRC="UC1">Mundhoehle</N>
    </NP>
  </PP>
</NP>"

```

The grammar used in syntactic parsing is organized in a modular way that allows to add or remove groups of rules. This is exploited when the sublanguage of a domain contains linguistic structures that are unusual or even ungrammatical in standard German.

###### Example 5

Excerpt from syntactic analysis

```
<PP CAS="AKK">
  <PRP CAS="AKK">durch</PRP>
    <NP TYPE="COMPLEX" RULE="NPC1" GEN="NTR" NUM="SG"
                                 CAS="AKK">
      <NP TYPE="FULL" RULE="NP1" CAS="AKK" NUM="SG"
                                 GEN="NTR">
        <N>Schaffen</N>
      </NP>
      <NP TYPE="FULL" RULE="NP2" CAS="GEN" NUM="SG"
                                  GEN="MAS">
         <DETD>des</DETD>
         <N>Zusammenhalts</N>
      </NP>
    </NP>
</PP>

```

## Semantic analysis

At the time of writing semantic analysis uses three methods:

### Semantic tagging

For semantic tagging we apply a semantic lexicon. This lexicon contains the semantic interpretation of a token and a case frame combined with the syntactic valence requirements. Similar to POS tagging the tokens are annotated with their meaning and a classification in semantic categories like e.g. concepts and relations. Again it is possible, that the classification of a token in isolation is not unique. Multiple classification can be resolved through the following analysis of the case frame and through its combination with the syntactic structure which includes the token.

### Analysis of case frames

By the case frame analysis of a token we obtain details about the type of recognized concepts (resolving multiple interpretations) and possible relations to other concepts. The results are tagged with XML tags. The following example describes the DTD for the annotation of the results of case frame analysis.

###### Example 6

DTD for the annotation by case frame analysis

```
    <!ELEMENT CONCEPTS (CONCEPT)*>

    <!ELEMENT CONCEPT (WORD, DESC, SLOTS?)>
    <!ATTLIST CONCEPT TYPE CDATA #REQUIRED>

    <!ELEMENT WORD (#PCDATA)>
    <!ELEMENT DESC (#PCDATA)>
    <!ELEMENT SLOTS (RELATION+)>

    <!ELEMENT RELATION (ASSIGN_TO, FORM, CONTENT)>
    <!ATTLIST RELATION TYPE CDATA #REQUIRED>

    <!ELEMENT ASSIGN_TO (#PCDATA)>
    <!ELEMENT FORM (#PCDATA)>
    <!ELEMENT CONTENT (#PCDATA)>

```

We use attributes to show the description of the concepts and we can annotate the relevant relations between the concepts through nested tags (e.g. the tag *SLOTS*).

###### Example 7

Excerpt from case frame analysis

```
 <CONCEPT TYPE=Prozess>
    <WORD>Fertigen</WORD>
    <DESC>Schaffung von etwas</DESC>
    <SLOTS>
        <RELATION>
        <RESULT FORM="N(gen, fak) P(akk, fak, von)">
                fester Koerper</RESULT>
        <SOURCE FORM="P(dat, fak, aus)">aus formlosem
                Stoff </SOURCE>
        <INSTRUMENT FORM="P(akk, fak, durch)">durch
                Schaffen des Zusammenhalts</INSTRUMENT>
        </RELATION>
    </SLOTS>
</CONCEPT>

```

The example above is part of the result of the analysis of the German phrase: Fertigen fester Koerper aus formlosem Stoff durch Schaffen des Zusammenhalts111In English: production of solid objects from formless matter by creating cohesion. The token Fertigen is classified as process with the relations source, result and instrument. The following phrases (noun phrases and preposition phrases) are checked to make sure that they are assignable to the relation requirements (semantic and syntactic) of the token Fertigen.

### Semantic interpretation of the syntactic structure

An other step to analyze the relations between tokens can be the interpretation of the syntactic structure of a phrase or sentences respectively. We exploit the syntactic structure of the sublanguage to extract the relation between several tokens. For example a typical phrase from an autopsy report: Leber dunkelrot.222In English: Liver dark red.

From semantic tagging we obtain the following information:

###### Example 8

results of semantic tagging

```
<CONCEPT TYPE="organ">Leber</CONCEPT>
<PROPERTY TYPE="color">dunkelrot</PROPERTY>
<XXX>.</XXX>

```

In this example we can extract the relation ”has-color” between the tokens Leber and dunkelrot. This is an example of a simple semantic relation. Other semantic relations can be described through more complex variations. In these cases we must consider linguistic structures like modifiers (e.g. *etwas*), negations (e.g. *nicht*), coordinations (e.g. *Beckengeruest unversehrt und fest gefuegt*) and noun groups (e.g. *Bauchteil der grossen Koerperschlagader*).

## Current state and future work

The XDOC document workbench is currently employed in a number of applications. These include:

* knowledge acquisition from technical documentation about casting technology
* extraction of company profiles from WWW pages
* analysis of autopsy protocols

The latter application is part of a joint project with the institute for forensic medicine of our university. The medical doctors there are interested in tools that help them to exploit their huge collection of several thousand autopsy protocols for their research interests. The confrontation with this corpus has stimulated experiments with ‘bootstrapping techniques’ for lexicon and ontology creation.

The core idea is the following:

When you are confronted with a new corpus from a new domain, try to find linguistic structures in the text that are easy to detect automatically and that allow to classify unknown terms in a robust manner both syntactically as well as on the knowledge level. Take the results from a run of these simple but robust heuristics as an initial version of a domain dependent lexicon and ontology. Exploit these initial resources to extend the processing to more complicated linguistic structures in order to detect and classify more terms of interest automatically.

An example: In the sublanguage of autopsy protocols (in German) a very telegrammatic style is dominant. Condensed and compact structures like the following are very frequent:

> *Harnblase leer.*
> *Harnleiter frei.*
> *Nierenoberflaeche glatt.*
> *Vorsteherdruese altersentsprechend.*

These structures can be abstracted syntactically as $`<`$Noun$`>`$$`<`$Adjective$`>`$$`<`$Fullstop$`>`$ and semantically as reporting a finding in the form $`<`$Anatomic-entity$`>`$ has $`<`$Attribute-value$`>`$ and they are easily detectable \[roesner02\].

In our experiments we have exploited this characteristic of the corpus extensively to automatically deduce an initial lexicon (with nouns and adjectives) and ontology (with concepts for anatomic regions or organs and their respective features and values). The feature values were further exploited to cluster the concept candidates into groups according to their feature values. In this way container like entities with feature values like ‘leer’ (empty) or ‘gefuellt’ (full) can be distinguished from e.g. entities of surface type with feature values like ‘glatt’ (smooth).

## Related Work

The work in XDOC has been inspired by a number of precursory projects:

In GATE \[GATESite, GATE\] the idea of piping simple modules in order to achieve complex functionality has been applied to NLP with such a rigid architecture for the first time. The project LT XML has been pioneering XML as a data format for linguistic processing.

Both GATE and LT XML \[ltxml99\] were employed for processing English texts. SMES \[Neumann97\] has been an attempt to develop a toolbox for message extraction from German texts. A disadvantage of SMES that is avoided in XDOC is the lack of a uniform encoding formalism, in other words, users are confronted with different encodings and formats in each module.

## System availability

Major components of XDOC are made publicly accessible for testing and experiments under the URL:

http://lima.cs.uni-magdeburg.de:8000/

## Summary

We have reported about the current state of the XDOC document suite. This collection of tools for the flexible and robust processing of documents in German is based on the use of XML as unifying formalism for encoding input and output data as well as process information. It is organized in modules with limited responsibilities that can easily be combined into pipelines to solve complex tasks. Strong emphasis is laid on a number of techniques to deal with lexical and conceptual gaps and to guarantee robust systems behaviour without the need for a priori investment in resource creation by users. When end users are first confronted with the system they typically are interested in quick progress in their application but should not be forced to be engaged e.g. in lexicon build up and grammar debugging, before being able to start with experiments. This is not to say that creation of specialized lexicons is unnecessary. There is a strong correlation between prior investment in resources and improved performance and higher quality of results. Our experience shows that initial results in experiments are a good motivation for subsequent efforts of users and investment in extended and improved linguistic resources but that a priori costs may be blocking the willingness of users to get really involved.


text:
# References

Noncommutative deformation and

a topological nature of singularity Koiter

Trinh V. K.

Department of physics of polymer and crystal, MGU

E-mail : khoa@polly.phys.msu.ru

1. From the point of view of phenomenology, till now there are many models of plastic deformation. Almost in all monographs , for example \[1-3,10\], the following models are considered: sliding, regular, current, singular, analytical and deformation plasticity. The determinant ratio in deformation are a major problem, which was submitted for discussion. The its most simple kind was presented in the form of Genki-Nadai:

$$\sigma _{ij}=2G_s\left(T\right)e_{ij}\text{at}T=\overline{T}$$
$`\left(1\right)`$
$$d\sigma _{ij}=2Ge_{ij}\text{at}T=\overline{T}$$

Where $`d\sigma _{ij},e_{ij}`$ \- deviatory of a stress and deformation, $`G_s`$\- current module of pure shift, $`G`$\- the elastic module of shift, $`T`$\- intensity tangent stress, $`\overline{T}`$\- the maximal value for a history of loading pre-eminent some value $`T_0`$. In singularity model it looks like

$$de_{ij}=\underset{n=1}{\overset{N}{}}df_n\frac{f_n}{\sigma _{ij}}$$
$`\left(2\right)`$
Where $`f_n`$\- the current limiting surface. Here not zero are those $`df_n`$, for which on the appropriate surface $`f_n`$ be performed loading, i.e.

$$\left(\frac{f_n}{\sigma _{ij}}\right)d\sigma _{ij}0$$

In model Iliusin it has the kind

$$\sigma _{ij}=\underset{n}{}A_n\frac{d^ne_{ij}}{de^n}:e=\sqrt{de_\alpha de_\alpha }$$
$`\left(3\right)`$
However, if we use concept of potential of deformation as in \[1-3,10-11\], on definition of function $`H`$, we shall receive following determining ratio

$$\sigma _{ij}=\frac{H}{e_{ij}}$$
$`\left(4\right)`$
In our model , in which we consider process of plastic deformation as process of phase transition, potential of deformation (or the free energy) plays role of function Hamilton or action of system. So, if we admit Hamilton as $`H=\alpha e^2+\beta e^4+\mathrm{},`$ $`\alpha ,\beta `$\- factors of elasticity, - deformation tensor, that as a first approximation we shall receive the formulas (1) (2) (3) (if $`H`$ is expressed through the metrics of spaces). Roughly speaking, form (4) is general determining ratio for process of plastic deformation. Here there are two problems: how the fields $`\sigma ,e`$ in the time and in space is distributed ? What nature of the singularity Koiter ? We try to solve these problems consistently. For this purpose we should construct more general model of deformation.

2. As is known, all experimental data have shown, that in a plastic status the stress and deformation is strongly fluctuation. Besides the diagram $`\sigma =f(e)`$ for pressure - deformation is similar to the diagram $`P=f(V)`$ for pressure - volume for a liquid . So, it is possible to consider a process of deformation as process of phase transition, or more generalized speaking, it is process of formation of structure . Density of distribution of probability transition to plasticity is solution of the equation Fokker-Plank. However, if using potential character of process of deformation it is possible instead of the equation Fokker-Plank to apply the equation Schrodinger. In this approach tensor of deformation is chosen in the roles of parameter of the order. In work  A.A. Iliusin has presented the concept of trajectory of deformation in space $`E_5`$. As a matter of fact the space $`E_5`$ is the one fiber in the deformation bundle  Further, because of fluctuation of the statuses of deformation, it is necessary for me to used ” the secondary calculation ” to receive a field deformations in a plastic status from the elasticity. Nevertheless it is necessary to notice, that a nature of process of plastic deformation is noncommutative nature.

So, the process of deformation $`P=\{e_{ij}(x,t)\}`$ is given. That is noncommutative deformation space. The deformation wave $`E_{ij}(x,t)`$ is noncommutative wave. For complete descriptions of distribution of this wave we have to use the space being tensor product $`C(R_4)M_n`$, and $`C(R_4)`$ \- algebra of smooth function determined on usual space-time, and $`M_n`$ \- algebra $`n\times n`$ of matrixes. In this general space the matrix function $`E(x,t)`$ looks like (for 2-dimention case): $`E(x,t)=\left(\begin{array}{cccc}a& b& & \\ c& d& & \end{array}\right)C(R_4)M_2`$. Now our task is construction of noncommutative deformation bundle on the $`C(R_4)M_2`$. This problem was performed in frameworks of noncommutative geometry Connes . There are two important objects in this geometry. It is associative algebra $`𝒜`$ and universal algebra of the differential forms $`\mathrm{\Omega }_D^1(𝒜)`$ on the algebra $`𝒜`$. Here space of deformation $`P=C(R_4)M_n`$ plays role of associative algebra $`𝒜`$, at the same time universal algebra $`\mathrm{\Omega }_D^1(𝒜)`$ looks like :

$$\mathrm{\Omega }_D^1\left(𝒜\right)=\mathrm{\Omega }_D^1\left(C\left(R_4\right)\right)\mathrm{\Omega }_D^1\left(M_n\right)$$

It is tensor product of algebra of the differential forms above algebra of smooth function and algebra of the differential forms above algebra $`n\times n`$ of matrixes. In turn this algebra $`\mathrm{\Omega }_D^1(𝒜)`$ can be present through the direct sum of a horizontal and vertical of part

$$\mathrm{\Omega }_D^1\left(𝒜\right)=\mathrm{\Omega }_H^1\mathrm{\Omega }_V^1$$

and

$$\mathrm{\Omega }_H^1=\mathrm{\Omega }_D^1\left(C\left(R_4\right)\right)M_n;\mathrm{\Omega }_V^1=C\left(R_4\right)\mathrm{\Omega }_D^1\left(M_n\right)$$

Thus, differential calculation on algebra $`𝒜`$, being in usual external differential, is chosen. More precisely, if through $`Der(𝒜)`$ we had designated Lie algebra of derivative on $`𝒜`$, we shall receive

$$Der\left(𝒜\right)=\left(Der\left(C\left(R_4\right)\right)1\right)\left(C\left(R_4\right)Der\left(M_n\right)\right)$$

Differential $`df`$ of an element $`f𝒜`$ also is divided for the sum

$$df=d_Hf+d_Vf$$

where $`d_hf`$ and $`d_Vf`$ belong to $`\mathrm{\Omega }_H^1`$ and $`\mathrm{\Omega }_V^1`$ respectively. If $`E_k,k\{1,2,\mathrm{},n^21\}`$ is basis of algebra $`M_n`$, then $`_k=ad(iE_k)`$ is a basis of algebra $`Der(M_n)`$ and $`[_k,_l]=C_{klm}_m`$. So, the basis $`\theta ^k`$ of algebras $`\mathrm{\Omega }_D^1(M_n)\mathrm{\Omega }_D^1(𝒜)`$ is determined as $`\theta ^k(_l)=\delta _l^k1`$. If through $`\theta ^i=(\theta ^\alpha ,\theta ^a)`$, where $`\alpha =0,1,2,3;a\{1,2,\mathrm{},n^21\}`$, the basis of algebra $`\mathrm{\Omega }_D^1(𝒜)`$ is designated, then the basis in $`Der(𝒜)`$ will be as $`e_i=(e_\alpha ,e_a)`$ and $`e_\alpha =e_\alpha ^\mu _\mu `$ are 4-generators for algebra $`C(R_4)`$, and $`e_\alpha =ad(E_\alpha )`$ is a basis of $`Der(M_n)`$. In this case $`\theta ^i`$ accepts a kind

$$\theta ^\alpha =\theta _\mu ^\alpha dx^\mu ,\theta ^a=E_bE^adE^b$$

Differential $`d_h`$ and $`d_V`$ will accept a kind

$$d_Hf=e_\alpha f\theta ^\alpha ,d_Vf=e_af\theta ^a,f𝒜$$

Besides there is an canonical element $`\theta \mathrm{\Omega }_D^1(Mn)`$ determined as $`\theta =E_k\theta ^k`$. Within the framework of geometry Connes, noncommutative bundle is right, or left $`𝒜`$ \- module, or is more exact, Hermit $`𝒜`$ \- module. On this module following connection is defined

$$\omega =A+\chi $$

and $`A`$\- it is an element of $`\mathrm{\Omega }_H^1`$ and $`\chi `$\- an element of $`\mathrm{\Omega }_V^1`$. In turn, the element $`\chi `$ is divided for the sum

$$\chi =\theta +\varphi $$

because of the fact, that $`\mathrm{\Omega }_D^1(M_n)\mathrm{\Omega }_V^1`$. Here $`\varphi `$ is a field Higgs. If through $`𝒰`$ it displays the group of gauge maps , then an element $`\omega ,A,\theta ,\varphi `$ will be transformed under action $`g𝒰`$ by a rule

$$\omega ^{}=g^1\omega g+g^1dg$$
$`\left(5\right)`$
$$A^{}=g^1+g^1d_Hg,\theta ^{}=\theta ,\varphi ^{}=g^1\varphi g$$

The intensity of a field of deformation is determined through the 2-form $`\mathrm{\Omega }`$

$$\mathrm{\Omega }=d\omega +\omega ^2$$

Or

$$\mathrm{\Omega }=F+D_H\varphi +\mathrm{\Omega }_V$$
$`\left(6\right)`$
And

$$F=d_HA+[A,A]=\frac{1}{2}F_{\alpha \beta }\theta ^\alpha \theta ^\beta $$

$$\mathrm{\Omega }_V=\frac{1}{2}\mathrm{\Omega }_{ab}\theta ^a\theta ^b,\mathrm{\Omega }_{ab}=[\varphi _a,\varphi _b]C_{ab}^c\varphi _c$$

$$D_H\varphi =d_h\varphi +A\varphi +\varphi A=\theta ^\alpha \left(e_\alpha +[A_\alpha ,\varphi ]\right)$$

$$\varphi =\varphi _\alpha \theta ^\alpha ,A=A_\alpha \theta ^\alpha $$

Action of a field will be as

$$L=\frac{1}{2}(\mathrm{\Omega }_{ij},\mathrm{\Omega }_{ij})$$
$`\left(7\right)`$
Thus, in frameworks model , we shall receive probability of transition to plasticity

$$Z=\mathrm{exp}\left(kL\right)$$

Where $`k`$\- it is factor of elasticity.

3. We discuss the problem put above. For simplicity we shall be limited to consideration only pure deformation also we shall not take into account the interaction with a field Higgs. Actually, all above mentioned models were considered only in such similar situations. In these models was assumed about existence of limiting surfaces deformations (or loading). It is equivalent to the requirement of finiteness of action of deformation fields . Then, for search of a configuration of a field of deformation, it is necessary to solve the equation, which is a condition extreme of action $`L`$,

$$\delta L\left[A_\mu \right]=0$$
$`\left(9\right)`$
However, we have also other condition for action

$$d^4xTr\left[(F_{\mu \nu }+\pm \stackrel{~}{F}_{\mu \nu })^2\right]$$

where $`\stackrel{~}{F}_{\mu \nu }`$ is dual with $`F_{\mu \nu }`$. Because of that $`Tr[F_{\mu \nu },F_{\mu \nu }]=Tr[\stackrel{~}{F}_{\mu \nu },\stackrel{~}{F}_{\mu \nu }]`$, we receive

$$d^4xTr[F_{\mu \nu },F_{\mu \nu }]d^4xTr[\stackrel{~}{F}_{\mu \nu },\stackrel{~}{F}_{\mu \nu }]$$

Is like as in work Belavin etc. , we are interested in importance of self-dual and self-antidual configuration of a field of deformation. Thus, one of a extreme condition of action shall write down as

$$\stackrel{~}{F}_{\mu \nu }=F_{\mu \nu }$$
$`\left(10\right)`$
Belavin etc. used this condition for giving the instanton for gauge field. So, instead of the equation (9) we shall solve more simple the equation (10). For it we use method in  for search of multiinstanton in special case, when $`n,r=1`$. Let the field of deformation $`A=A_{mu}dx^{mu}\mathrm{\Omega }_D^1(R_4)1`$, has a individual kind

$$A_\mu \left(x\right)=iB_{\mu \nu }_\nu \left(\mathrm{ln}\varphi \left(x\right)\right)$$
$`\left(11\right)`$
where $`\varphi (x)`$ \- it is a scalar function, and the matrix $`B_{\mu \nu }`$ looks like

$$B_{\mu \nu }=\frac{1}{2}\left\{\begin{array}{cccc}0& +\sigma _3& \sigma _2& \sigma _1\\ \sigma _3& 0& +\sigma _1& \sigma _2\\ +\sigma _2& \sigma _1& 0& \sigma _3\\ +\sigma _1& +\sigma _2& +\sigma _3& 0\end{array}\right\}$$

where $`\sigma _1,\sigma _2,\sigma _3`$\- it is matrixes Pauli. The gauge group is admitted as $`SU(2)`$. Setting $`A_\mu `$ in expression (6) that condition of self-dual (10) will be as

$$_\sigma _\sigma \left(\mathrm{ln}\varphi \right)+\left(_\sigma \mathrm{ln}\varphi \right)^2=0$$

or

$$\frac{^2\varphi }{\varphi }=0,\text{where}^2=_\sigma _\sigma $$
$`\left(12\right)`$
Using a method in  we receive instanton, being the solution (12)

$$\varphi \left(x\right)=1+\frac{\lambda _{1}^{}{}_{}{}^{2}}{x^2}$$
$`\left(13\right)`$
and multiinstanton

$$\varphi \left(x\right)=1+\underset{i=1}{\overset{N}{}}\frac{\lambda _{i}^{}{}_{}{}^{2}}{x_\mu a_{i\mu }^2}$$
$`\left(14\right)`$
where $`a_{i\mu },\lambda _i`$\- it is any valid constant. Setting $`\varphi (x)`$ from (12) and (13) in (11) we receive deformation wave as matrixes:

One-instanton

$$F_\mu \left(x\right)=2i\lambda ^2B_{\mu \nu }\frac{y_\nu }{y^2\left(y^2+\lambda ^2\right)}$$

Multi-instanton

$$A_\mu \left(x\right)=2iB_{\mu \nu }\left(\underset{i}{\overset{N}{}}\frac{\lambda _I^2y_{i\nu }}{y_{i}^{}{}_{}{}^{4}}\right)/\left(1+\underset{j}{\overset{N}{}}\frac{\lambda _i^2}{y_j^2}\right)$$

where $`(y_i)_\mu (xa_i)_\mu ,i,j=1,2,\mathrm{},N`$

4. What real nature of sectors Koiter on a limiting surface?. With the account of solution of instanton , this problem can be easily clear. Really, using these instanton, we receive topological number Pontriagin

$$Q=\frac{1}{16\pi ^2}d^4xTr\left[F_{\mu \nu }F_{\mu \nu }\right]$$

This number tell me, that space of deformation during plastic deformation will be is divided on $`Q`$ of independent sectors. Thus, limiting surface, constructed in the same space, too is divided on $`Q`$ of parts. It is sectors Koiter, described in works \[1-3,10-11\]. We see, that occurrence of singularity Koiter has the topological reasons. Now we can speak, that process of plastic deformation is a process of birth and distribution of solitons and like-solitons. It is necessary to notice, that the defect, from the point of view of fluctuation, is one of kinds of solitons. The soliton wave is considered as the concentrated located energy, which is distributed in space without dissipation with constant velocity. Thus, the development of process of plastic deformation occurs with constant velocity and it needlessly increases a load. It will allow to understand, why a part of plasticity of the diagram stress - deformation has the horizontal form.