© ESO 2019

1 Introduction

As part of the ESA Rosetta mission (Glassmeier et al. 2007a), the lander Philae touched down on comet 67P/Churyumov-Gerasimenko (67P) on November 12, 2014 at 15:34:04 UTC. Philae rebounded three times before it finally came to a stop in a region called Abydos. One of the primary objectives of the lander mission was to measure the surface properties (such as the compressive strength) of 67P in situ to determine what comets are made of and where they formed. While the growth of kilometer-sized planetesimals to larger planet-sized objects in the early solar system is well understood, it is still under debate how smaller objects formed out of sub-decimeter-sized dust aggregates (see Blum 2018 for a recent review). Because of their primitive nature, comets are the best candidates for planetesimals and possibly the sole small survivors of the planet-formation era (Davidsson et al. 2016), although they might have experienced a few collisions since formation (Morbidelli & Rickman 2015). The mechanical strength of the cometary material strongly depends on the formation history and evolution of the nucleus (Blum et al. 2014). Thus, measurements of mechanical properties can be used to study the origin and evolution of comets, and can thereby provide the missing link between protoplanetary dust and planets. Two different processes are currently discussed as the main driver behind comet nucleus formation: gravitational instability combined with a gentle collapse of pebble-sized dust aggregates (Johansen et al. 2007; Blum et al. 2017), or hierarchical agglomeration (Weidenschilling 1997; Davidsson et al. 2016).

In the first scenario, comets are believed to have formed in the young solar system by the gentle gravitational collapse of dust clouds, typically consisting of sub-decimeter-sized aggregates (Johansen et al. 2007). The nature of this process means that the nucleus mostly consists of intact dust aggregates (Blum et al. 2017) that have survived the comet formation owing to low impact velocities of less than 1 m s⁻¹ during collapse and have thus experienced relatively little compression (Wahlberg Jansson & Johansen 2014; Wahlberg Jansson et al. 2017). In contrast, hierarchical agglomeration is based on growth by high-velocity impacts (up to ~50 m s⁻¹) onto the cometesimal to form the cometary nucleus. Compared to a gentle gravitational collapse, this scenario leads to an increased mechanical strength and volume filling factor that is due to the higher compaction of the material (Blum 2018). Hence, in situ measurements of the mechanical properties like the compressive strength allows testing these formation hypotheses. The nature of 67P is most likely primordial (Davidsson et al. 2016), which makes it an ideal candidate for these types of studies.

Several different estimates for the compressive strength σ have been published. Spohn et al. (2015) derived a uniaxial lower limit of 2 MPa based on the results of a spot measurement with the Philae MUPUS hammer experiment. Biele et al. (2015) and Roll et al. (2016) estimated a surface compressive strength in the order of 1 kPa based on an analysis of the mechanics of Philae’s first touchdown. By analyzing the collapse of cliff overhangs observed from orbit, Groussin et al. (2015) inferred a compressive strength of between 30 Pa and 150 Pa.

In this paper we derive an upper limit for the surface compressive strength by using the Philae lander as an impact probe. The compressive strength is determined based on the energy balance and geometry of Philae’s contacts. While the lander initially touched down in an area covered by fallback dust, the other contacts (in the following the contacts between Philae and 67P are called first touchdown, collision, second touchdown, and third touchdown, in line with Auster et al. 2015 and Biele et al. 2015) took place in more consolidated areas without significant fallback. Deriving the compressive strength based on the energy balance requires reasonable knowledge of the contact duration and corresponding geometry. Reconstructions of the first touchdown (Biele et al. 2015; Roll & Witte 2016; Roll et al. 2016) showed that Philae pivoted several times between the individual landing feet. As the exact sequence of these events can only be estimated based on simulations, it is impossible to derive the contact geometry and corresponding contact duration. Additionally, in contrast to the subsequent surface contacts, the landing gear damper system was engaged during the initial touchdown, actively dissipating some of the kinetic energy. In this study we therefore did not revisit the first touchdown. The energy balance and compressive strength was only derived for the collision, second touchdown, and third touchdown.

