A&A 421, 771-774 (2004)
DOI: 10.1051/0004-6361:20047025

Numerical escape criteria for a symmetric four-body model

A. Széll ¹ - B. Steves ¹ - B. Érdi ²

1 - Glasgow Caledonian University, Cowcaddens Road, G4 0BA Glasgow, UK
2 - Eötvös Loránd University, Pázmány Péter Sétány, Budapest, Hungary

Received 7 January 2004 / Accepted 7 April 2004

Abstract
Numerical escape criteria is presented for the Caledonian Symmetric Four-Body problem (CSFBP). The numerical experiments show that escapes can be detected very early with the help of the method. Integrating a huge amount of orbits of the symmetric four-body system we found that for the equal mass case the double binary escape, and in the planetary case the single bodies escape are the most likely outcome of the disintegration of the system.

Key words: celestial mechanics

1 Introduction

The determination of whether a body escapes a system or not is a very difficult problem. There are no known analytical methods for the general n-body problem (n>2) which detect escapes. However, it is possible to determine numerical criteria for the detection of escapes. Numerically these criteria can be verified, but as of yet the validity of such criteria has not been proven analytically.

In the paper of Shebalin & Tippens (1996) a numerical escape criterion is given for the general problem of three bodies P₁, P₂, and P₃. They divided the system into three binary subsystems, ie. three sets of "composite binaries", or binaries P₁P₂, P₁P₃, and P₂P₃. They defined three energy-like parameters, corresponding to the energies of the three binary subsystems. They found that these parameters correlate with an escape of the system, or in other words these parameters can predict escape. Plotting the three energy-like parameters as functions of time shows that they fluctuate below and sometimes above zero. Shebalin & Tippens (1996) found that if all the parameters become greater than zero at the same time, then the system would break up in a finite time.

In our present work we use their idea in order to find escape criteria for a symmetric four-body problem called the Caledonian Symmetric Four Body Problem (CSFBP) first defined by Steves & Roy (1998, 2000, 2001). In their papers they investigated the hierarchical stability and evolution of the CSFBP using an analytical stability criterion which they showed depends solely on a parameter they called the Szebehely constant C₀, where C₀is a function only of the total energy and angular momentum of the system.

The hierarchical stability of the CSFBP was also investigated numerically by Sz $\acute{\rm e}$ll et al. (Sz $\acute{\rm e}$ll 2003), in which they confirmed the analytical criteria and explored in detail the dominant hierarchical configurations for different mass ratios $\mu$ between the bodies. The CSFBP has two dynamical symmetric pairs of mass m and Mwith mass ratio $\mu = \frac{m}{M}$. As part of their study, they followed the evolution of about 70 000 different CSFBP orbits having 23 different pairs of $(\mu, C_0)$ where $\mu = 1, 0.1,0.01$ and 0.001.

During these orbits many escapes of the system bodies occurred. It became a valuable saving of the computer time to be able to recognize those systems which would result in escapes long before they had been numerically integrated to break up point.

The aim of this paper is to derive such an early detection method for escapes in the CSFBP and to study the dominant characteristics of escapes in the symmetrical four body problem.

At first we define the four-body configuration. Then we introduce energy-like parameters that can be used to detect escapes. We apply these criteria to investigate the distribution of the different types of escape configurations in the symmetric four-body problem.

2 The symmetric four-body model

	Figure 1: The symmetric four-body configuration at any time t.
Open with DEXTER

Let us consider four mass-points P₁, P₂, P₃, P₄ in the two dimensional Euclidean space with position vectors $\vec{R}_i( \ne \vec{0} )$, and velocity vectors $\dot{\vec{R}_i}$, i=1,2,3,4. Let the origin be at the centre of mass C of the whole system and let it be at rest. Let M₁, M₂, M₃, M₄ be the masses of the mass-points. Let us introduce the following symmetry conditions

 

$$\vec{R}_1 = -\vec{R}_4 = \vec{r}_1, \quad \vec{R}_2 = -\vec{R}_3 = \vec{r}_2 ,$$

$$M_1 = M_4 = m , \quad M_2 = M_3 = M$$ (1)

for all t.

It is straightforward from Eq. (1) that

$$\begin{eqnarray*}\dot{\vec{R}}_1 = -\dot{\vec{R}}_4 = \dot{\vec{r}}_1, \quad \dot{\vec{R}}_2 = -\dot{\vec{R}}_3 = \dot{\vec{r}}_2 \end{eqnarray*}$$

