Turbulence-Induced Relative Velocity of Dust Particles III: The Probability Distribution

Turbulence-Induced Relative Velocity of Dust Particles III: The Probability Distribution

Liubin Pan Harvard-Smithsonian Center for Astrophysics, 60 Garden St., Cambridge, MA 02138; lpan@cfa.harvard.edu Paolo Padoan ICREA & ICC, University of Barcelona, Marti i Franquès 1, E-08028 Barcelona, Spain; ppadoan@icc.ub.edu John Scalo Department of Astronomy, University of Texas, Austin, TX 78712; parrot@astro.as.utexas.edu
Abstract

Motivated by its important role in the collisional growth of dust particles in protoplanetary disks, we investigate the probability distribution function (PDF) of the relative velocity of inertial particles suspended in turbulent flows. Using the simulation from our previous work, we compute the relative velocity PDF as a function of the friction timescales, and , of two particles of arbitrary sizes. The friction time of particles included in the simulation ranges from to , with and the Kolmogorov time and the Lagrangian correlation time of the flow, respectively. The relative velocity PDF is generically non-Gaussian, exhibiting fat tails. For a fixed value of , the PDF is the fattest for equal-size particles (), and becomes thinner at both and . Defining as the friction time ratio of the smaller particle to the larger one, we find that, at a given in , the PDF fatness first increases with the friction time, , of the larger particle, peaks at , and then decreases as increases further. For , the PDF shape becomes continuously thinner with increasing . The PDF is nearly Gaussian only if is sufficiently large (). These features are successfully explained by the Pan & Padoan model. Using our simulation data and some simplifying assumptions, we estimated the fractions of collisions resulting in sticking, bouncing, and fragmentation as a function of the dust size in protoplanerary disks, and argued that accounting for non-Gaussianity of the collision velocity may help further alleviate the bouncing barrier problem.

1. Introduction

This is the third paper of a series on turbulence-induced relative velocity of dust particles. The study is mainly motivated by the problem of dust particle growth and planetestimal formation in protoplanetary disks (e.g., Dullemond and Dominik 2005; Zsom et al. 2010, 2011; Birnstiel et al. 2011; Windmark et al. 2012a, 2012b; Garaud et al. 2013, Testi et al. 2014). In the first two papers of the series (Pan & Padoan 2013, Pan, Padoan & Scalo 2014; hereafter Paper I and Paper II, respectively), we conducted an extensive investigation of the root-mean-square or variance of the relative velocity of inertial particles in turbulent flows using both analytical and numerical approaches. In particular, we showed that the prediction of the Pan & Padoan (2010) model for the rms relative velocity is in satisfactory agreement with the simulation data, confirming the validity of its physical picture. In Paper I, we also explored the collision kernel for the case of equal-size particles, known as the monodisperse case, and analyzed the probability distribution function (PDF) of the turbulence-induced relative velocity as a function of the particle inertia. The main goal of the current paper is to study the relative velocity PDF in the bidisperse case for particles of arbitrarily different sizes.

The collision velocity of dust particles in protoplanetary disks plays an important role in coagulations models of particle growth, which is crucial for understanding the challenging problem of planetesimal formation. The collision velocity determines not only the collision rate, but also the collision outcome. As the size grows, the particles become less sticky, and, depending on the particle properties and the collision energy, the collisions may lead to bouncing or fragmentation (Blum & Wurm 2008, Güttler et al. 2010). Due to the stochastic nature of turbulence-induced relative velocity, the collision outcomes for dust particles with exactly the same sizes and properties may be completely different. It is thus not sufficient to use an average or mean collision velocity to predict the size evolution of dust particles. Instead, an accurate prediction would require the PDF of the collision velocity to evaluate the fractions of collisions resulting in sticking, bouncing or fragmentation.

The effect of the collision velocity PDF on particle coagulation has been considered by several recent studies (Windmark et al. 2012b, Garaud et al. 2013). These studies assumed a Maxwellian distribution for the 3D amplitude of the collision velocity, or equivalently a Gaussian distribution for each component111For simplicity, we will not distinguish “Maxwellian” for the 3D amplitude and “Gaussian” for one component, and refer to both as Gaussian.. Implementing the distribution into coagulation models, they found important differences in the prediction of the dust size evolution. In particular, they showed that accounting for the distribution of the collision velocity can soften the fragmentation barrier and overcome the bouncing barrier (Windmark et al. 2012b). The maximum particle size that can be reached by collisional growth appears to be significantly larger than the case using only a mean collision velocity. This alleviates the so-called meter-size barrier for planetesimal formation via collisional growth of dust particles.

The assumption of a Gaussian PDF made in the above-mentioned papers is not justified for the collision velocity induced by turbulent motions. For equal-size particles, turbulence-induced relative velocity has been found to be highly non-Gaussian, exhibiting extremely fat tails (Sundaram & Collins 1997, Wang et al. 2000, Cate et al. 2004, Gustavsson et al. 2008, Bec et al. 2010, 2011, de Jong 2010, Gustavsson et al. 2012, Gustavsson & Mehlig 2011, 2014, Lanotte et al. 2011, Gualtieri et al. 2012, Hubbard 2012, Salazar & Collins 2012). The earlier works showed that the PDF for equal-size particles may be fit by exponential or stretched exponential distributions (e.g., Sundaram & Collins 1997, Wang et al. 2000, Cate et al. 2004). A stretched exponent PDF with an index of -4/3 was predicted theoretically for inertial-range particles under the assumption of exactly Gaussian flow velocity and Kolmogorov scaling (Gustavsson et al. 2008, Paper I). Lanotte et al. (2011) found that the relative velocity for small particles of equal size shows a power-law PDF in the limit of small particle distance. A power-law PDF was predicted by recent theoretical models based on a Gaussian smooth velocity field with rapid temporal decorrelation (e.g., Gustavsson & Mehlig 2011, 2014).

