Subscriber Authentication Point
Free Access
 Issue A&A Volume 586, February 2016 A20 14 Planets and planetary systems https://doi.org/10.1051/0004-6361/201527533 21 January 2016

1. Introduction

Protoplanetary disks are the sites of planet formation. In these disks, the coagulation of microscopic dust particles, already present in the interstellar medium, into kilometer-sized planetesimals constitutes the first – and arguably least understood – step in the assembly of fully grown planets (Testi et al. 2014; Johansen et al. 2014). Initially, the dust aggregates, held together by surface forces, grow by sticking to each other in gentle, low-velocity collisions (e.g., Kempf et al. 1999). As a result, aggregates form very open, porous structures. As the aggregates gain mass and their relative velocities increase, collisions become more energetic, possibly leading to compaction and catastrophic fragmentation (Dominik & Tielens 1997; Blum & Wurm 2000; Ormel et al. 2007). A second hurdle is presented in the form of the radial drift barrier: particles with certain aerodynamic properties decouple from the gas and drift radially on short timescales (Whipple 1972; Weidenschilling 1977; Brauer et al. 2008).

In the inner few AU of the protoplanetary disk, dust grains consist mainly of silicates, and these aggregates bounce off each other in collisions or even disrupt completely upon impact at collision velocities above several m s-1 (Blum & Wurm 2008; Güttler et al. 2010). These collisional processes limit the growth beyond a centimeter or so in the inner disk (Brauer et al. 2008; Zsom et al. 2010; Windmark et al. 2012).

Outside the snow line, located typically at ~ 3 AU (Min et al. 2011), water ice becomes an important constituent of the dust grains. This is beneficial for growth because aggregates composed mostly of ice are capable of sticking at tens of m s-1 (Wada et al. 2009, 2013; Gundlach & Blum 2015). In addition, these icy particles maintain highly porous structures (Suyama et al. 2008, 2012), making them less likely to bounce in collisions (Wada et al. 2011; Seizinger & Kley 2013) and allowing them to out-grow the radial drift barrier in the inner ~ 10 AU of the protoplanetary nebula (Okuzumi et al. 2012; Kataoka et al. 2013a). However, the growth of these porous aggregates might be frustrated by high-velocity erosive collisions (Krijt et al. 2015) or sintering in certain regions of the disk (Sirono 2011).

Instead of coagulating directly, planetesimals can also be formed through particle concentration mechanisms (Johansen et al. 2014, and references therein). One promising mechanism is the streaming instability (SI; Youdin & Goodman 2005; Johansen et al. 2007; Bai & Stone 2010a,b), which can be triggered by a dense midplane layer of partially decoupled dust particles. Recently, Dra¸żkowska & Dullemond (2014) have defined a set of conditions for SI and compared them to dedicated models of compact coagulation. They found that in the inner disk, where the growth of silicates is limited by bouncing/fragmentation, particles cannot grow to Stokes numbers1 large enough (St ~ 10-2−1) for triggering SI. Outside the snow-line however, rapidly growing highly porous ice aggregates can grow to large Stokes numbers (Okuzumi et al. 2012) at which point their growth is possibly limited by erosive collisions (Krijt et al. 2015). The possibility of triggering SI through rapid porous coagulation has not yet been investigated, but it has been shown that erosion-limited porous growth can concentrate most of the dust mass in St ~ 1 particles (Krijt et al. 2015).

We set out to study the formation of the first generation of planetesimals. Giant planets have not yet formed, hence the protoplanetary disk is smooth. The focus is to understand the evolution of the mass-dominating particles in disks around Sun-like stars and understand how their evolution influences the dust surface density. Ultimately, the goal is to identify regions in both space and time where the first planetesimals can form, either through direct porous coagulation (e.g., Okuzumi et al. 2012) or through coagulation triggering SI (Dra¸żkowska & Dullemond 2014).

In order to answer these questions, we develop a global, semi-analytical, panoptic model that captures the evolution of the mass-dominating particles as they grow and drift radially in the protoplanetary disk (Sect. 3.1), including a detailed description of grain porosity (Sect. 3.2) and erosion (Sect. 3.3). After testing the method against two well-studied cases in Sect. 3.4, we use it to study rapid porous growth through the drift barrier in Sect. 4 and erosion-limited porous growth as a possible cause for SI (Sect. 5). The results are discussed in Sect. 6 and conclusions are presented in Sect. 7.

2. Nebula model and dust properties

We consider a turbulent disk of gas around a 1 M star. Neglecting the possible presence of pressure bumps, dead zones, etc., we concentrate on the outer regions where ices are an important part of the solid mass reservoir.

2.1. Nebula model

We adopt a truncated power-law for the gas surface density distribution, (1)The normalization constant will be determined by fixing the total mass of the gas disk MD, using (2)Typical power law exponents range between γ = 3/2, appropriate for the minimum mass solar nebula (MMSN; Hayashi 1981); and γ ≃ 1, consistent with observations (Andrews et al. 2009) and accretion disk theory (e.g., Armitage 2010). Disk masses between MD = 10-3M and 0.2 M and radii between 30−100 AU are consistent with observational constraints for disks in the Taurus star forming region (Andrews & Williams 2005; Andrews et al. 2013), though the drop in surface density in the outer edge is found to be much less dramatic than the drop-off in Eq. (1).

Irrespective of the disk mass, we adopt a temperature structure (3)appropriate for an optically thin disk and in agreement with observational constraints (Andrews & Williams 2005).

Most other quantities are derived, together with assumptions about the turbulence and vertical structure, from Eqs. (1) and (3). The isothermal gas sound speed is given by (4)with kB Boltzmann’s constant and mg = 2.34mp = 3.9 × 10-24 g with mp the mass of a proton. The Kepler frequency equals (5)Assuming an isothermal column, the gas density drops with increasing distance from the midplane z according to (6)with the vertical scale height of the gas hg = cs/ Ω. The turbulent viscosity is parametrized as following Shakura & Sunyaev (1973). Unless noted otherwise, we will use the disk parameters listed in Table 1.

Table 1

Benchmark disk model used throughout this work.

2.2. The dust content

Initially, the dust density follows the gas density through Σd/ Σg = Z0, with Z0 = 0.02 the vertically integrated metallicity. The individual dust particles are assumed to be compact spherical monomers with radii a and masses . We use ρ = 1.4 g cm-3, appropriate for icy particles. When dust particles grow, they are uniquely described by a mass m and volumetric filling factor φ (e.g., Krijt et al. 2015). Depending on the details of the growth process (i.e., the disk properties and collision physics), the filling factor (a measure for the porosity) can be as low as ~ 10-5 or very close to 1, resembling compact spheres.

2.2.1. Vertical and radial motions

The relative dust scale height is given by (Youdin & Lithwick 2007) (7)The dust scale height depends on the particle Stokes number Ωts, with the stopping time ts a function of the particle mass and porosity, and the local gas properties. Appendix A describes how the stopping time can be calculated.

The radial drift velocity is then given by (Weidenschilling 1977) (8)where vK = RΩ is the Keplerian orbital velocity and η can be written as (Nakagawa et al. 1986) (9)The drift timescale is defined as (10)and depends on the masses and porosities of the dust particles through their dimensionless stopping time Ωts (see Appendix A). Equation (8) neglects the backreaction of the dust particles on the gas and overestimates the drift velocity when the dust density becomes comparable to the gas density. In most models presented here such conditions are not reached and identifying regions and times where such collective effects become important is one of the goals of this paper.

2.2.2. Growth timescales in the single-size approximation

The model for porous coagulation is based on the semi-analytical model of Krijt et al. (2015, Sect. 5). At the heart of the semi-analytical approach lies the assumption that the local dust population can be approximated by a mono-disperse grain population, with a single characteristic mass and characteristic porosity. This assumption is valid when i) the full mass distribution has a clearly defined peak mass and porosity; and ii) the growth (and porosity-evolution, if dominated by collisions) of the peak-mass grains is mainly due to collisions with similar-size particles. These assumptions generally hold for populations resulting from porous coagulation (e.g., Ormel et al. 2007; Okuzumi et al. 2012).

