{'fraction_non_alphanumeric': 0.05311546335998764, 'fraction_numerical': 0.016378877428825277, 'mean_word_length': 4.168297065282491, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Perturbation theory for large Stokes number particles in random velocity fields

We derive a perturbative approach to study, in the large inertia limit, the dynamics of solid particles in a smooth, incompressible and finite-time correlated random velocity field. We carry on an expansion in powers of the inverse square root of the Stokes number, defined as the ratio of the relaxation time for the particle velocities and the correlation time of the velocity field. We describe in this limit the residual concentration fluctuations of the particle suspension, and determine the contribution to the collision statistics produced by clustering. For both concentration fluctuations and collision velocities, we analyze the differences with the compressible one-dimensional case.


Introduction
An important component in the transport of aerosols by turbulent flows is the tendency to form clusters, an ubiquitous phenomenon that has been observed e.g. inside clouds [1], and whose simplest instance is particle aggregation in a 1D (one dimensional) random force field [2]. In 1D, clustering is the combined result of particle slipping with respect to the random field, and the fact that the random forces pushing particles apart become smaller as the separation decreases [2,3]. In more than 1D, the picture is more complicated, and the process is accompanied by preferential concentration of the particles (supposed denser than the fluid) in the strain regions of the flow [4]. Thus, contrary to the intuition that would suggest a mixing behavior, a spatially homogeneous random velocity field will lead to the formation of clumps out of an initially uniform distribution [5]. This unmixing will take the form of concentration fluctuations, superimposed to a mean concentration field that remains uniform, and exeeding those due to discreteness in the distribution, de-scribed by Poisson statistics [6]. Additional segregation will be produced by inhomogeneity of the turbulent flow, such as in wall turbulence, resulting in a non-uniform mean concentration profile [7]. Gravity will affect the concentration fluctuations of heavier particles [8,9], decreasing the effective correlation time of the turbulent fluctuations sampled by the particles [10].
A practical motivation for the interest in clustering is clearly the possibility of enhanced binary collision, compared to a spatially homogeneous condition, and this could find application e.g. to rain formation [11,12]. It must be mentioned that there have been suggestions that clustering is of secondary importance in the collision dynamics [14,15]. It is not clear, in general, how clustering affects the relative velocity dynamics, and in this way also collisions [13].
Several models for the clustering dynamics have been proposed [5,16,17], however, although high Reynolds number turbulence is a multiscale flow, the bulk of the studies has been on the dynamics of inertial particles in smooth random fields (see e.g. [18,19]). There are good reasons for this. The dynamics of a sufficiently small (and sufficiently dense) spherical particle can be described by the relaxation time of its velocity relative to the fluid: the Stokes time τ S = 2/9 r 2 0 λ/ν 0 , where r 0 is the particle radius, λ is the ratio of the particle to fluid density and ν 0 is the kinematic viscosity of the fluid [20]. Now, most atmospheric aerosols are characterized by Stokes times shorter than the Kolmogorov times of the flows by which they are advected [11]. Furthermore, experimental data [21] and numerical simulations [4] both indicate that clustering is stronger for particles with τ S of the order of the Kolmogorov time, corresponding to a range of scales in which the velocity field could be approximated as smooth.
In the case of a smooth incompressible random field, clustering is peaked at S = τ S /τ E ∼ 1, where τ E is the correlation time of the field and S is called the Stokes number. For an incompressible velocity field, the limit S → 0 corresponds to passive scalar transport, and spatially homogeneous initial conditions will not lead to clustering. In the opposite limit S → ∞, one expects that the particles be scattered by the velocity fluctuations they cross in their motion [22] as if undergoing Brownian diffusion, resulting again in no clustering. One reason for being interested in this limit is that also particles in turbulent flows for which τ S lies in the inertial range (or above), see smaller vortices as if S ≫ 1. It must be said that the large S limit is more of interest for industrial application, since in cloud turbulence, the relative motion of heavier droplets is dominated by the different gravitational settling velocity of particles of different size [11] and the possibility of chaotic trapping [23,24].
We shall focus in this paper on the case of a strongly stirred fluid, such that, also for large-S particles, the largest component of the particlefluid relative motion is due to inertia and not gravity. We will try to understand, in particular, how the Brownian limit of [22] is achieved.
There is some evidence [25] that clustering destruction occurs at finite S, as the result of a crossover of the correlation dimension D 2 of the particle distribution, above the dimension of the space D. For S below this threshold, the probability density function (PDF) ρ(r) of finding a pair of particles at separation r diverges for r → 0 like r D2−D (see also [18,26]), while it should remain finite above theshold. Actually, residual (non-singular) concentration fluc-tuations remain present above threshold, decaying in 1D like ρ(r = 0) ∝ S −1/2 , and their origin lies in the fact that only particles moving at increasingly small relative velocities, as S → ∞, are able to stay close long enough to remain correlated [13].
The picture just described suggests the possibility of a perturbative expansion around S → ∞. This is formally identical to an expansion in powers of the Peclet number P e, for particle advection in the presence of strong molecular diffusion. In both the large S and the small P e regimes, the small quantity is the contribution to the relative particle motion produced at small separation, by the correlations in the random field. The perturbative approach will allow to carry on an analysis of the large-S inertial particle dynamics for D > 1, where the heuristic approaches like the one in [13] are made difficult by the complexity of geometry and incompressibility.
This paper is organized as follows. In Sec. 2, the model equations for the problem are derived. The perturbative approach is derived in Sec. 3. In Sec. 4, the approach is tested in 1D, comparing with numerical simulations and with the heuristic results in [13]. In Sec. 5, concentration fluctuations are analyzed in the case of a 3D incompressible random field. Section 6 is devoted to determination of the effect of concentration fluctuations on the relative particle velocity statistics. Section 7 contains the conclusions. Discussion of additional technical aspects is left to the appendices.

The stochastic model
Following the derivation in [13], we introduce a zero mean, smooth Gaussian random velocity field u(x, t), with correlation where g αβ (0) = δ αβ , F (0) = 1, and we can introduce correlation times and lengths τ E and r v obeying ∞ 0 dτ F (τ ) = τ E , and ∞ 0 drg αβ (rr) ∼ r v (r ≡ r/r fixed). We assume isotropy and homogeneity of the field in space and time.
From the field parameters r v , τ E and σ u , we can introduce the Kubo number K = σ u τ E /r v , which tells us whether the field is intrinsically short-or long-time correlated. The K → 0 limit would correspond to a random field with zero correlation time, like in the Kraichnan model [27]. The K → ∞ limit would correspond to a frozen field regime. In most of our analysis we shall assume K = O(1), as in actual turbulent flows.
A suspension of inertial particles is advected by the random field and their velocity v is taken to obey the Stokes equation: where we have allowed for a Brownian motion component η: η α (t)η β (0) = κτ −2 S δ αβ δ(t); κ is the molecular diffusivity of the particles.
We choose units such that σ u = τ S = 1; therefore: In the regime S ≫ 1, the particle displacement in a correlation time τ E becomes negligible with respect to r v and Eq. (2) could be approximated by a Langevin equation (the precise condition for a Langevin dynamics is K 2 ≪ S [13]). It is then possible to substitute into Eq. (2): with ξ(t) white noise: ξ α (t)ξ β (0) = δ αβ δ(t), ξ α η β = 0. We thus obtain for the particle velocity variance σ 2 v ∼ S −1 (1 + P e −1 ), where the Peclet number P e = 2σ 2 u τ E /κ parameterizes the relative strength of the random field and molecular contributions to particle diffusion.
Turning to the relative motion of particle pairs, let us introduce difference variables where 1 and 2 label members of a particle pair.
(The particles are assumed immaterial, so that they can cross without interaction). As in the one-particle case, we can approximate the equation for the relative motion of particles with the Langevin equation: (in this paper, summation over repeated indices is assumed), with a = (1 + P e) −1 P e and q = 2(Sa) −1/2 .
We see that a goes to zero with P e → 0, and, as expected, the relative particle dynamics described in Eqs. (3)(4) becomes uncorrelated. The zero molecular diffusion limit corresponds instead to P e → ∞ and a = 1. It should be mentioned that, in the Langevin equation limit, the dynamics becomes equivalent to that of a Kraichnan model and can be described in terms of the single parameter [3] To understand the connection with the Kraichnan model, notice that ǫ −1 τ S ∼ r 2 v /(σ 2 u τ E ) is the diffusion time of a single passive tracer across a distance r v and plays the role of effective correlation time of the field, so that ǫ becomes an effective Stokes number. More interestingly, for P e, ǫ ≫ 1, ǫ −1/2 τ S is the permanence time at separation < r v = (KS) −1 of a particle pair approaching at relative velocity σ v ∼ S −1/2 (it is easy to see that such particles cross at large enough speed to behave ballistically at scale r v ).
The limit S ≫ K ∼ 1, beyond allowing a Langevin equation based model, leads to decoupling of the difference variables (ν, r) from the center of mass variables 1 2 (x 1 + x 2 ) and . At equilibrium, this means that the two-particle PDF will be in the form, using spatial homogeneity: where Ω −1 , with Ω the volume of the domain for x, is just ρ((x 1 +x 2 )/2). (Where not ambiguous, we do not use subscripts to identify PDF's of different quantities). Multiplying by N 2 , with N the total number of particles in the domain and integrating over d 3 v 1 d 3 v 2 , we obtain the expression for the concentration correlation: wheren = N/Ω is the mean concentration. The quantity f (r) = Ωρ(r)−1, with ρ(r) = d 3 νρ(ν, r), gives the strength of the concentration fluctuations (we have exploited isotropy). The concentration variance can be expressed in terms of the function f by means of the relation (n − n) 2 =n 2 f (0).

Perturbation theory
For small P e, the two-particle PDF ρ(ν, can be determined solving Eqs. (3)(4) for ν(t) and r(t) perturbatively in a. We shall argue and verify numerically in the next section, that the same perturbative strategy remains valid for S ≫ K ∼ 1 also in the absence of molecular diffusion, i.e. for a = 1.
The perturbative expansion of ν(t) and r(t) is obtained substituting the Taylor expansion which is obtained from Eq. (4), into the solution of Eq. (3): and Taylor expanding again in a. The ground state of the expansion is therefore: while the first two corrections are: and b (2) where ∂ γ ≡ ∂/∂r γ . We focus here on the PDF for the particle separation r, leaving the analysis of the joint PDF ρ(ν, r) to Sec. 6. The PDF for r, given some distribution of initial conditions (r(−T ), ν(−T )), will be something in the form: Depending on the choice of initial conditions for the trajectories ending at r, ρ(r; t) will be in general a non-equilibrium PDF. Carrying out the integral over r (0) : Taylor expanding in a leads to a perturbative series ρ = ρ (0) + ρ (1) + . . ., with ρ (0) (r, t) ≡ ρ(r (0) (t) = r). In explicit form: and similarly for the other conditional averages. Sending T → ∞ and using lim T →∞ ρ(r (0) ) = ρ (0) (r) = Ω −1 , we obtain the perturbative expansion for the equilibrium PDF: Clearly, the dependence on the condition r (1) ( −T ) = 0 in Eq. (11) disappears when T → ∞, but the limit must be taken after the average. In order to take into account the simultaneous conditioning in the past and in the present in Eq. (11), the following expression for r (1) |r can be utilized: where the factor Ω is simply [ρ(r (0) (t)|r (1) (−∞) = 0)] −1 . Substituting Eqs. (8) and (9) into (13) and setting without lack of generality t = 0, we obtain: Using the functional integration by part formula [28], allows to treat the correlation between ξ β (τ ) and the other factors in the integral: γ (0) → −∂ γ and using Eq. (9), we obtain finally: Similar expressions can be derived for the higher orders in the expansion for ρ, and they will be required to take into account incompressibility in the D > 1 case.

Small Peclet vs. large Stokes
The two-particle distribution in Eq. (12) is expressed as a formal (regular) perturbation expansion in a. The real expansion parameter, however, is the ratior/r v , wherer ∼ r (1) is the contribution of the velocity correlation g αβ to the relative particle displacement. A relevant quantity is then the permanence time at r (0) < r v , i.e. the time during which g αβ (r (0) ) is significantly different from zero and can contribute to the integral in Eq. (14). If r < r v , the permanence time will be . Now, our perturbation expansion is carried on in powers of a for fixed q; in other words, a (and therefore P e) are considered small independently of q. Had K and S, in place of q, been chosen as the fixed quantities in our problem, we would have obtained instead, from Eq.
This suggests that the perturbation expansion of Eq. (12) could be translated, for large ǫ and P e generic, into a new expansion in powers of ǫ −1/2 . The validity of the new expansion rests on regularity of the solution ρ(r) in r = 0 (the estimate ∂ α , ∂ β ∼ r −1 v adopted implicitly this assumption). While this is guaranteed for small a by molecular diffusion, for a = 1, b αβ (0) = 0 and Eq. (3) becomes singular. Now, the correspondence between the small a and large ǫ cases guarantees finiteness of the terms in the series (it is possible to see that this is not verified in the case a ground state of ballistic particles at scale r v is adopted [13]). Hence, expanding in powers of ǫ −1/2 is a meaningful procedure. It is still possible that the perturbative approach is simply unable to catch singular behaviors [notice that Eq. (15) does not reveal break-down of the perturbation theory for a = 1 and ǫ ≪ 1, in spite of the singularity of ρ(r)]. The numerical evidence in [25], that singular behavior are absent for large ǫ, however, suggests that this is not the case.
Substituting into Eq. (17), we obtain the result wherer = r/q. Notice that q is the typical separation of a pair of particles at a time ∼ τ S = 1 after crossing (r played the role of outer scale in the matching asymptotics analysis in [13]). In other words, the correlation length of the concentration fluctuations is not expected to be of the order of r v , rather, of the larger scale q ∼ (ǫ/a) 1/2 r v .

The 1D case: small ǫ regime
Taking the limit ǫ → 0 under the condition K 2 /S ≪ 1, guarantees that the Langevin approach of Sec. 2 continues to be valid. (Notice that the above conditions imply K ≪ 1). This regime has been extensively studied in [3,19] as a stochastic model for the behavior of low inertia S ≪ 1 particles in realistic turbulent flows, for which K ∼ 1.
The ǫ → 0 limit is diffusive at scale r v for the variable r(t) and we can derive an exact expression for ρ(r), valid for all a < 1. This expression can then be used to test the perturbative results of the previous section. A diffusive limit means that the variable r(t) changes in a correlation time τ S much less than the spatial scale r v [on the contrary, for most pairs, the large ǫ limit is ballistic at scale r v and becomes diffusive only at the scale (ǫ/a) 1/2 r v ]. The drift and diffusivity in the effective equation for r are the leading contributions to [r(t) − r(0)]|r(0) /t and [r(t) − r(0)] 2 |r(0) /t, for t ≫ 1 and |r(t) − r(0)| ≪ r v . Actually, it is easy to prove from moment analysis of the Fokker-Planck equation for Eqs. (3), that the drift ṙ|r = ν|r is zero [13]. We only have to calculate the diffusivity; from Eq. (3), we can write, for t ≫ 1 and |r(t) − r(0)| ≪ r v : and [r(t) − r(0)] 2 |r(0) /t ≃ B(r(0)). The particle separation behaves therefore like a nonuniform Brownian motion:ṙ = b(r)ξ and from the associated Fokker-Planck equation we obtain the exact expression: Notice the quadratic divergence of the PDF at r → 0 for a = 1, that is approached at t → ∞ with an O(ǫ) rate [3]. Let us compare with Eq. (19). In the limit ǫ → 0, we can write G r (ǫ, t) = t + a/(4ǫ) and Eq. (19) becomes, after a few manipulations: but this is precisely the O(a) contribution to Eq. .

The 1D case: large ǫ regime
The prediction in [13] that f (0) ∼ ǫ −1/2 for ǫ → ∞ seems to be supported by Eqs. (15) and (19). However, the integral in Eq. (19) presents singularities for ǫ → ∞, and we must put some care in the analysis. Let us consider first the case r = 0, and write Eq. (19) in the form where ǫ −1/2 ≪ T ≪ 1; the dots indicate the remnant of the integral in Eq. (19), which is shown by inspection to be finite in the limit ǫ → ∞. The first integral in the formula above, instead, is dominated by the ballistic crossing time t ∼ ǫ −1/2 and is logarithmically divergent for ǫ → ∞. We obtain the leading order expression for f (0): The logarithmic divergence is eliminated and a pure power law is recovered, provided we coarse grain the function f at a fixed scale R: where w(r) is a smoothing function with w(r) > 0 and w(r)dr = 1. A simple analytic expression for f R is obtained choosing a Gaussian w(r) = (2π) −1/2 exp(−r 2 /2): These results confirm the heuristic prediction in [13]. The interesting point is the improved performance of perturbation theory at large ǫ, as illustrated in Fig. 1. The collapse of the rescaled profiles indicates a region in which the higher orders in the perturbation expansion are playing a negligible role. The perturbation expansion in a can therefore be converted into one in ǫ −1/2 , as claimed.
We focus next on the correlation profile f (r). We can obtain analytical expressions for largē r = r/q. We rewrite Eq. (19) in the form . We then substitute into Eq. (22) and integrate by parts, using the fact that R and all its derivatives are zero at both zero and infinity. The integrand in Eq. (22) takes the form:

Relative velocity statistics
The approach in Sec. III can be extended to the calculation of the joint PDF ρ(r, ν), with Eqs. (12,13) being replaced by and where y = (r, ν) and the vector indices run now from 1 to 6. Clearly, integrating Eq. (25) over ν leads to Eq. (12). From Eq. (25), we could in principle obtain information on the velocity components that contribute the most to clustering. Looking at Eq. (25), we see that the terms contributing to ρ (k) (r) for k > 0 are those containing only spatial derivatives (the others give zero after integration over ν). We can thus write: where Ω −1 ρ B (ν|r) contains the terms in the RHS of Eq. (25) involving velocities derivatives ∂ να ≡ ∂/∂y 3+α that do not contribute to Eq. (12); of course, d 3 ν ρ B (ν|r) = d 3 ν ρ cl (ν |r) = 0. Substituting into Eq. (25): We may thus interpret ρ cl (ν|r) as the cluster contribution to the velocity PDF ρ(ν|r) at separation r. It is important to notice that this contribution is not necessarily localized in the clusters: different moments with respect to ν of ρ cl (ν, r) = ρ cl (ν|r)ρ(r) do not come necessarily from the same spatial spots in the volume Ω.
To understand what happens, we derive from Eq. (3) the analog of the lowest order equation for the position (8), in the case of the velocity: Integrating Eq. (8) with initial conditions (ν = 0, r = 0), we see that the particles separate at r (0) = r v at a time ∼ ǫ −1/3 τ S , which is the lowest order estimate for τ exit (0|0). From Eq. (29), ǫ −1/6 σ v appears to be the relative velocity at t = τ exit (0|0), and is therefore the escape velocity out of the correlated region r r v for slow particle pairs, that do not behave ballistically in that region [13]. Notice that there are no singularities at ν = 0, and, in contrast with a perturbation theory with a ground state of ballistic particles, the individual terms in the expansion remain finite as ν → 0.
(30) and we recall that f (r) = Ωρ(r)−1 [see Eq. (7)]. Again, all conditional averages are intended in the sense of Eq. (11). Analogously to the analysis in Secs. IV and V, we see that the first non-zero contribution in the incompressible case arises at second order: while in the 1D case: (32) The calculation of the cluster contribution to the collision velocity variance at this point is a matter of lengthy but straightforward algebra. Leaving the calculation details to Appendix B, the result in 1D and in the incompressible 3D case are [see Eqs. (B2,B6)]: and (ν (0) ) 2 −1 ν 2 |r = 0 (2) with c and H positive O(1) constants. We may interpret these results by saying that concentration fluctuation produce collision hindering in the compressible 1D case, and to collision enhancement in the incompressible 3D case.
Notice the O(ǫ −1/2 ) behavior of ν 2 |r = 0 (2) cl in place of the expected O(ǫ −1 ) at second order. It is possible to show, however, that the O(ǫ −1/2 ) contributions from ρ cl and ρ B cancel, so that the total correction to the collision velocity variance is O(ǫ −1 ).

Conclusion
We have derived a perturbative approach for the two-particle statistics of a randomly advected inertial particle suspension, that is valid both in regimes of small Peclet number P e (strong molecular diffusion), and large Stokes number S (corresponding to high inertia). In both cases, one expands around a lowest order of independent Brownian particles. This is natural in the small P e regime; for large S [more precisely, for large ǫ, with ǫ defined in Eq. (6)], the mechanism is more subtle and is due to the fact that particle trajectories evolve on a characteristic scale ǫ 1/2 r v , much larger than the correlation length r v of the field u. The particle pair trajectories are dominated therefore by large separation, uncorrelated u contributions [22]. The perturbative approach we have derived provides an analytical description of how the Brownian particle limit of [22] is achieved, that goes beyond the qualitative considerations in [13]. The concentration fluctuation amplitudē n −2 (n −n) 2 , withn the mean concentration, appears to be O(ǫ −1/2 ) for compressible, and O(ǫ −1 ), for incompressible flows. The correlation length of the fluctuations is ∼ ǫ 1/2 r v , in contrast with the ǫ ≪ 1 case, in which, a power law at separations r < r v would occur.
The perturbative approach allows to identify a concentration fluctuation contribution to the statistics for the collision velocity ν, in those terms in the expansion for the joint PDF ρ(ν, r = 0), that lead to deviations from the uniform fluctuation-free regime in the separation PDF ρ(r). This goes beyond the observation that the expected relative velocity should decrease at small values of the separation r.
It has been noted in [29] that for ǫ ≪ 1, different collision velocities originate from par-ticle "jumps" starting at different initial separations. The collision velocity distribution is thus affected by the spatial structure of the clusters, which produce an effective collision hindering.
In the present large ǫ regime, no such direct physical association between clustering and collision dynamics exists. From an analysis of the collision velocity variance, we see that collision hindering occurs in 1D, while enhancement occurs in 3D if the flow is incompressible. The result in 1D is not unexpected, as clusters are in this case the result of particles slowing down relative to one another as they get closer. The role of incompressibility and 3D in leading to collision enhancement is less clear.
Additional information on the velocity statistics is contained in the velocity structure of the concentration fluctuations. From Eq. (28), we see that, for D < 3, the integral ρ(r) = d D νρ( ν, r) is dominated by small relative velocities, while for D = 3 all velocities contribute equally, down to the velocity scale ǫ −1/6 σ v , at which, particle pairs with r < r v cease to behave ballistically. (The Brownian regime at large ǫ is correlated in time at the scale of the Stokes time τ S , and one expects particle pairs with a relative velocity that is not too small, to behave ballistically at separation r < r v ).
At least in 1D, the predictions of perturbation theory begin to be valid only for rather large values of ǫ. The fact is that clustering at ǫ 1, and residual concentration fluctuations at ǫ ≫ 1, are very different in nature. At small ǫ, in first approximation, the particle separation r, for r ≪ r v evolves as a diffusion process with diffusivity ∝ ǫr 2 . Clustering arises therefore from trapping at r = 0 of the particle pairs. For large ǫ, most particle pairs at r < r v behave instead ballistically. In this regime, concentration fluctuations and associate particle velocity modifications, can be seen as corrections to ballistic motion at scale r v . They are at most only the remnants of the trapping behaviors that dominate the pair dynamics at ǫ 1.