A Hybrid-Driven Optimization Framework for Fixed-Wing UAV Maneuvering Flight Planning

: Performing autonomous maneuvering ﬂight planning and optimization remains a challenge for unmanned aerial vehicles (UAVs), especially for ﬁxed-wing UAVs due to its high maneuverability and model complexity. A novel hybrid-driven ﬁxed-wing UAV maneuver optimization framework, inspired by apprenticeship learning and nonlinear programing approaches, is proposed in this paper. The work consists of two main aspects: (1) Identifying the model parameters for a certain ﬁxed-wing UAV based on the demonstrated ﬂight data performed by human pilot. Then, the features of the maneuvers can be described by the positional/attitude/compound key-frames. Eventually, each of the maneuvers can be decomposed into several motion primitives. (2) Formu-lating the maneuver planning issue into a minimum-time optimization problem, a novel nonlinear programming algorithm was developed, which was unnecessary to determine the exact time for the UAV to pass by the key-frames. The simulation results illustrate the effectiveness of the proposed framework in several scenarios, as both the preservation of geometric features and the minimization of maneuver times were ensured.


Introduction
Recent times have witnessed a wide range of applications for unmanned aerial vehicles (UAVs) including the commercial, military, and research fields [1][2][3][4][5]. Most of the autonomous UAV flight missions are limited to cruise on a predefined path with steady flight states. However, in some scenarios, such as dog fights and high-speed obstacle avoidance, UAVs are required to perform fast and agile maneuver flights. During maneuvers, drastic changes in position and attitude will hinder UAVs from maintaining trim conditions. Therefore, developing maneuvering flight techniques is of great importance, and the current research mainly focuses on the quadrotor UAV [6][7][8][9]. It must be noted that fixed-wing UAVs have longer endurance and a larger payload capacity compared to those of the quadrotors. Yet the maneuvers of fixed-wing UAVs have not been thoroughly studied due to the fact of its sophisticated movement features [10]. Meanwhile, it is harder to conduct real flight tests also.
The research on autonomous maneuverable flight is assumed to be hierarchical including several topics such as maneuver decision-making, maneuver planning, and tracking control [11]. The planner first decides the category of the maneuver, then generates a specific trajectory for the tracking control. Among them, maneuver planning is essential for the maneuverability of UAVs. In this article, we were mainly concerned with maneuver planning issues. When validating the complex maneuver representation approaches, appropriate maneuver generation algorithms are also considered.
Generally, state-of-the-art maneuver generation algorithms incorporate data-driven [12] and model-driven approaches. One category of maneuver generation algorithms that have been successfully utilized in UAV maneuver generation is learning from demonstration [13]. By collecting the flight data taught by experts with a high-level planning layer, imitating learning methods can be applied to extract the maneuver features [14]. The resulting algorithm is practically capable of reproducing and generalizing the learned motions to some extent. However, due to the limitations of viewing distance, communication delay, and external disturbances, such as wind gusts, it is hard for the pilot to perform their flight skills optimally [11]. Therefore, optimality is essential in the learning procedure and few of the data-driven algorithms in the recent literature take the optimization into consideration.
In the robotic literature, model-driven planning has been widely studied and the typical approaches include polynomial interpolation, Dubins curve [10,15,16]. For some special scenarios, it is essential to plan the position and attitude trajectories simultaneously. Therefore, many researchers have also made explorations in SE(3) space [17][18][19]. Most of the methods are based on the differential flatness principle which can greatly simplify the calculation complexity. However, for fixed-wing UAV flight maneuvers, it has to be stressed that those approaches did not take the complicated physical models with complex aerodynamic features into consideration. Therefore, the non-differential flatness property of fixed-wing UAV [20] hinders the application of the aforementioned methods. Solving the optimal control problems is also an effective maneuver planning method and has attracted growing attention. In particular, both the direct optimization method and indirect optimization methods have been applied in UAV flights in recent years [21][22][23]. Through taking the physical model into consideration, dynamical feasibility is ensured. Furthermore, by combing various cost items, optimal maneuvers can be solved for different scenarios [24]. However, for complex maneuvers, it is always hard to model the optimization problem, which has a great impact on the calculation efficiency and trajectory quality. In visual tracking [25,26], maneuvering is also required, vision-and-language navigation [27][28][29][30] provides another novel perspective. This method combines vision, language, and action which can turn relatively general natural-language instructions into robot agent actions. It is generally used in indoor complex environments and also has certain lightening significance for fixed-wing tracking of dynamic targets or maneuvering in complex outdoor environments.
For the representation of complex maneuvers, one recent algorithm that has been successfully proposed is the idea of key-frames [31], which selects several specific points in 3D space as the key-frames and various tasks can be accomplished. In [32], both the position and yaw angle are taken into account in the key-frame and the minimum snap trajectory is generated using polynomials. Furthermore, the author formulates the trajectory optimization problem as a quadratic programming issue and then enables the quadrotors to pass through multiple circular hoops quickly. Ref. [33] emphasizes the Bang-Bang characteristics of the minimum time trajectory and compares the existing timeoptimal approaches, and through analysis, it is concluded that the polynomial method can only obtain non-optimal trajectories. On this basis, Foehn proposes complementary constraints [34] and solved both the time allocation and time-optimal trajectory planning problem elaborately. Through comparison and verification, it can be concluded that the trajectory designed by the algorithm is faster than that of professional pilots. For fixed-wing, a model is required, which can be obtained by identifying methods [35]. In addition, it focuses on the minimum time problem through setting a series of waypoints and utilizing the flight corridors and B-spline curves. Using numerical optimization methods to solve non-convex problems and cleverly designed initial guesses, the optimal trajectory is obtained. As for the fixed-wing UAV maneuver, it is hard to extract the trajectory features merely from position information. It is necessary to comprehensively consider the position and attitude features and design the corresponding key-frames for various maneuvers. Motion primitives also have a satisfactory effect on the decomposition of complex maneuvers [36]. Mueller and D'Andrea [15] proposed a framework for the efficient calculation of motion primitives. McGill University has presented a series of representative works based on the idea of motion primitives and maneuver optimization [37][38][39]. Dynamic motion primitives (DMPs) is a typical primitive which has been successfully utilized in the design of car driving motion libraries [40]. Inspired by the former works, we studied the dual quaternion-based dynamic motion primitives (DQ-DMPs) [41], which have the ability to learn and generalize maneuvers in SE(3) space. Nevertheless, it is hard for the dynamic motion primitive algorithm to ensure the kinodynamic feasibility.
Inspired by the above discussion, we propose a data-model-driven framework, which is shown in Figure 1, for fixed-wing UAV maneuvering optimization. The framework is based on a global model identified by teaching data, and uses an optimization method based on positional, attitude, and compound key-frames. Through numerical optimization, the time-optimal maneuver is finally obtained. Through comparison, the actions generated by this algorithm take less time than professional pilots. In complex actions, different primitives can be flexibly concatenated, which simplifies the generation of complex actions. The main contribution of this paper is listed as follows: (1) We proposed a novel data-driven approach for model identification and key-frames extraction using the learning from demonstration principles. Then, complex maneuvers are decomposed into multiple motion primitives; (2) Based on the motion primitives, the optimal maneuver generation issue is formulated into a time-optimal problem considering key-frames which the UAV must pass by. The connection method of different primitives was also considered in this paper for practicability; (3) The proposed framework was verified thoroughly in simulation experiments, and it was possible to deduce that this framework is applicable for flight maneuvers in reality.

