Citation: Sánchez-García, M.M.; Barderas, G.; Romero, P. Modelization of Low-Cost Maneuvers for an Areostationary Preliminary Mission Design. Math. Comput. Appl. 2023, 28, 105. https://doi.org/10.3390/ mca28060105 Academic Editor: Gianluigi Rozza Received: 20 September 2023 Revised: 15 October 2023 Accepted: 25 October 2023 Published: 27 October 2023 Copyright: © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). Mathematical and Computational Applications Article Modelization of Low-Cost Maneuvers for an Areostationary Preliminary Mission Design Marta M. Sánchez-García 1 , Gonzalo Barderas 1,2,* and Pilar Romero 1,2 1 U.D. Astronomía y Geodesia, Facultad de Matemáticas, Universidad Complutense de Madrid, Plaza de las Ciencias, 3, E-28040 Madrid, Spain; martsa08@ucm.es (M.M.S.-G.); pilar_romero@mat.ucm.es (P.R.) 2 Instituto de Matemática Interdisciplinar, Facultad de Matemáticas, Universidad Complutense de Madrid, Plaza de las Ciencias, 3, E-28040 Madrid, Spain * Correspondence: gonzalobm@mat.ucm.es Abstract: The aim of this paper is to analyze the determination of interplanetary trajectories from Earth to Mars to evaluate the cost of the required impulse magnitudes for an areostationary orbiter mission design. Such analysis is first conducted by solving the Lambert orbital boundary value problem and studying the launch and arrival conditions for various date combinations. Then, genetic algorithms are applied to investigate the minimum-energy transfer orbit. Afterwards, an iterative procedure is used to determine the heliocentric elliptic transfer orbit that matches at the entry point of Mars’s sphere of influence with an areocentric hyperbolic orbit imposing specific conditions on inclination and periapsis radius. Finally, the maneuvers needed to obtain an areostationary orbit are numerically computed for different objective condition values at the Mars entry point to evaluate an areostationary preliminary mission cost for further study and characterization. Results show that, for the dates of the minimum-energy Earth–Mars transfer trajectory, a low value for the maneuvers to achieve an areostationary orbit is obtained for an arrival hyperbola with the minimum possible inclination and a capture into an elliptical trajectory with a low periapsis radius and an apoapsis at the stationary orbit. For a 2026 mission with a TOF of 304 for the minimum-energy Earth–Mars transfer trajectory, for a capture with a periapsis of 300 km above the Mars surface the value achieved will be 2.083 km/s. Keywords: areostationary mission planning; Earth–Mars transfer trajectories; hyperbolic orbit matching; Lambert problem 1. Introduction Small relay satellites in areostationary orbit are considered the most efficient candidates to support the telecommunication needs in the 2020s [1–7]. Areostationary orbiters, like geostationary satellites for Earth [8,9], can provide continuous access at very high data rates to remotely supervise a significant population of probes and robotic missions on the Martian surface. The determination of transfer trajectories from Earth to Mars aimed at lowering costs in terms of impulses has become a key factor in mission planning, allowing for larger payloads to be transported at a minimum energy cost. In this work, we analyze the design of an interplanetary Earth–Mars transfer to reach the areostationary orbit with the minimum impulsive maneuvers cost. Several authors have studied the optimization of interplanetary trajectories: in [10,11] transfer trajectories to the Moon and Jupiter, respectively, passing close to a Lagrangian point, are considered; in [12], a method is developed to obtain approximate near-optimal low-thrust interplanetary transfers using solar electric propulsion spacecraft; in [13], the optimization is performed with a cost function with variable coefficients; in [14], launch constraints are imposed for the optimization. We derive the heliocentric elliptic transfer characterizing the launch windows using an heuristic optimization method for determining an optimal time of flight Math. Comput. Appl. 2023, 28, 105. https://doi.org/10.3390/mca28060105 https://www.mdpi.com/journal/mca https://doi.org/10.3390/mca28060105 https://doi.org/10.3390/mca28060105 https://creativecommons.org/ https://creativecommons.org/licenses/by/4.0/ https://creativecommons.org/licenses/by/4.0/ https://www.mdpi.com/journal/mca https://www.mdpi.com https://orcid.org/0009-0004-1017-1193 https://orcid.org/0000-0002-2156-0010 https://orcid.org/0000-0002-9657-1291 https://doi.org/10.3390/mca28060105 https://www.mdpi.com/journal/mca https://www.mdpi.com/article/10.3390/mca28060105?type=check_update&version=1 Math. Comput. Appl. 2023, 28, 105 2 of 16 (TOF) that minimizes the characteristic energies [15,16]. We will analyze the sensibility of this parameter in the optimization of impulsive maneuvers. The first step consists in solving the Lambert problem [17] for various combinations of departure and arrival dates. Departure characteristic energy and hyperbolic arrival velocity plots are usually examined to investigate possible transfer windows [18]. We use genetic algorithms [19] to simultaneously minimize these two key parameters within these launch windows, comparing their performance. Then, we match this interplanetary transfer with an entry hyperbola around Mars. The classic patched conic problem has been used to achieve a continuous trajectory composed of the trajectory between two planets and the planetocentric trajectory [20–24]. We use the iterative procedure [25] with imposed conditions on the periapsis distance, the arrival hyperbolic inclination, and a fixed radius for the Mars sphere of influence (SOI), and analyze the changes in the B-plane [26] due to the variations in the arrival asymptote direction. This iterative procedure enables the evaluation of these selected parameters in order to minimize fuel consumption for planning an areostationary mission. Once the fully matched trajectory to arrive at Mars is obtained, the maneuvers necessary to capture the orbiter and to place it in the areostationary orbit are analyzed. The paper is organized as follows. In Section 2, we describe the dynamical model of the minimum-energy launch window problem. The determination of the Earth–Mars transfer trajectory with imposed hyperbolic arrival trajectory conditions is presented in Section 3. Section 4 describes the maneuvers performed to capture the spacecraft into an areostationary orbit and the numerical simulations to evaluate these maneuvers for different conditions. Finally, in Section 5, we briefly summarize the main conclusions. 2. Minimum-Energy Launch Window for Earth–Mars Transfer Trajectories We first analyze the determination of interplanetary trajectories from Earth to Mars by minimizing the required energy at Earth departure and Mars arrival. We assume point mass gravitational forces for Earth and Mars within their respective spheres of influence and an unperturbed Keplerian orbit around the Sun. The key parameters commonly used [27,28] to analyze the Earth–Mars mission launch opportunities are the characteristic energy at departure from Earth, C3 = V2 ∞E , the hyper- bolic excess velocity to escape from Earth, V∞E , and Mars arrival hyperbolic excess velocity, V∞M . In order to obtain these two parameters, it is necessary to first solve the Lambert orbital boundary value problem for the heliocentric spacecraft position, rsS , r̈sS = −µS rsS r3 sS , (1) constrained by two points, P1 and P2, and an elapsed TOF, t2 = t1 + TOF, as illustrated in Figure 1, rsS(t1) = rE, (2) rsS(t2) = rM, (3) where rE and rM are the heliocentric position vectors for the Earth at t1 and Mars at t2, respectively. The solution to the Lambert problem results in an elliptic conic section connecting P1 and P2. We consider the shortcut solution that satisfies the boundary conditions, based on the iterative procedure, choosing the time transfer function introduced by Lancaster [29] as parameter for the iteration. The solution results in an elliptic conic section connecting P1 and P2, with departure and arrival velocities, VsS1 and VsS2 , at t1 and t2, respectively. Math. Comput. Appl. 2023, 28, 105 3 of 16 For each Earth departure and Mars arrival combination of dates, V∞E and V∞M change as VsS1 and VsS2 change according to V∞E = VsS1 −VES1 , (4) V∞M = VsS2 −VMS2 , (5) where VES1 and VMS2 are the heliocentric velocity vectors for the Earth at t1 and Mars at t2, respectively. Sun P : Mars at t 2 2 P : Earth at t 1 1 r r 1 2 Figure 1. Lambert orbital boundary value problem and the heliocentric elliptic trajectory. To analyze launch and arrival window opportunities, we first focus on data visual- ization [15,18] of the departure characteristic energy C3 and the hyperbolic arrival velocity V∞M for various combinations of departure and arrival dates. We use Matlab software avail- able from [30] to solve the Lambert problem, first obtaining a reduced launch window for the minimum-energy solution. The porkchop plots [31,32] shown in Figure 2 depict the contour lines of constant C3 in km2/s2 and V∞M in km/s for the 2019–2029 departure and the 2020–2030 arrival time frames. It is possible to observe that the launch and arrival windows that give the minimum values approximately repeat every Mars synodic period of about 780 days. In more detail, Figure 3a shows the departure characteristic energy and hyperbolic arrival velocity contour plots for the Earth–Mars transfer covering 27 months, from 1 July 2019 to 1 November 2021, and Figure 3b presents these plots for the departure window on July 2020 and the arrival window on February 2021. 1 3 13 13 14 14 14 15 15 15 15 15 16 16 16 16 20 20 20 20 20 20 25 25 25 25 2 5 25 2 5 25 25 30 30 3 0 3 0 3030 30 3 0 30 30 30 35 35 3 5 35 3 5 35 35 35 35 35 35 35 35 40 4 0 40 40 40 40 40 4 0 40 40 40 4 0 40 40 45 45 4 5 45 45 45 45 45 45 45 45 45 45 45 45 45 45 50 50 5 0 50 50 50 50 50 50 50 50 50 50 50 50 50 50 50 2.62.7 2.7 2.8 2.8 2.9 2.9 2.9 2.9 2.9 3 3 3 3 3 3 3. 5 3.5 3.5 3.5 3.5 3.5 3.5 44 4 4 4 4 4 4 5 5 5 5 55 5 5 5 5 66 6 66 6 6 66 66 7 7 7 7 7 7 77 77 77 8 8 8 8 8 8 8 8 8 8 8 8 8 8 -1500 -1250 -1000 -750 -500 -250 0 250 500 750 1000 1250 1500 Departure: days from 01-Jul-2024 -1500 -1250 -1000 -750 -500 -250 0 250 500 750 1000 1250 1500 A rr iv a l: d a y s f ro m 0 1 -J u l- 2 0 2 5 C 3 (km 2 /s 2 ) V (km/s) 8M Figure 2. Departure characteristic energy and hyperbolic arrival velocity contour plots for the Earth–Mars transfer spanning 10 years. Math. Comput. Appl. 2023, 28, 105 4 of 16 45 50 2. 7 44 5 6 7 7 8 8 8 -250 Departure: days from 01−Mar−2020 -250 -200 -150 -100 -50 0 50 100 150 200 250 A rr iv a l: d a ys f ro m 0 1 − M a r− 2 0 2 1 C 3 (km 2 /s ) 2 V M (km/s) 8 2.7 6 16 25 5 0 14 2 5 25 3 5 4 0 5 0 50 2.5 3 3 3.5 3 .5 3. 5 3 .5 4 4 5 6 7 5 20 40 35 -200 -150 -100 -50 0 50 100 150 200 250 C 3 (km 2 /s ) 2 V M (km/s) 8 −30 −20 −10 0 10 20 30 −30 −20 −10 0 10 20 30 20 2. 8 2. 9 3 3 3. 5 3. 5 14 15 16 2. 5 2 .62 .72. 82. 9 Departure: days from 15−Jul−2020 A rr iv a l: d a ys f ro m 1 5 − F e b − 2 0 2 1 a) b) Figure 3. Departure characteristic energy and hyperbolic arrival velocity contour plots for the Earth–Mars transfer (a) from July 2019 to November 2021 and (b) from July 2020 to February 2021. Now, we search for the solution minimizing the equation: C = C3 + V∞M . (6) The minimum C in Equation (6) tends to give lower values of the impulsive maneuvers required, first to obtain an Earth escape velocity, and after, at the Mars arriving hyperbolic orbit, to reduce the hyperbolic excess velocity to capture the probe. To minimize Equation (6), applied to the reduced windows previously estimated from porkchop plots, we now use a Matlab optimizer [33] that includes a library dedicated to genetic algorithms with different implementations of selection and crossover functions (see [19]). In order to analyze the accuracy of the genetic algorithms when applied to this problem, we compare the performance of the Remainder and the Stochastic Uniform functions as selection functions to select the individuals that contribute to the population at the next generation, and the Heuristic, the Scattered, and the Single point rules as crossover functions to combine two individuals to form the next generation for populations of 100, 500, and 2000 individuals. Table 1 summarizes the key results to compare the genetic algorithms performances. It can be concluded that the considered selection and crossover functions do not significantly change the results. CPU time depends on the population size, leading to equivalent results. Moreover, the results in Table 1 are in agreement with the trajectories defined by the different missions launched in the year 2020. The Mars 2020 mission (EEUU) [34] was launched on 30 July 2020, and its rover, the Perseverance rover, landed on 18 February 2021, with TOF = 203 days; The Tianwen-1 mission (China) [35] was launched on 23 July 2020 and arrived in Mars on 10 February 2021, with a TOF of 202 days. The Emirates Mars Mission (UAE) [36] was launched on 19 July 2020, and arrived at the orbit around Mars on 9 February 2021, with a TOF of 205 days. In order to analyze how a lower TOF impacts V∞E and V∞M , we compare the helio- centric Earth–Mars optimal transfer orbit with the resulting orbit when a TOF 31 days less than the optimal minimum-energy trajectory is imposed. We choose a population of 500 in- dividuals, the Stochastic Uniform function as the selection function, and the Heuristic as the crossover function. The resulting heliocentric elliptical orbits are illustrated in Figure 4, and their orbital elements are listed in Table 2. For the optimal minimum-energy orbit with a TOF of 197 days, departure on 20 July 2020 and arrival on 1 February 2021, we obtain values of V∞E = 3.6361 km/s and V∞M = 2.7682 km/s. When reducing the TOF to 166 days, de- parture and arrival dates change to 22 July 2020 and 4 January 2021, respectively, resulting in a more eccentric orbit, with V∞M = 3.5888 km/s significantly increased. Math. Comput. Appl. 2023, 28, 105 5 of 16 Table 1. Simulation scenarios in a launch window in 2020 to compare the genetic algorithm per- formances minimizing the cost function, C = C3 + V∞M , of the required C3 energy and the arrival velocity, V∞M . Pop. Crossover Selection CPU Departure TOF C3 V∞M C Time (s) Time (Days) (km2/s2) (km/s) Departure Date 20 July 2020 100 Heuristic remainder 10.86 01:19:08 196.9397 13.2216 2.7681 15.9897 stoch. unif. 10.36 01:05:27 196.9253 13.2212 2.7685 15.9897 Scattered remainder 10.30 01:04:01 196.9501 13.2213 2.7684 15.9897 stoch. unif. 10.49 01:01:17 196.9281 13.2212 2.7685 15.9897 Single pt remainder 10.14 01:06:02 196.9657 13.2218 2.7679 15.9897 stoch. unif. 10.35 01:14:40 196.9566 13.2217 2.7680 15.9897 500 Heuristic remainder 48.30 01:06:27 196.9288 13.2212 2.7685 15.9897 stoch. unif. 47.86 01:13:05 196.9420 13.2215 2.7682 15.9897 Scattered remainder 48.31 01:12:22 196.9420 13.2215 2.7682 15.9897 stoch. unif. 47.60 01:21:17 196.9537 13.2218 2.7679 15.9897 Single pt remainder 47.32 01:24:10 196.9407 13.2217 2.7681 15.9898 stoch. unif. 47.19 01:16:41 196.9450 13.2216 2.7681 15.9897 2000 Heuristic remainder 188.06 01:04:09 196.9393 13.2214 2.7683 15.9897 stoch. unif. 186.57 01:07:45 196.9332 13.2213 2.7684 15.9897 Scattered remainder 222.98 01:16:41 196.9450 13.2216 2.7681 15.9897 stoch. unif. 229.03 01:10:38 196.9361 13.2214 2.7683 15.9897 Single pt remainder 238.38 01:08:54 196.9343 13.2214 2.7683 15.9897 stoch. unif. 237.49 01:26:54 196.9480 13.2218 2.7679 15.9897 Table 2. Comparison of Earth–Mars optimal transfer heliocentric orbital parameters with respect to the mean ecliptic and equinox of J2000 for (a) the optimal launch and arrival dates in the window from 1 July 2019 to 1 November 2021 and (b) the launch and arrival dates in the window from 1 July 2019 to 1 November 2021 with the constraint of having 31 days of flight less than the optimal. Parameter (a) (b) Departure Date, t1 20 July 2020 01:13:05 22 July 2020 15:20:14 Arrival Date, t2 1 February 2021 4 January 2021 Arrival time 23:49:34 h 15:21:40 h Semimajor axis, asS (km) 198,312,598.97 202,972,264.04 Eccentricity, esS (unitless) 0.23346 0.25130 Inclination, isS (deg) 1.73626 0.72898 Ascending node long., ΩsS (deg) 297.4399 299.8067 Arg. of the perihelium, ωsS (deg) 359.6131 358.1402 True anomaly, vsS (deg) 0.4652 2.0420 V∞E (km/s) 3.6361 3.7803 V∞M (km/s) 2.7682 3.5888 Time of flight, TOF (days) 196.9420 166.0010 Total cost, C 15.9897 17.8795 For the analysis henceforth, we focus on the optimal orbit for a transference in the launch window of 2026. We consider the window with the earliest launch date on 1 March 2026 and latest arrival date on 1 November 2027. We consider the same parameters for the genetic algorithms as considered in the case of the launch window of 2020: a population of 500 individuals, the Stochastic Uniform function as the selection function, and the Heuristic as crossover function. The obtained trajectory has its launch date on 31 October 2026 and arrival date on 31 August 2027. Table 3 presents the heliocentric elliptic orbital elements for Math. Comput. Appl. 2023, 28, 105 6 of 16 this transference. The results, particularly the higher TOF value with respect to 2020, are in good agreement with the mission analysis referenced for year 2026 in [28]. Table 3. Elliptic heliocentric orbital elements for the optimal launch and arrival dates in the window from 1 March 2026 to 1 November 2027. Parameter Value Departure Date, t1 31 October 2026 05:42:13 h Arrival Date, t2 31 August 2027 16:47:12 h Semimajor axis, asS (km) 189,961,652.992134 Eccentricity, esS (unitless) 0.218496 Inclination, isS (deg) 0.8695 Periapsis argument, ωsS (deg) 4.1768 Right ascension node longitude, ΩsS (deg) 37.5815 True anomaly, vsS (deg) 197.9132 V∞E (km/s) 3.0311 V∞M (km/s) 2.5913 Time of flight, TOF (days) 304.4618 −0.5 0 0.5 1 1.5 −0.5 0 0.5 1 1.5 −0.035 0.035 Y axis (a u) (au) Z a x is (a u ) −0.5 0 0.5 1 1.5 −0.5 0 0.5 1 1.5 −0.035 0.035 Z a x is (a u ) 22−Jul−2020 20−Jul−2020 4−Jan−2021 X axis 1−Feb−2021 (au) X axis Y axis (a u) Figure 4. Minimum-energy optimal trajectory for the Earth–Mars transfer for departure date on 1 July 2019 and arrival date on 1 November 2021 (top) in comparison with a trajectory with 31 fewer days of flight (bottom). 3. Determination of Earth–Mars Trajectories with Hyperbolic Orbital Objective Values Once the launch and arrival dates for the optimal minimum-energy solution have been determined, we move on to deal with the determination of the orbit for the entry point at the Mars SOI. To this end, we implement a procedure for matching the Earth–Mars elliptic transfer orbit and the Mars arrival hyperbolic orbit, fixing the periapsis distance, rpsM , the arrival hyperbolic inclination, isM , and the radius for the SOI. We first present the determination of the arrival hyperbolic orbit and subsequent computation of the position and velocity at the matching point, which will be iterated with the Earth–Mars elliptic transfer afterwards. The goal is to analyze the changes in the B-plane [26] due to the variations in the arrival asymptote direction. The iterative procedure used to compute the entry point at the SOI, in which both the heliocentric elliptic transfer orbit for the obtained TOF and the areocentric hyperbolic arrival orbit match, is formulated as follows: r(i+1) sM = g ( f ( r(i)sM )) . (7) The function f provides the areocentric ecliptic transfer velocity, V∞M , at the Mars SOI entry point, rsM : f(r(i)sM) = V(i) ∞M , (8) obtaining first VsS2 , by solving Lambert’s problem (1) with a modified condition (3), rsS(t2) = rM + rsM , (9) Math. Comput. Appl. 2023, 28, 105 7 of 16 and then using Equation (5). The function g in (7) provides the position vector, rsM , in the areocentric ecliptic reference frame: g(V(i) ∞M) = r(i+1) sM . (10) This function gathers a set of expressions to obtain rsM at the SOI for the objective values of isM and rpsM , based on [22,24,25], fixing the radius of the SOI instead of the time that the spacecraft is inside the SOI, as in [22]. To this end, we first determine the areoequatorial coordinates (α, δ) of the arrival asymptote, given by−VsM (obtained by transforming V∞M to the areocentric areoequatorial reference frame), to obtain the parameter σ as σ = arcsin tan δ tan isM , (11) defining the minimum inclination as the value of ‖δ‖. For a direct orbit, the right ascension of the ascending node, ΩsM , can be computed in two different ways: ΩsM = α− σ (VsMz > 0), (12) ΩsM = α + σ + π (VsMz < 0). (13) The value of VsM determines the semimajor axis, asM . Then, the eccentricity of the hyperbola, esM , for a periapsis radius, rpsM , is fixed as esM = 1 + rpsM asM . (14) The true anomaly, vsM , of the spacecraft at the SOI is determined using the standard procedure (i.e., [37,38]). According to Figure 5, the unitary vectors defining the local reference system for the arrival asymptote, {uT , uB, uH}, are obtained as uT = −VsM ||VsM || , (15) uH =  sin isM sin ΩsM − sin isM cos ΩsM cos isM , (16) uB = uH × uT . (17) The plane containing uB, the planet center, and perpendicular to the arrival orbit is known as the B-plane, which is considered a fundamental tool in analyzing planetary arrivals [26]. With this iterative procedure, we update the B-plane according to the objective values imposed (see Figure 6). The vectors (15)–(17) are then transformed to a local reference system in the periapsis,{ uRp , uTp , uH } , according to uRp = − sin η uT + cos η uN , (18) uTp = − cos η uT − sin η uN , (19) where η = arcsin 1 esM . (20) Math. Comput. Appl. 2023, 28, 105 8 of 16 Then, the areocentric areoequatorial position RsM = (XsM , YsM , ZsM) vector compo- nents at the SOI are determined in terms of the base vectors { uRp , uTp , uH } , in a classical way, as RsM = rsM ( cos vsM uRp + sin vsM uTp ) . (21) Finally, the position vector RsM is transformed to obtain rsM in the areocentric ecliptic reference frame, as a result of the function g in (10). Figure 5. Local reference system vectors {~uT ,~uB} and { ~uRp ,~uTp } in the arrival hyperbolic orbit plane defined by the normal vector ~uH . Figure 6. Changes in the B-plane orientation due to the arrival asymptotic velocity variations. Math. Comput. Appl. 2023, 28, 105 9 of 16 We consider that the iterative method defined by (7) converges, for a given tolerance τ, when the following condition is achieved:∥∥∥r(i+1) sM − r(i)sM ∥∥∥ < τ. (22) For comparison with Table 3, in Table 4, we present the results for the elliptic transfer or- bit in the heliocentric ecliptic reference frame for the minimum-energy TOF (304.4618 days), the minimum arrival hyperbolic inclination (isMmin = 16◦.1167), rpsM = 20,428 km, and a SOI of rIM = 577,239 km. We note that V∞M changes by 0.015 km/s approximately. Also, the resulting areocentric hyperbolic arrival trajectory orbital elements in the areoequatorial reference frame are shown in Table 4. Table 4. Elliptic heliocentric orbital elements for (a) the Earth–Mars transfer orbit at the Mars sphere of influence entry point and final hyperbolic areocentric orbital elements and (b) launch date on 31 October 2026 and arrival date on 31 August 2027, for the objective values isM = isMmin = 16.◦1167 and rpsM = 20,428 km. Orbit Parameter Units Value (a) Departure date, t1 (UT) 31 October 2026 05:42:13 h Arrival date, t2 (UT) 31 August 2027 16:47:12 h Semimajor axis, asS (km) 189,905,238.422086 Eccentricity, esS (unitless) 0.218286 Inclination, isS (deg) 0.9311 Periapsis argument, ωsS (deg) 4.3094 Right ascension node longitude, ΩsS (deg) 37.5723 True anomaly, vsS (deg) 197.9302 V∞E (km/s) 3.0333 V∞M (km/s) 2.5763 Time of flight, TOF (days) 304.4618 (b) Semimajor axis, asM (km) 6600.229103 Eccentricity, esM (unitless) 4.0950441 Inclination, isM (deg) 16.1167 Periapsis argument, ωsM (deg) 194.1344 Right ascension node longitude, ΩsM (deg) 162.7116 True anomaly at t2, vsM (deg) 258.4533 Periapsis radius, rpsM (km) 20,467.9232 Arrival at periapsis date, tp (UT) 31 August 2027 20:03:55 4. Mars Arrival Maneuvers Evaluation for an Areostationary Mission We now conduct a preliminary evaluation of the total impulsive maneuvers ∆VM needed to capture the probe and to place it in an areostationary orbit: ∆VM = ∆Vc + ∆VP + ∆VA + ∆Vi, (23) where ∆Vc is the capture maneuver to avoid the probe leaving the SOI on a flyby trajec- tory; ∆VP and ∆VA are the two Hohmann transfer maneuvers, at the perigee and apogee, respectively, designed to obtain a transfer orbit from the capture orbit to the target are- ostationary orbit; and ∆Vi is the inclination correction maneuver to reach the final zero desired inclination. If we consider a ∆Vc maneuver at the periapsis of the hyperbolic orbit in order to obtain a circular capture orbit, its magnitude would be ∆Vc = √ V2 ∞M + 2 µM rpsM − √ µM rpsM . (24) Math. Comput. Appl. 2023, 28, 105 10 of 16 The hyperbolic periapsis distance, rpsM , is calculated as follows: rpsM = asM(esM − 1), (25) with asM as the hyperbolic semimajor axis and esM as the hyperbolic eccentricity. The two Hohmann transfer maneuvers required to achieve an orbit at the areosta- tionary semimajor axis, rA = 20,428 km, in the case that rpsM < rA, are calculated as ∆VP = √ µM rpsM (√ 2rA rpsM+rA − 1 ) , (26) ∆VA = √ µM rA ( 1− √ 2rpsM rpsM+rA ) . (27) Finally, the inclination maneuver is performed at the stationary distance in order to minimize its magnitude. For a maneuver at the node, this is obtained as ∆Vi = √ 2 µM rA − 2 µM rA cos isM , (28) where isM is the hyperbolic inclination. In the case that rpsM > rA, the inclination maneuver at the node would be the first to be performed at the rpsM distance: ∆Vi = √ 2 µM rpsM − 2 µM rpsM cos isM . (29) Then, the two Hohmann maneuvers would be performed according to the following equations: ∆VP = √ µM rpsM ( 1− √ 2rA rpsM+rA ) , (30) ∆VA = √ µM rA (√ 2rpsM rpsM+rA − 1 ) . (31) As can be seen from Equations (24) to (31), the periapsis distance, rpsM , and the orbital inclination of the arrival trajectory with respect to Mars, isM , are the key design parameters in order to minimize the required impulses. Corrections to change the approach asymptotic plane are not considered, and neither are the combinations of the capture maneuver with the first Hohmann maneuver or the inclination maneuver with the second Hohmann maneuver. Several numerical simulations are carried out to analyze the magnitude of these impulses due to the variations in the B-plane when different imposed values for TOF, isM , and rpsM are considered. We fix the launch window for 2026 according to the parameters given in Table 3 for the launch and arrival dates connecting Earth and Mars that minimize V∞M . In the first set of simulations, we consider different imposed objective values for the inclination, isM , and the periapsis radius, rpsM , to determine the entry point at the SOI, using the iterative procedure defined in (7) to evaluate the total impulse ∆VM in Equation (23). We consider different values for isM , ranging in the interval [ isMmin, 90◦ ) , and for rpsM , in the interval [15,000 km, 25,000 km], with the arrival velocity, V∞M , recalculated for the different entry points. Figure 7a,b represent the results obtained for V∞M with different values of rpsM and isM , respectively. For the objective conditions of isM = 16◦.1167, which corresponds to the isMmin for rpsM = rA = 20,428 km, a value of V∞M = 2.5763 km/s is obtained. As can be observed in Figure 7a for isMmin, V∞M decreases linearly when increasing rpsM from 2.5768 km/s to 2.5759 km/s. Figure 7b shows variations of about 2× 10−3 km/s in V∞M , for the areostationary distance, as the objective inclination increases from isMmin to 90◦. Table 5 summarizes the values for the C3, V∞M , and total cost C of Equation (6) for the minimum and maximum values of rpsM and isM considered. Math. Comput. Appl. 2023, 28, 105 11 of 16 15 00 0 16 00 0 17 00 0 18 00 0 19 00 0 20 00 0 21 00 0 22 00 0 23 00 0 24 00 0 25 00 0 r ps M (km) 2.576 2.5765 2.577 2.5775 2.578 V M ( k m /s ) (a) 10 20 30 40 50 60 70 80 90 i s M (deg) 2.5762 2.5764 2.5766 2.5768 2.577 2.5772 2.5774 2.5776 2.5778 2.578 2.5782 V M ( k m /s ) (b) Figure 7. V∞M values for (a) the isMmin for different periapsis radii, rpsM , ranging in the interval [15,000 km, 25,000 km] and (b) the objective values of rpsM = rA = 20,428 km for different inclination values, isM , ranging in the interval [ isMmin = 16.◦1167, 90◦ ] . Figure 8a,b represent the values of ∆Vc, ∆Vi, ∆VA, ∆VP, and the total ∆VM for dif- ferent values of rpsM and isM to achieve the circular areostationary orbit at rA with zero areoequatorial inclination. In the case that rpsM < rA, the maneuvers are computed using (24) and (26)–(28). In the other case, the maneuvers are computed by means of (24) and (29)–(31). Math. Comput. Appl. 2023, 28, 105 12 of 16 Table 5. Interplanetary transfer energy variations for the considered combinations of isM and the maximum and minimum rpsM and of rpsM = rA with the maximum and minimum isM . isM (deg) rpsM (km) C3 (km2/s2) V∞M (km/s) C isMmin = 16.1158 15,000 9.2014 2.5768 11.7782 isMmin = 16.1175 25,000 9.2008 2.5759 11.7767 isMmin = 16.1167 rA = 20, 428 9.2011 2.5763 11.7774 90 rA = 20, 428 9.2123 2.5781 11.7904 As can be seen in Figure 8a, with the above-described strategy, the values for the capture maneuver ∆Vc change, for isMmin, from 1.8246 km/s for rpsM = 15,000 km to 1.8631 km/s for rpsM = 25,000 km, depending on the different objective periapsis radius considered. Obviously, the Hohmann maneuver vanishes for a capture at rA, and its value varies from 0.2404 km/s for rpsM = 15,000 km to 0.1387 km/s for rpsM = 25,000. The in- clination maneuver, ∆Vi, to achieve zero inclination, from the isMmin for each rpsM , has a slight variation from 0.4062 km/s to 0.3673 km/s. The minimum for the total amount, ∆VM, has a value of 2.2493 km/s, which corresponds to rpsM = rA. This magnitude rises to 2.4712 km/s for rpsM = 15,000 km and to 2.3691 km/s for rpsM = 25,000 km. Figure 8b shows that, for rpsM = rA, with ∆VH = 0 km/s, ∆VM increases when different values for the objective hyperbolic inclination are considered. From the same minimum of ∆VM = 2.2493 km/s for isMmin = 16.◦1167, the total quantity of impulses rises to 3.8922 km/s, as illustrated in Figure 8a. In a second set of simulations, we consider a different strategy to evaluate the Oberth effect that, due to the potential energy for a capture near the Mars surface, allows one to reduce the required ∆VM. To this end, we consider a capture into an elliptical trajectory with a low periapsis radius and an apoapsis at the stationary orbit. Then, a circularization process is performed at the apoapsis. For a capture at a periapsis of rpsM = rM + 300 km into an ellipse with an apoapsis of rA, the value of the impulse required to obtain a circular orbit at the areostationary distance is decreased to 1.6764 km/s, with ∆Vc = 1.0294 km/s and ∆VA = 0.6470 km/s. The total ∆VM = 2.0834 km/s reduces 0.1659 km/s with respect to the previous strategy. Finally, we analyze the influence of the TOF in ∆VM for a possible launch date delay. We consider up to two weeks delay from 31 October 2026 for the launch date and up to two months after 31 August 2027 for the arrival date for the minimum isM . In the first case, a capture into a circular orbit of rpsM = rA is considered, as shown in Figure 9a. In the second case, a capture into an ellipse with rpsM = rM + 300 and an apoapsis of rA is considered, as shown in Figure 9b. It can be noticed, by comparing Figure 9a,b, that the second case is less sensitive to TOF variations. For case (a), a delay of 14 days at departure and 60 at arrival, increases the value of 2.25 km/s by 0.75 km/s, corresponding to the minimum-energy transfer TOF of 304 days, instead of an increase of 0.52 km/s from the value 2.08 km/s in case (b). In both cases, a launch delay of 14 days maintaining a TOF of 304 days would increase ∆VM by 0.25 km/s. Math. Comput. Appl. 2023, 28, 105 13 of 16 15 00 0 16 00 0 17 00 0 18 00 0 19 00 0 20 00 0 21 00 0 22 00 0 23 00 0 24 00 0 25 00 0 r ps M (km) 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 V ( k m /s ) V c V i V H V M (a) 10 20 30 40 50 60 70 80 90 i M (deg) 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 V ( k m /s ) V c V i V H V M (b) Figure 8. ∆Vc, ∆Vi, ∆VA, ∆VP, and the total ∆VM values for (a) the isMmin for different peri- apsis radii, rpsM , ranging in the interval [15,000 km, 25,000 km] and (b) the objective values of rpsM = rA = 20,428 km for different inclinations, isM , ranging in the interval [ isMmin = 16.◦1167, 90◦ ] . Math. Comput. Appl. 2023, 28, 105 14 of 16 2.3 2.35 2.35 2.4 2.4 2 .4 5 2.45 2.45 2.5 2.5 2.5 2.55 2.55 2.6 2.6 2.6 2.65 2.65 2.7 2.7 2.7 2.75 2.75 2.8 2.82.8 2.85 2.85 2.9 2.92.9 2.95 2.95 3 3 0 2 4 6 8 10 12 14 Days from Earth departure date 31-Oct-2026 0 10 20 30 40 50 60 D a y s f ro m M a rs a rr iv a l d a te 3 1 -A u g -2 0 2 7 (a) 2.1 2.15 2.2 2.2 2 .2 5 2.25 2.25 2 .3 2.3 2.3 2.35 2.35 2.35 2.4 2.4 2.4 2.45 2.45 2.5 2.5 2.5 2.55 2.55 2.6 2.6 0 2 4 6 8 10 12 14 Days from Earth departure date 31-Oct-2026 0 10 20 30 40 50 60 D a y s f ro m M a rs a rr iv a l d a te 3 1 -A u g -2 0 2 7 (b) Figure 9. ∆V values for a combination of dates ranging from 31 October 2026 to 13 November 2026 for the launch date and from 31 August 2027 to 30 October 2027 for the arrival date, for the following objective conditions: (a) minimum isM and rpsM = rA = 20,428 km; (b) minimum isM and a capture into an ellipse with a periapsis of rpsM = rM + 300 km and an apoapsis of rA. 5. Conclusions In this paper, we conducted a preliminary analysis of the impulsive maneuvers cost to transfer a spacecraft from Earth to Mars and to position it in an areostationary orbit. We first obtained the minimum-energy interplanetary transfer from Earth to Mars by applying genetic algorithms to select launch and arrival dates. Several simulations were carried out to analyze the performance of the genetic algorithms. Results show that differences in the Math. Comput. Appl. 2023, 28, 105 15 of 16 final energy cost and the time of flight are negligible, and the only significant change is the CPU time needed to converge, which is dependent on the population size. With the optimized launch and arrival dates, an iterative procedure was used to match the interplanetary trajectory obtained with the genetic algorithms, with an entry hyperbola defined by imposing objective conditions for the inclination and the periapsis radius of the orbit. Two different strategies were computed to evaluate the cost of the mission, ∆VM: The first includes a capture maneuver to a circular orbit at different periapsis radii, as well as Hohmann maneuvers and an inclination maneuver; The second includes a capture maneuver to an elliptic orbit with a low periapsis and an apoapsis at the areostationary orbit. Simulations with different imposed conditions on the entry hyperbola were con- ducted depending on the two key parameters: the hyperbolic inclination, isM , and the periapsis radius, rpsM . For a circular capture at the stationary radius, results show that for a 2026 mission with a TOF of 304 days for the minimum-energy Earth–Mars transfer trajectory, the values achieved are ∆Vc = 1.84, ∆VH = 0 and ∆Vi = 0.41, being the to- tal impulse ∆VM = 2.25 km/s for the minimum possible inclination isM = 16.◦1167 and rpsM = 20,428 km corresponding to an areostationary radius. For an elliptical capture with a periapsis of rpsM = rM + 300 km and an apoapsis of 20,428, the values achieved are ∆Vc = 1.03, ∆VA = 0.65, and ∆Vi = 0.41, the total impulse amounting to ∆VM = 2.08 km/s. A launch delay of two weeks would increase this minimum value of ∆VM by 0.25 km/s in both cases. The capture into an ellipse with a delay of 14 days at departure and 60 at arrival is less sensitive to TOF variations, increasing the total ∆VM by 0.52 km/s, instead of 0.75 for the direct circular capture. Author Contributions: All authors contributed equally to this work. All authors have read and agreed to the published version of the manuscript. Funding: This research received no external funding. Data Availability Statement: No new data were created or analyzed in this study. Data sharing is not applicable to this article. Conflicts of Interest: The authors declare no conflict of interest. Abbreviations The following abbreviations are used in this manuscript: CPU Central Processing Unit SOI Sphere Of Influence TOF Time Of Flight References 1. Edwards, C.; Arnold, B.; DePaula, R.; Kazz, G.; Lee, C.; Noreen, G. Relay communications strategies for Mars exploration through 2020. Acta Astronaut. 2006, 59, 310–318. [CrossRef] 2. Edwards, C.; DePaula, R. Key telecommunications technologies for increasing data return for future Mars exploration. Acta Astronaut. 2007, 61, 131–138. [CrossRef] 3. Jentsch, C.; Rathke, A.; Wallner, O. Interplanetary communication: A review of future missions. In Proceedings of the 2009 International Workshop on Satellite and Space Communications, Siena, Italy, 9–11 September 2009; pp. 291–294. [CrossRef] 4. Podnar, G.; Dolan, J.; Elfes, A. Telesupervised robotic systems and the human exploration of Mars. J. Cosmol. 2010, 12, 4058–4067. 5. Romero, P.; Pablos, B.; Barderas, G. Analysis of orbit determination from Earth-based tracking for relay satellites in a perturbed areostationary orbit. Acta Astronaut. 2017, 136, 434–442. [CrossRef] 6. Montabone, L.; Heavens, N.; Babuscia, A.; Barba, N.; Battalio, J.; Bertrand, T.; Edwards, C.; Guzewich, S.; Kahre, M.; Kass, D.; et al. Observing Mars from Areostationary orbit: Benefits and applications. In Proceedings of the Mars Exploration Program Analysis Group (MEPAG) #38, Virtual, 15–17 April 2020. 7. Montabone, L.; Heavens, N.; Alvarellos, J.L.; Lillis, R.; Aye, M.; Liuzzi, G.; Babuscia, A.; Mischna, M.A.; Barba, N.; Newman, C.E.; et al. Observing Mars from Areostationary Orbit: Benefits and Applications. Available online: https://doi.org/10.13140/RG.2.2. 21498.72643 (accessed on 20 September 2023). 8. Romero, P.; Gambi, J. Optimal control in the east/west station-keeping manoeuvres for geostationary satellites. Aerosp. Sci. Technol. 2004, 8, 729–734. [CrossRef] http://doi.org/10.1016/j.actaastro.2006.02.038 http://dx.doi.org/10.1016/j.actaastro.2007.01.016 http://dx.doi.org/10.1109/IWSSC.2009.5286356 http://dx.doi.org/10.1016/j.actaastro.2017.04.002 https://doi.org/10.13140/RG.2.2.21498.72643 https://doi.org/10.13140/RG.2.2.21498.72643 http://dx.doi.org/10.1016/j.ast.2004.08.004 Math. Comput. Appl. 2023, 28, 105 16 of 16 9. Romero, P.; Gambi, J.; Patiño, E. Stationkeeping manoeuvres for geostationary satellites using feedback control techniques. Aerosp. Sci. Technol. 2007, 11, 229–237. [CrossRef] 10. Prado, A.; Broucke, R. Transfer orbits in the Earth-Moon system using a regularized model. J. Guid. Control. Dyn. 1996, 19, 929–933. [CrossRef] 11. Broucke, R.; Prado, A. Jupiter swing-by trajectories passing near the Earth. Adv. Astronaut. Sci. 1993, 82, 1159–1176. 12. Kluever, C. Efficient Computation of Optimal Interplanetary Trajectories Using Solar Electric Propulsion. J. Guid. Control. Dyn. 2014, 38, 5. [CrossRef] 13. Li, X.; Qiao, D.; Chen, H. Interplanetary transfer optimization using cost function with variable coefficients. Astrodynamics 2019, 3, 173–188. [CrossRef] 14. Chen, L.; Li, J. Optimization of Earth-Mars transfer trajectories with launch constraints. Astrophys. Space Sci. 2022, 367, 12. [CrossRef] 15. Woolley, R.; Whetsel, C. On the nature of Earth-Mars Porkchop plots. Adv. Astronaut. Sci. 2013, 148, 413–426. 16. Sanchez-Garcia, M.M.; Barderas, G.; Romero, P. Analysis of the optimization for an Earth to Mars areostationary mission. In Proceedings of the European Planetary Science Congress, Virtual, 21 September–9 October 2020; Volume 14, p. EPSC2020-134. [CrossRef] 17. Gooding, R.H. A procedure for the solution of Lambert’s orbital boundary-value problem. Celest. Mech. Dyn. Astron. 1990, 48, 145–165. [CrossRef] 18. Conte, D. Survey of Earth-Mars trajectories using Lambert’s Problem and Applications. Ph.D. Thesis, The Pennsylvania State University, State College, PA, USA, 2014. 19. Matlab. Genetic Algorithms Options, 2022. Available online: https://es.mathworks.com/help/gads/genetic-algorithm-options. html (accessed on 20 September 2023). 20. Clarke, V. C., Jr.; Bollman, W.E.; Feitis, P.H.; Roth, R.Y. Design Parameters for Ballistic Interplanetary Trajectories, Part II: One-Way Transfers to Mercury and Jupiter; Technical Report 32–77. JPL; US Gov.: Washington, DC, USA, 1966. 21. Cornelisse, J.W. Trajectory analysis for interplanetary missions. ESA J. 1978, 2, 131–144. 22. Parvathi, S.P.; Ramanan, R.V. Direct Transfer Trajectory Design Options for Interplanetary Orbiter Missions using an Iterative Patched Conic Method. Adv. Space Res. 2016, 59, 1763–1774. [CrossRef] 23. Parvathi, S.P.; Ramanan, R.V. Direct interplanetary trajectory design with a precise V-infinity targeting technique. In Proceedings of the 2017 First International Conference on Recent Advances in Aerospace Engineering (ICRAAE), Coimbatore, India, 3–4 March 2017; pp. 1–6. [CrossRef] 24. Iwabuchi, M.; Satoh, S.; Yamada, K. Smooth and continuous interplanetary trajectory design of spacecraft using iterative patched-conic method. Acta Astronaut. 2021, 185, 58–69. [CrossRef] 25. Bond, V.R. Matched-Conic Solutions to Round-Trip Interplanetary Trajectory Problems That Insure State-Vector Continuity at All Boundaries; Technical Note D-4942; NASA: Washington, DC, USA, 1969. 26. Farnocchia, D.; Eggl, S.; Chodas, P.; Giorgini, J.; Chesley, S. Planetary encounter analysis on the B-plane: A comprehensive formulation. Celest. Mech. Dyn. Astron. 2019, 131, 36. [CrossRef] 27. George, L.E.; Kos, L.D. Interplanetary Mission Design Handbook: Earth-to-Mars Mission Opportunities and Mars-to-Earth Return Opportunities, 2009–2024; Technical Report TM-1998-208533; NASA: Washington, DC, USA, 1998. 28. Burke, L.M.; Falck, R.D.; McGuire, M.L. Interplanetary Mission Design Handbook: Earth-to-Mars Mission Opportunities 2026 to 2045; Technical Report TM-2010-216764; NASA: Washington, DC, USA, 2010. 29. Lancaster, E.R.; Blanchard, R.C. A Unified Form of Lambert’s Theorem; Technical Note D-5368; NASA: Washington, DC, USA, 1969. 30. Oldenhuis, R. Robust Solver for Lambert’s Orbital-Boundary Value Problem, 2017. Available online: http://es.mathworks. com/matlabcentral/fileexchange/26348-robust-solver-for-lambert-s-orbital-boundary-value-problem?? (accessed on 20 Septem- ber 2023). 31. Conte, D.; Di Carlo, M.; Ho, K.; Spencer, D.; Vasile, M. Earth-Mars transfers through Moon Distant Retrograde Orbits. Acta Astronaut. 2017, 143, 372–379. [CrossRef] 32. Conte, D.; Spencer, D. Mission Analysis for Earth to Mars-Phobos Distant Retrograde Orbits. Acta Astronaut. 2018, 151, 761–771. [CrossRef] 33. Trajectory Optimization Tool, v2.1.1. 2011. Available online: http://www.orbithangar.com/searchid.php?ID=5418 (accessed on 20 September 2023). 34. Farley, K.A.; Williford, K.H.; Stack, K.M.; Bhartia, R.; Chen, A.; de la Torre, M.; Hand, K.P.; Goreva, Y.; Herd, C.D.K.; Hueso, R.; et al. Mars 2020 Mission Overview. Space Sci. Rev. 2020, 216, 142. [CrossRef] 35. Jiang, X.; Yang, B.; Li, S. Overview of China’s 2020 Mars mission design and navigation. Astrodynamics 2018, 2, 1–11. [CrossRef] 36. Sharaf, O.; Amiri, S.; AlDhafri, S.; Withnell, P.; Brain, D. Sending Hope to Mars. Nat. Astron. 2020, 4, 722. [CrossRef] 37. Capderou, M. Satellites. Orbits and Missions; Springer: Paris, France, 2005. 38. Montenbruck, O.; Gill, E. Satellite Orbits: Models, Methods, and Applications; Springer: Berlin/Heidelberg, Germany, 2000. Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content. http://dx.doi.org/10.1016/j.ast.2006.08.003 http://dx.doi.org/10.2514/3.21720 http://dx.doi.org/10.2514/1.G000144 http://dx.doi.org/10.1007/s42064-018-0043-8 http://dx.doi.org/10.1007/s10509-021-04025-2 http://dx.doi.org/10.5194/epsc2020-134 http://dx.doi.org/10.1007/BF00049511 https://es.mathworks.com/help/gads/genetic-algorithm-options.html https://es.mathworks.com/help/gads/genetic-algorithm-options.html http://dx.doi.org/10.1016/j.asr.2017.01.023 http://dx.doi.org/10.1109/ICRAAE.2017.8297227 http://dx.doi.org/10.1016/j.actaastro.2021.04.021 http://dx.doi.org/10.1007/s10569-019-9914-4 http://es.mathworks.com/matlabcentral/fileexchange/26348-robust-solver-for-lambert-s-orbital-boundary-value-problem?? http://es.mathworks.com/matlabcentral/fileexchange/26348-robust-solver-for-lambert-s-orbital-boundary-value-problem?? http://dx.doi.org/10.1016/j.actaastro.2017.12.007 http://dx.doi.org/10.1016/j.actaastro.2018.06.049 http://www.orbithangar.com/searchid.php?ID=5418 http://dx.doi.org/10.1007/s11214-020-00762-y http://dx.doi.org/10.1007/s42064-017-0011-8 http://dx.doi.org/10.1038/s41550-020-1151-y Introduction Minimum-Energy Launch Window for Earth–Mars Transfer Trajectories Determination of Earth–Mars Trajectories with Hyperbolic Orbital Objective Values Mars Arrival Maneuvers Evaluation for an Areostationary Mission Conclusions References