Transfer, loss and physical processing of water in hit-and-run collisions of planetary embryos

Collisions between large, similar-sized bodies are believed to shape the final characteristics and composition of terrestrial planets. Their inventories of volatiles such as water are either delivered or at least significantly modified by such events. Besides the transition from accretion to erosion with increasing impact velocity, similar-sized collisions can also result in hit-and-run outcomes for sufficiently oblique impact angles and large enough projectile-to-target mass ratios. We study volatile transfer and loss focusing on hit-and-run encounters by means of smooth particle hydrodynamics simulations, including all main parameters: impact velocity, impact angle, mass ratio and also the total colliding mass. We find a broad range of overall water losses, up to 75% in the most energetic hit-and-run events, and confirm the much more severe consequences for the smaller body also for stripping of volatile layers. Transfer of water between projectile and target inventories is found to be mostly rather inefficient, and final water contents are dominated by pre-collision inventories reduced by impact losses, for similar pre-collision water mass fractions. Comparison with our numerical results shows that current collision outcome models are not accurate enough to reliably predict these composition changes in hit-and-run events. To also account for non-mechanical losses, we estimate the amount of collisionally vaporized water over a broad range of masses and find that these contributions are particularly important in collisions of ∼\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sim $$\end{document} Mars-sized bodies, with sufficiently high impact energies, but still relatively low gravity. Our results clearly indicate that the cumulative effect of several (hit-and-run) collisions can efficiently strip protoplanets of their volatile layers, especially the smaller body, as it might be common, e.g., for Earth-mass planets in systems with Super-Earths. An accurate model for stripping of volatiles that can be included in future planet formation simulations has to account for the peculiarities of hit-and-run events and track compositional changes in both large post-collision fragments.


