Morphological Diversity of Impact Craters on Asteroid (16) Psyche: Insight From Numerical Models

The asteroid (16) Psyche, target of NASA's “Psyche” mission, is thought to be one of the most massive exposed iron cores in the solar system. Earth‐based observations suggest that Psyche has a metal‐rich surface; however, its internal structure cannot be determined from ground‐based observations. Here we simulate impacts into a variety of possible target structures on Psyche and show the possible diversity in crater morphologies that the “Psyche” mission could encounter. If Psyche's interior is homogeneous, then the mission will find simple bowl‐shaped craters, with a depth‐diameter ratio diagnostic of rock or iron. Craters will be much deeper than those on other visited asteroids and possess much more spectacular rims if the surface is dominated by metallic iron. On the other hand, if Psyche has a layered structure, the spacecraft could find craters with more complex morphologies, such as concentric or flat‐floored craters. Furthermore, if ferrovolcanism occurred on Psyche, then the morphology of craters less than 2 km in diameter could be even more exotic. Based on three to four proposed large craters on Psyche's surface, model size‐frequency distributions suggest that if Psyche is indeed an exposed iron core, then the spacecraft will encounter a very old and evolved surface, that would be 4.5 Gyr old. For a rocky surface, then Psyche could be at least 3 Gyr old.


Introduction
Asteroid (16) Psyche is one of the most intriguing main-belt asteroids. Classified as an M-type asteroid (Tholen, 1984), early spectroscopic measurements showed consistency with iron-nickel meteorites observed in the laboratory, making Psyche a possible candidate parent body for iron meteorites (Bell et al., 1989;Britt et al., 1992;Cloutis et al., 1990). The metal-rich composition of Psyche is supported by radar albedos almost three times higher than those of S-type or C-type asteroids . Radar observations also suggest a high mean bulk density (3,600 + 600/−400 kg/m 3 ) of the upper 1 m , and a uniform metal content across Psyche's surface (Sanchez et al., 2017). The asteroid's thermal inertia measurements are some of the highest recorded, which also suggests a metal-rich surface (Matter et al., 2013).
Despite the compelling evidence for a metal-rich surface and interior composition, more recent density estimates (Drummond et al., 2018;Hanus et al., 2017;Shepard et al., 2017;Viikinkoski et al., 2018) and ground-based spectroscopic observations have challenged the interpretation that Psyche is an intact iron core remnant. Estimates of Psyche's bulk density range between 3.8 and 4.6 g/cm 3 (Drummond et al., 2018;Hanus et al., 2017;Shepard et al., 2017;Viikinkoski et al., 2018), which although much higher than the As the size and morphology of a crater on an asteroid's surface are determined by both the projectile's properties (e.g., mass, velocity), and the target's properties (e.g., strength, porosity, structure), Psyche's craters may help discriminate between alternative models of its internal structure and formation. Reconstructed 3-D shape models of Psyche Viikinkoski et al., 2018) have revealed several large quasi-circular depressions that could be interpreted as impact craters. Shepard et al. (2017) studied the southern hemisphere of the asteroid and identified two distinct units: a wider and shallower depression (D1), ≈70 km in diameter, and a smaller but deeper depression (D2), ≈50 km in diameter. Viikinkoski et al. (2018) studied the northern hemisphere of the asteroid, identifying two further units: a quasi-circular depression 80-100 km in diameter, and an irregularly shaped depression, 90 km wide and 10 km deep.
Asteroid (16) Psyche will be visited by NASA's "Psyche" mission in 2026 . The spacecraft will be equipped with a multispectral imager  capable of resolving craters on Psyche's surface larger than ≈100 m. Study of these craters will therefore help estimate the age of the asteroid surface Neukum et al., 1975), as well as determine the near-subsurface structure and composition. However, due to the unknown and possibly unusual internal structure and geological history of Psyche, the spacecraft may find a much more complex asteroid topography than expected, with impact crater morphologies that are very different from those observed on other visited asteroids.
In this paper, we present numerical simulations of impacts on Psyche analogs that comprise several possible internal asteroid structures. These include a homogeneous dunite mantle material, intact and porous exposed iron cores, and different mantle-core/iron-rock layering configurations. For each impact scenario we investigate the resulting crater size and morphology and use these results to derive scaling relationships. Our results reveal the diversity in possible crater morphologies that the "Psyche" mission could encounter and inform crater model estimates of Psyche's surface age.

Numerical Model
To determine the size and morphology of impact craters on Psyche we used the iSALE2D shock physics code (Wünnemann et al., 2006) to numerically simulate typical asteroid impacts into a variety of possible analog surfaces. iSALE2D is a multimaterial, multirheology extension of the SALE hydrocode (Amsden et al., 1980), specifically developed for simulating impact processes and similar to the older SALEB hydrocode (Ivanov & Artemieva, 2002;Ivanov et al., 1997). iSALE includes a strength model suitable for impacts into geologic targets (Collins et al., 2004) and a porosity compaction model: the − model Wünnemann et al., 2006). iSALE has been extensively validated against laboratory impact experiments (Raducan et al., 2019;Wünnemann et al., 2016), as well as benchmarked against other hydrocodes Pierazzo et al., 2008) for simulating crater size and morphology. RADUCAN ET AL.

