Magma ascent in planetesimals: control by grain size

Tim Lichtenberg; Tobias Keller; Richard F. Katz; Gregor J. Golabek; Taras V. Gerya

doi:10.1016/j.epsl.2018.11.034

Magma ascent in planetesimals: control by grain size

Lichtenberg, Tim;Keller, Tobias;Katz, Richard F.;Golabek, Gregor J.;Gerya, Taras V. 2018-02-06 00:00:00 Rocky planetesimals in the early solar system melted internally and evolved chemically due to radiogenic heating from Al. Here we quantify the parametric controls on magma genesis and transport using a coupled petrological and ﬂuid mechanical model of reactive two-phase ﬂow. We ﬁnd the mean grain size of silicate minerals to be a key control on magma ascent. For grain sizes & 1 mm, melt segregation produces distinct radial structure and chemical stratiﬁcation. This stratiﬁcation is most pronounced for bodies formed at around 1 Myr after formation of Ca,Al-rich inclusions. These ﬁndings suggest a link between the time and orbital location of planetesimal formation and their subsequent structural and chemical evolution. According to our models, the evolution of partially molten planetesimal interiors falls into two categories. In the magma ocean scenario, the whole interior of a planetesimal experiences nearly complete melting, which would result in turbulent convection and core-mantle differentiation by the rainfall mechanism. In the magma sill scenario, segregating melts gradually deplete the deep interior of the radiogenic heat source. In this case, magma may form melt-rich layers beneath a cool and stable lid, while core formation would proceed by percolation. Our ﬁndings suggest that grain sizes prevalent during the internal heating stage governed magma ascent in planetesimals. Regardless of whether evolution progresses toward a magma ocean or magma sill structure, our models predict that temperature inversions due to rapid Al redistribution are limited to bodies formed earlier than 1 Myr after CAIs. We ﬁnd that if grain size was . 1 mm during peak internal melting, only elevated solid–melt density contrasts (such as found for the reducing conditions in enstatite chondrite compositions) would allow substantial melt segregation to occur. Keywords: Planetary formation, Planetesimals, Magma ocean, Melt migration, Chemical differentiation, Achondrites 1. Introduction The interior evolution of early solar system planetesimals has broad implications for the formation of both rocky planets At the time of planet formation, the inner solar system was and main-belt asteroid populations, the most immediate rem- populated by rocky planetesimals that seeded today’s ter- nants of the accretion process. Meteorites, broken-up pieces restrial planets through dynamical accretion of many smaller of asteroids fallen to Earth, are currently our only source bodies (Johansen et al., 2015), and whose internal evolu- of direct evidence from the early solar system. Therefore, tion was governed by radiogenic heating from short-lived Al our understanding of planetary growth and evolution is fun- (Hevey and Sanders, 2006). For large planetesimal radii and damentally limited by our ability to reconstruct the thermo- sufﬁcient Al incorporated upon formation, the released en- chemical evolution of planetesimals as evidenced by mete- ergy led to volatile degassing (Castillo-Rogez and Young, orites. Achondritic meteorites, which are thought to originate 2017; Monteux et al., 2018) and signiﬁcant silicate melting, from differentiated planetesimals, show a remarkable diver- surpassing the rheological transition from solid-state creep sity and likely originate from more than 50–100 parent bodies to disaggregation and melt-dominated deformation at melt (Wasson, 1990). However, spectral properties of asteroids fractions & 0:4–0:6 (Costa et al., 2009). In comparison trans do not match this diversity, as most known asteroids with an with solid or partially molten interiors, which lose heat by con- achondritic surface are interpreted to represent the debris of duction and/or laminar convection, disaggregation results in only a few parent bodies (Burbine et al., 2017). This appar- signiﬁcantly increased heat ﬂux by turbulent convection and ent lack of achondritic asteroids is at odds with the available rapid metal-silicate differentiation by the raining out of iron meteorite record. droplets (Stevenson, 1990). A possible solution to this conundrum is that internally dif- ferentiated planetesimals can retain their primitive, chondritic surfaces if magma remains conﬁned to the interior instead Corresponding author, now at: Atmospheric, Oceanic and Planetary of being erupted by volcanism (Elkins-Tanton et al., 2011; Physics, University of Oxford, Parks Rd, Oxford OX1 3PU, United Kingdom, email: tim.lichtenberg@physics.ox.ac.uk. Weiss and Elkins-Tanton, 2013). Some paleomagnetic stud- Published in Earth and Planetary Science Letters January 8, 2019 arXiv:1802.02157v2 [astro-ph.EP] 7 Jan 2019 ies on CV and CM meteorites suggest previous dynamo ac- mordial rock, with parts of the latter potentially lost to the tivity consistent with this hypothesis (Carporzen et al., 2011; core by percolation before the onset of major silicate melt- Cournede et al., 2015). Furthermore, Rosetta spacecraft ing (Yoshino et al., 2003; Bagdassarov et al., 2009; Ceran- data indicates a carbonaceous or enstatite chondrite sur- tola et al., 2015; Ghanbarzadeh et al., 2017). Higher oxygen face for 21 Lutetia and an average density of 3400 kg m fugacity may therefore result in more Fe-rich silicate melts (Sierks et al., 2011; Pätzold et al., 2011), higher than known with reduced (or even inverted) density contrast relative to chondrites and consistent with past compaction and partial the host rock. Lower oxygen fugacity, in contrast, may pro- melting beneath a primitive, chondritic crust (Weiss et al., duce iron-poor, buoyant melts that ascend rapidly. 2012; Neumann et al., 2013). In this study, we seek to establish regime boundaries that Based on the available evidence, most current models pro- separate primary evolution scenarios of early solar system pose a magma ocean scenario, where high melt fractions planetesimals by assessing the effects of melt segregation ( & ) dominated the thermal and chemical evolution on thermal evolution and chemical differentiation. We focus trans of planetesimal interiors. For the purposes of this study, we on the melting and transport of the major lithophile phases characterize the magma ocean scenario as a planetesimal in primitive bodies and investigate the potential for melt ac- exhibiting a fully molten interior of a well-mixed composition cumulation and heat source redistribution. We employ a and an adiabatic temperature proﬁle located beneath a thin computational model of coupled ﬂuid dynamics and thermo- ( 10 km), unmolten, chemically primitive lid with a linear chemical evolution that combines multi-component petrolog- conductive thermal proﬁle. Recent modeling studies investi- ical reactions with a two-phase magma transport model. We gating this scenario have relied either on thermal modeling quantify the leading controls on melt segregation in planetes- with parameterized melting (e.g., Hevey and Sanders, 2006; imals using theoretical considerations and numerical calcu- Elkins-Tanton et al., 2011), or on one-phase convection mod- lations of idealized planetesimal evolution. Our results show els (e.g., Golabek et al., 2014; Lichtenberg et al., 2016a, that both the magma ocean and magma sill scenarios are 2018) that capture the collective ﬂow and thermo-chemical realized within a relevant parameter space. We will focus evolution of partially molten rock or partly crystalline magma. our discussion on the latter case, where melt segregation is However, two-phase theory of partially molten systems most important. (e.g., McKenzie, 1984) suggests that silicate melts may buoyantly ascend relative to the ambient rock matrix. De- 2. Melt segregation scaling pending on the compositional and rheological properties of silicate minerals and their melts, this segregation may have To gain a leading-order understanding of silicate melt as- delayed, or even precluded, the generation of a magma cent in Al-heated planetesimals, we ﬁrst consider the char- ocean structure. Ascending melts may instead have formed acteristic time scales of melt transport in partially molten melt-rich layers beneath the primitive lid, hereafter referred bodies. In two-phase theory (McKenzie, 1984), the charac- to as magma sills (Wilson and Keil, 2017). For the purposes teristic length scale of melt migration by porous ﬂow relative of this study, we deﬁne this magma sill scenario as a body to a permeable rock matrix is given by the compaction length with radial heterogeneities of melt fraction that differ from the ﬁducial magma ocean case. This scenario implies a poten- 0 0 = ; (1) 26 c tially signiﬁcant redistribution of Al, which is a moderately 0 0 incompatible element and preferentially partitions into sili- cate melts. The transfer of the major heat source into shallow with characteristic rock viscosity , melt viscosity , melt 0 0 magma sills might then result in a transient, inverted temper- fraction , and rock permeability ature proﬁle with a thermal history distinct from the magma 2 n 0 0 ocean scenario. To date, only few studies have quantitatively k = ; (2) b (1 ) investigated this effect. These were either based on melt 0 transport models with parameterized melt ascent velocities where a is the characteristic grain size, b a geometric factor, (Moskovitz and Gaidos, 2011; Wilson and Keil, 2012; Man- and m, n powerlaw exponents. The characteristic velocity of dler and Elkins-Tanton, 2013; Neumann et al., 2013, 2014, segregating melts is 2018), or focused on metal-silicate separation (Šrámek et al., k g 2012; Ghanbarzadeh et al., 2017). 0 0 0 w = ; (3) The efﬁciency of melt transport in planetesimals depends 0 0 on various parameters. The presence of primordial volatiles with the solid-melt density contrast, and g the surface 0 0 favors rapid segregation by increasing the buoyancy and low- gravity. In primordial, homogeneous planetesimals, gravity ering the viscosity of magmas. However, if volatiles are ex- increases with increasing distance r from the center, solved before the onset of silicate melting, Fu and Elkins- Tanton (2014) argue that the segregation rate of dry melt is g(r) = g r=R ; (4) 0 P mostly controlled by the oxygen fugacity and the degree of melting. The oxygen fugacity determines (or is determined where R is the total planetesimal radius. The ﬁrst silicate by) the relative abundance of FeO and Fe-FeS in the pri- melts in sufﬁciently large planetesimals form in an adiabatic 2 zone stretching from the center, where gravity is negligible, to below the upper conductive lid (Hevey and Sanders, 2006; Lichtenberg et al., 2016a). As melting progresses, the per- meability increases and melts in shallower regions of the planetesimal, where gravity is highest, begin to segregate from the residual rock. Melt segregation can alter the chem- ical and thermal structure of Al-heated planetesimals, be- cause early-formed melts are preferentially enriched in in- compatible elements including Al. We deﬁne the characteristic melt segregation time scale, , as segr = R =w : (5) segr P 0 To achieve a substantial redistribution of heat-producing Al, the rate of melt transport must exceed the rate of melt gener- ation. We thus deﬁne the heating time scale, , of a plan- heat etesimal at a given time t after Ca,Al-rich inclusions (CAIs) as = c T =H26 (t ) ; (6) heat p 0 0 Al Figure 1: Scaling analysis of melt segregation propensity, with melt seg- with the speciﬁc heat capacity of silicates, c = 1100 J 1 1 regation number R = log ( / ). (A) At t = 0 Myr (CAI formation) seg 10 heat mt kg K , the temperature difference between accretion and and with rising melt fraction , the migration velocity increases, and so the solidus temperature, T 1100 K, and the decay power of system is more likely to become segregated. At around & 0:4–0:6, the Al per unit mass, magma ocean regime is reached and the system would be dominated by turbulent convection. (B) For ﬁxed melt fraction = 0:4 and later times (= " # Al E26 weaker radiogenic heating) the melt segregation number rises further. Al t = 0 26 Al H26 (t ) = f e : (7) 0 Al Al Al Al Here, f is the chondritic abundance of aluminum per unit compared to the internal heating by Al in planetesimals Al h i Al 5 mass (Lodders, 2003), = 5:25 10 is the canonical (e.g., Hevey and Sanders, 2006; Lichtenberg et al., 2016a, Al 26 27 2018). However, it crucially depends on the dynamic evolu- ratio of Al to Al at CAI formation (Kita et al., 2013), E = Al tion of the melt fraction, which is controlled by the ﬂuid me- 3:12 MeV = 5 10 J is the decay energy, and = Al chanics of melt transport in a deforming rock matrix, and the t =ln(2) = 1:03 Myr is the mean lifetime of Al. 1=2; Al thermo-chemical evolution of the body. In particular, Figure Using these characteristic scales, we deﬁne the non- 1 highlights the importance of considering the evolution of dimensional melt segregation number, internal heating and melt segregation over time. The scal- ! ! k g c T heat 0 0 0 p 0 0 ing analysis does not yet capture any time-dependent ef- R = log = log ; (8) seg 10 10 R H (t ) fects of interest here, which include the potential accumu- segr P 0 Al 0 lation of magma beneath a primitive lid and the redistribution which quantiﬁes the propensity of a planetesimal to undergo of the heat source by transport of melt enriched in Al. In substantial melt segregation during the internal heating by order to assess these dynamic processes we require a time- Al as a function of the key model parameters. To antic- dependent evolution model, which we introduce in the next ipate the expected melt segregation regimes in planetesi- section. mals, we calculate R for a reasonable range of melt frac- seg tions ( 2 [0; 0:4]) below the rheological transition ( 0 trans 0:4–0:6, Costa et al., 2009), formation times (t 2 [0; 4] 3. Method 5 2 Myr), grain sizes (a 2 [10 ; 10 ] m) and melt-rock density contrasts ( 2 [120; 1200] kg m ). Figure 1 shows that 3.1. Melting and heat source partitioning a growing melt fraction increases the melt segregation num- ber through its effect on permeability. Larger grain sizes and Studies of primitive meteorites (McSween et al., 1991; higher density contrasts also signiﬁcantly enhance segrega- Dunn et al., 2010) and equilibrium condensation sequence tion, but the effect of the latter is less pronounced than the calculations (Ebel and Grossman, 2000) suggest that the former. In Figure 1b, holding melt fraction ﬁxed, the heat- main rock-forming mineral phases in solar system planetesi- ing rate decreases with later formation times, which again mals were olivine (dominantly forsterite, Mg SiO ), pyroxene 2 4 serves to favor melt segregation relative to melt production. (dominantly enstatite, MgSiO ), and feldspar (dominantly From this scaling analysis, we conclude that melt segre- anorthite, CaAl Si O ). Ignoring minor contributions from 2 2 8 gation can in principle occur on a time scale that is relevant CAIs and trace minerals, feldspar represents the major host 3 Table 1: Scaling quantities, deﬁnitions and parameter values introduced in the scaling analysis, R_DMC and two-phase ﬂow models. Varying model parameters are named in the text and ﬁgures. Parameters not listed here are as given in Table 1 in Keller and Katz (2016). Parameter Symbol Unit Value Geometric factor b non-dim. 100 Melt fraction exponent n non-dim. 3 Solid fraction exponent m non-dim. 2 Melt shear viscosity Pa s 1 1 5 Thermal expansivity K 3 10 1 1 Speciﬁc heat capacity c J kg K 1100 2 1 6 Thermal diffusivity m s 1:14 10 olv olv initial mass fraction c¯ wt % 50 pxn pxn initial mass fraction c¯ wt % 35 fsp fsp initial mass fraction c¯ wt % 15 olv olv melting point T K 2050 m;0 pxn pxn melting point T K 1500 m;0 fsp fsp melting point T K 1350 m;0 Entropy gain of fusion S J K 320 olv pxn fsp 1 1 Curvature coefﬁcients r , r , r J kg K 50, 20, 10 olv pxn fsp 1 Linear P-coefﬁcients B , B , B K GPa 60, 100, 120 Rock density kg m 3200 Reference rock viscosity Pa s 10 Shear viscosity cut-off Pa s 10 min Compaction viscosity cut-off Pa s 10 min Initial temperature T K 290 init Surface temperature T K 290 space Planetesimal radius R km 60 5 2 Grain size a m [10 , 10 ] Formation time t Myr [0, 4] form Melt-rock density contrast kg m [120, 1200] 26 i phase of Al in rocky planetesimals. In addition to the tim- Curvature coefﬁcients r adjust the temperature dependence ing of accretion and size of a planetesimal (e.g., Lichten- of partition coefﬁcients. The pressure-dependent pure- berg et al., 2016a, 2018), magma genesis depends on the component melting points are parameterized as relative abundance of these phases and the concentration i i i T (P) = T + B P ; (11) m m;0 of volatiles. However, to avoid further complexity, we will consider only dry melting here, which is justiﬁed if volatile with linear slopes B . degassing during accretion is efﬁcient (e.g., Monteux et al., At a given volume-averaged bulk composition (lever rule) 2018). We therefore formulate a model for melting and melt- i i i c¯ = (1 )c + c ; (12) s ` rock partitioning of these major mineral phases. eq We employ the R_DMC method of Keller and Katz (2016) we numerically determine the equilibrium melt fraction to calculate an idealized thermodynamic equilibrium at given that satisﬁes the partitioning coefﬁcients K by combining the temperature, pressure, and bulk composition, and linear ki- lever rules with the unity sum constraint on all components netic reaction rates for a multi-component system. We deﬁne n n X i X i c¯ c¯ three compositional pseudo-components and their mass- = : (13) eq i eq eq eq i i i =K + (1 ) + (1 )K concentrations in the solid (rock), c , and liquid (melt), c . s i=1 i=1 These capture the leading-order behavior of classes of min- eq From , we then calculate the equilibrium phase composi- erals grouped by similar melting points and partitioning be- tions for solid and melt, havior: olv (olivine-like, i = 1), pxn (pyroxene-like, i = 2) c¯ and fsp (feldspar-like, i = 3). The mass fraction of all three i;eq c = ; (14a) eq i eq components must sum to unity in both phases. =K + (1 ) Using a simpliﬁed form of ideal solid solution theory c¯ i;eq c = : (14b) (Rudge et al., 2011), we determine the component partition ` eq eq i + (1 )K coefﬁcients at given P; T -conditions, Dynamic changes in pressure, temperature or bulk com- " !# i;eq c L 1 1 position over time create disequilibrium. The mass transfer K = = exp ; (9) i;eq i i of component i from solid to liquid that drives the system r T T (P) c m back towards equilibrium is assumed to occur at linear ki- i;eq i;eq which are the ratios of solid, c , to liquid, c , compo- netic rates, nent concentrations at equilibrium. The latent heat of pure- i;eq i eq i = R c c ; (15) component fusion, with a rate factor R = = , which restores equilibrium over i i 0 L = S T ; (10) m;0 a time scale at reference density . The sum of all com- is given by the entropy gain of fusion, S , and the pure- ponent reaction rates is the total melting rate = . All component melting temperatures at zero pressure, T . parameter values are given in Table 1. m;0 4 pxn fsp of melting increases with temperature, more pxn is dissolved 90 ΔT = T - T = 150 K m m,0 m,0 into the melt. Finally, in absence of melt migration, the melt ΔT +/- 50 K composition would converge to the bulk composition as the ΔT +/- 100 K system approaches complete melting. At our chosen refer- pxn ence calibration T = 1500 K, the melt initially comprises m;0 80% fsp, but becomes relatively enriched in pxn by the time the melt fraction reaches . trans 3.2. Two-phase, multi-component ﬂuid model As partial melts segregate from their residual, the interior becomes gradually depleted of fsp and thus its heat source Al. To model this dynamic process, we couple the multi- fsp B component melting model above to the two-phase reactive transport model of Keller and Katz (2016). The ﬂuid mechan- pxn ics part of the model is based on McKenzie (1984). Here, we give a brief summary of the governing equations and consti- tutive relations. The physical model is derived from statements for the con- olv servation of phase and component mass, phase momentum, and total energy. We consider the model in a Cartesian co- ordinate system with gravity pointing down the vertical co- ordinate, g = gzˆ. The governing equations are formulated olv in the geodynamic limit (liquid viscosity rock viscosity ), using the extended Boussinesq approximation (densities = = taken equal and constant except when multi- s ` 0 plying gravity). The ﬂuid mechanics governing equations are pxn " # fsp r rv + (rv ) Ir v rP + g = rP ; s s s ` C 1500 1600 1700 1800 (16a) Temperature [K] r v = r v ; s D Figure 2: Temperature dependence of melt fraction (A), melt (B) and rock (16b) (C) composition, with varying melting point difference between fsp and pxn, pxn fsp T = T - T , which changes the partitioning coefﬁcient of fsp and m m v = (v v ) = (rP + g) ; m;0 D ` s ` the composition of the earliest melts. If the fsp melting point is close to the one of pxn, the initial melt composition is close to a 50–50 mixture. (16c) For higher melting point difference, the ﬁrst melts are dominated by fsp, P = (1 )(P P ) = r v : and thus the heating component ( Al) of the system preferentially follows C s ` s the dynamics of the earliest melts. When the absolute temperature of the (16d) system further rises and approaches the olv melting point, the composition converges towards the initial pure solid setting. They are posed in four independent variables, the dynamic pressures, P , P , and velocities v , v , of the solid and liq- s ` s ` uid phases. Two dependent variables, the Darcy segregation ﬂux, v , and the compaction pressure, P , express the me- We ﬁrst consider the aluminum partitioning behavior in a D C chanical interactions between the phases. If these vanish, closed, isobaric (P = 0) system under increasing tempera- the equations become identical with the Stokes system. As- ture. For consistency with earlier work (Tkalcec et al., 2013; suming a diffusion creep rheology with melt-weakening, the Golabek et al., 2014; Lichtenberg et al., 2016a; Monteux shear viscosity of the rock matrix is given by et al., 2018; Lichtenberg et al., 2018; Hunt et al., 2018), we fsp choose the lowest component melting point, T = 1350 K, m;0 E + PV a a to conform with previous estimates of the silicate solidus. To = A a exp ; (17) RT test different partitioning behaviors, we vary the relative tem- perature difference between the melting points for the Al- with prefactor A , activation energy, E , and activation vol- 0 a pxn fsp hosting fsp and the pxn components, T = T T , with m ume V as in Hirth and Kohlstedt (2003) and Mei et al. m;0 m;0 a pxn T 2 [1400, 1600] K. The resulting rock and melt composi- (2002). R is the universal gas constant, and 30 the melt- m;0 tions in the ternary system olv-pxn-fsp are shown in Figure weakening factor. Permeability is set by the Kozeny-Carman 2. The larger T between fsp and pxn, the more incipient relation (Eq. 2), with powerlaw exponents n = 3, m = 2 for melt will be enriched in fsp, and thus Al. As the degree the melt and solid fractions, respectively. The compaction Rock composition [wt %] Melt composition [wt %] Melt fraction [wt %] viscosity is set to = r =, with r 1 the shear to com- interior evolution of planetesimals of R & 50 km as these are paction viscosity ratio. dominated by a relatively large adiabatic interior and a thin To these equations we add thermo-chemical evolution (. 10 km) conductive lid, whereas the relative dimensions of equations for temperature, T (assuming thermal equilibrium these domains vary signiﬁcantly for planetesimal radii . 50 between phases), melt fraction, , and pseudo-component km (cf. Castillo-Rogez and Young, 2017; Lichtenberg et al., i i concentrations in the solid and liquid phases, c and c : 2018; Hunt et al., 2018). The computational domain includes 500 grid cells for a spatial resolution of 120 m. The sur- X i i DT L H (t) T face boundary is T = 290 K, similar to the temperature Al 2 space = r T + + + gw ; in the protoplanetary disk inside of the water snowline, while Dt c c c c 0 p 0 p p p i=1 the center boundary is insulating, @T=@z j = 0. We as- z=0 (18a) sume the accreted body is initially at ambient temperature, T = T . As noted above, gravity decreases linearly from = (1 )r v + = ; (18b) init space s 0 Dt the surface gravity down to zero at the center. i i i D c c ` ` Magma and rock composition are modeled in the three- = ; (18c) Dt 0 component compositional space of olv, pxn and fsp (Section i i i 3.1). We use component melting points as in Figure 2 (solid D c c s s = : (18d) lines). The solid-melt density contrast is varied as 2 [200, Dt (1 ) 1200] kg m between runs to reﬂect FeO content and thus Temperature evolves due to advection, thermal diffusion, density of the melt reﬂecting the planetesimal’s oxygen fu- latent heat exchange of reactions, heating by viscous dis- gacity (Fu and Elkins-Tanton, 2014; Wilson and Keil, 2017). sipation and internal heating H26 (t), and adiabatic de- Al Grain size a , which controls both the permeability and rock compression. Melt fraction evolves due to advection, com- viscosity, is held constant during calculations; we consider paction and reactions, and composition by advection and re- values a 2 [10 m, 1 cm], from chondrite matrix-like dust action. The material derivative of the two-phase mixture is to pallasite-like crystal sizes (cf. Figure 1). Heating is in- ¯ 26 D=Dt = (1 )D =Dt + D =Dt, with D =Dt = @=@t + v r, s ` s s duced solely by Al in the fsp component, whose redistribu- and D =Dt = @=@t + v r. is the thermal diffusivity, ` ` tion hence affects the local heating rate. The initial heating 7 1 w = [(1 )v + v ] zˆ the vertical bulk speed, and the z s ` rate is varied from H26 (0) 2 [1:52; 0:19] 10 W kg , re- Al thermal expansivity. ﬂecting planetesimal formation times t 2 [0; 3] t = form 1=2; Al The governing equations for the ﬂuid mechanics [0; 2:2] Myr after CAIs. To limit model complexity, we ignore (Eqs. 16a–16d) and thermo-chemistry (Eqs. 18a–18d) are the potential heat contribution from Fe (see, e.g., Tang and discretized using the ﬁnite difference method on a rectangu- Dauphas, 2015; Lichtenberg et al., 2016b). lar, staggered grid and solved by two coupled Newton-Krylov solvers. The simulation software uses parallel data struc- 4. Results tures and solvers provided by PETSC (Balay et al., 1997; Katz et al., 2007). Nonlinearities between the ﬂuid mechan- 4.1. Parameter study ics and thermo-chemical sub-problems are resolved using We ﬁnd that model results across the tested parameter a Picard ﬁxed-point iterative scheme. During every itera- range fall into three qualitative categories. Some models tion, Equation 13 is solved using Newton’s method to up- show no substantial melting or segregation; we will not fur- date the equilibrium melt fraction. The adopted model is ther discuss these undifferentiated models here. The time strictly valid only for melt transport by porous ﬂow below the evolution of the latter two categories shown in Figure 3 is disaggregation threshold. However, we cannot avoid mod- generally the same: rapid initial heating leads to substan- els producing regions with higher melt content. To ensure tial melt production, followed by some degree of segregation, that the equations do not produce numerically unstable so- before ending with slow cooling from the top down. One cat- lutions in these regions, we apply lower viscosity cut-offs to egory of models, which we identify as the magma ocean end- the shear viscosity ( = + ) and compaction viscosity num min member (red), evolve to where most of the interior is above ( = + =(1 )). The effect of this regularization is to num min the disaggregation threshold, whereas the other, the magma dampen the segregation velocity and compaction pressure at sill end-member (blue), result in a melt-depleted interior be- elevated . As a result, our model will produce stable solu- neath melt-rich sills (Fig. 3A–E). The latter clearly shows a tions above the rheological transition, but will underestimate thermal inversion (Fig. 3F–J) related to fsp-enrichment in the chemical mixing and heat transport by rapid convection in magma sills beneath the lid (Fig. 3K–O). the crystal-bearing suspensions that characterize this limit. The scaling analysis above predicts that grain size, a , and density contrast, , are pertinent controls on melt segrega- 3.3. Model setup and parameter space tion. Figure 4 shows the results of a detailed study of that We model magma genesis and transport along a 1D parameter space. To quantitatively analyze the results, we Cartesian column from the center to the surface of an initially introduce three metrics measuring the degree of homogeneous and isothermal body of 60 km radius. Plan- etesimals of this size are qualitatively representative of the melt segregation: = – , max ctr 6 1.5 Myr 2.0 Myr 4.0 Myr 7.0 Myr 10.0 Myr 0.9 0.7 0.5 0.3 0.1 A A B B C C D D E E 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 [wt] 0.9 0.7 0.5 0.3 0.1 F F G G H H I I J J 1450 1700 1950 1450 1700 1950 1450 1700 1950 1450 1700 1950 1450 1700 1950 [K] 0.9 0.7 0.5 0.3 0.1 K K L L M M N N O O 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 [wt] Figure 3: Time evolution of melt fraction (A–E), temperature (F–J) and fsp bulk composition (K–O) for two end-member models with radii R = 60 km and 3 4 formation time t = 1.25 t = 0.90 Myr after CAIs. Red lines show a magma ocean model with = 200 kg m and a = 10 m, blue lines a form 26 0 1=2; Al 3 3 magma sill model with = 800 kg m and a = 10 m. Upon progressive heating, the magma sill model builds up melt accumulates below the cold upper lid, depleting the center of the planetesimal of silicate melts. High melt fractions > := 0.5 (yellow areas) are only reached in the sub-lid sills. fsp crit enrichment in the sill structure leads to a temperature inversion of 400 K at peak melting. The magma ocean model, in contrast, shows a near-isothermal internal temperature and thus constant melt fraction structure in the interior. The fsp component shows notable deviations from the initial bulk composition only after t 2 Myr, when most of the Al is already decayed. A video showing the time evolution of the major thermo-chemical parameters and composition is linked to in the Supplementary Materials. temperature inversion: T = T – T , 4.2. Silicate differentiation max ctr fsp fsp fsp differentiation: c = c – c . fsp max ctr Here, () denotes the maximum value in the computational As the planetesimals cool and crystallize, magma sill end- max domain and () denotes the value at the base of the do- member cases exhibit a silicate differentiation trend towards ctr main, i.e., the planetesimal center. With these metrics, we compositional layering. In Figure 5, we compare a repre- quantify the deviation from an interior structure with near- sentative magma ocean with a magma sill model outcome. constant melt fraction, temperature, and bulk composition, Magma ocean models evolve towards a uniformly molten as it would be expected in the absence of signiﬁcant segre- interior, with melt fractions above the rheological transition gation. > and cool relatively unsegregated, such that the trans Figure 4 shows the three metrics across a range of a 2 bulk concentrations of the fsp, pxn and olv components re- 5 2 3 [10 , 10 ] m and 2 [200, 1200] kg m . We ﬁnd that main similar to the initial composition (Figure 5A). However, grain size strongly controls the interior evolution of the plan- in magma sill models, the melt-rich layers above a low-melt- etesimals. If grain size remains below a < 10 m, melt den- fraction interior crystallize into distinct layers enriched in fsp, sity contrasts of < 500 kg m are not sufﬁcient to drive signif- pxn and olv. This stratiﬁcation reﬂects the component melt- fsp pxn 3 olv icant segregation. For density contrasts > 1000 kg m and ing points (T < T < T , Figure 5B). The melt-depleted m m for larger grain sizes, melt segregation becomes signiﬁcant, central parts of the planetesimal are strongly depleted in fsp as evidenced by a step-increase in each of the three met- and somewhat less in pxn. In general, such compositional rics. However, we ﬁnd that since initially fsp-enriched melts layering forms during the cooling stage and therefore does migrate on time scales comparable or longer than t 26 ( not cause substantial temperature inversion. The densities 1=2; Al 0.72 Myr), the temperature inversion effect remains minor of the minerals represented by the pseudo-components sug- throughout (Figure 4, panels B & F). gest that such layering would be dynamically stable. Radius / = const. = const. Radial profiles 2 4 0.9 0.9 0.9 0.7 0.7 0.7 0.5 0.5 0.5 0.3 0.3 0.3 0.1 1 1 0.1 3 0.1 A B C D -5 -4 -3 -2 200 700 1200 0.2 0.5 0.8 0.2 0.5 0.8 log [m] [kg/m ] 0.9 0.9 0.7 0.7 90 0.5 0.5 1 3 0.3 0.3 1 2 0.1 0.1 E F G H -5 -4 -3 -2 200 700 1200 1400 1600 1400 1600 log [m] [K] [kg/m ] 0.9 0.9 0.7 0.7 0.7 0.5 0.5 0.5 0.3 0.3 0.3 0.1 1 1 3 0.1 0.1 I J K L -5 -4 -3 -2 200 700 1200 0.1 0.35 0.6 0.1 0.35 0.6 log [m] [kg/m ] [wt] Figure 4: Parameter study of the inﬂuence of grain size a and density contrast on melt segregation , temperature inversion T , and compositional stratiﬁcation c , for planetesimals with R = 60 km, t = 1.5 t . Panels (A,E,I) show the metric deviation for constant density contrast = 600 kg fsp P form 26 1=2; Al 3 4 3 3 m and varying grain size a , indicating a steep gradient between grain sizes of 10 m and 10 m. For these two values ﬁxed (blue: a = 10 m, green: 0 0 a = 10 m), panels (B,F,J) display the metric deviations for varying . Here, variations in density contrast are outweighed by those in grain size. Models 3 4 with a = 10 m feature notable melt segregation, temperature inversions, and compositional differentiation. Models with a = 10 m only do so towards the 0 0 high end of density contrasts, & 700 kg m . Panels (C/D, G/H, K/L) show the radial proﬁles for the end-member models of the variations from (A/B, E/F, I/J). In general, variations in grain size a outweigh effects from increasing density contrast . Magma sill structures only form for grain sizes a > 10 m. 0 0 Magma ocean Magma sill 4.3. Magma dynamics regimes A B Figure 6 shows the regimes of melt segregation and com- positional stratiﬁcation for different formation times, t , and form 0.8 melt segregation numbers, R . We quantify the boundaries seg fsp-enriched separating the three characteristic regimes as follows. pxn-enriched 0.6 (I) Magma ocean regime: > := 0.5, where the olv-enriched ctr crit whole interior melts to above the disaggregation thresh- old. 0.4 (II) Magma sill regime: > , where melt segregation crit fsp-depleted 0.2 generates a melt-rich layer beneath the lid and depletes the interior of melt. fsp (III) Undifferentiated regime: c < c = 0.15, where 0.2 0.5 0.8 0.2 0.5 0.8 fsp bulk;0 melting and melt segregation do not redistribute a sub- [wt] after 10 Myr stantial amount of fsp. Figure 5: Compositional stratiﬁcation after cooling and crystallization of In addition to these segregation and differentiation criteria, magma beneath the primordial lid for magma ocean (A) and magma sill - we show which models are most affected by substantial tem- type (B) models. Magma sill cases with intermediate temperatures and thus high concentrations of fsp in the upper layers produce this signature. perature inversions, T > 250 K. These inversions occur both for magma sill and magma ocean models and reﬂect rapid magma ascent on time scales shorter than t 26 . We 1=2; Al fsp [K] [wt] pxn olv Compositional stratification Temperature Melt fsp anomaly segregation pxn olv / / / Figure 6: Evolution of silicate melt segregation with formation time t versus melt segregation number R at a reference melt fraction of = 0:1. form seg 0 Colormap values are plotted for the time of peak melting for each model (circles). We identify three primary melt dynamics regimes. (I) Magma ocean models, where melting occurs more rapidly than melt migration, feature high melt fractions above the rheological transition in their center, > := 0.5. ctr crit Magma ocean-type evolution is preferred for early t and low R , i.e., small a and . (II) Magma sill models feature efﬁcient melt segregation, with form seg 0 additional compositional stratiﬁcation towards cooling-down of the planetesimals (cf. Figure 5). Magma sill -type evolution is preferred for intermediate t form 0.5–1.75 t 26 , and high R . (III) Undifferentiated models never show melt fractions > , and never experience substantial compositional seg ctr crit 1=2; Al redistribution. They are preferred for late formation times t & 2.0 t . In addition to these three regimes, we show the region of increasing form 26 1=2; Al temperature inversion in yellow. ﬁnd that the magma sill regime generally occurs at higher Sato, 1984; Fu and Elkins-Tanton, 2014) at a 10 m. segregation numbers—at larger grain sizes or elevated den- Moreover, we ﬁnd that even in the case of signiﬁcant magma sity contrasts—and formation times less than 1.5 Myr after redistribution into shallow sills, the segregation time scale is CAIs, with a peak at around 1 Myr after CAIs. Very early for- comparable to the evolutionary time scale of the protoplane- mation time, t 1 Myr, or lower melt segregation number, tary disk and thus the accretion time scale. Therefore, scal- form R . -1.5, favor magma ocean-type evolution. Formation ing analysis alone (see Section 2), does not adequately cap- seg later than 1.5 Myrs after CAIs results in limited melting and ture the time-dependent competition between melting, parti- undifferentiated planetesimals. Temperature inversions oc- tioning, and melt transport. Because of the protracted onset cur for t . 1 Myr and R & -1.5. of melt ascent during heat-up, Al-hosting phases release form seg at least parts of the Al decay energy in the deeper region of planetesimals, and substantial temperature inversions on 5. Discussion the order of a few hundred K are only observed for extreme cases with R & -1, or early formation times t . 1 Myr seg form 5.1. Parametric controls on magma segregation with R & -2. seg The models above present thermo-chemically coupled At t = 26 1 Myr after CAIs the models show a clus- form Al two-phase ﬂow calculations that resolve the partitioning of tering of magma sill cases. This peak is due to the competi- the major rock-forming components and their transport by tion between heating and segregation, as introduced in Sec- magma. Using this method, we show that magma segre- tion 2. For formation times earlier than 1 Myr, heating, and gation in planetesimals formed within 2 Myr after CAIs was thus melting, is dominant, and new melt is generated faster potentially signiﬁcant if magma ascent was rapid with re- than transported away. Around t 1 Myr, melt upwelling form spect to the rate of melt production controlled by Al de- velocities become fast enough to exceed melt generation. Fi- cay. Expressed in terms of the melt segregation number R , seg nally, at later times (t & 1.5 Myr after CAIs or 2 t 26 , form 1=2; Al magma sill models ( > 0:5) were produced for R -1.5. seg respectively) melt production remains low and minimal to no This requires that (i) the average grain size in these planetes- segregation occurs. imals was comparable to or larger than chondrules, on the order of a 10 m, or (ii) in reducing environments (iron- Melt–rock density contrast is thought to be controlled by wÃijstite IW-2.5, or 1200 kg m , respectively, Brett and the oxidation state of the body. The chemical composition of 9 accreting planetesimals is inherited from the oxidation state ﬂow regime in the magma, leading to a variety of possible in the midplane of the protoplanetary disk. Hence, the time scenarios (Solomatov, 2015). and place of formation is expected to inﬂuence the subse- As a comparison, grain sizes of meteorite classes span quent interior evolution as it relates to the effects of melt seg- a wide range, from m-sized dust to pallasite-like, cm- regation. For example, planetesimals accreted towards the sized phenocrysts (Hutchison, 2004). Chondrites, gener- inner disk may feature lower oxygen fugacities and therefore ally regarded to be the most pristine rock samples from higher , compared to planetesimals accreted further out the early solar system, display a bimodal size distribution, (Rubie et al., 2015). In our models, we ﬁnd that the effect split between m-sized dust (‘matrix’) and chondrules, which 4 3 of density contrast on magma segregation is of secondary show characteristic sizes of 10 m to 10 m (Friedrich importance. If grain sizes were too small to allow for a suf- et al., 2015), with drastically differing textures and min- ﬁcient permeability, variations in density contrast would not eral grain sizes. The ratio of chondrules to matrix varies have led to signiﬁcantly different outcomes. This ﬁnding is greatly between different chondrite groups, resulting in com- in contrast to previous studies (Moskovitz and Gaidos, 2011; plex mixtures and grain size distributions. Meteorites that Wilson and Keil, 2012; Šrámek et al., 2012; Fu and Elkins- likely underwent partial melting (‘primitive achondrites’), like Tanton, 2014; Wilson and Keil, 2017; Fu et al., 2017), which Acapulcoite-Lodranites, Winonites and Brachinites, display did not consider magnitude variations in grain size, or relied grain sizes around 10 m (Hutchison, 2004; Wilson and on grain growth by pure Ostwald ripening without taking into Keil, 2017). Ureilites, Aubrites, and Angrites, which come account mechanisms of growth inhibition and grain destruc- from bodies with more extensive, or even body-wide, silicate tion (Neumann et al., 2013, 2014). melting feature larger grains, on the order of 10 m. How- ever, these textures are the end result of million-year long 5.2. Implications for the role of grain size evolutionary processes, and may have undergone repeated destruction and reaccretion cycles that reset their thermal A main conclusion from our models is that the mean grain histories and textures (Lichtenberg et al., 2018). Therefore, size a is the dominant parameter controlling the magma based on the grain sizes observed in meteorites, we are transport rate inside planetesimals. There are two main strongly limited in assessing the grain size evolution in plan- mechanisms that affect grain size during planetesimal ac- etesimal interiors at the time of their ﬁrst melting. cretion and early evolution. First, grain sizes may differ de- Interpretation of our results in this context suggests that pending on the orbital location and mineralogical composi- the planetesimal interior evolution and the redistribution of tion (van Boekel et al., 2004) of the dust that agglomer- chemical and isotopic heterogeneities during planetary ac- ates to form the planetesimals. Coagulated dust aggre- cretion can be inﬂuenced by the planetesimal formation gates on the order of a 10 m are in the critical size mechanism, its accretion location, and the local composi- regime for fast radial drift towards the protostar, depending tional distribution of grains in the protoplanetary disk. Fur- on orbital location and evolutionary state of the protoplan- ther studies of planetesimal formation and mineralogical in- etary disk (Weidenschilling, 1977). However, these sizes ventory are needed to link their formation processes to their may facilitate planetesimal formation from local dust-gas in- subsequent dynamic evolution. The local grain-size distribu- teractions (Johansen et al., 2015; Birnstiel et al., 2016) and tion within planetesimals may be evolving during rapid grav- can trigger planetesimal formation at various orbital sepa- itational collapse (Wahlberg Jansson and Johansen, 2017) rations and times (Dra ¸ zk ˙ owska and Dullemond, 2018) with or more gradual growth (Kataoka et al., 2013). Also, ongoing varying dust particle distributions within newly-formed bod- accretion, for instance due to subsequent pebble accretion ies (Wahlberg Jansson and Johansen, 2017). Second, grain (Visser and Ormel, 2016), may inﬂuence whether magma sizes may evolve during the heating and melting in the plan- can reach the surface through fractures. Future experimental etesimal interior after accretion. In this process, the grain studies on grain size evolution of partially molten aggregates size evolves due to competing mechanisms for growth and spanning the meteoritic compositional space will be needed destruction (Rozel et al., 2011; Bercovici and Ricard, 2016). to advance our understanding of melt migration in early solar Grain growth by either normal grain growth before the ﬁrst system planetesimals. melts arise, or Ostwald ripening during partial melting in a solid–liquid aggregate, is driven by the reduction of interfa- 5.3. Implications for chemical differentiation cial energies, and competes with size reduction due to the Recently, it was shown that differentiation by percolation presence of secondary solid phases (Hiraga et al., 2010) of Fe,Ni-FeS liquids in primordial planetesimals may not be and mechanical work due to solid-state deformation in plan- complete and that at least some material remains trapped in etesimal interiors (Tkalcec et al., 2013). During later stages, the rock matrix (Bagdassarov et al., 2009; Cerantola et al., when the melt fraction reaches or exceeds the disaggrega- 2015). But once silicate melting has reached the disaggre- tion threshold, grain sizes are governed by the convective gation threshold, the remaining metal droplets will rain out rapidly towards the forming core (Lichtenberg et al., 2018). Note that in astrophysics literature grain size usually denotes the Therefore, even though the models in this study do not in- characteristic size of porous dust aggregates embedded in the disk ﬂow, clude a metal phase, they allow leading order predictions re- whereas we here use the term grain size, a , for the characteristic size of mineral grains in a granular, polymineralic rock aggregate. garding core formation. 10 In the magma ocean case, we expect core formation to be 5.4. Limitations & future directions rapid, with nearly complete loss of metals to the core. In the One of the limitations of our model is the use of a 1D case of a magma sill, a two-step process may occur. First, Cartesian geometry, which introduces systematic errors as a small proto-core may form from incomplete percolation. compared to a 1D spherical geometry assuming radial sym- Then, after the formation of the sill structure, the remain- metry. Among these errors, our model over-predicts the ing metal within this region may rain out and accumulate at heating-to-cooling ratio of planetesimals. Because we focus the interface between the melt-depleted deep interior and the on planetesimals of 60 km radius that show a nearly isother- magma sill zone. This emerging metal pool will either per- mal evolution in the deep interior (cf. Figure 3, Castillo- colate downwards or form diapirs sinking through the weak- Rogez and Young, 2017; Hunt et al., 2018), the heat-up ened partially molten rock. The thermo-mechanical evolu- phase is consistent with a radial model. However, geometric tion of such a two-step process needs to be tested by taking errors result in an under-prediction of the rate of heat loss into account metal phases in self-consistent multi-phase ﬂow through the surface. Similar problems apply to melt segre- models, in order to make robust predictions that can be com- gation velocities and melting rates. For example, the liquid pared to laboratory studies of the core formation process. mass conservation equation, in spherical coordinates with Neumann et al. (2018) recently suggested a multi-stage core radial symmetry, contains the term 2w =r, which is neglected formation scenario for the IVB parent body, which is quali- here (a similar term appears in the solid mass conserva- tatively consistent with the magma sill regime we propose tion equation). For radially outward melt migration, this term based on our models. has the effect of reducing @=@t; its absence therefore leads to an over-prediction of melt fraction. To quantify the er- The limited melt segregation in our undifferentiated mod- ror introduced by neglecting this term we compared its size, els may offer an explanation for the absence of remnant dif- computed a posteriori, with the ﬂux divergence @w =@z that ferentiated materials in the asteroid belt (Weiss et al., 2012; appears in the same equation. We made this comparison Mandler and Elkins-Tanton, 2013). Conversely, chemical across models of the magma ocean and magma sill regimes. stratiﬁcation resulting from melt segregation may help to ex- The results indicate that the geometrical error can reach a plain the paucity of olivine-rich deposits on 4 Vesta’s surface magnitude comparable to the divergence. However, through- (Clenet et al., 2014; Consolmagno et al., 2015; Raymond out most of the domain and once magma sill structures are et al., 2017). Furthermore, the magma sill models and the established, the term is negligible. resulting chemical stratiﬁcation we describe here are consis- In this work we consider a diffusion creep rheology only. A tent with earlier work by Neumann et al. (2013, 2014, “shal- more realistic rheology would be a composite of both diffu- low magma ocean”) and Mizzon (2015, “completely liquid sion and dislocation creep (e.g., Bercovici and Ricard, 2016; layer”), predicting a subsurface layer of accumulated silicate Mulyukova and Bercovici, 2018). At the highest temper- melts below a cold lid resulting from efﬁcient melt segrega- atures and grain sizes tested here, model behavior would tion (cf. discussion in Raymond et al., 2017; Neumann et al., likely fall into the dislocation creep regime, where the matrix 2018). viscosity becomes independent of grain size. In this case, the grain size sensitivity of the compaction length (Equation Finally, our results indicate that the importance of melt 1) would decrease. segregation in planetesimal interiors varied substantially and Furthermore, we use a constant melt viscosity of 1 Pa s affected the redistribution of heat-producing elements, such (Table 1). However, silicate melt viscosity varies with tem- as Al, during melt ascent. In the case of our magma sill perature and composition (e.g., Moskovitz and Gaidos, 2011; end-member models, we also expect other incompatible el- Mizzon, 2015). For the compositional space explored here, a ements to be preferentially segregated to shallow layers of temperature dependent viscosity of 1–100 Pa s may be con- a planetesimal. The crustal stripping paradigm of plane- sidered realistic. As we are interested in the consequences tary accretion assumes that frequent impacts during plan- of melt segregation in planetesimals here, we chose a rea- etary growth eroded and redistributed signiﬁcant amounts sonable lower limit. In addition, because the ratio of per- of shallow materials between planetesimals. The strongly meability to melt viscosity controls the upwelling timescales variable degree of melt segregation, and the resulting range (grain size squared versus linear melt viscosity, Equation 1), of variably differentiated major, trace, and isotopic composi- the order of magnitude variability in grain size outweighs po- tions of shallow planetesimal layers could result in compo- tential variations in melt viscosity. sitional differences between forming planets and chondritic meteorites (e.g., Carter et al., 2018). Our simulations indi- As a consequence of our limiting assumptions on the ge- cate that magma ascent governing material redistribution to ometry, rheology, and melt viscosity, the extent of the magma the planetesimal crust occurs on a Myr time scale, i.e., sill regime in Figure 6 may be overestimated. More efﬁ- comparable to the collisional time scale of planetary accre- cient cooling in a spherical geometry would lead to decreas- tion. This suggests that models quantifying compositional ing melt fractions and thus lower migration speeds (Figure deviations between planets and chondrites should take into 1), as would a stronger heat ﬂux from turbulent convection account the evolution of planetesimal interiors during plane- above the disaggregation threshold (Solomatov, 2015), and tary accretion. the geometric inﬂuence in spherical or higher-dimensional 11 geometry. A weaker dependence on grain size in the models the total retention of magma on the inside of planetesimals. with the highest temperatures and largest grain sizes would However, this does not affect our conclusions or any con- also decrease migration speeds. Therefore, while our mod- straints on temperature inversions of planetesimals unless els constrain the possible phase space of melt migration in this transport is faster or comparable to the magma segrega- early solar system planetesimals, more complex models not tion in the interior. bound to the above limitations would result in a reduced sta- bility ﬁeld where melt segregation becomes dominant than 6. Summary & conclusions shown in Figure 6. Our models do not feature a metal phase that would al- In this study we investigated magma genesis and redis- low the direct resolution of metal percolation (Ghanbarzadeh tribution in planetesimals during and shortly after the pro- et al., 2017; Neumann et al., 2018), and therefore our re- toplanetary disk phase. Using scaling relations we demon- sults only allow for qualitative inferences about possible core strated that the interior evolution of planetesimals sensitively formation scenarios. However, the complexities of textural depends on a variety of model parameters, with the grain equilibrium phase topologies are not yet fully understood (cf. size exerting the primary control on melt segregation. Based Rudge, 2018), in particular for systems with several immis- on average chondritic abundances of the most common min- cible liquid phases. For example, the wetting angles formed eral groups in meteorites, we calculated the composition for between metal liquids with silicate minerals in the presence rock–melt aggregates comprising idealized components us- of a silicate melt remains unclear, leaving open the debate ing a reactive, multi-component melting model. We quanti- around a possible percolation threshold for metal liquids at ﬁed the effects on Al partitioning and magma ascent with low melt fractions (Cerantola et al., 2015). Further work a coupled, two-phase ﬂow model. We deﬁned the melt seg- needs to be undertaken to better quantify these effects. regation number R as the ratio between the heating and seg Our choice of melting model in the form of a ternary ideal melt transport time scales in a planetesimal, which estab- solution limits the degree to which the model may represent lishes the leading order parametric control on the propensity natural melting behavior. For example, our model does not for magma redistribution during the heating stage of plan- reproduce the eutectic behavior expected for silicate compo- etesimal evolution. We predicted that the primary two con- sitions as the ones considered here (see discussion in Keller trols are the melt–rock density contrast and the mineral and Katz, 2016), nor does it include volatiles and incompat- grain size a . ible elements producing low-degree, incompatible-enriched Investigating the relative importance of model parameters melts at temperatures below the anhydrous solidus. Using for the evolution of planetesimals, we categorize model out- a more consistent petrological model that takes into account comes in terms of their melt segregation numbers R and seg the non-ideal thermodynamics of the full range of major ele- their formation times t . Using this scheme, we ﬁnd: form ments and mineral phases would likely lead to more complex relations between heating, melt production, and element par- Planetesimal melt migration behavior can be classiﬁed titioning (Keller and Suckale, 2018). Signiﬁcant differences in three general regimes: in aluminum partitioning, which is the focus here, are likely conﬁned to the onset of melting, where low-degree, enriched – The magma ocean regime with global interior melts may have important control on geochemical evolution. magma oceans, where rapid melting overwhelms Moreover, our dry models neglect volatile-driven eruptions, melt segregation. which were previously discussed as a catalyst for upward – The magma sill regime, where global interior transport via explosive volcanism (Fu et al., 2017; Wilson magma oceans are prevented by rapid magma and Keil, 2017). If substantial volatile quantities could be transport that concentrates melt in sills beneath retained, this mechanism would decrease the retention of the cool lid. magma in the upper layers of planetesimals and potential chemical stratiﬁcation upon crystallization of the silicate ma- – The undifferentiated regime with a low degree terial. However, the relatively low pressure at the planetes- of melting, minor melt segregation, and therefore imal sizes we consider disfavors a high volatile solubility in chemically primordial and largely undifferentiated silicate magma, and degassing would therefore be expected interiors. to be nearly complete in the earliest stages of melting and segregation (Monteux et al., 2018). Magma sill models crystallize to a compositionally strat- Finally, the melt in our magma ocean and magma sill mod- iﬁed structure, with shallow depth layers dominantly en- els cannot breach the cold surface layers, as the temperature riched in the low melting point components, feldspar and is too low for melt to exist in the porous medium. Our sim- pyroxene, and a paucity of high melting point compo- ulations cannot resolve potential melt transport through the nents, such as olivine (cf. Ghosh and McSween, 1998; upper lid by fracturing (Keller et al., 2013) or gravitational in- Mizzon, 2015). The crystallization sequence, and thus stability in the layered structure in Figure 5, which may bury the compositional stratiﬁcation, however, may be af- the primitive lid (Wilson and Keil, 2012), and lead to efﬁcient fected by convective motions beyond the disaggregation heat loss and magmatic resurfacing. This would decrease limit, which we do not model here. 12 Magma ocean and magma sill models show tempera- constructive comments by David Bercovici and an anony- ture inversions for high R and early t , where the mous referee, which helped to improve the manuscript. TL seg form temperatures in the shallow- to mid-mantle are higher was supported by ETH Zürich Research Grant ETH-17 13-1 than at the center of the planetesimal. These thermal in- and acknowledges partial ﬁnancial support from a MERAC versions, however, are restricted to formation times t travel grant of the Swiss Society for Astrophysics and As- form . 1 Myr and R & 1:5 for T 250 K. Therefore, the tronomy and from the National Center for Competence in seg majority of planetesimals likely underwent thermal evo- Research PlanetS supported by the Swiss National Science lutionary scenarios that can be qualitatively reproduced Foundation. TK acknowledges ﬁnancial support from Jenny with models that do not take into account melt segre- Suckale, Stanford University. TK and RFK received fund- gation and Al partitioning, depending on the level of ing from the European Research Council under the Euro- detail that needs to be assessed. pean Union’s Seventh Framework Programme (FP7/2007– 2013)/ERC grant agreement number 279925. The numeri- The magma sill regime can be achieved depending on cal simulations in this work were performed on the EULER a combination of a few key parameters: computing cluster of ETH Zürich, and were analyzed using the open source software environments MATPLOTLIB and – The formation time t controls the total amount form SEABORN . of energy available from Al and restricts the regime for melt segregation to 0.5 t 26 . t 1=2; Al form . 1.75 t 26 , with a peak at 1 Myr after CAIs, Supplementary Materials 1=2; Al when the rate of melt transport dominates over the A simulation video associated with this article can be found generation of new melts. attached to the publication and at this URL. The video shows – The grain size a controls the rate of melt trans- a comparison between magma ocean and magma sill end- port and thus whether a planetesimal with sufﬁ- member models. Magma ocean stages are indicated in yel- cient heating evolves towards a melt segregated low. The H3 annotation describes the heating value below structure. Below a . 0.1 mm no segregation which radioactive heating from Al is inefﬁcient. Shown are occurs; above a & 1 mm segregated structures various parameters for both models, from left to right: melt form. fraction , temperature T , radiogenic heating H, melt/liquid – The solid–melt density contrast is of secondary upwelling velocity w , composition fraction of liquid c , and liq liq importance, but can enhance melt segregation in composition fraction of solid c . The composition is broken- sol the regime transition from a = 0.1 mm to 1 mm. down into the deﬁned pseudo-components fsp, pxn, and olv. The model is started (= planetesimal formation time) at 0.9 To conclude, in this paper we have advanced the techni- Myr after CAIs. cal capabilities to simulate multi-phase and multi-component planetesimal evolution, gaining insights into features not ac- References cessible to single-phase ﬂuid dynamics models. However, unraveling more detailed evolutionary regimes of planetes- A’Hearn, M. F., 2017. The future of NASA’s missions. Nat. Astron. 1, 0095. imals will require a time-dependent treatment that includes Bagdassarov, N., Golabek, G. J., Solferino, G., Schmidt, M. W., 2009. Con- straints on the Fe-S melt connectivity in mantle silicates from electrical metal and volatile phases, which shape the structure and impedance measurements. Phys. Earth Planet. Inter. 177, 139–146. subsequent evolution of these bodies (Keller and Suckale, Balay, S., Gropp, W. D., McInnes, L. C., Smith, B. F., 1997. Efﬁcient man- 2018). Further work is required to understand planetesi- agement of parallelism in object oriented numerical software libraries. mal evolution and its connection to the meteoritic record and In: Arge, E., Bruaset, A. M., Langtangen, H. P. (Eds.), Modern Software Tools in Scientiﬁc Computing. Birkhäuser Press, pp. 163–202. rocky planet formation (Lichtenberg et al., 2018). In partic- Bercovici, D., Ricard, Y., Apr. 2016. Grain-damage hysteresis and plate tec- ular, robust scaling laws for the evolution of grain sizes in tonic states. Physics of the Earth and Planetary Interiors 253, 31–47. partially molten and heated systems relevant for planetesi- Birnstiel, T., Fang, M., Johansen, A., 2016. Dust Evolution and the Forma- mal settings are required to establish more detailed regimes tion of Planetesimals. Space Sci. Rev. 205, 41–75. Brett, R., Sato, M., 1984. Intrinsic oxygen fugacity measurements on seven for melt segregation. chondrites, a pallasite, and a tektite and the redox state of meteorite In the mid-term, future spacecraft missions (A’Hearn, parent bodies. Geochim. Cosmochim. Acta 48, 111–120. 2017) may be able to deliver further observational con- Burbine, T. H., Demeo, F. E., Rivkin, A. S., Reddy, V., 2017. Evidence for straints on asteroid-belt objects. In conjunction with self- differentiation among asteroid families. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: Early Differentiation and Consequences for consistent evolutionary models of metal–silicate and solid– Planets. Cambridge University Press, Cambridge, pp. 298–320. melt segregation, these can help to further decipher the in- Carporzen, L., Weiss, B. P., Elkins-Tanton, L. T., Shuster, D. L., Ebel, D., terior evolution of planetary bodies in the early solar system Gattacceca, J., 2011. Magnetic evidence for a partially differentiated car- bonaceous chondrite parent body. Proc. Natl. Acad. Sci. U.S.A. 108, and sharpen our understanding of terrestrial planet formation 6386–6389. in the solar system and elsewhere. Acknowledgements. The authors acknowledge stimulating 2 matplotlib.org (Hunter, 2007) discussions with Ian Sanders and Guy Consolmagno, and seaborn.pydata.org 13 Carter, P. J., Leinhardt, Z. M., Elliott, T., Stewart, S. T., Walter, M. J., Feb. Kataoka, A., Tanaka, H., Okuzumi, S., Wada, K., 2013. Fluffy dust forms icy 2018. Collisional stripping of planetary crusts. Earth and Planetary Sci- planetesimals by static compression. Astron. Astrophys. 557, L4. ence Letters 484, 276–286. Katz, R. F., Knepley, M., Smith, B., Spiegelman, M., Coon, E., 2007. Nu- Castillo-Rogez, J., Young, E. D., 2017. Origin and evolution of volatile-rich merical simulation of geodynamic processes with the portable extensible asteroids. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: toolkit for scientiﬁc computation. Phys. Earth Planet. Inter. 163, 52–68. Early Differentiation and Consequences for Planets. Cambridge Univer- Keller, T., Katz, R. F., 2016. The Role of Volatiles in Reactive Melt Transport sity Press, Cambridge, pp. 92–114. in the Asthenosphere. J. Petrol. 57, 1073–1108. Cerantola, V., Walte, N. P., Rubie, D. C., 2015. Deformation of a crystalline Keller, T., May, D. A., Kaus, B. J. P., 2013. Numerical modelling of magma olivine aggregate containing two immiscible liquids: Implications for early dynamics coupled to tectonic deformation of lithosphere and crust. Geo- core-mantle differentiation. Earth Planet. Sci. Lett. 417, 67–77. phys. J. Int. 195, 1406–1442. Clenet, H., Jutzi, M., Barrat, J.-A., Asphaug, E. I., Benz, W., Gillet, P., 2014. Keller, T., Suckale, J., 2018. A continuum model of multi-phase reactive A deep crust-mantle boundary in the asteroid 4 Vesta. Nature 511, 303– transports in igneous systems. arXiv:1809.00079. 306. Kita, N. T., Yin, Q.-Z., MacPherson, G. J., Ushikubo, T., Jacobsen, B., Na- 26 26 Consolmagno, G. J., Golabek, G. J., Turrini, D., Jutzi, M., Sirono, S., gashima, K., Kurahashi, E., Krot, A. N., Jacobsen, S. B., 2013. Al- Mg Svetsov, V., Tsiganis, K., 2015. Is Vesta an intact and pristine proto- isotope systematics of the ﬁrst solids in the early solar system. Meteorit. planet? Icarus 254, 190–201. Planet. Sci. 48, 1383–1400. Costa, A., Caricchi, L., Bagdassarov, N., 2009. A model for the rheology Lichtenberg, T., Golabek, G. J., Dullemond, C. P., Schönbächler, M., Gerya, of particle-bearing suspensions and partially molten rocks. Geochem. T. V., Meyer, M. R., 2018. Impact splash chondrule formation during plan- Geophys. Geosys. 10, Q03010. etesimal recycling. Icarus 302, 27–43. Cournede, C., Gattacceca, J., Gounelle, M., Rochette, P., Weiss, B. P., Lichtenberg, T., Golabek, G. J., Gerya, T. V., Meyer, M. R., 2016a. The Zanda, B., 2015. An early solar system magnetic ﬁeld recorded in CM effects of short-lived radionuclides and porosity on the early thermo- chondrites. Earth Planet. Sci. Lett. 410, 62–74. mechanical evolution of planetesimals. Icarus 274, 350–365. Dra ¸ zk ˙ owska, J., Dullemond, C. P., 2018. Planetesimal formation during pro- Lichtenberg, T., Parker, R. J., Meyer, M. R., 2016b. Isotopic enrichment of toplanetary disk buildup. Astron. Astrophys. 614, A62. forming planetary systems from supernova pollution. Mon. Not. R. As- Dunn, T. L., Cressey, G., McSween, Jr., H. Y. ., McCoy, T. J., 2010. Analysis tron. Soc. 462, 3979–3992. of ordinary chondrites using powder X-ray diffraction: 1. Modal mineral Lodders, K., 2003. Solar System Abundances and Condensation Tempera- abundances. Meteorit. Planet. Sci. 45, 123–134. tures of the Elements. Astrophys. J. 591, 1220–1247. Ebel, D. S., Grossman, L., 2000. Condensation in dust-enriched systems. Mandler, B. E., Elkins-Tanton, L. T., 2013. The origin of eucrites, diogenites, Geochim. Cosmochim. Acta 64, 339–366. and olivine diogenites: Magma ocean crystallization and shallow magma Elkins-Tanton, L. T., Weiss, B. P., Zuber, M. T., 2011. Chondrites as samples chamber processes on Vesta. Meteorit. Planet. Sci. 48, 2333–2349. of differentiated planetesimals. Earth Planet. Sci. Lett. 305, 1–10. McKenzie, D., 1984. The generation and compaction of partially molten Friedrich, J. M., Weisberg, M. K., Ebel, D. S., Biltz, A. E., Corbett, B. M., rock. J. Petrol. 5, 713–765. Iotzov, I. V., Khan, W. S., Wolman, M. D., 2015. Chondrule size and McSween, H. Y., Bennett, M. E., Jarosewich, E., 1991. The mineralogy of or- related physical properties: A compilation and evaluation of current data dinary chondrites and implications for asteroid spectrophotometry. Icarus across all meteorite groups. Chemie der Erde / Geochemistry 75, 419– 90, 107–116. 443. Mei, S., Bai, W., Hiraga, T., Kohlstedt, D. L., 2002. Inﬂuence of melt on the Fu, R. R., Elkins-Tanton, L. T., 2014. The fate of magmas in planetesimals creep behavior of olivine-basalt aggregates under hydrous conditions. and the retention of primitive chondritic crusts. Earth Planet. Sci. Lett. Earth Planet. Sci. Lett. 201, 491–507. 390, 128–137. Mizzon, H., 2015. The magmatic crust of vesta. Ph.D. thesis, Université de Fu, R. R., Young, E. D., Greenwood, R. C., Elkins-Tanton, L. T., 2017. Sil- Toulouse, Université Toulouse III-Paul Sabatier. icate melting and volatile loss during differentiation in planetesimals. In: Monteux, J., Golabek, G. J., Rubie, D. C., Tobie, G., Young, E. D., 2018. Wa- Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: Early Differentia- ter and the interior structure of terrestrial planets and icy bodies. Space tion and Consequences for Planets. Cambridge University Press, Cam- Sci. Rev. 214, 39. bridge, pp. 115–135. Moskovitz, N., Gaidos, E., 2011. Differentiation of planetesimals and the Ghanbarzadeh, S., Hesse, M. A., Prodanovic, ´ M., 2017. Percolative core thermal consequences of melt migration. Meteorit. Planet. Sci. 46, 903– formation in planetesimals enabled by hysteresis in metal connectivity. 918. Proc. Natl. Acad. Sci. 114, 13406–13411. Mulyukova, E., Bercovici, D., 2018. Collapse of passive margins by litho- Ghosh, A., McSween, H. Y., 1998. A Thermal Model for the Differentiation spheric damage and plunging grain size. Earth Planet. Sci. Lett. 484, of Asteroid 4 Vesta, Based on Radiogenic Heating. Icarus 134, 187–206. 341–352. Golabek, G. J., Bourdon, B., Gerya, T. V., 2014. Numerical models of Neumann, W., Breuer, D., Spohn, T., 2013. The thermo-chemical evolution the thermomechanical evolution of planetesimals: Application to the of Asteroid 21 Lutetia. Icarus 224, 126–143. acapulcoite-lodranite parent body. Meteorit. Planet. Sci. 49, 1083–1099. Neumann, W., Breuer, D., Spohn, T., 2014. Differentiation of Vesta: Implica- Hevey, P. J., Sanders, I. S., 2006. A model for planetesimal meltdown by tions for a shallow magma ocean. Earth Planet. Sci. Lett. 395, 267–280. Al and its implications for meteorite parent bodies. Meteorit. Planet. Neumann, W. O., Kruijer, T. S., Breuer, D., Kleine, T., 2018. Multi-Stage Sci. 41, 95–106. Core Formation in Planetesimals Revealed by Numerical Models and Hf- Hiraga, T., Tachibana, C., Ohashi, N., Sano, S., 2010. Grain growth system- W Chronometry of Iron Meteorites. J. Geophys. Res. Planets 123, 421– atics for forsterite enstatite aggregates: Effect of lithology on grain size 444. in the upper mantle. Earth Planet. Sci. Lett. 291, 10–20. Pätzold, M., Andert, T. P., Asmar, S. W., Anderson, J. D., Barriot, J.-P., Bird, Hirth, G., Kohlstedt, D., 2003. Rheology of the upper mantle and the mantle M. K., Häusler, B., Hahn, M., Tellmann, S., Sierks, H., Lamy, P., Weiss, wedge: A view from the experimentalists. Inside the subduction Factory, B. P., 2011. Asteroid 21 Lutetia: Low Mass, High Density. Science 334, 83–105. 491. Hunt, A. C., Cook, D. L., Lichtenberg, T., Reger, P. M., Ek, M., Golabek, Raymond, C. A., Russell, C. T., McSween, H. Y., 2017. Dawn at vesta: G. J., Schönbächler, M., 2018. Late metal-silicate separation on the IAB Paradigms and paradoxes. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), parent asteroid: Constraints from combined W and Pt isotopes and ther- Planetesimals: Early Differentiation and Consequences for Planets. mal modelling. Earth Planet. Sci. Lett. 482, 490–500. Cambridge University Press, Cambridge, pp. 321–339. Hunter, J. D., 2007. Matplotlib: A 2d graphics environment. Computing In Rozel, A., Ricard, Y., Bercovici, D., 2011. A thermodynamically self- Science & Engineering 9, 90–95. consistent damage equation for grain size evolution during dynamic re- Hutchison, R., 2004. Meteorites, Cambridge University Press. crystallization. Geophys. J. Int. 184, 719–728. Johansen, A., Mac Low, M.-M., Lacerda, P., Bizzarro, M., 2015. Growth Rubie, D. C., Jacobson, S. A., Morbidelli, A., O’Brien, D. P., Young, E. D., de of asteroids, planetary embryos, and Kuiper belt objects by chondrule Vries, J., Nimmo, F., Palme, H., Frost, D. J., 2015. Accretion and differ- accretion. Sci. Adv. 1, 1500109. entiation of the terrestrial planets with implications for the compositions 14 of early-formed Solar System bodies and accretion of water. Icarus 248, 89–108. Rudge, J. F., 2018. Textural equilibrium melt geometries around tetrakaidec- ahedral grains. Proc. Royal Soc. Lond. A 474, 20170639. Rudge, J. F., Bercovici, D., Spiegelman, M., 2011. Disequilibrium melting of a two phase multicomponent mantle. Geophys. J. Int. 184, 699–718. Sierks, H., Lamy, P., Barbieri, C., Koschny, D., Rickman, H., Rodrigo, R., A’Hearn, M. F., et al., 2011. Images of Asteroid 21 Lutetia: A Remnant Planetesimal from the Early Solar System. Science 334, 487. Solomatov, V. S., 2015. Magma oceans and primordial mantle differentia- tion. Treatise on Geophysics 2nd ed., pp. 81–104. Stevenson, D. J., 1990. Fluid dynamics of core formation. In: Newsom, H. E., Jones, J. H. (Eds.), Origin of the Earth. pp. 231–249. Tang, H., Dauphas, N., 2015. Low Fe Abundance in Semarkona and Sa- hara 99555. Astrophys. J. 802, 22–31. Tkalcec, B. J., Golabek, G. J., Brenker, F. E., 2013. Solid-state plastic de- formation in the dynamic interior of a differentiated asteroid. Nat. Geo. 6, 93–97. Šrámek, O., Milelli, L., Ricard, Y., Labrosse, S., 2012. Thermal evolution and differentiation of planetesimals and planetary embryos. Icarus 217, 339–354. van Boekel, R., Min, M., Leinert, C., Waters, L. B. F. M., Richichi, A., Ches- neau, O., Dominik, C., Jaffe, W., Dutrey, A., Graser, U., Henning, T., de Jong, J., Köhler, R., de Koter, A., Lopez, B., Malbet, F., Morel, S., Paresce, F., Perrin, G., Preibisch, T., Przygodda, F., Schöller, M., Wit- tkowski, M., 2004. The building blocks of planets within the ‘terrestrial’ region of protoplanetary disks. Nature 432, 479–482. Visser, R. G., Ormel, C. W., 2016. On the growth of pebble-accreting plan- etesimals. Astron. Astrophys. 586, A66. Wahlberg Jansson, K., Johansen, A., 2017. Radially resolved simulations of collapsing pebble clouds in protoplanetary discs. Mon. Not. R. Astron. Soc. 469, 149–157. Wasson, J. T., 1990. Ungrouped iron meteorites in Antarctica - Origin of anomalously high abundance. Science 249, 900–902. Weidenschilling, S. J., 1977. Aerodynamics of solid bodies in the solar neb- ula. Mon. Not. R. Astron. Soc. 180, 57–70. Weiss, B. P., Elkins-Tanton, L. T., 2013. Differentiated Planetesimals and the Parent Bodies of Chondrites. Annu. Rev. Earth Planet. Sci. 41, 529–560. Weiss, B. P., Elkins-Tanton, L. T., Antonietta Barucci, M., Sierks, H., Snod- grass, C., Vincent, J.-B., Marchi, S., Weissman, P. R., Pätzold, M., Richter, I., Fulchignoni, M., Binzel, R. P., Schulz, R., 2012. Possi- ble evidence for partial differentiation of asteroid Lutetia from Rosetta. Planet. Space Sci. 66, 137–146. Wilson, L., Keil, K., 2012. Volcanic activity on differentiated asteroids: A review and analysis. Chemie der Erde / Geochemistry 72, 289–321. Wilson, L., Keil, K., 2017. Arguments for the non-existence of magma oceans in asteroids. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Plan- etesimals: Early Differentiation and Consequences for Planets. Cam- bridge Planetary Science. Cambridge University Press, pp. 159–179. Yoshino, T., Walter, M. J., Katsura, T., 2003. Core formation in planetesimals triggered by permeable ﬂow. Nature 422, 154–157. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Astrophysics arXiv (Cornell University) http://www.deepdyve.com/lp/arxiv-cornell-university/magma-ascent-in-planetesimals-control-by-grain-size-H2lXM0UQOU

Loading next page...

References

References for this paper are not available at this time. We will be adding them shortly, thank you for your patience.

ISSN: 0012-821X
eISSN: ARCH-3330
DOI: 10.1016/j.epsl.2018.11.034
Publisher site: See Article on Publisher Site

Abstract

Rocky planetesimals in the early solar system melted internally and evolved chemically due to radiogenic heating from Al. Here we quantify the parametric controls on magma genesis and transport using a coupled petrological and ﬂuid mechanical model of reactive two-phase ﬂow. We ﬁnd the mean grain size of silicate minerals to be a key control on magma ascent. For grain sizes & 1 mm, melt segregation produces distinct radial structure and chemical stratiﬁcation. This stratiﬁcation is most pronounced for bodies formed at around 1 Myr after formation of Ca,Al-rich inclusions. These ﬁndings suggest a link between the time and orbital location of planetesimal formation and their subsequent structural and chemical evolution. According to our models, the evolution of partially molten planetesimal interiors falls into two categories. In the magma ocean scenario, the whole interior of a planetesimal experiences nearly complete melting, which would result in turbulent convection and core-mantle differentiation by the rainfall mechanism. In the magma sill scenario, segregating melts gradually deplete the deep interior of the radiogenic heat source. In this case, magma may form melt-rich layers beneath a cool and stable lid, while core formation would proceed by percolation. Our ﬁndings suggest that grain sizes prevalent during the internal heating stage governed magma ascent in planetesimals. Regardless of whether evolution progresses toward a magma ocean or magma sill structure, our models predict that temperature inversions due to rapid Al redistribution are limited to bodies formed earlier than 1 Myr after CAIs. We ﬁnd that if grain size was . 1 mm during peak internal melting, only elevated solid–melt density contrasts (such as found for the reducing conditions in enstatite chondrite compositions) would allow substantial melt segregation to occur. Keywords: Planetary formation, Planetesimals, Magma ocean, Melt migration, Chemical differentiation, Achondrites 1. Introduction The interior evolution of early solar system planetesimals has broad implications for the formation of both rocky planets At the time of planet formation, the inner solar system was and main-belt asteroid populations, the most immediate rem- populated by rocky planetesimals that seeded today’s ter- nants of the accretion process. Meteorites, broken-up pieces restrial planets through dynamical accretion of many smaller of asteroids fallen to Earth, are currently our only source bodies (Johansen et al., 2015), and whose internal evolu- of direct evidence from the early solar system. Therefore, tion was governed by radiogenic heating from short-lived Al our understanding of planetary growth and evolution is fun- (Hevey and Sanders, 2006). For large planetesimal radii and damentally limited by our ability to reconstruct the thermo- sufﬁcient Al incorporated upon formation, the released en- chemical evolution of planetesimals as evidenced by mete- ergy led to volatile degassing (Castillo-Rogez and Young, orites. Achondritic meteorites, which are thought to originate 2017; Monteux et al., 2018) and signiﬁcant silicate melting, from differentiated planetesimals, show a remarkable diver- surpassing the rheological transition from solid-state creep sity and likely originate from more than 50–100 parent bodies to disaggregation and melt-dominated deformation at melt (Wasson, 1990). However, spectral properties of asteroids fractions & 0:4–0:6 (Costa et al., 2009). In comparison trans do not match this diversity, as most known asteroids with an with solid or partially molten interiors, which lose heat by con- achondritic surface are interpreted to represent the debris of duction and/or laminar convection, disaggregation results in only a few parent bodies (Burbine et al., 2017). This appar- signiﬁcantly increased heat ﬂux by turbulent convection and ent lack of achondritic asteroids is at odds with the available rapid metal-silicate differentiation by the raining out of iron meteorite record. droplets (Stevenson, 1990). A possible solution to this conundrum is that internally dif- ferentiated planetesimals can retain their primitive, chondritic surfaces if magma remains conﬁned to the interior instead Corresponding author, now at: Atmospheric, Oceanic and Planetary of being erupted by volcanism (Elkins-Tanton et al., 2011; Physics, University of Oxford, Parks Rd, Oxford OX1 3PU, United Kingdom, email: tim.lichtenberg@physics.ox.ac.uk. Weiss and Elkins-Tanton, 2013). Some paleomagnetic stud- Published in Earth and Planetary Science Letters January 8, 2019 arXiv:1802.02157v2 [astro-ph.EP] 7 Jan 2019 ies on CV and CM meteorites suggest previous dynamo ac- mordial rock, with parts of the latter potentially lost to the tivity consistent with this hypothesis (Carporzen et al., 2011; core by percolation before the onset of major silicate melt- Cournede et al., 2015). Furthermore, Rosetta spacecraft ing (Yoshino et al., 2003; Bagdassarov et al., 2009; Ceran- data indicates a carbonaceous or enstatite chondrite sur- tola et al., 2015; Ghanbarzadeh et al., 2017). Higher oxygen face for 21 Lutetia and an average density of 3400 kg m fugacity may therefore result in more Fe-rich silicate melts (Sierks et al., 2011; Pätzold et al., 2011), higher than known with reduced (or even inverted) density contrast relative to chondrites and consistent with past compaction and partial the host rock. Lower oxygen fugacity, in contrast, may pro- melting beneath a primitive, chondritic crust (Weiss et al., duce iron-poor, buoyant melts that ascend rapidly. 2012; Neumann et al., 2013). In this study, we seek to establish regime boundaries that Based on the available evidence, most current models pro- separate primary evolution scenarios of early solar system pose a magma ocean scenario, where high melt fractions planetesimals by assessing the effects of melt segregation ( & ) dominated the thermal and chemical evolution on thermal evolution and chemical differentiation. We focus trans of planetesimal interiors. For the purposes of this study, we on the melting and transport of the major lithophile phases characterize the magma ocean scenario as a planetesimal in primitive bodies and investigate the potential for melt ac- exhibiting a fully molten interior of a well-mixed composition cumulation and heat source redistribution. We employ a and an adiabatic temperature proﬁle located beneath a thin computational model of coupled ﬂuid dynamics and thermo- ( 10 km), unmolten, chemically primitive lid with a linear chemical evolution that combines multi-component petrolog- conductive thermal proﬁle. Recent modeling studies investi- ical reactions with a two-phase magma transport model. We gating this scenario have relied either on thermal modeling quantify the leading controls on melt segregation in planetes- with parameterized melting (e.g., Hevey and Sanders, 2006; imals using theoretical considerations and numerical calcu- Elkins-Tanton et al., 2011), or on one-phase convection mod- lations of idealized planetesimal evolution. Our results show els (e.g., Golabek et al., 2014; Lichtenberg et al., 2016a, that both the magma ocean and magma sill scenarios are 2018) that capture the collective ﬂow and thermo-chemical realized within a relevant parameter space. We will focus evolution of partially molten rock or partly crystalline magma. our discussion on the latter case, where melt segregation is However, two-phase theory of partially molten systems most important. (e.g., McKenzie, 1984) suggests that silicate melts may buoyantly ascend relative to the ambient rock matrix. De- 2. Melt segregation scaling pending on the compositional and rheological properties of silicate minerals and their melts, this segregation may have To gain a leading-order understanding of silicate melt as- delayed, or even precluded, the generation of a magma cent in Al-heated planetesimals, we ﬁrst consider the char- ocean structure. Ascending melts may instead have formed acteristic time scales of melt transport in partially molten melt-rich layers beneath the primitive lid, hereafter referred bodies. In two-phase theory (McKenzie, 1984), the charac- to as magma sills (Wilson and Keil, 2017). For the purposes teristic length scale of melt migration by porous ﬂow relative of this study, we deﬁne this magma sill scenario as a body to a permeable rock matrix is given by the compaction length with radial heterogeneities of melt fraction that differ from the ﬁducial magma ocean case. This scenario implies a poten- 0 0 = ; (1) 26 c tially signiﬁcant redistribution of Al, which is a moderately 0 0 incompatible element and preferentially partitions into sili- cate melts. The transfer of the major heat source into shallow with characteristic rock viscosity , melt viscosity , melt 0 0 magma sills might then result in a transient, inverted temper- fraction , and rock permeability ature proﬁle with a thermal history distinct from the magma 2 n 0 0 ocean scenario. To date, only few studies have quantitatively k = ; (2) b (1 ) investigated this effect. These were either based on melt 0 transport models with parameterized melt ascent velocities where a is the characteristic grain size, b a geometric factor, (Moskovitz and Gaidos, 2011; Wilson and Keil, 2012; Man- and m, n powerlaw exponents. The characteristic velocity of dler and Elkins-Tanton, 2013; Neumann et al., 2013, 2014, segregating melts is 2018), or focused on metal-silicate separation (Šrámek et al., k g 2012; Ghanbarzadeh et al., 2017). 0 0 0 w = ; (3) The efﬁciency of melt transport in planetesimals depends 0 0 on various parameters. The presence of primordial volatiles with the solid-melt density contrast, and g the surface 0 0 favors rapid segregation by increasing the buoyancy and low- gravity. In primordial, homogeneous planetesimals, gravity ering the viscosity of magmas. However, if volatiles are ex- increases with increasing distance r from the center, solved before the onset of silicate melting, Fu and Elkins- Tanton (2014) argue that the segregation rate of dry melt is g(r) = g r=R ; (4) 0 P mostly controlled by the oxygen fugacity and the degree of melting. The oxygen fugacity determines (or is determined where R is the total planetesimal radius. The ﬁrst silicate by) the relative abundance of FeO and Fe-FeS in the pri- melts in sufﬁciently large planetesimals form in an adiabatic 2 zone stretching from the center, where gravity is negligible, to below the upper conductive lid (Hevey and Sanders, 2006; Lichtenberg et al., 2016a). As melting progresses, the per- meability increases and melts in shallower regions of the planetesimal, where gravity is highest, begin to segregate from the residual rock. Melt segregation can alter the chem- ical and thermal structure of Al-heated planetesimals, be- cause early-formed melts are preferentially enriched in in- compatible elements including Al. We deﬁne the characteristic melt segregation time scale, , as segr = R =w : (5) segr P 0 To achieve a substantial redistribution of heat-producing Al, the rate of melt transport must exceed the rate of melt gener- ation. We thus deﬁne the heating time scale, , of a plan- heat etesimal at a given time t after Ca,Al-rich inclusions (CAIs) as = c T =H26 (t ) ; (6) heat p 0 0 Al Figure 1: Scaling analysis of melt segregation propensity, with melt seg- with the speciﬁc heat capacity of silicates, c = 1100 J 1 1 regation number R = log ( / ). (A) At t = 0 Myr (CAI formation) seg 10 heat mt kg K , the temperature difference between accretion and and with rising melt fraction , the migration velocity increases, and so the solidus temperature, T 1100 K, and the decay power of system is more likely to become segregated. At around & 0:4–0:6, the Al per unit mass, magma ocean regime is reached and the system would be dominated by turbulent convection. (B) For ﬁxed melt fraction = 0:4 and later times (= " # Al E26 weaker radiogenic heating) the melt segregation number rises further. Al t = 0 26 Al H26 (t ) = f e : (7) 0 Al Al Al Al Here, f is the chondritic abundance of aluminum per unit compared to the internal heating by Al in planetesimals Al h i Al 5 mass (Lodders, 2003), = 5:25 10 is the canonical (e.g., Hevey and Sanders, 2006; Lichtenberg et al., 2016a, Al 26 27 2018). However, it crucially depends on the dynamic evolu- ratio of Al to Al at CAI formation (Kita et al., 2013), E = Al tion of the melt fraction, which is controlled by the ﬂuid me- 3:12 MeV = 5 10 J is the decay energy, and = Al chanics of melt transport in a deforming rock matrix, and the t =ln(2) = 1:03 Myr is the mean lifetime of Al. 1=2; Al thermo-chemical evolution of the body. In particular, Figure Using these characteristic scales, we deﬁne the non- 1 highlights the importance of considering the evolution of dimensional melt segregation number, internal heating and melt segregation over time. The scal- ! ! k g c T heat 0 0 0 p 0 0 ing analysis does not yet capture any time-dependent ef- R = log = log ; (8) seg 10 10 R H (t ) fects of interest here, which include the potential accumu- segr P 0 Al 0 lation of magma beneath a primitive lid and the redistribution which quantiﬁes the propensity of a planetesimal to undergo of the heat source by transport of melt enriched in Al. In substantial melt segregation during the internal heating by order to assess these dynamic processes we require a time- Al as a function of the key model parameters. To antic- dependent evolution model, which we introduce in the next ipate the expected melt segregation regimes in planetesi- section. mals, we calculate R for a reasonable range of melt frac- seg tions ( 2 [0; 0:4]) below the rheological transition ( 0 trans 0:4–0:6, Costa et al., 2009), formation times (t 2 [0; 4] 3. Method 5 2 Myr), grain sizes (a 2 [10 ; 10 ] m) and melt-rock density contrasts ( 2 [120; 1200] kg m ). Figure 1 shows that 3.1. Melting and heat source partitioning a growing melt fraction increases the melt segregation num- ber through its effect on permeability. Larger grain sizes and Studies of primitive meteorites (McSween et al., 1991; higher density contrasts also signiﬁcantly enhance segrega- Dunn et al., 2010) and equilibrium condensation sequence tion, but the effect of the latter is less pronounced than the calculations (Ebel and Grossman, 2000) suggest that the former. In Figure 1b, holding melt fraction ﬁxed, the heat- main rock-forming mineral phases in solar system planetesi- ing rate decreases with later formation times, which again mals were olivine (dominantly forsterite, Mg SiO ), pyroxene 2 4 serves to favor melt segregation relative to melt production. (dominantly enstatite, MgSiO ), and feldspar (dominantly From this scaling analysis, we conclude that melt segre- anorthite, CaAl Si O ). Ignoring minor contributions from 2 2 8 gation can in principle occur on a time scale that is relevant CAIs and trace minerals, feldspar represents the major host 3 Table 1: Scaling quantities, deﬁnitions and parameter values introduced in the scaling analysis, R_DMC and two-phase ﬂow models. Varying model parameters are named in the text and ﬁgures. Parameters not listed here are as given in Table 1 in Keller and Katz (2016). Parameter Symbol Unit Value Geometric factor b non-dim. 100 Melt fraction exponent n non-dim. 3 Solid fraction exponent m non-dim. 2 Melt shear viscosity Pa s 1 1 5 Thermal expansivity K 3 10 1 1 Speciﬁc heat capacity c J kg K 1100 2 1 6 Thermal diffusivity m s 1:14 10 olv olv initial mass fraction c¯ wt % 50 pxn pxn initial mass fraction c¯ wt % 35 fsp fsp initial mass fraction c¯ wt % 15 olv olv melting point T K 2050 m;0 pxn pxn melting point T K 1500 m;0 fsp fsp melting point T K 1350 m;0 Entropy gain of fusion S J K 320 olv pxn fsp 1 1 Curvature coefﬁcients r , r , r J kg K 50, 20, 10 olv pxn fsp 1 Linear P-coefﬁcients B , B , B K GPa 60, 100, 120 Rock density kg m 3200 Reference rock viscosity Pa s 10 Shear viscosity cut-off Pa s 10 min Compaction viscosity cut-off Pa s 10 min Initial temperature T K 290 init Surface temperature T K 290 space Planetesimal radius R km 60 5 2 Grain size a m [10 , 10 ] Formation time t Myr [0, 4] form Melt-rock density contrast kg m [120, 1200] 26 i phase of Al in rocky planetesimals. In addition to the tim- Curvature coefﬁcients r adjust the temperature dependence ing of accretion and size of a planetesimal (e.g., Lichten- of partition coefﬁcients. The pressure-dependent pure- berg et al., 2016a, 2018), magma genesis depends on the component melting points are parameterized as relative abundance of these phases and the concentration i i i T (P) = T + B P ; (11) m m;0 of volatiles. However, to avoid further complexity, we will consider only dry melting here, which is justiﬁed if volatile with linear slopes B . degassing during accretion is efﬁcient (e.g., Monteux et al., At a given volume-averaged bulk composition (lever rule) 2018). We therefore formulate a model for melting and melt- i i i c¯ = (1 )c + c ; (12) s ` rock partitioning of these major mineral phases. eq We employ the R_DMC method of Keller and Katz (2016) we numerically determine the equilibrium melt fraction to calculate an idealized thermodynamic equilibrium at given that satisﬁes the partitioning coefﬁcients K by combining the temperature, pressure, and bulk composition, and linear ki- lever rules with the unity sum constraint on all components netic reaction rates for a multi-component system. We deﬁne n n X i X i c¯ c¯ three compositional pseudo-components and their mass- = : (13) eq i eq eq eq i i i =K + (1 ) + (1 )K concentrations in the solid (rock), c , and liquid (melt), c . s i=1 i=1 These capture the leading-order behavior of classes of min- eq From , we then calculate the equilibrium phase composi- erals grouped by similar melting points and partitioning be- tions for solid and melt, havior: olv (olivine-like, i = 1), pxn (pyroxene-like, i = 2) c¯ and fsp (feldspar-like, i = 3). The mass fraction of all three i;eq c = ; (14a) eq i eq components must sum to unity in both phases. =K + (1 ) Using a simpliﬁed form of ideal solid solution theory c¯ i;eq c = : (14b) (Rudge et al., 2011), we determine the component partition ` eq eq i + (1 )K coefﬁcients at given P; T -conditions, Dynamic changes in pressure, temperature or bulk com- " !# i;eq c L 1 1 position over time create disequilibrium. The mass transfer K = = exp ; (9) i;eq i i of component i from solid to liquid that drives the system r T T (P) c m back towards equilibrium is assumed to occur at linear ki- i;eq i;eq which are the ratios of solid, c , to liquid, c , compo- netic rates, nent concentrations at equilibrium. The latent heat of pure- i;eq i eq i = R c c ; (15) component fusion, with a rate factor R = = , which restores equilibrium over i i 0 L = S T ; (10) m;0 a time scale at reference density . The sum of all com- is given by the entropy gain of fusion, S , and the pure- ponent reaction rates is the total melting rate = . All component melting temperatures at zero pressure, T . parameter values are given in Table 1. m;0 4 pxn fsp of melting increases with temperature, more pxn is dissolved 90 ΔT = T - T = 150 K m m,0 m,0 into the melt. Finally, in absence of melt migration, the melt ΔT +/- 50 K composition would converge to the bulk composition as the ΔT +/- 100 K system approaches complete melting. At our chosen refer- pxn ence calibration T = 1500 K, the melt initially comprises m;0 80% fsp, but becomes relatively enriched in pxn by the time the melt fraction reaches . trans 3.2. Two-phase, multi-component ﬂuid model As partial melts segregate from their residual, the interior becomes gradually depleted of fsp and thus its heat source Al. To model this dynamic process, we couple the multi- fsp B component melting model above to the two-phase reactive transport model of Keller and Katz (2016). The ﬂuid mechan- pxn ics part of the model is based on McKenzie (1984). Here, we give a brief summary of the governing equations and consti- tutive relations. The physical model is derived from statements for the con- olv servation of phase and component mass, phase momentum, and total energy. We consider the model in a Cartesian co- ordinate system with gravity pointing down the vertical co- ordinate, g = gzˆ. The governing equations are formulated olv in the geodynamic limit (liquid viscosity rock viscosity ), using the extended Boussinesq approximation (densities = = taken equal and constant except when multi- s ` 0 plying gravity). The ﬂuid mechanics governing equations are pxn " # fsp r rv + (rv ) Ir v rP + g = rP ; s s s ` C 1500 1600 1700 1800 (16a) Temperature [K] r v = r v ; s D Figure 2: Temperature dependence of melt fraction (A), melt (B) and rock (16b) (C) composition, with varying melting point difference between fsp and pxn, pxn fsp T = T - T , which changes the partitioning coefﬁcient of fsp and m m v = (v v ) = (rP + g) ; m;0 D ` s ` the composition of the earliest melts. If the fsp melting point is close to the one of pxn, the initial melt composition is close to a 50–50 mixture. (16c) For higher melting point difference, the ﬁrst melts are dominated by fsp, P = (1 )(P P ) = r v : and thus the heating component ( Al) of the system preferentially follows C s ` s the dynamics of the earliest melts. When the absolute temperature of the (16d) system further rises and approaches the olv melting point, the composition converges towards the initial pure solid setting. They are posed in four independent variables, the dynamic pressures, P , P , and velocities v , v , of the solid and liq- s ` s ` uid phases. Two dependent variables, the Darcy segregation ﬂux, v , and the compaction pressure, P , express the me- We ﬁrst consider the aluminum partitioning behavior in a D C chanical interactions between the phases. If these vanish, closed, isobaric (P = 0) system under increasing tempera- the equations become identical with the Stokes system. As- ture. For consistency with earlier work (Tkalcec et al., 2013; suming a diffusion creep rheology with melt-weakening, the Golabek et al., 2014; Lichtenberg et al., 2016a; Monteux shear viscosity of the rock matrix is given by et al., 2018; Lichtenberg et al., 2018; Hunt et al., 2018), we fsp choose the lowest component melting point, T = 1350 K, m;0 E + PV a a to conform with previous estimates of the silicate solidus. To = A a exp ; (17) RT test different partitioning behaviors, we vary the relative tem- perature difference between the melting points for the Al- with prefactor A , activation energy, E , and activation vol- 0 a pxn fsp hosting fsp and the pxn components, T = T T , with m ume V as in Hirth and Kohlstedt (2003) and Mei et al. m;0 m;0 a pxn T 2 [1400, 1600] K. The resulting rock and melt composi- (2002). R is the universal gas constant, and 30 the melt- m;0 tions in the ternary system olv-pxn-fsp are shown in Figure weakening factor. Permeability is set by the Kozeny-Carman 2. The larger T between fsp and pxn, the more incipient relation (Eq. 2), with powerlaw exponents n = 3, m = 2 for melt will be enriched in fsp, and thus Al. As the degree the melt and solid fractions, respectively. The compaction Rock composition [wt %] Melt composition [wt %] Melt fraction [wt %] viscosity is set to = r =, with r 1 the shear to com- interior evolution of planetesimals of R & 50 km as these are paction viscosity ratio. dominated by a relatively large adiabatic interior and a thin To these equations we add thermo-chemical evolution (. 10 km) conductive lid, whereas the relative dimensions of equations for temperature, T (assuming thermal equilibrium these domains vary signiﬁcantly for planetesimal radii . 50 between phases), melt fraction, , and pseudo-component km (cf. Castillo-Rogez and Young, 2017; Lichtenberg et al., i i concentrations in the solid and liquid phases, c and c : 2018; Hunt et al., 2018). The computational domain includes 500 grid cells for a spatial resolution of 120 m. The sur- X i i DT L H (t) T face boundary is T = 290 K, similar to the temperature Al 2 space = r T + + + gw ; in the protoplanetary disk inside of the water snowline, while Dt c c c c 0 p 0 p p p i=1 the center boundary is insulating, @T=@z j = 0. We as- z=0 (18a) sume the accreted body is initially at ambient temperature, T = T . As noted above, gravity decreases linearly from = (1 )r v + = ; (18b) init space s 0 Dt the surface gravity down to zero at the center. i i i D c c ` ` Magma and rock composition are modeled in the three- = ; (18c) Dt 0 component compositional space of olv, pxn and fsp (Section i i i 3.1). We use component melting points as in Figure 2 (solid D c c s s = : (18d) lines). The solid-melt density contrast is varied as 2 [200, Dt (1 ) 1200] kg m between runs to reﬂect FeO content and thus Temperature evolves due to advection, thermal diffusion, density of the melt reﬂecting the planetesimal’s oxygen fu- latent heat exchange of reactions, heating by viscous dis- gacity (Fu and Elkins-Tanton, 2014; Wilson and Keil, 2017). sipation and internal heating H26 (t), and adiabatic de- Al Grain size a , which controls both the permeability and rock compression. Melt fraction evolves due to advection, com- viscosity, is held constant during calculations; we consider paction and reactions, and composition by advection and re- values a 2 [10 m, 1 cm], from chondrite matrix-like dust action. The material derivative of the two-phase mixture is to pallasite-like crystal sizes (cf. Figure 1). Heating is in- ¯ 26 D=Dt = (1 )D =Dt + D =Dt, with D =Dt = @=@t + v r, s ` s s duced solely by Al in the fsp component, whose redistribu- and D =Dt = @=@t + v r. is the thermal diffusivity, ` ` tion hence affects the local heating rate. The initial heating 7 1 w = [(1 )v + v ] zˆ the vertical bulk speed, and the z s ` rate is varied from H26 (0) 2 [1:52; 0:19] 10 W kg , re- Al thermal expansivity. ﬂecting planetesimal formation times t 2 [0; 3] t = form 1=2; Al The governing equations for the ﬂuid mechanics [0; 2:2] Myr after CAIs. To limit model complexity, we ignore (Eqs. 16a–16d) and thermo-chemistry (Eqs. 18a–18d) are the potential heat contribution from Fe (see, e.g., Tang and discretized using the ﬁnite difference method on a rectangu- Dauphas, 2015; Lichtenberg et al., 2016b). lar, staggered grid and solved by two coupled Newton-Krylov solvers. The simulation software uses parallel data struc- 4. Results tures and solvers provided by PETSC (Balay et al., 1997; Katz et al., 2007). Nonlinearities between the ﬂuid mechan- 4.1. Parameter study ics and thermo-chemical sub-problems are resolved using We ﬁnd that model results across the tested parameter a Picard ﬁxed-point iterative scheme. During every itera- range fall into three qualitative categories. Some models tion, Equation 13 is solved using Newton’s method to up- show no substantial melting or segregation; we will not fur- date the equilibrium melt fraction. The adopted model is ther discuss these undifferentiated models here. The time strictly valid only for melt transport by porous ﬂow below the evolution of the latter two categories shown in Figure 3 is disaggregation threshold. However, we cannot avoid mod- generally the same: rapid initial heating leads to substan- els producing regions with higher melt content. To ensure tial melt production, followed by some degree of segregation, that the equations do not produce numerically unstable so- before ending with slow cooling from the top down. One cat- lutions in these regions, we apply lower viscosity cut-offs to egory of models, which we identify as the magma ocean end- the shear viscosity ( = + ) and compaction viscosity num min member (red), evolve to where most of the interior is above ( = + =(1 )). The effect of this regularization is to num min the disaggregation threshold, whereas the other, the magma dampen the segregation velocity and compaction pressure at sill end-member (blue), result in a melt-depleted interior be- elevated . As a result, our model will produce stable solu- neath melt-rich sills (Fig. 3A–E). The latter clearly shows a tions above the rheological transition, but will underestimate thermal inversion (Fig. 3F–J) related to fsp-enrichment in the chemical mixing and heat transport by rapid convection in magma sills beneath the lid (Fig. 3K–O). the crystal-bearing suspensions that characterize this limit. The scaling analysis above predicts that grain size, a , and density contrast, , are pertinent controls on melt segrega- 3.3. Model setup and parameter space tion. Figure 4 shows the results of a detailed study of that We model magma genesis and transport along a 1D parameter space. To quantitatively analyze the results, we Cartesian column from the center to the surface of an initially introduce three metrics measuring the degree of homogeneous and isothermal body of 60 km radius. Plan- etesimals of this size are qualitatively representative of the melt segregation: = – , max ctr 6 1.5 Myr 2.0 Myr 4.0 Myr 7.0 Myr 10.0 Myr 0.9 0.7 0.5 0.3 0.1 A A B B C C D D E E 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 [wt] 0.9 0.7 0.5 0.3 0.1 F F G G H H I I J J 1450 1700 1950 1450 1700 1950 1450 1700 1950 1450 1700 1950 1450 1700 1950 [K] 0.9 0.7 0.5 0.3 0.1 K K L L M M N N O O 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 0.2 0.5 0.8 [wt] Figure 3: Time evolution of melt fraction (A–E), temperature (F–J) and fsp bulk composition (K–O) for two end-member models with radii R = 60 km and 3 4 formation time t = 1.25 t = 0.90 Myr after CAIs. Red lines show a magma ocean model with = 200 kg m and a = 10 m, blue lines a form 26 0 1=2; Al 3 3 magma sill model with = 800 kg m and a = 10 m. Upon progressive heating, the magma sill model builds up melt accumulates below the cold upper lid, depleting the center of the planetesimal of silicate melts. High melt fractions > := 0.5 (yellow areas) are only reached in the sub-lid sills. fsp crit enrichment in the sill structure leads to a temperature inversion of 400 K at peak melting. The magma ocean model, in contrast, shows a near-isothermal internal temperature and thus constant melt fraction structure in the interior. The fsp component shows notable deviations from the initial bulk composition only after t 2 Myr, when most of the Al is already decayed. A video showing the time evolution of the major thermo-chemical parameters and composition is linked to in the Supplementary Materials. temperature inversion: T = T – T , 4.2. Silicate differentiation max ctr fsp fsp fsp differentiation: c = c – c . fsp max ctr Here, () denotes the maximum value in the computational As the planetesimals cool and crystallize, magma sill end- max domain and () denotes the value at the base of the do- member cases exhibit a silicate differentiation trend towards ctr main, i.e., the planetesimal center. With these metrics, we compositional layering. In Figure 5, we compare a repre- quantify the deviation from an interior structure with near- sentative magma ocean with a magma sill model outcome. constant melt fraction, temperature, and bulk composition, Magma ocean models evolve towards a uniformly molten as it would be expected in the absence of signiﬁcant segre- interior, with melt fractions above the rheological transition gation. > and cool relatively unsegregated, such that the trans Figure 4 shows the three metrics across a range of a 2 bulk concentrations of the fsp, pxn and olv components re- 5 2 3 [10 , 10 ] m and 2 [200, 1200] kg m . We ﬁnd that main similar to the initial composition (Figure 5A). However, grain size strongly controls the interior evolution of the plan- in magma sill models, the melt-rich layers above a low-melt- etesimals. If grain size remains below a < 10 m, melt den- fraction interior crystallize into distinct layers enriched in fsp, sity contrasts of < 500 kg m are not sufﬁcient to drive signif- pxn and olv. This stratiﬁcation reﬂects the component melt- fsp pxn 3 olv icant segregation. For density contrasts > 1000 kg m and ing points (T < T < T , Figure 5B). The melt-depleted m m for larger grain sizes, melt segregation becomes signiﬁcant, central parts of the planetesimal are strongly depleted in fsp as evidenced by a step-increase in each of the three met- and somewhat less in pxn. In general, such compositional rics. However, we ﬁnd that since initially fsp-enriched melts layering forms during the cooling stage and therefore does migrate on time scales comparable or longer than t 26 ( not cause substantial temperature inversion. The densities 1=2; Al 0.72 Myr), the temperature inversion effect remains minor of the minerals represented by the pseudo-components sug- throughout (Figure 4, panels B & F). gest that such layering would be dynamically stable. Radius / = const. = const. Radial profiles 2 4 0.9 0.9 0.9 0.7 0.7 0.7 0.5 0.5 0.5 0.3 0.3 0.3 0.1 1 1 0.1 3 0.1 A B C D -5 -4 -3 -2 200 700 1200 0.2 0.5 0.8 0.2 0.5 0.8 log [m] [kg/m ] 0.9 0.9 0.7 0.7 90 0.5 0.5 1 3 0.3 0.3 1 2 0.1 0.1 E F G H -5 -4 -3 -2 200 700 1200 1400 1600 1400 1600 log [m] [K] [kg/m ] 0.9 0.9 0.7 0.7 0.7 0.5 0.5 0.5 0.3 0.3 0.3 0.1 1 1 3 0.1 0.1 I J K L -5 -4 -3 -2 200 700 1200 0.1 0.35 0.6 0.1 0.35 0.6 log [m] [kg/m ] [wt] Figure 4: Parameter study of the inﬂuence of grain size a and density contrast on melt segregation , temperature inversion T , and compositional stratiﬁcation c , for planetesimals with R = 60 km, t = 1.5 t . Panels (A,E,I) show the metric deviation for constant density contrast = 600 kg fsp P form 26 1=2; Al 3 4 3 3 m and varying grain size a , indicating a steep gradient between grain sizes of 10 m and 10 m. For these two values ﬁxed (blue: a = 10 m, green: 0 0 a = 10 m), panels (B,F,J) display the metric deviations for varying . Here, variations in density contrast are outweighed by those in grain size. Models 3 4 with a = 10 m feature notable melt segregation, temperature inversions, and compositional differentiation. Models with a = 10 m only do so towards the 0 0 high end of density contrasts, & 700 kg m . Panels (C/D, G/H, K/L) show the radial proﬁles for the end-member models of the variations from (A/B, E/F, I/J). In general, variations in grain size a outweigh effects from increasing density contrast . Magma sill structures only form for grain sizes a > 10 m. 0 0 Magma ocean Magma sill 4.3. Magma dynamics regimes A B Figure 6 shows the regimes of melt segregation and com- positional stratiﬁcation for different formation times, t , and form 0.8 melt segregation numbers, R . We quantify the boundaries seg fsp-enriched separating the three characteristic regimes as follows. pxn-enriched 0.6 (I) Magma ocean regime: > := 0.5, where the olv-enriched ctr crit whole interior melts to above the disaggregation thresh- old. 0.4 (II) Magma sill regime: > , where melt segregation crit fsp-depleted 0.2 generates a melt-rich layer beneath the lid and depletes the interior of melt. fsp (III) Undifferentiated regime: c < c = 0.15, where 0.2 0.5 0.8 0.2 0.5 0.8 fsp bulk;0 melting and melt segregation do not redistribute a sub- [wt] after 10 Myr stantial amount of fsp. Figure 5: Compositional stratiﬁcation after cooling and crystallization of In addition to these segregation and differentiation criteria, magma beneath the primordial lid for magma ocean (A) and magma sill - we show which models are most affected by substantial tem- type (B) models. Magma sill cases with intermediate temperatures and thus high concentrations of fsp in the upper layers produce this signature. perature inversions, T > 250 K. These inversions occur both for magma sill and magma ocean models and reﬂect rapid magma ascent on time scales shorter than t 26 . We 1=2; Al fsp [K] [wt] pxn olv Compositional stratification Temperature Melt fsp anomaly segregation pxn olv / / / Figure 6: Evolution of silicate melt segregation with formation time t versus melt segregation number R at a reference melt fraction of = 0:1. form seg 0 Colormap values are plotted for the time of peak melting for each model (circles). We identify three primary melt dynamics regimes. (I) Magma ocean models, where melting occurs more rapidly than melt migration, feature high melt fractions above the rheological transition in their center, > := 0.5. ctr crit Magma ocean-type evolution is preferred for early t and low R , i.e., small a and . (II) Magma sill models feature efﬁcient melt segregation, with form seg 0 additional compositional stratiﬁcation towards cooling-down of the planetesimals (cf. Figure 5). Magma sill -type evolution is preferred for intermediate t form 0.5–1.75 t 26 , and high R . (III) Undifferentiated models never show melt fractions > , and never experience substantial compositional seg ctr crit 1=2; Al redistribution. They are preferred for late formation times t & 2.0 t . In addition to these three regimes, we show the region of increasing form 26 1=2; Al temperature inversion in yellow. ﬁnd that the magma sill regime generally occurs at higher Sato, 1984; Fu and Elkins-Tanton, 2014) at a 10 m. segregation numbers—at larger grain sizes or elevated den- Moreover, we ﬁnd that even in the case of signiﬁcant magma sity contrasts—and formation times less than 1.5 Myr after redistribution into shallow sills, the segregation time scale is CAIs, with a peak at around 1 Myr after CAIs. Very early for- comparable to the evolutionary time scale of the protoplane- mation time, t 1 Myr, or lower melt segregation number, tary disk and thus the accretion time scale. Therefore, scal- form R . -1.5, favor magma ocean-type evolution. Formation ing analysis alone (see Section 2), does not adequately cap- seg later than 1.5 Myrs after CAIs results in limited melting and ture the time-dependent competition between melting, parti- undifferentiated planetesimals. Temperature inversions oc- tioning, and melt transport. Because of the protracted onset cur for t . 1 Myr and R & -1.5. of melt ascent during heat-up, Al-hosting phases release form seg at least parts of the Al decay energy in the deeper region of planetesimals, and substantial temperature inversions on 5. Discussion the order of a few hundred K are only observed for extreme cases with R & -1, or early formation times t . 1 Myr seg form 5.1. Parametric controls on magma segregation with R & -2. seg The models above present thermo-chemically coupled At t = 26 1 Myr after CAIs the models show a clus- form Al two-phase ﬂow calculations that resolve the partitioning of tering of magma sill cases. This peak is due to the competi- the major rock-forming components and their transport by tion between heating and segregation, as introduced in Sec- magma. Using this method, we show that magma segre- tion 2. For formation times earlier than 1 Myr, heating, and gation in planetesimals formed within 2 Myr after CAIs was thus melting, is dominant, and new melt is generated faster potentially signiﬁcant if magma ascent was rapid with re- than transported away. Around t 1 Myr, melt upwelling form spect to the rate of melt production controlled by Al de- velocities become fast enough to exceed melt generation. Fi- cay. Expressed in terms of the melt segregation number R , seg nally, at later times (t & 1.5 Myr after CAIs or 2 t 26 , form 1=2; Al magma sill models ( > 0:5) were produced for R -1.5. seg respectively) melt production remains low and minimal to no This requires that (i) the average grain size in these planetes- segregation occurs. imals was comparable to or larger than chondrules, on the order of a 10 m, or (ii) in reducing environments (iron- Melt–rock density contrast is thought to be controlled by wÃijstite IW-2.5, or 1200 kg m , respectively, Brett and the oxidation state of the body. The chemical composition of 9 accreting planetesimals is inherited from the oxidation state ﬂow regime in the magma, leading to a variety of possible in the midplane of the protoplanetary disk. Hence, the time scenarios (Solomatov, 2015). and place of formation is expected to inﬂuence the subse- As a comparison, grain sizes of meteorite classes span quent interior evolution as it relates to the effects of melt seg- a wide range, from m-sized dust to pallasite-like, cm- regation. For example, planetesimals accreted towards the sized phenocrysts (Hutchison, 2004). Chondrites, gener- inner disk may feature lower oxygen fugacities and therefore ally regarded to be the most pristine rock samples from higher , compared to planetesimals accreted further out the early solar system, display a bimodal size distribution, (Rubie et al., 2015). In our models, we ﬁnd that the effect split between m-sized dust (‘matrix’) and chondrules, which 4 3 of density contrast on magma segregation is of secondary show characteristic sizes of 10 m to 10 m (Friedrich importance. If grain sizes were too small to allow for a suf- et al., 2015), with drastically differing textures and min- ﬁcient permeability, variations in density contrast would not eral grain sizes. The ratio of chondrules to matrix varies have led to signiﬁcantly different outcomes. This ﬁnding is greatly between different chondrite groups, resulting in com- in contrast to previous studies (Moskovitz and Gaidos, 2011; plex mixtures and grain size distributions. Meteorites that Wilson and Keil, 2012; Šrámek et al., 2012; Fu and Elkins- likely underwent partial melting (‘primitive achondrites’), like Tanton, 2014; Wilson and Keil, 2017; Fu et al., 2017), which Acapulcoite-Lodranites, Winonites and Brachinites, display did not consider magnitude variations in grain size, or relied grain sizes around 10 m (Hutchison, 2004; Wilson and on grain growth by pure Ostwald ripening without taking into Keil, 2017). Ureilites, Aubrites, and Angrites, which come account mechanisms of growth inhibition and grain destruc- from bodies with more extensive, or even body-wide, silicate tion (Neumann et al., 2013, 2014). melting feature larger grains, on the order of 10 m. How- ever, these textures are the end result of million-year long 5.2. Implications for the role of grain size evolutionary processes, and may have undergone repeated destruction and reaccretion cycles that reset their thermal A main conclusion from our models is that the mean grain histories and textures (Lichtenberg et al., 2018). Therefore, size a is the dominant parameter controlling the magma based on the grain sizes observed in meteorites, we are transport rate inside planetesimals. There are two main strongly limited in assessing the grain size evolution in plan- mechanisms that affect grain size during planetesimal ac- etesimal interiors at the time of their ﬁrst melting. cretion and early evolution. First, grain sizes may differ de- Interpretation of our results in this context suggests that pending on the orbital location and mineralogical composi- the planetesimal interior evolution and the redistribution of tion (van Boekel et al., 2004) of the dust that agglomer- chemical and isotopic heterogeneities during planetary ac- ates to form the planetesimals. Coagulated dust aggre- cretion can be inﬂuenced by the planetesimal formation gates on the order of a 10 m are in the critical size mechanism, its accretion location, and the local composi- regime for fast radial drift towards the protostar, depending tional distribution of grains in the protoplanetary disk. Fur- on orbital location and evolutionary state of the protoplan- ther studies of planetesimal formation and mineralogical in- etary disk (Weidenschilling, 1977). However, these sizes ventory are needed to link their formation processes to their may facilitate planetesimal formation from local dust-gas in- subsequent dynamic evolution. The local grain-size distribu- teractions (Johansen et al., 2015; Birnstiel et al., 2016) and tion within planetesimals may be evolving during rapid grav- can trigger planetesimal formation at various orbital sepa- itational collapse (Wahlberg Jansson and Johansen, 2017) rations and times (Dra ¸ zk ˙ owska and Dullemond, 2018) with or more gradual growth (Kataoka et al., 2013). Also, ongoing varying dust particle distributions within newly-formed bod- accretion, for instance due to subsequent pebble accretion ies (Wahlberg Jansson and Johansen, 2017). Second, grain (Visser and Ormel, 2016), may inﬂuence whether magma sizes may evolve during the heating and melting in the plan- can reach the surface through fractures. Future experimental etesimal interior after accretion. In this process, the grain studies on grain size evolution of partially molten aggregates size evolves due to competing mechanisms for growth and spanning the meteoritic compositional space will be needed destruction (Rozel et al., 2011; Bercovici and Ricard, 2016). to advance our understanding of melt migration in early solar Grain growth by either normal grain growth before the ﬁrst system planetesimals. melts arise, or Ostwald ripening during partial melting in a solid–liquid aggregate, is driven by the reduction of interfa- 5.3. Implications for chemical differentiation cial energies, and competes with size reduction due to the Recently, it was shown that differentiation by percolation presence of secondary solid phases (Hiraga et al., 2010) of Fe,Ni-FeS liquids in primordial planetesimals may not be and mechanical work due to solid-state deformation in plan- complete and that at least some material remains trapped in etesimal interiors (Tkalcec et al., 2013). During later stages, the rock matrix (Bagdassarov et al., 2009; Cerantola et al., when the melt fraction reaches or exceeds the disaggrega- 2015). But once silicate melting has reached the disaggre- tion threshold, grain sizes are governed by the convective gation threshold, the remaining metal droplets will rain out rapidly towards the forming core (Lichtenberg et al., 2018). Note that in astrophysics literature grain size usually denotes the Therefore, even though the models in this study do not in- characteristic size of porous dust aggregates embedded in the disk ﬂow, clude a metal phase, they allow leading order predictions re- whereas we here use the term grain size, a , for the characteristic size of mineral grains in a granular, polymineralic rock aggregate. garding core formation. 10 In the magma ocean case, we expect core formation to be 5.4. Limitations & future directions rapid, with nearly complete loss of metals to the core. In the One of the limitations of our model is the use of a 1D case of a magma sill, a two-step process may occur. First, Cartesian geometry, which introduces systematic errors as a small proto-core may form from incomplete percolation. compared to a 1D spherical geometry assuming radial sym- Then, after the formation of the sill structure, the remain- metry. Among these errors, our model over-predicts the ing metal within this region may rain out and accumulate at heating-to-cooling ratio of planetesimals. Because we focus the interface between the melt-depleted deep interior and the on planetesimals of 60 km radius that show a nearly isother- magma sill zone. This emerging metal pool will either per- mal evolution in the deep interior (cf. Figure 3, Castillo- colate downwards or form diapirs sinking through the weak- Rogez and Young, 2017; Hunt et al., 2018), the heat-up ened partially molten rock. The thermo-mechanical evolu- phase is consistent with a radial model. However, geometric tion of such a two-step process needs to be tested by taking errors result in an under-prediction of the rate of heat loss into account metal phases in self-consistent multi-phase ﬂow through the surface. Similar problems apply to melt segre- models, in order to make robust predictions that can be com- gation velocities and melting rates. For example, the liquid pared to laboratory studies of the core formation process. mass conservation equation, in spherical coordinates with Neumann et al. (2018) recently suggested a multi-stage core radial symmetry, contains the term 2w =r, which is neglected formation scenario for the IVB parent body, which is quali- here (a similar term appears in the solid mass conserva- tatively consistent with the magma sill regime we propose tion equation). For radially outward melt migration, this term based on our models. has the effect of reducing @=@t; its absence therefore leads to an over-prediction of melt fraction. To quantify the er- The limited melt segregation in our undifferentiated mod- ror introduced by neglecting this term we compared its size, els may offer an explanation for the absence of remnant dif- computed a posteriori, with the ﬂux divergence @w =@z that ferentiated materials in the asteroid belt (Weiss et al., 2012; appears in the same equation. We made this comparison Mandler and Elkins-Tanton, 2013). Conversely, chemical across models of the magma ocean and magma sill regimes. stratiﬁcation resulting from melt segregation may help to ex- The results indicate that the geometrical error can reach a plain the paucity of olivine-rich deposits on 4 Vesta’s surface magnitude comparable to the divergence. However, through- (Clenet et al., 2014; Consolmagno et al., 2015; Raymond out most of the domain and once magma sill structures are et al., 2017). Furthermore, the magma sill models and the established, the term is negligible. resulting chemical stratiﬁcation we describe here are consis- In this work we consider a diffusion creep rheology only. A tent with earlier work by Neumann et al. (2013, 2014, “shal- more realistic rheology would be a composite of both diffu- low magma ocean”) and Mizzon (2015, “completely liquid sion and dislocation creep (e.g., Bercovici and Ricard, 2016; layer”), predicting a subsurface layer of accumulated silicate Mulyukova and Bercovici, 2018). At the highest temper- melts below a cold lid resulting from efﬁcient melt segrega- atures and grain sizes tested here, model behavior would tion (cf. discussion in Raymond et al., 2017; Neumann et al., likely fall into the dislocation creep regime, where the matrix 2018). viscosity becomes independent of grain size. In this case, the grain size sensitivity of the compaction length (Equation Finally, our results indicate that the importance of melt 1) would decrease. segregation in planetesimal interiors varied substantially and Furthermore, we use a constant melt viscosity of 1 Pa s affected the redistribution of heat-producing elements, such (Table 1). However, silicate melt viscosity varies with tem- as Al, during melt ascent. In the case of our magma sill perature and composition (e.g., Moskovitz and Gaidos, 2011; end-member models, we also expect other incompatible el- Mizzon, 2015). For the compositional space explored here, a ements to be preferentially segregated to shallow layers of temperature dependent viscosity of 1–100 Pa s may be con- a planetesimal. The crustal stripping paradigm of plane- sidered realistic. As we are interested in the consequences tary accretion assumes that frequent impacts during plan- of melt segregation in planetesimals here, we chose a rea- etary growth eroded and redistributed signiﬁcant amounts sonable lower limit. In addition, because the ratio of per- of shallow materials between planetesimals. The strongly meability to melt viscosity controls the upwelling timescales variable degree of melt segregation, and the resulting range (grain size squared versus linear melt viscosity, Equation 1), of variably differentiated major, trace, and isotopic composi- the order of magnitude variability in grain size outweighs po- tions of shallow planetesimal layers could result in compo- tential variations in melt viscosity. sitional differences between forming planets and chondritic meteorites (e.g., Carter et al., 2018). Our simulations indi- As a consequence of our limiting assumptions on the ge- cate that magma ascent governing material redistribution to ometry, rheology, and melt viscosity, the extent of the magma the planetesimal crust occurs on a Myr time scale, i.e., sill regime in Figure 6 may be overestimated. More efﬁ- comparable to the collisional time scale of planetary accre- cient cooling in a spherical geometry would lead to decreas- tion. This suggests that models quantifying compositional ing melt fractions and thus lower migration speeds (Figure deviations between planets and chondrites should take into 1), as would a stronger heat ﬂux from turbulent convection account the evolution of planetesimal interiors during plane- above the disaggregation threshold (Solomatov, 2015), and tary accretion. the geometric inﬂuence in spherical or higher-dimensional 11 geometry. A weaker dependence on grain size in the models the total retention of magma on the inside of planetesimals. with the highest temperatures and largest grain sizes would However, this does not affect our conclusions or any con- also decrease migration speeds. Therefore, while our mod- straints on temperature inversions of planetesimals unless els constrain the possible phase space of melt migration in this transport is faster or comparable to the magma segrega- early solar system planetesimals, more complex models not tion in the interior. bound to the above limitations would result in a reduced sta- bility ﬁeld where melt segregation becomes dominant than 6. Summary & conclusions shown in Figure 6. Our models do not feature a metal phase that would al- In this study we investigated magma genesis and redis- low the direct resolution of metal percolation (Ghanbarzadeh tribution in planetesimals during and shortly after the pro- et al., 2017; Neumann et al., 2018), and therefore our re- toplanetary disk phase. Using scaling relations we demon- sults only allow for qualitative inferences about possible core strated that the interior evolution of planetesimals sensitively formation scenarios. However, the complexities of textural depends on a variety of model parameters, with the grain equilibrium phase topologies are not yet fully understood (cf. size exerting the primary control on melt segregation. Based Rudge, 2018), in particular for systems with several immis- on average chondritic abundances of the most common min- cible liquid phases. For example, the wetting angles formed eral groups in meteorites, we calculated the composition for between metal liquids with silicate minerals in the presence rock–melt aggregates comprising idealized components us- of a silicate melt remains unclear, leaving open the debate ing a reactive, multi-component melting model. We quanti- around a possible percolation threshold for metal liquids at ﬁed the effects on Al partitioning and magma ascent with low melt fractions (Cerantola et al., 2015). Further work a coupled, two-phase ﬂow model. We deﬁned the melt seg- needs to be undertaken to better quantify these effects. regation number R as the ratio between the heating and seg Our choice of melting model in the form of a ternary ideal melt transport time scales in a planetesimal, which estab- solution limits the degree to which the model may represent lishes the leading order parametric control on the propensity natural melting behavior. For example, our model does not for magma redistribution during the heating stage of plan- reproduce the eutectic behavior expected for silicate compo- etesimal evolution. We predicted that the primary two con- sitions as the ones considered here (see discussion in Keller trols are the melt–rock density contrast and the mineral and Katz, 2016), nor does it include volatiles and incompat- grain size a . ible elements producing low-degree, incompatible-enriched Investigating the relative importance of model parameters melts at temperatures below the anhydrous solidus. Using for the evolution of planetesimals, we categorize model out- a more consistent petrological model that takes into account comes in terms of their melt segregation numbers R and seg the non-ideal thermodynamics of the full range of major ele- their formation times t . Using this scheme, we ﬁnd: form ments and mineral phases would likely lead to more complex relations between heating, melt production, and element par- Planetesimal melt migration behavior can be classiﬁed titioning (Keller and Suckale, 2018). Signiﬁcant differences in three general regimes: in aluminum partitioning, which is the focus here, are likely conﬁned to the onset of melting, where low-degree, enriched – The magma ocean regime with global interior melts may have important control on geochemical evolution. magma oceans, where rapid melting overwhelms Moreover, our dry models neglect volatile-driven eruptions, melt segregation. which were previously discussed as a catalyst for upward – The magma sill regime, where global interior transport via explosive volcanism (Fu et al., 2017; Wilson magma oceans are prevented by rapid magma and Keil, 2017). If substantial volatile quantities could be transport that concentrates melt in sills beneath retained, this mechanism would decrease the retention of the cool lid. magma in the upper layers of planetesimals and potential chemical stratiﬁcation upon crystallization of the silicate ma- – The undifferentiated regime with a low degree terial. However, the relatively low pressure at the planetes- of melting, minor melt segregation, and therefore imal sizes we consider disfavors a high volatile solubility in chemically primordial and largely undifferentiated silicate magma, and degassing would therefore be expected interiors. to be nearly complete in the earliest stages of melting and segregation (Monteux et al., 2018). Magma sill models crystallize to a compositionally strat- Finally, the melt in our magma ocean and magma sill mod- iﬁed structure, with shallow depth layers dominantly en- els cannot breach the cold surface layers, as the temperature riched in the low melting point components, feldspar and is too low for melt to exist in the porous medium. Our sim- pyroxene, and a paucity of high melting point compo- ulations cannot resolve potential melt transport through the nents, such as olivine (cf. Ghosh and McSween, 1998; upper lid by fracturing (Keller et al., 2013) or gravitational in- Mizzon, 2015). The crystallization sequence, and thus stability in the layered structure in Figure 5, which may bury the compositional stratiﬁcation, however, may be af- the primitive lid (Wilson and Keil, 2012), and lead to efﬁcient fected by convective motions beyond the disaggregation heat loss and magmatic resurfacing. This would decrease limit, which we do not model here. 12 Magma ocean and magma sill models show tempera- constructive comments by David Bercovici and an anony- ture inversions for high R and early t , where the mous referee, which helped to improve the manuscript. TL seg form temperatures in the shallow- to mid-mantle are higher was supported by ETH Zürich Research Grant ETH-17 13-1 than at the center of the planetesimal. These thermal in- and acknowledges partial ﬁnancial support from a MERAC versions, however, are restricted to formation times t travel grant of the Swiss Society for Astrophysics and As- form . 1 Myr and R & 1:5 for T 250 K. Therefore, the tronomy and from the National Center for Competence in seg majority of planetesimals likely underwent thermal evo- Research PlanetS supported by the Swiss National Science lutionary scenarios that can be qualitatively reproduced Foundation. TK acknowledges ﬁnancial support from Jenny with models that do not take into account melt segre- Suckale, Stanford University. TK and RFK received fund- gation and Al partitioning, depending on the level of ing from the European Research Council under the Euro- detail that needs to be assessed. pean Union’s Seventh Framework Programme (FP7/2007– 2013)/ERC grant agreement number 279925. The numeri- The magma sill regime can be achieved depending on cal simulations in this work were performed on the EULER a combination of a few key parameters: computing cluster of ETH Zürich, and were analyzed using the open source software environments MATPLOTLIB and – The formation time t controls the total amount form SEABORN . of energy available from Al and restricts the regime for melt segregation to 0.5 t 26 . t 1=2; Al form . 1.75 t 26 , with a peak at 1 Myr after CAIs, Supplementary Materials 1=2; Al when the rate of melt transport dominates over the A simulation video associated with this article can be found generation of new melts. attached to the publication and at this URL. The video shows – The grain size a controls the rate of melt trans- a comparison between magma ocean and magma sill end- port and thus whether a planetesimal with sufﬁ- member models. Magma ocean stages are indicated in yel- cient heating evolves towards a melt segregated low. The H3 annotation describes the heating value below structure. Below a . 0.1 mm no segregation which radioactive heating from Al is inefﬁcient. Shown are occurs; above a & 1 mm segregated structures various parameters for both models, from left to right: melt form. fraction , temperature T , radiogenic heating H, melt/liquid – The solid–melt density contrast is of secondary upwelling velocity w , composition fraction of liquid c , and liq liq importance, but can enhance melt segregation in composition fraction of solid c . The composition is broken- sol the regime transition from a = 0.1 mm to 1 mm. down into the deﬁned pseudo-components fsp, pxn, and olv. The model is started (= planetesimal formation time) at 0.9 To conclude, in this paper we have advanced the techni- Myr after CAIs. cal capabilities to simulate multi-phase and multi-component planetesimal evolution, gaining insights into features not ac- References cessible to single-phase ﬂuid dynamics models. However, unraveling more detailed evolutionary regimes of planetes- A’Hearn, M. F., 2017. The future of NASA’s missions. Nat. Astron. 1, 0095. imals will require a time-dependent treatment that includes Bagdassarov, N., Golabek, G. J., Solferino, G., Schmidt, M. W., 2009. Con- straints on the Fe-S melt connectivity in mantle silicates from electrical metal and volatile phases, which shape the structure and impedance measurements. Phys. Earth Planet. Inter. 177, 139–146. subsequent evolution of these bodies (Keller and Suckale, Balay, S., Gropp, W. D., McInnes, L. C., Smith, B. F., 1997. Efﬁcient man- 2018). Further work is required to understand planetesi- agement of parallelism in object oriented numerical software libraries. mal evolution and its connection to the meteoritic record and In: Arge, E., Bruaset, A. M., Langtangen, H. P. (Eds.), Modern Software Tools in Scientiﬁc Computing. Birkhäuser Press, pp. 163–202. rocky planet formation (Lichtenberg et al., 2018). In partic- Bercovici, D., Ricard, Y., Apr. 2016. Grain-damage hysteresis and plate tec- ular, robust scaling laws for the evolution of grain sizes in tonic states. Physics of the Earth and Planetary Interiors 253, 31–47. partially molten and heated systems relevant for planetesi- Birnstiel, T., Fang, M., Johansen, A., 2016. Dust Evolution and the Forma- mal settings are required to establish more detailed regimes tion of Planetesimals. Space Sci. Rev. 205, 41–75. Brett, R., Sato, M., 1984. Intrinsic oxygen fugacity measurements on seven for melt segregation. chondrites, a pallasite, and a tektite and the redox state of meteorite In the mid-term, future spacecraft missions (A’Hearn, parent bodies. Geochim. Cosmochim. Acta 48, 111–120. 2017) may be able to deliver further observational con- Burbine, T. H., Demeo, F. E., Rivkin, A. S., Reddy, V., 2017. Evidence for straints on asteroid-belt objects. In conjunction with self- differentiation among asteroid families. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: Early Differentiation and Consequences for consistent evolutionary models of metal–silicate and solid– Planets. Cambridge University Press, Cambridge, pp. 298–320. melt segregation, these can help to further decipher the in- Carporzen, L., Weiss, B. P., Elkins-Tanton, L. T., Shuster, D. L., Ebel, D., terior evolution of planetary bodies in the early solar system Gattacceca, J., 2011. Magnetic evidence for a partially differentiated car- bonaceous chondrite parent body. Proc. Natl. Acad. Sci. U.S.A. 108, and sharpen our understanding of terrestrial planet formation 6386–6389. in the solar system and elsewhere. Acknowledgements. The authors acknowledge stimulating 2 matplotlib.org (Hunter, 2007) discussions with Ian Sanders and Guy Consolmagno, and seaborn.pydata.org 13 Carter, P. J., Leinhardt, Z. M., Elliott, T., Stewart, S. T., Walter, M. J., Feb. Kataoka, A., Tanaka, H., Okuzumi, S., Wada, K., 2013. Fluffy dust forms icy 2018. Collisional stripping of planetary crusts. Earth and Planetary Sci- planetesimals by static compression. Astron. Astrophys. 557, L4. ence Letters 484, 276–286. Katz, R. F., Knepley, M., Smith, B., Spiegelman, M., Coon, E., 2007. Nu- Castillo-Rogez, J., Young, E. D., 2017. Origin and evolution of volatile-rich merical simulation of geodynamic processes with the portable extensible asteroids. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: toolkit for scientiﬁc computation. Phys. Earth Planet. Inter. 163, 52–68. Early Differentiation and Consequences for Planets. Cambridge Univer- Keller, T., Katz, R. F., 2016. The Role of Volatiles in Reactive Melt Transport sity Press, Cambridge, pp. 92–114. in the Asthenosphere. J. Petrol. 57, 1073–1108. Cerantola, V., Walte, N. P., Rubie, D. C., 2015. Deformation of a crystalline Keller, T., May, D. A., Kaus, B. J. P., 2013. Numerical modelling of magma olivine aggregate containing two immiscible liquids: Implications for early dynamics coupled to tectonic deformation of lithosphere and crust. Geo- core-mantle differentiation. Earth Planet. Sci. Lett. 417, 67–77. phys. J. Int. 195, 1406–1442. Clenet, H., Jutzi, M., Barrat, J.-A., Asphaug, E. I., Benz, W., Gillet, P., 2014. Keller, T., Suckale, J., 2018. A continuum model of multi-phase reactive A deep crust-mantle boundary in the asteroid 4 Vesta. Nature 511, 303– transports in igneous systems. arXiv:1809.00079. 306. Kita, N. T., Yin, Q.-Z., MacPherson, G. J., Ushikubo, T., Jacobsen, B., Na- 26 26 Consolmagno, G. J., Golabek, G. J., Turrini, D., Jutzi, M., Sirono, S., gashima, K., Kurahashi, E., Krot, A. N., Jacobsen, S. B., 2013. Al- Mg Svetsov, V., Tsiganis, K., 2015. Is Vesta an intact and pristine proto- isotope systematics of the ﬁrst solids in the early solar system. Meteorit. planet? Icarus 254, 190–201. Planet. Sci. 48, 1383–1400. Costa, A., Caricchi, L., Bagdassarov, N., 2009. A model for the rheology Lichtenberg, T., Golabek, G. J., Dullemond, C. P., Schönbächler, M., Gerya, of particle-bearing suspensions and partially molten rocks. Geochem. T. V., Meyer, M. R., 2018. Impact splash chondrule formation during plan- Geophys. Geosys. 10, Q03010. etesimal recycling. Icarus 302, 27–43. Cournede, C., Gattacceca, J., Gounelle, M., Rochette, P., Weiss, B. P., Lichtenberg, T., Golabek, G. J., Gerya, T. V., Meyer, M. R., 2016a. The Zanda, B., 2015. An early solar system magnetic ﬁeld recorded in CM effects of short-lived radionuclides and porosity on the early thermo- chondrites. Earth Planet. Sci. Lett. 410, 62–74. mechanical evolution of planetesimals. Icarus 274, 350–365. Dra ¸ zk ˙ owska, J., Dullemond, C. P., 2018. Planetesimal formation during pro- Lichtenberg, T., Parker, R. J., Meyer, M. R., 2016b. Isotopic enrichment of toplanetary disk buildup. Astron. Astrophys. 614, A62. forming planetary systems from supernova pollution. Mon. Not. R. As- Dunn, T. L., Cressey, G., McSween, Jr., H. Y. ., McCoy, T. J., 2010. Analysis tron. Soc. 462, 3979–3992. of ordinary chondrites using powder X-ray diffraction: 1. Modal mineral Lodders, K., 2003. Solar System Abundances and Condensation Tempera- abundances. Meteorit. Planet. Sci. 45, 123–134. tures of the Elements. Astrophys. J. 591, 1220–1247. Ebel, D. S., Grossman, L., 2000. Condensation in dust-enriched systems. Mandler, B. E., Elkins-Tanton, L. T., 2013. The origin of eucrites, diogenites, Geochim. Cosmochim. Acta 64, 339–366. and olivine diogenites: Magma ocean crystallization and shallow magma Elkins-Tanton, L. T., Weiss, B. P., Zuber, M. T., 2011. Chondrites as samples chamber processes on Vesta. Meteorit. Planet. Sci. 48, 2333–2349. of differentiated planetesimals. Earth Planet. Sci. Lett. 305, 1–10. McKenzie, D., 1984. The generation and compaction of partially molten Friedrich, J. M., Weisberg, M. K., Ebel, D. S., Biltz, A. E., Corbett, B. M., rock. J. Petrol. 5, 713–765. Iotzov, I. V., Khan, W. S., Wolman, M. D., 2015. Chondrule size and McSween, H. Y., Bennett, M. E., Jarosewich, E., 1991. The mineralogy of or- related physical properties: A compilation and evaluation of current data dinary chondrites and implications for asteroid spectrophotometry. Icarus across all meteorite groups. Chemie der Erde / Geochemistry 75, 419– 90, 107–116. 443. Mei, S., Bai, W., Hiraga, T., Kohlstedt, D. L., 2002. Inﬂuence of melt on the Fu, R. R., Elkins-Tanton, L. T., 2014. The fate of magmas in planetesimals creep behavior of olivine-basalt aggregates under hydrous conditions. and the retention of primitive chondritic crusts. Earth Planet. Sci. Lett. Earth Planet. Sci. Lett. 201, 491–507. 390, 128–137. Mizzon, H., 2015. The magmatic crust of vesta. Ph.D. thesis, Université de Fu, R. R., Young, E. D., Greenwood, R. C., Elkins-Tanton, L. T., 2017. Sil- Toulouse, Université Toulouse III-Paul Sabatier. icate melting and volatile loss during differentiation in planetesimals. In: Monteux, J., Golabek, G. J., Rubie, D. C., Tobie, G., Young, E. D., 2018. Wa- Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Planetesimals: Early Differentia- ter and the interior structure of terrestrial planets and icy bodies. Space tion and Consequences for Planets. Cambridge University Press, Cam- Sci. Rev. 214, 39. bridge, pp. 115–135. Moskovitz, N., Gaidos, E., 2011. Differentiation of planetesimals and the Ghanbarzadeh, S., Hesse, M. A., Prodanovic, ´ M., 2017. Percolative core thermal consequences of melt migration. Meteorit. Planet. Sci. 46, 903– formation in planetesimals enabled by hysteresis in metal connectivity. 918. Proc. Natl. Acad. Sci. 114, 13406–13411. Mulyukova, E., Bercovici, D., 2018. Collapse of passive margins by litho- Ghosh, A., McSween, H. Y., 1998. A Thermal Model for the Differentiation spheric damage and plunging grain size. Earth Planet. Sci. Lett. 484, of Asteroid 4 Vesta, Based on Radiogenic Heating. Icarus 134, 187–206. 341–352. Golabek, G. J., Bourdon, B., Gerya, T. V., 2014. Numerical models of Neumann, W., Breuer, D., Spohn, T., 2013. The thermo-chemical evolution the thermomechanical evolution of planetesimals: Application to the of Asteroid 21 Lutetia. Icarus 224, 126–143. acapulcoite-lodranite parent body. Meteorit. Planet. Sci. 49, 1083–1099. Neumann, W., Breuer, D., Spohn, T., 2014. Differentiation of Vesta: Implica- Hevey, P. J., Sanders, I. S., 2006. A model for planetesimal meltdown by tions for a shallow magma ocean. Earth Planet. Sci. Lett. 395, 267–280. Al and its implications for meteorite parent bodies. Meteorit. Planet. Neumann, W. O., Kruijer, T. S., Breuer, D., Kleine, T., 2018. Multi-Stage Sci. 41, 95–106. Core Formation in Planetesimals Revealed by Numerical Models and Hf- Hiraga, T., Tachibana, C., Ohashi, N., Sano, S., 2010. Grain growth system- W Chronometry of Iron Meteorites. J. Geophys. Res. Planets 123, 421– atics for forsterite enstatite aggregates: Effect of lithology on grain size 444. in the upper mantle. Earth Planet. Sci. Lett. 291, 10–20. Pätzold, M., Andert, T. P., Asmar, S. W., Anderson, J. D., Barriot, J.-P., Bird, Hirth, G., Kohlstedt, D., 2003. Rheology of the upper mantle and the mantle M. K., Häusler, B., Hahn, M., Tellmann, S., Sierks, H., Lamy, P., Weiss, wedge: A view from the experimentalists. Inside the subduction Factory, B. P., 2011. Asteroid 21 Lutetia: Low Mass, High Density. Science 334, 83–105. 491. Hunt, A. C., Cook, D. L., Lichtenberg, T., Reger, P. M., Ek, M., Golabek, Raymond, C. A., Russell, C. T., McSween, H. Y., 2017. Dawn at vesta: G. J., Schönbächler, M., 2018. Late metal-silicate separation on the IAB Paradigms and paradoxes. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), parent asteroid: Constraints from combined W and Pt isotopes and ther- Planetesimals: Early Differentiation and Consequences for Planets. mal modelling. Earth Planet. Sci. Lett. 482, 490–500. Cambridge University Press, Cambridge, pp. 321–339. Hunter, J. D., 2007. Matplotlib: A 2d graphics environment. Computing In Rozel, A., Ricard, Y., Bercovici, D., 2011. A thermodynamically self- Science & Engineering 9, 90–95. consistent damage equation for grain size evolution during dynamic re- Hutchison, R., 2004. Meteorites, Cambridge University Press. crystallization. Geophys. J. Int. 184, 719–728. Johansen, A., Mac Low, M.-M., Lacerda, P., Bizzarro, M., 2015. Growth Rubie, D. C., Jacobson, S. A., Morbidelli, A., O’Brien, D. P., Young, E. D., de of asteroids, planetary embryos, and Kuiper belt objects by chondrule Vries, J., Nimmo, F., Palme, H., Frost, D. J., 2015. Accretion and differ- accretion. Sci. Adv. 1, 1500109. entiation of the terrestrial planets with implications for the compositions 14 of early-formed Solar System bodies and accretion of water. Icarus 248, 89–108. Rudge, J. F., 2018. Textural equilibrium melt geometries around tetrakaidec- ahedral grains. Proc. Royal Soc. Lond. A 474, 20170639. Rudge, J. F., Bercovici, D., Spiegelman, M., 2011. Disequilibrium melting of a two phase multicomponent mantle. Geophys. J. Int. 184, 699–718. Sierks, H., Lamy, P., Barbieri, C., Koschny, D., Rickman, H., Rodrigo, R., A’Hearn, M. F., et al., 2011. Images of Asteroid 21 Lutetia: A Remnant Planetesimal from the Early Solar System. Science 334, 487. Solomatov, V. S., 2015. Magma oceans and primordial mantle differentia- tion. Treatise on Geophysics 2nd ed., pp. 81–104. Stevenson, D. J., 1990. Fluid dynamics of core formation. In: Newsom, H. E., Jones, J. H. (Eds.), Origin of the Earth. pp. 231–249. Tang, H., Dauphas, N., 2015. Low Fe Abundance in Semarkona and Sa- hara 99555. Astrophys. J. 802, 22–31. Tkalcec, B. J., Golabek, G. J., Brenker, F. E., 2013. Solid-state plastic de- formation in the dynamic interior of a differentiated asteroid. Nat. Geo. 6, 93–97. Šrámek, O., Milelli, L., Ricard, Y., Labrosse, S., 2012. Thermal evolution and differentiation of planetesimals and planetary embryos. Icarus 217, 339–354. van Boekel, R., Min, M., Leinert, C., Waters, L. B. F. M., Richichi, A., Ches- neau, O., Dominik, C., Jaffe, W., Dutrey, A., Graser, U., Henning, T., de Jong, J., Köhler, R., de Koter, A., Lopez, B., Malbet, F., Morel, S., Paresce, F., Perrin, G., Preibisch, T., Przygodda, F., Schöller, M., Wit- tkowski, M., 2004. The building blocks of planets within the ‘terrestrial’ region of protoplanetary disks. Nature 432, 479–482. Visser, R. G., Ormel, C. W., 2016. On the growth of pebble-accreting plan- etesimals. Astron. Astrophys. 586, A66. Wahlberg Jansson, K., Johansen, A., 2017. Radially resolved simulations of collapsing pebble clouds in protoplanetary discs. Mon. Not. R. Astron. Soc. 469, 149–157. Wasson, J. T., 1990. Ungrouped iron meteorites in Antarctica - Origin of anomalously high abundance. Science 249, 900–902. Weidenschilling, S. J., 1977. Aerodynamics of solid bodies in the solar neb- ula. Mon. Not. R. Astron. Soc. 180, 57–70. Weiss, B. P., Elkins-Tanton, L. T., 2013. Differentiated Planetesimals and the Parent Bodies of Chondrites. Annu. Rev. Earth Planet. Sci. 41, 529–560. Weiss, B. P., Elkins-Tanton, L. T., Antonietta Barucci, M., Sierks, H., Snod- grass, C., Vincent, J.-B., Marchi, S., Weissman, P. R., Pätzold, M., Richter, I., Fulchignoni, M., Binzel, R. P., Schulz, R., 2012. Possi- ble evidence for partial differentiation of asteroid Lutetia from Rosetta. Planet. Space Sci. 66, 137–146. Wilson, L., Keil, K., 2012. Volcanic activity on differentiated asteroids: A review and analysis. Chemie der Erde / Geochemistry 72, 289–321. Wilson, L., Keil, K., 2017. Arguments for the non-existence of magma oceans in asteroids. In: Elkins-Tanton, L. T., Weiss, B. P. (Eds.), Plan- etesimals: Early Differentiation and Consequences for Planets. Cam- bridge Planetary Science. Cambridge University Press, pp. 159–179. Yoshino, T., Walter, M. J., Katsura, T., 2003. Core formation in planetesimals triggered by permeable ﬂow. Nature 422, 154–157.

Journal

Astrophysics – arXiv (Cornell University)

Published: Feb 6, 2018

References

The future of NASA's missions

A'Hearn
Constraints on the Fe–S melt connectivity in mantle silicates from electrical impedance measurements

Bagdassarov
Efficient management of parallelism in object oriented numerical software libraries

Balay
Grain-damage hysteresis and plate tectonic states

Bercovici
Dust evolution and the formation of planetesimals

Birnstiel
Intrinsic oxygen fugacity measurements on seven chondrites, a pallasite, and a tektite and the redox state of meteorite parent bodies

Brett
Evidence for differentiation among asteroid families

Burbine
Magnetic evidence for a partially differentiated carbonaceous chondrite parent body

Carporzen
Collisional stripping of planetary crusts

Carter
Origin and evolution of volatile-rich asteroids

Castillo-Rogez
Deformation of a crystalline olivine aggregate containing two immiscible liquids: implications for early core–mantle differentiation

Cerantola
A deep crust–mantle boundary in the Asteroid 4 Vesta

Clenet
Is Vesta an intact and pristine protoplanet?

Consolmagno
A model for the rheology of particle-bearing suspensions and partially molten rocks

Costa
An early solar system magnetic field recorded in CM chondrites

Cournede
Planetesimal formation during protoplanetary disk buildup

Drążkowska
Analysis of ordinary chondrites using powder X-ray diffraction: 1. Modal mineral abundances

Dunn
Condensation in dust-enriched systems

Ebel
Chondrites as samples of differentiated planetesimals

Elkins-Tanton
Chondrule size and related physical properties: a compilation and evaluation of current data across all meteorite groups

Friedrich
The fate of magmas in planetesimals and the retention of primitive chondritic crusts

Fu
Silicate melting and volatile loss during differentiation in planetesimals

Fu
Percolative core formation in planetesimals enabled by hysteresis in metal connectivity

Ghanbarzadeh
A thermal model for the differentiation of Asteroid 4 Vesta, based on radiogenic heating

Ghosh
Numerical models of the thermomechanical evolution of planetesimals: application to the acapulcoite-lodranite parent body

Golabek
A model for planetesimal meltdown by 26Al and its implications for meteorite parent bodies

Hevey
Grain growth systematics for forsterite ± enstatite aggregates: effect of lithology on grain size in the upper mantle

Hiraga
Rheology of the upper mantle and the mantle wedge: a view from the experimentalists

Hirth
Late metal–silicate separation on the IAB parent asteroid: constraints from combined W and Pt isotopes and thermal modelling

Hunt
Matplotlib: a 2d graphics environment

Hunter

Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Magma ascent in planetesimals: control by grain size

Magma ascent in planetesimals: control by grain size

Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Magma ascent in planetesimals: control by grain size

Magma ascent in planetesimals: control by grain size

References

Abstract

Journal

Recommended Articles

References

Our policy towards the use of cookies