Growing the seeds of pebble accretion through planetesimal accretion

We explore the growth of planetary embryos by planetesimal accretion up to and beyond the point at which pebble accretion becomes efﬁcient at the so-called Hill-transition mass. Both the transition mass and the characteristic mass of planetesimals that formed by the streaming instability increase with increasing distance from the star. We developed a model for the growth of a large planetesimal (embryo) embedded in a population of smaller planetesimals formed in a ﬁlament by the streaming instability. The model includes in a self-consistent way the collisional mass growth of the embryo, the fragmentation of the planetesimals, the velocity evolution of all involved bodies, and the viscous spreading of the ﬁlament. We ﬁnd that the embryo accretes all available material in the ﬁlament during the lifetime of the protoplanetary disc only in the inner regions of the disc. In contrast, we ﬁnd little or no growth in the outer parts of the disc beyond 5–10AU. Overall, our results demonstrate very long timescales for collisional growth of planetesimals in the regions of the protoplanetary disc in which giant planets form. This means that in order to form giant planets in cold orbits, pebble accretion must act directly on the largest bodies present in the initial mass function of planetesimals with little or no help from mutual collisions.


Introduction
The classic picture for the formation of planets requires growth over several orders of magnitude in size. Starting from micrometre-sized dust and ice grains, coagulation produces millimetre-sized pebbles. Collisional and dynamical processes, such as bouncing, fragmentation, and radial drift, limit the maximum particle size (Blum & Wurm 2008;Güttler et al. 2010;Zsom et al. 2010;Krijt et al. 2015), and the formation of larger bodies is effectively prevented. Porosity in combination with an increased stickiness of ice could bypass these barriers (Wada et al. 2009;Okuzumi et al. 2012;Kataoka et al. 2013). However, it requires that ice is indeed stickier than rocky material (Gundlach & Blum 2015;Arakawa & Krijt 2021;Schräpler et al. 2022), which might not necessarily be the case (Musiolik & Wurm 2019;Kimura et al. 2020), and that the initial grains are sub-micron in size.
An alternative mechanism that has been extensively studied since its discovery invokes the concentration of pebbles through streaming instability and the subsequent gravitational collapse of dense filament-like structures that converts about millimetresized pebbles directly to about kilometre-sized planetesimals (e.g. Youdin & Goodman 2005;Johansen et al. 2007Johansen et al. , 2014Simon et al. 2016;Schäfer et al. 2017;Abod et al. 2019). These planetesimals then would grow to planet-sized bodies through runaway and oligarchic growth.
Runaway growth occurs when a planetesimal that is slightly more massive than the rest of the bodies accretes more efficiently through gravitational focusing. The mass ratio of the two bodies then increases with time, and the more massive body quickly outgrows the other planetesimals. As the growing body becomes more massive, it starts to gravitationally stir the surrounding planetesimals, which reduces gravitational focusing, and runaway growth ceases (Ida & Makino 1993;Kokubo & Ida 1996;). Eventually, a number of planetary embryos form that grow in an oligarchic fashion by accreting the planetesimals in their respective feeding zones until they reach their isolation mass (Kokubo & Ida 1998, 2000. In the final assembly of planets, these bodies grow by collisions, and when they reach a threshold mass of ∼10 M ⊕ , they start to accrete gas from the surrounding nebula to form the terrestrial and giant planets. Planetesimal accretion has long been thought to be the main pathway of planet formation. However, accretion is only efficient when the planetesimals are small. For small planetesimals of the order of a few kilometres in size at most, gas drag damps eccentricities and inclinations, boosting the accretion rate. However, it is uncertain whether planetesimals actually formed small or if they were large to begin with, typically around 100 km in diameter (Morbidelli et al. 2009;Weidenschilling 2011;Johansen et al. 2015). Evidence for the latter case is seen not only in the size distribution of the asteroid belt (Bottke et al. 2005) and in the cold classical Kuiper belt objects (Kavelaars et al. 2021), which are most likely the unaltered remnants of the planetesimals that formed in the outer Solar System, but also in the absence of large craters on Pluto, which indicates a lack of bodies 1 to 2 km in diameter (Singer et al. 2019). Furthermore, numerical studies of planetesimal formation through the streaming instability indicate a large initial size (Youdin & Goodman 2005;Johansen et al. 2007). The fragmentation of dense pebble filaments into planetesimals results in an initial mass function (IMF) of plan-etesimals that is well described by a power law with exponential cut-off for bodies exceeding a characteristic mass. The characteristic mass translates into a characteristic size of ∼100 km at a heliocentric distance of the asteroid belt, and the largest bodies that form through this process are roughly the size of Ceres (∼10 −4 M ⊕ ) (Simon et al. 2016;Schäfer et al. 2017;Abod et al. 2019;Li et al. 2019).
Because of the long timescales for planetesimal accretion, an efficient formation of terrestrial planets and the cores of giant planets within the lifetime of the protoplanetary disc of typically only a few million years (Haisch et al. 2001) is problematic. Johansen & Bitsch (2019) explored the conditions for forming the cores of the giant planets through planetesimal accretion. Their model focused on the growth track of a single migrating protoplanet sweeping through a population of planetesimals. They found that their fiducial model with constraints from the Solar System of (i) a primordial population of planetesimals of a few hundred Earth masses, (ii) a characteristic planetesimal size of ∼100 km, and (iii) a weakly turbulent protoplanetary disc allow protoplanets to grow to only ∼0.1 M ⊕ within the disc lifetime of 3 Myr. Allowing for a massive disc of planetesimals of ∼1000 M ⊕ produces close-in giant planets, but fails to form cold giant planets, such as Jupiter or Saturn in the Solar System, unless the planetesimal size and turbulence strength are reduced at the same time. The conclusion of Johansen & Bitsch (2019) was that unless all three constraints are ignored, planetesimal accretion is insufficient to grow the cores of giant planets.
In contrast, fast growth can be achieved by the accretion of pebbles that are ubiquitous in the protoplanetary disc. Processes such as the streaming instability that explain planetesimal formation would convert between ∼10 % and up to 80 % of the pebble mass trapped in filaments into planetesimals (Abod et al. 2019). The remnant pebbles as well as newly forming pebbles in the outer disc, where growth timescales are longer, would then provide a mass reservoir for further growth through pebble accretion. Pebble accretion becomes efficient for sufficiently large embryos above the so-called transition mass (Lambrechts & Johansen 2012), at which the growth mode changes from slow Bondi accretion to fast Hill accretion. The transition mass is ∼2×10 −3 M ⊕ at 1 AU and ∼6×10 −3 M ⊕ at 10 AU (Ormel & Klahr 2010;Lambrechts & Johansen 2012). Such a body then accretes a large number of pebbles within a short timescale and can form the terrestrial planets and the cores of the giant planets, consistent with the lifetime of protoplanetary discs (Lambrechts & Johansen 2012Johansen & Lambrechts 2017;Johansen et al. 2021). However, a body of ∼10 −3 to 10 −2 M ⊕ needs to form in the first place, and planetesimal accretion could be the process for this.
In this paper, we investigate whether and under which conditions planetesimal accretion would lead to the formation of pebble-accreting embryos. Most simplified models that were developed to describe the growth of planets start with a narrow annulus of uniformly distributed planetesimals and an embryo in the centre. The growth of the embryo is followed until all planetesimals from the feeding zone are accreted (e.g. Thommes et al. 2003;Chambers 2006;Fortier et al. 2013). It is commonly assumed that there is only one planetesimal size, that the embryo mass follows from the transition from runaway to oligarchic growth, and that the feeding zone of the embryo is always populated with planetesimals. In our model, we wish to deviate in some aspects from this approach. While we also model the growth of an embryo embedded in a population of planetesimals, we employ a different initial situation where we (i) limit the available mass at a certain location by assuming that it is given by the mass budget of a streaming instability filament and (ii) use the streaming instability IMF to derive the characteristic planetesimal size and the size of the embryo. We furthermore deviate from the approach of Johansen & Bitsch (2019) by (i) ignoring migration, (ii) including fragmentation, and (iii) treating the growth rates and the eccentricity and inclination evolution of embryos, planetesimals, and fragments self-consistently. In this way, we study the growth from planetesimals to planetary embryos, while Johansen & Bitsch (2019) focused on the later growth stages in which migration is important. Liu et al. (2019) investigated the growth from planetesimals to embryos and beyond through planetesimal and pebble accretion at the water snowline by means of N-body simulations. They tested different initial conditions for the planetesimal population of (i) a mono-dispersed population of planetesimals, (ii) a poly-dispersed population with IMF from streaming instability simulations, and (iii) a two-component population emerging from runaway growth of planetesimals. They find that a monodispersed population of planetesimals with a size of 400 km fails to form planets because growth timescales are too long due to the rapid excitation of eccentricities and inclinations of the planetesimals. In the other two cases, however, the largest body that forms, either because of the IMF or as a result of runaway growth of 100 km sized planetesimals, grows to several Earth masses firstly by planetesimal accretion and later through pebble accretion when the embryo reaches a mass of ∼10 −3 to 10 −2 M ⊕ . Our work is complementary to their study because we explore the planetesimal accretion phase at various locations, whereas Liu et al. (2019) focused on a single site, the snowline at 2.7 AU.
The paper outline is as follows. In Sect. 2 we outline our model and our assumptions. In Sect. 3 we present the results of planetesimal accretion around a solar-like star, which is our fiducial model. In Sect. 4 we explore and discuss parameter variations of the fiducial model. Finally, in Sect. 5 we summarise and conclude the study.

Basic outline
We used a semi-analytic model to follow the growth of a planetary embryo at a fixed distance r from the central star (e.g. Chambers 2006). Three different types of bodies were considered. These are (i) an embryo, (ii) a population of planetesimals, and (iii) a population of fragments.
The embryo was treated as a single body with mass M em , radius R em , eccentricity e em , inclination i em , and surface density Σ em . For the surface density of the embryo, we simply assumed that the mass of the embryo was distributed uniformly in an annulus of area A em centred at r with a width of br h , where r h =r (M em /(3M )) 1/3 is the Hill radius (Chambers 2006). The value b=10 corresponds to the typical spacing of isolated embryos, which has been shown in N-body simulations to be ∼10 Hill radii (Kokubo & Ida 1998).
A single planetesimal has mass m p and radius R p , and the population has a surface density Σ p . We took the root-meansquare eccentricity and inclination to describe the orbits of the planetesimals. The fragments were treated in the same way as the planetesimals with a mass m fr and radius R fr for a single fragment, and a surface density Σ fr , root-mean-square eccentricity, and inclination for the population.

Mass and surface density evolution, fragmentation
The embryo grows by accreting planetesimals and fragments. The accretion rate of the embryo can be written aṡ where k stands for either planetesimals or fragments. Furthermore, we have the reduced mutual Hill radius h em,k = ((M em + m k ) / (3M )) 1/3 , the Keplerian frequency Ω K , and a dimensionless collision rate P col that is a function of the sizes of the colliding bodies, and their mutual eccentricities and inclinations (Inaba et al. 2001;Chambers 2006). We are now able to formulate the evolution equations for the surface densities. For the embryo surface density, we differentiate Eq. 1 with respect to time and obtaiṅ for the change of surface density due to the accretion of bodies from population k.
To derive the evolution of planetesimal and fragment surface densities, we first need to derive the number of fragments that is produced in a collision between planetesimals. To do this, we first calculated the number of collisions per unit time between planetesimals, where N s,p =Σ p /m p is the surface number density of planetesimals and P col is the collision rate between planetesimals (Inaba et al. 2001). Initially, there are no fragments. When planetesimals are excited to a high enough eccentricity, such that ev K ≈v esc , collisions become disruptive and fragments are produced. Typically, this results in a collisional cascade with a size distribution of fragments, but here we used a typical fragment size of 0.5 km to represent the fragment population. The fragmentation model of Kobayashi & Tanaka (2010) allowed us to determine the total mass of fragments ∆M fr that is produced in a disruptive collision between two planetesimals. The value of ∆M fr depends on the ratio of impact energy and materialdependent critical disruption energy of the planetesimals. By multiplying the collision rate of planetesimals with ∆M fr , we obtain the mass production rate of fragments, Our approach was to keep the total mass of solids constant, that is, we have the condition where the total mass of planetesimals M p is given by Σ p A p , where A p =2πr∆r is the area of the annulus of width ∆r that the planetesimals occupy, and likewise for the fragments. The initial width of the annulus is ηr (see below) for planetesimals and fragments, but the annuli widen diffusively with time because eccentricities and inclinations are excited by viscous stirring Tanaka et al. 2003).
With the assumption of constant total mass, we can now formulate the evolution equations for the surface densities given the mass accretion rate of the embryoṀ (k) em and the production rate of fragments by taking the time derivative of Eq. 6, which giveṡ Substituting the total mass changes with surface densities and areas, we obtain that the surface density of planetesimals reduces due to accretion of planetesimals by the embryo and due to fragmentation, The area ratio appears because of the mass conservation condition. Likewise, the surface density of fragments evolves according tȯ The set of Eqs. 3, 4, 5, 8, and 9 fully describes the mass growth of the embryo and the conversion of planetesimals into fragments. We emphasise that the embryo does not grow to the isolation mass in the classical sense by accreting all the material in the expanding feeding zone. Instead, growth is limited by the available mass contained in the annulus of width ∆r. We considered this approach suitable for studying the growth of an embryo in an isolated filament formed through streaming instability where a fixed amount of mass is converted into planetesimals in a confined narrow ring. Furthermore, while the embryo grows by the accretion of planetesimals and fragments, the mass distribution of planetesimals and fragments does not evolve over time owing to our choice of representing these populations by bodies of a characteristic mass.

