Research ArticleASTROPHYSICS

Growth of asteroids, planetary embryos, and Kuiper belt objects by chondrule accretion

See allHide authors and affiliations

Science Advances  17 Apr 2015:
Vol. 1, no. 3, e1500109
DOI: 10.1126/sciadv.1500109


Chondrules are millimeter-sized spherules that dominate primitive meteorites (chondrites) originating from the asteroid belt. The incorporation of chondrules into asteroidal bodies must be an important step in planet formation, but the mechanism is not understood. We show that the main growth of asteroids can result from gas drag–assisted accretion of chondrules. The largest planetesimals of a population with a characteristic radius of 100 km undergo runaway accretion of chondrules within ~3 My, forming planetary embryos up to Mars’s size along with smaller asteroids whose size distribution matches that of main belt asteroids. The aerodynamical accretion leads to size sorting of chondrules consistent with chondrites. Accretion of millimeter-sized chondrules and ice particles drives the growth of planetesimals beyond the ice line as well, but the growth time increases above the disc lifetime outside of 25 AU. The contribution of direct planetesimal accretion to the growth of both asteroids and Kuiper belt objects is minor. In contrast, planetesimal accretion and chondrule accretion play more equal roles in the formation of Moon-sized embryos in the terrestrial planet formation region. These embryos are isolated from each other and accrete planetesimals only at a low rate. However, the continued accretion of chondrules destabilizes the oligarchic configuration and leads to the formation of Mars-sized embryos and terrestrial planets by a combination of direct chondrule accretion and giant impacts.

  • minor planets
  • asteroids: general
  • methods: numerical
  • hydrodynamics
  • planets and satellites: formation
  • turbulence
  • protoplanetary disks


The formation of planetesimals and planetary embryos is an important step toward the assembly of planetary systems (1, 2). The asteroid belt contains planetesimals left over from the formation of the solar system and thus provides a record of the early stages of planet formation. Chondrite meteorites are fragments of asteroids that did not undergo large-scale melting and differentiation. Their constituent particles allow us to study the first stage of planet formation, in which micrometer-sized dust grains grow to asteroid-scale planetesimals. Chondrules dominate the mass budget of most chondrite meteorites; they are millimeter-sized glassy spherules formed by transient heating events in the protoplanetary disc (3). The small matrix grains that fill the space between these chondrules likely entered the meteorite parent bodies as accretionary rims on the chondrules (4).

Uranium-corrected Pb-Pb dates of chondrules show crystallization ages ranging from the earliest epoch of the solar protoplanetary disc—defined by the condensation of calcium-aluminum–rich inclusions (CAIs) 4567.30 ± 0.16 My ago—to about 3 My later (5). Moreover, chondrules from individual chondrites show variability in their 54Cr/52Cr and 50Ti/47Ti ratios (5, 6) that track genetic relationships between early formed solids and their respective reservoirs. Collectively, these observations indicate that chondrules from individual chondrite groups formed in different regions of the protoplanetary disc and were subsequently transported to the accretion regions of their respective parent bodies. Thus, the region of asteroid formation must have been dominated by chondrules over time scales comparable to the observed lifetimes of protoplanetary discs around young stars (7).

Meteorites provide a direct view of the primitive, chondritic material from which the planetesimals in the asteroid belt formed. In contrast, planetesimals originally located in the formation region of terrestrial planets accreted into planetary bodies, and hence little evidence exists with respect to the material from which those planetesimals formed. However, the abundance of chondrules in the asteroid belt suggests that chondrules could have been widespread in the terrestrial planet formation region as well. Indeed, some chondrules record 54Cr/52Cr and 50Ti/47Ti ratios indicating formation in the accretion regions of Earth and Mars (5, 6). Hence, understanding the role of chondrules in the formation of asteroids can provide critical insights into how the terrestrial planets formed closer to the Sun.

The mechanism by which large amounts of chondrules were incorporated into asteroids is poorly understood. Particles of chondrule sizes can concentrate near the smallest scales of the turbulent gas (8). The most extreme concentration events have been proposed to lead to asteroid formation by gravitational contraction of such intermittent chondrule clouds (9). However, it is still debated whether such high concentrations actually occur (10). The streaming instability (11) is an alternative mechanism that concentrates larger particles (12, 13), with typical sizes of 10 cm to 1 m when applied to the asteroid belt, because of an aerodynamical effect wherein particles pile up in dense filaments, which can reach densities of more than 1000 times the local gas density (14). The streaming instability nevertheless fails to concentrate chondrule-sized particles whose motion is strongly dampened by gas drag (15).


Here, we report the results of a model where asteroids and planetary embryos grow by accretion of chondrules onto planetesimal seeds formed by the streaming instability. Particles of large enough size to be concentrated by streaming instabilities could have been present in the early stages of the protoplanetary disc when planetesimals formed. Such particles include macrochondrules (16), chondrule aggregates (17), and ice-rich pebbles. All of these drift rapidly toward the Sun because of gas drag. Therefore, we assume that centimeter-sized particles were only present in the earliest stages of the protoplanetary disc, whereas smaller chondrules existed during most of the disc’s lifetime. Chondrule-sized particles may have been able to remain in the asteroid belt because of stellar outflows (18) or advection with the outward moving gas in the mid-plane layer of sedimented particles (15). We further discuss the conditions for planetesimal formation in the asteroid belt, particularly the formation of decimeter-sized particles, in the Supplementary Text.

Planetesimal formation simulations

We use high-resolution models of planetesimal formation through streaming instabilities to set the initial conditions for our model. The highest resolution reached before this work in such simulations is 1283 grid cells with 2.4 million superparticles representing the solids (14). We report on simulations with up to 5123 grid cells and 153.6 million superparticles representing the solids, which is 64 times more resolution elements than in previous simulations. Details of the simulation method are given in Materials and Methods.