for all time t.

We shall use the following notations:

• let $\mu = m/M$;
• let c₁₂, c₃₄, c₁₃ and c₂₄ denote the baricenter of the bodies P₁ and P₂, the baricenter of the bodies P₃ and P₄, the baricenter of the bodies P₁ and P₃ and the baricenter of the bodies P₂ and P₄, respectively.

The system has essentially four different kinds of escape configurations. They are

• "12'' escape: a double binary escape, where P₁ and P₂ form a binary and escape, while symmetrically P₃ and P₄ also form a binary and escape;

• "13'' escape: a double binary escape, where P₁ and P₃ form a binary and escape, while symmetrically P₂ and P₄ also form a binary and escape;

• "14'' escape: a two single bodies escape, where P₂ and P₃ form an inner binary and P₁ escapes, while symmetrically P₄ also escapes;

• "23'' escape: a two single bodies escape, where P₁ and P₄ form an inner binary and P₂ escapes, while symmetrically P₃ also escapes.

   
3 The escape criteria

	Figure 2: Typical time evolution of the E1, E2, E3, and E4 parameters in the case of a) "12'' type of escape (two binaries escaping) and b) "13'' type of escape (two binaries escaping).
Open with DEXTER

	Figure 3: Typical time evolution of the E1, E2, E3 and E4 parameters in the case of a) "14'' type of escape (two singles escaping, leaving a binary) and b) "23'' type of escape (two singles escaping, leaving a binary).
Open with DEXTER

In this section we define four new, energy-like parameters: E₁, E₂, E₃, and E₄, that can be used to detect escapes.

   
3.1 Detection of a "12'' escape using E₁

We recall, that in the Newtonian two-body problem the energy of the system can be written as

$$E = T + U = \frac{1}{2} m_1 \vert\vec v_1\vert^2 + \frac{1}{2} m_2 \vert\vec v_2\vert^2 - G \frac{m_1 m_2}{r_{12}^2} ,$$ (2)

where T is the kinetic energy, U is the potential energy, m₁ and m₂ are the masses, $\vec v_1$ and $\vec v_2$ are the velocities of the bodies, and r₁₂ is the distance between the bodies. G is the gravitational constant which with appropriate choice of time units can be set to 1. If we substitute $m_1 = m_2 = 1 + \mu$, $\vert\vec v_1\vert = \vert\vec v_2\vert = v_{c_{12}}$, and r₁₂ = r_c12c34, then Eq. (2) gives Eq. (3):

$$E_{1} = (1 + \mu) v_{c_{12}}^{2} - \frac{(1+\mu)^{2}}{r_{c_{12}c_{34}}} ,$$ (3)

where v_c12 is the absolute value of the velocity of the centre of mass c₁₂ of the bodies P₁ and P₂. Due to the symmetry, v_c34=v_c12, where v_c34 is the absolute value of the velocity of the centre of mass c₃₄ of the bodies P₃ and P₄. r_c12c34 is the distance between c₁₂ and c₃₄.

Thus E₁ can be interpreted as the formal energy expression of a two-body system that is formed by a point mass located at c₁₂ with mass $1+\mu$, and another point mass at c₃₄ with the same mass. But in the present case E₁ is not a constant with respect to time t. E₁ = E₁(t), since Eq. (3) is not an integral of motion of the CSFBP (the variations of v_c12 and r_c12c34 can be determined from the equations of motion).

Still, (3) can be used to predict the disruption of the system. If the energy E₁ at a given time was positive, then the two-body system corresponding to Eq. (3) at that time would break up. For increasing positive energies E₁, the points c₁₂ and c₃₄ would diverge from each other. According to our numerical experiments in the case of a "12'' type of escape, i.e when the system breaks into two binaries, with the P₁ P₂ and the P₃ P₄ binaries diverging from each other, the function E₁(t) behaves typically as shown in Fig. 2a. It can be seen that E₁(t) oscillates with positive minimum values that continue to increase.

3.2 Detection of a "13'' escape using E₂

Let E₂ be

$$E_{2} = (1 + \mu) v_{c_{13}}^{2} - \frac{(1+\mu)^{2}}{r_{c_{13}c_{24}}} ,$$ (4)

