This is the author prepublication version. This paper does not include the changes arising from the revision, formatting and publishing processes. The final version that should be used for referencing is: A. Gardi, M. Marino, S. Ramasamy, T. Kistan, R. Sabatini
, “
4Dimensional Trajectory Optimisation Algorithm for Air Traffic Management Systems
”,
IEEE/AIAA 35
th
Digital Avionics Systems Conference (DASC)
, Sacramento, CA, USA, 2016.
4Dimensional Trajectory Optimisation Algorithm for Air Traffic Management Systems
Alessandro Gardi, Matthew Marino, Subramaniam Ramasamy and Roberto Sabatini
RMIT University
–
School of Engineering Melbourne, VIC 3000, Australia
Trevor Kistan
THALES Australia
–
Air Traffic Management Melbourne, VIC 3000, Australia
Abstract
—
This paper presents Multi Objective Trajectory Optimization (MOTO) algorithms that were developed for integration in stateoftheart Air Traffic Management (ATM) and Air Traffic Flow Management (ATFM) systems. The MOTO algorithms are conceived for the automationassisted replanning of 4Dimensional Trajectories (4DT) when unforeseen perturbations arise at strategic and tactical online operational timeframes. The MOTO algorithms take into account updated weather and neighbouring traffic data, as well as the related forecasts from selected sources. Multiple userdefined operational, economic and environmental objectives can be integrated as necessary. Two different MOTO algorithms are developed for future implementation in ATM systems: an enroute variant and a Terminal Manoeuvring Area (TMA) variant. In particular, the automated optimal 4DT replanning algorithm for enroute airspace operations is restricted to constant flight level to avoid violating the current vertical airspace structure. As such, the complexity of the generated trajectories reduces to 2 dimensions plus time (2D+T), which are optimally represented in the present 2D ATM display formats. Departing traffic operations will also significantly benefit from MOTO4D by enabling steep/continuous climb operations with optimal throttle, reducing perceived noise and gaseous emissions.
Keywords
—
4Dimensional Trajectory; Air Traffic Management; Decision Support System; Trajectory Optimization
I.
I
NTRODUCTION
More effective and comprehensive implementations of flight trajectory optimisation techniques are being considered as a very promising pathway to enhance the efficiency and flexibility of online air traffic operations both for short and longhaul flights. The fundamental proposition is not new, but stateoftheart flight planning methods are still based on the optimised vertical planning techniques initially developed in the 1970s and on lateral path planning based on optimal wind routing [16]. The known limitations are associated to the small set of optimality criteria (currently only fuel and timecosts) and to the fact that the initial flight plan is the static entity assumed as a reference for every subsequent amendment. As a result, the limited initial optimality is progressively compromised when strategic or tactical Air Traffic Management (ATM) and Air Traffic Flow Management (ATFM) amendments are introduced. Suitably defined models can replicate the various operational and environmental aspects that depend on the flown aircraft trajectory, allowing for more comprehensive and realtime MultiObjective Trajectory Optimisation (MOTO) techniques to be integrated in airborne avionics and groundbased Communication, Navigation, Surveillance and ATM (CNS/ATM) Decision Support Systems (DSS). These novel DSS integrating realtime MOTO capabilities have the potential to support a further enhanced exploitation of airspace and airport capacities, which is emerging thanks to the advances in CNS technologies. This paper describes two dedicated MOTO functionalities for integration in stateoftheart ATM and ATFM DSS. The MOTO algorithms are conceived for the automationassisted replanning of 4Dimensional Trajectories (4DT) when unforeseen perturbations arise at strategic and tactical online operational timeframes. The MOTO algorithms take into account updated weather and neighbouring traffic data, as well as the related forecasts from selected sources. Multiple operational, economic and environmental objectives can be userdefined as necessary. The project specifically addresses current and shortterm future ATM operational paradigms and regulations and thus two different MOTO algorithms are developed: an enroute variant and a Terminal Maneuvering Area (TMA) variant. In particular, the automated optimal 4DT replanning algorithm for enroute airspace operations is restricted to constant flight level to avoid violating the current vertical airspace structure. As such, the complexity of the generated trajectories reduces to 2 dimensions plus time (2D+T), which are optimally represented in the present 2D ATM display formats. The Air Traffic Controller (ATCo) may amend the flight level of an optimized 2D+T trajectory in the traditional manner if necessary. The constant flight level limitation will theoretically produce a suboptimal flight trajectory, however, the computed trajectory will be more efficient than vectors as the atmospheric wind field can be exploited to maximize aircraft ground speed while reducing fuel burn and emissions. The operational, economic and environmental benefits from MOTO are maximized in the TMA variant as the full 4DT MOTO (MOTO4D) routines are exploited to increase the operational efficiency and reduce fuel burn, emissions and noise. In the current paradigm aircraft are handled on a "first come first serve" basis where multiple objectives can be user defined and applied to the optimization problem. However, an evolution to "just in time" 4DT based operations in the TMA is supported by the MOTO functionalities. As a consequence, TMA operations are allowed to evolve and become more efficient by minimising scheduling
The work presented in this article was supported by THALES Australia under RMIT University Contract ID 0200312837
delays, providing automated pathstretching, reducing pollutants/noise emissions and supporting continuous descent approaches. As a consequence, ATCo workload will be relieved thanks to automationassisted deconfliction and path planning, and this is expected to increase operational safety. Departing traffic operations will also significantly benefit from MOTO4D by enabling steep/continuous climb operations with optimal throttle, reducing perceived noise and gaseous emissions. II.
O
PTIMAL
C
ONTROL
F
ORMULATION
Trajectory optimisation studies methods to determine the best possible trajectory of a dynamical system in a finitedimensional manifold, in terms of specific objectives and adhering to given constraints and boundary conditions [7]. This definition corresponds to the definition of Optimal Control Problems (OCP), and consequently the most traditional and theoretically rigorous way to pose a Trajectory Optimisation Problem (TOP) is based on the optimal control theory. Most OCP solution methods are conventionally categorised as either
direct methods
if based on the transcription to a finite NonLinear Programming (NLP) problem, or
indirect methods
if theoretical derivations based on the calculus of variation are implemented to formulate a BoundaryValueProblem (BVP) [811]. An additional class of OCP solution strategies is represented by heuristic methods. The optimal control formulation of TOP is based on a scalar time
and on vectors of timedependent state variables
, timedependent control variables
and system parameters
. Based on these, the following components are defined: dynamic constraints, path constraints, boundary conditions and cost functions. These components characterise a wellposed OCP and guide the selection of an appropriate numerical solution method and multiobjective decision logic. Dynamic constraints describe the feasible motion of the system (i.e. the aircraft, in our case) within the TOP as in:
̇
(1)
All nondifferential constraints insisting on the system between the initial and final conditions are classified as
path constraints
, as they restrict the path of states and controls of the dynamical system. A generalised expression accounting for both equality and inequality constraints is [12]:
(2)
Boundary conditions define the values that state and control variables of the dynamical system shall have at the initial and final times. A generalised expression also accounting for relaxed boundary conditions is:
(
)
(
)
(3)
Performance indexes
quantify the achievement of a particular objective by means of suitably defined
cost functions
. The generic formulation of a performance index
that takes into account both integral and terminal costs was introduced by Bolza [1315], and is expressed as:
(
)(
)∫
(4)
The optimisation is classified as
singleobjective
when an individual performance index
J
is introduced and
multiobjective
when two or multiple performance indexes
J
i
are defined. Different objectives can be conflicting, that is the attainment of a better
would lead to a worse
}
. Hence, the optimisation in terms of two or more objectives generates a number of possible compromise choices, for which a tradeoff decision logic must be adopted to identify an individual solution. The first step of the theoretical derivations required by indirect methods to formulate a BVP problem consists in the Lagrangian relaxation as in:
((
)
)
(
)(
)
∫ ̇
}
(5)
Where
are the Lagrangian multipliers and
is the Hamiltonian function, defined as:
(6)
In addition to being mathematically complex, the application of all the analytical steps involved in indirect methods lead to very different BVP depending on the initial problem statement and this aspect limits the flexibility of indirect methods. Some notable examples are presented in [15]. Finally, since considerable nonlinearities are present in most aerospace TOP, the reduction of nontrivial cases into linear quadratic BVP is typically precluded, and therefore solutions have to be attempted by either iterative methods, which conventionally exhibit nonglobal convergence, or by heuristic solution techniques. Direct methods, on the other hand, prescribe the immediate parameterisation of states and controls in the case of direct collocation methods, or of controls only in the case of direct shooting. This involves adopting a basis of known linearly independent functions
with unknown coefficients
in the general form:
∑
(7)
In global collocation (pseudospectral) methods, the evaluation of state and control vectors is performed at discrete collocation points across the problem domain using suitably defined orthogonal (spectral) interpolating functions [11, 16]. Similarly to the nonlinear BVP arising in indirect solution approaches, the NLP problems encountered in the direct solution of TOP are typically solved by iterative algorithms, or through some kind of heuristics. The initial steps involve the adoption of an ndimensional Taylor series expansion of
F(
x
)
to the third term as in:
(
) (
)
(
)
(
)
(8)
where:
(9)
(10)
()
(11)
An iterative NLP solution can thus be formulated so that the search direction at step
k
based on the ndimensional Newton method is written:
(12)
Various factors have to be considered when developing computationally efficient NLP solution strategies and some further detail is given in [17, 18].
A.
Optimality criteria and dynamic constraints
Optimality criteria and constraints are introduced in the optimisation problem by means of suitable models. Three degrees of freedom (3DOF) pointmass flight dynamics models are currently preferred for TOP involving mediumlarge transport aircraft and for possible avionics and ATM system implementations, as the overall size of the resulting NLP is considerably lower than in the case of six degrees of freedom (6DOF). A fairly comprehensive 3DOF formulation used in our MOTO algorithms assumes variable aircraft mass, constant vertical gravity and the effects of winds and is therefore written as:
̇ ̇ ̇ ̇
̇
̇
̇
(13)
where the state vector consists of: longitudinal velocity
v
[m s
1
], flight path angle
γ
[rad]; track angle
χ
[rad]; geographic latitude
[rad]; geographic longitude
λ
[rad]; altitude
z
[m]; aircraft mass
m
[kg]; and the control vector includes: thrust force
T
[N]; load factor
N
[ ]; bank angle
μ
[rad]. Other variables and parameters include: aircraft weight
and aerodynamic drag
D
[N]; wind velocity
v
w
in its three scalar components [m s
1
]; gravitational acceleration
g
[m s
2
]; Earth radius
R
E
[m]; fuel flow
FF
[kg s
1
] and thrust angle of attack
[rad]. The aerodynamic drag is modelled as:
(14)
where
is the local air density retrieved from weather input data grid or a weather model,
S
is the reference wing surface,
and
are the parabolic drag coefficients typically available from aircraft performance databases such as
Eurocontrol’s Base of Aircraft Data
(BADA). The lift coefficient
can be calculated from:
(15)
The thrust force control variable is most frequently expressed as the product of the throttle coefficient
(defined as dimensionless and ranging between 0 and 1), and the maximum thrust
as in:
(16)
This allows a natural nondimensionalisation of the control variable. For turbofan aircraft, the following empirical expressions were adopted in the development of BADA, to determine the climb thrust and the fuel flow
, which operationally equates to the maximum thrust
in all flight phases excluding takeoff [19]:
(17)
(18)
where
is the throttle control,
is the geopotential pressure altitude in feet,
is the deviation from the standard atmosphere temperature in kelvin,
is the true airspeed.
are the empirical thrust and fuel flow coefficients, which are also supplied as part of BADA for a considerable number of currently operating aircraft [19]. In order to calculate pollutant emissions as a function of the fuel flow, the emission index
is introduced as per the following definition:
∫
(19)
where the generic Gaseous Pollutant (GP) should be replaced by the specific one being investigated. While carbon dioxide (CO
2
) emissions are characterised by an approximately
constant emission index of
, an empirical model for carbon monoxide (CO) and unburned hydrocarbons (HC) emission indexes (
) in [g/Kg] at mean sea level based on nonlinear fit of experimental data from the ICAO emissions databank is:
(20)
where the fitting parameters
accounting for the emissions of 165 currently operated civil turbofan engines are
for CO and
for HC [20]. The nitrogen oxides (NO
X
) emission index [g/Kg] based on the curve fitting of 177 currently operated civil aircraft engines is [20]:
(21)
III.
M
ULTI
O
BJECTIVE
O
PTIMALITY
In the aviation domain, single and biobjective optimisation techniques have been exploited for decades but they accounted only for flight timerelated costs and fuelrelated costs. These techniques have also been implemented in a number of current generation FMS in terms of the Cost Index (
CI
), that allowed an optimal selection of Calibrated Air Speed (CAS) / Mach number based on time and fuel costs only. Conflicting objectives arise when introducing multiple environmental, economic and operational criteria in our MOTO algorithms. Furthermore, the implementation of constraints that are either unfeasible or contrasting the attainment of better optimality also has to be addressed by adopting suitable multiobjective optimality techniques. These can either be expressed
a priori
(i.e. beforehand) or
a posteriori
(i.e. afterwards). A
priori methods analyse ways to articulate the preferences to identify a combined objective which is then supplied to a singleobjective TOP solution algorithm, and include weighted global criterion (including the simple weighted sum), weighted minmax, weighted product, exponential weighted criterion, lexicographic and physical programming methods [7, 21]. A posteriori methods allow the user (or suitably defined decision logics) to select an individual optimal solution from a Pareto front of tradeoff choices, and include normal boundary intersection, normal constraint and physical programming methods [7, 21]. In general, a specific performance objective can be defined for each route segment. This performance objective is a multiobjective generalisation of the CI. In general, the weightings can be varied dynamically among the different flight phases of the flight. Since computational times are a crucial aspect in online 4DT planning applications, an
a priori
articulation of preference involving the weighted sum of the various performance indexes
J
i
is employed to combine the multiple conflicting operational, economic and environmental objectives. For a more detailed discussion the reader is referred to [7, 9, 11]. IV.
MOTO
FOR
T
ERMINAL
M
ANOEUVRING
A
REA
The MOTO for Terminal Sequencing and Spacing (TSS) is based on the optimal control formulation presented in section 2. The optimal control formulation shows its best potential in this context as this is a truly 4D application where multiple variables and constraints are at play at the same time. The adoption of optimal controlbased techniques in the TSS context allows accurate continuous descent profiles to be natively generated, CO/HC emissions to be minimised, as well as optimised path stretching to be introduced as necessary to achieve the set constraints. The challenges associated with optimal controlbased MOTO implementations lie in the characteristics of the generated trajectory. More specifically, the mathematically optimal 4DT generated as output by these algorithms is a discretised Continuous/Piecewise Smooth (CPWS) curve, which in general may not be flyable by human pilots nor by conventional Automatic Flight Control Systems (AFCS), as it includes transition manoeuvres involving multiple simultaneous variations in the control inputs. Although continuous, the variations in control input profiles can also be particularly abrupt, exceeding what can be practically and safely feasible in the modelled aircraft. Moreover, the discretised CPWS consists of a very high number of 4D waypoints, which would have unacceptable impacts on the bandwidth of the datalink where they would be exchanged. Therefore, a postprocessing stage is introduced to segment the discretised CPWS trajectory in feasible flight legs, including straight and level flight, straight climbs and descents, level turns, and climbing/descending turns. The final result is a concisely described 4DT consisting of feasible flight segments. The MOTO algorithm implemented for TMA operations accounting for the combined objective is represented in Fig. 1.
Fig. 1.
Block diagram of the MOTO4D algorithm for TMA operations.
The currently employed solution technique is based on direct solution methods of the global orthogonal collocation family (pseudospectral), which are arguably the most effective solution methods currently available [7]. Additionally, to further enhance the algorithm stability and convergence
Mathematically Optimal 4DT (CPWS)O
PERATIONAL
S
MOOTHING
MOTO4D GFS Weather Field{
v
w
, p, RH, T
}(
, , z, t
)CombinedObjectiveTerminalConditions{
f
,
f
, z
f
,t
f
,
f
,
f
,
f
}InitialConditions{
0
,
0
, z
0
,t
0
, m
0
, v
0
,
0
,
0
,
0
}4DT GenerationGlobal orthogonalcollocationManoeuvre IdentificationOperationalCriteriaTransitionLikelihoodFinal Solution: Optimal and Feasible 4DTFeasible 4DT Optimiser LevelTurnModelStraightClimbModelStraight& LevelModelStraightDescentModelTurningClimbModelTurningDescentModel
performances, path constraints and boundary conditions are automatically strengthened on all state and control variables to restrict the search domain as much as feasible.The TMA MOTO algorithm considers the controlled time of arrival target defined by Arrival Management (AMAN) functionalities. This is used as the final time constraint by the TMA MOTO algorithm. The Estimated Time of Arrival (ETA) may be computed at multiple fixes along the flight path.
A.
Operational Smoothing
The 4DT intent data produced as output of the operational 4DT smoothing include 4D waypoints (latitude, longitude, altitude and time), flyby/flyover turn parameters, as well as performance criteria and restrictions. Based on the manoeuvre identification algorithm assessment on the CPWS trajectory, the lateral path is constructed in terms of segments (straight and turns), whose geometric characteristics depend on the required course change
and the predicted Ground Speed (GS) of the aircraft during turns. The ATM DSS computes turn radius
R
and True Air Speed
based on the selected altitude and takes into account the predicted wind speed
at that altitude. The bank angle
is selected based on the aircraft dynamics database, on the current altitude and on local/airline operational restrictions. The radius
R
of a turn can be conservatively calculated based on the maximum GS during the turn as:
(22)
where
g
is the gravity acceleration module. The turn arc length
is simply given by:
(23)
where the track change
is the difference in radians between the final and the initial ground tracks. The following prescriptions on
,
,
lead distance
(from the turn initiation to the 4DT waypoint) and
abeam distance
(between the 4DT waypoint and the point of the circular arc abeam it) are implemented to plan fixedradius flyby turns that comply with RTCA DO229D, DO236C and DO283B [22, 23]: 1.
The flyby transition is defined by equation (22) combined with the following equations [22]:
(24)
⃗
⃗

