Analytical investigations of quasi-circular frozen orbits in the Martian gravity field

Frozen orbits are always important foci of orbit design because of their valuable characteristics that their eccentricity and argument of pericentre remain constant on average. This study investigates quasi-circular frozen orbits and examines their basic nature analytically using two different methods. First, an analytical method based on Lagrangian formulations is applied to obtain constraint conditions for Martian frozen orbits. Second, Lie transforms are employed to locate these orbits accurately, and draw the contours of the Hamiltonian to show evolutions of the equilibria. Both methods are verified by numerical integrations in an 80\times80 Mars gravity field. The simulations demonstrate that these two analytical methods can provide accurate enough results. By comparison, the two methods are found well consistent with each other, and both discover four families of Martian frozen orbits: three families with small eccentricities and one family near the critical inclination. The results also show some valuable conclusions: for the majority of Martian frozen orbits, argument of pericentre are kept at 270 degrees because J3 has the same sign with J2; while for a minority of ones with low altitude and low inclination, argument of pericentre are possible to be kept at 90 degrees because of the effect of the higher degree odd zonals; for the critical inclinations cases, argument of pericentre can also be kept at 90 degrees. It is worthwhile to note that there exist some special frozen orbits with extremely small eccentricity, which could provide much convenience for reconnaissance. Finally, the stability of Martian frozen orbits is estimated based on the trace of the monodromy matrix. The analytical investigations can provide good initial conditions for numerical correction methods in the more complex models.

💡 Research Summary

The paper investigates quasi‑circular frozen orbits around Mars—orbits for which the eccentricity and argument of pericentre (ω) remain constant on average—using two complementary analytical techniques. First, the authors determine how many zonal harmonics must be retained for an accurate yet tractable analytical model. By comparing numerical integrations with the full GMM‑2B gravity field (degree and order 80), they find that including terms up to J₉ is sufficient: truncation at J₅ leads to large deviations (≈40° in ω after five years), while J₇ already yields small errors, and J₉ reduces them further to a few degrees. Consequently, all subsequent analytical work uses a zonal‑only model up to J₉.

The first analytical method is based on Lagrange’s planetary equations. The disturbing potential is expressed in terms of Kaula’s inclination‑ and eccentricity‑functions, averaged over the mean anomaly, and separated into even‑degree (secular) and odd‑degree (long‑period) contributions. The long‑period rates of eccentricity and ω caused by the odd‑degree zonals are derived (Eqs. 7‑8). Setting the averaged rates to zero yields the frozen‑orbit conditions. The ω‑equation (Eq. 12) shows that ω must be either 90° or 270° for the secular term to vanish. Substituting ω=90° or 270° into the eccentricity equation (Eq. 11) gives the required eccentricity as a function of semimajor axis a and inclination i. For a representative a = 3597 km, the authors plot log e versus i (Fig. 2) and identify four families of frozen orbits:

• Family 1S: small e, ω = 90°, away from the critical inclination;
• Families 2S and 4S: small e, ω = 270°, also away from the critical inclination;
• Family 3S: near the critical inclination (≈63.43°), ω = 90°.

A notable feature is a “boundary” orbit at i ≈ 12° where e becomes extremely small (≈10⁻⁴), offering an almost perfectly circular path that could be advantageous for high‑resolution reconnaissance. The authors also generate contour maps of e over a wide range of i (0–20°) and a (3397–4397 km) (Fig. 3), confirming that low‑inclination, low‑altitude frozen orbits can indeed have ω = 90° because higher‑order odd zonals (J₅, J₇, …) become significant.

The second analytical approach employs Lie‑transform perturbation theory. The full Hamiltonian, including zonal terms up to J₉, is normalized via successive Lie transformations, yielding a simplified (averaged) Hamiltonian that depends only on the slow variables (e, ω). By plotting level curves of this Hamiltonian in the (e, ω) phase plane, the authors locate equilibrium points that correspond to frozen orbits. The contour plots reproduce exactly the four families found with the Lagrangian method, confirming the consistency of the two techniques. Moreover, the Hamiltonian view clarifies how the equilibria evolve as i varies, especially the bifurcation that occurs near the critical inclination where the ω = 90° equilibrium appears.

To validate the analytical predictions, a numerical integration in the full 80×80 GMM‑2B field is performed for one year. Using the initial condition derived from the Lagrangian analysis (a = 3597 km, i = 50°, e ≈ 0.00746, ω = 270°), the eccentricity oscillates with an amplitude of ≈0.0025 around the nominal value, while ω oscillates within ±20°. This behavior satisfies the definition of a quasi‑frozen orbit and demonstrates that the analytical initial conditions are sufficiently accurate for practical mission design.

Stability is assessed by computing the trace of the monodromy matrix (the state‑transition matrix over one orbital period). For all inclinations between 0° and 90°, the trace remains within the bounds (|Tr| < 2) that indicate linear stability. Hence, the identified frozen orbits are not only stationary in the averaged sense but also dynamically stable over many revolutions.

The paper’s key contributions are:

1. Demonstrating that a zonal‑only model up to J₉ captures the essential dynamics of Martian frozen orbits, making analytical treatment tractable.
2. Providing closed‑form frozen‑orbit conditions that reveal why the majority of Martian frozen orbits have ω = 270° (the sign of J₃ matches J₂), while a minority, especially low‑altitude/low‑inclination cases, can maintain ω = 90° due to higher‑order odd zonals.
3. Identifying a new family (1S) with extremely low eccentricity, which could be exploited for high‑precision mapping or communication relay missions.
4. Validating two independent analytical methods (Lagrangian and Lie‑transform) against high‑fidelity numerical simulations, establishing confidence in the derived initial conditions for subsequent numerical correction or optimization.

Future work suggested includes incorporating non‑gravitational perturbations (atmospheric drag, solar radiation pressure), third‑body effects (Phobos, Deimos, Sun), and extending the stability analysis to the full nonlinear regime. Integrating these analytical seeds into differential correction or optimal‑control frameworks would enable the design of long‑lived, fuel‑efficient Martian missions that exploit the natural stability of frozen orbits.