Figure 1.
Hybrid-driven optimal maneuvering flight planning framework. We start with pilot demonstration and maneuver data collection and then perform model identification and maneuver key-frame and motion primitive analysis. Finally, based on the global model and key-frames, the optimal primitives are generated, and the corresponding primitives are concatenated into a complete maneuver.
In the next section, we discuss the basic knowledge. Sections 3 and 4 introduce data collection and acrobatic maneuver optimization, respectively. The experimental results are introduced in Section 5, and the conclusion in Section 6.

Global Fixed-Wing Model
The aircraft model plays an important role in our algorithm framework. In order to obtain feasible maneuvers, we need to establish a relatively accurate rigid-body model. Due to the singular problem of Euler angle during the big maneuver, this paper adopted the quaternion model [42] to calculate: For translational kinematics equation, the p i = [x, y, z] ∈ R 3 is the position of the inertial coordinate system, and v b = [u, v, w] is the speed of the body coordinate system. The two state quantities are connected by a conversion matrix, which is composed of q, where q = [q 0 , q 1 , q 2 , q 3 ] T ∈ SO(3) is a unit quaternion given q = 1. The rotation matrix is: In rotational kinematics equation, is angular rate, we can write: In the translational dynamic equation, T = [F T , 0, 0], F T is the thrust, and where F A represents the term of aerodynamic force and F g represents gravity-related items, F x , F y , F z are the aerodynamic forces in the x-, y-, and z-axes, respectively. The last one is the rotational dynamics equation, where τ = [M x , M y , M z ] T is the aerodynamic moment, and J is the moment of inertia which is composed by: There are many unknown coefficients in this model as well as the margin of state quantities. These quantities have a great impact on the performance of a UAV and the completion of the maneuver. In Section 3, we identify these unknown quantities in a data-driven way.

