A second-order solution to the two-point boundary value problem for rendezvous in eccentric orbits

A new second-order solution to the two-point boundary value problem for relative motion about orbital rendezvous in one orbit period is proposed. First, nonlinear differential equations to describe the relative motion between a chaser and a target are presented considering the second-order terms in the gravity. Then, by regarding the second-order terms as external accelerations, we establish second-order state transition equations. Moreover, the J2 perturbations effects can also be considered in the state transition equations. Last, the initial relative velocity to fulfill a rendezvous is determined by solving the state transition equations. Numerical simulations show that the new second-order state transition equations are accurate. The second-order solution to the two-point boundary value problem on eccentric orbits is valid even if the relative range is farther than 500 km.     

Approximate solutions of non-linear circular orbit relative motion in curvilinear coordinates

A compact, time-explicit, approximate solution of the highly non-linear relative motion in curvilinear coordinates is provided under the assumption of circular orbit for the chief spacecraft. The rather compact, three-dimensional solution is obtained by algebraic manipulation of the individual Keplerian motions in curvilinear, rather than Cartesian coordinates, and provides analytical expressions for the secular, constant and periodic terms of each coordinate as a function of the initial relative motion conditions or relative orbital elements. Numerical test cases are conducted to show that the approximate solution can be effectively employed to extend the classical linear Clohessy–Wiltshire solution to include non-linear relative motion without significant loss of accuracy up to a limit of 0.4–0.45 in eccentricity and 40–45\(^\circ \) in relative inclination for the follower. A very simple, quadratic extension of the classical Clohessy–Wiltshire solution in curvilinear coordinates is also presented.     

Second-order state transition for relative motion near perturbed,elliptic orbits

This paper develops a tensor and its inverse, for the analytical propagation of the position and velocity of a satellite, with respect to another, in an eccentric orbit. The tensor is useful for relative motion analysis where the separation distance between the two satellites is large. The use of nonsingular elements in the formulation ensures uniform validity even when the reference orbit is circular. Furthermore, when coupled with state transition matrices from existing works that account for perturbations due to Earth oblateness effects, its use can very accurately propagate relative states when oblateness effects and second-order nonlinearities from the differential gravitational field are of the same order of magnitude. The effectiveness of the tensor is illustrated with various examples.     

A satellite relative motion model including $$J_2$$ and $$J_3$$J3 via Vinti’s intermediary

Vinti’s potential is revisited for analytical propagation of the main satellite problem, this time in the context of relative motion. A particular version of Vinti’s spheroidal method is chosen that is valid for arbitrary elliptical orbits, encapsulating \(J_2\), \(J_3\), and generally a partial \(J_4\) in an orbit propagation theory without recourse to perturbation methods. As a child of Vinti’s solution, the proposed relative motion model inherits these properties. Furthermore, the problem is solved in oblate spheroidal elements, leading to large regions of validity for the linearization approximation. After offering several enhancements to Vinti’s solution, including boosts in accuracy and removal of some singularities, the proposed model is derived and subsequently reformulated so that Vinti’s solution is piecewise differentiable. While the model is valid for the critical inclination and nonsingular in the element space, singularities remain in the linear transformation from Earth-centered inertial coordinates to spheroidal elements when the eccentricity is zero or for nearly equatorial orbits. The new state transition matrix is evaluated against numerical solutions including the \(J_2\) through \(J_5\) terms for a wide range of chief orbits and separation distances. The solution is also compared with side-by-side simulations of the original Gim–Alfriend state transition matrix, which considers the \(J_2\) perturbation. Code for computing the resulting state transition matrix and associated reference frame and coordinate transformations is provided online as supplementary material.     

Satellite Relative Motion Using Differential Equinoctial Elements

For precise control, to minimize the fuel consumption, and to maximize the lifetime of satellite formations a precise analytic solution is needed for the relative motion of satellites. Based on the relationship between the relative states and the differential orbital elements, the state transition matrix for the linearized relative motion that includes the effects due to the reference orbit eccentricity and the gravitational perturbations is derived. This method is called the Geometric Method. To avoid any singularities at zero eccentricity and zero inclination, equinoctial variables are used to derive the relative motion state transition matrices for both mean and osculating elements. This approach can be extended easily to include other perturbing forces.     

A new family of multistep numerical integration methods based on Hermite interpolation

