Planet formation with envelope enrichment: new insights on planetary diversity
AAstronomy & Astrophysics manuscript no. venturini2016_arxiv c (cid:13)
ESO 2018June 27, 2018
Planet formation with envelope enrichment: new insights onplanetary diversity
Julia Venturini (cid:63) , Yann Alibert, and Willy Benz
Physikalisches Institut, Universität Bern, Sidlerstrasse 5, 3012 Bern, Switzerlande-mail: [email protected]
June 27, 2018
ABSTRACT
Aims.
We compute, for the first time, self-consistent models of planet growth including the e ff ect of envelope enrichment. The changeof envelope metallicity is assumed to be the result of planetesimal disruption or icy pebble sublimation. Methods.
We solve internal structure equations taking into account global energy conservation for the envelope to compute in-situplanetary growth. We consider di ff erent opacities and equations of state suited for a wide range of metallicities. Results.
We find that envelope enrichment speeds up the formation of gas giants. It also explains naturally the formation of low andintermediate mass objects with large fractions of H-He ( ∼
20 - 30 % in mass). High opacity models explain well the metallicity of thegiant planets of the solar system, whereas low opacity models are suited for forming small mass objects with thick H-He envelopesand gas giants with sub-solar envelope metallicities. We find good agreement between our models and the estimated water abundancefor WASP-43b. For HD 189733b, HD 209458b and WASP-12b we predict fractions of water larger than what is estimated fromobservations, by at least a factor ∼ Conclusions.
Envelope enrichment by icy planetesimals is the natural scenario to explain the formation of a large variety of objects,ranging from mini-Neptunes, to gas giants. We predict that the total and envelope metallicity decrease with planetary mass.
Key words.
Planet formation, exoplanet atmospheric composition.
1. Introduction
The core accretion model (Perri & Cameron 1974; Mizuno 1980;Bodenheimer & Pollack 1986; Pollack et al. 1996) is the mostaccepted scenario to explain the formation of a vast diversity ofplanets (e.g, Alibert et al. 2005; Mordasini et al. 2009; Guileraet al. 2011). The central idea of the core accretion model can besummarised as follows. First, a solid core must be formed fromthe accretion of planetesimals / pebbles. Once this core reachesapproximately a lunar mass, the core gravity is strong enoughto start to bind some gas from the protoplanetary disk. Thus,from this stage on, the protoplanet keeps growing by accret-ing both solids (planetesimals / pebbles) and gas (basically H-He). Planet formation models typically assume, for simplifica-tion, that solids and gas do not mix: all the solids deposit theirmass and energy at the top of the core, and the primordial H-Heis collected above, building the atmosphere (or envelope). Thisis, of course, a very strong and unrealistic simplification: bolidesthat traverse a gaseous atmosphere undergo thermal ablation andmechanical breakup. Hence, volatile material can vaporise andmix with the primordial H-He, changing the composition of theenvelope during the formation of a planet.If planetesimal / pebble disruption did not occur during for-mation, then the envelope metallicities of planets should berather sub-stellar, because the gas accreted into the planetsshould in principle be metal-poor compared to the central star(the metals condense to form planetesimals / pebbles). This is notwhat is observed in the solar system, where the giant plants showsome level of envelope enrichment (Irwin et al. 2014; Niemann (cid:63) Currently at the Institute for Computational Science, University ofZurich, Winterthurerstrasse 190, 8057 Zurich, Switzerland et al. 1998; Guillot & Gautier 2014). The alternative hypothe-sis to planetesimal / pebble dissolution for explaining envelopesenriched in heavy elements, is core erosion (Wilson & Militzer2012). From an energetic point of view, it is possible to mix partof the core upwards (Guillot et al. 2004). In addition, core mate-rial is miscible in metallic hydrogen (Wilson & Militzer 2012),which allows for the heavy elements (if they can be lifted up)not to settle to re-form a core. Regarding the mixing of an initialcore within the gaseous envelope, Vazan et al. (2015) showedthat this process is favoured if an initial compositional gradientexists in the interior of the planet. Thus, the mixing of heavyelements on the planetary envelope seems to be more likely ifthe heavy elements are not initially concentrated in well-definedcore. The formation of such a di ff use core requires planetesi-mal dissolution in the deep envelope during the formation of theplanet. Hence, even if core erosion could play a relevant role inmixing heavy elements in the planetary envelope, this processseems to demand as well that the envelope is initially enrichedby planetesimal / pebble dissolution.The problem of considering an envelope that is enriched withrespect to stellar values during the formation of a giant planethas been raised since the very early studies of planet formation.Already in 1986, Bodenheimer & Pollack (1986) mentioned,among other two problems that remained to be solved “the factthat the molecular weight of the envelope is expected to increasewith time as some of the icy planetesimals dissolve in it", andadded that this problem “could significantly change the accre-tion scenario". Indeed, Stevenson (1982) showed that the criticalcore mass (the mass required to trigger rapid accretion of gas)is reduced when the envelope mean molecular weight increases.Moreover, Wuchterl (1993) showed a direct dependence of the Article number, page 1 of 14 a r X i v : . [ a s t r o - ph . E P ] S e p & A proofs: manuscript no. venturini2016_arxiv critical core mass with the adiabatic gradient. The adiabatic gra-dient is expected to decrease when chemical reactions take place.This necessarily occurs when abundant elements such as H, C,O are present in the envelope. Therefore, the e ff ect of pollutingthe primordial envelopes reduces the critical core mass not onlydue to the increase of mean molecular weight but also via a re-duction of the adiabatic gradient (Hori & Ikoma 2011). Anothere ff ect that can reduce the adiabatic gradient is the condensationof species. We showed (Venturini et al. 2015) that if a planetforms in cold regions of the disk, such that the temperatures andpressures are low enough for certain species to condense in theatmosphere of the protoplanet, then this e ff ect leads to an evenlarger reduction of the critical core mass.Despite its expected importance for the formation of giantplanets, the e ff ect of envelope enrichment has so far never beenimplemented in a self-consistent way in any evolutionary calcu-lation of planet formation. This is of course a consequence of thedi ffi culty in modelling all the processes involved, which include:planetesimal thermal ablation and dynamical breakup, mass andenergy deposition in di ff erent envelope layers, and the inclusionof a self-consistent change in the envelope’s microphysics (opac-ities and equation of state) as the envelope metallicity evolves.The magnitude of the enrichment depends on the accretionrate of solids, and on their size and strength properties. Thesmaller and more porous the bolide, the easier it is to disrupt itand mix it with the envelope gas when crossing the atmosphere(Podolak et al. 1988). This tells us that the e ff ect of envelopeenrichment is more relevant for smaller bolides. Hence, whengrowing planets from cm-m size particles (the so-called pebbles ,Lambrechts & Johansen 2012; Lambrechts et al. 2014), includ-ing this e ff ect is necessary.In this work we compute, for the first time, the in-situ growthof a planet taking into account envelope enrichment by icy plan-etesimals / pebbles. Given the uncertainties in the initial size dis-tribution of planetesimals, we test the e ff ect of envelope enrich-ment corresponding to a wide range of particle sizes. We solveinternal structure equations using global energy conservation ar-guments, and taking into account the changes of envelope metal-licity in the opacities and equation of state. In Sect. 2 we explainthe numerical method, Sect. 3 is devoted to discussing our as-sumptions, in Sect. 4 and 5 we show and analyse our results, inSect. 6 and 7 we discuss our main results and predictions. Wesummarise our conclusions in Sect. 8.
2. Methodology
We use the numerical method developed in Venturini et al.(2015) to solve the equations of hydrostatic equilibrium,mass conservation and heat transport. We apply the usualSchwarzschild criterion to distinguish between radiative andconvective layers, adopting an adiabatic gradient for the lattercase.Regarding energy conservation, we use global energy con-servation arguments (Mordasini et al. 2012; Fortier et al. 2013;Piso & Youdin 2014) to find the total luminosity radiated awayby the protoplanet. In the following subsection we show how wecompute this total luminosity.
To explain how we compute the total luminosity ( L ), let us anal-yse the total energy of the system at time t and at time t + dt , asit is sketched in Fig. 1. Global energy conservation time t time t+dt Ldt E tot ( t ) = E g,env ( t ) + E i,env ( t ) + E g,core ( t )+ e gas,acc dm l vap dm z c p Tdm z E tot ( t + dt ) = E g,env ( t + dt ) + E i,env ( t + dt )+ E g,core ( t + dt ) L = dE g,core dt dE env dt + e gas,acc ˙ M gas (1 )( l vap + c p T ) ˙ M z P disk dVdt E tot ( t + dt ) E tot ( t ) = Ldt P disk dV Fig. 1: Sketch representing the total energy of the system at time t and t + dt . See main text for the explanation of the di ff erentterms.At a given time t , the total energy of the system is given bythe sum of : • the total energy of the envelope (gravitational plus internal): E g,env ( t ) + E i,env ( t ) • the gravitational potential energy of the core (we neglect heatsources coming from the core ): E g,core ( t ) • the total energy of the gas layer that will be accreted, repre-sented by the term e gas,acc dm , where e gas,acc is the total energyof the layer per unit mass, and dm the mass of gas that willbe accreted • the energy that needs to be subtracted from the envelope inorder to vaporise and thermalise the ices of the planetesi-mals that will be mixed in the envelope. The vaporisationis given by the term − l vap dm z , l vap being the latent heat ofvaporisation of the ice per unit mass, and dm z the mass ofice remaining in the envelope . The thermalisation term isgiven by c p ∆ T dm z , being c p the specific heat capacity of theice and ∆ T the change of temperature from 0 ◦ C to the meanenvelope temperature.At time t + dt , all the mass of the planetesimals and gas that wereaccreted in the time elapsed ( dt ) are now part of the protoplanet,so the total energy is just the gravitational potential of the coreand envelope, plus the internal energy of the envelope.Our protoplanet is not a closed system, it is embedded in aprotoplanetary disk, which pushes the outermost layer of the pro-toplanet, making work represented by − P disk dV . Also, the pro-toplanet cools, radiating energy away by photons ( − Ldt ). Hence, E tot ( t + dt ) − E tot ( t ) = − Ldt − P disk dV . (1) The planetesimals are assumed to have started with zero velocity atinfinity, therefore their total energy is zero. The release of their gravita-tional potential energy, which heats the envelope, appears as the term ofchange of gravitational potential energy of the core, as we show later inthe main text. The change of E g,core ( t ) with time has typically, minimum values of10 erg / s, whereas the change of E i,core ( t ) with time reaches values of10 erg / s (Lopez & Fortney 2014) and is therefore negligible. The cool-ing of the core should, on the other hand, be considered once the heatdue to planetesimal bombardment ( dE g,core / dt ) ceases, as is the case, forexample, in evolutionary simulations. We assume the ices to be water ice. We explain better this assumptionin Sect. 3.1 and Sect. 3.2.2.Article number, page 2 of 14enturini et al.: Planet formation with envelope enrichment
Regrouping terms, we find L = − dE g,core dt − dE env dt + e gas,acc ˙ M gas − P disk dVdt − (1 − β )( l vap + c p ∆ T ) ˙ M z . (2)The first term is just the change of gravitational potential en-ergy of the core, the second is the change of total energy of theenvelope, the next two are the surface terms already explained,and the last are the “ice heating" terms, also explained above. Wecall β the mass fraction of planetesimals that sinks to the core,that is why a factor 1 − β (what remains in the envelope) is neededin front of the “ice heating" terms.Now let us further develop the first two terms, which are al-ways at least two orders of magnitude larger than the others. Tocompute the first term (which we will call hereafter L core ), wemake use of the fact that the gravitational potential energy of asphere of uniform density is: E core = − GM R core . (3)Di ff erentiating with respect to time and using the fact thatthe core density is constant , one can readily find: L core = − dE g,core dt = GM core ˙ M core R core . (4)If we call ˙ M Z the total accretion rate of solids, since β isthe mass fraction of planetesimals / pebbles that sinks to the core,then ˙ M core = β ˙ M Z .We express the second term as L gas : L gas = − E env ( t + dt ) − E env ( t ) dt . (5)The problem with the computation of L gas is that the energyof the envelope at time t + dt is not known before obtaining thestructure at t + dt . Therefore, we must make a guess of E env ( t + dt )to be able to have L ( t + dt ) for computing the structure at t + dt .This guess is performed following Mordasini et al. (2012) andFortier et al. (2013). These authors assumed a similar functionalform for the total energy of the envelope as the one of the gravi-tational potential energy of the core of uniform density. Hence, E env = − k env M env g , (6)where g is a mean gravity, defined by g = G ( M core / R core + M P / R P ).For a given converged structure at time t , we know the totalenergy of the envelope at this time. Hence, Eq.(6) defines k env attime t , which is then used for the total energy of the envelopeat time t + dt . This assumption is quite good during most of thegrowth of the planet, since k env does not practically vary fromone timestep to another. Indeed, to test this we have run sim-ulations using a predictor-corrector scheme: we first computethe structure of the envelope using the k env from the previoustimestep, k env (t - dt), then compute the resulting k env of the ac-tual timestep, k env, corr , and finally recompute the structure of theplanet at time t using k env = k env (t - dt) + k env, corr ). We assume, throughout this work, a constant density for the core,an usual practice in most planet formation works. If the core were com-pressible, an extra term in L core should be included, reflecting the changeof core gravitational energy associated with the change of density. Nev-ertheless, since the core is in solid phase, that term is not really relevant. For the initial conditions, we assume that the gas in the disc -and therefore the planetary envelope- is made of hydrogen andhelium. This is justified by the fact that we are forming planetsbeyond the iceline, and thus, most of the metals have alreadycondensed into solids , which exist in the disc in the form ofembryos, planetesimals or pebbles.We start all simulations with M core = ⊕ . The coregrows at a given accretion rate ˙ M Z (see Sect. 3.3 for the dif-ferent schemes adopted to compute ˙ M Z ). When the core reachesa threshold value of M thresh , we assume that the core keeps grow-ing by an amount ∆ M core = β ˙ M Z dt while the amount ∆ M Z,env = (1 − β ) ˙ M Z dt remains uniformly mixed in the envelope. The de-pendence of M thresh on the planetesimals’ properties is discussedin Sect. 3.5.At the beginning of each timestep, an initial guess of enve-lope metallicity is assumed in order to compute the correspond-ing equilibrium envelope mass (the envelope mass depends onthe core mass and radius, total luminosity and envelope metal-licity). For the computation of M env we use the second numericalscheme described in Venturini et al. (2015): iteration on M env fora given core mass and radius, and Z env . The total luminosity isrecalculated at each iteration on M env , because L gas depends on M env , as explained in the above subsection. The envelope metal-licity is taken as uniform throughout the envelope and is definedas: Z env = M Z,env M env , (7) M Z,env being the mass of volatiles that remains mixed in the en-velope.Once convergence in M env is achieved, the envelope metal-licity is updated to the new M env found. Afterwards, an extraiteration on Z env is required in order to have the metallicity self-consistently defined (and thus, mass conservation of metals guar-anteed).
3. Assumptions β We recall that β is the mass fraction of the incoming planetesi-mals that is assumed to reach the core (hence, the fraction 1 − β remains homogeneously mixed in the envelope). For our nom-inal model, we take β = .
5. This choice is based on the re-fractory / volatile content of comets (Lambrechts et al. 2014; Thi-abaud et al. 2015). Of the material falling in, we assume thatonly ices mix in the envelope and that the refractory sinks intothe core. We assume for simplification that the ices are just wa-ter ice, because it is the main volatile species in comets (Kofmanet al. 2015), and because of self-consistency with our equation ofstate (see Sect. 3.2.2). The choice of β = . Strictly speaking, volatile species with condensation temperaturesbelow the one of water (e.g, CO) could exist in the gas phase beyondthe iceline. Nevertheless, the mass contribution of these species to thetotal gaseous disc is expected to be negligible (Thiabaud et al. 2015).Article number, page 3 of 14 & A proofs: manuscript no. venturini2016_arxiv
Opacities Crossover Mass [M ⊕ ]Semenov + Ferguson 17Mordasini + Freedman 6.9Table 1: Crossover mass for Z = Z (cid:12) , di ff erent opacities and ˙ M Z = − M ⊕ / yr.these ranges do not necessary cover all the pressure-temperaturevalues existent during formation (see Fig. 5 of Venturini et al.2015), which can extend from P ∼ − Pa to P ∼
10 GPa andfrom T ∼
100 K to T ∼ In Venturini et al. (2015) we assumed a total opacity resultingfrom contributions from dust and gas. We used dust opacitiesfrom Semenov et al. (2003) and a gas opacity from tables com-puted by Ferguson for arbitrary metallicity, based on the originalcalculations of Alexander & Ferguson (1994). The dust opac-ities from Semenov et al. (2003) are suited for protoplanetarydisks, and therefore do not consider grain growth and settling ina planetary envelope. This is taken into account in the analyticalopacities calculated by Mordasini (2014), which we adopt in thiswork. The dust opacities of Mordasini (2014) have the additionaladvantage of being independent on the accretion rate of metals(an equilibrium among the deposition of grains, grain growth andsettling was found in those calculations), making these opacitiesindependent of the envelope metallicity and therefore, suited forany value of Z env .Regarding gas opacities, new calculations by Freedman et al.(2014) are suited for planetary envelopes and for a wide rangeof metallicities ( 0 (cid:46) Z env (cid:46) . ≤ µ bar to 300 bars and75 K to 4000 K. These ranges do not reach the high tempera-tures and pressures that are typical at the surface of the core, butin practice, in that domain of high pressure and temperature theenvelopes are always convective. Therefore, even taking a con-stant extrapolation in these regions does not a ff ect the internalstructure computation.In any case, it is well known that opacities a ff ect drasticallythe timescale to form a giant planet (Pollack et al. 1996; Hu-bickyj et al. 2005; Ikoma et al. 2000). This is why we run com-parison tests using di ff erent extreme combinations of dust andgas opacities. On one side, we take the old combination of Se-menov + Ferguson, and on the other, Mordasini + Freedman.The results concerning the crossover mass are summarised inTable 1 for a solar composition envelope and a constant accre-tion rate of metals of 10 − M ⊕ / yr . The notorious decrease in thecrossover mass when using the new opacities is mainly due tothe fact that the dust opacities from Mordasini (2014) are muchlower than the ones from Semenov et al. (2003) (whose valuesare very similar to the interstellar ones). We repeat this compar-ison test taking into account envelope enrichment in Sect. 4.4.For all results shown before, we adopt the new opacities of Mor-dasini for the dust and Freedman for the gas. Initial solid surface density ( Σ ) 4 g / cm Planetesimal density ( ρ p ) 0.92 g / cm Planetesimal radius ( r p ) 100 kmTable 2: Parameters used for results where the Pollack-schemeis implemented and nominal opacities (Mordasini + Freedman)are adopted. When we use opacities of Semenov + Fergusonthe only parameter we change is the initial solid surface density,which is set at Σ =
10 g / cm . We assume in this work that the volatile content of the planetes-imals (or pebbles) that are destroyed while traversing the enve-lope is what remains well mixed in the envelope, thanks to icesublimation. We assume as well that this volatile-component isentirely made of water. This is justified mainly by two reasons.First, as we mentioned before, water is thought to be the mainvolatile molecule present in planetesimals. Second, there is noaccurate EOS available in the literature for an arbitrary mixtureof volatiles that covers the large ranges of temperature and pres-sure present in planetary interiors. The only EOS suited for thispurpose is that of water.Thus, we adopt an EOS for a mixture of H, He and H Owhich takes into account degeneracy due to free electrons. Forthe H-He component we implement the Saumon et al. (1995)EOS, and for the H O component we use an improved version ofANEOS (Thompson 1990). An important drawback in the stan-dard version of ANEOS is the assumption of the substance to bemonoatomic in the gas phase. We have corrected this for water(Benitez et al. 2016, in preparation) to include the proper degreesof freedom, following the approach of Melosh (2007). We haveimplemented this new version of ANEOS for our water compo-nent. The EOS of the mixture of H-He with water is obtained,as in Bara ff e et al. (2008), by means of the additive volume rule,which has been proven to yield adequate results for mixtures ofH, He and H O (Soubiran & Militzer 2015).
We adopt di ff erent models for the accretion rate of solids:1. Constant accretion rate.
This basic scheme allows us toanalyse the e ff ect of envelope enrichment in its simplestform. For nominal results, we adopt a ˙ M Z = − M ⊕ / yr.Other values of accretion rates are tested in Sect. 4.5.2. Accretion rate a la Pollack . We implement our enrichmentcode in a Pollack-scheme (Pollack et al. 1996) planetesimalaccretion rate (see Sect. 5 for an explanation). The purposeof this is to show the e ff ect of envelope enrichment in a morerealistic formation scenario (although still not very accurate,since the proper excitation of the planetesimals in the oli-garchic growth regime is not included as it is in Fortier et al.2007, 2013). The crucial parameters assumed for this sce-nario are summarised in Table 2. The di ff erent choices of Σ for the di ff erent sets of opacities are set to obtain simi-lar formation timescales, as is explained later in Sect.5. Forthe capture radius we use the prescription of Inaba & Ikoma(2003) (see as well, Guilera et al. 2011). Article number, page 4 of 14enturini et al.: Planet formation with envelope enrichment a T out
150 K P out / cm ρ core / cm Table 3: Boundary conditions used in the nominal model.
For the definition of the radius of the protoplanet, we use the oneproposed by Lissauer et al. (2009), which takes into account flowcirculation from 2D hydrodynamical simulations: R P = GM P ( c s + GM P / R H ) , (8) R P being the planet radius, M P the planet mass, c s the soundspeed of the gas and R H the Hill radius.For nominal results, we adopt boundary conditions suited forthe actual position of Jupiter. The values of these boundary con-ditions are given in Table 3. The results are very insensitive tothe boundary conditions (as long as the planet is beyond the ice-line and the accretion rate of solids is constant), this is why wedo not show results for other boundary conditions.Regarding the choice of core density, we set the conserva-tive number of 3.2 g / cm . This number should be larger once weconsider that just rocky material sinks to the core, but we havetested that increasing the core density has a very small impact onthe overall evolution. M thresh The amount of material that is deposited in the envelope whena planetesimal is disrupted, is a complicated function that de-pends upon many planetesimal properties (mass, density, ma-terial strength, etc). Therefore, if one wants to compute self-consistently the mass deposition, an accurate knowledge of plan-etesimal properties is needed. Since this is still unknown, weencode this lack of information in the parameter M thresh , whichrepresents the minimum mass of the core which can bind an en-velope massive enough to destroy completely the incoming plan-etesimals. For nominal results, we take M thresh = ⊕ . We testthe e ff ect of considering, as well, M thresh = . ,
1, and 4 M ⊕ in Sect. 4.2. Just to have an order of magnitude in mind, a coremass of 1 M ⊕ would have an envelope massive enough to disruptcompletely, before reaching the core, icy planetesimals of ∼ ⊕ , icy planetesimals of ∼ −
10 km, and a core of 4 M ⊕ , icy planetesimals of ∼ − . A core of 0.5 M ⊕ has an envelope massequivalent to 2 Earth’s atmospheres in our models, so pebbleswould surely sublimate but not km-size planetesimals.
4. Results for constant accretion rate of solids
Figure 2 shows the growth of a planet assuming ˙ M Z = − M ⊕ / yr for a standard case where all the solids go to the core(envelope composed of H and He), and for a case where enve-lope enrichment is taken into account. Throughout this work, we these values would be upper limits, planetesimals smaller that thesewould be fully disrupted in the envelope as well. m a ss [ M ⊕ ] time [Myr]M core M HHe M env Fig. 2: Solid lines: enriched case with M thresh = ⊕ and β = t cross , Pollack et al. 1996) as the characteristic timeto form a giant planet, we find t cross = t cross = M thresh = ⊕ , half of the met-als (water in our assumptions) remains uniformly mixed in theenvelope and the other half is deposited in the core (see changeof slope in M core in Fig. 2). Since the envelope is still quite thinat this M core (see Table 4), the envelope metallicity rises veryrapidly at the beginning of the enrichment. In 3 . × yrs, Z env reaches its maximum and then dilution slowly starts. There is animportant fact to remark: rapid gas accretion is not triggered im-mediately as a consequence of envelope enrichment. This is wellillustrated in Fig. 4, where we observe that the accretion rate ofH-He remains lower than that of metals until t ∼ τ HHe = M HHe / ˙ M HHe ) also as a functionof time for the case of M thresh = ⊕ . We see that immediatelyafter the onset of enrichment, the timescale to accrete H-He isvery short, but this timescale starts to increase until it reaches amaximum at t ≈ . τ HHe starts to decrease, HHe-accretion accelerates. So wedefine the onset of the runaway of gas as the time when τ HHe reaches its maximum, and we denote this time as t run .The time elapsed between the onset of enrichment ( t ≈ Di ff erent definitions for the onset of runaway of gas have been givenin the literature. Some define it when τ HHe drops to a fraction of its max-imum(Piso & Youdin 2014; Lee et al. 2014). Qualitatively, the crossovermass criterion (Pollack et al. 1996) works to infer fast accretion of gas,but the truth is that after τ HHe starts to decrease, runaway of gas is in-evitable (as long as there is gas left on the disk, of course).Article number, page 5 of 14 & A proofs: manuscript no. venturini2016_arxiv Z en v time [Myr] Fig. 3: Evolution of the envelope metallicity for the enrichedcase of Fig. 2. -8-7.5-7-6.5-6-5.5-5-4.5-4 2 2.5 3 3.5 4 4.5 l og ( d M / d t ) [ M ⊕ / y r ] time [Myr]HHe acc. ratesolid acc. rate τ HH e [ M y r ] time [Myr] Fig. 4: Top: accretion rate of H-He (red) vs. accretion rate ofsolids (blue) as a function of the core mass for the enriched case.Bottom: timescale to accrete H-He for the enriched case. Notethat for M core (cid:38) ⊕ , the accretion rate of H-He dominatesthe one of solids, and also the timescale to accrete H-He starts todecrease. τ HH e [ M y r ] GCR = M env /M core enriched casenon-enriched case Fig. 5: Timescale to accrete H-He as a function of the gas-to-core ratio, for the enriched case (purple) and the non-enrichedone (green).This scenario implies that a very high gas-to-core ratio (GCR = M env / M core ) is achievable before runaway gas accretion is trig-gered. Fig. 5 depicts this more clearly. It shows the timescale toaccrete H-He as a function of the GCR. For the non-enrichedcase, the maximum of τ HHe occurs at GCR ∼ ffi cult to explain in the framework of the standardcore accretion model, where envelopes are assumed to remainnon-enriched (Lee & Chiang 2016). Nevertheless, if envelopeenrichment is taken into account, much larger CGR ( ∼ .
8) canbe achieved before the onset of runaway (Fig. 5).In this context we may wonder why are giant planets formedfaster when envelope enrichment is taken into account. When en-velope enrichment sets in, the structure of the envelope changesradically: the mean molecular weight increases very rapidly dueto the increase in Z env . Since the boundary conditions are thesame, an increase of the mean molecular weight translates intoan increase of the density profile. Therefore, the self-gravity ofthe envelope is stronger, and the planet becomes more prone tocontract and accrete more gas. M thresh In this section we test the e ff ect of changing M thresh . As we statedin Sect. 3.5, this would correspond to consider di ff erent sizes(and / or composition, but mainly sizes) for the disrupted plan-etesimals / pebbles. Table 4 summarises some aspects of the sim-ulations for values of M thresh = ⊕ . The first rowshows the mass of the envelope when envelope enrichment starts.This information is important when trying to link the value of M thresh with the corresponding maximum size of planetesimalsto be fully disrupted.In Fig. 6 (top) we show the evolution of the gas-to-core ra-tio (GCR) and of the total mass fraction of H-He for the di ff er-ent M thresh . The black dots indicate the time of the onset of therunaway of gas (maximum in τ HHe ) and the numbers above, thecorresponding value of GCR at this time.It is interesting to note that the smaller M thresh , the larger thegas-to-core ratio that can be achieved before runaway of gas. Article number, page 6 of 14enturini et al.: Planet formation with envelope enrichment M thresh [M ⊕ ] 0.5 1 2 4M env,thresh [M ⊕ ] 2.6 × − × − × − r P [km] (cid:46) M cross [M ⊕ ] 1.15 1.82 2.83 4.58 t cross [Myr] 1.8 2.65 3.68 5.16 t run [Myr] 2.67 2.9 3.48 4.83 ∆ t enriched [Myr] 2.17 1.9 1.48 0.83Table 4: Values of envelope mass at the onset of enrichment(M env,thresh ), approximate size of icy particles dissolved in theenvelope for the corresponding M thresh , crossover mass ( M cross ),crossover time ( t cross ), time at which runaway of gas begins ( t run )and time during which the envelope is enriched but still not inrunaway of gas ( ∆ t enriched ), for di ff erent values of M thresh . ˙ M Z = − M ⊕ / yr. G CR time [Myr]solid: GCR = M env /M core dashed: total mass fraction of H-Heenriched case, M thresh = 0.5enriched case, M thresh = 1enriched case, M thresh = 2enriched case, M thresh = 4non-enriched case G CR M core [M ⊕ ]solid: GCR = M env /M core dashed: total mass fraction of H-Heenriched case, M thresh = 0.5enriched case, M thresh = 1enriched case, M thresh = 2enriched case, M thresh = 4non-enriched case Fig. 6:
Top : evolution of the gas-to-core ratio and total mass frac-tion of H-He of the planets for di ff erent M thresh (in Earth masses,indicated with di ff erent colours at the inset of the figure). Theblack dots on each curve indicate the time when the runawayof gas starts, and the numbers above, the corresponding GCR atthat time. Bottom: same as above but as function of the mass ofthe core. M thresh [M ⊕ ] 0.5 1 2 4 ∞ GCR 1.56 1.16 0.82 0.54 0.10f
HHe Z env M core [M ⊕ ] 1.59 1.95 2.75 4.45 5.31 M P [M ⊕ ] 4.10 4.30 5.00 6.84 5.79Table 5: Values of gas-to-core ratio (GCR), total mass fraction ofH-He (f HHe , i.e 1- Z tot ), Z env , core mass ( M core ) and total mass ofthe planet ( M P ) at the onset of the runaway of gas ( t = t run ) fordi ff erent M thresh ( M thresh = ∞ corresponds to the non-enrichedcase). ˙ M Z = − M ⊕ / yr.Also interestingly, the period during which enrichment takesplace but fast accretion of gas still does not occur, tends to belonger the lower M thresh (last row of Table 4). This means thatthe smaller the planetesimals, the more likely it is that at thetime the disk disappears, a small planet with high GCR remains.Another surprising fact is that the total mass fraction of H-He(which we denote f HHe ) at the onset of runaway of gas is always ∼
30% for the enriched cases (see Table 5), despite the value of M thresh . That value would be the maximum attainable for planetsthat do not become gas giants, smaller values of f HHe are alsopossible, as Fig. 6 depicts. This means that envelope enrichmentseems to be a natural way to explain the formation of low mass,low density planets (mini-Neptunes) and also the recently called“super-pu ff s": light planets that have a bulk composition of H-Heof presumably (cid:38)
20 % of H-He by mass (Lee & Chiang 2016;Lopez & Fortney 2014). We want to remark that the claim ofenvelope enrichment being a more natural scenario to explainthe formation of objects with ∼
20 % of H-He is not just due tothe fact that these values are reached before the onset of runawayof gas, but also, because these high amounts of H-He last longer when envelope enrichment is included, as Fig.6 shows.An interesting aspect in Table 5 is that the total mass ofthe planet at the onset of runaway of gas grows with increas-ing M thresh but decreases for M thresh = ∞ (the non-enriched case).This is related to the fast accretion of H-He that is triggered onceenrichment sets in. If we analyse the core masses of the di ff er-ent M thresh at the onset of runaway (Table 5), we note that for M thresh = ∞ the core is larger than for M thresh = ⊕ . However,because the mass fraction of H-He is much larger in the enrichedcases that in the non-enriched one, the total planetary mass of the M thresh = ⊕ case is larger than the M thresh = ∞ one.Finally, in Fig. 7, we show the total mass of solids ofthe planet as a function of the planetary mass acquired dur-ing growth. Note that with the low opacities of our nominalmodel (Mordasini + Freedman), even in the extreme case of non-enrichment, the maximum mass of solids that can be attainedduring growth is M Z ≈ ⊕ . Larger accretion rates can lead tolarger core masses (or total mass of solids, in general). We showthis in Sect. 4.5. β One could think that since in reality there should be a distribu-tion of planetesimal’s sizes, then there is not one exact value of M thresh for which the envelope starts to get enriched, but rathera range of M thresh . We have already tested the e ff ect of adopting Note that the “total metallicity" of the planet is Z tot =
1- f
HHe = ( M core + M Z,env ) / M P Article number, page 7 of 14 & A proofs: manuscript no. venturini2016_arxiv M Z [ M ⊕ ] M planet [M ⊕ ]enriched case, M thresh = 0.5enriched case, M thresh = 1enriched case, M thresh = 2enriched case, M thresh = 4non-enriched case Fig. 7: Total mass of solids ( M Z = M core + M Z,env ) as a functionof total planet mass ( M planet = M Z + M HHe ) during formation.Note that simulations typically stop when M Z becomes asymp-totically flat with planet mass. This flatness occurs because theaccretion rate of H-He at this stage is much higher than the oneof solids (usually by 2 orders of magnitude). We will use thisfact to extrapolate total metallicities in Fig. 12 (just by fixing thefinal value of M Z in the definition of Z tot ).di ff erent M thresh in the previous section, but considering a distri-bution of planetesimal sizes also means that β would not changeabruptly at a given M thresh , but it would transition for a range of M thresh .We perform, therefore, the following test. We assume β = M core ≤ ⊕ , β = . M core ≥ ⊕ , and a linear varia-tion of β in between. The growth of the planet is shown in Fig.8, as well as the evolution of envelope metallicity. Note that de-spite the smoother transition in β than in the previous cases, theevolution is very similar to the case M thresh = ⊕ . The enve-lope metallicity still grows quite abruptly, because the envelopeat M core = ⊕ is quite thin. Since the envelope at M core = ⊕ is thinner than at M core = ⊕ , in this case of smoothertransition of β , Z env reaches a higher maximum than in the case M thresh = ⊕ . But the overall formation time is practically thesame as the case M thresh = ⊕ . Therefore we conclude thatassuming a sharp transition of β at a given M core is not a relevantsimplification. As we mentioned in Sect. 3.2.1, the timescale to form giant plan-ets depends strongly on the choice of opacities. Despite the factthat in our nominal results we used the latest opacities published,it could be that the opacities do not follow exactly those low val-ues. For instance, the recondensation of upstreaming gas intograins would increase the opacities, and this is not taken intoaccount for the computation of Mordasini (2014) nor Freedmanet al. (2014).Just to test how choosing larger opacities would a ff ect ourresults, we run a simulation with the other extreme set of opac-ities that we introduced in Sect. 3.2.1: dust opacities from Se-menov et al. (2003) and gas opacities from Ferguson (based onAlexander & Ferguson (1994), but including arbitrary metallici-ties). We consider here the enriched case with ˙ M Z = − M ⊕ / yrand M thresh = ⊕ . m a ss [ M ⊕ ] time [Myr]M core M HHe M env Z en v time [Myr] 0.5 0.6 0.7 0.8 0.9 1 β Fig. 8: Evolution of masses and metallicity for a test case wherewe assume β varying linearly between 1 and 0.5 for core massesranging between 1 and 2 M ⊕ . This smooth change of β is shownin color-bar with the evolution of Z env .The results of the growth of the planet are shown in Fig. 9.The total mass of the planet is plotted with a color-bar that in-dicates the timescale to accrete hydrogen and helium ( τ HHe ). Wenote that the maximum of this timescale in this case is reachedafter crossover mass. It occurs for a total mass of ∼
12 M ⊕ , ofwhich M core ∼ ⊕ and M env ∼ . ⊕ . At that time, halfof the mass of the envelope is H-He and the other half, water.It means that when τ HHe reaches it maximum ( t ≈ Z tot ∼ M Z = − M ⊕ / yr we chose ˙ M Z = × − M ⊕ / yr, then t run ≈ t = t run , the total planetary mass is that of Neptune, with a totalmetallicity of Z tot ∼ / Uranus-like planet.
Article number, page 8 of 14enturini et al.: Planet formation with envelope enrichment m a ss [ M ⊕ ] time [Myr]M planet M core M env M HHe M z,env τ HH e [ M y r ] Fig. 9: Planet growth with Semenov + Ferguson opacities (seeSects.3.2.1,4.4). M thresh = ⊕ , β = . M Z = − M ⊕ / yr.˙ M Z [M ⊕ / yr] 10 − − × − M cross [M ⊕ ] 2.83 3.7 4.45 t cross [Myr] 3.68 0.54 0.14 t run [Myr] 3.48 0.56 0.16Table 6: Crossover mass, crossover time, and time of the onset ofthe runaway of gas for di ff erent accretion rates of solids. M thresh = ⊕ . For completeness, we test in this section the e ff ect of consideringother values of accretion rates of solids. We run simulations withenvelope enrichment and the microphysics of our nominal model(low opacities) but using ˙ M Z = − M ⊕ / yr and 5 × − M ⊕ / yrinstead of ˙ M Z × − M ⊕ / yr. We summarise results in Tables 6and 7.Increasing the accretion rate of solids increases the core lu-minosity. Thus, the larger the accretion rate of solids, the hot-ter the envelope, which prevents gas accretion. Therefore, by in-creasing ˙ M Z , we expect to obtain planets with larger planetarymasses before the onset of runaway of gas. Indeed, we note thatwith an accretion rate as high as 5 × − M ⊕ / yr (more typicalfor pebble accretion than for planetesimal accretion), when τ HHe reaches its maximum, the total mass of the planet is M P ≈ ⊕ , with a core of M core = ⊕ and total fraction of H-He of ∼ ∼ M Z [M ⊕ / yr] 10 − − × − M P [M ⊕ ] 5.0 7.8 11.3 M core [M ⊕ ] 2.75 3.77 5.1 M HHe [M ⊕ ] 1.5 2.24 3.26 τ HHe [yr] 7 . × . × . × Table 7: Total mass, core mass, mass of H-He and value of τ HHe at the onset of runaway of gas (t = t run ) for di ff erent accretionrates of solids. M thresh = ⊕ . In all the simulations presented in this work we have includedthe e ff ect of water condensation (as it should be). In Venturiniet al. (2015) we showed that water condensation decreases theadiabatic gradient, and this plays a relevant role in diminishingthe critical core masses to very low values, especially for veryhigh envelope metallicities ( Z env (cid:38) . ff ect of water condensation in the compu-tation of the adiabatic gradient. For doing this, we run simula-tions using CEA for the equation of state (Gordon et al. 1994)as in Venturini et al. (2015), because with this package it is pos-sible to switch o ff the e ff ect of water condensation in the com-putation of the adiabatic gradient. The results presented in thissection were, therefore, all performed with CEA (those wherewater condensation is included, and those were it is not).In the nominal case where we use Mordasini + Freedmanopacities, the di ff erence between including or not water conden-sation does not really a ff ect the evolution, because with theselow opacities the outer layers of the envelope are radiative, sothe structure in those layers is independent of the value of theadiabatic gradient.For higher opacities and / or accretion rates of solids, the en-velopes are more prone to be convective, so for those cases, wa-ter condensation plays an non-negligible role. For instance, con-sidering the nominal case but increasing the dust opacities bya factor of 300, and taking an accretion rate of solids of ˙ M Z = × − M ⊕ / yr, we find that the total formation time (when theplanet reaches a mass of ∼
40 M ⊕ ) is of 5 Myr for the case wherewater condensation is taken into account, and 6 Myr when it isnot. Concerning the maximum Z env attained, for the former caseit is 0.75, and for the latter, 0.85. This is consistent with whatwe found in Venturini et al. (2015): water condensation makesenvelopes thicker for the same core mass, which explains thefact of reaching a smaller envelope metallicity and accreting gasfaster than in the case where water condensation is neglected.
5. Results with an accretion rate of metals as in P96
In this section we study the e ff ect of envelope enrichment in theframework of a Pollack et al. (1996) accretion scheme. The maindi ff erence with fixing a given accretion rate of planetesimals isthat now the accretion rate of planetesimals depends on the totalmass of the protoplanet and that the availability of planetesimalsto be accreted depends on the initial amount at the neighbour-hood of the embryo (the initial surface density of solids).We implement first an accretion rate of solids a la Pollack with initial surface solids of Σ = / cm and nominal opacities(Mordasini + Freedman). It is important to note that the value of Σ is arbitrary and has been chosen so as to obtain a formationtimescale for the non-enriched case of 8 Myr, as in Pollack et al.(1996). However, in our calculations this value has no influence Article number, page 9 of 14 & A proofs: manuscript no. venturini2016_arxiv on the opacities even though it is likely that dust-to-gas ratioand total disk mass would change opacities. Recent results byMordasini (2014) and Ormel (2014) tend to show that a leastregarding the dust, an increase in the accretion rate of solids (or Σ in the context of the P96 accretion scheme) does not increasethe dust opacity, because the larger the flux of dust received bythe protoplanet, the faster the coagulation and settling within theenvelope.The results for the mentioned opacities and surface densityof solids are illustrated in the upper panel of Fig. 10. The en-riched case assumes M thresh = ⊕ and β changing at this coremass from 0 to 0.5. The choice of M thresh is to ensure that theenvelope mass is the appropriate to dissolve icy planetesimals of100 km radius (which is the choice of planetesimal size for theP96 scheme, see Tab.2). In this case, when phase 1 ends, the coremass is 2.9 M ⊕ , which means that M thresh is reached during phase2 of P96. We note that the decrease in formation time, comparedto the non-enriched case, is of a factor ∼ ff ect of using the combination of highopacities (Semenov + Ferguson) in this P96 accretion scheme(Fig. 10, bottom). To get the same overall formation time as be-fore, for the non-enriched case, we choose an initial surface den-sity of solids of
Σ =
10 g / cm (as in the baseline model of P96).Because of the high opacities and the high accretion rate of solidsachieved in phase 1, the envelope required to destroy icy plan-etesimals of 100 km size corresponds now to M thresh =
11 M ⊕ .We note that in this case, the formation timescale is reduced by afactor ∼ ff erencewith the low-opacity case is that now M thresh is achieved just be-fore phase 1 ends. This means that when envelope enrichmenttakes place during phase 1, there are feedback mechanisms act-ing that favour more rapid gas accretion. In order to understandthese mechanisms, we proceed to analyse the behaviour of thesolid accretion rate for the high-opacity scenario.The upper panel of Fig. 11 shows the change in the accretionrate of solids as a function of time for the enriched and non-enriched cases of the high-opacitiy simulation (Fig. 10, bottom).A bump in the accretion rate of solids occurs once the envelopebegins to be enriched (at M core =
11 M ⊕ ). If we observe the evo-lution of the capture radius (Fig. 11, bottom) we see that alsoat this moment, it grows considerably. While the outer bound-ary conditions are the same as in the non-enriched case, themean molecular weight of the envelope increases as accretionproceeds. Therefore, the density increases, so the planetesimalsare more e ffi ciently slowed down when crossing the atmosphere.This increases the capture radius which translates into a largeraccretion rate.Note that, in principle, the same e ff ect occurs if enrichmentsets in during phase 2. However, in that case, the accretion rateof solids is much smaller, and the availability of the protoplanetto increase the accretion rate of planetesimals is limited to theamount of planetesimals that can enter in the feeding zone ineach timestep. Hence, it makes sense that the increase of the ac-cretion rate of solids is smaller than when enrichment starts inphase 1.
6. Implications on the formation of different typesof planets
We have shown that by including the e ff ect of envelope en-richment during the growth of a planet, the timescale to forma gas giant is shorter than in the standard case of no envelope m a ss [ M ⊕ ] time [Myr]M core M HHe M env m a ss [ M ⊕ ] time [Myr]M core M HHe M env Fig. 10: Growth of the planet using a P96 accretion scheme, withthe assumption of accretion of icy planetesimals of 100 km size.For the enriched cases, M thresh was chosen, for self-consistency,to correspond to the breakup of icy planetesimals of 100 km size.The solid lines correspond to the enriched case, and the dashedlines to the non-enriched case (classical case where the envelopeis composed of H-He). Upper panel : nominal opacities (Mor-dasini + Freedman). The initial surface density of solids is
Σ = / cm . Lower panel : high opacities (Semenov + Ferguson). Theinitial surface density of solids is
Σ =
10 g / cm . Note that in thehigh-opacity case, when envelope enrichment is taken into ac-count, the overall formation time is reduced by a factor of ∼ ff erence on the initial surface density values between thetwo figures was chosen to get the same overall time formationfor the non-enriched cases.enrichment. In the classical context of planetesimal accretion,the formation of gas giants is challenging when the oligarchicgrowth regime is taken into account (Fortier et al. 2007, 2013)and even more challenging if planetesimal fragmentation is in-cluded (Guilera et al. 2014).Fortier et al. (2013) showed (for the classical scenario of H-He envelopes) that the formation of giant planets required smallplanetesimals ( ∼
100 m). They found that km-size planetesimals(and larger) were not su ffi ciently damped by gas drag, and there-fore, due to their large relative velocities, were not e ffi cientlyaccreted by the protoplanet. We have shown that envelope en-richment reduces formation timescales. Therefore, including thise ff ect in simulations that consider oligarchic growth might help Article number, page 10 of 14enturini et al.: Planet formation with envelope enrichment -7-6.5-6-5.5-5-4.5-4-3.5-3 0 5 10 15 20 l og ( d M Z / d t ) [ M ⊕ / y r ] M core [M ⊕ ]non-enriched caseenriched case -6.5-6-5.5-5-4.5-4 10.5 11 11.5 12 12.5 13 R c ap / R c o r e M core [M ⊕ ]non-enriched caseenriched-case Fig. 11:
Upper panel : ˙ M Z as a function of core mass for the en-riched and non-enriched cases of the high-opacities simulations(Fig. 10, bottom). Note that the planetesimal accretion rate ofthe enriched case increases when the core reaches the thresholdvalue of 11 M ⊕ . Lower panel : ratio of capture radius to core ra-dius as a function of core mass (same simulation as in the upperpanel).the formation of gas giants from km-size planetesimals. The rea-son for this is not related to envelope enrichment playing anyrole in the damping of planetesimals eccentricity, but to the factthat it accelerates the growth of the planetary envelope, produc-ing a positive feedback in the increase of the capture radius. Weshowed in Sect. 5 that envelope enrichment can reduce by a fac-tor of ∼ Another important result we obtained, is that if the disk disap-pears before the onset of runaway of gas, the type of planetswe can form with envelope enrichment is quite di ff erent thanin the non-enriched counterpart. With envelope enrichment weexpect failed giants to have large gas-to-core ratios ( ∼ ∼ (cid:38)
40 %). This result is very interesting, be-cause one ongoing problem with the classical picture of the coreaccretion model (without envelope enrichment) is the formationof low ( ∼ ⊕ ) and intermediate mass ( ∼
10 - 20 M ⊕ ) ob-jects with important contents of H-He and high total metallicity.This has been, for instance, a recurrent problem in the formationof Uranus and Neptune in the solar system, because even if thedisk disappears at the moment when the mass of the planet isin the appropriate range (15 -20 M ⊕ ), a core-dominated planetwith ∼
20% in mass of H-He is di ffi cult to obtain (Helled & Bo-denheimer 2014). The same regarding mini-Neptunes: the maxi-mum gas-to-core ratios expected from simulations for H-He en-velopes is typically (cid:46)
10% (Lee & Chiang 2016).Our nominal model is able to reproduce mini-Neptunes, butnot Neptunes. This is because of the low opacities. Gas opaci-ties from Freedman et al. (2014) include the e ff ect of a changingmetallicity (and therefore, increase with increasing Z env ), but thedust opacities of Mordasini (2014) are still very low. Therefore,with these opacities, crossover masses are low (even without en-velope enrichment, see Table 1). It is not possible with these lowopacities to obtain the static structure of a ∼
15 M ⊕ planet witha core of ∼
10 (no matter which value of envelope metallicity).Nevertheless, a Neptune-type planet can be formed with enve-lope enrichment if larger opacities are invoked, as we showed inSect. 4.4. This shows the importance of keep improving opacitymodels. The formation of clouds could be one mechanism forlarger opacities to exist in the envelope of protoplanets, and wehave shown that water clouds can be present in the atmospheresof protoplanets formed beyond the iceline (Venturini et al. 2015).So e ff ort on the direction of including this physics on the com-putation of opacities should continue. In Sect. 4.5 we showed that when considering large accretionrates, as the ones invoked in the context of pebble accretion, theformation of ice giants was di ffi cult because formation timeswere too short. In principle, the situation could be improvedconsidering larger opacities (as we showed in Sect. 4.4 forplanetesimal-like accretion rates). However, we note that whenwe run the non-enriched case with the high opacities of Semenov + Ferguson and a solid accretion rate of 5 × − M ⊕ / yr, run-away of gas is triggered at t run ≈ .
45 Myr. This tells us thatif we considered envelope enrichment with these high opacities,formation timescales of giant planets would still be too short (ofthe order of ∼ M iso (cid:38)
20 M ⊕ for distances fromthe star a (cid:38) Article number, page 11 of 14 & A proofs: manuscript no. venturini2016_arxiv planet is to the star, so they claim that Jupiter and Saturn reachedisolation mass, which caused the onset of runaway of gas (halt-ing solid accretion promotes gas accretion, Ikoma et al. 2000);whereas Uranus and Neptune did not.Indeed, Lambrechts et al. (2014) show that they can get thecorrect structure of the four giant planets of the solar system(in terms of total mass and total metallicity). However, they donot show any time evolution. Despite of the increase of criticalcore mass, if this mass is reached in 0.1 Myr, fast gas accretionfrom this stage on will be inevitable. Perhaps a way to avoid run-away of gas is to consider that planets accrete both pebbles andplanetesimals: when pebble isolation mass is reached, the planetcould remain accreting planetesimals that provide the necessaryluminosity to delay runaway of gas (just as in phase 2 of Pollacket al. 1996). This will be the subject of a future work.
7. Predictions and comparison with observations
We considered in this work enrichment by icy planetesimals.Therefore our predictions are relevant for planets formed be-yond the iceline. Planetesimal disruption inside the iceline couldtake place as well, but given that the composition would be moresilicate-iron rich, the thermal ablation of planetesimals would re-quire either very small planetesimals’ sizes, or thick envelopes.The last option would mean that enrichment would take placeprobably when the planet has already entered in the phase ofrunaway of gas. In addition, even if thermal ablation and me-chanical disruption occur, the fate of a mixture of silicates withhydrogen is just expected to be miscible for T (cid:38) ff ect of envelope enrichment would not op-erate during formation, meaning that critical core masses wouldhave the usual values of ∼
10 - 20 M ⊕ . This would imply thatforming giant planets from rocky planetesimals (i.e, inside theiceline) would be di ffi cult.Regarding formation beyond the iceline, our results showthat we can expect both gas giants and neptune-like planets toform via envelope enrichment, as explained above.We have shown that by including envelope enrichment by icyplanetesimals during formation, we expect the planets formedto present a decrease of envelope (and total) metallicity withincreasing total mass. . To illustrate this better, we show, inFig. 12, the output of planet formation for di ff erent choicesof opacities and accretion rates of planetesimals. In all cases,the envelope starts to get enriched when it reaches a mass of M env = × − M ⊕ , which corresponds to the fully disruptionof icy planetesimals of approximately, 1-10 km-size (Sect. 3.5).We have overplotted the total metallicities of the giant planetsof the solar system for reference. This figure suggest that ratherhigh opacities are preferable to explain the giants of the solarsystem. Of course, the behaviour of the curves depends some-how on the choice of parameters. Still, when we analysed, forinstance, the nominal model (Mordasini + Freedman opacities) The decrease of envelope metallicity with planetary mass waspointed out empirically by Kreidberg et al. (2014). The decrease of totalmetallicity with planetary mass is a natural outcome of core accretion(e.g., Alibert et al. 2005; Mordasini et al. 2009) for di ff erent M thresh , we saw that even the extreme case of non-enrichment could not lead to total high-Z content of more than ∼ ⊕ (Sect. 4). So it is clear that with those low opacities we can-not explain Jupiter, which has at least 15 M ⊕ of heavy elements(Bara ff e et al. 2014). This situation can be of course solved, atleast for Jupiter, by increasing the accretion rate of solids (Sect.4.5). However, for intermediate mass planets, an increase of thesolid accretion rate would make formation times so short, thatthey would likely become gas giants. So the combination of lowopacities and high accretion rates of solids could work for form-ing Jupiters, but forming Neptunes would be more challenging.We note that the dust opacities from Mordasini (2014) are an-alytical, and thus, many simplifications were considered to com-pute them. Perhaps the most important one is the assumption ofa predominant grain size. Mordasini (2014) mention that if plan-etesimal ablation is important in all layers of the envelope, thenthe constant supply of small grains could raise the dust opacities.It is likely that in our enriched scenario this plays a relevant role. Our formation model assumes that envelopes get enriched in wa-ter during formation, because icy planetesimals are expected tobe water rich. Water has indeed been detected in all the atmo-spheres of the giant planets of the solar system, although in verysmall amounts, even less than what is expected based on the de-tection of other volatile molecules, such as CH (Niemann et al.1998; Irwin et al. 2014). The problem with measuring water onthe outermost layers of the giant planets is that due to the lowtemperatures, this species is expected to have condensed deeperin, and therefore be present in the form of clouds.In this sense, transiting exoplanets o ff er a better opportunityto detect water in their atmospheres: since the equilibrium tem-peratures are much higher (these planets are much closer to theircentral star than the giants of the solar system, due to obser-vational bias), water is expected to be present in the form ofvapour. Indeed, of the 19 transiting planets whose spectra hasbeen measured, 10 of them present signatures of water vapour intheir atmospheres (the remaining 9 are thought to posses wateras well, but their signature to be obscured by clouds, Iyer et al.2015). From these 10 exoplanets, 9 have masses in the range of0.5 - 2 jovian masses, so of course, these planets are expectedto be hydrogen-helium rich. In other words, the amount of waterpresent is expected to be low, as our results suggest (envelopemetallicity should decrease with increasing planet mass).Two works report precise water abundances of hot Jupiter’satmospheres: Kreidberg et al. (2014) for WASP-43b and Mad-husudhan et al. (2014) for HD 189733b, HD 209458b andWASP-12b. We have converted the abundances they give of wa-ter mixing ratios into mass fraction of water ( Z env ) and plottedthe predictions of our models together with these estimations inFig. 13. We note that in the case of WASP-43 b, the predictionsof both high and low opacity models work surprisingly well toexplain the estimated abundance of water vapour. For the hotJupiters reported by Madhusudhan et al. (2014), we find that allour models predict a larger mass fraction of atmospheric water,at least by a factor of 2. Madhusudhan et al. (2014) claim aswell that the water abundances they find seem too low, and sug-gest that the presence of clouds might be obscuring some of themolecular features of the spectra (hypothesis supported as wellby Sing et al. 2016; Benneke 2015). Madhusudhan et al. (2014)suggest, alternatively, that the atmospheres of these hot-Jupiters See explanation for this conversion in Appendix A.Article number, page 12 of 14enturini et al.: Planet formation with envelope enrichment Z M planet [M ⊕ ] solid: Z tot dashed: Z env nominal modelopac = 400 Moropac = S+F JSNU
Fig. 12: Total metallicity (solid lines) and envelope metallicity(dashed lines) as a function of the planet mass for di ff erent mod-els of planet formation with envelope enrichment. The di ff erentcolours indicate di ff erent choices of opacities and accretion ratesof solids. The violet corresponds to our nominal model (opac-ity of Mordasini + Freedman, and ˙ M Z = − M ⊕ / yr). The or-ange uses opacities of Semenov + Ferguson, and ˙ M Z = × − M ⊕ / yr. The green curves were computed using the nominal opac-ities, but the dust opacities of Mordasini are enhanced by a fac-tor of 400. The accretion rate of solids is ˙ M Z = × − M ⊕ / yr. M thresh was chosen such as the corresponding envelope is mas-sive enough to destroy icy planetesimals of 1-10 km. The thicklines show the output of formation, and they all end at times of ∼ ff e et al. (2014) for Jupiter and Saturn). The circu-lar points indicate the onset of the runaway of gas. This meansthat masses at the right of these points, especially after the sim-ulations stop, are unlikely (the expected elapsed time betweenthe end of the simulations and the one required for M P ≈ ⊕ is ∼ yr ). This figure suggest that for the solar system,a combination of high opacities is preferred with respect to thelow opacities adopted in the nominal model.could bear more carbon-rich species than oxygen-rich. It is im-portant to remark that in our formation model, we assumed thatall the volatiles deposited in the planet’s envelope were just water(because of the equation of state). In this sense, it could perfectlybe that we overestimate the amount of water, since our model ne-glects the possibility of initial amounts of, e.g, CH , CO , CO.Another interesting case to link our results with observationsis the detection of water in HAT-P-11b (Fraine et al. 2014), theonly Neptune-mass exoplanet where the presence of water vaporhas been inferred from transmission spectroscopy. Fraine et al.(2014) report from retrieval models that the atmosphere of HAT-P-11b is expected to have a metallicity of 1 to 700 times solar,with no real preference for a specific value within this range(Benneke, priv. comm.). Unfortunately, this wide range couldimply an envelope metallicity from solar up to Z env ∼ .
9. HAT-P-11b has a total mass of M P ≈
26 M ⊕ , so our models would Z en v M planet [M ⊕ ] nominal modelopac = 400 Moropac = S+F WASP-43b HD HD WASP-12b
Fig. 13: Same as Fig. 12, but we show just the mass range corre-sponding to four hot Jupiters where estimations of atmosphericmass fraction of water ( Z env ) have been reported, and we overplotthese estimations.predict a Z env in the range of, approximately, 0.05 - 0.4 (see Fig.12); which would correspond to values of “metallicity times so-lar" between 1 and 100. More precise measurements would beneeded in order to determine if our models are able to predictthe metal atmospheric content of this planet or not.A last word regarding the detection of water in exoplanets:it is important to remark that even if the sample of transmissionspectra measured is still small, water shows to be common in theatmospheres of planets. This reinforces a formation scenario be-yond the iceline, with volatile-rich planetesimals being dissolvedin the atmospheres during formation. Even if in the future it isfound that statistically the atmospheres are water poor, this couldbe a hint of an initial planetesimal composition that is more car-bon rich (Madhusudhan et al. 2014). This is why it is importantto determine precisely other volatile abundances besides water.
8. Conclusions
We have performed the first self-consistent calculation of thegrowth of a planet including the e ff ect of envelope enrichmentdue to the dissolution of icy planetesimals / pebbles. We have im-plemented suited equations of state and opacities taking into ac-count di ff erent metallicities. Moreover, we have considered twodi ff erent sets of opacities in order to test the impact on our re-sults, which can be summarised as: • Envelope enrichment accelerates notably the formation ofgas giants. This is mainly a consequence of the increase ofthe mean molecular weight of the envelope. The thinner theenvelope (i.e, the smaller the planetesimal), the sooner enve-lope enrichment sets in and the shorter the timescale to forma giant planet. • When envelope enrichment is taken into account, low andintermediate mass planets (namely mini-Neptunes to Nep-tunes) can be formed, with total mass fractions of H-He upto 30%, this number being independent of the choice of opac-ities. • Low-opacities allow for the formation of mini-Neptunes,whereas high-opacities lead to the formation of Neptunes.
Article number, page 13 of 14 & A proofs: manuscript no. venturini2016_arxiv • High-opacities are preferable for explaining the total massand metallicity of the giant planets of the solar system. • We were able to quantify the amount of volatile materialremaining in the primordial atmospheres as a result of for-mation. These allowed us to compare our results with waterabundances inferred from the transmission spectra of tran-siting exoplanets. We find good agreement for WASP-43b,whereas for the other three cases, the measured water abun-dance is lower than what is predicted by our models. • We predict that envelope metallicity and total metallicityshould decrease with planetary mass.
Acknowledgments . We thank M. Ikoma for fruitful discussionsand F. Benitez for providing the EOS tables and revising thismanuscript. J.V. acknowledges financial support from the Swiss-based MERAC Foundation. This work has, in part, been carriedout within the framework of the National Centre for Competencein Research PlanetS supported by the Swiss National ScienceFoundation. The authors acknowledge financial support from theSNSF.
Appendix A: Conversion of water mixing ratios towater mass fractions
It is a common practice in papers that estimate abundance ofdi ff erent compounds to give them in terms of “volume mixingratios". In particular, for estimates on water abundances derivedfrom transmission spectroscopy measurements, the volume mix-ing ratio of water ( m r ) is defined as (Madhusudhan et al. 2014): m r = n (H O) n (H ) , (A.1)where n (H O) and n (H ) are the number density of watermolecules and hydrogen, respectively.Because of the molar masses of H O and H , in terms ofmass, the ratio of water to hydrogen reads: m H O m H = n ( H O ) n ( H ) = m r . (A.2)In the formation models presented here, we defined the en-velope metallicity as the mass fraction of water in the envelope.Thus, Z env = m H O m H + m He + m H O = m H O / m H + Y / X + m H O / m H = m r + ( Y / X ) (cid:12) + m r , (A.3)where we assumed the ratio of helium to hydrogen to be solar.This last equation is the one we use to estimate the values of Z env shown in Fig. 13 for the di ff erent hot-Jupiters whose watermixing ratios has been determined (Madhusudhan et al. 2014;Kreidberg et al. 2014). References
Alexander, D. R. & Ferguson, J. W. 1994, ApJ, 437, 879Alibert, Y., Mordasini, C., Benz, W., & Winisdoer ff er, C. 2005, A&A, 434, 343Bara ff e, I., Chabrier, G., & Barman, T. 2008, A&A, 482, 315 Bara ff e, I., Chabrier, G., Fortney, J., & Sotin, C. 2014, Protostars and Planets VI,763Benitez, F., Mordasini, C., & Venturini, J. 2016, in preparationBenneke, B. 2015, ArXiv e-printsBitsch, B., Lambrechts, M., & Johansen, A. 2015, A&A, 582, A112Bodenheimer, P. & Pollack, J. B. 1986, Icarus, 67, 391Fortier, A., Alibert, Y., Carron, F., Benz, W., & Dittkrist, K.-M. 2013, A&A, 549,A44Fortier, A., Benvenuto, O. G., & Brunini, A. 2007, A&A, 473, 311Fraine, J., Deming, D., Benneke, B., et al. 2014, Nature, 513, 526Freedman, R. S., Lustig-Yaeger, J., Fortney, J. J., et al. 2014, ApJS, 214, 25Gordon, S., McBride, B. J., & States., U. 1994, Computer program for calcu-lation of complex chemical equilibrium compositions and applications [mi-croform] / Sanford Gordon, Bonnie J. McBride (National Aeronautics andSpace Administration, O ffi ce of Management, Scientific and Technical Infor-mation Program ; National Technical Information Service, distributor [Cleve-land, Ohio] : [Springfield, Va), 2 v. :Guilera, O. M., de Elía, G. C., Brunini, A., & Santamaría, P. J. 2014, A&A, 565,A96Guilera, O. M., Fortier, A., Brunini, A., & Benvenuto, O. G. 2011, A&A, 532,A142Guillot, T. & Gautier, D. 2014, ArXiv e-printsGuillot, T., Stevenson, D. J., Hubbard, W. B., & Saumon, D. 2004, The interiorof Jupiter, ed. F. Bagenal, T. E. Dowling, & W. B. McKinnon, 35–57Helled, R., Anderson, J. D., Podolak, M., & Schubert, G. 2011, ApJ, 726, 15Helled, R. & Bodenheimer, P. 2014, ApJ, 789, 69Hori, Y. & Ikoma, M. 2011, MNRAS, 416, 1419Hubickyj, O., Bodenheimer, P., & Lissauer, J. J. 2005, Icarus, 179, 415Ikoma, M., Nakazawa, K., & Emori, H. 2000, ApJ, 537, 1013Inaba, S. & Ikoma, M. 2003, A&A, 410, 711Irwin, P. G. J., Lellouch, E., de Bergh, C., et al. 2014, Icarus, 227, 37Iyer, A. R., Swain, M. R., Zellem, R. T., et al. 2015, ArXiv e-printsKofman, W., Herique, A., Barbin, Y., et al. 2015, Science, 349Kreidberg, L., Bean, J. L., Désert, J.-M., et al. 2014, ApJ, 793, L27Lambrechts, M. & Johansen, A. 2012, A&A, 544, A32Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35Lee, E. J. & Chiang, E. 2016, ApJ, 817, 90Lee, E. J., Chiang, E., & Ormel, C. W. 2014, ApJ, 797, 95Lissauer, J. J., Hubickyj, O., D’Angelo, G., & Bodenheimer, P. 2009, Icarus,199, 338Lopez, E. D. & Fortney, J. J. 2014, ApJ, 792, 1Madhusudhan, N., Crouzet, N., McCullough, P. R., Deming, D., & Hedges, C.2014, ApJ, 791, L9Melosh, H. J. 2007, Meteoritics and Planetary Science, 42, 2079Mizuno, H. 1980, Progress of Theoretical Physics, 64, 544Mordasini, C. 2014, A&A, 572, A118Mordasini, C., Alibert, Y., & Benz, W. 2006, in Tenth Anniversary of 51 Peg-b:Status of and prospects for hot Jupiter studies, ed. L. Arnold, F. Bouchy, &C. Moutou, 84–86Mordasini, C., Alibert, Y., & Benz, W. 2009, A&A, 501, 1139Mordasini, C., Alibert, Y., Klahr, H., & Henning, T. 2012, A&A, 547, A111Niemann, H. B., Atreya, S. K., Carignan, G. R., et al. 1998, J. Geophys. Res.,103, 22831Ormel, C. W. 2014, ApJ, 789, L18Papaloizou, J. C. B. & Terquem, C. 1999, ApJ, 521, 823Perri, F. & Cameron, A. G. W. 1974, Icarus, 22, 416Piso, A.-M. A. & Youdin, A. N. 2014, ApJ, 786, 21Podolak, M., Pollack, J. B., & Reynolds, R. T. 1988, Icarus, 73, 163Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62Saumon, D., Chabrier, G., & van Horn, H. M. 1995, ApJS, 99, 713Semenov, D., Henning, T., Helling, C., Ilgner, M., & Sedlmayr, E. 2003, A&A,410, 611Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59Soubiran, F. & Militzer, B. 2015, ApJ, 806, 228Stevenson, D. J. 1982, Planet. Space Sci., 30, 755Thiabaud, A., Marboeuf, U., Alibert, Y., Leya, I., & Mezger, K. 2015, A&A,574, A138Thompson, S. 1990, Sandia National Laboratories.Vazan, A., Helled, R., Kovetz, A., & Podolak, M. 2015, ApJ, 803, 32Venturini, J., Alibert, Y., Benz, W., & Ikoma, M. 2015, A&A, 576, A114Wilson, H. F. & Militzer, B. 2012, Physical Review Letters, 108, 111101Wuchterl, G. 1993, Icarus, 106, 323ce of Management, Scientific and Technical Infor-mation Program ; National Technical Information Service, distributor [Cleve-land, Ohio] : [Springfield, Va), 2 v. :Guilera, O. M., de Elía, G. C., Brunini, A., & Santamaría, P. J. 2014, A&A, 565,A96Guilera, O. M., Fortier, A., Brunini, A., & Benvenuto, O. G. 2011, A&A, 532,A142Guillot, T. & Gautier, D. 2014, ArXiv e-printsGuillot, T., Stevenson, D. J., Hubbard, W. B., & Saumon, D. 2004, The interiorof Jupiter, ed. F. Bagenal, T. E. Dowling, & W. B. McKinnon, 35–57Helled, R., Anderson, J. D., Podolak, M., & Schubert, G. 2011, ApJ, 726, 15Helled, R. & Bodenheimer, P. 2014, ApJ, 789, 69Hori, Y. & Ikoma, M. 2011, MNRAS, 416, 1419Hubickyj, O., Bodenheimer, P., & Lissauer, J. J. 2005, Icarus, 179, 415Ikoma, M., Nakazawa, K., & Emori, H. 2000, ApJ, 537, 1013Inaba, S. & Ikoma, M. 2003, A&A, 410, 711Irwin, P. G. J., Lellouch, E., de Bergh, C., et al. 2014, Icarus, 227, 37Iyer, A. R., Swain, M. R., Zellem, R. T., et al. 2015, ArXiv e-printsKofman, W., Herique, A., Barbin, Y., et al. 2015, Science, 349Kreidberg, L., Bean, J. L., Désert, J.-M., et al. 2014, ApJ, 793, L27Lambrechts, M. & Johansen, A. 2012, A&A, 544, A32Lambrechts, M., Johansen, A., & Morbidelli, A. 2014, A&A, 572, A35Lee, E. J. & Chiang, E. 2016, ApJ, 817, 90Lee, E. J., Chiang, E., & Ormel, C. W. 2014, ApJ, 797, 95Lissauer, J. J., Hubickyj, O., D’Angelo, G., & Bodenheimer, P. 2009, Icarus,199, 338Lopez, E. D. & Fortney, J. J. 2014, ApJ, 792, 1Madhusudhan, N., Crouzet, N., McCullough, P. R., Deming, D., & Hedges, C.2014, ApJ, 791, L9Melosh, H. J. 2007, Meteoritics and Planetary Science, 42, 2079Mizuno, H. 1980, Progress of Theoretical Physics, 64, 544Mordasini, C. 2014, A&A, 572, A118Mordasini, C., Alibert, Y., & Benz, W. 2006, in Tenth Anniversary of 51 Peg-b:Status of and prospects for hot Jupiter studies, ed. L. Arnold, F. Bouchy, &C. Moutou, 84–86Mordasini, C., Alibert, Y., & Benz, W. 2009, A&A, 501, 1139Mordasini, C., Alibert, Y., Klahr, H., & Henning, T. 2012, A&A, 547, A111Niemann, H. B., Atreya, S. K., Carignan, G. R., et al. 1998, J. Geophys. Res.,103, 22831Ormel, C. W. 2014, ApJ, 789, L18Papaloizou, J. C. B. & Terquem, C. 1999, ApJ, 521, 823Perri, F. & Cameron, A. G. W. 1974, Icarus, 22, 416Piso, A.-M. A. & Youdin, A. N. 2014, ApJ, 786, 21Podolak, M., Pollack, J. B., & Reynolds, R. T. 1988, Icarus, 73, 163Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62Saumon, D., Chabrier, G., & van Horn, H. M. 1995, ApJS, 99, 713Semenov, D., Henning, T., Helling, C., Ilgner, M., & Sedlmayr, E. 2003, A&A,410, 611Sing, D. K., Fortney, J. J., Nikolov, N., et al. 2016, Nature, 529, 59Soubiran, F. & Militzer, B. 2015, ApJ, 806, 228Stevenson, D. J. 1982, Planet. Space Sci., 30, 755Thiabaud, A., Marboeuf, U., Alibert, Y., Leya, I., & Mezger, K. 2015, A&A,574, A138Thompson, S. 1990, Sandia National Laboratories.Vazan, A., Helled, R., Kovetz, A., & Podolak, M. 2015, ApJ, 803, 32Venturini, J., Alibert, Y., Benz, W., & Ikoma, M. 2015, A&A, 576, A114Wilson, H. F. & Militzer, B. 2012, Physical Review Letters, 108, 111101Wuchterl, G. 1993, Icarus, 106, 323