Resolution and Regridding
Each numerical impact simulation ran until the crater was completely formed and the ejecta in close proximity to the crater had landed. Some simulations required a high spatial resolution to numerically resolve the near-surface target structure. However, due to the low gravitational acceleration on Psyche (≈0.06 m/s 2 ), ejecta traveling below the escape velocity takes a long time to fall back on to the target. To reduce the computational expense of such simulations, the simulation domain had an initial spatial resolution of 40 cells per impactor radius (cppr) that was coarsened by a factor of 2 after predetermined amounts of time, ending up with a 10 cppr resolution, using the regridding method described in (Raducan et al., 2019).

Impactor Properties
The impacts simulated here considered a nonporous dunite impactor, modeled using the ANEOS equation of state (EOS) and strength model described by Collins et al. (2004), with input parameters summarized in Table 1. The impact velocity was kept at 5 km/s, which represents the average impact velocity in the Main Asteroid Belt (Farinella & Davis, 1992). Craters on Psyche will have been formed by impacts at a range of oblique angles to the target surface, with the most likely angle of 45 • . However, three-dimensional geometry calculations are required to simulate oblique impacts, which are computationally expensive. Moreover, previous studies have shown that the 50% of impacts that occur at angles steeper than 45 • produce craters with circular planforms with little departure from axial symmetry apart from very high speed ejecta (e.g., Elbeshausen et al., 2009). Hence, to perform a large number of simulations in a practical time, we consider vertical impacts only, modeled in two dimensions with an axially symmetric geometry.

Target Properties
To develop a target model for the modern asteroid Psyche, we consider a scenario that assumes that Psyche's parent body began as a Vesta-like asteroid (Ivanov & Melosh, 2013;Ivanov et al., 2010), ≈500 km in diameter, that underwent a long history of collisions before becoming the asteroid we observe today. The structure of a typical once molten, differentiated body has a layered configuration in which a thin basaltic crust and a thick, dense olivine-rich mantle covers an iron core (Bell et al., 1989;Gaffey et al., 1993;McCoy et al., 2006). We assume that all of the basaltic crust and the majority of the silicate mantle was eroded by hit-and-run collisions in Psyche's early history (Asphaug, 2018). We therefore consider Psyche analogs that are made entirely of iron or retained a small volume of remnant mantle.
In all scenarios that include mantle materials we use the same dunite material model to represent the mantle as that used to represent the impactor. The dunite material model uses an ANEOS-derived EOS table (Benz et al., 1989;Thompson & Lauson, 1972) and the pressure-dependent Rock strength model (Collins et al., 2004) in which shear strength reduces with the accumulation of plastic strain from an intact curve to a damaged one (see Table 1).
The nonporous iron material was modeled using an ANEOS-derived EOS table (Ivanov et al., 2010;Thompson & Lauson, 1972) and the Johnson-Cook strength model for iron (Johnson & Cook, 1983), with the input parameters of Armco iron (Johnson & Cook, 1985). The ANEOS-derived EOS table for iron we use was constructed in ca. 2000 by B. Ivanov (personal communication, October 2019) based on parameters by Thompson and Lauson (1972). It is distributed with the iSALE shock physics code (Collins et al., 2016) and has been used in several previous works (Ivanov et al., 2010;Lyons et al., 2019;Potter et al., 2012) to represent the metallic core of asteroids and planets. The EOS includes an approximation of the solid-solid phase transition at ≈11 GPa but does not include a melt transition and is restricted to a maximum pressure and temperature of 2,300 GPa and 86,700 K, respectively. While this tabular EOS is unsuitable for investigating the thermodynamic consequences of high-speed impacts (>20 km/s) it is adequate for the present work, which considers predominantly low-speed (5 km/s) impacts and focuses on strength-dominated crater formation. While ground-based observations of Psyche provide us with some constraints on the near-surface composition (e.g., iron and silicates), the exact structure is still unknown, and will remain so until the arrival of the Psyche spacecraft. To investigate the diversity of possible crater morphologies that might be encountered by the space mission, we modeled six different target scenarios ( Figure 1 and Table 2), that are described below. In all of our simulations, the surface gravity was kept constant, at 0.06 m/s 2 .
(a) Impacts into a homogeneous dunite half-space ( Figure 1a). This scenario considers impacts on Psyche in localized regions of residual thick mantle or before its silicate mantle was stripped off. (b) Impacts into a nonporous intact iron half-space ( Figure 1b). This scenario considers impacts on an exposed dense and intact core, after a long series of hit-and-run collisions (Asphaug et al., 2006) and subsequent impacts removed the entirety of the silicate mantle. (c) Impacts into a shattered, porous iron half-space. If Psyche's composition is almost exclusively Fe-Ni, the asteroid today is most probably not a pristine iron core, but rather a highly fractured body, with a bulk porosity of ≈40% . For this reason, in this impact scenario we assume that the entire mantle was removed by hit-and-run collisions, which also shattered the core, leaving behind a highly porous iron aggregate. At present, a good strength model for fractured, porous iron is not available and there are no impact experiments into such materials. However, several theoretical (Danninger et al., 1993;German, 1994) and experimental (Kutsch et al., 1997) studies have reported a reduction in strength due to porosity in iron-rich materials. Data from Kubicki (1995) for porous iron estimates that 40% porosity causes an 80% reduction in strength (Klar & Samal, 2007). In iSALE, the shear strength, as determined by the Johnson-Cook strength model (Johnson & Cook, 1983), can be reduced by a damage parameter D, which varies between intact (D = 0) and fully damaged with no strength (D = 1). To achieve a strength reduction of 80%, we assigned to the target an initial damage of D = 0.8. This had the effect of reducing the Johnson-Cook parameters A and B by 80% and reducing the strain and the strain rate used in the strength model by the same amount. For the simulations considered here we did not employ the tensile strength model in iSALE. The 40% bulk porosity was modeled using the − porosity model (Wünnemann et al., 2006). The − porosity model parameters for iron were determined by fits to Hugoniot data for porous iron  and took the following values: 0 = 1.7, e = 0, x = 1.15, = 0.97, and = 1.0. (d1), (d2) Impacts into an iron core covered by a layer of silicate regolith. Spectroscopic studies (Landsman et al., 2018;Neeley et al., 2014;Ockert-Bell et al., 2010;Sanchez et al., 2017) and radar observations  found evidence of fine-grained silicate regolith overlaying a metallic bedrock. These measurements suggest that Psyche might have retained some of its mantle; however, we do not have an estimate of this layer thickness to work with. Therefore, we arbitrarily chose to model two scenarios that consider either a 5 or a 10 km dunite layer overlaying an iron substrate. The dunite layer was modeled as in scenario (a), while the iron core was modeled as either an intact core (scenario b) or a shattered, porous core (scenario c). (e) Impacts into a dunite mantle covered by a thin iron layer. One hypothesis that may explain the observed elevated metal content on Psyche's surface is the presence of a thin iron layer covering the mantle as a result of ferrovolcanism (Abrahams & Nimmo, 2019;Johnson et al., 2019). During the cooling period of proto-Psyche, compressional stresses produced in the cooling crust were relieved by faults, allowing Fe-Ni rich material from the core to propagate to the surface in dikes. Abrahams and Nimmo (2019) proposed that repeated ferrovolcanic eruptions on a fully metallic asteroid would be able to raise enough molten iron material to the surface to cover a more dense iron shell. Johnson et al. (2019) showed that ferrovolcanic eruptions are also possible on iron asteroids covered by a thin mantle; however, an estimate of the amount of molten iron erupted was not given. Therefore, here we assumed that the impacts occurred after the iron melt cooled and formed a 50-m layer covering the mantle. The thickness of the iron layer was chosen based on a plausible estimate of the volume of iron that could be erupted in the absence of a mantle (Abrahams & Nimmo, 2019). The iron layer was modeled as nonporous iron (see scenario b), while the mantle substrate was modeled as in scenario (a). (f) Impacts into dunite-iron layered targets. This impact scenario assumed that the thin iron layer from scenario e) was submerged under a 50-m thick layer of mantle material, possibly due to ejecta from impacts after the emplacement of the iron lava by ferrovolcanism. The mantle layer was modeled in the same way as the mantle substrate material in (a).

