The origin of (47171) Lempo-like Kuiper belt triple systems during binary–binary interactions

ABSTRACT

The trans-Neptunian triple Lempo–Paha–Hiisi is composed of a tight inner binary with components of similar size and an outer companion about half their size orbiting 10 times further away. Large trans-Neptunian objects like Pluto also have multiple small moons, but Lempo’s structure is unique in the Solar system, and the place and timing of its origin is still a subject of debate. We propose a new formation mechanism able to form a large number of systems like Lempo–Paha–Hiisi, which involves binary–binary close encounters in the primordial planetesimal disc at 30–40 au. Some of these systems were then implanted in different populations of the trans-Neptunian region during Neptune’s outward migration. Our results strongly support that the 4:7 resonant multiple object Manwë–Thorondor was once a triple system similar to Lempo–Paha–Hiisi, but the orbit of the inner binary evolved by tides, becoming a contact binary. As with Lempo–Paha–Hiisi, it should have originated in the planetesimal disc below 30–40 au. Triple systems like Lempo–Paha–Hiisi or Manwë–Thorondor could not have formed in situ and the existence of this kind of system is not expected in the cold classical Kuiper belt.

Kuiper belt: general, Kuiper belt objects: individual: Lempo–Paha–Hiisi, Comets: general

1 INTRODUCTION

Lempo–Paha—Hiisi is a trans-Neptunian hierarchical triple system composed of a tight inner binary with components of similar size and an outer companion about half their size orbiting 10 times further away (Trujillo & Brown 2002; Benecchi et al. 2010). All large trans-Neptunian objects like Pluto have multiple small moons, but Lempo’s structure is unique in the Solar system. The place and timing of its origin are still a subject of debate (Nesvorný, Youdin & Richardson 2010; Correia 2018). In contrast to the Lempo system, all other known triples in the Solar system have their orbits almost regularly spaced, with one component much smaller than the others, with the most distant component being the largest (Johnston 2020). Unveiling the possible origin of the Lempo system is relevant, as the architecture of multiples holds clues to their formation story and the conditions prevailing in the primitive outer Solar system, but its origin is still unclear. Capture theories proposed so far failed to reproduce the orbital characteristics of observed binaries, especially the distribution of their orbital inclinations (Nesvorný et al. 2010; Brunini 2020). Also, triple formation requires multiple captures, a very unlikely event. Gravitational collapse of pebble clouds in a turbulent gas disc would be efficient in producing binaries and, in some particular conditions, can also produce triple systems. However, such triples do not seem to match the orbital structure of Lempo–Paha–Hiisi (Nesvorný et al. 2010). The non-detection of triple systems in the cold classical Kuiper belt, where the number of known binaries is much higher than in the resonant populations (Noll et al. 2020), argues against this formation mechanism of triple systems. The fragile dynamical stability of Lempo–Paha–Hiisi (Correia 2018) also casts doubt on the place and time of its origin, leading to speculations about a possible recent formation at the place in which it currently resides.

As Pluto is at 3:2 resonance with Neptune, all other bodies at this resonance are known as plutinos. They are found today at ∼39.4 au, but like the objects of other resonant populations in the trans-Neptunian region, the plutinos were born closer to the Sun, and were later gradually implanted in the zone where they reside today (Nesvorný & Vokrouhlický 2016; Malhotra 2019; Volk & Malhotra 2019) as Neptune migrated outwards to its present orbit. Many details of the implantation process are still unknown (Volk & Malhotra 2019), but reliable dynamical models agree that 1 out of 1000 of the original objects from this inner massive disc was successfully implanted in the plutino group (Nesvorný & Vokrouhlický 2016; Volk & Malhotra 2019), which still conserves ∼30 per cent of its original members (Tiscareno & Malhotra 2009).

