Abstract

Translate

Developing spacecraft for efficient aerocapture missions demands managing extreme aerothermal environments, precise controls, and atmospheric uncertainties. Successful designs must integrate vehicle airframe considerations with trajectory planning, adhering to launcher dimension constraints and ensuring robustness against atmospheric and insertion uncertainties. To advance robust multi-objective optimization in this field, a new framework is presented, designed to rapidly analyze and optimize non-thrusting, fixed angle-of-attack aerocapture-capable spacecraft and their trajectories. The framework employs a three-degree-of-freedom atmospheric flight dynamics model incorporating planet-specific characteristics. Aerothermal effects are approximated using established Sutton–Graves, Tauber–Sutton, and Stefan–Boltzmann relations. The framework computes the resulting post-atmospheric pass orbit using an orbital element determination algorithm to estimate fuel requirements for orbital corrective maneuvers. A novel algorithm that consolidates multiple objective functions into a unified cost function is presented and demonstrated to achieve superior optima with computational efficiency compared to traditional multi-objective optimization approaches. Numerical examples demonstrate the methodology’s effectiveness and computational cost at optimizing terrestrial and Martian aerocapture maneuvers for minimum fuel, heat loads, peak heat transfers, and an overall optimal trajectory, including volumetric considerations.

Full text

Turn on search term navigation

Translate

1. Introduction

A critical step in planetary exploration is the transition of a spacecraft from a hyperbolic trajectory to a closed orbit around a target planet. Conventional orbit insertion relies on high-thrust chemical propulsion systems, which require substantial propellant reserves, thereby limiting the available payload mass. Since the early 2000s, the adoption of electric propulsion for orbital insertion has grown significantly. Due to its extreme efficiency, its use for deep space missions and attitude control has rapidly increased. However, despite its benefits, electric propulsion has rarely been employed for full orbital capture. Instead, many missions have utilized electric propulsion primarily for station-keeping and proximity operations rather than achieving a bound orbit. To date, only three missions—SMART-1 [1], Dawn [2,3], and Psyche [4]—rely primarily on electric propulsion for orbital insertion, leveraging continuous low-thrust trajectories over extended durations.

Aerocapture offers a compelling alternative to conventional propulsive orbit insertion, utilizing atmospheric drag to decelerate the spacecraft and eliminating the need for large propellant reserves. This approach has the potential to significantly enhance payload capacity. Since its initial introduction [5], the feasibility and implementation of aerocapture have been extensively studied in the literature [6,7,8,9]. A comparative study assessing the relative risks of propulsive capture, aerobraking, and aerocapture concluded that the reliability of these three strategies differed only marginally [6]. This was largely because nearly all subsystems required for aerocapture had already demonstrated flight heritage and operational reliability in previous missions.

Recent studies demonstrate that NASA is technologically equipped to implement aerocapture Titan, Mars, and Venus [8,9]. These reports state that a dedicated flight demonstration mission is unnecessary, as several aeroassisted techniques—each presenting greater challenges than aerocapture—have been successfully executed. These include the skip-entry maneuver employed by the Soviet Zond [10], the Chinese Chang’e 5-T1 [11], the Apollo 4 and 6 missions [12], and the hypersonic aeromaneuvering task performed by the MSL [13], which all exemplify such achievements. Nonetheless, space exploration mission planners remain reluctant to adopt aerocapture due to flight-day uncertainties, including variations in atmospheric density and adverse meteorological conditions.

Recent optimization efforts in this field [14,15,16,17,18,19,20,21] have centered on the pursuit of efficient trans-atmospheric trajectories through meticulous manipulation of drag and lift forces to achieve an optimal aerocapture or reentry maneuver. Recent years have witnessed a push towards enhancing aerocapture technology through drag modulation (DM) techniques, particularly with the exploration of non-rigid aeroshell concepts such as HIAD and ADEPT [22]. These methods allow control over the spacecraft’s ballistic coefficient ( $β$ ) by manipulating either the reference area or the drag coefficient ( $C_{D}$ ). The manipulation of these quantities can be performed either via continuous strategies or discrete events. However, certain technical challenges remain, particularly those associated with jettison or deployment events during hypersonic flight [23]. In contrast, lift modulation, which has been successfully demonstrated in missions such as Apollo, MSL, and the Space Shuttle, relies on sophisticated avionics algorithms and control systems to enable continuous bank angle modulation. Despite the complexities associated with this technique, it remains better understood and established than DM techniques. Optimal aerocapture problems have often overlooked the complex engineering challenges in effecting control over these parameters. Furthermore, the influence of aerodynamic control on the resulting captured orbit is limited, and its nature is highly dependent upon the spacecraft’s state at the atmospheric interface $(A I)$ . Consequently, optimizing the atmospheric trajectory yields only incremental benefits to the overall mission, primarily impacting atmospheric phenomena such as heat loads.

Several optimization strategies for aerocapture trajectories avoid the continuous modulation of control variables during the atmospheric pass, instead relying on low-thrust maneuvers to guide the spacecraft to its desired operational orbit after the atmospheric phase. Giordano and Topputo [24] demonstrated this approach through ballistic capture. Their technique exploited gravitational interactions with multiple celestial bodies at the same time, allowing the Sun’s gravity to accelerate the spacecraft and naturally raise its periapsis, eliminating the need for a periapsis raise burn. Additionally, Fornari and Pontani [25] explored the impact of uncertainties in initial conditions and atmospheric density on aerocapture mission success, as well as the use of nonlinear orbit correction maneuvers to improve robustness. Their statistical analysis showed a high likelihood of three unfavorable outcomes—escape, near-capture, and destructive entry—necessitating distinct corrective maneuvers to ensure a successful trajectory. In a separate study, Gochenaur et al. [26] examined the feasibility of using lift vectoring for orbital plane rotations. While their focus was on assessing the maneuver’s viability rather than its practical implementation within an optimization framework, their work laid the foundation for future applications, as plane change maneuvers are known to be highly fuel-intensive.

However, because these studies assume specific vehicle geometries and initial conditions at the atmospheric interface, they are not well suited for preliminary mission design and planning applications. Given the multidisciplinary nature of such missions, an optimization framework must account for both the initial conditions—linked to the injection orbit—and vehicle sizing—determined by mission requirements—to provide planners with an optimal and viable preliminary solution. The AMAT framework [27], built on previous work by Girija [28], is a Python-based tool designed to aid in the design and sizing of re-entry vehicles for aeroassisted trajectories, making it particularly useful for preliminary mission planning. However, its primary focus is on trade-off analysis rather than design optimization. An example of a more suitable algorithm for preliminary design is the framework developed by Armellin and Lavagna [29], which enables a detailed analysis of the interplay between vehicle shape, trajectory control, and thermal protection system design. However, this framework assumes that the orbital plane achieved as the vehicle exits the atmosphere matches the target orbit plane, thereby yielding the target entry flight path angle (EFPA). However, that study overlooks uncertainties in EFPA and atmospheric density, meaning the entire interplanetary injection trajectory is designed without accounting for these critical uncertainties, which play a crucial role in determining mission success [25]. The approach proposed by Hanninen et al. [30] shares these same limitations.

Additionally, current strategies often employ multi-objective optimization (MOO) [14,15,16,17,29] for optimal re-entry trajectory problems, presenting numerical and decision-making challenges. MOO involves navigating a multi-dimensional objective space, resulting in trade-offs between objectives. The exploration of vast potential solutions and the non-convex and irregular Pareto front further complicates optimal solution searches. Selecting a single solution from a multitude of trade-offs involves considering subjective preferences, priorities, and uncertainty, requiring informed judgment, which can be challenging for algorithms due to ambiguous decision-making criteria.