Simple to Complex Crater Transition
Impact craters are widely distributed on asteroids, and most of them are simple craters; however, above a critical crater diameter, impact craters transition from a bowl-shaped internal structure, to a "complex" one, as they undergo late-stage gravity-driven modification. Complex craters generally have a smaller depth-to-diameter ratio, compared with simple craters in the same target material and can exhibit features such as central peaks, flat floors, or terraced walls.
Based on the observed trend of simple-to-complex transition diameter versus gravity on rocky planetary surfaces, including Vesta, the expected transition diameter on a body with Psyche's surface gravity is 66-95 km.
On the other hand, the transition diameter likely also depends on surface strength and density. If Psyche's near surface is predominantly metal, with an effective strength much higher than a silicate surface, then the transition diameter on Psyche could be significantly larger than 100 km. Given that only the largest observed quasi-circular depressions (Viikinkoski et al., 2018) approach these estimated transition diameters, it is unlikely that any craters found on the asteroid's surface will exhibit a classical complex morphology. Therefore, for simplicity, in this study we neglected the effects of acoustic fluidization that are often included in impact simulations to facilitate complex crater formation. We note that a prerequisite for acoustic fluidization to operate is that the target is pervasively fractured. We would therefore expect acoustic fluidization to operate on Psyche only for scenarios that involve a rocky mantle or prefractured, porous iron core. However, our understanding of complex crater formation is far from complete and any discrepancies between our simulated craters and those observed might plausibly be explained by our omission of acoustic fluidization.