In Paper I, we conducted a systematic study for the PDF of equal-size particles as a function of the particle inertia, and showed that the behavior of the PDF shape was successfully explained by the model of Pan & Padoan (2010, hereafter PP10). In this work, we further analyze the PDF for particles of different sizes, and show that non-Gaussianity is a generic feature of turbulence-induced relative velocity, which should be incorporated into coagulation models for an accurate prediction of the dust size evolution.

Despite extensive studies on the variance or rms of the relative velocity of different-size particles (e.g., Völk et al. 1980, Markiewicz et al. 1991, Zhou et al. 2001, Cuzzi & Hogan 2003, Ormel & Cuzzi 2007, Zaichik et al. 2006, 2008, Zaichik & Alipchenkov 2009), the PDF in the general bidisperse case has received little investigation (see, however, Johansen et al. 2007). A recent work by Hubbard (2013) explored the PDF for particles of different sizes using a synthetic or “model” flow velocity filed. However, the results obtained from such an approach are clearly not accurate because the artificial velocity field does not correctly account for the non-Gaussianity or intermittency of turbulent flows, which does leave an imprint on the relative velocity of inertial particles (Paper I). The commonly-adopted models for the rms relative velocity (Völk et al. 1980, Ormel & Cuzzi 2007) in the astronomy community are based on the responses of particles to turbulent eddies in Fourier space, and do not provide adequate insight for understanding or predicting the distribution of the relative velocity. A physical weakness of these models was pointed out and discussed in PP10 and Paper I. On the other hand, the PP10 model is based on the statistics of the flow velocity structures in real space, and we have shown in Paper I that the physical picture revealed by the model offers a satisfactory explanation for the PDF behavior in the case of equal-size particles. In Paper II, we further developed the model, and established physical connections between the particle relative velocity and the temporal and spatial velocity structures of the carrier flow. In this paper, we will continue to use the picture of PP10 to interpret the trend of the relative velocity PDF shape as a function of the particle friction times in the bisdisperse case.

We use the same simulation as in Papers I and II. It is carried out using the Pencil code (Brandenburg & Dobler 2002) with a periodic box. The simulated flow is driven and maintained by a large-scale force, , which is generated in Fourier space using all modes with wave length half box size. Each mode gives an independent contribution to the driving force. At each time step of the simulation, the direction of each mode is random, and the amplitudes of each mode in three spatial directions are independently drawn from a Gaussian distribution. This conventional method of driving produces a turbulent flow with the broadest inertial range at a given resolution and with the maximum degree of statistical isotropy. At steady state, the 1D rms flow velocity is , in units of the sound speed, corresponding to a (3D) rms Mach number of 0.085. The Taylor Reynolds number of the flow is estimated to be , and the regular Reynolds number is . The integral length, , of the flow is 1/6 of the box size, and the Komlogorov scale is times the size of a computational cell, so that . The Kolmogorov scale, , is computed from its definition, , using the viscosity, , and the average dissipation rate, , in the simulated flow. Using tracer particles, we estimated the Lagrangian correlation time of the flow, , which is relevant for the particle dynamics than the large-eddy turnover time, (). We found , where the Kolmogorov timescale, (), corresponds to the turnover time of the smallest eddies. The large eddy turnover time, , is estimated to . The Kolmogorov velocity scale, (), is related to the 1D rms velocity of the flow, , by .

The simulated flow evolved 14 species of inertial particles of different sizes, each containing 33.6 million particles. We integrated the trajectory of each particle by solving the momentum equation for the particle velocity,

 dvdt=1τp[u(X(t),t)−v] (1)

where is the flow velocity at the particle position, , at time . The friction timescale, , of the smallest particles in the simulation is , while that of the largest particles is Defining a Stokes number as , this range corresponds to . Spanning 4 orders of magnitude, this range of covers all the length scales of interest in the simulated turbulent flow. When integrating the particle trajectories, we interpolated the flow velocity inside computational cells using the triangular-shaped-cloud (TSC) method (Johansen and Youdin 2007). Our simulation run lasted (or ). At the end of the run, all the statistical measures reached a quasi steady-state and the dynamics of all particles was relaxed. For our statistical analysis, we use three well-separated snapshots toward the end of the run.

Following Paper II, we name two nearby particles under consideration as particles (1) and (2), and denote their friction times as and . The particle position, their velocity and the flow velocity at the particle position are denoted by , , and , respectively. We set the time at which the relative velocity is measured to be time zero, i.e., . The particle separation and the relative velocity at are denoted by () and (). The separation of the two particles as a function of time is , and . The spatial flow velocity difference across the particle distance is denoted by . We denote as the temporal flow velocity difference at a time lag of along the trajectories of particles (1) and (2), respectively. We also consider the relative velocity, , between a particle and the flow element at the particle location. For convenience, we denote the friction times (Stokes numbers) of the smaller and larger particles as and ( and ), where the subscripts “l” (or “”) and “h” (or “”) stand for low and high, respectively. We define a friction time or Stokes ratio as . By definition, . It is also convenient to define the ratio, , of the friction time to the Lagrangian correlation time.

In §2, we review the physical picture of the PP10 model for turbulence-induced relative velocity in the bidisperse case. In §3, we discuss several interesting limits, and, in particular, we present our simulation result for the PDF of the particle-flow relative velocity, . In §4, we compute the particle relative velocity PDF for all friction time (or Stokes number) pairs available in our simulation, and interpret the results using the physical picture of PP10. The implications of our results for dust particle collisions in protoplanetary disks are discussed in §5. We summarize the main conclusions of this study in §6.

2. The Physical Picture of the PP10 Model

In Papers I and II, we introduced and further developed the formulation of Pan & Padoan (2010, PP10) for the relative velocity of inertial particles induced by turbulent motions. We showed in Paper II that the model prediction for the root-mean-square relative velocity agrees well with our simulation results for particles of arbitrarily different sizes. Here, we briefly review the physical picture revealed by the PP10 model, and refer the reader to Papers I and II for further details. The physical picture will be used to qualitatively interpret the trend of the relative velocity PDF measured from our simulation data.

