Multiaircraft Optimal 4D Trajectory Planning Using Logical Constraints

,


Introduction
In air traffic management (ATM), the flight of several aircraft can be modeled as a hybrid dynamical system, which can be regarded as a set of interacting continuous dynamical systems.A number of frameworks have been proposed to model hybrid dynamical systems, in which, in general, differential equations describe the dynamics of each system, whereas logical constraints describe the behavior of the systems during the interactions among them and the interaction with the environment in which they operate.In the ATM context, logical constraints describe, for instance, policies to apply in conflict detection and resolution and operational constraints to be fulfilled during flight.The main operational constraints to be fulfilled during flight are separation constraints, keep-out constraints to avoid no-fly zones, and passage constraints through or by waypoints [1].
Given a set of aircraft, separation constraints between them can be expressed as follows: pairwise, they must keep a vertical distance greater than a minimum vertical safety distance or a horizontal distance greater than a minimum horizontal safety distance.The minimum horizontal separation distance can be fixed or variable.In the latter case, it can be established based on the turbulence generated by the preceding aircraft and the ability of the following aircraft to resist turbulence [2].Obstacles and no-fly zones are in general polyhedral regions of airspace.However, the corresponding keep-out constraints are usually introduced by bounding ellipsoids around obstacles.In some cases, this is a coarse approximation.Keep-out constraints from a polyhedral region of airspace can be expressed as follows: each aircraft must stay outside one of the half-spaces defined by the planes that supports the faces of the polyhedron.This method for modeling keep-out constraints from obstacles can be extended to model passage constraints through windows or waypoints in the airspace.In this manner, the multiphase modeling of the problem is avoided.This is an interesting possibility, since multiphase optimal control models imply the introduction of additional variables to the problem, such as the duration of the phases and additional constraints such as the linkage constraints between them to enforce continuity of the state variables between contiguous phases [3].
It is easy to see that all the constraints mentioned above are expressed in disjunctive form.Standard modeling techniques are able to tackle constraints in disjunctive form using binary variables.The four-dimensional (4D) trajectory planning problem for multiple aircraft with logical constraints in disjunctive form can be solved as an optimal control problem (OCP) for a hybrid dynamical system and a common approach for solving this class of problems is to formulate them as a mixed-integer programming problem.In [4], the optimal cooperative three-dimensional (3D) conflict resolution problem among multiple aircraft has been solved in which separation constraints among aircraft expressed in disjunctive form have been included in the model using continuous auxiliary variables.In [5], the optimal path planning problem for multiple unmanned aerial vehicles in the horizontal plane with collision avoidance has been studied, in which constraints for collision avoidance with rectangular obstacles expressed in disjunctive form are included in the model using continuous auxiliary variables.In [6], the trajectory optimization problem for multiple aircraft landing on a single runway in the presence of constraints on the airspace has been treated.The constraints considered are passage constraints through windows in the airspace, and optimal trajectories have been determined by solving a nonsequential constrained multiple-phase optimal control problem.
In this paper, the embedding technique proposed in [5] to model rectangular obstacle avoidance in the horizontal plane has been extended to model time-based and distancebased separation constraints and passage through waypoints constraints in trajectory optimization for multiple aircraft.The modeling of passage constraints through waypoints has been done by defining vertical walls in the airspace with a cuboidal window around the waypoint.In this way, introducing multiple phases in the model to enforce passage through waypoints is avoided.Moreover, the dimensions of the windows can be easily calibrated to induce a fly-by or a flythrough the waypoint.
This study can be classified into the category of continuous descent operations (CDO) [7].During CDO, aircraft descend from the cruise altitude to the final approach fix at or near idle thrust without level segments at low altitude minimizing the need for high thrust levels to remain at a constant altitude and reducing the environmental impact.Actually, the term CDO makes reference to the different techniques to maximize operational efficiency and, at the same time, fulfilling local airspace requirements and constraints.These operations are known as continuous descent arrivals, optimized profile descents (OPDs), tailored arrivals, 3D path arrival management, and continuous descent approaches (CDA).In particular, an OPD is a descent profile normally associated with a standard terminal arrival route (STAR) and designed to allow maximum use of a CDO.Planning CDOs is one of the functions of the so-called arrival managers (AMANs) whose purpose is to ensure an optimal sequencing and spacing of arrival traffic [8].
Most of the previous research on CDOs based on optimal control theory focused on the trajectory optimization of a single aircraft.In [9], a multiphase optimal control method based on the pseudospectral technique has been employed to optimize vertical trajectories for individual aircraft in CDAs.Since the lateral path is assumed to be given by a STAR procedure, this work focused on optimizing vertical profile only using time and fuel consumption as performance indices.All the phases are formulated based on operational constraints and flap/gear schedules.The initial along track distance is free.Hence, it is possible to calculate both the optimal top-of-descent (TOD) and CDA trajectory.The optimal trajectories have been computed for two aircraft types: a Boeing 737-500 and a Boeing 767-400.
In [10], the vertical trajectory optimization for the en route descent phase of an aircraft has been studied in the presence of both along-track and cross winds, which are both modeled as functions of altitude.Flight idle thrust was assumed during the entire descent phase.The problem is formulated as an optimal control problem.The flight range was specified from a point during the latter stages of the cruise to the meter fix.Calibrated air speed (CAS) and Mach constraints, which are the state path constraints, are considered, along with flight path angle constraints, and a maximum descent rate limit, which is a mixed input and state path constraint.The descent trajectory is optimized with respect to two cost functionals: fuel and emissions.The effects of wind speed, windshear, and cross-wind on the optimal trajectory have been analyzed using the models of two types of aircraft, Boeing 737-500 and Boeing 767-400.
Less research efforts have been devoted to combined optimization of trajectories of multiple aircraft and sequencing for approaching a terminal manoeuvring area (TMA) in which all aircraft follow CDAs, while satisfying the operational requirements.This fact motivated the study presented in this paper.
Two types of CDA exist depending on the lateral path followed, generally referred to as CDA under vectoring and advanced CDA.In the first case, the lateral path followed by the aircraft is assumed to be specified through instructions provided by the air traffic control (ATC).In the second case, the lateral path of the aircraft is predefined and is based on the STAR.
Each of the problems that has been solved to validate the method proposed in this paper can be stated as follows.Given the dynamic models of a set of aircraft, their initial and final states, a set of operational constraints, find the optimal trajectories that steer the aircraft from the initial to the final states, fulfilling all the constraints and optimizing an objective functional.

