Issue 
A&A
Volume 658, February 2022



Article Number  A118  
Number of page(s)  16  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/202142098  
Published online  09 February 2022 
The influence of gravity on granular impacts
II. A gravityscaled collision model for slow interactions
^{1}
Institut Supérieur de l’Aéronautique et de l’Espace,
Ave Edouard Belin,
31400
Toulouse cedex 4,
France
email: cecily.sunday@isaesupaero.fr
^{2}
Université Côte d’Azur, Observatoire de la Côte d’Azur, CNRS, Laboratoire Lagrange, Bd de l’Observatoire,
CS 34229,
06304
Nice cedex 4,
France
^{3}
Centre National d’Études Spatiales,
Ave Edouard Belin,
31400
Toulouse cedex 4,
France
Received:
27
August
2021
Accepted:
24
November
2021
Context. Slow interactions on small body surfaces occur both naturally and through human intervention. The resettling of grains and boulders following a cratering event, as well as observations made during small body missions, can provide clues regarding the material properties and the physical evolution of a surface. In order to analyze such events, it is necessary to understand how gravity influences granular behavior.
Aims. In this work, we study slow impacts into granular materials for different collision velocities and gravity levels. Our objectives are to develop a model that describes the penetration depth in terms of the dimensionless Froude number and to use this model to understand the relationship between collision behavior, collision velocity, and gravity.
Methods. We used the softsphere discrete element method to simulate impacts into glass beads under gravitational accelerations ranging from 9.81 m s^{−2} to 0.001 m s^{−2}. We quantified collision behavior using the peak acceleration, the penetration depth, and the collision duration of the projectile, and we compared the collision behavior for impacts within a Froude number range of 0–10.
Results. The measured penetration depth and collision duration for lowvelocity collisions are comparable when the impact parameters are scaled by the Froude number, and the presented model predicts the collision behavior well within the tested Froude number range. If the impact Froude number is low (0 < Fr < 1.5), the collision occurs in a regime that is dominated by a depthdependent quasistatic friction force. If the impact Froude number is high enough (1.5 < Fr < 10), the collision enters a second regime that is dominated by inertial drag.
Conclusions. The presented collision model can be used to constrain the properties of a granular surface material using the penetration depth measurement from a single impact event. If the projectile size, the collision velocity, the gravity level, and the final penetration depth are known and if the material density is estimated, then the internal friction angle of the material can be deduced.
Key words: minor planets, asteroids: general / planets and satellites: surfaces / methods: numerical / methods: miscellaneous
© C. Sunday et al. 2022
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1 Introduction
Two recent small body missions have involved direct interactions with asteroid surfaces: the JAXA Hayabusa2 mission (Watanabe et al. 2017) and the NASA OSIRISREx mission (Lauretta et al. 2017). The Hayabusa2 spacecraft arrived at asteroid (162173) Ryugu in June 2018 and deployed several modules (a lander and two minihoppers) to the asteroid’s surface (Van Wal et al. 2018; Ho et al. 2021). In addition, the spacecraft successfully collected a surface sample (Morota et al. 2020) and performed an impact experiment to observe the creation of an artificial crater (Arakawa et al. 2020). From December 2018 to June 2021, the OSIRISREx spacecraft surveyed and mapped the surface of asteroid (101955) Bennu. In October 2020, it performed a successful touchandgo maneuver to retrieve a sample from the asteroid’s surface (Bierhaus et al. 2018). Ryugu and Bennu, similar to most asteroids, are covered with boulders and a layer of loose granular material, referred to as regolith (Murdoch et al. 2015; Hestroffer et al. 2019). Granular materials display notoriously complex behaviors on Earth, and the fundamental contact laws used to describe particletoparticle behavior remain the same for different gravity levels. However, the electrostatic and cohesive forces that are often ignored on Earth can become nonnegligible in lowgravity (Scheeres et al. 2010; Hestroffer et al. 2019), and the bulk granular system can experience regime changes under different external triggers (Brucks et al. 2007; Murdoch et al. 2017). The nonintuitive nature of granular flow in lowgravity, paired with frequent unknowns regarding target body surface materials, makes it difficult to plan and analyze spacecraft operations similar to those listed above. It also makes it very challenging to deduce the origins of the various surface features observed on asteroids.
The ultimate objective of most small body missions is to understand the formation and evolution of these celestial objects. The Hayabusa2 and OSIRISREx spacecrafts are not the first to have visited small bodies, nor will they be the last to conduct in situ surface experiments. JAXA’s upcoming Martian Moons eXploration mission (MMX) will deploy a rover and will land on the surface of Phobos (Michel et al. 2022; Usui et al. 2020). The spacecraft will be equipped with a system to collect and return a sample to Earth in 2029. The NASA DART mission, planned for launch in late November 2021, will perform an impact experiment on the small moon called Dimorphos, which is part of the binary asteroid (65803) Didymos (Cheng et al. 2018). Following DART, ESA will launch the Hera mission to observe the outcomes of the impact and to track the dynamical and physical changes of the asteroid system (Michel et al. 2018). The Hera mission includes two CubeSats that may attempt to land on the surface of Dimorphos and return information on the surface response using accelerometer measurements (Goldberg et al. 2019; Ferrari et al. 2021). A basic understanding of regolith dynamics is required to design and operate the next generation of impactors, landers, and sampling systems. More critically, scientists need reasonably accurate models to help analyze observations from toolregolith interactions. With the right data, the act of impacting or landing on a granular surface can provide valuable insight regarding a surface’s friction, cohesion, and density properties, not to mention its response to external forces. This information is crucial for developing a robust timeline of the objects’ geological history, which is linked to that of the Solar System and its surface chronology.
In this study, we introduce a framework that can be used to analyze the outcome of slow impacts and interactions with granular materials, independent of gravity level. This model can therefore be applied to regolithcovered celestial bodies of any size. Several existing studies have investigated the link between lowvelocity impact dynamics and landing on small body surfaces. One popular area of research focuses on rebound behavior and bouncing dynamics in reducedgravity environments (Colwell & Taylor 1999; Colwell 2003; Çelik et al. 2021), while another investigates the mass transfer and ejecta patterns associated with slow collisions (Schwartz et al. 2014; Brisset et al. 2018, 2020). Several works have attempted to simulate and analyze specific lander and spacecraft interactions, such as those associated with the Hayabusa2 and the OSIRISREx missions (Maurel et al. 2018; Thuillet et al. 2018, 2020; Ballouz et al. 2021). Here, we investigate the general behavior of a spherical projectile impacting a bed of spherical grains at low speeds. In accordance with our previous experiments, we assume that the projectile comes to a rest without rebounding (Murdoch et al. 2017, 2021; Sunday et al. 2021).
The objective of this study is to develop a theoretical collision model that links a projectile’s final penetration depth to its impact velocity, the gravity level of the system, and the properties of the surface material. In Sect. 2.1, we present existing models and discuss why these solutions cannot be used to predict impact behavior for different gravitational accelerations. In Sects. 2.2 and 2.3, we introduce a revised collision model, and we reformulate the model in terms of the dimensionless Froude number, which takes into account both impact velocity and gravity. Then, we extend the work of Sunday et al. (2021) by performing impact simulations using the softsphere discrete element method (SSDEM) and the opensource code Chrono (Tasora et al. 2016; Sunday et al. 2020). In Sect. 3, we use the simulation results to show how collision behavior scales with impact velocity and gravity, and we assess the accuracy of the revised model. Finally, in Sect. 4, we demonstrate how collision behavior relates to the physical and mechanical properties of a surface material, and we discuss the implications for upcoming small body missions.
2 Method
2.1 Existing collision models
With the help of diverse experimental and numerical techniques, lowvelocity collisions into granular materials have been studied extensively for different surface materials and projectile shapes (Omidvar et al. 2014; Katsuragi 2016). While some of this research focuses on the refinement of empirical scaling relationships (Uehara et al. 2003; de Bruyn & Walsh 2004), a significant portion of the existing work attempts to characterize projectile behavior in terms of Newton’s second law of motion (see Omidvar et al. 2014 for more information): (1)
In Eq. (1), m is the mass of the projectile, z is the depth of the projectile, g is gravity, and f(z), b(z), and h(z) are the coefficients associated with the drag forces that act on the projectile as it moves through a granular medium. f(z) is a quasistatic resistance force, b(z) ż a viscous resistance force, and h(z) ż^{2} is an inertial resistance force. f(z) can either be constant or depthdependent depending on the model, but b(z) and h(z) are often considered to be constant and depthindependent (see Table 1).
By neglecting or modifying the form of the drag terms in Eq. (1), the general force law can be simplified to better represent impact behavior within specific ranges of collision velocities. For example, the viscous and inertial drag terms are generally neglected from the model when , where v_{c} is the projectile’s impact velocity and d is the average grain diameter (Albert et al. 1999). Above this threshold, however, the behavior of the projectile becomes velocitydependent. When the inertial term is neglected, Eq. (1) resembles the Bingham force model, where . If the viscous term is neglected instead of the inertial term, Eq. (1) takes on the form of the Poncelet model, where (Omidvar et al. 2014; Katsuragi 2016). Table 1 provides a summary of the most frequently cited models in reference to lowvelocity granular impacts (i.e., impact where v_{c} < 5−6 m s^{−1}).
The models proposed by de Bruyn & Walsh (2004) and Clark & Behringer (2013) seem appropriate for impacts where v_{c} ≈ 1–6 m s^{−1}, but they exhibit nonphysical behavior for dense materials and for lower collision velocities. More specifically, these models fail to capture the nonzero initial sinkage of a projectile when v_{c} = 0. Replacing the constant friction term f_{o} with a depthdependent f(z) term improves a model’s ability to predict collision behavior as v_{c} → 0. Katsuragi & Durian (2013), Lohse et al. (2004), and Kang et al. (2018) introduce models where f(z) varies linearly with depth, while Tsimring & Volfson (2005) assume that f(z) varies with the area of the projectile in contact with the grains. Alternatively, Brzinski III et al. (2013) assert that f(z) varies by the volume of the material displaced by the projectile. PachecoVázquez et al. (2011) introduce yet another form of f(z) that takes into accountthe effect of container’s size on the pressure in the system (i.e., the Janssen effect, described in Katsuragi 2016). The models listed in Table 1 succeed to various degrees when used to predict the penetration depth or the drag force for collisions where v_{c} < 6 m s^{−1}. With the exception of Brzinski III et al. (2013), however, these references do not address how their findings apply to different gravity levels, despite the fact that g appears as a variable in the expressions.
Several other groups have performed experiments using either fluidized granular beds or droptower setups in order to better understand the role of gravity within the context of select collision models. Goldman & Umbanhowar (2008) use a simple droptower to investigate how collision behavior scales with projectile size, surface density, and gravity while comparing their results against several of the models in Table 1. Costantino et al. (2011) use a fluidized granular bed to evaluate gravity scaling for the f(z) term in the Tsimring & Volfson (2005) model. Altshuler et al. (2014), in comparison, reference the f(z) term in the PachecoVázquez et al. (2011) model while studying settling behavior using a droptower setup. Murdoch et al. (2021) also use a droptower setup (described in Sunday et al. 2016), but they consider impact, not settling behavior, and they compare their experimental measurements to the model referenced in Clark & Behringer (2013).
The above studies highlight the limitations of applying the existing collision models to impacts in lowgravity environments. Even though the scaling relationships for gravity are still unclear, the shared and most notable issue with the current models is that the definition for a “lowvelocity” collision does not account for gravity. Interestingly, granular materials have been shown to fluidize more readily in lowgravity (Murdoch et al. 2017), so the model that best describes a 1 m s^{−1} collision on Earth is not likely to be valid for a 1 m s^{−1} collision on a small body.
In the remainder of this section, we describe a collision model that is applicable for all gravity levels. We begin by introducing a revised model that describes impacts up to 10 m s^{−1} under terrestrial gravity. Then, we reformulate the model to highlight the role that gravity plays in the collision process, and we define two distinct impact regimes that are characterized by both collision velocity and gravity. At the end of the section, we provide an overview of the simulation setup that was used to evaluate the revised model.
Frequently cited collision models that derive from the general force model, .
2.2 Revised collision model
There is strong evidence to support the hypothesis that the quasistatic friction force in Eq. (1) varies with depth when a projectile first impacts a surface, but that the friction force eventually saturates and becomes constant if the initial impact velocity of the projectile is high enough (Brzinski III et al. 2013; Tsimring & Volfson 2005). Here, we assume that f(z) becomes constant at a depth of z_{1} and that (2)
where f_{o} and h are constants. We then define two collision regimes based on z_{1} and the final penetration depth of the projectile, z_{stop}. In collision regime I, z_{stop} ≤ z_{1} and f(z) = f_{o}z∕z_{1} for the entire duration of the collision. In collision regime II, z_{stop} > z_{1}, and f(z) is depthdependent when z ≤ z_{1} but constant when z > z_{1}. Figure 1 illustrates the two phases of the collision during regime II.
Expressions for the penetration depth in collision regime I, z_{stop,I}, and the penetration depth in collision regime II, z_{stop,II}, can be obtained by recasting Eq. (1) in terms of kinetic energy K and solving the linear ordinary differential equation with the correct boundary conditions for each regime. Following the method detailed in Appendix A, (3)
Equation (3) is valid when z_{stop,I} ≤ z_{1} (collision regime I) and Eq. (4) is valid when z_{stop,II} > z_{1} (collision regime II). If v_{c} = 0, then z_{stop,I} = z_{o}, or the initial sinkage of the projectile. If z_{o} is a known value, then Eq. (3) can be rearranged to calculate z_{1} as follows: (5)
As will be demonstrated in Sect. 3.1, the above model succeeds in reproducing two important observation from previous works: (1) the projectile has a nonzero penetration depth when v_{c} = 0, and (2) the quasistatic friction force is constant for higher collision velocities. The key difference between Eq. (4) and existing solutions (i.e., those found in Clark & Behringer 2013), is that here, we consider that the kinetic energy of the projectile is when f(z) becomes constant, not . The velocity term v_{1} corresponds to the velocity of the projectile at depth z_{1}.
Fig. 1 Force models applied to collision regime II. The quasistatic drag force acting on the projectile f(z) varies linearly when 0 < z < z_{1} but is constant when z_{1} < z < z_{stop,II}. m is the mass of the projectile, g is gravity, z is the depth of the projectile, z_{1} is the depth where f(z) becomes constant, z_{stop,II} is the final penetration depth in collision regime II, and f_{o}, and h are constants. 
2.3 Collision model based on Froude scaling
When studying lowvelocity impacts, it is common to adopt the hypothesis that the f(z) term in Eq. (1) originates from hydrostatic pressure and Coulomb friction such that f(z) ∝ ρ_{g}gAz, where A is the cross sectional area of the projectile and ρ_{g} is the bulk density of the surface material (Tsimring & Volfson 2005; Katsuragi 2016). On the other hand, the inertial drag force, h(z)v^{2}, can be expressed as the momentum transfer between the projectile and the surface material, where h(z) ∝ ρ_{g}A (Katsuragi 2016). At first glance, f(z) has a clear dependence on gravity, while h(z)v^{2} does not.
With the above relationships, Katsuragi (2016) shows that the general force model can be broken down into five key dimensionless Πgroups using a, v, z, g, D, ρ_{p}, ρ_{g}, and μ. The ρ_{p} term is the density of the projectile, and μ is a coefficient related to the internal friction angle of the granular material. One of the dimensionless groups, v^{2} ∕gD, corresponds to the Froude number Fr. The Froude number is the ratio of the inertial to gravitational forces in a system and is a powerful dimensionless parameter for demonstrating how certain processes scale with gravity. For instance, Brucks et al. (2007) use Fr to study granular flow in a rotating drum for varied gravity levels, Housen & Holsapple (2003) use a term resembling Fr to develop crater size scaling relationships, Hilton & Tordesillas (2013) and Faug (2015) use Fr to study drag through granular materials, and Wright et al. (2020) reference Fr in the context of ricochets on asteroid surfaces.
The exact expression of Fr varies by application and is different in all of the above works. As in Wright et al. (2020), we define the impact Froude number as (6)
where D is the projectile diameter. The expressions in Sect. 2.2 can be reformulated in terms of Fr and by substituting Eq. (6) into Eqs. (3) and (4). Since f(z) is a function of gravity, we assume that , and we rewrite the penetration depth equations for regimes I and II as (7)
Here, , , , and h′ = hD. , , and are dimensionless parameters while and h′ have the units of mass. In the reformulated expressions, the explicit gravity term disappears.
The transition from collision regime I to collision regime II, denoted as Fr_{t}, occurs when . When Fr > Fr_{t}, the penetration depth exceeds z_{1} and the collision takes place in regime II. Making the appropriate substitutions and rearranging Eq. (7) gives the following expression for Fr_{t}: (9)
2.4 Numerical simulations and data analysis
To assess the validity of the proposed collision model, numerical simulations of a spherical projectile impacting glass beads at different collision velocities and gravity levels were conducted using the softsphere discrete element method (SSDEM) and the Multicore module in the opensource code Chrono (Tasora et al. 2016; Sunday et al. 2020). For all tests, a 1 kg, 10 cm diameter projectile was dropped into a cylindrical container filled with 1 ± 0.05 cm diameterbeads. The container had a nominal diameter of 45 cm and was filled to a height of 22 cm. The beads were initially mixed and settled under terrestrial gravity conditions. Then, for the lowgravity tests, the beads were permitted to relax and resettle under the desired gravity level. The same surface sample was used for all of the simulations, but the initial position of the projectile was randomized within a radius of two particle diameters at the beginning of each test to add variability to the impact configuration. The bead properties and the simulation parameters were selected by calibrating the simulation results against terrestrialgravity impact experiments. A detailed description of the code, the simulation setup, the boundary conditions, and the simulation parameters can be found in Sect. 2.2 of Sunday et al. (2021).
In order to quantify the projectile’s collision behavior for different gravity levels, four collision parameters were extracted from each simulated impact: (1) the projectile’s collision velocity v_{c}, (2) the projectile’s peak acceleration a_{peak}, (3) the projectile’s final penetration depth z_{stop}, and (4) the total duration of the collision t_{stop}. Figure 2 shows the raw data curves that were obtained for a typical impact simulation. In this example, the projectile hits the granular bed with a collision velocity of 1 m s^{−1}. The projectile’s acceleration a, velocity v, and penetration depth z as a function of time t are shown in Figs. 2a–c, respectively. The results for a test case where g = 9.81 m s^{−2} are shown in blue, and the results for a case where g = 0.1 m s^{−2} are shown in red. The definition of the peak acceleration, the collision duration, and the penetration depth parameters are indicated bythe dashed lines and the text annotations on the plot. The exact method and criteria that was used to identify a_{peak}, z_{stop}, and t_{stop} can be found in Sect. 2.3 of Sunday et al. (2021). Additional information regarding the filtering and treatment of the acceleration data for this particular study can be found in Appendix B. Next, in Sect. 3, we compare the collision parameters for over 150 impact simulations with collision velocities ranging from 0 to 10 m s^{−1}.
Fig. 2 Acceleration, velocity, and position profiles for example impact simulations where g = 9.81 m s^{−2} (blue) and g = 0.1 m s^{−2} (red). For both test cases, the collision velocity of the projectile is 1 m s^{−1}. The projectile’s approximate collision velocity v_{c}, peak acceleration a_{peak}, penetration depth z_{stop}, and collision duration t_{stop} are indicated by the circles, dashed lines, and text annotations on the plot. 
3 Results
3.1 Collision behavior by gravity level
In this section, we compare impact behavior for three different gravity levels: g = 9.81 and 10 m s^{−2} (terrestrial gravity), g = 0.1 m s^{−2}, and g = 0.001 m s^{−2}. The 0.1 m s^{−2} gravity case is lower than but comparable to the gravity found on asteroid (4) Vesta, and the 0.001 m s^{−2} gravity case is on the same order of magnitude as the gravities found on bodies such as asteroid (433) Eros and the moons of Mars, Phobos and Deimos (Murdoch et al. 2015; Murchie et al. 2015). Figure 3 shows the peak acceleration, the penetration depth, and the collision duration for the three gravity cases as a function of collision velocity. The measured peak acceleration in plot (a) has a quadratic dependence, but the data is shown on a loglog scale in order to better compare the full range of values.
As illustrated in Fig. 3, the penetration depth of the projectile and the duration of the collision clearly increase when gravity level decreases. This observation is consistent with the experimental results reported in Murdoch et al. (2021). Interestingly, it appears as though the peak acceleration has the same trend for all three gravity cases. It is difficult to draw conclusions from Fig. 3a, however, since there are so few overlapping test cases. The decision to test different collision velocities for different gravity levels was intentional, because the test cases were selected to cover a Froude number range of 0–10, not a collision velocity range of 0–10 m s^{−1}. Following Eq. (6), Fr = 10 for a 10 cm diameter projectile when g = 10 m s^{−2} and v_{c} = 10 m s^{−1}. For the same Froude number, the projectile will only have an impact velocity of 0.1 m s^{−1} when g = 0.001 m s^{−2}. The reason for comparing the collision behavior in terms of Froude number becomes apparent when inspecting Fig. 4.
In Fig. 4, the peak acceleration, the penetration depth, and the collision duration are expressed as dimensionless groups and are plotted with respect to the impact Froude number. Figure 4a is shown on a loglog scale for direct comparison with Fig. 3a. When z_{stop} is normalized by the projectile diameter D and plotted against Fr, the penetration depth collapses onto a single curve for all gravity cases. When t_{stop} is normalized by and plotted against Fr, the collision duration also collapses for all gravity levels. The normalized collision duration appears to be slightly lower for the g = 0.001 m s^{−2} case as compared to the other gravity cases. However, this offset is likely because the simulations for the g = 0.001 m s^{−2} case were terminated earlier than desired (see Appendix B). The only collision parameter that does not collapse by Fr is the peak acceleration measurement. Several possible explanations for this discrepancy are presented in Sect. 3.4.
The curves on Fig. 4 represent the model fit and model predictions produced by the collision model in Sect. 2.3. The constants , h′, and were determined by performing a nonlinear regression and leastsquares fit using Eq. (8) and the penetration depth data. The fit, illustrated by the solid curve in Fig. 4b, was obtained by considering data only within a Froude range of 3–10, because the collisions within this range are assumed to take place completely within regime II. For the given data set, we found that = 0.4, = 4.07 kg, h′ = 1.66 kg. These results indicate that f(z) becomes constant at a depth slightly less than the radius of the sphere. The significance of the and h′ terms will be discussed in Sect. 4.1. The transition from collision regime I to collision regime II is calculated using Eq. (9), and based on the given fit parameters, occurs when Fr_t ≈ 1.5.
The dotted line in Fig. 4b represents the model prediction for the penetration depth in regimes I and the beginning of regime II. In Fig. 4b, the prediction was obtained using the fit parameters and Eqs. (7) and (8). The dotted line in Fig. 4c shows the model prediction for the collision duration in both regimes, as determined using the fit parameters and numerical integration of Eq. (1). Though the model slightly underpredicts the collision duration, it captures the correct form of penetration behavior at the lowest Froude numbers and is a significant improvement of the model used for a comparable study in Murdoch et al. (2021). Unlike the model in Murdoch et al. (2021), the revised model results in a nonzero penetration depth when Fr = 0 (i.e., when v_{c} = 0) and shows that the collision duration increases as expected when v_{c} → 0 (Goldman & Umbanhowar 2008).
Fig. 3 Collision behavior for a 1 kg, 10 cm diameter spherical projectile as a function of collision velocity v_{c} for different gravity levels, where a_{peak} is the peak acceleration, z_{stop} is the penetration depth, and t_{stop} is the collision duration. The measured peak acceleration in plot (a) has a quadratic dependence, but the data is shown on a loglog scale in order to better compare the full range of values. 
Fig. 4 Collision behavior for a 1 kg, 10 cm diameter spherical projectile as a function of the Froude number for different gravity levels, where a_{peak} is the peak acceleration, z_{stop} is the penetration depth, t_{stop} is the collision duration, g is gravity, and D is the projectile diameter. The solid line in plot (b) represents the leastsquares fit to Eq. (8) when 3 < Fr < 10. The dotted line in plot (b) represents the model prediction based on Eqs. (7) and (8), and the dotted line in plot (c) represents the model prediction obtained through numerical integration of Eq. (1). Plot (a) is shown on a loglog scale for direct comparison with Fig. 3a. 
Fig. 5 Evolution of the quasistatic friction force f(z) and the inertial drag force h(z)v^{2} for a 1kg, 10 cm diameter spherical projectile impacting glass beads under terrestrial gravity. F_{d} is the magnitude of the drag force. The vertical dashed lines represent z_{h} and z_{stop}, where z_{h} is the depth where f(z) = h(z)v^{2} and z_{stop} is the final penetration depth of the projectile. 
3.2 Drag forces by collision regime
Murdoch et al. (2017) hypothesize that the inertial drag term dominates collision behavior at decreasing collision velocities as gravity level decreases. To test and better understand this concept, we compare the magnitude of the quasistatic friction force against the inertial drag force throughout each collision and calculate the percentage of the collision where the inertial drag term dominates the interaction. Since the friction drag f(z) depends on the depth of the projectile, the inertial drag component h(z)v^{2} is always superior to the friction component at the start of a collision. As the collision progresses, however, the friction component increases until a depth z_{1}, where it remains constant until the projectile comes to a rest, and the inertial drag component decreases to zero as the projectile’s velocity decreases. The two drag components intersect at some depth z_{h}, where z_{h} is dependent on f_{o}, h, and v_{c}. The ratio of z_{h} to z_{stop} gives the portion of the projectile’s total travel distance where the inertial drag term exceeds the friction term. Figure 5 provides example drag profiles for the cases where g = 9.81 m s^{−2} and v_{c} = 1.5 m s^{−1} (Fr = 1.5), v_{c} = 2.2 m s^{−1} (Fr = 2.2), and v_{c} = 3.5 m s^{−1} (Fr = 3.5). The vertical dashed lines identify z_{h} and z_{stop} with respect to f(z) and h(z)v^{2}. In Fig. 5, the profiles for the quasistatic drag force, f(z), were calculated using the relationship f(z) = f_{o} when z ≤ z_{1} and f(z) = f_{o}z∕z_{1} when z > z_{1}. From Sect. 3.1, f_{o} = 39.9 kg m s^{−2} and z_{1} = 0.04 m. The profiles for the inertial drag force, h(z)v^{2}, were calculated using the relationship h(z)v^{2} = hv^{2}. From Sect. 3.1, h = 16.6 kg m^{−1}. The velocity of the projectile v as a function of depth was found using Eq. (A.9) when z ≤ z_{1} and Eq. (A.16) when z > z_{1} (see Appendix A).
From z = 0 to z = z_{h}, the inertial drag term is superior to the friction term. From z = z_{h} to z = z_{stop}, the friction term is superior to the inertial drag term. As Fr increases, so does z_{h} and z_{stop}. The trend for z_{h}∕z_{stop} as a function of Fr is shown in Fig. 6. Here, the vertical dashed line indicates Fr_{t}, or the transition from collision regime I to collision regime II. The horizontal dashed line shows where the inertial drag term exceeds the friction term for more than 50% of the collision. In regime I, the friction drag dominates the collision, while in regime II, the inertial drag is more important. For the given material and fit parameters, z_{h} ∕z_{stop} ≈ 0.5 when Fr = Fr_{t} (see drag profile in Fig. 5a). Between Fr = 1.5 and Fr = 3.5, the impact occurs within regime II and the inertial drag term dominates the interaction, but the inertial term intersects the friction component when f(z) = f_{o}z∕z_{1} (see example in Fig. 5b). The change in curvature around Fr = 3.5 in Fig. 6 corresponds to the moment where the inertial drag term begins intersecting the friction term at f(z) = f_{o} (see example in Fig. 5c).
Through this analysis, we see that granular interactions occur within an inertial regime sooner (i.e., at lower collision velocities) on small bodies than they do on Earth. Here, the transition from collision regime I to collision regime II, or the transition from a friction dominated system to an inertial system, occurs when v_{c} = 1.56 m s^{−1} on Earth, but when v_{c} = 0.16 m s^{−1} at g = 0.1 m s^{−2}, and when v_{c} = 0.015 m s^{−1} at g = 0.001 m s^{−2}. It is importantto note that the transition Froude number, as well as the relative importance of the quasistatic and inertial drag terms, heavily depend on the properties of the granular surface material. For example, Katsuragi & Blum (2017) study impact dynamics into lowdensity dust agglomerates and find that friction drag dominates the majority of the collision, even for impact velocities of 3.2 m s^{−1}. The influence of different material properties on collision dynamics will be discussed further in Sect. 4.
Fig. 6 Ratio of a collision where the inertial drag force h(z)v^{2} is superior to the quasistatic friction force f(z) for a 1 kg, 10 cm diameter spherical projectile impacting glass beads at different Froude numbers. z_{h} is the depth where the f(z) = h(z)v^{2} and z_{stop} is the final penetration depth of the projectile. The vertical dashed line indicates Fr_{t}, or the transition from collision regime I to regime II, and the horizontal dashed line shows where the inertial drag term exceeds the friction term for more than 50% of the collision. 
Fig. 7 Snapshots of two different impact simulations when g = 10 m s^{−2} and v_{c} = 10 m s^{−1} (Fr = 10) and z = z_{stop}. L is the container diameter and D is the projectile diameter. In image a, the 1 kg, 10 cm diameter projectile is dropped into a 31.5 cm diameter container. In image b, the same projectile is dropped into a 45 cm diameter container. The particles are colored by their initial radial distance from the center of the container. 
3.3 Boundary effects by Froude regimes
Seguin et al. (2008) study the influence of container size on penetration depth for impacts up to 3 m s^{−1} and find that the lateral container walls affect collision behavior when the distance between the projectile and the container walls is less than one projectile diameter. If the projectile is dropped too close to a wall, its penetration depth will be less than if dropped onto an unconfined surface. Goldman & Umbanhowar (2008) report similar findings for a collision velocity range of 1–5 m s^{−1}. Based on the Froude analysis presented in the previous section, it seems logical that boundary effects should be considered as a function of both collision velocity and gravity level, not just collision velocity.
To explore the relationship between boundary effects, collision velocity, and gravity, we performed impact simulations using three different container sizes: L = 3.15D, L = 4.5D, and L = 6.0D. L and D are the container and projectile diameters, respectively. Figure 7 provides snapshots of two different impact simulations when Fr = 10 and z = z_{stop}. The particles are colored by their initial radial distance from the center of the container. When the projectile is dropped into the smaller container, it causes particles to lift across the entire surface of the granular bed, even near the container walls. When the projectile is dropped into the larger container, the particles at the outermost edge of the container remain more or less unaffected by the projectile’s motion. In fact, the particles within two grain diameters of the wall have an average velocity of about 2.7 cm s^{−1} at the end of the test with the 31.5 cm diameter container, but only 0.5 cm s^{−1} at the end of the test with the 45 cm diameter container. Both visually and analytically, it seems as though the size constraints proposed by Seguin et al. (2008) are not applicable for impacts within the upper range of the tested Froude numbers.
Figure 8 shows how the projectile’s penetration depth trends with the Froude number for different container sizes. The filled symbols correspond to the results for the smallest container size (L = 3.15D) and the open symbols correspond to the results for the largest container size (L = 4.5D). Due to the large computation time, simulations were not conducted for the smallest container size and the g = 0.001 m s^{−2} case when Fr > 2. The verticaldashed line in Fig. 8 indicates the transition from collision regime I to collision regime II, and the plot inset provides a zoomedin view of the data where Fr ≤ Fr_{t}. Two different “large container” sizes were tested for the g = 10 m s^{−2} gravity case,L = 4.5D and L = 6.0D. The penetration depth measurements for both of these cases are represented by the same open blue symbol, because the results are nearly identical for the two container sizes.
Expressing the collision behavior in terms of the Froude number helps demonstrate how boundary effects are sensitive to both collision velocity and gravity. Similar to Seguin et al. (2008) and Goldman & Umbanhowar (2008), we see that the projectile’s penetration depth is artificially low when the container diameter is too small. In collision regime I, a minimum separation of one projectile diameter between the projectile and container walls is sufficient to avoid influence from the walls. This regime includes impacts up to 1.5 m s^{−1} in terrestrial gravity, but only 0.16 m s^{−1} when g = 0.1 m s^{−2} and 0.016 m s^{−1} when g = 0.001 m s^{−2}. In comparison, a container diameter of at least 4.5D is required to avoid influence from the walls in collision regime II. We consider that the 4.5D size is sufficiently large because the penetration depth when L = 4.5D and L = 6.0D is the same for the terrestrial gravity case. As discussed in Sect. 3.2, collision regime I is characterized by friction drag while collision regime II is dominated by inertial drag. Since particles are mobilized more easily during impacts in regime II, it follows that a larger wallprojectile separation distance is required to avoid boundary effects.
Fig. 8 Penetration depth as a function of the Froude number for a 1 kg, 10 cm diameter projectile impacting glass beads under different gravity levels. z_{stop} is the penetration depth, v_{c} is the collision velocity, g is gravity, D is the projectile diameter, and L is the container diameter. The filled symbols correspond to the smallest container size (L = 3.15D) and the opensymbols correspond to the largest container sizes (L = 4.5D, or when g = 10 m s^{−2}, L = 4.5D and L = 6.0D). The vertical dashed line indicates Fr_{t}, or the transition from collision regime I to collision regime II, and the plot inset provides a zoomedin view of the data where Fr ≤ Fr_{t}. 
3.4 Peak acceleration measurements
In Sect. 3.1, we see that the penetration depth and the collision duration collapse when expressed as a function of the Froude number. However, the peak acceleration does not. To highlight the differences between the trends, Fig. 9 shows the acceleration, the velocity, and the position profiles for collisions at three different gravity levels when Fr = 6.3. The profiles are normalized by gravity and projectile diameter as shown. The inset of Fig. 9a provides a closer view of the acceleration behavior directly following the impact.
The variations in the peak acceleration measurements might occur for several reasons. The fluctuations in the acceleration profiles are due in large part to the coarse size of the grains and roughness of the surface layer (Goldman & Umbanhowar 2008). To introduce variability between simulation trials, the initial position of the projectile was randomized within a radius of two grain diameters before the start of each test. Some of the differences in the peak measurements can therefore be attributed to changes in the impact point between the projectile and the surface material. The filtering method that is used to process the data also changesthe peak measurements (see Appendix B). Lastly, it is worth noting that the bulk density of surface material varies slightly between the three gravity cases (ranging from approximately 1.36 g cm^{−3} when g = 0.001 m s^{−2} to 1.46 g cm^{−3} when g = 10 m s^{−2}). It is likely that the peak acceleration is affected by the surface density, but more work is required to truly understand the relationship between the peak value, the gravity level, and the surface properties (see Sects. 4.1 and 5).
The sensitivity of the acceleration profile to grain size, material properties, and processing method makes it difficult to extract meaningful information from the peak acceleration parameter. An additional challenge of using a_{peak} to analyze collision behavior is that a clear relationship between a_{peak} and v_{c} has yet to be determined. Murdoch et al. (2021) derive an expression for the peak acceleration as a function of collision velocity using a model where f(z) = f_{o} and h(z) = h. However, the solution assumes that the peak acceleration occurs right at the moment of impact. Figure 10 illustrates how the arrival of the peak changes with the Froude number, and by extension, the collision velocity and the gravity level.
At the higher end of the tested Froude range (Fr ≥ 6.3), the peak acceleration occurs within a few milliseconds of the initial impact. By comparison, the arrival of the peak is delayed with respect to the initial impact for midrange Froude values (Fr ≈ 3–4.5). At the lowest Froude numbers (Fr ≤ 1), the peak becomes almost indistinguishable from the other fluctuations in the profile. Goldman & Umbanhowar (2008) show how the general form of the acceleration profile also changes with surface material. All existing collision models, including this one, fail to reproduce the complex details of the acceleration profile. At the same time, however, they succeed in predicting macroscopic collision properties such as the penetration depth and the collision duration. For the above reasons, the penetration depth appears to be a more reliable metric for interpreting collision behavior than the peak acceleration.
Fig. 9 Example acceleration, velocity, and position profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different gravity levels when Fr = 6.3. The profiles are normalized by gravity g and projectile diameter D as shown. The inset of plot a highlights the fluctuations in the acceleration profiles directly following the impact. 
Fig. 10 Example acceleration, velocity, and position profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different Froude values when g = 10 m s^{−2}. The profiles are normalized by gravity g and projectile diameter D as shown. The inset of plot a highlights how the timing of the peak acceleration measurement changes with the impact Froude number. 
4 Discussion
4.1 Interpreting penetration depth measurements
An undeniable challenge for studying small body surfaces is the lack of opportunity for extensive in situ exploration. The surfaces of asteroids (25143) Itokawa, (162173) Ryugu, and (101955) Bennu have been examined using hopping rovers, impact experiments, and/or touchandgo sampling maneuvers (Kawaguchi et al. 2008; Lauretta et al. 2017; Watanabe et al. 2017). These operations have produced incredibly valuable, but limited data related to spacecraftregolith interactions. As a result, many in the field would like to know how to constrain the material properties of a surface using a single accelerometer or force profile.
As discussed in Sect. 2.3, the quasistatic friction and inertial drag terms in Eq. (1) are proportional to hydrostatic pressure and momentum transfer, respectively. This means that f(z) and h(z)v^{2} have physical links to the size and mass of the projectile and the bulk density and the friction angle of the surface material. If f(z) scales with Coulomb friction and hydrostatic pressure, then f_{o} = μρ_{g}gz_{1}A. If h(z)v^{2} scales with momentum transfer and is proportional to the friction coefficient μ, then h = μρ_{g}A (Katsuragi 2016). Empirical analysis of experimental data has shown that f(z) and h(z)v^{2} also scale by two dimensionless factors, α and η, and that , where ρ_{p} is the density of the projectile (Katsuragi & Durian 2013). The scaling for f_{o} and h then becomes (10)
Katsuragi & Durian (2013) find that α = 8.0 and η = 2.7 for a model where f(z) = kz and h(z) = h, where k and h are constants. In this study, Eqs. (10) and (11) can be used to calculate α and η. Using the fit parameters found in Sect. 3.1 and recalling that , h′ = hD, and , we find that α ≈ 19.6 and η ≈ 3.7. The revised model from Sect. 2.2 can then be reformulated as (12)
The scaling parameters α and η were determined using the mean bulk density of all simulated material beds (ρ_{g} ≈ 1.4 g cm^{−3}) along with the expression μ = tanϕ, where ϕ is the approximate repose angle for the 10 cm diameter glass beads (ϕ ≈ 22°) (Sunday et al. 2021).
The scaled model given by Eq. (12) is a preliminary solution because the scaling expressions for f_{o} and h are still being debated in literature. More work is required to verify how f_{o} and h scale as a function of ρ_{g}, ρ_{p}, ϕ, and D, and Eq. (12) needs to be tested for different impactor sizes and shapes as well as different surface materials. We introduce the preliminary model in this section because it allows us to demonstrate how one can infer certain material properties from a projectile’s penetration depth. Figure 11 illustrates how the penetration depth varies with the bulk density and the internal friction angle of the surface material. In this analysis, we consider impacts for three Froude values, Fr = 0, Fr = 3, and Fr = 10, and we vary either the friction angle of the material (solid black line in Fig. 11) or the bulk density of the material (dashed blue line in Fig. 11).
The case where Fr = 0 corresponds to the initial sinkage of the projectile, z_{o}. The sinkage varies minimally within the given range of the ϕ and ρ_{g} values. By contrast, the penetration depth can more than double for the same range of the ϕ and ρ_{g} values if Fr = 10. The higher the collision Froude number becomes, the more the surface properties influence the projectile’s behavior. Curves similar to those shownin Fig. 11 can be constructed for a given landing event and used to constrain the bulk density or the friction angle of a surface material, as long as the impactor size, gravity, impact velocity, and penetration depth are known.
Fig. 11 Penetration depth for a 1 kg, 10 cm diameter projectile into surface materials with different friction angles ϕ and bulk densities ρ_{g}. The solid black lines represent the case where ρ_{g} = 1.4 g cm^{−3} and ϕ is varied. The dashed blue lines represent the case where ϕ = 22° and ρ_{g} is varied. 
4.2 Application to small body missions
Over the past few years, the Hayabusa2 and OSIRISREx spacecrafts have been visiting and sampling the nearEarth asteroids Ryugu and Bennu, respectively. In the coming decade, the Martian Moons eXploration (MMX) mission will perform a Phobos sample return and a Deimos flyby, and the upcoming DART and Hera missions will contribute to the first asteroid deflection test and rendezvous with the binary asteroid Didymos and its secondary, Dimorphos. The bulk densities of Ryugu and Bennu are measured to be about 1.19 g cm^{−3} (Lauretta et al. 2019; Watanabe et al. 2019), while the densities of Phobos and Didymos are estimated to be higher and around 1.85 g cm^{−3} and 2.17 g cm^{−3}, respectively (Murchie et al. 2015; Naidu et al. 2020). The internal friction angles for the different surface materials are either assumed or are based on modeling and simulation results. The Environment Requirement Document (ERD) produced by the MMX rover team assumes that the regolith layer on Phobos could have a friction angle ranging anywhere 30–50° (S. Tardivel, priv. comm., 2021). Similarly, the friction angles on Ryugu, Bennu, and Didymos are thought to range from 30–40° (Barnouin et al. 2019; Zhang et al. 2021). Using the estimated properties for the different materials, we can predict how far a 1 kg, 10 cm diameter projectile might sink into the surface of each body. Figure 12 shows the estimated penetration depth for the sphere when dropped onto the surfaces of Bennu, Phobos, and Didymos at different Froude numbers. The case for a glass bead surface material is also provided for reference.
The shaded regions in Fig. 12 highlight the transition from collision regime I to collision regime II. In collision regime I (the gray shaded region on the plot), the penetration depth changes only slightly for the different material cases. The difference inthe penetration depth becomes much more pronounced for collisions in regime II (the blue shaded region on the plot). As expected, the highest penetration depth is associated with the lowest friction material (glass beads), and the lowest penetration depth occurs for impacts into the material with the highest bulk density and highest internal friction angle (Didymos and Phobos). It should be noted that the presented collision model has only been tested for collision regime I and when Fr ≤ 10 in collision regime II (the dark blue shaded region on the plot). Impact cases where Fr > 10 have not yet been tested (the light blue shaded region on the plot), but the prediction is thought to provide a reasonable or a lowerboundestimate for the penetration depth in this range. As the impact Froude number increases, the collision behavior might enter a third regime where the presented model is no longer valid (see Sect. 5 for additional discussion).
Table 2 provides an example of the impact Froude numbers that are associated with various in situ operations on past and upcoming small body missions. The Froude number is calculated using the estimated collision velocity v_{c} and the effective projectile diameter D′, where (Kang et al. 2018). The MASCOT, MMX, and Hera landers are cubic in shape (Scholten et al. 2019; Sedlmayr et al. 2020; Goldberg et al. 2019), and for the purpose of this study, are assumed to land on the largest face of the cube. The OSIRISREx and MMX samplers are taken to be closed cylinders, even though they are hollow in the center (Ballouz 2017; Usui et al. 2020).
The surface operations that occurred during the Hayabusa2 and OSIRISREx missions likely took place within the untested region of collision regime II (where Fr > 10). The future Hera CubeSat landing will also occur in the untested region, but the MMX rover will have an impact Froude number of approximately 10, and the MMX sampling operations could take place completely within collision regime I. The impact Froude numbers for these events give an indication of how much resistance the different systems will encounter during their respective surface interactions. The outcome of the actual operations, however, are still highly dependent on the given surface terrain. Ballouz et al. (2021) simulate the touchdown of the OSIRISREx sampling mechanism on Bennu and revise the model by Katsuragi & Durian (2013) for Bennulevel gravity. The authors investigate the influence of friction, packing fraction and cohesion on the sampling operation and observe that certain types of surface materials result in a rebound of the lander. The OSIRISREx TAGSAM penetrated deep into the surface of Bennu (Lauretta & Team 2021). On the other hand, the MASCOT lander appears to have hit a boulder when it touched down on the surface of Ryugu and consequently rebounded (Scholten et al. 2019). The model described in this work is appropriate for granular surfaces that can be treated as a near continuum. Additional testing is required to understand the extent to which the model applies to cohesive surfaces or surfaces with a large dispersion of grain sizes.
Fig. 12 Predicted penetration depth for a 1 kg, 10 cm diameter spherical projectile impacting different granular materials, including a bed of glass beads and surface materials similar to what might be found on Phobos, asteroid (101955) Bennu, and asteroid (65803) Didymos. z_{stop} is the penetration depth, D is the projectile diameter, v_{c} is the collision velocity, g is the gravity level, ρ_{g} is the bulk density of the material, and ϕ is the internal friction angle of the material. The shaded regions on the plot show the transition from collision regime I (gray) to collision regime II (blue). The presented model has been tested against numerical and laboratory experiments in this study for collisions in regime II where Fr ≤ 10 (dark blue shaded region), but not for collisions where Fr > 10 (light blue shaded region). 
Estimated Froude numbers associated with direct surface interactions during recent and upcoming small body missions.
5 Conclusions
In this study, we introduce a revised collision model that can be used to predict the final penetration depth of a projectile during slow granular impacts. The revised model resembles the Poncelet model because it includes both quasistatic and inertial drag components. The inertial drag force is quadratic with velocity, and we assume that the drag coefficient isconstant throughout the entire collision. Then, we assume that the quasistatic friction component is depthdependent during the initial part of a collision, but that the friction force becomes constant after some threshold depth z_{1}. If the projectile comes to a rest before reaching z_{1}, then the impact occurs within collision regime I. If the projectile comes to a rest after passing through z_{1}, then the impact occurs within collision regime II. The penetration depth in regimes I and II can be calculated as a function of impact velocity using Eqs. (3) and (4).
In order to highlight the influence that gravity has on collision behavior, we reformulate the expressions for the penetration depth in regimes I and II in terms of the dimensionless Froude number Fr, where . The new solutions, given by Eqs. (7) and (8), include three parameters that are dependent on the properties of the projectile and the surface material: f_{o}, h, and z_{1}. We treat f_{o}, h, and z_{1} as fit parameters, and we determine their values using simulation data for impacts under three different gravity levels and a Froude number ranging from 3 to 10. Then, we use the proposed model to predict the penetration depth for impacts where Fr ≤ 3 and the collision duration for impacts where 0 ≤Fr ≤ 10. By applying this method, we show that the revised model predicts collision behavior well for impacts where Fr < 10. In addition, we see that the penetration depth and collision duration measurements collapse onto a single curve when expressed as a function of the Froude number.
With the proper scaling, the presented framework can be used to analyze single landing or impact events from small body missions. It should be noted that it is difficult to use the peak acceleration measurement to assess collision behavior, because the relationship between the peak measurement, the surface material, and the collision velocity is still poorly understood. However, accelerometer or force sensor data can be used to obtain estimates for penetration depth and collision duration, which can in turn be used to determine certain surface properties, such as the internal friction angle or the bulk density of the material. More work is needed to verify the scaling relationships for the friction and inertial drag coefficients, and the expressions for f_{o} and h should be developed as a function of cohesion, grain size, grain size distribution, grain shape, internal friction angle, surface density, and projectile shape.
The presented model has not yet been tested for lower gravity levels or higher Froude numbers, because these test cases require larger containers, and consequently, higher computation times. However, collision behavior when g < 0.001 m s^{−2} is expected to follow the revised collision model and the Froude scaling presented in this work. For higher Froude number impacts, the projectile will likely enter a third collision regime where additional forms of energy dissipation must be taken into account (see Omidvar et al. 2014 for an overview of the different types of energy loss that occur for quasistatic through supersonic impacts). An objective of future work will be to assess the validity of the presented model for impacts where Fr > 10 and to identify the transition point between the second and third collision regimes.
Acknowledgements
C.S. is jointly funded by the Centre National d’Études Spatiales (CNES) and the Institut Supérieur de l’Aéronautique et de l’Espace (ISAE) under a PhD research grant. All authors acknowledge support from CNES. P.M. acknowledges funding from Academies of Excellence: Complex systems and space, environment, risk, and resilience, part of the IDEX JEDI of the Université Côte d’Azur, and Y.Z. acknowledges funding from the Doeblin Federation and the “Individual grants for young researchers” program, also part of the IDEX JEDI of the Université Côte d’Azur. P.M. and Y.Z. acknowledge funding from the European Union’s Horizon 2020 research and innovation program under grant agreement No. 870377 (project NEOMAPP). Simulations with Chrono were performed using HPC resources from CALMIP under grant allocation 2020P20020.
Appendix A Derivation of penetration depth expressions
We begin by considering the general force law (A.1)
where m is the mass of the projectile, g is gravity, f(z) is the quasistatic friction force, h(z) is the inertial drag term, and z is the depth of the projectile. As demonstrated in Clark & Behringer (2013), Eq. (A.1) can be recast in terms of kinetic energy K such that (A.2)
where and v is the velocity of the projectile. Eq. (A.2) is a linear ordinary differential equation with the following solution: (A.3)
In Eq. (A.3), K_{o} is the initialkinetic energy of the projectile, (A.4)
The projectile impacts the surface with a collision velocity of v_{c} and comes to a rest at z_{stop}. We assume that the inertial drag term is constant and that h(z) = h throughout the entire collision. The quasistatic friction force, however, takes on one of two forms depending on the projectile’s depth. At first, f(z) increases linearly with depth such that f(z) = kz. Then, f(z) reaches a maximum value f_{o} and remainsconstant until the end of the collision. The transition between the depthdependent and the constant friction states occurs at some depth z_{1}, and this transition separates what we refer to as collision regimes I and II.
If the collision velocity is small enough such that z_{stop} ≤ z_{1}, then the collision takes place completely within collision regime I. In this regime, f(z) = kz, where k = f_{o}∕z_{1}. f_{o} is a constant that depends the properties of both the surface material and the projectile. By evaluating Eqs. (A.4) and (A.5) from 0 to z with h(z) = h and f(z) = f_{o}∕z_{1}, K_{p} (z) and Φ(z) can be expressed as (A.6)
Substituting Eqs. (A.6) and (A.7) into Eq. (A.3) yields (A.8)
The above expression for velocity as a function of depth is valid when z ≤ z_{1}. Eq. (A.8) can also be simplified using the correct boundary conditions to obtain an expression for z_{stop,I}, or the penetration depth in collision regime I. In regime I, and K(z_{stop,I}) = 0, so (A.10)
If the collision velocity is large enough such that z_{stop} > z_{1}, then the collision occurs within regime II. In collision regime II, (A.11)
where z_{stop,II} is the penetration depth in regime II and k = f_{o}∕z_{1}. Since f(z) at the beginning of regime II has the same form as f(z) in regime I, the velocity of the projectile at depth z_{1}, or v_{1} can be determined using Eq. (A.8). During this phase of the collision, and , so (A.12)
The expression for K(z) during the second phase of the collision must be obtained by evaluating Eqs. (A.4) and (A.5) from z_{1} to z. When h(z) = h and f(z) = f_{o}, (A.13) (A.14)
The above expression for velocity as a function of depth is valid when z > z_{1}. The penetration depth in regime II, z_{stop,II}, is determined by applying the correct boundary conditions to Eq. (A.15). During the second phase of the collision, , where v_{1} is given by Eq. (A.12), and K(z_{stop,II}) = 0, so (A.17)
Appendix B Postprocessing of acceleration data
In this study, the projectile’s acceleration information is provided as a direct data output from the impact simulations. The data is reported at a high sampling frequency of 10 kHz. In Sunday et al. (2021), we calibrated the simulation parameters using terrestrial gravity impact experiments. In order to compare the peak acceleration measurements from the impacts experiments and simulations, we filtered the acceleration data using using a second order Butterworth filter with a cutoff frequency of 500 Hz (see Sunday et al. (2021) for more information). This filtering method resulted in a reasonable match between the different measurements because the impact simulation and experiment were conducted under the same gravity level. When comparing impact behavior under different gravity levels, however, the same processing technique cannot be used because the collision time scale varies with gravity. Rather than filtering the data in the time domain t as was donein Sunday et al. (2021), we filtered the data in this study using the normalized time domain t^{*}, where . Fig. B.1 illustrates why it is important to consider the collision time scale when comparing impacts under different gravity levels. The blue, red, and green lines on the plot represent the acceleration profiles for impacts where Fr = 6.3 and g = 10, 0.1 and 0.001 m s^{−2} respectively.Fig. B.1 (a) shows the acceleration profiles when the data is filtered in the time domain t, and Fig. B.1 (b) shows the acceleration profiles when the data is filtered in the normalized time domain t^{*}.
As seen in Fig. B.1, the peak acceleration measurements change drastically with the processing method. Based on the frequency of the fluctuations in the acceleration data, it seems more appropriate to filter the data in the normalized time domain. Though not shown, it should be noted that the filtering method does not influence the velocity or depth profiles. However, the collision duration, just like the peak acceleration, should be determined in the normalized domain rather than the standard domain. In Sunday et al. (2021), the end of the collision was identified as the moment when the projectile’s velocity fell below 1 mm s^{−1}. In this study, the end of the collision is identified as the moment when the normalized velocity v^{*} falls below 1x10^{−2}, where . If v^{*} = 1x10^{−2}, then v = 10, 1, and 0.1 mm s^{−1} when g = 10, 0.1 and 0.001 m s^{−2} respectively. The end criteria for the lowest gravity level is difficult to achieve within a reasonable computation time, so the impact simulations with the g = 0.001 m s^{−2} gravity case often terminated prematurely, when v ≈ 0.05 mm s^{−1}.
Fig. B.1 Example acceleration profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different gravity levels when Fr = 6.3. The acceleration data is normalized by gravity g, and the timeparameter t is normalized by g and the projectile diameter D. Plot (a) gives the acceleration profiles when the data is filtered in the time domain t, and plot (b) gives the profiles when the data is filtered in the normalized time domain t^{*} where . 
Appendix C List of symbols
Table C.1 lists of all of the symbols that are referenced in the main body of this text.
Description of symbols, with units expressed in length L, mass M, and time T.
References
 Albert, R., Pfeifer, M., Barabási, A.L., & Schiffer, P. 1999, Phys. Rev. lett., 82, 205 [NASA ADS] [CrossRef] [Google Scholar]
 Altshuler, E., Torres, H., GonzálezPita, A., et al. 2014, Geophys. Res. Lett., 41, 3032 [NASA ADS] [CrossRef] [Google Scholar]
 Ambroso, M., Kamien, R. D., & Durian, D. J. 2005, Phys. Rev. E, 72, 041305 [CrossRef] [Google Scholar]
 Arakawa, M., Saiki, T., Wada, K., et al. 2020, Science, 368, 67 [NASA ADS] [CrossRef] [Google Scholar]
 Ballouz, R. 2017, PhD thesis, University of Maryland, College Park, USA [Google Scholar]
 Ballouz, R.L., Walsh, K. J., Sánchez, P., et al. 2021, MNRAS, 507, 5087 [NASA ADS] [CrossRef] [Google Scholar]
 Barnouin, O., Daly, M., Palmer, E., et al. 2019, Nat. Geosci., 12, 247 [NASA ADS] [CrossRef] [Google Scholar]
 Bierhaus, E., Clark, B., Harris, J., et al. 2018, Space Sci. Rev., 214, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Brisset, J., Colwell, J., Dove, A., et al. 2018, Prog. Earth Planet. Sci., 5, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Brisset, J., Cox, C., Anderson, S., et al. 2020, A&A, 642, A198 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brucks, A., Arndt, T., Ottino, J. M., & Lueptow, R. M. 2007, Phys. Rev. E, 75, 032301 [NASA ADS] [CrossRef] [Google Scholar]
 Brzinski III, T. A., Mayor, P., & Durian, D. J. 2013, Phys. Rev. Lett., 111, 168002 [CrossRef] [Google Scholar]
 Çelik, O., Tsuda, Y., Yoshikawa, K., Kawakatsu, Y., et al. 2021, Adv. Space Res., 67, 436 [CrossRef] [Google Scholar]
 Cheng, A. F., Rivkin, A. S., Michel, P., et al. 2018, Planet. Space Sci., 157, 104 [CrossRef] [Google Scholar]
 Clark, A. H., & Behringer, R. P. 2013, Europhys. Lett., 101, 64001 [NASA ADS] [CrossRef] [Google Scholar]
 Colwell, J. E. 2003, Icarus, 164, 188 [Google Scholar]
 Colwell, J. E., & Taylor, M. 1999, Icarus, 138, 241 [NASA ADS] [CrossRef] [Google Scholar]
 Costantino, D., Bartell, J., Scheidler, K., & Schiffer, P. 2011, Phys. Rev. E, 83, 011305 [NASA ADS] [CrossRef] [Google Scholar]
 de Bruyn, J. R., & Walsh, A. M. 2004, Can. J. Phys., 82, 439 [NASA ADS] [CrossRef] [Google Scholar]
 Faug, T. 2015, Eur. Phys. J. E, 38, 1 [Google Scholar]
 Ferrari, F., Franzese, V., Pugliatti, M., Giordano, C., & Topputo, F. 2021, Adv. Space Res., 67, 2010 [NASA ADS] [CrossRef] [Google Scholar]
 Goldberg, H. R., Karatekin, Ö., Ritter, B., et al. 2019, in 33rd Annual AIAA/USU Conference on Small Satellites [Google Scholar]
 Goldman, D. I., & Umbanhowar, P. 2008, Phys. Rev. E, 77, 021308 [NASA ADS] [CrossRef] [Google Scholar]
 Hestroffer, D., Sánchez, P., Staron, L., et al. 2019, A&ARv, 27, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Hilton, J., & Tordesillas, A. 2013, Phys. Rev. E, 88, 062203 [NASA ADS] [CrossRef] [Google Scholar]
 Ho, T.M., Jaumann, R., Bibring, J.P., et al. 2021, Planet. Space Sci., 200, 105200 [NASA ADS] [CrossRef] [Google Scholar]
 Housen, K. R., & Holsapple, K. A. 2003, Icarus, 163, 102 [CrossRef] [Google Scholar]
 Jaumann, R., Schmitz, N., Ho, T.M., et al. 2019, Science, 365, 817 [Google Scholar]
 Kang, W., Feng, Y., Liu, C., & Blumenfeld, R. 2018, Nat. Commun., 9, 1 [Google Scholar]
 Katsuragi, H. 2016, Physics of Soft Impact and Cratering, 1st edn., Lecture Notes in Physics 910 (Japan: Springer) [Google Scholar]
 Katsuragi, H., & Blum, J. 2017, ApJ, 851, 23 [Google Scholar]
 Katsuragi, H., & Durian, D. J. 2013, Phys. Rev. E, 87, 052208 [NASA ADS] [CrossRef] [Google Scholar]
 Kawaguchi, J., Fujiwara, A., & Uesugi, T. 2008, Acta Astron., 62, 639 [NASA ADS] [CrossRef] [Google Scholar]
 Lauretta, D., & Team, O.R. T. 2021, Lunar Planet. Sci. Conf., 2548, 2097 [Google Scholar]
 Lauretta, D., BalramKnutson, S., Beshore, E., et al. 2017, Space Sci. Rev., 212, 925 [Google Scholar]
 Lauretta, D., DellaGiustina, D., Bennett, C., et al. 2019, Nature, 568, 55 [Google Scholar]
 Lohse, D., Rauhe, R., Bergmann, R., & Van Der Meer D. 2004, Nature, 432, 689 [NASA ADS] [CrossRef] [Google Scholar]
 Maurel, C., Michel, P., Biele, J., Ballouz, R.L., & Thuillet, F. 2018, Adv. Space Res., 62, 2099 [Google Scholar]
 Michel, P., Kueppers, M., Sierks, H., et al. 2018, Adv. Space Res., 62, 2261 [NASA ADS] [CrossRef] [Google Scholar]
 Michel, P., Ulamec, S., Boettger, U., et al. 2022, Earth, Planets and Space, 74, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Morota, T., Sugita, S., Cho, Y., et al. 2020, Science, 368, 654 [NASA ADS] [CrossRef] [Google Scholar]
 Murchie, S. L., Thomas, P. C., Rivkin, A. S., & Chabot, N. L. 2015, Phobos and Deimos (Tucson: The University of Arizona Press), 451 [Google Scholar]
 Murdoch, N., Sánchez, P., Schwartz, S. R., & Miyamoto, H. 2015, Asteroid Surface Geophysics (Tucson: The University of Arizona Press), 767 [Google Scholar]
 Murdoch, N., Avila Martinez, I., Sunday, C., et al. 2017, MNRAS, 468, 1259 [NASA ADS] [Google Scholar]
 Murdoch, N., Drilleau, M., Sunday, C., et al. 2021, MNRAS, 503, 3460 [NASA ADS] [CrossRef] [Google Scholar]
 Naidu, S., Benner, L., Brozovic, M., et al. 2020, Icarus, 348, 113777 [Google Scholar]
 Omidvar, M., Iskander, M., & Bless, S. 2014, Int. J. Impact Eng., 66, 60 [Google Scholar]
 PachecoVázquez, F., CaballeroRobledo, G., SolanoAltamirano, J., et al. 2011, Phys. Rev. Lett., 106, 218001 [CrossRef] [Google Scholar]
 Sawada, H., Kato, H., Satou, Y., et al. 2021, in 2021 IEEE Aerospace Conference 50100, IEEE, 1–8 [Google Scholar]
 Scheeres, D. J., Hartzell, C. M., Sánchez, P., & Swift, M. 2010, Icarus, 210, 968 [NASA ADS] [CrossRef] [Google Scholar]
 Scholten, F., Preusker, F., Elgner, S., et al. 2019, A&A, 632, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schwartz, S. R., Michel, P., Richardson, D. C., & Yano, H. 2014, Planet. Space Sci., 103, 174 [Google Scholar]
 Sedlmayr, H.J., Barthelmes, S., Bayer, R., et al. 2020, in 2020 IEEE Aerospace Conference, IEEE, 1–10 [Google Scholar]
 Seguin, A., Bertho, Y., & Gondret, P. 2008, Phys. Rev. E, 78, 010301 [NASA ADS] [CrossRef] [Google Scholar]
 Sunday, C., Murdoch, N., Cherrier, O., et al. 2016, Rev. Sci. Instrum., 87, 084504 [NASA ADS] [CrossRef] [Google Scholar]
 Sunday, C., Murdoch, N., Tardivel, S., Schwartz, S. R., & Michel, P. 2020, MNRAS, 498, 1062 [NASA ADS] [CrossRef] [Google Scholar]
 Sunday, C., Zhang, Y., Thuillet, F., et al. 2021, A&A, 656, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Tasora, A., Serban, R., Mazhar, H., et al. 2016, in High Performance Computing in Science and Engineering (Berlin: Springer International Publishing), 19 [Google Scholar]
 Thuillet, F., Michel, P., Maurel, C., et al. 2018, A&A, 615, A41 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Thuillet, F., Michel, P., Tachibana, S., Ballouz, R.L., & Schwartz, S. R. 2020, MNRAS, 491, 153 [NASA ADS] [CrossRef] [Google Scholar]
 Tsimring, L., & Volfson, D. 2005, Powders Grains, 2, 1215 [NASA ADS] [Google Scholar]
 Uehara, J., Ambroso, M., Ojha, R., & Durian, D. J. 2003, Phys. Rev. Lett., 90, 194301 [NASA ADS] [CrossRef] [Google Scholar]
 Usui, T., Bajo, K.i., Fujiya, W., et al. 2020, Space Sci. Rev., 216, 1 [CrossRef] [Google Scholar]
 Van Wal, S., Tsuda, Y., Yoshikawa, K., et al. 2018, J. Spacecraft Rockets, 55, 797 [NASA ADS] [CrossRef] [Google Scholar]
 Watanabe, S.i., Tsuda, Y., Yoshikawa, M., et al. 2017, Space Sci. Rev., 208, 3 [Google Scholar]
 Watanabe, S., Hirabayashi, M., Hirata, N., et al. 2019, Science, 364, 268 [NASA ADS] [Google Scholar]
 Wright, E., Quillen, A. C., South, J., et al. 2020, Icarus, 351, 113963 [Google Scholar]
 Zhang, Y., Michel, P., Richardson, D. C., et al. 2021, Icarus, 362, 114433 [Google Scholar]