Introduction
Collisions on vastly different size scales are ubiquitous during all stages of planet formation. Varying bulk densities, probably due to different ice mass fractions, of outer Solar system satellites and KBOs, for example, indicate a rich collisional history. Even though there is currently a lot of discussion about whether planetary embryos were formed by planetesimal accretion followed by oligarchic growth, or rather by rapid accretion of small pebbles (Lambrechts and Johansen 2012), it is well established that the late stages of planet formation-starting once the gas disk disappears and not enough background material for dynamical friction is left-are dominated by giant collisions between embryos, with sizes from hundreds to thousands of kilometers. This includes mixing of material that condensed at vastly different orbital radii to some degree (Raymond et al. 2014). Isotopic measurements as well as dynamical arguments point to material from-or at least isotopically similar to-the outer asteroid belt as the origin of the bulk of Earth's water. The delivery of volatiles such as water, in our Solar system as well as in extrasolar systems, is intimately connected to the long-term evolution of volatiles interacting with the environment in the form of gradual loss of atmospheric constituents (see, e.g., Odert et al. 2017) or sublimation of exposed surface ice (e.g., Schorghofer 2008), and to the collisional evolution of volatile-carrying bodies. In this study, we deal with the latter.
While impacts of small bodies onto much larger ones are usually mostly accretionary, it is well known that collisions between objects of comparable size can result in a much more diverse spectrum of outcomes. Leinhardt and Stewart (2012) developed a comprehensive analytical model to distinguish different collision outcome regimes. Low-velocity impacts result in accretion of at least some projectile material onto the target, while more energetic collisions often lead to erosion or even disruption. For sufficiently oblique impact angles (Genda et al. 2012), an additional outcome is possible-hit-and-run. These grazing encounters are characterized by two large post-collision fragments, originating from the target and the projectile. Therefore, hit-and-run is usually defined via an accretion efficiency (Asphaug 2009) close to zero, given by ξ = (M lf − M targ )/M proj , with the mass of the largest post-collision fragment M lf and of the (pre-collision) target and projectile. This means that the target body neither accretes a lot of projectile material nor it is substantially eroded (cf. Fig. 5). The hit-and-run regime covers a wider range of impact velocities the more oblique collisions are, but is also a strong function of the projectile-to-target mass ratio, with a higher probability for hit-and-run the more similar-sized impacting bodies are. Hit-and-run encounters are particularly interesting, but also complex, when it comes to transfer and loss of volatiles. While in the other outcome regimes there is at most one large fragment and usually some amount of orders-of-magnitude less massive debris, changes of the volatile fraction on both large survivors are of interest for hit-and-run, as well as transfer between the initial reservoirs on projectile and target in the course of the impact. The impact energy is partitioned between the colliding bodies; therefore, the smaller one of the colliding pair is more affected, resulting in large-scale deformation and stripping of its outer layers, particularly of volatiles. While central impacts are dominated by shocks, gravitational stresses become more important in oblique encounters, and the most grazing collisions are almost entirely controlled by tidal forces. Hit-and-run is also important for shaping the final spin of forming planets (Kokubo and Genda 2010) and is considered an efficient mantle stripping mechanism to explain for instance Mercury's high bulk density, or the formation of Earth's Moon (Reufer et al. 2012).
Current N-body simulations often include relatively simple models for collisional fragmentation, which were developed in recent years (Kokubo and Genda 2010;Genda et al. 2012;Marcus et al. 2010;Leinhardt and Stewart 2012;Stewart and Leinhardt 2012). These models are applied to compute the basic fragmentation outcome of similar-sized collisions (e.g., Chambers 2013;Quintana et al. 2016), and also in studies on compositional evolution and mantle stripping for modeling material fractionation (e.g., Carter et al. 2015, and references therein). A more direct, but computationally very demanding, approach are simulations that combine the long-term dynamical evolution and the hydrodynamics of individual collisions in hybrid codes (Genda et al. 2017). In addition to numerical (SPH) simulations, some analytical work on volatile losses by impacts is available, as, e.g., in Schlichting et al. (2015). While most of the mentioned studies preferentially treat refractory elements, a comprehensive model for volatile transfer and loss in similar-sized collisions is still not available. Therefore, current work on volatile mixing and water delivery usually assumes perfect conservation of volatile inventories in collisions (e.g., Raymond et al. 2004;Izidoro et al. 2013), even though such material is especially prone to loss processes, either mechanically, thermally, or in connection with the stellar environment. The reasons are not only its volatile behavior, but also its preferred location in the outer layers of sufficiently differentiated bodies. Studies on collisional water loss mostly distinguish merely material that is gravitationally bound to the largest (Canup and Pierazzo 2006), or the most significant (Maindl et al. , 2017 post-collision fragments, and usually do neither trace the fate of individual inventories on the two large hit-and-run fragments, nor volatile transfer between projectile and target during the collision. Often these studies do not consider dependencies on all important parameters and commonly fix especially the mass ratio of the colliding bodies. We try to close this gap in this work with a dedicated study on volatile transfer and loss in hit-and-run collisions. The obtained results and insights are a necessary basis for the development and inclusion of realistic models for volatile transport in future planet formation simulations.
We describe the applied numerical methods in Sect. 2 and our set of simulation scenarios in Sect. 3. The results are presented and described in Sect. 4, before we discuss them and conclude in Sect. 5. A description of our relaxation approach, as well as collision outcomes and results for water transfer and loss for all computed scenarios, is included in "Appendix".

Methods
We performed smooth particle hydrodynamics (SPH) simulations of collisions between similar-sized planetary embryos. Our SPH hydrocode utilizes GPU hardware to allow for high-resolution runs with still reasonable computing times and has been successfully applied to different aspects of collision and impact processes Dvorak et al. 2015;Haghighipour et al. 2016). In the following, we briefly summarize only its main features and refer to Schäfer et al. (2016) for a comprehensive description. The SPH code includes self-gravity, and in addition to hydrodynamic objects it is also capable of modeling full solid-body physics, with a von Mises plasticity model (Benz and Asphaug 1994) and a brittle failure/fragmentation model introduced to SPH by Benz and Asphaug (1995). This model is based on a Weibull distribution of flaws, with parameters from Nakamura et al. (2007) for basalt and from Lange et al. (1984) for water ice. The von Mises yield criterion does not consider the pressure dependence of shear strength and is therefore not the ideal choice for geologic materials, where shear strength is in general a complex function of pressure, and to a lesser degree also of temperature, strain, and even strain rate. These issues are discussed further in Sect. 4.4. To overcome the problems associated with different particle masses of the standard SPH method, we apply the modified SPH approach by Ott and Schnetter (2003) in all simulations.
The thermodynamical behavior is modeled by means of the nonlinear (Tillotson 1962) equation of state (eos), applicable over a wide range of physical conditions (see also Melosh 1989). The Tillotson eos has two distinct analytical forms, covering different regions of density and (specific) internal energy e. For compressed regions ( ≥ 0 ) and cold expanded states ( < 0 ; e < e iv -the energy of incipient vaporization), it reads with η = / 0 and μ = η − 1. Expanded and vaporized states ( < 0 ; e ≥ e cv -the energy of complete vaporization) are described by In the partial vaporization regime (e iv < e < e cv ), a weighted average of (1) and (2) is used to interpolate. For cold expanded states (e < e cv and < 0 ), a low-density pressure cutoff is applied by setting p = 0 for / 0 < 0.9 to avoid unphysical tension in states where the material rather forms droplets or fractures instead of remaining a continuum. The required material parameters for 0 , e 0 , e iv , e cv , A, B, a, b, α, β for basalt and water ice are from Benz and Asphaug (1999), and those for iron from Melosh (1989). To estimate the amount of vaporized material after a collision, we define vaporization as e ≥ e cv and < 0 , i.e., falling into region (2) of the Tillotson eos. For water, these values are 0 = 917 kg/m 3 and e cv = 3.04 MJ/kg. However, the water vapor fractions derived by this method should be considered rather qualitative estimates, while for a robust quantitative assessment other, more complex eos like ANEOS (see, e.g., Melosh 2007) are required. While ANEOS would provide a thermodynamically consistent treatment of phase changes, it comes with other drawbacks, like a requirement for higher resolution (Reinhardt and Stadel 2017). The simulated bodies were initialized with relaxed internal structures following selfconsistent, semi-analytically computed hydrostatic profiles, with internal energy values from adiabatic compression. The details of this relaxation approach are summarized in "Appendix A". In the initial configuration, the colliding bodies are placed apart at a distance of several times the sum of their radii to allow for buildup of tidal effects and are sent on analytical twobody orbits that lead to the desired collision parameters. All simulations were computed until a final state, characterized by a large distance of the major fragments (typically dozens of times the sum of projectile and target radius). The identification of the final fragments is then done as a post-processing step. First, clumps of spatially connected SPH particles are identified by a friends-of-friends algorithm. Then, we iteratively search for gravitationally bound aggregates of such clumps, where we first compute a seed, consisting of the most massive clump and all other clumps that are mutually gravitationally bound to it. This aggregate of clumps (its barycenter) serves as the starting point of an iterative procedure, where we go through the Fig. 1 Collision geometry in a frame centered on the target, with the conventional definition that the target is the larger body (left) and also from the perspective of the smaller one of the colliding pair (right). The impact angle α spans values up to 90 • and is 0 • for a head-on collision. The impact velocity vector is labeled v 0 , where the actual impact velocity is v 0 = |v 0 | list of clumps, check whether they are gravitationally bound to it and add/remove them if necessary from the aggregate of clumps. Once this procedure converges, we identify the result as the largest fragment and repeat it once again to identify the second largest fragment as well, by starting with the largest still unbound clump. (Clumps that are already included in the largest fragment are always left there.) Convergence toward a self-consistent final state is usually achieved after few iterations. A variety of similar methods is commonly applied to identify the largest collisional fragments, which either include a preceding friends-of-friends step (e.g., Genda et al. 2012;Benz and Asphaug 1999), or start identifying gravitationally bound mass directly on the individual particle level (as, e.g., in Movshovitz et al. 2016;Asphaug 2010;Marcus et al. 2009), or apply both approaches (e.g., Genda et al. 2015). In a hit-and-run collision, the two largest fragments usually originate from the target and the projectile, respectively.

Scenarios
Our scenarios comprise collisions between differentiated, embryo-sized bodies. Most of them have a 3-layered structure, consisting of an iron core (25% by mass), a silicate basalt mantle (65%) and a water (ice) shell (10%). Including volatile inventories in both colliding bodies allows us to track transfer from the projectile to the target and vice versa, as well as losses from both bodies, from the results of a single simulation run, by tracing the respective SPH particles (target vs. projectile water). Therefore, this eliminates the necessity to run an otherwise similar collision twice, once with a water-covered target hit by a dry projectile and once the other way around (see Sect. 4.3). Even with the composition fixed the parameter space is large, where the impact velocity in units of the mutual escape velocity 1 v/v esc , the impact angle α and the bodies' mass ratio γ = M proj /M targ are the most important ones. Figure 1 illustrates the definition of the impact velocity and angle, both at the moment of 1 The two-body escape velocity is given by v esc = √ 2 G M/r with the gravitational constant G, the combined mass M = m 1 + m 2 and the bodies' separation r . first contact, assuming spherical objects. While the total colliding mass M tot is not of crucial importance to some aspects of giant collisions (Asphaug 2010;Genda et al. 2012)-meaning that outcomes are scale invariant-recent results have shown that this holds only limited for transfer and loss of volatiles .
The range of parameters we chose for the simulation runs is focused on the region in parameter space occupied by hit-and-run collisions and is also typical for an active planet formation environment. We have not covered our whole parameter space uniformly, but started with its central region and performed exploratory simulations in various directions we deemed promising for gaining further insights (see Fig. 3). This was mainly done to reduce computational costs and the effort for post-processing in parameter regions where no new results were expected. The impact velocity varies between 1.5 and 5.5 × v esc . The impact angles include some head-on scenarios for comparison, but are otherwise rather grazing as is typical for hit-and-run events, and range from 30 • to 90 • . Mass ratios vary from 1:2 up to 1:50 according to the traditional definition that the target is always the larger body, but from the perspective of a single body which is hit by an arbitrarily large projectile the range of mass ratios encompasses values between 1:50 and 50:1 (cf. Fig. 4). The total mass is 10 23 kg (∼ Moon-mass) in most scenarios, representing typical planetary embryos, but was also varied between 10 22 (∼ 1/10 Moon-mass) and 10 25 kg (∼ Earth-mass), to study the dependence on M tot and thus impact energy, especially in the context of water vapor production. Most simulations were computed with 100k SPH particles, which was confirmed to be sufficient for our purposes by similar results of selected scenarios, which were run again with several higher resolutions up to 2.25 million particles (Sect. 4.6). In the majority of scenarios, the objects are purely hydrodynamical, which is common in this mass range and justified by the dominance of gravitational stresses over material strength. However, in order to clarify the possible influence of material strength we compare them also to results obtained with our solid-body material model (Sect. 4.4).