In the formulation of PP10, the relative velocity, , for particles of different sizes can be approximately written as the sum of two terms, , where and are named the generalized acceleration and shear contributions, respectively. The two contributions reduce to the corresponding terms in the formula of Saffman & Turner (1956) in the limit of small particles with much smaller than the Kolomogorov time of the flow, . The generalized acceleration term reflects different responses of particles of different sizes to the flow velocity. In Paper II, we showed that is physically associated with the temporal flow velocity difference, , on individual particle trajectories (see Fig. 1 of Paper II for an illustration). Based on a calculation of the contribution of the acceleration term to the rms relative velocity, we established in Paper II an approximate expression for ,

 wa≃1−f[(1+Ωl)]1/2ΔuT(τp,h) (2)

where is the friction time ratio and . The factor, , indicates that depends on the friction time difference, and it vanishes for equal-size particles. The equation connects to at a time lag, , equal to the friction time, , of the larger particle. The probability distribution of the acceleration contribution is thus related to the PDF of . In eq. (2), we ignored the possible difference in the statistics of and along the trajectories of two different particles, and here should be viewed as the average of and (Paper II). In the limit , the particle trajectories are close to the Lagrangian trajectories, and we have , where is the Lagrangian temporal flow velocity difference and is the local flow acceleration. Therefore, eq. (2) gives for small particles with , consistent with the Saffman-Turner formulation (see also Weidenschilling 1980).

Physically, the generalized shear term represents the particles’ memory of the spatial flow velocity difference, , across the distance, , of the two particles at given times in the past (see Fig. 1 of Paper I and Fig. 2 of Paper II for illustrations for the mono- and bi- disperse cases, respectively). In the case of equal-size particles, the generalized shear term is the only contribution to the relative velocity. As discussed in Paper I, the relative velocity of a particle pair of equal size can be approximated by , where , named the primary distance, is defined as the separation of the particle pair at a friction timescale ago, i.e., . This distance is of particular importance because it is the largest distance before the particle memory cutoff takes effect at . The estimate of depends on how the particle pair separates backward in time. Due to the particle inertia, an initial ballistic separation is expected, and, assuming that the duration of the ballistic behavior is , the primary distance is estimated as , where is the particle distance at the current time. Note that the primary distance of a particle pair depends on their relative velocity, , at , which has important implications for the relative velocity PDF (see Paper I and discussions in §3.2). The timescale is the correlation time of the flow velocity structures at the scale , which is essentially the turnover time of turbulent eddies of size . The factor is due to the cutoff by the flow “memory”, which in effect causes a reduction in the time range of the particle memory that can contribute to the relative velocity when (see Paper I).

Similar to the monodisperse case, Paper II proposed an approximate equation for for particles of different sizes. Based on eq. (37) in Paper II, we have222Here, we replaced and in eq. (37) of Paper II by and , respectively. In Paper II, denotes the average primary distance over all particle pairs at a given distance, , and the same average is implied in the timescale . Therefore, eq. (37) of Paper II is an approximation for the overall average shear contribution. Here, in order to understand the PDF of the relative velocity, we are interested in individual particle pairs, and thus used the primary distance and the timescale for a single particle pair to replace and . With this replacement, in eq. (3) corresponds to the generalized shear contribution for each individual pair.,

 ws≃Δu(rp)(tptp+τp,h)1/2. (3)

As discussed in Paper II, the primary distance for the bidisperse case is mainly determined by the smaller particle, and a rough estimate is , where is the 3D amplitude of the relative velocity of the particle pair including both the acceleration and shear contributions333Since the rate of the ballistic backward separation in the bidisperse case has a contribution from the acceleration term, the distribution of is complicated. Rigorously, for particles of different sizes, it depends not only on the sptatial flow velocity structures but also on the statistics of . On the other hand, the distribution of is simpler, as it is completely controlled by (see eq. (2)). Since depends on the amplitude of , eq. (3) is an implicit equation. As in the monodisperse case, the last factor in eq. (3) represents the effect that the flow memory cutoff may reduce the range of the larger particle’s memory (around the memory time of the smaller particle) that can contribute to the relative velocity. The direction of with respect to is stochastic due to the “turbulent” separation of the particle pair backward in time. In the small particle limit , we have , and is approximately equal to , meaning that the shear term follows the spatial flow velocity difference across the particle distance (Saffman & Turner 1956). For below the Kolmorgorov scale, the flow velocity is smooth and is linear with , i.e., , with the local velocity gradients. The distribution of in the small particle limit is thus related to the PDF of or equivalently the PDF of the velocity gradients.

There are several interesting qualitative differences between the generalized acceleration and shear terms. A fundamental difference is that the acceleration term is determined by the flow velocity along individual particle trajectories, while the shear term depends on the trajectories of the two particles relative to each other, as seen from the dependence of on . The shear contribution is thus more complicated. For example, for small particles, has a dependence, and it also depends on whether the particles are approaching or separating from each other. In contrast, the generalized acceleration term is independent of the particle distance or the relative motion of the two particles. Also, the direction of in eq. (2) is random with respect to , and thus provides equal contributions to the radial () and tangential () components of the relative velocity along and perpendicular to , respectively. On the other hand, the shear contributions to the two components may differ for small particles because the longitudinal () and transverse () velocity difference in a turbulent flow across a small distance are nonequal (Papers I and II).

3. Three Limiting Cases

In this section, we consider three special cases of interest, i.e., the particle-flow relative velocity, the relative velocity between equal-size particles, and a case where one of the particles is extremely large with . These limiting cases are helpful for the understanding of the general behavior of the relative velocity PDF in the bidisperse case.

3.1. The Particle-flow Relative Velocity

