Free Access
Volume 579, July 2015
Article Number A43
Number of page(s) 20
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/201425120
Published online 24 June 2015

Online material

Appendix A: A measure of particle clumping

Particle clumping is difficult to measure objectively. Previous authors have measured the peak particle density, but this is a very local measure that contains limited information about the particle behavior. In particular when small particles (τf ≤ 0.003) form filaments that are clearly visible in a spacetime diagram (Fig. 4), no single grid cell achieves a very high density. At the same time, large particles (τf> 3) have may have a few cells with very high density, but the density peaks are unstable and short lived.

Figure A.1 illustrates our approach. For each time step, we take the particle surface density Σsolid (top row). Particle overdensities are clearly visible on this plot. We then average Σsolid over 25 orbits (⟨ Σsolidt, second row). This serves to to smooth out overdensities that are unstable, short-lived, or have too much radial drift. We chose 25 orbits because it is half of the 50-orbit timescale for Z to increase. We considered using the radial speed to separate clumps that are drifting into the star. The appeal of radial speed is that it provides an instantaneous measure. However, radial speed does not capture the longevity of a particle clump. For example, as in a traffic jam, the pattern speed of the clump may be lower than the particle drift speed. At the same time, particles with low radial speed may produce short-lived clumps. Therefore, it is necessary to compare the location of a clump across time. Averaging Σsolid over time is a simple way to achieve this goal.

Finally, we sort the grid cells of ⟨ Σsolidt in order from highest to lowest density, and compute the cumulative distribution F of the sorted grid cells. If the particle distribution is uniform, then F forms a straight diagonal line G from (0, 0) to (1, 1). We use the difference FG as a large-scale measure of clumping. The Kolmogorov-Smirnov (Smirnov 1948; Wall & Jenkins 2003) test gives a standard way to compare F and G. It computes a p-value that measures probability that the data set that produced F follows the distribution Gwhere D is the maximum distance between F and G, and n is the number of independent observations. We caution the reader against interpreting these p-values as strict probabilities; we only use p as a convenient metric to measure clumping. With that in mind, we divide the p-values into four broad categories:

  • Clumping is unlikely:

    p ≥ 0.45

  • Clumping is somewhat likely:

    0.25 ≤ p< 0.45

  • Clumping is likely:

    0.10 ≤ p< 0.25

  • Clumping is very likely:

    p< 0.10.

Figure 8 shows the result of this analysis across our entire range of simulations. The figure is broadly consistent with the spacetime diagrams in Figs. 4 and 5.

Appendix A.1: Alternate measures of particle clumping

To confirm our results, we devised an alternate method to measure particle clumping, based as much as possible on very different principles. Our second method is to perform a chi-squared test on ρpzt, without sorting the grid cells. We use the reduced chi-squared statistic, (A.3)where k = n − 1 is the degrees of freedom, Oj is the jth observation (value of ρpzt on grid cell j), Ej = ⟨ ρp is the jth expected density for a uniform distribution, and σ2 is the expected variance between different measurements. If , that indicates that the difference between observations and the model is greater than can be explained by the variance σ2. This can be because the model is incorrect, or because the true variance between measurements is greater than σ2. We divide into four intervals:

  • Clumping is unlikely:

  • Clumping is somewhat likely:

  • Clumping is very likely:

That leaves σ2 as the only unspecified parameter. When the simulations start, σ2 is almost zero, giving unrealistically high values. As sedimentation proceeds, the value of σ2 increases. We opted for a σ2 that is typical for the last 25 orbits of the sedimentation phase. Because is sensitive to σ2, there is necessarily some amount of subjectivity in σ2. In the end, we chose σ2 so that the τf = 0.003 simulations wold give the same result as in Fig. 8. Therefore, the true test is whether the method remains consistent with the KS method across all the other particle sizes.

Figure A.2 shows the final results for the method. The similarity between this figure and Fig. 8 is striking. Considering how different the two techniques are, the agreement between Figs. 8 and A.2 indicates that our core results are robust.

