Brownian motion under annihilation dynamics
aa r X i v : . [ c ond - m a t . s t a t - m ec h ] J a n Brownian motion under annihilation dynamics
Mar´ıa Isabel Garc´ıa de Soria, Pablo Maynar,
2, 3 and Emmanuel Trizac Universit´e Paris-Sud, LPTMS, UMR 8626, Orsay Cedex, F-91405 and CNRS, Orsay, F-91405 Laboratoire de Physique Th´eorique (CNRS UMR 8627),Bˆatiment 210, Universit´e Paris-Sud, 91405 Orsay cedex, France F´ısica Te´orica, Universidad de Sevilla, Apartado de Correos 1065, E-41080, Sevilla, Spain (Dated: November 30, 2018)The behavior of a heavy tagged intruder immersed in a bath of particles evolving under ballisticannihilation dynamics is investigated. The Fokker-Planck equation for this system is derived andthe peculiarities of the corresponding diffusive behavior are worked out. In the long time limit, theintruder velocity distribution function approaches a Gaussian form, but with a different temperaturefrom its bath counterpart. As a consequence of the continuous decay of particles in the bath, themean squared displacement increases exponentially in the collision per particle time scale. Analyticalresults are finally successfully tested against Monte Carlo numerical simulations.
PACS numbers: 51.10.+y,05.20.Dd,82.20.Nk
I. INTRODUCTION
In recent years, there has been some interest for systems where particles annihilate ballistically [1, 2, 3, 4, 5, 6, 7, 8].In these studies, the model considered consists of an ensemble of hard particles which evolve freely until a binaryencounter, which leads either to the annihilation of the colliding partners with probability p , or to an elastic collisionwith probability 1 − p . For this probabilistic annihilation model, there are no collisional invariants, and numericalsimulations have shown that for a broad class of initial conditions, the system reaches an homogeneous state in whichall the time dependence of the one particle distribution function is encoded in the density and temperature (defined asthe second velocity moment of the distribution function) [4, 5]. This is the so-called “Homogeneous Decay State”. Sucha behavior resembles the one of granular fluids (see [9] and references therein) where, if the system is stable, it evolvesinto an homogeneous cooling state, in which all the time dependence is borne by the granular temperature (in this casethe density is conserved) [10]. For the annihilation model, the hydrodynamic equations have been derived using theChapmann-Enskog method [11] by the usual assumption of the existence of a “normal solution”, whose space and timedependence occurs only through the hydrodynamic fields [6]. Recently, the hydrodynamic equations linearized aroundthe homogeneous decay state have been derived relaxing such an assumption [12]. Nevertheless, it must be assumedthat there is scale separation, i.e that the spectrum of the linearized Boltzmann collision operator is such that theeigenvalues associated to the hydrodynamic excitations are separated from the faster “kinetic eigenvalues”. Althoughthis property is valid for elastic collisions [13], it has not been proven for the probabilistic ballistic annihilation modelin general, but only for Maxwell molecules [14] and for p smaller than a given threshold [12].The objective in this paper is to study the simplest transport process in this system in which we can rigorouslyprove that there is scale separation. We will consider a tagged particle in a fluid in the homogeneous decay state,but collisions between the tagged particle and the particles of the fluid will be always elastic. The equation for thetagged particle is the Boltzmann-Lorentz equation [13, 15] which depends on the one particle distribution functionof the bath. In the limit of asymptotically large relative mass for the tagged particle, this equation reduces to aFokker-Planck equation which depends on the time-dependent density and temperature of the bath. Due to thestructure of this equation, we can prove that there is scale separation and that, in the long-time limit, the velocitydistribution function of the tagged particle approaches a Gaussian distribution but with a temperature that differsfrom that of the bath. A similar breakdown of equipartition has been reported for a heavy particle in a granular bath[17, 18, 19, 20, 21], a problem that can be mapped onto an elastic situation [22], at variance with the situation underscrutiny here. We also study the diffusion of the heavy particle and identify the diffusion coefficient as a Green-Kuboformula in terms of the velocity autocorrelation function. Finally, we perform Monte Carlo numerical simulations inorder to test our theoretical results. II. FOKKER-PLANCK EQUATION
We consider a tagged particle of mass m and diameter σ immersed in a low-density gas. This gas is composed ofhard spheres or disks of mass m g and diameter σ g which move ballistically until one particle meets another one; suchbinary encounters lead to the annihilation of the colliding partners with probability p or to an elastic collision withprobability 1 − p [1, 2, 3, 4, 5, 8]. Collisions between the particles of the gas and the tagged particle are always elastic. A. From Boltzmann-Lorentz to Fokker-Planck
The evolution equation for the probability density F ( r , v , t ) of the tagged particle is the Boltzmann-Lorentz equation[13, 15] (cid:18) ∂∂t + v · ∇ (cid:19) F ( r , v , t ) = J [ r , v , t | F, f ] , (1)where the collision operator is given by J [ r , v , t | F, f ] = σ d − Z d v Z d ˆ σ Θ( g · ˆ σ )( g · ˆ σ ) { F ( r , v ∗ , t ) f ( r , v ∗ , t ) − F ( r , v , t ) f ( r , v , t ) } . (2)Here f ( r , v , t ) is the distribution function of the particles in the gas, d is the space dimension, g = v − v is therelative velocity, Θ is the Heaviside step function, ˆ σ is a unit vector pointing from the centre of the gas particle tothe centre of the tagged particle at contact, and σ = σ + σ g . The precollisional velocities v ∗ and v ∗ are given by v ∗ = v − g · ˆ σ ) ˆ σ , (3) v ∗ = v + 21 + ∆ ( g · ˆ σ ) ˆ σ , (4)with ∆ = m g /m the (gas/tagged particle) mass ratio.We shall consider that the gas is in the homogeneous decay state, so its distribution function has the scaling form[5] f H ( v , t ) = n g ( t ) v dg ( t ) χ H ( c ) , c = v v g ( t ) , (5)where n g ( t ) is the number density of the gas, v g ( t ) = (cid:16) T g ( t ) m g (cid:17) / is the thermal velocity of the particles in the gasand χ H is an isotropic function depending only on the modulus c = | c | of the rescaled velocity. It can be seen thatthe homogeneous density and temperature obey the following equations [6] ∂n g ( t ) ∂t = − pν g ( t ) ζ n n g ( t ) , (6) ∂T g ( t ) ∂t = − pν g ( t ) ζ T T g ( t ) , (7)where we have introduced the collision frequency of the corresponding hard sphere fluid in equilibrium (with sametemperature and density) ν g ( t ) = n g ( t ) T / g ( t ) σ d − g m / π d − ( d + 2)Γ( d/ . (8)Here the dimensionless decay rates ζ n and ζ T are functionals of the distribution function and are approximatelyknown in the first Sonine approximation [4, 7], see Appendix A. Equations (6) and (7) can be integrated to obtainthe following power laws for the decay of the density and temperature n g ( t ) = n g (0) [1 + ν g (0) p ( ζ n + ζ T / t ] − ζn ζn + ζT , (9) T g ( t ) = T g (0) [1 + ν g (0) p ( ζ n + ζ T / t ] − ζT ζn + ζT . (10)As a consequence, we get n g T / g ∝ t − , a simplified form of a scaling relation common to all ballistically controlledprocesses [16].We next study the evolution equation for the tagged particle in the limit of large relative mass for the taggedparticle. In the limit ∆ ≪
1, it is possible to expand the collision operator J [ r , v , t | F, f ] in powers of ∆. In AppendixB it is shown that the leading order is J [ r , v , t | F, f ] ≃ ∂∂ v · [ A ( v ) F ( r , v , t )] + 12 ∂∂ v ∂∂ v : [ N ( v ) F ( r , v , t )] , (11)where A ( v , t ) = γ ( t ) v , N ( v , t ) = 2¯ γ ( t ) I, (12)and I is the second order unit tensor. The definitions of γ and ¯ γ are respectively γ ( t ) = γ e [ n g ( t ) , T g ( t )] a ( p ) , (13)¯ γ ( t ) = γ e [ n g ( t ) , T g ( t )] a ( p ) b ( p ) T g ( t ) m , (14)The friction coefficient γ e ( t ) is the same as for elastic bodies, and appears here as a function of the time-dependentdensity n g ( t ) and temperature T g ( t ) γ e [ n g ( t ) , T g ( t )] = 4 π d − d Γ( d/
2) ∆ / n g ( t ) (cid:18) T g ( t ) m (cid:19) / σ d − , (15)with a ( p ) and b ( p ) functionals of the distribution function of the bath which depend only on the parameter pa ( p ) = Γ( d/ d + 1) / Z d c χ H ( c ) c , (16) b ( p ) = 2 d + 1 R d c χ H ( c ) c R d c χ H ( c ) c . (17)In Appendix B, it is shown that the two terms on the right hand side of equation (11) are both of order n g v g σ d − ∆,while the other contributions in the Kramers-Moyal expansion are at least of order n g v g σ d − ∆ / . In the sameAppendix, the expressions for a ( p ) and b ( p ) are evaluated to first order in a Sonine expansion.Taking into account the approximate expression for the collision operator, Eq. (11), it is possible to write theBoltzmann-Lorentz equation as a Fokker-Planck equation for asymptotically small ∆ (cid:20) ∂∂t + v · ∇ (cid:21) F ( r , v , t ) = γ e ( t ) a ( p ) ∂∂ v · (cid:20) v + b ( p ) T g ( t ) m ∂∂ v (cid:21) F ( r , v , t ) . (18)As in the inelastic case, the Einstein relation is violated due to the fact that the distribution function of the bathis not Maxwellian [18, 23, 24, 25], which in turn implies that b ( p ) = 1. On the other hand, if we suppose that thevelocity of the tagged particle obeys a Markov process and write the corresponding Fokker-Planck equation, in termsof the jump moments, lim ∆ t → h ∆ v i / ∆ t and lim ∆ t → h ∆ v i / ∆ t , we obtain exactly Eq. (18). Here h . . . i means averageover different noise (bath) realizations. B. Coarse grained fields and relevant scales
We now focus on the study of the hydrodynamic fields of the tagged particle with the aid of the Fokker-Planckequation, Eq. (18). We define the mean velocity and the temperature of the Brownian particle as u ( t ) = Z d r Z d vv F ( r , v , t ) , (19) d T ( t ) = Z d r Z d v m ( v − u ) F ( r , v , t ) . (20)Taking moments in the Fokker-Planck equation, we obtain (see Appendix C) ∂ u ( t ) ∂t = − γ e ( t ) a ( p ) u ( t ) , (21) ∂T ( t ) ∂t = − γ e ( t ) a ( p )[ T ( t ) − b ( p ) T g ( t )] . (22)As the function γ e is a known functional of the gas density n g and temperature T g , equations (21) and (22) can beintegrated, which yields u ( t ) = u (0) [1 + ν g (0) p ( ζ n + ζ T / t ] − a ( p ) ζT (2 ζn + ζT ) ǫ , (23)and T ( t ) = b ( p ) T g (0)1 − ǫ [1 + ν g (0) p ( ζ n + ζ T / t ] − ζT ζn + ζT + (cid:20) T (0) − b ( p ) T g (0)1 − ǫ (cid:21) [1 + ν g (0) p ( ζ n + ζ T / t ] − ζTǫ (2 ζn + ζT ) . (24)In the above equation, we have introduced the dimensionless coefficient ǫǫ = pζ T ν g ( t )2 a ( p ) γ e ( t ) = √ dζ T d + 2) a ( p ) (cid:18) σ g σ (cid:19) d − p ∆ , (25)that is not necessarily a small quantity.As can be seen in Eq. (24), the behavior of the temperature depends strongly on the value of ǫ . If ǫ < t →∞ T ( t ) T g ( t ) = b ( p )1 − ǫ , ǫ < . (26)On the other hand, if ǫ > ǫ is essentiallythe quotient between the cooling rate of the gas and the relaxation rate of the tagged particle’s temperature. If theformer is smaller than the latter, the tagged particle’s temperature is eventually slaved by T g due to the second termof the Eq. (22). In the reversed case, the tagged particle’s temperature evolves independently in the long time limitwith a cooling rate slower than that of the gas. It should be emphasized here that in the expansion made in theprevious section, we implicitly assumed that T /T g remains finite, because the coefficients A and N were expandedin powers of ∆ / TT g (see Appendix B). Hence, the Fokker-Planck equation for ǫ > TT g is small enough. In the following, we will only consider the case in which the Fokker-Planck equation isvalid for all times, i.e. the double limit∆ → , p → , ǫ < , ǫ ∝ (cid:18) σ g σ (cid:19) d − p ∆ , (27)where the requirement of small p stems from ∆ ≪ ǫ bounded from above. Hence, in order to be consistent withthis limit, we can substitute in the Fokker-Planck equation, Eq. (18), the values of the coefficients, a ( p ) and b ( p ) bytheir elastic limits lim p → a ( p ) = 1 , lim p → b ( p ) = 1 , (28)so that (cid:20) ∂∂t + v · ∇ (cid:21) F ( r , v , t ) = γ e ( t ) ∂∂ v · (cid:20) v + T g ( t ) m ∂∂ v (cid:21) F ( r , v , t ) . (29)This equation is formally identical to the one obtained for an elastic gas [15], except for the fact that the densityand temperature of the gas depend on time. This is a consequence of the elastic limit to which we are restricted. Ingeneral, the coefficients a ( p ) and b ( p ), Eqs. (16) and (17), differ from unity due to the non Maxwellian character ofthe distribution function of the bath, and Einstein relation is violated.In order to analyze this equation, it is convenient to introduce the dimensionless time scale, t ∗ , proportional to thenumber of collisions of the tagged particle t ∗ = (1 − ǫ ) Z t dt ′ γ e ( t ′ ) , (30)where ǫ is defined by substituting ζ T ( p ) /a ( p ) in (25) by its p → ǫ = √ (cid:18) σ g σ (cid:19) d − p ∆ . (31)The dimensionless time scale t ∗ is related to the real time t by t ∗ = γ e − ǫ ) ν g p (2 ζ n + ζ T ) log[1 + ν g (0) p ( ζ n + ζ T / t ] . (32)In this time scale, the evolution of the mean velocity and temperature of the tagged particle are particularly simple u ( t ∗ ) = u (0) e − t ∗ − ǫ , (33)and T ( t ∗ ) T g ( t ∗ ) = T (0) T g (0) e − t ∗ + 11 − ǫ (1 − e − t ∗ ) . (34)Such predictions will be compared against numerical simulations in section IV. Let us also introduce the scaleddistribution F ( r , v , t ) = 1 v dǫ ( t ) F ∗ ( r , v ∗ , t ∗ ) , v ∗ = v v ǫ ( t ) , (35)where v ǫ ( t ) = (cid:18) − ǫ (cid:19) / (cid:18) T g ( t ) m (cid:19) / . (36)The function v ǫ ( t ) is introduced because with these definitions we have v ǫ ( t ) → (cid:20) T ( t ) m (cid:21) / , (37)in the long time limit. In these variables the Fokker-Planck equation (29) reduces to (cid:18) ∂∂t ∗ + ℓ ( t ∗ ) v ∗ · ∇ (cid:19) F ∗ ( r , v ∗ , t ∗ ) = Λ F P ( v ∗ ) F ∗ ( r , v ∗ , t ∗ ) , (38)where we have introduced the standard homogeneous Fokker-Planck operatorΛ F P ( v ∗ ) = ∂∂ v ∗ · (cid:18) v ∗ + 12 ∂∂ v ∗ (cid:19) , (39)and the function proportional to the mean free path ℓ ( t ∗ ) = v ǫ ( t )(1 − ǫ ) γ e ( t ) = d Γ( d/ − ǫ ) / π d − ∆ − / (cid:18) σ g σ (cid:19) d − [ n g ( t ∗ ) σ d − g ] − . (40)Taking into account the definition of ǫ , Eq. (25), we can write explicitly ℓ as a function of the t ∗ variable as ℓ ( t ∗ ) = ℓ (0) e ǫ ∗ t ∗ , ǫ ∗ = 2 ζ n ǫ ζ T (1 − ǫ ) . (41)To sum up, we have obtained the evolution equation for the distribution function of a tagged particle in a bath ofparticles which annihilate, in the limit where the tagged particle is much heavier than the particles of the bath. Thereare some points in common with the elastic case, but also some important differences. The homogeneous operator,whose spectral properties are well-known [13, 15], is exactly the same, but the flux term is weighted by a functiondepending on time and that diverges in the long time limit. This will have important consequences in the study ofdiffusion as we will see in the following section. Moreover, the equation is not valid for all values of the probability p of annihilation in the bath and masses of the tagged particle but, as already mentioned, is limited to the double limitof Eq. (27), in which ǫ < III. LONG-TIME LIMIT SOLUTION OF THE FOKKER-PLANCK EQUATION
In this section we investigate the long time behavior of the solution of the Fokker-Planck equation, Eq. (38),starting with an arbitrary initial condition. The objective is to study if the tagged particle reaches some scaling statein the long time limit and also to analyze how the particle diffuses.
A. Evolution towards a scaling form
As the Fokker-Planck equation is linear, it is convenient to work in the Fourier space. The Fourier component ofthe tagged particle distribution function is defined as F k ( v ∗ , t ∗ ) = Z d r e − i k · r F ∗ ( r , v ∗ , t ∗ ) , (42)so that Eq. (38) yields F k ( v ∗ , t ∗ ) ∂∂t ∗ F k ( v ∗ , t ∗ ) = [Λ F P ( v ∗ ) − iℓ ( t ∗ ) k · v ∗ ] F k ( v ∗ , t ∗ ) . (43)The spectrum of the operator Λ F P ( v ∗ ) − iℓ ( t ∗ ) k · v ∗ is known [13, 15]. The eigenvalues are λ n ( k , t ) = −
12 [ kℓ ( t )] − d X j =1 n j , (44)where we have introduced the vector label n = ( n , . . . , n d ), with possible coordinate values n i = 0 , , , . . . ∞ .Hence, for any initial condition, all the k -Fourier components decay and only the k = remains. Moreover, as theeigenfunction associated to the vanishing eigenvalue is the Maxwellian distribution [13, 15] χ M ( v ∗ ) = 1 π d/ e − v ∗ , (45)we obtain F ( r , v , t ) → v dǫ ( t ) χ M ( v ∗ ) , (46)in the long time limit.As a consequence, the tagged particle distribution function approaches a scaling form similar to (5) for the gas, butwith a different temperature (see (26)). In this regime the cooling rates of the bath and of the tagged particle arethe same and the temperatures are proportional. The situation is similar to that of an elastic particle in a bath ofinelastic grains [18, 19]. Nevertheless, there is an important difference: in the inelastic case it has been proved thatthere exists an exact mapping with an elastic system. On the other hand, in our problem such a mapping fails dueto the flux term which explicitly depends on time. B. Characteristics of diffusive motion
Our objective is to study the evolution equation for the density of tagged particles, n ( r , t ) = R d v F ( r , v , t ) in a“macroscopic” scale, i.e. in a long time and length scale compared to the microscopic ones. The microscopic timescale is defined by the slowest kinetic modes of Λ F P , i.e. the modes with a single non vanishing component, labeledby n i = δ ij for a given value of j in [1 , d ]. The microscopic length scale is defined by the mean free path of the taggedparticle which is proportional to ℓ ( t ∗ ). The starting point will be the Fokker-Planck equation for F k , Eq. (43). As thegenerator of the dynamics, the operator Λ F P ( v ∗ ) − iℓ ( t ∗ ) k · v ∗ , does not commute with its time derivative, it is notpossible to write the general solution of equation (43) in terms of the initial condition in a simple way. Nevertheless,it is shown in Appendix D that in the hydrodynamic limit, i.e. kℓ ( t ∗ ) ≪ t ∗ ≫
1, aclosed equation for the Fourier component of the density, n k = R d v ∗ F k ( v ∗ , t ∗ ), is obtained ∂n k ( t ∗ ) ∂t ∗ = − D [ kℓ ( t ∗ )] n k ( t ∗ ) , (47)where the diffusion coefficient is D = 12(1 + ǫ ∗ ) . (48)This asymptotic behavior can be evaluated, taking advantage of the scale separation (i.e. the mode with n = isisolated from the other modes). From equation (47) we can derive the evolution equation for the density ∂n ( r , t ∗ ) ∂t ∗ = D ℓ ( t ∗ ) ∇ n ( r , t ∗ ) . (49)This is the equation we were looking for and it is only valid in the “macroscopic” time and length scale. If wetransform this equation to real time with the aid of formula (30) we obtain ∂n ( r , t ) ∂t = D ( t ) ∇ n ( r , t ) , D ( t ) = ζ T (1 − ǫ )[ ζ T + (2 ζ n − ζ T ) ǫ ] D e ( t ) , (50)where D e ( t ) = 2 v g ( t ) /γ e ( t ) is the same as the diffusion coefficient for elastic collisions except that it appears here asa function of the time-dependent gas temperature and density. As can be seen, for our system the diffusion coefficient D ( t ) is far from being a trivial generalization of the elastic diffusion coefficient.Now let us focus on the predictions of our diffusion equation. To this end, we introduce the mean square displacement h r ( t ∗ ) i = Z d r r n ( r , t ∗ ) . (51)If we consider an infinite system, we obtain from equation (49) ∂∂t ∗ h r ( t ∗ ) i = 2 dD ℓ ( t ∗ ) , (52)that will only be valid in the long time limit. It is straightforward to integrate equation (52), taking into account theexplicit formula for ℓ ( t ∗ ), Eq. (41). This gives h r ( t ∗ ) i = dD ℓ (0) e ε ∗ t ∗ − ε ∗ , (53)or in real time h r ( t ) i = dD ℓ (0) ǫ ∗ n [1 + ν g (0) p ( ζ n + ζ T / t ] ζn ζn + ζT − o . (54)As can be seen in Eq. (53), the diffusive behavior is completely different from its elastic or even inelastic counterparts,where it was found that the mean square displacement is proportional to the number of collision per particle [18, 19].The fact that the bath is loosing particles significantly affects this dynamics, and the mean square displacementincreases exponentially in the collision per particle scale. As we will see in the following section, simulation resultsagree well with our theoretical prediction. Roughly speaking, the exponent 4 ζ n / (2 ζ n + ζ T ) is close to 8 d/ (4 d + 1)(approximately 1.77 in two dimensions, and 1.84 in three dimensions). C. Diffusive behavior : an alternative derivation
In the remainder of this section, we show that, under plausible hypothesis, it is possible to re-derive “`a la Einstein”the formula for the mean square displacement, Eq. (53). This derivation has the merit of leading to a Green-Kubolike expression for the diffusion coefficient. We start by writing the mean-squared displacement as h r ( t ) i = Z t dt ′ Z t dt ′′ h V ( t ′ ) · V ( t ′′ ) i . (55)Here, the position and the velocity of the tagged particle, r ( t ) and V ( t ), are considered as a stochastic process and h . . . i denotes an ensemble average over different trajectories. Let us change the variables from t → t ∗ and let us alsointroduce the scaled velocity w ( t ∗ ) ≡ V v ǫ ( t ) , (56) τ T / T g m=20m g m=200m g Figure 1: (Color online) Evolution of the temperature ratio as a function of the number of collisions per particle τ for a systemwith p = 0 .
01 and two values of the tagged particle’s mass, m = 20 m g (∆ = 1 /
20) and m = 200 m g (∆ = 1 / where v ǫ ( t ) is defined in (36). With these definitions we have h r ( t ∗ ) i = 1(1 − ǫ ) Z t ∗ ds γ − e ( s ) Z t ∗ ds γ − e ( s ) h V ( s ) · V ( s ) i = ℓ (0) Z t ∗ ds Z t ∗ ds e ǫ ∗ ( s + s ) h w ( s ) · w ( s ) i , (57)where we have used the definitions of t ∗ and ℓ (0), Eq. (30) and Eq. (41). Now, if we assume that the tagged particleis in the scaled regime, i.e. the temperature is T g ( t ∗ ) / (1 − ǫ ) and that the correlation function h w ( s ) · w ( s ) i isa function of s − s , by integrating in the new variables S = ( s + s ) / s = s − s , the following relation isobtained h r ( t ∗ ) i e ǫ ∗ t ∗ ǫ ∗ dℓ (0) → d Z t ∗ ds h w ( s ) · w (0) i e − ǫ ∗ s , (58)in the long time limit. This formula is the generalization of the Einstein formula for the diffusion coefficient of a heavyparticle in a fluid in the homogeneous decay state. It relates the asymptotic behavior of the mean-square displacementwith the time integral of the autocorrelation function of the velocity weighted by the exponential e − ǫ ∗ s . So far wehave considered the Fokker-Planck equation as the equation for the one-time probability distribution function. Ifwe assume that the velocity of the tagged particle is a Markov process, then the Fokker-Planck equation is also theequation for the conditional probability and we can evaluate easily the correlation function h w ( s ) · w ( s ) i . Takinginto account Eqs. (33) and (36), we have h w ( s ) · w ( s ) i = d e −| s − s | . (59)By substituting this formula into Eq. (58) we re-derive Eq. (53) with the same diffusion coefficient D , that can bewritten in the Green-Kubo form D = 1 d Z ∞ ds h w ( s ) · w (0) i e − ǫ ∗ s . (60) IV. DIRECT SIMULATION MONTE CARLO RESULTS
The objective of this section is to put to the test the main results of the previous sections by means of the directsimulation Monte Carlo method (DSMC). More precisely, we will analyze the temperature and mean velocity evolution,together with the tagged particle diffusion. We have performed DSMC simulations of a system of N g hard disks ofmass m g and diameter σ g which annihilate with probability p or collide elastically with probability 1 − p everytime t* Φ p=0.01, m=38m g p=0.01, m=150m g p=0.1, m=20m g p=0.1, m=60m g Figure 2: (Color online) Scaling function Φ defined in the main text as a function of reduced time t ∗ for different systems with p = 0 . p = 0 .
01 and different values of the tagged particle’s mass. The dashed line is the theoretical prediction of Eq.(62). ∆ ( T / T g ) s t Figure 3: (Color online) Stationary values of the temperatures ratio for a system with p = 0 .
01 and different values of thetagged particle’s mass ∆ = m g /m . Symbols are for the Monte Carlo data and the dashed line shows the long time limit of Eq.(34). two particles meet. Bird’s algorithm [26] has been used. The parameters in all the simulations were m g = 1, σ g = 1, N g (0) = 10 and T g (0) = 1. We have considered only one tagged particle in each simulation, that collides elasticallywith the surrounding bath. The diameter of this particle has been set to unity ( σ = 1) and we have varied the valueof its mass m . The values of the probability of annihilation p of the particles in the bath have been p = 0 . p = 0 .
01, and the results have been averaged over 2 · and 4 · trajectories respectively. For a given value of p ,we have performed a series of simulations for different values of the mass of the tagged particle. Taking due accountof the constraint ǫ <
1, the value of the tagged particle’s mass must be smaller than 10 for p = 0 . for p = 0 . T /T g for a system with p = 0 .
01 and for two values of the tagged particle’smass, that is, m = 20 m g and m = 200 m g , as a function of the number of collisions per particle, τ , defined as τ = 12 Z t dt ′ ν g ( t ′ ) = 12(1 − ǫ ) ν g γ e t ∗ . (61)The initial value of temperature of the tagged particle is T (0) = 0 for all the trajectories, since at t = 0, the intruderhas a prescribed velocity. As we can see, in this scale, the evolution to the stationary value of the ratio of thetemperatures is slower as we increase the mass of the tagged particle. This implies that there are values of the taggedparticle’s mass for which the ratio of the temperatures will not reach its stationary value in the time of the simulation(for instance, the number of particles in the bath for τ = 200 is N g ≃ τ < v > / < v > m=15m g m=20m g m=60m g Figure 4: (Color online) Reduced fourth moment of the tagged particle velocity as a function of the dimensionless time τ for p = 0 . t*/(1- ε ) -5-4-3-2-10 l n [ u x ( t * ) / u x ( ) ] t*/(1- ε ) u x ( t * ) / u x ( ) m=15m g m=20m g m=38m g m=60m g m=80m g m=100m g Figure 5: (Color online) Mean velocity of the tagged particle as a function of t ∗ / (1 − ǫ ) for a system with p = 0 . t ∗ / (1 − ǫ ) for a system with p = 0 . The theoretical prediction in this time scale (dashed line) is obtained directly from Eq. (34), taking into account Eq.(61). The agreement between theory and simulations is good.Similar simulations were performed with different values of the tagged particle mass ( m ranging from 15 to 100 for p = 0 . p = 10 − ). Since T g (0) = 0, Eq. (34) predictsΦ ≡ (1 − ǫ ) TT g = 1 − e − t ∗ . (62)As can be observed in Fig. 2, all simulation data for Φ collapse onto a single curve. The time scale t ∗ can be calculatedfrom the scale τ defined in (61). In the same vein, we can obtain the stationary value of the temperature ratio, whichis plotted as a function of ∆ = m g /m in Fig. 3 ; this validates our theoretical analysis.In order to probe –at least partially– the Gaussian nature of the time dependent tagged particle velocity statistics,we have measured the reduced fourth moment 4 h v i / ( d ( d + 2) h v i ). As can be seen in Fig. 4, where we have plottedthe results for p = 0 . τ , the value ofthis quantity is in agreement with the Gaussian prediction (that is unity) within the statistical uncertainties.Consider next the mean velocity of the tagged particle, u ( t ∗ ). In order to study the decay of this quantity, we haveperformed a set of simulations starting with a component of the velocity field in the x direction when it is immersedin a bath in the homogeneous decay state, u x (0). In Fig. 5, we plot u x ( t ∗ ) /u x (0) as a function of t ∗ / (1 − ǫ ) for p = 0 . m . In this time scale, the data for all values of m collapse due to Eq. (33). In the inset,the same quantity is plotted on a logarithmic scale. If the theoretical prediction in Eq. (33) is verified, the aboveplot must lead to a straight line with slope ζ u = 1 (dashed line), where ζ u is the decay rate of the mean velocity. On1 ∆ ζ u Figure 6: Decay rate of the mean velocity as a function of ∆ for a system with p = 0 .
1. The symbols are from DSMC simulationsan the dashed line is the theoretical prediction. t* × × ×
01 and m = 60 m g as a function of t ∗ . Thedashed line is the theoretical prediction, Eq. (53). Inset : same quantity as a function of n g (0) /n g ( t ∗ ) . The dashed line isthe theoretical prediction of Eq. (63) where B follows from (65). the other hand, ζ u can be fitted on the logarithmic plot of the inset. Reporting the corresponding measures in Fig.6 against the mass ratio, it appears that the theoretical prediction ζ u = 1 is approached as ∆ →
0, as expected.Finally, the accuracy of the prediction for the diffusion equation has also been tested by measuring in DSMCsimulations the mean square displacement of the tagged particle. In contrast to the granular case phenomenologyand as a consequence of the continuous decay of particles in the bath, the mean squared displacement increasesexponentially in the t ∗ scale, see Eq. (53). In Fig. 7 we have plotted the time evolution of h r i in the scale t ∗ fora system with p = 0 .
01 and m = 60 m g . The dashed line is the theoretical prediction given by Eq. (53), and showsgood agreement with numerical data. The same quantity, written in terms of the bath density, reads h r ( t ∗ ) i = dD ℓ (0) e ǫ ∗ t ∗ − ǫ ∗ = B (cid:20) n g (0) n g ( t ∗ ) − (cid:21) , (63)where we have defined B = dD ℓ (0) ε ∗ , (64)and we have taken into account that n g ( t ∗ ) = n g (0) e − ǫ ∗ t ∗ . It then appears that the mean squared displacementincreases linearly with n g (0) /n g ( t ∗ ) . This is full agreement with the simulation results, see the inset of Fig 7. Sucha plot allows us to extract by linear fitting the coefficient B , which can then be compared against the prediction of2 ∆ B n g2 ( ) ∆ B n g2 ( ) Figure 8: Values of Bn g (0) as a function of ∆ for systems with p = 0 . p = 0 .
01 (right). The dashed line is for theprediction of Eq. (65).
Eq. (64), which explicitly reads B = − √ d π − d σ − d ) g Γ (cid:0) d (cid:1) n g (0) p ( −
16 + √ p ∆ )(16 + √ d − p ∆ ) , (65)where we have taken into account that σ = σ g . Such a comparison is worked out in Figure 8 and fully corroboratesthe theoretical analysis, with again an improved agreement when ∆ decreases.It would be also interesting to confirm our theoretical predictions with Molecular Dynamics simulations. In the lowdensity limit, it is expected to get similar results. In fact, some simulations were performed finding qualitatively thesame behavior but with much more statistical inaccuracies. V. CONCLUSIONS
In this paper, the diffusive behavior of a tagged intruder immersed in a gas of particles undergoing ballistic anni-hilation (i.e. which annihilate with probability p or scatter elastically otherwise), has been analyzed. The collisionsbetween the tagged particle and the surrounding gas are elastic. Some similarities are found between our systemand the elastic or inelastic case [15, 18], but, on the other hand, important differences arise as a consequence of thecontinuous decay of particle number in the system.We start from the Boltzmann-Lorentz equation for the distribution function of the tagged particle, which is valid,in principle, for arbitrary mass of the tagged particle. In the limit of a very massive tagged particle, a Fokker-Planckequation for the distribution function is derived by means of a systematic expansion in the mass ratio ∆. Our approachholds in the limit ∆ ≪
1, but we additionally have the more stringent condition that the parameter introduced inEq. (25), ǫ ∝ p ∆ , must be smaller than unity. Analysis of the Fokker-Planck equation leads to predictions for thetemperature ratio, the decay rate of the mean velocity of the tagged particle and the diffusion coefficient. As in theinelastic case [18], the theory predicts that the ratio between the temperatures of the tagged particle and the gas isconstant in the long time limit as a consequence of equilibrating cooling rates. When represented in the appropriatetime scale, which is proportional to the number of collisions experienced by the tagged particle, temperature ratioscollapse for all values of ∆ and p considered. Likewise, the mean intruder velocity (averaged over bath realizations)decays exponentially. The dynamics of the distribution function of the tagged particle is governed by a Fokker-Planckoperator which spectral properties are known. More specifically, as the eigenvalues of this operator are non-positive,the distribution of the tagged particle approaches a Gaussian in the long-time limit and admits a scaling form similarto the one for the distribution function for the particles in the gas but with a different temperature. At variance withthe situation of an elastic intruder in a bath of inelastic particles [22], there is apparently no mapping between ourproblem and a well chosen elastic system. A unique vanishing eigenvalue is responsible for the slow diffusive behaviorof the tagged particle density. The corresponding diffusion equation has been derived in the hydrodynamic limit, bymeans of a projector decomposition, which yields an explicit expression for the diffusion coefficient. From a differentpoint of view, the expression for the mean squared displacement has also been derived “`a la Einstein”. Followingthis route, the diffusion coefficient is expressed as a Green-Kubo formula in terms of a weighted time integral of thetagged particle velocity correlation function. This provides a more physical perspective on the results derived from the3projector method. As already mentioned, the mean squared displacement for this system does not increase linearlyin the collision per particle time scale, as is the case in the elastic and inelastic cases. This different behavior is dueto the time dependent bath density. In the elastic case, both the temperature and the density do not depend ontime. In an inelastic system [18, 19], the time dependence goes through the temperature and could be absorbed inthe collision per particle time scale, which turns out to be impossible in our system where the mean free path ℓ ( t ∗ ) isan increasing function of time.Finally, our analytical results have been tested by numerical simulations, and a very good agreement has beenreported for all the range of parameters considered. As expected, the agreement is all the better as ∆ is smaller.In summary, the work reported here provides an example of the accuracy of hydrodynamics to describe a system inwhich there are no conserved quantities in binary encounters (no collisional invariants). Acknowledgments
We acknowledge useful discussions with G. Schehr and A. Barrat. We would like to thank the Agence Nationale dela Recherche for financial support (grant ANR-05-JCJC-44482). M. I. G. S. and P. M. acknowledge financial supportfrom Becas de la Fundaci´on La Caixa y el Gobierno Franc´es. M. I. G. S. would like to thank the HPC-EUROPAproject (RII3-CT-2003-506079), with the support of the European Community Research Infrastructure Action, forfinancial support.
Appendix A: SOME USEFUL APPROXIMATIONS
For the sake of completeness, we provide here the approximate expressions for the density and temperature decayrates, that are relevant for explicit computation of several of the quantities discussed in the main text. They havebeen obtained from a truncated Sonine expansion (Sonine polynomials being particular types of Laguerre polynomials,particularly convenient for kinetic theory calculus) [4, 5, 6]. ζ n = d + 24 (cid:18) − a (cid:19) , (A1) ζ T = d + 28 d (cid:18) a d + 1116 (cid:19) , (A2)where a = 8(3 − √ p (4 d + 6 − √ p + 8 √ d − − p ) . (A3) Appendix B: FROM THE BOLTZMANN-LORENTZ EQUATION TO THE FOKKER-PLANCKEQUATION
In this Appendix we expand the collision operator, Eq. (2), in series of ∆. We start by multiplying the collisionoperator by a generic function H ( v ) and integrate in velocity space Z d v H ( v ) J [ r , v , t | F, f ]= σ d − Z d v Z d v H ( v ) Z d ˆ σ Θ( g · ˆ σ )( g · ˆ σ )[ F ( v ∗ ) f ( v ∗ ) − F ( v ) f ( v )] . (B1)The above expression can be written Z d v H ( v ) J [ r , v , t | F, f ]= σ d − Z d v Z d v F ( v ) f ( v ) Z d ˆ σ Θ( g · ˆ σ )( g · ˆ σ )[ H ( v − δ v ) − H ( v )] , (B2)4where we have introduced δ v = 2∆1 + ∆ ( g · ˆ σ ) ˆ σ , (B3)which is the increment of the tagged particle velocity due to collisions with a particle of the bath (it should beremembered that g = v − v ). Equation (B2) essentially tells us how the function H varies due to collisions. If weadmit that ∆ is small enough, we can expand H ( v − δ v ) around v in powers of δ v , keeping only the lower orders H ( v − δ v ) ≃ H ( v ) − (cid:20) ∂H ( v ) ∂ v (cid:21) · δ v + 12 (cid:20) ∂∂ v ∂∂ v H ( v ) (cid:21) : δ v δ v . (B4)If we introduce expansion (B4) in equation (B2) we obtain Z d v H ( v ) J [ r , v , t | F, f ] ≃ Z d v H ( v ) (cid:26) ∂∂ v · [ A ( v ) F ( v )] + 12 ∂∂ v ∂∂ v : [ N ( v ) F ( v )] (cid:27) , (B5)where we have introduced A ( v ) = 2∆ π d − σ d − (1 + ∆)Γ (cid:0) d +32 (cid:1) Z d v f ( v ) g g , (B6) N ( v ) = (cid:20) (cid:21) π d − σ d − Γ (cid:0) d +52 (cid:1) Z d v f ( v ) (cid:20) d + 32 d g I + 32 g (cid:18) gg − d g I (cid:19)(cid:21) . (B7)In the last expression I is the unit tensor. As H ( v ) is a generic function of v , we can compare equations (B1) and(B5), and we obtain that the collision operator can be written as J [ r , v , t | F, f ] ≃ ∂∂ v · [ A ( v ) F ( v )] + 12 ∂∂ v ∂∂ v : [ N ( v ) F ( v )] . (B8)We next specify A and N within the scaling form provided by the homogeneous decay state of the bath. This willlead us to identify the remaining ∆ dependence in these coefficients and to simplify the functional dependence in thetagged particle velocity. To this end, we introduce the dimensionless velocities c = v v ( t ) , c = v v g ( t ) , (B9)where v g ( t ) = (cid:16) T g ( t ) m g (cid:17) / and v ( t ) = (cid:16) T ( t ) m (cid:17) / , with T ( t ) the temperature of the tagged particle and T g ( t ) thetemperature of the suspending gas. The relative velocity g = v − v can be written as g = v g ( t ) (cid:20) T ( t ) T g ( t ) (cid:21) / ∆ / c − v g ( t ) c . (B10)A formal expansion in ( T /T g )∆ leads to A ( v , t ) = γ ( t ) v , N ( v , t ) = 2¯ γ ( t ) I . (B11)The definitions of γ and ¯ γ are respectively γ ( t ) = γ e [ n g ( t ) , T g ( t )] a ( p ) , (B12)¯ γ ( t ) = γ e [ n g ( t ) , T g ( t )] a ( p ) b ( p ) T g ( t ) m , (B13)where γ e ( t ) is the same friction coefficient as for an elastic system at the corresponding density and temperature γ e [ n g ( t ) , T g ( t )] = 4 π d − d Γ( d/
2) ∆ / n g ( t ) (cid:18) T g ( t ) m (cid:19) / σ d − , (B14)5and a , b are functionals of the distribution function of the bath which depend only on the parameter pa ( p ) = Γ( d/ d + 1) / Z d c χ H ( c ) c , (B15) b ( p ) = 2 d + 1 R d c χ H ( c ) c R d c χ H ( c ) c . (B16)These coefficients have been evaluated in the first Sonine approximation and depend very weakly on pa ( p ) = 8 − a ( p )8 , (B17) b ( p ) = 8 + 3 a ( p )8 − a ( p ) , (B18)where a , defined in (A3), is the gas velocity distribution kurtosis (a Gaussian ansatz would amount to setting a = 0).By dimensional analysis and taking into account the explicit formulas for A and N , equation (B11), we can see thatthe two terms we have considered in the expansion of the collision operator, Eq. (B8), are of order n g v g σ d − ∆, whilethe other terms in the Kramers-Moyal expansion are at least of order n g v g σ d − ∆ / . Hence, we can conclude thatthe leading order contribution in ∆ of the collision operator is actually the one written in (B8). Appendix C: EQUATIONS FOR THE VELOCITY AND TEMPERATURE OF THE TAGGED PARTICLE
In this Appendix we derive the equations for the mean velocity and temperature of the tagged particle. Takingmoments in the Fokker-Planck equation, Eq. (18), we obtain for the velocity ∂ u ( t ) ∂t = Z d r Z d vv (cid:26) − ( v · ∇ ) + γ e ( t ) a ( p ) ∂∂ v · v + a ( p ) b ( p ) T g m γ e ( t ) ∂ ∂ v (cid:27) F ( r , v , t ) . (C1)By integration we have ∂ u ( t ) ∂t = Z d r Z d vv γ e ( t ) a ( p ) ∂∂ v · ( v F ( r , v , t )) (C2)= − Z d r Z d v γ e ( t ) a ( p ) v F ( r , v , t ) (C3)= − γ e ( t ) a ( p ) u ( t ) . (C4)Consequently, the equation for the mean velocity is ∂ u ( t ) ∂t = − γ e ( t ) a ( p ) u ( t ) . (C5)Taking into account the definition of temperature, d T ( t ) = Z d r Z d v m [ v − u ( t )] F ( r , v , t )= Z d r Z d v m [ v − u ( t )] F ( r , v , t ) . (C6)we can write d ∂∂t T ( t ) = m (cid:20) ∂∂t Z d r Z d v v F ( r , v , t ) − u ( t ) · ∂ u ( t ) ∂t (cid:21) . (C7)In order to evaluate the first term on the right hand side, we make use of the Fokker-Planck equation: ∂∂t Z d r Z d v v F ( r , v , t ) = − γ e ( t ) a ( p ) Z d r Z d v v F ( r , v , t ) + 2 dm T g γ e ( t ) a ( p ) b ( p ) . (C8)Taking this formula and the equation for the velocity into account, we obtain d ∂T ( t ) ∂t = − dγ e ( t ) a ( p )[ T ( t ) − b ( p ) T g ( t )] . (C9)6 Appendix D: THE DIFFUSION EQUATION
In this Appendix we derive the diffusion equation for the tagged particle’s density. The starting point is theFokker-Planck equation (43) ∂∂t ∗ F k ( v ∗ , t ∗ ) = [Λ F P ( v ∗ ) − iℓ ( t ∗ ) k · v ∗ ] F k ( v ∗ , t ∗ ) , (D1)in which we introduce the two projectors P g ( v ∗ ) = h χ M ( v ∗ ) | g ( v ∗ ) i χ M ( v ∗ ) , (D2) P ⊥ g ( v ∗ ) = (1 − P ) g ( v ∗ ) . (D3)Here, we have introduced the maxwellian distribution, χ M ( v ∗ ), which is the eigenfunction of Λ F P associated with the0 eigenvalue and we have used the scalar product defined as h f ( v ∗ ) | g ( v ∗ ) i = Z d v ∗ χ − M ( v ∗ ) f † ( v ∗ ) g ( v ∗ ) , (D4) f † being the complex conjugate of f . In a next step, we decompose the function F k in P F k and P ⊥ F k , and write theequations for these two quantities (cid:20) ∂∂t ∗ + iℓ ( t ∗ ) P k · v ∗ − P Λ F P (cid:21)