Large Craters Into Homogeneous Targets
The first distinguishable surface features the 'Psyche" mission will encounter will be the largest impact craters on the asteroid's surface. To determine the likely morphology of these large craters on Psyche, as well as the relationship between crater size and impactor size, we first simulated impacts into homogeneous asteroid surfaces, and we considered three scenarios: (a) a nonporous dunite target; (b) a nonporous iron target, and (c) a 40% porous iron target. For each of these scenarios, we varied the dunite impactor radius between 500 m and 15 km and held the vertical impact speed constant at 5 km/s.
All the simulated impacts in these three scenarios formed simple craters: circular, bowl-shaped depressions with raised rims. Impacts into a nonporous dunite target formed wide and shallow craters, with very smooth rims. The crater radius varied with impactor size between 5.5 km, for a 1 km impactor, up to 98 km, for a 30 km impactor, while the crater depth varied between 2.8 and 53 km. The average crater depth to diameter ratio was 0.26, slightly larger than the average ratio on Vesta (Vincent et al., 2014). Figure 2a shows the crater profile resulting from an impact of a 10 km dunite projectile into a nonporous dunite target, which produced an 80 km wide and 21 km deep crater.
The much higher strength and density of iron compared to dunite implied that the craters formed by impacts into nonporous iron (scenario b) were much narrower and shallower than their scenario (a) counterparts. All simulated craters in scenario (b) displayed curled, overturned rims characteristic of laboratory impacts in metal targets (e.g., Libourel et al., 2019;Pierazzo et al., 2008), and their shape appeared to be independent of crater size. Simulated craters in nonporous iron varied in radius between 1.7 km for a 1 km impactor, to 50 km, for a 15 km impactor, while crater depth varied between 1.5 and 43 km. The average crater depth-to-diameter ratio was ≈0.41. Similar depth-to-diameter ratios were measured in craters formed by 10.1029/2020JE006466 recent laboratory impact experiments into metallic iron-rich targets (Marchi et al., 2019). In our simulations, a 10 km dunite projectile produced a crater 34 km wide and 13 km deep (Figure 2b).
For the same impact conditions, impacts into damaged, porous iron targets (scenario c) resulted in craters with sharp but not curled rims, slightly wider and deeper than the craters formed in the nonporous iron targets (scenario b). In this scenario, the crater radius was comparable with the crater depth, resulting in an average depth-to-diameter ratio of 0.65. Figure 2c shows the crater formed by a 10 km dunite projectile, which was 42 km wide and 27 km deep.

Large Craters in Layered Targets
To determine the morphology of craters formed on layered Psyche analogs, as well as the relationship between impactor size and crater dimensions, we simulated impacts into a 5 and a 10 km dunite mantle, covering both a solid iron core (scenario d1) and a porous iron core (scenario d2). The simulated impactors varied in radius between 500 m and 15 km; again the vertical impact speed was held constant at 5 km/s.
Previous experimental and numerical studies of impacts on layered targets showed that the greatest diversity in crater morphology caused by layering occurs when the ratio of upper layer thickness, h, to crater diameter, D, is in the range 0.08 < h/D < 0.3 (Quaide & Oberbeck, 1968;Raducan et al., 2020;Senft & Stewart, 2007). In these studies, which examined crater formation in weak regolith-like targets, underlain by a rocky substrate, the craters exhibited four distinct morphologies depending on the h/D ratio. For large ratios, h/D > 0.3, the upper layer is thick enough that the substrate does not influence the cratering process and the resulting crater is bowl-shaped. For 0.2 < h/D < 0.3, the crater forms mainly in the upper layer, creating flat-bottomed and central-mound craters. For h/D < 0.2, the impactor penetrates the upper layer and the crater forms in both the upper and the lower layers. The result is a concentric crater, which can have both an inner and an outer crater rim.
Although the layering in this work is different-a rocky layer over a denser, stronger metallic layer-we find approximately the same range in crater morphologies with h/D. Therefore, using the previous studies as a guide, we selected a mantle thickness h that would illustrate the full diversity of possible crater morphologies. Figure 3 shows the crater profiles for impacts of varying projectile radius into a 5 km dunite mantle ( Figure  3a) and a 10 km dunite mantle ( Figure 3b) covering a solid iron substrate (left side of each panel) and a porous iron substrate (right side of each panel). In all cases, for large h/D (>0.2), the crater forms in the dunite mantle alone and no core material is excavated. For smaller ratios, h/D < 0.2, the impactor's kinetic energy is largely dissipated in the dunite mantle layer, forming a wide crater in the dunite mantle and only a small shallow crater in the core. This is known as a concentric crater morphology (Quaide & Oberbeck, 1968).
Craters formed in a target of dunite covering a nonporous iron substrate (scenario d1) exhibit the full suite of two-layered target crater morphologies (e.g., bowl-shaped, flat-floor, and concentric craters), similar to craters formed in the regolith on the Moon (Quaide & Oberbeck, 1968). However, craters on dunite covering a porous iron substrate (scenario d2) seem to only create bowl-shaped and concentric craters. This result could reflect the limited number of impactor sizes considered here; however, it could also be due to the less extreme transition in material properties between the dunite and the porous iron layers. In their simulations of impacts in lunar layered targets, Prieur et al. (2018) found crater morphology transitioned directly from bowl-shaped to concentric if the impact velocity was higher than about 12 km/s or when the two layers did not have a high strength difference.