Results
The results of all simulated scenarios are summarized in Tables 1 and 2 in "Appendix B". In recent studies on collisional water loss, results are often presented as the escaping fraction of the colliding bodies initial water inventory, either w.r.t. the two largest fragments in case of a hit-and-run (e.g., Maindl et al. 2017), all significant 2 fragments , or considering only the largest one at all, regardless of the actual collision outcome (e.g., Canup and Pierazzo 2006). Especially the latter can result in misleading conclusions, since the smaller body in a hit-and-run encounter can carry large amounts of volatiles with it, which are then counted as lost, even though they are still part of a large body, possibly similar in size to the largest (Fig. 5). In order to connect to these studies and summarize results, we provide a similar plot in Fig. 2, where water loss refers to the fraction of water bound to neither of the two largest fragments after the collision 3 . However, we do not think that this is a good representation for studying transfer and loss of volatiles in hit-and-run collisions for two reasons, (1) because it does not contain any individual information on the two large postcollision fragments, and (2) because the amount of lost volatiles relative to the pre-collision inventory gives only an approximate figure of the fragments' actual post-collision water mass ig. 2 Water loss, defined as escaping the two largest fragments, as a function of impact velocity, for different impact angles α and mass ratios γ . Connected points and single large symbols are our numerical results, and not connected small ones are predictions of the model from Leinhardt and Stewart (2012, see Sect. 4.7). Note the differences in the y-axis range fractions 4 (wmf), since they can loose (or gain) other material as well in the process. To avoid these issues, we will discriminate the fate of the two largest hit-and-run fragments in the rest of this study and furthermore state results as changes in wmf-relative to overall fragment masses. We believe this measure is more instructive for tracking the compositional evolution of growing planets, than the change relative to the initial (pre-collision) water inventory. Figure 3 shows an illustration of the bulk of our results following these conventions.
While vaporization (mainly due to shocks) of significant amounts of refractory material happens only in the most energetic events, the situation is different for volatiles, owing to their much lower vaporization energy. Once vaporized, and perhaps heated to high temperatures, this material can escape either thermally, or due to interaction with the (stellar) environment, driven mainly by extreme ultraviolet (EUV) radiation of the young star. Odert et al. (2017) studied the loss of steam atmospheres, which is a complicated function of many factors, like the masses of the body and the atmosphere, its thermal state, the distance to the star and its activity, and the frequency of impacts (see also Schlichting et al. 2015). Due to the high uncertainties associated with this large number of influences (and with the determination of The pictograms illustrate the collision geometry, the sizes of projectile and target before the collision (light gray) and afterward (color), as well as the remaining post-collision wmf, which are initially 0.1. Circle sizes are proportional to mass 1/3 , assuming a constant bulk density. The enclosed rectangular region is further examined from a one-body perspective in Fig. 4, and the connecting line refers to the scenarios shown in Fig. 5 the vapor fraction, cf. Sect. 2), we do not directly include loss of vaporized material in the presented post-collision water contents, but rather estimate the amount of vaporized water for some selected scenarios in a separate Sect. (4.5) and discuss the implications thereof. Thus, the resulting water contents presented here include all material that is gravitationally bound to the respective fragments, independent of its actual physical state.