All Tables
Estimated Froude numbers associated with direct surface interactions during recent and upcoming small body missions.
All Figures
Fig. 1 Force models applied to collision regime II. The quasistatic drag force acting on the projectile f(z) varies linearly when 0 < z < z_{1} but is constant when z_{1} < z < z_{stop,II}. m is the mass of the projectile, g is gravity, z is the depth of the projectile, z_{1} is the depth where f(z) becomes constant, z_{stop,II} is the final penetration depth in collision regime II, and f_{o}, and h are constants. 

In the text 
Fig. 2 Acceleration, velocity, and position profiles for example impact simulations where g = 9.81 m s^{−2} (blue) and g = 0.1 m s^{−2} (red). For both test cases, the collision velocity of the projectile is 1 m s^{−1}. The projectile’s approximate collision velocity v_{c}, peak acceleration a_{peak}, penetration depth z_{stop}, and collision duration t_{stop} are indicated by the circles, dashed lines, and text annotations on the plot. 

In the text 
Fig. 3 Collision behavior for a 1 kg, 10 cm diameter spherical projectile as a function of collision velocity v_{c} for different gravity levels, where a_{peak} is the peak acceleration, z_{stop} is the penetration depth, and t_{stop} is the collision duration. The measured peak acceleration in plot (a) has a quadratic dependence, but the data is shown on a loglog scale in order to better compare the full range of values. 