Acrobatic Maneuver
Acrobatic maneuvers are generally summarized based on the pilot's actual flight experience and have strong practical significance. Acrobatic flights are a competitive sport or also a performance event. Therefore, remote control flight is also a way to achieve it. Figure 2 shows a typical fixed-wing maneuvering process. The International Federation of Aeronautics (FAI) is the main maker of acrobatics rules. FAI also provides a number of basic maneuvers for acrobatic events such as the Cuban eight [43]. This article mainly focuses on the realization of these basic maneuvers on UAVs. In contrast, manned pilots are restricted by physiological limits, remote control flight is safer and more flexible, but there are communication delays and the impact of visual distance. For autonomous drones, it removes the limitations of visual range and physiological limits and has the potential to achieve acrobatics.

Acrobatic Maneuver Data Collection
Maneuver attitude and other information can be obtained by collecting flight data. A maneuver has high requirements for data accuracy and frequency. With the improvement in sensor accuracy and miniaturization, a large amount of manual flight experience can be accurately recorded through data. UAV design, modeling, and flight testing can be realized through data collection. In actual flight, accelerometers, magnetometers, and other sensors can be used to record the aircraft's position, attitude, control inputs, and other data [43]. On the other hand, flight simulation technology is also an idea to solve the problem of flight maneuvers. We built a flight simulation system that can be used to collect flight data. As shown in Figure 3, the simulation system was built by the open-source flight control px4 and the flight simulation software X-Plane. The remote controller sends the control instructions to the flight controller. The control signal is processed by the internal program and then sent to the X-Plane simulator; this is a typical hardware-in-the-loop simulation system. For different maneuvers, we collected the position, attitude, and other information of the aircraft through expert teaching.

Model Parameter Identification
As mentioned in Section 1, the global fixed-wing model contains many unknown coefficients, mainly in aerodynamic forces and moments: where ρ is the air density, S re f , b re f , c re f are the reference wing area, reference wing span, and average aerodynamic chord length, respectively, C x , C y , C z , C l , C m , C n are the aerodynamic coefficient and moment coefficient, δ e , δ a , δ r are the deflection of the elevator, aileron and rudder, respectively. The calculation method of the angle of attack and the angle of sideslip are α = arctan(w/u) and According to [44], the aerodynamic coefficients could be expressed using the global aerodynamic model in Equations (7) and (8).
where C x 0 , C x α , C x α 2 , C y 0 , C y β , C y p , C y r , C y δa , C y δr , C z 0 , C z α , C z δe and C z α 2 represent the coefficients related to the aerodynamic force, and C l 0 , C l β , C l p , C l r , C l δa , C l δr , C m 0 , C m α , C m q , C m δe , C m α 2 , C m 0 , C m α , C m q , C m δe , C m α 2 , C n 0 , C n β , C n p , C n r , C n δa , C n δr are the coefficients related to the aerodynamic moments.
A least squares method was carried out to estimate the unknown coefficients in Equation (1). The optimization object was set as follows: subject to: It is worth noting that the optimization problem in Equation (9) could be solved utilizing the nonlinear least squares optimization methods such as the Levenberg-Marquardt algorithms. However, the divergence of the optimization problem or the convergence to the suboptimal solution might occur. Moreover, in order to diminish the side effects of over-fitting and the errors during the integration, we chose 10 s of flight data for the identification by trial and error. The results of system identification are listed in Section 5.