We first perform a convergence test of streaming instability simulations with no self-gravity between the particles. Figure 1 shows the maximum particle density versus time in these simulations. We run the initial 30 orbits with the full gas density. Here, the maximum particle density is very moderate, and no strong concentrations occur. Between 30 and 40 orbits, we remove 50% of the gas, effectively increasing the solids-to-gas ratio from Z = 0.01 to Z = 0.02. This increase in metallicity triggers strong particle concentration through the streaming instability of up to 10,000 times the local gas density. The maximum particle density increases about quadratically with the inverse size of the resolution element δx. The maximum particle density on different scales is shown in Fig. 2. Here, we have measured the maximum density in spheres of diameters from δx up to the full size of the simulation box. The results are well converged from 643 to 2563 at the scales that are present at all resolutions. The smaller scales made accessible at higher resolution yield increasingly higher densities. This is an effect of the filamentary structure of the overdense particle filaments forming by the streaming instability; higher resolution allows us to resolve finer structure and hence higher densities. Smaller planetesimals can thus form as the resolution is increased, and increasingly smaller-scale filaments cross the Roche density.

Fig. 1 The maximum particle density versus time for streaming instability simulations without self-gravity at resolutions 643 (black line), 1283 (yellow line), and 2563 (blue line).

Particle density is measured in units of the mid-plane gas density ρg and time in units of the orbital period Torb. The first 30 orbits are run with a solids-to-gas ratio of Z = 0.01, with only modest overdensities seen in the particle density. Half of the gas is then removed over the next 10 orbits, triggering strong particle concentration, of up to 10,000 times the local gas density at the highest resolution. Doubling the resolution leads to a more than a quadrupling of the maximum particle density.

Fig. 2 The maximum particle density versus the length scale (measured in gas scale heights, H).

We have taken spheres with diameters from one grid cell up to the largest scale of the simulation and noted the maximum value of the density over all simulation snapshots. The results are relatively converged on each scale. The increase in maximum particle density with increasing resolution is an effect of resolving ever smaller-scale filaments. Blue dotted lines mark the Roche density at three values of the gas column density at 2.5 AU relative to the minimum mass solar nebula (MMSN). The red dotted line indicates the characteristic length scale of the streaming instability. The black line shows a power law of slope −2, which shows that the maximum density follows approximately the inverse square of the length scale.

In Fig. 3, we show the size distributions of planetesimals forming in streaming instability simulations that include the self-gravity of the particles. The planetesimal size distribution is dominated in number by small planetesimals, but most of the mass is in the few largest bodies. We find that progressively smaller planetesimals form, alongside the large planetesimals, at higher resolution. The largest planetesimals, which contain the dominant mass of the population, decrease to about 100 km in radius for the column density of the Minimum Mass Solar Nebula at the location of the asteroid belt (19). The size distribution of asteroids above 60 km in radius has been shown to retain its primordial shape because the depletion of the asteroid belt happens in a size-independent fashion for those large bodies (20, 21). The differential size distribution resulting from streaming instabilities disagrees with the observed size distribution of the asteroids (dashed line in Fig. 3).

Fig. 3 Birth size distribution of planetesimals forming by the streaming instability in 25-cm-sized particles.

The differential number of planetesimals (per 1022 g) is calculated with respect to the nearest size neighbors in the simulation. Yellow, blue, and red circles indicate individual planetesimals forming in computer simulations at 1283, 2563, and 5123 grid cells, respectively. The size distribution of the highest-resolution simulation is fitted with a power law dN/dMMqM (red line; the fit is based on the cumulative size distribution shown in Fig. 4, but does not include the exponential tapering). Simulations with lower values of the particle column density Σp yield successively smaller radii for the largest planetesimals, down to below 100 km in radius for a column density similar to the Minimum Mass Solar Nebula model (Σp = 4.3 g/cm2 at r = 2.5 AU). Binary planetesimals (marked with circles) appear only in the highest-resolution simulation because binary survival requires sufficient resolution of the Hill radius. The differential size distribution of the asteroid belt (dashed line) shows characteristic bumps at 60 and 200 km. The planetesimal birth sizes from the simulations are clearly not in agreement with main belt asteroids.

Chondrule accretion on asteroids

We proceed to consider the subsequent evolution of the size distribution as the newly formed planetesimals accrete chondrules present in the gaseous surroundings. We assume that planetesimals are born within an ocean of small chondrules that do not directly participate in planetesimal formation. Accretion of these chondrules by the planetesimals over the following millions of years can lead to significant mass growth. We use initial size distributions inspired by the streaming instability simulations. However, chondrule accretion is not limited to this model; it would be efficient in the context of any planetesimal formation mechanism that produces planetesimals with characteristic sizes of 100 km (9, 22).

The gas in the protoplanetary disc is slightly pressure-supported in the radial direction and hence moves at a sub-Keplerian speed, with Δv marking the difference between the Keplerian speed and the actual gas speed. The speed difference is about 53 m/s in the optically thin Minimum Mass Solar Nebula model (19). The Bondi radius of a planetesimal with radius R and internal density ρ, as shown below,Embedded Image(1)measures the gravitational deflection radius of the planetesimal. Chondrules embedded in the sub-Keplerian gas flow are accreted from the Bondi radius (23, 24) when their friction time tf is comparable to the time to drift azimuthally across the Bondi radius, tB = RBv. An additional turbulent component to the gas motion can be ignored because this is expected to be much slower than the sub-Keplerian speed in the asteroid formation region. Peak accretion rates are obtained for tf/tB in the range from 0.5 to 10 (24). This happens at the orbital distance of the asteroid belt for particle sizes in the interval, as shown belowEmbedded Image(2)