The planetesimal disc between 15 and 30 au must have been massive enough (20–30 M_⊕) to be able to initiate Neptune’s migration, and to account for the present mass in some of the trans-Neptunian populations, like the scattered disc (Brasser & Morbidelli 2013). The fact that Neptune stopped at its present orbit indicates that the disc cannot be so massive beyond 30 au (Nesvorný et al. 2020).

Observational and theoretical evidence suggests that almost 100 per cent of the planetesimals in the outer Solar system were born as part of binary systems (Brunini & Zanardi 2016; Fraser et al. 2017). Several erosive mechanisms then dissolved most of them, giving rise to pairs of single objects (Brunini & Zanardi 2016; Nesvorný & Vokrouhlický 2019) and the fraction of binaries observed today. In this context, binary–binary close encounters should have been the norm during primordial stages, when the fraction of them was still high. Because of the huge number of parameters involved, the outcome of a binary–binary close encounter is complex to anticipate. There are, a priori, two possible situations where a triple might be formed during such a close encounter, as described below.

Binary–binary exchange-reaction (BBER) channel: A strong binary–binary interaction occurs when the closest approach distance of two binaries is of the order of the binary separations, and also when the encounter speed |$v$|_REL is small (of the order of or less than the escape velocity of the system). These interactions can cause the destruction of one of the binaries, with the hyperbolic ejection of one of its members. In favourable conditions, the other component could remain gravitationally bound to the remaining binary, forming a triple system. This mechanism has proved to be capable of producing triple stars in stellar clusters (Aarseth 2004).

Binary–binary collision (BBC) channel: During a close encounter, a component of one binary can impact one of the components of the other binary. If the impact is at low relative speed, such as one that can produce objects with bi-lobe-shaped structures (Jutzi & Benz 2017), both objects can merge, forming a single object. A fraction of the kinetic energy is dissipated in the inelastic collision, and if the total energy of the system becomes low enough, the resulting three bodies can remain gravitationally bound, forming a triple system.

To explore these two possibilities, we performed several numerical experiments of close encounters of a plausible synthetic population of binaries. We set the size and mutual orbits of their components according to the conditions prevailing in the region just outside Neptune’s primordial orbit (15 au), up to the edge of the massive protoplanetary disc (30 au). We did not perform simulations in the less-massive disc extension beyond 30 au, since it is still poorly constrained (Nesvorný et al. 2020). Nevertheless, we will make some considerations below and give special attention to the distribution of relative planetesimal velocities expected in that region.

In the next section, we describe the numerical simulations. In Section 3, we analyse the main results obtained, and the last section is devoted to the discussion of the results and some conclusions.

2 NUMERICAL SIMULATIONS OF BINARY–BINARY CLOSE ENCOUNTERS

In this section, we will describe the details of the numerical experiments that we have carried out and the initial conditions used for each one of the proposed formation channels.

For each numerical experiment, we have generated the size of the components and the initial mutual orbits of both binaries based on our present knowledge of binaries in the trans-Neptunian region (Brunini 2020; Johnston 2020; Noll et al. 2020). We have set the separation of the mutual orbits of the binaries in the interval 0.005 ≤ a_B/R_H ≤ 0.1, where a_B is the binary mutual semimajor axis and R_H is the Hill radius of the considered binary, which is given by

$$\begin{eqnarray*} R_\mathrm{H}=\mathrm{a}_\odot \left(\frac{ M_\mathrm{BIN} }{ 3 \, \mathrm{M}_\odot }\right)^{1/3}, \end{eqnarray*}$$

(1)