For a mono-disperse dust population, and assuming all collisions result in perfect sticking, the growth rate is given by (11)with hd the dust scale height, σcol the collisional cross section, and vrel the relative velocity between the (identical) dust particles. The growth timescale is then defined as (12)Appendix B describes how to calculate vrel between two grains with arbitrary properties.

In reality however, the dust size distribution will not be infinitely narrow and dust grains will grow by colliding with dust particles with a broad range of masses and porosities, see for example Okuzumi et al. (2012, Fig. 9) and Krijt et al. (2015, Fig. 10). In order to simulate the width of the distribution, instead of calculating the velocity between identical particles, we calculate the velocity between two particles with Stokes numbers Ωts and κΩts, with 0 ≤ κ ≤ 1 a numerical factor2. The two limiting cases refer to the relative velocity between identical particles (κ = 1) and between a particle and an extremely well-coupled small grain (κ = 0).

In Fig. 1, we compare Eq. (11) to the full Monte Carlo model of Krijt et al. (2015). For this comparison, radial drift is neglected and perfect sticking is assumed. For every particle mass, the corresponding porosity is calculated by assuming an initial hit-and-stick phase that is followed by gas-compaction and eventually self-gravity compaction (see Sect. 3.2 and Appendix C for more details). The disk model is identical to the one in Okuzumi et al. (2012) and Krijt et al. (2015) with γ = 3/2,Z0 = 0.01, Σg(5 AU) = 152 g cm-2, and α = 10-3.

Overall, the semi-analytical model captures the growth timescales very well. The largest discrepancies are seen in the mass-range where growth is dominated by turbulent velocities, but where the grains are still relatively well coupled to the gas: ts<tη. In that regime (see Appendix B and in particular Eq. (B.2)) the relative velocity scales with the difference in Stokes numbers and is very sensitive to the choice of κ. While a closer correspondence with the MC model could be obtained by having a mass-dependent κ, we will use a constant value of κ = 0.5 in the remainder of this work. The same value is used by Sato et al. (2015).

 Fig. 1 Growth timescales of highly porous aggregates in an MMSN-like disk at 5 AU with a turbulence of α = 10-3. The solid black line is obtained with Monte Carlo models that take into account the full particle size and porosity distributions (Krijt et al. 2015) and the colored lines show Eq. (11) for three different values of κ (see text). Open with DEXTER

 Fig. 2 Illustration of the batch method described in Sect. 3. A) Globally, the dust content of the protoplanetary disk is described by a collection of unique and independent batches. The trajectories that individual batches follow in space and time are dubbed lifelines. B) Each batch is composed of three separate legs. These legs can move closer together or further apart, but the mass enclosed between the inner two and outer two legs is always identical. For each leg, growth and radial drift are solved simultaneously using Eqs. (13) and (14), while the batch’s surface density and its first derivative (needed to calculate the growth timescales) are obtained using Eqs. (17) and (20). C) At every leg location, the local dust component is described by a monodisperse aggregate population with a single mass, size, and porosity. D) The aggregates themselves are composed of (sub)micron-sized monomers, whose properties profoundly affect the aggregate’s strength against fragmentation and compaction (see Appendix C). Open with DEXTER

2.2.3. Bouncing, fragmentation, and erosion

In reality, not every collision will result in sticking and particles can bounce off each other, erode one another, or even completely fragment (Blum & Wurm 2000, 2008). The regimes in which these processes occur are usually defined by relatively sharp threshold velocities. In realistic models, the transition velocities for these collisional outcomes to occur are complex functions of the particle size(s), porosity, and material properties (e.g., Güttler et al. 2010).

Though only few experimental investigations have been performed with icy grains, it is clear that they are much stickier than their silicate counterparts (Dominik & Tielens 1997; Wada et al. 2013; Gundlach & Blum 2015). Moreover, they are expected to grow into porous structures (Okuzumi et al. 2012), making them less likely to bounce (Wada et al. 2011; Seizinger & Kley 2013). Therefore, we do not include bouncing and fragmentation. Recently however, Krijt et al. (2015) argued that erosion by small particles can be an effective way of halting growth when radial drift starts to play a role. We will include the effects of erosion on the growth timescale in Sect. 5.

3. The batch method

Here we introduce a new numerical approach for modeling the evolution of solids in a protoplanetary disk on a global scale. The method differs from grid-based methods (Brauer et al. 2008; Birnstiel et al. 2010; Okuzumi et al. 2012) and Monte Carlo coagulation models (Ormel et al. 2007; Zsom & Dullemond 2008) in a fundamental way, as it employs a Lagrangian approach, in which dust particles are followed through the disk as they grow through collisions and drift radially through the disk. A similar Lagrangian way of thinking is found in particle-tracking models that simulate diffusion in evolving disks (e.g., Ciesla 2011; Hughes & Armitage 2012), but there are important differences between those methods and the approach presented here (see Sect. 6.1).

In this section, we first introduce the principles of the method (Sect. 3.1), describe how to include porosity and erosion (Sects. 3.2 and 3.3), and test the method against two well-studied cases (Sect. 3.4).

3.1. Description

In our approach we follow the evolution of a small groups of dust particles, referred to as a batch. The idea is to follow individual batches in time, while the constituent dust grains grow through collisions and drift inward as the result of radial drift. Specifically, we are interested in solving the Lagrangian (i.e., co-moving with the batch) derivatives for dust grains (13)and radial drift (14)where the drift- and growth timescales are given by Eqs. (10) and (12). By following a single batch in both time and space, information about the particle’s history is easily accessible. The key assumption of this method is that the dust particles inside a batch grow exclusively by colliding with similar particles. Similar not only because they have comparable masses and porosities, but also because they share the same growth- and drift histories. This assumption allows us to calculate the evolution of multiple batches independently, without necessarily knowing how the dust in the rest of the disk is evolving. The validity of this assumption will be discussed further in Sect. 6.

In order to capture the evolution of the dust surface density, we will give the batches a finite radial width. One could imagine following the front-side (slightly closer to the star) and back-side (slightly further from the star) of the batch separately. When they move closer because of small differences in the drift velocity at both sides, the surface density of dust increases; when they move away from each other, it decreases. However, this does not allow the surface density gradient to change. It transpires that we need three coordinates: one at the inner edge, one at the outer edge, and one at the center of mass. At the location of each leg, the dust is described by a mono-disperse distribution, the properties of which are initially identical, but will start to differ slightly in time because of the dependence on location in the details of the coagulation process. To create a complete picture of the dust evolution on a global scale, we will combine the results of many batches that start out at different locations. The approach is illustrated in Fig. 2.

An alternative approach would be to calculate batches with a single leg and obtain a single global function for the surface density by looking at the distribution of all the batches simultaneously. However, in such an approach, batches cannot overtake each other, as this gives rise to infinite surface densities. An important feature of our method, is that the local surface density – determined by the spacing of a batch’s three legs – is a property of the batch itself. This renders the batches truly independent, allowing one to calculate their evolution separately, but also making it possible for two batches with different constituent particles and different surface densities to occupy the same location. We will encounter such a situation in Sect. 4.

At a time t, the state of a batch is then fully described by the vector (15)where indices i, c, o correspond to the inside (or front), center of mass, and outside (or back) of the batch respectively, and by the batch’s dust surface density, assumed to be a power-law with index p(16)with both Σ0 and p functions of time. At t = 0, both Σ0 and p can be obtained from Eq. (1). Specifically, and p(t = 0) = γ. Later, when the dust is evolving, the surface density within the batch is acquired from the locations of the three ‘legs’ at the front, back, and center of mass of the batch. First, given that there is equal mass between ri and rc as there is between rc and ro, the slope can be found by solving (17)which can be approximated as (18)Second, the total mass inside the batch (19)is conserved. Since MB is known at t = 0, Eq. (19) can be inverted to obtain Σ0 at later times as a function of the location of the three legs3. (20)For the starting conditions, we will assume all dust particles start out as monomers, such that mi = mc = mo = m at t = 0. The initial locations of the three legs can be chosen freely, but batches should not become too wide, i.e., (rori) /rc ≪ 1.

