Äîêóìåíò âçÿò èç êýøà ïîèñêîâîé ìàøèíû. Àäðåñ îðèãèíàëüíîãî äîêóìåíòà : http://poly.phys.msu.ru/~bodrova/JStatMech.pdf Äàòà èçìåíåíèÿ: Mon Dec 14 16:28:32 2009 Äàòà èíäåêñèðîâàíèÿ: Sun Apr 10 22:43:57 2016 Êîäèðîâêà:

J

ournal of Statistical Mechanics: Theory and Experiment

An IOP and SISSA journal

A mo del of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Nikolay V Brilliantov1,2, Anna S Bo drova and P L Krapivsky3
1 2

2

University of Leicester, University Road, Leicester LE1 7RH, UK Moscow State University, Vorobievy Gory 1, 119899, Moscow, Russia 3 Department of Physics, Boston University, Boston, MA 02215, USA E-mail: nb144@leicester.ac.uk, b odrova@p olly.phys.msu.ru and paulk@bu.edu Received 7 Novemb er 2008 Accepted 11 February 2009 Published 16 June 2009 Online at stacks.iop.org/JSTAT/2009/P06011 doi:10.1088/1742-5468/2009/06/P06011

Abstract. A simple model of ballistic aggregation and fragmentation is
prop osed. The model is characterized by two energy thresholds, Eagg and Efrag , which demarcate different typ es of impacts: if the kinetic energy of the relative motion of a colliding pair is smaller than Eagg or larger than Efrag , particles resp ectively merge or break; otherwise they reb ound. We assume that particles are formed from monomers which cannot split any further and that in a collisioninduced fragmentation the larger particle splits into two fragments. We start from the Boltzmann equation for the mass­velocity distribution function and derive Smoluchowski-like equations for concentrations of particles of different mass. We analyze these equations analytically, solve them numerically and p erform Monte Carlo simulations. When aggregation and fragmentation energy thresholds do not dep end on the masses of the colliding particles, the model b ecomes analytically tractable. In this case we show the emergence of the two typ es of b ehavior: the regime of unlimited cluster growth arises when fragmentation is (relatively) weak and the relaxation towards a steady state occurs when fragmentation prevails. In a model with mass-dep endent Eagg and Efrag the evolution with a crossover from one of the regimes to another has b een detected.

Keywords: irreversible aggregation phenomena (theory), granular matter,
kinetic theory of gases and liquids, Boltzmann equation

c 2009 IOP Publishing Ltd and SISSA

1742-5468/09/P06011+18$30.00


A model of ballistic aggregation and fragmentation

Contents 1. Intro duction 2. The mo del 3. Theoretical analysis 3.1. Constant rates . . . . . . . . . . . . . . . . . . 3.1.1. Unlimited cluster growth, < 1. . . . . 3.1.2. Relaxation to a steady state, > 1. . . 3.2. Mass-independent energy thresholds . . . . . . 3.3. Dependence of the energy thresholds on masses 4. Numerical simulations 5. Conclusion References . . . . o . . . . f ..... ..... ..... ..... colliding ..... ..... ..... ..... particles . . . . . . . . . . . . . . . 2 3 6 .6 .7 .7 .9 . 12 14 15 17

J. Stat. Mech. (2009) P06011

