Transition trajectory generation with time-scaled trigonometric series control parameterization

This paper is concerned with shaping the time behavior of controls for generating transition trajectories. Parameterizing controls with trigonometric series enables infinite differentiability with respect to time, which is very appealing in real operations. In practice, one may desire the controls to have a specific distribution over time. However, to achieve this usually relies on judicious selection of the cost function. To this end, this paper introduces a time-scaling transformation for the trigonometric series. It yields an intuitive and flexible way to design the controls tailored to the desired shape of a transition maneuver. The process does not rely on additional design effort for the cost function. The proposed method is applied to generating a transition trajectory that converges to the desired condition earlier in comparison to the unscaled formulation. Numerical results demonstrate the effectiveness of the proposed method.


Introduction
Trajectory generation is an essential practice in aerospace operations [2,7,19]. While discontinuous controls are common and possibly welcome in some applications [4,16], smoothness is very often considered to be an important feature of a reliable reference trajectory, especially in aeronautical operations [6,17,18]. For this purpose, there are dedicated trajectory generation approaches using fundamental geometric elements, such as lines, arcs, and clothoids, capable of generating trajectories with only discontinuities in high-order derivatives at the connecting points [3,20,21] 1 are generated using numerical optimization. For aerospace applications, one can refer to [1,5]. However, many existing methods nowadays based on numerical approaches may not be able to guarantee the continuity in higher-order derivatives. This is because the controls at every discrete point are set to be independent decision variables.
A trigonometric series-based trajectory generation approach is developed in [10] and has been applied in various applications [11][12][13]. It parameterizes the controls as trigonometric series and uses numerical optimization methods including convex optimization to determine the unknown series coefficients. The controls being a time-dependent trigonometric series are of C ∞ , i.e., infinitely differentiable. Therefore, smoothness of trajectories can be achieved. By imposing equality constraints, the boundary time behavior, namely, the continuity to the adjacent flight segments can be ensured as required. However, except for the boundaries, to acquire certain time behavior of controls relies on a judicious choice of the cost function. In fact, the cost function traditionally plays the sole role in deciding the control behavior in trajectory optimization [19]. Therefore, we seek for an innovative means to approach the goal within the trigonometric seriesbased trajectory generation.
In this paper, we develop a simple but effective way to shape the time behavior of controls parameterized by trigonometric series through time-scaling. A user-defined time scaling function is introduced to replace the original argument of the series basis function, namely, time, with a new independent variable. This variable is designed to have the desired time behavior and will pass the characteristic to the parameterized controls. Hence, having a specific time behavior of controls is decoupled from the design process of the cost function. Still, the cost function may affect the control behavior. However, the proposed method shapes the behavior in an external manner without modifying on the cost function. This philosophy is a unique advancement to previous works on the trigonometric series and thus a contribution to the literature.
The implementation of the newly proposed method is demonstrated in the convex-optimization framework developed in [8]. In this paper, it is shown that the proposed method only requires minimal modification of the optimization framework in its original form in [10], leading to a direct port. We show an application of transition trajectory generation, where a transition trajectory connects two steady flight conditions. A smooth connection is desired, which corresponds well to the capability of the trigonometric series parameterization method. For the transition maneuver, the application in this paper demands the aircraft to be more stable, i.e., having slow-varying control commands, in the later phase of the transition. This time behavior will better prepare the aircraft for the following flight segment. To this end, a time-scaling function based on a logistic function is selected to achieve this feature within the proposed method. Numerical simulations are carried out to verify the design.
The rest of this paper is organized as follows: Section 2 introduces the trigonometric series parameterization and the time scaling that shapes the time behavior of it prior to the determination of the series coefficients. Section 3 shows the details of the implementation of the proposed approach to a convex optimization framework. Section 4 shows the simulation results of the proposed method in combination with the convex optimization being applied to a transition trajectory generation task. Section 5 presents the conclusion of this paper.

Time-scaled trigonometric series control parameterization
In [10], a control parameterization method is introduced. A time-varying control variable u is described as an N th-order trigonometric series as where sin π Here, s N (t) is a function of time t ∈ 0, T f . Naturally, u (t) is of C ∞ after the parameterization. Moreover, the parameterization specifies u as a linear function of the series coefficients c N . It is intuitive that for Eq. (1), the time behavior of u over the entire horizon 0, T f is solely decided by c N which is to be determined by optimization. This is generally not a disadvantage because of its simple structure. To achieve a particular time behavior of the control variable, one often relies on designing the cost function which affects the selection of c N . This is not always favorable as such an indirect process may lead to a Pareto objective and thus possibly causes the cost function to lose advantageous properties, e.g., convexity. In addition, it is not effortless to relate certain behaviors to mathematical expressions in a cost function.
To mitigate this, we introduce a time-scaling transformation to Eq. (1). The time-scaling transformation originates in [15] and yields a new optimization problem formulation over a new horizon. Meanwhile in this paper, the proposed time-scaling transformation only influences the trigonometric series parameterization in Eq. (1). It introduces a new variable that monotonically increases with respect to time: where p ∈ [0, P]. Here, v (t) > 0, ∀t ∈ 0, T f and it is to be designed by user in order to achieve desired behavior. It shall be noted that P does not need to be specified because the definite integral of v (t) over t ∈ 0, T f gives P. Next, we reformulate the control parameterization in (1) to be a function of p: The time behavior of the control u is no longer only affected by c N but preserves characteristics of p (t) which are defined already in the designing process of v (t). This can be inferred using the chain rule as It is a common measure for numerical optimization approaches that a time grid is given over the entire time horizon. The discrete time grid specifies the total number of the discrete points. Provided that T f is fixed, the value of s N at each discrete point is known a priori. That is to say, s N does not change during the optimization algorithm execution. Therefore, it can be evaluated only once before executing the optimization program. For some optimization approaches including the one utilized in this paper, this feature saves some computational effort and streamlines the formulation.
The only remaining issue is the design of v (t). Ideally, it should be rather flexible. In Sect. 4, the design is shown along with the actual need of the application. As a remark, p must be smooth in order to ensure that u is smooth.