Here, fgas represents the ratio of the actual gas column density to that of the Minimum Mass Solar Nebula (19). Chondrule accretion happens at a rate proportional to Embedded Image (or Embedded ImageR6). This runaway accretion can drive a very steep differential size distribution similar to what is observed for large asteroids in the asteroid belt.

We solve numerically for the temporal evolution of the sizes and orbits of planetesimals accreting chondrules and colliding with each other in an annulus extending from 2.4 to 2.6 AU (see Materials and Methods for details). The planetesimals have initial radii from 10 to 150 km, distributed in size along an exponentially tapered power law based on the cumulative size distribution of planetesimals forming in the streaming instability simulations (see Fig. 4), dN/dRR−2.8 exp[−(R/Rexp)4], with exponential tapering at Rexp = 100 km and a total mass of 0.04 Earth masses. Chondrules have diameters from 0.02 to 1.6 mm, in broad agreement with observed chondrule sizes (3). The growth of particles and the formation of chondrule precursors are discussed in the Supplementary Materials. We place 0.2 Earth masses of chondrules initially in the annulus while another 0.2 Earth masses is created continuously over a time span of 3 My. The chondrules are given a Gaussian density distribution around the mid-plane, with the scale height set through a nominal value of the turbulent viscosity of α = 10−4. The gas in the protoplanetary disc is depleted exponentially on a time scale of 3 My (7).

Fig. 4 Cumulative birth size distribution of planetesimals forming by the streaming instability in our highest-resolution simulation (5123).

The cumulative size distribution is fitted with an exponentially tapered power law N>M−0.6 exp[−(M/Mexp)(4/3)] (red line), with exponential tapering at Mexp = 1.2 × 1023 g (dotted line). Shallower or steeper power laws yield poorer fits to the populations of small and large planetesimals, respectively. We choose to fit the cumulative size distribution rather than the differential size distribution to avoid the noise inherently present in the latter. The fit can be translated to dN/dMM−1.6 exp[−(M/Mexp)(4/3)] (the power law part of this fit is indicated in Fig. 3) as well as to dN/dRR−2.8 exp[−(R/Rexp)4] in the differential size distribution per unit radius.

We show in the left panel of Fig. 5 the size distribution after 5 My of accreting chondrules. The steep drop in the initial number of planetesimals larger than 100 km in radius has become much shallower because of chondrule accretion, forming a characteristic bump at 60- to 70-km radius. The size distribution of smaller asteroids in today’s belt was likely filled later with fragments from collisional grinding (20). The steep size distribution of larger asteroids ends at around 200 km as even larger asteroids would need to accrete chondrules larger than a millimeter in size (see Eq. 2). At this point, the largest asteroids enter a more democratic growth phase wherein the accretional cross sections are reduced by gas friction within the Bondi radius (23, 24). The agreement with the observed size distribution of asteroids is excellent in the nominal model, except in the range around Ceres mass where it predicts too few objects by a factor of about 2 to 3. Changing the model parameters (using larger chondrules or smaller planetesimal birth sizes) yields poorer agreement with observed asteroid sizes. This implies that the observed size distribution of asteroids is directly determined by the size distribution of the chondrules and the birth sizes of the planetesimals.

Fig. 5 The size distribution of asteroids and embryos after accreting chondrules for 5 My (left panel) and selected masses and inclinations as a function of time (right panel).

The nominal model (black line) closely matches the steep size distribution of main belt asteroids (red line, representing the current asteroid belt multiplied by a depletion factor of 540) from 60 to 200 km in radius. The size distribution becomes shallower above 200 km; this is also seen in the asteroid belt, although the simulations underproduce Ceres-sized planetesimals by a factor of about 2 to 3. A simulation with no chondrules (blue line) produces no asteroids larger than 300 km. Inclusion of chondrules up to centimeter sizes (pink line in insert) gives a much too low production of Ceres-sized asteroids, whereas setting the exponential cutoff radius of planetesimals to 50 km (green line) leads to a poorer match to the bump at 60 km. The right plot shows that the formation of the first embryos after 2.5 My quenches chondrule accretion on the smaller asteroids by exciting their inclinations i (right axis). The dotted lines indicate the mass contribution from planetesimal-planetesimal collisions. Asteroids and embryos larger than 200 km in radius owe at least two-thirds of their mass to chondrule accretion.

The right panel of Fig. 5 shows that the major growth phase of asteroids with final radii ranging from 200 to 500 km occurs after 2.5 My. Beyond this time, the largest planetesimals in the population grow to become planetary embryos with sizes similar to the Moon and Mars. These growing embryos excite the inclinations of the smaller asteroids to i ~ 0.01, which disconnects the asteroids from the chondrules in the mid-plane layer, quenching their accretion of further chondrules. The beginning of embryo formation terminates the accretion of asteroids after 3 My, thus defining the final sizes of the asteroid belt population. The absence of such planetary embryos in today’s asteroid belt likely reflects a later depletion by gravitational perturbations from Jupiter (2527). Depletion of the asteroid belt is discussed further in the Supplementary Text. An important implication of chondrule accretion is that accretion of other planetesimals contributes only a minor fraction to the final masses of large asteroids and embryos (right panel of Fig. 5).

Size sorting of chondrules