Water transfer and loss
Due to the high dimensionality of the parameter space, most studies on volatile losses in similar-sized collisions hold one or more parameters constant or almost constant, where especially the mass ratio is often fixed to a single value. While we focus on hit-and-run collisions in this study, our scenarios cover a significant range of all basic parameters, and thus, we can comprehensively investigate their influence on post-collision volatile inventories. The general increase in volatile losses with impact velocity (Fig. 2) is not surprising, but the dependencies on impact angle and mass ratio are more interesting. There is a strong trend toward losses decreasing with increasing impact angle. While for head-on impacts water losses are considerable even for low v/v esc , and quickly rise above 50% for higher velocities, this figure changes drastically for more oblique collisions and for α 60 • water losses are small (< 10%) even for high-velocity collisions. For α = 45 • , they never rise above 50% for velocities up to 5 × v esc . These results for overall water loss are in good agreement with previous ones by Maindl et al. ( , 2017 who found also a significant decrease in water losses toward oblique impact angles, but studied only fixed mass ratios. There is also broad agreement to earlier results from Canup and Pierazzo (2006), when their definition of water loss as everything that is not bound to the largest fragment is taken into account. However, their treatment ignores that the impactor in a hit-and-run collision also escapes the target's gravity more or less intact and can still contain a large volatile reservoir which could be delivered to the same or nearby objects later. Reufer et al. (2013) focus on mantle stripping of the second largest fragment and consider variations in the four probably most important parameters (v/v esc , α, γ, M tot ) to some degree. Our results for mass ratios between 1:50 and 1:2 show a clear increase in combined water loss for more equally sized bodies, which seems to be more pronounced for lower impact angles (Fig. 2). These large differences can be explained by the fact that the specific energy of a collision Q R (see Sect. 4.5) is also a function 5 of γ and is roughly doubled when going from γ = 1:9 to 1:2 (cf. Fig. 7). In addition geometry comes into play, because the smaller the projectile, the more it is affected relative to the target, but small projectiles have only small water inventories to loose (for an initially equal wmf), while the target retains much of its volatile inventory. For collisions in the hit-and-run regime however, not only the combined volatile losses, but rather individual ones, as well as transfer between the colliding bodies are of interest. Figure 3 illustrates outcomes for the M tot = 10 23 kg scenarios and shows both bodies before (light gray circles and impact geometry) and after the collision (colored circles). Most interesting in this plot are the big differences between the larger and the smaller one of the colliding pair. The target barely looses neither mass nor large amounts of water, except for low impact angles and high velocities, and even then only if the projectile has almost the same size. Disruption of a larger object by the impact of a smaller one was found to be very difficult once bodies grow large enough (e.g., Asphaug 2010) and indeed in our scenario set target disruption happens only in high-velocity head-on collisions with large impactors. The projectile on the other hand is typically much more affected. Except for slow and oblique scenarios, it looses considerable amounts of mass, and especially of water, often above 50% and in some rather 5 For the simplifying assumption of equal bulk density , the specific collision energy (Eq. 4) can be expressed as which increases with γ (for all other parameters equal). ig. 4 One-body perspective as a function of mass ratio for six different parameter pairs of v/v esc and impact angle α (cf. Fig. 3). The sizes (gray before and color after the collision, ∝ mass 1/3 ) and color-coding indicate what happens to an individual body when hit by a projectile depicted by the light gray circles. The uppermost graph in each panel represents a collision of equally composed bodies (wmf = 0.1), the middle graph shows the outcome if only the target body contained water initially, and the lower one depicts the outcome for a dry target high-velocity impacts up to 90%. In addition, we also did some exploratory simulations with α = 75 • and 90 • , which can already be considered a tidal collision, but water losses in these cases are negligible. Our results also confirm those of earlier studies (e.g., Marcus et al. 2010), namely that the volatile fraction practically never increases in giant collisions between similarly composed objects, but always decreases for both bodies.
To explore the more subtle consequences of hit-and-run, like to what proportions a post-collision fragment's water inventory originated from projectile and target, we find it instructive to leave the target-projectile point of view and rather switch to a one-body perspective. In Fig. 5, which illustrates the effects of varying the main collision parameters, projectile and target water are highlighted in different colors. This visualizes compositional mixing (transfer) of projectile and target water (see also Fig. 6), as well as losses from the individual pre-collision inventories. Figure 4 illustrates the consequences for a single body when hit by a range of different impactors, smaller and larger ones. We denote the body of interest always as target and the impacting one as projectile, irrespective of which one is larger, and keep the definition γ = M proj /M targ . The six panels in Fig. 4 correspond to the six (v/v esc , α) parameter pairs enframed in Fig. 3 and exemplify the effects on a single body for common hit-and-run parameters. From the three graphs in each plot, the topmost one shows the post-collision mass and water content if both bodies initially have a wmf of 0.1. Not surprisingly water losses are small for γ < 1 and strongly increase up to around 50% for larger projectiles γ > 1, except for the calmest scenario (lower left panel). The transition from γ < 1 to γ > 1 is smooth, and no rapid increase in water losses-as one might expect-is found. The middle graphs in Fig. 4 represent the same collisions but accounting only for water initially on the target; hence, they represent the impact of an entirely dry projectile onto a 0.1 wmf target. The relatively small differences compared to the topmost graphs in most cases indicate only little transfer of water from projectile to target, even for large γ 1, where the projectile could in principle provide large amounts of water to the target due to its size. This is further emphasized in the bottom graphs, which correspondingly show the outcome if the target were initially dry. Interestingly, the target's wmf (originating solely in projectile water) increases with γ , peaks at relatively small, positive values of γ and decreases again for even larger projectiles (with even larger volatile contents). Except for rather central α = 30 • impacts the highest transferred water content is only around 1/10 of the (much larger) impactor's wmf. In the α = 30 • collisions, especially in the low-velocity v/v esc = 1.5 case, there is considerably more transfer of volatiles, and toward γ = 2:1 transfer from the projectile to the target becomes efficient and the combined wmf seems to even increase again. However, obviously this situation occurs only for the lowest velocity hitand-run encounters, and indeed, this scenario is on the edge to the graze-and-merge regime , indicated by post-collision velocities barely above v esc . For all larger mass ratios and thus lower collision energies, the outcome is not hit-and-run anymore (and therefore not plotted). Hence, it seems that transfer of volatiles between similarly composed bodies in hit-and-run events has a rather minor influence compared to collisional erosion.

Dependence on total mass
While the bulk of this study focuses on collisions of ∼ Moon-sized embryos (M tot = 10 23 kg) we also considered different masses between 10 22 and 10 25 kg, to study the influence of the total mass on the general collision outcome and especially on vaporization of volatile material (see Sect. 4.5). We already investigated the dependence of similar-sized collision outcomes on total mass in Burger and Schäfer (2017) in detail and will therefore limit  Simulation snapshots (cut views) to illustrate the mechanism of volatile transfer and loss. Water originally on the target/projectile in blue/white, and red and black particles represent basalt and iron, respectively. The left sequence shows a scenario with v/v esc = 2.5, α = 45 • and γ = 1:9, and the right image represents the same parameters except for γ = 1:2 (the same scenarios as in the two lower panels in Fig. 5). See the text for discussion ourselves to some important remarks in this paper. The results presented in Fig. 7 indicate that for the larger body the outcome is very similar over the whole investigated range of masses-their volatile inventory remains largely untouched (upper panels). The situation for the smaller body is substantially different (lower panels). For the three simulated mass ratios final wmf are only in the γ = 1:2 case relatively independent of total mass, but show large variations for γ = 1:9, ranging from 0.75 down to only 0.27, and for γ = 1:20 even from 0.79 down to basically zero. The main reasons for these differences are the transition from subsonic to supersonic collision velocities with increasing mass (for other parameters equal), and increasing gravitational compression toward higher masses, resulting in more compact objects with greater hydrostatic pressures to be partly released upon impact. This clearly shows that the total mass is a crucial parameter for stripping of volatiles from the smaller body in a hit-and-run encounter, and at least this aspect of similar-sized collisions can certainly not be assumed to be scale invariant (Asphaug 2010;Burger and Schäfer 2017).

Dependence on water amount and distribution
In order to check the influence of the chosen composition model, we performed additional simulation runs with entirely dry target bodies (with v/v esc = 2.5, γ = 1:9, M tot = 10 23 kg and α = 45 • and 60 • ; cf. Table 1) and compared them to runs with otherwise equal parameters and our usual composition model with both bodies covered in water (ice). Naturally the overall post-collision water contents are lower now, but the amount transferred from the projectile to the target and the amount lost from the projectile are expected to be approximately equal and independent of the target's precise composition. Our results confirm this with differences in wmf (and also in fragment masses and kinetics) around 20% or lower. A dry target consistently accretes more water (from the projectile) than a water-rich one and causes the projectile to loose a larger fraction of its initial volatile content. The projectile also looses more mass in general when colliding with a dry target instead of a water-rich one, where the difference is larger for the α = 45 • case than for 60 • . We suspect that this behavior is an interplay between the higher density of basalt compared to water, which results in a larger resistance against impacting material, and a slightly different collision geometry, because a dry target is more compact than a volatile-rich one. The latter means not least a higher specific impact energy in the dry-target runs (e.g., in the α = 45 • case about 5%, 1.27 vs. 1.33 MJ/kg), since all other parameters were kept constant. This is consistent with the higher losses in the dry-target run. A point closely related to the above is the influence of the extent of the water layer. We repeated the v/v esc = 2.5, α = 45 • , γ = 1:2, M tot = 10 23 kg scenario with wmf of 5 and 20%, in addition to the standard 10% run (see Table 1) and found only relatively small deviations from the expectation that relative water losses and transfers are similar. This means, for instance, that the post-collision water content of the target is roughly twice as large for the 10% water scenario than for the 5% one and again approximately doubled for the 20% run. Considering all three scenarios, deviations are around 25% at most, but if only the 5% and 10% water runs are compared the differences reduce to 10% and less, suggesting a probable convergence toward increasingly thinner volatile layers. The general trend and the expected underlying reasons are the same as in the dry-target case above. Bodies with thinner water ice shell accrete more (from the other body), but also loose more water from their initial inventory, probably again due to slightly different specific impact energies (2.96 vs. 2.86 vs. 2.70 MJ/kg). Maindl et al. (2017) also studied the influence of different initial volatile contents and also found only small variations (∼ 10%) in relative water losses, in concordance with our results. These authors additionally considered water SPH particles randomly distributed inside the whole projectile body. We did not include such uniform distributions because there is strong evidence that bodies in the mass range considered here are probably largely differentiated. It is also yet unclear how such water inclusions (e.g., in hydrated minerals) can be accurately modeled, and how related and perhaps important effects like degassing due to pressure unloading can be included in a consistent way (but see, e.g., Asphaug et al. 2006). However, at least for collisions in the hit-and-run regime Maindl et al. (2017) found no large differences between randomly distributed water and spherical shell configurations in terms of water losses. The generally low dependence of final (relative) wmf on the extent and distribution of volatile inventories is an important result, which helps to safely reduce the dimensionality of the parameter space in simulations on water delivery. We discuss this further in Sect. 5.

Influence of material strength
The justification for usually modeling sufficiently large bodies as strengthless fluids is motivated by highly dominating gravitational stresses over material strength. Jutzi (2015) studied collisions between similar-sized bodies and found that including a realistic rheology is still necessary for 100 km bodies, and recent results using the same SPH code and material model  have shown significant differences in fragment characteristics and also in water losses between solid-body simulations and (otherwise identical) purely hydrodynamic runs also for embryo-sized bodies. For a more comprehensive overview, we refer to the latter study and discuss the topic only briefly here, exemplified by two scenarios (see Table 1), which were computed again with the full solid-body model as outlined in Sect. 2. The results are generally in good agreement with their strengthless counterparts and differ by less than 20% for final wmf and losses, which is also consistent with the results of Burger and Schäfer (2017). Larger deviations of up to 50% are found only for the contribution of transferred volatiles (between projectile and target) to the final inventory. However, this is a more subtle process, and only a minor contribution to a fragment's post-collision water budget compared to the retained amount (cf. Fig. 4); therefore, these relatively large variations are not surprising. These contributions by transferred water are always smaller in the solid scenarios than in the respective strengthless runs, probably because tensile and shear strength act against removal and subsequent transfer. It is yet unclear to what extent and for which aspects of giant collisions material strength becomes important. The large pressures in the interiors of sufficiently massive bodies change the behavior of geologic materials toward high viscosities and ductile, plastic flows (Holsapple 2009), and in general their rheology exhibits a complex dependence on stress, strain, strain rate, temperature and pressure. In addition, even fully damaged rubble pile material is often not strengthless, but can support considerable shear stresses. This makes the here-used von Mises yield criterion not an ideal choice, as already indicated in Sect. 2, since it does not consider pressure-dependent shear strength. A more sophisticated model for the calculation of shear strength of geologic materials has been developed, e.g., by Collins et al. (2004), who use a pressure-dependent yield strength for intact rock and a Coulomb dry-friction law for completely fragmented material. We suggest and plan a dedicated future study to clarify the influence of different rheology models also beyond sizes of the 100 km studied by Jutzi (2015), including also volatile mate-rial like water ice. In addition to these complexities, the physical state of volatile inventories prior to the collision is also uncertain, where, for instance, a considerable fraction might be molten, which would put the application of material strength into perspective. The behavior of real objects probably lies somewhere between that of a strengthless fluid and a fully solid body, and thus, these two models may be considered rather limiting cases, where further studies are necessary to clarify the necessity for certain material models.

Water vapor production
Once impact energies rise above the vaporization energy of water (ice), large-scale vapor production can be expected. Here we treat only vaporization caused directly by the impact and the further development of heat-redistribution, like melting or vaporization of ice by an impact-heated mantle, is not included. To exemplarily estimate the post-collision water vapor fraction, we ran scenarios with fixed v/v esc = 2.5 and α = 45 • for several different total masses (and thus impact energies) between ∼ M Moon /10 and M Earth and for varying mass ratios of 1:2, 1:9 and 1:20 (see Table 1).
The results are summarized in Fig. 7, where the reduced-mass specific impact energy Q R of the individual scenarios is also indicated. It is given by (Stewart and Leinhardt 2009) with the reduced mass μ = M proj M targ /M tot and collision velocity v 0 (see Fig. 1). Not surprisingly the water vapor fraction strongly increases with M tot , and larger amounts of water vapor are only produced for energies greater than e vap . At this threshold, around 10% of the water inventories are already vaporized on both large fragments, and roughly equal throughout the three investigated mass ratios (note the different y-axis scales). If considered as a function of M tot the water vapor content is a function of γ as well, because Q R is also a function of γ (for all other parameters equal), as elaborated in Sect. 4.1 (see footnote 5), and there is certainly a strong correlation between the available energy and the amount of vaporization. The results in Fig. 7 indicate that water vapor production seems to scale indeed well with Q R , relatively independent of the mass ratio, well visible when comparing vapor fractions for specific Q R values (like indicated by the vertical lines). Also the individual vapor fractions on the two largest fragments are surprisingly similar, also for smaller mass ratios. From the remaining water on the two large fragments after a hit-and-run roughly up to 40% can be vaporized for Moon to Mars-sized embryos, and up to 80% in collisions involving Earth-mass objects. These numbers, however, depend on the mass ratio (for a fixed total mass) and may also deviate significantly for other values of impact velocity and angle. Note that due to the limitations of the Tillotson eos the computed vapor fractions are rather rough estimates, as a first step indicating possible directions for future work.

Dependence on resolution
The resolution in SPH simulations is defined by the particle number. Genda et al. (2015) investigated the influence of resolution on the critical specific impact energy for disruption (where the largest post-collision fragment has half the total mass) for gravity-dominated bodies and found a factor of two difference between 50k and 5 million SPH particles, due to different efficiencies of kinetic energy dissipation. They conclude that approximate convergence is reached only for their highest resolutions. Especially low-density regions, like impact ejecta, are affected by differing particle numbers, where the evolution in such regions   Table 1). a, b The wmf of the two largest fragments, c the transferred wmf originating only from the other body (see Sects. 4.1 and 4.3), and d water vapor fractions (see Sect. 4.5) of the two largest fragments after the collision is often dominated by resolution instead of physics (Reinhardt and Stadel 2017;Genda et al. 2015). Another closely connected issue is correctly resolving shock waves, where about 10 particles are required in one dimension (e.g., Genda et al. 2015). While our standard resolution (100k) is sufficient for projectile and target diameters being always clearly above this threshold, the thickness of the water envelope alone is usually below it. For example, a single body made of n part = 10 5 particles and a wmf of 0.1 has a diameter ∝ n 1/3 part of ∼ 57 particles, but a water shell thickness of only ∼ 4 particles, while the latter figure is well above 10 for our highest resolution of 2.25 million.
To check the reliability of our results, we ran three selected scenarios (see Table 1) also with increased resolutions of 300k, 1 million and one with 2.25 million particles (in addition to the standard 100k runs). The outcome for several main quantities of interest as a function of particle number is plotted in Fig. 8 for a collision between an ∼ Earth-sized and a ∼ Marssized object, indicating approximate resolution convergence for most quantities. The large masses of the involved bodies (M tot = 10 25 kg) result in a very energetic collision and thus in large-scale water vapor production (Sect. 4.5), ideal for clarifying the dependence on resolution for the purposes of our topic. In all three scenarios, the global outcomes (fragment masses and their kinetics) of runs with different resolutions show only minor deviations below 10%. Post-collision wmf of the largest fragment are even more accurate within only few %, and for the second largest fragment and the wmf transferred between projectile and target results are still within ∼ 25% (cf. Fig. 8). Note that the outcomes of the 100k runs seem to give an upper limit for the second largest fragment's wmf, meaning that resolution converged volatile stripping of the smaller of the colliding bodies is even more efficient than indicated by the 100k results. The situation for the water vapor fraction is twofold, with deviations of less than 10% for the largest fragment, and ∼ 20% for the second largest one. A more detailed analysis of the latter showed that these values take particularly long to converge after the actual collision itself, since ongoing accretion of debris continues to dissipate energy, leading to further vaporization even after the two main fragments have clearly separated and the material bridge between them has largely dissolved (cf. Fig. 5). For the 2.25 million particles run, this quantity has still not converged at the end of the simulation time and is therefore not included in Fig. 8. We conclude that simulations with probably even higher resolutions would have to be ran for particularly long (simulated) times to finally clarify this issue.
The deviations for the different quantities described above can also be considered a rough error estimate of the presented results w.r.t. resolution convergence. The results confirm that shock acceleration and energy dissipation are already fairly well resolved in the 100k runs, also for the relatively thin (in terms of particle layers) water shells at this resolution. Albeit true resolution convergence is hard to proof these tests indicate that our results are in general within ∼ 25% of this limit.