In this paper, a new family of explicit and implicit multistep methods is presented both for the error-controlled and uncontrolled modes. The main concept is to replace the Newton interpolation with the Hermite interpolation, where the Hermite polynomial is fitted to the function values and its derivatives. This idea is very useful in the numerical solution of problems (e.g., orbit propagation problem) where higher-order derivatives can easily be computed. In addition to the theoretical concept, the stability regions of the proposed methods are determined. The new methods are more stable than the well-known multistep numerical integrators (i.e., Adams–Bashforth and Adams–Bashforth–Moulton) in the explicit, implicit, and predictor–corrector forms. Using the second-order derivatives gives smaller error constants in the proposed method. The new integrators are numerically tested for a few examples, and the solutions are compared with those of the well-known multistep methods. Moreover, the CPU time and absolute integration error are compared in the satellite orbit propagation problem using various integration methods. The CHAMP mission, i.e., a German small-satellite mission for geoscientific and atmospheric research and applications, is considered as a case study for comparing the achievable accuracy of the proposed method with the existing method for solving the two-body problem.     

Nonlinear analytic solution for perturbed relative motion using differential equinoctial elements

This paper develops a nonlinear analytic solution for satellite relative motion in J₂-perturbed elliptic orbits by using the geometric method that can avoid directly solving the complex differential equations. The differential equinoctial elements (DEEs) are used to remove any singularities for zero-eccentricity or zero-inclination orbits. Based on the relationship between the relative states and the DEEs, state transition tensors (STTs) for transforming the osculating DEEs and propagating the mean DEEs have been derived. The formulation of these STTs has been split into a set of vector and matrix operations, which avoids directly expanding the complex second-order terms, and thus, the obtained STTs could be easy-to-understand and easy-to-code. Numerical results show that the proposed nonlinear solution is valid for zero-eccentricity and zero-inclination reference orbit and is more accurate than the previous linear or nonlinear methods for the long-term prediction of satellite relative motion.     

First-order analytic propagation of satellites in the exponential atmosphere of an oblate planet

The paper offers the fully analytic solution to the motion of a satellite orbiting under the influence of the two major perturbations, due to the oblateness and the atmospheric drag. The solution is presented in a time-explicit form, and takes into account an exponential distribution of the atmospheric density, an assumption that is reasonably close to reality. The approach involves two essential steps. The first one concerns a new approximate mathematical model that admits a closed-form solution with respect to a set of new variables. The second step is the determination of an infinitesimal contact transformation that allows to navigate between the new and the original variables. This contact transformation is obtained in exact form, and afterwards a Taylor series approximation is proposed in order to make all the computations explicit. The aforementioned transformation accommodates both perturbations, improving the accuracy of the orbit predictions by one order of magnitude with respect to the case when the atmospheric drag is absent from the transformation. Numerical simulations are performed for a low Earth orbit starting at an altitude of 350 km, and they show that the incorporation of drag terms into the contact transformation generates an error reduction by a factor of 7 in the position vector. The proposed method aims at improving the accuracy of analytic orbit propagation and transforming it into a viable alternative to the computationally intensive numerical methods.     

High order transfer maps for perturbed Keplerian motion

The paper presents a new semi-analytical technique for the propagation of near-Earth satellite motion. The approach uses differential algebra techniques to compute the high order expansion of the solution of the system’s ordinary differential equation for one orbital revolution, referred to as the transfer map. Once computed, a single high order transfer map (HOTM) can be reused to map an initial condition, or a set of initial conditions, forward in time for many revolutions. The only limiting factor is that the mapped objects must stay close to the reference orbit such that they remain within the region of validity of the HOTM. The performance of the method is assessed through a set of test cases in which both autonomous and non-autonomous perturbations are considered, including the case of continuously propelled trajectories.     

The theory of asynchronous relative motion I: time transformations and nonlinear corrections

Using alternative independent variables in lieu of time has important advantages when propagating the partial derivatives of the trajectory. This paper focuses on spacecraft relative motion, but the concepts presented here can be extended to any problem involving the variational equations of orbital motion. A usual approach for modeling the relative dynamics is to evaluate how the reference orbit changes when modifying the initial conditions slightly. But when the time is a mere dependent variable, changes in the initial conditions will result in changes in time as well: a time delay between the reference and the neighbor solution will appear. The theory of asynchronous relative motion shows how the time delay can be corrected to recover the physical sense of the solution and, more importantly, how this correction can be used to improve significantly the accuracy of the linear solutions to relative motion found in the literature. As an example, an improved version of the Clohessy-Wiltshire (CW) solution is presented explicitly. The correcting terms are extremely compact, and the solution proves more accurate than the second and even third order CW equations for long propagations. The application to the elliptic case is also discussed. The theory is not restricted to Keplerian orbits, as it holds under any perturbation. To prove this statement, two examples of realistic trajectories are presented: a pair of spacecraft orbiting the Earth and perturbed by a realistic force model; and two probes describing a quasi-periodic orbit in the Jupiter-Europa system subject to third-body perturbations. The numerical examples show that the new theory yields reductions in the propagation error of several orders of magnitude, both in position and velocity, when compared to the linear approach.     

