Dynamical Avenues for Mercury’s Origin. II. In Situ Formation in the Inner Terrestrial Disk

Matthew S. Clement; John E. Chambers

doi:10.3847/1538-3881/abfb6c

1. Introduction

The latter stages of terrestrial-planet formation in the solar system are thought to have unfolded as a series of giant impacts between a population of planetary-seed embryos engulfed in an ocean of smaller, D ∼ 1–100 km planetesimals (Wetherill & Stewart 1993; Weidenschilling et al. 1997; Kokubo & Ida 2002). Classic studies of this scenario (e.g., Chambers & Wetherill 1998; Chambers 2001; Raymond et al. 2004, 2006; O'Brien et al. 2006) tend to consider a total mass of terrestrial-forming material and surface density profile that are commensurate with the presumed minimum mass solar nebula (Weidenschilling 1977; Hayashi 1981) and generally consistent with modern models of dust evolution in protoplanetary disks (e.g., Birnstiel et al. 2012). While these types of models are successful at replicating numerous observed aspects of the modern terrestrial system, they systematically fail to generate several important features of the planet Mercury and its peculiar relationship to Venus. We provide a summary of the types of Mercury–Venus systems formed from a variety of common sets of initial conditions in past literature in Figure 1.

**Figure 1.** Summary of Mercury–Venus system mass and orbital period ratios in past dynamical studies. For the purposes of this plot, we simply define a Mercury–Venus system as any simulation finishing with exactly two planets inside of 0.85 au. The red star denotes the solar system value, and the different colored points correspond to simulation sets from the literature assuming different initial conditions as follows: black points—a set of 100 simulations from Clement et al. (2018) assuming the so-called classic initial conditions (∼5 M_⊕ of terrestrial-planet-forming material distributed uniformly between 0.5 and 4.0 au, e.g., Chambers & Wetherill 1998; Raymond et al. 2009). Here, Mercury forms directly within the massive disk in the same manner as Earth and Venus. Gray points—A batch of 100 simulations from Clement et al. (2019b) assuming identical classic initial conditions and also including a prescription (Section 2; Chambers 2013) accounting for the effects of collisional fragmentation. In successful realizations, Mercury forms small as the result of a series of imperfect accretion events. Gold points—A suite of 16 simulations reported in Chambers (2001) where the Mercury-forming region is modeled with a linearly increasing surface density of embryos and planetesimals between 0.3 and 0.7 au (hence the disk is "peaked"). Successful outcomes occur when Mercury forms directly from within this mass-depleted inner-disk component. Tan points—40 Mercury–Venus analog systems formed in 89 simulations in Lykawka & Ito (2019) considering similar inner-disk components to Chambers (2001). Blue points—The work of Izidoro et al. (2015) comprising 60 numerical simulations where the slope of the disk surface density profile was varied. Successful Mercury analogs are produced when an embryo is scattered out of the disk. Green points—A set of 125 simulations in Clement et al. (2019b) including a fragmentation prescription where the terrestrial planets form out of a narrow annulus of material (∼2 M_⊕) between 0.7 and 1.0 au as envisioned in Hansen (2009). Successful Mercury analogs are occasionally generated when an embryo is scattered out of the annulus. Magenta points—800 simulations from Clement et al. (2018) where the giant-planet instability (Tsiganis et al. 2005; Nesvorný & Morbidelli 2012) occurs within the first 10 Myr of the process of terrestrial-planet formation. Appropriate Mercury analogs are produced when the giant planets' resonant perturbations liberate an embryo from the planet-forming disk and strand it on a Mercury-like orbits. Pink points—A set of 600 simulations from Clement et al. (2019b) modeling the same early instability scenario and also incorporating the effects of collisional fragmentation. The initial conditions for the terrestrial disks in both instability batches plotted here (magenta and pink points) are identical to those of the classic disk models (5 M_⊕ of planet-forming material between 0.5 and 4.0 au; black and gray points).
Download figure:
Standard image High-resolution image

As the ultimate giant impact phase of terrestrial-planet formation is highly stochastic and unfolds over ∼10–100 Myr timescales (Touboul et al. 2007; Kleine et al. 2009), the computational cost of studying highly simplified systems of ≲1000 objects statistically can be substantial. Due to the small time step required to properly resolve the innermost (r < 0.5 au) region of planet-forming material, authors typically truncate the disk at ∼0.5 au in order to accelerate calculations (e.g., Fischer & Ciesla 2014; Izidoro et al. 2014; Jacobson & Morbidelli 2014; Clement et al. 2018). While a sharp interior edge in the distribution of solid material might seem loosely justified by pebble evaporation in the hot inner disk (Boley et al. 2014; we refer to the region of a ≲ 0.5 au as the inner disk throughout our manuscript), the preponderance of short-period planets in the exoplanet catalog suggests that planet formation in this locality is possible (Raymond et al. 2018). However, the majority of these so-called "hot" planets are, perhaps, best explained by gas-driven migration of giant-planet seeds (Izidoro et al. 2017; Bitsch et al. 2019; Lambrechts et al. 2019)—a process presumably halted in the solar system by the rapid formation of Jupiter's core (Kruijer et al. 2017)—or as the consequence of an intrinsic feature of the solar nebula such as a long-lived pressure maxima (Brasser & Mojzsis 2020). Nevertheless, the presumed truncation of the disk at 0.5 au is fairly arbitrary beyond the fact that it is necessary to prevent forming additional, ∼Earth-mass planets in the vicinity of Mercury's modern orbit (Chambers & Wetherill 1998). Indeed, authors investigating terrestrial-planet formation with various approaches often explicitly indicate their intention to ignore Mercury's formation (e.g., Raymond et al. 2009; Jacobson & Morbidelli 2014; Levison et al. 2015).

Early attempts to directly model Mercury's accretion considered an additional region of embryos and planetesimals with a linearly decreasing surface density profile between 0.7 and 0.3 au (Chambers 2001; O'Brien et al. 2006). To first order, these models succeeded in boosting the likelihood of forming a planet near Mercury's semimajor axis interior to three larger planets. However, the resulting Mercury analogs were too massive by around an order of magnitude (Mercury's modern mass, M_M = 0.055 M_⊕). Moreover, the dynamical separation between these planets and the systems' corresponding Venus analogs (in terms of their orbital period ratio: P_V/P_M) were more similar to the Earth–Venus spacing (P_E/P_V,ss = 1.6) than the modern Venus–Mercury period ratio (P_V/P_M,ss = 2.55; see Figure 1). Recent work by Lykawka & Ito (2017) and Lykawka & Ito (2019) investigating linear profiles between 0.2–0.5 au and 0.3–0.7 au, respectively, reached the same conclusion.

It is also possible that a primordial generation of close-in embryos or even planets existed in the vicinity of Mercury's modern orbit and were subsequently cleared out by some dynamical process. Volk & Gladman (2015) suggested that a system of tightly packed planets might have existed in the young solar system (similar to, for example, Kepler-11 or TRAPPIST-1: Lissauer et al. 2013; Gillon et al. 2017) and were eventually lost as a result of having formed in a quasi-stable configuration (Izidoro et al. 2017). While connecting the solar system to the Kepler catalog in this manner is compelling, it is unclear how the modern terrestrial planets might have survived such a scenario (see Raymond et al. 2016; Clement et al. 2019a; Bean et al. 2021; Lenz et al. 2020, for critiques of this hypothesis). Additionally, Raymond et al. (2016) proposed that the truncated initial conditions often assumed in the literature (i.e., the absence of planetesimals and embryos interior to 0.5 au) might be generated via planetesimal shepherding if Jupiter migrated through the terrestrial-forming region during the gas disk phase (as invoked in the so-called Grand Tack model: Walsh et al. 2011; Brasser et al. 2016; Walsh & Levison 2016). Alternatively, this classic setup might reflect the inefficiency of planetesimal accretion in Mercury's regime of the solar system due to the presence of a primordial Si snow line (Morbidelli et al. 2016). Thus, while it might be possible to roughly replicate the so-called "classic" initial conditions (and truncate the inner disk at ∼0.5 au: Chambers & Wetherill 1998; Raymond et al. 2004, 2009), the most plausible mechanisms for accomplishing this truncation are tied to the giant planets still unconstrained migration history (Pierens et al. 2014; Raymond et al. 2018; de Sousa et al. 2020) and the unknown properties of the protoplanetary disk. Moreover, the Mercury analogs generated from these classic initial conditions are quite scarce and are plagued by the same issues as those formed in studies assuming an additional inner-disk component (i.e., they are systematically overmassed and too close to Venus; Figure 1). In Clement et al. (2019a, henceforward Paper I), we demonstrated this with a systematic analysis of several different proposed terrestrial evolutionary schemes (Walsh et al. 2011; Raymond & Izidoro 2017a; Clement et al. 2018; see further discussion in Section 2). While it is possible at the ≲1% level for small (≲0.2M_⊕) objects with Mercury-like compositions to accrete inside of Venus' orbit, the solar system value of P_V/P_M is outside of the spectrum of simulation-generated outcomes.