where M_BIN is the combined mass of the corresponding binary, a_⊙ the semimajor axis of its heliocentric orbit, and M_⊙ the solar mass. Objects with larger a_B have little chance of surviving as binaries in the outer primordial disc, due to effects such as collisions with the sea of planetesimals and close encounters with Neptune during its outward migration phase (Nesvorný et al. 2011; Porter & Grundy 2012; Nesvorný, Vokrouhlický & Morbidelli 2013; Brunini & Zanardi 2016; Nesvorný & Vokrouhlický 2019). Also, smaller values of a_B are subject to strong orbital decay due to tidal friction, on time-scales much shorter than the duration of the pre-instability phase of the outer Solar system (Porter & Grundy 2012; Nesvorný et al. 2013; Brunini & Zanardi 2016). We adopted initial circular orbits and a random distribution of orbital inclinations for all the binaries. Each binary has a primary with a diameter chosen randomly in the interval D_p = 100–300 km, corresponding to the diameter of most trans-Neptunian binaries of equal size (Johnston 2020). Then, the diameter of the secondary is generated at random in the interval 0.5 D_p ≤ D_s ≤ 1.0 D_p (Johnston 2020). We used a bulk density of ρ = 0.5 g cm⁻³, in agreement with the value found for the Lempo system (Benecchi et al. 2010). In each close-encounter simulation, one binary starts with its centre of mass at the centre of the mutual Hill sphere, defined as

$$\begin{eqnarray*} R_\mathrm{H}^\mathrm{MUT}=[R_{\mathrm{H}_{1}}^3+R_{\mathrm{H}_{2}}^3]^{1/3}, \end{eqnarray*}$$

(2)

where the sub-indices 1 and 2 correspond to each of the binaries involved in the close encounter. We also suppose the same a_⊙ for both binaries. The centre of mass of the second binary starts at a random point on the surface of |$R_\mathrm{H}^\mathrm{MUT}$|⁠. The centre of mass starts approaching with a velocity vector oriented at random, and whose magnitude was also taken at random from given distribution functions. To compute the distribution of relative speeds we assumed the outer-disc model of the planetesimal disc between 15 and 30 au given by Nesvorný et al. (2011). There is a marked gradient of dynamical excitation in the disc with distance from the Sun. Therefore, as in Morbidelli & Rickman (2015), the disc was conveniently divided into three zones: Zone 1 with 15 < a < 20 au, Zone 2 with 20 < a < 25 au, and an outer Zone 3 with 25 < a < 30 au. Each zone has a particular intrinsic collisional probability P_i (i.e. the mean rate at which mutual encounters of bodies occur within a unit cross-section) and mean relative velocity 〈|$v$|_REL〉. We took previously computed P_i values for each zone (Morbidelli & Rickman 2015). However, as we have already discussed, the encounters that matter for triple formation are those produced at low relative velocity, so we must obtain their distribution. We proceeded as follows: Assuming that the dynamical state of the disc does not evolve much in the time previous to the instability phase of the outer planets, we have evolved a set of 10 000 test particles under the gravitational attraction of the Sun, neglecting self-interaction and assuming the same dynamical state given in Nesvorný et al. (2011). The orbits of the planetesimals are then only driven by the central attraction of the Sun, and their dynamical evolution can be obtained very quickly by two-body analytical equations. We believe that the inclusion of the gas and ice giants is not necessary for this task because, during the pre-instability phase, the outermost ice giant has a ∼ 12 au, whereas the outer planetesimal disc that we considered here starts at 15 au. In this situation, the viscous stirring of the disc is very fast (Nesvorný et al. 2011). It reaches dynamical equilibrium a few Myr after the dispersal of the nebular gas (Nesvorný et al. 2011), in such a way that the rms values of e and i do not evolve in time. Therefore, we assumed, as in Morbidelli & Rickman (2015), that the collision probability and relative velocity for each pair of simulated binaries were not necessary. Averaged values are enough for our study. We do not expect that the inclusion of the giant planets in the computation of P_i and v_REL could cause substantial changes in these average values. In addition, we remark that the accurate dynamical condition of the outer Solar system is still a subject of intense debate (see e.g. Volk & Malhotra 2019). We evolved the entire disc using a small step size, and each time two particles were closer one to another by 5 × 10⁵ km, which is representative of their mutual R_H, we calculated the relative speed. We computed their distributions for each zone after collecting 5 × 10⁵ values, which are shown in Fig. 1. One remarkable feature of these distributions is that the mean values are in excellent agreement with those given in Morbidelli & Rickman (2015). These are the distributions that we used to generate the random relative velocities for all the simulations.