In light of these multifaceted challenges—including the difficulty of continuously controlling lift and drag, uncertainties in EFPA, scattering of atmospheric densities, and the limitations of frameworks that impose fixed geometry or insertion trajectories—we propose a novel framework that supports the rapid preliminary design and optimization of aerocapture missions without the restrictive assumptions on vehicle geometry and insertion trajectory. Developed in MATLAB 2023a, the Determination of Aerocapture Successful Trajectories and Robust Optimization (D-ASTRO) algorithm is an efficient and robust tool designed for rapid analysis and optimization of generic aerocapture missions. Its unique capability lies in its ability to optimize the insertion trajectory and sizing of fixed angle-of-attack trajectories for aerocapture-capable spacecraft, while accounting for critical uncertainties. By requiring only minimal aeroshell information and mission parameters, D-ASTRO employs a nonlinear, single-objective constrained optimization strategy to identify optimal values for the flight path angle ( $γ$ ) at the AI and the fixed ballistic coefficient ( $β$ ) of the aeroshell.

The paper is structured as follows: Section 2 introduces the concept of the aerocapture corridor. The optimization problem, along with novel objective functions and performance metrics, is discussed in Section 3. Section 4 focuses on the pipeline of the algorithm and the novel methodology used to rapidly compute the aerocapture corridor boundaries. Followed by the modeling considerations in Section 5. Finally, Section 6 presents the algorithm’s numerical stability and convergence properties applied to a typical Martian aerocapture mission test case, with the computational considerations and performance of D-ASTRO explored in Section 7.

2. Aerocapture Corridor

Without employing lift or drag modulation, the flight path angle $γ$ at the AI and $β$ primarily govern the vehicle dynamics during the atmospheric pass, ultimately determining the success of the aerocapture maneuver. The aerocapture corridor defines the range of $γ$ and $β$ values that result in a successful maneuver, that is, any atmospheric pass that results in a bounded orbit with an apoapsis greater than a specified altitude. The lower boundary is reached when excessive drag causes the spacecraft to descend to the surface, a condition which is readily determined by checking when the altitude along the trajectory reaches zero. In contrast, the upper boundary is defined by the condition where insufficient energy is dissipated, leaving the vehicle in a hyperbolic orbit. This is determined by the eccentricity of the orbit after the atmospheric pass, where $e < 1$ indicates a closed orbit. Notably, trajectories exceeding the upper boundary may still be considered optimal if they significantly decelerate the spacecraft, although they necessitate a subsequent propulsive burn to reduce orbital energy, resulting in an aerobraking maneuver.

Figure 1 compares the aerocapture corridors for ballistic and lifting re-entry profiles of a spacecraft approaching Mars with a hyperbolic excess velocity, $V_{e}$ , of 3.5 ${kms}^{- 1}$ . The maneuver’s success is represented by green hashing in Figure 1. These contours illustrate the intricate relationship among $γ$ , $β$ , and the lift-to-drag ratio necessary for successful maneuver. A lifting re-entry profile provides a more favorable corridor, reducing sensitivity to precise insertion attitudes. While D-ASTRO supports various re-entry profiles, this study focused on non-zero $L / D$ profiles for Martian aerocapture, utilizing a fixed value of 0.2. This choice aligns with the design characteristics of heritage Martian re-entry vehicles [24,31].

Uncertainties in atmospheric properties significantly influence the width of the aerocapture corridor. Variations in atmospheric density lead to changes in drag forces, causing the spacecraft to dissipate energy differently than expected. Near the aerocapture boundaries, even minor density fluctuations can result in maneuver failure, leading to either a hyperbolic orbit or surface collision. Two corridors are computed to establish a robust corridor: one with reduced atmospheric density and another with increased density. The “robust” aerocapture corridor represents the intersection of these two, resulting in a narrowed corridor. These corridors are illustrated in Figure 2, depicting terrestrial and Martian corridors for various inbound $V_{e}$ , where the green and yellow regions represent the robust and non-robust corridors.

3. Optimal Aerocapture Maneuver via Parameter Optimization

Multiple approaches can be taken to optimize an aerocapture trajectory. However, we focused on the optimization of design parameters $γ$ and $β$ which yielded an optimal aerocapture, or aeroassisted, trajectory. Determining these two parameters allows mission planners to identify the most suitable interplanetary orbit, sizing of the spacecraft’s aeroshell, and structural design of the heat shield to ensure an optimal trajectory.

3.1. Formulation of Optimization Problem

The optimal aerocapture corridor problem is a constrained nonlinear optimization problem of the form:

(1) $\begin{matrix} \underset{β, γ_{A I}}{minimize} & f {(θ; β, γ_{A I})}_{i} i = 1, \dots, l \\ subject to & γ_{A I} - g_{u} (β) + Δ γ \leq 0, \\ g_{l} (β) - γ_{A I} + Δ γ \leq 0, \\ β > 0 \end{matrix}$

where f represents the objective function, and

g_{u}

and

g_{l}

are functions of

β

that define the upper and lower boundaries of the aerocapture corridor. The vector

θ

comprises selected performance indices.

Insertion inaccuracies arising from navigation and modeling errors can lead to significant deviations from the intended atmospheric insertion coordinates. Deviating from the intended $γ_{A I}$ leads to different atmospheric trajectories, impacting the maneuver’s optimality and sometimes its success. To account for this, an expected insertion inaccuracy is incorporated as $Δ γ$ , which narrows the aerocapture corridor to ensure mission success is robust to insertion inaccuracies.

Insertion inaccuracies become particularly critical when the optimal trajectory approaches the lower boundary of the aerocapture corridor, as even minor deviations can result in excessively deep atmospheric dives that demand prohibitively high thrust levels to ensure success. Conversely, inaccuracies in trajectories close to the upper boundary are easier to correct, as these lead to an aerobraking maneuver, although a propulsive capture maneuver is required.

3.2. Performance Indices

This section introduces a series of typically used performance metrics for re-entry trajectories.

3.2.1. Thermal Environment

The corresponding aerocapture trajectory is computed for every candidate point, and the ${\dot{q}}_{t o t}$ profile is available. Hence, ${\dot{q}}_{t o t, m a x}$ is computed by locating the maximum value in the profile and Q by using Equation (2):

(2) $Q = \int_{t_{0}}^{t_{t}} {\dot{q}}_{t o t} d t$

where

t_{t}

is the total duration of the atmospheric flight.

3.2.2. Volumetric Constraints

To address fairing constraints and current technological limitations, this optimization process aims to minimize the volume of the aeroshell by reducing the “equivalent axisymmetric aeroshell radius”, $R_{B}$ , computed by rearranging Equation (3). This metric is preferred because it can represent any aeroshell configuration, from axisymmetric to complex high lift-to-drag designs, as long as $C_{D}$ and $S_{r e f}$ are known.

(3) $β = \frac{m}{C_{D} S} = \frac{m}{C_{D} π R_{B}^{2}}$

3.2.3. Orbital Corrective Maneuver

$Δ V$ represents the fuel needed to adjust the post-atmospheric pass orbit to the predetermined target orbit. In reality, each candidate trajectory will have its unique optimal orbital corrective maneuver composed of complex transfers. To simplify comparisons between all candidate trajectories, trivial orbit transfers are used. These include two-impulse Hohmann transfers to raise the orbit’s periapsis and match the desired target orbital eccentricity. In addition, a plane change burn at the apoapsis is assumed to match the target orbital inclination. The $Δ V$ required can be computed using Equation (4):

(4) $\begin{matrix} Δ V_{T O T A L} = \sqrt{μ} [\sqrt{\frac{1}{a_{t}}} + \sqrt{\frac{2}{a_{i} (e_{i} + 1)} - \frac{2}{a_{i} (e_{i} + 1) + a_{t}}} - \sqrt{\frac{2}{a_{t}} - \frac{2}{(a_{i} (e_{i} + 1) + a_{t})}} \\ - \sqrt{\frac{2}{a_{i} (e_{i} + 1)} - \frac{1}{a_{i}}}] + 2 \sqrt{μ (\frac{2}{a_{i} (e_{i} + 1)} - \frac{1}{a_{i}})} sin (\frac{i_{t} - i_{i}}{2}) \end{matrix}$

where a is the semi-major axis of the orbit and subscripts “i” and “t” denote the quantity

(\cdot)

evaluated at the post-atmospheric pass orbit (determined at AE) and the target operational orbit, respectively.

3.3. Single-Objective Function Optimization