Now everything is in place to calculate the evolution of a single batch. For a given disk model, dust properties, and batch starting conditions, the evolution of a single batch, i.e., x(t), can be calculated by integrating Eqs. (13) and (14) simultaneously for the three legs, while using Eqs. 17 and 20 to calculate the surface density at their locations. To solve this initial value problem, we have made use of the scipy.integrate.ode package of python.

By connecting the properties4 of all batches (e.g., particle mass, surface density) at specific times, the disk-wide dust distribution can be obtained as a function of time. As the batches are completely independent, there are no restrictions on their distribution. For example, it is not necessary that the edges of one batch are connected to edges of the adjacent ones. In most of this work, we will use between 102 and 103 batches5, distributed logarithmically between 3 AU and rout.

 Fig. 3 Evolution of the dust surface density in time. A) No growth case. Particles always have a = 1 mm. The dashed lines show analytical results of Youdin & Shu (2002). B) Compact growth with perfect sticking. Dust particles start out as 0.1 μm monomers and grow on timescales given by Eq. (11), while staying compact at all times (φ = 1). The gas surface density is shown by the dotted line. The results are qualitatively different, with the no-growth case resulting in pile-ups (locally increasing Σd/ Σg) and the compact coagulation model showing an inside-out removal of dust (decreasing Σd/ Σg everywhere). Open with DEXTER

3.2. Treatment of porosity

In principle, the vector x of Eq. (15) can be expanded to include additional dust properties such as porosity, chemical composition, water fraction, etc., provided that a recipe for determining the corresponding time-derivative exists.

In the case of porosity, obtaining an expression for (Dφ/Dt) that covers fractal growth and the various different compaction mechanisms is not straightforward. Therefore, in this work, we will use a different method. Instead of adding the porosity to the vector in Eq. (15), we will construct a unique function that gives the filling factor as a function of the current particle mass and location, i.e., φ = φ(m,r). Appendix C discusses this function in more detail. Once φ(m,r) is available, we can use it whenever we need to know a particle’s porosity, for example when calculating the collisional cross section or the Stokes numbers.

3.3. Including erosion

The growth rate of Eq. (11) assumes all collisions result in sticking. However, as shown by Krijt et al. (2015), erosive collisions have the potential to effectively halt growth when the drift velocity exceeds the threshold velocity veros needed for erosion. While the method outlined in this section does not provide information about the population of small projectiles, we can nonetheless mimic the effect of efficient erosion by adjusting the growth timescale as follows (21)where is the relative velocity between the particle in question and a fiducial monomer grain. In this way, the growth timescale rapidly increases when the drift velocity exceeds the erosion threshold. Both catastrophic fragmentation and bouncing could be simulated in a similar fashion, but in that case the relevant vrel becomes the velocity between similar-size particles. In this work however, as we are limiting ourselves to the study of sticky ices in the outer disk, we only include erosion as a possible growth barrier.

3.4. Test cases

In the remainder of this Section, we illustrate and test the approach outlined above by comparing it against two well-known cases: the drift-only case (neglecting coagulation completely) and the compact-growth case (without fragmentation).

3.4.1. Test case I: drift-only

When the particle size is kept constant in time (i.e., no coagulation takes place), grains in the outer disk drift faster than particles in the inner disks, causing pileups and potentially planetesimal formation (Youdin & Shu 2002; Youdin & Chiang 2004). To study the drift-only case, we use the disk model of Table 1 and set = 0 everywhere. The dust particles are compact icy grains with a radius of 1 mm.

Figure 3A shows the resulting dust surface density evolution. Since, for a single grain size, Stokes numbers are largest in the outer part of the disk, particles further out will drift faster. As a result, the outer disk is slowly depleted of dust and pile-ups are created closer in (e.g., Youdin & Shu 2002; Youdin & Chiang 2004; Birnstiel & Andrews 2014). The pile-ups are caused by the slowing down of particles as they drift into the denser inner disk and are the result of the assumption of a fixed size of the particles. Below we will see that when grain sizes are limited by aerodynamical properties, no such pile-ups are created. Our model is in excellent agreement with the analytical prediction for this case Youdin & Shu (2002, Eq. (28)), shown by the overplotted dashed lines. The razor-sharp outer edge of the dust distributions in Fig. 3A is not actually resolved, but simply marks the location of the outermost batch at any given time. When an exponentially decaying disk is assumed instead of Eq. (1), the dust distributions at later times are smoother (e.g., Youdin & Shu 2002, Fig. 5). The outer edge of the dust distribution will be discussed further in Sect. 6.

 Fig. 4 Locations of dust batches in time for four standard cases. A) No growth, particles are always 1 mm and compact. B) Compact growth starting from 0.1 μm particles. C) Porous growth with perfect sticking. Successful planetesimal formation (Ωts> 103) through coagulation is indicated by the (•)-symbol. The lifeline of the batch that formed planetesimals furthest out is colored red. D) Porous growth followed by erosion above velocities of 20 m s-1. Dashed lines indicate the particles in the batch have 10-2< Ωts< 3 and background colors indicate the midplane dust-to-gas ratio of particles that have Ωts ≤ 3. Open with DEXTER

In Fig. 4A, we plot the lifelines of a selection of batches for the same simulation. In this representation, batches move up in time and to the left as they start to drift in. The places where lines move close together correspond to the pile-ups in Fig. 3A. An important observation is that while batches get close to each other, the lifelines do not cross, justifying the assumption that they evolve independently. Dashed lines indicate the Stokes numbers of the dust particles inside the batch are 10-2< Ωts< 3. Since particles are always mm-sized in the no-growth scenario, Stokes numbers are highest in the outer disk, with the transition Ωts = 10-2 around 10 AU. In addition, the background colors indicate the midplane dust-to-gas ratio, calculated as (22)where the fraction of scale-heights on the right-hand-side is a function of the particle Stokes number (see Eq. (7)). Initially, because of settling, the highest midplane dust-to-gas ratios are found exclusively in the outer disk. Later, when pile-ups are formed, high dust-to-gas ratios are reached at smaller radii as well. After ~ 105 yr, all batches have drifted inside the snow-line and the solid content of the disk has vanished.

3.4.2. Test case II: compact growth

The next test includes growth, but assumes that the particles stay compact (φ = 1) at all times. For this test, the dust particles start out as monomers with a = 0.1 μm and their growth timescale is calculated using Eq. (11). Disk-wide numerical simulations of compact coagulation (with or without fragmentation) and drift have been performed by many authors, including Brauer et al. (2008), Birnstiel et al. (2010), Okuzumi et al. (2012), Birnstiel & Andrews (2014).

Figures 3B and 4B show the results of the panoptic model for compact growth in the disk of Table 1. The behavior of the dust component is fundamentally different from the no-growth case in a number of ways. We focus on Fig. 4B first. Since particles start out as well-coupled sub-micron monomers, Ωts ≪ 10-2 everywhere in the disk and radial drift is negligible. However, Brownian motion and the presence of turbulence causes the particles to collide and grow. As the grains gain mass, they slowly decouple from the gas and their Stokes number increases. At some point, the drift timescale becomes comparable to the growth timescale and the particles will start to drift inward. Because the timescales for growth are shortest at small radii, grains that start out in the inner disk reach high Stokes numbers much earlier. This causes an inside-out clearing of the dust component, very different from the no-growth scenario. As illustrated in Fig. 3B, the inside-out removal of dust does not allow the formation of pile-ups and the dust surface density decreases in time at every location.