(25)
(26)
(27)
2.
The geometry of the Fixed Radius Transition (FRT) is defined by the track change
and the radius
. The Lead Distance
and the Abeam Distance
are defined based on the radius and the track change as per the following equations:
(28)
(29)
3.
When transitioning from one airway to another, if both require a FRT at the common waypoint, the smaller of the two radii applicable shall be selected. In the case one of the two airways does not involve a FRT, the FRT of the other shall be implemented. The criteria for implementing the definitions above in different flight phases and at various altitudes are given in Table 27 of RTCA DO229D [23]. V.
MOTO
FOR
E
NROUTE
O
PERATIONS
Even neglecting current operational restrictions, the enroute context typically involves limited or no variations in cruise altitude. Additionally, the route distance and time involved are substantially greater than the transients associated with the aircraft flight dynamics model, which can therefore be neglected since the entire phase takes place in quasisteadystate conditions. As a result, optimal control based techniques have limited or no advantages in this context. As a result, the aviation research community is investigating alternative formulations typically based on geometric steadystate trajectory approximations. The mainstream philosophies pursued involve either the generation of a grid of discrete 2D points across the cruise space and the implementation of suitable tree search algorithms to identify the optimal path, or the implementation of evolutionary algorithms or more general NLP solvers to optimise a discrete geometric trajectory. Predefined flight level (FL) conditions are typically considered to reduce the number of optimisation parameters and also to meet current operational paradigms that prescribe semicircular cruise levels. The MOTO algorithm implemented for enroute operations adopts a polygonal trajectory model, simplified aircraft performance equations based on the steadystate approximation of eq. 1321 and NLP solvers. The variables considered include the latitude and longitude of each and every trajectory waypoint, the airspeed flown (which dictates the throttle setting and therefore the pollutant emissions according to the steadystate approximation of equations 1321. VI.
TMA
MOTO
C
ASE
S
TUDY
The sequencing of dense arrival traffic towards a single final approach path was extensively assumed as a representative case study of online tactical TMA operations. The results of one exemplary simulation run are depicted in Fig. 2. The TMA MOTO considers the selected arrival sequence. Longitudinal separation is enforced at the merge point to ensure safe separation upon landing based on waketurbulence category and to prevent separation infringements in the approach phase itself. The assumed minimum longitudinal separation is 4 nautical miles on the approach path for medium category aircraft approaching at 140 knots, therefore the generated time slots are characterized by a 90~160 seconds separation depending on the waketurbulence categories of two consecutive traffics. Fig. 3 depicts the computed 4DT in the AMAN schedule display format. Waypoints and lines depicted