Comparison with collision outcome models
The most basic collision outcome model is perfect merging with conservation of linear momentum and the assumption of 100% retention of volatiles. This is certainly not realistic in most situations. Leinhardt and Stewart (2012) and Stewart and Leinhardt (2012) developed a comprehensive analytical model that distinguishes several collision outcome regimes, including partial accretion, erosion and hit-and-run. It has been included in N-body simulations and used in several studies on planet formation (e.g., Dwyer et al. 2015;Bonsor et al. 2015;Leinhardt et al. 2015). It is based on scaling laws to determine the mass of the largest collisional fragment M lf and also applies corrections for varying mass ratios as well as the reduced interacting mass in oblique collisions. To also track basic changes in bulk composition, they suggest to include a simple model for mantle stripping introduced by Marcus et al. (2010), which considers two idealized cases for the core mass: (1) M core = min(M lf , M core,proj + M core,targ ), i.e., mantle is only added to the largest fragment once both cores are used up, and (2) either M core = M core,targ + min(M core,proj , M lf − M targ ) on accretion (i.e., if M lf > M targ ), or M core = min(M core,targ , M lf ) for target erosion (i.e., if M lf < M targ ). This means model (2) first adds projectile core material to the whole target in case of accretion, but removes material from the target (starting with the mantle) for erosion. To compute the actual change in composition, it is suggested to use the average of (1) and (2). For the second large survivor in hit-and-run events, they suggest to consider the reverse impact, where the projectile is hit by a hypothetical body consisting only of the part of the target that geometrically overlaps with the projectile, and to apply the above model to this (reversed) collision situation to determine compositional changes of the second largest fragment. Even though this model was not directly developed for volatile transfer and loss in hit-and-run collisions, but rather for the less subtle effects of major bulk compositional changes, we compare it to our numerical results. It turned out that it gives only a very crude estimation of the changes in water contents for scenarios like ours. While predictions for the target body are mostly at least in broad agreement with numerical results, those for the projectile are often very far off. Therefore, we have not included these predictions in the onebody perspective results (Fig. 4), but only in the combined view in Fig. 2, as disconnected and smaller but otherwise equal symbols for comparison. The predicted combined water losses are mostly much larger than the numerical results, mainly because the reverse impact is predicted much more destructive than it actually is, leaving the second largest fragment often entirely water-stripped in this model. We suspect that this is at least partly due to the geometry of a hit-and-run event, where the interacting fraction of the impacting body still grazes by the impacted one, which is likely to be less destructive than a more central collision with the same impacting mass. In addition, this model is based on absolute masses of the cores (which is everything except the water layer for our purposes) and the (predicted) largest fragment; thus, uncertainties increase the lower the amount of the material considered for stripping is. Also, taking a closer look at the components (1) and (2) and the analysis and impression of hit-and-run collisions (cf. Fig. 6), it seems that outcomes are typically rather close or even beyond predictions of component (2) instead of an average of (1) and (2), for both accretion of projectile material by the target and also for target erosion. Note also that the model can not be reasonably used for volatile transfer, for instance from a water-rich projectile to a dry target.