In the text 
Fig. 4 Collision behavior for a 1 kg, 10 cm diameter spherical projectile as a function of the Froude number for different gravity levels, where a_{peak} is the peak acceleration, z_{stop} is the penetration depth, t_{stop} is the collision duration, g is gravity, and D is the projectile diameter. The solid line in plot (b) represents the leastsquares fit to Eq. (8) when 3 < Fr < 10. The dotted line in plot (b) represents the model prediction based on Eqs. (7) and (8), and the dotted line in plot (c) represents the model prediction obtained through numerical integration of Eq. (1). Plot (a) is shown on a loglog scale for direct comparison with Fig. 3a. 

In the text 
Fig. 5 Evolution of the quasistatic friction force f(z) and the inertial drag force h(z)v^{2} for a 1kg, 10 cm diameter spherical projectile impacting glass beads under terrestrial gravity. F_{d} is the magnitude of the drag force. The vertical dashed lines represent z_{h} and z_{stop}, where z_{h} is the depth where f(z) = h(z)v^{2} and z_{stop} is the final penetration depth of the projectile. 

In the text 
Fig. 6 Ratio of a collision where the inertial drag force h(z)v^{2} is superior to the quasistatic friction force f(z) for a 1 kg, 10 cm diameter spherical projectile impacting glass beads at different Froude numbers. z_{h} is the depth where the f(z) = h(z)v^{2} and z_{stop} is the final penetration depth of the projectile. The vertical dashed line indicates Fr_{t}, or the transition from collision regime I to regime II, and the horizontal dashed line shows where the inertial drag term exceeds the friction term for more than 50% of the collision. 

