The threat of Centaurs for terrestrial planets and their orbital evolution as impactors

Abstract

Centaurs are solar system objects with orbits situated among those of the giant planets. Originating mainly from the Trans-Neptunian objects, Centaurs represent one of the sources of Near-Earth Objects. Thus, it is crucial to understand their orbital evolution which in some cases might end in collision with terrestrial planets and produce catastrophic events. We study the orbital evolution of the Centaurs toward the inner solar system and estimate the number of close encounters and impacts with the terrestrial planets after the Late Heavy Bombardment, assuming a steady-state population of Centaurs. We also estimate the possible crater sizes. We compute the approximate amount of water released on the Earth, which is about 10⁻⁵ of the total water present now. We also found subregions of the Centaurs where the possible impactors originate from. While crater sizes could extend up to hundreds of kilometres in diameter, given the presently known population of Centaurs the majority of the craters would be less than ∼10 km in all the planets. For an average impactor size of ∼12 km in diameter, the average impact frequency since the Late Heavy Bombardment is one every ∼1.9 Gyr for the Earth and 2.1 Gyr for Venus. For smaller bodies (e.g. >1 km), the impact frequency is one every 14.4 Myr for the Earth, 13.1 Myr for Venus, and 46.3 for Mars in the recent solar system. Only 53 per cent of the Centaurs can enter into the terrestrial planet region and ∼7 per cent can interact with terrestrial planets.

celestial mechanics, minor planets, asteroids: general

1 INTRODUCTION

Centaurs are objects in the solar system, with their orbits situated between those of the giant planets. There is the lack of general agreement among various definitions; thus, in this paper we will consider Centaurs as objects with a semimajor axis ranging from 5.5 to 30 au.

They are generally thought as a by-product or an evolutionary product of the Trans-Neptunian objects, and thus should have generally lower density than the main belt bodies. This is because objects coming into the inner solar system from beyond Neptune can become comets more often than asteroids in the main belt and as such should be consistent with volatile material. Centaurs have a short lifetime compared to the main belt asteroids (MBAs), and are typically removed from the solar system on time-scales of only a few millions of years (Di Sisto, Brunini & de Elía 2010; Galiazzo, Wiegert & Aljbaae 2016). These bodies are one of the main sources of Near-Earth Objects (NEOs). Thus, it is crucial to understand their orbital evolution, since it might end in a collision with one of the terrestrial planets (TPs) and lead to catastrophic events.

Some former comets that had performed many passages in the inner solar system could also have exhausted part of their water content and thus can be considered as hydrated asteroids, very likely of type C, D, or P (Barucci, Fulchignoni & Lazzarin 1996). All this material can be delivered to the planets via impacts and become important for their atmospheres and have astrobiological implications. This work aims to (i) constrain the rate of impacts on the TPs for the entire Centaur population in the recent solar system (after 3.8 Gyr ago; Hartmann et al. 2012; Kirchoff et al. 2013), meaning after the Late Heavy Bombardment; (ii) find the typical orbits of the Centaur impactors; and (iii) determine from which sub-region of the Centaurs the impactors come from. Some past works have already attempted to compute this rate, but none of them considered the entire Centaurs population, e.g. Napier (2015) and Napier et al. (2015). The paper is structured as follows. Section 2 describes the population sample, the tools to analyse the sample, and how the analysis was made. The results are given in Section 3: the source regions of the impactors are described in Section 3.1. Sections 3.2 and 3.3 show the statistical analysis and the results on close encounters (CEs) and the impact probability, interactions with the TPs, as well as a brief overview of interactions with the giant planets. In addition, we estimate the contribution of the Centaurs to the NEOs, potential hazardous objects (PHOs), and short period comets (SPCs). The conclusions are given in Section 4.

2 MODEL