Implementation in model predictive convex programming
In general, the proposed time-scaled trigonometric series parameterization works with various optimization approaches, such as successive convex optimization approaches [11,13], because it is only a linear constraint. In this paper, to show the difference (or the similarity in a positive sense) from the formulation in [10], we again consider the model predictive convex programming (MPCP) [8].
Only the most salient steps of the MPCP in [8] are presented in this section.
First, we write the trigonometric series parameterization in a discrete time form for the control vector where and the superscript indicates the index of the components. The total number of discrete points is denoted by N t . Because of the linearity, the incremental change of u k , denoted as du k can be expressed by where dc N is the incremental change of the coefficient vector c N . It is important to notice that, as mentioned in the previous section, regardless of whether S N is a function of t as in Eq. (1) or p as in Eq. (7), as long as v and T f are specified, [S N ] k is fully defined for all k = 1, 2, . . . , N t − 1. Hence, the time-scaling of the parameterization and the functional dependency make no difference to Eq. (12) and thus have no impact on the following derivations. For introducing the MPCP, a discrete-time form of the system dynamics is given as where F is continuously differentiable, and x and u are the state and the control vector, respectively. The MPCP is an iterative approach that updates the state and control histories as where the superscript p denotes the variables in the previous iteration. Using a Taylor series expansion and neglecting higher-order terms, Eq. (13) can be expressed as which yields Eq. (15) can be further expanded: where Here, the initial condition x 1 is considered to be fixed, i.e., d x 1 = 0. Hence, Eq. (16) can be written as where k = 1, 2, . . . , N t − 1. This equality is derived in [14], connecting the state increments d x to the control corrections du. Using a static Lagrangian multiplier, Eq. (21) can be used to solve a class of optimization problems and it is referred to as the model predictive static programming. For incorporating more constraints, especially path constraints, a convex optimization problem can be formulated. Substituting Eq. (12) into Eq. (21) yields While S N is different due to the time scaling, Eq. (22) is identical to the one in [10]. This shows that the proposed method is compatible with the previous works, and a very intuitive port can be realized.
To formulate a convex optimization problem, we select the cost function as As one might notice, this quadratic cost function simply minimizes the incremental changes of the coefficients. It can be interpreted as minimizing the change of controls between iterations. It does not address the time behavior of the controls. Constraints can be reformulated in the following way where the elements of the lower triangular matrix B are B k, j in Eq. (21). Notice that du is the only unknown variable in Eq. (24) and it can be related to dc using Eq. (12). Therefore, constraints are ultimately functions of dc. A convex optimization problem can be formulated to iteratively search for suitable dc that satisfies all constraints. More details can be found in [10].
With a dedicated solver, the order of the trigonometric series does not significantly influence the computational efficiency for solving the convex optimization problems. Moreover, for a maneuver that is not complex, a low order suffices. In this work, we select N = 3. On the other hand, the total number of the discrete points, N t , does affect the computational efficiency. However, this setting is decided by the required accuracy of the trajectory rather than the proposed control parameterization. Over a 20 s time horizon, N t = 101 is considered to give a sufficiently accurate result for the application introduced in the next section.

Transition flight trajectory generation
In this section, we adopt the application considered in [9] to demonstrate the effectiveness of the proposed method and in Sect. 3, while also comparing the results to the original form in [10].
The scenario is to generate a transition flight trajectory connecting two steady flight conditions. The objective is to make the maneuver relatively milder at the later phase of the transition, which prepares the aircraft for the next flight segment.