For turbulence-driven growth in the Epstein drag regime (as is the case here), a pure drift-growth balance leads to a dust surface density profile Σd ∝ (Σgr-2Ω-2)1/2 (e.g., Birnstiel et al. 2012). For γ = 3/2, this results in Σdr-1. The surface density curves of Fig. 3B are in excellent agreement with this prediction and with grid-based numerical models of similar disks (e.g., Okuzumi et al. 2012, Fig. 3).

4. Porous growth with perfect sticking

While compact growth typically results in an inside-out removal of the disk’s solids, highly porous growth could potentially lead to planetesimal formation by allowing fluffy icy grains to quickly grow through the radial drift barrier (Okuzumi et al. 2012). In this Section, we will study such highly porous growth using our Lagrangian approach and quantify in which regions of the protoplanetary disk such rapid growth is feasible.

The set-up is the same as in Sect. 3.4.2, but the dust particle porosity φ(m,r) is now calculated assuming hit-and-stick growth followed by compaction by gas-pressure and eventually self-gravity (Sect. 3.2 and Appendix C). When the aggregates inside a batch grow to Stokes numbers Ωts> 103, we classify them as planetesimals6 and no longer follow their evolution.

Figure 4C shows the result of the porous growth case assuming 0.1 μm icy monomers and perfect sticking. The ()-symbols indicate successful planetesimal formation. Our results confirm the findings of Okuzumi et al. (2012), Kataoka et al. (2013a) and Krijt et al. (2015), indicating that growth through the radial drift barrier is possible in the region just outside the snow-line. Indeed, dust can grow from sub-micron to planetesimal sizes within ~ 104 yr. In the outer disk, the behavior is similar to that in Fig. 4B: dust grains grow to Ωts ~ 1 and drift inward. We note, however, that in the porous case these Ωts = 1 grains have masses and porosities that differ from the compact grains by many orders of magnitude (see Fig. 3 of Johansen et al. 2014).

4.1. Planetesimals and pebbles

Using Fig. 4C, three distinct regions can be identified inside the protoplanetary disk. First, we identify the batch that formed the outermost planetesimals. The lifeline of this batch is shown in red in Fig. 4C. Dust particles that are located inside of this batch at t = 0 will eventually grow into planetesimals on a timescale of ~ 104 yr.

Material located outside of the red batch at t = 0 will grow to Ωts ~ 1 and drift in. We call these particles pebbles7. Around 104 yr, material from the outer disk will start to drift into the planetesimal belt located just behind the snow-line. Since our model at this point does not include any interaction between batches, the pebbles fly through the planetesimals unaffected. In reality, the planetesimals can potentially accrete the pebbles very efficiently (Ormel & Klahr 2010; Guillot et al. 2014). In Fig. 4C, material that starts out outside ~ 10 AU will form pebbles. For the disk profile of Eq. (1), the fraction of the disk’s total mass located outside r equals (23)Thus, the pebbles formed outside of 10 AU together hold ~ 2/3 of the total dust mass and have the potential to boost protoplanet growth after the formation of planetesimals. We will discuss pebble accretion further in Sect. 6.

 Fig. 5 Influence of total disk mass on the location where planetesimals and pebbles form (Sect. 4.1). Material that starts out left of the solid blue line successfully grows through the drift barrier, eventually forming planetesimals a bit closer in. The red outermost planetesimal of Fig. 4C has been plotted for the MD/M⊙ = 10-2. Material that starts out in the gray region cannot coagulate fast enough, experience rapid radial drift, and at some point drift past the already-formed planetesimals in the inner disk. For reference, the current locations of the four Jovian planets are indicated by capital letters. Open with DEXTER

4.2. Influence of disk mass

In Fig. 5 we combine simulations performed with different disk masses (the other properties, i.e., turbulence, disk outer radius, metallicity, are identical to those in Table 1) and show how disk mass influences the location where planetesimals and pebbles originate. The figure should be read as follows: For a given disk mass, material that starts to the left of the solid blue line will form planetesimals in the ‘planetesimal zone’ marked by the shaded region. Material that starts out in the gray area will grow and drift to enter the previously formed planetesimal zone at Stokes numbers 10-2< Ωts< 3. The red lifeline of Fig. 4C is also plotted at the corresponding disk mass of 10-2M.

From the figure we can identify two trends. First, more massive disks form planetesimal at larger radii. This is to be expected, as more massive disks have a higher dust surface density, making growth through the radial drift barrier easier (e.g., Okuzumi et al. 2012, Sect. 4). For disk masses comparable to the MMSN, planetesimal formation around the current location of Jupiter is a realistic possibility. Second, the boundary between the regions where planetesimals and pebbles originate, moves out for increasing disk mass. This means that for more massive disks, a larger fraction of the dust content can form planetesimals directly, leaving a smaller fraction available to be accreted in the form of pebbles later. Nonetheless, for the most massive disks considered here the pebble origin region still contains ~ 1/2 of all the dust.

5. Porous growth with erosion

So far we have assumed that all collisions result in perfect sticking. And while collisions between icy same-size aggregates in the outer disk are indeed likely to occur below the critical fragmentation threshold, Krijt et al. (2015) have shown that high-velocity collisions with large mass ratios can result in efficient erosion and possibly frustrate the further growth of aggregates with Stokes numbers Ωts ~ 1. Here, we mimic the effect of efficient erosion by adjusting the growth timescale according to Eq. (21), using an erosion threshold velocity of 20 m s-1 (Gundlach & Blum 2015).

The results of the simulation with erosion are shown in Fig. 4D. The largest difference with the perfect sticking scenario of Fig. 4C is that there are no planetesimals being formed in the inner disk. Instead, when particles grow to Ωts ~ 1, their growth is frustrated by erosion and they drift inward.

 Fig. 6 Similar to Fig. 4D (porous growth + erosion), but with A) an increased metallicity of Z0 = 0.06; B) a reduced turbulence strength of α = 3 × 10-4; and C) a gas temperature 50% lower than that in Eq. (3). All three modifications result in a higher maximum dust-to-gas ratio compared to Fig. 4D. Open with DEXTER

5.1. Streaming instability

Apart from coagulating directly, planetesimals can also be formed through SI (e.g., Youdin & Goodman 2005; Johansen et al. 2007; Bai & Stone 2010a,b). In this section we investigate how porous coagulation followed by erosion can lead to conditions suitable for triggering SI.

We make use of the work of Dra¸żkowska & Dullemond (2014), who, based on the work of Johansen et al. (2007, 2009) and Bai & Stone (2010a,b) defined three conditions for triggering SI. For a mono-disperse particle distribution, these conditions are equivalent to:

• i)

The Stokes number of the mass dominating particles needs to beclose to unity. Typically, one needs 10-2 ≤ Ωts ≤ 3.

• ii)

The midplane dust-to-gas ratio of these particles needs to exceed, or be close to, unity.

• iii)

The vertically integrated metallicity should be at least a few times Solar, i.e., Z ≳ 0.02−0.03.

See also Fig. 8 of Carrera et al. (2015) for points i) and iii). Lastly, SI needs time to develop and will not be triggered if the growth timescales of the particles are too short. In other words, tgrow> Ω-1. The first two conditions are related to efficient momentum transfer between dust particles and the gas. Particles with much higher Stokes numbers do not effectively interact with the gas, while particles with with much smaller stopping times do not result in strong clumping. Condition iii) is related to suppressing midplane turbulence. For low metallicities, the strength of this turbulence drops sharply (Bai & Stone 2010a).

Concentrating on conditions i) and ii) and armed with the semi-analytical model of Sect. 3, we can identify regions in space and time where these conditions are met. In Fig. 4, dashed lines indicate the particles in the batch meet condition i), while the background colors show the midplane dust-to-gas ratio (Eq. (22)), important for condition ii).