The relative velocity, , between an inertial particle and the instantaneous local flow velocity is a special bidisperse case with one of the particles, say particle (2), being a tracer particle with , and at the same position as particle (1). We refer to this limiting case with as Limit I. In this limit, only the generalized acceleration term contributes444 The generalized shear term is determined by the flow velocity differences the two particles saw within their friction times. Due to the zero memory time of the tracer particle, the shear term depends on the flow velocity difference across the particle distance at time zero. It thus vanishes if the tracer particle is at the same position as the inertial particle at time zero. This can also be also be seen from eq. (3), as the primary distance in this case., and it thus gives useful information for the distribution of for particles of different sizes. As discussed in Paper II, for a particle with a friction time , can be approximately estimated by the temporal flow velocity difference, , on the particle trajectory at . This can be seen from eq. (2) with , , and .

Using our simulation data, we computed the PDF of . The flow velocity at the position of each particle was obtained by the same interpolation method (TSC) used in the simulation. The computed relative velocity is at zero particle-flow distance. In the left panel of Fig. 1, we plot the PDF, , of one component, , of the particle-flow relative velocity for 6 particles with , 0.39, 1.55, 6.21, 24.9, and . Note that here is not the 3D amplitude of . Each curve shows the PDF averaged over the three components of along the base directions of the simulation grid. The PDF of is symmetric for all particles, as expected from isotropy for a 1-point statistical quantity. The width or rms of the PDF increases with , corresponding to the increase of with the time lag (see Paper II). The dotted black line is the Gaussian fit to the case.

To see the PDF shape more clearly, the right panel of Fig. 1 plots the normalized PDF, , with unit variance, where , with the rms of . For clarity, we briefly introduce our terminology for the description of the PDF shape. Following Paper I, we use “fat” or “thin” to specifically describe the shape of the PDF, while the extension or width of the PDF, corresponding to the rms, will be described as “broad” or “narrow”. Conventionally, the fatness of a PDF is used as a measure of the departure from a Gaussian distribution, and in particular, it refers to higher probabilities than a Gaussian distribution at the tail part of the PDF. The fatness of a PDF may be quantified, e.g., by its kurtosis, and a fatter PDF has a larger kurtosis. Typically, a larger kurtosis corresponds to a fatter tail part, and at the same time a sharper or thinner shape at the inner part of the PDF. Despite a sharper inner part, the convention is to describe the PDF with a larger kurtosis as fatter based on the tail part. In the right panel of Fig. 1, we see that, for small particles, the PDF is highly non-Gaussian with fat tails. For example, for , the kurtosis of the PDF is . As increases, the PDF shape becomes thinner and thinner. The bottom line in the right panel corresponds to the largest particles () in our simulation, and the PDF for these particles is Gaussian (the black dotted line). In fact, becomes close to Gaussian at (or ).

From the PP10 picture, the PDF of is controlled by the distribution of , which, however, is unknown. For very small particles with , one may approximate by the temporal velocity difference, , along Lagrangian trajectories of tracers. On the other hand, due to their slow motions, the flow velocity seen by large particles with may be better described as Eulerian, and is likely close to the Eulerian temporal flow velocity difference at fixed points. In Appendix A, we compute the PDFs of and in our simulated flow and show that both PDFs are non-Gaussian at small time lags . The PDF tails of both and have a thinning trend with increasing and reach Gaussianity at (see, e.g., Mordant et al. 2001, Chevillard et al. 2005). The non-Gaussianity of and and their trend with is a well-known phenomenon in turbulent flows, usually referred to as intermittency (e.g., Frisch 1995). The behavior of the PDFs of and explains the non-Gaussianity of and the thinning of with increasing . In this sense, the particle-flow relative velocity “inherits” non-Gaussianity from the turbulent flow. From eq. (2), the PDF trend of suggests that the distribution of the acceleration contribution becomes thinner and thinner with increasing .

A comparison of the normalized PDF of with that of (Appendix A) at a time lag shows that, for , the tails of are thinner than . The physical reason is that inertial particles do not exactly follow the fluid elements. In particular, inertial particles tend to be pushed out of regions of high vorticity by a centrifugal force arising from the rotation (e.g., Pan et al. 2011). Therefore, they do not experience the flow velocity in strong vortex tubes at small scales, which are intense structures responsible for the high non-Gaussianity (or intermittency) of in turbulent flows (Bec et al. 2006). This effect contributes to make the PDF of thinner. For particles with , the PDF shape for and are close to each other, and this is because both the dynamics of these larger particles and the temporal flow velocity structures at larger timescales are less sensitive to the intense, coherent vortices at small scales.

For all particles, the innermost part of the PDF of shows a smooth Gaussian-like shape, suggesting that the PDF of trajectory temporal velocity difference, , also takes a similar shape in the inner part. This is supported by the observation in Appendix A that the central parts of the PDFs of and have a more-or-less Gaussian shape. The fat non-Gaussian PDF tails of can be fit by stretched exponential functions defined as,

 Pse(x)=α2βΓ(1/α)exp[−(|x|β)α], (4)

where is the Gamma function (Paper I). The variance of is given by . The parameter controls the fatness of the PDF, and a smaller corresponds to fatter tails. As expected, the value for the tails of increases with increasing or . For the eight lines shown in the right panel of Fig. 1, the best-fit for the PDF tails are 0.7, 0.78, 0.95, 1.23, 1.33, 1.7, 1.95 and 2 for , 0.39, 1.55, 6.21, 12.4, 24.9, 99.4, and 795, respectively.

The black dashed line in Fig. 1 is the stretched exponential function with that best fits the PDF tail of particles. The friction time of these particles is close to the Lagrangian correlation time of the flow. In Paper I, we found that this stretched exponential can also well fit the relative velocity PDF of equal-size particles with and (see also §3.2). This similarity in the shape of the PDF tails for and the monodisperse relative velocity of particles is of particular interest. Its implication will be discussed in §4.2.1.

Finally, we point our that, if the flow velocity were exactly Gaussian, then the temporal flow velocity difference would be Gaussian, and our model picture indicates that the PDF of the particle-flow relative velocity would be Gaussian for particles of any sizes. In that case, the relative velocity PDF between different particles approaches Gaussian once the generalized acceleration contribution dominates. This suggests that using a Gaussian model flow velocity (see, e.g., Hubbard 2013) to study the particle relative velocity PDF would underestimate the degree of non-Gaussianity.

