© ESO, 2015

1. Introduction

In a co-orbital system, two or more bodies share a mean orbit around a central body in a 1:1 orbital resonance. The possible orbits have been described by Lagrange in his Analytical Mechanics in 1788. We speak of tadpole orbits when one body librates around the Lagrange points L₄ or L₅ of another. Libration around L₃, L₄, and L₅ constitutes a horseshoe orbit (Brown 1911). Many trojan asteroids are in co-orbital trajectories with planets, but only in the Saturnian system are satellites known to be in co-orbital motion. While there are many examples of bodies with much smaller co-orbital partners, the Saturn moons Janus and Epimetheus are a unique system of two bodies with similar mass in horseshoe orbit. This co-orbital motion was confirmed in 1980 by the Voyager I probe (Synnott et al. 1981). Harrington & Seidelmann (1981) showed the stability of this system through numerical integrations including the oblateness of Saturn. Dermott & Murray (1981) developed an analytical study estimating the radial and angular amplitudes of each horseshoe orbit. A theory of motion was proposed two years later by Yoder et al. (1983), who explored the phase space of co-orbital systems with comparable masses, identifying the equilibrium points and showing the possible trajectory types (tadpole and horseshoe). Considering a system with up to nine equal-mass co-orbital bodies, Salo & Yoder (1988) found the location of the equilibrium points and their stability. Following a Hamiltonian approach, Sicardy & Dubois (2003) demonstrated the connection between co-orbital systems of comparable masses with that of the restricted case.

The formation of co-orbital systems has been the subject of intensive studies. In 1979, Yoder (1979) proposed three possible mechanisms for the formation of the Jovian trojans:

• 1.

  early capture during the accretion of Jupiter;

• 2.

  the continuous capture of asteroids; and

• 3.

  transfer of Jovian satellites to the Lagrange points.

In the first scenario, gas drag plays an important role. Peale (1993) and Murray (1994) showed that, in the presence of a gas, L₄ and L₅ can still be stable for proto-Jupiter. In 1995, Kary & Lissauer (1995) showed that a protoplanet on an eccentric orbit can indeed trap planetesimals in a 1:1 resonance, decaying from a large-libration horseshoe orbit to a tadpole orbit around L₅. Fleming & Hamilton (2000) explored the effect due to mass accretion and orbital migration on the dynamics of co-orbital systems. They showed how these perturbations affect the radial and angular amplitudes of the trajectories. Kortenkamp (2005) later proposed a mechanism to efficiently capture planetesimals on a satellite orbit after gravitational scattering by the protoplanet, making this a more likely fate than trapping in a co-orbital. Chanut et al. (2008) determined the size limit of trapping planetesimals via nebular gas drag to be 0.2 m per AU of orbit of the protoplanet. They (Chanut et al. 2013) later extended their approach to a non-uniform gas with a cavity in the vicinity of a heavy planet, but found little to no co-orbitals for eccentricities comparable to Janus and Epimetheus. All these scenarios result in co-orbital configurations with a large secondary body and a small partner, that is, with a low mass ratio. This almost excludes trapping as a mechanism for the formation of the Janus-Epimetheus system. Furthermore, Hadjidemetriou & Voyatzis (2011) showed that a system of two equal-mass bodies in a 1:1 resonance tends to evolve under gas drag into a closely bound binary satellite system. Another take on the subject by Cresswell & Nelson (2006) was the study of the evolution of a swarm of similarly sized planetesimals in a gas disk. They showed that the formation of mean-motion resonances, including 1:1 resonances, is commonplace, but that the bodies are still prone to be lost to accretion by the central body. Giuppone et al. (2012) followed a similar approach to Chanut et al. (2013), but studying a hypothetical extrasolar system. They showed that a pair of inwardly migrating planets can enter into a 1:1 resonance at the edge of a cavity in a gas disk. Whether their result is applicable to Janus and Epimetheus is not known to us.

