How to propagate the orbit of a satellite analytically using Lagrange planetary equations/VOPs and compare them with numerically propagated solution?

Question

I have propagated the orbit of a satellite from TLE elements numerically using RK4 method. But I'm bit confused when it comes to the propagation of satellite orbits using analytical approach. From what I have come across Lagrange planetary equations/VOPs (Variation of parameters) are in term of differential equation form, so if I tried to integrate them using numerical method such as RK4 and convert them into state vector, it will become numerical or semi-analytical propagation. Could someone clarify this confusion? And please provide some clarity over difference between analytical and semi-analytical approaches in terms of propagator implementation.

Answer 1

The question is about orbit-integration problems using Lagrange planetary equations and methods of variation of parameters (VP), especially for (earth-)satellite orbits.

It's worth noting first the independence in practice of three kinds of elements forming part of such numerical problems.

These elements are (a) an underlying theory of the motion and the accelerations to be integrated, (b) a choice amongst different parameter-types, how to represent the numerical data for the integration, and (c) a choice among generic methods or computation tools for the integration.

Lagrange's planetary equations, and methods of variation of parameters, are the focus of the present question, and they form part of aspect (b), how to choose parameter-types and arrange the data most suitably for integration. This aspect is somewhat generic, in that one and the same VP method can handle the numerical output from any theory of orbital motion. Irrespective of the underlying theory, its output can be put into a form that matches the chosen VP method, which remains a method that advances or propagates the orbit numerically step by step.

The categories 'analytical' or 'semi-analytical', also mentioned in the question, do not concern this method of arranging or parameterizing the integration. In the context of satellite motion, these categories are about the nature of different kinds of theories of satellite motion. The outputs of all such theories can be inputs to the chosen numerical integration arrangement, whether VP or not; they need have no connection to the numerical integration method employed.

Thus the analytical or semi-analytical character of an underlying theory of satellite motion has no necessary effect in terms of implementing an orbit propagation method. (An example of a paper discussing different kinds, analytical and semi-analytical, of satellite theory is Setty, Cefola et al (2013) (Paper AAS 13-769 from 2013 AAS/AIAA conference).)

The rest of this answer concerns aspect (b), arrangements and parameterization of the integration process.

Traditionally, methods of determining the orbit and motion of an orbiting body in the presence of perturbing forces have been divided into 'general perturbations' and 'special perturbations'.

'General perturbations' represent a fully analytical approach, often a very arduous one. The equations of motion including perturbing accelerations derived from the relevant central force and disturbing functions are expressed analytically, in forms that involve heavy use of trig. functions. The solution for the integration is then developed in expanded trig.-series. Where this approach is feasible at all, it is usually a great help nowadays to use dedicated software for the manipulation of trig. series.

'Special perturbations' involve an at least partly numerical approach, and come in several varieties. The equations of motion may be expressed e.g. in terms of rectangular coordinates, and an integration developed purely numerically from initial values of the positions and velocities with suitable parameters for the masses and other sources of acceleration.

Methods of this class are used by the Jet Propulsion Laboratory to produce the integrated solar-system ephemerides used as basis for NASA space navigation, for much of the data in official astronomical almanacs, &c. An example of such equations of motion and outline plan of integration is discussed in Newhall et al (1983), Astron & Astrophys 125, 150-167.

Or, an attempt may be made to avoid the work of integrating the well-known component of the motion due to the central force, by integrating the initial values of a suitably chosen set of orbital parameters ('variation of parameters'). Here there can be an element of analysis involved in deriving the perturbational rates of change of the orbital parameters (which would be constant in an ideal unperturbed system). Some analytical work can also be involved in setting up the interconversions between positions and velocities and the corresponding values of the orbital parameters and of their perturbational rates of change.

In this connection, sometimes the older traditional forms of orbital elements turn out to be numerically unsuited to a VP approach for progressing the orbit and the motion. (See this recent Q & A for a discussion of some of the issues, especially about the choice of orbital parameter-sets to work with: Why is the eccentricity vector used to describe near-circular orbits? .)

If you decide to implement a VP method in your own problem, you might find good value in A E Roy's 'Orbital Motion' (either 4th ed. 2004/2005, or earlier editions almost equally useful on this topic, many sources of the text are shown on a google search). Roy is a very helpful resource on general and special perturbations, and describes several VP methods. For other perspectives on these problems, the works of R H Battin (e.g. "Introduction to the Mathematics and Methods of Astrodynamics" 1999), or older works of S Herrick, or R Baker, among others, are also very instructive. Also very useful is Bate, Mueller & White (1971) "Fundamentals of Astrodynamics".

It's not always easy to verify whether a VP approach is overall advantageous in terms of the work of setup and integration-time, compared with integration of the entire equations of motion.

The flow of a VP method for propagating an orbit and the motion of the orbiting body can consist of:- converting input 'state' data for position and motion of an orbiting body into the form of the chosen parameters/elements; obtaining from the chosen theory of motion and from the 'state' data the contemporaneous perturbing accelerations on the body; applying the relationships given by the Lagrange equations (or equivalent relations) to obtain (from the parameters and the accelerations) the contemporaneous rates of change of the parameters; using the parameters and their rates of change as inputs for integration; converting the resulting integrals (i.e. updated parameters) into updated 'state' data; and iterating the process as often as required.

Although the Lagrange equations and other parameter relations look intricate and some can certainly lead into involved computations, it is possible to express some parts of the work surprisingly simply using vector equations related to the Lagrange relationships. For each body there are normally reckoned to be six orbital elements, and VP integration schemes can call for at least six and possibly seven variables per body. If three of those variables are chosen as the components of the body's specific angular momentum vector $(H)$, and $r$, $v$, and $F_p$ are its vectors of position, velocity and the overall perturbing acceleration on it, then

$(H)$ = ($r$ x $v$),

and its rate of change

$(Hdot)$ = ($r$ x $F_p$),

which sets up three of the variable parameters and their rates of change with little computational effort.

Answer 2

In addition to the Lagrange planetary equations, there are lots of other SomeFamousPerson (no derogation intended) planetary equations. Examples include the Gauss planetary equations, the Delaunay planetary equations, the Hill planetary equations, and so on. These alternatives attempt to address two key shortcomings in the Lagrange planetary equations: The Lagrange planetary equations have singularities at zero eccentricity and at zero inclination.

One can use an integrator of choice to propagate SomeFamousPerson's planetary equations, and then convert the propagated state to Cartesian. One key advantage of integrating SomeFamousPerson's planetary equations is that doing so inherently accounts for geometry. A geometric integrator will almost inevitably fare better than a non-geometric integrator such as using RK4 to propagate Cartesian position and velocity. Integrating SomeFamousPerson's planetary equations can enable using a significantly larger step size before running into numerical issues.

A variation of parameters (VOP) approach can go one step further. For example, it's VOP that shows that there exist specific inclinations combined with specific semi-major axis lengths that result in sun synchronous orbits. Those VOP analyses also result in time-varying derivatives. Handling the secular effects separately from the non-secular derivatives can enable using an even larger step size before running into numerical issues.