Trajectory optimization of a reentry vehicle based on artificial emotion memory optimization

FU Shengnan, WANG Liang, and XIA Qunli

1. School of Mechatronical Engineering, Beijing Institute of Technology, Beijing 100081, China;

2. Beijing Aerospace Automatic Control Institute, Beijing 100854, China;

3. School of Aerospace Engineering, Beijing Institute of Technology, Beijing 100081, China

Abstract: The trajectory optimization of an unpowered reentry vehicle via artificial emotion memory optimization (AEMO) is discussed. Firstly, reentry dynamics are established based on multiple constraints and parameterized control variables with finite dimensions are designed. If the constraint is not satisfied, a distance measure and an adaptive penalty function are used to address this scenario. Secondly, AEMO is introduced to solve the trajectory optimization problem. Based on the theories of biology and cognition, the trial solutions based on emotional memory are established. Three search strategies are designed for realizing the random search of trial solutions and for avoiding becoming trapped in a local minimum. The states of the trial solutions are determined according to the rules of memory enhancement and forgetting. As the iterations proceed, the trial solutions with poor quality will gradually be forgotten. Therefore,the number of trial solutions is decreased, and the convergence of the algorithm is accelerated. Finally, a numerical simulation is conducted, and the results demonstrate that the path and terminal constraints are satisfied and the method can realize satisfactory performance.

Keywords: trajectory optimization, adaptive penalty function,artificial emotion memory optimization (AEMO), multiple constraint.

1. Introduction

The reentry trajectory optimization problem [1-5] is an optimal control problem, which is limited by path and terminal constraints and aiming at maximizing the performance index. Thus, an optimization method should have a superior fast convergence and a global optimization performance.

In recent years, the intelligent optimization method has been extensively studied as a type of direct method. Janin et al. [6] applied the genetic algorithm (GA) to test the problems and space flight problems. Miyamoto et al. [7]proposed a dynamic distributed GA for the trajectory generation to maintain the diversity for its optimization solutions. Naresh et al. [8] used GA for generating a trajectory for the maximum range under various in-flight and impact constraints. Lu et al. [9] proposed the improved particle swarm optimization (IPSO) with dynamic parameters, which outperforms the traditional particle swarm optimization (PSO) in global search and converges faster. For the trajectory optimization of Boeing aircraft in vertical flight profile, Zhu et al. [10] proposed an IPSO algorithm with object-oriented performance database. Zheng et al. [11] proposed an innovative parameter-adaptive strategy for ant colony optimization (ACO)algorithms for tracking prespecified paths. Kiyak et al.[12] applied ACO to a single-stage solid missile design problem involving six degrees of freedom (6-DOF) flight trajectory modeling.

These intelligent optimization methods provide new approaches for solving large-scale nonlinear problems,but there are still some shortcomings in some ways[13-16]. If the number of the population is too large, the optimization algorithm may take a long time to reach the convergence. In addition, the autonomy of the optimization algorithm needs to be strengthened. For ACO, the foraging behavior of ants is acted with blindness. For PSO, the autonomy of the behavior for particles is relatively weak. For the global optimization, the intelligent optimization methods are easy to fall into the local extreme point when dealing with multi-peak problems.

Some studies have introduced the memory theory into the optimization algorithm. Duan et al. [17] used an extended memory to store the historical information of particles for PSO. Huang et al. [18] proposed a novel forecasting method inspired by memory. However, the memory theory in these studies is quite different from the real human memory. By introducing the theory of cognitive psychology [19], Huang et al. [20,21] proposed the artificial memory optimization (AMO). Sun et al. [22] introduced the emotion theory into AMO, and proposed artificial emotion memory optimization (AEMO). By introducing different types of emotional memory cells [23-25],the short-term memory cell will be forgotten, and the long-term memory will have a stronger impact on the optimization process. Therefore, the convergence process of the optimization algorithm can be accelerated.

For the problem that optimization time of the previous intelligent optimization methods is too long, AEMO introduces emotional memory cells and enables poor quality memory cells gradually forgotten to accelerate the convergence velocity. For the optimization methods may not be so autonomous, AEMO makes the optimization direction close to the solution with a better quality, by introducing long-term memory cells. For the problem that the optimization methods are easy to fall into the local extreme point, AEMO can retain or introduce a few cells with poor quality to keep the diversity of the population when approaching the optimal solution. Thus, the global optimization can still be carried out.

