Issue 
A&A
Volume 625, May 2019



Article Number  A138  
Number of page(s)  12  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/201834168  
Published online  28 May 2019 
Planetesimal fragmentation and giant planet formation
II. Dependencies with planetesimal relative velocities and compositions
^{1}
Facultad de Ciencias Astronómicas y Geofísicas, UNLP,
Paseo del Bosque s/n, 1900FWA La Plata, Argentina
email: irina@fcaglp.unlp.edu.ar
^{2}
CONICET,
Godoy Cruz 2290, CABA, Argentina
^{3}
Instituto de Astrofísica de La Plata, CCT La Plata, CONICETUNLP,
Paseo del Bosque s/n, 1900FWA La Plata, Argentina
^{4}
Instituto de Astrofísica, Pontificia Universidad Catolica de Chile,
Santiago, Chile
^{5}
Instituto Argentino de Radioastronomía, CCTLa Plata, CONICETCICPBA,
CC N. 5 (1894), Villa Elisa, Argentina
Received:
31
August
2018
Accepted:
21
March
2019
Context. Most planet formation models that incorporate planetesimal fragmentation consider a catastrophic impact energy threshold for basalts at a constant velocity of 3 km s^{−1} throughout the process of the formation of the planets. However, as planets grow, the relative velocities of the surrounding planetesimals increase from velocities of the order of meters per second to a few kilometers per second. In addition, beyond the ice line where giant planets are formed, planetesimals are expected to be composed roughly of 50% ices.
Aims. We aim to study the role of planetesimal fragmentation on giant planet formation considering the planetesimal catastrophic impact energy threshold as a function of the planetesimal relative velocities and compositions.
Methods. We improved our model of planetesimal fragmentation incorporating a functional form of the catastrophic impact energy threshold with the planetesimal relative velocities and compositions. We also improved in our model the accretion of small fragments produced by the fragmentation of planetesimals during the collisional cascade considering specific pebble accretion rates.
Results. We find that a more accurate and realistic model for the calculation of the catastrophic impact energy threshold tends to slow down the formation of massive cores. Only for reduced grain opacity values at the envelope of the planet is the crossover mass achieved before the disk timescale dissipation.
Conclusions. While planetesimal fragmentation favors the quick formation of massive cores of 5–10 M_{⊕} the crossover mass could be inhibited by planetesimal fragmentation. However, grain opacity reduction or pollution by the accreted planetesimals together with planetesimal fragmentation could explain the formation of giant planets with lowmass cores.
Key words: planets and satellites: formation / planets and satellites: gaseous planets / methods: numerical
© ESO 2019
1 Introduction
The standard model of giant planet formation is the coreaccretion mechanism (Mizuno 1980; Stevenson 1982; Pollack et al. 1996; Guilera et al. 2010; Helled et al. 2014). According to this model, the planet forms from the accretion of planetesimals onto a solid core until it has enough mass to start accreting gas from the protoplanetary disk. When the crossover mass is reached, the mass of the envelope equals the core mass (at ~10 M_{⊕}), and the planet starts to accrete large quantities of gas in a short period of time; this process is known as runaway gas accretion. At some point, due to mechanisms that are not yet well understood, the accretion of gas onto the planet is limited. Finally, the planet evolves, cooling and contracting at constant mass.
An alternative model has been developed in the past few years, also based on the accretion of solid material, for the formation of giant planets. This model is based on the accretion of smallsized particles (called pebbles) for the formation of the planet core. Unlike planetesimals, pebbles can be accreted by the full Hill sphere making their accretion rates significantly larger than planetesimal accretion rates (Ormel & Klahr 2010; Lambrechts & Johansen 2012).
In the standard core accretion mechanism an important issue is the formation of a massive core before the dissipation of the gaseous component of the protoplanetary disk. In this framework, the process of accretion of planetesimals plays a fundamental role. As the planetary embryo grows, the surrounding planetesimals increase their relative velocities due to the gravitational stirring produced by the embryo, leading to disruptive collisions among them. This effect produces a cascade of collision fragments. Inaba et al. (2003), Kobayashi et al. (2011), Ormel & Kobayashi (2012), and Guilera et al. (2014) found that excessive fragmentation may cause the smaller fragments produced by these collisions to drift inwards through gas drag and a significant amount of mass could be lost stalling the oligarchic growth. Moreover, when the core reaches several Earth masses the relative planetesimal velocities become high enough to produce supercatastrophic collisions among planetesimals (Guilera et al. 2014; Chambers 2014). In this case, Guilera et al. (2014) considered that the planetesimals involved in such a collision are pulverized and the mass of such bodies is lost. On the other hand, Chambers (2014) found that small planetesimals collide frequently and produce rapid embryo growth. Unlike Guilera et al. (2014), Chambers (2014) considered that the mass that would go into fragments smaller than the minimum size considered in the model (referred to as pebbles of second generation) is assigned to become an equivalent mass of the minimum size adopted in their model, since they assume that tiny fragments quickly coagulate into new pebbles.
Despite the different planetesimal fragmentation models, during the process of planetary formation, fragmentation of planetesimals is an important effect that may inhibit or favor the formation of giant planet cores. In a collision between two planetesimals, the catastrophic impact energy threshold required to fragment and disperse 50% of the target mass (also known as specific impact energy) is an important function that has to be defined in models of planetesimal fragmentation. This function depends on many factors of the collision, such as the impact velocity, the planetesimal type of material and its porosity, the sizes of the planetesimals, and other factors that may affect the outcome of the collision (see Jutzi et al. 2015). For example, Benz & Asphaug (1999) found that a target of a fixed size composed of ice material has a catastrophicimpact energy threshold lower than if it were composed of basalt. Benz (2000) also found that lowvelocity collisions are more efficient in destroying and dispersing basaltic bodies than collisions at higher velocities. From smoothed particle hydrodynamics (SPH) simulations, Jutzi et al. (2010) studied the porosity of the planetesimals in a collision and found that, in the strength regime, porous targets are more difficult to disrupt than nonporous ones, while in the gravity regime the outcome of the collision for porous targets depends on gravity and porosity, unlike nonporous ones that are also controlled by the strength. Recently, Beitz et al. (2011) and Bukhari Syed et al. (2017) showed in lowvelocity impact experiments that the strength of compacted dust aggregates is much weaker than that ofporous rocks. This outcome is analogous to Jutzi et al. (2010). It is important to remark here that most models of planet formation that incorporate the planetesimal fragmentation process use the prescription of Benz &Asphaug (1999) for the catastrophic impact energy threshold for nonporous basalts at impact velocities of 3 km s^{−1} (Kobayashi et al. 2011; Ormel & Kobayashi 2012; Guilera et al. 2014; Chambers 2014). However, Guilera et al. (2014) showed thatduring the formation of massive cores the relative planetesimal velocities are incremented from centimeters per second and meters per second to kilometers per second as planet grows.
In our previous work (Guilera et al. 2014), we incorporated a model of planetesimal fragmentation into our global model of giant planet formation (Guilera et al. 2010). We found that the formation of massive cores is only possible from an initial population of big planetesimals and massive disks. In other situations, planetesimal fragmentation tends to inhibit the formation of massive cores before the dissipation of the disk. In this new work, we incorporate into our planetesimal fragmentation model the dependence of the catastrophic impact energy threshold with material composition and planetesimal relative velocities. We also incorporate different velocity regimes for planetesimal relative velocities (Keplerian shear and dispersiondominated regimes) into the planetesimal fragmentation model, and pebble accretion rates for the small particles produced by the collision process into the global model of giant planet formation. The aim of this work is to analyze the impact on the formation of massive cores of the newly incorporated phenomena, especially planetesimals’ relative velocities and planetesimal compositions.
This paper is organized as follows: in Sect. 2, we describe the improvements of our model of giant planet formation and the improvements of the model of planetesimal fragmentation; in Sect. 3, we present the results of our simulations for the formation of a giant planet, and finally, in Sect. 5 we present the summary and conclusions of our results.
2 Our model of giant planet formation
In a series of previous works (Guilera et al. 2010, 2014) we developed a numerical model that describes the formation of giant planets immersed in a protoplanetary disk that evolves in time. In our model, the protoplanetary disk is represented by a gaseous and a solid components. Planets grow by the simultaneous accretion of solids and gas. The solid component of the disk evolves by planet accretion, radial drift due to nebular drag, and collisional evolution, while the gaseous component evolves by an exponential decay (see Appendix A.1 for technical details on the model).
In this work we improve our planetesimal fragmentation model developed in Guilera et al. (2014), incorporating a dependence of , the catastrophic impact energy threshold, with the planetesimal relative velocities and the planetesimal composition considering a mix of nonporous ices and, as a proxy for the rocky material of planetesimals, nonporous basalts. We also incorporate pebble accretion rates for the small fragments that are products of the planetesimal fragmentation process (hereafter, pebbles of second generation), and different velocity regime models for the calculation of low and high relative planetesimal velocities. The improvements in our model are presented in the next subsections.
2.1 Improvements in the solid accretion rates
In the past few years, a new model of solid particles accretion has been proposed as an alternative to the core accretion model for the formation of solid massive cores. Such a model proposes that massive cores, precursors of giant planet cores, can be quickly formed by 100–1000 km embryos accreting centimetersized particles known as pebbles. Ormel & Klahr (2010) and Lambrechts & Johansen (2012) found that pebbles, particles with Stokes numbers below unity (St ≲ 1), are strongly coupled to the gas and can be very efficiently accreted by hundredthousand kilometersized embryos. While planetesimals can be accreted by a fraction of the Hill radius of the growing planet, pebbles can be accreted by the full Hill radius, making pebble accretion rates significantly larger than planetesimal accretion rates.
However, the dominant initial size of the accreted solids (centimetersized particles named pebbles or 1–100 kmsized bodies named planetesimals) is unknown (Helled et al. 2014; Johansen et al. 2014). At the beginning of our simulations, the solid component of the protoplanetary disk is composed of planetesimals of a radius of 100 km with the distribution given by the framework of the oligarchic growth (Kokubo & Ida 1998, 2000, 2002; Ida & Makino 1993). Thus, we assume that there are no initial pebbles in our model. However, pebbles appear as a result of the planetesimal collisional evolution, that is, second generation pebbles.
Following Guilera (2016) and Guilera & Sándor (2017), we improve the solid accretion rates by incorporating the pebble accretion rate for pebbles of second generation in addition to the planetesimal accretion. For planetesimals we adopt the accretion rates given by Inaba et al. (2001) and for pebbles of second generation we use the pebble accretion rates given by Lambrechts et al. (2014) with a reduction factor that takes into account the scale height of the pebbles in comparison with the Hill radius of the planet (see Appendix A.2).
2.2 Improvements to our planetesimal fragmentation model
Our planetesimal fragmentation model, which is based on the Boulder code (Morbidelli et al. 2009 and Supplementary Material), describes the collisional evolution of the planetesimal population. A brief description of the planetesimal fragmentation model and technical details are presented in Appendix B.
2.2.1 Velocities and probabilities of collision regimes
Following Morbidelli et al. (2009), the impact rate between targets and projectiles estimated in our model is given by (1)
where V_{rel} is the relative velocity between targets and projectiles, α is a coefficient that depends on V_{rel}, H is the mutual scale height, a and δa are the mean semimajor axis and the width of the annulus that contains the targets and projectiles, respectively, e_{P} is the eccentricity of the projectiles, E is a coefficient that considers the deviation of the gravitational focusing corresponding to the twobody problem at low relative velocities, b is also a function of V_{rel}, and V_{esc} is the mutual escape velocity. Finally, R_{T} and R_{P} are the radius of the targets and projectiles, respectively (see Morbidelli et al. 2009, for details).
It is important to note, that the relative velocity among planetesimals used in the Boulder code is the dispersion velocity (the velocity corresponding to the dispersion regime). However, when the planetesimal relative velocity tends to zero, the impact rate between them is undetermined. Weidenschilling (2011) noticed this problem and incorporated the relative velocity corresponding to the Keplerian shear. Following Weidenschilling (2011), we incorporate in our planetesimal fragmentation model different velocity regimes and probabilities of collision for lower relative velocities given by Greenberg et al. (1991) to calculate collision rates more accurately. We consider three different regimes and their transitions according to Greenberg et al. (1991): regime A: dominance by random motion (impact rate given by Eq. (1)); regime B: dominance by Keplerian shear motion; regime C: Keplerian shear dominance in a very thin disk (see Appendix B.1 for details).
2.2.2 Catastrophic impact energy threshold
The catastrophic impact energy threshold per unit target mass is the energy needed to fragment and disperse half of the target mass in an impact, that is, the threshold for catastrophic disruption. This quantity plays an important role in the collisional evolution of the planetesimal population. As we mentioned in Sect. 1, depends on many factors of the collision; particularly important are the relative velocity among planetesimals, which determines the impact velocity of the collision, and the planetesimals’ composition.
Figure 1 represents as a function of the radius of a nonporous basalt target for three different values of the impact velocity (~ 25 m s^{−1}, 3 and 5 km s^{−1}), while Fig. 2 represents this threshold for a nonporous icy target for two values of the impact velocity (0.5 and 3 km s^{−1}) (see Appendix B.2). For basalt targets, we can see that for a fixed target’s radius, the smaller the impact velocity, the smaller the . This phenomenon can have an important effect during the formation of massive cores due to the fact that, initially, planetesimal relative velocities are low and are then incremented as the planet grows. We remark again that most works of planetary formation that include planetesimal fragmentation use a fixed value of for nonporous basalt targets at an impactvelocity of 3 km s^{−1} (Ormel & Kobayashi 2012; Chambers 2014). Given the functional dependency of with target’s radius for different impact velocities given by Benz & Asphaug (1999) and Benz (2000), for a fixed target’s radius, we implemented an interpolation between the impact velocity dependencies to get an improved value of as a function of the impact velocity. We also include an interpolation using the different curves of of Benz & Asphaug (1999) for different velocities for ices. We note that if the impact velocities are greater or lower than the velocities corresponding to the upper and lower curves of adopted, we do not extrapolate . In that case, we adopt corresponding to the maximum or minimum velocity used in Benz & Asphaug (1999) and Benz (2000).
On the other hand, according to Lodders (2003), beyond the iceline the amount of solid mass in the protoplanetary disk increases by a factor of two meaning that 50% of the material behind the iceline underwent condensation in the primitive solar system. The icetorock ratio derived from transNeptunian objects, comets, and irregular satellites of giant planets (McDonnell et al. 1987; Stern et al. 1997; Johnson & Lunine 2005) confirms this. As we are interested in the formation of giant planets behind the iceline, we assume that planetesimals are composed of ices and basalts following Lodders (2003). Finally, to implement the dependence of with the impact velocity and the target’s composition, we first interpolate between the curves of as a function of the target’s radius for the different impact velocities, obtaining the values of for each pure material (basalts and ices). Then, we make a linear combination between them, depending on the percentage of basalts and ices that we define for the solids (in this work we considered three cases, planetesimals composed purely by basalts, composed purely by ices, and composed 50% by basalts and 50% by ices). We remark here again that in all the cases is calculated by Eq. (B.8), but using the effective radius (r_{eff} = 3(M_{T} + M_{P})∕4πρ) instead of the target’s radius.
It should be noted that Leinhardt & Stewart (2012) also obtain a derivation of a general catastrophic disruption law. However, their work is focused on the gravitational regime, and for bodies of different porosities, while in our model we adopt compatible laws valid in the gravitational as well as in the strength regime for nonporous bodies.
Fig. 1 Catastrophic disruption thresholds for basalt targets at impact velocities of 20–30 m s^{−1} (Benz 2000), 3 and 5 km s^{−1} (Benz & Asphaug 1999). The blue points correspond to the discrete data extracted from Benz (2000). 
Fig. 2 Catastrophic disruption thresholds for icy targets at impact velocities of 0.5 and 3 km s^{−1} (Benz & Asphaug 1999). 
3 Results
We aim to study the impact of the improvements of our planetesimal fragmentation model discussed above on the formation of a giant planet. We carried out a set of different simulations including one phenomenon at a time. Our simulations start at the beginning of the oligarchic growth with a moonsized embryo located at 5 au from the central star. Initially, the embryo is immersed in an homogeneous singlesized population of planetesimals of 100 km radius and the disk is ten times more massive than the minimum mass solar nebula (MMSN, Hayashi 1981). The simulations stop when the mass of the envelope equals the core mass, that is, when the crossover mass is achieved (in this case we consider that the planet ends its formation in a short period of time after the gaseous runaway growth starts), or at 6 Myr when we consider that the gas of the disk is dissipated.
3.1 Baseline model
In Guilera et al. (2014) we included a first approach to the modeling of planetesimal fragmentation in our model of giant planet formation (Guilera et al. 2010), and we studied, for different initial planetesimal sizes and disk masses, how the collisional evolution of the planetesimal population modified the planetary formation process. We remark that in Guilera et al. (2014) the accretion rates of particles with Stokes number less than the unity were calculated using the prescription derived by Inaba et al. (2001) in the lowvelocity regime. We found that only for initial large planetesimals (r_{P} = 100 km) and massive disks, and only if most of the mass loss in collisions was distributed in larger fragments (see Eq. (B.3)) did planetesimal fragmentation favor the relatively rapid formation of a massive core (larger than 10 M_{⊕}). If smaller planetesimals are considered, the planetesimal fragmentation process inhibits the formation of massive cores. We remark that these results are in agreement with previous works that adopt a similar hypothesis for the model of planetesimal fragmentation (Inaba et al. 2003; Kobayashi et al. 2011; Ormel & Kobayashi 2012). Thus, we choose as our baseline case for this work the most favorable simulation in Guilera et al. (2014), which corresponds to the case of initial planetesimals of 100 km radius and a 10 MMSN initial disk mass.
In Fig. 3, we show the time evolution of the core mass and envelope mass of the planet, for the baseline model and for the formation of the planet without considering the planetesimal fragmentation process. We can see that the inclusion of planetesimal fragmentation reduces by ~7% the time to reach the crossover mass. We can also observe that the planet reaches the crossover mass at a lower core mass (18% lower). However, the fragmentation model incorporated in Guilera et al. (2014) considered the catastrophic disruption threshold given by Benz & Asphaug (1999) for basalts at impact velocities of 3 km s^{−1}. This is clearly a simplification. On one hand, as we are studying the formation of a giant planet behind the iceline, it should be taken into account that planetesimals are composed of rocks and ices (Lodders 2003). Also, during the process of planet formation, the planetesimal relative velocities are increased due to the gravitational perturbations produced by the growing planet. In the left panel of Fig. 4, we show the increase of eccentricities and inclinations for planetesimals of 100 km radius, and in the right panel their relative velocities, near the planet while it grows. We can see that eccentricities and inclinations are quickly increased near the planet. This leads to an increase of the planetesimal relative velocities from velocities of some meters per second to velocities of about 5 km s^{−1} (right panel of Fig. 4). These results show that clearly the assumption of a given for a constant impact velocity is a simplification in the model of planetesimal fragmentation, and thus, it motivated us to incorporate the following dependencies in the catastrophic disruption threshold.
Fig. 3 Core masses (solid lines) and envelope masses (dashed lines) as function of time. Here, as in Figs. 5, 8, 9, 12, and 13, the red thick lines correspond to the baseline model, while the gray thin lines correspond to the formation of the planet without planetesimal fragmentation. 
3.2 Dependencies of Q with the planetesimal impact velocities and compositions
As we mentioned above, in general models of planet formation that include planetesimal fragmentation consider the catastrophic impact energy threshold for a fixed velocity of basaltic planetesimals. In this section, we explore how the dependencies of the catastrophicdisruption threshold with the planetesimal compositions and with the planetesimal impact relative velocities impact affect the formation of the planet. To do this, we first calculate the formation of the planet considering the same initial parameters as in the baseline model (an embryo located at 5 au in a disk ten times more massive than the MMSN and with initial planetesimals of 100 km radius) but now we implement the dependence of with different impact velocities following the approaches described in Sect. 2.2.2. The solid accretion rates are calculated following Eq. (A.6). In Fig. 5, we show the time evolution of the core mass and envelope mass of the planet for the model wherein planetesimal fragmentation is not considered, for the baseline model, and for the case where we consider basaltic planetesimals with as a function of the relative velocities among planetesimals. We can see that the crossover mass, for the last case, is achieved in a longer time than for the baseline model. Moreover, the formation of the giant planet takes ~ 25% more time compared to the case wherein planetesimal fragmentation is not considered. In Fig. 6, we plot the solid accretion rates for the different sizes of planetesimals and small particles for both models, the baseline model, and the model that considers as a function of the relative velocities among planetesimals. We can see that the accretion of planetesimals smaller than 100 km and small particles generated by the collisional evolution becomes effective at a longer time for the model where is a function of the relative velocities. We can also see that in both models, the total solid accretion rate is always dominated by plantesimalsof 100 km. However, while in the baseline model planetesimals between 1 and 25 km are the products of the fragmentation process that most contribute to the total solid accretion rate (being the planetesimals of ~ 25 km the most important), in the modelwhere is a functionof the relative velocities, the most important contribution from the fragments is given by planetesimals between 10 and ~ 16 km. Moreover, in this last case, planetesimals of 1 km make a negligible contribution to the total solid accretion rate. This is due to the fact that these bodies are near the minimum of for low impact velocities, which is for this size one order of magnitude lower with respect to the at impact velocities of 3 km s^{−1}. We can see from Fig. 7 that relative velocities for planetesimals of 1 m reach 20–30 m s^{−1} at 1 Myr (when the mass of the core is ~1 M_{⊕}), and 1 kmsized planetesimals exceed thesevalues rapidly; at 10^{−2} Myr (when the mass of the core is only ~0.2 M_{⊕}) the relative velocity is over 200 m s^{−1}. It is important to remark that, in both cases, their velocities are always below 3 km s^{−1}.
On the other hand, in Fig. 8 we plot the time evolution of the core mass and envelope mass of the planet for a model that considers planetesimals purely composed of ices and adopting for ices at impact velocities of 3 km s^{−1}, and for a model wherein planetesimals are composed of 50% of basalts and 50% of ices, and is calculated by a linear combinationbetween for basalts at 3 km s^{−1} and for ices at 3 km s^{−1}. We also show, for completeness and comparison, the baseline model and the case where planetesimal fragmentation is not considered. We can see that when we considered planetesimals purely composed of ices, the crossover mass is never achieved. This is due to the fact that for the same fixed impact velocity (in this case 3 km s^{−1}), is lower for ices thanfor basalts. Thus, the collisional evolution of the planetesimal population is different, small fragments are disrupted more efficiently, and the accretion rates of such fragments do not compensate for the fragmentation of big bodies. Besides, if we compare from Figs. 1 and 2 we can see that for ices at 3 km s^{−1} is lower than for basalts at 20–30 m s^{−1} for planetesimals of ~100 km. Therefore, this difference can explain why in the model that considers as a function of the relative velocities for basalts the giant planet core could form during the disk lifetime and, for the model that adopts for ices at 3 km s^{−1}, the crossover mass was not achieved in that time. When planetesimals are composed of 50% of basalts and 50% of ices, and is calculated as a linear combination of the corresponding for basalts and ices at impact velocities of 3 km s^{−1}, respectively, the crossover mass is achieved in less than 6 Myr with a core mass of 16 M_{⊕}. However, the planet achieves the crossover mass after a longer time, by ~12%, compared to the case without planetesimal fragmentation.
All these results show that the collisional evolution of the planetesimal population and the accretion of the fragments produced in such collisional evolution play an important role in the formation of a giant planet. We point out that in this subsection, the accretion rates of particles with Stokes number less than unity was also computed as in the baseline model, using the prescription derived by Inaba et al. (2001) in the lowvelocity regime.
Fig. 4 Left panel: time evolution of the radial profiles of the eccentricity (solid lines) and inclination (dashed lines) for 100 kmsized planetesimals. Right panel: time evolution of the relative velocity between planetesimals with 100 km radius. In brackets we show the total mass of the planet at the given times. The black curve corresponds to the time (and mass) at which the planet reaches the crossover mass. 
Fig. 5 Coremasses (solid lines) and envelope masses (dashed lines) as a function of the time. Gray lines: model without planetesimal fragmentation. Red lines: baseline model. Pink lines: model in which is a functionof relative velocities. 
Fig. 6 Top panel: solid accretion rates for basaltic planetesimals and small particles of various radii, products of the collisional evolution of planetesimals, for the baseline model. Bottom panel: same as top panel but for the model that considers as a function of the planetesimal relative velocities. 
Fig. 7 Time evolution of the relative velocity between basaltic planetesimals of the same size for the model that considers as a function of the planetesimal relative velocities. Top panel: 1 msized planetesimals. Bottom panel: 1 kmsized planetesimals. In brackets we show the total mass of the planet at the given times. 
Fig. 8 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model with no planetesimal fragmentation. Red lines: baseline model. Pink lines: model with for ices at impact velocities of 3 km s^{−1}. Blue lines: model where is calculated by a linear combination between for basalts at 3 km s^{−1} and for ices at 3 km s^{−1}. 
3.3 Accretion of small fragments: pebble accretion of second generation
As we mentioned before, in out model we improved the solid accretion rates, incorporating in this work the accretion rates of pebbles (Eq. (A.7)) derived by Lambrechts et al. (2014) for particles with Stokes number less than or equal to unity. In Guilera et al. (2014), the accretion rates of such particles were calculated using the prescription derived by Inaba et al. (2001) in the lowvelocity regime (Eq. (A.6)). In that work, we showed that the accretion rates for particles smaller than ~ 1 m could become greater than the pebble accretion rates (see Guilera et al. 2014, for the details). Thus, in order to calculate in a more accurate way the accretion of small particles, we implemented the pebble accretion rates mentioned above.
In Fig. 9, we show the time evolution of the core mass and the envelope mass for the case where planetesimals are purely composed of basalts adopting the corresponding at impact velocities of 3 km s^{−1}, and where small particles with Stokes number less or equal to unity are accreted as pebbles using Eq. (A.7), while in the baseline model the accretion rate for planetesimals (Eq. (A.6)) is applied also for small particles. We also plot for comparison the baseline model and the case wherein planetesimal fragmentation is not considered. We can see that a more accurate treatment for the accretion of small particles delays the formation of the giant planet by about 30%, despite the fact that pebbles do not mostly contribute to the total solid accretion rate (Eq. (A.7)).
This effect is discussed in Guilera et al. (2014) where the rates of small fragments using the prescriptions of Inaba et al. (2001) and Lambrechts & Johansen (2012) are compared. Following that analysis we here obtain the accretion rates of Inaba et al. (2001) in terms of the pebble accretion rates of Lambrechts et al. (2014): (2)
where is the enhanced radius due to the planet’s gaseous envelope, and , where Ṁ_{H} for 0.1 ≤ St < 1 and Ṁ_{H2} for St < 0.1 are the pebble accretion rate dM_{C}∕dt of Eq. (A.7). We can see from Fig. 15 of Guilera et al. (2014) that when the planet’s core mass becomes larger than ~0.2 M_{⊕} planetesimal accretion rates are higher than pebble accretion rates for fragments with r_{P} ≲ 1 m. In our work this is the case for particles with 0.1 ≤ St < 1, but for fragments with St < 0.1, the difference between the planetesimal and pebble accretion rates is even higher due to the factor . Moreover, we introduced a factor β = min(1, R_{H}∕H_{P}) to include a reduction in the pebble accretion rates if the scale height of the pebbles becomes greater than the Hill radius of the planet; this is not taken into account in the planetesimal accretion rate for small fragments in the baseline model. In Fig. 10 we show R_{H}∕H_{P} as a function of time, where we can observe that objects between 1 and 10 cm have values of R_{H} ∕H_{P} below 1. The quantified differences are shown in Fig. 11 where we present the comparison of the solid accretion rates for small particles with St ≲ 1 of the baseline case and the simulation including pebble accretion. We can see from Fig. 11 that the accretion rates in the baseline model are higher than in the simulation that takes into account the accretion of pebbles as explained before.
Fig. 9 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model with no planetesimal fragmentation. Red lines: baseline model. Pink lines: model with for basalts at impact velocities of 3 km s^{−1} and where small particles with Stokes numbers less than or equal to unity are accreted as pebbles of second generation. 
Fig. 10 Ratio between the Hill radius and the scale height of small particles, with Stokes number less than or equal to unity, for the baseline model (solid lines) and for the simulation including pebble accretion (dashed lines) as a function of time. 
Fig. 11 Solid accretion rates for small particles, product of the collisional evolution of planetesimals, for the baseline model (solidlines) and for the simulation including pebble accretion (dashed lines) as a function of time. Small particles of these sizes have Stokes number less than or equal to unity. 
3.4 Global model
In this section, we compare the results for the formation of the giant planet obtained with the baseline model and with the global model. The global model includes all the improvements in the calculation of , that is, the dependencies of with the planetesimal impact velocities and their compositions described in Sect. 3.2, and the accretion of second generation pebbles described in Sect. 3.3.
In the global model where we include all the new effects on and a more accurate and realistic treatment for the accretion of the small particles, the core does not achieve the crossover mass within the 6 Myr as shown in Fig. 12. We point out that this model also includes the different regimes for the calculation of the relative velocities (and the impact rates) discussed in Sect. 2.2.1. However, as the planet quickly excites the relative velocities of the near planetesimals, they are in general in the dispersion regime, and thus the Keplerian shear regime does not play a relevant role. Finally, we remark that despite the fact that the planet does not achieve the crossover mass, a solid core of a few Earth masses, up to five Earth masses, is formed more quickly compared to the case where planetesimal fragmentation is not considered. Ikoma et al. (2000) and Hubickyj et al. (2005) showed that a reduction in the grain opacity in the planet envelope, as well as the pollution of the envelope (due to evaporated materials of icy planetesimals in the envelope, see Hori & Ikoma 2011), could reduce the mass of the core at which the planet reaches the crossover mass, and thus reduce the formation time. We will explore this possibility in the next section.
Fig. 12 Coremasses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model without planetesimal fragmentation. Red lines: baseline model. Pink lines: model with all the improvements in the calculation of (the dependencies of with the planetesimal impact velocities and compositions) and the pebble accretion rates. 
3.5 Reduced grain opacities
In the previous sections we studied the process of giant planet formation considering fragmentation of planetesimals including, for the catastrophic disruption threshold, the dependencies with relative velocity and composition of the colliding bodies. We also include the accretion rates for small fragments called pebbles. We analyzed every improvement separately and finally all together and we found that, in all the cases, the formation of the giant planet core is slower.
Motivated by this result we incorporate a physical phenomenon that could act in an opposite way, accelerating the formation of a giant planet, in this case a reduction in the grain opacities of the planet envelope. The grain opacities in the planet envelope play an important role in the formation of a giant planet. Ikoma et al. (2000) found that if the grain opacities are small enough, a smallsized giant planet core can capture significant amounts of gas. Later, Hubickyj et al. (2005) found that if the grain opacities are assumed to be 2% of the interstellar medium value, the planet can reach the crossover mass with a core mass between 5 and 10 M_{⊕}. More recently, Hori & Ikoma (2010) showed that a core of only ~2 M_{⊕} is able to capture enough gas to form a giant planet if the accreting envelope is grainfree. Thus, following these works, especially the work of Hubickyj et al. (2005), we compute a set of simulations with and without planetesimal fragmentation, reducing the grain opacities in the envelope to up to 2% of the interstellar medium value. For the case where we consider the planetesimal fragmentation, we adopt the global model described above.
In Fig. 13, we present the time evolution of the core mass and the envelope mass of the planet for the models that consider reduced grain opacities with and without planetesimal fragmentation and the baseline case. When we calculate the formation of the giant planet without considering the planetesimal fragmentation process, the time at which the planet achieves the crossover mass is reduced by more than 50% when the grain opacities in the planet envelope are reduced. However, the crossover mass remains practically similar. These results are in very good agreement with the previous results of Hubickyj et al. (2005). A reduction in the envelope opacity allows the planet to release more efficiently the heat generated by the accretion of the planetesimals, and as a consequence, the gas accretion becomes more efficient. We can see this effect in Fig. 14 where we plot the envelope mass as a function of the core mass for the two models without planetesimal fragmentation. At a fixed core mass, the model with the reduced grain opacities has a greater value for the envelope mass.
Finally, we note that when planetesimal fragmentation is considered, the planet achieves the crossover mass at a time (1.66 Myr) and with a core mass (16.44 M_{⊕}) that are lower than for the case where planetesimal fragmentation is not considered (~7 and 24% lower in time and core mass, respectively). We associate this result to the fact that the amount of envelope mass in the planet at lowmass cores plays an important role, significantly enhancing the capture radius of the planet.
Fig. 13 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model without planetesimal fragmentation. Black lines: model without planetesimal fragmentation but considering grain opacities reduction (G.O.R). Red lines: baseline model. Pink lines: model taking into account all the improvements in the calculation of (the dependencies of with the planetesimal impact velocities and compositions), pebble accretion rates, and considering reduced grain opacities. 
Fig. 14 Envelope mass as a function of the core mass for the two cases where planetesimal fragmentation is not considered. We plot the core mass until 10 M_{⊕} in order to show more clearly the differences at lowmass core values. The black line represents the case where grain opacities are not reduced, while the red line represents the case where grain opacities are reduced in the planet envelope. 
4 Comparison with previous work
A similar model of giant planet formation was developed by Chambers (2014). They carried out numerical simulations of oligarchic growth including pebble accretion and planetesimal fragmentation. One of the major differences between their model and ours is that they consider that all the collisions involving targets smaller than their minimum bin size are assumed to coagulate and grow rapidly when they approach such minimum size, while we assume that the mass distributed under cm is lost. In our fragmentation model most of the mass is in the largest fragments. As discussed in Guilera et al. (2014), the results are analogous when the minimum bin size is reduced to lower values. Also, Chambers (2014) uses a fixed value of (for basalts at 3 km s^{−1}) and includes an additional factor to reduce 10 and 100 times the strength for all the sizes of the bodies to see how this could change the growth of the embryo. He found the most significant differences in the case where is reduced by a factor of 100, changing the growth track of the embryo (shown in Fig. 20 of his paper). Even though we cannot compare directly this result with ours since our models are different, we agree that the relative importance of pebble accretion for giant planet formation will depend strongly on the rate at which pebbles are generated from planetesimalplanetesimal collisions, depending in particular on the impact parameters.
Anothermodel of giant planet formation that includes pebble accretion is that of Alibert et al. (2018). They developed a model of the formation of Jupiter that could explain the constraints that the cosmochemical evidence gives on the existence of two main reservoirs of small bodies that remain separated during the early solar system (Kruijer et al. 2017). They conclude that Jupiter formed in a threestep process. First, Jupiter’s core grew by the accretion of pebbles. Second, the pebble accretion stopped and the core continued growing by the accretion of small planetesimals. Third, the core was massive enough for the runaway gasaccretion to start. It is important to highlight that our initial conditions correspond to the beginning of the oligarchic growth where almost all the solid material is contained in planetesimals, in contrast to Alibert et al. (2018) where they start their simulations with a first generation of pebbles. Moreover, they do not include a fragmentation model since they estimate, for the start of their second stage, the amount of fragments that could have been produced during their first stage of pebble accretion. For this estimation, they used a fixed value of for 100 kmsized planetesimals. Different initial conditions for giant planet formation including a first generation of pebbles and a variation of the initial size of planetesimals will be studied in future works.
5 Summary and conclusions
We studied the formation of a giant planet considering the collisional evolution of the planetesimal population, including the dependence of the catastrophic impact energy threshold with composition and relative velocities of the colliding planetesimals. We also included the pebble accretion rates for small particles produced by the collisional cascade, which are pebbles of second generation strongly coupled to the gas, to calculate in a more accurate way than in our baseline case. We showed that as the planet grows, the relative velocities among the planetesimals near the planet are increased from velocities of the order of meters per second to velocities up to ~ 4.5 km s^{−1} before the planet reaches the crossover mass.
We analyzed the improvements incorporated in the calculation of and the inclusion of the pebble accretion rates one at each time, and finally all the improvements together. We point out that, for each improvement in , as the collisional cascade increases, the timescale to reach the crossover mass is delayed. We also found that in all the cases, the formation of the core is less efficient compared with the case where planetesimal fragmentation is not considered, and with our baseline model.
When we include all the improvements at once there is a sum of effects that inhibit the formation of the giant planet core. On one hand, pebble accretion rates are lower than the planetesimal accretion rates for the small particles adopted in the baseline model. This effect contributes to slow down the growth of the giant planet core. On the other hand, when we include the dependence of with impact velocity and composition we first interpolate between the curves of for a given velocity obtaining for each pure material. Then, we make a linear combination between these values where we assume that planetesimals are half basalts and half ices. For low impact velocities basalts are weaker than ices (except for ~100 kmsized bodies) and for higher velocities, ice is weaker than basalts. These effects are averaged since the composition of planetesimals is 50% basalts and 50% ices. Moreover, it is important to remark that the catastrophic disruption thresholds for ices at impact velocities of 0.5 and 3 km s^{−1} are lower for ~100 kmsized planetesimals than for basalts at any of the three impact velocities (20–30 m s^{−1}, 3 and 5 km s^{−1}) for that size of bodies. Regardless of whether the planet arrives at the crossover mass or not within the disk lifetime, in all the cases, planetesimal fragmentation favors the relative rapid growth of the giant planet core up to a few Earth masses compared to the case where planetesimal fragmentation is not considered.
This last result motivated us to explore alternatives for the formation of a giant planet with less massive cores. In this way, we followed the work of Hubickyj et al. (2005). These authors showed that if the grain opacity of the planet envelope is reduced to 2% of the interstellar medium value, the planet can reach the crossover mass with a core mass between 5 and 10 M_{⊕}. Thus we computed two new simulations, reducing the grain opacity of the envelope up to 2% of the interstellar medium value, considering planetesimal fragmentation (including all the improvements in our model) and without planetesimal fragmentation. For this last case, we found that the planet reached the crossover mass at a time significantly lower (by a factor two) with respect to the case without planetesimal fragmentation but not considering reduced grain opacities. However, the mass of the core at which the crossover mass is achieved remained practically similar. For the case where planetesimal fragmentation and reduced grain opacities are included, the planet reached the crossover time at a lower time and with a less massive core compared to the case where reduced grain opacities are included but without planetesimal fragmentation. We associated this result to the fact that if the opacity of the planet envelope is reduced, the planet releases more efficiently the luminosity through the envelope, and then, the gas accretion becomes more efficient. Thus the planet accretes significantly more gas for a lowmass core, enhancing the radius of capture of the planet and increasing the solid accretion rate.
Hori & Ikoma (2011) studied the formation of a giant planet where the gas envelope is polluted by icy planetesimals that are dissolved in the envelope. They found that the increase in the molecular weight and the reduction in the adiabatic temperature gradient produced by the pollution of planetesimals significantly reduce the mass of the core at which the planet achieves the crossover mass. More recently, Venturini et al. (2015, 2016) found similar results, computing for the first time selfconsistent models of giant planet formation that include the effect of envelope enrichment by the pollution of planetesimals dissolved as they enter the planet envelope. The increment in the envelope molecular weight and the reduction in adiabatic temperature gradient mean that the envelope layers are compressed more efficiently and the gas accretion is significantly increased. We will incorporate this phenomenon in a future work.
Acknowledgements
We would like to thank John Chambers and the anonymous referee for the comments and suggestions that helped us to improve this paper. This work was supported by the PIP 11220150100699CO grant from Consejo Nacional de Investigaciones Científicas y Técnicas, and by PIDT 11/G144 grant from Universidad Nacional de La Plata, Argentina. O.M.G is also supported by the PICT 20160023 from ANPCyT, Argentina. O.M.G. also acknowledges being hosted as an invited researcher by IAPUC.
Appendix A Accretion model
A.1 Evolution of the envelope and planetesimal radial migration
The evolution of the envelope is calculated solving the standard equations of the stellar evolution theory. The following equations correspond to the conservation of mass, the hydrostatic equilibrium, the energetic balance, and the energy transport: (A.1) (A.2) (A.3) (A.4)
where G is the universal gravitational constant, ρ is the density of the envelope, S is the entropy per unit mass, ɛ_{pl} is the energy release rate due to the accretion of planetesimals, and is the dimensionless temperature gradient, which depends on whether the energy is carried by radiation or by convection (see Fortier et al. 2009; Guilera et al. 2010, for details).
In our simulations we consider 2500 radial bins logarithmically equally spaced in a protoplanetary disk defined between 0.4 and 20 AU. The particles migrate inward in the protoplanetary disk due to gas drag, where the radial migration velocity for each drag regime is given by (A.5)
being a the semimajor axis, η the ratio of the gas velocity to the local Keplerian velocity, St = t_{stop}ω_{k} the Stokes number, with t_{stop} the stopping time depending on the drag regime (Epstein, Stokes or quadratic regime) and ω_{k} the Keplerian frequency. The rest of the model is explained in detail in Guilera et al. (2010, 2014).
A.2 Solid accretion rates
For planetesimals, we use the accretion rates given by Inaba et al. (2001), (A.6)
where M_{C} is the core’s mass, R_{H} is the planet’s Hill radius, Σ_{P} is the surface density of solids at the location of the planet, Ω_{P} is the Keplerian frequency atthe planet location, and P_{coll} is the collision probability, which is a function of the core radius R_{C}, the Hill radius of the planet, and the relative velocity between the planetesimals and the planet v_{rel}, thus P_{coll} = P_{coll}(R_{C}, R_{H}, v_{rel}). In fact, as we also consider the drag force that planetesimals experience on entering the planetary envelope (following Inaba & Ikoma 2003), the collision probability is a function of the enhanced radius instead of R_{C} (Guilera et al. 2014).
For the pebbles of second generation, we use the pebble accretion rates given by Lambrechts et al. (2014), (A.7)
In Eq. (A.7), we introduce the factor β = min(1, R_{H}∕H_{p}) in order to take into account a reduction in the pebble accretion rates if the scale height of the pebbles, H_{p}, becomes greater than the Hill radius of the planet. This can happen if the vertical turbulent dispersion of small particles becomes significant. The scale height of the solids at a given distance from the central star is given by (Youdin & Lithwick 2007) (A.8)
being H_{g} the scale height of thedisk, and α the dimensionless Shakura & Sunyaev viscosityparameter (Shakura & Sunyaev 1973); we adopt α =10^{−3}.
Appendix B Fragmentation model
According to our model, a collision between a target of mass M_{T} and a projectile of mass M_{P} results in a remnant of mass M_{R} given by (B.1)
where Q is the collisional energy per unit target mass and is the catastrophic impact energy threshold per unit target mass required to fragment and disperse half of the target mass. Usually, is a function of the target’s radius. However, it is important to remark here that in our model (following Morbidelli et al. 2009) has to be calculated with an effective radius r_{eff} = 3(M_{T} + M_{P})∕4πρ, being ρ the planetesimal density. Besides this, in our model the mass loss in the collision, which we define as (M_{T} + M_{P} − M_{R}), is distributed between the minimum mass bin considered and the mass bin corresponding to the biggest fragment M_{F} given by (B.2)
As in Guilera et al. (2014), we note again that for some supercatastrophic collisions (which occur when M_{R} ≪ M_{T} + M_{P}), M_{F} > M_{R}. For such collisions we assume that M_{F} = 0.5M_{R}.
Unlike Morbidelli et al. (2009), the fragments are distributed following a powerlaw distribution given by (Kobayashi et al. 2011; Ormel & Kobayashi 2012) (B.3)
meaning that most of the mass is distributed in the larger fragments. In Guilera et al. (2014), we found that massive core formation is favored when the exponent of the powerlaw mass distribution is lower than 2, thus in this work we analyze the most favorable scenarios, taking the value of 5/3 for the exponent. To calculate the growth of the planet, we adopt a discrete planetesimal size distribution using 36 size (or mass) bins logarithmic equally spaced between 1 cm and 100 km, where initially all the solid mass is in the nonporous planetesimals of 100 km radius. However, when we calculate the planetesimal fragmentation process, we extrapolate the planetesimal size distribution two orders of magnitude below the minimum size of the main model to avoid the accumulation of spurious mass in the smaller fragments. Therefore, only the mass ejected from the collision distributed between the mass of the larger fragment and the minimum size considered ( cm) is taken into account to calculate later the solid accretions rates; that is, we assume that the mass distributed below 1 cm is lost.
The feeding zone of the embryo extends to four Hill radii at either side of the embryo. Adopting 2500 radial bins along the protoplanetary disk guarantees that there are at least ten radial bins between R_{P} − 4R_{H} and R_{P} + 4R_{H} at the beginning of the simulation. We define the width of the fragmentation zone as twice the feeding zone, that is, eight times the embryo Hill radius at both sides of the embryo. The excitation of eccentricities and inclinations decay with the distance to the embryo, especially outside the feeding zone. Then, our definition of the fragmentation zone guarantees that collisions are well determined within this zone. The amount of bins increases as the core mass grows, for example, in the baseline case (see Sect. 3.1), when the planet reaches the crossover mass, the fragmentation zone has ~ 600 radial bins.
B.1 Velocities and probabilities of collision regimes
Following Greenberg et al. (1991) we adopt three different regimes and their transitions, regime A: dominance by random motion; regime B: dominance by Keplerian shear motion; regime C: Keplerian shear dominance in a very thin disk. For the dispersion regime (dominance by random motion), where Keplerian behavior is unimportant, we adopt the impact rate given by Eq. (1).
The transition between regime A and regime B is given by (B.4)
where is the Hill radius of the target, and a_{P}, e_{P} and a_{T}, e_{T} are the semimajor axis and eccentricities of the projectile and the target, respectively. Then, for larger values of e and a the system is dominated by random motion (regime A) and for smaller values of these orbital parameters the system is in the Keplerian shear regime (regime B). If a_{p}i_{p} < R_{G}, being R_{G} the target’s gravitational diameter and i_{p} the inclination of the projectile’s orbit, the particles are in regime C, wherein still dominates the Keplerian shear motion but the system is two dimensional.
Free parameters that determine for different types of materials and different impact velocities from Benz & Asphaug (1999) and Benz (2000).
For regimes B and C the relative velocity is given by (B.5)
where μ = M_{p}∕M_{⊙} and . The impact rate of regime B follows (B.6)
where σ is the projectile’s surface density and ω is the Keplerian frequency using a = (a_{P} + a_{T})∕2.
Finally, the impact rate of regime C is given by (B.7)
B.2 Catastrophic impact energy threshold
From SHP simulations, Benz & Asphaug (1999) found that can be expressed by the functional form (B.8)
being Q_{0}, B, a, and b parameters that depend on the properties of the material and on the impact velocity over the target, and ρ the density of the nonporous planetesimals (in this work we adopted ρ = 1.5 g cm^{−3}). They performed simulations for nonporous basalts at 3 and 5 km s^{−1}, and for nonporous ices at 0.5 and 3 km s^{−1}. Later, Benz (2000), performed new SPH simulations to calculate for nonporous basalts at low impact velocities (between 20 and 30 m s^{−1}) finding that targets impacted at such low impact velocities are weaker than targets impacted at greater velocities. Despite the fact that Benz (2000) did not provide a functional form for , we fit the results of their simulations adopting the same functional form proposed by Benz & Asphaug (1999). Table B.1 summarizes the values of the free parameters that determine the values of used in this work for different types of materials and different impact velocities.
References
 Alibert, Y., Venturini, J., Helled, R., et al. 2018, Nat. Astron., 2, 873 [NASA ADS] [CrossRef] [Google Scholar]
 Beitz, E., Güttler, C., Blum, J., et al. 2011, ApJ, 736, 34 [NASA ADS] [CrossRef] [Google Scholar]
 Benz, W. 2000, Space Sci. Rev., 92, 279 [NASA ADS] [CrossRef] [Google Scholar]
 Benz, W., & Asphaug, E. 1999, Icarus, 142, 5 [NASA ADS] [CrossRef] [Google Scholar]
 Bukhari Syed, M., Blum, J., Wahlberg Jansson, K., & Johansen, A. 2017, ApJ, 834, 145 [NASA ADS] [CrossRef] [Google Scholar]
 Chambers, J. E. 2014, Icarus, 233, 83 [NASA ADS] [CrossRef] [Google Scholar]
 Fortier, A., Benvenuto, O. G., & Brunini, A. 2009, A&A, 500, 1249 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Greenberg, R., Bottke, W. F., Carusi, A., & Valsecchi, G. B. 1991, Icarus, 94, 98 [NASA ADS] [CrossRef] [Google Scholar]
 Guilera, O. M. 2016, Boletin de la Asociacion Argentina de Astronomia La Plata Argentina, 58, 316 [NASA ADS] [Google Scholar]
 Guilera, O. M., & Sándor, Z. 2017, A&A, 604, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Guilera, O. M., Brunini, A., & Benvenuto, O. G. 2010, A&A, 521, A50 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Guilera, O. M., de Elía, G. C., Brunini, A., & Santamaría, P. J. 2014, A&A, 565, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hayashi, C. 1981, Prog. Theor. Phys. Suppl., 70, 35 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Helled, R., Bodenheimer, P., Podolak, M., et al. 2014, Protostars and Planets VI, (Tucson: University of Arizona Press), 643 [Google Scholar]
 Hori, Y., & Ikoma, M. 2010, ApJ, 714, 1343 [NASA ADS] [CrossRef] [Google Scholar]
 Hori, Y., & Ikoma, M. 2011, MNRAS, 416, 1419 [NASA ADS] [CrossRef] [Google Scholar]
 Hubickyj, O., Bodenheimer, P., & Lissauer, J. J. 2005, Icarus, 179, 415 [NASA ADS] [CrossRef] [Google Scholar]
 Ida, S., & Makino, J. 1993, Icarus, 106, 210 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Ikoma, M., Nakazawa, K., & Emori, H. 2000, ApJ, 537, 1013 [NASA ADS] [CrossRef] [Google Scholar]
 Inaba, S., & Ikoma, M. 2003, A&A, 410, 711 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Inaba, S., Tanaka, H., Nakazawa, K., Wetherill, G. W., & Kokubo, E. 2001, Icarus, 149, 235 [NASA ADS] [CrossRef] [Google Scholar]
 Inaba, S., Wetherill, G. W., & Ikoma, M. 2003, Icarus, 166, 46 [NASA ADS] [CrossRef] [Google Scholar]
 Johansen, A., Blum, J., Tanaka, H., et al. 2014, Protostars and Planets VI, (Tucson: University of Arizona Press), 547 [Google Scholar]
 Johnson, T. V., & Lunine, J. 2005, in Lunar and Planetary Inst. Technical Report, Vol. 36, eds. S. Mackwell, & E. Stansbery, 36th Annual Lunar and Planetary Science Conference [Google Scholar]
 Jutzi, M., Michel, P., Benz, W., & Richardson, D. C. 2010, Icarus, 207, 54 [NASA ADS] [CrossRef] [Google Scholar]
 Jutzi, M., Holsapple, K., Wünneman, K., & Michel, P. 2015, ArXiv eprints [arXiv:1502.01844] [Google Scholar]
 Kobayashi, H., Tanaka, H., & Krivov, A. V. 2011, ApJ, 738, 35 [NASA ADS] [CrossRef] [Google Scholar]
 Kokubo, E., & Ida, S. 1998, Icarus, 131, 171 [NASA ADS] [CrossRef] [Google Scholar]
 Kokubo, E., & Ida, S. 2000, Icarus, 143, 15 [NASA ADS] [CrossRef] [Google Scholar]
 Kokubo, E., & Ida, S. 2002, ApJ, 581, 666 [NASA ADS] [CrossRef] [Google Scholar]
 Kruijer, T. S., Burkhardt, C., Budde, G., & Kleine, T. 2017, Proc. Natl. Acad. Sci., 114, 6712 [NASA ADS] [Google Scholar]
 Lambrechts, M., & Johansen, A. 2012, A&A, 544, A32 [Google Scholar]
 Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Leinhardt, Z. M., & Stewart, S. T. 2012, ApJ, 745, 79 [NASA ADS] [CrossRef] [Google Scholar]
 Lodders, K. 2003, ApJ, 591, 1220 [NASA ADS] [CrossRef] [Google Scholar]
 McDonnell, J. A. M., Evans, G. C., Evans, S. T., et al. 1987, A&A, 187, 719 [NASA ADS] [Google Scholar]
 Mizuno, H. 1980, Prog. Theor. Phys., 64, 544 [Google Scholar]
 Morbidelli, A., Bottke, W. F., Nesvorný, D., & Levison, H. F. 2009, Icarus, 204, 558 [NASA ADS] [CrossRef] [Google Scholar]
 Ormel, C. W., & Klahr, H. H. 2010, A&A, 520, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ormel, C. W., & Kobayashi, H. 2012, ApJ, 747, 115 [NASA ADS] [CrossRef] [Google Scholar]
 Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62 [NASA ADS] [CrossRef] [Google Scholar]
 Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
 Stern, S. A., McKinnon, W. B., & Lunine, J. I. 1997, On the Origin of Pluto, Charon, and the PlutoCharon Binary, eds. S. A. Stern & D. J. Tholen, (Tucson: University of Arizona Press), 605 [Google Scholar]
 Stevenson, D. J. 1982, Planet. Space Sci., 30, 755 [NASA ADS] [CrossRef] [Google Scholar]
 Venturini, J., Alibert, Y., Benz, W., & Ikoma, M. 2015, A&A, 576, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Venturini, J., Alibert, Y., & Benz, W. 2016, A&A, 596, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Weidenschilling, S. J. 2011, Icarus, 214, 671 [NASA ADS] [CrossRef] [Google Scholar]
 Youdin, A. N., & Lithwick, Y. 2007, Icarus, 192, 588 [NASA ADS] [CrossRef] [Google Scholar]