Chondrules in chondrites appear strongly size-sorted (8, 28), with the average chondrule diameter varying among the ordinary chondrites (28) from 0.3 mm (H chondrites) up to 0.5 mm (L and LL chondrites). Carbonaceous chondrites exhibit a larger range in chondrule sizes, from 0.1 to 1 mm. In Fig. 6, we show that the mean diameter of accreted chondrules in our model lies in a range similar to that observed for chondrites. The size distribution of accreted chondrules is very narrow, also in agreement with the observed size distribution in the ordinary chondrites (28). The chondrule accretion process, which leads to aerodynamical sorting of chondrules, may thus represent the underlying physical mechanism for the size sorting of chondrules observed in chondrite meteorites. Although aerodynamical sorting has been previously proposed to arise from the gas flow around asteroids (29), our simulations show that accretion from the full Bondi radius is much more efficient in achieving aerodynamical sorting of chondrules.

Fig. 6 Mean chondrule sizes (〈d〉, upper panel) as a function of layer radius R and size distribution width (σφ, lower panel) as a function of mean chondrule size.

Yellow lines in the upper panel indicate the chondrule size evolution in individual asteroids and embryos, whereas red lines indicate mean accreted chondrule sizes at different times. The accreted chondrule size increases approximately linearly with planetesimal size at t = 0 My. Asteroids stirred by the growing embryos over the next 2 My accrete increasingly larger chondrules because asteroid velocities align with the sub-Keplerian chondrule flow at aphelion. Finally, asteroids accrete surface layers of mainly very small chondrules, down to below 0.1 mm in diameter, at late times when their large inclinations decouple the asteroid orbits from the large chondrules in the mid-plane. The width of the chondrule size distribution in the lower panel is given in terms of σφ, the base-2-logarithmic half-width of the cumulative mass distribution of chondrules (σφ = 1 implies that two-thirds of the chondrules have diameters within a factor 21 = 2 from the mean). Dots indicate chondrule layers inside asteroids and embryos. Chondrules are aerodynamically sorted by the accretion process to values of σφ in agreement with the chondrules found in ordinary chondrites (hashed region). Unfiltered accretion from the background size distribution of chondrules (blue dot, size distribution of unsedimented particles; yellow dot, size distribution in the mid-plane) yields specific pairs of 〈d〉 and σφ that are not consistent with constraints from ordinary chondrites.

Terrestrial planet formation with chondrules

Our identification of chondrule accretion as driving planetesimal growth in the asteroid belt implies that chondrules also could play an important role in terrestrial planet formation. To test this hypothesis, we have performed a numerical integration of the evolution of planetesimal orbits and sizes at 1 AU (Fig. 7). The left panel shows the size distribution of the planetesimals at 0, 1, 3, and 5 My. Growth to super-Ceres–sized planetesimals is rapid and happens within 1 My. This growth is driven mainly by planetesimal accretion, in stark contrast to the situation at 2.5 AU. This is seen in the right panel of Fig. 7, where the mass of the largest body is shown as a function of time (full line) together with the mass contribution from chondrules (dashed line) and from planetesimals (dash-dotted line). However, chondrule accretion becomes dominant after 1.5 My and drives further growth up to Mars-mass embryos after 4 My. The largest body experiences a giant impact after 4 My, after which it continues to grow by chondrule accretion toward 0.9 Earth masses.

Fig. 7 Growth of embryos and terrestrial planets at 1 AU.

Left panel shows the size distribution at four different times, and right panel shows the mass of the most massive body in the simulation as a function of time. The growth up to 1000-km-sized embryos is mainly driven by planetesimal accretion because chondrule-sized pebbles are tightly coupled to the gas and hence hard to accrete. However, chondrule accretion gradually comes to dominate the accretion as the embryos grow. The largest body reaches Mars size after 3 My, with more than 90% contribution from chondrule accretion. A giant impact occurs just before 4 My, wherein the largest body accretes the third largest body in the population. The continued accretion of chondrules leads to the formation of an Earth-mass body after 5 My. A simulation with no chondrules evolves very differently: a large number of embryos form with masses just below the isolation mass of Miso ≈ 0.01 ME.

Chondrule accretion can thus promote the growth of the largest embryos from Moon masses toward Mars masses and finally Earth masses. The dominance of planetesimal accretion in the initial growth toward Moon masses occurs because chondrules couple tightly to the gas at 1 AU. Here, the gas density is more than a factor of 10 higher than at 2.5 AU. This prevents sedimentation, such that all chondrule sizes are well mixed with the gas, and it reduces the cross section for accreting chondrules because tightly coupled particles cannot be accreted from the full Bondi radius. This situation changes as gas dissipates exponentially over 3 My, increasing the friction times and thus the accretion efficiency of the chondrules. Furthermore, the increasingly large embryos obtain high accretion radii for chondrules and hence high chondrule accretion rates, despite the relatively low degree of sedimentation of chondrules present at 1 AU orbits.

Chondrules also play a critical role in terrestrial planet formation in a different way, namely, by breaking the isolation mass configuration. Growth by pure planetesimal accretion tends to end in oligarchic growth wherein the largest embryos are isolated from each other by about 10 Hill radii (30). We have implemented the effect of this isolation tendency by identifying the group of isolated bodies as those that can fit their combined reach of 10 Hill radii into the annulus of 0.2 AU in diameter (31, 32). Isolated bodies are not allowed to accrete each other and only affect each other dynamically via distant viscous stirring in the eccentricity. Inclination perturbations, as well as dynamical friction in the eccentricity, are ignored between isolated bodies.

We mark the isolation mass of about 0.01 ME in the right plot of Fig. 7 (blue line). A simulation performed with no chondrule accretion (red line) shows the expected growth of the largest embryos to just below the isolation mass. The growth curve for pure planetesimal accretion follows closely the growth curve for the full simulation until about 10−4 Earth masses. At this point, chondrule accretion becomes significant and the two curves diverge. Chondrule accretion can destabilize a set of isolated bodies because their continued growth by chondrule accretion drives the least massive of the isolated bodies out of isolation. Hence, the giant impacts experienced by the largest embryo between 1.5 and 4 My are all driven by chondrule accretion because these impacts happen beyond the isolation mass.