We forward integrate 319 Centaurs (data¹ from JPL Small-Body Database Search Engine: http://ssd.jpl.nasa.gov/sbdb query.cgi, see Table 1 for their osculating elements) for 30 Myr, which is a span of time larger than the average lifetime of Centaurs (typically removed from the solar system on time-scales of only a few million years; Dones, Levison & Duncan 1996; Levison & Duncan 1997; Tiscareno & Malhotra 2003; Horner, Evans & Bailey 2004a; Di Sisto & Brunini 2007; Bailey & Malhotra 2009; Galiazzo 2013), or until the body has undergone an impact (see Table 2), or reaches a hyperbolic orbit and escapes. For each Centaur, we compute 14 clones (the total of 4785 orbital evolutions, considering real initial orbits plus each clone orbits), distributed as described in Horner et al. (2004a) and Galiazzo et al. (2016) where the osculating elements are a = a₀ ± 0.005, e = e₀ ± 0.003, and i = i₀ ± 0.01, here a₀, e₀, and i₀ are the initial elements. For the Sun, we consider an impact when the body reaches a distance² of 0.00465 au relative to its centre. Some comets were observed to survive in relatively close proximity to the Sun, such as C/2011 W3 (Lovejoy) that survived at a distance of 0.0055 au (Sekanina & Chodas 2012).

We use the Lie-integrator (Hanslmeier & Dvorak 1984), an orbital numerical integrator with an adaptive step-size (Eggl & Dvorak 2010; Bancelin, Hestroffer & Thuillot 2012), also utilized in several previous works dealing with CEs of asteroids with the planets (Galiazzo 2013; Galiazzo, Bazso & Dvorak 2014; Galiazzo & Schwarz 2014). In this work, the accuracy parameter of the Lie-integrator is set to 10⁻¹³. The output stepsize for the numerical integration is 2 kyr, and a simplified solar system (SSS) is considered for the orbital propagation. In the SSS, all planets from Venus to Neptune are considered, with the mass of Mercury added to the mass of the Sun and assuming the Earth and the Moon are one body in their barycentric point³; physical destruction and non-gravitational effects in the orbital propagation were neglected. We assume a CE occurs when a Centaur is within a distance from the centre of the planet equal to r_CE = 0.01 au from the perturbing body. This is because we aim to establish a region that can adequately cover the Hill sphere of the largest TP (the Earth). However, for statistical analysis of the CEs we separate TPs from giant planets (GPs). For TPs, we consider the Hill radius of the Hill sphere as r_CE⁴ and, for GPs we scale to the largest Hill sphere, the one of Neptune, taking into consideration the maximum value we have set for CEs, 0.01 au (see Table 2).

We consider the population comprised of observed Centaurs for the statistical computations and the probability results will pertain to a Centaur with the average size in the sample considered. Then, to mitigate bias, we only consider the unbiased (computed through a Monte Carlo simulation on the size distribution of the observed Centaurs, see Sheppard et al. 2000) population comprised of known Centaurs with D > 1 km. Based on the results from Sheppard et al. (2000), there are ∼10⁷ Centaurs, and, in the work of Horner et al. (2004a), km-size Centaurs should be 44 300 (which later, see Section 3.1, will be the favourable value).

We calculated the percentage of Centaurs that encounter the planets, |$\rm CE_{Ct,{{\ \rm per\ cent}}}$|⁠. The following was computed: the number of CEs for all the Centaurs CE, the rate of encounters per million years, CE_T (My⁻¹), the average time for CEs with a planet, |$\bar{T}$| (My), the time it takes for one respective body to experience a CE, |$T_{\rm CE_{Ct}}$|⁠, and the time required to undergo a CE per planet, T_CE (My):

\begin{eqnarray} {\rm CE}_{{\rm Ct,{{\ \rm per\ cent}}}}&=&N_{\rm CE}/15/N_{T,\rm CE}\cdot 100,\nonumber\\ {\rm CE}&=&{\rm CE}_{\rm C}/15,\nonumber \\ {\rm CE}_T&=&{\rm CE}/T_{l,R},\nonumber\\ T_{{\rm CE_{Ct}}}&=&{\rm (CE_{Ct,{{\ \rm per\ cent}}}} N_{T,{\rm CE}}/100/T_{l,R})^{-1},\nonumber\\ T_{\rm CE}&=&({\rm CE}/T_{l,R})^{-1}, \end{eqnarray}

(1)

where N_CE, N_T,CE, CE_C, and T_l,R are the total number of clones encountering the planet, the total number of observed Centaurs, the total number of CEs for all the clones, and the average lifetime of a Centaur, respectively. Then, the following quantities were computed: the impact probability (see Section 3.3), and the average, maximum and minimum impact velocity for each TP, in order to estimate the size of the craters, which are related to the impact velocity. The impact velocity (⁠|$v$|_imp; Collins, Melosh & Marcus 2005; Galiazzo, Dvorak, Bazso 2014) was computed from the velocity at infinity (before the atmospheric entry, at the proximity of the Hill sphere border, |$v$|_∞), taking into consideration the density, the diameter of the asteroid, and the impact angle (D_A, ρ_A, and θ_imp), the surface pressure and the gravitational acceleration of the given planet (P_S and g_P):

\begin{eqnarray*} v_{\rm imp}= \frac{\rho _{\rm A} D_{\rm A}}{\frac{\rho _{\rm A} D_{\rm A}+3 P_{\rm S}}{2g_{\rm P} \sin \theta _{\rm imp}}} v_\infty. \end{eqnarray*}

(2)

The mass of the body was obtained by assuming a spherical shape and the average density⁵ for Centaurs with a diameter less than 400 km: (940 kg m⁻³). Taking into account the average impact angle of 45° (Le Feuvre & Wieczorek 2008), we determined the impact velocity using the velocity at infinitum (as described earlier). These data, along with asteroid diameters (Table 3), are needed to compute the minimum, average and maximum impact crater diameters on TPs. The asteroid diameter was obtained from the Johnstone’s archive (where in the case of multiple values, we selected the most recent one). If the diameter was not available from the archive, we computed it via the absolute magnitude and the average albedo (taken from Johnstone’s archive) for Centaurs ρ_V = 0.078 using the standard equation of Tedesco et al. (1992):

\begin{eqnarray*} D=\frac{1329}{\rho _V} 10^{\frac{-H_V}{5}}, \end{eqnarray*}

(3)

where H is the absolute magnitude of the Centaur.

The computations were done for two specific regions associated with Centaurs and for the Centaurs population as a whole. The two Centaur regions were divided based on the relative semimajor axis of the given minor bodies, because found by Galiazzo et al. (2016) there are two distinct average sizes borderlining J4:1 mean motion resonance (at 13.11 au). Beyond this resonance, the largest Centaurs (with the diameter larger than ∼100 km) are present in the outer part of the Centaurs’ region. We define Region 1 (R1) as the region beyond J4:1 (⁠|$13.11 \, \mathrm{au} \lesssim a \le 30 \, \mathrm{au}$|⁠), and Region 2 (R2) as |$5.5 \, \mathrm{au} \le a \lesssim 13.11 \, \mathrm{au}$|⁠. Region 1 contains 160 known Centaurs and Region 2, 159. For the whole population, we summed up the probabilities and obtained the average of the average of the diameter in Regions 1 and 2. Thus, the probability is associated with a Centaur having that average value of diameter.

We compute the amount of water released (W_M, in Earth masses per million years), assuming that the amount of water released by a single active Centaur, as a cometary body, is P_C = 15 per cent [this is consistent with the amount of water on the surface of the nucleus of comet 67P G-C, which is taken as a proxy for the whole body (Fulle et al. 2016) of its own mass], and P_NC = 10 per cent from other non-active Centaurs (asteroids), assuming their content to be similar to the most hydrated meteorite (Morbidelli et al. 2000). Guilberte-Lepoutre (2012) estimated that only ∼11 per cent (⁠|$P_{\rm C_A}$|⁠, percentage of active comets) of the Centaurs would become active (comets).

We used the average diameter (mass) per region and considering the rate of impacts from 0.1 Gyr ago until now and from the Late Heavy Bombardment (LHB) until now, the probability was then multiplied by 1.49. Given that the population of Centaurs (such as their progenitors, the TNOs) from now until the LHB was on average 149 per cent the present one and that in the last 0.1 Gyr it did not significantly change (only ∼1 per cent)⁶ (I_r = 1/K_r Myr, where K_r is the impact rate in million years per each planet), it follows:

\begin{eqnarray} W_{\rm M}&=&W_{\rm M_C}+W_{\rm M_{NC}},\nonumber\\ W_{\rm M_C}&=&\frac{4\pi \rho r^3_{\rm c,region}}{3} I_r (P_{\rm C} P_{\rm C_A}),\nonumber\\ W_{\rm M_{NC}}&=&\frac{4\pi \rho r^3_{\rm c,region}}{3} I_r [P_{\rm NC}(1- P_{\rm C_A})],\nonumber\\ W_{\rm M}&=&\frac{4\pi \rho r^3_{\rm c,region}}{3} I_r [P_{\rm C} P_{\rm C_A}+P_{\rm NC}(1- P_{\rm C_A})]. \end{eqnarray}

(4)

Here |$W_{\rm M_C}$| and |$W_{\rm M_{NC}}$| is the amount of water released per million year, respectively, by active and non-active Centaurs.

Considering that we are interested in establishing (roughly) the amount of water released since the end of the LHB, this result was multiplied by the average number of bodies expected from that time until now. This results in 149 per cent of the present population⁷ (see Fig. 1 for the expected population from the LHB to the present).

Figure 1.

Decay of the Centaurs (from TNOs obtained in the model of Brasser & Morbidelli 2013). This is a cumulative number of surviving minor bodies versus time in million years. The arrow shows the end of the LHB era, when the analysis of the Centaurs was initiated. The data comprising of the population decay of the TNOs were kindly provided by A. Morbidelli.

For the statistical computation of impacts, we also took into consideration the (average) amount of material lost per region and for the whole Centaurs population. Here, we take into account the law proposed by Öpik (1976), and used also by Napier (2015), to determine the disruption (D_kill) in one dimension⁸:

\begin{eqnarray*} D_{\rm kill}=\frac{3.8}{\sqrt{\bar{q}}} \frac{\bar{T}_{\rm life}}{\bar{T}_{\rm Rev}} \, \mathrm{[m]}. \end{eqnarray*}

(5)

The expression above is a function of the perihelion (⁠|$\bar{q}$|⁠) and apparition (perihelion passage). The latter represents the ratio between the average lifetime (⁠|$\bar{T^{\prime }}_{\rm life}$|⁠) and the orbital revolution (⁠|$\bar{T}_{\rm Rev}$|⁠).

We compute D_kill when a given Centaur is supposed to be active (as a comet, or if it can become a comet) – this takes place when a Centaur remains inside the radius |$a_c=a(1-e^2)\lt 12 \, \mathrm{au}$| (see Guilberte-Lepoutre 2012). This limit is appropriate because it is valid for Region 2, beyond Saturn (further away from the usual snowline, roughly |${\sim } 5 \, \mathrm{au}$|⁠, as the model of Cyr, Sears & Lunine 1998 asserts), and it has been determined to occur for real comets: the furthest active comet seen to-date is C/2017 K2 (PANSTARRS), observed at a distance of ∼24 au (Jewitt et al. 2017). We consider the average values (when |$a_c\lt 12 \, \mathrm{au}$|⁠) for all the clones of

• the average perihelion, |$\bar{q}$|⁠;

• the average period of revolution (from the averaged semimajor axis): |$\bar{T}_{\rm Rev} = \bar{a}^{1.5}\,y;$|

• the average lifetime, |$\bar{T^{\prime }}_{\rm life}$|⁠.

The impact crater dimensions can also be estimated using the analytical expression (see Collins, Melosh & Marcus 2005):

\begin{eqnarray*} D_{\rm tc}=1.161 \left(\frac{\rho _{\rm i}}{\rho _{\rm t}} \right)^{\frac{1}{3}} L^{0.78} v^{0.44}_{\rm i} g^{-0.22}_{\rm P} {\rm sin}^{\frac{1}{3}}{\theta }. \end{eqnarray*}

(6)

Here, D_tc is the transient crater diameter, ρ_i and ρ_t are the impactor (projectile) and target density, respectively, L is the impactor diameter, |$v$|_i is the impact velocity, g_P is the gravitational acceleration for a given planet, and θ is the impact angle (here assumed to be 45 deg). The final diameter (D_f) for simple (bowl shaped) craters is given by D_f ∼ 1.25D_tc (equation 22 in Collins et al. 2005). For complex craters, we use the expression that relates the transient crater diameter to the diameter associated with the simple-to-complex transition (D_c) on a planetary body (see equation 27 of Collins et al. 2005):

\begin{eqnarray*} D_{\rm fr}=1.17\frac{D^{1.13}_{\rm tr}}{D^{0.13}_{\rm C}}. \end{eqnarray*}

(7)

For the sake of completeness, we also performed numerical simulations of impact crater formation using iSALE-2D, a multimaterial, multirheology shock physics code (Melosh, Ryan & Asphaug 1992; Ivanov, Deniem & Neukum 1997; Collins, Melosh & Ivanov 2004; Wunnemann, Collins & Melosh 2006), an extension of the SALE hydrocode (Amsden, Ruppen & Hirt 1980). To accommodate crater formation by impacts with a wide range of projectile sizes, in all our simulations we apply a two-layer model, with the crust and the mantle being represented with the equation of state (EOS) for granite (Pierazzo, Vickery & Melosh 1997) and dunite (Benz, Cameron & Melosh 1989), respectively. These are considered analogues for planetary crust and mantle (e.g. Yue et al. 2013). Since the projectile is considered to be of cometary origin, we applied the parameters appropriate for ice (see Bray et al. 2014), including the Tillotson EOS (Tilllotson 1962; Ivanov et al. 2002). Due to axial symmetry, all impacts are vertical, with the impact velocity corresponding to the vertical component of the velocity calculated earlier in this section. The crustal thickness (d_crust) and thermal gradient (dT/d|$z$|⁠) for the three planets are as follows: Venus, d_crust = 25 km, dT/d|$z$| = 10 K km⁻¹ (Grimm & Solomon 1988; James, Zuber & Phillips 2013); Mars d_crust = 50 km, dT/d|$z$| = 8 K km⁻¹ (Zuber et al. 2000; Neumann et al. 2004); and the Earth d_crust = 32 km, dT/d|$z$| = 25 K km⁻¹ (e.g. Christensen & Mooney 1995 and references therein). The grid resolution was varied according to the impact size, with the span of 10–30 cells per projectile radius. Such choice in resolution is adequate for validating impact crater dimensions, while reducing computational expense. To prevent numerical instabilities, the high-resolution zone was set up such that it extends sufficiently far from the crater rim horizontally and below the transient crater floor vertically. The simulations were run as long as it was required for the crater collapse to seize. In simulation time, this could be anywhere from a few hundred seconds up to several hours (basin size impacts). The crater diameter was measured from rim to rim, following the methodology described in Silber & Johnson (2017). Considering that the crater size obtained through simulations is consistent with that obtained using equations (6) and (7), it was sufficient to preform numerical simulations for a representative sample of impacts of varying sizes, in order to capture a dynamic range of impact craters on the TEs.

3 RESULTS

From numerical orbital integrations, we find that Centaur clones undergo CEs with all the planets; however, real impacts take place only with GPs, but not with TPs. Thus, to estimate the impact probability with TPs, we extrapolate it via fitting the number of CEs per relative distance to the planet inside its Hill’s sphere. We also compute the probability for Centaurs with D > 1 km (as mentioned in Section 2, these represent the current population).

3.1 Source regions for terrestrial planet encounters in the Centaurs region

The entire Centaurs region is very chaotic, apart from some constrained orbital sub-regions. The source sites associated with terrestrial planet encounters (TPEs) can be seen in the brightest sections of Fig. 2. The sub-regions that produce more TPEs are all of Region 2 (specifically, in the sub-region with semimajor axis, |$5.5 \, \mathrm{au} \le a \lesssim 8 \, \mathrm{au}$| and with eccentricity, 0.09 ≤ e ≲ 0.75). There is also a minor contribution from Region 1 in the following two phase space sub-regions:

• |$a\gtrsim 20 \, \mathrm{au}$| and e ≲ 0.63,

• |$26 \, \mathrm{au} \le a \lesssim 28 \, \mathrm{au}$| and 0.8 ≲ e ≲ 1.

Figure 2.

Source regions in the a − e and a − i space. The colour code is related to the density of the relative points in each region, Gaussian scaled. Note that the colour bar scale is different in each panel. The full colour figure is available in the online version of this paper.

It is interesting to point out the inclination, where there is initially a predominant region for 8° ≲ i ≲ 40°, subsequently decreasing to |$a\sim 24 \, \mathrm{au}$| at i ∼ 30°. There is a minor contribution, although significant compared to other sub-regions in R1, in the sub-region for retrograde orbits at i ∼ 140° at |$22 \, \mathrm{au}$| (see Fig. 2, bottom panel). However, the most predominant inclination is i ∼ 20°, similar to that of the Jupiter-family comets. Most of Centaurs are the sources of Jupiter family comets (e.g. Fig. 3), where a Centaur becomes a Jupiter-family comet at 24 kyr [its Tisserand parameter in respect to Jupiter (T_J) is between 2 and 3, which by definition, corresponds to Jupiter family comets]. Thus, all the sub-regions are also the major sources of possible impactors. The average lifetime for different regions varies. In some, the scattering from the average value is very large, and some also have very short lifespan (see i.e. Fig. 3), of the order of only thousands of years. We find the following lifetimes for Region 1, Region 2, and the whole region (all the Centaurs), respectively, T_life,1 = 10.158 Myr, T_life,2 = 1.824 Myr, and T_life,_TOT ∼ 6.041 Myr. This value is much higher than that found by Napier (2015); however, they used only Chiron-like orbits and performed orbital evolutions for only 1 Myr (in our case, the average Chiron lifetime is 3.884 Myr) as a representation of the entire Centaurs population in fact Chiron is also in R1, a region where the Centaurs have an average larger life than the ones in R2. Of all the Centaurs population, 6.2 per cent collide with the Sun, 82.5 per cent escape the outer solar system in hyperbolic orbits, about 0.9 per cent have impacts with planets [Napier (2015) computed 86 per cent, in accordance with our result] and the rest survive 30 Myr. Region 2 is the region that experiences more impacts with the Sun, 56.4 per cent of the total number of impacts with the Sun of the all Centaurs. This region also has relatively more impacts with GPs, and it dominates (as source of impactors again) for TPs (the latter one ‘statistically’), the ratio⁹ is 1.6 and 20.7, respectively (see Section 3.3).

Figure 3.

A short-life Centaur (coming from a source region of TPEs) that becomes a Jupiter family comet and then collides with the Sun. Plotted are the semimajor axis (a), perihelion (q), eccentricity (e), sinus of inclination (i), and Tisserand parameter in respect to the time. Horizontal lines at values 2 and 3 of y-axis represented the borders for Tisserand parameters (T_J ≥ 3 for typical asteroidal orbits and 2 ≤ T_J ≤ 3 for Jupiter family comet orbits).

3.2 Close encounters

CEs with Centaurs occur for all the planets, as shown in Table 4, where we can gauge that the Earth undergoes most encounters.

Table 1.

The range of osculating elements for the Centaurs sample considered in this work. The data are divided into regions (Region 1, outer region, and Region 2, inner region). Listed are (in this order) semimajor axis, eccentricity, inclination, perihelion, and aphelion (minimum and maximum for the observed Centaurs in the range defined for the respective region).

Region	a (au)	e	i	q (au)	Q (au)
1	13.15–29.92	0.009–0.942	2.0°–170.7°	1.42–27.38	15.67–56.43
2	5.50–13.02	0.022–0.867	2.0°–175.2°	1.37–12.31	6.37–22.99

Table 2.

CE and impact radii per planet.

Planets	rCE (10−3 au)	rimp
Venus	6.71	0.04045
Eartha	9.88	0.04259
Mars	6.66	0.02266
Jupiter	4.38	0.46733
Saturn	5.33	0.38926
Uranus	5.84	0.16953
Neptune	10.00	0.16459

^aMeans that the ‘Earth’ is in reality the Earth and Moon together in the barycentric position.

We found that about 53 per cent of the Centaurs enter the TP region, and of these, 57 per cent are from Region 2 and the rest are from Region 1.

The Centaurs, in general, arrive to the planets up to a couple million years (Table 4) comparable to that of the MBA becoming NEOs (see table III of Bottke et al. 2002).

This demonstrates that the contribution of the Centaurs to the NEOs population is very important. The Centaurs are a relevant source of NEOs, considering that these objects have typical lifetimes consistent with a mixture of orbital evolutions of (former) MBAs and Centaurs. The latter usually evolve much faster as NEOs, only about 0.2 Myr, and consequently the Centaur-NEOs have shorter lifetimes. The whole population of NEOs survives for T_N ∼ 10 My (Gladman, Michel & Froeschlé 2000). We also found that 43.4 per cent of the Centaurs become NEOs during their lifetime. According to the data available from HORIZON JPL (https://ssd.jpl.nasa.gov/sbdb_query.cgi, 977 NEOs with a diameter larger than 1 km¹⁰ are known (and observed) now. It is thought that about 10 per cent of the population still remains undiscovered (https://cneos.jpl.nasa.gov/stats/ and Tricarico 2017). From these values, we can compute approximately the present contribution of Centaurs to km-size NEOs. We find that the observed population of Centaurs should produce 60 NEOs from Region 1 and 78 NEOs from Region 2, respectively, with the average dynamical time of 96 kyr for both (see Table 5). Because the predicted numbers of NEOs from Bottke et al. (2002) are extrapolated, using a source zone with these borders: a < 4.2 au, e < 1.0, and i < 90° (a zone called ETR, Extended Target Region, from Bottke et al. 2002). We cannot compare our results with the standard debiased population of NEOs (because Centaurs were not considered). Therefore to compute the NEO contribution of Centaurs of kilometre size, we make a rough estimate, considering NEOs in the ETR region only. The number of (ETR-) NEOs known now is 239. As for Centaur-NEOs (CNEOs) found in the ‘Bottke’s region’, there are 20 from R2 and 11 from R1 with a dynamical lifetime of 21.3 and 35.8 kyr, respectively. Considering the entire lifetime of these objects, the frequencies of such episodes are 31/6.038 CNEOs per Myr, about 5 every Myr, or in terms of the average lifetime for the whole NEOs population, 50 CNEOs per 10 Myr. Nevertheless, the probability to see them together with the other NEOs is low because the whole population of CNEOs in the ETR should live on average about 28.5 kyr. Thus, the debiased population present in the ETR region should be currently 50 × 44 300 × 0.0285/10 CNEOs¹¹ (we scale the probability to the average dynamical time to that of NEOs) and 239 + 239 × 0.01 = 253 (10 per cent still unknown), the contribution of CNEOs to the current population of NEOs should be ∼4 per cent. This is consistent with Morbidelli et al. (2002), who assert that the contribution from the ‘comets’ should be up to 10 per cent of all the NEOs. However, a debiased estimate of NEOs for large eccentricities is also needed. As shown in Fig. 4, no NEOs are observed at e > 0.96 and beyond Jupiter very few. These results suggest that there should be more NEOs to be discovered than currently predicted. The unknown 1 km size NEOs (10 per cent of the current estimation) represent only a minimum.

Figure 4.

(a–e) and (a–i) phase space of real NEOs and all the possible orbits of the Centaur-NEOs (all clones).

We found that the Centaurs have the highest probability of encountering the Earth, followed Venus and then Mars, if we count per single body. However, if we count the total number of CEs, Venus has a slightly larger probability of CE. The probability of an Earth encounter for a body larger than 1 km is 7.5 per cent (Table 4).

The rate of arrival of large bodies (from Region 1, where the average size is D ∼ 50 km) in TP region is about 27.58 |$\, \mathrm{Myr}^{-1}$|⁠, which translates to 1 every 0.04 Myr. For the Earth-crossing orbits (CEs inside the Hill sphere), this rate is |$1.3 \, \mathrm{Myr}^{-1}$| (1 every ∼0.8 Myr). Napier (2015), for ≥100 km Centaurs, found 1.9 Chiron|$^{-1} \, \mathrm{Myr}^{-1}$| (1 every 0.5 Myr). However, that study considered any Centaur with a perihelion smaller than the semimajor axis of the Earth as an Earth-crosser. In our work, we found similar results (but more frequent), suggesting that Chiron-clones can have 185 Earth-crossing epochs in all, meaning that 185/15/3884 Chiron⁻¹ Myr⁻¹ or 1 every ∼0.3 Myr.

For R2, the rate of arrival is about 19 times greater than that for R1, with a similar rate for Venus (21) and Mars (17). There is a general trend of increasing rate with the decreasing distance to the Sun. In principle, Venus undergoes more CEs than the Earth, and the Earth experiences more encounters than Mars (Table 4), and each planet at least 1 every 50 kyr for the biased population. For Centaurs with D > 1 km, there is a CE approximately every 3–5 yr (3.5 yr for the Earth).

Napier (2015) found that integrating the orbits of 100 Chiron clones, leads to the mean total dynamical lifetime on an Earth-crossing configuration of ∼17.4 kyr, with a highly skewed distribution, and the median lifetime of 6.7 kyr. We found 3.767 Myr for Region 1, and 1.440 Myr for Region 2. This implies that Chiron should evolve relatively quickly in the planetary region and sooner than an average Centaur (for Chiron in particular, we found a mean total dynamical lifetime on an Earth-crossing configuration of ∼1010 kyr, an order of magnitude greater than Napier (2015)). Fig. 5 shows the rate of CEs of the clones of the known Centaurs every 1050 kyr. It is evident that in the first 3.5 million years (and in particular up to 1 Myr) the evolution is the most important factor for interactions with the Earth. After that, and for the next ∼10 Myr, Centaurs play an important role for R1. This implies that in a window of about 25 Myr, there is a possibility of high rate of CEs, which is higher in two intervals of time at the beginning and at ∼11 Myr from the initial orbital evolution. Again, this is due to Centaurs from R1; the average time to enter the planetary region from R1 is much longer than R2, 4414 and 722 kyr, respectively. The median values are 1418 and 100 kyr, respectively, for R1 and R2, meaning that more than the half of large Centaurs from R1 arrive after the first peak. The converse is also true for R2, where more than half of Centaurs entering the planetary region arrive in the interval time of the first bin (Fig. 5). An exponential decay can be seen in the entire population. We consider an exponential fitting function (the proportional rate growth):

\begin{eqnarray*} N_{\rm CE}=N_{\rm CE,0}-\frac{a}{k}\big (1-e^{kt}\big) \end{eqnarray*}

(8)

and we found N_{CE, 0} = 1109.53 ± 77.13, a = 1.6086 ± 0.2385, and k = 0.001461 ± 0.000115.

Figure 5.

Rate of CEs of the Centaurs clones with the Earth per kyr. Every bin (one every 1050 kyr) considered includes the contribution from the clones from the two different regions. Upper plot: the period of time when there are CEs. Bottom plot: magnification of the time interval around the second peak.

Some real (PHOs) are on very eccentric and large semimajor axis cometary orbits, thus they might belong to former comets, due to a break-up event (as also found in past studies, e.g. Horner, Evans & Bailey 2004a, Horner et al. 2004b; Galiazzo et al. 2016). In this work, we find those objects with a ≳ 2.6 au and e > 0.63 (see Fig. 6), but in particular with a > a_j for Region 2 (i.e. asteroid 1999 XS₃₅ that could be a former comet or a piece of it coming from Region 1, considering that its size is small ∼1 km, having an H = 17.2). We also found three sets of Centaur Earth encounters with the following orbits (Table 6) in the ‘main belt orbital plane’:

• a ∼ 2.7 au, e ∼ 0.64, and i ∼ 9°, i.e. asteroids 2016 UF₁₀₁, 2008 VE, 2016 YJ₄, 2016 EV₂₇, and 2008 DV₂₂, all asteroids smaller than 100 m in diameter (H ≲ 21.8)

• a ∼ 3.2 au, e ∼ 0.76, and i ∼ 15°. Asteroids closest to the average values are 2012 DG₆₁, 2008 OB₉. In terms of comets with similar orbits, those with the closest osculating elements are 3D/Biela, 15P/Finlay, 2007 T2 (Kowalski), 222 P/LINEAR, 364P/PANSTARRS, P/2009 WX₅₁ (Catalina), and of 141P/Machholz 2 (and 2-A and 2-D).

• a ∼ 3.4 au, e ∼ 0.77, and i ∼ 18°. In practice, this range includes the range of set (2), but with a higher trend in semimajor axis. Some minor bodies close to this orbital range are (those that are more likely to come from Region 2), e.g. 2012 GS₅ (a small asteroid) and the comets 73P/Schwassmann–Wachmann group and 103P/Hartley 2.

Figure 6.

Centaurs’ Earth encounters comparison with real NEOs which can cross the Earth’s orbit in the a − e and a − i phase space. Centaurs-R1 and -R2 represent the orbit of the clones found in the regions in the plots. The Earth-crosser asteroids (Earth-crossers ast.) and comets (Earth-crossers com.) are the orbits of the real observed population. In the panel which describes the a − i phase space, the dotted square and the dotted circles represent the real observed asteroids and the real observed comets, respectively. Clones of Centaurs are also shown as in the other panels. The full colour figure is available in an online version of this paper.

The asteroids with similar orbits to set 1 and 2 might have originated from former comets coming from Region 1. Additionally, comets in set 2 come from Region 1. For set 3, there are bodies coming from Region 2.

Typically, the Centaurs coming from Region 1 are on the edge of the Earth-crossing bodies in the a−e phase space, and they match the orbit of the Near-Earth-Comets (Fig. 6).

In the (a–i) space, there are also preferred regions at certain values. It is interesting to point out a gap beyond |$11^{+2}_{-1}$| au, in particular, at a specific value in retrograde orbits with the following inclinations range 85° ≲ i ≲ 110°. Unlike our clones, no currently known real NEO is seen at this range, corroborating these results. However, further studies, including observations, are needed to confirm the presence of this gap. From the numerical orbital evolution standpoint, it might be possible to cover this region by a larger number of clones. From the observational standpoint, this might be a bias due to difficulties of finding asteroids in this region. Nevertheless, we cannot exclude the possibility that this could be a real selection effect due to peculiar gravitational effect or due to the orbital evolution of the asteroids arriving in particular orbits into the planetary region. Integrations for regions at semimajor axis larger than 20 au do not show Centaur Earth encounters at low inclinations (below 70º) like the real ones. This could be a bias deflection due to the integrator or an observational gap. We also note that on average, Centaurs coming from Region 1, thus the largest ones, remain at higher inclinations (see Fig. 6 and Table 6).

The leakage from the MBAs appears inadequate for NEOs larger than 3–4 km and it is also an order of magnitude less than that cited in the previous studies (Menichella, Paolicchi & Farinella 1996; Minton & Malhotra 2010). The Centaurs complete this lack of NEOs in this aforementioned size range. Thus, the impact rate of the bodies of this size is complete with the contribution of the Centaurs added to the NEOs coming only from the main belt.

From Table 4, we also found that the contribution to the PHOs is at least of 1 object per Myr, and, applying the debiased model of Horner et al. (2004a), there should be at least 139 (=44300/318) km-size objects per Myr. However, the probability to see these objects is very low, considering that their dynamical time as NEOs, determined earlier, is 25 kyr. Thus, performing the same computation as done for all the NEOs population, leads to only one former Centaur or none, among the observed PHOs with a MOID¹² < 0.01 au (we consider the 11 observed PHOs inside 0.01 au because our numerical integration was limited to this planetocentric distance). In agreement of this finding, among the observed PHOs, there is only one comet, C/1983 H1 (IRAS-Araki-Alcock), which has also very similar orbits (a = 98.02, e = 0.9899, and i = 73.25) to our computed CNEOs (see Fig. 6).

Following the same approach and considering that only ∼11 per cent of the Centaurs can become active, we can compute the contribution to the short periods comets, or objects with a period less than 200 yr and a_C < 19 au. We find that from R2 and R1, 2335 and 1743 clones, respectively, have a_c in the ‘activity range’ with a dynamical lifetime of ∼581 kyr (R1) and ∼406 kyr (R2). Given that we do not have an estimate of a unbiased population, we can at least compare the observed bodies: from our estimate, 0.2(=(2335/158 + 1743/160)/15 × 0.11∼1 SPCs should be produced by the observed population of Centaurs. From the Horizon JPL catalogue, we find 96 short period km-size comets, thus considering that the number of comets missing to get the unbiased value is larger for the Centaurs than the SPC, the contribution of the Centaurs to the SPC is ≳1 per cent.

These minor bodies arrive on average at cometary speeds (Table 7), which are generally higher than those for Centaurs coming from Region 1, due to their more eccentric orbits (Table 6).

3.3 Impacts

Before computing the impact probability, we calculate the disruption rate for our sample using the Öpik relation (Öpik 1976) with which we can find the new average diameter per region, considering only the active comets (11 per cent). With this calculation, ∼62.5 per cent are in Region 1 [10 out of 16 given the results from Guilberte-Lepoutre (2012, table III)]. We use the averaged version of the Öpik’s relation:

\begin{eqnarray*} r_{\rm lost}=3.8/\sqrt{q}/T_{\rm rev}\qquad \mathrm{(cm)}. \end{eqnarray*}

(9)

Using the average perihelion (q) and orbital revolution period around the Sun (T_rev) of the Centaurs that enter in the terrestrial region, r_lost is the radius that loses the Centaur through its activity as a comet. The average period of revolution is |$\bar{T}=\bar{a}^{1.5}$|⁠, where a is the semimajor axis (in au taken through all its evolution, from the beginning of the integration to the end of its dynamical life) of the Centaurs that enter the terrestrial region.

We determined that the average dismantle per planet in Region 1, Region 2 and the total region is overall negligible,¹³ having a maximum average disruption of about 200 m in diameter (Table 8). This result suggests that the surviving bodies of the known sample coming into the inner solar system are generally not affected by the activity: the main reason for disruption by the activity before encountering a TP could be due for some isolated cases, such as comets with a small nucleus (D ≲ 1 km), thus Centaurs with D < <1 km are not taken into account for statistical computations.

Table 9 gives a summary of the impact probability per million years, impact velocity and diameter of possible craters. The probability of impacts for TPs was found by correlating the number of CEs to the relative distance for all the binning values available (see Fig. 7), and fitting these data by a parabolic equation:

\begin{eqnarray*} y=ax^2 \end{eqnarray*}

(10)

Figure 7.

Fit of the data for CEs of Centaurs from Region 2 with the Earth in its relative Hill sphere.

Table 3.

Size of Centaurs per region in kilometers and relative masses in kilograms.

Region	Dmin	|$\bar{D}$|	Dmax
	Mmin	|$\bar{M}$|	Mmax
1	1.6	51.5	300.2
	1.769 × 1013	3.134 × 1018	1.066 × 1020
2	0.5	12.6	72.0
	3.982 × 1011	7.396 × 1017	1.471 × 1018

Table 4.

CE statistics for the clones in Regions 1, 2, and the whole population of Centaurs. ‘R1’, ‘R2’, ‘TOT’, and ‘D₁’ stand, respectively, for Region 1, Region 2, the entire region (all the observed Centaurs), and all the Centaurs with D > 1 km (debiased population with 10⁷ or 44 300 and current rate: the average rate since the LHB has to be scaled by a factor of 1.49). CE_Ct%, CE, CE_T represent, respectively, the percentage of Centaurs (clones) having CEs inside the Hill’s sphere of the planet encountered, the total number of CEs inside the planet relative Hill’s sphere, the rate (over the lifetime of a Centaur in the region) of CEs (inside the planet Hill’s sphere) per time in million years, T is the average time of a CE with its minimum and maximum. |$T_{\rm CE_{Ct}}$|⁠, T_CE are the Centaur’s rate arrival in a planet crossing like orbit and rate of CEs in million years, respectively.

Region/planet	CECt%	CE	CET (Myr−1)	T (Myr)	|$T_{\rm CE_{Ct}}$| (Myr)	TCE (Myr)
R1 Venus	3.9	≈13	∼1	4.470 (0.016–25.010)	∼1.655	∼0.769
R2	10.4	≈51	∼28	1.644 (0.003–26.097)	∼0.115	∼0.036
D1	7.2	∼2022989  ∼8962	∼335043  ∼1484	∼3.057	∼8.423 × 10−6	∼2.99 × 10−6  ∼6.74 × 10−4
R1 Earth	4.3	≈13	∼1	3.767 (0.023–24.942)	∼1.493	∼0.801
R2	10.7	≈42	∼23	1.440 (0.004–26.201)	∼0.107	∼0.043
D1	7.5	∼1719958  ∼7619	∼284856  ∼1262	∼2.604	∼8.070 × 10−6	∼3.51 × 10−6  ∼7.92 × 10−4
R1 Mars	4.1	≈33	∼1	3.555 (0.017–24.925)	∼1.586	∼0.962
R2	10.0	≈33	∼18	1.395 (0.170–13.892)	∼0.113	∼0.056
D1	7.1	∼1343783  ∼5953	∼222554  ∼986	∼2.375	∼8.548 × 10−6	∼4.49 × 10−6  ∼1.01 × 10−3

Table 5.

Typical orbits and lifetime of Centaur-NEOs. Each second row per zone represents the minimum value for the respective osculating elements.

Zone	|$\bar{a}$| (au)	|$\bar{e}$|	|$\bar{i}$|	|$\bar{T}$| (ky)
All	2606.2	0.96	92.5	96
	2.5	0.51	0.3
ETR	3.4	0.78	32.9	∼15
	2.5	0.51	0.3

Table 6.

Statistics (average and standard deviation) of the orbital parameters for the clones in Regions 1 and 2, which undergo CEs with the TPs. ‘R1’ and ‘R2’ stand for Regions 1 and 2, respectively. For the Centaur clones having the Tisserand parameter range in the given region during the CE, Column CE is the number, CE%R and CE%T are the percentage in respect to the total and to the relative region.

Region (planet)	CE	CE%R	CE%T	|$\bar{a}\pm \sigma _a$|	|$\bar{e}\pm \sigma _e$|	|$\bar{i}\pm \sigma _i$|
R1 Venus TJ > 3	–	–	–	–	–	–
2 < TJ < 3	63	32	7	3.16 ± 0.40	0.820 ± 0.068	17.4 ± 18.1
TJ < 2	137	69	14	13.07 ± 19.89	0.92 ± 0.07	102.6 ± 55.7
R2 TJ > 3	–	–	–	–	–	–
2 < TJ < 3	185	24	19	3.22 ± 0.46	0.82 ± 0.06	18.7 ± 16.0
TJ < 2	583	76	60	9.80 ± 14.10	0.91 ± 0.06	114.1 ± 52.0
R1 Earth TJ > 3	6	3	1	2.70 ± 0.03	0.64 ± 0.01	8.98 ± 3.10
2 < TJ < 3	92	48	11	3.23 ± 0.51	0.76 ± 0.09	14.6 ± 11.4
TJ < 2	94	49	12	15.27 ± 19.05	0.93 ± 0.06	104.9 ± 54.7
R2 TJ > 3	–	–	–	–	–	–
2 < TJ < 3	176	28	22	3.42 ± 0.63	0.77 ± 0.09	17.5 ± 14.5
TJ < 2	455	72	56	10.19 ± 12.05	0.89 ± 0.08	116.6 ± 52.0
R1 Mars TJ > 3	–	–	–	–	–	–
2 < TJ < 3	63	41	10	3.66 ± 0.78	0.69 ± 0.11	13.0 ± 11.6
TJ < 2	91	59	14	18.73 ± 28.55	0.88 ± 0.09	111.8 ± 50.9
R2 TJ > 3	9	2	1	2.91 ± 0.16	0.53 ± 0.07	10.6 ± 3.8
2 < TJ < 3	155	32	24	3.98 ± 1.03	0.72 ± 0.11	19.3 ± 14.8
TJ < 2	325	66	51	10.80 ± 9.31	0.85 ± 0.09	117.0 ± 49.9

Table 7.

Velocity at infinitum for a CE with each planet in kilometers per second.

Planet	|$v$|min	|$\bar{v}$|	|$v$|max
Venusa	12.4	48.9	84.8
Eartha	9.4	38.8	71.7
Marsa	7.2	30.7	59.1
Venusb	11.2	46.3	84.2
Earthb	9.8	40.9	71.5
Marsb	5.9	27.6	60.3

^aand ^bmean Regions 1 and 2, respectively.

Table 8.

Average semimajor axis (⁠|$\bar{a}$|⁠), perihelion (⁠|$\bar{q}$|⁠), orbital revolution period (⁠|$\bar{T}_{\rm rev}$|⁠), the lifetime of a body when it is in a cometary orbit (when it can be active) a_c < 12 au (⁠|$\bar{T}_{\rm Opik}$|⁠). All values are computed per relative region. These values provide the remaining diameter of the clones, where we report the values of the average disruption in one dimension |$\bar{D}_{\rm kill}$| and the average value found for diameters after disruption, |$\bar{D}_O$|⁠.

Planet	|$\bar{q}$| (au)	|$\bar{a}$| (au)	|$\bar{T}_{\rm rev}$| (Myr)	|$\bar{T}_{\rm Opik}$|	|$\bar{D}_{\rm kill}$| (km)
R1	3.66	784.55	21975	724	0.065
R2	3.79	325.16	5863	677	0.225

Table 9.

Impact probability (IP) per million years and impact rate (K) in Myr for each planet. ‘I’, ‘II’ (for the observed population and thus for a body of their average diameter), and ‘C’ stand for, respectively, Region 1, Region 2, the entire region (all the Centaurs) for the entire real population [unbiased with 10⁷ bodies (‘S’) and with 44 300 (‘H’) population with a diameter larger than 1 km]. Out means more than the solar system age (no possible impacts after the LHB). W_Earth is the water release in Earth masses from the LHB until now. For the water released by the unbiased population, we considered an asteroid with D = 5  km, which is the mode value of all the diameters.

Planets	IPI (10−8ky−1)	IPII (10−8ky−1)	IPC(10−2ky−1) – S	IPC (10−5ky−1) – H
	KI (My)	KII (My)	KC (ky) – S	KC (ky) – H
			W (WC, Earth10−6) – S	W (WC, Earth10−6) – H
Venus	2.67 ± 0.01  3.98a ± 0.02	47.26 ± 0.01  70.42a ± 0.01	1.6 ± 0.1  2.3a ± 0.1	7.0 ± 0.1  10.3a ± 0.1
	out  outa	2116 ± 2  1420a ± 1	64 ± 1  43a ± 1	out  965a ± 1
			∼3.40  ∼192.4a	∼0.02  ∼0.9a
Earth	2.26 ± 0.01  3.37a ± 0.01	52.69 ± 0.01  78.5a ± 0.01	1.7 ± 0.1 2.6a ± 0.1	7.7 ± 0.1  11.4a ± 0.1
	out  outa	1898 ± 1  1274a ± 1	58 ± 1  39a ± 1	out  877a ± 1
			∼3.7  ∼211.8a	∼0.02  ∼0.9a
Mars	0.62 ± 0.02  0.09a ± 0.01	14.90 ± 0.01  22.19a ± 0.01	0.5 ± 0.1  0.7a ± 0.1	2.2 ± 0.1  3.2a ± 0.1
	out  outa	out  outa	205 ± 1  138a ± 1	out  3104a ± 1
			∼1.1  ∼59.8a	∼0.01  ∼0.3a
Jupiter			1112 ± 45  1657a ± 67	4.9 ± 0.2  7.3a ± 0.2
			0.090 ± 0.002  0.060a ± 0.002	20.03 ± 0.42  13.63a ± 0.28
			∼2406.4  ∼136247.6a	∼10.7  ∼603.5a
Saturn			278 ± 11  414a ± 17	1.2 ± 0.1  1.8a ± 0.1
			0.360 ± 0.009  0.242a ± 0.006	81.21 ± 1.70  54.51a ± 1.14
			∼601.6  ∼34061.9a	∼2.7  ∼150.9a
Uranus			–	–
Neptune			35 ± 1  52a ± 2	0.2 ± 0.1  0.2a ± 0.1
			2.878 ± 0.071  1.932a ± 0.047	64.97 ± 13.56  436.05a ± 9.34
			∼75.2  ∼4257.7a	∼0.3  ∼18.9a

^aMeans a period from now to the LHB, without ^a means from now until 0.1 Gyr ago.

Here, y is the number of CEs per given distance (x) and a is the relative fitting parameter. This parameter was subsequently adapted to the number of clones (N_c = 15) per Centaur and computed per million years (considering the average lifetime |$\bar{T}_{\rm life}$| of the Centaurs), such that the impact rate per million year is

\begin{eqnarray*} \begin{split} P=\frac{y}{15 \bar{T}_{\rm life}},\\ P=\frac{ax^2}{15 \bar{T}_{\rm life}}. \end{split} \end{eqnarray*}

(11)

The error associated with the measurement was computed with the partial derivatives:

\begin{eqnarray*} \Delta P &=& \sqrt{\bigg (\frac{\partial P}{\partial a}\sigma _a \bigg)^2 +\bigg (\frac{\partial P}{\partial T_{\rm life}} \sigma _{T_{\rm life}}\bigg)^2}\nonumber\\ &=&\frac{1}{15}\sqrt{\bigg (\frac{x^2}{T_{\rm life}} \sigma _y \bigg)^2+ \bigg (\frac{ax^2}{T^2_{\rm life}} \sigma _{T_{\rm life}} \bigg)^2},\nonumber\\ \sigma _y &=&\frac{\partial y}{\partial a}\sigma _a = x^2_0 \sigma _a. \end{eqnarray*}

(12)

The error for the average lifetime and for the fitting constant, a, are their standard deviation of the mean. x₀ is the radius of the planet engaged.

Thus, the final probability of impacts per million years is |$P=a*x^2_{\rm imp}/15/\bar{T}_{\rm life}\pm \Delta P$|(My⁻¹).

For GPs, we only count the number of real impacts (found in numerical orbital integrations) per given time, because we are more interested in TPs in this work. Table 9 shows the forecasted number of impacts per region of the TPs, the frequency and recurrence time expected per regions and global region (for Centaurs with D > 1 km). Centaurs from Region 1 do not seem to present any threat for any of the 3 TPs in the recent solar system (computations are for an impact in more than 3.8 Gyr).

We find that no impact should happen from bodies from Region 1. In Region 1, even though there are bodies on average ∼51 km in size, the possibility of an impact is even more unlikely. For example, Napier (2015) found that a 100 km comet striking the Earth would release ∼1000x the energy required to form the 150 km Chicxulub crater, and would presumably remove the entire biosphere. The impact crater statistics show this scenario to had been an unlikely event over the 600 Myr time-scale of Phanerozoic, and for other results even longer than ∼4 Gyr, assuming the current orbits of the presently known¹⁴ planets in our solar system.

As for Region 2, for an average Centaur of size 12.1  km, an impact on the Earth and Venus happen every 1.3 and 1.4 Gyr, respectively. This implies that one or two catastrophic impacts on these planets after the LHB could have been possible. These two impacts should create craters of size 133  km on Venus and 119  km on the Earth. Although it is not feasible to confirm whether such an impact had occurred on the Earth, it is interesting to entertain the possibility. Of all confirmed impact structures on the Earth, the Sudbury crater falls within the size range (age 1850 ± 3 Ma and D ∼ 130 km). The Vredefort crater might have been created by a Centaur (Johnson & Bowling 2014); however, if this is correct, then it follows that the Earth is right now within the interval of time susceptible to the collision with a large impactor. Thus, an extinction event could hypothetically occur at present and in a time window of about 60 kyr. A large crater like that would not come only from the supposed E-belt (Bottke et al. 2012) objects but also from Centaurs. Additionally, a few of the Centaurs seem to become active. Out of nearly 150 Centaurs detected to date, the comet-like activity has been reported in only 16 of them, with the current fraction of active Centaurs being f ∼ 10⁻¹. Therefore, many of Centaurs could potentially preserve most of their mass if they do not approach too close to the Sun.

The rate of impacts for all the Centaurs having a smaller average size (roughly ∼5  km, the mode value of the presently known Centaurs), Venus, Earth, and Mars experience an impact in a time less than 64, 58, and 205 kyr, respectively. Thus, there is at least one impact in this limit. For an impactor of 1 km in size, the resulting craters on the Earth would range from 6.9 (see Fig. 8) to 18.9  km in diameter.

Figure 8.

The time series of the impact crater as a result of a Centaur from R1 colliding with the Earth. The initial time at 0 s is when the collision with the surface starts. It takes ∼3 min for the crater collapse to finish. The impactor diameter is 1 km. Considering that the simulations are axisymmetric, the impact velocity was scaled to the vertical component only. The non-scaled impact velocity is |$v$|_i = 9.4  km s⁻¹(at 45° and scaled to the vertical component is |$v$|_i = 6.6  km s⁻¹).

The transport of water is negligible compared to the total present amount of water (∼5 × 10⁻⁴M_Earth; Lecuyer, Gillet & Robert 1998; Morbidelli et al. 2000). For the Earth, it should be of the order of at least 10⁻⁴ times the current amount of water, which implies that the water present now comes mostly from other sources or minor bodies older than 3.8 Gyr. However, this amount of water is still significant compared to the amount of with the water in the atmosphere, which is about¹⁵ 8.54 × 10⁻⁹ the total amount of water on the Earth.

We provide a summary of the craters (Table 10) and include three sets of speeds (average, maximum, and minimum) for each impactor size.

Table 10.

Crater dimensions (D_C) for each TP in kilometers, impactor dimension (D_i) in kilometres, and impact velocity (⁠|$v$|_i) in kilometres per second at an impact angle of 45° (minimum, average, maximum).

Region	Di,min|$v$|i,min	|$\bar{D_{\rm i}}$||$\bar{v_{\rm i}}$|	Di, max|$v$|i,max	|$D_{i,D_{1km}}$||$\bar{v_i}$|
Planet	DC,min|$v$|i,min	|$\bar{D_{\rm C}}$||$\bar{v_{\rm {\rm i}}}$|	DC,max|$v$|i,max	DC,1 km|$\bar{v_i}$| (DC,1 km,min − DC,1 km,max)
Region 1	1.6	51.2	300.2	1.0
Venus	12.6 12.4	484.9 48.9	3225.6 84.8	15.1 48.9 (8.1–21.1)
Earth	10.5 9.3	418.5 38.8	2844.1 71.7	13.0 38.8 (6.8–18.6)
Mars	10.7 7.2	436.7 30.7	3005.2 59.1	14 30.7 (7.6–19.7)
Region 2	0.5	12.1	72.0	1.0
Venus	1.8 11.2	133.1 46.3	912.9 84.2	14.8 46.3 (4.2–10.7)
Earth	3.9 9.8	118.7 40.8	610.7 71.5	13.2 40.8 (7.7–21.0)
Mars	3.0 5.9	118.8 27.6	862.5 60.2	13.5 27.6 (7.0–19.9)

While it would be speculative to place any constraints on possible Earth craters (based on the confirmed impact structures record) that might have been formed by Centaurs, such possibility should not be entirely dismissed. A large range in crater sizes could be a consequence of very high impact velocities of the Centaurs to the planets, i.e. they arrive to the Earth with an average speed of ∼40  km s⁻¹, which is at least 1.5 times the average typical velocity of a V-type NEA asteroid impacting the Earth, |$v$|_imp = 20.0  km s⁻¹(Galiazzo et al. 2016).

4 CONCLUSIONS

In this paper, we presented our results from modelling the total Centaur population to investigate the threat of impacts with TPs. The Centaur population analysed in our work was sub-divided into the biased population (BP), the known observed Centaurs, and the unbiased population (UP), the estimated real population of Centaurs larger than 1 km in size. In particular, two regions of Centaurs are emphasized and divided by the mean motion resonance J4:1 at 13.11 au.

Region 2 is the source of Centaurs that undergo more encounters (consequently more impacts) with TPs, particularly in these orbits: 5.5 au ≤ a ≲ 8 au and with eccentricity 0.09 ≤ e ≲ 0.75 and i ≲ 40°. A minor, but relevant contribution, is in Region 1, further divided into two sub-regions: (a) a ≲ 20 au, e ≲ 0.63 and i ≲ 35°, and (b) 26 au ≤ a ≲ 28 au and with eccentricity 0.8 ≲ e ≲ 1 and 15° ≲ a ≲ 30°.

There is also a considerable contribution to TPEs from Centaurs with retrograde orbits, i ∼ 140° and a ∼ 22 au. However, the dominant inclination for Centaur sources of TPEs is i ∼ 20°, similar to Jupiter comets. Centaurs can become Jupiter comets in only ∼20 kyr.

In general, the lifetime of Centaurs is 6.0 Myr. This period of time is comparable to the replenishment of NEOs by the MBAs via the strongest mean motion resonances (ν₆ and J3:1). This shows that the orbital and physical distribution of the NEOs is strictly dependent on the Centaurs and MBAs together. Thus, distinguishing their orbital evolution can help to understand and characterize different groups of the NEOs and their interaction with other bodies in the solar system. The lifetime for Centaurs in Region 2 is substantially shorter than that of Centaurs in Region 1; these lifetimes for R2 and R1 are 1.8 and 10.3 Myr, respectively. This work shows that in order to study the orbital evolution of the Centaurs, it is not enough to study only one Centaur with its relative clones, as it was performed in previous studies (i.e. Napier (2015), even if such studies could provide some constraints on the orbital evolution. Centaurs like Chiron can be used as a proxy for its own relative region; in particular, Chiron can be used as a proxy for R1.

The majority of the Centaurs escape (82.4 per cent), 6.2 per cent collide with the Sun (of which 56.4 per cent are from Region 2) and 0.9 per cent undergo impacts with planets. For the latter value, Region 2 dominates in a ratio 13:1 over Region 1.

53 per cent of the Centaurs enter the TP region, and of these, 57 per cent come from Region 2. They arrive on a time-scale of millions years, comparable to the one relative to the main belt bodies in becoming NEOs, which are a mixture of the main belt objects and Centaurs. They survive as NEOs 96 kyr and their contribution to this population is ∼4 per cent. This study also confirms that the preferable unbiased population size of km-size Centaurs is 44 300, found by Horner et al. (2004a), but not consistent with the higher value found by Sheppard et al. (2000). Our results indicate that the contribution of the Centaurs to the PHOs is about 9 per cent, and for the SPC, it is at least ∼1 per cent, with an average dynamical lifetime (as an SPC) ranging from 0.4 for R1 to 0.6 for R2. CEs occur with all the planets, mainly with Jupiter. Among the TPs, Venus undergoes the most encounters. For the biased Centaurs population from Region 1, about one Centaur encounters the TPs every million years, as opposed to Region 2, where 28, 23, and 18 encounters take place every 1 million years, meaning one encounter every 36, 43, and 56 kyr, respectively, for Venus, the Earth, and Mars. For the unbiased Centaurs population, there is one encounter at a more frequent scale, every 0.7, 0.8, and 1.0 kyr for Venus, the Earth, and Mars, respectively. About 7–8 per cent of all the Centaurs (the unbiased population) could encounter TPs (the Earth 7.5 per cent). For large bodies (D > 50 km), encounters occur at a rate of about 1 every 40 kyr, and because there are about 11 per cent active Centaurs (11 per cent of all Centaurs become comets), it means that the probability to have comets arriving in the TP region every 40 kyr is 11 per cent: therefore, there could be comets with a nucleus larger than 50 km, roughly 1 every ∼0.4 Myr.

For all the Centaurs, we find that in a time frame of 25 Myr, there are two time intervals with a high rate of CEs. The intervals between the two maxima are 8 and 17 Myr. The rate of arrivals for Centaurs coming from Region 2 is 19 times greater than that from Region 1.

We find that PHOs with a ≳ 2.6 au and e > 0.63 are potentially former Centaurs, in particular those with a > a_J which should come from Region 1. One example is the asteroid 1999 XS₃₅, which could be a former nucleus of a comet or a part of it (diameter  = D > 1 km). In addition, comet C/1983 H1 (IRAS-Araki-Alcock) is likely a former Centaur. Centaurs that undergo CE with the Earth have three sets of main orbits, with peculiar orbits that mostly reside in the main belt:

• a ∼ 2.7 au, e ∼ 0.64, and i ∼ 9°, i.e. real small asteroids in range (D < 1000 m) are 2016 UF₁₀₁, 2008 VE, etc.

• a ∼ 3.2 au, e ∼ 0.76, and i ∼ 15°. Here, we find real asteroids and comets in this range, like i.e. 2012 DG₆₁ and 222 P/LINEAR.

• a ∼ 3.4 au, e ∼ 0.77, and i ∼ 18°, e.g. like 2012 GS₅ and 74P/Schwassmann–Wachmann group.

The Centaurs in CEs with Venus and Mars orbits are mainly in the outer main belt at high eccentricities (e ≳ 0.8 for Venus and e ≳ 0.4 for Mars) and i ≳ 10° and even retrograde orbits. Usually, the Centaurs coming from Region 1 are on the edge of the a–e phase space for the Earth-crossing bodies. The sub-region with inclinations, 85° ≲ i ≲ 110°, and a > 11 au is unlikely a source of terrestrial planetary encounters. Objects with a > 20 au and i ≲ 65° do not come from Centaurs.

Concerning impacts, the impact rate for Centaurs complete the gap found in the main belt objects with sizes larger than ∼3 km.

Using the Öpik relation, we find that the largest Centaur that can be destroyed completely after activity as a comet is 200 m in diameter. Thus, this result is irrelevant for the rates of impacts computed for the biased (no minor body of the sample is less than 400 m, and the minimum size is ∼0.47 km) and unbiased populations (D > 1 km).

Our results suggest that Region 1 is not a source of catastrophic collisions in the recent solar system (after the LHB). Instead, Region 2 is the likely source of potentially dangerous objects: for a Centaur with an average size of 12.6 km, an impact with the Earth and Venus could happen every 1.3 and 1.4 Gyr, respectively. This means that one or two extinction events could have happened on the Earth due to the Centaurs (apart from the main belt objects) in the last 3 billion years. These impacts can be as large as ∼130 km on Venus and ∼120 km on the Earth, considering an average impact velocity of 46.9 and 40.9  km s⁻¹, respectively. Although it is not possible to make any firm conclusions, we note that some Earth craters fit these size ranges and are consistent with the timing, i.e. the Sudbury Crater (age 1.85 Gyr and 130 km in diameter) and the Vredefort crater (age 2.023 Gyr and D = 120 km). If these craters had been created by a former Centaur, it follows that we might be presently living in a time epoch when a catastrophic impact event capable of extinction might be possible.

For Centaurs larger than 1 km, the unbiased population, there should be on average one impact every ∼9.7, ∼8.8, and ∼31.0 Myr for Venus, the Earth, and Mars, respectively, since the LHB. In the recent solar system, however, it should happen less frequently, approximately every 14.4, 13.1, and 46.3 Myr, respectively, for Venus, the Earth, and Mars. On the Earth, the craters should be at least 7 km in diameter for these kind of impacts.

In terms of water release to the Earth, we compute that Centaurs could release a value of roughly 10⁻⁵ of the total water present now since the LHB and, at the moment they should release to our planet 2 × 10⁻¹⁰ (current) Earth-water-like content every million year and for Mars at least less than half of this value. It is interesting to note that this result is also applicable to Mars, where the water delivery is 3 × 10⁻⁵ the total water present on the Earth now, the amount consistent with the amount for a ‘dry’ planet (and only in the recent solar system).

Our work suggests that the Centaurs are very important for the evolution of a planet, and could have profound implications geologically and also astrobiologically. Finally, we would like to invite the scientific community to observe more this population of objects given that it is of vital impact to our planet (and also to Mars), because we need the tools to better discriminate their size distribution (and thus a better estimate of their impact rate) and their physical characteristics, the amount of water, and other important material they might bring to the TPs.

ACKNOWLEDGEMENTS

MG wishes to acknowledge the support by the Austrian FWF project J-3588-N27, and thanks Y. Cavecchi for suggestions to improve the computations; S. Aljbaae running some of our numerical orbital integrations on a cluster at FEGP-UNESP (Brazil); and Western University for providing computational facilities (the TITAN supercluster and the SHARCNET). EAS gratefully acknowledges the Natural Sciences and Engineering Research Council of Canada (NSERC) Postdoctoral Fellowship programme for supporting this project. RD acknowledges the support by the Austrian FWF project S11603-N16. We gratefully acknowledge the developers of iSALE-2D (http://www.isale-code.de), the simulation code used in our research, including G. Collins, K. Wunnemann, D. Elbeshausen, B. Ivanov, and J. Melosh. The authors also thank the referee Mario Melita for valuable advice allowing us to improve this paper.

Footnotes

1

Updated at 2015 August 24.

2

695700 km, the volumetric mean radius of the Sun, see https://nssdc.gsfc.nasa.gov/planetary/factsheet/sunfact.html.

3

For the Hill sphere of the Earth, we considered the sum of the masses of the Earth and the Moon.

4

A CE depends entirely on the relative gravitational attraction of each planet, and thus the CE region must be at least proportional to the relative Hill sphere; see also Galiazzo, Bazso & Dvorak (2013).

5

We used the value given by the Johnstone’s archive, http://www.johnstonsarchive.net/astro/tnodiam.html.

6

Data from Morbidelli from the work Brasser & Morbidelli (2013).

7

This is an extrapolation from the results for the TNOs population, which is the source region of the Centaurs (Brasser & Morbidelli 2013; Galiazzo, Silber & Bancelin 2016).

8

For the purpose of our study this is sufficient, considering that this computation provides a rough estimate on which comet size might highly affect (completely or almost completely destroy) an active Centaur by solar radiation or tidal effect due to CEs.

9

The ratio between the number of impacts per clone per region for the giant planets and the impact probabilities for TPs.

10

If the diameter is not known, we selected the objects absolute magnitude H < 17.75.

11

We use the debiased value found by Horner et al. (2004a), otherwise the number of CNEOs would be unreal, unrealistically larger than the present population of NEOs.

12

Minimum orbit intersection distance.

13

We consider the average ‘Öpik’ disruption, and use the average parameters to compute it.

14

Escluding, i.e. a possible, planet IX, or more after Neptune Gomes, Matese & Lissauer (2006).

15

1 per cent of the total mass of the atmosphere of the Earth: 5.1 × 10¹⁸ kg (https://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html).