Focussing on the porous growth + erosion case (Fig. 4D), we see that conditions for SI are not reached. While every batch eventually forms particles with high Stokes numbers (evidenced by the dashed lines), the achieved midplane densities are about a factor 10 too low. We can see a front of moderate midplane densities moving outward in time, starting in the inner disk at 103 yr and reaching 100 AU after almost 106 yr.

This behavior can be understood in the following way. Essentially, the midplane dust-to-gas ratio is the combination of the vertically integrated metallicity and the degree of settling (24)where Z = (Σd/ Σg) and we have used that hd/hg ~ (α/ Ωts)1/2 for high Stokes numbers (see Eq. (7)). Second, since the removal of dust occurs from the inside-out, the surface density – at a given r – will decrease, very much like in Fig. 3B, and pile-ups like the ones seen in the no-growth case do not occur. As a result, the highest vertically integrated metallicity is achieved at the very beginning and ZZ0. The highest midplane densities are then achieved in the interval when growth toward Ωts> 10-2 has occurred, but before enough time has past for radial drift to reduce the local dust surface density. The diagonal shape of high (ρd/ρg)z = 0 in Fig. 4D occurs because growth takes longer in the outer disk.

5.2. Disk metallicity, turbulence, and temprature

Inspired by Eq. (24), we briefly discuss three ways that can help to create conditions suitable for SI. All three are shown in Fig. 6.

Increasing the initial dust content: the most straightforward way to reach higher midplane dust-to-gas ratios is to start out with a higher metallicity from the beginning. Figure 6A shows that for Z0 = 0.06, a mass loading close to unity can be achieved in the inner ~ 10 AU.

Decreasing the turbulence strength: the turbulence strength influences the degree of settling and for a given metallicity and particle size, the midplane dust-to-gas ratio is inversely proportional to . For example, reaching a mass loading of 1 in a column with Z = 0.01 and Ωts = 0.1 particles requires α ~ 10-5. However, the turbulence cannot be decreased to arbitrarily small values. The presence of marginally decoupled grains will give rise to Kelvin-Helmholtz instability (Weidenschilling 1980, 1995; Johansen et al. 2006) and SI (e.g., Takeuchi et al. 2012), resulting in turbulence that can be parametrized with for example Eq. (8) of Dra¸żkowska & Dullemond (2014). For the disk models used in this work, the strength of this turbulence amounts to αmp ≃ 5 × 10-4 and ~ 10-3 at 5 and 50 AU, respectively, for particles with Stokes numbers around unity.

Decreasing the temperature: the temperature structure used so far (Eq. (3)) is based on an optically thin disk. However, midplane temperatures in disks might be significantly lower, especially if the disk is optically thick (e.g., Andrews et al. 2009). In contrast, the presence of viscous heating can affect the temperature profile significantly (Bitsch et al. 2015). Lowering the temperature affects many things. Most important for the calculations presented here is that a lower temperature decreases the pressure gradient η (see Eq. (9)), causing the radial drift to be slower. This makes erosion less efficient (lower ηvK) and also gives grains more time to grow, allowing them to reach larger sizes. Reducing the temperature by 50% indeed results in higher midplane dust-to-gas ratios (Fig. 6C), but the effect is not large enough to trigger SI.

6. Discussion

6.1. Batch approach

The numerical approach introduced in this work (Sect. 3), in which individual batches of dust are followed as they grow and drift in the protoplanetary nebula, provides an intuitive and flexible way to calculate (porous) coagulation and the evolution of the global surface density simultaneously, connecting the properties of the (sub)micron monomers to the growing aggregates and the global redistribution of solids (Fig. 2). The flexibility and speed of this approach allow us to study the impact of different coagulation models (Fig. 4) and or disk properties (Fig. 6) quickly, while preserving the essential characteristics of the growth process.

At this point, the method has three main drawbacks. First, it traces only the mass-dominating particles and does not provide information about the number distribution for smaller masses. If the distribution can be assumed to be in growth/fragmentation equilibrium, the complete mass distribution may be reconstructed (e.g., Birnstiel et al. 2011), though this has not yet been attempted for porous growth or for a steady state between growth and erosion.

Second, batches do not interact with each other. In most situations encountered here this is justified, as evidenced by the fact that the different lifelines do not cross in Fig. 4. In some cases however, this assumption does not hold. One example, in which pebbles from the outer disk drift past planetesimals in the inner disk, is found in Sect. 4 and will be discussed below. Another situation where batches come uncomfortably close is in the outer edge of an exponential disk. Birnstiel & Andrews (2014) have shown that in exponentially decaying disks, a pile-up will be created at the outer edge (e.g., their Figs. 3 and 4). In theory, these pile-ups host all dust grains that were originally located further out. Coagulation inside such a pile-up would be hard to model with independent dust batches.

Third, the approach does not include radial mixing as a result of turbulent diffusion. With one of the underlying assumptions of the model being that dust grains that start out together evolve in a similar manner, it appears not to be straightforward to add these processes. Radial mixing and diffusion are accurately captured by particle tracking models that simulate individual trajectories for a large number of individual grains (e.g., Ciesla 2011; Hughes & Armitage 2012). Such models differ from the method presented here in a number of ways. First, in our approach, the dust surface density is solved from the evolution of a batch itself, rather than being obtained by looking at the dynamical evolution of all the particles in the simulation. Second, where we focus on the effects of the details of the coagulation process in a static disk, these works primarily look at evolving gas disks and either neglect particle coagulation (Ciesla 2011) or treat it in a simplified way; for example by injecting particles of increasing size at a later stages (Hughes & Armitage 2012). Despite these differences, the main conclusions of Hughes & Armitage (2012) and this work are similar: both do not find significant enhancements of the vertically-integrated dust-to-gas ratio.

6.2. Rapid porous growth

For porous icy aggregates, fragmentation occurs only at very high velocities (Wada et al. 2013) and the resulting fluffy aggregates can swiftly grow through the radial drift barrier (Fig. 4C). Rapid growth is harder to achieve in the outer disk, since timescales are longer and dust and gas spatial densities are lower. As a result, planetesimals form in a zone just behind the snow-line. Figure 5 shows how the size of this region depends on disk mass. Comparing the planetesimal formation zone to the current locations of the Solar System giant planets8, we find that while it is hard to form planetesimals all the way out to the current location Uranus and Neptune, it is very well possible to form planetesimals in the location where Jupiter is now. One could envision a scenario where enough planetesimals accumulate to trigger the formation of Jupiter early on. At that point, the assumptions in our model break down and the presence of Jupiter will steer the evolution of the disk and planet formation therein (e.g., Pollack et al. 1996; Pinilla et al. 2012; Kobayashi et al. 2012).

6.3. Pebble accretion

After being studied by Ormel & Klahr (2010), Lambrechts & Johansen (2012), pebble accretion has received a lot of attention as a robust way of growing planetesimals and planetary embryos (Lambrechts & Johansen 2014; Kretke & Levison 2014; Levison et al. 2015). The method developed in this work is ideally suited for studying pebble accretion. First, it provides the location and formation time of the planetesimal population (e.g., Fig. 4C), but it also accurately captures the properties of the pebbles that drift in later and provides information about their history (e.g., where the pebbles originated from, how long they spent in which part of the disk, etc.). The next step will be to include interaction between different batches when they overlap. Once could imagine syphoning mass from the batch with pebbles to the batch with planetesimals with a certain efficiency. This accretion efficiency depends sensitively on the properties of the local gas, the size of the planetesimals, and the aerodynamic properties of the pebbles (e.g., Guillot et al. 2014; Visser & Ormel 2015), information that is all readily available in our approach.

6.4. Erosion and streaming instability