Laughlin & Chambers (2002) showed that stable 1:1 resonances of equal-sized planets around a central body exist. They hypothesized that at the L₅ point of a planet in the accretion disk of a protostar, another planet could accrete. Beaugé et al. (2007) investigated formation by accretion of a swarm of planetesimals in a tadpole orbit around L₄ or L₅ of a giant extrasolar planet. They showed that an Earth-type body can form, but with a very low efficiency and thus a low final mass. Chiang & Lithwick (2005) considered this a likely process for the formation of the Neptune trojans. Izidoro et al. (2010) studied the mechanism for the Saturnian system and showed that formation of co-orbital satellite pairs with a relative mass of 10^-13 to 10^-9 is possible. This again rules out the formation of a co-orbital pair with similar masses.

While much research has been done on the subject of co-orbital systems, there is as yet no satisfactory explanation for the origin of the Janus-Epimetheus system in particular. Our approach here is based on the collision of two bodies at or close to the current orbit of Epimetheus. This is consistent with the equal densities (see Table 1) and compositions of Janus and Epimetheus. Moreover, results by Charnoz et al. (2010) show a recent formation of the Saturn moonlets from its rings through series of subsequent collisions. Crida & Charnoz (2012) later extended this scenario to all regular moons of Saturn. In this context, low-velocity collisions are frequent and are almost inevitable, especially for the predecessors of Janus and Epimemetheus. We assume the disruption into four main fragments that later collide again constructively to form a pair of bodies in a co-orbital configuration. We establish the statistical likelihood of this outcome in 8000 numerical simulations. We used Chambers’ Mercury (Chambers 1999) N-Body integrator with a Burlisch-Stoer algorithm.

In Sect. 2, we explain our model in detail and introduce the methods used to obtain our results. Section 3 shows the setup and results of our simulations. Finally, we discuss the viability of our model considering cases with six and eight fragments, and propose opportunities for future studies in Sect. 4.

2. Model and method

Fig. 1 Sketch of typical co-orbital formation scenario from our simulations.

We investigated whether the Janus/Epimetheus co-orbital system might have formed by disruption of a large body in the mean orbit of Epimetheus. Figures 1a, b illustrates the scenario that we believe can lead to co-orbital configurations. We assumed that four large fragments are formed by a collision event. While the initial collision is not subject of this study, we note that a low collision velocity and low density (compare Table 1) of the colliding bodies is a requirement. The collision is neither completely destructive nor is it constructive. After the collision event, fragments spread out on orbits around Saturn. Fragments can be lost by ejection from the system, collision with Saturn, or collision with each other. In the last case, fragments can unite to form a larger body, given a sufficiently low collision velocity (see Kortenkamp & Wetherill 2000).

Our simulations start right after the initial collision event, with the four fragments α, β, γ, and δ on the vertices of a square, moving in opposite directions with equal speed (see Fig. 2 for details). We always chose equal masses for bodies opposite on the square. The length a₀ of the sides of the square was chosen to be a₀ = 2R₁ + d₀ , where R₁ is the radius of the largest body and d₀ = 100m. The parameters varied in this study include the masses m_i of the fragments (given as fractions of the mass of Janus and Epimetheus, see Tables 1 and 3), the orientation of the square in the orbital plane denoted by the angle ϕ (going from 0 to π because of symmetry), and the initial speed v_init (see Fig. 2) of the fragments in the coordinate system moving with v₀ relative to Saturn. The position and velocity { r₀,v₀ } of the centre of the square are equal to the state vector of Epimetheus at an arbitrarily chosen point of its orbit. This is the same for all simulations.