This paper considered two formulations of optimization Problem 1. The first strategy is the most widely employed multi-objective optimization (MOO), where $i \neq 1$ , and a novel single-objective optimization formulation is presented, where $i = 1$ .

MOO is the preferred approach for re-entry trajectory optimization, addressing problems with multiple conflicting criteria that cannot be combined into a single objective. However, MOO methods face challenges, including non-unique optimal solutions and subjective decision-making. Additionally, achieving both convergence and diversity—how close and spread the candidate solutions are relative to the Pareto front—is complex. Consequently, a new optimization problem may need to be formulated to address these issues and improve the performance of the MOO process.

This paper opted to embed the weighting function in a single-objective cost function to enable the use of traditional gradient descent (GD) methods to solve the problem. The proposed objective function utilized the weighted least squares (WLS) approach and is presented below.

(5) $\begin{matrix} f (θ; x) = \sum_{i} w (i) {[\tilde{θ} (i) - {\tilde{θ}}_{t a r g e t} (i)]}^{2} \end{matrix}$

where

w

is the weight vector that includes the biases of each performance index, and the definition of

\tilde{θ}

follows.

To combine multiple performance metrics into a single-objective function, feature scaling was applied prior to computing the weighted least squares, ensuring all metrics were on a comparable scale. For instance, Q and $Δ V$ differ by five to seven orders of magnitude. Unlike the conventional use of feature scaling—typically employed to expedite the convergence of gradient descent algorithms—this normalization method was specifically implemented to preserve the relative weightings of the performance indices. By rescaling these indices to values between zero and one, via Equation (6), each metric’s contribution to the objective function was consistent with their weighting.

(6) $\tilde{θ} (i) = \frac{θ (i) - \underset{̲}{θ} (i)}{\bar{θ} (i) - \underset{̲}{θ} (i)}$

where

\bar{(\cdot)}

and

\underset{̲}{(\cdot)}

represent the upper and lower bounds of property

(\cdot)

. To perform accurate normalization when no constraints are imposed on the performance metrics, it is essential to understand the range of possible values each feature can take. However, computing this range can be challenging, and an approximate procedure for this is explained in Section 4. When constraints are applied to the performance metrics, the upper and lower bounds of the constraint can replace

\bar{θ} (i)

and

\underset{̲}{θ} (i)

. If strict enforcement of these constraints is required—though this is typically unrealistic during the preliminary design phase—barrier functions should be introduced.

4. The D-ASTRO Framework

D-ASTRO distinguishes itself through its pre-optimization phase, where the two-dimensional space spanned by $β$ and $γ_{A I}$ is studied to apply appropriate feature normalization and to implement a novel method of handling the nonlinear constraints introduced by $g_{l}$ and $g_{u}$ . The strategy used to find the optimal sizing of the spacecraft and corresponding optimal interplanetary orbit is defined as follows:

(1). A user-specified $β$ range is sampled with seven linearly spaced points, including the upper and lower limits.

(2). At each $β$ , a binary search determines the lower and upper corridor boundary $γ$ values, $γ_{l}$ and $γ_{u}$ , respectively.

(3). Normalization quantities $\bar{{\dot{q}}_{m a x}}$ , $\underset{̲}{{\dot{q}}_{m a x}}$ , $\bar{Q}$ , $\underset{̲}{Q}$ , $\bar{R_{B}}$ , $\underset{̲}{R_{B}}$ , $\bar{Δ V}$ , and $\underset{̲}{Δ V}$ can be computed:

(a). For unconstrained performance metrics, these values have already been computed when $γ_{l}$ and $γ_{u}$ were devised. Section 4.1 explains how this is performed.

(b). For constrained performance metrics, the upper and lower bounds of the constraint can replace $\bar{(\cdot)}$ and $\underset{̲}{(\cdot)}$ .

(4). Sixth-order polynomials are fitted to $γ_{l}$ and $γ_{u}$ with $β$ to obtain a one-to-one relationship. This process is explored in Section 4.2.

(5). After computing normalization quantities and polynomial fits, MATLAB’s fmincon optimizer is used, employing sequential quadratic programming with an active set strategy, as described by Gill et al. [32].

The Martian robust aerocapture corridor and polynomial fits for an injection orbit with $V_{e} = 3.5$ ${kms}^{- 1}$ are shown in Figure 3, with the normalization values listed in Table 1.

4.1. Computation of Normalization Quantities

As discussed in Section 3.3, computing the extremal values of the performance metrics in $θ$ is crucial for casting the optimization problem into an SOO form. In the following, we investigate which $γ$ and $β$ values may achieve these targets.

For volumetric constraints, the bounds are easily computed based on $β$ , with maximum and minimum volumes corresponding to the lower and upper limits of $β$ , respectively. Thermal metrics, such as ${\dot{q}}_{m a x}$ , increase with deeper atmospheric entries, so ${\dot{q}}_{m a x}$ reaches its maximum at the lower corridor boundary and upper $β$ limit, while the minimum occurs at the upper corridor boundary and lower $β$ limit. The same applies to Q, which is maximized for deeper entries. Determining the bounds of $Δ V$ is more complex, as post-atmospheric trajectory eccentricities can span from 0 to 1 within the aerocapture corridor. If circularization and plane-change maneuvers are neglected, $\underset{̲}{Δ V} \approx 0$ . Hence, we assign $\underset{̲}{Δ V} = 0$ , while $\bar{Δ V}$ depends on the operational orbit and requires an assessment across the $β$ ’s.

This analysis focuses on identifying the maximum and minimum values within the corridor. While the optimizer can explore $β$ and $γ$ values beyond those used for normalization, this discrepancy does not significantly impact the process. Normalization ensures that values remain on a similar scale, with values exceeding the estimated maximum normalized above 1 and smaller values becoming negative. However, the nonlinear constraints $g_{l}$ and $g_{u}$ enforce candidate points within the corridor, mitigating this issue.

4.2. Approximating Aerocapture Corridor Boundaries via Polynomial Fitting

The primary challenge in aerocapture optimization lies in handling the nonlinear constraints introduced by $g_{u}$ and $g_{l}$ in Problem 1. Current optimization techniques necessitate computing $g_{u}$ and $g_{l}$ at each internal candidate point, significantly increasing execution time despite efficient search algorithms. This is quantitatively characterized in Section 7.

Multiple Martian aerocapture corridors with various hyperbolic approach orbits were analyzed. Figure 2 illustrates two examined corridors. From this figure, it is evident that the upper and lower corridor boundaries exhibited an exponential trend over the considered $β$ range. Several candidate functions were considered to rapidly map $β$ to $g_{u}$ and $g_{l}$ , including exponential and polynomial curves.

(7) $f (x) = \sum_{i = 1}^{n} a_{i} e^{b_{i} x}$

Initial testing revealed the limitations of a primary exponential curve ( $n = 1$ in Equation (7)), leading to the adoption of a second-degree exponential fit ( $n = 2$ ), significantly enhancing predictive accuracy with error bounded by a $\pm 0.8 %$ Martian aerocapture. However, when terrestrial aerocapture was explored, the lower corridor boundary experienced a non-exponential-like trend, due to differing atmospheric structures of the planets at low altitudes. Hence, to expand the applicability of this technique to multiple celestial bodies, sixth-degree polynomial curves were examined due to the potential flat nature of some boundaries. As illustrated by Figure 4, the polynomial curves resulted in excellent predictive capabilities, exceeding the former’s accuracy. These observations were consistent for terrestrial aerocapture missions.

Exploration into the minimum number of points required for accurate fitting demonstrated that the accuracy of the fits was very sensitive to the selected fitting points. Randomly selecting seven points within the interest $β$ range demonstrated poor fit accuracy, as illustrated by Figure A1. However, when the fitting points were evenly distributed across the range, both fits exhibited increased accuracy with fewer fitting points. Specifically, fitting the sixth-degree polynomial with seven linearly spaced points proved to capture boundaries accurately, resulting in superior approximations in flat boundaries. Figure 4 illustrates the error between the actual boundary and the proposed fits, comparing the second-degree exponential fit with the sixth-degree polynomial fit, both fitted with 7 linearly spaced points.

