Abstract

In the field of aerospace, solving the boundary problem associated with the parachute-capsule system remains a big challenge. The conventional Monte Carlo method proves inadequate for acquiring comprehensive boundary information. To address this issue, this paper introduces a novel tube prediction scheme by leveraging the natural geometric characteristics of the reachable tube and employing a multiobjective optimization strategy. Initially, a multibody dynamic model with nine degrees of freedom was established and verified by the airdrop test data to ensure the accuracy and reliability of the model. Subsequently, the Sobol sensitivity analysis method was employed to assess uncertain factors that affect the deceleration phase of the reentry capsule. These factors are then utilized to determine the optimization parameters for the multiobjective optimization model. Ultimately, the multiobjective evolutionary algorithm based on decomposition was employed to solve the multiobjective optimization model, and the geometric boundary of the tube corresponds to the Pareto front of the multiobjective optimization. The proposed methodology was validated through a simulation experiment utilizing the Chang’e-5 reentry capsule as an engineering case. The experimental results unequivocally demonstrate the superior accuracy of our approach in predicting the boundary of the reachable tube compared to the Monte Carlo method. This research serves as a valuable reference for calculating reachable tubes in practical engineering scenarios and can be effectively applied to spacecraft search and rescue operations during the reentry phase.

1. Introduction

With the rapid advancement of aerospace technology and an increasing number of deep-space exploration activities in China, research on the recovery and landing technology of reentry capsules has become significantly important. This recovery methodology prioritizes the safety of the landing sites. Predicting the recovery trajectory ensures that all items land in a secure area, away from observers and buildings, during the recovery process. Anticipating the landing positions of the reentry capsule components is crucial for ensuring the safety of the recovery process, the hardware, and individuals in the vicinity of the landing sites [1].

As a commonly employed pneumatic deceleration device, parachutes have found extensive applications in the recovery of spacecraft reentry capsules because of their lightweight and efficient characteristics. Parachutes are commonly used for deceleration, descent, stabilization, flight termination, and landing [2, 3]. However, the unpredictability associated with a parachute system remains relatively high owing to factors such as weight, additional mass, drag coefficient, and wind, making accurate predictions challenging. This complexity poses significant challenges for predicting reachable trajectories.

In recent years, the prediction of the landing points of reentry capsules has reached an immature stage. The conventional approach involves generating scattered information about the airdrop point through multiple calculations using the Monte Carlo method, often leading to challenges such as time-consuming computations and complex data processing. For instance, the consolidated airdrop tool (CAT; the USA) calculates the scattered area of landing locations for airdrops, and Zhou et al. employed Monte Carlo simulations to predict the descent trajectory of controllable parafoils [4]. Cao and Wei incorporated interfering factors into the flight dynamic model and employed the Monte Carlo method to forecast the flight trajectory and landing point distribution of the parachute system [5]. Guo et al. conducted predictions of the flight trajectory for the Shenzhou-7 spacecraft by employing various wind profiles [6]. Spencer and Braun conducted a Monte Carlo study on the landing point dispersion of the Mars Pathfinder spacecraft [7]. NASA developed the CPAS Sasquatch software for forecasting airdrop landing positions [8] to enhance the accuracy of predictions and ensure range safety. In 2013, they released an improved version of the Sasquatch tool. The principle behind this approach is to calculate the recovery trajectories for various groups, factoring in only the uncertain elements. The landing positions of all the generated trajectories were then combined to establish the boundary of the scattered region. A novel prediction method for reachable tubes based on optimal control was introduced in the literature [9]. This method established a polar coordinate system by placing the origin in a reachable set. Then, it defines a series of directional vectors originating from the origin and traversing the entire state space at specific angle intervals. The farthest reachable points are calculated for each direction. This approach, derived from the aircraft flight envelope calculation method, offers a new perspective for predicting the recovery trajectories. Furthermore, this approach can be coupled with optimization techniques to address the reachable region under uncertain disturbances.