To determine the total kinetic energy balance (translational and rotational energy) for each of the surface contacts, it is necessary toknow the incoming and outgoing velocities and the rotation rate of the lander during flight. While the velocities can be derived from the corresponding trajectories, the rotational energy was determined based on the rotation frequencies provided by Heinisch et al. (2017a). In the following section, the locations of the individual touchdowns are discussed and used to derive the required velocities based on assumed ballistic flight parabolas, taking into account previous reconstructions by Auster et al. (2015), Biele et al. (2015), and Heinisch et al. (2017a). In a next step, attitude information gained from two-point magnetic field observations by the ROMAP (Auster et al. 2007) and RPC-MAG (Glassmeier et al. 2007b) magnetometers based on the work of Heinisch et al. (2016, 2017a) and Richter et al. (2016) is used in conjunction with a digital terrain model of 67P (Preusker et al. 2017) to determine the geometry of the contacts. This is necessary to calculate the contact pressure (force per contact area). The compressive strength is finally calculated based on the total energy balance, contact geometry, and duration of the surface contacts.

Next, we combine attitude information for the last part of the trajectory up until third touchdown with scratches visible in the images of the lander at Abydos. As these features were most likely created by Philae, they can be used to constrain the surface compressive strength as well. Finally, the implications of these results are discussed.

2 Trajectory and attitude reconstruction

As the first and final landing sites were photographed from orbit (see Fig. 1), the respective positions were determined based on these images. For the first touchdown at Agilkia, the Cartesian (x, y, z) coordinates (2.12, −0.96, 0.50) km given by Biele et al. (2015) were used, while for Abydos (final location of Philae), the estimates (2.45, −0.07, −0.34) km of Heinisch et al. (2017a) were used. All coordinates are given in the 67P body-fixed CHEOPS coordinate system (SPICE: 67P/C-G_CK). The NASA NAIF SPICE system (Acton 1996) was used to facilitate the coordinate calculations.

The trajectory between first touchdown and collision was constrained by the flight time, surface structure, and an image of Philae and its shadow in flight, taken approximately 10 min (15:45:02 UTC) after first touchdown (see Fig. 2a). Based on this image, Philae’s approximate position at 15:45:02 UTC was (2.27, −0.88, 0.36) km. By combining this information with the coordinates of first touchdown, the direction of flight was determined (see Fig. 2). The flight time t₁ = (2756 ± 3) s was derived based on the time between the initial touchdown and the collision. Assuming a free ballistic trajectory (similar to Auster et al. 2015) with a mean local gravitational acceleration of g ≈ 1.6 × 10⁻⁴ m s⁻² (based on Biele et al. 2015), a possible area (approximately 120 m × 90 m) for the collision was identified using the latest SHAP7 digital terrain model (Preusker et al. 2017). This area was selected based on the flight direction, taking into account the preliminary collision coordinates proposed by Heinisch et al. (2017a), the flight time, and uncertainties in the digital terrain model. As the aim of this analysis is to estimate the velocities and kinetic energy for the rebound, it is not necessary to determine the exact contact locations or definitive trajectories. To simplify computation and be independent of a specific version of the terrain model, the area of interest was approximated by three planes placed parallel to the local surface to reflect the different slopes of the Hatmehit crater rim.

Possible locations for second touchdown are constrained by the estimated position of the collision (taking uncertainties in the location into account) and the known final landing site. A cliff-like structure visible in the upper right-hand corner of Fig. 2 further restricts the trajectory, as Philae would not have been able to clear the top of this cliff (relative height ~60 m) on a ballistic flight coming from the collision site. Hence the site of second touchdown must be between this cliff and the location of the collision. This yields an estimated flight direction relative to third touchdown given by (0.15, −0.95, 0.27) km, as previously described by Heinisch et al. (2017a). Based on the flight time between second and third touchdown of only t₃ = (215 ± 3) s, the locations of the two touchdowns must be close together. A rough initial estimate for the horizontal velocity between the collision and second touchdown can be derived from the corresponding flight time t₂ = (3926 ± 3), preliminarily assuming the locations of second and third touchdown to be identical. When we combine this estimated velocity with the flight directions, we find the area highlighted in red in Fig. 2b. The center of this region is approximately at (2.45, −0.07, −0.34) km. Even though the available information is not enough to determine a definitive touchdown site, the derived area is sufficient for the purpose of this study.

