Chaos in Planar, Circular, Restricted Three-Body Problem

In this article we analyze the motion of a test particle of a planar, circular, restricted three-body problem in resonance, using the Kustaanheimo-Stiefel formalism. We show that a good qualitative description of the motion can be reduced to three simple equations for semi-major axis, eccentricity and resonance angle. Studying these equations reveals the onset of chaos, and sheds a new light on its weak nature. The 7:4 resonance is used as an example.


Introduction
The purpose of this article is to study the motion of an effectively massless "particle" (an asteroid) orbiting a massive "primary" (the sun), and being perturbed by yet another small (relative to the primary) "satellite" (Jupiter) whose motion is assumed to be perfectly circular and in the same plane as the massless-particle's orbit.This is known as the planar, circular, restricted three-body problem which has been studied extensively in the past by various numerical, semi-analytical and (exceptionally, and with many restrictions) analytical techniques.In this article, we build a fully general, analytical (even though iterative) solution which, furthermore, does not rely on the averaging principle (which would eliminate all chaotic solutions-the very focus of our study).
Our main objective is to trace the source of chaos, which arises when the two orbital periods are in (or close to) a resonance (i.e.their ratio is a simple fraction, such as 3:1).It is known [1] that, under these circumstances, chaos occurs quite readily when the perturbing orbit has a non-zero eccentricity or inclination (mathematically, one can show that terms proportional to either of these two parameters are responsible for its creation); it is a lot more difficult to incite chaos (and find its exact source) in a planar, circular situation.From this point of view, the seemingly simplest case of the restricted three-body problem is also the most challenging!Furthermore, only some resonances lead to chaos in this case; we have selected the 7:4 resonance (the massless particle completing seven orbits, while the perturber finishes four) as their typical representative.It is hoped that this example will be a catalyst for a future, more systematic study of the phenomenon.
The technique we employ is based on Kustaanheimo-Stiefel (K-S for short) formulation of the perturbed Kepler problem, which has several advantages over other techniques: it is fully analytic; it does not rely on any averaging; it can be made arbitrarily accurate (even though, in case of a chaotic motion, no solution can accurately follow the actual orbit indefinitely); it can also be made very simple (when we sacrifice high accuracy, yet maintain a good quantitative description of the motion: this is the approach taken by this article).
Vast literature exists on the subject of chaos in the context of the restricted three-body problem in resonance (a good, even though slightly dated, review of various techniques and results can be found in [2]).Most of it relies on numerical or, at best, semi-analytic exploration of the resulting motion (being of little relevance to our goal).The remaining articles usually concentrate on the elliptic version of the problem (the easy-to-explain chaos mentioned earlier), and are "tailored" to deal with only one specific resonance at a time.They also rely heavily on the averaging principle [3], thus eliminating the chaosinducing terms we are seeking to identify!
In the following three sections we  present formulas for building a solution to a perturbed Kepler problem when the direction of perturbing force (yet unspecified) is always in the plane of the massless particle's orbit;  computing the perturbing force due to another orbiting body of relatively small mass that has a circular orbit and is in 7:4 resonance with the perturbed particle, and expanding this force in terms of orbital elements of the latter;  finding and solving a set of differential equations for the orbital elements, and demonstrating the onset of chaos in some solutions.

