EFFECTS OF DYNAMICAL EVOLUTION OF GIANT PLANETS ON THE DELIVERY OF ATMOPHILE ELEMENTS DURING TERRESTRIAL PLANET FORMATION

Soko Matsumura; Ramon Brasser; Shigeru Ida

doi:10.3847/0004-637X/818/1/15

1. INTRODUCTION

Recent observations of extrasolar planetary systems have shown that even similar-sized planets have a wide range of densities, which indicates a variety of planetary compositions ranging from a water world to an iron core (e.g., Howe et al. 2014; Lopez & Fortney 2014). The observations of planetary atmospheres revealed the importance of the C/O ratio in estimating the compositions and formation processes of planets (e.g., Madhusudhan et al. 2011). Furthermore, the advance in observations of protostellar disks started revealing the distribution of refractory and volatile molecules across the disk, as well as the locations of the condensation fronts of some volatiles (e.g., Mathews et al. 2013; Qi et al. 2013). Future observations will place further constraints on the compositions of protostellar disks and planets, and help us understand the solar system as one of the planetary systems.

The compositions of planets may be largely determined by the primordial compositions of a protostellar disk. However, there are also good indications that the dynamical evolution of giant planets could significantly affect the compositions of terrestrial planets. Two major pathways for giant planets to influence the formation and survival of rocky planets are migration in the disk and dynamical instabilities among them. The migration of a giant planet through the inner disk could potentially be hazardous to the formation of terrestrial planets (Armitage 2003). Some numerical studies showed that small bodies interior to a migrating giant planet could survive either by being shepherded toward the star or by being scattered outward (e.g., Fogg & Nelson 2005; Zhou et al. 2005), although the survival rate of such bodies depends on the migration rate of the massive body (e.g., Ward & Hahn 1995; Tanaka & Ida 1999; Edgar & Artymowicz 2004). Raymond et al. (2006) showed that such migration could mix the disk materials radially and create a variety of planetary systems with both dry and water-rich planets.

The dynamical instability of giant planets, on the other hand, could be a major threat during the final assembly of terrestrial planets, which is expected to occur after the disk's dispersal. Indeed, the presence of dynamically active giant planets could drive building blocks of rocky planets to merge with the central star due to the sweeping of secular resonances, or to be ejected out of the system due to scattering by giant planets (e.g., Veras & Armitage 2005, 2006; Raymond et al. 2011, 2012; Matsumura et al. 2013). Raymond et al. (2012) showed that some outer disk materials could be scattered into the inner disk by giant planets, creating a water-rich planet. Thus, both the migration and dynamical instability of giant planets could help move the building blocks of terrestrial planets from their initial locations.