In this paper, AEMO is improved and used to solve the trajectory optimization of a reentry vehicle. Compared with the method described in [22], there are two improvements in this paper: (i) A distance measure and adaptive normalized penalty function are introduced to redefine the objective function. The dynamic adjustment strategy is designed for memory thresholds and forgetting thresholds.This design enables memory cells to be uniformly distributed in instantaneous, short-term and long-term memory systems, and controls the number of memory cells in the forgetting system. The adjustment strategy ensures the diversity of memory cells in the final stage of convergence.(ii) A reasonable and feasible trajectory optimization model of the reentry vehicle is established, and AEMO is applied to the trajectory solution.

The remainder of this paper is organized as follows:The reentry dynamics and constraints are described in Section 2. The optimization model is established in Section 3. The AEMO algorithm is presented in Section 4,and the numerical results are discussed in Section 5. Finally, the most important conclusions are summarized in Section 6.

2. Problem formulation

2.1 Reentry dynamics

This article assumes that the sideslip angle of the unpowered hypersonic reentry vehicle is maintained at0°throughout the flight. Based on a spherical and irrational earth model, 3-dimensional equations [26] of motion with respect to time can be established as follows in polar coordinates:

There are six state variables inX=(V,θ,φ,r,γ,ψ)T:Vis the velocity, θ is the longitude, φ is the latitude,ris the radial distance from the Earth’s center to the vehicle’s center, γ is the flight path angle and ψ is the heading angle. In addition, σ is the bank angle, andmis the mass of the vehicle. The aerodynamic lift accelerationLand the aerodynamic drag accelerationDare expressed as

whereq=ρV2/2 is the dynamic pressure,Srefis the reference area of the vehicle,CLis the aerodynamic lift coefficient,CDis the aerodynamic drag coefficient, α is the angle of attack,Ma=V/cis the Mach number, andcis the velocity of sound.

The atmospheric density ρ in an exponential form is expressed as follows:

where ρ0is the atmospheric density at sea level,h=r-r0is the altitude,r0is the equatorial radius of the Earth andhsis the atmospheric scale height.

The acceleration of gravitygis expressed as follows:

where μ is the gravitational constant.

2.2 Constraints

Reentry flight is a highly dynamic process with strict constraints, which ensures the reliability of the structure and thermal protection. The path constraints mainly include constraints on the heating rate, the dynamic pressure and the aerodynamic acceleration [27]:

wherekQis the thermal model coefficient,is the heating rate, andnis the aerodynamic acceleration.

Terminal constraints are affected by the flight mission,which include constraints on the altitude, velocity and position:

wherehfandVfare expected terminal values andtfis the terminal time.

3. Trajectory optimization design

3.1 Performance index

Optimization objectives can be designed according to the requirements. To realize the minimum flight range, the objective function can be described as the performance index as follows:

wherergois the flight range of the vehicle.

3.2 Reference profiles

For the simplified trajectory model, the lateral motion can be ignored and the control variable is only the angle of attack. Because the velocity of a reentry vehicle decreases monotonously, the control variable can be designed as follows:

wherenis the number of nodes, and [α0,α1,···,αn-1] corresponds to the monotonically decreasing [V0,V1,···,Vn-1].

Since the control nodes may cause the state variables change unsteadily, in order to make the control characteristics of the vehicle stable, the angle of attack profile can be described by the cubic spline interpolation function for each interpolation segment:

The optimal method needs to solve the angle of attack at each velocity node. After obtaining the trial solution,the data is used to generate the real-time angle of attack by the cubic spline interpolation. Thus, the final fitness of each memory cell can be calculated by integrating the dynamic equations.

3.3 Optimization model

A general optimization model can be described [28] as

where u is the control variable, S is the search field and Rnis ann-dimensional Euclidean space.f(X,u) is the objective function, andgi(u) is the inequality constraint,whereIis the number of inequality constraints.

Corresponding to (17), control variables u=α are determined via the optimization method. Meanwhile, the path constraints in (11)-(13) and the terminal constraints in (14) and (15) must be satisfied. In the search field, u is subject to (17) and (18).