3.2. Equal-size Particles

The relative velocity PDF of equal-size particles has been discussed in details in Paper I. For convenience, we refer to this monodisperse case with as Limit II. In Paper I, we found that the relative velocity PDF of identical particles is extremely non-Gaussian. We identified two sources of non-Gaussianity, i.e., the imprint of the non-Gaussianity (or intermittency) of the turbulent flow and an intrinsic contribution arising from the particle dynamics. Our simulation data showed that the fatness of the monodisperse PDF first increases with , reaches a peak at , and then continuously decreases as increases above 1. In the following, we briefly review the explanation for this behavior based on the PP10 picture.

For equal-size particles, only the shear term contributes, and the relative velocity, , can be estimated by eq. (3) for . For a particle pair with , the factor in eq. (3) is of oder unity, and, qualitatively, can be approximated by the spatial flow velocity difference, , across the primary distance, (Paper I). At a given , the primary distance, , of particle pairs with is close to . Therefore, the inner part of the PDF with follows at a distance of . The distribution of at small scales is already non-Gaussian (e.g., Frisch 1995), and it leaves a non-Gaussian imprint on the PDF of . Particle pairs lying at higher PDF tails separate faster backward in time and have larger , so their relative velocity “remembers” at larger scales . Here, denotes a component of at a separation . Since the PDF of is wider at larger , this leads to an enhancement at the tail part of the PDF. The effect, named as a self-amplification in Paper I, corresponds to the sling effect or caustic formation that occurs in flow regions with large velocity gradients, (Falkovich et al. 2002, Wilkinson & Mehlig 2005, Wilkinson et al. 2006, Falkovich & Pumir 2007, Bewley et al. 2013). As increases, the self-amplification proceeds from the tails toward smaller , while the innermost part of the PDF remains unchanged at first. As a result, the overall shape of the PDF becomes fatter.

The fattening trend ends at . As increases above 1, the amplification moves deeper toward the central part, and the PDF at the same value of samples at larger scales. Since the PDF shape of is progressively thinner with increasing (see Appendix B of Paper I), this causes the overall fatness of the PDF of to decrease continuously at . The PDF shape is thus the fattest at .

For particles with close to , the PDF tails of were found to be well described by a 4/3 stretched exponential distribution (see eq. (4)). Such a distribution were obtained in Paper I from a phenomenological argument based on the PP10 picture. Assuming exactly Gaussian flow velocity and Kolmogorov scaling for with the Kolmogorov constant, the probability of finding a particle pair with a relative velocity of is estimated as , which is a 4/3 stretched exponential at (see also Gustavsson et al. 2008). Using the same argument, one can show that non-Gaussianity is an intrinsic feature of the particle dynamics: Even if the flow velocity were exactly Gaussian, the particle relative velocity would be non-Gaussian for any scaling behavior of with (unless it is constant with as in the case of a white noise).

Finally, for , the flow correlation or memory time is typically shorter than the particle memory, and the last factor in eq. (3) reduces the PDF width by a factor of . In the limit , the particle motions are similar to Brownian motion, because even the largest eddies in the flow act as random noise when viewed over the friction time of these large particles. Therefore, the relative velocity PDF is expected to approach Gaussian as . A Gaussian shape is almost reached for the largest particles with in the simulation of Paper I.

The comparison between Limits I and II at equal is useful to understand the PDF trend in the general bidisperse case (see §4). For a small particle with , the PDFs of its relative velocity with respect to the local flow velocity and to a nearby particle of the same size are determined by the distributions of the flow acceleration, , and the velocity gradients, , respectively (see §2). In Appendix B, we compute the normalized PDFs of and in our simulated flow, and show that the distribution of is generally fatter, consistent with previous studies (e.g., Ishihara et al. 2007). This suggests that, in the limit, the PDF of has a fatter shape than the relative velocity of equal-size particles. As increases toward , the PDF for equal-size particles fattens significantly due to the tail amplification (see discussion above), while the PDF for the particle-flow relative velocity becomes thinner continuously (see Fig. 1). The opposite trends of the two PDFs for have an important consequence. The PDF for equal-size particles is found to be already fatter than the distribution of at the smallest Stokes number, , in our simulation, and, at , the former has a much stronger non-Gaussianity than the latter. For above , the PDF fatness for both cases decreases with increasing , and we will show in §4.1.1 that the relative velocity PDF of equal-size particles remains to be significantly fatter than at all . Thus, for all particle species in our simulation, the PDF for Limit II has a higher degree of non-Gaussianity than for Limit I.

3.3. The τp2→∞ Limit and the One-particle Velocity

If one of the two particles, say particle (2), is extremely large with , its velocity is tiny, and the relative velocity with respect to the smaller particle (1) is essentially the one-particle velocity, , of particle (1). We refer to this case as Limit III, corresponding to . Paper I found that the PDF of the one-particle velocity, , for particles of any size is Gaussian. As discussed in Paper I, this is easy to see in the small and large particle limits. Small particles closely follow the local flow velocity, and the PDF of is approximately given by the 1-point PDF of the flow velocity, , which is known to be approximately Gaussian. On the other hand, the motions of large particles with are similar to Brownian motion, because, for these particles, turbulent eddies of all sizes act as random noise when viewed over a timescale of . Therefore, the one-particle velocity PDF in the large particle limit is also Gaussian. The Gaussianity of the one-particle velocity for particles of any size suggests that, for any given , the bidisperse relative velocity PDF approaches Gaussian as .

4. The PDF of the Particle Relative Velocity

We analyze the probability distribution function of the particle relative velocity for all Stokes number pairs in our simulation. To compute the relative velocity statistics, we search pairs of particles from all species at given distances, . We will consider , , and . It would be useful to measure the PDF at even smaller because the size of dust particles in protoplanetary disks is much smaller than the Kolmogorov scale of protoplanetary turbulence. However, the number of particle pairs available in the simulation decreases with decreasing , and, at , the number is too limited to provide sufficient statistics. We will thus show simulation results for , and discuss the dependence (or convergence) of the PDF in this range.