Problem formulation
A fixed-wing unmanned aerial vehicle point-mass model is considered: where x, y, h are coordinates of the aerial vehicle, and χ , γ , and V are the flight path course angle, the flight path climb angle, and the speed, respectively. μ is the bank angle. m is the mass of aerial vehicle. g is the gravitational acceleration. The thrust and aerodynamic forces in Eq. (25) are expressed by where T max is the maximum thrust, δ T is the thrust lever position, ρ is the air density, S is the reference area, and k i is the induced drag factor. C L and C D 0 are the lift coefficient and the zero-lift drag coefficient, respectively. The state, control, and output vectors of the dynamics in Eq. (25) are given as The simulation scenario is to generate a transition trajectory connecting a coordinated turn with a steady descent. The equality constraints are given as (32) On the right-hand side of the equality constraints in Eq. (32) are all constants. The path constraints are specified by

Time-scaling function
The goal of designing the time-scaling function is to properly distribute the control authority for the transition maneuver. In accordance with the need of flying more stably control-wise in the later phase of the transition, we select the following function: It is a logistic function and p can be obtained using Eqs. (4) and (5). The time histories of v and p are shown in Fig. 1.
Because the problem at hand is a fixed-final-time trajectory optimization problem, the variable p and the series basis functions are fully known and they do not change with the optimization runs. As can be seen from Fig. 1, v is monotonically decreasing as time evolves, where v| t=0 = 0.77 and v| t=20 = 0.22. Here, the slope is gentle, as an aggressive maneuver during the earlier stage, or technically any stage, is not desired. Referring to Eq. (8), the control derivatives are relatively more slowly varying during the later phase of the transition.

Simulation results
The simulation results are presented in this section. The quadratic programming problems are solved with a dedicated solver-OSQP [22]. The formulation in [10] is referred to as the original method, while the time-scaled trigonometric series is called the proposed method.
The flight trajectories generated by the original method and the proposed method are shown in Fig. 2. The trajectories are similar but the maneuver is less aggressive for the proposed method in the later phase. This can be confirmed by the state histories given in Fig. 3. The trajectories generated by both methods lead the aircraft to zero course angle and −5 deg climb angle. The speed is maintained at 50 m/s as desired. The main occurrence of the state variations has been shifted to the earlier part of the maneuver. It is also visible that Although the proposed method parameterizes the controls, the state trajectories are the real system responses to be manipulated. Speaking of the controls, their histories can be seen in Fig. 4. The state histories are all smooth and satisfy every constraint. It is worth noting that boundary conditions, i.e., the trim conditions, can be easily enforced for both the original and the proposed methods. These are important conditions that ensure continuous connections to the steady flights. For the proposed method, slightly more control authorities are used in the earlier phase of the transition maneuver. The maximum utilization increases by 6.6%, 5.0%, and 0%, for C L , μ, and δ T , respectively. On the other hand, the controls approach faster to being steady-state and the desired terminal conditions in the later phase. This means that the aircraft at a relatively earlier time flies practically at a trim condition using the proposed method. Inspecting the time histories of the control derivatives in Fig. 5 reveals more details. The trade-off for getting steadier controls at the later phase is committing marginal but still higher change rate of the controls in the earlier phase. The incremental changes of the maximum absolute control rates are 0.01 s −1 , −0.05 deg/s, and 0.02 s −1 , forĊ L ,μ, andδ T , respectively. Neverthe-  Fig. 6 Time histories of the normalized cumulative sums of the control derivatives less, the control derivatives are smooth and well within the bounds. In fact, the mathematical property of the trigonometric functions guarantees the existence and continuity of arbitrary derivatives. The derivatives are set to be zero at both ends of the transition trajectory, which is achieved by both methods. Next, we introduce an indicator to quantify completeness of the transition over time. This gives a very intuitive view of the shifted distribution of the control derivatives using the proposed method. The indicator is expressed by The indicator stands for the normalized cumulative sum of the control derivative over time, which is depicted in Fig. 6. The vertical lines identify where the indicator reaches 95%. The indicators of the proposed method approach the condition earlier than those of the original method. Specifically for C L , the indicator reaches 95% 2.6 s earlier for the proposed method, which is already 13% of the 20 s time horizon for the transition maneuver. The differences of μ and δ T are found to be 4% and 8%, respectively. As the new independent variable v has been designed to be of a slow-varying nature, the results of the proposed method and the original method are not vastly different, and yet we can observe the effectiveness of the proposed method. In fact, the difference can be more visible if a very aggressive time-scaling function had been selected. However, it is beyond the purpose of this application that requires an overall smooth transition. It is to be noted that there can be different ways of designing the time-scaling function and the time-scaled trigonometric series. The shown application only represents one specific concern that has been satisfactorily addressed.

Conclusions
For shaping the time behavior of the controls in generating trajectories, this paper develops a new time-scaled trigonometric series control parameterization method. The proposed method can incorporate a user-defined time-scaling function that alters the fundamental time behavior of the basis functions such that the controls inherit the desired characteristic. It has been shown that in the model predictive convex optimization framework, the modification required for using the newly developed method is minor. For the transition maneuver scenario in this paper, the aircraft is required to maneuver more gently in the later stage of the transition trajectory, for which the proposed method introduces a smooth, decreasing time-scaling function. Simulations are carried out to show the effectiveness of the proposed method and the comparison to the original method demonstrates that the proposed method accomplished the task well.

Conflict of interest
The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.