The prediction of the recovery-reachable tube boundary for reentry capsules plays a pivotal role in airdrop mission design and planning systems. Currently, predicting the recovery of the reachable tube boundary of the reentry capsule involves generating scattered information about the landing point through multiple simulation trials and subsequently determining the reachable boundary. However, this method does not consider the interactions among the various influencing factors. To enhance prediction efficiency, ensure data accuracy, reduce calculation time, and ultimately improve the performance of the planning system, this paper proposes a novel method for predicting the reachable tube boundary of a parachute recovery system. This method utilizes the concept of multiobjective optimization by transforming the recovery reachable boundary at a fixed height into a multiobjective optimization problem to overcome the limitations of the existing approach.

The remainder of this paper is organized as follows. Section 2 introduces the nine-degree-of-freedom nonlinear model of the Chang’e-5 (CE-5) reentry capsule and verifies its accuracy and effectiveness. Section 3 employs the Sobol method to analyze the sensitivity of uncertain perturbation parameters and determine the parameters with the most significant impact on the output. In Section 4, an accessibility theory is introduced, and a prediction method for capsule reentry flight is proposed using multiobjective optimization. Section 5 elaborates on the simulation calculations and compares and discusses the prediction results. Finally, Section 6 concludes the study.

2. Models

2.1. Model Hypothesis

Before simulating a dynamic system on a computer, it is crucial to develop a mathematical model that accurately describes the physics of the system. However, owing to the complexity and unknown real-world effects of the system, simplification of the model is often necessary to make it more computationally tractable. Although this simplification can limit the effectiveness of the model, a balance between accuracy and computational feasibility is necessary. In this study, the main simplification involved considering the entire parachute and reentry capsule as rigid bodies with six degrees of freedom. Subsequently, a nonlinear mathematical model with nine degrees of freedom (9-DOF) was established to accurately depict the dynamic characteristics of the recovery system. In this model, factors such as aerodynamic force, aerodynamic moment, and gravity are transformed into a body-fixed coordinate system. Additionally, the effect of unsteady drag on the parachute in air was considered as an additional mass. To summarize, the model is based on several simplified assumptions as follows: (1)The parachute system is assumed to be in a stable descent state, with both the parachute and reentry capsule treated as rigid bodies connected at the point at the bottom of the capsule(2)The parachute is allowed to yaw, pitch, or roll relative to the capsule(3)The fixed reference frame of the Earth is considered an inertial frame(4)The aerodynamic effects of the reentry capsule are neglected

2.2. Model Coordinate System

In the 9-DOF model of the parachute-capsule system, an earth-fixed coordinate system was established, along with corresponding fixed coordinate systems on both the parachute and reentry capsule, as shown in Figure 1. Their details are provided below. (1)Earth-Fixed Coordinate System. The origin is the projection point of the reentry capsule on the ground when the recovery program starts. The axis is oriented in the downward direction of a vertical gravity vector(2)Parachute-Fixed Coordinate System. The origin is at the hinge point. The axis is aligned along the symmetry axis of the parachute(3)Capsule-Fixed Coordinate System. The origin is also at the hinge point. The axis points down along the symmetry axis of the reentry capsule

In the parachute-fixed coordinate system, the velocity of the system at the hinge point is denoted as , and the rotational angular velocities of the main parachute and reentry capsule are recorded as and , respectively. The transformation matrix from to is described by the pitch angle , yaw angle , and roll angle , and the corresponding transformation matrix is denoted as . Similarly, Euler angles , , and describe the transformation matrix from to , which is represented as . The transformation matrix from to is denoted as .

2.3. Multibody Dynamic Model

The parachute-capsule system is characterized by assigning a mass to the parachute and to the reentry capsule. The generalized mass and inertia of the parachute are defined as follows. where , , and represent the moment of inertia of the parachute at the origin and , , , and represent the additional mass of the parachute [10].

The aerodynamic force and moment acting on the parachute are denoted as and , respectively, whereas those acting on the reentry capsule are denoted as and , respectively. The gravitational acceleration is represented as . The dynamic equation of the main parachute-reentry capsule two-body system is established at the hinge point and expressed in matrix form as follows: where is a identity matrix. , , and represent the accelerations of , , and , respectively. and represent the vectors extending from the hinge point to the centroid of the reentry capsule and parachute, respectively. and are antisymmetric matrices constructed based on and , respectively, and they are defined as follows:

By combining (2)–(4) with the kinematic equation, we obtain where denotes the derivative of , which denotes the displacement at the hinge point , and

Finally, the 9-DOF dynamic equations of the parachute-capsule system were obtained as follows: where

