# Tangential velocity constraint for orbital maneuvers with Theory of Functional Connections

A. K. de Almeida, T. Vaillant, V. M. de Oliveira, D. Barbosa, D. Maia, S.

Aljbaae, B. Coelho, M. Bergano, J. Pandeirada, A. F. B. A. Prado, et al.

► To cite this version:

A. K. de Almeida, T. Vaillant, V. M. de Oliveira, D. Barbosa, D. Maia, et al. Tangential velocity constraint for orbital maneuvers with Theory of Functional Connections. Scientific Reports, 2024, 14, 10.1038/s41598-024-57986-z. insu-04822210

HAL Id: insu-04822210

https://insu.hal.science/insu-04822210v1

Submitted on 6 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Tangential velocity constraint for orbital maneuvers with Theory of Functional Connections

A. K. de Almeida Jr.<fnref n="1" />,✉, T. Vaillant<fnref n="1" />,², V. M. de Oliveira²,³, D. Barbosa⁴, D. Maia<fnref n="1" />, S. Aljbaae⁵,⁶, B. Coelho²,⁷, M. Bergano<fnref n="1" />,⁷, J. Pandeirada⁴,⁸, A. F. B. A. Prado⁵, A. Guerman⁹ & A. C. M. Correia¹⁰,²

Maneuvering a spacecraft in the cislunar space is a complex problem, since it is highly perturbed by the gravitational influence of both the Earth and the Moon, and possibly also the Sun. Trajectories minimizing the needed fuel are generally preferred in order to decrease the mass of the payload. A classical method to constrain maneuvers is mathematically modeling them using the Two Point Boundary Value Problem (TPBVP), defining spacecraft positions at the start and end of the trajectory. Solutions to this problem can then be obtained with optimization techniques like the nonlinear least squares conjugated with the Theory of Functional Connections (TFC) to embed the constraints, which recently became an effective method for deducing orbit transfers. In this paper, we propose a tangential velocity (TV) type of constraints to design orbital maneuvers. We show that the technique presented in this paper can be used to transfer a spacecraft (e.g. from the Earth to the Moon) and perform gravity assist maneuvers (e.g. a swing-by with the Moon). In comparison with the TPBVP, solving the TV constraints via TFC offers several advantages, leading to a significant reduction in computational time. Hence, it proves to be an efficient technique to design these maneuvers.

In order to perform planetary or interplanetary maneuvers within the solar system, a spacecraft changes its velocity. Such a variation can be obtained, for instances, with a gravity assist maneuver after a close encounter with a massive body or via specific impulses due to the burn of propellant preloaded in the spacecraft. In the case of a spacecraft initially orbiting the Earth, a velocity impulse can be applied in a specific direction to guide the spacecraft toward its target. This is characterized by the flight-path angle, which is the angle between the local horizontal (the line perpendicular to the radius vector) and the velocity vector<fnref n="1" />. When this angle is equal to zero, the velocity is tangential to the local orbit. A tangential velocity represents the direction of the initial (or final) velocity after (or before) the application of the impulse of a trajectory that minimizes the costs to transfer a spacecraft between circular orbits<fnref n="1" />. This result also holds for transfers between elliptical orbits, with the point of burn located at either its perigee or apogee². A tangential velocity is the solution of the optimization procedure within a two-body problem approach³, although different solutions can take advantage of symmetries for some specific transfers between co-planar elliptical orbits. In a model more complex than the two-body problem, tangential burns turned out also from numerical evidences to be able to minimize the costs involved in impulsive maneuvers, although there is no analytical proof. For instance, this result can be seen in⁴ for Earth-to-Moon or Moon-to-Earth transfers. Furthermore, a tangential initial velocity with respect to the point of burn at the initial circular orbit represents the maximum efficiency to transfer a spacecraft from the Earth to the Moon⁵. This tangential velocity assumption is used to solve the transfer problem (and assess costs) in various works, such as⁶,⁷, employing the shooting method.

The Theory of Functional Connections (TFC) is a recently developed mathematical framework to perform linear interpolation⁸,⁹. TFC derives functionals that contain a free function and satisfy linear constraints, regardless

<footnote n="1">¹CICGE, DGAOT, FCUP, Vila Nova de Gaia, Portugal. ²CFisUC, Departamento de Física, Universidade de Coimbra, 3004-516 Coimbra, Portugal. ³Instituto de Matemática e Estatística, Universidade de São Paulo, 05508-090 São Paulo, SP, Brazil. ⁴Instituto de Telecomunicações, Universidade de Aveiro, 3810-193 Aveiro, Portugal. ⁵Postgraduate Division, National Institute for Space Research (INPE), São Paulo, Brazil. ⁶Make The Way, R. Elvira Ferraz, 250 - FL Office 305 e 306 - Vila Olímpia, São Paulo, SP 04545-015, Brazil. ⁷ATLAR Innovation, Ed. Multiusos, Rua Rangel de Lima, 3320-229 Pampilhosa da Serra, Portugal. ⁸Instituto Superior Técnico, Avenida Rovisco Pais 1, 1049-001 Lisboa, Portugal. ⁹University of Beira Interior, Covilhã, Portugal. ¹⁰IMCCE, UMR8028 CNRS, Observatoire de Paris, PSL Université, 77 Avenue Denfert-Rochereau, 75014 Paris, France. ¹¹email: allan.junior@inpe.br</footnote>
of the form of the free function. These functionals enable the search for solutions in a subspace of functions that analytically satisfy the linear constraints. TFC can be applied to transform a constrained problem into an unconstrained one, thereby converting differential equations subject to constraints into unconstrained ones¹⁰. For example, these differential equations can represent the equations of motion of a spacecraft¹¹. The procedure searches for a solution (finds the trajectory) in a subspace that exactly satisfies the specified linear constraints. This is important, because TFC embeds the type of constraints proposed in this work into the dynamical model. A nonlinear least squares method¹² is then employed to enhance efficiency of a numerical procedure to solve the problem, and, hence, to design the proposed orbital maneuvers.

The Two Point Boundary Value Problem (TPBVP) in an astrodynamics context consists of searching for a solution that enables a spacecraft to travel from point A to point B. Thus, only the positions of the spacecraft are given, and the velocities at these points are unknown and need to be determined. One technique for solving this problem is the shooting method, which reduces the problem to a sequence of initial value problems. The method relies on the integration of the initial conditions multiple times to match the boundary conditions. In astrodynamics, this procedure can leverage solutions of the Lambert's problem as an initial guess¹³, can be combined with a direct transcription¹⁴ and multiple shooting approach⁶, among others in the literature, like the sequential gradient-restoration algorithm¹⁵. The TPBVP can also be efficiently solved using TFC¹¹. For instance, spacecraft transfers from Earth to the Moon¹⁶ or between other bodies in the Solar system¹⁷ have been evaluated using TFC via the TPBVP constraints. In this paper, it is proposed the use of a set of constraints different from the TPBVP, called TV (tangential velocity) type constraints. In these constraints, the final position is free to be anywhere on a circumference around the Moon, rather than a single specific point. The approach combines analytical and numerical techniques to efficiently search for solutions satisfying the TV type constraints.

The primary objective of this paper is to introduce the TV type constraints to design orbital maneuvers and demonstrate their efficiency - via the TFC procedure - to transfer a spacecraft from Earth to the Moon and to perform a gravity assist maneuver with the Moon. The TFC allows us to improve efficiency by analytically embedding the TV type constraints into the equations of motion of the spacecraft to be maneuvered. The proposed method is compared to the TPBVP types of solution also solved via a TFC procedure. Although the TV type constraints are coupled and nonlinear in rectangular coordinates, as will be shown later, these constraints can become linear after a proper change of coordinates, as was demonstrated in¹⁸. Thus, a set of polar coordinates centered at the Moon is adopted in this paper in order to linearize the proposed constraints.

In Section “Dynamics modelling”, the reference frames and the proper (polar) variables used in this work are defined and the equations of motion are derived. Section “TV type constraints” introduces the TV type constraints proposed in this paper, outlines their application scope with TFC, and highlights their advantages in comparison to the TPBVP type of constraints. In Sect. 4, a practical application of the TV type constraints is provided by transferring a spacecraft from a circular orbit around Earth to a circular orbit around the Moon, aiming to determine the minimum ΔV required for this maneuver. The advantages of using TV type constraints are also discussed for this maneuver. Section “Moon’s gravity assist maneuver” offers a second application of the TV type constraints involving a gravity assist with the Moon. The influence of the parameters on the efficiency of the maneuver is investigated. Finally, Section “Conclusions” summarizes the main conclusions drawn from this work.

## Dynamics modelling