Satellite onboard orbit propagation using Deprit’s radial intermediary

Short-term satellite onboard orbit propagation is required when GPS position measurements are unavailable due to an obstruction or a malfunction. In this paper, it is shown that natural intermediary orbits of the main problem provide a useful alternative for the implementation of short-term onboard orbit propagators instead of direct numerical integration. Among these intermediaries, Deprit’s radial intermediary (DRI), obtained by the elimination of the parallax transformation, shows clear merits in terms of computational efficiency and accuracy. Indeed, this proposed analytical solution is free from elliptic integrals, as opposed to other intermediaries, thus speeding the evaluation of corresponding expressions. The only remaining equation to be solved by iterations is the Kepler equation, which in most of cases does not impact the total computation time. A comprehensive performance evaluation using Monte-Carlo simulations is performed for various orbital inclinations, showing that the analytical solution based on DRI outperforms a Dormand–Prince fixed-step Runge–Kutta integrator as the inclination grows.     

Three-dimensional subsatellite motion under air drag and oblateness perturbations

The three-dimensional relative motion of a subsatellite with respect to a reference station in an elliptical orbit is studied. A general theory based on the variation of the relative elements, i.e. the instantaneous differences between the orbital parameters of the subsatellite and those of the station, is formulated in order to incorporate arbitrary perturbing forces acting on both satellites. The loss of precision inherent in the subtraction of almost identical quantities is avoided by the consistent use of difference variables. In the absence of perturbations exact analytical representations can be obtained for the relative state parameters. The influences of air drag and Earth's oblateness on the relative motion trajectories are investigated and illustrated graphically for a number of cases.     

Variational and symplectic integrators for satellite relative orbit propagation including drag

Orbit propagation algorithms for satellite relative motion relying on Runge–Kutta integrators are non-symplectic—a situation that leads to incorrect global behavior and degraded accuracy. Thus, attempts have been made to apply symplectic methods to integrate satellite relative motion. However, so far all these symplectic propagation schemes have not taken into account the effect of atmospheric drag. In this paper, drag-generalized symplectic and variational algorithms for satellite relative orbit propagation are developed in different reference frames, and numerical simulations with and without the effect of atmospheric drag are presented. It is also shown that high-order versions of the newly-developed variational and symplectic propagators are more accurate and are significantly faster than Runge–Kutta-based integrators, even in the presence of atmospheric drag.     

Element sets for high-order Poincaré mapping of perturbed Keplerian motion

The propagation and Poincaré mapping of perturbed Keplerian motion is a key topic in Celestial Mechanics and Astrodynamics, e.g., to study the stability of orbits or design bounded relative trajectories. The high-order transfer map (HOTM) method enables efficient mapping of perturbed Keplerian orbits using the high-order Taylor expansion of a Poincaré or stroboscopic map. The HOTM is only accurate close to the expansion point and therefore the number of revolutions for which the map is accurate tends to be limited. The proper selection of coordinates is of key importance for improving the performance of the HOTM method. In this paper, we investigate the use of different element sets for expressing the high-order map in order to find the coordinates that perform best in terms of accuracy. A new set of elements is introduced that enables extremely accurate mapping of the state, even for high eccentricities and higher-order zonal perturbations. Finally, the high-order map is shown to be very useful for the determination and study of fixed points and center manifolds of Poincaré maps.     

The determination of the satellite orbit of Mariner 9

This paper presents a comprehensive analysis of the Mars orbital phase of the Mariner 9 trajectory as determined from Earth based radio data. Both the method and accuracy of the orbit determination process are reviewed. Analysis is presented to show the effects of Mars gravity model and node in the plane of the sky errors on the accuracy of orbit determination. In addition the long term evolution of the orbit from insertion through the first 500 revolutions is presented, and decomposed into effects from the Mars garvity field,n-body perturbations, and solar radiation pressure. Since the orbit period is nearly commensurable with the Mars rotational period, the orbit experiences significant resonance perturbations. The primary perturbation is in-track with a maximum amplitude of 1000 km and a wavelength of 39 spacecraft revolutions.This paper was presented at the AIAA/AAS Astrodynamics Conference, Palo Alto, California, September 11 and 12, 1972. At this time Mariner 9 operations were still underway. The operational life of Mariner 9 ended October 27, 1972, when the supply of nitrogen gas, used for attitude stabilization, was depleted. This paper represents one phase of research carried out at the Jet Propulsion Laboratory, California Institute of Technology, Pasadena, California, under NASA Contract No. NAS 7-100.     