Based on the derived touchdown positions, we calculated the velocities and corresponding kinetic energies for each possible part of the trajectory. Combined with the rotation rates determined from magnetic field observations, the total energy balance of the reboundwas derived. To facilitate this, the possible areas for the individual touchdowns were discretized with a spatial resolution of1 m in each direction and the trajectories for each of these points were calculated. To determine the horizontal and vertical velocity components during the ballistic flights, a simple horizontal motion with constant velocity v_h and a vertical motion upwards with a constant initial velocity, combined with a free fall were used. Because of the difference in relative height (Δh), the outgoing vertical velocity (v_vO) after contact can differ from the incoming vertical velocity (v_vI) of the following contact. The horizontal and vertical components of the velocity were calculated using $$v_h=\frac{d}{t_f},$$(1) $$v_{vO}=g{t}_f/2+\Delta h/t_f,$$(2)

where d denotes the horizontal distance and Δh the height difference between the touchdowns. The incoming vertical velocity was calculated based on the outgoing velocity of the previous surface contact: $$v_{vI}=v_{vO}-g{t}_f.$$(3)

The flight time t_f was determined based on the time of the individual contacts to be t₁ = (2756 ± 3) s for the flight between first touchdown and collision, t₂ = (3926 ± 3) s for the flight between collision and second touchdown, and t₃ = (215 ± 3) s for the last part of free flight after second touchdown. The rotational energy was calculated using $$E_{\text{rot}}=2{\text{π}}^2{f}_{\text{rot}}^2{I}_{Pl},$$(4)

with f_rot being the rotation frequency and I_PL denoting the moment of inertia, which is ~16.6 kg m² for rotation around the main body-axis (lander z-axis), ~ 15.9 kg m² for rotation around the front landing gear axis (lander x-axis), and ~14.6 kg m² along the axis parallel to the lander balcony (y-axis; Biele et al. 2015; Roll et al. 2016). During collision, the rotation pattern of the lander changed significantly, from a simple rotation along the body-axis (perpendicular to the landing gear plane) with a frequency of f₁ = (77.0 ± 1.0) mHz to a precession with a frequency of f₂ = (42.0 ± 1.0) mHz. After second touchdown, Philae had lost most of the rotational energy, spinning only with about f₃ = (2.0 ± 0.5) mHz. The flight times, the rotation rates, contact times, contact geometry, and lander attitude used in the following sections were all taken from Heinisch et al. (2016, 2017a). All these values are also summarized in Table 1. The resulting total velocity and energy distributions for collision and second touchdown are depicted in Fig. 3.

During the collision, no significant acceleration along the body-axis (Philae z-axis) perpendicular to the landing gear was detected, in contrast to the other touchdowns. This constrains the direction of the contact force to the landing gear plane. As the acceleration vector was approximately orthogonal to the lander body-axis and in line with the landing gear during contact, Philae must have struck the comet with the side structure of the foot assembly and not with the soles of the feet. This is the major difference between this collision and the other touchdowns and the reason that Auster et al. (2015) introduced the different naming scheme. The approximation of the direction of acceleration was made possible by the fact that the ROMAP magnetometer was mounted at the end of a spring-supported boom. An assembly like this causes the magnetometer to swing from its nominal position as soon as the lander is accelerated perpendicular to the hinge axis. As the natural frequency of the boom assembly is significantly above the excitation caused by the surface contacts, the boom movement correlates to the impact acceleration. As the background magnetic field is known (Heinisch et al. 2017b), this movement can be tracked using the magnetic field measurements, turning ROMAP into a simple accelerometer.