Volatile transfer and loss
Our results show that volatile loss and transfer in hit-and-run collisions is a function of four main parameters, the impact velocity, angle, the colliding bodies' mass ratio, and also of the total colliding mass. The results for combined water loss in Fig. 2 clearly illustrate that losses increase with impact velocity and mass ratio and decrease with impact angle. We find significant values in most parts of the parameter space, except for the most grazing impacts (α 60 • ), lowest mass ratios ( 1:50) and impact velocities close to v esc . In a common protoplanet encounter (e.g.,  with α = 45 • , γ = 1:9 and v/v esc between 1.5 and 2.5, the water mass fraction (wmf) of the target decreases only by less than ∼ 10%, while the projectile looses between 5% and up to 40-50% of its initial wmf in this velocity range (Fig. 4). However, it is important to keep in mind that this figure can change over a wide range (up as well as down) for sufficiently different scenarios. Figure 7 shows that losses also increase strongly with the colliding bodies' mass, which is mainly due to strong shocks once impact speeds become supersonic. For our scenarios, the speed of sound is around 3 km/s. For v/v esc = 2.5, impact velocities approach this value for M tot between 10 22 and 10 23 kg. Our results for overall water loss are in good agreement with previous studies. The reason why head-on and also low-obliquity hit-and-run impacts strip the most volatiles is particularly associated with the large interacting mass in such events compared to more grazing encounters. The cut views in Fig. 6 visualize this for a γ = 1:9 and 1:2 scenario. The smaller the projectile compared to the target the more of it is directly affected and mechanically disrupted due to enormous shear stresses. Additionally shocks and gravitational (tidal) stresses enhance disruption. For γ = 1:9 in Fig. 6 (left), the projectile is ripped apart down to its core, while for γ = 1:2 and the same α (right) both bodies stay more or less intact. During the main interaction phase, the projectile plows through the target (and the target through the projectile) and both shear away the outer layers of their opponent and accelerate most of it away, while also dragging some material with them. This is how the bulk of the volatiles accreted from the other body are transferred. It is well visible from Fig. 6 that the water shell is rather unaffected where no direct mechanical interaction with the impactor took place (even though strong shocks can blow parts of this off too). It is due to this behavior that approximately equal fractions of the overall volatile inventory are lost or transferred even for different extents of the water shell (cf. Sect. 4.3). Marcus et al. (2010) found that volatile fractions practically never increase in giant collisions of similarly composed Earth to Super-Earth-mass bodies, which we can confirm for all our scenarios, and this also holds for the two large fragments in a hit-and-run individually. Note that this does not hold anymore if pre-collision wmf are sufficiently different. The sim-ulations with dry targets (Sect. 4.3) are an example thereof. While embryo-sized and larger objects become increasingly difficult to disrupt by smaller impactors (Asphaug et al. 2006), which we can confirm also for their volatile inventories, hit-and-run collisions are highly transformative for the smaller one of the colliding pair (e.g., Asphaug 2010). We find that this holds especially also for stripping of outer volatile layers, like the water shells in our scenarios. The inefficiency of transfer of projectile volatiles to the target as long as pre-collision mass fractions are similar (illustrated in Fig. 4), means that hit-and-run post-collision inventories in these cases are dominated by the bodies' pre-collision inventories reduced by impact losses. These impact losses are an increasing function of mass ratio for the whole range of simulated values between 1:50 and 50:1. Interestingly, the transferred wmf on the contrary exhibits a maximum value for a given body but differently massive projectiles and decreases again for even larger impactors despite their greater (absolute) volatile inventory. This peak is reached for projectiles a couple of times more massive than the target (Fig. 4), and even in these cases the target can accrete a wmf of only up to 10% of the projectile's pre-collision value, except for slow, low-obliquity encounters. However, transferred volatiles can still significantly modify isotopic fingerprints, while deeper, more protected regions are not or only little affected.