Figure 1.

The relative encounter velocity of planetesimals in the different zones of the disc between 15 and 30 au previous to the start of Neptune’s outward migration. The gradient in the excitation of the disc is due to the shorter orbital periods in the inner zones and to the distribution of spatial density of planetesimals. The relative velocities reflect this gradient.

In addition, the fraction of planetesimals present in each zone is shown in Fig. 2.

Figure 2.

The planetesimal disc starts at t = 0 with a surface density such that in each zone there are the same number of planetesimals. After reaching a dynamical equilibrium, the number of planetesimals redistributes in the way shown in the figure. These fractions were used to compute the number of triple systems formed by both proposed mechanisms.

We adopted the size distribution of planetesimals in the disc that is required to fit the observed size distribution of Jupiter trojans (Wong & Brown 2015; Yoshida & Terai 2017). Dynamical models suggest that Jupiter trojans probably formed in a massive disc between 20 and 30 au, and were subsequently implanted in their present location during the planetary instability phase of the outer planets (Leinhardt, Marcus & Stewart 2010; Leinhardt & Stewart 2012). The size distribution of Jupiter trojans is well characterized in our range of interest, from D ≃ 100–300 km (Wong & Brown 2015; Yoshida & Terai 2017). They have a cumulative size distribution characterized by a power law of the form N > D ∝ D^−γ, with γ ≃ 2.1 below D = 100 km, whereas at larger diameters, the size distribution is much steeper, with γ ≃ 6. As the probability that a planetesimal from the outer disc ends up into the Jupiter trojan clouds is of the order of 5 × 10⁻⁷ (Leinhardt et al. 2010; Leinhardt & Stewart 2012), ≃6 × 10⁹ planetesimals with D > 10 km have to be present in the primordial disc to match the present population of Jupiter trojans.

We also performed simulations exploring the possibility that the Lempo system originated from some of the proposed formation channels once the plutinos were established at their present location. For these simulations we used the intrinsic probability of collision and the distribution of relative velocities given by Dell’Oro et al. (2001).

The numerical integration of close encounters was carried out using an adaptive Bulirsch and Stoer routine (Press et al. 1992). We used an initial step size of 0.01 d, but if objects of different binaries approach closer than 10 times their combined radii, the step size was reduced so as not to miss possible collisions.

In the following, we will explain the details of the simulations for each formation channel.

2.1 Binary–binary collisions

In this case, we assumed that the collision should be inelastic and that both colliding objects remain stuck together after the collision. Detailed hydrodynamical simulations show that a projectile merges with a target of similar size when the impact velocity is ≤|$v$|_ESC, the mutual escape velocity (Leinhardt & Stewart 2012). Even for relative velocities, |$v$|_REL above |$v$|_ESC, merging occurs because the two bodies hit, losing kinetic energy, separate as nearly intact bodies, and then merge after a second collision (Leinhardt et al. 2010).

In a binary–binary collision, the impact velocity is difficult to anticipate, because the projectile and target are components of different binaries with their own orbital velocity. A body orbiting in a circular orbit has an orbital velocity of |$v_\mathrm{ESC}/\sqrt{2}$|⁠, and therefore, only due to the orbital motion, components of two different binaries may have a relative velocity at the moment of impact which may be of the same order as the mutual escape velocity. Taking these arguments into consideration, and due to the great computational expense of the simulations, we restrict ourselves to integration of only those encounters where |$v$|_REL ≤ |$v$|_ESC. We do not lose generality, because we also register as close encounters all those with |$v$|_REL > |$v$|_ESC: we compute the relative velocity; if |$v$|_REL < |$v$|_ESC then we proceed to the numerical integration of the close encounter; if not, we register it as a trial and compute a new relative velocity. Therefore, the only difference is that close encounters with |$v$|_REL > |$v$|_ESC are not numerically integrated. These high-velocity close encounters would not have given merging of bodies, but catastrophic collisions, which would not produce triple systems. As the mutual escape velocity is approximately given by