Of increasing interest in the recent literature, Mercury's massive iron core (70%–80% of its total mass: Hauck et al. 2013; Nittler et al. 2017) and depleted bulk volatile content (note, however, that MESSENGER inferred an unexpectedly high surface volatile inventory for the planet: Nittler et al. 2011; Peplowski et al. 2011) have been interpreted to imply that its diminutive size is partially the result of an energetic, mantle-stripping collision (Benz et al. 1988, 2007; Asphaug & Reufer 2014). Indeed, hydrodynamical simulations of such a scenario consistently succeed at matching Mercury's mass and core mass fraction (CMF; Chau et al. 2018). Moreover, recent work demonstrated that this disrupted mantle material might be easily removed via interactions with the young Sun's intense solar wind (Spalding & Adams 2020). However, in Paper I we found that the preferred impact geometries and velocities for such an event (e.g., Asphaug & Reufer 2014; Chau et al. 2018; Jackson et al. 2018) are highly unlikely occurrences in conventional terrestrial-planet formation models, though it is possible that Mercury's CMF was further altered by accreting predominantly collisionally altered planetesimals in its late veneer (Hyodo et al. 2021). While Venus might seem to be a logical target for a Mercury-forming impact, in Paper I we were unable to reproduce the solar system value of P_V/P_M in simulations designed to model a proto-Mercury–Venus erosive impact. Moreover, the fact that Venus lacks a natural satellite and internally generated magnetic dynamo has been interpreted to suggest that its growth was not interrupted by such an energetic collision (Jacobson et al. 2017). An intriguing alternative to the impact hypothesis might be tidal stripping via repeated close encounters with Venus (Deng 2020). While such interactions do occur in N-body simulations of terrestrial-planet formation (Fang & Deng 2020), it is not clear how Mercury's orbit would be adequately recircularized or its semimajor axis driven away from Venus such that P_V/P_M,ss is matched after such a series of encounters. Thus, in spite of substantial effort and investigation, Mercury's enigmatic origin remains, arguably, the most significant outstanding problem in the terrestrial-planet formation literature (Figure 1).