All Tables
Free parameters that determine for different types of materials and different impact velocities from Benz & Asphaug (1999) and Benz (2000).
All Figures
Fig. 1 Catastrophic disruption thresholds for basalt targets at impact velocities of 20–30 m s^{−1} (Benz 2000), 3 and 5 km s^{−1} (Benz & Asphaug 1999). The blue points correspond to the discrete data extracted from Benz (2000). 

In the text 
Fig. 2 Catastrophic disruption thresholds for icy targets at impact velocities of 0.5 and 3 km s^{−1} (Benz & Asphaug 1999). 

In the text 
Fig. 3 Core masses (solid lines) and envelope masses (dashed lines) as function of time. Here, as in Figs. 5, 8, 9, 12, and 13, the red thick lines correspond to the baseline model, while the gray thin lines correspond to the formation of the planet without planetesimal fragmentation. 

In the text 
Fig. 4 Left panel: time evolution of the radial profiles of the eccentricity (solid lines) and inclination (dashed lines) for 100 kmsized planetesimals. Right panel: time evolution of the relative velocity between planetesimals with 100 km radius. In brackets we show the total mass of the planet at the given times. The black curve corresponds to the time (and mass) at which the planet reaches the crossover mass. 

In the text 
Fig. 5 Coremasses (solid lines) and envelope masses (dashed lines) as a function of the time. Gray lines: model without planetesimal fragmentation. Red lines: baseline model. Pink lines: model in which is a functionof relative velocities. 