We used the package Mercury to numerically integrate the system using a Burlisch-Stoer algorithm (Chambers 1999). We included the moons Mimas and Tethys into our simulations, since they exert the strongest influence on the Janus-Epimetheus system, but neglected other objects. We also considered the oblateness of Saturn through J₂, J₄, and J₆ (Table 2). To obtain the state vectors of Mimas, Tethys, and Epimetheus (centre of the square from which the fragments start), we used the data given by Jacobson et al. (2008) and followed the approach presented by Renner & Sicardy (2006) to transform the geometric elements. We took the mean distance between the system Janus-Epimetheus and Saturn as the unit of length R₀. The simulation code was set to record the state vectors of all bodies every 1 / 100 day. We need this high time resolution to study the velocities and angles of collisional events between fragments. Mercury assumes a perfectly inelastic collision where the colliding bodies aggregate into a new body of the added masses of the colliding bodies, and the momentum is preserved. We therefore needed to analyse the collisions to ensure that this assumption is realistic. The outcome of a collision is strongly dependent on the collision velocity and the physical characteristics (density, strength, etc.) of each body. Kortenkamp & Wetherill (2000) established an escape velocity term ${\mathit{v}}_{\max }^{\mathrm{2}}\mathrm{=}\frac{\mathrm{2}\mathit{G}\mathrm{(}{\mathit{m}}_{\mathrm{1}}\mathrm{+}{\mathit{m}}_{\mathrm{2}}\mathrm{)}}{{\mathit{R}}_{\mathrm{1}}\mathrm{+}{\mathit{R}}_{\mathrm{2}}}$(1)for the encounter of rocky bodies, where m_i is the mass of body i, R_i is its radius, and G is the gravitational constant. For collision velocities below v_max, a collision can be assumed to be constructive. For a collision of Janus and Epimetheus, v_max ≈ 58 m/s. However, their study involved bodies with a density about five times as high as that of Janus and Epimetheus. Since in a collision event between fluffy bodies (such as Janus and Epimetheus, whose density is lower than 700 kg / m³) kinetic energy is lost to deformation, higher velocity collision events can still be constructive. High-velocity collisions can also result in the catastrophic disruption of one or both colliding bodies. Multiple studies on the disruption threshold (Colwell 1994; Benz & Asphaug 1999; Leinhardt & Stewart 2009) have been conducted. Using Fig. 11 from Leinhardt & Stewart (2009) and a target radius of 6 × 10⁴ m (the rough size of Epimetheus), we obtain a disruption threshold of ${\mathit{Q}}_{\mathrm{D}}^{\mathrm{\ast }}\mathrm{\approx }\mathrm{6.5}\mathrm{\times }{\mathrm{10}}^{\mathrm{3}}$J/kg. This gives a cutoff velocity ${\mathit{v}}_{\mathrm{D}}\mathrm{=}\sqrt{\mathrm{2}{\mathit{Q}}_{\mathrm{D}}^{\mathrm{\ast }}\frac{{\mathit{m}}_{\mathrm{T}}}{{\mathit{m}}_{\mathrm{P}}}}$(2)(where m_T is the mass of the target and m_P the mass of the projectile) of approximately 60 m/s for a collision with Janus. Therefore we used the first estimate of 58 m/s as a cutoff, keeping in mind that the actual threshold might still vary and is in fact dependent on which of the fragments collide.

Fig. 2 Initial configuration of the simulations. ${\begin{array}{c}{\mathrm{˜}}_{}\\ {v}\end{array}}_{\mathrm{i}}$ are given in the reference frame moving with v0. The initial velocity of each fragment in the simulation frame of reference is given by vi = ${\begin{array}{c}{}_{\mathrm{˜}}\\ {v}\end{array}}_{\mathrm{i}}\mathrm{+}{{v}}_{\mathrm{0}}$. The direction of ${\begin{array}{c}{\mathrm{˜}}_{}\\ {v}\end{array}}_{\mathrm{i}}$ is given by the vector from the centre of the square to the respective vertex. The collision fragments α, β, γ, and δ are positioned on the vertices of a square centred on r0, the position vector with respect to Saturn. The angle ϕ gives the rotation of the square around the vector r0 × v0 (the normal of the orbital plane). For ϕ = 0α, γ, and Saturn are connected by a straight line, with γ between α and Saturn. The fragments are always set on the vertices in a counter-clockwise sense, as shown in the sketch. The dimensions here are not to scale.

Because of the large number of simulations, we needed an automatic filtering algorithm to find the systems that evolved into a co-orbital. We used a four-step filtering process.

• 1.

  We selected the systems in which at least two fragments havesurvived until the end of the simulation (20–30 yr). Fragments canbe lost through ejection, collision with Saturn, or collision witheach other.

• 2.

  From the remaining systems, we selected the systems with exactly one pair of fragments that survived until the end of the simulation and overlap in their range of semi-major axes for the time after the last collision in the simulation.

• 3.

  From these, we selected the systems in which the mean semi-major axis of the surviving fragments is smaller than 1.2 R₀ and the minima and maxima of the semimajor axes are no more than 0.02 R₀ apart.