Crater Sizes on Psyche
For the 5 km/s vertical impacts we consider, the crater radius and crater depth for a given impactor radius a can differ dramatically depending on the target material ( Figure 4). An additional complication in determining crater size arises from the ambiguity in defining the crater rim. For homogeneous targets a single rim is produced and the only ambiguity is the altitude at which to measure the crater diameter. In this work, crater radius and depth were measured at the preimpact level. The difference between the rim-to-rim diameter and the preimpact level diameter was ≈1.1 for impacts into iron and porous iron and ≈1.3 for impacts into dunite. However, the large craters formed in layered targets have concentric morphologies, which exhibit characteristic double rims. In these cases, it is ambiguous which rim should be interpreted as the topographic crater rim (for comparison with observations): the outer crater rim, formed in the rocky layer, or the inner crater rim, formed in the iron substrate. For example, for h/D > 0.1, the outer crater rim is much more prominent than the inner crater rim (Figure 3), while for h/a < 0.1 the inner crater rim becomes more prominent than the outer crater rim. For this reason, Figure 4 shows both the outer and the inner rim diameters (dashed lines) as well as the highest crater rim out of the two (solid line). Figure 4 shows crater radius as a function of impactor radius for 5 km/s vertical impacts into each of the Psyche target models considered in this work. Three of the quasi-circular depressions identified by Shepard et al. (2017) and Viikinkoski et al. (2018), which are denoted by D1, D2, and D3 are plotted for comparison. While the radius of the impactor responsible for each of the observed craters will depend on the unknown, unique combination of impactor density, shape, angle and speed, our results for one possible combination (vertical impact of a dunite sphere at 5 km/s) provide an indication of the sensitivity of impactor size to the nature of the cratered target material. For example, if Psyche's surface is rocky, the vertical impact simulations that correspond best to the observed crater diameters are those with impactor radii between 2 and 7 km (Figure 4a). For a solid iron surface and the same impactor parameters, impactor radii between 8 and 15 km produce the same range in crater diameter (Figure 4a). For a layered structure of rock above solid iron, and assuming that we measure the most prominent crater rim, then the three craters could have been formed by impactors with radii between 2 and 15 km, depending in part on the thickness of the rocky mantle (Figures 4a and 4 d). Similarly, if Psyche is a porous, damaged iron core, then the three craters could have been created by impactors with radii between 5 and 12 km (Figure 4b). And for a layered structure of rock above porous iron, the implied impactor radii range between 2 and 12 km (Figure 4b and 4 e).  Crater radius, R, as a function of impactor radius, for impacts into homogeneous targets (dunite, iron, porous iron) and different layering configurations: (a) into a 5 km mantle covering a solid iron substrate (5 km layer); (b) into a 5 km mantle covering a porous iron substrate (5 km layer + porosity); (d) into a 10 km mantle covering a solid iron substrate (10 km layer), and (e) into a 10 km mantle covering a porous iron substrate (10 km layer + porosity). The dashed lines show the inner and the outer crater radii for the craters formed on layered targets. The radii and depths of the crater like depressions identified by Shepard et al. (2017), D1 and D2, and the radius of the Meroe unit (Viikinkoski et al., 2018), denoted by D3, are plotted for comparison. Panels (c) and (f) show the crater depth to diameter ratios. For the craters into layered targets, the hollow symbols show the ratio for both the inner and outer crater rims, while the filled symbols show the ratios for the topographic rim.

Morphology of Small Craters on Psyche
The morphology of a small impact crater is more susceptible to small-scale variations in the surface and subsurface properties of an asteroid than to the large-scale variation. As a result the morphology of small craters might differ substantially from those of large craters in the presence of small-scale target heterogeneity. Ferrovolcanism has the potential to cause substantial small-scale variations in near-surface strength, density, and surface topography on Psyche (Abrahams & Nimmo, 2019;Johnson et al., 2019). As discussed previously, Johnson et al. (2019) showed that repeated ferrovolcanic eruptions can produce a global iron layer on Psyche's surface. To determine the effect of such a layer on small-crater morphology, we conducted iSALE simulations of vertical impacts of 10 to 100 m radius dunite spheres impacting at 5 km/s into a 50 m thick iron layer over a dunite substrate (scenario c). If ferrovolcanic eruptions occurred early in Psyche's evolution it is likely that ejecta blankets produced by subsequent large impacts redistributed rocky material on top of the iron layer. Hence, we also considered a target scenario with a 50 m thick dunite layer over a 50 m thick iron layer covering a dunite substrate (scenario f).  Small, a < 10 m impactors create a bowl-shaped crater when impacting directly in the iron layer. In addition, the strong shock wave transmitted from the iron layer to the rock beneath is sufficient to generate excavation flow and a transient subsurface cavity in the rock beneath it, even though the impactor does not penetrate the iron layer. For a ≥ 25 m impactors, the transient crater exhibits hanging iron rims, that then collapse, causing dunite and iron material to be intimately mixed at the bottom. Similarly, for the a = 50 and a = 100 m impactors, if the impact occurs directly on the iron layer then the resulting crater undergoes major rim collapse owing to the extra weight of the iron in the uplifted rim.
On the other hand, if the iron layer is covered by a rocky layer that is thicker than about one impactor diameter (right-hand side of each panel in Figure 5), then target deformation is restricted to the upper layer and no iron material gets excavated, resulting in a shallow, flat-floored crater. For the a = 25 m impactor, the hanging crater rim mixes with the substrate rocky material. Interestingly, for the a = 50 m impactor, if the impact occurs over a rocky layer rather than directly on the iron layer, then the uplifted iron layer in the transient crater rim remains largely intact and is strong enough to resist complete subsequent collapse. For

10.1029/2020JE006466
impacts >100 m, in scenario (f), the iron is excavated and repositioned outside the crater, leaving the inside of the crater iron-free.