In Krijt et al. (2015), it was seen that the growth of Ωts ~ 1 particles stagnates (i.e., → 0) soon after the drift velocity exceeds the erosion threshold velocity. The growth stagnated because a balance was reached between growth through similar-size collisions and erosion by small grains (Krijt et al. 2015, Fig. 10). Inspired by these results, we adopted Eq. (21) to simulate efficient erosion. While instructive, this treatment of erosion is very rudimentary, but since our monodisperse approach does not hold information about the population of small grains, a self-consistent treatment of erosion is not possible at this point. Unfortunately at the moment, no global codes are available that can self-consistently treat the porosity evolution and radial drift of the full mass distribution while taking into account destructive and erosive collisions.

In Sect. 5.1, we investigated when porous growth might lead to SI, making use of the conditions as defined by Dra¸żkowska & Dullemond (2014). Our simulations indicate that, in a smooth disk, reaching conditions for SI is not straightforward. While porous growth (possibly followed by erosion) naturally leads to a population of aggregates with 10-2 ≤ Ωts ≤ 3 in a large portion of the protoplanetary disk (Figs. 4C and D), the creation of a dense midplane layer of solids is problematic. The main reason for this is that growth followed by drift inevitably results in an inside-out removal of the solid content and does not lead to pile-ups (see also Fig. 3). These results are in line with Dra¸żkowska & Dullemond (2014), who, focussing on compact growth, concluded that SI can only follow if the bouncing barrier can be overcome, the local metallicity can somehow be enhanced, and the turbulence is sufficiently weak. A number of ways to increase the maximum midplane dust-to-gas ratio are discussed in Sect. 5.2 and depicted in Fig. 6. Such variations in the disk properties, perhaps working together, could result in conditions suitable for SI in a subset of protoplanetary disks.

6.5. Sintering

One intriguing possibility of creating regions of enhanced dust density is by having dust sintering in specific regions of the disk. Sintering is expected to lower the fragmentation threshold velocity in regions around the ice-lines of major volatile species (Sirono 2011), effectively creating alternating regions of fragmentation-limited and drift-limited growth. Recently, Okuzumi et al. (2015) have argued that the existence of such sintering regions can result in pile-ups of material, offering a possible explanation for the peculiar shape of the HL Tau protoplanetary disk (Brogan et al. 2015). While it appears sintering can play an important role in the radial redistribution of solids in some systems, there are still many uncertainties. For example, it is not yet clear how sintering affects the porosity evolution described in Appendix C or how exactly it alters the expected collisional outcomes.

7. Conclusions

We have developed a novel Lagrangian approach for calculating the evolution of the mass-dominating dust aggregates as they grow and drift in a protoplanetary disk. The method, summarized in Fig. 2, allows the calculation of the global evolution of the dust surface density on Myr timescales while preserving the essential characteristics of the porous growth process and can be used to study planetesimal formation and pebble delivery.

After testing the new approach against two well-known cases (Fig. 3), we use it to study the formation of the first generation of planetesimals – those that can form in a smooth disk structure – in disks around Sun-like stars. When fragmentation and erosion are inefficient, we find that:

• 1.

Planetesimals can coagulate very rapidly, within ~ 104 yr, around the current location of Jupiter. While they end up in the region just behind the snow-line, the planetesimals include material from a broader region that extends out to ~ 10 AU for an MMSN-like disk (Fig. 4C).

• 2.

For more massive disks, both the region where planetesimals eventually form and the region where they originate from move outward (Fig. 5). Thus, in these massive disks, a larger fraction of the dust content can directly form planetesimals, leaving less material to be accreted as pebbles later on.

These scenarios rely on the assumption that fragmentation and erosion are relatively unimportant at collision speeds up to several tens of m s-1. Alternatively, when erosion balances growth around Ωts ~ 1 (Krijt et al. 2015), further coagulation is not possible but conditions necessary for streaming instability (SI) might be reached. In these cases, we find that:

• 3.

While porous growth limited by drift-induced erosion is an effective way of creating aggregates with 10-2 ≤ Ωts ≤ 3 in a large region of the protoplanetary disk (Figs. 4D and 6), conditions needed for SI are generally not reached.

• 4.

The most stringent condition is creating and maintaining a dense midplane layer of solids. In a smooth gas disk, rapid porous growth followed by erosion leads to an inside-out clearing of the dust disk. In such a scenario, no pile-ups are created and the metallicity decreases.

• 5.

The highest midplane densities are reached in the inner disk first and then move out toward the outer parts. A reduced turbulence level and lower gas temperature increase the maximum midplane dust-to-gas ratios slightly, but the biggest impact comes from increasing the initial metallicity (Fig. 6).

Future improvements to the method, in particular the addition of interaction between batches, will help build a coherent picture of the planet formation process.

1

A measure for the degree of coupling between the dust particle and the surrounding gas.

2

See Sato et al. (2015) for further discussion on this approach in the case of compact growth.

3

When p = 2, Eqs. (18)–(20) break down. In that case, they should be replaced by rc = (riro)1/2 and .

4

For this purpose, we use the dust properties at the batch center of mass rc.

5

While the total number of batches used does not influence the evolution of the individual batches, a higher number does result in a better sampling of the disk at a given time. We found that using between 102 and 103 batches was sufficient to provide smooth curves in Fig. 3.

6

At 5 AU, these Stokes numbers correspond to aggregates with a mass of about ~ 1014 g and a porosity of ~ 1%.

7

While some batches starting outside of the red batch form planetesimals, these grains have at some point drifted past the already-formed planetesimals in the inner disk.

8

In the Nice model however, the giant planets of the Solar System were originally located much closer in, roughly between 5 and 17 AU, migrating out at later times (Tsiganis et al. 2005; Morbidelli et al. 2005; Gomes et al. 2005).

Acknowledgments

Dust studies at Leiden Observatory are supported through the Spinoza Premie of the Dutch science agency, NWO. The authors would like to thank C. P. Dullemond, A. Johansen, J. Dra¸żkowska, T. Birnstiel and S. Okuzumi for encouraging discussions and the anonymous reviewer for comments that helped improve the manuscript.

References