Optimal Acrobatic Maneuver Design and Generation
In this section, we propose a maneuver optimization algorithm based on teaching data and accurate models. We first introduce the two basic concepts of the algorithm in maneuver design: key-frames and motion primitives. Then different kinds of key-frames are listed and the entire optimization problem for maneuver is formulated. Finally, we conducted a detailed analysis of the solution of the proposed optimization problem.

Key-Frames and Motion Primitives
The concept of key-frames is often used in computer animation and simultaneous localization and mapping (SLAM) to represent frames that are decisive over a period of time. As mentioned in Section 1, this concept is also used to represent the necessary waypoints in motion planning. We introduced it into maneuvers.
As we can see in Figure 4, the same maneuver has the same typical characteristics of position and attitude changes. The key-frames are used to indicate the position and attitude that play a key role in a maneuver. The momentary state is a short-term key-frame, and the continuous state is a long-term key-frame.
Motion primitives explain the execution of complex motions based on action units. In maneuvers, simple maneuvers can be regarded as a single motion primitive without segmentation. For complex maneuvers, since it contains a variety of different pose change segments, it will be difficult to calculate if viewed as a whole, and different optimization goals should be considered for different segments, so it is necessary to divide it into multiple motion primitives. Motion primitives explain the execution of complex motions based on action In maneuvers, simple maneuvers can be regarded as a single motion primitive wi segmentation. For complex maneuvers, since it contains a variety of different pose ch segments, it will be difficult to calculate if viewed as a whole, and different optimiz goals should be considered for different segments, so it is necessary to divide i multiple motion primitives.