• 4.

  We plot the semi-major axes of the remaining fragments against time and inspected them visually to check for periodic crossings. We thus determined false positives of the automatic filtering process and removed them.

3. Simulations and results

Fig. 3 Cumulative distribution of collision velocities for co-orbital systems in batch 3.

In each simulation, we sampled the parameters v_init and ϕ in equally spaced steps. We started our investigation with a batch of 4000 simulations, covering 20 different initial speeds v_init and 200 different angles ϕ (see Table 3–1). The initial speeds were between 10 and 100 times as high as the escape velocity of Janus. We found that around 1.5% of the simulations evolved into co-orbital systems. From the analysis of these systems, we were able to establish a common formation scenario as depicted in Fig. 1: after the initial collision event, the fragments first spread out on orbits around Saturn. Then, they re-collide pairwise, most commonly α with γ and β with δ. The velocities of these collision events are very high (1000 m/s and above), far beyond where a perfectly inelastic joining, as is assumed by Mercury, is realistic. But the moons created in this way quickly (in just a few orbits) enter into co-orbital horseshoe motion after the second moon has been formed. No co-orbital systems are the result of ejection of fragments or collision with Saturn.

Then, we investigated whether lower initial velocities (Table 3–2) result in lower collision velocities. We reduced the parameter range of v_init by a full order of magnitude. This indeed reduced the collision velocities by nearly an order of magnitude, but not to the point where simple merging of bodies was a realistic outcome for any of the collisions. We note, however, that co-orbital formation increased to 5.7%.

Our third batch of simulations was divided into three parts: first, we simulated with the same masses as before, but at even lower initial velocities (Table 3–3a). As expected, the collision velocities decrease even more (see Fig. 3 for a comparison of collision velocities). The portion of co-orbitals increases to 8.4%. Another result is that the frequency of co-orbital formation decreases for very low initial velocities (below ~0.015 R₀/ d or around half of the Janus escape velocity). The highest observed frequency of co-orbitals is around 0.025 R₀/ d.

Fig. 4 Frequency of co-orbitals as a function of the initial speed vinit of the fragments for batch 3.

Fig. 5 Period of the horseshoe orbits. Note the similar distribution for all masses versus selected masses. The masses of the bodies have almost no influence on the period.

For the second and third part of this batch, we decreased the masses of the fragments by a half. One part was ran at initial velocities like in the previous batches (Table 3–3b) and the other at lower ones (Table 3–3c). By comparing simulations 3a and 3c, we see that the mass of the fragments has little effect on the distribution of generated co-orbital systems. Only the decrease in occurrence of co-orbitals for the lowest initial velocities (see Fig. 4) is less prominent for lower mass. The decrease at these initial velocities indicates a lower limit for the formation of co-orbitals. Supplementary simulations indicate that for even lower velocities v_init< 2 × 10^-5R₀/ d, formation of co-orbitals is quite unlikely (~0.7%) compared to the other simulations. A cut-off would probably depend on the mass of the fragments, assuming it is related to the escape velocity of the fragments. This could explain the higher formation probability for lower-mass fragments at low initial velocities.

Fig. 6 Four classes of co-orbitals by mass ratio and $\begin{array}{c}{\mathrm{˜}}_{}\\ \mathit{A}\end{array}$ (see Table 4) The variables are inversely proportional for horseshoe orbits (Dermott & Murray 1981, Eq. (1)). The red circle indicates the only pair in a tadpole orbit we observed in this batch of simulations.

We thus have two constraints on the initial conditions if we wish to determine a range likely to spawn co-orbitals: the collision velocities are highly dependent on the initial velocities. Low collision velocities are required for the fragments to collide constructively, thus low initial velocities. At the same time, the formation frequency of co-orbitals decreases after a certain point, so the initial velocity cannot be arbitrarily low. The existence and width of a sweet spot for the formation of co-orbitals is key to decide the viability of this model.

Fig. 7 Example orbits from all four classes. Radial axis is the semi-major axis relative to the mean semi-major axis. Angular is the mean longitude in a frame rotating with the mean motion associated with the mean semi-major axis. To view the evolution of each of these cases, see the associated movies available online.