For each particle (1), we locate particles (2) in a distance shell . The shell thickness, , is set to , , and for , , and , respectively. For with , the number of particle pairs for any two species is , about sufficient for the PDF measurement. We examine the PDFs of the radial and tangential components, as well as the 3D amplitude, of the relative velocity. The 3D amplitude is of practical interest as it determines the total collision energy and hence the collision outcome. The PDFs of the radial and tangential components, on the other hand, are helpful for a thorough understanding of the underlying physics. In particular, the features of the radial and tangential PDFs provide very useful tests for the PP10 picture. To obtain the radial and tangential components, we set up a local coordinate system for each selected pair. In terms of the grid base vectors, , , and , the local system is taken to be , and . The Euler angles are defined as , , , and , with , and the components of in the grid coordinate. We decompose as , , and . The tangential PDF, , is set to the average of and of . The analysis is the same as in Paper I for equal-size particles. In §4.1, we fix the Stokes number, , of particles (1), and examine the trend of the relative velocity PDF as a function of . In §4.2, we will show the PDFs at fixed Stokes number ratios, . We will pay particular attention to the non-Gaussianity of the PDF, and discuss how the PDF shape changes with the Stokes numbers. As a reminder, in our terminology, a “fat” (or “thin” ) PDF shape refers to higher (or lower) probabilities specifically at the tail part of the PDF.

4.1. The PDF at Fixed St1

If one of the Stokes numbers is fixed at , the three interesting limits, i.e., , , and , discussed in §3 are useful to confine the trend of the PDF shape as a function of . These cases were named as limits I, II, and III, respectively. The shape of the PDF for is expected to lie in between Limit I and Limit II, while corresponds to a range between Limits II and III.

4.1.1 The PDFs of the radial and tangential relative speeds

In Fig. 2, we plot the PDFs, , of the radial relative velocity for fixed at . The relative velocity is measured at . The left and right panels show results for and , respectively. In both panels, the solid black line is the PDF for identical particles with , corresponding to Limit II. The monodisperse PDF at is highly non-Gaussian, and the physical origin of the non-Gaussianity has been explained in §3.2 using the PP10 picture. This PDF is also negatively skewed. The left and right wings of the PDF correspond to approaching and separating particle pairs. The asymmetry of the two wings indicates faster relative velocity for approaching pairs. As discussed in Paper I, there are two reasons for this asymmetry. First, the radial relative speed PDF for small equal-size particles inherits a negative skewness from the PDF of the longitudinal flow velocity difference or gradients of the flow (see Appendix B and Paper I). Second, approaching particle pairs have a larger separation in the near past than separating pairs. This suggests that, for small particles, the primary distance, , of approaching pairs is larger than separating ones, which enhances the asymmetry as increases at small . As keeps increasing, the particle distance () at a friction time ago is less dependent of the particle separation in the near past, and the asymmetry decreases (see Figs. 4 and 5).

The black dashed line in the left panel plots the PDF of the radial component, , of the relative velocity between particles and the flow element, corresponding to Limit I with and . For a consistent comparison with the particle-particle case, here the PDF of is measured at the same distance, . It is thus not the same as the particle-flow relative velocity PDFs shown in Fig. 1 at zero distance. To compute at a finite distance, we used the TSC interpolation to obtain the flow velocities at a separation of from the position of each particle in the three base directions of the simulation grid. We then averaged the PDFs of measured from the three directions. Unlike the particle-flow relative velocity at zero distance, here the PDF of at for is asymmetric. This is because, at a finite distance, the generalized shear term has a nonzero contribution (see eq. 3) to the relative velocity, which gives rise to a noticeable negative skewness in the PDF of for small particles with .

As expected, the bidisperse PDFs in the left panel of Fig. 2 lie in between Limit II (black solid line) and Limit I (dashed solid line). Since a negative skewness exists in both limits, the PDF is asymmetric for all . At the smallest (; the blue line) shown in the figure, the radial PDF, , approaches Limit I. As decreases below , the central part of the PDF widens, and the rms of the PDF increases (Paper II). This corresponds to the increase in the contribution of the generalized acceleration term, . As discussed in §2, is related to the temporal flow velocity difference, , on the particle trajectory. Applying eq. (2) to the Stokes number pairs in the left panel of Fig. 2, we have , where for . Clearly, as decreases, decreases and the contribution of increases. On the other hand, it can be seen from eq. (3) that the shear contribution decreases with decreasing (see Paper II).

The black solid line for Limit II has a fatter overall shape than the black dashed line for Limit II (see the discussion at the end of §3.2), suggesting that the distribution of the generalized acceleration contribution is thinner than the shear contribution. As increases with decreasing , the overall fatness of the relative velocity PDF decreases. Interestingly, despite the increase in the rms width, the probability at the far tails of the relative velocity PDF becomes smaller at smaller , indicating a lower probability of finding particle pairs that collide with extremely high velocity. In general, whether the probability at the far tails increases or decreases with decreasing is determined by the competition of two effects. First, if the PDF shape is given, an increase in the rms tends to give a larger probability at the far tails. On the other hand, if the rms is fixed, the thinning of the PDF shape would lead to a lower probability at high tails. It appears that the effect of thinning PDF shape wins in the case of and .

The right panel of Fig. 2 shows the radial relative speed PDF for . As increases, the PDF becomes wider, and its shape becomes thinner. This is again due to the increase in the contribution from the acceleration term. Using eq. (2) here, we have , where . As increases, both the time lag in and the factor increase, leading to the increase of the rms width of the PDF. The PDF shape of is expected to become thinner with increasing . This is based on the observation of the thinning trend of the Eulerian and Lagrangian temporal velocity differences with increasing time lag (see Appendix A). The argument was used earlier to explain the thinning trend of the particle-flow relative velocity with particle inertia (see §3.1). The thinning of the distribution of with increasing makes the PDF of the relative velocity thinner. At , the PDF is close to Gaussian, consistent with the expectation for Limit III (§2.3). The black dashed line in the figure is the Gaussian fit to the case. We note that the PDF becomes symmetric for .