$$\begin{eqnarray*} v_\mathrm{ESC} \simeq \sqrt{\frac{2 G (M_\mathrm{TARG}+M_\mathrm{PROJ})}{R_\mathrm{TARG}+R_\mathrm{PROJ}}} \sim 0.01 \!-\!0.1 \,\, \mathrm{km} \,\, \mathrm{s}^{-1}, \end{eqnarray*}$$

(3)

the close encounters that matter are a very small fraction of the complete number of encounters because the distribution functions of relative velocities fall very fast with decreasing relative velocity. The number of trials varies for each zone because the relative velocity distributions are different.

The numerical integration stops if: (1) both binaries leave the Hill sphere intact; (2) one or both binaries are disrupted, either if e > 1 or if the semimajor axis exceeds 0.5R_H, which is a stability limit when considering the Hill stability criterion (Markellos & Roy 1981); and (3) a collision occurs; in this case, we must check if the collision would be potentially successful at forming a triple system. We define a successful collision as one where the total energy of the three remaining bodies after the collision is negative (Heggie 2005), and also by additional criteria based on the two-body orbits that we computed for the possible formed triple: we first compute the orbit of the two closer objects as a binary. If the orbit fulfils the condition of having a_B < 0.5R_H (a_B being the semimajor axis of the mutual orbit) then we compute the orbit of the third object around the centre of mass of the closer binary. This orbit must fulfil the same condition; in addition, we impose the requirement that neither orbit crosses. If all these conditions are satisfied, then we consider this system as a potential triple system. We stopped the simulation for each zone (1, 2, and 3) when 1000 successful collisions were recorded. For each zone, we performed ∼1.5 × 10⁶ integrations, and in total, for the entire disc, ∼1.2 × 10¹¹ trials were counted. For each collision, we computed the specific impact energy Q and the critical impact energy Q_D for catastrophic disruption (Benz & Asphaug 1999), and with them, the mass of the largest remnant was computed. In all cases, this mass was |$\gt 98{{\ \rm per\ cent}}$| of the original target mass. In all the 1000 successful cases the orbit of the three bodies was further integrated for 100 periods of the inner binary. As the model does not include important effects such as tidal friction and perturbations due to the non-spherical shape of the objects, we do not go further with this integration. As shown by Correia (2018), triple systems like Lempo’s have very unstable dynamics. Therefore it was not a surprise that most systems dissolved during this second integration. We are left with a small fraction of ∼100 stable triple systems.

2.2 Binary–binary exchange reaction

We have repeated almost the same simulations, but looking for possible close encounters where one binary is dissolved. For a binary to be dissolved by a close encounter, two basic conditions must be fulfilled: the relative velocity should be smaller than the orbital velocity of the binary, and the encounter’s closest approach distance should be comparable to or less than the binary separation. Both factors are important as the change in orbital energy after a close encounter increases with decreasing closest approach distance and relative velocity (Chauvineau & Farinella 1995; Funato et al. 2004). The impact parameter cannot be constrained because at such a low relative velocity, even for 100 km planetesimals, the gravitational enhancement of the cross-section, given by

$$\begin{eqnarray*} \sigma = \sigma _0 \left[ 1+\left(\frac{ v_\mathrm{ESC}}{v_\mathrm{REL}}\right)^2\right], \end{eqnarray*}$$

(4)