Velocity evolution
The velocity distributions, that is, the eccentricities and inclinations, of the bodies evolve through viscous stirring and dynamical friction. To take this into account, we used the rate equations from Ohtsuki et al. (2002) for the root-mean-square eccentricities and inclinations. We included viscous stirring and dynamical friction between all populations with the exception that the single embryo did not interact with itself. Gas drag dampens the orbits of the bodies, and we included damping for the embryo, the planetesimals, and the fragments (Adachi et al. 1976;Inaba et al. 2001). We did not include turbulent stirring of planetesimals through the disc gas.

Protoplanetary disc
We used the self-similar solution for a viscously evolving α-disc (Shakura & Sunyaev 1973;Lynden-Bell & Pringle 1974). The disc is heated by stellar irradiation with a temperature profile of (Ida et al. 2016). The viscously heated part of the disc would initially extend to ∼5 AU for our fiducial parameter choices (Ida et al. 2016). However, we neglected viscous heating here (i) because recent work indicates that irradiation rather than viscous heating might be the relevant heat source for protoplanetary discs (Mori et al. 2019(Mori et al. , 2021, and (ii) because we verified by running a model with a viscous temperature profile where it applies that the choice of temperature profile has negligible impact on the planetesimal accretion process studied here. The filament mass (see Sect. 2.6) and the initial radii of planetesimals and embryos (see Sect. 2.7) vary only within a factor of unity, which does not affect the results.
Article number, page 3 of 12 A&A proofs: manuscript no. 44333corr The viscosity in the α-disc model is ν=αc 2 s Ω −1 K , where we used α=10 −2 . The value of α determines the viscous evolution timescale of the disc and the accretion of the gas onto the star. The value is consistent with what is determined from observations of protoplanetary discs (Hartmann et al. 1998). For the temperature profile used here, the viscosity is a power law in radial distance, ν∝r γ , with an exponent of γ=15/14.
The surface density profile of the self-similar solution for a power-law viscosity is where t vis = r 2 1 / 3(2 − γ) 2 ν 1 is the characteristic time for the viscous evolution, τ= (t/t vis ) + 1, r 1 is the characteristic radius of the disc, and ν 1 =ν(r 1 ) is the viscosity at distance r 1 . For the initial mass accretion rate, we usedṀ ,0 =10 −7 M yr −1 at 0.5 Myr, which corresponds to a typical class II object (Hartmann et al. 1998(Hartmann et al. , 2016. The total disc mass was set to be 10 % of the stellar mass, and by integrating the surface density at t=0 from the inner edge of the disc, which we set to 0.1 AU, to infinity, we can determine the characteristic disc radius, which is ∼72 AU in our fiducial case.

Formation time of planetesimals and embryo
Our model started at time t 0 when planetesimals are expected to have formed. To estimate t 0 , we assumed that the planetesimal accretion phase takes place in the class-II phase of the disc, ∼0.5 Myr after star formation (Evans et al. 2009;Williams & Cieza 2011).
To form planetesimals, dust first grows by coagulation to pebbles. The growth timescales t grow =R/Ṙ of the dust is where g ≈0.5 is a coagulation efficiency and Z is the solid-togas ratio of the dust in the disc (Birnstiel et al. 2012;Lambrechts & Johansen 2014). The time for grains of radius R 0 to grow to pebbles of radius R peb is found to be (Lambrechts & Johansen 2014). Dust growth is limited by radial drift and fragmentation (Birnstiel et al. 2012), and we used the minimum of both for the pebbles. The fragmentation-limited Stokes number is where we set the collision velocity of similar-sized pebbles driven by turbulence (Ormel & Cuzzi 2007) equal to the fragmentation threshold velocity v frag . The value of v frag ranges from ∼1 m s −1 for silicate pebbles (Blum & Wurm 2008) to ∼10 m s −1 for icy pebbles (Gundlach & Blum 2015). However, more recent studies have shown that ice might not be as sticky as previously thought (e.g. Musiolik & Wurm 2019; Kimura et al. 2020) and that the tensile strength of ice aggregates is comparable to the tensile strength of silicates (Gundlach et al. 2018). We therefore used a common fragmentation threshold velocity of v frag =1 m s −1 10 −1 10 0 10 1 10 2 r (AU) for fragmentation-limited growth in our study. The turbulent collision velocity depends on the midplane turbulence α t , which is different from the gas α that drives the viscous evolution of the protoplanetary disc. The value of α t obtained from observations of the dust component of protoplanetary discs ranges from ∼10 −5 to a few times 10 −4 (Pinte et al. 2016;Villenave et al. 2022). We decided to use a value of α t =10 −4 in our study.
In the drift-limited case, the Stokes number is which is obtained from setting the growth timescale of the pebbles (Eq. 12) equal to the drift timescale t drift =r/v r , with v r =2ηv K τ s being the radial velocity of the pebbles. The pebbles are subsequently concentrated by the streaming instability and form planetesimals through the gravitational collapse of dense particle filaments. Depending on the size of the pebbles that form, the streaming instability takes about ten (for large pebbles, τ s ≈0.3) to a few thousand (for small pebbles, τ s ≈10 −3 −10 −2 ) local orbital periods to create dense filaments that subsequently fragment gravitationally into planetesimals (Yang & Johansen 2014;Yang et al. 2017;Li et al. 2018Li et al. , 2019. We calculated the pebble growth time at distance r according to Eq. 13 and the streaming instability timescale as t SI ∼500 2πΩ −1 K and added both to the 0.5 Myr to obtain the initial time t 0 for our simulations.

Filament mass
The initial conditions of the planetesimal population are derived in the framework of planetesimal formation through the streaming instability. We assumed that the streaming instability forms a dense filament of pebbles which fragments into planetesimals. The typical radial width of a filament is ∼ηr, where η is related to the pressure gradient of the disc gas and r is the distance from the star (Yang & Johansen 2014;Liu et al. 2019;Gerbig et al. 2020). This length scale can be thought of as the length scale We show the mass of planetesimals (solid) and embryos (dashed) at initial time t 0 , which ranges from ∼0.5 Myr in the inner disc to ∼1.5 Myr in the outer disc. The planetesimal size is the characteristic size, which follows from the streaming instability IMF, and the embryo size is the single most massive body. The initial sizes for different stellar masses are colour-coded.
over which the Keplerian flow adjusts to the gas flow. We set the solid-to-gas ratio in the filament to Z fil =0.1, and the mass contained in one filament is therefore M fil =2πrΣ gas Z fil ηr (Liu et al. 2019). Because not all pebbles are converted into planetesimals, we introduced a planetesimal formation efficiency p eff . The total mass of planetesimals is then p eff M fil (Liu et al. 2019). For an optimistic upper limit on embryo growth, we assumed p eff =1 throughout the study. Figure 1 shows the mass of the filaments as a function of distance for different stellar masses.

Initial masses of planetesimals and embryo
The initial sizes of planetesimals and embryos follow from the IMF found in streaming instability simulations. The IMF is a power law with an exponential cut-off above a characteristic mass, (Schäfer et al. 2017). The slope of the power law was p≈0.6 and the steepness of the cut-off was q≈0.4 (Simon et al. 2016;Schäfer et al. 2017). The IMF is top heavy, which means that most of the mass is in the large bodies of characteristic mass. We used the characteristic mass above which the IMF drops exponentially as a proxy for the planetesimal mass, which is where γ=4πGρ g Ω −2 K and h is the aspect ratio of the disc (Liu et al. 2020). To determine the mass of the embryo, we calculated the mass of the single most massive body that forms from the IMF. To do this, we set N > (m)=1 in Eq. 16 and solved for m. The value of N tot was found by noting that the total number of bodies are determined by the smallest bodies. We set the minimum mass to m min = 10 −3 m p and calculated N tot =M fil /m min . Figure 2 shows the initial mass and size of the embryo and the planetesimals as a function of distance for different stellar masses.

Diffusion of planetesimals and fragments
We included the diffusive widening of the planetesimal and fragment rings due to viscous stirring Tanaka et al. 2003). This reduced the surface densities with time which impacts both the accretion and the stirring rates. We assumed that the initial width of ηr increases with time as √ 2Dt, where D is the diffusion coefficient and t is the time, as it is characteristic for a random walk (Liu et al. 2019). The diffusion coefficient is related to the viscous stirring rates of eccentricity and inclination Tanaka et al. 2003).

Results
In the fiducial model, we used a solar-mass central star with a viscously evolving disc heated by solar irradiation. We placed filaments at various distances ranging from 0.1 to 100 AU and simulated the growth of the embryo for 10 Myr. Table 1 summarises the model parameters. We later varied the stellar mass, the solid-to-gas ratio of the filament, and other parameters to explore how they affect the growth of the embryo. time t 0 of the simulation, which ranged from ∼0.5 Myr in the inner disc to ∼1.5 Myr in the outer disc. The planetesimal size increased with distance, and the initial mass of the embryo was typically a factor 10 2 to 10 3 higher than the planetesimal mass, as shown in Fig. 2. Inside ∼2 AU, the growth timescale is short and embryos reach their final mass by accreting all the available mass in the filament within 10 Myr. Farther out, growth slows down and ceases to almost zero for 20 AU.

Growth of the embryo
The reasons for the rapid growth inside ∼2 AU are the high surface density of planetesimals, which results in a short growth timescale, and the excitation of the eccentricities of the planetesimals, which results in fragmentation and the boost of growth through accretion of fragments. This is visible in Fig. 4, which shows the time evolution of the surface densities of the planetesimals and the fragments and the eccentricity evolution of the planetesimals in the middle and bottom panels for various distances. The eccentricity at which the collision speed exceeds the escape speed of the planetesimal is approximately given by where we set the random speed of the planetesimals ∼ev K equal to their escape speed. The bottom panel in Fig. 4 shows that within ∼2 AU, the embryo excites the planetesimals above the threshold in short times ( 10 4 yr). As a consequence, the embryo efficiently accretes the small fragments (which we assumed here to have a constant radius of 0.5 km). However, in the outer disc, this effect is negligible because planetesimal eccentricities are not excited enough to result in fragmentation. For example, at 5 AU and for 200 km planetesimals (see Fig. 2), this requires eccentricities 2×10 −2 , which are reached only after ∼1 Myr (after t 0 ) and later. Therefore, fragmentation, if at all, sets in late and embryo growth is not boosted by fragment accretion as is the case in the inner disc. At even larger distances, fragmentation plays no role because the stirring of planetesimals by the embryo and by self-stirring is not sufficient to reach e esc . Therefore, the long accretion timescale limits the growth in the outer disc. Figure 5 shows the eccentricities of embryos, planetesimals, and fragments as a function of distance for different times. Initially, embryos and planetesimals have the same eccentricity and inclination such that e=η/2 and i/e=1/2 because we assumed that just after formation there has not been enough time for dynamical friction to result in a mass dependent eccentricity. The initial eccentricity of the fragments was set to 10 % of the escape speed of a planetesimal.

Planetesimals
In the inner disc, the eccentricity of the planetesimals is determined by the equilibrium of viscous stirring by the embryo and gas drag because the gas density is sufficiently high. The damping timescale for gas drag is where m p and R p are mass and radius of the planetesimal, ρ gas is the gas density, and C D is the drag coefficient (Adachi et al. 1976;Inaba et al. 2001). The timescale on which viscous stirring of planetesimals by an embryo of mass M em and surface density Σ em excites eccentricities is given by (Ida & Makino 1993). When viscous stirring by the embryo and gas drag on the planetesimal are in equilibrium, the eccentricity can be calculated by setting the damping timescale equal to the stirring timescale, Eccentricities of embryos, planetesimals, and fragments at different distances for different times. The different panels show the eccentricities of the embryo (top), the planetesimals (middle), and the fragments (bottom). The eccentricities are shown for different times (relative to t 0 ) as coloured lines. In the inner disc, the bodies follow the equilibrium eccentricity given by the balance of viscous stirring through the embryo and friction with the gas. In the outer disc, eccentricities are determined by viscous stirring with minor damping from gas drag. (Kokubo & Ida 1998) and enters via the embryo surface density Σ em =M em /(2πrbr h ). Figure 2 shows that the mass of planetesimals and embryos scales with distance approximately as m p ∝r 3/2 . Because of the distance dependence of m p and ρ gas , T drag ∝r 47/14 (for our viscous α-disc) increases strongly with distance and the equilibrium eccentricity increases with distance approximately as e eq ∝r 61/70 , which is the slope of e with r, as shown in Fig. 5 inside of ∼0.8 AU.
In the outer disc, damping by gas drag becomes inefficient because of the lower gas density and because of the large planetesimals. Therefore, viscous stirring by the embryo is the process that determines the planetesimal eccentricity. The viscous stirring of planetesimals by the embryo is expressed as (Ida & Makino 1993). Inserting Eq. 20, we can integrate Eq. 22, which gives where T vis,0 is the viscous stirring timescale for initial planetesimal eccentricity e 0 . Figure 4 shows that e∝t 1/4 for large distances, that is, outside of ∼2 AU. Evaluating the scaling with distance, we find that e∝r 1/4 , when taking the r dependency of the embryo mass (initial mass because there is some growth up to ∼20 AU) and of the initial eccentricity (∝r 4/7 because η∝(c s /v K ) 2 ) into account. However, Fig. 5 shows a more flat scaling with r. The reason is that even though the eccentricity evolution is determined by viscous stirring because the damping timescale is too long to result in equilibrium eccentricities, gas drag still damps the eccentricities.

Fragments
The eccentricities of the fragments within ∼1 AU are given by the equilibrium eccentricity. Because fragments are smaller, they are more strongly damped by the gas and hence acquire lower eccentricities. The ratio of fragment size to planetesimal size is ∼10 −2 close to the star, which translates into a ratio of e eq to ∼0.4, as shown in Fig. 5. In the outer disc, eccentricities are excited by viscous stirring and damped by gas drag, where equilibrium values might be reached at late times.

Embryo
The eccentricity of the embryo is more complex, as shown in Fig. 5, because it is determined through the interplay of viscous stirring through the planetesimals, dynamical friction from planetesimals and fragments, and damping through gas drag. The mass growth further complicates the picture, and simple scaling arguments as provided for planetesimals and fragments no longer suffice. However, qualitatively, the embryo keeps a low eccentricity ( 10 −3 ) throughout the simulation. Close to the star, the high gas and planetesimal surface density in combination with the fast growth circularises the orbit. In the outer disc, where no growth occurs, the embryo eccentricity remains close to the initial value experiencing some damping through gas drag and dynamical friction.

Varying the stellar mass
The stellar mass affects the mass accretion rateṀ , the luminosity L , and the density and temperature structure of the disc.
Here, we investigate the growth of embryos around different stellar masses, ranging from 0.1 M to 1 M . The mass accretion rate of low-mass stars is also lower. Manara et al. (2012) provided a fit for the mass accretion rate as a function of stellar age and mass. For our chosen initial time of 0.5 Myr, we finḋ M ∝M . Hartmann et al. (2016) reported that the mass accretion rate correlates with stellar mass asṀ ∝M 2.1 . The linear scaling of Manara et al. (2012), withṀ ∝M , results in smaller and more rapidly evolving discs than the quadratic scaling. We ran simulations with both relations to scale the initial mass accretion rate of 10 −7 M yr −1 for M =1 M to lower stellar masses. The luminosity scales with mass as L ∝M 1..2 for stellar ages 10 Myr (Liu et al. 2020), and hence we set the slope of the L −M -relation to an intermediate value of 1.5. Figure 6 shows the final embryo mass at 10 Myr as a function of distance for different stellar masses forṀ ∝M 2.1 . We find that the maximum distance out to which embryos accrete all available mass of the filament scales with stellar mass, ranging from ∼0.3 AU for a M =0.1 M to ∼1 AU for a M =1 M . Farther out, accretion becomes less efficient and ceases for distances 10 AU for M =0.1 M and 30 AU for M =1 M . In comparison to the quadratic scaling, the final embryo masses forṀ ∝M are shown in Fig. 7. The general finding is the same as for Fig. 6, but the embryos are more massive for all stellar masses <1 M . The reason for this is that Σ gas ∝Ṁ . For the shallower scaling,  the discs around the lower mass stars have higher surface densities and hence the masses of the filaments are higher, which consequently leads to higher masses of embryos and more mass that is available for the embryos to accrete. Figure 8 shows the outcome of our model for different values of the filament solid-to-gas ratio. This Z fil provides a minimum value for how much pebble mass is turned into planetesimals (assuming that the planetesimal formation efficiency is 100 %, which we used in our model). We varied Z fil from the canonical value of Z fil =0.01 to Z fil =0.5, which means that the mass in planetesimals would be 50 % of the gas mass at distance r. We find that varying the value of Z fil does not change the general picture of efficient growth for distances 1 to 2 AU. The final embryo masses vary according to the value of Z fil simply because the filaments are more massive.

Reducing the initial embryo mass
In the fiducial run and variations thereof, we used the single most massive body from the streaming instability IMF as the embryo. In this case, the embryo is typically a factor 10 3 more massive than the planetesimals (see Fig. 2). We ran a model in which we reduced the embryo mass by a factor of 10 while keeping the total mass of the filament fixed. The mass ratio of embryo to planetesimals is hence ∼10 to 100. Figure 9 shows that reducing the initial mass of the embryo does not change the final outcome for 3 AU, where the embryo accretes nearly all the available mass. Outside ∼3 AU, however, accretion efficiency decreases and for 20 AU, the embryo does not grow significantly. etesimal and fragment rings, or extended the simulation from 10 Myr to 1 Gyr. Disabling fragmentation results in longer growth timescales. The final mass of the embryo after 10 Myr, however, is not strongly affected. Inside of ∼1 AU, we find the same final mass while between ∼1 and ∼5 AU, the final mass is lower by less than a factor of 2 at most. Outside ∼5 AU, we find the same final mass as in the fiducial case because fragmentation does not play a role.

Fragmentation, eccentricities, and diffusive widening
The eccentricity of the fragments affect the growth behaviour more strongly. Fixing the fragments on orbits with zero eccentricity and inclination (i.e. assuming that gas drag is very efficient) allows the embryo to accrete fragments at a constant rate in the low-velocity regime (Inaba et al. 2001;Chambers 2006), while in the fiducial case, the accretion rate decreases as frag-ments are excited by viscous stirring through the embryo and the planetesimals. As a consequence, we find that the embryos are more massive than in the fiducial case out to distances of ∼20 AU. We also ran a model in which we set the eccentricity and inclination of the embryo to zero. In this case, we did not find any significant deviation from the fiducial run. We conclude that a circular and planar embryo orbit as used in other studies (e.g. Chambers 2006) is a valid approximation because the eccentricity of the embryo is ∼10 −3 e p , e fr because of dynamical friction and gas drag (Fig. 5).
Lastly, we considered the diffusive widening of the planetesimal and fragment rings. Disabling diffusion has a strong impact on embryo growth. Within 10 Myr and out to ∼5 AU, the embryo accretes all the filament mass, which allows growth to up to ∼0.2 M ⊕ . On much longer timescales than 10 Myr, embryos would grow up to ∼1 M ⊕ out to ∼20 AU. This is shown in Fig. 10, in which we show the final mass after 1 Gyr for the fiducial and the no-diffusion case for comparison. Without diffusion, embryos grow up to ∼1 M ⊕ , whereas with diffusion, the mass even after 1 Gyr is ∼0.1 M ⊕ at most. The reason for the strongly enhanced growth without diffusive widening is that the surface densities of planetesimals and fragments decrease only through accretion. However, the increase in eccentricities and inclinations through viscous stirring causes the bodies to occupy a larger volume, which additionally reduces the surface density and hence reduces the accretion rate of the embryo, which is proportional to the surface density of the accreted bodies.

Implications for pebble accretion
The growth of embryos by planetesimal accretion in the filaments formed through streaming instability is efficient only in the inner part of the protoplanetary disc. In the inner disc, the collision timescale is short enough and fragmentation of planetesimals is efficient enough for an embryo to accrete all the available material. At larger distances and especially outside ∼5 to 10 AU, planetesimal accretion is highly inefficient, even though the available material in the filament increases. The larger sizes of the planetesimals, the excitation of planetesimal eccentricities, and the lack of fragmentation prevents embryos from growing massive within the lifetime of the disc of 10 Myr.
The inefficient growth by planetesimal accretion does not necessarily imply that planets cannot form at all. We did not consider pebble accretion in our model because we focused on the accretion of the filament material turned into planetesimals. However, embryos might still be able to reach masses for which pebble accretion becomes an highly efficient growth process. Pebble accretion becomes important when the friction time of the pebbles is shorter than the time in which they would pass by the embryo (Ormel 2017). This condition leads to the onset mass for pebble accretion (24) (Visser & Ormel 2016;Ormel 2017). Above the so-called transition mass, pebble accretion becomes very efficient. The transition mass marks the change from drift-driven (Bondi) accretion to shear-driven (Hill) accretion of pebbles (Lambrechts & Johansen 2012;Johansen & Lambrechts 2017;Ormel 2017). This means that in the latter case, pebbles from the entire Hill sphere are accreted by the embryo. The transition mass can be found by equating the Bondi radius r B =GM em /(ηv K ) 2 and the Hill radius, and it reads (25) (Ormel 2017). Pebble accretion stops when the embryo reaches the pebble isolation mass. At this mass, the embryo carves a gap in the gas disc that creates a pressure bump outside its orbit, which stops pebbles from drifting inward and being accreted by the embryo. The pebble isolation mass is which is derived from fits to hydrodynamic simulations of pebble accretion (Bitsch et al. 2018).
We now compare the embryo masses to the characteristic masses for the pebble accretion given above. Figure 11 shows a map in which we highlight the different regimes of pebble accretion. For an embryo to accrete pebbles efficiently, the mass needs to be above the transition mass. Below the transition mass and above the onset mass, embryos would still be able to accrete pebbles, but on the less efficient Bondi branch. Figure 11 shows that the initial embryo mass is below the transition mass for all stellar masses; even though the difference is small for M =1 M . For M =0.1 M , the initial embryo mass is even below the onset mass for distances 2 AU. For M =0.3 M , the initial embryo mass is below the onset mass outside of ∼50 AU. Figure 11 also shows that the maximum mass an embryo can reach through planetesimal accretion in a filament (the filament mass M fil ) is above the transition mass (except for the M =0.1 M outside ∼20 AU), but well below the pebble isolation mass, which means that enough mass in planetesimals would be available for embryos to grow into the pebble-accreting regime. However, we find that embryos growing through accretion of planetesimal would reach the transition mass only out to a distance of ∼20 AU for a solar-like central star within the lifetime of a protoplanetary disc. For stars of lower mass, this distance shifts significantly inwards to 1 AU for M =0.1 M . Therefore, we conclude that planetesimal accretion might be a channel for forming the seeds for pebble accretion out to ∼20 AU. Farther out, where planetesimal accretion becomes negligible, pebble accretion even though on the slow Bondi branch would be the only growth channel. Our result agrees with Liu et al. (2019), who investigated the growth of planetesimals to planets at a single site, namely the water snowline at 2.7 AU, through planetesimal and pebble accretion using N-body simulations. In their work, embryos would also grow to masses of 10 −3 to 10 −2 M ⊕ through planetesimal accretion, after which pebble accretion would take over. The final embryo masses in our model at 2.7 AU (Fig. 3 and the bottom right panel of 11) show comparable masses.

Limitations of the model
Our model describes the growth of an embryo at a fixed location. We neglected the migration of the embryo, the planetesimals, and the fragments for several reasons: (i) For the embryo masses considered here ( 0.1 M ⊕ ), the migration timescales are some million years and longer (Tanaka et al. 2002;Cresswell & Nelson 2008;Ida et al. 2020), (ii) in addition to a more complicated model, radial drift of planetesimals and fragments would only reduce the accretion efficiency due to an additional drain of available material, and (iii) gas-drag-induced radial drift peaks for metre-sized bodies, but planetesimals, fragments, and embryos have sizes of one kilometre to several hundreds of kilometres, resulting in slow radial drift. Therefore, the non-migrating bodies provide an upper limit on the final masses. Any migration would result in bodies with lower mass.
In our study, we did not assume that streaming instability forms filaments at special locations in the disc, such as snowlines, where a pressure bump would naturally lead to an increased solid-to-gas ratio due to pile-up of pebbles that would trigger the streaming instability (Drążkowska & Alibert 2017;Schoonenberg et al. 2018). Instead, we studied what would happen if filaments occurred at any location Lenz et al. 2019). A consequence of this would be a significant reservoir of planetesimals that are not accreted. This reservoir can nevertheless interact with the planets that might form by pebble accretion, leading to scattering and populating of the Kuiper belt, scattered disc, and Oort cloud, thus providing the bodies for comets and Kuiper belt objects (Brasser & Morbidelli 2013). The fact that outside of 10 AU neither significant growth nor fragmentation occurs might imply that the size distribution of the planetesimals remains also largely unchanged. The coldclassical Kuiper belt might be a remnant of this. In our model, the total mass of planetesimals in the whole disc is ∼333 M ⊕ . The mass contained in the asteroid belt region between 2 AU and 4 AU is ∼10 M ⊕ , and in the region of the primordial disc region between 15 AU and 30 AU it is ∼50 M ⊕ , both of which are orders of magnitude higher than the current mass in asteroids and in the Kuiper belt, which are estimated to be 4×10 −4 M ⊕ and ∼10 −2 M ⊕ , respectively (DeMeo & Carry 2013;Fraser et al. 2014). Therefore, efficient depletion becomes necessary, such as the giant planet instability in the Nice model, which could have been responsible for sculpturing the outer Solar System and depleting the asteroid belt by scattering of planetesimals and ejection of planetesimals from the Solar System Tsiganis et al. 2005;Morbidelli 2010;Brasser & Morbidelli 2013). On the other hand, we assumed that 100 % of the filament mass is converted into planetesimals, which gives an upper limit on the available mass. The planetesimal formation efficiencies in streaming instability simulations are not well constrained and can vary significantly from 10 % to as high as ∼80 % (Abod et al. 2019). Converting fewer pebbles into planetesimals will reduce the available mass significantly, and additionally, lifting the assumption that filaments form throughout the disc reduces the amount of planetesimals even further.
We furthermore neglected interactions of filaments with each other and that especially the large embryos in the outer disc whose Hill radii might exceed the typical spacing of the filaments would be able to accrete from neighbouring filaments. However, the relevance of this might be low because even though there is a huge mass reservoir of planetesimals of several Earth masses, the embryos accrete almost no planetesimals.

Conclusion
We modelled the planetesimal accretion phase that follows the birth of planetesimals. To do this, we developed a model for the growth of a large planetesimal (embryo) embedded in a population of smaller planetesimals of characteristic size. The model  Fig. 11. Mapping final regimes of embryo masses to pebble accretion. The figure shows the masses of the embryos at 10 Myr for four different stellar masses as a function of distance. The initial masses (dashed lines) and the filament masses (dotted lines) are shown for comparison. Where the mass is below the onset mass M on (Eq. 24) or above the pebble isolation mass M iso (Eq. 26), pebble accretion would be absent (red area). Above the onset mass but below the transition mass M tr (Eq: 25), embryos would accrete pebbles on the Bondi branch (yellow area). Efficient pebble accretion would be possible for embryo masses above M tr (green).
included mass growth of the embryo, the fragmentation of planetesimals, and the velocity evolution of all involved bodies in a self-consistent fashion. We represented the planetesimal size distribution at birth with bodies of characteristic masses, the planetesimals and the embryo. Our growth model hence described an oligarchic-like growth. Fragmentation assumed a representative fragment size. We found that embryos accrete the available material efficiently only in the inner disc where a combination of high planetesimal surface density and fragmentation ensures short growth timescales for the embryo. On the other hand, we find little to no growth in the outer parts of the disc beyond ∼5 to 10 AU on a 10 Myr timescale. The embryos typically reached masses in the range ∼10 −3 to 10 −1 M ⊕ . When we compared the embryo masses to the transition mass for pebble accretion, we found that embryos would be able to grow into the pebble-accreting regime through planetesimal accretion out to ∼20 AU. Pebble accretion on the less efficient Bondi branch might help embryos to reach the transition mass also beyond ∼20 AU.