2
International Journal of Aerospace Engineering In particular, the optimal trajectories of multiple aircraft in converging arrival routes are computed taking into account time separation constraints, distance separation constraints, and their optimized profile descent along a STAR lateral profile.The problem has been solved using optimal control techniques.In particular, the OCP is transcribed using a Hermite-Simpson collocation method [11].The resulting nonlinear programming (NLP) problem has been solved using the NLP solver IPOPT [12].
The paper is structured as follows.In Section 2, the general optimal control problem for multiple dynamical systems is stated and the direct collocation approach for its resolution is described.In Section 3, the aircraft equation of motion and the flight envelope constraints are stated.In Section 4, the general approach to model logical constraints is presented, which is then particularized to model timebased and distance-based separation constraints between aircraft, obstacle avoidance constraints, and waypoint constraints.In Section 5, the results of the application of the proposed method to solve several trajectory optimization problems for multiple aircraft with logical constraints in disjunctive form are reported and discussed.Finally, in Section 6, some conclusions are drawn.

Optimal Control Approach
2.1.Statement of the Optimal Control Problem.The multiaircraft flight planning problem considered in this paper can be regarded as a multitrajectory optimization problem in which the motion of each aircraft has been modeled as a differential algebraic dynamic system for p = 1, 2, ⋯, N p , where f p describes the right-hand side of the differential equation and g p describes the algebraic constraints where Since this multiaircraft flight planning problem also involves operative performances and flight envelope conditions for multiple aircraft, as well as the optimization of a specified performance index, the multitrajectory optimization problem can be formulated as an OCP of a set of dynamic systems in which the goal is to find the trajectories and the corresponding control inputs that steer the states of the systems between two configurations, satisfying a set of constraints on the state and/or control variables while minimizing an objective functional.
Therefore, the optimal control problem considered in this work can be stated as follows: where The objective function is given in Bolza form.It is expressed as a combination of a Mayer term and a Lagrange term Functions provides the terminal conditions at the final time t F , and it is assumed to be twice differentiable.The system must also satisfy algebraic path constraints within the time interval ½t I , t F given by the vector function with lower bound ϕ l ∈ ℝ n ϕ and upper bound ϕ u ∈ ℝ n ϕ .Function ϕ is assumed to be twice differentiable.
In the objective function (4), the Lagrange term represents a running cost, whereas the Mayer terms represent a terminal cost.A usual Lagrange objective function is to minimize the total amount of fuel consumed during the flight.A typical Mayer objective function is to minimize the duration of the flight.Equations (4b) and (4c) represent the differential-algebraic equation system that governs the motion of the dynamical system, e.g., the aircraft.Equation (4d) models the physical performance limitations of the dynamical system, typically expressed as upper and lower bounds on both states and control variables.Equations (4e) and (4f) denote the boundary (initial and final, respectively) conditions of the aircraft.Note that Equations (4c) and (4d) will also include the logical constraints that model conflict detection and resolution, and operational constraints as described in Section 4, which are of special interest for the problem studied in this paper.
Hence, the optimal control problem (4a), (4b), (4c), (4d), (4e) and (4f) consists in finding an admissible control u * ðtÞ such that the set of aircraft follow an admissible trajectory x * ðtÞ between the initial and final state that minimizes the performance index JðxðtÞ, uðtÞ, s, tÞ.The final time, t F , may be fixed or free.