Preliminaries
We start by reviewing the analytical, K-S based solution to a planar perturbed Kepler problem, which (unlike its three-dimensional analog), requites only complex (not quaternion) algebra.The following is a summary of the necessary formulas, presented in its general form (i.e. able to deal with any small perturbing force  f ).
To solve the corresponding set of two equations, namely where r   r  and is the semi-major axis of the corresponding solution (the bar denotes complex conjugation).This is the celebrated K-S transformation (reduced to two dimensions).a A prospective (trial) solution can be parametrized in the following way

 
O  function of (taking care of the orbit's small deviation from a perfect ellipse).At this point, we are yet to find expressions for the time derivatives of the orbital elements, and for the exact form of (thus constructing the actual solution).
  Using the proposed transformation ( 2) and ( 3), one can show that (1) is equivalent to where prime implies differentiation with respect to .s To solve the last equation, we start with the unperturbed solution (assuming constant orbital elements), evaluate the left hand side of (4) in the particle's Kepler frame (i.e. a non-inertial frame in which the new x axis points towards the aphelion; this is achieved by multiplying the equation by exp i 2 ) and call the resulting expression .When is an explicit function of , we must also make the following replacement: where o is a "typical" value of (in our case set to the resonant value of 1-see the next section), and Once (the left hand side of Equation ( 40)) has been computed, it needs to be converted to ( 0 r being the unperturbed ), based on which we can find r  -accurate values of and of the orbital elements' s-derivatives, by means of the following procedure.
,  First, needs to be expanded in powers of (we are assuming a periodic perturbing force, which makes this expansion discrete), thus: where is an integer, and m (these can be computed with surprising ease by a few lines of Mathematica code).
To build an 2  -accurate, 3   -accurate, etc. solution, we simply repeat the whole procedure in an iterative manner (the exact criteria for convergence of this procedure have not been worked out yet, but there is is enough computational evidence to feel more than safe with our 0.001

 
).The ensuing fast convergence enables us to construct an analytic solution to any level of accuracy (and be able to readily verify its correctness).To complete the solution, we divide , a   and   by   -thus replacing the physically meaningless independent variable s by the so called eccentric anomaly -and integrate the resulting differential equations.For full details, the reader should consult [4].
Quite often (as in this article), instead of high accuracy, we are interested in a solution which captures all qualitative features of the particle's motion, using the smallest number of terms (to gain a valuable conceptual understanding of the problem).Luckily, in almost all cases (critical inclination of the main satellite problem being one of the few exception), this can be achieved with only a handful of terms of the  -accurate solution.

Perturbing Force for 7:4 Resonance
We now apply the formalism of the previous section to a special case of the restricted three-body problem (solving for the motion of the "massless" particle).
We assume that, by a choice of units, the primary's (sun's) gravitational mass  is equal to 1, the perturbing body's (Jupiter's) gravitational mass is  and its (circular) orbit is (to give a concrete example, we will concentrate on this specific resonance; one can investigate other resonances using the same approach-no special tailoring required).We investigate the motion of a particle (asteroid) of a negligible mass, whose semimajor axis has values very close to 1 (we can thus take ), i.e. whose motion is in or near the corresponding resonance.In addition to this, we assume that both orbits are in the same plane (the so called planar, circular, restricted problem of three bodies; see [5]).
The perturbing force (on the massless particle) is   Utilizing (6), we can express in the Kepler frame of the massless particle as (based on the denominator of  ), and where  is the so called resonance variable.It is a multiple of the angle from the small body's aphelion to its conjunction with the larger orbiting body.
The perturbing force (excluding the  factor) must then be expanded, first in  (an inevitable step of this procedure), getting (in the actual computation we use terms up to an including 4 ,  but they are too long to be quoted in full), and then in both and (as generalized Fourier expansion), utilizing where we are displaying only the coefficients (in terms of the expansion) needed to compute time derivatives of the orbital elements (our only concern).We have also expanded the result in powers of , quoting the leading term only (equivalent to setting i     > 2 ; i none of these contribute to the solution of the next section.

Resulting Solution
It is now fairly routine, by utilizing formulas of Section 2, to build the corresponding solution (at this point, we bring back additional terms of the  -expansion).For the purpose of this paper, we perform only one step of the iterative procedure, resulting in an  -accurate solu- tion (whose typical accuracy within the next orbits is in the 50 3 10  to 2 10  range-anyone interested in constructing a more accurate analytic solution using this method should consult [6]).Furthermore, due to high sensitivity to initial conditions, any solution will sooner or later diverge from the exact answer, and become practically useless for a long-range forecast.Luckily this is not our aim here-we want to find a solution which would capture all qualitative features of the corresponding motion, with the least number of terms-so let us proceed doing just that.
As the crudest approximation, we can apply the averaging "principle" (a euphemism for yet another, often crude, approximation) and discard all terms proportional to any non-zero power of Furthermore, for each orbital-element derivative, we will keep only the leading term of the corresponding .  -expansion.This yields:   Even though this set of equations has a surprising variety of solutions, these are always very regular (no trace of chaos, which we know some solutions should possess).Furthermore, the agreement with the "exact" (highly accurate numerical) solution is quite poor-the averaging principle is thus totally inadequate in this case.This is due to to the resonance being of the third order (the same kind of approximation would be sufficiently accurate for first-order resonances, such as 2:1, 3:2, etc.).
To correct the situation, we proceed to add the leading (in terms of the  expansion)  -dependent term to each right-hand side of ( 16), getting This time, the resulting solution is always (for relatively small  values) sufficiently close to the exact solution to be considered its adequate representation, with all the qualitative details (see Figure 1, where the solid, dotted and dashed lines represent the exact, the fully  -accurate and the current, truncated solution, respectively; note that the solid and dotted lines practically overlap).
It is worth mentioning that we use our definition of orbital elements (adjusted for the orbit's distortion) for all three solutions (the osculating orbital elements have a much more irregular behaviour than what we see in Figure 1).The exact solution has been obtained by highly accurate numerical integration of (1), converting the results to the corresponding variation in  .The accurate solution corresponds to a numerical solution of an untruncated (i.e. using all terms of the  expansion) version of (13), and is nearly indistinguishabe from the exact solution.The "current" solution refers to the result of numerical integration of (13); this is the only solution which we continue to investigate further.
It is important to note that, unlike (12), the new set of Equations ( 13) already captures one of the most important feature of the exact solution, namely the occasional onset of chaos.This can be conveniently explored by the following analog of the Poincaré surfaces of section: Using the numerical solution of (13), we collect the and a  values as found at the end of each anomalistic cycle (i.e. at 0, 2π, 4π,    ), plotting these in a twodimensional graph.As a result, we discover an amazing variety of both open and closed, connected and disconnected curves on one hand (see Figure 2, for an example), or (using a slightly different set of initial conditions) a scattered cloud of disconnected point on the other (see Figure 3), indicating the presence of chaotic motion.Incidentally, when we apply the same technique to different solutions of (12), we always end up with a simple, (near) straight-line segment.

Conclusion
We have thus illustrated that all qualitative features of the motion of a small body, in the 7:4 resonance with a bigger one, are encapsulated in three simple differential Equations (13).A similar set of equations can be easily constructed for any other resonance (a routine exercisea systematic study of these is yet to be carried out).If desired, one can also add the effect of the large-body's eccentricity and the two orbits' inclination (a threedimensional version of the K-S technique will be required to deal with the latter).We believe that this technique opens a new level of understanding of resonances related to the three-body restricted problem, and of the resulting weak chaos.

.
The last identity serves to establish the following convention: from now on, in all expressions involving an extra summation (over the original ) is implied.
find and the orbital elements' s-derivatives, using ,  Copyright © 2013 SciRes.AM leading (trailing) ellipsis indicates terms proportional to i Z with 1