could be very large. Instead, for the relative velocity, we can proceed as in the collision simulations, restricting ourselves to starting the integration only for those close encounters where |$v$|_REL ≤ |$v$|_ESC. These cases occur for only a very small fraction of flybys. As in the previous channel, for the statistics, we counted all the initial conditions as a trial. In the simulations where both binaries survived the encounter, the integration was followed until they leave the Hill sphere. As the simulations for this channel are much more expensive than those exploring the BBC triple-formation channel, we decided to stop the simulations when 100 successful trials were recorded for each zone. This is when three objects remain bounded (negative energy and non-crossing orbits, as in the previous case) and one component leaves the Hill sphere alone. When this happens, the triple system is integrated additionally for 100 orbits of the closest system. Those that exhibit stable behaviour (checked by visual inspection) are then recorded as successful trials. In total, we have done ∼1.7 × 10¹¹ exchange-reaction simulations for the three zones, but only ∼2 × 10⁶ integrations were performed.

3 RESULTS

The main results of our simulations are shown in Table 1.

Table 1.

For each zone, we show its main characteristics: Rate of close encounters within the unit cross-section, mean relative encounter speed, and the number of planetesimals in the size range of interest. We also show the average number of close encounters needed to form one stable triple system, the number of stable triple systems found in the simulations, and the total number of stable triple systems expected to originate from each one of the mechanisms in T = 10⁸ yr for each formation mechanism.

Population	Pi	<|$v$|REL>	N	BBC			BBER
	km−2 yr−1	km s−1	50 ≤ D ≤ 300 km	Ne/Nsts	sim	total	Ne/Nsts	sim	total
15–20 au	1.8 × 10−20	0.78	∼5.5 × 107	2.6 × 109	32	47 642	1.8 × 1010	9	6991
20–25 au	8.9 × 10−21	0.44	∼7.9 × 107	8.3 × 108	35	225 391	9.3 × 108	11	201 764
25–30 au	7.3 × 10−21	0.24	∼7.0 × 107	3.6 × 108	33	1525 014	9.1 × 107	13	2028 925
Plutinos	4.4 × 10−22	1.44	∼5.6 × 105	1.3 × 1012	35	∼10−5	5.6 × 1011	10	∼10−5

In total, we have performed ∼1.7 × 10¹¹ close-encounter experiments for the BBER channel and ∼1.2 × 10¹¹ for the BBC channel.

Two examples of stable triple systems that we have found are shown in Figs 3(a) and (b). The main characteristics of them are depicted in Figs 4(a)–(d). We also show the corresponding values for the Lempo–Paha–Hiisi system.

Figure 3.

Upper panel: example of the process of formation of a triple system by collision of the components of two interacting binaries. The collision always gives rise to the larger component of the triple system. Lower panel: one component, which is usually the smaller one, is ejected during a close flyby, and a triple system is formed.

Fig. 4 shows that the orbital architecture and the size of the components of Lempo–Paha–Hiisi are reproduced quite well by the systems formed by both proposed mechanisms. In the BBC formation channel, one of the inner components tends to be larger than in the BBER case. This is because this component is formed by the merging of two objects. We also observe that, in some systems formed by the BBER channel, the external component is further away than in the systems formed by the BBC channel. Qualitatively, we explain this difference by the fact that the energy lost in an inelastic collision is less than the energy carried away by an escaping object. Therefore, the systems formed by the BBER channel are left with less orbital energy than those formed by the BBC channel. The relative inclinations of the orbital planes of the three components are quite random, reflecting the way in which we generated the initial condition. Several systems have almost coplanar orbits, like Lempo–Paha–Hiisi.

Figure 4.

Panels (a) (BBC) and (b) (BBER) show the mean separation between the components of the triple systems formed through each of the formation channels versus their diameters. The size of the circle is proportional to the diameter of the component. Open light-brown and light-green circles correspond to the inner binary components. Open cyan circles correspond to the outer component. The vertical green line indicates the centre of mass of the inner binary, from which the separation of the outer binary is measured. Panels (c) (BBC) and (d) (BBER) show the orbital distribution of the obtained systems for each formation channel. The colour of the circle is indicative of the relative inclination of the orbital plane of the outer binary with respect to the orbital plane of the inner binary. In all panels, the corresponding values of Lempo–Paha–Hiisi are also shown by blue open circles.