4. AEMO

The process of memorization includes identification, recording and recall of information, and memory can be divided into three phases: instantaneous memory, shortterm memory and long-term memory. The instantaneous memory identifies and filters external information; the short-term memory is a buffer that is applied before the information enters into the long-term memory; and the long-term memory is the relatively stable result of longterm stimulus, which enables the algorithm to recall the past information. The relationships are illustrated in Fig. 1.

Fig. 2 Process of the algorithm

Regarding emotional intelligence, individuals generate emotions according to their own states and external stimuli. If a stimulus benefits the individual, it will generate positive emotions; otherwise, it will generate negative emotions. Emotional memory is the continuation of emotion in time and space. Behavior is decided by emotions, which are stored in the memory system.

4.1 Memory cells

The artificial emotional memory system is composed of a set of memory cells {M1,M2,···,MN}, and these cells consist of trial solutions, memory residual values (MRVs)and other features. Trial solutions {X1,X2,···,XN} are selected from the solution space randomly. Each trial solution corresponds to a memory cell, and this relationship remains constant throughout the optimization process.

Under the stimulation of external information, memory cellMi(i=1,2,···,N) is activated, and emotional memory can be generated, which can be classified into two types:

(i) Positive emotional memory: The quality of the trial solutionXi(i=1,2,···,N) has improved, namely, the value of the objective function is smaller than that prior to the stimulus. Thus, the memory cellMiwill be strengthened.al solutionXi(i=1,2,···,N) has decreased, namely, the

(ii) Negative emotional memory: The quality of the trivalue of the objective function is larger than that prior to the stimulus. Thus, the memory cellMiwill be weakened.

The memory system consists of an instantaneous memory system, a short-term memory system and a longterm memory system. The system to which a memory cellMibelongs is decided based on memory and forgetting thresholds. When receiving external information for the first time, the memory cell enters the instantaneous memory system. Once MRV reaches the memory thresholdMSof the short-term memory system, the memory cell enters the short-term memory system. Similarly, once MRV reaches the memory thresholdMLof the long-term memory system, the memory cell enters the long-term memory system. Three forgetting thresholds are designed for these memory systems. If MRV is below the corresponding forgetting threshold (FI,FSandFLfor the instantaneous, short-term and long-term memory systems,respectively), the memory cellMiwill be forgotten and will no longer be stimulated by the external information.Memory cells are constantly stimulated by the external information; hence, the satisfactory memory cells will gradually approach the globally optimal solution, and the poorly performing memory cells will be gradually forgotten. To store all the memory information, a memory cell is designed as follows:

whereXiis a trial solution,miis the MRV (if the emotional memory is positive, Δmi＞0 ; otherwise, Δmi=0),tiis the number of iterations,siis the memory state of the memory cell (si=1 if the memory cell is in the instantaneous memory system,si=2 if the memory cell is in the short-term memory system, andsi=3 if the memory cell is in the long-term memory system),fiis the forgetting state of the memory cell (fi=1 if the memory cell has been forgotten; otherwise,fi=0 ) , andkiis the behavior type (ki=1,ki=2 , andki=3 correspond to behaviors 1,2, and 3).

4.2 Penalty function

When a behavior is executed, new external information will be input into the optimization algorithm, and a new trial solution will be generated, which may not satisfy all constraints. In this paper, a distance measure and an adaptive penalty function are used to address this scenario.

Since the magnitudes of the objective function and constraints differ, it is necessary to normalize them. The normalized objective function is expressed as

wherefmaxandfminare the maximum and minimum values of the performance indices for all trial solutions. The degree of violation for all constraints is defined as

whereci(X,u)=max(0,gi(X,u)),i=1,2,···,Iandci,maxis the maximum value of theith constraint over all trial solutions.

The distance is defined as follows:

whererfis the ratio of the number of feasible solutions(whereXsatisfies all constraints) to the total number of trial solutions.

The adaptive penalty function is defined as follows:

where

and

The new objective function is defined as

Compared to the traditional penalty function, this method has the following features: (i) If the multiple trial solutions are feasible, the trial solution with the smallest value of the objective function and the lowest degree of violation for all constraints is dominant among all trial solutions. (ii) If all trial solutions are feasible, the trial solution with the smallest value of the objective function is dominant. (iii) The values of the objective function and constraints are normalized to avoid complex design of the penalty function and to increase generality of this method.