1. Andrews, S. M., & Williams, J. P. 2005, ApJ, 631, 1134 [NASA ADS] [CrossRef] [Google Scholar]
2. Andrews, S. M., Wilner, D. J., Hughes, A. M., Qi, C., & Dullemond, C. P. 2009, ApJ, 700, 1502 [NASA ADS] [CrossRef] [Google Scholar]
3. Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129 [NASA ADS] [CrossRef] [Google Scholar]
4. Armitage, P. J. 2010, Astrophysics of Planet Formation, ed. P. J. Armitage (Cambridge University Press) [Google Scholar]
5. Bai, X.-N., & Stone, J. M. 2010a, ApJ, 722, 1437 [NASA ADS] [CrossRef] [Google Scholar]
6. Bai, X.-N., & Stone, J. M. 2010b, ApJ, 722, L220 [NASA ADS] [CrossRef] [Google Scholar]
7. Birnstiel, T., & Andrews, S. M. 2014, ApJ, 780, 153 [NASA ADS] [CrossRef] [Google Scholar]
8. Birnstiel, T., Dullemond, C. P., & Brauer, F. 2010, A&A, 513, A79 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
9. Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2011, A&A, 525, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
10. Birnstiel, T., Klahr, H., & Ercolano, B. 2012, A&A, 539, A148 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
11. Bitsch, B., Johansen, A., Lambrechts, M., & Morbidelli, A. 2015, A&A, 575, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
12. Blum, J., & Wurm, G. 2000, Icarus, 143, 138 [NASA ADS] [CrossRef] [Google Scholar]
13. Blum, J., & Wurm, G. 2008, ARA&A, 46, 21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
14. Brauer, F., Dullemond, C. P., & Henning, T. 2008, A&A, 480, 859 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
15. Brogan, C. L., Pérez, L. M., ALMA Partnership, et al. 2015, ApJ, 808, L3 [NASA ADS] [CrossRef] [Google Scholar]
16. Carrera, D., Johansen, A., & Davies, M. B. 2015, A&A, 579, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
17. Ciesla, F. J. 2011, ApJ, 740, 9 [NASA ADS] [CrossRef] [Google Scholar]
18. Dominik, C., & Tielens, A. G. G. M. 1997, ApJ, 480, 647 [NASA ADS] [CrossRef] [Google Scholar]
19. Drążkowska, J., & Dullemond, C. P. 2014, A&A, 572, A78 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
20. Gomes, R., Levison, H. F., Tsiganis, K., & Morbidelli, A. 2005, Nature, 435, 466 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
21. Guillot, T., Ida, S., & Ormel, C. W. 2014, A&A, 572, A72 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
22. Gundlach, B., & Blum, J. 2015, ApJ, 798, 34 [NASA ADS] [CrossRef] [Google Scholar]
23. Gundlach, B., Kilias, S., Beitz, E., & Blum, J. 2011, Icarus, 214, 717 [NASA ADS] [CrossRef] [Google Scholar]
24. Güttler, C., Blum, J., Zsom, A., Ormel, C. W., & Dullemond, C. P. 2010, A&A, 513, A56 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
25. Hayashi, C. 1981, Progress of Theoretical Physics Suppl., 70, 35 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
26. Heim, L.-O., Blum, J., Preuss, M., & Butt, H.-J. 1999, Phys. Rev. Lett., 83, 3328 [NASA ADS] [CrossRef] [Google Scholar]
27. Hughes, A. L. H., & Armitage, P. J. 2012, MNRAS, 423, 389 [NASA ADS] [CrossRef] [Google Scholar]
28. Johansen, A., Henning, T., & Klahr, H. 2006, ApJ, 643, 1219 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
29. Johansen, A., Oishi, J. S., Mac Low, M.-M., et al. 2007, Nature, 448, 1022 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
30. Johansen, A., Youdin, A., & Mac Low, M.-M. 2009, ApJ, 704, L75 [NASA ADS] [CrossRef] [Google Scholar]
31. Johansen, A., Blum, J., Tanaka, H., et al. 2014, Protostars and Planets VI, 547 [Google Scholar]
32. Kataoka, A., Tanaka, H., Okuzumi, S., & Wada, K. 2013a, A&A, 557, L4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
33. Kataoka, A., Tanaka, H., Okuzumi, S., & Wada, K. 2013b, A&A, 554, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
34. Kempf, S., Pfalzner, S., & Henning, T. K. 1999, Icarus, 141, 388 [NASA ADS] [CrossRef] [Google Scholar]
35. Kobayashi, H., Ormel, C. W., & Ida, S. 2012, ApJ, 756, 70 [NASA ADS] [CrossRef] [Google Scholar]
36. Kretke, K. A., & Levison, H. F. 2014, AJ, 148, 109 [NASA ADS] [CrossRef] [Google Scholar]
37. Krijt, S., Dominik, C., & Tielens, A. G. G. M. 2014, J. Phys. D Appl. Phys., 47, 175302 [NASA ADS] [CrossRef] [Google Scholar]
38. Krijt, S., Ormel, C. W., Dominik, C., & Tielens, A. G. G. M. 2015, A&A, 574, A83 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
39. Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
40. Lambrechts, M., & Johansen, A. 2014, A&A, 572, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
41. Levison, H. F., Kretke, K. A., & Duncan, M. J. 2015, Nature, 524, 322 [NASA ADS] [CrossRef] [Google Scholar]
42. Min, M., Dullemond, C. P., Kama, M., & Dominik, C. 2011, Icarus, 212, 416 [NASA ADS] [CrossRef] [Google Scholar]
43. Morbidelli, A., Levison, H. F., Tsiganis, K., & Gomes, R. 2005, Nature, 435, 462 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
44. Nakagawa, Y., Sekiya, M., & Hayashi, C. 1986, Icarus, 67, 375 [NASA ADS] [CrossRef] [Google Scholar]
45. Okuzumi, S., Tanaka, H., & Sakagami, M.-A. 2009, ApJ, 707, 1247 [NASA ADS] [CrossRef] [Google Scholar]
46. Okuzumi, S., Tanaka, H., Kobayashi, H., & Wada, K. 2012, ApJ, 752, 106 [NASA ADS] [CrossRef] [Google Scholar]
47. Okuzumi, S., Momose, M., Sirono, S.-I., Kobayashi, H., & Tanaka, H. 2015, ApJ, submitted [arXiv:1510.03556] [Google Scholar]
48. Ormel, C. W., & Cuzzi, J. N. 2007, A&A, 466, 413 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
49. Ormel, C. W., & Klahr, H. H. 2010, A&A, 520, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
50. Ormel, C. W., Spaans, M., & Tielens, A. G. G. M. 2007, A&A, 461, 215 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
51. Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
52. Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62 [NASA ADS] [CrossRef] [Google Scholar]
53. Sato, T., Okuzumi, S., & Ida, S. 2015, A&A, submitted [arXiv:1512.02414] [Google Scholar]
54. Seizinger, A., & Kley, W. 2013, A&A, 551, A65 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
55. Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
56. Sirono, S.-I. 2011, ApJ, 735, 131 [NASA ADS] [CrossRef] [Google Scholar]
57. Suyama, T., Wada, K., & Tanaka, H. 2008, ApJ, 684, 1310 [NASA ADS] [CrossRef] [Google Scholar]
58. Suyama, T., Wada, K., Tanaka, H., & Okuzumi, S. 2012, ApJ, 753, 115 [NASA ADS] [CrossRef] [Google Scholar]
59. Takeuchi, T., Muto, T., Okuzumi, S., Ishitsu, N., & Ida, S. 2012, ApJ, 744, 101 [NASA ADS] [CrossRef] [Google Scholar]
60. Testi, L., Birnstiel, T., Ricci, L., et al. 2014, Protostars and Planets VI, 339 [Google Scholar]
61. Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Nature, 435, 459 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
62. Visser, R. G., & Ormel, C. W. 2015, A&A, accepted [arXiv:1511.03903] [Google Scholar]
63. Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2009, ApJ, 702, 1490 [NASA ADS] [CrossRef] [Google Scholar]
64. Wada, K., Tanaka, H., Suyama, T., Kimura, H., & Yamamoto, T. 2011, ApJ, 737, 36 [NASA ADS] [CrossRef] [Google Scholar]
65. Wada, K., Tanaka, H., Okuzumi, S., et al. 2013, A&A, 559, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
66. Weidenschilling, S. J. 1977, MNRAS, 180, 57 [NASA ADS] [CrossRef] [Google Scholar]
67. Weidenschilling, S. J. 1980, Icarus, 44, 172 [NASA ADS] [CrossRef] [Google Scholar]
68. Weidenschilling, S. J. 1995, Icarus, 116, 433 [NASA ADS] [CrossRef] [Google Scholar]
69. Whipple, F. L. 1972, in From Plasma to Planet, ed. A. Elvius (New-York: Wiley-Interscience), 211 [Google Scholar]
70. Windmark, F., Birnstiel, T., Güttler, C., et al. 2012, A&A, 540, A73 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
71. Youdin, A. N., & Chiang, E. I. 2004, ApJ, 601, 1109 [NASA ADS] [CrossRef] [Google Scholar]
72. Youdin, A. N., & Goodman, J. 2005, ApJ, 620, 459 [NASA ADS] [CrossRef] [Google Scholar]
73. Youdin, A. N., & Lithwick, Y. 2007, Icarus, 192, 588 [NASA ADS] [CrossRef] [Google Scholar]
74. Youdin, A. N., & Shu, F. H. 2002, ApJ, 580, 494 [NASA ADS] [CrossRef] [Google Scholar]
75. Zsom, A., & Dullemond, C. P. 2008, A&A, 489, 931 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
76. Zsom, A., Ormel, C. W., Güttler, C., Blum, J., & Dullemond, C. P. 2010, A&A, 513, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

Appendix A: Particle stopping time