An inertial frame of reference with coordinates $(x_t, y_t)$ centered at the Sun is defined. The Earth and the Moon rotate in circular orbits around their common barycenter, which is the center of a rotating frame of reference with coordinates $(x_b, y_b)$ at a distance $R_s$ from the Sun, where $R_s$ is the position of the Sun with respect to the barycenter. This frame rotates with the same constant angular speed $\omega$ of the Earth and the Moon with respect to the inertial frame. Additionally, its center undergoes rotation around the Sun. Both frames are co-planar, as shown in Fig. 1 with respect to the inertial frame.

The equations of motion of a satellite in the rotating frame are given by

$$ \frac{d^2 \mathbf{r}_b}{dt^2} + 2 \omega \times \frac{d\mathbf{r}_b}{dt} + \omega \times (\omega \times \mathbf{r}_b) = - \frac{\mu_e}{r_e^3} \mathbf{r}_e - \frac{\mu_m}{r^3} \mathbf{r} - \frac{\mu_s}{r_s^3} \mathbf{r}_s - \frac{\mu_s}{R_s^3} \mathbf{R}_s, \quad (1) $$

where $\mathbf{r}_b$ is the position of the body with respect to the barycenter in the rotating frame; $\omega = (0, 0, \omega)$ is the angular velocity of the rotating frame; $\mu_e, \mu_m$, and $\mu_s$ are the gravitational parameters of the Earth, Moon, and Sun, respectively; $\mathbf{r}_e, \mathbf{r}, \mathbf{r}_s$ and $\mathbf{R}_s$ are the position of the satellite with respect to the Earth, Moon, and Sun, respectively; and $\mathbf{R}_s = \mathbf{r}_b - \mathbf{r}_s$ is the position of the Sun in the rotating frame of reference. The vectors and bodies involved in the motion can be seen in Fig. 2 from the perspective of the rotating frame of reference.

The motion of the satellite must be described relative to the Moon in order to satisfy linearity conditions in the constraints. The position of the satellite can be written as a function of the position relative to the Moon using the transformation

$$ \mathbf{r}_b = \mathbf{r} + \mathbf{R}_m, \quad (2) $$

where $\mathbf{R}_m = d_2\hat{\mathbf{x}}_b$ is the position of the Moon relative to the barycenter, with $d_2 = R\mu_e/(\mu_e + \mu_m)$, where R is the distance between the Earth and the Moon. Thus, the equations of motion with respect to the Moon can be derived by replacing Eq. (2) into Eq. (1), resulting in

$$ \frac{d^2 \mathbf{r}}{dt^2} + 2 \omega \times \frac{d\mathbf{r}}{dt} + \omega \times (\omega \times \mathbf{r}) - \omega^2 \mathbf{R}_m = - \frac{\mu_e}{r_e^3} \mathbf{r}_e - \frac{\mu_m}{r^3} \mathbf{r} - \frac{\mu_s}{r_s^3} \mathbf{r}_s - \frac{\mu_s}{R_s^3} \mathbf{R}_s. \quad (3) $$
<figure><img src="image_3.png" /><figcaption>**Figure 1.** The rotating frame of reference $(x_b, y_b)$, with the Earth and the Moon, with respect to the inertial frame of reference $(x_I, y_I)$ centered at the Sun.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** The four bodies and main vectors in the rotating frame perspective.</figcaption></figure>

The components of $r$ can be written in rectangular coordinates as the pair $(x, y)$, defining a system of coordinates centered at the Moon, as illustrated in Fig. 3. Here, $r$ represents the position of the satellite with respect to the Moon. This clarification is essential, because the forthcoming constraints will be linearly defined with respect to the Moon, not the barycenter. The coordinate $r$ is defined as the distance between the satellite and the Moon, while the coordinate $\theta$ is the angle between the position vector $r$ and the y-axis, with positive values in the counterclockwise direction (refer to Fig. 3).

Thus, the polar coordinates are defined based in the following transformation from rectangular coordinates

$$ \begin{cases} \hat{r} = -\sin\theta \hat{x} + \cos\theta \hat{y}, \\ \hat{\theta} = -\cos\theta \hat{x} - \sin\theta \hat{y}. \end{cases} \qquad (4) $$
<figure><img src="image_5.png" /><figcaption>**Figure 3.** The systems of coordinates rectangular ($x_b, y_b$) centered at the barycenter, rectangular $(x, y)$ centered at the Moon, and polar $(r, \theta)$ centered at the Moon.</figcaption></figure>

Since $r$ and $\theta$ are variables dependent on time, the position, velocity, and acceleration are given by

$$ \begin{cases} r = r \hat{r}, \\ \dot{r} = \dot{r} \hat{r} + r \dot{\theta} \hat{\theta}, \\ \ddot{r} = (\ddot{r} - r \dot{\theta}^2) \hat{r} + \left(\frac{1}{r} \frac{d}{dt}(r^2 \dot{\theta})\right) \hat{\theta}. \end{cases} \quad (5) $$

In accordance with the polar coordinates $(r, \theta)$ defined above, the equations of motion of the satellite given by Eq. (3) become

$$ \begin{aligned} & \ddot{r} - r \dot{\theta}^2 - r\omega^2 - 2r\omega\dot{\theta} + \omega^2 d_2 \sin\theta + \frac{\mu_m}{r^2} + \frac{\mu_e(r - R \sin\theta)}{(r^2 - 2rR \sin\theta + R^2)^{3/2}} - P_r = 0, \\ & 2\dot{r}(\dot{\theta} + \omega) + r\ddot{\theta} + \omega^2 d_2 \cos\theta - \frac{\mu_e R \cos\theta}{(r^2 - 2rR \sin\theta + R^2)^{3/2}} - P_\theta = 0, \end{aligned} \quad (6) $$

where the independent variable is $t$ and the dependent variables are $r(t)$ and $\theta(t)$. Here, $P_r$ and $P_\theta$ are the components of the perturbation of the Sun on the satellite in, respectively, $\hat{r}$ and $\hat{\theta}$ directions, as expressed by

$$ P_r \hat{r} + P_\theta \hat{\theta} = -\frac{\mu_s}{r_s^3} \mathbf{r}_s - \frac{\mu_s}{R_s^3} \mathbf{R}_s, \quad (7) $$

where

$$ \mathbf{r}_s = (r + R_s \cos \theta_s \sin \theta - d_2 \sin \theta - R_s \sin \theta_s \cos \theta) \hat{r} + (R_s \cos \theta_s \cos \theta - d_2 \cos \theta + R_s \sin \theta_s \sin \theta) \hat{\theta}, \quad (8) $$

with magnitude $r_s$ and

$$ \mathbf{R}_s = (-R_s \cos \theta_s \sin \theta + d_2 \sin \theta + R_s \sin \theta_s \cos \theta) \hat{r} + (-R_s \cos \theta_s \cos \theta + d_2 \cos \theta - R_s \sin \theta_s \sin \theta) \hat{\theta}, \quad (9) $$

where $R_s = 1$ AU and $\theta_s$ is the angular position of the Sun relative to the rotating frame centered at the barycenter. This angle is given by $\theta_s = \omega_s t + \gamma$, where $\gamma$ and $\omega_s$ are the initial angular position and velocity of the Sun with respect to the rotating frame, respectively. The values of the parameters used in this work are shown in Table 1. For comparison purposes, these values are the same as those used in$^{6,16,19,20}$.

## TV type constraints

The constraints introduced in this paper dictate that the spacecraft travels from a point A to any location on a circle of radius $r_p$ around the Moon. Additionally, the final velocity must be such that its component in the radial direction is zero. In this section, these constraints are precisely formulated in mathematical terms and integrated into the equations of motion using TFC. Subsequently, it is demonstrated that these constraints serve to address

<table><caption>**Table 1.** The values of the parameters used in this research.</caption><tr><td>μ<sub>e</sub></td><td>3.975837768911438 × 10<sup>14</sup> m<sup>3</sup>/s<sup>2</sup></td></tr><tr><td>μ<sub>m</sub></td><td>4.890329364450684 × 10<sup>12</sup> m<sup>3</sup>/s<sup>2</sup></td></tr><tr><td>μ<sub>s</sub></td><td>1.3237395128595653 × 10<sup>20</sup> m<sup>3</sup>/s<sup>2</sup></td></tr><tr><td>ω</td><td>2.66186135 × 10<sup>-6</sup> s<sup>-1</sup></td></tr><tr><td>ω<sub>s</sub></td><td>-2.462743433827215 × 10<sup>-6</sup> s<sup>-1</sup></td></tr><tr><td>R</td><td>3.84405000 × 10<sup>8</sup> m</td></tr></table>
two distinct problems: assessing the associated costs of an Earth-to-Moon transfer and evaluating the gains (e.g. specific energy) obtained by utilizing the Moon for a gravity assist maneuver, propelling the spacecraft to higher speeds in the trans-lunar space.

## Definition

