Circumlunar free-return trajectories: a patched-conics approach

A spacecraft on a free-return trajectory departs a primary body and is sent back to its vicinity by a gravity assist from a secondary body; in other words, after orbit injection, no additional propulsive maneuver is required for the return trip. Free-return trajectories between the Earth and the Moon, in particular, have been of great practical importance for both unmanned and manned missions since the beginning of the Space Age.

In this post I discuss the construction of a circumlunar free-return trajectory by a patched-conics approach. For the Earth-Moon case, this type of analysis must account for the fact that a significant portion of the trajectory is spent inside the Sphere of Influence (SOI) of the Moon, whose extension cannot be neglected as in the case of heliocentric interplanetary missions.

I will focus in particular on the method described in Section 9.3 of Battin (1999), which was originally published in a technical report from the early Apollo era [2]. Besides its historical significance, this method affords greater insight and understanding of the problem than a more modern, computationally intensive approach. As the description in the textbook omits some key details that are important for a practical implementation of the algorithm, a secondary goal of this post is to leave a record of my personal efforts in filling those gaps.

Constraints and general strategy

A circumlunar free-return trajectory is completely determined by six independent parameters (equivalently to specifying the state vector at injection). In practice, the requirements of a real-life mission, such as Apollo, usually impose constraints at both the near-Earth end points (i.e., at departure and re-entry), as well as on the cislunar side. Our goal is to construct a circumlunar free-return trajectory for which the following quantities have pre-assigned values (cf. [1]):

The perilune distance, or minimum passing distance from the Moon, $r_m$;
The time of arrival at the lunar SOI, $t_A$;
The orbital inclination of the outbound trajectory, $i_L$;
The orbital inclination of the inbound trajectory, $i_R$;
The perigee distance of the outbound trajectory, $r_L$;
The perigee distance of the inbound trajectory, $r_R$.

Conditions (1) and (2) are directly related to mission objectives and total duration constraints (e.g., life support requirements). Conditions (3) and (4) on the orbital inclination of the geocentric arcs (with respect to the equatorial plane of the Earth) are constrained by the latitude of the launch point (e.g., $28.3^\circ$ for Cape Canaveral) and of the re-entry point. Condition (5) is especially useful under the simplifying assumption that the translunar injection maneuver occurs tangentially from a parking circular orbit of radius $r_L$, while condition (6) on the (in vacuo) perigee distance of the inbound arc controls the flight-path angle at re-entry.

The general strategy comprises two independent tasks. First, two geocentric arcs, one outbound (Earth to Moon), one inbound (Moon to Earth) are constructed, satisfying conditions (2) to (6), and terminating at the lunar SOI with velocity vectors relative to the Moon directed radially, i.e., aligned with the center of the Moon. In the second step, the two geocentric arcs are offset in the plane of the selenocentric orbit within the SOI, to satisfy condition (1). The simplicity and elegance of this method cannot be overstated.

Establishing the geocentric arcs

The procedure for establishing the geocentric arcs involves a series of iterations, but it is relative straightforward; for reference, here is a brief sketch of the key logical steps:

Assume a value for the time of flight $t_{\rm FL}$ from injection until arrival at the lunar SOI at time $t_A$.
Construct an outbound arc which satisfies the constraints (3) and (5), and whose relative velocity with respect to the Moon at arrival on the SOI at $t_A$ is directed toward the center of the Moon.
With the state vector at $t_A$ and the value of perilune distance from constraint (1), it is possible to estimate the time $t_S$ spent inside the lunar SOI.
Assume a value of the time of flight $t_{\rm FR}$ of the return trip, i.e., from the departure from lunar SOI at time $t_D = t_A+t_S$ until re-entry at Earth.
Construct an inbound arc which satisfies the constraints (4) and (6), and whose relative velocity with respect to the Moon at departure from the SOI at $t_D$ is directed radially away from the center of the Moon.
Iterate Steps 4-5 on $t_{\rm FR}$ until the relative velocity vectors at $t_A$ and $t_D$, on the outbound and inbound legs, respectively, have the same magnitudes.
Iterate Steps 3-4-5-6 on $t_{\rm FL}$ until the magnitudes and turning angle of the relative velocity vectors at $t_A$ and $t_D$ are consistent with a perilune distance formally satisfying constraint (1).

The description of this part of the algorithm in [1] is very thorough, and most of the mathematical details required in a practical implementation are provided, so I will not repeat it here. The following remarks cover a few points which were left implied in the text, but are of practical importance.

Steps 2 and 5 both involve a systematic search of the lunar SOI for a suitable geocentric terminal position vector ${\bf r}_T$ before starting a lower level iteration that brings the impact parameter with respect to the Moon close to zero. No details are provided in [1] about how the systematic search is to be performed. I found it convenient to parametrize the surface of the lunar SOI in terms of longitude and latitude angles $(\lambda_1, L_1)$ relative to a reference frame having the reference plane coinciding with the lunar orbital plane, and the principal axis pointing in the Earth-Moon direction, away from the Moon (referred to lunar state vectors calculated at instants $t_A$ and $t_D$ for Step 2 and 5, respectively). In this way, it is easy to ensure that a circumlunar trajectory is obtained, i.e., one for which the point of minimum lunar distance occurs behind the Moon (indeed, other types of free-return trajectories exist: see, for instance, [3]). To this end, the preliminary search for ${\bf r}_T$ should be limited to $\lambda_1 \in [\frac{\pi}{2}, \pi]$ and $\lambda_1 \in [\pi, \frac{3\pi}{2}]$ for the outbound and inbound arcs, respectively (see the figure).
Battin (1999) stresses the fact that we seek to construct elliptical geocentric arcs, and discusses the conditions for the existence of such solutions for given ${\bf r}_T$ and time of flight. Specifically, the geocentric sweep angle $\vartheta$, i.e., the difference in true anomaly between perigee and the lunar SOI contact point, is limited by the value $\vartheta_p$, which corresponds to a parabolic arc. The textbook does not mention explicitly that the sign difference of the sweep angle for the outbound and inbound arcs implies that the bracketing values for $\vartheta$ are $(\vartheta_p, \pi)$, and $(\pi, 2\pi - \vartheta_p)$ for the outbound and the inbound arcs, respectively. An independent condition must also be imposed due to the orbital inclination requirements, which in general should read $|\tan L| \leq |\tan i_L|$ (and not "$i_L > L$" as in the book), where $L$ is the (geocentric equatorial) latitude of ${\bf r}_T$, and analogously for the return arc.