In the text 
Fig. 6 Top panel: solid accretion rates for basaltic planetesimals and small particles of various radii, products of the collisional evolution of planetesimals, for the baseline model. Bottom panel: same as top panel but for the model that considers as a function of the planetesimal relative velocities. 

In the text 
Fig. 7 Time evolution of the relative velocity between basaltic planetesimals of the same size for the model that considers as a function of the planetesimal relative velocities. Top panel: 1 msized planetesimals. Bottom panel: 1 kmsized planetesimals. In brackets we show the total mass of the planet at the given times. 

In the text 
Fig. 8 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model with no planetesimal fragmentation. Red lines: baseline model. Pink lines: model with for ices at impact velocities of 3 km s^{−1}. Blue lines: model where is calculated by a linear combination between for basalts at 3 km s^{−1} and for ices at 3 km s^{−1}. 

In the text 
Fig. 9 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model with no planetesimal fragmentation. Red lines: baseline model. Pink lines: model with for basalts at impact velocities of 3 km s^{−1} and where small particles with Stokes numbers less than or equal to unity are accreted as pebbles of second generation. 

In the text 
Fig. 10 Ratio between the Hill radius and the scale height of small particles, with Stokes number less than or equal to unity, for the baseline model (solid lines) and for the simulation including pebble accretion (dashed lines) as a function of time. 

In the text 
Fig. 11 Solid accretion rates for small particles, product of the collisional evolution of planetesimals, for the baseline model (solidlines) and for the simulation including pebble accretion (dashed lines) as a function of time. Small particles of these sizes have Stokes number less than or equal to unity. 

In the text 
Fig. 12 Coremasses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model without planetesimal fragmentation. Red lines: baseline model. Pink lines: model with all the improvements in the calculation of (the dependencies of with the planetesimal impact velocities and compositions) and the pebble accretion rates. 

In the text 
Fig. 13 Core masses (solid lines) and envelope masses (dashed lines) as a function of time. Gray lines: model without planetesimal fragmentation. Black lines: model without planetesimal fragmentation but considering grain opacities reduction (G.O.R). Red lines: baseline model. Pink lines: model taking into account all the improvements in the calculation of (the dependencies of with the planetesimal impact velocities and compositions), pebble accretion rates, and considering reduced grain opacities. 

In the text 
Fig. 14 Envelope mass as a function of the core mass for the two cases where planetesimal fragmentation is not considered. We plot the core mass until 10 M_{⊕} in order to show more clearly the differences at lowmass core values. The black line represents the case where grain opacities are not reduced, while the red line represents the case where grain opacities are reduced in the planet envelope. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.