Vaporization of volatile inventories
Losses due to escape of collisionally vaporized material have not been included so far in studies on volatile loss in giant collisions. Our results indicate that the produced vapor fraction is similar on both large fragments, even for smaller mass ratios (Fig. 7). To assess how much of this vaporized material eventually escapes, remains in the object's atmosphere or perhaps recondenses, a lot of additional factors come into play, among them the mass of body and atmosphere, their thermal state (e.g., further vaporization by an impact-heated mantle), and external factors as the distance to the host star, the stellar EUV flux and the frequency of impacts. Recent modeling by Odert et al. (2017) that takes most of these factors into account suggests that embryos in the terrestrial planet region up to Mars size quickly loose their water vapor atmosphere (∼ Myrs), either thermally or driven by the large EUV fluxes of young stars. Due to the large number of unknowns in this figure and the only approximate treatment of vaporization with the Tillotson eos, we cannot make quantitative conclusions on these additional volatile losses. However, our results ( Fig. 7) indicate that vaporization losses are probably most important in our intermediate mass scenarios around M tot = 10 24 kg (∼ Mars-sized bodies), where vapor fractions are already considerable, but the (hit-and-run) fragments' masses, particularly of the smaller body, are still low enough for large-scale escape. According to our estimates, fragments in this mass range could loose another ∼ 10 to 40% of their remaining volatile inventory. Therefore, it seems that collisional volatile stripping is dominated by direct losses due to mechanical/gravitational stresses and shock acceleration, but vaporization and subsequent loss can still enhance this further, not as the dominating contribution, but also not negligibly. Continuing from these first steps, a dedicated future study including an improved treatment of phase changes by a thermodynamically consistent eos like ANEOS (see, e.g., Melosh 2007) and higher resolution to resolve also low-density plumes in detail (Reinhardt and Stadel 2017) could provide improved estimates on this topic.