Since the PDF shape becomes thinner for both above or below , the fatness of the PDF peaks at equal-size particles. The PDF shapes in both limits I and III are thinner than Limit II. A computation of the skewness of the radial PDFs shows that it decreases as the Stokes number difference increases. This is because the generalized acceleration is independent of the relative motions of the two particles and thus provides a symmetric contribution to the two wings of the radial relative velocity PDF (see discussion in §2).

In Fig. 3, we show the PDF of the tangential relative speed, at a distance of for . The conditioning on a negative radial relative speed () indicates that only particle pairs approaching each other are counted. Since only approaching particles may lead to collisions, it is of practical interest to separate them out from the pairs moving away from each other. Fig. 3 is plot in the same way as Fig. 2. The trend of as a function of is similar to that of . Unlike the radial component, the two wings of the tangential PDF are symmetric, as expected from statistical isotropy. We find that, for any , the left wing of coincides with that of . The coincidence in the monodisperse case was found in Paper I for all particles with in our simulation.555Because the longitudinal () and transverse () flow velocity differences in a turbulent flow are not equal, as , the coincidence of the radial and tangential PDFs for approaching particle pairs of equal size are not trivial. As discussed in Paper I, the equalization of and is essentially caused by the randomization of the direction of the primary separation, , with respect to , due to the deviation of the particle trajectories from the flow elements and the turbulent separation of particle pairs backward in time (see §3.2). For the bidisperse case, the generalized acceleration term gives equal contributions to and (§2), and it thus enhances the equalization of the radial and tangential PDF.

The central part of the monodipserse PDF is very sharp with a cusp-like shape666The formation of the sharp cusp in the monodisperse PDF can be explained by the picture described in §3.2. For small particles, the faster backward separation of particle pairs at higher PDF tails causes a self-amplification of the tails, while the innermost part remains unchanged. As the amplification proceeds toward the center of the PDF with increasing , the range of the unaffected central part is narrower, and the inner part would appear sharper relative to the outer parts. As increases further, the innermost part of the PDF would eventually be affected, and the cusp then shrinks and finally disappears for . (see also Fig. 4 for the case, where the cusp in the monodisperse PDF is even sharper). In the bidisperse case, the central cusp smoothes rapidly as the Stokes number difference increases. This is caused by the contribution of the acceleration term. As argued in §3.1, the PDF of the temporal trajectory velocity difference, , is expected to take a smooth, Gaussian-like shape in the innermost part. Then, since (eq. (2)), the bidisperse PDF would have a smoother central part, when the acceleration contribution increases.

Fig. 4 plots the PDF of the radial relative velocity for . The figure is plot in the same way as Fig. 2 for the radial PDF for . As found in Paper I, the fatness of the PDF shape for identical particles decreases with for , and the tails of the monodisperse PDF (the black solid lines) for are thinner than the case shown in Fig. 2 (see explanation in §3.2). At , the skewness for the monodisperse radial PDF disappears. An explanation for the recovery of the symmetry at large was discussed earlier in this section.

Similar to the case, the rms width of the PDF increases as decreases below or increases above , while the overall shape of the PDF becomes thinner. Interestingly, in the left panel, the high PDF tails almost coincide for all below . In particular, the far tails in the dashed black line for Limit I (i.e., the particle-flow relative velocity) happen to be close to those in the black solid line for equal-size particles. As mentioned earlier, two effects, i.e., the increase of the rms and the thinning of the overall PDF shape with decreasing , determine the trend of the probability at the far tails. The coincidence of the far tails suggests that, at , these two effects roughly cancel out. In the right panel, as increases above , the PDF keeps widening, and the PDF shape approaches Gaussian (the black dashed line) at .

Fig. 5 plots the PDF of the radial relative velocity for . In the left panel, we see that the PDF width increases with decreasing . The widening of the PDF width wins over the thinning trends of the overall PDF shape with decreasing , and thus the probability at the far tails keeps increasing with decreasing . In the right panel, the PDF becomes thinner as increases above , and finally approaches Gaussian, as expected for Limit III. The rms width of all PDFs shown in this panel is almost the same.

To summarize, for a given , the tail shape of the relative velocity PDF is the fattest for the monodisperse case with . As increase above or decreases below , the generalized acceleration contribution increases, and the overall PDF shape becomes thinner.

An important conclusion of theoretical interest is that the distribution of the generalized shear contribution, , is fatter than the acceleration term, . This is seen from a comparison of Limits I (black dashed lines) and II (black solid lines) in Figs. 2, 4 and 5, which are dominated by the generalized acceleration and shear terms, respectively. Since and are related to the temporal () and spatial () flow velocity differences, respectively, one may understand their distributions by considering the statistics of spatial and temporal velocity structures in turbulence. We first note that the spatial velocity structures in turbulent flows do not have a higher degree of non-Gaussianity than the temporal structures. In fact, the Lagrangian temporal structures are known to be more intermittent than the Eulerian spatial structures (see Appendix A). Therefore, the finding that the distribution of is fatter than cannot be interpreted by a simple comparison of the temporal and spatial intermittency. A key to understand the fatter distribution of is that and sample the temporal and spatial flow velocity differences in different ways. Unlike the PDF of , which roughly samples at a single time lag , is not controlled by the spatial flow velocity difference, , at a single length scale. In fact, different parts of the distribution of sample the flow velocity at different scales. For example, in the case of equal-size particles with a friction time of , the PDF at relative velocities below/above the rms value, i.e., at , depends on the flow velocity at scales , respectively. Considering that the PDF width of increases with increasing , this implies that the PDF shape of is significantly fatter than the PDF of at a single scale . This fattening effect is responsible for why the distribution of has a higher degree of non-Gaussianity than , even though the spatial flow velocity structures are not more intermittent than the temporal ones.

4.1.2 The PDF of the 3D amplitude