The TV (Tangential Velocity) type constraints characterize a trajectory starting from an initial point (specified by coordinates $r_0$ and $\theta_0$) at time zero and extending to any point on a circle of radius $r_f$ around the center of the system of coordinates (the center of the Moon) at time $T$. Furthermore, the velocity at this specified time is such that its radial component is zero. The mathematical definition of the TV type constraints is shown in Table 2, and a visual representation of these constraints concerning the main celestial bodies is depicted in Fig. 4. It is noteworthy that the final position is not uniquely determined, allowing for any point on a circle around the Moon. This means that, if a solution exists (and a convergence is reached by the method), it will converge to any point at a distance $r_f$ from the center of the Moon such that the other constraint given by a null radial component velocity is also satisfied.

## Using TFC to embed the TV type constraints into the equations of motion

The TFC is extensively developed to address uncoupled and linear constraints on dependent variables⁸. TFC proves versatile for accommodating any number of constraints of the form

<table><caption>**Table 2.** The constraints of the Two Point Boundary Value Problem (TPBVP) and the Tangential Velocity (TV) proposed in this paper.</caption><thead><tr><th rowspan="2">Type</th><th rowspan="2">Time</th><th colspan="2">Position</th><th colspan="2">Velocity</th></tr><tr><th>Radial</th><th>Angular</th><th>Radial</th><th>Angular</th></tr></thead><tbody><tr><td rowspan="2">TPBVP</td><td>t = 0</td><td>r = r<sub>0</sub></td><td>θ = θ<sub>0</sub></td><td>none</td><td>none</td></tr><tr><td>t = T</td><td>r = r<sub>f</sub></td><td>θ = θ<sub>f</sub></td><td>none</td><td>none</td></tr><tr><td rowspan="2">TV</td><td>t = 0</td><td>r = r<sub>0</sub></td><td>θ = θ<sub>0</sub></td><td>none</td><td>none</td></tr><tr><td>t = T</td><td>r = r<sub>f</sub></td><td>none</td><td>dot = 0</td><td>none</td></tr></tbody></table>

<figure><img src="image_6.png" /><figcaption>**Figure 4.** (a) The TPBVP involves a transfer from point A to point B without constraints on the velocity. (b) The TV type constraints represent a transfer from point A to any location at a distance $r_f$ from the center of the Moon, with tangential velocity relative to the dashed circular final orbit.</figcaption></figure>
$$ \left. \frac{d^k x}{dt^k} \right|_{t_n} = x_n^{(k)}, \qquad (10) $$

where k is an integer, $t_n$ represents a specific time t value, and $x_n^{(k)}$ denotes the particular value of the k-th derivative of x evaluated at time $t_n$. For instance, the combination $n = 0$ and $k = 0, 1$ in Eq. (10) characterizes the initial value problem (for the time $t_0$), while $n = 0, 1$ and $k = 0$ in Eq. (10) characterizes the TPBVP (for the times $t_0$ and $t_1$). TFC can also handle relative constraints of the type

$$ \left. \frac{d^i x}{dt^i} \right|_{t_j} = \left. \frac{d^k x}{dt^k} \right|_{t_j}, \qquad (11) $$

where, if $i = k$, then $j \neq l$, and if $j = l$, then $i \neq k$. Finally, TFC can manage integral type constraints for certain cases of linear equations$^{21}$.

The TV type constraints at time T represent a fixed final radius around the Moon with a null tangential velocity at this point. In rectangular coordinates, these constraints are expressed as