Water delivery to terrestrial planets
What does this mean for water (volatile) delivery in early planetary systems? With a focus on the target planet for water delivery the prime example is proto-Earth, with strong isotopic evidence that the bulk of its water inventory originated from the outer main belt (or the same source as this material). In a scenario where the majority of these volatiles is delivered by only a relatively small number of large impactors (Morbidelli et al. 2000), hit-and-run can play a decisive role, because it is a frequent outcome in such similar-sized collisions (∼ 50%; Kokubo and Genda 2010;Genda et al. 2017). The chances for growing planets like proto-Earth to be hit by a larger body are slim; therefore, the majority of impactors are likely smaller, and the generally low transfer efficiency, combined with additional losses as from vaporization, strongly limits the amount of accreted volatiles. This implies a substantial difference compared to the still mostly used assumption of perfect merging of water inventories in (N-body) planet formation simulations. It depends on the dynamical environment and resulting impact velocities, how efficient volatiles can be accreted from a hit-and-run impactor, but in such a scenario they might be high since the impactors were scattered over large orbital distances. For collisions not in the hit-and-run regime, the results for head-on scenarios in Fig. 2 indicate that volatile losses for rather central impacts are also high, but this is only the combined value (i.e., for identical pre-collision wmf of both bodies). A more detailed study on the transfer efficiency in collisions beyond hit-and-run (less oblique, either accretionary or erosive) could help to reduce these uncertainties. The consequences of hitand-run for objects just below the top of the size distribution when colliding with the largest are often large-scale stripping of their volatiles. Our results confirm that this suggests a tendency toward dry Earth-mass planets in extrasolar systems with Super-Earths (e.g., Asphaug 2010).
In a young planetary system that is rather viewed as an ensemble of interacting embryos/protoplanets, the collisional evolution is also strongly driven by the distribution of impact velocities. Even for modest values of v/v esc combined water losses in hit-and-run events can go up to ∼ 40% (Fig. 2). It may not be unlikely that a specific embryo suffers multiple hit-and-run events, often as the smaller opponent, stripping most or all of its volatile content, before it is finally incorporated into a growing planet. If bodies experience a series of such encounters, either as the larger or the smaller of the colliding pair, combined with additional losses of vaporized material and other environmental effects, it is probably a rather conservative estimate to place final water contents after the cumulative effects of a tens-of-Myrs giant collision phase at ∼ 1/10 of what perfect merging simulations suggest. A combination of these loss processes is a natural explanation for the too high water contents that typically result from perfect merging simulations. Currently available collision outcome models cannot reliably predict volatile losses and changes in wmf, as shown by a comparison with our numerical results (Sect. 4.7 and Fig. 2), particularly not for the second largest hit-and-run fragments. These models were developed rather for large-scale changes in bulk composition and not for the more subtle effects concerning (surface) volatile inventories. Hybrid codes that simulate both, the long-term dynamical evolution as well as the hydrodynamics of individual collisions, are computationally very expensive and can thus be ran only in low resolution currently, as recently done by Genda et al. (2017) with 10k SPH particles per collision. This is probably a too low resolution for studying the fate of surface volatiles. Realistic modeling of volatile delivery in planet formation simulations must begin with a profound understanding of transfer and loss mechanisms in individual collisions, to eventually include models developed for this particular aspect of collision outcomes into terrestrial planet formation computations. The high dimensionality of the collision parameter space makes the development of these models-e.g., in the form of scaling laws-a formidable task, which we plan to address as the next logical step in a future study. The similarity in (relative) transfer and loss results (Sect. 4.3), independent of the actual extent of the volatile shell (on either body), might proof very helpful for the development of such a model, because it likely eliminates the necessity to include absolute volatile amounts and their distribution as additional parameters.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.

Appendix A: Semi-analytical relaxation of initial conditions
Simulations of giant collisions during the late stages of planet formation comprise large, selfgravitating bodies, which naturally exhibit a certain internal structure representing hydrostatic equilibrium. The initial conditions for such computations should be in a relaxed configuration (in terms of the applied numerical model) to resemble reality. The common approach is to apply some dynamical settling (numerical relaxation) prior to the actual simulation run, which, however, results in considerable amount of additional computing time. To overcome the drawbacks of numerical relaxation, we use a semi-analytical approach instead, by computing hydrostatic equilibrium structures for multi-material spherical bodies, before setting up the initial particle configuration according to them. This is done in a self-consistent manner, i.e., consistent with the physical model implemented in the simulations themselves, to ensure a configuration already very close to equilibrium. Figure 9a illustrates this by comparing radial surface speeds for a numerical relaxation run (i.e., starting from homogeneous-density bodies) and a run comprising a body of equal mass and composition-but set up following our semi-analytical method. While numerical relaxation naturally leads to strong radial oscillations (due to initial gravitational collapse and subsequent rebound), the results of the semi-analytical approach show only some minor residual fluctuations, mainly due to material boundary effects, inherent to SPH (see also Reinhardt and Stadel 2017;Woolfson 2007). Maximum particle velocities falling below some threshold (relative to the collision velocity) are an often considered criteria for a sufficiently relaxed body (e.g., Hosono et al. 2016, and references therein). In the example in Fig. 9, the peak velocity at the onset of numerical relaxation is above 1 km/s, where the surface escape velocity is about 4.8 km/s, while the maximum value after semi-analytical relaxation is around 50 m/s and drops quickly, mostly still during the initial approach phase of projectile and target. This is only about 1% of v esc (and therefore of a typical collision velocity), where this value is also an often used criterion for a sufficiently relaxed configuration (as, e.g., in Hosono et al. 2016;Reufer et al. 2012). The validity of the semi-analytically relaxed initial conditions is also evident from Fig. 9b, which illustrates that the calculated profiles are-except for the typical boundary effects-very close to equilibrium configurations. The semi-analytical approach works almost instantaneously, and practical experience has shown no significant differences to numerically relaxed (but otherwise equal) runs. Thus, it is a convenient alternative for producing equilibrated initial conditions in giant collision simulations. The algorithm's implementation in C is freely available to the community upon request 6 . In the following, we briefly outline the calculation of the hydrostatic structure and the adiabatic internal energy profile.

A.1. Hydrostatic structure calculation
For the given geometry, i.e., spherical symmetry, the continuum mechanics problem reduces to simple hydrostatics, even if full solid-body physics is considered in principle. The Euler representation of the set of equations necessary for describing a spherical, hydrostatic configuration is with the pressure p, density and integrated (enclosed) mass m, all functions of the radial coordinate r , and the gravitational constant G. An eos p = p ( , e), along with a convenient treatment of the (specific) internal energy e, like, e.g., e = e ( ), close the system of equations for p(r ), (r ), m(r ) and e(r ). The treatment of e is basically open to many different approaches; here, its final value is calculated for a given state of adiabatic compression, determined by density alone, neglecting other non-mechanical contributions (see below). Even though it is possible to directly use these equations, transforming them to their Lagrangian form (with m as the independent variable instead of r ) is very advantageous when it comes to numerical stability and handiness. This is due to the already initially well- Since these are differential equations for r and p, respectively, it is necessary to also calculate for a given p (and e), each time the differential equations' numerical solution is advanced by one step (in m). This can be achieved by inverting the eos, usually given as p = p ( , e), to obtain = (p, e), which is in principle possible as long as the eos is bijective, even though it might not be possible to invert it analytically (as for the Tillotson eos). Inserting e ( ) into the eos results in the self-consistency problem = p, e( ) . We solve this problem by a simple fixed-point iteration procedure, following the steps p = p( , e) −→ = ( p, e) = p, e( ) −→ (n+1) = p, e (n) + (n−1) 2 .
Collecting the pieces together finally allows to calculate (m), p(m) and r (m) from the differential equations (6) and from the iteration (7), and subsequently e(m) = e (m) . How this can be carried out practically can be understood with the aid of Fig. 10, depicting a body with a core-shell structure. It is clear that the numerical integration has to start at the outer boundary and be advanced inwards, since initially the properties at the center are usually unknown, while conditions at the outer boundary are specified by the overall mass M, a given surface pressure (typically p = 0) and the density at the surface. A suitable initial guess for the body's radius R serves as the starting point for the first inwards integration down to the core. In the following, an iterative approach varies R until the found solution (between surface and center) converges toward a consistent internal structure, meaning in particular r (m = 0) ∼ = 0 (where the sign of r (m = 0) can be used to test the current R and adjust it for the next iteration). We found a simple fourth-order Runge-Kutta integrator sufficient for this task.