4.3 Memory model

The state of a memory cell is determined by MRV and the thresholds of the memory systems, as discussed in Subsection 4.1. To record the state of a memory cell, we denote the instantaneous memory state, the short-term memory state and the long-term memory state asI,S,andL, respectively, which are expressed as follows:

(i) Forgetting model

As the iterative generation increases, the memory will decay and the MRV will decrease continuously. The forgetting model consists of forgetting of memory and fading of recall according to the Ebbinghaus forgetting curve[29], which is defined as

wheremi(t) is the memory value of the memory for timet, Δtis the time increment, andaand λ are accommodation coefficients (a＞0 andaI,aSandaLcorrespond to the states of the instantaneous, short-term and long-term memories; λ is defined similarly).

Intuitively, the decay rate is directly proportional to MRV and inversely proportional to the time. Thus, the forgetting model must satisfy the following limits:b

(ii) Updating MRV

For a positive emotion, the increment of MRV Δmi(t)is proportional to the degree of the external stimulus.Considering the increment of the external stimulus and the magnitude of the stimulus, the update model can be defined as

wherehis an accommodation coefficient andh＞0.

The new MRV can be calculated according to the forgetting model and the update model:

Substituting (31) and (32) into (33) yields:

For simplicity, we design a forgetting velocity coefficient β=e-aΔt+λe-bΔt( 0＜β＜1). Thus, (34) can be reformulated as

In order to ensure the diversity of memory cells, the memory thresholds are dynamically adjusted. By adjusting the instantaneous memory threshold and the shortterm memory threshold, the long-term memory cells can also be adjusted. The adjustment strategy is designed as follows:

End for

The forgetting threshold also needs to be adjusted, and the strategy is designed as follows:

4.4 Behavior operator

At timet-1 , the trial solutions areRandomly selectL(L≥1) memory cells from the longthe long-term memory is smaller thanL, selectL(L≥1)term memory system. If the number of memory cells in memory cells for which MRV is smaller than that of the current memory cell from all the trial solutions, and form a new set of trial solutionsWe design the following three behaviors for updating

4.5 Behavior selection

4.6 Time complexity analysis

In Fig. 2, we assume that the maximum number for the iteration stopping condition is max_iter, the population size of memory cells isN, and the dimension of the memory cell isn. In each iteration, the time complexity of the algorithm corresponding to Fig.2 is as follows:

Step 1 determines whether the memory cell enters the forgetting system, and the time complexity isO(N).

Step 2 determines which behavior will be selected, and the time complexity isO(N).

Step 3 executes the behavior, and the time complexity isO(Nn).

Step 4 determines which kind of the emotional memory will be produced, and the time complexity isO(N).

Step 5 updates the memory cell, and the time complexity isO(N).

Step 6 updates the probability, and the time complexity isO(N).

Thus, the time complexity of the algorithm isO(Nnmax_iter) . Since max_iter,N, andnare finite constants, the algorithm has the fastest speed for the time complexity.

5. Simulation

The proposed algorithm is verified and compared with other optimization algorithms. Then the algorithm is applied to the trajectory optimization to verify its applicability.

5.1 Algorithm verification

In this section, the Griewank function is selected to test the performance of the algorithm. The geometric characteristics of the function is shown in Fig. 3, and described as follows:

Fig. 3 Geometric characteristics

The Griewank function is a typical multimodal function, which has multiple suboptimal solutions and a global optimal solution. It is used to test the optimization accuracy of the algorithm. The global best value is min(f(X*))=f(0,0,···,0)=0.

The Griewank function runs 30 times for the intelligent optimization algorithms of AEMO, differential evolution (DE), PSO and GA, and the average values are taken for statistical analysis.

In the optimization, the maximum number of iterations is 500, and the number of optimization variables is 20.For the AEMO algorithm, the number of memory cells is set as 200. The number of memory cells that are selected for the long-term memory system isL=3, the selection probability ise=0.3 , and the variation factor is φ=0.03.The forgetting velocity coefficients of the instantaneous,short-term and long-term memory systems are βI=0.99,βS=0.95 and βS=0.90, respectively. The memory thresholds of the short-term and long-term memory systems areMS=0.025 andML=0.15. The forgetting thresholds of the instantaneous, short-term and long-term memory systems areFS=0.005 ,FS=0.003 andFL=0.0005 . The accommodation coefficient ish=5.