5. Modelling

5.1. Atmospheric Flight Dynamics

At the preliminary design stages of an aerocapture-capable spacecraft, key geometric parameters such as the inertia matrix, center of mass, and center of pressure are undefined, limiting the analysis of rotational kinematics during re-entry. Consequently, D-ASTRO employs a three-degree-of-freedom (three-DOF) model. The equations describing flight over a spherical planet in a planet-fixed reference frame are as follows [33]:

(8) $\begin{matrix} \dot{δ} = \frac{V cos (γ) cos (χ)}{r cos (λ)} - Ω \\ \dot{λ} = \frac{V cos (γ) sin (χ)}{r} \\ \dot{r} = V sin (γ) \\ \dot{V} = - \frac{μ sin (γ)}{r^{2}} - \frac{D}{m} + Ω^{2} r cos (λ) [sin (γ) cos (λ) - cos (γ) sin (χ) sin (λ)] \\ \begin{matrix} \dot{γ} = (\frac{V}{r} - \frac{μ}{V r^{2}}) cos (γ) + \frac{L}{V m} + 2 Ω cos (χ) cos (λ) + \\ \frac{Ω^{2} r cos (λ)}{V} [cos (γ) cos (λ) + sin (γ) sin (χ) sin (λ)] \end{matrix} \\ \dot{χ} = \frac{- V}{r} cos (γ) cos (χ) tan (λ) + 2 Ω [tan (γ) sin (χ) cos (λ) - sin (λ)] \end{matrix}$

where r is the radial distance from the planet’s center to the vehicle, V is the planetocentric velocity,

γ

is the flight path angle,

δ

and

λ

are the longitude and latitude, respectively, and

χ

is the heading angle. Additionally, m denotes the spacecraft’s mass, L and D are the lift and drag forces, respectively, and

μ

and

Ω

represent the planet’s gravitational parameter and rotational period of the planet, respectively. The atmosphere is assumed to be at rest with respect to the planet.

Let the trajectory state in Equation (8) be $x = {[r, V, γ, δ, λ, χ]}^{T}$ . The initial condition of $x$ is assumed to be user-specified and coincident with the atmospheric interface (AI).

5.2. Aerodynamic Model

Accurately approximating the aerodynamic coefficients, $C_{L}$ and $C_{D}$ , is essential for properly simulating the trans-atmospheric trajectory. These coefficients depend on the angle of attack, which, in reentry applications, varies with the Mach number. However, given the high velocity of this maneuver, the spacecraft operates in a fully hypersonic regime, where the Mach independence principle allows these coefficients to be treated as constant [34].

Nevertheless, at high altitudes, gas rarefaction effects can cause deviations from the Mach independence principle. To address this, D-ASTRO enables users to define $C_{L}$ and $C_{D}$ as functions of Mach number $(M)$ , Reynolds number $(R e)$ , and Knudsen number $(K n)$ , enhancing accuracy, as outlined in Equations (9) and (10):

(9) $\begin{matrix} D = \frac{1}{2} ρ V^{2} S_{r e f} C_{D} (K n, M, R e) \end{matrix}$

(10) $\begin{matrix} L = \frac{1}{2} ρ V^{2} S_{r e f} C_{L} (K n, M, R e) or L = D (\frac{L}{D}) \end{matrix}$

where

ρ

is the freestream density, and

S_{r e f}

is the reference area of the aeroshell.

5.3. Aerothermal Model

To compute the heat flux, ${\dot{q}}_{m a x}$ , and heat load, Q, of the trajectory in a simplistic manner, the thermal model focuses on stagnation-point aerothermal heating, as this is the expected peak heat transfer region of the aeroshell. Thus, the Sutton–Graves relation, Equation (11), approximates the convective stagnation-point heating [35]. Furthermore, D-ASTRO incorporates the Tauber–Sutton relations, Equation (12), to account for radiative heat transfer, ${\dot{q}}_{r a d}$ , in cases where re-entry speeds are moderately high [36]:

(11) $\begin{matrix} {\dot{q}}_{c o n v e c t i v e} = k {(\frac{ρ}{R_{N}})}^{0.5} V^{3} \end{matrix}$

(12) $\begin{matrix} {\dot{q}}_{r a d i a t i v e} = C {R_{N}}^{a} ρ^{m} f (V) \end{matrix}$

where

R_{N}

is the nose radius, k, C, a, and m are atmospheric-specific constants, and

f (V)

is an empirically determined function that is also planet-specific.

In cases where wall temperature is a concern, it is common practice to employ the radiative equilibrium approach to approximate the equilibrium wall temperature $T_{W}$ [37]. This technique assumes that the surface temperature rapidly rises until it reaches a value at which the surface rejects all the convective heat flux, and sometimes radiation heat fluxes due to high-temperature shock layers, by black body re-radiation, leading to an energy balance between total heating, ${\dot{q}}_{t o t} = {\dot{q}}_{c o n v e c t i v e} + {\dot{q}}_{r a d i a t i v e}$ , and re-radiation cooling, ${\dot{q}}_{r e - r a d}$ , calculated using the Stefan–Boltzmann law, Equation (13). $T_{W}$ is approximated by Equation (14):

(13) ${\dot{q}}_{r e - r a d i a t i o n} = ϵ σ {(T_{W} - T_{\infty})}^{4}$

(14) $T_{W} = T_{\infty} + \sqrt[4]{\frac{{\dot{q}}_{t o t}}{ϵ σ}}$

where

ϵ

is the emissivity, and

σ

is the Boltzmann constant.

5.4. Martian Atmospheric Model

A reliable atmospheric model is essential for accurately calculating the vehicle’s aerodynamic forces and aerothermal heating during the aerocapture maneuver. This report investigated the Martian atmosphere, highlighting the specific model below.

Mars-GRAM [38] was employed to generate atmospheric data for Martian aerocapture maneuvers. A simplified table was generated from the extensive data available in the Mars-GRAM database, encompassing atmospheric variations concerning seasons, latitude and longitude positions, and numerous other variables. This table is based on a yearly average of equatorial properties at zero longitude. Additionally, a statistical analysis was conducted to obtain an accurate variation in atmospheric properties for $\pm 3 σ$ , which encompasses 99.73% of cases.

5.5. Computation of Keplerian Orbital Elements

Determining orbital elements is crucial for characterizing a spacecraft’s orbit after an atmospheric pass. The resulting orbit determines the amount of fuel required to obtain the target operational orbit via corrective maneuvers. D-ASTRO employs an algorithm for computing the Keplerian orbital elements from the Cartesian state of the spacecraft [39].

6. Results and Analysis

A Martian aerocapture mission was considered to assess the robustness and convergence properties of the proposed optimization strategy. The Martian atmosphere was modeled using data from Mars GRAM 2010 [38], an engineering-level atmospheric tool developed for mission design applications. D-ASTRO requires multiple input parameters, including spacecraft preliminary properties, insertion conditions, and target operational orbit. These input quantities are provided in Table 2 and Table 3.

The optimization algorithm was subjected to various initial conditions to assess its convergence properties. The range of considered initial conditions comprised extreme and random values both within and outside the aerocapture corridor, with some having limited physical meaning. Table 4 provides detailed information on the initial conditions used ( $β_{0}$ and $γ_{0}$ ).

6.1. Dependencies on Solver Initial Conditions

A bias vector $w$ assigning equal importance to all performance indices (denoted as $w = [1, 1, 1, 1]$ ) was employed for this examination. In Figure 5, level lines of the cost function are depicted within and around the aerocapture corridor. The region delineated by magenta solid lines represents the aerocapture corridor, defined by the novel polynomial fit approach discussed in Section 4. The figure illustrates that the novel objective function formulation was well suited for gradient descent (GD) algorithms, as it had no visible sharp discontinuities. Notably, areas where aerocapture failed, such as surface collisions in the lower left corner of the figure, incurred high costs, expediting solution convergence.