$$ \left\{ \begin{array}{l} \sqrt{x^2 + y^2}|_{t=T} = r_f, \\ (\dot{x} x + \dot{y} y)|_{t=T} = 0. \end{array} \right. \qquad (12) $$

In rectangular coordinates, these constraints are coupled and nonlinear. Consequently, applying TFC directly to embed the TV type constraints in these coordinates is not feasible. However, a proper change of coordinates can be implemented to transform coupled and nonlinear constraints into uncoupled and linear ones$^{18}$. The coordinates $r, \theta$ defined in Section “Dynamics modelling” are adopted in this paper because the TV type constraints are linear on these variables. At time T, they become

$$ \left\{ \begin{array}{l} r|_{t=T} = r_f, \\ \dot{r}|_{t=T} = 0. \end{array} \right. \qquad (13) $$

As a result, the constraints are now uncoupled and linear, making it possible to use TFC to embed them into the equations of motion, as described next.

The constraints shown in Table 2 are linear on the dependent variables, allowing them to be incorporated into the equations of motion (Eq. 6) using TFC in polar coordinates. Any solution of the resulting equations exactly satisfies these constraints. Following the methodology in$^{8}$, the constrained expression for the dependent variable $r(t)$ is derived from a linear combination of a free function and a given number of support functions, which are linearly independent with regarding to each other and to the free functions, and whose number is equal to the number of constraints. Since there are three constraints for the radial coordinate and its derivatives, the constrained expression can be derived from

$$ r(t) = g_r(t) + \eta_0(t, g(t))s_0(t) + \eta_1(t, g(t))s_1(t) + \eta_2(t, g(t))s_2(t), \qquad (14) $$

where $g_r(t)$ is the free function associated to the r coordinate and $\eta_k(t, g(t))$ and $s_k(t)$, for $k = 0, 1, 2$, are coefficients and given support functions, respectively. For the purpose of this work, the coefficients $\eta_k(t, g(t))$ are constants in time (e.g.$^{8}$). To simplify, we also choose the support functions as elements of the canonical basis of polynomials of order 2 ($s_0(t) = 1, s_1(t) = t, s_2(t) = t^2$), according to$^{9}$. The TV type constraints shown in Table 2 for the $r(t)$ coordinates are given by

$$ \left\{ \begin{array}{l} r(0) = r_0, \\ r(T) = r_f, \\ \dot{r}(T) = 0. \end{array} \right. \qquad (15) $$

A set of three equations are generated by substituting $r(t)$ from Eq. (14) into Eq. (15). This set of equations is solved for the three coefficients $\eta_k$ for $k = 0, 1, 2$. This operation ensues the derivation of the following constrained expression for $r(t)$

$$ r(t) = g_r(t) + (-r_f t(t - 2T) + r_0(t - T)^2 - (t - T)^2 g_r(0) + t(t - 2T)g_r(T) + tT(-t + T)g'_r(T))/T^2, \qquad (16) $$

where $g'_r(T)$ represents the derivative of the free function $g_r(t)$ with respect to $t$ evaluated at $t = T$. Analogously, the constrained expression for the $\theta$ variable can be derived as

$$ \theta(t) = g_{\theta}(t) + \theta_{0} - g_{\theta}(0), \qquad (17) $$

where $g_{\theta}(t)$ is the free function for the $\theta$ coordinate. The free functions are written as

$$ \left\{ \begin{array}{l} g_r(t) = \sum_{j=3}^{m} \xi_j h_j(\tau), \\ g_{\theta}(t) = \sum_{j=1}^{m} \zeta_j h_j(\tau), \end{array} \right. \qquad (18) $$
where $\xi_j$ and $\zeta_j$ are unknown coefficients, $h_j$ are the orthogonal Chebyshev polynomials (e.g.$^{22}$), and $m$ is its highest order truncation. It is important to note that the orthogonal Chebyshev polynomials are valid only in the interval $-1 \le \tau \le 1$, hence the change of the independent variable from $t$ to $\tau$ given by $\tau = 2t/T - 1$ must be performed in order to satisfy this time interval. Note that the sum starts at $j = 3$ and $j = 1$ for $g_r$ and $g_\theta$, respectively, because the terms of the polynomials must be linearly independent with the support functions, and there are three support functions (for the three constraints given by Eq. (15)) for $r$ (and its derivative) and one only support function (for the constraint $\theta(0) = \theta_0$) for $\theta$.

The free functions shown in Eq. (18) are used in the expressions for the $r$ and $\theta$ variables shown in Eqs. (16) and (17). After this, the dependent variables $r$ and $\theta$ in Eq. (6) are replaced by the respective ones given by Eqs. (16) and (17). This procedure ensues in a set of two unconstrained differential equations. This set is discretized through the collocation points method using the Chebyshev-Gauss-Lobatto nodes$^{23}$ distributed between 0 and $T$ as follows

$$t_k = \left(1 - \cos\left(\frac{k\pi}{N}\right)\right)^{\frac{T}{2}}, \quad \text{for } k = 0, 1, 2, \dots, N, \qquad (19)$$

where $(N + 1)$ is the number of equations generated by this procedure for each equation in Eq. (6). These equations are then solved for the unknown coefficients $\xi_j$ for $j = 3, 4, \dots, m$ and $\zeta_j$ for $j = 1, 2, \dots, m$ using the nonlinear least squares optimization method$^{24}$ to minimize the sum of the squares of the left sides of Eq. (6) under the distribution given by Eq. (19). Further details on the TFC procedure can be seen in$^{8,24}$.

## Comparison with the TPBVP

The problem represented by the TPBVP is to transfer a spacecraft from point A around the Earth to point B around the Moon in a given transfer time $T$, as seen in Fig. 4a. The constraints of the type TPBVP are shown in Table 2. Although two positions are specified at two times, there is no constraint on the velocity.

In this paper, it is proposed the use of the TV type constraints for this problem. These constraints are identical to the traditional TPBVP at the initial time ($t = 0$), but at the final time ($t = T$), a constraint is applied in the radial component of the final velocity $\dot{r}(t)|_{t=T}$, instead of applying the constraint at the angle $\theta(T)$. A comparison between the TV and TPBVP types of constraints is shown in Table 2. The main distinction between the TV and TPBVP types is that the final velocity and position are perpendicular in the TV type, i.e. $(r \cdot \dot{r})|_{t=T} = 0$. Since the system of coordinates is centered at the Moon, the final velocity is tangential to a final circular orbit of radius $r|_{t=T} = r_f$ to which the spacecraft is transferred.

As was mentioned earlier, for maximum efficiency of the burn, tangential burn occurs only at apogee or peri-gee on elliptical orbits$^{1,2}$. This result is valid in a two-body problem local approximation$^3$. In a three body-problem context, the impulse applied in the same direction as the velocity relative to the rotating frame maximizes the change in the Jacobi constant$^5$. In a restricted four body-problem context, a tangential impulse is assumed in the numerical procedure given by the shooting method to evaluate the costs for Earth-to-Moon transfers in$^6$, which is adopted also using lunar gravity assist$^7$. A tangential burn applied at perigee or apogee is represented by the tangential velocity constraint (included in the TV type constraints) and it is analytically integrated into the equations of motions using TFC. Hence, any solution numerically obtained using the procedures explained in this section must satisfy this constraint, which is not present in the TPBVP. Besides that, the tangential velocity in the two-body context means that this point is the periapsis of the orbit around the Moon, so there is no risk of collision of the spacecraft with the Moon. In fact, the motion close to the Moon for short periods of time can be approximated by the two-body problem. The TV type constraints can be visualized in Fig. 4b with the dashed circular initial and final orbits of the transfer. Note that the point B is not defined in the final orbit because there is no constraint in the final angle $\theta(T)$.

## Earth-to-Moon transfer using the TV type constraints

The task in this section is to transfer a spacecraft from an initial circular orbit around the Earth to another circular orbit around the Moon within a specified time transfer duration $T$ using bi-impulsive maneuvers. A first impulse $\Delta V_1$ is applied at the initial orbit around the Earth to initiate the transfer, and another impulse $\Delta V_2$ is applied at the final time $t = T$ to circularize the orbit around the Moon. The spacecraft is transferred under the gravitational influence of the Earth, Moon, and Sun. This problem can be formulated as a TPBVP. Several methods are available to solve it, such as the patched restricted three-body problem$^{25,26}$, gradient shooting method$^{13}$, direct transcription and multiple shooting method$^6$, and Jacobi integral variational equation$^{27}$. The bi-circular restricted four-body problem has been commonly employed for such transfers in numerous works$^{6,7,20,28-33}$. Furthermore, the TPBVP is solved using TFC for the Earth-to-Moon transfer problem based on the bi-circular restricted four-body problem in$^{16}$. In this paper, an alternative approach is proposed by adopting the TV type constraints outlined in Section “TV type constraints”, instead of the TPBVP, to address the Earth-to-Moon transfer problem. The efficiency of the solution is enhanced by the TFC procedure described in Section “Using TFC to embed the TV type constraints into the equations of motion”. The influence of the Sun can be neglected for preliminary works on Earth-to-Moon transfers$^{16}$. However, the influence of the Sun will be considered in this section for the purpose of comparison with previous works.

## Advantages of the TV type constraints

The initial position must be specified for both the TPBVP and the TV types of constraints. On the other hand, the angle of the final position must be provided only when using TPBVP, because the final angle is unconstrained
in the TV type constraints. The TV type constraints show advantages in comparison with the TPBVP type in evaluating the costs for transfers between the Earth and the Moon. They are:

* The value of $\theta(T)$ satisfying $\dot{r}|_{t=T} = 0$ (to minimize costs) is directly provided by the solution (by the convergence of the optimization method) when using TV type constraints. On the other side, for the TPBVP constraints, an additional optimization procedure must be implemented in order to find the angle of the position satisfying $\dot{r}|_{t=T} = 0$. Hence, the numerical solution provided by using the TV type is faster than TPBVP type. A quantitative comparison will be shown in Section “Computational efficiency”.

* There is no need to use patched conics to solve the transfer problem. In patched conics, the motion is cut out in different sections, and each of them corresponds to a two-body approximation, where the main body is the one which has the highest gravitational influence on the spacecraft.

* The use of integrators (such as Runge-Kutta, etc.) in the code is unnecessary, as required by shooting methods to solve TPBVP type constraints. There is no need to iterate the initial conditions multiple times to match boundary conditions.

* The solution ensures no collision with the Moon. Assuming that the velocity of the satellite in the rotating frame is approximately equal to its relative velocity with the Moon in an inertial frame centered at the Moon (in a two-body problem context, as shown in³⁴), the condition $\dot{r}|_{t=T} = 0$ implies that $r(T)$ is locally the periapsis (or apoapsis) of the orbit. In a two-body problem context, the spacecraft has zero radial velocity only at the periapsis or apoapsis. The position $r(T)$ cannot be the apoapsis, since the spacecraft comes from a region out of the sphere of influence of the Moon.

* Since the risk of collision with the Moon is null, there is no need to choose a large $r_f$, contrary to TPBVP.

## Results for a bi-impulsive maneuver to solve the Earth-to-Moon transfer

This subsection presents the costs associated to the use of bi-impulsive maneuvers to transfer a spacecraft from an initial circular orbit of 167 km altitude around the Earth to a final orbit of 100 km altitude around the Moon (considering a radius of the Earth and Moon of 6378 km and 1738 km, respectively). The total cost $\Delta V = \Delta V_1 + \Delta V_2$ is the sum of the impulse $\Delta V_1$ applied at the initial orbit around the Earth and the one applied at the final orbit around the Moon $\Delta V_2$. For comparison purposes, these costs are evaluated in this paper using both the TPBVP and TV types of constraints defined in Section “TV type constraints”. Note that the solutions obtained using these constraints are deduced for a given initial point A, which can be anywhere in the initial orbit around the Earth (see point A in Fig. 4). Hence, an additional procedure must be implemented in order to find the best point A related to the lowest cost $\Delta V$, as described next:

* Simulations are launched for each position $(r_0, \theta_0)$ of point A in the initial orbit of 167 km altitude around the Earth (Fig. 4) spanned with a first rough step.

* For each position $(r_0, \theta_0)$, the cost $\Delta V$ is computed from the simulation.

* The position $(r_0, \theta_0)$ presenting the minimum value of $\Delta V$, denoted as $(r_0, \theta_0)_{min1}$, is determined.

This procedure is repeated, but this time in a smaller interval of the orbit centered on $(r_0, \theta_0)_{min1}$ with a smaller spanning step, in order to determine a more accurate position $(r_0, \theta_0)_{min2}$ minimizing the value of $\Delta V$, and so on. The procedure is then repeated several times by successively decreasing the step until a new iteration of the procedure no longer modifies the position $(r_0, \theta_0)$ minimizing $\Delta V$ within a given accuracy. Simultaneously, the same procedure is applied for the initial position of the Sun, indicated by the angle $\gamma$. The TV type constraints do not require an additional procedure to find the best point B (with the lowest cost $\Delta V$) to apply the second impulse for the final circular orbit. However, if the TPBVP type is adopted, an additional procedure must also be implemented to find the best point B on the final orbit.

The $\Delta V$ costs evaluated using the TV type constraints are depicted in Fig. 5a as a function of the flight time $T$, where $3 \le T \le 6.8$ days with a step of 0.01 day. These costs consider the gravitational influence of the Sun. The lowest cost found in this range is 3945.6619 m/s, associated with a flight time $T = 4.59$ days. This cost aligns with other works in the literature that also consider the solar perturbation for this Earth-to-Moon transfer. For instance, the minimum cost of 3949.53 m/s associated with $T = 4.59446$ days is shown in²⁰, the minimum cost of 3944.8 m/s associated with $T = 4.6$ days in⁶, and the minimum cost of 3944.83 m/s associated with $T = 4.625$ days in¹⁶. In comparison, a Hohmann-like transfer requires 5 days with a total cost of 3954 m/s⁶.

The same $\Delta V$ costs shown in Fig. 5a were also evaluated using the TPBVP type constraints, but this time with the aid of the additional procedure to find the position of the point B in the final orbit (the value of $\theta_f$) that minimizes the $\Delta V$. The differences between the $\Delta V$ evaluated using the TV and TPBVP types of constraints are shown in Fig. 5b. Note that these differences are very small. The closer to the minimum $\Delta V$, the smaller this difference. Such a difference around $T \approx 4.59$ days could not be properly evaluated in this narrow region because it is smaller than the machine level error involved in the numerical procedures. This is an important result, since the very small differences between the TV and the TPBVP types of constraints shown in Fig. 5b mean that the TV type is almost as accurate as the TPBVP type to minimize the $\Delta V$ costs.

## Computational efficiency

In this subsection, a comparison between the traditional TPBVP and the proposed TV types of constraints is done concerning the numerical efficiency of the methods. The Python language with the assistance of the TFC module³⁵ is chosen to numerically evaluate the nonlinear least squares method using singular-value decomposition in a just-in-time (JIT) compiler³⁶. This procedure is adopted to solve both the TPBVP and TV types of
<figure><img src="image_7.png" /><figcaption>**Figure 5.** The costs of Earth-to-Moon transfers evaluated using the TV type constraints as function of flight time *T* are shown in (a). The difference between these costs evaluated using the TV and the TPBVP types of constraints is shown in (b).</figcaption></figure>

constraints by embedding the respective constraints into the equations of motion (Eq. 6). In order to determine the coefficients $\xi_j$ and $\zeta_j$ of the free functions (Eq. 18), the TFC method proceeds iteratively with a nonlinear least squares procedure to best fit the trajectory with the dynamical equations. An initial guess is then necessary. The first initial guess adopted is a trajectory with a constant velocity between a point with a distance of 167 km from the surface of the Earth and a point with an altitude of 100 km from the surface of the Moon. After a first iteration, a new set of coefficients is obtained and it is used as an initial guess for the procedure explained in Section “Using TFC to embed the TV type constraints into the equations of motion”. The procedure is then solved for a flight time $T = 1$ day, and the solution obtained is used as an initial guess for the following flight time, increased by the step used in this work, i.e., 0.01 day.

The average times to find the minimum costs for a given flight time *T* are shown in Table 3. They are evaluated from an ensemble of 100 cases each. The TV type constraints show a reduced computational time of 31.39% in comparison with the TPBVP type of constraint for the 4BP case (considering the perturbation of the Sun), and also a reduced computational time of 45.35% for the 3BP case. Note that the TFC method to obtain the solution is very fast - a solution is usually obtained in less than one second for a given set of constraints. The major computational cost is the additional optimization procedure, which must be implemented in the code to find the position of the point A of the initial orbit, the initial best position of the Sun (the angle $\gamma$ in the case of the 4BP), and/or the position B of the final orbit around the Moon in the case of the TPBVP. This last optimization procedure (to find the final position around the Moon) is not needed in case the TV type constraints are used. This is the main advantage of the TV type in comparison with the TPBVP type constraints.

A residual is defined as the sum of the magnitude of the left side of Eq. (6) allocated for the $N + 1$ different times $t_k$ given by Eq. (19). The maximum value of the order of magnitude of the residual divided by N is of the order of $10^{-11}$ (m/s²) for all cases among the TV or TPBVP for the 3BP or 4BP.

## Moon's gravity assist maneuver

A gravity assist maneuver can be completely described by the following three parameters³⁴, which are also shown in Fig. 6:

* the periapsis distance $r_p$,

* the magnitude of the velocity at the periapsis $v_p$,

<table><caption>**Table 3.** Average time to find the minimum $\Delta V$ cost for a given flight time. The number of points discretized in the time domain and the number of basis functions are such that $N = 400$ and $m = 396$. The transfer using the TPBVP and the TV types of constraints are both solved using the TFC procedure.</caption><thead><tr><th rowspan="2">Type of constraints</th><th colspan="2">Computational time (s)</th></tr><tr><th>4BP</th><th>3BP</th></tr></thead><tbody><tr><td>TPBVP</td><td>49.6575</td><td>38.4789</td></tr><tr><td>TV</td><td>34.0725</td><td>21.0278</td></tr></tbody></table>
<figure><img src="image_8.png" /><figcaption>Figure 6. Parameters of a gravity assist with the Moon (following$^{34}$). P is the periapsis.</figcaption></figure>

• the close approach angle $\theta_p$.

The efficiency of the gravity assist can also be quantified by the velocity gain which corresponds to the increase of velocity of the spacecraft with respect to the Earth-Moon barycenter due to the gravity assist maneuver with the Moon. The gain in velocity is defined as

$$DV_b = \|V_f - V_i\|, \qquad (20)$$

where $V_i$ and $V_f$ represent the velocities of the spacecraft at the sphere of influence of the Moon in an inertial frame of reference centered at the barycenter of the Earth-Moon system before and after the gravity assist maneuver, respectively. Assuming that the motion of the spacecraft close to the Moon can be approximated by the two-body problem dynamics, this gain is$^{34}$

$$DV_b = 2v_{\infty} \sin \delta, \qquad (21)$$

where the angle $\delta$ (Fig. 6) is defined by

$$\sin \delta = \frac{1}{1 + \frac{r_p v_{\infty}^2}{\mu_m}}, \qquad (22)$$

and $v_{\infty}$ is the norm of the relative velocity of the spacecraft with respect to the Moon when approaching it. The constant energy of the two-body problem is used to connect this velocity to the velocity at the periapsis, according to$^{34}$

$$v_{\infty}^2 = v_p^2 - 2 \frac{\mu_m}{r_p}. \qquad (23)$$

The gain $DV_b$ shows the change in velocity due to the gravity assist, but cannot indicate if the velocity increases or decreases. Hence, another gain in velocity can be defined by

$$DV_g = \|V_f\| - \|V_i\|. \qquad (24)$$

It measures the gain in the norm of the velocity with respect to the Earth, and can be deduced from the expression of the velocity given in$^{34}$ with

$$\|V_f\| = \sqrt{v_{\infty}^2 + V_2^2 - 2v_{\infty}V_2 \sin(\theta_p + \delta)}, \qquad (25)$$

$$\|V_i\| = \sqrt{v_{\infty}^2 + V_2^2 - 2v_{\infty}V_2 \sin(\theta_p - \delta)}, \qquad (26)$$

where $V_2 = d_2\omega$ is the velocity of the Moon with respect to the barycenter. The specific energy gain is given by$^{34}$

$$DE = -2V_2v_{\infty} \cos \theta_p \sin \delta. \qquad (27)$$
Note that the notation used by³⁴ to characterize the gravity assist maneuver is $γ$, which is the angle between the periapsis and the x axis. In this paper, $θ_p$ is the angle between the periapsis and the y axis. This choice enables a direct relation with $θ$ defined in Section “Dynamics modelling” and shown in Fig. 3 for the case where the position is such that $θ(T) = θ_p$. Hence, the relation between the angle adopted in this paper and the angle used in³⁴ is $θ_p = γ - 90°$. The relation $sin γ = cos θ_p$ is used to write Eq. (27) obtained from the respective energy gain shown in³⁴. Note that from this equation, the specific energy decreases for $-90° < θ_f < 90°$, and increases for $90° < θ_f < 270°$, with a maximum at $180°$.

## Using the TV type constraints to perform a gravity assist maneuver

The first parameter of interest in a gravity assist maneuver is the distance at the perigee $r_p$. In this maneuver, we propose that $r_p$ be given by the final position of the transfer. The connection between the perigee and the TV type constraint is

$$r_p = r_f. \qquad (28)$$

The maneuver is then evaluated as a function of this distance, which must be provided. The following parameter of the transfer is the magnitude of the velocity at periapsis $v_p$. In a two-body problem, the radial component of the velocity at the periapsis is null, hence, this magnitude is evaluated only with the tangential component $v_p = r_p\dot{θ}_p$. Since $r_p$ is provided, the parameter of the gravity assist maneuver becomes $\dot{θ}_p$, instead of $v_p$. The connection between this parameter and the TV type constraints is

$$\dot{θ}_p = \dot{θ}_f. \qquad (29)$$

Hence, the velocity at the periapsis is $v_p = r_f\dot{θ}_f$. The last parameter of the gravity assist maneuver is $θ_p$. The connection between this parameter and the TV type constraints is

$$θ_p = θ_f. \qquad (30)$$

## Advantages

The constraints of the type TV show the following advantages in comparison with the type TPBVP:

* In the procedure proposed in this manuscript to perform the gravity assist maneuver, the parameters $\dot{θ}_p$ and $θ_p$ are obtained by the solution (through the convergence) of the optimization method.

* It is not necessary to use patched conics to solve the transfer problem.

* There is no need to use optimization procedures to satisfy the main constraint of the gravity assist maneuver, which is $\dot{r}|_{t=0} = 0$ for a given flight time, because it is embedded into the equations using TFC.

* Once again, since $r_f$ is the periapsis ($r_f = r_p$), the solution does not show any collision with the Moon.

* Since a collision with the Moon is not possible, there is no need to choose a large $r_f$ (in comparison with methods to solve the TPBVP) in order to avoid solutions that collide with the Moon.

* The periapsis distance $r_p$, which is one of the three parameters of the gravity assist, is specified in the constraint. It is then easier to study its influence on the trajectory, and to determine what initial conditions must be chosen in order to have a gravity assist with a given value of periapsis distance. In the TPBVP, the periapsis distance can only be deduced a posteriori from the obtained solutions as the two other parameters (the magnitude of the velocity at the periapsis $v_p$, the close approach angle $θ_p$). With the TV type constraints, the periapsis is imposed and the two other parameters are obtained from the solution.

## Methodology

In the gravity assist maneuver, because it is not necessary to put the spacecraft in an orbit around the Moon, only one impulsive maneuver is applied, contrary to the transfer maneuver where two impulses are used. The spacecraft starts from a circular orbit at an altitude of 167 km around the Earth, and a unique impulse characterized by $ΔV$ is then applied so that the spacecraft heads towards the Moon for a gravity assist. In Section “Earth-to-Moon transfer using the TV type constraints”, the influence of the Sun was considered for the purpose of comparison with previous works on Earth-to-Moon transfers, although its influence can be neglected. However, for the lunar gravity assist maneuver, the influence of the Sun will be not considered in order to simplify the problem. The terms due to its perturbation are set to zero ($P_r = 0$ and $P_θ = 0$) in the equations of motion given by Eq. (6), and the problem is reduced to the CR3BP (with the inertial frame centered at the Earth-Moon barycenter for simplicity). The influence of the Sun is relatively weak for a short time (less than 10 days) orbit transfer from the Earth to the Moon. For instance, the influence of the Sun can reduce the costs of an Earth-to-Moon transfer by about 5 m/s for short time travels¹⁶. Although this value may be significant for a mission, it represents a small percentage of the total cost. Thus, it is expected that the perturbation of the Sun does not change the general behavior of the gain for a gravity assist maneuver, as was noted, for instance, in³⁷.

Like in Section “Results for a bi-impulsive maneuver to solve the Earth-to-Moon transfer”, an additional optimization procedure is performed in order to find the most efficient trajectory, where, from the initial impulse $ΔV$, the final velocity after the gravity assist $V_f = \|V_f\|$ is as high as possible. Here, the optimization with TFC gives a solution for the trajectory between the starting point A and the periapsis, and then all three parameters of the gravity assist maneuver are known (the periapsis distance $r_p$, the magnitude of the velocity at the periapsis $v_p$, the close approach angle $θ_p$). From these parameters, the velocity $V_f$ is then estimated from Eq. (25), which was obtained from a two-body approximation. In order to obtain an efficient gravity assist, the procedure then
minimizes $\Delta V/V_f$, which represents the relative velocity gain. The optimization procedure used for this is the following:

* The initial circular orbit of 167 km of altitude around the Earth is spanned with a rough first step, and simulations are carried out for each of the initial positions $(r_0, \theta_0)$ obtained.

* For each position $(r_0, \theta_0)$, the ratio $\Delta V/V_f$ is computed from the simulation.

* The position $(r_0, \theta_0)$, that provides the lowest value of $\Delta V/V_f$, labeled $(r_0, \theta_0)_{min}$, is then determined.

This operation is then repeated several times, as it was done in Section “Results for a bi-impulsive maneuver to solve the Earth-to-Moon transfer”, in order to refine the position $(r_0, \theta_0)_{min}$, always working with smaller intervals around this value and smaller steps. It is important to note that the specific position of the periapsis represented by $\theta_p$ in Fig. 6 is obtained by the convergence of the procedure. The use of the TV type constraints allows the code to converge to a value of $\theta_p$ that minimizes the objective function for a given flight time. This is very advantageous for a gravity assist maneuver, since the $\theta_p$ value is not known a priori.

The procedure used for the gravity assist is then identical to the one used for the orbit transfer in Section “Results for a bi-impulsive maneuver to solve the Earth-to-Moon transfer”. The only difference is the quantity to minimize. In the same way, the computational efficiency performance is identical to the one presented in Section “Computational efficiency”.

## Results

As was explained in Section “Advantages”, one of the advantages of the TV type constraints is the ability to choose the periapsis distance $r_p$. This makes it easier to observe its influence on the gravity assist since it is not deduced a posteriori from the solution. For clarity, a new variable related to the altitude of the perigee is defined as $a_p = r_p - r_m$, where $r_m = 1738$ km is the radius of the Moon. The procedure explained in the last subsection is performed for several values of the altitude of the perigee, given by $a_p = 50, 100, 150, 200, 500, 1000, 2000, 5000, 10000$ km, imposed by the constraints for which the flight time varies along $1 \le T \le 9.99$ days with a step of 0.01 day. An example of trajectory obtained with this procedure is shown in Fig. 7 for $T = 4.58$ days and $a_p = 100$ km using the frame of reference centered at the barycenter. The velocity with respect to the inertial frame is in red, and the dashed black line represents the escape velocity from the Earth-Moon system defined as $\sqrt{2(\mu_e + \mu_m)/r_b}$, where $r_b$ is the distance between the spacecraft and the barycenter.

Therefore, for each value of $a_p$, the evolution of the different parameters and quantity can be obtained. For instance, the evolution of the non-dimensional quantity given by the ratio $\Delta V/V_f$ and minimized for the lunar

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Trajectory obtained with the TV type constraints in TFC for a flight time $T = 4.58$ days and a periapsis altitude of $a_p = 100$ km in the inertial frame of reference associated with the Earth-Moon barycenter. In the inset, the red curve represents the magnitude of the velocity of the spacecraft, while the dashed black curve is the escape velocity from the Earth-Moon system. As the motion of the Earth is very close to the Earth-Moon barycenter, only the initial position of the Earth is indicated.</figcaption></figure>
gravity assist (Section “Methodology”) is represented in Fig. 8. For, $a_p \le 2000$ km, this quantity is lower than 6, but this ratio increases for higher $a_p$, in particular for large times of flight. In Fig. 9 the evolution of $V_f$ with respect to $\Delta V$ is represented for the different values of $T$. From Fig. 9, it can be deduced the specific impulse $\Delta V$ necessary for a given flight time in order to obtain a requested velocity $V_f$ after the gravity assist. As expected, as the periapsis altitude $a_p$ decreases, the velocity increases. The largest values of the velocity $V_f$ are obtained for short flight times ($T \le 3$ days), but they require a higher impulse $\Delta V$. For longer flight times ($T > 3$ days), a given value of $\Delta V$, and, consequently, a given quantity of propellant can lead to two different values of $V_f$, the highest one being given by the shortest flight time among the two possibilities. This leads to the presence of a minimum value of $\Delta V$ and a corresponding flight time $T$ for each value of the periapsis altitude, whose values are given in Table 4. Depending on the periapsis altitude, the minimum $\Delta V$ varies between 3131.4447 and 3134.6159 m/s, and the corresponding $T$ between 4.45 and 4.58 days.

The velocity gains, $DV_b$, defined in Eq. (21), and $DV_g$, defined in Eq. (24), and the energy gain, $DE$, defined in Eq. (27), are also shown in Figs. 10, 11, 12, for each value of $a_p$, respectively. From Figs. 11 and 12, it can be verified that the gravity assist allows the spacecraft to gain energy and then to increase its velocity. The evolution of $DE$ presents a maximum different for each periapsis altitude $a_p$, whose value and corresponding flight time $T$ are indicated in Table 5. With respect to $a_p$, the maximum $DE$ varies between 0.4843 and 1.6717 (km/s)$^2$, for a corresponding $T$ between 2.05 and 2.90 days. However, it can be observed that, for higher values of $T$, $DE$ can be negative for $a_p \ge 2000$ km. These negative values $DE$ should correspond to the observed negative values of $DV_g$ and then lead to a decrease in velocity. Nevertheless, this decrease in velocity cannot be indicated by the velocity gain, $DV_b$, which only enables to observe the absolute variation of velocity. On the other hand, the peak of the velocity gain, $DV_b$, matches with a reasonable approximation with the peak of the energy gain, $DE$, for $a_p \le 2000$ km. In Fig. 13, the evolution of the angle at the periapsis $\theta_p$ is represented with respect to the flight time. The maximum gain in energy is for $\theta_p = 180$ deg, according to Eq. (27), which coincides with flight times of about 2.2 days and values of $a_p$ up to 200 km. The angle can be below 90 deg for large flight times if $a_p \ge 2000$ km.

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Ratio $\Delta V/V_f$ for a gravity assist maneuver around the Moon, evaluated using the TV type constraints as a function of time for different periapsis altitudes $a_p$. The initial point is at an altitude of 167 km around the Earth.</figcaption></figure>

<table><caption>**Table 4.** Minimum value of the velocity impulse $\Delta V$ required and corresponding flight time for each value of the periapsis altitude.</caption><thead><tr><td>r<sub>p</sub> (km)</td><td>ΔV (m/s)</td><td>T (day)</td></tr></thead><tbody><tr><td>50</td><td>3134.6159</td><td>4.58</td></tr><tr><td>100</td><td>3134.5947</td><td>4.58</td></tr><tr><td>150</td><td>3134.5736</td><td>4.57</td></tr><tr><td>200</td><td>3134.5527</td><td>4.57</td></tr><tr><td>500</td><td>3134.4310</td><td>4.57</td></tr><tr><td>1000</td><td>3134.2389</td><td>4.56</td></tr><tr><td>2000</td><td>3133.8823</td><td>4.55</td></tr><tr><td>5000</td><td>3132.9201</td><td>4.51</td></tr><tr><td>10000</td><td>3131.4447</td><td>4.45</td></tr></tbody></table>
<table><caption>**Table 5.** Maximum value of the energy gain *DE* obtained and flight time corresponding for each value of the periapsis altitude.</caption><thead><tr><td>r<sub>p</sub> (km)</td><td>DE (km/s)<sup>2</sup></td><td>T (day)</td></tr></thead><tbody><tr><td>50</td><td>1.6717</td><td>2.05</td></tr><tr><td>100</td><td>1.6487</td><td>2.07</td></tr><tr><td>150</td><td>1.6266</td><td>2.09</td></tr><tr><td>200</td><td>1.6053</td><td>2.11</td></tr><tr><td>500</td><td>1.4909</td><td>2.22</td></tr><tr><td>1000</td><td>1.3389</td><td>2.37</td></tr><tr><td>2000</td><td>1.1199</td><td>2.57</td></tr><tr><td>5000</td><td>0.7553</td><td>2.82</td></tr><tr><td>10000</td><td>0.4843</td><td>2.90</td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Velocity *V*<sub>f</sub> after the gravity assist maneuver with respect to the initial impulse Δ*V*, and for different periapsis altitudes *a*<sub>p</sub>. For each value of *a*<sub>p</sub>, the curve is obtained with different values of the flight time *T*. The integer values of *T* in days are indicated for each curve with a special symbol.</figcaption></figure>

This was expected because, when the gravity assist maneuver leads to a decrease in energy, as this is the case for these conditions (Fig. 12), $\theta_p$ should be below 90 deg. In Fig. 14, the velocity gain $DV_g$ is represented as a function of the initial impulse $\Delta V$. It can be noted that $DV_g$ presents a maximum around the flight time of 4 days, whose exact value depends on the periapsis altitude. It can then be observed that, for values of the velocity impulse close to the minimum (see Table 4), it is possible to slow down the spacecraft at the minimum or to speed it up at the maximum, depending on the objectives of the mission, by choosing suitable values for the flight time and periapsis altitude. Therefore, the minimum quantity of propellant with only a small variation can be used for all the possible gains in velocity.

## Conclusions

This paper proposed the use of the tangential velocity (TV) type of constraints in the Theory of Functional Connections (TFC) in order to find possible trajectories of orbital maneuvers. The TFC method enables to include linear constraints, as initial or final conditions into a functional interpolation of the solution. Therefore, any solution given by TFC will always analytically satisfy the constraints. In the example of an orbit transfer from a circular orbit around the Earth to a circular orbit around the Moon, the constraints given to TFC in the case of the Two Point Boundary Value Problem (TPBVP) are the initial and final positions. However, for such a problem, the choice of a tangential impulse enables to minimize the cost. Therefore, imposing a tangential velocity constraint in TFC enables a direct search for solution solutions minimizing the transfer cost. In this paper, it was shown that this can be only done in polar coordinates in order to keep the constraints linear.
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Velocity gain $DV_b$ of a gravity assist maneuver around the Moon, evaluated using the TV type constraints as a function of time for different periapsis altitudes $a_p$. The initial point is at an altitude of 167 km around the Earth.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 11.** Velocity gain $DV_g$ of a gravity assist maneuver around the Moon, evaluated using the TV type constraints as a function of time for different periapsis altitudes $a_p$. The initial point is at an altitude of 167 km around the Earth.</figcaption></figure>

As a result, using the TV type constraints in TFC significantly reduces the computational time. Indeed, for the example of the orbit transfer from the Earth to the Moon, if the final velocity is imposed as tangential, only one optimization procedure is necessary to find the best position of the spacecraft on its initial orbit, as the final position is directly given by the solution. In the TPBVP, two optimization procedures are necessary, one for the initial position, and another for the final one. Using the TV type constraints in TFC also avoids the use of patched conics approximation and integrators like Runge-Kutta required by shooting methods. Finally, this also avoids the risk of collision with the final body, because the final position given by the TFC solution corresponds to the periapsis of the trajectory around the final body, thus making it possible to choose a small radius for the final position.

The use of the TV type constraints in TFC is applied in this paper in the cislunar space for two maneuvers: an orbit transfer from a circular orbit around the Earth to a circular orbit around the Moon and a lunar gravity assist maneuver. For the orbit transfer, where the objective is to find the trajectory minimizing the quantity of propellant required, it has been shown that the TV type constraints yield identical results to those of the TPBVP with a significant decrease in computational time. Concerning the lunar gravity assist, the TV type constraints
<figure><img src="image_14.png" /><figcaption>**Figure 12.** Energy gain *DE* of a gravity assist maneuver around the Moon, evaluated using the TV type constraints as a function of time for different periapsis altitudes *a*<sub>p</sub>. The initial point is at an altitude of 167 km around the Earth.</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 13.** Periapsis angle θ<sub>p</sub> of a gravity assist maneuver around the Moon, evaluated using the TV type constraints as a function of time for different periapsis altitudes a<sub>p</sub>. The initial point is at an altitude of 167 km around the Earth.</figcaption></figure>

enable a rapid determination of the minimum initial impulse with respect to the required final velocity (after the encounter) depending on the different parameters of the gravity assist. Indeed, the TV type constraints enable to impose the altitude of the periapsis of the gravity assist, which is one of the determinant parameters, being the two other parameters given by the solution. It can be noted that, for both applications considered here, the computation performed with TFC is identical, and that the unique difference lies in the quantities to minimize.

The choice of the cislunar space was only made in order to work in a well-known and studied framework, where the obtained results can be compared with previous reference works. However, the TV type constraints are not limited to the cislunar space, and could be used in any spatial environment. In the same way, the use of the TV type constraints is not restricted to the two maneuvers considered here, orbit transfer and gravity assist, but could be applied in any spatial environment or maneuver requiring a tangential velocity constraint. This study is also here reduced to the framework of the planar case. However, the use of the TV constraints could also be applied in a more realistic model with three dimensions. For instance, with a spherical coordinates system (r,θ,φ) centered at the Moon, the TV constraints would be still valid. However, it would be then necessary to add the new initial constraint φ(0) = φ₀. The solutions could then be obtained by TFC subject to convergence
<figure><img src="image_16.png" /><figcaption>**Figure 14.** Velocity gain $DV_g$ after the gravity assist maneuver with respect to the initial impulse $\Delta V$, and for different periapsis altitudes $a_p$. For each value of $a_p$, the curve is obtained with different values of the flight time $T$. The integer values of $T$ in days are indicated for each curve with a special symbol.</figcaption></figure>

with the constraints $r(0) = r_0$, $r(T) = r_f$, $\theta(0) = \theta_0$, $\phi(0) = \phi_0$, $\dot{r}(T) = 0$. As TFC turned out recently in several studies as an effective and promising method to design orbital maneuvers, and as a tangential velocity in orbital maneuver corresponds to a optimal solution related to the lowest propellant consumption, the use of the TV type constraints in TFC could then have many other useful applications in the future.

## Data availability

The datasets analysed during the current study can be generated by following the procedures explained in this paper. They can also be obtained from the corresponding author on reasonable request.

Received: 8 January 2024; Accepted: 24 March 2024

Published online: 29 March 2024

## References

1. Vallado, D. A. *Fundamentals of Astrodynamics and Applications* (Springer, 2007), 3 edn.

2. Lawden, D. F. Impulsive transfer between elliptical orbits. In Leitmann, G. (ed.) *Optimization Techniques*, vol. 5 of *Mathematics in Science and Engineering*, 323–351, https://doi.org/10.1016/S0076-5392(08)62099-8 (Elsevier, 1962).

3. Marchal, C., Synthèse des résultats analytiques sur les transferts optimaux entre orbites Képlériennes (Durée indifférente). In De Veubeke, B. F. (ed.) *Advanced Problems and Methods for Space Flight Optimization*, 91–156, https://doi.org/10.1016/B978-0-08-013290-7.50010-8 (Pergamon, 1969).

4. Miele, A. & Mancuso, S. Optimal trajectories for earth-moon-earth flight. *Acta Astronaut.* **49**, 59–71. https://doi.org/10.1016/S0094-5765(01)00007-8 (2001).

5. Pernicka, H. J., Scarberry, D. P., Marsh, S. M. & Sweetser, T. H. Search for low $\Delta V$ Earth-to-Moon trajectories. *J. Astronaut. Sci.* **43** (1995).

6. Topputo, F. On optimal two-impulse Earth-Moon transfers in a four-body model. *Celest. Mech. Dyn. Astron.* **117**, 279–313 (2013).

7. Qi, Y. & Xu, S. Optimal earth-moon transfers using lunar gravity assist in the restricted four-body problem. *Acta Astronaut.* **134**, 106–120. https://doi.org/10.1016/j.actaastro.2017.02.002 (2017).

8. Mortari, D. The theory of connections: connecting points. *Mathematics* **5**, 57. https://doi.org/10.3390/math5040057 (2017).

9. Leake, C., Johnston, H. & Mortari, D. *The Theory of Functional Connections: A Functional Interpolation*. Morrisville NC, Framework with Applications. (Lulu, 2022).

10. Mortari, D. & Leake, C. The multivariate theory of connections. *MDPI Math.7*. https://doi.org/10.3390/math7030296 (2019).

11. Furfaro, R. & Mortari, D. Least-squares solution of a class of optimal space guidance problems via theory of connections. *Acta Astronaut.* **168**, 92–103. https://doi.org/10.1016/j.actaastro.2019.05.050 (2020).

12. Mortari, D., Johnston, H. & Smith, L. High accuracy least-squares solutions of nonlinear differential equations. *J. Comput. Appl. Math.* **352**, 293–307. https://doi.org/10.1016/j.cam.2018.12.007 (2019).

13. Prado, A. & Broucke, R. Transfer orbits in restricted problem. *J. Guid. Control. Dyn.* **18**, 593. https://doi.org/10.2514/3.21428 (1995).

14. Enright, P. J. & Conway, B. A. Discrete approximations to optimal trajectories using direct transcription and nonlinear programming. *J. Guid. Control. Dyn.* **15**, 994–1002. https://doi.org/10.2514/3.20934 (1992).

15. Miele, A., Huang, H. Y. & Heideman, J. C. Sequential gradient-restoration algorithm for the minimization of constrained function- ordinary and conjugate gradient versions. *J Optim. Theory Appl.* **4**, https://doi.org/10.1007/BF00927947 (1969).

16. de Almeida Junior, A. K., Johnston, H., Leake, C. & Mortari, D. Fast 2-impulse non-keplerian orbit transfer using the theory of functional connections. *Eur. Phys. J. Plus* **136**, https://doi.org/10.1140/epjp/s13360-021-01151-2 (2021).
17. de Almeida Junior, A. K. & de Almeida Prado, A. F. B. Comparisons between the circular restricted three-body and bi-circular four body problems for transfers between the two smaller primaries. *Sci. Rep.* **12**, 1–19, https://doi.org/10.1038/s41598-022-08046-x (2022).

18. de Almeida Jr., A. K., Prado, A. F. & Mortari, D., Orbit transfer using theory of functional connections via change of variables. *Eur. Phys. J. Spec. Top*. https://doi.org/10.1140/epjs/s11734-023-01013-1 (2023).

19. Simó, C., Gómez, G., Jorba, À. & Masdemont, J. The bicircular model near the triangular libration points of the RTBP. In *From Newton to Chaos*, 343–370 (Springer, 1995).

20. Yagasaki, K. Sun-perturbed earth-to-moon transfers with low energy and moderate flight time. *Celest. Mech. Dyn. Astron.* **90**, 197. https://doi.org/10.1007/s10569-004-0406-8 (2004).

21. De Florio, M., Schiassi, E., D'Ambrosio, A., Mortari, D. & Furfaro, R. Theory of functional connections applied to linear ODEs subject to integral constraints and linear ordinary integro-differential equations. *Math. Comput. Appl.* **26**. https://doi.org/10.3390/mca26030065 (2021).

22. Abramowitz, M. & Stegun, I. A. *Handbook of Mathematical Functions* (Dover Publications, 1972).

23. Lanczos, C. *Applied analysis* (Courier Corporation, 1988).

24. Mortari, D. & Leake, C. The multivariate theory of connections. *Mathematics* **7**, 296. https://doi.org/10.3390/math7030296 (2019).

25. da Silva Fernandes, S. & Maranhão Porto Marinho, C. Optimal two-impulse trajectories with moderate flight time for earth-moon missions. *Math. Problems Eng.* **2012** (2012).

26. Peng, L., Wang, Y., Dai, G., Chang, Y. & Chen, F. Optimization of the Earth-Moon low energy transfer with differential evolution based on uniform design. In *IEEE Congress on Evolutionary Computation*, 1–8 (2010).

27. Filho, L. A. G. & da Silva Fernandes, S. A method based on Jacobi Integral variational equation for computing Earth-Moon trajectories in the four-body problem. *Acta Astronautica* **165**, 312–330. https://doi.org/10.1016/j.actaastro.2019.09.021 (2019).

28. Yamakawa, H., Kawaguchi, J., Ishii, N. & Matsuo, H. A numerical study of gravitational capture orbit in the earth-moon system. *Spaceflight Mech.* **1992**, 1113–1132 (1992).

29. Da Silva Fernandes, S. & Marinho, C. Sun influence on two-impulsive Earth-to-Moon transfers. In *Proceedings of the 22nd International Symposium on Space Flight Dynamics* (2011).

30. Moore, A., Ober-Blöbaum, S. & Marsden, J. E. Trajectory design combining invariant manifolds with discrete mechanics and optimal control. *J. Guid. Control. Dyn.* **35**, 1507–1525. https://doi.org/10.2514/1.55426 (2012).

31. Oshima, K., Topputo, F., Campagnola, S. & Yanao, T. Analysis of medium-energy transfers to the Moon. *Celest. Mech. Dyn. Astron*. https://doi.org/10.1007/s10569-016-9727-7 (2017).

32. Onozaki, K., Yoshimura, H. & Ross, S. D. Tube dynamics and low energy Earth-Moon transfers in the 4-body system. *Adv. Space Res.* https://doi.org/10.1016/j.asr.2017.07.046 (2017).

33. Oshima, K., Topputo, F. & Yanao, T. Low energy transfers to the Moon with long transfer time. *Celestial Mech. Dyn. Astron.* **131**. https://doi.org/10.1007/s10569-019-9883-7 (2019).

34. Broucke, R. A. *The Celestial Mechanics of Gravity Assist*. https://doi.org/10.2514/6.1988-4220 (1988).

35. Leake, C. & Johnston, H. TFC: A Functional Interpolation Framework (2021). https://github.com/leakec/tfc.

36. Bradbury, J. et al. JAX: composable transformations of Python+NumPy programs (2018).

37. Negri, R. B., Sukhanov, A. & Prado, A. F. B. D. A. Lunar gravity assists using patched-conics approximation, three and four body problems. *Adv. Space Res.* **64**, 42–63. https://doi.org/10.1016/j.asr.2019.04.002 (2019).

## Acknowledgements

DM, AKAJ and TV were partially supported by the project Centro de Investigação em Ciências Geo-Espaciais, reference UIDB/00190/2020, funded by COMPETE 2020 and FCT, Portugal. The team acknowledges further support from ENGAGE-SKA Research Infrastructure, ref. POCI-01-0145-FEDER-022217, funded by COMPETE 2020 and FCT, Portugal; IT team members acknowledge support from Projecto Lab. Associado UID/EEA/50008/2019. We acknowledge support by the European Commission H2020 Programme under the grant agreement 2-3SST2018-20, and support from Center for Mechanical and Aerospace Science and Technologies C-MAST, funded by FCT-Fundação para a Ciência e a Tecnologia, Portugal through projects https://doi.org/10.54499/UIDB/00151/2020 and https://doi.org/10.54499/UIDP/00151/2020. This work is also supported by CFisUC (UIDB/04564/2020 and UIDP/04564/2020). This work is supported by the Fundação para a Ciência e Tecnologia (FCT), Ph.D. Grant No.2022.12341.BDANA. AKAJ and TV acknowledge partial support from ATLAR internal project n. 01012024. VMO is partially supported by the São Paulo Research Foundation (FAPESP, Brazil), under Grants No. 2021/11306-0 and 2022/12785-1.

## Author contributions

AKAJ and TV conceived the idea, generated the figures, and wrote the first version of the manuscript. All authors analyzed the results and reviewed the manuscript.

## Competing interests

The authors declare no competing interests.

## Additional information

Correspondence and requests for materials should be addressed to A.K.A.

Reprints and permissions information is available at www.nature.com/reprints.

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
<figure><img src="image_17.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International

License, which permits use, sharing, adaptation, distribution and reproduction in any medium or
format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the
Creative Commons licence, and indicate if changes were made. The images or other third party material in this
article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the
material. If material is not included in the article's Creative Commons licence and your intended use is not
permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from
the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2024