The attitude just after second touchdown implies that the balcony was pointing toward the comet and at least one foot, the tip of the ROMAP boom, and possibly part of the balcony were in contact with the surface. This was confirmed by a deflection of the ROMAP boom away from the landing gear and the fact that the ROMAP plasma sensor was covered by dust (Heinisch et al. 2017a). The position of the entrances of the plasma monitor means that the sensor head must have penetrated the surface at least ~ 5 cm for the senor entrances to become covered by dust. This contact caused an upward deflection of the boom, as it is only supported by springs. As the plasma sensors housekeeping data remained nominal for the entire time, the sensor could not have sustained damage during contact. Especially the delicate charged retarding grids at the sensor entrance were not damaged. Considering an impact velocity of ≈ 0.22 m s⁻¹, this is already a qualitative indication that the surface material has to be soft.

Fig. 1 Overview of the four contact sites on the upper lobe of comet 67P. All contacts took place during Philae’s descent on November 12, 2014. The zoomed-in inlays show the differences in morphology between the four contacts. The possible areas of contact are circled in white. The flight direction is illustrated by the red arcs.

Fig. 2 Panela: possible area for the collision event overlaid on the digital terrain model (SHAP7) and the estimated trajectory (constrained by red lines, direction indicated by white arrow) from first touchdown, including the position of Philae in flight as seen in OSIRIS camera images. Panel b: possible area for second touchdown overlaid on an OSIRIS image. The incoming flight direction is indicated by a red arrow. Philae is visible to the right at its final locations (circled).

3 Compressive stress

The rebound of Philae offers the possibility of determining the compressive strength σ, which can be described as the threshold between plastic and elastic deformation. If the applied load is below the compressive strength limit of the material for a given area, the contact will be virtually elastic, with negligible loss of kinetic energy leading to rebounding. This would cause a strong but transient acceleration. Otherwise, the contact becomes inelastic as a result of plastic deformation or reordering of the material (e.g., Schräpler et al. 2015). As an inelastic contact is linked to a penetration of the material, the acceleration is weaker but spread out over a longer interval. Thus, an upper estimate can be derived by examining the load on the cometary material applied during contact and change in kinetic energy linked to it. Granular material behaves differently under load than regular solids (Biele et al. 2009), for example, exhibiting slip-stick interaction (Cole & Peters 2007). Hence, standard elastoplastic models (e.g., Li et al. 2001; Yigit et al. 2011)cannot be used. As there is no general description of the behavior of granular material under load (Cole & Peters 2007), a constant deceleration and hence a constant breaking force during contact is assumed. Even though this is not an exact description of real granular cometary material, it generally leads to an overestimation of the compressive strength, which is of no concern as we are only interested in an upper bound. Furthermore, additional energy losses inside the lander structure are neglected (similar on Biele et al. 2015). With these assumptions, an upper bound for the compressive strength of the material can be derived based on the total loss of energy using $$\text{σ}=\frac{F}{A}=\frac{a_d\cdot m_{Pl}}{A}=\frac{\left(\sqrt{\frac{2E_0}{m_{Pl}}}-\sqrt{\frac{2E_1}{m_{Pl}}}\right)m_{Pl}}{A{\text{τ}}_c}.$$(5)

Here F is the breaking force linked to the deceleration a_d. E₀ and E₁ denote the total energy before and after impact, respectively, τ_c is the contact duration, and A is the area of contact between the landing gear and the surface. The contact duration τ_c was determined from the movement of the ROMAP magnetometer boom caused by the acceleration during contact. As the boom is only supported by two springs, it can be deflected if an external acceleration is applied. This movement creates a characteristic signature in the observed magnetic field, which can be used to derive the duration of the acceleration (Heinisch et al. 2017a). The area of contact A was determined based on the known geometry of the landing gear, particularly the feet, combined with the angle of attack of the lander relative to the surface (derived from lander attitude) just before contact. The parameter m_Pl = 97.63 kg (Biele et al. 2015) is the mass of Philae. As we are primarily interested in an upper bound for σ, it is assumed that the change in total kinetic (E = E_rot + E_kin) energy is only linked to the surface interaction. Hence the deceleration is determined based on the total loss of energy during contact and not just on the translational velocity.