The importance of chondrule accretion increases if chondrules in the terrestrial planet formation region are larger, namely, up to 1 cm (macrochondrules). We show the results of expanding the size distribution of chondrules up to 1 cm in Fig. 8. Here, the initial growth has equal contribution from chondrules and planetesimals. The growth toward embryos and terrestrial planets is much more rapid than in Fig. 7; a planet with the mass of the Earth already forms after 1.5 My.

Fig. 8 Growth of embryos and terrestrial planets at 1 AU, with chondrule sizes up to 1 cm.

The increase in chondrule size compared to the previous figure enhances the chondrule accretion rate substantially. The initial growth to 1000-km-sized embryos now has approximately equal contribution from planetesimal accretion and chondrule accretion. An Earth-mass terrestrial planet forms already after 2 My, driven by a combination of chondrule accretion and giant impacts.

The final size of embryos and planets forming by combined planetesimal accretion and chondrule accretion depends on the exact timing of disc dissipation. If the inner disc already dissipates after 3 My, then the terrestrial planet formation region will be dominated by a number of Mars-sized embryos. These bodies can go on to partake in a traditional terrestrial planet formation scenario of gradual orbital perturbations and giant impacts over the next 30 to 100 My (33). In contrast, later disc dissipation or the presence of large chondrules in the terrestrial planet formation region allows for the direct formation of planetary bodies within a few million years. This is far more rapid than traditional terrestrial planet formation scenarios that do not consider chondrule accretion.

Pebble accretion in the Kuiper belt

The size distribution of trans-Neptunian objects has a characteristic bump around 50 km in radius (34), similar to main belt asteroids, followed by a steep decline toward larger sizes (20). We therefore explore here whether accretion of chondrules and icy pebbles could be responsible for the observed size distribution of Kuiper belt objects as well.

The asteroid belt displays a value of q ≈ 4.5 in the size distribution dN/dRRq of asteroids larger than 60 km in radius, whereas the hot population of the classical Kuiper belt and the Neptune Trojans have a much higher value of q ≈ 5…6 for large planetesimals (35). The cold population of the classical Kuiper belt has an even steeper decline, with q ≈ 8…9. The planetesimals in the cold population are believed to have formed in situ in the outermost regions of the solar system (36), whereas the hot population and scattered disc objects formed outside the birth location of Neptune and were subsequently pushed outward by Neptune’s migration to its current orbit (37, 38). The planetesimals that formed between the orbits of Jupiter and Neptune may today be represented by the C-type asteroids, scattered to the asteroid belt during Jupiter’s migration (27). To probe planetesimal growth by chondrule accretion (or icy particles of a similar size range) beyond the asteroid belt, we therefore consider two characteristic distances, 10 and 25 AU.

Planetesimal growth at 10 AU. The region around 10 AU represents the conditions in which the gas giants formed. The formation of Jupiter and Saturn must have had an enormous effect on the planetesimals in that region, with some ejected to the Oort cloud and others potentially mixed into the outer asteroid belt during Jupiter’s and Saturn’s migration (27). To test how the size distribution of planetesimals evolves at 10 AU, we set up an annulus of width 0.2 AU. The planetesimals are again given sizes from 10 to 150 km in radius, with a steep exponential tapering above 100 km. The internal density of the planetesimals is set to 2 g/cm3, similar to that of Ceres, a typical representative of the icy asteroids in the outer part of the main belt. The pebbles have radii from 0.01 to 1 mm. These pebbles could represent a mixture of chondrules transported outward from the inner solar system and icy pebbles formed beyond of the ice line. Chondrules may also have formed in situ in the outer parts of the protoplanetary disc. Some chondrules present in carbonaceous chondrites have 50Ti/47Ti ratios comparable to those of the most pristine chondrites, namely, the CI carbonaceous chondrites (6). This suggests that at least a fraction of chondrules in carbonaceous chondrites formed from pristine, thermally unprocessed precursor material typical of the accretion regions of CI chondrites. Given that these chondrites are believed to have accreted beyond the ice line (27), chondrule formation appears to not have been limited to the asteroidal belt, but also occurred in the outer solar system.

The results at 10 AU are shown in Fig. 9. Planetesimals of Ceres size grow within about 2 My, followed by the usual phase of embryo growth also seen in Fig. 5. The largest embryo reaches a mass comparable to the Earth’s before 3 My. The mid-plane density of pebbles is much lower at 10 AU than at 2.5 AU, but this is more than counteracted by increased sedimentation of particles in the dilute gas. The size distribution is very steep, much steeper than in the asteroid belt, and better resembles the steepness of the size distribution of the trans-Neptunian populations (35). Finally, it is interesting to note that, at these orbital distances, planetesimal-planetesimal collisions are even less important than at 2.5 AU (right panel of Fig. 9). The Earth-mass protoplanet has only a few parts in a thousand contributions from planetesimal collisions.

Fig. 9 Growth of icy planetesimals at 10 AU.

The growth rate by pebble accretion is as high as in the asteroid belt because the lower column density of pebbles is counteracted by the increased sedimentation in the more dilute gas. The annulus of 0.2-AU width produces in the end an Earth-sized protoplanet and a single Moon-sized embryo. Pebble accretion overwhelmingly dominates the growth (right panel). The icy protoplanet that forms has only a few parts in a thousand mass contribution from collisions. Large Ceres-sized planetesimals have a contribution from collisions of less than 5%, whereas the 200-km planetesimal owes about 1/10 of its growth from 130 km to planetesimal collisions. Note how the Ceres-sized planetesimal (blue line in the right panel) got a head start for efficient accretion of pebbles by experiencing a significant collision after 700,000 years.