Direct Collocation Transcription of the Optimal Control
Problem.A direct numerical method has been employed to transcribe the OCP into a NLP problem.More specifically, a Hermite-Simpson direct collocation method [11] has been used.The time interval ½t I , t F has been subdivided into N subintervals of equal length, whose endpoints are with t 0 = t I and t N = t F .In each subinterval ½t i , t i+1 , i = 0, ⋯, N − 1, the Hermite-Simpson numerical integration scheme has been used.
The set of constraints of the resulting NLP problem includes the Hermite-Simpson system constraints that correspond to the differential constraint (4b) and the discretized versions of the other constraints of the optimal control problem.They include the algebraic constraints (4c), the state and control envelope constraints (4d), and the boundary conditions (4e) and (4f).The unknowns of the NLP problem are the values of the state and the control variables at the endpoints of each subinterval ½t i , t i+1 , i = 0, ⋯, N − 1.
To solve the resulting NLP problem, the open source IPOPT solver [12] has been employed.It implements an interior point line search filter method and it is able to handle properly large-scale sparse nonconvex NLP problems, with a large number of equality and inequality constraints.Source and binary files are available at the Computational Infrastructure for Operations Research (COIN-OR) website (https://www.coin-or.org/).

Aircraft Model Description
Following [13], a common three-degree-of-freedom dynamic model has been used which describes the point variable-mass motion of the aircraft over a spherical Earth model.In particular, a symmetric flight has been considered.Thus, it has been assumed that there is no sideslip and all forces lie in the plane of symmetry of aircraft.

Equations of Motion.
The following equations of motion of the aircraft have been considered: The three dynamic equations in (14) are expressed in an aircraft-attached reference frame ðx w , y w , z w Þ, and the three kinematic equations are expressed in a ground-based reference frame ðx e , y e , z e Þ as shown in Figure 1.
The states of the system ( 14) are V, χ, γ, λ e , θ e , h e , and m.Thus, xðtÞ = ðVðtÞ, χðtÞ, γðtÞ, λ e ðtÞ, θ e ðtÞ, h e ðtÞ, mðtÞÞ.The state variables V,χ, and γ are the true airspeed, the heading angle, and the flight path angle, respectively.The state variables λ e , θ e , and h e are the aircraft 3D position, the longitude, the latitude, and the altitude, respectively.The components of the aircraft position vector in two dimensions (x e , y e ) are Note that differential equations in ( 14) take the form of (4b) of the continuous optimal control problem stated in Section 2.1.

Flight Envelope Constraints.
Flight envelope constraints are derived from the geometry of the aircraft, structural limitations, engine power, and aerodynamic characteristics.The performance limitation model and the parameters have been obtained from the Base of Aircraft Data (BADA), version 3.6 [14]: In (15), h M0 is the maximum reachable altitude and h u ðtÞ is the maximum operative altitude at a given mass (it increases as fuel is burned).MðtÞ is the Mach number and M M 0 is the maximum operating Mach number.C v is the minimum speed coefficient, V s ðtÞ is the stall speed, V M 0 is the maximum operating CAS, and a n and a l are, respectively, the maximum normal and longitudinal accelerations for civilian aircraft.Finally, T min ðtÞ and T max ðtÞ correspond to the minimum and maximum available thrust, respectively, and μ corresponds to the maximum bank angle due to structural limitations.
Note that inequality constraints in (15) take the form of (4d) of the continuous optimal control problem stated in Section 2.1.

Logical Constraints Modeling
In this section, the approach proposed in [5] has been followed in which an extension of the embedding optimal control technique stated in [15] and developed in [16] was proposed.The embedding technique introduced in [15,16], to transform hybrid optimal control problems into traditional smooth optimal control problems, in which the discrete aspect of the system arisen only from switches in the dynamic equations, was adapted in [5] to deal with the logical (discrete) components which also might appear as constraints.
It was shown in [17] that every Boolean expression can be transformed into conjunctive normal form (CNF). Thus, it has been assumed that any logical constraint considered in this study can be written as a CNF expression where is in disjunctive form, in which only one of the propositions must be satisfied and the proposition P j i is either X j i or ¬X 5 International Journal of Aerospace Engineering where g j i : ℝ n x ⟶ ℝ is assumed to be a C 1 function.In order to include the logical constraint (16) in a smooth continuous optimal control problem formulation, it must be converted into a set of equality or inequality constraints in which binary variables are not employed.In this way, the combinatorial complexity of integer programming is eluded.
Transcribing the conjunction in ( 16) is straightforward since it is equivalent to the following expression: Thus, taking into account (17), the logical expression ( 16) can be represented as follows: To transcribe these disjunctions into a set of inequality constraints, a continuous variable β j i ∈ ½0, 1 is defined and related to each P j i in (18).Thus, (20) can be expressed as follows: It is immediate to check that if β j i = 0 in the first term in (21), then constraint g j i ðxðtÞÞ ≤ 0 is not fulfilled.On the other hand, if 0 < β j i ≤ 1, then β j i • g j i ðxðtÞÞ ≤ 0 is in fact g j i ðxðtÞÞ ≤ 0, and thus, constraint g j i ðxðtÞÞ ≤ 0 is enforced.Finally, the last term in (21) guarantees that at least one of the propositions P j i holds.Note that, as expected, equality and inequality constraints in (21) take the form of (4c) and (4d), respectively, of the traditional continuous optimal control problem stated in Section 2.1.In the following subsections, the application of this technique to several instances of interest in the ATM context will be presented in detail.

Collision Avoidance Constraints
4.1.1.Time-Based Separation between Aircraft.One of the problems of most interest in ATM is the conflict detection and resolution problem [18].Two different instances of this problem will be studied in this subsection and in the next one.In a first instance, a collision avoidance model among different aircraft along routes converging at the same waypoint has been considered, in which a safety time-based separation is guaranteed at the merging waypoint.
Let t p F and t q F be the unfixed time at the merging waypoint of aircraft p and q, respectively, and let d t be the safety time difference between two consecutive aircraft.Since in terms of the discretization (13), t p F = t p N and t q F = t q N , then multiaircraft time-based separation constraints can be expressed as follows: where condition q > p prevents unnecessary duplication of constraints.Constraints (22) can be rewritten as follows: Following the technique described above to tackle constraints in disjunctive form, if we define new variables δ 1 p N ,q N , δ 2 p N ,q N ∈ ½0, 1 satisfying condition Equation ( 23) can be transformed into The last constraint in (25) ensures that at least one of the constraints in (23) is fulfilled, that is, the safety timebased separation between aircraft p and q is guaranteed.
4.1.2.Distance-Based Separation between Aircraft.In the second instance, to model collision avoidance among different aircraft, at each endpoint of the subintervals of the discretization (13) and for every couple of aircraft, a safety distancebased separation is guaranteed in the horizontal or vertical directions, that is, λ − θ or h directions.
Let ðλ p i , θ p i , h p i Þ and ðλ q i , θ q i , h q i Þ be the positions of aircraft p and q at the endpoint t i of the discretization, respectively.If the safety distance-based separation in the λ − θ and h directions is denoted by d λθ and d h , respectively, then the distance-based separation constraints can be expressed as follows: where the haversine formula has been considered for the distance in the λ − θ direction with 6 International Journal of Aerospace Engineering As before, condition q > p prevents unnecessary duplication of constraints, and the set of constraints (26) can be rewritten as follows: Again, following the technique described above, new variables ν are introduced, and Equation (28) can be transformed into The last constraint in (30) ensures that at least one of the constraints in (28) is fulfilled, that is, aircraft p and q are guaranteed to be a safe distance apart.Moreover, since this embedded logical constraint approach is quite general, any other distance separation model described in terms of Equation ( 17) can be considered.4.1.3.Obstacle Avoidance.In the third instance, an obstacle avoidance problem has been considered.Following [5], in which a set of stationary rectangular obstacles in a twodimensional environment were assumed, stationary obstacles on the three-dimensional space have been enveloped by cuboids.Cuboids give computational advantages with respect to other models, such as ellipsoids, since they give rise to simple bound constraints instead of quadratic constraints, and mathematical programming solvers are able to deal with bound constraints more efficiently than with nonlinear constraints.Additionally, this simple approach allows us to specify each cuboid by giving only the coordinates of two opposite corners as shown in Figure 2.
Let ðλ l , θ l , h l Þ and ðλ u , θ u , h u Þ be the positions of opposite corners of the cuboid.Cuboid avoidance for a single aircraft involves that at every endpoint t i of the subintervals of the discretization( 13), the position of the aircraft ðλ i , θ i , h i Þ must remain outside it.In terms of logical constraints, this condition can be expressed as follows: Note that constraints (31) are not enforced at the initial and final points of the discretization, t 0 and t N , since in practice both positions are fixed.Moreover, due to the nature of the collocation technique described in Section 2.2, constraints are not enforced between the endpoints of the discretized differential system of equations.Therefore, small intrusions inside the cuboids might occur.A smaller discretization size implies a larger intrusion, and vice versa.In practice, this drawback can be easily overcome by considering a sufficient safety margin in the modeling of the cuboids.
Following the same technique used in the previous sections, new variables are introduced, and Equation (31) can be transformed into 7 International Journal of Aerospace Engineering Again, as mentioned above, the last constraint in (33) ensures that at least one of the constraints in (31) is fulfilled, which means that aircraft flies outside the cuboid.Also note that the set of constraints in (33) can be straightforward extended to multiobstacle and multiaircraft settings.As already pointed out in Section 4.1.2,since this embedded logical constraint approach is quite general, any other obstacle model described in terms of Equation ( 17) can be considered.

Waypoint Constraints.
Finally, the modeling of an aircraft flying through a waypoint has been considered.The model has been based on the use of cuboids similar to those considered in Section 4.1.3.In this case, a cuboid centered at each waypoint has been defined.However, unlike the obstacle avoidance constraints, in this setting, the aircraft has forced to pass through the cuboid instead of avoiding it.
Let ðλ W l , θ W l , h W l Þ and ðλ W u , θ W u , h W u Þ be the positions of opposite corners of a cuboid surrounding a single waypoint.Flying by this waypoint (that is, passing through the corresponding cuboid) involves that at every endpoint t i of the subintervals of the discretization (13), the position of the aircraft ðλ i , θ i , h i Þ must remain inside it.In terms of logical constraints, this condition can be expressed as follows: Note that, for one hand, constraints (34) are enforced at every point of the discretization except for the initial and final points, t 0 and t N , since there is no a priori knowledge about when the aircraft is going to fly by the waypoint.On the other hand, these constraints obviously make sense only when the aircraft is closed enough to the waypoint.
To overcome this drawback, a second auxiliary cuboid is considered to modelled free flight mode of the aircraft.Let λ min , θ min , h min , λ max , θ max , and h max be the minimum and maximum values of the state variables λ, θ, and h, respectively.In terms of logical constraints, the free flight mode condition can be expressed as follows: Then, the transcription into a logical disjunction which allows to select along the whole trajectory between flying by the waypoint mode (WM) or free flight mode (FM), namely, WM∨FM, can be expressed as follows: where WM i and FM i denote if at discretization instant i, the aircraft is in waypoint mode or free flight mode, respectively.Once again, following the technique described above, if we define new variables κ 1 i , κ 2 i ∈ ½0, 1 satisfying condition Equation (36) can be transformed into The last constraint in (38) ensures that at least one of the conditions in (36) is fulfilled, which means that at each 8 International Journal of Aerospace Engineering discretization instant i, the aircraft flies in FM or WM.Note that, depending on the performance index (4a) considered in the optimal control problem, two related undesired issues could potentially arise.On the one hand, the optimal solution could provide a trajectory in which the aircraft flies in FM for each discretization instant i.Therefore, in order to force the aircraft to actually fly by the waypoint, a penalty term must also be added to the numerical transcription of the performance index (4a).The use of this penalty term, which is set forth to encode the desired control objectives, implies that the numerical solution of the optimal control problem must combine the usual collocation technique describe in Section 2.2 with a penalty function methodology.In particular, in this work, the well-known continuation method has been implemented following a similar approach to [19].In this specific case of waypoint modeling, the penalty term takes the form where c 1 and c 2 are suitable constants determined by the continuation method.In the context of a minimization performance index such as (4a), a large enough value c 2 such that c 2 ≫ c 1 > 0, which penalizes the free flight mode, guarantees that the aircraft actually flies by the waypoint.
On the other hand, the optimal solution could provide a trajectory in which the aircraft flies by the waypoint more than once.This situation can be easily avoid introducing in the model the following simple constraint: where c 3 ∈ f1, 2, ⋯, N − 1g is a suitable constant.In the context of the problem considered in this work including the penalty term (39), a value of c 3 between 2 and 4 is enough to avoid this potential undesired issue.Note that this waypoint constraint setting can be straightforward extended to multiwaypoint and multiaircraft settings.Moreover, once again, since this embedded logical constraint approach is quite general, any other waypoint constraint described in terms of Equation ( 17) can be considered.

Numerical Results
To show the effectiveness of the methodology describe in Section 4, the following three numerical experiments have been carried out: All of these numerical experiments involve Airbus A-320 BADA 3.6 aircraft models in which the performance index is the sum of the duration of the flights of the three aircraft.In some of them, constraints derived from current flight regulation have been introduced, namely, time-based separation and distance-based separation operational constraints have been imposed.In particular, a minimum horizontal distance separation between aircraft of 5000 m and vertical separation between aircraft of 1000 m has been considered, whereas the considered minimum time separation between aircraft has been 200 s.These specific values have been chosen taking as a reference the aviation regulation [20] in which, in general, aircraft have to be separated by at least 3 NM or 3 min in the TMA.
The numerical experiments have been conducted on an Intel Core i7 2.8 GHz CPU with 16 GB RAM.

Experiment 1: Minimum-Time Continuous Descent with
Time-Based Separation Constraints.In this experiment, a CDA under vectoring has been considered, that is, the lateral path followed by the aircraft has been assumed to be specified through instructions provided by the ATC.In particular, the boundary conditions of the state variables have been selected from the chart of the Adolfo Suárez Madrid-Barajas (LEMD/MAD) TMA shown in Figure 3.The initial position of Aircraft 1, Aircraft 2, and Aircraft 3 has been supposed to be coincident with the ROLDO, SOTUK, and MORAL waypoints, respectively.Their common final position has been assumed to be the LALPI waypoint.The three aircraft have been supposed to go directly from their initial positions to the final waypoint, that is, no constraints on the lateral path have been considered.Neither timebased nor distance-based separation constraints have been included in the model.The initial mass of the three aircraft has been assumed equal to the maximum landing weight of the aircraft.The specific boundary conditions of the state variables are given in Table 1 and the final mass and time given by the solution of the corresponding OCP are reported for each aircraft in the first row of Table 2.
In Figure 4, the 3D view of the paths obtained in the solution is represented in thick lines, whereas in Figures 5  and 6, the horizontal and vertical profiles are represented in solid lines.
The computation time to find the solution has been 15:8098 s.The final time of the three aircraft is 1539 s, 1419 s, and 1350 s, respectively, corresponding to 22-25 min of flight.The difference between the final time of aircraft 2 and 3 is 69 s.The difference between the final time of aircraft 1 and 2 is 120 s.The difference between the final time of aircraft 1 and 3 is 189 s.These differences in time among the three aircraft at LALPI waypoint are below the required time separation of 200 s established above.Therefore, a new experiment has been conducted, which includes time-based separation logical constraint of at least 200 s between aircraft 9 International Journal of Aerospace Engineering as described in Section 4.1.1.The final mass and time given by the solution are reported for each aircraft in the second row of Table 2.
In Figure 4, the 3D view of the paths is represented in thin lines, and in Figures 5 and 6, the horizontal and vertical profiles are represented in dashed lines.It can be seen that     3, the information about the positions of the aircraft at the conflict region is reported in Table 4, and the final mass and time given by the solution of the corresponding OCP are reported for each aircraft in the first row of Table 5.
In Figure 8, the 3D view of the paths obtained in the solution is represented in thick lines, whereas Figures 9 and 10, the horizontal and vertical profiles are represented in solid lines.The computation time to find the solution has been 65:6747 s, and the final time of three aircraft is 1838 s.
To avoid the conflict, a new experiment has been conducted, which includes distance-based separation logical constraints of at least 5000 m (2:7 NM) in the horizontal profile and at least 1000 m (0:5 NM) in the vertical profile as described in Section 4.1.2.The final mass and time obtained in the solution are reported for each aircraft in the second row of Table 5.In Figure 8, the 3D view of the paths is represented in thin lines, whereas in Figures 9 and 10, the horizontal and vertical profiles are represented in dashed lines.
It can be seen from Figures 9 and 11 that the conflict is avoided by activating the logical constraints associated with the horizontal profile.In particular, from Figure 11, it can     12, in which the solution without distance-based separation constraints is represented in solid lines, whereas the solution with distance-based separation constraints is represented in dashed lines.It can be seen that introducing distance-based separation constraints does not significantly change the mass consumption of the aircraft.

Experiment 3: Minimum-Time STAR-Based Continuous
Descent along Converging Routes.In this experiment, a STAR-based CDA has been considered, that is, the lateral path followed by the aircraft has been assumed to be specified in a navigation chart.In particular, as in Experiment 1, the boundary conditions of the state variables have been selected from the chart of the Adolfo Suárez Madrid-Barajas (LEMD/MAD) TMA shown in Figure 3.The initial position of Aircraft 1, Aircraft 2, and Aircraft 3 is supposed to be coincident with the ROLDO, SOTUK, and MORAL waypoints, respectively.Aircraft are constrained to pass through TODNO and RESBI waypoints and their common final position is assumed to be the LALPI waypoint.For the setting of the cuboids centered at the given waypoints, the modeling defined in Equation (34) has been considered.In particular, the cuboid centered at TODNO waypoint has been defined by the two corners ð39:560 ∘ ,−4:24 ∘ , 5400 mÞ and ð39:640 ∘ ,−4:160 ∘ , 6000 mÞ, whereas the cuboid centered at RESBI waypoint has been defined by the two corners ð40:400 ∘ ,−4:150 ∘ , 4200 mÞ and ð40:480 ∘ ,−4:070 ∘ , 4800 mÞ.
The initial mass of the three aircraft has been assumed equal to the maximum landing weight of the aircraft.The specific boundary conditions of the state variables are the same as in Experiment 1 and are given in Table 1, and the final mass and time given by the solution of the corresponding OCP are reported for each aircraft in the second row of Table 6.In spite of the fact that the three aircraft have to pass through two waypoints, neither distance-based conflict nor time-based conflict arises during the flight except for the last part of the time discretization.In particular, Aircraft Therefore, a new experiment has been conducted, which, besides waypoint constraints, also includes time-based separation logical constraint of at least 200 s between aircraft as described in Section 4.1.1.The final mass and time given by the solution are reported for each aircraft in the third row of Table 6.In Figure 13, the 3D view of the paths is represented in thin lines, whereas in Figures 14 and 15, the horizontal and vertical profiles are represented in dashed lines.
It can be seen that, as expected, there are significant changes in both the horizontal and vertical profiles.The final time of the three aircraft is 2109 s, 1619 s, and 1909 s, respectively.The difference between the final time of aircraft 1 and 3 is 200 s, whereas the difference between the final time of aircraft 3 and 2 is 290 s.Notice that, unlike Experiment 1, aircraft 3 and 2 do not maintain the minimum required time separation of 200 s.This is due to the fact that the performance index includes in this case, besides the minimization of the duration of the flights, the penalty term associated with the waypoints constraints.Therefore, the saturation of the time-based constraints may not happen in this setting.In this case, the computation time to find the solution has been 49:0766 s.
The mass consumption is depicted in Figure 16, in which the solution without logical constraints is represented in solid lines, whereas the solution with waypoints constraints is represented in dashed lines.It can be seen that introducing waypoint constraints increases the mass consumption of Aircraft 1 and Aircraft 3, whereas that of Aircraft 2 remains almost unchanged.

Conclusions
In this paper, the trajectory planning problem for multiple aircraft has been studied in which logical constraints in disjunctive form are included in the model.The logical constraints in disjunctive form have been transformed into inequality and equality constraints which involves only  14 International Journal of Aerospace Engineering continuous auxiliary variables.In this way, the optimal control problem with logical constraints has been converted into a smooth optimal control problem which has been solved using standard techniques.This approach has been applied to the computation of the optimized profile descent of multiple aircraft in converging and intersecting arrival routes within the Adolfo Suárez Madrid-Barajas (LEMD/MAD) TMA.The results show the effectiveness of the proposed technique.

4
International Journal of Aerospace Engineering approximated as x e = λ e • ðR + h e Þ • cos θ e and y e = θ e • ðR + h e Þ.Finally, the state variable m is the aircraft mass.The control inputs are the bank angle μ, the engine thrust T, and the lift coefficient C L .Thus, uðtÞ = ðTðtÞ, μðtÞ, C L ðtÞÞ.The parameter R is the radius of Earth and η is the speeddependent fuel efficiency coefficient.Lift, L = C L Sq, and drag, D = C D Sq, are the components of the aerodynamic force.The parameter S is the reference wing surface area and q = 1/2ρV 2 is the dynamic pressure.A parabolic drag polar C D = C D0 + KC 2 L and an International Standard Atmosphere (ISA) model are assumed.The lift coefficient C L is a known function of the angle of attack α and the Mach number.

Figure 2 :
Figure 2: Two opposite corners of a cuboid in the 3D space.
Minimum-time continuous descent of three aircraft along routes converging at the same waypoint with time-based separation constraint among them (ii) Experiment 2. Minimum-time continuous descent of three aircraft along intersecting routes with distance-based separation constraint (iii) Experiment 3. Minimum-time STAR-based continuous descent of three aircraft along converging routes.The considered STAR procedure includes sequencing the aircraft at a merging point and passing through two other waypoints

5. 2 .
Experiment 2: Minimum-Time Continuous Descent of Three Aircraft along Intersecting Routes with Distance-Based Separation Constraint.This experiment has been designed to test the effectiveness of the logical constraint formalism described in Section 4.1.2,where three aircraft start at and have to reach three different positions flying along intersecting routes in which a mid-air conflict arises, in the sense that the experiment is intentionally designed in such a way that they reach almost the same point at the same time.Neither time-based nor distance-based separation constraints have been included in the model.The initial mass of the three aircraft has been assumed equal to the maximum

Figure 8 :Figure 9 :
Figure 8: Experiment 2. 3D view of the paths with (thin lines) and without (thick lines) distance-based separation constraints.

Figure 13 :Figure 15 :
Figure 13: Experiment 3. 3D view of the paths with (thin lines) and without (thick lines) waypoints and time-based separation constraints.

Table 1 :
Boundary conditions for Experiment 1 and Experiment 3.

Table 4 :
Positions of the aircraft at conflict region (at time 980 s) without distance constraint for Experiment 2.

Table 6 :
Figure 12: Experiment 2. Mass consumption with (DC, dashed lines) and without (FF, solid lines) distance-based separation constraints.Results of Experiment 3. Aircraft 3 reach the LALPI waypoint with a time separation of 6 s only.The computation time to find the solution has been 41:7981 s.