Crater Scaling on Psyche
Our simulation results demonstrate that depending on the target material and structure, the crater radius produced by a given set of impactor parameters can vary considerably. To predict the radius of the crater formed by an impact into a homogeneous target, other than the scenarios considered here, we next interpret our simulation results in the framework of the widely used group scaling relationships (e.g., Holsapple, 1993).
The most general form of the scaling relationship for crater radius, R, as a function of impactor properties (assuming vertical impact of a spherical impactor): velocity, U; radius, a; density, ; mass, m and target properties: density, ; strength, Y , and gravitational acceleration, g, is given by (e.g., Holsapple, 1993): (1) where , , K R1 , and K R2 are empirically determined constants. The velocity scaling exponent, , is constrained between theoretical limits of = 1∕3 if the crater size scales with the momentum of the impactor, and = 2∕3 if the crater size scales with the impact energy (Holsapple & Schmidt, 1987).
To apply Equation 1 to a planetary surface we must determine values for the constants , , K R1 , and K R2 as well as the target properties , g, and Y . A particular challenge is the selection of the appropriate measure of target strength Y , as this is poorly defined and unlikely to be well characterized by a single value (Holsapple, 2009). For this reason, it is common to eliminate one of the material-specific constants K R2 by subsuming it into the definition of Y to define an effective "cratering" strengthȲ (e.g., Holsapple, 1993;Marchi et al., 2019;Prieur et al., 2017). In principle, this allowsȲ to be determined empirically for a given target material (e.g., Marchi et al., 2019); however, in practice this approach requires an independent method to determine K R1 , which is typically only possible with numerical simulations (e.g., Prieur et al., 2017).
The relationship described by Equation 1 is dominated by different terms depending on whether target strength Y or the gravitational stress ∼ ga dominates the resistance to crater growth. For small craters or materials with very high strength, the dominant resisting force to crater growth is the strength, and the crater formation is said to occur in the "strength" regime. In this case, the term 2 = ga U 2 , can be neglected, and crater radius can be expressed as a function of the strength-scale size, 3 where R2 . For large craters or weak target materials, the gravitational stress dominates over strength and the crater is said to form in the "gravity" regime. In such cases, 3 can be neglected and crater radius can be expressed in terms of the gravity-scaled size, 2 = ga U 2 : To determine suitable constants for a rocky, nonporous-iron or porous-iron Psyche surface, we performed additional impact simulations of a 10 km impactor at impact velocities varying between 1 and 25 km/s, to span a large range of 2 and 3 . For impacts into a dunite target our simulations spanned both the gravity and the strength regimes, thus we determined the scaling parameters by fitting Equation 1 to the simulation data ( Figure 6a). As we did not vary the impactor-target density ratio we adopted a fixed value of = 0.4, suggested by previous experimental studies (Housen & Holsapple, 2003;Schmidt, 1980) and varied only , K R1 , and K R2 . To define the target strength Y we used the zero-pressure pressure shear strength (cohesion) from our dunite shear strength model Y d = 10 kPa. The best fit gave K R1 = 0.78, K R2 = 0.60, and = 0.39. For the nonporous and porous iron targets our simulations did not reach the gravity regime and were mostly in the strength regime (rather than the transition regime). In this case, it was only possible to determine the constants H R2 and by fitting Equation 2 to our simulation results (Figures 6a and 6b). For these fits we also assumed = 0.4. For the porous iron targets, the cohesion of the target material accounting for a damaged of 0 where is the equivalent plastic strain, . is the strain rate, T is temperature, T ref is a reference temperature, and T m is the melt temperature. A, B, C, N, and M are constants. Although both temperature and strain can be high near the impact point, neither of these terms is likely to vary significantly with impactor size. On the other hand, it is well documented that iron's yield strength rises significantly at high strain rates and recent laboratory experiments into iron targets (Marchi et al., 2019) have noted the importance of accounting for the effect of strain rate when extrapolating from laboratory-scale crater measurements to planetary scales, or vice-versa. For example, laboratory measurements of the quasistatic yield strength of iron meteorites give values between 200 and 400 MPa (Petrovic, 2001), however at high strain rates, for example, . ≫ 10 3 s −1 , the equivalent strength can almost double (Marchi et al., 2019). We therefore define Y using a simplified version of the Johnson-Cook strength model that includes both the reference shear strength A and the influence of strain rate: We note that the effect of strain rate was only considered for the solid iron targets and not for the predamaged porous iron targets owing to the reduced rate-dependence in this case.
In our simulations we adopt a value of A = 175 MPa for iron (Johnson & Cook, 1983). C is a material specific constant determined from lab experiments and here we use C = 0.06 (Johnson & Cook, 1983). To define the strain rate . for a given impact scenario, we approximate it as the reciprocal of the crater growth time, T g Based on the average strain in the deformed target material, we adopted a strain of one to define the strain rate.
The crater growth time in the strength regime can be expressed as a function of crater radius and target properties using another group scaling relationship (Holsapple & Housen, 2007) which can be rewritten in terms of impactor and target properties by substitution of Equation 2 in Equation 7 to give where C 3 is a constant. By fitting Equation 8 to our simulation data we found C 3 ≈ 0.14.
To determine the scaled crater size for a given impact scenario in iron the procedure is therefore as follows: (a) estimate the crater growth time T g using Equation 8 with Y = A; (b) define the rate-dependent effective strength of the target using Equation 5 with . = 1∕T g ; and finally (c) use this rate-dependent strength in Equation 2 to determine the scaled-crater size. We note that, strictly speaking, steps (a) and (b) should be solved iteratively as T g (Y ) and Y (T g ); however, in practice we find that using Y = A in Equation 8 is sufficient to determine the rate-dependent Y to within a few percent and avoids the need for iteration. We used this procedure to determine Y when fitting Equation 2 to our simulation results for nonporous and porous iron. The best fit parameters for iron were H R2 = 0.48 and = 0.50.
To test the scale-dependence of our strength-regime scaling relationship for nonporous iron, we also applied it to data from recent laboratory-scale impact experiments in metallic iron-rich target materials (Libourel et al., 2019;Marchi et al., 2019). Figure 6c shows the scaled crater radius as a function of impact velocity for iSALE simulations of 10 km dunite impacts into iron targets at 1, 5, 10, 15, 20, and 25 km/s, together with the best fit scaling relationship (Equation 2) to the simulation data. The dotted line on Figure 6c shows the same scaling relationship adjusted to the higher strain rates of laboratory-scale (1 cm diameter) impacts. In other words, the dotted line represents Equation 2 with Y ( . ) defined using Equation 5 as before, but with the crater growth timescale (Equation 8) calculated for an impactor mass, m, appropriate for a cm-scale, laboratory impactor rather than a 10 km wide asteroid. The excellent agreement between the scaling relationship and the experimental data from impacts into Gibeon and iron ingot targets from Libourel et al. (2019) and Marchi et al. (2019) suggests that the scaling relationship is broadly applicable for impacts on an exposed iron target over a large range in impactor size.