In the text 
Fig. 7 Snapshots of two different impact simulations when g = 10 m s^{−2} and v_{c} = 10 m s^{−1} (Fr = 10) and z = z_{stop}. L is the container diameter and D is the projectile diameter. In image a, the 1 kg, 10 cm diameter projectile is dropped into a 31.5 cm diameter container. In image b, the same projectile is dropped into a 45 cm diameter container. The particles are colored by their initial radial distance from the center of the container. 

In the text 
Fig. 8 Penetration depth as a function of the Froude number for a 1 kg, 10 cm diameter projectile impacting glass beads under different gravity levels. z_{stop} is the penetration depth, v_{c} is the collision velocity, g is gravity, D is the projectile diameter, and L is the container diameter. The filled symbols correspond to the smallest container size (L = 3.15D) and the opensymbols correspond to the largest container sizes (L = 4.5D, or when g = 10 m s^{−2}, L = 4.5D and L = 6.0D). The vertical dashed line indicates Fr_{t}, or the transition from collision regime I to collision regime II, and the plot inset provides a zoomedin view of the data where Fr ≤ Fr_{t}. 

In the text 
Fig. 9 Example acceleration, velocity, and position profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different gravity levels when Fr = 6.3. The profiles are normalized by gravity g and projectile diameter D as shown. The inset of plot a highlights the fluctuations in the acceleration profiles directly following the impact. 