For the DE algorithm, the weighting factor is 0.5, and the mutation probability is 0.5. For the PSO algorithm,the cognitive coefficient is 1.2, and the social coefficient is 1.8.

The iteration stopping rule is designed as follows:

(i) Reach the maximum number of iterations;

(ii) The minimum objective function value remains unchanged in 100 iterations;

(iii) The minimum objective function reaches the convergence threshold 1.0e-3.

Run the intelligent optimization algorithms of AEMO,DE, PSO and GA independently, and the optimization results are shown in Table 1. The similarity of the optimal solution refers to the proportion of the trial solutions similar to the optimal solution in the whole population when the optimal problem reaches the convergence at the end. It can be seen that AEMO has a great advantage in convergence velocity and accuracy compared with other algorithms and the diversity of the population is well maintained at the end compared with other optimization algorithms. To compare with the unimproved AEMO, the distribution of memory cells at the end is added in Table 2.It can be seen that this design enables memory cells to be uniformly distributed in instantaneous, short-term and long-term memory systems and controls the number of memory cells in the forgetting system. Thus, the diversity of the trial solutions can be maintained.

Table 1 Contrast results 1

Table 2 Contrast results 2

Fig. 4 shows the optimal solution convergence process of the optimization. It can be seen that the convergence velocity of the AEMO algorithm is significantly faster than other intelligence algorithms.

Fig. 4 Convergence process of the optimal solution

In order to verify the applicability of AEMO, other functions are selected. In this simulation, the maximum number of iterations is 30 000. Other parameters remain unchanged. The optimization results are shown in Table 3.

Table 3 Optimization results

For each optimization function, other intelligence algorithms are applied, and the convergence time cannot exceed the convergence time of AEMO in Table 3. Results are recorded in Table 4.

Table 4 Contrast results for different functions

The iteration stopping rules are designed as follows:

(i) Reach the maximum number of iterations;

(ii) The minimum objective function reaches the convergence threshold 1.0e-3;

(iii) Reach the convergence time for AEMO (If other optimization algorithms meet the iteration stopping rules(i) or (ii) within the convergence time of AEMO, the global optimal solution and the convergence time are recorded as solution (time)).

In order to evaluate the performance of the applied algorithms, the scoring criteria can be designed as follows:

(i) For the same convergence time, the smaller the error between the final global optimal solution and the theoretical global optimal solution is, the higher the score is;

(ii) For different convergence time, the smaller the convergence time is, the higher the score is.

Based on the results in Table 3 and Table 4, evaluation results of optimization algorithms can be concluded in Table 5. According to the average scoring results, the performance ranking of applied optimization algorithms can be obtained:

Table 5 Evaluation results of optimization algorithms

In conclusion, AEMO is suitable for different kinds of optimization situations, with strong advantages in the convergence accuracy and the convergence time.

5.2 Trajectory simulation

The structural and aerodynamic parameters correspond to the common aero vehicle (CAV) [30], which is a hypersonic vehicle that was designed by the Lockheed-Martin company in 1998. The reference area is 0.483 9 m2, and the mass is 907.2 kg. To verify the algorithm, Matlab is used to model the problem [31]. The control variables are obtained via AEMO and introduced into the dynamic equations, and a trajectory can be obtained via numerical integration with the classical Runge-Kutta method [32].The integration step size is set as 1 s. The initial state variables are designed asV0=5000 m/s, θ0=-1°,φ0=a ndh0=60km.0°,r0=R0+h0,γ0=0°, ψ0=0°, whereR0=6378.136 km

In the optimization, the maximum number of iterations is 300. The number of memory cells is set as large as 50 in this paper. The optimization process ends when the maximum number of iterations is reached or the change of the objective function value is less than 0.1% for more than 50 iterations. The initial control variables are generated randomly within the range that is specified in Subsection 3.2, but the corresponding trajectories may not satisfy all constraints. Therefore, before the optimization process, we design a certain number of random trial solutions. In this paper, we design the number as 500. For the trial solutions satisfying all constraints, we select 50 trial solutions as the initial data.