Oblique Impacts
Our numerical simulations of impacts on Psyche considered only vertical impacts and the crater scaling relationships above do not account for impact angle. However, impacts on Psyche, as on any planetary body, will occur at a range of impact angles, with a most probable angle of 45 • to the target plane. Previous simulations and experiments of oblique impacts on metallic targets suggest that crater volume and diameter are expected to decrease with impact angle and that a good approach to include the effect of impact angle in crater scaling relationships is to replace the impact velocity in Equations 1 and 2 with the vertical component of the impact velocity (U sin ;Chapman & McKinnon, 1986;Davison et al., 2011;Elbeshausen et al., 2009).
It is also worth noting that while most craters are expected to be circular in planform, in very oblique angle impacts the crater planform is expected to be elliptical, with elongation in the along-trajectory direction.  showed that the transition from a circular to an elliptical crater occurs at a threshold impact angle, c , which depends on the vertical-impact cratering efficiency (the ratio of crater diameter to impactor diameter in a vertical impact scenario). For craters on Psyche, this cratering efficiency is mainly determined by the target strength. The elliptical crater transition angle, and hence the proportion of craters that show significant ellipticity, is therefore likely to be substantially larger for an iron Psyche than for a rocky Psyche. Based on the cratering efficiency in our vertical impact simulations, elliptical craters are expected to occur at impact angles lower than c ≈ 23 • for impacts into iron, while for porous iron the expected threshold angle is c ≈ 20 • . These threshold angles correspond to elliptical crater populations of 15% and 12% of all craters, respectively. Similar critical impact angles ( c ≈ 20 • ) have been observed for oblique impacts into aluminium . For the 5 km/s impacts into dunite considered here, which span the transition from the strength to the gravity regime, the elliptical crater threshold impact angle is expected to vary with the size of the impactor, between c ≈ 13 • for 500 m impactors and 17 • for 10 km impactors. In this case, elliptical craters would constitute less than 5-9% of all craters.