In summary, this dynamic model can be used to simulate and analyze the dynamic characteristics of the main parachute-capsule system, such as attitude, velocity, and acceleration.

2.4. Model Verification

The recovery process of the CE-5 capsule is depicted in Figure 2. The recovery program was initiated when the reentry capsule altitude reached 11 km. The system control assembly issued instructions to eject the parachute container cover, deploy the parachute pack, and extend the forward bay cover parachutes (FBCPs) and drogues until they separated. Under the influence of drogues, the reentry capsule decelerated and maintained a steady attitude. After the drogues had been deployed for a certain period, they detached from the reentry capsule, and the main parachute pack was deployed. The main parachute was initially filled in a closed state. After a brief delay, the closure was released, allowing the main parachute to be fully inflated. The reentry capsule safely descended to the ground at a controlled speed using the main parachute.

In this section, we validate the theoretical model using an airdrop experiment. An airdrop test was conducted using a full-scale replica of a reentry capsule. By launching it at an altitude of 5124.50 m, we obtained the time-dependent changes in the altitude and velocity of the reentry capsule [11]. We employed the dynamic model established in Section 2.3 to simulate and calculate the outer ballistic data of characteristic points within the parachute cabin system under identical initial conditions. These results were then compared with the data from the airdrop test.

Figure 3 shows a comparison between the simulated and measured values of altitude and velocity during the airdrop test. The model-based simulation results closely matched the airdrop test results, confirming the correctness and effectiveness of the established 9-DOF dynamic model.

3. Model Parameter Sensitivity Analysis

Parameter sensitivity analysis is a valuable method for evaluating the effects of different parameters on the output. This enabled us to address the uncertainty in the output of each variable in the entire system [12]. The Sobol method [13] involves the Monte Carlo sampling of parameters within a feasible space, simulation of the sampled values, and generation of a large number of output results. During the sampling process, each parameter is associated with a distribution variance, and the ratio of the variance of a parameter to the variance of the output result represents the first-order sensitivity. Assume that the model/system is denoted by , where denotes the output and denotes a set of influence parameters. The first-order sensitivity analysis coefficient reflects the independent effect of a single parameter on the variance of the model output. The first-order sensitivity of is defined as where indicates the set of all parameters except , and denote expectation and variance operators, respectively. represents the expectation operation performed on parameters , given the condition of . It yields the average expected value of when is known.

Moreover, the global sensitivity of parameter is defined as follows:

The difference between the global and first-order sensitivity coefficients indicates the extent to which higher-order parameter interactions affect the model. The closer this difference is to zero, the less influence the parameter exerts on the output.

In practice, , , and are estimated using the following sampling procedures. (1)Generate an sample matrix, where each row represents a sample point of 2D parameters. The value for each parameter is generated based on its probability distribution. The -th and -th columns correspond to the -th parameter such that they are generated based on the probability distribution of the -th parameter(2)The first column of the sample matrix is used as matrix , and the remaining columns are used as matrix (3)Construct matrices , where the -th column of is equal to the -th column of , and the remaining columns are from (4)The calculated outputs of the samples in the matrices , , and are denoted by , , and , respectively.(5), , and are estimated as follows:

In our system, 11 uncertain parameters affect the recovery of the reachable tube of the reentry capsule. These parameters included the mass, drag area, atmospheric density, wind direction, wind speed, wind direction angle, initial velocity, reentry inclination angle, reentry deviation angle, height, and longitude and latitude deviations. The longitude and latitude deviations affect only the translation of the recovery reachable tube. Instead of considering them as uncertain parameters, we analyzed the offset of the parachute opening point. We treated 10 of these uncertainties as inputs, and the distance between the final and unbiased landing points (The unbiased landing point denotes the landing point without considering any uncertain parameters.) was considered as the output. In this paper, we discuss how these uncertain parameters can influence the boundary prediction of the reachable tube of a reentry capsule.

The detailed parameter sampling settings for the sensitivity analysis used in our method are listed in Table 1, where denotes the base value for the -th parameter.