For the collision event, the total energy E₀ before contact was (7.02 ± 0.37) J, which decreased to E₁ = (2.83 ± 0.22) J afterward. Based on the reconstructed attitude (Heinisch et al. 2017a), Philae struck the surface with only one foot. In contrast to first touchdown, the foot hit the comet not with the soles, but with the side of the foot structure, which leads to a contact area of A = (225 ± 40) cm². This area was estimated based on the CAD model of the foot assembly (depicted in Fig. 4) and a plane angled relative to the foot to represent the attitude of the lander. Even though the exact angle between the local surface and the foot is unknown, angular deviations only lead to small changes in the area of contact because of the specifics of the foot design. The contact time τ_c = (1.5 ± 1.0) s was determined from the acceleration derived from ROMAP observations, as explained above. This leads to an upper limit for the compressive strength of σ = (399 ± 393) Pa.

The same approach was used for second touchdown. The kinetic energy E₀ before contact was (2.97 ± 0.15) J, which is slightly higher than directly after the collision because of the difference in height between thearea of the collision and second touchdown. Afterward, the energy dropped to E₁ (0.10 ± 0.02) J. The different contact geometry (at least one foot and magnetometer boom) leads to a higher contact area of (350 ± 75) cm² and a longer contact duration of τ_c = (4 ± 1) s compared to the collision. This results in a compressive strength σ = (147 ± 77) Pa. These results are summarized in Table 1.

Images of Philae at the final Abydos landing site (Fig. 5) revealed several scratch features (S_A, S_B, and S_C). Scratches S_A and S_B are located to the left of the lander and have a length of ~1.8 m, a width of ~11 cm, and are ~42 cm apart. The length was derived from the number of pixels while accounting for projection effects and using the known dimensions of the lander for validation. It is possible that the lengths of S_A and S_B are slightly underestimated, as the boulder is only visible up to the edge of what is depicted in Fig. 5. The imaging geometry means that some parts could be hidden behind smaller boulders in the foreground. Based on the size of these structures, a possible length error of 0.2 m was determined. The scratches are aligned with the incoming flight path of the lander coming from second touchdown (see Fig. 2 for reference). As all high-resolution images of this area were taken almost at the same time with similar geometry and lighting, it is impossible to reconstruct the exact 3D topography. Therefore the directionality of the scratch marks cannot be used to extrapolate the exact position of second touchdown. Based on the available local terrain model, it is only possible to estimate the approximate flight direction after second touchdown. No other parallel structures at scales comparable to the scratch marks are visible in the images of the landing area (Bibring et al. 2015; Poulet et al. 2016), which suggests that the scratches were made by Philae. From the attitude information derived from the magnetometer measurements (Heinisch et al. 2017a), we can infer that two of the landing legs were pointing toward the surface shortly before third touchdown. This makes it likely that the lander skidded across the surface for the final ~ 1.8 m, causing scratch marks S_A and S_B.

A third feature (scratch S_C) with a length of approximately 30 cm is visible above Philae next to one of the feet. We assume that scratch S_C was created while Philae came to a stop after it cleared the cliff (bottom left-hand corner in Fig. 5) and tilted while hitting the wall it now rests on. As the direction of scratch S_C is parallel to the lander z-axis, the deceleration caused an observable movement of the ROMAP boom, which was visible in the magnetic field observations from 17:31:16 ± 1 s until 17:31:26 ± 1 s UTC.