Constructive-analytical solution of the problem of the secular evolution of polar satellite orbits

The well-known twice-averaged Hill problem is considered by taking into account the oblateness of the central body. This problem has several integrable cases that have been studied qualitatively by many scientists, beginning with M.L. Lidov and Y. Kozai. However, no rigorous analytical solution can be obtained in these cases due to the complexity of the integrals. This paper is devoted to studying the case where the equatorial plane of the central body coincides with the plane of its orbital motion relative to the perturbing body, while the satellite itself moves in a polar orbit. A more detailed qualitative study is performed, and an approximate constructive-analytical solution of the evolution system in the form of explicit time dependences of the eccentricity and pericenter argument of the satellite orbit is proposed. The methodical accuracy for the polar orbits of lunar satellites has been estimated by comparison with the numerical solution of the system.     

Coordinate Additional Perturbations to Mars Orbiters and Choice of Corresponding Coordinate System

Similar to the study of the related problems of Earth satellites, in the research of the motion of Mars orbiter especially for low-orbit satellites, it is more appropriate to choose an epoch Mars-centered and Mars-equator reference system, which indeed is called the Mars-centered celestial coordinate system. In this system, the xy-plane and the direction of the x-axis correspond to the mean equator and mean equinox. Similar to the precession and nutation of the Earth, the wiggling of instantaneous Mars equator causes the coordinate additional perturbations in this Mars coordinate system. The paper quotes a method which is similar to the one used in dealing with the coordinate additional perturbations of Earth. According to this method, based on the IAU2000 Mars orientation model and under the precondition of a certain accuracy, we are able to figure out the precession part of the change of Mars gravitation. This lays the foundation for further study of its influence on the Mars orbiter's orbit of precession and the solution of the corresponding coordinate additional perturbations. The obtained analytical solution is easy to use. Compared with the numerical solution with higher accuracy, the result shows that the accuracy of this analytical solution could satisfy the general requirements in use. Therefore, our result verifies that a unified coordinate system, the Mars-centered celestial system in which J2000.0 is chosen as its current initial epoch, could be applied to deal with the relative problems of Mars orbiters, especially for low-orbit satellites. It is different from the method we previously used in dealing with the corresponding problems of Earth satellites, where we adopted the instantaneous equator and epoch (J1950.0) mean equinox as xy-plane and the direction of x -axis. In contrast, the coordinate transformation brings heavy workload and certain inconvenience in relative former works in which the prior system is used. If adopting the unified coordinate system, the transformation could be simply avoided and the computation load could be decreased significantly.     

The Determination of an Intermediate Perturbed Orbit from Two Position Vectors

Based on the theory of intermediate orbits developed earlier by the author of this paper, a new approach is proposed to the solution of the problem of finding the orbit of a celestial body with the use of two position vectors of this body and the corresponding time interval. This approach makes it possible to take into account the main part of perturbations. The orbit is constructed, the motion along which is a combination of two motions: the uniform motion along a straight line of a fictitious attracting center, whose mass varies according to the first Meshchersky law, and the motion around this center. The latter is described by the equations of the Gylden–Meshchersky problem. The parameters of the constructed orbit are chosen so that their limiting values at any reference epoch determine a superosculating intermediate orbit with third-order tangency. The accuracy of approximation of the perturbed motion by the orbits calculated by the classical Gauss method and the new method is illustrated by an example of the motion of the unusual minor planet 1566 Icarus. Comparison of the results obtained shows that the new method has obvious advantages over the Gauss method. These advantages are especially prominent in cases where the angular distances between the reference positions are small.     

Generation of a secular system in the analytical theory for the motion of the moon

In order to generate an analytical theory of the motion of the Moon by considering planetary perturbations, a procedure of general planetary theory (GPT) is used. In this case, the Moon is considered as an addition planet to the eight principal planets. Therefore, according to the GPT procedure, the theory of the Moon’s orbital motion can be presented in the form of series with respect to the evolution of eccentric and oblique variables with quasi-periodic coefficients, which are the functions of mean longitudes for principal planets and the Moon. The relationship between evolution variables and the time is determined by a trigonometric solution for the independent secular system that describes the secular motion of a perigee and the Moon node by considering secular planetary inequalities. Principal planetary coordinates required for generating the theory of the motion of the Moon includes only Keplerian terms, the intermediate orbit, and the linear theory with respect to eccentricities and inclinations in the first order relative to the masses. All analytical calculations are performed by means of the specialized echeloned Poisson Series Processor EPSP.