Abstract
Particle swarm optimization (PSO) is a population-based stochastic optimization technique in a smooth search space. However, in a category of trajectory optimization problem with arbitrary final time and multiple control variables, the smoothness of variables cannot be satisfied since the linear interpolation is widely used. In the paper, a novel Legendre cooperative PSO (LCPSO) is proposed by introducing Legendre orthogonal polynomials instead of the linear interpolation. An additional control variable is introduced to transcribe the original optimal problem with arbitrary final time to the fixed one. Then, a practical fast one-dimensional interval search algorithm is designed to optimize the additional control variable. Furthermore, to improve the convergence and prevent explosion of the LCPSO, a theorem on how to determine the boundaries of the coefficient of polynomials is given and proven. Finally, in the numeral simulations, compared with the ordinary PSO and other typical intelligent optimization algorithms GA and DE, the proposed LCPSO has traits of lower dimension, faster speed of convergence, and higher accuracy, while providing smoother control variables.
1. Introduction
Swarm intelligence is a collective dynamic behavior of distributed, self-organized systems, natural or artificial, employed in work on artificial intelligence. It introduces many simple agents with very general rules to achieve an “intelligent” global optimal behavior. Swarm intelligence-based techniques can be used in a number of applications on optimization. The US military is investigating the swarming techniques to control unmanned vehicles. The European Space Agency is thinking about an orbital swarm for self-assembly and interferometry. NASA is investigating the use of swarm technology for planetary mapping.
In particular, trajectory optimization problem is one of the most important tasks in the preliminary design of the next generation of high speed vehicles, such as NASA’s X-43 unmanned hypersonic vehicle (HV), and has a great effect on the choice of conceptual design [1]. It is a daunting work to get the solution of nonlinear optimal control problem with the arbitrary final time and multiconstraints.
Trajectory optimization with multiconstraints has been linked with some stochastic search algorithms. The typical approaches include genetic algorithms (GA) [2], differential evolution (DE) [3], and particle swarm optimization (PSO) [4–6]. Besides, several novel bionic optimization algorithms spring up in these years and show remarkable efficiency in the industrial domain, such as the honeybee mating optimization [7], harmony search algorithm [8], and ants swarm optimization [9]; however, they are rarely used in the aerospace field because of their excessive novelty. The PSO algorithm, as an optimization algorithm based on swarm intelligence, similar to GA and DE for their origination from population-based heuristic search mechanisms, has recently become more popular due to its simplicity and effectiveness. It is widely used because of its simple principle, small number of parameters to be adjusted, and easy realization. Considering that the dynamic parameter optimization needs to be transformed into the static format, it is validated to evaluate the trajectory parameters optimization using a spline function [10].
Under the guidance of thought of the hybrid algorithm [11], a novel method of trajectory optimization is proposed to improve the global convergence ability in this paper, the Legendre cooperative PSO (LCPSO) method, which is a kind of cooperative PSO based on Legendre orthogonal polynomials.
In LCPSO, the following improvements make it promising in solving complex problems. Firstly, the search space is divided into certain subspaces and different swarms are arranged to optimize the different parts of the search space. In this way, the optimized scale for each swarm can be reduced directly. The algorithm is suitable for the problems with larger scale and higher dimension.
Secondly, for the optimal control problem by PSO, the discretization is necessary before solving the parameter optimization problem. Herein, the Legendre orthogonal polynomial approximation can achieve higher smoothness of control variables with lower dimensions. Additionally, a theorem on how to find the precise range of optimal parameters is given and proven, as well as how to find the boundaries of the coefficient of polynomials. Therefore, the proposed LCPSO is expected to realize the optimization with higher accuracy and efficiency.
Thirdly, in order to solve the optimal control problems with arbitrary final time rather than the fixed one, the proposed Legendre orthogonal polynomial approximation method introduces an additional control variable to transcribe the original optimal problem to the one with fixed final time. Then, a traditional one-dimensional search method based on the interval analysis is proposed to optimize the additional control variable. This way, the specific optimal problem with single boundary can be solved.
Finally, we use two typical trajectory optimization problems to illustrate efficiency of the proposed LCPSO algorithm. One is the ascent trajectory optimization of X-43 hypersonic vehicle, and the other one is the classic optimal orbit transfer problem. The simulation results demonstrate the advantages of the LCPSO in terms of solution accuracy and convergence rate by comparing with some traditional intelligent optimization algorithms.
The organization of the remainder of this paper is as follows. Section 2 formulates the optimal control problem by 3-DOF mass point dynamics. Section 3 describes the algorithm of novel Legendre cooperative PSO and some properties. Section 4 presents the numerical simulations on the performance of the proposed optimal algorithms. Section 5 draws conclusions to the paper.
2. Problem Description
HV is one of the main workhorses for most of the nations of the world for scientific studies and military and commercial applications [1]. Efforts are ongoing in the 21st century to enhance flexibility and reliability and reduce the overall cost of such systems. Here, without loss of generality, we provide the proposed algorithm to solve a trajectory optimal problem of HV. Its 3-DOF dynamics over a spherical earth are described by the following motion equations [12, 13]:
The 3-DOF dynamics described in (1) are dimensionless equations of motion with six state variables and two control variables. The real variables are normalized as follows:where is the radial distance from the earth’s center to the vehicle, is the longitude, is the latitude, is the velocity of the vehicle, is the flight path angle, is the velocity azimuth angle measured clockwise from the north, the control variables are the angle of attack and the fuel throttle opening , is the mass of the vehicle, and the terms and are the aerodynamic drag and lift acceleration, which are defined byHere, is the reference area of the vehicle. The terms of and are the coefficients of drag and lift, which are also functions of and .
The control variables and are approximated by Legendre orthogonal polynomial functions, respectively. Assuming is the Legendre polynomial of degree , the continuous function of variable can be approximated as
This way, the dynamic parameter optimization problem is transformed into a static one before the PSO algorithm performed.
The evaluation of the performance index starts when the scramjet is launched. The vehicle releases from the boost phase and then enters into the later ascent and cruise phase. Here, we define the switching time as . By analyzing the dynamic characteristics of HV, the fuel throttle opening shall be minimized during the later ascent stage.
In this problem, the optimal time of ascent stage is free but satisfying some constraints. Thus, an additional control variable is introduced to transcribe the proposed problem to a fixed final time problem. The switching step is available and then the ascent trajectory can be divided into two phases. In the first phase, only the fuel throttle opening needs to be optimized to guarantee that it is a minimized constant value at the beginning of the second phase. In order to simplify the equations and problem, we give a reasonable assumption that the switching step satisfies , , with a fixed time interval , and the system’s continuous state does not jump at each switching point. Therefore, the arbitrary terminal time optimal control problem is turned into a two-point boundary value problem with the fixed terminal time. The optimal time would be obtained as long as all the terminal conditions are fulfilled.
The design space of HV ascent trajectory optimization using LCPSO can be defined by the following vector:
Here, the left two parts in the vector represent the two groups of coefficients of the orthogonal polynomials which approximate the control variables and , while the last parameter represents the ratio between the switching step and the supposed final time. Hence, the problem space can be divided into three subspaces , , and . Then, in LCPSO, we provide two groups of particles to optimize the subspaces and , respectively, and use a one-dimensional search method based on interval analysis to optimize the variable .
To seek the solutions of both the continuous input and the switching step according to the given initial state , we define a performance function as follows:In this case, the above nonlinear optimal problem can be described as a basic optimal control problem of nonlinear Bolza type.
Considering that the performance function here is selected to minimize the fuel expendable ratio (FER) through the control variables and , normally, we use the following relationship to describe it directly:
Furthermore, the ascent trajectory terminate conditions can be represented by the following inequalities:where , , and are three state variables of the terminal point, while the tolerance is given by the upper bound and the lower bounds , , and .
The terminal conditions are added to the performance function to ensure that all the constraints are satisfied. Then, the performance function of LCPSO is defined as follows:Here, in (9), the first part represents the optimal index of minimum FER described in (7). The three terms , , and represent the terminal constraints index. is the radial distance from the center of the earth to the vehicle and represents the velocity at the final time, while the flight path angle index is .
Considering that the proposed LCPSO provides two groups of particles based on Legendre orthogonal polynomial approximation, the assumption of the fixed time interval can be transformed into a closed interval of Legendre orthogonal polynomials to facilitate subsequent processing. Meanwhile, the additional switching step can make the optimal interval of more accurate. After that, the PSO algorithm is employed to solve the optimal control problem about in the fixed interval and the next interval , . Furthermore, the 4th-order Runge-Kutta numeric integration is employed to evaluate the performance function, while the optimal time can be obtained as long as all the terminal conditions are satisfied.
3. Algorithm Formulation
3.1. LCPSO
The proposed LCPSO is designed as a double iteration scheme: the inner iteration of dual cooperative PSO (CPSO) [14] and the outer iteration of interval analysis. The mathematical formulation of the CPSO algorithm is redesigned in (10) and (11).with , , , ,where is the coordinate of particle in dimension at time as subinterval at interval iteration step .
The pseudocode of CPSO algorithm is reported in Pseudocode 1.
| ||||||||||||||||||||||||||||||||||||||||||||||||||
The mathematical formulation of the one-dimensional search method based on interval analysis is designed as
The one-dimensional search algorithm based on the interval analysis is reported in Pseudocode 2.
The midpoint and the initial width of the switching step interval are defined aswhere is the ratio between the switching step and the assumption final time in the interval iteration and the subinterval . The switch step interval is divided into subintervals at each interval iteration step .
The method starts with the initial switch step interval, which is split into multiple subdivisions. Then, those subdivisions are either sent to the solution list which are considered later, or removed from further test list by certain cut-off condition. The above process is repeated by choosing a new switch step interval until no switch step could be considered or a global optimal point is found.
3.2. Theorem of Boundaries Selection of LCPSO
The Legendre orthogonal polynomials can be generated using Gram-Schmidt orthonormalization [15] in the interval with the function of the weight . represents a Legendre polynomial with degree . We provide the following transformation in (15) of the independent variable of 3-DOF dynamics when to guarantee that :This way, the continuous function of variable can be approximated as
Then, using Legendre polynomial approximation, the original optimal control problem can be transformed into a parameter optimization problem: to find the optimal coefficients so as to minimize the fitness in (9) and satisfy the terminal constraints in (8).
Obviously, the appropriate ranges of the coefficients can improve the convergence of the optimization procedure and prevent explosion. In our formal studies, a theorem to determine the coefficient of the Legendre orthogonal polynomial is proposed in [16]. However, the statement and proof of the theorem were incomplete and less rigorous. The complete expression of the theorem of how to determine the boundaries of the coefficients of the orthogonal polynomials is proposed and proven.
Theorem 1. Assume that represents an th degree Legendre polynomial. The continuous function of variable can be approximated as in (16). And the range of belongs to the interval with the variable presented in (15). Then, the ranges of optimal coefficients can be calculated byorwhere
Proof. According to the assumption, we have the following equations:According to the orthogonal relationship of the Legendre polynomials, we havein which is called Kronecker Delta and satisfies Substituting (23) and (24) into (22), we haveAssume that and are continuous in and satisfy ; then, we haveexcept the finite numbers of points.
Since the range of is and from (20) and (21), we haveThen, the following inequations should be achieved:Hence,According to (26), (27), and (29), we can obtainMoreover, according to (25) and (30), we have the effective range of :Therefore, the upper bound and lower bound of coefficients can be defined:with .
The closed form for the orthogonal polynomials is given in (19). Assume thatThen,Therefore,Substituting (35) into (32), we haveor When , importing and into (32), we haveThe proof is achieved.
Consider that the control value always works in a symmetric feasible interval. Assuming that and , to simplify the results in practical application, we can obtain
4. Simulation Results
4.1. Ascent Trajectory Optimization of the X-43 Vehicle
Firstly, we use the model of the X-43 vehicle [13] as the test case. Particularly, the vector of the initial states is , while the terminal states are . The tolerance of state variables in the terminal point is . Some other experiment parameters are as follows: the number of LCPSO iterations that occurred , the population size of a single particle swarm , the dimension of the particle , the number of swarm iterations of the cooperated particle , the acceleration coefficients , and the range of the inertia weight and . In addition, the weights of the performance function are , and the parameters of interval analysis iteration are .
In Table 1, during the optimal processing, five subintervals are introduced and the optimal costs vary with the switch steps in each iteration step of interval analysis, and the results converged at the global optimal point only through three steps of iteration. It takes about 580 seconds under the hardware environment of Intel Core i3 3220 CPU and DDR4G 1600 MHz memory and the software environment of Windows 7 operating system with VC++ 6.0 compiling environment.
The optimal flight path angle profile is shown in Figure 1, with the corresponding optimal control variables profile in Figures 2 and 3.
In Table 2, according to the ratios of the feasible results, the best and the mean results of 200 times simulations using interval analysis of LCPSO with different coefficients’ ranges are presented. Apparently, the proposed method with the coefficients’ ranges selection has the best performance, while all state variables meet the state terminal constraints with higher precision. Furthermore, the feasible results distribution with standard range and half-range selection is shown in Figure 4, while the feasible results distributions with standard, two, and five times range selection are displayed in Figure 5. From the comparison, we can find that if the optimal solutions existed in this reductive area, the reduction of the search range to its half will probably decrease the search complexity. However, if the performance index function corresponding to the reductive area is of bad quality, the optimized solution is likely to be lost. On the other hand, the multifold coefficient will provide the unwanted search space, which will increase the complexity and reduce efficiency of the optimization.
4.2. Earth-Mars Transfer Problem
The orbit transfer problem is another kind of hotspot issues on trajectory optimization research [17]. Here, we use the proposed LCPSO method to solve the Earth-Mars trajectory transfer problem and compare it with the traditional DE and GA method.
The Earth-Mars transfer orbit can be divided into three stages: geocentric escape section, heliocentric transition section, and areosynchronous capture section. Each segment has different constraints. In the geocentric escape section, the initial state constraints are and the terminal state constraint is . On the other hand, the constraints required in the third segment are the terminal radius and the terminal speeds and , which satisfy . And the initial radius is the radius of the sphere of influence for Mars, .
The LCPSO algorithm parameters are set as follows: the population size and the maximum number of iterations . Meanwhile, the parameters of DE and GA algorithm are set as follows: the population size , the maximum number of iterations , the cross probability , and the mutation probability . The optimization variables contain the thrust angle of the aircraft, the initial and terminal orbital angles, and the time of the whole three-segment orbits flight. Currently, the range of the thrust angle and the initial/terminal orbital angles of the vehicle is set in . The problem is to transfer the orbit within the shortest time to minimize the fuel consumption, so the fitness function can be described as , where , , and are the transfer periods of three sections, respectively. Figure 6 shows the best fitness process using LCPSO, DE, and GA algorithms. It is shown that the LCPSO has a faster speed of convergence and higher accuracy than DE and GA algorithms.
The final resulting candidate solutions by LCPSO, DE, and GA are shown in Figures 7–9. Because the intelligent algorithm obtains the feasible solution to meet the accuracy requirements, the state curves of any two algorithms are not necessarily completely consistent. But both of them satisfy the constraints and accuracy requirements. For instance, as shown in Figure 7, although the final radius of the geocentric orbit is about 145 by LCPSO, 145.1 by GA, and 143 by DE, the orbital radius and the velocity are all within the permitted constrained ranges. And in this case, the optimization result of the final geocentric orbital radius by LCPSO is closer to the constraint requirements. Hence, the LCPSO exhibits better efficiencies on the constraints requirements than the other two algorithms. Moreover, the results of the whole transfer time by different optimization algorithms are shown in Table 3. It can be found that, in the case of satisfying constraints, the optimal trajectory by LCPSO will spend the shortest flight time among the three algorithms, which has the best performance.
5. Conclusion
In this paper, a novel interval analysis based Legendre cooperative PSO algorithm is proposed and applied to solve the trajectory optimization problems. The Legendre orthogonal polynomial approximation is synthesized with the dual cooperative particle swarms to transfer the arbitrary final time optimal problem into a two-point boundary value problem with fixed terminal one. Then, a fast one-dimensional interval search method is provided in each selected interval to reduce the search space of the particles in the iterations. Furthermore, a theorem that determines the range of the parameters of the Legendre polynomial is investigated, and the problem can be solved to get closer to the global optimal solution.
Lastly, in numeral simulation, the results demonstrate that the LCPSO algorithm can solve the unsmooth trajectory optimization problem of X-43 effectively and obtain a smooth control variable. And the appropriate range of parameter values will significantly reduce the complexity of the optimization search. Moreover, in the solution to the orbit transfer problem, the comparisons with the existing GA and DE algorithms represent the notion that the proposed LCPSO method has better performance in the speed of convergence, final accuracy, and constraints satisfaction.
The flaw of the proposed LCPSO algorithms is that the parameters of the Legendre polynomial are sensitive in adjustment. Nevertheless, in the trajectory optimization problem of HV, the variety of parameters is a relatively gentle process because of the dynamic characters of HV, and the parameters are adjusted less frequently in iterations. Normally, the parameters will remain the same after they are adjusted once before. The optimization between parameter sensitivity and optimal solution will be studied to make the method more practical in engineering in our future work.
Nomenclature and Abbreviations
| CPSO: | Cooperative particle swarm optimization |
| DE: | Differential evolution |
| DOF: | Degree of freedom |
| FER: | Fuel expendable ratio |
| HV: | Hypersonic vehicle |
| GA: | Genetic algorithm |
| LCPSO: | Legendre cooperative particle swarm optimization |
| PSO: | Particle swarm optimization. |
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This work was supported in part by the National Natural Science Foundation of China (Grants nos. 61573161 and 61473124).