The width of the scratches and the size of Philae’s leg structure indicates that only the soles of the feet have been in contact with the ground. This constrains the maximum penetration depth to approximately ~ 5 cm (see Fig. 4 for the dimensions of the foot), which is half the diameter of the soles. Otherwise, the significantly larger upper foot structure (see Fig. 5 for the total scale of the landing gear) would have been in contact with the surface, causing substantial drag and considerably larger scratch marks in the range of ~ 25 cm. A width-to-depth ratio of ~2 is also consistent with other depressions found on 67P (El-Maarry et al. 2015) and can be considered realistic for the surface material. It is impossible to derive the actual depth of the scratches from the images because of the diffuse lighting conditions. From the characteristics of the scratches, it is possible to estimate the contact pressure σ analogous to Eq. (5) based on the mechanical work done during creation of the scratches: $$\text{σ}=\frac{F}{A}=\frac{\Delta E}{A{l}_s}.$$(6)

Here A is again the area of contact, F is the breaking force, ΔE denotes the mechanical work done while creating the scratches, and l_s is the length of the scratches. While the material in front of the foot is compressed, additional energy is required to overcome the shear and tensile strength (see, e.g., Akono & Ulm 2011 for scratch modeling) at the boundary of the cometary material and the outer rim of the foot. As the aim of this study is to derive an upper limit for the compressive strength, it was assumed that the entire energy lost during the creation of the scratches was due to compression, hence leading to the highest possible compression force and preventing underestimation of the compressive strength.

To determine ΔE for scratches S_A and S_B, the remaining velocity afterward was estimated as half of the mean initial velocity (again assuming constant deceleration) as a conservativefirst-order approximation, based on the impact drag force analysis performed by Katsuragi & Blum (2017). The overall ratio of the individual values means that minor changes in the velocities would not cause major changes in the resulting compressive strength, $$v_{\text{out}}=\frac{l_s}{2{\text{ τ}}_{AB}},$$(7)

with τ_AB = (36 ± 5) s and l_s = (1.8 ± 0.2) m, resulting in v_out = (0.03 ± 0.006) m s⁻¹. This yields a remaining kinetic energy of the lander of (0.04 ± 0.01) J, which combined with an incoming kinetic energy of (0.10 ± 0.02) J (see above) translates into ΔE = (0.06 ± 0.03) J. With a contact area of A_AB = (41 ± 10) cm² estimated based on the foot geometry and penetration depth, this yields an upper limit for the compressive strength of σ = (8 ± 7) Pa. Because of the change in attitude and single foot contact while scratch S_C was created, the contact area was smaller, with A_C = (20 ± 10) cm². With the remaining energy of (0.04 ± 0.01) J lost during contact and a length of the scratch of l_s = (0.3 ± 0.1) m, the resulting limit for the compressive strength is σ = (73 ± 70) Pa. These results are summarized in Table 2.

Even though the force necessary to plow through the material depends on several material properties, the resulting pressure is dominated by the compressive properties of the material. The compressive strength is expected to be the highest by afactor of roughly 10, especially due to the low gravitational compression of the material (e.g., Groussin et al. 2015).

Fig. 3 Histogram of the incoming and outgoing total velocities for the collision (panel a) and second touchdown (panel b) calculated based on each possible trajectory using a spatial resolution of 1 m. Panels c and d: difference in energy before and after collision and second touchdown, respectively, which is equal to the work done during theses contacts.

Fig. 4 Illustration of one of Philae’s feet based on a rendering of the 3D lander CAD model. The approximate dimensions of the soles and the upper feet structure are given in red.

Fig. 5 OSIRIS camera image of the Philae lander at the final Abydos landing site. The right panel is annotated to show the scratch marks “A”, “B” and “C”, the approximated incoming flight direction is depicted by a white arrow. The image resolution is approximately 7 mm px−1.

4 Results and conclusions