The particle stopping time is a function of an aggregate’s mass m and size a, and the properties of the local gas. Depending on the aggregate size in relation to a gas molecule mean free path λmfp, the stopping time is given by the Epstein or Stokes drag regime through (A.1)with the mean thermal velocity of the gas molecules. The mean free path depends on the gas density and is given by λmfp = mg/ (σmolρg), with σmol = 2 × 10-15 cm2 the molecular cross section. For porous aggregates, the cross section A equals the orientation-averaged projected cross-section Okuzumi et al. (2009).

Equation (A.1) is valid when the particle Reynolds number Rep = 4avdg/ (vthλmfp) < 1, with vdg the relative velocity between the gas and the dust particle. For the largest bodies however, this condition is often not met. In these cases, it is useful to write the stopping time as (A.2)with CD the drag coefficient. Following Weidenschilling (1977), we use (A.3)For large Reynolds numbers, the stopping time depends the velocity relative to the gas, and one has to iterate to obtain ts.

Appendix B: Particle relative velocity

The relative velocity between particle 1 and particle 2 is obtained by adding various velocity sources quadratically. We take into account relative velocities arising from Brownian motion, turbulence, radial drift, and azimuthal drift. The Brownian motion relative velocity is given by (B.1)and depends only on particle masses and temperature.

The turbulence-induced relative velocity between two particles with stopping times ts,1 and ts,2ts,1 has three regimes (Ormel & Cuzzi 2007) (B.2)where δvg = α1/2cs is the mean random velocity of the largest turbulent eddies, and tη = Ret1/2tL is the turn-over time of the smallest eddies. The turbulence Reynolds number is Ret=, with the molecular viscosity νmol = vthλmfp/ 2.

The relative velocity from radial drift just equals Δvr = | vdriftts,1)−vdriftts,2) |, with the drift velocity given by Eq. (8). The azimuthal relative velocity is obtained in a similar way, as Δvφ = | vφts,1)−vφts,2) | with (B.3)Finally, the total relative velocity is given by (B.4)The particle stopping times and relative velocities are calculated in the midplane of the gas disk, as this is where most of the coagulation occurs.

Appendix C: Particle porosity evolution

 Fig. C.1 Particle porosity as a function of aggregate mass at various locations in the benchmark disk model (see Table 1), assuming 0.1 μm monomers. When calculating φ(m,r), we assume the aggregates grow though hit-and-stick collisions, followed by gas- and self-gravity compaction (see text). Gas compaction is more efficient at smaller radii. Open with DEXTER

Initially, particles grow through low-energy hit-and-stick collisions. During this growth phase, the fractal dimension is ~ 2 (Kempf et al. 1999), and the porosity is given by (C.1)where m is the monomer mass. The fractal growth regime ends when collisions become energetic enough for compaction (Dominik & Tielens 1997) or when gas ram pressure compaction becomes effective (Kataoka et al. 2013b).

For low internal densities, Kataoka et al. (2013b) found that the external pressure a dust aggregate can just withstand equals (C.2)An important parameter is the critical rolling energy Eroll (Dominik & Tielens 1997). Based on experimental investigations (Heim et al. 1999; Gundlach et al. 2011) and theoretical work (Krijt et al. 2014) we obtain (C.3)The pressure of Eq. (C.2) can then be compared to the pressure arising form the surrounding gas and from self-gravity (Kataoka et al. 2013a) (C.4)with \begin{lxirformule}$G$\end{lxirformule} the gravitational constant.

In our semi-analytical model, the porosity is determined by fractal growth followed by gas- and eventually self-gravity compaction. Specifically, we make use of Eq. (C.1), until the external pressures become too large for the aggregate to withstand. At that point, we obtain the porosity by setting Eq. (C.2) equal to Eq. (C.4). Of these processes, only gas pressure compaction varies with disk location, being more efficient in at smaller disk radii.

As particles only move inward in the cases considered in this work (Fig. 4), we can construct a single function φ(m,r) that is independent of the dust particle’s growth- and drift history. More specifically, for our assumptions (i.e., fractal growth followed by gas- and self-gravity compaction), φ(m,r1) ≥ φ(m,r2>r1). Figure C.1 shows porosity as a function of mass for three different radii, assuming 0.1 μm monomers and the disk properties of Table 1.

All Tables

Table 1

Benchmark disk model used throughout this work.

All Figures

 Fig. 1 Growth timescales of highly porous aggregates in an MMSN-like disk at 5 AU with a turbulence of α = 10-3. The solid black line is obtained with Monte Carlo models that take into account the full particle size and porosity distributions (Krijt et al. 2015) and the colored lines show Eq. (11) for three different values of κ (see text). Open with DEXTER In the text
 Fig. 2 Illustration of the batch method described in Sect. 3. A) Globally, the dust content of the protoplanetary disk is described by a collection of unique and independent batches. The trajectories that individual batches follow in space and time are dubbed lifelines. B) Each batch is composed of three separate legs. These legs can move closer together or further apart, but the mass enclosed between the inner two and outer two legs is always identical. For each leg, growth and radial drift are solved simultaneously using Eqs. (13) and (14), while the batch’s surface density and its first derivative (needed to calculate the growth timescales) are obtained using Eqs. (17) and (20). C) At every leg location, the local dust component is described by a monodisperse aggregate population with a single mass, size, and porosity. D) The aggregates themselves are composed of (sub)micron-sized monomers, whose properties profoundly affect the aggregate’s strength against fragmentation and compaction (see Appendix C). Open with DEXTER In the text
 Fig. 3 Evolution of the dust surface density in time. A) No growth case. Particles always have a = 1 mm. The dashed lines show analytical results of Youdin & Shu (2002). B) Compact growth with perfect sticking. Dust particles start out as 0.1 μm monomers and grow on timescales given by Eq. (11), while staying compact at all times (φ = 1). The gas surface density is shown by the dotted line. The results are qualitatively different, with the no-growth case resulting in pile-ups (locally increasing Σd/ Σg) and the compact coagulation model showing an inside-out removal of dust (decreasing Σd/ Σg everywhere). Open with DEXTER In the text
 Fig. 4 Locations of dust batches in time for four standard cases. A) No growth, particles are always 1 mm and compact. B) Compact growth starting from 0.1 μm particles. C) Porous growth with perfect sticking. Successful planetesimal formation (Ωts> 103) through coagulation is indicated by the (•)-symbol. The lifeline of the batch that formed planetesimals furthest out is colored red. D) Porous growth followed by erosion above velocities of 20 m s-1. Dashed lines indicate the particles in the batch have 10-2< Ωts< 3 and background colors indicate the midplane dust-to-gas ratio of particles that have Ωts ≤ 3. Open with DEXTER In the text
 Fig. 5 Influence of total disk mass on the location where planetesimals and pebbles form (Sect. 4.1). Material that starts out left of the solid blue line successfully grows through the drift barrier, eventually forming planetesimals a bit closer in. The red outermost planetesimal of Fig. 4C has been plotted for the MD/M⊙ = 10-2. Material that starts out in the gray region cannot coagulate fast enough, experience rapid radial drift, and at some point drift past the already-formed planetesimals in the inner disk. For reference, the current locations of the four Jovian planets are indicated by capital letters. Open with DEXTER In the text
 Fig. 6 Similar to Fig. 4D (porous growth + erosion), but with A) an increased metallicity of Z0 = 0.06; B) a reduced turbulence strength of α = 3 × 10-4; and C) a gas temperature 50% lower than that in Eq. (3). All three modifications result in a higher maximum dust-to-gas ratio compared to Fig. 4D. Open with DEXTER In the text
 Fig. C.1 Particle porosity as a function of aggregate mass at various locations in the benchmark disk model (see Table 1), assuming 0.1 μm monomers. When calculating φ(m,r), we assume the aggregates grow though hit-and-stick collisions, followed by gas- and self-gravity compaction (see text). Gas compaction is more efficient at smaller radii. Open with DEXTER In the text

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.