The number of stable triple systems generated per close encounter within the binary mutual Hill sphere (the Hill radius R_H is a measure of the distance at which the gravitational attraction of the binary and that of the Sun becomes comparable) depends strongly on the zone of the outer disc where the encounter happens. This is because the relative velocity among planetesimals decreases with increasing distance from the Sun. In the outer edge of the disc, near 30 au, 1 out of ∼3.4 × 10⁸ close encounters with impact parameter ≤R_H produces one stable triple system by BBER, whereas for BBC this number is ∼4.5 × 10⁸. As the integration was only carried out when |$v$|_REL ≤ |$v$|_ESC, the efficiencies in the generation of triple systems shown in Table 1 reflect mainly the difficulty of obtaining such a low relative velocity for the different zones. For example,

$$\begin{eqnarray*} \frac{N(v_\mathrm{REL}\le v_\mathrm{ESC})_\mathrm{ZONE\ 3}}{N(v_\mathrm{REL}\le v_\mathrm{ESC})_\mathrm{ZONE\ 1}}\sim 22, \end{eqnarray*}$$

(5)

in reasonable agreement with the efficiencies shown in Table 1. In total, ∼3300 triple systems were created in all our numerical simulations, but only 133 of them showed stable orbits.

4 CONCLUSIONS

We are now in position to compute the total number of triple systems created in the outer massive disc by these two channels, during the time T from just after the dispersal of the nebular gas up to the start of the instability phase of the outer planets. T is constrained by several observational arguments, including the survival of the Jupiter trojan binary Patroclus–Menoetius, to 100 Myr, but values up to 600 Myr cannot be ruled out (Bottke et al. 2012; Morbidelli et al. 2012; Marchi et al. 2013; Nesvorný et al. 2018). The size distribution of planetesimals in the outer disc, which is required to fit the observed size distribution of Jupiter trojans (Wong & Brown 2015; Yoshida & Terai 2017), implies the existence of N_p ∼ 2 × 10⁸ planetesimals with diameters between 50 ≤ D ≤ 300 km (Nesvorný & Vokrouhlický 2019), which is the range of diameters that results in triple systems like Lempo–Paha–Hiisi by BBC or BBER. We know that a fraction f_b of them were part of binary systems. Again, this fraction must be f_b ∼ 0.18–1 to account for the existence of the right number and characteristics of binaries in the Jupiter trojan clouds (Marchi et al. 2013), and the presence of a population of blue binaries in the cold classical Kuiper belt (Fraser et al. 2017). The number of created triple systems also depends on the average rate P_i at which close encounters within the unit cross-section occur in different zones of the disc (Morbidelli & Rickman 2015). The formation channels are not mutually exclusive, and there is no obvious reason for both not to work simultaneously. Considering all these factors, the number of triple systems that were formed during T = 100 Myr in the massive disc was

$$\begin{eqnarray*} N = (N_\mathrm{p} f_\mathrm{b})^2 \,\, R_\mathrm{H}^2 \,\, P_\mathrm{i} \,\, T \,\, N_\mathrm{sts}/N_\mathrm{e} \approx 4 \times 10^6 \times (f_\mathrm{b}/1)^2, \end{eqnarray*}$$

(6)

where N_e is the total number of simulated close encounters within a distance R_H and N_sts is the number of stable triple systems generated in these simulations. We expect N × 10⁻³ × 0.3 × f_b ≈ 200–1200 triple systems like Lempo–Paha–Hiisi residing today in the plutino population, depending on the value of f_b. The other two factors are the efficiency of implantation of planetesimals from the massive disc to the plutino population (Volk & Malhotra 2019) and the fraction of plutinos that have survived from the implantation epoch up to the present (Tiscareno & Malhotra 2009). An additional factor <1 must be applied, considering that a fraction of the triple systems formed is dynamically unstable on long time-scales (Correia 2018). However, this factor is hard to obtain because the number of parameters to explore is huge. An exhaustive exploration of the stability of orbits of trans-Neptunian triples would be important, giving us the remaining unknown factor and providing stronger constraints for the primordial environment where they form. Correia (2018) performed 100 simulations of the Lempo–Paha–Hiisi triple system and included tidal friction in the dynamical model, which in most cases is an orbital stabilizing factor. By only changing slightly the initial orbital eccentricity of the inner pair, Correia (2018) found that in 42 cases the system is circularized, while the remaining 58 systems become unstable. This result shows that the dynamical behaviour of triple systems is very complex and that the remaining unknown factor is surely much less than 0.5.