The in situ results we presented for the compressive strength span a range between 7 Pa and 399 Pa, which agrees well with the compressive strength of 30 Pa to 150 Pa estimated from remote observations by Groussin et al. (2015). Deviations between the upper limitsfor the compressive strength values derived for the individual contacts can most likely be attributed to variations in local material. While the comet is homogeneous on larger sales, local differences in insolation, the dust-to-ice ratio (due to different outgassing rates) or volume filling factor are expected and lead to local deviations in the compressive strength. Possible uncertainties in the area of contact between the lander and the surface and errors in the estimated trajectories can also contribute to the differences, but these deviations should be covered by the conservative error range. Owing to the simplified mechanical model and the limited range of available measurements, the resulting compressive strength is only a local upper limit, without taking higher-order mechanical characteristics of granular material into account.

Based on the mechanical interaction during the surface contacts, a compressive pressure below 100 Pa is already enough to cause significant compaction of the surface material. This compaction led to an increase in compressive strength, which caused Philae to rebound at ~160 Pa. This kind ofmaterial makes comets like 67P an ideal target for lander missions as the surface is strong enough to support the weight of a probe without overpenetration, but is also soft enough to dissipate the kinetic energy during touchdown.

The range for the compressive strength derived above is slightly below the value of 1 kPa estimated by Biele et al. (2015) from first touchdown, which can largely be attributed to the difference in dust coverage. In contrast to the first landing site Agilkia, neither the cliff near the Hatmehit crater nor Abydos have significant surface dust layers (see Fig. 1; Schröder et al. 2017; Keller et al. 2015; Thomas et al. 2015). Hence, the mechanics of the subsequent contacts were not governed by surface dust layers and therefore provide insight into the properties of the more consolidated surface regions. This is also the reason why the scratch marks at Abydos are still visible two years after Philae landed, and it shows that regions that might look consolidated can in reality be extremely soft. Based on the hammering of the MUPUS instrument, Spohn et al. (2015) postulated a compressive strength of 2 MPa for the Abydos site, which is significantly above the independently derived results of Biele et al. (2015) and Groussin et al. (2015) or the range presented above. A possible explanation for this discrepancy is an overestimation of the compressive strength caused by uncertainties in the MUPUS deployment, especially as the placement of the hammer could not be confirmed independently by camera images.

The results of the scratch mark analysis are also inconsistent with a solidified surface ice layer with a thickness between 10 and 50 cm as proposed by Knapmeyer et al. (2018). It is also questionable if such an ice layer is consistent with the fact that the ROMAP plasma sensor was covered by dust during second touchdown (Heinisch et al. 2017a). This requires a relatively soft upper layer, as the ROMAP instrument was not damaged, but capable of breaking up the upper surface. As the results of Knapmeyer et al. (2018) rely on the correct deployment of the MUPUS hammer, uncertainties in the MUPUS operation might have caused an overprediction of the thickness of a possible surface ice layer.

Several laboratory experiments have been performed to determine the mechanical properties of possible building blocks of comets. In particular, Güttler et al. (2009), Schräpler et al. (2015), Lorek et al. (2016), and Katsuragi & Blum (2017) investigated how the compressibility of ice layers and different dust structures created by random ballistic deposition as well as pebble sub-structures (with scales of centimeters, millimeters, and 0.1 mm) depend on the volume filling factor. The respective results are summarized in Fig. 6. Remote observations lead to multiple independent estimates for the volume filling factor. Kofman et al. (2015) derived values of 0.15–0.25, Pätzold et al. (2016) estimated 0.25–0.30, and a study by Fulle et al. (2016) resulted in 0.21–0.37. The corresponding range is depicted in Fig. 6 by the horizontal bar shaded in blue. As thevolume filling factor is known, the laboratory measurements can be used to determine the most likely material of the surface. For the observed volume filling factor, dust and ice layers have a compressive strength above 10³ Pa. A surface made up of such dust or ice layers is therefore inconsistent with the results of Groussin et al. (2015) or the compressive strength derived as part of this work. In contrast, a surface consisting of layers of dust and ice aggregates could explain all observations, except for the MPa range derived by Spohn et al. (2015). This is strong evidence for the presence of aggregate layers on the surface of 67P, that is, a surface composed of sub-decimeter-sized aggregates, as inferred from the gravitational collapse scenario. This conclusion can be drawn without precise knowledge of the actual value of the compressive strength, as the upper limit derived as part of this work and by Groussin et al. (2015), even considering all possible errors, is well below 10³ Pa for different parts of 67P. The presence of such aggregate layers can only be explained by a formation dominated by relatively low-velocity collisions. Hierarchical agglomeration would cause a significantly higher compressive strength as a result of impact-compaction of the pebbles and possible break-up of aggregates (Blum 2018). The derived mechanical properties combined with the available laboratory models therefore suggest that the gentle collapse of aggregate ensembles played a major role in the formation ofcomets like 67P.