thumbnail Fig. A.1

Particle density distribution for three particle sizes. Particle size is measured in the stopping time τf, where τf = 0.001, 0.003 and 3 correspond to the suspension regime, streaming regime and radial drift regime respectively. Row a) is the Surface density of solids Σsolid as a function of the radial coordinate x, evaluated at two points in time 25 orbits apart (first blue, then green). Row b) averages Σsolid over the 25-orbit interval ⟨ Σsolidt. Row c) shows ⟨ Σsolidt with the grid cells sorted from highest to lowest density. A steep curve indicates stable particle clumps. Row d) shows the cumulative distribution of ⟨ Σsolidt in row c) (solid red curve) along with the cumulative distribution of a perfectly uniform particle distribution (dashed blue line). The distance between these two curves serves as a measure of particle clumping. In these plots the particle concentrations Z = ⟨ Σsolid ⟩ / ⟨ Σtotal are (left to right) 0.095, 0.095 and 0.05.

Open with DEXTER

thumbnail Fig. A.2

As in Fig. 8, but using the method to estimate the likelihood of clumping instead of the KS-derived method. The figure marks the region of the particle size vs. concentration phase space where the streaming instability is active. Particle size is measured in the stopping time τf and the particle concentration is Z = Σsolid/ Σtotal. The colors mark the probability that particle clumps can form, where red is “unlikely”, magenta is “somewhat likely”, blue is “likely”, and green is “very likely”. When different simulations give different results, the two extreme values are shown. The symbols indicate the number of simulations available. The circle, cross and triangle indicate, respectively, one, two and three simulations. Regions without symbols are extrapolations. To facilitate the comparison with Fig. 8, the boundary marked by the black lines has been copied from Fig. 8 without modification.

Open with DEXTER

Appendix B: Convergence tests

Appendix B.1: Three-dimensional box

thumbnail Fig. B.1

Spacetime diagram of two 3D runs with particle size τf = 0.03 and Δ = 0.05. The color indicates the column density Σsolid/ ⟨ Σsolid using the same color bar as Figs. 4 and 5. Both runs begin with Z = 0.01 and have Z increased only for a short interval. In one run (left), Z grows to 0.02 and in the other (right) Z grows to 0.03. Clumps are visible for Z ≥ 0.02, consistent with the 2D runs.

Open with DEXTER

thumbnail Fig. B.2

A snapshot of run 3D.Z2 taken at t = 150 orbits (Z = 0.02). Image a) is the view from the side (xz plane). Image b) is the view from above (xy plane). The color scale marks the column solid density normalized to the initial gas density (Σsolid/ Σgas,0). The two images have a different color scale, since the column density in the azimuthal direction is much higher. Note that the particle structure is nearly axisymmetric.

Open with DEXTER

thumbnail Fig. B.3

A sequence of vertical snapshots (x-z plane) showing the formation of distinct particle clumps in one of our 2D simulations (left) and our two 3D simulations (middle, right). For each run, the images are taken 10 orbits apart. The time and Z values are shown in each plot. The times were chosen to show the formation of the first distinct clumps, and to give a similar peak density on the final image (ρp,max = 1.8,1.5,2.1left to right). The color indicates the column density using the same color scale as in Fig. B.2a. The other 2D runs give similar results to the one shown. The 2D run is consistent with the 3D runs. In particular, the 3D runs do not show any evidence of additional turbulence stirring compared to the 2D runs.

Open with DEXTER

thumbnail Fig. B.4

Maximum particle density (ρp,max) for one of the 2D simulations and the two 3D simulations. The top plot compares run 3D.Z2 (blue) and a 2D run (red). The bottom plot compares run 3D.Z3 (blue) and the same 2D run (red). In each plot the purple line shows ρp,max for the 3D run after averaging along the azimuthal direction. The vertical dashed lines mark the places where the 3D runs transition from Z = 0.01 to Z = 0.019 (for 3D.Z2) or Z = 0.029 (for 3D.Z3). Similarly, the vertical solid lines mark the progressive increase of Z during the course of the 2D runs.