Planetesimal growth at 25 AU. In the outer regions of the protoplanetary disc, we adopt a model where the surface density profile is proportional to the inverse of the orbital radius. This is in agreement with surface density profiles of observed protoplanetary discs (39, 40). Hence, we set the total column density of particles at r = 25 AU to Σp = 0.54 g/cm2. We set the temperature at 25 AU according to full radiative transfer models of protoplanetary discs. This is important because the optically thin assumption in the Minimum Mass Solar Nebula model of Hayashi (19) overestimates the temperature in the outer disc. The growth rate of a planetesimal with mass Mp by pebble accretion is given byEmbedded Image(3)

Here, fB denotes the ratio of the actual accretion radius to the Bondi radius (24), Δv is the sub-Keplerian speed difference, and ρp is the mid-plane particle density. The sub-Keplerian speed is given byEmbedded Image(4)

The speed difference scales with the sound speed as Embedded Image because H/r contains another cs. Hence, the accretion rate scales with the sound speed to the negative sixth power, and the sound speed must be set carefully to yield realistic results. We set the sound speed at 25 AU through the aspect ratio H/r = 0.05 on the basis of radiative transfer models of protoplanetary accretion discs (41). This is much smaller than the standard value in the Minimum Mass Solar Nebula (H/r = 0.072 at r = 25 AU) and leads to much higher accretion rates in the outer nebula.

We consider two models for planetesimal growth at 25 AU: first, a low density model where the internal density of the planetesimals is set to 0.5 g/cm3, similar to comets and some binary Kuiper belt objects (42); second, a high density model that has instead an internal density of 2 g/cm3, more in agreement with the icy dwarf planets Pluto and Ceres. The low density model experiences much lower pebble accretion rates because of the reduced gravity of the planetesimals. We find that a turbulent stirring of α = 10−4, which we used for the asteroid belt, results in growth times for pebble accretion of more than 10 My. Hence, we set the turbulent stirring to a low value of α = 10−6 in the low density model. The high density model has α = 10−5. We discuss in the next section how such low values could arise physically.

We show the resulting size distributions and growth curve of the largest body at 25 AU in Fig. 10. The size distributions in the left panel of the figure both display a steep power law beyond 100 km in radius. Planetesimals up to 300 km in radius accrete significant amounts of pebbles. Beyond 300 km, a runaway growth sets in wherein the largest body detaches from the remaining population (right panel of Fig. 10). This runaway growth is possible because of the weak gas friction at 25 AU. Hence, runaway pebble accretion is not slowed down by reaching the strong coupling branch but can instead continue to very large sizes (see the Supplementary Materials and figs. S6 and S7).

Fig. 10 Planetesimal growth at 25 AU.

Two models are considered: a low density model where the internal density is set to ρ = 0.5 g/cm3, similar to comets and binary Kuiper belt objects, and a high density model where the internal density is set to ρ = 2 g/cm3, similar to the dwarf planet Ceres. The low density model has turbulent stirring α = 10−6, whereas the high density model has α = 10−5. Both models display ordered growth up to 300-km radii, with a steep size distribution beyond 100-km sizes. This is followed by a runaway growth of a single, massive body. The right panel shows the size of the largest body as a function of time as well as the speed relative to a circular orbit. The runaway growth is facilitated by a steep decline in the eccentricity of the orbit because the high pebble accretion rate damps the eccentricity.

Altogether, it is clear that pebble accretion at 25 AU is much slower than at 2.5 AU, and actually requires very low turbulent activity for significant growth. However, planetesimal growth in such wide orbits apparently also has the potential to result in runaway accretion up to ice giant sizes. The size distribution of planetesimals from 100 to 300 km in radius is very steep, in good qualitative agreement with what is observed in the present Kuiper belt (35). Our results show that 25 AU can be considered the very limit to where pebble accretion is significant, and then only for weak turbulence. Beyond this radius, planetesimals maintain their birth sizes (and size distribution), and hence, the very steep size distribution of large planetesimals in the cold component of the classical Kuiper belt may reflect the actual exponential cutoff in the underlying birth size distribution.

In the solar system, Neptune stopped its outward migration at 30 AU. Therefore, the primordial Kuiper belt cannot have extended much beyond that distance because otherwise Neptune would have continued to migrate outward by scattering planetesimals (38). However, the presence of the cold population, which appears detached from Neptune, indicates that the Kuiper belt could simply have transitioned to a much lower surface density around 30 AU. The planetesimals there did not grow beyond their birth sizes, resulting in a jump in the planetesimal column density between the regions inside 25 AU that experienced growth by pebble accretion and the regions outside that did not.


Variation of turbulence strength

An important and relatively unconstrained parameter in our chondrule accretion model for the formation of asteroids and planetary embryos is the strength of the turbulence in the protoplanetary disc because the turbulent diffusion coefficient sets the scale height of the particle layer and the mid-plane density of the chondrules (1). The asteroid formation region is located in a region of the protoplanetary disc where the degree of ionization is believed to be too low to sustain turbulence driven by the magnetorotational instability. The ionization degree in the surface layers may nevertheless be high enough to drive accretion there. The turbulent motion in these layers will stir the otherwise laminar mid-plane. Effective turbulent viscosities between α = 10−5 and α = 10−4 have been measured in the mid-plane in computer simulations of dead zone stirring (43).