Fig. 6 Dependence between volume filling factor and compressive pressure for dust aggregate, dust, and water-ice layers measured inthe laboratory (SC15: Schräpler et al. 2015; GU09: Güttler et al. 2009; LO16: Lorek et al. 2016). The experimental results for dust and ice layers have been extrapolated for lower volume filling factors (dashed lines). The derived compression range for the surface contacts of Philae are depicted by the vertical blue area. Measurements from Kofman et al. (2015; K), Pätzold et al. (2016; F), and Fulle et al. (2016; P) have been used to determine the possible range for the volume filling factor.

Acknowledgements

Rosetta is an ESA mission with contributions from its member states and NASA. We are indebted to the whole Rosetta Mission Team, LCC, SONC, SGS and RMOC for their outstanding efforts making this mission possible. The contributions of the ROMAP and RPC-MAG teams were financially supported by the German Ministerium für Wirtschaft und Energie and the Deutsches Zentrum für Luft- und Raumfahrt (DLR) under contract 50QP1401. Jürgen Blum acknowledges funding through DFG project BL298/24-1 and DLR project 50WM1536. The SPICE kernels used for this work are available as ROS-E/M/A/C-SPICE-6-V1.0 through PSA and PDS. We thank the anonymous referee for a very constructive review.

References

All Tables

All Figures

	Fig. 1 Overview of the four contact sites on the upper lobe of comet 67P. All contacts took place during Philae’s descent on November 12, 2014. The zoomed-in inlays show the differences in morphology between the four contacts. The possible areas of contact are circled in white. The flight direction is illustrated by the red arcs.
In the text

Fig. 2 Panela: possible area for the collision event overlaid on the digital terrain model (SHAP7) and the estimated trajectory (constrained by red lines, direction indicated by white arrow) from first touchdown, including the position of Philae in flight as seen in OSIRIS camera images. Panel b: possible area for second touchdown overlaid on an OSIRIS image. The incoming flight direction is indicated by a red arrow. Philae is visible to the right at its final locations (circled).

In the text

	Fig. 3 Histogram of the incoming and outgoing total velocities for the collision (panel a) and second touchdown (panel b) calculated based on each possible trajectory using a spatial resolution of 1 m. Panels c and d: difference in energy before and after collision and second touchdown, respectively, which is equal to the work done during theses contacts.
In the text

	Fig. 4 Illustration of one of Philae’s feet based on a rendering of the 3D lander CAD model. The approximate dimensions of the soles and the upper feet structure are given in red.
In the text

	Fig. 5 OSIRIS camera image of the Philae lander at the final Abydos landing site. The right panel is annotated to show the scratch marks “A”, “B” and “C”, the approximated incoming flight direction is depicted by a white arrow. The image resolution is approximately 7 mm px−1.
In the text

Fig. 6 Dependence between volume filling factor and compressive pressure for dust aggregate, dust, and water-ice layers measured inthe laboratory (SC15: Schräpler et al. 2015; GU09: Güttler et al. 2009; LO16: Lorek et al. 2016). The experimental results for dust and ice layers have been extrapolated for lower volume filling factors (dashed lines). The derived compression range for the surface contacts of Philae are depicted by the vertical blue area. Measurements from Kofman et al. (2015; K), Pätzold et al. (2016; F), and Fulle et al. (2016; P) have been used to determine the possible range for the volume filling factor.

In the text