Maneuver Optimization
In maneuver or maneuver primitives, a drone is required to complete a ce position and attitude change in the shortest time; inspired by [24,[31][32][33][34][35], we formu this problem as a keyframe-based time optimization problem. A schematic diagram o problem is shown in Figure 5.

Maneuver Optimization
In maneuver or maneuver primitives, a drone is required to complete a certain position and attitude change in the shortest time; inspired by [24,[31][32][33][34][35], we formulated this problem as a keyframe-based time optimization problem. A schematic diagram of the problem is shown in Figure 5. Motion primitives explain the execution of complex motions based on action In maneuvers, simple maneuvers can be regarded as a single motion primitive w segmentation. For complex maneuvers, since it contains a variety of different pose c segments, it will be difficult to calculate if viewed as a whole, and different optimi goals should be considered for different segments, so it is necessary to divide multiple motion primitives.

Maneuver Optimization
In maneuver or maneuver primitives, a drone is required to complete a c position and attitude change in the shortest time; inspired by [24,[31][32][33][34][35], we form this problem as a keyframe-based time optimization problem. A schematic diagram problem is shown in Figure 5.  First, we set the state quantity to x dyn = [p, q, v, ω], input u = [δ e , δ a , δ r , F T ] T which needs to meet the constraints of the kinematics and dynamics model in Equation (1). The state quantity must first satisfy the start and end constraints and must be given or set free. The control input and some state quantities need to meet the upper and lower limits of UAV performance such as u min ≤ u k ≤ u max and ω min ≤ ω k ≤ ω max . We selected the direct multiple shooting method to solve the optimization problem. Suppose the total time is t N , discretize it into N segments, dt = t N /N, and the index is k. Meanwhile, the state quantity needs to be discretized. In order to reduce the calculation error, the 4th order Runge-Kutta was used for numerical integration.
where k 1 , k 2 , k 3 , k 4 are integral terms composed of x dyn,k , u k , and dt. Suppose the number of key-frames is set to M, indexed by i, and M-dimensional process state variables λ and M-dimensional process change variables µ are introduced to record the completion of key-frames, which meet the constraints: As is shown in Figure 5, process state variables λ saves the state of event completion. λ is 1 at the beginning and 0 at the end. µ is a process change variable, which can be used to record the occurrence of instantaneous events. Through the constraint that µ i multiplied by the state quantity part is equal to 0, it is required that when the event does not occur, the state quantity part is not 0 and µ i is 0; when the event occurs, the state quantity part is 0 and µ i is 1, then the corresponding λ i changes from 1 to 0 permanently. Under this framework, λ i can be used to represent the time period scale, which is used to record the requirements that need to be met during the process of two instantaneous events.
Corresponding to the sphere in Figure 5, the slack variable needs to be added due to the influence of discrete calculation, and different types of key-frames are introduced below: (1) Positional key-frame (KF-P) Short-term position constraints: where i ∈ [1, M], k ∈ [1, N + 1], and in incomplete position constraints, only a part of the position variables are constrained; take flight altitude as an example: During the flight, some state variables will be constrained for a long time such as altitude maintenance. Process state variables λ i divides t N into M + 1 fragments; take altitude maintenance as an example, we can set long-term position constraints as: (2) Attitude key-frame (KF-A) Short-term attitude constraints: If the requirements for posture are not so strict, we can set incomplete attitude constraint by Appendix A, for example, only constrain a certain angle: In the calculation process, the range of pitch angle is generally [−π/2, π/2], and the range of roll and yaw is [−π, π]. The above key-frames of roll angle and yaw angle are not one-to-one correspondence. The calculation will produce singularities, the constraint can be changed to the following equation, which will increase a certain amount of calculation.
(3) Compound key-frame (KF-C) In some special scenarios, there are requirements for the position and attitude of the drone. We proposed a compound key-frame, if the same µ i is used, the position and attitude can be constrained at the same time.
In addition to these constraints, if we want to ensure the order of key-frames, we need to set the following constraints.
Based on the above statements, our optimization variables are integrated as x opt = [t N , x 0 , . . . , x N ], where: The most basic goal of the designed maneuver is to minimize the time and require the maneuver to be completed in the shortest time. In order to improve the execution effect, we added the minimum energy term, and every motion primitive need to be adjusted according to its characteristics.
Algorithm 1 flow is summarized as follows: Algorithm 1. Fixed-wing UAV optimal maneuver generation.

1:
Input: L pieces of maneuver trajectory obtained by demonstration 2: Output: optimal trajectory of single maneuver ξ (divided into primitives ξ j ) 3: for ∀l ∈ [1, L] do 4: split the teaching trajectory into multiple maneuvers ξ 5: extract the key-frames p i , q i of each maneuver 6: decompose ξ into motion primitives ξ j according to certain principles 7: end for 8: for ∀j ∈ [1, J] do 9: if ξ j exists, or similar ξ j exists 10: use the related ξ j directly 11: else 12: constructing ξ j as keyframe-based optimization problems 13: end for 24: process ξ j and concatenate ξ j to ξ 25: return primitives ξ j concatenated maneuver data ξ

Experiments and Discussion
In this section, the identification results of the UAV model were obtained based on the flight data and three types of maneuvers are studied in this section. As mentioned in Section 4, the corresponding key-frames and motion primitives could be extracted for a specific maneuver. Therefore, we propose a technique to construct different minimum-time maneuver problem through adjusting the parameters of the boundary and key-frames. This algorithm could help to find the optimal and physical-realizable flight maneuvers.

Model Identification Results
The main parameters of the UAV are: mass m = 3.24 kg, air density ρ = 1.225 kg · m −3 , reference area S re f = 0.56 m 2 , wing span b re f = 1.83 m, mean chord length c re f = 0.30 m, moment of inertia J xx = 0.22 kg · m 2 , J yy = 0.31 kg · m 2 , J zz = 0.48 kg · m 2 , J xy = J xz = J yz = 0 kg · m 2 , and the gravity acceleration was assumed to be g = 9.8 m · s −2 .
According to Section 3, the structure of the aerodynamic coefficients were defined in Equations (7) and (8). The estimation of the unknown parameters was performed from the flight using the least square approach. During the flight demonstration, the pilot performed several types of maneuvers. Both of the control surfaces' deflections and the output of the system are recorded in time domain. Furthermore, to satisfy the kinematic and dynamic constraints, the upper and lower limits of the control inputs and system states are listed in Table 1.

Optimization Simulation Setup
In this section, we investigate several types of maneuvers and conducted simulation experiments separately. First, we evaluate the loop maneuver which contains a single primitive. Meanwhile, only the positional key-frames were employed for this motion. Then, another classical maneuver named the Immelmann turn was taken into consideration. In this motion, two parts of the primitives, which contains the positional and attitude keyframes, respectively, are concatenated. Thirdly, the half Cuban eight, which was similar to the former one, was also evaluated, and the compound key-frames were proposed. Furthermore, the concatenation of two half Cuban eight results in a whole Cuban eight, which is a sophisticated motion containing multiple positional key-frames and compound key-frames.
As is shown in Figure 1, maneuver optimization is an important part of the framework, which can be referred to in detail in Algorithm 1. In this paper, all the maneuver optimization problems were carried out on CasADi [45] with Ipopt [46] optimization approach. The initial flight status and parameters of the optimization are listed prior to the results. For ease of presentation, we chose the original point at the initial position before the aerobatic flight.

Loop Maneuver
The Loop is a maneuver mainly performed in the vertical plane. At the beginning of the motion, the UAV keeps trim flight and starts to pitch up. After the pitch angle finishes a 360-degree turn, the UAV returns to the beginning position and maintains the trim condition.
It is worthy to mention that, even though the human pilot is well trained, the geometrical size and shape of different demonstrated trajectories are not completely consistent. Therefore, it is almost impossible for the human pilot to realize an optimal maneuver. However, the non-optimal trajectories still share similar features. Through sufficient trialand-error, we found that the dominant feature of the loop was the position of key-frames. Nevertheless, the number/values of the positional key-frames are essential to obtain a reasonable circular trajectory. The set of the positional key-frames and the initial flight conditions are listed in Table 2. As for the optimization, we set the number of interval points N = 210, the positional tolerance d tol_p = 0.4 m, the parameters of the Equation (29) are set as σ = 1, τ = 0.1, and the simulation results are illustrated in Figure 6. As shown in Figure 6, the drone passed through all the positional key-frames in a short period of time (t N = 3.22 s). Once the drone meets the tolerance of each key-frames, the corresponding λ i changes from 1 to 0. Furthermore, the UAV returns to the initial position after finishing the whole maneuver, which is difficult for a human pilot. Even though only one primitive is considered in this scenario, the half loop maneuver can be extracted from the results naturally.   As shown in Figure 6, the drone passed through all the positional key-frames i short period of time ( = 3.22s N t ). Once the drone meets the tolerance of each key-fram the corresponding i λ changes from 1 to 0. Furthermore, the UAV returns to the ini position after finishing the whole maneuver, which is difficult for a human pilot. E though only one primitive is considered in this scenario, the half loop maneuver can extracted from the results naturally.

The Immelmann Maneuver
The Immelmann maneuver is also known as upside-down half-roll. During maneuver, the UAV first performs a half loop from trim flight. As soon as the airc reaches the top of the circular trajectory, it spins around the x-axis and executes a 1 roll. Finally, the UAV resumes to trim flight with an opposite direction compared to initial condition.
We decompose the Immelmann into two concatenated primitives: half Loop a 180  roll for a better description. For the first primitive, there are two ways to obtain near-optimal form: extracting from the Loop maneuver directly and modeling this mot into a two-point boundary value problem (BVP).
Modeling the optimization problem of 180  roll is significantly different from former cases. Rolling with high angular rates will not only result in a large displacem in the forward direction, but also lead to deviations both in height and heading an

The Immelmann Maneuver
The Immelmann maneuver is also known as upside-down half-roll. During the maneuver, the UAV first performs a half loop from trim flight. As soon as the aircraft reaches the top of the circular trajectory, it spins around the x-axis and executes a 180 • roll. Finally, the UAV resumes to trim flight with an opposite direction compared to the initial condition.
We decompose the Immelmann into two concatenated primitives: half Loop and 180 • roll for a better description. For the first primitive, there are two ways to obtain its near-optimal form: extracting from the Loop maneuver directly and modeling this motion into a two-point boundary value problem (BVP).
Modeling the optimization problem of 180 • roll is significantly different from the former cases. Rolling with high angular rates will not only result in a large displacement in the forward direction, but also lead to deviations both in height and heading angle. Firstly, we set a short-term attitude key-frame (KF-A) which contained a 90 • roll to ensure the motion completion. Then, in addition to the time factor, the displacement in the forward direction is also considered. The optimal solution is expressed as follows: Meanwhile, the final states (position, attitude, velocities, and angular rates) of the half loop maneuver was set as the initial condition of the rolling primitive. Furthermore, the y-axis and z-axis constraints were also added to the final states of the rolling primitive. The initial condition of the second primitive is listed in Table 3. We set N = 100, d tol_q = 0.04, the results can be obtained by calculation and the result of concatenated primitives for the whole Immelman maneuver is shown in Figure 7.

Key-Frame Value
Short-term angle key-frame θ = −π/2  As is shown in Figure 7, the entire Immelmann maneuver takes only 2.6155 s, and each state quantity is well connected. It can be seen that the entire roll primitive takes only 0.9631 s, the drone smoothly passes the 90 degree roll angle key-frame and flies forward 26.1 m, finally flies out horizontally. The entire Immelmann maneuver takes only 2.6155 s, and each state quantity is well connected.
It is worth noting that the concatenated primitives might be intrinsically sub-optimal. However, the sub-optimal solutions are sufficient for practical applications and capable of generalizing more maneuvers that never demonstrated.

Half Cuban Eight and Cuban Eight Maneuver
The Cuban eight is a sophisticate maneuver consisting of two 3/4 Loops followed by 180 • rolling. From the view of the ground, the UAV's trajectory is a vertical figure of eight. Due to the symmetric property, the left/right part of the motion, named half Cuban eight, can be evaluated first. To facilitate the calculation, the half Cuban eight was decomposed into a half loop (the same as the Immelmann turn's) and a special rolling primitive.
For the rolling primitive, we set a compound key-frame (KF-C), which requires the UAV to pass through the center of the Cuban eight with a 90 • roll angle. Due to the inertia of rotation, the UAV will continue to rotate around the x-axis. Furthermore, the angular limitations are considered to avoid large angle deviations. The initial flight conditions are listed in Table 4. The two key-frames can be formulated as follows: In the rolling primitive, we set N = 100, d tol_p = 0.1m, d tol_q = 0.04, the coefficients in Equation (29) are set as σ = 1, τ = 0.1 for both primitives, and we can connect the two primitives into a half Cuban eight, and the results are shown in Figure 8.  The results of the optimization of the half Cuban eight are illustrated in Figure 9. It takes 1.4118 s for the rolling primitive and 3.0642 s for the half Cuban eight. It can be seen that the UAV accurately crossed the center point with a 90 • roll angle, and the roll angle was always less than 5 • in the subsequent primitive. Nevertheless, the velocity at the end of the half Cuban eight maneuver is much larger than the one at the start point. This phenomenon will hinder the completion of the concatenate loop. Therefore, we considered limiting the speed to [20,0,0] of the drone at the end of the first rolling primitive, although this will increase the entire maneuver time to some extent.
Even though the two identical half-Cuban eights cannot be joined directly, there are still many similarities between them. Modeling a mirroring problem from the original problem and using data symmetry facilitate to obtain the optimization results. As shown in Figure 9, with the limitation on the terminal velocity, the duration time of the rolling primitive increased by 0.3898 s, and the engine thrust maintained at 0 N. Basically, the sub-optimal Cuban eight maneuver is completed, which takes 6.4856 s in total, and the primitives are well connected.

Benchmark Comparisons
To further illustrate the superiority of our approach, a comparison with three other methods is conducted. The first strategy is through collecting the maneuvering trajectories performed by the human pilot. The second one is the Gaussian pseudo-spectral method proposed in Morales's work [24]. Moreover, the 3-DOF model introduced [23] was also utilized in our approach, formulated in Appendix B. The simulation results are depicted in the Figure 10.
It can be seen from Figure 10; Figure 11 that the trajectories obtained from manual flights are inferior to the other methods. Even though Morales's method's performance is slightly better than the manual flight, it cannot solve the period constraints. The flight time is significantly longer than those of our proposed method and the 3-DOF model. The method proposed in this paper was based on the global model, and motion segmentation is needed increasing computational costs. It is worth noting that our approach can flexibly connect those segments, resulting in various maneuvers. In the result of the mass point model, since our model was simplified, and there was no need to segment the trajectories, the maneuvering time optimized was slightly shorter than our method. Therefore, the optimization results based on the 3-DOF model can be utilized as an initial guess in our approach or for benchmark comparison. Since the 3-DOF model was simplified by ignoring the kinodynamic constraints, it cannot be directly used for maneuvering planning.   In order to conduct a more thorough analysis, we performed statistics on the running time of several methods. All the experiments were executed on an Intel Core i7-4710MQ CPU with 8 GB RAM that runs at 2.50 GHz.
It can be seen from Table 5 that the minimal-time maneuver optimization is timeconsuming. Our work decomposes the optimization into the calculation of the 3-DOF model (initial guess) and the re-optimization of the global model. Compared with Morales's work, the optimized solution can be obtained faster.

Analysis of Maneuver Optimization Algorithm
It worth noting that the optimization process of maneuvers is quite time-consuming due to the inherent highly non-convex property. There are several methods that can accelerate the calculations and avoid the infeasible problem in Ipopt solver.
(1) Global Fixed-wing model Our optimization algorithm requires an accurate model in order to obtain the extreme motion of the UAV. Although the simplification of the model was conducive to simplifying the calculation, it was difficult to get the most realistic movement conditions. In our multiple maneuvering experiments, the performance of the drone was pushed to the limit, which was in line with the minimum time movement characteristics.
(2) Key-frames and primitives design To the best knowledge of the authors, it's unfeasible to propose a unified approach of key-frame extraction and motion decomposition for various maneuvers. Choosing the number of the key-frames and primitives requires a balance between computational complexity and geometric accuracy. Furthermore, most of the maneuvers can be decomposed into several longitudinal and lateral primitives.
(3) Initial guess The initial guess is essential for the nonlinear programming problem. In this paper, a reasonable initial guess was obtained either from the optimization results based on the point-mass model ( . p i = a, with input u = a, a ≤ a max ) or simply from the demonstration data using interpolation.
(4) Optimization constraints The small slacking variables can be added into the nonlinear constraints to accelerate the calculation and avoid the infeasible solutions. For instance, the complimentary constraint can be adjusted into: In addition, excessive slacking variables can be added into nonlinear constraints to obtain an initial solution, which can be considered as the initial guess in the next step of optimization iteratively.

Conclusions
In this paper, a hybrid-driven flight maneuver optimization framework was proposed. Based on the demonstrated maneuvers performed by an experienced pilot, the parameters of the UAV model along with the features of maneuvers were extracted. Then, the idea of key-frames and maneuver primitives were proposed which can design and concatenate different maneuvers more flexibly. Based on the above work, we formulated maneuver planning into a time-optimal problem, the feasibility of the framework was fully verified through designing the positional, attitude, and compound key-frames in the loop, Immelmann, half Cuban eight, and Cuban eight maneuvers. Through comparison, it can be concluded that the designed action was faster than professional pilots, and the results provide verifications on the optimality of the solutions and the effectiveness of the proposed approach.
Future works will concentrate on realizing our algorithm on a UAV platform experimentally.  Acknowledgments: Thanks to the research environment provided by the Institute of Unmanned Systems, the National University of Defense Technology, and the Nanjing Institute of Telecommunications Technology, National University of Defense Technology, for its administrative support.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix B
Aircraft Point-mass model: where κ is the aerodynamic roll angle, D, C, L are the components of the aerodynamic forces, described in the aerodynamic reference system, calculated as follows: −1/2 · ρV 2 a S re f C x −1/2 · ρV 2 a S re f C y −1/2 · ρV 2 a S re f C z    (A5)