Open with DEXTER

Our simulations are axisymmetric, they neglect the effect of the Kelvin-Helmholtz Instability (KHI). Whether the KHI is present is dictated by the Richardson number (Chandrasekhar 1961), where u is the gas speed, g = Ω2z is the vertical gravitational acceleration, and ρ is the effective fluid density (treating gas and solids as a single fluid). Bai & Stone (2010a) have shown that the streaming instability is able to keep the Richardson number above the critical value needed for the KHI, so that the KHI is absent in this problem. To confirm this, we performed two 3D simulations with τf = 0.03 and Δ = 0.05. The experiment performed in our 3D simulation is slightly different from that of our 2D simulations. Instead of increasing Z gradually, we only test the behavior of the system at three discrete values: Z = 0.01, Z = 0.02 and Z = 0.03. The reason for this is that 3D simulations are very computationally expensive, especially for small particles, so that a long-term simulation is prohibitive. Our 3D runs have three phases:

  • Phase I: each run starts with Z = 0.01, and Z is held constant for 30 orbits.

  • Phase II: the value of Z is increased quickly. In run “3D.Z2” we increase Z to Z = 0.02, and in run “3D.Z3” we increase Z to Z = 0.03.

  • Phase III: the value of Z is once again held constant for the remainder of the run.

Figure B.1 shows the spacetime diagrams for both runs. The figure shows that both runs produce visible particle clumps for Z ≥ 0.02 and no clumps for Z = 0.01. It is also clear that clumps form more readily with Z = 0.03. These results are consistent with our 2D simulations (Fig. 4). For particles smaller than τf = 0.03, the small integration steps make 3D simulations prohibitive. However, as explained in Sect. 2.6, the long response time of very small particles (τf ≤ 0.0064) means that our simulations will overestimate Z; so our results can be considered robust upper bounds for this size range.

Figure B.2 gives a side view (xz axis) and a top view (xy axis) of the 3D.Z2 run at t = 150 orbits. The clumps and their filamentary structure are very prominent. An important feature of the 3D runs is that they form visible clumps more easily than the 2D runs. Figure B.3 shows the formation of clumps in the two 3D runs and one of the 2D runs. This figure shows that the 3D runs are consistent with the 2D run. In particular, the clumps are readily visible at Z = 0.02 and the 3D runs show no evidence of additional stirring from the KHI.

Figure B.4 shows the maximum particle density (ρp,max) for the two 3D simulations and one of the 2D simulations. One salient feature of the figure is that in the early phase of the simulations, before clumping occurs, the 3D runs consistently show a peak density 23 times greater than the 2D runs. This occurs because particles can also accumulate along the azimuthal direction. The figure also shows ρp,max for the 3D runs after averaging along the azimuthal direction (purple lines). This averaging removes the apparent discrepancy between the 2D and 3D runs. In the later stages the 2D runs reach higher peak densities because they have very high Z values.

Appendix B.2: Resolution

We performed two simulations at 256 × 256 grid resolution with τf = 0.003 particles. Figure B.5 shows the spacetime diagrams for both sets of simulations. Since the behavior of the 256 × 256 runs is almost identical to the corresponding 128 runs, we conclude that, in terms of resolution, our simulations are fully converged.

thumbnail Fig. B.5

Spacetime diagrams showing the solid surface density Σsolid as a function of the radial coordinate x and simulation time τf = 0.003 and τf = 0.010 particles. The surface density is shown by color, using the same color scale as in Figs. 4 and 5. The figure shows both 128 × 128 simulations and the corresponding 256 × 256 simulations. In all simulations Δ = 0.05. The higher-resolution runs produce results consistent with the 128 × 128 runs.

Open with DEXTER

© ESO, 2015

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.