We examined the results of simulation batch 3c more closely. In many but not all cases, the fragments collided in such a way that the final mass of the bodies was equal to Janus and Epimetheus. This is the case for 36 of the 87 co-orbital systems we obtained in this batch, or 3.6% of all simulations in the batch (see Table 4). The collision velocities reach up to about 450 m/s. The collisions occur at angles smaller than 5° for nearly all cases. Assuming a limit of 58 m/s for constructive collision, we find that 15.8% of all collision events in this batch remain (see Fig. 3). By multiplying the original success rate of 3.6% by 15.8% squared (since both collision events need to be below the threshold), we obtain a corrected success rate of 0.09% to generate Janus and Epimetheus or 0.2% to generate any co-orbital.

We compared the libration period of the horseshoe orbits to the period of the horseshoe orbit of Janus and Epimetheus (see Fig. 5). Janus and Epimetheus cross orbits roughly every four years. This period is highly variable for different simulations and shows no particular connection to the mass of the moons. Instead, there seems to be a general preference for short periods.

We see more influence of the mass of the fragments on the trajectories by examining the relative range of the semi-major axes. With our setup of four fragments, two each of identical mass, there are four possible combinations of a final co-orbital pair, listed in Table 4.

We now consider the ratio of radial amplitudes $\begin{array}{c}\mathrm{˜}\\ \mathit{A}\end{array}\mathrm{\equiv }\frac{\max \mathrm{\left(}{\mathit{a}}_{\mathrm{1}}\mathrm{\right)}\mathrm{-}\min \mathrm{\left(}{\mathit{a}}_{\mathrm{1}}\mathrm{\right)}}{\max \mathrm{\left(}{\mathit{a}}_{\mathrm{2}}\mathrm{\right)}\mathrm{-}\min \mathrm{\left(}{\mathit{a}}_{\mathrm{2}}\mathrm{\right)}}\mathit{,}$(3)where a₁ is the set of semi-major axes of the heavier body of the final co-orbital pair after the last collision event in the simulation, and a₂ is the same set for the lighter body. For the Janus-Epimetheus system, this is $\begin{array}{c}\mathrm{˜}\\ {\mathit{A}}_{\mathrm{J}\mathit{/}\mathrm{E}}\end{array}\mathrm{=}\mathrm{0.2102.}$(4)By plotting $\begin{array}{c}{\mathrm{˜}}_{}\\ \mathit{A}\end{array}$ against the ratio of masses (see Fig. 6), we recover the expected inversely linear correlation between the two. Each class has similar ratio of radial amplitudes. There is one simulation with significant deviance (see Fig. 6, red circle) from this law. Closer investigation reveals that in this case, the surviving bodies enter into a tadpole orbit (see Fig. 8). We have found no other case of tadpole orbits in our simulations.

Fig. 8 Single case of tadpole orbit. The fragments have the mass of Janus and Epimetheus. The plot shows the semi-major axis relative to the mean semi-major axis plotted against the mean longitude in a frame rotating with the mean motion of the mean semi-major axis. (Online movie.)

4. Conclusion and outlook

We here investigated the possibility of a formation process of the Janus-Epimetheus system starting from a collision event that spawns four large fragments. Additional debris was neglected. We used numerical integration of the N-body problem to obtain the evolution of 8000 initial configurations. We selected the configurations that led to the formation of a co-orbital system. By analysing these, we established a common scenario that leads to the formation of co-orbital systems: after the fragments spread out on orbit, they re-collide and form a pair of co-orbiting bodies.

We investigated the circumstances of these secondary collisions and showed that the initial velocity has a high influence on the collision velocities. At the same time, we found that lower initial velocities yield a higher probability of co-orbital formation. This reinforces our assertion that for our assumed initial condition to be achieved, the prior collision event must have occurred at a low velocity. It is unclear, however, whether the collision velocities observed in these simulations are low enough that a perfectly inelastic collision and joining of the bodies is a good approximation. We assume that relative velocities below 58 m/s allow the fragments to collide constructively. Previous studies on the collisions of celestial bodies are not fully applicable to our system. Most works dealt either with collisions of rocky bodies and/or of a small projectile hitting a large target. Therefore it is necessary to study the collision behaviour of large low-density fluffy objects such as Janus and Epimetheus in more detail. Studies in this area would help to ultimately decide the viability of our formation model. We also found that there is a lower cut-off below which the formation of co-orbitals decreases significantly. We attribute this to the inability of the fragments to escape their combined gravity when the kinetic energy is too low. This imposes two constraints on the initial velocity: it must be low enough to obtain reasonably low collision velocities. At the same time, it must be high enough to allow the fragments to escape their mutual gravitational pull. Our results indicate that this sweet spot is quite narrow, but indeed exists.

