How primordial is the structure of comet 67P/C-G? Combined collisional and dynamical models suggest a late formation
AAstronomy & Astrophysics manuscript no. Jutzi_etal_final_corrected c (cid:13)
ESO 2018May 25, 2018
How primordial is the structure of comet 67P/C-G?
Combined collisional and dynamical models suggest a late formation
M. Jutzi , W. Benz , A. Toliou , , A. Morbidelli , R. Brasser Physics Institute, University of Bern, NCCR PlanetS, Sidlerstrasse 5, 3012 Bern, Switzerlande-mail: [email protected]; [email protected] Department of Physics, Aristotle University of Thessaloniki, GR-54124 Thessaloniki, Greecee-mail: [email protected] Laboratoire Lagrange, Université Côte d’Azur, CNRS, Observatoire de la Côte d’Azur, Nice, France Earth Life Science Institute, Tokyo Institute of Technology, Meguro-ku, Tokyo, 152-8550, JapanReceived – ; accepted –
ABSTRACT
Context.
There is an active debate about whether the properties of comets as observed today are primordial or, alternatively, if theyare a result of collisional evolution or other processes.
Aims.
We investigate the e ff ects of collisions on a comet with a structure like 67P / Churyumov-Gerasimenko (hereafter 67P / C-G). Wedevelop scaling laws for the critical specific impact energies Q reshape required for a significant shape alteration. These are then used insimulations of the combined dynamical and collisional evolution of comets in order to study the survival probability of a primordiallyformed object with a shape like 67P / C-G. Although the focus of this work is on a structure of this kind, the analysis is also performedfor more generic bi-lobe shapes, for which we define the critical specific energy Q bil . The simulation outcomes are also analyzed interms of impact heating and the evolution of the porosity. Methods.
The e ff ects of impacts on comet 67P / C-G are studied using a state-of-the-art smooth particle hydrodynamics (hereafterSPH) shock physics code. In the 3D simulations, a publicly available shape model of 67P / C-G is applied and a range of impactconditions and material properties are investigated. The resulting critical specific impact energy Q reshape (as well as Q bil for generic bi-lobe shapes) defines a minimal projectile size which is used to compute the number of shape-changing collisions in a set of dynamicalsimulations. These simulations follow the dispersion of the trans-Neptunian disk during the giant planet instability, the formationof a scattered disk, and produce 87 objects that penetrate into the inner solar system with orbits consistent with the observed JFCpopulation. The collisional evolution before the giant planet instability is not considered here. Hence, our study is conservative in itsestimation of the number of collisions. Results.
We find that in any scenario considered here, comet 67P / C-G would have experienced a significant number of shape-changingcollisions, if it formed primordially. This is also the case for generic bi-lobe shapes. Our study also shows that impact heating is verylocalized and that collisionally processed bodies can still have a high porosity.
Conclusions.
Our study indicates that the observed bi-lobe structure of comet 67P / C-G may not be primordial, but might haveoriginated in a rather recent event, possibly within the last 1 Gy. This may be the case for any kilometer-sized two-componentcometary nuclei.
Key words.
Comets: general – Comets: individual: 67P / C-G – Kuiper belt: general – Planets and satellites: formation
1. Introduction
Comets or their precursors formed in the outer planet region dur-ing the initial stages of solar system formation. They may havebeen assembled by hierarchical accretion (e.g. Weidenschilling1997; Windmark et al. 2012b,a; Kataoka et al. 2013) or, alterna-tively, were born big in gravitational instabilities (e.g. Youdin& Goodman 2005; Johansen et al. 2007; Cuzzi et al. 2010;Morbidelli et al. 2009), thereby bypassing the primary accre-tion phase entirely. Whether the cometary nuclei structures asobserved today are pristine and preserve a record of their orig-inal accumulation, or are a result of later collisional or otherprocesses is much debated (e.g. Weissman et al. 2004; Mummaet al. 1993; Sierks et al. 2015; Rickman et al. 2015; Morbidelli &Rickman 2015; Jutzi & Asphaug 2015; Davidsson et al. 2016).The shape, density, composition, and other properties of comet67P / Churyumov-Gerasimenko (67P / C-G) have been determined in detail by the European Space Agency’s Rosetta rendezvousmission (e.g. Sierks et al. 2015; Hässig et al. 2015; Capaccioniet al. 2015). Based on this data, it has been suggested that itsstructure is pristine and was formed in the early stages of the so-lar system (Massironi et al. 2015), possibly by low velocity ac-cretionary collisions (Jutzi & Asphaug 2015). What is less clearis whether or not a structure like comet 67P / C-G would havebeen able to survive until today.The collisional evolution of an object of the size of comet67P / C-G was studied by Morbidelli & Rickman (2015) usingdynamical models of the evolution of the early solar system. Inthe "standard model", as defined by the so-called Nice model(Tsiganis et al. 2005), cometary nuclei, or their precursors, orig-inated from an initial trans-planetary disk of icy planetesimalswith a lifetime of a few hundred Myr. In this concept, the trans-planetary disk formed in the infant stages of the solar system
Article number, page 1 of 13 a r X i v : . [ a s t r o - ph . E P ] N ov & A proofs: manuscript no. Jutzi_etal_final_corrected beyond the original orbits of all giant planets, which were ini-tially closer to the Sun. This disk may have given rise to both theScattered Disk and the Oort cloud (Brasser & Morbidelli 2013)and thus, it may once have been the repository for all the cometsobserved today. According to the standard assumption, the dis-persal of the disk coincided with the beginning of the so-calledLate Heavy Bombardment (Gomes et al. 2005; Morbidelli et al.2012), and had a lifetime of about 450 Myr before it was dynam-ically dispersed.As shown in Morbidelli & Rickman (2015), it is clear thatin this standard model, an object of the size of comet 67P / C-Gwould have experienced a high number of catastrophic collisionsand thus could not have survived. However, it was also shownthat in the (hypothetical) case that the dispersal of the disk oc-curred early, right after gas removal, the collisional evolution ofkm-size bodies ending in the Scattered Disk would have beenless severe, and a fraction of these objects may have escaped allcatastrophic collisions. We also note that in alternative models(e.g. Davidsson et al. 2016), the total number of comets is con-sidered to be lower than previously thought. Therefore, the fateof cometary-sized objects appears to depend upon the details ofthe dynamical scenario considered.However, whether or not an object like comet 67P / C-Gwould have been able to survive until today does not only dependupon its dynamical evolution but even more so on the "strength"of the body. In other words, it is crucial to know the critical spe-cific impact energy at which the shape and structure of such anobject are altered significantly. Previous studies of the collisionalevolution of comet 67P / C-G (Morbidelli & Rickman 2015) usedscaling laws for catastrophic disruption energies that are basedon idealized spherical, solid icy bodies (Benz & Asphaug 1999),which may not represent well the properties of a highly porouscometary nuclei. It is well known that the impacts in highlyporous material, given its dissipative properties, can lead to verydi ff erent results compared to impacts involving solid materials(e.g. Housen & Holsapple 2003; Jutzi et al. 2008). Furthermore,complex shapes such as the one of 67P / C-G may already be sub-stantially a ff ected by relatively low energy, sub-catastrophic im-pacts.It has been suggested recently that rotational fission and re-configuration may be a dominant structural evolution processfor short-period comet nuclei having a two-component struc-ture with a volume ratio larger than ∼ / C-G) as observed today would then bethe result of the last merger in this cycle. In this context, it isimportant to also study the survival probability of more generaltwo-component structures, as such structures are required for thefission-merging cycle to begin.In this paper, we consider both the dynamical evolution andthe response to impacts of objects with a 67P / C-G-like shape aswell as generic bi-lobe structures. This combined approach al-lows us to compute the expected number of shape-changing col-lisions for such objects, as well as the related survival probabilityand possible formation age.In the first part of the paper, we describe our modeling ap-proach to study the e ff ects of impacts on comet 67P / C-G andgeneric bi-lobe shapes. In section 2, we determine the specificenergies Q reshape required to change a 67P / C-G-like shape sig-nificantly, as well as the corresponding Q bil for reshaping genericbi-lobe objects. The catastrophic disruption threshold Q ∗ D forbodies of 67P / C-G size, with the same properties, is computed as well here. Using the result of this modeling, we develop scalinglaws for Q reshape , Q bil and Q ∗ D . Finally, the simulation outcomesare analyzed in terms of impact heating and the evolution of theporosity.In the second part of the paper, we first describe the details ofthe dynamical simulations used in this study and discuss the dif-ferences and the improvements with respect to the previous workby Morbidelli & Rickman (2015) (section 3). We then computethe average number of shape-changing collisions for a body witha 67P / C-G-like shape as well as for generic bi-lobe shapes, usingthe corresponding scaling laws ( Q reshape and Q bil ). In section 4,the uncertainties of our model as well as alternative models arediscussed, followed by conclusions in section 5.A scenario of the late formation of 67P / C-G-like (two-lobe)shapes by a new type of sub-catastrophic impacts is presentedin a companion paper (Jutzi&Benz, 2016 submitted; hereafterPaper II).
2. The effects of impacts on bi-lobe structures
Here, in a suite of 3D smooth particle hydrodynamics (SPH)code calculations, we compute the specific impact energy Q reshape required to significantly change the shape of comet67P / C-G as well as of generic bi-lobe structures. The catas-trophic disruption threshold Q ∗ D for spherical objects of the samemass is computed as well. We consider a range of material(strength) properties and various impact conditions. The simu-lation outcomes are also analysed in terms of impact heating andthe evolution of the porosity. Cometary nuclei come apart easily due to tides (Asphaug &Benz 1994) and other gentle stresses (Boehnhardt 2004). Lab-oratory materials analysis (Skorov & Blum 2012), observationsof comet disruptions by tides (Asphaug & Benz 1994) or frag-mentation through dynamic sublimation pressure (Stecklo ff et al.2015), suggest a bulk strength of <
10 - 100 Pa for these weaklyconsolidated bodies. On the other hand, a high compressivestrength of surface layers on comet 67P / C-G (Biele et al. 2015)was found at 0.1-1 m scales. For our analysis of the overall sta-bility, this kind of small scale ( < ∼
10 m) strength is not relevant,as we are interested in the bulk properties. In our modeling, wethus consider bulk tensile strengths of up to 1 kPa. The corre-sponding values of cohesion and compressive strength are ∼ anorder of magnitude higher (see section 2.2).The low bulk densities of comets indicate substantial poros-ity; in the case of comet 67P / C-G it is about 75% (e.g. Pätzoldet al. 2016). In our modeling approach (section 2.2) it is im-plicitly assumed that porosity is at small scales and the body ishomogenous. In the case of comet 67P / C-G, recent gravity fieldanalysis (Pätzold et al. 2016) indicate that the interior of the nu-cleus is homogeneous (down to scales of ∼ The modeling approach used here has already been successfullyapplied in a previous study to the regime of cometesimal colli-sions (Jutzi & Asphaug 2015). We use a parallel smooth parti-cle hydrodynamics (SPH) impact code (Benz & Asphaug 1995;Ny ff eler 2004; Jutzi et al. 2008; Jutzi 2015) which includes self- Article number, page 2 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G?
Table 1.
Material parameters used in the simulations. Crush curve parameters P e and P s (Jutzi et al. 2008), density of matrix material ρ s , initialbulk density ρ , density of the compacted material ρ compact , initial distention α , bulk modulus A , friction coe ffi cient µ , cohesion Y , average tensilestrength Y T . type P e (Pa) P s (Pa) ρ s (kg / m ) ρ (kg / m ) ρ compact (kg / m ) α A (Pa) µ Y (Pa) Y T (Pa)low strength 10
910 440 1980 4.5 2.67 × medium str. (nominal) 10
910 440 1980 4.5 2.67 × high strength 10
910 440 1980 4.5 2.67 × P o r o s i t y Pressure (Pa)Pressure-porosity relation (high strength)Pressure-porosity relation (medium strength)Pressure-porosity relation (low strength)Dust agglomerates (3D); Guetller et al., 2009Ice pebbles (quasi 3D); Lorek et al., 2016
Fig. 1.
Pressure-porosity relations (crush curve) used in the simulationsfor the three di ff erent sets of parameters (low, medium, high strength)as defined in Table 1. Also shown are the results from laboratory exper-iments dust agglomerates (Güttler et al. 2009) and ice pebbles (Loreket al. 2016). gravity as well as material strength models. To avoid numericalrotational instabilities, the scheme suggested by Speith (2006) isused.In our modeling, we include an initial cohesion Y > Y T varying from ∼
10 to ∼ Y T =
100 Pa asthe nominal case. To model fractured, granular material, a pres-sure dependent shear strength (friction) is included by using astandard Drucker-Prager yield criterion (Jutzi 2015). As shownin Jutzi (2015) and Jutzi et al. (2015), granular flow problems arewell reproduced using this method.Porosity is modeled using a P-alpha model (Jutzi et al. 2008)with a simple quadratic crush curve defined by the parameters P e , P s , ρ , ρ s and α . We further introduce the density of thecompacted material ρ compact = / m , which is used definethe initial distention α = ρ compact /ρ = . − /α ∼ ρ s in this model is a parameter determin-ing the form of the crush curve and does not correspond to thedensity of the fully compacted material). As an estimate of thecompressive strength σ c = P s / c e for a porous aggregate body can be very low, of the order of c e ∼ / s. To take this into account, we apply a reduced bulkmodulus (leading term in the Tillotson EOS; see Table 1). Theapproach has the additional major advantage that the time-stepsbecome large enough to carry out the simulations over many dy-namical timescales. Di ff erent values of c e =
10 - 100 m / s areinvestigated.The relevant material parameters used in the simulations areindicated in Table 1. To setup the target, we apply a publicly available shape model ofcomet 67P / C-G , which defines the surface of the body. Threedi ff erent sets of material parameters as indicated in Table 1 areused, corresponding to di ff erent target strength.To determine Q reshape for 67P / C-G-like shapes, we investi-gate a range of impact energies using a range of impactor sizesof R p = ρ ∼
440 kg / m . We only consider impactsinto the smaller of the the lobes of comet 67P / C-G. Two di ff er-ent impact geometries are investigated (Figure 2).To determine the critical shape-changing impact energy Q bil in the case of more general bi-lobe structures, we set up a targetconsisting of two overlapping ellipsoids (Figure 6). Each ellip-soid has an axis ratio of 0.6. The volume ratio between the twocomponents is ∼ . M t = × kg. For these targets, we only use the nominal set of strengthparameters (Table 1) and an impactor size of R p =
100 m.The simulations are carried out using a moderately high res-olution of ∼ × SPH particles.
In addition to Q reshape and Q bil , we also investigate the criticalspecific energy for catastrophic disruption Q ∗ D of spherical bod-ies of the same mass and material properties as in the model ofcomet 67P / C-G. To do this we consider 3 di ff erent size ratios ofprojectile and target (1:2; 1:4; 1:8) and varying impact velocities.The impact angle is fixed to 45 ◦ . The results of our modeling of impacts on 67P / C-G are displayedin Figures 2-5. We find that this particular structure, with two http: // sci.esa.int / rosetta / / Article number, page 3 of 13 & A proofs: manuscript no. Jutzi_etal_final_corrected v (m/s) 10 30 50 80 100Q (J/Kg) 0.01 0.08 0.23 0.59 0.92 R proj = 100 m Q crit ~ 0.1 - 0.3 J/Kg Tensile Strength: 10 Pa Geometry 1 c e = 10 m/sGeometry 1 c e = 100 m/sGeometry 2 c e = 100 m/sGeometry 1 c e = 100 m/srotating Fig. 2.
Shape-changing collisions on comet 67P / C-G. We use SPH to simulate impacts of a R p =
100 m projectile on the smaller of the two lobes ofcomet 67P / C-G. The minimal specific energy required to cause a significant change of the comet’s shape by such impacts, Q reshape , is estimated fordi ff erent impact geometries and rotation axis. The material strength is the same in all cases shown here ( Y T =
10 Pa). The e ff ect of the material’ssound speed is investigated as well (top row; in this case, a bulk modulus of A = × Pa instead of the nominal A = × Pa is used).Plotted is a surface of constant density which represents the surface of the comet; shown in red are regions on the surface with materials whoseprescribed tensile strength was exceeded. As a rough average, the minimal specific energy required to cause a significant shape change is estimatedas Q reshape ∼ ± / kg, as indicated by the horizontal yellow line. lobes connected by a neck, is significantly altered even by rela-tively low energy impacts. For a fixed set of material parameters(i.e. constant strength), the di ff erent impact geometries and ro-tation states considered here lead to slightly di ff erent outcomes(Figure 2), but there are no major, order of magnitude, di ff er-ences between the various runs.As it can be observed in Figure 3, higher material strengthlead to higher specific impact energy required to reach a certaindegree of change in the overall shape.There is no unique way to define the critical shape-changingspecific impact energy from these results, but rough estimatesare possible. Based on visual inspection, we define Q reshape forthe di ff erent strength as: Q reshape ∼ ± / kg for Y T = Q reshape ∼ ± / kg for Y T =
100 Pa; Q reshape ∼ ± / kg for Y T = R p =
100 m projectile. For the simulations with the larger projectileswe obtain Q reshape ∼ ± / kg ( R p =
200 m; Figure 4) and Q reshape ∼ ± / kg ( R p =
300 m; Figure 5), using thenominal strength of Y T =
100 Pa. These values are plotted inFigure 7 and compared to the catastrophic disruption threshold,as discussed below. We note that impacts into the larger lobe maylead to slightly di ff erent values for Q reshape , but we do not expectorder of magnitude di ff erences.The results of our modeling of impacts on generic bi-lobeshapes (using nominal strength properties) are displayed in Fig-ure 6. The estimated minimal specific impact energies for re- shaping are Q bil ∼ ± / kg, which is slightly higher than inthe case of the 67P / C-G-like shape with the same strength ( Q bil [ Y T =
100 Pa] ∼ corresponds to Q reshape for the Y T = The results of our modeling of catastrophic disruptions of spher-ical bodies with the same mass and material properties as in themodel of comet 67P / C-G are shown in Figure 7. We define thespecific impact energy as Q = . µ r V / ( M t + M p ) (1)where µ r = M p M t / ( M t + M p ) is the reduced mass, M p is themass of the projectile and V the impact velocity. For the oblique(45 ◦ ) impacts considered here, we also take into account thatonly a part of the mass of the colliding bodies is interacting(Leinhardt & Stewart 2012), and compute the Q ∗ D values of theequivalent head-on collisions.As expected, the energy threshold for catastrophic disruption Q ∗ D >> Q reshape , by ∼ two orders of magnitude.As found in previous studies (e.g. Jutzi 2015), in the disrup-tion regime, the results for Q ∗ D are almost independent of thematerial (tensile) strength.Our values of Q ∗ D for di ff erent impact velocities (Figure 7)agree well with scaling law predictions (Housen & Holsapple Article number, page 4 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G? v (m/s) 10 50 100 150 200 300 400Q (J/Kg) 0.01 0.23 0.92 2.07 3.69 8.29 14.74 R proj = 100 m, Geometry 1, cs =100 m/s P a100 P a1 k P a
10 Pa: Q crit ~ 0.10 - 0.30 J/Kg (v ~ 50m/s) 100 Pa: Q crit ~ 0.50 - 1.50 J/Kg (v ~ 100 m/s)1 kPa: Q crit ~ 1.00 - 3.00 J/Kg (v ~ 150 m/s)
Fig. 3.
Same as Fig. 2 but for di ff erent material strength Y T of the target. c e =
100 m / s in all cases. The critical specific energies are: Q reshape ∼ ± / kg for Y T =
10 Pa (corresponds to average in Figure 2); Q reshape ∼ ± / kg for Y T =
100 Pa; Q reshape ∼ ± / kg for Y T = v (m/s) 5 10 20 30Q (J/Kg) 0.02 0.07 0.29 0.66 R = 200: Q crit ~ 0.2 - 0.4 J/Kg (v ~ 20 m/s)R = 300: Q crit ~ 0.075 - 0.0225 J/Kg (v ~ 7.5 m/s)
100 Pa, R = 200
Geometry 1, cs =100 m/s v (m/s) 5 10 20 30Q (J/Kg) 0.06 0.25 1.00 2.24
100 Pa, R = 300
Fig. 4.
Same as Fig. 3 but for R p =
200 m ( Y T =
100 Pa). Q reshape ∼ ± / kg. v (m/s) 5 10 20 30Q (J/Kg) 0.02 0.07 0.29 0.66 R = 200: Q crit ~ 0.2 - 0.4 J/Kg (v ~ 20 m/s)R = 300: Q crit ~ 0.075 - 0.0225 J/Kg (v ~ 7.5 m/s)
100 Pa, R = 200
Geometry 1, cs =100 m/s v (m/s) 5 10 20 30Q (J/Kg) 0.06 0.25 1.00 2.24
100 Pa, R = 300
Fig. 5.
Same as Fig. 3 but for R p =
300 m ( Y T =
100 Pa). Q reshape ∼ ± / kg. v (m/s) 50 100 200Q (J/Kg) 0.23 0.92 3.69 R proj = 100 m Bi-ellipsoid
10 Pa: Q crit ~ 0.1 - 0.3 J/Kg (v ~ 50 m/s)100 Pa: Q crit ~ 1.0 - 3.0 J/Kg (v ~ 150 m/s)1 kPa: Q crit ~ 4.0 - 12 J/Kg (v ~ 300 m/s)
100 Pa
Fig. 6.
Results of impacts on generic bi-lobe shapes with nominalstrength properties ( Y T =
100 Pa) for two di ff erent impact geometries. R p =
100 m. The minimal specific energy required to cause a significantshape change is estimated as Q bil ∼ ± / kg. µ = Q ∗ D values for the weak, highly porous bodies consideredhere are slightly higher than the specific energies Q ∗ D , ice found for solid bodies made of strong ice (Benz & Asphaug 1999) (Fig-ure 8). This result reflects the dissipative properties of material -1
10 100 1000 10 -3 -2 -1 Q c r i t ( J / K g ) M a x . g l oba l t e m pe r a t u r e i n c r ea s e d T m a x ( K ) v (m/s) Q *D (this study), Y T = 10 Pa (fit)Q *D (this study), Y T = 100 Pa (fit)Q *D (this study), Y T = 1000 Pa (fit)Q *D (this study), u = 0.42 (scaling)Q reshape (this study), Y T = 10 Pa (scaling)Q reshape (this study), Y T = 100 Pa (scaling)Q reshape (this study), Y T = 1000 Pa (scaling) Fig. 7.
Critical specific impact energies Q crit . The energy thresholds forshape-changing impacts on a 67P / C-G-like shape ( Q reshape ) for di ff er-ent material strength are shown, as well as the catastrophic disruptionenergies Q ∗ D for various impact velocities. We note that the Q bil valuesfound for shape-changing collisions on generic bi-lobe shapes overlapthe results for Q reshape with Y T = dT shownon the right y-axis is estimated by assuming that all kinetic impact en-ergy is converted into internal energy: dT = Q crit / c p where a constantheat capacity c p =
100 J / kg / K is used. porosity and is consistent with previous studies (e.g. Jutzi et al.2010).Also shown in Figure 8 is the value of Q ∗ D suggested by Lein-hardt & Stewart (2009) for weak ice as well as Q ∗ D predictedfrom scaling laws for collisions between gravity-dominated bod-ies (Leinhardt & Stewart 2012). In these studies, the e ff ects ofmaterial porosity were not taken into account. Article number, page 5 of 13 & A proofs: manuscript no. Jutzi_etal_final_corrected -1
1 10 100 1000 10 -3 -2 -1 Q , Q c r i t ( J / K g ) M a x . g l oba l t e m pe r a t u r e i n c r ea s e d T m a x ( K ) v (m/s) Q for bodies of similar size colliding at 40 m/s (Davidsson et al., 2016)Q for bodies of similar size colliding at 1 m/s (Jutzi&Asphaug, 2015)Q *D (Benz&Asphaug, 1999; solid ice, ρ = 440 kg/m )Q *D (Leinhardt&Stewart, 2009; weak ice, ρ = 440 kg/m )Q *D (Leinhardt&Stewart 2012; "small bodies", c* = 5, µ = 0.37)Q *D (this study; µ = 0.42)Q reshape (this study; µ = 0.42, Y T = 10 Pa)Q reshape (this study; µ = 0.42, Y T = 100 Pa)Q reshape (this study; µ = 0.42, Y T = 1000 Pa) Fig. 8.
Comparison of critical specific impact energies Q crit . The scalinglaws shown Figure 7 are compared here with Q ∗ D values found in previ-ous studies (Benz & Asphaug 1999; Leinhardt & Stewart 2009, 2012).Also displayed are the specific energies Q of collisions involving bodiesof similar size (mass ratio of 1:2) for the bi-lobe forming collisions instudy by Jutzi & Asphaug (2015) with very low velocities ( V ∼ V esc ∼ / s) as well as for collisions with a velocity of V =
40 m / s, corre-sponding to the average random velocity in the primordial disk duringthe first 25 Myr in the model by Davidsson et al. (2016). Finally, we also display in Figure 8 the specific energies Q involved in collisions of bodies of similar size (mass ratio 1:2)in the bi-lobe forming low-velocity regime investigated by Jutzi& Asphaug (2015). As expected, those low-velocity ( V ∼ V esc )accretionary collisions have specific impact energies far belowthe disruption threshold. For reference, we also show the specificenergy for collisions with much higher velocities (v =
40 m / s),which correspond to the average random velocity in the initialprimordial disk in the model by Davidsson et al. (2016). For amass ratio of 1:2, the specific impact energies are even aboveenergy threshold for catastrophic disruptions Q ∗ D . The results obtained in the previous section allow us to derivea Q ∗ D scaling law for porous cometary nuclei, which is a func-tion of impact velocity V and target size R (Housen & Holsapple1990): Q ∗ D = aR µ V − µ (2)where µ and a are scaling parameters.For Q reshape and Q bil , we use a fixed target size R = µ = Q ∗ D , Q reshape and Q bil are given in Table 2. The e ff ects of the impacts considered in this study (shape-changing impacts as well as catastrophic disruptions) are ana-lyzed in terms of impact heating and porosity evolution (below).First, in order to get an idea of the maximal global heating, we Table 2.
Parameters (SI units) for the scaling law Q crit = aR µ V − µ ,where R is the target radius and V the impact velocity. The scalingfor shape-changing impacts on 67P / C-G ( Q reshape ) and for impacts ongeneric bi-lobe shapes ( Q bil ) only hold for a fixed size ( R = Scaling µ aQ ∗ D Q reshape (10 Pa) 0.42 9.0e-7 Q reshape (100 Pa; nominal) 0.42 2.5e-6 Q reshape (1000 Pa) 0.42 3.8e-6 Q bil (nominal) 0.42 3.8e-6 -4 -3 -2 -1 -3 -2 -1 T e m pe r a t u r e i n c r ea s e d T ( K ) Mass fractionY T = 100 Pa, R p = 100 m, v = 200 m/sY T = 100 Pa, R p = 300 m, v = 20 m/s Fig. 9.
Cumulative post-impact temperature increase dT for two specificcases of shape-changing collisions, as indicated in the plot. -4 -3 -2 -1 -3 -2 -1 T e m pe r a t u r e i n c r ea s e d T ( K ) Mass fractionCatastrophic disruption, v = 600 m/sCatastrophic disruption, v = 200 m/sCatastrophic disruption, v = 100 m/sCatastrophic disruption, v = 10 m/s
Fig. 10.
Fraction of material in the final body that experienced a tem-perature increase larger than a certain value dT in catastrophic dis-ruptions with di ff erent velocities V . The mass of the largest remnant M lr / M t ∼ simply convert the total specific impact energy into a global tem-perature increase dT = Q crit / c p where a constant heat capac-ity c p =
100 J / kg / K is used. The value of c p is a rough massweighted average of the heat capacity of ice (Klinger 1981) andsilicates (Robie & Hemingway 1982) at low temperatures ( ∼ dT values cor-responding to collisions with a given specific impact energy. Article number, page 6 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G?
From this simple estimation, it is already obvious that im-pacts with energies comparable to Q reshape , the maximal globaltemperature increase must be limited to small values ( dT << dT distribution for a fewspecific cases of the shape-changing (section 2.4.1) as well ascatastrophic collisions (section 2.4.2). In the later, we only con-sider the material which ends up in the largest remnant ( ∼ ∼
50% of the original target mass)is not a ff ected much be the collision (Figure 10) and the heatingis limited ( <
1% of the mass is heated by dT > / s). Porosity is changed by impacts in multiple ways. First, materialis compacted due to the pressure wave generated by the impact.On the other hand, material is ejected and the process of reac-cumulation of the gravitationally bound material can give rise toadditional macroporosity. Our porosity model computes the de-gree of compaction (change of the distention variable). In orderto specify the increase of the macroporosity, we treat each SPHparticle individually according to its ejection / reaccumulationhistory. Particles which are lifted o ff the surface or are ejectedand reaccumulated experience a density decrease, resulting in anincrease of porosity. We assume that reaccumulated material canlead to the addition of macroporosity of maximal 40%, a typicalporosity of rubble-pile asteroids (Fujiwara et al. 2006). To com-pute the total final porosity φ total resulting from compaction andreaccumulation, we use the relation φ total = − /α total (3)and define the distention α total = min( ρ compact /ρ min , α max ) (4)where ρ min is the minimal density reached by the SPH parti-cle and ρ compact = / m . For this calculation we considerall particles which are gravitationally bound to the main body(largest remnant). The upper limit of the distention is given by α max = α α v (5)where α v is the distention value corresponding to 40% macrop-orosity, α v = (1 − φ v ) − with φ v = .
4, and α = ff the sur-face / ejected by the impact. Due to the addition of macroporos-ity resulting from reaccumulation, the final average porosity isabout the same as the initial porosity (Figure 11).In the catastrophic disruptions, most of the material whichundergoes collisional induced compaction does not remain on P o r o s i t y Mass fractionY T = 100 Pa, R p = 100 m, v = 200 m/s, including macro-porosityY T = 100 Pa, R p = 100 m, v = 200 m/s, compaction onlyY T = 100 Pa, R p = 300 m, v = 20 m/s, including macro-porosityY T = 100 Pa, R p = 300 m, v = 20 m/s, compaction onlyInitial porosity Fig. 11.
Post-impact porosity distribution for two specific cases ofshape-changing collisions, as indicated in the plot. The porosity calcu-lation takes into account compaction as well as the increase of macrop-orosity. For comparison, the porosity distributions resulting from com-paction only are shown as well. The final average porosity (compactionplus addition of macroporosity by reaccumulation) is 78.8% ( R p = R p =
300 m), respectively, while the initial porosity was77.8%. P o r o s i t y Mass fractionCatastrophic disruption, v = 800 m/s, including macro-porosityCatastrophic disruption, v = 800 m/s, compaction onlyCatastrophic disruption, v = 200 m/s, including macro-porosityCatastrophic disruption, v = 200 m/s, compaction onlyCatastrophic disruption, v = 100 m/s, including macro-porosityCatastrophic disruption, v = 100 m/s, compaction onlyCatastrophic disruption, v = 10 m/s, including macro-porosityCatastrophic disruption, v = 10 m/s, compaction onlyInitial porosity
Fig. 12.
Same as Figure 11 but for catastrophic collisions. M lr / Mtot ∼ M lr is considered.The final average porosity (compaction plus addition of macroporosityby reaccumulation) is 83.3% ( V =
10 m / s), 83.3% ( V =
100 m / s), 83.3%( V =
200 m / s), 82.4% ( V =
800 m / s), respectively, while the initialporosity was 77.8%. the main body (largest remnant). As a result, only ∼
10% of thematerial in the final main body has experienced significant com-paction. On the other hand reaccumulation plays a major role inthis collision regime, resulting in a significant increase of macro-porosity. The final porosity is therefore even slightly higher thanthe initial porosity (Figure 12).In Paper II, the interior porosity distribution of bi-lobe struc-tures resulting from sub-catastrophic collisions are compared toobservations of comet 67P / C-G.
Article number, page 7 of 13 & A proofs: manuscript no. Jutzi_etal_final_corrected
3. The combined dynamical and collisionalevolution of comet 67P/C-G
We follow the approach described in Morbidelli & Rickman(2015) in order to combine the dynamical evolution of the plan-etesimals precursors of Jupiter family comets with their colli-sional evolution. We do not repeat here a detailed description ofthe procedure, but we report on the di ff erences and the improve-ments in the new implementation.These are of three kinds. First, we consider here only thedynamical dispersal of the original trans-Neptunian disk of plan-etesimals, which generates the Scattered Disk (the current sourcereservoir of JFCs). Thus, we neglect the phase ranging from thetime when the gas was removed from the protoplanetary diskto the time when the giant planets developed a dynamical insta-bility that dispersed the planetesimal disk (Tsiganis et al. 2005;Gomes et al. 2005). This choice is made because Morbidelli &Rickman (2015) already showed that in the standard model, acomet the size of 67P / C-G has no chance to survive intact duringthis phase, if protracted for ∼
400 My. On the other hand the de-bate on the timing of the giant planet instability is still open (seefor instance Kaib & Chambers 2016; Toliou et al. 2016), so itmight be possible that the aforementioned phase is short. Thereis no doubt, however, that the dispersal of the trans-Neptuniandisk occurred and that this formed the Scattered Disk. In thiscase, Morbidelli & Rickman (2015) showed that during this pro-cess the number of catastrophic collisions for planetesimals thesize of 67P / C-G is ∼
1, so there might be some objects escapingbreak-up events. Thus, in this work we focus on this case, usingimproved assessments on the specific energies for catastrophicbreak-up and for reshaping, described in the previous sections.The second improvement over Morbidelli & Rickman (2015)concerns the dynamical simulations. Morbidelli & Rickman(2015) used the simulation of Gomes et al. (2005), which cov-ered only the first 350 My after the giant planet instability. Thisis when most of the action happens, but the subsequent 3.5-4.0 Gy cannot be neglected. Morbidelli & Rickman (2015) as-sumed that over this remaining time the orbital distribution ofthe Scattered Disk does not evolve any more, but its populationdecays exponentially down to 1% of the original population after4 Gy. The 1% fraction comes from previous studies of the longterm evolution of the Scattered Disk (Duncan & Levison 1997).Here we use the simulations presented in Brasser & Morbidelli(2013), which constitute a much more coherent set. Brasser &Morbidelli (2013) studied the dispersal of the trans-Neptunianplanetesimal disk during the giant planet instability using a largenumber of particles (1,080,000; including clones). At the end ofthe instability, they drove the giant planets towards their exactcurrent orbits, so to avoid artefacts in the subsequent long-termevolution of the Scattered Disk. The evolution of the ScatteredDisk was followed for 4 Gy. Because the number of active parti-cles decays over time, the test particles have been cloned 3 times,at 0.2, 1.0 and 3.5 Gy. In the final 0.5 Gy simulation, the parti-cles leaving the Scattered Disk to penetrate into the inner solarsystem as JFCs have been followed, in order to compare theirorbital distribution with that of the observed comets. This finalstep is crucial to demonstrate that the Scattered Disk generatedfrom the dispersal of the trans-Neptunian disk by the giant planetinstability is a valid source of JFCs.The third improvement over Morbidelli & Rickman (2015)is that the collisional evolution is followed only for the parti-cles that eventually become JFCs in the final 0.5 Gy simulation. number of disruptive collisions number of JFC q=-2.5q=-3.0q=-3.5
Fig. 13.
The number of expected catastrophic collisions N disrupt duringthe formation and evolution of the Scattered Disk for the particles thateventually become JFCs in the final 0.5 Gy simulation. N disrupt is com-puted using the scaling parameters for our new Q ∗ D . The symbols depictdi ff erent values for the exponent of the di ff erential size distribution q ,as labeled in the plot. These are 87 particles. We think that, potentially, this is an im-portant improvement. The particles that penetrate the inner solarsystem at the present time might have had specific orbital histo-ries relative to the other particles that either became JFCs muchearlier or are still in the Scattered Disk today. Averaging the col-lisional histories of these three categories of particles may not besignificant to address the specific case of 67P / C-G, which clearlybecame JFC only in recent time.Like in Morbidelli & Rickman (2015) the number of colli-sions su ff ered by each considered body is computed assumingthat the initial disk particles represent a population of 2 × planetesimals with diameter D > . ff erential sizedistribution characterized by an exponent q . The minimum pro-jectile size is determined by the scaling law (equation 2) for thecritical specific energy, with parameters given in Table 2. As forthe exponent q , in agreement with Morbidelli & Rickman (2015)and previous studies of the comet size distribution, we considerhere the cases with q = − . , − . − .
5. However, in themeantime the New Horizons mission to Pluto and Charon hasmeasured the crater size frequency distribution, allowing the as-sessment of the size distribution of the trans-Neptunian objectslarger than a few km in diameter (Singer et al. 2015). The pre-liminary results suggest q = − .
3. Thus, we consider the re-sults for q = − . − . × with D > . ffi cient sources of the LPC andJFC fluxes that we observe today. We note that based on the most recent results it has been suggestedthat there may be a deficit of small objects (Singer et al. 2016); seediscussion in section 4.Article number, page 8 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G? number of reshaping collisions ( σ T =10 Pa) number of JFC q=-2.5q=-3.0q=-3.5 number of reshaping collisions ( σ T =100 Pa) number of JFC q=-2.5q=-3.0q=-3.5 number of reshaping collisions ( σ T =1 kPa) number of JFC q=-2.5q=-3.0q=-3.5 Fig. 14.
Same as Figure 13, but shown is the number of shape-changingcollisions on a 67P / C-G-like body N reshape , computed using scaling theparameters for Q reshape for di ff erent strengths. We note that the numberof shape-changing collisions N bil in the case of a generic bi-lobe shapewith nominal strength properties is the same as N reshape for Y T = The number of events for each particle surviving in the ScatteredDisk at the end of the disk dispersal simulation is shown in Figs.13 - 14. The results for the various types of collisions, using thecorresponding scaling laws ( Q ∗ D and Q reshape ), are plotted. We C u m u l a t i v e f r a c t i on Number of collisions
Disruption, q=-2.5Disruption, q=-3.0Disruption, q=-3.5Shape-change, 1000 Pa, q=-2.5Shape-change, 1000 Pa, q=-3.0Shape-change, 1000 Pa, q=-3.5 Shape-change, 100 Pa, q=-2.5Shape-change, 100 Pa, q=-3.0Shape-change, 100 Pa, q=-3.5Shape-change, 10 Pa, q=-2.5Shape-change, 10 Pa, q=-3.0Shape-change, 10 Pa, q=-3.5
Fig. 15.
Cumulative fraction of particles (that eventually become JFCs)as a function of the number of collisions. This is an alternative represen-tation of the results already shown in Figs. 13 and 14. The solid linescorrespond to the Q ∗ D scaling; the dotted lines were computed using Q reshape (for di ff erent strength values). -4 -3 -2 -1 c u m u l a t i v e f r a c t i on probability to miss all disruptive/reshaping collisionsq=-2.5 (disruptive)q=-3.0 (disruptive)q=-3.5 (disruptive)q=-2.5 ( σ T = 1000 Pa)q=-2.5 ( σ T = 100 Pa)q=-2.5 ( σ T = 10 Pa) Fig. 16.
Cumulative fraction of particles (that eventually become JFCs)as a function of the probability P (0) to escape all collisions. The dif-ferent line styles refer to di ff erent exponents for the di ff erential sizedistribution q , as labeled on the plot. The three curves on the right cor-respond to the Q ∗ D scaling; the three curves on the left correspond to Q reshape (with di ff erent strength values 10 Pa, 100 Pa and 1 kPa fromleft to right). For q = − . q = − .
5, the probability to miss allreshaping collisions is P (0) << − and the corresponding curves arenot displayed here. note that the results for Q bil (impacts on generic bi-lobe shapeusing nominal material properties) are the same as in the caseof Q reshape with Y T = Q ∗ D scaling law used here leadsto disruption energies which are higher than the ones by Benz& Asphaug (1999) (which were used in the previous study). Asdiscussed in section 2.4.2, this can be explained by the highlydissipative properties of porous materials, which are taken intoaccount in the new Q ∗ D . As Figure 13 shows, for shallow sizedistributions, it is possible in principle that a significant fractionof the 67P / C-G sized objects escaped all catastrophic collisions.
Article number, page 9 of 13 & A proofs: manuscript no. Jutzi_etal_final_corrected -5 -4 -3 -2 -1
0 5x10 Mean number of expected collisions ( σ T = 10Pa) Time (years) q=-2.5q=-3.0q=-3.5 -5 -4 -3 -2 -1
0 5x10 Mean number of expected collisions ( σ T = 100Pa) Time (years) q=-2.5q=-3.0q=-3.5 -5 -4 -3 -2 -1
0 5x10 Mean number of expected collisions ( σ T = 1000Pa) Time (years) q=-2.5q=-3.0q=-3.5
Fig. 17.
Mean number of reshaping collisions N reshape expected for67P / C-G-like objects as a function of time for di ff erent strengths, asindicated on the y -axis. We note that the number of shape-changing col-lisions N bil in the case of a generic bi-lobe shape with nominal strengthproperties is the same as N reshape for Y T = t = t = × years is now. On the other hand, the number of shape-changing collisions(Figure 14), requiring a much smaller impact energy ( Q reshape ),is substantially larger than the number of catastrophic events.As expected, the weaker the strength the larger the number ofreshaping collisions taking place. Also, the steeper the size dis- -5 -4 -3 -2 -1
0 5x10 Mean number of expected collisions
Time (years) q=-2.5q=-3.0q=-3.5 -5 -4 -3 -2 -1
0 5x10 Mean number of expected collisions
Time (years) q=-2.5q=-3.0q=-3.5
Fig. 18.
Top: Mean number of potential 67P / C-G-forming catastrophiccollisions of a parent body of R = Q ∗ D ) as a func-tion of time t (defined as in Figure 17). Bottom: Same, but for the sce-nario of 67P / C-G formation by low energy sub-catastrophic collisions. tribution (larger q ), the larger the number of collisions happen-ing. However, in any case, even for the largest strength (1000Pa) and the shallowest slope ( q = -2.5), the number of reshapingcollisions largely exceeds 1 for all comets.The results are summarized in Figure 15 which shows the cu-mulative fraction of particles as a function of the number of col-lisions. In Figure 16, the number of collisions N coll is convertedinto a probability to avoid all collisions P (0) = exp( − N coll ) andthe normalized cumulative distribution of the P (0) values is plot-ted. The average number of collisions and the related probabili-ties are given in Table 3.It is also interesting to look at the number of collisions as afunction of time (Figs. 17 and 18) as this in principle allows us todetermine the time at which on average the last event of a certaintype took place.For size distributions with q ≤ − .
0, the last shape-changingevent (on average) would have taken place in the last 1 Gy (Fig-ure 17), suggesting that the structure of comet 67P / C-G musthave formed in a recent period.In Figure 18, we plot the average number of events as a func-tion of time for two potential formation scenarios. In the firstscenario, it is assumed that the structure of 67P / C-G formed asa result of a catastrophic break-up of a parent body of R = Q ∗ D scaling. In the second case, we consider im- Article number, page 10 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G?
Table 3.
Average number of shape-changing collision on a 67P / C-G-like object ( N reshape ), shape-changing collisions on a generic bi-lobe body( N bil ) and catastrophic collisions ( N disrupt ). The corresponding probability P (0) to avoid all collisions is given in parenthesis. Type q = -2.5 q = -3 q = -3.5 N disrupt N reshape (10 Pa) 4.92 (7.3E-3) 35.1 (6.0E-16) 258 (7.3E-113) N reshape (100 Pa) 3.06 (4.7E-2) 18.1 (1.4E-8) 112 (2.4E-49) N reshape (1000 Pa) 2.53 (8.0E-2) 13.8 (1.0E-6) 79.6 (2.7E-35) N bil (nominal) 2.53 (8.0E-2) 13.8 (1.0E-6) 79.6 (2.7E-35)pact energies corresponding to a scenario of 67P / C-G formationby low energy sub-catastrophic collisions, as proposed in PaperII. Clearly, the number of events in the later case are substantiallylarger. This suggests that it may be a more probable formationmechanism than the catastrophic break-up scenario (see a moredetailed discussion on this topic in Paper II).
4. Uncertainties and alternative models
In this section we discuss some aspects of the robustness anduncertainties of our modeling approach and alternative models.
The values for the specific catastrophic disruption energies Q ∗ D are well defined and follow the expected scaling (Figure 7). Thecritical specific impact energies for reshaping are not as well de-fined and do depend on the material properties. However, weexplore a reasonably large range of material properties and alsoapply large error bars to the results in this case. In any case, thereis no doubt that Q reshape << Q ∗ D and consequently, there must bemany more shape-changing events than catastrophic disruptions. A crucial quantity in the dynamical model is the initial numberof comets. The assumption of the existence of 2 × cometsis in line with estimates of the current Scattered Disk and Oortcloud populations and numerical estimates of the fractions of theprimoridal disk that end up in these populations. Both could bewrong, in principle. However the fractions of the primordial diskpopulation implanted in the Scattered Disk and Oort cloud thatwe use (from Brasser & Morbidelli (2013)) are not very di ff er-ent from those found in quite di ff erent dynamical models (Doneset al. (2004) for the Oort cloud to Duncan and Levison et al.(2008), for the Scattered Disk). Therefore, they seem to be ro-bust.The number of comets used in our model are based also ona flux of Jupiter family comets which is assumed to be currentlyin a steady state. If this is not the case, the Scattered Disk couldbe less (or more) populated than predicted by the model. How-ever, we find this unprobable for the following reason. The cur-rent estimates for the populations in the Scattered Disk and theOort cloud are consistent with these two reservoirs being gen-erated from the same parent disk (Brasser & Morbidelli 2013).Thus, if the Jupiter family comet flux is now - say 10 × - themean flux (so to argue for a Scattered Disk 10 × less populated),the same should apply for the flux of long period comets. Butthe fluxes out of Scattered Disk and Oort cloud follow di ff erentprocesses: for the Scattered Disk, this is resonant di ff usion andscattering from Neptune; for the Oort cloud it is stellar perturba- tions. Therefore, it seems unlikely that both fluxes increased bythe same amount relative to the mean values.Another crucial quantity in our modeling is the slope of thesize distribution q , which determines the number of projectiles ofa given size and thus the number of impacts with energies abovethe critical value. There is an ongoing debate about the form ofthe size distribution in the Scattered Disk population. We arguethat the observations of the crater size distributions in the Plutosystem by the New Horizons mission provides one of the bestavailable constraints. The cratering of Pluto and Charon is dom-inated by the hot population (Greenstreet et al. 2015). All modelsagree that the hot population and the Scattered Disk populationare the same population in terms of physical properties and ori-gin. In fact, the collisional evolution of the hot population is notmore severe than that of the Scattered Disk. Both su ff ered mostcollisions during the dispersal of the primitive disk (or before,if the dispersal was late). It is true that comets have a shallowerdistribution (Snodgrass et al. 2011) as well as have the craters onthe Jovian satellites (Bierhaus 2006; Bierhaus et al. 2009). Butthis is probably because small comets disintegrate very quickly.On the satellites of Saturn, the crater size distribution is similarto the one expected from a projectile population with a size dis-tribution like that of the main asteroid belt (e.g. Plescia & Boyce1982; Neukum et al. 2005, 2006), i.e. it is the same as measuredby New Horizons on Pluto and Charon.We note that based on the most recent analysis of the NewHorizons data, it has been suggested (Singer et al. 2016) thatthe size distribution for small ( < q (cid:39)− .
5) than at large scales. However, this result is still preliminarywith uncertainties to be clarified. As discussed above, the TNOsize distribution looks very similar to the size distribution of theasteroid belt, which is a result of a collisional equilibrium (below ∼
100 km). This suggests that the size distribution of TNOs is ina collisional equilibrium as well. A change of slope below 2 kmwould produce waves in the TNO size distribution above 2 km.This is not observed, which may indicate that the change of slopeis not as pronounced.To check the e ff ects of a varying slope on our results, we per-formed additional calculations using q = − . > q s for small ( < q s (cid:38) -2, which meansthat if indeed q s = -1.5, a 67P / C-G-like shape would survive.However, we reiterate that this calculation considers a conserva-tive scenario without any collisional evolution in the primordialdisk.
Alternative models to the standard model such as suggested byDavidsson et al. (2016) predict a much smaller collisional evo-lution and are consistent with the idea of comets being primi-tive unprocessed objects, formed primordially. However, these
Article number, page 11 of 13 & A proofs: manuscript no. Jutzi_etal_final_corrected models require the number of objects in the Scattered Disk tobe orders of magnitude smaller. We note that there is no directobservational measure of the Scattered Disk population and allestimates are indirect and pass through models, so such a smallnumber can in principal not be excluded.In is not clear, however, how bi-lobe structures wouldform / survive in these models. Previous studies indicate that theprimordial formation of bi-lobed shapes, such as the one ofcomet 67P / C-G, by direct merging requires extremely low col-lision velocities of V / V esc ∼ V = / s during the first 25 Myrs. For kilometer-sized bodies this im-plies a ratio V / V esc ∼
40! In fact, the corresponding specific im-pact energies are larger than the catastrophic disruption thresh-old (Figure 8). Our results show that even relative velocities ofa few m / s are destructive and lead to reshaping (Figures 3- 5).Therefore, it is unlikely that primordial bi-lobe structures wouldsurvive this phase, and at the same time their formation by colli-sional merging is implausible due to the high relative velocities.
5. Summary and conclusions
We have estimated the number of shape-changing collisions foran object with a shape like comet 67P / C-G, considering a dy-namical evolution path typical for a Jupiter family comet, usinga "standard model" of the early solar system dynamics.First, we computed the e ff ects of impacts on comet 67P / C-Gusing a state-of-the-art shock physics code, investigating rangeof impact conditions and material properties. We found that theshape of comet 67P / C-G, with two lobes connected by a neck,can be destroyed easily, even by impacts with a low specific im-pact energy. From these results, scaling laws for the specific en-ergy required for a significant shape alteration ( Q reshape ) weredeveloped. For more general applications, the critical specificenergies to alter the shape of generic bi-lobe objects ( Q bil ) wasinvestigated as well.These scaling laws for Q reshape and Q bil were then used toanalyze the dynamical evolution of a 67P / C-G-like object andgeneric bi-lobe shapes in terms of shape-changing collisions.We considered a conservative scenario without any collisionalevolution before the dynamical instability of the giant planets.Rather, we tracked the collisions during the dispersion of thetrans-Neptunian disk caused by the giant planet instability, theformation of a scattered disk of objects and the penetration oftens of objects into the inner solar system. To do this we useda set of simulations (Brasser & Morbidelli 2013) that producesorbits consistent with the observed JFC population.We find that even in this conservative scenario, comet 67P / C-G would have experienced a significant number of shape-changing collisions, assuming that its structure formed primor-dially. For size distributions with q ≤ − .
0, the last reshapingevent (on average) would have taken place in the last 1 Gy. Thepreliminary results of the New Horizons missions concerningthe crater size-frequency distribution on Pluto and Charon sug-gest that the current trans-Neptunian population (i.e. includingthe Scattered Disk) has a di ff erential power-law size distributionwith an exponent q (cid:39) − . < ∼ ffi cient to destroy their two-componentstructure. This strongly suggests that the two-component bodywhich is required to exist at the beginning of the fission-mergingcycle cannot be primordial.Thus, according to our model, comets are not primordial inthe sense that their current shape and structure did not form inthe initial stages of the formation of the solar system. Rather,they evolve through the e ff ect of collisions and the final shape isa result of the last major reshaping impact, possibly within thelast 1 Gy. A scenario of a late formation of 67P / C-G-like two-component structures is presented in Paper II.It is clear that the results presented here are based on the as-sumption that the standard model of dynamical evolution is cor-rect. Although some of its parameters are debated, as discussedin section 4, we believe that the model is robust. We note that itis so far the only model which produces the correct number ofobjects in the inner solar system with orbits consistent with theobserved JFC population.Our results clearly show that if this standard model of solarsystem dynamics is correct, it means that the cometary nucleias they are observed today must be collisionally processed ob-jects. Therefore, the remaining important question is whether ornot such collisionally processed bodies can still have primitiveproperties (i.e. high porosity, presence of supervolatiles). If thisis not the case, then the standard model must be wrong. Thiswould mean for instance that either the primordial disk was dy-namically cold and contained a much lower number of objects,as proposed by Davidsson et al. (2016) or that there is a lack ofsmall comets, implying an abrupt change in the slope of the sizedistribution.However, the analysis of the outcomes of the detailed im-pact modeling carried out here (for shape-changing impactsand catastrophic disruptions) suggest that collisionally processedcometary nuclei can still have a high porosity, and could haveretained their volatiles, since there is no significant large-scaleheating. Therefore, they may still look primitive, meaning thatthe standard model is consistent with the observations of comet67P / C-G. This question is investigated further in Paper II andalso in an ongoing study of bi-lobe formation in large-scalecatastrophic disruptions (Schwartz et al., 2016, in prep).Primordial or not, the structure of comet 67P / C-G is an im-portant probe of the dynamical history of small bodies.
Acknowledgements.
M.J. and W.B. acknowledge support from the Swiss NCCRPlanetS. A.T. wishes to thank OCA for their kind hospitality during her staythere. We thank the referees B. Davidsson and J. A. Fernandez for their thoroughreview which helped improve the paper substantially.
References
Asphaug, E. & Benz, W. 1994, Nature, 370, 120Benz, W. & Asphaug, E. 1995, Computer Physics Communications, 87, 253Benz, W. & Asphaug, E. 1999, ICARUS, 142, 5Biele, J., Ulamec, S., Maibaum, M., et al. 2015, Science, 349, 9816Bierhaus, E. B. 2006, in LPI Contributions, Vol. 1320, Workshop on SurfaceAges and Histories: Issues in Planetary Chronology, 14–15Bierhaus, E. B., Zahnle, K., & Chapman, C. R. 2009, Europa’s Crater Distri-butions and Surface Ages, ed. R. T. Pappalardo, W. B. McKinnon, & K. K.Khurana, 161Boehnhardt, H. 2004, Comets II, University of Arizona Press., 301
Article number, page 12 of 13. Jutzi , W. Benz, A. Toliou, A. Morbidelli, R. Brasser: How primordial is the structure of comet 67P / C-G?
Brasser, R. & Morbidelli, A. 2013, Icarus, 225, 40Capaccioni, F., Coradini, A., Filacchione, G., et al. 2015, Science, 347, aaa0628Cuzzi, J. N., Hogan, R. C., & Bottke, W. F. 2010, ICARUS, 208, 518Davidsson, B. J. R., Sierks, H., Güttler, C., et al. 2016, A&A, 592, A63Dones, L., Weissman, P. R., Levison, H. F., & Duncan, M. J. 2004, Oort cloudformation and dynamics, ed. M. C. Festou, H. U. Keller, & H. A. Weaver,153–174Duncan, M. J. & Levison, H. F. 1997, Science, 276, 1670Fujiwara, A., Kawaguchi, J., Yeomans, D. K., et al. 2006, Science, 312, 1330Gomes, R., Levison, H. F., Tsiganis, K., & Morbidelli, A. 2005, Nature, 435, 466Greenstreet, S., Gladman, B., & McKinnon, W. B. 2015, Icarus, 258, 267Güttler, C., Krause, M., Geretshauser, R. J., Speith, R., & Blum, J. 2009, ApJ,701, 130Hässig, M., Altwegg, K., Balsiger, H., & etal. 2015, Science, 347Hirabayashi, M., Scheeres, D. J., Chesley, S. R., et al. 2016, Nature, 534, 352Housen, K. R. & Holsapple, K. A. 1990, Icarus (ISSN 0019-1035), 84, 226Housen, K. R. & Holsapple, K. A. 2003, Icarus, 163, 102Johansen, A., Oishi, J. S., Mac Low, M.-M., et al. 2007, Nature, 448, 1022Jutzi, M. 2015, Planetary and Space Science, 107, 3Jutzi, M. & Asphaug, E. 2015, Science, 348, 1Jutzi, M., Benz, W., & Michel, P. 2008, Icarus, 198, 242Jutzi, M., Holsapple, K. A., Wünnemann, K., & Michel, P. 2015, Asteroids IV, 1Jutzi, M., Michel, P., Benz, W., & Richardson, D. C. 2010, Icarus, 207, 54Kaib, N. A. & Chambers, J. E. 2016, MNRAS, 455, 3561Kataoka, A., Tanaka, H., Okuzumi, S., & Wada, K. 2013, A&A, 557, L4Klinger, J. 1981, Icarus, 47, 320Kofman, W., Herique, A., Barbin, Y., et al. 2015, Science, 349Leinhardt, Z. M. & Stewart, S. T. 2009, Icarus, 199, 542Leinhardt, Z. M. & Stewart, S. T. 2012, The Astrophysical Journal, 745, 79Levison, H. F., Morbidelli, A., Van Laerhoven, C., Gomes, R., & Tsiganis, K.2008, Icarus, 196, 258Lorek, S., Gundlach, B., Lacerda, P., & Blum, J. 2016, A&A, 587, A128Massironi, M., Simioni, E., Marzari, F., et al. 2015, Nature, 526, 402Melosh, H. J. 1989, Impact cratering: A geologic processMorbidelli, A., Bottke, W. F., Nesvorný, D., & Levison, H. F. 2009, ICARUS,204, 558Morbidelli, A., Marchi, S., Bottke, W. F., & Kring, D. A. 2012, Earth and Plane-tary Science Letters, 355, 144Morbidelli, A. & Rickman, H. 2015, A&A, 583, A43Mumma, M. J., Weissman, P. R., & Stern, S. A. 1993, In: Protostars and planetsIII (A93-42937 17-90), 1177Neukum, G., Wagner, R., Wolf, U., & Denk, T. 2006, in European PlanetaryScience Congress 2006, 610Neukum, G., Wagner, R. J., Denk, T., Porco, C. C., & Cassini Iss Team. 2005, inLunar and Planetary Inst. Technical Report, Vol. 36, 36th Annual Lunar andPlanetary Science Conference, ed. S. Mackwell & E. StansberyNy ff eler, B. 2004, PhD thesis, University of BernPätzold, M., Andert, T., Hahn, M., et al. 2016, Nature, 530, 63Plescia, J. B. & Boyce, J. M. 1982, Nature, 295, 285Rickman, H., Marchi, S., A’Hearn, M. F., et al. 2015, A&A, 583, A44Robie, A. & Hemingway, B. 1982, American Mineralogist, 67, 470Rotundi, A., Sierks, H., Della Corte, V., et al. 2015, Science, 347, aaa3905Sierks, H., Barbieri, C., Lamy, P. L., et al. 2015, Science, 347, aaa1044Singer, K. N., McKinnon, W. B., Greenstreet, S., et al. 2016, in AAS / Division forPlanetary Sciences Meeting Abstracts, Vol. 48, AAS / Division for PlanetarySciences Meeting AbstractsSinger, K. N., Schenk, P. M., Robbins, S. J., et al. 2015, in AAS / Division forPlanetary Sciences Meeting Abstracts, Vol. 47, AAS / Division for PlanetarySciences Meeting Abstracts, 102.02Skorov, Y. & Blum, J. 2012, ICARUS, 1Snodgrass, C., Fitzsimmons, A., Lowry, S. C., & Weissman, P. 2011, MNRAS,414, 458Speith, R. 2006, Habilitation, University of TübingenStecklo ff , J. K., Johnson, B. C., Bowling, T., et al. 2015, Icarus, 258, 430Toliou, A., Morbidelli, A., & Tsiganis, K. 2016, A&A, 592, A72Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Nature, 435, 459Weidenschilling, S. J. 1997, ICARUS, 127, 290Weissman, P. R., Asphaug, E., & Lowry, S. C. 2004, Comets II, 337Windmark, F., Birnstiel, T., Güttler, C., et al. 2012a, A&A, 540, A73Windmark, F., Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2012b, A&A,544, L16Youdin, A. N. & Goodman, J. 2005, ApJ, 620, 459, J. K., Johnson, B. C., Bowling, T., et al. 2015, Icarus, 258, 430Toliou, A., Morbidelli, A., & Tsiganis, K. 2016, A&A, 592, A72Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Nature, 435, 459Weidenschilling, S. J. 1997, ICARUS, 127, 290Weissman, P. R., Asphaug, E., & Lowry, S. C. 2004, Comets II, 337Windmark, F., Birnstiel, T., Güttler, C., et al. 2012a, A&A, 540, A73Windmark, F., Birnstiel, T., Ormel, C. W., & Dullemond, C. P. 2012b, A&A,544, L16Youdin, A. N. & Goodman, J. 2005, ApJ, 620, 459