where v_c13 is the absolute value of the velocity of the centre of mass c₁₃ of the bodies P₁ and P₃. Due to the symmetry, v_c24=v_c13, where v_c24 is the absolute value of the velocity of the centre of mass c₂₄ of the bodies P₂ and P₄. r_c13c24 is the distance between c₁₃ and c₂₄.

Note, that E₂ is very similar to E₁ in nature, only the binaries studied are now P₁ P₃ and P₂ P₄, instead of the previous P₁ P₂ and P₃ P₄ cases. The arguments given in Sect. 3.1 can also be applied here. The E₂(t) curve can therefore be used to detect "13'' type of escapes.

In Fig. 2b the E₂(t) function can be seen for a "13'' type of escape. It can be seen that the minimum values of E₂ are positive and continue to increase. When this behaviour is found the eventually system breaks up exhibiting a "13'' escape.

   
3.3 Detection of a "14'' escape using E₃

Let E₃ be

$$E_{3} = \sqrt{8 \mu + \mu^2} v_{a}^{2} - \frac{8 \mu + \mu^2}... ...mu v_{1}^{2} - \frac{4\mu}{r_{1}} - \frac{\mu^{2}}{2r_{1}} ,$$ (5)

where v₁ is the absolute value of the velocity of the body P₁, and due to the symmetry, v₄=v₁ is the absolute value of the velocity of the body P₄, and r₁ is the distance of the body P₁ from the centre of mass of the four-body system. In order to gain the formal energy formula of the two-body problem we introduce $v_a = \sqrt[4]{\frac{\mu}{8+\mu}} v_1$. If we substitute $m_1 = m_2 = \sqrt{8 \mu + \mu^2}$, $\vert\vec v_1\vert = \vert\vec v_2\vert = v_a$ and r₁₂ = 2 r₁ into Eq. (2), then it gives Eq. (5).

Thus E₃ can be interpreted as the formal energy of a two-body system that is formed by a point mass located at P₁ with mass $\sqrt{8 \mu + \mu^2}$, and another point mass at P₄ with the same mass. But again in this case E₃ is not a constant with respect to time t. E₃ = E₃(t), since v_a and r₁ can be determined from the equations of motion to which Eq. (5) is not an integral.

If the bodies P₁ and P₄ were placed very far from each other, then E₃ would tend to be a constant, since the system could be approximated by a two-body system formed by P₁ and P₄ with masses $\sqrt{8 \mu + \mu^2}$. This is the reason why E₃ is defined as the energy-like variable belonging to the "14'' type of escape.

Now we can repeat the train of thought we used for the "12'' type hierarchy case. If the energy E₃(t) at a given time was positive, then the two-body system corresponding to Eq. (5) at that time would break up. For increasing positive energy P₁ and P₄ would diverge from each other.

According to our numerical experiments in the case of a "14'' type of escape, the E₃(t) function behaves typically as shown in Fig. 3a. It can be seen that the E₃(t) function is positive, first increasing and then almost constant indicating escape has occurred.

3.4 Detection of a "23'' escape using E₄

Let E₄ be

$$E_{4} = \sqrt{8 \mu + 1} v_{b}^{2} - \frac{8 \mu + 1}{2 r_{2}} = v_{2}^{2} - \frac{4\mu}{r_{2}} - \frac{1}{2r_{2}} ,$$ (6)

where v₂ is the absolute value of the velocity of the body P₂, and due to the symmetry, v₃=v₂ is the absolute value of the velocity of the body P₃, and r₂ is the distance of the body P₂ from the centre of mass of the four-body system. In order to gain the formal energy formula of the two-body problem we introduce $v_b = \sqrt[4]{\frac{1}{8\mu+1}} v_2$.

Note, that E₄ is very similar to E₃ in nature, only the central binary pair is now P₁ and P₄, instead of the previous pair P₂ and P₃. The arguments given in Sect. 3.3 can also be applied here. The E₄(t) curve can therefore be used to detect "23'' type of escapes.

In Fig. 3b the E₄(t) function can be seen for a "23'' type of escape. It can be seen that the E₄(t) function is positive and almost constant in time after a certain limit.

4 Numerical integrations

We examined in the aggregate about 15 000 orbits which ended in escapes, for four different values of $\mu$, for $\mu = 1$, 0.1, 0.01, and 0.001, in order to investigate the distribution of the different kinds of escape configurations. For greater $\mu$values the system can be a model for quadruple stellar systems, while for smaller $\mu$ values the system can be a model for planetary systems.