Models of protoplanetary accretion discs that also include the effect of ambipolar diffusion in the surface layers find that, even there, the growth of the magnetorotational instability is quenched by the lack of coupling between electrons and neutrals (44). The magnetic field may enter a configuration where gas (and angular momentum) is expelled upward, whereas gas connected to the magnetic field lines closer to the mid-plane is accreted toward the star. This disc wind accretion leaves the mid-plane completely laminar. The sedimentation of dust will nevertheless cause a mild stirring of the mid-plane (45). The stirring caused by streaming and Kelvin-Helmholtz instabilities in the dense mid-plane layer of chondrules will be much milder than the stirring from the active layers, with a typical value of α = 10−6 for chondrule-sized particles.

The nominal value for the turbulent diffusion coefficient in our model is α = 10−4, corresponding to conditions in a dead zone stirred by active surface layers. To quantify the effect of the turbulent stirring of chondrules on our results, we have run simulations with a lower (α = 2 × 10−6) and a higher (α = 10−3) value of the turbulent diffusion coefficient. The size of the largest planetesimal in the asteroid belt is shown as a function of time in Fig. 11. Clearly, a lower value of the turbulent diffusion coefficient results in a much higher growth rate of the planetesimals because chondrules are allowed to sediment much more strongly to the mid-plane. For α = 10−3, the formation of the first embryos is delayed to after 5 My, which is longer than the typical lifetimes of protoplanetary discs around young stars.

Fig. 11 The maximum planetesimal radius in the asteroid belt versus time for three different values of the turbulent viscosity α.

Here, α = 2 × 10−6 represents the strength of turbulence caused by streaming instabilities and Kelvin-Helmholtz instabilities in a sedimented mid-plane layer of chondrules, α = 10−4 represents the turbulence strength in a dead zone stirred by active surface layers, and α = 10−3 represents the turbulence strength caused directly by the magnetorotational instability. Turbulent stirring of chondrules sets the scale height and mid-plane density of the chondrule layer and hence dictates the planetesimal growth rate. The formation time of the first embryo depends strongly on the degree of stirring.

Layered accretion in the asteroid belt

Our model provides a direct connection between chondrules and asteroid growth. Whereas early accreted asteroids and their chondrules would melt from the energy released by the decay of the short-lived 26Al radionuclide, layers accreted later than 2 My may preserve their pristine, undifferentiated nature (46). Indeed, the existence of differentiated asteroids overlaid with chondritic crusts has been proposed to explain the systematic magnetizations of chondrules in the Allende meteorite (47, 48), imprinted from a geodynamo operating in the liquid central region of the parent body. The asteroid Lutetia has been suggested to be partially differentiated because of its high mean density and primitive surface (49).

The characteristic maximum size of planetesimals formed through streaming instabilities is 50 to 150 km in radius for column densities comparable to the Minimum Mass Solar Nebula (Fig. 3). Planetesimals much smaller than 100 km in radius are relatively inefficient at accreting chondrules and hence accrete only a thin veneer of chondrules on top of the primordially formed planetesimal (Figs. 5 and 6). Planetesimals forming before 1 to 2 My after CAIs will likely differentiate, as was the case for the parent bodies of the differentiated meteorites (50, 51). Hence, we predict that asteroids currently residing in the knee at the size distribution at 60-km radii represent primordial, differentiated bodies overlaid with a veneer of chondrules of thickness up to a few dozen kilometers.

Asteroids larger than 100 km in radius, on the other hand, grow mainly by accretion of chondrules. The initial planetesimal seeds will be covered with extended layers of narrowly size-sorted chondrules. The Allende meteorite, whose systematic magnetization makes the parent body a prime candidate for layered accretion, has been proposed to originate from the Eos family in the asteroid belt (52). This family is known to display significant spectral variations between its members (52), indicating partial differentiation and internal heterogeneity of the 110-km-radius parent body. Another example of an asteroid family with large internal variation in the spectral properties is the Eunomia family (48), believed to originate from a parent body of 150 km in radius. Intermediate-sized asteroids of 100 to 200 km in radius appear to be the best candidates for producing such families with large internal variation. Bodies of this size range have comparable mass fractions in the differentiated planetesimal seed and in the accreted chondrule layers.

Many other asteroid families display internal albedo distributions that are much narrower than the spread of albedos across families (53), in stark contrast to the heterogeneous Eos and Eunomia families discussed above. The inference of the internal structure of the parent bodies of such apparently homogeneous families is nevertheless complicated by the identification of interlopers in (and exclusion from) the asteroid families. A family member whose albedo is different from the family mean could in fact be mistaken for an interloper. However, the identification of asteroid families that appear to originate from internally homogeneous asteroids is not in conflict with our model. Our results show that many asteroids in the 50 to 100 km radius range will have experienced little chondrule accretion and represent primordial asteroid seeds. Additionally, very large asteroids of radii larger than 300 to 400 km must have undergone significant accretion of chondrules between 2 and 3 My after the formation of the seed planetesimal. Families formed mainly from the chondrule layers of such large asteroids would also appear relatively homogeneous. Fragments as large as 100 km in radius could consist of purely chondritic material, and families produced from those fragments would in turn appear entirely homogeneous.

Implications for terrestrial planet formation

The largest planetesimals in our simulations of the asteroid belt and the terrestrial planet formation region grow to sizes of several thousand kilometers, forming the embryos whose mutual collisions drive the subsequent terrestrial planet formation stage. With only 11% of Earth’s mass, Mars may be one of these remaining embryos. Mars is inferred to have accreted 2 ± 1 My after the first condensations in the solar system (54), in good agreement with the formation time scale of planetary embryos by chondrule accretion in our simulations. This implies that chondrule accretion is a driving mechanism ensuring the rapid growth of planetary embryos in the terrestrial planet formation region. Thus, the production of chondrules may be a key process promoting the assembly of rocky planets.


Planetesimal formation simulations