Size-Frequency Distribution of Craters on Psyche
To determine an estimate of the size-frequency distribution of craters on Psyche, we used the CHESS (Collisional Histories in the Early Solar System) Monte Carlo model (Davison et al., 2013(Davison et al., , 2014. CHESS uses the results from collisional and dynamical evolution models, coupled with the group scaling relationships to estimate the impact history of a given body. The size-frequency distribution of the population of potential impacting asteroids varied though time and was derived from collisional evolution models (O'Brien, 2009;O'Brien et al., 2006O'Brien et al., , 2007, and the intrinsic collision probability was taken from those models for the first 100 Myr and then set to a fixed value of 2.86 × 10 −18 km −2 yr −1 (Bottke et al., 1994). The impact velocity was randomly selected from a Maxwellian distribution with a mean of 5.3 km/s. The impact angle was chosen randomly, with P(> ) = cos 2 (Gilbert, 1893;Shoemaker, 1962), where P(> ) is the probability of the impact angle being greater than .
To apply the CHESS model to Psyche, the target asteroid diameter was set to 225 km and the model was run for several different Psyche surface ages between 1 and 4.5 Gyr. In each case, 10 5 target asteroids were considered to produce a statistically meaningful result. CHESS incorporates the group scaling relationships defined above (Equations 1-3), modified to include the effect of impact angle to the horizontal by replacing the velocity term with the vertical component of the velocity (U sin ; Chapman & McKinnon, 1986;Davison et al., 2011;Elbeshausen et al., 2009). To determine how the size of the craters produced by the population of impacting asteroids might differ depending on the target material on Psyche's surface we applied our derived scaling exponents for dunite, iron, and porous iron. For impacts into iron and porous iron we employed the procedure described above for determining the appropriate strain rate-dependent strength for use in the group scaling relationships.
Because almost no data exists on the cratering record on Psyche, the aim of this analysis was not to determine a specific age for Psyche, but rather to determine age constraints for different target scenarios and determine whether the possible few observed craters are consistent with the expected crater population on a rocky or an iron-rich surface. Figure 7 shows the crater size-frequency distribution (SFD) for a Psyche-sized asteroid that experienced the same impact history as the main belt population, over the last 2 Gyr (Figure 7a), the last 3 Gyr (Figure 7b), the last 4 Gyr (Figure 7c), and 4.5 Gyr (Figure 7d). The observed crater-like depressions identified by Shepard et al. (2017) and Viikinkoski et al. (2018) are plotted on top of the SFDs. Horizontal error bars represent the uncertainty in the size of those craters, while the vertical error bars represent the uncertainty in the number of craters on Psyche. This uncertainty is calculated as where N is the cumulative number of craters observed on Psyche. Here we have included the tentative fourth ≈ 90 km depression identified by Viikinkoski et al. (2018) in the observed crater population; the uncertainty of this potential extra crater is smaller than the √ N uncertainty.
The group scaling relationships predict the crater diameter for a given impact measured at the preimpact surface, while from observations crater diameter is usually measured from rim to rim. The ratio between the rim-diameter and the preimpact level diameter D r /D for our Psyche impact simulations in rock was ≈1.3, which is consistent with the ratio for experimental craters in sand and weak soils (Housen & Holsapple, 2011). We used this ratio in our simulations to convert from preimpact level diameter to rim-to-rim diameter for impacts into rocky mantle. For impacts into iron and porous iron we adopt D r /D≈1.1, based on our Psyche impact simulation results, to convert from the preimpact level to rim-height crater diameters (Figure 7).
For a dunite target scenario, if the proposed large craters on Psyche are confirmed they would suggest a surface age that is at least ≈2-3 Gyr old, and up to about 4 Gyr. If, on the other hand, Psyche's surface is dominated by iron or porous iron, the same size craters are consistent with model crater SFDs of a surface that is ≈4.5 Gyr old. This suggests that any mantle stripping events and large impacts must have occurred very early in solar system formation.
Our results suggest that for a very old Psyche, having an iron-rich surface is not only possible, but even likely. However, these surface age estimates are highly uncertain owing to the large uncertainties in the observed 10.1029/2020JE006466 crater sizes. More precise measurements of the craters, together with a better understanding of the surface composition on Psyche will yield a much more accurate estimation of surface age.

Limitations and Future Work
This study used a simplistic and poorly constrained model to describe the strength of porous, fractured iron. Future studies should investigate alternative strength models, that could include more sophisticated strength models for metals that include damage and or a strength model for granular materials that might better represent an iron rubble pile.
Another limitation of this study is that in all simulations we assumed a planar target even though the diameters of the largest craters are comparable to the radius of Psyche. The formation of the largest craters on Psyche will likely be influenced by the target curvature and future simulations should investigate how this might effect crater scaling and morphology. These craters might be wider and shallower than a similar impact into a planar target (Billingham & Collins, 2010).

Conclusions
The asteroid Psyche is one of the most intriguing main-belt asteroids. Earth-based observations suggest that the asteroid has an iron-rich surface; however, the low bulk density estimates and the detection of hydrated minerals on the surface have questioned the interpretation that Psyche has a predominantly metallic composition. To aid efforts to infer Pysche's subsurface properties from its craters, here we conducted numerical simulations of impacts into Psyche analogs that considered several possible target surface materials and structures: a homogeneous rocky (dunite) target; an exposed, dense iron core; a porous, fractured exposed iron core; as well as layered structures with a rocky mantle covering a dense or porous iron core, and a thin near-surface iron layer within or atop the rocky mantle.
Our simulations reveal a wide diversity of potential crater morphologies. The classical "complex" craters observed on large planets, moons, and even large asteroids are not expected owing to the high strength and low gravity of Psyche's surface. However, because of the large difference in strength between rock and iron, the presence of near-surface iron-rock layering or other forms of heterogeneity is likely to result in a range of crater morphologies, many of which may be unique or witnessed for the first time. If Psyche's near-surface structure is homogeneous, then the spacecraft is expected to find simple, bowl-shaped craters, although they will be much deeper than those on other visited asteroids and possess much more spectacular rims if the surface is dominated by metallic iron. On the other hand, if Psyche has a layered structure, the mission could also find shallow craters with flat floors or concentric craters with at least two prominent rims.
If the surface strength is more representative of iron than rock, then the cratering efficiency on Psyche will be substantially lower than on more familiar rocky asteroids. Consequently, the threshold impact angle (relative the horizontal) to produce elliptical craters will be greater, resulting in a higher proportion of elliptical craters on Psyche compared with other asteroids.
Small craters, more susceptible to small-scale variations in the target structure, have the potential to exhibit even more exotic morphologies. If ferrovolcanism occurred on Psyche and produced thin surficial or near-surface lenses of iron, then craters less than about 1 km in diameter in such areas might display a range of complex morphologies, including central mounds of iron or over-hanging metallic rims covering a larger cavity in the rocky mantle below. Ferrovolcanism could have also brought iron-rich material to the surface, and iron mixed with rock would reside on the surface around small craters. Moreover, as ferrovolcanism requires a relatively thin or absent mantle, direct excavation of the iron core by large impacts offers another mechanism for iron exhumation, in addition to volcanism.
Several quasi-circular depressions with diameters larger than 50 km have been identified in 3-D reconstructions of Psyche and interpreted as possible impact craters. If the observed depressions on Psyche are impact craters, then model crater size-frequency distributions based on our results suggest that the Psyche spacecraft will most likely find an old asteroid surface that is at least 3 Gyr old. Moreover, the presence of three or four large D > 50 km craters on Psyche seems to be more compatible with a surface material that is rocky than made of iron. If Psyche's surface was predominantly nonporous iron we would expect fewer than three D > 50 km craters, even if the asteroid's surface was as old as the solar system. However, more accurate measurements of the large craters on Psyche are needed in order to better estimate the surface age or to infer the surface material, which might not be available until the arrival of the spacecraft at Psyche.