The results depicted in Figure 4 highlight that certain parameters, specifically wind direction, wind speed, atmospheric density, altitude of parachute deployment, parachute-drag area, and eastward offset, exert a significantly greater influence than the others. Consequently, only these six uncertain parameters were considered as optimization variables in the subsequent analysis. This strategic approach not only reduces the complexity arising from uncertain parameters but also substantially enhances the optimization efficiency.

4. Reentry Capsules Reachable Tube Boundary Prediction Model and Algorithm

The concept of reentry capsule reachable tubes is based on accessibility theory. Therefore, this section briefly introduces the accessibility theory for nonlinear dynamic systems, followed by multiobjective optimization models for predicting the reentry capsule reachable tube boundary. Finally, a multiobjective evolutionary algorithm based on decomposition (MOEA/D) is introduced to solve the formulated multiobjective optimization problem.

4.1. Accessibility Theory

Accessibility analysis is commonly employed in automatic control and computer science [14]. It is primarily used to address a class of complex systems in which continuous dynamic and discrete events coexist and interact, essentially addressing complex nonlinear systems [15]. Reachability analysis generally encompasses forward and backward reachable sets. In this context, we focus on discussing an approximate method for calculating the forward reachable set. Throughout the following discussion, the term “reachable set” specifically refers to the forward reachable set. The fundamental approach involves employing a constructed model and iteratively advancing through each segment based on a given initial set and specified time step. The continuous dynamics of each period were approximated, and the set of reachable states within a specified time was computed. A pictorial illustration of the reachable set is shown in Figure 5.

For the following nonlinear dynamic system where denotes the time, represents the input of the power system, denotes the input uncertain parameters, is the current height, and is the system output (e.g., represents the position coordinates in our system). The reachable set of the initial set is the set of all end states that can be reached along the system trajectory from the initial states in that are within the allowable control input [9]. The reachable tube is formally defined as follows: where denotes the reachable set at time and denotes the number of status variables.

Based on the concept of a reachable tube, we can similarly define a reentry capsule reachable tube that comprises all possible positions of the reentry capsule at any allowable height. In the following section, we formulate the prediction of the reentry capsule reachable tube boundary as a multiobjective optimization problem and then use a multiobjective evolutionary algorithm to address it.

4.2. Multiobjective Optimization Model of Reentry Capsule Reachable Tube Boundary Prediction

We assume that the position of the reentry capsule at height when no uncertain parameters are considered (i.e., all uncertain parameters are set to 0). The boundary of the reachable tube, denoted as , at a height of , consists of the farthest positions from in all directions at height . The boundary of the reachable tube at all heights, given by , forms the complete boundary of the reentry capsule reachable tube during the recovery process, where denotes the height of the reentry capsule at the recovery start-up point. Therefore, as shown in Figure 6, the boundary point at the height along the direction can be obtained using the following single-objective optimization model: where denotes the decision vector, including the six uncertain parameters (namely, wind direction, wind speed, atmospheric density, altitude of parachute deployment, parachute-drag area, and eastward offset) and denotes their possible values at height .

Because of the computationally expensive nature of the system output calculation, this single-objective model has the following drawbacks. (1)For a specific height , it needs to provide a direction to determine one boundary point, and the total number of directions is infinite(2)The optimal solution only corresponds to one boundary point at the height along direction . Hence, the optimization algorithm must be run multiple times to approximate the complete boundary at the height

To address these issues of the single-objective optimization model, we developed the following multiobjective optimization model:

Assume the position of the equilibrium point without any disturbance at height as the origin and establish the coordinate system at height . In the established coordinate system, the model (16) is a multimodal multiobjective optimization problem [16, 17], which means that different decision vector values can result in the same objective function value. The complex nature of the problem poses challenges for the optimization algorithm in determining the entire Pareto front (PF) that can cover the entire boundary. Therefore, we address this issue by separating it (16) into four independent multiobjective optimization subproblems as follows: (i)The MOP in quadrant 1 is as follows:(ii)The MOP in quadrant 2 is as follows:(iii)The MOP in quadrant 3 is as follows:(iv)The MOP in quadrant 4 is as follows:

Based on the separation, we have four MOP subproblems in the coordinate system . However, since we remove the direction of the single-objective optimization model in the multiobjective optimization model, the PFs of all these four MOPs may not effectively cover the entire boundary. Therefore, we rotate the coordinate system clockwise by 45° and get another four MOP subproblems in the rotated coordinate system. Finally, we get eight MOP subproblems at the height , and the combination of the PFs of all these eight MOP subproblems forms the final reachable boundary at height . A pictorial illustration of the PFs in the original coordinate system and the corresponding rotated coordinate system are shown in Figure 7.