In the text 
Fig. 10 Example acceleration, velocity, and position profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different Froude values when g = 10 m s^{−2}. The profiles are normalized by gravity g and projectile diameter D as shown. The inset of plot a highlights how the timing of the peak acceleration measurement changes with the impact Froude number. 

In the text 
Fig. 11 Penetration depth for a 1 kg, 10 cm diameter projectile into surface materials with different friction angles ϕ and bulk densities ρ_{g}. The solid black lines represent the case where ρ_{g} = 1.4 g cm^{−3} and ϕ is varied. The dashed blue lines represent the case where ϕ = 22° and ρ_{g} is varied. 

In the text 
Fig. 12 Predicted penetration depth for a 1 kg, 10 cm diameter spherical projectile impacting different granular materials, including a bed of glass beads and surface materials similar to what might be found on Phobos, asteroid (101955) Bennu, and asteroid (65803) Didymos. z_{stop} is the penetration depth, D is the projectile diameter, v_{c} is the collision velocity, g is the gravity level, ρ_{g} is the bulk density of the material, and ϕ is the internal friction angle of the material. The shaded regions on the plot show the transition from collision regime I (gray) to collision regime II (blue). The presented model has been tested against numerical and laboratory experiments in this study for collisions in regime II where Fr ≤ 10 (dark blue shaded region), but not for collisions where Fr > 10 (light blue shaded region). 

In the text 
Fig. B.1 Example acceleration profiles for a 1 kg, 10 cm diameter projectile impacting glass beads at three different gravity levels when Fr = 6.3. The acceleration data is normalized by gravity g, and the timeparameter t is normalized by g and the projectile diameter D. Plot (a) gives the acceleration profiles when the data is filtered in the time domain t, and plot (b) gives the profiles when the data is filtered in the normalized time domain t^{*} where . 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.