Another limiting feature of our physical model is the number of initial fragments adopted in the simulations. To strengthen the scenario proposed here, we performed two sets of simulations, one considering six, the other eight initial fragments. Considering initial velocities as those of simulation 3c, 100 simulations were performed for each case (six and eight fragments). The results showed an even higher rate of co-orbital generation (18–20%). Therefore, the scenario proposed here is somewhat robust.

Furthermore, we are aware that our initial conditions are not very natural. Instead, we chose a simple configuration to constrain the otherwise quite large number of free parameters. The main result of this study is therefore that the proposed formation model is possible. We are confident to find co-orbital formation for a wider range of initial conditions building upon this result. Ideally, a future work would include the initial collision as part of the simulation and use a better approximation for the secondary collision events.

Online material

Movie of Fig. 7 Access here

Movie of Fig. 7 Access here

Movie of Fig. 7 Access here

Movie of Fig. 7 Access here

Movie of Fig. 8 Access here

Acknowledgments

We thank Rafael Sfair for providing file storage for our results on short notice. Othon C. Winter thanks FABESP (proc. 2011/08171-3) and CNPq. The constructive comments and suggestions of the referee significantly improved this paper.

References

All Tables

All Figures

	Fig. 1 Sketch of typical co-orbital formation scenario from our simulations.
In the text

Fig. 2 Initial configuration of the simulations. ${\begin{array}{c}{\mathrm{˜}}_{}\\ {v}\end{array}}_{\mathrm{i}}$ are given in the reference frame moving with v0. The initial velocity of each fragment in the simulation frame of reference is given by vi = ${\begin{array}{c}{}_{\mathrm{˜}}\\ {v}\end{array}}_{\mathrm{i}}\mathrm{+}{{v}}_{\mathrm{0}}$. The direction of ${\begin{array}{c}{\mathrm{˜}}_{}\\ {v}\end{array}}_{\mathrm{i}}$ is given by the vector from the centre of the square to the respective vertex. The collision fragments α, β, γ, and δ are positioned on the vertices of a square centred on r0, the position vector with respect to Saturn. The angle ϕ gives the rotation of the square around the vector r0 × v0 (the normal of the orbital plane). For ϕ = 0α, γ, and Saturn are connected by a straight line, with γ between α and Saturn. The fragments are always set on the vertices in a counter-clockwise sense, as shown in the sketch. The dimensions here are not to scale.

In the text

	Fig. 3 Cumulative distribution of collision velocities for co-orbital systems in batch 3.
In the text

	Fig. 4 Frequency of co-orbitals as a function of the initial speed vinit of the fragments for batch 3.
In the text

	Fig. 5 Period of the horseshoe orbits. Note the similar distribution for all masses versus selected masses. The masses of the bodies have almost no influence on the period.
In the text

	Fig. 6 Four classes of co-orbitals by mass ratio and $\begin{array}{c}{\mathrm{˜}}_{}\\ \mathit{A}\end{array}$ (see Table 4) The variables are inversely proportional for horseshoe orbits (Dermott & Murray 1981, Eq. (1)). The red circle indicates the only pair in a tadpole orbit we observed in this batch of simulations.
In the text

	Fig. 7 Example orbits from all four classes. Radial axis is the semi-major axis relative to the mean semi-major axis. Angular is the mean longitude in a frame rotating with the mean motion associated with the mean semi-major axis. To view the evolution of each of these cases, see the associated movies available online.
In the text

	Fig. 8 Single case of tadpole orbit. The fragments have the mass of Janus and Epimetheus. The plot shows the semi-major axis relative to the mean semi-major axis plotted against the mean longitude in a frame rotating with the mean motion of the mean semi-major axis. (Online movie.)
In the text