In Fig. 6, we show the PDF of the 3D amplitude, , of the relative velocity of approaching particles with . In both panels, is fixed at , and the left and right panels plot the PDFs for and , respectively. We normalized to its rms value, i.e., , where is the 3D rms relative velocity of approaching particle pairs with negative . As discussed in Paper I, is related to the particle collisional energy and its PDF plays an important role in determining the collision outcome.

In both panels of Fig. 6, the solid and dashed black curves show the PDFs of equal-size particles with at and , respectively. The black dotted lines correspond to the normalized PDF for the amplitude of a 3D Gaussian vector. As shown in Paper I, the monodipserse PDF is extremely non-Gaussian, with large probabilities distributed at both very small and large relative speeds. The large probabilities at small () correspond to the sharp cusps at the central parts of the radial and tangential PDFs (see Figs. 2 and 3). The solid and dashed color curves are bidisperse PDFs at and , respectively. As moves away from , the degree of non-Gaussianity in the PDF decreases, consistent with the results in §4.1.1. When the Stokes numbers differ by a factor of 2, there is a rapid decrease in the probability at the left tail. This is because, due to the acceleration contribution, the central parts of the radial and tangential PDFs become smooth and the sharp central cusps disappear (see §4.1.1). As the Stokes number difference increases, the right PDF tail of become thinner, and we also see that the peak of the PDF moves to the right toward , meaning that more probabilities are distributed around the rms relative velocity. In the left panel, as decreases to , the PDF shape approaches that of the particle-flow relative velocity (brown lines), which is still fatter than Gaussian (the black dotted line). In the right panel, the PDF of approaches Gaussian in the limit of large . At , the PDF coincides with the Gaussian distribution. The trend of the PDF shape as a function of is again consistent with the expectation that it lies in between the three limits, , and . Fig. 7 shows the simulation result for , which is similar to Fig. 6 for the case.

In the monodisperse case, the PDF converges with decreasing already at , for particles with (see black solid and dashed lines in Fig. 7). However, for smaller particles of equal-size with , the PDF shape has an dependence at (see Fig. 6). The convergence for the PDF of identical particles of small size is very challenging to reach, e.g., in the study of Lanotte et al. (2011), the convergence is barely reached at for particles with . In the bidisperse case, the convergence is easier to achieve because the contribution from the generalized acceleration term is independent (see §2 and Paper II). All the bidisperse PDFs for different particles shown in Figs. 6 and Fig. 7 already converge at . In fact, we find that, except for the two smallest particles with and in our simulation (see Fig. 14 in §4.2.4), the convergence of the PDF shape is reached at for all non-equal Stokes pairs.

4.2. The Relative Velocity PDF at Fixed Stokes Ratios

In this subsection, we examine the relative velocity PDF for particle pairs with fixed Stokes number ratios, , and show how the PDF changes with . As a reminder, and are the Stokes numbers of the smaller and larger particles, respectively. With a fixed , it is easier to compare the bidisperse PDF with that of equal-size particles () discussed in Paper I. At a given , there are two interesting limits, and . For , both particles become tracers, and the particle relative velocity PDF approaches the PDF, , of the flow velocity difference, , across the particle distance, . In the opposite limit , the relative velocity is essentially the 1-particle velocity of the smaller particle, and its PDF approaches Gaussian (§3.3). A particularly interesting result we find is that, if the friction time of the larger particle is close to the Lagrangian correlation time of the flow, the PDF tails can be approximately described by a 4/3 stretched exponential for any value of .

4.2.1 The PDFs of the radial and tangential relative speeds

In Fig. 8, we show the PDF, , of the radial relative velocity at for a Stokes ratio of . The left and right panels plot results for and , respectively. Each PDF is normalized to its value at the central peak (i.e., at ). In the left panel, the black short-dashed line is the PDF of the radial relative velocity between tracer particles (i.e., ) at , and the long-dashed line is the stretched exponential function with (see eq. (4)) that best fits the PDF tails for . The dashed line in the right panel is the Gaussian fit for and . The figure is plot in the same way as Fig. 10 of Paper I for the monodisperse case (i.e., ). The PDF width first increases with , corresponding to the increase of with in the generalized acceleration term (eq. (2)) and the increase of with in the shear term (eq. (3)). But for large and , the and factors in eqs. (2) and (3) take effect, leading to the decrease of the PDF width at (the right panel). The behavior of the width or rms of the PDF at fixed as a function of has been studied and explained in Paper II in the context of the PP10 picture. The left and right wings are asymmetric due to the shear contribution. The asymmetry first increases with as increases to 0.39, then decreases at larger . This is similar to the case of equal-size particles (; see Paper I) and consistent with the expected behavior of the shear contribution (see a physical discussion in §4.1.1). The PDF reaches symmetry at .

A comparison of the left panel of Fig. 8 with that of Fig. 10 in Paper I for identical particles reveals an interesting difference. As explained in §3.2 (see Paper I for details), in the monodisperse case, the innermost part of the PDF closely follows the flow velocity difference PDF (the short-dashed line for tracers in Fig. 10 of Paper I). On the other hand, the central part of the bidisperse PDF at any in the left panel of Fig. 8 is wider than the short-dashed line for tracers. This is because the generalized acceleration term, which is absent for equal-size particles, contributes to the central part of the bidisperse PDF. Also, at the same , the PDF shape for is thinner than for the case of equal-size particles with (see §4.1.1 and also §4.2.3 below).

In Fig. 9, we show the PDF of the tangential relative velocity, , for approaching particle pairs with . The figure is plot in the same way as Fig. 8 for the radial case. The two wings of the tangential PDF are symmetric for any . The tangential PDF as a function of shows a similar trend as the case of the radial relative speed (Fig. 8). We find again that the left wings of the radial PDF and the tangential PDF, , of approaching pairs coincide (see §4.1.1). In Appendix C, we consider the tangential PDF, , of separating particle pairs (with ), and compare it with . Similar to the asymmetric wings of the radial PDF, for small particles of similar sizes there is a difference between and for approaching and separating pairs.