# How can we determine the most appropriate central body for orbit propagation?

I have recently been experimenting with expanding my numerical propagator to satellites that orbit central bodies other than Earth. I started with the Moon (using GRGM1200B gravity model), and found that, even though in principle it should be possible to integrate the trajectory using GCRF (Geocentric ICRF), in reality the results were massively better when using Lunar-centered ICRF.

I can always select the appropriate central body if I know beforehand what a given satellite is orbiting, but I would like to ask, is there any way to determine the best central body without a priori knowledge? Or said in another way, is it possible to determine the optimal central body for propagation, given only a state vector of position and velocity?

For simplification, let's say the choice only needs to be made amongst the bodies for which JPL DE ephemerides provide Chebyshev coefficients (i.e., Sun, Moon, main planets and Pluto).

Several years ago I had an intern investigate this very problem. A good intern task is an interesting but nonessential problem. This qualified as such; he even managed to turn that internship work as the basis of his masters thesis. Note: I obtained the author's permission to link to the thesis, so I've added a link to his thesis at the end of this answer.

What he found was that

• Using single precision floating point arithmetic is a bad idea for modeling an Earth to Moon trajectory or a Moon to Earth trajectory.
• Using double precision floating point arithmetic can be a good idea if one switches from geocentric to Moon-centric (or vice versa) at a "reasonable" place.
• Using extended precision extends the concept of "reasonableness".

A "reasonable" place to switch from Earth-centered to Moon-centered (or vice versa) was a place that resulted in an acceptably small loss in accuracy. Switching from Earth-centered to Moon-centered while the vehicle was still in low Earth orbit was not "reasonable", nor was waiting to switch until after performing insertion into a low lunar orbit. (Exception: Almost any place, or no place at all qualified as "reasonable" with 100 decimal place precision arithmetic.) With doubles, there was a fairly broad band of what qualified as "reasonable".

The Apollo program used the Laplace sphere of influence:

$$r_{\text{SOI}} \approx a \left(\frac m M\right)^{2/5}$$

where $$a$$ is the semi-major axis length, $$m$$ is the mass of the smaller object (the Moon in the case of a transfer from the Earth to the Moon), and $$M$$ is the mass of the larger object (the Earth in this case). This turned out to be "reasonable". Then again, so did the Hill sphere:

$$r_{\text{Hill}} \approx a \left(\frac m {3M}\right)^{1/3}$$

For the Earth-Moon system, with $$a\approx 385000\,\text{km}$$ and $$m/M\approx 0.0123$$ (a handy number to remember), the Laplace sphere of influence radius is about 66300 km while the Hill sphere radius is about 61600 km.

The optimal place to switch was in between the sphere of influence and the Hill sphere, but a bit closer to the sphere of influence. With double precision floating point arithmetic, the error in using either the sphere of influence or the Hill sphere was rather small compared to the error that resulted from using this optimal transition point. You won't get fired for choosing either the sphere of influence or the Hill sphere as the transition point.

Reference:

M. Vautier, Effect of Coordinate Switching on Simulation Accuracy of Translunar Trajectories, August 2008

• Very interesting results and project indeed! Thinking about these results, I am now also wondering not just at what point the change of center of coordinates should be performed, but also at what point the gravitational attraction by the central body cannot just be considered as a point attraction and instead needs to be modeled with spherical harmonics of higher degrees... I guess that will be a separate question though!
– Rafa
Feb 17, 2022 at 7:40
• I've asked that former intern whether I can link to his thesis. Doing so would put his name out in public, which might or might not be a good thing. So I'm waiting. Feb 19, 2022 at 6:58
• That would be amazing, thanks a lot for asking him! I would be extremely interested in reading it, if he is OK with it of coursen
– Rafa
Feb 19, 2022 at 9:18
• @Rafa Approval granted, so I'm adding a link to Mana's thesis. Feb 21, 2022 at 14:35
• I'm going through it at the moment and this is excellent! Thanks a lot again to you and Mana for providing this. I will certainly cite it in my propagator!
– Rafa
Feb 22, 2022 at 2:09

A common approach is to calculate the sphere of influence of the celestial objects whose gravity you're accounting for in your propagation.

$$r_{SOI} \simeq a \left( \frac M m \right)^{\frac 2 5}$$

Where $$a$$ is the semi major axis of the smaller object compared to the larger one, and $$M$$ and $$m$$ are the masses of the larger and smaller objects, respectively.

In practice, the exact moment at which the switching of the central body happens does not typically matter a whole lot as long as it is done. Usually, astrodynamics engineers choose one point in the trajectory where the central body is updated, and continue with that for a while.

Keep in mind that every frame transformation in the propagation is going to lead to a small discontinuity in the orbit itself due to the maximum precision of the floating point values on a computer and the precision of the ephemerides used. As such, it is strongly recommended to not repeatedly change the central body in a short duration.

Finally, in terms of physics, the reason that change is important is due to the third body formulation:

$$\mathbf{r_{j}} = \mathbf{r}_t - \mathbf{r_{ij}}$$

$$\dot{\mathbf{v}}_{t'} = \dot{\mathbf{v}}_{t'} - \mu_i \left( \frac{\mathbf{r_j}}{|\mathbf{r_j}|^3} - \frac{\mathbf{r_{ij}}}{|\mathbf{r_{ij}}|^3} \right)$$

If those third body effects $$r_j$$ lead to significantly greater perturbation than the central body, the propagation will become unstable.

• This is excellent, thanks a lot! I will go in detail through it and play with the concepts :) Exactly, instability is what I saw when using GCRF to propagate Moon-orbiting satellites (the orbit oscillated in a strange way), which went away simply by moving the center of the frame to the Moon. Nice to see everything falling into place slowly!
– Rafa
Feb 17, 2022 at 1:52
• That formulation of third body effects implies that acceleration toward the central body should be calculated first and then third body perturbations be added to that. Using finite precision arithmetic (e.g., double precision arithmetic), it is better from a numerical precision perspective to perform the summation from smallest in magnitude to largest, which typically means that the gravitational acceleration to the central body is the last to enter the sum. Feb 22, 2022 at 9:02
• Also note that in An Introduction to the Mathematics and Methods of Astrodynamics, Battin provides an alternative formulation of the third body effect that can sometimes be more precise due to precision loss in the subtraction. By the way despite the name ("An Introduction To ..."), this is a graduate level text. Feb 22, 2022 at 9:07
• Concerning the limitations of finite precision arithmetic, is it still valid today with 64 bit floating points? On a 64 bit machine, adding smaller and smaller numbers in ascending or descending order changes the final value only at the 15th decimal: play.rust-lang.org/… . I get the same results in Python on a 64 bit processor as well. Feb 22, 2022 at 21:48