1. Intro duction Collision-induced aggregation and fragmentation are ubiquitous pro cesses underlying numerous natural phenomena. For a gentle collision with a small relative velo city, colliding particles can merge; a violent collision with a large relative velo city can cause fragmentation. For intermediate relative velo cities, particles usually rebound. These collisions may still be irreversible--the kinetic energy could be lost in inelastic collisions. Important examples of such systems are dust agglomerates in the Earth atmosphere or in interstellar dust clouds and proto-planetary discs [1]­[4]. Another example is dynamic ephemeral bo dies in planetary rings; see e.g. [5]­[7]. A comprehensive description of the aggregation and fragmentation kinetics in such systems is very complicated. Therefore it is desirable to develop idealized mo dels that involve three kinds of collisions in the simplest possible way. The understanding of the ballistic-controlled reactions is still quite incomplete [8]. Ballistic aggregation has attracted most attention (see [3, 4], [9]­[14] and references therein) and a few studies were also devoted to ballistic fragmentation (see [15]­[18]). The situation with aggregation and fragmentation operating simultaneously has been analyzed only for a very special case in which all particles have the same relative velo city and the after-collision fragment mass distribution obeys a power law [7]. Moreover, studies of pure ballistic fragmentation are usually based on the assumption that all particles may split independently of their mass and relative velo city between the colliding grains [16, 17]. In reality, the type of an impact strongly depends on the relative velo city [4]; furthermore, the agglomerates are comprised of primary particles (`grains') that cannot split into smaller fragments [5, 6]. The fragmentation mo del with the splitting probability depending on energy has been studied in [19]; this mo del, however, do es not consider ballistic impacts of many particles, but rather an abstract pro cess of a successive fragmentation of one bo dy with a random distribution of the bulk energy between fragments.
doi:10.1088/1742-5468/2009/06/P06011 2


A model of ballistic aggregation and fragmentation

In this paper we propose a mo del of ballistic aggregation and fragmentation which accounts for three kinds of collisions depending on masses and relative velo city of a colliding pair. In section 2 we intro duce the mo del, write the Boltzmann equation for the joint mass­velo city distribution function, and deduce from the Boltzmann equation the rate equations for concentrations of various mass species. Section 3 is devoted to the theoretical analysis of rate equations. Numerical verification of theoretical results and simulation results in situations intractable theoretically is given in section 4. Section 5 concludes the paper. 2. The mo del Consider a system comprised of primary particles (monomers) of mass m1 and radius r1 , which aggregate to form clusters of 2, 3,... ,k ,... monomers with masses mk = km1 . In some applications (e.g. in mo deling of dynamic ephemeral bo dies [5, 6]) it is appropriate to consider clusters as ob jects with fractal dimension D ; for compact clusters D = 3. The characteristic radius of an agglomerate containing k monomers scales with mass as rk r1 k 1/D . We assume that when the kinetic energy of two colliding clusters in the center-of-mass reference frame (the `relative kinetic energy' in short) is less than Eagg , they merge. In this case a particle of mass (i + j )m1 is formed. If the relative kinetic energy is larger than Eagg , but smaller than Efrag , the colliding particles rebound without any change of their properties. Finally, if the relative kinetic energy exceeds Efrag , one of the particles (we assume that the larger one ) splits into two fragments. We denote as pi,k-i the probability that a particle of mass k splits into particles of masses i and k - i. Obviously, i pi,k-i = 1 and pi,k-i = 0 if k i. We restrict ourselves to dilute and spatially uniform systems. Let fi f (vi ,t) be the mass­velo city distribution function which gives the concentration of particles of mass mi with the velo city vi at time t. The mass­velo city distribution function evolves according to the Boltzmann equation agg fra re fk (vk ,t) = Ik + Ik b + Ik g , t (1)

J. Stat. Mech. (2009) P06011

agg fra re where Ik , Ik b and Ik g are respectively the collision integrals describing collisions leading to aggregation, rebound, and fragmentation. The first integral reads agg Ik (vk ) = 1 2 i+j =k 2 ij

dvi

dvj

de (-vij · e ) |vij · e|

â fi (vi ) fj (vj ) (Eagg - Eij ) (mk vk - mi vi - mj vj ) -
j 2 kj

dvj

de (-vkj · e) |vkj · e| fk (vk ) fj (vj ) (Eagg - Ekj ) . (2)

Here ij = r1 (i1/D + j 1/D ) is the sum of radii of the two clusters, while mk = mi + mj and mk vk = vi mi + mj vj , due to the conservation of mass and momentum. We also intro duce the relative velo city, vij = vi - vj , the reduced mass, ij = mi mj /(mi + mj ), and the 2 relative kinetic energy, Eij = 1 ij vij . The unit vector e specifies the direction of the 2 inter-center vector at the collision instant. The factors in the integrand in equation (2) 2 have their usual meaning (see e.g. [20]): ij |vij · e| defines the volume of the collision
doi:10.1088/1742-5468/2009/06/P06011 3


A model of ballistic aggregation and fragmentation

cylinder, (-vij · e) selects only approaching particles and (Eagg - Eij ) guarantees that the relative kinetic energy do es not exceed Eagg to cause the aggregation. The first sum in the right-hand side of equation (2) refers to collisions where a cluster of mass k is formed from smaller clusters of masses i and j , while the second sum describes the collisions of k -clusters with all other aggregates. For collisions leading to fragmentation we have
fra Ik g (vk ) = j ij

p

k,j -k

2 1 - 1 i,j ij 2

dvj

dvi

de (-vij · e)

â |vij · e| fj (vj ) fi (vi )(Eij - Efrag )(vi , vj , vk ) -
ik

J. Stat. Mech. (2009) P06011

1 - 1 i, 2

k

2 ki

dvi

de (-vki · e) |vki · e| (3)

â fk (vk ) fi (vi ) (Eki - Efrag ) ,

where mj = mk + mj -k and we use the abbreviation, (vi , vj , vk ) = (mj vj + mi vi - mk vk - mk-j vk-j + mi vi ) for the factor which guarantees the momentum conservation at the collision. The after-collisional velo cities vk-j and vi are determined by a particular fragmentation mo del. The first sum in equation (3) describes the collision of particles of mass i and j (j k , j i) with the relative kinetic energy above the fragmentation threshold Efrag . The larger particle, i.e. the particle of mass j , splits with the probability pk,j -k into two particles of mass k and j - k , thereby giving rise to a particle of mass k . The second sum describes the opposite pro cess, when particles of mass k break in collisions with smaller particles. In the present study we do not need an explicit expression for the velo cities vk-j and vi of the fragments. We also do not need an expression for the re collision integral Ik b ; it has the usual form (see e.g. [20]) with a slight mo dification to account for the requirement that the relative kinetic energy Eij belongs to the interval (Eagg < Eij < Efrag ). Thus we have a mixture of particles of different masses and each species generally has its own temperature. For this (granular) mixture we write ni = dvi fi (vi ), N=
i

ni ,

(4)

where ni is the number density (concentration) of particles of mass i and N is the total number density. Using the mean kinetic energy of different species one can also define the partial granular temperatures Ti for clusters of mass i and effective temperature T of the mixture [21]. We assume that the distribution function fi (vi ,t) may be written as [4, 14, 21] fi (vi ,t) = ni (t) i (c i ), 3 v0,i (t) ci vi , v0,i (5)

2 where v0,i (t) = 2Ti (t)/mi is the thermal velo city and (ci ) the reduced distribution function. For the force-free granular mixtures the velo city distribution functions of the

doi:10.1088/1742-5468/2009/06/P06011

4


A model of ballistic aggregation and fragmentation

components are not far from the Maxwellian distribution [21], which reads, in terms of the reduced velo city c = v/v0 , M (c ) =
-3/2

exp(-c2 ).

(6)

The equipartition between different components may, however, break down, in the sense that the partial temperatures Ti are not equal and differ from the effective temperature T [21]. Here we ignore the deviation from the Maxwellian distribution4 and possible violation of the equipartition and use the approximation i (ci ) M (ci ) and Ti T for all i. We also assume that the temperature of the system do es not depend on time. This is formally inconsistent within the Boltzmann equation (1), yet in many applications the temperature is approximately constant on an astronomical timescale, e.g. this happens in planetary rings where the viscous heating due to the shearing mo de of the particle orbital motion keeps the granular temperature constant [22, 4]. The consistent approach would be to mo dify the Boltzmann equation to take into account gradients of the lo cal hydro dynamic velo city, which will result in additional terms in the velo city distributions fi , proportional to these gradients. If we assume that such gradients are very small but still sufficient to support constant temperature due to viscous heating, we can neglect the small corrections to the distribution functions and approximate them with a gradient-free form (5). Integrating equation (1) over vk we obtain the equations for the zero-order moments of the velo city distribution functions fk , that is, for the concentrations nk . Taking into account that collisions resulting in rebounds do not change the concentrations of different species and using (2)­(3) together with (4)­(6) we arrive at rate equations d 1 nk = dt 2
j

J. Stat. Mech. (2009) P06011

Ci,j ni nj - nk
i+j =k i=1 k

Ck,ini +
j =k +1 i=1

1 Ai,j ni nj 1 - i, 2

j

p

k,j -k

- nk (1 - 1k )
i=1

Ak,ini 1 - 1 i, 2

k

,

(7)

with rates given by C
i,j 2 = 2ij

2T ij 2T ij

1/2

1 - 1+
1/2

Eagg T

e-

E

agg

/T

(8)

2 Ai,j = 2ij

e-

E

frag

/T

.

It is useful to verify that the above kinetic equation (7) fulfills the condition of mass conservation, k km1 nk = M = const, where M is the total mass density. The probability of splitting pik depends on geometric and mechanical properties of the aggregates and generally it is quite complicated. For concreteness we fo cus on splitting into (almost) equal fragments. Namely, we assume that a particle of mass 2km1 splits
Note that for a slightly modified model, where only in a small criterion, particles merge and only in a small fraction of collisions split, the velocity distribution is close to the Maxwellian, like modified model would lead to the same kinetic equation (9), but
4

fract that in a with

ion of collisions that fulfill the aggregation fulfill the fragmentation criterion, particles granular mixture. At the same time this the renormalized timescale.

doi:10.1088/1742-5468/2009/06/P06011

5


A model of ballistic aggregation and fragmentation

into two equal halves, while a particle of mass (2k +1)m1 splits into particles of mass km1 and (k +1)m1 . For this choice of the splitting probability, the kinetic equation reads d 1 nk = dt 2
k

Ci,j ni nj - nk
i+j =k 2k i=1

Ck,ini - nk
i=1

Ak,ini (1 - k,i /2)
2k +1

+2
i=1 2k -1

A2k,i n2k ni (1 - 2k,i /2) +
i=1

A2k /2) .