Based on Table 5 and Figure 5, the problem’s minimum was near $[13.20, - 9 . 76^{\circ}]$ . While Case 2 indicated a slightly different optimal point (explored further in Section 6.2, the algorithm consistently converged to the global minimum, even with initial conditions lacking physical significance (e.g., Cases 1, 3, and 4). The strong convergence properties of D-ASTRO were further demonstrated by the similarity in trajectory profiles shown in Figure 6.

6.2. Minimum Fuel Mission Design

The adoption of an aerocapture maneuver is often motivated by a desire to minimize fuel requirements to achieve a desired operational orbit. The minimum fuel trajectory for the proposed case study could be identified by selecting $w = [0, 1, 0, 0]$ . The same initial conditions as in the previous section were studied to assess convergence properties. Additionally, various target semi-major axis values, specifically $3500$ km, $4621$ km, and $5621$ km, were analyzed to explore the behavior of the optimal solution.

As illustrated by the $Δ V$ ( ${kms}^{- 1}$ ) cost contours in Figure 7 and the information in Table 6, this optimization problem yielded optima that were scattered over a broad range. However, it is crucial not to interpret this as a sign of poor convergence of D-ASTRO’s optimizer. Instead, understanding the underlying physics and formulation of the optimization problem is essential to justify this apparent discrepancy.

Figure 7 illustrates the presence of a low-cost valley in the $Δ V$ cost function, directing the global minimum along a specific direction. The valley-like structure of the $Δ V$ cost function characterized a typical aerocapture corridor, with the trough’s vertical position determined by the target a. Smaller values required deeper atmospheric entries, shifting the valley downward, while larger values led to shallower entries, moving it upward. The valley’s convex structure arose because one side required increasing orbital energy due to excessive energy depletion, while the other side required decreasing energy due to insufficient depletion. The valley in $Δ V$ did not span the entire $β$ range; the increase in cost at low $β$ was linked to plane-change maneuvers. Figure 8 illustrates the $Δ V$ contours without considering these maneuvers, revealing the anticipated valley that spans the entire $β$ range. The evolution of plane changes $Δ V$ $(Δ V_{p c})$ resulted from the forces acting on the spacecraft during atmospheric entry, which did not align with the orbital plane, thereby altering the inclination of the orbit’s plane.

To enhance convergence and mitigate the impact on other cases, the cost function for $Δ V$ must be revised. Comparing $Δ V$ contours in Figure 7a and Figure 8 shows that including $Δ V_{p c}$ significantly increased fuel requirement. Applying feature standardization to the total $Δ V = Δ V p c + Δ V_{H o h m a n n}$ flattened the valley as $Δ V_{p c} > > Δ V_{H o h m a n n}$ , causing $Δ V_{p c}$ to dominate the problem. To address this, a potential solution is to decompose the $Δ V$ into two components and apply separate standardizations.

Figure 9 depicts the optimal trajectory profiles. The extensive variation in both $β$ and $γ$ resulted in a multitude of trajectories. It is crucial to emphasize that despite discrepancies in the h-v profiles shown in Figure 9a, they all converged to a similar velocity (V) at the atmospheric exit, thereby leading to approximately comparable values of $Δ V$ . As anticipated, thermal quantities exhibited significant differences across all trajectories, reflecting the lack of assigned importance to these features in the optimization process.

6.3. Minimum Heat Load Mission Design

Ensuring mission success during re-entry requires equipping the spacecraft with a suitably robust thermal protective system (TPS). By minimizing the heat load experienced during the atmospheric trajectory, the TPS can be made more lightweight, allowing for an increased payload mass fraction.

As outlined by Equation (2), calculating stagnation-point heat loads requires complete trajectory knowledge. In practice, researchers often use Equation (15) for simplified heat load estimations during Entry, Descent, and Landing (EDL), based only on parameters like $β$ and $γ$ at AI. Utilizing the Sutton–Graves model for ${\dot{q}}_{c o n v}$ and an exponential atmospheric model, Allen and Eggers [40] formulated the complete expression, suggesting that heat load could be reduced by minimizing $β$ and maximizing $sin γ$ . This approach is later explored as a substitute for the current Q cost function for aerocapture missions.

(15) $Q \propto \sqrt{\frac{β}{sin γ}}$

Figure 10 illustrates the exploration of multiple initial conditions aimed at determining the optimal trajectory with the minimum heat load, revealing partially anticipated trends. Table 7 details the outcomes of the explored cases. All cases suggested $β$ values near the lower limit of feasible $β$ $(3 {kgm}^{- 2})$ , consistent with Equation (15). However, only Case 1 aligned with the minimum of Q predicted by Equation (15), with $γ$ at the lower boundary of feasibility, while other cases showed $γ$ near the upper boundary.

To investigate the disparity between the anticipated minimum of Equation (15) and the complete trajectory analysis, both heat loads were normalized with respect to their corresponding minima and maxima and are illustrated in Figure 11. This figure reveals differences in the contours of Q, particularly in the slopes with respect to $γ$ $(\partial Q / \partial γ$ ). The approximate expression in Equation (15) incorrectly suggested increments in Q with $γ$ $(\partial Q / \partial γ > 0)$ , whereas the complete trajectory analysis showed the opposite trend $(\partial Q / \partial γ < 0)$ . This discrepancy stemmed from the assumption of an EDL trajectory made by Allen and Eggers [40].

Optimizing for the minimum Q highlighted the importance of the optimizer’s step tolerance in D-ASTRO’s robustness. Initially, fmincon’s default convergence criteria yielded reliable and consistent results, as demonstrated in Figure 5 and Figure 7. However, minimizing Q initially showed significant variations in $γ_{f}$ due to the flat nature of the normalized Q ( $\tilde{Q}$ ) in the $γ$ direction, with $\tilde{Q}$ varying only slightly at the leftmost region of the corridor, of the order of $O (10^{- 4})$ . This was mitigated by reducing the function and step tolerance of the solver, causing D-ASTRO to identify the global minimum consistently, and in Case 1, a local minimum coincident with the feasibility boundary (see Figure A4). However, local minima could be avoided by selecting a reasonable initial condition. The decrement tolerance only benefited the convergence properties of the minimal Q trajectory, hindering the other problems with no changes in optimal trajectory for the other test cases. The resulting optimal trajectory profiles are shown in Figure 12, where the difference between the global and local minima trajectories can be seen.

6.4. Minimum Peak Heat Flux Trajectory

Similar to the parameter Q discussed earlier, the peak heat transfer rate ${\dot{q}}_{m a x}$ plays a crucial role in determining the design specifications of the TPS. This metric underwent an independent optimization process to validate the underlying physics principles and the cost function employed.

The minimization of peak heat transfer demonstrated strong convergence properties, as shown in Figure 13 and Table 8, where multiple initial points converged to the same optima. The variation in convergence quality between Q and ${\dot{q}}_{m a x}$ arose from the latter’s greater sensitivity to $γ$ , illustrated by curved contours in Figure 13 versus the quasi-vertical contours of Figure 10. Figure 14 shows indistinguishable trajectory profiles, highlighting the optimizer’s robust convergence.

The results presented in Table 8 align with established expectations that smaller values of $β$ lead to reduced ${\dot{q}}_{m a x}$ . However, the practical implementation of this solution is limited, as $β = 3.00 {kgm}^{- 2}$ corresponds to an unfeasibly small value for a realistic mission, since unprecedented aeroshell surface areas (A) would be required to achieve that $β$ . However, introducing volumetric considerations in the optimization (non-zero first entry of $w$ ) would lead to a more realistic $β$ . The resulting optimal points of combining thermal loads and volumetric considerations is explored in the following Section, where $w = [1, 0, 1, 0]$ and $[5, 0, 1, 0]$ were considered.

6.5. Direct Comparison with MOO

This section compares D-ASTRO with a basic MOO algorithm using MATLAB’s fgoalattain optimizer. Multiple initial conditions were again tested to evaluate the solver’s convergence dynamics and stability. Figure 15 shows the convergence dynamics of both approaches, revealing significant disagreement in optimal solutions. Table 9 presents an additional comparison, incorporating data from Table 5.

From Figure 15, it is clear that the MOO optimizer experienced similar poor convergence in the leftmost region as D-ASTRO when only the minimal heat load trajectory was considered, yielding one of two possible solutions: a local or global minimum. The fact that the MOO strategy suffered from the same convergence issues as D-ASTRO, even when subjected to a wide range of initial conditions, suggested that the subjective decision-making inherent in the MOO strategy was not responsible for the disagreement in solutions. Ill conditioning of the optimization problem was not the reason behind that discrepancy, as the MOO strategy was conducted with and without D-ASTRO’s normalization technique and yielded identical results.

Further analysis revealed that the MOO strategy prioritized the optimization of Q and ${\dot{q}}_{m a x}$ over $Δ V$ and $R_{B}$ . As illustrated by Figure 10 and Figure 13, the heat load and peak heat rate cost functions exhibited a strong convex behavior in the $β$ direction, facilitating the rapid convergence of GD algorithms to the leftmost regions. In contrast, the volumetric consideration cost function experienced an inverse proportionality to $β$ , leading to rapid convergence to the rightmost regions. Finally, the $Δ V$ cost function experienced a convex-like behavior in the $γ$ direction, causing the optimality valley explored in Section 6.2.

The dynamics of the individual metrics cost contours indicated three convergence behaviors: one toward the left region of the corridor, another toward the rightmost regions, and a third toward the corridor’s center. We conclude that the MOO optimizer was attempting to globally optimize Q and ${\dot{q}}_{m a x}$ , at the expense of $R_{B}$ , to later find the optimal $Δ V$ in regions where Q and ${\dot{q}}_{m a x}$ were quasi-globally optimal, as $Δ V$ had a relatively low cost over the entire $β$ range. The resulting optimal trajectories are depicted in Figure 16.

For completeness, all optimal trajectory profiles presented in this paper are illustrated in Figure 17 to aid the comparison between the different optima.

7. Computational Performance and Considerations of D-ASTRO

This section examines the computational performance and considerations within the D-ASTRO framework. The analysis is divided into two parts: the first addresses the computational costs and considerations of the binary search algorithm used to compute the robust aerocapture corridor; the second evaluates the performance of the optimization process.

7.1. Dependency of Corridor Fits on Binary Search Tolerance

To balance accuracy and computational efficiency, a study was conducted to analyze how the accuracy of the computed corridor boundaries varied with different tolerance levels. To evaluate this, we computed the robust coarse aerocapture corridor depicted in Figure 3 using varying binary search tolerances. The overall accuracy of the resulting polynomial fits was assessed by calculating the Euclidean L-2 norm between actual corridor boundaries, computed from a finely discretized aerocapture corridor (570 points), and those computed using the fits. Figure 18 illustrates the evolution of polynomial fit accuracy for both lower and upper corridor boundaries over the considered $β$ range. Additionally, Figure 19 presents the cumulative accuracy via the Euclidean L-2 norm approach and the execution time of the corridor, examining a Martian aerocapture maneuver. This figure shows that the Euclidean norm reached a plateau at a tolerance of approximately $2 \times 10^{- 5}$ in both scenarios, indicating that more accurate fits were not achievable with a further decrement in tolerance. Hence, for expedited solutions of D-ASTRO, a binary search tolerance of $2 \times 10^{- 5}$ was selected as the optimal balance between computational efficiency and accuracy. For the selected tolerance level, 31 trajectory simulations were required to compute $g_{l}$ and $g_{u}$ , with the number of simulations being unaffected by the chosen $β$ .

7.2. Execution Time of Optimizers

To evaluate the optimizers’ performance, 10 iterations of the optimization process were performed using a consistent initial guess ( $β_{0} = 40.60 {kgm}^{- 2}$ and $γ_{0} = - 10.113 °$ ) to minimize execution timing outliers. Various weightings ( $w$ ) were also considered to further explore the convergence properties of MOO and D-ASTRO, highlighting the challenges of the multi-objective approach in optimizing these problems. Table 10 presents the execution time, number of function calls, time per function call, and the suggested optima for both methods. This table indicates a significant dependency of execution time on the weight vector and optimization strategy employed. However, the time required per function call was roughly consistent across all cases, with the MOO function requiring marginally less time than D-ASTRO’s. This difference arises from the distinct optimizers used, as fgoalattain and fmincon are fundamentally different algorithms.

In general, D-ASTRO converged more quickly to the global optimum compared to the MOO strategy due to a reduced number of function evaluations, except in the case where $w = [1, 1, 1, 1]$ . Under that condition, D-ASTRO required more functions evaluations, and hence time, but yielded a solution that reflected a more balanced trade-off among all metrics, in contrast to the MOO strategy, which prioritized Q and ${\dot{q}}_{m a x}$ . The same weighting where provided to D-ASTRO and fgoalattain to obtain comparable results. Solutions with different weightings can be seen in Table 10.

Comparing rows 3 and 4 reveals that different weightings did not impact the solution provided by MOO, whereas they clearly affected D-ASTRO’s solution. This behavior was consistent across multiple initial guesses and weightings. The only modification that enabled the MOO solution to depend on the weight vector was performing feature regularization. These findings underscore the superior capabilities of the D-ASTRO strategy over the multi-objective approach.

As described in Section 4, determining the robust aerocapture corridor is essential to compute the scaling quantities required for normalizing the performance metrics and performing the single-objective optimization. Therefore, the time required to compute the robust corridor must be added to D-ASTRO’s optimizer execution times listed in Table 10 to obtain the total execution time. From Figure 19 and recalling that a binary search tolerance of $2 \times 10^{- 5}$ was used, the overhead time required for computing the robust corridor was approximately 2.2 s.

To quantify the computational benefits of using polynomial fits over a binary search for computing point feasibility, D-ASTRO and MOO strategies were executed with a binary search algorithm instead of polynomial fits to evaluate $g_{l}$ and $g_{u}$ . Table 11 reports results of performing the same optimization as in Table 10, including that modification. The significant difference in execution times was evident, with one to two orders of magnitude increments. Although there was a slight disagreement in the optima suggested by D-ASTRO when using polynomial fits versus a binary search, this discrepancy resulted in minimal changes to performance metrics. Table 12 shows the percentage difference for the first three rows. This slight difference was likely introduced by the accuracy of the constraint gradients employed by fmincon, as slight variations in this quantity could lead to different optimization paths.

8. Conclusions

This paper presented a novel framework (D-ASTRO) for the rapid optimization of aerocapture maneuvers for spacecraft mission design. D-ASTRO efficiently identifies trajectories that optimize a spacecraft’s flight path angle (γ) at the atmospheric interface and its ballistic coefficient (β), while meeting user-specified target orbits and design constraints. Computational efficiency was significantly improved by the use of sixth-degree polynomials to model the aerocapture corridor boundaries. By eliminating the need for a double binary search at each candidate point, execution time per function call was reduced from 300–400 ms to around 10 ms.

Results from Martian aerocapture test missions demonstrated D-ASTRO’s capabilities, and the feasibility of the proposed solutions were validated. It was shown that D-ASTRO reliably converged to the global optimum, consistently achieving lower costs than multi-objective optimal trajectories, even under challenging initial conditions.

A significant limitation of D-ASTRO is its reliance on a prescribed lift-to-drag ratio and hyperbolic excess velocity. Future iterations should explore optimizing the hyperbolic excess velocity of the approach orbit and treating the lift-to-drag ratio as a variable, broadening the algorithm’s utility. Future research may involve coupling D-ASTRO with an interplanetary trajectory design tool to ensure coherence between aerocapture and interplanetary trajectories by accurately determining the atmospheric interface. This integration would also enable the use of more advanced atmospheric models, further improving the accuracy of the framework.

Author Contributions

Conceptualization, S.U.A.; methodology, S.U.A.; software, S.U.A.; validation, P.B. and S.U.A.; formal analysis, P.B. and S.U.A.; investigation, S.U.A.; resources, P.B.; data curation, S.U.A.; writing—original draft preparation, S.U.A.; writing—review and editing, P.B. and S.U.A.; visualization, S.U.A.; supervision, P.B.; project administration, P.B. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Data and Algorithm is available in the following GitHub repository.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

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.

Figures and Tables

Figure 1 Example of Martian aerocapture corridors with post-atmospheric pass orbital eccentricity contours: (a) ballistic re-entry profile, $L / D = 0.0$ ; (b) lifted re-entry profile, $L / D = 0.2$ .

Figure 2 Martian and terrestrial aerocapture corridors for multiple hyperbolic excess velocities, $V_{e}$ . The polynomial fits for low, mean, and high-density corridor boundaries are shown: (a) Martian aerocapture corridor with $V_{e} = 3.0$ ; (b) Martian aerocapture corridor with $V_{e} = 6.0$ .

Figure 3 Coarse robust Martian aerocapture corridor for $V_{e} = 3.5$ ${kms}^{- 1}$ used to compute normalization quantities listed in Table 1 and sixth-degree polynomial fits.

Figure 4 Percentage error between finely discretized corridor limits (570 points) and computed values using exponential and polynomial constructed from 7 points. (a) Mars, $V_{e}$ = 1.0 ${kms}^{- 1}$ . (b) Mars, $V_{e}$ = 3.0 ${kms}^{- 1}$ . (c) Mars, $V_{e}$ = 6.0 ${kms}^{- 1}$ .

Figure 5 Convergence property of D-ASTRO for overall optimal mission design, $w = [1, 1, 1, 1]$ , with magenta lines delineating the boundaries of the aerocapture corridor. (a) Convergence of $y_{0}$ to $y_{f}$ . (b) Zoom into global minimum region.

Figure 6 Aerocapture trajectory profiles resulting from overall optimal mission design, $w = [1, 1, 1, 1]$ . (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 7 Convergence property of D-ASTRO for minimum fuel mission design, $w = [0, 1, 0, 0]$ , with magenta lines delineating the boundaries of the aerocapture corridor. (a) Convergence of $y_{0}$ to $y_{f}$ . (b) Zoom into global minimum region.

Figure 8 $Δ V$ contours for minimal fuel trajectory neglecting plane-change maneuvers, $Δ V_{p c} = 0$ , and $a = 4620.5$ km.

Figure 9 Aerocapture trajectory profiles resulting from minimal fuel mission design, $w = [0, 1, 0, 0]$ . (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 10 Convergence property of D-ASTRO for minimum heat load mission design, $w = [0, 0, 1, 0]$ , with magenta lines delineating the boundaries of the aerocapture corridor. (a) Convergence of $y_{0}$ to $y_{f}$ . (b) Zoom into global minimum region.

Figure 11 Normalized heat load, $\tilde{Q}$ . (a) Complete trajectory analysis. (b) Predicted by Equation (15).

Figure 12 Aerocapture trajectory profiles resulting from minimal heat load mission design, $w = [0, 0, 1, 0]$ . (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 13 Convergence property of D-ASTRO for minimum peak heating rate mission design, $w = [0, 0, 0, 1]$ , with magenta lines delineating the boundaries of the aerocapture corridor. (a) Convergence of $y_{0}$ to $y_{f}$ . (b) Zoom into global minimum region.

Figure 14 Aerocapture trajectory profiles resulting from minimal peak heating rate mission design, $w = [0, 0, 0, 1]$ . (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 15 Convergence property of MOO strategy for overall optimal mission design (blue) compared with D-ASTRO (green), $w = [1, 1, 1, 1]$ , with magenta lines delineating the aerocapture boundaries.

Figure 16 Aerocapture trajectory profiles resulting from overall optimal mission design using MOO strategy, $w = [1, 1, 1, 1]$ . (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 17 Aerocapture trajectory profiles resulting from all optimal mission designs presented in this study. Initial conditions correspond to those that result in the lowest cost. (a) Altitude vs. velocity. (b) Q and $\dot{q}$ vs. time.

Figure 18 Accuracy evolution of sixth-degree polynomial curves as a function of binary search tolerance used for computing fitting points. (a) Corridor boundaries. (b) Lower-bound fit. (c) Upper-bound fit.

Figure 19 Evolution of Euclidean L-2 norm of upper- and lower-boundary fits and execution time of robust coarse corridor for a Martian aerocapture with $V_{e} = 3.5 {kms}^{- 1}$ .

Table 1

Minimum and maximum values of considered metrics for a Martian aerocapture corridor for $V_{e} = 3.5 {kms}^{- 1}$ .

Variable	Minimum	Maximum
$R_{B}$ $(m)$	1.15	5.15
$Δ V$ $({kms}^{- 1})$	0.00	3.56
Q $({MJm}^{- 2})$	3.42	41.53
${\dot{q}}_{m a x}$ $({Wcm}^{- 2})$	3.45	46.32

Table 2

Martian aerocapture simulation parameters.

Parameter	Value		Parameter	Value
Planetary parameters			Insertion trajectory parameters
Radius of Mars, $R_{P}$	3390	km	Hyperbolic excess velocity, $V_{e}$	3.5	${kms}^{- 1}$
Gravitational parameter, $μ$	$4.28 \times 10^{4}$	${km}^{3} s^{- 2}$	Atmospheric interface altitude, $h_{A I}$	125	$km$
Angular frequency of Mars, $Ω$	$7.08 \times 10^{- 5}$	${rads}^{- 1}$	Initial Radius, $R_{0}$	$R_{P} + h_{A I}$	$km$
Vehicle parameters			Initial velocity, $V_{0}$	$\sqrt{\frac{2 μ}{R_{0}} + V_{e}^{2}}$	${kms}^{- 1}$
Mass of vehicles, m	400	$kg$	Initial longitude, $δ_{0}$	0.5798	$\deg$
Drag coefficient, $C_{D}$	1.6		Initial latitude, $λ_{0}$	34.49	$\deg$
Lift-to-drag ratio, $L / D$	0.2		Initial heading angle, $χ_{0}$	−18.24	$\deg$
Nose-radius-to-body-radius ratio, $R_{N} / R_{B}$	0.5		Insertion flight path angle accuracy, $Δ γ$	$\pm 0.2$	$\deg$

Table 3

Target operational orbit.

Target parameter	Semi-major axis, a	Inclination, i	Eccentricity, e
Value	4621 $km$	70 $\deg$	0.05

Table 4

Initial conditions used.

Case	Initial Conditions
Case	$β_{0}$ ( ${Kgm}^{- 2}$ )	$γ_{0}$ (deg)
1	59.0	−60.000
2	40.6	−10.113
3	59.0	−0.057
4	3.0	−0.057
5	39.6	−10.852
6	37.8	−10.502

Table 5

Convergence properties of optimization algorithm with different initial conditions.

Case	Cost	Optimized Parameters		Performance Metrics
Case	Cost	$β_{f}$ ( ${Kgm}^{- 2}$ )	$γ_{f}$ ( $\deg$ )	$R_{B}$ ( $m$ )	$Δ V$ ( ${kms}^{- 1}$ )	$Q$ ( ${MJm}^{- 2}$ )	${\dot{q}}_{\max}$ ( ${Wcm}^{- 2}$ )
1	0.5280	13.03	−9.756	2.472	0.777	12.608	12.805
2	0.5280	13.69	−9.784	2.411	0.776	13.094	13.301
3	0.5303	13.08	−9.758	2.467	0.778	12.649	12.846
4	0.5283	13.10	−9.759	2.465	0.778	12.665	12.862
5	0.5280	13.17	−9.762	2.458	0.778	12.714	12.913
6	0.5283	13.16	−9.762	2.459	0.778	12.704	12.902

Table 6

Optimal results for minimum fuel trajectory.

Case	Cost	$β_{f}$ ( ${Kgm}^{- 2}$ )	$γ_{f}$ ( $\deg$ )	$Δ V$ ( ${kms}^{- 1}$ )
1	0.1092	40.85	−10.393	0.7163
2	0.1064	48.46	−10.487	0.7070
3	0.1122	34.17	−10.294	0.7259
4	0.1392	7.58	−9.445	0.8087
5	0.1078	44.40	−10.439	0.7117
6	0.1105	37.80	−10.350	0.7205

Table 7

Optimal results for minimum heat load trajectory.

Case	Cost	$β_{f}$ ( ${Kgm}^{- 2}$ )	$γ_{f}$ ( $\deg$ )	Q ( ${MJm}^{- 2}$ )
1	$4.47 \times 10^{- 4}$	3.00	−9.6289	4.221
2	$3.63 \times 10^{- 4}$	3.00	−8.7822	4.141
3	$3.62 \times 10^{- 4}$	3.00	−8.7803	4.140
4	$3.67 \times 10^{- 4}$	3.00	−8.7934	4.145
5	$3.64 \times 10^{- 4}$	3.00	−8.7862	4.142
6	$3.66 \times 10^{- 4}$	3.00	−8.7908	4.144

Table 8

Optimal results for minimal ${\dot{q}}_{m a x}$ trajectory.

Case	Cost	$β_{f}$ ( ${Kgm}^{- 2}$ )	$γ_{f}$ (deg)	${\dot{q}}_{\max}$ ( ${Wcm}^{- 2}$ )	$Δ γ$ (deg)
1	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000
2	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000
3	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000
4	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000
5	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000
6	$7.92 \times 10^{- 3}$	3.0	−8.780	4.0754	0.2000

Table 9

Convergence properties of MOO algorithm and D-ASTRO for overall optimal trajectory.

Case	Strategy	$β_{f}$ ( ${Kgm}^{- 2}$ )	$γ_{f}$ ( $\deg$ )	$R_{B}$ (m)	$Δ V$ ( ${kms}^{- 1}$ )	Q ( ${MJcm}^{- 2}$ )	${\dot{q}}_{\max}$ ( ${Wcm}^{- 2}$ )
1	MOO	3.000	−9.526	5.150	1.969	4.225	4.723
1	D-ASTRO	13.03	−9.756	2.472	0.779	12.608	12.805
2	MOO	3.007	−9.635	5.145	2.062	4.228	4.818
2	D-ASTRO	13.69	−9.784	2.411	0.776	13.094	13.301
3	MOO	3.000	−9.635	5.150	2.063	4.221	4.810
3	D-ASTRO	13.08	−9.758	2.467	0.778	12.649	12.846
4	MOO	3.000	−8.780	5.150	0.884	4.140	4.075
4	D-ASTRO	13.10	−9.759	2.465	0.778	12.665	12.862
5	MOO	3.000	−8.780	5.150	0.884	4.140	4.075
5	D-ASTRO	13.17	−9.762	2.458	0.778	12.714	12.913
6	MOO	3.000	−9.590	5.150	2.025	4.222	4.774
6	D-ASTRO	13.16	−9.762	2.459	0.778	12.704	12.902

Table 10

Computational performance and suggested optima of D-ASTRO and MOO strategies with $β_{0} = 40.60 {kgm}^{- 2}$ and $γ_{0} = - 10.113 °$ for the Martian test case.

Weights $(w)$	Strategy	Execution Time (s)	Function Counts	Time per Function Call (ms)	$β_{f}$	$γ_{f}$ $(^{\circ})$
$[1, 1, 1, 1]$	D-ASTRO	1.889	1820	10.381	13.23	−9.765
$[1, 1, 1, 1]$	MOO	0.381	360	9.979	3.03	−9.638
$[1, 1, 0, 0]$	D-ASTRO	2.512	2010	12.496	43.97	−10.434
$[1, 1, 0, 0]$	MOO	3.066	3020	10.151	40.61	−10.113
$[1, 0, 1, 0]$	D-ASTRO	0.564	500	11.286	15.58	−9.653
$[1, 0, 1, 0]$	MOO	2.788	3010	9.263	12.39	−9.853
$[5, 0, 1, 0]$	D-ASTRO	0.506	470	10.773	37.84	−10.101
$[5, 0, 1, 0]$	MOO	2.927	3010	9.723	12.39	−9.853
$[10, 5, 1, 0]$	D-ASTRO	2.033	1630	12.473	41.49	−10.401
$[10, 5, 1, 0]$	MOO	2.767	3010	9.192	12.40	−9.163

Table 11

Computational performance and suggested optima of D-ASTRO and MOO strategies without employing polynomial fits to evaluate $g_{l}$ and $g_{u}$ . $β_{0} = 40.60 {kgm}^{- 2}$ and $γ_{0} = - 10.113 °$ used as initial guess for Martian test case.

Weights $(w)$	Strategy	Execution Time (s)	Function Counts	Time per Function Call (ms)	$β_{f}$	$γ_{f}$ $(^{\circ})$
$[1, 1, 1, 1]$	D-ASTRO	65.1	202	322	12.22	−9.719
$[1, 1, 1, 1]$	MOO	90.2	304	297	3.00	−9.634
$[1, 1, 0, 0]$	D-ASTRO	48.9	162	302	40.60	−10.390
$[1, 1, 0, 0]$	MOO	93.1	304	306	21.78	−8.501
$[1, 0, 1, 0]$	D-ASTRO	71.0	202	3528	16.00	−9.657
$[1, 0, 1, 0]$	MOO	109	301	363	12.39	−9.853

Table 12

Percentage difference in performance metrics when fits and binary search were used to assess the feasibility of candidate points.

Weights $(w)$	Strategy	Percentage Difference (%)
Weights $(w)$	Strategy	$R_{B}$	$Δ V$	$Q$	${\dot{q}}_{\max}$
$[1, 1, 1, 1]$	D-ASTRO	4.036	0.566	−5.781	−5.798
$[1, 1, 1, 1]$	MOO	0.454	0.097	−0.682	−0.673
$[1, 1, 0, 0]$	D-ASTRO	4.066	0.604	−5.842	−5.804
$[1, 1, 0, 0]$	MOO	36.545	92.352	−58.701	−62.593
$[1, 0, 1, 0]$	D-ASTRO	−1.327	0.526	1.895	1.672
$[1, 0, 1, 0]$	MOO	0.000	−0.046	−0.001	−0.005

Appendix A. Sensitivity of Boundary Fit on Seeding Points

Figure A1 Sensitivity of exponential and polynomial fit to seeding points.

Appendix B. Minimum Heat Load Mission Design: Optimization with Default Fmincon Convergence Criteria

Figure A2 Convergence property of D-ASTRO for minimal heat load mission design, $w = [0, 0, 1, 0]$ , with defaultfmincon convergence criteria. (a) Convergence of $y_{0}$ to $y_{f}$ . (b) Zoom into global minimum region.

Figure A3 Aerocapture trajectory profiles resulting from minimal heat load mission design, $w = [0, 0, 1, 0]$ , with defaultfmincon convergence criteria. (a) Altitude vs. velocity. (b) Altitude and $\dot{q}$ vs. time.

Appendix C. Minimum Heat Load Trajectory: Local Minima vs. Global Minima at β = 3

Figure A4 Local vs. global minima in minimum heat load trajectory, blue cross, and green dot, respectively.

© 2025 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/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.

Framework for the Multi-Objective Design Optimization of Aerocapture Missions

Content area

Abstract

Full text

1. Introduction

2. Aerocapture Corridor

3. Optimal Aerocapture Maneuver via Parameter Optimization

3.1. Formulation of Optimization Problem

3.2. Performance Indices

3.2.1. Thermal Environment

3.2.2. Volumetric Constraints

3.2.3. Orbital Corrective Maneuver

3.3. Single-Objective Function Optimization

4. The D-ASTRO Framework

4.1. Computation of Normalization Quantities

4.2. Approximating Aerocapture Corridor Boundaries via Polynomial Fitting

5. Modelling

5.1. Atmospheric Flight Dynamics

5.2. Aerodynamic Model

5.3. Aerothermal Model

5.4. Martian Atmospheric Model

5.5. Computation of Keplerian Orbital Elements

6. Results and Analysis

6.1. Dependencies on Solver Initial Conditions

6.2. Minimum Fuel Mission Design

6.3. Minimum Heat Load Mission Design

6.4. Minimum Peak Heat Flux Trajectory

6.5. Direct Comparison with MOO

7. Computational Performance and Considerations of D-ASTRO

7.1. Dependency of Corridor Fits on Binary Search Tolerance

7.2. Execution Time of Optimizers

8. Conclusions