This manuscript is part of a series of papers reexamining the Mercury problem from the ground up. The goal of our initial few investigations is to characterize functional models for the genesis of the precise Mercury–Venus system. In future work, we intend to more robustly develop these potentially viable scenarios within the larger context of planet formation and dynamical evolution in the solar system. In this paper, we revisit the scenario where Mercury accretes directly from an additional, mass-depleted, inner-disk component of material (Chambers 2001; O'Brien et al. 2006; Lykawka & Ito 2017). Specifically, we search for disk structures that are capable of generating the Venus–Mercury period and mass ratios with a first-order investigation of the applicable parameter space (disk surface density profile, total mass, extent, and embryo–planetesimal mass distribution). Our study is motivated by the work of Izidoro et al. (2014, 2015), where the authors varied the a ≳ 1.0 au region of the disk's mass and structure to find the parameters that best reproduced Mars (an extrapolation of the disk profile envisioned by Chambers & Cassen 2002).

2. Methods

We perform 210 N-body simulations (Table 1) of terrestrial-planet formation with the Mercury6 Hybrid integrator (Chambers 1999). All of our simulations are integrated for 200 Myr¹ and include algorithms designed to model the effects of collisional fragmentation (the fragmentation scheme is described in greater detail in Paper I; see also Leinhardt & Stewart 2012; Stewart & Leinhardt 2012; Chambers 2013; Clement et al. 2019b). Each computation includes Jupiter and Saturn (Levison & Agnor 2003) on their presumed pre-instability orbits (i.e., in a 3:2 MMR with a_J = 5.6 au: Tsiganis et al. 2005; Nesvorný & Morbidelli 2012; Deienno et al. 2017; Clement et al. 2021), use a time step of 1.6 days, remove objects that make perihelion passages within 0.02 au of the central body (Chambers 2001), consider objects ejected at a heliocentric distance of 15 au, and incorporate a prescription to account for the effects of relativity. As in Paper I, the minimum fragment mass (MFM) in our simulations is set to 0.0055 M_⊕ (10% of Mercury's modern mass: Wallace et al. 2017; Clement et al. 2019b).

Table 1. Summary of Initial Conditions for Our Various Simulation Sets

N_sim	a_in (au)	a_out (au)	M_tot (M_⊕)	α	R
10	0.2	0.6	0.05	0.5	4
10	0.2	0.6	0.05	1.0	4
10	0.2	0.6	0.05	1.5	4
10	0.2	0.6	0.05	2.0	4
10	0.2	0.6	0.1	0.5	4
10	0.2	0.6	0.1	1.0	4
10	0.2	0.6	0.1	1.5	4
10	0.2	0.6	0.1	2.0	4
10	0.2	0.6	0.5	0.5	4
10	0.2	0.6	0.5	1.0	4
10	0.2	0.6	0.5	1.5	4
10	0.2	0.6	0.5	2.0	4

10	0.35	0.75	0.1	0.5	1
10	0.35	0.75	0.25	0.5	1
10	0.35	0.75	0.5	0.5	1
10	0.35	0.75	0.1	0.5	4
10	0.35	0.75	0.25	0.5	4
10	0.35	0.75	0.5	0.5	4
10	0.35	0.75	0.1	0.5	8
10	0.35	0.75	0.25	0.5	8
10	0.35	0.75	0.5	0.5	8

Note. The columns are as follows: (1) the number of simulations in each set, (2) the inner-disk component's inner edge, (3) outer edge, (4) total mass, (5) surface density profile power law, and (6) ratio of total embryo to planetesimal mass (R). In all simulations, the inner disk comprises 20 embryos and 200 planetesimals, while the outer disk extends from the inner disk's outer edge to 1.0 (top set of simulations) or 1.1 au (bottom set), and contains 40 embryos and 400 planetesimals with M_tot = 2.0 M_⊕ and R = 4.

Download table as: ASCII Typeset image

2.1. Terrestrial Disk Structure

Our work follows the example of Izidoro et al. (2014, 2015) by experimenting with various total masses and surface density profiles for an additional, mass-depleted region of material. Izidoro et al. (2015) studied the Mars- and asteroid-belt-forming regions with a depleted outer disk. In our case, we apply the same logic to the Mercury problem by modeling a depleted inner-disk component possessing a surface density profile that increases with radial distance (Lykawka & Ito 2017):

$\begin{eqnarray}&&{\rm{\Sigma }}(r)=\left\{\begin{array}{l}{{\rm{\Sigma }}}_{\mathrm{in}}{\left(\displaystyle \frac{r}{1\mathrm{au}}\right)}^{+\alpha }:\quad \mathrm{inner}\,\mathrm{disk}\\ {{\rm{\Sigma }}}_{\mathrm{out}}{\left(\displaystyle \frac{r}{1\mathrm{au}}\right)}^{-3/2}:\quad \mathrm{outer}\,\mathrm{disk}\end{array}\right.\end{eqnarray} \tag{ 1 }$

where Σ_in and Σ_out are the disk's surface density at 1 au and are calibrated to achieve the desired total mass for the respective disk components. Thus, our work investigates a scenario where the terrestrial-forming disk's surface density profile did not possess a sharp inner edge (discussed further in Section 2.2, below). In all of our simulations, the inner-disk component is modeled with 20 equal-mass embryos and 200 equal-mass planetesimals. The precise initial masses for the respective particles depend on the parameters M_tot (the total inner-disk mass) and R (the ratio of total embryo to planetesimal mass: $R={M}_{\mathrm{tot},\mathrm{emb}}/{M}_{\mathrm{tot},\mathrm{pln}}$ ), which vary in our different simulations according to the values provided in Table 1. In general, our Mercury-forming embryos range in mass from 0.002 to 0.02 M_⊕, and the planetesimals in the region are assigned masses between 5.0 × 10⁻⁵ and 1.25 × 10⁻³ M_⊕. Each of the 220 inner-disk particles interacts gravitationally and can experience fragmenting collisions with all other objects in the simulation (note that particles initially smaller than the MFM cannot fragment until they grow larger than the MFM). Thus, the main difference between embryos and planetesimals in the inner disk is their mass. In 120 simulations, we investigate inner disks that extend from 0.2 to 0.6 au and possess an embryo–planetesimal ratio of R = 4 (based on the results of high-resolution simulations of embryo growth; e.g., Walsh & Levison 2019; Clement et al. 2020) and test a range of values for α and M_tot (Table 1). Our selection of boundary locations for our depleted inner-disk regions in these simulations is loosely based on Mercury's modern orbit and previous works considering inner-disk components spanning the region of 0.3–0.7 au that find Mercury analogs systematically form too close to Venus with such disk parameters (Chambers 2001; O'Brien et al. 2006; Lykawka & Ito 2017). Based on the most successful disk parameters from this first set of simulations (see further discussion in Section 3) we perform 90 additional simulations that test an inner-disk component extending from 0.35 to 0.75 au. We also vary the prescribed value of R in this follow-on suite of computations.

We structure our outer disks (Equation (1)) in a manner that maximizes the probability of forming Venus, Earth, and Mars analogs with the correct masses and semimajor axes. Thus, we take the so-called "annulus" (Agnor et al. 1999; Morishima et al. 2008; Hansen 2009) disk conditions with a truncated outer edge (1.0 au for simulations with 0.2–0.6 au inner disks and 1.1 au for those considering 0.35–0.75 au inner disks). In all of our simulations, the outer disk's total mass is 2.0 M_⊕ (R = 4) and is composed of 40 equal-mass embryos (M_emb = 0.04 M_⊕) and 400 equal-mass planetesimals (M_pln = 0.0025 M_⊕). Embryos in the outer disk interact gravitationally and can experience fragmentation events with all other simulation particles. Conversely, the 400 planetesimals cannot collide with or feel the gravitational effects of one another. When one of these planetesimals undergoes a fragmenting collision with an embryo or inner-disk planetesimal, the resulting new fragment particles are treated as embryos. Therefore, embryos and planetesimals in the outer disk differ in both their masses and numerical treatment within the integration.

Our initial conditions are roughly analogous to those of the low-mass asteroid belt model (Izidoro et al. 2014, 2015; Levison et al. 2015; Raymond & Izidoro 2017a, 2017b) or Grand Tack hypothesis (Walsh et al. 2011; Jacobson & Morbidelli 2014; O'Brien et al. 2014; Brasser et al. 2016; Walsh & Levison 2016), and do not consider the possibility that late giant-planet migration sculpted the outer terrestrial disk and Mars-forming regions (e.g., Lykawka & Ito 2013; Bromley & Kenyon 2017; Clement et al. 2018, 2019b, 2019c). For recent reviews of the various proposed terrestrial evolutionary scenarios, we direct the reader to Izidoro & Raymond (2018) and Raymond et al. (2018).

2.2. Iron Enrichment of the Mercury-forming Planetesimals

While the main motivation for our study is to reproduce the Venus–Mercury period and mass ratios in N-body terrestrial-planet formation simulations, it is worthwhile to discuss our scenario in the context of the various hypotheses that aim to explain how Mercury acquired its massive iron core. It is certainly possible that a random series of erosive impacts (Chambers 2013; Clement et al. 2019a) ensue in our mass-depleted inner disk such that Mercury finishes in a mantle-depleted state. However, our initial conditions and proposed scenario of direct, in situ formation are perhaps more consistent with ideas suggesting that the Mercury-forming planetesimals were already iron enriched prior to the giant impact phase (for a review of the differences between the various giant impact hypotheses and more "orderly" explanations for Mercury's origin in the context of MESSENGER's findings, see Ebel & Stewart 2017). MESSENGER determined that the inventories of less volatile, lithophile elements (Si, Ca, Al, and Mg) in Mercury's crust are not abnormal when compared to the other terrestrial planets and chondritic compositions (Weider et al. 2015). These results broadly refute ideas that silicates in the Mercury-forming region were evaporated by the intense solar activity or that the most refractory elements preferentially condense and accrete (e.g., Morgan & Anders 1980) as moderately volatile elements like Mg would be lost as well. However, there are three promising scenarios for iron enrichment of the material in the inner terrestrial disk that are still potentially viable.

2.2.1. Dynamical Fractionation

Wurm et al. (2013) argued that Mercury's peculiar composition might be explained by the photophoretic effect: a process through which small (e.g., ∼millimeter scale) particles in the gaseous disk can migrate rather substantially as the result of nonisotropic solar radiation. The magnitude of this additional force on particles is proportional to the regions' thermal gradient and the size of the particles themselves (Krauss & Wurm 2005). Thus, it is possible for particles to be size sorted through this process (Loesche et al. 2016). Because the midplane of protoplanetary disks is thought to be highly opaque, it is unclear whether this mechanism plays a significant role in altering the chemistry of the material near Mercury's modern orbit (e.g., Cuzzi et al. 2008). Nevertheless, several authors have attempted to blend dynamical and chemical models in this manner. Notably, Moriarty et al. (2014) used an analytical disk chemistry model that accounted for sequential elemental condensation during planetesimal formation in conjunction with N-body simulations to show that carbon-enriched planets can form at a variety of radial distances in systems with supersolar carbon abundances. However, the authors used the same inner-disk profile as Chambers (2001) and were thus unable to form reasonable Mercury analogs. Additionally, Pignatale et al. (2016) used a two-dimensional condensation model to show that two enstatite-rich regimes develop in the inner terrestrial disk, with the innermost region extending in as far as Mercury's modern semimajor axis. While more sophisticated chemical and dynamical models are still required to fully understand whether fractionation played a significant role in altering the composition of the Mercury-forming region, such a scenario is interesting in that it does not invoke a low-probability, violent dynamical event to explain Mercury's iron content.

2.2.2. Magnetic Aggregation

Ferromagnetic fluids of suspended magnetized particles that behave like dipoles tend to generate chain-like structures, and the viscosity of these chains is directly related to the strength of the applied magnetic field. As simulations (e.g., Dudorov & Khaibrakhmanov 2014) and observations (e.g., FU Ori: Donati et al. 2005) of protoplanetary disks indicate magnetic field strengths in the Mercury-forming regions as high as 10–100 mT, Kruss & Wurm (2018) proposed that magnetic interactions would preferentially generate large chains of iron-rich aggregates in the innermost regions of the disk (see Hubbard 2014 for a similar idea). In turn, these larger chains (Kruss & Wurm 2020) would presumably be more likely to become incorporated into the planetesimal precursors to Mercury formed via gravitational collapse (Youdin & Goodman 2005; Johansen et al. 2015; Simon et al. 2016). However, it is still unclear whether such a scenario is viable for the solar system.

2.2.3. A Carbon-rich Inner Disk

Ebel & Alexander (2011) showed that mixtures of presolar interplanetary dust particles in high-temperature, carbon-enriched, oxygen-depleted environments form condensates with Fe/Si ratios as high as half that of the value presumed for Mercury's bulk composition. In particular, the effect can be quite substantial if the dust concentrated at the disk midplane in the Mercury-forming region and is highly enriched in anhydrous chondritic interplanetary dust particles.² If isolated from the gas, these condensates could explain the iron enrichment of the precursors to the Mercury-forming planetesimals. However, it is still unknown whether graphite-rich silicates were present in the innermost regions of the terrestrial disk (e.g., Vander Kaaden & McCubbin 2015; Peplowski et al. 2016). While this process might explain the Fe/Si content of Mercury's core and mantle (i.e., it is bulk composition), some mechanical processes would still be required to reconcile the physical size of its core (i.e., its bulk structure and high CMF). Thus, while we scrutinize our fully formed Mercury analogs' final CMFs in the subsequent sections assuming a chondritic distribution of planet-forming material, we note that this is not a strict constraint if the depleted, inner-disk component planetesimals already possessed high Fe/Si ratios.

2.3. Success Criteria

In contrast to many classic studies of terrestrial-planet formation, we focus our investigation almost exclusively on our simulations' ability to generate Mercury–Venus analog systems. Thus, we do not present an analysis of our simulations' success in terms of traditional metrics related to Mars' mass, the planets' formation timescales, or the terrestrial-planet systems' dynamical excitation and orbital spacing (e.g., the angular momentum deficit and radial mass concentration statistics employed throughout the literature; Laskar 1997; Chambers 2001; Raymond et al. 2009) as we use initial conditions that are already well studied and verified to be reasonably successful in this manner (e.g., Hansen 2009; Walsh & Levison 2016; Clement et al. 2019b; Lykawka 2020). Instead, we scrutinize each system against four constraints (Table 2) that are designed to systematically compare simulated Mercury–Venus analogs with the actual pair of planets. Additionally, we introduce a fifth constraint that evaluates the general structure of the entire terrestrial system in order to ascertain whether successfully forming Mercury–Venus pairs and Venus–Earth–Mars systems are mutually exclusive results or not.

Table 2. Summary of Success Criteria for Our Simulations

Success Criterion	Parameter	Actual Value	Accepted Value
A	Mercury–Venus system	N/A	N_pln = 2 for a < 0.85 au; M_M < M_V
B	M_V/M_M	14.75	>5.0
C	P_V/P_M	2.55	>1.75; q_V > Q_M
D	CMF	∼0.7–0.8	>0.5
E	Terrestrial system	N/A	M_E,M_V > 0.6 M_⊕;
			0.05 < M_Ma < 0.3 M_⊕;
			a_E, a_V < 1.3 au; a_Ma > 1.3 au

Note. The columns are as follows: (1) the particular criterion's identifier, (2) the parameter scrutinized, (3) the solar system value, and (4) the value required for a system to be successful.

Download table as: ASCII Typeset image

2.3.1. Mercury–Venus Analog Systems

Criterion A separates systems that successfully form a Mercury–Venus duo of planets from those that do not. Specifically, A requires that a simulation finish with exactly two planets (defined here as any object with M > 0.01 M_⊕) possessing semimajor axes interior to 0.85 au (equidistant between Earth and Venus' semimajor axes in the modern solar system; see similar classification schemes in Clement et al. 2018; Lykawka & Ito 2019) with the innermost planet less massive than the outer one. Moreover, simulations that fail in this regard are discarded for the majority of the discussion and analysis sections of our manuscript. Thus, only systems that form exactly two objects with M > 0.01 M_⊕, a < 0.85 au, and the innermost planet less massive than the next planet are considered, regardless of the presence of additional leftover embryos and planetesimals less massive than 0.01 M_⊕ in the region.

2.3.2. Mass Ratio

Criterion B analyzes the mass ratio of the Mercury–Venus pair as defined by criterion A. While many authors scrutinize terrestrial-planet formation simulations by requiring that each respective planet analog finish below or above an established minimum or maximum mass tolerance (typically ≲0.1–0.3 M_⊕ for Mercury and Mars and ≳0.5–0.8 M_⊕ for Earth and Venus; Raymond et al. 2009; Clement et al. 2018; Lykawka & Ito 2019), this biases "successful" systems toward more massive Mercury analogs and diminutive Venus's. However, it is worth noting the obvious challenge of properly resolving Mercury's formation with initial embryo masses close to the modern mass of the planet. For instance, the embryos in our systems testing the most massive inner-disk components already possess 36% of Mercury's current mass at time zero. To avoid overconstraining our simulations, we require that criterion B satisfying simulations acquire Venus–Mercury mass ratios of at least 5.0 in order to permit a small number of embryo accretion events on Mercury. Indeed, the combination of criterion A and B adequately imposes an upper mass limit on our Mercury analogs as the most massive planet to satisfy both constraints in any of our simulations is 0.18 M_⊕.

2.3.3. Period Ratio

Criterion C considers the dynamical spacing of the Mercury–Venus system. In the solar system, Mercury's evolution is dynamically coupled to that of Venus by strong mutual nodal forcing perturbations (Nobili et al. 1989). Additionally, Mercury's proximity to the ν₅ secular resonance with Jupiter's perihelia precession drives chaotic orbital evolution in the inner solar system (Laskar 1997; Batygin et al. 2015) and has been shown to result in possible collisional trajectories between Venus and Mercury over the expected life of the solar system (Laskar & Gastineau 2009). While we leave the assembly of the inner solar system's secular architecture (e.g., Brasser et al. 2009) to future work, a major motivation of our study is the inability of contemporary terrestrial-planet formation simulations to replicate the modern value of P_V/P_M = 2.55 (Clement et al. 2019a; Lykawka & Ito 2019). Thus, criterion C requires that the final Venus–Mercury period ratio exceed 1.75. Additionally, in Paper I we found that simulations that successfully reproduced the inner two planets' masses (criterion B) often possessed Mercury and Venus analogs on crossing orbits. Therefore, we stipulate that criterion C satisfying runs finish with Venus' perihelion (q_V) beyond Mercury's aphelion (Q_M).

2.3.4. Core Mass Fraction

As discussed in Section 2.2, Mercury's CMF may not represent a strict constraint for our models if the planets' high mean density is not the result of mantle removal during a giant impact (Ebel & Stewart 2017). However, for consistency we still compute the final CMF of each Mercury analog formed in our simulations using the same procedure described in Paper I. Criterion D stipulates that successful Mercury analogs attain a final CMF of at least 0.5 (recall that Mercury's modern CMF is ∼0.7–0.8; Hauck et al. 2013; Nittler et al. 2017). Each embryo and planetesimal in our simulations is initialized as fully differentiated with a CMF of 0.3 (30% of the total mass in an iron core and 70% represented by a silicate-rich layer of mantle material; motivated by the Earth's CMF). Fragments are produced during the integration by dividing the mass of ejected material (determined by utilizing relationships from Leinhardt & Stewart 2012; Stewart & Leinhardt 2012) into a number of equal-mass fragments with masses greater than the MFM. The particles are then ejected in uniformly spaced directions within the collisional plane and with velocities ∼5% greater than the mutual two-body escape velocity. When such an event occurs, the fragments are first generated from the mantle of the projectile, followed by its core, the target object's mantle, and finally the target's core material. We experimented with alternative methodologies such as assigning each fragment a uniform mixture of the net ejected material. However, we determined that the particular choice of algorithm does not qualitatively alter the final results.

2.3.5. Score

In addition to evaluating each set of 10 simulations based on their rates of success for each individual criterion (A–D), we assign each batch of runs a numerical score that combines the respective success rates. To accomplish this, we linearly interpret between the solar system and simulation values of M_V/M_M, P_V/P_M, and CMF_M to convert each system's result into a number between 0.0 and 1.0. Period or mass ratios less than 1.0 (i.e., those that fail criterion A) and CMFs less than 0.3 receive scores of 0.0. Similarly, results exceeding the solar system value (i.e., M_V/M_M > 14.75, P_V/P_M > 2.55, and CMF >0.7) earn scores of 1.0. As an example, consider a Mercury–Venus system that finishes with M_V/M_M = 10.0, P_V/P_M = 1.5, and CMF_M = 0.3. The system's mass score would be (10.0–1.0)/(14.75–1.0) = 0.65, its orbital period score would be (1.5–1.0)/(2.55–1.0) = 0.32, and the corresponding CMF score would be 0.0. For simplicity, we then add the 30 scores for each batch of 10 simulations (Table 1) and multiply by 10/3 to report the final cumulative scores as a percentage.

2.3.6. Inner Solar System Analogs

Finally, we investigate the correlation between the ability of our systems to form accurate Mercury analogs and their success in terms of the other three terrestrial planets. While there are numerous ways in which we might constrain our Venus–Earth–Mars systems, we focus on the planets' respective masses and semimajor axes to avoid overconstraining our simulations. Moreover, our selected initial disk configuration of a narrow annulus is highly successful at replicating the low orbital eccentricities and inclinations of the terrestrial planets (Hansen 2009; Walsh & Levison 2016; Clement et al. 2019b; Lykawka 2020). Therefore, we adopt a similar classification scheme to the one proposed in Clement et al. (2018). Criterion E requires that a system form exactly two planets with m > 0.6 M_⊕ and a < 1.3 au, and exactly one Mars analog with 0.05 < m < 0.3 M_⊕ and a > 1.3 au (approximately equal to Mars' modern perihelion). This ensures that Earth and Venus are more massive than Mars by at least a factor of 2. Therefore, a system satisfying both criteria A and E necessarily contains analogs of all four terrestrial planets (albeit not necessarily an adequate Mercury analog in terms of criteria B and C). While this scheme obviously classifies some systems that are not exact solar system analogs as successful, we find it to be an adequate prescription for analyzing how our results for Mercury depend on the properties of the other three planets. As with our CMF calculation, we experimented with several alternative classification schemes and found that the overall trends are not particularly dependent on the precise methodology utilized.

3. Results

We tabulate the number of simulations from each of our different batches (Table 1) that satisfy our various success criteria (Table 2) in Table 3. Additionally, we denote the number of Venus–Earth–Mars analog systems (comprising a subset of 26% of all our simulations; criterion E) that meet each respective constraint in parentheses. In general Mercury–Venus analog systems (criterion A) are fairly pervasive in the majority of our simulation batches, and several sets of initial conditions yield 7 systems (of 10 total simulations) that are successful in this manner. However, consistent with previous studies that explored narrower regions of parameter space (Chambers 2001; O'Brien et al. 2006; Lykawka & Ito 2017), we find it quite difficult to reproduce the modern Venus–Mercury mass ratio (criterion B). We note that the final value of M_V/M_M is only moderately sensitive to the particular choice of initial inner-disk mass and radial slope (M_tot and α; see Figures 2 and 4). Interestingly, more massive disks tend to provide improved results when paired with a shallower inner-disk slope, while the combination of a lower value of M_tot and a larger α also yields reasonable success rates. This conclusion is supported by the cumulative scores (last column of Table 3) for the respective simulation batches. We elaborate further on these trends in Sections 3.1 and 3.2.

**Figure 2.** Cumulative distribution of Mercury–Venus system properties (top panel: mass ratio; bottom panel: period ratio) in criterion A–satisfying simulations. The different line colors (blue, gold, red, and gray) represent the four different total inner-disk masses (M_tot) employed in our simulations (Table 1). Analogs originating from different inner-disk radial ranges (0.2–0.6 au vs. 0.35–0.75 au) are plotted in solid and dashed lines, respectively. The vertical gray dashed lines represent the delimiting values for our success criteria B and C, while the solid gray line denotes the respective modern solar system values. The data plotted with a thick, gray, transparent line represent results from integrations presented in Paper I (Clement et al. 2019a, 2019b) that did not include a mass-depleted, inner-disk component (the "annulus" initial conditions of Hansen 2009).
Download figure:
Standard image High-resolution image

Table 3. Number of Simulations Satisfying Our Various Success Criteria (Table 2)

N_sim	a_in (au)	a_out (au)	M_tot (M_⊕)	α	R	A (A+E)	B (B+E)	C (C+E)	D (D+E)	Score
						Mercury–Venus Pair	M_V/M_M	P_V/P_M	CMF
10	0.2	0.6	0.05	0.5	4	3 (1)	1 (0)	3 (1)	1 (1)	17.9
10	0.2	0.6	0.05	1.0	4	5 (0)	0 (0)	2 (0)	2 (0)	25.2
10	0.2	0.6	0.05	1.5	4	6 (0)	2 (0)	4 (0)	2 (0)	35.8
10	0.2	0.6	0.05	2.0	4	6 (1)	1 (0)	4 (1)	2 (0)	34.9
10	0.2	0.6	0.1	0.5	4	5 (1)	0 (0)	4 (0)	1 (0)	27.2
10	0.2	0.6	0.1	1.0	4	5 (3)	1 (1)	4 (2)	0 (0)	29.1
10	0.2	0.6	0.1	1.5	4	5 (0)	1 (0)	3 (0)	4 (0)	32.5
10	0.2	0.6	0.1	2.0	4	3 (1)	0 (0)	3 (1)	2 (1)	15.2
10	0.2	0.6	0.5	0.5	4	0 (0)	⋯	⋯	⋯	0.0
10	0.2	0.6	0.5	1.0	4	0 (0)	⋯	⋯	⋯	0.0
10	0.2	0.6	0.5	1.5	4	0 (0)	⋯	⋯	⋯	0.0
10	0.2	0.6	0.5	2.0	4	4 (0)	0 (0)	4 (0)	2 (0)	25.7

10	0.35	0.75	0.1	0.5	1	7 (0)	2 (0)	5 (0)	2 (0)	39.8
10	0.35	0.75	0.25	0.5	1	7 (3)	3 (1)	4 (1)	1 (0)	48.0
10	0.35	0.75	0.5	0.5	1	4 (1)	1 (0)	3 (1)	2 (0)	24.7
10	0.35	0.75	0.1	0.5	4	6 (0)	0 (0)	1 (0)	4 (0)	39.1
10	0.35	0.75	0.25	0.5	4	3 (0)	1 (0)	3 (0)	1 (0)	20.0
10	0.35	0.75	0.5	0.5	4	6 (0)	0 (0)	6 (0)	2 (0)	41.6
10	0.35	0.75	0.1	0.5	8	4 (2)	1 (1)	3 (2)	3 (1)	23.7
10	0.35	0.75	0.25	0.5	8	6 (1)	1 (1)	3 (0)	3 (1)	37.5
10	0.35	0.75	0.5	0.5	8	4 (1)	0 (0)	2 (0)	3 (0)	24.1

Note. The first six columns are the same as those of Table 1. The subsequent four columns give the results for each of the success criteria: (A) the ability to form a Mercury–Venus pair with a < 0.85 au and M_M < M_V, (B) a mass ratio (M_V/M_M) greater than 5.0 for the inner two planets, (C) a period ratio (P_V/P_M) greater than 1.75, and (D) a CMF for Mercury of at least 0.5. The final column provides the simulation set's score (see Section 2.3.5 for a detailed explanation). The values in parentheses are the number of systems with successful Venus, Earth, and Mars analogs (criterion E) satisfying each respective success criteria. Note that the final four criteria can only be satisfied if criterion A is met.

Download table as: ASCII Typeset image

Our simulations indicate that criterion C (P_V/P_M) is easier to match than B, although the solar system value lies at the extreme of the distribution of possible outcomes for nearly all of our various disk structures. While adequate results are produced from nearly the full range of our tested parameter space, we find that satisfactory Venus–Mercury period ratios slightly correlate with low to moderate total planetesimal masses (R = 1 or 4). We explore these trends further in Section 3.3.

Finally, we note that our scenario provides a rather effective means of altering Mercury's final CMF, a consequence we attribute to the relatively high mutual collision velocities in the inner disk. Indeed, all but one criterion A–satisfying Mercury analog in our current investigation finish with CMFs that are altered from the initial value of 0.3 to some degree. However, we find no significant trends or commonalities between our various simulations that are successful with regard to criterion D and thus conclude that these systems succeed purely by happenstance (i.e., a high-CMF object coincidentally survives the simulation as the Mercury analog, see also Chambers 2013; Clement et al. 2019a). Thus, while none of our simulations simultaneously satisfy all five of our success metrics, we assess this to be the result of an overmultiplication of constraints, rather than strong mutual exclusivities between the individual criteria. However, several simulations perform adequately when scrutinized against four of our five constraints. We present a subset of these inner solar system analogs in Section 3.5.

3.1. Heavy versus Light Disks

We begin our analysis by scrutinizing the dependence of our results on the particular selection of inner-disk mass (M_tot). Figure 2 plots the cumulative fraction of final M_V/M_M (top panel) and P_V/P_M (bottom panel) values for the various values of M_tot tested in our simulations. It is important to note that we only plot and discuss systems that satisfy criterion A (for a measure of each system's performance that is negatively influenced by the failure of A, consult the score metric provided in Table 3). In this section (as well as Sections 3.2, 3.3, and 3.4) we focus on the relative success of these Mercury–Venus pairs formed in different disk structures in terms of our other three constraints. In Sections 3.5 and 3.6 we analyze the properties of the individual systems that are successful when measured against multiple metrics. As A stipulates that a system form exactly two planets with M_M < M_V interior to 0.85 au, we remind the reader that a nonnegligible fraction of Mercury-like planets are omitted from our analysis. For instance, a small number of systems form more than one Mercury analog interior to Venus. Conversely, some simulations yield two undermassed Venus analogs (note that neither class of simulation is considered in the majority of our analyses; specifically, only Figure 5 compiles data from all 210 simulations). However, the majority of the systems that fail criterion A are those that begin with M_tot = 0.5, a_in = 0.2, and a_out = 0.6. In these scenarios, the system begins with an excessive amount of mass concentrated interior to 0.85 au. Typically, such simulations finish with Earth, Venus, and an overly massive Mercury analog all with a < 0.85 au. This effect is slightly lessened in our set testing the steepest inner-disk mass profile (α = 2.0) as the more compact mass distribution tends to help Venus form closer to its modern semimajor axis. However, the Mercury analogs in this set still tend to be overmassed.

To prevent Mercury, Venus, and Earth from all growing in the vicinity of Venus' modern orbit, we performed an additional set of 90 simulations where the inner-disk component extends from 0.35 to 0.75 au, and the outer section stretches between 0.75 and 1.1 au. Figure 3 demonstrates how this change in annulus boundaries affects the range of initial semimajor axes of embryos and planetesimals incorporated into each final planet (referred to in the subsequent text as each planet's "feeding zone"; e.g., Kaib & Cowan 2015). While it is difficult to establish clear trends with only 10 simulations for each combination of varied initial conditions, it is clear from Table 3 that our most successful sets of simulations are overwhelmingly those that investigate inner disks stretching between 0.35 and 0.75 au.

**Figure 3.** Feeding zones of analog planets formed in our study depicted by the mass-weighted cumulative distribution of initial semimajor axes of all chondritic (i.e., not collisional fragments) objects incorporated into the final planets. Analog planets are defined as follows—Mercury: 0.01 < m < 0.2 M_⊕, a < 0.55 au; Venus: m > 0.6 M_⊕, a < 0.85 au; Earth: m > 0.6 M_⊕, 0.85 < a < 1.3 au; and Mars: 0.05 < m < 0.3 M_⊕, a > 1.3 au formed in our study. Simulations considering inner disks extending from 0.2 to 0.6 au (Table 1) are plotted in solid lines, and runs where the inner disk was distributed between 0.35 and 0.75 au are plotted in dashed lines. Mercury analogs formed in simulations testing the most extreme combinations of M_tot and α are plotted separately in gray and pink for comparison. In the same manner, Venus analogs from simulations testing the lowest value of M_tot are plotted in magenta.
Download figure:
Standard image High-resolution image

In spite of the tendency of our more massive disks to fail criterion A when the inner-disk boundaries are set to 0.2–0.6 au, we still assess moderate to larger values of M_tot to be generally more successful than smaller disk masses. While appropriate Venus–Mercury mass ratios are only produced in a small fraction of simulations, regardless of initial disk configuration (top panel of Figure 2), slightly more massive disks (M_tot ≃ 0.1) tend to more consistently yield larger values of P_V/P_M. This is a consequence of the total disk mass initially being less radially concentrated when the inner component is more massive. In the opposite scenario, the total disk mass is more condensed toward the outer disk, thus making it easier for Venus to accrete embryos and planetesimals in the Mercury-forming region (Figure 3). When this is the case, Mercury tends to remain dynamically coupled to Venus, in a sense forming as a by-product of Venus' formation. Contrarily, when the total inner-disk mass is larger, Mercury is more likely to form independently at a smaller heliocentric distance within the inner disk. This is also evident from the respective simulation sets' success rates for criteria A, B, and C. While the M_tot = 0.05 M_⊕ batch boasts success rates that are similar to those of the M_tot = 0.1 M_⊕ sets, on closer inspection we find that these "successful" simulations formed from less massive disks are better characterized as Venus–Earth analog systems inside of 0.85 au. Thus, these systems satisfy criterion A because an Earth analog and a less massive Venus analog both finish inside of 0.85 au. Increasing the value of M_tot to just 0.1 M_⊕ triples the total number of systems forming four terrestrial planets. For these reasons, our supplementary set of 90 simulations only consider inner-disk masses of 0.1, 0.25, and 0.5 M_⊕.

Analyzing our supplementary set of simulations' (0.35–0.75 au disks) success rates in terms of criteria B and C, it is clear that low-moderate values of M_tot tend to more consistently yield satisfactory results. Indeed, our most successful batch of simulations in terms of our score metric considers disks with M_tot = 0.25 M_⊕, α = 0.5, and R = 1. This is not particularly surprising, given the trends discussed above. An upper limit on the range of viable values of M_tot is necessary to prevent Mercury from growing too large, while excessively low total masses tend to increase the chances of Venus accreting the totality of material in the inner disk. Thus, future investigations modeling Mercury's formation from a mass-depleted interior component of the material should focus on moderate total masses (0.1–0.25 M_⊕) and truncate the outer disks' inner edge around Venus' modern semimajor axis.

3.2. Steep versus Shallow Disks

Figure 4 depicts the cumulative fraction of final M_V/M_M (top panel) and P_V/P_M (bottom panel) ratios for the various sets of our simulations testing different inner-disk surface density profiles (α; note that only our simulations investigating 0.2–0.6 au inner disks varied this parameter). In general, steeper slopes tend to be more successful at producing more realistic Venus–Mercury period ratios. Conversely, we were unable to decipher any general conclusive correlations between the presumed value of α and the planets' final mass ratios. Intriguingly, steeper slopes (α = 1.5 or 2.0) tend to be more successful for lower initial inner-disk masses (M_tot), while shallower slopes are advantageous in more massive disks. This is directly reflected in the scores for our M_tot = 0.05 M_⊕ (17.9, 25.2, 35.8, and 34.9 in order of increasing α) and 0.1 M_⊕ (27.2, 29.1, 32.5, and 15.2 in order of increasing α) batches investigating R = 4. Indeed, only two simulations considering inner-disk masses of 0.05 M_⊕ simultaneously satisfy criteria A, B, and C, and both are in the α = 1.5 batch. As discussed in Section 3.1, our simulations investigating lower values of M_tot systematically struggle to form Mercury analogs as Venus tends to accrete the majority of the inner-disk material. This broadening of Venus' feeding zone into the Mercury-forming region is depicted in Figure 3 by the difference between the solid magenta and gold lines. However, when the inner-disk surface density is sufficiently steep, the higher relative concentration of material in the Venus-forming region tends to restrict the planets' feeding zone. This, in turn, occasionally allows Mercury to form as a stranded embryo in the inner-disk component. The clear distinction between the solid gray (M_tot = 0.05; α = 2.0) and black (all simulations) lines in Figure 3 demonstrates Mercury's tendency to form in isolation in this manner in our integrations considering M_tot = 0.05 M_⊕ and α = 2.0. When the initial inner-disk mass is larger, Mercury forms more efficiently and in dynamical isolation from Venus with a shallower inner-disk surface density profile. This is depicted in Figure 3 by the difference between the solid pink (M_tot = 0.5;α = 0.5) and black (all simulations) lines. However, as these trends are rather weak, and we do not vary the initial value of α in our simulations investigating 0.35–0.75 au inner-disk components, we conclude that the inner disk's surface density profile only mildly affects Mercury's formation.

3.3. The Effect of the Total Planetesimal Mass (R)

We varied the total embryo–planetesimal mass ratio ( $R={M}_{\mathrm{tot},\mathrm{emb}}/{M}_{\mathrm{tot},\mathrm{pln}}$ ) in our additional batch of 90 simulations investigating inner disks with embryos and planetesimals distributed between 0.35 and 0.75 au. Naively, one might expect the additional dynamical friction (e.g., O'Brien et al. 2006; Raymond et al. 2006; Jacobson & Morbidelli 2014; Lykawka & Ito 2019) from a more massive swarm of planetesimals to restrict the feeding zones of each planet, but we find this effect to be almost negligible. However, we note that our runs investigating the most extreme bimodal mass distribution (R = 8) tend to struggle to satisfy criterion C (Table 3) and adequately reproduce the Venus–Mercury period ratio. We assess this to be a direct consequence of the systematically weaker dynamical interactions between the growing embryos and smaller planetesimals (though it is difficult to decouple this trend from the effects of varying M_tot in the final batch scores). In the top panel of Figure 5 we show the mean eccentricity of all embryos in the Mercury-forming region (a < 0.55 au) in our three batches of simulations investigating different values of R. While the difference between the R = 1 and 4 sets is minor, the embryos in the R = 8 runs possess substantially hotter eccentricity distributions for the majority of Mercury's growth timescale (the average time to reach 90% of its final mass: ≲50 Myr for all of our different simulation sets). As the embryos attain higher eccentricities in this manner, they interact more strongly with the growing Venus and are thus more likely to eventually be incorporated into Venus. Contrarily, when the embryo's eccentricity distribution is colder, the innermost embryos are dynamically sequestered from Venus in a manner that allows them to combine and accrete to form a Mercury analog in isolation. Through this process, the final planets tend to possess larger radial offsets from Venus, as demonstrated in the bottom panel of Figure 5.

**Figure 5.** Effect of the total embryo–planetesimal mass ratio ( $R={M}_{\mathrm{tot},\mathrm{emb}}/{M}_{\mathrm{tot},\mathrm{pln}}$ ) on the evolution of Mercury and Venus in the first 60 Myr of our simulations. The different line colors represent the three different values of R utilized in our study (Table 1). The top panel depicts the time evolution of the average eccentricity of all embryos in the Mercury-forming region (a < 0.55 au). The bottom panel shows the time evolution of the orbital spacing (Δa) of the two largest embryos with a < 0.85 au (ideally Mercury and Venus). As not all of our simulations form a Mercury–Venus pair (criterion A), the earlier times plotted in this figure average values from all 90 of our simulations investigating 0.35–0.75 au inner disks, while the latter times only depict the 47 simulations that meet criterion A (Table 3).
Download figure:
Standard image High-resolution image

It is worth pointing out that the more extreme bimodal distributions of material we find to be less successful are, perhaps, more consistent with studies of runaway growth (Kokubo & Ida 1996) during the gas disk phase. Indeed, recent high-resolution investigations of planetesimal collisional evolution increasingly find a strong radial dependence on the efficiency of runaway growth (Carter et al. 2015; Walsh & Levison 2019; Clement et al. 2020; Woo et al. 2021). Thus, we acknowledge that our R = 1 disk might not be an accurate representation of the potential disk conditions in the Mercury-forming region. However, as our R = 4 disks still yield reasonable results (Figure 5 and Table 3) we do not assess this to be a serious shortcoming of our model. Moreover, as the precise orbital and size distributions of the r < 0.5 au section of the solar system's terrestrial-forming disk remain unconstrained, we plan to further investigate the feasibility of our proposed initial conditions with high-resolution disk models in future work. It seems reasonable to speculate that a population of planetesimals and embryos drifting inward via Type I migration might be stranded in the Mercury-forming region during the disk's photoevaporation phase. Through this process, it might be possible for a low-mass inner component of the terrestrial disk to develop an overabundance of planetesimals relative to larger embryos.

3.4. Eroding Mercury's Mantle

While we argue that Mercury's final CMF might not represent a strict constraint for our models (Section 2.2), the collisional environment in our systems' inner disks provides an efficient mechanism for altering the planets' bulk composition. Indeed, ∼40% of the Mercury analogs in our criterion A–satisfying systems possess final CMFs in excess of 0.5 (criterion D). Figure 6 plots the cumulative distribution of Mercury-analog CMFs compared to a similar set of simulations taking the same annulus initial conditions (without an inner-disk component: Hansen 2009) from Paper I. It is not surprising that the Mercury-like planets generated in these simulations tend to have CMFs that concentrate around 0.3, as the majority of these planets form when an embryo is scattered from the annulus onto a relatively isolated orbit. Thus, the final analogs are far less compositionally processed than those produced in this manuscript from an inner disk of embryos and planetesimals. It is clear, however, that the iron contents of our contemporary simulations' Mercury analogs (independent of their success in terms of criterion B or C) are systematically enriched compared to that of their initial embryo precursors. Indeed, 100% of our final Mercury analogs experience at least one imperfect collision over the duration of our simulations. Curiously, while the ejected mantle material is predominantly removed from our systems via mergers with the Sun, it does so rather indirectly. As our terrestrial-forming disks extend over a more expansive radial range than those in Paper I, collisional fragments are far more likely to be transferred to different embryos rather than be absorbed by the target embryo or central body. Indeed, fragments produced in Paper I are ∼320% more likely to be removed from the simulation via collision with the Sun and ∼14% less likely to be accreted by an embryo other than the original target particle. However, embryos in our simulations considering an inner disk are nearly twice as likely to be removed via merger with the central star. Thus, while the ultimate fate of the mantle material removed from our Mercury analogs is the Sun, the material tends to be transferred to different embryos within the inner disk before eventually being expelled from the system.

**Figure 6.** Cumulative distribution of Mercury-analog CMFs in systems satisfying criterion A for all simulations presented in this work (black line), those satisfying both A and B (blue line), analogs successful in terms of both A and C (gold line), and comparison integrations from Paper I (Clement et al. 2019a, 2019b) that did not include a mass-depleted, inner-disk component (the "annulus" initial conditions of Hansen 2009). Recall that the fragmentation algorithm, MFM setting, and CMF calculation utilized here are identical to the methodology employed in Paper I.
Download figure:
Standard image High-resolution image

3.5. Preferred Disk Structure

Our simulations lead us to favor an inner-disk component of moderate total mass, shallow slope, and low to moderate R. Additionally, we find that mass depletion interior to 0.75 au (rather than 0.6 au) tends to boost the probability of forming a successful Mercury–Venus system. A sample of eight systems that simultaneously satisfy criteria A, B, and C is plotted in Figure 7. An example of a successful evolution of one of these systems (the one depicted in panel six of Figure 7) is plotted in Figure 8. The system begins with M_tot = 0.1 M_⊕, α = 1.0, and R = 4. Venus and Earth grow rapidly from seed embryos near the center of the disk at a_V,o = 0.69 and a_E,o = 0.87 au. By t = 2.8 Myr Venus attains half its ultimate mass. Similarly, Earth grows to 50% of its final size in 4.8 Myr. However, the evolution of the two larger planets subsequently bifurcate as Venus continues to rapidly accrete embryos and planetesimals from both the outer and inner disk in a manner such that it attains 80% of its eventual mass at t = 12 Myr. Conversely, Earth slowly grows to ∼65% of its ultimate size before experiencing a final giant impact with a 0.19 M_⊕ protoplanet originally seeded at 0.97 au at t = 25.4 Myr (we note that these divergent accretion histories are a reasonable example of the scenario proposed by Jacobson et al. 2017, which aims to explain Venus' lack of an internally generated magnetic field and natural satellite).

**Figure 7.** Semimajor axis/eccentricity plots of eight systems forming successful Mercury–Venus pairs (in terms of criteria A, B, and C) from our investigation. The size of each point corresponds to the mass of the particle. The solar system's terrestrial planets are plotted in the top panel for comparison.
Download figure:
Standard image High-resolution image

Curiously, the Mercury analog in Figure 8 originates when a mantle-only fragment at r = 0.45 au produced in a previous event collides with an embryo initially positioned at a = 0.56 au at t = 10.0 Myr. The collision yields three remnant particles with roughly equal masses of ∼0.005 M_⊕. Interestingly, one of these remnants is composed entirely of mantle material, another is totally derived from the core of the target embryo, and the final fragment (that goes on to become the Mercury analog) is constructed from a combination of the mantle-only projectile, the embryo's mantle, and part of the embryo's core. However, as the Mercury analog only possesses ∼6% of its ultimate mass at this point, the random assignment of material by our algorithm (Section 2.3.4) has little effect on the planets' ultimate composition. The fragment continues to slowly accrete planetesimals and embryos, thereby increasing its total mass to 0.042 M_⊕ (half its final value) by t = 20.0 Myr. Over the subsequent 7 Myr, the eventual Mercury analog experiences a series of 7 giant impacts with embryos in the inner disk which cumulatively boost its mass to 0.067 M_⊕. The planets' mass remains largely unchanged throughout the next ∼30 Myr aside from a few impacts with diminutive planetesimals and three minor fragmenting collisions, which erode a small amount of mantle material. At t = 86.6 Myr Mercury experiences its final giant impact with a 0.009 M_⊕ embryo. Thus, the Mercury analog in our example simulation has the longest accretion timescale of the four final terrestrial planets. Through this complex formative epoch, the planet attains a final mass of 0.083 M_⊕ (M_V/M_M = 13.5), semimajor axis of 0.33 au (P_V/P_M = 2.33), and a CMF of 0.32.

It is worth mentioning that this simulation's Mercury analog's orbit is not nearly as dynamically excited (e = 0.013; i = 2 fdg 0) as in the actual solar system (e_M = 0.21; i_M = 7 fdg 0). This is also the case for the vast majority of the other Mercury-like planets generated in our study. Figure 9 plots the orbits of our various criterion A–satisfying Mercury analogs in a/e, a/i, and a/m space, with the systems that are successful in terms of both A and B isolated in the right panel of the plot. It is clear from this figure and Figure 7 that even the Mercury analogs with masses most akin to that of the real planet tend to possess dynamically cold orbits. However, lower primordial eccentricity and inclination are likely advantageous as Mercury's orbit is easily excited during the giant-planet instability (Tsiganis et al. 2005; Nesvorný & Morbidelli 2012). Moreover, Roig et al. (2016) found that Mercury's precise orbit is reasonably explained by dynamical perturbations from Jupiter's stepwise semimajor axis evolution on an initially circular, coplanar Mercury (the so-called "jump"; Brasser et al. 2009; Nesvorný et al. 2013). Thus, we conclude that the tendency of our scenario to yield dynamically cold Mercury analogs is not particularly concerning as our simulations do not consider giant-planet migration. However, it is unclear whether an early dynamical instability (transpiring in conjunction with the giant impact phase of terrestrial-planet formation as proposed in Clement et al. 2018) would be compatible with our scenario and adequately excite Mercury's eccentricity and inclination.

3.6. Comparison with Paper I

The left panel of Figure 9 compares the orbits of the criterion A–satisfying Mercury analogs generated in this manuscript to a set presented in Paper I (note the fragmentation algorithm and MFM setting used in this study are the same as in Paper I). As this panel compiles many systems that are unsuccessful in terms of matching Mercury's mass, we isolate the 15 simulations (Table 3) that satisfy both criteria A and B in the right panel of Figure 9. In Paper I we noted that a narrow annulus (Hansen 2009) yielded larger values of P_V/P_M and M_V/M_M (though not as large as those of the real system) than the classic, extended disk model of Chambers & Wetherill (1998) perturbed by the giant-planet instability as envisioned in Clement et al. (2018). As we have reiterated throughout our manuscript, our new simulations' ability to form adequate Mercury–Venus analogs represents a marked improvement from the annulus models of Paper I. While the properties of the planets themselves are not particularly different (Figures 2, 4, and 9), our new models substantially improve the rate of generating Mercury–Venus systems of appropriate mass and orbital spacing. Specifically, 20% of our most successful set of initial conditions (0.35–0.75 au disk, M_tot = 0.1 M_⊕, α = 0.5, R = 1) and 10% of all our 0.35–0.75 au disk simulations simultaneously satisfy criteria A, B, and C. Conversely, only 5.6% of the annulus models from Paper I were successful in this regard.

In Figure 10 we revisit the parameter space of Venus–Mercury period and mass ratios depicted in Figure 1 by comparing the criterion A–satisfying Mercury analogs generated in this paper, with those from past dynamical modeling efforts (Chambers 2001; Izidoro et al. 2015; Clement et al. 2018, 2019b). It is clear from the distribution of simulation outcomes that, while our proposed scenario substantially improves the efficiency of forming Mercury–Venus systems in broad strokes, the precise configuration remains a low-probability outcome in each dynamical model. Encouragingly, the five or six most successful systems in Figure 10 are derived from our present study. Thus, our current manuscript serves as a proof of concept and motivation for future study of mass-depleted inner-disk components as a potentially viable means of consistently reproducing Mercury-like planets without invoking a low-probability giant impact. Follow-on investigation should further probe the successful parameter space of disk structures uncovered in this paper and work to incorporate and test the scenario's feasibility within the various proposed terrestrial evolutionary scenarios (e.g., Walsh et al. 2011; Levison et al. 2015; Bromley & Kenyon 2017; Raymond & Izidoro 2017a; Clement et al. 2018).

4. Conclusions

In this manuscript we presented a dynamical analysis of a scenario where Mercury forms directly within a mass-depleted, inner component of the terrestrial disk. Our study follows the mold of previous investigations into the formation of Mars and the asteroid belt that varied the mass and radial concentration of planet-forming material to determine the parameters that best reproduced the modern system (Chambers & Cassen 2002; Izidoro et al. 2014, 2015). In general, we find our scenario markedly improves upon previous efforts to form Mercury in terms of its ability to consistently generate analogs of the modern Mercury–Venus system. Additionally, we varied several disk parameters to determine which inner-disk structures most efficiently formed adequate Mercury analogs. Through this process, we conclude that our results are most sensitive to the total mass and radial extent of the inner disk, while the ratio of total embryo to planetesimal mass and the surface density profile (or slope) only mildly affect the final statistics. Specifically, we determine that the planets' precise orbits are best reproduced when the mass depletion extends to just outside of Venus' modern orbit. Furthermore, we find that moderate to larger total inner-disk masses (0.1–0.25 M_⊕) increase the likelihood of Mercury forming in dynamical isolation from Venus. When the inner-disk mass is too small, Mercury often survives as a stranded embryo and is typically too close to Venus. Finally, we note that our scenario provides a potent collisional environment that is capable of moderately modifying Mercury's composition through erosive collisions. In successful simulations, Mercury's mantle is slowly eroded throughout the growth process and the final planet possesses a massive iron core, in reasonable agreement with that of the actual planet. The ejected fragments of mantle material combine with other embryos in the simulation that disproportionately tend to be excited onto orbits that collide with the Sun.

In spite of all efforts made, Mercury's precise orbit and mass remain difficult to explain in our scenario. Specifically, the Venus–Mercury mass and orbital period ratios lie at the extreme of the distribution of outcomes generated in our simulations. Nevertheless, our results are promising in terms of their consistent ability to generate Mercury–Venus pairs, and several of our simulations generate remarkably accurate inner solar system analogs. Thus, future investigations must thoroughly investigate the parameter space of initial conditions determined to be plausible in this manuscript and simultaneously incorporate more detailed giant-planet evolution and migration models.

We thank Conel Alexander and an anonymous reviewer for insightful comments and input that greatly improved the manuscript. We are also grateful to André Izidoro, Seth Jacobson, and Patryk Lykawka for graciously and promptly sharing their data from Izidoro et al. (2015) Jacobson & Morbidelli (2014), and Lykawka & Ito (2019) for the production of Figures 1 and 10. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1548562. Specifically, it used the Comet system at the San Diego Supercomputing Center (SDSC) and the Bridges system, which is supported by NSF award number ACI-1445606, at the Pittsburgh Supercomputing Center (PSC; Nystrom et al. 2015). Additional computation for the work described in this paper was supported by Carnegie Science's Scientific Computing Committee for High-Performance Computing (hpc.carnegiescience.edu).

Dynamical Avenues for Mercury's Origin. II. In Situ Formation in the Inner Terrestrial Disk

Author e-mails

Author affiliations

ORCID iDs

Dates

DOI

Abstract

1. Introduction

2. Methods

2.1. Terrestrial Disk Structure

2.2. Iron Enrichment of the Mercury-forming Planetesimals

2.2.1. Dynamical Fractionation

2.2.2. Magnetic Aggregation

2.2.3. A Carbon-rich Inner Disk

2.3. Success Criteria

2.3.1. Mercury–Venus Analog Systems

2.3.2. Mass Ratio

2.3.3. Period Ratio

2.3.4. Core Mass Fraction

2.3.5. Score

2.3.6. Inner Solar System Analogs

3. Results

3.1. Heavy versus Light Disks

3.2. Steep versus Shallow Disks

3.3. The Effect of the Total Planetesimal Mass (R)

3.4. Eroding Mercury's Mantle

3.5. Preferred Disk Structure

3.6. Comparison with Paper I

4. Conclusions

Footnotes