The 30–40 au low-mass disc extension is still poorly constrained (Nesvorný et al. 2020). We can only make some speculations, considering a mass fraction f_m (∼ 10 per cent) of the mass in the 15–30 au disc and a dynamical state similar to that of the outer zone of the massive disc below 30 au (e.g. relative velocities and the size distribution of planetesimals). P_i should be smaller because of the larger volume occupied by the planetesimals and the longer orbital periods. We suppose that P_i ∼ 4 × 10⁻²² km⁻² yr⁻¹, similar to the value in the cold classical Kuiper belt (Dell’Oro et al. 2001). With these values, we would have obtained ∼3 × 10⁴ triple systems. We do not yet know the implantation efficiency of this population (Nesvorný et al. 2020) nor the binary fraction there. Nevertheless, even if the outer-disc extension has fewer planetesimals than the massive disc below 30 au, they have more chances of ending up as part of one of the implanted populations in the Kuiper belt (Nesvorný et al. 2020). It is thus plausible that triple systems like Lempo–Paha–Hiisi could come from the 30–40 au region, a possibility supported by the very red colours of its three components (Benecchi et al. 2010).

We conclude that both formation mechanisms proposed here can explain the existence of triple systems like Lempo–Paha–Hiisi in the 3:2 resonant population and that a few more of them could still be found.

Manwë–Thorondor is a binary located at the 4:7 Neptune resonance, near 43.8 au (Nesvorný et al. 2011). Thorondor has an estimated diameter of 92 ± 14 km and orbits at a distance of 6674 ± 41 km from the primary Manwë, which is a barbell-shaped object composed of two elongated lobes of equal size (Nesvorný et al. 2011). Most of the triple systems formed by the mechanisms proposed here are hierarchical systems, composed of a tight inner binary with components of similar size and an outer distant companion. The tidal friction time-scale of the inner binary, assuming plausible values of the tidal parameters for icy bodies (Goldreich & Sari 2009), could be very short. This supports that Manwë–Thorondor was formed in the disc between 15 and 40 au by one of the formation channels previously proposed, and the tight inner binary evolved by tides very quickly, becoming a contact binary. All these processes happened when the system was still in its natal environment.

Additional numerical experiments of triple formation in the plutino region (at ∼39.4 au) show that, given the small number of objects in this population (∼10⁻³ of all the objects originally in the 15–30 au massive disc) and the high relative velocities among them (Dell’Oro et al. 2001), we cannot expect triple systems formed by either of both mechanisms there, or in any other resonant population. According to the same arguments, it is very unlikely that triple systems like Lempo–Paha–Hiisi or Manwë–Thorondor, formed by BBC or BBER, could be found in the cold classical Kuiper belt, even assuming that the initial mass in this region was ∼100–1000 times larger than its current mass, ∼3 × 10⁻⁴M_⊕ (Fraser et al. 2014). Nevertheless, we cannot rule out the existence of triple systems in the scattered-disc population.

An additional consequence of the mechanisms proposed here is that one of the components of triple systems formed by binary–binary collisions should have a bi-lobed shape, as a result of the low relative velocity impact giving rise to it (Jutzi & Benz 2017). This fact would provide a possible direct confirmation of the plausibility of the formation mechanisms proposed in this paper.

ACKNOWLEDGEMENTS

We acknowledge the financial support by CONICET and UNPA.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the corresponding author.