+1,i n2k +1 ni

(1 - 2k

+1,i

/2)

+
i=1

A2k

-1,i n2k -1 ni

(1 - 2k

J. Stat. Mech. (2009) P06011

-1,i

(9)

In the following sections we study this mo del theoretically and numerically. 3. Theoretical analysis To understand the qualitative behavior it is instructive to start with the simplest mo del which allows an analytical treatment.
3.1. Constant rates

Consider first the mo del with constant rates Ai,j and Ci,j . Without loss of generality we can cho ose these rates to be C
i,j

= 2,

Ai,j = 2.

(10)

The parameter quantifies the relative intensity of fragmentation with respect to aggregation. Fragmentation prevails when > 1 while aggregation wins in the opposite case of < 1. If = 1 two pro cesses are in a balance. Even in this simple case we still ought to analyze a cumbersome system of infinitely many equations. To gain insight it is useful to consider the evolution of the total density N= nk . (In many problems involving aggregation and fragmentation this quantity satisfies a simple equation that do es not contain other densities.) Summing up all of equation (9) we obtain N = -(1 - )N 2 - n2 1 (11)

which do es indeed have a neat form, although it additionally involves the density of monomers. This density evolves according to n1 = -2n1 N +2 n2 (2n1 + n2 )+ n3 (n1 + n2 )+ 1 n2 . 23 Although we do not have a closed system conclusions. Equation (11) indicates that two aggregation prevails, the system continues to larger clusters; when > 1, one expects the analyze these situations in more detail.
doi:10.1088/1742-5468/2009/06/P06011

(12)

we can already reach some qualitative regimes are possible. If < 1, i.e. when evolve leading to formation of larger and system to reach a steady state. We now
6


A model of ballistic aggregation and fragmentation

3.1.1. Unlimited cluster growth, < 1. In this case larger and larger Since the total mass is conserved, one expects the concentration of small N when t 1 and therefore one can omit t decrease. Therefore, n1 the right-hand side of (11). Similarly one can keep only the first term side of (12). This leads to the simplified equations n1 -2n1 N N -(1 - )N 2 ,

clusters will arise. clusters to rapidly he second term on on the right-hand (13)

which are solved to yield the large time behavior: N 1 , (1 - )t
2/(1-)

(14) (15)

J. Stat. Mech. (2009) P06011

n1 t-

.

Further, one anticipates that the density distribution approaches the scaling form k x= z (16) nk t-2z (x), t in the scaling limit t , k , with the scaled mass x = (Here z is the dynamic exponent characterizing the average mass: scaling form agrees with mass conservation: knk dxx(x) is independent. The exponent z can be found from the known asymptotic behavior writing N=
k 1

k/tz kept finite. k tz .) The manifestly time of N (t). Indeed,

nk

t-

z 0



dx (x) t-

z

and matching this with already known asymptotic behavior (14) we conclude that z = 1. If 1 and combine this asymptotic with z = 1 and we further assume that (x) x for x -2 - the scaling ansatz (16) we obtain n1 t t . Matching with (15) we get = 2/(1 - ). Therefore 1 nk 2 t when k nk t-
2/(1-)

k t

2/(1-)

(17)

t. Obviously, the above equation implies the asymptotic time dependence and the mass dependence nk k 2/(1-) for x = k/t 1.

3.1.2. Relaxation to a steady state, > 1.

For > 1 the system evolves to a steady state with constant concentration of clusters. In this case nk = N = 0 and equation (11) yields n1 = N 1 - -1 . (18) 1 and the governing

The densities nk rapidly decay with k . Therefore n2k nk for k equation (9) for the stationary concentrations simplifies to
k -1

ni nk
i=1

-i

- 2(1 + )nk N = 0,

(19)
7

doi:10.1088/1742-5468/2009/06/P06011


A model of ballistic aggregation and fragmentation

where we have approximated a finite sum up to k 1 by an infinite sum and ignore the terms containing n2k ,n2k±1 . The above equation is supposed to be valid for large k ; it is certainly invalid for k = 1 when the right-hand side do es not vanish. The qualitative form of the large k asymptotic behavior is determined by the mathematical structure of (19). To extract this asymptotic let us consider the simplest version when equation (19) is valid for all k 2. Specifically, let us probe the mo del
k -1

ni nk
i=1

-i

- 2(1 + )nk N = -(1 + 2)Nk,1 ,

(20)

J. Stat. Mech. (2009) P06011

where the amplitude (1 + 2) was chosen to set N = 1. (For mo del (20), this choice merely sets the overall amplitude.) The infinite system (20) forms a recurrence and therefore it is solvable. Intro ducing the generating function N (z ) =
k 1

nk z

k

(21)

we recast (20) into a quadratic equation N 2 - 2(1 + )N +(1+2)z = 0 which is solved to yield N = (1 + ) - Expanding N (z ) we arrive at 1 1+ 1- nk = (1 + ) 4
k 2

(22)

(1 + )2 - (1 + 2)z. (k - 1/2) . (k +1)

(23)

(24)

From this solution one gets n1 = ( +1/2)/( + 1), which of course directly follows from equation (20) as well. For large k , equation (24) simplifies to nk 1+ k 4
-3/2

1 1- (1 + )

k 2

.

(25)

We consider other tractable versions where equation (19) is exact above a certain threshold, k +1, while for k = 1,... , we use the same mo dification as in equation (20) for k = 1. In this case instead of (22) one gets N 2 - 2(1 + )N + P (z ) = 0 with P (z ) = A1 z + ··· + A z . The ro ot of P (z ) = (1 + )2 closest to the origin is positive (one can show that it exceeds unity, z > 1) and non-degenerate. Expanding the generating - function N = 1 + - (1 + )2 - P (z ) leads to the asymptotic nk k -3/2 z k . The above argument favors the asymptotic behavior nk Ak
-3/2 -k

e

.

(26)

This asymptotic form is universal and only the parameters A, depend on the specificity of the mo del, that is on the parameter . It is impossible to determine A, since mo dels like (20) are uncontrolled approximations. Let us still use such mo dels and cho ose the simplest one which obeys the exact relation of equation (18). The mo del (20) is inappropriate as it fails to satisfy (18):
doi:10.1088/1742-5468/2009/06/P06011 8


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 1. Evolution of cluster concentrations nk in the case of constant rates with = 0.5. The lines corresp ond to the numerical solution of 1000 differential equations; symb ols are the results of MC simulation with 100 000 monomers for the monodisp erse initial conditions. Note that while n1 (t) always decays, the nk (t) initially increase and then decay to zero.

( + 1/2)/( + 1) > N = 1)