In this paper we explore how the dynamical evolution of giant planets affects the composition of terrestrial planets by investigating solar system formation. There are two major models for the formation of rocky planets in the solar system: a classical formation model and a so-called Grand Tack model. In the classical model, giant planets are assumed to have present masses and current orbital parameters. Chambers (2001) studied such a scenario and reproduced various features of the solar system planets, including the formation timescale (within 200 Myr), the number of planets (typically three to four), the fraction of mass in the largest planet compared to the others, and wide orbital spacings compared to giant planets. However, this model disagreed with the observed systems in some ways: eccentricities and inclinations were higher than the terrestrial planets in the solar system, and masses of Mercury and Mars analogues were higher than those planets. For the former point, the dynamical friction due to a large number of small planetesimals was proposed to be a potential solution to high eccentricities and inclinations of simulated rocky planets (e.g., Chambers 2001; O'Brien et al. 2006; Raymond et al. 2006). For the latter point, Hansen (2009) showed that the mass distribution of the solar system could be reproduced when the disk materials were concentrated in a very narrow region near Venus and Earth orbits. However, this work did not explain how such a high density region could be created.

Walsh et al. (2011) proposed a so-called Grand Tack model, which could at least partially explain the creation of a high density region in the inner disk. In the first step of their model, Jupiter is assumed to form before Saturn, clear the gas in an annulus whose width is comparable to its Hill radius, and undergo inward Type II migration (e.g., Ward 1997). As Jupiter migrates, it scatters some material out of the system and shepherds others down to the inner disk, creating an enhanced density region for terrestrial planet formation. Once Saturn grows to a critical mass, it is assumed to open a partial gap, migrate rapidly, catch up with Jupiter, and get trapped in a mean-motion resonance. In this process, these two giant planets share a gap in the disk together, which effectively means that Jupiter receives an angular momentum from the inner disk and Saturn loses an angular momentum to the outer disk while they exchange an angular momentum with each other. Since Jupiter is more massive than Saturn, the pair gain angular momentum overall and migrate outward (e.g., Masset 2001; Morbidelli & Crida 2007; Pierens & Nelson 2008; Pierens & Raymond 2011). In other words, the planets reverse their migration and they "tack." Once the giant planets move away and the gas disk dissipates, terrestrial planets could form via collisions between embryos and planetesimals, free from close encounters with giant planets. Their work successfully reproduced the mass distribution of rocky planets in the solar system including Mars analogues and also explained the compositional differences across the asteroid belt. However, it is unclear whether each step of this scenario would work correctly. For example, the particular configuration and outward migration favored by Grand Tack is only supported for a specific set of initial conditions (D'Angelo & Marzari 2012) and is not universal (Zhang & Zhou 2010). There may also be other pathways to produce such a high density region (e.g., Izidoro et al. 2014). Furthermore, the Grand Tack model needs to explain how the innermost disk (<0.5 AU) is cleaned up.

Here, we will numerically study both classical and Grand Tack formation models and explore the differences in the elemental abundances of planets. In Section 2, we describe chemical models as well as numerical methods and the initial conditions. In Section 3, we present results for different formation models and compare their results. In Section 4, we discuss and summarize our work.

Throughout this paper, we use the following terminologies for various elements. The standard classification scheme identifies the major rock-forming elements (e.g., Mg, Fe, Si, and Ni), and categorizes the other elements into four groups according to their volatility (e.g., McDonough 2003; Davis 2006). The refractory elements (e.g., Ca, Al, and Ti in Table 2) are the ones with equilibrium condensation temperatures higher than the major elements, and the critical temperature is >1335 K for gas of solar composition at a pressure of 10⁻⁴ atm. The moderately volatile elements (e.g., Cr, Li, Na, K, and P) are the ones with 1335–665 K for the same condition, and the highly volatile elements (e.g., Cd and Pb, not included in Table 2) are the ones with <665 K. The atmophile elements (e.g., H, C, N, and O in Table 2, as well as noble gases) condense at and below the condensation temperature of water (180 K at 10⁻⁴ atm).

2. NUMERICAL METHODS AND INITIAL CONDITIONS

2.1. Dynamical Model

All the simulations in this paper start with small planetesimals, nearly Mars-mass embryos as well as Jupiter- and Saturn-like planets. We do not include Uranus and Neptune in any of the simulations because they do not have as direct and significant effects as Jupiter and Saturn on terrestrial planet formation.

In all Grand Tack simulations, we do not take into account the effect of Saturn's mass growth because Walsh et al. (2011) showed that the overall results were largely determined by the tack location and the initial location of Jupiter, and that other effects—including the presence of Uranus and Neptune, the radial evolution of Saturn, and the migration timescale—did not change the results much. We assume that Jupiter has a current mass and is initially placed on a near-circular orbit at 3.5 AU, and that a fully grown Saturn is placed in the 2:3 resonance with Jupiter at 4.5 AU. Grand Tack models have been explored by assuming the tack location of 1.5 AU (e.g., Walsh et al. 2011; Jacobson et al. 2014). However, as noted by Walsh et al. (2011), the tack location does not need to be 1.5 AU, and this should be constrained by the observed properties of the solar system planets. Here, we test two sets of Grand Tack simulations with different tack locations: the first at 1.5 AU (hereafter GT1.5AU) and the second at 2.0 AU (hereafter GT2.0AU).

The migration of the gas giants was mimicked through fictitious forces as detailed in Walsh et al. (2011). In GT1.5AU, during the first 0.1 Myr, Jupiter and Saturn migrate from the initial locations (3.5 AU and 4.5 AU) to the tack locations (1.5 AU and 2.0 AU), respectively. For the next 0.4 Myr, Jupiter and Saturn migrate out to 5.2 and 6.9 AU, and end up in the appropriate initial conditions for the NICE model (Gomes et al. 2005; Morbidelli et al. 2005; Tsiganis et al. 2005). While the gas disk is around, the planetary embryos experience tidal damping of their eccentricities and inclinations from the gas disk, and the planetesimals suffer that from gas drag. After 0.5 Myr, the disk is assumed to be dispersed completely and bodies just experience gravitational interactions with one another. For gas drag purposes, we assume a planetesimal size is 50 km. The migration in GT2.0AU is set up in a similar manner to GT1.5AU.

For the classical formation model, we follow the approaches by Chambers (2001) and O'Brien et al. (2006) and assume that Jupiter and Saturn have current masses, radii, and orbital properties. We will call these simulations Eccentric Jupiter and Saturn (hereafter EJS) runs, following O'Brien et al. (2006).

We employ another set of simulations, where Jupiter and Saturn are on circular and nearly aligned orbits, and call them Circular Jupiter and Saturn (hereafter CJS) runs. Specifically, Jupiter and Saturn initially have semimajor axes of a_J = 5.45 AU and a_S = 8.18 AU, eccentricities of e_J = e_S = 0, and inclinations of i_J = 0 and i_S = 0 fdg 5. These conditions mimic the Nice Model that leads to the late heavy bombardment in the time comparable to the observational estimate (∼700 Myr; Gomes et al. 2005). In both EJS and CJS simulations, Jupiter and Saturn do not migrate in the disk.

To compare the effects of different formation models on the compositions of planets, we use a similar distribution of embryos and planetesimals for all cases. All simulations in this paper assume the equal-mass embryos, which agrees better with a pebble-accretion scenario (Morbidelli et al. 2015b), rather than a more traditional oligarchic-growth scenario (Kokubo & Ida 1998). We will discuss the differences between these two scenarios in a follow-up paper (R. Brasser et al. 2016, in preparation).

We also assume that the total mass ratio of embryos and planetesimals in the inner disk (0.3–3 AU) is 1:1. More specifically, the total inner disk mass is 3.74 M_⊕, with half the mass in 20 embryos with a mass comparable to Mars (9.34 × 10⁻² M_⊕) and the other half in 1571 planetesimals. Such a bimodal distribution was shown to produce a number of planets similar to the solar system (3.5 ± 0.5), compared to more uniform mass distributions (2.9 ± 0.6, Chambers 2001). We do not explore different total mass ratios for embryos and planetesimals in this work. However, Jacobson et al. (2014) studied such effects and reported that the larger mass ratios tend to lead to higher excitation and less concentration of planets compared to terrestrial planets in the solar system.

All simulations in this paper also include the outer planetesimals, unless noted otherwise. These consist of 500–650 planetesimals with an individual mass of 1.2 × 10⁻⁴ M_⊕, following Rubie et al. (2015). The outer planetesimals are distributed over 5–9 AU for GT, 3–8 AU for CJS, and 3–9.5 AU for EJS models. The eccentricities and inclinations of both embryos and planetesimals are randomly chosen from a uniform distribution between 0–0.01 and 0°–0 fdg 5 as in Chambers (2001) and O'Brien et al. (2006). The other angular orbital elements were randomly chosen from 0°–360°.

The system consisting of gas giants, planetary embryos, and planetesimals is integrated with the SyMBA software package (Duncan et al. 1998), with a time step of 0.02 year for 200 Myr. We compute the mutual gravity between gas giants and embryos, but the planetesimals are treated as test particles that are only perturbed by other bodies and do not affect others. The approximation is used to keep the CPU time within reasonable limits, and is justified because Jupiter destroys the disk beyond 1 AU in 0.1 Myr. Planets and planetesimals are removed once they are farther than 100 AU from the Sun (whether bound or unbound) or when they collide with a planet or venture closer than 0.1 AU from the Sun. We run eight simulations without and eight simulations with outer planetesimals for both GT1.5AU and GT2.0AU, whereas we run six simulations with outer planetesimals for both CJS and EJS.

2.2. Chemistry Model

Knowing the primordial compositions of a protostellar disk is a starting point for estimating the compositions of planets in the context of planet formation. The primordial compositions of the protostellar disk in the solar system can be inferred from the solar photosphere abundances because a protostellar disk is formed around a protostar as a molecular core gravitationally collapses, and thus should be formed from the same materials as the star. Such a collapse and the following disk accretion can heat the disk up intensely and reset the chemistry that occurred in the molecular cloud core. In this picture, the elements condense as the disk cools by following the condensation sequence. The assumption is supported by the elemental abundance pattern of meteorites and bulk terrestrial planets (e.g., Davis 2006; Pontoppidan et al. 2014). The abundances of these bodies are depleted with respect to the most primitive bodies, CI chondrites, which have abundances that are comparable to the solar photosphere except for the most volatile elements (e.g., Lodders 2003), and the abundances overall decrease with the 50% condensation temperatures of these elements (e.g., Cassen 1996; Davis 2006). The equilibrium condensation in an evolving protostellar disk has been shown to explain such a depletion pattern well (Cassen 1996, 2001). The assumption also predicts that refractory minerals such as Ca-, Al-, and Ti-oxides would be the first to condense, which agrees with the observation that Calcium–Aluminum Inclusions (CAIs) are the oldest objects in the solar system (e.g., Ciesla & Charnley 2006; Pontoppidan et al. 2014).

However, in particular in the outer region of a disk, the temperature is expected to never reach the sublimation temperatures of some of the most volatile species (e.g., Messenger et al. 2006). Thus, in such regions, presolar grains from a molecular cloud core could survive. The assumption is supported by the widespread presence of presolar grains (e.g., Davis 2006; Pontoppidan et al. 2014). The picture of the primordial compositions of a proto-solar disk is further complicated because the path to equilibrium condensation can be altered by various processes, such as transport of chemical species and transient heating events (Ciesla & Charnley 2006; Ebel 2006). Overall, most regions of the disk show evidence of both the inheritance from the cloud core and the disk processing after resetting the previous chemical features, whereas the latter plays a more important role in the inner disk (<20 AU, Pontoppidan et al. 2014). In this paper, we are interested in the disk compositions and thus the compositions of planet building blocks within 10 AU, as shown in Section 2.1. Therefore, we ignore the effects of presolar grains and other disk processing effects, such as shocks and lightning, and instead focus on the equilibrium condensation model.

Bond et al. (2010) followed the idea by Cassen (1996) and applied the chemical equilibrium model to planet formation simulations by O'Brien et al. (2006). For given temperature and pressure in a certain region of a disk, they determined which gaseous species condense out of the disk by using a commercial software HSC chemistry, and used them as compositions of planetesimals and embryos in that region. They found that the abundances of refractory and moderately volatile elements in planets agree well with the estimated values of terrestrial planets in the solar system. A notable disagreement between their model and the actual planets is a lack of accretion of the atmophile elements C and N. The abundances of these elements in the Earth are very small (see Table 6), but they are two of six major biogenic elements (H, C, N, O, S, and P) and thus are crucial for the origin of life. In Bond et al. (2010), only C and N did not accrete onto planets because of their low condensation temperatures; major species such as CO and NH₃ were not included as solid species in their calculations. Marboeuf et al. (2014b) developed a model that calculates the compositions of ices in a protostellar disk. We take into account the effects of icy species including C and N by following their work as described below.

Another assumption made in Bond et al. (2010) is that the composition of the building blocks of planets (e.g., planetesimals and planetary embryos) are determined at a single instant in the disk's evolution. Moriarty et al. (2014) recently improved this assumption and took account of the sequential condensation effect by converting a fraction of a disk into planetesimals, evolving the disk, and recalculating the equilibrium chemical composition at every time step. For the solar system C/O value (0.54), they found that the planetary compositions are similar for both equilibrium and sequential condensation assumptions. On the other hand, for C-rich systems with C/O > 0.8, they showed that the C-rich planetesimals and planets form over a wider orbital radii. Because we focus on the solar system in this study, we determine the initial compositions of building blocks at a single instant, as in Bond et al. (2010) and Elser et al. (2012). This also helps us clearly see a difference between the iceline evolution and the dynamical effects.

For refractory and moderately volatile elements, we follow previous studies and use the HSC Chemistry that determines the equilibrium composition by iteratively minimizing the system's Gibbs free energy (e.g., an isothermal-isobaric system in equilibrium—ΔT = 0 and ΔP = 0). The underlying assumption here is that the disk evolves slower than the time required to reach the chemical equilibrium, so that the temperature and the pressure are approximately constant. This assumption is likely valid, except for the icy midplane beyond the iceline, because the chemical reaction timescale is expected to be short (∼10⁴ year or shorter; Woitke et al. 2009) compared to the disk's evolution timescale.

As in the previous studies, we start with the solar photosphere abundances of 16 elements in Table 1, and calculate their equilibrium condensation at each radius of a disk model. This assumption is reasonable because the Sun and planets are formed out of the same molecular cloud core, and because the stellar pollution effect is expected to be limited. The initial compositions of planetesimals and embryos are determined for temperatures and pressures corresponding to a Chambers' disk model (Chambers 2009). We use 12 different initial disk ages— ${t}_{{\rm{init,disk}}}\;=\;0$ , 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1 Myr. Note that t_init,disk does not indicate the time evolution of a disk, and simply determines the initial temperature and pressure profiles of the disk. We assume that all planetesimals and embryos form at t_init,disk and determine their initial compositions at this instant. Table 2 summarizes the condensation temperatures for solid species that condense for t_init,disk = 0.2 Myr.

Table 1. Elemental Abundances used for HSC chemistry

Element Name	Abundance (mol)
H	1.00 × 10¹²
He	8.51 × 10¹⁰
C	2.45 × 10⁸
N	6.03 × 10⁷
O	4.57 × 10⁸
Na	1.48 × 10⁶
Mg	3.39 × 10⁷
Al	2.34 × 10⁶
Si	3.24 × 10⁷
P	2.29 × 10⁵
S	1.38 × 10⁷
Ca	2.04 × 10⁶
Ti	7.94 × 10⁴
Cr	4.37 × 10⁵
Fe	2.82 × 10⁷
Ni	1.70 × 10⁶

Note. Column 1: name of the element, and column 2: elemental abundance in moles. These values are scaled with respect to H to use in HSC chemistry. All values are taken from Bond et al. (2010).

Download table as: ASCII Typeset image

Table 2. Condensation Temperatures for Solid Species included in HSC calculations

Solid Name	Temperature (K)
${\mathrm{CaAl}}_{12}{{\rm{O}}}_{19}$	1790.5
$\mathrm{CaS}$	1633.5
CaTiO₃	1633.5
${\mathrm{Ca}}_{2}{\mathrm{Al}}_{2}{\mathrm{SiO}}_{7}$	1452.9
${\mathrm{CaMgSi}}_{2}{{\rm{O}}}_{6}$	1343.3
Fe	1343.3
MgSiO₃	1319.3
${\mathrm{Mg}}_{2}{\mathrm{SiO}}_{4}$	1319.3
Ni	1295.7
${\mathrm{Fe}}_{3}{\rm{P}}$	1191.1
${\mathrm{Cr}}_{2}{\mathrm{FeO}}_{4}$	1162.8
${\mathrm{MgAl}}_{2}{{\rm{O}}}_{4}$	1135.1
${\mathrm{Ti}}_{2}{{\rm{O}}}_{3}$	1135.1
${\mathrm{Ca}}_{3}{({\mathrm{PO}}_{4})}_{2}$	876.4
${\mathrm{NaAlSi}}_{3}{{\rm{O}}}_{8}$	850.5
$\mathrm{FeS}$	723.0
FeSiO₃	723.0
${\mathrm{Fe}}_{2}{\mathrm{SiO}}_{4}$	410.7
${\mathrm{Fe}}_{3}{{\rm{O}}}_{4}$	305.9
${\mathrm{Al}}_{2}{{\rm{O}}}_{3}$	227.8
${\mathrm{Mg}}_{3}{\mathrm{Si}}_{2}{{\rm{O}}}_{5}{(\mathrm{OH})}_{4}$	227.8
${{\rm{H}}}_{2}{\rm{O}}$	175.9

Note. Column 1: names of the solid species, and column 2: condensation temperatures (K) for t_init,disk = 0.2 Myr.

Download table as: ASCII Typeset image

For atmophile elements, we follow Marboeuf et al. (2014b) and take into account both condensation and clathrate hydrate formation on the surface of refractory grains (also see Marboeuf et al. 2014a; Thiabaud et al. 2014). A clathrate hydrate is a "cage" compound where a gas molecule such as ${{\rm{H}}}_{2}{\rm{S}}$ , CH₄, CO, or N₂ is trapped inside the water ice lattice. The clathrate formation could occur at higher temperatures than pure condensation, and thus may help incorporate ices into building blocks of planets. They assume solar abundances and consider condensation and clathrate formation of most abundant atmophile species in ISM and comets of the solar system: ${{\rm{H}}}_{2}{\rm{O}}$ , CO, CO₂, ${\mathrm{CH}}_{3}\mathrm{OH}$ , CH₄, N₂, NH₃, and ${{\rm{H}}}_{2}{\rm{S}}$ . In their model, these species either condense or form clathrates when the equilibrium pressures become lower than or equal to their partial pressures (i.e., total pressure of the disk times the abundance of the species). They have also taken into account the effects of refractory organics. However, they found that the ice/rock ratio of observed comets is more consistent with that of icy planetesimals estimated without these organics. Thus, we focus on the cases without refractory organics here.

Figure 1 shows the location changes of icelines for these molecules in a Chambers' disk model. The iceline for serpentine ( ${\mathrm{Mg}}_{3}{\mathrm{Si}}_{2}{{\rm{O}}}_{5}{(\mathrm{OH})}_{4}$ ) is also plotted for comparison. The closest iceline is for water, but icelines for carbon and nitrogen such as ${\mathrm{CH}}_{3}\mathrm{OH}$ and NH₃ also move near the rocky planet formation region as the disk evolves and the temperature decreases. The icelines for CO and N₂ and their clathrates are beyond 10 AU all the time. The figure indicates that, as the disk evolves, the building blocks of terrestrial planets are more likely to include atmophile species such as C and N. Note that in a Chambers' disk model the outer disk is dominated by stellar irradiation rather than by viscous heating, and thus the outer disk temperature stays the same for constant stellar properties as the disk pressure decreases. The irradiation dominated region also expands into the inner region as the disk evolves (see Figure 1 of Chambers 2009). This is why some of the icelines (e.g., CH₄ clathrates) move inward and then outward with time. Also note that the standard planet formation expects a factor of a few increase in the surface mass density beyond the water iceline, and thus the mass of embryos (e.g., Kokubo & Ida 2000; Morbidelli et al. 2015b). However, we do not change the mass of embryos across the water iceline for the ease of the comparison of models, so that all the embryos have the same mass in all the disk models we test.

3. RESULTS

3.1. Dynamical Results

Before presenting the elemental abundances of planets generated in our simulations, we first discuss the dynamical outcomes of Grand Tack, CJS, and EJS models. These models have been explored by other authors (e.g., O'Brien et al. 2006, 2014; Walsh et al. 2011; Jacobson et al. 2014), and our results overall agree with theirs. However, we discuss the outcomes of these models and compare them with one another, because they give us clues to understand the compositional differences of planets, which we discuss in Section 3.2.

3.1.1. Overall Results

All the models generate a few planets per system in ∼10–100 Myr, which is consistent with the Earth's formation timescale as estimated from radioactive elements (10–150 Myr, e.g., Halliday & Kleine 2006; Kleine et al. 2009). The overall agreement with properties of the solar system is shown in Figure 2, and the average properties are summarized in Table 3.

**Figure 2.** The comparison of results of GT1.5AU and GT2.0AU (top panels) and CJS and EJS (bottom panels). The left panels show the mass distributions of planets and right panels show the separation of the neighboring pairs in the mutual Hill radius as a function of the average eccentricities of the pairs. The red stars represent corresponding results for Mercury, Venus, Earth, and Mars.
Download figure:
Standard image High-resolution image

Table 3. Outcomes of Runs for Embryos

Run Name	e_ave	i_ave (degree)	${K}_{{\rm{ave}}}({R}_{{\rm{mut,}}{\rm{Hill}}})$	S	M	E	CE	CJ	CS
GT1.5AU	0.024	1.26	28.0	19.7	0.31	50.3	16.3	12.5	0.94
GT2.0AU	0.033	1.60	28.1	24.7	0	40.0	29.4	5.3	0.63
CJS	0.045	2.65	34.8	28.8	0	2.5	68.8	0	0
EJS	0.042	6.39	45.0	19.0	19.0	9.0	53.0	0	0
solar system	0.096	5.01	43.3	...	...	...	...	...	...

Note. Column 1: name of a set of runs; column 2–4: the average eccentricity, inclination, and separation of a neighboring pair in the mutual Hill radius, respectively; column 5–10: percentages of fates of embryos: survived, merged with the central star, ejected out of the system, collided with another embryo, collided with Jupiter, and collided with Saturn, respectively. The information for the solar system is added for comparison.

Download table as: ASCII Typeset image

In all of these models, the mass distributions have a peak near Venus and Earth locations. However, for the CJS model, the mass distribution is much broader than the other models, which indicates that the circular giant planet assumption is inconsistent with the terrestrial planet formation in the solar system. The mass distributions from GT1.5AU and GT2.0AU are broadly similar to each other. For the further tack location (2.0 AU), we find that the mass distribution becomes slightly broader than the default case and that the location of the peak also shifts slightly outward. The production of Mars analogues is more common for GT2.0AU. We investigate the effects of tack locations on the planetary mass distribution in our future paper (R. Brasser et al. 2016, in preparation).

A general trend for the eccentricity and inclination distributions is similar for all the models as well. For GT1.5AU and GT2.0AU, the average eccentricities are 0.024 and 0.033, while the average inclinations are 1 fdg 3 and 1 fdg 5, respectively. These are lower than the values for the terrestrial planets in the solar system (e_ave ∼ 0.096 and i_ave = 5 fdg 01), and thus these planets have more circular and aligned orbits. The eccentricities and inclinations of the terrestrial planets could increase as giant planets go through the NICE-model phase of resonance crossing and migration (Brasser et al. 2013). Both CJS and EJS models also have lower average eccentricities than the solar system, and lower or comparable average inclinations (see Table 3).

The models also reproduce the overall trend of the separation of neighboring pairs as a function of the average eccentricities of the pairs, where the separation increases with the average eccentricity. To calculate the separation K in a mutual Hill radius, we define the mutual Hill radius as the average of the Hill radii of two planets. For the solar system, the separation of the neighboring pairs is K = 63.4, 26.3, and 40.1 from inward to outward, and the average value is K_ave = 43.3 (see Table 3). For GT models, the final planetary systems tend to be more compact than the solar system, and the average separation is K_ave ∼ 28 for both GT1.5AU and GT2.0AU. The average separation for CJS is larger than GT (K_ave ∼ 34), and that for EJS (K_ave ∼ 45) is comparable to the solar system value. The larger average separations for CJS and EJS models compared to the GT models are consistent with their slightly higher eccentricities and inclinations.

3.1.2. Dynamical Fates and the Mixing of the Inner Disk

Although the distributions of planetary masses are similar for GT and EJS models (see Figure 2), the dynamical fates of embryos are markedly different for these models, as seen in Figure 3. The break down of the various fates of embryos is summarized in Table 3. In GT models, most survived embryos (purple) were initially within the tacking radius of Jupiter (1.5 AU or 2.0 AU), while in EJS and CJS models survived planets come from a wider range of orbital radii. A difference in the initial locations of the survived embryos indicates a potential importance of tacking locations in determining their compositions. In Section 3.2, we show that these differences do not lead to significant changes in the elemental composition of planets. In these models, embryos colliding with other embryos (red) come from similar regions to survived ones.

**Figure 3.** The comparison of the dynamical outcomes of embryos of GT1.5AU (left panels) and EJS models (right panels). The top panels show these outcomes as a function of the initial semimajor axis of an embryo, while the bottom panels are plotted against the minimum semimajor axis reached during the simulation.
Download figure:
Standard image High-resolution image

A characteristic difference between GT and EJS models is seen among the fates of the other embryos that are removed from the systems. In GT simulations, some embryos collide with Jupiter and a few collide with Saturn, but the vast majority are ejected out of the system (see Table 3). The peak in the number of embryos that collided with Jupiter is seen near the minimum semimajor axis comparable to the Jupiter's tack location, which indicates that most collisions occur near this radius. Both ejection and collision in this scenario can be explained as a result of close encounters of embryos with the giant planets. On the other hand, in EJS simulations, most other embryos merge with the central star, and only a few are ejected out of the system. The merger with the central star in EJS is likely due to secular resonances induced by Jupiter and Saturn, which gradually increase the eccentricities of embryos by keeping their semimajor axes constant (e.g., Matsumura et al. 2013). The lack of collisions with Jupiter and Saturn for EJS makes sense because embryos in this model (and CJS) are far away and scattering is more effective further from the star.

The corresponding result for the CJS model is similar to that of EJS, and survived embryos come from all over the inner disk. An interesting difference between EJS and CJS models is that most embryos in the CJS model collided with each other. No embryos collide with Jupiter or Saturn, or merge with the central star, which indicates that planet formation proceeds without significant gravitational perturbations from the giant planets.

The dynamical fates of inner planetesimals for the GT1.5AU and EJS models are compared in Figure 4, and the breakdowns for all the models are summarized in Tables 4 and 5. Similar to embryos, a merger with the central star (∼41%) is a more common outcome than ejection (∼16%) for inner planetesimals in EJS model, but ejection is more common in GT models (∼26% compared to a few percent). A notable difference from the trend of embryos is that the planetesimals that collide with embryos are from all over the inner disk for GT, while they tend to come from the innermost disk <1.5 AU for EJS. Because the minimum semimajor axes of collided planetesimals are ≲1 AU for GT1.5AU, it is expected that these planetesimals mostly collided with embryos that eventually became planets. These trends indicate that the inner planetesimals are well mixed in the GT model. On the other hand, because the planetesimals that ended up in planets appear to come most commonly from the innermost region in EJS, it is expected that the planetesimals are less well mixed in EJS compared to GT. Indeed, as shown in Table 5, most planetesimals stay in the initial region for EJS, while ∼50% of planetesimals from 1 to 2 AU and ∼20% from 2 to 3 AU arrive within 1 AU for GT models.

Table 4. Outcomes of runs for planetesimals

Run Name	${{\rm{S}}}_{\mathrm{in}}$	${{\rm{M}}}_{\mathrm{in}}$	${{\rm{E}}}_{\mathrm{in}}$	${\mathrm{CE}}_{\mathrm{in}}$	${\mathrm{CJ}}_{\mathrm{in}}$	${\mathrm{CS}}_{\mathrm{in}}$	${{\rm{S}}}_{\mathrm{out}}$	${{\rm{M}}}_{\mathrm{out}}$	${{\rm{E}}}_{\mathrm{out}}$	${\mathrm{CE}}_{\mathrm{out}}$	${\mathrm{CJ}}_{\mathrm{out}}$	${\mathrm{CS}}_{\mathrm{out}}$
GT1.5AU	1.5	1.9	25.9	65.5	4.8	0.4	1.6	0.2	89.9	0.9	4.0	3.3
GT2.0AU	2.0	4.6	25.9	61.8	3.5	0.3	1.0	0.8	91.8	0.8	3.3	2.3
CJS	3.0	3.3	33.3	59.3	0.9	0.08	5.1	0.6	85.2	1.4	6.9	0.8
EJS	2.2	41.3	15.9	40.2	0.3	0.04	2.0	1.5	93.2	0.05	2.8	0.4

Note. Column 1: name of a set of runs; column 2–7: mass percentages of inner planetesimals that survived, merged with the central star, ejected out of the system, collided with another embryo, collided with Jupiter, and collided with Saturn, respectively; column 8–13: the corresponding mass percentages of outer planetesimals.

Download table as: ASCII Typeset image

Table 5. Initial and Final Locations of Planetesimals

Run Name	a_in (AU)	a_fin < 1 (AU)	a_fin = 1–2 (AU)	a_fin = 2–3 (AU)	a_fin ≥ 3 (AU)
GT1.5AU	<1	77.7	14.3	6.0	2.0
	1−2	58.3	20.2	7.0	14.5
	2−3	23.5	7.3	6.8	62.4
	≥3	0.7	0.4	2.9	96.0
GT2.0AU	<1	76.1	15.3	2.4	6.3
	1−2	48.1	27.5	4.8	19.6
	2−3	21.7	11.8	4.2	62.2
	≥3	0.7	0.5	1.1	97.8
CJS	<1	70.5	18.3	2.1	9.1
	1−2	26.3	35.3	7.4	31.1
	2−3	9.5	18.3	10.4	61.8
	≥3	0.5	0.9	1.2	97.4
EJS	<1	73.9	16.4	6.2	3.5
	1−2	26.2	32.4	31.6	9.8
	2−3	2.1	9.5	51.1	37.2
	≥3	0	0.03	0.87	99.1

Note. Column 1: name of a set of runs; column 3–6: mass percentages of planetesimals with the initial semimajor axis a_in < 1 AU, 1–2 AU, 2–3 AU, and ≥3 AU that end up with the final semimajor axis a_fin < 1 AU, 1–2 AU, 2–3 AU, and ≥3 AU.

Download table as: ASCII Typeset image

The distribution for CJS looks similar to GT, and planetesimals from all over the inner disk collide with embryos. However, the percentages of planetesimals that end up within 1 AU are much smaller than GT, and ∼26% and ∼10% for planetesimals from 1 to 2 AU and 2 to 3 AU, respectively. Thus, it appears that the radial mixing of materials is most efficient in GT models, less efficient in CJS, and least efficient in EJS. In Section 3.2, we show that this leads to differences in compositions.

The dominant fate of outer planetesimals for all the models is ejection out of the system. However, percentages of outer planetesimals that reach the inner disk are about one order of magnitude lower in EJS compared to the other models for each region (see Table 5). Thus, the scattering of outer planetesimals appears to be much less efficient in EJS compared to the others. The scattering sources of outer planetesimals are different for GT and CJS models, and we will discuss this below.

3.1.3. Collisional Speed and the Scattering of Outer Planetesimals

The collisions between two bodies could lead to agglomeration or destruction, depending on the magnitude of the encounter speed relative to the escape speed. Marcus et al. (2009) studied collisions between terrestrial-size planets (1−10M_E) via smoothed particle hydrodynamics simulations and showed that the erosive outcome is more common when the impact speed is larger than twice the escape speed (v_imp/v_esc > 2), although the outcome strongly depends on the impact angle.

In Figures 5, 6, we show the ratios of the relative encounter speeds of colliding bodies to mutual escape speeds as a function of simulation time for GT1.5AU, CJS, and EJS models. Here, the mutual escape speed is defined as ${v}_{{\rm{esc}}}\;=\;\sqrt{2G({M}_{1}+{M}_{2})/({R}_{1}+{R}_{2})}$ . For the GT1.5AU model, all embryo–embryo collisions but one have a speed ratio of less than 1.5, which indicates the agglomeration rather than the erosion. For embryo–planetesimal collisions, the ratio increases with time, and some collisions may be destructive. There is no dependence of this ratio on the initial locations of planetesimals, because both inner and outer planetesimals (defined as planetesimals within and beyond 3 AU here) show similar distributions.

**Figure 5.** The ratio of a relative encounter speed to a mutual escape speed as a function of time for embryo–embryo collisions for the GT1.5AU (top), CJS (middle), and EJS (bottom) models.
Download figure:
Standard image High-resolution image

**Figure 6.** The ratio of a relative encounter speed to a mutual escape speed as a function of time for embryo–planetesimal collisions for the GT1.5AU (top), CJS (middle), and EJS (bottom) models.
Download figure:
Standard image High-resolution image

The outcomes for the CJS model are overall similar to the GT model. We find that all embryo–embryo collisions have a speed ratio less than 1.4, which indicates agglomeration. The overall shape of the speed ratio for embryo–planetesimal collisions is similar to the GT case as well, with some of them being outer planetesimals. This indicates that, even without a giant planet migration, the CJS case has an efficient inward scattering of outer planetesimals. The critical difference, however, is the timing of the scattering. For the GT case, the outer planetesimals are scattered inward and collide with embryos throughout the simulations, from 0.1 to 100 Myr. On the other hand, for the CJS case, most inward scattering occurs after 10 Myr. This difference arises because, in the CJS model, the scattering is done by growing terrestrial planets rather than by giant planets. Indeed, Figure 2 shows that more massive terrestrial planets tend to be at large orbital radii (i.e., near the asteroid belt) for the CJS model compared to the others.

The EJS model shows a similar trend for the embryo–embryo collision to the other cases, but the relative velocity varies more than GT or CJS models. This may be due to their slightly higher eccentricities and inclinations compared to the other models (see Table 3). The distribution of embryo–planetesimal collisions is similar to that of the GT and CJS models. However, in the EJS model there are almost no outer planetesimals that collide with embryos. The scattering is inefficient in EJS since massive rocky planets do not exist in the asteroid belt region, as in CJS, because they tend to be removed by secular effects (see Figures 2 and 3).

An interesting difference, as shown in Figure 5, is that the major impact starts and finishes earlier in GT models (∼0.1–10 Myr) compared to CJS and EJS models (∼1–100 Myr). The early occurrence of giant impacts is due to the migration of giant planets. However, as we discuss below, the timing of the last giant impact depends on the initial mass ratios of planetesimals and embryos.

In summary, we find that the scattering of outer planetesimals is inefficient in EJS and efficient in CJS and GT. In GT, the giant planets scatter planetesimals from the outer disk to the inner disk; whereas in CJS the terrestrial planets formed near the asteroid belt are responsible for the scattering. As we show in Section 3.2, the difference in the scattering efficiency leads to the difference in abundances of highly volatile species.

3.1.4. Moon-forming Timescale

For completeness, we estimate the Moon-forming timescale in our simulations. Jacobson et al. (2014) calculated the late accreted mass of a planet, which is the total mass of planetesimals accreted onto an embryo after the last giant impact between embryos, and found that the late accreted mass decreases with the time of the last giant impact. They argued that the correlation could be used to estimate the timing of the Moon-forming event by comparing the estimated late accreted mass with that of Earth, and proposed that the Moon-forming event could have happened long after 40 Myr and around 95 ± 32 Myr. Figure 7 is the corresponding figure for our simulations. The estimated time is ∼109 Myr for Grand Tack simulations, and is consistent with their estimate. Our simulations did not produce many planets that experienced the last impact at such late times, compared to Jacobson et al. (2014). This is because our simulations assume the same total mass for embryos and planetesimals, while Jacobson et al. (2014) includes various mass ratios. As they showed in the supplementary information, the time of the last giant impact increases when the dynamical friction effect is reduced due to the smaller total mass in planetesimals compared to embryos.

**Figure 7.** The late accreted mass normalized by the planetary mass as a function of the last impact time between embryos. The late accreted mass is the mass accreted onto an embryo after the last giant impact between embryos. The horizontal solid and dashed lines are the best estimate and 1σ uncertainty of the late accreted mass estimated from the highly siderophile element abundances in the mantle: 4.8 ± 1.6 × 10^− 3 M_⊕ (Jacobson et al. 2014). The intersection of the best estimate and the polynomial fit to the combined data (104.34 Myr) is shown as a vertical red line. The intersection is 109.04 Myr for only the results from GT1.5AU and GT2.0AU.
Download figure:
Standard image High-resolution image

For CJS and EJS simulations, we find a similar trend to GT runs. By taking account of all the simulations, the estimated time of the Moon-forming event is ∼104 Myr. The trend that CJS and EJS both give slightly earlier times than GT is the opposite of Jacobson et al. (2014), because their corresponding simulations point to a comparable, but overall longer, last giant impact time for the Moon formation. It is unclear whether this is a statistically significant effect. In R. Brasser et al. (2016, in preparation), we will study various initial conditions for the GT model by taking into account type I migration and the size distribution of embryos to investigate the late accretion problem further.

3.2. Elemental Abundances of Planets

Overall, the dynamical results imply that (1) the radial mixing of materials is most efficient in GT and less so in CJS and EJS models, and (2) icy outer planetesimals could be scattered into the inner disk region and thus may be able to contribute to the delivery of atmophile species to rocky planets in GT and CJS models, but not in EJS model. In this subsection, we discuss how these differences affect the elemental abundances of planets. All the simulated planets will be compared with the composition of the Earth. This is partly because the composition of the Earth is the best known of all terrestrial planets, and the estimates for the other planets, except Mercury, are similar to those of the Earth. The composition of the bulk Earth used in this study is taken from McDonough (2003) and shown in Table 6. We discuss the abundances of refractory and moderately volatile elements, water, and more volatile C and N elements, separately.

Table 6. Elemental Abundances used for HSC Chemistry

Element Name	Earth	H	LL	CI	CM	EH	EH
H	260 ppm	...	...	2.0	1.4	...	...
C	730 ppm	0.11	0.12	3.2	2.2	0.40	0.36
N	25 ppm	48 ppm	70 ppm	0.15	0.152	...	...
O	29.7	35.7	40.0	46.0	43.2	28.0	31.0
Na	0.18	0.64	0.70	0.49	0.41	0.68	0.58
Mg	15.4	14.0	15.3	9.7	11.7	10.6	14.1
Al	1.59	1.13	1.19	0.86	1.18	0.81	1.05
Si	16.1	16.9	18.9	10.5	12.9	16.7	18.6
P	0.121	0.108	0.085	0.102	0.090	0.200	0.117
S	0.635	2.0	2.3	5.9	3.3	5.8	3.3
Ca	1.71	1.25	1.30	0.92	1.27	0.85	1.01
Ti	813 ppm	600 ppm	620 ppm	420 ppm	580 ppm	450 ppm	580 ppm
Cr	0.466	0.366	0.374	0.265	0.305	0.315	0.305
Fe	31.9	27.5	18.5	18.2	21.0	29.0	22.0
Ni	1.82	1.60	1.02	1.07	1.20	1.75	1.30

Note. Column 1: name of the element; and columns 2–8: elemental concentration (ppm = 1.e − 4%) for Earth and various chondrites, where the values are in % unless stated otherwise. The values for Earth are taken from Table 3 of McDonough (2003) and through a private communication with Bill McDonough in 2015. The values for chondrites are taken from Wasson & Kallemeyn (1988).

Download table as: ASCII Typeset image

3.2.1. Refractory and Moderately Volatile Elements

Figure 8 shows the elemental abundance comparisons for all the planets formed in GT1.5AU. The abundances of refractory and moderately volatile elements are listed roughly in increasing order of volatility, and normalized by the corresponding abundances of the Earth for four different initial disk ages of ${t}_{{\rm{init,disk}}}\;=\;0$ , 0.1, 0.2, and 1 Myr. The abundance ratio for an element X is defined as follows.

$\begin{eqnarray}&&{N}_{X}\;=\;\displaystyle \frac{{\rm{wt.\%}}\;{\rm{of}}\;{\rm{X}}}{{\rm{wt.\%}}\;{\rm{of}}\;{\rm{Si}}}.\end{eqnarray} \tag{ 1 }$

In our model, we find the best agreement for t_init,disk = 0.2 Myr. The abundances of refractory and moderately volatile elements agree with those of the Earth very well, except for some of the most volatile elements in this figure ( ${\rm{N}}{\rm{a}}$ and S). This trend is consistent with previous works (e.g., Bond et al. 2010; Elser et al. 2012; Moriarty et al. 2014). The agreement could be improved by taking account of a volatile loss, but such a study is beyond the scope of this paper.

**Figure 8.** The abundances of refractory and moderately volatile elements normalized by those of the Earth for all survived embryos in GT1.5AU. The abundance ratio N_X is taken with respect to Si, thus the ratio for Si is always one. Different panels correspond to different initial disk ages.
Download figure:
Standard image High-resolution image

Figure 9 shows the overall agreement of the abundances of refractory and moderately volatile elements of planets from GT1.5AU with those of Earth for different initial disk ages (i.e., different initial compositions of embryos and planetesimals). The error in the abundance ratio with respect to the Earth is calculated for each element, except for Na and S, and the average error was determined for each planet. Na and S were not included here because they tend to give the largest errors, as shown in Figure 8. The solid line indicates the median value of the calculated average errors for each t_init,disk. The overall average error becomes the smallest for the initial disk age of 0.2 Myr, and increases on both sides. For earlier times, the agreement is better for embryos that were initially further out, while for later times, the agreement is better for closer-in embryos. The figure indicates that, if the initial disk was cold enough ( ${t}_{{\rm{init,disk}}}\geqslant 0.3$ Myr), there is little difference in abundances of refractory and moderately volatile elements for planets over a wide range (∼0.5–2.0 AU). The corresponding figures for GT2.0AU, CJS, and EJS models are all very similar to Figure 9. Therefore, we conclude that we cannot distinguish different formation models from refractory and moderately volatile elements alone.

3.2.2. Hydrous Species

To investigate the abundance differences of more volatile species, we estimate the amount of water accreted onto a planet from hydrous species (water ice and serpentine ${\mathrm{Mg}}_{3}{\mathrm{Si}}_{2}{{\rm{O}}}_{5}{(\mathrm{OH})}_{4}$ ). Following Bond et al. (2010), we assume that all the hydrogen accreted as hydrous species is converted to water. The results are plotted in Figure 10 for GT1.5AU and EJS models. The solid lines indicate the median values of a water mass fraction per planet normalized by that of the Earth. In both cases, the estimated amount of water increases with the initial disk age. The critical disk age to get water comparable to the Earth is ∼0.2 Myr and ∼0.3 Myr for GT and EJS, respectively. The results for the CJS model are similar to the GT model. As seen in Figure 1, the water iceline moves within 3 AU for ${t}_{{\rm{init,disk}}}\gtrsim 0.3$ Myr while the serpentine iceline is slightly inward to water and moves into the inner disk for ${t}_{{\rm{init,disk}}}\gtrsim 0.2$ Myr. Because the hydrous species are included in some of the embryos and inner planetesimals beyond this time, the water delivery in EJS is largely explained by the inward migration of the icelines. This is consistent with the lack of scattering of outer planetesimals into the inner disk in EJS model (see Figure 6). Although the accreted fraction of water may seem fine for EJS, the model is not favorable to explain the water delivery. This is because the requirement of the icelines to be inside the inner disk region implies the formation of larger, icy cores, which are different from those expected for terrestrial planets.

**Figure 10.** A fraction of the water mass of a planet is normalized by that of the Earth for GT (blue) and EJS (orange) models, by assuming that the Earth's water mass is 2.34 × 10⁻⁴ M_E. The solid lines show the median values of the water content for these cases. For t_init,disk = 0.2 and 0.3 Myr, planets in GT are accreting about one order of magnitude more water compared to those in EJS for the same initial disk age.
Download figure:
Standard image High-resolution image

On the other hand, for GT and CJS, outer planetesimals are clearly contributing to the delivery of water, because the water accretion becomes significant earlier than EJS. Our results for the EJS and CJS cases are in agreement with the results indicated in Bond et al. (2010). The comparable water accretion for GT and CJS indicates that both the migration of giants and the presence of terrestrial planets in the asteroid belt region can help deliver water via scattering of outer planetesimals.

We also find that there is a correlation between the water mass fraction and final semimajor axis. Figure 11 shows these results for four models by using the initial disk age of ${t}_{{\rm{init,disk}}}\;=\;0.4$ Myr. The results are qualitatively similar for other initial disk ages. For the GT1.5AU and GT2.0AU cases, we find that the amount of accreted water is nearly constant independent of the initial semimajor axes of the survived embryos. This agrees with our expectation from Figure 4 that the materials are well mixed in the GT model. For the EJS and CJS cases, on the other hand, the amount of water is weakly dependent on the initial semimajor axis of an embryo, and the inner embryos tend to have a smaller amount of water. For ${t}_{{\rm{init,disk}}}\;=\;0.4$ Myr, the water iceline would be ∼2.4 AU and thus embryos and planetesimals beyond this radius are expected to initially contain water. The trend in water mass confirms that these water-rich materials did not reach the innermost disk as efficiently as GT models in CJS and EJS models, and that the materials are less well mixed in these scenarios.

**Figure 11.** A fraction of the water mass of a planet normalized by that of the Earth is plotted as a function of the final semimajor axis of a survived embryo for GT1.5AU (top left), GT2.0AU (top right), CJS (bottom left), and EJS (bottom right) models. These are calculated for the Chamber's disk model at ${t}_{{\rm{init,disk}}}\;=\;0.4$ Myr. The vertical dashed lines indicate the location of the water iceline for this disk age (∼2.4 AU).
Download figure:
Standard image High-resolution image

3.2.3. Carbon and Nitrogen

We also check the abundances of even more volatile species: carbon and nitrogen. Figure 12 shows how much C and N are delivered to terrestrial planets in GT1.5AU. For the case of carbon, we find that an abundance ratio comparable to Earth or larger is delivered for the initial disk age of ≥0.6 Myr. Because the iceline for ${\mathrm{CH}}_{3}\mathrm{OH}$ enters within 3 AU for ${t}_{{\rm{init,disk}}}\geqslant 0.4$ Myr, this indicates that the C is incorporated into planets as a result of the iceline evolution. For younger initial disks, the icelines of C species are further than 3 AU, and therefore the delivery of C should be done solely by outer planetesimals. Because the abundances of C for ${t}_{{\rm{init,disk}}}\lt 0.4$ Myr are ∼2 orders of magnitude smaller than that of the Earth, we find that our default outer planetesimal disk is not massive enough to explain the abundance of C in the Earth.

On the other hand, the accreted amount of N is nearly independent of the initial disk age, indicating that all the N delivery is due to the scattering of outer planetesimals. This is consistent with the expectation from Figure 1; because the NH₃ iceline never enters within 3 AU, N is never incorporated initially within embryos and inner planetesimals. However, the accreted amount of N is about an order of magnitude too low compared to the Earth in our model.

The results of the CJS model are similar to the GT results and are consistent with the expectation from Figure 6, where the efficient incorporation of outer planetesimals into embryos occurs for GT and CJS models. The corresponding results for EJS are shown in Figure 13. Because there is almost no scattering of outer planetesimals to the inner disk, the accretion of C occurs after the iceline enters 3 AU ( ${t}_{{\rm{init,disk}}}\geqslant 0.4$ Myr) and the accretion of N occurs very rarely. We confirm that the trend is similar to GT models without outer planetesimal disks (i.e., no scattering of these planetesimals into the inner disk region). However, the variations of C and N abundances for the same disk age are much larger in the EJS model due to a similar radial dependence to Figure 11. Also, the median value of C becomes comparable to the Earth's amount for ${t}_{{\rm{init,disk}}}\geqslant 0.6$ Myr, while that occurs for ${t}_{{\rm{init,disk}}}\geqslant 0.8$ Myr in EJS. This is also a consequence of the efficient radial mixing in GT models compared to EJS.

**Figure 13.** The abundance ratio of C and N of a planet normalized by that of the Earth for EJS. The solid lines show the median values of the C and N contents. The carbon comparable to Earth or larger is delivered for the initial disk age of ${t}_{{\rm{init,disk}}}\geqslant 0.8$ Myr. The median line for nitrogen does not show up because the accretion is inefficient.
Download figure:
Standard image High-resolution image

**Figure 13.** The abundance ratio of C and N of a planet normalized by that of the Earth for EJS. The solid lines show the median values of the C and N contents. The carbon comparable to Earth or larger is delivered for the initial disk age of ${t}_{{\rm{init,disk}}}\geqslant 0.8$ Myr. The median line for nitrogen does not show up because the accretion is inefficient.
Download figure:
Standard image High-resolution image

It may be possible to accrete larger amounts of C and N via more efficient gas drag on the outer planetesimals. For the GT simulations shown so far, we have assumed a planetesimal size of 50 km for gas drag purposes. However, smaller planetesimals would migrate much faster, and thus icy materials could be incorporated into planets more efficiently. We assume three different planetesimal sizes (1, 10, and 100 km) for gas drag, and run eight simulations for each size until the gas disk dissipates. We find that the mass percentages of outer planetesimals that arrive within 1.5 AU increase with decreasing planetesimal radii: 2.9%, 1.9%, 0.75%, and 0.45% for 1, 10, 50, and 100 km, respectively. Thus, if planetesimals are 1 km in size, about four times more icy materials could be delivered to rocky planets. However, as shown in Figure 12, the amount of C is lower by ∼2 orders of magnitude compared to Earth for default cases. Therefore, a more realistic size distribution of planetesimals is unlikely to resolve this issue.

To check whether it is possible to explain the delivery of C and N solely by the scattering of outer planetesimals, we assume that the outer planetesimal mass is a factor of ∼20 larger and is comparable to the inner one (2.3 × 10⁻³ M_E). This means that the outer planetesimal disk mass changes from ∼0.06 M_E to ∼1.2 M_E. In this case, the accreted amounts of C and N increase and become comparable to those of Earth, nearly independent of the initial disk age. This change in the mass of outer planetesimals also affects the amount of water accreted on planets. However, the effect is largely limited to the initial disk's age below 0.2 Myr, where the accreted amount of water becomes comparable to the amount of the Earth's water. The accreted amounts of C and N show a similar dependence on a semimajor axis to the case of water in Figure 11.

The abundance of atmophile elements could potentially distinguish GT, CJS, and EJS models. For the EJS model, the delivery of atmophile species via outer planetesimals is inefficient. Therefore, the main building blocks of terrestrial planets (i.e., embryos and inner planetesimals) would need to include these species due to the iceline evolution. However, such a scenario could lead to the larger errors for the abundances of more refractory species, as well as the existence of larger, icy planetary embryos in the inner disk. For GT and CJS models, the delivery of atmophile species could be aided by outer planetesimals scattered into the planet formation region. Distinguishing these two models is harder, but the timing of accretion may be later for the CJS model than GT (see Figure 6), and the abundance may be dependent on the distance from the Sun for CJS and independent for GT (see Figure 11).

In summary, to satisfy the constraints from both the refractory elements and the amount of water, the initial disk age of ∼0.2 Myr is most preferable. This in turn implies that, if the C and N are incorporated into rocky planets during their formation, the scattering of the outer planetesimals, rather than the iceline evolution, would be responsible for the delivery of these species.

3.3. Elemental Abundances of Planetesimals

Finally, we discuss the distributions and compositions of survived planetesimals. Figure 14 shows the final distribution of planetesimals for the GT1.5AU case. Most survived planetesimals are near the asteroid belt location (∼1.8–3.5 AU), and inner planetesimals (orange), which are initially within 3 AU, are located closer to the star compared to outer planetesimals scattered into the inner disk (blue). The average eccentricity is 0.288 and the average inclination is 17 fdg 7, and thus is more excited compared to the observed asteroids. These trends are mostly consistent with Walsh et al. (2011).

To compare the compositions of survived planetesimals with those of chondrites, we choose planetesimals with semimajor axes within 1.8–3.5 AU, eccentricities ≤0.5, and inclinations ≤40° at the end of the simulations, which cover roughly the same ranges of orbital parameters to asteroids. Figures 15 and 16 show the elemental abundances of the survived inner and outer planetesimals, respectively, compared to H, LL, CI, CM, EH, and EL chondrites. Because the inner planetesimals are more rocky than the outer ones, it is expected that they have more similar properties to S-type chondrites (H and LL) rather than C-type (CI and CM) or E-type (EH and EL) chondrites. These figures show that inner planetesimals have a good agreement with H chondrites, except for S, C, and N. Note that, because we ignore the interactions of planetesimals, their compositions do not change during the simulations due to collisions. Also, for t_init,disk = 0.2 Myr, the C and N icelines are not within the inner disk. The abundances of O of inner planetesimals tend to be depleted compared to CI and CM. Compared to LL, EH, and EL chondrites, the survived inner planetesimals tend to have overestimated refractory elements.

**Figure 15.** The abundances of elements in inner planetesimals normalized by those of H (top left panels), LL (top right), CI (middle left), CM (middle right), EH (bottom left), and EL (bottom right) chondrites. The initial disk age of 0.2 Myr is assumed here. Similar to Figure 8, the abundance ratio N_X is taken with respect to Si, thus the ratio for Si is always one.
Download figure:
Standard image High-resolution image

**Figure 16.** The abundances of elements in outer planetesimals normalized by those of H (top left panels), LL (top right), CI (middle left), CM (middle right), EH (bottom left), and EL (bottom right) chondrites. The initial disk age of 0.2 Myr is assumed here. Similar to Figure 8, the abundance ratio N_X is taken with respect to Si, thus the ratio for Si is always one.
Download figure:
Standard image High-resolution image

On the other hand, the outer planetesimals show a good agreement with CI chondrites except for N, while the abundances of O, S, C, and N tend to be enhanced compared to H, LL, EH, and EL chondrites. Therefore, our results show that inner planetesimals are more consistent with S-type chondrites while outer planetesimals are more consistent with C-type chondrites. These agree with the expectation from the Grand Tack scenario (Walsh et al. 2011).

4. DISCUSSIONS AND CONCLUSIONS

We studied how the planetary compositions change depending on the dynamical effects of giant planets by testing three popular setups of terrestrial planet formation in the solar system (GT, CJS, and EJS). The abundances of refractory and moderately volatile elements in terrestrial planets are nearly independent of orbital radii for ${t}_{{\rm{init,disk}}}\geqslant 0.2$ Myr for all the models (see Figure 9 and the discussion in text). Therefore, we cannot distinguish different formation models from these elements. The mass fractions of water and other more volatile species (e.g., C and N), on the other hand, can be used to distinguish different models.

We found that there are two key mechanisms to distinguish different models: the radial mixing of the inner planetesimal disk and the scattering of outer, icy planetesimals into the inner disk. The former tends to homogenize the planetary compositions over the inner disk, while the latter helps deliver the atmophile species, which only condense in the cold region of a disk, to terrestrial planets. The dynamical simulations show that the mixing is efficient in GT models but less so in CJS and least efficient in EJS (see Table 5). On the other hand, the scattering of outer planetesimals is efficient both in GT and CJS models, while it is inefficient in EJS (see Figure 6). These effects together lead to the following trends for each formation model.

In the classical formation model (EJS), neither the scattering nor the mixing are efficient. This makes it very difficult to deliver atmophile species to terrestrial planets. Thus, if terrestrial planets are formed in the presence of giant planets with moderately eccentric orbits, we expect that such planets would be dry. It is still possible for these planets to receive atmophile species as the disk cools and the icelines move inward. However, in such a case, the formed planets may not be terrestrial-size but rather like mini-Neptunes, because their cores may be born beyond the iceline and thus grow larger than typical rocky planets. Due to the poor radial mixing, planets in such a system may have a larger fraction of atmophile species for larger orbital radii.

In the GT model, both the scattering and the mixing are efficient. The atmophile elements can be efficiently delivered to terrestrial planets, and the mass fractions of these elements would be at least initially nearly independent of orbital radii. For the solar system case, the GT model can reproduce the overall abundances of most refractory and moderately volatile elements, the amount of water, and the abundances of C and N consistently, if we assume a massive enough outer planetesimal disk with a total mass $\gtrsim 1.2{M}_{E}$ . This disk condition further produces the asteroid analogues that have comparable compositions to S-type and C-type chondrites, where the S-type analogues tend to locate interior to the C-type ones, as shown in Walsh et al. (2011).

The CJS model is a "hybrid" of these two models, where the scattering is efficient but the mixing is not as efficient as the GT model. The model would allow terrestrial planets to form much further than the EJS and GT models, and the growing rocky planets on such orbits (rather than giant planets) could efficiently scatter the icy planetesimals into the inner disk. We calculated Tisserand's parameter with respect to Jupiter in our simulations, and found that the dynamics of planetesimals within ∼3 AU are not controlled by Jupiter. The atomophile elements would be delivered to terrestrial planets, but the mass fractions of these elements may increase with orbital radii, which is different from the GT model (see Figure 11). The CJS model is inappropriate for describing the formation path of the solar system because the delivery of water and other atmophile species relies on the presence of terrestrial planets in the asteroid belt region that do not exist now. The model has a prediction for some of the closely packed extrasolar planetary systems: the delivery of atmophile species to rocky planets may be possible even in systems that have no giant planets or those without strong perturbations from them. A recent study suggested that the inner planets of systems with no giant planets are in fact more water-rich because the radial drift of icy materials would not be stopped by giants (Morbidelli et al. 2015a).

Our model relies on a critical assumption that the compositions of building blocks of terrestrial planets are determined by the condensation sequence. However, the compositions of protostellar disks would be altered by many other mechanisms, such as transient heating events and transport of chemical species, and also affected by the presence of grains inherited from the molecular cloud core. Moreover, the loss of volatiles during the impacts (e.g., Okeefe & Ahrens 1982) or the erosion of embryos and planetesimals (Carter et al. 2015) could also affect the final compositions of planets. The former may help explain the overabundance of Na and S, which are common in this type of study (our work and, e.g., Bond et al. 2010; Elser et al. 2012; Moriarty et al. 2014). However, that would also affect the abundance of more volatile species, such as C and N. The inclusion of these effects will be important for the future study.

In this work, we did not study the compositions of gas giants. It has been proposed that the C/O ratio of a planetary atmosphere is important in characterizing the compositions of an exoplanet (Madhusudhan et al. 2012), and potentially distinguishing different formation models of hot Jupiters (Madhusudhan et al. 2014). Recent studies highlighted the importance of the locations of icelines (Öberg et al. 2011; Ali-Dib et al. 2014), as well as those of disk evolution and planet migration in determining the C/O ratio (Moriarty et al. 2014; Thiabaud et al. 2015). The future study of planet formation will need to incorporate chemical processes and dynamical evolution to better understand the solar system and extrasolar planetary systems.

We thank Kevin Walsh for providing us with SyMBA, including a Grand Tack model. We thank the anonymous referee for the careful reading and comments. S.M. thanks Prof. Bill McDonough, Prof. Qing-Zhu Yin, and Prof. Yuri Aikawa for useful discussions. S.M. is thankful for the support provided by the Kavli Institute for Theoretical Physics in the United States and the Earth-Life Science Institute in Japan, where some parts of this work were done, and a support from the Northern Research Partnership. R.B. is supported by the Astrobiology Centre Project of National Institutes of Natural Sciences (NINS) (Grant Number AB271017).

EFFECTS OF DYNAMICAL EVOLUTION OF GIANT PLANETS ON THE DELIVERY OF ATMOPHILE ELEMENTS DURING TERRESTRIAL PLANET FORMATION

Article metrics

Permissions

Author e-mails

Author affiliations

Author notes

Dates

ABSTRACT

1. INTRODUCTION