For the height , we can obtain its reachable boundary by solving the eight MOP subproblems. The reachable tube boundary of the reentry capsule during the entire recovery process was obtained by combining all reachable boundaries within the height range of . In this study, we discretized the height using and used a combination of the boundaries at all discrete heights to approximate the entire boundary of the reachable tube. A pictorial illustration of the discretization process is shown in Figure 8.

4.3. Multiobjective Evolutionary Algorithm Based on Decomposition

Multiobjective evolutionary algorithms (MOEAs) have been a popular research topic in various fields [18, 19] for several reasons: (1) They can handle a wide range of multiobjective optimization problems. (2) They can generate a set of Pareto optimal solutions in a single run. Current mainstream MOEAs can be divided into three categories: dominance-[20], decomposition- [21], and indicator-based MOEAs [22]. In this study, because complex nonlinear relationships exist between the decision vectors and objectives in the formulated multiobjective optimization models, we adopted a specific MOEA, namely, MOEA/D [21, 23], to solve all the MOP subproblems. MOEA/D first decomposes an MOP into a number of single-objective optimization subproblems using a set of uniformly distributed weight vectors and then collaboratively solves all these single-objective subproblems. The two main reasons for using MOEA/D in this paper are as follows [24, 25]: (1)The PFs of all MOP subproblems were continuous and regular(2)The two objectives of all MOP subproblems have the same magnitude

The details of the MOEA/D algorithm can be found in the literature [21].

5. Simulation Experiments

5.1. Experimental Settings

In our simulation experiments, the search space for the uncertain parameters remained the same as that listed in Table 1, and five discrete heights were considered (, , , , and ). For the MOEA/D, the Tchebycheff method was selected as the decomposition method, and the population size and maximum number of iterations were set as 200 and 100, respectively. The other parameters of MOEA/D were the same as those used in the original study [21]. Each MOP subproblem at a fixed height was independently solved using MOEA/D.

5.2. Experimental Results

In this study, we utilized the commonly used Monte Carlo target shooting [26] and our proposed methods to predict the reachable boundary of the reentry capsule at five discrete heights. The experimental results are shown in Figure 9. The following observations were derived from the results: (1)The Monte Carlo target-shooting method cannot provide a completely reachable boundary at each height, whereas our method can(2)The boundary determined by our method fully encompasses the region predicted by the Monte Carlo target shooting method, indicating that our method has higher accuracy in predicting the reachable boundary than the Monte Carlo target shooting method

Finally, by applying the interpolation method to the results obtained at these five heights, the reachable boundaries at any height can be determined. The resulting reentry capsule-reachable tube obtained using the proposed method is illustrated in Figure 10. Figure 10 displays the primary structure of the reachable tube, and its accuracy can be readily enhanced by incorporating data from a greater number of heights.

6. Conclusion

This paper introduces a novel approach by modeling the reentry capsule reachable tube boundary prediction problem as a multiobjective optimization problem and solving it by employing a decomposition-based multiobjective evolutionary algorithm. Its performance was validated through comparison with the commonly used Monte Carlo target shooting method in the context of estimating the reachable tube of the Chang’e-5 reentry capsule. The experimental results demonstrate the effectiveness of both the model and optimization algorithm.

In terms of algorithms, the proposed method involves solving numerous multiobjective optimization problems to determine a reachable tube. There were positive correlations between these multiobjective optimization problems. Therefore, we plan to employ evolutionary multiobjective multitask optimization algorithms in the future to enhance both the solution accuracy and efficiency [2729]. In addition, because of the high computational requirements of the system, we intend to consider data-driven evolutionary algorithms [30, 31] to deal with it.

Our research provides a new way for the recovery prediction of other types of parachute-capsule systems, which can be applied to the recovery of general spacecraft in activities such as deep space exploration, improve the efficiency and accuracy of search and rescue, and avoid economic losses.

Data Availability

No underlying data was collected or produced in this study.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was supported by the National Natural Science Foundation of China under grants 62206120 and 11772353.