1 - -1 . Mo difying equation (19) at k = 1, 2 yields (we still set (27)

ni nj - 2(1 + )nk = -qk,1 - (2 +1 - q )k,2 ,
i+j =k

with

q = 2 (1 + ) 1 -

-1

ensuring the validity of (18). The same approach as before gives (26) with A and . In particular, 1 - -1 + - 1 - -1 , (28) = ln (1 + ) where we have used the shorthand notation = 1 + 2 - 2(1 + ) 1 - -1 . For = 2 (which we have studied numerically) one gets = 0.495 156 ... . (29)

This is an uncontrolled approximation, of course. Interestingly, the result is rather close to the numerically obtained value 0.465.
3.2. Mass-indep endent energy thresholds

We now turn to the analysis of the situation when aggregation and fragmentation energy thresholds Eagg and Efrag are constant. In this case the total density of clusters evolves
doi:10.1088/1742-5468/2009/06/P06011 9


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 2. Evolution of the total numb er of clusters N for the same system as in figure 1. The solid line corresp onds to the numerical solution of the differential equations, symb ols to the MC simulation and the dotted line shows the theoretical 1. prediction, equation (14), N (t) t-1 for t

Figure 3. The long time limit b ehavior of same system as in figure 1. In accordance equations (15) and (17), the cluster concentr the same slop e t-2/(1-) t-4 , shown by the

cluster concentrations nk for the with the theoretical predictions, 1 with ations nk (t) decay for t solid line.

according to


N = - 1 (1 - ) 2
i=1 j =1

Ci,j ni nj - 1 C1,1 n2 , 1 2

(30)

where -1 = eEfrag /T (1 - (1 + Eagg /T )e-Eagg /T ). This again implies the existence of the two opposite evolution regimes: for > 1 the relaxation to a steady state is expected, while for < 1 we have the regime of the unlimited cluster growth.
doi:10.1088/1742-5468/2009/06/P06011 10


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 4. The cluster mass distribution at different time instants for the same system as in figure 1. The initial cluster distribution at t = 10 (long-dashed line) drastically differs from that in the scaling regime, t 1. The dotted­dashed and dashed lines show resp ectively the cluster mass distribution for t = 500 and t = 1000. The solid line shows the theoretical prediction, equation (17), nk k , with = 2/(1 - ) = 2.

The rates Cij = C (i, j ) and Aij = A(i, j ) differ by a constant factor ; moreover, they are homogeneous functions of their arguments: A(ai, aj ) = a A(i, j ), C (ai, aj ) = a C (i, j ), (31) with the exponent 1 2 -, (32) D2 which follows from the relation for a particle's mass mi = im1 , the cross-section of the 2 collision cylinder, ij (i1/D + i1/D )2 , and equations (8). Plugging the scaling ansatz (16) into equation (9), taking into account that for k 1 the summation may be approximated by integration, and exploiting the homogeneity of the rate kernels, equation (31), we obtain (see e.g. [23, 25, 24] for analysis of similar integro-differential equations for the re-scaled mass distribution) = dy (y )[(Cx,y + Ax,y )(x) t(3- )z 0 - 1 Cy,x-y (x - y ) - 4A2x,y (2x)]. 2 From equation (33) we find the scaling exponent -1 3 2 1 = - . z= 1- 2D t2z
+1

z

(2(x)+ x (x)) =

1



(33)

(34)

Although the scaling theory do es not allow one to determine the scaling function (x), one can find the total concentration of clusters from equation (34): N (t) t-z t-2D/(3D-4) . (35) Correspondingly, the average clusters mass grows as k = M/N tz .
doi:10.1088/1742-5468/2009/06/P06011 11


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 5. Evolution of cluster con coefficients with = 2. After a steady state. Lines corresp ond equations, symb ols to the results monodisp erse initial conditions.

centrations nk (t) for the case of constant kinetic certain p eriod of time the system relaxes to a to the numerical solution of 1000 differential of MC simulation (100 000 monomers) for the

3.3. Dep endence of the energy thresholds on masses of colliding particles

In the preceding analysis we have assumed that Eagg and Efrag do not depend on the mass of colliding particles so that ij = Aij /Cij = is constant. In reality, however, such dependence do es exist, implying that ij is a function of i and j . Still, if ij > 1 or ij < 1 for all i and j , the qualitative behavior of a system is similar to that for the case of constant Eagg and Efrag : for ij > 1 a relaxation to a steady state is expected, while for ij < 1 an unlimited cluster growth is observed. The most interesting behavior is expected when ij - 1 changes its sign with increasing cluster masses i and j . In this case one anticipates a crossover from one type of evolution to another. To cho ose a realistic dependence of Eagg and Efrag on the masses of colliding particles, one needs more details of the collision pro cess. We shall use the threshold energy for ballistic aggregation that takes into account surface adhesion [26]. In this case
0 Eagg (i, j ) = Eagg

ij i+j

4/3

,

(36)

0 where Eagg is expressed in terms of the monomer radius, particle surface tension, the Young mo dulus and the Poisson ratio of the particle material (see [26] for the explicit 0 expression for Eagg ). For the energy of fragmentation we assume that it is equal to the energy required to create an additional surface, which may be roughly estimated as twice the area of the equatorial cross-section of the larger particle (recall that the mo del assumes, that the larger particle in a collision pair breaks down). Hence we adopt the following mass

doi:10.1088/1742-5468/2009/06/P06011

12


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 6. The steady state distribution of cluster mass for the same system as in figure 5. The solid line corresp onds to the numerical solution, the dotted line to the theoretical prediction, equation (26), nk = Ak-3/2 e-k . The constant = 0.465, obtained by fitting, is very close to the theoretical value of = 0.495, equation (29).

Figure 7. Evolution of cluster concentrations nk (t) for the case of ballistic kinetic coefficients with constant aggregation and fragmentation energies Eagg /T = 0.9, Efrag /T = 3 and < 1. The cluster dimension is D = 3.

dependence for Efrag :
0 Efrag (ij ) = Efrag (ij i + ji j )2 ,

(37)

0 2 where ij = 1 if i > j , ij = 0 if i < j and ii = 1/2; Efrag = 2s r1 , with s being the surface tension.

doi:10.1088/1742-5468/2009/06/P06011

13


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 8. Evolution of cluster concentrations nk (t) in the case of ballistic kinetic coefficients with constant aggregation and fragmentation energies Eagg /T = 0.3, Efrag /T = 3 and > 1. The cluster dimension is D = 3. Like for the case of constant kinetic coefficients, the system relaxes to a steady state. Lines corresp ond to the numerical solution of 1000 equations, symb ols to the results of MC simulation (100 000 monomers) for the monodisp erse initial conditions.

4. Numerical simulations In our numerical studies we apply two different approaches--the solutions of the system of differential equations and the direct mo deling of random aggregation and fragmentation pro cesses (with the corresponding rates Ci,j and Ai,j ) by means of a Monte Carlo (MC) metho d. In the former case we use 1000 equations and in the later one 100 000 monomers (we always used the mono disperse initial conditions). The approach based on the solution of differential equation has an obvious deficiency as one must approximate an infinite system of equations with a finite one. The MC approach is more time-consuming, but it has an advantage of directly imitating the physical pro cesses in which particles are involved. To mo del the fragmentation and aggregation kinetics with MC we use the standard Gillespie algorithm [27, 28] (see [29] for the application of this algorithm to the aggregation and fragmentation pro cesses). The results presented in figures 1­11 confirm our theoretical predictions qualitatively and quantitatively. For constant rates two opposite types of evolution have indeed been observed: the relaxation to a steady state for dominating fragmentation ( > 1) and the unlimited cluster growth when aggregation prevails ( < 1). The two numerical approaches (the solution of the differential equations and MC) yield very close results. In figure 1 the evolution of the concentration of clusters of different mass is shown for the < 1 regime when the cluster growth continues ad infinitum. Note that all concentrations nk (t), except for n1 (t) which always decays, initially increase and then decay to zero. Figure 2 shows that the decay of cluster density N (t) agrees well with the theoretical prediction (14). Figures 3, 4 show respectively the asymptotic evolution
doi:10.1088/1742-5468/2009/06/P06011 14


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 9. Steady state cluster mass distribution of nk for the same system as in figure 8. Like for the case of constant kinetic coefficients with > 1, the distribution nk has a nearly exp onential form, equation (26), nk = Ak-3/2 e-k . The constant = 0.315 is obtained by fitting.

of cluster concentrations and the distribution of the cluster mass for x = k/t 1. The theoretical predictions, equations (15) and (17), are in a go o d agreement with the simulations. Relaxation to a steady state in the case when fragmentation dominates ( > 1) is illustrated in figure 5, while figure 6 demonstrates the corresponding stationary cluster mass distribution. Note that the numerical simulations confirm the theoretical form of the steady state cluster mass distribution. Qualitatively similar behavior is observed for the case of the ballistic kinetic co efficients, equations (8), with the constant aggregation and fragmentation energies. Again, for < 1, as for constant kinetic co efficients, clusters grow unlimitedly, figure 7, while for > 1 the system relaxes to a steady state, figure 8. The cluster mass distribution in the steady state may be anew very well fitted with the nearly exponential form, equation (26); see figure 9. In figure 10 the prediction (35) of the scaling theory is compared with the numerical data. Again we see that the agreement between the theory and simulations is rather satisfactory. Finally figure 11 illustrates evolution of the system with the ballistic co efficients that depend on the cluster mass in accordance with equations (36)and (37). It is interesting to note that the system tends initially to a quasi-steady state, as previously for the case of > 1, but then a crossover to a different evolution regime corresponding to < 1 takes place. In the latter regime all cluster concentrations decay with a similar slope, close to t-1 , still to be explained theoretically. 5. Conclusion We analyzed the dynamics of a system where particles move ballistically and undergo collisions which can lead to decrease or increase of the number of particles. The precise outcome depends on the kinetic energy Ekin in the center-of-mass reference frame. We
doi:10.1088/1742-5468/2009/06/P06011 15


A model of ballistic aggregation and fragmentation

J. Stat. Mech. (2009) P06011

Figure 10. Evolution of total numb er of clusters N (t) for the case of ballistic kinetic coefficients with constant aggregation and fragmentation energies Eagg /T = 0.9, Efrag /T = 3 and < 1 for different cluster dimensions D. Solid lines corresp ond to the numerical solution of 1000 differential equations; dotted lines show the prediction of the scaling theory, equation (35).

Figure 11. Evolution of the cluster concentrations nk (t) for the case of ballistic coefficients Cij , Aij with the mass-dep endent aggregation and fragmentation 0 energies Eagg (i, j ) and Efrag (i, j ), given by equations (36) and (37) with Eagg /T = 0 0.1, Efrag /T = 0.6. Lines corresp ond to the numerical solution of 1000 differential equations, symb ols to the results of MC simulation (100 000 monomers) for the monodisp erse initial conditions. The cluster dimension is D = 3. Note that the system tends initially to a steady state, clearly seen for the first few cluster masses. Later its evolution alters to the regime corresp onding to the unlimited cluster growth. In this regime the cluster concentrations decay with a slop e close to t-1 , shown in the figure by the dotted line. doi:10.1088/1742-5468/2009/06/P06011 16


A model of ballistic aggregation and fragmentation

proposed a simple mo del with two threshold energies, Eagg and Efrag , which define the type of an impact: for Ekin < Eagg the colliding particles merge, for Eagg < Ekin < Efrag they rebound, and for Efrag < Ekin one of the particles (the larger one) splits upon the collision. We assume that the aggregates are composed of 1, 2,... ,k ,... monomers and split into two equal (for an even number of monomers in the cluster) or almost equal (for an o dd number of monomers) pieces. The monomers are assumed to be stable, that is, they do not split further. For this mo del we wrote the Boltzmann kinetic equation for the mass­velo city distribution function of the aggregates and derived rate equations for the time evolution of the cluster concentrations. The bal listic rates were obtained in terms of the aggregation and fragmentation energy thresholds Eagg and Efrag , masses of the colliding particles and the temperature of the system (which was assumed to be constant). The Maxwellian velo city distribution for all species in the system was also assumed. We analyzed theoretically and studied numerically the rate equations. In the numerical studies we used two different metho ds--the solution of the system of differential equations and Monte Carlo mo deling. The two numerical metho ds yielded very close results. We started with the simplest case of constant rates and observed two opposite evolution regimes--the regime of unlimited cluster growth and of the relaxation to a steady state; we described both these cases analytically. For the regime of the unlimited cluster growth we obtained the asymptotic time dependence for the cluster concentrations and for their mass distribution. For the relaxation regime, which corresponds to the prevailing fragmentation, we derived the asymptotic behavior of the stationary mass distribution. In the evolving regime, the cluster concentrations decay as a power law in time; the stationary mass distribution has a nearly exponential form. Theoretical predictions are in a go o d agreement with numerical results. We also studied the case of mass-dependent rates arising in the situation when aggregation and fragmentation energy thresholds are constant. We observed that the behavior of the system is qualitatively similar to that of the system with the constant rates. Surprisingly, we detected that the steady state cluster mass distribution has also a nearly exponential form. We developed a scaling theory for the asymptotic large time behavior of the cluster concentrations and checked it numerically for different fractal dimensions of the aggregates. The numerical data agree well with the results of our theory. Finally, we explored numerically the case of the ballistic kinetic co efficients with the aggregation and fragmentation energies depending on the mass of colliding particles. For the aggregation energy threshold we use the result available in the literature for a collision of particles with surface adhesion. For the fragmentation energy threshold we adopted a mo del where Efrag is proportional to the surface energy of the maximal cross-section of the larger particle in the colliding pair. For this mo del the dependence on mass of Efrag is much stronger than that of Eagg . As the result, the evolution of the system, where the fragmentation initially prevails and drives it to a steady state, alters at later time when the unlimited cluster growth eventually wins and then it continues ad infinitum. References
[1] Chokshi A, Tielens A G G and Hollenbach D, 1993 Astrophys. J. 407 806 [2] Dominik C and Tielens A G G, 1997 Astrophys. J. 480 647

J. Stat. Mech. (2009) P06011

doi:10.1088/1742-5468/2009/06/P06011

17


A model of ballistic aggregation and fragmentation [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] Ossenkopf V, 1993 Astron. Astrophys. 280 617 Spahn F, Alb ers N, Sremcevic M and Thornton C, 2004 Europhys. Lett. 67 545 Greenb erg R, Davis D R, Weidenschilling S J and Chapman C R, 1983 Bull Am. Astron. Soc. 15 812 Weidenschilling S J, Chapman C R, Davis D R and Greenb erg R, 1984 Planetary Rings ed R Greenb erg and A Brahic (Tucson, AZ: Arizona University Press) pp 367­416 Longaretti P-Y, 1989 Icarus 81 51 Ben-Naim E, Krapivsky P, Leyvraz F and Redner S, 1994 J. Phys. Chem. 98 7284 Carnevale G F, Pomeau Y and Young W R, 1998 Phys. Rev. Lett. 64 2913 Trizac E and Hansen J-P, 1995 Phys. Rev. Lett. 74 4114 Fracheb ourg L, 1999 Phys. Rev. Lett. 82 1502 Fracheb ourg L, Martin Ph A and Piasecki J, 2000 Physica 279 69 Trizac E and Krapivsky P L, 2003 Phys. Rev. Lett. 91 218302 Brilliantov N V and Spahn F, 2006 Math. Comput. Simul. 72 93 Cheng Z and Redner S, 1990 J. Phys. A: Math. Gen. 23 1233 Krapivsky P and Ben-Naim E, 2003 Phys. Rev. E 68 021102 Pagonabarraga I and Trizac E, 2003 Granular Gas Dynamics ed T Poeschel and N Brilliantov (Berlin: Springer) Hidalgo R C and Pagonabarraga I, 2008 Phys. Rev. E 77 061305 Marsili M and Zhang Yi-C, 1996 Phys. Rev. Lett. 77 3577 Brilliantov N V and P¨ hel T, 2004 Kinetic Theory of Granular Gases (Oxford: Oxford University Press) osc Garzo V and Dufty J W, 1999 Phys. Rev. E 59 5895 Greenb erg R and Brahic A (ed), 1984 Planetary Rings (Tucson, AZ: Arizona University Press) van Dongen P G J and Ernst M H, 1985 Phys. Rev. Lett. 54 1396 Leyvraz F, 2003 Phys. Rep. 383 95 Cheng Z and Redner S, 1988 Phys. Rev. Lett. 60 2450 Brilliantov N V, Alb ers N, Spahn F and Poeschel T, 2007 Phys. Rev. E 76 051302 Gillespie D T, 1976 J. Comput. Phys. 22 403 Feistel R, 1977 Wiss. Z. Univ. Rostock 26 663 Poeschel T, Brilliantov N and Frommel C, 2003 Biophys. J. 85 3460

J. Stat. Mech. (2009) P06011

doi:10.1088/1742-5468/2009/06/P06011

18