Other parameters for AEMO must also be designed.The number of memory cells that are selected for the longterm memory system isL=3. The forgetting velocity coefficients of the instantaneous, short-term and longterm memory systems are 0.99, 0.95 and 0.90, respectively. The memory thresholds of the short-term and longterm memory systems are 1 and 6, respectively. The forgetting thresholds of the instantaneous, short-term and long-term memory systems are 0.01, 0.005 and 0.001, respectively.

The overall constraints are listed in Table 6.

Table 6 Constraints of trajectory optimization

The heating rate, the dynamic pressure and the aerodynamic acceleration that are obtained via AEMO can satisfy the path constraints strictly, and the errors of the terminal velocity and altitude are within thresholds.

In order to make the simulation results have less occasionality, optimization methods are simulated for 30 times, and average values are taken for analysis. The traditional DE and PSO algorithms are compared with AEMO to analyze the performance of the proposed method in reentry trajectory optimization. The parameters refer to Subsection 5.1.

The simulation results of different optimization algorithms for the trajectory optimization are demonstrated in Table 7 and Fig. 5. If the two trial solutions are the same, they are considered to be similar. Table 7 shows that the similarity of the trial solutions and the similarity of the optimal solution of AEMO are small when the optimization stops, which means that the population still keeps diversity at the end. Thus, the optimization process does not converge to the local extremum value. The uptime and number of iterations also show the advantages of AEMO. Since the maximum number of iterations may not be reached, the fitness of 50 iterations is shown in Fig. 5. Table 7 and Fig. 5 show that the premature convergence of PSO easily results in early ending of the iteration. The convergence speed of the DE algorithm is fast in the initial stage, but slower in the final stage. Compared with other algorithms, the result shows that AEMO has obvious advantages in solving trajectory optimization problems.

Table 7 Contrast for trajectory optimization

Fig. 5 Convergence process of the optimal solution for trajectory optimization

Select the typical optimization results to get the comparison of trajectory curves. In the initial phase of the reentry trajectory, the path constraint depends on the heating rate constraint due to the high velocity. In the middle and terminal phases, with the decrease of the altitude and the increase of the atmospheric density, the path constraint begins to depend on the dynamic pressure and the aerodynamic acceleration constraints gradually in Fig. 6. Due to the optimization requirement of the minimum range, the terminal trajectory oscillates greatly. The trajectories obtained by the three methods are similar to each other to a certain extent. Fig. 7 shows the curve of the flight path angle.

Fig. 6 Optimal reentry trajectory

Fig. 7 Curve of the flight path angle

Fig. 8 shows that the velocity decreases monotonically and satisfies the terminal constraint.

Fig. 8 Curve of the velocity

Fig. 9 presents the changes of the attack angle. The overall variation of control variables obtained by different methods are similar, especially the profiles obtained by the AEMO and DE. The cubic spline interpolation makes the control variable change smoothly.

Fig. 9 Profile of the attack angle

Fig. 10-Fig. 12 show that the heating rate, the dynamic pressure and the aerodynamic acceleration satisfy the path constraints. Although the objective function is to obtain the minimum range trajectory, the aerodynamic heat,the dynamic pressure and the aerodynamic acceleration still do not reach the maximum values, which is due to the constraints of terminal height and velocity. The satisfactory performance of AEMO in trajectory optimization is demonstrated in this section.

Fig. 10 Curve of the heating rate

Fig. 12 Curve of the aerodynamic acceleration

6. Conclusions

In this paper, the application of AEMO in trajectory optimization is studied. An adaptive penalty function is de-signed for handling complex constraints and enhancing the efficiency and versatility. Based on the parametric design of the control variables, the emotional memory model from biology and cognition is introduced. By applying rules of memory to the optimization of trial solutions, invalid solutions can be effectively removed to accelerate the convergence of AEMO. The simulation results for the reentry of a hypersonic vehicle demonstrate that the proposed method can generate a minimum range 3-DOF optimal trajectory under path and terminal constraints. This method does not require the estimation of initial values of the control variables according to engineering experience, and it enhances the operability of trajectory optimization.