The planetesimal formation simulations were performed with the Pencil Code (which can be freely downloaded, including modifications done for this work, at A new sink particle algorithm developed for this project is presented in the Supplementary Materials.

Our simulations of planetesimal formation through streaming instabilities are performed in the local shearing box frame. The simulation box orbits with the local Keplerian frequency Ω at an arbitrary radial distance r from the central star. We fix the box size to a cube of side length L = 0.2 H, where H is the gas scale height. The box size is chosen to capture the relevant wavelengths of the streaming instability (55). The gas is initialized with a Gaussian density profile around the mid-plane, in reaction to the component of stellar gravity directed toward the mid-plane. The temperature is assumed to be a constant set through the sound speed cs, which is kept fixed during the simulation. Particles are given positions according to a Gaussian density distribution, with a scale height of 1% of the gas scale height. We assume that particles in the entire column have sedimented into the simulation domain, giving a mean solids-to-gas ratio of 〈ρp〉/ρ0 ≈ 0.25 for the chosen metallicity Z = 0.01. The particles have uniform sizes given by the Stokes number St = Ωtf = 0.3, where tf is the friction time of the particles, corresponding to about 25-cm-sized particles at 2.5 AU. It has previously been shown that the consideration of a range of particle sizes gives similar results to the monodisperse case (12).

The strength of the particle self-gravity in a scale-free local box simulation is set by the nondimensional parameterEmbedded Image(5)

Here, G is the gravity constant, ρ0 is the gas density in the mid-plane, and Ω is the Keplerian frequency at the chosen radial distance from the star. The scaling with semi-major axis is based on the properties of the Minimum Mass Solar Nebula (19). The parameter Γ multiplied by the nondimensional particle density ρp0 enters the Poisson equation for self-gravity, ∇2Φ = 4πGρp, when lengths are measured in units of the scale height and times in units of the inverse Keplerian frequency Ω−1. The parameter also converts a given mass density of particles in a grid cell to an actual mass because the choice of Γ defines the mid-plane gas density ρ0. Thus, both the dynamics and resulting planetesimal masses are strongly influenced by the value of Γ (this is also evident from the results presented in Fig. 3).

Chondrule accretion simulations

We have developed a statistical code that evolves the masses and orbits of planetesimals as they accrete chondrules and collide mutually. Details of the code algorithms and test calculations as well as comparisons to other published results are presented in Supplementary Materials and Methods.


Supplementary material for this article is available at

Materials and Methods


Fig. S1. The nondimensionless planetesimal masses μ = ρpx)3 (here ρp is the particle density represented by the planetesimal in a grid cell and δx is the size of the cell) and corresponding contracted radii as a function of time after self-gravity is turned on.

Fig. S2. The accretion radius Racc (that is, the impact parameter required for accretion, given in units of the Bondi radius RB) versus the particle friction time tf (in units of the Bondi time-scale tB), for planetesimal radii R between 10−4 RB and 104 RB.

Fig. S3. The evolution of eccentricity e and inclination i of 1000 planetesimals with mass M = 1024 g located at 1 AU with a surface density of 10 g/cm2, for three values of the time-step.

Fig. S4. The evolution of eccentricity e and inclination i of 800 planetesimals with mass M = 1024 g located at 1 AU with a surface density of 10 g/cm2.

Fig. S5. Phase points covering planetesimal orbits with inclinations i = 0, i = 0.0016, i = 0.0032 and i = 0.0064.

Fig. S6. The effect of eccentricity and inclination on the chondrule accretion rate, at an orbital distance of 2.5 AU.

Fig. S7. The accretion rate of planetesimals at 25 AU, as a function of the planetesimal size.

Fig. S8. The stratification integral (that is, the mean particle density over the accretion cross section) as a function of the accretion radius of the planetesimal.

Fig. S9. The number of remaining bodies versus time for the constant kernel test.

Fig. S10. Damping of eccentricities and inclinations of 10,000 1-cm–sized planetesimals located between 30 and 35 AU, emulating a test problem defined in Morbidelli et al. (21).

Fig. S11. Cumulative size distribution after 3 Myr of coagulation within a population of planetesimals of initial diameters 100 km (50 km in radius), as shown by the cross.

Fig. S12. Formation of planetesimals from centimeter-sized particles in a 2D shearing sheet simulation with the same spatial extent as the streaming instability simulations presented in the main text.

References (5677)

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: We are grateful to G. Stewart, A. Morbidelli, J. Cuzzi, and M. Lambrechts for stimulating discussions. Funding: A.J. is grateful for the financial support from the European Research Council (ERC Starting Grant 278675-PEBBLE2PLANET), the Knut and Alice Wallenberg Foundation, and the Swedish Research Council (grant 2010-3710). M.-M.M.L. was partly supported by NASA Origins of Solar Systems grant #NNX14AJ56G, NSF grant AST10-09802, and the Alexander von Humboldt-Stiftung. M.B. acknowledges funding from the Danish National Research Foundation (grant number DNRF97) and the European Research Council (ERC Consolidator Grant 616027-STARDUST2ASTEROIDS). This work was granted access to the HPC resources of Research Centre Jülich, Irish Centre for High-End Computing and Rechenzentrum Garching, made available within the Distributed European Computing Initiative by the PRACE-2IP (PRACE project “PLANETESIM”), receiving funding from the European Community’s Seventh Framework Programme (FP7/2007-2013) under grant agreement no. RI-283493. Author contributions: A.J. developed the chondrule accretion code and performed the numerical simulations. The manuscript was written by A.J. and M.B., and M.-M.M.L. and P.L. contributed to the discussions. Competing interests: The authors declare that they have no competing financial interests.
View Abstract

Stay Connected to Science Advances

Navigate This Article