## 5 Model Validation and Numerical Simulations  
The proposed procedure was tested in a set of circle-to-circle orbit, two-dimensional transfers that model some typical heliocentric (in which the primary body is the Sun and \(\mu=\mu_{\odot}=132,712,439,935.5\,\text{km}^{3}/\text{s}^{2}\)) and geocentric (in which the primary body is Earth and \(\mu=\mu_{\odot}=398,600\,\text{km}^{3}/\text{s}^{2}\)) mission scenarios. In particular, regarding a potential heliocentric mission scenario, the initial parking orbit radius was set equal to \(1\,\text{AU}\). In fact, the case of \(r_{0}=1\,\text{AU}\) models the situation in which the spacecraft begins the heliocentric phase of a typical interplanetary transfer after escaping from Earth using a parabolic (escape) orbit. In this context, three different (circular) target obits were considered; that is, three different values of \(r_{f}\) were used in the numerical simulations, namely (i) \(r_{f}=r_{\mathbb{Q}}=0.723\,\)AU, which models an Earth-Venus transfer; (ii) \(r_{f}=r_{\mathbb{Q}^{\prime}}=1.524\,\)AU, which models an Earth-Mars transfer; and (iii) \(r_{f}=r_{\mathbb{Q}^{\prime}}=5.203\,\)AU, which models an Earth-Jupiter transfer. Accordingly, in this selected set of heliocentric mission scenarios, the dimensionless value of the target orbit radius is \(\bar{r}_{f}=r_{f}/r_{0}=\{\bar{r}_{\mathbb{Q}},\bar{r}_{\mathbb{Q}^{\prime}},\bar{r}_{\mathbb{Q}^{\prime}}\}=\{0.723,1.524,5.203\}\). As regards the geocentric mission scenario, a classical (two-dimensional) LEO-GEO transfer is considered, in which the radius of the circular parking orbit is \(r_{0}=r_{\text{LEO}}=6578\,\)km (i.e., a circular LEO of 200 km of altitude is assumed) and the target orbit has a radius \(r_{f}=r_{\text{GEO}}=42,\)164 km. Therefore, in the geocentric case, the dimensionless value of the target orbit radius is \(\bar{r}_{f}=r_{\text{GEO}}/r_{\text{LEO}}\triangleq\bar{r}_{\mathbb{G}}\simeq 6.41\). The number \(n\) of completed revolutions around the primary body, in the four proposed mission scenarios, can be calculated using Equation (38) as a function of the value of \(\bar{a}_{m}\). The function \(n=n(\bar{a}_{m})\) is shown in Figure 3 when \(\bar{a}_{m}<0.05\) and \(\bar{r}_{f}=\{\bar{r}_{\mathbb{Q}},\bar{r}_{\mathbb{Q}^{\prime}},\bar{r}_{ \mathbb{Q}^{\prime}}\}\). In particular, the stepped curves in Figure 3 allow us to quickly estimate the maximum value of \(\bar{a}_{m}\) for which the proposed analytical approximation can be used to initialize the associated TPBVP. In fact, according to Figure 3, the condition \(n\geq 2\) is met when (i) \(\bar{a}_{m}\leq 0.018\) in the Earth-Venus scenario; (ii) \(\bar{a}_{m}\leq 0.011\) in the Earth-Mars scenario; (iii) \(\bar{a}_{m}\leq 0.018\) in the Earth-Jupiter scenario; and (iv) \(\bar{a}_{m}\leq 0.019\) in the LEO-GEO scenario. As a result, the validity of the proposed method was tested by considering, for each mission scenario, about 20 values of \(\bar{a}_{m}\) ranging between 0.001 and 0.02 (with a step of 0.001). Figure 3: Number \(n\) of completed revolutions around the primary body in the four analyzed mission scenarios. The value of \(n\) is calculated using Equation (38). For each value of \(\tilde{a}_{m}\), the optimal (i.e., the rapid) transfer trajectory was calculated by solving the associated TPBVP using a single shooting method [20], which was initialized according to the expressions given by Equations (26), (29), and (33). In that context, the equations of motion (6) and the Euler-Lagrange equations (15), (19), and (20) were integrated in double precision using the Adams-Bashforth method [23] with an absolute and relative tolerance of \(1\times 10^{-10}\), while the tolerance on the convergence of the TPBVP was set to \(1\times 10^{-8}\).