During the numerical studies we chose the initial conditions to satisfy the symmetry conditions of the CSFBP and initially the bodies were collinear.

We detected the escapes with the help of our escape criteria. However, the integration was not stopped when an escape was detected. It was always executed for 4 000 000 time steps. Since the step-size was set to 0.001, the result was a numerical integration time of 4 000 time units. Since the gravitational constant G and the mass of the more massive bodies were set to one, one time unit was about 368 years (Széll 2003). Thus the total integration time was about 1 500 000 years. We found that in each case, when the escape criteria observed an escape, at the end of the integration time two pairs or two of the bodies were always located considerably far from the baricentre of the system.

The results of the integrations can be seen in Table 1. Physically there is no difference between the "12'' and "13'' type of escapes. They are both double binary escapers. In the case of "14'' and "23'' type of escapes there are two single bodies escapers. In the $\mu = 1$ equal masses case they are similar. For $\mu < 1$ there are two possibilities: the two smaller ("14'' type) or the two greater ("23'' type) bodies can escape. The numerical integrations showed that always the two smaller bodies escape from the system. Thus in Table 1 we tabulate only the two different kind of outcomes: double binary escapes when the system falls apart into two escaping binaries, and two single escapers when the two smaller single bodies escape ("14'' type of escapes).

Table 1: Percentage of the total number of escape configuration.

The table shows that as the parameter $\mu$ decreases, the single body escapes become dominant. For $\mu = 1$ the most likely escape is the double binary escape. Thus four-star systems break up most likely into two binary systems. In the $\mu = 0.1$ case the escape of 2 single bodies escape is dominant, but the double binary escapes are still considerable. For small values of $\mu$ only the single escapers are considerable. Thus for planetary systems, the escape of the planets is the most likely outcome.

5 Conclusions

In this paper we presented numerical escape criteria for the Caledonian Symmetric Four-Body Problem. By approximating the symmetrical four body problem on the verge of break up as a two body problem, we were able to derive four energy like parameters E₁, E₂, E₃ and E₄ which were related to each of the four types of possible escapes. A study of the energy like parameters for 15 000 orbits integrated to escapes showed that the energy like parameters exhibit distinctive behaviour when the system is about to break up. this behaviour can be used as an escape criteria to predict escapes of the CSFBP. For example when either E₃ or E₄ becomes a positive value increasing to a constant, this indicates that a "14'' or "23'' type of escape respectively will occur. When either E₁ or E₂ has a consecutive minimum values which are positive and increasing and the energy parameters have some negative value, then this indicate that a "12'' or "13'' type of escape respectively will occur. The CSFBP symmetry conditions were employed in the derivation of energy parameters, but it is likely that the same train of thought can be applied for arbitrary systems. This will be a subject of future research.

Numerous orbits were investigated with the help of the numerical escape criteria. We found that for stellar systems the most likely escape configuration is the double binary escape, i.e. the system falls apart into two binary-star configurations. For small values of $\mu$, the system can be a model for two stars two planets systems. The integrations show that the most likely escape is the escape of the two planets.

Acknowledgements

This work was partly supported by the British-Hungarian Scientific and Technological Foundation under grant number GB-53/01.

References

• McMillan, S. L. W., & Hut, P. 1996, ApJ, 467, 348 [NASA ADS] [CrossRef]
• Shebalin, J. V., & Tippens, A. L. 1996, A&A, 309, 459 [NASA ADS] (In the text)
• Széll, A. 2003, Ph.D. Thesis (In the text)
• Valtonen, M., & Mikkola, S. 1991, ARA&A, 29, 9 [NASA ADS]
• Marchal, C., Yoshida, D., & Sun, Y.-S. 1984, Celest. Mech., 34, 65 [NASA ADS] [MathSciNet]
• Roy, A. E., & Steves, B. A. 2000, Cel. Mech. Dyn. Astr., 78, 299 [NASA ADS] [MathSciNet]
• Steves, B. A., & Roy, A. E. 1998, Planet Space Sci., 46, 1465 [NASA ADS] [CrossRef] (In the text)
• Steves, B. A., & Roy, A. E. 2001, in The Restless Universe: Application of Gravitational N-Body Dynamics to Planetary, Stellar and Galactic Systems, ed. B. A. Steves, & A. J. Maciejewski (Bristol: IOP Publishing), 301 (In the text)