Satisfying the perilune distance condition

Battin's textbook provides only a few pithy remarks on the final steps of the algorithm, so from now on I will discuss mostly my own interpretation and re-derivation of the missing steps. Needless to say, the following may contain mistakes, and likely cannot match the accuracy and elegance of the original version.

Using the same notation of the book, the position and velocity vectors relative to the Moon at $t_A$ on the outbound arc constructed in the first part of the algorithm are ${\bf r}_{\rm TML}$ and ${\bf v}_{\rm TML}$, while their counterparts at $t_D$ on the inbound arc have the subscript "TMR". The book states that the vectors ${\bf r}_{\rm TML}$ and ${\bf r}_{\rm TMR}$ are to be offset "in the plane determined by the relative velocity vectors"; I take this sentence to mean that the two terminal relative position vectors should be changed as:

${\bf r}_{\rm TML} = {\bf r}_{\rm TML} + r_a\, {\bf u}_{\rm a, L}$,

with

${\bf u}_{\rm a, L} = \dfrac{{\bf r}_{\rm TML} \times ({\bf v}_{\rm TML} \times {\bf v}_{\rm TMR})}{| {\bf r}_{\rm TML} \times ({\bf v}_{\rm TML} \times {\bf v}_{\rm TMR}) |}$,

and similar for the return arc. With the new ${\bf r}_{\rm TML}$, but without changing ${\bf v}_{\rm TML}$, I determine the perilune distance using elementary two-body problem relations, and then iterate on the value of $r_a$ until constraint (1) on $r_m$ is satisfied. After the iteration on $r_a$ is completed, I recalculate the time spent in the lunar SOI, $t_S$, for overall consistency.

Displacing the terminal position vectors on the lunar SOI, however, changes the geocentric arcs in such a way that they do not satisfy the near-Earth constraints anymore. After testing many possible alternatives, I settled on the following approximate correction procedure. I run the complete algorithm a first time, and take note of the mismatch in the final values of the perigee distances, say $\delta r_L$ and $\delta r_R$. I then recalculate the trajectory from the beginning, but with the geocentric perigees constraints modified to match $r_L - \delta r_L$ and $r_R - \delta r_R$. The recalculated trajectory deviates from the original target values of $r_L$ and $r_R$ by a few km, or less, and by less than a degree in $i_L$ and $i_R$ (the $t_A$ and $r_m$ conditions are satisfied to a much greater precision by construction).

This final correction can also be turned into a full top-level iteration cycle to achieve greater accuracy (a cryptic sentence in the book may be a hint that such an iteration is part of the original method). By successively adjusting $\delta r_L$ and $\delta r_R$, three iterations of this outer cycle are sufficient to bring $r_L$ and $r_R$ within a tenth of a mile, or better, from their nominal values. My converged trajectory, calculated for the same target values of the constraints as in [1], has the following parameters:

$r_m = 1180.0$ miles (target value: $1180$ miles)
$ i_L = 28.8^\circ$ ($28.3^\circ$)
$ i_R = 35.6^\circ$ ($35.0^\circ$)
$ r_L = 4076.9$ miles ($4077$ miles)
$ r_R = 4008.0$ miles ($4008$ miles)

In my converged trajectory a moderate inconsistency exists between the state vectors at the two ends ends of the cislunar arc, $({\bf r}_{\rm TML}, {\bf v}_{\rm TML})$ and $({\bf r}_{\rm TMR}, {\bf v}_{\rm TMR})$. Namely, propagating $({\bf r}_{\rm TML}, {\bf v}_{\rm TML})$ forward in time by $t_S$, the difference in the terminal position vector with respect to ${\bf r}_{\rm TMR}$ is about $300$ km (to be compared with their absolute magnitude, which is approximately the radius of the lunar SOI, or $66300$ km).

The figures above show my converged circumlunar free-return trajectory projected in the orbital plane of the Moon and with the principal axis in the Earth-Moon direction (calculated at time $t_A+t_S/2$). They should be compared with Figures 9.11 and 9.12 in Battin's textbook. From this comparison, one can note a difference of about one hour in the time labels of my trajectory with respect to those shown in the book.

References

[1] Battin, R. H. (1999). An Introduction to the Mathematics and Methods of Astrodynamics, Revised Edition. (1999).

[2] Battin, R. H. and Miller, J. S. (1962). Circulunar Trajectory Calculations, MIT Instrumentation Laboratory Report R-353

[3] Schwaniger, A. J. (1963). Trajectories in the Earth-Moon Space with Symmetrical Free Return Properties. Technical Note D-1833. Huntsville, Alabama: NASA / Marshall Space Flight Center

orbital insights